Zadziwiająco skutecznym sposobem na oba przedstawione wyżej kłopoty jest losowanie istotne (Importance Sampling, w skrócie IS). Przypuśćmy, że umiemy losować z rozkładu o gęstości q. Zauważmy, że

gdzie w⁢x=p⁢x/q⁢x. Piszemy tu wzory dla całek ale oczywiście dla sum jest tak samo. Niech X1,…,Xn będą niezależnymi zmiennymi losowymi o jednakowym rozkładzie g,

gdzie

traktujemy jako wagi wylosowanych punktów Xi. Z tego co wyżej powiedzieliśmy wynika, że Eq⁢θ⌃n=θ oraz θ⌃n→θ prawie na pewno, gdy n→∞. Mamy więc estymator nieobciążony i zgodny. Milcząco założyliśmy, że q⁢Xi>0 w każdym wylosowanym punkcie, czyli, że x:p⁢x>0⊆x:q⁢x>0. Z tym zwykle nie ma wielkiego kłopotu. Ponadto musimy założyć, że umiemy obliczać funkcję w w każdym wylosowanym punkcie. Jeśli znamy tylko p′⁢x=Z⁢p⁢x a nie znamy stałej Z to jest kłopot. Możemy tylko obliczyć Wi′=p′⁢Xi/q⁢Xi=Z⁢Wi. Używamy zatem nieco innej postaci estymatora IS, mianowicie

Ten estymator wygląda dziwnie, bo w mianowniku występuje estymator jedynki ale chodzi o to, że we wzorze (8.2) w przeciwieństwie do (8.1) możemy zastąpić Wi przez Wi′, bo nieznany czynnik Z skraca się. Możemy jeszcze zapisać nasz estymator w zgrabnej formie

gdzie W~i=Wi/∑jWj są unormowanymi wagami. Trzeba jednak pamiętać, że to ,,unormowanie” polega na podzieleniu wag przez zmienną losową. W rezultacie estymator IS⁢2 jest obciążony. Jest jednak mocno zgodny, bo licznik we wzorze (8.2) po podzieleniu przez n zmierza do θ a mianownik do 1, prawie na pewno. Estymator IS⁢2 jest znacznie częściej używany niż IS⁢1. Oprócz uniezależnienia się od stałej normującej, jest jeszcze inny powód. Okazuje się, że (pomimo obciążenia) estymator IS⁢2 może być bardziej efektywny. Sens tego pojęcia wyjaśnimy w następnym podrozdziale. Zarówno dobór przykładów jak i ogólny tok rozważań jest zapożyczony z monografii Liu [15].

Naturalne jest żądanie, aby estymator był mocno zgodny, θ⌃n→θ prawie na pewno przy n→∞. Znaczy to, że zbliżymy się dowolnie blisko aproksymowanej wielkości, gdy tylko dostatecznie długo przedłużamy symulacje. Chcielibyśmy jednak wiedzieć coś konkretniejszego, oszacować jak szybko błąd aproksymacji θ⌃n-θ maleje do zera i jakie n wystarczy do osiągnięcia wystarczającej dokładności. Przedstawimy teraz najprostsze i najczęściej stosowane w praktyce podejście, oparte na tzw. asymptotycznej wariancji i konstrukcji asymptotycznych przedziałów ufności. W typowej sytuacji estymator θ⌃n ma następującą własność, nazywaną asymptotyczną normalnością:

w sensie zbieżności według rozkładu. W konsekwencji, dla ustalonej liczby a>0,

gdzie Φ jest dystrybuantą rozkładu N⁢0,1. Często nie znamy asymptotycznej wariancji σ2 ale umiemy skonstruować jej zgodny estymator σ⌃n2. Dla ustalonej ,,małej” dodatniej liczby α łatwo dobrać kwantyl rozkładu normalnego z=z1-α/2 tak, żeby 2⁢Φ⁢z-1=1-α. Typowo α=0.05 i z=z0.975=1.9600≈2. Otrzymujemy

Ten wzór interpretuje się w praktyce tak: estymator θ⌃n ma błąd nie przekraczający z⁢σ⌃n/n z prawdopodobieństwem około 1-α. W żargonie statystycznym 1-α jest nazywane asymptotycznym poziomem ufności.

Chociaż opisame powyżej podejście ma swoje słabe strony, to prowadzi do prostego kryterium porównywania estymatorów. Przypuśćmy, że mamy dwa estymatory θ⌃nI i θ⌃nII, oba asymptotycznie normalne, o asymptotycznej wariancji σI2 i σII6 odpowiednio. Przypuśćmy dalej, że dla obliczenia pierwszego z tych estymatorów generujemy nI2 punktów, zaś dla drugiego nI2. Błędy obu estymatorów na tym mym poziomie istotności są ograniczone przez, odpowiednio, z⁢σI/nI i z⁢σII/nII. Przyrównując te wyrażenia do siebie dochodzimy do wniosku, że oba estymatory osiągają podobną dokładność, jeśli nI/nII=σII2/σI2. Liczbę

nazywamy względną efektywnością (asymptotyczną). Czasami dobrze jest wybrać za ,,naturalny punkt odniesienia” estymator CMC i zdefiniować efektywność estymatora θ⌃n o asymptotycznej wariancji σ2 jako

Mówi się też, że jeśli wygenerujemy próbkę n punktów i obliczymy θ⌃n to ,,efektywna liczność próbki” (ESS, czyli effective sample size) jest n/eff⁢θ⌃n. Tyle bowiem należałoby wygenerować punktów stosując CMC, żeby osiągnąć podobną dokładność.

Asymptotyczna normalność prostego estymatora CMC wynika wprost z Centralnego Twierdzenia Granicznego (CTG). Istotnie, jeśli generujemy niezależnie X1,…,Xn o jednakowym rozkładzie π, to

gdzie

Zupełnie podobnie, dla losowania istotnego w formie (8.1) otrzymujemy asymptotyczną normalność i przy tym

Z tego wzorku widać, że estymator może mieć wariancję zero jeśli q⁢x∝f⁢x⁢p⁢x. Niestety, żeby obliczyć q⁢x potrzebna jest znajomość współczynnika proporcjonalności, który jest równy…θ. Nie wszystko jednak stracone. Pozostaje ważna reguła heurystyczna:

Gęstość q⁢x należy tak dobierać, aby jej ,,kształt” był zbliżony do funkcji f⁢x⁢p⁢x.

Dla drugiej wersji losowania istotnego, (8.2), asymptotyczna normalność wynika z następujących rozważań:

gdzie

Wyrażenie w kwadratowym nawiasie może być ujemne, jeśli jest duża dodatnia korelacja zmiennych w⁢X i w⁢X⁢f⁢X. W tej sytuacji estymator IS2 jest lepszy od IS1. Okazuje się więc rzecz na pozór paradoksalna: dzielenie przez estymator jedynki może poprawić estymator. Poza tym oba estymatory IS1 i IS2 mogą mieć mniejszą (asymptotyczną) wariancję, niż CMC. Jeśli efektywność jest większa niż 100%, to używa się czasem określenia ,,estymator superefektywny”.

Jeśli interesuje nas obliczenie wartości oczekiwanej dla wielu różnych funkcji f, to warto wprowadzić uniwersalny, niezależny of f wskaźnik efektywności. Niezłym takim wskaźnikiem jest Varq⁢w⁢X. Istotnie, ,,odchylenie χ2” pomiędzy gęstością instrumentalną p i docelową p, jest określone poniższym wzorem:

Niestety, to ,,odchylenie” nie jest odległością, bo nie spełnia warunku symetrii. Niemniej widać, że małe wartości Varq⁢w⁢X świadczą o bliskości q i p. Zauważmy, jeszcze, że Varq⁢w⁢X jest wariancją asymptotyczną estymatora IS1 dla funkcji f⁢x=1.

Zgodnym estymatorem wielkości χ2 jest

Dzielenie przez n-1 (a nie n) wynika z tradycji. Dzielenie przez W¯2 powoduje, że estymator nie wymaga znajomości czynnika normującego gęstość p. Jest to w istocie estymator typu IS2. Estymator typu IS1 jest równy po prostu ∑1-Wi2/n, ale można go stosować tylko gdy znamy czynnik normujący.

Wprowadźmy teraz bardziej formalnie pojęcie ważonej próbki. Niech X,W będzie parą zmiennych losowych o wartościach w X×[0,∞[. Jeśli rozkład brzegowy X ma gęstość q zaś docelowa gęstość jest postaci p⁢x=p′⁢x/Z to żądamy aby

dla prawie wszystkich x. Mówimy wtedy, że próbka jest dopasowana do rozkładu p. Idea jest taka jak w losowaniu istotnym IS2: dopasowanie gwarantuje, że dla każdej funkcji f mamy Eq⁢f⁢X⁢W=Z⁢Ep⁢f⁢X. Nie jest konieczne, aby W było deterministyczną funkcją X. Pozostawiamy sobie możliwość generowania wag losowo.

Stała normująca Z jest na ogół nieznana, a więc waga W jest określona z dokładnością do proporcjonalności. W poprzednich podrozdziałach dla podkreślenia opatrywaliśmy takie nieunormowane wagi znaczkiem ,,prim” ale teraz to pomijamy. Gdy generujemy ciąg ważonych zmiennych losowych Xi,Wi,…,Xn,Wn, to oczywiście wymaga się aby E(Wi|Xi=x)=Zp(x)/q(x), gdzie stała Z musi być jednakowa dla wszystkich i=1,…,n. Jeśli założy się niezależność par Xi,Wi, to zgodność i asymptotyczną normalność estymatora θ⌃nIS⁢2 danego równaniem (8.2) wykazuje się tak samo jak poprzednio. Asymptotyczna wariancja jest równa Z-2⁢Varq⁢W⁢f⁢X-θ; czynnik Z-2 pojawia się dlatego, że operujemy nieunormowanymi wagami.

Jest jednak dużo ciekawych zastosowań, w których punkty Xi są zależne. Wtedy nawet zgodność estymatora θ⌃nIS⁢2 nie jest automatyczna. Asymptotyczna normalność może zachodzić z graniczną wariancją zupełnie inną niż w przypadku niezależnym lub nie zachodzić wcale.

Interesujące jest powiązanie idei eliminacji z ważeniem próbek. Czysta metoda eliminacji pracuje w sytuacji, gdy gęstość ,,instrumentalna” majoryzuje funkcję proporcjonalną do gęstości docelowej, czyli Z⁢p⁢x=p′⁢x⩽q⁢x. Jeśli ten warunek nie jest spełniony, to można naprawić odpowiednio ważąc wylosowaną próbkę. Dokładniej, algorytm ważonej eliminacji (Rejection Control, w skrócie RC) jest następujący.

Zauważmy jeszcze, jak należy modyfikować wagi w RC, jeśli na wejściu mamy próbkę ważoną, dopasowaną do p. Jeżeli punkt Y, pochodzący z rozkładu q, ma na wejściu wagę WY=p′⁢Y/q⁢Y, to na wyjściu zaakceptowany punkt X=Y ma rozkład ∝p′⁢x∧q⁢x i powinien mieć wagę WX=p′⁢Y/p′⁢Y∧q⁢Y. Widać stąd, że

gdzie a(⋅) jest napisanym wyżej prawdopodobieństwem akceptacji w RC.