k-means Clustering

Context

k-means splits the rows of a table into k groups, a number you choose, with no column to predict: every listing joins the nearest of k centers, every center moves to the mean of the listings that joined it, and the two steps repeat until no listing changes group.

A listings team uses it to describe thousands of listings by a few typical ones, to show a seller the listings most like their own, to place each new listing in a group the moment it is posted, and to pull out the listings far from every group for a check by hand before any model is fitted on them.

Its failures hand back groups without a word. On raw columns the price, whose spread runs to hundreds of thousands, decides the groups alone, one start at the library's defaults stops at a grouping ten starts improve on, the number k-means reports about its fit falls with every group added, and a group of mostly houses gets read as the houses. Lesson 3.2.1 named k-means among the methods whose distances let scale decide what counts, Lesson 3.2.3 built the pipeline that fills and scales, and this lesson fits the module's first method with no target at all.

Goal

• Implement k-means' assign and move steps with NumPy, run them until no listing changes cluster, and check the result against KMeans started from the same listings.
• Distinguish a run stopped by KMeans' tolerance from one resting at a local minimum, and explain precisely what n_init=10 does about both.
• Build KMeans into one pipeline behind SimpleImputer and StandardScaler, place new listings that arrive with gaps, and diagnose the unscaled run that splits the listings by price alone.
• Compute inertia by hand, trace the one-cluster inertia of 18,000 to the standardized columns, and explain why inertia falls at every k and its smallest value cannot choose k.
• Trace centers back to listing units with inverse_transform, check each against its listings' mean, and compare clusters with property_type through pd.crosstab without naming a cluster after a type.

Prerequisites

None of this is required, and the lesson stands on its own. Both are here if you want to watch k-means run on a handful of points, or read scikit-learn's own account of it, before you start.

• YouTube - "StatQuest: K-means clustering" by StatQuest with Josh Starmer. Points on a line and on a plane grouped by picking starting points at random, assigning every point to the nearest and averaging, then rerunning from other starts and keeping the grouping with the least variation.
• Official docs - scikit-learn 1.9 User Guide, "Clustering", the section "K-means". The two steps, inertia as the within-cluster sum of squares, why a run can end in a local minimum, and the k-means++ start.

Table of Contents

• Learning Without a Target
• How k-means Assigns and Moves
• KMeans in a Pipeline
• Inertia, and Why Its Smallest Value Cannot Choose k
• Reading the Centers

Theory

Nothing to install: KMeans, from sklearn.cluster, arrived with Lesson 3.1.1's install. Every block below adds to one script, clusters.py, saved in ~/dev/bootcamp beside listings.csv and run with python3 clusters.py, each block under the one before it.

Learning Without a Target

What. Unsupervised learning (fitting with no target column, to find structure in the features alone) hands fit the features and nothing else. Clustering, one kind of it, splits the rows into clusters (groups of rows that sit near one another), and k-means, the method Lesson 3.2.1 named, describes each cluster by a center (a point with one value per column that stands for the cluster), so that every listing belongs to the nearest of k centers and every center is the mean of its listings. k, the number of clusters, is a hyperparameter you set before the fit. Why. Every model since Lesson 3.1.1 predicted a column whose true values were known, so rows held out from its fit could score it against them. Some questions a listings site asks have no such column, what kinds of listings it holds, how many of each and which look like no others, so no answer marks a grouping right or wrong and no held-out score against an answer keeps one honest. When. When the question is about the listings themselves rather than one of their columns: describing the stock, finding the listings most like a new one, picking out the ones that sit far from all the rest.

Every block works on all 3,000 listings and the six numeric columns, the price among them, because here the price is one more fact about a listing rather than an answer, and with no answer to check them against, no rows are held out. Lesson 3.2.3's imputer and Lesson 3.2.1's scaler fill and standardize them in one pipeline, and the third subtitle shows what each step is for:

import numpy as np
import pandas as pd
from sklearn.cluster import KMeans
from sklearn.impute import SimpleImputer
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

listings = pd.read_csv("listings.csv")
cols = ["area_sqm", "rooms", "bathrooms", "year_built", "distance_km", "price"]
prep = make_pipeline(SimpleImputer(strategy="median"), StandardScaler())
Z = prep.fit_transform(listings[cols])
print(Z.shape)                     # (3000, 6)
print(np.round(Z.std(axis=0), 3))  # [1. 1. 1. 1. 1. 1.]

Six columns, each with a standard deviation of 1. KMeans takes n_clusters, the k, and random_state, because it starts from a random choice, the next subtitle's subject:

km = KMeans(n_clusters=2, random_state=42)
print(km.fit(Z) is km)            # True
print(km.labels_[:8])             # [0 1 0 1 1 1 1 0]
print(km.labels_.shape)           # (3000,)
print(km.cluster_centers_.shape)  # (2, 6)

fit received Z alone and returned the estimator, as every fit since Lesson 3.1.1 has. labels_ holds one cluster number per listing, and cluster_centers_ one row per center with one value per column, in the standardized units of Z. Listing 0, the 175.7 square meter listing that opens the file, against both centers with Lesson 2.3.3's distance:

first = Z[0]              # listing 0, standardized
dist = np.linalg.norm(km.cluster_centers_ - first, axis=1)
print(np.round(dist, 3))  # [1.461 4.289]
print(km.labels_[0])      # 0

It sits 1.461 from center 0 and 4.289 from center 1, so it belongs to cluster 0: a cluster is the set of listings that share a nearest center, and its number only says which center that is. A target handed to KMeans changes nothing:

with_y = KMeans(n_clusters=2, random_state=42)
with_y.fit(Z, listings["price"])  # a target, passed anyway
print(np.array_equal(with_y.labels_, km.labels_))

