R/calculate_and_plot_cpue.R

In FESta: Fishing Effort Standardization

  #=========================================================
  # Global variable declarations for R CMD check
  # ========================================================
utils::globalVariables(
  c(
    "xmin",
    "xmax",
    "ymin",
    "ymax"

  )
)
#' Calculate Nominal CPUE and Plot with Standardized CPUE
#'
#' Calculates Nominal CPUE from raw catch and effort data, produces a
#' summary table, and generates a dual-axis plot showing total catch (bars),
#' Nominal CPUE, and standardized CPUE.
#'
#' @param year Numeric vector of year values (one per row).
#' @param catch Numeric vector of catch values (one per row).
#' @param effort Numeric vector of effort values (one per row).
#' @param std_cpue Numeric vector of standardized CPUE, one value per year
#'   (same length as unique years, in chronological order).
#'
#' @return A data frame (invisibly) with columns: Year, Total_Catch,
#'   Nominal_CPUE, Standardized_CPUE. The summary table is printed and
#'   the plot is rendered.
#' @keywords internal
#' @examples
#' \dontrun{
#' calculate_and_plot_cpue(
#'   year     = mydata$Year,
#'   catch    = mydata$Catch,
#'   effort   = mydata$Effort,
#'   std_cpue = glm_result$standardized_index$Standardized_CPUE
#' )
#' }
#'@noRd
calculate_and_plot_cpue <- function(year, catch, total_catch = NULL, effort, std_cpue, nom_cpue = NULL,log_transform=NULL,aic=NULL,bic=NULL) {

  # ---- Validation ----
  n <- length(year)

  if (any(effort <= 0, na.rm = TRUE))
    stop("All 'effort' values must be positive.")

  # ---- Step 1: Build row-level data frame ----
  df <- data.frame(
    Year   = as.integer(year),
    Catch  = as.numeric(catch),
    Effort = as.numeric(effort)
  )
  yearly_efforts <- stats::aggregate(Effort ~ Year, data = df, FUN = sum)
  names(yearly_efforts) <- c("Year", "Total_Effort")
  # ---- Step 2: Compute yearly total catch ----
  if (is.null(total_catch)) {
    yearly_totals <- stats::aggregate(Catch ~ Year, data = df, FUN = sum)
  } else {
    yearly_totals <- data.frame(
      Year        = yearly_efforts$Year,
      Total_Catch = total_catch
    )
  }
  names(yearly_totals) <- c("Year", "Total_Catch")

  # ---- Step 3: Build base summary table ----
  summary_df <- merge(yearly_totals, yearly_efforts, by = "Year", all = TRUE)
  summary_df <- summary_df[order(summary_df$Year), ]

  # ---- Step 4: Handle Nominal CPUE ----
  if (is.null(nom_cpue)) {
    summary_df$Nominal_CPUE<-summary_df$Total_Catch/summary_df$Total_Effort
    # df <- merge(df, yearly_totals, by = "Year", all.x = TRUE)
    # df <- merge(df, yearly_efforts, by = "Year", all.x = TRUE)

    #df$Nominal_CPUE <- (df$Catch / df$Total_Catch) * df$Total_Effort


    # norm_by_year <- stats::aggregate(
    #   Nominal_CPUE ~ Year,
    #   data = df,
    #   FUN = mean
    # )

    # summary_df <- merge(
    #   summary_df,
    #   norm_by_year,
    #
    #   by = "Year",
    #   all = TRUE
    # )

  } else {
    if(length(nom_cpue) != nrow(summary_df)) {
      stop(sprintf(
        "'nom_cpue' must contain %d values, one per unique year. Got %d.",
        nrow(summary_df), length(nom_cpue)
      ))
    }
    summary_df$Nominal_CPUE <- as.numeric(nom_cpue)

  }
  # if(is.null(log_transform)==FALSE){
  # if (log_transform) {
  #   summary_df$Nominal_CPUE <- exp(summary_df$Nominal_CPUE)
  # }else{
  #   summary_df$Nominal_CPUE <- as.numeric(nom_cpue)
  # }
  # }
  if (isTRUE(log_transform)) {
    summary_df$Nominal_CPUE <- exp(summary_df$Nominal_CPUE)
  }
  # ---- Step 5: Attach standardized CPUE ----
  if (length(std_cpue) != nrow(summary_df)) {
    stop(sprintf(
      "'std_cpue' must contain %d values, one per unique year. Got %d.",
      nrow(summary_df), length(std_cpue)
    ))
  }
  summary_df$Standardized_CPUE <- as.numeric(std_cpue)
  # ---- Step 5b: Attach model selection criteria ----
  if (!is.null(aic) && !is.null(bic)){
    lm=list(AIC=round(as.numeric(aic), 4),BIC=round(as.numeric(bic), 4))
    print(lm)
  }
  # ---- Step 6: Print table ----
  cat("\n", strrep("=", 60), "\n", sep = "")
  cat("  CPUE Summary Table\n")
  cat(strrep("=", 60), "\n\n", sep = "")
  row.names(summary_df) <- NULL
  summary_df<-round(summary_df, digits = 4L)

  # ---- Step 7: Plot ----

  col_catch <- "orange"
  col_norm  <- "#6a4c93"
  col_std   <- "skyblue3"

  yr <- summary_df$Year
  tc <- summary_df$Total_Catch
  nc <- summary_df$Nominal_CPUE
  sc <- summary_df$Standardized_CPUE

  bad_rows <- summary_df[!is.finite(summary_df$Total_Catch) |
                           !is.finite(summary_df$Nominal_CPUE)   |
                           !is.finite(summary_df$Standardized_CPUE), ]
  if (nrow(bad_rows) > 0) {
    print(bad_rows)
    stop("[calculate_and_plot_cpue] Non-finite value(s) found in the row(s) above ",
         "(Total_Catch / Nominal_CPUE / Standardized_CPUE). Cannot build plot.")
  }

  if(isTRUE(log_transform)==TRUE){
    if (any(nc <= 0, na.rm = TRUE)) {
      stop("[calculate_and_plot_cpue] Nominal_CPUE contains non-positive value(s); ",
           "cannot take log() for plotting. Check Total_Catch/Total_Effort computation ",
           "(a year with zero total catch would produce this).")
    }
    summary_df$Nominal_CPUE<-log(nc)
    nc<-log(nc)
  }
  # ── x-axis breaks ─────────────────────────────────────────────────────
  all_yr  <- sort(unique(as.integer(yr)))
  yr_brks <- all_yr[seq(1L, length(all_yr), by = 2)]

  # ── Internal panel builder ─────────────────────────────────────────────

  .build_cpue_panel <- function(yr, tc, cpue, cpue_colour,
                                cpue_label, cpue_linetype, cpue_shape,
                                panel_title,log_transform = NULL) {
    if(isTRUE(log_transform)==TRUE){
      cpue_name<-"log(CPUE)"
    }else{
      cpue_name<-"CPUE"
    }

    all_yr_local <- sort(unique(as.integer(yr)))

    cpue_min  <- min(cpue, na.rm = TRUE)
    cpue_max  <- max(cpue, na.rm = TRUE)
    catch_max <- max(tc,   na.rm = TRUE)

    cpue_breaks  <- pretty(c(cpue_min, cpue_max), n = 6)
    y_range      <- range(cpue_breaks)
    y_pad        <- 0.15 * diff(y_range)
    y_limit_low  <- y_range[1] - y_pad
    y_limit_high <- y_range[2] + y_pad

    catch_breaks <- pretty(c(0, catch_max), n = 6)
    catch_breaks <- catch_breaks[catch_breaks >= 0]
    catch_limit  <- max(catch_breaks) * 1.15

    catch_scale_f <- (y_limit_high - y_limit_low) / catch_limit
    tc_scaled     <- tc * catch_scale_f + y_limit_low

    # geom_bar(stat="identity") ALWAYS baselines bars at y = 0 internally,
    # regardless of tc_scaled's value or the axis limits. When the CPUE
    # axis doesn't span across 0 (e.g. all-negative log-CPUE, or all-positive
    # but far from 0), that hidden 0 baseline falls outside
    # scale_y_continuous(limits=...), gets censored to NA, and ggplot2 drops
    # the whole bar ("Removed N rows" warning, no bars rendered). geom_rect()
    # lets us anchor the bar's bottom explicitly at y_limit_low -- the real
    # bottom of the visible axis -- instead of an implicit, possibly
    # off-screen, zero.
    half_width <- 0.3
    plot_df <- data.frame(
      Year   = yr,
      xmin   = yr - half_width,
      xmax   = yr + half_width,
      ymin   = y_limit_low,
      ymax   = tc_scaled,
      tc_scaled = tc_scaled,
      CPUE   = cpue
    )

    min_year <- min(all_yr_local)
    max_year <- max(all_yr_local)
    yr_brks  <- seq(min_year, max_year + 1, by = 2)

    ggplot2::ggplot(plot_df, ggplot2::aes(x = Year)) +
      ggplot2::geom_rect(
        ggplot2::aes(xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax,
                     fill = "Total Catch"),
        alpha = 0.6
      ) +
      ggplot2::geom_line(
        ggplot2::aes(y = CPUE, colour = cpue_label),
        linewidth = 1.1, linetype = cpue_linetype
      ) +
      ggplot2::geom_point(
        ggplot2::aes(y = CPUE, colour = cpue_label),
        shape = cpue_shape, fill = "white", size = 2.5, stroke = 1.3
      ) +
      ggplot2::scale_y_continuous(
        name   = cpue_name,
        labels = scales::number_format(accuracy = 0.01),
        expand = ggplot2::expansion(mult = c(0, 0)),
        limits = c(y_limit_low, y_limit_high),
        breaks = cpue_breaks,
        sec.axis = ggplot2::sec_axis(
          ~ (. - y_limit_low) / catch_scale_f,
          name   = "Total Catch",
          breaks = catch_breaks,
          labels = scales::comma(catch_breaks)
        )
      ) +
      ggplot2::scale_x_continuous(
        breaks = yr_brks,
        labels = as.character(yr_brks),
        limits = c(min_year - 0.5, max_year + 1.5),
        expand = ggplot2::expansion(mult = c(0, 0))
      ) +
      ggplot2::scale_colour_manual(
        name   = NULL,
        values = stats::setNames(cpue_colour, cpue_label),
        guide  = ggplot2::guide_legend(
          order = 2,
          override.aes = list(linetype = cpue_linetype, shape = cpue_shape, fill = "white")
        )
      ) +
      ggplot2::scale_fill_manual(
        name = NULL,
        values = c("Total Catch" = col_catch),
        guide = ggplot2::guide_legend(order = 1)
      ) +
      ggplot2::labs(title = panel_title, x = "Year") +
      ggplot2::theme_minimal(base_size = 11) +
      ggplot2::theme(
        plot.title         = ggplot2::element_text(face = "bold", hjust = 0.5, size = 12),
        axis.title.y       = ggplot2::element_text(colour = cpue_colour, face = "bold", size = 10),
        axis.text.y        = ggplot2::element_text(colour = "black", size = 9),
        axis.title.y.right = ggplot2::element_text(colour = col_catch, face = "bold", size = 10,angle=90),
        axis.text.y.right  = ggplot2::element_text(colour = "black", size = 9),
        axis.text.x        = ggplot2::element_text(colour = "black", size = 9, angle = 45, hjust = 1),
        axis.title.x       = ggplot2::element_text(size = 10),
        legend.position    = "bottom",
        legend.key.width   = grid::unit(1.4, "cm"),
        legend.text        = ggplot2::element_text(size = 9),
        panel.grid.minor   = ggplot2::element_blank(),
        panel.grid.major   = ggplot2::element_line(colour = "grey92", linewidth = 0.35),
        plot.margin        = ggplot2::margin(8, 8, 6, 8)
      )
  }

  # ── Build two panels ──────────────────────────────────────────────────
  p1 <- .build_cpue_panel(
    yr            = yr,
    tc            = tc,
    cpue          = nc,
    cpue_colour   = col_norm,
    cpue_label    = "Nominal CPUE",
    cpue_linetype = "solid",
    cpue_shape    = 21,
    panel_title   = "Nominal CPUE vs Total Catch",
    log_transform = log_transform
  )

  p2 <- .build_cpue_panel(
    yr            = yr,
    tc            = tc,
    cpue          = sc,
    cpue_colour   = col_std,
    cpue_label    = "Standardized CPUE",
    cpue_linetype = "solid",
    cpue_shape    = 21,
    panel_title   = "Standardized CPUE vs Total Catch",
    log_transform = log_transform
  )

  # ── Combine into one frame ────────────────────────────────────────────
  if (requireNamespace("patchwork", quietly = TRUE)) {
    pw <- patchwork::wrap_plots(p1, p2, ncol = 2) +
      patchwork::plot_annotation(
        title = "CPUE vs Total Catch",
        theme = ggplot2::theme(
          plot.title = ggplot2::element_text(face = "bold",
                                             hjust = 0.5,
                                             size = 14)
        )
      )
    print(pw)
  } else {
    grobs <- lapply(list(p1, p2), ggplot2::ggplotGrob)
    arr   <- gridExtra::arrangeGrob(
      grobs = grobs, ncol = 2,
      top   = grid::textGrob(
        "CPUE vs Total Catch",
        gp = grid::gpar(fontface = "bold", fontsize = 13)
      )
    )
    grid::grid.newpage()
    grid::grid.draw(arr)
  }

  # ---- Return ----

  return(summary_df[,-3])
}