Markov chain Monte Carlo

In statistics, Markov chain Monte Carlo (MCMC) methods are a class of algorithms for sampling from a probability distribution based on constructing a Markov chain that has the desired distribution as its equilibrium distribution. The state of the chain after a number of steps is then used as a sample of the desired distribution. The quality of the sample improves as a function of the number of steps.

Convergence of the Metropolis-Hastings algorithm. MCMC attempts to approximate the blue distribution with the orange distribution

Random walk Monte Carlo methods make up a large subclass of MCMC methods.

Application domains

• MCMC methods are primarily used for calculating numerical approximations of multi-dimensional integrals, for example in Bayesian statistics, computational physics, computational biology and computational linguistics.^[1][2]
• In Bayesian statistics, the recent development of MCMC methods has been a key step in making it possible to compute large hierarchical models that require integrations over hundreds or even thousands of unknown parameters.^[3]
• They are also used for generating samples that gradually populate the rare failure region in rare event sampling.

Classification

Random walk Monte Carlo methods

Multi-dimensional integrals

When an MCMC method is used for approximating a multi-dimensional integral, an ensemble of "walkers" move around randomly. At each point where a walker steps, the integrand value at that point is counted towards the integral. The walker then may make a number of tentative steps around the area, looking for a place with a reasonably high contribution to the integral to move into next.

Random walk Monte Carlo methods are a kind of random simulation or Monte Carlo method. However, whereas the random samples of the integrand used in a conventional Monte Carlo integration are statistically independent, those used in MCMC methods are correlated. A Markov chain is constructed in such a way as to have the integrand as its equilibrium distribution.

Examples

Examples of random walk Monte Carlo methods include the following:

• Metropolis–Hastings algorithm: This method generates a random walk using a proposal density and a method for rejecting some of the proposed moves.
• Gibbs sampling: This method requires all the conditional distributions of the target distribution to be sampled exactly. When drawing from the full-conditional distributions is not straightforward other samplers-within-Gibbs are used (e.g., see ^[4][5][6]). Gibbs sampling is popular partly because it does not require any 'tuning'.
• Slice sampling: This method depends on the principle that one can sample from a distribution by sampling uniformly from the region under the plot of its density function. It alternates uniform sampling in the vertical direction with uniform sampling from the horizontal 'slice' defined by the current vertical position.
• Multiple-try Metropolis: This method is a variation of the Metropolis–Hastings algorithm that allows multiple trials at each point. By making it possible to take larger steps at each iteration, it helps address the curse of dimensionality.^[7][8]
• Reversible-jump: This method is a variant of the Metropolis–Hastings algorithm that allows proposals that change the dimensionality of the space.^[9] MCMC methods that change dimensionality have long been used in statistical physics applications, where for some problems a distribution that is a grand canonical ensemble is used (e.g., when the number of molecules in a box is variable). But the reversible-jump variant is useful when doing MCMC or Gibbs sampling over nonparametric Bayesian models such as those involving the Dirichlet process or Chinese restaurant process, where the number of mixing components/clusters/etc. is automatically inferred from the data.

Other MCMC methods

Interacting MCMC methodologies are a class of mean field particle methods for obtaining random samples from a sequence of probability distributions with an increasing level of sampling complexity.^[10] These probabilistic models include path space state models with increasing time horizon, posterior distributions w.r.t. sequence of partial observations, increasing constraint level sets for conditional distributions, decreasing temperature schedules associated with some Boltzmann-Gibbs distributions, and many others. In principle, any MCMC sampler can be turned into an interacting MCMC sampler. These interacting MCMC samplers can be interpreted as a way to run in parallel a sequence of MCMC samplers. For instance, interacting simulated annealing algorithms are based on independent Metropolis-Hastings moves interacting sequentially with a selection-resampling type mechanism. In contrast to traditional MCMC methods, the precision parameter of this class of interacting MCMC samplers is only related to the number of interacting MCMC samplers. These advanced particle methodologies belong to the class of Feynman-Kac particle models,^[11][12] also called Sequential Monte Carlo or particle filter methods in Bayesian inference and signal processing communities.^[13] Interacting MCMC methods can also be interpreted as a mutation-selection genetic particle algorithm with MCMC mutations.

Markov Chain quasi-Monte Carlo (MCQMC)^[14][15] The advantage of low-discrepancy sequences in lieu of random numbers for simple independent Monte Carlo sampling is well known.^[16] This procedure, known as Quasi-Monte Carlo method (QMC),^[17] yields an integration error that decays at a superior rate to that obtained by IID sampling, by the Koksma-Hlawka inequality. Empirically it allows to reduce both estimation error and convergence time by an order of magnitude.

Reducing correlation

More sophisticated methods use various ways of reducing the correlation between successive samples. These algorithms may be harder to implement, but they usually exhibit faster convergence (i.e. fewer steps for an accurate result).

Examples

Examples of non-random walk MCMC methods include the following:

• Hybrid Monte Carlo (HMC): Tries to avoid random walk behaviour by introducing an auxiliary momentum vector and implementing Hamiltonian dynamics, so the potential energy function is the target density. The momentum samples are discarded after sampling. The end result of Hybrid Monte Carlo is that proposals move across the sample space in larger steps; they are therefore less correlated and converge to the target distribution more rapidly.
• Some variations on slice sampling also avoid random walks.^[18]
• Langevin MCMC and other methods that rely on the gradient (and possibly second derivative) of the log posterior avoid random walks by making proposals that are more likely to be in the direction of higher probability density.^[19]

Convergence

Usually it is not hard to construct a Markov chain with the desired properties. The more difficult problem is to determine how many steps are needed to converge to the stationary distribution within an acceptable error. A good chain will have rapid mixing: the stationary distribution is reached quickly starting from an arbitrary position.

Typically, MCMC sampling can only approximate the target distribution, as there is always some residual effect of the starting position. More sophisticated MCMC-based algorithms such as coupling from the past can produce exact samples, at the cost of additional computation and an unbounded (though finite in expectation) running time.

Many random walk Monte Carlo methods move around the equilibrium distribution in relatively small steps, with no tendency for the steps to proceed in the same direction. These methods are easy to implement and analyze, but unfortunately it can take a long time for the walker to explore all of the space. The walker will often double back and cover ground already covered.

See also

• Bayesian inference
• Bayesian network
• Coupling from the past
• Gibbs sampling
• Quasi-Monte Carlo method
• Hybrid Monte Carlo
• Metropolis-adjusted Langevin algorithm
• Metropolis–Hastings algorithm
• Multiple-try Metropolis
• Particle filter
• Reversible-jump
• Slice sampling

Notes

1. See Gill 2008.
2. See Robert & Casella 2004.
3. Banerjee, Sudipto; Carlin, Bradley P.; Gelfand, Alan P. Hierarchical Modeling and Analysis for Spatial Data (Second ed.). CRC Press. p. xix. ISBN 978-1-4398-1917-3.
4. Gilks, W. R.; Wild, P. (1992-01-01). "Adaptive Rejection Sampling for Gibbs Sampling". Journal of the Royal Statistical Society. Series C (Applied Statistics) 41 (2): 337–348. doi:10.2307/2347565. JSTOR 2347565.
5. Gilks, W. R.; Best, N. G.; Tan, K. K. C. (1995-01-01). "Adaptive Rejection Metropolis Sampling within Gibbs Sampling". Journal of the Royal Statistical Society. Series C (Applied Statistics) 44 (4): 455–472. doi:10.2307/2986138. JSTOR 2986138.
6. Martino, L.; Read, J.; Luengo, D. (2015-06-01). "Independent Doubly Adaptive Rejection Metropolis Sampling Within Gibbs Sampling". IEEE Transactions on Signal Processing 63 (12): 3123–3138. doi:10.1109/TSP.2015.2420537. ISSN 1053-587X.
7. Liu, Jun S.; Liang, Faming; Wong, Wing Hung (2000-03-01). "The Multiple-Try Method and Local Optimization in Metropolis Sampling". Journal of the American Statistical Association 95 (449): 121–134. doi:10.1080/01621459.2000.10473908. ISSN 0162-1459.
8. Martino, Luca; Read, Jesse (2013-07-11). "On the flexibility of the design of multiple try Metropolis schemes". Computational Statistics 28 (6): 2797–2823. doi:10.1007/s00180-013-0429-2. ISSN 0943-4062.
9. See Green 1995.
10. Del Moral, Pierre (2013). Mean field simulation for Monte Carlo integration. Chapman & Hall/CRC Press. p. 626. Monographs on Statistics & Applied Probability
11. Del Moral, Pierre (2004). Feynman-Kac formulae. Genealogical and interacting particle approximations. Springer. p. 575. Series: Probability and Applications
12. Del Moral, Pierre; Miclo, Laurent (2000). Branching and Interacting Particle Systems Approximations of Feynman-Kac Formulae with Applications to Non-Linear Filtering. (PDF). Lecture Notes in Mathematics 1729. pp. 1–145. doi:10.1007/bfb0103798.
13. "Sequential Monte Carlo samplers - P. Del Moral - A. Doucet - A. Jasra - 2006 - Journal of the Royal Statistical Society: Series B (Statistical Methodology) - Wiley Online Library". onlinelibrary.wiley.com. Retrieved 2015-06-11.
14. Chen, S., Josef Dick, and Art B. Owen. "Consistency of Markov chain quasi-Monte Carlo on continuous state spaces." The Annals of Statistics 39.2 (2011): 673-701.
15. Tribble, Seth D. Markov chain Monte Carlo algorithms using completely uniformly distributed driving sequences. Diss. Stanford University, 2007.
16. Papageorgiou, Anargyros, and J. F. Traub. "Beating Monte Carlo." Risk 9.6 (1996): 63-65.
17. Sobol, Ilya M. "On quasi-monte carlo integrations." Mathematics and Computers in Simulation 47.2 (1998): 103-112.
18. See Neal 2003.
19. See Stramer 1999.

