En estadística , el algoritmo Monte Carlo de cadena de Markov ( MCMC ) es una clase de algoritmos que se utilizan para extraer muestras de una distribución de probabilidad . Dada una distribución de probabilidad, se puede construir una cadena de Markov cuya distribución de elementos se aproxime a ella, es decir, la distribución de equilibrio de la cadena de Markov coincide con la distribución objetivo. Cuantos más pasos se incluyan, más se ajustará la distribución de la muestra a la distribución deseada real.
Los métodos de Monte Carlo de cadenas de Markov se utilizan para estudiar distribuciones de probabilidad que son demasiado complejas o demasiado dimensionales para estudiarlas únicamente con técnicas analíticas. Existen varios algoritmos para construir dichas cadenas de Markov, incluido el algoritmo Metropolis-Hastings .
En las estadísticas bayesianas, los métodos de Monte Carlo de cadenas de Markov se utilizan normalmente para calcular momentos e intervalos creíbles de distribuciones de probabilidad posteriores . El uso de métodos MCMC permite calcular grandes modelos jerárquicos que requieren integraciones sobre cientos o miles de parámetros desconocidos. [5]
Los métodos de Monte Carlo de cadena de Markov crean muestras de una variable aleatoria continua , con una densidad de probabilidad proporcional a una función conocida. Estas muestras se pueden utilizar para evaluar una integral sobre esa variable, como su valor esperado o varianza .
En la práctica, se suele desarrollar un conjunto de cadenas, a partir de un conjunto de puntos elegidos arbitrariamente y suficientemente distantes entre sí. Estas cadenas son procesos estocásticos de "caminantes" que se desplazan aleatoriamente según un algoritmo que busca lugares con una contribución razonablemente alta a la integral para pasar a continuación, asignándoles mayores probabilidades.
Aunque los métodos MCMC se crearon para abordar problemas multidimensionales mejor que los algoritmos genéricos de Monte Carlo, cuando el número de dimensiones aumenta, también tienden a sufrir la maldición de la dimensionalidad : las regiones de mayor probabilidad tienden a estirarse y perderse en un volumen creciente de espacio que contribuye poco a la integral. Una forma de abordar este problema podría ser acortar los pasos del caminante, de modo que no intente continuamente salir de la región de mayor probabilidad, aunque de esta manera el proceso estaría altamente autocorrelacionado y sería costoso (es decir, se requerirían muchos pasos para obtener un resultado preciso). Métodos más sofisticados como el Hamiltoniano Monte Carlo y el algoritmo de Wang y Landau utilizan varias formas de reducir esta autocorrelación, al tiempo que logran mantener el proceso en las regiones que dan una mayor contribución a la integral. Estos algoritmos generalmente se basan en una teoría más complicada y son más difíciles de implementar, pero generalmente convergen más rápido.
Ejemplos
Paseo aleatorio
Algoritmo Metropolis-Hastings : Este método genera una cadena de Markov utilizando una densidad de propuestas para los nuevos pasos y un método para rechazar algunos de los movimientos propuestos. En realidad, se trata de un marco general que incluye como casos especiales el primer y más simple algoritmo MCMC (algoritmo Metropolis) y muchas alternativas más recientes que se enumeran a continuación.
Muestreo de Gibbs : cuando la distribución objetivo es multidimensional, el algoritmo de muestreo de Gibbs [6] actualiza cada coordenada a partir de su distribución condicional completa dadas otras coordenadas. El muestreo de Gibbs puede considerarse un caso especial del algoritmo Metropolis-Hastings con una tasa de aceptación uniformemente igual a 1. Cuando la extracción de las distribuciones condicionales completas no es sencilla, se utilizan otros muestreadores dentro de Gibbs (por ejemplo, consulte [7] [8] ). El muestreo de Gibbs es popular en parte porque no requiere ningún "ajuste". La estructura del algoritmo del muestreo de Gibbs se parece mucho a la de la inferencia variacional de ascenso de coordenadas en que ambos algoritmos utilizan las distribuciones condicionales completas en el procedimiento de actualización. [9]
Algoritmo de Langevin ajustado a Metropolis y otros métodos que se basan en el gradiente (y posiblemente la segunda derivada) de la densidad objetivo logarítmica para proponer pasos que tienen más probabilidades de estar en la dirección de una mayor densidad de probabilidad. [10]
Hamiltoniano (o híbrido) Monte Carlo (HMC): intenta evitar el comportamiento de caminata aleatoria introduciendo un vector de momento auxiliar e implementando la dinámica hamiltoniana , de modo que la función de energía potencial sea la densidad objetivo. Las muestras de momento se descartan después del muestreo. El resultado del Monte Carlo híbrido es que las propuestas se mueven a través del espacio muestral en pasos más grandes; por lo tanto, están menos correlacionadas y convergen a la distribución objetivo más rápidamente.
Metrópolis pseudomarginal-Hastings: este método reemplaza la evaluación de la densidad de la distribución objetivo con una estimación imparcial y es útil cuando la densidad objetivo no está disponible analíticamente, por ejemplo, modelos de variables latentes .
Muestreo por cortes : este método se basa en el principio de que se puede tomar una muestra de una distribución tomando muestras de manera uniforme de la región que se encuentra debajo del gráfico de su función de densidad. Alterna el muestreo uniforme en dirección vertical con el muestreo uniforme de la "porción" horizontal definida por la posición vertical actual.
Metropolis de múltiples intentos : este método es una variación del algoritmo Metropolis-Hastings que permite múltiples intentos en cada punto. Al permitir realizar pasos más grandes en cada iteración, ayuda a abordar el problema de la dimensionalidad.
Salto reversible : este método es una variante del algoritmo Metropolis-Hastings que permite propuestas que cambian la dimensionalidad del espacio. [11] Los métodos de Monte Carlo de cadena de Markov que cambian la dimensionalidad se han utilizado durante mucho tiempo en aplicaciones de física estadística , donde para algunos problemas se utiliza una distribución que es un gran conjunto canónico (por ejemplo, cuando el número de moléculas en una caja es variable). Pero la variante de salto reversible es útil cuando se hace un muestreo de Monte Carlo de cadena de Markov o de Gibbs sobre modelos bayesianos no paramétricos como los que involucran el proceso de Dirichlet o el proceso del restaurante chino , donde el número de componentes/grupos/etc. de mezcla se infiere automáticamente a partir de los datos.
Métodos de partículas interactuantes
Las metodologías MCMC interactuantes son una clase de métodos de partículas de campo medio para obtener muestras aleatorias de una secuencia de distribuciones de probabilidad con un nivel creciente de complejidad de muestreo. [12] Estos modelos probabilísticos incluyen modelos de estados del espacio de trayectoria con horizonte temporal creciente, distribuciones posteriores con respecto a la secuencia de observaciones parciales, conjuntos de niveles de restricción crecientes para distribuciones condicionales, programas de temperatura decrecientes asociados con algunas distribuciones de Boltzmann-Gibbs y muchos otros. En principio, cualquier muestreador de Monte Carlo de cadena de Markov se puede convertir en un muestreador de Monte Carlo de cadena de Markov interactuante. Estos muestreadores de Monte Carlo de cadena de Markov interactuantes se pueden interpretar como una forma de ejecutar en paralelo una secuencia de muestreadores de Monte Carlo de cadena de Markov. Por ejemplo, los algoritmos de recocido simulado interactuantes se basan en movimientos Metropolis-Hastings independientes que interactúan secuencialmente con un mecanismo de tipo selección-remuestreo. A diferencia de los métodos tradicionales de Monte Carlo de cadena de Markov, el parámetro de precisión de esta clase de muestreadores de Monte Carlo de cadena de Markov interactuantes solo está relacionado con el número de muestreadores de Monte Carlo de cadena de Markov interactuantes. Estas metodologías de partículas avanzadas pertenecen a la clase de modelos de partículas de Feynman-Kac, [13] [14] también llamados métodos de Monte Carlo secuencial o de filtro de partículas en las comunidades de inferencia bayesiana y procesamiento de señales . [15] Los métodos de Monte Carlo de cadena de Markov interactuantes también se pueden interpretar como un algoritmo de partículas genéticas de selección de mutaciones con mutaciones de Monte Carlo de cadena de Markov.
Cuasi-Monte Carlo
El método cuasi-Monte Carlo es un análogo del método Monte Carlo normal que utiliza secuencias de baja discrepancia en lugar de números aleatorios. [16] [17] Produce un error de integración que decae más rápido que el del muestreo aleatorio verdadero, como se cuantifica mediante la desigualdad de Koksma-Hlawka . Empíricamente, permite la reducción tanto del error de estimación como del tiempo de convergencia en un orden de magnitud. [ cita requerida ] Los métodos cuasi-Monte Carlo de cadena de Markov [18] [19] como el método Array-RQMC combinan la simulación aleatoria cuasi-Monte Carlo y de cadena de Markov simulando cadenas simultáneamente de una manera que se aproxima mejor a la distribución real de la cadena que con el MCMC ordinario. [20] En experimentos empíricos, la varianza del promedio de una función del estado a veces converge a la tasa o incluso más rápido, en lugar de la tasa de Monte Carlo. [21]
Convergencia
Por lo general, no es difícil construir una cadena de Markov con las propiedades deseadas. El problema más difícil es determinar cuántos pasos se necesitan para converger a la distribución estacionaria dentro de un error aceptable. [22] Una buena cadena tendrá una mezcla rápida : la distribución estacionaria se alcanza rápidamente a partir de una posición arbitraria. Un método empírico estándar para evaluar la convergencia es ejecutar varias cadenas de Markov simuladas independientes y verificar que la relación entre las varianzas intercadena e intracadena para todos los parámetros muestreados sea cercana a 1. [22] [23]
Muchos métodos de Monte Carlo de paseo aleatorio se mueven alrededor de la distribución de equilibrio en pasos relativamente pequeños, sin que los pasos tiendan a avanzar en la misma dirección. Estos métodos son fáciles de implementar y analizar, pero desafortunadamente puede llevar mucho tiempo para que el caminante explore todo el espacio. El caminante a menudo retrocederá y cubrirá el terreno ya cubierto.
Se puede encontrar una consideración más amplia de la convergencia en el teorema del límite central de la cadena de Markov . Véase [24] para una discusión de la teoría relacionada con la convergencia y la estacionariedad del algoritmo Metropolis-Hastings.
Software
Varios programas de software proporcionan capacidades de muestreo MCMC, por ejemplo:
El software Monte Carlo paralelo ParaMonte está disponible en varios lenguajes de programación, incluidos C , C++ , Fortran , MATLAB y Python .
Paquetes que utilizan dialectos del lenguaje del modelo BUGS :
^ Kasim, MF; Bott, AFA; Tzeferacos, P.; Lamb, DQ; Gregori, G.; Vinko, SM (septiembre de 2019). "Recuperación de campos de radiografía de protones sin perfiles de fuente". Physical Review E . 100 (3): 033208. arXiv : 1905.12934 . Bibcode :2019PhRvE.100c3208K. doi :10.1103/PhysRevE.100.033208. PMID 31639953. S2CID 170078861.
^ Gupta, Ankur; Rawlings, James B. (abril de 2014). "Comparación de métodos de estimación de parámetros en modelos cinéticos químicos estocásticos: ejemplos en biología de sistemas". AIChE Journal . 60 (4): 1253–1268. Bibcode :2014AIChE..60.1253G. doi :10.1002/aic.14409. PMC 4946376 . PMID 27429455.
^ Véase Gill 2008.
^ Véase Robert y Casella 2004.
^ Banerjee, Sudipto; Carlin, Bradley P.; Gelfand, Alan P. (12 de septiembre de 2014). Modelado y análisis jerárquico de datos espaciales (segunda edición). CRC Press. pág. xix. ISBN978-1-4398-1917-3.
^ Geman, Stuart; Geman, Donald (noviembre de 1984). "Relajación estocástica, distribuciones de Gibbs y restauración bayesiana de imágenes". IEEE Transactions on Pattern Analysis and Machine Intelligence . PAMI-6 (6): 721–741. doi :10.1109/TPAMI.1984.4767596. ISSN 0162-8828. PMID 22499653. S2CID 5837272.
^ Gilks, WR; Wild, P. (1 de enero de 1992). "Muestreo de rechazo adaptativo para el muestreo de Gibbs". Revista de la Royal Statistical Society. Serie C (Estadística aplicada) . 41 (2): 337–348. doi :10.2307/2347565. JSTOR 2347565.
^ Gilks, WR; Best, NG ; Tan, KKC (1 de enero de 1995). "Muestreo de metrópolis con rechazo adaptativo dentro del muestreo de Gibbs". Revista de la Royal Statistical Society. Serie C (Estadística aplicada) . 44 (4): 455–472. doi :10.2307/2986138. JSTOR 2986138.
^ Lee, Se Yoon (2021). "Inferencia variacional mediante el muestreador de Gibbs y el ascenso de coordenadas: una revisión de la teoría de conjuntos". Comunicaciones en estadística: teoría y métodos . 51 (6): 1–21. arXiv : 2008.01006 . doi :10.1080/03610926.2021.1921214. S2CID 220935477.
^ Véase Stramer 1999.
^ Véase Green 1995.
^ Del Moral, Pierre (2013). Simulación de campo medio para integración de Monte Carlo. Chapman & Hall/CRC Press. pág. 626.
^ Del Moral, Pierre (2004). Fórmulas de Feynman-Kac. Aproximaciones genealógicas y de partículas interactuantes. Springer. pág. 575.
^ Del Moral, Pierre; Micló, Laurent (2000). "Aproximaciones de sistemas de partículas ramificadas e interactivas de fórmulas de Feynman-Kac con aplicaciones al filtrado no lineal". En Jacques Azéma; Michel Ledoux; Michel Émery; Marc Yor (eds.). Seminario de Probabilités XXXIV (PDF) . Apuntes de conferencias de matemáticas. vol. 1729, págs. 1-145. doi :10.1007/bfb0103798. ISBN978-3-540-67314-9.
^ Del Moral, Pierre (2006). "Muestreadores secuenciales de Monte Carlo". Revista de la Royal Statistical Society. Serie B (Metodología estadística) . 68 (3): 411–436. arXiv : cond-mat/0212648 . doi :10.1111/j.1467-9868.2006.00553.x. S2CID 12074789.
^ Sobol, Ilya M (1998). "Sobre integraciones cuasi-monte carlo". Matemáticas y computadoras en simulación . 47 (2): 103–112. doi :10.1016/s0378-4754(98)00096-2.
^ Chen, S.; Dick, Josef; Owen, Art B. (2011). "Consistencia de cadenas de Markov cuasi-Monte Carlo en espacios de estados continuos". Anales de Estadística . 39 (2): 673–701. arXiv : 1105.1896 . doi : 10.1214/10-AOS831 .
^ Tribble, Seth D. (2007). Algoritmos de Monte Carlo con cadenas de Markov que utilizan secuencias de conducción completamente uniformemente distribuidas (Tesis doctoral). Universidad de Stanford. ProQuest 304808879.
^ L'Ecuyer, P.; Lécot, C.; Tuffin, B. (2008). "Un método de simulación cuasi-Monte Carlo aleatorio para cadenas de Markov" (PDF) . Investigación de operaciones . 56 (4): 958–975. doi :10.1287/opre.1080.0556.
^ L'Ecuyer, P.; Munger, D.; Lécot, C.; Tuffin, B. (2018). "Métodos de ordenamiento y tasas de convergencia para Array-RQMC: algunas comparaciones empíricas". Matemáticas y computadoras en simulación . 143 : 191–201. doi :10.1016/j.matcom.2016.07.010.
^ ab Gelman, A.; Rubin, DB (1992). "Inferencia a partir de simulación iterativa utilizando múltiples secuencias (con discusión)" (PDF) . Ciencia estadística . 7 (4): 457–511. Bibcode :1992StaSc...7..457G. doi : 10.1214/ss/1177011136 .
^ Cowles, MK; Carlin, BP (1996). "Diagnóstico de convergencia de Monte Carlo de cadena de Markov: una revisión comparativa". Revista de la Asociación Estadounidense de Estadística . 91 (434): 883–904. CiteSeerX 10.1.1.53.3445 . doi :10.1080/01621459.1996.10476956.
^ Hill, SD; Spall, JC (2019). "Estacionariedad y convergencia del algoritmo Metropolis-Hastings: perspectivas sobre aspectos teóricos". Revista IEEE Control Systems . 39 (1): 56–67. doi :10.1109/MCS.2018.2876959. S2CID 58672766.
^ Foreman-Mackey, Daniel; Hogg, David W.; Lang, Dustin; Goodman, Jonathan (25 de noviembre de 2013). "emcee: The MCMC Hammer". Publicaciones de la Sociedad Astronómica del Pacífico . 125 (925): 306–312. arXiv : 1202.3665 . Código Bibliográfico :2013PASP..125..306F. doi :10.1086/670067. S2CID 88518555.
^ Phan, Du; Pradhan, Neeraj; Jankowiak, Martin (24 de diciembre de 2019). "Efectos componibles para programación probabilística flexible y acelerada en NumPyro". arXiv : 1912.11554 [stat.ML].
Fuentes
Christophe Andrieu, Nando De Freitas, Arnaud Doucet y Michael I. Jordan Introducción a MCMC para el aprendizaje automático, 2003
Asmussen, Søren; Glynn, Peter W. (2007). Simulación estocástica: algoritmos y análisis . Modelado estocástico y probabilidad aplicada. Vol. 57. Springer.
Atzberger, P. "Introducción a los métodos de Montecarlo" (PDF) .
Casella, George; George, Edward I. (1992). "Explicación del muestreador de Gibbs". The American Statistician . 46 (3): 167–174. CiteSeerX 10.1.1.554.3993 . doi :10.2307/2685208. JSTOR 2685208.
Chib, Siddhartha ; Greenberg, Edward (1995). "Entendiendo el algoritmo Metropolis-Hastings". The American Statistician . 49 (4): 327–335. doi :10.1080/00031305.1995.10476177. JSTOR 2684568.
Gill, Jeff (2008). Métodos bayesianos: un enfoque de las ciencias sociales y del comportamiento (2.ª ed.). Chapman y Hall /CRC. ISBN978-1-58488-562-7.
Green, PJ (1995). "Cálculo de Monte Carlo de cadena de Markov con salto reversible y determinación de modelo bayesiano". Biometrika . 82 (4): 711–732. CiteSeerX 10.1.1.407.8942 . doi :10.1093/biomet/82.4.711.
Neal, Radford M. (2003). "Muestreo por cortes". Anales de estadística . 31 (3): 705–767. doi : 10.1214/aos/1056562461 . JSTOR 3448413.
Neal, Radford M. (1993). "Inferencia probabilística utilizando métodos de Monte Carlo con cadenas de Markov".
Robert, Christian P.; Casella, G. (2004). Métodos estadísticos de Monte Carlo (2ª ed.). Saltador. ISBN978-0-387-21239-5.
Smith, RL (1984). "Procedimientos de Monte Carlo eficientes para generar puntos distribuidos uniformemente sobre regiones delimitadas". Investigación de operaciones . 32 (6): 1296–1308. doi :10.1287/opre.32.6.1296. hdl : 2027.42/7681 .
Spall, JC (abril de 2003). "Estimación mediante el método Monte Carlo de cadenas de Markov". Revista IEEE Control Systems . 23 (2): 34–45. doi :10.1109/mcs.2003.1188770.
Stramer, O.; Tweedie, R. (1999). "Modelos de tipo Langevin II: Candidatos autodirigidos para algoritmos MCMC". Metodología y computación en probabilidad aplicada . 1 (3): 307–328. doi :10.1023/A:1010090512027. S2CID 1512689.
Lectura adicional
Diaconis, Persi (abril de 2009). "La revolución de Monte Carlo de la cadena de Markov" (PDF) . Bull. Amer. Math. Soc. 46 (2): 179–205. doi : 10.1090/s0273-0979-08-01238-x . S 0273-0979(08)01238-X.
Richey, Matthew (mayo de 2010). "La evolución de los métodos de Monte Carlo de cadenas de Markov" (PDF) . The American Mathematical Monthly . 117 (5): 383–413. CiteSeerX 10.1.1.295.4478 . doi :10.4169/000298910x485923. S2CID 13630404.