[Math] Cholesky for Non-Positive Definite Matrices

cholesky decompositiondiagonalizationlinear algebramatricesreference-request

I am trying to approximate a NPD matrix with the nearest PD form and compute its Cholesky decomposition.

I know that the usual method is to perform an eigenvalue decomposition, zero out the negative eigenvalues, and recompute the eigenvalue decomposition. However, I don't want to do any eigenvalue decompositions, and I want to deal with this issue within the Cholesky algorithm.

For those of you who are familiar with the Cholesky algorithm, there is a step where you compute the diagonal of the lower diagonal matrix as the square root of a value (let's say x). If x<0 then, this means the matrix is NPD.

A simple way to deal with this could be to set x = 0 or x = 10^-10, just to work around this problem. However, it lacks the mathematical rigor.

My question is, has any of you heard such a method, or something similar to the above method, which is either published or somewhat proven?

Thanks!

Best Answer

You should be a bit more precise what you mean by NPD. My guess is: a symmetric/Hermitian (so, indefinite) matrix.

There is a Cholesky factorization for positive semidefinite matrices in a paper by N.J.Higham, "Analysis of the Cholesky Decomposition of a Semi-definite Matrix". I don't know of any variants that would work on indefinite matrices and find the closest positive (semi)definite matrix, but read this paper and see if you can work something out.

You might also be interested in Bunch-Parlett's symmetric indefinite decomposition described in their classic paper "Direct Methods for Solving Symmetric Indefinite Systems of Linear Equations".

To the best of my knowledge, the closest approximations like this have come to focus relatively recently and there are less demanding problems that are not yet solved, so I doubt that you'll find something already done for what you're interested in.

Edit

Of course, if you use the eigenvalue decomposition, you can get the closest semidefinite matrix. However, note that there is no need to multiply the obtained factors. Let's say you have $A = U \Lambda U^T$, where $U$ is orthogonal and $\Lambda$ is diagonal. Obtain nonnegative $\Lambda'$ by replacing negative elements in $\Lambda$ by zero and get the closest semidefinite aproximation of $A$: $A' = U \Lambda' U^T$.

Now, take QR of $(U\sqrt{\Lambda'})^T$. This will give you orthogonal $Q$ and upper triangular $R$ such that $QR = (U\sqrt{\Lambda'})^T$. Now, you have: $$A' = (U\sqrt{\Lambda'}) (U\sqrt{\Lambda'})^T = (QR)^T (QR) = R^T R,$$ which is a Cholesky(-like) decomposition you want.

I don't know why you want to avoid the eigenvalue decomposition. This looks quite clean to me.