Skip to main content

Section 3.3 Régression et maximum à posteriori

Jusqu’a présent on a considérer des problèmes de régression ou il existe une unique solution. Pour cela il faut que le matrice \(A^t A\) soit inversible dans l’équation normale. A priori si on a suffisamment d’échantillons (on y reviendra) on se retrouvera dans ce cas là. On va maintenant voir qu’en grande dimension \(d \gt \gt 1\) on peut se trouver fréquemment dans un cadre ou on perd l’unicité de la solution. On va introduire cela et montrer comment on peut modifier le problème de régression pour aider à la résolution du problème.

Subsection 3.3.1 Grande dimension et surapprentissage

Subsubsection 3.3.1.1 Cas linéaire

Les résultats centraux de la méthode des moindres carrés Lemme 3.2 Lemme 3.3 nous assure l’unicité de la solution si le \(n \gt k\text{.}\) Lorsqu’on applique ce résultat à la régression linéaire de la Proposition 3.19 on a donc unicité de la solution pour \(n > d+1\) avec \(n\) le nombre d’exemples et \(d\) la dimension.
Il peut arriver de regarder des problèmes en très grande dimension \(d \gt \gt 1\) (notamment dans des problèmes de réduction de dimension). Dans ce cas, il est possible que le nombre d’exemples soit inférieur à la dimension du problème. Dans ce cas on perd l’unicité de la solution du problème de régression. En effet la matrice \(A^t A\) est de rang maximal \(n\) alors que matrice est de dimension \(d+1\text{.}\) On a donc une infinité de solution.
Cette perte d’unicité permet au modèle d’interpolation exactement les \(n\) exemples. Cependant il reste des degrés de liberté au modèle qui lui permettent de capturer une infinité de solutions parfois très éloignées de la bonne. Lorsqu’on appliquera le modèle à de nouvelles données, on risque de faire une erreur très importante si la solution capturée est loin de la bonne.
Figure 3.39. Erreur de généralisation pour une régression linéaire en fonction de la dimension et sur-apprentissage. On peut retrouver le notebook ici 1 .
Sur la figure Figure 3.39 on voit que plus la dimension \(d\) se rapproche de \(n \) plus l’erreur sur le jeu d’entraînement diminue. Pour \(d+1 \ge n\) l’erreur devient nulle sur le jeux d’entraînement. Par contre l’erreur de généralisation sur de nouvelles données s’envole totalement.
Définition 3.40.
Lorsqu’on un modèle approche bien les données d’entraînement, mais très mal de nouvelles de données on par de surapprentissage.

Subsubsection 3.3.1.2 Cas Polynomial

Précédemment, sur la méthode de régression linéaire on a vu émerger un 1er problème de perte d’unicité de la solution et de surapprentissage lorsque le nombre de paramètres du modèles (linéaire par rapport à la dimension \(d\)) est plus grand que le nombre d’échantillons. Ce problème existe évidemment de la même façon en régression polynomiale lorsque le le nombre de paramètres ( qui croît exponentiellement avec la dimension) est plus grand la aussi que le nombre d’échantillons.
Cependant le surapprentissage peut advenir même dans un cadre ou le problème de régression a une unique solution. Pour cela on va utiliser une régression polynomiale en dimension un.
Figure 3.42. Comparaison de régression polynomial avec des polynômes de degré 1, 3, 5 et 18 avec 70 points d’entrainement (en haut) et 20 points d’entrainement (en bas). On peut retrouver le notebook ici 2 .
Sur la figure Figure 3.42 on voit lorsqu’on prend suffisamment de points avec peu de bruit que l’ensemble des modèles polynomiaux à part le polynôme de degré un marche bien. Lorsqu’on baisse le nombre de points et/ou augmente le bruit, on voit que le polynôme de degré 18 donne une reconstruction très fortement oscillante et est beaucoup moins satisfaisant que la reconstruction donné par le polynôme de degré 5. Pour autant on a encore plus de points dans le jeu d’entraînement que de degré de liberté dans le modèle. On a donc pas de perte d’unicité de la régression.
Pour autant on voit que le modèle avec le polynôme de degré 18 passe bien par les points de l’échantillon et donc l’erreur au sens des moindres carrés est bien minimisée. Le fait qu’on approche très bien l’échantillon, mais très mal la fonction est justement du surapprentissage qu’on avait introduit précédemment. Si on approche quasiment parfaitement l’échantillon, on risque d’apprendre aussi le bruit. On peut aussi exacerber les défauts du modèle. Ici on utilise un polynôme de degré très élevé qui naturellement peu oscillé fortement et c’est ce qu’il va faire pour approcher au mieux les points.
Remarque 3.43.
De façon général, utiliser un modèle trop complexe (comportant trop de paramètres ou ayant tendance à trop osciller) peut mener a du sur-apprentissage si on n’a pas assez de données.

