x <- c(19.4,19.4,19.4,19.5,19.5,19.5,28.7,28.7,28.7,21.8,21.8,21.8,20,20,20,22.6,22.6,22.6,22.4,22.4,22.4,18.6,18.6,18.6,16.2,16.2,16.2,13.8,13.8,13.8,24.1,24.1,24.1,19.9,19.9,19.9,22.3,22.3,22.3,20.9,20.9,20.9,23.5,23.5,23.5,23.1,23.1,23.1,25.7,25.7,25.7,19.7,19.7,19.7,23.1,23.1,23.1,20.7,20.7,20.7,18,18,18,16.9,16.9,16.9,24.4,24.4,24.4,15.5,15.5,15.5,18.4,18.4,18.4,16.2,16.2,16.2,20.6,20.6,20.6,19.8,19.8,19.8) par2 = '12' par1 = '50' par2 <- '12' par1 <- '50' #'GNU S' R Code compiled by R2WASP v. 1.2.291 () #Author: root #To cite this work: Wessa P., (2012), Blocked Bootstrap Plot for Central Tendency (v1.0.4) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_bootstrapplot.wasp/ #Source of accompanying publication: Office for Research, Development, and Education # par1 <- as.numeric(par1) par2 <- as.numeric(par2) if (par1 < 10) par1 = 10 if (par1 > 5000) par1 = 5000 if (par2 < 3) par2 = 3 if (par2 > length(x)) par2 = length(x) library(lattice) library(boot) boot.stat <- function(s) { s.mean <- mean(s) s.median <- median(s) s.midrange <- (max(s) + min(s)) / 2 c(s.mean, s.median, s.midrange) } (r <- tsboot(x, boot.stat, R=par1, l=12, sim='fixed')) postscript(file="/var/wessaorg/rcomp/tmp/1dlqm1387740125.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) plot(r$t[,1],type='p',ylab='simulated values',main='Simulation of Mean') grid() dev.off() postscript(file="/var/wessaorg/rcomp/tmp/2sgvt1387740125.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) plot(r$t[,2],type='p',ylab='simulated values',main='Simulation of Median') grid() dev.off() postscript(file="/var/wessaorg/rcomp/tmp/3y8su1387740125.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) plot(r$t[,3],type='p',ylab='simulated values',main='Simulation of Midrange') grid() dev.off() postscript(file="/var/wessaorg/rcomp/tmp/4dp1o1387740125.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) densityplot(~r$t[,1],col='black',main='Density Plot',xlab='mean') dev.off() postscript(file="/var/wessaorg/rcomp/tmp/50jus1387740125.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) densityplot(~r$t[,2],col='black',main='Density Plot',xlab='median') dev.off() postscript(file="/var/wessaorg/rcomp/tmp/6yhey1387740125.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) densityplot(~r$t[,3],col='black',main='Density Plot',xlab='midrange') dev.off() z <- data.frame(cbind(r$t[,1],r$t[,2],r$t[,3])) colnames(z) <- list('mean','median','midrange') postscript(file="/var/wessaorg/rcomp/tmp/7gf671387740125.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) boxplot(z,notch=TRUE,ylab='simulated values',main='Bootstrap Simulation - Central Tendency') grid() dev.off() #Note: the /var/wessaorg/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab load(file="/var/wessaorg/rcomp/createtable") a<-table.start() a<-table.row.start(a) a<-table.element(a,'Estimation Results of Blocked Bootstrap',6,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'statistic',header=TRUE) a<-table.element(a,'Q1',header=TRUE) a<-table.element(a,'Estimate',header=TRUE) a<-table.element(a,'Q3',header=TRUE) a<-table.element(a,'S.D.',header=TRUE) a<-table.element(a,'IQR',header=TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'mean',header=TRUE) q1 <- quantile(r$t[,1],0.25)[[1]] q3 <- quantile(r$t[,1],0.75)[[1]] a<-table.element(a,q1) a<-table.element(a,r$t0[1]) a<-table.element(a,q3) a<-table.element(a,sqrt(var(r$t[,1]))) a<-table.element(a,q3-q1) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'median',header=TRUE) q1 <- quantile(r$t[,2],0.25)[[1]] q3 <- quantile(r$t[,2],0.75)[[1]] a<-table.element(a,q1) a<-table.element(a,r$t0[2]) a<-table.element(a,q3) a<-table.element(a,sqrt(var(r$t[,2]))) a<-table.element(a,q3-q1) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'midrange',header=TRUE) q1 <- quantile(r$t[,3],0.25)[[1]] q3 <- quantile(r$t[,3],0.75)[[1]] a<-table.element(a,q1) a<-table.element(a,r$t0[3]) a<-table.element(a,q3) a<-table.element(a,sqrt(var(r$t[,3]))) a<-table.element(a,q3-q1) a<-table.row.end(a) a<-table.end(a) table.save(a,file="/var/wessaorg/rcomp/tmp/8yd3o1387740125.tab") try(system("convert tmp/1dlqm1387740125.ps tmp/1dlqm1387740125.png",intern=TRUE)) try(system("convert tmp/2sgvt1387740125.ps tmp/2sgvt1387740125.png",intern=TRUE)) try(system("convert tmp/3y8su1387740125.ps tmp/3y8su1387740125.png",intern=TRUE)) try(system("convert tmp/4dp1o1387740125.ps tmp/4dp1o1387740125.png",intern=TRUE)) try(system("convert tmp/50jus1387740125.ps tmp/50jus1387740125.png",intern=TRUE)) try(system("convert tmp/6yhey1387740125.ps tmp/6yhey1387740125.png",intern=TRUE)) try(system("convert tmp/7gf671387740125.ps tmp/7gf671387740125.png",intern=TRUE))