k-Means vs. Spectral Clustering on Two Circles

Math 123: Mathematical Aspects of Data Analysis — Fall 2026

Published

October 20, 2026

Two concentric rings are the standard example of clusters that k-means cannot find: its clusters are always convex, and no straight line separates one ring from the other. Spectral clustering builds a similarity graph from the points instead, and on that graph the two rings come apart. This page runs both methods through scikit-learn on the data from Lecture 14 and compares the results.

import warnings

import numpy as np
import plotly.graph_objects as go
from plotly.subplots import make_subplots
from scipy.sparse.csgraph import connected_components
from sklearn.cluster import KMeans, SpectralClustering
from sklearn.datasets import make_circles
from sklearn.manifold import spectral_embedding
from sklearn.metrics import adjusted_rand_score

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)"
COLOR = ["#2a78d6", "#eb6834"]          # outer ring / cluster 0, inner ring / cluster 1
RING_NAME = ["outer ring", "inner ring"]


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


def equal_aspect(fig, n_panels):
    """Make one unit the same length on both axes of every subplot."""
    for i in range(1, n_panels + 1):
        suffix = "" if i == 1 else str(i)
        fig.layout[f"yaxis{suffix}"].scaleanchor = f"x{suffix}"
        fig.layout[f"xaxis{suffix}"].constrain = "domain"

The data

make_circles puts half the points on a circle of radius $1$ and half on a circle of radius factor, then adds Gaussian noise. It also returns the true ring of each point: y = 0 for the outer ring and y = 1 for the inner one. Neither method gets to see y; we only use it to score the results.

X, y = make_circles(n_samples=400, factor=0.4, noise=0.05, random_state=0)
print("X:", X.shape, " points per ring:", np.bincount(y))

X: (400, 2)  points per ring: [200 200]

To score a clustering we use two numbers. The agreement is the fraction of points whose cluster matches their ring, after renaming the clusters to fit the rings as well as possible. The adjusted Rand index (ARI) from sklearn.metrics is $1$ for a perfect match and about $0$ for a random labelling, whatever the names of the clusters.

def align(labels, y):
    """Rename the two clusters so that cluster r matches ring r as well as possible."""
    return labels if np.mean(labels == y) >= 0.5 else 1 - labels


def agreement(labels, y):
    return np.mean(align(labels, y) == y)


def ari(y, labels):
    """Adjusted Rand index, rounded to two digits for figure titles (no "-0.00")."""
    return round(adjusted_rand_score(y, labels), 2) + 0.0

k-means

Every scikit-learn estimator works the same way: set the parameters in the constructor, call .fit(X), and read the results from the attributes that end in an underscore. For KMeans, n_init=10 runs Lloyd’s algorithm from ten k-means++ seedings and keeps the run with the lowest cost.

kmeans = KMeans(n_clusters=2, n_init=10, random_state=0)
kmeans.fit(X)

print("labels_          ", kmeans.labels_[:12], "...")
print("cluster_centers_ ", kmeans.cluster_centers_.round(3).tolist())
print("inertia_ (cost)  ", round(kmeans.inertia_, 3))
print("n_iter_          ", kmeans.n_iter_)
print(f"agreement {agreement(kmeans.labels_, y):.1%},  ARI {adjusted_rand_score(y, kmeans.labels_):.3f}")

labels_           [0 1 0 0 1 1 0 1 0 0 1 0] ...
cluster_centers_  [[0.27, -0.339], [-0.299, 0.358]]
inertia_ (cost)   153.458
n_iter_           2
agreement 51.2%,  ARI -0.002

A fitted KMeans also has .predict, which assigns new points to the nearest center. This is how we shade its two cells in the figure below.

print(kmeans.predict([[0.0, 0.0], [0.9, 0.0], [-0.9, 0.0]]))

[0 0 1]

The agreement is about $50\%$ and the ARI is about $0$: no better than a coin flip. Both centers land near the middle, and the cut between them goes straight through both rings. More restarts do not help, because the best convex split is still a straight cut.

Spectral clustering

SpectralClustering has the same interface. With affinity="rbf" it builds the Gaussian similarity graph with weights $w_{ij}=\exp(-\gamma\|\bfx_i-\bfx_j\|_2^2)$, so the lecture’s bandwidth $\varepsilon$ corresponds to gamma $=1/(2\varepsilon)$. We use $\varepsilon=0.01$, as in Lecture 13.

eps = 0.01
spectral = SpectralClustering(n_clusters=2, affinity="rbf", gamma=1 / (2 * eps), random_state=0)
spectral.fit(X)

