Intuitive Anwendung von Monte-Carlo-Methoden mit Markov-Ketten

Ursprünglicher Autor: Rahul Agarwal
  • Übersetzung

Ist es einfach Ich habe es versucht


Alexey Kuzmin, Leiter der Entwicklung und Arbeit mit den Daten „DomKlik“ Lehrer Richtung Daten Wissenschaft in Netologii, übersetzten Artikel Rahul Agarwal, wie die Monte Carlo Markov - Ketten Methoden laufen für Probleme mit einem großen Zustandsraum zu lösen.

Jeder, der mit Data Science in Verbindung steht, hat mindestens einmal von Monte-Carlo-Methoden mit Markov-Ketten (MCMCs) gehört. Manchmal wird das Thema beim Studium der Bayes'schen Statistik angesprochen, manchmal bei der Arbeit mit Werkzeugen wie Prophet.

Aber MCMC ist schwer zu verstehen. Jedes Mal, wenn ich über diese Methoden las, bemerkte ich, dass sich die Essenz von MCMC in den tiefen Schichten des mathematischen Rauschens verbirgt und es schwierig ist, hinter diesem Rauschen etwas zu erkennen. Ich musste viele Stunden aufwenden, um dieses Konzept zu verstehen.

In diesem Artikel wird versucht, Monte-Carlo-Methoden mit Markov-Ketten zu erklären, damit klar wird, wofür sie verwendet werden. In meinem nächsten Beitrag werde ich mich auf einige weitere Möglichkeiten konzentrieren, diese Methoden zu verwenden.

Also fangen wir an. MCMC besteht aus zwei Begriffen: Monte-Carlo- und Markov-Ketten. Sprechen wir über jeden von ihnen.

Monte Carlo




Im einfachsten Sinne können Monte-Carlo-Methoden als einfache Simulationen definiert werden.

Monte Carlo Methoden haben ihren Namen vom Monte Carlo Casino in Monaco. Bei vielen Kartenspielen müssen Sie die Gewinnwahrscheinlichkeit des Dealers kennen. Manchmal kann die Berechnung dieser Wahrscheinlichkeit mathematisch kompliziert oder schwierig zu lösen sein. Wir können jedoch immer eine Computersimulation ausführen, um das gesamte Spiel viele Male zu spielen, und die Wahrscheinlichkeit als die Anzahl der Siege geteilt durch die Anzahl der gespielten Spiele betrachten.
Dies ist alles, was Sie über Monte-Carlo-Methoden wissen müssen. Ja, es ist nur eine einfache Modellierungstechnik mit einem ausgefallenen Namen.

Markov-Ketten




Da der Begriff MCMC aus zwei Teilen besteht, müssen Sie noch verstehen, was Markov-Ketten sind . Aber bevor wir zu den Markov-Ketten übergehen, wollen wir uns ein wenig mit den Markov-Eigenschaften befassen.

Angenommen, es gibt ein System von M-möglichen Zuständen, und Sie wechseln von einem Zustand in einen anderen. Lass dich noch nicht verwirren. Ein spezielles Beispiel für ein solches System ist das Wetter, das sich von heiß zu kalt zu mäßig ändert. Ein weiteres Beispiel ist der Aktienmarkt, der von einem bärischen in einen bullischen und stagnierenden Zustand übergeht.

Die Markov-Eigenschaft legt nahe, dass für einen gegebenen Prozess, der sich zu einem bestimmten Zeitpunkt im Zustand X n befindet  , die Wahrscheinlichkeit X n + 1 ist= k (wobei k einer der M-Zustände ist, in die der Prozess gehen kann) hängt nur davon ab, wie dieser Zustand im Moment ist. Und nicht darüber, wie es seinen gegenwärtigen Zustand erreicht hat.
In mathematischer Sprache können wir dies in der Form der folgenden Formel schreiben: Der

Klarheit halber ist Ihnen die Abfolge der Bedingungen, die der Markt angenommen hat, um „bullisch“ zu werden, egal. Die Wahrscheinlichkeit, dass der nächste Zustand „bärisch“ ist, wird nur dadurch bestimmt, dass sich der Markt derzeit in einem „bullischen“ Zustand befindet. Das macht auch in der Praxis Sinn.

Ein Prozess mit einer Markov-Eigenschaft wird als Markov-Prozess bezeichnet. Warum ist die Markov-Kette wichtig? Aufgrund seiner stationären Verteilung.

