Appariement statistique par imputation fractionnaire 6. Études par simulationAppariement statistique par imputation fractionnaire 6. Études par simulation

• Précédent
• Suivant

Pour mettre à l’essai notre théorie, nous présentons deux études par simulation limitées. La première porte sur la combinaison de deux enquêtes indépendantes avec observation partielle afin d’effectuer une analyse conjointe. La deuxième porte sur la définition de modèles d’erreur de mesure avec calage externe.

6.1 Première simulation

Pour comparer les méthodes proposées avec les méthodes actuelles, nous avons généré 5 000 échantillons Monte Carlo comprenant $\left(x_i\mathrm{,}{y}_{1i}\mathrm{,}{y}_{2i}\right)$ et de taille $n\mathrm{<mo>=</mo>\; 400},$ où

$$\left(\begin{array}{c}{y}_{1i}\\ x_i\end{array}\right)\sim N\left(\left[\begin{array}{l}2\hfill \\ 3\hfill \end{array}\right]\mathrm{,}\left[\begin{array}{cc}1& 0\text{,}7\\ 0\text{,}7& 1\end{array}\right]\right)\mathrm{,}$$

$$y_{2i}\mathrm{<mo>=</mo>}{\beta }_0+{\beta }_1{y}_{1i}+e_i\mathrm{,}(6.1)$$

$e_i\sim N\left(\mathrm{0,}{\sigma }^2\right),$ et $\beta \mathrm{<mo>=</mo>}({\beta }_0\mathrm{,}{\beta }_1\mathrm{,}{\sigma }^2{)}^{\prime }\mathrm{<mo>=</mo>}(\mathrm{1,1,1}{)}^{\prime }.$ Soulignons que dans ce scénario, on trouve $f\left(y_2|x\mathrm{,}{y}_1\right)\mathrm{<mo>=</mo>}f\left(y_2|y_1\right)$ ; la variable $x$ joue donc le rôle de variable instrumentale pour $y_1.$

Au lieu d’observer $\left(x_i\mathrm{,}{y}_{1i}\mathrm{,}{y}_{2i}\right)$ conjointement, on présume que seuls $\left(y_1\mathrm{,}x\right)$ sont observés dans l’échantillon A, et que seuls $\left(y_2\mathrm{,}x\right)$ sont observés dans l’échantillon B; l’échantillon A est obtenu par sélection des $n_a\mathrm{<mo>=</mo>\; 400}$ premiers éléments de l’échantillon initial, et l’échantillon B, par sélection des $n_b\mathrm{<mo>=</mo>\; 400}$ éléments qui restent. On veut estimer quatre paramètres : trois paramètres de régression ${\beta }_0\mathrm{,}{\beta }_1\mathrm{,}{\sigma }^2$ et $\pi \mathrm{<mo>=</mo>}P\left(y_1\mathrm{<mo><</mo>\; 2,}{y}_2\mathrm{<mo><</mo>\; 3}\right),$ la proportion de $y_1\mathrm{<mo><</mo>\; 2}$ et de $y_2\mathrm{<mo><</mo>\; 3.}$ Quatre méthodes sont envisagées pour estimer ces paramètres :

1. Estimation de l’échantillon complet (EEC) : Utiliser toutes les observations de $\left(y_{1i}\mathrm{,}{y}_{2i}\right)$ de l’échantillon B.
2. Imputation par régression stochastique (IRS) : Utiliser la régression de $y_1$ sur $x$ dans l’échantillon A pour obtenir $\left({\widehat{\alpha }}_0\mathrm{,}{\widehat{\alpha }}_1\mathrm{,}{\widehat{\sigma }}_1^2\right),$ où le modèle de régression est $y_1\mathrm{<mo>=</mo>}{\alpha }_0+{\alpha }_{1}x+e_1$ avec $e_1\sim \left(\mathrm{0,}{\sigma }_1^2\right).$ Pour chaque $i\in B,$ $m\mathrm{<mo>=</mo>\; 10}$ valeurs imputées sont générées par $y_{1i}^{\mathrm{<mo>*</mo>}\left(j\right)}\mathrm{<mo>=</mo>}{\widehat{\alpha }}_0+{\widehat{\alpha }}_1{x}_i+e_i^{\mathrm{<mo>*</mo>}\left(j\right)}$ où $e_i^{\mathrm{<mo>*</mo>}\left(j\right)}\sim N\left(\mathrm{0,}{\widehat{\sigma }}_1^2\right).$
3. Imputation fractionnaire paramétrique (IFP) avec $m\mathrm{<mo>=</mo>\; 10}$ selon l’hypothèse de variable instrumentale.
4. Imputation fractionnaire hot deck (IFHD) avec $m\mathrm{<mo>=</mo>\; 10}$ selon l’hypothèse de variable instrumentale.

