The Project | About Us | Contribute | Donations | License

HOME

The entire module up to this point has been supervised learning: every historical delivery came with its label (delayed, delivery_minutes) and the model learned to predict it. But there are Rutalia questions with no possible label: which zones of the city behave alike? What types of delivery actually exist? Nobody has marked the answer in the data because nobody knows it: the structure has to be discovered, not imitated. That is unsupervised learning, and its central technique is clustering: grouping examples so that those in a group resemble each other more than they resemble those in other groups. In this lesson we will implement k-means from scratch (Lloyd's algorithm), look at hierarchical clustering and DBSCAN, close a circle left open since module 3 (cutting Kruskal's MST is clustering) and apply everything to a real Rutalia problem: grouping zones by demand pattern to redesign the routes — which will bring us back, elegantly, to the TSP of module 2.

Learning without labels: what changes
k-means: Lloyd's algorithm
Initialization and k-means++
Choosing k: the elbow method and the silhouette
Limitations of k-means
Hierarchical clustering: the dendrogram
Single-linkage = MST: closing the circle from module 3
DBSCAN: density-based clustering
Rutalia application: zones by demand pattern (and back to the TSP)
Comparison table

Learning without labels: what changes

Without y, three fundamental things change:

There is no correct answer. In classification, the confusion matrix settled the argument. Here "the true clustering" does not exist: different groupings can be valid for different purposes.
Evaluation is internal. We measure geometric properties (compactness, separation) and, above all, usefulness to the business.
Distance is everything. "Resembling" means being close in feature space — so the scaling from 05-01 goes from important to critical: a badly scaled feature defines the clusters all by itself.

k-means: Lloyd's algorithm

k-means looks for k centroids (mean points) that minimize the inertia: the sum of squared distances from each point to its nearest centroid. Minimizing it exactly is NP-hard (an old acquaintance from 02-02), so an excellent iterative heuristic is used, Lloyd's algorithm:

flowchart TD
    A["1. Pick k initial centroids"] --> B["2. ASSIGN: each point,\nto its nearest centroid"]
    B --> C["3. UPDATE: each centroid,\nto the mean of its points"]
    C --> D{"Did any assignment change?"}
    D -->|yes| B
    D -->|no| E["Converged: final clusters"]

Each step carries a guaranteed-improvement logic: assigning each point to its nearest centroid cannot increase the inertia, and moving each centroid to the mean of its points cannot either (the mean is the point that minimizes the sum of squared distances — the same fact that appeared in least squares in 05-03). The inertia decreases or stays put at every iteration and the possible assignments are finite: the algorithm always terminates. What it does not guarantee is terminating at the global optimum — like the greedy algorithms of 02-02, it converges to whichever local minimum happens to be nearby.

def kmeans_lloyd(X, k, seed=42, max_iter=100):
    """k-means from scratch. Returns assignments, centroids and inertia."""
    rng = np.random.default_rng(seed)
    centroids = X[rng.choice(len(X), k, replace=False)]    # k random points
    for _ in range(max_iter):
        # ASSIGN: distance matrix (n × k) via broadcasting
        dists = np.linalg.norm(X[:, None, :] - centroids[None, :, :], axis=2)
        assign = dists.argmin(axis=1)              # nearest centroid
        # UPDATE: each centroid to the mean of its points
        new = np.array([X[assign == c].mean(axis=0) if (assign == c).any()
                        else centroids[c] for c in range(k)])
        if np.allclose(new, centroids):            # no changes: converged
            break
        centroids = new
    inertia = float(((X - centroids[assign]) ** 2).sum())
    return assign, centroids, inertia

Complexity per iteration: O(n·k·d) — linear in everything, which is why k-means scales to Rutalia's millions of records. The number of iterations to convergence is usually small in practice.

Initialization and k-means++

Lloyd's Achilles heel: the result depends on the initial centroids. Two unlucky starts (e.g., two centroids landing in the same natural group) converge to different local minima, some of them bad. Two complementary remedies:

Multiple restarts: run the algorithm several times with different seeds and keep the lowest inertia (sklearn: n_init=10 by default).
k-means++: pick the initial centroids already spread apart. The first one at random; each subsequent one is drawn with probability proportional to D(x)², the squared distance from the point to its nearest already-chosen centroid — points far from everything chosen so far get more tickets in the raffle.

def kmeans_pp_init(X, k, rng):
    """k-means++ initialization: well-separated centroids with high probability."""
    centroids = [X[rng.integers(len(X))]]          # the first one, uniform
    for _ in range(k - 1):
        d2 = np.min([((X - c) ** 2).sum(axis=1) for c in centroids], axis=0)
        probs = d2 / d2.sum()                      # ∝ squared distance
        centroids.append(X[rng.choice(len(X), p=probs)])
    return np.array(centroids)

k-means++ doesn't just work well empirically: it guarantees, in expectation, an inertia within an O(log k) factor of the optimum — a probabilistic approximation bound, a relative of the guarantees we discussed in 02-02. It is sklearn's default (init="k-means++").

Choosing k: the elbow method and the silhouette

In supervised learning, hyperparameters were chosen with cross-validation against the label. Without a label, we need internal criteria:

Elbow method: plot the inertia against k. Inertia always drops as k grows (with k=n it would be 0: each point its own cluster — clustering's version of overfitting), but it drops in large leaps while k is below the natural number of groups and in small steps afterwards. The "elbow" of the curve suggests the natural k. It is visual and somewhat subjective.
Silhouette coefficient: for each point, compare a (mean distance to the members of its own cluster) with b (mean distance to the members of the nearest foreign cluster): s = (b − a) / max(a, b) ∈ [−1, 1]. Close to 1: well grouped; close to 0: on the border; negative: probably in the wrong cluster. The dataset's mean silhouette is maximized over several candidate k values — objective and automatable (sklearn.metrics.silhouette_score), at the price of costing O(n²).

for k in range(2, 9):
    _, _, inertia = kmeans_lloyd(X_s, k)
    print(f"k={k}: inertia={inertia:10.1f}")   # look for the elbow in the sequence

Limitations of k-means

k-means has its geometry written into it: by assigning via distance to the centroid, it assumes clusters that are convex, roughly spherical and of similar size. It fails predictably when reality doesn't cooperate:

Situation	What k-means does	What would be needed
Elongated or nested clusters (two rings, two half moons)	Cuts them in half with straight borders	DBSCAN or single-linkage hierarchical
Very different densities	The dense cluster "steals" points from the sparse one	DBSCAN
Outliers	They drag the centroids (the mean is not robust)	DBSCAN (marks them as noise)
Unknown k	You must supply it	Hierarchical (decide afterwards) or DBSCAN (infers it)

These limitations motivate exactly the next two algorithms.

Hierarchical clustering: the dendrogram

Agglomerative clustering does not fix k in advance: it builds all the granularities at once.

Start with n clusters (each point its own).
Merge the two closest clusters.
Repeat until only one remains.

The merge history forms a tree, the dendrogram: the height of each merge is the distance at which it happened. Cutting the dendrogram at a height yields a partition; lower cuts give more clusters. You choose the granularity after seeing the structure — the luxury k-means doesn't offer.

What remains is defining "distance between clusters" (the linkage), and the choice changes the algorithm's character:

Linkage	Distance between clusters A and B	Tendency
Single	minimum over pairs (a∈A, b∈B)	Chains: follows elongated structures (and creates unwanted "bridges")
Complete	maximum over pairs	Compact, round clusters
Average	mean over pairs	Compromise
Ward	smallest inertia increase when merging	The most k-means-like; the usual default

In scipy: scipy.cluster.hierarchy.linkage(X, method="ward") and dendrogram(...) to draw it; fcluster cuts at the chosen height. Cost O(n² log n) and memory O(n²): perfect for grouping Rutalia's 9 zones or a few thousand customers, unworkable for millions of rows.

Single-linkage = MST: closing the circle from module 3

In the final exercise of 03-04 we mentioned, almost as a curiosity, that cutting the k−1 most expensive edges of Kruskal's minimum spanning tree splits the graph into k "natural" groups. Time to close that circle: that was single-linkage hierarchical clustering, exactly.

The correspondence is precise. Kruskal (03-04) processes the edges from lowest to highest weight, joining components with union-find (01-04). Single-linkage merges at each step the two clusters with the closest pair of points between them. They are the same algorithm: Kruskal's sequence of unions over the complete distance graph is single-linkage's sequence of merges, and the dendrogram heights are the weights of the MST edges in the order Kruskal accepts them. That is why cutting the k−1 most expensive edges of the MST is equivalent to cutting the dendrogram to obtain k clusters.

def clusters_via_mst(X, k):
    """Single-linkage via Kruskal: cuts the k-1 most expensive MST edges."""
    n = len(X)
    edges = sorted((np.linalg.norm(X[i] - X[j]), i, j)
                   for i in range(n) for j in range(i + 1, n))
    parent = list(range(n))                   # union-find from 01-04
    def find(a):
        while parent[a] != a:
            parent[a] = parent[parent[a]]     # path compression
            a = parent[a]
        return a
    accepted = []
    for w, i, j in edges:                     # Kruskal, exactly as in 03-04
        ri, rj = find(i), find(j)
        if ri != rj:
            parent[ri] = rj
            accepted.append((w, i, j))
    # Full MST: n-1 edges. Removing the k-1 most expensive = k components.
    parent = list(range(n))
    for w, i, j in accepted[:-(k - 1)]:       # we skip the last k-1 (most expensive)
        parent[find(i)] = find(j)
    roots = {find(i) for i in range(n)}
    return np.array([sorted(roots).index(find(i)) for i in range(n)])

A graph algorithm from module 3, with the data structure from module 1, turns out to be an unsupervised learning algorithm from module 5. The borders between "classical algorithmics" and "machine learning" are thinner than they seem: ML is, to a large extent, algorithmics applied to data.

DBSCAN: density-based clustering

DBSCAN abandons centroids and hierarchies: a cluster is a dense region of points, whatever its shape. Two parameters: eps (neighborhood radius) and min_samples (how many neighbors make a point "dense").

Core point: has ≥ min_samples neighbors within distance ≤ eps.
A cluster is a maximal set of cores connected by neighborhood chains (density reachability), plus the border points attached to them.
Whatever no core reaches is noise: it gets no cluster, labeled −1.

Algorithmically it is an old friend: a BFS-style exploration (03-02) over the "density graph" — from each unvisited core, expand the neighborhood; if a neighbor is a core, keep expanding through it. Cost O(n²) naive, O(n log n) with spatial indexes.

When it beats k-means:

Arbitrary shapes: it follows the density, imposing no spheres — it solves the rings and half moons where k-means fails.
No k required: the number of clusters emerges from the data.
Explicit noise: outliers contaminate no cluster; at Rutalia, anomalous deliveries (wrong address, serious incident) get flagged for free — clustering and anomaly detection in one.

Its weaknesses: choosing eps takes finesse (rule of thumb: plot the distance to the k-th neighbor and look for the elbow), and it struggles when clusters have very disparate densities. In sklearn: DBSCAN(eps=0.5, min_samples=5).fit_predict(X_s).

Rutalia application: zones by demand pattern (and back to the TSP)

Let's put the pieces together on a real problem. Rutalia wants to redesign its routes: which zones behave alike and could share a delivery strategy? Using the canonical dataset, we build each zone's demand profile: a vector of features aggregated per zone.

data = generate_dataset()                      # canonical, seed 42 (05-01)
profiles, names = [], []
for z in ZONES:
    m = data["zone"] == z
    profiles.append([
        data["delivery_minutes"][m].mean(),    # mean duration
        data["delayed"][m].mean(),             # delay rate
        (data["departure_hour"][m] >= 17).mean(),  # evening delivery fraction
        data["weight_kg"][m].mean(),           # mean weight
    ])
    names.append(z)
P = np.array(profiles)
P_s = (P - P.mean(axis=0)) / P.std(axis=0)     # standardize: critical

assign, cents, _ = kmeans_lloyd(P_s, k=3)
for c in range(3):
    print(f"cluster {c}: {[names[i] for i in range(9) if assign[i] == c]}")

With the generator's hidden truth, the groups that emerge are interpretable: the congested zones (CEN, MER, HOS) cluster together through their high durations and delay rates; the free-flowing ones (IND, PAR, UNI) form another profile; the intermediate ones (ALM, RIO, EST), the third. Nobody labeled the zones: the structure was in the data and the algorithm made it visible. (With 9 points, hierarchical clustering with a dendrogram would be equally appropriate and more illustrative; k-means wins when grouping thousands of customers or millions of deliveries.)

And here the course's grand circle closes: with the clusters in hand, operations assigns couriers per group of homogeneous zones... and planning each courier's route within their group is, once again, the TSP of module 2 — that canonical instance with optimum 35.22 km — over the module 3 graph, with the A* of module 4 for the legs and this module's models predicting the times that feed the weights. Machine learning does not replace classical algorithmics: it feeds it better data and better groupings, and that final redesign decision is made by the people in operations with the analyses in hand.

Comparison table

Criterion	k-means	Hierarchical (agglomerative)	DBSCAN
Requires k?	Yes, up front	No (cut the dendrogram afterwards)	No (emerges from density)
Cluster shapes	Spherical, convex	Depends on linkage (single: elongated)	Arbitrary
Noise handling	No (outliers drag centroids)	No	Yes, label −1
Cost	O(n·k·d) per iteration	O(n² log n), memory O(n²)	O(n²), O(n log n) with index
Scales to millions	Yes	No	Acceptable with spatial indexes
Deterministic	No (initialization)	Yes	Yes (given eps/min_samples)
Delicate parameters	k, initialization	linkage, cut height	eps, min_samples
Use it when…	Compact groups, massive data	You want the hierarchy / little data	Odd shapes, outliers, unknown k

Common Mistakes and Tips

Clustering without standardizing. In the zone profile, delivery_minutes (~20-60) would crush delayed (~0-0.4): the clusters would come from a single feature. With no label to give the mistake away, it can go unnoticed — always standardize.
A single k-means run. Lloyd converges to whatever local minimum comes up. Use k-means++ and several restarts (n_init), and compare inertias.
Taking the elbow's k as revealed truth. The elbow suggests; the silhouette quantifies; the business decides. If k=3 and k=4 have similar silhouettes but operations works with 4 delivery teams, the answer is 4.
Interpreting clusters as real categories. A cluster is a geometric regularity, not a truth about the world. Always validate with domain knowledge before acting on it.
Using DBSCAN with very disparate densities. A single eps cannot serve a dense city center and a sparse industrial estate at the same time; consider hierarchical clustering or variants (HDBSCAN).
Tip: plot whenever you can (project to 2D if needed). In clustering, one look at a colored scatter is worth more than any internal metric.

Exercises

Lloyd under the magnifying glass. Modify kmeans_lloyd to store the inertia at every iteration and run it on the zone profiles (k=3) with 5 different seeds. Check that (a) the inertia never rises within a run, and (b) different seeds may end at different final inertias. Which of the two properties did we prove and which is only empirical?
The elbow and the silhouette. On a standardized sample of 500 deliveries from the canonical dataset (features: distance_km, departure_hour, weight_kg), compute the k-means inertia for k=2..8 and the mean silhouette (sklearn.metrics.silhouette_score) for each k. Do the elbow and the silhouette maximum agree? What would you do if they don't?
Kruskal as a clusterer. Apply clusters_via_mst with k=3 to the standardized zone profiles and compare the partitions with k-means and with sklearn's AgglomerativeClustering(n_clusters=3, linkage="single"). Verify that MST and single-linkage match exactly and explain why k-means may differ.

Solutions

Exercise 1:

# Inside the kmeans_lloyd loop, after assigning:
#   history.append(float(((X - centroids[assign]) ** 2).sum()))
for s in range(5):
    assign, _, inertia = kmeans_lloyd(P_s, 3, seed=s)
    print(f"seed {s}: final inertia = {inertia:.3f}")

(a) is a theorem: each assignment step and each update step cannot increase the inertia (we argued this when presenting Lloyd), so the history is monotonically non-increasing in every run. (b) is empirical: convergence is to a local minimum, and which one you hit depends on the start — you will see seeds ending at different inertias. Hence k-means++ and the restarts.

Exercise 2:

from sklearn.cluster import KMeans
from sklearn.metrics import silhouette_score

Xm = np.column_stack([data["distance_km"], data["departure_hour"],
                      data["weight_kg"]])[:500]
Xm = (Xm - Xm.mean(axis=0)) / Xm.std(axis=0)
for k in range(2, 9):
    km = KMeans(n_clusters=k, n_init=10, random_state=42).fit(Xm)
    print(f"k={k}: inertia={km.inertia_:8.1f}  "
          f"silhouette={silhouette_score(Xm, km.labels_):.3f}")

With individual delivery data (fairly continuous, with no strongly marked groups), it is normal for the elbow to be soft and the silhouette modest with its maximum at a small k: an honest signal that the group structure is weak. If the elbow and the silhouette disagree, neither one rules: you examine the candidates (visualization, mean profiles of each cluster) and interpretability and business use decide. "There are no crisp clusters" is also a valid outcome of the analysis.

Exercise 3:

from sklearn.cluster import AgglomerativeClustering

mst = clusters_via_mst(P_s, 3)
single = AgglomerativeClustering(n_clusters=3, linkage="single").fit_predict(P_s)
km = kmeans_lloyd(P_s, 3)[0]
print(mst, single, km, sep="\n")

MST and single-linkage produce the same partition (the labels may be permuted — compare the grouping, not the numbers): they are the same algorithm, as we argued, with the merge sequence dictated by the MST edges. k-means may differ because it optimizes something else (inertia with respect to centroids, favoring round, balanced groups), while single-linkage follows chains of proximity even when they form elongated groups or groups of disparate sizes.

Conclusion

With clustering we complete the module's map: unsupervised learning, where no label is imitated but structure is discovered. You saw k-means from the inside (Lloyd and its monotone inertia, k-means++ with its probabilistic guarantee, the elbow and the silhouette for choosing k, and its spherical assumptions), hierarchical clustering with its multiscale dendrogram, and DBSCAN finding free-form clusters and flagging the noise — with the revelation that single-linkage is exactly module 3's Kruskal cut at its most expensive edges: classical algorithmics and machine learning are the same continent. And with it the entire module 5 closes: we started from a paradigm shift — letting the rules be learned from the data — and traveled it in full: the honest evaluation protocol (05-01), classifying (05-02), predicting quantities with gradient descent as the engine (05-03), composing nonlinearity with networks and backpropagation (05-04) and discovering structure without labels (05-05), always on the same 2,000 deliveries of the canonical Rutalia dataset. You now hold the course's complete arsenal: analysis and data structures (module 1), optimization (module 2), graphs (module 3), searching and sorting (module 4) and machine learning (module 5). In module 6 we stop sharpening tools and step out into the world: real case studies where these algorithms combine — optimization in industry, graphs in social networks, searching and sorting at scale, and this module's ML put to work in production — starting in 06-01 with optimization in industry.

Clustering Algorithms

Contents

Learning without labels: what changes

k-means: Lloyd's algorithm

Initialization and k-means++

Choosing k: the elbow method and the silhouette

Limitations of k-means

Hierarchical clustering: the dendrogram

Single-linkage = MST: closing the circle from module 3

DBSCAN: density-based clustering

Rutalia application: zones by demand pattern (and back to the TSP)

Comparison table

Common Mistakes and Tips

Exercises

Solutions

Conclusion

Advanced Algorithms

Module 1: Introduction to Advanced Algorithms

Module 2: Optimization Algorithms

Module 3: Graph Algorithms

Module 4: Search and Sorting Algorithms

Module 5: Machine Learning Algorithms

Module 6: Case Studies and Applications

Module 7: Final Project