Random effects meta-analysis for correlated test statistics

Icahn School of Medicine at Mount Sinai, New York

Standard approaches to meta-analysis assumes that effect sizes are statistically independent. Here we provide methods for fixed and random effects meta-analysis when the correlation between effect sizes are known.

Fixed effects meta-analysis

LS() implements fixed effect meta-analysis for correlated test statistics using method of Lin and Sullivan (2009). By default, correlation is set to identity matrix to for independent test statistics.

Random effects meta-analysis

RE2C() implements random effect meta-analysis for correlated test statistics that jointly tests deviation of the mean from zero as well as effect size heterogenity. This method uses the RE2 method of Han and Eskin (2011), or RE2 for correlated test statistics from Han et al. (2016). By default, correlation is set to identity matrix to for independent test statistics. (In addition, this function computes the two step RE2C method of Lee, Eskin, and Han (2017) to further test for heterogenity in effect size after applying a fixed effect test.)

• stat1: statistic testing effect mean

• stat2: statistic testing effect heterogeneity

• RE2Cp: RE2 p-value accounting for correlelation between tests. (This is the p-value appropriate for most questions)

• RE2Cp.twoStep: two step RE2C test after fixed effect test. Only evaluated if twoStep==TRUE. (not typically used)

• QE: test statistic for the test of (residual) heterogeneity

• QEp: p-value for the test of (residual) heterogeneity

• Isq: I^2 statistic

  QE, QEp and Isq are only evaluted if correlation is diagonal

Examples

library(remaCor)
library(metafor)
library(mvtnorm)
library(clusterGeneration )

# sample size
n = 30

# number of response variables
m = 2

# Error covariance
Sigma = genPositiveDefMat(m)$Sigma

# regression parameters
beta = matrix(0, 1, m)

# covariates
X = matrix(rnorm(n), ncol=1)

# Simulate response variables
Y = X %*% beta + rmvnorm(n, sigma = Sigma)

# Multivariate regression
fit = lm(Y ~ X)

# Correlation between residuals
C = cor(residuals(fit))

# Extract effect sizes and standard errors from model fit
df = lapply(coef(summary(fit)), function(a) 
  data.frame(beta = a["X", 1], se = a["X", 2]))
df = do.call(rbind, df)

# Standard fixed effects meta-analysis
# of independent effects with metafor pacakge
rma( df$beta, sei=df$se, method="FE")

## 
## Fixed-Effects Model (k = 2)
## 
## I^2 (total heterogeneity / total variability):   62.83%
## H^2 (total variability / sampling variability):  2.69
## 
## Test for Heterogeneity:
## Q(df = 1) = 2.6901, p-val = 0.1010
## 
## Model Results:
## 
## estimate      se     zval    pval    ci.lb   ci.ub    
##  -0.0304  0.2675  -0.1137  0.9095  -0.5547  0.4939    
## 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# Standard random effects meta-analysis
# of independent effects with metafor pacakge
rma( df$beta, sei=df$se, method="REML")

## 
## Random-Effects Model (k = 2; tau^2 estimator: REML)
## 
## tau^2 (estimated amount of total heterogeneity): 0.2449 (SE = 0.5513)
## tau (square root of estimated tau^2 value):      0.4949
## I^2 (total heterogeneity / total variability):   62.83%
## H^2 (total variability / sampling variability):  2.69
## 
## Test for Heterogeneity:
## Q(df = 1) = 2.6901, p-val = 0.1010
## 
## Model Results:
## 
## estimate      se     zval    pval    ci.lb   ci.ub    
##  -0.0611  0.4411  -0.1384  0.8899  -0.9256  0.8035    
## 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# Run fixed effects meta-analysis, assume identity correlation  
# Use Lin-Sullivan method
LS( df$beta, df$se)

##          beta        se         p
## 1 -0.03041698 0.2675163 0.9094745

# Run fixed effects meta-analysis, accounting for correlation  
# Use Lin-Sullivan method
LS( df$beta, df$se, C)

##          beta        se         p
## 1 -0.01968027 0.2905851 0.9460035

# Run random effects meta-analysis, assume identity correlation  
RE2C( df$beta, df$se)

