Diagnosing MCMC Convergence With Rank-Normalized R-Hat

Why Compare Independent Chains?

A QSP calibration can produce a long sequence of parameter values that looks stable while still representing only one region of the posterior. For example, several chains may appear to fluctuate around different levels, or they may have similar centers but different spreads. Comparing independent chains helps reveal these problems: if the chains have not settled into comparable distributions, their pooled draws may not reliably represent the target posterior.

R-hat, also called the potential scale reduction factor, summarizes how much the estimated variation across chains exceeds the variation within chains. Rank-normalized split R-hat improves the basic comparison by splitting chains and transforming pooled draws to normal scores based on their ranks. A folded version adds sensitivity to differences in spread. It is a diagnostic, not a proof that a Markov chain has converged.

Definition: Given $m$ original chains with retained draws, split R-hat divides each chain into two equal-length halves, producing $M=2m$ split chains of $n$ draws each. After pooling the $S=Mn$ split-chain draws, rank normalization assigns each draw its pooled rank $r$ and transforms it to $$ z=\Phi^{-1}\left(\frac{r-3/8}{S+1/4}\right), $$ where $\Phi^{-1}$ is the standard normal quantile function. Tied values receive their average rank. Rank-normalized split R-hat is the ordinary split-chain R-hat calculated on these $z$ values.

The rank transform depends on the ordering of draws rather than their original numerical spacing. This makes it useful when posterior draws are skewed or have heavy tails. It does not make the chains independent or repair poor exploration; those remain requirements of the sampling procedure and diagnostic setup.

Calculating Split R-Hat

Suppose the transformed draws in split chain $j$ have mean $\bar{z}_j$ and sample variance $s_j^2$. The within-chain variance $W$ is the average of those variances. The between-chain variance $B$ scales the sample variance of the split-chain means by the half-chain length $n$. The estimated marginal variance $\widehat{\mathrm{var}}^+$ combines these quantities.

\bar{z}_{\cdot}=\frac{1}{M}\sum_{j=1}^{M}\bar{z}_j,\qquad W=\frac{1}{M}\sum_{j=1}^{M}s_j^2,\qquad B=\frac{n}{M-1}\sum_{j=1}^{M}(\bar{z}_j-\bar{z}_{\cdot})^2

\widehat{\mathrm{var}}^+=\frac{n-1}{n}W+\frac{B}{n},\qquad \hat{R}=\sqrt{\frac{\widehat{\mathrm{var}}^+}{W}}

The same calculation applies to folded rank-normalized draws. For each original-scale pooled draw $\theta$, first calculate its absolute distance from the pooled median, $|\theta-\operatorname{median}(\theta)|$. Rank-normalize those distances, then calculate split R-hat again. The reported rank-normalized split R-hat is the larger of the ordinary rank-normalized value and the folded rank-normalized value.

Assumptions and setup: Use multiple independently initialized chains targeting the same posterior distribution, and exclude warm-up draws. Split chains must have equal length; if the retained length is odd, an end draw can be discarded so the two halves match. The R-hat comparison is meaningful only if chains are sampling the same target and the retained draws are produced by a valid sampling procedure. Starting chains in dispersed, plausible regions can make failures easier to detect. A value near one indicates similar within- and between-chain variation, but does not establish that all relevant posterior regions have been visited.

The formula has a simple justification. If the split chains explore the same distribution, their means should differ only by ordinary sampling variation. In that case, variation between chain means should be comparable to the uncertainty expected from within-chain variation, and $\widehat{\mathrm{var}}^+$ should be close to $W$. If chain means are more dispersed than expected, $B$ rises, increasing the ratio. This is a diagnostic comparison, not an exact test with a universal pass/fail probability.

Worked Rank Calculation

Worked Example: Rank Scores When Chains Occupy Different Regions

Consider a deliberately stark simulated illustration with two original chains, each retaining eight draws. Splitting gives $M=4$ chains of $n=4$, and $S=16$ pooled draws. Suppose their pooled ranks occur in consecutive groups: the first split chain contains ranks 1–4, the second 5–8, the third 9–12, and the fourth 13–16. This ordering illustrates chains separated in location, not a realistic claim about a QSP study.

For rank 1, the Blom probability is $(1-3/8)/(16+1/4)=0.625/16.25=0.03846$, giving a normal score of about $-1.769$. Applying the same formula to all ranks gives the approximate score sequence:

$-1.769,-1.282,-0.988,-0.761,-0.569,-0.396,-0.233,-0.077,0.077,0.233,0.396,0.569,0.761,0.988,1.282,1.769$.

The four split-chain means are approximately $-1.200,-0.3188,0.3188,1.200$. Their grand mean is zero. The first and fourth chains each have sample variance $0.56815/3=0.1894$; the second and third each have sample variance $0.13439/3=0.04480$. Thus:

W=\frac{2(0.1894)+2(0.04480)}{4}=0.1171

The sample variance of the four means is approximately $(1.200^2+0.3188^2+0.3188^2+1.200^2)/3=1.0277$, so $B=4(1.0277)=4.111$. Therefore:

\widehat{\mathrm{var}}^+=\frac{3}{4}(0.1171)+\frac{4.111}{4}=1.116,\qquad \hat{R}=\sqrt{\frac{1.116}{0.1171}}\approx 3.087

This very large value reflects the non-overlapping rank groups: between-chain differences dominate within-chain variation. The calculation demonstrates the diagnostic; it does not substitute for inspecting the actual traces and model context.

