[PYTHON] Gaußscher Prozess mit pymc3

Ich werde versuchen, anhand des Gaußschen Prozesses anhand der auf dieser Website enthaltenen Daten zu Shampoo-Verkäufen zu schätzen.

https://machinelearningmastery.com/time-series-datasets-for-machine-learning/

Der Gaußsche Prozess ist ein typisches nichtparametrisches Regressionsmodell. Es hat den Vorteil, dass es gegenüber nichtlinearen Daten stark ist und die Schätzungsunsicherheit bestätigen kann.

import os
import pandas as pd
shampoo_df = pd.read_csv('../data/shampoo.csv')
shampoo_df
shampoo_df.plot.line(x='Month', y='Sales',c="k");

Im Folgenden wird die Grenzwahrscheinlichkeit berechnet, um eine angemessene Wahrscheinlichkeit und ein vorheriges Verteilungspaar zu erhalten.

with pm.Model() as shampoo_model:
    ρ=pm.HalfCauchy('ρ', 1) #Je größer die Längenskala ist, desto näher scheinen die Bewegungen zu sein.<= Selbstkorreliert
    η = pm.HalfCauchy('η', 1) #Je größer die Längenskala ist, desto näher scheinen die Bewegungen zu sein.<= Selbstkorreliert
    
    M = pm.gp.mean.Linear(coeffs=(shampoo_df.index/shampoo_df.Sales).mean()) #Ursprünglicher Wert
    K= (η**2) * pm.gp.cov.ExpQuad(1, ρ)

    σ = pm.HalfCauchy('σ', 1)
    gp = pm.gp.Marginal(mean_func=M, cov_func=K) 
    gp.marginal_likelihood("recruits", X=shampoo_df.index.values.reshape(-1,1), 
                           y=shampoo_df.Sales.values, noise=σ)
    trace = pm.sample(1000)

pm.traceplot(trace);

Die Schätzung der meisten Wahrscheinlichkeiten erzeugt ein gut angepasstes Modell durch Maximieren der Wahrscheinlichkeitsfunktion p (X | θ) in Bezug auf den Parameter θ, verwendet jedoch nicht die vorherige Verteilung p (θ). Andererseits scheint der Fall, in dem eine vorherige Verteilung vorliegt, als maximale posteriore Wahrscheinlichkeitsschätzung bezeichnet zu werden.

Weitere Informationen finden Sie in den folgenden Büchern. Atsushi Suyama. Startup-Serie zum maschinellen Lernen Einführung in das maschinelle Lernen durch Bayesianische Inferenz (japanische Ausgabe) (Kindle Locations 2154-2156). Kindle Edition.

Erhalten Sie die optimalen Parameter, indem Sie die maximale hintere Wahrscheinlichkeit schätzen.

with shampoo_model:
    marginal_post = pm.find_MAP()

{'ρ_log__': array(2.88298779), 'η_log__': array(5.63830367), 'σ_log__': array(4.13876033), 'ρ': array(17.86757799), 'η': array(280.98567051), 'σ': array(62.72501496)}

X_pred=np.array([r for r in range(X_pred.max()+10)]).reshape(-1,1)

with shampoo_model:
    shampoo_pred = gp.conditional("shampoo_", X_pred)
    samples=pm.sample_ppc([marginal_post], vars=[shampoo_pred] ,samples=100)



fig = plt.figure(figsize=(10,5)); 
ax = fig.gca()

from pymc3.gp.util import plot_gp_dist
plot_gp_dist(ax, samples["shampoo_"], X_pred);

ax.plot(shampoo_df.index, shampoo_df.Sales, 'dodgerblue', ms=5, alpha=0.5);
ax.plot(shampoo_df.index, shampoo_df.Sales, 'ok', ms=5, alpha=0.5);
plt.show()

In pm.sample_ppc wird eine Stichprobe aus der posterioren Vorhersageverteilung basierend auf den Parametern in marginal_post durchgeführt.