True

KMeans.fit accepts a y only so that it can be called like every estimator's fit, and never reads it: the price counts only as one of the six columns of Z. Leaving out n_clusters does not let the data choose k either:

default = KMeans(random_state=42).fit(Z)
print(np.bincount(default.labels_))  # [289 461 425 216 596 317 423 273]

np.bincount counts the listings that carry each label, label 0 first: eight clusters, because n_clusters defaults to 8, for a table nothing says holds eight groups.

DON'T:

• DON'T pass the price as y and expect clusters built around it. fit(Z, listings["price"]) returned every label fit(Z) did, because KMeans.fit never reads its y.
• DON'T leave n_clusters out. KMeans(random_state=42) split the listings into its default of 8 clusters, whatever the table holds.

How k-means Assigns and Moves

What. k-means reaches its centers by repeating two steps from k starting points: assign every listing to its nearest center by Lesson 2.3.3's distance, then move every center to the mean of the listings assigned to it, until an assignment changes no listing. What the steps lower is inertia (the sum of every listing's squared distance to its own center), and where they come to rest is a local minimum (a grouping neither step can improve, though another grouping may be lower). Why. Every grouping KMeans returns comes out of this loop, and two choices made around it, where it starts and when it stops, decide which grouping that is, so a result cannot be read without knowing both. When. Whenever a KMeans result is read, compared with another or reported, and before a single run is trusted.

The loop by hand, named lloyd after the algorithm setting under which KMeans runs it, started from the first three listings:

def lloyd(Z, centers):
    k = len(centers)
    labels = np.full(len(Z), -1)           # a label no center has
    rounds = 0
    while True:
        rounds += 1
        # assign: each listing to its nearest center
        gaps = Z[:, np.newaxis, :] - centers[np.newaxis, :, :]
        nearest = np.linalg.norm(gaps, axis=2).argmin(axis=1)
        if np.array_equal(nearest, labels):
            return labels, centers, rounds
        labels = nearest
        # move: each center to the mean of its listings
        centers = np.array([Z[labels == j].mean(axis=0) for j in range(k)])


labels, centers, rounds = lloyd(Z, Z[:3])  # the first three listings
# each listing's own center beside it, then the squared distances
inertia = ((Z - centers[labels]) ** 2).sum()
print(rounds, round(inertia, 1))           # 45 9050.1
print(np.bincount(labels))                 # [ 861  943 1196]

The broadcast is Lesson 2.3.3's table of distances from every listing to every center, (3000, 3) after the norm, and centers[labels], Lesson 2.2.2's fancy indexing, sets each listing's own center beside it. The loop ran 45 rounds, the last an assignment that changed no listing, and came to rest at an inertia of 9,050.1.

Neither step can raise the inertia. Assigning gives every listing its nearest center, so no squared distance grows. Moving puts each center at the mean of its listings, and Lesson 3.1.3 showed the mean is the one value with the smallest sum of squared distances to a set of numbers, column by column. Each round lowers the inertia or leaves it, there are only so many ways to split 3,000 listings three ways, and once an assignment changes nothing, the move returns the same centers and the loop has nowhere left to go.

KMeans takes a start of your own as init, an array with one row per center, and n_init=1 runs it once:

same = KMeans(n_clusters=3, init=Z[:3], n_init=1).fit(Z)
print(round(same.inertia_, 1))         # 9050.1
print(same.n_iter_)                    # 41
print((same.labels_ == labels).sum())  # 2998
own = [Z[same.labels_ == j].mean(axis=0) for j in range(3)]
# each center against the mean of its own cluster's listings
print(np.allclose(own, same.cluster_centers_))

False

The same inertia to one decimal, after 41 rounds where the loop took 45, with 2 listings in another cluster and centers not all at their clusters' means. KMeans does not wait for an assignment that changes nothing. After each round it adds up the squared distances its centers moved, and stops once that sum is no larger than tol times the columns' average variance, 1e-4 on Z at the default tol of 1e-4. It then assigns every listing to the nearest of those centers and moves no center again, so two listings the loop went on to move stay put. A tolerance a hundred times smaller:

to_end = KMeans(n_clusters=3, init=Z[:3], n_init=1, tol=1e-6).fit(Z)
print(to_end.n_iter_)                    # 45
print((to_end.labels_ == labels).sum())  # 3000
own = [Z[to_end.labels_ == j].mean(axis=0) for j in range(3)]
print(np.allclose(own, to_end.cluster_centers_))

True

At tol=1e-6 the centers never moved that little before an assignment changed nothing, so KMeans ran the loop's 45 rounds, agreed with it on all 3,000 listings, and left every center at its cluster's mean.

Where the loop comes to rest depends on where it starts. Four centers, started once from the first four listings and once from listings 32 to 35:

for start in [0, 32]:
    four = Z[start:start + 4]  # four listings in a row
    found, where, steps = lloyd(Z, four)
    sum_sq = ((Z - where[found]) ** 2).sum()
    print(start, steps, round(sum_sq, 1))

0 18 7742.8
32 18 8017.5

Both runs ended the way the loop ends, on an assignment that changed no listing, so each is a local minimum: neither step improves either grouping, and the start from listings 32 to 35 settled at 8,017.5, against 7,742.8 from listings 0 to 3.

At its defaults, KMeans picks its own start with k-means++ (its starting rule: a first center drawn at random from the listings, then each further center picked from a few listings drawn with a chance that grows with their squared distance to the nearest center already chosen), so random_state fixes the draw, and n_init sets how many starts it runs:

