# Bootstrapping clusters in R

I am running a negative binomial regression of clinic counts in each county in the entire country (~3k counties). I'd like to at least partially account for the non-independence of neighboring counties by bootstrapping the confidence intervals in a "clustered" fashion--e.g. draw an entire state's (50 states total) worth of data at once. This has become standard practice, for better or for worse, in the econometric literature.

I could write the code to do this myself, but the boot package seems like it should have the ability to do this somehow, and in general I prefer tested, general solutions to one-off hacks. Is there a way to coerce the boot package to do a clustered bootstrap?

I tried the strata argument, but that randomizes within strata rather than randomizing which cluster gets taken, as the following code confirms:

dat <- data.frame( cluster=rep(letters[1:5],each=10), x=runif(5*10), stringsAsFactors=TRUE )
boot.stat <- function(dat,idx) {
print(dat[idx,]$cluster) print(table(dat[idx,]$cluster))
mean(dat[idx,]$x) } boot( data=dat, statistic=boot.stat, strata=dat$cluster,
stype="i",
R=5
)

• Why would an unconditional bootstrap be inappropriate for your purposes? It would generate a random number of observations from each cluster and there would be no permutation of cluster labels from bootstrap samples. – AdamO May 8 '13 at 19:19
• @AdamO Because it destroys the within-cluster correlation. Maintaining that correlation is the logic behind sampling the entire cluster. – Ari B. Friedman May 8 '13 at 19:22
• Does it? If values were constant within each cluster, they would still be constant within each cluster for any unconditional bootstrap resampling of the data. I think it would maintain cluster level correlation. That wouldn't be the case with permutation testing, of course. – AdamO May 8 '13 at 21:02
• @AdamO I think it does but I have been crazy recently. I'll try to come up with a simulation demonstrating soon, but regardless it's the method I'd like to replicate. – Ari B. Friedman May 9 '13 at 19:56
• @AdamO if I'm not mistaken it's the same logic behind time series bootstraps--drawing blocks to account for serial correlation. And Ari, that link is dead, do you by any chance recall the paper you had linked? I could go for some academic references on this topic. – MichaelChirico Sep 27 '15 at 11:57

If I understand you correctly you want to estimate a statistic per state and that average that statistic to get a bootstrapped estimation of the overall statistic.

Stratified sampling does something different. It ensures that the label is samples representatively in each sample. I do not think that is what you want to do.

You could do this manually without being hacky. Using the dplyr, tidyr and purrr package from the tidyverse this becomes transparant and clean code.

library(tidyr)
library(dplyr)
library(purrr)

dat <- data.frame(cluster=rep(letters[1:5],each=10),
x=runif(5*10), stringsAsFactors=TRUE)

boot.stat2 <- function(df) {
mean(df\$x)
}

dat %>%
nest(x) %>%
mutate(stat = map_dbl(data, boot.stat2))