przewodnikPoPakiecieR3

138 Wybrane procedury statystyczne

Za autora metody Monte Carlo uznaje się (>o)ukiego matematyka Stanisława Ularnu. który współpracując w Los Alamos z Enrico Fermim, Jolm von Ncumannem i Nicholasem Metropotisem, użył tej metody do modelowania reakcji łańcuchowej przebiegającej podczas wybuchu bomby wodorowej.

3.2 Liczby losowe

W wielu zagadnieniach, szczególnie tych związanych z modelowaniem metoda^S Monte Carlo bardzo złożonych zjawisk rzeczywistych (finansowych, biologicznych.^ fizycznych, chemicznych, meteorologicznych itp.), potrzebne są narzędzia pozywsjj łające emulować losowość, tak aby odzwierciedlić stochastyczną naturę badaneUj zjawiska. Dlatego też w wielu pakietach statystycznych oraz bibliotekach różnych fS* języków programowania dostępne są funkcje do generowania liczb psęudolo8owy<$jł|p czyli liczb o właściwościach podobnych do liczb losowych.

Takie funkcje dostępne są również w R. W tym podrozdziale opiszemy podstawo- * we zagadnienia związane z generatorami liczb losowych a w następnym podrozdziąf|||| przedstawimy funkcje do generowania liczb losowych z różnych popularnych rozkłą-^Sr & dów do oraz wyznaczania charakterystyk tych rozkładów.

3.2.1 Generatory liczb losowych

Generowane przez komputer liczby nazywane są pseudolosowymi, ponieważ mają emulować losowość, ale są wyznaczane w deterministyczny, choć często bardzo skont plikowany, sposób. Z uwagi na nieustannie rosnące zapotrzebowanie na dobre liczby £■ pseudoloso^e w literaturze opisanych jest wiele metod ich generowania a nad nowymi metodami wciąż pracuje wielu badaczy.

W R mamy możliwość wyboru jednego z już zaimplementowanych algorytmów: generujących liczby losowe lub napisania własnej funkcji do generowania tych liczb..;*4 . Aktualnie w pakiecie base zaimplementowanych jest sześć algorytmów generowania liczb losowych. Zmienić generator liczb losowych można wywołując fimkcję RNGkind(}' «(; oraz podając za jej argumenty nazwy dwóch wybranych generatorów, odpowiedni nio dla rozkładu jednostajnego i dla rozkładu normalnego. Wynikiem tej funkr^f cji są nazwy uprzednio używanych generatorów. Do generowania liczb z rozktójfł: du jednostajnego dostępne są następujące nazwy algorytmów: "Wichmann-Hili%||j "Marsaglia-Multicarry", "Super-Duper", "Mersenne-Twister", "Knuth-TAOCP^li)^{j "Knuth-TA0CP-2002", “user-supplied" (ostatnia pozycja umożliwia użytkowniko-V wi własnoręczne wskazanie funkcji generującej liczby losowe). Domyślnie stosowany^ jest najpopularniejszy "Mersenne-Twister".

Liczby losowe dla innych rozkładów można otrzymać z liczb losowych dla ro?(|ę| kładu jednostajnego wykorzystując np. metodę odwrotnej dystrybuanty, nie ma^ więc potrzeby określania generatora dla każdego rozkładu oddzielnie. Wyjątkiem....:, jest rozkład normalny, dla którego mamy możliwość wybrania metody generowa- ':-nia liczb losowych. Dostępne metody oznaczone są nazwami: “Kinderman-Ramage“,;jt| "Buggy Kinderaan-Ramage", "Ahrens-Dieter", "Box-Muller", "Inversion\^ "user-supplied". Domyślnie wykorzystywana jest metoda odwrotnej dystrybuantyjg; ("Inversion"). Ponieważ jednak dystrybuanta rozkładu normalnego nie ma prostej® postaci analitycznej, wykorzystywane są zaawansowane i wymagające oblic zeniowraw aproksymacje numeryczne (opis tych aproksymacji znaleźć można w opisie funkcji ąnormO). Innym popularnym algorytmem generowania liczb z rozkładu normulnegdjli znacznie szybszym od metody odwrotnej dystrybuanty, jest algorytm Boxa-Muller» . ("Box-Muller"). Więcej informacji o tym i innych generatorach przeczytać ino?-na w dokumentacji do funkcji RNGkind(base). Osoby zainteresowane tym tematem ; znajdą więcej informacji o generatorach w książce (31|. Poniżej przedstawiamy przy--:.| kład użycia funkcji RNGkindO do zmiany generatora.

.....SB

Liczby losowe

130

K-f.f uzbieramy generator “Super-Dup ^'■'.jednostajnego i metodę "Boz-Muller" do generowania liczb fy normalnego, wynikiem funkcji RNGkindO są poprzednio sto l£V generatory

dla generowania liczb z rozkładu z rozkładu tosowane

>(RHGkind("Super" , "Bor")) "Mersenne-Twister" "I™

Xnversion"

&-•# wybierać generatory można podając tylko pierwsze znaki nazw V’fRNGRind("Super" , "Ahr"))

j-j] “Super-Duper" "Box-Mulłer"

f'• Generator to funkcja deterministyczna. Do losowaniu kolejnych liczb wykorzystuje tzw. ziarno (ang. seed), całkowicie determinujące wartości kolejnych liczb pseu-

djjlosowych.

j y Dla ustalonego generatora i ziarna generowane będą identyczne liczby ‘ ’ losowe bez względu na system operacyjny, nazwę komputera, rasę użyt-

/ '• kownika czy temperaturę w pokoju.

3& Tą właściwość programowych generatorów liczb pseudolosowych (istnieją też sprzętowe generatory, ale to osobny temat) można łatwo wykorzystać! Sterując wyborem ziarna umożliwiamy otrzymywanie identycznych ciągów liczb losowych na różnych komputerach. W ten sposób możemy powtarzać wyniki symulacji, odtwarzać te wyniki na innych komputerach lub kontynuować obliczenia przerwane w wyniku wystąpienia jakiegoś błędu.

; Informacje o ziarnie i generatorze liczb losowych przechowywana jest w wektorze . Random. seed. Pierwszym elementem tego wektora jest informacja, który generator jest wykorzystywany, pozostałe elementy przechowują ziarno generatora. Dla generatora Mersenne-Twister ziarno jest wektorem 625 liczb całkowitych, dla pozostałych 'generatorów ziarno ma inną, zwykle krótszą, długość. Poniżej przykład zapisywania i odtwarzania ziarna generatora.

> .Random.seedfl :7]

59254738

[13 403 624 37667055 -1124359580 -1088152043

786696638

> tt zapisujemy ziarno generatora

> save.seed <- .Random.seed

V # losujemy 10 liczb losowych ‘> rnorm(lO)

[13 -0.393 1.175 0.394-0.732 0.930-2.766 0.289 0.341-0.559 1.133

> # odtwarzamy stan generatora liczb losowych

> Random.seed <- save.seed

> * Deja vu ?

> rnorm(lO)

CU -0.393 1.175 0.394 -0.732 0.930 -2.766 0.289 0.341 -0.559 1.133

Pamiętanie ziarna złożonego z-625 liczb nic jest specjalnie wygodne. Do łatwej inicjacji ziarna zalecane jest korzystanie z funkcji set.seed(base). Argumentem tej funkcji jest jedna liczba, która jest zamieniana na ziarno odpowiedniej długości. Poniżej przykład wywołania tej funkcji.