The Maximum Variance Perspective

Page 1001 fixed the code length $M$ and left the basis open. §10.2 closes it, from one assumption:

If we interpret information content in the data as how “space filling” the dataset is, then we can describe the information contained in the data by looking at the spread of the data.

That is a modelling choice, not a theorem — and it is worth noticing that it is one, because §10.3 will minimise something apparently different and arrive at exactly the same $\mathbf{B}$.

What you’ll learn

• Equations 10.7–10.9b: how “maximise the variance of the code” collapses to $\mathbf{b}_1^\top\mathbf{S}\mathbf{b}_1$ — one Chapter 6 identity and one line of algebra. Measured, all three readings agree to $8.5\times10^{-14}$.
• Why the constraint $\lVert\mathbf{b}_1\rVert^2 = 1$ is not cosmetic: without it the objective has no maximum at all.
• Equations 10.11–10.14: the Lagrangian, and the moment $\mathbf{S}\mathbf{b}_1 = \lambda_1\mathbf{b}_1$ appears — “the Lagrange multiplier plays the role of the corresponding eigenvalue.”
• Equation 10.15, $V_1 = \lambda_1$ — measured for the first six components, largest gap $8.5\times10^{-14}$. The variance you keep is the eigenvalue; there is nothing else to compute.
• Equation 10.17’s deflation and Equation 10.21’s claim that $\mathbf{S}$ and $\hat{\mathbf{S}}$ share eigenvectors — verified by deflating four times and recovering each $\mathbf{b}_m$ to under $2.1\times10^{-6}$ degrees.
• Equations 10.24 and 10.25: $V_M = \sum_{m\le M}\lambda_m$ and $J_M = \sum_{j>M}\lambda_j$, and the fact that the measured reconstruction error matches $J_M$ to $2.8\times10^{-13}$ — a result §10.3 will spend a whole section deriving.
• Why a random search is hopeless in $64$ dimensions: $200{,}000$ random unit vectors reach $37.18\%$ of $\lambda_1$, while power iteration reaches it in $500$ cheap steps.

Intuition: the objective is a function on a circle

Strip §10.2 to two dimensions and the whole argument fits in one picture. A unit vector in the plane is just an angle. Feed it to $\mathbf{b}^\top\mathbf{S}\mathbf{b}$ and you get a number. Sweep the angle and you get a curve — and the curve has a maximum somewhere.

figure Equation 10.10 is a search over directions, and the answer is an eigenvector matplotlib

The top-left curve is the objective. Its maximum is exactly the largest eigenvalue and its minimum exactly the smallest — for a symmetric S those are the only two stationary points on the circle. The bottom row is the same three numbers as spread: the amber projection is the fattest histogram, the red the thinnest.

Two things the picture makes obvious that the algebra hides:

The curve is smooth and has exactly two stationary points on the half-circle — one maximum, one minimum. Those are the two eigenvectors. In $D$ dimensions there are $D$ of them, which is why one eigendecomposition answers the question for every $M$ at once.

You cannot find the peak by sampling. In two dimensions you could. In $64$ dimensions the surface of the unit sphere is overwhelmingly far from any particular direction: measured below, $200{,}000$ random unit vectors get to $37.18\%$ of the maximum.

§10.2.1 From “spread” to an eigenvalue problem

The variance of the first coordinate of the code, using the i.i.d. assumption:

$$V_1 := \mathbb{V}[z_1] = \frac{1}{N}\sum_{n=1}^{N} z_{1n}^2, \qquad z_{1n} = \mathbf{b}_1^\top\mathbf{x}_n \qquad \text{(10.7), (10.8)}$$

Substituting and using the symmetry of the dot product, $\mathbf{b}_1^\top\mathbf{x}_n = \mathbf{x}_n^\top\mathbf{b}_1$:

$$V_1 = \frac{1}{N}\sum_{n=1}^{N}\mathbf{b}_1^\top\mathbf{x}_n\mathbf{x}_n^\top\mathbf{b}_1 = \mathbf{b}_1^\top\!\left(\frac{1}{N}\sum_{n=1}^{N}\mathbf{x}_n\mathbf{x}_n^\top\right)\!\mathbf{b}_1 = \mathbf{b}_1^\top\mathbf{S}\mathbf{b}_1 \qquad \text{(10.9a), (10.9b)}$$

The only move is pulling $\mathbf{b}_1$ outside the sum, which is legal because it does not depend on $n$. Measured on the digit-”8” data, three routes to the same number:

route	value
$\frac{1}{N}\sum_n z_{1n}^2$, Equation 10.7	$151.568412$
$\mathbf{b}_1^\top\mathbf{S}\mathbf{b}_1$, Equation 10.9b	$151.568412$
np.var of the code	$151.568412$
$\lambda_1$, Equation 10.15	$151.568412$
largest gap between any two	$8.5\times10^{-14}$

(The code’s mean comes out at $2.5\times10^{-16}$, so Equation 10.7’s “mean of the squares” really is the variance, as the book’s centred-data remark promises.)

Why Equation 10.6 lets the whole section assume mean zero

$\mathbb{V}_z[\mathbf{z}] = \mathbb{V}_x[\mathbf{B}^\top(\mathbf{x}-\boldsymbol\mu)] = \mathbb{V}_x[\mathbf{B}^\top\mathbf{x} - \mathbf{B}^\top\boldsymbol\mu] = \mathbb{V}_x[\mathbf{B}^\top\mathbf{x}]$ — a variance is blind to a shift, so subtracting the mean changes nothing about the objective.

This is “without loss of generality” in the strict sense and it is worth separating from page 1001’s measurement. Centring does not change what the objective means; it changes what $\frac{1}{N}\sum\mathbf{x}_n\mathbf{x}_n^\top$ is. If you centre, Equation 10.1 computes the covariance and maximising $\mathbf{b}^\top\mathbf{S}\mathbf{b}$ maximises the variance. If you don’t, it computes something else, and you maximise that instead — which is how the leading direction ended up $0.2851^\circ$ from the mean image.

The constraint is load-bearing

“Arbitrarily increasing the magnitude of the vector $\mathbf{b}_1$ increases $V_1$” — doubling $\mathbf{b}_1$ quadruples $\mathbf{b}_1^\top\mathbf{S}\mathbf{b}_1$. So the unconstrained problem has a supremum of $+\infty$ and no maximiser. Restricting to $\lVert\mathbf{b}_1\rVert^2 = 1$ gives

$$\max_{\mathbf{b}_1}\ \mathbf{b}_1^\top\mathbf{S}\mathbf{b}_1 \quad\text{subject to}\quad \lVert\mathbf{b}_1\rVert^2 = 1 \qquad \text{(10.10)}$$

and Chapter 7’s Lagrangian

$$\mathfrak{L}(\mathbf{b}_1, \lambda) = \mathbf{b}_1^\top\mathbf{S}\mathbf{b}_1 + \lambda_1\!\left(1 - \mathbf{b}_1^\top\mathbf{b}_1\right) \qquad \text{(10.11)}$$ $$\frac{\partial\mathfrak{L}}{\partial\mathbf{b}_1} = 2\mathbf{b}_1^\top\mathbf{S} - 2\lambda_1\mathbf{b}_1^\top, \qquad \frac{\partial\mathfrak{L}}{\partial\lambda_1} = 1 - \mathbf{b}_1^\top\mathbf{b}_1 \qquad \text{(10.12)}$$

Setting both to $\mathbf{0}$:

$$\mathbf{S}\mathbf{b}_1 = \lambda_1\mathbf{b}_1, \qquad \mathbf{b}_1^\top\mathbf{b}_1 = 1 \qquad \text{(10.13), (10.14)}$$

“By comparing this with the definition of an eigenvalue decomposition, we see that $\mathbf{b}_1$ is an eigenvector of the data covariance matrix $\mathbf{S}$, and the Lagrange multiplier $\lambda_1$ plays the role of the corresponding eigenvalue.” Measured — the gradient $2\mathbf{S}\mathbf{b} - 2\lambda_1\mathbf{b}$ at three unit vectors:

unit vector, with $\lambda = \lambda_1$	$\lVert$gradient$\rVert$
$\mathbf{b}_1$	$\mathbf{0.000000}$
$\mathbf{b}_2$	$127.709219$
a random direction (seed $0$)	$285.800195$

Equation 10.15, which is the whole point

Because $\mathbf{S}\mathbf{b}_1 = \lambda_1\mathbf{b}_1$ and $\mathbf{b}_1^\top\mathbf{b}_1 = 1$:

$$V_1 = \mathbf{b}_1^\top\mathbf{S}\mathbf{b}_1 = \lambda_1\mathbf{b}_1^\top\mathbf{b}_1 = \lambda_1 \qquad \text{(10.15)}$$

