This is one issue. However, there are more issues at hand here, since regplot(res, mod = ~ Continent) isn't gong to work either, as this will trigger the error "Can only specify a single variable via argument 'mod'." The reason is that the 'mod' argument is not meant to take a formula as input. However, when you do, then the length of 'mod' will be 2 (length(~ blah)), which triggers the error.
@Emanuel: While I don't know what you really intend to do here, I suspect you would like to show all levels of 'Continent' on the x-axis. This is not so easily accomplished with the regplot() function, which is designed for visualizing quantitative/continuous moderators. For example:
dat <- escalc(measure="RR", ai=tpos, bi=tneg,
ci=cpos, di=cneg, data=dat.bcg)
res <- rma(yi, vi, mods = ~ ablat, data=dat)
res
regplot(res, mod="ablat")
A dichotomous moderator is also easily handled:
dat$random <- ifelse(dat$alloc == "random", 1, 0)
res <- rma(yi, vi, mods = ~ random, data=dat)
res
regplot(res, mod="random")
We might want to make this look a bit nicer, maybe like this:
regplot(res, mod="random", xlab="Method of Treatment Allocation", xaxt="n")
axis(side=1, at=c(0,1), labels=c("Non-Random", "Random"))
However, things become more tricky with a factor variable that has 3 or more levels:
res <- rma(yi, vi, mods = ~ alloc, data=dat)
res
regplot(res, mod="alloc")
This will not work, since there are two 'alloc' dummy variables, but regplot() is designed to place a single variable on the x-axis. One could do:
par(mfrow=c(2,1))
regplot(res, mod="allocrandom")
regplot(res, mod="allocsystematic")
par(mfrow=c(1,1))
but this is showing the difference between random and not random (consisting of 'systematic' and 'alternate') in the first plot and the difference between systematic and not systematic (consisting of 'random' and 'alternate') in the second plot. Probably not how most people would want to visualize this. Instead, I suspect most would want to show three 'columns' of points, corresponding to the three levels, with lines connecting the fitted/predicted values for these levels.
I had not quite considered that some may want to do something like this with this function. I am not sure how easy it would be to add this kind of functionality directly to regplot() given some of the internal intricacies. However, with a small trick, we can still accomplish this. A model with a categorical moderator with p levels can be represented as a polynomial regression model to the degree p-1 where the first term is the linear one (which we want to place on the x-axis). This, combined with the possibility to pass predicted values to regplot() via the 'pred' argument, we can place all levels on the x-axis as follows:
dat$anum <- as.numeric(factor(dat$alloc))
res <- rma(yi, vi, mods = ~ poly(anum, degree=2, raw=TRUE), data=dat)
res
pred <- predict(res, newmods=unname(poly(1:3, degree=2, raw=TRUE)))
pred
regplot(res, mod=2, pred=pred, xvals=c(1:3), xlim=c(1,3), xlab="Allocation Method", xaxt="n")
axis(side=1, at=1:3, labels=levels(factor(dat$alloc)))
Best,
Wolfgang
-----Original Message-----
From: R-sig-meta-analysis <r-sig-meta-analysis-bounces at r-project.org> On Behalf
Of Michael Dewey via R-sig-meta-analysis
Sent: Friday, September 6, 2024 11:04
To: R Special Interest Group for Meta-Analysis <r-sig-meta-analysis at r-
project.org>
Cc: Michael Dewey <lists at dewey.myzen.co.uk>
Subject: Re: [R-meta] Metafor - regplot() for a categorical (and interaction)
variable
Dear Emanuel
I think you will find that the parameter is named mods in rma.uni but
mod in regplot.
Michael
On 06/09/2024 00:56, Emanuel Schembri via R-sig-meta-analysis wrote:
Hi
I am trying to plot a bubble plot using the replot() function in
Metafor. However, I cannot make it work for a categorical moderator
(and an interaction between a categorical and numerical variable),
while it does work when inputting a numerical or date.
Any help would be greatly appreciated.
Regards,
Emanuel
res <- rma(yi, vi, mods = ~ Continent, data= df)
regplot(res, mods = ~ Continent)
Mixed-Effects Model (k = 24; tau^2 estimator: REML)
tau^2 (estimated amount of residual heterogeneity): 3.5128 (SE = 1.1529)
tau (square root of estimated tau^2 value): 1.8743
I^2 (residual heterogeneity / unaccounted variability): 99.79%
H^2 (unaccounted variability / sampling variability): 477.76
R^2 (amount of heterogeneity accounted for): 8.05%
Test for Residual Heterogeneity:
QE(df = 21) = 1150.4471, p-val < .0001
Test of Moderators (coefficients 2:3):
QM(df = 2) = 5.7609, p-val = 0.0561
Model Results:
estimate se zval pval ci.lb ci.ub
intrcpt 2.2738 1.3346 1.7037 0.0884 -0.3420 4.8897 .
ContinentAsia 4.3152 1.8165 2.3756 0.0175 0.7549 7.8755 *
ContinentEurope 2.6747 1.4053 1.9032 0.0570 -0.0797 5.4291 .
---
Signif. codes: 0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1
regplot(res, mods = ~ Continent)
Error in regplot.rma(res, mods = ~Continent) :
Must specify 'mod' argument for models with multiple predictors.