Fondations : estimateur, MSE, CRLB

Cette page pose les outils pour juger un estimateur, sans encore dire comment en construire un. La section 1 définit ce qu'est un estimateur et montre, sur un exemple, qu'il s'agit d'une variable aléatoire. La section 2 en tire les trois quantités qui le caractérisent, biais, variance et MSE. La section 3 montre comment les estimer par simulation lorsque le calcul analytique n'aboutit pas.

1. Problématique

⚠ Bascule de notation

Le chapitre 1 écrivait $X$ (majuscule) la variable aléatoire et $x$ (minuscule) sa réalisation, une distinction utile pour apprendre à voir le hasard. À partir de maintenant, nous suivons la convention de Kay : une seule casse, $x[n]$, pour la variable aléatoire et sa réalisation ; le contexte tranche. De même $\hat{\theta }$ désignera aussi bien l'estimateur (une fonction des données, donc une v.a.) que l'estimation qu'il produit sur un jeu de données particulier (un nombre).

Modèle d'estimateur

Considérons que nous disposons de données $\mathbf{x}=[x[0],\dots ,x[N-1]{]}^T$, modélisées par une loi $p(\mathbf{x};\theta )$ à un paramètre $\theta$ inconnu. Un estimateur de $\theta$ est une fonction des données censée approcher la valeur de $\theta$. Mathématiquement, un estimateur peut s'exprimer sous la forme d'une fonction :

$$\hat{\theta }=g(\mathbf{x})$$

où

• $g$ est une fonction des données, choisie par nous,
• $\hat{\theta }$ est lui-même une variable aléatoire, puisque fonction de données aléatoires : il possède une loi, une espérance et une variance.

Confusion

Attention à ne pas confondre la valeur vraie $\theta$ et son estimation $\hat{\theta }$.

Exemple : estimation de la borne supérieure

Une machine produit des pièces dont la cote varie uniformément entre $0$ et une valeur maximale $\theta$ qui lui est propre ; nous mesurons $N$ pièces et cherchons $\theta$. Le modèle est donc $N$ mesures i.i.d. $x[n]\sim \mathcal{U}(0,\theta )$, une loi uniforme continue (fiche uniforme) dont le seul paramètre inconnu est aussi la borne supérieure du support.

Répétons trois fois la même expérience, huit mesures à chaque fois, et plaçons les valeurs obtenues sur un axe.

Figure 2.1a, Trois répétitions indépendantes de la même expérience, $N=8$ mesures tirées selon $\mathcal{U}(0,\theta )$ avec $\theta =10$ (script figures/ch2/01-echantillon-uniforme.py).

Les points se répartissent sans privilégier aucune zone, c'est le propre de la loi uniforme, et aucun ne dépasse $\theta$, qui borne le support. Chaque série suggère donc une réponse : $\theta$ se trouve quelque part au-dessus de la plus grande valeur observée. Mais cette plus grande valeur change d'une série à l'autre, $9,1$, $9,7$, puis seulement $6,6$, et elle reste toujours en deçà de $\theta$.

Trois candidats se présentent pour estimer $\theta$ :

$$\begin{array}{rl}{\hat{\theta }}_1& =\underset{n}{max}x[n],\\ {\hat{\theta }}_2& =\underset{n}{max}x[n]-\underset{n}{min}x[n],\\ {\hat{\theta }}_3& =2\overline{x}=\frac{2}{N}\sum _{n}x[n].\end{array}$$

${\hat{\theta }}_1$ mise sur l'observation la plus informative : aucune mesure ne peut dépasser $\theta$, donc la plus grande mesure observée doit en être proche. ${\hat{\theta }}_2$ retranche la plus petite mesure, dans l'idée de corriger un possible décalage. ${\hat{\theta }}_3$ mobilise toutes les mesures à la fois via la moyenne empirique (chapitre 1) : $E(\overline{x})=\theta /2$, donc ${\hat{\theta }}_3$ vise juste en espérance. Trois fonctions légitimes des données, donc trois variables aléatoires. Laquelle choisir ? Il faut un critère de comparaison.

Illustration

Fixons $\theta =10$ et tirons une seule série de 200 mesures ; à chaque $N$, recalculons les trois estimateurs sur les $N$ premières mesures de cette même série.

Figure 2.1b, Les trois estimateurs recalculés sur les $N$ premières mesures d'un même tirage ($\theta =10$), comparés à la valeur vraie (script figures/ch2/01-fondations.py).

Les trois courbes se rapprochent de $\theta =10$ quand $N$ grandit, mais à des vitesses très différentes, ${\hat{\theta }}_1$ colle à la cible presque immédiatement, ${\hat{\theta }}_3$ tâtonne plus longtemps. Ce graphe n'est qu'une seule réalisation pour une meme valeur de $\theta$ : un autre tirage déplacerait les trois courbes, sans changer leur allure générale. Comparer des estimateurs sur un tirage unique est donc fragile. Il faut un critère qui résume leur comportement sur l'ensemble des tirages possibles.

