本文以經典的協同過濾為切入點，重點介紹了被工業界廣泛使用的矩陣分解演算法，從理論與實踐兩個維度介紹了該演算法的原理，通俗易懂，希望能夠給大家帶來一些啟發。 ...
作者：vivo 互聯網伺服器團隊-Tang Shutao
現如今推薦無處不在，例如抖音、淘寶、京東App均能見到推薦系統的身影，其背後涉及許多的技術。本文以經典的協同過濾為切入點，重點介紹了被工業界廣泛使用的矩陣分解演算法，從理論與實踐兩個維度介紹了該演算法的原理，通俗易懂，希望能夠給大家帶來一些啟發。筆者認為要徹底搞懂一篇論文，最好的方式就是動手復現它，復現的過程你會遇到各種各樣的疑惑、理論細節。
一、 背景
1.1 引言
在信息爆炸的二十一世紀，人們很容易淹沒在知識的海洋中，在該場景下搜索引擎可以幫助我們迅速找到我們想要查找的內容。在電商場景，如今的社會物質極大豐富，商品琳琅滿目，種類繁多。消費者很容易挑花眼，即用戶將會面臨信息過載的問題。為瞭解決該問題，推薦引擎應運而生。例如我們打開淘寶App，JD app，B站視頻app，每一個場景下都有推薦的模塊。那麼此時有一個幼兒園小朋友突然問你，為什麼JD給你推薦這本《程式員頸椎康復指南》？你可能會回答，因為我的職業是程式員。接著小朋友又問，為什麼《Spark大數據分析》這本書排在第6個推薦位，而《Scala編程》排在第2位？這時你可能無法回答這個問題。
為了回答該問題，我們設想下麵的場景：
在JD的電商系統中，存在著用戶和商品兩種角色，並且我們假設用戶都會對自己購買的商品打一個0-5之間的分數，分數越高代表越喜歡該商品。
基於此假設，我們將上面的問題轉化為用戶對《程式員頸椎康復指南》，《Spark大數據分析》，《Scala編程》這三本書打分的話，用戶會打多少分（用戶之前未購買過這3本書）。因此物品在頁面的先後順序就等價於預測用戶對這些物品的評分，並且根據這些評分進行排序的問題。
為了便於預測用戶對物品的評分問題，我們將所有三元組(User, Item, Rating)，即用戶User給自己購買的商品Item的評分為Rating，組織為如下的矩陣形式：
其中，表格包含\(m\)個用戶和\(n\)個物品，將表格定義為\({\bf{R}}_{m \times n}\)評分矩陣_，其中的元素\(r_{u,i})\)表示第\(u\)個用戶對第\(i\)個物品的評分。
例如，在上面的表格中，用戶user-1購買了物品 item-1， item-3， item-4，並且分別給出了4，2，5的評分。最終，我們將原問題轉化為預測白色空格處的數值。
1.2 協同過濾
協同過濾，簡單來說是利用與用戶興趣相投、擁有共同經驗之群體的喜好來推薦給用戶感興趣的物品。興趣相投使用數學語言來表達就是相似度 (人與人，物與物)。因此，根據相似度的對象，協同過濾可以分為基於用戶的協同過濾和基於物品的協同過濾。
以評分矩陣為例，以行方向觀測評分矩陣，每一行代表每個用戶的向量表示，例如用戶user-1的向量為 4, 0, 2, 5, 0, 0。以列方向觀測評分矩陣，每一列表示每個物品的向量表示，例如物品item-1的向量為4, 3, 0, 0, 5。
基於向量表示，相似度的計算有多種公式，例如餘弦相似度，歐氏距離，皮爾森。這裡我們以餘弦相似度為例，它是我們中學學過的向量夾角 (中學只涉及2維和3維) 的高維推廣，餘弦相似度公式很容易理解和使用。給定兩個向量\(\mathbf{A}=\{a_1, \cdots, a_n\}\)和\(\mathbf{B}=\{b_1, \cdots, b_n\}\)，其夾角定義如下：
\(\cos(\theta)=\frac{\bf{A}\cdot \bf{B}}{{|\bf{A}}|{|\bf{B}}|}=\frac{a_1 b_1 + \cdots + a_n b_n}{\sqrt{a_1^2+\cdots a_n^2}\sqrt{b_1^2 + \cdots b_n^2}}\)
例如，我們計算user-3和user-4的餘弦相似度，二者對應的向量分別為 0, 2, 0, 3, 0, 4，0, 3, 3, 5, 4, 0
\(\text{cos_sim}(u_3, u_4)=\frac{2 \times 3 + 3 \times 5}{\sqrt{2^2+3^2+4^2}\sqrt{3^2+3^2+5^2+4^2}} \approx 0.507685\)
向量夾角的餘弦值越接近1代表兩個物品方向越接近平行，也就是越相似，反之越接近-1代表兩個物品方向越接近反向，表示兩個物品相似度接近相反，接近0，表示向量接近垂直/正交，兩個物品幾乎無關聯。顯然，這和人的直覺完全一致。
例如，我們在視頻App中經常能看到"相關推薦"模塊，其背後用到的原理之一就是相似度計算，下圖展示了一個具體的例子。
我們用《血族第一部》在向量庫 (存儲向量的資料庫，該系統能夠根據輸入向量，用相似度公式在庫中進行檢索，找出TopN的候選向量) 裡面進行相似度檢索，找到了前7部高相似度的影片，值得註意的是第一部是自己本身，相似度為1.0，其他三部是《血族》的其他3部同系列作品。
1.2.1 基於用戶的協同過濾 (UserCF)
基於用戶的協同過濾分為兩步
找出用戶相似度TopN的若幹用戶。
根據TopN用戶評分的物品，形成候選物品集合，利用加權平均預估用戶u對每個候選物品的評分。
例如，由用戶u的相似用戶{u1, u3, u5, u9}可得候選物品為
\(\{i_1, i_2, i_3, i_4, i_5, i_6, i_7\}\)
我們現在預測用戶u對物品i1的評分，由於物品在兩個用戶{u1, u5}的購買記錄里，因此用戶u對物品i1的預測評分為：
\({r}_{u,i_1} = \frac{\text{sim}{(u,u_1)}\times {r}_{u_1,i_1}+\text{sim}{(u,u_5)}\times {r}_{u_5,i_1}}{\text{sim}{(u,u_1)}+\text{sim}{(u,u_5)}}\)
其中\(\text{sim}{(u,u_1)}\)表示用戶\(u\)與用戶\(u_1\)的相似度。
在推薦時，根據用戶u對所有候選物品的預測分進行排序，取TopM的候選物品推薦給用戶u即可。
1.2.2 基於物品的協同過濾 (ItemCF)
基於物品的協同過濾分為兩步：
在用戶u購買的物品集合中，選取與每一個物品TopN相似的物品。
TopN相似物品形成候選物品集合，利用加權平均預估用戶u對每個候選物品的評分。
例如，我們預測用戶u對物品i3的評分，由於物品i3與物品{i6, i1, i9}均相似，因此用戶u對物品i3的預測評分為：
\({r}_{u,i_3} = \frac{\text{sim}{(i_6,i_3)}\times {r}_{u,i_6}+\text{sim}{(i_1,i_3)}\times {r}_{u,i_1}+\text{sim}{(i_9,i_3)}\times {r}_{u,i_9}}{\text{sim}{(i_6,i_3)}+\text{sim}{(i_1,i_3)}+\text{sim}{(i_9,i_3)}}\)
其中\(\text{sim}{(i_6,i_3)}\)表示物品\(i_6\)與物品\(i_3\)的相似度，其他符號同理。
1.2.3 UserCF與ItemCF的比較
我們對ItemCF和UserCF做如下總結：
UserCF主要用於給用戶推薦那些與之有共同興趣愛好的用戶喜歡的物品，其推薦結果著重於反映和用戶興趣相似的小群體的熱點，更社會化一些，反映了用戶所在的小型興趣群體中物品的熱門程度。在實際應用中，UserCF通常被應用於用於新聞推薦。ItemCF給用戶推薦那些和他之前喜歡的物品類似的物品，即ItemCF的推薦結果著重於維繫用戶的歷史興趣，推薦更加個性化，反應用戶自己的興趣。在實際應用中，圖書、電影平臺使用ItemCF，比如豆瓣、亞馬遜、Netflix等。除了基於用戶和基於物品的協同過濾，還有一類基於模型的協同過濾演算法，如上圖所示。此外基於用戶和基於物品的協同過濾又可以歸類為基於鄰域 (K-Nearest Neighbor, KNN) 的演算法，本質都是在找"TopN鄰居"，然後利用鄰居和相似度進行預測。
二、矩陣分解
經典的協同過濾演算法本身存在一些缺點，其中最明顯的就是稀疏性問題。我們知道評分矩陣是一個大型稀疏矩陣，導致在計算相似度時，兩個向量的點積等於0 (以餘弦相似度為例)。為了更直觀的理解這一點，我們舉例如下：
我們從評分矩陣中抽取item1 - item4的向量，並且利用餘弦相似度計算它們之間的相似度
rom sklearn.metrics.pairwise import cosine_similarity
a = [
[ 0, 0, 0, 3, 2, 0, 3.5, 0, 1 ],
[ 0, 1, 0, 0, 0, 0, 0, 0, 0 ],
[ 0, 0, 1, 0, 0, 0, 0, 0, 0 ],
[4.1, 3.8, 4.6, 3.8, 4.4, 3, 4, 0, 3.6]
]
cosine_similarity(a)
# array([[1. , 0. , 0. , 0.66209271],
# [0. , 1. , 0. , 0.34101639],
# [0. , 0. , 1. , 0.41280932],
# [0.66209271, 0.34101639, 0.41280932, 1. ]])
我們從評分矩陣中抽取item1 - item4的向量，並且利用餘弦相似度計算它們之間的相似度。
通過相似度矩陣，我們可以看到物品item-1, item-2, item-3的之間的相似度均為0，而且與item-1, item-2, item-3最相似的物品都是item-4，因此在以ItemCF為基礎的推薦場景中item-4將會被推薦給用戶。
但是，物品item-4與物品item-1, item-2, item-3最相似的原因是item-4是一件熱門商品，購買的用戶多，而物品item-1, item-2, item-3的相似度均為0的原因僅僅是它們的特征向量非常稀疏，缺乏相似度計算的直接數據。
綜上，我們可以看到經典的基於用戶/物品的協同過濾演算法有天然的缺陷，無法處理稀疏場景。為瞭解決該問題，矩陣分解被提出。
2.1 顯示反饋
我們將用戶對物品的評分行為定義為顯示反饋。基於顯示反饋的矩陣分解是將評分矩陣\({\bf{R}}_{m \times n}\)用兩個矩陣\({\bf{X}}_{m \times k}\)和\({\bf{Y}}_{n \times k}\)的乘積近似表示，其數學表示如下：
\({\bf{R}}_{m \times n} \approx {\bf{X}}_{m \times k}\left({\bf{Y}}_{n \times k}\right)^{\text T}\)
其中，\(k \ll m/n\)表示隱性因數，以用戶側來理解，\(k=2\)表示的就是用戶的年齡和性別兩個屬性。此外有個很好的比喻就是物理學的三棱鏡，白光在三棱鏡的作用下被分解為7種顏色的光，在矩陣分解演算法中，分解的作用就類似於"三棱鏡"，如下圖所示，因此，矩陣分解也被稱為隱語義模型。矩陣分解將系統的自由度從\(\mathcal{O}(mn)\)降到了\(\mathcal{O}((m+n)k)\)，從而實現了降維的目的。
為了求解矩陣\({\bf{X}}_{m \times k}\)和\({\bf{Y}}_{n \times k}\)，需要最小化平方誤差損失函數，來儘可能地使得兩個矩陣的乘積逼近評分矩陣\({\bf{R}}_{m \times n}\)，即
\(\min\limits_{{\bf{x}}^*,{\bf{y}}^*} L({\bf{X}},{\bf{Y}})=\min\limits_{{\bf{x}}^*,{\bf{y}}^*}\sum\limits_{r_{u,i} \text{ is known}}(r_{u,i}-{\bf{x}}_u^{\text T}{\bf{y}}_i)^2+\lambda \left( \sum\limits_{u}{\bf{x}}_u^{\text T}{\bf{x}}_u+\sum\limits_{i}{\bf{y}}_i^{\text T}{\bf{y}}_i\right)\)
其中，\(\lambda \left( \sum\limits_{u}{\bf{x}}_u^{\text T}{\bf{x}}_u+\sum\limits_{i}{\bf{y}}_i^{\text T}{\bf{y}}_i\right)\)為懲罰項，\(\lambda\)為懲罰繫數/正則化繫數，\(\mathbf{x}_u\)表示第\(u\)個用戶的\(k\)維特征向量，\(\mathbf{y}_i\)表示第\(i\)個物品的\(k\)維特征向量。
\({\bf{x}}_u = \begin{pmatrix} x_{u,1} \ \vdots \ x_{u,k} \ \end{pmatrix} \qquad {\bf{y}}_i = \begin{pmatrix} y_{i,1} \ \vdots \ y_{i,k} \ \end{pmatrix}\)
全體用戶的特征向量構成了用戶矩陣\({\bf{X}}_{m \times k}\)_，全體物品的特征向量構成了物品矩陣\({\bf{Y}}_{n \times k}\)。
\({\bf{X}}_{m \times k}= \begin{pmatrix} {\bf{x}}_{1}^{\text T} \ \vdots \ {\bf{x}}_{m}^{\text T} \ \end{pmatrix} \qquad {\bf{Y}}_{n \times k}= \begin{pmatrix} {\bf{y}}_{1}^{\text T} \ \vdots \ {\bf{y}}_{n}^{\text T} \ \end{pmatrix}\)
我們訓練模型的時候，就只需要訓練用戶矩陣中的\(m \times k\)參數和物品矩陣中的\(n \times k\)個參數。因此，協同過濾就成功轉化成了一個優化問題。
2.2 預測評分
通過模型訓練 (即求解模型繫數的過程)，我們得到用戶矩陣\({\bf{X}}_{m \times k}\)和物品矩陣\({\bf{Y}}_{n \times k}\)，全部用戶對全部物品的評分預測可以通過\({\bf{X}}_{m \times k}\left({\bf{Y}}_{n \times k}\right)^{\text T}\)獲得。如下圖所示。
得到全部的評分預測後，我們就可以對每個物品進行擇優推薦。需要註意的是，用戶矩陣和物品矩陣的乘積，得到的評分預估值，與用戶的實際評分不是全等關係，而是近似相等的關係。如上圖中兩個矩陣粉色部分，用戶實際評分和預估評分都是近似的，有一定的誤差。
2.3 理論推導
矩陣分解ALS的理論推導網上也有不少，但是很多推導不是那麼嚴謹，在操作向量導數時有的步驟甚至是錯誤的。有的博主對損失函數的求和項理解出現錯誤，例如
\(\sum\limits_{\color{red}{u=1}}^{\color{red} m}\sum\limits_{\color{red}{i=1}}^{\color{red} n}(r_{u,i}-{\bf{x}}_u^{\text T}{\bf{y}}_i)^2\)
但是評分矩陣是稀疏的，求和並不會貫穿整個用戶集和物品集。正確的寫法應該是
\(\sum\limits_{\color{red}{(u,i) \text{ is known}}}(r_{u,i}-{\bf{x}}_u^{\text T}{\bf{y}}_i)^2\)
其中，\({(u,i) \text{ is known}}\)表示已知的評分項。
我們在本節給出詳細的、正確的推導過程，一是當做數學小練習，其次也是對演算法有更深層的理解，便於閱讀Spark ALS的源碼。
將\({(u,i) \text{ is known}}\)使用數學語言描述，矩陣分解的損失函數定義如下：
\(L({\bf{X}},{\bf{Y}})=\sum\limits_{\color{red}{(u,i) \in K}}(r_{u,i}-{\bf{x}}_u^{\text T}{\bf{y}}_i)^2+\lambda \left( \sum\limits_{u}{\bf{x}}_u^{\text T}{\bf{x}}_u+\sum\limits_{i}{\bf{y}}_i^{\text T}{\bf{y}}_i\right)\)
其中\(K\)為評分矩陣中已知的\((u, i)\)集合。例如下麵的評分矩陣對應的\(K\)為
\({\bf{R}}_{4 \times 4} = \begin{pmatrix} 0 & r{1,2} & r{1,3} & 0 \ r{2,1} & 0 & r{2,3} & 0 \ 0 & r{3,2} & 0 & r{3,4} \ 0 & r{4,2} & r{4,3} & r_{4,4} \end{pmatrix} \ \Rightarrow \color{red}{K = {(1,2), (1,3), (2,1), (2,3), (3,2), (3,4), (4,2), (4,3), (4,4)}}\)
求解上述損失函數存在兩種典型的優化方法，分別為
- 交替最小二乘 (Alternating Least Squares, ALS)
- 隨機梯度下降 (Stochastic Gradient Descent, SGD)
交替最小二乘，指的是固定其中一個變數，利用最小二乘求解另一個變數，以此交替進行，直至收斂或者到達最大迭代次數，這也是“交替”一詞的由來。
隨機梯度下降，是優化理論中最常用的一種方式，通過計算梯度，然後更新待求的變數。
在矩陣分解演算法中，Spark最終選擇了ALS作為官方的唯一實現，原因是ALS很容易實現並行化，任務之間沒有依賴。
下麵我們動手推導一下整個計算過程，在機器學習理論中，微分的單位一般在向量維度，很少去對向量的分量為偏微分推導。
首先我們固定物品矩陣\({\bf{Y}}\)，將物品矩陣\({\bf{Y}}\)看成常量。不失一般性，我們定義用戶\(u\)評分過的物品集合為\(I_u\)，利用損失函數對向量\(\mathbf{x}_u\)求偏導，並且令導數等於0可得：
\(\displaystyle \frac{\partial L}{\partial {\bf{x}}u}=-2\sum\limits{i \in Iu}(r{u,i}-{\bf{x}}_u^{\text T}{\bf{y}}_i)\frac{\partial {(\bf{x}}_u^{\text T}{\bf{y}}_i)}{\partial {\bf{x}}_u}+2\lambda \frac{\partial {(\bf{x}}_u^{\text T}{\bf{x}}_u)}{\partial {\bf{x}}_u}=0, \quad u=1, \cdots, m \ \begin{split} & \quad \Rightarrow \sum\limits{i \in I_u}(r{u,i}-{\bf{x}}_u^{\text T}{\bf{y}}_i){\bf{y}}_i^{\text T}=\lambda {\bf{x}}_u^{\text T} \ & \quad \Rightarrow \sum\limits{i \in I_u}r{u,i}{\bf{y}}i^{\text T}-\sum\limits{i \in I_u}{\bf{x}}_u^{\text T}{\bf{y}}_i{\bf{y}}_i^{\text T}=\lambda {\bf{x}}_u^{\text T} \ & \quad \Rightarrow \sum\limits{i \in I_u}{\bf{x}}_u^{\text T}{\bf{y}}_i{\bf{y}}_i^{\text T}+\lambda {\bf{x}}_u^{\text T}=\sum\limits{i \in Iu}r{u,i}{\bf{y}}_i^{\text T} \end{split}\)
因為向量\(\mathbf{x}_u\)與求和符號\(\sum\limits_{i \in I_u}\)無關，所有將其移出求和符號，因為\({\bf{x}}_u^{\text T}{\bf{y}}_i{\bf{y}}_i^{\text T}\)是矩陣相乘 (不滿足交換性)，因此\(\mathbf{x}_u\)在左邊\(\begin{split} & {\bf{x}}u^{\text T}\sum\limits{i \in Iu}{\bf{y}}_i{\bf{y}}_i^{\text T}+\lambda {\bf{x}}_u^{\text T}=\sum\limits{i \in Iu}r{u,i}{\bf{y}}_i^{\text T} \ \end{split}\) 等式兩邊取轉置，我們有
\(\begin{split} & \quad \Rightarrow \left({\bf{x}}u^{\text T}\sum\limits{i \in Iu}{\bf{y}}_i{\bf{y}}_i^{\text T}+\lambda {\bf{x}}_u^{\text T}\right)^{\text T}=\left(\sum\limits{i \in Iu}r{u,i}{\bf{y}}_i^{\text T}\right)^{\text T} \ & \quad \Rightarrow \left(\sum\limits{i \in I_u}{\bf{y}}_i{\bf{y}}_i^{\text T}\right){\bf{x}}_u+\lambda {\bf{x}}_u=\sum\limits{i \in Iu}r{u,i}{\bf{y}}_i \ & \quad \Rightarrow \left(\sum\limits{i \in I_u}{\bf{y}}_i{\bf{y}}_i^{\text T}+\lambda {\bf{I}}_k\right){\bf{x}}_u=\sum\limits{i \in Iu}r{u,i}{\bf{y}}_i \end{split}\)
為了化簡\(\sum\limits_{i \in I_u}{\bf{y}}_i{\bf{y}}_i^{\text T}\)與\(\sum\limits_{i \in I_u}r_{u,i}{\bf{y}}_i\)，我們將\(I_u\)展開。
假設\(I_u=\{i_{c_1}, \cdots, i_{c_N}\}\), 其中\(N\)表示用戶\(u\)評分過的物品數量，\(i_{c_i}\)表示第\(c_i\)個物品對應的索引/序號，藉助於\(I_u\)，我們有
\(\sum\limits_{i \in I_u}{\bf{y}}_i{\bf{y}}_i^{\text T}= \begin{pmatrix} {\bf{y}}_{c_1}, \cdots,{\bf{y}}_{c_N} \end{pmatrix} \begin{pmatrix} {\bf{y}}_{c_1}^{\text T} \\ \vdots \\{\bf{y}}_{c_N}^{\text T} \end{pmatrix}={\bf{Y}}_{I_u}^{\text T}{\bf{Y}}_{I_u} \\ \sum\limits_{i \in I_u}r_{u,i}{\bf{y}}_i= \begin{pmatrix}{\bf{y}}_{c_1}, \cdots,{\bf{y}}_{c_N} \end{pmatrix} \begin{pmatrix} r_{u,c_1} \\ \vdots \\ r_{u,c_N} \end{pmatrix}={\bf{Y}}_{I_u}^{\text T}{\bf{R}}_{u,I_u}^{\text T}\)
其中，\({\bf{Y}}_{I_u}\)為以\(I_u=\{i_{c_1}, \cdots i_{c_N}\}\)為行號在物品矩陣\({\bf{Y}}\)中選取的\(N\)個行向量形成的子矩陣
\({\bf{R}}_{u,I_u}\)為以\(I_u=\{i_{c_1}, \cdots i_{c_N}\}\)為索引，在評分矩陣\({\bf{R}}\)的第\(u\)行的行向量中選取的\(N\)個元素，形成的子行向量
因此，我們有
\(\left({\bf{Y}}{I_u}^{\text T}{\bf{Y}}{Iu}+\lambda {\bf{I}}_k\right){\bf{x}}_u={\bf{Y}}{Iu}^{\text T}{\bf{R}}{u,I_u}^{\text T} \ \quad \Rightarrow {\bf{x}}u=\left({\bf{Y}}{Iu}^{\text T}{\bf{Y}}{Iu}+\lambda {\bf{I}}_k\right)^{-1}{\bf{Y}}{Iu}^{\text T}{\bf{R}}{u,I_u}^{\text T}\)
網上的博客，許多博主給出類似下麵形式的結論不是很嚴謹，主要是損失函數的理解不到位導致的。
\({\bf{x}}_u=\left({\bf{\color{red} Y}}^{\text T}{\bf{\color{red} Y}}+\lambda {\bf{I}}_k\right)^{-1}{\bf{\color{red} Y}}^{\text T}{\bf{\color{red} R}}_{u}^{\text T}\)
同理，我們定義物品\(i\)被評分的用戶集合為\(U_i=\{u_{d_1}, \cdots u_{d_M}\}\)
根據對稱性可得
\({\bf{y}}_i=\left({\bf{X}}_{U_i}^{\text T}{\bf{X}}_{U_i}+\lambda {\bf{I}}_k\right)^{-1}{\bf{X}}_{U_i}^{\text T}{\bf{R}}_{i,U_i}\)
其中，\({\bf{X}}_{U_i}\)為以\(U_i=\{u_{d_1}, \cdots, u_{d_M}\}\)為行號在用戶矩陣\({\bf{X}}\)中選取的\(M\)個行向量形成的子矩陣
\({\bf{R}}_{i,U_i}\)為以\(U_i=\{u_{d_1}, \cdots, u_{d_M}\}\)為索引，在評分矩陣\({\bf{R}}\)的第\(i\)列的列向量中選取的\(M\)個元素，形成的子列向量
此外，\(\mathbf{I}_k\)為單位矩陣
如果讀者感覺上述的推導還是很抽象，我們也給一個具體實例來體會一下中間過程
\(\begin{pmatrix} 0 & r_{1,2} & r_{1,3} & 0 \\ r_{2,1} & 0 & r_{2,3} & 0 \\ 0 & r_{3,2} & 0 & r_{1,3} \\ 0 & r_{2,2} & r_{2,3} & r_{2,4} \\ \end{pmatrix} \approx \begin{pmatrix} x_{1,1} & x_{1,2} \\ x_{2,1} & x_{2,2} \\ x_{3,1} & x_{3,2} \\ x_{4,1} & x_{4,2} \end{pmatrix} \begin{pmatrix} y_{1,1} & y_{1,2} \\ y_{2,1} & y_{2,2} \\ y_{3,1} & y_{3,2} \\ y_{4,1} & y_{4,2} \end{pmatrix}^{\text T} \\ \Rightarrow {\bf{R}} \approx {\bf{X}} {\bf{Y}}^{\text T}\)
註意到損失函數是一個標量，這裡我們只展開涉及到\(x_{1,1}, x_{1,2}\)的項，如下所示
\(\sum\limits_{\color{red}{(u,i) \text{ is known}}}(r_{u,i} - {\bf{x}}_u^{\text{T}}{\bf{y}}_i)^2 =(\color{blue}{r_{1,2}} - \color{red}{x_{1,1}}y_{2,1} - \color{red}{x_{1,2}}y_{2,2})^2 + (\color{blue}{r_{1,3}} - \color{red}{x_{1,1}}y_{3,1} - \color{red}{x_{1,2}}y_{3,2})^2 + \cdots\)
讓損失函數對\(x_{1,1}, x_{1,2}\)分別求偏導數可以得到
\(\frac{\partial{L}}{\partial{\color{red}{x_{1,1}}}}=2(\color{blue}{r_{1,2}} - \color{red}{x_{1,1}}y_{2,1}-\color{red}{x_{1,2}}y_{2,2})(-y_{2,1}) + 2(\color{blue}{r_{1,3}} - \color{red}{x_{1,1}}y_{3,1}-\color{red}{x_{1,2}}y_{3,2})(-y_{3,1}) = 0 \\ \frac{\partial{L}}{\partial{\color{red}{x_{1,2}}}}=2(\color{blue}{r_{1,2}} - \color{red}{x_{1,1}}y_{2,2}-\color{red}{x_{1,2}}y_{2,2})(-y_{2,2}) + 2(\color{blue}{r_{1,3}} - \color{red}{x_{1,1}}y_{3,1}-\color{red}{x_{1,2}}y_{3,2})(-y_{3,2}) = 0\)
寫成矩陣形式可得
利用我們上述的規則，很容易檢驗我們導出的結論。
總結來說，ALS的整個演算法過程只有兩步，涉及2個迴圈，如下圖所示：
演算法使用RMSE（root-mean-square error）評估誤差。
\(rmse = \sqrt{\frac{1}{|K|}\sum\limits_{(u,i) \in K}(r_{u,i}-{\bf{x}}_u^{\text T}{\bf{y}}_i)^2}\)
當RMSE值變化很小時或者到達最大迭代步驟時，滿足收斂條件，停止迭代。
“Talk is cheap. Show me the code.” 作為小練習，我們給出上述偽代碼的Python實現。
import numpy as np
from scipy.linalg import solve as linear_solve
# 評分矩陣 5 x 6
R = np.array([[4, 0, 2, 5, 0, 0], [3, 2, 1, 0, 0, 3], [0, 2, 0, 3, 0, 4], [0, 3, 3,5, 4, 0], [5, 0, 3, 4, 0, 0]])
m = 5 # 用戶數
n = 6 # 物品數
k = 3 # 隱向量的維度
_lambda = 0.01 # 正則化繫數
# 隨機初始化用戶矩陣, 物品矩陣
X = np.random.rand(m, k)
Y = np.random.rand(n, k)
# 每個用戶打分的物品集合
X_idx_dict = {1: [1, 3, 4], 2: [1, 2, 3, 6], 3: [2, 4, 6], 4: [2, 3, 4, 5], 5: [1, 3, 4]}
# 每個物品被打分的用戶集合
Y_idx_dict = {1: [1, 2, 5], 2: [2, 3, 4], 3: [1, 2, 4, 5], 4: [1, 3, 4, 5], 5: [4], 6: [2, 3]}
# 迭代10次
for iter in range(10):
for u in range(1, m+1):
Iu = np.array(X_idx_dict[u])
YIu = Y[Iu-1]
YIuT = YIu.T
RuIu = R[u-1, Iu-1]
xu = linear_solve(YIuT.dot(YIu) + _lambda * np.eye(k), YIuT.dot(RuIu))
X[u-1] = xu
for i in range(1, n+1):
Ui = np.array(Y_idx_dict[i])
XUi = X[Ui-1]
XUiT = XUi.T
RiUi = R.T[i-1, Ui-1]
yi = linear_solve(XUiT.dot(XUi) + _lambda * np.eye(k), XUiT.dot(RiUi))
Y[i-1] = yi
最終，我們列印用戶矩陣，物品矩陣，預測的評分矩陣如下，可以看到預測的評分矩陣非常逼近原始評分矩陣。
# X
array([[1.30678487, 2.03300876, 3.70447639],
[4.96150381, 1.03500693, 1.62261161],
[6.37691007, 2.4290095 , 1.03465981],
[0.41680155, 3.31805612, 3.24755801],
[1.26803845, 3.57580564, 2.08450113]])
# Y
array([[ 0.24891282, 1.07434519, 0.40258993],
[ 0.12832662, 0.17923216, 0.72376732],
[-0.00149517, 0.77412863, 0.12191856],
[ 0.12398438, 0.46163336, 1.05188691],
[ 0.07668894, 0.61050204, 0.59753081],
[ 0.53437855, 0.20862131, 0.08185176]])
# X.dot(Y.T) 預測評分
array([[4.00081359, 3.2132548 , 2.02350084, 4.9972158 , 3.55491072, 1.42566466],
[3.00018371, 1.99659282, 0.99163666, 2.79974661, 1.98192672, 3.00005934],
[4.61343295, 2.00253692, 1.99697545, 3.00029418, 2.59019481, 3.99911584],
[4.97591903, 2.99866546, 2.96391664, 4.99946603, 3.99816006, 1.18076534],
[4.99647978, 2.31231627, 3.02037696, 4.0005876 , 3.5258348 , 1.59422188]])
# 原始評分矩陣
array([[4, 0, 2, 5, 0, 0],
[3, 2, 1, 0, 0, 3],
[0, 2, 0, 3, 0, 4],
[0, 3, 3, 5, 4, 0],
[5, 0, 3, 4, 0, 0]])
三、Spark ALS應用
Spark的內部實現並不是我們上面所列的演算法，但是核心原理是完全一樣的，Spark實現的是上述偽代碼的分散式版本，具體演算法參考Large-scale Parallel Collaborative Filtering for the Netflix Prize。其次，查閱Spark的官方文檔，我們也註意到，Spark使用的懲罰函數與我們上文的有細微的差別。
\(\lambda \left( \sum\limits_{u}{\color{red}{n_u}\bf{x}}_u^{\text T}{\bf{x}}_u+\sum\limits_{i}{\color{red}{n_i}\bf{y}}_i^{\text T}{\bf{y}}_i\right)\)
其中\(n_u, n_i\)別表示用戶\(u\)打分的物品數量和物品\(i\)被打分的用戶數量。即
\(\begin{cases} n_u = |I_u| \\ n_i = |U_i| \\ \end{cases}\)
本小節通過兩個案例來瞭解Spark ALS的具體使用，以及在面對互聯網實際工程場景下的應用。
3.1 Demo案例
以第一節給出的數據為例，將三元組(User, Item, Rating)組織為als-demo-data.csv，該demo數據集涉及5個用戶和6個物品。
userId,itemId,rating
1,1,4
1,3,2
1,4,5
2,1,3
2,2,2
2,3,1
2,6,3
3,2,2
3,4,3
3,6,4
4,2,3
4,3,3
4,4,5
4,5,4
5,1,5
5,3,3
5,4,4
使用Spark的ALS類使用非常簡單，只需將三元組(User, Item, Rating)數據輸入模型進行訓練。
import org.apache.spark.sql.SparkSession
import org.apache.spark.ml.recommendation.ALS
val spark = SparkSession.builder().appName("als-demo").master("local[*]").getOrCreate()
val rating = spark.read
.options(Map("inferSchema" -> "true", "delimiter" -> ",", "header" -> "true"))
.csv("./data/als-demo-data.csv")
// 展示前5條評分記錄
rating.show(5)
val als = new ALS()
.setMaxIter(10) // 迭代次數，用於最小二乘交替迭代的次數
.setRank(3) // 隱向量的維度
.setRegParam(0.01) // 懲罰繫數
.setUserCol("userId") // user_id
.setItemCol("itemId") // item_id
.setRatingCol("rating") // 評分列
val model = als.fit(rating) // 訓練模型
// 列印用戶向量和物品向量
model.userFactors.show(truncate = false)
model.itemFactors.show(truncate = false)
// 給所有用戶推薦2個物品
model.recommendForAllUsers(2).show()
上述代碼在控制台輸出結果如下：
+------+------+------+
|userId|itemId|rating|
+------+------+------+
| 1| 1| 4|
| 1| 3| 2|
| 1| 4| 5|
| 2| 1| 3|
| 2| 2| 2|
+------+------+------+
only showing top 5 rows
+---+------------------------------------+
|id |features |
+---+------------------------------------+
|1 |[-0.17339179, 1.3144133, 0.04453602]|
|2 |[-0.3189066, 1.0291641, 0.12700711] |
|3 |[-0.6425665, 1.2283803, 0.26179287] |
|4 |[0.5160747, 0.81320006, -0.57953185]|
|5 |[0.645193, 0.26639006, 0.68648624] |
+---+------------------------------------+
+---+-----------------------------------+
|id |features |
+---+-----------------------------------+
|1 |[2.609607, 3.2668495, 3.554771] |
|2 |[0.85432494, 2.3137972, -1.1198239]|
|3 |[3.280517, 1.9563107, 0.51483333] |
|4 |[3.7446978, 4.259611, 0.6640027] |
|5 |[1.6036265, 2.5602736, -1.8897828] |
|6 |[-1.2651576, 2.4723763, 0.51556784]|
+---+-----------------------------------+
+------+--------------------------------+
|userId|recommendations |
+------+--------------------------------+
|1 |[[4, 4.9791617], [1, 3.9998217]]| // 對應物品的序號和預測評分
|2 |[[4, 3.273963], [6, 3.0134287]] |
|3 |[[6, 3.9849386], [1, 3.2667015]]|
|4 |[[4, 5.011649], [5, 4.004795]] |
|5 |[[1, 4.994258], [4, 4.0065994]] |
+------+--------------------------------+
我們使用numpy來驗證Spark的結果，並且用Excel可視化評分矩陣。
import numpy as np
X = np.array([[-0.17339179, 1.3144133, 0.04453602],
[-0.3189066, 1.0291641, 0.12700711],
[-0.6425665, 1.2283803, 0.26179287],
[0.5160747, 0.81320006, -0.57953185],
[0.645193, 0.26639006, 0.68648624]])
Y = np.array([[2.609607, 3.2668495, 3.554771],
[0.85432494, 2.3137972, -1.1198239],
[3.280517, 1.9563107, 0.51483333],
[3.7446978, 4.259611, 0.6640027],
[1.6036265, 2.5602736, -1.8897828],
[-1.2651576, 2.4723763, 0.51556784]])
R_predict = X.dot(Y.T)
R_predict
輸出預測的評分矩陣如下：
array([[3.99982136, 2.84328038, 2.02551472, 4.97916153, 3.0030386, 3.49205357],
[2.98138452, 1.96660155, 1.03257371, 3.27396294, 1.88351875, 3.01342882],
[3.26670123, 2.0001004 , 0.42992289, 3.00003605, 1.61982132, 3.98493822],
[1.94325135, 2.97144913, 2.98550149, 5.011649 , 4.00479503, 1.05883274],
[4.99425778, 0.39883335, 2.99113433, 4.00659955, 0.41937014, 0.19627587]])
從Excel可視化的評分矩陣可以觀察到預測的評分矩陣非常逼近原始的評分矩陣，以user-3為例，Spark推薦的物品是item-6和item-1, [6, 3.9849386, 1, 3.2667015]，這和Excel展示的預測評分矩陣完全一致。
從Spark函數recommendForAllUsers()給出的結果來看，Spark內部並沒有去除用戶已經購買的物品。
3.2 工程應用
在互聯網場景，用戶數 \(m\)(千萬~億級別) 和物品數 \(n\)(10萬~100萬級別) 規模很大，App的埋點數據一般會保存在HDFS中，以互聯網的長視頻場景為例，用戶的埋點信息最終聚合為用戶行為表 t_user_behavior。
行為表包含用戶的imei，物品的content-id，但是沒有直接的用戶評分，實踐中我們的解決方案是利用用戶的其他行為進行加權得出用戶對物品的評分。即
rating = w1 * play_time (播放時長) + w2 * finsh_play_cnt (完成的播放次數) + w3 * praise_cnt (點贊次數) + w4 * share_cnt (分享次數) + 其他適合於你業務邏輯的指標
其中, wi為每個指標對應的權重。
如下的代碼塊演示了工程實踐中對大規模用戶和商品場景進行推薦的流程。
import org.apache.spark.ml.feature.{IndexToString, StringIndexer}
// 從hive載入數據，並利用權重公式計算用戶對物品的評分
val rating_df = spark.sql("select imei, content_id, 權重公式計算評分 as rating from t_user_behavior group by imei, content_id")
// 將imei和content_id轉換為序號，Spark ALS入參要求userId, itemId為整數
// 使用org.apache.spark.ml.feature.StringIndexer
val imeiIndexer = new StringIndexer().setInputCol("imei").setOutputCol("userId").fit(rating_df)
val contentIndexer = new StringIndexer().setInputCol("content_id").setOutputCol("itemId").fit(rating_df)
val ratings = contentIndexer.transform(imeiIndexer.transform(rating_df))
// 其他code，類似於上述demo
val model = als.fit(ratings)
// 給每個用戶推薦100個物品
val _userRecs = model.recommendForAllUsers(100)
// 將userId, itemId轉換為原來的imei和content_id
val imeiConverter = new IndexToString().setInputCol("userId").setOutputCol("imei").setLabels(imeiIndexer.labels)
val contentConverter = new IndexToString().setInputCol("itemId").setOutputCol("content_id").setLabels(contentIndexer.labels)
val userRecs = imeiConverter.transform(_userRecs)
// 離線保存供線上調用
userRecs.foreachPartition {
// contentConverter 將itemId轉換為content_id
// 保存redis邏輯
}
值得註意的是，上述的工程場景還有一種解決方案，即隱式反饋。用戶給商品評分很單一，在實際的場景中，用戶未必會給物品打分，但是大量的用戶行為，同樣能夠間接反映用戶的喜好，比如用戶的購買記錄、搜索關鍵字，加入購物車，單曲迴圈播放同一首歌。我們將這些間接用戶行為稱之為隱式反饋，以區別於評分對應的顯式反饋。胡一凡等人在論文Collaborative filtering for implicit feedback datasets中針對隱式反饋場景提出了ALS-WR模型 (ALS with Weighted-λ-Regularization)，並且Spark官方也實現了該模型，我們將在以後的文章中介紹該模型。
四、總結
本文從推薦的場景出發，引出了協同過濾這一經典的推薦演算法，並且由此講解了被Spark唯一實現和維護的矩陣分解演算法，詳細推導了顯示反饋下矩陣分解的理論原理，並且給出了Python版本的單機實現，能夠讓讀者更好的理解矩陣這一演算法，最後我們以demo和工程實踐兩個實例講解了Spark ALS的使用，能夠讓沒有接觸過推薦演算法的同學有個直觀的理解，用理論與實踐的形式明白矩陣分解這一推薦演算法背後的原理。
參考文獻：
- 王喆, 深度學習推薦系統
- Hu, Yifan, Yehuda Koren, and Chris Volinsky. "Collaborative filtering for implicit feedback datasets." 2008 Eighth IEEE International Conference on Data Mining. IEEE, 2008.
- Zhou, Yunhong, et al. "Large-scale parallel collaborative filtering for the Netflix prize." International conference on algorithmic applications in management. Springer, Berlin, Heidelberg, 2008.
