矩阵分解是推荐系统常用的手段,经常用来做用户偏好预测.在当下的推荐系统中,我们得到用户对于物品的评分矩阵往往是非常稀疏的,一个有m个用 户,n个商品的网站,它所收集到的m*n用户评分矩阵R可能只有不到万分之一的数据非零.矩阵分解算法常用来构造出多个矩阵, 用这些矩阵相乘的结果R’来拟合原来的评分矩阵R，目标是使得到的矩阵R’在R的非零元素那些位置上的值尽量接近R中的元素，同时对于R中非零值进行补 全.我们定义R和R’之间的距离,把它作为优化的目标,那么矩阵分解就变成了最优化问题,而这类最优化问题常用梯度下降的方法来求出一个局部最优解。最近 我们就对于这个问题进行了一下预研，发现用分布式图计算框架来解决这类问题比较方便流畅,尤其是在矩阵非常稀疏的时候..
1:SGD(Stochastic Gradient Descent)随机梯度下降的公式推导
设评分矩阵R中的每个不为0的元素为R_ui,所有的不为0的元素对应的下标对(u,i)组成的集合为S,我们认为存在一个K维的隐空间,用户偏好在隐空间可以表示为一个m k的矩阵P,物品偏好在隐空间表示为一个k n的矩阵Q. P是m k的矩阵,每一行p_u表示用户u在K维隐空间上的偏好, Q是k n的矩阵,每一列表示是的是物品i在K维空间上的特性. 第u个用户的偏好可以记为p_u ,是一个1 k的向量,第i个商品的偏好记为q_i,是一个k 1的向量,用户u对于物品i的打分预估值就是p_u*q_i,整个预测的误差可以记为
(1)式可以理解为对于所有有预测的数据记录,评分预测值和实际值的差值平方和,通常用这个做为最初的损失函数(loss function).我们希望任何一个p_u, q_i的值都不能过大,这样容易出现过拟合的现象,所以我们给(1)加上一个正则化项,于是(1)式变为 
其中 
分别是p_u,q_i 的二范数,之所以选择二范数是因为我们认为用户偏好，物品特性这些维度上的分值分布是符合正态分布的，这种情况下常使用二范数，即向量所有元素的平方和. （对于L0,L1,L2的深入探讨详见参考资料1）.我们的目标就是找到矩阵P和矩阵Q,使得(2)的值最小.因为(2)是左边是个连续函数,所以极值出 现在梯度为0的时候,所以我们可以通过求梯度的方法得到使得L最小的Q和P
(3)式就是L对于p_u的梯度R(u)表示所有用户u评价过的物品,从(3)式中可以看出,一个用户的在K维空间的向量,只和用户自身的以及他 评价过的物品的K维空间上的向量有关.(3)的结果是一个K维的向量,每个位置上的元素对应的是用户u在K维空间的向量对应的梯度值 类似地,我们可以得到
通过对于(3) (4)两个式子的对比我们可以发现它们的形式是一致的,唯一不同的就是p,q呼唤了位置,一个是在所有的用户点上通过相邻物品点的值进行计算 (i∈R(u)),一个是在所有的物品点上通过相邻用户点的值进行计算(u∈R(i)) .和(3)式的计算结果类似,(4)式也是一个K维的向量,每个位置上的元素对应的是物品i在K维空间的值对应的梯度值. 如果我们把学习率记为 
整个迭代的过程可以看作在一个K维超空间的曲面上寻找极值点的过程,每一步都寻找一个比之前更低的一个位置,如下图所示,图中的θ1,θ2可以认 为相当于P和Q,差别就是一个是值,一个是矩阵,但是他们都是某个凸函数的输入值.从另一个角度来说,一个值也可以看成是一个1 * 1的矩阵,如果P和Q退化为一个1 * 1的矩阵的话,那就是θ1,θ2了.
上面的说法直接理解可能比较抽象,我们用一个例子来说明.图1说明了梯度下降就是在一个超空间的曲面上寻找一个点，这个点对应的J(θ),即损失 函数(loss function)最小，步长η就是上图每条线段的长度，-∂L/∂Pu表示了接下来往哪个方向走。图一表示的是一个梯度下降的过程得到了全局最优解。但 是这个情况还是比较少见的，更多的是像图二中的情况，得到一个局部最优解。为了避免出现局部最优解比较差的情况，我们通常可以多次随机地初始化Q和P矩 阵，相当于是选择不同的初始点开始进行梯度下降，选择一个损失函数(loss function)最小的那一次对应的P,Q值作为梯度下降的结果。
2:SGD的一个例子说明
下图是我目前得到的一个评分文件,3列的含义分别是UID:User ID,IID:Item ID,score:用户评分.可以看到一共有3个用户,4个物品. 
他们可以构成一个3 * 4的评分矩阵矩阵.我现在取k=2,要把它们分解成为一个3 2的P矩阵和一个2 4的Q矩阵.
首先初始化P和Q矩阵,一般都用符合正态分布的随机数X~N(0,1)来填充矩阵.
现在我计算u=1时的∂L/∂Pu,取λ的值为1 首先计算u=1,i=1这个评分记录带来的梯度分量,即u=1,i=1时的 
的值,这时 
然后计算u=1,i=2这个评分带来的梯度分量: 
所以,u=1时的梯度为
刚才提到的迭代公式: p_u=p_u-η ∂L/∂Pu,其中η表示学习率,就是平时我们说的梯度下降时的步长,取η=0.05 于是有 
可见,通过这一次对于p_1的梯度下降, p_1对应的向量从随机产生的 
我们来验证一下这个变化是否是的预估的值更加准确了 原来对于R_11的估计误差为 
新的对于R_11的估计误差为 
确实比原来有提升.按照上述步骤分别更新所有P矩阵中的p_i,然后用更新好的P矩阵,用公式(4)中的方法更新Q,再用更新好的Q第二次更新P.这种迭代方法最终可以使得|R-P Q|的值缩小,当|R-P Q|在下一轮迭代的过程中不再变小时,我们就认为在当前的学习率η下,P,Q已经得到了局部最优解.
3:用图计算框架来解SGD
通过刚才的演示你可以发现,每次某个用户的K维向量需要更新的时候,只有这个用户评分过的商品对应的K维向量会参加运算, 和其他的元素无关.而且整个评分矩阵在现实应用中往往是非常稀疏的.这种计算过程和数据分布的特点使得这类问题用图计算框架来解决比起传统的 mapreduce过程更加流畅高效. 对于表1的评分矩阵我们可以用一个二部图来表示它,每个点分别表示一个用户或者一个物品,用户对于物品的评分就是用户点和物品点之间边的属性,这个二部图 如图三. 
图三:根据评分矩阵生成二部图
图四:每个点的K维向量初始化后的效果演示.
图五:SGD第一步数据在图计算框架中的传输演示.
正如上文说的,我们要初始化Q和P矩阵,实际上就是在每个点上初始化一个K维的数组,这个数组表示了这个点的在K维空间的映射结果. . 用图计算框架来做SGD遵循以下步骤来进行迭代: 第零步:初始化,生成可以表达Q矩阵,P矩阵和评分矩阵的图数据结构(如图四),令L=0 第一步,每个用户节点接受边上的值以及每个和这个用户评过分的物品的K维向量,如图五所示. 第二步,每个用户点根据自己每条边上传来的值计算条边和边上的邻居造成的梯度分量,将这些梯度分量相加,得到这轮迭代的梯度值. 第三步,每个用户点将自己得到的梯度值,根据公式p_u=p_u-η ∂L/∂Pu更新自己的K维向量,也就是P矩阵得到了更新. 第四步,每个物品点接受边上的评分值以及每个和对于这个物品进行过评分的所有用户的K维向量. 第五步, 每个物品点根据自己每条边上传来的值计算条边和边上的邻居造成的梯度分量,将这些梯度分量相加,得到这轮迭代的梯度值. 第六步, 每个物品点将自己得到的梯度值,根据公式q_i=q_i-η ∂L/∂qi更新自己的K维向量,也就是Q矩阵得到了更新. 第七步,根据 
检测当前的L值是否比上次的L小,如果是,回到第一步继续迭代.如果否则表示迭代介绍,按顺序遍历图中的用户和物品节点,读取每个节点的K维向量,组合后输出矩阵Q和P.
图六:用SGD做矩阵分解流程图
4:通过逐步减小参与计算的数据量加速
在上面的介绍中我们可以看到整个的迭代都是矩阵Q和矩阵P不断变化,使得P*Q的结果接近R的过程。对于某些用户和物品节点而言,这些点的K维向量很快就收敛到了我们的目标值,即对于这些用户节点而言,
对于某些点,误差值在迭代开始没几轮的时候就达到了目标范围,这个某些点提前结束迭代的功能对于物品节点也成立.这些点可以不用参加之后的迭代,在基于spark的图计算框架graphx提供了这样功能的函数：
mapReduceTriplets(sendMsg, mergeMsg, Some((newVerts, activeDir))
这个函数是GraphX中最核心的一个函数，表示对于图中的所有以两个点一条边组成的三元组Triplets进行一次数据传输与计算，具体过程 为：sendMsg，相当于map，应用于每一个Triplet上，生成一个或者多个消息，消息以Triplet关联的两个顶点中的任意一个或两个为目标 顶点；mergeMsg，相当于reduce，应用于每一个Vertex上，将发送给每一个顶点的消息合并起来。（对于内存式图计算框架graphx更详 细的介绍可以见参考资料2），前面的两个参数别是在”think like vertex”的情况下每个点对于邻居节点发出的信息sendMsg和对于自己的邻居发给自己的信息的处理函数mergeMsg .而第三个参数确定了满足一定条件的点参与这一次的图的迭代. GraphX针对这种应用场景进行了较为深层的优化 没有更新的顶点在下一轮迭代时不需要向邻居重新发送消息。因此，mapReduceTriplets遍历边时，如果一条边邻居点值在上一轮迭代时没有更 新，则直接跳过，避免了大量无用的计算和通信。 这样做的好处是，随着收敛的点越来越多,每轮迭代参与计算的点越来越少,每次迭代需要的通讯和计算开销都变小了,每次迭代的时间也逐步减小,如下图所示:
图七:在Netflix数据集上用SGD做矩阵分解耗时图 整个的算法效果可以从上图中看到的，随着迭代次数增加，整个每次迭代参与的用户和物品点变少，耗时迅速减少，整个算法很快就收敛了，这正是使用图计算框 架，尤其是经过底层优化后的图计算框架（比如graphx）的优势所在。更详细的介绍巨型图的存储以及基于分布式图存储，图计算框架并行计算时的通信方式 的优化可以看参考资料3，4。
5:结论
对于很多的迭代问题，尤其是涉及到矩阵分解的问题，使用图计算框架如graphx比起使用传统的mapreduce不论是在编码效率还是执行速度 上都有着有明显的优势。图计算框架不仅仅是用来解决“图”，即所谓网络科学的问题，它在传统的大规模机器学习领域也有这广泛的应用场景。除了文中演示的矩 阵分解问题，概率图模型（如PLSA，LDA）用图计算框架也可以完成。据说google的图计算框架pregel承担了google20%的计算任务， 可惜它没有开源，不能一看究竟。开源图计算框架的性能提升，应用场景扩展有着很多潜力供将来深入挖掘。
参考资料
1： csdn博文：机器学习中的范数规则化之（一）L0、L1与L2范数
2： 程序员杂志2014年9月刊： 快刀初试：Spark GraphX在淘宝的实践
3： 淘宝技术部博客：关于图计算和graphx的一些思考
4： GraphLab:新的面向机器学习的并行框架
