To zamieszanie między procedurami ładowania początkowego a procedurami Monte Carlo powtarza się, więc być może jest to tak dobre miejsce, jak każde inne, aby je rozwiązać. (Przykłady R
kodu mogą również pomóc w odrabianiu lekcji.)
Rozważ tę implementację bootstrap w R
:
boot <- function(x, t) { # Exact bootstrap of procedure t on data x
n <- length(x) # Must lie between 2 and 7 inclusive.
if (n > 7) {
stop("Sample size exceeds 7; use an approximate method instead.")
}
p <- c(n, 1:(n-1))
a <- rep(x, n^(n-1))
dim(a) <- rep(n, n)
y <- as.vector(a)
while (n > 1) {
n <- n-1
a <- aperm(a, p)
y <- cbind(as.vector(a), y)
}
apply(y, 1, t)
}
Szybkie spojrzenie potwierdzi, że jest to obliczenie deterministyczne : nie są generowane ani używane żadne losowe wartości. (Zostawię szczegóły jego wewnętrznego funkcjonowania zainteresowanym czytelnikom, aby sami się zorientowali).
Argumenty boot
to partia danych liczbowych w tablicy x
i odwołanie t
do funkcji (którą można zastosować do tablic dokładnie tak samo x
), aby zwrócić pojedynczą wartość liczbową; innymi słowy, t
jest statystyką . Generuje wszystkie możliwe próbki z zamianą x
i stosuje się t
do każdej z nich, tworząc w ten sposób jedną liczbę dla każdej takiej próbki: to jest bootstrap w pigułce. Wartość powrotna jest tablicą reprezentujący dokładnej dystrybucji startowej z t
dla próbki x
.
Jako mały przykład , załóżmy średnią dla próbki x
= c(1,3)
:
> boot(c(1,3), mean)
> [1] 1 2 2 3
Istnieją rzeczywiście cztery możliwe próbki wielkości z zamiennikiem z ; mianowicie , , i . generuje je wszystkie (w podanej kolejności) i stosuje się do każdego z nich. W tym przypadku oblicza średnią, a te okazują się odpowiednio , , i , jak pokazano na wyjściu.( 1 , 3 ) ( 1 , 1 ) ( 1 , 3 ) ( 3 , 1 ) ( 3 , 3 ) 1 2 2 32)( 1 , 3 )( 1 , 1 )( 1 , 3 )( 3 , 1 )( 3 , 3 )boot
t
t
12)2)3)
To, dokąd się wybierasz, zależy od tego, jak chcesz użyć bootstrap. Pełne informacje o bootstrap jest zawarte w tej tablicy wyjściowej, więc jest to zwykle dobry pomysł, aby go wyświetlić. Oto przykład, w którym standardowe odchylenie jest z próbki :( 1 , 3 , 3 , 4 , 7 )
hist(boot(c(1,3,3,4,7), sd))
Teraz jesteśmy przygotowani do rozmowy o symulacji Monte Carlo. Załóżmy, powiedzmy, że zamierzamy załadować 95% górny limit ufności dla SD z próbki , używając górnego 95. percentyla jego rozkładu ładowania. Jakie właściwości miałaby ta procedura? Jednym ze sposobów na to jest założenie, że próbka została uzyskana losowo z, powiedzmy, jednolitego rozkładu. (Aplikacja często wskazuje, jakie może być rozsądne założenie dystrybucyjne; tutaj arbitralnie wybrałem takie, które jest proste do obliczeń, ale nie jest łatwe do rozwiązania analitycznego.) Możemy symulować, co się stanie, pobierając taką próbkę i obliczając UCL:5
> set.seed(17)
> quantile(boot(runif(5, min=0, max=10), sd), .95)[1]
95%
3.835870
Wynik dla tej konkretnej próby losowej wynosi 3,83587. Jest to określone: gdybyś zadzwonił boot
ponownie z tym samym zestawem danych, odpowiedź byłaby dokładnie taka sama. Ale jak odpowiedź może się zmienić w przypadku różnych losowych próbek? Dowiedz się, powtarzając ten proces kilka razy i rysując histogram wyników:
> boot.sd <- replicate(100, quantile(boot(runif(5, min=0, max=10), sd), .95)[1])
> hist(boot.sd)
Gdybyśmy mieli wykonać inny zestaw symulacji, losowe losowania byłyby inne, tworząc (nieco) inny histogram - ale nie bardzo różny od tego. Możemy go używać z pewną pewnością, aby zrozumieć, jak działa bootstrap UCL SD. Dla porównania, zauważ, że odchylenie standardowe rozkładu równomiernego (obejmującego zakres od do jak tu określono) wynosi . Jak można się spodziewać, że każdy UCL warty swojej soli, większość (trzy czwarte lub 0,75) wartości na histogramie przekracza to:10 10 / √01010 / 12--√≈ 2,887
> length(boot.sd[boot.sd >= 10/sqrt(12)]) / length(boot.sd)
[1] 0.75
Ale to nie jest blisko nominalnego 95%, które podaliśmy (i mieliśmy nadzieję)! To jedna z wartości symulacji: porównuje nasze nadzieje z tym, co naprawdę się dzieje. (Dlaczego ta rozbieżność? Myślę, że dzieje się tak, ponieważ ładowanie SD nie działa dobrze z naprawdę małymi próbkami.)
Przejrzeć
Statystyki bootstrap są koncepcyjnie takie same jak każda inna statystyka, jak średnia lub odchylenie standardowe; obliczenia zajmują dużo czasu. (Zobacz komunikat ostrzegawczy w boot
kodzie!)
Symulacja Monte-Carlo może być przydatna do badania, jak zmienia się statystyka bootstrap ze względu na losowość w uzyskiwaniu próbek. Zmienność zaobserwowana w takiej symulacji jest spowodowana zmianami w próbkach, a nie zmianami paska ładującego.
(Nie pokazano tutaj) Ponieważ statystyki ładowania początkowego mogą wymagać wielu obliczeń (najwyraźniej do obliczeń dla próbek o rozmiarze ), wygodne jest przybliżenie rozkładu ładowania początkowego. Zwykle odbywa się to poprzez utworzenie programu „czarnej skrzynki” w celu uzyskania losowo jednej wartości z prawdziwej dystrybucji bootstrap i wielokrotne wywoływanie tego programu. Łączna produkcja zbliżona jest do dokładnego rozkładu. Przybliżenie może się różnić ze względu na losowość w czarnej skrzynce - ale ta odmiana jest artefaktem procedury zbliżenia. Nie jest (koncepcyjnie) nierozerwalnie związany z samą procedurą ładowania początkowego. nnnn