Why is the inverse of a symmetric matrix not symmetric?!

Hi all, As far as I know, the inverse of symmetric matrix is always symmetric. However, I have a symmetric covariance matrix, call it C, and when I invert it (below), the solution, invC, is not symmetric!

    >> invC = inv(C);   % (inefficient I know, but it should still work...)

    >> isequal(invC,invC')

    ans = 0

Has anyone had this issue? Can this be due to rounding errors? My matrix is 1810x1810 with many entries like 0.0055, etc.

Thanks in advance!

Answer 1

[Roger Stafford]

Yes, it's roundoff error. Instead of 'isequal' which demands exact equality, try displaying the difference invC-invC' to see if the differences fall within the range of what you would regard as reasonable round off errors. With a matrix which is close to being singular these can be surprisingly large sometimes.

[Ted]

Thanks for the tip.

The errors were small. I've now forced my "inverse" to be symmetric:

>> invCtrans = invC';

>> invC = triu(invC) + tril(invCtrans) - diag(diag(invC);

What I'm really trying to do now is use the chol() command for cholesky factorization. I know I need a symmetric positive definite matrix (spd), and I've checked by using eigs(invC,10,0), which tells me the 10 smallest eigenvalues are all positive. I'm still getting the following though:

>> U = chol(invC)

error using chol: matrix must be positive definite

... Any thoughts? Doesn't my eigenvalues test show that invC is spd?

[Ted]

Never mind. I was getting the 10 eigenvalues with smallest magnitude, rather than the "most-negative."

Answer 2

[Youssef Khmou]

hi,

Try to use a tolerance criterion :

 C=symdec(100,100);
 C=C/max(C(:))
 I1=inv(C);
 I2=inv(C');
 norm(I1-I2) % its not zeros but saturated to zero (1e-n , n>20 )

[Ted]

The "symdec" command doesn't help me. I don't have the Robust Control Toolbox...