svd降维 python案例_SVD（奇异值分解）Python实现

注：在《SVD(异值分解)小结》中分享了SVD原理，但其中只是利用了numpy.linalg.svd函数应用了它，并没有提到如何自己编写代码实现它，在这里，我再分享一下如何自已写一个SVD函数。但是这里会利用到SVD的原理，如何大家还不明白它的原理，可以去看看《SVD(异值分解)小结》，或者自行百度/google。

1、SVD算法实现

1.1 SVD原理简单回顾

有一个\(m \times n\)的实数矩阵\(A\)，我们可以将它分解成如下的形式

A = U\Sigma V^T

\tag{1-1}

其中\(U\)和\(V\)均为单位正交阵，即有\(UU^T=I\)和\(VV^T=I\)，\(U\)称为左奇异矩阵，\(V\)称为右奇异矩阵，\(\Sigma\)仅在主对角线上有值，我们称它为奇异值，其它元素均为0。上面矩阵的维度分别为\(U \in \mathbf{R}^{m\times m},\ \Sigma \in \mathbf{R}^{m\times n}\),\(\ V \in \mathbf{R}^{n\times n}\)。

正常求上面的\(U,V,\Sigma\)不便于求，我们可以利用如下性质

AA^T=U\Sigma V^TV\Sigma^TU^T=U\Sigma \Sigma^TU^T

\tag{1-2}

A^TA=V\Sigma^TU^TU\Sigma V^T=V\Sigma^T\Sigma V^T

\tag{1-3}

1.2 SVD算法

据1.1小节，对式(1-3)和式(1-4)做特征值分解，即可得到奇异值分解的结果。但是样分开求存在一定的问题，由于做特征值分解的时候，特征向量的正负号并不影响结果，比如，我们利用式(1-3)和(1-4)做特征值分解

AA^T\mathbf{u}_i = \sigma_i \mathbf{u}_i\quad \text{or} \quad AA^T(-\mathbf{u}_i) = \sigma_i (-\mathbf{u}_i)\\

A^TA\mathbf{v}_i = \sigma_i \mathbf{v}_i\quad \text{or} \quad A^TA(-\mathbf{v}_i) = \sigma_i (-\mathbf{v}_i)

如果在计算过程取，取上面的\(\mathbf{u}_i\)组成左奇异矩阵\(U\)，取\(-\mathbf{v}_i\)组成右奇异矩阵\(V\)，此时\(A\ne U\Sigma V^T\)。因此求\(\mathbf{v}_i\)时，要根据\(\mathbf{u}_i\)来求，这样才能保证\(A= U\Sigma V^T\)。因此，我们可以得出如下1.1计算SVD的算法。它主要是先做特性值分解，再根据特征值分解得到的左奇异矩阵\(U\)间接地求出部分的右奇异矩阵\(V'\in \mathbf{R}^{m\times n}\)。

算法1.1：SVD

输入：样本数据

输出：左奇异矩阵，奇异值矩阵，右奇异矩阵

计算特征值：特征值分解\(AA^T\)，其中\(A \in \mathbf{R}^{m\times n}\)为原始样本数据

AA^T=U\Sigma \Sigma^TU^T

得到左奇异矩阵\(U \in \mathbf{R}^{m \times m}\)和奇异值矩阵\(\Sigma' \in \mathbf{R}^{m \times m}\)

间接求部分右奇异矩阵：求\(V' \in \mathbf{R}^{m \times n}\)

利用\(A=U\Sigma'V'\)可得

V' = (U\Sigma')^{-1}A = (\Sigma')^{-1}U^TA

\tag{1-4}

返回\(U,\ \Sigma',\ V'\)，分别为左奇异矩阵，奇异值矩阵，右奇异矩阵。

注：这里得到的\(\Sigma'\)和\(V'\)与式(1-2)所得到的\(\Sigma,\ V\)有区别，它们的维度不一样。\(\Sigma'\)是只取了前\(m\)个奇异值形成的对角方阵，即\(\Sigma' \in \mathbf{R}^{m \times m}\)；\(V'\)不是一个方阵，它只取了\(V \in \mathbf{R}^{m \times n}\)的前\(m\)行(假设\(m < n\))，即有\(V' = V(:m,\cdot)\)。这样一来，我们同样有类似式(1-1)的数学关系成立，即

A = U\Sigma' (V')^T\tag{1-5}

我们可以利用此关系重建原始数据。

2、SVD的Python实现

以下代码的运行环境为python3.6+jupyter5.4。

2.1 SVD实现过程

读取数据

这里面的数据集大家随便找一个数据就好，如果有需要我的数据集，可以下在面留言。

import numpy as np

import pandas as pd

from scipy.io import loadmat

# 读取数据，使用自己数据集的路径。

train_data_mat = loadmat("../data/train_data2.mat")