(some) generalized eigenvalues in numpy
Hi, in econometrics/statistics it's relatively common to solve a special generalized eigenvalue problem, where the matrices are real, symmetric, and positive definite. It bothered me to depend on scipy just because of that. Also, scipy.linalg.eig is much more general and returns complex types which is inconvenient for sorting etc. So I've written the workaround function below, but it's probably not very efficient (at least I hope the algebra is right).
I'd be happy to hear your comments, suggestions; first with respect to the function itself, but also whether you agree that the implemented functionality would be useful enough to include in official numpy (in a more professional way than I have done).
cheers, Sven
def geneigsympos(A, B): """ Solves symmetricpositivedefinite generalized eigenvalue problem Az=lBz.
Takes two realvalued symmetric matrices A and B (B must also be positivedefinite) and returns the corresponding generalized (but also realvalued) eigenvalues and eigenvectors. Return format: as in scipy.linalg.eig, tuple (l, Z); l is directly taken from eigh output (a 1dim array of length A.shape[0] ?), and Z is an array or matrix (depending on type of input A) with the corresponding eigenvectors in columns (hopefully). Eigenvalues are ordered ascending (thanks to eigh). """ from numpy import asmatrix, asarray, linalg #fixme: input checks on the matrices LI = asmatrix(linalg.cholesky(B)).I C = LI * asmatrix(A) * LI.T evals, evecs = linalg.eigh(C) if type(A) == type(asarray(A)): output = "array" # A was passed as numpyarray else: output = "matrix" #but the evecs need to be transformed back: evecs = LI.T * asmatrix(evecs) if output == "array": return evals, asarray(evecs) else: return asmatrix(evals), evecs
participants (1)

Sven Schreiber