Aby zasymulować rozkład normalny z zestawu zmiennych jednorodnych, istnieje kilka technik:
Algorytm Boxa-Mullera , w którym jeden próbkuje dwa niezależne jednolite zmienia się na i przekształca je w dwa niezależne standardowe rozkłady normalne poprzez:
metoda CDF , w której można zrównać normalny cdf do wariantu Uniform:
i wyprowadzić
Moje pytanie brzmi: która jest bardziej wydajna obliczeniowo? Myślałem, że to ta druga metoda - ale większość artykułów, które czytam, używa Box-Mullera - dlaczego?
Dodatkowe informacje:
Odwrotność normalnego CDF jest znana i podana przez:
Stąd:
normal-distribution
simulation
uniform
użytkownik2350366
źródło
źródło
Odpowiedzi:
Z czysto probabilistycznego punktu widzenia oba podejścia są poprawne, a zatem równoważne. Z punktu widzenia algorytmu porównanie musi uwzględniać zarówno precyzję, jak i koszt obliczeniowy.
Box-Muller opiera się na jednolitym generatorze i kosztuje mniej więcej tyle samo, co ten jednolity generator. Jak wspomniano w moim komentarzu, możesz uciec bez połączeń sinusoidalnych lub cosinusowych, jeśli nie bez logarytmu:
Ogólny algorytm inwersji wymaga na przykład wywołania odwrotnego normalnego pliku cdf
qnorm(runif(N))
w R, co może być bardziej kosztowne niż powyższe i, co ważniejsze, może zawieść w ogonach pod względem precyzji, chyba że funkcja kwantylu jest dobrze zakodowana.Aby podążać za komentarzami Whubera , porównanie
rnorm(N)
iqnorm(runif(N))
przewaga odwrotnego cdf, zarówno w czasie wykonania:i pod względem dopasowania do ogona:
Po komentarzu Radforda Neala na moim blogu , chcę zauważyć, że domyślna wartość
rnorm
w R wykorzystuje metodę inwersji, stąd powyższe porównanie odzwierciedla interfejs, a nie samą metodę symulacji! Aby zacytować dokumentację R dotyczącą RNG:źródło
R 3.0.2
rowSums
qnorm(runif(N))
InverseCDF[NormalDistribution[], #] &
qnorm(runif(N))
jest nawet o 20% szybszy niżrnorm(N)
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)
Dostajęmean 11.38363 median 11.18718
inwersję imean 13.00401 median 12.48802
dla Box-Mullera