30分钟学会SVD矩阵分解

SVD(Singular Value Decomposition)奇异值分解分解是机器学习中最重要的矩阵分解方法。

它能够将一个任意形状的矩阵分解成一个正交矩阵和一个对角矩阵以及另一个正交矩阵的乘积。

SVD分解具有非常深刻的几何含义。矩阵实际上对应着一种线性变换，一个矩阵作用到一个向量上，会得到一个新的向量。任何一个矩阵的操作效果可以分解成一次旋转，一次拉伸和维度改变，以及另外一次旋转三者作用效果的合成。

SVD分解通常用于数据压缩和数据降维。用于数据降维时，既可以对列降维，也可以对行降维，其中对列的降维等价于PCA的降维。

不仅如此，SVD算法还可以用于在声音和图像处理中剥离背景信号，在推荐算法中也经常出现它的身影。

一，SVD矩阵分解简介

SVD分解将任意矩阵分解成一个正交矩阵和一个对角矩阵以及另一个正交矩阵的乘积。

对角矩阵的对角元称为矩阵的奇异值，可以证明，奇异值总是大于等于0的。

当对角矩阵的奇异值按从大到小排列时，SVD分解是唯一的。

SVD分解有着非常深刻的几何含义。

矩阵实际上是对应着一种线性变换。一个矩阵作用到一个向量上，会得到一个新的向量。任何一个矩阵的操作效果可以分解成一次旋转，一次拉伸和维度改变，以及另外一次旋转三者作用效果的合成。

注意正交矩阵和作用到向量后是不会改变向量长度的，所以对应着旋转变换。

二，SVD分解的数学推演

三，SVD分解和数据压缩

假设 m = 10000,n = 8000,原来存储矩阵A需要存储8000万个数字，如果经过奇异值分解发现前100个奇异值贡献了99%的奇异值之和，于是可以近似只保留这100个奇异值及对应的左右奇异向量，那么只需要保留100+10000×100+100×8000= 180.01万个数字，只有原来的不到2.3%。

# 下面的范例示范SVD分解用于图片数据压缩。
%matplotlib inline 
%config InlineBackend.figure_format = 'svg'
import numpy as np 
from matplotlib import pyplot as plt
from skimage import data

def compressBySVD(img,r):
    u,s,vt = np.linalg.svd(img)
    ur = u[:,0:r]
    sr = s[0:r]
    vtr = vt[0:r,:]
    return (ur,sr,vtr)

def rebuildFromSVD(ur,sr,vtr):
    img = ur@np.diag(sr)@vtr
    return(img)


img = data.camera()/255.0

plt.figure(figsize=(10,8)) 
for i,r in enumerate([5,10,20,30,40,50,100,200],start = 1):
    ur,sr,vtr = compressBySVD(img,r)
    compress_ratio = (np.product(ur.shape) + len(sr) + 
                      np.product(vtr.shape))/np.product(img.shape)
    img_rebuild = rebuildFromSVD(ur,sr,vtr)
    ax=plt.subplot(3,3,i)
    ax.imshow(img_rebuild,cmap = "gray")
    ax.set_title("r=%d"%r+", compress_ratio=%.2f"%compress_ratio)
    ax.set_xticks([])
    ax.set_yticks([]) 

ax = plt.subplot(3,3,9)
ax.imshow(img,cmap = "gray")
ax.set_title("r = 512, original image")
ax.set_xticks([])
ax.set_yticks([]) 

plt.show()

四，SVD分解和PCA降维

PCA降维可以看成是SVD分解的一个应用。PCA降维使用的变换矩阵恰好是SVD分解的右奇异矩阵。

实际上，由于SVD分解存在着无需通过计算特征值和特征向量的可并行的数值迭代计算算法，sklearn的PCA降维算法正是通过SVD分解计算的。

# 演示SVD用于PCA降维的计算

%matplotlib inline 
%config InlineBackend.figure_format = 'svg'
import numpy as np 
from sklearn.decomposition import PCA

from matplotlib import pyplot as plt
from skimage import data

X = np.array([[-1.0, -3, -2], [-2, -1, -3], [-3, -2, -5], [2, 1, 3], [6, 1, 3], [2, 2, 3]])

pca = PCA(n_components= 2)
X_new = pca.fit_transform(X)
print("ndecomposition by pca:")
print("singular value:")
print(pca.singular_values_)
print("X_new:")
print(X_new)

print("ndecomposition by svd:")
U,S,Vt = np.linalg.svd(X-X.mean(axis = 0))
print("singular value:n",S[:2])
print("X_new:")
print(np.dot(X-X.mean(axis = 0),np.transpose(Vt)[:,0:2]))

# 注：降维结果中正负号的差异是因为PCA调整了SVD分解后的U和Vt符号以保持各列最大值取正

输出如下：

decomposition by pca:
singular value:
[11.31375337  2.89544001]
X_new:
[[ 3.23378083  1.06346839]
 [ 3.88607412 -0.50763321]
 [ 6.25267378  0.08479886]
 [-3.50509914 -0.96584476]
 [-6.02398361  1.89494314]
 [-3.84344598 -1.56973242]]

decomposition by svd:
singular value:
 [11.31375337  2.89544001]
X_new:
[[-3.23378083 -1.06346839]
 [-3.88607412  0.50763321]
 [-6.25267378 -0.08479886]
 [ 3.50509914  0.96584476]
 [ 6.02398361 -1.89494314]
 [ 3.84344598  1.56973242]]