Le tableau 6.1 présente les moyennes et les variances Monte Carlo des estimateurs ponctuels des quatre paramètres d’intérêt. La méthode IRS comporte des biais importants pour tous les paramètres étudiés, parce qu’elle repose sur l’hypothèse d’indépendance conditionnelle. Les méthodes IFP et IFHD fournissent des estimateurs presque sans biais pour tous les paramètres. Les estimateurs obtenus à l’aide de la méthode IFHD sont légèrement plus efficients que ceux obtenus par la méthode IFP, parce que la démarche en deux étapes de la méthode IFHD utilise l’ensemble complet des répondants à la première étape. La variance asymptotique théorique de ${\widehat{\beta }}_1$ calculée à partir de la méthode IFP est

$$V\left({\widehat{\beta }}_1\right)\doteq \frac{1}{{\left(0\text{,}7\right)}^2}\frac{1}{400}2\left(1-\frac{0\text{,}{7}^2}{2}\right)+\frac{1}{{\left(0\text{,}7\right)}^2}\frac{1}{400}\left(1-0\text{,}{7}^2\right)\doteq 0\text{,}0103$$

ce qui correspond au résultat de simulation présenté dans le tableau 6.1. En plus de l’estimation ponctuelle, on calcule aussi les estimateurs de variance pour les méthodes IFP et IFHD. Les estimateurs de variance montrent de faibles biais relatifs (moins de 5 % en valeur absolue) pour tous les paramètres. Les résultats de l’estimation de la variance ne sont pas présentés ici par souci de concision.

La méthode proposée repose sur l’hypothèse de variable instrumentale. Pour déterminer la sensibilité de la méthode d’imputation fractionnaire proposée aux violations de l’hypothèse de variable instrumentale, nous avons réalisé une étude par simulation supplémentaire. Au lieu de générer $y_{2i}$ à partir de (6.1), on utilise

$$y_{2i}\mathrm{<mo>=</mo>\; 0}\text{,}5+y_{1i}+\rho \left(x_i-3\right)+e_i\mathrm{,}(6.2)$$

où $e_i\sim N\left(\mathrm{0,1}\right)$ et $\rho$ peuvent prendre des valeurs non nulles. Nous avons utilisé trois valeurs de $\rho ,$ $\rho \in \left\{0;0\text{,}1;0\text{,}2\right\},$ pour l’analyse de sensibilité et nous avons employé la même procédure d’IFP et d’IFHD fondée sur l’hypothèse que $x$ est une variable instrumentale pour $y_1.$ Cette hypothèse est satisfaite pour $\rho \mathrm{<mo>=</mo>\; 0},$ mais elle est légèrement enfreinte pour $\rho \mathrm{<mo>=</mo>\; 0}\text{,}1$ ou $\rho \mathrm{<mo>=</mo>\; 0}\text{,}2.$ À l’aide des données obtenues par imputation fractionnaire dans l’échantillon B, nous avons estimé trois paramètres, ${\theta }_1\mathrm{<mo>=</mo>}E\left(Y_1\right),$ ${\theta }_2$ la pente de la régression simple de $y_2$ sur $y_1,$ et ${\theta }_3\mathrm{<mo>=</mo>}P\left(y_1\mathrm{<mo><</mo>\; 2,}{y}_2\mathrm{<mo><</mo>\; 3}\right),$ la proportion de $y_1\mathrm{<mo><</mo>\; 2}$ et $y_2\mathrm{<mo><</mo>\; 3.}$ Le tableau 6.2 présente les moyennes et les variances Monte Carlo des estimateurs ponctuels pour les trois paramètres en vertu des trois différents modèles. Dans le tableau 6.2, on constate que les valeurs absolues de l’écart entre l’estimateur obtenu par imputation fractionnaire et l’estimateur obtenu à partir de l’échantillon complet augmente avec la valeur de $\rho ,$ ce qui est normal puisque l’hypothèse de variable instrumentale est plus gravement enfreinte pour les valeurs élevées de $\rho ,$ mais les écarts sont relativement faibles dans tous les cas. Plus particulièrement, l’estimateur de ${\theta }_1$ n’est pas touché par la dérogation à l’hypothèse de variable instrumentale, parce que l’estimateur par imputation obtenu en vertu du modèle d’imputation inexact fournit quand même un estimateur non biaisé pour la moyenne de population, à condition que le modèle d’imputation par régression contienne un terme d’ordonnée à l’origine (Kim et Rao 2012). Ainsi, cette analyse de sensibilité limitée indique que la méthode proposée semble produire des estimations comparables lorsque l’hypothèse de variable instrumentale est légèrement enfreinte.

