Clustering

import numpy as np
from sklearn.cluster import KMeans
from sklearn.metrics import silhouette_score

rng = np.random.default_rng(0)
noise = rng.uniform(0, 1, (600, 2))            # no structure whatsoever
labels = KMeans(n_clusters=4, n_init=10, random_state=0).fit_predict(noise)

np.bincount(labels)                             # four tidy, balanced groups
silhouette_score(noise, labels)                 # about 0.36 — not obviously bad

import numpy as np
from sklearn.cluster import KMeans
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

km = KMeans(n_clusters=3, n_init=10, init="k-means++", random_state=0)
pipe = make_pipeline(StandardScaler(), km).fit(X)

labels = km.labels_                # one integer per row
km.cluster_centers_                 # shape (3, n_features), in scaled space
km.inertia_                         # the objective it minimised
km.n_iter_                          # how many assign/update rounds it took

# Predicting is just 'which centroid is nearest', so it extends to new rows.
pipe.predict(X_new)

import numpy as np
from sklearn.cluster import KMeans

inertias = []
ks = range(1, 11)
for k in ks:
    km = KMeans(n_clusters=k, n_init=10, random_state=0).fit(X_scaled)
    inertias.append(km.inertia_)

# The bend, found numerically: the largest second difference.
drops = -np.diff(inertias)                 # how much each extra cluster bought
bend = int(np.argmax(-np.diff(drops)) + 2)
print(f"inertias {[round(v) for v in inertias]}  elbow at k={bend}")

import numpy as np
from sklearn.cluster import KMeans
from sklearn.metrics import silhouette_score, silhouette_samples

for k in range(2, 7):
    labels = KMeans(n_clusters=k, n_init=10, random_state=0).fit_predict(X_scaled)
    print(k, round(silhouette_score(X_scaled, labels), 3))

# Per-point, to find the cluster that is dragging the average down.
labels = KMeans(n_clusters=4, n_init=10, random_state=0).fit_predict(X_scaled)
sil = silhouette_samples(X_scaled, labels)
for c in range(4):
    print(f"cluster {c}: mean silhouette {sil[labels == c].mean():.3f}")

import numpy as np
from sklearn.cluster import KMeans
from sklearn.datasets import make_blobs, make_moons

# 1. Elongated: k-means cuts across the long axis, not between the bands.
X, _ = make_blobs(n_samples=600, centers=2, cluster_std=0.6, random_state=0)
X_long = X @ np.array([[4.0, 0.0], [0.0, 0.25]])

# 2. Non-convex: two interleaving crescents. No pair of straight boundaries works.
X_moons, y_moons = make_moons(n_samples=600, noise=0.06, random_state=0)
KMeans(n_clusters=2, n_init=10, random_state=0).fit_predict(X_moons)

# 3. Unequal size: the big cluster gets two centroids, the small one shares.
X_uneq, _ = make_blobs(n_samples=[1000, 60], centers=[[0, 0], [6, 0]],
                       cluster_std=[1.6, 0.3], random_state=0)

from scipy.cluster.hierarchy import dendrogram, fcluster, linkage
from sklearn.cluster import AgglomerativeClustering

# scipy: the full tree, then cut it wherever you like.
Z = linkage(X_scaled, method="ward")        # shape (n-1, 4): merge records
labels_3 = fcluster(Z, t=3, criterion="maxclust")
labels_by_dist = fcluster(Z, t=7.5, criterion="distance")   # cut at a height

# scikit-learn: one cut, chosen up front, and it fits in a Pipeline.
agg = AgglomerativeClustering(n_clusters=3, linkage="ward").fit(X_scaled)
agg.labels_

# Or state the merge distance and let it decide how many clusters that implies.
AgglomerativeClustering(n_clusters=None, distance_threshold=7.5,
                        linkage="ward").fit(X_scaled).n_clusters_

import numpy as np
from sklearn.cluster import DBSCAN
from sklearn.neighbors import NearestNeighbors
from sklearn.preprocessing import StandardScaler

Xs = StandardScaler().fit_transform(X)

# Choosing eps: plot the sorted distance to the min_samples-th neighbour and
# take the value at the knee.
d, _ = NearestNeighbors(n_neighbors=5).fit(Xs).kneighbors(Xs)
kdist = np.sort(d[:, -1])
eps = float(np.quantile(kdist, 0.95))

db = DBSCAN(eps=eps, min_samples=5).fit(Xs)
labels = db.labels_
n_clusters = len(set(labels)) - (1 if -1 in labels else 0)
n_noise = int((labels == -1).sum())
print(f"{n_clusters} clusters, {n_noise} noise points")

import numpy as np
from sklearn.cluster import KMeans
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

# income in the tens of thousands, age in years, score in 0-1.
raw = KMeans(n_clusters=3, n_init=10, random_state=0).fit(X_raw)
np.round(raw.cluster_centers_, 2)      # the centres differ only in column 0

scaled = make_pipeline(
    StandardScaler(),
    KMeans(n_clusters=3, n_init=10, random_state=0),
).fit(X_raw)
np.round(scaled[-1].cluster_centers_, 2)   # every column separates the centres

from sklearn.metrics import (adjusted_mutual_info_score, adjusted_rand_score,
                             calinski_harabasz_score, davies_bouldin_score,
                             silhouette_score)

for name, labels in (("kmeans", km_labels), ("dbscan", db_labels)):
    print(name,
          round(silhouette_score(X_scaled, labels), 3),
          round(calinski_harabasz_score(X_scaled, labels), 1),
          round(davies_bouldin_score(X_scaled, labels), 3))

# When a ground truth does exist (a benchmark, or a held-out human labelling):
adjusted_rand_score(y_true, km_labels)          # 0 for a random labelling, 1 for exact
adjusted_mutual_info_score(y_true, km_labels)   # insensitive to how many clusters

from sklearn.cluster import KMeans
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import GridSearchCV
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler

# KMeans as a transformer: transform() returns distances to every centroid.
pipe = Pipeline([
    ("scale", StandardScaler()),
    ("rbf", KMeans(n_init=10, random_state=0)),   # n_clusters columns out
    ("clf", LogisticRegression(max_iter=2000)),
])
search = GridSearchCV(pipe, {"rbf__n_clusters": [5, 20, 50, 100]},
                      cv=5, scoring="roc_auc").fit(X_train, y_train)
search.best_params_