rmetaformultilevel-analysis

Metafor updated degrees of freedom


Edit: changed code to include test = "t"

I'm hoping to better understand how the updated dev version of Metafor 2.5-101 will help me to adjust my degrees of freedom in a multi-level model to provide some protection against type 1 error.

My understanding of this has come from the Nakagawa preprint "Methods for testing publication bias in ecological and evolutionary meta-analyses" https://ecoevorxiv.org/k7pmz/ and their "Supplemental_Impleentation_Example.Rmd" file, following along with their line 133-142:

Before moving on to some useful corrections, users should be aware that the most up-to-date version of metafor (version 2.5-101) does now provide users with some protection against Type I errors. Instead of using the number of effect sizes in the calculation of the degrees of freedom we can instead make use of the total numbers of papers instead. We show in our simulations that a "papers-1" degrees of freedom can be fairly good. This can be implemented as follows after installing the development version of metafor (see "R Packages Required" above):

mod_multilevel_pdf = rma.mv(yi = yi, V = vi, mods = ~1, 
                                    random=list(~1|study_id,~1|obs), 
                                    data=data, test="t", dfs = "contain")
summary(mod_multilevel_pdf)

We can see here that the df for the model has changes from 149 to 29, and the p-value has been adjusted accordingly.

So my understanding here is that the model now shows df as 29 (the original no. of papers (30) -1, instead of the no. of papers x no. of effects (30 papers with 5 effects each (150) -1)

Adapting this to my code, where I have n=18 papers and total of n=24 effects, I would expect using the above code would adjust my df to 17 (the original no. of papers (18) -1), however I still have df as 23 (total no. of effects (24) -1).

The output using the df code:

mod_multilevel_pdf = rma.mv(yi = yi, V = vi, mods = ~1, 
                            random=list(~1|study_id,~1|es_id), 
                            data=dat, test="t", dfs = "contain")
summary(mod_multilevel_pdf)

Is:

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

  logLik  Deviance       AIC       BIC      AICc 
-30.2270   60.4540   66.4540   69.8604   67.7171   

Variance Components:

            estim    sqrt  nlvls  fixed    factor 
sigma^2.1  0.6783  0.8236     18     no  study_id 
sigma^2.2  0.1416  0.3763     24     no     es_id 

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

Model Results:

estimate      se     tval  df    pval    ci.lb   ci.ub 
 -0.3508  0.2219  -1.5809  17  0.1323  -0.8190  0.1174    

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

Quite stumped on this one! Any help would be majorly appreciated.


Solution

  • You neither have df=17 nor df=23, since you did not specify that you want a t-test. With test="t", dfs = "contain", you will get the expected t-test with df=17.