Puisqu'un modèle ne peut pas capturer les caractéristiques des données, un modèle mixte est utilisé lorsque vous souhaitez combiner plusieurs modèles. Cette fois, nous allons mettre en œuvre un mélange de modèles de régression linéaire.

Mélange de modèles de régression linéaire

Représentant de manière probabiliste un modèle de régression linéaire normal avec entrée $ x $ et sortie $ t $,

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

On dirait. Cependant, $ \ phi (\ cdot) $ est le vecteur de caractéristiques, $ {\ bf w} $ est le coefficient de pondération et $ \ beta $ est le paramètre de précision. Dans le modèle de régression linéaire mixte, nous introduisons le coefficient de mélange $ \ pi_k $ (où $ \ sum_k \ pi_k = 1 $) et additionnons les K modèles de régression linéaire pour l'exprimer comme suit.

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

Comme d'habitude, le paramètre $ {\ bf W lorsque les données d'entraînement $ \ {x_n, t_n \} _ {n = 1} ^ N $ (ci-après dénommé $ {\ bf t, x} $) est donnée. }, \ pi, \ beta $ est le plus susceptible d'être estimé. La fonction de vraisemblance logarithmique à ce moment est

\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)

Ce sera. Nous introduisons maintenant une variable latente $ z_ {nk} $ qui représente à partir de quel modèle les données individuelles ont été générées. Cette variable latente prend une valeur de 0 ou 1, et si $ z_ {nk} = 1 $, cela signifie que le nième point de données est généré à partir du kième modèle, $ 1, \ cdots, k-1. , K + 1, \ cdots, K $ valent 0. Maintenant, vous pouvez également écrire la fonction de vraisemblance de type journal pour les paramètres étant donné les données complètes.

De cette manière, l'algorithme EM peut être utilisé pour l'estimation la plus probable du modèle de régression linéaire mixte. À l'étape E, trouvez $ \ mathbb {E} [z_ {nk}] $ pour chaque point de données et modèle, et à l'étape M, après coup sur la variable latente $ z $ de la fonction de vraisemblance logarithmique étant donné les données complètes. Maximise la valeur attendue de la distribution pour le paramètre.

code

import Cette fois, en plus des habituels numpy et matplotlib, nous importerons également la bibliothèque python standard, itertools. (Je n'en ai pas vraiment besoin)

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

Matrice de planification

Créez une matrice de planification dans laquelle les vecteurs de ligne représentent chaque vecteur d'entités.

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()

Modèle mixte conditionnel

C'est le code pour faire le modèle de régression linéaire mixte, qui est la partie principale.

class ConditionalMixtureModel(object):

    def __init__(self, n_models=2):
        #Spécifiez le nombre de modèles
        self.n_models = n_models

    #Méthode pour estimer le plus probable
    def fit(self, X, t, n_iter=100):
        #Initialisation des paramètres à estimer (coefficient de poids, coefficient de mélange, variance)
        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

        #Répétez l'étape EM un nombre spécifié de fois
        for i in xrange(n_iter):
            #E taux de charge de pas E[z_nk]Est calculé pour chaque donnée et chaque modèle
            resp = self._expectation(X, t)

            #Maximiser les paramètres des étapes M
            self._maximization(X, t, resp)

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

    #Gamma de taux de charge d'étape E_nk=E[z_nk]Calculer
    def _expectation(self, X, t):
        #Formule PRML(14.37)
        responsibility = self.weight * self._gauss(X.dot(self.coef), t[:, None])
        responsibility /= np.sum(responsibility, axis=1, keepdims=True)
        return responsibility

    #Maximiser les paramètres des étapes M
    def _maximization(self, X, t, resp):
        #Formule PRML(14.38)Calcul du coefficient de mélange
        self.weight = np.mean(resp, axis=0)

        for m in xrange(self.n_models):
            #Formule PRML(14.42)Calculez le coefficient pour chaque modèle
            self.coef[:, m] = np.linalg.inv((X.T * resp[:, m]).dot(X)).dot(X.T * resp[:, m]).dot(t)

        #Formule PRML(14.44)Calcul de la variance
        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)

Données d'entraînement

y=\|x\|Préparez les données d'entraînement selon.x\ge0Quandx\le0Il est possible d'ajuster un modèle de régression linéaire dans chaque plage de, mais il est préférable de préparer deux modèles pour s'adapter aux deux plages.

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

Fonction principale

def main():
    #Création de données d'entraînement
    x_train, y_train = create_toy_data()

    #Définition du vecteur de caractéristiques (ordre 1)
    feature = PolynomialFeatures(degree=1)

    #Conversion en matrice de planification
    X_train = feature.transform(x_train)

    #Préparation du modèle mixte conditionnel (2 modèles)
    model = ConditionalMixtureModel(n_models=2)

    #Effectuer l'estimation du paramètre le plus probable
    model.fit(X_train, y_train)

    #Résultats illustrés
    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()

Code entier

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()

résultat

Fonction non linéaire en combinant plusieurs modèles linéairesy=\|x\|L'ajustement aux données selon est terminé.

À la fin

Lorsque le modèle mixte conditionnel implémenté cette fois est utilisé pour la prédiction, la ligne s'étend même là où il n'y a pas de points de données, comme indiqué dans le résultat ci-dessus. Le modèle expert mixte, qui est un développement ultérieur du modèle mixte conditionnel, semble résoudre ce problème en faisant du coefficient de mélange $ \ pi $ une fonction de l'entrée $ x $.

PRML Chapitre 14 Implémentation Python de modèle mixte conditionnel