6.2  Deuxième simulation

Dans la deuxième étude par simulation, on examine une variable de réponse binaire $y_{2i},$ où

$$y_{2i}\sim \text{Bernoulli}\left(p_i\right)\mathrm{,}(6.3)$$

$$\text{logit}\left(p_i\right)\mathrm{<mo>=</mo>}{\gamma }_0+{\gamma }_1{y}_{1i}\mathrm{,}$$

et $y_{1i}\sim N\left({\mu }_1\mathrm{,}{\sigma }_1^2\right).$ Dans l’échantillon principal, désigné par la lettre $B,$ au lieu d’observer $\left(y_{1i}\mathrm{,}{y}_{2i}\right),$ on observe $\left(x_i\mathrm{,}{y}_{2i}\right),$ où

$$x_i\mathrm{<mo>=</mo>}{\beta }_0+{\beta }_1{y}_{1i}+u_i\mathrm{,}(6.4)$$

et $u_i\sim N\left(\mathrm{0,}{\sigma }^2{\left|y_{1i}\right|}^{2\alpha }\right).$ On observe $\left(x_i\mathrm{,}{y}_{1i}\right),$ $i\mathrm{<mo>=</mo>\; 1,}\dots \mathrm{,}{n}_A$ dans un échantillon de calage, désigné par la lettre A. Aux fins de la simulation, $n_A\mathrm{<mo>=</mo>}{n}_B\mathrm{<mo>=</mo>\; 800},$ ${\gamma }_0\mathrm{<mo>=</mo>\; 1},$ ${\gamma }_1\mathrm{<mo>=</mo>\; 1},$ ${\beta }_0\mathrm{<mo>=</mo>\; 0},$ ${\beta }_1\mathrm{<mo>=</mo>\; 1},$ ${\sigma }^2\mathrm{<mo>=</mo>\; 0}\text{,}25,$ $\alpha \mathrm{<mo>=</mo>\; 0}\text{,}4,$ ${\mu }_1\mathrm{<mo>=</mo>\; 0},$ et ${\sigma }_1^2\mathrm{<mo>=</mo>\; 1.}$ Le but principal est d’estimer ${\gamma }_1$ et de mettre à l’essai l’hypothèse nulle que ${\gamma }_1\mathrm{<mo>=</mo>\; 1.}$ La taille de l’échantillon Monte Carlo (MC) est 1 000.

On compare les estimateurs IFP de ${\gamma }_1$ aux trois autres estimateurs. Comme la distribution conditionnelle de $y_{1i}$ sachant $x_i$ n’est pas standard, on utilise les poids de (5.3) pour réaliser l’IFP, où la distribution proposée ${\widehat{h}}_a\left(y_{1i}|x_i\right)$ est une estimation de la distribution marginale de $y_{1i}$ fondée sur les données de l’échantillon A. On considère les quatre estimateurs suivants :

1. IFP : Pour l’IFP, la distribution proposée pour générer $y_{1i}^{\mathrm{<mo>*</mo>}\left(j\right)}$ est une distribution normale de moyenne ${\widehat{\mu }}_1$ et de variance ${\widehat{\sigma }}_1^2,$ où ${\widehat{\mu }}_1$ et ${\widehat{\sigma }}_1^2$ sont les estimations du maximum de vraisemblance de ${\mu }_1$ et ${\sigma }_1^2,$ respectivement, en fonction de l’échantillon A. Les poids fractionnaires définis en (5.3) se présentent sous la forme

$$w_{ij}^{*}\propto {\widehat{p}}_{ij}^{y_{2i}}{\left(1-{\widehat{p}}_{ij}\right)}^{1-y_{2i}}{\widehat{f}}_a\left(y_{1i}^{\mathrm{<mo>*</mo>}\left(j\right)}|x_i\right)\mathrm{,}(6.5)$$

• où ${\widehat{p}}_{ij}\mathrm{<mo>=</mo>}{\left\{1+\exp \left(-{\widehat{\gamma }}_0-{\widehat{\gamma }}_1{y}_{1i}^{\mathrm{<mo>*</mo>}\left(j\right)}\right)\right\}}^{-1}$ et ${\widehat{f}}_a\left(y_{1i}|x_i\right)$ est l’estimation de $f\left(y_{1i}|x_i\right)$ fondée sur l’estimation du maximum de vraisemblance à partir des données de l’échantillon A. La taille de la classe d’imputation est $m\mathrm{<mo>=</mo>\; 800.}$