Subsubsection 3.3.1.3 Séparation des données de données et arrêt précoce

On a vu que le surapprentissage est un phénomène ou un apprentissage très précis sur les données du jeu d’apprentissage empêche d’obtenir des bonnes propriétés de généralisation sur de nouvelles données. Afin de déterminer ce surapprentissage on procède de la façon suivante:
  1. On prend notre jeu de données \(\mathcal{X}\) et on le divise en deux jeux de données \(\mathcal{X}_{train}\) et \(\mathcal{X}_{test}\)
  2. On procède à l’apprentissage avec le jeu de données d’entraînement \(\mathcal{X}_{train}\text{.}\)
  3. On vérifie les propriétés de généralisation du modèle avec le jeu de données test \(\mathcal{X}_{train}\text{.}\)
Le phénomène d’apprentissage d’un point de vue des erreurs ressemble souvent à la courbe Figure 3.44.
Figure 3.44. Exemple de courbe d’apprentissage
En général les deux erreurs descendent au début et au bout d’un moment l’erreur sur le test décroche et commence a remonter alors que celle associée au jeu d’entraînement continue de diminuer.
On propose l’exemple suivant: on propose une régression polynomiale avec un polynôme de degré \(30\) et un échantillon de taille \(n=30\) sur une fonction sinus. On va utiliser un processus itératif (type gradient) et stopper a plusieurs moments. On obtient le résultat de la figure:
Figure 3.45. Arrêt précoce pour une régression polynomiale avec un polynôme de degré 30.
On remarque, que lorsque le nombre d’itérations est importante (donc on a bien convergé sur le jeu d’apprentissage), on retrouve le phénomène de surapprentissage. On voit que finalement en arrêtant après 12 itérations (seconde image) on obtient un modèle très précis. On obtient la méthode d’arrêt précoce qui consiste à:
  1. On fait quelques itérations d’une méthode d’apprentissage itératif (gradient stochastique, etc).
  2. On calcule l’erreur sur les jeux de données "test" et "train".
  3. Si l’erreur sur le jeu de données test recommence a remonter on arrête sinon on repart à l’étape une.
Cette méthode est très simple, mais pose des problèmes. Le principal est la détection du moment où on doit arrêter. En effet la courbe d’erreur n’est pas forcement monotone, on peut facilement avoir des oscillations notamment avec les méthodes stochastiques. On ne peut donc pas arrêter dès que l’erreur sur le jeu test remonte. On doit donc souvent arrêter lorsque l’erreur remonte de façon continue pendant quelques itérations et donc garder un peu d’effet du surapprentissage.

Subsection 3.3.2 Régression MAP et régularisation

