Applications: Image Processing and Statistics

Principal component analysis applies the singular value decomposition and the Spectral Theorem to multivariate data: a collection of objects, each described by several measurements, so each datum is a vector in $R^{p}$ . PCA finds new coordinates in which the measurements are uncorrelated and sorted by decreasing variance, then discards the coordinates with little variance.

For example, a chemical plant might test each of $300$ plastic samples on eight properties: melting point, density, tensile strength, and so on. Each report is a vector in $R^{8}$ , and the set forms an $8 \times 300$ matrix of observations. A satellite image measured in three spectral bands assigns each pixel a vector in $R^{3}$ ; a two-megapixel image is then a $3 \times (4, 000, 000)$ matrix. The dimension that PCA reduces is the number of measurements per object, not the number of objects.

Mean, deviation, and covariance

Let $[X_{1} \dots X_{N}]$ be a $p \times N$ matrix of observations. The sample mean is the average vector

M = \frac{1}{N} (X_{1} + \dots + X_{N}),

the center of the data cloud. Subtracting it recenters every observation: $\hat{X}_{k} = X_{k} - M$ , and the matrix $B = [\hat{X}_{1} \dots \hat{X}_{N}]$ is said to be in mean-deviation form, with zero mean.

Because $S$ has the form $\frac{1}{N - 1} B B^{⊤}$ , it is symmetric and positive semidefinite: $x^{⊤} S x = \frac{1}{N - 1} ∥ B^{⊤} x ∥^{2} \geq 0$ . Its entries carry the spread and coupling of the measurements:

Variances on the diagonal. The entry $s_{j j}$ is the variance of the $j$ -th measurement — how widely the $j$ -th coordinate ranges over the data.
Covariances off the diagonal. The entry $s_{ij}$ ( $i \neq = j$ ) is the covariance of measurements $i$ and $j$ . When $s_{ij} = 0$ the two are uncorrelated.
Total variance is the trace $tr (S) = s_{11} + \dots + s_{pp}$ , the sum of the individual variances.

Analysis is easiest when the covariance matrix is diagonal, meaning all measurements are mutually uncorrelated; PCA produces exactly this.

Principal component analysis

PCA looks for an orthogonal change of variable that decorrelates the data. With the observations in mean-deviation form, find an orthogonal $P = [u_{1} \dots u_{p}]$ giving new coordinates $X = P Y$ , equivalently $Y = P^{⊤} X$ , so that the new variables $y_{1}, \dots, y_{p}$ are uncorrelated and ordered by decreasing variance. The covariance matrix of the transformed data is $P^{⊤} S P$ , and making it diagonal is orthogonal diagonalization of $S$ .

Since $S$ is symmetric, the Spectral Theorem supplies the diagonalizing $P$ . The eigenvalue $λ_{j}$ is the variance of the new variable $y_{j}$ , and the $j$ -th principal component defines that variable as a weighted combination of the originals: if $u_{1}$ has entries $c_{1}, \dots, c_{p}$ , then

y_{1} = u_{1}^{⊤} X = c_{1} x_{1} + c_{2} x_{2} + \dots + c_{p} x_{p} .

The weights are the entries of the eigenvector, and $y_{1}$ has the largest variance any unit-weighted combination can achieve. This last claim is the extreme-values-on-the-unit-sphere theorem of the constrained optimization lesson: for a unit vector $u$ , the variance of $y = u^{⊤} X$ is $u^{⊤} S u$ , whose maximum over unit $u$ is the top eigenvalue $λ_{1}$ , attained at $u_{1}$ . Each later component maximizes variance subject to being uncorrelated with the ones before.

The first principal component u₁ lies along the long axis of the data cloud — the direction of greatest variance — and u₂ is perpendicular; the mean M is the center from which deviations are measured.

Multispectral imaging

Reducing dimension

An orthogonal change of variable preserves total variance, because $P$ does not change lengths: $tr (D) = tr (P^{⊤} S P) = tr (S)$ . So the variance is redistributed across the new variables,

new variables λ_{1} + \dots + λ_{p} = tr (S) = original variables s_{11} + \dots + s_{pp},

and the fraction of total variance captured by $y_{j}$ is $λ_{j} / tr (S)$ . When a few eigenvalues dominate, the corresponding components describe most of the data.

For the Railroad Valley image, the total variance is $tr (S) = 7614.23 + 427.63 + 98.10 = 8139.96$ , split as

Component	Eigenvalue	Fraction of variance
First	$7614.23$	$93.5%$
Second	$427.63$	$5.3%$
Third	$98.10$	$1.2%$

The first component alone holds $93.5%$ of the variation, so a single grayscale image displays almost all the information the three bands carry. The third component is nearly constant across the scene; the data effectively lie in a plane, and largely along one line within it. Three-dimensional data has been compressed to essentially one dimension with a $1.2%$ loss.

The variance spectrum of the covariance matrix falls off sharply; the cumulative curve shows the first component already explains most of the total, so truncating after it discards little.

The SVD route

Diagonalizing $S$ directly is not the preferred computation. Forming $S = \frac{1}{N - 1} B B^{⊤}$ squares the data, magnifying errors, the same hazard noted for $A^{⊤} A$ in the SVD lesson. Instead, run the SVD on the data matrix itself. Set

A = \frac{1}{N - 1} B^{⊤},

an $N \times p$ matrix. Then $A^{⊤} A = \frac{1}{N - 1} B B^{⊤} = S$ , so:

The squares of the singular values of $A$ are the eigenvalues of $S$ , hence the variances $λ_{j}$ of the principal components.
The right singular vectors of $A$ are the eigenvectors of $A^{⊤} A = S$ , hence the principal components themselves.

Computing the SVD of $A$ iteratively is faster and more accurate than an eigenvalue decomposition of $S$ , and the gap widens with dimension. For a hyperspectral image with $p = 224$ bands, PCA runs in seconds this way. The numerical SVD lesson develops the algorithm.

Low-rank approximation and compression

The rank-one expansion $A = \sum_{i = 1}^{r} σ_{i} u_{i} v_{i}^{⊤}$ orders the pieces of $A$ by singular value, so truncating after $k$ terms,

A_{k} = σ_{1} u_{1} v_{1}^{⊤} + \dots + σ_{k} u_{k} v_{k}^{⊤},

gives the best rank- $k$ approximation of $A$ : no other rank- $k$ matrix is closer in either the operator or Frobenius norm. Image compression uses this truncation directly. A grayscale image is an $m \times n$ matrix of pixel intensities; storing it in full takes $mn$ numbers. Storing $A_{k}$ takes only $k (m + n + 1)$ : the $k$ singular values, plus $k$ vectors of length $m$ and $k$ of length $n$ . When the singular values decay quickly, a small $k$ reconstructs the image with little visible loss at a fraction of the storage.

A rank-k SVD truncation stores k(m+n+1) numbers instead of the full mn; when the singular values decay fast, small k reconstructs the image with little loss.

The choice of $k$ is read off the singular value spectrum: keep enough terms to capture a target fraction of the total energy $\sum σ_{i}^{2}$ , and drop the tail where the singular values are negligible. The same truncation denoises, since noise usually spreads across the small singular values that get discarded.

The ratio $σ_{1} / σ_{r}$ of largest to smallest retained singular value is the condition number, a measure of how close the matrix sits to rank-deficiency and how sensitive its inverse problems are. PCA itself is equivalent to orthogonal regression: the line through the mean along $u_{1}$ minimizes the sum of squared perpendicular distances to the data.

In the u₁ coordinate the first principal component is the horizontal axis; PCA places it to minimize the total squared perpendicular distance — the vertical drops — from the data to the line.