Comment générer des points uniformément répartis sur la surface de la sphère unité 3-d?

Je me demande comment générer des points uniformément répartis sur la surface de la sphère d'unité 3D? Aussi, après avoir généré ces points, quel est le meilleur moyen de visualiser et de vérifier s’ils sont vraiment uniformes sur la surface ?$x^2+y^2+z^2=1$

Une méthode standard consiste à générer trois normales normales et à en construire un vecteur unitaire. C'est-à-dire que lorsque et , alors est uniforme distribué sur la sphère. Cette méthode fonctionne bien pour sphères de dimension, aussi.$X_i \sim N(0,1)$$\lambda^2 = X_1^2 + X_2^2 + X_3^2$$(X_1/\lambda, X_2/\lambda, X_3/\lambda)$$d$

En 3D, vous pouvez utiliser l'échantillonnage de rejet: dessinez partir d'une distribution uniforme jusqu'à ce que la longueur de soit inférieure ou égale à 1, puis - comme avec la méthode précédente - normaliser le vecteur en unité de longueur. Le nombre d'essais prévus par point sphérique est égal à = 1,91. Dans les dimensions supérieures, le nombre d'essais escompté devient si important que cela devient rapidement irréalisable. [ - 1 , 1 ] ( X 1 , X 2 , X 3 ) 2 3 / ( 4 π / 3 $X_i$$[-1,1]$$(X_1, X_2, X_3)$$2^3/(4 \pi / 3)$

Il y a plusieurs façons de vérifier l'uniformité . La fonction K de Ripley est un moyen astucieux, même s’il est quelque peu intensif en calcul . Le nombre attendu de points à l'intérieur d'une distance (euclidienne 3D) de n'importe quel emplacement de la sphère est proportionnel à l'aire de la sphère à une distance égale à . En calculant toutes les distances entre points, vous pouvez comparer les données à cet idéal.ρ π ρ $\rho$$\rho$$\pi\rho^2$

Les principes généraux de construction des graphiques statistiques suggèrent un bon moyen de faire la comparaison: tracer les résidus stabilisés par la variance rapport à où est le plus petit des distances mutuelles et . L'intrigue devrait être proche de zéro. (Cette approche est non conventionnelle.)i = 1 , 2 , ... , n ( n - 1 ) / 2 = m d [ i ] i e e i = 2 $e_i(d_{[i]} - e_i)$$i = 1, 2, \ldots, n(n-1)/2=m$$d_{[i]}$$i^\text{th}$$e_i = 2\sqrt{i/m}$

Voici une image de 100 tirages indépendants d'une distribution sphérique uniforme obtenue avec la première méthode:

Voici l'intrigue de diagnostic des distances:

L'échelle y suggère que ces valeurs sont toutes proches de zéro.

Voici l'accumulation de 100 de ces graphiques pour suggérer quels écarts de taille pourraient en réalité être des indicateurs significatifs de non-uniformité:

(Ces parcelles ressemblent énormément à des ponts browniens ... il y a peut-être d'intéressantes découvertes théoriques qui se cachent ici.)

Enfin, voici l’intrigue de diagnostic pour un ensemble de 100 points aléatoires uniformes plus 41 autres points uniformément répartis dans l’hémisphère supérieur uniquement:

Par rapport à la distribution uniforme, il montre une diminution significative des distances moyennes entre les points jusqu'à un hémisphère. Cela en soi n'a pas de sens, mais l'information utile ici est que quelque chose n'est pas uniforme à l'échelle d'un hémisphère. En effet, ce graphique détecte facilement qu'un hémisphère a une densité différente de celle de l'autre. (Un test chi-carré plus simple ferait cela avec plus de puissance si vous saviez à l'avance quel hémisphère il faut tester sur une infinité d'infinis.)

Voici un code R assez 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)     

Il est très simple de voir d'après la construction que et donc mais si cela doit être testé, alors$x^2+y^2 = 1- z^2$$x^2+y^2+z^2=1$

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

et facile à tester que et sont uniformément répartis sur ( est évidemment) avec$x$$y$$[-1,1]$$z$

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

Clairement, étant donné une valeur de , et sont uniformément répartis autour d'un cercle de rayon et ceci peut être testé en regardant la distribution de l'arctangente de leur rapport. Mais puisque a la même distribution marginale que et que , une déclaration similaire est vraie pour toute paire, et cela aussi peut être testé. $z$$x$$y$$\sqrt{1-z^2}$$z$$x$$y$

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

Si elles ne sont toujours pas convaincues, les prochaines étapes consisteront à examiner une rotation arbitraire à trois dimensions ou le nombre de points correspondant à un angle solide donné, mais cela commence à devenir plus compliqué et, à mon avis, inutile.

Si vous voulez échantillonner des points uniformément répartis sur la sphère 3D (c'est-à-dire la surface d'une boule 3D), utilisez un simple rejet, ou utilisez la méthode de Marsaglia (Ann. Math. Statist., 43 (1972), p. 645– 646). Pour les petites dimensions, le taux de rejet est assez faible.

Si vous souhaitez générer des points aléatoires à partir de sphères et de billes de dimensions supérieures, cela dépend de l'objectif et de l'échelle de la simulation. Si vous ne souhaitez pas effectuer de grandes simulations, utilisez la méthode de Muller (Commun. ACM, 2 (1959), p. 19–20) ou sa version "à boule" (voir le document de Harman & Lacko cité ci-dessus). C'est:

pour obtenir un échantillon uniformément distribué sur une n-sphère (surface) 1) générer X à partir de la distribution normale standard n-dimensionnelle 2) diviser chaque composant de X par la norme euclidienne de X

pour obtenir un échantillon uniformément réparti sur une n-boule (intérieur) 1) générer X à partir de (n + 2) -dimensional distribution normale standard 2) diviser chaque composant de X par la norme euclidienne de X et ne prendre que les n premiers composants

Si vous souhaitez effectuer de grandes simulations, vous devez alors rechercher des méthodes plus spécialisées. Sur demande, je peux vous envoyer le document de Harman et Lacko sur les méthodes de distribution conditionnelle, qui fournit la classification et les généralisations de certains algorithmes mentionnés dans cette discussion. Le contact est disponible sur mon site Web (http://www.iam.fmph.uniba.sk/ospm/Lacko)

Si vous voulez vérifier si vos points sont vraiment uniformes à la surface ou à l’intérieur d’une balle, regardez les marginaux (tous devraient le même, en raison de l’invariance de rotation, la norme au carré d’un échantillon projeté est distribuée bêta).

J'ai eu un problème similaire (n-sphère) au cours de ma thèse et l'un des «experts» locaux a suggéré l'échantillonnage de rejet à partir d'un n-cube! Ceci, bien sûr, aurait pris l'âge de l'univers puisque je regardais n dans l'ordre des centaines.

L'algorithme que j'ai fini par utiliser est très simple et publié dans:

WP Petersen et A. Bernasconic Échantillonnage uniforme dans une n-sphère: méthode, méthode isotrope, Rapport technique, TR-97-06, Centre suisse de l'informatique scientifique

J'ai aussi cet article dans ma bibliographie que je n'ai pas encore regardé. Vous pouvez peut-être le trouver utile.

Harman, R. & Lacko, V. Algorithmes de décomposition pour l'échantillonnage uniforme de sphères et de billes Journal of Multivariate Analysis, 2010$n$$n$

J'ai déjà eu ce problème auparavant, et voici une alternative que j'ai trouvée,

En ce qui concerne la distribution elle-même, la formule que j’ai trouvée qui fonctionne décemment consiste à utiliser des coordonnées polaires (j’utilise en fait une variation des coordonnées de pales développées), puis à les convertir en coordonnées cartésiennes.

Le rayon est bien sûr le rayon de la sphère sur laquelle vous tracez. Ensuite, vous avez la deuxième valeur d'angle sur le plan plat, suivie de la troisième valeur qui est l'angle au-dessus ou au-dessous de ce plan.

Pour obtenir une distribution décente, supposons que U soit un nombre aléatoire uniformément réparti, r représente le rayon, a est la deuxième coordonnée polaire et b est la troisième coordonnée polaire,

a = U * 360 b = U + U-1 puis converti en cartésien via x = r * sin (b) sin (a) z = r sin (b) cos (a) y = r sin (b)

J'ai récemment trouvé ce qui suit qui est meilleur mathématiquement parlant, a = 2 (pi) * U b = cos ^ -1 (2U-1)

Pas très différent de ma formule originale en fait, bien que le mien soit en degrés vs radians.

Cette version récente est censée être utilisée pour les hypersphères, bien qu'aucune mention n'ait été faite sur la manière de l'obtenir.

Bien que je vérifie visuellement l'uniformité par la méthode plutôt peu coûteuse de créer des cartes pour Homeworld 2, puis de "jouer" ces cartes. En fait, comme les cartes sont créées avec des scripts Lua, vous pouvez intégrer votre formule directement dans la carte et ainsi vérifier plusieurs échantillonnages sans jamais quitter la partie. Peut-être pas scientifique, mais est une bonne méthode pour voir visuellement les résultats.

Voici le pseudocode:

1. $v \sim MultiVariateGaussian(\mu,\sigma I)$
2. $v = \frac{v}{ \| v \| }$

En pytorch:

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

Je ne comprends pas assez bien cela, mais Whuber m'a dit ceci:

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

est également correct, c’est-à-dire que l’échantillonnage d’une normale univariée pour chaque coordonnée.

Ma meilleure hypothèse serait de générer d'abord un ensemble de points uniformément répartis dans un espace à 2 dimensions, puis de projeter ces points sur la surface d'une sphère en utilisant une sorte de projection.

Vous devrez probablement mélanger et faire correspondre la façon dont vous générez les points avec la façon dont vous les mappez. En ce qui concerne la génération de génération de points 2D, je pense que les séquences brouillées à faible divergence seraient un bon point de départ (c’est-à-dire une séquence Sobol brouillée) car elles produisent généralement des points qui ne sont pas "groupés". Je ne suis pas aussi sûr du type de cartographie à utiliser, mais Woflram a fait apparaître la projection Gnonomic ... alors peut-être que cela pourrait fonctionner?

MATLAB a une implémentation décente de séquences à faible discordance que vous pouvez générer en utilisant q = sobolset(2)et brouiller en utilisant q = scramble(q). Il existe également une boîte à outils de mappage dans MATLAB avec un ensemble de fonctions de projection différentes que vous pouvez utiliser si vous ne souhaitez pas coder vous-même le mappage et les graphiques.