DBSCAN vs HDBSCAN for Spatial Clustering

DBSCAN and HDBSCAN both find dense groups of points without being told how many clusters to expect, but they differ in one decisive way: DBSCAN uses a single fixed density threshold, while HDBSCAN adapts across densities. This guide compares them on spatial point data and shows when each is right — it is for anyone clustering crime incidents, wildlife sightings, or IoT sensor hits stored as a GeoPandas point layer. It sits under Spatial Clustering Algorithms in Spatial Analysis & Advanced Query Techniques, and once your points carry labels you can render the groups directly with clustering map markers in Folium.

Why This Approach / What Goes Wrong

DBSCAN classifies every point as one of three things. A core point has at least min_samples neighbours within eps; a border point falls inside a core point's eps neighbourhood but is not itself dense enough to be core; everything else is noise, labelled -1. Clusters grow by chaining core points that are density-reachable from one another. This makes DBSCAN fast, deterministic, and easy to reason about — but the whole model hinges on eps, and eps is a single global distance. On data that mixes a dense downtown with sparse rural tracts, no single eps fits: set it small enough to resolve the rural clusters and the downtown fragments into dozens of pieces; set it large enough to hold the downtown together and the rural points dissolve into noise or merge across streets that should separate them.

HDBSCAN removes eps altogether. It reweights the space by a mutual reachability distance (each point's distance inflated by how isolated its neighbourhood is), builds a minimum spanning tree over that metric, and condenses it into a hierarchy of clusters that appear and disappear as a density threshold sweeps from tight to loose. Instead of cutting that hierarchy at one fixed level, it extracts the clusters that persist over the widest range of thresholds — the most stable ones. The practical payoff is that a dense cluster and a sparse cluster can both be selected at the density that suits each, so mixed-density data that defeats DBSCAN comes out cleanly. You pay for it with runtime and with a probabilistic output: probabilities_ gives each point a soft membership strength, and noise points score near zero.

The mistake that sinks both algorithms is geographic coordinates. Every distance here is Euclidean by default, and Euclidean distance on longitude/latitude degrees is meaningless — a degree of longitude is ~111 km at the equator but shrinks toward the poles, so an eps in "degrees" silently distorts with latitude. Project to a metre-based CRS first (see Coordinate Reference System Transformations) so eps is an honest distance, and never cluster in Web Mercator (EPSG:3857), whose scale factor stretches with latitude just as badly. If you must cluster in lon/lat, pass metric="haversine" on coordinates converted to radians — never raw degrees with the default Euclidean metric.

Prerequisites

• geopandas>=0.14 — reads the point layer and manages the CRS
• scikit-learn>=1.4 — ships both DBSCAN and HDBSCAN under sklearn.cluster (HDBSCAN was merged into scikit-learn in 1.3, so no separate hdbscan package is needed)
• numpy>=1.26 — coordinate array handling and np.radians

conda install -c conda-forge "geopandas=0.14.*" "scikit-learn=1.4.*" "numpy=1.26.*"

Install from conda-forge rather than mixing pip wheels — GeoPandas binds the GDAL/PROJ C stack, and letting conda resolve one consistent build avoids the PROJ-path mismatches that produce silent CRS failures.

Step-by-Step Implementation

1. Load points and project to a metric CRS so distances are in metres, not degrees.

import geopandas as gpd
import numpy as np

# crime_incidents: point events across a metro region, stored in EPSG:4326
crime_incidents = gpd.read_file("crime_incidents.gpkg").to_crs(epsg=25832)
assert crime_incidents.crs.is_projected, "Cluster in a metric CRS, not degrees"

coords_m = np.column_stack([crime_incidents.geometry.x, crime_incidents.geometry.y])

EPSG:25832 (ETRS89 / UTM zone 32N) is a metre-based grid for central Europe; pick the UTM zone that matches your data, or derive it with crime_incidents.estimate_utm_crs(). The zone-selection mechanics live in choosing a UTM zone automatically in Python.

2. DBSCAN with an explicit metric eps (e.g. 250 m, ≥5 points to form a core).

from sklearn.cluster import DBSCAN

db = DBSCAN(eps=250, min_samples=5).fit(coords_m)   # eps is in metres
crime_incidents["dbscan"] = db.labels_              # -1 = noise

Because coords_m is projected, eps=250 means a literal 250-metre radius. That is the single lever DBSCAN gives you, and its value is the entire result — halving it can turn one cluster into ten.

Choosing it is not guesswork, and "try a few values and look at the map" is the slow way. Sort every point's distance to its k-th nearest neighbour (with k = min_samples) and plot the sorted curve: it stays flat while you are inside clusters and turns sharply upward at the distance where points stop having dense neighbourhoods. That knee is the eps the data is asking for. The curve is cheap — one k-nearest-neighbour query per point over the same tree DBSCAN would build anyway.

import numpy as np
from sklearn.neighbors import NearestNeighbors

MIN_SAMPLES = 5

# Distance from every incident to its 5th nearest neighbour, in metres
nn = NearestNeighbors(n_neighbors=MIN_SAMPLES).fit(coords_m)
distances, _ = nn.kneighbors(coords_m)
kth = np.sort(distances[:, -1])

# The knee sits where the curve's slope jumps; the percentile band brackets it
for q in (90, 95, 97, 99):
    print(f"{q}th percentile of 5-NN distance: {np.percentile(kth, q):>7.0f} m")
# 90th percentile of 5-NN distance:     186 m
# 95th percentile of 5-NN distance:     241 m
# 97th percentile of 5-NN distance:     318 m
# 99th percentile of 5-NN distance:     994 m

Read that output as a decision, not a number: the curve is flat through the 95th percentile and then climbs steeply, so anything between roughly 190 m and 250 m clusters the same dense core, and the run above at eps=250 sits at the top of that stable band. Past the 97th percentile the value is being set by the sparse tail — choose eps=994 and the isolated 1% starts chaining clusters together across the whole city. A min_samples change moves the whole curve, so redo the plot whenever you change it; the two parameters are not independent, and the usual starting point is min_samples = 2 × dimensions, which for planar coordinates means 4 or 5.

3. HDBSCAN with no eps — only a minimum cluster size.

from sklearn.cluster import HDBSCAN

hdb = HDBSCAN(min_cluster_size=15).fit(coords_m)
crime_incidents["hdbscan"] = hdb.labels_              # -1 = noise
crime_incidents["hdbscan_prob"] = hdb.probabilities_  # 0..1 membership strength

min_cluster_size is the smallest group HDBSCAN will accept as one cluster — a semantic knob ("what counts as a hotspot?") rather than a distance. Keep probabilities_: it lets you keep only confident members (e.g. prob > 0.5) and discard fringe points without re-running the fit.

Three further parameters decide most real outcomes, and skipping them is why people conclude HDBSCAN "does not work" on their data. min_samples — separate from min_cluster_size and defaulting to it — controls how conservative the density estimate is: raise it and more points are declared noise, which is the direct lever when a run labels too much of a sparse region as belonging to a dense group. cluster_selection_method chooses how the hierarchy is cut: the default "eom" (excess of mass) prefers a few large, stable clusters, while "leaf" takes the finest branches and returns many small homogeneous ones — the right choice when the question is "where are the individual hotspots?" rather than "what are the broad zones?". cluster_selection_epsilon is the hybrid escape hatch: it sets a distance floor below which clusters are merged rather than split, which is how you tell HDBSCAN that two hotspots 40 metres apart are one place regardless of what the hierarchy says.

from sklearn.cluster import HDBSCAN

# Finer-grained hotspots, but nothing closer than 50 m is ever split in two
hdb_fine = HDBSCAN(
    min_cluster_size=15,
    min_samples=8,                    # stricter density -> more honest noise
    cluster_selection_method="leaf",  # many small clusters instead of a few broad ones
    cluster_selection_epsilon=50.0,   # METRES, because coords_m is projected
).fit(coords_m)

crime_incidents["hotspot"] = hdb_fine.labels_
print("hotspots:", hdb_fine.labels_.max() + 1)   # hotspots: 141

Note that cluster_selection_epsilon is in the units of the coordinate array, so it inherits the same metric-CRS requirement as DBSCAN's eps — it is the one HDBSCAN parameter that will silently misbehave on degrees.

4. Cluster lon/lat directly with the Haversine metric when reprojecting is undesirable — pass radians, not degrees.

from sklearn.cluster import DBSCAN

# Work straight from EPSG:4326 without projecting: Haversine expects RADIANS.
ll = crime_incidents.to_crs(epsg=4326)
coords_rad = np.radians(np.column_stack([ll.geometry.y, ll.geometry.x]))  # lat, lon

eps_km = 0.25
earth_radius_km = 6371.0088
db_hav = DBSCAN(eps=eps_km / earth_radius_km, min_samples=5,
                metric="haversine").fit(coords_rad)
crime_incidents["dbscan_haversine"] = db_hav.labels_

The Haversine eps is an angular distance in radians, so divide your target ground distance by Earth's radius. Full coordinate-system handling with PyProj covers when a true projection beats the great-circle approximation (large extents, high latitudes).

5. Summarize how many groups each algorithm found and compare.

def summarize(labels):
    n_clusters = len(set(labels)) - (1 if -1 in labels else 0)
    n_noise = int((labels == -1).sum())
    return n_clusters, n_noise

print("DBSCAN :", summarize(crime_incidents["dbscan"]))   # DBSCAN : (38, 5120)
print("HDBSCAN:", summarize(crime_incidents["hdbscan"]))  # HDBSCAN: (52, 3987)

On mixed-density incident data HDBSCAN typically recovers more clusters with less noise, because it resolves sparse groups a single global eps would have thrown away.

6. Know where each one stops. Both algorithms are index-backed and near O(n log n) in time on well-behaved planar data, so runtime is rarely what breaks first — memory is, and it breaks for different reasons on each side. DBSCAN materialises neighbourhoods: with a small eps that is a handful of indices per point, but the count grows with the area of the neighbourhood, so doubling eps roughly quadruples what has to be held. A run that is comfortable at 250 m on two million points can exhaust memory at 2 km on the same array, and the traceback names an internal radius_neighbors call rather than your parameter. HDBSCAN's cost is structural instead of parametric: it builds a minimum spanning tree over the whole dataset, which is why it is consistently several times slower than DBSCAN and why its memory is set by the point count rather than by any setting you chose.

The practical envelope on a normal workstation: below a million points either algorithm runs in seconds to a couple of minutes and the choice is purely about which result is more truthful. Between one and roughly ten million, DBSCAN with a tuned eps and n_jobs=-1 stays comfortable while HDBSCAN starts to need real patience and headroom. Above that, neither is the right tool applied to the full array at once — cluster a spatially coherent subset, or tile the extent and merge, as covered under Spatial Clustering Algorithms. Fit on a random subsample to explore parameters, then run once on the full set: the k-distance curve and the stable eps band are properties of the density, and a 10% sample reproduces both closely enough to tune against.

Version differences matter here too. scikit-learn's in-tree HDBSCAN (1.3 and later) gives you labels_ and probabilities_ and little else. The standalone hdbscan package adds the pieces you may actually need in production: approximate_predict() to assign new points to an existing clustering without refitting — the difference between a nightly rebuild and a streaming pipeline — plus GLOSH outlier scores and a density-based validity index for judging a result without ground truth. Both are worth knowing about before you commit to the scikit-learn-only path; DBSCAN, by contrast, has no such split and behaves identically across recent releases.

Verification

Confirm clustering ran in metric space and that labels are internally consistent.

assert crime_incidents.crs.is_projected, "Cluster in a metric CRS, not degrees"

# Every point is labelled (cluster id or -1 noise), none left unassigned
assert crime_incidents["dbscan"].notna().all()
assert crime_incidents["hdbscan"].notna().all()

# HDBSCAN noise points carry ~0 membership probability
noise_mask = crime_incidents["hdbscan"] == -1
assert crime_incidents.loc[noise_mask, "hdbscan_prob"].max() < 1e-6
print("Both clusterings labelled all", len(crime_incidents), "points")
# Both clusterings labelled all 12000 points

For a spatial sanity check, dissolve every cluster to its convex hull and confirm the hulls do not wildly overlap — heavy overlap usually means eps (DBSCAN) or min_cluster_size (HDBSCAN) is too permissive.

Edge Cases & Debugging

• Everything is one giant cluster. eps too large for DBSCAN; lower it, or switch to HDBSCAN if density genuinely varies across the extent.
• Everything is noise. eps too small or min_samples too high; relax them, or drop min_cluster_size for HDBSCAN.
• Clusters make no geographic sense. You clustered in degrees or in EPSG:3857 — reproject to a UTM/national metric CRS first so eps is a real distance.
• HDBSCAN slower than expected. It is heavier than DBSCAN on large N; subsample for exploration, then run once on the full set.
• Need lat/lon directly. Use metric="haversine" with coordinates in radians (np.radians) and an angular eps, not Euclidean on degrees.
• Unstable cluster ids across runs. Cluster numbering is arbitrary and not comparable between fits; join on geometry or a stable attribute, never on the raw label integer.
• MemoryError inside radius_neighbors. DBSCAN's eps is too large for the point density, so the neighbourhood lists no longer fit. Lower eps to the k-distance knee, or switch to HDBSCAN, whose memory does not scale with a radius.
• Two clusters that a human would call one place. Duplicate or near-duplicate coordinates — several incidents geocoded to the same building centroid — inflate local density and split a group. De-duplicate on rounded coordinates before fitting, or set cluster_selection_epsilon to the distance below which you consider points the same place.
• A single point becomes its own cluster. Only possible with min_samples=1, which makes every point a core point and defeats the noise concept entirely. Keep min_samples at 4 or above for planar coordinates.
• HDBSCAN and DBSCAN disagree completely on the same data. Usually correct behaviour, not a bug: the disagreement is the mixed-density signal. Compare the k-distance curve against the run — if it has no clear knee, DBSCAN has no honest eps and HDBSCAN's answer is the one to trust.
• Labels change when you re-read the file. The row order changed (a different driver, a different sort), and label integers follow discovery order. The partition is stable even though the numbering is not; compare clusterings by membership, not by id.

Frequently Asked Questions

Which one should I default to on a dataset I have not seen before? Run HDBSCAN first. It needs one semantic parameter instead of a distance you have not calibrated yet, and its output tells you something DBSCAN's cannot: whether the data has one density scale or several. If the resulting clusters all turn out to sit at a similar scale, switch to DBSCAN for the production run — it is faster, its single parameter is auditable, and reviewers can reason about a radius in metres.

Is HDBSCAN always better on mixed-density data? It is better at recovering groups a single radius would miss, which is not the same as being right for your question. If the analysis has a fixed operational meaning — "premises within 250 m of each other count as one incident site" — that radius is part of the definition, and DBSCAN encodes it exactly while HDBSCAN will happily merge or split against it. An adaptive algorithm is the wrong choice when the threshold is a policy rather than a discovery.

How many noise points is too many? There is no universal figure, but noise carries information: a run labelling 60% of points as noise on data you believe is clustered means the density parameters describe a different phenomenon than the one present. Check the k-distance curve first. Genuinely diffuse point processes — random calls for service, evenly spread sensors — should come back almost entirely as noise, and forcing clusters out of them by relaxing parameters produces groups that will not reproduce on next month's data.

Can I cluster on attributes as well as location? Yes, by appending scaled attribute columns to the coordinate array — but the scaling is the whole problem. Euclidean distance over a matrix of metres and, say, decibels is meaningless until you decide how many metres one decibel is worth, and that conversion factor silently determines the result. Do it deliberately with an explicit multiplier per column, never with a blanket standardisation that assigns weights by whatever the variances happen to be.

Do I need to worry about determinism? DBSCAN's cluster membership is deterministic apart from border points, which may attach to whichever core cluster reaches them first and so can flip with row order. HDBSCAN is deterministic given the same input array. Neither guarantees stable label integers between runs, so persist a join key derived from geometry — a centroid or convex hull of the grouped points — rather than the label itself.

Should I use the Haversine metric or reproject? Reproject whenever the extent fits a single UTM zone or national grid: the projected run is faster, and every parameter reads in metres. Reach for Haversine on genuinely global or high-latitude extents where no single projection is honest across the data — and remember it costs you a ball_tree implementation and radian units in exchange.