Binary aggregate outcomes

metaGLMM authors

Data contract

For a binomial analysis, y is the observed proportion and ni is the number of trials. With two treatment arms per study, use re_group to identify the study and trt to identify the 0/1 treatment loading. This shares one random treatment effect within each study.

binary_dat <- data.frame(
  study = factor(rep(paste0("Study ", 1:6), each = 2)),
  treatment = rep(c(0, 1), 6),
  events = c(20, 60, 25, 20, 30, 25, 15, 55, 35, 15, 25, 70),
  ni = rep(c(100, 90, 110, 95, 105, 120), each = 2)
)
binary_dat$y <- binary_dat$events / binary_dat$ni
binary_dat$vi <- 1 / (binary_dat$ni * binary_dat$y * (1 - binary_dat$y))
binary_dat
#>      study treatment events  ni         y         vi
#> 1  Study 1         0     20 100 0.2000000 0.06250000
#> 2  Study 1         1     60 100 0.6000000 0.04166667
#> 3  Study 2         0     25  90 0.2777778 0.05538462
#> 4  Study 2         1     20  90 0.2222222 0.06428571
#> 5  Study 3         0     30 110 0.2727273 0.04583333
#> 6  Study 3         1     25 110 0.2272727 0.05176471
#> 7  Study 4         0     15  95 0.1578947 0.07916667
#> 8  Study 4         1     55  95 0.5789474 0.04318182
#> 9  Study 5         0     35 105 0.3333333 0.04285714
#> 10 Study 5         1     15 105 0.1428571 0.07777778
#> 11 Study 6         0     25 120 0.2083333 0.05052632
#> 12 Study 6         1     70 120 0.5833333 0.03428571

The finite vi values above are useful for displaying the arm-level contrasts later. The binomial likelihood uses ni to represent the number of trials.

Fit the grouped model

binary_fit <- metaGLMM(
  y ~ treatment,
  data = binary_dat,
  vi = binary_dat$vi,
  ni = binary_dat$ni,
  tau2 = NA,
  family = binomial(link = "logit"),
  tau2_var = TRUE,
  re_group = binary_dat$study,
  trt = "treatment",
  fast = TRUE,
  ghq_Q = 40L
)

summary(binary_fit)
#> Aggregate-data generalized linear mixed-effects meta-analysis
#> 
#> Call:
#> metaGLMM(formula = y ~ treatment, data = binary_dat, vi = binary_dat$vi, 
#>     ni = binary_dat$ni, tau2 = NA, family = binomial(link = "logit"), 
#>     tau2_var = TRUE, re_group = binary_dat$study, trt = "treatment", 
#>     fast = TRUE, ghq_Q = 40L)
#> 
#> Family:binomial(logit)
#> Random-effect structure: grouped_treatment_slope 
#> Integration: ghq 
#> 
#> Fixed effects:
#>             Estimate Std. Error z value Pr(>|z|)    
#> (Intercept) -1.14209    0.09378 -12.179   <2e-16 ***
#> treatment    0.50537    0.36480   1.385    0.166    
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Heterogeneity:
#>   tau^2: 0.6835 
#>   tau:   0.8268 
#> 
#> Fit diagnostics:
#>   Optimizer: L-BFGS-B 
#>   Convergence code: 0 
#>   Boundary tau^2: no 
#>   Positive-definite Hessian: yes
coef(binary_fit)
#> (Intercept)   treatment 
#>   -1.142094    0.505369
confint(binary_fit, parm = "treatment", method = "wald")
#>                lower   upper
#> treatment -0.2096316 1.22037
stopifnot(is.finite(binary_fit$tau), binary_fit$tau > 0)

The treatment effects vary deliberately across studies (including effects in both directions), making the grouped random slope and its prediction interval visible in this small example.

The treatment coefficient is on the log-odds-ratio scale. A two-level factor can be used in place of the numeric treatment column; the second factor level is normalized to one. Logical treatment indicators are also accepted.

Study-level contrasts and a forest plot

The paired two-arm structure identifies one treatment contrast per study. forest() obtains these study estimates and variances through as_metafor_data() and selects the fixed treatment coefficient automatically:

binary_contrasts <- as_metafor_data(binary_fit)
binary_contrasts
#>           yi         vi    slab
#> 1  1.7917595 0.10416667 Study 1
#> 2 -0.2972515 0.11967033 Study 2
#> 3 -0.2429462 0.09759804 Study 3
#> 4  1.9924302 0.12234848 Study 4
#> 5 -1.0986123 0.12063492 Study 5
#> 6  1.6714733 0.08481203 Study 6

forest(binary_fit,
       xlab = "Log odds ratio",
       ci_methods = c("Wald", "profile", "SBC"))

Here yi is the within-study treated-minus-control log-odds contrast and vi is its sampling variance. The treatment-arm rows are therefore combined into one comparable effect estimate per study before plotting.

Explicit estimate, vi, labels, and parm arguments remain available for custom displays. If a study has a zero cell, its automatically derived display variance may be infinite; supply a finite, pre-specified correction before requesting a forest plot.

mirror server hosted at Truenetwork, Russian Federation.