4.1 Inférence en régression Poisson

À partir d’un jeu de données, nous allons voir quelques méthodes d’estimation de modèles de données de comptage.

Nous déutons avec le modèle de régression le plus simple pour les données de comptage, soit la régression Poisson.

Pour les MLE de la Poisson, nous allons voir que:

  1. La convergence des estimateurs nécessite une juste spécificiation de la moyenne conditionnelle; il n’est pas nécessaire que la variable aléatoire soit Poisson;

  2. Une inférence valide basée sur les écarts-types du MLE nécessite une juste spécificiation de la moyenne et de la variance conditionnelles. L’équidispersion est donc nécessaire, mais il n’est pas nécessaire que la variable aléatoire soit Poisson. Pour fin de compréhension, on pourrait comparer l’hypothèse de l’équidispersion de la Poisson à celle de l’homoscédasticité de la normale pour la régression linéaire.

  3. Une inférence valide basée sur les écarts-types du MLE peut être faite en utilisant un MLE modifié de manière appropriée s’il n’y a pas d’équidispersion, mais la moyenne conditionnelle doit être correctement spécifiée;

  4. Des estimateurs plus efficaces que le MLE de la Poisson peuvent être utilisés si les données ne sont pas équidispersée.

La juste spécificiation de la moyenne conditionnelle, de la variance conditionnelle ou de la densité signifie que la forme de la moyenne et que les variables explicatives de la moyennes soient les mêmes que dans le processus de création des données.

4.1.1 Poisson MLE

Définition 4.1 La fonction de vraisemblance, notée \(L(\theta|\mathbf{k})\), est une fonction de \(\theta\) et dépend des observations \(\mathbf{k}\). Elle s’exprime comme:

\[\begin{eqnarray*} L(\theta|\mathbf{k}) = \prod_{j=1}^n f_{\theta}(k_j) \end{eqnarray*}\]

où \(f_{\theta}(k_j)\) est la fonction de densité (ou de probabilité) associée à la distribution dont on cherche à estimer les paramètres.

Intuitivement, il faut considérer \(L(\theta|\mathbf{k})\) comme la chance d’observer les valeurs \(k_1, k_2, ..., k_n\) de l’échantillon, pour une valeur déterminée \(\theta\) du paramètre.

Définition 4.2 L’estimateur du maximum de vraisemblance \(\hat{\theta}\) de \(\theta\) est obtenu en maximisant la probabilité d’observer les valeurs \(k_1, k_2, ..., k_n\):

\[\hat{\theta} = arg \max_{\theta \in \Theta} L(\theta|\mathbf{k}) \]

En pratique, il est souvent plus pratique de passer au logarithme avant l’étape de maximisation. On définit ainsi le log-vraisemblance comme:

\[\begin{eqnarray*} \ell(\theta|\mathbf{k}) = \ln L(\theta|\mathbf{k}) \end{eqnarray*}\]

Donc, le MLE, \(\hat{\theta}_{MLE}\), est la solution à:

\[\begin{eqnarray*} \frac{\delta \ell(\theta|\mathbf{k})} {\delta \theta} = 0 \text{ , un vecteur } q \times 1. \end{eqnarray*}\]

4.1.1.1 Propriétés du MLE

Sous certaines conditions normalement respectées, on peut montrer que l’estimateur est convergent, et que \(\hat{\theta}_{MLE}\) est sans biais:

\[E[\hat{\theta}_{MLE}] =\theta_0\]

avec \(\theta_0\) représentant la vraie valeur.

4.1.1.2 Distribution asymptotique

Le maximum de vraisemblance a la propriété asymptotique suivante:

\[\begin{eqnarray*} \sqrt{n} (\hat{\theta}_{MLE} - \theta_0) &\sim_.& Normal(0,\mathbf{A}^{-1}) \end{eqnarray*}\]

où la matrice \(\mathbf{A}\) de dimension \(q \times q\) est définie comme:

