Skip to content

[R-meta] multi-model prediction with rma.mv ?

3 messages · Kelly Gravuer, Viechtbauer Wolfgang (STAT)

Original
Kelly Gravuer
# 
Hello all,
Using the rma.mv() function in metafor, I'm running models with two
random effects and four moderator variables. The format of the model
is:

m1 <- rma.mv(yi ~ mod1 + mod2 + mod3 + mod4, vi, random = list(~ 1 |
rand1, ~1 | rand2), tdist=TRUE, method="ML", data=d)

I would like to do two things: (1) understand the relationship of each
of the four moderator variables to the effect size; (2) predict an
effect size with confidence interval (CI) for some new suites of
values of the moderators.

I have found multi-model inference using glmulti() to be a very useful
approach to my first objective, using the methods described here:
http://www.metafor-project.org/doku.php/tips:model_selection_with_glmulti.
I would now like to use this multi-model suite to make predictions
with CIs.

I realize I can run the model as specified above (including all four
moderators as main effects) and then use predict.rma() to generate
estimated effect sizes + CIs.  I realize I can also calculate an
effect size estimate that is informed by the multiple model fits using
the coefficient estimates in the multi-model coefficient table.
However, I don't know how to get CIs for these effect size estimates
that are informed by the multiple model fits.

I think one way to do this would be to write a predict function that
links rma.mv with glmulti.  I don't feel up to this task myself, but
wanted to inquire whether anyone on this list has written such a
function?  Or alternatively, does anyone have other suggestions for
how to use the information from the multiple model fits to inform the
predicted CIs?

Thanks so much for any thoughts,
Kelly
Viechtbauer Wolfgang (STAT)
# 
Dear Kelly,

If I understand you correctly, you want to obtain multimodel predicted values and corresponding CIs.

See, for example, Anderson (2008; Model Based Inference in the Life Sciences: A Primer on Evidence), Chapter 5.

I don't see a way to get this to work 'out of the box', but you can make this work with predict() from metafor and glmulti with a bit of code. To illustrate, I will use the same example that you linked to (i.e., http://www.metafor-project.org/doku.php/tips:model_selection_with_glmulti). So, to start:

library(metafor)
library(glmulti)

dat <- get(data(dat.bangertdrowns2004))
dat <- dat[!apply(dat[,c("length", "wic", "feedback", "info", "pers", "imag", "meta")], 1, anyNA),]

rma.glmulti <- function(formula, data, ...)
   rma(formula, vi, data=data, method="ML", ...)

res <- glmulti(yi ~ length + wic + feedback + info + pers + imag + meta, data=dat,
               level=1, fitfunction=rma.glmulti, crit="aicc", confsetsize=128)

### get model weights
weights <- weightable(res)$weights

### vector with the values for which you want to make a prediction
x <- c("length"=15, "wic"=1, "feedback"=1, "info"=0, "pers"=0, "imag"=1, "meta"=1)

### get predicted values based on all models

preds <- list()

for (j in 1:length(weights)) {

   model <- res at objects[[j]]
   vars <- names(coef(model))[-1]

   ### need special handling for the intercept-only model
   if (length(vars) == 0) {
      preds[[j]] <- predict(model)
   } else {
      preds[[j]] <- predict(model, newmods=x[vars])
   }

}

### multimodel prediction
yhat <- sum(weights * sapply(preds, function(x) x$pred))
yhat

### multimodel SE (based on the unconditional variance)
se <- sqrt(sum(weights * sapply(preds, function(x) x$se^2 + (x$pred - yhat)^2)))
se

### multimodel 95% CI
yhat + c(-1,1)*qnorm(.975)*se

That's basically it. You should be able to adapt this to your case.

Best,
Wolfgang

-----Original Message-----
From: R-sig-meta-analysis [mailto:r-sig-meta-analysis-bounces at r-project.org] On Behalf Of Kelly Gravuer
Sent: Tuesday, 02 January, 2018 19:43
To: r-sig-meta-analysis at r-project.org
Subject: [R-meta] multi-model prediction with rma.mv ?

Hello all,
Using the rma.mv() function in metafor, I'm running models with two
random effects and four moderator variables. The format of the model
is:

m1 <- rma.mv(yi ~ mod1 + mod2 + mod3 + mod4, vi, random = list(~ 1 |
rand1, ~1 | rand2), tdist=TRUE, method="ML", data=d)

I would like to do two things: (1) understand the relationship of each
of the four moderator variables to the effect size; (2) predict an
effect size with confidence interval (CI) for some new suites of
values of the moderators.

I have found multi-model inference using glmulti() to be a very useful
approach to my first objective, using the methods described here:
http://www.metafor-project.org/doku.php/tips:model_selection_with_glmulti.
I would now like to use this multi-model suite to make predictions
with CIs.

I realize I can run the model as specified above (including all four
moderators as main effects) and then use predict.rma() to generate
estimated effect sizes + CIs.  I realize I can also calculate an
effect size estimate that is informed by the multiple model fits using
the coefficient estimates in the multi-model coefficient table.
However, I don't know how to get CIs for these effect size estimates
that are informed by the multiple model fits.

I think one way to do this would be to write a predict function that
links rma.mv with glmulti.  I don't feel up to this task myself, but
wanted to inquire whether anyone on this list has written such a
function?  Or alternatively, does anyone have other suggestions for
how to use the information from the multiple model fits to inform the
predicted CIs?

Thanks so much for any thoughts,
Kelly
Kelly Gravuer
# 
Dear Wolfgang,
Thank you so much - this was exactly what I needed!
Best,
Kelly

On Tue, Jan 2, 2018 at 12:46 PM, Viechtbauer Wolfgang (SP) <
wolfgang.viechtbauer at maastrichtuniversity.nl> wrote: