PythonMastery
intermediate 22 min read · lesson 9 of 9 in Data Science & ML

Clustering: Unsupervised Pattern Discovery

1 · The lesson

read

In supervised learning, the data has labels, so you know what to predict. In unsupervised learning there are no labels. Your job is to find structure in the raw data: which customers behave alike, which transactions look unusual, which documents cover similar topics. Clustering is the main tool when someone says "we have a database, so tell us what's in it".

This lesson covers the three canonical algorithms (K-Means, DBSCAN, Hierarchical), the dimensionality reduction methods you almost always need first (PCA, UMAP), and the hardest evaluation question: "do these clusters mean anything to the business?"

Runtime note: scikit-learn and matplotlib run right here. Press Run on any block. umap-learn and hdbscan do not have WebAssembly builds, so the few snippets using those need pip install umap-learn hdbscan on your own machine. Their expected output is in comments.


1. Unsupervised: No Labels, Just Geometry

You have X of shape (n_samples, n_features). No y. Two big families:

FamilyOutputExamples
ClusteringAssignment of each row to a groupK-Means, DBSCAN, hierarchical, GMM
Dimensionality reductionLower-dimensional embedding of each rowPCA, t-SNE, UMAP, autoencoders

In practice you use them together: reduce dimensions first (so distances are meaningful), then cluster the embedding. If you skip the reduction on high-dimensional data, you hit the curse of dimensionality: every point is roughly the same distance from every other point, so distances stop being useful.


2. K-Means: The Default

The algorithm in five steps:

1. Pick K random points as initial centroids.
2. Assign every data point to its nearest centroid.
3. Recompute each centroid as the mean of its assigned points.
4. Repeat 2-3 until assignments stop changing.
5. Output: K centroid coordinates + cluster label for each row.

python
import numpy as np
from sklearn.cluster import KMeans
from sklearn.preprocessing import StandardScaler

# Customer data: (income, age, spend)
X = np.array([
    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],
    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],
    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],
])

# Step 1: always scale before clustering based on distance
X_scaled = StandardScaler().fit_transform(X)

# Step 2: fit
km = KMeans(n_clusters=3, n_init=10, random_state=42)
km.fit(X_scaled)

