Example of using a bootstrap t-test for a non-normal distribution.
Here, take a sample of size n=20 from
\(X\sim \mbox{Exp}(\lambda)\) with
\(\mu = E[X] = 1/\lambda\) and \(\sigma^2 =
\mbox{Var}[X]=1/\lambda^2\)
and treat this as the given data. Then find the bootstrap distribution
from the data and use that to evaluate the p-value for the two-sided
hypothesis test.
Generate a sample of size n:
n <- 20
lambda <- 0.5
mu <- 1/lambda
sigma <- 1/lambda
set.seed(11) # sets starting point for pseudorandom generator
data <- rexp(n,lambda)
data_mean <- mean(data)
data_sd <- sd(data)
data_max <- max(data)
data_min <- min(data)
theory mean = 2
data mean = 2.9010401
theory std dev = 2
data std dev = 3.3027627
Plot the histogram of the sample and scaled version of the pdf
hist(data,xlim=c(0,10))
xfine<-seq(0,10,length=101)
scale = 30
lines(xfine,scale*dexp(xfine,lambda),col='blue',type='l',xlab='x',ylab='pdf')
abline(v=mu,col='blue')
abline(v=data_mean,col='red',lty='dashed')
Note: Unless n is large, the mean of the data (red dashed) is not necessarily close to the theoretical mean (blue solid).
Exact result for p-value from Gamma
distribution
Note, the Exponential distribution \(X \sim
\mbox{Exp}(\lambda)\) is a special case of the Gamma distribution
\(G \sim \mbox{Gamma}(r,\lambda)\) with
\(r=1\) (see Appendix B of
textook).
Since the result for the sampling mean for the Gamma distribution is
known:
\(\bar{G} \sim \mbox{Gamma}(nr,n
\lambda)\),
we have the exact result for the sampling distribution of the
Exponential distribution as
\(\bar{X} \sim
\mbox{Gamma}(n,n\lambda)\).
Thus the exact two-sided p value can be calculated from
pexact_minus = pgamma(data_mean,n,n*lambda)
pexact_plus=1 - pexact_minus
pexact2 = 2*(min(pexact_minus,pexact_plus))
pexact_minus = 0.9675011
pexact_plus = 0.0324989
pexact2 = 0.0649978
xfine <- seq(0,5,length=101)
plot(xfine,dgamma(xfine,n,n*lambda),col='blue',type='l')
abline(v=data_mean,col='red',lty='dashed')
Bootstrap distribution for T-statistic
We now use the bootstrap approximation for the T-statistic
\(T = \frac{\bar{X} -
\mu}{S/\sqrt{n}}\),
where the bootstrap distribution is based on resampling of the data with
replacement and using
\(T^* = \frac{\bar{X^*} -
\bar{x}}{S^*/\sqrt{n}}\)
where \(\bar{X^*}\) is the mean and
\(S^*\) is the standard deviation of
each resample. By computing \(T^*\)
using \(\bar{x}\) instead of \(\mu\), we obtain the shape of the bootstrap
distribution relative to the mean of the data.
Then, for determining the p-value, we use
\(t =
\frac{\bar{x}-\mu}{s/\sqrt{n}}\)
to determine how large the observed \(\bar{x}\) is relative to the assumed mean
of the null distribution \(\mu\).
N <- 10**4 # number of samples
tboot <- numeric(N) # array for bootstrap dist
# loop to create sample, calculate t statistic,
for (i in seq(1,N)) {
resample <- sample(data,n,replace=TRUE)
xbar <- mean(resample)
s <- sd(resample)
tboot[i] <- (xbar-data_mean)/(s/sqrt(n)) # gives shape relative to data_mean
}
t_observed <- (data_mean - mu)/(data_sd/sqrt(n))
hist(tboot,breaks=20)
abline(v=t_observed,col='red',lty='dashed')
Note the skewness of the bootstrap T distribution is reversed from the sampling distribution of \(\bar{X}\). This is because of the variability of the standard deviation of the resamples. Because all the data is positive, resamples with many small values are biased to have a smaller standard deviation, leading to larger values of \(T^* = \frac{\bar{X^*} - \bar{x}}{S^*/\sqrt{n}}\)
The p-value can be estimated from the bootstrap distribution from \(P(\bar{X} \le \bar{x}) \approx P(T^* \le t)\).
pt_minus <- (sum(tboot <= t_observed) + 1)/(N+1)
pt_plus <- (sum(tboot >= t_observed)+1)/(N+1)
pt2 <- 2*min(pt_minus,pt_plus) # two-sided p value
Standard T-test
Here it is assumed that the data is from a normal distribution with unknown mean and variance.
t <- (data_mean - mu)/(data_sd/sqrt(n))
pstd_minus = pt(t,df=n-1)
pstd_plus = 1 - pstd_minus
pstd2 = 2*min(pstd_minus,pstd_plus)
Comparison of methods
Method | P minus | P plus | P two-sided |
---|---|---|---|
Xbar exact | 0.9675011 | 0.0324989 | 0.0649978 |
T bootstrap | 0.9361064 | 0.0639936 | 0.1279872 |
standard T | 0.8813153 | 0.1186847 | 0.2373694 |
According to Sec 8.2 of the text, the error for bootstrap t-test scales like 1/n for n large, whereas the error for the standard T-test scales like \(1/\sqrt{n}\) for skewed distributions. Both approximations become exact as \(n\rightarrow \infty\) but the T bootstrap converges faster.