Machine Learning

Imperial College London

Nuri Cingillioglu

Part 6: Unsupervised Learning

The art of latent

Where are we?

In supervised learning, every training sample came with an explicit teacher target:

  • Regression: Predict continuous target (e.g. server power, house price).
  • Classification: Predict discrete class (e.g. benign vs malignant).

Human labels are expensive, noisy, and rare. Unlabelled data is everywhere: millions of unlabelled images, streaming sensor telemetry, genomes, and raw web text.

What about self-supervised learning?

Enter Unsupervised Learning

In unsupervised learning, we are handed only feature vectors with no labels :

There are no ground-truth answers and no loss to minimise.

What Can We Discover?

  1. Natural Groupings: Do the data points naturally fall into discrete families or clusters?
  2. Compact Representations: Can we capture 95% of the information with a few primary axes?
  3. Complex Manifolds: Can we unroll high-dimensional data onto a 2D map to visualise it?
  4. Density Estimation: What probability distribution generated the data?

What we'll cover

Two clean pairs of Linear Baseline Non-linear / Probabilistic Generalisation:

Grouping Without Labels

Imagine you are given customer purchase records, hospital patient metrics, or galaxy spectral emissions.

You are asked: "Are there distinct types or groups here?"

The goal:

Partition unlabelled observations into disjoint groups such that:

  • Points in the same group are as close (similar) as possible.
  • Points in different groups are as far apart as possible.

Mathematical Formulation of k-Means

Let each cluster be represented by a prototype centre:

For each data point , define a categorical cluster assignment:

where denotes that data point is assigned to cluster :

The objective function

We want to minimise the Within-Cluster Sum of Squares (WCSS), also called Inertia or Distortion:

In my words, this means bundle the points up. The more bundled the better.

The Great British Bundle Off.

Surprise! Finding the global minimum of is NP-hard.

Lloyd's Coordinate Descent Algorithm

Stuart Lloyd (1957) proposed an elegant alternating algorithm:

Step 1: Assignment Step

Hold centroids fixed. Assign each point to its nearest centroid:

Step 2: Update Step

Hold assignments fixed. Recompute each centroid as the mean of its assigned points:

Iterations Progress in 2D

Notice how the centroids shift towards the centres of mass at each step until the objective stops decreasing.

We got clusters Steve

Convergence & Stopping Criteria

Monotonic Convergence: Each step is guaranteed to decrease (or leave unchanged) the objective . Because there are only finitely many assignments, Lloyd's algorithm is guaranteed to converge in a finite number of iterations!

Practical Stopping Criteria:

  1. Centroids stop moving: for all .
  2. Assignments stabilise: for all (no points switch clusters).
  3. Objective threshold: Relative decrease threshold .
  4. Iteration budget: A fixed maximum number of iterations is reached (e.g. max_iter = 300).

The Initialisation Trap & k-Means++

Because k-Means is iterative, it converges to a local minimum. Bad initial seeds can lead to terrible clusterings.

Arthur & Vassilvitskii (2007) introduced k-Means++:

  1. Pick the first centroid uniformly at random from the data.
  2. For each point , compute its distance to the closest chosen centroid:

  3. Pick the next centroid with probability proportional to squared distance:

  4. Repeat until centroids are chosen, then run standard Lloyd's.

k-Means++ Sampling Geometry

Points far from existing centroids have high probability of becoming the next seed, spreading initial centres across clusters:

Selecting : The Elbow Method

  • As increases, always decreases (at , ).
  • Look for the "Elbow": the point where adding another cluster gives diminishing returns.
  • Supplement with the Silhouette Score measuring cluster cohesion vs separation.

Silhouette Analysis

For each sample in cluster :

  • Mean intra-cluster distance (cohesion):

  • Mean nearest-cluster distance (separation):

The Silhouette coefficient of sample is:

Silhouette Analysis

  • : Point is tightly inside its own cluster and well separated from neighbours.
  • : Point sits right on the boundary between two clusters.
  • : Point is closer to a neighbouring cluster than its own (probably misclustered).

Where Does k-Means Struggle?

k-Means makes three assumptions:

  1. Spherical Clusters: Assumes isotropic variance in all directions.
  2. Equal Cluster Sizes & Densities: Struggles when clusters have different spreads.
  3. Hard Boundaries: Forces every point to belong 100% to a single cluster, with zero uncertainty.

