Równowaga Kwasowo-Zasadowa

Odpowiadanie na pytanie poprzez rzucenie kośćmi

Metoda Monte Carlo zastępuje obliczenie, którego nie można wykonać za pomocą eksperymentu, który można powtórzyć. Zapisz wartość, którą chcesz oszacować jako wartość oczekiwaną, wyciągnij losowe próbki z odpowiedniej rozkładu i oblicz średnią. Prawo wielkich liczb gwarantuje, że średnia zbiega się do prawdziwej wartości; twierdzenie o granicy centralnej mówi Ci, jak szybko.

Pomysł narodził się w Los Alamos w latach 40. XX wieku. Stanislaw Ulam, powracając do zdrowia po chorobie, próbował obliczyć prawdopodobieństwo, że gra solitaire zakończy się danym sposobem — i zdał sobie sprawę, że łatwiej jest rozdać sto kart i zliczać niż enumerate kombinatorykę. John von Neumann natychmiast dostrzegł, że ta sama sztuczka dotyczy transportu neutronów, a Nicholas Metropolis podsunął nazwę kodową, nawiązującą do kasyna w Monako, gdzie wuj Ulama lubił grać.

π z kwadratu i koła

Klasyczna demonstracja, oraz ta uruchamiana w symulacji na tej stronie, polega na oszacowaniu π. Weź kwadrat jednostkowy i ćwiartkę koła wpisane w niego. Losowo rozłóż punkty w kwadracie. Prawdopodobieństwo, że punkt trafi do środka ćwiartki koła, jest równe stosunkowi pól, czyli π/4. Zlicz trafienia, podziel przez rzuty, pomnóż przez cztery.