Was ist stationäre Verteilung?


Ich werde versuchen, die stationäre Verteilung zu erklären, indem ich sie für das folgende Beispiel berechne. Angenommen, Sie haben einen Markov-Prozess für die Börse, wie unten gezeigt.

Sie haben eine Übergangswahrscheinlichkeitsmatrix, die die Wahrscheinlichkeit eines Übergangs von Zustand X i nach X j bestimmt .

Matrix der Übergangswahrscheinlichkeiten, Q

In der gegebenen Matrix der Übergangswahrscheinlichkeiten Q ist die Wahrscheinlichkeit, dass der nächste Zustand "bull" ist, bei dem gegenwärtigen Zustand "bull" = 0,9; die Wahrscheinlichkeit, dass der nächste Zustand "bärisch" ist, wenn der aktuelle Zustand "bull" ist = 0,075. Usw.

Fangen wir mit einem bestimmten Zustand an. Unser Staat wird durch den Vektor [Bulle, Bär, Stagnation] bestimmt. Wenn wir mit dem bärischen Zustand beginnen, ist der Vektor wie folgt: [0,1,0]. Wir können die Wahrscheinlichkeitsverteilung für den nächsten Zustand berechnen, indem wir den aktuellen Zustandsvektor mit der Übergangswahrscheinlichkeitsmatrix multiplizieren.

Beachten Sie, dass die Wahrscheinlichkeiten eine Summe von 1 ergeben.

Die folgende Verteilung von Zuständen kann durch die Formel gefunden werden:


Und so weiter. Am Ende erreichen Sie einen stationären Zustand, in dem sich der Zustand stabilisiert:


Für die oben beschriebene Übergangswahrscheinlichkeitsmatrix Q lautet die stationäre Verteilung s wie folgt:

Sie erhalten die stationäre Verteilung mit folgendem Code:

Q = 
np.matrix([[0.9,0.075,0.025],[0.15,0.8,0.05],[0.25,0.25,0.5]])
init_s = np.matrix([[0, 1 , 0]])
epsilon =1
while epsilon>10e-9:
   next_s = np.dot(init_s,Q)
   epsilon = np.sqrt(np.sum(np.square(next_s - init_s)))
   init_s = next_s
print(init_s)
------------------------------------------------------------------
matrix([[0.62499998, 0.31250002, 0.0625    ]])

Sie können auch von jedem anderen Zustand aus starten - erreichen Sie die gleiche stationäre Verteilung. Ändern Sie den Ausgangszustand im Code, wenn Sie dies überprüfen möchten.

Jetzt können wir die Frage beantworten, warum die stationäre Verteilung so wichtig ist.

Die stationäre Verteilung ist wichtig, weil damit zufällig bestimmt werden kann, mit welcher Wahrscheinlichkeit sich ein System in einem bestimmten Zustand befindet.

In unserem Beispiel können wir sagen, dass sich der Markt in 62,5% der Fälle in einem „bullischen“ Zustand befindet, 31,25% in einem „bearischen“ Zustand und 6,25% in einer Stagnation.

Intuitiv kann man dies als zufälliges Umherwandern der Kette sehen.


Zufälliger Spaziergang

Sie befinden sich an einem bestimmten Punkt und wählen den nächsten Zustand aus, wobei Sie die Wahrscheinlichkeitsverteilung des nächsten Zustands unter Berücksichtigung des aktuellen Zustands beobachten. Wir können einige Knoten öfter als andere besuchen, basierend auf den Wahrscheinlichkeiten dieser Knoten.

So löste Google das Suchproblem zu Beginn des Internets. Das Problem bestand darin, die Seiten nach ihrer Wichtigkeit zu sortieren. Google hat das Problem mithilfe des Pagerank-Algorithmus gelöst. Der Google Pagerank-Algorithmus sollte den Status als Seite und die Wahrscheinlichkeit einer Seite in einer stationären Verteilung als relative Bedeutung betrachten.

Nun wenden wir uns direkt der Betrachtung von MCMC-Methoden zu. 

Was sind Monte-Carlo-Methoden mit Markov-Ketten (MCMC)?


