# How to compute SVD of a huge sparse matrix?

What is the best way to compute singular value decomposition (SVD) of a very large positive matrix (65M x 3.4M) where data is extremely sparse?

Less than 0.1% of the matrix is non zero. I need a way that:

• will fit into memory (I know that online methods exists)
• will be computed in a reasonable time: 3,4 days
• will be accurate enough however accuracy is not my main concern and I would like to be able to control how much resources I put into it.

It would be great to have a Haskell, Python, C# etc. library which implements it. I am not using mathlab or R but if necessary I can go with R.

• How much memory do you have? 0.1% of 65M*3.4M is still 221e9 non zero values. If you use 4 bytes per value, that is still more than 55 gb assuming no overhead, so the sparsity still doesn't solve the problem... Do you need to load the whole set into memory at once? Oct 26, 2012 at 18:02
• I should have been more precise. No more than 250-500mb with 32-bit integer. Probably much less, but the dimensionalilty is the problem as I understand it. I have a 16GB machine. Oct 26, 2012 at 18:05
• This webpage links to a Python library which implements "a fast, incremental, low-memory, large-matrix SVD algorithm": en.wikipedia.org/wiki/Latent_semantic_analysis Oct 26, 2012 at 19:56
• May 23, 2017 at 19:16

If it fits into memory, construct a sparse matrix in R using the Matrix package, and try irlba for the SVD. You can specify how many singular vectors you want in the result, which is another way to limit the computation.

That's a pretty big matrix, but I've had very good results with this method in the past. irlba is pretty state-of-the-art. It uses the implicitly restarted Lanczos bi-diagonalization algorithm.

It can chew through the netflix prize dataset (480,189 rows by 17,770 columns, 100,480,507 non-zero entries) in milliseconds. You dataset is ~ 200,000 times bigger than the Netflix dataset, so it take significantly longer than that. It might be reasonable to expect that it could do the computation in a couple of days.

• the data matrix fits into memory, will irlba will handle the decomposition in a memory efficient way as well? Oct 26, 2012 at 17:56
• @Sonia: irlba is very memory efficient: it computes an approximate solution, you can limit the number of singular vectors, and it was designed to work on sparse matrices. As far as I know, it's as fast as you're going to get for computing partial SVDs.
– Zach
Oct 26, 2012 at 18:05
• @Sonia: Good luck!
– Zach
Oct 26, 2012 at 18:14
• Gave it a try - out of memory... I will compute a triangle block form before running it. Nov 3, 2012 at 13:13
• @Sonia do you have it stored as a sparse Matrix? Try limiting the number of singular values you compute... maybe just look at the top 10?
– Zach
Nov 3, 2012 at 18:43

If you're willing to have a low-rank approximation (as you would with Lanczos-type algorithms and a limited number of singular vectors), an alternative is stochastic SVD. You get similar accuracy and computational effort to things like irlba, but a much, much simpler implementation -- which is relevant if none of the available implementations handle precisely your situation.

If you start with an $$N\times M$$ matrix and you want a rank $$k$$ approximation, you need to be able to multiply your matrix by a dense matrix with, say, $$M\times (k+10)$$ entries and then do an ordinary SVD on the resulting $$N\times (k+10)$$ matrix. As with Lanczos-type algorithms, the computation only uses the matrix for multiplication (it's a "matrix-free" algorithm), so it will take advantage of sparsity and other structure. (My use case was in genetics and the matrix was the product of a sparse matrix and a projection orthogonal to a low-rank subspace)