Check for positive definiteness or positive semide

Question

I want to check if a matrix is positive definite or positive semidefinite using Python.

How can I do that? Is there a dedicated function in SciPy for that or in other modules?

Answer 1

One good solution is to calculate all the minors of determinants and check they are all non negatives.

Answer 2

an easier method is to calculate the determinants of the minors for this matrx.

Answer 3

Check whether the whole eigenvalues of a symmetric matrix A are non-negative is time-consuming if A is very large, while the module scipy.sparse.linalg.arpack provides a good solution since one can customize the returned eigenvalues by specifying parameters.(see Scipy.sparse.linalg.arpack for more information)

As we know if both ends of the spectrum of A are non-negative, then the rest eigenvalues must also be non-negative. So we can do like this:

from scipy.sparse.linalg import arpack
def isPSD(A, tol = 1e-8):
    vals, vecs = arpack.eigsh(A, k = 2, which = 'BE') # return the ends of spectrum of A
    return np.all(vals > -tol)

By this we only need to calculate two eigenvalues to check PSD, I think it's very useful for large A

Answer 4

I assume you already know your matrix is symmetric.

A good test for positive definiteness (actually the standard one !) is to try to compute its Cholesky factorization. It succeeds iff your matrix is positive definite.

This is the most direct way, since it needs O(n^3) operations (with a small constant), and you would need at least n matrix-vector multiplications to test "directly".

Answer 5

Cholesky decomposition is a good option if you're working with positive definite (PD) matrices.

However, it throws the following error on positive semi-definite (PSD) matrix, say,

A = np.zeros((3,3)) // the all-zero matrix is a PSD matrix
np.linalg.cholesky(A)
LinAlgError: Matrix is not positive definite -
Cholesky decomposition cannot be computed

For PSD matrices, you can use scipy/numpy's eigh() to check that all eigenvalues are non-negative.

>> E,V = scipy.linalg.eigh(np.zeros((3,3)))
>> E
array([ 0.,  0.,  0.])

However, you will most probably encounter numerical stability issues. To overcome those, you can use the following function.

def isPSD(A, tol=1e-8):
  E,V = scipy.linalg.eigh(A)
  return np.all(E > -tol)

Which returns True on matrices that are approximately PSD up to a given tolerance.