k-Means on Handwritten Digits

Math 123: Mathematical Aspects of Data Analysis — Fall 2026

Published

October 13, 2026

The scikit-learn digits dataset is $1797$ handwritten digits, each an $8\times8$ grayscale image with pixel values $0,\ldots,16$. Flattening each image gives a point $\mathbf{x}_i\in\mathbb{R}^{64}$. We run Lloyd’s algorithm for the k-means problem with $k=10$ on these points, without ever looking at the true labels, and use PCA only to watch what happens.

import numpy as np
import plotly.graph_objects as go
from plotly.subplots import make_subplots
from sklearn.datasets import load_digits
from sklearn.decomposition import PCA

digits = load_digits()
X = digits.data                    # (1797, 64), one flattened image per row
y = digits.target                  # true digit, used only to evaluate at the end
n, d = X.shape
k = 10

print(f"{n} images, {d} pixels each")

1797 images, 64 pixels each

Figure style

# A muted ink colour and a translucent grid read correctly against both the
# light and the dark version of this site, so the figures need no theme switch.
INK = "#7c7c78"
GRID = "rgba(128, 128, 128, 0.28)"
# Dark bold text with a white halo stays legible over the points in both themes.
LABEL_FONT = dict(color="#1a1a19", size=13, weight="bold",
                  shadow="0 0 3px white, 0 0 3px white, 0 0 2px white")
CLUSTER_COLOR = [
    "#2a78d6", "#eb6834", "#1baf7a", "#eda100", "#e87ba4",
    "#008300", "#6250d6", "#e34948", "#0f8b8d", "#8c564b",
]


def style(fig, height=420, title=None, showlegend=None):
    """Apply the shared layout to a plot with axes."""
    fig.update_layout(
        template="simple_white",
        paper_bgcolor="rgba(0,0,0,0)",
        plot_bgcolor="rgba(0,0,0,0)",
        font=dict(color=INK, size=13),
        title=dict(text=title, x=0, xanchor="left", font=dict(size=15)),
        margin=dict(l=70, r=30, t=60 if title else 30, b=55),
        height=height,
        showlegend=showlegend,
        hovermode="closest",
    )
    fig.update_xaxes(gridcolor=GRID, linecolor=GRID, zerolinecolor=GRID)
    fig.update_yaxes(gridcolor=GRID, linecolor=GRID, zerolinecolor=GRID)
    return fig

Lloyd’s algorithm in 64 dimensions

This is the lloyd function from the “k-means in practice” box after Lloyd’s algorithm, changed only to record the assignment, the centers, and the cost after every update step. We start from $k$ distinct data points chosen at random.

def lloyd_history(X, centers, max_iter=100):
    centers = np.array(centers, dtype=float)
    labels = None
    history = []
    for _ in range(max_iter):
        dists = ((X[:, None, :] - centers[None, :, :]) ** 2).sum(axis=2)  # (n, k)
        new_labels = dists.argmin(axis=1)          # assignment step
        if labels is not None and np.array_equal(new_labels, labels):
            break
        labels = new_labels
        for l in range(len(centers)):              # update step
            if np.any(labels == l):
                centers[l] = X[labels == l].mean(axis=0)
        cost = ((X - centers[labels]) ** 2).sum()
        history.append((labels.copy(), centers.copy(), cost))
    return history


rng = np.random.default_rng(0)
initial_centers = X[rng.choice(n, size=k, replace=False)]
history = lloyd_history(X, initial_centers)

for t, (_, _, cost) in enumerate(history, start=1):
    print(f"J_{t:<2d} = {cost:,.0f}")

J_1  = 1,516,554
J_2  = 1,337,221
J_3  = 1,254,299
J_4  = 1,237,883
J_5  = 1,228,740
J_6  = 1,218,411
J_7  = 1,203,291
J_8  = 1,182,637
J_9  = 1,170,242
J_10 = 1,167,157
J_11 = 1,165,931
J_12 = 1,165,552
J_13 = 1,165,322
J_14 = 1,165,160

The cost decreases at every step and the algorithm stops once the assignment repeats, exactly as the termination theorem promises.

Watching it in the PCA plane

We cannot plot points in $\mathbb{R}^{64}$, so we project onto the top $2$ principal directions, fit once on the data, and send the centers through the same projection. Since projection is linear (an affine map once centering is included), the projected center of a cluster is exactly the mean of its projected points, so each cross sits at the center of mass of its color in the plot.

