6. Études par simulation

Alina Matei et M. Giovanna Ranalli

Précédent | Suivant

Nous évaluons la performance de l’estimateur présenté à la section 5 au moyen d’une simulation Monte Carlo sous deux scénarios différents. Le premier utilise un ensemble de données réelles comme population et des variables d’intérêt qui sont toutes binaires, tandis que le second utilise des données de population simulées avec des variables d’intérêt continues. Les résultats pour le premier scénario sont présentés à la section 6.1, tandis que ceux pour le deuxième sont présentés à la section 6.2.

Sous les deux scénarios, nous utilisons l’échantillonnage aléatoire simple sans remise et considérons les estimateurs suivants :

• $\text{HT=}{\displaystyle {\sum }_{k\in s}{y}_{kj}/{\pi }_k}:$  l’estimateur de Horvitz-Thompson dans le cas d’une réponse complète est calculé comme valeur de référence en l’absence de non-réponse.

• ${\widehat{Y}}_{j,\text{naïf}}:$  l’estimateur naïf donné en (5.1); aucune mesure explicite n’est prise pour corriger les non-réponses totale et partielle. Notons que, sous échantillonnage aléatoire simple sans remise, il se réduit à ${\widehat{Y}}_{j,\text{naïf}}=N{\displaystyle {\sum }_{k\in r_j}{y}_{kj}/n_{r_j}},$  où $n_{r_j}$  est la taille de l’ensemble $r_j\mathrm{,}$  et il est identique à l’estimateur de Horvitz-Thompson ajusté pour la non-réponse totale sous l’hypothèse de probabilités de réponse uniformes estimées par $n_{r_j}/n\mathrm{.}$

• ${\widehat{Y}}_{j,pq}:$  l’estimateur sous échantillonnage à trois phases proposé à la section 5, équation (5.2).

• ${\widehat{Y}}_{j,pq,\text{réelles}}:$  l’estimateur sous échantillonnage à trois phases qui utilise les valeurs réelles des probabilités de réponse $p_k$  et $q_{kj}$  est également calculé aux fins de comparaison avec ${\widehat{Y}}_{j,pq}$  pour comprendre l’effet de l’estimation des probabilités de réponse.

Les simulations sont exécutées en R version 2.15, en utilisant le module R « ltm » (Rizopoulos 2006) pour ajuster les modèles à traits latents. Les mesures de performance suivantes sont calculées pour chaque estimateur, ci-après désigné génériquement par $\widehat{Y}$  où le suffixe $j$  est abandonné pour simplifier la notation $(Y$  désigne le total de population):

• le biais Monte Carlo

  $$\text{B=}{E}_{\text{sim}}\left(\widehat{Y}\right)-Y\mathrm{,}$$

  où $E_{\text{sim}}\left(\widehat{Y}\right)={\displaystyle {\sum }_{i=1}^M{\widehat{Y}}_i}/M,{\widehat{Y}}_i$  est la valeur de l’estimateur $\widehat{Y}$  à la $i^{\text{e}}$  exécution de la simulation et $M$  est le nombre total d’exécutions de la simulation;

• le biais relatif

  $$\text{BR=}\frac{\text{B}}{Y};$$

• l’écart-type Monte Carlo

  $$\sqrt{\text{VAR}}=\sqrt{\frac{1}{M-1}{\displaystyle \sum _{i=1}^M}{\left({\widehat{Y}}_i-E_{\text{sim}}\left(\widehat{Y}\right)\right)}^2};$$

• l’erreur quadratique moyenne Monte Carlo

  $${\text{EQM=B}}^2+\text{VAR}\mathrm{.}$$

6.1 Scénario de simulation 1

Nous considérons un jeu de données formé de quatre variables binaires extraites de la British Social Attitudes Survey de 1986 et qui concernent l’attitude à l’égard de l’avortement. Les données sont disponibles dans le module R « ltm » (Rizopoulos 2006). $N=\text{379}$  personnes ont répondu aux questions suivantes après qu’on leur ait demandé si la loi devrait permettre l’avortement dans les circonstances présentées par chaque item :

