22 Statistique Bayésienne

L’objectif du point de vue bayésien est de quantifier l’incertitude sur le paramètre à estimer en le traitant comme une variable aléatoire. Au lieu de considérer un paramètre \(\theta\) fixé, on introduit une loi de probabilité \(\pi\) sur \(\Theta\) (la loi a priori), puis on met à jour cette connaissance à partir des observations via la formule de Bayes.

Dans cette section, on se place dans un cadre général : on observe une variable aléatoire \(X\) qui dépend du paramètre \(\theta\) ; la loi de \(X\) conditionnellement à \(\theta\) a pour densité (ou masse) \(p_\theta(x)\) par rapport à une mesure de référence.

22.1 Définitions

Définition 22.1 (Loi a priori) Une loi a priori est une mesure de probabilité \(\pi\) sur l’ensemble des paramètres \(\Theta\).

Cette loi mesure notre connaissance sur le paramètre \(\theta\) avant de voir les données. Typiquement, si l’on est sûr que le « vrai » paramètre est égal à \(0\), alors la loi a priori est une masse de Dirac en \(0\). Si l’on est sûrs que le « vrai » paramètre est dans un ensemble \(A \subset \Theta\), alors la loi a priori aura toute sa masse sur \(A\).

Étant donné une rélisation \(\theta\) de \(\pi\), on génère un échantillon \(X \sim p_\theta\). À ce stade, il y a deux niveaux d’aléas : 

l’aléa sur le paramètre \(\theta\),
puis l’aléa sur la variable \(X\) conditionnellement à \(\theta\).

Définition 22.2 (Loi marginale, loi a posteriori) La loi de \(X\) s’appelle la loi marginale de \(X\), et on la note abusivement \(p(X)\). Elle est donnée par \[p(X) = \int_\Theta p_\theta(X) \pi(\theta) d\theta.\]

La loi de \(\theta\) sachant \(X\) s’appelle la loi a posteriori et on la note abusivement \(p(\theta \mid X)\). La formule de Bayes dit que

\[p(\theta \mid X) = \frac{p_\theta(X) \pi(\theta)}{p(X)}.\]

La loi a posteriori donne en quelque sorte la distribution de probabilité sur le paramètre \(\theta\) qui est la plus compatible avec les données observées \(X\). Comme la loi marginale \(p(X)\) ne dépend pas du tout du paramètre, on écrit souvent l’équation ci-dessus avec une constante de proportionnalité, \[p(\theta \mid x) \propto p_\theta(x) \pi(\theta), \] le symbole \(\propto\) signifiant « est proportionnel à ». Cette formule s’interprète de la façon classique : \[\mathsf{a~posteriori} \propto \mathsf{vraisemblance} \times \mathsf{a~priori}.\]

La loi a posteriori permet d’estimer \(\theta\) de plusieurs façons possibles : moyenne, médiane, mode, etc.

Constante de proportionnalité. Nous y reviendrons plus tard, mais certains calculs (comme l’estimateur MAP) ne nécessitent pas la connaissance de la constante de proportionnalité, et ne nécessitent même pas que \(\pi\) soit une mesure de probabilité.

Calcul de l’a posteriori. En revanche, le calcul de la distribution a posteriori nécessite, lui, de connaître la loi marginale \(p(X)\). Cette loi marginale est égale à une intégrale, \[p(x) = \int_\Theta p_\theta(x) \pi(\theta) d\theta.\] Cette intégrale peut être très difficile à calculer, et dans ce cas les statisticiens ont recours à des approximations numériques, via des méthodes de Monte-Carlo.

Remarque. Dans la construction bayésienne décrite ci-dessus, on doit considérer \(p_\theta\) comme une mesure de probabilité qui est elle-même aléatoire, puisque \(\theta\) l’est. Or, cela nécessiterait d’avoir défini ce qu’est une mesure aléatoire, c’est-à-dire d’avoir doté l’espace des mesures de probabilités d’une tribu, et de vérifier la mesurabilité de \(\theta \mapsto p_\theta\). Le même problème se pose quand on parle de loi conditionnelle : si l’espérance conditionnelle est définie rigoureusement comme une variable aléaoire, alors la loi conditionnelle est elle-même définie comme une variable aléatoire sur l’espace des mesures.

