跟人工智能学线性代数

第 14 章 奇异值分解 SVD ​

学习目标 ​

  • 理解对角化的局限,以及为什么 SVD 适用于任意矩阵
  • 掌握 SVD 的定义
  • 从 的特征值推导奇异值,理解 、、 的来源
  • 理解 SVD 的几何意义:旋转 → 缩放 → 旋转
  • 掌握截断 SVD 与低秩近似
  • 用 SVD 做图片压缩,并用 NumPy 与 PyTorch 计算

14.1 对角化不够用 ​

第 12 章的对角化 有两个限制:

  1. 只有方阵才能对角化;
  2. 方阵还要有 个线性无关的特征向量(可能缺)。

现实中,数据矩阵几乎都是非方阵(样本数 ≠ 特征数)。能不能给任意矩阵一个「对角分解」?

出路藏在 里。回顾第 6 章练习: 一定是对称矩阵;对称矩阵的特征值全是实数,而且一定有 个两两正交的特征向量(第 12 章预告过)。所以 一定能对角化,而且用正交矩阵就能做到。SVD 就是把这个事实包装成对 本身的分解。

14.2 SVD 的定义 ​

奇异值分解:任意 矩阵 都可以写成

其中:

  • 是 正交矩阵(列叫左奇异向量);
  • 是 正交矩阵(列叫右奇异向量);
  • 是 「对角」矩阵,对角线上的数是奇异值

( 是秩,其余奇异值为 0)。

正交矩阵 满足 ,几何上就是旋转或镜像:它不改变长度与夹角,只改变方向。所以 SVD 说的是:任何线性变换 = 先旋转()→ 沿坐标轴缩放()→ 再旋转()。

SVD 几何意义

例:对角矩阵的 SVD 就是它自己。

python
import matplotlib.pyplot as plt
from matplotlib import font_manager
for f in font_manager.findSystemFonts():
    if any(k in f for k in ("NotoSansCJK", "NotoSansSC", "wqy", "SimHei", "msyh", "PingFang")):
        font_manager.fontManager.addfont(f)
plt.rcParams["font.sans-serif"] = ["Noto Sans CJK SC", "WenQuanYi Zen Hei", "SimHei", "Microsoft YaHei", "PingFang SC"]
plt.rcParams["axes.unicode_minus"] = False

import numpy as np

A = np.array([[2.0, 1.0],
              [1.0, 2.0]])
U, s, Vt = np.linalg.svd(A)

theta = np.linspace(0, 2 * np.pi, 200)
circle = np.column_stack([np.cos(theta), np.sin(theta)])
ellipse = circle @ A.T

fig, ax = plt.subplots(figsize=(6, 6))
ax.axhline(0, color="black", linewidth=1)
ax.axvline(0, color="black", linewidth=1)
ax.grid(True, linestyle=":", alpha=0.5)
ax.set_xlim(-3.5, 3.5)
ax.set_ylim(-3.5, 3.5)
ax.set_aspect("equal")

ax.plot(circle[:, 0], circle[:, 1], color="tab:blue", linewidth=2, label="单位圆")
ax.plot(ellipse[:, 0], ellipse[:, 1], color="tab:red", linewidth=2, label="A 变换后:椭圆")

# 半轴:奇异值方向
for i in range(2):
    v = Vt[i]
    u = U[:, i]
    ax.plot([0, s[i] * u[0]], [0, s[i] * u[1]], color="tab:green",
            linestyle="--", linewidth=1.5)
    ax.annotate(f"σ{i+1} = {s[i]:.0f}", s[i] * u, textcoords="offset points",
                xytext=(8, 6), color="tab:green", fontsize=12)

ax.set_title("SVD:单位圆被旋转、拉伸成椭圆")
ax.legend()
fig.savefig("book/public/figs/fig14-svd_ellipse.png", dpi=120, bbox_inches="tight")
plt.close(fig)

14.3 推导:奇异值来自 的特征值 ​

第 1 步: 是 对称矩阵,可正交对角化:

的列 是 的特征向量, 对角线是特征值 (对称半正定矩阵特征值非负)。定义

即奇异值 = 特征值的算术平方根。

第 2 步:构造 的列。对 ,令

验证 两两正交:

因为 。所以 是正交的单位向量组,补足成整个 维空间的正交基,得到 。

第 3 步:验证分解成立。按定义,

两边「除以」(即乘 ),得到 ——对矩阵的每一列都成立,所以整个等式成立。

完整例子 ​

特征值 ,奇异值 、。,特征向量 、。构造:

验证: ✓

python
import numpy as np

A = np.array([[0.0, 2.0],
              [3.0, 0.0]])
U, s, Vt = np.linalg.svd(A)
print("U:\n", U)
print("奇异值 s:", s)
print("Vt:\n", Vt)

# 手动验证 AᵀA 的特征值
AtA = A.T @ A
print("AᵀA 特征值:", np.linalg.eigvalsh(AtA))
print("奇异值 = 特征值开方:", np.sqrt(np.linalg.eigvalsh(AtA)))

