Principal Component Analysis (PCA): Theory and Practice

Datasets often contain dozens or thousands of correlated features. Many of those dimensions are redundant: height in centimetres and height in inches; the pixel values of neighbouring pixels. Principal Component Analysis (Pearson 1901, Hotelling 1933) finds a new coordinate system aligned with the directions along which the data varies most, letting us compress, visualise and denoise data. It is the most widely used dimensionality-reduction technique in science.

Two equivalent goals#

Let the data matrix $\mathbf{X} \in \mathbb{R}^{n \times d}$ be centred (each column has mean zero). We seek a unit vector $\mathbf{u}$ onto which to project the data.

1. Maximise variance#

The variance of the projections $\mathbf{X}\mathbf{u}$ is $\mathbf{u}^\top\mathbf{S}\mathbf{u}$, where $\mathbf{S} = \frac{1}{n-1}\mathbf{X}^\top\mathbf{X}$ is the sample covariance matrix. Maximise it subject to $\|\mathbf{u}\| = 1$ using a Lagrange multiplier:

So $\mathbf{u}$ must be an eigenvector of the covariance matrix, and the variance along it equals the eigenvalue $\lambda$. The first principal component is the eigenvector with the largest eigenvalue; the second is the next eigenvector (orthogonal to the first), and so on.

2. Minimise reconstruction error#

Project onto a $k$-dimensional subspace spanned by orthonormal $\mathbf{U}_k$ and reconstruct: $\hat{\mathbf{x}} = \mathbf{U}_k\mathbf{U}_k^\top\mathbf{x}$. The subspace minimising the average squared reconstruction error $\frac{1}{n}\sum_i\|\mathbf{x}_i - \hat{\mathbf{x}}_i\|^2$ is spanned by the top $k$ eigenvectors. Since total variance = retained variance + reconstruction error, maximising one minimises the other. (This is the Eckart–Young theorem in disguise.)

Computing PCA with the SVD#

In practice we never form $\mathbf{X}^\top\mathbf{X}$ explicitly. Take the SVD of the centred data:

• The columns of $\mathbf{V}$ are the principal directions (loadings).
• The eigenvalues are $\lambda_j = \sigma_j^2/(n - 1)$.
• The projected coordinates (scores) are $\mathbf{X}\mathbf{V}_k = \mathbf{U}_k\boldsymbol{\Sigma}_k$.

import numpy as np
from sklearn.datasets import load_digits

X, y = load_digits(return_X_y=True)              # 1797 images, 64 pixels
Xc = X - X.mean(0)
U, s, Vt = np.linalg.svd(Xc, full_matrices=False)
explained = s**2 / (s**2).sum()
cum = np.cumsum(explained)
print("variance explained by first 2 PCs:", explained[:2].round(3))
print("components needed for 90% variance:", int(np.searchsorted(cum, 0.90) + 1))

Z = Xc @ Vt[:2].T                                 # 2-D scores for plotting
X_rec = (Xc @ Vt[:20].T) @ Vt[:20] + X.mean(0)    # reconstruct from 20 components
print("relative reconstruction error (20 PCs):",
      round(np.linalg.norm(X - X_rec) / np.linalg.norm(X), 3))

The 64-dimensional digits need only around 20–30 components to retain 90% of the variance.

Choosing the number of components#

• Cumulative explained variance — keep enough components for, say, 90–95%.
• Scree plot — plot eigenvalues and look for an elbow.
• Downstream performance — choose $k$ by cross-validating the model that uses the PCA features.
• Kaiser rule (for standardised data) — keep components with eigenvalue > 1; crude but common in social science.

Standardise or not?#

PCA is sensitive to scale: a feature measured in large units dominates the variance. If features have different units, standardise them (PCA on the correlation matrix). If all features share meaningful units (pixel intensities, gene expression on the same scale), centring alone may be appropriate.

Interpreting components#

Each component is a weighted combination of original features. Inspecting the loadings can reveal latent factors — in a survey, one component might load on all income-related questions ("economic status"). But be careful: components are defined by variance, not meaning, and their signs are arbitrary.

Uses#

• Visualisation — project to 2-D or 3-D.
• Compression and speed — fewer features for downstream models.
• Noise reduction — discard low-variance components that mostly contain noise.
• Decorrelation / whitening — transform features to be uncorrelated with unit variance.
• Eigenfaces — the classic face-recognition representation.
• Population genetics — top PCs of genotype data mirror geography.
• Analysing neural network representations and embeddings.

Variants#

• Incremental PCA — for data that does not fit in memory.
• Randomised PCA — fast approximate SVD for large matrices.
• Kernel PCA — non-linear PCA via kernels.
• Sparse PCA — loadings with many zeros for interpretability.
• Probabilistic PCA / Factor analysis — generative latent-variable versions.