今天小编给大家分享一下Python如何实现两种稀疏矩阵的最小二乘法的相关知识点,内容详细,逻辑清晰,相信大部分人都还太了解这方面的知识,所以分享这篇文章给大家参考一下,希望大家阅读完这篇文章后有所收获,下面我们一起来了解一下吧。
最小二乘法
scipy.sparse.linalg实现了两种稀疏矩阵最小二乘法lsqr和lsmr,前者是经典算法,后者来自斯坦福优化实验室,据称可以比lsqr更快收敛。
这两个函数可以求解Ax=b,或arg minx ∥Ax−b∥2,或arg minx ∥Ax−b∥2 +d2∥x−x0∥2,其中A必须是方阵或三角阵,可以有任意秩。
通过设置容忍度at ,bt,可以控制算法精度,记r=b-Ax 为残差向量,如果Ax=b是相容的,lsqr在∥r∥⩽at∗∥A∥⋅∥x∥+bt∥b∥时终止;否则将在∥ATr∥⩽at∥A∥⋅∥r∥。
如果两个容忍度都是10−6 ,最终的∥r∥将有6位精度。
lsmr
的参数如下
lsmr(A, b, damp=0.0, atol=1e-06, btol=1e-06, conlim=100000000.0, maxiter=None, show=False, x0=None)
参数解释:
A 可谓稀疏矩阵、数组以及线性算子
b 为数组
damp 阻尼系数,默认为0
atol, btol 截止容忍度,是lsqr迭代的停止条件,即at ,bt 。
conlim 另一个截止条件,对于最小二乘问题,conlim应该小于108,如果Ax=b是相容的,则conlim最大可以设到1012
iter_limint 迭代次数
show 如果为True,则打印运算过程
calc_var 是否估计(A.T@A + damp**2*I)^{-1}的对角线
x0 阻尼系数相关
lsqr和lsmr相比,没有maxiter参数,但多了iter_lim, calc_va参数。
上述参数中,damp为阻尼系数,当其不为0时,记作δ,待解决的最小二乘问题变为
返回值
lsmr的返回值依次为:
x 即Ax=b中的x
istop 程序结束运行的原因
itn 迭代次数
normr ∥b−Ax∥
normar ∥AT (b−Ax)∥
norma ∥A∥
conda A的条件数
normx ∥x∥
lsqr的返回值为
x 即Ax=b中的x
istop 程序结束运行的原因
itn 迭代次数
r1norm
anorm 估计的Frobenius范数Aˉ
acond Aˉ的条件数
arnorm ∥ATr−δ2(x−x0)∥
xnorm ∥x∥
var (ATA)−1
二者的返回值较多,而且除了前四个之外,剩下的意义不同,调用时且须注意。
测试
下面对这两种算法进行验证,第一步就得先有一个稀疏矩阵
import numpy as npfrom scipy.sparse import csr_arraynp.random.seed(42) # 设置随机数状态mat = np.random.rand(500,500)mat[mat<0.9] = 0csr = csr_array(mat)
然后用这个稀疏矩阵乘以一个x,得到b
xs = np.arange(500)b = mat @ xs
接下来对这两个最小二乘函数进行测试
from scipy.sparse.linalg import lsmr, lsqrimport matplotlib.pyplot as pltmx = lsmr(csr, b)[0]qx = lsqr(csr, b)[0]plt.plot(xs, lw=0.5)plt.plot(mx, lw=0, marker='*', label="lsmr")plt.plot(qx, lw=0, marker='.', label="lsqr")plt.legend()plt.show()
为了对比清晰,对图像进行放大,可以说二者不分胜负
接下来比较二者的效率,500 × 500 500\times500500×500这个尺寸显然已经不合适了,用2000×2000
from timeit import timeitnp.random.seed(42) # 设置随机数状态mat = np.random.rand(500,500)mat[mat<0.9] = 0csr = csr_array(mat)timeit(lambda : lsmr(csr, b), number=10)timeit(lambda : lsqr(csr, b), number=10)
测试结果如下
>>> timeit(lambda : lsqr(csr, b), number=10)
0.5240591000001587
>>> timeit(lambda : lsmr(csr, b), number=10)
0.6156221000019286
看来lsmr并没有更快,看来斯坦福也不靠谱(滑稽)。
以上就是“Python如何实现两种稀疏矩阵的最小二乘法”这篇文章的所有内容,感谢各位的阅读!相信大家阅读完这篇文章都有很大的收获,小编每天都会为大家更新不同的知识,如果还想学习更多的知识,请关注编程网行业资讯频道。