Monte Carlo & Randomized Methods

randomized SVD

A giant data matrix — millions of rows and columns — often has hidden simplicity: almost all of its meaningful content lies in a handful of dominant directions, and the rest is noise. Computing the full singular value decomposition to find those directions is enormously expensive. Randomized SVD is a surprisingly powerful shortcut: it uses a dose of randomness to quickly find the few directions that matter, skipping the rest, and delivers a near-best low-rank approximation in a fraction of the time.

The idea is to first SKETCH the big matrix down to a small one using random projection. To approximate the top-k structure of an m-by-n matrix A: draw a random n-by-(k+p) matrix Omega (a few extra columns p for safety), form the sketch Y = A * Omega — multiplying A by random vectors mixes its columns and, with high probability, the columns of Y span almost the same space as A's top-k left singular vectors. Orthonormalize Y (a thin QR factorization) to get an m-by-(k+p) matrix Q whose columns form an orthonormal basis for that captured subspace. Then project A down: B = Q^T * A is now a tiny (k+p)-by-n matrix, and you run an ordinary, cheap SVD on B and lift the result back through Q. The cost drops from the full SVD's O(m*n*min(m,n)) to roughly O(m*n*k) plus small dense factorizations — a huge saving when k is much less than m and n.

Randomized SVD is a flagship of randomized numerical linear algebra and is the practical engine behind large-scale PCA, recommender systems, and data compression. The honesty: the result is APPROXIMATE, with an error governed by how fast A's singular values decay — it is excellent when the spectrum drops off quickly (a genuinely low-rank-plus-noise matrix) and weaker when the tail is flat. The oversampling parameter p and one or two 'power iterations' (replacing A*Omega with (A*A^T)^q * A * Omega to sharpen the dominant directions) control the accuracy. The error is RANDOM — it holds with high probability, not with certainty — but the failure probability is made vanishingly small, and the Eckart-Young theorem still tells you the best possible rank-k error you are aiming near.

A 100,000-by-50,000 ratings matrix is well approximated by rank 50. A full SVD is infeasible, but randomized SVD draws a 50,000-by-60 random Omega, forms Y = A*Omega, orthonormalizes to Q, computes B = Q^T*A (only 60-by-50,000), and SVDs that small B. The whole job runs in seconds and recovers the top 50 singular directions to high accuracy.

Random projection sketches the matrix small, then a cheap SVD recovers its top directions.

Accuracy depends on the singular values DECAYING fast: randomized SVD shines on low-rank-plus-noise matrices and weakens on flat spectra. The error is probabilistic — extremely reliable, but not guaranteed; oversampling and a power iteration or two buy more accuracy when the tail is slow.

Also called
randomized low-rank approximationRSVDrandom projection SVD隨機SVD隨機投影低秩逼近