# Patch data
par(mfrow=c(2,2))     # Sets up multiple plots 2x2

library(bootstrap)

# First we look at the data
#str(patch)
attach(patch)   #attaches dataframe so we can just use the variable names
plot(subject,placebo,ylim=c(6000,30000))
points(subject,oldpatch,pch=4,col="red")
points(subject,newpatch,pch=18,col="blue")
abline(h=mean(oldpatch),col="red")
abline(h=mean(newpatch),col="blue")
# Set up the dataframe for bootstrapping

patchdata<-data.frame(new_old=patch$y,old_placebo=patch$z)
plot(patchdata, xlim=c(-4000,10000),ylim=c(-4000,10000))


# With bivariate data we define the basic data set as the case indices
# ie 1:8 in this example.  The bootstrap command will then
# generate bootstrap samples from 1:8 and input these samples
# to the statistic.  So the function theta defining the
# statistic is written as a function of `ind' - the case indices.

theta <- function(ind){
      y <- patchdata[ind,"new_old"]
      z <- patchdata[ind,"old_placebo"]
      mean(y)/mean(z) }
patch.boot <- bootstrap(1:8, 400, theta)
hist(patch.boot$thetastar)
abline(v=c(-0.2, 0.2), col="red2")  
abline(v=theta(1:8) , col="blue")
abline(v=mean(patch.boot$thetastar),col="purple")
theta(1:8) 
show(mean(patch.boot$thetastar) - theta(1:8))
show(sd(patch.boot$thetastar))