Técnica de estadística computacional
En el análisis numérico y la estadística computacional , el muestreo por rechazo es una técnica básica que se utiliza para generar observaciones a partir de una distribución . También se denomina comúnmente método de aceptación-rechazo o "algoritmo de aceptación-rechazo" y es un tipo de método de simulación exacta. El método funciona para cualquier distribución con una densidad .
El muestreo por rechazo se basa en la observación de que para muestrear una variable aleatoria en una dimensión, se puede realizar un muestreo aleatorio uniforme del gráfico cartesiano bidimensional y mantener las muestras en la región bajo el gráfico de su función de densidad. [1] [2] [3] Nótese que esta propiedad se puede extender a funciones de N dimensiones.
Descripción
Para visualizar la motivación detrás del muestreo de rechazo, imagine que grafica la función de densidad de probabilidad (PDF) de una variable aleatoria en un tablero rectangular grande y lanza dardos hacia él. Suponga que los dardos están distribuidos uniformemente alrededor del tablero. Ahora retire todos los dardos que estén fuera del área bajo la curva. Los dardos restantes se distribuirán uniformemente dentro del área bajo la curva y las posiciones de estos dardos se distribuirán de acuerdo con la densidad de la variable aleatoria. Esto se debe a que hay más espacio para que los dardos caigan donde la curva es más alta y, por lo tanto, la densidad de probabilidad es mayor.
La visualización que se acaba de describir es equivalente a una forma particular de muestreo de rechazo donde la "distribución de propuestas" es uniforme. Por lo tanto, su gráfico es un rectángulo. La forma general de muestreo de rechazo supone que el tablero no es necesariamente rectangular, sino que tiene una forma acorde a la densidad de alguna distribución de propuestas (no necesariamente normalizada a ) de la que sabemos cómo tomar muestras (por ejemplo, utilizando el muestreo por inversión ). Su forma debe ser al menos tan alta en cada punto como la distribución de la que queremos tomar muestras, de modo que la primera encierre completamente a la segunda. De lo contrario, habría partes del área curva de la que queremos tomar muestras que nunca podrían alcanzarse.
El muestreo de rechazo funciona de la siguiente manera:
- Muestrear un punto en el eje de la distribución propuesta.
- Dibuje una línea vertical en esta posición, hasta el valor y máximo de la función de densidad de probabilidad de la distribución propuesta.
- Muestrear uniformemente a lo largo de esta línea desde 0 hasta el máximo de la función de densidad de probabilidad. Si el valor muestreado es mayor que el valor de la distribución deseada en esta línea vertical, rechazar el valor ‑ y volver al paso 1; de lo contrario, el valor ‑ es una muestra de la distribución deseada.
Este algoritmo se puede utilizar para tomar muestras del área bajo cualquier curva, independientemente de si la función se integra a 1. De hecho, escalar una función mediante una constante no tiene efecto sobre las posiciones muestreadas . Por lo tanto, el algoritmo se puede utilizar para tomar muestras de una distribución cuya constante de normalización es desconocida, lo que es común en las estadísticas computacionales .
Teoría
El método de muestreo de rechazo genera valores de muestreo a partir de una distribución objetivo con una función de densidad de probabilidad arbitraria utilizando una distribución propuesta con densidad de probabilidad . La idea es que se puede generar un valor de muestra a partir de muestreando en lugar de y aceptando la muestra de con probabilidad , repitiendo las extracciones de hasta que se acepte un valor. aquí hay un límite finito y constante en la razón de verosimilitud , que satisface sobre el soporte de ; en otras palabras, M debe satisfacer para todos los valores de . Tenga en cuenta que esto requiere que el soporte de debe incluir el soporte de —en otras palabras, siempre que .
La validación de este método es el principio de la envolvente: al simular el par , se produce una simulación uniforme sobre el subgrafo de . Aceptando solo pares tales que entonces se producen pares distribuidos uniformemente sobre el subgrafo de y, por lo tanto, marginalmente, una simulación de
Esto significa que, con suficientes réplicas, el algoritmo genera una muestra de la distribución deseada . Existen varias extensiones de este algoritmo, como el algoritmo Metropolis .
Este método se relaciona con el campo general de las técnicas de Monte Carlo , incluidos los algoritmos de Monte Carlo de cadena de Markov que también utilizan una distribución proxy para lograr la simulación a partir de la distribución objetivo . Forma la base de algoritmos como el algoritmo Metropolis .
La probabilidad de aceptación incondicional es la proporción de muestras propuestas que se aceptan, que es donde , y el valor de cada tiempo se genera bajo la función de densidad de la distribución de la propuesta .
El número de muestras necesarias para obtener un valor aceptado sigue una distribución geométrica con probabilidad , que tiene media . Intuitivamente, es el número esperado de iteraciones que se necesitan, como medida de la complejidad computacional del algoritmo.
Reescribe la ecuación anterior,
Observa que , debido a la fórmula anterior, donde es una probabilidad que solo puede tomar valores en el intervalo . Cuando se elige más cerca de uno, la probabilidad de aceptación incondicional es mayor cuanto menos varía esa razón, ya que es el límite superior para la razón de verosimilitud . En la práctica, se prefiere un valor más cercano a 1, ya que implica menos muestras rechazadas, en promedio, y por lo tanto menos iteraciones del algoritmo. En este sentido, uno prefiere tener lo más pequeño posible (mientras sigue siendo satisfactorio , lo que sugiere que generalmente debería parecerse de alguna manera. Sin embargo, observa que no puede ser igual a 1: esto implicaría que , es decir, que las distribuciones objetivo y propuesta son en realidad la misma distribución.
El muestreo por rechazo se utiliza con mayor frecuencia en los casos en que la forma de muestreo dificulta el muestreo. Una única iteración del algoritmo de rechazo requiere tomar muestras de la distribución propuesta, extraer de una distribución uniforme y evaluar la expresión. Por lo tanto, el muestreo por rechazo es más eficiente que cualquier otro método siempre que M veces el costo de estas operaciones (que es el costo esperado de obtener una muestra con el muestreo por rechazo) sea menor que el costo de obtener una muestra utilizando el otro método.
Algoritmo
El algoritmo, que fue utilizado por John von Neumann [4] y se remonta a Buffon y su aguja , [5] obtiene una muestra de una distribución con densidad utilizando muestras de una distribución con densidad de la siguiente manera:
- Obtenga una muestra de la distribución y una muestra de (la distribución uniforme sobre el intervalo unitario).
- Comprueba si .
- Si esto es así, aceptar como muestra extraída de ;
- En caso contrario, rechace el valor y vuelva al paso de muestreo.
El algoritmo tomará un promedio de iteraciones para obtener una muestra. [6]
Ventajas sobre el muestreo mediante métodos ingenuos
El muestreo por rechazo puede ser mucho más eficiente en comparación con los métodos ingenuos en algunas situaciones. Por ejemplo, dado un problema como el muestreo condicional en un conjunto dado , es decir, , a veces se puede simular fácilmente, utilizando los métodos ingenuos (por ejemplo, mediante el muestreo por transformada inversa ):
- Pruebe de forma independiente y acepte aquellos resultados satisfactorios.
- Salida: (ver también truncamiento (estadísticas) )
El problema es que este muestreo puede ser difícil e ineficiente si . El número esperado de iteraciones sería , lo que podría ser cercano al infinito. Además, incluso cuando se aplica el método de muestreo de rechazo, siempre es difícil optimizar el límite para la razón de verosimilitud. En la mayoría de los casos, es grande y la tasa de rechazo es alta, por lo que el algoritmo puede ser muy ineficiente. La familia exponencial natural (si existe), también conocida como inclinación exponencial, proporciona una clase de distribuciones de propuestas que pueden reducir la complejidad computacional, el valor de y acelerar los cálculos (ver ejemplos: trabajar con familias exponenciales naturales).
Muestreo de rechazo mediante inclinación exponencial
Dada una variable aleatoria , es la distribución objetivo. Supongamos, para simplificar, que la función de densidad se puede escribir explícitamente como . Elija la propuesta como
donde y . Claramente, , es de una familia exponencial natural . Además, la razón de verosimilitud es
Nótese que esto implica que de hecho es una función de generación de cumulantes , es decir,
- .
Es fácil derivar la función de generación de cumulantes de la propuesta y, por lo tanto, los cumulantes de la propuesta.
Como ejemplo simple, supongamos que bajo , , con . El objetivo es tomar una muestra de , donde . El análisis se realiza de la siguiente manera:
- Elija la forma de distribución propuesta , con función generadora de cumulantes como
- ,
- lo que implica además que es una distribución normal .
- Decide cuál es el pozo elegido para la distribución de la propuesta. En esta configuración, la forma intuitiva de elegir es establecer
- ,
- Es decir la distribución propuesta es por lo tanto .
- Escriba explícitamente el objetivo, la propuesta y la relación de probabilidad.
- Derive el límite para la razón de verosimilitud , que es una función decreciente para , por lo tanto
- Criterio de muestreo de rechazo: para , si
se mantiene, se acepta el valor de ; si no, se continúa muestreando nuevo y nuevo hasta su aceptación.
Para el ejemplo anterior, como medida de la eficiencia, el número esperado de iteraciones del método de muestreo de rechazo basado en la familia exponencial natural es del orden de , es decir , mientras que bajo el método ingenuo, el número esperado de iteraciones es , que es mucho más ineficiente.
En general, la inclinación exponencial de una clase paramétrica de distribución propuesta resuelve los problemas de optimización de manera conveniente, con sus propiedades útiles que caracterizan directamente la distribución de la propuesta. Para este tipo de problema, para simular condicionalmente en , entre la clase de distribuciones simples, el truco es usar la familia exponencial natural, que ayuda a obtener cierto control sobre la complejidad y a acelerar considerablemente el cálculo. De hecho, existen razones matemáticas profundas para usar la familia exponencial natural.
Desventajas
El muestreo de rechazo requiere conocer la distribución objetivo (específicamente, la capacidad de evaluar la PDF objetivo en cualquier punto).
El muestreo por rechazo puede dar lugar a que se tomen muchas muestras no deseadas si la función que se está muestreando está muy concentrada en una región determinada, por ejemplo, una función que tiene un pico en alguna ubicación. Para muchas distribuciones, este problema se puede resolver utilizando una extensión adaptativa (véase muestreo por rechazo adaptativo), o con un cambio apropiado de variables con el método de la razón de uniformes . Además, a medida que las dimensiones del problema se hacen más grandes, la razón del volumen incrustado con respecto a las "esquinas" del volumen de incrustación tiende a cero, por lo que pueden producirse muchos rechazos antes de que se genere una muestra útil, lo que hace que el algoritmo sea ineficiente y poco práctico. Véase la maldición de la dimensionalidad . En dimensiones altas, es necesario utilizar un enfoque diferente, normalmente un método de Monte Carlo de cadena de Markov como el muestreo de Metropolis o el muestreo de Gibbs . (Sin embargo, el muestreo de Gibbs, que descompone un problema de muestreo multidimensional en una serie de muestras de baja dimensión, puede utilizar el muestreo por rechazo como uno de sus pasos).
Muestreo de rechazo adaptativo
Para muchas distribuciones, es difícil encontrar una distribución propuesta que incluya la distribución dada sin desperdiciar mucho espacio. Una extensión del muestreo de rechazo que se puede utilizar para superar esta dificultad y muestrear de manera eficiente de una amplia variedad de distribuciones (siempre que tengan funciones de densidad logarítmicamente cóncavas , lo que de hecho es el caso de la mayoría de las distribuciones comunes, incluso aquellas cuyas funciones de densidad no son cóncavas en sí mismas) se conoce como muestreo de rechazo adaptativo (ARS) .
Hay tres ideas básicas en esta técnica, introducida finalmente por Gilks en 1992: [7]
- Si resulta de ayuda, defina su distribución de envolvente en el espacio logarítmico (por ejemplo, probabilidad logarítmica o densidad logarítmica). Es decir, trabaje con en lugar de hacerlo directamente.
- A menudo, las distribuciones que tienen funciones de densidad algebraicamente desordenadas tienen funciones de densidad logarítmica razonablemente más simples (es decir, cuando son desordenadas, pueden ser más fáciles de trabajar o, al menos, más cercanas a la linealidad por partes).
- En lugar de una única función de densidad de envolvente uniforme, utilice una función de densidad lineal por partes como envolvente.
- Cada vez que tenga que rechazar una muestra, puede utilizar el valor de la que evaluó para mejorar la aproximación por partes . Esto reduce la posibilidad de que su próximo intento sea rechazado. Asintóticamente, la probabilidad de tener que rechazar su muestra debería converger a cero y, en la práctica, a menudo muy rápidamente.
- Como se propone, cada vez que elegimos un punto que es rechazado, estrechamos la envolvente con otro segmento de línea que es tangente a la curva en el punto con la misma coordenada x que el punto elegido.
- Un modelo lineal por partes de la distribución logarítmica propuesta da como resultado un conjunto de distribuciones exponenciales por partes (es decir, segmentos de una o más distribuciones exponenciales, unidos de extremo a extremo). Las distribuciones exponenciales se comportan bien y se comprenden bien. El logaritmo de una distribución exponencial es una línea recta y, por lo tanto, este método implica esencialmente encerrar el logaritmo de la densidad en una serie de segmentos de línea. Esta es la fuente de la restricción de cóncava logarítmica: si una distribución es cóncava logarítmica, entonces su logaritmo es cóncavo (con forma de U invertida), lo que significa que un segmento de línea tangente a la curva siempre pasará por encima de la curva.
- Si no se trabaja en el espacio logarítmico, también se puede muestrear una función de densidad lineal por partes mediante distribuciones triangulares [8]
- Podemos aprovechar aún más el requisito de concavidad (logaritmo) para evitar potencialmente el costo de evaluar cuándo se acepta su muestra.
- Así como podemos construir un límite superior lineal por partes (la función "envolvente") usando los valores que tuvimos que evaluar en la cadena actual de rechazos, también podemos construir un límite inferior lineal por partes (la función "compresión") usando estos valores también.
- Antes de evaluar (lo potencialmente costoso) para ver si su muestra será aceptada, es posible que ya sepamos si será aceptada comparándola con la función de compresión (o en este caso) (idealmente más barata) que tenga disponible.
- Este paso de compresión es opcional, incluso cuando lo sugiere Gilks. En el mejor de los casos, le ahorra una sola evaluación adicional de su densidad objetivo (complicada y/o costosa). Sin embargo, presumiblemente para funciones de densidad particularmente costosas (y suponiendo la rápida convergencia de la tasa de rechazo hacia cero) esto puede marcar una diferencia considerable en el tiempo de ejecución final.
El método consiste básicamente en determinar sucesivamente una envolvente de segmentos de línea recta que se aproxima cada vez mejor al logaritmo sin dejar de estar por encima de la curva, comenzando con un número fijo de segmentos (posiblemente sólo una línea tangente). El muestreo a partir de una variable aleatoria exponencial truncada es sencillo. Sólo hay que tomar el logaritmo de una variable aleatoria uniforme (con el intervalo apropiado y el truncamiento correspondiente).
Desafortunadamente, ARS solo se puede aplicar para muestrear densidades de objetivos logarítmicos-cóncavos. Por esta razón, se han propuesto varias extensiones de ARS en la literatura para abordar distribuciones de objetivos no logarítmicos-cóncavos. [9] [10] [11] Además, se han diseñado diferentes combinaciones de ARS y el método Metropolis-Hastings para obtener un muestreador universal que construya densidades de propuesta autoajustables (es decir, una propuesta construida y adaptada automáticamente al objetivo). Esta clase de métodos a menudo se denominan algoritmos de muestreo de Metropolis de rechazo adaptativo (ARMS) . [12] [13] Las técnicas adaptativas resultantes siempre se pueden aplicar, pero las muestras generadas están correlacionadas en este caso (aunque la correlación desaparece rápidamente a cero a medida que aumenta el número de iteraciones).
Véase también
Referencias
- ^ Casella, George; Robert, Christian P.; Wells, Martin T. (2004). Esquemas de muestreo generalizados de aceptación-rechazo . Instituto de Estadística Matemática. págs. 342–347. doi :10.1214/lnms/1196285403. ISBN. 9780940600614.
- ^ Neal, Radford M. (2003). "Muestreo por cortes". Anales de estadística . 31 (3): 705–767. doi : 10.1214/aos/1056562461 . MR 1994729. Zbl 1051.65007.
- ^ Bishop, Christopher (2006). "11.4: Muestreo por sectores". Reconocimiento de patrones y aprendizaje automático . Springer . ISBN 978-0-387-31073-2.
- ^ Forsythe, George E. (1972). "Método de comparación de von Neumann para muestreo aleatorio de distribuciones normales y otras". Matemáticas de la computación . 26 (120): 817–826. doi :10.2307/2005864. ISSN 0025-5718. JSTOR 2005864.
- ^ Legault, Geoffrey; Melbourne, Brett A. (1 de marzo de 2019). "Contabilización del cambio ambiental en modelos poblacionales estocásticos de tiempo continuo". Ecología teórica . 12 (1): 31–48. doi :10.1007/s12080-018-0386-z. ISSN 1874-1746.
- ^ Thomopoulos, Nick T. (19 de diciembre de 2012). Fundamentos de la simulación de Monte Carlo: métodos estadísticos para la construcción de modelos de simulación (edición de 2013). Nueva York, NY Heidelberg: Springer. ISBN 978-1-4614-6021-3.
- ^ Gilks, WR; Wild, P. (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.
- ^ Thomas, DB; Luk, W. (2007). "Generación de números aleatorios no uniformes mediante aproximaciones lineales por partes". IET Computers & Digital Techniques . 1 (4): 312–321. doi :10.1049/iet-cdt:20060188.
- ^ Hörmann, Wolfgang (1995-06-01). "Una técnica de rechazo para el muestreo de distribuciones T-cóncavas". ACM Trans. Math. Softw . 21 (2): 182–193. CiteSeerX 10.1.1.56.6055 . doi :10.1145/203082.203089. ISSN 0098-3500.
- ^ Evans, M.; Swartz, T. (1998-12-01). "Generación de variables aleatorias utilizando propiedades de concavidad de densidades transformadas". Journal of Computational and Graphical Statistics . 7 (4): 514–528. CiteSeerX 10.1.1.53.9001 . doi :10.2307/1390680. JSTOR 1390680.
- ^ Görür, Dilan; Teh, Yee Whye (1 de enero de 2011). "Muestreo de rechazo adaptativo cóncavo-convexo". Revista de estadística computacional y gráfica . 20 (3): 670–691. doi :10.1198/jcgs.2011.09058. ISSN 1061-8600.
- ^ 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.
- ^ Meyer, Renate; Cai, Bo; Perron, François (15 de marzo de 2008). "Muestreo de Metropolis de rechazo adaptativo utilizando polinomios de interpolación de Lagrange de grado 2". Computational Statistics & Data Analysis . 52 (7): 3408–3423. doi :10.1016/j.csda.2008.01.005.
Lectura adicional
- Robert, CP; Casella, G. (2004). Métodos estadísticos de Monte Carlo (segunda edición). Nueva York: Springer-Verlag.