metafor, aggregate() data gives a different estimate than the multilevel model

36 views Asked by At

I'm following this guide to aggregate data at the study level for a meta-analysis i'm conducting.

Viechtbauer (at least i think it was him who wrote the article) says that "When aggregating the estimates in this manner, it turns out that a meta-analysis of these pooled estimates yields exactly the same results as those obtained from the multilevel model above"

Which is not what happens in my case.

I think the problem is that in his case, he was using measure = "SMD" in escalc(), while i'm using "SMCR" for my meta-analysis.

df <- escalc(measure="SMCR", 
             m1i=me_post_int, 
             m2i=me_pre_int, 
             sd1i=sd_pre_int, 
             ni=num_int,
             ri=r_int_imp,
             data=df, 
             append=TRUE,  
             replace = F,
             var.names = c("yi_int", "vi_int")) 

df <- escalc(measure="SMCR", 
             m1i  = me_post_c,
             m2i  = me_pre_c, 
             sd1i = sd_pre_c, 
             ni   = num_c, 
             ri   = r_c_imp, 
             data = df, 
             append=TRUE, 
             replace = F,
             var.names = c("yi_c", "vi_c")) 

#reversing ES for outcomes with a negative direction
df$yi_int <- ifelse(df$direction == "neg", yes = df$yi_int*(-1), no = df$yi_int)
df$yi_c <- ifelse(df$direction == "neg", yes = df$yi_c*(-1), no = df$yi_c)


df$yi <- df$yi_int - df$yi_c

df$vi <- df$vi_int + df$vi_c

This is the code and output of the multilevel model, where df_pain is a subset of df:

pain_3lv<- rma.mv(yi, vi, random = ~ 1 | study_id/Unique_ID_ES, 
                  method="REML", data=df_pain)
pain_3lv

Multivariate Meta-Analysis Model (k = 50; method: REML)

Variance Components:

            estim    sqrt  nlvls  fixed                 factor 
sigma^2.1  0.2641  0.5139     14     no               study_id 
sigma^2.2  0.0764  0.2764     50     no  study_id/Unique_ID_ES 

Test for Heterogeneity:
Q(df = 49) = 271.0640, p-val < .0001

Model Results:

estimate      se    zval    pval   ci.lb   ci.ub     
  0.4517  0.1581  2.8568  0.0043  0.1418  0.7616  ** 

---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

And this is the code and output of the aggregated model:

agg <- aggregate(df_pain, cluster=study_id, V=vcov(pain_3lv, type="obs"),  addk=TRUE)

res <- rma(yi, vi, method="EE", data=agg, slab = study_id)
res

Equal-Effects Model (k = 14)

I^2 (total heterogeneity / total variability):   79.70%
H^2 (total variability / sampling variability):  4.93

Test for Heterogeneity:
Q(df = 13) = 64.0479, p-val < .0001

Model Results:

estimate      se    zval    pval   ci.lb   ci.ub      
  0.3016  0.0475  6.3448  <.0001  0.2084  0.3948  *** 

---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

The difference in the method use to measure the effect sizes (SMD vs SMCR) is the only difference i can see in the code, i have otherwise copy-pasted his code. I used the "EE", equal effect model as prompted in the guide, but even if i try to use the same REML model, i still get a different estimate from the original one.

0

There are 0 answers