Skip to content

[R-meta] Plotting rma.mv as forest plots for individual moderator levels

2 messages · Julian Packheiser, Wolfgang Viechtbauer

#
Dear members of the mailings list,

Currently, we are conducting a multivariate multilevel meta-analysis with multiple outcomes coming from the same study. I use both study_id and effectsize_id for the calculation of V as well as as random effects in the model. We also added the laboratory in which the study was conducted as a random effect. The code thus looks like this:

V <- vcalc(vi, cluster=study_id, obs=esid, data=mydata, rho=0.6)

# random-effects model
res <- rma.mv(yi, V, mods = ~ outcome-1, random = ~ 1 | lab/study_id/esid, data=mydata)

My first question is, should we use a sensitivity approach for rho values as they are unknown?

My second question is about forest plots. Ideally, I would like to plot the results for each moderator separately. If I use forest.rma, I get all results in one plot. An overall result for this meta-analysis is however not meaningful since some outcomes are positively correlated and some are negatively correlated with the variable of interest.
I assume that I cannot simply do a forest plot of a rma.mv result for each outcome individually as the results obviously differ from a multivariate approach. Could I simply re-arrange my outcomes and only plot them ordered by outcome? I assume the plotted values are yi.f and vi.f for the means and variances, but I could not find values for the fitted polygons.

Apologies if this question has already been asked, but I could not find it in the archives.

Best
Julian

---
Julian Packheiser
Postdoc in the Social Brain Lab
Netherlands Institute for Neuroscience
1 day later
#
Dear Julian,

Please see below for my responses.

Best,
Wolfgang
That's always a good idea when 'guestimating' a value. You could also consider using cluster-robust inference methods as a follow-up.
With 'moderators', I assume you mean 'outcome' (as in the model above).

I would use forest() (forest.default to be precise) to draw a forest plot per outcome and then use addpoly() to add the estimated pooled effect from the model above to each plot. For example:

library(metafor)

dat <- dat.berkey1998
V <- vcalc(vi=1, cluster=author, rvars=c(v1i, v2i), data=dat)
res <- rma.mv(yi, V, mods = ~ outcome - 1, random = ~ outcome | trial, struct="UN", data=dat)
res

with(dat, forest(yi, vi, subset=outcome=="AL", slab=paste0(author, year), header=TRUE, ylim=c(-1.5,8)))
addpoly(predict(res, newmods=c(1,0)), row=-1)

with(dat, forest(yi, vi, subset=outcome=="PD", slab=paste0(author, year), header=TRUE, ylim=c(-1.5,8)))
addpoly(predict(res, newmods=c(0,1)), row=-1)

I would recommend to add a note to these plots that the pooled estimate is based on the multivariate/multilevel model (as indeed pooling the estimates within each plot would lead to a slightly different result).