【统计学习方法死磕系列Ⅱ】奇异值分解

一、引言：矩阵分解的”金牌标准”

如果说线性代数中有一个分解是”万能”的，那一定是奇异值分解（Singular Value Decomposition, SVD）。与特征值分解（EVD）不同，SVD 对任意矩阵（不要求方阵、不要求对称）都适用，这使得它成为连接线性代数与数据科学的桥梁。

SVD 的应用极其广泛：PCA 降维、图像压缩、推荐系统、潜在语义分析（LSA）、伪逆计算、矩阵补全……几乎任何一个涉及矩阵的数据分析任务，背后都有 SVD 的影子。

本文将深入 SVD 的数学定义、几何意义、计算方法、核心定理和应用场景，帮你彻底掌握这一基础而强大的工具。

二、SVD 的定义与几何意义

2.1 定义

对于任意矩阵 $A \in \mathbb{R}^{m \times n}$，其奇异值分解为：

$$ A = U \Sigma V^T $$

其中：

$U \in \mathbb{R}^{m \times m}$ 是正交矩阵（$U^T U = I$），其列向量称为左奇异向量
$V \in \mathbb{R}^{n \times n}$ 是正交矩阵（$V^T V = I$），其列向量称为右奇异向量
$\Sigma \in \mathbb{R}^{m \times n}$ 是对角矩阵（严格说是矩形对角矩阵）：

$$ \Sigma = \begin{bmatrix} \sigma_1 & 0 & \cdots & 0 & \cdots & 0 \\ 0 & \sigma_2 & \cdots & 0 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots & \cdots & 0 \\ 0 & 0 & \cdots & \sigma_r & \cdots & 0 \\ \vdots & \vdots & \cdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & 0 & \cdots & 0 \end{bmatrix}_{m \times n} $$

其中 $\sigma_1 \geq \sigma_2 \geq \dots \geq \sigma_r > 0$，$r = \text{rank}(A)$。

2.2 几何理解

SVD 有一个非常优美的几何解释：任何线性变换 $A: \mathbb{R}^n \to \mathbb{R}^m$ 可以分解为三个基本变换的复合：

旋转/反射（$V^T$）：将输入空间的标准正交基旋转到右奇异向量方向
缩放（$\Sigma$）：沿各坐标轴缩放，缩放因子为奇异值 $\sigma_i$
旋转/反射（$U$）：将结果旋转到左奇异向量方向

换句话说：矩阵 $A$ 将 $\mathbb{R}^n$ 中的单位球映射为 $\mathbb{R}^m$ 中的一个超椭球。超椭球的各半轴长度就是奇异值 $\sigma_i$，各半轴方向就是左奇异向量 $u_i$。

2.3 紧凑形式

SVD 的紧凑形式（compact SVD）保留了秩信息：

$$ A = U_r \Sigma_r V_r^T = \sum_{i=1}^{r} \sigma_i u_i v_i^T $$

其中 $U_r$ 取 $U$ 的前 $r$ 列，$V_r$ 取 $V$ 的前 $r$ 列，$\Sigma_r = \text{diag}(\sigma_1, \dots, \sigma_r)$。

这个秩-1 分解的和形式揭示了 SVD 的核心洞察：任何矩阵都是若干个秩-1 矩阵的加权和。每个 $\sigma_i u_i v_i^T$ 是一个秩-1 矩阵，$\sigma_i$ 是其”重要性”权重。

三、SVD 与特征值分解（EVD）的关系

3.1 从 EVD 到 SVD

对于对称半正定矩阵，SVD 与 EVD 完全一致。实际上，SVD 可以通过 $A^T A$ 和 $AA^T$ 的特征分解导出：

$A^T A = V \Sigma^T U^T U \Sigma V^T = V \Sigma^2 V^T$

即 $V$ 的列是 $A^T A$ 的特征向量，$\sigma_i^2$ 是 $A^T A$ 的特征值。
$AA^T = U \Sigma V^T V \Sigma^T U^T = U \Sigma \Sigma^T U^T$

