矩阵分解的含义与应用

一、LU分解

1. 基本概念

• 含义：将矩阵分解为下三角矩阵（L）和上三角矩阵（U）的乘积
• 公式：A = L × U
• 扩展：A = P × L × U（包含置换矩阵P）

2. 几何意义

• L矩阵：记录行变换操作（消元过程）
• U矩阵：消元后的上三角形式
• P矩阵：记录行交换操作（用于主元选取）

python

import numpy as np
from scipy.linalg import lu

A = np.array([[2, 5, 8], [5, 3, 1], [7, 9, 4]])
P, L, U = lu(A)

print("原矩阵A:\n", A)
print("置换矩阵P:\n", P)
print("下三角矩阵L:\n", L)
print("上三角矩阵U:\n", U)
print("验证 P×L×U = A:\n", P @ L @ U)

3. 应用场景

• 解线性方程组：Ax = b → LUx = b
  • 先解 Ly = Pb（前向替代）
  • 再解 Ux = y（后向替代）
• 计算行列式：det(A) = det(P) × ∏Uᵢᵢ
• 求逆矩阵：A⁻¹ = U⁻¹L⁻¹Pᵀ

4. 计算优势

python

# LU分解后求解方程组更高效
import numpy as np
from scipy.linalg import lu_factor, lu_solve

A = np.array([[2, 5, 8], [5, 3, 1], [7, 9, 4]])
b = np.array([1, 2, 3])

# 传统方法（每次都要计算）
x1 = np.linalg.solve(A, b)

# LU分解方法（分解一次，多次使用）
lu, piv = lu_factor(A)
x2 = lu_solve((lu, piv), b)

二、QR分解

1. 基本概念

• 含义：将矩阵分解为正交矩阵（Q）和上三角矩阵（R）的乘积
• 公式：A = Q × R
• 性质：QᵀQ = I（正交性）

2. 几何意义

• Q矩阵：一组标准正交基，是A列空间的正交基
• R矩阵：A在Q基下的坐标表示
• 过程：通过Gram-Schmidt正交化得到

python

import numpy as np

A = np.array([[1, 2], [3, 4], [5, 6]])
Q, R = np.linalg.qr(A)

print("原矩阵A (3×2):\n", A)
print("正交矩阵Q (3×2):\n", Q)
print("Q的列正交性 QᵀQ:\n", Q.T @ Q)  # 近似单位矩阵
print("上三角矩阵R (2×2):\n", R)
print("验证 Q×R = A:\n", Q @ R)

3. 应用场景

• 最小二乘法：QR分解是求解最小二乘问题的稳定方法python# 解最小二乘问题 Ax ≈ b A = np.array([[1, 1], [1, 2], [1, 3]]) b = np.array([1, 3, 2]) Q, R = np.linalg.qr(A) # 解 Rx = Qᵀb x = np.linalg.solve(R, Q.T @ b) print(“最小二乘解:”, x)
• 特征值计算：QR算法用于计算特征值
• 正交投影：Q的列张成A的列空间

4. 计算优势

• 数值稳定性：比正规方程 (AᵀA)x = Aᵀb 更稳定
• 处理秩亏矩阵：可以检测和处理列相关性

三、SVD分解（最重要！）

1. 基本概念

• 含义：将矩阵分解为三个特殊矩阵的乘积
• 公式：A = U × Σ × Vᵀ
  • U：m×m 正交矩阵（左奇异向量）
  • Σ：m×n 对角矩阵（奇异值，降序排列）
  • Vᵀ：n×n 正交矩阵的转置（右奇异向量）

2. 几何意义

python

import numpy as np

A = np.array([[3, 1], [1, 3]])
U, S, Vh = np.linalg.svd(A)

print("原矩阵A:\n", A)
print("左奇异向量U:\n", U)
print("奇异值S:", S)
print("右奇异向量Vᵀ:\n", Vh)

# 重构矩阵
Sigma = np.zeros(A.shape)
np.fill_diagonal(Sigma, S)
A_reconstructed = U @ Sigma @ Vh
print("重构的A:\n", A_reconstructed)

几何解释：

1. Vᵀ：在输入空间（ℝⁿ）中进行旋转/反射
2. Σ：沿坐标轴进行缩放（奇异值决定缩放因子）
3. U：在输出空间（ℝᵐ）中进行旋转/反射

3. SVD的四个基本子空间

• 行空间：V的前r列（r=秩）
• 零空间：V的后n-r列
• 列空间：U的前r列
• 左零空间：U的后m-r列

4. 应用场景（非常广泛！）

(1) 低秩近似（数据压缩）

python

# 图像压缩示例
import numpy as np
import matplotlib.pyplot as plt

# 创建测试矩阵（模拟图像）
image = np.random.rand(100, 100)

U, S, Vh = np.linalg.svd(image)

# 不同秩的近似
ranks = [5, 20, 50]
fig, axes = plt.subplots(1, 4, figsize=(15, 4))

