Spectral Clustering and the Relaxation of Normalized Cut

Spectral clustering does not partition points directly. It first converts the data into a weighted graph whose edge weights encode local similarity, then partitions the graph by minimizing a cost that counts the weight of edges running between groups. The discrete optimization is relaxed to an eigenvector problem, and the final labels come from a centroid method applied to a small number of eigenvectors. The mechanism is worth understanding because the whole method’s behaviour, including its sensitivity to a scale parameter, follows from the two steps that precede the clustering.

From points to a weighted graph

Given n observations, form an n by n symmetric matrix W with non-negative entries. The common construction is the Gaussian kernel, w_ij = exp(-||x_i - x_j||^2 / (2 sigma^2)), which is positive everywhere and decays with distance; a second construction connects each point to its k nearest neighbours and sets other entries to zero, which produces a sparse matrix but usually an asymmetric one, so it is symmetrized by taking the maximum or the average of w_ij and w_ji. A third connects points within a fixed radius.

The matrix is not a distance matrix, and its entries are not comparable to distances: the kernel converts nearness into a weight, so two clusters separated by a gap larger than a few multiples of sigma still receive small but nonzero edge weights. The scale parameter sigma therefore controls how far similarity reaches, and it has no analogue in a centroid method. The degree of point i is the sum of its row, d_i = sum_j w_ij, and the diagonal matrix D holding those degrees is the other half of the construction. The volume of a set A is vol(A) = sum_{i in A} d_i, the total weight attached to its members.

Two cut costs, one trivial solution

The cut between two groups is the total weight of the edges crossing them, cut(A, B) = sum_{i in A, j in B} w_ij. Minimizing that quantity over all bipartitions is useless, because isolating a single point usually produces the smallest cut: one point has few edges to the rest, so separating it costs almost nothing. The fix is to normalize by the size of the groups. Ratio cut divides the crossing weight by the number of vertices in each part; the normalized cut of Shi and Malik divides it by the volume, giving

$$Ncut(A, B) = \frac{cut(A, B)}{vol(A)} + \frac{cut(A, B)}{vol(B)}$$

Each term is the fraction of the weight attached to one group that leaves it, so the objective penalizes cutting a group into pieces whose members are strongly attached to their own side. Both objectives are NP-hard to optimize over discrete partitions, so both are relaxed.

The relaxation and the Laplacian

Write the partition as a vector with two distinct values on the two parts. The cut cost becomes a quadratic form in that vector through the graph Laplacian L = D - W, since

$$x^T L x = \frac{1}{2} \sum_{i,j} w_{ij} (x_i - x_j)^2$$

which counts each crossing edge once per unit of squared difference across it. Constraints on the vector’s norms and orthogonality turn the minimization into a Rayleigh quotient problem whose solution is an eigenvector. Ratio cut leads to the second smallest eigenvector of L, the Fiedler vector, used as a real-valued indicator whose sign splits the graph. Normalized cut leads to the generalized problem L x = lambda D x, or equivalently to the symmetric form D^{-1/2} L D^{-1/2} y = lambda y after the substitution y = D^{1/2} x. For k groups, the embedder takes the k eigenvectors belonging to the smallest eigenvalues.

Three Laplacians appear in the literature and their spectra are related. The unnormalized L = D - W has eigenvalue zero with multiplicity equal to the number of connected components of the graph, and the corresponding eigenvectors can be chosen as the indicator vectors of those components. The symmetric normalized Laplacian L_sym = I - D^{-1/2} W D^{-1/2} has the same eigenvalues as the random-walk Laplacian L_rw = I - D^{-1} W, because the two are related by a similarity transform, but their eigenvectors differ by a factor of D^{-1/2}. Eigenvalues of L_sym and L_rw lie in the interval from 0 to 2, and eigenvalue zero occurs once per connected component.

The second smallest eigenvalue measures how hard the graph is to split. Cheeger-type bounds for the normalized Laplacian sandwich the graph’s conductance h, the minimum over vertex subsets of the ratio of crossing weight to volume, between lambda_2 / 2 and sqrt(2 lambda_2). A small lambda_2 guarantees a sparse cut exists and a large one forbids it, so the eigenvalue is a certificate rather than a heuristic.

The rounding step and why it is k-means again

