"Régression linéaire" et "Version probabiliste de la régression linéaire" en Python "Régression linéaire de Bayes"

motivation

--Je veux écrire un algorithme d'apprentissage automatique à partir de zéro en Python ――Il existe une Page de révision technique qui a une bonne réputation lorsque vous recherchez sur Google "l'apprentissage automatique", et elle sert également d'étude. Bien qu'il s'agisse d'une ancienne page, elle est recommandée à ceux qui ne l'ont pas lue car l'explication était facile à comprendre.

Flux de la recherche d'une régression linéaire

Fonction de base

――Dans "régression linéaire", définissez d'abord ce que l'on appelle une "fonction de base"

** Base polygonale **

\phi_i(x)=x^i\hspace{1.0em}(i=0,\cdots,M-1)

** Base gaussienne **

\phi_i(x)=\exp\left\{-\frac{(x-c_i)^2}{2s^2}\right\}\hspace{1.0em}(i=0,\cdots,M-1)

Comment trouver une régression linéaire

f(x)= \sum^{M-1}_{i=0}w_i\phi_i(x) = w^T\phi(x) \hspace{1.0em}(Équation 1)

--F (x) peut être obtenu en déterminant le paramètre $ w_i $, il devrait donc y avoir un moyen de déterminer $ w_i $ de manière appropriée.

E(w)=\frac{1}{2}\sum^N_{n=1}(f(x_n)-t_n)^2 
=\frac{1}{2}\sum^N_{n=1}\left(\sum^{M-1}_{i=0}w_{i}x^{i}_{n}-t_n\right)^2 \hspace{1.0em}(Équation 2)

--C'est compliqué à première vue, mais en se rappelant que $ x_n $ et $ t_n $ sont des constantes (points déjà donnés), cela peut être considéré comme une fonction de $ w_i $, juste une expression quadratique de $ w_i $. Est --Par conséquent, vous pouvez trouver $ w_i $ qui minimise la fonction d'erreur $ E (w) $, et vous pouvez également décider de la fonction f (x) que vous voulez trouver.

Exemple de données

Étant donné un graphique comme celui ci-dessous, je veux trouver sa régression linéaire

%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt

X = np.array([0.02, 0.12, 0.19, 0.27, 0.42, 0.51, 0.64, 0.84, 0.88, 0.99])
t = np.array([0.05, 0.87, 0.94, 0.92, 0.54, -0.11, -0.78, -0.79, -0.89, -0.04])

plt.xlim([0,1.0])
plt.ylim([-1.5,1.5])
plt.plot(X, t, 'o')
plt.show()

output_3_0.png

Trouvez la relation entre X et t par régression linéaire

Trouver une régression linéaire à l'aide d'une base polynomiale

Base de polygone

Défini par la formule suivante (republication)

\phi_i(x)=x^i\hspace{1.0em}(i=0,\cdots,M-1)

--La base polynomiale a l'avantage que la solution peut être obtenue dans le format familier de "polypole", mais elle a la caractéristique qu'elle influence fortement l'estimation en un point éloigné.

--Cette fois, nous allons essayer de résoudre la fonction caractéristique φ en utilisant la base polynomiale de la fonction cubique. Par conséquent, l'équation de régression $ f (x) $ à obtenir est la suivante.

f(x)=w_{1}+w_{2}x+w_{3}x^2+w_{4}x^3
\frac{\partial E(w)}{\partial w_m}=\sum^N_{n=1}\phi_m(x_n)\left(\sum^M_{j=1}w_j\phi_j(x_n)-t_n\right)=0 \hspace{1.0em}(m=1,\cdots,M) \hspace{1.0em}(Équation 3)

--Si ces équations simultanées sont arrangées à l'aide d'une matrice, la solution $ w $ peut être obtenue par l'équation suivante.

w = (\phi^T\phi)^{-1}\phi^{T}t \hspace{1.0em}(Équation 4)

Cependant, Φ est une matrice définie comme suit.

	\phi=\left(
	\begin{array}{ccc} 
    	\phi_{1}(x_1) & \phi_{2}(x_1) & \cdots & \phi_{M}(x_1)  \\ 
    	\vdots & \vdots& & \vdots      \\ 
        \phi_{1}(x_N) & \phi_{2}(x_N) & \cdots & \phi_{M}(x_N)\\ 
	\end{array}  
	\right) \hspace{1.0em}(Équation 5)

