【推荐系统入门到项目实战】（五）：SVD矩阵分解
- 🌸个人主页：JOJO数据科学
- 📝个人介绍：统计学top3高校统计学硕士在读
- 💌如果文章对你有帮助，欢迎✌
关注、👍点赞、✌收藏、👍订阅专栏- ✨本文收录于【推荐系统入门到项目实战】本系列主要分享一些学习推荐系统领域的方法和代码实现。
引言
之前我们介绍了矩阵分解ALS算法，并介绍了几个案例，下面我们来看看另一种使用广泛的矩阵分解方法——SVD，及其在推荐系统上的应用。
老规矩，我们首先来回顾一下推荐算法的常见方法框架。主要分为两大类
- 1.基于内容的推荐
- 2.基于协同过滤的推荐
而基于协同过滤的推荐是推荐系统的主流思想之一。其中矩阵分解是一个重要的模块。本文我们来讨论一下SVD矩阵分解的原理与实现。
矩阵分解前置知识（MF）
首先，矩阵分解是将一个原始的大的用户矩阵进行分解成多个矩阵的乘积。在推荐系统中，矩阵分解属于隐语义模型一个大的分支。之前我们讨论的ALS算法和推荐系统，我们是以一个电影评分数据集为例，将其分解为user矩阵和item矩阵。然后我们使用ALS算法将其估计。如下图所示。
之前我们是假设评分矩阵可以分解成User矩阵和Item矩阵。这一章我们来讨论一些背后的数学知识和直觉。
现在我们来看看SVD（奇异值分解）是怎么做的。首先，我们来回顾一下一些线性代数的基本知识：
矩阵的特征分解：
特征分解，是将矩阵分解为特征值和特征向量表示的矩阵之积的方法，也称为谱分解
N 维非零向量 v 是 N×N 的矩阵 A 的特征向量，当且仅当下式成立
$$ Av=\lambda{v} $$
$\lambda$为特征值（标量），v为特征值对应的特征向量。特征向量被施以线性变换 A 只会使向量伸长或缩短，而方向保持不变。
接下来我们求解这个特征向量呢?在线性代数中，我们知道，求解特征方程：$|A-\lambda I|=0$即可。
$p(\lambda):=|A-\lambda I|=0$称作矩阵的特征多项式
这个特征多项式是关于$\lambda$的N次多项式，用N个解
$$ p(\lambda)=(\lambda-\lambda_1)^{n_1}...(\lambda-\lambda_k)^n_k=0 $$
其中$\sum_{i=1}^{k}n_i=N$,对于每一个特征值，我们都有$(A-\lambda_iI)v=0$
可能大家对于这些公式看起来觉得有些枯燥，下面我们来通过一个具体的例子来理解。
我们要求矩阵A的特征值和对应的特征向量，即求解下述方程
得到：
求解的$\lambda_1=1,\lambda_2=\lambda_3=0$。
当$\lambda=1$
简化得到：
然后有：
即：
令$x_1=1$,则有特征向量为：
同理，当$\lambda_2=\lambda_3=0$，计算可得特征矩阵
现在我们知道了如何计算特征值和特征向量，为什么我们要讲这些呢？是因为当$A$是一个$N\times N$的方针时，那么有A的特征分解如下：
- U 的列向量是 A 的特征向量
- Λ 是对角矩阵，元素是特征向量的特征值
如下所示，对于这个方阵，我们可以很轻松的得到矩阵分解后的结果。
具体的计算步骤大家不用担心，可以直接在python中轻松实现，代码如下：
import numpy as np
A = np.array([[5,3],[1,1]])
lam,U=np.linalg.eig(A)#特征分解
inv = np.linalg.inv(U)#求逆矩阵
print(A)
print('特征值:',lam)
print('特征向量',U)
print('特征向量的逆',inv)[[5 3]
[1 1]]
特征值: [5.64575131 0.35424869]
特征向量 [[ 0.97760877 -0.54247681]
[ 0.21043072 0.84007078]]
特征向量的逆 [[ 0.89807344 0.5799321 ]
[-0.2249599 1.04510774]]现在，我们已经解决了对于对称的矩阵分解问题。但是，在实际中很多矩阵都是 非对称的矩阵,A不是方阵，维度为m*n。这种情况就需要使用我们的SVD（奇异值分解）方法。
SVD矩阵分解
基本理论
对于非对称矩阵，我们可以通过下面的方式将其变换为对称矩阵：$AA^T$和$A^TA$
将A和A的转置矩阵进行相乘，得到对称方阵,然后我们有
此时$\Lambda_1,\Lambda_2$均为对角矩阵，具有相同的非零特征值。
假设这些特征值为$\sigma_1,\sigma_2,...\sigma_k$。k不超过m和n，也就是k<=min(m,n)
此时矩阵A的特征值
我们可以得到为奇异值分解
P为左奇异矩阵，m*m维
Q为右奇异矩阵，n*n维。如下图所示
Λ对角线上的非零元素为特征值λ1, λ2, ... , λk
在推荐系统中
左奇异矩阵：User矩阵
右奇异矩阵：Item矩阵
接下来我们来看一个具体的例子
奇异值分解：
$λ_1$为特征值，$p_1$为左奇异矩阵的特征向量,$q_1$为右奇异矩阵的特征向量
具体的数学推导大家不用太过于纠结，我们一样在python中实现。
from scipy.linalg import svd
import numpy as np
from scipy.linalg import svd
A = np.array([[1,2],
[1,1],
[0,0]])
p,s,q = svd(A,full_matrices=False)
print('P=', p)
print('S=', s)
print('Q=', q)
P= [[-0.85065081 -0.52573111]
[-0.52573111 0.85065081]
[ 0. 0. ]]
S= [2.61803399 0.38196601]
Q= [[-0.52573111 -0.85065081]
[ 0.85065081 -0.52573111]]Excellent！现在给我们一个任意的矩阵，我们都可以使用SVD将其进行分解，那么我们如何理解SVD的作用呢？
对于奇异值,它跟我们特征分解中的特征值类似，在奇异值矩阵中也是按照从大到小排列，并且往往前10%甚至1%的奇异值的和就占了全部的奇异值之和的90%以上。也就是说，我们也可以用最大的k个的奇异值和对应的奇异向量来近似描述矩阵。如下图所示，这样大大的减少了计算成本。直观上也就是用较少的特征反应了绝大部分信息，因此在本质上SVD是一个降纬的作用，类似我们的主成分分析。
SVD的应用
现在我们来看看对于特征值矩阵，我们如果只包括某部分特征值，结果会怎样？
矩阵A：大小为1440*1080的图片
- Step1，将图片转换为矩阵
- Step2，对矩阵进行奇异值分解，得到p,s,q
- Step3，包括特征值矩阵中的K个最大特征值，其余特征值设置为0
- Step4，通过p,s',q得到新的矩阵A'，对比A'与A的差别
import numpy as np
from scipy.linalg import svd
from PIL import Image
import matplotlib.pyplot as plt
# 取前k个特征，对图像进行还原
def get_image_feature(s, k):
# 对于S，只保留前K个特征值
s_temp = np.zeros(s.shape[0])
s_temp[0:k] = s[0:k]
s = s_temp * np.identity(s.shape[0])
# 用新的s_temp，以及p,q重构A
temp = np.dot(p,s)
temp = np.dot(temp,q)
plt.imshow(temp, cmap=plt.cm.gray, interpolation='nearest')
plt.show()
print(A-temp)
# 加载256色图片
image = Image.open('./256.bmp')
A = np.array(image)
# 显示原图像
plt.imshow(A, cmap=plt.cm.gray, interpolation='nearest')
plt.show()
# 对图像矩阵A进行奇异值分解，得到p,s,q
p,s,q = svd(A, full_matrices=False)
# 取前k个特征，对图像进行还原
get_image_feature(s, 5)
get_image_feature(s, 50)
get_image_feature(s, 500)
可以看出，当k=50的时候，效果就很不错了， 此时，我们只需要保存 (1440+1+1080)*50=126050个元素，占比126050/(1440*1080)=8%。
因此，我们只用少量的信息（比如10%），可以还原大部分图像信息（比如99%）
将user-item评分问题，转化为SVD矩阵分解
如果我们想要看user2对item3评分，可以得到
$4=-0.46*16.47*(-0.29)+(-0.30)*6.21*(-0.38)+(-0.65)*4.40*(-0.13)+0.28*2.90*0.87+0.02*1.58*(-0.03)$
A中各元素=user行向量，item列向量，奇异值的加权内积
实际上，我们发现user矩阵的最后一列是没有用到的，而且我们还可以使用更少的特征，比如特征个数=2
得到近似解A'
传统SVD在推荐系统中的应用
我们可以通过k来对矩阵降维
第i个用户对第j个物品的评分
完整的SVD，可以将M无损的分解成为三个矩阵
为了简化矩阵分解，我们可以使用k，远小于min(m,n)，对矩阵M近似还原。乍一看，感觉没什么问题，但是传统SVD在使用上是有局限：
- SVD分解要求矩阵是稠密的 => 矩阵中的元素不能有缺失。
- 因此实际上传统SVD更适合做降维
SVD的要求是稠密的，但是一般的推荐问题中，数据都是很稀疏的。 因此我们需要使用一个近似的方法FunkSVD
FunkSVD
FunkSVD的基本思想如下：
- 避开稀疏问题，而且只用两个矩阵进行相乘，如下所示
$$ M_{m\times n}\approx P_{m\times k}^TQ_{k\times n} $$
- 损失函数=P和Q矩阵乘积得到的评分，与实际用户评分之差
- 让损失函数最小化 => 最优化问题。
然后我们就有以下
$$ min\sum_{p,q}(m_{ij}-q_j^Tp_i)^2+\lambda(||p_i||^2+||q_j||^2) $$
到这里是不是有一种很熟悉的感觉，这里与我们之介绍的ALS矩阵分解是一样的，只不过当时我们使用ALS方法进行优化，在这里我们使用梯度下降法(SGD)
Step1，通过梯度下降法，求解P和Q使得损失函数最小化
Step2，通过P和Q将矩阵补全
Step3，针对某个用户i，查找之前值缺失的位置，按照补全值从大到小进行推荐
具体代码实现如下：
数据集：MovieLens下载地址：https://www.kaggle.com/jneupane12/movielens/download
主要使用的文件：ratings.csv
格式：userId, movieId, rating, timestamp
记录了用户在某个时间对某个movieId的打分情况
在这里我们使用的数据依然是MovieLens数据集，大家可以自行下载。下面看具体代码，要注意的是，我们要设置algo = SVD(biased=False)
!pip install surprise
from surprise import SVD
from surprise import Dataset
from surprise.model_selection import cross_validate
from surprise import Reader
from surprise import accuracy
from surprise.model_selection import KFold
import pandas as pd
import time
time1=time.time()
# 数据读取
reader = Reader(line_format='user item rating timestamp', sep=',', skip_lines=1)
data = Dataset.load_from_file('/content/drive/MyDrive/02名企/4.1 SVD矩阵分解与基于内容的推荐/L4-code/L4/MovieLens/ratings.csv', reader=reader)
train_set = data.build_full_trainset()
# 使用funkSVD
algo = SVD(biased=False)#这里默认是True，使用的是biasSVD
# 定义K折交叉验证迭代器，K=3
kf = KFold(n_splits=3)
for trainset, testset in kf.split(data):
# 训练并预测
algo.fit(trainset)
predictions = algo.test(testset)
# 计算RMSE
accuracy.rmse(predictions, verbose=True)
uid = str(196)
iid = str(302)
# 输出uid对iid的预测结果
pred = algo.predict(uid, iid, r_ui=4, verbose=True)RMSE: 0.8724
RMSE: 0.8712
RMSE: 0.8728
user: 196 item: 302 r_ui = 4.00 est = 4.33 {'was_impossible': False}biasSVD算法
BiasSVD算法原理：
用户有自己的偏好(Bias)，比如乐观的用户打分偏高
商品也有自己的偏好(Bias)，比如质量好的商品，打分偏高
将与个性化无关的部分，设置为偏好(Bias)部分
$$ min\sum_{i,j}(m_{ij}-\mu-b_i-b_j-q_j^Tp_i)^2+\lambda(||b_i||^2+||b_j||^2+||p_i||^2+||q_j||^2) $$
$u$:表示表示所有数据的平均值
$b_i$:表示用户的评分倾向
$b_j$:表示商品的评分倾向
$p_iq_j^T$:表示个性化推荐部分
这里和我们之前的baseline方法很类似，但是在这里，我们增加了一个个性化推荐部分。
下面来看代码实现：
biasSVD的实现和funkSVD基本一致，但是这个时候我们直接设置algo = SVD()
!pip install surprise
from surprise import SVD
from surprise import Dataset
from surprise.model_selection import cross_validate
from surprise import Reader
from surprise import accuracy
from surprise.model_selection import KFold
import pandas as pd
import time
time1=time.time()
# 数据读取
reader = Reader(line_format='user item rating timestamp', sep=',', skip_lines=1)
data = Dataset.load_from_file('/content/drive/MyDrive/02名企/4.1 SVD矩阵分解与基于内容的推荐/L4-code/L4/MovieLens/ratings.csv', reader=reader)
train_set = data.build_full_trainset()
# 使用biasSVD
algo = SVD()
# 定义K折交叉验证迭代器，K=3
kf = KFold(n_splits=3)
for trainset, testset in kf.split(data):
# 训练并预测
algo.fit(trainset)
predictions = algo.test(testset)
# 计算RMSE
accuracy.rmse(predictions, verbose=True)
uid = str(196)
iid = str(302)
# 输出uid对iid的预测结果
pred = algo.predict(uid, iid, r_ui=4, verbose=True)RMSE: 0.8437
RMSE: 0.8466
RMSE: 0.8432
user: 196 item: 302 r_ui = 4.00 est = 4.09 {'was_impossible': False}从这个结果来看，BiasSVD比FunkSVD效果是要好一点的，因为我们增加了要估计的参数。
SVD++
在bias的基础上加入了隐性评分
SVD++算法原理：
在BiasSVD算法基础上进行了改进，考虑用户的隐式反馈
隐式反馈：没有具体的评分，但可能有点击，浏览等行为
对于某一个用户i，假设他的隐式反馈item集合为I(i)
用户i对商品j对应的隐式反馈修正值为$c_{ij}$
用户i所有的隐式反馈修正值之和为$\sum_{s\in N(i)}c_{sj}$
优化目标函数：
在考虑用户隐式反馈的情况下，最终得到P和Q
代码实现如下：
需要注意的是SVD++模型需要导入SVDpp，其余用法一致。SVDpp跑的时间要久一些
#!pip install surprise
from surprise import SVDpp
from surprise import Dataset
from surprise.model_selection import cross_validate
from surprise import Reader
from surprise import accuracy
from surprise.model_selection import KFold
import pandas as pd
import time
time1=time.time()
# 数据读取
reader = Reader(line_format='user item rating timestamp', sep=',', skip_lines=1)
data = Dataset.load_from_file('/content/drive/MyDrive/02名企/4.1 SVD矩阵分解与基于内容的推荐/L4-code/L4/MovieLens/ratings.csv', reader=reader)
train_set = data.build_full_trainset()
# 使用svd++
algo = SVDpp()
# 定义K折交叉验证迭代器，K=3
kf = KFold(n_splits=3)
for trainset, testset in kf.split(data):
# 训练并预测
algo.fit(trainset)
predictions = algo.test(testset)
# 计算RMSE
accuracy.rmse(predictions, verbose=True)
uid = str(196)
iid = str(302)
# 输出uid对iid的预测结果
pred = algo.predict(uid, iid, r_ui=4, verbose=True)RMSE: 0.8427
RMSE: 0.8416
RMSE: 0.8412
user: 196 item: 302 r_ui = 4.00 est = 4.08 {'was_impossible': False}可以看出，SVD++相对于BiasSVD效果差不多。但是需要计算的时间要长。
总结
奇异值分解（SVD）可以对矩阵进行无损分解
在实际中，我们可以抽取前K个特征，对矩阵进行降维，使得只使用较少的特征就可以包含大部分信息 （10%的特征包含99%的信息）
在评分预测中我们使用funkSVD，只根据实际评分误差进行目标最优化，这里其实就是我们上一章所用的矩阵分解方法，只不过在这里我们使用的优化方法不再是ALS而是SGD（随机梯度下降）
在funkSVD的基础上，加入用户/商品偏好 => BiasSVD，这里可以类比上一章介绍的baseline方法。
在BiasSVD的基础上，考虑用户的隐式反馈 => SVD++。
现在，我们基本把矩阵分解的基本介绍完成了，但是矩阵分解（MF）存在不足，因为它只考虑了user和item特征，对于其他特征的利用我们需要使用新的工具，后续会介绍因子分解机(FM)和DeepFM模型。