即 $U$ 的列是 $AA^T$ 的特征向量，$\sigma_i^2$ 也是 $AA^T$ 的非零特征值。

因此，奇异值是 $A^T A$（或 $AA^T$）特征值的算术平方根：

$$ \sigma_i = \sqrt{\lambda_i(A^T A)} $$

3.2 EVD 与 SVD 对比

性质	EVD	SVD
适用矩阵	方阵	任意矩阵
正交性要求	对称矩阵才有正交特征向量	U 和 V 总是正交的
数值稳定性	对非正规矩阵不稳定	总是数值稳定
特征值/奇异值	可为复数（非对称矩阵）	总是非负实数
与 PCA 关系	PCA 通过 EVD 计算	PCA 通过 SVD 计算更稳定

3.3 构造 SVD 示例

以下是一个手工计算 SVD 的数值例子：

设 $A = \begin{bmatrix} 1 & 0 \ 0 & 1 \ 1 & 1 \end{bmatrix}$。

步骤 1：计算 $A^T A$

$$ A^T A = \begin{bmatrix} 1 & 0 & 1 \\ 0 & 1 & 1 \end{bmatrix} \begin{bmatrix} 1 & 0 \\ 0 & 1 \\ 1 & 1 \end{bmatrix} = \begin{bmatrix} 2 & 1 \\ 1 & 2 \end{bmatrix} $$

步骤 2：对 $A^T A$ 做特征分解

特征方程：$(2-\lambda)^2 - 1 = 0 \implies \lambda^2 - 4\lambda + 3 = 0 \implies \lambda_1 = 3, \lambda_2 = 1$

奇异值：$\sigma_1 = \sqrt{3}, \sigma_2 = 1$

特征向量（归一化）：

$$ \lambda_1 = 3: \quad v_1 = \frac{1}{\sqrt{2}}\begin{bmatrix} 1 \\ 1 \end{bmatrix} $$
$$ \lambda_2 = 1: \quad v_2 = \frac{1}{\sqrt{2}}\begin{bmatrix} 1 \\ -1 \end{bmatrix} $$

步骤 3：计算左奇异向量

$$ u_1 = \frac{1}{\sigma_1} A v_1 = \frac{1}{\sqrt{3}} \begin{bmatrix} 1 & 0 \\ 0 & 1 \\ 1 & 1 \end{bmatrix} \frac{1}{\sqrt{2}}\begin{bmatrix} 1 \\ 1 \end{bmatrix} = \frac{1}{\sqrt{6}}\begin{bmatrix} 1 \\ 1 \\ 2 \end{bmatrix} $$

$$ u_2 = \frac{1}{\sigma_2} A v_2 = 1 \cdot \begin{bmatrix} 1 & 0 \\ 0 & 1 \\ 1 & 1 \end{bmatrix} \frac{1}{\sqrt{2}}\begin{bmatrix} 1 \\ -1 \end{bmatrix} = \frac{1}{\sqrt{2}}\begin{bmatrix} 1 \\ -1 \\ 0 \end{bmatrix} $$

$u_3$ 需要与前两个向量正交：$u_3 = \frac{1}{\sqrt{3}}\begin{bmatrix} 1 \ 1 \ -1 \end{bmatrix}$

最终：

$$ A = \begin{bmatrix} 1/\sqrt{6} & 1/\sqrt{2} & 1/\sqrt{3} \\ 1/\sqrt{6} & -1/\sqrt{2} & 1/\sqrt{3} \\ 2/\sqrt{6} & 0 & -1/\sqrt{3} \end{bmatrix} \begin{bmatrix} \sqrt{3} & 0 \\ 0 & 1 \\ 0 & 0 \end{bmatrix} \begin{bmatrix} 1/\sqrt{2} & 1/\sqrt{2} \\ 1/\sqrt{2} & -1/\sqrt{2} \end{bmatrix}^T $$

四、Eckart-Young 定理：SVD 的最优低秩近似

4.1 定理陈述

Eckart-Young 定理（1936）：设 $A \in \mathbb{R}^{m \times n}$ 的 SVD 为 $A = U \Sigma V^T$。对任意 $k \leq \text{rank}(A)$，令：

