Matrix Factorization推荐系统矩阵分解详解

实验和代码实现

https://github.com/Myolive-Lin/RecSys--deep-learning-recommendation-system

引言

推荐系统已经成为现代互联网的重要组成部分，特别是在电商、社交媒体、视频平台等领域，帮助用户发现他们可能感兴趣的商品或内容。推荐系统的核心目标是根据用户的历史行为或偏好，预测他们对未见过物品的评分或兴趣程度。矩阵分解作为推荐系统中的一种经典技术，在大规模数据上有着高效的表现。本文将详细介绍矩阵分解的基本原理、应用场景及其在推荐系统中的实现方法。（注其他方法可以在往期中进行查看）

---

矩阵分解的基本原理

矩阵分解（Matrix Factorization）是将一个大的稀疏矩阵（例如用户-物品评分矩阵）分解成两个低秩矩阵的过程。这些低秩矩阵表示了用户和物品之间的隐含特征。通过矩阵分解，我们可以更高效地进行评分预测。

假设我们有一个用户-物品评分矩阵 $R$，其中行表示用户，列表示物品，矩阵中的每个元素 $r_{ui}$ 表示用户 $u$ 对物品 $i$ 的评分。如果评分矩阵是稀疏的，即很多评分为空，我们的目标是填补这些空缺的评分，预测用户对未评过的物品的兴趣。

矩阵分解的核心思想是将评分矩阵 $R$ 分解成两个低秩矩阵：

• $P$：用户隐特征矩阵，形状为 $m\times k$，其中 $m$ 是用户数量，$k$ 是隐特征的维度。
• $Q$：物品隐特征矩阵，形状为 $k\times n$，其中 $n$ 是物品数量，$k$ 是隐特征的维度。

我们的目标是将评分矩阵 $R$ 近似为： $$R\approx P\times Q$$

其中 $P\times Q$ 是对评分矩阵的预测。每一行 $P_u$ 表示用户 $u$ 的隐特征向量，每一列 $Q_i$ 表示物品 $i$ 的隐特征向量。

评分矩阵示例

假设我们有如下的用户-物品评分矩阵 $R$，其中用户对物品的评分是一个稀疏矩阵，一些评分缺失：

$$R=\left[\begin{array}{cccc}5& 3& 0& 1\\ 4& 0& 0& 1\\ 1& 1& 0& 5\\ 1& 0& 0& 4\\ 0& 1& 5& 4\end{array}\right]$$

其中，$r_{ui}$ 表示用户 $u$ 对物品 $i$ 的评分。例如，用户 1 对物品 1 的评分为 5，用户 2 对物品 2 的评分为 0（即缺失）。

用户矩阵 $P$ 和物品矩阵 $Q$

在矩阵分解中，我们将评分矩阵 $R$ 分解为两个矩阵：用户矩阵 $P$ 和物品矩阵 $Q$。用户矩阵 $P$ 表示用户对潜在特征的偏好，而物品矩阵 $Q$ 表示物品在这些潜在特征上的表现。

假设我们将矩阵分解为 $k=2$ 个潜在特征，用户矩阵 $P$ 和物品矩阵 $Q$ 如下：

用户矩阵 $P$:

$$P=\left[\begin{array}{cc}2.31045989& 0.74082556\\ 1.8393161& 0.67679475\\ 0.29166742& 2.17961771\\ 0.32361523& 1.75546914\\ 0.4216019& 1.77386313\end{array}\right]$$

物品矩阵 $Q$:

$$Q=\left[\begin{array}{cccc}2.10425048& 1.20169726& 1.97073538& -0.31480952\\ 0.1800648& 0.29058456& 2.34189001& 2.33241682\end{array}\right]$$

预测评分矩阵 $\hat{R}$

通过矩阵分解，预测评分矩阵 $\hat{R}$ 可以通过 $P$ 和 $Q$ 的乘积得到：

$$\hat{R}=PQ$$

这里，$P$ 是用户矩阵，$Q$ 是物品矩阵。我们通过计算 $P$ 和 $Q$ 的矩阵乘积来得到预测的评分矩阵：

