DBSCAN, and when density beats centroids

Plot the two features and the picture is obvious to anyone in the room: a dense blob near the origin, a long crescent curving around it, and a hundred points scattered far out with nothing near them. Ask k-means for three clusters and it returns three wedges. It cuts the crescent in half, hands one half to the blob, and drags the scattered points into whichever centre is closest — which pulls that centre a long way out to reach them, so the wedge boundaries move too.

Every point got a label. Two of the three groups are artefacts of the method.

Dense regions, not nearest centres

DBSCAN — Ester, Kriegel, Sander and Xu’s 1996 algorithm — asks a different question. Instead of “which centre is this point nearest to”, it asks “is this point in a crowd, and does that crowd join up with another one”.

Two parameters define a crowd. eps is a radius. min_samples is how many points have to sit inside that radius, the point itself included. A point that clears the bar is a core point. Any point inside a core point’s radius is reachable from it. A cluster is a maximal set of points connected through chains of core points, and it swallows the non-core points sitting on its fringe as border points. Anything left over is labelled noise, and sklearn gives it the label -1.

Two clusters found by chaining neighbours, and thirteen points that join nothing. Eight points sit inside the marked radius, so at min_samples = 6 that point is a core point.

The run below uses the public dataset: one row per subscription, four columns counting what that customer did in their first thirty days.

import numpy as np
import pandas as pd
from sklearn.cluster import DBSCAN
from sklearn.preprocessing import StandardScaler

subs = pd.read_csv("https://dataacademy.ai/data/subscriptions.csv",
                   parse_dates=["signed_up_at"])
events = pd.read_csv("https://dataacademy.ai/data/events.csv",
                     parse_dates=["occurred_at"])

ev = events.merge(subs[["subscription_id", "customer_id", "signed_up_at"]],
                  on="customer_id")
month1 = ev[(ev.occurred_at - ev.signed_up_at).dt.days.between(0, 29)]
X = month1.pivot_table(index="subscription_id", columns="event_type",
                       values="event_id", aggfunc="count", fill_value=0)

Xs = StandardScaler().fit_transform(X)
labels = DBSCAN(eps=0.35, min_samples=10).fit_predict(Xs)

n_clusters = len(set(labels) - {-1})
print(n_clusters, (labels == -1).mean().round(3))

10 0.312

Ten clusters, and 31.2% of the 1,198 subscriptions attached to none of them. Hold that result; the last section works out what it means, and it is not what the count suggests.

Three consequences follow, and they are the whole reason to reach for this method.

Clusters can be any shape. Membership spreads through chains of neighbours, so a crescent, a ring or an S-curve comes out intact. k-means cannot produce a non-convex cluster at all; the straight boundaries between centres forbid it.

The number of clusters is an output. You do not supply k, which changes the conversation with stakeholders: nobody has to defend the choice of four segments, because nobody made it.

Points are allowed to belong to nothing. That single change removes the mechanism by which one distant outlier distorts a whole partition.

eps is the parameter that matters

min_samples is the easy one. It sets how many neighbours make a crowd, and so how aggressively the run treats sparse regions as noise. A common starting point is twice the number of dimensions, with a floor around 4 or 5, and higher when the data is noisy.

eps is the hard one, because it is a distance in the units of the scaled feature space and nobody has intuition for that. Too small and everything is noise. Too large and the whole dataset becomes one cluster. The useful range is often narrow.

The practical way to find it is the k-distance plot. For each point, measure the distance to its min_samples-th nearest neighbour, sort those distances, and look at where the curve turns up.

from sklearn.neighbors import NearestNeighbors

min_samples = 10
nn = NearestNeighbors(n_neighbors=min_samples).fit(Xs)
dist, _ = nn.kneighbors(Xs)          # column 0 is the point itself, at 0.0
kdist = np.sort(dist[:, -1])

print(np.quantile(kdist, [0.5, 0.9, 0.95, 0.99]).round(3))
print((kdist == 0).mean().round(3), len(X.drop_duplicates()), len(X))

[0.114 1.006 1.409 2.537]
0.33 546 1198

The curve is flat and low across the points inside dense regions, then bends sharply upward across the points that are not. The elbow is where “close” stops meaning anything, and eps belongs just before it. Read the quantiles next to it: an elbow at the 95th percentile means you have already decided that roughly 5% of the data comes back as noise.

The four quantiles above are a warning, not an elbow. The median point’s tenth neighbour sits 0.114 away and the 90th percentile point’s sits at 1.006 — an order of magnitude apart with nothing in between to aim at. A third of the rows have a tenth neighbour at distance exactly zero, because the four columns are small whole numbers and only 546 of the 1,198 rows are distinct. There is no value of eps here that separates dense from sparse; there are ties and then there is everything else.

When the features have physical meaning, skip the elbow and use the meaning. For GPS points clustered with the haversine metric, eps is a real distance and the right value is however far apart two pings can be and still count as the same stop.

One eps cannot fit two densities