$$ A_k = \sum_{i=1}^{k} \sigma_i u_i v_i^T $$

即保留最大的 $k$ 个奇异值对应的秩-1 分量。则 $A_k$ 是 $A$ 在所有秩不超过 $k$ 的矩阵中的 Frobenius 范数最优近似：

$$ \|A - A_k\|_F = \min_{\text{rank}(B) \leq k} \|A - B\|_F = \sqrt{\sum_{i=k+1}^{r} \sigma_i^2} $$

同样，对谱范数也成立：

$$ \|A - A_k\|_2 = \min_{\text{rank}(B) \leq k} \|A - B\|_2 = \sigma_{k+1} $$

4.2 意义

Eckart-Young 定理是 SVD 最具实用价值的定理之一。它告诉我们：

SVD 提供了最优的低秩近似——如果你想用秩-k 矩阵近似 $A$，没有什么比截断 SVD 更好的了
奇异值衰减越快，低秩近似越精确——这解释了为什么 PCA（基于 SVD）能有效降维
误差可以直接从丢弃的奇异值计算——你不需要实际构造 $A_k$ 就能知道近似的质量

4.3 数值示例

用 SVD 压缩一张图像的示例：

import numpy as np
from PIL import Image

def svd_compress(image_path, k):
    # Load and convert to grayscale
    img = Image.open(image_path).convert('L')
    A = np.array(img, dtype=float)  # m x n matrix
    
    # SVD
    U, s, Vt = np.linalg.svd(A, full_matrices=False)
    
    # Keep top k singular values
    A_k = U[:, :k] @ np.diag(s[:k]) @ Vt[:k, :]
    
    # Clip to valid range
    A_k = np.clip(A_k, 0, 255).astype(np.uint8)
    
    # Compression ratio
    original_size = A.size
    compressed_size = k * (A.shape[0] + A.shape[1] + 1)
    ratio = compressed_size / original_size
    
    return A_k, ratio

# Example: compress a 500x500 image with k=50
A_k, ratio = svd_compress('lena.png', k=50)
# ratio ≈ 50*(500+500+1) / 250000 ≈ 0.20 (80% compression)

对于一张 500x500 的图像（250,000 个像素），用 k=50 的截断 SVD 只需要存储 $50 \cdot (500 + 500 + 1) = 50,050$ 个数，压缩率约 20%。

五、SVD 的数值计算

5.1 标准算法

完整的 SVD 计算分两步：先二对角化，再迭代。

Golub-Reinsch 算法：

1. Bidiagonalization: 
   Apply Householder reflections to reduce A to bidiagonal form B
   A = U1 B V1^T (via Golub-Kahan bidiagonalization)

2. Diagonalization of B:
   Use implicit QR iteration (Golub-Kahan SVD step) 
   to reduce B to diagonal form Σ
   B = U2 Σ V2^T

3. Combine:
   U = U1 U2, V = V1 V2, A = U Σ V^T

复杂度：$O(mn^2 + n^3)$（假设 $m \geq n$）。

5.2 随机化 SVD（Randomized SVD）

对于大规模矩阵，标准 SVD 的计算成本可能无法接受。随机化 SVD 是一种近似算法，效率极高：

Algorithm: Randomized SVD
----------------------------------------------
Input:  A ∈ R^{m×n}, target rank k, oversampling p
Output: Approximate truncated SVD: U_k, Σ_k, V_k

1. Generate random matrix Ω ∈ R^{n×(k+p)} 
   (e.g., standard Gaussian entries)

2. Compute sketch: Y = A Ω  (size m × (k+p))

3. Orthogonalize: Q = orth(Y) via QR decomposition

4. Form small matrix: B = Q^T A  (size (k+p) × n)

5. Compute SVD of small matrix: B = U_B Σ V^T
   (standard SVD on (k+p) × n matrix, cheap when k+p ≪ m,n)

6. Approximate left singular vectors: U_k = Q (U_B)_{:,1:k}

7. Approximate singular values: Σ_k = (Σ)_{1:k, 1:k}