1. La femme décide toute seule qu’elle ne souhaite pas garder le bébé.
2. Le couple est d’accord qu’il ne souhaite pas avoir un enfant.
3. La femme n’est pas mariée et ne souhaite pas épouser l’homme.
4. Le couple n’a pas les moyens d’avoir un autre enfant.

La variable d’intérêt $y_j$  est choisie comme étant la deuxième $\left(j=2\right)$  avec un total $Y_j\text{=225}$  dans la population.

Les données sont analysées par Bartholomew et coll. (2002) à titre d’exemple de situation où l’on peut trouver une variable latente qui mesure l’attitude à l’égard de l’avortement. Au niveau de la population, nous calculons la variable latente (notée ici ${\theta }_k^a)$  en utilisant le modèle (4.2) sur les données ${\left\{y_{k\ell }\right\}}_{k=1,\dots \mathrm{,}N;\ell =1,\dots \mathrm{,4}}.$  La corrélation entre les valeurs de $y_{k\ell }$  et ${\theta }_k^a$  est approximativement égale à 0,85, pour $\ell =1,\dots \mathrm{,4.}$  Ensuite, nous avons fixé ${\theta }_k={\widehat{\theta }}_k^a\mathrm{,}$  pour tout $k=1,\dots \mathrm{,}N\mathrm{.}$

Au niveau de la population, les probabilités de réponse totale sont générées en utilisant la structure de réponse suivante :

$$p_k=1/\left(1+\exp \left(-\left(\text{0,7}+y_{k2}+{\theta }_k+\text{0,2}{\epsilon }_k\right)\right)\right)\mathrm{,}(6.1)$$

avec ${\epsilon }_k\sim U\left(\mathrm{0,1}\right),$  pour simuler une non-réponse non ignorable. La moyenne de population de $p_k$  est égale à environ 0,74.

Pour générer les probabilités de réponse partielle au niveau de la population, le modèle utilisé est le suivant :

$$q_{k\ell }=1/\left(1+\exp \left(-\left(b_{\ell }{\theta }_k+a_{\ell }+y_{k\ell }\right)\right)\right)\mathrm{,}\text{pour }\ell =1,\dots \mathrm{,4,}(6.2)$$

où $b_{\ell }=3,$  pour $\ell =1,\dots \mathrm{,4},$  tandis que $a_{\ell }$  prend différentes valeurs en fonction de $\ell ;$  en particulier, $a_1=1,a_2=0,a_3=-0\text{,5}$  et $a_4=1.$  Le taux nominal de non-réponse partielle pour les quatre items dans la population est de 35 %, 42 %, 47 % et 31 %, respectivement.

Nous tirons $M=\text{10 000}$  échantillons aléatoires simples sans remise à partir de la population en utilisant deux tailles d’échantillon : $n=\text{50}$  et $n=\text{100}\mathrm{.}$  Dans chaque échantillon $s\mathrm{,}$  les unités sont classées comme étant des répondants conformément à un échantillonnage de Poisson, en utilisant les probabilités $p_k$  calculées comme dans l’équation (6.1) et résultant en l’ensemble $r.$  Alors, sachant $r,$  nous construisons la matrice ${\left\{x_{k\ell }\right\}}_{k\in r;\ell =1,\mathrm{...},4},$  où les valeurs de $x_{k\ell }$  sont tirées selon un échantillonnage de Poisson avec probabilités $q_{k\ell }$  définies en (6.2). Dans chaque simulation, le modèle (4.2) et l’ensemble de répondants $r$  sont utilisés pour calculer la variable ${\widehat{\theta }}_k$  pour tout $k\in s$  comme il est décrit à la section 4.4. Le modèle (4.4) est ajusté pour obtenir ${\widehat{p}}_k\mathrm{.}$  Le taux moyen de non-réponse partielle sur les simulations pour les quatre items est de 26 %, 33 %, 38 % et 23 %, respectivement. L’estimateur de variance jackknife a été calculé comme il est décrit à la section 5 en utilisant la fonction gencalib() dans le module R « sampling » (Tillé et Matei 2012) et la distance logistique (Deville, Särndal et Sautory 1993).

