4.3 Inférence en régression binomiale négative
Au lieu de supposer une distribution Poisson et de transformer ses propriétés afin d’admettre de la surdispersion, une alternative est de simplement supposer une distribution de comptage qui admet la possibilité que les données soient surdispersées.
La famille de distribution de la binomaile négative est une possibilité.
4.3.1 Binomiale négative 2
Une première option est la NB2. Tel que nous l’avons vu, la fonction de probabilité de la NB2 peut s’exprimer comme:
\[\begin{eqnarray*} \Pr(Y_i = y_i) = \frac{\Gamma( y_i + \alpha^{-1})}{\Gamma( y_i+1) \Gamma(\alpha^{-1})}\Big(\frac{\lambda_i}{\alpha^{-1}+\lambda_i}\Big)^{ y_i} \Big(\frac{\alpha^{-1}}{\alpha^{-1}+\lambda_i}\Big)^{\alpha^{-1}}. \label{NB2} \end{eqnarray*}\]
L’utilisation de \(\alpha^{-1}\) dans la fonction de probabilité (au lieu de \(\alpha\), tel que nous l’avions vu) permet de mieux voir que le modèle converge vers une Poisson pour \(\alpha \rightarrow 0\). Cette reparamétrisation nous sera utile plus tard.
4.3.1.1 Fonction de logvraisemblance
Les estimateurs de MLE pour la NB2 peuvent s’obtenir en dérivant la fonction de logvraisemblance. Nous avons donc:
\[\begin{eqnarray*} L(\beta, \alpha) &=& \prod_{i=1}^n \Pr(Y_i = y_i) = \prod_{i=1}^n\Big( \frac{\Gamma( y_i + \alpha^{-1})} {\Gamma(\alpha^{-1})} \Big) \frac{1}{ y_i!} \Big(\frac{\lambda_i}{\alpha^{-1}+\lambda_i}\Big)^{ y_i} \Big(\frac{\alpha^{-1}}{\alpha^{-1}+\lambda_i}\Big)^{\alpha^{-1}} \end{eqnarray*}\]
Il peut être difficile de travailler avec les fonctions gamma. Ainsi, sachant qu’on peut montrer que si \(y\) est entier:
\[\frac{\Gamma(y+a)}{\Gamma(a)} = \prod_{j=0}^{y-1} (j+a),\]
la fonction de vraisemblance devient:
\[\begin{eqnarray*} L(\beta, \alpha) &=& \prod_{i=1}^n\Big(\prod_{j=0}^{ y_i -1}(j + \alpha^{-1})\Big) \frac{1}{ y_i!} \Big(\frac{\lambda_i}{\alpha^{-1}+\lambda_i}\Big)^{ y_i} \Big(\frac{\alpha^{-1}}{\alpha^{-1}+\lambda_i}\Big)^{\alpha^{-1}} \end{eqnarray*}\]
La fonction de logvraisemblance, qu’on doit dériver pour obtenir les MLE de \(\beta\) et \(\alpha\) s’exprime comme:
\[\begin{eqnarray*} \ell(\beta, \alpha) &=& \ln(L(\beta, \alpha)) \\ &=& \sum_{i=1}^n \bigg[\Big(\sum_{j=0}^{ y_i -1}\ln(j + \alpha^{-1})\Big) -\ln(y_i!) - ( y_i + \alpha^{-1}) \ln(\alpha^{-1}+\lambda_i) + y_i \ln(\lambda_i) +\alpha^{-1}\ln(\alpha^{-1}) \bigg] \end{eqnarray*}\]
4.3.1.2 MLE - \(\beta\)
En prenant un lien logarithmique, avec \(\lambda_i = \exp(\mathbf{X}_{i}' \mathbf{\beta})\), nous cherchons tout d’abord à estimer \(\beta\).
\[\begin{eqnarray*} \frac{\delta \ell}{\delta \beta} &=& \frac{\delta}{\delta \beta} \sum_{i=1}^n \bigg[\Big(\sum_{j=0}^{ y_i -1}\ln(j + \alpha^{-1})\Big) -\ln(z_i!) - ( y_i + \alpha^{-1}) \ln(\alpha^{-1}+\lambda_i) + y_i \ln(\lambda_i) +\alpha^{-1}\ln(\alpha^{-1}) \bigg]\\ &=& \sum_{i=1}^n \left( - \frac{\delta}{\delta \beta}( y_i + \alpha^{-1}) \ln(\alpha^{-1}+\lambda_i) + \frac{\delta}{\delta \beta} y_i \ln(\lambda_i) \right) \\ &=& \sum_{i=1}^n \left( - \frac{y_i + \alpha^{-1}}{\alpha^{-1}+\lambda_i} + \frac{y_i}{\lambda_i} \right) \lambda_i \mathbf{X}_{i} \\ \end{eqnarray*}\]
En reparamétrisant, nous obtenons la condition de premier ordre à satisfaire pour obtenir \(\hat{\beta}_{MLE}\):
\[\begin{eqnarray*} \sum_{i=1}^n \Big(\frac{ y_i - \lambda_i}{1 + \lambda_i \alpha}\Big) \mathbf{X}_{i} = 0 \end{eqnarray*}\]
4.3.1.3 MLE - \(\alpha\)
La fonction de probabilité de la NB2 contient aussi le paramètre \(\alpha\), que nous devons aussi estimer par maximum de vraisemblance. Nous avons ainsi:
\[\begin{eqnarray*} \frac{\delta \ell}{\delta \alpha} &=& \frac{\delta}{\delta \alpha} \sum_{i=1}^n \bigg[\Big(\sum_{j=0}^{ y_i -1}\ln(j + \alpha^{-1})\Big) -\ln(z_i!) - ( y_i + \alpha^{-1}) \ln(\alpha^{-1}+\lambda_i) + y_i \ln(\lambda_i) +\alpha^{-1}\ln(\alpha^{-1}) \bigg]\\ &=& \sum_{i=1}^n \frac{\delta}{\delta \alpha} \Big(\sum_{j=0}^{ y_i -1}\ln(j + \alpha^{-1})\Big) - \frac{\delta}{\delta \alpha} ( y_i + \alpha^{-1}) \ln(\alpha^{-1}+\lambda_i) + \frac{\delta}{\delta \alpha} \alpha^{-1}\ln(\alpha^{-1}) \\ &=& \frac{-1}{\alpha^2} \sum_{i=1}^n \Big(\sum_{j=0}^{ y_i -1}\frac{1}{j + \alpha^{-1}} \Big) - \frac{\delta}{\delta \alpha} ( y_i + \alpha^{-1}) \ln(\alpha^{-1}+\lambda_i) + \frac{\delta}{\delta \alpha} \alpha^{-1}\ln(\alpha^{-1}) \end{eqnarray*}\]
Après quelques manipulations mathématiques, nous trouvons la condition à satisfaire pour obtenir \(\hat{\alpha}_{MLE}\):
\[\begin{eqnarray*} \sum_{i=1}^n \left[ \frac{1}{\alpha^2} \left( \ln(1+ \alpha \lambda_i) - \sum_{j=0}^{ y_i - 1} \frac{1}{j + \alpha^{-1}} \right) + \Big(\frac{ y_i - \lambda_i}{\alpha(1 + \lambda_i \alpha)}\Big) \right] = 0 \end{eqnarray*}\]
Ainsi, les estimateurs du maximum de vraisemblance pour la NB2, \(\hat{\mathbf{\beta}},\hat{a}\), sont obtenus lorsque les deux conditions de premier ordre sont respectées simultanément.
4.3.1.4 Variance des estimateurs
Puisque l’estimateur du maximum de vraisemblance est asymptotiquement sans biais et atteint la borne inférieure de Cramer-Rao pour la variance, nous avons la distribution asymptotique suivante :
\[ \left[ \begin{array}{c} \hat{\mathbf{\beta}} \\ \hat{\alpha} \end{array} \right] \stackrel{.}{\sim} N \left(\left[ \begin{array}{c} \mathbf{\beta} \\ \alpha \end{array} \right], \left[ \begin{array}{cc} Var_{MLE}[\hat{\mathbf{\beta}}] & Cov_{MLE}[\hat{\mathbf{\beta}}, \hat{\alpha}] \\ Cov_{MLE}[\hat{\mathbf{\beta}}, \hat{\alpha}] & Var_{MLE}[\hat{\alpha}] \\ \end{array} \right]\right) \]
où la matrice variance-covariance a été trouvée en trouvant la deuxième dérivée du logvraisemblance par rapport à \(\beta\) pour avoir:
\[\begin{eqnarray*} Var_{MLE}[\hat{\beta}] &=& \bigg[ \sum_{i=1}^n \mathbf{X}_i \mathbf{X}_i' \frac{\lambda_i} {1 + \lambda_i \alpha} \bigg]^{-1}, \label{Varbetamle} \end{eqnarray*}\]
et la deuxième dérivée du logvraisemblance par rapport à \(\alpha\) pour avoir:
\[\begin{eqnarray*} Var_{MLE}[\hat{\alpha}] &=& \bigg[ \sum_{i=1}^n \frac{1}{\alpha^4} \left( \ln(1+\alpha \lambda_i) - \Big[\sum_{j=0}^{y_i -1} \frac{1}{\alpha^{-1}} \Big] \right)^2 - \frac{\lambda_i}{\alpha^2(1+\alpha \lambda_i)}\bigg]^{-1}. \end{eqnarray*}\]
La dérivée du logvraisemblance par rapport à \(\beta\), et ensuite par rapport à \(\alpha\) (ou vice-versa) génère:
\[\begin{equation*} Cov_{MLE}[\hat{\mathbf{\beta}}, \hat{\alpha}] = 0, \end{equation*}\]
car l’espérance de la dérivée de \(\frac{\delta \ell}{\delta \beta}\) par rapport à \(\alpha\) est égal à:
\[\begin{equation*} E\bigg[\frac{\delta^2 \ell}{\delta \mathbf{\beta} \delta \alpha}\bigg] = E\bigg[ - \sum_{i=1}^n \frac{Y_i - \lambda_i}{(1+\lambda_i \alpha)^2} \lambda_i \mathbf{X}_i \mathbf{X}_i' \bigg] = 0 \end{equation*}\]
puisque l’espérance de \(Y_i\) est égale à \(\lambda_i\).
Ceci fait en sorte qu’il devient beaucoup plus simple de calculer la matrice variance-covariance car:
\[ \left[ \begin{array}{cc} E\left[ \frac{\delta^2 \ell}{\delta \beta \delta \beta'}\right] & 0 \\ 0 & E\left[ \frac{\delta^2 \ell}{\delta \alpha^2}\right] \\ \end{array} \right]^{-1} = \left[ \begin{array}{cc} \left[E\left[ \frac{\delta^2 \ell}{\delta \beta \delta \beta'}\right] \right]^{-1} & 0 \\ 0 & \left[E\left[ \frac{\delta^2 \ell}{\delta \alpha^2}\right] \right]^{-1}\\ \end{array} \right] \]
4.3.2 Binomiale négative 1
La seconde option est la NB1. Tel que nous l’avons vu, la fonction de probabilité de la NB1 peut s’exprimer comme:
\[\begin{eqnarray*} \Pr(Y_i = y_i) = \frac{\Gamma(\alpha^{-1} \lambda_i + y_i)}{\Gamma(\alpha^{-1} \lambda_i )\Gamma(y_i+1)} \left(\frac{\alpha^{-1}}{\alpha^{-1} + 1} \right)^{\alpha^{-1} \lambda_i} \left(\frac{1}{\alpha^{-1} + 1} \right)^{y_i} \end{eqnarray*}\]
Tout comme pour la NB2, l’utilisation de \(\alpha^{-1}\) dans la fonction de probabilité (au lieu de \(\alpha\), tel que nous l’avions vu) permet de mieux voir que le modèle converge vers une Poisson pour \(\alpha \rightarrow 0\). Cette reparamétrisation nous sera utile plus tard.
Des développements similaires à ce que nous avons fait avec la NB2 peuvent être faits avec la NB1. C’est d’ailleurs demandé pour le devoir #1.
En procédant exactement de la même manière, vous pourrez trouver les conditions de premier ordre pour obtenir les paramètres MLE de la NB1. En dérivant deux fois le logvraisemblance, un estimateur de la variance pourra être obtenu. Attention, contrairement à la NB2, les MLE de \(\hat{\mathbf{\beta}}\) et \(\hat{\alpha}\) ne sont pas indépendants et la covariance entre les estimateurs n’est pas égale à 0.