scipy.sparse.linalg.lsqr¶
- scipy.sparse.linalg.lsqr(A, b, damp=0.0, atol=1e-06, btol=1e-06, conlim=100000000.0, iter_lim=None, show=False, calc_var=False, x0=None)[源代码]¶
找出大型稀疏线性方程组的最小二乘解。
函数求解
Ax = b
或min ||Ax - b||^2
或min ||Ax - b||^2 + d^2 ||x||^2
。矩阵A可以是正方形或矩形(超定或欠定),并且可以具有任何秩。
1. Unsymmetric equations -- solve A*x = b 2. Linear least squares -- solve A*x = b in the least-squares sense 3. Damped least squares -- solve ( A )*x = ( b ) ( damp*I ) ( 0 ) in the least-squares sense
- 参数
- A{稀疏矩阵,ndarray,LinearOperator}
m×n矩阵的表示形式。或者,
A
可以是线性运算符,它可以产生Ax
和A^T x
使用,例如,scipy.sparse.linalg.LinearOperator
。- b类似阵列,形状(m,)
右侧向量
b
。- damp浮动
阻尼系数。默认值为0。
- ATOL,BTOL浮动,可选
止动公差。
lsqr
继续迭代,直到取决于ATOL和BTOL的某个向后误差估计小于某个量。让我们r = b - Ax
是当前近似解的残差向量x
。如果Ax = b
似乎是一致的,lsqr
在下列情况下终止norm(r) <= atol * norm(A) * norm(x) + btol * norm(b)
。否则,lsqr
在下列情况下终止norm(A^H r) <= atol * norm(A) * norm(r)
。如果两个公差均为1.0E-6(默认值),则最终norm(r)
应该精确到6位左右。(决赛)x
通常会有较少的正确数字,具体取决于cond(A)
和Lambda的大小。)如果 atol 或 btol 为None,则将使用默认值1.0E-6。理想情况下,它们应该是对以下条目中的相对误差的估计A
和b
分别为。例如,如果的条目A
有7个正确的数字,设置atol = 1e-7
。这可以防止算法在输入数据的不确定性之外进行不必要的工作。- conlim浮动,可选
又一个停车容差。如果估计为
cond(A)
超过 conlim 。对于兼容的系统Ax = b
, conlim 可能高达1.0E+12(比方说)。对于最小二乘问题,余弦应小于1.0E+8,最大精度可通过设置atol = btol = conlim = zero
,但是迭代次数可能会过多。默认值为1e8。- iter_lim整型,可选
明确限制迭代次数(出于安全考虑)。
- show布尔值,可选
显示迭代日志。默认值为False。
- calc_var布尔值,可选
是否估计对角线
(A'A + damp^2*I)^{{-1}}
。- x0ARRAY_LIKE,Shape(n,),可选
如果没有使用零,则初始猜测x。默认值为None。
1.0.0 新版功能.
- 退货
- xNdarray of Floor(浮子线)
最终的解决方案。
- istop集成
给出终止的原因。1表示x是Ax=b的近似解。2表示x近似解决最小二乘问题。
- itn集成
终止时的迭代编号。
- r1norm浮动
norm(r)
,在哪里r = b - Ax
。- r2norm浮动
sqrt( norm(r)^2 + damp^2 * norm(x)^2 )
。等于 r1norm 如果damp == 0
。- anorm浮动
的Frobenius范数的估计
Abar = [[A]; [damp*I]]
。- acond浮动
估计
cond(Abar)
。- arnorm浮动
估计
norm(A'*r - damp^2*x)
。- xnorm浮动
norm(x)
- varNdarray of Floor(浮子线)
如果
calc_var
为真,则估计所有对角线(A'A)^{{-1}}
(如果damp == 0
)或更一般地说(A'A + damp^2*I)^{{-1}}
。如果A具有完整的列排名或damp > 0
。(不确定var是什么意思,如果rank(A) < n
和damp = 0.
)
注意事项
LSQR使用迭代方法近似求解。达到一定精度所需的迭代次数很大程度上取决于问题的规模。因此,在可能的情况下,应避免A的行或列的缩放不佳。
例如,在问题1中,解决方案不会因行缩放而改变。如果A的某一行与A的其他行相比非常小或很大,则(A,b)的相应行应该放大或缩小。
在问题1和问题2中,解决方案x在列缩放之后很容易恢复。除非知道更好的信息,否则A的非零列应该进行缩放,以使它们都具有相同的欧几里得范数(例如,1.0)。
在问题3中,如果阻尼不为零,则没有重新缩放的自由。但是,只有在注意到A的比例之后,才应指定阻尼值。
参数DAMP旨在通过防止真正的解非常大来帮助使病态系统正则化。另一种对正则化的帮助是由参数aSecond提供的,该参数可用于在计算的解变得非常大之前终止迭代。
如果一些初步估计
x0
是已知的,并且如果damp == 0
,可以按以下步骤进行:计算残差向量
r0 = b - A*x0
。使用LSQR求解系统
A*dx = r0
。添加校正DX以获得最终解
x = x0 + dx
。
这就要求
x0
在调用LSQR之前和之后均可用。为了判断好处,假设LSQR使用K1次迭代来求解A x = b and k2 iterations to solve A dx=R0。如果X0是“好的”,则范数(R0)将小于范数(B)。如果每个系统使用相同的停止公差ATOL和BTOL,则K1和K2将相似,但最终解X0+DX应该更精确。减少总工作量的唯一方法是对第二个系统使用更大的停止容差。如果某个值bto1适用于A x = b, the larger value btol 范数(B)/范数(R0)应适用于A*DX=R0。预处理是减少迭代次数的另一种方式。如果有可能解决一个相关的系统
M*x = b
有效地,在M以某种有用的方式逼近A的情况下(例如,M-A具有较低的秩或者其元素相对于A的元素较小),LSQR可能在系统上更快地收敛A*M(inverse)*z = b
之后,可以通过求解M*x=z来恢复x。如果A是对称的,则不应使用LSQR!
替代方法有对称共轭梯度法(CG)和/或SYMMLQ。SYMMLQ是对称CG的实现,它适用于任何对称A,并且收敛速度比LSQR更快。如果A是正定的,那么还有其他的对称CG实现,其每次迭代所需的工作量比SYMMLQ略少(但需要相同的迭代次数)。
参考文献
- 1
C.C.Paige和M.A.Saunders(1982a)。“LSQR:稀疏线性方程和稀疏最小二乘的算法”,“美国医学杂志”8(1),43-71。
- 2
C.C.Paige和M.A.Saunders(1982b)。“算法583.LSQR:稀疏线性方程和最小二乘问题”,ACM TOMS 8(2),195-209。
- 3
M.A.Saunders(1995)。“使用LSQR和Craig求解稀疏矩形系统”,第35位,588-604。
示例
>>> from scipy.sparse import csc_matrix >>> from scipy.sparse.linalg import lsqr >>> A = csc_matrix([[1., 0.], [1., 1.], [0., 1.]], dtype=float)
第一个示例具有平凡的解决方案 [0, 0]
>>> b = np.array([0., 0., 0.], dtype=float) >>> x, istop, itn, normr = lsqr(A, b)[:4] >>> istop 0 >>> x array([ 0., 0.])
停止码 istop=0 RETURN表示找到了一个零向量作为解决方案。退回的解决方案 x 确实包含了 [0., 0.] 。下一个示例有一个非常重要的解决方案:
>>> b = np.array([1., 0., -1.], dtype=float) >>> x, istop, itn, r1norm = lsqr(A, b)[:4] >>> istop 1 >>> x array([ 1., -1.]) >>> itn 1 >>> r1norm 4.440892098500627e-16
如图所示 istop=1 ,
lsqr
找到了符合容差限制的解决方案。给定的解决方案 [1., -1.] 显然是解了这个方程式。其余返回值包括有关迭代次数的信息 (itn=1 )和所解方程的左侧和右侧的剩余差。最后一个示例演示了方程没有解的情况下的行为:>>> b = np.array([1., 0.01, -1.], dtype=float) >>> x, istop, itn, r1norm = lsqr(A, b)[:4] >>> istop 2 >>> x array([ 1.00333333, -0.99666667]) >>> A.dot(x)-b array([ 0.00333333, -0.00333333, 0.00333333]) >>> r1norm 0.005773502691896255
istop 指示系统不一致,因此 x 而是相应的最小二乘问题的近似解。 r1norm 包含找到的最小残差的范数。