低秩矩阵的计算

Question

鉴于我有一个等级3的矩阵3x3 A.我想创建一个秩2的矩阵,它在$ {l} _ {2} $/Frobenius范数中最接近A.我们称之为矩阵F.

很容易通过SVD实现,即如果$ A = US {V} ^ {H} $由SVD分解$ F = U\hat {S} {V} ^ {H} $.其中$\hat {S} $与$ S $相同,最后一个奇异值归零.

问题是,是否有一个计算密集度较低的方法来创建F但使用SVD分解？

谢谢.

Answer 1

如果你知道矩阵是秩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

因此,通过在V方向上减去C的那个分量,我们得到一个秩2矩阵.

类似地,我们可以将V转换为原始问题空间,并使用它将矩阵B(我们的原始矩阵)转换为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等,但最后,实际计算相对微不足道.我只需要......

1. 3户主变革.(存储它们以供以后使用.请注意,Householder转换只是对矩阵的一级更新.)
2. 计算特征多项式.
3. 使用您喜欢的方法计算最小根用于三次多项式.
4. 恢复相应的特征向量.高斯消除就足够了.
5. 使用我们最初做的户主变换将特征向量引回到原始空间.
6. 排名第一的矩阵更新.

所有这些计算步骤对于任何大小矩阵都是快速有效的,只要我们知道实际等级为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

但是,如果你真的想自己做这项工作,你就可以做到.

1. 直接从3x3矩阵A'*A计算特征多项式.
2. 计算该三次多项式的根.选择最小的根.
3. 生成相应的特征向量.
4. 对A进行一级更新,杀死位于该特征向量方向的A部分.

这是否比对svd的简单调用更简单或更少计算效率,然后是一级更新？我完全不确定3x3的努力是否值得.一个3x3的DVD确实非常快速的计算.