M50: Clustering with k-Means; scaling effects

CSCI 1109 — Practical Data Science

Frank Rudzicz

Learning outcomes

By the end of this module, you should be able to:

  1. Explain in plain language what clustering is and how k-means forms clusters from unlabeled points.
  2. Write the k-means objective function and walk through the assignment–update iteration.
  3. Implement k-means in scikit-learn, both with and without feature scaling, and compare the resulting clusters.
  4. Diagnose when feature scales and outliers are distorting k-means results, using scatterplots, cluster centers, and simple diagnostics.
  5. Reflect on the risks of clustering people (e.g., customers, patients) and articulate simple safeguards and evaluation strategies.

From labels to nothing

So far, our models have had labels:

  • “Does this patient have cancer?” (classification)
  • “How many days until readmision?” (regression)

Now we flip the question:

What if we only have features, no labels, and we want to discover structure?

  • Unsupervised learning:
    • Input: feature vectors \(x_i\).
    • No “right answers” \(y_i\) for training.
    • Goal: find patterns that might be useful for:
      • understanding the data;
      • designing interventions;
      • building downstream models.

Clustering is a core unsupervised task:

  • Group points into “similar” subsets.
  • Here: focus on k-means clustering and how scaling changes the story.

Conceptual scaffold

Key definitions

Clustering.

  • Given points \(x_1, \dots, x_n \in \mathbb{R}^d\), find groups (“clusters”) such that:
    • points within a cluster are similar/close;
    • points in different clusters are less similar/farther.

\(k\)-means clustering.

  • Choose a number of clusters \(k\) (hyperparameter).

  • Algorithm finds:

    • Cluster centers (or centroids) \(\mu_1, \dots, \mu_k\).
    • An assignment of each point \(x_i\) to exactly one cluster \(z_i \in \{1, \dots, k\}\).

  • Objective (informally):

    Put points with nearby neighbours, and place each centroid in the “middle” of its cluster.

Centroid.

  • The mean vector of the points currently in a cluster: \[ \mu_j = \frac{1}{|\mathcal{C}_j|} \sum_{i \in \mathcal{C}_j} x_i. \]

Within-cluster sum of squares (WCSS; ‘distortion’).

  • Measures total squared distance from each point to its cluster centre: \[ \text{WCSS} = \sum_{i=1}^n \left\| x_i - \mu_{z_i} \right\|^2. \]
  • \(k\)-means tries to minimize WCSS.

What k-means is doing

Math lens

The k-means objective

Formally, in \(k\)-means we choose:

  • assignments \(z_1, \dots, z_n \in \{1,\dots,k\}\),
  • centroids \(\mu_1, \dots, \mu_k \in \mathbb{R}^d\),

to minimize the within-cluster sum of squares:

\[ \min_{\mu_1,\dots,\mu_k,\,z_1,\dots,z_n} \sum_{i=1}^n \left\| x_i - \mu_{z_i} \right\|^2 \]

  • This is a hard optimization problem (non-convex, combinatorial).
  • We don’t search all assignments; instead we use a greedy, alternating algorithm.

Lloyd’s algorithm (standard \(k\)-means).

  1. Initialize centroids \(\mu_1,\dots,\mu_k\) (e.g., random points, or k-means++ seeding).
  2. Assignment step.
    • For each point \(x_i\), assign it to the nearest centroid: \[ z_i \leftarrow {\arg\min}_{j \in \{1,\dots,k\}} \, \|x_i - \mu_j\|^2. \]
  3. Update step.
    • For each cluster \(j\), recompute its centroid as the mean: \[ \mu_j \leftarrow \frac{1}{|\mathcal{C}_j|} \sum_{i \in \mathcal{C}_j} x_i. \]
  4. Repeat 2–3 until assignments ‘stop changing’ (or we hit a max number of iterations).

Lloyd’s algorithm (standard \(k\)-means).

Non-trivial but important fact.

  • Each assignment step and each update step never increases WCSS:
    • assignment: we pick the closest centroid;
    • update: mean is the point that minimizes squared distances within that cluster.
  • So WCSS decreases or stays the same each iteration.
  • Because there are only finitely many assignments, the algorithm converges, but:
    • it may end in a local minimum;
    • different initializations can give different solutions.
  • In sklearn’s 🔗KMeans, the attribute n_init controls the number of iterations we run, from which we choose the best result.

Silhouette score

how separated are the clusters?

So far we’ve looked at k-means geometry and some sanity checks.

The silhouette score is an internal metric that tries to answer:

For each point, am I closer to my own cluster than to any other cluster?

