Singular Value Decomposition

SVD definition, low-rank approximation, connection to PCA

You work at a streaming service. Your rating matrix is 1 million users by 10,000 movies. Most entries are missing — users have rated a tiny fraction of the catalog. Your job is to predict the missing entries so you can recommend movies users have not seen yet.

The naive approach: fill in missing entries with the average rating. This ignores everything you know about which users are similar and which movies are similar. What you want is to find latent structure — users who like action movies, comedies, prestige dramas — and represent both users and movies in that latent space.

SVD gives you this. Factor the matrix $M \approx U Sigma V^T$. $U$ is a 1M × k matrix — each user represented as a $k$-dimensional latent preference vector. $V$ is a 10,000 × k matrix — each movie as a $k$-dimensional latent attribute vector. $Sigma$ is $k \times k$ diagonal — the importance of each latent factor. To predict user $i$'s rating of movie $j$: compute $u_i^T Sigma v_j$. The dot product measures how well user $i$'s preferences align with movie $j$'s attributes across all $k$ latent dimensions. This was the foundation of most Netflix Prize solutions.

Why not eigendecomposition instead? Eigendecomposition requires a square symmetric matrix. Your rating matrix is 1M × 10,000 — rectangular. SVD works on any matrix. The singular values $sigma_i = \sqrt{lambda_i(M^T M)}$ — they are the square roots of the eigenvalues of $M^T M$. Left singular vectors $U$ are eigenvectors of $MM^T$; right singular vectors $V$ are eigenvectors of $M^T M$. SVD is strictly more general.

Two SVD variants matter here, and the quiz below tests the difference. The full SVD of an m×n matrix keeps every dimension: U is m×m orthogonal, Σ is m×n (padded with zero rows or columns beyond the rank), V is n×n orthogonal — U and V are square because their extra columns span the null spaces of M and Mᵀ, not just the row/column space. The compact SVD drops that padding: for a rank-r matrix, U shrinks to m×r, Σ to r×r diagonal holding only the r nonzero singular values, V to n×r — this still reconstructs M exactly, since rank(M)=r, it just discards the null-space columns that contributed nothing. The truncated (rank-k) SVD used above for the ratings matrix goes further: keep only the top k<r singular vectors — U becomes m×k, Σ becomes k×k, V becomes n×k — which is lossy by design, since only the k dominant latent factors matter for prediction.

The Eckart-Young theorem proves that keeping only the top-$k$ singular vectors gives the best possible rank-$k$ approximation: $M_k = sum_{i=1}^k sigma_i u_i v_i^T$. No other rank-$k$ matrix is closer to $M$ in either Frobenius or spectral norm. This is the mathematical guarantee behind every truncated SVD application — compression, denoising, dimensionality reduction.

Rank itself is a hard quantity to optimise directly: minimising rank(M) subject to matching a set of observed entries is NP-hard, because rank counts how many singular values are nonzero — the same L0-style combinatorial difficulty as minimising the number of nonzero entries in a vector. The standard fix mirrors the vector case: just as the L1 norm (sum of absolute values) is the convex relaxation used to encourage a sparse vector, the nuclear norm ‖M‖* = Σσᵢ (the sum of the singular values) is the convex relaxation used to encourage a low-rank matrix. Minimising the nuclear norm instead of rank is a tractable convex problem, and because it penalises the sum of singular values rather than merely their count, it pushes small singular values toward zero — the matrix version of L1 shrinking small coefficients to zero. This is the trick behind matrix completion algorithms (recovering a full matrix from a few observed entries, e.g. filling in the rest of that ratings matrix).

NOT this. Most people treat SVD and PCA as synonyms. They are not. SVD is a matrix factorization — a numerical decomposition that works on any matrix. PCA is a statistical method for finding directions of maximum variance in a dataset. PCA uses SVD as its computational engine: on mean-centered data, the right singular vectors of the data matrix equal the eigenvectors of the covariance matrix, and the two methods produce the same result. But SVD is more general (works on rectangular matrices, handles non-statistical applications like pseudo-inverses), and numerically more stable than computing the covariance matrix $X^T X$ explicitly, which squares the condition number (condition number κ = σ_max/σ_min, the ratio of the largest to smallest singular value — the higher it is, the more small numerical errors get amplified when you invert or solve with that matrix).

Key points

• Use SVD when you need a low-rank approximation, a pseudo-inverse, or intrinsic dimensionality. Call `np.linalg.svd(X, full_matrices=False)` and inspect the singular value spectrum. A sharp elbow in the scree plot tells you the effective rank of your data. Singular values below the elbow are noise — including them adds variance without signal to any downstream model.
• The production trap: computing PCA via eigendecomposition of the covariance matrix $X^T X$. This squares the condition number: if $κ(X) = 100$, then $κ(X^T X) = 10,000$. Numerical errors get amplified by $κ^2$. Always compute PCA via SVD of the data matrix directly. sklearn's `PCA` does this by default. If you wrote your own PCA from scratch using `np.linalg.eig` on $X^T X$, replace it.
• The diagnostic: plot the singular value spectrum and check the condition number $κ = sigma_{max}/sigma_{min}$. A condition number above $10^6$ means the matrix is nearly singular and any computation involving its inverse (least squares, linear regression) will be numerically unstable. The pseudoinverse $A^+ = V Sigma^+ U^T$ handles this by zeroing near-zero singular values rather than inverting them — sklearn's LinearRegression uses this by default.

Recap

• SVD $M \approx U\Sigma V^T$ factors any matrix — rectangular OK, unlike eigendecomposition (square only).
• Latent factors: U = user vectors, V = item vectors, Σ = factor importance; predict via $u_i^T\Sigma v_j$ (Netflix Prize).
• Singular values $\sigma_i=\sqrt{\lambda_i(M^TM)}$; U eigenvectors of MMᵀ, V of MᵀM.
• Eckart-Young: top-k singular vectors give the best rank-k approximation — the guarantee behind truncated SVD.
• SVD ≠ PCA: SVD is a matrix factorisation; PCA is variance-finding that *uses* SVD as its engine.
• Compute PCA via SVD of X, not eig of XᵀX — forming XᵀX squares the condition number ($\kappa^2$).
• Condition number $\kappa=\sigma_{max}/\sigma_{min}$; κ > 1e6 → nearly singular, inverse-based computations unstable.

Check your understanding