An End-to-End MCMC Diagnostic Workshop

From Chain Output to a Defensible Summary

MCMC diagnostics are most useful when they answer a specific question: does this output support the posterior summary we plan to report, to the precision we need? No single number answers that. Acceptance rates describe proposal behavior; trace plots reveal patterns over time; rank-normalized R-hat checks agreement among chains; effective sample size and Monte Carlo standard error address the simulation precision of particular summaries. These checks can expose problems, but they cannot prove that a chain has explored every relevant region of a target distribution.

Let the target distribution be \(\pi\) on a state space \(S\), and suppose \(C\) chains have been run. Write \(X_{c,n}\) for the state at retained iteration \(n\) of chain \(c\), where \(c=1,\ldots,C\) and \(n=1,\ldots,N\). Each chain is a stochastic process on a probability space \((\Omega,\mathcal{F},\mathbb{P})\), with discrete time \(n\). If \(\mathcal{F}_{c,n}\) records the output of chain \(c\) through iteration \(n\), then its transition rule uses the current state (and, for an algorithm with auxiliary variables, its specified update mechanism). Diagnostics are applied to retained output after warmup; they do not make an inadequate warmup adequate.

For rank-normalized split R-hat, each retained chain is divided into two equal-length halves, discarding an unmatched final draw if necessary. The half-chains are treated as separate chains. Pool their values, assign ranks (using average ranks for ties), and transform rank \(r\) among \(L\) values to \(\Phi^{-1}((r-3/8)/(L+1/4))\), where \(\Phi\) is the standard normal distribution function. Compute the split R-hat on these transformed values. A folded version first ranks the absolute deviations from the pooled median; the reported rank-normalized R-hat is the larger of the ordinary and folded versions. This construction is more sensitive than the original variance-based diagnostic to some differences in location, scale and tails. It remains a diagnostic, not a convergence certificate.

For \(m\) chains of equal split length \(n\), the usual variance calculation on any transformed values \(Y_{c,j}\) uses chain means \(\bar Y_c\), within-chain sample variances \(s_c^2\), their average \(W=m^{-1}\sum_c s_c^2\), and \(B=n(m-1)^{-1}\sum_c(\bar Y_c-\bar Y_{\cdot})^2\), where \(\bar Y_{\cdot}\) is the mean of the chain means. The variance estimate is \(\widehat V=(n-1)W/n+B/n\), and the split R-hat is \(\sqrt{\widehat V/W}\) when \(W>0\). If \(W=0\), this ratio is undefined; constant chains must be investigated, not assigned a reassuring finite value.

Proof. The pooled sample variance of the \(mn\) transformed observations is the sum of squared deviations from \(\bar Y_{\cdot}\), divided by \(mn-1\). Within each chain, the sum of squared deviations from \(\bar Y_c\) is \((n-1)s_c^2\), so summing over chains gives \(m(n-1)W\). The remaining between-chain sum is \(n\sum_c(\bar Y_c-\bar Y_{\cdot})^2=(m-1)B\). Thus the pooled sum of squares is \(m(n-1)W+(m-1)B\). Dividing by \(mn-1\) gives the pooled variance. The conventional \(\widehat V\) above is a weighted within- and between-chain variance estimate, rather than that pooled variance; its two weights \((n-1)/n\) and \(1/n\) sum to one. The formula therefore gives more weight to within-chain variation as \(n\) grows while retaining a between-chain component. The ratio compares this estimate with \(W\), which proves the stated calculation, not that the chains have converged. \(\square\)

Precision for a Particular Quantity

ESS depends on the function being estimated. A chain can estimate a smooth posterior mean adequately while providing little information about a rare tail event. Under a Markov chain central limit theorem for \(f(X)\), with finite positive posterior variance \(\sigma_f^2=\operatorname{Var}_{\pi}(f)\) and finite long-run variance \(\sigma^2_{\mathrm{asym}}\), the ESS for \(f\) is \(N\sigma_f^2/\sigma^2_{\mathrm{asym}}\). The corresponding MCSE is estimated by \(\sqrt{\widehat{\sigma}^2_{\mathrm{asym}}/N}\), equivalently \(\widehat{\sigma}_f/\sqrt{\mathrm{ESS}_f}\). These estimates require a suitable estimate of long-run variance and adequate chain behavior; positive recurrence alone does not guarantee a central limit theorem or a reliable finite-run MCSE.

Proof. Under the stated central limit theorem, the variance of the sample average of \(f\), to first order in \(N\), is \(\sigma^2_{\mathrm{asym}}/N\). By the definition of ESS, \(\mathrm{ESS}_f=N\sigma_f^2/\sigma^2_{\mathrm{asym}}\). Rearranging gives \(\sigma^2_{\mathrm{asym}}/N=\sigma_f^2/\mathrm{ESS}_f\). Taking square roots yields \(\sigma_{\mathrm{asym}}/\sqrt N=\sigma_f/\sqrt{\mathrm{ESS}_f}\), the relation used for the plug-in MCSE. This is a variance relation under the assumed limit theorem; it does not certify that an estimated long-run variance is accurate. \(\square\)

There is no universal ESS threshold that makes every summary adequate. Choose a tolerance based on the use of the result. For example, a posterior mean used to distinguish effects of practical importance may need a much smaller MCSE than a preliminary exploratory estimate. For a tail probability or quantile, inspect the relevant tail diagnostics and its own uncertainty, rather than borrowing the ESS for the mean.

Worked Diagnostic Workflow

The examples treat rank-normalized R-hat and ESS entries as supplied diagnostic output calculated from the stated retained chains. The arithmetic below checks proposal rates and MCSEs, then evaluates the evidence. In practice, also inspect the actual trace plots: numerical summaries cannot reveal every pattern visible in the time series.

Reading the Evidence Together

A practical assessment proceeds from the chain behavior to the particular reported quantity. First verify that the chains target the intended distribution and that warmup and retained draws have been handled consistently. Next inspect traces for drift, sticking, abrupt transitions and chain-to-chain separation; running means and autocorrelation plots can add context. Then review acceptance rates in light of the update method. A low or high rate may indicate a tuning problem for some Metropolis proposals, but there is no single ideal rate for all algorithms, dimensions or targets. Finally examine rank-normalized split R-hat, bulk and tail ESS, and the MCSE for each intended summary.

R-hat values near one are commonly sought, with \(1.01\) often used as a practical warning threshold, not a theorem. ESS should be judged against the precision required for the estimand. MCSE describes simulation uncertainty, not posterior uncertainty: a narrow MCSE does not mean the posterior itself is concentrated, and a wide posterior can be estimated precisely. Conversely, a small estimated MCSE is not trustworthy if the chains failed to explore important target regions or the long-run variance estimate is unreliable.

Common mistake. Suppose someone sees the \(0.20\) acceptance rate in the separated-chain example and concludes the chains mix well. That conclusion confuses proposal acceptance with exploration of the target. The arithmetic \(2{,}000/10{,}000=0.20\) is correct; the inference is not. The correction is to weigh acceptance behavior alongside trace plots and between-chain diagnostics, then check precision for the specific estimand. Likewise, do not report \(0.02\) as the MCSE for every summary just because it is the mean’s MCSE: the tail indicator has its own ESS and uncertainty.