8. Return: U_k, Σ_k, V_{:,1:k}

随机化 SVD 将复杂度从 $O(mn^2)$ 降为 $O(mn\log k)$，在 k 远小于 m,n 时优势显著。

def randomized_svd(A, k, p=5):
    """Randomized SVD approximation"""
    m, n = A.shape
    # Step 1: Random projection
    Omega = np.random.randn(n, k + p)
    # Step 2: Sketch
    Y = A @ Omega
    # Step 3: Orthogonalize
    Q, _ = np.linalg.qr(Y)
    # Step 4: Small matrix
    B = Q.T @ A
    # Step 5: SVD of small matrix
    U_B, s, Vt = np.linalg.svd(B, full_matrices=False)
    # Step 6: Recover left singular vectors
    U = Q @ U_B[:, :k]
    return U, s[:k], Vt[:k, :]

六、SVD 的核心应用

6.1 伪逆（Moore-Penrose Pseudoinverse）

SVD 为矩阵伪逆提供了最清晰的计算方式。对于 $A = U \Sigma V^T$，其伪逆为：

$$ A^+ = V \Sigma^+ U^T $$

其中 $\Sigma^+$ 是将 $\Sigma$ 的非零奇异值取倒数后转置：

$$ \Sigma^+_{ii} = \begin{cases} 1/\sigma_i & \text{if } \sigma_i > 0 \\ 0 & \text{otherwise} \end{cases} $$

最小二乘中的角色：对于线性方程 $Ax = b$，最小范数最小二乘解为：

$$ x_{\text{LS}} = A^+ b = V \Sigma^+ U^T b $$

当 $A$ 列满秩时，这等价于正规方程的解 $x = (A^T A)^{-1} A^T b$。但当 $A$ 不满秩时，伪逆给出的是唯一的最小范数解，而正规方程无唯一解。

6.2 应用对比表

应用	如何使用 SVD	关键参数
PCA 降维	$X = U\Sigma V^T$，主成分 = $U_k$	$k$（保留的主成分数）
图像压缩	$A_k = \sum_{i=1}^{k} \sigma_i u_i v_i^T$	$k$（压缩等级）
潜在语义分析（LSA）	词-文档矩阵的截断 SVD	$k$（主题数）
推荐系统	评分矩阵的截断 SVD	$k$（潜在因子数）
伪逆求解	$A^+ = V \Sigma^+ U^T$	$\tau$（奇异值截断阈值）
矩阵补全	迭代软阈值 SVD	$\lambda$（正则化参数）
总最小二乘	增广矩阵 $[A, b]$ 的最小奇异值对应右奇异向量	—
信号去噪	保留大奇异值分量，丢弃小奇异值分量	$k$（信号秩）

6.3 推荐系统：SVD 在 Netflix Prize 中的应用

在协同过滤推荐中，用户-物品评分矩阵 $R \in \mathbb{R}^{m \times n}$ 通常是稀疏的（大多数用户只评价过少数物品）。SVD 的应用思路是：

将评分矩阵分解为 $R \approx U_k \Sigma_k V_k^T$（或在 FunkSVD 中直接用 $P Q^T$）
$k$ 是潜在因子的数量
用户 $i$ 对物品 $j$ 的预测评分为：$\hat{r}_{ij} = (U_k \Sigma_k^{1/2})_i \cdot (\Sigma_k^{1/2} V_k^T)_j$

FunkSVD（Simon Funk 在 Netflix Prize 中实现的版本）通过 SGD 直接优化：

$$ \min_{P, Q} \sum_{(i,j) \in \Omega} (r_{ij} - p_i^T q_j)^2 + \lambda(\|p_i\|^2 + \|q_j\|^2) $$

其中 $\Omega$ 是已知评分集合，$p_i$ 是用户 $i$ 的潜在因子向量，$q_j$ 是物品 $j$ 的潜在因子向量。

6.4 条件数与数值稳定性

SVD 直接给出了矩阵的条件数（condition number）：

$$ \kappa(A) = \frac{\sigma_{\max}}{\sigma_{\min}} $$

