From 2082f4efa3923fe083f8f29a36c6cad427de1013 Mon Sep 17 00:00:00 2001
From: corentinchoisy <corentin.choisy@proton.me>
Date: Fri, 2 Feb 2024 16:30:57 +0100
Subject: [PATCH] Updated README and descriptive analysis

---
 README.md                |  7 ++++
 RProject/desc_analysis.R | 91 ++++++++++++++++++++++++++++++++++++++++
 RProject/pcm.R           |  3 +-
 3 files changed, 100 insertions(+), 1 deletion(-)

diff --git a/README.md b/README.md
index 9843793..1b3e8ab 100644
--- a/README.md
+++ b/README.md
@@ -31,7 +31,14 @@ Ce dépôt contient l'ensemble du code pour les simulations basées sur simIRT.
 
 ## Conventions de nomenclature
 
+### Jeux de données initiaux
+
 **XX_N** - Scénarios de référence sans confusion / N réplications  
 **1XX_N** - TBA / N réplications  
 **2XX_N** - TBA / N réplications  
 **3XX_N** - TBA / N réplications
+
+### Jeux de données analysés
+
+**XX_N_nodif** - Analyse du scénario XX_N par PCM __sans__ prise en compte du DIF
+**XX_N** - Analyse du scénario XX_N par PCM __avec__ prise en compte du DIF
diff --git a/RProject/desc_analysis.R b/RProject/desc_analysis.R
index 97eb70c..9093daa 100644
--- a/RProject/desc_analysis.R
+++ b/RProject/desc_analysis.R
@@ -240,6 +240,60 @@ points(y=res.null5$h0.rejected.p,x=rep(1,nrow(res.null5)),col='#053305',pch=3)
 
 par(mfrow=c(1,1))
 
+##############################################################################
+#----------------------------------------------------------------------------#
+############################ BOXPLOTS H1 SCENARIOS ###########################
+#----------------------------------------------------------------------------#
+##############################################################################
+
+## Proportion of rejected h0 per dif value in h1 scenarios by DIF size // eff.size positive
+
+res.null <- res.dat[res.dat$eff.size>0,]
+boxplot(h0.rejected.p~dif.size,data=res.null,col=c(3,2,4,2,3),xlab='DIF size',ylab='H0 rejection proportion in target scenario')
+res.null0 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==0,]
+points(y=res.null0$h0.rejected.p,x=rep(3,nrow(res.null0)),col='darkblue',pch=3)
+res.null3 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==-0.3,]
+points(y=res.null3$h0.rejected.p,x=rep(2,nrow(res.null3)),col='#590b0c',pch=3)
+res.null5 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==-0.5,]
+points(y=res.null5$h0.rejected.p,x=rep(1,nrow(res.null5)),col='#053305',pch=3)
+res.null3 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==0.3,]
+points(y=res.null3$h0.rejected.p,x=rep(4,nrow(res.null3)),col='#590b0c',pch=3)
+res.null5 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==0.5,]
+points(y=res.null5$h0.rejected.p,x=rep(5,nrow(res.null5)),col='#053305',pch=3)
+
+
+############# By N 
+
+####### N=100
+
+res.null <- res.dat[res.dat$eff.size>0 & res.dat$N==100,]
+boxplot(h0.rejected.p~dif.size,data=res.null,col=c(3,2,4,2,3),xlab='DIF size',ylab='H0 rejection proportion in target scenario')
+res.null0 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==0 & res.dat$N==100,]
+points(y=res.null0$h0.rejected.p,x=rep(3,nrow(res.null0)),col='darkblue',pch=3)
+res.null3 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==-0.3 & res.dat$N==100,]
+points(y=res.null3$h0.rejected.p,x=rep(2,nrow(res.null3)),col='#590b0c',pch=3)
+res.null5 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==-0.5 & res.dat$N==100,]
+points(y=res.null5$h0.rejected.p,x=rep(1,nrow(res.null5)),col='#053305',pch=3)
+res.null3 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==0.3 & res.dat$N==100,]
+points(y=res.null3$h0.rejected.p,x=rep(4,nrow(res.null3)),col='#590b0c',pch=3)
+res.null5 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==0.5 & res.dat$N==100,]
+points(y=res.null5$h0.rejected.p,x=rep(5,nrow(res.null5)),col='#053305',pch=3)
+
+####### N=300 // DIF à 0.5 - QUELS SONT LES SCENARIOS EN HAUT
+
+res.null <- res.dat[res.dat$eff.size>0 & res.dat$N==300,]
+boxplot(h0.rejected.p~dif.size,data=res.null,col=c(3,2,4,2,3),xlab='DIF size',ylab='H0 rejection proportion in target scenario')
+res.null0 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==0 & res.dat$N==300,]
+points(y=res.null0$h0.rejected.p,x=rep(3,nrow(res.null0)),col='darkblue',pch=3)
+res.null3 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==-0.3 & res.dat$N==300,]
+points(y=res.null3$h0.rejected.p,x=rep(2,nrow(res.null3)),col='#590b0c',pch=3)
+res.null5 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==-0.5 & res.dat$N==300,]
+points(y=res.null5$h0.rejected.p,x=rep(1,nrow(res.null5)),col='#053305',pch=3)
+res.null3 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==0.3 & res.dat$N==300,]
+points(y=res.null3$h0.rejected.p,x=rep(4,nrow(res.null3)),col='#590b0c',pch=3)
+res.null5 <- res.dat[res.dat$eff.size>0 & res.dat$dif.size==0.5 & res.dat$N==300,]
+points(y=res.null5$h0.rejected.p,x=rep(5,nrow(res.null5)),col='#053305',pch=3)
+
 
 ##############################################################################
 #----------------------------------------------------------------------------#
@@ -354,3 +408,40 @@ res.null3 <- res.dat.temp[res.dat.temp$dif.size==-0.3,]
 points(y=res.null3$true.value.in.ci.p,x=rep(2,nrow(res.null3)),col='#053305',pch=3)
 res.null3 <- res.dat.temp[res.dat.temp$dif.size==0.3,]
 points(y=res.null3$true.value.in.ci.p,x=rep(3,nrow(res.null3)),col='#053305',pch=3)
+
+
+##############################################################################
+#----------------------------------------------------------------------------#
+########################## BETA SIGN CHANGE BOXPLOTS #########################
+#----------------------------------------------------------------------------#
+##############################################################################
+
+# Overall
+
+boxplot(beta.same.sign.truebeta.p~dif.size,data=res.dat,col=c(2,3),xlab='DIF size',
+        ylab='Proportion of estimates with same sign as true value in target scenario',main='DIF on 3 items',ylim=c(0,1))
+res.null3 <- res.dat[res.dat$dif.size==-0.5,]
+points(y=res.null3$beta.same.sign.truebeta.p,x=rep(1,nrow(res.null3)),col='#590b0c',pch=3)
+res.null3 <- res.dat[res.dat$dif.size==0.5,]
+points(y=res.null3$beta.same.sign.truebeta.p,x=rep(5,nrow(res.null3)),col='#590b0c',pch=3)
+res.null3 <- res.dat[res.dat$dif.size==-0.3,]
+points(y=res.null3$beta.same.sign.truebeta.p,x=rep(2,nrow(res.null3)),col='#053305',pch=3)
+res.null3 <- res.dat[res.dat$dif.size==0.3,]
+points(y=res.null3$beta.same.sign.truebeta.p,x=rep(4,nrow(res.null3)),col='#053305',pch=3)
+res.null3 <- res.dat[res.dat$dif.size==0,]
+points(y=res.null3$beta.same.sign.truebeta.p,x=rep(3,nrow(res.null3)),col='gray',pch=3)
+
+# Overall // H0 rejected 
+
+boxplot(beta.same.sign.truebeta.signif.p~dif.size,data=res.dat,col=c(2,3),xlab='DIF size',
+        ylab='Proportion of estimates with same sign as true value in target scenario',main='When H0 rejected',ylim=c(0,1))
+res.null3 <- res.dat[res.dat$dif.size==-0.5,]
+points(y=res.null3$beta.same.sign.truebeta.signif.p,x=rep(1,nrow(res.null3)),col='#590b0c',pch=3)
+res.null3 <- res.dat[res.dat$dif.size==0.5,]
+points(y=res.null3$beta.same.sign.truebeta.signif.p,x=rep(5,nrow(res.null3)),col='#590b0c',pch=3)
+res.null3 <- res.dat[res.dat$dif.size==-0.3,]
+points(y=res.null3$beta.same.sign.truebeta.signif.p,x=rep(2,nrow(res.null3)),col='#053305',pch=3)
+res.null3 <- res.dat[res.dat$dif.size==0.3,]
+points(y=res.null3$beta.same.sign.truebeta.signif.p,x=rep(4,nrow(res.null3)),col='#053305',pch=3)
+res.null3 <- res.dat[res.dat$dif.size==0,]
+points(y=res.null3$beta.same.sign.truebeta.signif.p,x=rep(3,nrow(res.null3)),col='gray',pch=3)
\ No newline at end of file
diff --git a/RProject/pcm.R b/RProject/pcm.R
index 3d4a279..d67760a 100644
--- a/RProject/pcm.R
+++ b/RProject/pcm.R
@@ -2134,7 +2134,8 @@ compile_simulation <- function(scenario) {
                   m.high.ci.beta=mean(s$high.ci.beta),
                   true.value.in.ci.p=mean(s$true.value.in.ci),
                   h0.rejected.p=mean(s$h0.rejected),
-                  beta.same.sign.truebeta.p=mean(s$beta.same.sign.truebeta,na.rm=T))
+                  beta.same.sign.truebeta.p=mean(s$beta.same.sign.truebeta,na.rm=T),
+                  beta.same.sign.truebeta.signif.p=mean(s[s$h0.rejected==1,]$beta.same.sign.truebeta,na.rm=T))
   d <- cbind(b,a,z)
   d$prop.
   return(d)