Avantages de Box-Muller par rapport à la méthode CDF inverse pour simuler une distribution normale?

Afin de simuler une distribution normale à partir d'un ensemble de variables uniformes, il existe plusieurs techniques:

L'algorithme de Box-Muller , dans lequel on échantillonne deux variables uniformes indépendantes sur et les transforme en deux distributions normales standard indépendantes via: $(0,1)$
$Z_{0} = \sqrt{- 2 ln U_{1}} cos (2 π U_{0}) Z_{1} = \sqrt{- 2 ln U_{1}} sin (2 π U_{0})$ $Z_0 = \sqrt{-2\text{ln}U_1}\text{cos}(2\pi U_0)\\ Z_1 = \sqrt{-2\text{ln}U_1}\text{sin}(2\pi U_0)$
la méthode CDF , où l'on peut assimiler le cdf normal à une variable uniforme: et dériver $(F(Z))$
$F (Z) = U$ $F(Z) = U$ $Z = F^{- 1} (U)$ $Z = F^{-1}(U)$

Ma question est la suivante: lequel est le plus efficace sur le plan informatique? Je pense que c'est la dernière méthode - mais la plupart des articles que je lis utilisent Box-Muller - pourquoi?

Information additionnelle:

L'inverse du CDF normal est connu et donné par:

F^{- 1} (Z) = \sqrt{2} \erf^{- 1} (2 Z - 1), Z \in (0, 1) .

$F^{-1}(Z)\; =\; \sqrt2\;\operatorname{erf}^{-1}(2Z - 1), \quad Z\in(0,1).$

D'où:

Z = F^{- 1} (U) = \sqrt{2} \erf^{- 1} (2 U - 1), U \in (0, 1) .

$Z = F^{-1}(U)\; =\; \sqrt2\;\operatorname{erf}^{-1}(2U - 1), \quad U\in(0,1).$

normal-distribution simulation uniform user2350366
la source

Quel est l'inverse du cdf normal? Il ne peut pas être calculé analytiquement, uniquement si le CDF d'origine est approximé avec une fonction linéaire par morceaux.

Artem Sobolev

Les deux ne sont-ils pas étroitement liés? Box Muller, je pense, est un cas particulier pour la génération à 2 variables.

ttnphns

Salut Barmaley, j'ai ajouté plus d'informations ci-dessus. Le CDF inverse a une expression - cependant le doit être calculé par calcul - ce qui pourrait être la raison pour laquelle la boîte Muller est préférée? J'ai supposé que serait calculé dans les tables de recherche, tout comme les valeurs de et

sont? Donc, pas beaucoup plus cher en calcul? Je peux me tromper cependant.

{erf}^{- 1}

$\text{erf}^{-1}$

{erf}^{1}

$\text{erf}^1$

sin

$\text{sin}$

cosine

$\text{cosine}$

user2350366

Il existe des versions de Box-Muller sans péché ni cosinus.

Xi'an

@Dilip Pour les applications de très faible précision, telles que l'infographie, le sinus et le cosinus peuvent en effet être optimisés en utilisant des tables de recherche appropriées. Pour les applications statistiques, cependant, une telle optimisation n'est jamais utilisée. En fin de compte, il n'est pas vraiment plus difficile de calculer

que

, mais sur les systèmes informatiques modernes, les fonctions élémentaires liées à

compris les fonctions trig - ont tendance à être optimisées (

étaient des instructions de base sur Intel) 8087 chip!), Alors que erf n'est pas disponible ou a été codé à un niveau supérieur (= plus lent).

{erf}^{- 1}

$\text{erf}^{-1}$

\log

$\log$

sqrt

$\text{sqrt}$

\exp

$\exp$

\cos

$\cos$

\log

$\log$

whuber

Réponses:

D'un point de vue purement probabiliste, les deux approches sont correctes et donc équivalentes. D'un point de vue algorithmique, la comparaison doit considérer à la fois la précision et le coût de calcul.

Box-Muller s'appuie sur un générateur uniforme et coûte environ le même prix que ce générateur uniforme. Comme mentionné dans mon commentaire, vous pouvez vous en sortir sans appels sinus ou cosinus, sinon sans le logarithme:

générer jusqu'à $U_{1}, U_{2} \overset{iid}{\sim} U (- 1, 1)$ $U_1,U_2\stackrel{\text{iid}}{\sim}\mathcal{U}(-1,1)$ $S=U_1^2+U_2^2\le 1$
prendre et définissez $Z=\sqrt{-2\log(S)/S}$ $X_{1} = Z U_{1}, X_{2} = Z U_{2}$ $X_1=ZU_1\,,\ X_2=Z U_2$

L'algorithme d'inversion générique nécessite l'appel au cdf normal inverse, par exemple qnorm(runif(N))dans R, qui peut être plus coûteux que ce qui précède et plus important peut échouer dans les queues en termes de précision, à moins que la fonction quantile ne soit bien codée.

Pour suivre les commentaires de whuber , la comparaison rnorm(N)et qnorm(runif(N))l'avantage de l'inverse du cdf, tant en temps d'exécution:

> system.time(qnorm(runif(10^8)))
sutilisateur     système      écoulé
 10.137           0.120      10.251 
> system.time(rnorm(10^8))
utilisateur     système      écoulé
 13.417           0.060      13.472` `

et en termes d'ajustement dans la queue: entrez la description de l'image ici

Suite à un commentaire de Radford Neal sur mon blog , je tiens à souligner que la valeur rnormpar défaut dans R utilise la méthode d'inversion, d'où que la comparaison ci-dessus se reflète sur l'interface et non sur la méthode de simulation elle-même! Pour citer la documentation R sur RNG:

‘normal.kind’ can be ‘"Kinderman-Ramage"’, ‘"Buggy
 Kinderman-Ramage"’ (not for ‘set.seed’), ‘"Ahrens-Dieter"’,
 ‘"Box-Muller"’, ‘"Inversion"’ (the default), or ‘"user-supplied"’.
 (For inversion, see the reference in ‘qnorm’.)  The
 Kinderman-Ramage generator used in versions prior to 1.7.1 (now
 called ‘"Buggy"’) had several approximation errors and should only
 be used for reproduction of old results.  The ‘"Box-Muller"’
 generator is stateful as pairs of normals are generated and
 returned sequentially.  The state is reset whenever it is selected
 (even if it is the current normal generator) and when ‘kind’ is
 changed.

Xi'an
la source

\log

$\log$

\sqrt{}

$\sqrt{}$

Φ^{- 1}

$\Phi^{-1}$

Φ^{- 1}

$\Phi^{-1}$

X_{1}

$X_1$

X_{2}

$X_2$

U_{i}

$U_i$

- 1

$-1$

1

$1$

0

$0$

1

$1$

R 3.0.2rowSums

S

$S$ qnorm(runif(N))InverseCDF[NormalDistribution[], #] &

Je suis d'accord, qnorm(runif(N))est même 20% plus rapide quernorm(N)

Xi'an

Φ^{- 1}

$\Phi^{-1}$

\sin

$\sin$

\cos

$\cos$

À titre de comparaison, en utilisant un i7-3740QM @ 2.7Ghz et R 3.12, pour les appels suivants:

RNGkind(kind = NULL, normal.kind = 'Inversion');At <- microbenchmark(A <- rnorm(1e5, 0, 1), times = 100L);RNGkind(kind = NULL, normal.kind = 'Box-Muller');Bt <- microbenchmark(B <- rnorm(1e5, 0, 1), times = 100L)

Je reçois mean 11.38363 median 11.18718pour l'inversion et mean 13.00401 median 12.48802pour Box-Muller

Avraham