Few changes in R code

RProject/Scripts/Analysis/aggregation.R

@@ -215,6 +215,41 @@ compile_simulation <- function(scenario) {
   }
   N <- ifelse(substr(scenario,start=nchar(scenario)-1,stop=nchar(scenario))=="50","50",substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario)))
   zz <- ifelse(N=="50",substr(scenario,start=0,stop=nchar(scenario)-3),substr(scenario,start=0,stop=nchar(scenario)-4)) 
+  zz.type <- substr(zz,start=nchar(zz),stop=nchar(zz))
+  
+  if (unique(is.na(s$dif.size))) {
+    s$dif.size <- 0
+  }
+  if(zz=="10B") {
+    s$dif.size <- 0.3
+  }
+  if (substr(zz,1,1)%in%c("6","8")) {
+    s$nb.dif <- 1
+  }
+  if (substr(zz,1,2)%in%c("10","12","14",'16')) {
+    s$nb.dif <- 2
+  }
+  if (substr(zz,1,1)%in%c("18","20")) {
+    s$nb.dif <- 3
+  }
+  
+  if (zz.type!="A") {
+  if (zz.type=="B") {
+    s$eff.size <- 0.2
+  } else if (zz.type=="C" & unique(s$dif.size)==0) {
+    s$eff.size <- 0.4
+  } else if (zz.type=="C" & unique(s$dif.size)!=0) {
+    s$eff.size <- 0.2
+  } else if (zz.type=="D" & unique(s$dif.size)!=0) {
+    s$eff.size <- 0.4
+  } else if (zz.type=="E" & unique(s$dif.size)!=0) {
+    s$eff.size <- 0.4
+  }
+  } else {
+    s$eff.size <- 0
+    s$dif.size <- -1*s$dif.size
+  }
+  
   b <- data.frame(scenario=zz,
                   scenario.type=substr(zz,start=nchar(zz),stop=nchar(zz)),
                   N=N,
@@ -229,6 +264,7 @@ compile_simulation <- function(scenario) {
                   se.analytical.beta=mean(s$se.beta),
                   m.low.ci.beta=mean(s$low.ci.beta),
                   m.high.ci.beta=mean(s$high.ci.beta),
+                  bias=mean(s$beta-s$eff.size),
                   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),
@@ -247,18 +283,8 @@ for (x in results[seq(2,length(results))]) {
   res.dat <- bind_rows(res.dat,y)
 }
 
-res.dat[res.dat$scenario.type=='A','dif.size'] <- -res.dat[res.dat$scenario.type=='A','dif.size']
-res.dat[is.na(res.dat$dif.size),'dif.size'] <- 0
-res.dat[res.dat$scenario=="10B",]$dif.size <- 0.3 
-res.dat[substr(res.dat$scenario,1,1)%in%c("6","8"),'nb.dif'] <- 1
-res.dat[substr(res.dat$scenario,1,2)%in%seq(10,16,2),'nb.dif'] <- 2
-res.dat[substr(res.dat$scenario,1,2)%in%seq(18,20,2),'nb.dif'] <- 3
 res.dat[res.dat$N==50,"dif.size"] <- res.dat[which(res.dat$N==50)+1,"dif.size"]
-res.dat[res.dat$scenario.type=="B",]$eff.size <- 0.2
+
-res.dat[res.dat$scenario.type=="C" & res.dat$dif.size==0,]$eff.size <- 0.4
-res.dat[res.dat$scenario.type=="C" & res.dat$dif.size!=0,]$eff.size <- 0.2
-res.dat[res.dat$scenario.type=="D" & res.dat$dif.size!=0,]$eff.size <- 0.4
-res.dat[res.dat$scenario.type=="E" & res.dat$dif.size!=0,]$eff.size <- 0.4
 
 is.nan.data.frame <- function(x) {
   do.call(cbind, lapply(x, is.nan))
@@ -266,7 +292,6 @@ is.nan.data.frame <- function(x) {
 
 
 res.dat[is.nan(res.dat)] <- NA
-res.dat$bias <- res.dat$eff.size-res.dat$m.beta
 
 
 ##############################################################################

RProject/Scripts/Analysis/functions/resali.R

@@ -30,6 +30,7 @@ resali <- function(df=NULL,items=NULL,group=NULL,verbose=T) {
     pval[i] <- res.anova[[i]][1,"Pr(>F)"]
     fval[i] <- res.anova[[i]][1,'F value']
   }
+  print(res.anova)
   if (verbose) {
     cat('DONE\n')
     cat('-----------------------------------------------------------\n')

RProject/Scripts/Analysis/resali_analysis.R

@@ -4,7 +4,7 @@
 #----------------------------------------------------------------------------#
 ##############################################################################
 
-source(paste0(getwd(),"/functions/resali.R"))
+source("/home/corentin/Documents/These/Recherche/Simulations/RProject/Scripts/Analysis/functions/resali.R")
 
 generate_resali <- function(scenario=NULL,grp=NULL) {
   scen <- as.numeric(gsub("[A,B,C,D,E,F,G,_]","",substr(scenario,0,3)))

RProject/Scripts/Analysis/residif_analysis.R

@@ -0,0 +1,279 @@
+##############################################################################
+#----------------------------------------------------------------------------#
+################################### RESIDIF ##################################
+#----------------------------------------------------------------------------#
+##############################################################################
+
+generate_residif <- function(scenario=NULL,grp=NULL) {
+  scen <- as.numeric(gsub("[A,B,C,D,E,F,G,_]","",substr(scenario,0,3)))
+  if (substr(scenario,start=nchar(scenario)-1,stop=nchar(scenario))=="50") {
+    N <- 50
+  }
+  if (substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario))=="100") {
+    N <- 100
+  }
+  if (substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario))=="200") {
+    N <- 200
+  }
+  if (substr(scenario,start=nchar(scenario)-2,stop=nchar(scenario))=="300") {
+    N <- 300
+  }
+  if (scen<5) {
+    dat <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Data/NoDIF/N',N,'/scenario_',scenario,'.csv'))
+  }
+  if (scen>=5) {
+    dat <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Data/DIF/N',N,'/scenario_',scenario,'.csv'))
+  }
+  if (scen%in%c(3,4,13:20)) {
+    res <- residif(df=dat[dat$replication==1,],items = paste0("item",1:7),grp="TT",verbose=FALSE)
+    df_res <- data.frame(dif.detect.1=ifelse(length(res$dif.items)>=1,res$dif.items[1],NA),
+                         dif.detect.2=ifelse(length(res$dif.items)>=2,res$dif.items[2],NA),
+                         dif.detect.3=ifelse(length(res$dif.items)>=3,res$dif.items[3],NA),
+                         dif.detect.4=ifelse(length(res$dif.items)>=4,res$dif.items[4],NA),
+                         dif.detect.5=ifelse(length(res$dif.items)>=5,res$dif.items[5],NA),
+                         dif.detect.6=ifelse(length(res$dif.items)>=6,res$dif.items[6],NA),
+                         dif.detect.7=ifelse(length(res$dif.items)>=7,res$dif.items[7],NA),
+                         dif.detect.unif.1=ifelse(length(res$uniform)>=1,res$uniform[1],NA),
+                         dif.detect.unif.2=ifelse(length(res$uniform)>=2,res$uniform[2],NA),
+                         dif.detect.unif.3=ifelse(length(res$uniform)>=3,res$uniform[3],NA),
+                         dif.detect.unif.4=ifelse(length(res$uniform)>=4,res$uniform[4],NA),
+                         dif.detect.unif.5=ifelse(length(res$uniform)>=5,res$uniform[5],NA),
+                         dif.detect.unif.6=ifelse(length(res$uniform)>=6,res$uniform[6],NA),
+                         dif.detect.unif.7=ifelse(length(res$uniform)>=7,res$uniform[7],NA),
+                         N=N,
+                         nbdif=ifelse(scen<=4,0,ifelse(scen<=16,2,3)),
+                         true.dif.1=ifelse(scen<=4,NA,unique(dat[dat$replication==1,]$dif1)),
+                         true.dif.2=ifelse(scen<=4,NA,unique(dat[dat$replication==1,]$dif2)),
+                         true.dif.3=ifelse(scen<=16,NA,unique(dat[dat$replication==1,]$dif3))
+    )
+    for (k in 2:1000) {
+      if (k%%100==0) {
+        cat(paste0('N=',k,'/1000\n'))
+      }
+      res <- residif(df=dat[dat$replication==k,],items = paste0("item",1:7),grp="TT",verbose=FALSE)
+      df_res2 <- data.frame(dif.detect.1=ifelse(length(res$dif.items)>=1,res$dif.items[1],NA),
+                            dif.detect.2=ifelse(length(res$dif.items)>=2,res$dif.items[2],NA),
+                            dif.detect.3=ifelse(length(res$dif.items)>=3,res$dif.items[3],NA),
+                            dif.detect.4=ifelse(length(res$dif.items)>=4,res$dif.items[4],NA),
+                            dif.detect.5=ifelse(length(res$dif.items)>=5,res$dif.items[5],NA),
+                            dif.detect.6=ifelse(length(res$dif.items)>=6,res$dif.items[6],NA),
+                            dif.detect.7=ifelse(length(res$dif.items)>=7,res$dif.items[7],NA),
+                            dif.detect.unif.1=ifelse(length(res$uniform)>=1,res$uniform[1],NA),
+                            dif.detect.unif.2=ifelse(length(res$uniform)>=2,res$uniform[2],NA),
+                            dif.detect.unif.3=ifelse(length(res$uniform)>=3,res$uniform[3],NA),
+                            dif.detect.unif.4=ifelse(length(res$uniform)>=4,res$uniform[4],NA),
+                            dif.detect.unif.5=ifelse(length(res$uniform)>=5,res$uniform[5],NA),
+                            dif.detect.unif.6=ifelse(length(res$uniform)>=6,res$uniform[6],NA),
+                            dif.detect.unif.7=ifelse(length(res$uniform)>=7,res$uniform[7],NA),
+                            N=N,
+                            nbdif=ifelse(scen<=4,0,ifelse(scen<=16,2,3)),
+                            true.dif.1=ifelse(scen<=4,NA,unique(dat[dat$replication==k,]$dif1)),
+                            true.dif.2=ifelse(scen<=4,NA,unique(dat[dat$replication==k,]$dif2)),
+                            true.dif.3=ifelse(scen<=16,NA,unique(dat[dat$replication==k,]$dif3)))
+      df_res <- rbind(df_res,df_res2)
+    }
+  }
+  else if (scen%in%c(1,2,5:12)) {
+    res <- residif(df=dat[dat$replication==1,],items = paste0("item",1:4),grp="TT",verbose=FALSE)
+    df_res <- data.frame(dif.detect.1=ifelse(length(res$dif.items)>=1,res$dif.items[1],NA),
+                         dif.detect.2=ifelse(length(res$dif.items)>=2,res$dif.items[2],NA),
+                         dif.detect.3=ifelse(length(res$dif.items)>=3,res$dif.items[3],NA),
+                         dif.detect.4=ifelse(length(res$dif.items)>=4,res$dif.items[4],NA),
+                         dif.detect.unif.1=ifelse(length(res$uniform)>=1,res$uniform[1],NA),
+                         dif.detect.unif.2=ifelse(length(res$uniform)>=2,res$uniform[2],NA),
+                         dif.detect.unif.3=ifelse(length(res$uniform)>=3,res$uniform[3],NA),
+                         dif.detect.unif.4=ifelse(length(res$uniform)>=4,res$uniform[4],NA),
+                         N=N,
+                         nbdif=ifelse(scen<=4,0,ifelse(scen<=8,1,2)),
+                         true.dif.1=ifelse(scen<=4,NA,unique(dat[dat$replication==1,]$dif1)),
+                         true.dif.2=ifelse(scen<=8,NA,unique(dat[dat$replication==1,]$dif2))
+    )
+    for (k in 2:1000) {
+      if (k%%100==0) {
+        cat(paste0('N=',k,'/1000\n'))
+      }
+      res <- residif(df=dat[dat$replication==k,],items = paste0("item",1:4),grp="TT",verbose=FALSE)
+      df_res2 <- data.frame(dif.detect.1=ifelse(length(res$dif.items)>=1,res$dif.items[1],NA),
+                            dif.detect.2=ifelse(length(res$dif.items)>=2,res$dif.items[2],NA),
+                            dif.detect.3=ifelse(length(res$dif.items)>=3,res$dif.items[3],NA),
+                            dif.detect.4=ifelse(length(res$dif.items)>=4,res$dif.items[4],NA),
+                            dif.detect.unif.1=ifelse(length(res$uniform)>=1,res$uniform[1],NA),
+                            dif.detect.unif.2=ifelse(length(res$uniform)>=2,res$uniform[2],NA),
+                            dif.detect.unif.3=ifelse(length(res$uniform)>=3,res$uniform[3],NA),
+                            dif.detect.unif.4=ifelse(length(res$uniform)>=4,res$uniform[4],NA),
+                            N=N,
+                            nbdif=ifelse(scen<=4,0,ifelse(scen<=8,1,2)),
+                            true.dif.1=ifelse(scen<=4,NA,unique(dat[dat$replication==k,]$dif1)),
+                            true.dif.2=ifelse(scen<=8,NA,unique(dat[dat$replication==k,]$dif2)))
+      df_res <- rbind(df_res,df_res2)
+    }
+  }
+  return(df_res)
+}
+
+
+
+results <- c(sapply(c(2,4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
+
+results2 <- c(sapply(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
+
+results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
+
+results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
+
+results <- sort(results)
+
+results2 <- sort(results2)
+
+results <- c(results,results2)
+
+for (r in c(results[73:138])) {
+  cat(paste0(r,"\n"))
+  cat(paste0("-------------------------------------------","\n"))
+  write.csv(generate_residif(r,"TT"),paste0("/home/corentin/Documents/These/Recherche/residif/detection/",r,".csv"))
+  cat(paste0("-------------------------------------------","\n"))
+}
+
+##############################################################################
+#----------------------------------------------------------------------------#
+################################### NEWDATA ##################################
+#----------------------------------------------------------------------------#
+##############################################################################
+
+## Liste des scenarios
+
+results <- c(sapply(c(2,4),function(x) paste0(x,c('A','B','C'))),sapply(c(6,8),function(x) paste0(x,c('A','B','C','D','E'))))
+
+results2 <- c(sapply(seq(10,20,2),function(x) paste0(x,c('A','B','C','D','E'))))
+
+results <- c(sapply(c(50,100,300),function(x) paste0(results,'_',x)))
+
+results2 <- c(sapply(c(50,100,300),function(x) paste0(results2,'_',x)))
+
+results <- sort(results)
+
+results2 <- sort(results2)
+
+results <- c(results,results2)
+
+## Importer l'analyse resali pour chaque scenario
+
+for (r in results) {
+  cat('--------------------------------------------------------------------------\n')
+  cat(paste0(r,"\n"))
+  cat('--------------------------------------------------------------------------\n')
+  #### Importer les datas
+  scen <- as.numeric(gsub("[A,B,C,D,E,F,G,_]","",substr(r,0,3)))
+  if (substr(r,start=nchar(r)-1,stop=nchar(r))=="50") {
+    N <- 50
+  }
+  if (substr(r,start=nchar(r)-2,stop=nchar(r))=="100") {
+    N <- 100
+  }
+  if (substr(r,start=nchar(r)-2,stop=nchar(r))=="200") {
+    N <- 200
+  }
+  if (substr(r,start=nchar(r)-2,stop=nchar(r))=="300") {
+    N <- 300
+  }
+  if (scen<5) {
+    datt <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Data/NoDIF/N',N,'/scenario_',r,'.csv'))
+  }
+  if (scen>=5) {
+    datt <- read.csv(paste0('/home/corentin/Documents/These/Recherche/Simulations/Data/DIF/N',N,'/scenario_',r,'.csv'))
+  }
+  #### Importer l'analyse
+  analyse <- read.csv(paste0('/home/corentin/Documents/These/Recherche/residif/detection/',r,".csv"))
+  #### Pour chaque replication
+  for (k in 1:1000) {
+    if (k%%100==0) {
+      cat(paste0("N = ",k," / 1000\n"))
+    }
+    datt[datt$replication==k,"dif.detect.1"] <- analyse[analyse$X==k,"dif.detect.1"]
+    datt[datt$replication==k,"dif.detect.2"] <- analyse[analyse$X==k,"dif.detect.2"]
+    datt[datt$replication==k,"dif.detect.3"] <- analyse[analyse$X==k,"dif.detect.3"]
+    datt[datt$replication==k,"dif.detect.4"] <- analyse[analyse$X==k,"dif.detect.4"]
+    datt[datt$replication==k,"dif.detect.unif.1"] <- analyse[analyse$X==k,"dif.detect.unif.1"]
+    datt[datt$replication==k,"dif.detect.unif.2"] <- analyse[analyse$X==k,"dif.detect.unif.2"]
+    datt[datt$replication==k,"dif.detect.unif.3"] <- analyse[analyse$X==k,"dif.detect.unif.3"]
+    datt[datt$replication==k,"dif.detect.unif.4"] <- analyse[analyse$X==k,"dif.detect.unif.4"]
+    if (scen==3 | scen==4 | scen>=13) {
+      datt[datt$replication==k,"dif.detect.5"] <- analyse[analyse$X==k,"dif.detect.5"]
+      datt[datt$replication==k,"dif.detect.6"] <- analyse[analyse$X==k,"dif.detect.6"]
+      datt[datt$replication==k,"dif.detect.7"] <- analyse[analyse$X==k,"dif.detect.7"]
+      datt[datt$replication==k,"dif.detect.unif.5"] <- analyse[analyse$X==k,"dif.detect.unif.5"]
+      datt[datt$replication==k,"dif.detect.unif.6"] <- analyse[analyse$X==k,"dif.detect.unif.6"]
+      datt[datt$replication==k,"dif.detect.unif.7"] <- analyse[analyse$X==k,"dif.detect.unif.7"]
+    }
+  }
+  datt[is.na(datt$dif.detect.1),"dif.detect.1"] <- ""
+  datt[is.na(datt$dif.detect.2),"dif.detect.2"] <- ""
+  datt[is.na(datt$dif.detect.3),"dif.detect.3"] <- ""
+  datt[is.na(datt$dif.detect.4),"dif.detect.4"] <- ""
+  datt[is.na(datt$dif.detect.unif.1),"dif.detect.unif.1"] <- ""
+  datt[is.na(datt$dif.detect.unif.2),"dif.detect.unif.2"] <- ""
+  datt[is.na(datt$dif.detect.unif.3),"dif.detect.unif.3"] <- ""
+  datt[is.na(datt$dif.detect.unif.4),"dif.detect.unif.4"] <- ""
+  if (scen==3 | scen==4 | scen>=13) {
+    datt[is.na(datt$dif.detect.5),"dif.detect.5"] <- ""
+    datt[is.na(datt$dif.detect.6),"dif.detect.6"] <- ""
+    datt[is.na(datt$dif.detect.7),"dif.detect.7"] <- ""
+    datt[is.na(datt$dif.detect.unif.5),"dif.detect.unif.5"] <- ""
+    datt[is.na(datt$dif.detect.unif.6),"dif.detect.unif.6"] <- ""
+    datt[is.na(datt$dif.detect.unif.7),"dif.detect.unif.7"] <- ""
+  }
+  write.csv(datt,paste0("/home/corentin/Documents/These/Recherche/residif/detection_data/",r,".csv"))
+}
+
+
+
+
+
+
+
+
+
+for (r in results) {
+  cat('--------------------------------------------------------------------------\n')
+  cat(paste0(r,"\n"))
+  cat('--------------------------------------------------------------------------\n')
+  #### Importer les datas
+  scen <- as.numeric(gsub("[A,B,C,D,E,F,G,_]","",substr(r,0,3)))
+  if (substr(r,start=nchar(r)-1,stop=nchar(r))=="50") {
+    N <- 50
+  }
+  if (substr(r,start=nchar(r)-2,stop=nchar(r))=="100") {
+    N <- 100
+  }
+  if (substr(r,start=nchar(r)-2,stop=nchar(r))=="200") {
+    N <- 200
+  }
+  if (substr(r,start=nchar(r)-2,stop=nchar(r))=="300") {
+    N <- 300
+  }
+  #### Importer l'analyse
+  analyse <- read.csv(paste0("/home/corentin/Documents/These/Recherche/residif/detection_data/",r,".csv"))
+  analyse[is.na(analyse)] <- ""
+  names(analyse)[names(analyse)=="dif.detect.1"] <- "dif_detect_1"
+  names(analyse)[names(analyse)=="dif.detect.2"] <- "dif_detect_2"
+  names(analyse)[names(analyse)=="dif.detect.3"] <- "dif_detect_3"
+  names(analyse)[names(analyse)=="dif.detect.4"] <- "dif_detect_4"
+  names(analyse)[names(analyse)=="dif.detect.unif.1"] <- "dif_detect_unif_1"
+  names(analyse)[names(analyse)=="dif.detect.unif.2"] <- "dif_detect_unif_2"
+  names(analyse)[names(analyse)=="dif.detect.unif.3"] <- "dif_detect_unif_3"
+  names(analyse)[names(analyse)=="dif.detect.unif.4"] <- "dif_detect_unif_4"
+  
+  
+  if (scen==3 | scen==4 | scen>=13) {
+    names(analyse)[names(analyse)=="dif.detect.5"] <- "dif_detect_5"
+    names(analyse)[names(analyse)=="dif.detect.6"] <- "dif_detect_6"
+    names(analyse)[names(analyse)=="dif.detect.7"] <- "dif_detect_7"
+    names(analyse)[names(analyse)=="dif.detect.unif.5"] <- "dif_detect_unif_5"
+    names(analyse)[names(analyse)=="dif.detect.unif.6"] <- "dif_detect_unif_6"
+    names(analyse)[names(analyse)=="dif.detect.unif.7"] <- "dif_detect_unif_7"
+    
+  }
+  analyse <- analyse[,!names(analyse) %in% c("X","X.1","X.2")]
+  write.csv(analyse,paste0("/home/corentin/Documents/These/Recherche/residif/detection_data/",r,".csv"))
+}
+

RProject/Scripts/Tools/generate_item_difficulties.R

@@ -15,7 +15,7 @@ generate_diff_irt <- function(J=7,M=4) {
   }
   difficulties = matrix(c(0), J,M-1)
   rownames(difficulties)=paste("item",1:J)
-  colnames(difficulties)=paste("Moda", 1:(M-1))
+  colnames(difficulties)=paste("delta", 1:(M-1))
   for (j in 1:J){
     difficulties[j,1] = qnorm(j/(J+1))
   }