Il faut bien comprendre que, pour 99.9% des applications de la statistique bayésienne, ces considérations techniques sont parfaitement inutiles et pourraient être contournées d’une façon ou d’une autre. Les difficultés ne se poseront que dans des cas très particuliers, comme par exemple si on étudie des objets plus complexes comme les mouvements browniens.

Exemple 22.1 (Bernoulli) Plaçons-nous dans le cadre d’un modèle de Bernoulli \(\mathrm{Ber}(\theta)\), avec paramètre \(\theta \in [0,1]\). On utilise un a priori uniforme sur \([0,1]\), c’est-à-dire \(\pi(\theta) = 1\). Si l’on observe \(X_1, \dots, X_n\) i.i.d. de loi \(\mathrm{Ber}(\theta)\), la vraisemblance des observations est \(\theta^{s}(1-\theta)^{n-s}\), où \(s = x_1 + \dotsc + x_n\) est le nombre de succès. La loi marginale \(p(x_1,\dotsc, x_n)\) est donnée par \[\int_0^1 u^{s}(1-u)^{n-s} du\] et on reconnaît la fonction bêta, \(B(s+1, n-s+1)\). La loi a posteriori \(p(\theta \mid x_1,\dotsc, x_n)\) est donc donnée par \[\frac{\theta^{s}(1-\theta)^{n-s}}{B(s+1, n-s+1)}.\] Conditionnellement aux observations, le paramètre \(\theta\) est donc une variable aléatoire qui suit une loi bêta \(B(s+1, n-s+1)\).

22.2 Estimateurs a posteriori

Une fois que l’on a calculé la loi a posteriori \(p(\theta \mid X)\), on peut en déduire des quantités d’intérêt. Par exemple, dans l’exemple de Bernoulli que nous avons développé plus haut, on peut calculer le \(\theta\) qui est “le plus probable” conditionnellement aux observations, c’est-à-dire le mode de la loi a posteriori. Cet estimateur est appelé estimateur MAP (Maximum A Posteriori). Pour calculer l’estimateur MAP dans le modèle donné dans l’exemple plus haut (Exemple 22.1), il suffit de maximiser en \(\theta\) la vraisemblance, et on trouve \[\hat{\theta}_{\text{MAP}} = \frac{s}{n}. \tag{22.1}\]

Dans le cas général, le théorème suivant décrit plusieurs quantités centrales en pour utiliser la loi a posteriori afin d’estimer \(\theta\).

Dans le théorème, on suppose pour simplifier que tous ces objets sont bien définis et existent (en pratiquent, cela revient à dire que la loin a posteriori est \(L^2\) pour la MMSE ou \(L^1\) pour la MAVE).

Théorème 22.1 \(~\)

L’estimateur qui maximise la densité a posteriori est le mode a posteriori (MAP), \[\hat{\theta}_{\text{MAP}} = \arg \max_{\theta} p(\theta \mid X).\]
L’estimateur qui minimise l’erreur quadratique moyenne (MMSE, mean square error) est la moyenne a posteriori, \[\hat{\theta}_{\text{MMSE}} = \mathbb{E}[\theta \mid X].\]
En dimension 1, l’estimateur qui minimise la distance en valeur absolue (MAVE, minimum absolute value error) est la médiane a posteriori,

\[\hat{\theta}_{\text{MAVE}} = \text{médiane}(p(\theta \mid X)).\]

Preuve. Le premier point est une évidence. Pour le second, il s’agit de vérifier que la moyenne d’une variable aléatoire \(X\) est le point qui minimise l’erreur quadratique moyenne, c’est-à-dire que \[\mathbb{E}[X] = \arg \min_{x} \mathbb{E}[|X-x|^2].\] En écrivant la condition du premier ordre et en dérivant sous l’espérance, on voit immédiatement que l’unique point critique est donné par \(x = \mathbb{E}[X]\).

