PCA is the most widely used dimensionality reduction technique: a linear rotation of the data onto the directions of greatest variance. This page covers the intuition, the eigen/SVD math, how to choose the number of components, and how to use it correctly in a modeling pipeline.
Principal Component Analysis finds a new orthogonal coordinate system for the data in which the first axis captures the most variance, the second the most remaining variance, and so on. Keeping the first k axes gives the best rank-k linear approximation of the centered data in the least-squares sense. It is the default linear method for dimensionality reduction, denoising, decorrelating features, and getting a first look at high-dimensional data.
Center the data matrix X (n × p) by subtracting column means. The sample covariance is Σ = XTX / (n - 1). Its eigenvectors are the principal directions and its eigenvalues λj the variance along each. The first component solves:
w1 = argmax||w||=1 wTΣw
In practice PCA is computed from the singular value decomposition X = U S VT, which is more numerically stable than forming Σ. The rows of VT are the components, λj = Sj2 / (n - 1), and the projected scores are Z = X Vk = Uk Sk. For large data, randomized SVD computes only the top k components at a fraction of the cost.
Common criteria: a cumulative explained-variance target, the "elbow" in the scree plot, or simply the k that maximizes downstream model performance under cross-validation. scikit-learn accepts a float, e.g. PCA(n_components=0.95), to keep enough components for 95% of variance.
import numpy as np
from sklearn.datasets import load_breast_cancer
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
from sklearn.pipeline import make_pipeline
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import cross_val_score
X, y = load_breast_cancer(return_X_y=True) # 569 x 30, correlated features
pca = make_pipeline(StandardScaler(), PCA(n_components=0.95)).fit(X)
p = pca.named_steps["pca"]
print("components kept:", p.n_components_)
print("variance ratio:", p.explained_variance_ratio_.round(3)[:5])
print("top loadings PC1:", np.argsort(-np.abs(p.components_[0]))[:5])
# Reconstruction error as an anomaly score
Z = pca.transform(X)
X_hat = pca.inverse_transform(Z)
err = ((StandardScaler().fit_transform(X) - pca.named_steps["standardscaler"].transform(X_hat)) ** 2).sum(1)
# PCA inside a model pipeline: fit only on training folds
clf = make_pipeline(StandardScaler(), PCA(n_components=10), LogisticRegression(max_iter=1000))
print("CV accuracy:", cross_val_score(clf, X, y, cv=5).mean().round(3))Putting PCA inside a Pipeline matters: fitting it on the full dataset before cross-validation leaks test-fold information into the components.
TruncatedSVD, which does not center and so keeps the matrix sparse (this is LSA on TF-IDF).Take five points in two features, already centred: (2, 1.8), (1, 1.1), (0, 0), (−1, −0.9), (−2, −2). The sample covariance matrix, using n − 1 = 4, is approximately
Σ = [[2.50, 2.40], [2.40, 2.32]]
Its eigenvalues are about 4.81 and 0.006. The first eigenvector points almost exactly along the diagonal, roughly (0.72, 0.69), and the second is perpendicular to it. PC1 therefore explains 4.81 / 4.815 ≈ 99.9% of the variance: the two features are almost the same measurement, and one number per point (its projection onto the diagonal) captures nearly everything. The projected score for the first point is 0.72·2 + 0.69·1.8 ≈ 2.69.
This is the typical real-world case in miniature. Sensor channels, survey items that ask the same question twice, or pixel neighbours are strongly correlated, so a few components capture most of the variation and the rest is mostly noise.
Projecting onto the top k components and mapping back gives the best rank-k reconstruction in squared error (the Eckart–Young theorem). The residual, ||x − x̂||2, equals the variance in the discarded directions for that point. Points that follow the dominant correlation structure reconstruct well; points that break it reconstruct poorly. That residual is a classic score in process monitoring (often called the Q or SPE statistic) and in anomaly detection, alongside Hotelling's T2, which measures how far a point lies inside the retained subspace.
PCA is exactly a linear autoencoder with tied weights: a linear network trained to minimise reconstruction error spans the same subspace as the top components. Nonlinear autoencoders generalise the idea to curved manifolds at the cost of training and tuning. Dropping the trailing components also denoises: in image patches or spectra, the small-eigenvalue directions are dominated by noise, so reconstructing from the top components smooths it out.
In scikit-learn, pca.inverse_transform(pca.transform(X)) gives the reconstruction; for whitening (unit variance in every component, useful before some algorithms such as ICA) set whiten=True, but be aware that whitening also amplifies the noisy low-variance components you keep.
| Method | Difference from PCA | Use when |
|---|---|---|
| TruncatedSVD / LSA | No centring; works on sparse matrices | TF-IDF text, very wide sparse data |
| Incremental PCA | Fits in mini-batches | Data does not fit in memory |
| Kernel PCA | PCA in an implicit feature space via a kernel | Moderate n with nonlinear structure |
| Sparse PCA | L1 penalty gives loadings with many zeros | You need interpretable components |
| Factor Analysis | Models per-feature noise variance separately | Latent-variable interpretation (psychometrics) |
| ICA | Seeks statistically independent, not just uncorrelated, components | Source separation (EEG, audio) |
| LDA (linear discriminant) | Supervised: maximises between-class separation | Labels exist and separation matters |
| t-SNE / UMAP | Nonlinear, neighbour-preserving | 2-D visualisation |