# Bootstrapping mean difference: standard error versus quantiles

I am implementing a bootstrap procedure (in R) to calculate the confidence interval of a difference of two means.

I have little experience with bootstrap but I am aware of two methods:

• calculate the standard deviation of the bootstrap sample and construct the CI accordingly
• get the CI from the quantiles of the bootstrap sample

Here is some commented R code for 500 replications:

x=mtcars$$mpg by=ifelse(mtcars$$vs==0, "vshaped", "straight")
R=500
ref = "vshaped"
conf_level=0.95

#effect: difference in means
effect = mean(x[by!=ref], na.rm = TRUE) - mean(x[by==ref], na.rm = TRUE)

#bootstrap effect:
beffect = numeric(length=R)                   #allocate a vector of length 500
for (i in 1:R) {                              #loop over 500 times
ib = sample(1:length(x), replace = TRUE)  #get sample row numbers
xi = x[ib]                                #select samples in x
gi = by[ib]                               #select samples in by
#allocate the effect
beffect[i] = mean(xi[gi!=ref], na.rm = TRUE) - mean(xi[gi==ref], na.rm = TRUE)
}

#method 1: standard error
sd.effect = sd(beffect, na.rm=TRUE)
effect + qnorm(c((1-conf_level)/2, 1-(1-conf_level)/2))*sd.effect
#> [1]  4.773336 11.107617

#method 2: quantiles
quantile(beffect, c(0.025,0.975))
#>      2.5%     97.5%
#>  4.827105 10.952341


Created on 2021-03-22 by the reprex package (v1.0.0)

As you can see, the results are slightly different, which is somehow expected.

However, I could not find any resources about the difference between these procedures.

Is one better suited in some cases? What hypothesis do they imply each?

Your bootstrapped beffect is drawn from some normal distribution:

hist(beffect)
qqnorm(beffect)
qqline(beffect)


So your question is: I have a large number of samples from a normal distribution, how to determine the quantiles of the underlying distribution. In any case, you can only estimate these quantiles, never determine them. When faced with the task to produce the best possible estimates the most important thing to do is investigate many more cars. Man these are 32 cars form 1974! If you want to know about cars, sample more cars and maybe more up-to-date ones. That is the most important aspect: get as much data as you can. Then the next step is to draw more then 500 bootstrap samples Sample a million while you are sitting at your desk or ten million while you pour a cup of coffee. That will bring the two estimates close together. Investing the time into these two steps will make far more difference then thinking to long about whether or not a parametric or a non-parametric bootstrap will bring a tiny advantage in precision.

Whilst typing this I ran you code with R = 1e7 in the background. The result was

> effect + qnorm(c((1-conf_level)/2, 1-(1-conf_level)/2))*sd.effect
[1]  4.655487 11.225466
>
> quantile(beffect, c(0.025,0.975))
2.5%     97.5%
4.701619 11.266667


Compare the precision of the lower and upper borders with the width of the confidence interval and you will find, that most of the time it really does not matter. Remove or add just one car and see the influence that makes in comparison.

• Thanks, this is very interesting! So basically, both methods are the same when R is high? But when it's not, what hypothesis should I consider? – Dan Chaltiel Mar 22 at 14:18
• They will only be the same to the extent that the bootstrapped samples are normally distributed. This may not be the case with some statistics (eg., medians) and/or with small sample sizes. Using the search terms "nonparametric vs parametric bootstrap confidence intervals" and "bias-corrected and accelerated bootstrap (BCa) confidence intervals" you'll find some research on the topig, eg. frontiersin.org/articles/10.3389/fpsyg.2019.02215/full . For practical purposes it may be a good idea to use the boot package to compute all three types of confidence ... – Bernhard Mar 22 at 14:52
• ... intervals and see, if they come to the same conclusion. If they do, it really does not matter. If they don't, you are in a special situation that needs closer inspection in detail. My advice is to start reading here: blog.methodsconsultants.com/posts/… (Understanding Bootstrap Confidence Interval Output from the R boot Package, Jeremy Albright, posted on Sep 30, 2019) – Bernhard Mar 22 at 14:54
• Thanks, I'll look into this. I have only one bootstrap to do and I don't want to add a dependency to my package. Moreover, I find the boot package very tedious to use but that's a matter of taste I guess. – Dan Chaltiel Mar 22 at 14:56
• It is generally thoughtfull to avoid unneccessary dependencies. What I personally like most about bootstrapping is not having to rely on distribution assumptions so it feels more natural for me to use quantiles of the bootstrapped samples. That is a feeling, not a decision made in full awareness of all the available literature. – Bernhard Mar 23 at 7:10