let inside = 0; for (let i = 0; i < n; i++) { const x = Math.random(), y = Math.random(); if (x * x + y * y <= 1) inside++; // nie potrzeba sqrt } const piEstimate = 4 * inside / n; Nic w tym nie jest specyficzne dla okręgów. Każde całkowanie to pole, a każde pole to prawdopodobieństwo, że punkt równomiernie rozłożony wpada pod krzywą — więc ten sam trzy-liniowy pętlę jest uniwersalnym integratorem. Jedynym, co się zmienia z problemu na problem, jest to, czego próbujemy pobrać i ile liczymy.

let inside = 0;
for (let i = 0; i < n; i++) {
  const x = Math.random(), y = Math.random();
  if (x * x + y * y <= 1) inside++;      // no sqrt needed
}
const piEstimate = 4 * inside / n;

Prawo 1/√N i dlaczego jest zarówno straszne, jak i wspaniałe

Każde rzut jest próbą Bernoulliego z prawdopodobieństwem sukcesu p = π/4 ≈ 0.7854, więc jego wariancja wynosi p(1 − p) ≈ 0.1685. Estymator jest średnią z N takich prób, a zatem jego odchylenie standardowe wynosi √(p(1−p)/N), a szacunek π — czterokrotność tego — ma typowy błąd:

σ(π̂) = 4 · sqrt(p(1 − p) / N) ≈ 1.64 / sqrt(N) N = 10⁴ → typowy błąd ≈ 0.016 (około 2 poprawne cyfry) N = 10⁶ → typowy błąd ≈ 0.0016 (około 3 poprawne cyfry) N = 10⁸ → typowy błąd ≈ 0.00016 (około 4 poprawne cyfry) To jest właściwość definiująca Monte Carlo: błąd spada w miarę 1/√N, więc zakup kolejnego ułamka dziesiętnego kosztuje setki razy więcej próbek. Jako sposób na obliczenie π jest to bezsensowne — szeregi podobne do Machina dają więcej cyfr w kilkudziesięciu wyrazach niż ta pętla w miliardzie rzutów. Jednakże nie ma w tym wzorze jednej rzeczy: wymiaru. Deterministyczna metoda kwadratury na siatce potrzebuje mᵈ punktów, aby umieścić m punktów wzdłuż każdego z d osi, a jej dokładność zawodzi wraz ze wzrostem d — przekleństwo wymiarowości. 1/√N Monte Carlo nie zależy od d w ogóle. W dwóch wymiarach jest to żart; w dwudziestu lub stu dwudziestu wymiarach integralnego transportu światła, wypłaty finansowej lub funkcji podziału mechaniki statystycznej, często jest jedynym działającym rozwiązaniem. Każde śledzenie ścieżek, każdy kalkulator opcji i każdy próbnik posteriora Bayesa opiera się na tym jednym spostrzeżeniu.

σ(π̂) = 4 · sqrt(p(1 − p) / N) ≈ 1.64 / sqrt(N)

N = 10⁴    → typical error ≈ 0.016     (about 2 correct digits)
N = 10⁶    → typical error ≈ 0.0016    (about 3 correct digits)
N = 10⁸    → typical error ≈ 0.00016   (about 4 correct digits)

Redukcja wariancji: jedyna dźwignia, którą masz

Ponieważ wykładnik w 1/√N jest stały, jedynym sposobem na uzyskanie lepszej odpowiedzi za tę samą cenę jest zmniejszenie stałej przed nią — wariancji tego, co średnim.

Przykłady czterech standardowych technik:

Sampling ważony próbuj, gdzie funkcja integracyjna jest duża, a następnie dziel przez gęstość próbkowania (wagę). Klasycznym przypadkiem jest próbkowanie źródła światła bezpośrednio zamiast celowania wiązkami w większość czarnego półkuli.

Sampling warstwowy podziel domenę na K warstw i próbuj każdej, zapobiegając skupiskom i przerwom czystej losowości.

Wariaty przeciwstawne dla każdego przykładu u, użyj również 1 − u. Odwrócony korelacja wymaza część błędu.

Wariaty kontrolowane odejmij skorelowaną ilość, której dokładna średnia jest już znana, i dodaj tę średnią z powrotem.

importance sampling   sample where the integrand is large, then divide
                      by the sampling density (the weight). The classic
                      case: sampling a light source directly instead of
                      firing rays into a mostly black hemisphere.

stratified sampling   split the domain into K strata and sample each,
                      forbidding the clumps and gaps of pure randomness.

antithetic variates   for every sample u, also use 1 − u. The negative
                      correlation cancels part of the error.

control variates      subtract a correlated quantity whose exact mean
                      you already know, and add that mean back.

Kiedy nie można bezpośrednio próbować z rozkładu

Wszystko powyżej zakłada, że można pobierać niezależne próbki z rozkładu docelowego. W statystyce bayesowskiej i fizyce statystycznej zazwyczaj tego nie można: zna się tylko gęstość do niewiadomego stałego normalizującego. Odpowiedź to Monte Carlo o zmiennym stopniu – zbuduj losowy spacer, którego długoterminowa rozkład stacjonarny jest celem, a następnie średnia wartość wzdłuż spaceru. Zasada Metropolis–Hastings jest zaskakująco prosta: zaproponuj ruch i zaakceptuj go z prawdopodobieństwem

α = min(1, π(x_nowe) / π(x_przestepne)) // symetryczny proponowany zaakceptuj → łańcuch się porusza rzeknieć → łańcuch pozostaje w tym miejscu i ponownie liczy się x_przestepne stały normalizujący wygasa w stosunku, który jest całym punktem. Koszt polega na tym, że kolejne próbki są skorelowane, więc N MCMC kroków przenosi mniej informacji niż N niezależnych próbek, a należy odrzucić początkowy okres wypalania, podczas którego łańcuch zapomina, gdzie zaczął się.

α = min(1, π(x_new) / π(x_old))       // symmetric proposal

accept  → the chain moves
reject  → the chain stays where it is, and you count x_old AGAIN

Losy numery same w sobie

Wynik Monte Carlo jest dobry tylko do tak dobrego generatora. Trzy praktyczne zasady. Po pierwsze, użyj dobrze przetestowanego generatora liczb pseudolosowych — nowoczesnego opartego na liczniku lub z rodziny xorshift, a nie domowego liniowego konwencjonalnego, którego dolne bity słynnie są nieszandomowe i których punkty słynnie padają na hiperplane. Po drugie, zainicjalizuj go i zarejestruj wartość inicjalizacji: uruchomienie Monte Carlo, które nie może być odtworzone, nie może być debugowane. Trzeci, obserwuj okres i korelacje między strumieniami, jeśli równolegle — dwa pracownicy przypadkowo dzielący się strumieniem to cicha, katastrofalna usterka, której żaden test nie wykryje, ponieważ odpowiedź nadal wygląda na prawdopodobną.

Użycie Math.random() w przeglądarce jest w porządku dla wizualizacji takiej jak ta (nie jest inicjalizowane i nie jest kryptograficznie silne, ale jest szybkie i jego jakość statystyczna jest wystarczająca); do uruchomienia naukowego użyj generatora, który można zainicjalizować i odtworzyć.

Frequently asked questions

Wypróbuj na żywo

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