SVD(Singular Value Decomposition)奇异值分解的理解

  **SVD(Singular Value Decomposition)奇异值分解**
在很多情况下，数据中的一小段携带了数据集中的大部分信息，而其他信息要么是噪声，要么就是毫不相关的信息。矩阵分解可将原始矩阵表示成新的易于处理的形式，新形式是两个或多个矩阵的乘积。不同的矩阵分解技术具有不同的性质，其中有些更适合于某个应用，有些则更适合于其他应用。最常见的一种矩阵分解技术就是SVD。

1、简介
奇异值分解（singular Value Decomposition），简称SVD，可以用来简化数据，去除噪声和冗余，提高算法的结果。它是线性代数中矩阵分解的方法。假如有一个矩阵A，对它进行奇异值分解，可以得到三个矩阵：A=UΣVT，这三个矩阵的大小用图形来表示

图1 矩阵的分解图
由公式可以看出如果原始矩阵A是m行n列，那么分解成的三个矩阵U，Σ，V依次是m行m列，m行n列，n行n列。矩阵sigma(即图1 U和V中间的矩阵)除了对角元素，其他元素都为0，并且对角元素是从大到小排列的，前面的元素比较大，后面的很多元素接近0。这些对角元素就是奇异值（Singular Value），他们对应了原始数据集矩阵A的奇异值。奇异值就是矩阵AAT（A的转置）特征值的平方根。sigma中有n个奇异值，但是由于排在后面的很多接近0，意味着原始数据集中仅有r个重要特征，其余特征则都是噪声或冗余特征，所以我们可以仅保留比较大的r个奇异值，于是矩阵的大小又可表示如下：

图2 去除不重要特征的矩阵分解
当然在实际应用中，我们仅需保留着三个比较小的矩阵，就能表示A，不仅节省存储量，在计算的时候更是减少了计算量。从生物信息学到金融学等在内的很多应用中，SVD都是提取信息的强大工具。
SVD的应用领域包括但不限于以下：
 信息检索（隐性语义索引）（Latent Semantic Indexing，LSI）
 图像压缩
 推荐系统
 金融领域
2、利用Python实现SVD
Python3 numpy包里有linalg的线性代数工具箱，其有包括处理求范式、逆、行列式、伪逆、秩、qr分解、svd分解等等函数。下面用一个矩阵分解的例子来说明。
案例基于Python3 需要安装numpy包，这里是用pip install numpy来安装的，当然也可以用其它方法。具体代码及注释如下。
from numpy import
#import numpy as np
from numpy import linalg as la
def loadExData() :
return [[1, 1, 1, 0, 0],
[2, 2, 2, 0, 0],
[1, 1, 1, 0, 0],
[5, 5, 5, 0, 0],
[1, 1, 0, 2, 2],
[0, 0, 0, 3, 3],
[0, 0, 0, 1, 1]]
data=loadExData()
u,sigma,vt=la.svd(data)#利用线代库特征提取,sigma为奇异值
print (sigma)#输出sigma的目的是为了看下可以取哪些，舍弃哪些，即取出前面重要的数据
#通过测试发现前3个数据比其他的值大很多，后两个值在不同机器上结果可能会稍有差异，但数量级差不多，于是保留前3
‘----------由于sigma是以向量的形式存储，故需要将其转为矩阵，转为矩阵的方法有以下两种---------’
sig3=mat([[sigma[0],0,0],
[0,sigma[1],0],
[0,0,sigma[2]]])
#也可以利用numpy的diag函数，直接sig3=diag(sigma)[:3,:3]即可，diag将行向量转为矩阵，值放在对角线上,并取前面3行3列
‘---------下面是重构原始矩阵，只要前3行3列-----------’
#重构原始矩阵,取u的3列，sigma的3行3列，vt的3行相乘得到的结果近似原始矩阵。Data m×n≈Um×3Σ3×3VT3×n
print (u[:,:3]sig3vt[:3,:])

