The Gap Statistic for Choosing Cluster Number

Why Compare Dispersion With a Reference?

Within-cluster dispersion usually decreases as the number of clusters increases. A small dispersion by itself therefore does not establish that a clustering has found meaningful structure: a sufficiently large number of clusters can make even unstructured data look compact. The gap statistic addresses this issue by asking how much smaller the observed dispersion is than the dispersion expected under a specified reference model with no cluster structure.

For single-cell RNA sequencing, the comparison is meaningful only after defining the representation being clustered. For example, the analysis might use cells represented by a fixed set of normalized, transformed, and scaled features or by fixed coordinates in a reduced-dimensional space. The observed cells and every reference data set must be treated in a comparable space, with the same distance, clustering method, and procedure for choosing a partition at each $k$.

Definition: For a partition into $k$ clusters $C_1,\ldots,C_k$, let $W_k$ denote its within-cluster dispersion. For Euclidean $k$-means, a standard choice is the sum of squared distances from each point to its assigned cluster centroid. Let $W^*_{kb}$ be the dispersion at $k$ for reference data set $b$, where $b=1,\ldots,B$. The gap statistic is $$ \operatorname{Gap}(k)=\frac{1}{B}\sum_{b=1}^{B}\log(W^*_{kb})-\log(W_k). $$

The logarithms compare dispersions multiplicatively: a reference dispersion twice as large as the observed dispersion gives a log difference of $\log(2)$, regardless of their absolute scale. Since the same fixed distance scale is used for both, multiplying every dispersion by a common positive constant does not change the gap. A positive gap means that the observed dispersion is smaller than the average log dispersion under the chosen reference model; it does not, by itself, establish biological validity.

Dispersion and the Reference Model

For $n$ points $x_1,\ldots,x_n$ and a cluster $C_r$ of size $n_r$, the usual $k$-means dispersion can also be written as a sum of pairwise squared distances:

W_k=\sum_{r=1}^{k}\frac{1}{2n_r} \sum_{i,j\in C_r}\lVert x_i-x_j\rVert^2 =\sum_{r=1}^{k}\sum_{i\in C_r}\lVert x_i-\bar{x}_r\rVert^2,

where $\bar{x}_r$ is the centroid of cluster $r$. This identity follows by expanding the pairwise squared distances: the terms involving each point’s squared norm contribute $\sum_{i\in C_r}\lVert x_i\rVert^2$, while the cross-products contribute $n_r\lVert\bar{x}_r\rVert^2$. Their difference is exactly $\sum_{i\in C_r}\lVert x_i-\bar{x}_r\rVert^2$. Thus either expression gives the same dispersion for this definition.

A reference model is a distribution used to generate data sets that represent the absence of the structure being assessed while retaining relevant features of the analysis. A common construction generates points uniformly within the observed data’s bounding region; a rotated bounding region can be used when the feature axes are strongly correlated. These are choices, not universal truths. In single-cell data, a uniform reference may fail to preserve dependence, manifold structure, or technical variation. The resulting gap is always relative to the reference model actually chosen.

The comparison requires several conditions. The observed and reference data must have the same number of points and the same features or coordinates; the same distance and definition of $W_k$ must be used; and the clustering procedure must be applied comparably at every candidate $k$. If the clustering algorithm is randomized, use a controlled, sufficiently reliable optimization procedure for both observed and reference data. Otherwise, variation in optimization can be mistaken for variation under the reference model. These conditions make the numerical comparison interpretable, but they do not guarantee that the reference distribution is scientifically appropriate.

Worked Example: Computing Dispersion and a Gap

Suppose four hypothetical cells have one-dimensional coordinates $0,2,8,10$, and a two-cluster $k$-means partition assigns $\{0,2\}$ and $\{8,10\}$. The centroids are $1$ and $9$. The observed dispersion is $$ W_2=(0-1)^2+(2-1)^2+(8-9)^2+(10-9)^2=1+1+1+1=4. $$ The two clusters have equal size, and the same point representation and squared Euclidean distance will be used for the reference data, so this dispersion is comparable with reference dispersions.

For illustration, suppose $B=3$ reference data sets produce log dispersions $2.0,\ 2.2,\ 2.1$ at $k=2$. Their mean is $(2.0+2.2+2.1)/3=2.1$, while $\log(W_2)=\log(4)=1.3863$. Therefore, $$ \operatorname{Gap}(2)=2.1-1.3863=0.7137. $$ The positive value says the observed dispersion is smaller than the reference average on the log scale. It does not mean the probability that the partition is biologically correct is $0.7137$, nor does it test a null hypothesis by itself.

Uncertainty and the Selection Rule

The reference data sets approximate the distribution of log dispersion under the specified null model. Let $L_{kb}=\log(W^*_{kb})$, and let $sd_k$ be the sample standard deviation of the $B$ values $L_{k1},\ldots,L_{kB}$. The Monte Carlo standard error of the estimated gap, conditional on the observed dispersion, is $sd_k/\sqrt{B}$: the observed $\log(W_k)$ is fixed, and the reference contribution is an average of $B$ simulated values.

The conventional quantity used in the gap selection rule is different:

s_k=sd_k\sqrt{1+\frac{1}{B}}.

This quantity is not the Monte Carlo standard error of the estimated gap. Its factor has a prediction-uncertainty interpretation: under independent draws with variance $\sigma_k^2$, the variance of one future reference log dispersion minus the mean of $B$ reference log dispersions is $\sigma_k^2+\sigma_k^2/B$. Replacing $\sigma_k$ with $sd_k$ gives the conventional $s_k$ used for the selection rule. This rule compares adjacent $k$ values and chooses the first $k$ for which the observed gap is close enough to the next value, allowing for this uncertainty quantity.

Selection rule: Choose the smallest $k$ such that $$ \operatorname{Gap}(k)\geq \operatorname{Gap}(k+1)-s_{k+1}. $$ The comparison uses $s_{k+1}$, the uncertainty quantity for the next candidate cluster count. It is a practical one-standard-error rule, not a confidence interval for the true number of biological cell types.

Worked Example: Applying the One-Standard-Error Rule

Suppose the following hypothetical results were obtained using the same reference construction and clustering procedure at each $k$.

k	Gap(k)	s(k)
1	0.20	Not used in the first comparison
2	0.55	0.12
3	0.72	0.10
4	0.68	0.09
5	0.62	0.08

At $k=1$, the rule requires $0.20\geq 0.55-0.12=0.43$, which is false. At $k=2$, it requires $0.55\geq 0.72-0.10=0.62$, also false. At $k=3$, it requires $0.72\geq 0.68-0.09=0.59$, which is true. The rule therefore selects the smallest qualifying value, $k=3$. The gap itself is slightly larger at $k=3$ than at $k=4$; the rule’s point is to stop at the first $k$ for which the next increase is not compelling relative to the prescribed uncertainty quantity.

Worked Example: Separating Monte Carlo Error From the Rule Quantity

Return to the three reference log dispersions $2.0,\ 2.2,\ 2.1$ from the first example. Their mean is $2.1$. The sample standard deviation is $$ sd_2=\sqrt{\frac{(2.0-2.1)^2+(2.2-2.1)^2+(2.1-2.1)^2}{3-1}} =\sqrt{\frac{0.02}{2}}=0.1000. $$

The Monte Carlo standard error of the estimated gap is $0.1000/\sqrt{3}=0.05774$. The conventional quantity for the selection rule is instead $$ s_2=0.1000\sqrt{1+\frac{1}{3}}=0.1155. $$ These numbers answer different questions and should not be labeled interchangeably. With only three simulations, both are estimates based on a small reference sample; increasing $B$ can reduce Monte Carlo error and stabilize the estimated reference mean and standard deviation.

Reference-Model Sensitivity and Diagnostics

A gap statistic should be interpreted together with a description of how its reference data were generated. The following checks help determine whether the calculation is answering the intended question:

Check the scale and representation. Confirm that observed and reference points use the same features, transformations, distance, and dispersion definition. A mismatch makes the log-dispersion difference difficult to interpret.
Inspect reference simulations. Check for failed clustering runs, zero or non-finite dispersions, and implausible generated values. Since the logarithm is used, $W^*_{kb}$ must be positive.
Check Monte Carlo precision. A small $B$ can make the reference average and $sd_k$ noisy. Report $B$, and distinguish $sd_k/\sqrt{B}$ from $s_k$.
Consider whether the null is plausible. If the reference destroys important correlation or geometry, a large gap may reflect that mismatch rather than convincing cluster structure in the observed data.
Check the clustering optimization. Different optimization quality between observed and reference data can distort $W_k$. Use comparable settings and enough starts to reduce avoidable optimization effects.

