30 votos

Simular una distribución uniforme en un disco

Intentaba simular la inyección de puntos aleatorios dentro de un círculo, de manera que cualquier parte del círculo tuviera la misma probabilidad de tener un defecto. Esperaba que el recuento por área de la distribución resultante siguiera una distribución de Poisson si dividía el círculo en rectángulos de igual área.

Como sólo requiere colocar puntos dentro de un área circular, he inyectado dos distribuciones aleatorias uniformes en coordenadas polares: $R$ (radio) y $\theta$ (ángulo polar).

Pero después de hacer esta inyección, obtengo claramente más puntos en el centro del círculo en comparación con el borde.

enter image description here

¿Cuál sería la forma correcta de realizar esta inyección a través del círculo de manera que los puntos se distribuyan aleatoriamente a lo largo del mismo?

43voto

jldugger Puntos 7490

Se desea que la proporción de puntos sea uniformemente proporcional a zona en lugar de la distancia al origen. Dado que el área es proporcional a la distancia al cuadrado, genere áreas aleatorias uniformes y tomar sus raíces cuadradas; escalar los resultados como se desee. Combínalo con un ángulo polar uniforme.

Esto es rápido y sencillo de codificar, eficiente en la ejecución (especialmente en una plataforma paralela), y genera exactamente el número de puntos prescrito.

Ejemplo

Esto funciona R para ilustrar el algoritmo.

n <- 1e4
rho <- sqrt(runif(n))
theta <- runif(n, 0, 2*pi)
x <- rho * cos(theta)
y <- rho * sin(theta)
plot(x, y, pch=19, cex=0.6, col="#00000020")

enter image description here

7voto

David Puntos 41

Muestreo de rechazo se puede utilizar. Esto significa que podemos muestrear a partir de una distribución uniforme 2D, y seleccionar muestras que satisfagan la condición de disco.

He aquí un ejemplo.

x=runif(1e4,-1,1)
y=runif(1e4,-1,1)

d=data.frame(x=x,y=y)
disc_sample=d[d$x^2+d$y^2<1,]
plot(disc_sample)

enter image description here

7voto

Aksakal Puntos 11351

Te daré una respuesta general de n dimensiones que funciona también para el caso de dos dimensiones, por supuesto. En tres dimensiones un análogo de un disco es el volumen de una bola sólida (esfera).

Hay dos enfoques que voy a discutir. Uno de ellos lo llamaría "preciso" y obtendrás una solución completa con ella en R. La segunda la llamo heurística Y es sólo la idea, no se proporciona ninguna solución completa.

"Solución "precisa

Mi solución se basa en Trabajos de Marsaglia y Muller . Básicamente, sucede que el vector gaussiano normalizado a su norma daría los puntos uniformemente distribuidos en una hiperesfera d-dimensional:

enter image description here

Es lo mismo que los puntos distribuidos uniformemente en un círculo en dos dimensiones. Para extenderlo a toda la superficie de un disco, hay que escalarlos aún más por el radio. El cuadrado de un radio procede de una distribución uniforme en dos dimensiones o elevado al poder $d$ en las dimensiones d. Entonces, se eleva a la potencia $1/d$ un número aleatorio uniforme para obtener el radio correctamente distribuido. Aquí tienes un código completo en R para dos dimensiones, que puedes extender fácilmente a cualquier número de dimensiones:

n <- 1e4
rho <- sqrt(runif(n))
# d - # of dimensions of hyperdisk
d = 2
r = matrix(rnorm(n*d),nrow=n,ncol=d)
x = r/rep(sqrt(rowSums(r^2))/rho,1)
plot(x[,1], x[,2], pch=19, cex=0.6, col="#00000020")

enter image description here

Aquí hay un fragmento de código para el caso 3d, es decir, una bola sólida:

library(scatterplot3d)
n <- 1e3
# d - # of dimensions of hyperdisk

d=3
rho <- (runif(n))^(1/d)
r = matrix(rnorm(n*d),nrow=n,ncol=d)
x = r/rep(sqrt(rowSums(r^2))/rho,1)

scatterplot3d(x[,1], x[,2], x[,3])

enter image description here

Enfoque heurístico

Este enfoque se basa en un hecho no tan obvio: la proporción del volumen de la hiperesfera unitaria sobre el volumen de un hipercubo unitario que la encierra se reduce a cero cuando el número de dimensiones aumenta hasta el infinito. Esto puede verse fácilmente en la expresión para un volumen de una hiperesfera : $$V_n(R) = \frac{\pi^\frac{n}{2}}{\Gamma\left(\frac{n}{2} + 1\right)}R^n$$ Aquí se puede ver cómo el coeficiente delante de $R^n$ disminuye rápidamente a cero. Esta es otra manifestación del fenómeno que está ligado a lo que se conoce como maldición de la dimensionalidad en el aprendizaje automático.

¿Por qué es esto relevante para nuestro problema? Supongamos que queremos generar d números aleatorios uniformes, que serían los puntos aleatorios dentro del hipercubo d-dimensional. A continuación, se aplica el muestreo de rechazo para elegir los puntos dentro de la hiperesfera (también conocida como n-ball): $\sum_{i=1}^d x_i^2<R^2$ . El problema es que para un número elevado de dimensiones d, ¡casi todos los puntos estarán fuera de la esfera! Acabarás desechando la gran mayoría de tus muestras.

La solución que propongo es utilizar el muestreo de rechazo con sobremuestreo de los puntos cercanos al centro. Resulta que si se observara una de las coordenadas cartesianas de la muestra aleatoria uniforme del interior de la bola, su distribución convergería a una gaussiana con varianza $\frac 1 {\sqrt{d+2}}$ . Así, en lugar de elegir puntos uniformemente del cubo, muestreamos las coordenadas cartesianas utilizando la gaussiana, y luego aplicamos el muestreo de rechazo sobre ellas. De esta manera no estaríamos desperdiciando tantas variantes aleatorias generadas. Esto sería una forma de técnica de muestreo de importancia.

2voto

Q_Li Puntos 58

He aquí una solución alternativa en R :

n <- 1e4
## r <- seq(0, 1, by=1/1000)
r <- runif(n)
rho <- sample(r, size=n, replace=T, prob=r)
theta <- runif(n, 0, 2*pi)
x <- rho * cos(theta)
y <- rho * sin(theta)
plot(x, y, pch=19, cex=0.6, col="#00000020")

enter image description here

i-Ciencias.com

I-Ciencias es una comunidad de estudiantes y amantes de la ciencia en la que puedes resolver tus problemas y dudas.
Puedes consultar las preguntas de otros usuarios, hacer tus propias preguntas o resolver las de los demás.

Powered by:

X