¿Cómo generar puntos distribuidos uniformemente en la superficie de la esfera de la unidad tridimensional?

68

Me pregunto cómo generar puntos uniformemente distribuidos en la superficie de la esfera de la unidad tridimensional. Además, después de generar esos puntos, ¿cuál es la mejor manera de visualizar y verificar si son realmente uniformes en la superficie ?x2+y2+z2=1

Qiang Li
fuente
Si por uniforme quieres decir "regular", no hay forma de hacerlo fuera de = 2, 4, 6, 8, 12, 20.n
Marcos
1
qué hay de malo con la muestra de un MultiVariateGaussian y ese vector simplemente lo normaliza: v = MultivariateNormal(torch.zeros(10000), torch.eye(10000))y luego v = v/v.norm(10000)
Pinocho el

Respuestas:

72

Un método estándar es generar tres normales estándar y construir un vector unitario a partir de ellas. Es decir, cuando y , entonces es uniforme distribuido en la esfera. Este método también funciona bien para esferas dimensionales.λ 2 = X 2 1 + X 2 2 + X 2 3 ( X 1 / λ , X 2 / λ , X 3 / λ ) dXiN(0,1)λ2=X12+X22+X32(X1/λ,X2/λ,X3/λ)d

En 3D puede usar el muestreo de rechazo: dibuje desde una distribución uniforme hasta que la longitud de sea ​​menor o igual a 1, luego, al igual que con el método anterior, normalizar el vector a la unidad de longitud. El número esperado de ensayos por punto esférico es igual a = 1.91. En dimensiones más altas, el número esperado de ensayos aumenta tanto que rápidamente se vuelve impracticable. [ - 1 , 1 ] ( X 1 , X 2 , X 3 ) 2 3 / ( 4 π / 3 )Xi[1,1](X1,X2,X3)23/(4π/3)

Hay muchas formas de verificar la uniformidad . Una forma ordenada, aunque algo computacionalmente intensiva, es con la función K de Ripley . El número esperado de puntos dentro de la distancia (3D Euclidiana) de cualquier ubicación en la esfera es proporcional al área de la esfera dentro de la distancia , que es igual a . Al calcular todas las distancias entre puntos, puede comparar los datos con este ideal.ρ π ρ 2ρρπρ2

Los principios generales de la construcción de gráficos estadísticos sugieren que una buena manera de hacer la comparación es trazar los residuos estabilizados por varianza contra donde es el más pequeño de las distancias mutuas y . La trama debe estar cerca de cero. (Este enfoque no es convencional).i = 1 , 2 , , n ( n - 1 ) / 2 = m d [ i ] i th e i = 2 ei(d[i]ei)i=1,2,,n(n1)/2=md[i]ithei=2i/m

Aquí hay una imagen de 100 dibujos independientes de una distribución esférica uniforme obtenida con el primer método:

100 puntos esféricos uniformes

Aquí está la gráfica de diagnóstico de las distancias:

Trama de diagnóstico

La escala y sugiere que estos valores están todos cerca de cero.

Aquí está la acumulación de 100 de estos gráficos para sugerir qué desviaciones de tamaño podrían ser indicadores significativos de falta de uniformidad:

Valores simulados

(Estas parcelas se parecen mucho a los puentes brownianos ... puede haber algunos descubrimientos teóricos interesantes que acechan aquí).

Finalmente, aquí está el diagrama de diagnóstico para un conjunto de 100 puntos aleatorios uniformes más otros 41 puntos distribuidos uniformemente solo en el hemisferio superior:

Valores no uniformes simulados

En relación con la distribución uniforme, muestra una disminución significativa en las distancias entre puntos promedio a un rango de un hemisferio. Eso en sí mismo no tiene sentido, pero la información útil aquí es que algo no es uniforme en la escala de un hemisferio. En efecto, este gráfico detecta fácilmente que un hemisferio tiene una densidad diferente que el otro. (Una prueba de chi-cuadrado más simple haría esto con más potencia si supiera de antemano qué hemisferio probar de los infinitos posibles).

whuber
fuente
@whuber: muy bien! muchas gracias por tu publicación " se distribuye uniformemente en la esfera". ¿Dónde puedo encontrar referencias sobre su prueba, o es simplemente demostrable? (X1/λ,X2/λ,X3/λ)
Qiang Li
23
@Qiang, aquí está la esencia de la prueba: Sea donde denota la matriz de identidad . Luego, para cualquier matriz ortogonal , . Por lo tanto, la distribución de es invariante bajo rotaciones. Deje y observe que para cualquier ortogonal . Dado que es invariante a las rotaciones, también lo es , y dado que casi seguro, entonces debe estar distribuido uniformemente en la esfera. XN(0,In)Inn×nQQXN(0,In)XY=X/X2YQ=QX/QX2=QX/X2QXYY2=1
cardenal
3
@ Mike No, porque una distribución uniforme de la latitud no proporciona una distribución uniforme en la esfera. (La mayor parte de la superficie de la esfera se encuentra en latitudes más bajas cerca del ecuador, lejos de los polos. En su lugar, necesita una distribución uniforme de .)ϕcos(ϕ)
whuber
1
@Ahsan Debido a que las matrices ortogonales forman un grupo transitivo de transformaciones de la esfera que preservan el área, la distribución es uniforme sobre el subconjunto de la esfera de la forma : pero esta es toda la esfera. X/||X||2
whuber
1
@Cesar "Distribución uniforme" (en la esfera).
whuber
19

Aquí hay un código R bastante simple

n     <- 100000                  # large enough for meaningful tests
z     <- 2*runif(n) - 1          # uniform on [-1, 1]
theta <- 2*pi*runif(n) - pi      # uniform on [-pi, pi]
x     <- sin(theta)*sqrt(1-z^2)  # based on angle
y     <- cos(theta)*sqrt(1-z^2)     

Es muy simple ver en la construcción que y por lo tanto pero si necesita ser probado entoncesx2+y2=1z2x2+y2+z2=1

mean(x^2+y^2+z^2)  # should be 1
var(x^2+y^2+z^2)   # should be 0

y fácil de probar que cada uno de e están distribuidos uniformemente en ( obviamente es) conxy[1,1]z

plot.ecdf(x)  # should be uniform on [-1, 1]
plot.ecdf(y)
plot.ecdf(z)

Claramente, dado un valor de , e se distribuyen uniformemente alrededor de un círculo de radio y esto se puede probar observando la distribución del arcotangente de su relación. Pero dado que tiene la misma distribución marginal que y como , una afirmación similar es cierta para cualquier par, y esto también se puede probar. zxy1z2zxy

plot.ecdf(atan2(x,y)) # should be uniform on [-pi, pi]
plot.ecdf(atan2(y,z))
plot.ecdf(atan2(z,x))

Si aún no está convencido, los siguientes pasos serían observar una rotación arbitraria en 3-D o cuántos puntos cayeron dentro de un ángulo sólido dado, pero eso comienza a complicarse más, y creo que es innecesario.

Enrique
fuente
Me pregunto si su método de generar puntos (x, y, z) es esencialmente el mismo que el método de Whuber.
Qiang Li
3
No, no lo es: whuber usa tres números aleatorios mientras yo uso dos. El mío es un caso especial de "generar un punto en con una densidad adecuada [proporcional a ] y luego bajar una dimensión". Aquí convenientemente ya que esto es formalmente una esfera de 2 . [1,1](1z2)n/21n=2
Henry
3
O, más generalmente, genere puntos uniformes en el mapa usando cualquier proyección de área igual (la suya es una de área igual cilíndrica), luego proyecte hacia atrás. (+1)
whuber
@whuber: De hecho. Fuera de tema, pero para cualquier persona interesada, tengo una selección interactiva de proyecciones de mapas del mundo aquí , algunas de las cuales son de igual área
Henry
2
Este es más o menos el enfoque estándar utilizado en gráficos por computadora, basado en el teorema de la caja de sombreros de Arquímedes: mathworld.wolfram.com/ArchimedesHat-BoxTheorem.html
Edward KMETT
10

