Adding title, axis labels and facet wrap to brms marginal_effects

[GabriellaS-K]

Hi,

I am using marginal_effects to plot my results, and there are a number of different predictors, which give me a number of different graphs. When I run it in R, my plots are separated out, and I see Hit <Return> to see next plot: in the console. Is there a way to add all the plots onto one, so I can visualise them together? It took a long time to run the model, and then the plot-so if there is a way to do this that doesn’t involve re-running, that would be great!

code:
mod1plot<-conditional_effects(mod1_all, categorical=TRUE)

Once they’re all combined, I’d like to add a title as well.

Thank you for your help!

[torkar]

In the manual for conditional_effects() you can see how one can save the plot and pick them out one by one.

library(ggplot2)
foo <- conditional_effects(mod1_all, categorical=TRUE)
plot(foo, plot = FALSE)[[1]] 

The last line you can save in a variable, i.e., p1 <- plot(foo, plot = FALSE)[[1]], and do that for each plot you want to combine into one big figure. I use patchwork for this: https://patchwork.data-imaginist.com

[GabriellaS-K]

@torkar

https://patchwork.data-imaginist.com

Thanks!
This works well:
foo ← conditional_effects(mod1_all, categorical=TRUE)
plot(foo, plot = FALSE)[[1]]

But I get stuck after, I can’t seem to combine them using patchwork, as ggplot won’t accept adding it to a ggplot object.

[torkar]

So this doesn’t work?

library(ggplot2)
library(patchwork)
foo ← conditional_effects(mod1_all, categorical=TRUE)
p1 <- plot(foo, plot = FALSE)[[1]]
p2 <- plot(foo, plot = FALSE)[[2]]
p1 + p2

I’m not 100% sure but I think I did something like that the last time I needed to. If that doesn’t work get back again and then I’ll check tomorrow when I’m at the laptop.

[scellus]

Hmm, conditional_effects seems to give a list of objects, of class brms_conditional_effects, and plotting them makes them ggplot2 plot objects. Sure you can make a list of them, etc., but the plot method that runs through the list and plots them all is gone, as it is for the class of conditional effects.

[GabriellaS-K]

No that doesn’t work, gives the error “Error: Can’t add p1 to a ggplot object.”

[joels]

It looks like conditional_effects() returns a list containing several plots (four plots in the example below). These can be laid out together using wrap_plots from the patchwork package:

library(brms)
library(tidyverse)
library(grid)
library(patchwork)

fit <- brm(count ~ zAge + zBase * Trt + (1 | patient),
           data = epilepsy, family = poisson()) 

# Get the plots generated by conditional_effects()
#  p1 is a list containing four elements. Each element is a plot.
p1 = plot(conditional_effects(fit), plot=FALSE)

# Lay out all four plots and add a title
wrap_plots(p1) + 
  plot_annotation(title="My Title", 
                  theme=theme(plot.title=element_text(hjust=0.5)))

wrap_plots(p1) is equivalent to “adding” the individual plots like this: p1[[1]] + p1[[2]] + p1[[3]] + p1[[4]]. However, wrap_plots is easier if you already have your plots in a list, as we do here.

Note that all of the plots have the same y-axis title. It would be nice to have a single title instead. Unfortunately, patchwork’s plot_annotation function doesn’t seem to have a way to add a general y-axis title. Here’s a more complicated approach where we create a separate object for the y-axis title and combine it with the four plots:

# Create a graphical object (grob) for the y-axis title
ytitle = textGrob("Count", rot=90)

# Remove y-axis title from individual plots
p2 = map(p1, ~.x + labs(y=NULL)) 

# Combine the four plots (inner wrap_plots) and then 
#  add the y-axis title (outer wrap_plots)
wrap_plots(list(ytitle, wrap_plots(p2))) + 
  plot_layout(ncol=2, widths=c(1,100)) +
  plot_annotation(title="My Title", 
                  theme=theme(plot.title=element_text(hjust=0.5)))

[torkar]

And here’s more info:

David Issa Mattos | Plot Panels in R (from brms models)

[GabriellaS-K]

This worked perfectly!! Thank you so much!

-a smaller question- do you know how to change the legend values? It is so easy in ggplot to rename (1,2,3 mean very little for example, I would recode to negative, positive, neutral).