Soft Clustering

A patient whose blood markers are on the borderline between "Healthy" and "Pre-diabetic".

  • In k-Means, they are forced into one cluster ( or with 100% certainty).
  • In reality, we want a probability expressing confidence:

where is the probability of assigning data point to cluster .

Gaussian Mixture Model (GMM)

For every data point , imagine nature generating the sample in two steps:

  1. Pick a cluster :

  2. Sample the observation :

Marginal Density (The Resulting Data Distribution)

By the law of total probability:

where the multivariate Gaussian density is:

This is a bit ugly, especially the multivariate Gaussian but don't panic.

GMM Visually

is a latent (unobserved) variable; is the observed data:

The Marginal Likelihood

For i.i.d. observations , the log-likelihood is:

Notice the sum inside the logarithm! We cannot distribute the log. Setting derivatives to zero does not yield a closed-form solution. All cluster parameters are coupled together!

Plan B quit and build an AI startup without GMMs, wait for GPT-X to solve it. 😉

Expectation-Maximisation (EM)

Dempster, Laird, and Rubin (1977) proposed the EM algorithm to solve this:

1. E-Step (Expectation / Soft Responsibilities)

Bayes' theorem to compute the posterior probability (responsibility ) that cluster produced point :

Here is a number, you plug in and get a number.

2. M-Step (Maximisation / Update Parameters)

Re-estimate parameters using the soft responsibilities as weights:

GMM EM: A Concrete 1D Walkthrough

Imagine we record for 4 students preparing for exams:

Two natural groups exist (Casual Revision vs Intensive Cramming).
We fit Gaussians with initial priors and initial means .

E-Step: Evaluating Responsibilities

Given Gaussian likelihood :

Student (Hours) Density Density Responsibility Responsibility
0.2420 0.0000 1.000 0.000
0.3989 0.0001 0.999 0.001
0.0000 0.2420 0.000 1.000
0.0000 0.0540 0.000 1.000

M-Step Parameter Updates

Effective Cluster Sizes

Updated Priors & Means

Updated Variances (General Case)

In a single EM step, the means snap to the cluster averages and variance shrinks to capture the tight spread!

GMM EM Progress in 2D

The elliptical contours adapt their orientation and scale to fit the data:

k-Means is GMM

What happens to GMM if we fix equal priors and isotropic spherical covariance , and let the variance shrink to zero?

As :

  • The component closest to dominates exponentially.
  • The soft responsibility collapses to a 0 or 1 hard step function:

There is a connection here with LLM temperatures, can you spot it?

Hard vs Soft Clustering

Covariance Matrix Structures

We can control the geometric complexity of the Gaussian ellipses through :

  • Spherical: Isotropic circles (1 variance parameter per cluster).
  • Diagonal: Axis-aligned ellipses (no feature correlations, parameters per cluster).
  • Tied: All clusters share the exact same orientation and shape ( parameters).
  • Full: Freely oriented, unconstrained ellipses ( parameters).

Model Selection: Why Likelihood Fails

In k-Means, we used inertia and the elbow heuristic to choose . Can we choose and the covariance structure simply by maximising training log-likelihood?

  • Adding more Gaussian components () or choosing more flexible covariance matrices ('full' vs 'spherical') always increases (or maintains) the training log-likelihood.

  • In the extreme limit , every single point can have its own zero-variance Gaussian spike, sending likelihood to .

Bayesian Information Criterion (BIC)

Introduced by Gideon Schwarz (1978), BIC balances data likelihood against parameter count:

where:

  • is the maximised likelihood of the fitted model.
  • is the total number of free parameters estimated by EM.
  • is the number of observed data points.

Regularisation:

  • Goodness-of-Fit Error (): Decreases as the model fits the data closer.
  • Complexity Penalty (): Increases with more clusters and parameters.

We select the model that minimises BIC!

Counting Parameters

Why does the choice of covariance type drastically alter BIC?

For a mixture of Gaussians in dimensions:

  • Mixing weights : parameters (since ).
  • Cluster centres : parameters.
  • Covariances :
    • Spherical: parameters.
    • Diagonal: parameters.
    • Tied: shared parameters.
    • Full: parameters.