pca = PCA(n_components=2).fit(X)
Z = pca.transform(X)                                       # (1797, 2)
print("explained variance ratio:", np.round(pca.explained_variance_ratio_, 3))

explained variance ratio: [0.149 0.136]

Plot each iteration, with a slider

frames = [(None, initial_centers, None)] + history
fig = go.Figure()
for t, (labels, centers, cost) in enumerate(frames):
    C = pca.transform(centers)
    colors = INK if labels is None else [CLUSTER_COLOR[l] for l in labels]
    fig.add_trace(go.Scatter(
        x=Z[:, 0], y=Z[:, 1], mode="markers", visible=(t == 0), showlegend=False,
        marker=dict(color=colors, size=5, opacity=0.35),
        customdata=np.c_[y, labels if labels is not None else np.full(n, -1)],
        hovertemplate="true digit %{customdata[0]}<br>cluster %{customdata[1]}<extra></extra>",
    ))
    fig.add_trace(go.Scatter(
        x=C[:, 0], y=C[:, 1], mode="markers+text", visible=(t == 0), showlegend=False,
        text=[f"c{l}" for l in range(k)], textposition="top center",
        textfont=LABEL_FONT,
        marker=dict(symbol="x", size=14, color=CLUSTER_COLOR,
                    line=dict(color="white", width=1.5)),
        hovertemplate="center %{text}<extra></extra>",
    ))

steps = []
for t, (_, _, cost) in enumerate(frames):
    visible = [False] * len(fig.data)
    visible[2 * t] = visible[2 * t + 1] = True
    title = ("Step 0: random initial centers" if cost is None
             else f"Step {t}: cost J = {cost:,.0f}")
    steps.append(dict(method="update", label=str(t),
                      args=[{"visible": visible}, {"title.text": title}]))

fig.update_layout(sliders=[dict(active=0, steps=steps, pad=dict(t=50),
                                currentvalue=dict(prefix="step "))])
fig.update_xaxes(title="z1 (PC1)", constrain="domain")       # shrink the plot, not pad the range
fig.update_yaxes(title="z2 (PC2)", scaleanchor="x", scaleratio=1)  # one unit is the same length on both axes
style(fig, height=600, title="Step 0: random initial centers")
fig.show()

Figure 1: Lloyd’s algorithm on the digits, viewed in the top-2 principal subspace. Drag the slider to step through the iterations; step 0 shows the random initial centers.

Two things are worth noticing.

Most of the motion happens in the first few steps; after that, only points near the boundaries between clusters change color.
The colors overlap in the plot even though k-means clusters are convex. The nearest-center regions live in $\mathbb{R}^{64}$, and the top two principal directions keep only about $29\%$ of the variance, so clusters separated in the full space land on top of each other in this shadow.

What the clusters found

Each center is the mean of a set of images, so it is itself an $8\times8$ image. Only now do we use the true labels, to see which digit dominates each cluster.

labels, centers, cost = history[-1]
for l in range(k):
    counts = np.bincount(y[labels == l], minlength=10)
    print(f"c{l}: {counts.sum():4d} images, {counts.max() / counts.sum():.0%} are the digit {counts.argmax()}")



def majority_fraction(labels):
    """Fraction of images whose cluster's most common digit is their own."""
    return sum(np.bincount(y[labels == l]).max() for l in np.unique(labels)) / n


frac_64 = majority_fraction(labels)
print(f"\nfraction of images whose cluster's majority digit is their own: {frac_64:.1%}")

c0:  177 images, 87% are the digit 3
c1:  147 images, 93% are the digit 5
c2:  211 images, 83% are the digit 7
c3:  249 images, 56% are the digit 9
c4:  169 images, 98% are the digit 4
c5:   87 images, 62% are the digit 1
c6:  174 images, 84% are the digit 2
c7:  178 images, 99% are the digit 0
c8:  182 images, 97% are the digit 6
c9:  223 images, 45% are the digit 1

fraction of images whose cluster's majority digit is their own: 79.2%

Plot the centers as images

fig = make_subplots(rows=2, cols=5, horizontal_spacing=0.02, vertical_spacing=0.12,
                    subplot_titles=[f"c{l}" for l in range(k)])
