最小二乘法：正规方程与 svd

Question

我尝试编写自己的线性回归代码，遵循正常方程beta = inv(X'X)X'Y。lstsq然而，平方误差比中的函数大得多numpy.linalg。有人可以向我解释为什么 SVD 方法（lstsq 使用的）比正规方程更准确吗？谢谢

Answer 1

我怀疑您的数据矩阵的条件数X'X很高。尝试计算此类矩阵的数值逆可能会导致较大的错误。显式计算逆矩阵通常是一个坏主意（例如，请参见http://www.johndcook.com/blog/2010/01/19/dont-invert-that-matrix/或http://epubs .siam.org/doi/abs/10.1137/1.9780898718027.ch14）。

您可以使用 检查条件编号numpy.linalg.cond。

这是一个例子。首先创建X并Y：

In [186]: X = np.random.randn(500, 30)

In [187]: Y = np.linspace(0, 1, len(X))

对于这个 random X，条件数并不大：

In [188]: np.linalg.cond(X.T.dot(X))
Out[188]: 2.4456380658308148

正规方程 并lstsq给出相同的结果（根据numpy.allclose使用该函数的默认参数时）：

In [189]: betan = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(Y)

In [190]: betal, res, rnk, s = np.linalg.lstsq(X, Y)

In [191]: np.allclose(betan, betal)
Out[191]: True

现在X通过使两列几乎相同来进行调整。这使得X'X几乎是奇异的，并给它一个很大的条件数：

In [192]: X[:,0] = X[:,1] + 1e-8*np.random.randn(len(X))

In [193]: np.linalg.cond(X.T.dot(X))
Out[193]: 3954529794300611.5

现在正规方程给出的结果与 不同lstsq：

In [194]: betan = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(Y)

In [195]: betal, res, rnk, s = np.linalg.lstsq(X, Y)

In [196]: np.allclose(betan, betal)
Out[196]: False