Grundlegendes zu MCMC und dem Metropolis-Hastings-Algorithmus

In den letzten Tagen habe ich versucht zu verstehen, wie Markov Chain Monte Carlo (MCMC) funktioniert. Insbesondere habe ich versucht, den Metropolis-Hastings-Algorithmus zu verstehen und zu implementieren. Bisher denke ich, dass ich den Algorithmus allgemein verstehe, aber es gibt ein paar Dinge, die mir noch nicht klar sind. Ich möchte MCMC verwenden, um einige Modelle an Daten anzupassen. Aus diesem Grund werde ich mein Verständnis des Metropolis-Hastings-Algorithmus zum Anpassen einer geraden Linie an einige beobachtete Daten D beschreiben :$f(x)=ax$$D$

1) Machen Sie eine erste Vermutung für . Setzen Sie dies a als unser aktuelles a ( a 0 ). Füge auch ein am Ende der Markov-Kette hinzu ( C ).$a$$a$$a$$a_0$$a$$C$

2) Wiederholen Sie die folgenden Schritte mehrmals.

3) Auswertung gegenwärtiges Likelihood ( ) gegeben ein 0 und D .${\cal L_0}$$a_0$$D$

4) Schlagen Sie ein neues ( a 1 ) durch Abtasten aus einer Normalverteilung mit μ = a 0 und σ = s t e p s i z e vor . Fürs Erste s t e p s i z $a$$a_1$$\mu=a_0$$\sigma=stepsize$$stepsize$ konstant.

5) Berechnen neuer Likelihood ( ) gegeben ein 1 und D .${\cal L_1}$$a_1$$D$

6) Wenn größer als L 0 ist , akzeptiere eine 1 als die neue a ${\cal L_1}$${\cal L_0}$$a_1$$a_0$ und füge sie am Ende von $C$ und fahren Sie mit Schritt 2 fort.

7) Wenn kleiner als L ${\cal L_1}$${\cal L_0}$ erzeugen Sie eine Zahl ( $U$ ) im Bereich [0,1] aus einer gleichmäßigen Verteilung

8) Wenn kleiner ist als die Differenz zwischen den beiden Wahrscheinlichkeiten ( L 1 - L 0 ), akzeptiere eine 1 als die neue a 0 und füge sie am Ende von $U$${\cal L_1}$${\cal L_0}$$a_1$$a_0$$C$ und fahren Sie mit Schritt 2 fort.

9) Wenn größer ist als die Differenz zwischen den beiden Wahrscheinlichkeiten ( L 1 - L 0 ), fügen Sie am Ende von C die a 0 hinzu und verwenden Sie weiterhin dieselbe a $U$${\cal L_1}$${\cal L_0}$$a_0$$C$$a_0$ , und fahren Sie mit Schritt 2 fort.

10) Ende der Wiederholung.

11) Entfernen Sie einige Elemente vom Anfang von $C$ (Einbrennphase).

12) Nehmen Sie nun den Durchschnitt der Werte in . Dieser Durchschnitt ist der geschätzte $C$$a$ .

Jetzt habe ich einige Fragen zu den oben genannten Schritten:

• Wie konstruiere ich die Wahrscheinlichkeitsfunktion für $f(x)=ax$ aber auch für eine beliebige Funktion?
• Ist dies eine korrekte Implementierung des Metropolis-Hastings-Algorithmus?
• Wie kann die Auswahl der Zufallszahlengenerierungsmethode in Schritt 7 die Ergebnisse ändern?
• Wie wird sich dieser Algorithmus ändern, wenn ich mehrere Modellparameter habe? Wenn ich zum Beispiel das Modell .$f(x)=ax+b$

Anmerkungen / Credits: Die Hauptstruktur des oben beschriebenen Algorithmus basiert auf dem Code eines MPIA Python-Workshops.

Es scheint einige Missverständnisse darüber zu geben, was der Metropolis-Hastings (MH) -Algorithmus in Ihrer Beschreibung des Algorithmus ist.

Zunächst muss man verstehen, dass MH ein Abtastalgorithmus ist. Wie angegeben in Wikipedia angegeben

In der Statistik und in der statistischen Physik ist der Metropolis-Hastings-Algorithmus eine Markov-Chain-Monte-Carlo-Methode (MCMC), mit der eine Folge von Zufallsstichproben aus einer Wahrscheinlichkeitsverteilung erhalten wird, für die eine direkte Stichprobe schwierig ist.

$Q(\cdot\vert\cdot)$$f(\cdot)$ , kann der MH-Algorithmus wie folgt implementiert werden:

1. $x_0$ .
2. $x^{\star}$$Q(\cdot\vert x_0)$ .
3. Berechnen Sie das Verhältnis $\alpha=f(x^{\star})/f(x_0)$.
4. Akzeptieren $x^{\star}$ als eine Realisierung von $f$ mit wahrscheinlichkeit $\alpha$.
5. Nehmen $x^{\star}$ als neuen Ausgangszustand und fahren Sie mit der Probenahme fort, bis Sie die gewünschte Probengröße erhalten.

Sobald Sie die Probe erhalten haben, müssen Sie sie noch brennen und verdünnen : Da der Sampler asymptotisch arbeitet, müssen Sie den ersten entfernen$N$ Samples (Burn-In), und da die Samples abhängig sind, müssen Sie jeweils ein Subsampling durchführen $k$ Iterationen (Ausdünnung).

Ein Beispiel in R finden Sie unter folgendem Link:

http://www.mas.ncl.ac.uk/~ndjw1/teaching/sim/metrop/metrop.html

Diese Methode wird hauptsächlich in der Bayes'schen Statistik für die Stichprobe aus der posterioren Verteilung der Modellparameter verwendet.

Das Beispiel, das Sie verwenden, erscheint mir angesichts dessen unklar $f(x)=ax$ ist keine Dichte, es sei denn, Sie beschränken $x$auf einem begrenzten Satz. Mein Eindruck ist, dass Sie daran interessiert sind, eine gerade Linie an eine Reihe von Punkten anzupassen, für die ich empfehlen würde, die Verwendung des Metropolis-Hastings-Algorithmus im Kontext der linearen Regression zu überprüfen. Der folgende Link enthält einige Ideen, wie MH in diesem Kontext verwendet werden kann (Beispiel 6.8):

Robert & Casella (2010), Einführung in Monte-Carlo-Methoden mit R , Ch. 6, "Metropolis-Hastings-Algorithmen"

Auf dieser Site gibt es auch viele Fragen mit Hinweisen auf interessante Referenzen, die sich mit der Bedeutung der Wahrscheinlichkeitsfunktion befassen.

Ein weiterer interessanter Zeiger ist das R-Paket mcmc, das den MH-Algorithmus mit Gaußschen Vorschlägen im Befehl implementiert metrop().