‘------下面是如何取出前3个奇异值的方法，先算出奇异值的平方再求和为总能量，看前几维的数据超过总能量的90%-----’
‘’’
sig2=sigma**2
print (sig2)
print(sum(sig2))
print(sum(sig2)*0.9)
print(sum(sig2[:3]))

‘’’

运行结果如下：

通过上方截图发现sigma只有一行，由于矩阵除了对角元素其他均为0，因此这种仅返回对角元素的方式能够节省空间，这就是由Numpy的内部机制产生的。以后一旦看到Sigma就知道其为一个矩阵。
同时也发现Simga中前三个数值比其他的值大多了，所以我们的原始数据集A可以用如下结果近似：A m×n≈Um×3Σ3×3VT3×n
通过运行结果可以发现得到的结果和原始数据矩阵没有什么差别。
问题：如何知道保留前三个奇异值的呢？
策略1：保留矩阵中90%的能量信息。为了计算总能量信息，将所有的奇异值求其平方和，将奇异值的平方和累加到总值的90%为止。
策略2：当数据有上万的奇异值时，那么就保留前面的2000或3000个。
由于该例数据较少，可以肉眼看出舍弃哪些奇异值。
3、用Python3实现基于SVD的图像压缩
通过可视化的方式，可以让我们很容易就能看到SVD对数据近似的效果。用一个手写的数字0图像来进行举例说明，首先在txt文档中存有一个图像大小是32x32=1024像素的0,1数字，试图对图像进行压缩，希望可以节省空间。方法是可以使用SVD来对数据降维，从而实现图像的压缩。下面就会看到利用SVD的手写数字图像的压缩过程了。下面的程序包含了数字的读入和压缩代码。要了解最后的压缩效果，对压缩后的图像进行了重构。
代码如下：
#coding=utf-8
from numpy import *

图像压缩函数

用于打印矩阵，

由于矩阵含有浮点数，因此必须定义浅色和深色。这里通过一个阈值来界定。

该函数遍历所有的矩阵元素，当元素大于阈值时打印1，否则打印0

#printMat()函数用于打印矩阵
def printMat(inMat, thresh=0.8):
for i in range(32):
for k in range(32):
if float(inMat[i,k]) > thresh:
print(1,end=""),
else: print(0,end=""),
print("")
#imgCompress()实现图像压缩 ,它允许基于任意给定的奇异值数目来重构图像
def imgCompress(numSV=3, thresh=0.8):
# 构建一个列表myl
myl = []
# 打开文本文件，以数值方式读入字符
for line in open(‘0_5.txt’,“r”).readlines():
newRow = []
for i in range(32):
newRow.append(int(line[i]))
myl.append(newRow)
myMat = mat(myl)
print (“original matrix**”)
#调用函数printMat
printMat(myMat, thresh)
# 对原始图像进行SVD分解并重构图像，通过将Sigma重构成SigRecon来实现
U,Sigma,VT =linalg.svd(myMat)
# Sigma是一个对角矩阵，需要建立一个全0矩阵，然后将前面的那些奇异值填充到对角线上。
SigRecon = mat(zeros((numSV, numSV)))
for k in range(numSV):#construct diagonal matrix from vector
SigRecon[k,k] = Sigma[k]
#通过截断的U和VT矩阵，用SigRecon得到重构后的矩阵
reconMat = U[:,:numSV]SigReconVT[:numSV,:]
print (“reconstructed matrix using %d singular values**” % numSV)
printMat(reconMat, thresh)
def main():
imgCompress(2)
if name == ‘main’:
main()

运行的结果如下：

数字0的原始矩阵 数字0重构后的矩阵
可知，只需要两个奇异值就能相当精确地对图像实现重构，U和VT都是32×2的矩阵，有两个奇异值，也即只需要32×2+2+32×2=130个0、1进行存储。原来需要1024个0、1，几乎获得10倍的压缩比。