On a vu précédemment qu’en grande dimension (lorsque le nombre d’exemples n’est pas suffisant par rapport à la dimension on souffre de surapprentissage. Dans ce cas le modèle approche très bien les exemples, mais généralise très mal. On a vu à travers deux exemples distincts que le modèle avait trop de liberté pour approcher le jeu de données. On va voir ici on comment ajouter des a priori sur la forme des fonctions cibles dans le problème afin de contraindre plus le problème et de limiter le surapprentissage. Pour cela on va se placer dans le cadre bayésien.
Les statistiques Bayesiennes sont une branche des statistiques fondée sur des probabilités exprimant un degré de croyance en l’événement que l’ont souhaite regarder. Cette croyance peut être basée sur des connaissances a priori modélisées à travers les probabilités conditionnelles. Le théorème de Bayes étant central pour les probabilités conditionnelles il a donné le nom à cette branche. L’apprentissage bayésien et une variante de l’apprentissage ou l’on peut aussi se donner une connaissance a priori sur le modèle à travers une loi de probabilité sur les paramètres. On voit ici donc la différence avec les méthodes d’apprentissage précédentes. Avant on cherchait un jeu de paramètre optimal en un certain sens.
On va commencer par introduire une série de définition et de notions qui nous seront utiles par la suite. On se donne toujours des données observées notées \(x\text{,}\) les données cibles \(y\) et les paramètres du modèle \(\theta\text{.}\)

Définition 3.46.

On appelle la distribution d’échantillonnage (ou aussi la fonction de vraisemblance):
\begin{equation*} p(\mathcal{Y} \mid \mathcal{X}, \theta) \end{equation*}
si on a des données observées et cibles.

Définition 3.47.

On appelle la distribution à priori:
\begin{equation*} p(\theta \mid \alpha) \end{equation*}
avec \(\alpha\) les hyperparamètres.
Le coeur de l’apprentissage Bayésien est dans l’existence de cette distribution a priori \(p(\theta \mid \alpha)\text{.}\) Dans le cas particulier on cette distribution est une loi uniforme on se retrouve dans le cadre de l’apprentissage standard.

Définition 3.48.

On appelle la vraisemblance marginale:
\begin{equation*} p(\mathcal{Y} \mid \mathcal{X}, \alpha) = \int_{\theta} p(\mathcal{Y} \mid \mathcal{X},\theta)p(\theta \mid \alpha) d\theta \end{equation*}
si on a des données observées et cibles.

Définition 3.49.

On appelle la distribution a posteriori:
\begin{equation*} p(\theta \mid \mathcal{X},\mathcal{Y}, \alpha) = \frac{p(\mathcal{Y}\mid \mathcal{X},\theta,\alpha)p(\theta \mid \alpha)}{p(\mathcal{Y}\mid, \mathcal{X}, \alpha)} \end{equation*}
si on a des données observées et cibles.
jusqu’à présent on a défini les distributions liées à la phase d’apprentissage puisqu’elle font intervenir les jeux de données. On va maintenant introduire les notions permettant d’utiliser les modèles bayésiens pour la prédiction.

Définition 3.50.

On appelle la distribution a posteriori prédictive:
\begin{equation*} p(y \mid \mathcal{Y}, \mathcal{X}, \alpha) = \int p(y\mid x, \theta)p(\theta \mid \mathcal{Y}, \mathcal{X},\alpha) d \theta \end{equation*}
si on a des données observées et cibles.
Se pose la question de comment construire le modèle. L’idée est toujours d’appliquer une sorte de maximum de vraisemblance, mais en tenant compte de l’a priori. On va donc introduire une généralisation de l’estimateur du maximum de vraisemblance dans le cas bayésien. On parle d’estimateur du maximum a posteriori.

Définition 3.51.

Soit une distribution d’échantillonnage \(f(x\mid \theta)\text{.}\) Soit \((x_1,...x_n)\) un échantillon de vraisemblance
\begin{equation*} \mathcal{L}(x_1,...x_n \mid \theta) \end{equation*}
On appelle estimateur du maximum a posteriori la solution de
\begin{equation*} \hat{\theta}_{map}=\operatorname{argmax}_{\theta} \mathcal{L}(\theta \mid x_1,...x_n, \alpha)= \operatorname{argmax}_{\theta} \mathcal{L}( x_1,...x_n \mid \theta, \alpha)p(\theta \mid \alpha) \end{equation*}
avec \(p(\theta \mid \alpha)\) distribution a priori.
Dans la définition de cet estimateur on néglige souvent \(\alpha\) qui est un paramètre fixe. On voit que dans le cadre Bayésien le maximum de vraisemblance revient à maximiser la distribution d’échantillonnage alors que l’estimateur MAP maximise la distribution a posteriori. Les deux approches sont identiques lorsqu’on considère que \(p(\theta\mid \alpha)\) est une loi uniforme.

Subsubsection 3.3.2.1 Régression Ridge et a priori Gaussien

On a vu jusqu’à présent, qu’en grande dimension la trop grande liberté laissée au modèle mène au surapprentissage et problème de généralisation. Le cadre bayésien du MAP permet de se donner un a priori sur le modèle. On va ici se donner un premier a priori sur le modèle et montrer qu’il permet de régler le problème d’unicité de la solution et dans une certaine mesure du surapprentissage.
On considère ici un modèle linéaire du type (3.7). Un a priori raisonnable (notamment si on normalise les données) et lorsqu’on est en grande dimension c’est que poids \(\theta\) soit petits. En effet si \(y\) et les \(x_i\) dans le modèle sont en \(O(1)\) on voit que pour \(d\) grand une solution naturelle (pas la seule) est d’avoir des petits poids \(\theta>\text{.}\) Cela nous permet de choisir un premier a priori.
Définition 3.52.
On appelle a priori Gaussien le distribution à priori:
\begin{equation} p(\theta \mid \alpha) = \mathcal{N}(0,\frac{1}{\alpha})\tag{3.29} \end{equation}
avec \(\alpha\) un hyper paramètre de la loi.
Preuve.
Afin de déterminer les paramètres \(\theta=(\omega,b)\) On va donc calculer l’estimateur du maximum a posteriori pour le modèle linéaire (3.8)- (3.9) avec a priori gaussien (3.29) On cherche donc a résoudre:
\begin{equation*} \operatorname{argmax}_{\theta} \mathcal{L}( x_1,...x_n \mid \theta, \alpha)p(\theta \mid \alpha) \end{equation*}
Cela revient donc à maximiser
\begin{equation*} \mathcal{L}(y_1,....y_n,\theta\alpha)p(\theta \mid \alpha)= \left(\Pi_{i=1}^n \mathcal{N}(y_i \mid (\omega,x_i) + b,\sigma^2)\right)\mathcal{N}(0,\frac{1}{\alpha}) \end{equation*}
Comme dans le cas du maximum de vraisemblance on passe au log ce qui revient à
\begin{equation*} \operatorname{log}\left(\mathcal{L}(y_1,....y_n,\theta,\alpha)p(\theta\mid \alpha)\right)= - \frac{1}{2\sigma^2}\sum_{i=1}^n ((\omega,x_i) + b -y_i)^2 + \frac{N}{2} \operatorname{ln} (\sigma^{-2}) - N\operatorname{ln} (2\pi) - \frac{\alpha}{2}\left(\omega,\omega) + b^2\right) + \frac{N}{2} \operatorname{ln} \alpha \end{equation*}
Résoudre ce problème revient donc au final a minimiser la fonctionnelle
\begin{equation} \mathcal{J}(\theta)= \frac{1}{2}\sum_{i=1}^n ((\omega,x_i) + b -y_i)^2 + \frac{\sigma^2\alpha}{2}\left(\omega,\omega) + b^2\right)\tag{3.31} \end{equation}
En prenant \(\lambda=\sigma^2\alpha\text{,}\) on conclut.
Preuve.
On ajoute un terme strictement convexe aux problèmes des moindres carrés. Par conséquent la preuve de convexité et coercivité sont similaire à celle de Lemme 3.2. On va donc maintenant calculer le gradient
\begin{equation*} \nabla_{\theta} \mathcal{J}(\theta)=0 \end{equation*}
On calcul le gradient
\begin{align*} \mathcal{J}(\theta+h) \amp =\langle A(\theta+h)-b ,A(\theta+h)-b\rangle + \lambda\langle \theta+h ,\theta+h\rangle\\ \amp=\mathcal{J}(\theta)+2\langle A \theta -b,A h\rangle + \lambda 2 \langle \theta,h\rangle + o(h)\\ \amp =\mathcal{J}(\theta)+\langle 2A^t(A\theta-b)+2\theta,h\rangle+o(h) \end{align*}
Donc
\begin{equation*} \nabla \mathcal{J}_r(x)h=\lim\limits_{\varepsilon\rightarrow 0}\frac{\mathcal{J}_r(\theta+\varepsilon h)-\mathcal{J}(\theta)}{\varepsilon}=\langle 2A^t(A \theta-b),h\rangle \end{equation*}
et ainsi
\begin{equation*} \nabla_{\theta} \mathcal{J}(\theta) =0 \Longleftrightarrow 2A^t(A\theta+\lambda\theta - b)=0 \end{equation*}
On passe à la limite
\begin{equation*} \nabla_{\theta} \mathcal{J}(\theta)h=\lim\limits_{\varepsilon\rightarrow 0}\ frac{\mathcal{J}(\theta+\varepsilon h)-\mathcal{J}(\theta)}{\varepsilon}= \langle 2A^T(A \theta-b) + 2\lambda I_d \theta,h\rangle \end{equation*}
On obtient donc une équation normale modifiée:
\begin{equation} (A^tA+\lambda I_d)\theta =A^t b\tag{3.32} \end{equation}
La matrice est symétrique. On va vérifier qu’elle est définie positive. Soit \(\theta \in \mathbb{R}^{d+1,1}\) privé du vecteur nul. On calcul
\begin{equation*} \langle (A^tA+\lambda I_d)\theta ,\theta\rangle=\langle A^tA \theta\rangle+\lambda\langle\theta,\theta\rangle =\Norm A\theta\Norm_2^2 +\lambda \Norm \theta\Norm_2^2>0 \end{equation*}
La matrice est donc inverible ce qu conclu la preuve.
On peut voir ce problème de minimisation comme le problème de régression linéaire auquel on ajouté une contrainte faible sur la norme des poids. Cela revient à chercher de poids de pas trop grande taille. On voit qu’il s’agit de la traduction en termes d’optimisation de notre a priori Gaussien qui distribue les poids autour d’une moyenne nulle. Cela permet de retrouver l’unicité du problème. De plus si \(\lambda\) est pas trop grand on ne perturbe pas trop le problème initial. Évidemment on peut espérer améliorer les résultats que si la solution est peu éloignée de notre a priori.
Figure 3.55. On effectue une régression Ridge en dimension \(d=50\) avec des poids \(\theta\) générés par une loi \(\mathcal{N}(0,\frac{}{\alpha})\text{,}\) avec 100 données d’entrainement. On trace les erreurs en fonction du coefficient de régularisation \(\lambda\text{.}\)
Figure 3.56. On effectue une régression Ridge en dimension \(d=250\) avec des poids \(\theta\) générés par une loi \(\mathcal{N}(0,\frac{1}{\alpha})\text{,}\) avec 100 données d’entrainement. On trace les erreurs en fonction du coefficient de régularisation \(\lambda\text{.}\)
Les figures - montre que la régularisation Ridge que quel soit la dimension et \(\alpha\) augmente l’erreur sur le jeu de données d’entraînement. Sur le jeu de données tests (pour la généralisation) si on est loin de l’hypothèse des poids petits (\(\alpha\) grand) on peut déjà gagner un peu en grande dimension. Lorsque \(\alpha\) est grand on obtient en grande et moyenne dimension on obtient des erreurs nettement plus faible que pour (\(\alpha\) petit) avec la régression Ridge. Cela correspond bien au fait que \(\alpha\) grand correspond à l’a priori Gaussien fait par la régression Ridge.

