计算低秩矩阵

Question

计算低秩矩阵

performancealgorithmmathmatrix

4

假设我有一个排名为3的3x3矩阵A。我想创建一个最接近A的秩为2的矩阵，其$ {l}_{2} $ / Frobenius范数最小。我们称这个矩阵为F。

可以通过SVD轻松实现，即如果$ A = U S {V}^{H} $通过SVD分解，则$ F = U \hat{S} {V}^{H} $。其中$ \hat{S} $与$ S $相同，只是最后一个奇异值为零。

问题是，是否有一种更少计算量的方法来创建F，而不使用SVD分解？

谢谢。

- Royi

3个回答

0

你可能已经考虑过了，但如果A是正常的，SVD可以通过特征值分解来计算。这相当于解决特征多项式，由于对于秩为3的矩阵而言，它是一个三次多项式，因此根据众所周知的三次公式找到根。

总的来说，SVD似乎也必须可归约为解决秩为3的矩阵的三次方程，但我不记得读过任何相关内容。快速搜索发现this piece of code声称可以用这种方式解决秩为3的SVD问题。不幸的是，没有附带的论文。如果使用此代码，请使用非正定矩阵进行测试。

编辑：第二次阅读时，似乎作者也在使用特征分解。在非正定矩阵上可能不太好，但我希望在这里被证明是错误的。

- kkm

0

您可以使用幂迭代来找到与最大奇异值对应的奇异向量。一旦找到了它，您可以使用修改后的幂迭代来找到第二大奇异值的向量；在每次迭代之后，您需要从部分解向量投影中减去最大奇异值的向量。

这是一种找到所有奇异向量的不好方法，但对于前两个奇异向量来说，应该可以胜任。

- Michael J. Barber

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- user85109 · Accepted Answer

如果你已知矩阵的秩为3，那么只需要进行3次Householder变换就足以将一个nxm矩阵化简为一个3xm的矩阵。接着，我们可以很容易地将其转换为一个特征值问题。计算特征多项式即可。例如，考虑以下矩阵（我会使用MATLAB来完成计算）：

>> A = rand(10,3);
>> B = A*rand(3,6);

显然，B会有3个等级，但如果你不相信我，等级可以证实这一点。

>> rank(B)
ans =
     3

>> size(B)
ans =
    10     6

B矩阵是一个10x6的矩阵，它的秩为3。SVD也证实了这一点。

>> svd(B)
ans =
          6.95958965358531
          1.05904552889497
         0.505730235622534
      7.37626877572817e-16
      2.36322066654691e-16
      1.06396598411356e-16

我有点懒得写Householder变换的代码（虽然我已经写了一些，但是……），所以我将使用QR分解来帮助我。

>> [Q,R] = qr(B,0);
>> C = Q(:,1:3)'*B
C =
   -2.0815   -1.7098   -3.7897   -1.6186   -3.6038   -3.0809
    0.0000    0.91329   0.78347   0.44597  -0.072369   0.54196
    0.0000    0.0000   -0.2285   -0.43721  -0.85949  -0.41072

在这里进行的乘法显示了我们在三个Householder变换后所看到的内容。如我们所预期的那样，它是上三角形的。接下来，计算特征多项式。我将在此使用自己的工具以符号方式进行操作，但计算只需进行一些代数运算即可。

>> sympoly lambda
>> det(C*C' - lambda*eye(3))
ans =
    13.8942 + 66.9996*lambda + 49.8132*lambda^2 + lambda^3

>> P = det(C*C' - lambda*eye(3))
P =
    13.8942 - 66.9996*lambda + 49.8132*lambda^2 - lambda^3

什么是P的根，即C*C'的特征值？

>> r = roots(P)
r =
       48.436
       1.1216
      0.25576

我们知道这里的特征值必须是奇异值的平方，因此一个好的测试是看看我们是否恢复了svd找到的奇异值。因此，再次扩展显示格式，我们可以看到它非常好地完成了这个任务。

>> sqrt(roots(P))
ans =
          6.95958965358531
          1.05904552889497
         0.505730235622533

当数学运作起来时，它可以很有趣。我们可以用这个信息做什么呢？如果我知道一个特定的特征值，我就可以计算相应的特征向量。基本上，我们需要解决线性3x3齐次方程组。

(C*C' - eye(3)*r(3)) * X = 0

我在这里懒得写代码来找到解决方案，但高斯消元法可以解决问题。

>> V = null((C*C' - eye(3)*r(3)))
V =
        -0.171504758161731
        -0.389921448437349
         0.904736084157367

所以我们有V，是C*C'的特征向量。通过以下调用svd，我们可以确信它是一个特征向量。

>> svd(C - V*(V'*C))
ans =
           6.9595896535853
          1.05904552889497
      2.88098729108798e-16

因此，通过从C在V方向上减去该分量，我们得到一个秩为2的矩阵。

同样地，我们可以将V转换回原始问题空间，并使用它来通过对B进行一次秩为1的更新来转换我们的原始矩阵B。

>> U = Q(:,1:3)*V;
>> D = B - U*(U'*B);

D的排名是什么？

>> svd(D)
ans =
          6.95958965358531
          1.05904552889497
      2.62044567948618e-15
      3.18063391331806e-16
      2.16520878207897e-16
      1.56387805987859e-16

>> rank(D)
ans =
     2

正如您所看到的，即使我进行了大量的数学计算，多次调用svd、QR、rank等函数，最终实际计算仍然相对简单。我只需要：

进行3次Householder变换（将其存储以备后用。请注意，Householder变换只是矩阵的秩1更新）
计算特征多项式
使用自己喜欢的立方多项式求根方法计算最小根
使用高斯消元恢复相应的特征向量
使用我们最初进行的Householder变换将该特征向量回退到原始空间
对矩阵进行秩1更新。

只要我们知道实际秩为3，这些计算步骤对于任何大小的矩阵都将快速高效。它甚至不值得写一篇论文。

编辑：

由于问题已经修改为矩阵本身的大小只有3x3，因此计算甚至更加简单。无论如何，我将保留帖子的第一部分，因为它描述了适用于任何秩为3的一般矩阵的完全有效的解决方案。

如果您的目标是消除3x3矩阵的最后一个奇异值，则在3x3上进行svd似乎是相当高效的。但是，通过间接手段生成最后一个奇异值也会导致一些精度损失。例如，在此处比较使用svd计算的最小奇异值和使用特征值计算的最小奇异值。因此，您可能会看到一些小误差，因为形成A'*A将导致某些精度损失。这种损失的程度可能取决于A的条件数。

>> A = rand(3);
>> sqrt(eig(A'*A))
ans =
        0.0138948003095069
         0.080275195586312
          1.50987693453097
>> svd(A)
ans =
          1.50987693453097
        0.0802751955863124
        0.0138948003095054

然而，如果你真的想自己完成这项工作，你可以这样做。

直接从3x3矩阵A'*A中计算特征多项式。
计算该三次多项式的根。选择最小的根。
生成相应的特征向量。
在A上进行一次秩一更新，杀死位于特征向量方向上的A部分。

这个方法跟简单地调用SVD然后进行秩一更新相比是否更简单或者计算效率更高呢？对于3x3的情况，我不确定是否值得这样做。因为3x3 SVD 的计算非常快速。