$$\hat{R}=\left[\begin{array}{cccc}4.995183& 2.991746& 6.288237& 1.000559\\ 3.992249& 2.406967& 5.209784& 0.999533\\ 1.006214& 0.983859& 5.679224& 4.991957\\ 0.997066& 0.899000& 4.748876& 3.992609\\ 1.206566& 1.022095& 4.985058& 4.004664\end{array}\right]$$

这就是根据用户矩阵和物品矩阵的乘积得到的预测评分矩阵,通过这个矩阵，就可以预测用户对未评价物品的评分

矩阵分解中的损失函数

矩阵分解的目标是最小化预测评分与实际评分之间的误差。为了衡量预测评分与实际评分之间的差异，我们使用均方误差（Mean Squared Error，MSE）作为损失函数：

$$\begin{array}{r}L=\sum _{(u,i)\in \text{non-zero entries}}(r_{ui}-q{{}_i}^T{p}_u{)}^2+\lambda (||q_i|{|}^2+||p_u|{|}^2)\end{array}$$

其中：

• $r_{ui}$ 是用户 $u$ 对物品 $i$ 的实际评分。
• $q_i^T p_u $ 是用户 $u$ 对物品 $i$ 的预测评分。其中$p_u$是用户u在用户矩阵P中的对应行向量，$q_i$是物品i在物品矩阵Q中的对应列向量。
• $\lambda$ 是正则化参数，防止模型过拟合。

正则化项 $\lambda \left(||p_u|{|}^2+||q_i|{|}^2\right)$ 是用来惩罚较大的特征值，使得模型能够避免过拟合。通过调整 $\lambda$ 的值，可以控制正则化的强度。

梯度下降法优化矩阵分解

矩阵分解的目标是通过优化损失函数来得到 $P$ 和 $Q$。常用的优化方法是梯度下降法（Gradient Descent）。梯度下降法通过计算损失函数的梯度并沿着梯度方向更新参数。

损失函数关于 $P_u$ 和 $Q_i$ 的梯度分别为：

$$\begin{array}{r}\frac{\partial L}{\partial P_u}=-2\cdot (r_{ui}-q_i^T{p}_u)\cdot q_i+2\lambda P_u\end{array}$$

$$\begin{array}{r}\frac{\partial L}{\partial q_i}=-2\cdot (r_{ui}-q_i^T{p}_u)\cdot p_u+2\lambda q_i\end{array}$$

然后通过更新规则：

$$\begin{array}{r}{P}_u:=P_u-\eta \cdot \frac{\partial L}{\partial P_u}\\ Q_i:=Q_i-\eta \cdot \frac{\partial L}{\partial Q_i}\end{array}$$

其中，$\eta$ 是学习率，控制更新步长。

通过迭代优化这些隐特征矩阵 $P$ 和 $Q$，我们能够逐渐减少预测误差，最终得到较为准确的评分预测。

Code

#使用for循环方法
def matrix_factorization_original(Matrix, k, learning_rate=0.01, reg_param=0.1, iterations=1000, seed=None):
    """
    Matrix: 评分矩阵
    k: 隐特征维度
    learning_rate: 学习率
    iterations: 迭代次数
    seed: 随机数种子
    """
    
    if seed is not None:
        np.random.seed(seed)

    rows, cols = Matrix.shape

    # Initialize user matrix and item matrix
    user_matrix = np.random.rand(rows, k).astype(np.float32)
    item_matrix = np.random.rand(k, cols).astype(np.float32)
    
    pbar = tqdm(range(iterations), desc="Training Progress", ncols=100, unit="iter")

    for iter in pbar:
        total_loss = 0  # 用来计算每次的总误差
        reg_loss = 0  # 用来计算每次的正则化损失

        for i in range(rows):
            for j in range(cols):
                if Matrix[i, j] > 0:
                    error = Matrix[i, j] - np.dot(item_matrix[:, j], user_matrix[i, :])

                    # 计算物品矩阵的梯度
                    q_i = item_matrix[:, j].copy()
                    q_grad = -2 * error * user_matrix[i, :] + 2 * reg_param * item_matrix[:, j]
                    item_matrix[:, j] -= learning_rate * q_grad

                    # 计算用户矩阵的梯度
                    p_grad = -2 * error * q_i + 2 * reg_param * user_matrix[i, :]
                    user_matrix[i, :] -= learning_rate * p_grad

                    # 计算总损失（错误项）
                    total_loss += (error ** 2)

                    # 计算正则化损失
                    reg_loss += reg_param * (np.sum(user_matrix[i, :]**2) + np.sum(item_matrix[:, j]**2))
        total_loss += reg_loss      
        pbar.set_postfix(error=total_loss)
         
    return user_matrix, item_matrix

