How expensive is it to compute the eigenvalues of a matrix?

Question

Most of the algorithms for eigen value computations scale to big-Oh(n^3), where n is the row/col dimension of the (symmetric and square) matrix.

For knowing the time complexity of the best algorithm till date you would have to refer to the latest research papers in Scientific Computing/Numerical Methods.

But even if you assume the worse case, you would still need at least 1000^3 operations for a 1000x1000 matrix.

R uses the LAPACK routine's (DSYEVR, DGEEV, ZHEEV and ZGEEV) implementation by default. However you could specify the EISPACK=TRUE as a parameter to use a EISPACK's RS, RG, CH and CG routines.

The most popular and good open source packages for eigenvalue computation are LAPACK and EISPACK.

With big matrices you usually don't want all the eigenvalues. You just want the top few to do (say) a dimension reduction.

The canonical algorithm is the Arnoldi-Lanczos iterative algorithm implemented in ARPACK:

www.caam.rice.edu/software/ARPACK/

There is a matlab interface in eigs:

http://www.mathworks.com/access/helpdesk/help/techdoc/index.html?/access/helpdesk/help/techdoc/ref/eigs.html

eigs(A,k) and eigs(A,B,k) return the k largest magnitude eigenvalues.

And there is now an R interface as well:

http://igraph.sourceforge.net/doc-0.5/R/arpack.html

I assume it helps if the matrix is sparse?

Yes, there are algorithms, that perform well on sparse matrices.

See for example: http://www.cise.ufl.edu/research/sparse/

How long might it take in practice if I have a 1000x1000 matrix?

MATLAB (based on LAPACK) computes on a dual-core 1.83 GHz machine all eigenvalues of a 1000x1000 random in roughly 5 seconds. When the matrix is symmetric, the computation can be done significantly faster and requires only about 1 second.

I would take a look at Eigenvalue algorithms, which link to a number of different methods. They'll all have different characteristics, and hopefully one will be suitable for your purposes.

You can use the GuessCompx package from CRAN to estimate the empirical complexity of your eigenvalues computation and predict the full running time (although it's still small in your example). You need a little helper function because the fitting process only subsets the rows, so you must make the matrix square:

library(GuessCompx)
m = matrix(rnorm(1e6), ncol=1000, nrow=1000)
# custom function  to subset the increasing-size matrix to a square one:
eigen. = function(m) eigen(as.matrix(m[, 1:nrow(m)]))
CompEst(m, eigen.)
#### $`TIME COMPLEXITY RESULTS`
#### $`TIME COMPLEXITY RESULTS`$best.model
#### [1] "CUBIC"
#### $`TIME COMPLEXITY RESULTS`$computation.time.on.full.dataset
#### [1] "5.23S"
#### $`TIME COMPLEXITY RESULTS`$p.value.model.significance
#### [1] 1.784406e-34

You get a cubic complexity for time, and a Nlog(N) complexity for memory usage of the R base eigen() function. It takes 5.2 secs and 37Mb to run the whole computation.

Apache Mahout is an open-source framework built on map-reduce (i.e. it works for really really big matrices). Note that for a lot of matrix stuff the question isn't "whats the big-o runtime" but rather "how parallelizable is it?" Mahout says they use Lanczos, which can essentially be run in parallel on as many processors as you care to give it.