Trouver une régression linéaire à base polynomiale (fonction du troisième ordre) (implémentation du contenu ci-dessus)

Écrivez le contenu ci-dessus tel qu'il est dans le code et trouvez la régression linéaire avec la base polynomiale (fonction du troisième ordre) (* Les termes de normalisation ne sont pas pris en compte ici)

#Lorsqu'une base polynomiale est utilisée pour la fonction caractéristique φ
def phi(x):
    return [1, x, x**2, x**3]

PHI = np.array([phi(x) for x in X])
#Calcul de la partie ci-dessus (équation 4). np.dot:Trouvez le produit intérieur np.linalg.inv:Trouvez l'inverse
w = np.dot(np.linalg.inv(np.dot(PHI.T, PHI)), np.dot(PHI.T, t))
#Np pour trouver la solution d'équations simultanées.linalg.Vous pouvez également utiliser la fonction de résolution. C'est plus rapide. np.linalg.solve(A, b): A^{-1}renvoie b
#w = np.linalg.solve(np.dot(PHI.T, PHI), np.dot(PHI.T, t))

xlist = np.arange(0, 1, 0.01) #Le point x de l'équation de régression linéaire. 0.Tracer finement en 01 unités et sortie sous forme de ligne
ylist = [np.dot(w, phi(x)) for x in xlist] #Équation de régression linéaire y

plt.plot(xlist, ylist) #Diagramme d'équation de régression linéaire
plt.plot(X, t, 'o')
plt.show()

output_8_0.png

#Sortie w obtenue par le calcul ci-dessus
w
array([ -0.14051447,  11.51076413, -33.6161329 ,  22.33919877])

En d'autres termes, à la suite de la régression linéaire, l'équation de régression $ f (x) = -0,14 + 11,5x-33,6x ^ 2 + 22,3x ^ 3 $ a été obtenue.

Trouver une régression linéaire en utilisant une base gaussienne

Base gaussienne

La base gaussienne est une fonction en forme de cloche (pas une distribution) définie par l'équation suivante

\phi_i(x)=\exp\left\{-\frac{(x-c_i)^2}{2s^2}\right\}\hspace{1.0em}(i=0,\cdots,M-1)\hspace{1.0em}(Équation 6)

Variable de base gaussienne s

Variable de base gaussienne c_i

Trouver une régression linéaire sur une base gaussienne (implémentation de ce qui précède)

Tout d'abord, sortons simplement un diagramme dans lequel la fonction caractéristique φ utilisée dans la régression linéaire ci-dessus est réécrite de "base polypoly" à "base de Gauss".

%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt

X = np.array([0.02, 0.12, 0.19, 0.27, 0.42, 0.51, 0.64, 0.84, 0.88, 0.99])
t = np.array([0.05, 0.87, 0.94, 0.92, 0.54, -0.11, -0.78, -0.79, -0.89, -0.04])

#Lorsqu'une base gaussienne est utilisée pour la fonction caractéristique φ
#phi renvoie un vecteur à 12 dimensions en ajoutant 1 pour représenter le terme constant en plus des 11 fonctions de base gaussiennes.
#Base gaussienne (0 à 1).0 à 0.11 points par incréments de 1)+Renvoie un total de 12 points (dimensions), un point du terme constant
def phi(x): 
    s = 0.1 #Base gaussienne "largeur"
    #Écrire l'équation 6
    return np.append(1, np.exp(-(x - np.arange(0, 1 + s, s)) ** 2 / (2 * s * s)))

PHI = np.array([phi(x) for x in X])
w = np.linalg.solve(np.dot(PHI.T, PHI), np.dot(PHI.T, t))

xlist = np.arange(0, 1, 0.01)
ylist = [np.dot(w, phi(x)) for x in xlist]

plt.plot(xlist, ylist)
plt.plot(X, t, 'o')
plt.show()

output_14_0.png

La base gaussienne à 12 dimensions est complètement surapprise car elle a plus de liberté et d'expressivité que la base polymorphe à 4 dimensions.

Régression linéaire bayésienne qui probabilise la régression linéaire

――La régression linéaire stochastique est basée sur l'idée que "du bruit aléatoire est ajouté aux données observées".

