The Singular Value Decomposition

Eigenvalues demanded a square, cooperative matrix. The singular value decomposition demands nothing. Every matrix, any shape, real or rectangular, factors into a rotation, a set of axis-aligned stretches, and another rotation. It is the closest thing linear algebra has to a universal normal form.

The factorization

For any real $A \in \mathbb{R}^{m\times n}$ there exist orthogonal matrices $U$ and $V$ and a diagonal $\Sigma$ so that

$$A = U\,\Sigma\,V^{\top}, \qquad U^{\top}U = I_m,\quad V^{\top}V = I_n,$$ $$\Sigma = \operatorname{diag}(\sigma_1,\dots,\sigma_r), \qquad \sigma_1 \ge \sigma_2 \ge \dots \ge \sigma_r > 0 .$$

The $\sigma_i$ are the singular values, the columns $u_i$ of $U$ the left singular vectors, and the columns $v_i$ of $V$ the right singular vectors. The defining relation is $A v_i = \sigma_i u_i$: each right singular vector is sent to a scaled left singular vector. Equivalently $A$ is a sum of rank-one pieces,

$$A = \sum_{i=1}^{r} \sigma_i\, u_i v_i^{\top} .$$

The $v_i$ are the orthonormal eigenvectors of $A^{\top}A$, the $u_i$ those of $AA^{\top}$, and $\sigma_i^{2}$ their shared eigenvalues.¹

The geometry: sphere to ellipsoid

Read $A = U\Sigma V^{\top}$ right to left on a unit vector $x$. First $V^{\top}$ rotates $x$ into the singular basis. Then $\Sigma$ stretches coordinate $i$ by $\sigma_i$. Then $U$ rotates the result out. Writing $x = \sum_i (v_i^{\top}x)\, v_i$,

$$A x = \sum_{i=1}^{r} \sigma_i\, (v_i^{\top} x)\, u_i .$$

So the image of the unit sphere is an ellipsoid whose semi-axes are the vectors $\sigma_i u_i$: the directions come from $U$, the lengths are the singular values. The rank $r$ is exactly the number of nonzero $\sigma_i$, the number of axes the ellipsoid actually has.

A maps the unit circle to an ellipseinteractive

domain: unit circle

image: ellipse

v₁ , σ₁u₁v₂ , σ₂u₂

A = [[2.0, 0.0], [1.0, 1.0]]

A₁₁ = 2.0A₁₂ = 0.0A₂₁ = 1.0A₂₂ = 1.0

σ₁

2.29

σ₂

0.87

σ₁·σ₂

2.00

|det A|

2.00

The circle stretches into an ellipse with semi-axes σ₁ and σ₂. Their product σ₁σ₂ equals |det A|, the area scale factor.

Rank, approximation, and the pseudoinverse

Truncating the sum after $k$ terms gives the best rank-$k$ approximation of $A$,

$$A_k = \sum_{i=1}^{k} \sigma_i\, u_i v_i^{\top}, \qquad \lVert A - A_k \rVert_2 = \sigma_{k+1}, \qquad \lVert A - A_k \rVert_F^{2} = \sum_{i>k} \sigma_i^{2} .$$

No rank-$k$ matrix comes closer in either the spectral or the Frobenius norm. This is the Eckart-Young theorem, and it is why the SVD underlies image compression, denoising, and latent-factor models.²

When $A$ is not invertible the SVD still supplies the right inverse. Inverting each nonzero stretch gives the pseudoinverse

$$A^{+} = V\,\Sigma^{+}\,U^{\top}, \qquad \Sigma^{+} = \operatorname{diag}\!\left(\tfrac{1}{\sigma_1},\dots,\tfrac{1}{\sigma_r}\right)^{\!\top}\!,$$

and $x = A^{+}b$ is the least-squares solution of $Ax = b$ of smallest norm.

PCA is the SVD of centered data

Center a data matrix $X \in \mathbb{R}^{n\times d}$ by subtracting the column means, then take $X = U\Sigma V^{\top}$. The right singular vectors $v_i$ are the principal directions, the eigenvectors of the sample covariance $\tfrac{1}{n-1}X^{\top}X$, and the variance captured along $v_i$ is $\sigma_i^{2}/(n-1)$. Principal component analysis is not a separate algorithm: it is the SVD of the centered data, and the leading singular values say how much of the spread lives in how few directions.¹

A 2×2 by hand

Take

$$A = \begin{bmatrix} 2 & 0 \\ 1 & 1 \end{bmatrix}, \qquad A^{\top}A = \begin{bmatrix} 5 & 1 \\ 1 & 1 \end{bmatrix} .$$

The singular values are the square roots of the eigenvalues of $A^{\top}A$. For a symmetric $2\times2$ matrix the eigenvalues follow from the trace and determinant. Here $\operatorname{tr} = 6$ and $\det = 5\cdot 1 - 1\cdot 1 = 4$, so

$$\lambda = \frac{6 \pm \sqrt{6^{2} - 4\cdot 4}}{2} = 3 \pm \sqrt{5}, \qquad \lambda_1 \approx 5.236,\quad \lambda_2 \approx 0.764 .$$

Taking square roots,

$$\sigma_1 = \sqrt{3+\sqrt5} \approx 2.288, \qquad \sigma_2 = \sqrt{3-\sqrt5} \approx 0.874 .$$

Two checks. The product of the singular values is the area scale factor: $\sigma_1\sigma_2 = \sqrt{\det(A^{\top}A)} = \sqrt{4} = 2$, and indeed $\det A = 2\cdot 1 - 0\cdot 1 = 2$. The sum of their squares is the squared Frobenius norm: $\sigma_1^{2} + \sigma_2^{2} = 6 = 2^{2} + 0^{2} + 1^{2} + 1^{2}$. Both hold, so the ellipse this $A$ draws has semi-axes near $2.29$ and $0.87$, enclosing an area $\pi\,\sigma_1\sigma_2 = 2\pi$ times that of the unit circle.

Footnotes

1. Trefethen & Bau, Numerical Linear Algebra (1997). ↩ ↩²

2. Strang, Introduction to Linear Algebra (5th ed., 2016). ↩