[PYTHON] Versuchen Sie es mit einer probabilistischen Programmiersprache (Pyro).

Die probabilistische Programmiersprache wird für die statistische Modellierung nach Bayes verwendet. Ich denke, dass die Methode der Verwendung von Stan, BUGS und JAGS in der R-Sprache üblich ist, da es viele Bücher gibt. Zu den probabilistischen Programmierbibliotheken von Python gehören PyMC, PyStan und Edward (jetzt TensorFlow Probability), eine Erweiterung von TensorFlow. Pyro wurde 2018 veröffentlicht. Pyro ist eine probabilistische Programmierbibliothek, die auf PyTorch basiert und von Uber entwickelt wurde. Vor kurzem hat das Unternehmen auch NumPyro basierend auf Jax veröffentlicht. Pyros Artikel: http://jmlr.org/papers/v20/18-403.html NumPyros Artikel: https://arxiv.org/abs/1912.11554

Pyro und NumPyro können Bayesian Deep Learning problemlos implementieren, und ich denke, sie werden in Zukunft angewendet. (Implementierter Artikel)

Dieses Mal werde ich die Grundlagen aus Pyros offiziellem Tutorial ausprobieren. Ich habe auch auf den folgenden Blog-Artikel verwiesen. https://www.hellocybernetics.tech/entry/2019/12/08/033824 https://www.hellocybernetics.tech/entry/2019/12/08/193220

Installieren Sie Pyro

Es wird davon ausgegangen, dass Python 3.4 oder höher ist und PyTorch installiert ist. Die Umgebung des Autors war Ubuntu 16.04, CUDA10, Anaconda, Python3.6, Pytorch 1.2.0.

pip install pyro-ppl

Pyro 1.3.1 ist installiert. Außerdem wurde PyTorch automatisch auf 1.4.0 aktualisiert.

Implementierung (1): Definition und Stichprobe der Wahrscheinlichkeitsverteilung

Samen importieren und setzen. Für die Reproduzierbarkeit der Probenahme ist ein fester Samen erforderlich. Wenn Sie pyro.set_rng_seed ausführen, wird anscheinend auch der Startwert in pytorch festgelegt.

import torch
import pyro
pyro.set_rng_seed(101)

pyro.distributions kann ähnlich wie torch.distributions sein. Dieses Mal werden wir das Wetter und die Temperatur mit der Bernoulli-Verteilung und der Normalverteilung definieren.

def weather():
    cloudy = pyro.sample('cloudy', pyro.distributions.Bernoulli(0.3))
    cloudy = 'cloudy' if cloudy.item() == 1.0 else 'sunny'
    mean_temp = {'cloudy': 55.0, 'sunny': 75.0}[cloudy]
    scale_temp = {'cloudy': 10.0, 'sunny': 15.0}[cloudy]
    temp = pyro.sample('temp', pyro.distributions.Normal(mean_temp, scale_temp))
    return cloudy, temp.item()

Stapeln Sie dies, um den Verkauf von Eis zu definieren. Mit anderen Worten, es ist ein Wahrscheinlichkeitsmodell der Hierarchie.

def ice_cream_sales():
    cloudy, temp = weather()
    expected_sales = 200. if cloudy == 'sunny' and temp > 80.0 else 50.
    ice_cream = pyro.sample('ice_cream', pyro.distributions.Normal(expected_sales, 10.0))
    return ice_cream

Lassen Sie uns den Verkauf von Eis probieren und in einem Histogramm anzeigen.

import matplotlib.pyplot as plt
N = 1000 #Anzahl von Beispielen
x_list = []
for _ in range(N):
    x = ice_cream_sales()
    x_list.append(x)
plt.hist(x_list)

Eine bimodale Verteilung wurde erhalten.

Vektorisierung mit der pyro.plate-Methode