Le tableau 6.1 donne les résultats pour $n=\text{50}$  et $n=\text{100}.$  Comme prévu, les estimateurs $\text{HT}$  et ${\widehat{Y}}_{j,pq,\text{réelles}}$  ont un biais presque nul, tandis que le second présente une EQM relativement plus grande qui est due uniquement à la plus petite taille d’échantillon. L’estimateur naïf donne un biais négatif très important. Cela tient au fait que les unités dont la valeur de $y_j$  est nulle sont moins susceptibles de répondre et que le total est clairement sous-estimé. L’estimateur ${\widehat{Y}}_{j,pq}$  présente un bien plus petit biais que l’estimateur naïf. Notons que la performance de l’estimateur proposé est dictée principalement par le biais absolu, de sorte qu’elle ne diffère pas particulièrement lorsqu’on augmente la taille de l’échantillon, mis à part une diminution de la variance. Si nous comparons ${\widehat{Y}}_{j,pq,réell\text{es}}$  et ${\widehat{Y}}_{j,pq},$  nous notons que ${\widehat{Y}}_{j,pq}$  souffre encore d’un certain biais qui provient de la spécification incorrecte de la structure de réponse (nous ne tenons pas compte des valeurs des variables d’intérêt).

Pour l’estimateur proposé, l’estimateur de variance jackknife a également été testé en examinant la couverture empirique d’un intervalle de confiance à 95 % calculé pour chaque réplique comme ${\widehat{Y}}_{j,pq}\pm \text{1,96}\sqrt{{\widehat{V}}_r}.$  Pour $n=\text{50}\mathrm{,}$  la valeur moyenne de $\sqrt{{\widehat{V}}_r}$  sur l’ensemble des simulations était de 54,8, tandis que pour $n=\text{100}\mathrm{,}$  elle était de 53,3, avec un taux de couverture de l’IC à 95 % de 94,6 % et de 96,3 %, respectivement. L’estimateur par rééchantillonnage surestime l’écart-type Monte Carlo donné pour ${\widehat{Y}}_{j,pq}$  au tableau 6.1 dans les deux cas, mais possède de bons taux de couverture.

Tableau 6.1
Résultats des simulations sous le scénario 1 – Ensemble de données sur l’avortement
Sommaire du tableau
Le tableau montre les résultats de Résultats des simulations sous le scénario 1 – Ensemble de données sur l’avortement. Les données sont présentées selon Estimateur (titres de rangée) et B, $\sqrt{\text{VAR}}$, EQM et BR %(figurant comme en-tête de colonne).

Estimateur	B	$\sqrt{\text{VAR}}$	EQM	BR %
$n=\text{50}$
$\text{HT}$	0,05	24,5	600,5	$<$ 0,1
${\widehat{Y}}_{j,\text{naïf}}$	-126,5	19,4	16 378,6	-56,2
${\widehat{Y}}_{j,pq}$	20,6	32,4	1 474,1	9,1
${\widehat{Y}}_{j,pq,\text{réelles}}$	0,02	35,0	1 225,0	$<$ 0,1
$n=\text{100}$
$\text{HT}$	-0,06	16,0	255,5	$<$ 0,1
${\widehat{Y}}_{j,\text{naïf}}$	-126,9	13,5	16 284,1	-56,4
${\widehat{Y}}_{j,pq}$	17,9	21,9	802,2	8,0
${\widehat{Y}}_{j,pq,\text{réelles}}$	-0,1	23,7	559,9	$<$ 0,1

Pour étudier la performance du modèle à traits latents au niveau de la population et la corrélation entre la variable d’intérêt et la variable latente estimée, nous avons appliqué la procédure décrite plus haut en utilisant $q_{k\ell }$  définie en (6.2) pour construire la matrice ${\left\{x_{k\ell }\right\}}_{k=1,\dots \mathrm{,}N;\ell =1,\mathrm{...},4}$  pour toutes les unités de la population. Nous avons ajusté le modèle (4.2) au niveau de la population et calculé la variable ${\theta }_k$  pour tout $k=1,\dots \mathrm{,}N.$  Le coefficient alpha de Cronbach prend la valeur de 0,83 ce qui indique une bonne cohérence interne des items. Le coefficient de corrélation entre la variable d’intérêt et la variable latente estimée prend la valeur de 0,76, ce qui indique que l’information auxiliaire latente possède un fort pouvoir de prédiction de $y_{k2}\mathrm{,}$  comme il l’a été prôné pour le modèle de Cassel et coll. (1983). L’inspection des marges de tableau à double entrée pour la matrice $\left\{x_{k\ell }\right\}$  donne les résidus ${\left(O-E\right)}^2/E,$  compris entre 0,03 et 0,23. De même, les marges de tableau à triple entrée pour la matrice $\left\{x_{k\ell }\right\}$  donnent des résidus entre 0 et 1,19. Cela indique que nous n’avons aucune raison de rejeter ici le modèle à un facteur latent (4.2) (voir Bartholomew et coll. 2002, page 186).

