Table of Contents
Stijnen et al. (2010)
The Methods and Data
Dichotomous and event count data (based on which one can calculate effect size or outcome measures, such as odds ratios, incidence rate ratios, proportions, and incidence rates) are often assumed to arise from binomial and Poisson distributed data. Stijnen et al. (2010) describe meta-analytic models that are directly based on these distibutions. The collection of models are essentially special cases of generalized linear (mixed-effects) models (i.e., mixed-effects logistic and Poisson regression models). For 2×2 table data, a mixed-effects conditional logistic model (based on the non-central hypergeometric distribution) can also be used. Mixed-effects logistic regression models for dichotomous data are often referred to as "binomial-normal models" in the meta-analytic literature. Analogously, mixed-effects Poisson regression models for event count data could be referred to as "Poisson-normal models". The results obtained from such models can be contrasted with those obtained from the standard (inverse-variance) approach (sometimes called the "normal-normal model"). The various models described in the article are illustrated with datasets from two meta-analyses comparing the risk of catheter-related bloodstream infection (CRBSI) when using anti-infective-treated versus standard catheters (Niel-Weise et al., 2007, 2008).
Binomial-Normal Model for the Meta-Analysis of Proportions
Niel-Weise et al. (2007) conducted a meta-analysis comparing the CRBSI risk when using anti-infective-treated versus standard catheters in the acute care setting. These data can be used to illustrate the binomial-normal model for the meta-analysis of proportions (i.e., infection risks) within the standard catheter and anti-infective catheter groups. The data can be extracted from Figure 2 (p. 2063) in Niel-Weise et al. (2007) and can be loaded with:
library(metafor) dat <- dat.nielweise2007 dat
(I copy the dataset into dat
, which is a bit shorter and therefore easier to type further below). The contents of the dataset are:
study author year ai n1i ci n2i 1 1 Bach 1996 0 116 3 117 2 2 George 1997 1 44 3 35 3 3 Maki 1997 2 208 9 195 4 4 Raad 1997 0 130 7 136 5 5 Heard 1998 5 151 6 157 6 6 Collin 1999 1 98 4 139 7 7 Hannan 1999 1 174 3 177 8 8 Marik 1999 1 74 2 39 9 9 Pierce 2000 1 97 19 103 10 10 Sheng 2000 1 113 2 122 11 11 Chatzinikolaou 2003 0 66 7 64 12 12 Corral 2003 0 70 1 58 13 13 Brun-Buisson 2004 3 188 5 175 14 14 Leon 2004 6 187 11 180 15 15 Yucel 2004 0 118 0 105 16 16 Moretti 2005 0 252 1 262 17 17 Rupp 2005 1 345 3 362 18 18 Osma 2006 4 64 1 69
Variables n1i
and n2i
denote the number of patients receiving an anti-infective or standard catheter, respectively, and ai
and ci
the number of CRBSIs in the respective groups. We can now fit a "normal-normal model" to the data for the two groups separately. Here, we use the logit transformed proportions (i.e., log odds) for the meta-analysis. First, for the standard catheter groups, this can be done with:
res <- rma(measure="PLO", xi=ci, ni=n2i, data=dat) print(res, digits=3)
Random-Effects Model (k = 18; tau^2 estimator: REML) tau^2 (estimated amount of total heterogeneity): 0.663 (SE = 0.338) tau (square root of estimated tau^2 value): 0.814 I^2 (total heterogeneity / total variability): 74.01% H^2 (total variability / sampling variability): 3.85 Test for Heterogeneity: Q(df = 17) = 72.166, p-val < .001 Model Results: estimate se zval pval ci.lb ci.ub -3.302 0.238 -13.883 <.001 -3.768 -2.836 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
The estimated average log odds (i.e., $\hat{\mu} = -3.302$) can be back-transformed using the inverse logit transformation with:
predict(res, transf=transf.ilogit, digits=3)
pred ci.lb ci.ub pi.lb pi.ub 0.036 0.023 0.055 0.007 0.163
As pointed out in the article, the proper interpretation of the back-transformed value (i.e., $0.036$) is that it reflects the median infection risk (although it would often be interpreted as the average risk).
Next, for the anti-infective catheter groups, we obtain:
res <- rma(measure="PLO", xi=ai, ni=n1i, data=dat) print(res, digits=3)
Random-Effects Model (k = 18; tau^2 estimator: REML) tau^2 (estimated amount of total heterogeneity): 0.393 (SE = 0.374) tau (square root of estimated tau^2 value): 0.627 I^2 (total heterogeneity / total variability): 37.60% H^2 (total variability / sampling variability): 1.60 Test for Heterogeneity: Q(df = 17) = 23.671, p-val = 0.129 Model Results: estimate se zval pval ci.lb ci.ub -4.260 0.259 -16.457 <.001 -4.768 -3.753 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
And again, using the back-transformation:
predict(res, transf=transf.ilogit, digits=3)
pred ci.lb ci.ub pi.lb pi.ub 0.014 0.008 0.023 0.004 0.051
Note that REML estimation is the default for the rma()
function. Also, calculation of the log odds (and the corresponding sampling variance) is problematic when a group has zero events, as is the case in some of the studies. The standard continuity correction (adding 1/2 to the number of patients with and without an event) is automatically applied inside the function. These results match what is reported by Stijnen et al. (2010) in Table II (p. 3050).
Instead of the normal-normal model, we can use a "binomial-normal model", which in this case is just a logistic regression model with a random intercept. This model can be fitted to the data for the standard catheter groups with:
res <- rma.glmm(measure="PLO", xi=ci, ni=n2i, data=dat) print(res, digits=3)
Random-Effects Model (k = 18; tau^2 estimator: ML) tau^2 (estimated amount of total heterogeneity): 0.812 tau (square root of estimated tau^2 value): 0.901 I^2 (total heterogeneity / total variability): 77.73% H^2 (total variability / sampling variability): 4.49 Tests for Heterogeneity: Wld(df = 17) = 69.103, p-val < .001 LRT(df = 17) = 85.284, p-val < .001 Model Results: estimate se zval pval ci.lb ci.ub -3.496 0.257 -13.605 <.001 -4.000 -2.993 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
predict(res, transf=transf.ilogit, digits=3)
pred ci.lb ci.ub pi.lb pi.ub 0.029 0.018 0.048 0.005 0.160
And again, for the anti-infective catheter groups:
res <- rma.glmm(measure="PLO", xi=ai, ni=n1i, data=dat) print(res, digits=3)
Random-Effects Model (k = 18; tau^2 estimator: ML) tau^2 (estimated amount of total heterogeneity): 0.827 tau (square root of estimated tau^2 value): 0.909 I^2 (total heterogeneity / total variability): 55.91% H^2 (total variability / sampling variability): 2.27 Tests for Heterogeneity: Wld(df = 17) = 16.145, p-val = 0.514 LRT(df = 17) = 37.990, p-val = 0.002 Model Results: estimate se zval pval ci.lb ci.ub -4.812 0.355 -13.537 <.001 -5.509 -4.115 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
predict(res, transf=transf.ilogit, digits=3)
pred ci.lb ci.ub pi.lb pi.ub 0.008 0.004 0.016 0.001 0.052
By default, a random/mixed-effects model using ML estimation is fitted by rma.glmm()
as was done in the article. Again, these results match what is reported in the article.
Hypergeometric-Normal Model for the Meta-Analysis of Odds Ratios
Next, the authors describe a "hypergeometric-normal model" for the meta-analysis of odds ratios. This is essentially a mixed-effects conditional logistic model. The model results from conditioning on the total number of cases within each study, leading to the non-central hypergeometric distribution for the 2×2 tables. Since fitting this model can be difficult and computationally expensive, one can approximate the exact likelihood by a binomial distribution, which works well when the number of cases in each study is small relative to the group sizes (as is the case here). Again, the results obtained from this model can be compared with the normal-normal model approach.
In fact, again, we start with the normal-normal model, using the log odds ratio as the outcome measure. This model can be fitted with:
res <- rma(measure="OR", ai=ai, n1i=n1i, ci=ci, n2i=n2i, data=dat, drop00=TRUE) print(res, digits=3)
Random-Effects Model (k = 17; tau^2 estimator: REML) tau^2 (estimated amount of total heterogeneity): 0.036 (SE = 0.299) tau (square root of estimated tau^2 value): 0.189 I^2 (total heterogeneity / total variability): 3.46% H^2 (total variability / sampling variability): 1.04 Test for Heterogeneity: Q(df = 16) = 15.812, p-val = 0.466 Model Results: estimate se zval pval ci.lb ci.ub -0.980 0.243 -4.027 <.001 -1.458 -0.503 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
The estimated average log odds ratio (i.e., $\hat{\mu} = -0.980$) can be back-transformed through exponentiation with:
predict(res, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 0.38 0.23 0.60 0.21 0.69
Note that there is one study with zero events in both arms. Such studies can be explicitly dropped from the analysis by setting drop00=TRUE
(careful: this is not the default in the rma()
function).
We can fit the hypergeometric-normal model with:
res <- rma.glmm(measure="OR", ai=ai, n1i=n1i, ci=ci, n2i=n2i, data=dat, model="CM.EL") print(res, digits=3)
Model fitting may take a few moments to complete. The results are:
Random-Effects Model (k = 17; tau^2 estimator: ML) Model Type: Conditional Model with Exact Likelihood tau^2 (estimated amount of total heterogeneity): 0.693 (SE = 0.640) tau (square root of estimated tau^2 value): 0.833 I^2 (total heterogeneity / total variability): 41.13% H^2 (total variability / sampling variability): 1.70 Tests for Heterogeneity: Wld(df = 16) = 11.866, p-val = 0.753 LRT(df = 16) = 28.609, p-val = 0.027 Model Results: estimate se zval pval ci.lb ci.ub -1.353 0.351 -3.855 <.001 -2.041 -0.665 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
And we can back-transform the estimated log odds ratio with:
predict(res, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 0.26 0.13 0.51 0.04 1.52
Since the likelihood for studies with zero events is flat, such studies are automatically dropped by the rma.glmm()
function (i.e., drop00=TRUE
by default for the rma.glmm()
function).
Finally, using the approximation to the exact likelihood, we can fit the same model with:
res <- rma.glmm(measure="OR", ai=ai, n1i=n1i, ci=ci, n2i=n2i, data=dat, model="CM.AL") print(res, digits=3)
Random-Effects Model (k = 17; tau^2 estimator: ML) Model Type: Conditional Model with Approximate Likelihood tau^2 (estimated amount of total heterogeneity): 0.601 tau (square root of estimated tau^2 value): 0.775 I^2 (total heterogeneity / total variability): 37.70% H^2 (total variability / sampling variability): 1.61 Tests for Heterogeneity: Wld(df = 16) = 11.115, p-val = 0.802 LRT(df = 16) = 27.392, p-val = 0.037 Model Results: estimate se zval pval ci.lb ci.ub -1.303 0.339 -3.847 <.001 -1.966 -0.639 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
predict(res, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 0.27 0.14 0.53 0.05 1.43
All of these results match what is reported in Table III (p. 3052).
Poisson-Normal Model for the Meta-Analysis of Incidence Rates
Niel-Weise et al. (2008) conducted a meta-analysis comparing the CRBSI risk when using anti-infective-treated versus standard catheters for total parenteral nutrition or chemotherapy. These data can be used to illustrate Poisson-normal models for the meta-analysis of incidence rates and incidence rate ratios. The data can be extracted from Table III (p. 119) in Niel-Weise et al. (2008) and can be loaded with:
dat <- dat.nielweise2008 dat
(I copy the dataset into 'dat', which is a bit shorter and therefore easier to type further below). The contents of the dataset are:
study authors year x1i t1i x2i t2i 1 1 Bong et al. 2003 7 1344 11 1988 2 2 Ciresi et al. 1996 8 1600 8 1461 3 3 Hanna et al. 2004 3 12012 14 10962 4 4 Harter et al. 2002 6 1536 10 1503 5 5 Jaeger et al. 2001 1 370 1 483 6 6 Jaeger et al. 2005 1 729 8 913 7 7 Logghe et al. 1997 17 6760 15 6840 8 8 Ostendorf et al. 2005 3 1107 7 1015 9 9 Pemberton et al. 1996 2 320 3 440
Here, the number of infections (x1i
and x2i
) and the total number of catheter days (t1i
and t2i
) are given for patients in the anti-infective-treated and standard catheter groups, respectively. Note that these data differ very slighty (for a few studies) from the data reported by Stijnen et al. (2010) in the appendix of their article (p. 3062)1) but these differences are inconsequential for the results.
Following the authors, we first divide the total number of catheter days by 1000, so that the estimated (average) incidence rates reflect the expected number of events per 1000 days:
dat$t1i <- dat$t1i/1000 dat$t2i <- dat$t2i/1000
Next, the normal-normal model can be fitted using the log incidence rate as the outcome measure. For patients in the standard catheter groups, this can be done with:
res <- rma(measure="IRLN", xi=x2i, ti=t2i, data=dat) print(res, digits=3)
Random-Effects Model (k = 9; tau^2 estimator: REML) tau^2 (estimated amount of total heterogeneity): 0.370 (SE = 0.259) tau (square root of estimated tau^2 value): 0.608 I^2 (total heterogeneity / total variability): 75.38% H^2 (total variability / sampling variability): 4.06 Test for Heterogeneity: Q(df = 8) = 36.384, p-val < .001 Model Results: estimate se zval pval ci.lb ci.ub 1.468 0.243 6.051 <.001 0.992 1.943 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
and the estimated average log incidence rate (i.e., $\hat{\mu} = 1.468$) can be back-transformed through exponentiation:
predict(res, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 4.34 2.70 6.98 1.20 15.66
The same model can be fitted to patients in the anti-infective-treated catheter groups with:
res <- rma(measure="IRLN", xi=x1i, ti=t1i, data=dat) print(res, digits=3)
Random-Effects Model (k = 9; tau^2 estimator: REML) tau^2 (estimated amount of total heterogeneity): 0.639 (SE = 0.467) tau (square root of estimated tau^2 value): 0.800 I^2 (total heterogeneity / total variability): 75.41% H^2 (total variability / sampling variability): 4.07 Test for Heterogeneity: Q(df = 8) = 25.427, p-val = 0.001 Model Results: estimate se zval pval ci.lb ci.ub 0.981 0.326 3.009 0.003 0.342 1.620 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
predict(res, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 2.67 1.41 5.05 0.49 14.49
Now we fit the Poisson-normal model to the same data. For patients in the standard catheter groups, this can be done with:
res <- rma.glmm(measure="IRLN", xi=x2i, ti=t2i, data=dat) print(res, digits=3)
Random-Effects Model (k = 9; tau^2 estimator: ML) tau^2 (estimated amount of total heterogeneity): 0.317 tau (square root of estimated tau^2 value): 0.563 I^2 (total heterogeneity / total variability): 72.38% H^2 (total variability / sampling variability): 3.62 Tests for Heterogeneity: Wld(df = 8) = 36.384, p-val < .001 LRT(df = 8) = 38.330, p-val < .001 Model Results: estimate se zval pval ci.lb ci.ub 1.401 0.231 6.063 <.001 0.948 1.853 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
predict(res, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 4.06 2.58 6.38 1.23 13.37
For patients in the anti-infective-treated catheter groups, we obtain:
res <- rma.glmm(measure="IRLN", xi=x1i, ti=t1i, data=dat) print(res, digits=3)
Random-Effects Model (k = 9; tau^2 estimator: ML) tau^2 (estimated amount of total heterogeneity): 0.654 tau (square root of estimated tau^2 value): 0.809 I^2 (total heterogeneity / total variability): 75.84% H^2 (total variability / sampling variability): 4.14 Tests for Heterogeneity: Wld(df = 8) = 25.427, p-val = 0.001 LRT(df = 8) = 44.488, p-val < .001 Model Results: estimate se zval pval ci.lb ci.ub 0.849 0.330 2.572 0.010 0.202 1.497 * --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
predict(res, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 2.34 1.22 4.47 0.42 12.96
These findings match what is reported in Table VII (p. 3054).
Binomial-Normal Model for the Meta-Analysis of Incidence Rate Ratios
Next, we consider the mixed-effects conditional Poisson regression model described by the authors. Through the conditioning on the total number of events within each study, we actually obtain a binomial-normal model for the data within each study. But first, we start again with the normal-normal model, using the log incidence rate ratio as the outcome measure:
res <- rma(measure="IRR", x1i=x1i, t1i=t1i, x2i=x2i, t2i=t2i, data=dat) print(res, digits=3)
Random-Effects Model (k = 9; tau^2 estimator: REML) tau^2 (estimated amount of total heterogeneity): 0.094 (SE = 0.212) tau (square root of estimated tau^2 value): 0.306 I^2 (total heterogeneity / total variability): 20.88% H^2 (total variability / sampling variability): 1.26 Test for Heterogeneity: Q(df = 8) = 9.698, p-val = 0.287 Model Results: estimate se zval pval ci.lb ci.ub -0.396 0.227 -1.747 0.081 -0.841 0.048 . --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
And the estimated average log incidence rate ratio (i.e., $\hat{\mu} = -0.396$) can again be back-transformed through exponentiation:
predict(res, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 0.67 0.43 1.05 0.32 1.42
The mixed-effects conditional Poisson regression model (i.e., the binomial-normal model) can be fitted with:
res <- rma.glmm(measure="IRR", x1i=x1i, t1i=t1i, x2i=x2i, t2i=t2i, data=dat, model="CM.EL") print(res, digits=3)
Random-Effects Model (k = 9; tau^2 estimator: ML) Model Type: Conditional Model with Exact Likelihood tau^2 (estimated amount of total heterogeneity): 0.123 tau (square root of estimated tau^2 value): 0.350 I^2 (total heterogeneity / total variability): 25.68% H^2 (total variability / sampling variability): 1.35 Tests for Heterogeneity: Wld(df = 8) = 9.698, p-val = 0.287 LRT(df = 8) = 11.602, p-val = 0.170 Model Results: estimate se zval pval ci.lb ci.ub -0.476 0.238 -2.004 0.045 -0.942 -0.010 * --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
predict(res, transf=exp, digits=2)
pred ci.lb ci.ub pi.lb pi.ub 0.62 0.39 0.99 0.27 1.42
These results match up with Table VIII (p. 3055) in the article.
References
Niel-Weise, B. S., Stijnen, T., & van den Broek, P. J. (2007). Anti-infective-treated central venous catheters: A systematic review of randomized controlled trials. Intensive Care Medicine, 33(12), 2058–2068.
Niel-Weise, B. S., Stijnen, T., & van den Broek, P. J. (2008). Anti-infective-treated central venous catheters for total parenteral nutrition or chemotherapy: A systematic review. Journal of Hospital Infection, 69(2), 114–123.
Stijnen, T., Hamza, T. H., & Ozdemir, P. (2010). Random effects meta-analysis of event outcome in the framework of the generalized linear mixed model with applications in sparse data. Statistics in Medicine, 29(29), 3046–3067.