Worked Examples: Near-One and Folded R-Hat

Worked Example: Eight Split Chains With Similar Locations

Suppose four independent hypothetical chains are split, producing $M=8$ half-chains, each of length $n=500$. On the rank-normalized values, assume the average within-chain variance is $W=1.04$. The eight split-chain means are $-0.04,-0.02,0.02,0.04$, each occurring twice. Their grand mean is zero.

The sum of squared deviations is $2(0.04^2+0.02^2+0.02^2+0.04^2)=0.008$. The sample variance of the eight means is $0.008/7=0.001143$, giving $B=500(0.001143)=0.5714$. Hence:

\widehat{\mathrm{var}}^+=\frac{499}{500}(1.04)+\frac{0.5714}{500} =1.03906,\qquad \hat{R}=\sqrt{\frac{1.03906}{1.04}}\approx 0.9995

This value is essentially one, indicating little location disagreement by this summary. A value slightly below one can occur from the finite-sample variance calculation; it is not evidence of “better than converged” chains. The conclusion is limited to this diagnostic and these draws.

Worked Example: Folded R-Hat Detects Differences in Spread

Now consider a different simulated calibration in which the ordinary rank-normalized split R-hat is $1.006$, but the folded rank-normalized values have $W=0.80$. There are again eight split chains of length $n=100$. Suppose their folded-score means are $-0.60,-0.20,0.20,0.60$, each occurring twice. The grand mean is zero.

The squared deviations sum to $2(0.60^2+0.20^2+0.20^2+0.60^2)=1.6$. The sample variance of the eight means is $1.6/7=0.2286$, so $B=100(0.2286)=22.86$. Then:

\widehat{\mathrm{var}}^+=\frac{99}{100}(0.80)+\frac{22.86}{100} =1.02057,\qquad \hat{R}_{\mathrm{folded}}=\sqrt{\frac{1.02057}{0.80}}\approx 1.1295

The final rank-normalized split R-hat is the larger value, $1.1295$. The ordinary rank comparison did not flag a substantial location difference, but the folded comparison indicates that the chains differ in their distances from the pooled median, consistent with a scale mismatch. This warrants investigation rather than pooling the chains as if the discrepancy were harmless.

Diagnostics, Alternatives, and Limitations

R-hat is most informative when used alongside trace plots and knowledge of the model. Plot each chain separately for important parameters and derived QSP quantities. Look for persistent separation, long periods without movement, abrupt transitions, and different apparent spreads. Calculate diagnostics for parameters that matter to the scientific conclusions, including relevant derived quantities; a favorable value for one parameter does not guarantee favorable values for the rest.

A practical convention is to investigate rank-normalized split R-hat values above about $1.01$, especially when the parameter affects decisions or mechanistic interpretation. This is a useful warning threshold, not a theorem that guarantees adequate sampling below it. Short chains, rare posterior regions, and insufficiently dispersed starting values can make a diagnostic look reassuring when exploration is still inadequate.

Compared with unsplit, unranked R-hat, the split rank-normalized version is more sensitive to within-chain nonstationarity and less dependent on the original scale and tail shape. Folding adds a check for disagreement in spread as well as location. Trace plots provide a visual account of time ordering that a single number cannot show; conversely, visual inspection alone can miss modest chain differences. These checks complement each other.

Do not interpret near one as proof. Chains can agree while all miss the same region, particularly if initialization and exploration are poor.
Do not average away a warning. A problematic parameter or folded value deserves investigation even if most other R-hat values are near one.
Do not treat the threshold as universal. The consequences of a discrepancy depend on the posterior quantity and the intended use of the calibration.
Do not ignore implementation issues. Check that warm-up was excluded, chains target the same posterior, and parameter transformations and recorded quantities are consistent.

A common mistake is to say, “The trace plot looks flat, so the chains converged,” or “R-hat is close to one, so the posterior is fully explored.” Flatness can conceal chains stuck in different regions, and R-hat compares the chains that were actually run rather than proving that they found every important region. The correction is to use independent chains, examine split rank-normalized and folded diagnostics, inspect chain-specific traces, and treat disagreement as a reason to investigate sampling and model geometry.

Key takeaway: Rank-normalized split R-hat compares within- and between-chain variation after splitting and pooling ranks. The folded version checks spread differences. A value near one is reassuring only as one part of a broader convergence assessment.

Check Your Understanding

Try each problem, then compare your reasoning with the solutions.

Four original chains each have an even number of retained draws. How many split chains result, and what information does splitting help reveal?
Eight split chains have length $n=100$, $W=1.00$, and a sample variance of split-chain means equal to $0.02$. Calculate $B$, $\widehat{\mathrm{var}}^+$, and split R-hat.
Why can folded rank-normalized R-hat be high when ordinary rank-normalized R-hat is near one?
What can and cannot be concluded if all monitored R-hat values are below $1.01$?

Solutions. (1) Splitting produces $M=8$ half-chains. It helps reveal disagreement between early and late portions of a chain, which an unsplit summary can obscure. (2) $B=100(0.02)=2$. Then $\widehat{\mathrm{var}}^+=(99/100)(1.00)+2/100=1.01$, and $\hat{R}=\sqrt{1.01/1.00}\approx1.005$. (3) Folding ranks absolute distances from the pooled median, so it is sensitive to chains with different spreads even when their locations are similar. (4) The values indicate little measured between-chain disagreement for the monitored quantities, but do not prove that all posterior regions were explored or that unmonitored quantities are reliable. Chain traces and model context still matter.