The relaxation returns real-valued vectors, and labels require integers. The standard recipe takes the k smallest eigenvectors of the normalized Laplacian, stacks them as the columns of an n by k matrix, scales each row to unit Euclidean norm, and applies a centroid method to the rows. The row normalization matters: eigenvectors of L_sym are supported on the vertices with weights that vary with degree, so rows of small degree have small norm; normalizing to the sphere gives every vertex comparable influence, and it makes the rows of a graph with k perfect components coincide with the coordinate axes up to scale. In that idealized case the centroid step recovers the components exactly, which is the precise sense in which the method detects non-convex groups that no Voronoi partition can describe, unlike the cells discussed in Lloyd’s Algorithm and the Local Minima of k-Means.

The rounding step is where the method inherits its initialization dependence. The embedding is a deterministic function of the affinity matrix apart from eigenvector signs, but the final partition comes from a local search over centroid configurations, so different seeds can give different labels for the same embedding. The same caveat that applies to the agglomerative alternative, namely that greedy decisions are not globally optimal, applies here one level deeper.

Affinity construction decides the answer

The scale parameter dominates the outcome. A Gaussian affinity with a large sigma makes all pairs nearly equally similar, the graph approaches a uniform complete graph, and the leading eigenvectors carry little structural information, so the method returns one diffuse cluster or an arbitrary split. A small sigma makes the graph effectively a disjoint union of tiny components, with many points having no meaningful connection at all; then a row of the adjacency matrix is zero, the degree is zero, and the normalization divides by zero. Local scaling replaces the global sigma with a per-point scale, using the distance to each point’s own neighbours, so a point in a sparse region gets a wide kernel and a point in a dense region gets a narrow one.

Nearest-neighbour graphs add a second pitfall: an affinity graph that is disconnected into more components than the requested k produces those components regardless of k, because each component contributes its own eigenvalue zero. Disconnected graphs are frequent with sparse data, and the practical response is to raise the number of neighbours until the graph is connected, or to accept the connectivity structure as the answer and choose k accordingly.

Choosing k follows the eigengap: the sorted eigenvalues are examined and the largest jump between consecutive values is taken as a boundary between the structural eigenvalues and the rest. The heuristic fails on graphs whose degrees vary widely, and on graphs where the structural eigenvalues are close together, so it is a suggestion rather than a test. It is a different kind of criterion from the ones examined in Cluster Validity, Stability, and Choosing the Number of Clusters, since it reads the spectrum of the graph rather than the geometry of a produced partition.

import numpy as np
from scipy.linalg import eigh
from sklearn.cluster import KMeans

def spectral(X, k, sigma):
    d2 = ((X[:, None, :] - X[None, :, :]) ** 2).sum(-1)
    W = np.exp(-d2 / (2.0 * sigma ** 2))
    np.fill_diagonal(W, 0.0)                       # no self edges
    d = W.sum(1)
    Dm = np.diag(1.0 / np.sqrt(np.maximum(d, 1e-12)))
    Lsym = np.eye(len(X)) - Dm @ W @ Dm            # symmetric normalized Laplacian
    vals, vecs = eigh(Lsym)                        # dense: fine for modest n
    Y = vecs[:, :k]
    Y /= np.maximum(np.linalg.norm(Y, axis=1, keepdims=True), 1e-12)
    return KMeans(n_clusters=k, n_init=10).fit_predict(Y), vals[:k]

Cost, approximation, and where the method fits

The dense path costs O(n^3) time and O(n^2) memory for the factorization, which puts a hard ceiling on the number of points, since the affinity matrix alone holds n^2 entries. Two approximations change the scaling. A sparse affinity built from k nearest neighbours has O(n k) non-zeros, and a restarted Lanczos or block eigensolver extracts the few smallest eigenpairs of the sparse matrix without factoring it; whether that is faster in practice depends on the smallest spectral gap it has to resolve, which is exactly the quantity the method cares about. The Nyström extension samples a subset of columns, solves a smaller eigenproblem, and extends the eigenvectors to the remaining points, with accuracy governed by how well the sampled columns span the leading subspace.

Spectral clustering is best understood as a centroid method operating in a space whose inner product is the affinity, since the relaxation is the continuous analogue of kernel k-means with the graph’s kernel. The embedding it produces is also usable on its own: the eigenvectors with the smallest eigenvalues are the smoothest functions on the graph, and the ones with large eigenvalues oscillate across individual edges, which is the same separation of scales that makes the linear projection in Principal Component Analysis Through the Singular Value Decomposition a useful first look at the data. The neighbourhood-graph machinery reappears in nonlinear embedding, where the graph is used for a different objective rather than for a cut.