# Compute bootstrap distribution of sample medians
m<-1000
med.boot<-numeric(m)
for (i in 1:m) med.boot[i]<-median(sample(sids,replace=T))
hist(med.boot,nclass=25,prob=T)
summary(med.boot)

#this was run three times storing results in x1, x2, and x3
plot(density(x1,n=200),type='l')
lines(density(x2,n=200),lty=2)
lines(density(x3,n=200),lty=2)
title('3 Bootstrap Distributions for the Median')

m<-1000
mean.boot<-numeric(m)
for (i in 1:m) mean.boot[i]<-mean(sample(sids,replace=T))
summary(mean.boot)

#this was run 3 times, storing the result in y1, y2, y3
plot(density(y1,n=200),type='l')
lines(density(y2,n=200),lty=2)
lines(density(y3,n=200),lty=2)
title('3 Bootstrap Distributions for the Mean')

#just for fun, plot medians against means
for (i in 1:m) {x<-sample(sids,replace=T)
med.boot[i]<-median(x)
mean.boot[i]<-mean(x)}
plot(mean.boot,med.boot,main='Medians vs Means, Bootstrap Values')
lines(lowess(mean.boot,med.boot),lwd=3)

#compute approximate 95% CI based on percentiles of the bootstrap
distributions
quantile(med.boot,c(.025,.975))
quantile(mean.boot,c(.025,.975))