Dimensionality Reduction (PCA, t-SNE)

The Curse of Dimensionality

As the number of features (dimensions) grows, data becomes sparse: the volume of the space increases exponentially, so a fixed number of data points covers less and less of the space. Consequences:

• Distance-based methods (KNN, k-means, SVM-RBF) lose discriminative power — all pairwise distances converge
• More parameters to estimate → more risk of overfitting
• Visualization becomes impossible above 3 dimensions
• Computation cost scales poorly

Dimensionality reduction projects data from high-D to low-D, ideally preserving the structure that matters for downstream tasks.

Two Camps

Principal Component Analysis (PCA)

PCA finds the directions of maximum variance in the data — the principal components — and projects data onto those directions. It is a linear, unsupervised method.

Intuition

Imagine a cloud of 3D data points shaped like a flat oval disk. PCA finds the plane of the disk (first 2 PCs) and the thin direction (third PC). Projecting onto the plane loses only the thin direction, retaining most variance.

Algorithm

1. Center the data: X̃ = X − μ (subtract column means)
2. Compute covariance matrix: C = (1/(n−1)) X̃ᵀX̃ (shape: d×d)
3. Eigendecompose C: C = V Λ Vᵀ (V: eigenvectors, Λ: diagonal eigenvalues)
4. Sort eigenvectors by eigenvalue (descending)
5. Project: Z = X̃ V_k where V_k contains the top-k eigenvectors

Equivalently, via SVD: X̃ = U S Vᵀ → principal components = columns of V, scores = U·S.

Explained Variance Ratio

Explained variance ratio for component i:

EVR_i = λᵢ / Σⱼ λⱼ

Choose k such that cumulative EVR ≥ 90% (or 95%, domain-dependent).

from sklearn.decomposition import PCA
import numpy as np
import matplotlib.pyplot as plt

pca = PCA()
pca.fit(X_train)

# Scree plot
plt.plot(np.cumsum(pca.explained_variance_ratio_))
plt.xlabel("Number of components")
plt.ylabel("Cumulative explained variance")
plt.axhline(0.95, color="red", linestyle="--")
plt.show()

# Choose k=50
pca_k = PCA(n_components=50, random_state=42)
X_train_pca = pca_k.fit_transform(X_train)
X_test_pca  = pca_k.transform(X_test)   # use same transform!

print(f"Explained variance: {pca_k.explained_variance_ratio_.sum():.2%}")

Properties of PCA

• Components are orthogonal (uncorrelated)
• First component has maximum variance; each subsequent component has maximum remaining variance orthogonal to previous ones
• PCA is sensitive to feature scale → always standardize before PCA
• PCA can whiten data (divide by std of each component): useful for neural nets

When PCA Works Well vs. Poorly

Works well: redundant features, strong linear correlations, noise reduction.

Fails: when important structure is nonlinear (spirals, clusters not aligned with variance directions), when dimensions with high variance are not the discriminative ones for classification.

Randomized PCA

For large matrices (millions of rows), exact SVD is slow. Randomized PCA uses a sketch:

pca_rand = PCA(n_components=50, svd_solver="randomized", random_state=42)
X_reduced = pca_rand.fit_transform(X)   # much faster for large n

Incremental PCA

For datasets too large to fit in RAM:

from sklearn.decomposition import IncrementalPCA

ipca = IncrementalPCA(n_components=50, batch_size=256)
for batch in np.array_split(X, 100):
    ipca.partial_fit(batch)
X_reduced = ipca.transform(X)

Linear Discriminant Analysis (LDA)

LDA is a supervised linear dimensionality reduction method. Instead of maximizing variance (PCA), it maximizes the ratio of between-class to within-class scatter — finding projections that separate classes.

Objective: maximize

J(w) = (wᵀ S_B w) / (wᵀ S_W w)

• S_B = between-class scatter matrix: Σₖ nₖ (μₖ − μ)(μₖ − μ)ᵀ
• S_W = within-class scatter matrix: Σₖ Σᵢ∈Cₖ (xᵢ − μₖ)(xᵢ − μₖ)ᵀ

LDA produces at most K−1 discriminant components (K = number of classes).

from sklearn.discriminant_analysis import LinearDiscriminantAnalysis

lda = LinearDiscriminantAnalysis(n_components=2)
X_lda = lda.fit_transform(X_train, y_train)   # supervised — needs y

Use: classification pre-processing, better class separation than PCA when labels available.

t-SNE (t-Distributed Stochastic Neighbor Embedding)

