Linear Algebra (source: Gilbert Strang notes) · 31 / 37

Lecture 29: Singular Value Decomposition

This summary is based on the lecture's transcript on MIT OCW.

Strang calls the SVD "the final and best factorization of a matrix". It applies to every matrix, including rectangular and singular ones:

\[\boxed{A = U\Sigma V^\top}\]

  • \(U\) is orthogonal.
  • \(\Sigma\) is diagonal.
  • \(V\) is orthogonal.

The new ingredient is two different orthogonal matrices. Strang sees the SVD as bringing the whole course together.


1. The special case we already know

Symmetric positive definite \(A\). The spectral theorem gives

\[A = Q\Lambda Q^\top\]

  • The eigenvectors are orthonormal, so the eigenvector matrix \(S\) becomes an orthogonal \(Q\).
  • The eigenvalues are positive.

This already is an SVD, with \(U = V = Q\) and \(\Sigma = \Lambda\).

Why general matrices need more. For a general \(A\), the eigenvector matrix is not orthogonal, and \(A\) may not even be square. So eigenvectors alone won't work. We need two orthogonal matrices, \(U\) and \(V\).


2. What the SVD means geometrically

The four-subspace picture. \(A\) takes vectors in the row space to vectors in the column space.

The goal. Find an orthonormal basis \(v_1, \dots, v_r\) of the row space that \(A\) maps to an orthogonal set in the column space:

\[Av_i = \sigma_iu_i, \qquad i = 1, \dots, r\]

  • The \(u_i\) are orthonormal in the column space.
  • Each \(\sigma_i > 0\) is a stretching factor, called a singular value. Strang uses \(\sigma\) rather than \(\lambda\).

Why this is special. Any orthonormal basis of the row space is easy to find, for example by Gram–Schmidt. But a random orthonormal basis will generally not map to perpendicular vectors. The SVD finds the one basis that does.

Nullspaces are no trouble.

  • Complete the \(v\)'s with an orthonormal basis \(v_{r+1}, \dots, v_n\) of \(N(A)\). These give the full orthonormal basis of \(\mathbb{R}^n\).
  • Complete the \(u\)'s with an orthonormal basis \(u_{r+1}, \dots, u_m\) of \(N(A^\top)\). These give the full orthonormal basis of \(\mathbb{R}^m\).
  • The corresponding diagonal entries of \(\Sigma\) are zero.

Matrix form

Stack the equations \(Av_i = \sigma_iu_i\) as columns:

\[A\begin{bmatrix} v_1 & \cdots & v_r \end{bmatrix} = \begin{bmatrix} u_1 & \cdots & u_r \end{bmatrix}\begin{bmatrix} \sigma_1 & & \\ & \ddots & \\ & & \sigma_r \end{bmatrix}\]

That is, \(AV = U\Sigma\).

Solve for \(A\). \(V\) is square and orthogonal, so \(V^{-1} = V^\top\):

\[A = U\Sigma V^\top\]


3. How to find \(U\), \(\Sigma\), \(V\)

The trick: eliminate one orthogonal matrix at a time. Form \(A^\top A\) and use \(U^\top U = I\):

\[A^\top A = V\Sigma^\top\underbrace{U^\top U}_{I}\Sigma V^\top = V\,\Sigma^\top\Sigma\,V^\top\]

This is exactly the spectral decomposition \(Q\Lambda Q^\top\) of the symmetric, positive semidefinite matrix \(A^\top A\). Therefore:

  • The columns of \(V\) are the eigenvectors of \(A^\top A\).
  • \(\sigma_i^2\) are the eigenvalues of \(A^\top A\), so \(\sigma_i = \sqrt{\lambda_i(A^\top A)}\), taking the positive square root.

The same trick in the other order:

\[AA^\top = U\,\Sigma\Sigma^\top\,U^\top\]

  • The columns of \(U\) are the eigenvectors of \(AA^\top\).
  • The eigenvalues of \(AA^\top\) are the same \(\sigma_i^2\).

Why the eigenvalues match. \(AB\) and \(BA\) always have the same nonzero eigenvalues.


