Source-linked AI summary
MCMC using Hamiltonian dynamics
Radford M. Neal
TL;DR
The paper addresses inefficient random-walk exploration and the practical difficulty of tuning Hamiltonian Monte Carlo. It reviews HMC and several extensions, showing favorable dimensional scaling while identifying important limitations for multimodal sampling.
Problem
Hamiltonian Monte Carlo requires suitable choices of leapfrog stepsize and trajectory length, making it harder to tune than simple Metropolis methods.
Method
The paper reviews Hamiltonian Monte Carlo, its volume-preserving dynamics, and extensions including windows, fast approximations, tempering, and trajectory shortcuts.
Results
Computation time typically grows as d^5/4 for HMC to reach a nearly independent point, versus d^2 for random-walk Metropolis.
Takeaways & Limitations
HMC can avoid inefficient random-walk exploration and scale more favorably than random-walk Metropolis as dimensionality increases.
Takeaways & Limitations
Tempered HMC may fail to move well between modes when one is high and narrow and the other low and broad, especially in high dimensions.
Abstract
from arXiv · showhide
Hamiltonian dynamics can be used to produce distant proposals for the Metropolis algorithm, thereby avoiding the slow exploration of the state space that results from the diffusive behaviour of simple random-walk proposals. Though originating in physics, Hamiltonian dynamics can be applied to most problems with continuous state spaces by simply introducing fictitious "momentum" variables. A key to its usefulness is that Hamiltonian dynamics preserves volume, and its trajectories can thus be used to define complex mappings without the need to account for a hard-to-compute Jacobian factor - a property that can be exactly maintained even when the dynamics is approximated by discretizing time. In this review, I discuss theoretical and practical aspects of Hamiltonian Monte Carlo, and present some of its variations, including using windows of states for deciding on acceptance or rejection, computing trajectories using fast approximations, tempering during the course of a trajectory to handle isolated modes, and short-cut methods that prevent useless trajectories from taking much computation time.
1 Introduction
The introduction traces MCMC from molecular-state simulation to deterministic molecular dynamics and identifies the 1987 work that united these approaches as “Hybrid Monte Carlo.”
- Origins: MCMC originated in Metropolis et al.’s 1953 simulation of states for idealized molecules.
- Origins: Alder and Wainwright introduced deterministic molecular simulation in 1959, with molecular motion following Newton’s laws.
- Hybrid Monte Carlo: In 1987, Duane, Kennedy, Pendleton, and Roweth united MCMC with molecular dynamics and called their method “Hybrid Monte Carlo”.
2 MCMC USING HAMILTONIAN DYNAMICS
This section introduces Hamiltonian Monte Carlo by combining position variables with auxiliary Gaussian momentum variables and discretized Hamiltonian dynamics. It then develops practical tuning guidance, scaling results, and variations such as windowed acceptance decisions and approximate trajectories.
- Hamiltonian dynamics: Hamiltonian dynamics must be discretized for computer implementation, typically using the simple leapfrog scheme.The leapfrog integrator provides the standard computational approximation to the continuous dynamics.
- Basic method: HMC defines a Hamiltonian from the target distribution, augments position variables with typically independent Gaussian momentum variables, and alternates simple updates.The method constructs an MCMC procedure from the desired sampling distribution using auxiliary momentum variables.
- Practical and theoretical issues: The review addresses tuning the leapfrog stepsize and number of steps, alongside theoretical results on how HMC scales with dimensionality.These practical and theoretical issues follow the presentation of basic HMC.
- HMC variations: It also presents HMC variations that use windows of states to increase acceptance rates and approximate trajectories for statistical problems.Windowing considers states at the beginning and end of a trajectory, while approximate computation can use subsets of the data.
2 Hamiltonian dynamics
Hamiltonian dynamics augments position variables with momenta and evolves them through a Hamiltonian over a 2d-dimensional state space. Its reversibility, volume preservation, and exact preservation under suitable discretizations make it useful for Metropolis updates without Jacobian corrections.
- Hamiltonian formulation: Hamiltonian dynamics evolves a d-dimensional position q and d-dimensional momentum p using a Hamiltonian H(q, p), yielding a 2d-dimensional state space.In MCMC, q represents variables of interest, potential energy is minus their log density, and momenta are introduced artificially.
- Hamiltonian formulation: For separable Hamiltonians, potential energy is minus the log target density, while quadratic kinetic energy corresponds to a zero-mean Gaussian momentum with covariance M.M is typically a diagonal symmetric positive-definite mass matrix or a scalar multiple of the identity.
- Dynamical properties: Hamiltonian trajectories are reversible and preserve the Hamiltonian, while volume preservation removes the need to include a Jacobian determinant in Metropolis acceptance probabilities.For arbitrary non-Hamiltonian dynamics, computing the Jacobian determinant may be infeasible.
- Dynamical properties: Hamiltonian dynamics preserves volume exactly because det(B_s) = 1, and reversibility, volume preservation, and symplecticness can remain exact after approximation.The volume-preservation result follows by taking the limit of bounded local determinant changes along nonsingular trajectories.
- Discretization: The leapfrog method is volume-preserving and reversible, producing stable trajectories where Euler’s method diverges, while its global error is order ε^2 versus order ε for Euler.With ε = 0.3, leapfrog results are indistinguishable from the true trajectory at plot scale; with ε = 1.2, approximation error is visible but stability remains.
3 MCMC from Hamiltonian dynamics
HMC augments continuous target variables with fictitious momentum and uses Hamiltonian-dynamics proposals within a Metropolis algorithm. Momentum resampling and the Metropolis update preserve the canonical joint distribution, while systematic trajectories can explore more efficiently than random-walk proposals.
- Hamiltonian Monte Carlo construction: HMC represents target variables as position q, introduces independent momentum p, and alternates momentum resampling with a Hamiltonian-dynamics Metropolis update.The target density is encoded through potential energy U(q), while kinetic energy K(p) specifies the momentum distribution.
- Hamiltonian Monte Carlo construction: HMC applies to continuous distributions whose density and partial derivatives of its log can be evaluated, initially assuming the density is everywhere non-zero.The density need only be available up to an unknown normalizing constant.
- Correctness: Both momentum resampling and the Metropolis update leave the canonical joint distribution invariant, so the complete HMC iteration does as well.This invariance follows because the two component steps each preserve the same canonical distribution.
- Correctness: HMC is typically ergodic, but exact periodicity produced by the leapfrog trajectory can cause ergodicity to fail.When ergodic, the chain asymptotically converges to its unique invariant distribution.
- Efficiency: Systematic motion during HMC trajectories produces larger state changes than comparable random-walk Metropolis iterations, with movement growing roughly linearly for n steps before reaching distributional width.The efficiency advantage is roughly tied to the ratio of standard deviations in the least confined and most confined directions.
- Efficiency: For 100 variables, autocorrelation was clearly higher for random-walk Metropolis than for HMC when computation time per iteration was made comparable.The comparison counted 150 random-walk Metropolis updates as one iteration.
4 HMC in practice and theory
HMC’s practical performance depends critically on choosing momentum geometry, leapfrog stepsize, and trajectory length, while its computational cost scales more favorably with dimension than random-walk Metropolis. Theoretical analysis explains both tuning constraints and this dimensionality advantage.
- Momentum and variable transformations: HMC performance is unchanged by invertible linear transformations of position and momentum variables when the transformed dynamics replicates the original dynamics.Orthogonal transformations preserve performance directly; rotationally symmetric momentum with M = mI also leaves the kinetic energy and momentum distribution unchanged.
- Momentum and variable transformations: Matrix-based kinetic energy can compensate for correlations, but its computational cost limits usefulness in high dimensions; diagonal covariance captures scales but not correlations.Choosing optimal diagonal scales is difficult when correlations are nonzero, though conditional-standard-deviation approximations may be possible.
- Tuning HMC: Selecting leapfrog stepsize ε and step count L is a major practical impediment, because unstable or overly large ε can make H grow exponentially and acceptance probability extremely small.HMC is more sensitive than random-walk Metropolis to an excessively large tuning value, with sharp degradation when ε exceeds the stability limit; stability can vary across state-space regions.
- Tuning HMC: Trajectory length must be long enough to move far along weakly constrained directions, but excessive length can undermine systematic exploration.The choice of trajectory length is crucial for avoiding random-walk-like exploration in distributions with directions of widely differing constraints.
- Scaling with dimensionality: HMC avoids inefficient random-walk exploration and, for moderate-to-high dimensions, computation time increases more slowly with dimension for a given accuracy.Random-walk Metropolis requires at least order d updates to explore the distribution, reinforcing HMC’s scaling advantage.
- Scaling with dimensionality: ε must scale as d^-1/4 to maintain reasonable acceptance, making leapfrog updates scale as d1/4 and computation time typically as d5/4, versus d2 for random-walk Metropolis.For fixed ε and trajectory length, leapfrog energy-error mean and variance grow linearly with d; maintaining acceptance requires keeping mean total energy error near one.
5 Extensions and variations on HMC
This section surveys HMC modifications that improve efficiency or extend applicability, including alternative discretizations, approximate dynamics, momentum refreshment, windowed acceptance, tempering, and computational shortcuts. These variations preserve correctness under specified constructions while addressing scaling, tuning, or trajectory-cost limitations.
- Alternative discretizations: Composing exact dynamics for component Hamiltonians yields a valid discretization that converges to the target dynamics with global error O(ε) or less.The construction assumes the component Hamiltonians are twice differentiable and reproduces exact dynamics as ε approaches zero.
- Approximate dynamics: When U0 is cheaper than U1, using M substeps for U0 permits larger ε, fewer trajectory steps, and fewer evaluations of U1.The effective U0 stepsize is ε/M; practical gains are smaller because terminal Hamiltonian error exceeds that of standard leapfrog, though neural-network models may benefit.
- One-step HMC and LMC: For LMC, maintaining acceptance as dimension d grows requires ε to scale as d^-1/6, while reaching a nearly independent point takes computation proportional to d^1/3.The result follows from average Hamiltonian error of order dε^6 and random-walk exploration.
- One-step HMC and LMC: The optimal LMC acceptance target is 0.57 when its behavior resembles sampling replicated independent variables.This value is obtained numerically and agrees with a formal result by Roberts and Rosenthal.
- Partial momentum refreshment: Partial momentum refreshment preserves the canonical distribution but gives only a small constant-factor improvement, with no better dimensional scaling, for high-dimensional Gaussian targets.Best performance uses long trajectories; moderate trajectories with α > 0 may become more effective when combined with windowed HMC, which can reduce rejection rates.
- Windowed HMC: Windowed HMC preserves the correct distribution, improves performance by a factor of two or more in anisotropic Gaussian examples, and avoids a stated ergodicity problem when W > 1.Averaging densities across trajectory windows smooths Hamiltonian oscillations caused by inexact simulation.
- Tuning and computational shortcuts: HMC performance depends critically on tuning ε and L, while early failure checks can limit wasted computation when an inappropriate stepsize is detected.The optimal trajectory length εL depends on the distribution’s global extent, and the shortcut behaves like standard HMC when grouped tests do not fail.
- Tempering and further extensions: Tempered trajectories retain a unit Jacobian determinant, so their endpoints can be accepted as proposals without adding a Jacobian factor.A proposed unexplored extension would split the exact potential into an approximation and its error, using splitting techniques to exploit tractability or low computational cost.