print(KMeans().get_params()["n_init"])  # auto
single = KMeans(n_clusters=3, random_state=42).fit(Z)
print(round(single.inertia_, 1))        # 9439.0
print(single.n_iter_)                   # 12
once = KMeans(n_clusters=3, random_state=42, n_init=1).fit(Z)
print(np.array_equal(once.labels_, single.labels_))

True

"auto" means one start for k-means++, and n_init=1 gave the same labels. That start stopped at 9,439.0 after 12 rounds, far above the loop's 9,050.1. The same start at the smaller tolerance says whether 9,439.0 is a local minimum, here beside n_init=10, which runs ten starts and keeps the lowest inertia:

slow = KMeans(n_clusters=3, random_state=42, tol=1e-6).fit(Z)
print(round(slow.inertia_, 1), slow.n_iter_)
ten = KMeans(n_clusters=3, random_state=42, n_init=10).fit(Z)
print(round(ten.inertia_, 1))

9050.1 47
9050.1

It is not one. With tol=1e-6 the same start ran 47 rounds to 9,050.1, so the default run had stopped in a round that happened to move its centers by less than the tolerance, while they still had far to go. Ten starts kept 9,050.1 too: a start that stops early or settles in a poor local minimum loses to any other start that ends lower. The rest of the Theory fits KMeans with n_init=10.

DON'T:

• DON'T trust the default single run. Seed 42's one start stopped at 9,439.0 at k=3, where ten starts keep 9,050.1.
• DON'T expect KMeans to stop where the loop stops. At its default tolerance it stopped the first three listings' run after 41 rounds, with 2 listings in other clusters and its centers not all at their clusters' means, and seed 42's after 12 rounds, which tol=1e-6 ran to 47.
• DON'T take the loop's resting point as the lowest grouping. Started from listings 32 to 35 at k=4 it settled at 8,017.5, a local minimum, where listings 0 to 3 led to 7,742.8.

KMeans in a Pipeline

What. KMeans refuses a gap and measures distance in each column's own units, so it goes last in a pipeline behind the imputer and the scaler, built with Lesson 3.2.3's make_pipeline. The pipeline's fit and its fit_predict (fit, then return the labels, in one call) take the raw columns and no y, and its predict places a new listing in a cluster, its gaps filled and its values scaled with the numbers the fit learned. Why. The fill values and the scaler's means and spreads decide where every listing sits before any distance is measured, so they belong to the clustering as much as the centers do, and a listing placed later has to pass through exactly the same numbers. When. Whenever k-means meets a table with gaps or with columns in different units, which on listings.csv is every time.

The raw columns, handed to KMeans directly:

try:
    KMeans(n_clusters=3, random_state=42, n_init=10).fit(listings[cols])
except ValueError as err:
    print(str(err).splitlines()[0])  # Input X contains NaN.

The message goes on to say that KMeans does not accept missing values encoded as NaN, and suggests an imputer in a pipeline or dropping the rows. The pipeline:

pipe = make_pipeline(
    SimpleImputer(strategy="median"),
    StandardScaler(),
    KMeans(n_clusters=3, random_state=42, n_init=10),
)
scaled = pipe.fit_predict(listings[cols])  # raw columns, no y
print([name for name, step in pipe.steps])
print(round(pipe["kmeans"].inertia_, 1))   # 9050.1
print(np.bincount(scaled))                 # [ 942 1197  861]

['simpleimputer', 'standardscaler', 'kmeans']

make_pipeline named the steps after their classes, as in Lesson 3.2.3, and the pipeline reached the 9,050.1 that ten starts reached on Z. It places listings it has never seen, gaps and all:

new = pd.DataFrame({
    "area_sqm": [np.nan, 160.0],     # the first seller left it blank
    "rooms": [2, 5],
    "bathrooms": [1, 3],
    "year_built": [1975.0, np.nan],  # the second seller left it blank
    "distance_km": [2.5, 9.0],
    "price": [310000, 690000],
})
print(pipe.predict(new))             # [0 2]

The missing area and year were filled with the imputer's medians from the 3,000 listings, both listings were scaled with the scaler's means and spreads, and each went to its nearest center.

The same three clusters without the scaler, on the filled columns in their own units:

unscaled = make_pipeline(
    SimpleImputer(strategy="median"),
    KMeans(n_clusters=3, random_state=42, n_init=10),
)
raw = unscaled.fit_predict(listings[cols])
print(listings.groupby(raw)["price"].agg(["min", "max", "size"]))
print(listings.groupby(raw)["area_sqm"].agg(["min", "max"]))

      min      max  size
0  344000   613000  1095
1   60000   343000  1544
2  614000  1930000   361
    min    max
0  34.7  237.4
1  24.0  171.2
2  51.4  278.5

Three price bands that do not overlap, 60,000 to 343,000, 344,000 to 613,000 and 614,000 to 1,930,000, while each band's areas run from small to large. The spreads the scaler learned, the standard deviations of the filled columns, say why:

spreads = pd.Series(prep["standardscaler"].scale_, index=cols)
print(spreads.round(2))
# each column's share of the summed variances
share = spreads ** 2 / (spreads ** 2).sum()
print(share.round(6))

area_sqm           43.18
rooms               1.71
bathrooms           0.81
year_built         22.12
distance_km         3.04
price          195267.73
dtype: float64
area_sqm       0.0
rooms          0.0
bathrooms      0.0
year_built     0.0
distance_km    0.0
price          1.0
dtype: float64

A squared distance adds up two listings' squared gaps column by column, and averaged over every pair of listings, each column adds in proportion to its variance, the square of its spread. The price's spread is 195,267.73 against 43.18 square meters and 22.12 years, so its share of the summed variances rounds to 1.0 at six decimals, and three clusters on these columns come out as three ranges of price. The scaled run is the same fit with that cause taken away, every spread set to 1:

print(listings.groupby(scaled)["price"].agg(["min", "max"]))

      min      max
0   60000  1233000
1   60000   850000
2  228000  1930000

Each scaled cluster's prices overlap the others', cluster 0's running from 60,000 to 1,233,000, because all six columns now count. Neither run is the listings' true grouping: standardizing decides that every column counts alike, which suits a question about the listings, and a team that wants price bands can draw them with Lesson 2.6.1's pd.cut.

DON'T:

• DON'T hand KMeans a table with gaps. The raw columns stopped it with Input X contains NaN., and the imputer in front of it fills the 125 missing areas and the 140 missing years.
• DON'T cluster columns measured in different units without scaling them. The price's spread, 195,267.73 against 43.18 square meters, made three clusters three price bands, 60,000 to 343,000, 344,000 to 613,000 and 614,000 to 1,930,000.

Inertia, and Why Its Smallest Value Cannot Choose k

What. inertia_ is the sum the loop lowers, every listing's squared distance to its own center added up: Lesson 3.1.3's squared error with each listing's center standing in for its prediction, which makes inertia the loss k-means minimizes. With one center, the center is the mean listing, and the inertia is the whole squared spread of the table, which on standardized columns is the number of listings times the number of columns. Why. Inertia is the one measure of fit KMeans reports, and what it measures, how near the listings sit to their centers, improves with every center added, whether or not the table holds that many groups. When. To compare starts at one k on one table, the comparison n_init makes with it, and never to choose k by its smallest value or to compare two tables by their raw inertias.

one = KMeans(n_clusters=1, random_state=42, n_init=10).fit(Z)
print(round(one.inertia_, 3))                # 18000.0
print(np.allclose(one.cluster_centers_, 0))  # True
print(np.round((Z ** 2).sum(axis=0), 3))
by_hand = ((Z - ten.cluster_centers_[ten.labels_]) ** 2).sum()
print(round(by_hand, 1), round(ten.inertia_, 1))

[3000. 3000. 3000. 3000. 3000. 3000.]
9050.1 9050.1

With one center, the loop puts it at the mean listing, 0 in every column of Z, so the inertia is the sum of every squared standardized value: each column's 3,000 squares add up to 3,000 times its variance, which standardizing set to 1, and six columns make 18,000. The last line sums the ten-start fit's squared distances by hand, and inertia_ holds the same 9,050.1. One fit for every k from 1 to 8, each with ten starts:

for k in range(1, 9):
    fitted = KMeans(n_clusters=k, random_state=42, n_init=10).fit(Z)
    print(k, round(fitted.inertia_, 1))

1 18000.0
2 11057.0
3 9050.1
4 7742.9
5 6831.9
6 6179.9
7 5799.2
8 5470.8

Inertia fell at every step, and the lowest inertia k + 1 centers can reach is never above the lowest k centers can: add a center on any one listing of the best k-grouping, and no listing ends farther from its nearest center while that one listing's distance becomes 0. Taken to the end, 3,000 centers, one on each listing, make every distance 0, Lesson 3.1.2's memorizer again: a perfect score on the rows the fit was made from, earned by remembering them.

So a smaller inertia never says that a k is right, as a training score never says that a model is, Lesson 3.6.1's point. With no target to check held-out rows against, k is chosen by whether the clusters hold when the start changes and how cleanly they separate, Lesson 3.8.3's two questions, and by what the clusters are for.

DON'T:

• DON'T choose k by the smallest inertia. It fell at every k, from 18,000.0 at one cluster to 5,470.8 at eight, and a center on every listing brings it to 0.
• DON'T compare inertias across tables. On the five columns without the price, three clusters reach 7,226.4, below the six columns' 9,050.1 because one column fewer adds up, and against each table's one-cluster inertia, 15,000.0 and 18,000.0, the two leave 0.482 and 0.503 of the spread.

Reading the Centers

What. cluster_centers_ holds each center in the standardized units KMeans was fitted on. The scaler's inverse_transform, Lesson 3.2.1's undo, maps them back to square meters, rooms, bathrooms, years, kilometers and currency, and each center of a loop run to its end is the mean of its listings' values as the imputer left them. Why. A cluster is acted on through what it holds. A center in standard deviations tells a seller nothing, and the same center in square meters and currency describes its cluster in one row. When. Every time a clustering is reported, and before any cluster is given a name.

Two clusters, read in listing units. The fit takes tol=1e-6, the second subtitle's tolerance, so that every center is exactly the mean of its listings:

reader = make_pipeline(
    SimpleImputer(strategy="median"),
    StandardScaler(),
    KMeans(n_clusters=2, random_state=42, n_init=10, tol=1e-6),
)
groups = reader.fit_predict(listings[cols])
print(np.round(reader["kmeans"].cluster_centers_, 2))

[[ 1.1   1.1   1.03  0.2   0.43  0.91]
 [-0.56 -0.56 -0.52 -0.1  -0.22 -0.46]]

In standard deviations, cluster 0 sits above the average listing in every column and cluster 1 below it. In listing units:

in_units = reader["standardscaler"].inverse_transform(
    reader["kmeans"].cluster_centers_
)
described = pd.DataFrame(in_units, columns=cols)
described["price"] = described["price"].round().astype(int)
described["listings"] = np.bincount(groups)
print(described.round(1))

   area_sqm  rooms  bathrooms  year_built  distance_km   price  listings
0     146.0    5.5        2.5      1994.8          6.8  558041      1015
1      74.4    2.6        1.3      1988.0          4.8  290114      1985

