* For data without auxiliary variables ## Program to run Mplus for the simulation ## write the model gen.model<-function(datafile){ ## generate mplus scripts for data analysis ## Mplus input scripts mplus.null<-'model.inp' ## for null model cat('TITLE: Model\n',file=mplus.null); cat('DATA:\n',file=mplus.null, append=T); cat(' FILE=', datafile, ';\n',file=mplus.null, append=T); cat('VARIABLE: \n',file=mplus.null, append=T); cat(' NAMES ARE X M Y;\n', file=mplus.null, append=T); cat(' USEVARIABLES ARE X M Y;\n', file=mplus.null, append=T); cat(' MISSING = ALL(99999); \n',file=mplus.null, append=T); cat('ANALYSIS: \n',file=mplus.null, append=T); cat(' COVERAGE=.01; \n',file=mplus.null, append=T); cat(' bootstrap=1000; \n',file=mplus.null, append=T); cat('MODEL: \n',file=mplus.null, append=T); cat(' M ON X ;\n',file=mplus.null, append=T); cat(' Y ON M X; \n',file=mplus.null, append=T); cat('MODEL INDIRECT: \n',file=mplus.null, append=T); cat(' Y IND X; \n',file=mplus.null, append=T); cat('OUTPUT:\n', file=mplus.null, append=T); cat(' CINTERVAL (BCBOOTSTRAP); \n', file=mplus.null, append=T); } ## function to process output to get the parameter estimates and the bc intervals proc.res<-function(){ out.file<-"model.out" output <- suppressWarnings(readLines(out.file)) lnumber <- grep(" Total indirect", output) est<- output[lnumber[2]] #as.numeric(strsplit(output[lnumber[2]], " ")[[1]][c(16,22,29)]) est } ## run the simulation mis.pat<-c('mar', 'mcar', 'mnar') rate1<-c(1.5, 4, 9) rate2<-c(.1,.2,.4) for (pat in mis.pat){ if (pat == 'mar'){ rate <- rate1 }else{ rate <- rate2 } for (r in rate){ resfile <- paste(pat ,'-',r,'.txt', sep='') for (i in 1:1000){ datafile <- paste('C:\\zzy\\MNARmediation\\newdata\\', pat ,'\\data-',i,'-',r,'.txt', sep='') gen.model(datafile) system("c:/programs/mplus/mplus model.inp",show.output.on.console=F) temp<-proc.res() cat(c(i, temp), file=resfile, append=T) cat("\n", file=resfile,append=T) } } } * For data with auxiliary variables ## Program to run Mplus for the simulation ## write the model gen.model<-function(datafile){ ## generate mplus scripts for data analysis ## Mplus input scripts mplus.null<-'model.inp' ## for null model cat('TITLE: Model\n',file=mplus.null); cat('DATA:\n',file=mplus.null, append=T); cat(' FILE=', datafile, ';\n',file=mplus.null, append=T); cat('VARIABLE: \n',file=mplus.null, append=T); cat(' NAMES ARE X M Y U V;\n', file=mplus.null, append=T); cat(' USEVARIABLES ARE X M Y U V;\n', file=mplus.null, append=T); cat(' MISSING = ALL(99999); \n',file=mplus.null, append=T); cat('ANALYSIS: \n',file=mplus.null, append=T); cat(' COVERAGE=.01; \n',file=mplus.null, append=T); cat(' bootstrap=1000; \n',file=mplus.null, append=T); cat('MODEL: \n',file=mplus.null, append=T); cat(' M ON X ;\n',file=mplus.null, append=T); cat(' Y ON M X; \n',file=mplus.null, append=T); cat(' U with X M Y; \n',file=mplus.null, append=T); cat(' V with X M Y; \n',file=mplus.null, append=T); cat(' U with V; \n',file=mplus.null, append=T); cat('MODEL INDIRECT: \n',file=mplus.null, append=T); cat(' Y IND X; \n',file=mplus.null, append=T); cat('OUTPUT:\n', file=mplus.null, append=T); cat(' CINTERVAL (BCBOOTSTRAP); \n', file=mplus.null, append=T); } ## function to process output to get the parameter estimates and the bc intervals proc.res<-function(){ out.file<-"model.out" output <- suppressWarnings(readLines(out.file)) lnumber <- grep(" Total indirect", output) est<- output[lnumber[2]] #as.numeric(strsplit(output[lnumber[2]], " ")[[1]][c(16,22,29)]) est } ## run the simulation mis.pat<-c('mcar0','mcar1','mcar2','mar0','mar1','mar2','mnar0','mnar1','mnar2','mnar0-1','mnar1-1','mnar2-1') rate1<-c(1.5, 4, 9) rate2<-c(.1,.2,.4) for (pat in mis.pat){ if (substr(pat,1,3)=="mar"){ rate <- rate1 }else{ rate <- rate2 } for (r in rate){ resfile <- paste(pat ,'-',r,'.txt', sep='') for (i in 1:1000){ datafile <- paste('C:\\zzy\\MNARmediation\\newdata\\', pat ,'\\data-',i,'-',r,'.txt', sep='') gen.model(datafile) system("c:/programs/mplus/mplus model.inp",show.output.on.console=F) temp<-proc.res() cat(c(i, temp), file=resfile, append=T) cat("\n", file=resfile,append=T) } } }