El método polar de Marsaglia [1] es un método de muestreo de números pseudoaleatorios para generar un par de variables aleatorias normales estándar independientes . [2] Si bien es superior a la transformada Box-Muller , [3] el algoritmo Ziggurat es aún más eficiente. [4]
Las variables aleatorias normales estándar se utilizan con frecuencia en informática , estadística computacional y, en particular, en aplicaciones del método de Monte Carlo .
El método polar funciona eligiendo puntos aleatorios ( x , y ) en el cuadrado −1 < x <1, −1 < y <1 hasta
y luego devolviendo el par requerido de variables aleatorias normales como
o equivalente,
dónde y representan el coseno y el seno del ángulo que forma el vector ( x , y ) con el eje x .
Bases teóricas
La teoría subyacente se puede resumir de la siguiente manera:
Si u se distribuye uniformemente en el intervalo 0 ≤ u <1, entonces el punto (cos (2π u ), sin (2π u )) se distribuye uniformemente en la circunferencia unitaria x 2 + y 2 = 1, y multiplicar ese punto por una variable aleatoria independiente ρ cuya distribución es
producirá un punto
cuyas coordenadas se distribuyen conjuntamente como dos variables aleatorias normales estándar independientes.
Historia
Esta idea se remonta a Laplace , a quien Gauss atribuye haber encontrado lo anterior.
tomando la raíz cuadrada de
La transformación a coordenadas polares hace evidente que θ se distribuye uniformemente (densidad constante) de 0 a 2π, y que la distancia radial r tiene densidad
( r 2 tiene la distribución de chi cuadrado apropiada ).
Este método de producir un par de variables normales estándar independientes proyectando radialmente un punto aleatorio en la circunferencia de la unidad a una distancia dada por la raíz cuadrada de una variante de chi-cuadrado-2 se llama método polar para generar un par de variables aleatorias normales ,
Consideraciones prácticas
Una aplicación directa de esta idea,
se llama la transformada de Box-Muller , en la que la variante chi generalmente se genera como
pero esa transformación requiere funciones de logaritmo, raíz cuadrada, seno y coseno. En algunos procesadores, el coseno y el seno del mismo argumento se pueden calcular en paralelo usando una sola instrucción. [5] En particular, para las máquinas basadas en Intel, se puede utilizar la instrucción del ensamblador fsincos o la instrucción expi (disponible, por ejemplo, en D), para calcular
y simplemente separe las partes real e imaginaria.
Nota: Para calcular explícitamente la forma polar compleja, utilice las siguientes sustituciones en la forma general,
Dejar y Luego
Por el contrario, el método polar aquí elimina la necesidad de calcular un coseno y un seno. En cambio, mediante la resolución de un punto en el círculo unidad, estas dos funciones pueden ser reemplazados con las x y Y coordenadas normalizadas a laradio. En particular, un punto aleatorio ( x , y ) dentro del círculo unitario se proyecta sobre la circunferencia unitaria estableciendo y formando el punto
que es un procedimiento más rápido que calcular el coseno y el seno. Algunos investigadores sostienen que la instrucción if condicional (para rechazar un punto fuera del círculo unitario) puede hacer que los programas sean más lentos en los procesadores modernos equipados con predicción de ramificaciones y canalizaciones. [6] Además, este procedimiento requiere aproximadamente un 27% más de evaluaciones del generador de números aleatorios subyacente (solo de los puntos generados se encuentran dentro del círculo unitario).
A continuación, ese punto aleatorio de la circunferencia se proyecta radialmente a la distancia aleatoria requerida por medio de
utilizando la misma s porque esa s es independiente del punto aleatorio de la circunferencia y está distribuida uniformemente de 0 a 1.
Implementación
Implementación simple en Java usando la desviación estándar y media:
doble repuesto estático privado ; private static boolean hasSpare = false ; público estático sincronizado doble generateGaussian ( doble media , doble stdDev ) { if ( hasSpare ) { hasSpare = false ; devuelve repuesto * stdDev + mean ; } else { doble u , v , s ; haz { u = Math . aleatorio () * 2 - 1 ; v = Matemáticas . aleatorio () * 2 - 1 ; s = u * u + v * v ; } mientras ( s > = 1 || s == 0 ); s = Matemáticas . sqrt ( - 2.0 * Math . log ( s ) / s ); repuesto = v * s ; hasSpare = true ; return mean + stdDev * u * s ; } }
Una implementación no segura para subprocesos en C ++ usando la desviación estándar y media:
double generateGaussian ( doble media , doble stdDev ) { reserva doble estática ; static bool hasSpare = false ; if ( hasSpare ) { hasSpare = false ; devuelve repuesto * stdDev + mean ; } else { doble u , v , s ; hacer { u = ( rand () / (( doble ) RAND_MAX )) * 2.0 - 1.0 ; v = ( rand () / (( doble ) RAND_MAX )) * 2,0 - 1,0 ; s = u * u + v * v ; } while ( s > = 1.0 || s == 0.0 ); s = sqrt ( -2,0 * log ( s ) / s ); repuesto = v * s ; hasSpare = true ; return mean + stdDev * u * s ; } }
La implementación de C ++ 11 GNU GCC libstdc ++ de std :: normal_distribution usa el método polar de Marsaglia, como se cita en este documento .
Referencias
- ↑ Marsaglia, G .; Bray, TA (1964). "Un método conveniente para generar variables normales". Revisión SIAM . 6 (3): 260–264. doi : 10.1137 / 1006063 . JSTOR 2027592 .
- ^ Peter E. Kloeden Eckhard Platen Henri Schurz, Solución numérica de SDE a través de experimentos informáticos, Springer, 1994.
- ^ Glasserman, Paul (2004), Métodos Monte Carlo en Ingeniería Financiera , Aplicaciones de las Matemáticas: Modelado Estocástico y Probabilidad Aplicada, 53 , Springer, p. 66, ISBN 9780387004518.
- ^ Thomas, David B .; Luk, Wayne; Leong, Philip HW; Villasenor, John D. (2007). "Generadores de números aleatorios gaussianos". Encuestas de computación ACM . 39 (4): 11 – es. CiteSeerX 10.1.1.127.5773 . doi : 10.1145 / 1287620.1287622 .
- ^ Kanter, David. "Arquitectura de gráficos Ivy Bridge de Intel" . Tecnología del mundo real . Consultado el 8 de abril de 2013 .
- ^ Este efecto se puede intensificar en una GPU que genera muchas variantes en paralelo, donde un rechazo en un procesador puede ralentizar muchos otros procesadores. Consulte la sección 7 de Thomas, David B .; Howes, Lee W .; Luk, Wayne (2009), "Una comparación de CPU, GPU, FPGA y matrices de procesadores masivamente paralelos para la generación de números aleatorios", en Chow, Paul; Cheung, Peter YK (eds.), Proceedings of the ACM / SIGDA 17th International Symposium on Field Programmable Gate Arrays, FPGA 2009, Monterey, California, EE. UU., 22-24 de febrero de 2009 , Association for Computing Machinery , págs. 63-72 , CiteSeerX 10.1.1.149.6066 , doi : 10.1145 / 1508128.1508139 , ISBN 9781605584102.