From 8bffc896b997da5852355e5abf1765ed1dcb964a Mon Sep 17 00:00:00 2001 From: corentinchoisy Date: Tue, 2 Apr 2024 14:38:53 +0200 Subject: [PATCH] Updated RESALI with new bonferroni correction --- RProject/resali.R | 217 +++++++++++++++++++++++++++++++++++++++++++--- 1 file changed, 203 insertions(+), 14 deletions(-) diff --git a/RProject/resali.R b/RProject/resali.R index bced126..968a53c 100644 --- a/RProject/resali.R +++ b/RProject/resali.R @@ -2,9 +2,9 @@ library(TAM) resali <- function(df=NULL,items=NULL,group=NULL,verbose=T) { if (verbose) { - cat('-----------------------------------------------------------\n') - cat('COMPUTING INITIAL PCM\n') - cat('-----------------------------------------------------------\n') + cat('-----------------------------------------------------------\n') + cat('COMPUTING INITIAL PCM\n') + cat('-----------------------------------------------------------\n') } nbitems <- length(items) nbitems_o <- nbitems @@ -14,28 +14,29 @@ resali <- function(df=NULL,items=NULL,group=NULL,verbose=T) { pcm_initial <- TAM::tam.mml(resp=resp,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F) dat <- df dat$score <- rowSums(dat[,items_n]) - dat$score_q5 <- cut(dat$score,quantile(dat$score,seq(0,1,0.2)),labels=1:5,include.lowest=T) + nqt <- ifelse(length(unique(quantile(dat$score,seq(0,1,0.2))))==6,5,length(unique(quantile(dat$score,seq(0,1,0.2))))-1) + dat$score_q5 <- cut(dat$score,unique(quantile(dat$score,seq(0,1,1/nqt))),labels=1:nqt,include.lowest=T) res.anova <- rep(NA,nbitems) pval <- rep(NA,nbitems) fval <- rep(NA,nbitems) for (i in items) { - dat[,paste0('res_',i)] <- rowMeans(predict(pcm_initial)$stand.resid[,,i]) + dat[,paste0('res_',i)] <- IRT.residuals(pcm_initial)$stand_residuals[,i] res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat)) pval[i] <- res.anova[[i]][1,"Pr(>F)"] fval[i] <- res.anova[[i]][1,'F value'] } if (verbose) { - cat('DONE\n') - cat('-----------------------------------------------------------\n') + cat('DONE\n') + cat('-----------------------------------------------------------\n') } res.items <- c() res.uniform <- c() k <- 1 - while(any(pval<0.05/nbitems_o)) { + while(any(pval<0.05/(nbitems_o*3))) { k <- k+1 if (verbose) { - cat(paste('COMPUTING STEP',k,'\n')) - cat('-----------------------------------------------------------\n') + cat(paste('COMPUTING STEP',k,'\n')) + cat('-----------------------------------------------------------\n') } res.item <- gsub("[a-z]", "",colnames(resp)[which.max(fval)]) res.items <- c(res.items,res.item) @@ -52,19 +53,103 @@ resali <- function(df=NULL,items=NULL,group=NULL,verbose=T) { pval <- rep(NA,nbitems) fval <- rep(NA,nbitems) for (i in 1:nbitems) { - dat[,paste0('res_',i)] <- rowMeans(predict(pcm_while)$stand.resid[,,i]) + dat[,paste0('res_',i)] <- IRT.residuals(pcm_while)$stand_residuals[,i] res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat)) pval[i] <- res.anova[[i]][1,"Pr(>F)"] fval[i] <- res.anova[[i]][1,'F value'] } if (verbose) { + cat('DONE\n') + cat('-----------------------------------------------------------\n') + } + } + if (verbose) { + cat("DETECTED DIF ITEMS\n") + cat('-----------------------------------------------------------\n') + } + if (length(res.items>0)) { + results <- data.frame(dif.items=res.items, + uniform=1*res.uniform) + return(results) + } + else { + if (verbose) { + cat("No DIF was detected\n") + } + return(NULL) + } +} + + + + + +resali_BF2 <- function(df=NULL,items=NULL,group=NULL,verbose=T) { + if (verbose) { + cat('-----------------------------------------------------------\n') + cat('COMPUTING INITIAL PCM\n') + cat('-----------------------------------------------------------\n') + } + nbitems <- length(items) + nbitems_o <- nbitems + items_n <- paste0('item',items) + resp <- df[,items_n] + grp <- df[,group] + pcm_initial <- TAM::tam.mml(resp=resp,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F) + dat <- df + dat$score <- rowSums(dat[,items_n]) + nqt <- ifelse(length(unique(quantile(dat$score,seq(0,1,0.2))))==6,5,length(unique(quantile(dat$score,seq(0,1,0.2))))-1) + dat$score_q5 <- cut(dat$score,unique(quantile(dat$score,seq(0,1,1/nqt))),labels=1:nqt,include.lowest=T) + res.anova <- rep(NA,nbitems) + pval <- rep(NA,nbitems) + fval <- rep(NA,nbitems) + for (i in items) { + dat[,paste0('res_',i)] <- rowMeans(predict(pcm_initial)$stand.resid[,,i]) + res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat)) + pval[i] <- res.anova[[i]][1,"Pr(>F)"] + fval[i] <- res.anova[[i]][1,'F value'] + } + if (verbose) { cat('DONE\n') cat('-----------------------------------------------------------\n') + } + res.items <- c() + res.uniform <- c() + k <- 1 + while(any(pval<0.05/(nbitems_o-k+1))) { + k <- k+1 + if (verbose) { + cat(paste('COMPUTING STEP',k,'\n')) + cat('-----------------------------------------------------------\n') + } + res.item <- gsub("[a-z]", "",colnames(resp)[which.max(fval)]) + res.items <- c(res.items,res.item) + res.uni <- res.anova[[which.max(fval)]][3,"Pr(>F)"]>0.05 + res.uniform <- c(res.uniform,res.uni) + items_n <- c(items_n[items_n!=paste0('item',res.item)],paste0("item",res.item,c("noTT","TT"))) + dat[,paste0("item",res.item,'TT')] <- dat[dat$TT==1,paste0('item',res.item)] + dat[,paste0("item",res.item,'noTT')] <- dat[dat$TT==0,paste0('item',res.item)] + resp <- dat[,items_n] + grp <- dat[,group] + pcm_while <- TAM::tam.mml(resp=resp,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F) + nbitems <- length(items_n) + res.anova <- rep(NA,nbitems) + pval <- rep(NA,nbitems) + fval <- rep(NA,nbitems) + for (i in 1:nbitems) { + dat[,paste0('res_',i)] <- rowMeans(predict(pcm_while)$stand.resid[,,i]) + res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat)) + pval[i] <- res.anova[[i]][1,"Pr(>F)"] + fval[i] <- res.anova[[i]][1,'F value'] + } + if (verbose) { + cat('DONE\n') + cat('-----------------------------------------------------------\n') } } if (verbose) { - cat("DETECTED DIF ITEMS\n") - cat('-----------------------------------------------------------\n') + cat("DETECTED DIF ITEMS\n") + cat('-----------------------------------------------------------\n') } if (length(res.items>0)) { results <- data.frame(dif.items=res.items, @@ -72,7 +157,111 @@ resali <- function(df=NULL,items=NULL,group=NULL,verbose=T) { return(results) } else { - cat("No DIF was detected") + if (verbose) { + cat("No DIF was detected\n") + } + return(NULL) + } +} + + + + + + + + +resali_LRT <- function(df=NULL,items=NULL,group=NULL,verbose=T) { + if (verbose) { + cat('-----------------------------------------------------------\n') + cat('COMPUTING INITIAL PCM\n') + cat('-----------------------------------------------------------\n') + } + nbitems <- length(items) + nbitems_o <- nbitems + items_n <- paste0('item',items) + resp <- df[,items_n] + grp <- df[,group] + items_n_alt <- paste0(items_n,c("_TT","_noTT")) + for (i in items_n) { + df[df$TT==0,paste0(i,"_noTT")] <- df[df$TT==0,i] + df[df$TT==1,paste0(i,"_TT")] <- df[df$TT==1,i] + } + resp_alt <- df[,items_n_alt] + pcm_initial <- TAM::tam.mml(resp=resp,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F) + pcm_alt <- TAM::tam.mml(resp=resp_alt,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F) + dat <- df + dat$score <- rowSums(dat[,items_n]) + nqt <- ifelse(length(unique(quantile(dat$score,seq(0,1,0.2))))==6,5,length(unique(quantile(dat$score,seq(0,1,0.2))))-1) + dat$score_q5 <- cut(dat$score,unique(quantile(dat$score,seq(0,1,1/nqt))),labels=1:nqt,include.lowest=T) + res.anova <- rep(NA,nbitems) + pval <- rep(NA,nbitems) + fval <- rep(NA,nbitems) + for (i in items) { + dat[,paste0('res_',i)] <- rowMeans(predict(pcm_initial)$stand.resid[,,i]) + res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat)) + pval[i] <- res.anova[[i]][1,"Pr(>F)"] + fval[i] <- res.anova[[i]][1,'F value'] + } + if (verbose) { + cat('DONE\n') + cat('-----------------------------------------------------------\n') + } + res.items <- c() + res.uniform <- c() + k <- 1 + if (anova(pcm_initial,pcm_alt)$p[1]<0.05) { + while(any(pval<0.05/nbitems_o)) { + k <- k+1 + if (verbose) { + cat(paste('COMPUTING STEP',k,'\n')) + cat('-----------------------------------------------------------\n') + } + res.item <- gsub("[a-z]", "",colnames(resp)[which.max(fval)]) + res.items <- c(res.items,res.item) + res.uni <- res.anova[[which.max(fval)]][3,"Pr(>F)"]>0.05 + res.uniform <- c(res.uniform,res.uni) + items_n <- c(items_n[items_n!=paste0('item',res.item)],paste0("item",res.item,c("noTT","TT"))) + dat[,paste0("item",res.item,'TT')] <- dat[dat$TT==1,paste0('item',res.item)] + dat[,paste0("item",res.item,'noTT')] <- dat[dat$TT==0,paste0('item',res.item)] + resp <- dat[,items_n] + grp <- dat[,group] + pcm_while <- TAM::tam.mml(resp=resp,Y=grp,irtmodel = "PCM",est.variance = T,verbose=F) + nbitems <- length(items_n) + res.anova <- rep(NA,nbitems) + pval <- rep(NA,nbitems) + fval <- rep(NA,nbitems) + for (i in 1:nbitems) { + dat[,paste0('res_',i)] <- rowMeans(predict(pcm_while)$stand.resid[,,i]) + res.anova[i] <- summary(aov(dat[,paste0('res_',i)]~TT*score_q5,data=dat)) + pval[i] <- res.anova[[i]][1,"Pr(>F)"] + fval[i] <- res.anova[[i]][1,'F value'] + } + if (verbose) { + cat('DONE\n') + cat('-----------------------------------------------------------\n') + } + } + if (verbose) { + cat("DETECTED DIF ITEMS\n") + cat('-----------------------------------------------------------\n') + } + if (length(res.items>0)) { + results <- data.frame(dif.items=res.items, + uniform=1*res.uniform) + return(results) + } + else { + if (verbose) { + cat("No DIF was detected\n") + } + return(NULL) + } + } + else { + if (verbose) { + cat("No DIF was detected\n") + } return(NULL) } } \ No newline at end of file