for l in range(k):
    fig.add_trace(go.Heatmap(z=centers[l].reshape(8, 8), colorscale="Greys",
                             zmin=0, zmax=16, showscale=False,
                             hovertemplate="%{z:.1f}<extra></extra>"),
                  row=l // 5 + 1, col=l % 5 + 1)
fig.update_xaxes(visible=False)
fig.update_yaxes(visible=False, autorange="reversed", scaleanchor="x")
fig.update_annotations(font=dict(size=13))
style(fig, height=380)
fig.show()

Figure 2: The ten final centers, displayed as $8\times8$ images.

Most centers are recognizably an “average” handwritten digit, so k-means rediscovers much of the digit structure with no labels at all. It is not perfect. With this initialization the 1s are split between c5 and c9, and c9 is about half 1s and half 8s, which is why its center looks like a thin blend of the two. Most of the remaining 8s land with the 9s in c3, so no cluster belongs to the 8s. Rerunning from a different seed (or with k-means++ and several restarts, as sklearn.cluster.KMeans does) can give a different local minimum.

Clustering the shadow instead

scikit-learn’s own k-means digits example does things in the other order: it projects to the PCA plane first and then runs k-means on the $2$-dimensional scores $Z$.

from sklearn.cluster import KMeans

km_2d = KMeans(n_clusters=k, n_init=4, random_state=0).fit(Z)
frac_2d = majority_fraction(km_2d.labels_)

print(f"k-means in R^64, viewed in 2D:  {frac_64:.1%} of images in a cluster whose majority digit is their own")
print(f"k-means on the 2D PCA scores:  {frac_2d:.1%}")

k-means in R^64, viewed in 2D:  79.2% of images in a cluster whose majority digit is their own
k-means on the 2D PCA scores:  59.4%

Plot the 2D clusters and their cells

grid_x, grid_y = np.meshgrid(np.linspace(Z[:, 0].min() - 1, Z[:, 0].max() + 1, 400),
                             np.linspace(Z[:, 1].min() - 1, Z[:, 1].max() + 1, 400))
cells = km_2d.predict(np.c_[grid_x.ravel(), grid_y.ravel()]).reshape(grid_x.shape)
discrete = [[edge, c] for l, c in enumerate(CLUSTER_COLOR) for edge in (l / k, (l + 1) / k)]
C2 = km_2d.cluster_centers_
center_digit = [np.bincount(y[km_2d.labels_ == l]).argmax() for l in range(k)]

fig = go.Figure()
fig.add_trace(go.Heatmap(x=grid_x[0], y=grid_y[:, 0], z=cells, colorscale=discrete,
                         zmin=-0.5, zmax=k - 0.5, opacity=0.3, showscale=False,
                         hoverinfo="skip"))
fig.add_trace(go.Scatter(x=Z[:, 0], y=Z[:, 1], mode="markers", showlegend=False,
                         marker=dict(color=INK, size=4, opacity=0.4),
                         customdata=y, hovertemplate="true digit %{customdata}<extra></extra>"))
fig.add_trace(go.Scatter(x=C2[:, 0], y=C2[:, 1], mode="markers+text", showlegend=False,
                         text=[f"mostly {dg}" for dg in center_digit], textposition="top center",
                         textfont=LABEL_FONT,
                         marker=dict(symbol="x", size=14, color=CLUSTER_COLOR,
                                     line=dict(color="white", width=1.5)),
                         hovertemplate="%{text}<extra></extra>"))
fig.update_xaxes(title="z1 (PC1)", constrain="domain", range=[grid_x.min(), grid_x.max()])
fig.update_yaxes(title="z2 (PC2)", scaleanchor="x", scaleratio=1,
                 range=[grid_y.min(), grid_y.max()])
style(fig, height=600, title=f"k-means on the 2D PCA scores: {frac_2d:.1%} correct")
fig.show()

Figure 3: k-means run directly on the top-2 PCA scores. Shaded cells are the nearest-center regions in the plane; each cross is labeled with its cluster’s most common digit.

The picture is tidier, since in the plane the clusters are the convex cells, but the clustering is worse: about $20$ percentage points fewer images share their cluster’s majority digit. Projecting first throws away the roughly $71\%$ of the variance that PCA leaves out, and with it much of what distinguishes one digit from another. A clean $2$D picture is a statement about the shadow, not about the data.

--- title: "k-Means on Handwritten Digits" subtitle: "Math 123: Mathematical Aspects of Data Analysis — Fall 2026" date: 2026-10-13 number-sections: false jupyter: python3 format: html: toc: true toc-depth: 2 code-line-numbers: true code-tools: true out-width: "100%" column: page execute: echo: true draft: false --- The scikit-learn [digits dataset](https://scikit-learn.org/stable/modules/generated/sklearn.datasets.load_digits.html) is $1797$ handwritten digits, each an $8\times8$ grayscale image with pixel values $0,\ldots,16$. Flattening each image gives a point $\mathbf{x}_i\in\mathbb{R}^{64}$. We run [Lloyd's algorithm](../../lectures/lecture-12.qmd#alg-lloyd) for the [k-means problem](../../lectures/lecture-12.qmd#def-kmeans) with $k=10$ on these points, without ever looking at the true labels, and use PCA only to *watch* what happens. ```{python} import numpy as np import plotly.graph_objects as go from plotly.subplots import make_subplots from sklearn.datasets import load_digits from sklearn.decomposition import PCA digits = load_digits() X = digits.data # (1797, 64), one flattened image per row y = digits.target # true digit, used only to evaluate at the end n, d = X.shape k = 10 print(f"{n} images, {d} pixels each") ``` ```{python} #| code-fold: true #| code-summary: "Figure style" # A muted ink colour and a translucent grid read correctly against both the # light and the dark version of this site, so the figures need no theme switch. INK = "#7c7c78" GRID = "rgba(128, 128, 128, 0.28)" # Dark bold text with a white halo stays legible over the points in both themes. LABEL_FONT = dict(color="#1a1a19", size=13, weight="bold", shadow="0 0 3px white, 0 0 3px white, 0 0 2px white") CLUSTER_COLOR = [ "#2a78d6", "#eb6834", "#1baf7a", "#eda100", "#e87ba4", "#008300", "#6250d6", "#e34948", "#0f8b8d", "#8c564b", ] def style(fig, height=420, title=None, showlegend=None): """Apply the shared layout to a plot with axes.""" fig.update_layout( template="simple_white", paper_bgcolor="rgba(0,0,0,0)", plot_bgcolor="rgba(0,0,0,0)", font=dict(color=INK, size=13), title=dict(text=title, x=0, xanchor="left", font=dict(size=15)), margin=dict(l=70, r=30, t=60 if title else 30, b=55), height=height, showlegend=showlegend, hovermode="closest", ) fig.update_xaxes(gridcolor=GRID, linecolor=GRID, zerolinecolor=GRID) fig.update_yaxes(gridcolor=GRID, linecolor=GRID, zerolinecolor=GRID) return fig ``` ## Lloyd's algorithm in 64 dimensions {#sec-lloyd} This is the `lloyd` function from the "k-means in practice" box after [Lloyd's algorithm](../../lectures/lecture-12.qmd#alg-lloyd), changed only to record the assignment, the centers, and the cost after every update step. We start from $k$ distinct data points chosen at random. ```{python} def lloyd_history(X, centers, max_iter=100): centers = np.array(centers, dtype=float) labels = None history = [] for _ in range(max_iter): dists = ((X[:, None, :] - centers[None, :, :]) ** 2).sum(axis=2) # (n, k) new_labels = dists.argmin(axis=1) # assignment step if labels is not None and np.array_equal(new_labels, labels): break labels = new_labels for l in range(len(centers)): # update step if np.any(labels == l): centers[l] = X[labels == l].mean(axis=0) cost = ((X - centers[labels]) ** 2).sum() history.append((labels.copy(), centers.copy(), cost)) return history rng = np.random.default_rng(0) initial_centers = X[rng.choice(n, size=k, replace=False)] history = lloyd_history(X, initial_centers) for t, (_, _, cost) in enumerate(history, start=1): print(f"J_{t:<2d} = {cost:,.0f}") ``` The cost decreases at every step and the algorithm stops once the assignment repeats, exactly as [the termination theorem](../../lectures/lecture-12.qmd#thm-lloyd) promises. ## Watching it in the PCA plane {#sec-watch} We cannot plot points in $\mathbb{R}^{64}$, so we project onto the top $2$ [principal directions](../../lectures/lecture-06.qmd#def-pca), fit once on the data, and send the centers through the same projection. Since projection is linear (an affine map once centering is included), the projected center of a cluster is exactly the mean of its projected points, so each cross sits at the center of mass of its color in the plot. ```{python} pca = PCA(n_components=2).fit(X) Z = pca.transform(X) # (1797, 2) print("explained variance ratio:", np.round(pca.explained_variance_ratio_, 3)) ``` ```{python} #| label: fig-lloyd-steps #| fig-cap: "Lloyd's algorithm on the digits, viewed in the top-2 principal subspace. Drag the slider to step through the iterations; step 0 shows the random initial centers." #| fig-alt: "A scatter plot of the 1797 digit images in their first two principal-component coordinates, with ten cluster centers drawn as labeled crosses. A slider steps through the iterations of Lloyd's algorithm: at step 0 the points are gray and the centers are ten random images; at later steps each point is colored by its current cluster and the centers move toward the middle of their clusters." #| code-fold: true #| code-summary: "Plot each iteration, with a slider" frames = [(None, initial_centers, None)] + history fig = go.Figure() for t, (labels, centers, cost) in enumerate(frames): C = pca.transform(centers) colors = INK if labels is None else [CLUSTER_COLOR[l] for l in labels] fig.add_trace(go.Scatter( x=Z[:, 0], y=Z[:, 1], mode="markers", visible=(t == 0), showlegend=False, marker=dict(color=colors, size=5, opacity=0.35), customdata=np.c_[y, labels if labels is not None else np.full(n, -1)], hovertemplate="true digit %{customdata[0]}<br>cluster %{customdata[1]}<extra></extra>", )) fig.add_trace(go.Scatter( x=C[:, 0], y=C[:, 1], mode="markers+text", visible=(t == 0), showlegend=False, text=[f"c{l}" for l in range(k)], textposition="top center", textfont=LABEL_FONT, marker=dict(symbol="x", size=14, color=CLUSTER_COLOR, line=dict(color="white", width=1.5)), hovertemplate="center %{text}<extra></extra>", )) steps = [] for t, (_, _, cost) in enumerate(frames): visible = [False] * len(fig.data) visible[2 * t] = visible[2 * t + 1] = True title = ("Step 0: random initial centers" if cost is None else f"Step {t}: cost J = {cost:,.0f}") steps.append(dict(method="update", label=str(t), args=[{"visible": visible}, {"title.text": title}])) fig.update_layout(sliders=[dict(active=0, steps=steps, pad=dict(t=50), currentvalue=dict(prefix="step "))]) fig.update_xaxes(title="z1 (PC1)", constrain="domain") # shrink the plot, not pad the range fig.update_yaxes(title="z2 (PC2)", scaleanchor="x", scaleratio=1) # one unit is the same length on both axes style(fig, height=600, title="Step 0: random initial centers") fig.show() ``` Two things are worth noticing. - Most of the motion happens in the first few steps; after that, only points near the boundaries between clusters change color. - The colors overlap in the plot even though k-means clusters are [convex](../../lectures/lecture-12.qmd#prp-voronoi-convex). The nearest-center regions live in $\mathbb{R}^{64}$, and the top two principal directions keep only about $29\%$ of the variance, so clusters separated in the full space land on top of each other in this shadow. ## What the clusters found {#sec-centers} Each center is the mean of a set of images, so it is itself an $8\times8$ image. Only now do we use the true labels, to see which digit dominates each cluster. ```{python} labels, centers, cost = history[-1] for l in range(k): counts = np.bincount(y[labels == l], minlength=10) print(f"c{l}: {counts.sum():4d} images, {counts.max() / counts.sum():.0%} are the digit {counts.argmax()}") def majority_fraction(labels): """Fraction of images whose cluster's most common digit is their own.""" return sum(np.bincount(y[labels == l]).max() for l in np.unique(labels)) / n frac_64 = majority_fraction(labels) print(f"\nfraction of images whose cluster's majority digit is their own: {frac_64:.1%}") ``` ```{python} #| label: fig-centers #| fig-cap: "The ten final centers, displayed as $8\\times8$ images." #| fig-alt: "Ten small grayscale images, the final k-means centers. Most are blurry but recognizable digits: 3, 5, 7, 9, 4, 1, 2, 0, and 6. The last one is a thin blend that could be read as a 1 or an 8." #| code-fold: true #| code-summary: "Plot the centers as images" fig = make_subplots(rows=2, cols=5, horizontal_spacing=0.02, vertical_spacing=0.12, subplot_titles=[f"c{l}" for l in range(k)]) for l in range(k): fig.add_trace(go.Heatmap(z=centers[l].reshape(8, 8), colorscale="Greys", zmin=0, zmax=16, showscale=False, hovertemplate="%{z:.1f}<extra></extra>"), row=l // 5 + 1, col=l % 5 + 1) fig.update_xaxes(visible=False) fig.update_yaxes(visible=False, autorange="reversed", scaleanchor="x") fig.update_annotations(font=dict(size=13)) style(fig, height=380) fig.show() ``` Most centers are recognizably an "average" handwritten digit, so k-means rediscovers much of the digit structure with no labels at all. It is not perfect. With this initialization the 1s are split between `c5` and `c9`, and `c9` is about half 1s and half 8s, which is why its center looks like a thin blend of the two. Most of the remaining 8s land with the 9s in `c3`, so no cluster belongs to the 8s. Rerunning from a different seed (or with k-means++ and several restarts, as `sklearn.cluster.KMeans` does) can give a different local minimum. ## Clustering the shadow instead {#sec-shadow} scikit-learn's own [k-means digits example](https://scikit-learn.org/stable/auto_examples/cluster/plot_kmeans_digits.html) does things in the other order: it projects to the PCA plane *first* and then runs k-means on the $2$-dimensional scores $Z$. ```{python} from sklearn.cluster import KMeans km_2d = KMeans(n_clusters=k, n_init=4, random_state=0).fit(Z) frac_2d = majority_fraction(km_2d.labels_) print(f"k-means in R^64, viewed in 2D: {frac_64:.1%} of images in a cluster whose majority digit is their own") print(f"k-means on the 2D PCA scores: {frac_2d:.1%}") ``` ```{python} #| label: fig-shadow #| fig-cap: "k-means run directly on the top-2 PCA scores. Shaded cells are the nearest-center regions in the plane; each cross is labeled with its cluster's most common digit." #| fig-alt: "The same scatter plot of the digits in the PCA plane, now tiled by ten convex shaded cells, one per cluster center. The cells partition the plane cleanly, but two cells are labeled with the digit 8 and no cell is labeled 5." #| code-fold: true #| code-summary: "Plot the 2D clusters and their cells" grid_x, grid_y = np.meshgrid(np.linspace(Z[:, 0].min() - 1, Z[:, 0].max() + 1, 400), np.linspace(Z[:, 1].min() - 1, Z[:, 1].max() + 1, 400)) cells = km_2d.predict(np.c_[grid_x.ravel(), grid_y.ravel()]).reshape(grid_x.shape) discrete = [[edge, c] for l, c in enumerate(CLUSTER_COLOR) for edge in (l / k, (l + 1) / k)] C2 = km_2d.cluster_centers_ center_digit = [np.bincount(y[km_2d.labels_ == l]).argmax() for l in range(k)] fig = go.Figure() fig.add_trace(go.Heatmap(x=grid_x[0], y=grid_y[:, 0], z=cells, colorscale=discrete, zmin=-0.5, zmax=k - 0.5, opacity=0.3, showscale=False, hoverinfo="skip")) fig.add_trace(go.Scatter(x=Z[:, 0], y=Z[:, 1], mode="markers", showlegend=False, marker=dict(color=INK, size=4, opacity=0.4), customdata=y, hovertemplate="true digit %{customdata}<extra></extra>")) fig.add_trace(go.Scatter(x=C2[:, 0], y=C2[:, 1], mode="markers+text", showlegend=False, text=[f"mostly {dg}" for dg in center_digit], textposition="top center", textfont=LABEL_FONT, marker=dict(symbol="x", size=14, color=CLUSTER_COLOR, line=dict(color="white", width=1.5)), hovertemplate="%{text}<extra></extra>")) fig.update_xaxes(title="z1 (PC1)", constrain="domain", range=[grid_x.min(), grid_x.max()]) fig.update_yaxes(title="z2 (PC2)", scaleanchor="x", scaleratio=1, range=[grid_y.min(), grid_y.max()]) style(fig, height=600, title=f"k-means on the 2D PCA scores: {frac_2d:.1%} correct") fig.show() ``` The picture is tidier, since in the plane the clusters *are* the convex cells, but the clustering is worse: about $20$ percentage points fewer images share their cluster's majority digit. Projecting first throws away the roughly $71\%$ of the variance that PCA leaves out, and with it much of what distinguishes one digit from another. A clean $2$D picture is a statement about the shadow, not about the data.