p(t|w,x) = N(t|\mu,\beta^{-1})=\frac{1}{Z}\exp\left\{-\frac{1}{2}\beta(t-\mu)^2\right\} \hspace{1.0em}(Équation 7)\\
cependant,\mu=f(x)=\sum^M_{i=1}w_i\phi_i(x), Z=\sqrt{2\pi\beta^{-1}}

--Lorsque le point $ (x_1, t_1) $ est affecté à l'équation 7 ci-dessus, la fonction $ p (t_1 | w, x_1) $ de $ w $ obtenue par celle-ci est appelée "fonction de probabilité".

――La probabilité peut être obtenue en multipliant la vraisemblance à chaque point, et $ w $ qui la maximise peut être obtenue.

\ln{p(T|w,X)}=-\frac{1}{2}\beta\sum^N_{n=1}\left(t_n-\sum^{M}_{i=1}w_{i}\phi_i(x_n)\right)^2+C \hspace{1.0em}(Équation 8)

Trouver la distribution a posteriori p (t | w, x) à partir de la distribution a priori p (w)

-Distribution ex postp(t|w,x)Surtout que x est une valeur d'entrée, traitez-la comme une constantep(t|w,x)Aussip(t|w)Noté comme --Mettre à jour la distribution de $ w $ vers la distribution postérieure $ p (w | t) $ en utilisant $ p (w) $ et les points de données (x, t). Cette distribution a posteriori est la «confiance de la réponse» que vous voulez trouver. --La formule de mise à jour de la distribution postérieure peut être obtenue à partir de la formule bayésienne et du théorème de multiplication / addition.

