Hamiltonian Monte Carlo
In computational physics and statistics, the Hamiltonian Monte Carlo algorithm (also known as hybrid Monte Carlo), is a Markov chain Monte Carlo method for obtaining a sequence of random samples which converge to being distributed according to a target probability distribution for which direct sampling is difficult. This sequence can be used to estimate integrals with respect to the target distribution (expected values).
Hamiltonian Monte Carlo corresponds to an instance of the Metropolis–Hastings algorithm, with a Hamiltonian dynamics evolution simulated using a time-reversible and volume-preserving numerical integrator (typically the leapfrog integrator) to propose a move to a new point in the state space. Compared to using a Gaussian random walk proposal distribution in the Metropolis–Hastings algorithm, Hamiltonian Monte Carlo reduces the correlation between successive sampled states by proposing moves to distant states which maintain a high probability of acceptance due to the approximate energy conserving properties of the simulated Hamiltonian dynamic when using a symplectic integrator. The reduced correlation means fewer Markov chain samples are needed to approximate integrals with respect to the target probability distribution for a given Monte Carlo error. The algorithm was originally proposed by Simon Duane, Anthony Kennedy, Brian Pendleton and Duncan Roweth in 1987 for calculations in lattice quantum chromodynamics.
Suppose the target distribution to sample is and a chain of samples is required. The Hamilton's equations are
where is the potential energy. The potential energy for a target is given as
which comes from the Boltzmann's factor.
The algorithm requires a positive integer for number of leap frog steps and a positive number for the step size . Suppose the chain is at . Let . First, a random Gaussian momentum is drawn from . Next, the particle will run under Hamiltonian dynamics for time , this is done by solving the Hamilton's equations numerically using the leap frog algorithm. The position and momentum vectors after time using the leap frog algorithm are
These equations are to be applied to and times to obtain and .
The leap frog algorithm is an approximate solution to the motion of non-interacting classical particles. If exact, the solution will never change the initial randomly-generated energy distribution, as energy is conserved for each particle in the presence of a classical potential energy field. In order to reach a thermodynamic equilibrium distribution, particles must have some sort of interaction with, for example, a surrounding heat bath, so that the entire system can take on different energies with probabilities according to the Boltzmann distribution.
One way to move the system towards a thermodynamic equilibrium distribution is to change the state of the particles using the Metropolis–Hastings. So first, one applies the leap frog step, then a Metropolis-Hastings step.
The transition from to is
This is repeated to obtain .
No U-Turn Sampler
The No U-Turn Sampler (NUTS) is an extension by controlling automatically. Tuning is critical. For example in the one dimensional case, the potential is which corresponds to the potential of a simple harmonic oscillator. For too large, the particle will oscillate and this waste computational time. For too small, the particle will behave like a random walk.
Loosely, NUTS runs the Hamiltonian dynamics both forwards and backwards in time randomly until a U-Turn condition is satisfied. When that happens, a random point from the path is chosen for the MCMC sample and the process is repeated from that new point.
In detail, a binary tree is constructed to trace the path of the leap frog steps. To produce a MCMC sample, an iterative procedure is conducted. A slice variable is sampled. Let and be the position and momentum of the forward particle respectively. Similarly, and for the backward particle. In each iteration, the binary tree selects at random uniformly to move the forward particle forwards in time or the backward particle backwards in time. Also for each iteration, the number of leap frog steps increase by a factor of 2. For example, in the first iteration, the forward particle moves forwards in time using 1 leap frog step. In the next iteration, the backward particle moves backwards in time using 2 leap frog steps.
The iterative procedure continues until the U-Turn condition is met, that is
or when the Hamiltonian becomes inaccurate
where, for example, .
Once the U-Turn condition is met, the next MCMC sample, , is obtained by sampling uniformly the leap frog path traced out by the binary tree which satisfies
This is usually satisfied if the remaining HMC parameters are sensible.
- Dynamic Monte Carlo method
- Software for Monte Carlo molecular modeling
- Metropolis-adjusted Langevin algorithm
- Duane, Simon; Kennedy, Anthony D.; Pendleton, Brian J.; Roweth, Duncan (3 September 1987). "Hybrid Monte Carlo". Physics Letters B. 195 (2): 216–222. Bibcode:1987PhLB..195..216D. doi:10.1016/0370-2693(87)91197-X.
- Hoffman, Matthew D; Gelman, Andrew (2014). "The No-U-turn sampler: adaptively setting path lengths in Hamiltonian Monte Carlo". Journal of Machine Learning Research. 15 (1): 1593-1623.
- Neal, Radford M (2011). "MCMC Using Hamiltonian Dynamics" (PDF). In Steve Brooks; Andrew Gelman; Galin L. Jones; Xiao-Li Meng (eds.). Handbook of Markov Chain Monte Carlo. Chapman and Hall/CRC. ISBN 9781420079418.
- Betancourt, Michael (2018). "A Conceptual Introduction to Hamiltonian Monte Carlo". arXiv:1701.02434. Bibcode:2017arXiv170102434B. Cite journal requires
- Liu, Jun S. (2004). Monte Carlo Strategies in Scientific Computing. Springer Series in Statistics, Springer. pp. 189-203. ISBN 978-0-387-76369-9.
- Betancourt, Michael. "Efficient Bayesian inference with Hamiltonian Monte Carlo". MLSS Iceland 2014 – via YouTube.
- Hamiltonian Monte Carlo from scratch
- Optimization and Monte Carlo Methods