SVD的计算和低秩近似

系列第三篇：从“为什么存在”走向“怎样算、怎样取舍”

一、从前两篇走到这里

第一篇把矩阵看成一台“变形机器”：任意三维方形矩阵$A$都会把输入空间中的单位球变成输出空间中的椭球。椭球的主轴来自一组特殊的输入方向$v_i$；它们在输入端彼此正交，经过$A$之后仍然彼此正交。每个方向被缩放为

$$A{v}_i={\sigma }_i{u}_i,$$

其中，$v_i$是右奇异向量，${\sigma }_i$是缩放倍数，$u_i$是变换后的单位向量。把所有方向并排放在一起，就得到

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

顺着矩阵乘法从右往左看，$V^T$先把输入方向对齐到主轴，$\Sigma$沿各条主轴分别拉伸，$U$再把它们摆到输出空间中的位置。任何线性变换，都可以拆成“换一组正交坐标轴、沿轴缩放、再换一组正交坐标轴”。

第二篇继续追问：这样的方向为什么一定存在，又该从哪里找到？从单位球面上“寻找被$A$拉伸最长的方向”出发，问题被转换成：

$$\underset{\Vert v\Vert =1}{\max }\Vert Av{\Vert }^2=\underset{\Vert v\Vert =1}{\max }{v}^T{A}^{T}Av.$$

拉格朗日乘子法把这个极值问题化成了特征方程

$$A^{T}Av=\lambda v.$$

由于$A^{T}A$是对称半正定矩阵，它拥有一组标准正交特征向量，并且所有特征值都非负。由此可以依次构造 SVD 的三个组成部分：

$$v_i:A^{T}A\text{ 的单位特征向量},{\sigma }_i=\sqrt{{\lambda }_i},u_i=\frac{A{v}_i}{{\sigma }_i}.$$

至此，前两篇建立的几何直觉与存在性证明已经落到具体的构造方法上。本文将把这一构造整理成计算步骤，再用一个具体矩阵走完整个过程；最后把$U\Sigma V^T$展开为一层层秩为$1$的“信息薄片”，由此引出低秩近似。

二、把推导收成计算配方

1. 给定矩阵$A\in {\mathbb{R}}^{m\times n}$，第二篇的推导可以整理为下面这套构造规则：

要求的量	计算方法	延续前文的几何含义
右奇异向量$v_i$	求$A^{T}A$的单位特征向量	输入空间中那组“正交映成正交”的方向
奇异值${\sigma }_i$	对相应特征值开平方：${\sigma }_i=\sqrt{{\lambda }_i}$	对应方向的缩放倍数，也是椭球主轴的半长
左奇异向量$u_i$	对${\sigma }_i>0$，令$u_i=A{v}_i/{\sigma }_i$	输入主方向经过$A$后，在输出空间中的单位方向

将$A^{T}A$的特征值按从大到小排列，

$${\lambda }_1\ge {\lambda }_2\ge \cdots \ge 0,$$

对应的奇异值也自然满足

$${\sigma }_1\ge {\sigma }_2\ge \cdots \ge 0.$$

排序并不是形式上的整理。到了低秩近似中，排在前面的方向还会成为最先被保留的信息层。

完整 SVD 的标准形式为：

$$A_{m\times n}=U_{m\times m}{\Sigma }_{m\times n}{V}_{n\times n}^T.$$

$U$与$V$都是正交矩阵：它们的列分别构成输出空间和输入空间的一组标准正交基。因此

$$U^{T}U=U{U}^T=I_m,V^{T}V=V{V}^T=I_n.$$

$\Sigma$是一个$m\times n$的“矩形对角矩阵”，奇异值依次放在主对角线上，其余位置均为零。SVD 给出$\min (m,n)$个对角位置，其中真正非零的奇异值有

$$r=\mathrm{rank}(A)$$

个。若${\sigma }_i=0$，则$A{v}_i=0$，对应的$v_i$落在$A$的零空间中。

对非零奇异值，我们可以直接用$u_i=A{v}_i/{\sigma }_i$构造左奇异向量；如果完整的$U$还缺少列，就在输出空间中寻找正交补，把已有的$u_i$补成一组完整的标准正交基。完整 SVD 中$U$和$V$都是满秩的正交方阵；矩阵$A$是否满秩，则由非零奇异值的个数决定。

