东北网站建设公司,空投网站建设,网站建设业务员主动话术,网站建制作公司文章目录基本原理scipy实现测试基本原理
当AAA是方阵时#xff0c;可以很容易地进行特征分解#xff1a;AWΣW−1AW\Sigma W^{-1}AWΣW−1#xff0c;其中Σ\SigmaΣ是AAA的特征值组成的对角矩阵。如果WWW由标准正交基组成#xff0c;则W−1WTW^{-1}W^TW−1WT#xff0c;…
文章目录基本原理scipy实现测试基本原理
当AAA是方阵时可以很容易地进行特征分解AWΣW−1AW\Sigma W^{-1}AWΣW−1其中Σ\SigmaΣ是AAA的特征值组成的对角矩阵。如果WWW由标准正交基组成则W−1WTW^{-1}W^TW−1WT特征分解可进一步写成WTΣWW^T\Sigma WWTΣW。
然而当AAA不是方阵时情况大不一样了但仍然可以将AAA表示成AUΣVTAU\Sigma V^TAUΣVT的形式其中Σ\SigmaΣ也是对角矩阵对角线上的每个元素被称作奇异值。
奇异值的求解过程和特征值息息相关因为把AAA变成方阵很简单只要乘以转置就行。故令LAATLAA^TLAATRATARA^TARATA则L,RL, RL,R都可以求特征值λi\lambda_iλi和特征向量其中LLL的特征向量为AAA的左奇异向量RRR的特征向量为右奇异向量。对应的奇异值σiλi\sigma_i\sqrt{\lambda_i}σiλi。
scipy实现
scipy.sparse.linalg中实现了稀疏矩阵奇异值分解算法其参数列表如下
svds(A, k6, ncvNone, tol0, whichLM, v0None, maxiterNone, return_singular_vectorsTrue, solverarpack, random_stateNone, optionsNone)各参数含义如下
A 待分解矩阵k 奇异值个数必须在[k,kmax][k, k_{\max}][k,kmax]之间 当solverpropack时kmaxmin(M,N)k_{max}\min(M,N)kmaxmin(M,N)否则kmaxmin(M,N)−1k_{max}\min(M,N)-1kmaxmin(M,N)−1ncv solverarpack时此为Lanczos向量个数否则此项忽略。tol 奇异值容忍度为0表示达到机器的精度which 为LM时选取最大的奇异值SM则选取最小奇异值v0 迭代初值maxiter 迭代次数return_singular_vectors 可选4个值 True 返回奇异向量False 不返回奇异向量u: 如果M N只计算左奇异向量vh: 如果M N只计算右奇异向量如果 solverpropack这个选项将忽略矩阵维度 solver 可选arpack, propack, lobpcg但比较吊诡的是似乎并没有关于这三者区别的文档random_state 设置随机数状态optionsdict 求解器参数
其返回值有三
u 即UUUs 即奇异值数组也就是Σ\SigmaΣ的对角线vh 即VTV^TVT
测试
下面对奇异值分解做个测试
import numpy as np
from scipy.linalg import svd
from scipy.sparse import csc_array
from scipy.sparse.linalg import svds
np.random.seed(42) # 设置随机数状态
mat np.random.rand(500,800)
mat[mat0.9] 0
csc csc_array(mat)
u1, s1, vh1 svds(csc, k10)
u2, s2, vh2 svd(mat)结果是svds得到的结果和svd的前十个值完全相同只是排序不一样但也无关紧要。
下面测试一下二者的时间由于在Windows下用不了propack所以svds计算的奇异值数最多只能是M−1M-1M−1也就是499所以只能测试这个和svd返回500个奇异值的结果相比对结果如下 from timeit import timeittimeit(lambda : svds(csc, k499), number10)
3.651770199999987timeit(lambda : svd(mat), number10)
0.47201400000005833可见稀疏矩阵在计算上的确是比不上规整的矩阵。