Si desea muestrear puntos distribuidos uniformemente en la esfera 3D (es decir, la superficie de una bola 3D), utilice un rechazo simple o el método de Marsaglia (Ann. Math. Statist., 43 (1972), págs. 645– 646). Para dimensiones bajas, la relación de rechazo es bastante baja.

Si desea generar puntos aleatorios a partir de esferas y bolas de dimensiones superiores, entonces depende del propósito y la escala de la simulación. Si no desea realizar simulaciones grandes, utilice el método de Muller (Commun. ACM, 2 (1959), págs. 19–20) o su versión "ball" (consulte el documento de Harman & Lacko citado anteriormente). Es decir:

para obtener una muestra distribuida uniformemente en una esfera n (superficie) 1) generar X a partir de la distribución normal estándar n-dimensional 2) dividir cada componente de X por la norma euclidiana de X

para obtener una muestra distribuida uniformemente en una bola n (interior) 1) generar X a partir de la distribución normal estándar (n + 2) -dimensional 2) dividir cada componente de X por la norma euclidiana de X y tomar solo los primeros n componentes

Si desea realizar simulaciones grandes, debe investigar métodos más especializados. Si lo solicita, puedo enviarle el documento de Harman y Lacko sobre métodos de distribución condicional, que proporciona la clasificación y generalizaciones de algunos algoritmos mencionados en esta discusión. El contacto está disponible en mi sitio web (http://www.iam.fmph.uniba.sk/ospm/Lacko)

Si desea verificar, si sus puntos son verdaderamente uniformes en la superficie o en el interior de una pelota, observe los márgenes (todos deberían ser iguales, debido a la invariancia rotacional, la norma al cuadrado de una muestra proyectada es beta distribuida).

V Lacko
fuente
qué hay de malo con la muestra de un MultiVariateGaussian y ese vector simplemente lo normaliza: v = MultivariateNormal(torch.zeros(10000), torch.eye(10000))y luego v = v/v.norm(10000)
Pinocho el
8

¡Tuve un problema similar (n-esfera) durante mi doctorado y uno de los 'expertos' locales sugirió el muestreo de rechazo de un n-cubo! Esto, por supuesto, habría tomado la edad del universo ya que estaba viendo n en el orden de los cientos.

El algoritmo que terminé usando es muy simple y publicado en:

WP Petersen y A. Bernasconic Muestreo uniforme de una esfera n: Informe técnico del método isotrópico, TR-97-06, Centro suizo de informática científica

También tengo este artículo en mi bibliografía que no he visto. Lo podrías encontrar útil.

Harman, R. y Lacko, V. Sobre algoritmos de descomposición para muestreo uniforme de esferas y bolas Journal of Multivariate Analysis, 2010nn

emakalic
fuente
¿Es posible publicar los enlaces donde puedo encontrar el texto completo de estas referencias? Gracias.
Qiang Li
No tengo el documento sobre mí, pero esta página parece describir el algoritmo (y varios otros) mlahanas.de/Math/nsphere.htm
emakalic
3
Según tengo entendido, (del documento de Petersen y Bernasconic) para una bola d-dimensional, uno puede generar el radio elevando una variante U (0,1) a la potencia (1 / d), y el último ángulo como un Variante U (0,2 ). Los ángulos intermedios se pueden obtener como , donde es . Para mí esto suena bastante simple. Lo que me pregunto es esto: si uso una secuencia cuasialeatoria para mis uniformes, ¿obtendré también la amabilidad en la pelota? πC.asin(uk)C1
πΓ(k2+0.5)Γ(k2+1)
Mohit
3

He tenido este problema antes, y aquí hay una alternativa que encontré,

En cuanto a la distribución en sí, la fórmula que encontré que funciona decentemente es usar coordenadas polares (en realidad uso una variación de coordenadas poler que se desarrolló), luego convertir a coordenadas cartesianas.

El radio es, por supuesto, el radio de la esfera en la que está trazando. Luego tiene el segundo valor para el ángulo en el plano, seguido del tercer valor, que es el ángulo por encima o por debajo de ese plano.

Para obtener una distribución decente, suponga que U es un número aleatorio distribuido uniformemente, r es radio, a es la segunda coordenada polar yb es la tercera coordenada polar,

a = U * 360 b = U + U-1 y luego convertir a cartesiano a través de x = r * sin (b) sin (a) z = r sin (b) cos (a) y = r sin (b)

Recientemente encontré lo siguiente que es mejor matemáticamente hablando, a = 2 (pi) * U b = cos ^ -1 (2U-1)

En realidad, no es muy diferente de mi fórmula original, aunque la mía es en grados frente a radianes.

Esta versión reciente supuestamente puede usarse para hiperesferas, aunque no se mencionó cómo lograrla.

Aunque compruebo visualmente la uniformidad mediante el método bastante barato de hacer mapas para Homeworld 2 y luego "jugar" esos mapas. De hecho, debido a que los mapas están hechos con scripts lua, puedes construir tu fórmula directamente en el mapa y así verificar múltiples muestreos sin salir del juego. Quizás no sea científico, pero es un buen método para ver visualmente los resultados.

TheDerpyAlicorn
fuente
2

Aquí está el pseudocódigo:

  1. vMultiVariateGaussian(μ,σI)
  2. v=vv

En pytorch:

v = MultivariateNormal(torch.zeros(10000), torch.eye(10000))
v = v/v.norm(2)

No entiendo esto lo suficientemente bien, pero Whuber me dijo que:

v = torch.normal(torch.zeros(10000), torch.eye(10000))
v = v/v.norm(2)

también es correcto, es decir, el muestreo de una normal univariada para cada coordenada.

Pinocho
fuente
0

Mi mejor conjetura sería generar primero un conjunto de puntos distribuidos uniformemente en un espacio bidimensional y luego proyectar esos puntos en la superficie de una esfera utilizando algún tipo de proyección.

Probablemente tendrá que mezclar y combinar la forma en que genera los puntos con la forma en que los asigna. En términos de la generación de generación de puntos 2D, creo que las secuencias codificadas de baja discrepancia serían un buen lugar para comenzar (es decir, una secuencia codificada de Sobol) ya que generalmente produce puntos que no están "agrupados". No estoy tan seguro sobre qué tipo de mapeo usar, pero Woflram apareció la proyección gnonómica ... ¿entonces tal vez eso podría funcionar?

MATLAB tiene una implementación decente de secuencias de baja discrepancia que puede generar usando q = sobolset(2)y codificar usando q = scramble(q). También hay una caja de herramientas de mapeo en MATLAB con un montón de diferentes funciones de proyección que puede usar en caso de que no quiera codificar el mapeo y los gráficos usted mismo.

Berk U.
fuente
1
¿Puede alguna de estas proyecciones preservar la uniformidad de la aleatoriedad? Nuevamente, ¿cómo puedo verificar si la distribución final de estos puntos está realmente uniformemente distribuida en la superficie de la esfera? Gracias.
Qiang Li
Lo siento, solo estaba hablando hipotéticamente ... Creo que las funciones de mapeo en MATLAB le permitirán verificar eso ya que tienen algunas visualizaciones incrustadas en ellas. Si no, también encontré un buen sitio web que habla sobre cómo generar puntos distribuidos uniformemente en una esfera en 3D usando cosas como ángulos aleatorios, etc. También tienen un código C allí. Echa un vistazo
Berk U.
3
Los puntos aleatorios uniformes en una proyección gnomónica no serán uniformes en la esfera, porque el gnomónico no es un área igual. La proyección propuesta por Henry, -> (de longitud-latitud a un rectángulo en ), es igual a área. ( λ , sin ( ϕ ) ) R 2(λ,ϕ)(λ,sin(ϕ))R2
whuber