1,015 listings around 146.0 square meters, 5.5 rooms and 2.5 bathrooms at 558,041, and 1,985 around 74.4 square meters, 2.6 rooms and 1.3 bathrooms at 290,114. No listing has 5.5 rooms: a center is an average of listings, and need not match any one of them. Each center against the listings it holds:

seen = reader["simpleimputer"].transform(listings[cols])
means = pd.DataFrame(seen, columns=cols).groupby(groups).mean()
print(np.allclose(means, in_units))  # True
measured = listings.groupby(groups)["area_sqm"].mean()
print(measured.round(1).tolist())    # [148.5, 73.7]
print(listings.groupby(groups)["price"].median().tolist())

[520000.0, 276000.0]

Each center is the mean of its cluster's listings as the fit saw them, gaps filled. Cluster 0's center sits at 146.0 square meters while its measured areas average 148.5, because its listings with no area count at the median of 89.2 square meters, and in cluster 1 the same fill lifts 73.7 to 74.4. The price has no gaps, so a center's price is its cluster's mean price, 558,041 in cluster 0, whose median price is 520,000: a mean of prices, not a typical price and not a price per anything.

pd.crosstab (a table that counts the rows holding each pair of values from two columns, one row per value of the first and one column per value of the second) sets the clusters against the property types, a column the fit never saw, with rownames naming its rows:

print(pd.crosstab(groups, listings["property_type"], rownames=["cluster"]))

property_type  apartment  house  townhouse
cluster                                   
0                     57    663        295
1                   1535    148        302

1,535 of cluster 1's 1,985 listings are apartments, and 958 of cluster 0's 1,015 are houses and townhouses. That describes what the clusters hold, and it stops there. Cluster 0 is not the houses, since it holds 57 apartments and cluster 1 holds 450 houses and townhouses, and neither cluster is a kind of listing the market defines: k-means makes as many clusters as n_clusters asks for, as the first subtitle's eight showed, wherever its loop comes to rest. Why cluster 0's listings are both large and expensive, the clustering does not say, and this lesson claims no cause for it.

DON'T:

• DON'T name a cluster after a property type. Cluster 0 holds 663 houses and 295 townhouses and also 57 apartments, while 148 houses and 302 townhouses sit in cluster 1.
• DON'T read a center as a listing, or its price as a typical price. Cluster 0's center has 5.5 rooms, which no listing has, and its 558,041 is a mean where the cluster's median price is 520,000.
• DON'T read a center's area as its listings' measured average. Cluster 0's center sits at 146.0 square meters and its measured areas average 148.5, because the filled areas count at the median.

Lab Exercises

Lab Exercise 1: The Loop on Size and Price

In size_loop.py in ~/dev/bootcamp, load listings.csv and build S from four columns, the size columns area_sqm, rooms and bathrooms and the price, filled with their medians and standardized as the lesson built Z.

Fit KMeans(n_clusters=3, random_state=42) on S once with no target and once with listings["price"] passed as the target, and print whether every listing got the same label from both.

Write two functions: assign(S, centers), which returns each listing's nearest center, and move(S, labels, k), which returns the mean of each cluster's listings. Starting from the listings at positions 0, 1 and 2, alternate the two until an assignment changes no listing, recording the inertia after every move. Print the number of rounds, the last assignment counted, the first and last inertia recorded, to one decimal, whether no inertia rose above the one before it, and each cluster's size.

Finally fit KMeans from the same three listings with n_init=1, once at its default tol of 1e-4 and once with tol=1e-6. For each, print on one line the inertia rounded to one decimal, n_iter_, how many listings it labels as your loop does, and whether its centers are the means of its own clusters. End with a comment saying which fit stopped where your loop stopped, and why the other did not.

Lab Exercise 2: Four Clusters in Listing Units

In four_clusters.py, load listings.csv and keep the lesson's six columns in cols.

Show KMeans(n_clusters=4, random_state=42) refusing the raw columns by printing the first line of its message. Then build the lesson's pipeline with KMeans(n_clusters=4, random_state=42, n_init=10, tol=1e-6) as its last step, fit it on the raw columns with fit_predict, and print how many listings each cluster holds and the inertia rounded to one decimal.

Print the four centers in listing units, the price in whole units and the other columns to one decimal, check that every center is the mean of its cluster's filled values, and count each cluster's listings by condition with pd.crosstab.

Then take the table the pipeline's first two steps return, check its one-cluster inertia against 3,000 times its number of columns, and print each n_init=10 fit's inertia for k from 1 to 6 as a share of it, to three decimals.

Finally fit the pipeline again without its scaler, at the same k, and print each cluster's lowest and highest price and its size. End with a comment saying what the condition counts tell you about the scaled clusters, and what the unscaled clusters describe.

import numpy as np
import pandas as pd
from sklearn.cluster import KMeans
from sklearn.impute import SimpleImputer
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

listings = pd.read_csv("listings.csv")
size_price = ["area_sqm", "rooms", "bathrooms", "price"]
prep = make_pipeline(SimpleImputer(strategy="median"), StandardScaler())
S = prep.fit_transform(listings[size_price])

plain = KMeans(n_clusters=3, random_state=42).fit(S)
priced = KMeans(n_clusters=3, random_state=42).fit(S, listings["price"])
# True, the target was never read
print(np.array_equal(plain.labels_, priced.labels_))

def assign(S, centers):
    gaps = S[:, np.newaxis, :] - centers[np.newaxis, :, :]
    return np.linalg.norm(gaps, axis=2).argmin(axis=1)


def move(S, labels, k):
    return np.array([S[labels == j].mean(axis=0) for j in range(k)])


