基于协同过滤的推荐引擎（理论部分）
记得原来和朋友猜测过网易云的推荐是怎么实现的，大概的猜测有两种：一种是看你听过的和收藏过的音乐，再看和你一样听过这些音乐的人他们喜欢听什么音乐，把他喜欢的你没听过的音乐推荐给你；另一种是看他听过的音乐或者收藏的音乐中大部分是什么类型，然后把那个类型的音乐推荐给他。当然这些都只是随便猜测。但是能发现一个问题，第二种想法很依赖于推荐的东西本身的属性，比如一个音乐要打几个类型的标签，属性的粒度会对推荐的准确性产生较大影响。今天看了协同过滤后发现其实整个算法大概和第一种的思想差不多，它最大的特点就是忽略了推荐的东西本身的属性，而是根据其他用户对它的喜好程度进行推荐的。
我在网上随便找了个电影评分数据集——Movielens电影评分数据，选择recommended for education and development里面的ml-latest-small.zip，打算边学边用这个数据集做测试并看看效果。写在一篇太长了，下一篇发实战。也可以用其他数据集，这个博客有几个推荐引擎测试数据汇总。
整个算法学习是通过图灵程序设计丛书的《Machine Learning in Action》学习的。
什么是协同过滤
协同过滤（collaborative filtering）是通过将用户和其他用户的数据进行对比来实现推荐的算法。
协同过滤流程图
协同过滤流程图.jpg
图用Grafio 3画的，看着好有成就感哈哈哈哈哈哈哈:-):-):-)，但是上传了以后这个分辨率好像低了T.T。
相似度
相似度计算就是看两个物品（或用户）有多相似，拿电影来说，可能会比较类型、导演、地区等等，但是在协同过滤里，不关心这些属性，严格地按照许多用户的观点来计算相似度。比如下面的电影和用户评分矩阵：
电影_用户矩阵.png
相似度计算
欧氏距离（euclidean metric）
欧氏距离指在m维空间中两个点之间的真实距离，或者向量的自然长度（即该点到原点的距离）。在二维和三维空间中的欧氏距离就是两点之间的实际距离，就是那个“根号下横坐标差的平方加纵坐标差的平方”。
假设我们看《海洋奇缘》和《一条狗的使命》的相似度，《海洋奇缘》列向量是（4， 3，1），《一条狗的使命》列向量是（4，3，2），做差再平方最后开根号。
一条狗与海洋奇缘的欧氏距离.png
Python实现代码：
import numpy as np
def eulid_sim(colA, colB):
return 1.0/(1.0 + np.linalg.norm(colA - colB))
代码解析：
numpy的线性代数（Linear algebra）库linalg里有一个norm函数，用于求范数（normal form），如果不指定范数的阶数，就默认指2范数，就是向量各个元素平方和开根号，比如向量A=(4，2，2)，向量A的2范数||A||=根号下(4^2+2^2+2^2)（简书不支持LaTeX真是不方便……），所以np.linalg.norm(colA - colB))就是向量A和B的欧式距离了。
1.0/(1.0 + 欧式距离)的作用是使相似度的值在0到1之间变化，越相似，相似度的值越大，距离为0时，相似度为1。
皮尔逊相关系数
皮尔逊相关系数.jpg
资料参考这里
- 皮尔逊相关系数
皮尔逊相关系数可以用来度量两个向量之间的相似度，比欧氏距离好的一点是它对用户评级不敏感，比如某个狂躁者对所有电影评分都是5，一个忧郁者对所有电影评分都是1，皮尔逊相关系数会认为这两个向量相等。看最后一个公式，对比两个向量的余弦公式，长得挺像，据说皮尔逊系数是两组向量的余弦。 - Z分数
z分数（z-score）,也叫标准分数（standard score）是一个数与平均数的差再除以标准差的过程。
z分数可以回答这样一个问题：一个给定分数（指评分）距离平均数多少个标准差?在平均数之上的分数会得到一个正的标准分数，在平均数之下的分数会得到一个负的标准分数。
z分数是一种可以看出某分数在分布中相对位置的方法。z分数能够真实的反应一个分数距离平均数的相对标准距离。如果我们把每一个分数都转换成z分数，那么每一个z分数会以标准差为单位表示一个具体分数到平均数的距离或离差。 - Python代码实现
def pearson_sim(colA, colB):
if len(colA) < 3:
return 1.0
return 0.5 + 0.5 * np.corrcoef(colA, colB, rowvar=0)[0][1]
- 代码解析
len(colA) < 3是检查是否有3个或更多的点，如果不存在，则返回1，两向量完全相关。
corrcoef(colA, colB, rowvar=0)返回的是变量的相关系数矩阵，第[0][1]个元素是相关系数，rowvar=0代表列是variables。API在这里。
0.5 + 0.5 *皮尔逊相关系数目的也是将取值范围归一化到0～1之间，皮尔逊相关系数的取值范围是-1～1，所以用0.5+0.5*系数的方式归一化。
余弦相似度
余弦相似度就是计算两个向量夹角的余弦值，如果夹角为90度，则相似度为0；如果方向相同，相似度为1。因为余弦值的范围也是-1～1，所以需要用同样的方法进行归一化。
前面说过，||A||代表向量的2范数。
Python代码
def cos_sim(colA, colB):
num = float(colA.T * colB)
denom = np.linalg.norm(colA) * np.linalg.norm(colB)
return 0.5 + 0.5 * (num / denom)
相似度选择
计算两个电影之间的距离，是基于物品（item-based）的相似度，计算用户的距离，是基于用户（user-based）的相似度。到底使用哪种相似度，取决于用户和物品的数量。基于物品的相似度会随着物品增加而增加，基于用户的相似度会随着用户的增加而增加。如果用户很多，则倾向于物品相似度计算方法。对于大部分推荐引擎而言，用户数目往往大于物品数目，所以一般用物品相似度。
构建评分估计函数（就是预估用户会给要推荐的电影打多少分）
对于用户没看过的电影，我们要预测他会为这些电影打多少分，然后把排名前N个的电影推荐给他。
一般评分估计
流程图：
一般评分预估算法流程图.PNG
# 计算某个物品和所有其他物品的相似度，进行累加，连评分也累加，最后用累加的总评分／总相似度得到预测该用户对新物品的评分
# data_mat:物品-用户矩阵
# user:用户编号
# item:要预测评分的物品编号
# sim_meas:相似度计算方法
def stand_est(data_mat, user, item, sim_meas):
n = np.shape(data_mat)[1] # 取第1轴的元素个数，在这里也就是列数
sim_total = 0.0
rat_sim_total = 0.0
# 遍历整行的每个元素
for j in range(n):
user_rating = data_mat[user, j] # 取一个评分
if user_rating == 0: # 如果用户没有评分，就跳过这个物品
continue
# 找出要预测评分的物品列和当前取的物品j列里评分都不为0的下标（也就是所有评过这两个物品的用户对这两个物品的评分）
overlap = np.nonzero(np.logical_and(data_mat[:, item].A > 0, data_mat[:, j].A >0))[0]
# 如果预测评分的物品列（假设叫列向量A）和当前取的物品j列（假设叫列向量B）没有都非零的项（也就是说两组向量里要么A评分B没评分，要么B评分A没评分），则认为相似度为0，否则，计算相似度
if len(overlap) == 0:
similarity = 0
else:
similarity = sim_meas(data_mat[overlap, item], data_mat[overlap, j]) # 注意overlap是一个array，所以这里还是传的两个列向量，两个元素中都没有0的列向量
print('the %d and %d similarity is %f' % (item, j, similarity))
sim_total += similarity # 累加相似度
rat_sim_total += similarity * user_rating # 累加相似度*评分，
if sim_total == 0:
return 0
else:
return rat_sim_total / sim_total # 总评分／总相似度，除以总相似度是为了归一化，将评分归到相似度的范围（比如0～5）
代码解析：
这里比较难理解的就是overlap一句，data_ma[:,item]代表取矩阵中编号为item的那一列，.A操作是将返回值变为ndarray，data_ma[:,item].A>0会产生一个shape相同的布尔型矩阵，根据是否大于零置True或False，logical_and方法对两个布尔矩阵求逻辑与，nonzero方法找出逻辑与后非零值的下标。
整个过程的作用就是从两个物品列中晒出两物品都被评分的行的下标，用于相似度计算。
为了搞清楚overlap一句的作用，我做了如下测试：
data_mat = np.mat([[1, 2, 0, 0, 0],[6, 7, 8, 1, 10]])
a = data_mat[:,0]
b = data_mat[:,3]
print(a) # [[1] [6]]
print(a.A) # [[1] [6]] a.A和a长得一样，有什么差别呢？打印类型看看
print(type(a)) # <class 'numpy.matrixlib.defmatrix.matrix'>
print(type(a.A)) # <type 'numpy.ndarray'> 查看API，是这样写的：return 'self' as an'ndarray' object
print(type(a.A1)) # <type 'numpy.ndarray'> 查看API，是这样写的：return 'self' as a flattened 'ndarray'
print(a.A > 0) # [[ True][ True]]
print(b.A > 0) # [[ False][ True]]
print(np.logical_and(a.A > 0, b.A > 0)) # [[False][ True]] 每个值做逻辑与运算
print(np.nonzero(np.logical_and(a.A > 0, b.A > 0))) # (array([1]), array([0])) np.nonzero，return the indices of the elements that are non-zero.第0个元素是第0轴的下标，第1个元素是第1轴的下标
print(np.nonzero(np.logical_and(a.A > 0, b.A > 0))[0]) # [1] 因为是二维单列向量，取第0轴下标就行
SVD评分估计
流程图：
svd评分预测流程图.png
# 用svd将矩阵变换到低维空间，再给出预估评分
# data_mat:物品-用户矩阵
# user:用户编号
# item:物品编号
# sim_meas:相似度计算方法
def svd_est(data_mat, user, item, sim_meas):
n = np.shape(data_mat)[1]
sim_total = 0.0
rat_sim_total = 0.0
u, sigma, vt = np.linalg.svd(data_mat) # 对原数据矩阵做svd操作
sig4 = np.mat(np.eye(4) * sigma[:4]) # sigma[:4]是取sigma矩阵的前四个，python为了节省空接，sigma矩阵存成了行向量，所以通过eye(4)将其变回对角矩阵
x_formed_items = data_mat.T * u[:,:4] * sig4.I # 利用u矩阵将其转换到低维空间,I操作是矩阵的逆.
for j in range(n):
user_rating = data_mat[user, j]
if user_rating == 0 or j == item:
continue
similarity = sim_meas(x_formed_items[item, :].T, x_formed_items[j,:].T) #行向量转成列向量
print('the %d and %d similarity is %f' % (item, j, similarity))
sim_total += similarity
rat_sim_total += similarity * user_rating
if sim_total == 0:
return 0
else:
return rat_sim_total / sim_total
代码解析：
sig4取4个奇异值不是必须的，需要根据原数据矩阵看，要找的是包含原矩阵90%能量的个数的奇异值。x_formed_items一句，假设原数据矩阵的shape是（m,n)，则u[:,4]的shape是（m，4），sig4的shape是（4，4），逆也是（4，4），想成后得到的x_formed_items的shape是（n，4），原来的n变成了行，我们求相似度传的是一个个的列向量，所以转置。
推荐
流程图：
推荐流程图.png
# data_mat:数据矩阵
# user:用户编号
# N: 要返回的前N个要推荐的items
# sim_meas: 相似度方法
# est_method:评分预估法
def recommend(data_mat, user, N, sim_meas, est_method):
un_rated_items = np.nonzero(data_mat[user, :].A == 0)[1] # 上面说过，nonzero的第1个元素是第1轴的下标，这里也就是列下标
if len(un_rated_items) == 0:
return 'you rated everything'
item_scores = []
for item in un_rated_items:
estimate_score = est_method(data_mat, user, item, sim_meas)
item_scores.append((item, estimate_score))
return sorted(item_scores, key=lambda jj:jj[1], reverse=True)[:N] # 排序，key=lambda jj:jj[1]表示按每个元素的下标为1的参数从大到小排序，取前N个
评价
由于推荐引擎建好后既没有预测的目标值，也没有用户来调查他们对推荐的满意程度，所以常常将某些已知的评分值去掉，然后对它们进行预测，计算预测值和真实值之间的差异。
通常用于推荐引擎评价的指标是最小均方根误差（Root Mean Squared Error，RMSE），先计算均方误差的平均值，在取平方根。如果评级在1～5星级，我们得到的RMSE为1，则和真实评价差了1个星级。