Thanks again

[GabriellaS-K]

and combine, so as not to have 4 of the same legends, but just one. I have a feeling this will be a little e more complicated

[GabriellaS-K]

Just reading the link you sent-i think all the answers are there! Will come back if there are any issues. Thanks for this!

[GabriellaS-K]

There are issues! I can remove the legend, change colour of the writing, but I cannot change the values of the legend, it just adds a second legend (see attached). I would like 1 legend for all my graphs (sentiment: positive, negative, neutral) instead of one on each.

[GabriellaS-K]

Now using this code (after having fixed each graph individually)

a <- plot(conditional_effects(mod1, effects="age",categorical = T))[[1]] + 
  coord_cartesian(ylim = c(0,.82)) +
  labs(x="Age",
       title = "Age",
       fill = "Age",
       colour = "Age")+
  theme(axis.title.x=element_blank()) +
    theme(axis.title.y = element_blank())+
  theme(legend.position = "none")+
  theme(axis.text.x = element_text(color="black", 
                           size=4))

The above is done four times, for graph a,b,c,d.
Next, these two seem to work:

Plot <- (a+ b+
c+d) +
    plot_layout(nrow=2)+theme(legend.position = 'bottom')+ scale_colour_discrete(name="Sentiment",
                         labels=c("Negative", "Neutral", "Positive"))


m1allplot<-wrap_plots(list(ytitle, wrap_plots(Plot))) + 
  plot_layout(ncol=2, widths=c(1,100)) +
plot_annotation(title="Effects", 
                   theme=theme(plot.title=element_text(hjust=0.5)))

This makes the attached set of graphs (grey boxes over results as it is sensitive data)

[joels]

Each element of the list returned by conditional_effects is a ggplot, so you can use any ggplot function on each plot. For example:

# Create list of plots
p = plot(conditional_effects(mod1, effects="age",categorical = T))

# Format the first plot in the list
p[[1]] = p[[1]] + 
  coord_cartesian(ylim = c(0,.82)) +
  labs(x="Age",
       title = "Age",
       fill = "Age",
       colour = "Age")+
  theme(axis.title.x=element_blank(),
        axis.title.y = element_blank(),
        legend.position = "none",
        axis.text.x = element_text(color="black", size=4))

If you want to do the same thing to more than one plot, you can use the map function (from the tidyr package).

Format plots 1 and 3 (elements 1 and 3 of p):

p[c(1,3)] = map(p[c(1,3)],  
                ~ .x + coord_cartesian(ylim = c(0,.82)) +
                  theme(axis.title.x=element_blank(),
                        axis.title.y = element_blank(),
                        legend.position = "none",
                        axis.text.x = element_text(color="black", size=4))

Format all plots:

p = map(p,  
        ~ .x + coord_cartesian(ylim = c(0,.82)) +
            theme(axis.title.x=element_blank(),
                  axis.title.y = element_blank(),
                  legend.position = "none",
                  axis.text.x = element_text(color="black", size=4))

And you can, of course, apply common formatting to several or all plots, and custom formatting to any one of them. Note also in the examples above, that (1) all theme changes can be done in a single call to the theme function, and (2) we only need to run plot(conditional_effects(mod1, effects="age",categorical = T)) once. This gives us the list of plots p that we can then modify as desired.

Regarding the legend: It looks like conditional_effects generates plot “D” (in your example) with both colour and fill aesthetics. To keep those legends combined as a single legend, they both need to have the same name and the same labels. So, in your example, it would probably be:

# Create list of plots
p = plot(conditional_effects(mod1, effects="age",categorical = T))

p[[4]] = p[[4]] +
      scale_colour_discrete(name="Sentiment", 
                            labels=c("Negative", "Neutral", "Positive")) +
      scale_fill_discrete(name="Sentiment",
                          labels=c("Negative", "Neutral", "Positive"))

Also, regarding axis titles, in your example labs(x="Age") sets the x-axis title to “Age”, but then theme(axis.title.x=element_blank()) removes the x-axis title.

[GabriellaS-K]

Thank you so so much for such a comprehensive answer!! It all makes sense, and my graphs look great! Thanks