PRML Chapter 7 Implémentation de Python Vector Machine associée pour les problèmes de régression
Cette fois, j'ai implémenté une machine vectorielle associée. Ce n'est pas si célèbre par rapport à la machine à vecteurs de support, mais il semble avoir des avantages tels que la valeur de sortie étant une probabilité. Je pense que beaucoup de gens l'utiliseront s'il est dans scikit-learn (une bibliothèque qui implémente des méthodes d'apprentissage automatique), mais je ne l'ai pas implémenté car Microsoft a un brevet.
Machine vectorielle associée
Le chapitre 7 est "Machines du noyau avec des solutions éparses", mais la discussion sur les machines vectorielles associées est généralement valable sans utiliser la méthode du noyau. Paramètres tels que la régression linéaire
p({\bf w}|{\bf \alpha}) = \prod_i\mathcal{N}(w_i|0,\alpha_i^{-1})
** Introduisez des super paramètres pour chaque paramètre **. Ensuite, le super paramètre $ {\ bf \ alpha} $ est estimé de la même manière que Chapter 3 Evidence Approximation. Le $ \ alpha $ estimé a de nombreux composants qui deviennent infinis. Cela signifie que la variance du paramètre $ w $ $ \ alpha ^ {-1} $ sera proche de 0, et le paramètre $ w $ sera une distribution raide avec la moyenne de la distribution antérieure initialement définie de 0. La valeur du paramètre $ w $ devient clairsemée et les fonctionnalités à haute pertinence sont automatiquement sélectionnées (détermination automatique de la pertinence).
Machine vectorielle associée pour les problèmes de régression
La fonction de vraisemblance pour le paramètre de poids $ {\ bf w} $ lorsque N données d'entraînement $ \ {x_n, t_n \} _ {n = 1} ^ N $ sont observées est
\begin{align}
p({\bf t}|{\bf\Phi},{\bf w},\beta) &= \prod_{n=1}^N p(t_n|\phi(x_n),{\bf w}, \beta)\\
&= \prod_{n=1}^N \mathcal{N}(t_n|{\bf w}^{\rm T}\phi(x_n), \beta^{-1})
\end{align}
Cependant, $ {\ bf t} = (t_1, \ dots, t_N) ^ {\ rm T} $ et $ \ phi (\ cdot) $ utilisaient les fonctions gaussiennes suivantes centrées sur les points de données d'entraînement. Vecteur de caractéristiques,
\phi(x) =
\begin{bmatrix}
\phi_1(x)\\
\vdots\\
\phi_N(x)
\end{bmatrix}
=
\begin{bmatrix}
a\exp(-b(x - x_1)^2)\\
\vdots\\
a\exp(-b(x-x_N)^2)
\end{bmatrix}、
$ {\ bf \ Phi} $ est une matrice de planification de $ N \ fois N $ dont les éléments sont $ \ Phi_ {ni} = \ phi_i (x_n) $. La distribution postérieure du paramètre de poids $ {\ bf w} $ est la suivante du théorème de Bayes, en utilisant $ p ({\ bf w} | {\ bf \ alpha}) $ comme distribution antérieure.
p({\bf w}|{\bf t},{\bf\Phi},{\bf\alpha},\beta) = \mathcal{N}({\bf w}|{\bf m},{\bf \Sigma})
Cependant, la moyenne et la covariance sont
\begin{align}
{\bf m} &= \beta{\bf\Sigma}{\bf\Phi}^{\rm T}{\bf t}\\
{\bf\Sigma} &= \left({\bf A} + \beta{\bf\Phi}^{\rm T}{\bf\Phi}\right)^{-1}
\end{align}
La matrice $ {\ bf A} $ est une matrice diagonale avec des éléments $ A_ {ii} = \ alpha_i $.
Sur la base de la discussion jusqu'à présent, l'estimation la plus probable des super paramètres $ {\ bf \ alpha}, \ beta $ est effectuée. La fonction de preuve logarithmique est $ {\ bf C} = \ beta ^ {-1} {\ bf I} + {\ bf \ Phi} {\ bf A} ^ {-1} {\ bf \ Phi} ^ {\ Comme rm T} $
\begin{align}
\ln p({\bf t}|{\bf\Phi},{\bf\alpha},\beta) &= \ln\mathcal{N}({\bf t}|{\bf 0},{\bf C})\\
&= -{1\over2}\left\{N\ln(2\pi) + \ln|{\bf C}| + {\bf t}^{\rm T}{\bf C}^{-1}{\bf t}\right\}
\end{align}
Sera. Ceci est différencié pour $ {\ bf \ alpha}, \ beta $ pour obtenir l'équation de mise à jour des super-paramètres. Comme $ \ gamma_i = 1- \ alpha_i \ Sigma_ {ii} $
\begin{align}
\alpha_i^{new} &= {\gamma_i\over m_i^2}\\
\beta^{new} &= {N - \sum_i\gamma_i\over||{\bf t} - {\bf\Phi}{\bf m}||^2}
\end{align}
En utilisant les nouveaux super paramètres $ {\ bf \ alpha} ^ {new}, \ beta ^ {new} $ obtenus de cette manière, calculez à nouveau la distribution postérieure du paramètre de poids $ {\ bf w} $, puis super Répétez la mise à jour des paramètres.
Vecteur connexe
Cette fois, puisque la fonction Gauss centrée sur les points de données d'apprentissage est utilisée pour le vecteur de caractéristiques, la valeur du paramètre $ w $ peut être considérée comme une valeur indiquant dans quelle mesure les données d'apprentissage contribuent à la prédiction. ** Lorsque la détermination automatique de pertinence est utilisée par la machine vectorielle associée, de nombreux paramètres $ w $ deviennent 0, et les données d'apprentissage correspondant à ceux qui ne le sont pas sont appelées vecteur associé **. Lors du calcul de la distribution de prédiction, n'utiliser que les caractéristiques correspondant aux vecteurs associés ne devrait pas changer beaucoup le résultat de la prédiction.
la mise en oeuvre
Bibliothèque
Encore une fois, la partie algorithme est codée uniquement avec numpy.
import matplotlib.pyplot as plt
import numpy as np
Régression vectorielle associée
#Une classe qui effectue une régression vectorielle associée
class RelevanceVectorRegression(object):
#Initialisation des super paramètres
def __init__(self, alpha=1., beta=1.):
self.alpha = alpha
self.beta = beta
#Calcul des vecteurs de caractéristiques à l'aide du noyau gaussien
def _kernel(self, x, y):
return np.exp(-10 * (x - y) ** 2)
#Super paramètre alpha,estimation bêta
def fit(self, x, t, iter_max=1000):
self.x = x
self.t = t
N = len(x)
#Matrice de planification
Phi = self._kernel(*np.meshgrid(x, x))
self.alphas = np.zeros(N) + self.alpha
for _ in xrange(iter_max):
params = np.hstack([self.alphas, self.beta])
#Equation PRML de covariance pour la distribution postérieure du paramètre de poids w(7.83)
self.precision = np.diag(self.alphas) + self.beta * Phi.T.dot(Phi)
self.covariance = np.linalg.inv(self.precision)
#Équation PRML moyenne pour la distribution postérieure du paramètre de poids w(7.82)
self.mean = self.beta * self.covariance.dot(Phi.T).dot(t)
#Expression PRML de validité des paramètres(7.89)
gamma = 1 - self.alphas * np.diag(self.covariance)
#Expression PRML de mise à jour des super paramètres(7.87)
self.alphas = gamma / np.square(self.mean)
#10 pour que 0% ne se produise pas^10 alpha sur 10^Définir sur 10
self.alphas = np.clip(self.alphas, 0, 1e10)
#Expression PRML de mise à jour des super paramètres(7.88)
self.beta = (N - np.sum(gamma)) / np.sum((t - Phi.dot(self.mean)) ** 2)
#Si la quantité de mise à jour des paramètres est faible, elle se termine
if np.allclose(params, np.hstack([self.alphas, self.beta])):
break
else:
#Renvoie l'instruction suivante si elle ne se termine pas même après la mise à jour du nombre de fois spécifié
print "paramters may not have converged"
#Calculer la distribution de prédiction postérieure pour l'entrée x
def predict_dist(self, x):
K = self._kernel(*np.meshgrid(x, self.x, indexing='ij'))
#Formule PRML moyenne pour la distribution de prédiction postérieure(7.90)
mean = K.dot(self.mean)
#Dispersion PRML formule de distribution de prédiction postérieure(7.91)
var = 1 / self.beta + np.sum(K.dot(self.covariance) * K, axis=1)
#Renvoie la moyenne et l'écart type de la distribution de prédiction postérieure
return mean, np.sqrt(var)
Code entier
relevance_vector_regression.py
import matplotlib.pyplot as plt
import numpy as np
class RelevanceVectorRegression(object):
def __init__(self, alpha=1., beta=1.):
self.alpha = alpha
self.beta = beta
def _kernel(self, x, y):
return np.exp(-10 * (x - y) ** 2)
def fit(self, x, t, iter_max=1000):
self.x = x
self.t = t
N = len(x)
Phi = self._kernel(*np.meshgrid(x, x))
self.alphas = np.zeros(N) + self.alpha
for _ in xrange(iter_max):
params = np.hstack([self.alphas, self.beta])
self.precision = np.diag(self.alphas) + self.beta * Phi.T.dot(Phi)
self.covariance = np.linalg.inv(self.precision)
self.mean = self.beta * self.covariance.dot(Phi.T).dot(t)
gamma = 1 - self.alphas * np.diag(self.covariance)
self.alphas = gamma / np.square(self.mean)
self.alphas = np.clip(self.alphas, 0, 1e10)
self.beta = (N - np.sum(gamma)) / np.sum((t - Phi.dot(self.mean)) ** 2)
if np.allclose(params, np.hstack([self.alphas, self.beta])):
break
else:
print "paramters may not have converged"
def predict_dist(self, x):
K = self._kernel(*np.meshgrid(x, self.x, indexing='ij'))
mean = K.dot(self.mean)
var = 1 / self.beta + np.sum(K.dot(self.covariance) * K, axis=1)
return mean, np.sqrt(var)
def create_toy_data(func, low=0., high=1., n=10, std=0.1):
x = np.random.uniform(low, high, n)
t = func(x) + np.random.normal(scale=std, size=n)
return x, t
def main():
def func(x):
return np.sin(2 * np.pi * x)
x, t = create_toy_data(func, n=10)
plt.scatter(x, t, color="blue", alpha=0.5, label="observation")
regression = RelevanceVectorRegression()
regression.fit(x, t)
relevance_vector = np.abs(regression.mean) > 0.1
x_test = np.linspace(0, 1, 100)
plt.scatter(x[relevance_vector], t[relevance_vector], color="green", s=100, marker="D", label="relevance vector")
plt.plot(x_test, func(x_test), color="blue", label="sin($2\pi x$)")
y, y_std = regression.predict_dist(x_test)
plt.plot(x_test, y, color="red", label="predict_mean")
plt.fill_between(x_test, y - y_std, y + y_std, color="pink", alpha=0.5, label="predict_std")
plt.legend()
plt.show()
if __name__ == '__main__':
main()
résultat
En raison du retrait des points bleu et vert en tant que données d'apprentissage avec la machine vectorielle associée, le résultat est comme indiqué dans la figure ci-dessous (reproduction de PRML Figure 7.9). Les points verts sont représentés comme des vecteurs associés.
À la fin
Il semble qu'il puisse faire des prédictions plus rapides tout en ayant les mêmes performances de généralisation que la machine à vecteurs de support, je voudrais donc profiter de cette opportunité pour utiliser non seulement la machine à vecteurs de support, mais également les machines vectorielles associées. C'est aussi un bel avantage de pouvoir traiter la sortie de manière bayésienne et d'évaluer la distribution. Cependant, il est dommage que ce soit contre l'intuition humaine que la variance prédictive devienne faible là où il n'y a pas de points de données d'entraînement. Cette fois, nous avons résolu le problème de régression en utilisant la machine vectorielle associée, mais nous pouvons bien sûr l'étendre au problème de classification en utilisant la fonction sigmoïde logistique. Si j'ai une chance, je voudrais également mettre en œuvre cela.
Recommended Posts