How to choose resample size when drawing without replacement?

Question

Say I have some second-order statistic $m(x)$ where $x$ is a data vector of length $n$. Let's also assume that the limiting distribution of $x$ is gaussian-ish, but generally unknown, so that the assumptions that enable one to derive the usual error analysis expressions do not hold. In this case, if I want to get an estimate in the uncertainty in the measure of $m$, I will have to simulate it, using a bootstrap or something. So, I generate 1,000 unique realizations of $x$, $x_i$ (1 $\leq$ $i$ $\leq$ 1000), and use the distribution of all $m(x_i)$ to get an idea of the error in $m$

Now, since $m$ is second-order, it is preferable to draw without replacement when generating the 1,000 resamples of $x$. This is all fine, and the bootstrapping routines I've implemented work well. Here is my problem:

• I have to choose a size of the resample
• If I choose that size to be $n$, then all $x_i$ will be identical since I'm sampling without replacement
• So, the resample size must obviously be smaller than the size of $x$

Problem is, if I choose the size of each $x_i$ to be high, say $0.9n$, then my error is going to be very small. If I choose the size to be small, $0.1n$, then the error can blow up. So, I can effectively make the error in $m$ whatever I'd like, which obviously isn't right...

What do I do at this point, while maintaining integrity?!

Answer 1

I'm posting my own answer here because I think it is probably correct, but feel free to please input ideas.

@MossMurderer brought my attention, in the comments to the original question, an interesting proof. That is, ${n \choose k}$ is maximum at $k=n/2$. This means that if I choose my resample size $k$ to be $n/2$, I will be maximizing the number of unique draws $x_i$.

This is desirable for a very simple reason: to preform the bootstrap, I compute the statistic $m$ for each realization $x_i$. Let's call the results of that calculation a vector $M = [m(x_1),... m(x_{1000})]$. The estimated error in the statistic $m$ is then something like $\sigma(M)$. Of course, you only ever use the standard deviation $\sigma$ as a statistic if you can safely assume your limiting distribution to be gaussian.

If we ensure that we maximize the number of unique $x_i$ by setting $k = n/2$, then we have sampled the $M$ space as well as we possibly can, and therefore $M$ will be as close to gaussian as is possible. In this way, the measure $\sigma(M)$ will be most reliable when the length $k$ of each $x_i$ is $n/2$.

I apologize if this is difficult to read - please edit if you feel you are more eloquent than I.

Answer 2

I included a bit of R code as a quick way to test intuitions about the bootstrap and variance estimation. For this data set, the bootstrap will tend to slightly underestimate the variance, consistent with the intuition that sample overlap (caused by drawing with replacement) will reduce the variance. I added an example with N / 2, to go along with your suggestion of maximizing the number of unique samples.

But the difference is very small, and with the bootstrap you don't have to justify the sample size you used for the resampling procedure, which takes less words and less citations to accomplish almost the same thing.

If you really don't want to draw with replacement and are keen on the N / 2 idea, you could take a look at Antal & Tillé's paper on doubled half bootstrap.

library("tidyverse")
library("data.table")

# get the iris data set
dt.iris <- iris %>% data.table

B = 2000
n_sample = nrow(dt.iris)

# look at the true variance for the dataset
var(dt.iris$Sepal.Length)  # 0.6856935
sd(dt.iris$Sepal.Length)  # 0.8280661

set.seed(123)

# bootstraped variance
dt.boot.results <- data.table(i = 1:B, boot_result_var = as.double(NA), boot_result_sd = as.double(NA))
for (i in 1:B) {
  cur.sample <- dt.iris %>% sample_n(n_sample, replace = TRUE)
  set(dt.boot.results, i = i, j = "boot_result_var", value = var(cur.sample$Sepal.Length))
  set(dt.boot.results, i = i, j = "boot_result_sd", value = sd(cur.sample$Sepal.Length))
}
summary(dt.boot.results$boot_result_var)  # mean var: 0.6801
summary(dt.boot.results$boot_result_sd)  # mean sd: 0.8237 

# resampled sd (half sample)
dt.res.results <- data.table(i = 1:B, res_result = as.double(NA), res_result_sd = as.double(NA))
for (i in 1:B) {
  cur.sample <- dt.iris %>% sample_n(n_sample / 2, replace = FALSE)
  set(dt.res.results, i = i, j = "res_result_var", value = var(cur.sample$Sepal.Length))
  set(dt.res.results, i = i, j = "res_result_sd", value = sd(cur.sample$Sepal.Length))
}
summary(dt.res.results$res_result_var)  # mean var: 0.6895
summary(dt.res.results$res_result_sd)  # mean sd: 0.8293