$m$	$\mathbf{b}_m^\top\mathbf{S}\mathbf{b}_m$	$\lambda_m$	gap
$1$	$151.56841162$	$151.56841162$	$8.5\times10^{-14}$
$2$	$87.71380214$	$87.71380214$	$5.7\times10^{-14}$
$3$	$75.28494722$	$75.28494722$	$1.4\times10^{-14}$
$4$	$54.05513090$	$54.05513090$	$7.1\times10^{-15}$
$5$	$52.41321109$	$52.41321109$	$3.6\times10^{-14}$
$6$	$42.57742529$	$42.57742529$	$0.0$

The eigendecomposition does not just tell you where to look — it tells you what you will find there. No projection needs to be formed to know how much variance a component keeps. The book’s margin note adds the name: $\sqrt{\lambda_1}$ is the loading of $\mathbf{b}_1$, the standard deviation accounted for by $\mathrm{span}[\mathbf{b}_1]$.

Why nobody searches for the maximum, and what they do instead

On the digit-”8” covariance, $200{,}000$ random unit vectors in $\mathbb{R}^{64}$:

method	best $\mathbf{b}^\top\mathbf{S}\mathbf{b}$ found	share of $\lambda_1$
$200{,}000$ random unit vectors	$56.351776$	$37.18\%$
power iteration, $500$ steps	$151.568412$	$100\%$

The second lands $1.2\times10^{-6}$ degrees from $\mathbf{b}_1$ using nothing but repeated $\mathbf{b} \leftarrow \mathbf{S}\mathbf{b}/\lVert\mathbf{S}\mathbf{b}\rVert$ — which is §10.4.2’s Equation 10.52, and the original PageRank algorithm.

§10.2.2 The second component, and why there is nothing new to do

For $\mathbf{b}_m$ with $m > 1$, the book’s device is to remove what the earlier components already explain:

$$\hat{\mathbf{X}} := \mathbf{X} - \sum_{i=1}^{m-1}\mathbf{b}_i\mathbf{b}_i^\top\mathbf{X} = \mathbf{X} - \mathbf{B}_{m-1}\mathbf{X} \qquad \text{(10.17)}$$

then maximise $V_m = \mathbf{b}_m^\top\hat{\mathbf{S}}\mathbf{b}_m$ against the deflated covariance, which is the §10.2.1 problem again. Run literally — deflate, re-decompose, take the top eigenvector, repeat:

$m$	top eigenvalue of $\hat{\mathbf{S}}$	$\lambda_m$ of $\mathbf{S}$	angle to $\mathbf{b}_m$
$1$	$151.568412$	$151.568412$	$1.7\times10^{-6}$ deg
$2$	$87.713802$	$87.713802$	$2.1\times10^{-6}$ deg
$3$	$75.284947$	$75.284947$	$1.2\times10^{-6}$ deg
$4$	$54.055131$	$54.055131$	$0.0$ deg
$5$	$52.413211$	$52.413211$	$0.0$ deg

The book proves this must happen. Split the eigenvectors of $\mathbf{S}$ by whether they have been used:

$$\mathbf{B}_{m-1}\mathbf{b}_i = \mathbf{b}_i \ \ \text{if } i < m, \qquad \mathbf{B}_{m-1}\mathbf{b}_i = \mathbf{0} \ \ \text{if } i \geq m \qquad \text{(10.20)}$$

Expanding $\hat{\mathbf{S}} = \mathbf{S} - \mathbf{S}\mathbf{B}_{m-1} - \mathbf{B}_{m-1}\mathbf{S} + \mathbf{B}_{m-1}\mathbf{S}\mathbf{B}_{m-1}$ and applying each case gives

$$\hat{\mathbf{S}}\mathbf{b}_m = \mathbf{S}\mathbf{b}_m = \lambda_m\mathbf{b}_m \qquad \text{(10.21)}$$ $$\hat{\mathbf{S}}\mathbf{b}_i = \mathbf{0} = 0\,\mathbf{b}_i \quad\text{for } i < m \qquad \text{(10.22)}$$

figure Deflating removes one eigenvalue and leaves every other one alone matplotlib

Left: after removing m components, the deflated spectrum is the original one with its first m entries deleted — the curves overlay exactly, measured to 7.8e-14. Right: the three used-up directions become null vectors of S-hat, with residual norms 1.0e-14, 5.9e-15 and 7.1e-15 against roughly 50 for the untouched ones.

The deflation is a proof, not a procedure