References

• Christophe Andrieu, Nando De Freitas and Arnaud Doucet, An Introduction to MCMC for Machine Learning, 2003
• Asmussen, Søren; Glynn, Peter W. (2007). Stochastic Simulation: Algorithms and Analysis. Stochastic Modelling and Applied Probability 57. Springer.
• Atzberger, P. "An Introduction to Monte-Carlo Methods" (PDF).
• Berg, Bernd A. (2004). Markov Chain Monte Carlo Simulations and Their Statistical Analysis. World Scientific.
• Bolstad, William M. (2010). Understanding Computational Bayesian Statistics. Wiley. ISBN 0-470-04609-0.
• Casella, George; George, Edward I. (1992). "Explaining the Gibbs sampler". The American Statistician 46: 167–174. doi:10.2307/2685208. (Basic summary and many references.)
• Gelfand, A.E.; Smith, A.F.M. (1990). "Sampling-Based Approaches to Calculating Marginal Densities". Journal of the American Statistical Association 85: 398–409. doi:10.1080/01621459.1990.10476213.
• Gelman, Andrew; Carlin, John B.; Stern, Hal S.; Rubin, Donald B. (1995). Bayesian Data Analysis (1st ed.). Chapman and Hall. (See Chapter 11.)
• Geman, S.; Geman, D. (1984). "Stochastic Relaxation, Gibbs Distributions, and the Bayesian Restoration of Images". IEEE Transactions on Pattern Analysis and Machine Intelligence 6: 721–741.
• Gilks, W.R.; Richardson, S.; Spiegelhalter, D.J. (1996). Markov Chain Monte Carlo in Practice. Chapman and Hall/CRC.
• Gill, Jeff (2008). Bayesian methods: a social and behavioral sciences approach (2nd ed.). Chapman and Hall/CRC. ISBN 1-58488-562-9.
• Green, P.J. (1995). "Reversible-jump Markov chain Monte Carlo computation and Bayesian model determination". Biometrika 82 (4): 711–732. doi:10.1093/biomet/82.4.711.
• Neal, Radford M. (2003). "Slice Sampling". Annals of Statistics 31 (3): 705–767. doi:10.1214/aos/1056562461. JSTOR 3448413.
• Neal, Radford M. (1993). "Probabilistic Inference Using Markov Chain Monte Carlo Methods".
• Robert, Christian P.; Casella, G. (2004). Monte Carlo Statistical Methods (2nd ed.). Springer. ISBN 0-387-21239-6.
• Rubinstein, R.Y.; Kroese, D.P. (2007). Simulation and the Monte Carlo Method (2nd ed.). Wiley. ISBN 978-0-470-17794-5.
• Smith, R.L. (1984). "Efficient Monte Carlo Procedures for Generating Points Uniformly Distributed Over Bounded Regions". Operations Research 32: 1296–1308. doi:10.1287/opre.32.6.1296.
• Spall, J.C. (April 2003). "Estimation via Markov Chain Monte Carlo". IEEE Control Systems Magazine 23 (2): 34–45. doi:10.1109/mcs.2003.1188770.
• Stramer, O.; Tweedie, R. (1999). "Langevin-Type Models II: Self-Targeting Candidates for MCMC Algorithms". Methodology and Computing in Applied Probability 1 (3): 307–328. doi:10.1023/A:1010090512027.

Further reading

• Diaconis, Persi (April 2009). "The Markov chain Monte Carlo revolution" (PDF). Bull. Amer. Math. Soc. 46 (2): 179–205. doi:10.1090/s0273-0979-08-01238-x. S 0273-0979(08)01238-X.
• Press, W.H.; Teukolsky, S.A.; Vetterling, W.T.; Flannery, B.P. (2007), "Section 15.8. Markov Chain Monte Carlo", Numerical Recipes: The Art of Scientific Computing (3rd ed.), Cambridge University Press, ISBN 978-0-521-88068-8
• Richey, Matthew (May 2010). "The Evolution of Markov Chain Monte Carlo Methods" (PDF). The American Mathematical Monthly 117 (5): 383–413. doi:10.4169/000298910x485923.

External links

• MCMC sampling and other methods in a basic overview, by Alexander Mantzaris (original link - now broken)
• Interacting Metropolis-Hastings methodologies
• Visual demonstration of MCMC sampling methods (Java applet), by Laird Breyer
• A Toy Example of MCMC sampling, by Zhiyuan Weng
• MCL - a cluster algorithm for graphs, by Stijn van Dongen
• PyMC - Python module implementing Bayesian statistical models and fitting algorithms, including Markov chain Monte Carlo.
• IA2RMS is a Matlab code of the Independent Doubly Adaptive Rejection Metropolis Sampling method^[1] for drawing from the full-conditional densities within a Gibbs sampler.