bootmnsdcv <- function(x,iter) { x <- x[!is.na(x)] meanx <- mean(x) ; print(c("The observed mean value = ",meanx)) sdx <- sd(x) ; print(c("The observed sd value = ",sdx)) cvx <- 100*sdx/meanx ; print(c("The observed cv value = ",cvx)) bootdistmean <- numeric(iter) bootdistsd <- numeric(iter) bootdistcv <- numeric(iter) for (i in 1:iter) { bootsampx <- sample(x,replace=TRUE) bootdistmean[i] <- mean(bootsampx) bootdistsd[i] <- sd(bootsampx) bootdistcv[i] <- 100*bootdistsd[i]/bootdistmean[i] } biasmean <- mean(bootdistmean) - meanx biassd <- mean(bootdistsd) - sdx biascv <- mean(bootdistcv) - cvx print(c(biasmean,biassd,biascv)) msemean <- var(bootdistmean) + biasmean^2 msesd <- var(bootdistsd) + biassd^2 msecv <- var(bootdistcv) + biascv^2 print("Bootstrap MSE estimates for the mean, sd, and cv") print(c(msemean,msesd,msecv)) } # these values should be close to those on page 251 table811 <- c(7,11,15,16,20,22,24,25,29,33,34,37,41,42,49,57,66,71,84,90) bootmnsdcv(table811,5000)