R/plot_raincloud.R

In DUToolkit: Visualizing and Quantifying Decision Uncertainty

#' Generate raincloud plot
#'
#' @description The probability densities of the highest (or lowest if the
#' threshold is a minimum) projected outcome across simulation runs are
#' plotted for each policy alternative alongside corresponding box plots,
#' which indicate the mean and 50th percentile range. These plots are
#' presented collectively on a single graph to facilitate visual comparison
#' of the policy alternatives. The decision threshold is shown directly on
#' the plot as a vertical line to provide a clear reference point for
#' interpreting the outputs.
#'
#' @param max_min_values_list A list generated by [get_max_min_values()] that
#' must contain at least two elements.
#' @param D A single threshold value.
#'
#' @return A raincloud plot for peak (or minimum) values
#' @export
#'
#' @examples
#' tmin <- "2021-01-01"
#' tmax <- "2021-04-10"
#' D <- 750
#'
#' peak_values_list <- get_max_min_values(
#'   psa_data,
#'   tmin = tmin,
#'   tmax = tmax,
#'   Dt_max = TRUE
#' )
#'
#' plot_raincloud(
#'   peak_values_list,
#'   D = D
#' )
plot_raincloud <- function(max_min_values_list, D) {
  if (!inherits(max_min_values_list, "list") ||
    length(max_min_values_list) <= 1) {
    rlang::abort("The first argument must be a list with at least two elements.",
      class = "data_type"
    )
  }
  raincloud_df <- c()
  for (i in 1:length(max_min_values_list)) {
    if (!is.null(names(max_min_values_list)) && all(names(max_min_values_list) != "")) {
      max_min_values_list[[i]]$strategy <- paste(names(max_min_values_list)[i])
    } else {
      max_min_values_list[[i]]$strategy <- paste0("Intervention_", i)
    }
    raincloud_df <- rbind(raincloud_df, max_min_values_list[[i]])
  }

  raincloud_plot <- ggplot2::ggplot(
    data = raincloud_df,
    ggplot2::aes(
      x = factor(!!rlang::sym("strategy")),
      y = !!rlang::sym("outcome"),
      fill = factor(!!rlang::sym("strategy"))
    )
  ) +
    ggplot2::scale_fill_manual(values = c(
      "#6C568CFF", "#9386A6FF",
      "#BFCDD9FF", "#7F8C72FF", "#607345FF",
      "#6C568CFF", "#9386A6FF",
      "#BFCDD9FF", "#7F8C72FF", "#607345FF",
      "#6C568CFF", "#9386A6FF",
      "#BFCDD9FF", "#7F8C72FF", "#607345FF",
      "#6C568CFF", "#9386A6FF",
      "#BFCDD9FF", "#7F8C72FF", "#607345FF"
    )) +
    # scale_fill_paletteer_d("calecopal::lupinus") +
    # scale_fill_grey()+
    ggdist::stat_halfeye(
      point_colour = NA,
      justification = -0.1,
      # set slab interval to 95% data range
      .width = 0,
      # aes(thickness = after_stat(pdf*n)),
      scale = 0.7, slab_color = "black", slab_linewidth = 0.5
    ) +
    ggplot2::geom_boxplot(
      width = .1,
      outlier.shape = NA,
      alpha = 0.2
    ) +
    ggplot2::geom_hline(yintercept = D, color = "red", linetype = "dashed", linewidth = 0.75) +
    # stat_dots(side = "left",justification = 1.03) +
    ggplot2::theme_classic() +
    ggplot2::labs(y = "Outcome at peak") +
    ggplot2::theme(legend.position = "none", axis.title.y = ggplot2::element_blank()) +
    ggplot2::coord_flip()

  raincloud_plot <- gen_stand_descr(raincloud_plot,
                                    link = "Raincloud plots")
  return(raincloud_plot)
}