Table of Contents
Berkey et al. (1995)
The Methods and Data
Berkey et al. (1995) describe the meta-analytic random- and mixed-effects models and provide the equation for the empirical Bayes estimator for the amount of (residual) heterogeneity (p. 398). The models and methods are illustrated with the BCG vaccine dataset (Colditz et al., 1994). The data are provided in Table 1 in the article and can be loaded with:
library(metafor) dat.bcg
The contents of the dataset are:
trial author year tpos tneg cpos cneg ablat alloc 1 1 Aronson 1948 4 119 11 128 44 random 2 2 Ferguson & Simes 1949 6 300 29 274 55 random 3 3 Rosenthal et al 1960 3 228 11 209 42 random 4 4 Hart & Sutherland 1977 62 13536 248 12619 52 random 5 5 Frimodt-Moller et al 1973 33 5036 47 5761 13 alternate 6 6 Stein & Aronson 1953 180 1361 372 1079 44 alternate 7 7 Vandiviere et al 1973 8 2537 10 619 19 random 8 8 TPT Madras 1980 505 87886 499 87892 13 random 9 9 Coetzee & Berjak 1968 29 7470 45 7232 27 random 10 10 Rosenthal et al 1961 17 1699 65 1600 42 systematic 11 11 Comstock et al 1974 186 50448 141 27197 18 systematic 12 12 Comstock & Webster 1969 5 2493 3 2338 33 systematic 13 13 Comstock et al 1976 27 16886 29 17825 33 systematic
Next, we can calculate the log risk ratios and corresponding sampling variances with:
dat <- escalc(measure="RR", ai=tpos, bi=tneg, ci=cpos, di=cneg, data=dat.bcg)
The metafor package calculates the sampling variances of the log risk ratios as described by Berkey et al. (1995) on page 399 (right under the 2×2 table). However, note that Berkey et al. (1995) actually use 'smoothed' estimates of the sampling variances for their analyses (to reduce the correlation between the log risk ratios and corresponding sampling variances). We can compute these estimates with:
dat$vi <- with(dat, mean(tneg/tpos)/(tneg+tpos) + mean(cneg/cpos)/(cneg+cpos))
Random-Effects Model
Now we can fit a random-effects model to these data (using the empirical Bayes estimator for the amount of heterogeneity) with:
res.RE <- rma(yi, vi, data=dat, method="EB") res.RE
Random-Effects Model (k = 13; tau^2 estimator: EB) tau^2 (estimated amount of total heterogeneity): 0.2682 (SE = 0.1801) tau (square root of estimated tau^2 value): 0.5178 I^2 (total heterogeneity / total variability): 87.49% H^2 (total variability / sampling variability): 7.99 Test for Heterogeneity: Q(df = 12) = 85.8625, p-val < .0001 Model Results: estimate se zval pval ci.lb ci.ub -0.5429 0.1842 -2.9474 0.0032 -0.9040 -0.1819 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
These results match exactly what Berkey et al. (1995) report on page 408: The amount of heterogeneity ('between-trial variance') is estimated to be $\hat{\tau}^2 = 0.268$ and the pooled estimate is $\hat{\mu} = -0.5429$ with a standard error of $SE[\hat{\mu}] = 0.1842$.
Mixed-Effects Meta-Regression Model
A mixed-effects meta-regression model with absolute latitude as moderators (centered at 33.46 degrees) can be fitted with:
res.ME <- rma(yi, vi, mods=~I(ablat-33.46), data=dat, method="EB") res.ME
Mixed-Effects Model (k = 13; tau^2 estimator: EB) tau^2 (estimated amount of residual heterogeneity): 0.1572 (SE = 0.1262) tau (square root of estimated tau^2 value): 0.3965 I^2 (residual heterogeneity / unaccounted variability): 77.42% H^2 (unaccounted variability / sampling variability): 4.43 R^2 (amount of heterogeneity accounted for): 41.38% Test for Residual Heterogeneity: QE(df = 11) = 35.8827, p-val = 0.0002 Test of Moderators (coefficient 2): QM(df = 1) = 5.9018, p-val = 0.0151 Model Results: estimate se zval pval ci.lb ci.ub intrcpt -0.6303 0.1591 -3.9630 <.0001 -0.9421 -0.3186 *** I(ablat - 33.46) -0.0268 0.0110 -2.4294 0.0151 -0.0484 -0.0052 * --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Again, the results match the findings from Berkey et al. (1995): The residual amount of heterogeneity is now $\hat{\tau}^2 = 0.157$ and the estimated model is $log(RR) = -0.6303 - 0.0268 (x - 33.46)$, where $x$ is the distance from the equator (in degrees latitude). The standard errors of the model coefficients are $SE[b_0] = 0.1591$ and $SE[b_1] = 0.0110$.
The amount of variance (heterogeneity) accounted for by the absolute latitude moderator is provided in the output above. It can also be obtained with:1)
anova(res.RE, res.ME)
df AIC BIC AICc logLik LRT pval QE tau^2 R^2 Full 3 29.3011 30.9960 31.9678 -11.6506 35.8827 0.1572 Reduced 2 32.2578 33.3877 33.4578 -14.1289 4.9566 0.0260 85.8625 0.2682 41.3844%
The value below R^2
indicates that approximately 41% of the heterogeneity has been accounted for.
The predicted average risk ratios at 33.46 and 42 degrees reported by Berkey et al. (1995) can be computed with:
predict(res.ME, newmods=c(33.46,42)-33.46, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 1 0.53 0.39 0.73 0.23 1.23 2 0.42 0.28 0.64 0.18 1.02
Fixed-Effects Meta-Regression Model
The results from a fixed-effects meta-regression model with absolute latitude as the predictor can be obtained with:
res.FE <- rma(yi, vi, mods=~I(ablat-33.46), data=dat, method="FE") res.FE
Fixed-Effects with Moderators Model (k = 13) I^2 (residual heterogeneity / unaccounted variability): 69.34% H^2 (unaccounted variability / sampling variability): 3.26 Test for Residual Heterogeneity: QE(df = 11) = 35.8827, p-val = 0.0002 Test of Moderators (coefficient 2): QM(df = 1) = 49.9798, p-val < .0001 Model Results: estimate se zval pval ci.lb ci.ub intrcpt -0.5949 0.0696 -8.5462 <.0001 -0.7314 -0.4585 *** I(ablat - 33.46) -0.0282 0.0040 -7.0696 <.0001 -0.0360 -0.0204 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
The predicted risk ratios at 33.46 and 42 degrees are now:
predict(res.FE, newmods=c(33.46,42)-33.46, transf=exp, digits=2)
pred ci.lb ci.ub 1 0.55 0.48 0.63 2 0.43 0.36 0.52
References
Berkey, C. S., Hoaglin, D. C., Mosteller, F., & Colditz, G. A. (1995). A random-effects regression model for meta-analysis. Statistics in Medicine, 14(4), 395–411.
Colditz, G. A., Brewer, T. F., Berkey, C. S., Wilson, M. E., Burdick, E., Fineberg, H. V., et al. (1994). Efficacy of BCG vaccine in the prevention of tuberculosis: Meta-analysis of the published literature. Journal of the American Medical Association, 271(9), 698–702.