3.3.4 奇异值分解-SVD

奇异值分解是线性代数中一种重要的矩阵分解，在信号处理、统计学等领域都有重要应用。假设M是一个m×n阶矩阵，存在一个分解使得：M =UΣV*。其中U是m×m阶酉矩阵；Σ是半正定m×n阶对角矩阵；而V*，即V的共轭转置，是n×n阶酉矩阵。这样的分解就称作M的奇异值分解。Σ对角线上的元素为M的奇异值。通常奇异值按照从大到小的顺序排列。

奇异值的数学描述读起来有些难懂，让我们通过一个实例说明奇异值分解的用途。下面的例子对一幅灰度图像进行奇异值分解，然后从三个分解矩阵中提取奇异值较大的部分数据还原原始图像。首先读入一幅图像，并通过其红绿蓝三个通道计算灰度图像img，图像的宽为375个像素、高为505个像素。

    r, g, b = np.rollaxis(pl.imread("vinci_target.png"), 2).astype(float)
    img = 0.2989 *
 r + 0.5870 *
 g + 0.1140 *
 b
    img.shape
    (505, 375)

调用scipy.linalg.svd( )对图像矩阵进行奇异值分解，得到三个分解部分：

●U：对应于公式中的U。

●s：对应于公式中的Σ，由于它是一个对角矩阵，只有对角线上的元素非零，因此s是一个一维数组，保存对角线上的非零元素。

●Vh：对应于公式中的V *。

下面的程序查看这三个数组的形状：

    U, s, Vh = linalg.svd(img)
    U.shape       s.shape      Vh.shape 
    ----------     -------     ----------
    (505, 505)     (375,)      (375, 375)

s中的每个值与Vh的行向量以及U中的列向量对应，默认按照从大到小的顺序排列，它表示与其对应的向量的重要性。由图3-9可知s中的奇异值大小差别很大，注意Y轴是对数坐标系。

图3-9 按从大到小排列的奇异值

    pl.semilogy(s, lw=3)

下面的composite( )选择U和Vh中的前n个向量重新合成矩阵，当用上所有向量时，重新合成的矩阵和原始矩阵相同：

    def composite(U, s, Vh, n):
        return np.dot(U[:, :n], s[:n, np.newaxis] *
 Vh[:n, :])
    
    print np.allclose(img, composite(U, s, Vh, len(s)))
    True

下面演示选择前10个、20个以及50个向量合成时的效果，如图3-10所示，可以看到使用的向量越多，结果就越接近原始图像：

图3-10 原始图像以及使用10个、20个、50个向量合成的图像（从左到右）

    img10 = composite(U, s, Vh, 10)
    img20 = composite(U, s, Vh, 20)
    img50 = composite(U, s, Vh, 50)
    
    %array_image img; img10; img20; img50