# We fit an exponential model to the air-conditioning data and use
# that for a parametric bootstrap. Then we look at plots of the
# resampled means.
air.rg <- function(data, mle)
rexp(length(data), 1/mle)
data(aircondit)
air.boot <- boot(aircondit$hours, mean, R=999, sim="parametric",
ran.gen=air.rg, mle=mean(aircondit$hours))
plot(air.boot)
# In the difference of means example for the last two series of the
# gravity data
data(gravity)
grav1 <- gravity[as.numeric(gravity[,2])>=7,]
grav.fun <- function(dat, w)
{ strata <- tapply(dat[, 2], as.numeric(dat[, 2]))
d <- dat[, 1]
ns <- tabulate(strata)
w <- w/tapply(w, strata, sum)[strata]
mns <- tapply(d * w, strata, sum)
mn2 <- tapply(d * d * w, strata, sum)
s2hat <- sum((mn2 - mns^2)/ns)
c(mns[2]-mns[1],s2hat)
}
grav.boot <- boot(grav1, grav.fun, R=499, stype="w", strata=grav1[,2])
plot(grav.boot)
# now suppose we want to look at the studentized differences.
grav.z <- (grav.boot$t[,1]-grav.boot$t0[1])/sqrt(grav.boot$t[,2])
plot(grav.boot,t=grav.z,t0=0)
# In this example we look at the one of the partial correlations for the
# head dimensions in the dataset frets.
pcorr <- function( x )
{
# Function to find the correlations and partial correlations between
# the four measurements.
v <- cor(x);
v.d <- diag(var(x));
iv <- solve(v);
iv.d <- sqrt(diag(iv));
iv <- - diag(1/iv.d) %*% iv %*% diag(1/iv.d);
q <- NULL;
n <- nrow(v);
for (i in 1:(n-1))
q <- rbind( q, c(v[i,1:i],iv[i,(i+1):n]) );
q <- rbind( q, v[n,] );
diag(q) <- round(diag(q));
q
}
frets.fun <- function( data, i )
{ d <- data[i,];
v <- pcorr( d );
c(v[1,],v[2,],v[3,],v[4,])
}
data(frets)
frets.boot <- boot(log(as.matrix(frets)), frets.fun, R=999)
plot(frets.boot, index=7, jack=TRUE, stinf=FALSE, useJ=FALSE)
Run the code above in your browser using DataLab