Subsubsection 3.3.2.2 Régression Lasso

On considère ici un modèle linéaire du type (3.7). On considère toujours des données normalisées. Si \(y\) et les \(x_i\) dans le modèle sont en \(O(1)\) on voit que pour \(d\) grand une autre solution est de considérer qu’un nombre important de poids sont nuls. Il s’agit d’une hypothèse de Parcimonie. On va voir qu’on peut proposer un autre a priori modélisant cela.
Définition 3.57.
On appelle a priori de Laplace le distribution à priori:
\begin{equation} p(\theta \mid \alpha) = \mathcal{L}(0,\frac{1}{\alpha})\tag{3.33} \end{equation}
avec \(\alpha\) un hyper paramètre de la loi et \(\mathcal{L}\) une loi de Laplace.
Preuve.
Admis
On ne donnera pas de résultat sur l’unicité de la solution du problème Lasso. Il n’est pas a priori évident de voir pourquoi le terme en norme \(l_1\) permet de favoriser la parcimonie. On Va le montrer sur un cas simple.
Preuve.
On développe \(\frac12\parallel b-A\theta\parallel_2^2\) ce qui donne
\begin{equation*} \frac12\parallel b-A\theta\parallel_2^2 = \frac{1}{2}(b,b)-(b,A\theta)+(\theta, A^t A\theta) = \frac{1}{2}(b,b)-(b,A\theta)+(\theta,\theta) \end{equation*}
puisque \(b\) ne dépend pas de \(\theta\text{,}\) le problème de régression Lasso revient à minimiser
\begin{equation*} \operatorname{argmin}_{\theta} -(b,A \theta)+(\theta,\theta) + \lambda \sum_{i=1}^{d+1}\mid\theta_i\mid \end{equation*}
On pose \(\theta^{mc}\) la solution du problème des moindres carrés sans régularisation:
\begin{equation*} \theta^{mc}= (A^tA)^{-1}A^tb_y=A^tb \end{equation*}
Par le produit scalaire on obtient
\begin{equation*} \operatorname{argmin}_{\theta} -(A^tb,\theta)+(\theta,\theta) + \lambda \sum_{i=1}^{d+1}\mid\theta_i\mid \end{equation*}
Le problème de minimisation se réécrit donc
\begin{equation*} \operatorname{argmin}_{\theta} -(\theta^{mc},\theta)+(\theta,\theta) + \lambda \sum_{i=1}^{m+1}\mid\theta_i\mid \end{equation*}
On remarque que si \(\theta_{j}^{mc} \gt 0\) alors \(\theta_{j}\geq 0\) et si \(\theta_{j}^{mc}\lt 0\) alors \(\theta_{j}\leq 0\text{.}\) La fonction n’étant pas dérivable en zéro on sépare les cas:
  1. \(\theta_{j}^{mc} \gt 0\text{:}\) dans ce cas la fonction minimisée est
    \begin{equation*} f(\theta_j)= -\theta_{j}^{mc}\theta_j+\frac12\theta_j^2+\lambda \theta_j \end{equation*}
    donc en résolvant \(f^{'}(\theta_j)=0\) implique
    \begin{equation*} \theta_j= \theta_{j}^{mc}-\lambda \end{equation*}
    Comme \(\theta_j \geq 0\) on a \(\theta_j= max(0,\theta_{j}^{jmc}-\lambda)\text{.}\)
  2. \(\theta_{j}^{mc} \lt 0\text{:}\) Par le même raisonnement on obtient
    \begin{equation*} \alpha_j= max(0,-\theta_{j}^{mc}-\lambda) \end{equation*}
    donc à la fin
    \begin{equation*} \theta_j = sign(\theta_{j}^{mc})max(0,\mid\theta_{j}^{mc}\mid-\lambda) \end{equation*}
    On voit donc que plus \(\lambda\) est large plus le nombre de coefficients nuls sera important.
Figure 3.60. Régression Ridge/Lasso en dimension \(d\text{.}\) Erreur de généralisation Ridge (haut gauche). Norme des poids (Haut droite). Erreur de généralisation Lasso (bas gauche). Proportion de poids non nuls (bas droite). \(n=100\text{.}\) Les données sont générés avec le modèle
\begin{equation*} y=x_2-2x_4+ 1.5x_{d-2} \end{equation*}
la figure montre que la régularisation Lasso permet de bien capturer la solution en cas de parcimonie. Ici la solution de référence à trois éléments non nuls.

Subsubsection 3.3.2.3 Dropout

Pour finir on va introduire une autre façon d’éviter le surapprentissage appelée le dropout. On se place dans le cadre non régularisé que l’on résout avec une méthode itérative de type gradient stochastique.
La meilleure façon de "régulariser" un modèle de taille fixe est de faire une moyenne pondérée des prédictions de tous les jeux de paramètre possible. On peut essayer de faire cela pour des petits modèles. C’est évidemment impossible dans le cas des grands modèles. Il s’agit d’abandonner, de désactiver certains paramètres du modèle temporairement pendant les itérations du gradient. En changeant régulièrement les paramètres gelés, on retrouve à la fin un effet de moyenne. On se rappelle que le surapprentissage risque d’arriver lorsque \(d\) le nombre de paramètres est supérieur à \(n\text{.}\) On voit que le dropout consiste à baisser le nombre de paramètres libres à chaque itération. Pour la régression:
  1. Si on gèle une proportion \(P\) de paramètres à chaque itération, c’est comme si on était sur un problème de dimension \((1-P) (d+1)\text{.}\) On peut donc se retrouver dans des cas ou \((1-P) (d+1) \lt n\) et donc pas de surapprentissage.
  2. Chaque paramètre à une probabilité \(p\) (loi de Bernoulli) d’être temporairement abandonné.
  3. Pour la régression il suffit de mettre à zéro une sortie pour l’abandonner.
On peut le formaliser sous la forme suivante.
Définition 3.61.
On définit \(R\) une matrice telle que
\begin{equation*} R_{ij} = \mathcal{B}(p) \end{equation*}
avec \(\mathcal{B}\) une loi de Bernouilli. La dropout s écrit donc
\begin{equation*} R * \mathcal{X} \end{equation*}
avec \(*\) le produit terme à terme. Le problème aux moindres carrés avec Dropout s’écrit:
\begin{equation} \theta^* =\operatorname{argmin}_{\theta} \frac12\Norm b-R *A\theta\Norm_2^2\tag{3.35} \end{equation}
Maintenant on va montrer que cette méthode peut être vue comme une méthode de régularisation.
Preuve.
On part de la fonction coût:
\begin{equation*} \Norm b- R* A\theta \Norm_2^2= \langle b,b\rangle- 2 \langle R*A \theta,b \rangle +\langle R*A \theta,R* A \theta \rangle \end{equation*}
On prend l’espérance pour obtenir:
\begin{equation*} \langle b,b \rangle - 2 \langle \mathbb{E}[R*A]\theta,b \rangle + \langle \theta,\mathbb{E}[(R*A)^t (R*A)]\theta \rangle \end{equation*}
On commence par la première espérance:
\begin{equation*} (\mathbb{E}[R* A])_{ij}=\mathbb{E}_{R}[(R* A)_{ij}]= (A)_{ij}\mathbb{E}_{R}[R_{ij}]= p (A)_{ij} \end{equation*}
On va maintenant calculer la seconde:
\begin{equation*} ((R*A)^t(R*A))_{ij}=\sum_{k=1}^n R_{ki}R_{jk}(A)_{ki}(A)_{kj} \end{equation*}
On en déduit (car les lois sont indépendantes) que:
\begin{equation*} (\mathbb{E}[(R* A)^t (R* A))_{ij}= p^2 (A^t A)_{ij}, \quad \forall i\neq j \end{equation*}
et
\begin{equation*} (\mathbb{E}[(R* A)^tR* A])_{ij}= p (A^tA)_{ii}, \quad \forall i= j \end{equation*}
Maintenant on revient au calcul global, on utilise la 1ère espérance pour avoir:
\begin{equation*} \Norm b- R*A\theta \Norm_2^2 =\Norm b-pA\theta\Norm_2^2 -p^2 \langle A_x\theta ,A_x \theta \rangle + \langle \theta,\mathbb{E}[(R*A)^t(R*A)]\theta \rangle \end{equation*}
On obtient donc
\begin{equation*} \Norm b- R*A\theta \Norm_2^2 =\Norm b-pA\theta\Norm_2^2+ \langle \theta,(\mathbb{E}[(R* A)^t(R* A)]-p^2(A^tA))\theta \rangle \end{equation*}
En utilisant les calculs sur l’espérance \(\mathbb{E}[(R* A)^t(R* A)]\) on voit que
\begin{equation*} (\mathbb{E}[(R* A)^t(R* A)]-p^2(A^t A))=p(1-p)diag(A^t A) \end{equation*}
ce qui permet de conclure.
La méthode de Dropout revient a une méthode de régulation Ridge ou les poids associés aux plus grosses entrées sont les plus régularisés à travers la matrice \(\Gamma\text{.}\)

Subsection 3.3.3 Régression Bayesienne

Jusqu’a présent on a toujours construit le meilleur au modèle dans un certain sens (celui qui maximise la vraisemblance ou le maximum a posteriori). en choisissant les meilleurs paramètres \(\theta\text{.}\) L’idée de la régression linéaire Bayesienne vient du constat que peut plusieurs modèles sont très bon (plusieurs choix de \(\theta\)) et que certains un peu moins bons que le \(\theta_{opt}\text{,}\) mais peut être seront plus robustes au surapprentissage. Pour cela on va donc construire la la loi des paramètres et pas juste le meilleur.
Principe: Par rapport à l’approche précédente cela va revenir à ne pas chercher un jeu unique de poids \(\theta=(\omega,b)\) qui maximise la loi de probabilité a posteriori
\begin{equation*} p(\theta \mid x) \end{equation*}
mais la loi en elle-même.
Comme précédemment on se donne:
  1. Loi a priori \(p(\theta \mid \alpha) = \mathcal{N}(0,\frac{1}{\alpha})\text{.}\)
  2. Loi d’échantillonnage:
    \begin{equation*} p_{\theta}(y \mid x, \theta)=\mathcal{N}((\omega,x) + b,\sigma^2) \end{equation*}
    et donc
    \begin{equation*} p_{\theta}(\mathcal{Y} \mid \mathcal{X}, \theta, b)=\Pi_{i=1}^n\mathcal{N}(y_i\mid (\theta,x_i) + b,\sigma^2) \end{equation*}
    Cela revient a dire \(b_y = A_x \omega +\epsilon\) avec \(b\in \mathbb{R}^{n}\) le vecteur des données cibles et \(A \in \mathcal{M}_{n,d+1}(\mathbb{R})\) et \(\epsilon\) un variable aléatoire de loi \(\mathcal{N}(0,\sigma^2)\text{.}\)
L’objectif est donc de calculer la loi a posteriori.

Preuve.

Cette preuve est basée sur les propriétés des vecteurs Gaussiens. On remarque que les lois d’échantillonnage et a priori sont Gaussiennes. On écrit le vecteur aléatoire
\begin{equation*} z=\begin{pmatrix} \theta \\ y \end{pmatrix} \end{equation*}
avec \(\theta= (\omega, b)^t\) les paramètres du modèle linéaire. On peut écrire ce vecteur aléatoire sous la forme
\begin{equation*} \begin{pmatrix} \theta \\ y \end{pmatrix}= \begin{pmatrix} I_d \amp 0\\ A_x \amp I_d \end{pmatrix} \begin{pmatrix} \theta \\ \epsilon \end{pmatrix} \end{equation*}
On sait par hypothèse que \(\epsilon\) et \(\theta\) sont des variables aléatoires Gaussiennes donc le vecteur des deux est aussi Gaussien. Ensuite on applique une transformation linéaire. Or une transformation linéaire d’un vecteur aussi est encore un vecteur Gaussien. Notre vecteur est donc Gaussien. En utilisant la loi d’une transformation linéaire d’un vecteur Gaussien on obtient que
\begin{equation*} \mathbb{E}[z]= \begin{pmatrix} \mathbb{E}[\theta] \\ \mathbb{E}[y] \end{pmatrix}= \begin{pmatrix} \mathbb{E}[\theta] \\ \mathbb{E}[A \theta + \epsilon] \end{pmatrix}= \begin{pmatrix} \mathbb{E}[\theta] \\ A\mathbb{E}[\theta] + \mathbb{E}[\epsilon] \end{pmatrix}= \begin{pmatrix} 0 \\ 0 \end{pmatrix} \end{equation*}
et
\begin{equation*} \Sigma_z= \begin{pmatrix} \operatorname{cov}(\theta,\theta) \amp \operatorname{cov}(\theta,A \theta+\epsilon) \\ \operatorname{cov}(A\theta+\epsilon,\theta) \amp \operatorname{cov}(A\theta+\epsilon,A\theta+\epsilon) \end{pmatrix} \end{equation*}
Puisqu’on sait que les paramètres \(\theta\) suivent une Gausssienne \(\mathcal{N}(0,\frac{1}{\alpha})\) (loi a priori) on à
\begin{equation*} \operatorname{cov}(\theta,\theta) =\frac{1}{\alpha} I_d \end{equation*}
On passe ensuite au terme
\begin{equation*} \operatorname{cov}(A\theta+\epsilon) = \operatorname{cov}(\theta,A\theta)+\operatorname{cov}(\theta,\epsilon)= \operatorname{cov}(\theta,A\theta) \end{equation*}
car le bruit est indépendant des paramètres du modèle. Puisque la covariance est une application bilinéaire on a
\begin{equation*} \operatorname{cov}(\theta,A\theta+\epsilon) = \operatorname{cov}(\theta,A\theta) = \frac{A^t}{\alpha} \end{equation*}
Le troisième terme est le transposé. Il faut maintenant estimer le dernier terme:
\begin{equation*} \operatorname{cov}(A\theta+\epsilon,A\theta+\epsilon)= \operatorname{Var}(A\theta) + \operatorname{Var}(\epsilon) = \frac{A A^t}{\alpha} + \sigma^2I_d \end{equation*}
On a donc la matrice final de covariance pour la loi jointe
\begin{equation*} \Sigma_z= \frac{1}{\alpha}\begin{pmatrix} I_d \amp A^t \\ x_1 \amp AA^t +\alpha\sigma^2 I_d \end{pmatrix} \end{equation*}
On conclut en appliquant la formule qui permet à partir de la densité jointe de construire la densité conditionnelle pour un vecteur Gaussien.
Maintenant que la loi a posteriori et donc notre modèle est construit on souhaitera faire une prédiction. On va donc, pour finir, voir comment utiliser le modèle pour faire une prédiction. Pour faire une prédiction il faut calculer la loi a posteriori de prédiction:
\begin{equation*} p(x \mid \mathcal{X}, \alpha) = \int p(x\mid \theta)p(\theta \mid \mathcal{X},\alpha) d \theta \end{equation*}

Preuve.

Soit \(x\) une donnée observable. On nomme \(x_1\) le vecteur augmenté d’une ligne contenant un \(1\text{.}\) On se rappelle qu’on a choisi la loi d’échantillonnage suivante
\begin{equation*} y = (\theta,x_1)+\epsilon \end{equation*}
Donc
\begin{equation*} \mathbb{E}[y \mid \mathcal{X}, \alpha]= \mathbb{E}[\theta^t x_1+\epsilon \mid \mathcal{X}, \alpha] = \mathbb{E}[\theta \mid \mathcal{X}, \alpha] + \mathbb{E}[\epsilon]= \mu_{pos}^t \end{equation*}
On utilise le même argument pour calculer la variance de la combinaison affine de vecteurs aléatoires. Cela permet de conclure.
En général pour la prédiction on prend la moyenne de la loi a posteriori de prédiction et la variance nous permet de construire des intervalles de confiances pour nos prédictions. On peut évidemment échantillonner directement la loi. On remarque que
\begin{equation*} \mu_{pos}^t= \left(\frac{A_x A_x^t}{\alpha}+\sigma^2 I_d\right)^{-t} \frac{b_x}{\alpha} b_y \end{equation*}
ce qui donne
\begin{equation*} \mu_{pos}^t= \left(\frac{A_x^t A_x}{\alpha}+\sigma^2 I_d\right)^{-t} \frac{b_x}{\alpha} b_y \end{equation*}
On voit que l’on retrouve la solution de la régression Ridge déterministe de paramètre \(\lambda=\lambda^2\alpha\text{.}\) Par conséquent on voit que la régression Bayesienne donne en moyenne la régression Ridge, mais nous donne accès à la variance de l’estimateur.