Quant au troisième point, il s’agit de vérifier que la médiane d’une variable aléatoire \(X\) est le point qui minimise la distance en valeur absolue \(\mathbb{E}[|X-t|]\). Soit \(m\) une médiane, dont on rappelle qu’elle vérifie \[\mathbb{P}(X \leqslant m) \geqslant \frac{1}{2} \quad \text{et} \quad \mathbb{P}(X \geqslant m) \geqslant \frac{1}{2}.\]

Pour n’importe quel \(t\), on va étudier \(\delta = |X-t| - |X-m|\). Sans perte de généralité, on peut suppose que \(t\) est plus grand que \(m\), l’autre sens étant similaire. Si \(x\) est plus grand que \(t\), alors \(\delta\) vaut \(m-t\). Si \(x\) est plus petit que \(m\), alors \(\delta\) vaut \(t-m\). Et si \(x\) est entre les deux, alors \(\delta\) vaut \(2x - t-m\). Dans tous les cas, \(|x-t| - |x-m|\) est plus grand que
\[(t-m)(\mathbf{1}_{X \leqslant m} - \mathbf{1}_{X > m}).\] En particulier, \[\mathbb{E}[|X-t|] - \mathbb{E}[|X-m|]\] est plus grand que \[(t-m)(\mathbb{P}(X \leqslant m) - \mathbb{P}(X > m)).\] Le second terme du produit est aussi égal à \(2\mathbb{P}(X \leqslant m) - 1\), donc par définition de la médiane il est plus grand que 0, ce qui montre que la médiane est un minimiseur de \(\mathbb{E}[|X-t|]\).

Normalisation de l’a priori. Le calcul de l’etimateur MAP ne nécessite pas la connaissance de la loi marginale \(p(X)\), ni même la connaissance totale de la loi a priori \(\pi\). En fait, on peut connaître \(\pi\) seulement à une constante de proportionnalité près : en effet, supposons que \(\pi\) ne soit pas normalisée et que \(\int \pi = c\), pas forcément égal à 1, et notons \(\pi_0 = \pi / c\) la loi normalisée. Il est clair que \[\arg \max_{\theta} p_\theta(x)\pi(\theta) = \arg \max_{\theta} p_\theta(x) \pi_0(\theta).\]

Cette remarque n’est pas si anodine : elle implique que le calcul de l’estimateur MAP ne nécessite pas que \(\pi\) soit une mesure normalisée… et en réalité, ne nécessite même pas que ce soit une mesure de probabilité. En effet, rien ne nous empêche de choisir un « a priori impropre », c’est-à-dire non intégrable. C’est typiquement le cas sur \(\mathbb{R}\) : prenons l’exemple des lois gaussiennes \(\mathscr{N}(\mu, 1)\). Si l’on n’a aucune connaissance a priori sur la localisation de \(\mu\), c’est que \(\mu\) pourrait se trouver n’importe où dans \(\mathbb{R}\). Il paraît donc naturel de choisir un a priori « uniforme sur \(\mathbb{R}\), à savoir \(\pi(\mu) = 1\). Dans ce cas, l’estimateur MAP devient \[\arg \max_{\mu} p_\mu(x)\] c’est-à-dire l’estimateur du maximum de vraisemblance.

Exemple 22.2 (Bernoulli) Dans le cadre du modèle \(\mathrm{Ber}(\theta)\) avec a priori uniforme, on a vu que la loi a posteriori est \(B(s+1, n-s+1)\). Dans ce cadre, l’estimateur MAP est égal à la moyenne empirique (cf Équation 22.1).

Or, l’estimateur MMSE n’est pas égal à l’estimateur MAP. En effet, le théorème ci-dessus dit que \(\hat{\theta}_{\text{MMSE}} = \mathbb{E}[\theta \mid X]\), et conditionnellement à \(X\), \(\theta\) suit une loi bêta. La formule pour l’espérance d’une loi bêta dit alors que \[\hat{\theta}_{\text{MMSE}} = \frac{s+1}{n+2}.\]

On pourrait aussi calculer l’estimateur MAVE ; cependant, il n’y a pas de formule simple pour la médiane d’une loi bêta.