Show that FE model has nominal coverage for conditional inferences (fixed theta scenario).

fe_coverage_fixed_theta

rm(list=ls())

library(metafor)

iters <- 1000

### number of studies/effects
k <- 10

### SD of the true effects
tau <- 0.2

### generate heterogeneous true effects whose average is 0
thetai <- rnorm(k, 0, tau)
thetai <- thetai - mean(thetai)

cov.u <- rep(NA, iters)

for (j in 1:iters) {

   ### simulate estimates with heterogeneous true effects
   vi <- runif(k, .001, 1)
   yi <- rnorm(k, thetai, sqrt(vi))

   ### estimate unweighted mean of thetai values
   res <- rma(yi, vi, method="FE", weighted=FALSE)

   ### check for coverage of the unweighted mean
   theta.mean <- mean(thetai) # 0 by construction
   cov.u[j] <- (res$ci.lb <= theta.mean) && (res$ci.ub >= theta.mean)

}

### print coverage rate
print(mean(cov.u))