Plot all plots created in a function with the minimum y axis range necessary to show confidence interval (r)

Question

Suppose I have code to generate and save plots of multiple gam() models such as the following:

library(mgcv)
library(gratia)
library(tidyverse)

iterate = function(z){

    model1 <- gam(mpg ~ cyl + s(disp, k = 25) + wt, 
        data = mtcars %>% dplyr::filter(carb == z), method = "REML")

    gam_plot <- draw(model1, ci_col = "red", smooth_col = "red")
    assign(paste0("gam_plot", z), gam_plot, envir = .GlobalEnv)
}

z_list = unique(mtcars$gear)

sapply(X = z_list,FUN = iterate)

After generating and saving all of these plots, I want to plot all the plots together using ggarrange(), but because the output is the same, I need to use particular limits on the y-axis that both are not too narrow that the shaded confidence interval disappears, but also not so wide that the plots become unintelligible. What is the most efficient way to set the scale_y_continuous of all plots to be the narrowest possible that contains all of the plots without stopping the confidence intervals from being shaded? Does this need to be calculated through a function prior to this function used to create and plot the plots, or can it all be done in one function run?

Answer 1

One possible option would be to create your lists of plots first. In a second step get the ranges of the y scale for which I use layer_scales and compute the range of the "ranges". Finally apply the combined range to each of your plots before combining them:

library(mgcv)
#> Loading required package: nlme
#> This is mgcv 1.9-1. For overview type 'help("mgcv-package")'.
library(gratia)
library(tidyverse)

iterate <- function(z) {
  model1 <- gam(mpg ~ cyl +
    s(disp, k = 3) + wt, data = mtcars %>%
    dplyr::filter(gear == z), method = "REML")


  draw(model1,
    ci_col = "red",
    smooth_col = "red"
  )
}

z_list <- unique(mtcars$gear)

plot_list <- lapply(z_list,iterate)
#> Warning in newton(lsp = lsp, X = G$X, y = G$y, Eb = G$Eb, UrS = G$UrS, L = G$L,
#> : Fitting terminated with step failure - check results carefully

range_y <- plot_list |> 
  lapply(
    \(x) layer_scales(x)$y$range$range
  ) |> 
  unlist() |> 
  range()

plot_list |> 
  lapply(\(x) {
    x +
      scale_y_continuous(limits = range_y)
  }) |> 
  ggpubr::ggarrange(plotlist = _, nrow = 1)