条件数衡量了矩阵对数值误差的敏感程度。$\kappa(A)$ 越大，$Ax = b$ 的求解对 $b$ 的微小扰动越敏感。

L2 正则化与 SVD 的关系：岭回归（Ridge regression）等价于对奇异值做如下修正：

$$ \hat{\beta}_{\text{ridge}} = V (\Sigma^2 + \lambda I)^{-1} \Sigma U^T y $$

对比 OLS：$\hat{\beta}_{\text{OLS}} = V \Sigma^{-1} U^T y$

岭回归通过将每个奇异值 $\sigma_i$ 替换为 $\sigma_i / (\sigma_i^2 + \lambda)$ 来抑制小奇异值的放大效应，从而改善了条件数。

七、SVD 的 Python 实现

7.1 numpy/scipy 用法

import numpy as np
from scipy.linalg import svd

# Full SVD
U, s, Vt = np.linalg.svd(A, full_matrices=True)
# U: (m, m), s: (min(m,n),), Vt: (n, n)

# Economy SVD (recommended for m >> n or n >> m)
U, s, Vt = np.linalg.svd(A, full_matrices=False)
# U: (m, min(m,n)), s: (min(m,n),), Vt: (min(m,n), n)

# Truncated SVD
from sklearn.decomposition import TruncatedSVD
tsvd = TruncatedSVD(n_components=k)
A_k = tsvd.fit_transform(A)  # m x k (U_k @ Σ_k)

# Randomized SVD
from sklearn.utils.extmath import randomized_svd
U, s, Vt = randomized_svd(A, n_components=k, n_iter=5)

7.2 矩阵秩估计

def estimate_rank(A, energy_threshold=0.99):
    """
    Estimate effective rank via singular value energy.
    """
    _, s, _ = np.linalg.svd(A, full_matrices=False)
    cumulative_energy = np.cumsum(s) / np.sum(s)
    rank = np.searchsorted(cumulative_energy, energy_threshold) + 1
    return rank

6.5 SVD 与矩阵补全（Matrix Completion）

矩阵补全是推荐系统中广泛遇到的问题：用户-物品评分矩阵 $R \in \mathbb{R}^{m \times n}$ 中只观测到少量条目，如何预测缺失的评分？

将矩阵补全形式化为以下优化问题：

$$ \min_{X} \text{rank}(X), \quad \text{s.t. } X_{ij} = R_{ij}, \forall (i,j) \in \Omega $$

其中 $\Omega$ 是已观测到的位置集合。由于秩函数是 NP-hard 的，通常使用核范数（nuclear norm，即奇异值之和）作为秩的凸替代：

$$ \min_{X} \|X\|_*, \quad \text{s.t. } X_{ij} = R_{ij}, \forall (i,j) \in \Omega $$

软阈值 SVD（Soft-Impute） 是求解这一问题的经典迭代算法：

Algorithm: Soft-Impute (Mazumder et al., 2010)
----------------------------------------------
Input: Incomplete matrix R (with missing entries as 0), 
       regularization λ, tolerance ε
Output: Completed matrix X

Initialize X^(0) = R (zeros for missing entries)

Repeat:
    1. Compute SVD of current estimate:
       X^(t) = U Σ V^T
    
    2. Apply soft-thresholding to singular values:
       Σ_λ = diag(max(σ_i - λ, 0))
    
    3. Reconstruct:
       X_new = U Σ_λ V^T
    
    4. Fill in observed entries:
       For (i,j) in Ω: X^(t+1)_{ij} = R_{ij}
       For (i,j) not in Ω: X^(t+1)_{ij} = X_new_{ij}
    
    5. If ||X^(t+1) - X^(t)||_F < ε: break

软阈值 SVD 的直观理解：小的奇异值被”压缩”到零（去噪），大的奇异值被保留（信号），从而得到一个低秩的重构矩阵。参数 $\lambda$ 控制正则化强度——$\lambda$ 越大，秩越低。

这个框架与 Robust PCA（RPCA）有紧密联系。RPCA 将数据矩阵分解为低秩部分 $L$（结构信息）和稀疏部分 $S$（离群噪声）：

