Chapter 9 — Unsupervised Learning

Up to now, almost every dataset we touched came with labels — the "right answer" written next to each example. Unsupervised learning is what we do when those answers are missing. We hand the algorithm a pile of unlabeled examples and ask it to find structure on its own: groups, a simpler description, the shape of the distribution, or a guess about what a "typical" example looks like.

In this chapter you will learn:

• How to estimate the shape of a data distribution (density estimation), with and without assuming a bell curve
• How to group similar examples with k-means, Gaussian mixture models, hierarchical clustering, and DBSCAN
• How to choose the number of clusters with the elbow method
• How to squeeze many features into a few with Principal Component Analysis (PCA)
• How recommender systems guess missing ratings with collaborative filtering

9.1 Density Estimation

Density estimation asks: "I have a sample of data points — what curve did they probably come from?" The answer is an estimate of the probability density function (pdf), the function f whose value at any point x tells you how likely data is to appear near x.

Two families of approaches exist:

• Parametric: you assume a specific shape — most commonly the bell-shaped Gaussian (normal) distribution — and just fit its parameters: the mean μ and the variance σ² (or the covariance matrix Σ in higher dimensions). Fast and simple, but if the real data isn't bell-shaped, the fit is poor.
• Non-parametric: you make no strong shape assumption. Kernel density estimation (KDE) is the classic example: you drop a little smooth "bump" (a kernel, usually a tiny Gaussian) on top of every data point and add them all up. The result is a smooth curve that follows whatever shape the data actually has.

A bandwidth b controls how wide each bump is. Too small and every point becomes its own spike (overfitting); too large and everything melts into one fat blob (underfitting). It's the same bias–variance trade-off we keep meeting.

Everyday analogy: parametric estimation is like insisting "the crowd must be one normal blob, let me find its center." KDE is like pouring a handful of sand on every person in the crowd and looking at the resulting dune — it reveals whatever shape the crowd really has, even if it's two separate groups.

%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import make_blobs
from sklearn.neighbors import KernelDensity

plt.rcParams["figure.figsize"] = (7, 4)
plt.rcParams["axes.grid"] = True

# Build a 1-D dataset that is clearly NOT a single bell curve:
# two separate blobs stitched together (a bimodal distribution).
np.random.seed(42)
data = np.concatenate([
    np.random.normal(loc=-3.0, scale=0.8, size=120),
    np.random.normal(loc=3.0, scale=0.8, size=120),
])

# --- Parametric fit: assume ONE Gaussian, estimate mean & variance ---
mu_hat = data.mean()
sigma2_hat = data.var()

# --- Non-parametric fit: Kernel Density Estimation ---
kde = KernelDensity(bandwidth=0.5, kernel="gaussian").fit(data[:, None])

# Grid of x values where we will evaluate both estimates
x_grid = np.linspace(-7, 7, 400)

# Gaussian pdf:  1/sqrt(2*pi*sigma^2) * exp(-(x-mu)^2 / (2*sigma^2))
gaussian_pdf = (1.0 / np.sqrt(2 * np.pi * sigma2_hat)) * \
               np.exp(-(x_grid - mu_hat) ** 2 / (2 * sigma2_hat))

# KDE pdf: score_samples returns log-density, so we exponentiate
kde_pdf = np.exp(kde.score_samples(x_grid[:, None]))

# Plot the histogram of the data plus both density estimates
plt.hist(data, bins=30, density=True, color="lightgray", edgecolor="white",
         label="data histogram")
plt.plot(x_grid, gaussian_pdf, color="tab:red", lw=2,
         label=f"Gaussian fit (mu={mu_hat:.1f})")
plt.plot(x_grid, kde_pdf, color="tab:blue", lw=2, label="KDE (bandwidth=0.5)")
plt.title("Density Estimation: one Gaussian vs KDE")
plt.xlabel("value")
plt.ylabel("density")
plt.legend()
plt.show()

Reading the plot

• The gray histogram is the raw data — notice the two peaks around −3 and +3.
• The red Gaussian tries to describe all of it with a single bell centered at 0. It spreads wide to "cover" both bumps, but it dips where the data actually peaks and rises where the data is actually empty. That is the cost of a wrong assumption.
• The blue KDE places a small bump on each point and adds them up, so it naturally reproduces the two peaks. No shape assumption needed — just a bandwidth choice.

When you have no reason to believe the data is bell-shaped, KDE is the safer bet. When you do expect a single bell (many measurement errors are roughly Gaussian), the parametric fit is cheaper and good enough.

9.2 Clustering

Clustering is the unsupervised version of classification. We want to assign each example a cluster label — but nobody gave us the labels, and we don't even know how many groups there are. The algorithm has to discover the grouping from the geometry of the data alone.

This makes clustering hard to evaluate: with no true labels, "is this clustering good?" becomes a judgment call. Different algorithms make different geometric assumptions, so the "best" one depends on the unknown shape of your data. Let's meet the four most useful families.

9.2.1 K-Means

K-means is the workhorse of clustering. You choose k (the number of clusters) and the algorithm does this:

1. Initialize: place k points called centroids somewhere in the feature space (often randomly).
2. Assign: give every example the label of its nearest centroid (Euclidean distance).
3. Update: move each centroid to the mean (average position) of all the examples now assigned to it.
4. Repeat steps 2–3 until the assignments stop changing.

The quantity it minimizes is the sum of squared distances from each point to its centroid — scikit-learn calls this inertia. Lower inertia means tighter clusters.

$$\text{inertia} = \sum_{i=1}^{N} \min_{j}\ \lVert x_i - c_j \rVert^2$$

Two catches: you must pick k yourself, and because the start is random, two runs can give different answers (scikit-learn runs several starts and keeps the best by default).

Everyday analogy: k-means is like dropping k postmen in a city, having each house join its nearest postman, then moving each postman to the middle of their new round, and repeating until the postmen settle.

from sklearn.cluster import KMeans

# Three natural groups in 2-D (we pretend we don't know the labels)
X, y_true = make_blobs(n_samples=300, centers=3, cluster_std=0.9, random_state=7)

# Fit k-means with k=3. n_init=10 tries 10 random starts and keeps the best.
km = KMeans(n_clusters=3, n_init=10, random_state=7)
labels = km.fit_predict(X)

plt.scatter(X[:, 0], X[:, 1], c=labels, cmap="tab10", s=20, edgecolor="white")
# Mark the learned centroids with big red X's
plt.scatter(km.cluster_centers_[:, 0], km.cluster_centers_[:, 1],
            marker="X", s=200, color="red", edgecolor="black", label="centroids")
plt.title("K-Means clustering (k=3)")
plt.xlabel("feature 1")
plt.ylabel("feature 2")
plt.legend()
plt.show()

C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=2.
  warnings.warn(

Reading the plot

The three colored clouds are the clusters k-means found, and the big red X marks are the final centroids — each sits at the average of its cluster's points. Because k-means measures distance with straight lines, it produces round (hyperspherical) clusters. If your real groups are stretched, nested, or crescent-shaped, k-means will force them into circles and get it wrong — we'll see this clearly when we meet DBSCAN.

# Try k from 1 to 8 and record the inertia (within-cluster sum of squares)
ks = range(1, 9)
inertias = []
for k in ks:
    km = KMeans(n_clusters=k, n_init=10, random_state=7)
    km.fit(X)
    inertias.append(km.inertia_)

plt.plot(list(ks), inertias, marker="o", color="tab:purple")
plt.axvline(3, color="tab:red", linestyle="--", label="k=3 (the 'elbow')")
plt.title("Elbow method: inertia vs number of clusters")
plt.xlabel("number of clusters k")
plt.ylabel("inertia (sum of squared distances to centroid)")
plt.legend()
plt.show()

C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=2.
  warnings.warn(
C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=2.
  warnings.warn(
C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=2.
  warnings.warn(
C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=2.
  warnings.warn(
C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=2.
  warnings.warn(
C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=2.
  warnings.warn(
C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=2.
  warnings.warn(
C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=2.
  warnings.warn(

Reading the elbow

Inertia always decreases as k grows — add more centroids and points get closer to one. So we don't look for the minimum; we look for the kink (the "elbow") where the drop suddenly flattens. Here inertia falls fast from k=1 to k=3, then barely improves after. That kink at k=3 matches the true number of blobs. The elbow method is subjective (the kink isn't always obvious), but it's a quick, practical first guess. The book also mentions more rigorous options: prediction strength and the gap statistic.

9.2.2 Gaussian Mixture Models (EM)

K-means gives a hard assignment: each point belongs to exactly one cluster. A Gaussian Mixture Model (GMM) is softer. It models the data as a weighted sum of k Gaussian blobs:

$$f(x) = \sum_{j=1}^{k} \pi_j\, \mathcal{N}(x \mid \mu_j, \Sigma_j)$$

Each cluster j has a mean μⱼ, a covariance Σⱼ (which lets the blob be an ellipse of any size, stretch, and rotation), and a weight πⱼ (how big that cluster is overall). The parameters are learned with the Expectation-Maximization (EM) algorithm, which alternates two steps:

• E-step: given the current Gaussians, compute for each point the probability it came from each cluster (a soft assignment).
• M-step: given those probabilities, update each Gaussian's mean, covariance, and weight (using probability-weighted averages).

Repeat until the parameters stop moving. It's k-means' cousin: k-means is the special case where every point is 100% in one cluster and the blobs are identical round spheres. GMM's ellipses and soft probabilities let it handle overlapping clusters that k-means would chop arbitrarily.

from sklearn.mixture import GaussianMixture
from matplotlib.patches import Ellipse

# Two overlapping blobs so soft assignment actually matters
Xg, _ = make_blobs(n_samples=220, centers=[(-2, 0), (2, 0)],
                   cluster_std=1.6, random_state=3)

gmm = GaussianMixture(n_components=2, covariance_type="full",
                      random_state=3, n_init=1, max_iter=100)
gmm.fit(Xg)
labels_g = gmm.predict(Xg)

plt.scatter(Xg[:, 0], Xg[:, 1], c=labels_g, cmap="coolwarm", s=20, edgecolor="white")

# Draw an ellipse for each component from its mean & covariance
ax = plt.gca()
for j in range(gmm.n_components):
    mean = gmm.means_[j]
    cov = gmm.covariances_[j]
    # Eigen-decomposition gives the ellipse's axes & rotation
    eigvals, eigvecs = np.linalg.eigh(cov)
    order = eigvals.argsort()[::-1]
    eigvals, eigvecs = eigvals[order], eigvecs[:, order]
    angle = np.degrees(np.arctan2(eigvecs[1, 0], eigvecs[0, 0]))
    # 2-sigma ellipse: full width/height = 2 * (2 * sigma) = 4 * sqrt(eigenvalue)
    width, height = 4 * np.sqrt(eigvals)
    ell = Ellipse(xy=mean, width=width, height=height, angle=angle,
                  edgecolor="black", facecolor="none", lw=2, linestyle="--")
    ax.add_patch(ell)

plt.title("Gaussian Mixture Model: soft clusters with 2-sigma ellipses")
plt.xlabel("feature 1")
plt.ylabel("feature 2")
plt.show()

C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=1.
  warnings.warn(

Reading the plot

The two colors are the most-likely cluster for each point, and the dashed ellipses show each Gaussian's shape (a 2σ contour). Notice the ellipses overlap in the middle — points there genuinely could belong to either cluster. With gmm.predict_proba(X) you'd get the actual soft probabilities (e.g. 0.7 / 0.3) instead of a forced single label. Because the covariance matrix can stretch and rotate, GMM captures elongated, tilted clusters that round k-means would misplace.

9.2.3 Hierarchical Clustering

Hierarchical clustering builds a tree (a dendrogram) of clusters. The common agglomerative (bottom-up) version starts with every point as its own cluster, then repeatedly merges the two closest clusters until only one giant cluster remains. You don't have to pick k up front: you build the whole tree once, then "cut" it at the height that gives you the number of clusters you want.

A dendrogram visualizes these merges: each horizontal line is a merge, and its height shows how dissimilar the two clusters were when they joined. Tall merges mean "these groups are far apart" — a natural place to cut.

Everyday analogy: like a family tree of species, built by joining the most similar organisms first. Cut the tree high to get a few big kingdoms; cut it low to get many small families.

from scipy.cluster.hierarchy import linkage, dendrogram

# Tiny dataset so the dendrogram stays readable
Xh, _ = make_blobs(n_samples=20, centers=3, cluster_std=0.8, random_state=11)

# linkage matrix: 'ward' merges the pair that least increases total variance
Z = linkage(Xh, method="ward")

dendrogram(Z, color_threshold=4.0)
plt.title("Hierarchical clustering dendrogram (20 points)")
plt.xlabel("example index")
plt.ylabel("merge distance")
plt.show()

Reading the dendrogram

The three colored branches near the bottom are the three tight clusters forming first. They then join each other higher up (around distance 4–6) — that big jump in height is the visual cue that three is a natural number of clusters. Cutting the tree with a horizontal line at height ≈ 4 would yield exactly three clusters. The dendrogram is a nice bonus: it shows the whole hierarchy, not just one flat answer.

9.2.4 DBSCAN

DBSCAN is density-based. Instead of asking for the number of clusters, you give it two numbers:

• ε (eps): the radius to look around each point.
• min_samples: how many points must sit within that radius for a spot to be "crowded."

It grows each cluster from a crowded point outward, absorbing neighbors, then their neighbors, and so on, until it hits a region that's too sparse. Points in sparse regions are labeled outliers (cluster −1). The payoff: DBSCAN finds clusters of arbitrary shape — rings, crescents, blobs — because it follows density, not distance to a center.

The catch: picking ε is fiddly, and one fixed ε struggles when clusters have very different densities (the book recommends HDBSCAN to fix that — it keeps DBSCAN's strengths but only needs min_samples).

Let's put DBSCAN head-to-head with k-means on two interleaved moons, a shape k-means fundamentally cannot handle.

from sklearn.cluster import DBSCAN
from sklearn.datasets import make_moons

# Two interleaved half-moons -- a classic non-convex shape
Xm, _ = make_moons(n_samples=250, noise=0.07, random_state=9)

# K-means insists on 2 round clusters -- it will slice each moon in half
km_m = KMeans(n_clusters=2, n_init=10, random_state=9)
labels_km = km_m.fit_predict(Xm)

# DBSCAN follows the density of each moon
db = DBSCAN(eps=0.2, min_samples=5)
labels_db = db.fit_predict(Xm)

fig, axes = plt.subplots(1, 2, figsize=(11, 4))

axes[0].scatter(Xm[:, 0], Xm[:, 1], c=labels_km, cmap="coolwarm", s=20, edgecolor="white")
axes[0].set_title("K-Means on two moons (FAILS)")
axes[0].set_xlabel("feature 1")
axes[0].set_ylabel("feature 2")

# Mark DBSCAN outliers (-1) in black
outlier = labels_db == -1
axes[1].scatter(Xm[~outlier, 0], Xm[~outlier, 1], c=labels_db[~outlier],
                cmap="coolwarm", s=20, edgecolor="white")
axes[1].scatter(Xm[outlier, 0], Xm[outlier, 1], color="black", s=30,
                marker="x", label="outlier")
axes[1].set_title("DBSCAN on two moons (follows the curves)")
axes[1].set_xlabel("feature 1")
axes[1].set_ylabel("feature 2")
axes[1].legend()

plt.tight_layout()
plt.show()

C:\Users\DELL\anaconda3\Lib\site-packages\sklearn\cluster\_kmeans.py:1419: UserWarning: KMeans is known to have a memory leak on Windows with MKL, when there are less chunks than available threads. You can avoid it by setting the environment variable OMP_NUM_THREADS=1.
  warnings.warn(

Reading the comparison

• Left (k-means): forced to produce two round regions, it draws a straight cut right through the middle of each moon. Each "cluster" contains half of the top moon and half of the bottom one — geometrically wrong even though k=2 was "correct."
• Right (DBSCAN): it simply walks along each dense crescent, labeling the whole top moon as one cluster and the whole bottom moon as another. A few sparse points (black ×) are flagged as outliers rather than forced into a group.

The lesson: match the algorithm's assumption to the data's shape. Round blobs → k-means. Overlapping tilted blobs → GMM. Weird shapes / unknown count → DBSCAN/HDBSCAN.

Algorithm	Assumes	Picks k?	Cluster shape	Outliers?
K-Means	round, similar-size blobs	yes (you set it)	spherical	no
GMM (EM)	Gaussian ellipses	yes (you set it)	elliptical	no
Hierarchical	a merge tree	cut later	flexible	no
DBSCAN	dense regions	no	arbitrary	yes

9.3 Dimensionality Reduction

Sometimes you have too many features — hundreds or thousands — and most are redundant, noisy, or correlated with each other. Dimensionality reduction rewrites each example using fewer new features while keeping as much information as possible. Two classic payoffs:

• Visualization: humans can read at most 3D plots, so we squash high-dimensional data to 2D/3D to "see" it.
• Simpler, cleaner models: fewer features mean faster training, less overfitting, and often better interpretability.

The book names three workhorses: Principal Component Analysis (PCA), UMAP, and autoencoders. We'll focus on PCA (the oldest and most transparent) and sketch the others.

9.3.1 Principal Component Analysis (PCA)

PCA finds a new set of axes — the principal components — that are just rotations of the original features. The first component points in the direction of the greatest variance (spread) in the data; the second is perpendicular to it and points in the direction of the next greatest variance; and so on. Each component comes with a number: how much variance it explains.

To reduce to D_new dimensions, you keep the top D_new components and project the data onto them. Because the early components capture most of the spread, dropping the later ones usually loses little. You can also reconstruct an approximation of the original point by going back up — the lost dimensions are exactly the low-variance ones, which mostly held noise.

Everyday analogy: photographing a 3D object from the angle that shows the most detail. The first photo captures the longest silhouette; a second photo from the side captures what's left. Throw away the blurry head-on shot and you've kept the essence in two pictures.

from sklearn.decomposition import PCA

# A synthetic 6-dimensional dataset with strong structure:
# 3 tight blobs live roughly in a 2-D plane, so 2 components should
# capture most of the variance.
Xp, yp = make_blobs(n_samples=200, centers=3, n_features=6,
                    cluster_std=0.5, random_state=5)

pca = PCA(n_components=2)
Xp_2d = pca.fit_transform(Xp)

fig, axes = plt.subplots(1, 2, figsize=(12, 4))

# 2-D projection colored by the true (hidden) blob id
axes[0].scatter(Xp_2d[:, 0], Xp_2d[:, 1], c=yp, cmap="coolwarm",
                s=25, edgecolor="white")
axes[0].set_title("PCA: 6D data projected to 2D")
axes[0].set_xlabel("1st principal component")
axes[0].set_ylabel("2nd principal component")

# Explained variance ratio for ALL 6 components
pca_full = PCA(n_components=6).fit(Xp)
axes[1].bar(range(1, 7), pca_full.explained_variance_ratio_,
            color="tab:purple", label="per component")
axes[1].plot(range(1, 7), np.cumsum(pca_full.explained_variance_ratio_),
             marker="o", color="tab:red", label="cumulative")
axes[1].set_title("Explained variance ratio (6 features)")
axes[1].set_xlabel("principal component")
axes[1].set_ylabel("fraction of variance explained")
axes[1].legend()

plt.tight_layout()
plt.show()

Reading the PCA plots

• Left: even though each point really has 6 numbers, PCA found two new axes that separate the three blobs almost perfectly. We can now plot 6-dimensional data on a flat page.
• Right: the first two purple bars are tall and the red cumulative line jumps to near 1.0 by component 2 — meaning two components capture almost all the variance. The remaining four components explain almost nothing (mostly noise), so dropping them costs us little.

# Reconstruct the 6-D points from only 2 components and compare to originals.
pca2 = PCA(n_components=2).fit(Xp)
Xp_compressed = pca2.transform(Xp)               # 2-D code
Xp_reconstructed = pca2.inverse_transform(Xp_compressed)  # back to 6-D (approx)

print("Original      point 0:", np.round(Xp[0], 2))
print("Reconstructed point 0:", np.round(Xp_reconstructed[0], 2))
print()
print("Original      point 1:", np.round(Xp[1], 2))
print("Reconstructed point 1:", np.round(Xp_reconstructed[1], 2))
print()

# Mean squared reconstruction error across all points and features
err = np.mean((Xp - Xp_reconstructed) ** 2)
print(f"Mean squared reconstruction error: {err:.3f}")

Original      point 0: [-0.52 -7.49  7.42 -4.5  -2.24 -4.28]
Reconstructed point 0: [-0.68 -7.05  7.68 -4.93 -1.92 -4.04]

Original      point 1: [-0.35 -6.64  7.64 -4.6  -2.01 -4.27]
Reconstructed point 1: [-0.72 -6.88  7.51 -4.78 -1.91 -3.95]

Mean squared reconstruction error: 0.159

Reading the reconstruction

The reconstructed 6-number vectors are close to the originals but not exact — PCA kept the two "important" directions and threw away the four low-variance ones. The mean squared reconstruction error is small precisely because the discarded components held little variance. This is the core PCA trade: give up a little accuracy for a big drop in size.

9.3.2 UMAP and Autoencoders (conceptual)

PCA is fast and linear, but it can only capture straight-axis structure. Two more powerful alternatives:

• UMAP (and its cousin t-SNE) are non-linear methods built for visualization. They define a similarity between nearby points in the high-dimensional space, then place points in a 2D/3D space so those similarities are preserved as well as possible (the book writes this as a cross-entropy between the high- and low-dimensional similarity graphs). The result: clusters of weird shape often separate cleanly on a 2D plot — better than PCA for eyeballing. Slower than PCA, faster than an autoencoder.

• Autoencoders are tiny neural nets trained to copy their input to their output through a narrow bottleneck layer. Because the bottleneck has fewer neurons than the input, the network is forced to find a compact code. That bottleneck output is your reduced representation; the decoder half reconstructs the original. We met autoencoders in Chapter 7 — they also double as outlier detectors (an outlier reconstructs badly).

Reach for PCA when you want speed and interpretability, UMAP/t-SNE when you want an honest-to-the-eye 2D picture, and an autoencoder when the data is highly non-linear and you need a learned, reusable encoder.

9.4 Collaborative Filtering

How does a streaming service know you'll like a show you've never seen? Collaborative filtering is the classic recommender idea: use the ratings of many users to predict the rating of one user. The data lives in a user–item matrix R, where R[u, i] is the rating user u gave to item i — and most entries are missing (nobody has rated everything).

Two simple flavors:

• User-based: find users whose taste is most similar to yours (by comparing rating vectors with cosine similarity), then predict your missing rating as a weighted average of what those similar users gave the item.
• Item-based: find items most similar to the target item (again by cosine similarity over who-rated-what), then predict your rating as a weighted average of your own ratings on those similar items.

A more powerful approach is matrix factorization: approximate the big sparse matrix R as the product of two small matrices — a user matrix and an item matrix. Each user and item gets a short "taste vector," and a predicted rating is just the dot product of the two vectors. This is what won the Netflix Prize and still powers many real recommenders.

Everyday analogy: asking friends who share your taste to recommend a movie. The more a friend's past likes overlap with yours, the more you trust their suggestion — that weighting is exactly cosine similarity at work.

import pandas as pd
from sklearn.metrics.pairwise import cosine_similarity

# Tiny ratings matrix: 5 users x 5 movies, 0 means "not rated"
ratings = np.array([
    [5, 4, 0, 1, 0],   # user 0
    [4, 5, 0, 1, 0],   # user 1
    [1, 0, 5, 4, 4],   # user 2
    [0, 1, 4, 5, 5],   # user 3
    [5, 4, 1, 0, 0],   # user 4
])
users  = ["U0", "U1", "U2", "U3", "U4"]
movies = ["M0", "M1", "M2", "M3", "M4"]

df = pd.DataFrame(ratings, index=users, columns=movies)
print("Ratings matrix (0 = unrated):")
print(df)
print()

# Heatmap of the raw ratings
plt.imshow(ratings, cmap="YlGnBu", aspect="auto")
plt.colorbar(label="rating")
plt.xticks(range(len(movies)), movies)
plt.yticks(range(len(users)), users)
plt.title("User-Movie ratings matrix (0 = unrated)")
plt.xlabel("movie")
plt.ylabel("user")
plt.show()

# --- Predict U0's missing rating for M4 using USER-based CF ---
# For this small demo we treat unrated (0) as literal 0 in the cosine
# computation -- a simplification that real systems avoid.
sim = cosine_similarity(ratings.astype(float))

print("Cosine similarity of U0 to each user:")
for u, s in zip(users, sim[0]):
    print(f"  {u}: {s:.2f}")

# Predict U0 -> M4: similarity-weighted average of OTHER users' M4 ratings
target_user, target_movie = 0, 4
others = [u for u in range(len(users)) if u != target_user]
weights = sim[target_user, others]
other_ratings = ratings[others, target_movie]
rated = other_ratings > 0   # keep only users who actually rated M4
pred = np.average(other_ratings[rated], weights=weights[rated])
print(f"\nPredicted rating for U0 -> M4: {pred:.2f}")
print(f"(based on {int(rated.sum())} users who rated M4, similarity-weighted)")

Ratings matrix (0 = unrated):
    M0  M1  M2  M3  M4
U0   5   4   0   1   0
U1   4   5   0   1   0
U2   1   0   5   4   4
U3   0   1   4   5   5
U4   5   4   1   0   0

Cosine similarity of U0 to each user:
  U0: 1.00
  U1: 0.98
  U2: 0.18
  U3: 0.17
  U4: 0.98

Predicted rating for U0 -> M4: 4.48
(based on 2 users who rated M4, similarity-weighted)

Reading the demo

• The heatmap shows two clear taste groups: U0/U1/U4 love M0/M1 and rate M3 low; U2/U3 love M2/M3/M4 and rate M0 low. The zeros are the gaps a recommender must fill.
• U0 never rated M4. Their strongest matches are U1 and U4 (cosine similarity ≈ 0.98 — almost identical taste), while U2 and U3 are only weakly similar (≈ 0.18).
• To predict U0 → M4 we take a similarity-weighted average of the ratings other users gave M4. But here's the catch: U0's true taste-mates (U1, U4) didn't rate M4 either, so the prediction must lean on the weakly-similar U2/U3, who both gave M4 high marks (4 and 5). The result is a high predicted rating (~4.5) — yet we should distrust it, because it rests on users who barely share U0's taste.

This is the sparsity problem, and it's why real recommenders prefer matrix factorization: instead of comparing raw rating rows, it learns a short latent "taste vector" for every user and item, so a prediction can borrow strength across all users and items at once — even when a specific taste-mate hasn't rated the target item.

Key Takeaways

• Unsupervised learning finds structure in unlabeled data — there's no "right answer" to score against, so judging quality is harder than in supervised learning.
• Density estimation models the data's pdf: a single Gaussian is cheap but assumes a bell shape; KDE follows any shape by summing little bumps, with a bandwidth controlling smoothness.
• K-means iterates assign → update centroids and minimizes inertia; it makes round clusters and needs you to pick k (the elbow method is a practical first guess).
• Gaussian Mixture Models (fit by EM) give soft probabilities and ellipse-shaped clusters — ideal for overlapping data.
• Hierarchical clustering builds a merge dendrogram you cut at the height that gives the cluster count you want.
• DBSCAN is density-based: no k needed, handles arbitrary shapes, and flags outliers — but its ε parameter is fiddly (HDBSCAN helps).
• PCA rotates data onto axes of maximum variance, letting you project high-dimensional data to 2D/3D and reconstruct it with small error; UMAP and autoencoders extend the idea to non-linear structure.
• Collaborative filtering predicts missing user–item ratings via cosine similarity (user-/item-based) or matrix factorization, powering recommender systems.

What's Next

Next we move beyond the standard supervised/unsupervised split in Chapter 10 — Other Forms of Learning, covering semi-supervised learning, active learning, transfer learning, one-shot learning, and more.

Exercises

These exercises cover unsupervised learning from Chapter 9: density estimation, the four clustering families (k-means, GMM, hierarchical, DBSCAN), dimensionality reduction (PCA, UMAP, autoencoders), and collaborative filtering. Try each before reading the hint.

1. (Conceptual) Contrast parametric (single Gaussian) and non-parametric (KDE) density estimation. When is each the better choice, and what does the KDE bandwidth control? Hint: a Gaussian is cheap but assumes a bell shape; KDE follows any shape by summing little bumps; bandwidth too small ⇒ spikes (overfit), too large ⇒ one fat blob (underfit).
2. (Conceptual) Describe the k-means loop and what inertia measures. Why does k-means produce round (spherical) clusters? Hint: assign each point to its nearest centroid, then move each centroid to its cluster's mean; inertia is the sum of squared distances to centroids; Euclidean distance ⇒ spherical clusters.
3. (Conceptual) Why can't you just pick the k with the lowest inertia, and how does the elbow method work? Hint: inertia always decreases as k grows (add a centroid and points get closer); instead look for the "kink" where the drop suddenly flattens.
4. (Conceptual) How does a GMM differ from k-means, and what do the EM E-step and M-step do? Hint: GMM is a weighted sum of Gaussians giving soft probabilities with elliptical clusters; the E-step computes each point's responsibility per cluster, the M-step updates means/covariances/weights.
5. (Conceptual) What advantage does a hierarchical clustering dendrogram offer over k-means? Hint: you build the merge tree once and can cut it at any height to get a different cluster count, and it reveals the whole hierarchy of merges.
6. (Conceptual) How does DBSCAN differ from k-means, and what are its two parameters? Hint: DBSCAN needs no k, follows dense regions to find clusters of arbitrary shape, and flags sparse points as outliers; its parameters are eps (radius) and min_samples (crowd threshold).
7. (Conceptual) Explain PCA: what are principal components, how do you choose how many to keep, and where does reconstruction error come from? Hint: components are rotated axes of maximum variance; keep the top components that capture most variance; reconstruction error comes from discarding the low-variance (mostly-noise) directions.
8. (Conceptual) Contrast user-based and item-based collaborative filtering, and explain why matrix factorization beats raw cosine similarity on sparse data. Hint: user-based compares similar users, item-based compares similar items (both by cosine similarity); matrix factorization learns latent taste vectors that borrow strength across all users/items even when a specific taste-mate hasn't rated the target.

Hands-On Coding Problems

9. (Coding) Density estimation: make a bimodal sample (two Gaussian blobs at −3 and +3), fit a single Gaussian (mean/std) and a KDE, and plot the histogram, the Gaussian curve, and the KDE together. Hint: scipy.stats.gaussian_kde (or sklearn.neighbors.KernelDensity) for the smooth curve; compute the single Gaussian as norm.pdf with the sample mean and std.
10. (Coding) K-means on make_blobs(centers=3, random_state=0): fit KMeans(n_clusters=3, n_init=10, random_state=0), plot the points colored by cluster with the centroids marked as red "X", and print the inertia. Hint: kmeans.cluster_centers_ and kmeans.inertia_.
11. (Coding) Elbow method: loop k from 1 to 8, fit KMeans(n_clusters=k, n_init=10, random_state=0) each time, record the inertia, and plot inertia vs k; mark the elbow. Hint: inertia always drops — spot where it flattens.
12. (Coding) GMM soft assignment: fit GaussianMixture(n_components=2, random_state=0) on two overlapping blobs and print predict_proba for the first few points (each row should sum to 1). Hint: gmm.predict_proba(X); optionally plot the clusters.
13. (Coding) DBSCAN vs k-means on make_moons(noise=0.05, random_state=0): fit DBSCAN(eps=0.2, min_samples=5) and KMeans(n_clusters=2), and plot both results side by side; note that DBSCAN follows the crescents while k-means cuts straight through. Hint: DBSCAN labels outliers as −1.
14. (Coding) PCA on a 6-feature make_blobs dataset: fit PCA(n_components=2), plot the 2D projection colored by blob, print explained_variance_ratio_, and reconstruct one point with inverse_transform. Hint: pca.fit_transform(X) then pca.inverse_transform(X_2d).
15. (Coding) User-based collaborative filtering: given a small 5×5 ratings matrix (0 = unrated), compute cosine similarities between users and predict the missing R[0,4] as a similarity-weighted average of other users' ratings for item 4. Hint: sklearn.metrics.pairwise.cosine_similarity; only average over users who actually rated item 4.

# Exercise 9: density estimation -- single Gaussian vs KDE
import numpy as np
import matplotlib.pyplot as plt

rng = np.random.default_rng(0)
data = np.concatenate([rng.normal(-3, 1.0, 200), rng.normal(3, 1.0, 200)])

# Single Gaussian fit (mean, std)
mu, sigma = data.mean(), data.std()

# TODO: build a KDE over `data` (scipy.stats.gaussian_kde or sklearn.neighbors.KernelDensity)
# TODO: plot the histogram, the Gaussian pdf (1/(sigma*sqrt(2pi))*exp(-(x-mu)^2/(2 sigma^2)), and the KDE
xs = np.linspace(-8, 8, 300)
plt.hist(data, bins=40, density=True, color="lightgray", alpha=0.6)
# gaussian_pdf = ...  # TODO
# kde_curve = ...     # TODO
# plt.plot(xs, gaussian_pdf, "r-")
# plt.plot(xs, kde_curve, "b-")
plt.show()

# Exercise 10: k-means with centroids and inertia
import matplotlib.pyplot as plt
from sklearn.cluster import KMeans
from sklearn.datasets import make_blobs

X, _ = make_blobs(n_samples=300, centers=3, cluster_std=1.0, random_state=0)

# TODO: fit KMeans(n_clusters=3, n_init=10, random_state=0)
km = None  # placeholder

# TODO: scatter X colored by km.labels_, mark km.cluster_centers_ with red "X"
# TODO: print km.inertia_

# Exercise 11: elbow method for choosing k
import matplotlib.pyplot as plt
from sklearn.cluster import KMeans
from sklearn.datasets import make_blobs

X, _ = make_blobs(n_samples=300, centers=3, cluster_std=1.0, random_state=0)

ks = range(1, 9)
inertias = []
for k in ks:
    # TODO: fit KMeans(n_clusters=k, n_init=10, random_state=0) and append .inertia_
    inertias.append(0.0)  # placeholder

# TODO: plot inertias vs ks and mark the elbow
plt.show()

# Exercise 12: GMM soft probabilities
import numpy as np
from sklearn.mixture import GaussianMixture
from sklearn.datasets import make_blobs

X, _ = make_blobs(n_samples=300, centers=[[-3, 0], [3, 0]],
                  cluster_std=1.8, random_state=0)  # overlapping blobs

# TODO: fit GaussianMixture(n_components=2, random_state=0)
gmm = None  # placeholder

# TODO: print gmm.predict_proba(X[:5]) -- each row should sum to ~1
print("soft probabilities:", None)

# Exercise 13: DBSCAN vs k-means on the moons
import matplotlib.pyplot as plt
from sklearn.cluster import KMeans, DBSCAN
from sklearn.datasets import make_moons

X, _ = make_moons(n_samples=300, noise=0.05, random_state=0)

# TODO: fit KMeans(n_clusters=2, n_init=10, random_state=0) and DBSCAN(eps=0.2, min_samples=5)
km = None   # placeholder
db = None   # placeholder

# TODO: plot the two results side by side (km.labels_ vs db.labels_);
#       note that DBSCAN labels outliers as -1
plt.show()

# Exercise 14: PCA projection, variance, and reconstruction
import matplotlib.pyplot as plt
from sklearn.decomposition import PCA
from sklearn.datasets import make_blobs

X, y = make_blobs(n_samples=300, n_features=6, centers=3, random_state=0)

# TODO: fit PCA(n_components=2) and transform X to 2D
pca = None       # placeholder
X_2d = X[:, :2]  # placeholder: pca.fit_transform(X)

# TODO: scatter X_2d colored by y; print pca.explained_variance_ratio_
# TODO: reconstruct the first point with pca.inverse_transform(X_2d[:1]) and print it

# Exercise 15: user-based collaborative filtering
import numpy as np
from sklearn.metrics.pairwise import cosine_similarity

# Rows = users U0..U4, columns = items M0..M4; 0 means "unrated"
R = np.array([
    [5, 4, 0, 1, 0],   # U0 has NOT rated M4
    [5, 5, 0, 0, 0],
    [0, 0, 5, 4, 5],
    [0, 1, 4, 5, 4],
    [5, 4, 0, 1, 0],
])
target_user, target_item = 0, 4

# TODO: compute cosine_similarity between users (use R as the vectors)
sim = None  # placeholder

# TODO: predict R[target_user, target_item] as a similarity-weighted average of
#       OTHER users' ratings for target_item (only those who rated it, i.e. R[u, target_item] != 0)
prediction = 0.0  # placeholder

print("user similarities:", None)
print("predicted rating U0 -> M4:", prediction)