Publié le :
Statistique CREATION_INTERNE

Simulation de variables Bernoulli corrélées et analyse Monte Carlo

Ce code est également disponible en : Deutsch English Español
En attente de validation
Le programme génère d'abord un jeu de données de 20 000 clusters contenant chacun 5 unités élémentaires avec une probabilité de succès de 0.6 et une corrélation intra-cluster de 0.3. Il transpose les données pour vérifier la structure de corrélation. Ensuite, il effectue des agrégations par cluster et analyse la variance. Enfin, une simulation de Monte Carlo est réalisée pour estimer le taux d'erreur de type I sur des données sur-dispersées.
Analyse des données

Type : CREATION_INTERNE


Toutes les données sont générées algorithmiquement via des étapes DATA utilisant la fonction aléatoire 'uniform'.

1 Bloc de code
DATA STEP Data
Explication :
Génération de 20 000 observations simulées structurées en clusters (m=5) avec une corrélation intra-classe définie (rho2=0.3).
Copié !
1DATA correlated_bernoullis;
2 n = 20000; *--- Number of clusters;
3 m = 5; *--- Number of elemental units within the cluster;
4 pi = 0.6; *--- Probability of success of each elemental unit;
5 rho2 = 0.3; *--- Intra-cluster correlation;
6 seed = 16670;
7 rho = sqrt( rho2 );
8 DO subjid = 1 to n;
9 yy = 0; *--- Variable yy plays the role of y of eq. (1.7);
10 u = uniform( seed );
11 IF u < pi THEN yy = 1;
12 DO i=1 to m;
13 y = 0;
14 u = uniform( seed );
15 IF u < rho THEN y = yy;
16 ELSE DO;
17 uu = uniform( seed );
18 IF uu < pi THEN y = 1;
19 END;
20 OUTPUT;
21 END;
22 END;
23 keep subjid i y;
24RUN;
2 Bloc de code
PROC TRANSPOSE Data
Explication :
Pivot des données pour transformer les observations intra-cluster en colonnes distinctes (_1 à _5) afin de calculer les corrélations.
Copié !
1PROC TRANSPOSE DATA=correlated_bernoullis out=new;
2 BY subjid;
3 id i;
4 var y;
5RUN;
3 Bloc de code
PROC CORR
Explication :
Calcul de la matrice de corrélation entre les unités du cluster pour vérifier si la simulation respecte les paramètres définis.
Copié !
1ods html;
2PROC CORR DATA=new;
3 var _1 - _5;
4RUN;
5ods html close;
4 Bloc de code
PROC MEANS Data
Explication :
Agrégation des données : calcul de la somme des succès 'y' par cluster (subjid).
Copié !
1PROC MEANS DATA=correlated_bernoullis sum noprint;
2 class subjid;
3 var y;
4 OUTPUT out=out1 sum=t;
5RUN;
6 
7DATA out1;
8 SET out1;
9 IF _type_ = 1;
10 drop _type _freq_;
11RUN;
5 Bloc de code
PROC MEANS
Explication :
Analyse descriptive (moyenne et variance) de la variable agrégée 't'.
Copié !
1ods html;
2PROC MEANS DATA=out1 mean var maxdec=2;
3 var t;
4RUN;
5ods html close;
6 Bloc de code
DATA STEP Data
Explication :
Préparation de la simulation Monte Carlo : ajout d'une variable de réplication 'reps' pour diviser les données en 1000 sous-ensembles.
Copié !
1DATA correlated_bernoullis; *--- This data set was created in Program 1.1;
2 SET correlated_bernoullis;
3 reps = ceil( subjid / 20 ); *--- Creates 1000 reps of n=20 clusters of size m=5 each;
4RUN;
7 Bloc de code
PROC MEANS Data
Explication :
Calcul de la proportion estimée (Pi_hat) pour chaque réplication.
Copié !
1PROC MEANS DATA=correlated_bernoullis mean noprint;
2 class reps;
3 var y;
4 OUTPUT out=test mean=Pi_hat;
5RUN;
8 Bloc de code
DATA STEP Data
Explication :
Calcul du score Z et détermination de la significativité (P-value) pour tester l'hypothèse nulle (Pi=0.6) sur chaque réplication.
Copié !
1DATA test;
2SET test;
3 IF _type_ = 1;
4 z = (Pi_hat - 0.6) / sqrt( Pi_hat*(1-Pi_hat) / 100 );
5 P_val = 0;
6 IF z > 1.6449 THEN P_val = 1;
7RUN;
9 Bloc de code
PROC MEANS
Explication :
Estimation finale du taux d'erreur de type I en calculant la moyenne des rejets de l'hypothèse nulle sur l'ensemble des simulations.
Copié !
1title "Monte Carlo Estimate of Type I Error Rate";
2PROC MEANS DATA=test n mean maxdec=3;
3 var P_val;
4RUN;
Ce matériel est fourni "tel quel" par We Are Cas. Il n'y a aucune garantie, expresse ou implicite, quant à la qualité marchande ou à l'adéquation à un usage particulier concernant le matériel ou le code contenu dans les présentes. We Are Cas n'est pas responsable des erreurs dans ce matériel tel qu'il existe maintenant ou existera, et We Are Cas ne fournit pas de support technique pour celui-ci.