2. Estimateur naïf : Un estimateur naïf est l’estimateur de la pente dans la régression logistique de $y_{2i}$ sur $x_i$ pour $i\in B.$
3. Estimateur bayésien : On utilise l’approche de Guo et Little (2011) pour définir un estimateur bayésien. Le modèle de notre simulation diffère de celui de Guo et Little (2011) par le fait que la réponse d’intérêt est binaire. On établit un échantillonnage de Gibbs à l’aide du programme JAGS (Plummer 2003), en précisant les distributions a priori diffuses appropriées pour les paramètres du modèle. Soit

$${\theta }_1\mathrm{<mo>=</mo>}\left(\text{log}\left({\sigma }^2\right)\mathrm{,}\text{log}\left({\sigma }_1^2\right)\mathrm{,}{\mu }_1\mathrm{,}{\beta }_0\mathrm{,}{\beta }_1\mathrm{,}{\gamma }_0\mathrm{,}{\gamma }_1\right);$$

• on suppose a priori que ${\theta }_1\sim N\left({\mathrm{0,10}}^6{I}_7\right),$ où $I_7$ est une matrice identité de dimensions  $7\times 7$ et où la notation $N\left(\mathrm{0,}V\right)$ désigne une distribution normale de moyenne 0 et de matrice de covariances  $V.$ La distribution a priori pour la puissance $\alpha$ est uniforme sur l’intervalle  $\left[-\mathrm{5,5}\right].$
• Pour évaluer la convergence, on examine les tracés des courbes et les facteurs de réduction d’échelle possibles définis par Gelman, Carlin, Stern et Rubin (2003) pour 10 ensembles de données simulées préliminaires. On produit trois chaînes MCMC, chacune d’une longueur de 1 500, à partir de valeurs initiales aléatoires, et on rejette les 500 premières itérations, considérées comme faisant partie du rodage. Les facteurs de réduction d’échelle possibles dans les 10 ensembles de données simulées vont de 1,0 à 1,1, et les tracés des courbes indiquent que les chaînes se combinent bien. Pour réduire le temps de calcul, on utilise 1 000 itérations d’une seule chaîne pour la simulation principale, après rejet des 500 premières itérations de rodage.

4. Un estimateur par calage par régression pondérée (CRP). Cet estimateur par CRP est une modification de l’estimateur par calage par régression pondérée défini par Guo et Little (2011) pour une variable de réponse binaire. On calcule l’estimateur par calage par régression pondérée comme suit :

• (i)   À l’aide des moindres carrés ordinaires (MCO), effectuer une régression de $y_{1i}$ sur $x_i$ pour l’échantillon de calage.
• (ii)  Effectuer une régression du logarithme des résidus quadratiques obtenus à l’étape (i) sur le logarithme de $x_i^2$ pour l’échantillon de calage. Soit $\widehat{\lambda }$ la pente estimée de la régression.
• (iii) À l’aide des moindres carrés pondérés (MCP) avec le poids ${\left|x_i\right|}^{2\widehat{\lambda }}\text{},$ effectuer la régression de $y_{1i}$ sur $x_i$ pour l’échantillon de calage. Soient ${\widehat{\eta }}_0$ et ${\widehat{\eta }}_1$ l’ordonnée à l’origine et la pente estimées, respectivement, de la régression des MCP.
• (iv) Pour chaque unité $i$ de l’échantillon principal, posons ${\widehat{y}}_{1i}\mathrm{<mo>=</mo>}{\widehat{\eta }}_0+{\widehat{\eta }}_1{x}_i.$
• (v)  L’estimation de $\left({\gamma }_0\mathrm{,}{\gamma }_1\right)$ est obtenue à partir de la régression logistique de $y_{2i}$ sur ${\widehat{y}}_{1i}$ dans l’échantillon principal.

Le tableau 6.3 indique le biais, la variance et l’EQM Monte Carlo des quatre estimateurs de ${\gamma }_1.$ L’estimateur naïf a un biais négatif parce que $x_i$ est une version contaminée de $y_{1i}.$ L’estimateur par IFP est supérieur à l’estimateur bayésien et à l’estimateur par CRP.

On calcule une estimation de la variance des estimateurs par IFP de ${\gamma }_1$ à l’aide de l’expression de la variance fondée sur l’approximation linéaire. On définit le biais relatif MC comme étant le ratio de la différence entre la moyenne MC de l’estimateur de variance et la variance MC de l’estimateur à la variance MC de l’estimateur. Le biais relatif MC des estimateurs de variance pour l’IFP est négligeable (moins de 2 % en valeur absolue).

• Précédent
• Suivant