$$ \min_{L,S} \|L\|_* + \lambda\|S\|_1, \quad \text{s.t. } X = L + S $$

这同样可以通过迭代软阈值 SVD 来求解（交替方向乘子法，ADMM）。

Python 快速示例：

from sklearn.decomposition import TruncatedSVD
import numpy as np

def soft_impute(R, mask, k=20, lambda_=1.0, max_iter=100):
    """Soft-thresholded SVD for matrix completion"""
    X = R.copy()
    X[~mask] = 0
    
    for _ in range(max_iter):
        # SVD with truncation
        svd = TruncatedSVD(n_components=k)
        U = svd.fit_transform(X)
        s = svd.singular_values_
        Vt = svd.components_
        
        # Soft-threshold
        s_thresh = np.maximum(s - lambda_, 0)
        
        # Reconstruct low-rank estimate
        X_lowrank = U @ np.diag(s_thresh) @ Vt
        
        # Keep observed values, fill missing with low-rank estimate
        X_new = R.copy()
        X_new[~mask] = X_lowrank[~mask]
        
        if np.linalg.norm(X_new - X, 'fro') < 1e-4:
            break
        
        X = X_new
    
    return X

八、面试高频问题

Q1：SVD 和 EVD 有什么本质区别？为什么 SVD 更通用？

本质区别有三个层面：

适用范围：EVD 只适用于方阵（$n \times n$），SVD 适用于任意形状矩阵（$m \times n$）
正交性：EVD 保证特征向量正交的充分条件是矩阵为正规矩阵（$AA^T = A^T A$），对称矩阵是充分条件。SVD 的左右奇异向量 $U$ 和 $V$ 总是正交的，无需任何条件
特征值的实数性：EVD 的特征值在非对称矩阵中可能是复数，SVD 的奇异值永远是非负实数

SVD 更通用的数学根源是：对于任意 $A$，$A^T A$ 和 $AA^T$ 总是对称半正定的，因此总能正交对角化。这种”通过构造对称矩阵绕开限制”的技巧是 SVD 存在的关键。

Q2：Eckart-Young 定理有什么实际应用？

图像压缩：保留前 $k$ 个奇异值分量，将 $m \times n$ 矩阵的存储从 $mn$ 降为 $k(m+n+1)$
数据去噪：假设”信号”对应大奇异值分量，”噪声”对应小奇异值分量，去噪就是丢弃小奇异值
PCA：PCA 选择 $k$ 个主成分就是在应用 Eckart-Young 定理——这 $k$ 个主成分给出了数据在 Frobenius 范数下的最优秩-k 重构
潜在语义分析（LSA）：词-文档矩阵的秩-k 近似能揭示隐含的”主题”结构

Q3：SVD 的计算复杂度是多少？大规模数据怎么办？

标准的 Golub-Reinsch SVD 复杂度为 $O(mn^2 + n^3)$（假设 $m \geq n$）。

大规模数据的处理方法：

Economy/Thin SVD：只计算前 $\min(m,n)$ 个奇异向量，不计算 $U$ 的额外列
Truncated SVD：只计算前 $k$ 个奇异值和向量
随机化 SVD：通过随机投影将 $m \times n$ 问题转化为 $m \times (k+p)$ 问题，复杂度降至 $O(mn \log k)$
增量 SVD：逐批处理数据，适合流式场景

Q4：如何用 SVD 解释 PCA？

数据矩阵 $X \in \mathbb{R}^{n \times m}$（$n$ 特征，$m$ 样本）的 SVD 为 $X = U \Sigma V^T$。

协方差矩阵 $C = \frac{1}{m-1} XX^T = \frac{1}{m-1} U \Sigma^2 U^T$
所以 $U$ 的列就是主成分方向，$\frac{\sigma_i^2}{m-1}$ 就是第 $i$ 个主成分的方差
降维后的数据：$Z = U_k^T X = \Sigma_k V_k^T$

SVD 实现 PCA 的优势是避免构造协方差矩阵——直接对 $X$ 做 SVD，数值更稳定，且不需要计算 $n \times n$ 的矩阵。

