PRML Kapitel 14 Bedingte gemischte Modell-Python-Implementierung

Da ein Modell die Eigenschaften von Daten nicht erfassen kann, wird ein gemischtes Modell verwendet, wenn Sie mehrere Modelle kombinieren möchten. Dieses Mal werden wir eine Mischung aus linearen Regressionsmodellen implementieren.

Mischen linearer Regressionsmodelle

Probabilistische Darstellung eines normalen linearen Regressionsmodells mit Eingabe $ x $ und Ausgabe $ t $,

p(t|x,{\bf w},\beta) = \mathcal{N}(t|{\bf w}^\top\phi(x),\beta^{-1})

Es sieht aus wie. $ \ Phi (\ cdot) $ ist jedoch der Merkmalsvektor, $ {\ bf w} $ ist der Gewichtungskoeffizient und $ \ beta $ ist der Präzisionsparameter. Im gemischten linearen Regressionsmodell führen wir den Mischungskoeffizienten $ \ pi_k $ (wobei $ \ sum_k \ pi_k = 1 $) ein und addieren die K linearen Regressionsmodelle, um sie wie folgt auszudrücken.

p(t|x,{\bf W},\beta,\pi) = \sum_{k=1}^K\pi_k\mathcal{N}(t|{\bf w}_k^\top\phi(x),\beta^{-1})

Wie üblich wird der Parameter $ {\ bf W angegeben, wenn die Trainingsdaten $ \ {x_n, t_n \} _ {n = 1} ^ N $ (im Folgenden als $ {\ bf t, x} $ bezeichnet) angegeben werden. }, \ pi, \ beta $ wird höchstwahrscheinlich geschätzt. Die logarithmische Wahrscheinlichkeitsfunktion zu diesem Zeitpunkt ist

\ln p({\bf t}|{\bf x},{\bf W},\pi,\beta) = \sum_{n=1}^N\ln\left(\sum_{k=1}^K\pi_k\mathcal{N}(t_n|{\bf w}_k^\top\phi(x_n),\beta^{-1})\right)

Es wird sein. Nun führen wir eine latente Variable $ z_ {nk} $ ein, die darstellt, aus welchem Modell die einzelnen Daten generiert wurden. Diese latente Variable hat den Wert 0 oder 1, und wenn $ z_ {nk} = 1 $ ist, bedeutet dies, dass der n-te Datenpunkt aus dem k-ten Modell $ 1, \ cdots, k-1 generiert wird. , K + 1, \ cdots, K $ sind 0. Jetzt können Sie auch die logarithmische Wahrscheinlichkeitsfunktion für die Parameter bei vollständigen Daten aufschreiben.

Auf diese Weise kann der EM-Algorithmus für die wahrscheinlichste Schätzung des gemischten linearen Regressionsmodells verwendet werden. Suchen Sie in Schritt E $ \ mathbb {E} [z_ {nk}] $ für jeden Datenpunkt und jedes Modell und in Schritt M nach der Tatsache über die latente Variable $ z $ der logarithmischen Wahrscheinlichkeitsfunktion unter Berücksichtigung der vollständigen Daten. Maximiert den erwarteten Wert der Verteilung für den Parameter.

Code

import Dieses Mal werden wir zusätzlich zu den üblichen numpy und matplotlib auch die Standard-Python-Bibliothek itertools importieren. (Ich brauche es nicht wirklich)

import itertools
import matplotlib.pyplot as plt
import numpy as np

Planungsmatrix

Erstellen Sie eine Planungsmatrix, in der die Zeilenvektoren jeden Merkmalsvektor darstellen.

class PolynomialFeatures(object):

    def __init__(self, degree=2):
        assert isinstance(degree, int)
        self.degree = degree

    def transform(self, x):
        if x.ndim == 1:
            x = x[:, None]
        x_t = x.transpose()
        features = [np.ones(len(x))]
        for degree in range(1, self.degree + 1):
            for items in itertools.combinations_with_replacement(x_t, degree):
                features.append(reduce(lambda x, y: x * y, items))
        return np.asarray(features).transpose()

Bedingtes gemischtes Modell

Dies ist der Code für das gemischte lineare Regressionsmodell, das den Hauptteil darstellt.

class ConditionalMixtureModel(object):

    def __init__(self, n_models=2):
        #Geben Sie die Anzahl der Modelle an
        self.n_models = n_models

    #Methode zur Schätzung der wahrscheinlichsten
    def fit(self, X, t, n_iter=100):
        #Initialisierung der zu schätzenden Parameter (Gewichtskoeffizient, Mischungskoeffizient, Varianz)
        coef = np.linalg.pinv(X).dot(t)
        self.coef = np.vstack((
            coef + np.random.normal(size=len(coef)),
            coef + np.random.normal(size=len(coef)))).T
        self.var = np.mean(np.square(X.dot(coef) - t))
        self.weight = np.ones(self.n_models) / self.n_models

        #Wiederholen Sie den EM-Schritt eine bestimmte Anzahl von Malen
        for i in xrange(n_iter):
            #E Stufenbelastungsrate E.[z_nk]Wird für jede Daten und jedes Modell berechnet
            resp = self._expectation(X, t)

            #Maximieren Sie ungefähr M Schrittparameter
            self._maximization(X, t, resp)

    #Gaußsche Funktion
    def _gauss(self, x, y):
        return np.exp(-(x - y) ** 2 / (2 * self.var))

    #E Schritt Belastungsrate Gamma_nk=E[z_nk]Berechnung
    def _expectation(self, X, t):
        #PRML-Formel(14.37)
        responsibility = self.weight * self._gauss(X.dot(self.coef), t[:, None])
        responsibility /= np.sum(responsibility, axis=1, keepdims=True)
        return responsibility

    #Maximieren Sie ungefähr M Schrittparameter
    def _maximization(self, X, t, resp):
        #PRML-Formel(14.38)Berechnung des Mischungskoeffizienten
        self.weight = np.mean(resp, axis=0)

        for m in xrange(self.n_models):
            #PRML-Formel(14.42)Berechnen Sie den Koeffizienten für jedes Modell
            self.coef[:, m] = np.linalg.inv((X.T * resp[:, m]).dot(X)).dot(X.T * resp[:, m]).dot(t)

        #PRML-Formel(14.44)Berechnung der Varianz
        self.var = np.mean(
            np.sum(resp * np.square(t[:, None] - X.dot(self.coef)), axis=1))

    def predict(self, X):
        return X.dot(self.coef)