axes[0].imshow(image, cmap='gray')
axes[0].set_title("原图 (秩100)")

for i, rank in enumerate(ranks, 1):
    # 使用前k个奇异值重构
    U_k = U[:, :rank]
    S_k = np.diag(S[:rank])
    Vh_k = Vh[:rank, :]
    approx = U_k @ S_k @ Vh_k
    
    axes[i].imshow(approx, cmap='gray')
    axes[i].set_title(f"秩{rank}近似")
    compression_ratio = (rank * (100 + 100 + 1)) / (100 * 100)
    axes[i].set_xlabel(f"压缩比: {compression_ratio:.2%}")

(2) 推荐系统（协同过滤）

python

# 用户-物品评分矩阵
ratings = np.array([
    [5, 3, 0, 1],
    [4, 0, 0, 1],
    [1, 1, 0, 5],
    [1, 0, 0, 4],
    [0, 1, 5, 4]
])

U, S, Vh = np.linalg.svd(ratings, full_matrices=False)

# 使用前2个奇异值
k = 2
U_k = U[:, :k]
S_k = np.diag(S[:k])
Vh_k = Vh[:k, :]

# 预测缺失评分
predicted = U_k @ S_k @ Vh_k
print("原始评分:\n", ratings)
print("预测评分:\n", predicted.round(2))

(3) 主成分分析（PCA）

python

# PCA本质上是中心化数据的SVD
data = np.random.randn(100, 5)  # 100个样本，5个特征

# 中心化
data_centered = data - data.mean(axis=0)

# SVD（相当于PCA）
U, S, Vh = np.linalg.svd(data_centered, full_matrices=False)

# 主成分（Vh的前k行）
principal_components = Vh[:2, :]  # 前两个主成分
print("前两个主成分方向:\n", principal_components)

# 投影到主成分空间
scores = data_centered @ principal_components.T

(4) 求解最小二乘问题

python

# SVD求解病态方程组
A = np.array([[1, 1], [1, 1.0001], [1, 1.0002]])
b = np.array([2, 2.0001, 1.9999])

U, S, Vh = np.linalg.svd(A, full_matrices=False)

# 使用SVD求伪逆
S_inv = np.diag(1/S)
A_pinv = Vh.T @ S_inv @ U.T
x = A_pinv @ b

print("SVD解:", x)
print("直接求逆解:", np.linalg.lstsq(A, b, rcond=None)[0])

(5) 自然语言处理（潜在语义分析）

python

# 词-文档矩阵
term_doc_matrix = np.array([
    [2, 0, 1, 0, 0],  # 文档1：词频
    [1, 1, 0, 0, 1],  # 文档2
    [0, 0, 2, 1, 0],  # 文档3
    [1, 0, 0, 2, 1]   # 文档4
])

U, S, Vh = np.linalg.svd(term_doc_matrix)

# 潜在语义空间
latent_topics = 2
U_k = U[:, :latent_topics]    # 词在潜在主题上的分布
V_k = Vh[:latent_topics, :].T  # 文档在潜在主题上的分布

5. SVD的特殊性质

python

# 1. 奇异值与特征值的关系
A = np.random.rand(5, 3)
U, S, Vh = np.linalg.svd(A)

# 奇异值的平方是AᵀA的特征值
ATA = A.T @ A
eigenvalues = np.linalg.eigvals(ATA)
print("奇异值平方:", S**2)
print("AᵀA的特征值:", np.sort(eigenvalues)[::-1])

# 2. 矩阵的秩
rank = np.sum(S > 1e-10)  # 非零奇异值的个数
print("矩阵的秩:", rank)

# 3. 矩阵的范数
frobenius_norm = np.sqrt(np.sum(S**2))  # Frobenius范数
spectral_norm = S[0]                    # 2-范数（最大奇异值）

6. 经济型SVD

python

# 当矩阵不是满秩时，使用经济型SVD节省空间
A = np.random.rand(100, 10)  # 100×10矩阵，最大秩为10

# full_matrices=False 只计算必要的奇异向量
U, S, Vh = np.linalg.svd(A, full_matrices=False)

print("U形状:", U.shape)  # (100, 10) 而不是 (100, 100)
print("S长度:", len(S))   # 10
print("Vh形状:", Vh.shape) # (10, 10)

三大分解对比表

实际选择建议

1. 解线性方程组：优先使用LU分解（方阵）或QR分解（非方阵）
2. 最小二乘问题：使用QR分解或SVD分解
3. 数据压缩/降维：使用SVD分解
4. 病态问题/数值稳定：使用SVD分解
5. 特征值计算：QR算法（基于QR分解）

关键总结

• LU分解：高斯消元的矩阵形式，适合解方程
• QR分解：正交化过程，适合最小二乘
• SVD分解：最强大的分解，揭示了矩阵的本质结构，是数据科学中的”瑞士军刀”

SVD是数据科学中最重要的矩阵分解，它不仅能解决数值计算问题，还能揭示数据的潜在结构，广泛应用于机器学习、信号处理、图像压缩等领域。