PCA in Practice

Running PCA on real data requires a handful of practical steps that textbook treatments sometimes skip. First, center the data by subtracting the mean: $\tilde{\mathbf{x}}_n = \mathbf{x}_n - \bar{\mathbf{x}}$. Without centering, the first 'principal component' will simply point toward the mean, not toward the direction of real variation. Second, usually standardize (divide each feature by its standard deviation) when features live on different scales — otherwise PCA is dominated by whichever feature happens to have the largest range.

How many components to keep? The scree plot shows eigenvalues in decreasing order. Look for an 'elbow' where eigenvalues flatten — components beyond the elbow are mostly noise. A more principled alternative is the cumulative explained variance: choose $k$ so that $\sum_{i=1}^{k} \lambda_i / \sum_j \lambda_j \geq 0.95$ (or whatever threshold the downstream task requires). For purely visualization purposes, $k = 2$ or $3$ is fixed by necessity.

Computationally, avoid forming $S = X^\top X$ explicitly when $D$ is large — it is $D \times D$ and may not fit in memory. Instead, use SVD of $X$ directly: $X = U \Sigma V^\top$, where $V$ columns are principal components and $\Sigma^2 / N$ are the eigenvalues of $S$. For huge datasets, randomized SVD (Halko et al.) computes the top-$k$ components in $O(NDk)$ time — orders of magnitude faster than full SVD.

PCA has important limitations. It is linear — it cannot capture curved manifolds (use kernel PCA, Isomap, UMAP, or autoencoders for that). It assumes variance = signal, which fails when noise has large variance (consider ICA instead). It is not scale-invariant — standardize carefully. And it is a lossy, unsupervised method; discarded components may contain exactly the features your downstream task needs. Always validate PCA choices against the downstream metric, not just explained variance.

Despite these caveats, PCA is the first tool to try for any high-dimensional dataset. Canonical applications: face recognition (eigenfaces), finance (principal portfolios / risk factors), genomics (population structure), NLP (latent semantic analysis is PCA on word-document matrices), image compression (JPEG's DCT is cousin to PCA). The machine-learning lesson is deeper than the algorithm: covariance matrices encode structure, and eigen-decomposition reveals it.

Worked example — rank-1 reconstruction error: A $2\times 2$ image $X = \begin{pmatrix} 2 & 2 \\ 2 & 2 \end{pmatrix}$ has SVD $\Sigma = \text{diag}(4, 0)$, so the rank-1 approximation is exact — zero reconstruction error. Perturb to $X = \begin{pmatrix} 2 & 2 \\ 2 & 3 \end{pmatrix}$: singular values are $\sigma_1 \approx 4.56$, $\sigma_2 \approx 0.44$. Keeping only $\sigma_1$ gives Frobenius reconstruction error $\sigma_2 = 0.44$, while the total 'energy' is $\sqrt{\sigma_1^2 + \sigma_2^2} \approx 4.58$ — the top component explains $\sigma_1^2/(\sigma_1^2 + \sigma_2^2) \approx 99.1\%$ of variance.

Worked example — choosing $k$ for 95% variance: Eigenvalues $[10, 6, 2, 1, 1]$, total $20$. Cumulative ratios: $k=1 \to 0.5$; $k=2 \to 0.8$; $k=3 \to 0.9$; $k=4 \to 0.95$. So $k=4$ is the smallest number of components to retain 95% of the variance. The scree plot would show a clear elbow after the second component — a judgment call between aggressive compression ($k=2$, 80%) and conservative preservation ($k=4$, 95%).