利用矩阵的计算方法

上面的式子可以推广到利用矩阵进行计算

$$\begin{array}{r}E=(R-UV)\end{array}$$

更新用户矩阵 (U)：
$$\begin{array}{r}U=U+\eta \cdot \left(2\cdot E\cdot V^T-2\lambda U\right)\end{array}$$

更新物品矩阵 (V)：
$$\begin{array}{r}V=V+\eta \cdot \left(2\cdot E^T\cdot U-2\lambda V\right)\end{array}$$

其中

• $(U\in {\mathbb{R}}^{m\times k})$: 用户矩阵，每一行表示一个用户的特征向量 $(p_u)$。
• $(V\in {\mathbb{R}}^{k\times n}$: 物品矩阵，每一列表示一个物品的特征向量 $(q_i)$。

Code

def matrix_factorization(matrix,k, learning_rate = 0.01, reg_param=0.1,iterations = 1000, seed = 42):
    """
    matrix: 评分矩阵
    k: 隐特征维度
    learning_rate: 学习率
    iterations: 迭代次数
    seed: 随机数种子
    """
    
    if seed:
        np.random.seed(seed)

    rows,cols = matrix.shape

    # Initialize user matrix and item matrix
    user_matrix = np.random.rand(rows, k)
    item_matrix = np.random.rand(k,cols)

    #Only use non-zero elements in the matrix to optimize
    mask = (matrix > 0).astype(np.float32) 
    
    pbar = tqdm(range(iterations),desc= "Training",ncols = 100, unit = 'iter')
    for i in pbar:
        Predicted = np.dot(user_matrix,item_matrix)
        #这里值计算了非零元素的损失
        error = (matrix - Predicted) * mask

        #Gradient calculation and update

        user_grad = -2 *np.dot(error,item_matrix.T) + 2 * reg_param * user_matrix
        user_matrix -= learning_rate * (user_grad)

        item_grad = -2 * np.dot(user_matrix.T, error) + 2 * reg_param * item_matrix
        item_matrix -= learning_rate * item_grad

        #Calculate the error
        mse = np.sum(error ** 2) 

        # 计算正则化项
        reg_term = reg_param * (np.sum(user_matrix ** 2) + np.sum(item_matrix ** 2))

        # 计算带正则化的代价函数
        cost = mse + reg_term
        pbar.set_postfix(MSE=mse, Cost=cost)  # 显示MSE和代价函数
        
    return user_matrix,item_matrix

矩阵分解中的偏差项

在实际应用中，除了用户和物品的隐特征矩阵外，我们还需要考虑用户和物品的偏差。不同用户的评分习惯不同，可能导致评分系统的偏差。例如，某些用户可能普遍给出较高的评分，而另一些用户可能给出较低的评分。此外，不同物品的评分标准也可能存在差异。 为了消除这些偏差，通常在矩阵分解模型中加入以下偏差项：

• 全局偏差（$u$）：所有评分的平均值。
• 用户偏差（$b_u$）：用户 $u$ 对所有物品的评分平均值。
• 物品偏差（$b_i$）：物品 $i$ 的平均评分。 因此，评分的预测公式可以更新为：

$$\begin{array}{r}{r}_{ui}=u+b_u+b_i+q_i^T{p}_u\end{array}$$

---

损失函数
$$
L = \sum_{(u,i) \in \text{non-zero entries}} \left[ (r_{ui} - (u + b_u + b_i + p_u^T q_i))^2 \right] \

• \lambda \left(||p_u||^2 + ||q_i||^2 + b_u^2 + b_i^2 \right)
  $$

对用户矩阵 $p_u$求导：
$$\frac{\partial L}{\partial p_u}=-2\sum _{i\in \text{non-zero entries}}(r_{ui}-(u+b_u+b_i+p_u^T{q}_i))q_i+2\lambda p_u$$
对物品矩阵 $q_i$ 求偏导：
$$\frac{\partial L}{\partial q_i}=-2\sum _{u\in \text{non-zero entries}}(r_{ui}-(u+b_u+b_i+p_u^T{q}_i))p_u+2\lambda q_i$$
对用户偏差 $b_u$ 求偏导：
$$\frac{\partial L}{\partial b_u}=-2\sum _{i\in \text{non-zero entries}}(r_{ui}-(u+b_u+b_i+p_u^T{q}_i))+2\lambda b_u$$
对物品偏差 $b_i$求偏导：
$$\frac{\partial L}{\partial b_i}=-2\sum _{u\in \text{non-zero entries}}(r_{ui}-(u+b_u+b_i+p_u^T{q}_i))+2\lambda b_i$$

• 是用户 u 对物品 i 的评分。
• $\mu$ 是全局偏差，$b_u$ 是用户偏差，$b_i$ 是物品偏差。
• $p_u$ 和 $q_i$ 分别是用户 u 和物品 i 的隐特征向量。
• λ 是正则化系数，用于控制过拟合。

Code

import numpy as np
from tqdm import tqdm

def matrix_factorization_by_bias(matrix, k, learning_rate=0.01, reg_param=0.1, iterations=1000, seed=42):
    """
    matrix: 评分矩阵
    k: 隐特征维度
    learning_rate: 学习率
    reg_param: 正则化参数
    iterations: 迭代次数
    seed: 随机数种子
    """
    np.random.seed(seed)
    rows, cols = matrix.shape
    
    # 初始化全局偏差
    u = np.mean(matrix[matrix > 0])  # 所有评分的全局平均分
    
    # 物品和用户的偏差初始化
    b_i = np.mean(matrix, axis=0)  # 物品偏差系数
    b_u = np.mean(matrix, axis=1)  # 用户偏差系数
    
    # 随机初始化隐因子矩阵
    user_matrix = np.random.rand(rows, k)
    item_matrix = np.random.rand(k, cols)

    # 迭代训练
    pbar = tqdm(range(iterations), desc="Training", ncols=100, unit='iter')
    
    for iter in pbar:
        total_loss = 0

        # 遍历每个评分
        for i in range(rows):
            for j in range(cols):
                if matrix[i, j] > 0:  # 只对非零评分进行优化
                    Prediction = u + b_i[j] + b_u[i] + np.dot(user_matrix[i, :], item_matrix[:, j])
                    error = matrix[i, j] - Prediction

                    # 计算梯度
                    user_grad = -2 * error * item_matrix[:, j] + 2 * reg_param * user_matrix[i, :]
                    item_grad = -2 * error * user_matrix[i, :] + 2 * reg_param * item_matrix[:, j]

                    # 更新隐因子矩阵
                    user_matrix[i, :] -= learning_rate * user_grad
                    item_matrix[:, j] -= learning_rate * item_grad

                    # 更新物品和用户的偏差
                    b_i[j] -= learning_rate * (-2 * error + 2 * reg_param * b_i[j])
                    b_u[i] -= learning_rate * (-2 * error + 2 * reg_param * b_u[i])

                    # 累加误差
                    total_loss += error ** 2

        # 更新进度条
        pbar.set_postfix({"Loss": total_loss})

    return u, b_i, b_u, user_matrix, item_matrix

总结

矩阵分解通过将评分矩阵 $R$ 分解为两个低秩矩阵 $P$ 和 $Q$，利用隐因子表示用户和物品之间的潜在特征，从而帮助我们预测用户对未评分物品的兴趣。这一方法能够有效地填补评分矩阵中的缺失值，从而实现个性化推荐。

Reference

王喆 《深度学习推荐系统》