One practical sensitivity check is to recompute gaps under scientifically defensible reference constructions. Consider this hypothetical calculation at $k=3$, with observed $\log(W_3)=2.00$ and $B=4$. Under one reference construction, the log dispersions are $2.40,\ 2.50,\ 2.30,\ 2.40$, with mean $2.40$, so the gap is $2.40-2.00=0.40$. The sample standard deviation is $\sqrt{0.02/3}=0.08165$, giving Monte Carlo standard error $0.08165/\sqrt{4}=0.04083$ and rule quantity $0.08165\sqrt{1+1/4}=0.09129$. Under a second construction, suppose the log dispersions are $2.15,\ 2.20,\ 2.10,\ 2.15$. Their mean is $2.15$, so the gap is $0.15$; the sample standard deviation is $\sqrt{0.005/3}=0.04082$, the Monte Carlo standard error is $0.02041$, and the rule quantity is $0.04564$. The different gaps show that the reference choice matters. Neither result is automatically preferable; the constructions must be justified by the scientific question and data representation.

What the Gap Statistic Does Not Establish

The gap statistic is one internal validity criterion, not a universal measure of clustering quality. Silhouette widths focus on how separated points are from neighboring clusters; Calinski–Harabasz and Davies–Bouldin indices summarize different relationships between within- and between-cluster dispersion. Agreement with reference labels addresses a different question, and requires labels whose meaning and reliability are understood. These measures can disagree because they target different properties. The gap statistic’s distinctive feature is comparison with a specified reference distribution, not proof that the selected partition captures discrete biological types.

In single-cell data, a selected $k$ can be sensitive to the chosen features, distance, normalization, reference geometry, and clustering algorithm. Continuous biological gradients, technical effects, unequal cluster sizes, or rare populations can all complicate interpretation. A gap-based choice should therefore be reported as a choice supported relative to a stated reference model and dispersion definition. It does not demonstrate that the data contain exactly $k$ natural populations, and it does not establish differential expression or causal biological mechanisms.

Common Mistakes

Choosing the largest gap automatically. The rule selects the smallest $k$ meeting an adjacent-$k$ condition; it is not simply the $k$ with the maximum gap.
Calling $sd_k/\sqrt{B}$ the rule quantity. This is the Monte Carlo standard error of the estimated gap conditional on observed dispersion. The conventional selection quantity is $sd_k\sqrt{1+1/B}$.
Treating the reference model as neutral. The gap is conditional on the null construction. A poorly matched null can make the comparison misleading.
Interpreting a positive gap as biological confirmation. It only says observed dispersion is below the reference average on the log scale; biological interpretation requires additional evidence.
Changing the analysis between observed and reference data. Inconsistent preprocessing or optimization prevents a fair dispersion comparison.

Key takeaway: Calculate $\operatorname{Gap}(k)$ as mean reference log dispersion minus observed log dispersion, and use $sd_{k+1}\sqrt{1+1/B}$ in the stated adjacent-$k$ selection rule. Keep this rule quantity distinct from the Monte Carlo standard error $sd_k/\sqrt{B}$, and interpret the result relative to the reference model.

Check Your Understanding

Answer each question, then compare with the solutions below.

For a fixed $k$, reference log dispersions are $1.8,\ 2.0,\ 2.2,\ 2.0$, and the observed log dispersion is $1.5$. Calculate the gap.
For the values in Question 1, calculate the sample standard deviation, the Monte Carlo standard error of the gap, and the conventional selection-rule quantity.
If $\operatorname{Gap}(2)=0.40$, $\operatorname{Gap}(3)=0.50$, and $s_3=0.12$, does $k=2$ satisfy the selection rule?
Why can the same observed data produce different gap statistics under two reference constructions?
A colleague says, “The gap is positive, so the probability that these clusters are real is high.” Identify the error and give a correct interpretation.

Solutions. (1) The reference mean is $(1.8+2.0+2.2+2.0)/4=2.0$, so the gap is $2.0-1.5=0.5000$. (2) The sum of squared deviations from 2.0 is $0.04+0+0.04+0=0.08$, so $sd_k=\sqrt{0.08/3}=0.1633$. The Monte Carlo standard error is $0.1633/\sqrt{4}=0.08165$. The selection-rule quantity is $0.1633\sqrt{1+1/4}=0.1826$. (3) The right side is $0.50-0.12=0.38$, and $0.40\geq0.38$, so $k=2$ satisfies the rule. (4) The reference construction defines what dispersion is expected without the targeted structure; changing it changes the simulated log dispersions and therefore the reference mean and gap. (5) A positive gap is not a probability that clusters are real. It means observed dispersion is below average reference log dispersion under the chosen null model; it is not a posterior probability or stand-alone biological validation.