2

我尝试使用以下方法在 scipy 中解决线性最小二乘问题 Ax = b：

x = numpy.linalg.inv(A.T.dot(A)).dot(A.T).dot(b) #Usually not recommended

和

x = numpy.linalg.lstsq(A, b)

两者都给出几乎相同的结果。我还尝试使用 QR 算法手动执行此操作，即：

Qmat, Rmat = la.qr(A)

bpr = dot(Qmat.T,b)
n=len(bpr)
x = np.zeros(n)
for i in xrange(n-1, -1,-1):
    x[i] = bpr[i]
    for j in xrange(i+1, n):
        x[i] -= Rmat[i, j]*x[j]
    x[i] /= Rmat[i,i]

然而，这种方法给出了非常不准确的结果（大约 1e-2 的错误）。我在代码或数学上犯了 n00b 错误吗？或者，是方法的问题，还是 scipy 本身的问题？

我的 numpy 版本是 1.6.1（来自http://www.lfd.uci.edu/~gohlke/pythonlibs/的 mkl 编译版本），在 x86_64 上使用 Python 2.7.3。

4

2 回答 2

2

如果您使用这些二进制文件，则 QR 分解由英特尔 MKL 计算，并且可能是正确的。

对我来说，上面的代码1e-12为随机矩阵提供了正确结果内的解决方案。你用什么矩阵来测试它，你如何测量误差？

在某些情况下，最小二乘问题是病态的。例如，对于具有大空空间的矩阵，舍入误差会影响结果。考虑 rank-1 矩阵：

np.random.seed(1234)
v = np.random.rand(100)
A = v[:,None] * v[None,:]
b = np.random.randn(100)

x = scipy.linalg.lstsq(A, b)[0]
print(np.linalg.norm(A.dot(x) - b))
# -> 9.63612833771

# xp obtained using your above code
print(np.linalg.norm(A.dot(xp) - b))
# -> 3262.61161684

与 lstsq 中使用的更仔细编写的 LAPACK 例程相比，您自制的三角求解更容易出现舍入误差，因此它的准确度会有所降低。

于 2013-04-18T08:32:39.163 回答
0

您也可以尝试使用截断的特征分解。这意味着使用前 k 个特征值。我使用以下代码对最小二乘回归进行正则化，y=kc

在此处输入图像描述

其中u是特征向量，lambda是特征值

线性模型的核矩阵：

k=np.dot(X,X.T)

然后特征分解：

 w, v = scipy.linalg.eigh(k, eigvals=(lo, hi))

然后

temp= np.dot(np.dot(v,np.diag(1.0/w)),v.T)
c=np.dot(temp,y)

同样对于正则化，您应该0.001在内核矩阵的对角线上添加一个小值（如 ），否则您将拥有负特征值，这会阻止您的内核矩阵不是正定的。

于 2013-04-21T21:29:55.590 回答