For a given point \(i\):

  • Let \(a(i)\) be the average distance from \(i\) to all other points in the same cluster.
  • Let \(b(i)\) be the smallest average distance from \(i\) to any other cluster.

Then we compare \(a(i)\) vs \(b(i)\):

  • If \(a(i)\) is much smaller than \(b(i)\):
    \(\Rightarrow\) \(i\) is tightly inside its own cluster, far from others.
  • If \(a(i) \approx b(i)\):
    \(\Rightarrow\) \(i\) is near a boundary between clusters.
  • If \(a(i) > b(i)\):
    \(\Rightarrow\) \(i\) is (on average) closer to another cluster than its own.

a(i) ≈ 0.55, b(i) ≈ 3.81

Silhouette score: formula & interpretation

For each point \(i\), the silhouette value is:

\[ \boxed{s(i) = \frac{b(i) - a(i)}{\max\{a(i), b(i)\}}} \]

where:

  • \(a(i)\) = average distance from \(i\) to its own cluster

  • \(b(i)\) = smallest average distance from \(i\) to points in any other cluster

  • \(s(i)\) is always between \(-1\) and \(+1\).

We often look at the average silhouette score over all points:

\[ \bar{s} = \frac{1}{n} \sum_{i=1}^n s(i) \]

  • Rough rule-of-thumb:
    • \(\bar{s} < 0.1\): very weak / overlapping clusters
    • \(\bar{s} \approx 0.2\)–\(0.5\): reasonable structure
    • \(\bar{s} > 0.5\): clearly separated clusters (for simple cases)

In code (via 🔗sklearn.metrics.silhouette_score), we pass in:

  • the feature matrix (usually scaled), and
  • the cluster labels.

❌ Anti-example: wrong tool for the job

Scenario. You have data for hospital patients:

  • age (0–100),
  • systolic_bp (80–220),
  • home_postal_code (string),
  • diagnosis_code (categorical),
  • lab_risk_score (continuous).

You feed all columns directly into k-means using Euclidean distance on raw columns.

What goes wrong?

  • k-means assumes a numeric vector space; categorical strings like diagnosis_code and postal_code are not meaningful in this geometry.
  • Huge-scale features (e.g., lab score vs age) can dominate distances.
  • Clusters might be driven by weird numeric encodings (e.g., postal codes mapped to arbitrary integers).

Better approach?

  • Restrict to meaningful numeric features (or encode categories properly).
  • Scale features (standardization) so units don’t dominate.
  • Check whether “roughly spherical” clusters in feature space are plausible for your use case.

This anti-example highlights:

  • k-means is not a magical “find groups” button;
  • You must align the geometry (features & scaling) with the domain question.

Worked example 1:
Toy customer segments

Setup and first look

We’ll simulate a small “customer” dataset:

  • annual_spend (in $),
  • visits_per_month (in counts),

and try to discover 3 customer segments using k-means.

Code 
import numpy as np
import pandas as pd
from sklearn.cluster import KMeans
from sklearn.preprocessing import StandardScaler
from sklearn.pipeline import Pipeline
import matplotlib.pyplot as plt

rng = np.random.default_rng(1109)

# Simulate three rough groups
n_per_cluster = 80

group1 = np.column_stack([
    rng.normal(500, 50, n_per_cluster),    # low spend
    rng.normal(2, 0.5, n_per_cluster)      # few visits
])

group2 = np.column_stack([
    rng.normal(1500, 100, n_per_cluster),  # medium spend
    rng.normal(4, 0.7, n_per_cluster)      # moderate visits
])

group3 = np.column_stack([
    rng.normal(3000, 200, n_per_cluster),  # high spend
    rng.normal(6, 1.0, n_per_cluster)      # many visits
])

data = np.vstack([group1, group2, group3])
customers = pd.DataFrame(data, columns=["annual_spend", "visits_per_month"])
customers.head()
annual_spend visits_per_month
0 491.491400 1.839378
1 524.150273 2.645292
2 510.448798 1.900917
3 498.516362 1.163426
4 409.563898 1.749586

\(k\)-means without scaling

Code 
kmeans_raw = KMeans(n_clusters=3, n_init=10, random_state=1109)
customers["cluster_raw"] = kmeans_raw.fit_predict(customers[["annual_spend", "visits_per_month"]])

customers.groupby("cluster_raw")[["annual_spend", "visits_per_month"]].mean().round(1)
annual_spend visits_per_month
cluster_raw
0 1495.0 4.0
1 2996.6 5.8
2 500.4 1.9

Now plot:

Code 
fig, ax = plt.subplots()
scatter = ax.scatter(
    customers["annual_spend"],
    customers["visits_per_month"],
    c=customers["cluster_raw"]
)
centers = kmeans_raw.cluster_centers_
ax.scatter(centers[:, 0], centers[:, 1], marker="X", s=200, edgecolor="black")
ax.set_xlabel("Annual spend ($)")
ax.set_ylabel("Visits per month")
ax.set_title("k-means on raw features")
plt.show()

Note

This looks okay because our features were designed nicely.

k-means with distorted scaling

Suppose we accidentally record:

  • annual_spend in thousands of dollars;
  • visits_per_month unchanged.

Now annual_spend has much smaller numerical range.

Code 
customers_bad = customers.copy()
customers_bad["annual_spend_k"] = customers_bad["annual_spend"] / 1000

kmeans_bad = KMeans(n_clusters=3, n_init=10, random_state=1109)
customers_bad["cluster_bad"] = kmeans_bad.fit_predict(
    customers_bad[["annual_spend_k", "visits_per_month"]]
)

customers_bad.groupby("cluster_bad")[["annual_spend_k", "visits_per_month"]].mean().round(2)
annual_spend_k visits_per_month
cluster_bad
0 2.95 6.00
1 0.54 1.97
2 1.64 4.05
Code 
fig, ax = plt.subplots()

scatter = ax.scatter(
    customers_bad["annual_spend_k"],
    customers_bad["visits_per_month"],
    c=customers_bad["cluster_bad"]
)
centers_bad = kmeans_bad.cluster_centers_
ax.scatter(centers_bad[:, 0], centers_bad[:, 1], marker="X", s=200, edgecolor="black")
ax.set_xlabel("Annual spend (thousands of $)")
ax.set_xlim(0,9)
ax.set_ylabel("Visits per month")
ax.set_title("k-means with distorted scaling")
plt.show()

What often happens conceptually

  • If one feature’s range is tiny relative to another, distances are dominated by the larger-scale feature.
  • Here, distances are mostly about visits_per_month; annual_spend_k contributes much less.
  • Clusters may become horizontal bands (by visits), not meaningful spend–visit groups.

This is why we need feature scaling.

Fixing the geometry with StandardScaler

Code 
pipe_scaled = Pipeline([
    ("scaler", StandardScaler()),
    ("kmeans", KMeans(n_clusters=3, n_init=10, random_state=1109))
])

customers["cluster_scaled"] = pipe_scaled.fit_predict(
    customers[["annual_spend", "visits_per_month"]]
)

customers.groupby("cluster_scaled")[["annual_spend", "visits_per_month"]].mean().round(1)
annual_spend visits_per_month
cluster_scaled
0 1495.0 4.0
1 2996.6 5.8
2 500.4 1.9
Code 
fig, ax = plt.subplots()
scatter = ax.scatter(
    customers["annual_spend"],
    customers["visits_per_month"],
    c=customers["cluster_scaled"]
)
centers_scaled = pipe_scaled.named_steps["kmeans"].cluster_centers_
# Need to inverse-transform centers back to original scale for plotting:
centers_scaled_orig = pipe_scaled.named_steps["scaler"].inverse_transform(centers_scaled)
ax.scatter(centers_scaled_orig[:, 0], centers_scaled_orig[:, 1],
           marker="X", s=200, edgecolor="black")
ax.set_xlabel("Annual spend ($)")
ax.set_ylabel("Visits per month")
ax.set_title("k-means with standardization")
plt.show()

Takeaways

  • \(k\)-means is sensitive to feature scales because it uses Euclidean distance.
  • Scaling (e.g., standardization) can:
    • restore a balanced geometry;
    • produce more interpretable, stable clusters.
  • In practice: use pipelines with StandardScaler for k-means on mixed-scale numeric features.

Worked example 2:
Das ist eine grosse pinguin

Setup: a small “discover the species” task

We’ll use the classic palmer penguins dataset (via seaborn):

  • Numeric features (e.g., flipper_length_mm, body_mass_g).
  • True species labels (Adelie, Chinstrap, 🔗Gentoo) — but we’ll pretend we don’t see them at clustering time.
Code 
import seaborn as sns
from sklearn.metrics import silhouette_score

penguins = sns.load_dataset("penguins").dropna()

features = ["flipper_length_mm", "body_mass_g"]
X = penguins[features].to_numpy()
penguins[features].head()
flipper_length_mm body_mass_g
0 181.0 3750.0
1 186.0 3800.0
2 195.0 3250.0
4 193.0 3450.0
5 190.0 3650.0

We’ll cluster into \(k = 3\) groups and see how they line up with species.

\(k\)-means on penguins (with scaling)

