Capstone: PCA on real handwritten digits

Ten thousand images of 784 pixels each, squeezed to two numbers apiece. Two components keep a small share of the variance, and the digits already sort themselves.

11 min Prerequisites: Module 3 (SVD and PCA)

Symbols in this lessonℝⁿ real vector / matrix space x̄ bar ᵀ transpose rank rank

1. The idea

The third capstone runs PCA on real data: the 10,000 test images of MNIST, the handwritten digits that Core Course chapter 14 learns to read. Each image is 28 × 28 = 784 pixels, so each one is a point in ℝ⁷⁸⁴. PCA finds the directions in that space along which the images vary most.

The recipe is Module 3's: centre every pixel, take the SVD of the centred data, and read off the components. The first two together keep about 18% of the pixel variance, far from all of it. Yet project the images onto them and the digits already form separate clouds: the 1s, mostly a single upright stroke, sit far from the round 0s.

That is the usual job PCA does: not keeping everything, but finding a few directions that show the structure. The same steps are used to take a quick look at embeddings, to compress features and to remove noise.

2. Try it: digits as points

Each mark is one test image, placed by its first two principal components: 60 each of the digits 0, 1, 4 and 7. Turn digits on and off to see which overlap.

3. The math

Each image is a vector

\mathbf{x}^{(i)} \in \mathbb{R}^{784}, \qquad X \in \mathbb{R}^{10000 \times 784}

One row per image, one column per pixel.

Centre every pixel

\tilde{X} = X - \bar{\mathbf{x}}

x̄ is the average image.

The SVD gives the components

\tilde{X} = U\Sigma V^\top, \qquad \text{component } k = \text{row } k \text{ of } V^\top

Each component is itself a 784-pixel pattern.

Share of the variance

\frac{\sigma_k^2}{\sum_{j}\sigma_j^2}

From the singular values.

Two numbers per image

Z = \tilde{X}\,V_{:,\,1:2} \in \mathbb{R}^{10000 \times 2}

What the widget plots.

4. The code (python)

mm7_pca.py

# MNIST is by Yann LeCun, Corinna Cortes and Christopher Burges, licensed CC BY-SA 3.0.
# This downloads the two test files (about 1.7 MB) from our copy and checks they are the right ones.
import gzip
import hashlib
import os
import urllib.request
import numpy as np

URL = "https://epochzero.org/datasets/mnist/"
FILES = {"t10k-images-idx3-ubyte.gz": "8d422c7b0a1c1c79245a5bcf07fe86e33eeafee792b84584aec276f5a2dbc4e6",
         "t10k-labels-idx1-ubyte.gz": "f7ae60f92e00ec6debd23a6088c31dbd2371eca3ffa0defaefb259924204aec6"}
for name, sha in FILES.items():
    if not os.path.exists(name):
        urllib.request.urlretrieve(URL + name, name)
    with open(name, "rb") as f:
        assert hashlib.sha256(f.read()).hexdigest() == sha, name + " is not the file we expect"

with gzip.open("t10k-images-idx3-ubyte.gz") as f:
    X = np.frombuffer(f.read(), np.uint8, offset=16).reshape(-1, 784) / 255.0
with gzip.open("t10k-labels-idx1-ubyte.gz") as f:
    y = np.frombuffer(f.read(), np.uint8, offset=8)
print(X.shape, y.shape)                          # (10000, 784) (10000,)

Xc = X - X.mean(axis=0)                          # centre every pixel
U, s, Vt = np.linalg.svd(Xc, full_matrices=False)
evr = s**2 / np.sum(s**2)                        # share of the variance, per component
print(evr[:2].round(3), round(evr[:50].sum(), 3))   # [0.1   0.075] 0.832
Z = Xc @ Vt[:2].T                                # every image as two numbers
print(Z.shape)                                   # (10000, 2)
print({d: Z[y == d].mean(axis=0).round(2).tolist() for d in (0, 1, 4, 7)})   # {0: [-3.92, 0.45], 1: [3.31, 1.75], 4: [0.42, -2.32], 7: [1.29, -2.11]}

5. You’ll see this in

6. Check yourself

Answer all 4 questions correctly to complete the chapter · 0 / 4 done

Q1/4 Seen as a vector, where does one 28 × 28 MNIST image live?

Q2/4 What does PC1's share of the variance tell you?

Q3/4 Why centre the pixels before the SVD?

Q4/4 What shape is Z for 10,000 images and 2 components?

Finished this chapter?

Weekly newsletter

What changed in ML architectures this week

New models worth knowing, benchmark moves, and one architecture explained in plain words. Free, once a week, unsubscribe anytime.

Learn machine learning from scratch