Near positive definite (NearPD) of a symmetric matrix in Julia

General Usage Performance
1

Hi all,

Is there a package in Julia which can be used to find a near positive definite matrix to a symmetric matrix?
In R there is a function nearPD that I can use or the function LaplacesDemon::as.positive.definite.

Thanks

4 Likes
2

ProximalOperators exposes the IndPSD type, whose proximal operator is the projection onto the positive semidefinite cone (with the Frobenius norm taken as distance).

More info here.

Is this what you’re looking for?

Edit: an example snippet

using ProximalOperators

# X is your symmetric matrix
P, _ = prox(IndPSD(), X)
9 Likes
3

@lostella This perfect! Exactly what I was looking for! This is awesome. Thank you very much!
I will need to incorporate that in a inplace function since I will use this in a MCMC loop.

Also, inv (inverse of a matrix) sometimes fails for a matrix borderline invertible. In R, there is a function as.inverse in the LaplacesDemon package which quickly replaces the smallest eigen values by an epsilon.
Are you aware of any function that can do that in Julia?

Thanks

4

You should be able to use the in-place version

prox!(P, IndPSD(), X)

which will write directly on a preallocated P.

I don’t know of anything built-in that takes care of the singularity issue: you could either modify the code from ProximalOperators to account for that (instead of projecting eigenvalues to nonnegative values, you project them to a strictly positive region) or just add epsilon*I to P.

1 Like
5

Great! Indeed adding epsilon*I to P seems to do the trick. This works just perfect for me now!
I greatly appreciated your help.
Thank you!

1 Like
6

Feel free to mark my answer as a solution, so that other people can easily spot it. And of course if a better solution is given, feel free to mark that instead :wink:

7

indeed your answer works for me! How do I do mark your answer as the solution?
I am new to this.

8

This is an old thread, to which Google directed me when I was looking for something equivalent to R’s nearPD. But these days I think a simpler solution is:

9

These implementations differ on whether they constrain the resulting matrix to have 1’s along the diagonal, as required for a correlation matrix. General positive definite matrices such as covariance matrices need not have unit diagonal.

From what I can tell:

NearestCorrelationMatrix.jl seems to require unit diagonal, without an option to turn that constraint off.

ProximalOperators.jl’s IndPSD has no option to enforce unit diagonal.

R’s nearPD gives a selectable option.

If the unit diagonal constraint isn’t required and the input is symmetric (as for OP), basically one merely needs to clip the eigenvalues at zero (or epsilon).

1 Like
10

For what it’s worth, here is my own working solution (probably suboptimal in performance) following the excellent sources (hope it’s available) of https://nhigham.com/2021/01/26/what-is-the-nearest-positive-semidefinite-matrix/:

# Needs LinearAlgebra
function nearestposdef(A; minimumeigenvalue = 1e-6)

    Evalues, Evector = eigen(A)

    newA = Evector * Diagonal(max.(Evalues, minimumeigenvalue)) * Evector'

    Symmetric(newA)

end
1 Like