Lassen Sie mich vor der Beantwortung von MCMC eine Frage stellen. Wir kennen uns mit der Betaverteilung aus. Wir kennen seine Wahrscheinlichkeitsdichtefunktion. Aber können wir aus dieser Distribution eine Stichprobe ziehen? Können Sie einen Weg finden, dies zu tun?


Denken Sie ...

MCMC gibt Ihnen die Möglichkeit, aus einer beliebigen Wahrscheinlichkeitsverteilung zu wählen. Dies ist besonders wichtig, wenn Sie eine Auswahl aus der posterioren Verteilung treffen müssen.

Die Abbildung zeigt den Bayes'schen Satz

: Sie müssen beispielsweise eine Stichprobe aus der posterioren Verteilung erstellen. Aber ist es einfach, die hintere Komponente zusammen mit der Normalisierungskonstante (Evidenz) zu berechnen? In den meisten Fällen können Sie sie in Form eines Wahrscheinlichkeitsprodukts und einer A-priori-Wahrscheinlichkeit finden. Die Normierungskonstante (p (D)) kann jedoch nicht berechnet werden. Warum? Schauen wir uns das genauer an.

Angenommen, H nimmt nur 3 Werte an:

p (D) = p (H = H1) .p (D | H = H1) + p (H = H2) .p (D | H = H2) + p (H = H3 ) .p (D | H = H3)

In diesem Fall ist p (D) einfach zu berechnen. Was aber, wenn der Wert von H stetig ist? Wäre es möglich, dies so einfach zu berechnen, insbesondere wenn H unendliche Werte annehmen würde? Dazu müsste ein komplexes Integral gelöst werden.

Wir wollen eine zufällige Auswahl aus der posterioren Verteilung treffen, wollen aber auch p (D) als Konstante betrachten.

Wikipedia schreibt :

Monte-Carlo-Methoden mit Markov-Ketten sind eine Klasse von Algorithmen zur Abtastung aus einer Wahrscheinlichkeitsverteilung, die auf dem Aufbau einer Markov-Kette basiert, die als stationäre Verteilung die gewünschte Form hat. Der Kettenzustand nach einer Reihe von Schritten wird dann als Auswahl aus der gewünschten Verteilung verwendet. Die Abtastqualität verbessert sich mit zunehmender Anzahl von Schritten.

Schauen wir uns ein Beispiel an. Angenommen, Sie benötigen ein Beispiel aus der Betaverteilung . Seine Dichte:


wobei C die Normalisierungskonstante ist. Eigentlich ist dies eine Funktion von α und β, aber ich möchte zeigen, dass es für eine Stichprobe aus der Beta-Verteilung nicht benötigt wird, also nehmen wir es als Konstante.

Das Problem der Betaverteilung ist sehr schwierig, wenn nicht fast unlösbar. In der Realität müssen Sie möglicherweise mit komplexeren Verteilungsfunktionen arbeiten, und manchmal kennen Sie die Normalisierungskonstanten nicht.

MCMC-Methoden vereinfachen das Leben, indem sie Algorithmen bereitstellen, mit denen eine Markov-Kette mit einer Beta-Verteilung als stationäre Verteilung erstellt werden kann, vorausgesetzt, wir können aus einer (relativ einfachen) gleichmäßigen Verteilung wählen.

Wenn wir mit einem zufälligen Zustand beginnen und mehrmals auf der Grundlage eines Algorithmus in den nächsten Zustand übergehen, werden wir schließlich eine Markov-Kette mit einer Beta-Verteilung als stationäre Verteilung erstellen. Und die Zustände, in denen wir uns seit langem befinden, können als Beispiel aus der Beta-Distribution verwendet werden.

Einer dieser MCMC-Algorithmen ist der Metropolis-Hastings-Algorithmus.

Metropolis-Hastings-Algorithmus




Intuition:


Wozu also?

Intuitiv möchten wir ein Stück Oberfläche (unsere Markov-Kette) so entlanglaufen, dass die Zeit, die wir an jedem Ort verbringen, proportional zur Höhe der Oberfläche an diesem Ort ist (die gewünschte Wahrscheinlichkeitsdichte, aus der wir eine Auswahl treffen möchten).

Zum Beispiel möchten wir auf einem Hügel mit einer Höhe von 100 Metern doppelt so viel Zeit verbringen wie auf einem benachbarten Hügel mit einer Höhe von 50 Metern. Es ist gut, dass wir dies auch dann tun können, wenn wir die absoluten Höhen der Punkte auf der Oberfläche nicht kennen: Alles, was Sie wissen müssen, sind die relativen Höhen. Wenn zum Beispiel die Spitze des Hügels A zweimal höher ist als die Spitze des Hügels B, möchten wir doppelt so viel Zeit in A verbringen wie in B.

