How do I successfully make a contour-enhanced funnel plot with a legend?

I am trying to create a contour-enhanced funnel plot in R for a meta-analysis. A contour-enhanced funnel plot also shows the areas of the plot corresponding to different p-values (see Peters et al, 2008, Contour-enhanced meta-analysis funnel plots help distinguish publication bias from other causes of asymmetry).

I found this very useful guide to meta-analytic procedure online (https://bookdown.org/MathiasHarrer/Doing_Meta_Analysis_in_R/smallstudyeffects.html) suggesting the following code to create such a plot:

funnel(m.hksj, xlab="Hedges' g", 
       contour = c(.95,.975,.99),
       col.contour=c("darkblue","blue","lightblue"))+
legend(1.4, 0, c("p < 0.05", "p<0.025", "< 0.01"),bty = "n",
       fill=c("darkblue","blue","lightblue"))

When executing this code, the plot is produced, but without legend. The error says that there is a non-numeric argument for the binary operator. I also make sure to use specifically target the funnel() function from the meta package, as it says it in the guide, as metafor is loaded as well in my environment. Does anybody know how to fix this and create a nice legend?

Also, I tried the following code suggested by the creators of the metafor package (http://www.metafor-project.org/doku.php/plots:contour_enhanced_funnel_plot):

funnel(res, level=c(.90, .95, .99),shade=c("white", "gray55", "gray75"), refline=0, legend=TRUE)

There is a lot of warnings regarding the level argument, and it does not even produce a working funnel plot, let alone with a legend. In this case, I assumed that this is caused by the fact that I submit not a rma object, but a meta object from the meta package. On the other hand, then it should work with the first code above, but it does not.

Here is an example based on the code provided in the guide mentioned above

data = data.frame("Author" = c("Jones", "Goldman", "Townsend", "Martin", "Rose"),
                  "TE" = c(0.23, 0.56, 0.78, 0.23, 0.33),
                  "seTE" = c(0.324, 0.235, 0.394, 0.275, 0.348),
                  "subgroup" = c("one", "one", "two", "two", "three"))

m.data <- metagen(TE,
                  seTE,
                  data = data,
                  studlab = paste(Author),
                  comb.fixed = FALSE,
                  comb.random = TRUE,
                  method.tau = "SJ",
                  hakn = TRUE,
                  prediction = TRUE,
                  sm = "SMD")

funnel(m.data, xlab="Hedges' g", 
       contour = c(.95,.975,.99),
       col.contour=c("darkblue","blue","lightblue"))+
legend(1.4, 0, c("p < 0.05", "p<0.025", "< 0.01"),bty = "n",
       fill=c("darkblue","blue","lightblue"))

Solution

You're getting the error because of the plus sign, and the legend isn't printing because the coordinates are probably outside of the range of the plot.

Try removing the plus sign and changing the coordinates of the legend, like this:

funnel(m.data, xlab="Hedges' g", 
           contour = c(.95,.975,.99),
           col.contour=c("darkblue","blue","lightblue"))
legend("topright", c("p < 0.05", "p<0.025", "< 0.01"),bty = "n",
           fill=c("darkblue","blue","lightblue"))