t-SNE is a nonlinear method designed primarily for visualization (2D or 3D). It preserves local structure — nearby points in high-D stay nearby in low-D.

How It Works

1. Compute pairwise similarities in high-D using Gaussian kernel: p_{ij} ∝ exp(−‖xᵢ−xⱼ‖²/2σ²)
2. In low-D space, use Student t-distribution (heavier tails) for similarity: q_{ij} ∝ (1+‖yᵢ−yⱼ‖²)⁻¹
3. Minimize KL divergence between p and q via gradient descent

The t-distribution prevents crowding: distant points can spread out freely.

Key Hyperparameter: Perplexity

Perplexity (5–50, default 30): controls the effective number of neighbors. Low perplexity = focus on very local structure; high perplexity = more global.

from sklearn.manifold import TSNE

tsne = TSNE(
    n_components=2,
    perplexity=30,
    learning_rate="auto",
    n_iter=1000,
    random_state=42
)
X_2d = tsne.fit_transform(X)   # no transform() — must refit for new data

import matplotlib.pyplot as plt
plt.scatter(X_2d[:, 0], X_2d[:, 1], c=y, cmap="tab10", alpha=0.5, s=5)
plt.title("t-SNE visualization")
plt.show()

t-SNE Caveats

• Cannot transform new points — must rerun for any new data
• Does NOT preserve global distances — cluster sizes and inter-cluster distances in the plot are meaningless
• Stochastic — two runs with different seeds give different plots
• Slow: O(n²) memory/time (BH-t-SNE is O(n log n) via approximation)
• Multiple scales: run at several perplexity values; structure consistent across perplexity is real

UMAP (Uniform Manifold Approximation and Projection)

UMAP is faster than t-SNE and preserves both local and global structure better. Based on Riemannian geometry and topological data analysis.

Key differences vs. t-SNE:

import umap   # pip install umap-learn

reducer = umap.UMAP(
    n_components=2,
    n_neighbors=15,     # local neighborhood size (like perplexity)
    min_dist=0.1,       # minimum distance between embedded points
    metric="euclidean",
    random_state=42
)
X_umap = reducer.fit_transform(X)   # or fit then transform separately

# Transform new points (unlike t-SNE):
X_new_umap = reducer.transform(X_new)

Autoencoders (Neural Dimensionality Reduction)

An autoencoder learns to compress data to a low-dimensional latent code and reconstruct it:

Input x → Encoder f_θ(x) → Latent z → Decoder g_φ(z) → Reconstruction x̂

Loss: reconstruction error = ‖x − x̂‖²

The bottleneck (latent code dimension) forces the encoder to learn a compact representation.

import torch
import torch.nn as nn

class Autoencoder(nn.Module):
    def __init__(self, input_dim, latent_dim):
        super().__init__()
        self.encoder = nn.Sequential(
            nn.Linear(input_dim, 256), nn.ReLU(),
            nn.Linear(256, 64),        nn.ReLU(),
            nn.Linear(64, latent_dim)
        )
        self.decoder = nn.Sequential(
            nn.Linear(latent_dim, 64), nn.ReLU(),
            nn.Linear(64, 256),        nn.ReLU(),
            nn.Linear(256, input_dim)
        )

    def forward(self, x):
        z = self.encoder(x)
        return self.decoder(z), z

Variants: - Variational Autoencoder (VAE): latent space is a distribution N(μ, σ²); generates new samples - Denoising Autoencoder: reconstruct clean data from corrupted input (robust features) - Sparse Autoencoder: enforce sparsity in latent code (interpretable features)

Truncated SVD / LSA

For sparse matrices (like TF-IDF text vectors), use Truncated SVD (sklearn's TruncatedSVD) instead of PCA — it works directly on sparse format without centering:

from sklearn.decomposition import TruncatedSVD

svd = TruncatedSVD(n_components=100, random_state=42)
X_lsa = svd.fit_transform(X_tfidf)  # X_tfidf is sparse

This is the basis of Latent Semantic Analysis (LSA) in NLP.

Choosing a Dimensionality Reduction Method

Non-Negative Matrix Factorization (NMF)

Factorizes a non-negative matrix V ≈ WH where W and H are also non-negative. Unlike PCA, components are additive parts (e.g., face parts, topics). Common in topic modeling and image decomposition.

from sklearn.decomposition import NMF

nmf = NMF(n_components=20, random_state=42)
W = nmf.fit_transform(X)   # document-topic matrix
H = nmf.components_         # topic-word matrix