Python奇异值分解

当$A$是方阵时，可以很容易地进行特征分解：$A=W\Sigma W^{-1}$，其中$\Sigma$是$A$的特征值组成的对角矩阵。如果$W$由标准正交基组成，则$W^{-1}=W^T$，特征分解可进一步写成$W^T\Sigma W$。

然而，当$A$不是方阵时，情况大不一样了，但仍然可以将$A$表示成$A=U\Sigma V^T$的形式，其中$\Sigma$也是对角矩阵，对角线上的每个元素被称作奇异值。

奇异值的求解过程和特征值息息相关，因为把$A$变成方阵很简单，只要乘以转置就行。故令$L=A{A}^T$，$R=A^{T}A$，则$L,R$都可以求特征值${\lambda }_i$和特征向量，其中$L$的特征向量为$A$的左奇异向量，$R$的特征向量为右奇异向量。对应的奇异值${\sigma }_i=\sqrt{{\lambda }_i}$。

numpy.linalg中提供了奇异值分解函数svd，参数为

svd(a, full_matrices=True, compute_uv=True, hermitian=False)

其中

• a 待分解矩阵，维度为(M, N)
• full_matrices 若为True，则U, Vh分别为(M, M)和(N, N)；否则分别为(M, K), (K, N)，K为M, N中较小的那个
• compute_uv 如果为False则不计算U, Vh
• hermitian 为True时，表示处理的是实对称矩阵

scipy.linalg中也提供了奇异值分解函数svd，其参数为

svd(a, full_matrices=True, compute_uv=True, overwrite_a=False, check_finite=True, lapack_driver='gesdd')

其中与numpy.linalg相同的参数，其意义也相同，不相同的部分，各参数含义如下

• overwrite_a 如果为True，则直接对a进行修改
• check_finite 如果为True，则进行有限性检查
• lapack_driver SVD分解的方法，有两个选择
  • 'gesdd' 效率更高
  • 'gesvd' 此为Matlab和R中使用的方法

其返回值即$U,\Sigma ,V^T$。

scipy.linalg还提供了两个和SVD相关的函数，svdvals(a)用于求a的奇异值；diagsvd(s, M, N)通过s, M, N，创建一个$\Sigma$矩阵。

下面测试一下svd

import numpy as np
import scipy.linalg as sl
a = np.random.rand(5,5)
u1, s1, vh1 = sl.svd(a)
u2, s2, vh2 = np.linalg.svd(a)
print(s1)
# [2.63698545 0.94063722 0.36159198 0.21052102 0.19014115]
print(s1-s2)
# [ 0.0  0.0  1.11022302e-16 -2.77555756e-17 0.0]

numpy和scipy的结果是几乎相同的，下面测试一下不同方法进行奇异值分解的时间

from timeit import timeit
a = np.random.rand(1000,1000)
timeit(lambda:sl.svd(a), number=10)
# 1.870287900000001
timeit(lambda:np.linalg.svd(a), number=10)
# 13.355788999999998
timeit(lambda:sl.svd(a, lapack_driver='gesvd'), number=10)
# 3.873418600000001