I have a dense real symmetric square matrix. The dimension is about 1000x1000. I need to compute the first principal component and wonder what the best algorithm to do this might be.
It seems that MATLAB uses the Arnoldi/Lanczos algorithms (for eigs). But from reading about them I'm not sure whether they have any advantages over simple power iteration, since my matrix is not sparse and I'm only interested in the first eigenvector.
Any recommendations what's the fastest algorithm in this case?
The fastest method will likely depend upon the spectrum and normality of your matrix, but in all cases Krylov algorithms should be strictly better than power iteration. G.W. Stewart has a nice discussion of this issue in Chapter 4, Section 3 of Matrix Algorithms, Volume II: Eigensystems:
The power method is based on the observation that if $A$ has a dominant eigenpair then under mild restrictions on $u$ the vectors $A^k u$ produce increasingly accurate approximations to the dominant eigenvector. However, at each step the power method considers only the single vector $A^k u$, which amounts to throwing away information contained in the previously generated vectors. It turns out that this information is valuable..."
and he goes on to show that, for $100 \times 100$ diagonal matrix with the $i$'th diagonal value set to $.95^i$ (counting from $i=0$), after 25 iterations the Krylov subspace captures the dominant eigenvector eight orders of magnitude better than power iteration.
Power iteration is the simplest, but as mentioned above it would likely converge very slowly if the matrix is very non-normal. You get a "hump" phenomenon where the sequence appears to diverge for many iterations before asymptotic behavior kicks in.
Since your matrix is symmetric you could consider RQI iterations, which in the symmetric case yields cubic convergence: http://en.wikipedia.org/wiki/Rayleigh_quotient_iteration .
What makes Arnoldi or Lanczos iterations very nice (at least in my opinion, but I do not research numerical linear algebra) is that they are very versatile. It is usually possible to control which eigenvalues they give you, and how many you get. This is especially true in the symmetric case (and even better if your matrix is definite). For symmetric problems they are very robust. As a black box they work well, but they are also very receptive to new problem information, such as the ability to solve systems involving the matrix.