因此计算的顺序是：

1. 从$A^{T}A$求出按特征值排序的$v_i$；
2. 令${\sigma }_i=\sqrt{{\lambda }_i}$；
3. 对每个${\sigma }_i>0$，用$u_i=A{v}_i/{\sigma }_i$顺着$v_i$构造$u_i$；
4. 对剩余方向使用正交补。

这样每一组三元组$(u_i,{\sigma }_i,v_i)$都由关系$A{v}_i={\sigma }_i{u}_i$直接拴在一起，符号会自动匹配。

三、一个完整的计算示例

考虑矩阵

$$A=\left[\begin{array}{cc}3& 2\\ 2& 3\\ 2& -2\end{array}\right].$$

它把二维输入空间映射到三维输出空间，所以

$$A\in {\mathbb{R}}^{3\times 2},U\in {\mathbb{R}}^{3\times 3},\Sigma \in {\mathbb{R}}^{3\times 2},V\in {\mathbb{R}}^{2\times 2}.$$

（一）计算 A^TA

右奇异向量在二维输入空间中，因此先构造$2\times 2$的对称矩阵$A^{T}A$：

$$\begin{array}{rl}{A}^{T}A& =\left[\begin{array}{ccc}3& 2& 2\\ 2& 3& -2\end{array}\right]\left[\begin{array}{cc}3& 2\\ 2& 3\\ 2& -2\end{array}\right]=\left[\begin{array}{cc}17& 8\\ 8& 17\end{array}\right]\end{array}$$

（二）求右奇异向量与奇异值

特征方程为

$$\begin{array}{rl}\det (A^{T}A-\lambda I)& =\left|\begin{array}{cc}17-\lambda & 8\\ 8& 17-\lambda \end{array}\right|\\ & =(17-\lambda {)}^2-64\\ & ={\lambda }^2-34\lambda +225=0.\end{array}$$

解得

$${\lambda }_1=25,{\lambda }_2=9.$$

因此奇异值为

$${\sigma }_1=\sqrt{25}=5,{\sigma }_2=\sqrt{9}=3.$$

这两个数正是$A$沿两条输入主方向的缩放倍数。由于都不为零，矩阵$A$的秩为$2$。换成第一篇的几何语言，它没有把二维输入平面中的任何主方向彻底压扁。

接着求单位特征向量。

当${\lambda }_1=25$时，解

$$(A^{T}A-25I)v=0$$

可取

$$v_1=\frac{1}{\sqrt{2}}\left[\begin{array}{c}1\\ 1\end{array}\right].$$

当${\lambda }_2=9$时，解

$$(A^{T}A-9I)v=0$$

可取

$$v_2=\frac{1}{\sqrt{2}}\left[\begin{array}{c}1\\ -1\end{array}\right].$$

两者互相垂直，组成右奇异矩阵

$$V=\left[\begin{array}{cc}\frac{1}{\sqrt{2}}& \frac{1}{\sqrt{2}}\\ \frac{1}{\sqrt{2}}& -\frac{1}{\sqrt{2}}\end{array}\right],V^T=\left[\begin{array}{cc}\frac{1}{\sqrt{2}}& \frac{1}{\sqrt{2}}\\ \frac{1}{\sqrt{2}}& -\frac{1}{\sqrt{2}}\end{array}\right].$$

这里$V^T=V$，只是这个具体例子的巧合，并不是 SVD 的一般性质。

按照奇异值的位置构造

$$\Sigma =\left[\begin{array}{cc}5& 0\\ 0& 3\\ 0& 0\end{array}\right].$$

（三）顺着 V 构造 U

第一组输出方向为

$$\begin{array}{rl}{u}_1=\frac{A{v}_1}{{\sigma }_1}& =\frac{1}{5}\left[\begin{array}{c}5/\sqrt{2}\\ 5/\sqrt{2}\\ 0\end{array}\right]=\left[\begin{array}{c}1/\sqrt{2}\\ 1/\sqrt{2}\\ 0\end{array}\right]\end{array}$$

它满足

$$A{v}_1=5u_1$$

这句话与第一篇中的几何图像完全对应：输入方向$v_1$被$A$拉伸$5$倍，落到三维输出空间中的方向$u_1$上。

第二组输出方向为

$$\begin{array}{rl}{u}_2& =\frac{A{v}_2}{{\sigma }_2}=\frac{1}{3}\left[\begin{array}{c}1/\sqrt{2}\\ -1/\sqrt{2}\\ 4/\sqrt{2}\end{array}\right]=\left[\begin{array}{c}1/(3\sqrt{2})\\ -1/(3\sqrt{2})\\ 4/(3\sqrt{2})\end{array}\right]\end{array}$$

同样有

$$A{v}_2=3u_2.$$

$u_1$、$u_2$已经是三维输出空间中的两个正交单位向量。完整 SVD 需要把$U$补成$3\times 3$的正交矩阵，所以还要找一个单位向量$u_3$，使它同时垂直于$u_1$和$u_2$。

设

$$u_3=\left[\begin{array}{c}x\\ y\\ z\end{array}\right]$$

“同时垂直”意味着它与$u_1$、$u_2$的内积都等于零。先由$u_1^T{u}_3=0$得

$$\begin{array}{rl}{u}_1^T{u}_3& =\left[\begin{array}{ccc}1/\sqrt{2}& 1/\sqrt{2}& 0\end{array}\right]\left[\begin{array}{c}x\\ y\\ z\end{array}\right]\\ & =\frac{x+y}{\sqrt{2}}=0,\end{array}$$

所以

$$y=-x.$$

再由$u_2^T{u}_3=0$得

$$\begin{array}{rl}{u}_2^T{u}_3& =\left[\begin{array}{ccc}1/(3\sqrt{2})& -1/(3\sqrt{2})& 4/(3\sqrt{2})\end{array}\right]\left[\begin{array}{c}x\\ y\\ z\end{array}\right]\\ & =\frac{x-y+4z}{3\sqrt{2}}=0.\end{array}$$

代入$y=-x$：

$$2x+4z=0⟹z=-\frac{x}{2}.$$

因此，所有同时垂直于$u_1$、$u_2$的向量都可以写成

$$u_3=x\left[\begin{array}{c}1\\ -1\\ -1/2\end{array}\right].$$

最后要求$u_3$是单位向量，即$\Vert u_3\Vert =1$：

$$\begin{array}{rl}\Vert u_3{\Vert }^2& =x^2+(-x{)}^2+{\left(-\frac{x}{2}\right)}^2\\ & =\frac{9}{4}{x}^2=1.\end{array}$$

由此得到$x=\pm 2/3$。取$x=2/3$，便有

$$u_3=\left[\begin{array}{c}2/3\\ -2/3\\ -1/3\end{array}\right].$$

于是有：

$$U=\left[\begin{array}{ccc}1/\sqrt{2}& 1/(3\sqrt{2})& 2/3\\ 1/\sqrt{2}& -1/(3\sqrt{2})& -2/3\\ 0& 4/(3\sqrt{2})& -1/3\end{array}\right].$$

（四）拼回并验证分解

现在三部分已经齐全：

$$U\approx \left[\begin{array}{ccc}0.7071& 0.2357& 0.6667\\ 0.7071& -0.2357& -0.6667\\ 0& 0.9428& -0.3333\end{array}\right],\Sigma =\left[\begin{array}{cc}5& 0\\ 0& 3\\ 0& 0\end{array}\right],V^T\approx \left[\begin{array}{cc}0.7071& 0.7071\\ 0.7071& -0.7071\end{array}\right].$$

将它们相乘，得到

$$U\Sigma V^T=\left[\begin{array}{cc}3& 2\\ 2& 3\\ 2& -2\end{array}\right]=A.$$

手算过程由此闭合。从$A^{T}A$找到输入主方向，对特征值开根号得到缩放倍数，再把每条输入方向送过$A$，便能构造出与之配对的输出方向。第二篇的存在性推导由此转化为一套完整的计算路径。

四、从矩阵乘积到信息薄片

前三步“$V^T$对齐、$\Sigma$缩放、$U$摆放”解释了$A$如何作用于空间。若把矩阵乘积按列展开，同一个 SVD 还会显露出另一层结构：

$$A=U\Sigma V^T=\sum _{i=1}^r{\sigma }_i{u}_i{v}_i^T.$$

这个求和式把第一篇开头的“信息分解”写成了精确的代数形式：

$$A=\underset{\text{第一层，权重最大}}{\underbrace{{\sigma }_1{u}_1{v}_1^T}}+\underset{\text{第二层}}{\underbrace{{\sigma }_2{u}_2{v}_2^T}}+\cdots +\underset{\text{第 }r\text{ 层}}{\underbrace{{\sigma }_r{u}_r{v}_r^T}}.$$

其中，每个外积$u_i{v}_i^T$都是一个秩$1$矩阵；乘上权重${\sigma }_i$后，就形成一张最薄的“信息薄片”。$v_i$决定这一层从输入空间的哪个方向读取信息，$u_i$决定它写入输出空间的哪个方向，${\sigma }_i$则决定这一层的强弱。作用在向量$x$上时，矩阵的作用可以写成

$$Ax=\sum _{i=1}^r{\sigma }_i(v_i^{T}x)u_i.$$

这里，$v_i^{T}x$先读出$x$在第$i$个输入方向上的分量，${\sigma }_i$对它进行缩放，$u_i$再把结果写入输出空间。复杂的矩阵变换由此被拆成若干条互不干扰的信息通道。

接下来需要回答一个关键问题：这些薄片究竟有多大，又该怎样比较它们？为此，把矩阵看成由所有元素排成的长向量，并定义 Frobenius 范数与 Frobenius 内积：

$$\Vert A{\Vert }_F=\sqrt{\sum _{i=1}^m\sum _{j=1}^n{a}_{ij}^2},⟨X,Y{⟩}_F=\sum _{i=1}^m\sum _{j=1}^n{x}_{ij}{y}_{ij}.$$

在这个定义下，不同薄片两两正交：

$${⟨u_i{v}_i^T,u_j{v}_j^T⟩}_F=(u_i^T{u}_j)(v_i^T{v}_j)=0,i\ne j.$$

每张加权薄片的 Frobenius 范数则恰好等于相应的奇异值：

$$\Vert {\sigma }_i{u}_i{v}_i^T{\Vert }_F={\sigma }_i.$$

因此，按${\sigma }_1\ge {\sigma }_2\ge \cdots \ge {\sigma }_r>0$排列奇异值，也就等于按薄片的 Frobenius 范数从大到小排列。由于这些薄片彼此正交，矩阵整体的范数平方还可以像勾股定理一样逐层相加：

$$\Vert A{\Vert }_F^2=\Vert \sum _{i=1}^r{\sigma }_i{u}_i{v}_i^T{\Vert }_F^2=\sum _{i=1}^r{\sigma }_i^2.$$

这说明${\sigma }_i^2$是第$i$张薄片对矩阵总能量的贡献。既然各层已经按贡献从大到小排好，一个自然的问题随之出现：如果只保留最前面的几层，能否用更简单的矩阵抓住$A$的主要结构？

五、为什么可以只保留前几层

上一节已经把矩阵拆成按贡献从大到小排列的正交薄片。低秩近似顺着这个结构做一件很直接的事：保留前面的强模式，舍弃后面的弱模式，用较低秩的矩阵逼近原矩阵。

（一）完整、紧致与截断 SVD

在真正舍弃信息之前，需要先区分三种容易混淆的写法。

完整 SVD 会把$U$和$V$的正交基全部补齐；紧致 SVD 只保留与非零奇异值配对的$r$组方向，去掉纯粹用于补基的零方向。二者虽然尺寸不同，却都能精确重构$A$。

在前面的例子中，$A\in {\mathbb{R}}^{3\times 2}$且$r=2$。完整 SVD 中的$U$有三列，但第三列只用于补全输出空间的正交基。删去这一列得到

$$A=U_2\left[\begin{array}{cc}5& 0\\ 0& 3\end{array}\right]V^T,$$

这仍然是等式，因此属于紧致 SVD，而不是低秩近似。

只有继续取$k<r$，舍弃一部分非零奇异值，才得到截断 SVD：

$$\boxed{A_k=\sum _{i=1}^k{\sigma }_i{u}_i{v}_i^T=U_k{\Sigma }_k{V}_k^T,A\approx A_k}$$

$A_k$的秩不超过$k$。它保留了前$k$张权重最大的薄片，而被舍弃的部分为

