Vorteile von Box-Muller gegenüber der inversen CDF-Methode zur Simulation der Normalverteilung?

15

Um eine Normalverteilung aus einer Reihe einheitlicher Variablen zu simulieren, gibt es verschiedene Techniken:

  1. Der Box-Muller-Algorithmus , bei dem zwei unabhängige Uniformvariablen abgetastet werden, variiert auf und transformiert sie in zwei unabhängige Standardnormalverteilungen über: Z 0 = (0,1)

    Z0=2lnU1cos(2πU0)Z1=2lnU1sin(2πU0)
  2. die CDF-Methode , bei der man das normale cdf einer Uniformvariablen gleichsetzen kann : und daraus F ( Z ) = U Z = F - 1 ( U )(F(Z))

    F(Z)=U
    Z=F1(U)

Meine Frage ist: Was ist rechnerisch effizienter? Ich würde denken, es ist die letztere Methode - aber die meisten Zeitungen, die ich lese, verwenden Box-Müller - warum?

Zusätzliche Information:

Die Umkehrung der normalen CDF ist bekannt und gegeben durch:

F1(Z)=2erf1(2Z1),Z(0,1).

Also:

Z=F1(U)=2erf1(2U1),U(0,1).
user2350366
quelle
1
Was ist die Umkehrung des normalen cdf? Sie kann nur dann analytisch berechnet werden, wenn die ursprüngliche CDF mit der stückweisen linearen Funktion approximiert wird.
Artem Sobolev
Sind die beiden nicht eng miteinander verwandt? Ich glaube, Box Muller ist ein besonderer Fall für die Erzeugung von 2 Variablen.
TTNPHNS
Hallo Barmaley, ich habe oben einige weitere Informationen hinzugefügt. Die inverse CDF hat einen Ausdruck - jedoch muss der rechnerisch berechnet werden - weshalb wird möglicherweise Box Müller bevorzugt? Ich nahm an, dass in Nachschlagetabellen berechnet wird, ähnlich wie die Werte von und ? Also nicht viel rechenintensiver? Ich kann mich jedoch irren. erf 1 sin cosineerf1erf1sincosine
user2350366
2
Es gibt Versionen von Box-Muller ohne Sünde und Kosinus.
Xi'an,
2
@Dilip Für Anwendungen mit sehr geringer Genauigkeit wie Computergrafiken können Sinus und Cosinus in der Tat durch Verwendung geeigneter Nachschlagetabellen optimiert werden. Für statistische Anwendungen wird eine solche Optimierung jedoch niemals verwendet. Letztendlich ist es nicht wirklich schwieriger, als oder zu berechnen , aber auf modernen Computersystemen sind elementare Funktionen, die sich auf beziehen - einschließlich der Triggerfunktionen - in der Regel vorhanden optimiert ( und waren grundlegende Anweisungen auf dem Intel 8087-Chip!), während erf entweder nicht verfügbar ist oder auf einer höheren (= langsameren) Ebene codiert wurde. log sqrt exp cos logerf1logsqrtexpcoslog
whuber

Antworten:

15

Aus rein probabilistischer Sicht sind beide Ansätze richtig und daher gleichwertig. Aus algorithmischer Sicht muss der Vergleich sowohl die Genauigkeit als auch die Rechenkosten berücksichtigen.

Box-Muller setzt auf einen einheitlichen Generator und kostet ungefähr das Gleiche wie dieser einheitliche Generator. Wie in meinem Kommentar erwähnt, können Sie ohne Sinus- oder Cosinus-Aufrufe davonkommen, wenn nicht ohne den Logarithmus:

  • generiere bis
    U1,U2iidU(1,1)
    S=U12+U221
  • nimm und definiereZ=2log(S)/S
    X1=ZU1, X2=ZU2

Der generische Inversionsalgorithmus erfordert beispielsweise den Aufruf der inversen normalen cdf qnorm(runif(N)) in R, die teurer als die oben genannten sein und, was noch wichtiger ist, in Bezug auf die Genauigkeit in den Schwänzen versagen kann, sofern die Quantilfunktion nicht gut codiert ist.

Um den Kommentaren von whuber zu folgen , ist der Vergleich von rnorm(N)und qnorm(runif(N))der Vorteil des inversen cdf, sowohl in Bezug auf die Ausführungszeit:

> 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` `

und in Bezug auf die Passform im Schwanz: Bildbeschreibung hier eingeben

Nach einem Kommentar von Radford Neal in meinem Blog möchte ich darauf hinweisen, dass der Standardwert rnormin R die Inversionsmethode verwendet, sodass der obige Vergleich sich auf die Schnittstelle und nicht auf die Simulationsmethode selbst auswirkt! So zitieren Sie die R-Dokumentation zu 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
quelle
3
logΦ1Φ1X1X2Ui1101
whuber
2
R 3.0.2rowSumsSqnorm(runif(N))InverseCDF[NormalDistribution[], #] &
1
Ich stimme zu, qnorm(runif(N))ist sogar 20% schneller alsrnorm(N)
Xi'an
3
Φ1sincos
1
Zum Vergleich mit einem i7-3740QM bei 2,7 GHz und R 3,12 für die folgenden Anrufe: 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)Ich bekomme mean 11.38363 median 11.18718für die Inversion und mean 13.00401 median 12.48802für Box-Muller
Avraham