\[\begin{eqnarray*} \mathbf{A} &=& -\lim_{n \rightarrow \infty} \frac{1}{n} E\left[\sum_{i=1}^n \frac{\delta^2\ln f(k_i;\theta) }{\delta \theta \delta \theta'} \Bigg|_{\theta_0} \right] \end{eqnarray*}\]

ou encore, comme nous l’avons vu plus tôt, peut utiliser l’égalité \(\mathbf{A} = \mathbf{B}\) et:

\[\begin{eqnarray*} \mathbf{B} &=& \lim_{n \rightarrow \infty} \frac{1}{n} E\left[\sum_{i=1}^n \frac{\delta\ln f(k_i;\theta) }{\delta \theta} \frac{\delta\ln f(k_i;\theta) }{\delta \theta'} \Bigg|_{\theta_0} \right] \end{eqnarray*}\]

4.1.1.3 En pratique

Avec des données, on cherche à évaluer les résultats précédents. Plusieurs options sont possibles

  1. L’estimateur hessien évalue directement \(A\) avec \(\hat{\theta}\), sans prendre l’espérance;

\[\begin{eqnarray*} \mathbf{A} &=& -\frac{1}{n} \sum_{i=1}^n \frac{\delta^2\ln f(k_i;\theta) }{\delta \theta \delta \theta'} \Bigg|_{\hat{\theta}} \end{eqnarray*}\]

  1. L’estimateur de l’information de Fisher évalue \(A\) en prenant l’espérance, évaluée à \(\hat{\theta}\):

\[\begin{eqnarray*} \mathbf{A} &=& -\frac{1}{n} E\left[\sum_{i=1}^n \frac{\delta^2\ln f(k_i;\theta) }{\delta \theta \delta \theta'} \Bigg|_{\hat{\theta}} \right] \end{eqnarray*}\]

  1. L’estimateur du produit-croisé utilise \(B\) avec \(\hat{\theta}\), sans prendre l’espérance;

\[\begin{eqnarray*} \mathbf{B} &=& \frac{1}{n} \sum_{i=1}^n \frac{\delta\ln f(k_i;\theta) }{\delta \theta} \frac{\delta\ln f(k_i;\theta) }{\delta \theta'} \Bigg|_{\hat{\theta}} \end{eqnarray*}\]

  1. Une dernière forme, combinant \(A\) et \(B\) est appelé l’estimateur sandwich est de forme \(\mathbf{A}^{-1} \mathbf{B} \mathbf{A}^{-1}\). Nous y reviendrons.

Exemple 4.1 On suppose que \(Y_i \sim Poisson(\lambda_i)\), avec \(\lambda_i = \exp(\mathbf{X_i}' \mathbf{\beta})\). Développez les différentes manières de calculer la variance asymptotique des MLE de \(\beta\).

L’exemple est similaire à ce que nous avions vu au chapitre 2, mais simplement un peu plus formel.

On sait que:

\[\begin{eqnarray*} \frac{\delta \ell(\beta)}{\delta \beta} = \sum_{i=1}^n \mathbf{X}_{i} \left(y_i - \lambda_i \right) = 0. \end{eqnarray*}\]

Ainsi, la matrice \(\mathbf{A}\) se développe comme:

\[\begin{eqnarray*} \mathbf{A} = \lim_{n \rightarrow \infty} \frac{1}{n} \left[\frac{\delta^2 \ell(\beta)}{\delta \beta \delta \beta'} \right]= &=& \lim_{n \rightarrow \infty} \frac{1}{n} \sum_{i=1}^n \mathbf{X}_i \mathbf{X}_i^T \lambda_i \end{eqnarray*}\]

Pour le calcul de \(\mathbf{B}\) nous avons:

\[\begin{eqnarray*} \mathbf{B} = \lim_{n \rightarrow \infty} \frac{1}{n} \left[\left(\frac{\delta \ell(\mathbf{\beta})}{\delta \mathbf{\beta}}\right) \left(\frac{\delta \ell(\mathbf{\beta})}{\delta \mathbf{\beta}}\right)^T \right] &=& \lim_{n \rightarrow \infty} \frac{1}{n} \sum_{i=1}^n \left(\mathbf{X}_{i} \left(y_i - \lambda_i \right) \right) \left(\mathbf{X}_{i} \left(y_i - \lambda_i \right)\right)^T \\ &=& \lim_{n \rightarrow \infty} \frac{1}{n} \sum_{i=1}^n \mathbf{X}_{i} \mathbf{X}_{i}^T \left(y_i - \lambda_i \right)^2 \end{eqnarray*}\]

  1. L’estimateur hessien évalue directement \(A\) avec \(\hat{\theta}\), sans prendre l’espérance;

\[\begin{eqnarray*} \mathbf{A} &=& \frac{1}{n} \sum_{i=1}^n \mathbf{X}_i \mathbf{X}_i^T \hat{\lambda}_i \end{eqnarray*}\]

  1. Il n’y aucun \(y_i\) dans l’expression de \(\mathbf{A}\), ainsi prendre l’espérance de \(\mathbf{A}\) n’a aucune influence. L’estimateur de l’information de Fisher est donc le même que l’estimateur hessien.

\[\begin{eqnarray*} \mathbf{A} &=& \frac{1}{n} E[\sum_{i=1}^n \mathbf{X}_i \mathbf{X}_i^T \hat{\lambda}_i] = \frac{1}{n} \sum_{i=1}^n \mathbf{X}_i \mathbf{X}_i^T \hat{\lambda}_i \end{eqnarray*}\]

  1. L’estimateur du produit-croisé utilise \(B\) avec \(\hat{\theta}\), sans prendre l’espérance; \[\begin{eqnarray*} \mathbf{B} &=& \frac{1}{n} \sum_{i=1}^n \mathbf{X}_{i} \mathbf{X}_{i}^T \left(y_i - \lambda_i \right)^2 \end{eqnarray*}\] Sous l’hypothèse d’une Poisson, nous aurons \(E[\mathbf{B}] = \mathbf{A}\) car \(E[\left(y_i - \lambda_i \right)^2] = Var[Y_i] = \lambda_i\). Par contre, si nous utilisons \(\mathbf{B}\) avec \(\hat{\theta}\), sans prendre l’espérance, nous aurons un estimateur de la variance différent.

  2. Pour l’estimateur sandwich utilisant \(\mathbf{A}\) et \(\mathbf{B}\), nous y reviendrons.

4.1.2 Poisson Quasi-MLE

Un défaut majeur de l’estimateur par maximum de vraisemblance (MLE) est que nous devons supposer une forme précise de la fonction de densité (ou de probabilité). Avec des données réelles, dans notre cas, il serait fortement surprenant que les données soient exactement poissonniennes.

Deux articles importants:

  • White, H. (1982), Maximum Likelihood Estimation of Misspecified Models, Econometrica, 50, 1–25.
  • Huber, P.J. (1967), The Behavior of Maximum Likelihood Estimates under Nonstandard Conditions, in L. LeCam and J. Neyman, eds., Proceedings of the Fifth Berkeley Symposium on Mathematical Statistics and Probability, 221–234, Berkeley, University of California Press

ont étudié les propriétés du MLE lorsque la fonction de densité (ou de probabilité) est incorrectement définie. Dans une telle situation, l’estimateur obtenu est appelé quasi-MLE (QMLE) ou encore pseudo-MLE. Nous garderons QMLE pour le cours.

Généralement, une mauvaise spécification de la fonction de densité (ou de probabilité) de la variable aléatoire sous-jacente mène à un MLE qui sera incorrect. Les distributions de la famille exponentielle linéaire, dont fait partie la Poisson, sont une exception importante de cette conclusion. Par contre, la variance des estimateurs ne sera pas similaire à ce que nous venons de voir.

4.1.2.1 Construction du QMLE

Pour construire un QMLE, on doit réaliser que l’hypothèse que les données proviennent d’une Poisson est plus forte que nécessaire lorsque nous utilisons l’équation:

\[\begin{eqnarray*} \sum_{i=1}^n (y_i - \exp(\mathbf{X}_i' \mathbf{\beta})) \mathbf{X}_{i} = \mathbf{0}. \end{eqnarray*}\]

Pour les distributions de la famille exponentielle linéaire, il peut être montré que l’estimateur \(\mathbf{\beta}\) est convergent en ne supposant qu’une forme correcte de la moyenne. Intuitivement, on peut expliquer cette situation en voyant l’évidence que \(E[y_i|\mathbf{X}_{i}] = \exp(\mathbf{X_i}' \mathbf{\beta})\) et donc:

\[\begin{eqnarray*} E[ (y_i - \exp(\mathbf{X_i}' \mathbf{\beta})) \mathbf{X}_{i}] = \mathbf{0}. \end{eqnarray*}\]

Étant donnée cette robustesse, on peut ainsi toujours utiliser l’équation de premier ordre (provenant de la fonction de probabilité d’une Poisson avec un lien logarithmique), même si nous avons l’intuition que les données ne proviennent pas d’une Poisson. Dans ce cas, on peut montrer que:

\[\hat{\beta}_{QMLE} \sim_d N[\beta, Var[\hat{\beta}_{QMLE} ]],\]

avec

\[Var[\hat{\beta}_{QMLE} ] = \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1} \Big[ \sum_{i=1}^n \omega_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big] \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1}\]

où \(\omega_i = Var[y_i|\mathbf{X}_{i}]\). Il s’agit de l’estimateur sandwich vu plus tôt.

Ainsi, pour l’estimation par QMLE, il suffit de spécifier la forme de variance, \(\omega_i\), pour définir la variance de l’estimateur \(\hat{\beta}_{QMLE}\). Nous avons quelques possibilités.

4.1.2.2 Fonction de variance Poissonnienne;

Si \(Y_i\) est Poisson, alors \(\omega_i = \mu_i\) et donc

\[\begin{align*} Var_{QMLE}[\hat{\beta}] &= \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1} \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big] \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1}\\ &= \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1} \\ &= Var_{MLE}[\hat{\beta} ] \end{align*}\]

4.1.2.3 Fonction de variance de type NB1;

La forme de cette variance est \(\omega_i = \phi \mu_i\), et donc:

\[\begin{align*} Var_{NB1}[\hat{\beta}] &= \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1} \Big[ \sum_{i=1}^n \phi \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big] \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1}\\ &= \phi \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1} \\ &= \phi Var_{MLE}[\hat{\beta}] \end{align*}\]

Ainsi, une manière simple de traiter la surdispersion (ou la sous-dispersion) est d’utiliser le MLE d’une Poisson, et d’ensuite multiplier les variance des estimateurs pat \(\phi\).

Exemple 4.2 Proposez un estimateur de \(\phi\)

On sait que:

\[\omega_i = Var[Y_i|\mathbf{X}_{i}] = E[(Y_i - \mu_i)^2] = \phi \mu_i \]

et donc

\[\phi_i = E[(Y_i - \mu_i)^2 / \mu_i] \]

Ainsi, l’estimateur classique de \(\phi\) est:

\[\hat{\phi} = \frac{1}{n-k} \sum_{i=1}^n \frac{(y_i - \hat{\mu}_i)^2}{\hat{\mu}_i}\]

où on a divisé par \(n-k\) au lieu de \(n\) afin de corriger pour le nombre de degrés de liberté.

Un autre estimateur pertinent est

\[\hat{\phi} = \frac{\sum_{i=1}^n(y_i - \hat{\mu}_i)^2}{\sum_{i=1}^n \hat{\mu}_i}\] Il existe une multitude d’estimateur pour \(\phi\). Certains, comme nous le verrons, sont plus intéressants que d’autres dans un contexte d’assurance.

La correction pour tenir compte \(\phi\) est triviale et est souvent en option dans les logiciels d’estimation classique comme R. Nous verrons un exemple.

4.1.2.4 Fonction de variance de type NB2;

La forme de cette variance est \(\omega_i = \mu_i + \alpha \mu_i^2\), et donc:

\[\begin{align*} Var_{NB2}[\hat{\beta}] &= \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1} \Big[ \sum_{i=1}^n (\mu_i + \alpha \mu_i^2) \mathbf{X}_{i} \mathbf{X}_{i}^T \Big] \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1} \end{align*}\]

Exemple 4.2 Proposez un estimateur de \(\alpha\)

On sait que:

\[Var[y_i|\mathbf{X}_{i}] = E[(y_i - \mu_i)^2] = \mu_i + \alpha \mu_i^2\]

et donc

\[\alpha \mu_i^2 = E[(y_i - \mu_i)^2 - \mu_i] \]

et

\[\alpha = E[(y_i - \mu_i)^2 - \mu_i)/mu_i^2] \]

Ainsi, l’estimateur classique de \(\alpha\) est:

\[\hat{\alpha} = \frac{1}{n-k} \sum_{i=1}^n \frac{(y_i - \hat{\mu}_i)^2 - \hat{\mu}_i}{\hat{\mu}_i^2}\]

Mais encore une fois, plusieurs estimateurs sont possibles. Un autre estimateur est:

\[\hat{\alpha} = \frac{\sum_{i=1}^n (y_i - \hat{\mu}_i)^2 - \hat{\mu}_i}{\sum_{i=1}^n \hat{\mu}_i^2}\]

Classiquement, pour les cas de variance de forms NB1 ou NB2, on ne calcule par les variances des estimateurs de \(\phi\) et de \(\alpha\). Toutefois, il serait possible de les calculer en utilisant diverses méthodes statistiques (simulations, bootstrap, etc.).

4.1.2.5 Fonction de variance non-spécifiée

La variance de l’estimateur \(\beta\) peut être trouvée même si aucune forme pour \(\omega_i\) n’est proposée. Un estimateur robuste de la variance est:

\[Var_{UNS}[\hat{\beta}] = \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1} \Big[ \sum_{i=1}^n (y_i - \mu_i)^2 \mathbf{X}_{i} \mathbf{X}_{i}^T \Big] \Big[ \sum_{i=1}^n \mu_i \mathbf{X}_{i} \mathbf{X}_{i}^T \Big]^{-1},\]

évalué à \(\hat{\mu}_i\).