Source-linked AI summary
Data analysis recipes: Using Markov Chain Monte Carlo
David W. Hogg, Daniel Foreman-Mackey
TL;DR
MCMC is widely useful for probabilistic inference, but its practical use requires guidance on method choice, tuning, convergence, troubleshooting, and reporting. This paper provides a pedagogical user guide centered on sampling difficult posterior densities, emphasizing autocorrelation time for assessing sampling reliability and recommending representative posterior samples for reporting.
Problem
MCMC practice involves decisions about method choice, tuning, convergence, troubleshooting, and how to report posterior results from a sampling chain.
Method
The paper gives a practical overview of MCMC, including Metropolis-Hastings, initialization, tuning, convergence diagnosis, troubleshooting, and use of chain output.
Results
Autocorrelation time connects sampling efficiency and variance estimates, while randomly chosen posterior samples are recommended over reporting a single best sample.
Takeaways & Limitations
MCMC output is most useful when treated as a basis for sampling-based integrals and when posterior samples retain information needed for subsequent analyses.
Takeaways & Limitations
Autocorrelation time is difficult to estimate reliably, and MCMC convergence depends on assumptions including sufficient run length and connected posterior support.
Abstract
from arXiv · showhide
Markov Chain Monte Carlo (MCMC) methods for sampling probability density functions (combined with abundant computational resources) have transformed the sciences, especially in performing probabilistic inferences, or fitting models to data. In this primarily pedagogical contribution, we give a brief overview of the most basic MCMC method and some practical advice for the use of MCMC in real inference problems. We give advice on method choice, tuning for performance, methods for initialization, tests of convergence, troubleshooting, and use of the chain output to produce or report parameter estimates with associated uncertainties. We argue that autocorrelation time is the most important test for convergence, as it directly connects to the uncertainty on the sampling estimate of any quantity of interest. We emphasize that sampling is a method for doing integrals; this guides our thinking about how MCMC output is best used.
1. WHEN DO YOU NEED MCMC?
MCMC is useful for sampling posterior probability densities when their normalization is difficult to compute, but it should not be confused with an optimizer or parameter-space search algorithm. The paper presents MCMC as a sampling tool and introduces practical guidance for applying and diagnosing it in data analysis.
- 1. WHEN DO YOU NEED MCMC?: MCMC samples probability distributions or densities using ratios of values, without requiring a fully normalized analytic description.This normalization-insensitive property makes it suitable for posterior densities whose evidence or marginal likelihood is difficult to calculate.
- 1. WHEN DO YOU NEED MCMC?: In probabilistic inference, the posterior p(θ | D) combines the likelihood p(D | θ) with the prior p(θ), up to the normalization constant Z.Z is also called the evidence, marginal likelihood, Bayes integral, or prior predictive probability.
- 1. WHEN DO YOU NEED MCMC?: MCMC is appropriate when a model supplies a likelihood and prior for a potentially large parameter vector, while the posterior is hard to sample directly.The method does not require computing Z, and simple forms need neither derivatives nor integrals.
- 1. WHEN DO YOU NEED MCMC?: MCMC should not be used primarily to search parameter space or optimize the likelihood or posterior; those goals call for search algorithms or optimizers.Its intended role is sampling ill-normalized or otherwise hard-to-sample probability densities.
- 1. WHEN DO YOU NEED MCMC?: The paper is organized as a practical user manual covering Metropolis-Hastings, initialization, tuning, method choice, likelihood and prior construction, diagnosis, and troubleshooting.It emphasizes how to use MCMC in real data-analysis problems rather than presenting detailed mathematical properties of many methods.
2. WHAT IS A SAMPLING?
A sampling represents a probability density through draws whose empirical averages approximate integrals and expectation values. MCMC remains useful even when the target density is known only up to an unknown normalization constant, and samples can be marginalized or transformed to answer inference questions.
- What is a sampling?: A sampling is a set of draws whose large-sample histogram resembles the probability density from which they were drawn.The resemblance is understood up to total normalization.
- What is a sampling?: A good sampling is defined by how accurately its sample averages approximate expectation-value integrals, with accuracy improving as K increases.The required K depends on the expectations and accuracies of interest.
- What is a sampling?: MCMC can estimate expectations from an unnormalized function f(θ), without evaluating either the numerator or denominator integrals or learning the normalization constant Z.The function must be non-negative and have finite integral Z.
- What is a sampling?: Histograms reconstruct a pdf by averaging samples within bins, so even pdf reconstruction remains an integral-based use of sampling.The histogram is normalized so its step-function representation integrates to unity.
- What is a sampling?: Means, medians, and quantiles are preferred summaries because they involve integrals, whereas the mode does not.The paper therefore advises reporting mean, median, and quantiles rather than mode.
- What is a sampling?: Dropping nuisance-parameter dimensions from full samples produces samples from the marginalized pdf in the subspace of interest.The same idea generalizes from subspaces to arbitrary functions of the parameters.
3. METROPOLIS–HASTINGS MCMC
Metropolis–Hastings MCMC generates samples by proposing moves and accepting or repeating states, producing a biased Markov random walk whose stationary distribution is proportional to the target density. Correct implementation depends on the proposal distribution, coordinate consistency, numerical safeguards, and eventual convergence assessment.
- Algorithm: Metropolis–Hastings requires a target-density function f(θ) and a proposal pdf q(θ′ | θ) that generates candidate positions.In inference, f(θ) is typically the prior multiplied by the likelihood at the observed data.
- Algorithm: Each iteration draws a uniform random number, accepts a proposed position when the acceptance rule permits it, and otherwise repeats the previous sample.Iterating this procedure K times produces K samples.
- Why it works: The resulting chain is a biased Markov random walk: each step depends only on the current state, while acceptance ratios bias residence toward regions favored by the target function.The associated Markov process is intended to have a unique stationary distribution proportional to f(θ).
- Correctness conditions: A symmetric, properly normalized proposal satisfies detailed balance; asymmetric proposals require a corrected Metropolis–Hastings acceptance rule.The authors discourage violating detailed balance casually because both sampling and evaluating q(θ′ | θ) must then be handled correctly.
- Implementation: Proposals must use the same coordinates as the priors unless the acceptance probability includes the appropriate Jacobian.For example, proposing in ln θ changes a prior intended to be flat in θ if the Jacobian is omitted.
- Implementation: Log-density acceptance tests reduce underflow and overflow risks, but zero density can produce negative infinities and NaNs that code must catch.The specific failure case occurs when both proposed and current log densities are negative infinity.
4. LIKELIHOODS AND PRIORS
In data analysis, MCMC samples parameter values from a posterior density represented by a prior times a likelihood, often computed in log form. Using a likelihood alone does not sample the likelihood; it produces a posterior proportional to it.
- Posterior sampling: MCMC samples model parameters θ from a probability density proportional to an unnormalized positive function f(θ).The unknown normalization factor does not prevent sampling.
- Posterior construction: The ideal input combines a prior and likelihood, with the log-density computed as ln_f() = ln_prior() + ln_likelihood().Working in logarithms is usually advisable for implementation.
- Prior implementation: Bounded proper priors should return -Inf outside their parameter bounds, preventing likelihood evaluation for invalid parameter vectors.The pseudo-code checks prior validity before computing the likelihood and assumes the prior is easier to evaluate.
- Likelihoods versus posteriors: Using a likelihood or log-likelihood as the MCMC input does not sample the likelihood; it samples a posterior probability directly proportional to that function.A likelihood describes probability for data given parameters, whereas MCMC samples a density for parameters.
- Prior validity: Proper priors have finite normalization, although an improper prior can still yield a proper posterior; testing with a nearly constant likelihood can expose sampler and setup problems.The paper presents proper priors as good Bayesian practice while noting that posterior propriety is the relevant exception.
5. AUTOCORRELATION & CONVERGENCE
MCMC convergence cannot be established absolutely in general, so practitioners rely on heuristics and diagnostics. Integrated autocorrelation time connects chain dependence to sampling precision and is the fundamental basis for comparing sampler performance.
- Convergence limits: Absolute convergence is generally untestable because verifying complete posterior coverage would require solving an intractable global-search problem.Separated high-probability modes can cause a sampler to explore one mode while effectively never reaching another.
- Heuristic diagnostics: Practical convergence assessment uses heuristics: inspect chain traversal, compare substantial chain subsets, and apply diagnostics such as Gelman–Rubin.The Gelman–Rubin diagnostic compares within-chain variance with variance across chains.
- Autocorrelation time: τint is the number of steps required to produce an independent sample and determines the variance of sampling-based integral estimates.MCMC sampling is framed as Monte Carlo integration, with autocorrelation affecting the resulting error.
- Autocorrelation time: Lower integrated autocorrelation time means fewer function evaluations and less runtime for accurate sampling-based integral estimates.An autocorrelation time of unity is optimal, but is generally attainable only for analytically tractable sampling problems.
- Practical limitations: Autocorrelation time is difficult to estimate because it is a two-point statistic, making it harder to estimate than one-point statistics.The authors therefore present it as important for variance estimates while cautioning that its estimation is notoriously difficult.
- Sampler comparison: Burn-in time and acceptance ratio are less fundamental performance criteria than the precision of inferences achieved for a given computational cost.Burn-in depends strongly on initialization, tuning, and dynamics, whereas autocorrelation time governs precision for reasonably converged chains.
6. TUNING
MCMC tuning should target efficient exploration, with autocorrelation time as the principled objective and acceptance fraction or jump distance as practical proxies. Tuning belongs in burn-in, while diagnostics include autocorrelation estimates across window sizes and chain segments.
- The proposal distribution controls the kinds of steps a walker takes, with a D-dimensional Gaussian proposal offering D(D+1)/2 covariance parameters to set.
- The optimal proposal scale minimizes autocorrelation time, although estimating this second-order statistic reliably requires substantial data.
- Acceptance fraction provides a simple tuning proxy: nearly all accepted moves indicate steps that are too small, while nearly none indicate steps that are too large.
- Tuning must occur only during burn-in because adapting proposals from chain history violates the Markov property and undermines its guarantees.
- Expected Squared Jump Distance is more directly related to autocorrelation time than acceptance fraction and measures mean squared movement per step.
- Problem 14 treats minimizing autocorrelation time as preferable to targeting a fixed acceptance fraction, while Figure 7 displays both diagnostics against proposal scale.
7. INITIALIZATION AND BURN-IN
Initialization affects burn-in and can be improved by starting near a typical posterior location, while separated modes require multiple initializations and chain comparisons. MCMC is not designed for exhaustive parameter-space search or reliable optimization.
- Initializing near a typical or otherwise good posterior location can minimize burn-in, although starting exactly at the optimum is not generally ideal.
- Burn-in is the discarded beginning of a chain when initialization may be non-typical; it should be removed before inference.
- Badly multimodal posteriors should be tested with multiple chains started at different locations and compared for consistent inferences.
- If different initializations produce different posterior inferences, combining the resulting samples is not straightforward and may require long traversal or specialized samplers.
- MCMC samples generally need not lie near the posterior optimum, especially in high dimensions, so MCMC is not an optimizer or exhaustive search algorithm.
- Optimization before sampling is proposed as a comparison for reducing the dependence of burn-in on a distant starting position.
8. RESULTS, ERROR BARS, AND FIGURES
MCMC outputs should be reported as integral-based summaries and representative posterior samples rather than a single best sample. Means, medians, quantiles, histograms, and sample examples convey uncertainty, but marginal summaries can fail to form a jointly good-fitting parameter vector.
- MCMC results should emphasize expectations, medians, quantiles, and one- or multi-dimensional histograms because sampling is a method for computing integrals.
- For a scalar function of sampled parameters, the mean is an expectation estimate and the median is the midpoint of the ordered sampled values.
- A default one-dimensional report uses the posterior median with the 16th and 84th quantiles as lower and upper one-sigma error bars.
- The posterior mean can lie outside a sensible one-sigma interval when the posterior is strongly skewed.
- One-dimensional parameter summaries are means or medians of the corresponding marginalized posterior distributions.
- Coordinate-wise means or medians need not jointly form a good-fitting parameter vector in curved or banana-shaped posteriors.
- Randomly chosen example samples are recommended because they retain probabilistic information while also representing good-fitting models, unlike the best sample.
- Posterior samples intended for reuse should include prior-density values or an executable prior specification.
9. TROUBLESHOOTING AND ADVICE
Troubleshooting MCMC requires testing the code, priors, likelihood, initialization, tuning, and convergence diagnostics separately. Functional tests, visual checks, and reparameterization help identify failures and improve sampling behavior.
- Functional testing: End-to-end functional tests should sample distributions with known moments to verify that the sampler produces expected output.Testing a known Gaussian with non-trivial covariance provides one concrete example.
- Prior and likelihood checks: Sampling the prior with a flat likelihood tests both the sampler and whether coded priors are proper and correctly implemented.Replacing the log-likelihood with a function that always returns zero isolates prior behavior.
- Prior and likelihood checks: Prior samples that diverge or disagree with expectations indicate improper priors, coding errors, or mistaken assumptions about parameter transformations and boundaries.Practical priors can differ from intended priors because of logarithmic parameterizations, inversions, limits, or censoring.
- Prior and likelihood checks: Likelihood slices should vary smoothly across each parameter; noisy or jagged behavior can reveal numerical problems or bugs in the likelihood implementation.Internal numerical integration, truncated expansions, and adaptive approximations can create small-scale irregularities.
- Tuning and chain diagnostics: Low acceptance fractions call for smaller proposal steps, whereas high acceptance fractions should respond to larger proposal variance; otherwise, suspect a bug.The proposal distribution’s variance controls the step width in the described diagnostic.
- Tuning and chain diagnostics: Long horizontal chain segments indicate sticking, while initialization-dependent results indicate non-convergence and often problematic multimodality.Multiple modes may require broader-exploration methods such as nested sampling or simulated tempering.
- Parameterization: Reparameterization can reduce artificial multimodality from angular boundaries and correlations, but transformed priors must include the appropriate Jacobian.For amplitude-angle to Cartesian transformations, the cited Jacobian is 1/A.
10. MORE SOPHISTICATED SAMPLING METHODS
No single MCMC strategy works for every difficult target distribution, so method choice should reflect dimensionality, parameter structure, and available derivative information. The paper surveys ensemble, Gibbs, and Hamiltonian approaches while emphasizing their trade-offs and scope limits.
- Method choice: There is no general, problem-independent solution to slow or unreliable MCMC, because different strategies succeed or fail for different functions f(θ).The section offers field notes and pointers rather than a comprehensive algorithmic treatment.
- Ensemble methods: Ensemble methods use multiple walkers to generate tuned proposals when parameters are sufficiently equivalent or naturally form a vector-space representation.Their main target is the proposal-tuning challenge in such parameter spaces.
- Method choice: Most Metropolis–Hastings and ensemble methods are poor choices at very large parameter counts, where autocorrelation times can become extremely large.The passage notes that even dozens of parameters may be enough to create this difficulty.
- Gibbs sampling: Gibbs sampling is worth considering for many parameters with different scopes, including mixtures of global, local, and conditionally linear parameters.Local updates can exploit simpler calculations and, in hierarchical models, may be parallelized across data subsets.
- Hamiltonian methods: Hamiltonian methods perform well in high dimensions by using gradient information to inform proposals, but they require analytic derivatives for speed.Automatic differentiation systems can generate code for gradients of functions that can be written computationally.