##        stat1      stat2     RE2Cp RE2Cp.twoStep       QE       QEp      Isq
## 1 0.01292801 0.09177158 0.7780907            NA 2.690102 0.1009734 62.82668

# Run random effects meta-analysis, accounting for correlation 
RE2C( df$beta, df$se, C)

##         stat1     stat2     RE2Cp RE2Cp.twoStep QE QEp Isq
## 1 0.004586858 0.4210307 0.5708562            NA NA  NA  NA

##         stat1     stat2     RE2Cp RE2Cp.twoStep QE QEp Isq
## 1 0.004586858 0.4210307 0.5708562            NA NA  NA  NA

Session info

sessionInfo()

## R version 4.3.0 (2023-04-21)
## Platform: x86_64-apple-darwin22.4.0 (64-bit)
## Running under: macOS 14.2.1
## 
## Matrix products: default
## BLAS:   /Users/gabrielhoffman/prog/R-4.3.0/lib/libRblas.dylib 
## LAPACK: /usr/local/Cellar/r/4.3.0_1/lib/R/lib/libRlapack.dylib;  LAPACK version 3.11.0
## 
## locale:
## [1] C/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## time zone: America/New_York
## tzcode source: internal
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
## [1] clusterGeneration_1.3.7 MASS_7.3-60             mvtnorm_1.2-2          
## [4] metafor_4.2-0           numDeriv_2016.8-1.1     metadat_1.2-0          
## [7] Matrix_1.5-4.1          remaCor_0.0.18          ggplot2_3.4.4          
## 
## loaded via a namespace (and not attached):
##  [1] sass_0.4.6       utf8_1.2.3       generics_0.1.3   EnvStats_2.7.0  
##  [5] stringi_1.7.12   lattice_0.21-8   digest_0.6.33    magrittr_2.0.3  
##  [9] evaluate_0.21    grid_4.3.0       fastmap_1.1.1    plyr_1.8.8      
## [13] jsonlite_1.8.5   fansi_1.0.4      scales_1.2.1     codetools_0.2-19
## [17] jquerylib_0.1.4  Rdpack_2.4       cli_3.6.1        rlang_1.1.1     
## [21] rbibutils_2.2.13 munsell_0.5.0    withr_2.5.0      cachem_1.0.8    
## [25] yaml_2.3.7       tools_4.3.0      reshape2_1.4.4   dplyr_1.1.2     
## [29] colorspace_2.1-0 mathjaxr_1.6-0   vctrs_0.6.3      R6_2.5.1        
## [33] lifecycle_1.0.3  stringr_1.5.0    pkgconfig_2.0.3  pillar_1.9.0    
## [37] bslib_0.4.2      gtable_0.3.3     glue_1.6.2       Rcpp_1.0.11     
## [41] xfun_0.39        tibble_3.2.1     tidyselect_1.2.0 knitr_1.43      
## [45] htmltools_0.5.5  nlme_3.1-162     rmarkdown_2.22   compiler_4.3.0

References

Han, Buhm, Dat Duong, Jae Hoon Sul, Paul IW de Bakker, Eleazar Eskin, and Soumya Raychaudhuri. 2016. “A General Framework for Meta-Analyzing Dependent Studies with Overlapping Subjects in Association Mapping.” Human Molecular Genetics 25 (9): 1857–66. https://doi.org/10.1093/hmg/ddw049.

Han, Buhm, and Eleazar Eskin. 2011. “Random-Effects Model Aimed at Discovering Associations in Meta-Analysis of Genome-Wide Association Studies.” The American Journal of Human Genetics 88 (5): 586–98. https://doi.org/10.1016/j.ajhg.2011.04.014.

Lee, CH, Eleazar Eskin, and Buhm Han. 2017. “Increasing the Power of Meta-Analysis of Genome-Wide Association Studies to Detect Heterogeneous Effects.” Bioinformatics 33 (14): i379–88. https://doi.org/10.1093/bioinformatics/btx242.

Lin, Dan-Yu, and Patrick F Sullivan. 2009. “Meta-Analysis of Genome-Wide Association Studies with Overlapping Subjects.” The American Journal of Human Genetics 85 (6): 862–72. https://doi.org/10.1016/j.ajhg.2009.11.001.