Estymacja Pi Monte Carlo: Co Losowe Punkty Mogą i Nie Mogą Powiedzieć

Oszacowanie powierzchni, którą nie da się zintegrować, poprzez rzucanie dartami

Narysuj okrąg jednostkowy wpisany w kwadrat 2×2. Pole tego okręgu wynosi π; pole kwadratu to 4. Rozsyp punkty równomiernie losowo wewnątrz kwadratu i policz, jaki ułamek ląduje wewnątrz okręgu – geometrycznie ten ułamek powinien zbiegać się do π/4. Pomnóż obserwowany ułamek przez 4 i uzyskasz oszacowanie π zbudowane wyłącznie z liczenia, bez kalkulusu. Jest to kanoniczny przykład integracji Monte Carlo: przekształcanie prawdopodobieństwa geometrycznego w odpowiedź numeryczną poprzez powtarzające się losowe próbkowanie, nazwane na cześć kasyna Monte Carlo, gdzie fizycy z Los Alamos w latach 40-tych używali tej samej sztuczki do oszacowania całek dyfuzji neutronów, które nie miały zamkniętej postaci.

let inside = 0, total = 0;
function samplePoint() {
  const x = Math.random() * 2 - 1;   // uniform in [-1, 1]
  const y = Math.random() * 2 - 1;
  total++;
  if (x * x + y * y <= 1) inside++;   // inside the unit circle?
  return 4 * inside / total;           // running estimate of pi
}

Dlaczego zbiega się jak 1/√n, a nie szybciej

Każdy próbek jest próbą Bernoulliego: trafia on do wnętrza okręgu z prawdopodobieństwem p = π/4 ≈ 0.7854, lub znajduje się poza nim z prawdopodobieństwem 1 − p. Obserwowana część trafnych po n próbach jest bezstronnym oszacowaniem p, a zgodnie z twierdzeniem o rozkładzie normalnym jej odchylenie standardowe maleje do σ/√n, gdzie σ² = p(1 − p) jest ustalone przez geometrię. W praktyce oznacza to, że błąd spada tylko o około 3.16× (√10) za każdym razem, gdy pobierasz dziesięć razy więcej próbek – aby uzyskać jeden dodatkowy poprawny cyfrę π, potrzebujesz w przybliżeniu 100× większej liczby punktów. Ta szybkość 1/√n jest definiującą, nieuniknioną cechą prostych metod Monte Carlo i dlatego wygląda ona niewiele, obok deterministycznych szeregów dla π (Leibniz, formuły Machina-podobne lub algorytm Chudnowskiego), które zbiegają się wykładniczo szybko wraz z liczbą wyrazów.

Gdzie naprawdę wygrywa Monte Carlo

Metoda tablicy wędrownych jest złą metodą do obliczania π, a to tylko dlatego, że jest to specyficzne. Wiele klasycznych szeregów przewyższa ją po kilku terminach. Jej prawdziwą wartość ujawnia się w odwrotnej sytuacji: w wysokowymiarowych całkach, gdzie siatki kwadratury deterministycznej stają się niemożliwe do wykonania obliczeniowo. Metoda z siatką punktową o 100 punktach na osi wymaga 100^d oceny dla całki wymiaru d – przy d = 10 jest to 10^20 punktów, co jest bezcelowe. Błąd Monte Carlo wciąż maleje wraz z 1/√n niezależnie od wymiaru, więc milion losowych próbek daje używalne oszacowanie przy d = 10, gdzie żadna rozsądna siatka nie mogła tego zrobić. To właśnie dlatego metody Monte Carlo dominują w finansach obliczeniowych (ceny opcji zależnych od ścieżki), fizyce statystycznej (model Isinga, QCD na macierzy) i całkach transportu światła wewnątrz fizycznie opartego renderowania.

Redukcja wariancji: uzyskanie większej dokładności na próbkę

Ponieważ błąd jest zdeterminowany przez σ/√n, a n jest kosztowne, większość praktycznych ulepszeń metod Monte Carlo koncentruje się na zmniejszaniu σ zamiast zwiększania n. Próbkowanie stratyfikowane dzieli kwadrat na siatkę i gwarantuje co najmniej jeden punkt próbkowania na komórce, unikając nieprzychylnego skupiania się, które może czasami generować czyste, równomierne próbkowanie losowe. Ważne próbkowanie pobiera więcej próbek tam, gdzie funkcja całki wnosi największy wkład (nieistotne dla kwadratu jednostkowego, ale krytyczne, gdy funkcja do zintegrowania jest ostro skośna). Quasi-Monte Carlo zastępuje prawdziwą losowość sekwencjami o niskiej dyspersji, takimi jak Sobol lub Halton, które wypełniają przestrzeń bardziej równomiernie niż punkty i.i.d. losowe i mogą podnosić szybkość zbieżności do prawie 1/n dla gładkich całek — asymptoticznie znacznie lepiej niż 1/√n, kosztem utraty prostych progów błędów statystycznych.

Sprawdzenie rozsądności: standardowy błąd w skrócie

Przy p = π/4, σ = √(p(1−p)) ≈ 0.4105, więc standardowy błąd oszacowania π (które jest 4-krotnym szacunkiem udziału trafień) po n próbach wynosi około 4 × 0.4105 / √n ≈ 1.642/√n. Przy n = 10 000 przewiduje się błąd rzędu ±0,016 — co jest zgodne z tym, co zobaczymy zmieniającego się na żywym liczniku nawet po tysiącach punktów, oraz pokazuje, że "więcej próbek" przynosi dokładność powoli, a nie za darmo.

Frequently asked questions

Wypróbuj na żywo

Wszystko powyżej działa bezpośrednio w Twojej przeglądarce — otwórz Monte Carlo Pi Estimation i zmieniaj parametry podczas działania. Nic nie jest instalowane ani przesyłane na serwer, cały model działa w jednej karcie.