# stdmean.s: standardized mean difference model # + cumulative meta analysis # d<-scan("c:/MetaAnalysis/Rcode/normal.txt",list(n1=0,m1=0,s1=0,n0=0,m0=0,s0=0)) k<-length(d$n1) # # parameters for figures xposi <- 3 # x-axis for titles yup <- 11 # upper limit of y-axis ylow <- -4 # lower limit of y-axis xup <- 5 # upper limit of x-axis xlow <- -5 # lower limit of x-axis label<-rep("",k) label[1]<-"Edinburgh" label[2]<-"Orpington-Mild" label[3]<-"Orpington-Moderate" label[4]<-"Orpington-Severe" label[5]<-"Montreal-Home" label[6]<-"Montreal-Transfer" label[7]<-"Newcastle 1993" label[8]<-"Umea 1985" label[9]<-"Uppsala 1982" # par(mar=c(5,4,4,3)) dd<-d ad<-d$m1-d$m0 s<-( (d$n1-1)*d$s1^2 + (d$n0-1)*d$s0^2 )/(d$n1+d$n0-2) #s<-sqrt((1/d$n1+1/d$n0)*s) : modifed in June 2004 s<-sqrt(s) std<-ad/s se<-sqrt(1/d$n1+1/d$n0+std^2/2/(d$n1+d$n0) ) low<-std-1.96*se upp<-std+1.96*se tn<-d$n1+d$n0 # ---------- individual graph ------------------ id<-k:1 ind<-1:k res<-0.3 plot(std,id,ylim=c(ylow, yup), xlim=c(xlow, xup),yaxt="n",pch=" ", ylab="Citation",xlab="Standardized difference in means") title(main=" Standardized mean difference model ( circles are cumulative MA ) ") symbols(std,id+res,squares=sqrt(tn), add=T,inches=0.20) abline(v=0) for (i in 1:k){ j<-k-i+1 x<-c(low[i],upp[i]) y<-c(j,j)+0.3 lines(x,y,type="l") text(xposi ,i, label[j]) } # ----- cumulative graph ---- hh<-1:k xh<-hh xhl<-hh xhu<-hh for (ir in 1:k){ d$n1<-dd$n1[ind<=ir] d$m1<-dd$m1[ind<=ir] d$s1<-dd$s1[ind<=ir] d$n0<-dd$n0[ind<=ir] d$m0<-dd$m0[ind<=ir] d$s0<-dd$s0[ind<=ir] ad<-d$m1-d$m0 s<-( (d$n1-1)*d$s1^2 + (d$n0-1)*d$s0^2 )/(d$n1+d$n0-2) s<-sqrt(s) std<-ad/s se<-sqrt(1/d$n1+1/d$n0+std^2/2/(d$n1+d$n0) ) tn<-sum(d$n1+d$n0) hh[ir]<-tn # ----- fixed effects ---- w<-1/se/se sw<-sum(w) stdm<-sum(std*w)/sw stdl<-stdm-1.96*sqrt(1/sw) stdu<- stdm+1.96*sqrt(1/sw) q1<-sum(w*(std-stdm)^2) df1<-k-1 pval1<- 1-pchisq(q1,df1) # I2 <- (q1-df1)/q1 * 100 # q2<-stdm^2*sw df2<-1 pval2<- 1-pchisq(q2,df2) # ----- random-effects ---- tau2<-(q1-(k-1))/(sw-sum(w*w)/sw) tau2<-max(0,tau2) wx<-1/(tau2+se*se) swx<-sum(wx) stdmx<-sum(std*wx)/swx stdxl<- stdmx-1.96*sqrt(1/swx) stdxu<- stdmx+1.96*sqrt(1/swx) qx2<-stdmx^2*swx pvalx2<- 1-pchisq(qx2,df2) # -- graph of cumulative meta ---- j<-k-ir+1 x<-c(stdxl,stdxu) xh[ir]<- stdmx xhl[ir]<-stdxl xhu[ir]<-stdxu y<-c(j,j) lines(x,y, lty=3) } # -------- graph ---------------- symbols(xh,id,circles=sqrt(hh), add=T,inches=0.10) x<-c(stdl,stdu) y<-c(-1, -1) lines(x,y,type="b") x<-c(stdxl,stdxu) y<-c(-2, -2) lines(x,y,type="b") abline(v=c(stdm, stdmx), lty=2) text(xposi,-1, "Combined : fixed") text(xposi,-2, "Combined : random") # ------ output variables ------- out <- round( c(stdm, stdl, stdu), 4) outq1<- round( c(q1, df1, pval1), 6) outq2<- round( c(q2, df2, pval2), 6) outR<- round( c(stdmx, stdxl, stdxu), 4) outq2R<- round( c(qx2, df2, pvalx2), 6)