深入理解SVD

一、定义

    任何一个实矩阵  均可分解为三个子矩阵的乘积：
                                                                                             

    其中， 和  是正交矩阵，是  阶矩阵， 是  阶矩阵， 是  对角矩阵，对角线上的元素为矩阵  的奇异值，任何一个实矩阵的SVD分解必存在，但不唯一。

 

二、紧凑SVD与截断SVD

    设矩阵的秩为 ，，则称  为  的紧凑SVD。其中， 表示  的前  列， 表示  的前  行和前  列构成的对角矩阵， 表示  的前  行。

    如果只取最大的  （）个奇异值对应的部分，就得到矩阵的截断SVD，即 ， 表示  的前  列， 表示  的前  行和前  列构成的对角矩阵， 表示  的前  行。

    截断奇异值经常用于数据压缩，特别是对于图片的有损压缩。

    

三、几何解释

    这部分比较有意思，着重说明一下。线性变化

                                                                                           

    相当于将一个  维空间的向量  映射到一个  维空间中（矩阵乘法的公式）。再拆分一下即

                                                                                       

    上式相继对  做了如下变换： 首先， 在正交矩阵 的作用下，进行了一定的旋转，但是长度保持不变（正交变换不改变向量长度和向量之间的夹角）；

    然后， 对角线上的数值依次对旋转后的向量的各个维度进行缩放；最后，在正交矩阵  的作用下，被缩放的向量再次进行了旋转。

    说明： ，由于正交矩阵的逆是正交的，故  也是正交的。

 

四、重要性质

    一个关于SVD的重要性质如下：

                                                                  

    两边同时右乘 ，可得

                                                                                     

    展开即得

                                                                     

    其中， 是一个对角阵，主对角线上的元素是  的特征值，  的列向量是  的特征向量。同理有

                                                                                      

                                                                    

    其中， 是一个对角阵，主对角线上的元素是  的特征值，  的列向量是  的特征向量。

 

五、重新表示SVD

    关于公式 ，若令 ， 表示  位于第 k 行第 k 列的值， 则其等价形式表示如下：

                                                                                

    其中， 表示  的第  列， 表示  的第  行。可以看出上式是  个  的矩阵之和，这些矩阵之和还原了原始矩阵 。若  非常大，则对应的  矩阵  的值普遍较大，蕴含了原始矩阵中较多的信息；若  非常小，则可能是原始矩阵中的一些噪音。下面以图像去噪为例说明上式的作用。

---

        

---

    从左到右，分别是原图（左边）、保留90%主成分的图像（中间）和保留80%主成分的图像（右边）。可以看到 右图消除了原图中的大部分噪音，保留了主要图像信息。

    结合第四部分， 可看做是维度为  的  个样本的协方差矩阵，由图片的行表示特征， 可看做是维度为  的  个样本的协方差矩阵，由图片的列表示特征。无论是哪种情况，噪音点都不可能成为主成分，因为噪音点对应的维度方差几乎为零。所以，图像中的一些散点被剔除的特别彻底。下面是去噪代码：

import numpy as np
from numpy import linalg
from PIL import Image

img = Image.open("C:\\Users\\Paul\\Desktop\\test_img.png")
# 转为灰度图，并将灰度图转换为ndarray
gray_img = img.convert('L')
img_array = np.array(gray_img)

# m: 左奇异矩阵，sig: 奇异值对角阵，v: 右奇异矩阵
m,sig,v = np.linalg.svd(img_array)
sig = np.eye(sig.shape[0]) * sig

# 计算奇异值能量分布
sum_sig_list = []
s = 0.0
for i in range(sig.shape[0]):
    s += sig[i][i]
    sum_sig_list.append(s)

# 选择最优的奇异值个数, ratio控制阈值
i = 0
ratio = 0.7
while i < sig.shape[0]:
    r = sum_sig_list[i] / sum_sig_list[-1]
    if r >= ratio:
        print("%.2f : %.2f" % (sum_sig_list[i], sum_sig_list[-1]))
        i += 1
        break
    i += 1
print("%d, %.2f" % (i, sum_sig_list[i]))

sub_m = m[:,:i]       # 左奇异矩阵的前i列
sub_sig = sig[:i,:i]  # 对角阵的前i行和前i列组成的对角阵
sub_v = v[:i,:]       # 右奇异矩阵的前i行

# 压缩图片
compr_img_array = np.dot(np.dot(sub_m, sub_sig), sub_v)
compr_img = Image.fromarray(compr_img_array).convert("RGB")
compr_img.save("C:\\Users\\Paul\\Desktop\\test_img_70.png", "PNG")