Wenn Sie die pyro.plate-Methode verwenden, können Sie anscheinend ohne Wiederholung der for-Anweisung durch Vektorisierung eine Stichprobe erstellen. Dies bedeutet, dass die Abtastung parallel von derselben unabhängigen Verteilung durchgeführt wird. Ich habe bisher versucht, den Code mit pyro.plate neu zu schreiben. Modelldefinition

def model(N):
    with pyro.plate("plate", N):
        cloudy = pyro.sample('cloudy', pyro.distributions.Bernoulli(0.3))
        mean_temp = 75.0 - 20.0 * cloudy
        scale_temp = 15.0 - 5.0 * cloudy
        temp = pyro.sample('temp', pyro.distributions.Normal(mean_temp, scale_temp))
        expected_sales = 50.0 + 150.0 * (1 - cloudy) * (temp > 80.0)
        ice_cream = pyro.sample('ice_cream', pyro.distributions.Normal(expected_sales, 10.0))
    return ice_cream

Zeigen Sie das Histogramm an.

import matplotlib.pyplot as plt
N = 1000 #Anzahl von Beispielen
x_list = model(N)
plt.hist(x_list)

Ich habe ein ähnliches Histogramm. Ist der geringfügige Unterschied der Effekt der parallelen Abtastung aufgrund der Vektorisierung?

Als wichtiges Merkmal von pyro wird der Variablenname durch das erste Argument angegeben, z. B. pyro.sample ("hoge", dist). Dies scheint notwendig zu sein, um zu definieren, dass es sich um eine globale stochastische Variable handelt. (Ich bin mir über die Details nicht sicher.)

Implementierung (2): Bayesianische Inferenz

Im vorherigen Abschnitt haben wir herausgefunden, wie Daten aus einem Wahrscheinlichkeitsmodell generiert werden. In diesem Abschnitt werden wir versuchen, aus den Daten ein stochastisches Modell abzuleiten.

import torch
import pyro
pyro.set_rng_seed(101)
import matplotlib.pyplot as plt

Bereiten Sie als Daten insgesamt 10 Daten vor, 1 Wert x 6 und 0 Wert x 4.

data = []
for _ in range(6):
    data.append(torch.tensor(1.0))
for _ in range(4):
    data.append(torch.tensor(0.0))
data = torch.tensor(data) #Machen Sie es zu einem Tensor für die Handhabung mit Pyro
plt.hist(data)

Das Histogramm sieht wie oben aus.

Modellieren Sie eine Wahrscheinlichkeitsverteilung für diese Daten. Dieses Mal betrachten wir die Bernoulli-Verteilung mit dem Parameter "latent_fairness", der der Beta-Verteilung folgt.

import pyro.distributions as dist
def model(data):
    alpha0 = torch.tensor(10.0)
    beta0 = torch.tensor(10.0)
    f = pyro.sample("latent_fairness", dist.Beta(alpha0, beta0))
    with pyro.plate("data_plate"):
        pyro.sample("obs", dist.Bernoulli(f), obs=data)

Schliessen Sie Parameter, die zu diesem Modell passen. Es gibt verschiedene Inferenzmethoden, aber Variant Bayes und MCMC sind die wichtigsten.

Variante Bayes

Unter den Inferenzmethoden scheint Pyro die vollständigsten Funktionen von Variant Bayes zu haben. Dies liegt wahrscheinlich daran, dass die automatischen Differenzierungs- und Optimierungsmethoden von PyTorch, die die Grundlage von Pyro bilden, für die SVI (Probabilistic Variant Inference) von Variant Bayes verwendet werden.

Bereiten Sie zunächst eine variable Approximationsverteilung vor. Sie können es selbst definieren oder automatisch mit der Klasse pyro.infer.autoguide.guides.

guide = pyro.infer.autoguide.guides.AutoDelta(model)

Führt eine Differentialinferenz durch.