2. Métriques de performance : Biais, variance et MSE

Comment résumer, en un seul nombre, la qualité d'un estimateur sur l'ensemble de ses tirages possibles ? La mesure la plus naturelle est une distance moyenne à la cible, l'erreur quadratique moyenne, l'espérance du carré de l'écart entre l'estimateur et la vraie valeur.

Définition

Définition, Erreur quadratique moyenne (MSE)

$$\mathrm{MSE}(\hat{\theta })=E[(\hat{\theta }-\theta {)}^2].$$

où

• $\hat{\theta }=g(\mathbf{x})$ est l'estimateur évalué,
• $\theta$ est la valeur vraie, inconnue en pratique.

WARNING

La définition de $\mathrm{MSE}(\hat{\theta })$ considère ici que le paramètre $\theta$ est un paramètre déterministe mais inconnu (approche fréquentiste). Dans certaines parties du cours, nous assignerons une loi de probabilité à $\theta$ (approche Bayésienne)

Cette quantité peut se décomposer en deux termes de natures différentes.

Définition, Biais et variance

Pour un estimateur $\hat{\theta }$ de $\theta$ :

$$\mathrm{MSE}(\hat{\theta })=\mathrm{var}(\hat{\theta })+b(\hat{\theta }{)}^2$$

où

• $b(\hat{\theta })=E[\hat{\theta }]-\theta$ correspond au biais de l'estimateur,
• $\mathrm{var}(\hat{\theta })=E[(\hat{\theta }-E[\hat{\theta }]{)}^2]$ correspond à la variance de l'estimateur.

Preuve

Posons $m=E[\hat{\theta }]$ et insérons-le :

$$\hat{\theta }-\theta =(\hat{\theta }-m)+(m-\theta ).$$

Le second terme, $m-\theta$, est une constante (il ne dépend pas des données). En développant le carré et en prenant l'espérance, le terme croisé s'annule car $E[\hat{\theta }-m]=0$ :

$$\mathrm{MSE}(\hat{\theta })=\underset{\text{dispersion autour de }m}{\underbrace{E[(\hat{\theta }-m{)}^2]}}+2(m-\theta )\underset{=0}{\underbrace{E[\hat{\theta }-m]}}+\underset{\text{écart entre }m\text{ et }\theta }{\underbrace{(m-\theta {)}^2}}.◼$$

Le biais mesure une erreur systématique (nous visons à côté, toujours du même côté) ; la variance mesure une erreur de fluctuation (nous dispersons autour de notre propre cible). La MSE les réunit en un seul nombre, mais un nombre ne dit pas lequel des deux défauts pèse le plus.

Propriétés

En pratique, nous n'aurons pas de technique sur étagère pour choisir l'estimateur optimal $g(.)$. Néanmoins, dans la classe des estimateurs sans biais, il sera possible de comparer la variance de n'importe quel estimateur avec une borne présentant la variance optimale. Lorsqu'un estimateur sans biais atteint cette borne, il deviendra de facto optimal.

Théorème, Borne de Cramér-Rao (CRLB)

Pour tout estimateur sans biais $\hat{\theta }=g(\mathbf{x})$ de $\theta$ (sous conditions de régularité) :

$$\mathrm{var}(\hat{\theta })\ge \frac{1}{I(\theta )}.$$

où $I(\theta )$ correspond à l'information de Fisher et est définie par :

$$I(\theta )=-E\left[\frac{{\partial }^2\ell (\theta )}{\partial {\theta }^2}\right].$$

et $\ell (\theta )=\ln p(\mathbf{x};\theta )$ correspond à la log-vraisemblance.

Illustration

Cette décomposition se lit aussi sur un schéma classique : imaginez trois séries de tirs sur une cible, chaque tir étant une estimation $\hat{\theta }$ obtenue sur un tirage différent, et le centre marquant la valeur vraie $\theta$.

Figure 2.1c, Trois cibles, de gauche à droite : tirs groupés mais décalés du centre (biais fort, variance faible) ; tirs dispersés mais centrés en moyenne (sans biais, variance forte) ; tirs groupés et proches du centre malgré un léger décalage (biais faible, MSE la plus faible des trois) (script figures/ch2/01-fondations.py).

Le biais déplace le nuage de tirs loin du centre ; la variance l'étale autour de sa propre moyenne. Les cibles de gauche et du centre n'illustrent chacune qu'un seul défaut, et leur MSE reste élevée, celle de gauche à cause du biais, celle du centre à cause de la variance. La cible de droite cumule un peu des deux, un biais non nul et une variance non nulle, mais tous deux petits : sa MSE est la plus faible des trois. Un peu de biais, bien maîtrisé, peut valoir mieux qu'aucun biais du tout.

