BrevFeed
Sign in

36. Singular Value Decomposition (SVD)

📖 Reading · 9 min
💡 Every code box below is live — edit it and hit Run.

Every decomposition so far in this module needed A square and symmetric. This lesson removes both restrictions at once: the Singular Value Decomposition exists for any matrix, of any shape, and it's arguably the single most useful factorization in applied linear algebra — §24.1's PCA and §22.2's condition numbers both build directly on it.

Construction, from the Spectral Theorem

For any m\times n matrix A, consider A^TA — always symmetric ((A^TA)^T=A^TA^T{}^T=A^TA, §16.2) and always PSD (§21.2's B^TB fact, with B=A). By §21.0, A^TA=V\Sigma^2V^T for orthogonal V and \Sigma^2=\operatorname{diag}(\sigma_1^2,\dots,\sigma_n^2) with every \sigma_i\ge0 (square roots of A^TA's nonnegative eigenvalues, ordered largest to smallest). The \sigma_i are the singular values of A.

For each \sigma_i>0, define \vec u_i=\frac1{\sigma_i}A\vec v_i (where \vec v_i is column i of V) — these turn out orthonormal too, and extend to a full orthonormal basis U of \mathbb{R}^m. The result:

A=U\Sigma V^T

U (m\times m) and V (n\times n) both orthogonal, \Sigma (m\times n) diagonal (in the rectangular sense: nonzero only on the leading diagonal, zero-padded otherwise) with \sigma_1\ge\sigma_2\ge\cdots\ge0.

SVD generalizes the Spectral Theorem

If A is symmetric PSD, its SVD and its Spectral Theorem decomposition coincide — U=V=Q and \Sigma=\Lambda — since A^TA=A^2=Q\Lambda^2Q^T directly gives \sigma_i=\lambda_i in that case. SVD is a genuine generalization: it works for non-symmetric, non-square, even singular matrices, where eigenvalues alone (possibly complex, possibly undefined for non-square A) aren't available at all.

Reading off everything else from SVD

  • Rank: \operatorname{rank}(A)= number of nonzero singular values — a numerically robust rank test (§22.2), since "how many singular values are essentially zero" tolerates rounding noise far better than pivot-counting in RREF (§16.1) ever could.
  • Column space, null space: the \vec u_i for nonzero \sigma_i span \operatorname{Col}(A); the \vec v_i for zero \sigma_i span \operatorname{Null}(A) — every one of §17.4's fundamental subspaces, read directly off U and V.
  • Best low-rank approximation: truncating the sum A=\sum_i\sigma_i\vec u_i\vec v_i^T to its k largest terms gives the provably best rank-k approximation of A (the Eckart-Young theorem, not proved here) — the mathematical foundation of §24.1's PCA, and of image/data compression generally.
  • Pseudoinverse: A^+=V\Sigma^+U^T (invert each nonzero \sigma_i, leave the zeros as zero) generalizes A^{-1} to non-square or singular A, solving §20.5's least-squares problem even when A^TA isn't invertible.

Doing it in Python

import numpy as np

A = np.array([[3., 2., 2.], [2., 3., -2.]])   # 2x3, not square
U, s, Vt = np.linalg.svd(A)

print("singular values:", [round(v, 4) for v in s.tolist()])
print("U shape:", U.shape, " Vt shape:", Vt.shape)

Sigma = np.zeros_like(A)
np.fill_diagonal(Sigma, s)
reconstructed = U @ Sigma @ Vt
print("U Sigma V^T == A:", bool(np.allclose(reconstructed, A)))
print("rank (nonzero singular values):", int(np.sum(s > 1e-10)))

singular values: [5.0, 3.0]
U shape: (2, 2)  Vt shape: (3, 3)
U Sigma V^T == A: True
rank (nonzero singular values): 2

A best rank-1 approximation, via truncated SVD:

import numpy as np

A = np.array([[4., 0.], [3., -5.], [0., 4.]])
U, s, Vt = np.linalg.svd(A)

