[R] Factor moderators in metafor

[Viechtbauer Wolfgang (STAT)]

There is a simple explanation:

1) The command:

factor(ftype)

does not actually turn 'ftype' permanently into a factor, since you are not re-assigning it back to the object 'ftype'. You have to use:

ftype <- factor(ftype)

2) If you want to use the formula interface for specifying moderators, you have to use mods = ~ <formula>, so in other words:

rma(KL, VL, mods = ~ ftype)

after you have made 'ftype' a factor (see 1).

Or you can simply use:

rma(KL, VL, mods = ~ factor(ftype))

which does the conversion of 'ftype' into a factor within the model formula.

Best,

Wolfgang

--   
Wolfgang Viechtbauer, Ph.D., Statistician   
Department of Psychiatry and Psychology   
School for Mental Health and Neuroscience   
Faculty of Health, Medicine, and Life Sciences   
Maastricht University, P.O. Box 616 (VIJV1)   
6200 MD Maastricht, The Netherlands   
+31 (43) 388-4170 | http://www.wvbauer.com   

> -----Original Message-----
> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
> On Behalf Of John Hodgson
> Sent: Thursday, August 09, 2012 12:56
> To: r-help at r-project.org
> Subject: [R] Factor moderators in metafor
> 
> I'm puzzled by the behaviour of factors in rma models, see example and
> comments below. I'm sure there's a simple explanation but can't see it...
> 
> Thanks for any input
> 
> John Hodgson
> 
> 
> ------------------------------------- code/selected output ---------------
> --
> 
> 
> library(metafor)
> 
> ##    Set up data (from Lenters et al  A Meta-analysis of Asbestos and
> Lung
> Cancer...
> ##    Environmental Health Perspectives • volume 119 | number 11 |
> November
> 2011)
> 
> KL = c(0.02905, 0.06929, -0.1523, 1.6441, 0.1215, 0.3975, 1.0566, 0.1257,
> 0.2277, 0.06791, 0.08164, 0.2526, 0.07577, 0.03266, 0.1141, 0.1836,
> 1.8276,
> 0.4149, 15.4974)
> SE = c(0.006633, 0.09335, 0.08909, 0.4297, 0.07858, 0.1753, 0.3679,
> 0.1837,
> 0.2172, 0.2775, 0.4201, 0.1976, 0.7688, 0.06507, 0.06239, 0.09061, 0.9509,
> 0.2181, 7.331)
> 
> VL = SE*SE
> 
> amph =   c(0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
> mix =    c(0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1)
> ftype =  c(0, 0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2)
> 
> 
> factor(amph)
> factor(ftype)
> factor(mix)
> 
> ##  Fit ftype...
> 
> > rma(KL,VL,mods=ftype)
> 
> 
> Mixed-Effects Model (k = 19; tau^2 estimator: REML)
> 
> tau^2 (estimate of residual amount of heterogeneity): 0.0111 (SE = 0.0095)
> tau (sqrt of the estimate of residual heterogeneity): 0.1054
> 
> Test for Residual Heterogeneity:
> QE(df = 17) = 43.0937, p-val = 0.0005
> 
> Test of Moderators (coefficient(s) 2):
> QM(df = 1) = 1.1069, p-val = 0.2928
> 
> Model Results:
> 
>          estimate      se    zval    pval    ci.lb   ci.ub
> intrcpt    0.0811  0.0606  1.3380  0.1809  -0.0377  0.2000
> mods       0.0473  0.0449  1.0521  0.2928  -0.0408  0.1353
> 
> 
> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> 
> 
> ##   why only one df for the 2-level factor?
> ##   in other words, why isn't the above model the same as the
> following...
> 
> 
> >
> >
> >
> > rma(KL,VL,mods=cbind(amph,mix))
> 
> Mixed-Effects Model (k = 19; tau^2 estimator: REML)
> 
> tau^2 (estimate of residual amount of heterogeneity): 0.0030 (SE = 0.0046)
> tau (sqrt of the estimate of residual heterogeneity): 0.0549
> 
> Test for Residual Heterogeneity:
> QE(df = 16) = 37.5762, p-val = 0.0017
> 
> Test of Moderators (coefficient(s) 2,3):
> QM(df = 2) = 6.9220, p-val = 0.0314
> 
> Model Results:
> 
>          estimate      se    zval    pval    ci.lb   ci.ub
> intrcpt    0.0380  0.0402  0.9448  0.3447  -0.0408  0.1169
> amph       0.2879  0.1163  2.4754  0.0133   0.0599  0.5158  *
> mix        0.0888  0.0625  1.4199  0.1556  -0.0338  0.2114
> 
> ---
> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1