Exemple : estimation de la borne supérieure

En utilisant les propriétés de l'espérance mathématique, il est parfois possible d'obtenir analytiquement l'expression de la MSE, du biais et de la variance d'un estimateur.

Pour les trois candidats du §1 ($x[n]$ i.i.d. $\mathcal{U}(0,\theta )$), le tableau suivant présente les propriétés de chaque estimateur.

	${\hat{\theta }}_1=max$	${\hat{\theta }}_2=max-min$	${\hat{\theta }}_3=2\overline{x}$
Biais	$-{\displaystyle \frac{\theta }{N+1}}$	$-{\displaystyle \frac{2\theta }{N+1}}$	$0$
Variance	${\displaystyle \frac{N{\theta }^2}{(N+1{)}^2(N+2)}}$	${\displaystyle \frac{2(N-1){\theta }^2}{(N+1{)}^2(N+2)}}$	${\displaystyle \frac{{\theta }^2}{3N}}$
MSE	${\displaystyle \frac{2{\theta }^2}{(N+1)(N+2)}}$	${\displaystyle \frac{6{\theta }^2}{(N+1)(N+2)}}$	${\displaystyle \frac{{\theta }^2}{3N}}$

Calcul des moments de ${\hat{\theta }}_1$ et ${\hat{\theta }}_3$

Sous indépendance, $P({\hat{\theta }}_1\le x)=P(\text{tous les }x[n]\le x)=F_X(x{)}^N=(x/\theta {)}^N$ sur $[0,\theta ]$, la CDF du maximum s'obtient en élevant celle d'une seule mesure à la puissance $N$ (fonction de répartition). En dérivant, la densité est $f_{(N)}(x)=N{x}^{N-1}/{\theta }^N$, d'où

$$E({\hat{\theta }}_1)={\int }_0^{\theta }x\frac{N{x}^{N-1}}{{\theta }^N}\mathrm{d}x=\frac{N}{N+1}\theta ,E({\hat{\theta }}_1^2)=\frac{N}{N+2}{\theta }^2,$$

et $\mathrm{var}({\hat{\theta }}_1)=E({\hat{\theta }}_1^2)-E({\hat{\theta }}_1{)}^2=N{\theta }^2/((N+1{)}^2(N+2))$ par König-Huygens. Pour ${\hat{\theta }}_3$, le corollaire sur la loi de la moyenne empirique donne $\mathrm{var}(\overline{x})={\sigma }^2/N$ avec ${\sigma }^2={\theta }^2/12$ (fiche uniforme), d'où $\mathrm{var}(2\overline{x})=4{\sigma }^2/N={\theta }^2/(3N)$. Les moments de ${\hat{\theta }}_2$ s'obtiennent de façon analogue, à partir de la loi jointe du couple (min, max), nous les admettons ici.

Un seul de ces estimateurs est sans biais, ${\hat{\theta }}_3$, et pourtant ce n'est pas lui qui présente la MSE la plus faible. En effet, $\mathrm{MSE}({\hat{\theta }}_1)$ et $\mathrm{MSE}({\hat{\theta }}_2)$ décroissent en $1/N^2$, bien plus vite que $\mathrm{MSE}({\hat{\theta }}_3)$, qui décroît seulement en $1/N$. Dès $N=3$, le maximum biaisé a donc une MSE plus petite que la moyenne sans biais, et l'écart se creuse ensuite. Entre les deux estimateurs biaisés, ${\hat{\theta }}_1$ l'emporte partout sur ${\hat{\theta }}_2$, d'un facteur $3$ exactement.

Illustration

Figure 2.1d, Biais, variance et MSE des trois estimateurs en fonction de $N$ ($\theta =10$), $N\ge 2$ ; variance et MSE en échelle log-log (script figures/ch2/01-fondations.py).

En échelle log-log, une décroissance en $1/N^k$ devient une droite de pente $-k$ : les MSE de ${\hat{\theta }}_1$ et ${\hat{\theta }}_2$ y sont deux fois plus pentues que celle de ${\hat{\theta }}_3$, l'écart en $1/N^2$ contre $1/N$ se lisant directement. Le biais reste en échelle linéaire, ses valeurs étant négatives et nulles pour ${\hat{\theta }}_3$.

3. Simulation de Monte-Carlo

Les trois définitions du §2 posent une difficulté pratique immédiate. Elles reposent toutes sur une espérance, rarement calculable analytiquement dès que l'estimateur s'écarte des quelques cas d'école. La simulation de Monte-Carlo contourne cet obstacle : elle génère $B$ réalisations de la variable aléatoire $\hat{\theta }$, puis remplace chaque espérance par une moyenne empirique sur ces $B$ valeurs.

Définition

Définition, Estimation Monte-Carlo du biais, de la variance et de la MSE

Soient ${\hat{\theta }}^{(1)},\dots ,{\hat{\theta }}^{(B)}$ les estimations obtenues sur $B$ répétitions indépendantes de l'expérience, et $\bar{\hat{\theta }}=\frac{1}{B}\sum _b{\hat{\theta }}^{(b)}$ leur moyenne empirique. Le biais, la variance et la MSE s'estiment par

$$\hat{b}=\bar{\hat{\theta }}-\theta ,$$$$\widehat{\mathrm{var}}=\frac{1}{B}\sum _{b=1}^B({\hat{\theta }}^{(b)}-\bar{\hat{\theta }}{)}^2,$$$$\widehat{\mathrm{MSE}}=\frac{1}{B}\sum _{b=1}^B({\hat{\theta }}^{(b)}-\theta {)}^2.$$

Ces trois expressions se justifient par la loi des grands nombres : une moyenne empirique se rapproche de l'espérance qu'elle estime à mesure que $B$ grandit. La qualité de l'estimation dépend donc directement de $B$. Une valeur trop faible laisse subsister des fluctuations, une valeur trop élevée allonge le temps de simulation.

Implémentation

La mise en œuvre suit toujours les trois mêmes étapes, quelle que soit la loi des données et quel que soit l'estimateur étudié :

1. simuler un jeu de $N$ mesures selon le modèle, avec une valeur de $\theta$ que nous fixons ;
2. estimer $\hat{\theta }$ sur ce jeu, en appliquant l'estimateur à évaluer ;
3. répéter $B$ fois, puis moyenner selon les trois formules ci-dessus.

Seules les deux premières étapes dépendent du problème traité. Elles s'isolent donc dans deux fonctions, simulate et g, que le reste du code ne fait qu'appeler :

import numpy as np
rng = np.random.default_rng(2026)

theta, N, B = 10.0, 8, 100_000

def simulate(theta, N):
    """Draw N measurements from the model, here x[n] ~ U(0, theta)."""
    return rng.uniform(0, theta, size=N)

def g(x):
    """The estimator under study: from data to an estimate of theta."""
    return x.max()

# steps 1 and 2, repeated B times: B estimates of theta
theta_hat = np.empty(B)
for b in range(B):
    x = simulate(theta, N)      # one repetition of the experiment
    theta_hat[b] = g(x)         # the corresponding estimate

# step 3: each expectation becomes an average over these B values
bias = theta_hat.mean() - theta
variance = theta_hat.var()
mse = ((theta_hat - theta) ** 2).mean()

Changer d'estimateur revient à changer le corps de g, changer de modèle à changer celui de simulate : le reste du code est inchangé. Lorsque l'estimateur s'exprime de façon vectorielle, comme ici avec max, il est plus rapide de remplacer la boucle par un unique tirage dans un tableau $B\times N$, puis d'appliquer l'estimateur ligne par ligne.

Exemple

Appliquons cette méthode à ${\hat{\theta }}_1=\underset{n}{max}x[n]$ pour $N$ variable, en confrontant sa MSE théorique à son estimation Monte-Carlo pour deux valeurs de $B$.

Figure 2.2, MSE de ${\hat{\theta }}_1=\underset{n}{max}x[n]$ en fonction de $N$ ($\theta =10$), valeur théorique $2{\theta }^2/((N+1)(N+2))$ et estimations Monte-Carlo pour $B=10$ et $B=100$, en échelle log-log (script figures/ch2/01-monte-carlo-convergence.py).

Les deux estimations suivent la décroissance théorique sans jamais s'en écarter systématiquement : elles oscillent autour de la courbe, tantôt au-dessus, tantôt en dessous. C'est la signature d'un estimateur sans biais de la MSE, et c'est ce qui rend la méthode fiable. Leur dispersion, en revanche, diffère nettement. À $B=10$, l'écart à la valeur théorique atteint $104\mathrm{\%}$ et vaut $46\mathrm{\%}$ en médiane : dix répétitions ne permettent pas de conclure. À $B=100$, il tombe à $17\mathrm{\%}$ en médiane. Multiplier $B$ par dix a donc divisé l'écart typique par $2,7$, proche du facteur $\sqrt{10}\approx 3,2$ attendu, la précision d'une moyenne empirique progressant en $1/\sqrt{B}$. Atteindre une décimale sûre exige de ce fait des $B$ de l'ordre de ${10}^4$, ce qui reste bon marché : le coût est du temps de calcul, jamais des mesures supplémentaires.

Cette technique dépasse largement l'exemple : dès qu'un estimateur est trop compliqué pour un calcul analytique, ce qui est la règle plutôt que l'exception, la simulation de Monte-Carlo reste disponible, et c'est alors le seul accès à son risque.