# 重建 A = U Σ Vᵀ
Sigma = np.diag(s)
print("重建 A:\n", U @ Sigma @ Vt)

运行输出:

text
U:
 [[ 0. -1.]
 [-1.  0.]]
奇异值 s: [3. 2.]
Vt:
 [[-1. -0.]
 [-0. -1.]]
AᵀA 特征值: [4. 9.]
奇异值 = 特征值开方: [2. 3.]
重建 A:
 [[0. 2.]
 [3. 0.]]

注意 eigvalsh 返回的特征值按升序排列,所以打印出来是 ,对应 ;而 svd 的奇异值按降序 。另一个细节:库返回的 、 与 14.3 手算的选择符号不同( 两列都取了反方向, 同理)。这没关系——同时翻转 的某列与 的对应列,乘积不变,两种都是合法的 SVD。

14.4 几何意义:单位圆 → 椭圆 ​

SVD 的三个因子对应三步几何操作:

  1. :旋转,把标准基转到右奇异向量的方向;
  2. :沿坐标轴拉伸 倍;
  3. :再旋转到左奇异向量的方向。

因此,单位圆经过 变成椭圆,半轴长就是奇异值;单位球变成椭球。第一个奇异值 是 的最大拉伸倍数。

python
import numpy as np

A = np.array([[2.0, 1.0],
              [1.0, 2.0]])
U, s, Vt = np.linalg.svd(A)
print("奇异值(椭圆半轴):", s)
print("最大拉伸倍数 σ1 =", s[0])
print("验证:||A x|| 最大不超过 σ1 * ||x||")

运行输出:

text
奇异值(椭圆半轴): [3. 1.]
最大拉伸倍数 σ1 = 2.9999999999999996
验证:||A x|| 最大不超过 σ1 * ||x||

14.5 低秩近似:丢掉小奇异值 ​

SVD 还提供「压缩」的数学依据。把分解按列展开:

每个 都是一个秩 1 矩阵(一列 × 一行)。 被拆成一串秩 1 矩阵之和,贡献大小由 决定。

截断 SVD(低秩近似):只保留前 项:

Eckart-Young 定理(不证明):在所有秩不超过 的矩阵中, 是离 最近的(按 Frobenius 范数)。也就是说,丢掉最小的奇异值,是信息损失最少的一种压缩方式——这正是 PCA、图像压缩、推荐系统背后的原理。

14.6 实战:用 SVD 压缩图片 ​

先生成一张合成图(避免依赖网络下载),再用 SVD 压缩:

python
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import font_manager
for f in font_manager.findSystemFonts():
    if any(k in f for k in ("NotoSansCJK", "NotoSansSC", "wqy", "SimHei", "msyh", "PingFang")):
        font_manager.fontManager.addfont(f)
plt.rcParams["font.sans-serif"] = ["Noto Sans CJK SC", "WenQuanYi Zen Hei", "SimHei", "Microsoft YaHei", "PingFang SC"]
plt.rcParams["axes.unicode_minus"] = False

import torch

# 合成 128×128 灰度图:渐变 + 两个圆
yy, xx = np.mgrid[0:128, 0:128]
img = 0.5 * np.sin(xx / 8) + 0.5 * np.cos(yy / 10)
img += 1.5 * np.exp(-((xx - 50) ** 2 + (yy - 60) ** 2) / 400)
img += 1.5 * np.exp(-((xx - 85) ** 2 + (yy - 35) ** 2) / 500)
img -= img.min()
img /= img.max()

# 1. NumPy SVD
U, s, Vt = np.linalg.svd(img, full_matrices=False)
print("奇异值个数:", len(s))
print("前 5 个奇异值:", np.round(s[:5], 3))

# 2. PyTorch SVD 对比
img_t = torch.tensor(img)
U_t, s_t, Vt_t = torch.linalg.svd(img_t)
print("torch 前 5 个奇异值:", np.round(s_t[:5].numpy(), 3))

# 3. 截断重建
k = 8
img_k = (U[:, :k] * s[:k]) @ Vt[:k, :]   # 等价于 U_k Σ_k Vt_k

fig, axes = plt.subplots(1, 3, figsize=(13, 4))
axes[0].imshow(img, cmap="gray")
axes[0].set_title(f"原图(128×128={img.size} 个像素)")
axes[1].imshow(img_k, cmap="gray")
axes[1].set_title(f"保留前 {k} 个奇异值")
axes[2].semilogy(s)
axes[2].set_title("奇异值(对数坐标)")
axes[2].set_xlabel("i")
axes[2].set_ylabel("σi")
axes[2].grid(True, linestyle=":", alpha=0.5)
for ax in axes[:2]:
    ax.axis("off")
fig.savefig("book/public/figs/fig14-svd_compress.png", dpi=120, bbox_inches="tight")
plt.close(fig)