Code 
pipe_penguins = Pipeline([
    ("scaler", StandardScaler()),
    ("kmeans", KMeans(n_clusters=3, n_init=20, random_state=1109))
])

penguins["cluster"] = pipe_penguins.fit_predict(penguins[features])

penguins[["species", "cluster"]].head()
species cluster
0 Adelie 0
1 Adelie 0
2 Adelie 0
4 Adelie 0
5 Adelie 0
Code 
ct = pd.crosstab(penguins["cluster"], penguins["species"])
ct
species Adelie Chinstrap Gentoo
cluster
0 116 45 0
1 30 23 34
2 0 0 85

✏️ Think

  • Did clustering rediscover the labels well?
  • Do you remember that this is not a form of supervised learning?

\(k\)-means on penguins (with scaling)

Code 
fig, ax = plt.subplots()
scatter = ax.scatter(
    penguins["flipper_length_mm"],
    penguins["body_mass_g"],
    c=penguins["cluster"]
)
ax.set_xlabel("Flipper length (mm)")
ax.set_ylabel("Body mass (g)")
ax.set_title("Penguins clustered by k-means (k=3)")
plt.show()

\(k\)-means on penguins (with scaling)

We can also compute a silhouette score as a crude internal metric.

Code 
# Silhouette score: ranges roughly from -1 (bad) to 1 (good)
X_scaled_penguins = pipe_penguins.named_steps["scaler"].transform(penguins[features])
labels_penguins = penguins["cluster"].to_numpy()
sil_penguins = silhouette_score(X_scaled_penguins, labels_penguins)
sil_penguins
0.4760097141728845

Interpretation.

  • Silhouette score near 0.3–0.5: clusters are somewhat separated but not perfect.
  • Scores near 0 or negative would suggest overlapping or messy clusters.

Tip

Silhouette is an internal metric; whether these clusters are useful depends on the domain question.

Summary

Evaluating \(k\)-means: what “good” looks like

Unlike supervised learning, we often don’t have labels. Evaluation strategies:

  1. Internal metrics.
    • Inertia / WCSS (lower is better, but always decreases with more clusters).
    • Silhouette score (higher is better; balances cohesion & separation).
  2. Domain sanity checks.
    • Are cluster centers in realistic ranges?
    • Do clusters correspond to interpretable groups (e.g., “small-light”, “large-heavy” penguins)?
  3. Stability.
    • Run k-means with different seeds or subsamples; do we get similar clusters?
  4. Downstream usefulness.
    • Use clusters as features in a later supervised task.
    • Check whether they improve performance or understanding.

Choosing \(k\)

  • Common heuristic: elbow method.
    • Plot WCSS vs. \(k\); look for a point where improvement slows.
  • But: no universally “correct” \(k\); it’s part technical, part judgement.
Code 
import numpy as np
import matplotlib.pyplot as plt
from sklearn.cluster import KMeans

# -----------------------------
# Simulate a small 2D dataset
# -----------------------------
rng = np.random.default_rng(1109)

n_per_cluster = 80
centers = np.array([
    [0.0, 0.0],
    [4.0, 0.5],
    [0.5, 4.0]
])

points = []
for cx, cy in centers:
    pts = rng.normal(loc=[cx, cy], scale=[0.6, 0.6], size=(n_per_cluster, 2))
    points.append(pts)

X = np.vstack(points)

# -----------------------------
# Compute WCSS for k = 1..8
# -----------------------------
Ks = range(1, 9)
wcss = []

for k in Ks:
    model = KMeans(n_clusters=k, n_init=20, random_state=1109)
    model.fit(X)
    wcss.append(model.inertia_)  # inertia_ is the WCSS

# -----------------------------
# Plot elbow curve
# -----------------------------
fig, ax = plt.subplots(figsize=(6, 4))

ax.plot(list(Ks), wcss, marker="o")
ax.set_xlabel("k (number of clusters)")
ax.set_ylabel("Within-cluster sum of squares (WCSS)")
ax.set_title("Elbow method for choosing k")
ax.grid(alpha=0.3)

# Optional: annotate a visible elbow around k=3
elbow_k = 3
elbow_wcss = wcss[elbow_k - 1]
ax.scatter([elbow_k], [elbow_wcss], s=80, edgecolor="black")
ax.annotate(
    "possible elbow",
    xy=(elbow_k, elbow_wcss),
    xytext=(elbow_k + 0.4, elbow_wcss * 1.1),
    arrowprops=dict(arrowstyle="->", lw=1),
    fontsize=9,
)

plt.tight_layout()
plt.show()

Scaling: practical checklist for \(k\)-means

