Improving Precision of Eigenvectors with Large Eigenvalues

Question

bil 2024-6-24

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2131286-improving-precision-of-eigenvectors-with-large-eigenvalues

评论： Christine Tobler 2024-6-26

Hi all,

This is a bit of a generic question, but I was hoping someone could provide some insight on how I could improve the precision of eigenvectors using "projection techniques", like in this post, where a matrix will have 1 or 2 very large eigenvalues, and the remaining eigenvalues are much smaller (by several orders of magnitude). They have code written in R, which I am not too familiar with, but is the idea generically that I should subtract out the overlap of smaller eigenvectors with those of larger eigenvectors to improve their precision? When I say precision, what I mean is that applying eig to a matrix M will generate a set of eigenvectors, but once I check M*v - λ*v, the resultant array will deviate from 0, i.e. applying M to the "eigenvector" has resulted in a linear combination of multiple other eigenvectors.

Any guidance is appreciated.

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Christine Tobler 2024-6-24

2
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2131286-improving-precision-of-eigenvectors-with-large-eigenvalues#answer_1476266

The linked post is about a symmetric matrix, is this also your case? In that case (if issymmetric returns true for your matrix), the eigenvectors returned by EIG will be orthogonal up to numerical round-off (an exact 0 is not possible in numerical computation, practically speaking).

The proposed solution doesn't improve the eigenvectors, but instead applies a projection to the computed residual vectors, to project out any components along the eigenvectors of the larger eigenvalues. Whether this is useful will depend on your application - that is, do you need to do computations with that residual?

11 个评论
显示 9更早的评论隐藏 9更早的评论

Christine Tobler 2024-6-24

在 MATLAB Online 中打开

You're right - they don't improve the residual, although it turns out they also don't improve the eigenvectors. I had the R-code autotranslated to MATLAB (don't know much R myself), and here's what it seems to do:

% Matrix copied over from

% https://stats.stackexchange.com/questions/539436/get-accurate-eigenvectors-when-eigenvalues-are-minuscule

M = [

0.002307414, -1.318435973, 8.63e-11, -3.33e-11, 0.486238053, 0.41182779, -0.169525454, 0.624275558,

-1.318435973, 753.3425186, -4.93e-08, 1.90e-08, -277.8320732, -235.3147146, 96.86532753, -356.7054684,

8.63e-11, -4.93e-08, 3.23e-18, -1.24e-18, 1.82e-08, 1.54e-08, -6.34e-09, 2.34e-08,

-3.33e-11, 1.90e-08, -1.24e-18, 4.79e-19, -7.01e-09, -5.94e-09, 2.44e-09, -9.00e-09,

0.486238053, -277.8320732, 1.82e-08, -7.01e-09, 102.4642297, 86.78386443, -35.72384951, 131.5526701,

0.41182779, -235.3147146, 1.54e-08, -5.94e-09, 86.78386443, 73.50310589, -30.25693671, 111.4208258,

-0.169525454, 96.86532753, -6.34e-09, 2.44e-09, -35.72384951, -30.25693671, 12.45501408, -45.8654479,

0.624275558, -356.7054684, 2.34e-08, -9.00e-09, 131.5526701, 111.4208258, -45.8654479, 168.8989909

];

% Extract the size, verify M is exactly symmetric, and call eig

n = size(M, 1);

issymmetric(M)

ans = logical

1

[U, D] = eig(M);

% We can see the residual is accurate to machine precision. Also, the

% eigenvalues are orthogonal to each other up to machine precision.

residual = norm(M*U - U*D)

residual = 4.7680e-13

orthogonality = norm(U'*U - eye(n))

orthogonality = 9.7243e-16

% Extract eigenvalues into a vector, logical map for the two "large"

% eigenvalues.

lambda = diag(D);

largeEigenvalue = lambda >= max(lambda) * 1e-12;

for ii=1:500

% Make a random vector dx that should only be comprised of components

% along the small eigenvectors:

randVec = randn(8, 1);

randVec(largeEigenvalue) = 0; % set to zero for large eigenvalues

dx = U*randVec;

dx = dx / norm(dx);

normOrig(ii) = norm(M*dx);

% Correct this random vector dx by projecting away components along the

% large eigenvectors that have been introduced by round-off error in

% U*randVec operation.

dxCorrected = dx - U(:, largeEigenvalue) * (U(:, largeEigenvalue)' * dx);

dxCorrected = dxCorrected / norm(dxCorrected);

normCorrected(ii) = norm(M*dxCorrected);

end

semilogy(normOrig, 'x'); hold on; semilogy(normCorrected, 'o')

legend('orig', 'corrected')

So in this case, the accuracy of the original and the corrected version are about the same, and about match the accuracy of the corrected result from the blog post you linked.

Perhaps the eigenvalues and eigenvectors computed there are less accurate? Or I may have mixed something up - as I said, I don't know R very well and am relying on some automatic translation.

Torsten 2024-6-26

But I suspect this might not be the most accurate and the most straightforward method is, as you said, to just take the reciprocal of the eigenvalues of M to get the eigenvalues of M^(-1), and the eigenvectors are the same for both matrices.

Yes, at least I cannot think of any advantage to work with the inverse for your case.

Christine Tobler 2024-6-26

Closing the loop, I agree with Torsten that computing the eigenvalues and eigenvectors of the original matrix and then inverting the eigenvalues will be more accurate than calling eig on the inverse.

请先登录，再进行评论。