【机器学习】这次终于彻底理解了奇异值分解(SVD)原理及应用

思想年代 2022-02-09

展开全文

奇异值分解(Singular Value Decomposition，以下简称SVD)是在机器学习领域广泛应用的算法，有相当多的应用与奇异值都可以扯上关系，它不光可以用于降维算法中的特征分解，比如做feature reduction的PCA，做数据压缩（以图像压缩为代表）的算法还可以用于推荐系统。以及自然语言处理等领域，如做搜索引擎语义层次检索的LSI（Latent Semantic Indexing）。本文就对SVD的原理做一个总结。

特征值分解(EVD）

在开始学习奇异值分解（SVD）分解之前，先学习下特征值分解（EVD）。当然，理解SVD之前，还是需要一些基础知识储备，如果你对如下基础知识基本了解，可以跳过本节内容。

方阵

方块矩阵，也称方阵、方矩阵或正方矩阵，是行数及列数皆相同的矩阵。由矩阵组成的集合，连同矩阵加法和矩阵乘法，构成环。除了，此环并不是交换环。

单位矩阵

的对角线全是1而其他位置全是0，对所有矩阵及矩阵都有及。例如，若：

单位矩阵是方块矩阵环的单位元。

方块矩阵环的可逆元称为可逆矩阵或非奇异方阵。矩阵是可逆当且仅当存在矩阵使得

此时称为的逆矩阵，并记作。

特征向量与特征值

线性代数中对特征值和特征向量的定义：
设是阶方阵，如果存在和维非零向量，使，则称为方阵的一个特征值，为方阵对应于或属于特征值的一个特征向量。

从定义可以看出，对特征向量进行变换的实质是将特征向量进行缩放，缩放因子为特征值。

比如说下面的一个矩阵：

它其实对应的线性变换是下面的形式：

上面的矩阵是对称的，所以这个变换是一个对x，y轴的方向一个拉伸变换（每一个对角线上的元素将会对一个维度进行拉伸变换，当值>1时，是拉长，当值<1时时缩短）。

因此，特征向量的代数上含义是：
将矩阵乘法转换为数乘操作；

特征向量的几何含义是：
特征向量通过方阵变换只进行伸缩，而保持特征向量的方向不变。

特征值分解可以得到特征值与特征向量，特征值表示的是这个特征到底有多重要，而特征向量表示这个特征是什么，可以将每一个特征向量理解为一个线性的子空间，我们可以利用这些线性的子空间干很多的事情。

不过，特征值分解也有很多的局限，比如说变换的矩阵必须是方阵。

矩阵分解

现有特征值和特征向量：

其中是一个的实对称矩阵，是一个维向量，是矩阵的一个特征值，而是矩阵的特征值所对应的特征向量。

通过求解齐次方程来计算矩阵的特征值及特征向量：

求出特征值和特征向量就可以将矩阵特征分解。如果求出了矩阵的个特征值，以及这个特征值所对应的特征向量，，如果这个特征向量线性无关，那么矩阵就可以用下式的特征分解表示：

其中是这个特征向量所张成的维矩阵，而为这个特征值为主对角线的维矩阵。

一般我们会把的这个特征向量标准化，即满足，或者说，此时的个特征向量为标准正交基，满足，即也就是说为酉矩阵。

这样特征分解表达式可以写成：

矩阵分解了，相应的，其对应的映射也分解为三个映射。现在假设有向量，用将其变换到的列空间中，那么首先由先对做变换。

这个假设是向量原来的坐标，那么经过第一个变换之后，就可以把向量变成。

紧接着，在新的坐标系表示下，由中间那个对角矩阵对新的向量坐标换，其结果就是将向量往各个轴方向拉伸或压缩：

从上图可以看到，如果不是满秩的话，那么就是说对角阵的对角线上元素存在0，这时候就会导致维度退化，这样就可以降维了看没看到，这样就会使映射后的向量落入维空间的子空间中。

最后一个变换就是对拉伸或压缩后的向量做变换，由于和是互为逆矩阵，所以变换是变换的逆变换。

因此，从对称阵的分解对应的映射分解来分析一个矩阵的变换特点是非常直观的。假设对称阵特征值全为1那么显然它就是单位阵，如果对称阵的特征值有个别是0其他全是1，那么它就是一个正交投影矩阵，它将维向量投影到它的列空间中。

注意到要进行特征分解，矩阵必须为方阵。那么如果不是方阵，即行和列不相同时，此时需要借助SVD对矩阵进行分解。

奇异值分解(SVD)

本节进入我们的正题奇异值分解（SVD）。

奇异值

特征值及特征值分解都是针对方阵而言，现实世界中，我们看到的大部分矩阵不是方阵，比如每道数据有个点，一共采集了道数据，这样就形成了一个的矩阵，那么怎样才能像方阵一样提取出它的特征，以及特征的重要性。

奇异值分解就是来干这个事情的。奇异值相当于方阵中的特征值，奇异值分解相当于方阵中的特征值分解。

假设是一个的矩阵，那么得到的是一个的方阵（里面的向量是正交的，里面的向量称为左奇异向量），是一个的矩阵（除了对角线的元素都是0，对角线上的元素称为奇异值），是一个的矩阵，里面的向量也是正交的，里面的向量称为右奇异向量）。

那么奇异值和特征值是怎么对应起来的呢？首先，我们将一个矩阵的转置，将会得到一个方阵，我们用这个方阵求特征值可以得到：

这里得到的，就是我们上面的右奇异向量。此外我们还可以得到：

这里的就是上面说的奇异值，就是上面说的左奇异向量。奇异值跟特征值类似，在矩阵中也是从大到小排列，而且的减少特别的快，在很多情况下，前10%甚至1%的奇异值的和就占了全部的奇异值之和的99%以上了。也就是说，我们也可以用前大的奇异值来近似描述矩阵，这里定义一下部分奇异值分解：

右边的三个矩阵相乘的结果将会是一个接近于的矩阵，在这儿，越接近于，则相乘的结果越接近于。而这三个矩阵的面积之和（在存储观点来说，矩阵面积越小，存储量就越小）要远远小于原始的矩阵，我们如果想要压缩空间来表示原矩阵，我们存下这里的三个矩阵：就好了。

SVD 算法的原理推导

SVD也是对矩阵进行分解，但是和特征分解不同，SVD并不要求要分解的矩阵为方阵，是一种适用于任意矩阵的分解方法。假设我们的矩阵是一个的矩阵，那么我们定义矩阵的SVD为：

其中是一个的矩阵，是一个的矩阵，除了主对角线上的元素以外全为0，主对角线上的每个元素都称为奇异值，是一个的矩阵。和都是酉矩阵，即满足，，为单位矩阵。

下图可以很形象的看出上面SVD的定义：

那么我们如何求出SVD分解后的这三个矩阵呢？

如果我们将的转置和做矩阵乘法，那么会得到的一个方阵。既然是方阵，那么我们就可以进行特征分解，得到的特征值和特征向量满足下式：

这样我们就可以得到矩阵的个特征值和对应的个特征向量了。将的所有特征向量张成一个的矩阵，就是我们SVD公式里面的V矩阵了。一般我们将中的每个特征向量叫做的右奇异向量。

如果我们将和的转置做矩阵乘法，那么会得到的一个方阵。既然是方阵，那么我们就可以进行特征分解，得到的特征值和特征向量满足下式：

这样我们就可以得到矩阵的个特征值和对应的个特征向量了。将的所有特征向量张成一个的矩阵，就是我们SVD公式里面的矩阵了。一般我们将中的每个特征向量叫做的左奇异向量。

和我们都求出来了，现在就剩下奇异值矩阵没有求出了。由于除了对角线上是奇异值其他位置都是0，那我们只需要求出每个奇异值就可以了。

我们注意到:

这样我们可以求出我们的每个奇异值，进而求出奇异值矩阵。

上面还有一个问题没有讲，就是我们说的特征向量组成的就是我们SVD中的矩阵，而的特征向量组成的就是我们SVD中的矩阵，这有什么根据吗？这个其实很容易证明，我们以矩阵的证明为例。

上式证明使用了：。可以看出的特征向量组成的的确就是我们SVD中的V矩阵。类似的方法可以得到的特征向量组成的就是我们SVD中的矩阵。

进一步我们还可以看出我们的特征值矩阵等于奇异值矩阵的平方，也就是说特征值和奇异值满足如下关系：

这样也就是说，我们可以不用来计算奇异值，也可以通过求出的特征值取平方根来求奇异值。

这里我们用一个简单的例子来说明矩阵是如何进行奇异值分解的。我们的矩阵A定义为：

SVD计算举例

我们首先求出和

进而求出的特征值和特征向量：

接着求的特征值和特征向量：

利用求奇异值：

当然，我们也可以用直接求出奇异值为和 1.

最终得到的奇异值分解为：

特征分解和奇异值分解的区别

特征值分解和奇异值分解的目的都是一样，就是提取出一个矩阵最重要的特征。

所有的矩阵都可以进行奇异值分解，但只有方阵才可以进行特征值分解。当所给的矩阵是对称的方阵，，二者的结果是相同的。也就是说对称矩阵的特征值分解是所有奇异值分解的一个特例。但是二者还是存在一些小的差异，奇异值分解需要对奇异值进行从大到小的排序，而且全部是大于等于零。

对于如上的奇异值分解公式，我们可以看到奇异值分解同时包含了旋转，缩放和投影三种作用，并且和都起到了对进行旋转的作用，而起到了对缩放的作用，特征值分解只有缩放的效果。

在应用上，奇异值分解主要用于数据矩阵，而特征值分解主要用于方型的相关矩阵。

SVD的性质

上面几节我们对SVD的定义和计算做了详细的描述，似乎看不出我们费这么大的力气做SVD有什么好处。那么SVD有什么重要的性质值得我们注意呢？

对于奇异值，它跟我们特征分解中的特征值类似，在奇异值矩阵中也是按照从大到小排列，而且奇异值的减少特别的快，在很多情况下，前10%甚至1%的奇异值的和就占了全部的奇异值之和的99%以上的比例。也就是说，我们也可以用最大的k个的奇异值和对应的左右奇异向量来近似描述矩阵。也就是说：

其中要比小很多，也就是一个大的矩阵A可以用三个小的矩阵来表示。如下图所示，现在我们的矩阵只需要灰色的部分的三个小矩阵就可以近似描述了。

由于这个重要的性质，SVD可以用于PCA降维，来做数据压缩和去噪。也可以用于推荐算法，将用户和喜好对应的矩阵做特征分解，进而得到隐含的用户需求来做推荐。同时也可以用于NLP中的算法，比如潜在语义索引（LSI）。

SVD的应用

PCA降维中的应用

下面我们就对SVD用于PCA降维做一个介绍。

在主成分分析（PCA）原理总结中，我们讲到要用PCA降维，需要找到样本协方差矩阵的最大的个特征向量，然后用这最大的d个特征向量张成的矩阵来做低维投影降维。可以看出，在这个过程中需要先求出协方差矩阵，当样本数多样本特征数也多的时候，这个计算量是很大的。

注意到我们的SVD也可以得到协方差矩阵最大的个特征向量张成的矩阵，但是SVD有个好处，有一些SVD的实现算法可以不求先求出协方差矩阵，也能求出我们的右奇异矩阵。也就是说，我们的PCA算法可以不用做特征分解，而是做SVD来完成。这个方法在样本量很大的时候很有效。实际上，scikit-learn的PCA算法的背后真正的实现就是用的SVD，而不是我们我们认为的暴力特征分解。

另一方面，注意到PCA仅仅使用了我们SVD的右奇异矩阵，没有使用左奇异矩阵，那么左奇异矩阵有什么用呢？

假设我们的样本是的矩阵，如果我们通过SVD找到了矩阵最大的个特征向量张成的维矩阵，则我们如果进行如下处理：

可以得到一个的矩阵，这个矩阵和我们原来的维样本矩阵X相比，行数从减到了，可见对行数进行了压缩。也就是说，左奇异矩阵可以用于行数的压缩。相对的，右奇异矩阵可以用于列数即特征维度的压缩，也就是我们的PCA降维。

SVD作为一个很基本的算法，在很多机器学习算法中都有它的身影，特别是在现在的大数据时代，由于SVD可以实现并行化，因此更是大展身手。SVD的原理不难，只要有基本的线性代数知识就可以理解，实现也很简单因此值得仔细的研究。当然，SVD的缺点是分解出的矩阵解释性往往不强，有点黑盒子的味道，不过这不影响它的使用。

数据压缩的应用

奇异值分解可用于有效地表示数据。例如，假设我们希望传输以下图像，该图像由个黑色或白色像素的阵列组成。

由于该图像中只有三种类型的列，如下所示，应该可以用更紧凑的形式表示数据。我们将图像表示为一个矩阵，其中每个条目要么是 0，代表黑色像素，要么是 1，代表白色。因此，矩阵中有 375 个条目。

如果我们对执行奇异值分解，我们会发现只有三个非零奇异值。

因此，矩阵可以表示为

此时有三个向量，每个向量有 15 个条目，三个向量，每个向量有 25 个条目，以及三个奇异值。这时我们可以仅使用 123 个数字而不是出现在矩阵中的 375 个数字来表示矩阵。通过这种方式，奇异值分解删除了矩阵中的冗余，并提供了一种消除冗余后的格式。

为什么只有三个非零奇异值？请记住，非零奇异值的数量等于矩阵的秩。在这种情况下，我们看到矩阵中有三个线性独立的列，这意味着秩将为 3。

降噪的应用(noise reduction)

前面例子的奇异值都不为零，或者都还算比较大，下面我们来探索一下拥有零或者非常小的奇异值的情况。通常来讲，大的奇异值对应的部分会包含更多的信息。比如，我们有一张扫描，带有噪声的图像，如下图所示

我们可以以同样的方式进行：使用矩阵表示数据并执行奇异值分解。我们找到以下奇异值：

显然，前三个奇异值是最重要的，因此我们假设其他奇异值是由于图像中的噪声并进行近似

这导致以下改进的图像。

数据分析中的应用

每当我们收集数据时，噪音也会出现：无论仪器有多好，测量中总会有一些误差。如果我们还记得大奇异值指向矩阵中的重要特征这一主题，那么一旦收集到数据，使用奇异值分解来研究数据似乎很自然。

例如，假设我们收集一些数据，如下所示：

我们可以将数据放入矩阵中：

-1.03	0.74	-0.02	0.51	-1.31	0.99	0.69	-0.12	-0.72	1.11
-2.23	1.61	-0.02	0.88	-2.39	2.02	1.62	-0.35	-1.67	2.46

并执行奇异值分解。我们找到奇异值

由于一个奇异值比另一个大得多，可以安全地假设的小值是由于数据中的噪声，并且这个奇异值理想地为零。在这种情况下，矩阵的秩为 1 意味着所有数据都位于定义的线上。

这个简短的例子是属于主成分分析的领域，这是一组使用奇异值来检测数据中的依赖性和冗余性的技术。

以类似的方式，奇异值分解可用于检测数据中的分组，这解释了为什么奇异值分解被用于尝试改进 Netflix 的电影推荐系统。你观看过的电影的评分允许程序将你分类为与你的评分相似的其他人。可以通过选择你小组中其他人评价很高的电影来提出建议。

在图像压缩中的应用

一个图像矩阵，我们总可以将它分解为以下形式，通过选取不同个数的 Σ 中的奇异值，就可以实现图像的压缩：

如果只想实现图像压缩，我们可以使用python numpy 库中的 linalg.svd 对图像矩阵进行分解，进而提取前K个奇异值便能实现SVD图像压缩的效果，下面我们看一下代码:

#_*_ coding:utf-8_*_
import numpy as np
import cv2
 
img = cv2.imread('harden.jpg')
print('origin image shape is ', img.shape)
# 表示 RGB 中各有一个矩阵，都为300*532
#  origin image shape is  (300, 532, 3)
 
 
def svd_compression(img, k):
    res_image = np.zeros_like(img)
    for i in range(img.shape[2]):
        # 进行奇异值分解, 从svd函数中得到的奇异值sigma 是从大到小排列的
        U, Sigma, VT = np.linalg.svd(img[:,:,i])
        res_image[:, :, i] = U[:,:k].dot(np.diag(Sigma[:k])).dot(VT[:k,:])
 
    return res_image
 
 
# 保留前 k 个奇异值
res1 = svd_compression(img, k=300)
res2 = svd_compression(img, k=200)
res3 = svd_compression(img, k=100)
res4 = svd_compression(img, k=50)
 
row11 = np.hstack((res1, res2))
row22 = np.hstack((res3, res4))
res = np.vstack((row11, row22))
 
cv2.imshow('img', res)
cv2.waitKey(0)
cv2.destroyAllWindows()

我们这里分别提取了前300， 200， 100， 50 的奇异值，图如下：

可以看到，当我们取到前面300个奇异值来重构图片时，基本与原图看不出来差别，甚至两百的都是都还OK。

所以从上面的压缩结果来看，奇异值可以被看作是一个矩阵的代表值，或者说奇异值能够代表这个矩阵的信息，当奇异值越大时，它代表的信息越多。因此我们取前面若干个最大的奇异值，就可以基本还原出数据本身。

Numpy实现SVD

在Python中的Numpy模块中，已经实现了矩阵的奇异值分解。以下为示例的应用代码：

import numpy as np
#generate a random 3*4 matrix 
A =  np.random.randint(5, size=(3, 4))
#parameter full_matrices: control the size of P and Q
#d returns as numpy.ndarray, not matrix 
P,d,Q = np.linalg.svd(A, full_matrices=True)
print('A:',A)
print('P:',P)
#D return as diagonal 3*4 matrix
D = np.zeros(12).reshape(3,4)
for i in range(len(d)):
    D[i][i] = d[i]
print('D:',D)
print('Q:',Q)
#check if P*D*Q == A
print('P*D*Q:',np.dot(P,np.dot(D,Q)))

A: 
[[4 1 0 1]
 [2 1 4 4]
 [3 1 1 0]]
P: 
[[-0.47560601  0.6510767   0.59152181]
 [-0.79274458 -0.60868324  0.03256919]
 [-0.38125445  0.45343561 -0.80563093]]
D: 
[[7.07739778 0.         
  0.         0.        ]
 [0.         3.8600754  
  0.         0.        ]
 [0.         0.         
  1.00511622 0.        ]]
Q: 
[[-0.65443214 -0.23308073 
  -0.50191227 -0.51524366]
 [ 0.71170815  0.12845062 
  -0.51327944 -0.46207809]
 [ 0.01425989 -0.18061586 
  -0.67191651  0.71812448]
 [ 0.25492496 -0.94686415  
   0.18208926 -0.0728357 ]]
P*D*Q: 
[[4.00000000e+00 1.00000000e+00 
  1.11355742e-15 1.00000000e+00]
 [2.00000000e+00 1.00000000e+00 
  4.00000000e+00 4.00000000e+00]
 [3.00000000e+00 1.00000000e+00 
  1.00000000e+00 1.03872001e-15]]