On 19 Jun 2018, at 17:15, Viechtbauer, Wolfgang (SP) <wolfgang.viechtbauer at maastrichtuniversity.nl> wrote:
Hi C?sar,
I do not know of an automated way of doing this. I looked at the 'pdp' package, but it might take some effort to make it work together with metafor. However, doing this manually shouldn't be too complicated. Here is an example using a relatively simple model with two predictors:
library(metafor)
dat <- get(data(dat.bcg))
dat <- escalc(measure="RR", ai=tpos, bi=tneg, ci=cpos, di=cneg, data=dat.bcg)
res <- rma(yi, vi, mods = ~ ablat + year, data=dat)
res
xs <- seq(min(dat$ablat), max(dat$ablat), by=1)
pred <- sapply(xs, function(x) mean(predict(res, newmods = cbind(x, dat$year))$pred))
plot(xs, pred, type="o")
### same as:
pred <- predict(res, newmods = cbind(xs, mean(dat$year)))$pred
lines(xs, pred, col="red")
xs <- seq(min(dat$year), max(dat$year), by=1)
pred <- sapply(xs, function(x) mean(predict(res, newmods = cbind(dat$ablat, x))$pred))
plot(xs, pred, type="o")
### same as:
pred <- predict(res, newmods = cbind(mean(dat$ablat), xs))$pred
lines(xs, pred, col="red")
In your case:
M <- rma(yi, vi, dat, mods= ~ A*B + C + D)
xs <- seq(min(dat$C), max(dat$C), by=1) # or use an appropriate 'by' value
pred <- sapply(xs, function(x) mean(predict(res, newmods = cbind(dat$A, dat$B, dat$A*dat$B, x, dat$D))$pred))
plot(xs, pred, type="o")
xs <- seq(min(dat$D), max(dat$D), by=1) # or use an appropriate 'by' value
pred <- sapply(xs, function(x) mean(predict(res, newmods = cbind(dat$A, dat$B, dat$A*dat$B, dat$C, x))$pred))
plot(xs, pred, type="o")
Best,
Wolfgang
-----Original Message-----
From: R-sig-meta-analysis [mailto:r-sig-meta-analysis-bounces at r-project.org] On Behalf Of Cesar Terrer Moreno
Sent: Tuesday, 19 June, 2018 14:52
To: r-sig-meta-analysis at r-project.org
Subject: [R-meta] Partial dependence plots
Dear all,
I have a mixed-effects model of the form:
M <- rma(yi, vi, dat, mods= ~ A*B + C + D)
How can I produce partial dependence plots of e.g. C and D (predicted effect sizes of C and D as a function of the value of each predictor variable). The idea is to show that once the interaction A*B is taken into account, C and D explain very little of the overall effect size.
Can this be coded into a function to make partial dependence plots for all variables?
Thanks
C?sar