Model Selection with BIC in Action

  • Log-Likelihood (A): Monotonically rises as — blind to overfitting.
  • BIC Score (B): Penalises extra parameters, producing a clear minimum at (full).

Remember likelihood

To evaluate how well a fitted GMM explains data :

  1. Per-Cluster Density: For point , evaluate each Gaussian component:

  2. Mixture Density (Per Sample): Weight each Gaussian by its mixing proportion :

  3. Total Log-Likelihood (All Samples): Sum the logs assuming i.i.d. observations:

Toy Walkthrough: BIC Calculation

Consider fitting GMMs with full covariance to samples in dimensions ().
From our parameter formulas: .

Model Fit Error () Penalty () BIC Decision
Underfitting
Optimal (Minimum)
Overfitting (Penalised)

Fitting a GMM in scikit-learn

import numpy as np
from sklearn.mixture import GaussianMixture

# 1. Instantiate GMM with chosen K and covariance structure
gmm = GaussianMixture(
    n_components=3,
    covariance_type='full',   # 'spherical', 'diag', 'tied', or 'full'
    n_init=10,                # Run 10 EM restarts to avoid local minima
    random_state=42
)

# 2. Fit parameters via EM
gmm.fit(X)

# 3. Model inspection
print("Mixing weights:", gmm.weights_)          # pi_k
print("Cluster means:\n", gmm.means_)           # mu_k
print("Log-likelihood:", gmm.score(X) * len(X)) # ln p(X)
print("Model BIC:", gmm.bic(X))                 # -2 ln L + p ln N

# 4. Inference on new points
labels = gmm.predict(X_new)           # Hard assignment: argmax_k gamma_nk
probs  = gmm.predict_proba(X_new)     # Soft responsibilities: gamma_nk in [0, 1]

GMM Model Selection in Python

We can automate finding the optimal pair in a few lines of code:

lowest_bic = np.inf
best_gmm = None

# Grid search across candidate cluster counts and covariance structures
for cov_type in ['spherical', 'diag', 'tied', 'full']:
    for k in range(1, 8):
        candidate = GaussianMixture(n_components=k, covariance_type=cov_type, n_init=5, random_state=42)
        candidate.fit(X)
        bic = candidate.bic(X)
        if bic < lowest_bic:
            lowest_bic, best_gmm = bic, candidate

print(f"Optimal Model: K = {best_gmm.n_components}, Covariance = '{best_gmm.covariance_type}'")
print(f"Minimum BIC:   {lowest_bic:.2f}")

Interactive Clustering Demo

Source: https://scikit-learn.org/stable/modules/clustering.html

Dimensionality: The Redundancy Problem

Now we turn from grouping data points to finding compact representations.

Suppose we record 100 features for every Imperial student: revision hours, mock exam score, problem sheet grades, library visits, tutorial attendance...

  • Correlation: E.g. Revision Hours and Mock Exam Score vary together.
  • Useless: Most features are utterly useless, e.g. Study Room

Core Question: Can we find a compact representation of data?

Principal Component Analysis (PCA)

Karl Pearson (1901)

"On lines and planes of closest fit to systems of points in space."

The Intuitive Idea:

  • If data points are spread out widely along a direction, that direction carries real information.
  • If points barely vary along a direction, that direction carries mostly redundant noise.
  • Find the direction of maximum variance.

Step 1: Zero Mean (Centring the Data)

Given samples :

  1. Compute the sample mean vector:
  2. Centre every point by subtracting the mean:

  3. Stack the centred observations as rows into the design matrix :

Centring Data & Origin Shift

  • Centred Data (B): Coordinate origin translated directly to the sample mean .
  • Every observation is now a displacement vector emanating from the centre of mass.

Deriving the Covariance Matrix

Multiplying by computes all pairwise column dot products simultaneously:

  • is symmetric () and positive semi-definite ().
  • Diagonal entries ; off-diagonal entries .

This tells us how individual features increase or decrease together.

The Maximum Variance Principle

Let be a unit direction vector ().

The scalar projection of every point onto is:

The sample variance of these projected coordinates along is:

Maximising Variance

We want to find the direction that maximises the projected variance:

We constrain because we only care about the direction of the axis, so we restrict to lie on the unit sphere ().

In other words, Steve, give me the direction that maximises the variance.

Visualising Projected Variance

  • Max Variance Direction (): Projections are widely spread out, preserving maximal information about the data.

Eigenvalue Magic

  • Objective function:
  • Constraint function:

At the optimum, any move along the sphere must not increase the objective, so the gradient can only point straight out from the sphere, parallel to .

This is constrained optimisation, we don't really care so don't worry about the math too much. It is called Lagrange multipliers and is out of scope.

Geometric Intuition: Gradient Alignment

Eigenvalue Decomposition

This is the standard eigenvalue equation! The projected variance is simply:

must be an eigenvector of , and the eigenvalue is the projected variance! To maximise variance, we simply choose the eigenvector with the largest eigenvalue .

Subsequent Principal Components

What about the second direction ?

We maximise variance subject to:

  1. Unit length:
  2. Orthogonal to :
  • is the eigenvector with the second largest eigenvalue .
  • In general, the -th principal component is the -th eigenvector of , capturing variance .

PCA by Hand: A 2D Toy Walkthrough

Consider 4 students with 2 features ():

Student Study Hours () Mock Exam Score ()
1 2 1
2 4 5
3 6 5
4 8 9

Step 1: Means

Toy Walkthrough: Covariance Matrix

Step 2: Centred Data

Step 3: Covariance Matrix

Toy Walkthrough: Eigenvalues

Step 4: Characteristic Equation

Using the quadratic formula:

Notice: .

The first principal component alone captures:

Principal Axes & Orthogonal Projections

Singular Value Decomposition (SVD)

In practice, PCA is computed via the (SVD) of centred data :

where and are orthogonal (), and contains singular values .

  1. Numerical Stability:
    Forming is problematic. SVD operates directly on .
  2. Efficiency When :
    If we have samples and features (e.g. genomics), is a massive matrix (). Truncated SVD computes the top components without ever forming the matrix.

Why Does SVD Solve PCA?

Substitute the SVD of into the definition of the covariance matrix :

This becomes the exact spectral eigendecomposition .

  • The right vectors are the eigenvectors of .
  • The covariance eigenvalues are directly: .

Selecting

  • Elbow Rule: Retain components up to the "elbow" where eigenvalue drop-off flattens.
  • Threshold Rule: Retain components to explain 90% or 95% of cumulative variance.

Interactive PCA Demo

Rotate the projection axis to see the projected variance bar peak at the first eigenvector:

MNIST PCA

import matplotlib.pyplot as plt
from sklearn.datasets import fetch_openml
from sklearn.decomposition import PCA

# 1. Fetch MNIST (70,000 samples, 784 features)
X, y = fetch_openml('mnist_784', version=1, return_X_y=True, as_frame=False)

# 2. Fit PCA and project down to 2 components
pca = PCA(n_components=2)
X_pca = pca.fit_transform(X / 255.0)  # scale pixel values to [0, 1]

# 3. Plot the projection
plt.figure(figsize=(8, 6))
scatter = plt.scatter(
    X_pca[:, 0], X_pca[:, 1], c=y.astype(int), cmap='tab10', s=1, alpha=0.5
)
plt.colorbar(scatter, label='Digit Class')
plt.xlabel(f'PC1 ({pca.explained_variance_ratio_[0]:.1%} variance)')
plt.ylabel(f'PC2 ({pca.explained_variance_ratio_[1]:.1%} variance)')
plt.title('MNIST Projected onto First 2 Principal Components')
plt.tight_layout()
plt.show()

MNIST PCA

The Dual Perspective: Reconstruction Error

Pearson's original 1901 formulation pursued an identical goal from a different angle: finding the subspace of minimum reconstruction error!

Given centred data , find an orthonormal linear bottleneck basis matrix () minimising the average squared reconstruction error:

Pythagoras strikes back

For any data point , its squared distance from the origin is fixed and constant, that's where the data point lives.

The Zero-Sum Game

Averaging the Pythagorean identity over all observations:

Total variance (remember mean is 0) is a fixed property of the dataset. Therefore, every bit of variance captured along the subspace directly subtracts from the reconstruction error:

Maximising the shadow along the axis forces the perpendicular drop to be as short as possible.

Cool...

Where Does PCA Fail?

PCA is fast, global, and computed in closed form via SVD. But it has one fundamental limitation:

Linear
PCA can only project onto flat linear hyperplanes.

If data lives on a curved, non-linear manifold (e.g. an S-curve, a Swiss roll, or natural images), a linear projection fails...

So what?

In PCA, our bottleneck projection is strictly linear:

What if we replace this linear projection with deep neural networks?

It's happening again, we are going from some linear transformation to a non-linear one in the name of AI!

From PCA to Autoencoders

Instead of a linear matrix , decompose dimensionality reduction into:

The Encoder Network :
Compresses high-dimensional input into a low-dimensional bottleneck code ():

The Decoder Network :
Reconstructs the original input from the latent code:

The Objective Function:

The Autoencoder Architecture

An Autoencoder learns an identity mapping through an informational bottleneck:

  • Trained end-to-end via stochastic gradient descent (SGD) and backpropagation.

Linear AE is PCA

What happens if the encoder and decoder are strictly linear (no activation functions)?

where and .

  • linear algebra (spectral decomposition) and connectionist learning (gradient descent) converge to the exact same geometric solution...

Ouch... Where is my AI? AGI, SI?

Linear AE is PCA

I don't know what I was expecting when I asked AI to generate a plot, they are equal so we have 2 overlapping lines with some fancy callout boxes 🫠

Non-linear Autoencoders: Manifold Learning

By introducing non-linear activations (), the network can parametrise curved, non-linear manifolds:

  • PCA Fails: Folded branches collide and overlap in linear projection .
  • Autoencoder Succeeds: A 1D latent code cleanly unrolls the non-linear manifold, discovering the true intrinsic coordinate!

Bottleneck Capacity

What controls what an autoencoder learns? The bottleneck dimension relative to :

Undercomplete Autoencoder ():

  • Latent dimension is smaller than input dimension.
  • The network is forced to compress: it must discard noise and capture the dominant manifold structure.

Overcomplete Autoencoder ():

  • What if ? Without restrictions, the network can cheat by learning the trivial identity function:

  • You can enforce regularisation like sparsity Sparse Auto Encoders.

Representation Learning

I want to take this moment to highlight few connecting points:

  1. We said Neural Networks can learn intermediate representations, e.g. to solve XOR.
  2. These representations are hopefully some compact version of the input data.
  3. They are learned through some goal and optimisation: reconstruction in the case of autoencoders, next-token prediction in the case of language models.

In the grand scheme, it is about compressing and processing information in a compact and meaningful way.

Denoising Autoencoders

Introduced by Pascal Vincent et al. (2008):

  1. Take a clean data point .
  2. Inject stochastic noise to create a corrupted version :

  3. Feed corrupted into the encoder, but force the decoder to reconstruct the original clean :

Diffusion models generate images by chaining this idea.

Interactive Autoencoder Demo

Method Target Output Optimisation Key Strength Main Limitation
k-Means Hard clusters Lloyd's coordinate descent Extremely fast, scales to millions of points Assumes spherical, equal-variance clusters
GMM Soft probabilities Expectation-Maximisation (EM) Models uncertainty Sensitive to initialisation
PCA Continuous linear axes Exact SVD / Eigendecomposition Closed-form, global optimum, captures max variance Restricted to flat linear subspaces
Autoencoders Continuous non-linear coordinates Gradient descent (Backprop) on MSE Parametric, unrolls curved manifolds, reusable Non-convex loss, requires architecture tuning

The Practitioner's Workflow

                         UNLABELLED DATA X in R^{N x d}
                                       |
                 ┌─────────────────────┴─────────────────────┐
                 ▼                                           ▼
          DISCRETE CLUSTERING                         CONTINUOUS EMBEDDINGS
      (Which family does this belong to?)         (What are the coordinates?)
                 │                                           │
         Are clusters spherical?                     Is the data linear?
          /                 \                         /               \
        YES                  NO                     YES                NO
         |                   |                       |                  |
    k-Means++           GMM with EM                 PCA          Deep Autoencoder
   (Elbow / Sil.)     (BIC / Covariance)        (Scree Plot)      (MLP / CNN / DAE)

Where is my artificial brain?

Don't underestimate representation learning as a core concept. I think we faired better than expected but still:

AI ML Numbers game

Next Part: We'll let machines take actions and evolve.