Jak należy podejść do problemu Euler Projektu 213 („Flea Circus”)?

Chciałbym rozwiązać Project Euler 213, ale nie wiem od czego zacząć, ponieważ jestem laikiem w dziedzinie statystyki, zauważ, że wymagana jest dokładna odpowiedź, aby metoda Monte Carlo nie zadziałała. Czy mógłbyś polecić mi kilka tematów statystycznych do przeczytania? Proszę nie zamieszczać rozwiązania tutaj.

Flea Circus

Siatka kwadratów 30 × 30 zawiera 900 pcheł, początkowo jedną pchłę na kwadrat. Kiedy dzwonek jest dzwoniony, każda pchła losowo skacze na sąsiedni kwadrat (zwykle 4 możliwości, z wyjątkiem pcheł na krawędzi siatki lub w rogach).

Jaka jest oczekiwana liczba niezajętych kwadratów po 50 dzwonkach dzwonu? Podaj odpowiedź w zaokrągleniu do sześciu miejsc po przecinku.

Masz rację; Monte Carlo jest niewykonalne. (W naiwnej symulacji - tzn. Takiej, która dokładnie odtwarza sytuację problemową bez żadnych uproszczeń - każda iteracja wymagałaby 900 ruchów pcheł. Szacunkowy procent proporcji pustych komórek wynosi , co sugeruje wariancję Monte - Oszacowanie Carlo po takich iteracjach wynosi około Aby określić odpowiedź z do sześciu miejsc po przecinku, należy ją oszacować z dokładnością do 5.E -7 i, aby osiągnąć poziom ufności 95 +% (powiedzmy), musiałbyś w przybliżeniu zmniejszyć o połowę tę precyzję do 2,5E-7. Rozwiązanie dajeN 1 / N 1 / e ( 1 - 1 / e ) = 0,2325 … / N $1/e$$N$$1/N 1/e (1 - 1/e) = 0.2325\ldots /N$N>4E$\sqrt(0.2325/N) \lt 2.5E-7$$N > 4E12$około To byłoby około 3,6E15 ruchów pcheł, z których każdy wymagałby kilku kliknięć procesora. Z jednym dostępnym nowoczesnym procesorem będziesz potrzebował całego roku (wysoce wydajnego) przetwarzania. I nieco niepoprawnie i nadmiernie optymistycznie założyłem, że odpowiedź jest podana w postaci proporcji zamiast liczby: jako liczba będzie potrzebowała jeszcze trzech znaczących liczb, co pociągnie za sobą milion-krotny wzrost obliczeń ... Czy możesz długo czekać?)

Jeśli chodzi o rozwiązanie analityczne, dostępne są pewne uproszczenia. (Można ich również użyć do skrócenia obliczeń Monte Carlo). Oczekiwana liczba pustych komórek jest sumą prawdopodobieństw pustki we wszystkich komórkach. Aby to znaleźć, możesz obliczyć rozkład prawdopodobieństwa liczby zajętości każdej komórki. Rozkłady te uzyskuje się poprzez zsumowanie (niezależnego!) Wkładu każdej pchły. Zmniejsza to problem ze znalezieniem liczby ścieżek o długości 50 wzdłuż siatki 30 na 30 między dowolną parą komórek na tej siatce (jedna jest początkiem pcheł, a druga komórką, dla której chcesz obliczyć prawdopodobieństwo obłożenie pcheł).

Czy nie możesz iterować po prawdopodobieństwie zajęcia komórek dla każdej pcheł. Oznacza to, że pchła k jest początkowo w komórce (i (k), j (k)) z prawdopodobieństwem 1. Po 1 iteracji ma prawdopodobieństwo 1/4 w każdej z 4 sąsiednich komórek (zakładając, że nie jest na krawędzi ani w róg). Następnie, w następnej iteracji, każda z tych ćwiartek zostaje „rozmazana” z kolei. Po 50 iteracjach masz macierz prawdopodobieństw okupacji dla pcheł k. Powtórz ponad 900 pcheł (jeśli skorzystasz z symetrii, zmniejszy to prawie o współczynnik 8) i dodaj prawdopodobieństwa (nie musisz przechowywać ich wszystkich naraz, tylko macierz bieżącej pchły (hmm, chyba że jesteś bardzo sprytne, możesz potrzebować dodatkowej działającej macierzy) i bieżącej sumy macierzy). Wydaje mi się, że istnieje wiele sposobów na przyspieszenie tego tu i tam.

Nie wymaga to żadnej symulacji. Wymaga to jednak sporo obliczeń; nie powinno być bardzo trudno obliczyć rozmiar symulacji wymagany do uzyskania odpowiedzi z nieco lepszą dokładnością niż 6 dp z dużym prawdopodobieństwem i dowiedzieć się, które podejście będzie szybsze. Oczekuję, że takie podejście przewyższy symulację o pewien margines.