Q5：岭回归与 SVD 之间的关系是什么？

岭回归的解为：

$$ \hat{\beta}_{\text{ridge}} = (X^T X + \lambda I)^{-1} X^T y $$

代入 $X = U \Sigma V^T$：

$$ \begin{aligned} \hat{\beta}_{\text{ridge}} &= (V \Sigma^2 V^T + \lambda V V^T)^{-1} V \Sigma U^T y \\ &= V(\Sigma^2 + \lambda I)^{-1} V^T V \Sigma U^T y \\ &= V(\Sigma^2 + \lambda I)^{-1} \Sigma U^T y \end{aligned} $$

第 $i$ 个分量乘以 $\frac{\sigma_i}{\sigma_i^2 + \lambda}$（OLS 是 $1/\sigma_i$）。岭回归对小的 $\sigma_i$（噪声方向）大幅衰减权重，从而防止过拟合。这称为收缩估计（shrinkage estimator）。

Q6：SVD 在协同过滤中有哪些具体变体？Bias-SVD 和 SVD++ 是什么？

标准的 SVD 协同过滤模型 $R \approx U_k \Sigma_k V_k^T$ 可以简化为 $R \approx P Q^T$（将 $\Sigma_k$ 吸收进 $P$ 和 $Q$）。

Bias-SVD 在基础模型上加入了用户偏好偏置 $b_u$、物品偏好偏置 $b_i$ 和全局均值 $\mu$：

$$ \hat{r}_{ui} = \mu + b_u + b_i + p_u^T q_i $$

这解决了”有些用户整体评分偏高、有些物品整体评分偏低”的问题。

**SVD++**（Koren, 2008）进一步融合了隐式反馈：

$$ \hat{r}_{ui} = \mu + b_u + b_i + q_i^T \left(p_u + |N(u)|^{-1/2} \sum_{j \in N(u)} y_j\right) $$

其中 $N(u)$ 是用户 $u$ 有过隐式反馈（浏览、点击等）的物品集合，$y_j$ 是物品 $j$ 的隐式因子向量。SVD++ 的核心见解是：用户不仅由”他们是谁”（$p_u$）定义，也由”他们接触过什么”（$\sum y_j$）定义。

SVD++ 在 Netflix Prize 中表现极为出色，是获胜方案的关键组成部分。

模型	优化目标	参数	特点
基础 SVD	$\sum (r_{ui} - p_u^T q_i)^2$	$p_u, q_i$	仅显式评分
Bias-SVD	$\sum (r_{ui} - \mu - b_u - b_i - p_u^T q_i)^2$	$p_u, q_i, b_u, b_i, \mu$	考虑评分偏差
SVD++	加隐式反馈项	再加上 $y_j$	显式+隐式反馈融合

九、总结

奇异值分解是线性代数中”没有之一”的最重要的矩阵分解。它以优雅的 $U \Sigma V^T$ 形式统一了特征分解和矩阵近似，并凭借 Eckart-Young 定理为低秩近似提供了严格的最优性保证。

从 PCA 到推荐系统，从图像压缩到自然语言处理，SVD 的身影无处不在。理解和熟练运用 SVD，是每个数据科学家和机器学习工程师的基本功。

在实际编码中，记住两个原则：

能用 full_matrices=False 就不要用完整 SVD（省内存）
当 $k \ll \min(m,n)$ 时使用随机化 SVD（省时间）

参考文献：

Golub, G. H., & Reinsch, C. (1970). Singular value decomposition and least squares solutions. Numerische Mathematik, 14(5), 403-420.
Eckart, C., & Young, G. (1936). The approximation of one matrix by another of lower rank. Psychometrika, 1(3), 211-218.
Halko, N., Martinsson, P. G., & Tropp, J. A. (2011). Finding structure with randomness: Probabilistic algorithms for constructing approximate matrix decompositions. SIAM Review, 53(2), 217-288.
李航. (2019). 《统计学习方法》（第2版）. 第15章.
Strang, G. (2016). Introduction to Linear Algebra (5th ed.). Chapter 7.