from pyro.optim import Adam
from pyro.infer import SVI, Trace_ELBO
adam_params = {"lr": 0.0005, "betas": (0.90, 0.999)}
optimizer = Adam(adam_params)
pyro.clear_param_store()
svi = SVI(model, guide, optimizer, loss=Trace_ELBO())
n_steps = 1000 #Anzahl der Schritte
losses, f  = [], [] #Zur Aufnahme für eine spätere Handlung
for step in range(n_steps):
    loss = svi.step(data)
    losses.append(loss)
    f.append(pyro.param("AutoDelta.latent_fairness").item())

Jetzt können Sie schließen. Zeichnen Sie den Fortschritt.

plt.subplot(2,1,1)
plt.plot(losses)
plt.xlabel("step")
plt.ylabel("loss")
plt.subplot(2,1,2)
plt.plot(f)
plt.xlabel("step")
plt.ylabel("latent_fairness")

Der Parameter "latent_fairness" ist auf ungefähr 0,5357 konvergiert. Der endgültige Wert kann unten angezeigt werden.

for name in pyro.get_param_store():
    print("{}: {}".format(name, pyro.param(name)))

AutoDelta.latent_fairness: 0.535700261592865

Zeigen Sie die vorhergesagte Wahrscheinlichkeitsverteilung aus dem abgeleiteten Modell an. Nehmen Sie 1000 Proben.

from pyro.infer import Predictive
predictive_dist = Predictive(model=model, guide=guide,num_samples=1000, return_sites=["obs"])
plt.hist(predictive_dist.get_samples(None)["obs"].view(1,-1))

Die Verteilung ähnelt dem Histogramm der Daten. (Ist es etwas anders, weil die ursprüngliche Anzahl der Daten nur 10 beträgt?)

MCMC MCMC ist eine Brute-Force-Technik und braucht Zeit.

from pyro.infer import　NUTS, MCMC
nuts_kernel = NUTS(model)
mcmc = MCMC(nuts_kernel, num_samples=5000, warmup_steps=1000)
mcmc.run(data)

Jetzt können Sie schließen. Siehe das Ergebnis.

mcmc.summary()

MCMC kann die vorhergesagte Wahrscheinlichkeitsverteilung ausdrücken, indem eine Stichprobe von Parametern an posterior_samples übergeben wird.

f_mcmc = mcmc.get_samples()["latent_fairness"]
from pyro.infer import Predictive
predictive_dist = Predictive(model=model, posterior_samples={"latent_fairness": f_mcmc}, return_sites=["obs"])
plt.hist(predictive_dist.get_samples(None)["obs"].view(1,-1))

Dies hat auch eine Verteilung ähnlich dem Histogramm der Daten.

Andere Inferenzmethoden

Es gibt MAP-Schätzung und Laplace-Approximation. Die MAP-Schätzung ist eine spezielle Form der Bayes'schen Inferenz, und da es sich um eine Punktschätzung handelt, ist eine Wahrscheinlichkeitsverteilung nicht erforderlich. Es ist keine Schlussfolgerung, sondern eine Schätzung, da keine Wahrscheinlichkeitsverteilung erforderlich ist. (Referenzartikel) *** Die oben durchgeführte Variation Bayes entspricht der MAP-Schätzung, da die Delta-Verteilung als Variations-Approximationsverteilung ausgewählt wird. ** ** ** Die Laplace-Approximation scheint mit der AutoLaplace-Approximation des Autoguides möglich zu sein.

Anwendungen

Es scheint, dass Sie verschiedene Dinge tun können, wie das gemischte Gaußsche Modell, den Gaußschen Prozess, die Bayes'sche Optimierung, NMF, den Kalman-Filter, das versteckte Markov-Modell und so weiter.

Am Ende

Ich fange gerade erst mit Pyro an, aber irgendwie habe ich herausgefunden, wie man es benutzt. Ich denke, es ist eine wunderbare Technik, eine Theorie implementieren zu können, die mit mathematischen Formeln so präzise schwierig ist.