我们知道一个矩阵可以通过SVD分解来得到一种特定的表示形式，而低秩近似是其中一个很重要的应用，这次我们通过简单的程序通过图像来实现这个过程。

SVD

首先就是先介绍什么是SVD，也就是奇异值分解，以下粘贴自维基百科：

假设 $M$ 是一个 $m \times n$ 阶矩阵，其中的元素全部属于域 $K$（实数域或复数域）。如此则存在一个分解使得：

$$M = U \Sigma V^\ast$$

其中 $U$ 是 $m \times m$ 阶酉矩阵；$\Sigma$ 是 $m \times n$ 阶非负实数对角矩阵；而 $V^\ast$ 即 $V$ 的共轭转置，是 $n \times n$ 阶酉矩阵。这样的分解就称作 $M$ 的奇异值分解。$\Sigma$ 对角线上的元素 $\sigma_i$ 即为 $M$ 的奇异值。

若将对角线元素相同但排列顺序不同的 $\Sigma$ 视为等价，$\Sigma$ 由 $M$ 唯一确定。（虽然 $U$ 和 $V$ 仍然不能确定。）

这里通过通俗一点的说法就是一个矩阵 $M \in \mathbb{R}^{m \times n}$ 可以看成是一个从 $\mathbb{R}^n$ 到 $\mathbb{R}^m$ 的映射过程，而SVD则将这个过程分成了几个矩阵映射并直观地展示出来，SVD还可以表示成：

$$M = \sum_{i=1}^r \sigma_i u_i v_i^\ast$$

其中 $v_i^\ast$ 表示第 $i$ 个右奇异向量的转置，$u_i$ 表示第 $i$ 个左奇异向量，$\sigma_i$ 表示第 $i$ 个奇异值。这里可以看出奇异值分解的本质是将所有不同输入方向和输出方向进行聚合后乘以一个系数（可以认为是能量的大小），然后再进行求和，而能量的大小则反映了矩阵在对应输出方向和输入方向变换的一个权重大小，我们可以把 $\sigma_i u_i v_i^\ast$ 看成是一个通道，而矩阵则是由多个通道加和组成的一个总变换。

奇异值的数量和矩阵的秩有关，矩阵的秩表示了矩阵有多少个线性无关的方向，也就是最多能对输入做几个方向的变换。

如果以矩阵乘法的角度上看的话，我们知道 $MX=B$ 本质上是将 $M$ 空间中的 $X$ 表示映射回标准坐标系中，而奇异值分解中我们会看到，矩阵的这一性质可以看成是先将输入 $X$ 从标准坐标系中映射到 $V$ 空间中，然后在 $V$ 空间中的不同方向进行不同程度的拉伸，最后再将其视为 $U$ 空间坐标，最后映射回标准坐标系中的过程。本质上是先聚合后分类的一个过程。

从数据的角度上看，对于一个矩阵来说大部分的能量，也就是权重都聚集在奇异值较大的方向上，因为 $U$ 和 $V$ 都是正交矩阵，所以 $\sum_{i=1}^r \sigma_i u_i v_i^\ast$ 的数值大小主要是通过奇异值来决定的，如果我们尝试将整个矩阵进行奇异值的截断会怎么样呢？

程序实现

我们这里选取一张图片

amakano.jpg

从主观上看，我们如果只保留前K个奇异值，则整张图片的能量会发生衰减，我们在这里可以先简单理解能量为数据在某个通道上数据的方差大小（并不完全相等），可以表示数据的离散程度，而人眼对整张图片的辨识是体现在图片的色彩丰富度上的，如果数据在某个通道上方差越大，则其中不同强度的像素点越多，所以整张图像对于人眼来说特征则越丰富，反之人眼则看不出区别，所以奇异值大的通道上往往会聚集比较大的辨识特征，而较小的通道则对辨识图片几乎没有贡献

然后我们可以写一个简单的代码看看变化

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
from PIL import Image
from torchvision.transforms import v2
import torch

img=Image.open("amakano.jpg").convert("RGB")


def low_rank_approximation(output_path:str,img,K:int,only_one_color=False):
    input=v2.functional.to_image(img).float()
    if not only_one_color:
        U,S,Vh=torch.linalg.svd(input,full_matrices=False)
        print(S.shape)
        S_diag = torch.diag_embed(S[:, :K])
        X_approx = U[:, :, :K] @ S_diag @ Vh[:, :K, :]
        X_approx = X_approx.clamp(0, 255).to(torch.uint8)
        img_out = v2.functional.to_pil_image(X_approx)
        img_out.save(output_path)
    else:
        for i in range(3):
            one_color_input=input[i,:,:]
            U,S,Vh=torch.linalg.svd(one_color_input,full_matrices=False)
            print(S.shape)
            S_diag = torch.diag_embed(S[:K])
            X_approx = U[:, :K] @ S_diag @ Vh[:K, :]
            X_approx = X_approx.clamp(0, 255).to(torch.uint8)
            img_out = v2.functional.to_pil_image(X_approx)
            img_out.save(str(i)+'_'+output_path)
    
low_rank_approximation("amakano_low_rank_approximated.jpg",img,60)

运行结果：

1
torch.Size([3, 1434])

也就是说整张图片有1434个奇异值

保留了前60个奇异值（从大到小）

可以看到如果只有前60个奇异值的情况下，我们仍然可以辨识出图像的大部分特征（也许）

SVD 奇异值渐进截断失真过程

如果我们观察整个变化过程的话会发现其实保留了前200个奇异值的情况下整张图片几乎没有任何变化，也就是说1434个奇异值中，对整张图片的特征识别起作用的其实只占了大概百分之14左右

那么为什么会产生这样的现象呢，笔者这里只能给出一个勉强地解释。

像类似于自然世界中的图像或者画师手绘的图像如果转化成矩阵本身就是近似低秩的，因为图像并非随机噪点，而是充满了强烈的几何规律与空间相关性（如大面积平涂色块与连续轮廓），这种结构上的高度有序性，必然导致奇异值分布极其不均匀，能量高度集中在少数几个低频主方向上，如果奇异值分布是均匀的，也就是都差不多大，那么整个图片在任何通道上携带的能量都是差不多的，通道与通道之间没有对比，最终整张图象会呈现出杂乱无章的现象。