$$A-A_k=\sum _{i=k+1}^r{\sigma }_i{u}_i{v}_i^T.$$

由于保留部分和舍弃部分都由正交薄片组成，截断效果可以从两个互补角度衡量。一个看“留下了多少”，即能量保留比例

$${\eta }_k=\frac{{\sum }_{i=1}^k{\sigma }_i^2}{{\sum }_{i=1}^r{\sigma }_i^2}.$$

另一个看“丢掉了多少”，即 Frobenius 重构误差

$$\Vert A-A_k{\Vert }_F=\sqrt{\sum _{i=k+1}^r{\sigma }_i^2}.$$

这两个指标来自同一个正交分解：${\eta }_k$越接近$1$，尾部奇异值的平方和就越小，$A_k$也越接近$A$。

回到前面的例子，矩阵只有两张薄片：

$$A=5u_1{v}_1^T+3u_2{v}_2^T.$$

取$k=1$时，只保留第一张薄片，得到

$$\begin{array}{rl}{A}_1& =5u_1{v}_1^T\\ & =\left[\begin{array}{cc}5/2& 5/2\\ 5/2& 5/2\\ 0& 0\end{array}\right].\end{array}$$

第二张薄片正是被舍弃的全部，因此

$$\Vert A-A_1{\Vert }_F={\sigma }_2=3,$$

而保留的能量比例为

$${\eta }_1=\frac{{\sigma }_1^2}{{\sigma }_1^2+{\sigma }_2^2}=\frac{25}{34}\approx 73.53\%.$$

这个例子把截断的取舍具体化了：用秩$1$矩阵代替原来的秩$2$矩阵，换来更简单的表示，同时承担大小为$3$的 Frobenius 误差。

（二）压缩和去噪

截断 SVD 的价值来自同一件事：它用少数强模式代替原矩阵，同时舍弃较弱的尾部模式。

在压缩中，只需存储$U_k$、$k$个奇异值和$V_k$。它们一共包含

$$mk+k+nk=k(m+n+1)$$

个标量。只有当$k(m+n+1)<mn$时，参数数量才少于直接存储$A$所需的$mn$个标量；实际文件大小还会受到量化精度、编码方式和格式开销的影响。

在去噪中，如果噪声主要分散在小奇异值对应的弱模式里，截断就会在保留主体结构的同时滤去部分噪声。

下面的动图展示了随着$k$增大，主要轮廓和细节如何逐层恢复。

$k$的选择因此是一种取舍：较小的$k$带来更强的压缩或平滑效果，但也会增大重构误差；较大的$k$保留更多细节，却需要更多参数。能量比例${\eta }_k$可以提供定量参考，但不能把“小奇异值”直接等同于“噪声”，因为真实而细微的结构也可能落在尾部模式中。

（三）最优低秩近似

现在还剩最后一个问题：在所有秩不超过$k$的候选矩阵中，为什么偏偏保留前$k$张薄片？

Eckart–Young 定理给出了答案。在所有秩不超过$k$的矩阵$B$中，截断 SVD 得到的$A_k$满足

$$\Vert A-A_k{\Vert }_F=\underset{\mathrm{rank}(B)\le k}{\min }\Vert A-B{\Vert }_F=\sqrt{{\sigma }_{k+1}^2+\cdots +{\sigma }_r^2}.$$

因此，按奇异值从大到小截断并非经验规则，而是 Frobenius 范数下固定秩的最优选择。若尾部奇异值迅速衰减，少数薄片就足以逼近原矩阵；若奇异值下降缓慢，信息分散在许多方向上，任何低秩表示都会承担较大的误差。

回顾三篇文章，第一篇用几何图像解释矩阵“做了什么”，第二篇从极值问题说明奇异方向“为什么存在”，本篇则把这些方向算出来，并进一步说明怎样按奇异值取舍信息。从$A^{T}A$的特征向量，到$U\Sigma V^T$的构造，再到正交秩$1$薄片和最优截断，SVD 的几何直觉、计算方法与低秩近似由此连成了一条完整的主线。

---

系列文章目录

1、 深入理解主成分分析（PCA）
2、SVD 的核心思想及几何理解
3、SVD 公式推导
4、SVD的计算和低秩近似