在 Python 中用 QR 分解求解超定系统

Question

我正在尝试使用 QR 分解和 linalg.solve 求解一个超定系统，但我得到的错误是

LinAlgError: 数组的最后两个维度必须是正方形。

当 R 数组不是正方形时会发生这种情况，对吗？代码看起来像这样

import numpy as np
import math as ma

A = np.random.rand(2,3)
b = np.random.rand(2,1) 
Q, R = np.linalg.qr(A)
Qb = np.matmul(Q.T,b)
x_qr = np.linalg.solve(R,Qb)

有没有办法以更有效的方式为任意 A 维度编写此代码？如果没有，我该如何使这个代码片段起作用？

Answer 1

原因确实是矩阵R不是方阵，可能是系统多定的缘故。您可以尝试 np.linalg.lstsq，找到最小化平方误差的解决方案（如果存在，应该产生精确的解决方案）。

import numpy as np

A = np.random.rand(2, 3)
b = np.random.rand(2, 1) 
x_qr = np.linalg.lstsq(A, b)[0]

Answer 2

如numpy.linalg.solve的文档所示：

Computes the “exact” solution, x, of the well-determined, i.e., full rank, linear matrix equation ax = b.

您的方程组 underdetermined 没有多定。请注意，其中有 3 个变量和 2 个方程，因此方程比未知数少。

还要注意它还提到在 numpy.linalg.solve(a,b) 中，a 必须是一个 MxM 矩阵。这背后的原因是求解方程组Ax=b涉及计算A的逆，只有方阵是可逆的。

在这些情况下，一种常见的方法是采用 Moore-Penrose 伪逆，它将计算系统的最佳拟合（最小二乘）解。因此，与其尝试求解确切的解决方案，不如使用 numpy.linalg.lstsq:

x_qr = np.linalg.lstsq(R,Qb)

Answer 3

您需要使用标志模式='reduced'调用QR。默认的 Q R 矩阵返回为 M x M 和 M x N，因此如果 M 大于 N，则矩阵 R 将是非方阵。如果您选择简化（经济）模式，您的矩阵将为 M x N 和 N x N，在这种情况下，求解例程可以正常工作。

但是，对于超定系统，您也有 equations/unknowns 倒退。您的代码段应该是

import numpy as np 

A = np.random.rand(3,2)
b = np.random.rand(3,1) 
Q, R = np.linalg.qr(A, mode='reduced')
#print(Q.shape, R.shape)
Qb = np.matmul(Q.T,b)
x_qr = np.linalg.solve(R,Qb)

正如其他贡献者所指出的，您也可以直接调用 lstsq，但有时直接使用 Q 和 R 会更方便（例如，如果您还打算计算投影矩阵）。

在 Python 中用 QR 分解求解超定系统

Solve overdetermined system with QR decomposition in Python

python

numpy

least-squares

qr-decomposition