Chapter 3 Analysing model variability

3.0.1 Using Parametric and Non-Parametric Bootstrapping to Estimate Confidence Intervals in Meta-Analyses

Bootstrapping is a statistical technique used to estimate the sampling distribution of an estimator by repeatedly resampling from the observed data. In meta-analysis, bootstrapping is commonly applied to construct more robust confidence intervals (CIs) for parameters such as the overall effect size (μ) and the between-study variance (τ2^2τ2). This chapter will demonstrate how to implement both parametric and non-parametric bootstrapping methods using the boot package in R, emphasizing best practices for obtaining reliable results.

  • Parametric Bootstrapping: Assumes a specific distributional form for the data (e.g., normal distribution of residuals). The data are generated based on parameter estimates from the fitted model.

  • Non-Parametric Bootstrapping: Does not rely on specific distributional assumptions. The bootstrap samples are created directly by resampling the original data with replacement.

3.0.2 Parametric Bootstrapping: Step-by-Step Example

In parametric bootstrapping, we generate new datasets by simulating values from a specified distribution using the parameter estimates from the fitted model. This process requires defining two functions:

  1. A function to compute the statistics of interest (e.g., the mean effect size μand the between-study variance τ2^2τ2) based on the bootstrap data.

  2. A function to generate the bootstrap datasets.

Let’s illustrate this approach with a random-effects model:

Defining the Function

The function boot.func() calculates the effect size and variance components based on each bootstrap dataset:

# Load necessary libraries
library(metafor)
library(boot)


# 1. Fit the initial random-effects model
initial_model <- rma(yi, vi, data=dat)

# Extract estimated parameters for later use
mu_estimate <- coef(initial_model)
tau2_estimate <- initial_model$tau2

# 2. Define the Statistic Function for Bootstrapping
boot.func <- function(data.boot) {
  # Fit the random-effects model to the bootstrap data
  res <- try(suppressWarnings(rma(yi, vi, data=data.boot)), silent=TRUE)
  
  # Return NA if the model did not converge
  if (inherits(res, "try-error")) {
    return(rep(NA, 4))  # Return a vector of NAs
  } else {
    # Extract the estimated effect size (mu), its variance, tau², and its variance
    return(c(coef(res), diag(vcov(res)), res$tau2, res$se.tau2^2))
  }
}

# 3. Define the Data Generation Function for Bootstrapping
data.gen <- function(dat, mle) {
  # Generate effect sizes based on the estimated mu and tau²
  data.frame(yi = rnorm(nrow(dat), mle$mu, sqrt(mle$tau2 + dat$vi)), vi = dat$vi)
}

# 4. Running the Parametric Bootstrap
set.seed(1234)  # For reproducibility
res.boot <- boot::boot(dat, 
                        statistic = boot.func, 
                        R = 100, 
                        sim = "parametric", 
                        ran.gen = data.gen, 
                        mle = list(mu = mu_estimate, tau2 = tau2_estimate))

# Check results
print(res.boot)
## 
## PARAMETRIC BOOTSTRAP
## 
## 
## Call:
## boot::boot(data = dat, statistic = boot.func, R = 100, sim = "parametric", 
##     ran.gen = data.gen, mle = list(mu = mu_estimate, tau2 = tau2_estimate))
## 
## 
## Bootstrap Statistics :
##         original        bias     std. error
## t1* 2.552978e-01  1.086396e-03 1.908634e-02
## t2* 3.922816e-04 -9.255823e-06 5.149671e-05
## t3* 2.620568e-02 -7.892906e-04 4.629509e-03
## t4* 2.819226e-05 -8.527037e-07 7.784439e-06

Extracting Confidence Intervals

After running the bootstrap, we can calculate the confidence intervals for the mean effect size (μ) and the between-study variance (τ2^2τ2) using different bootstrap methods (normal, basic, studentized, percentile):

# Confidence intervals for mu
boot.ci(res.boot, type=c("norm", "basic", "stud", "perc"), index=1:2)
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 100 bootstrap replicates
## 
## CALL : 
## boot.ci(boot.out = res.boot, type = c("norm", "basic", "stud", 
##     "perc"), index = 1:2)
## 
## Intervals : 
## Level      Normal              Basic         
## 95%   ( 0.2168,  0.2916 )   ( 0.2137,  0.2933 )  
## 
## Level    Studentized          Percentile     
## 95%   ( 0.2123,  0.2956 )   ( 0.2173,  0.2969 )  
## Calculations and Intervals on Original Scale
## Some basic intervals may be unstable
## Some studentized intervals may be unstable
## Some percentile intervals may be unstable
# Confidence intervals for tau²
boot.ci(res.boot, type=c("norm", "basic", "stud", "perc"), index=3:4)
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 100 bootstrap replicates
## 
## CALL : 
## boot.ci(boot.out = res.boot, type = c("norm", "basic", "stud", 
##     "perc"), index = 3:4)
## 
## Intervals : 
## Level      Normal              Basic         
## 95%   ( 0.0179,  0.0361 )   ( 0.0160,  0.0373 )  
## 
## Level    Studentized          Percentile     
## 95%   ( 0.0184,  0.0432 )   ( 0.0152,  0.0364 )  
## Calculations and Intervals on Original Scale
## Some basic intervals may be unstable
## Some studentized intervals may be unstable
## Some percentile intervals may be unstable

These commands compute CIs based on various bootstrap methods, allowing comparison of the interval estimates.