print(km.labels_)                      # [0 0 0 1 1 1 2 2 2]
print(km.cluster_centers_)             # 3 centroids in scaled space
print(km.inertia_)                     # within-cluster sum of squared distances
print(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster

Key points:

  • n_init=10: K-Means is sensitive to where the first centroids are placed, so sklearn runs 10 random inits and keeps the best. (sklearn ≥ 1.4 defaults to "auto", which with k-means++ runs only 1 init, so pass n_init=10 yourself when you want several.)
  • inertia_: the sum of squared distances from each point to its centroid. Lower means "tighter" clusters. It almost always decreases as K grows, which is why you need the elbow method below.
  • cluster_centers_ are in scaled space. Invert the scaler if you want to report them in the original units.

3. Choosing K: Elbow and Silhouette

K-Means makes you choose K before it runs. Two checks help you choose:

Elbow method: plot inertia_ against K. The "elbow", where the curve bends, is your candidate K.

python
import matplotlib.pyplot as plt

inertias = []
ks = range(1, min(11, len(X_scaled)))     # K can't exceed the number of rows (9 in this toy set)
for k in ks:
    km = KMeans(n_clusters=k, n_init=10, random_state=42).fit(X_scaled)
    inertias.append(km.inertia_)

plt.plot(ks, inertias, "o-")
plt.xlabel("K"); plt.ylabel("inertia")
plt.show()
# Look for the "elbow", where more clusters stop helping much
+ setup added so this can run · re-runs 1 earlier block on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_scaled = _AutoMock('X_scaled')
def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

Silhouette score: for each point, (b - a) / max(a, b) where a is mean distance to its own cluster and b is mean distance to the nearest other cluster. Ranges from -1 (wrong cluster) to +1 (well-separated). Average over all points.

python
from sklearn.metrics import silhouette_score

for k in range(2, 8):
    km = KMeans(n_clusters=k, n_init=10, random_state=42).fit(X_scaled)
    s = silhouette_score(X_scaled, km.labels_)
    print(f"K={k}  silhouette={s:.3f}")
# Pick the K with the highest silhouette
+ setup added so this can run · re-runs 1 earlier block on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_scaled = _AutoMock('X_scaled')
def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

The elbow method is faster, but the silhouette score has a clearer mathematical basis. In practice, run both. When they agree, you can be more confident in your K.


4. Where K-Means Fails

K-Means assumes clusters are spherical and roughly equal-sized. It gives poor results on:

  • Non-spherical shapes: two concentric rings, or two crescents (the classic make_moons dataset).
  • Very different cluster sizes: one giant blob and one tiny one.
  • Different densities: a sparse outer cluster and dense inner cluster.
python
from sklearn.datasets import make_moons
X_moons, _ = make_moons(n_samples=300, noise=0.05, random_state=42)
km = KMeans(n_clusters=2, n_init=10, random_state=42).fit(X_moons)
# K-Means cuts the moons vertically down the middle, which is completely wrong.
+ setup added so this can run · re-runs 1 earlier block on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

When the clusters are not spherical, use DBSCAN.


5. DBSCAN: Density-Based, Finds Any Shape

DBSCAN does not need K in advance. Instead it groups points that are densely packed and labels sparse outliers as noise.

Two hyperparameters:

  • eps: the neighbourhood radius. Points within eps of each other are "neighbours".
  • min_samples: a "core point" must have at least this many neighbours.

A cluster is a connected component of core points plus the non-core points within eps of any core. Points in neither category get label -1 (noise).

python
from sklearn.cluster import DBSCAN
from sklearn.preprocessing import StandardScaler

X_scaled = StandardScaler().fit_transform(X_moons)
db = DBSCAN(eps=0.3, min_samples=5).fit(X_scaled)
print(db.labels_)                       # 0, 1, ..., or -1 for noise
print(set(db.labels_.tolist()))         # {0, 1}  two clusters and no noise points at this eps
+ setup added so this can run · re-runs 2 earlier blocks on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_moons = _AutoMock('X_moons')
def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
    from sklearn.datasets import make_moons
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
    'from sklearn.datasets import make_moons\nX_moons, _ = make_moons(n_samples=300, noise=0.05, random_state=42)\nkm = KMeans(n_clusters=2, n_init=10, random_state=42).fit(X_moons)\n# K-Means cuts the moons vertically down the middle, which is completely wrong.',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

When DBSCAN wins:

  • Arbitrary cluster shapes (moons, rings, blobs of any geometry).
  • You want noise detection without extra work: outliers get label = -1.
  • You don't know K up front.

When it loses:

  • Very different densities: one eps can't fit dense and sparse clusters together (HDBSCAN solves this).
  • Very high dimensions: eps becomes meaningless (reduce dimensions first).

Tuning eps: compute the k-distance graph (sklearn.neighbors.NearestNeighbors), sort, plot. The knee is your eps. min_samples = 2 × n_features is a common rule of thumb.


6. Hierarchical (Agglomerative) Clustering

Build a tree of clusters by repeatedly merging the two closest clusters until everything is one big cluster. The output is a dendrogram (a tree diagram). Cut it at any height to get a clustering.

python
from sklearn.cluster import AgglomerativeClustering
from scipy.cluster.hierarchy import linkage, dendrogram
import matplotlib.pyplot as plt

ac = AgglomerativeClustering(n_clusters=3, linkage="ward").fit(X_scaled)
print(ac.labels_)

# Visualise the merge history
Z = linkage(X_scaled, method="ward")
dendrogram(Z)
plt.show()
+ setup added so this can run · re-runs 3 earlier blocks on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_scaled = _AutoMock('X_scaled')
def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')
X_moons = _AutoMock('X_moons')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
    from sklearn.datasets import make_moons
    from sklearn.cluster import DBSCAN
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
    'from sklearn.datasets import make_moons\nX_moons, _ = make_moons(n_samples=300, noise=0.05, random_state=42)\nkm = KMeans(n_clusters=2, n_init=10, random_state=42).fit(X_moons)\n# K-Means cuts the moons vertically down the middle, which is completely wrong.',
    'from sklearn.cluster import DBSCAN\nfrom sklearn.preprocessing import StandardScaler\n\nX_scaled = StandardScaler().fit_transform(X_moons)\ndb = DBSCAN(eps=0.3, min_samples=5).fit(X_scaled)\nprint(db.labels_)                       # 0, 1, ..., or -1 for noise\nprint(set(db.labels_.tolist()))         # {0, 1}  two clusters and no noise points at this eps',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

linkage choices:

  • "ward": minimise the increase in variance at each merge. Best general default.
  • "average": mean pairwise distance between clusters.
  • "complete": max pairwise distance.
  • "single": min pairwise distance. It tends to produce long, thin "chained" clusters.

Hierarchical is interpretable (you can see the merge order) but it needs O(n²) memory, so it doesn't scale well past ~10⁴ rows. For larger data, use K-Means or DBSCAN.


7. Feature Scaling Matters (a Lot)

In distance-based methods (K-Means, DBSCAN, hierarchical), features with large numeric ranges dominate the distance. If income is in 0–200,000 and age is in 18–80, the distance is almost entirely income. Always scale before clustering.

See the feature engineering lesson: StandardScaler for Gaussian-like data, MinMaxScaler for bounded data, RobustScaler when outliers exist.

python
from sklearn.preprocessing import StandardScaler
X_scaled = StandardScaler().fit_transform(X)
+ setup added so this can run · re-runs 1 earlier block on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X = _AutoMock('X')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

This one line often decides whether the clusters mean something or are useless.


8. Dimensionality Reduction: Why and How

In high dimensions, every pair of points is roughly equidistant (the curse of dimensionality) and visualisation is impossible. Reduce first:

MethodLinear?Use caseCaveat
PCAYesCompression, preprocessing for clustering, baselineCaptures variance, not necessarily structure
t-SNENo2D visualisation onlyDistorts global structure, slow, don't feed downstream
UMAPNoVisualisation + preprocessingFaster than t-SNE, preserves more global structure
AutoencodersNo (NN)Very high-D, large dataNeeds a GPU and a lot of tuning

PCA

python
from sklearn.decomposition import PCA

pca = PCA(n_components=2).fit(X_scaled)
X_2d = pca.transform(X_scaled)
print(pca.explained_variance_ratio_.round(2))    # [0.74 0.25]  so 2D keeps about 99% of the variance (the toy data has only 3 features)
print(pca.explained_variance_ratio_.sum())
+ setup added so this can run · re-runs 4 earlier blocks on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_scaled = _AutoMock('X_scaled')
def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')
X_moons = _AutoMock('X_moons')
X = _AutoMock('X')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
    from sklearn.datasets import make_moons
    from sklearn.cluster import DBSCAN
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
    'from sklearn.datasets import make_moons\nX_moons, _ = make_moons(n_samples=300, noise=0.05, random_state=42)\nkm = KMeans(n_clusters=2, n_init=10, random_state=42).fit(X_moons)\n# K-Means cuts the moons vertically down the middle, which is completely wrong.',
    'from sklearn.cluster import DBSCAN\nfrom sklearn.preprocessing import StandardScaler\n\nX_scaled = StandardScaler().fit_transform(X_moons)\ndb = DBSCAN(eps=0.3, min_samples=5).fit(X_scaled)\nprint(db.labels_)                       # 0, 1, ..., or -1 for noise\nprint(set(db.labels_.tolist()))         # {0, 1}  two clusters and no noise points at this eps',
    'from sklearn.preprocessing import StandardScaler\nX_scaled = StandardScaler().fit_transform(X)',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

explained_variance_ratio_ tells you how much variance each component captures. Rule of thumb: keep enough components for ~95% cumulative variance.

PCA is linear: it finds orthogonal directions of maximum variance. It is fast and deterministic (same input, same result), so it is almost always worth trying first.

t-SNE: Visualisation Only

python
from sklearn.manifold import TSNE

X_2d = TSNE(n_components=2, perplexity=min(30, len(X_scaled) - 1),   # perplexity must be < number of rows
            random_state=42).fit_transform(X_scaled)
+ setup added so this can run · re-runs 4 earlier blocks on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_scaled = _AutoMock('X_scaled')
def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')
X_moons = _AutoMock('X_moons')
X = _AutoMock('X')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
    from sklearn.datasets import make_moons
    from sklearn.cluster import DBSCAN
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
    'from sklearn.datasets import make_moons\nX_moons, _ = make_moons(n_samples=300, noise=0.05, random_state=42)\nkm = KMeans(n_clusters=2, n_init=10, random_state=42).fit(X_moons)\n# K-Means cuts the moons vertically down the middle, which is completely wrong.',
    'from sklearn.cluster import DBSCAN\nfrom sklearn.preprocessing import StandardScaler\n\nX_scaled = StandardScaler().fit_transform(X_moons)\ndb = DBSCAN(eps=0.3, min_samples=5).fit(X_scaled)\nprint(db.labels_)                       # 0, 1, ..., or -1 for noise\nprint(set(db.labels_.tolist()))         # {0, 1}  two clusters and no noise points at this eps',
    'from sklearn.preprocessing import StandardScaler\nX_scaled = StandardScaler().fit_transform(X)',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

t-SNE keeps local neighbourhoods well, but it distorts global distances: the gap between two t-SNE clusters means nothing. Never use t-SNE output as features for a downstream model. It is only a plotting tool.

UMAP

python
import umap

reducer = umap.UMAP(n_components=2, n_neighbors=15, min_dist=0.1, random_state=42)
X_2d = reducer.fit_transform(X_scaled)
+ setup added so this can run · re-runs 4 earlier blocks on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_scaled = _AutoMock('X_scaled')
def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')
X_moons = _AutoMock('X_moons')
X = _AutoMock('X')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
    from sklearn.datasets import make_moons
    from sklearn.cluster import DBSCAN
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
    'from sklearn.datasets import make_moons\nX_moons, _ = make_moons(n_samples=300, noise=0.05, random_state=42)\nkm = KMeans(n_clusters=2, n_init=10, random_state=42).fit(X_moons)\n# K-Means cuts the moons vertically down the middle, which is completely wrong.',
    'from sklearn.cluster import DBSCAN\nfrom sklearn.preprocessing import StandardScaler\n\nX_scaled = StandardScaler().fit_transform(X_moons)\ndb = DBSCAN(eps=0.3, min_samples=5).fit(X_scaled)\nprint(db.labels_)                       # 0, 1, ..., or -1 for noise\nprint(set(db.labels_.tolist()))         # {0, 1}  two clusters and no noise points at this eps',
    'from sklearn.preprocessing import StandardScaler\nX_scaled = StandardScaler().fit_transform(X)',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

UMAP is faster than t-SNE, keeps more global structure, and can be used as input to downstream models (with care). In many modern data science notebooks it has replaced t-SNE.

The standard pipeline for high-dim data: scale → UMAP to 5–10 dims → cluster. PCA-then-cluster is the conservative alternative.


9. Evaluating Clustering

Without labels, you have three options:

1. Internal metrics (cluster-shape only):

python
from sklearn.metrics import silhouette_score, davies_bouldin_score

print(silhouette_score(X_scaled, labels))        # higher is better, in [-1, 1]
print(davies_bouldin_score(X_scaled, labels))    # lower is better
+ setup added so this can run · defines X_scaled, labels
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_scaled = _AutoMock('X_scaled')
labels = _AutoMock('labels')

2. Stability: re-cluster a bootstrap sample and measure how well the labels agree (adjusted Rand index). If the clusters are stable, they probably reflect real structure, not chance.

3. Business validation: in the end, this is the one that matters. Hand the clusters to a domain expert. Do they correspond to recognisable customer types? Do the cluster centroids tell a coherent story? Do the actions you'd take differ per cluster?

A clustering with excellent scores is worthless if nobody in the business can use it. A clustering with average scores is very valuable if it finds a clear "high-value but at-risk" segment that marketing can target.


10. Real-World Uses

  • Customer segmentation: group users by behaviour, target each segment differently.
  • Anomaly detection: cluster, then flag points far from any centroid (or DBSCAN noise points).
  • Document grouping: embed with sentence-transformers, then UMAP + HDBSCAN. The basis of topic modelling.
  • Image compression: cluster pixel colours with K-Means, replace each pixel with its centroid colour. This is the classic K=16 colour quantisation.
  • Recommendation pre-filtering: cluster items, then pick candidates from the same cluster as a fast first step before ranking.

11. Soft Clustering: Gaussian Mixture Models

K-Means gives every point a hard label. GMM gives every point a probability distribution over clusters, for example "67% cluster A, 33% cluster B". This is useful when the boundaries are not clear (a customer who buys both luxury and budget products).

python
from sklearn.mixture import GaussianMixture

gmm = GaussianMixture(n_components=3, random_state=42).fit(X_scaled)
print(gmm.predict(X_scaled))           # hard labels
print(gmm.predict_proba(X_scaled))     # soft probabilities (rows sum to 1)
+ setup added so this can run · re-runs 4 earlier blocks on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_scaled = _AutoMock('X_scaled')
def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')
X_moons = _AutoMock('X_moons')
X = _AutoMock('X')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
    from sklearn.datasets import make_moons
    from sklearn.cluster import DBSCAN
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
    'from sklearn.datasets import make_moons\nX_moons, _ = make_moons(n_samples=300, noise=0.05, random_state=42)\nkm = KMeans(n_clusters=2, n_init=10, random_state=42).fit(X_moons)\n# K-Means cuts the moons vertically down the middle, which is completely wrong.',
    'from sklearn.cluster import DBSCAN\nfrom sklearn.preprocessing import StandardScaler\n\nX_scaled = StandardScaler().fit_transform(X_moons)\ndb = DBSCAN(eps=0.3, min_samples=5).fit(X_scaled)\nprint(db.labels_)                       # 0, 1, ..., or -1 for noise\nprint(set(db.labels_.tolist()))         # {0, 1}  two clusters and no noise points at this eps',
    'from sklearn.preprocessing import StandardScaler\nX_scaled = StandardScaler().fit_transform(X)',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

GMM is more flexible than K-Means (clusters can be ellipses, not just spheres) and it tells you how uncertain each assignment is. The cost is more parameters and slower fitting.


12. HDBSCAN: DBSCAN's Modern Successor

DBSCAN uses one eps for the whole dataset, so it fails when the density varies. HDBSCAN runs DBSCAN at all eps scales and extracts the most stable clusters.

python
import hdbscan

clusterer = hdbscan.HDBSCAN(min_cluster_size=10, min_samples=5)
labels = clusterer.fit_predict(X_scaled)
print(labels)                           # cluster IDs, -1 for noise
print(clusterer.probabilities_)         # how strongly each point belongs to its cluster
+ setup added so this can run · re-runs 4 earlier blocks on this page first
# Lightweight mock for objects whose attributes/methods aren't critical
class _AutoMock:
    def __init__(self, name='mock'): self._name = name
    def __getattr__(self, k): return _AutoMock(self._name + '.' + k)
    def __call__(self, *a, **kw):
        print('-> ' + self._name + '() called')
        return _AutoMock(self._name + '()')
    def __repr__(self): return '<mock ' + self._name + '>'
    def __str__(self): return '<mock ' + self._name + '>'
    def __bool__(self): return True
    def __iter__(self): return iter([])
    def __len__(self): return 0
    def __getitem__(self, k): return _AutoMock(self._name + '[...]')
    def __setitem__(self, k, v): pass
    def __enter__(self): return self
    def __exit__(self, *a): return False
    async def __aenter__(self): return self
    async def __aexit__(self, *a): return False
    def __add__(self, o): return self
    def __radd__(self, o): return self
    def __sub__(self, o): return self
    def __mul__(self, o): return self
    def __rmul__(self, o): return self
    def __truediv__(self, o): return self
    def __eq__(self, o): return isinstance(o, _AutoMock)
    def __hash__(self): return hash(self._name)
    def __lt__(self, o): return True
    def __le__(self, o): return True
    def __gt__(self, o): return False
    def __ge__(self, o): return False
    def __mro_entries__(self, bases): return (object,)

X_scaled = _AutoMock('X_scaled')
def KMeans(*_a, **_kw):
    print('-> KMeans() called')
    return _AutoMock('KMeans()')
X_moons = _AutoMock('X_moons')
X = _AutoMock('X')

# Earlier blocks on this page that this one builds on, re-run quietly first.
try:
    import numpy as np
    from sklearn.cluster import KMeans
    from sklearn.preprocessing import StandardScaler
    from sklearn.datasets import make_moons
    from sklearn.cluster import DBSCAN
except Exception:
    pass
_pm_earlier = [
    'import numpy as np\nfrom sklearn.cluster import KMeans\nfrom sklearn.preprocessing import StandardScaler\n\n# Customer data: (income, age, spend)\nX = np.array([\n    [ 30_000, 22,  500], [ 35_000, 25,  600], [ 28_000, 21,  450],\n    [ 80_000, 45, 2_500], [ 85_000, 50, 2_800], [ 75_000, 42, 2_300],\n    [150_000, 35, 8_000], [160_000, 38, 9_000], [145_000, 33, 7_500],\n])\n\n# Step 1: always scale before clustering based on distance\nX_scaled = StandardScaler().fit_transform(X)\n\n# Step 2: fit\nkm = KMeans(n_clusters=3, n_init=10, random_state=42)\nkm.fit(X_scaled)\n\nprint(km.labels_)                      # [0 0 0 1 1 1 2 2 2]\nprint(km.cluster_centers_)             # 3 centroids in scaled space\nprint(km.inertia_)                     # within-cluster sum of squared distances\nprint(km.predict([[0.5, 0.5, 0.5]]))   # assign a new (scaled) point to a cluster',
    'from sklearn.datasets import make_moons\nX_moons, _ = make_moons(n_samples=300, noise=0.05, random_state=42)\nkm = KMeans(n_clusters=2, n_init=10, random_state=42).fit(X_moons)\n# K-Means cuts the moons vertically down the middle, which is completely wrong.',
    'from sklearn.cluster import DBSCAN\nfrom sklearn.preprocessing import StandardScaler\n\nX_scaled = StandardScaler().fit_transform(X_moons)\ndb = DBSCAN(eps=0.3, min_samples=5).fit(X_scaled)\nprint(db.labels_)                       # 0, 1, ..., or -1 for noise\nprint(set(db.labels_.tolist()))         # {0, 1}  two clusters and no noise points at this eps',
    'from sklearn.preprocessing import StandardScaler\nX_scaled = StandardScaler().fit_transform(X)',
]
def _pm_run_earlier(_srcs, _g):
    import sys as _s, io as _io, warnings as _w
    for _src in _srcs:
        _o, _e = _s.stdout, _s.stderr
        _s.stdout = _s.stderr = _io.StringIO()
        try:
            with _w.catch_warnings():
                _w.simplefilter('ignore')
                exec(compile(_src, '<earlier block>', 'exec'), _g)
        except (Exception, SystemExit):
            pass
        finally:
            _s.stdout, _s.stderr = _o, _e
    if 'matplotlib.pyplot' in _s.modules:
        _s.modules['matplotlib.pyplot'].close('all')   # their plots belong to their own blocks
_pm_run_earlier(_pm_earlier, globals())
del _pm_earlier, _pm_run_earlier

If your data has clusters of clearly different densities and you would otherwise use DBSCAN, try HDBSCAN first. Document-clustering and topic-modelling tools often use it by default (BERTopic, for example).


Common Mistakes

1. Not scaling features. A feature on 0–1,000,000 dominates a feature on 0–1 in any Euclidean distance. Always apply StandardScaler (or a similar scaler) before clustering. This is the most common beginner error.

2. Choosing K by intuition. "I think there are about 5 segments" is a guess, not a method. Use the elbow method and silhouette score, and prefer the K where both agree.

3. Treating cluster IDs as ordinal. Cluster 2 is not "between" cluster 1 and cluster 3. The labels are arbitrary. Don't sort, regress against, or compute means of cluster IDs. One-hot encode them if they are inputs to a downstream model.

4. Clustering raw text or images. A bag-of-words vector or a raw pixel array gives very poor clusters. Embed first: sentence-transformers for text, a pretrained CNN/ViT for images. Then cluster the embedding. UMAP + HDBSCAN on embeddings is the usual modern approach.

5. Using t-SNE coordinates downstream. t-SNE distorts global structure, so the (x, y) coordinates are not features. Use t-SNE for plots, and use PCA, UMAP, or the raw scaled features for downstream models.

6. Reading too much into K-Means in high dimensions. Above ~50 features, Euclidean distance is unreliable. PCA/UMAP to 5–20 dimensions first, then cluster.

7. Reporting "X% silhouette" as if it's accuracy. Silhouette in [-1, 1] is not a percentage, and 0.4 is good for many real datasets. Compare across K values, not against an absolute threshold.


🎯 Your Turn: Customer Segmentation with K-Means

Cluster the inline (income, age, spend) customer data. Choose K via the elbow method on inertia_, then return:

1. The chosen K.
2. The cluster centroids in original (unscaled) units.
3. The count of customers in each cluster.

Requirements:

  • StandardScaler before clustering.
  • Try K = 2 through 7. Find the first K where the relative drop in inertia (from K-1 to K) is below 30%, and pick the K just before it.
  • Use KMeans(n_init=10, random_state=42) for reproducibility.

Skeleton:

python
import numpy as np
from sklearn.cluster import KMeans
from sklearn.preprocessing import StandardScaler

# (income, age, monthly_spend)
X = np.array([
    [ 28_000, 22,  300], [ 32_000, 24,  400], [ 25_000, 21,  250],
    [ 30_000, 23,  350], [ 27_000, 22,  280], [ 33_000, 25,  420],
    [ 70_000, 40, 1_800], [ 75_000, 42, 2_000], [ 68_000, 38, 1_700],
    [ 72_000, 41, 1_900], [ 78_000, 44, 2_100], [ 71_000, 39, 1_750],
    [140_000, 35, 6_000], [155_000, 38, 6_500], [148_000, 36, 6_200],
    [160_000, 40, 7_000], [145_000, 34, 5_800], [152_000, 37, 6_400],
    [ 95_000, 55, 3_200], [ 98_000, 58, 3_400], [102_000, 60, 3_500],
    [ 92_000, 53, 3_000], [ 99_000, 57, 3_300], [105_000, 62, 3_600],
])

def segment(X):
    # TODO 1: scale X with StandardScaler
    # TODO 2: fit KMeans for K from 2 to 7, store the inertias
    # TODO 3: pick the first K where the relative drop in inertia is below 30%
    # TODO 4: refit KMeans with chosen K
    # TODO 5: invert the scaler on cluster_centers_ to report in original units
    # TODO 6: return (K, centroids_original_units, counts_per_cluster)
    ...

print(segment(X))
Hint 1: Relative inertia drop For each consecutive pair of K values, compute (inertias[k-1] - inertias[k]) / inertias[k-1]. If this fraction is below 0.30 at some K for the first time, then K-1 is your elbow: adding the K-th cluster no longer helps much. Add a check so you don't go below K=2.
Hint 2: Unscaling centroids StandardScaler stores mean_ and scale_. To invert: centroids_original = scaler.inverse_transform(km.cluster_centers_). For counts, np.unique(km.labels_, return_counts=True) gives both unique labels and their counts.
Show full solution
python
import numpy as np
from sklearn.cluster import KMeans
from sklearn.preprocessing import StandardScaler

X = np.array([
    [ 28_000, 22,  300], [ 32_000, 24,  400], [ 25_000, 21,  250],
    [ 30_000, 23,  350], [ 27_000, 22,  280], [ 33_000, 25,  420],
    [ 70_000, 40, 1_800], [ 75_000, 42, 2_000], [ 68_000, 38, 1_700],
    [ 72_000, 41, 1_900], [ 78_000, 44, 2_100], [ 71_000, 39, 1_750],
    [140_000, 35, 6_000], [155_000, 38, 6_500], [148_000, 36, 6_200],
    [160_000, 40, 7_000], [145_000, 34, 5_800], [152_000, 37, 6_400],
    [ 95_000, 55, 3_200], [ 98_000, 58, 3_400], [102_000, 60, 3_500],
    [ 92_000, 53, 3_000], [ 99_000, 57, 3_300], [105_000, 62, 3_600],
])


def segment(X):
    scaler = StandardScaler().fit(X)
    X_scaled = scaler.transform(X)

    inertias = {}
    for k in range(2, 8):
        km = KMeans(n_clusters=k, n_init=10, random_state=42).fit(X_scaled)
        inertias[k] = km.inertia_

    # Find the elbow: first K where relative drop falls below 30%
    chosen_k = 2
    prev = inertias[2]
    for k in range(3, 8):
        rel_drop = (prev - inertias[k]) / prev
        if rel_drop < 0.30:
            chosen_k = k - 1
            break
        prev = inertias[k]
        chosen_k = k

    final_km = KMeans(n_clusters=chosen_k, n_init=10, random_state=42).fit(X_scaled)
    centroids_original = scaler.inverse_transform(final_km.cluster_centers_)

    unique, counts = np.unique(final_km.labels_, return_counts=True)
    counts_by_cluster = dict(zip(unique.tolist(), counts.tolist()))

    return chosen_k, centroids_original, counts_by_cluster


k, centroids, counts = segment(X)
print(f"Chosen K: {k}")
print(f"Centroids (income, age, spend):\n{centroids.round(0).astype(int)}")
print(f"Counts per cluster: {counts}")

# Output (the order of the clusters can differ between sklearn versions):
# Chosen K: 4
# Centroids (income, age, spend):
# [[ 72333     41   1875]     # mid income, mid spend
#  [150000     37   6317]     # high income, high spend
#  [ 29167     23    333]     # young, low income, low spend
#  [ 98500     58   3333]]    # older, comfortable
# Counts per cluster: {0: 6, 1: 6, 2: 6, 3: 6}

The four clusters make sense as customer groups: young/low-spend, mid/professional, older/comfortable, and high-income/high-spend. This is the real test of a clustering: does the marketing team agree with the groups when you show them the centroids?

Production extensions worth knowing:

  • Add descriptive cluster names after looking at the centroids. Don't give a stakeholder cluster_id = 3. Give them "affluent_seniors".
  • For 50+ features, do PCA or UMAP first, then K-Means on the embedding.
  • Replace K-Means with GMM if you want soft assignments ("60% premium, 40% budget").
  • Replace K-Means with HDBSCAN if the densities vary or you want a "doesn't fit anywhere" noise bucket for outlier-flagging.

What You Learned

  • Unsupervised learning has no labels. It finds structure in raw data. Clustering and dimensionality reduction are the two big families.
  • K-Means is the default: fast and simple, but it assumes spherical, equal-sized clusters. Pick K via elbow (on inertia_) and silhouette.
  • DBSCAN finds arbitrary shapes and labels noise. Tune eps and min_samples. HDBSCAN is its modern successor for varying-density data.
  • Hierarchical clustering produces a dendrogram. It is easy to interpret but needs O(n²) memory.
  • Always scale features before clustering. Otherwise the features with the largest values dominate the distance.
  • PCA for linear dimensionality reduction, t-SNE for 2D plots only (never downstream), and UMAP for both visualisation and preprocessing.
  • The standard high-dim pipeline: scale → UMAP/PCA → cluster.
  • Evaluate with silhouette / Davies-Bouldin but trust business validation most. Clusters that don't lead to any action are worthless.
  • GMM gives soft cluster probabilities when boundaries are fuzzy.
  • Cluster IDs are categorical. Never treat them as ordinal numbers.

You have now covered the standard supervised and unsupervised methods. Next come the specialist tracks: time-series forecasting, NLP, and putting models into production. See the Data Science track for what's next.

Found this useful? Share it with a friend who’s learning Python.