6.2  Scénario de simulation 2

Nous générons $\left\{y_{k1}\mathrm{,}\dots \mathrm{,}{y}_{k6}\mathrm{,}{\theta }_k\right\}$  pour $k=1,\dots \mathrm{,}N=\text{2 000}$  en utilisant une loi normale multivariée de moyenne 1. Le degré de corrélation entre $y_{\ell }$  et $y_{{\ell }^{\prime }}$  est 0,8, avec $\ell \mathrm{,}{\ell }^{\prime }=1,\dots \mathrm{,6,}\ell \ne {\ell }^{\prime }.$  Nous posons que la variable d’intérêt est $y_6$  et considérons divers degrés de corrélation entre ses valeurs et celles prises par ${\theta }_k,$  à savoir 0,3, 0,5, 0,8. Les valeurs de ${\theta }_k$  sont ensuite centrées et réduites afin qu’elles soient de moyenne 0 et de variance 1.

Les probabilités de réponse sont obtenues en calculant d’abord

$$p_k^{\circ }=1/\left[1+\exp \left(-\left(0,5+y_{k1}+{\theta }_k\right)\right)\right]\mathrm{,}\text{pour }k=1,\dots \mathrm{,}N\mathrm{,}(6.3)$$

puis en les rééchelonnant afin qu’elles prennent des valeurs comprises entre 0,1 et 0,9, avec une moyenne de population approximativement égale à 0,7.

Les probabilités de réponse partielle sont générées en calculant d’abord :

$$q_{k\ell }^{\circ }=1/\left(1+\exp \left(-\left(b_{\ell }{\theta }_k+a_{\ell }+y_{k\ell }\right)\right)\right)\mathrm{,}\text{pour }k=1,\dots \mathrm{,}N\text{ et }\ell =1,\dots \mathrm{,6,}(6.4)$$

où ${\left\{a_{\ell }\right\}}_{\ell =1,\dots \mathrm{,6}}=\left\{1;0;\text{-0,5;}1;0;\text{-0,5}\right\}$  et ${\left\{b_{\ell }\right\}}_{\ell =1,\dots \mathrm{,6}}=\left\{1;1;1;\text{1,5;}\text{1,5;}\text{1,5}\right\}\mathrm{,}$  puis en rééchelonnant les valeurs pour qu’elles soient comprises entre 0,1 et 0,95.

Nous tirons $M=\text{10 000}$  échantillons par échantillonnage aléatoire simple sans remise de taille $n=\text{200}.$  Pour chaque échantillon $s,$  un ensemble de réponses $r$  est créé en réalisant un échantillonnage de Poisson de paramètre $p_k$  défini en (6.3). Chaque élément de la matrice ${\left\{x_{k\ell }\right\}}_{k\in r\mathrm{,}\ell =1,\dots \mathrm{,6}}$  est généré en utilisant l’échantillonnage de Poisson de paramètre $q_{k\ell }$  défini en (6.4). Les taux de non-réponse partielle sur l’ensemble des simulations prennent approximativement les valeurs de 18 %, 28 %, 35 %, 19 %, 29 % et 34 %, pour $\ell =1,\dots \mathrm{,6},$  respectivement. Pour chaque simulation, le modèle (4.2) est utilisé pour calculer la variable ${\widehat{\theta }}_k$  pour tout $k\in s.$  Le modèle (4.4) est alors ajusté pour obtenir ${\widehat{p}}_k\mathrm{.}$

Tableau 6.2
Résultats des simulations sous le scénario 2 – Données continues simulées
Sommaire du tableau
Le tableau montre les résultats de Résultats des simulations sous le scénario 2 – Données continues simulées. Les données sont présentées selon Estimateur (titres de rangée) et B, $\sqrt{\text{VAR}}$ , EQM et BR % (figurant comme en-tête de colonne).

Estimateur	B	$\sqrt{\text{VAR}}$	EQM	BR %
Coefficient de corrélation de 0,3
$\text{HT}$	-0,7	131,6	17 331,2	$\approx$ -0,0
${\widehat{Y}}_{j\mathrm{,}\text{naïf}}$	825,6	177,1	713 039,3	41,0
${\widehat{Y}}_{j\mathrm{,}pq}$	-227,4	188,0	87 033,0	-11,3
${\widehat{Y}}_{j,pq,\text{réelles}}$	48,4	231,8	56 073,2	2,4
Coefficient de corrélation de 0,5
$\text{HT}$	0,1	135,0	18 220,5	$\approx$ 0,0
${\widehat{Y}}_{j\mathrm{,}\text{naïf}}$	972,6	176,2	977 009,5	50,7
${\widehat{Y}}_{j\mathrm{,}pq}$	-180,0	175,5	63 552,0	-9,4
${\widehat{Y}}_{j,pq,\text{réelles}}$	74,8	212,7	50 844,0	3,9
Coefficient de corrélation de 0,8
$\text{HT}$	-0,1	134,1	17 992,0	$\approx$ -0,0
${\widehat{Y}}_{j\mathrm{,}\text{naïf}}$	1 154,6	168,1	1 361 388,1	57,7
${\widehat{Y}}_{j\mathrm{,}pq}$	-184,8	164,4	61 173,0	-9,2
${\widehat{Y}}_{j,pq,\text{réelles}}$	100,6	196,2	48 597,9	5,0

Le tableau 6.2 donne la performance des estimateurs pour les trois valeurs prises par le coefficient de corrélation nominal entre $y_{k1}$  et ${\theta }_k\text{:0,3},\text{0,5}$ et $\text{0,8}\text{.}$ L’estimateur proposé est toujours capable de réduire le biais comparativement à l’estimateur naïf, même quand la corrélation entre la variable d’intérêt et la variable latente devient plus faible. Le biais relatif prend des valeurs acceptables dans la plupart des cas. Le biais mérite d’être examiné de plus près. Dans tous les cas, l’estimateur naïf surestime fortement le total. Cela n’est pas étonnant, parce que les valeurs de $p_k,q_{k6},{\theta }_k$  et $y_{k6}$  vont toutes dans la même direction. Par conséquent, dans notre échantillon de répondants, nous sommes plus susceptibles de trouver des valeurs relativement grandes de $y_6,$  ce qui donne lieu à une surestimation pour l’estimateur naïf. Par ailleurs, ${\widehat{Y}}_{j,pq}$  sous-estime le total, parce qu’il est fondé uniquement sur les unités observées de $r_j$  qui ont des valeurs relativement grandes pour $y_6,$  mais aussi des valeurs relativement grandes pour $p_k$  et $q_{k6}$  et, par conséquent, qui ont à la fin un faible poids.

La matrice des valeurs de population ${\left\{x_{k\ell }\right\}}_{k=1,\dots \mathrm{,}\text{2 000}\mathrm{,}\ell =1,\dots \mathrm{,6}}$  est construite de la même façon qu’à la section 6.1 pour valider les hypothèses qui sous-tendent le modèle 2PL. Le coefficient alpha de Cronbach prend approximativement la valeur de 0,5 pour le coefficient de corrélation égal à 0,3, de 0,6 pour le coefficient de corrélation égal à 0,5, et de 0,7 pour le coefficient de corrélation égal à 0,8; l’association par paire entre les six items révèle des valeurs $p$ plus petites que 0,01. L’inspection des marges à double entrée et à triple entrée de la matrice $\left\{x_{k\ell }\right\}$  donne des résidus ${\left(O-E\right)}^2/E$  qui prennent tous des valeurs inférieures à 4. Par conséquent, le modèle à un facteur latent peut être accepté et les items semblent tous mesurer le même facteur latent.

Précédent | Suivant