Subsampling to determine a standard error, how does it work?

Question

I need to calculate the standard error on a complicated dataset (> 1700 records) which uses genetic matching. Using bootstrap results in very high computation time (because of the genetic matching).

My professor gives subsampling as an alternative:

An alternative which is less computational expensive is subsampling. You can perform this by selecting for example 100 smaller subsamples (consisting of like 500 measured subjects) from the big dataset and analysing these.

The empirical variance of these 100 estimators will be too high because you analyzed smaller datasets. But because the variance is inversely proportional to the number of measured subjects one can calculate the variance of the estimator in the complete dataset.

I seem to get stuck in the calculation of the inversely proportion. I tried this on a different (simpler) dataset in R:

library(boot)
x <- rnorm(5000,mean=0,sd=1)
##then the SE for the mean is 1/sqrt(5000)
se <- 1/sqrt(5000); se

mean.fn <- function(data, index){
  return(mean(data[index]))
}

# bootstrapping works good
boot(x,mean.fn,R=100)

## subsampling
res <- c()
count <- 100
size <- 500
for(i in 1:count){
  indices <- sample(x = length(x),size=size)
  res[i] <- mean.fn(x,indices)
}

emp_se <- var(res)
emp_var <- emp_se*sqrt(count)
emp_se; emp_var; se

Somehow emp_var seems to be a good estimator for se, but when I change the count to 200 this gets worse. How can I turn var(res) into an estimator for se?

Answer 1

There are a few problems in your code.

First of all, why use var(res) if you are interested in the standard deviation of the mean? Instead just use sd(res) (or equivalently sqrt(var(res))). Also your Professor of course meant that the errors scales according to sample and subsample size, not the number of your Monte Carlo replicates !!

So in your case you would get a reasonable estimator by:

mean_sd <- sd(res)*sqrt(size/5000)