Trainingsdaten

y=\|x\|Bereiten Sie die Trainingsdaten gemäß vor.x\ge0Wannx\le0Sie können mit einem linearen Regressionsmodell in jeden Bereich von passen, es ist jedoch besser, zwei Modelle so vorzubereiten, dass sie in beide Bereiche passen.

def create_toy_data(sample_size=20, std=0.1):
    x = np.linspace(-1, 1, sample_size)
    y = (x > 0.).astype(np.float) * 2 - 1
    y = x * y
    y += np.random.normal(scale=std, size=sample_size)
    return x, y

Hauptfunktion

def main():
    #Erstellung von Trainingsdaten
    x_train, y_train = create_toy_data()

    #Definition des Merkmalsvektors (Reihenfolge 1)
    feature = PolynomialFeatures(degree=1)

    #Umstellung auf Planungsmatrix
    X_train = feature.transform(x_train)

    #Vorbereitung eines bedingten gemischten Modells (2 Modelle)
    model = ConditionalMixtureModel(n_models=2)

    #Führen Sie die wahrscheinlichste Parameterschätzung durch
    model.fit(X_train, y_train)

    #Illustrierte Ergebnisse
    x = np.linspace(-1, 1, 100)
    X = feature.transform(x)
    Y = model.predict(X)
    plt.scatter(x_train, y_train, s=50, facecolor="none", edgecolor="g")
    plt.plot(x, Y[:, 0], c="b")
    plt.plot(x, Y[:, 1], c='r')
    plt.show()

Ganzer Code

import itertools
import matplotlib.pyplot as plt
import numpy as np


class PolynomialFeatures(object):

    def __init__(self, degree=2):
        assert isinstance(degree, int)
        self.degree = degree

    def transform(self, x):
        if x.ndim == 1:
            x = x[:, None]
        x_t = x.transpose()
        features = [np.ones(len(x))]
        for degree in range(1, self.degree + 1):
            for items in itertools.combinations_with_replacement(x_t, degree):
                features.append(reduce(lambda x, y: x * y, items))
        return np.asarray(features).transpose()


class ConditionalMixtureModel(object):

    def __init__(self, n_models=2):
        self.n_models = n_models

    def fit(self, X, t, n_iter=100):
        coef = np.linalg.pinv(X).dot(t)
        self.coef = np.vstack((
            coef + np.random.normal(size=len(coef)),
            coef + np.random.normal(size=len(coef)))).T
        self.var = np.mean(np.square(X.dot(coef) - t))
        self.weight = np.ones(self.n_models) / self.n_models

        for i in xrange(n_iter):
            resp = self._expectation(X, t)
            self._maximization(X, t, resp)

    def _gauss(self, x, y):
        return np.exp(-(x - y) ** 2 / (2 * self.var))

    def _expectation(self, X, t):
        responsibility = self.weight * self._gauss(X.dot(self.coef), t[:, None])
        responsibility /= np.sum(responsibility, axis=1, keepdims=True)
        return responsibility

    def _maximization(self, X, t, resp):
        self.weight = np.mean(resp, axis=0)
        for m in xrange(self.n_models):
            self.coef[:, m] = np.linalg.inv((X.T * resp[:, m]).dot(X)).dot(X.T * resp[:, m]).dot(t)
        self.var = np.mean(
            np.sum(resp * np.square(t[:, None] - X.dot(self.coef)), axis=1))

    def predict(self, X):
        return X.dot(self.coef)


def create_toy_data(sample_size=20, std=0.1):
    x = np.linspace(-1, 1, sample_size)
    y = (x > 0.).astype(np.float) * 2 - 1
    y = x * y
    y += np.random.normal(scale=std, size=sample_size)
    return x, y


def main():
    x_train, y_train = create_toy_data()
    feature = PolynomialFeatures(degree=1)
    X_train = feature.transform(x_train)
    model = ConditionalMixtureModel(n_models=2)
    model.fit(X_train, y_train)
    x = np.linspace(-1, 1, 100)
    X = feature.transform(x)
    Y = model.predict(X)
    plt.scatter(x_train, y_train, s=50, facecolor="none", edgecolor="g")
    plt.plot(x, Y[:, 0], c="b")
    plt.plot(x, Y[:, 1], c='r')
    plt.show()


if __name__ == '__main__':
    main()

Ergebnis

Nichtlineare Funktion durch Kombination mehrerer linearer Modelley=\|x\|Die Anpassung an die Daten gemäß ist abgeschlossen.

Am Ende

Wenn das diesmal implementierte bedingte gemischte Modell zur Vorhersage verwendet wird, erstreckt sich die Linie auch dort, wo keine Datenpunkte vorhanden sind, wie im obigen Ergebnis gezeigt. Das gemischte Expertenmodell, das eine Weiterentwicklung des bedingten gemischten Modells darstellt, scheint dieses Problem zu lösen, indem der Mischungskoeffizient $ \ pi $ eine Funktion der Eingabe $ x $ ist.