How to find the nearest positive definite matix

Hi, I have a quesiton, I want to use the chol function in matlab but the function requires a strict positive definite matrix as input. How can I find the nearest positive definite matrix?

Thanks

Answer 1

[Bruno Luong]

First step of preparation is to "symmetrize" the matrix

A = 0.5*(A + A')

Why? Because MATLAB engine selects the algorithm depending on the matrix form, it can fail to detect the matrix is hermitian (symmetric) if there is numerical error and isequal(A, A') returns FALSE, even if A is very close to A'.

Next, if it's not enough find the smallest eigen value lbdmin of A (use for example EIGS function), then change A by adding

A + max(-lbdmin,0)*speye(size(A)) % for sparse matrix, to keep it sparse

for dense matrix you might prefer not bother of speye.

A + max(-lbdmin,0)*eye(size(A))

This is the theory, in practice you might add some margin to overcome floating point truncation by boosting with some factor.

boost = 2; % close to 1 is better
A + boost*max(-lbdmin,0)*speye(size(A));

NOTE: This is not the nearest matrix (the nearest is to project negative eigen space to 0 and untouch the positive one, see John's answer), but convenient to get SDP matrix. Nearest SPD of sparse matrix is likely a dense matrix, which might not be desirable for large-side sparse matrix.

[Ron van de Sand]

Hello,

I ran in into a similar problem. Can you suggest any referene (Paper, Book etc.) regarding this topic?

Answer 2

[John D'Errico]

You can always just download my nearestSPD function from the file exchange. (Note that it does come with a reference to how it works.)

https://www.mathworks.com/matlabcentral/fileexchange/42885-nearestspd

From some work by N. Higham: "The nearest symmetric positive semidefinite matrix in the Frobenius norm to an arbitrary real matrix A is shown to be (B + H)/2, where H is the symmetric polar factor of B=(A + A')/2."

http://www.sciencedirect.com/science/article/pii/0024379588902236

nearestSPD uses the above idea, then verifies the result is indeed SPD using chol to verify it has suceeded. If there is still an issue, it applies a further tweak as needed.

[Ron van de Sand]

Thnak you for your quick response! I just needed some more information about the topic as I'm starting to dive into the topic. Thank you.

[John D'Errico]

Then just read the Higham paper.