The weakness is visible from the definition. There is one radius for the whole dataset, so DBSCAN assumes every cluster is dense in the same way.

Real data often is not. A dataset with one tight cluster and one diffuse one has no correct eps: small enough to keep the tight cluster separate and the diffuse one shatters into noise; large enough to hold the diffuse one together and the tight clusters merge. Chaining makes that merge easy — a thin bridge of points between two dense blobs connects them into one, because connectivity needs a chain, not a consensus.

HDBSCAN is the answer, and it is not a patch. It runs the density argument across all values of eps at once, builds a hierarchy using mutual reachability distance, then keeps the clusters that persist over the widest range of thresholds. Different clusters in the same run may live at different densities, and the parameter you supply becomes min_cluster_size — the smallest group you are willing to call a cluster, which is a question about the business rather than about the feature space.

from sklearn.cluster import HDBSCAN     # scikit-learn 1.3 and later

h = HDBSCAN(min_cluster_size=25).fit_predict(Xs)
print(len(set(h) - {-1}), (h == -1).mean().round(3))

13 0.433

If clusters at mixed densities are plausible in your data, start here rather than tuning eps for an afternoon. It is not a rescue, though: on the counts above it returns thirteen clusters and leaves 43.3% unassigned, which is more clusters and more noise than DBSCAN gave. Sweeping every eps at once cannot manufacture a density gap that the data does not contain.

Distance is the whole model

Both methods stand or fall on the distance function, and so does the parameter you spent all that effort choosing. eps is measured in whatever units the columns happen to have, so a revenue column in euros and an age column in years are not comparable, and the euro column decides every neighbourhood by itself. Standardise, or scale deliberately to encode which features matter — the same mechanism that makes scaling the real decision in k-means applies here with a sharper edge, because a bad scale does not merely distort the clusters, it makes a single global radius meaningless.

Pick the metric to match the data too. Cosine distance for text embeddings, haversine for coordinates, Gower or an explicit distance matrix for mixed types. DBSCAN accepts metric="precomputed", which is often the cleanest route.

Dimensionality is the last trap. As dimensions grow, distances concentrate, the k-distance curve loses its elbow, and no eps separates anything. Reduce first, and accept that the clusters you then find are clusters in the reduced space and have to be described there.

Noise is the point

The -1 label is what you are buying. k-means is obliged to place every point somewhere, so it has no vocabulary for “this customer resembles nobody”, and the rare cases end up quietly diluting whichever segment is nearest. DBSCAN says it plainly, which makes the noise set worth reading on its own: it is where the fraud, the test accounts, the bots and the data errors collect.

That honesty cuts both ways. A run that labels 60% of the data as noise is not a failed run — it is a result, and it means one of two things. Either eps is too small or min_samples too large, which the k-distance plot will show as a threshold sitting well below the elbow. Or there are no dense regions at that scale, and the data is a continuum with a few concentrations in it. The second is common in customer data, where nearly every variable is a smooth gradient.

Both apply to the run at the top, and the way to see it is to print what the ten clusters contain rather than how many there are.

print(X.assign(cluster=labels).query("cluster >= 0")
       .groupby("cluster").agg(n=("view", "size"), view=("view", "mean"),
                               login=("login", "mean"), export=("export", "mean"),
                               purchase=("purchase", "mean")).round(2))

           n  view  login  export  purchase
cluster
0        196  3.79   1.47     0.0       1.0
1        209  3.30   1.35     0.0       0.0
2        129  4.39   1.64     1.0       1.0
3         25  7.60   3.60     1.0       2.0
4        134  4.56   1.78     1.0       0.0
5         54  5.78   2.07     2.0       1.0
6         20  5.00   1.75     0.0       2.0
7         17  2.00   1.18     2.0       0.0
8         26  7.62   2.88     2.0       0.0
9         14  6.57   2.79     3.0       1.0

Every mean in the last two columns is a whole number, and that is not rounding: each cluster holds one exact pair of values. Cluster 1 is not a segment of low-engagement customers, it is 209 of the 218 rows with export = 0 and purchase = 0. Cluster 0 is export = 0, purchase = 1. Clusters 7 and 8 are both export = 2, purchase = 0, cut in half by view. The algorithm found the integer grid, because on columns whose useful range is 0 to 3 the densest regions are the lattice points themselves. No amount of eps tuning turns that into customer types. Bin such columns deliberately and call the result a crosstab, or cluster on something continuous.

Density-based labels still need the usual scrutiny, and border points are the place to apply it: they sit in no core neighbourhood of their own, so their assignment can flip between runs on resampled data, which makes the resampling stability check as necessary here as anywhere else. Report the noise fraction, the cluster sizes and the parameters together, because all three are part of the finding.

So make it a reporting habit: publish the noise fraction and both parameters next to the labels, every time. A reader who knows what share of points went unassigned, and at what eps, can judge the segmentation and argue with it. A reader handed four coloured wedges and four names has nothing to argue with.