Before running \(k\)-means on a dataset:

  1. Feature types.
    • Are features numeric and roughly continuous?
    • For categories, consider:
      • one-hot encoding;
      • or a different clustering method.
  2. Feature scales.
    • Are units very different (e.g., dollars vs counts vs probabilities)?
    • Use StandardScaler or MinMaxScaler in a pipeline.
  3. Outliers.
    • Extreme points can pull centroids away
  4. Number of clusters \(k\).
    • Use prior knowledge if possible (“we expect ~4 customer segments”).
    • Explore a small range of k (e.g., 2–6) with elbow/silhouette and domain input.
  5. Initialization.
    • Use k-means++ (scikit-learn’s default) and multiple runs (n_init).

Ethics & risks: clustering people

Clustering often gets applied to people:

  • “high-value customers”, “at-risk students”, “high-risk patients”.

Risks:

  • Essentialism & stereotyping.
    • Clusters can be over-interpreted as “natural kinds”.
  • Unequal treatment.
    • Interventions (discounts, extra support, more scrutiny) may be targeted differently across clusters.
  • Data & feature bias.
    • If features encode existing inequities (income, postcode, etc.), clusters can bake in those patterns.

Lightweight safeguards:

  • Always ask: “What is this clustering used for?”
    • Exploration? Targeted help? Pricing? Policing?
  • Inspect cluster compositions:
    • Are some demographic groups concentrated in “undesirable” clusters?
  • Prefer opt-in, supportive uses (e.g., extra tutoring outreach) over punitive or exclusionary ones.
  • Document:
    • feature choices,
    • scaling decisions,
    • and known limitations of the clusters.

Curiosity: Where else doth \(k\)-means ply its trade?

Beyond toy customer segments and penguins, k-means underlies:

  • Image compression / colour quantization.
    • Cluster pixel colours into a small palette; store palette + assignments.
  • Vector quantization in signal processing.
    • Represent many high-dimensional vectors via a smaller “codebook” of centroids.
  • Initialization for other algorithms.
    • E.g., seeding Gaussian mixture models.
  • As a building block in:
    • text clustering,
    • recommendation systems,
    • anomaly detection (points far from all centroids).
    • modern speech recognition (Hsu et al. 2021; Chen et al. 2022)

In later courses, you may see:

  • density-based clustering (DBSCAN) (Ester et al. 1996),
  • spectral clustering,
  • hierarchical clustering,

which address some limitations of \(k\)-means (non-spherical clusters, varying density, etc.).

Quick-checks

  1. Conceptual MCQ.
    You run k-means on a dataset with two features: age (0–100) and income (0–200,000). You do no scaling. Which is most likely?

    A. Clusters mostly reflect differences in age.
    B. Clusters mostly reflect differences in income.
    C. Both features influence clusters equally.
    D. k-means fails to converge.

  2. Short answer.
    In one sentence, explain what the k-means algorithm is optimizing.

  3. Calculation-ish.
    Suppose a cluster currently contains the 1D points \(\{2, 4, 6, 8\}\). What is the centroid in 1D, and why is this the best choice for k-means?

  4. Diagnosis MCQ.
    You run k-means with \(k=4\) on standardized features and get a silhouette score near 0.02. Scatterplots show heavily overlapping clusters. What is the best immediate interpretation?

    A. There is no structure; the data are purely random noise.
    B. k=4 is perfect; low silhouette is expected after scaling.
    C. The clusters are weak/overlapping; you might need a different k or a different algorithm.
    D. The model overfit; you must regularize more.

References

Chen, Sanyuan, Chengyi Wang, Zhengyang Chen, Yu Wu, Shujie Liu, Zhuo Chen, Jinyu Li, et al. 2022. “WavLM: Large-Scale Self-Supervised Pre-Training for Full Stack Speech Processing.” IEEE Journal of Selected Topics in Signal Processing 16 (6): 1505–18. https://doi.org/10.1109/JSTSP.2022.3188113.
Ester, Martin, Hans-Peter Kriegel, Jörg Sander, and Xiaowei Xu. 1996. “A Density-Based Algorithm for Discovering Clusters in Large Spatial Databases with Noise.” In Proceedings of the Second International Conference on Knowledge Discovery and Data Mining (KDD-96), 226–31. AAAI Press.
Hsu, Wei-Ning, Benjamin Bolte, Yao-Hung Hubert Tsai, Kushal Lakhotia, Ruslan Salakhutdinov, and Abdelrahman Mohamed. 2021. “HuBERT: Self-Supervised Speech Representation Learning by Masked Prediction of Hidden Units.” IEEE/ACM Transactions on Audio, Speech, and Language Processing 29: 3451–60. https://doi.org/10.1109/TASLP.2021.3122291.