rank1 = s[0] * np.outer(U[:, 0], Vt[0, :])
print("A =")
for row in A.tolist():
    print([round(v, 4) for v in row])
print("best rank-1 approximation =")
for row in rank1.tolist():
    print([round(v, 4) for v in row])
print("approximation error (Frobenius norm):", round(float(np.linalg.norm(A - rank1)), 4))

A =
[4.0, 0.0]
[3.0, -5.0]
[0.0, 4.0]
best rank-1 approximation =
[1.0588, -1.7647]
[3.0, -5.0]
[-1.7647, 2.9412]
approximation error (Frobenius norm): 4.0

Worked example

Find the SVD of A=\begin{pmatrix}2&0\\0&0\end{pmatrix} by inspection.

A already acts diagonally: A\vec e_1=(2,0)=2\vec e_1, A\vec e_2=(0,0)=\vec0. So A scales the first axis by 2 and collapses the second entirely — reading singular values directly off this action: \sigma_1=2,\sigma_2=0.

\boxed{U=I,\quad\Sigma=\begin{pmatrix}2&0\\0&0\end{pmatrix},\quad V=I}

Sanity check. \operatorname{rank}(A)=1 (only one nonzero singular value) — confirmed directly: A's second row and second column are both entirely zero, so \operatorname{Col}(A)=\operatorname{span}\{(1,0)\} and \operatorname{Null}(A)=\operatorname{span}\{(0,1)\}, matching "\vec u_1=(1,0) spans the column space, \vec v_2=(0,1) spans the null space" from the "Reading off everything else" section exactly.

Your turn

1. A matrix has singular values 5,3,0. What is its rank?

2. True or false: every matrix, including non-square ones, has an SVD.

3. For a symmetric matrix with eigenvalues -2,4, what are its singular values? (Careful — singular values are never negative.)

Solutions

1. Rank 2 — the count of nonzero singular values, ignoring the zero.

2. True. This is exactly the point of the "Construction" section above: A^TA is symmetric and PSD for any A (§21.2), so the Spectral Theorem applies to it regardless of A's own shape — SVD never requires A to be square, unlike every earlier decomposition in this module.

3. Singular values 2,4. Singular values are \sigma_i=\sqrt{\lambda_i(A^TA)}, and for symmetric A, A^TA=A^2 has eigenvalues \lambda_i^2 — so \sigma_i=\sqrt{\lambda_i^2}=|\lambda_i|. For \lambda=-2: \sigma=|-2|=2. For \lambda=4: \sigma=|4|=4. (This is exactly why "SVD generalizes the Spectral Theorem" above only claimed the two decompositions coincide for PSD matrices — a negative eigenvalue still gives a nonnegative singular value, just not an identical one.)

Check yourself in code

Find the singular values and rank of A=\begin{pmatrix}1&2\\2&4\\0&0\end{pmatrix} (note: column 2 is 2\times column 1).

Print exactly this:

singular values: [5.0, 0.0]
rank: 1

import numpy as np

A = np.array([[1., 2.], [2., 4.], [0., 0.]])
U, s, Vt = np.linalg.svd(A)
print("singular values:", [round(v, 4) or 0.0 for v in s.tolist()])
# print the rank: the count of singular values above a small tolerance

import numpy as np

A = np.array([[1., 2.], [2., 4.], [0., 0.]])
U, s, Vt = np.linalg.svd(A)
print("singular values:", [round(v, 4) or 0.0 for v in s.tolist()])

rank = int(np.sum(s > 1e-10))
print("rank:", rank)

SVD, A=U\Sigma V^T, exists for any matrix by applying the Spectral Theorem to the always-symmetric-PSD A^TA. Singular values generalize eigenvalues (coinciding for PSD matrices, always \ge0 otherwise), rank is the count of nonzero ones, and truncating the SVD sum gives the provably best low-rank approximation — the engine behind §24.1's PCA.

Next, closing this module's theory arc: the most general statement of when a matrix admits an orthonormal eigenbasis at all — normal matrices, and the fully general Spectral Theorem.

← Previous Next →