print("labels_          ", spectral.labels_[:12], "...")
print("affinity_matrix_ ", spectral.affinity_matrix_.shape, "(the weight matrix W)")
print(f"agreement {agreement(spectral.labels_, y):.1%},  ARI {adjusted_rand_score(y, spectral.labels_):.3f}")

labels_           [1 0 0 1 0 1 1 0 1 0 1 1] ...
affinity_matrix_  (400, 400) (the weight matrix W)
agreement 100.0%,  ARI 1.000

There is no cluster_centers_ and no .predict. Spectral clustering labels the vertices of one particular graph, so a new point would mean a new graph. For the same reason, fit_predict is the usual way to call it.

labels_nn = SpectralClustering(n_clusters=2, affinity="nearest_neighbors", n_neighbors=10,
                               random_state=0).fit_predict(X)
print(f"nearest-neighbour graph: agreement {agreement(labels_nn, y):.1%}")

nearest-neighbour graph: agreement 100.0%

/home/lotus/Repositories/tufts-math123-fall2026/.venv/lib/python3.13/site-packages/sklearn/manifold/_spectral_embedding.py:324: UserWarning: Graph is not fully connected, spectral embedding may not work as expected.
  warnings.warn(

That run, with the nearest-neighbour graph instead of the Gaussian one, prints the warning “Graph is not fully connected”. In this case the warning is good news. With $10$ neighbours, no edge crosses between the rings, so the graph has exactly two connected components, and those components are the two rings.

with warnings.catch_warnings():
    warnings.simplefilter("ignore", UserWarning)
    W_nn = SpectralClustering(n_clusters=2, affinity="nearest_neighbors", n_neighbors=10,
                              random_state=0).fit(X).affinity_matrix_
n_comp, comp = connected_components(W_nn)
print(f"{n_comp} components; points per ring in each component:",
      [np.bincount(y[comp == c], minlength=2).tolist() for c in range(n_comp)])

2 components; points per ring in each component: [[200, 0], [0, 200]]

Side by side

Plot the comparison

results = {"k-means": align(kmeans.labels_, y), "spectral clustering": align(spectral.labels_, y)}
lim = 1.45
gx, gy = np.meshgrid(np.linspace(-lim, lim, 200), np.linspace(-lim, lim, 200))
centers = kmeans.cluster_centers_
cells = kmeans.predict(np.c_[gx.ravel(), gy.ravel()])
if not np.array_equal(results["k-means"], kmeans.labels_):   # align swapped the names
    centers, cells = centers[::-1], 1 - cells
two_colors = [[0, COLOR[0]], [0.5, COLOR[0]], [0.5, COLOR[1]], [1, COLOR[1]]]

titles = [f"<b>{name}</b>: agreement {np.mean(lab == y):.0%}, ARI {ari(y, lab):.2f}"
          for name, lab in results.items()]
fig = make_subplots(rows=1, cols=2, horizontal_spacing=0.06, subplot_titles=titles)

# Shaded nearest-center cells and the cut between the two centers (k-means panel only).
fig.add_trace(go.Heatmap(x=gx[0], y=gy[:, 0], z=cells.reshape(gx.shape), colorscale=two_colors,
                         zmin=0, zmax=1, opacity=0.13, showscale=False, hoverinfo="skip"),
              row=1, col=1)
# The cut is the perpendicular bisector: points x with (mu1 - mu0) . x = (|mu1|^2 - |mu0|^2) / 2.
normal = centers[1] - centers[0]
midpoint = centers.mean(axis=0)
direction = np.array([-normal[1], normal[0]]) / np.linalg.norm(normal)
cut = midpoint + np.outer([-3, 3], direction)
fig.add_trace(go.Scatter(x=cut[:, 0], y=cut[:, 1], mode="lines", showlegend=False,
                         line=dict(color=INK, width=2, dash="dash"), hoverinfo="skip"),
              row=1, col=1)

# Points: for each panel, one trace per cluster and one trace per ring.
# The buttons switch between the two sets of traces.
point_traces = {"cluster": [], "ring": []}
for col, (name, lab) in enumerate(results.items(), start=1):
    hover = [f"{RING_NAME[r]}, cluster {c}" for r, c in zip(y, lab)]
    for mode, groups in [("cluster", lab), ("ring", y)]:
        for g in (0, 1):
            mask = groups == g
            legend_name = f"cluster {g}" if mode == "cluster" else RING_NAME[g]
            fig.add_trace(go.Scatter(
                x=X[mask, 0], y=X[mask, 1], mode="markers",
                name=legend_name, legendgroup=f"{mode}{g}", showlegend=(col == 1),
                visible=(mode == "cluster"),
                marker=dict(size=7, color=COLOR[g], line=dict(color="white", width=0.8)),
                text=np.array(hover)[mask],
                hovertemplate="%{text}<br>(%{x:.2f}, %{y:.2f})<extra></extra>"),
                row=1, col=col)
            point_traces[mode].append(len(fig.data) - 1)

fig.add_trace(go.Scatter(x=centers[:, 0], y=centers[:, 1], mode="markers", name="k-means centers",
                         marker=dict(symbol="x", size=15, color=COLOR,
                                     line=dict(color="white", width=1.5)),
                         hovertemplate="center (%{x:.2f}, %{y:.2f})<extra></extra>"),
              row=1, col=1)


def visibility(mode):
    vis = [True] * len(fig.data)
    for m, idx in point_traces.items():
        for i in idx:
            vis[i] = (m == mode)
    return vis


fig.update_layout(updatemenus=[dict(
    type="buttons", direction="right", x=0, xanchor="left", y=1.2, yanchor="bottom",
    showactive=True, bgcolor="rgba(128, 128, 128, 0.12)", bordercolor=GRID, font=dict(color=INK),
    buttons=[dict(label="Color by cluster found", method="restyle",
                  args=[{"visible": visibility("cluster")}]),
             dict(label="Color by true ring", method="restyle",
                  args=[{"visible": visibility("ring")}])])])
fig.update_xaxes(range=[-lim, lim], title="first coordinate")
fig.update_yaxes(range=[-lim, lim])
fig.update_yaxes(title="second coordinate", col=1)
equal_aspect(fig, 2)
style(fig, height=520, showlegend=True)
fig.update_layout(margin=dict(t=110, b=100),
                  legend=dict(orientation="h", x=0, y=-0.16, xanchor="left", yanchor="top"))
fig.show()

Figure 1: k-means (left) and spectral clustering with $\varepsilon=0.01$ (right) on the same $400$ points. Colors are the clusters found, renamed to match the rings where possible. On the left, the shading shows the two nearest-center cells, crosses are the centers, and the dashed line is the straight cut between them. Use the buttons to recolor both panels by the true rings, and hover over a point to see its ring and cluster.

Switching the buttons leaves the right panel unchanged, because spectral clustering found the rings exactly. On the left, the cut goes through both rings, so each k-means cluster holds half of the outer ring and half of the inner ring.

Why the graph helps

Spectral clustering runs k-means too, but on the spectral embedding rather than on the raw coordinates. sklearn.manifold.spectral_embedding computes that embedding from the weight matrix that the fit stored. With $k=2$ the embedding is one number $\phi_2(i)$ per point; below, we plot it against the angle of the point around the origin.

phi2 = spectral_embedding(spectral.affinity_matrix_, n_components=1, drop_first=True,
                          random_state=0)[:, 0]
angle = np.degrees(np.arctan2(X[:, 1], X[:, 0]))
for r in (0, 1):
    print(f"{RING_NAME[r]}: phi_2 = {phi2[y == r].mean():+.5f} ± {phi2[y == r].std():.1e}")

outer ring: phi_2 = +0.02474 ± 4.9e-06
inner ring: phi_2 = -0.00915 ± 3.5e-07

Plot the embedding

fig = go.Figure()
for r in (0, 1):
    m = y == r
    fig.add_trace(go.Scatter(x=angle[m], y=phi2[m], mode="markers", name=RING_NAME[r],
                             marker=dict(size=7, color=COLOR[r], line=dict(color="white", width=0.8)),
                             hovertemplate="angle %{x:.0f}°<br>φ₂ = %{y:.5f}<extra>"
                                           + RING_NAME[r] + "</extra>"))
    fig.add_annotation(x=180, y=phi2[m].mean(), text=RING_NAME[r], showarrow=False,
                       xanchor="right", yanchor="bottom", yshift=6, font=dict(color=INK, size=13))
pad = 0.15 * np.ptp(phi2)
fig.update_xaxes(title="angle of the point around the origin (degrees)", range=[-185, 185],
                 tickvals=[-180, -90, 0, 90, 180])
fig.update_yaxes(title="φ₂", range=[phi2.min() - pad, phi2.max() + pad])
style(fig, height=360, showlegend=True)
fig.update_layout(legend=dict(x=0.01, y=0.5, xanchor="left", yanchor="middle"))
fig.show()

Figure 2: The one-dimensional spectral embedding $\phi_2$ for $\varepsilon=0.01$, plotted against the angle of each point around the origin. Each ring collapses to almost a single value of $\phi_2$, so k-means in the embedding only has to separate two points.

With $\varepsilon=0.01$, the edges between the rings are so weak that the graph is nearly disconnected. The embedding then looks like the ideal case, where each component collapses to a single point. The rings, which are not convex in $\Rbb^2$, become two points in the embedding, and two points are easy for k-means to separate.

The bandwidth matters

The graph is the real modelling choice, as Lecture 13 notes (see also the notes on choosing hyperparameters). Refit spectral clustering for a range of bandwidths $\varepsilon$ and drag the slider to see the clusters change.

eps_values = [0.003, 0.01, 0.02, 0.03, 0.04, 0.05, 0.07, 0.1]
sweep = {e: align(SpectralClustering(n_clusters=2, affinity="rbf", gamma=1 / (2 * e),
                                     random_state=0).fit_predict(X), y)
         for e in eps_values}
for e, lab in sweep.items():
    print(f"eps = {e:<5}  agreement {np.mean(lab == y):6.1%}   ARI {adjusted_rand_score(y, lab):.3f}")

eps = 0.003  agreement 100.0%   ARI 1.000
eps = 0.01   agreement 100.0%   ARI 1.000
eps = 0.02   agreement 100.0%   ARI 1.000
eps = 0.03   agreement 100.0%   ARI 1.000
eps = 0.04   agreement 100.0%   ARI 1.000
eps = 0.05   agreement  70.8%   ARI 0.171
eps = 0.07   agreement  54.5%   ARI 0.006
eps = 0.1    agreement  51.7%   ARI -0.001

Plot the bandwidth slider

start = eps_values.index(0.01)


def colors(lab):
    return [COLOR[c] for c in lab]


def hover(lab):
    return [f"{RING_NAME[r]}, cluster {c}" for r, c in zip(y, lab)]


def title(e):
    lab = sweep[e]
    return (f"ε = {e}  (gamma = {1 / (2 * e):g}): agreement {np.mean(lab == y):.0%}, "
            f"ARI {ari(y, lab):.2f}")


fig = go.Figure(go.Scatter(
    x=X[:, 0], y=X[:, 1], mode="markers", showlegend=False,
    marker=dict(size=7, color=colors(sweep[eps_values[start]]), line=dict(color="white", width=0.8)),
    text=hover(sweep[eps_values[start]]),
    hovertemplate="%{text}<br>(%{x:.2f}, %{y:.2f})<extra></extra>"))
for g in (0, 1):                                   # legend entries for the two colors
    fig.add_trace(go.Scatter(x=[None], y=[None], mode="markers", name=f"cluster {g}",
                             marker=dict(size=9, color=COLOR[g])))
steps = [dict(label=f"{e}", method="update",
              args=[{"marker.color": [colors(sweep[e])], "text": [hover(sweep[e])]},
                    {"title.text": title(e)}, [0]])
         for e in eps_values]
fig.update_layout(sliders=[dict(active=start, steps=steps, x=0, len=1, pad=dict(t=50),
                                currentvalue=dict(prefix="bandwidth ε = ", font=dict(color=INK)),
                                font=dict(color=INK), bordercolor=GRID)])
lim = 1.45
fig.update_xaxes(range=[-lim, lim], title="first coordinate", constrain="domain")
fig.update_yaxes(range=[-lim, lim], title="second coordinate", scaleanchor="x")
style(fig, height=600, title=title(eps_values[start]), showlegend=True)
fig.update_layout(legend=dict(x=1, y=1, xanchor="right", yanchor="top"))
fig.show()

Figure 3: Spectral clustering with the Gaussian graph for several bandwidths $\varepsilon$; drag the slider to change $\varepsilon$. For small $\varepsilon$ the clusters are the two rings. As $\varepsilon$ grows, the edges between the rings get stronger and the clustering falls apart.

The rings are recovered for every $\varepsilon$ up to about $0.04$. Past that, the Gaussian kernel is wide enough to put real weight on edges between the rings, so there is no longer a bottleneck between them for the embedding to find. k-means fails for a different reason: no choice of parameters gives it a cluster shaped like a ring.

--- title: "k-Means vs. Spectral Clustering on Two Circles" subtitle: "Math 123: Mathematical Aspects of Data Analysis — Fall 2026" date: 2026-10-20 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 --- Two concentric rings are the standard example of clusters that [k-means](../../lectures/lecture-12.qmd#def-kmeans) cannot find: its clusters are always [convex](../../lectures/lecture-12.qmd#prp-voronoi-convex), and no straight line separates one ring from the other. [Spectral clustering](../../lectures/lecture-14.qmd#alg-spectral-clustering) builds a [similarity graph](../../lectures/lecture-13.qmd#from-point-clouds-to-graphs) from the points instead, and on that graph the two rings come apart. This page runs both methods through scikit-learn on the data from [Lecture 14](../../lectures/lecture-14.qmd#alg-spectral-clustering) and compares the results. ```{python} import warnings import numpy as np import plotly.graph_objects as go from plotly.subplots import make_subplots from scipy.sparse.csgraph import connected_components from sklearn.cluster import KMeans, SpectralClustering from sklearn.datasets import make_circles from sklearn.manifold import spectral_embedding from sklearn.metrics import adjusted_rand_score ``` ```{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)" COLOR = ["#2a78d6", "#eb6834"] # outer ring / cluster 0, inner ring / cluster 1 RING_NAME = ["outer ring", "inner ring"] 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 def equal_aspect(fig, n_panels): """Make one unit the same length on both axes of every subplot.""" for i in range(1, n_panels + 1): suffix = "" if i == 1 else str(i) fig.layout[f"yaxis{suffix}"].scaleanchor = f"x{suffix}" fig.layout[f"xaxis{suffix}"].constrain = "domain" ``` ## The data {#sec-data} [`make_circles`](https://scikit-learn.org/stable/modules/generated/sklearn.datasets.make_circles.html) puts half the points on a circle of radius $1$ and half on a circle of radius `factor`, then adds Gaussian noise. It also returns the true ring of each point: `y = 0` for the outer ring and `y = 1` for the inner one. Neither method gets to see `y`; we only use it to score the results. ```{python} X, y = make_circles(n_samples=400, factor=0.4, noise=0.05, random_state=0) print("X:", X.shape, " points per ring:", np.bincount(y)) ``` To score a clustering we use two numbers. The *agreement* is the fraction of points whose cluster matches their ring, after renaming the clusters to fit the rings as well as possible. The [adjusted Rand index](https://scikit-learn.org/stable/modules/generated/sklearn.metrics.adjusted_rand_score.html) (ARI) from `sklearn.metrics` is $1$ for a perfect match and about $0$ for a random labelling, whatever the names of the clusters. ```{python} def align(labels, y): """Rename the two clusters so that cluster r matches ring r as well as possible.""" return labels if np.mean(labels == y) >= 0.5 else 1 - labels def agreement(labels, y): return np.mean(align(labels, y) == y) def ari(y, labels): """Adjusted Rand index, rounded to two digits for figure titles (no "-0.00").""" return round(adjusted_rand_score(y, labels), 2) + 0.0 ``` ## k-means {#sec-kmeans} Every scikit-learn estimator works the same way: set the parameters in the constructor, call `.fit(X)`, and read the results from the attributes that end in an underscore. For [`KMeans`](https://scikit-learn.org/stable/modules/generated/sklearn.cluster.KMeans.html), `n_init=10` runs [Lloyd's algorithm](../../lectures/lecture-12.qmd#alg-lloyd) from ten [k-means++](../../lectures/lecture-12.qmd#alg-kmeans-pp) seedings and keeps the run with the lowest cost. ```{python} kmeans = KMeans(n_clusters=2, n_init=10, random_state=0) kmeans.fit(X) print("labels_ ", kmeans.labels_[:12], "...") print("cluster_centers_ ", kmeans.cluster_centers_.round(3).tolist()) print("inertia_ (cost) ", round(kmeans.inertia_, 3)) print("n_iter_ ", kmeans.n_iter_) print(f"agreement {agreement(kmeans.labels_, y):.1%}, ARI {adjusted_rand_score(y, kmeans.labels_):.3f}") ``` A fitted `KMeans` also has `.predict`, which assigns *new* points to the nearest center. This is how we shade its two cells in the figure below. ```{python} print(kmeans.predict([[0.0, 0.0], [0.9, 0.0], [-0.9, 0.0]])) ``` The agreement is about $50\%$ and the ARI is about $0$: no better than a coin flip. Both centers land near the middle, and the cut between them goes straight through both rings. More restarts do not help, because the [best convex split](../../lectures/lecture-12.qmd#prp-voronoi-convex) is still a straight cut. ## Spectral clustering {#sec-spectral} [`SpectralClustering`](https://scikit-learn.org/stable/modules/generated/sklearn.cluster.SpectralClustering.html) has the same interface. With `affinity="rbf"` it builds the [Gaussian similarity graph](../../lectures/lecture-13.qmd#from-point-clouds-to-graphs) with weights $w_{ij}=\exp(-\gamma\|\bfx_i-\bfx_j\|_2^2)$, so the lecture's bandwidth $\varepsilon$ corresponds to `gamma` $=1/(2\varepsilon)$. We use $\varepsilon=0.01$, as in Lecture 13. ```{python} eps = 0.01 spectral = SpectralClustering(n_clusters=2, affinity="rbf", gamma=1 / (2 * eps), random_state=0) spectral.fit(X) print("labels_ ", spectral.labels_[:12], "...") print("affinity_matrix_ ", spectral.affinity_matrix_.shape, "(the weight matrix W)") print(f"agreement {agreement(spectral.labels_, y):.1%}, ARI {adjusted_rand_score(y, spectral.labels_):.3f}") ``` There is no `cluster_centers_` and no `.predict`. Spectral clustering labels the vertices of one particular graph, so a new point would mean a new graph. For the same reason, `fit_predict` is the usual way to call it. ```{python} labels_nn = SpectralClustering(n_clusters=2, affinity="nearest_neighbors", n_neighbors=10, random_state=0).fit_predict(X) print(f"nearest-neighbour graph: agreement {agreement(labels_nn, y):.1%}") ``` That run, with the [nearest-neighbour graph](../../lectures/lecture-13.qmd#from-point-clouds-to-graphs) instead of the Gaussian one, prints the warning "Graph is not fully connected". In this case the warning is good news. With $10$ neighbours, no edge crosses between the rings, so the graph has exactly two connected components, and those components are the two rings. ```{python} with warnings.catch_warnings(): warnings.simplefilter("ignore", UserWarning) W_nn = SpectralClustering(n_clusters=2, affinity="nearest_neighbors", n_neighbors=10, random_state=0).fit(X).affinity_matrix_ n_comp, comp = connected_components(W_nn) print(f"{n_comp} components; points per ring in each component:", [np.bincount(y[comp == c], minlength=2).tolist() for c in range(n_comp)]) ``` ## Side by side {#sec-compare} ```{python} #| label: fig-compare #| fig-cap: "k-means (left) and spectral clustering with $\\varepsilon=0.01$ (right) on the same $400$ points. Colors are the clusters found, renamed to match the rings where possible. On the left, the shading shows the two nearest-center cells, crosses are the centers, and the dashed line is the straight cut between them. Use the buttons to recolor both panels by the true rings, and hover over a point to see its ring and cluster." #| fig-alt: "Two scatter plots of the same two noisy concentric rings. On the left, k-means splits the picture with a straight diagonal dashed line through the middle; each half is one color and contains arcs of both rings, and the two centers sit just inside the inner ring on either side of the line. On the right, spectral clustering colors the whole outer ring blue and the whole inner ring orange." #| code-fold: true #| code-summary: "Plot the comparison" results = {"k-means": align(kmeans.labels_, y), "spectral clustering": align(spectral.labels_, y)} lim = 1.45 gx, gy = np.meshgrid(np.linspace(-lim, lim, 200), np.linspace(-lim, lim, 200)) centers = kmeans.cluster_centers_ cells = kmeans.predict(np.c_[gx.ravel(), gy.ravel()]) if not np.array_equal(results["k-means"], kmeans.labels_): # align swapped the names centers, cells = centers[::-1], 1 - cells two_colors = [[0, COLOR[0]], [0.5, COLOR[0]], [0.5, COLOR[1]], [1, COLOR[1]]] titles = [f"{name}: agreement {np.mean(lab == y):.0%}, ARI {ari(y, lab):.2f}" for name, lab in results.items()] fig = make_subplots(rows=1, cols=2, horizontal_spacing=0.06, subplot_titles=titles) # Shaded nearest-center cells and the cut between the two centers (k-means panel only). fig.add_trace(go.Heatmap(x=gx[0], y=gy[:, 0], z=cells.reshape(gx.shape), colorscale=two_colors, zmin=0, zmax=1, opacity=0.13, showscale=False, hoverinfo="skip"), row=1, col=1) # The cut is the perpendicular bisector: points x with (mu1 - mu0) . x = (|mu1|^2 - |mu0|^2) / 2. normal = centers[1] - centers[0] midpoint = centers.mean(axis=0) direction = np.array([-normal[1], normal[0]]) / np.linalg.norm(normal) cut = midpoint + np.outer([-3, 3], direction) fig.add_trace(go.Scatter(x=cut[:, 0], y=cut[:, 1], mode="lines", showlegend=False, line=dict(color=INK, width=2, dash="dash"), hoverinfo="skip"), row=1, col=1) # Points: for each panel, one trace per cluster and one trace per ring. # The buttons switch between the two sets of traces. point_traces = {"cluster": [], "ring": []} for col, (name, lab) in enumerate(results.items(), start=1): hover = [f"{RING_NAME[r]}, cluster {c}" for r, c in zip(y, lab)] for mode, groups in [("cluster", lab), ("ring", y)]: for g in (0, 1): mask = groups == g legend_name = f"cluster {g}" if mode == "cluster" else RING_NAME[g] fig.add_trace(go.Scatter( x=X[mask, 0], y=X[mask, 1], mode="markers", name=legend_name, legendgroup=f"{mode}{g}", showlegend=(col == 1), visible=(mode == "cluster"), marker=dict(size=7, color=COLOR[g], line=dict(color="white", width=0.8)), text=np.array(hover)[mask], hovertemplate="%{text} (%{x:.2f}, %{y:.2f})<extra></extra>"), row=1, col=col) point_traces[mode].append(len(fig.data) - 1) fig.add_trace(go.Scatter(x=centers[:, 0], y=centers[:, 1], mode="markers", name="k-means centers", marker=dict(symbol="x", size=15, color=COLOR, line=dict(color="white", width=1.5)), hovertemplate="center (%{x:.2f}, %{y:.2f})<extra></extra>"), row=1, col=1) def visibility(mode): vis = [True] * len(fig.data) for m, idx in point_traces.items(): for i in idx: vis[i] = (m == mode) return vis fig.update_layout(updatemenus=[dict( type="buttons", direction="right", x=0, xanchor="left", y=1.2, yanchor="bottom", showactive=True, bgcolor="rgba(128, 128, 128, 0.12)", bordercolor=GRID, font=dict(color=INK), buttons=[dict(label="Color by cluster found", method="restyle", args=[{"visible": visibility("cluster")}]), dict(label="Color by true ring", method="restyle", args=[{"visible": visibility("ring")}])])]) fig.update_xaxes(range=[-lim, lim], title="first coordinate") fig.update_yaxes(range=[-lim, lim]) fig.update_yaxes(title="second coordinate", col=1) equal_aspect(fig, 2) style(fig, height=520, showlegend=True) fig.update_layout(margin=dict(t=110, b=100), legend=dict(orientation="h", x=0, y=-0.16, xanchor="left", yanchor="top")) fig.show() ``` Switching the buttons leaves the right panel unchanged, because spectral clustering found the rings exactly. On the left, the cut goes through both rings, so each k-means cluster holds half of the outer ring and half of the inner ring. ## Why the graph helps {#sec-embedding} Spectral clustering runs k-means too, but on the [spectral embedding](../../lectures/lecture-14.qmd#def-spectral-embedding) rather than on the raw coordinates. `sklearn.manifold.spectral_embedding` computes that embedding from the weight matrix that the fit stored. With $k=2$ the embedding is one number $\phi_2(i)$ per point; below, we plot it against the angle of the point around the origin. ```{python} phi2 = spectral_embedding(spectral.affinity_matrix_, n_components=1, drop_first=True, random_state=0)[:, 0] angle = np.degrees(np.arctan2(X[:, 1], X[:, 0])) for r in (0, 1): print(f"{RING_NAME[r]}: phi_2 = {phi2[y == r].mean():+.5f} ± {phi2[y == r].std():.1e}") ``` ```{python} #| label: fig-embedding #| fig-cap: "The one-dimensional spectral embedding $\\phi_2$ for $\\varepsilon=0.01$, plotted against the angle of each point around the origin. Each ring collapses to almost a single value of $\\phi_2$, so k-means in the embedding only has to separate two points." #| fig-alt: "A scatter plot with angle from minus 180 to 180 degrees on the horizontal axis and phi 2 on the vertical axis. The 200 outer-ring points form a flat blue line at one height, and the 200 inner-ring points form a flat orange line at another height, with a wide empty gap between them." #| code-fold: true #| code-summary: "Plot the embedding" fig = go.Figure() for r in (0, 1): m = y == r fig.add_trace(go.Scatter(x=angle[m], y=phi2[m], mode="markers", name=RING_NAME[r], marker=dict(size=7, color=COLOR[r], line=dict(color="white", width=0.8)), hovertemplate="angle %{x:.0f}° φ₂ = %{y:.5f}<extra>" + RING_NAME[r] + "</extra>")) fig.add_annotation(x=180, y=phi2[m].mean(), text=RING_NAME[r], showarrow=False, xanchor="right", yanchor="bottom", yshift=6, font=dict(color=INK, size=13)) pad = 0.15 * np.ptp(phi2) fig.update_xaxes(title="angle of the point around the origin (degrees)", range=[-185, 185], tickvals=[-180, -90, 0, 90, 180]) fig.update_yaxes(title="φ₂", range=[phi2.min() - pad, phi2.max() + pad]) style(fig, height=360, showlegend=True) fig.update_layout(legend=dict(x=0.01, y=0.5, xanchor="left", yanchor="middle")) fig.show() ``` With $\varepsilon=0.01$, the edges between the rings are so weak that the graph is nearly disconnected. The embedding then looks like the [ideal case](../../lectures/lecture-14.qmd#prp-ideal-embedding), where each component collapses to a single point. The rings, which are not convex in $\Rbb^2$, become two points in the embedding, and two points are easy for k-means to separate. ## The bandwidth matters {#sec-bandwidth} The graph is the real modelling choice, as Lecture 13 notes (see also the notes on [choosing hyperparameters](../../notes/choosing-hyperparameters.qmd)). Refit spectral clustering for a range of bandwidths $\varepsilon$ and drag the slider to see the clusters change. ```{python} eps_values = [0.003, 0.01, 0.02, 0.03, 0.04, 0.05, 0.07, 0.1] sweep = {e: align(SpectralClustering(n_clusters=2, affinity="rbf", gamma=1 / (2 * e), random_state=0).fit_predict(X), y) for e in eps_values} for e, lab in sweep.items(): print(f"eps = {e:<5} agreement {np.mean(lab == y):6.1%} ARI {adjusted_rand_score(y, lab):.3f}") ``` ```{python} #| label: fig-bandwidth #| fig-cap: "Spectral clustering with the Gaussian graph for several bandwidths $\\varepsilon$; drag the slider to change $\\varepsilon$. For small $\\varepsilon$ the clusters are the two rings. As $\\varepsilon$ grows, the edges between the rings get stronger and the clustering falls apart." #| fig-alt: "A scatter plot of the two rings with a slider underneath labelled with bandwidth values from 0.003 to 0.1. At the smallest bandwidths the outer ring is blue and the inner ring orange. At the largest bandwidths the colors no longer follow the rings." #| code-fold: true #| code-summary: "Plot the bandwidth slider" start = eps_values.index(0.01) def colors(lab): return [COLOR[c] for c in lab] def hover(lab): return [f"{RING_NAME[r]}, cluster {c}" for r, c in zip(y, lab)] def title(e): lab = sweep[e] return (f"ε = {e} (gamma = {1 / (2 * e):g}): agreement {np.mean(lab == y):.0%}, " f"ARI {ari(y, lab):.2f}") fig = go.Figure(go.Scatter( x=X[:, 0], y=X[:, 1], mode="markers", showlegend=False, marker=dict(size=7, color=colors(sweep[eps_values[start]]), line=dict(color="white", width=0.8)), text=hover(sweep[eps_values[start]]), hovertemplate="%{text} (%{x:.2f}, %{y:.2f})<extra></extra>")) for g in (0, 1): # legend entries for the two colors fig.add_trace(go.Scatter(x=[None], y=[None], mode="markers", name=f"cluster {g}", marker=dict(size=9, color=COLOR[g]))) steps = [dict(label=f"{e}", method="update", args=[{"marker.color": [colors(sweep[e])], "text": [hover(sweep[e])]}, {"title.text": title(e)}, [0]]) for e in eps_values] fig.update_layout(sliders=[dict(active=start, steps=steps, x=0, len=1, pad=dict(t=50), currentvalue=dict(prefix="bandwidth ε = ", font=dict(color=INK)), font=dict(color=INK), bordercolor=GRID)]) lim = 1.45 fig.update_xaxes(range=[-lim, lim], title="first coordinate", constrain="domain") fig.update_yaxes(range=[-lim, lim], title="second coordinate", scaleanchor="x") style(fig, height=600, title=title(eps_values[start]), showlegend=True) fig.update_layout(legend=dict(x=1, y=1, xanchor="right", yanchor="top")) fig.show() ``` The rings are recovered for every $\varepsilon$ up to about $0.04$. Past that, the Gaussian kernel is wide enough to put real weight on edges between the rings, so there is no longer a [bottleneck](../../lectures/lecture-14.qmd#why-the-embedding-reveals-clusters) between them for the embedding to find. k-means fails for a different reason: no choice of parameters gives it a cluster shaped like a ring.