用NumPy进行SVD - 结果的解释

Question

我正在尝试学习奇异值分解（SVD）。我找到了这个包含示例的YouTube讲座。然而，当我在numpy中尝试这个示例时，我得到了“有点”不同的结果。在这个示例中，输入矩阵是：

A = [ [1,1,1,0,0], [3,3,3,0,0], [4,4,4,0,0], [5,5,5,0,0], [0,2,0,4,4], [0,0,0,5,5], [0,1,0,2,2] ]
A = np.asarray(A)
print(A)

[[1 1 1 0 0]
 [3 3 3 0 0]
 [4 4 4 0 0]
 [5 5 5 0 0]
 [0 2 0 4 4]
 [0 0 0 5 5]
 [0 1 0 2 2]]

这个矩阵的秩是3 (np.linalg.matrix_rank(A))。讲座指出奇异值的数量就是矩阵的秩，在这个例子中，Sigma矩阵S的大小确实为3=3。然而，当我执行以下操作时：

U, S, V = np.linalg.svd(A)

矩阵S包含5个值。另一方面，前3个值与示例中的值匹配，其余2个基本上为0。我是否可以假定由于SVD算法背后的数值算法和计算机上实数的有限表示等原因，可以获得更多的奇异值而不是排名 - 或者类似的东西？

- Christian

---

看起来像是数值精度不准确的问题，因为np.allclose(A,np.dot(U * S, V))仍然返回True。 - DYZ

Answer 1

如此页面所述，numpy在内部使用LAPACK例程_gesdd来获取SVD分解。现在，如果您查看_gesdd文档，它提到：

为了找到一般矩阵A的SVD，请调用LAPACK例程?gebrd或?gbbrd将A通过酉（正交）变换减少为一个双对角矩阵B：A = QBPH。然后调用?bdsqr以形成双对角矩阵的SVD：B = U1ΣV1H。

因此，这里涉及到两个步骤：

• 通过正交变换（Householder转换）进行双对角化
• 使用隐式零移位QR算法获取双对角矩阵的SVD。

QR算法是一种迭代算法，意味着您不会得到“精确”的答案，但每次迭代都会得到越来越好的近似值，并在变化的值低于阈值时停止，因此从这个意义上讲它是“近似”的。

因此，除了由于实数的有限机器表示而产生的数值精度问题外，即使我们具有无限的表示容量，由于算法的迭代性质，我们也会在有限的时间内得到“近似”的结果（如果我们运行算法）。

- Deepak Saini

---

感谢详细的解释。只是补充一下，当我将所有三个矩阵“缩小”到预期的形状，例如 U = np.delete(U, np.s_[rank:], 1) 等等，并基于缩小后的矩阵重构 A 时，np.allclose() 也返回 True。 - Christian