I am a new user using R metafor package to develop a forest plot for meta-analysis project. I am struggling a bit to develop three level forest plot with subgroups (I managed well to develop two-level plots with subgroups which is simple and easy to perform).

I use following commands to develop three level plots: but end up receiving this error:

"Error in verbose > 2 :

comparison (6) is possible only for atomic and list types"

"Error in eval(mf.yi, data, enclos = sys.frame(sys.parent())) :

object 'yi' not found"

"Error in addpoly(DDD.r, row = 59.5, mlab = mlabfun("RE Model for Subgroup", :

object 'DDD.r' not found"

Here are the script:

DDDnew$var <- (DDDnew$se_log_rate)^2

DDD<- rma.mv(yi=log_rate,

v= var,

slab = study,

data = DDDnew,

random = ~ 1 | study/log_rate,

test = "t",

method = "REML")

(Error in verbose > 2 :

comparison (6) is possible only for atomic and list types)

mlabfun <- function(text, DDD) {

list(bquote(paste(.(text),

" (Q = ", .(formatC(DDD$QE, digits=2, format="f")),

", df = ", .(DDD$k - DDD$p),

", p ", .(metafor:::.pval(DDD$QEp, digits=2, showeq=TRUE, sep=" ")), "; ",

I^2, " = ", .(formatC(DDD$I2, digits=1, format="f")), "%, ",

tau^2, " = ", .(formatC(DDD$tau2, digits=2, format="f")), ")")))}

tiff(file="~/R plot/Subgroup-Income.tiff", units="in", width=15, height=12, res=500, compression = "lzw")

units="in", width=15, height=10, res=500, compression = "lzw")

forest(DDD, xlim=c(-16, 4.6), at=log(c(.05, .25, 1, 4)), atransf=exp,

ilab=cbind(DDDnew$antclass,DDDnew$country), ilab.xpos=c(-9.5,-6),

cex=.75, ylim=c(0.5, 64),

order=DDDnew$income, rows=c(1:9,14:61),

showweights = TRUE,

xlab = "Rate ratio (axis in log scale)",

digits = 2,

col = "blue",

border = "black",

mlab=mlabfun("RE Model for All Studies",DDD),

header="Author(s) and Year")

op <- par(cex=0.75, font=2)

text(c(-9.5,-6), 63, c("Antibiotic", "Country"))

op <- par(font=4)

text(-16, c(62,10), c("HICs", "LMICs"), pos=4)

### fit random-effects model in the two subgroups

DDD.s <- rma.mv(yi, v, subset=(income==c("HIC")), data=DDDnew)

(Error in eval(mf.yi, data, enclos = sys.frame(sys.parent())) :

object 'yi' not found)

DDD.r <- rma.mv(yi, v, subset=(income==c("LMIC")), data=DDDnew)

(Error in eval(mf.yi, data, enclos = sys.frame(sys.parent())) :

object 'yi' not found)

addpoly(DDD.s, row=7.5, mlab=mlabfun("RE Model for Subgroup", DDD.s))

(Error in addpoly(DDD.s, row = 7.5, mlab = mlabfun("RE Model for Subgroup", :

object 'DDD.s' not found)

addpoly(DDD.r, row=59.5, mlab=mlabfun("RE Model for Subgroup", DDD.r))

(Error in addpoly(DDD.r, row = 59.5, mlab = mlabfun("RE Model for Subgroup", :

object 'DDD.r' not found)

DDD <- rma.mv(DDDnew$log_rate,DDDnew$se_log_rate, mods = ~ income, data=DDDnew)

text(-16, -1.8, pos=4, cex=0.75, bquote(paste("Test for Subgroup Differences: ",

Q[M], " = ", .(formatC(DDD$QM, digits=2, format="f")), ", df = ", .(DDD$p - 1),

", p = ", .(formatC(DDD$QMp, digits=2, format="f")))))

dev.off()

I would be very provide more info if something is not clear. Your advice and guidance would be much appreciated.