Every eigenvector of $\mathbf{S}$ is still an eigenvector of $\hat{\mathbf{S}}$ — the used ones with eigenvalue $0$, the rest with the eigenvalue they had. So the $m$-th round of the sequential search returns the $m$-th eigenvector of the original $\mathbf{S}$, which one call to np.linalg.eigh already handed you.

Implementing Equation 10.17 literally costs $m$ eigendecompositions to obtain what one produces. Its job in the book is to justify the greedy argument — that taking the top $M$ eigenvectors is optimal, and not merely a sequence of locally optimal choices.

With Equation 10.21 in hand, $V_m = \mathbf{b}_m^\top\mathbf{S}\mathbf{b}_m = \lambda_m$ for every $m$ (Equation 10.23), and summing:

$$V_M = \sum_{m=1}^{M}\lambda_m, \qquad J_M := \sum_{j=M+1}^{D}\lambda_j = V_D - V_M \qquad \text{(10.24), (10.25)}$$

Example 10.2, and the number the next section has to derive

figure The spectrum of the digit-8 covariance, and what it costs to keep most of it matplotlib

The variance is heavily concentrated: one component of sixty-four holds 20.45 percent of it, and eighteen hold ninety percent. The cliff after index fifty-two is page 1001's twelve dead border pixels — twelve directions along which every image is identical, so their eigenvalues fall to machine zero.

$M$	$V_M$ (Eq 10.24)	$J_M$ (Eq 10.25)	share kept	measured $\frac{1}{N}\sum\lVert\mathbf{x}_n-\tilde{\mathbf{x}}_n\rVert^2$
$1$	$151.568412$	$589.590460$	$20.45\%$	$589.590460$
$2$	$239.282214$	$501.876658$	$32.28\%$	$501.876658$
$5$	$421.035503$	$320.123369$	$56.81\%$	$320.123369$
$10$	$574.968069$	$166.190803$	$77.58\%$	$166.190803$
$20$	$681.807067$	$59.351805$	$91.99\%$	$59.351805$
$30$	$722.255016$	$18.903856$	$97.45\%$	$18.903856$
$52$	$741.158872$	$0.000000$	$100.00\%$	$0.000000$

The last column was never mentioned in §10.2. It is the average squared reconstruction error, formed by actually building $\tilde{\mathbf{X}} = \mathbf{B}\mathbf{B}^\top\mathbf{X}$ and measuring — and across every $M$ from $1$ to $63$ it agrees with $J_M$ to $2.8\times10^{-13}$.

That equality is §10.3’s entire result, arrived at here for free. §10.2 maximised the variance kept; the error paid came out as the leftover eigenvalues without anyone asking. §10.3 turns the observation into a derivation, and Equation 10.44 will state it.

variance retained	components needed	share of the $64$ dimensions
$50\%$	$M = 5$	$7.8\%$
$80\%$	$M = 12$	$18.8\%$
$90\%$	$M = 18$	$28.1\%$
$95\%$	$M = 25$	$39.1\%$
$99\%$	$M = 36$	$56.2\%$

Reading the flat part of the spectrum honestly

The book’s Example 10.2 says “only a few of them have a value that differs significantly from 0”, and at MNIST’s $D = 784$ that is a fair description. At $D = 64$ it is less dramatic: $\lambda_1$ is $20.45\%$ of the total and $\lambda_1/\lambda_2 = 1.73$, so there is no single dominant direction. Half the variance still costs $5$ components and $99\%$ costs $36$.

Two lessons that survive both scales. First, the twelve exact zeros are free — they are directions on which the data has no spread at all, and discarding them loses literally nothing. Second, everything after that is a judgement call about how much loss is acceptable, and §10.2 gives you the exact price of every option in one sorted list.

pch.quizTag Does the variance argument hold together?

1. Why does Equation 10.10 constrain b to have unit norm?

   A To make the arithmetic tidier B Because without it the objective is unbounded — doubling b quadruples b-transpose S b, so there is no maximiser at all C Because eigenvectors are always unit vectors D To guarantee S is positive definite
   pch.quizShowAnswer

   B — Because without it the objective is unbounded — doubling b quadruples b-transpose S b, so there is no maximiser at all — The book says it plainly: arbitrarily increasing the magnitude of b increases V1. The supremum of the unconstrained problem is infinity and no vector attains it. The constraint turns a meaningless problem into the one Chapter 7's Lagrangian solves, and the multiplier that appears turns out to be the eigenvalue.