centers = S[:3]               # listings 0, 1 and 2
labels = np.full(len(S), -1)  # a label no center has
history = []                  # inertia after every move
rounds = 0
while True:
    rounds += 1
    nearest = assign(S, centers)
    if np.array_equal(nearest, labels):
        break                 # no listing changed cluster
    labels = nearest
    centers = move(S, labels, 3)
    history.append(((S - centers[labels]) ** 2).sum())
print(rounds, round(history[0], 1), round(history[-1], 1))
# never rising means the list is already sorted from largest down
print(history == sorted(history, reverse=True))
print(np.bincount(labels))

15 4264.6 3772.0
True
[ 443 1516 1041]

for tol in [1e-4, 1e-6]:           # the default, then smaller
    fit = KMeans(n_clusters=3, init=S[:3], n_init=1, tol=tol).fit(S)
    agree = (fit.labels_ == labels).sum()
    own = move(S, fit.labels_, 3)  # its own clusters' means
    is_mean = np.allclose(own, fit.cluster_centers_)
    print(round(fit.inertia_, 1), fit.n_iter_, agree, is_mean)
# tol=1e-6 stopped where the loop did, after 15 rounds, every center
# its cluster's mean. The default stopped after 11, once its centers
# moved less than its tolerance, and never moved them again.

3772.1 11 2995 False
3772.0 15 3000 True

import numpy as np
import pandas as pd
from sklearn.cluster import KMeans
from sklearn.impute import SimpleImputer
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

listings = pd.read_csv("listings.csv")
cols = ["area_sqm", "rooms", "bathrooms", "year_built", "distance_km", "price"]
try:
    KMeans(n_clusters=4, random_state=42).fit(listings[cols])
except ValueError as err:
    print(str(err).splitlines()[0])       # Input X contains NaN.

pipe = make_pipeline(
    SimpleImputer(strategy="median"),
    StandardScaler(),
    KMeans(n_clusters=4, random_state=42, n_init=10, tol=1e-6),
)
clusters = pipe.fit_predict(listings[cols])
print(np.bincount(clusters))              # [1015  418  899  668]
print(round(pipe["kmeans"].inertia_, 1))  # 7742.8

in_units = pipe["standardscaler"].inverse_transform(
    pipe["kmeans"].cluster_centers_
)
centers = pd.DataFrame(in_units, columns=cols)
centers["price"] = centers["price"].round().astype(int)
print(centers.round(1))
filled = pipe["simpleimputer"].transform(listings[cols])
means = pd.DataFrame(filled, columns=cols).groupby(clusters).mean()
print(np.allclose(means, in_units))  # True
print(pd.crosstab(clusters, listings["condition"], rownames=["cluster"]))

   area_sqm  rooms  bathrooms  year_built  distance_km   price
0      67.5    2.3        1.1      2000.9          6.3  263648
1     174.2    6.5        3.1      1993.7          6.7  670756
2     117.1    4.4        2.0      1995.1          6.6  451143
3      73.7    2.6        1.3      1965.7          2.1  282533
condition  excellent  fair  good  poor
cluster                               
0                166   333   400   116
1                 87   113   186    32
2                144   271   416    68
3                114   193   293    68

Z = pipe[:-1].transform(listings[cols])  # filled, then scaled
base = KMeans(n_clusters=1, random_state=42).fit(Z).inertia_
print(round(base, 3), Z.shape[0] * Z.shape[1])
for k in range(1, 7):
    fitted = KMeans(n_clusters=k, random_state=42, n_init=10).fit(Z)
    print(k, round(fitted.inertia_ / base, 3))

unscaled = make_pipeline(
    SimpleImputer(strategy="median"),
    KMeans(n_clusters=4, random_state=42, n_init=10, tol=1e-6),
)
raw = unscaled.fit_predict(listings[cols])
print(listings.groupby(raw)["price"].agg(["min", "max", "size"]))
# Every scaled cluster holds listings in all four conditions, from 32
# to 416 of them, so no cluster is a condition. The unscaled clusters
# are four price bands, 60,000 to 291,000, 292,000 to 475,000, 476,000
# to 743,000 and 747,000 to 1,930,000: they describe the price alone.

18000.0 18000
1 1.0
2 0.614
3 0.503
4 0.43
5 0.38
6 0.343
      min      max  size
0  292000   475000  1131
1  747000  1930000   163
2  476000   743000   577
3   60000   291000  1129

Theory Assessment

Q1: What fit Took, and What Each Step Does

The lesson fitted KMeans(n_clusters=2, random_state=42) on Z with no target, then again with listings["price"] as the target, and every label came back the same. Explain why, what labels_ and cluster_centers_ hold, and why listing 0 belongs to cluster 0. Then explain why one round of the loop, an assignment and then a move, can never raise the inertia, and why the loop stops once an assignment changes no listing.

Hint: For the first part, ask what KMeans.fit does with a y, and what a cluster's number refers to. For the second, ask what each step holds fixed while it changes the other, and what Lesson 3.1.3 showed the mean does to a sum of squared distances.

Q2: Four Places to Stop

Four runs stopped in four places. The loop from the first three listings came to rest after 45 rounds at 9,050.1. KMeans from the same three listings stopped after 41 rounds with 2 listings in other clusters. Seed 42's default single run stopped after 12 rounds at 9,439.0. At k=4 the loop from listings 32 to 35 came to rest at 8,017.5, where listings 0 to 3 led to 7,742.8. For each of the last three, say whether the tolerance stopped it or it rested at a local minimum, and which run tells you. Then say what n_init=10 does about both.

Hint: A stop forced by the tolerance moves when the tolerance shrinks or the loop runs on. A local minimum moves under neither, and only a different start does better.

Q3: Price Bands, 18,000 and a Center With 5.5 Rooms

Explain, from the spreads the scaler learned, why the pipeline without its scaler returned three price bands that do not overlap while the scaled pipeline's price ranges overlapped. Then explain why the one-cluster inertia on Z is 18,000, and why inertia fell at every k from 1 to 8. Finally say what cluster 0's center of 146.0 square meters, 5.5 rooms and 558,041 describes, and why the crosstab does not make cluster 0 the houses.

Hint: Compare each column's share of the summed variances before and after standardizing. For 18,000, ask where the one center sits on Z and what each column's squared values add up to. For the center, ask what it is the mean of.

Q4: The Same Listings, Priced in Millions

A teammate reruns the lesson's unscaled clustering and reports that its three price bands prove price is what separates the listings on the site. A second teammate keeps prices in millions, as the finance team does, and runs the same comparison. Read this script without running it.

import pandas as pd
from sklearn.cluster import KMeans
from sklearn.impute import SimpleImputer
from sklearn.pipeline import make_pipeline

listings = pd.read_csv("listings.csv")
cols = ["area_sqm", "rooms", "bathrooms", "year_built", "distance_km", "price"]
in_millions = listings[cols].assign(price=listings["price"] / 1_000_000)

for table in [listings[cols], in_millions]:
    unscaled = make_pipeline(
        SimpleImputer(strategy="median"),
        KMeans(n_clusters=3, random_state=42, n_init=10),
    )
    groups = unscaled.fit_predict(table)
    print(f"{unscaled['kmeans'].inertia_:.3e}")
    print(listings.groupby(groups)["price"].agg(["min", "max"]))

The run on the first table is the lesson's unscaled clustering. Say what the second table's prints look like and why. Then answer the first teammate, a third teammate who reads the second inertia as proof of tighter clusters, and the finance team, who ask whether the lesson's scaled pipeline would group the listings differently with prices in millions.

Hint: Work out the price's spread in millions beside the area's in square meters, and ask which column the distances now listen to. Then ask what units each inertia is in, and what standardizing returns for a column whatever unit it was written in.

2.191e+13
      min      max
0  344000   613000
1   60000   343000
2  614000  1930000
2.414e+06
      min      max
0  101000  1930000
1   60000  1255000
2  262000  1619000

import numpy as np
from sklearn.preprocessing import StandardScaler

print(listings.groupby(groups)["area_sqm"].agg(["min", "max"]))
filled = SimpleImputer(strategy="median").fit_transform(in_millions)
print(np.round(filled.std(axis=0), 3))
scaled = make_pipeline(
    SimpleImputer(strategy="median"),
    StandardScaler(),
    KMeans(n_clusters=3, random_state=42, n_init=10),
)
in_currency = scaled.fit_predict(listings[cols])
in_millions_too = scaled.fit_predict(in_millions)
print(np.array_equal(in_currency, in_millions_too))

     min    max
0   78.4  139.4
1   24.0   95.1
2  139.0  278.5
[43.176  1.708  0.815 22.12   3.042  0.195]
True

Assignment

Assignment: Typical Listings, and the Ones No Cluster Fits

The listings team wants the lesson's three clusters, from the third subtitle's pipeline, turned into material a seller can read and a reviewer can act on. No listing has 5.7 rooms, so for each cluster they want the real listing nearest its center, among the listings with all six values measured, because a seller sees a listing as it was posted. A reviewer wants the five listings farthest from their own cluster's center, to check by hand. And each new listing should be placed by the fitted pipeline, with its distance to its cluster's center, and flagged when that distance is larger than 99 in 100 of the current listings' distances.

One call does most of the work. A fitted KMeans, and a pipeline that ends in one, has transform(X), which returns each row's distance to every center in the centers' standardized units, not squared: on the listings, a (3000, 3) array whose row i, column j is listing i's distance to center j.

Write typical.py and show, in this order:

• How many listings each cluster holds, and the three centers in listing units, the price in whole units and the rest to one decimal.
• That every listing's smallest distance is to its own cluster's center, and that the squares of those distances add up to inertia_.
• Each cluster's typical listing, its listing_id and six values, with its distance to three decimals.
• The five listings farthest from their own center, their listing_id and six values, then their clusters and their distances to two decimals.
• The 99th percentile of the distances to own centers, to three decimals, and how many listings lie beyond it.
• Three new listings, placed and flagged: one of 68 square meters with 2 rooms and 1 bathroom, built in 1985, 3.5 kilometers out, at 265,000, one of 150 square meters with 5 rooms and 2 bathrooms, 7.0 kilometers out, at 540,000, whose seller left the year blank, and one of 95 square meters with 3 rooms and 1 bathroom, built in 2020, 1.0 kilometer out, at 1,450,000.

End with a comment saying how a typical listing differs from its center and why the team should show the listing, and what the five farthest listings have in common, described without claiming a cause.

Hints:

• Every listing's own center is its nearest one, so dist.min(axis=1) is each listing's distance to its own center, and dist.argmin(axis=1) should give back the labels.
• Lesson 2.2.2's np.where(condition, values, np.inf) keeps a distance only where the condition holds, so its argmin never picks a listing outside the cluster or with a gap, and np.argsort reversed with [::-1] puts the farthest first.

import numpy as np
import pandas as pd
from sklearn.cluster import KMeans
from sklearn.impute import SimpleImputer
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

listings = pd.read_csv("listings.csv")
cols = ["area_sqm", "rooms", "bathrooms", "year_built", "distance_km", "price"]
pipe = make_pipeline(
    SimpleImputer(strategy="median"),
    StandardScaler(),
    KMeans(n_clusters=3, random_state=42, n_init=10),
)
clusters = pipe.fit_predict(listings[cols])
print(np.bincount(clusters))  # [ 942 1197  861]
in_units = pipe["standardscaler"].inverse_transform(
    pipe["kmeans"].cluster_centers_
)
centers = pd.DataFrame(in_units, columns=cols)
centers["price"] = centers["price"].round().astype(int)
print(centers.round(1))

   area_sqm  rooms  bathrooms  year_built  distance_km   price
0      72.6    2.6        1.3      1975.9          2.2  294206
1      80.9    2.8        1.4      1999.0          7.3  302048
2     151.7    5.7        2.7      1994.0          6.6  585215

dist = pipe.transform(listings[cols])  # (3000, 3), not squared
own = dist.min(axis=1)                 # distance to own center
print(np.array_equal(dist.argmin(axis=1), clusters))
print(round((own ** 2).sum(), 1), round(pipe["kmeans"].inertia_, 1))

True
9050.1 9050.1

complete = listings[cols].notna().all(axis=1).to_numpy()
show = ["listing_id"] + cols
typical = []
for c in range(3):
    # own distances in cluster c with no gap, inf everywhere else
    pool = np.where(complete & (clusters == c), own, np.inf)
    typical.append(pool.argmin())
print(listings.iloc[typical][show])
print(own[typical].round(3))
far = np.argsort(own)[::-1][:5]  # farthest first
print(listings.iloc[far][show])
print(clusters[far], own[far].round(2))

      listing_id  area_sqm  rooms  bathrooms  year_built  distance_km   price
58         30059      68.3      3          1      1974.0          2.1  275000
1167       31168      76.7      3          1      2003.0          7.6  302000
994        30995     149.9      6          3      1989.0          6.0  604000
[0.461 0.502 0.564]
      listing_id  area_sqm  rooms  bathrooms  year_built  distance_km    price
534        30535     117.1      4          2         NaN          5.5  1930000
26         30027     220.4      9          4      2017.0          4.7  1619000
2963       32964     273.6     10          5      1978.0         10.9  1059000
124        30125      96.0      3          1      1985.0          0.3  1410000
2950       32951     208.4      6          3      1996.0          5.6  1544000
[2 2 2 2 2] [7.06 6.21 5.57 5.53 5.12]

cutoff = np.percentile(own, 99)
print(round(cutoff, 3), (own > cutoff).sum())
new = pd.DataFrame({
    "area_sqm": [68.0, 150.0, 95.0],
    "rooms": [2, 5, 3],
    "bathrooms": [1, 2, 1],
    "year_built": [1985.0, np.nan, 2020.0],
    "distance_km": [3.5, 7.0, 1.0],
    "price": [265000, 540000, 1450000],
})
placed = pipe.predict(new)
new_own = pipe.transform(new).min(axis=1)
print(placed, new_own.round(2), new_own > cutoff)
# A typical listing is one a seller can open: cluster 2's has 6 rooms
# and 3 bathrooms at 604,000, where its center has 5.7 rooms, 2.7
# bathrooms and the cluster's mean price, 585,215. The five farthest
# listings all sit in cluster 2, priced from 1,059,000 up to 1,930,000,
# the highest price in the file, and the first was posted with no
# year. That says no center describes them well, and nothing about why.
# The third new listing, 5.72 from its center, goes to the reviewer.

3.505 30
[0 2 2] [0.79 0.94 5.72] [False False  True]

Additional Resources

Optional, and nothing later depends on either. Worth opening if you want a clustering method that needs no k, or the shapes of data k-means splits badly.

• YouTube - "Clustering with DBSCAN, Clearly Explained!!!" by StatQuest with Josh Starmer. A clustering method that grows clusters from points packed closely together, needs no number of clusters, and leaves points far from every cluster out of all of them.
• Official docs - scikit-learn 1.9 example, "Demonstration of k-means assumptions". Four made-up datasets where k-means returns clusters that look wrong, too few clusters, stretched blobs, unequal spreads and unequal sizes, each with the remedy the page suggests.

Summary

k-means turned out to be two earlier lessons joined in a loop. The assignment is Lesson 2.3.3's nearest point and the move is Lesson 3.1.3's mean, so k-means fits the way LinearRegression does, by lowering a sum of squares. What it lacks is the column holding the answer, and with it went the honest score that checks held-out rows against their answers, leaving one number about the fit itself that improves with every center added.

What took the answer's place were the choices around the loop. Where it starts was a random draw, when it stops was a tolerance that halted one run at 9,439.0 on centers still drifting toward 9,050.1, and which columns count was settled by their units until a scaler put them on one footing. Each changed the grouping without a word, and each was caught only by a run that took it away.

What comes back is a description. A center is an average that no listing need match, a cluster is the set of listings nearest one such average, and a cluster of mostly houses still holds apartments. Whether the clusters survive another start and how cleanly they separate are the questions Lesson 3.8.3 asks, which is as near as a method without a target comes to being checked.

Cheat Sheet

Part A - Reference Table

Part B - Compact Code Block

pipe = make_pipeline(
    SimpleImputer(strategy="median"),
    StandardScaler(),
    KMeans(n_clusters=3, random_state=42, n_init=10),
)
labels = pipe.fit_predict(listings[cols])  # no y, a label per row
pipe["kmeans"].inertia_                    # squared distances, summed
pipe.predict(new)                          # place rows, no refit
pipe.transform(new)                        # distance to every center
C = pipe["kmeans"].cluster_centers_        # standardized units
# back to square meters, years and currency
pipe["standardscaler"].inverse_transform(C)
# count each cluster's rows per property type
pd.crosstab(labels, listings["property_type"])