Unsupervised Learning: Clustering, PCA and Anomaly Detection

Unsupervised methods find structure without labels. In interviews you will be asked to implement k-means, explain PCA through the covariance matrix and SVD, compare clustering algorithms, and say how you would evaluate something that has no ground truth. This chapter gives you working code and the reasoning to defend each choice.

1. Clustering: what it is and when it is useful

Clustering groups similar items. Typical uses: customer segmentation, grouping documents, finding patterns in sensor data, compressing data (vector quantisation), creating features (cluster IDs) and exploring data. There is no single "correct" clustering: the result depends on the distance measure, the features, and the algorithm's assumptions.

2. k-means

Objective: choose $k$ centroids to minimise the within-cluster sum of squared distances (inertia):

${\sum }_{j=1}^k{\sum }_{x\in C_j}\Vert x-{\mu }_j{\Vert }^2$

Algorithm (Lloyd's):

1. Initialise $k$ centroids.
2. Assign each point to its nearest centroid.
3. Update each centroid to the mean of its assigned points.
4. Repeat until assignments stop changing.

Each step cannot increase the objective, so it converges, but only to a local minimum. Therefore run several initialisations and keep the best, and use k-means++ initialisation (pick each new centroid with probability proportional to squared distance from the nearest chosen one), which spreads centroids and gives provable guarantees.

import numpy as np

def lloyd(X, k, iters=100, seed=0):
    rng = np.random.default_rng(seed)
    centers = X[rng.choice(len(X), k, replace=False)]
    for _ in range(iters):
        d = ((X[:, None, :] - centers[None, :, :]) ** 2).sum(-1)       # squared distances, shape (n, k)
        labels = d.argmin(1)
        new = np.array([X[labels == j].mean(0) if (labels == j).any() else centers[j] for j in range(k)])
        if np.allclose(new, centers):
            break
        centers = new
    inertia = ((X - centers[labels]) ** 2).sum()
    return labels, centers, inertia

def kmeans(X, k, n_init=10):
    runs = [lloyd(X, k, seed=s) for s in range(n_init)]       # restarts guard against a poor local minimum
    return min(runs, key=lambda r: r[2])

rng = np.random.default_rng(1)
blobs = np.vstack([rng.normal(c, 0.4, (100, 2)) for c in ([0, 0], [5, 5], [0, 5])])
labels, centers, inertia = kmeans(blobs, 3)

true_centers = np.array([[0, 0], [5, 5], [0, 5]], dtype=float)
for tc in true_centers:
    assert np.min(np.linalg.norm(centers - tc, axis=1)) < 0.3          # every true cluster is found
assert len(set(labels[:100])) == 1 and len(set(labels[100:200])) == 1  # clean separation

Complexity: $O(nkd)$ per iteration. It scales well, and mini-batch k-means scales further.

Assumptions and weaknesses: clusters that are roughly spherical, similar in size and density; sensitive to scale (standardise features) and outliers; must choose $k$ in advance; fails on elongated, ring-shaped or very unequal clusters.

Choosing $k$

• Elbow method: plot inertia against $k$ and look for the bend. Inertia always decreases with $k$, and the elbow is often ambiguous.
• Silhouette score: for each point, $s=\frac{b-a}{\max (a,b)}$ where $a$ is the mean distance to its own cluster and $b$ the mean distance to the nearest other cluster. Ranges from -1 to 1; higher is better.
• Gap statistic and stability under resampling.
• Business meaning: the number of segments a team can actually act on.

import numpy as np

def silhouette(X, labels):
    n = len(X); s = np.zeros(n)
    D = np.linalg.norm(X[:, None] - X[None], axis=-1)
    for i in range(n):
        same = (labels == labels[i]); same[i] = False
        a = D[i, same].mean() if same.any() else 0
        b = min(D[i, labels == c].mean() for c in set(labels) if c != labels[i])
        s[i] = (b - a) / max(a, b)
    return s.mean()

rng = np.random.default_rng(1)
X = np.vstack([rng.normal(c, 0.4, (60, 2)) for c in ([0, 0], [5, 5], [0, 5])])
lab3 = kmeans(X, 3)[0]; lab2 = kmeans(X, 2)[0]
assert silhouette(X, lab3) > silhouette(X, lab2) and silhouette(X, lab3) > 0.7

3. Other clustering methods

k-means is a special case of a GMM with equal spherical covariances and hard assignments. The EM algorithm alternates between estimating soft responsibilities (E-step) and re-estimating parameters (M-step), and it increases the likelihood every iteration.

4. Principal component analysis (PCA)

PCA finds orthogonal directions of maximum variance and projects data onto the top few. It is used for compression, visualisation, denoising and decorrelating features.

Derivation sketch: centre the data; the first principal component is the unit vector $v$ maximising $\mathrm{Var}(Xv)=v^{\top }\Sigma v$, where $\Sigma =\frac{1}{n-1}{X}^{\top }X$. The maximiser is the eigenvector of $\Sigma$ with the largest eigenvalue; the eigenvalue equals the variance captured. The next components are the next eigenvectors, orthogonal to the earlier ones.

Practical computation: use the SVD $X=US{V}^{\top }$. The rows of $V^{\top }$ are the principal directions, and the variance along component $i$ is $s_i^2/(n-1)$.

import numpy as np

rng = np.random.default_rng(0)
# 2D data stretched along the direction (1, 1)
z = rng.normal(0, 3, 500)
X = np.column_stack([z + rng.normal(0, 0.3, 500), z + rng.normal(0, 0.3, 500)])
Xc = X - X.mean(0)

U, S, Vt = np.linalg.svd(Xc, full_matrices=False)
var = S ** 2 / (len(X) - 1)
ratio = var / var.sum()

assert abs(abs(Vt[0] @ np.array([1, 1]) / np.sqrt(2)) - 1) < 0.02      # first component points along (1,1)
assert ratio[0] > 0.98                                                  # almost all variance is on it

# the eigen-decomposition of the covariance matrix gives the same variances
eigvals = np.sort(np.linalg.eigvalsh(np.cov(Xc.T)))[::-1]
assert np.allclose(eigvals, var)

# reconstruction from one component loses little
proj = Xc @ Vt[0]
recon = np.outer(proj, Vt[0])
assert np.mean((Xc - recon) ** 2) < 0.1

Using PCA well

• Always centre the data, and usually standardise features first; otherwise large-scale features dominate.
• Choose the number of components by explained variance (such as 90 to 95%), a scree plot, or downstream performance.
• PCA is linear and unsupervised: the highest-variance directions are not necessarily the most predictive of a target.
• It makes features less interpretable.
• PCA fitted on training data only; transform the others with the same components.

Beyond PCA

• t-SNE and UMAP: nonlinear embeddings for visualisation. They preserve local neighbourhoods, but distances between far-apart clusters and cluster sizes are not meaningful, and results depend on hyperparameters. Do not use them as features without care.
• Autoencoders: neural nonlinear compression; a linear autoencoder spans the same subspace as PCA.
• Matrix factorisation / SVD: the basis of classic recommenders and topic models.

5. Anomaly detection

Find rare points that do not fit the pattern: fraud, faults, intrusions.

Evaluate using whatever labels exist (precision at a review budget), human review of the top flagged cases, and the cost of alerts. In production, alert fatigue is the real constraint.

import numpy as np

rng = np.random.default_rng(0)
normal = rng.normal(0, 1, (500, 2))
outliers = np.array([[8.0, 8.0], [-7.0, 9.0]])
X = np.vstack([normal, outliers])

mu = normal.mean(0)
cov = np.cov(normal.T)
inv = np.linalg.inv(cov)
d = np.sqrt(np.einsum("ij,jk,ik->i", X - mu, inv, X - mu))         # Mahalanobis distance
flagged = np.where(d > 4.5)[0]
assert set(flagged) >= {500, 501}                                   # both planted anomalies are caught
assert len(flagged) <= 5                                            # almost no false alarms

6. Evaluating without labels

• Internal metrics: silhouette, Davies-Bouldin, inertia (compare within the same $k$ and method).
• Stability: do the clusters persist under resampling or different seeds?
• External validation if some labels exist: adjusted Rand index, normalised mutual information.
• Business validation: can you name and act on each cluster? Do clusters differ on metrics that matter?
• Downstream use: does adding the cluster feature improve a supervised model?

7. Common mistakes

• Not scaling features before k-means or PCA.
• Running k-means once and trusting a local minimum.
• Choosing $k$ only from the elbow without business or stability checks.
• Using k-means on categorical data (use k-modes or a different representation).
• Treating t-SNE or UMAP plots as quantitative evidence.
• PCA fitted on test data, or components selected using the target by accident.
• Reading PCA components as causal "factors".

8. Practice questions

1. Implement k-means. Why does it converge and why only to a local minimum?
2. How do you choose $k$? What are the limits of the elbow method?
3. Explain PCA via eigenvectors of the covariance matrix and via SVD.
4. When would you use DBSCAN rather than k-means?
5. Relate k-means to a Gaussian mixture model fitted with EM.
6. Why must you standardise features before PCA?
7. How would you detect anomalies in a stream of transaction data? How would you evaluate it?
8. What can and can't you conclude from a t-SNE plot?