Randomized SVD

Why Does This Exist?

Truncated SVD on a million-by-million matrix costs $O(n^3)$ and materialises dense factors that dwarf memory. Yet the goal is only the top $k$ directions, and randomness finds them cheaply. Randomized SVD multiplies $A$ by a thin random sketch $\Omega$ (say $n \times (k + 10)$), orthogonalises the sample $Y = A\Omega$ into $Q$, compresses to the small core $B = Q^T A$, and decomposes $B$ instead. Cost drops to $O(mn \log k)$ with structured sketches, and the error sits within a small factor of the optimal rank-$k$ error with overwhelming probability.

This powers PCA on data that never fits in RAM, fast LSA on term-document matrices, and near-real-time recommender updates. Exact foundations are in singular value decomposition; the compression payoff in PCA derivation.

The one idea you need first

A random projection preserves geometry: $Y = A\Omega$ captures the range of $A$'s top directions because random vectors have components along everything, amplified where $A$ stretches most. Oversampling (a few extra columns beyond $k$) plus one or two power iterations ($A A^T$ applied to the sketch) sharpens the capture on flat spectra.

Think of It Like This

How It Actually Works

The four stages for target rank $k$ with oversampling $p$:

• Sketch: draw Gaussian $\Omega$ ($n \times (k+p)$) and form $Y = A\Omega$. Optionally sharpen: $Y = (AA^T)^q A\Omega$ with $q = 1$ or $2$.
• Orthogonalise: thin QR, $Y = QR$. $Q$'s columns span the sampled range.
• Compress: $B = Q^T A$, small ($(k+p) \times n$). The giant problem now fits in memory.
• Solve small: SVD of $B = \tilde{U}\Sigma V^T$, lift $U = Q\tilde{U}$. Truncate to $k$.

Error bounds: within $(1 + \epsilon)$ of optimal with probability $1 - \delta$, where oversampling and power iterations trade cost against both. Fixed seeds make it reproducible; without them, two runs return different-but-equally-valid bases.

Worked example

Approximate the rank-one $A = \begin{bmatrix} 3 & 4 \\ 6 & 8 \end{bmatrix}$. Its Gram matrix $A^T A = \begin{bmatrix} 45 & 60 \\ 60 & 80 \end{bmatrix}$ has trace 125 and determinant $3600 - 3600 = 0$, so singular values $\sqrt{125} \approx 11.1803$ and $0$. Sketch with the deterministic $\Omega = [1, 0]^T$:

One column captures the entire range (the matrix is rank one, so any in-range sketch suffices). Normalise to $Q = [1, 2]^T/\sqrt{5}$, compress $B = Q^T A = [15/\sqrt{5}, 20/\sqrt{5}] = [3\sqrt{5}, 4\sqrt{5}]$, and $\lVert B \rVert = \sqrt{45 + 80} = \sqrt{125} \approx 11.1803$. The sketch route recovers $\sigma_1$ exactly, with zero dense $m \times n$ factor ever formed.

Code

The point is the check: confirm the sketched singular value matches the exact one.

1import numpy as np2
3
4A = np.array([[3.0, 4.0], [6.0, 8.0]])5exact = np.linalg.svd(A, compute_uv=False)6Y = A @ np.array([[1.0], [0.0]])  # one-column sketch7Q, _ = np.linalg.qr(Y)8B = Q.T @ A9sketched = np.linalg.svd(B, compute_uv=False)10print([round(float(v), 4) for v in exact])     # -> [11.1803, 0.0]11print(round(float(sketched[0]), 4))            # -> 11.1803

Watch Out For

The Quick Version

• Randomized SVD sketches $Y = A\Omega$, orthogonalises to $Q$, solves small $B = Q^T A$, lifts back.
• The rank-one demo recovers $\sigma_1 = 11.1803$ through a single sketch column.
• Cost falls from $O(n^3)$ to sketch-plus-small-solve; error stays near optimal with high probability.
• Oversample 5 to 10 and power-iterate once on flat spectra.
• Seed every sketch: different valid bases per run break caches and tests.

What to Read Next

• Singular Value Decomposition
• PCA Derivation
• Matrix Multiplication