# 4. 误差与压缩率
err = np.linalg.norm(img - img_k) / np.linalg.norm(img)
print(f"k={k}: 相对误差 = {err:.3f}")
print(f"存储:原图 {img.size} 个数,压缩后 U、s、Vt 共 {U.shape[0]*k + k + k*img.shape[1]} 个数")

运行输出:

text
奇异值个数: 128
前 5 个奇异值: [55.807  7.482  6.493  2.56   0.   ]
torch 前 5 个奇异值: [55.807  7.482  6.493  2.56   0.   ]
k=8: 相对误差 = 0.000
存储:原图 16384 个数,压缩后 U、s、Vt 共 2056 个数

128×128 = 16384 个像素,只用 8 个奇异值(2056 个数,约八分之一)就能几乎无损重建(误差 0.000)——因为这张合成图由几个可分离的函数(正弦、余弦、二维高斯)相加而成,秩不超过 4,信息本来就集中在前几个奇异值里。真实照片不是可分离的,误差会明显更大;把 img 换成一张真实灰度图,就能看到「奇异值下降较慢,需要更多 才能保真」的典型情形。

14.7 AI 视角:SVD 与主成分分析 ​

PCA(主成分分析)本质上就是数据中心化后的 SVD:

  1. 把数据矩阵每一列减去均值(中心化);
  2. 对中心化矩阵做 SVD;
  3. 右奇异向量 的列就是主方向,奇异值平方正比于该方向上的方差,投影 就是降维后的数据。

SVD 还用于:推荐系统的隐因子分解(Netflix 奖)、文本的潜在语义分析(LSA)、神经网络的模型压缩、图像与视频压缩。它被称为「线性代数里最重要的分解」——因为任何矩阵都能分解,而且分解出的每一块都有明确的几何与统计含义。

动手实践 ​

  1. 写出 的 SVD(提示:它就是对角矩阵本身),并用 np.linalg.svd 验证。
  2. 对 手算 的特征值与奇异值,并与 14.3 的例子对照。
  3. 用 SVD 求 的最大拉伸倍数,并用随机向量数值验证 。
  4. 把 14.6 的 改成 2、16、32,观察重建图的误差变化,画出误差随 的曲线。
  5. 用 torch.linalg.svd 计算一个 3×5 矩阵的奇异值,观察奇异值个数与矩阵形状的关系(提示:min(m, n) 个)。

常见错误 ​

错误写法/理解原因
把 SVD 写成 第三项是 ,不是 ; 一般不等于
认为奇异值必须全部非零奇异值个数是 ,其中非零个数 = 秩
把奇异值当作特征值奇异值是 特征值的平方根,不是 的特征值(方阵时
截断时直接 U[:, :k] @ s[:k] @ Vt[:k] 尺寸不匹配需要用广播 (U[:, :k] * s[:k]) @ Vt[:k, :]
认为 PCA 不需要中心化不中心化时第一个主方向会被均值带偏;PCA 先减均值再做 SVD

章末练习 ​

基础

  1. 写出 的奇异值(按降序)。
  2. 用 np.linalg.svd 求 的 、、,验证重建。
  3. 解释为什么 一定可正交对角化。

提高

  1. 对 用 NumPy 求 SVD,手算 的特征值,验证奇异值关系。
  2. 用秩 1 近似 近似 ,计算相对误差。

挑战

  1. 证明 的特征值非负:对特征向量 ,有 。补全推导。
  2. 用 SVD 对 14.6 的图做「彩色」版本:把 RGB 三个通道分别 SVD 压缩再合成,比较效果。

章末自测 ​

  1. SVD 分解的形式是?
    • A.
    • B.
    • C.
    • D.
  2. 奇异值来自哪个矩阵的特征值?
    • A.
    • B.
    • C.
    • D.
  3. 正交矩阵满足?
    • A.
    • B.
    • C.
    • D.
  4. 矩阵的奇异值个数是?
    • A.
    • B.
    • C.
    • D.
  5. 奇异值 的几何意义是?
    • A. 矩阵的最大拉伸倍数
    • B. 矩阵的行列式
    • C. 矩阵的秩
    • D. 矩阵的迹
  6. 单位圆经过 变成?
    • A. 正方形
    • B. 椭圆,半轴为奇异值
    • C. 直线
    • D. 圆形
  7. 截断 SVD 保留前 项,得到的是?
    • A. 精确的
    • B. 秩不超过 的最佳近似(按 Frobenius 范数)
    • C. 的特征分解
    • D. 随机矩阵
  8. PCA 的核心步骤是?
    • A. 中心化后做 SVD
    • B. 直接做矩阵乘法
    • C. 求行列式
    • D. 转置
  9. 秩 1 矩阵的形状是?
    • A. 一列 × 一行()
    • B. 一行 × 一列
    • C. 对角矩阵
    • D. 单位矩阵
  10. 的奇异值是?
    • A.
    • B.
    • C.
    • D.