# How to diagonalize a large sparse symmetric matrix, to get the eigen values and eigenvectors?

The matrix could be as large as $2500\times 2500$, what is the best algorithm to do that, is there some algorithm that is easy to write a program, is there any convenient packages for that?

There is a survey of decomposition algorithms in the first table of this NIPS paper.

It lists modern algorithms (with links to known implementations), including the stochastic decomposition of Halko et al., arguably the state-of-the-art method today.

You ask for convenient programming packages but don't state your platform or language of choice. Assuming it's:

• Python:
• use scipy for in-core decompositions (input must fit in RAM)
• use gensim for both in-core and out-of-core sparse decompositions (also supports incremental decomposition updates)
• Java:
• Mahout has several scalable decomposition algos
• LingPipe (in-core) supports missing input values
• C++
• redsvd (in-core) very clean and elegant, efficient implementation
• MATLAB
• pca.m by Mark Tygert, one of the co-authors of the stochastic method.

Your problem isn't particularly big though, so I guess any existing package (using the iterative Lanczos algorithm) will do fine, eigen decompositions have been around for a while.

• Welcome to the site! I'm sure you will receive enough upvotes from this question to post as many links as you want in the future. – Andy W Jul 19 '11 at 12:57
• +1 Great answer: it offers many solutions, provides links to them all, and adds some advice about newer algorithms. – whuber Jul 19 '11 at 13:42

Take a look at A Survey of Software for Sparse Eigenvalue Problems by Hernández et al.

I don't know much about eigenvalues or what they are applicable to, but R seems to have a built in function for this purpose named eigen(). Calculating the eigenvalues & eigenvectors for a 2500 * 2500 matrix took ~ 1 minute on my machine.

> sampData <- runif(6250000, 0, 2)
> x <- matrix(sampData, ncol = 2500, byrow = TRUE)
> system.time(eigen(x))
user  system elapsed
79.74    2.90   65.69


This question has also come up on Stack Overflow.

2500x2500 is not such a large problem. Even without leveraging the sparsity the SVD implementation of scipy.linalg is able to decompose it in less than a minute. See my answer to a related questions for more details.

For larger problems you will need to use the sparsity explicitly. The gensim project my help you for middle size problems that fit on a single computer but not in RAM and the mahout implementation is able to deal with sparse matrices that don't even fit on a single hard-drive.