4. Example 1: an invertible matrix (with Strang's sign problem fixed)

\[A = \begin{bmatrix} 4 & 4 \\ -3 & 3 \end{bmatrix}\]

\(A\) is invertible, so it has rank 2.

Step 1: compute \(A^\top A\)

\[A^\top A = \begin{bmatrix} 25 & 7 \\ 7 & 25 \end{bmatrix}\]

Step 2: eigenvectors and eigenvalues of \(A^\top A\)

  • \(\lambda_1 = 32\), with eigenvector \(v_1 = \tfrac{1}{\sqrt2}(1, 1)\).
  • \(\lambda_2 = 18\), with eigenvector \(v_2 = \tfrac{1}{\sqrt2}(1, -1)\).

Singular values: \(\sigma_1 = \sqrt{32}\) and \(\sigma_2 = \sqrt{18}\).

Step 3: Strang's attempt at \(U\)

He computed \(AA^\top\):

\[AA^\top = \begin{bmatrix} 32 & 0 \\ 0 & 18 \end{bmatrix}\]
This happens to be diagonal. Its eigenvalues are again 32 and 18, as expected.

He took \(U = I\). Multiplying out gave \(\begin{bmatrix} 4 & 4 \\ 3 & -3 \end{bmatrix}\). The second row has the wrong sign. Strang couldn't find the error on the spot and left it for the next lecture.

The fix

What went wrong. Each eigenvector of \(AA^\top\) is determined only up to sign. The \(u\)'s cannot be chosen independently of the \(v\)'s; they must satisfy \(Av_i = \sigma_iu_i\).

The safe way to get \(U\):

\[u_i = \frac{Av_i}{\sigma_i}\]

For this example:

  • \(Av_1 = \tfrac{1}{\sqrt2}(8, 0)\), so \(u_1 = (1, 0)\).
  • \(Av_2 = \tfrac{1}{\sqrt2}(0, -6)\), so \(u_2 = (0, -1)\).
\[A = \underbrace{\begin{bmatrix} 1 & 0 \\ 0 & -1 \end{bmatrix}}_{U}\begin{bmatrix} \sqrt{32} & 0 \\ 0 & \sqrt{18} \end{bmatrix}\underbrace{\frac{1}{\sqrt2}\begin{bmatrix} 1 & 1 \\ 1 & -1 \end{bmatrix}}_{V^\top}\]

Check. \(\Sigma V^\top = \begin{bmatrix} 4 & 4 \\ 3 & -3 \end{bmatrix}\). Multiplying by \(U\) flips the sign of row 2, which gives exactly \(A\).


5. Example 2: a rank-1 matrix

\[A = \begin{bmatrix} 4 & 3 \\ 8 & 6 \end{bmatrix}\]

The four subspaces

Subspace Contents
Row space multiples of \((4, 3)\)
Nullspace the perpendicular line
Column space multiples of \((1, 2)\)
Left nullspace the perpendicular line

The bases

Each subspace is a line, so each basis vector is just the unit vector along it.

  • \(v_1 = (0.8,\ 0.6)\), from the row space.
  • \(v_2 = (0.6,\ -0.8)\), from the nullspace.
  • \(u_1 = \tfrac{1}{\sqrt5}(1, 2)\), from the column space.
  • \(u_2 = \tfrac{1}{\sqrt5}(2, -1)\), from the left nullspace.

The singular value

\[A^\top A = \begin{bmatrix} 80 & 60 \\ 60 & 45 \end{bmatrix}\]

This has rank 1, so its eigenvalues are \(0\) and the trace, \(125\). Therefore \(\sigma_1 = \sqrt{125} = 5\sqrt5\), and there is no second nonzero singular value.

The decomposition

\[A = \frac{1}{\sqrt5}\begin{bmatrix} 1 & 2 \\ 2 & -1 \end{bmatrix}\begin{bmatrix} \sqrt{125} & 0 \\ 0 & 0 \end{bmatrix}\begin{bmatrix} 0.8 & 0.6 \\ 0.6 & -0.8 \end{bmatrix}\]

Check.

\[Av_1 = (3.2 + 1.8,\ 6.4 + 3.6) = (5, 10) = 5\sqrt5 \cdot \tfrac{1}{\sqrt5}(1, 2) = \sigma_1u_1\]

The nullspace vectors cost nothing. \(v_2\) and \(u_2\) are multiplied by the zero singular value, so they contribute nothing to the product.


6. The big picture: the right bases for the four subspaces

Vectors Orthonormal basis for Dimension
\(v_1, \dots, v_r\) row space \(C(A^\top)\) \(r\)
\(u_1, \dots, u_r\) column space \(C(A)\) \(r\)
\(v_{r+1}, \dots, v_n\) nullspace \(N(A)\) \(n - r\)
\(u_{r+1}, \dots, u_m\) left nullspace \(N(A^\top)\) \(m - r\)

Two properties at once. The bases are orthonormal, which Gram–Schmidt alone could give. They also diagonalize \(A\): each \(Av_i\) points along \(u_i\). That second property requires eigenvalues.

No coupling. In these bases, \(A\) just scales each direction by its own \(\sigma_i\).


Big ideas from Lecture 29

  1. \(A = U\Sigma V^\top\) for every matrix, with orthogonal \(U\), \(V\) and a nonnegative diagonal \(\Sigma\).
  2. Meaning: \(Av_i = \sigma_iu_i\). An orthonormal basis of the row space maps to an orthogonal basis of the column space.
  3. Computation:
    • \(V\) = eigenvectors of \(A^\top A\).
    • \(U\) = eigenvectors of \(AA^\top\).
    • \(\sigma_i = \sqrt{\lambda_i(A^\top A)}\).
  4. Sign consistency: compute \(u_i = Av_i/\sigma_i\) instead of choosing the \(u\)'s independently. This resolves Strang's on-the-board sign problem.
  5. The SVD gives the best bases for all four fundamental subspaces at once.

Connection to your first question. Center a data matrix \(X\) (\(n\) samples by \(d\) features). Then the right singular vectors of \(X\) are the PCA principal directions, and the eigenvalues of the covariance \(C = \frac1nX^\top X\) are \(\sigma_i^2/n\). This is why PCA is computed via the SVD in practice, as mentioned in that answer.

Connection to interpretability. SVD of weight matrices, such as OV and QK circuit matrices, is a standard tool. Keeping the top \(k\) terms \(\sum_{i\le k}\sigma_iu_iv_i^\top\) gives the best rank-\(k\) approximation (Eckart–Young).