Es gibt komplexere Schemata, um neue Orte und Regeln für ihre Annahme vorzuschlagen, aber die Hauptidee lautet wie folgt:

  1. Wählen Sie einen neuen "vorgeschlagenen" Ort.
  2. Finden Sie heraus, wie viel höher oder niedriger dieser Ort im Vergleich zum aktuellen ist.
  3. Mit einer Wahrscheinlichkeit proportional zur Höhe der Orte an Ort und Stelle bleiben oder an einen neuen Ort ziehen.

Der Zweck der MCMC besteht darin, aus einer Wahrscheinlichkeitsverteilung auszuwählen, ohne dass die genaue Höhe an einem beliebigen Punkt bekannt sein muss (C muss nicht bekannt sein).
Wenn der „Wandervorgang“ korrekt eingerichtet ist, können Sie sicherstellen, dass diese Proportionalität (zwischen der aufgewendeten Zeit und der Verteilungshöhe) erreicht wird
.

Algorithmus:


Definieren und beschreiben wir die Aufgabe nun in formelleren Begriffen. Sei s = (s1, s2, ..., sM) die gewünschte stationäre Verteilung. Wir wollen eine Markov-Kette mit einer solchen stationären Verteilung erstellen. Wir beginnen mit einer beliebigen Markov-Kette mit M-Zuständen mit der Übergangsmatrix P, so dass pij die Wahrscheinlichkeit des Übergangs von Zustand i nach j darstellt.

Intuitiv wissen wir, wie man die Markov-Kette durchstreift, aber die Markov-Kette hat nicht die erforderliche stationäre Verteilung. Diese Kette hat eine stationäre Verteilung (die wir nicht brauchen). Unser Ziel ist es, die Art und Weise zu ändern, in der wir uns in der Markov-Kette bewegen, damit die Kette die gewünschte stationäre Verteilung aufweist.

Um dies zu tun:

  1. Beginnen Sie mit einem zufälligen Ausgangszustand i.
  2. Wählen Sie einen neuen angenommenen Zustand zufällig aus, indem Sie die Übergangswahrscheinlichkeiten in der i-ten Zeile der Übergangsmatrix P betrachten.
  3. Berechnen Sie ein Maß namens Entscheidungswahrscheinlichkeit, das wie folgt definiert ist: aij = min (sj.pji / si.pij, 1).
  4. Wirf nun eine Münze, die mit der Wahrscheinlichkeit aij auf der Oberfläche des Adlers landet. Wenn ein Adler fällt, nehmen Sie das Angebot an, gehen Sie zum nächsten Status über, andernfalls lehnen Sie das Angebot ab, dh bleiben Sie im aktuellen Status.
  5. Wiederholen Sie viele Male.

Nach einer Vielzahl von Tests konvergiert diese Kette und hat eine stationäre Verteilung s. Dann können wir Kettenzustände als Beispiel für jede Verteilung verwenden.

Wenn Sie dies tun, um die Beta-Verteilung abzutasten, müssen Sie die Wahrscheinlichkeitsdichte nur verwenden, um nach der Wahrscheinlichkeit zu suchen, eine Entscheidung zu treffen. Teilen Sie dazu sj durch si (dh die Normalisierungskonstante C wird aufgehoben).

Beta-Auswahl




Kommen wir nun zum Problem des Samplings aus der Betaverteilung.

Eine Beta-Verteilung ist eine kontinuierliche Verteilung für [0,1] und kann für [0,1] unendlich viele Werte haben. Angenommen, eine beliebige Markov-Kette P mit unendlichen Zuständen auf [0,1] hat eine Übergangsmatrix P mit pij = pji = alle Elemente in der Matrix.

Wir brauchen keine Matrix P, wie wir später sehen werden, aber ich möchte, dass die Beschreibung des Problems dem von uns vorgeschlagenen Algorithmus so nahe wie möglich kommt.

  • Beginnen Sie mit einem zufälligen Ausgangszustand i, der sich aus einer Gleichverteilung auf (0,1) ergibt.
  • Wählen Sie einen neuen angenommenen Zustand zufällig aus, indem Sie die Übergangswahrscheinlichkeiten in der i-ten Zeile der Übergangsmatrix P betrachten. Nehmen wir an, wir haben einen anderen Zustand Unif (0,1) als angenommenen Zustand j gewählt.
  • Berechnen Sie die Kennzahl, die als Entscheidungswahrscheinlichkeit bezeichnet wird:


Was vereinfacht:

Da pji = pij, und wo

  • Wirf jetzt eine Münze. Mit großer Wahrscheinlichkeit wird ein Adler fallen. Wenn ein Adler fällt, sollten Sie das Angebot annehmen, dh in den nächsten Zustand wechseln. Ansonsten lohnt es sich, das Angebot abzulehnen, dh im gleichen Zustand zu bleiben.
  • Wiederholen Sie den Test viele Male.

Code:


Es ist Zeit, von der Theorie zur Praxis überzugehen. Lassen Sie uns unser Beta-Beispiel in Python schreiben.

impo
rt
rand
om
# Lets define our Beta Function to generate s for any particular state. We don't care for the normalizing constant here.
def beta_s(w,a,b):
	return w**(a-1)*(1-w)**(b-1)
# This Function returns True if the coin with probability P of heads comes heads when flipped.
def random_coin(p):
	unif = random.uniform(0,1)
	if unif>=p:
    	return False
	else:
    	return True
# This Function runs the MCMC chain for Beta Distribution.
def beta_mcmc(N_hops,a,b):
	states = []
	cur = random.uniform(0,1)
	for i in range(0,N_hops):
    	states.append(cur)
    	next = random.uniform(0,1)
    	ap = min(beta_s(next,a,b)/beta_s(cur,a,b),1) # Calculate the acceptance probability
    	if random_coin(ap):
        	cur = next
	return states[-1000:] # Returns the last 100 states of the chain

Vergleichen Sie die Ergebnisse mit der tatsächlichen Betaverteilung. 

impo
rt
num
py
as
np
import pylab as pl
import scipy.special as ss
%matplotlib inline
pl.rcParams['figure.figsize'] = (17.0, 4.0)
# Actual Beta PDF.
def beta(a, b, i):
	e1 = ss.gamma(a + b)
	e2 = ss.gamma(a)
	e3 = ss.gamma(b)
	e4 = i ** (a - 1)
	e5 = (1 - i) ** (b - 1)
	return (e1/(e2*e3)) * e4 * e5
# Create a function to plot Actual Beta PDF with the Beta Sampled from MCMC Chain.
def plot_beta(a, b):
	Ly = []
	Lx = []
	i_list = np.mgrid[0:1:100j]
	for i in i_list:
    	Lx.append(i)
        Ly.append(beta(a, b, i))
    pl.plot(Lx, Ly, label="Real Distribution: a="+str(a)+", b="+str(b))
    pl.hist(beta_mcmc(100000,a,b),normed=True,bins =25, histtype='step',label="Simulated_MCMC: a="+str(a)+", b="+str(b))
	pl.legend()
	pl.show()
plot_beta(0.1, 0.1)
plot_beta(1, 1)
plot_beta(2, 3)




Wie Sie sehen, sind die Werte der Betaverteilung sehr ähnlich. Somit hat das MCMC-Netzwerk einen stationären Zustand erreicht.Im

obigen Code haben wir einen Beta-Sampler erstellt, aber das gleiche Konzept gilt für alle anderen Distributionen, aus denen wir eine Auswahl treffen möchten.

Schlussfolgerungen




Es war ein großartiger Beitrag. Herzlichen Glückwunsch, wenn Sie es bis zum Ende gelesen haben.

Im Wesentlichen können MCMC-Methoden komplex sein, bieten uns jedoch eine große Flexibilität. Sie können eine beliebige Verteilungsfunktion auswählen, indem Sie diese über die MCMC auswählen. Typischerweise werden diese Methoden verwendet, um Proben aus posterioren Verteilungen zu entnehmen.

Sie können MCMC auch verwenden, um Probleme mit einem großen Zustandsraum zu lösen. Zum Beispiel bei einem Rucksackproblem oder zur Entschlüsselung. Ich werde versuchen, Ihnen im nächsten Beitrag weitere interessante Beispiele zu liefern . Bleib dran.

Von den Redakteuren



Jetzt auch beliebt: