Is estimating confidence intervals (CI) with different sample sizes in each bootstrap valid?

Question

I am trying to estimate a confidence interval using bootstrapping. As R data.frame my data looks like

library(data.table)
df <- data.table(compound= c(rep("ala", 5), rep("beta", 3), rep("phe", 8)),
             obs = c(rep(FALSE, 7), rep(TRUE, 9)))

The statistic I am interested in is the percentage of TRUE values compared to the number of observations (9/16*100 = 56% for my example data). In my confidence interval I would like to account for the fact that my compounds were selected at random from a large number of compounds. Hence I would have intuitively done something like that (as written in R):

boot::boot.ci(boot::boot(data.frame(var = df$compound),
                         function(data, indices, stat_tab = df){
                       comp_samp <- data[indices,]

                       fin_tab <- 
                       lapply(as.list(comp_samp), function(x, stat_tab_l = stat_tab ){
                         stat_tab_l[x == compound]
                         })

                       fin_tab <- rbindlist(fin_tab )

                       round(nrow(fin_tab[obs == TRUE])/nrow(fin_tab )*100,1)
                     },
                     R = 1000),
          index=1,
          type='basic')$basic 

Is that a valid thing to do? I am a bit confused since my compounds can lead to different numbers of observations (rows in df) which means that in the different bootstrap samples I will have different numbers of observations when sampling by compound. In case it is not valid, why is that and is there a better way to estimate the CI in my scenario? Thank you

Answer 1

For a TRUE/FALSE outcome variable you should use logistic regression instead, and evaluate all the compounds at once in a single model. If you are primarily interested in the set of compounds you evaluated the model could be something like the following fixed-effects model (in R):

glm(obs ~ compound, family = binomial)

Here compound would be a multi-level categorical variable. One of the compounds would be specified as the reference; the intercept would be the log-odds of obs=TRUE for that compound. The regression coefficients for the other compounds would be the differences from that reference in log-odds. The standard errors reported for the intercept and regression coefficients provide (with some calculation) confidence intervals for the individual compounds. You would use standard post-hoc tests based on those coefficients and standard errors to examine differences among compounds.

If you instead want to model sampling of these specific compounds from a larger universe of compounds, you could consider instead a random-effects model. In R:

glmer(obs ~ (1|compound), family = binomial)

Then the intercept is an overall intercept for all the compounds, and the individual compounds in your sample are modeled with a Gaussian distribution of intercepts around that value. The model will report the variance among compounds around the intercept. But you won't get confidence intervals for the individual compounds this way. The results will of course depend on the particular compounds in your sample, and the quality of generalization would depend on the representativeness of your sample.

You certainly could in either case use bootstrapping in addition, which could be a good check on the quality of the model. Bootstrap from all of the cases. There will be different representation of the compounds among the bootstrap samples, but the total sample size (which is what matters) will be same for all. For the fixed-effects model, try modeling on a large number of bootstrap samples and see how well the bootstrap-derived models perform on the full original data set. For random-effects modeling you could see how stable the reported individual random effects were from bootstrap sample to bootstrap sample.