2. Equation 10.15 says V1 equals lambda-1. What is the practical consequence?

   A You still have to project the data to find out how much variance a component keeps B The eigendecomposition already tells you the variance each component keeps — no projection needed C It only holds for the first component D It holds only if the data is Gaussian
   pch.quizShowAnswer

   B — The eigendecomposition already tells you the variance each component keeps — no projection needed — Measured for the first six components of the digit-8 covariance, b-transpose S b matches lambda to at worst 8.5e-14. One sorted list of eigenvalues prices every possible choice of M: V_M is the running sum and J_M is what is left, and the measured reconstruction error matches J_M to 2.8e-13.

3. What does Equation 10.21 establish about the deflated covariance S-hat?

   A S-hat has entirely new eigenvectors that must be recomputed each round B S and S-hat share eigenvectors; the used-up ones simply acquire eigenvalue zero C S-hat is no longer symmetric D S-hat has the same eigenvalues as S
   pch.quizShowAnswer

   B — S and S-hat share eigenvectors; the used-up ones simply acquire eigenvalue zero — Measured: deflating three times leaves S-hat b1, b2 and b3 with norms 1.0e-14, 5.9e-15 and 7.1e-15 — they are null vectors now — while the remaining eigenvalues match those of S to 7.8e-14. So the sequential search returns the eigenvectors of the ORIGINAL S, which one eigendecomposition already produced. Equation 10.17 is a proof device.

4. Why is random search a bad way to find the direction of maximum variance?

   A It is fine; it just takes a few hundred samples B In 64 dimensions the unit sphere is too large — 200,000 random unit vectors reached only 37.18 percent of lambda-1 C Random vectors are never unit vectors D The objective is discontinuous
   pch.quizShowAnswer

   B — In 64 dimensions the unit sphere is too large — 200,000 random unit vectors reached only 37.18 percent of lambda-1 — Measured on the digit-8 covariance: the best of 200,000 random unit vectors scored 56.351776 against lambda-1 of 151.568412. Power iteration — repeatedly multiply by S and renormalise, Equation 10.52 — reaches the exact answer in 500 steps, landing 1.2e-06 degrees from b1.

Exercises

Exercise 1 – Three ways to say the same variance

pch.exerciseTag pch.exerciseTitle

Exercise 2 – The Lagrangian’s gradient, and where it vanishes

pch.exerciseTag pch.exerciseTitle

Exercise 3 – Is lambda-1 really the maximum?

pch.exerciseTag pch.exerciseTitle

Exercise 4 – Deflate, and watch the eigenvalues not move

pch.exerciseTag pch.exerciseTitle

Exercise 5 – Equation 10.25 predicts an error nobody measured yet

pch.exerciseTag pch.exerciseTitle

Recall card

• Section 10.2 assumes that information is spread, and derives PCA by maximising the variance of the code. That is a modelling choice, and Section 10.3 reaches the same answer from a different one.
• Equation 10.9b turns the objective into b-transpose S b in one line, using only the symmetry of the dot product. Measured three ways on the digit-8 data, all agreeing to 8.5e-14.
• The unit-norm constraint is load-bearing: doubling b quadruples the objective, so without it there is no maximum at all.
• The Lagrangian’s stationarity condition is the eigenvalue equation. Measured, the gradient is exactly zero at the leading eigenvector and 127.709219 at the second.
• Equation 10.15: the variance kept equals the eigenvalue. So one sorted eigenvalue list prices every possible choice of M, with no projection formed.
• Two hundred thousand random unit vectors in 64 dimensions reach 37.18 percent of the largest eigenvalue. Power iteration reaches it exactly in 500 steps.
• Equation 10.17’s deflation is a proof device, not an algorithm. S and the deflated S-hat share eigenvectors; the used-up ones simply acquire eigenvalue zero, measured at about 1e-14.
• So the greedy sequential choice is genuinely optimal, not merely locally optimal — the m-th round returns the m-th eigenvector of the original S.
• Equations 10.24 and 10.25 split the total variance into what is kept and what is lost, and the two always sum to the trace.
• The lost part equals the measured reconstruction error, to 2.8e-13, even though Section 10.2 never mentioned reconstruction. Section 10.3 derives that.
• On 174 images of the digit eight, half the variance costs 5 components and ninety-nine percent costs 36. Twelve eigenvalues are machine zero — the dead border pixels, free to discard.

Next: The Projection Perspective — the same $\mathbf{B}$, derived by minimising the reconstruction error instead, starting with the coordinates rather than the basis.