Group Abstract Group Abstract

Message Boards Message Boards

0
|
474 Views
|
5 Replies
|
5 Total Likes
View groups...
Share
Share this post:

Unexpected output from matrix inversion

Posted 1 month ago

Hi,
I hope you are all doing well. I am trying to find the inverse of this matrix, but as you can see, it is giving me problems, and these values are definitely not correct. I tried both PseudoInverse and Inverse but they gave me the same results, If any of you have an idea where the problem may come from it will help me a lot. Thank you all in advance.

POSTED BY: Leo Murphy
5 Replies
Posted 1 month ago

The "matrice" in your example is ill - conditioned; thus, the computation of its inverse matrix is unstable (e.g., a small change in a matrix causes huge changes in the inverse matrix). I believe you already know the followings but just a reminder. If a matrix, X, is ill - conditioned, the inverse of the matrix is a "meaningless bad approximation". The condition number of the "matrice" in your example, (i.e., ratio of max singular value divided by min singular value), is gigantic (i.e., 2.6*10^10). Orthogonal matrix should have singular value of one.

Ill-conditioning is also known as multicollinearity in regression, confounding (or alias structure) in design of experiments, rank-deficient in matrix algebra. Note that the inverse of an orthogonal matrix is the same as its transpose (i.e., the U and V are orthogonal matrices, U'U=I, V'V=I.

X=U*S*V'.

I always compute Singular Value Decomposition (SVD) if a matrix computation is involved because SVD is the most reliable way to compute an inverse matrix. I checked Mathematica commands such as PseudoInverse, Inverse, LinearSolve, LeastSquares for ill-conditioning in the past and these commands do not provide a warning sometimes. Matlab has also the same issue.

{u, s, v} = SingularValueDecomposition[matrice] ListPlot[Diagonal[s], PlotRange -> {All, {0, 60}}] condNumber=Max[Diagonal[s]]/Min[Diagonal[s]] (* 2.6*10^10 *)

POSTED BY: Sangdon Lee
Posted 1 month ago

Thank you for your informative response Like you said, even in this case, using Least Squares and LinearSolve didn't help solve the problem. I tried Tikhonov regularization

alpha = 10^(-6);
Ainv = Inverse[matrice + alpha IdentityMatrix[Length[matrice]]];

It seems to give better results, but it's still not as precise as I wanted

POSTED BY: Leo Murphy

I’m not seeing any error or warning message regarding the matrix computation (as distinct from its display). Is there something not shown that implies the inversion gave a bad result?

POSTED BY: Daniel Lichtblau
Posted 1 month ago

Yes there are no warnings or errors but the values of the inverted matrix are definitely incorrect (10^6, 10^8......)

POSTED BY: Leo Murphy

I'm with Daniel. The following check returns zero, which is hard to improve upon: Norm[Inverse[matrice].matrice-IdentityMatrix[Length@matrice]]. If matrice has an eigenvalue with a modulus as small as $10^{-8}$, the inverse will have one as large as $10^8$. And the inverse will have to have some entries of a similar order of magnitude.

The following visualizes the magnitudes of the eigenvalues:

ListLogPlot[Abs[Eigenvalues[matrice]]]

If the small eigenvalues are the result of roundoff and are supposed to be zero, then the matrix is not invertible. If you want the pseudoinverse, then increase the Tolerance to a level you find acceptable.

PseudoInverse[matrice, Tolerance -> 1/1000]

You can use SingularValueList[matrice] and SingularValueList[matrice, Tolerance -> 1/1000] to examine what is discarded.

POSTED BY: Michael Rogers
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard