The Clustering Problem

The Supervised Learning course tried seven model families in twenty lessons on the same fictional subscriber table and wrote its budget next to every number. All of those measurements shared one foundation that was never questioned — the label. The baseline predicted the most frequent label, the measure compared the predicted label with the true one, the tuning search read the validation labels. Remove the label column and all three mechanisms stop working at once.

This lesson starts there. There is a table without labels and the goal is to gather similar subscribers together. The trouble is not finding a method — producing a partition is easy, and random assignment produces one too. The trouble is what the partition is measured against. Before the measure is built, there is an even more basic decision: the definition of “similar”.

• CL1. A cluster, in this course, is a subgroup of points; the database cluster of the Databases curriculum was a group of servers, while the dataset is all the rows.
• CL2. The dataset is the fictional subscriber table from the Introduction to Machine Learning, Data Preparation and Feature Engineering, and Supervised Learning courses: seed 20260218, 1,260 subscribers, split 756/252/252. The procedure works on the 756 subscribers of the training set.
• CL3. The suspicious column is not generated in this course. No method sees the label; the label will be used only once, in the last lesson of the topic, when external validity is measured.
• CL4. The columns are six numeric fields: mean consumption, volatility, per-person consumption, household size, satisfaction and period count. Region indicators are left out, because the “distance” between two regions is always the same constant in binary columns.
• CL5. A partition is an arrangement that assigns every point to exactly one cluster. Its measure is the within-cluster sum: the sum of every point’s distance to the centroid of its own cluster, divided by the number of points; smaller is better. The numbers printed in the application areas lesson of the Introduction to Machine Learning course are not repeated; they came from a different set of columns.
• CL6. Three baselines are used: random assignment, single cluster and slicing by a single column. The number for random assignment is the mean of five seeds; a single seed hides its fluctuation.
• CL7. The number of clusters is held at 4, and so is the procedure: four starting points are chosen with farthest-point seeding and two refinement rounds are run. The number of clusters as a hyperparameter and the starting point moving the result belong to the next lesson; here only the similarity definition changes.
• CL8. The numbers of different definitions cannot be compared with each other, because each measures its own distance; only the difference against its own baseline is comparable.

Similarity Definition Is a Choice

A clustering procedure does not learn what “close” means from the data; it takes it from outside. The first definition leaves the columns on the raw scale and sums squared differences: because the range of mean consumption is about eight and a half times as wide as the range of household size, the sum is quietly a weighting. The second turns every column into a standard score, so the six columns get equal weight. The third turns every value into the percentage of its rank within the column and sums absolute differences; extreme values are squeezed and only the ordering remains.

The definition determines not only the distance but also what the center is: the center that minimizes the squared distance is the mean vector, and the one that minimizes the absolute difference is the median vector. The three definitions do not pick the same nearest neighbor even for a single subscriber.

# clustering.py — MODEL: the FICTIONAL subscriber table of the Introduction to Machine Learning, Data Preparation and Feature Engineering and Supervised Learning courses is rebuilt with the same seed;
# the suspicious column is not generated in this course, no procedure sees the label.
import math
import statistics

SEED, RAW, M32 = 20260218, 1400, 0xFFFFFFFF
REGION = [(0.28, 21), (0.22, 17), (0.18, 26), (0.14, 14), (0.18, 23)]


def rng(t):
    x = ((t ^ (t >> 16)) * 2246822507) & M32
    x = ((x ^ (x >> 13)) * 3266489909) & M32
    s = [(x ^ (x >> 16)) & M32]

    def next_value():
        s[0] = (s[0] * 1664525 + 1013904223) & M32
        return s[0] / 4294967296
    return next_value


def pick(u, w):
    t = 0.0
    for i, x in enumerate(w):
        t += x
        if u < t:
            return i
    return len(w) - 1


DATA = []
for i in range(RAW):
    r = rng(SEED + i)
    base = REGION[pick(r(), [x[0] for x in REGION])][1]
    household = 1 + pick(r(), [0.06, 0.24, 0.30, 0.24, 0.11, 0.05])
    rating = 1 + pick(r(), [0.08, 0.14, 0.27, 0.34, 0.17])
    if r() < 0.046:
        continue
    r, v = rng(SEED + 17001 + i), []
    for d in range(pick(r(), [0.05, 0.12, 0.21, 0.62])):
        v.append(0.0 if r() < 0.038 else math.floor(
            base * math.exp((r() + r() + r() - 1.5) * 0.62)
            * (1 - d * 0.05) * 100 + 0.5) / 100)
        r()
    if v:
        x = {"mean_consumption": round(sum(v) / len(v), 2), "volatility": round(max(v) - min(v), 2),
             "household": household, "satisfaction": rating, "period": len(v)}
        x["per_person"] = round(x["mean_consumption"] / household, 3)
        DATA.append(x)

FIELDS = ["mean_consumption", "volatility", "per_person", "household", "satisfaction", "period"]


def shuffle(data, seed):
    r, s = rng(seed), list(range(len(data)))
    for i in range(len(s) - 1, 0, -1):
        j = int(r() * (i + 1))
        s[i], s[j] = s[j], s[i]
    return [data[i] for i in s]


TRAIN = shuffle(DATA, SEED + 90000)[:756]
N = len(TRAIN)
RAW_VALUES = [[x[a] for a in FIELDS] for x in TRAIN]
MEANS = [statistics.fmean(c) for c in zip(*RAW_VALUES)]
STD_DEVS = [statistics.pstdev(c) for c in zip(*RAW_VALUES)]
STD_VALUES = [[(v[j] - MEANS[j]) / STD_DEVS[j] for j in range(len(FIELDS))] for v in RAW_VALUES]


def rank_percent(column):                   # value turns into the percentage of its mean rank in the column
    s, place = sorted(range(len(column)), key=lambda i: column[i]), {}
    for p, i in enumerate(s):
        place.setdefault(column[i], []).append(p)
    mid = {v: statistics.fmean(p) / (len(s) - 1) for v, p in place.items()}
    return [mid[x] for x in column]


RANK_COLUMNS = [rank_percent(c) for c in zip(*RAW_VALUES)]
RANK_VALUES = [[RANK_COLUMNS[j][i] for j in range(len(FIELDS))] for i in range(N)]


squared = lambda a, b: sum((p - q) * (p - q) for p, q in zip(a, b))
absolute = lambda a, b: sum(abs(p - q) for p, q in zip(a, b))
mean_centroid = lambda g: [statistics.fmean(c) for c in zip(*g)]
median_centroid = lambda g: [statistics.median(c) for c in zip(*g)]

DEFINITIONS = [("raw scale, squared distance", RAW_VALUES, squared, mean_centroid),
               ("standard score, squared distance", STD_VALUES, squared, mean_centroid),
               ("rank percentage, absolute difference", RANK_VALUES, absolute, median_centroid)]

print(f"clustered subscribers {N}, columns {len(FIELDS)}, no label column")
print(f"{'':<38}{'neighbor':>10}" + "".join(f"{a:>18}" for a in FIELDS))
print(f"{'subscriber 48 (itself)':<38}{'--':>10}" + "".join(f"{TRAIN[48][a]:>18.3f}" for a in FIELDS))
for name, V, dist, _ in DEFINITIONS:
    j = min((x for x in range(N) if x != 48), key=lambda x: (dist(V[48], V[x]), x))
    print(f"{name:<38}{j:>10}" + "".join(f"{TRAIN[j][a]:>18.3f}" for a in FIELDS))

clustered subscribers 756, columns 6, no label column
                                        neighbor  mean_consumption        volatility        per_person         household      satisfaction            period
subscriber 48 (itself)                        --            16.090            10.560             3.218             5.000             4.000             3.000
raw scale, squared distance                  711            16.280            10.600             4.070             4.000             4.000             2.000
standard score, squared distance              52            13.540            11.990             2.708             5.000             4.000             3.000
rank percentage, absolute difference         241            17.400             9.830             2.900             6.000             4.000             3.000

The three rows give three different neighbors for the same subscriber, and the reason is visible in each row. Subscriber 711, chosen on the raw scale, is only 0.19 away from subscriber 48 in mean consumption; the household and period columns do not match, but their differences add at most 2 to the sum. Subscriber 52, chosen with the standard score, is 2.55 away in consumption, but the household, satisfaction and period columns all match: a perfect match in one column offsets a difference of one standard deviation.

None of the three neighbors is wrong. All three fit “nearest”; only for three different definitions.

Three Baselines

The within-cluster sum cannot be interpreted alone; whether the number is large or small cannot be known without a comparison point from the same measure. Three baselines answer three questions. Random assignment produces the same number of clusters without any structure: “how much better is this partition than chance”. A single cluster keeps the points together: “what did splitting gain over not splitting at all”. Slicing by a single column splits the points into four equal slices by mean consumption: “what did looking at all the columns gain over a solution that looks at a single column”. The third is the hardest baseline: it uses the data even though it does not learn.

CLUSTERS = 4


def within_sum(V, dist, centroid, assign):  # within-cluster sum, per point
    t = 0.0
    for c in set(assign):
        g = [V[i] for i in range(len(V)) if assign[i] == c]
        t += sum(dist(v, centroid(g)) for v in g)
    return t / len(V)


def random_assignment(n, k, seed):          # BASELINE 1: uses no structure
    r = rng(seed)
    return [int(r() * k) for _ in range(n)]


def by_column(V, j, k):                     # BASELINE 3: splits a single column into k equal slices
    s, a = sorted(range(len(V)), key=lambda i: (V[i][j], i)), [0] * len(V)
    for p, i in enumerate(s):
        a[i] = min(k - 1, p * k // len(V))
    return a


def partition(V, dist, centroid, k, rounds=2):  # MODEL: farthest-point seeding + two rounds
    starts = [min(range(len(V)), key=lambda i: (dist(V[i], centroid(V)), i))]
    while len(starts) < k:
        starts.append(max(range(len(V)),
                          key=lambda i: (min(dist(V[i], V[c]) for c in starts), -i)))
    M = [V[c] for c in starts]
    for _ in range(rounds + 1):
        a = [min(range(k), key=lambda c: dist(v, M[c])) for v in V]
        M = [centroid([V[i] for i in range(len(V)) if a[i] == c]) for c in range(k)]
    return a


_, V2, D2, C2 = DEFINITIONS[1]
RND = [within_sum(V2, D2, C2, random_assignment(N, CLUSTERS, SEED + 400 + t)) for t in range(5)]
A2 = partition(V2, D2, C2, CLUSTERS)
single, rnd = within_sum(V2, D2, C2, [0] * N), statistics.fmean(RND)
sliced = within_sum(V2, D2, C2, by_column(V2, 0, CLUSTERS))
method = within_sum(V2, D2, C2, A2)
print(f"with the standard score, squared distance definition, k = 4\n"
      f"  baseline - single cluster          {single:8.4f}\n"
      f"  baseline - random assignment       {rnd:8.4f}   (five seeds: {min(RND):.4f} - {max(RND):.4f})\n"
      f"  baseline - mean_consumption slice  {sliced:8.4f}\n"
      f"  method                             {method:8.4f}"
      f"   cluster sizes {sorted(A2.count(c) for c in range(CLUSTERS))}\n"
      f"  difference: against random {rnd - method:+.4f}, against slice {sliced - method:+.4f}")

with the standard score, squared distance definition, k = 4
  baseline - single cluster            6.0000
  baseline - random assignment         5.9790   (five seeds: 5.9763 - 5.9819)
  baseline - mean_consumption slice    4.9044
  method                               3.8910   cluster sizes [61, 67, 117, 511]
  difference: against random +2.0880, against slice +1.0134

The number for the single cluster is exactly 6.0000, that is, the number of columns. This is no coincidence: a standard score carries unit variability, so the per-point sum of the distance to a single center equals the number of columns. The measure’s upper bound is known by calculation, and every partition is read by how far down it pulls this 6.0000.

The second row is the most important number of this course. Random assignment gives 5.9790, only 0.0210 below the single cluster. Four random groups gain almost nothing over not splitting, because the center of a random cluster falls on the overall center of the data. The fluctuation across the five seeds is 0.0056, negligible beside the measured differences.

The third row changes the table. Slicing into four equal slices by mean consumption alone gives 4.9044. The method’s number is 3.8910 and the difference is 2.0880 or 1.0134, depending on which baseline is written. Had the slice not been written, the method would have credited itself with the 1.0134 that looking at all six columns gained. The choice of baseline is a measurement decision; when it is not written with its reasoning, the reported difference is not a finding but a choice.

The Same Cluster, Three Partitions

The third block measures what happens when the procedure and the number of clusters are held fixed and only the similarity definition changes. The disagreement of two partitions is counted over point pairs, without matching the clusters: 756 subscribers have 285,390 pairs, and for every pair each partition says either “in the same cluster” or “in different clusters”. The disagreement is the ratio of the pairs on which they answer differently.

PARTITIONS = []
print(f"{'definition':<38}{'single cluster':>16}{'random':>16}{'slice':>16}"
      f"{'method':>16}{'share over random':>21}")
for name, V, dist, centroid in DEFINITIONS:
    a = partition(V, dist, centroid, CLUSTERS)
    PARTITIONS.append((name, a))
    single = within_sum(V, dist, centroid, [0] * N)
    rnd = statistics.fmean(within_sum(V, dist, centroid, random_assignment(N, CLUSTERS, SEED + 400 + t))
                           for t in range(5))
    cut = within_sum(V, dist, centroid, by_column(V, 0, CLUSTERS))
    y = within_sum(V, dist, centroid, a)
    print(f"{name:<38}{single:>16.4f}{rnd:>16.4f}{cut:>16.4f}{y:>16.4f}{(rnd - y) / rnd:>21.3f}")


def disagreement(a, b):                     # on how many point pairs the two partitions differ
    differing = sum((a[i] == a[j]) != (b[i] == b[j])
                    for i in range(len(a)) for j in range(i + 1, len(a)))
    return 2 * differing / (len(a) * (len(a) - 1))


print("\npair disagreement between partitions (ratio of differing point pairs)")
for i in range(3):
    for j in range(i + 1, 3):
        print(f"  {PARTITIONS[i][0]:<38} - {PARTITIONS[j][0]:<38} {disagreement(PARTITIONS[i][1], PARTITIONS[j][1]):.4f}")
print("\nmeans of the largest cluster")
for name, a in PARTITIONS:
    g = [TRAIN[i] for i in range(N) if a[i] == max(a, key=a.count)]
    print(f"  {name:<38} n {len(g):>4}" + "".join(
        f"  {x}={statistics.fmean(k[x] for k in g):.2f}"
        for x in ("mean_consumption", "volatility")))

definition                              single cluster          random           slice          method    share over random
raw scale, squared distance                   130.6679        130.2948         92.4562         61.8604                0.525
standard score, squared distance                6.0000          5.9790          4.9044          3.8910                0.349
rank percentage, absolute difference            1.3983          1.3907          1.1475          0.9809                0.295

pair disagreement between partitions (ratio of differing point pairs)
  raw scale, squared distance            - standard score, squared distance       0.4192
  raw scale, squared distance            - rank percentage, absolute difference   0.4801
  standard score, squared distance       - rank percentage, absolute difference   0.3750

means of the largest cluster
  raw scale, squared distance            n  509  mean_consumption=17.40  volatility=6.86
  standard score, squared distance       n  511  mean_consumption=19.18  volatility=11.97
  rank percentage, absolute difference   n  325  mean_consumption=21.50  volatility=12.66

One rule applies when reading the first table: no comparison across rows. 61.8604 and 0.9809 are not two values of the same measure; one is on the order of cubic meters squared, the other a sum of percentage differences. The only comparable column is the share on the right: each definition’s gain over its own random baseline, 0.525, 0.349 and 0.295 in order. The share is highest on the raw scale, but that does not make the definition better — there almost all of the total scatter comes from a single column, so splitting is the same as cutting that column. The share measures not the success of the definition but how one-sidedly it looks at the data.

The second table is the real result of the lesson. The three partitions part on a share of point pairs between 0.3750 and 0.4801. The 0.4801 between the raw scale and the rank percentage means the two partitions give opposite answers on 137,002 of the 285,390 pairs. These are three runs on the same data, with the same procedure and number of clusters; only the definition of distance changes.

The third table shows where the difference sits. On the raw scale the largest cluster has 509 subscribers and a mean volatility of 6.86 — mid-consumption, low-volatility subscribers together. With the standard score the cluster is of similar size but its mean volatility is 11.97: the volatility column now has as much say as consumption. With the rank percentage the largest cluster drops to 325, because extreme values, once turned into ranks, do not push points apart.

No measure can say that one of these three partitions is the correct one. Each definition finds its own partition best on its own measure, and no higher court can be set up between definitions, because setting one up means choosing a distance definition once more. In the supervised course the label closed this question; here it stays open.

Summary

• When the label is removed, all three components of the measure must be built: the baseline, the method’s number and the difference. In clustering the baseline is not given ready-made; it is chosen and justified.
• Under the standard score definition a single cluster gives 6.0000 and random assignment 5.9790; a partition that does not learn moves the measure by only 0.0210 per point, and its fluctuation across five seeds is 0.0056.
• The method’s number is 3.8910; its difference is 2.0880 against the random baseline and 1.0134 against the slice baseline. Which baseline is written moves the reported gain nearly twofold.
• When three similarity definitions split the same 756 subscribers with the same procedure, three partitions come out that disagree on a share of point pairs between 0.3750 and 0.4801; the raw scale and the rank percentage give opposite answers on 137,002 of 285,390 pairs.
• The numbers of the definitions cannot be compared; the only thing that can be read is the proportional gain over their own baseline, and that share measures not the correctness of the definition but how one-sidedly it looks at the data.

Next Step

In this lesson the procedure was deliberately frozen: the starting points came from a single rule, two refinement rounds were run, the number of clusters was held at 4. All three are decisions. The next lesson resolves two of them: how many different partitions come out at the same number of clusters when the iteration continues until it stops and the starting points are chosen at random, and how large the gap is between the best and the worst run. After that the number of clusters is treated as a hyperparameter; because the within-cluster sum falls in only one direction as the number of clusters grows, it cannot point to the best value, and the “elbow” formed by the successive differences is counted as a decision point.