# Bootstrapping the t-Test

my problem is of practical nature regarding bootstrapping the t-Test.

My approach is to code a function which resamples the data vector and calculates the t-statistic for each new (resampled) vector. Then I want that the type 1 error of this function converges to something around 5% by using:

mean(replicate(2000, my_fun(rnorm(10),...)))


in R (my_fun is the bootstrapped t-Test returning TRUE if alternative is rejected). My underlying calculations of the t-Statistics are correct (I checked them manually and also tried coding them in C, C++, R). Following problem occurs: If I use classic non-parametric resampling such as:

x.sim = matrix(sample(x, replace = T, size = nboot * length(x)), nrow = nboot)


where x is the original data vector fed and nboot the number of resamplings, my type 1 error is 0 (irrespective of number of simulations and resamplings), so it never accepts alternative hypothesis. However, if I change the resampling to:

x.sim = sample(x, replace = T, size = nboot * length(x)) * matrix(sample(c(1, -1), replace = T, size = nboot * length(x)), nrow = nboot)


which is basically a non-centered wild bootstrap it works and converges to 5%.

I really wonder why the original non-parametric does not work in my case. I checked other possible error sources in my codes but didn't found any. The calculations of means, variances and test statistics are correct (checked manually, correct up to e-15). So it seems that the reason is actually the resampling of the values itself and I don't get why. If you have any ideas or suggestions for improvements, please answer. :)

Also: I am aware of packages for R such as boot, however it is not my intention to use them but to understand the coding of such bootstrapping itself (e.g. for courses at university).

Edit: Here is a MWE. (Here x.sim is named x.perm)

  t.sim = function(x, nboot = 1000, alpha = 0.05, mu = 0){
n = length(x)
x.perm = matrix(sample(x, size = nboot * n, replace = T), ncol = n, nrow = nboot, byrow = T) # Resampled Bootstrap Matrix

x.bar = rowMeans(x.perm)  # Means of Bootstrap Matrix
x.var = (rowSums(x.perm^2) - n * x.bar^2)/(n - 1) # Variance of the bootstrap matrix
#Bootstrapped t-statistics
T.star = (x.bar - mu)/sqrt(x.var / n)

#T-Value of original data
T.true = (mean(x) - mu)/(sqrt(var(x)/n))

#Compute Critical values
c.crit  = quantile(T.star, c(alpha/2, 1-alpha/2))
#Check whether the true t-Statistic lies outside those values
perm  = (T.true < c.crit[1] | T.true > c.crit[2])
return(perm)
}
# Check how often it returns true (thus rejects the null)
mean(replicate(1000, t.sim(rnorm(10))))
[1] 0


Please note that the only modification which achieved the correct convergence of type 1 error rate was the one with using rademacher weights. Different calculations in different languages didn't do the work.

• p values are based on sampling distributions under the null hypothesis. Your approach seems to be working with sampling distributions under the alternative hypothesis. Commented Apr 5, 2020 at 6:17
• My_fun uses mu_0 = 0 and rnorm(10) generates data points with mean = 0 and var = 1. So if you could explain where I am in the alternative and why one approach works and one not I'd be thankful. Identifying the specific problem is exactly where I'm stuck. Commented Apr 5, 2020 at 11:13
• Clearly something is wrong with your code--but without the details of what you're doing, all we can offer is speculation. Please include enough details to support your contention that the Type I error is zero.
– whuber
Commented Apr 5, 2020 at 14:51
• Edited. Provided a MWE. Commented Apr 5, 2020 at 15:01

I found the mistake, it lies in the bootstrapped distribution. In the above code the dist. is calculated by

T.star = (x.bar - mu)/sqrt(x.var / n)


Which is equal to $$T^{*} = \sqrt{n}\frac{\bar{X.}^{*} - \mu_0}{\hat{\sigma}^{*}}$$

However, the computation should be $$T^{*} = \sqrt{n}\frac{\bar{X.}^{*} - \mathbb{E}[\bar{X.}^{*}\mid X]}{\hat{\sigma}^{*}}$$ So, since $$\mathbb{E}[\bar{X.}^{*}\mid X] = \bar{x}$$, changing my code to

T.star = (x.bar - mean(x))/sqrt(x.var / n)


does the job and convergence is reached in type-1-error rate.