Chociaż nie sprzeciwiam się praktycznej niemożności (lub niepraktyczności) rozwiązania tego problemu z Monte Carlo z dokładnością do 6 miejsc po przecinku wskazanej przez whubera , sądzę , że można uzyskać rozdzielczość z sześciocyfrową dokładnością.

Po pierwsze, zgodnie z Glen_b , cząstki są wymienialne w trybie stacjonarnym, a zatem wystarczające (jak w wystarczającym stopniu ) jest monitorowanie zajętości różnych komórek, ponieważ stanowi to również proces Markowa. Rozkład zajętości w następnym kroku jest zakończony, określony przez zajętości w bieżącym czasie t . Napisanie macierzy przejścia K jest zdecydowanie niepraktyczne, ale symulacja przejścia jest prosta.$t+1$$t$$K$

Po drugie, jak zauważył shabbychef można śledzić proces użytkowanie na 450 nieparzyste (parzyste) lub kwadratów, która pozostaje na nieparzystych kwadratów kiedy tylko rozważa nawet razy, czyli kwadratu Markowa macierzy .$K^2$

Po trzecie, oryginalny problem uważa tylko częstotliwości zerowej po 50 przejściach Markowa. Uwzględniając fakt, że punkt początkowy ma bardzo duże znaczenie dla stacjonarnego rozkładu prawdopodobieństwa łańcucha Markowa ( X ( t ) ) oraz fakt, że koncentrują się na pojedynczej średnią dla wszystkich komórek, p 0 = $\hat{p}_0$$50$$(\mathbf{X}^{(t)})$możemy uznać, że realizacja łańcucha(X(t))w czasiet=50jest realizacją ze stacjonarnego rozkładu prawdopodobieństwa. Zapewnia to znaczną redukcję kosztów obliczeniowych, ponieważ możemy symulować bezpośrednio z tego rozkładu stacjonarnegoπ, który jest rozkładem wielomianowym z prawdopodobieństwami proporcjonalnymi do 2, 3 i 4 w parzystym narożniku, innymi komórkami na krawędzi i komórkami wewnętrznymi odpowiednio.

pot=rep(c(rep(c(0,1),15),rep(c(1,0),15)),15)*c(2,
    rep(3,28),2,rep(c(3,rep(4,28),3),28),2,rep(3,28),2)
pot=pot/sum(pot)
sum((1-pot)^450)-450
[1] 166.1069

$166.11$

Jak skomentował Whuber , szacunki należy pomnożyć przez 2, aby poprawnie odpowiedzieć na pytanie, stąd ostateczna wartość 332,2137,

Podejście analityczne może być nużące i nie zastanawiałem się nad zawiłościami, ale oto podejście, które warto rozważyć. Ponieważ interesuje Cię oczekiwana liczba komórek pustych po 50 pierścieniach, musisz zdefiniować łańcuch markowa nad „liczbą pcheł w komórce”, a nie pozycję pcheł (patrz odpowiedź Glen_b, która modeluje pozycję pchła jako łańcuch markowa. Jak zauważył Andy w komentarzach do tej odpowiedzi, takie podejście może nie osiągnąć tego, czego chcesz.)

W szczególności pozwól:

$n_{ij}(t)$$i$$j$

Następnie łańcuch markowa zaczyna się od następującego stanu:

$n_{ij}(0) =1$$i$$j$

Ponieważ pchły przenoszą się do jednej z czterech sąsiednich komórek, stan komórki zmienia się w zależności od liczby pcheł znajdujących się w komórce docelowej i liczby pcheł w czterech sąsiednich komórkach oraz prawdopodobieństwa, że ​​zostaną przeniesione do tej komórki. Za pomocą tej obserwacji możesz zapisać prawdopodobieństwo przejścia stanu dla każdej komórki w zależności od stanu tej komórki i stanu sąsiednich komórek.

Jeśli chcesz, mogę rozszerzyć odpowiedź, ale to wraz z podstawowym wprowadzeniem do łańcuchów Markowa powinno zacząć.

jeśli masz zamiar iść drogą numeryczną, prosta obserwacja: problem wydaje się podlegać czerwono-czarnej parzystości (pchła na czerwonym kwadracie zawsze przesuwa się na czarny kwadrat i odwrotnie). Może to pomóc zmniejszyć rozmiar problemu o połowę (rozważ tylko dwa ruchy na raz i, na przykład, patrz tylko na pchły na czerwonych kwadratach).

Podejrzewam, że pewna znajomość łańcuchów Markowa w dyskretnym czasie może okazać się przydatna.