p(w|t)=p(t|w)p(w)/p(t) \\
p(t)=\int p(t,w)dw=\int p(t|w)p(w)dw
p(w|t)=\frac{1}{Z'}\exp\left\{-\frac{1}{2}(w-m)^T\sum^{-1}(w-m)\right\} \hspace{1.0em}(Équation 9)\\
pourtant\\
m=\beta t\sum\phi(x)\\
\sum^{-1}=\alpha I+\beta\phi(x)\phi(x)^T
N(w|\mu,\sum)=\frac{1}{Z}\exp\left\{(-(w-\mu)^T\sum^{-1}(w-\mu)\right\}
p(w|t)=N(w|m,\sum)
p(w|t,x)=N(w|m_N,\sum_N) \hspace{1.0em}(Équation 10)\\
pourtant\\
m_N=\beta \sum_N\phi^{T}t\\
\sum^{-1}_N=\alpha I+\beta\phi^T\phi

-En conséquence, d'abord les points de données(x_1,t_1)Pour sa distribution postérieurep(w|t_1,x_1)Et considérez-le à nouveau comme une distribution antérieure, un autre point de données(x_2,t_2)La distribution postérieure suivante utilisantp(w|t_1,t_2,x_1,x_2)Peut être obtenu en répétant l'opération de recherche ...p(w|t_1,…,t_N,x_1,…,x_N)Aura la même distribution que

Régularisation de la régression linéaire bayésienne

Implémentation de la régression linéaire bayésienne

Comparons la régression linéaire bayésienne (ligne bleue) et la régression linéaire normale (ligne verte) sur le graphique.

%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt

X = np.array([0.02, 0.12, 0.19, 0.27, 0.42, 0.51, 0.64, 0.84, 0.88, 0.99])
t = np.array([0.05, 0.87, 0.94, 0.92, 0.54, -0.11, -0.78, -0.79, -0.89, -0.04])

#Lorsqu'une base gaussienne est utilisée pour la fonction caractéristique φ
def phi(x): 
    s = 0.1
    return np.append(1, np.exp(-(x - np.arange(0, 1 + s, s)) ** 2 / (2 * s * s)))

PHI = np.array([phi(x) for x in X])

#(Idem que ci-dessus) Résolution d'une régression linéaire en utilisant une base gaussienne → Superapprentissage
w = np.linalg.solve(np.dot(PHI.T, PHI), np.dot(PHI.T, t))

#Résolution de la régression linéaire bayésienne en utilisant une base gaussienne → Éviter le surapprentissage
#Une transcription de l'équation 10 ci-dessus
alpha = 0.1 #Stockage temporaire
beta = 9.0  #Stockage temporaire
Sigma_N = np.linalg.inv(alpha * np.identity(PHI.shape[1]) + beta * np.dot(PHI.T, PHI)) #np.identity(PHI.shape[1])Est une matrice unitaire à 12 dimensions spécifiée dans la base gaussienne
mu_N = beta * np.dot(Sigma_N, np.dot(PHI.T, t))

xlist = np.arange(0, 1, 0.01) 
plt.plot(xlist, [np.dot(w, phi(x)) for x in xlist], 'g')     #Tracez la solution de la régression linéaire ordinaire
plt.plot(xlist, [np.dot(mu_N, phi(x)) for x in xlist], 'b')  #Tracer la solution de la régression linéaire bayésienne
plt.plot(X, t, 'o') #Tracer des points d'échantillonnage
plt.show()

output_20_0.png

Quelle était la covariance Σ_N?

--Covariance $ \ sum_N $ représentait la "confiance de probabilité" du point de données --La "dispersion" d'une distribution unidimensionnelle représente la "dispersion" des données

#Matrice de covariance Sigma_Affichage de N facile à lire
print "\n".join(' '.join("% .2f" % x for x in y) for y in Sigma_N)
 2.94 -2.03 -0.92 -1.13 -1.28 -1.10 -1.21 -1.14 -1.23 -1.06 -0.96 -2.08
-2.03  2.33 -0.70  1.95  0.13  1.02  0.85  0.65  0.98  0.70  0.65  1.44
-0.92 -0.70  2.52 -1.86  1.97 -0.29  0.42  0.59  0.13  0.40  0.32  0.63
-1.13  1.95 -1.86  3.02 -1.66  1.50  0.17  0.29  0.73  0.33  0.36  0.82
-1.28  0.13  1.97 -1.66  2.82 -1.11  1.39  0.22  0.55  0.49  0.40  0.92
-1.10  1.02 -0.29  1.50 -1.11  2.39 -1.35  1.72 -0.29  0.53  0.46  0.69
-1.21  0.85  0.42  0.17  1.39 -1.35  2.94 -2.06  2.39 -0.02  0.25  1.01
-1.14  0.65  0.59  0.29  0.22  1.72 -2.06  4.05 -2.72  1.43  0.37  0.67
-1.23  0.98  0.13  0.73  0.55 -0.29  2.39 -2.72  3.96 -1.41  1.23  0.59
-1.06  0.70  0.40  0.33  0.49  0.53 -0.02  1.43 -1.41  3.30 -2.27  2.05
-0.96  0.65  0.32  0.36  0.40  0.46  0.25  0.37  1.23 -2.27  3.14 -0.86
-2.08  1.44  0.63  0.82  0.92  0.69  1.01  0.67  0.59  2.05 -0.86  2.45

--Il s'avère que la matrice de covariance est une matrice diagonale

Visualisez la distribution prévue

-(Idée) distribution postérieurep(w|X)Est la distribution du paramètre w, doncp(y|w,x)Une "fonction" qui renvoie la distribution de y lorsque x est donné.p(y|X,x)(En fait, il n'est pas possible de substituer la distribution à la distribution, c'est donc juste une image) --La distribution de la variable cible qui reflète la distribution postérieure est appelée distribution prédite.

#Fonction de densité de probabilité de distribution normale
def normal_dist_pdf(x, mean, var): 
    return np.exp(-(x-mean) ** 2 / (2 * var)) / np.sqrt(2 * np.pi * var)

#Format secondaire( x^Calculer T A x)
def quad_form(A, x):
    return np.dot(x, np.dot(A, x))

xlist = np.arange(0, 1, 0.01)
tlist = np.arange(-1.5, 1.5, 0.01)
z = np.array([normal_dist_pdf(tlist, np.dot(mu_N, phi(x)),1 / beta + quad_form(Sigma_N, phi(x))) for x in xlist]).T

plt.contourf(xlist, tlist, z, 5, cmap=plt.cm.binary)
plt.plot(xlist, [np.dot(mu_N, phi(x)) for x in xlist], 'r')
plt.plot(X, t, 'go')
plt.show()

output_26_0.png

--La partie sombre a une valeur élevée de la fonction de densité de probabilité, c'est-à-dire la partie où la fonction estimée est susceptible d'y passer.

À propos de l'hyper paramètre α de la distribution antérieure

--Les paramètres contenus dans la distribution a priori basilienne $ p (w) = N (w | 0, \ alpha ^ -1I) $ sont également appelés hyperparamètres. ―― Plus α est grand, plus la dispersion est petite, c'est-à-dire plus la connaissance préalable que w est une valeur proche de 0 est forte. --Si vous résolvez la régression linéaire bayésienne dans cet état, la force d'estimer w à une valeur proche de 0 est forte, il est donc facile de supprimer le résultat du soi-disant surapprentissage, mais vous arrivez à la vraie solution. Peut nécessiter beaucoup de données ――Au contraire, si α est petit, la force pour maintenir w devient faible.

Ancienne histoire

--Plusieurs pages de "Commençons avec l'apprentissage automatique" http://gihyo.jp/dev/serial/01/machine-learning

Recommended Posts

"Régression linéaire" et "Version probabiliste de la régression linéaire" en Python "Régression linéaire de Bayes"
Implémenté en Python PRML Chapitre 3 Régression linéaire bayésienne
Régression linéaire en ligne en Python
Implémentation python de la classe de régression linéaire bayésienne
Régression linéaire en Python (statmodels, scikit-learn, PyMC3)
Valeurs authentiques et vecteurs propres: Algèbre linéaire en Python <7>
Régression linéaire en ligne en Python (estimation robuste)
Reproduction sur plaque de régression linéaire bayésienne (PRML §3.3)
Indépendance et base linéaires: Algèbre linéaire en Python <6>
Liste des solveurs et modélisateurs de conception linéaire (LP) disponibles en Python
Projet Euler # 1 "Multiple de 3 et 5" en Python
J'ai essayé d'implémenter la régression linéaire bayésienne par échantillonnage de Gibbs en python
Matrice unitaire et matrice inverse: Algèbre linéaire en Python <4>
Produit intérieur et vecteur: Algèbre linéaire en Python <2>
Calcul matriciel et équations linéaires: Algèbre linéaire en Python <3>
Explication de la distance d'édition et de l'implémentation en Python
Recherche linéaire en Python
Analyse de régression avec Python
J'ai comparé la vitesse des expressions régulières en Ruby, Python et Perl (version 2013)
Traitement pleine largeur et demi-largeur des données CSV en Python
Calcul de l'écart type et du coefficient de corrélation en Python
Différence entre Ruby et Python en termes de variables
[python] Calcul des mois et des années de différence de date / heure
Coursera Machine Learning Challenge en Python: ex1 (régression linéaire)
Implémentation d'estimation bayésienne de variante du modèle de sujet en python
Exemple d'obtention du nom du module et du nom de la classe en Python
Récapitulatif du traitement de la date en Python (datetime et dateutil)
Gestion des versions de Node, Ruby et Python avec anyenv
Remarques sur l'utilisation de StatsModels qui peuvent utiliser la régression linéaire et GLM en python
Expression de régression multiple en Python
Jugement d'équivalence d'objet en Python
Mise à niveau de python Anaconda
Vérifiez la version OpenSSL de python 2.6
Pile et file d'attente en Python
Analyse de régression simple avec Python
[Python] Régression linéaire avec scicit-learn
Unittest et CI en Python
Implémentation du tri rapide en Python
Installation source et installation de Python
Avoir le graphique d'équation de la fonction linéaire dessiné en Python
Ordre de référence des variables de classe et des variables d'instance dans "self. Variables de classe" en Python
Comparaison de l'utilisation des fonctions d'ordre supérieur dans Python 2 et 3
Première physique computationnelle: mécanique quantique et algèbre linéaire avec python.
[Python] Forces et faiblesses de DataFrame en termes de temps requis
[Astuces] Problèmes et solutions dans le développement de python + kivy
[Python] Version Taple du menu déroulant de la préfecture
Manipulation des pixels d'image en Python
L'histoire de Python et l'histoire de NaN
Paquets qui gèrent le MIDI avec Python midi et pretty_midi
Comptez bien le nombre de caractères thaïlandais et arabes en Python
Différence entre list () et [] en Python
Première analyse de régression simple en Python
Différence entre == et est en python
Installer SciPy et matplotlib (Python)
Afficher les photos en Python et html
Algorithme de tri et implémentation en Python
Python: Application de l'apprentissage supervisé (retour)
Diviser timedelta dans la série Python 2.7