3.0.3 Non-Parametric Bootstrapping: Step-by-Step Example

Non-parametric bootstrapping involves generating new datasets by resampling the original data with replacement. We only need to define a single function for this purpose:

Defining the Function

# Load necessary libraries
library(metafor)
library(boot)


# 1. Fit the initial random-effects model
initial_model <- rma(yi, vi, data=dat)

# Extract estimated parameters for later use
mu_estimate <- coef(initial_model)
tau2_estimate <- initial_model$tau2

# 2. Define the Statistic Function for Bootstrapping
boot.func <- function(data.boot, indices) {
   # Resample the data based on the given indices
   sel <- data.boot[indices, ]
   
   # Fit the random-effects model to the resampled data
   res <- try(suppressWarnings(rma(yi, vi, data=sel)), silent=TRUE)
   
   # Return NA if the model did not converge
   if (inherits(res, "try-error")) {
      return(rep(NA, 4))  # Return a vector of NAs
   } else {
      # Extract the estimated effect size (mu), its variance, tau², and its variance
      return(c(coef(res), diag(vcov(res)), res$tau2, res$se.tau2^2))
   }
}

# 3. Define the Data Generation Function for Bootstrapping
data.gen <- function(dat, mle) {
   # Generate effect sizes based on the estimated mu and tau²
   data.frame(yi = rnorm(nrow(dat), mle$mu, sqrt(mle$tau2 + dat$vi)), vi = dat$vi)
}

# 4. Running the Parametric Bootstrap
set.seed(1234)  # For reproducibility
res.boot <- boot::boot(dat, 
                        statistic = boot.func, 
                        R = 100, 
                        sim = "parametric", 
                        ran.gen = data.gen, 
                        mle = list(mu = mu_estimate, tau2 = tau2_estimate))

# Check results
print(res.boot)
## 
## PARAMETRIC BOOTSTRAP
## 
## 
## Call:
## boot::boot(data = dat, statistic = boot.func, R = 100, sim = "parametric", 
##     ran.gen = data.gen, mle = list(mu = mu_estimate, tau2 = tau2_estimate))
## 
## 
## Bootstrap Statistics :
##         original        bias     std. error
## t1* 2.552978e-01  1.086396e-03 1.908634e-02
## t2* 3.922816e-04 -9.255823e-06 5.149671e-05
## t3* 2.620568e-02 -7.892906e-04 4.629509e-03
## t4* 2.819226e-05 -8.527037e-07 7.784439e-06

Extracting Confidence Intervals

After running the bootstrap, we can calculate the confidence intervals for the mean effect size (μ) and the between-study variance (τ2^2τ2) using different bootstrap methods (normal, basic, studentized, percentile):

# Confidence intervals for mu
boot.ci(res.boot, type=c("norm", "basic", "stud", "perc"), index=1:2)
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 100 bootstrap replicates
## 
## CALL : 
## boot.ci(boot.out = res.boot, type = c("norm", "basic", "stud", 
##     "perc"), index = 1:2)
## 
## Intervals : 
## Level      Normal              Basic         
## 95%   ( 0.2168,  0.2916 )   ( 0.2137,  0.2933 )  
## 
## Level    Studentized          Percentile     
## 95%   ( 0.2123,  0.2956 )   ( 0.2173,  0.2969 )  
## Calculations and Intervals on Original Scale
## Some basic intervals may be unstable
## Some studentized intervals may be unstable
## Some percentile intervals may be unstable
# Confidence intervals for tau²
boot.ci(res.boot, type=c("norm", "basic", "stud", "perc"), index=3:4)
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 100 bootstrap replicates
## 
## CALL : 
## boot.ci(boot.out = res.boot, type = c("norm", "basic", "stud", 
##     "perc"), index = 3:4)
## 
## Intervals : 
## Level      Normal              Basic         
## 95%   ( 0.0179,  0.0361 )   ( 0.0160,  0.0373 )  
## 
## Level    Studentized          Percentile     
## 95%   ( 0.0184,  0.0432 )   ( 0.0152,  0.0364 )  
## Calculations and Intervals on Original Scale
## Some basic intervals may be unstable
## Some studentized intervals may be unstable
## Some percentile intervals may be unstable

Comparing Bootstrap Methods

The choice between parametric and non-parametric bootstrapping depends on the underlying assumptions and the sample size:

  • Parametric Bootstrapping is more appropriate when we have strong assumptions about the distribution of effect sizes.

  • Non-Parametric Bootstrapping is recommended when the sample size is relatively small or when we want to avoid making strict assumptions.

Visualization of Bootstrap Distributions

To visualize the bootstrap distributions, we can create kernel density plots to inspect the variability and shape of the bootstrap samples:

# Visualize the bootstrap distribution for mu
plot(density(res.boot$t[,1]), main="Bootstrap Distribution of Mu")
abline(v=quantile(res.boot$t[,1], probs=c(0.025, 0.975)), col="red")

3.0.4 Conclusion

Bootstrapping provides a flexible and powerful method for estimating confidence intervals in meta-analysis. However, it is important to consider potential issues such as non-convergence of models during the bootstrap process, as well as the assumptions underlying each method. When applying bootstrap methods, carefully evaluate the coverage of different interval types and compare the results to standard methods (e.g., Wald-type CIs or Knapp-Hartung adjustments).

In the next section, we will delve deeper into prediction intervals, which provide a complementary measure of uncertainty, reflecting the variability in effect sizes for new studies.