Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Generated biomass plot #7

Merged
merged 6 commits into from
Sep 27, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -4,4 +4,5 @@ export(add_theme)
export(plot_landings)
export(plot_recruitment)
export(plot_spawning_biomass)
export(plot_total_biomass)
export(table_indices)
120 changes: 120 additions & 0 deletions R/plot_biomass.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,120 @@
#' Plot Total Biomass
#'
#' @inheritParamsplot_recruitment
#' @param show_warnings Option to suppress warnings
#' @param units units for biomass
#' @param scaled TRUE/FALSE; indicate whether the output values for biomass and recruitment are scaled
#' @param scale_amount indicate the exact amount of scale (i.e. 1000)
#' @param ref_line choose with reference point to plot a reference line and use
#' in relative totb calculations
#' @param end_year input the end year of the stock assessment data (not including
#' projections). This parameter will be deprecated once the output converter is fully developed.
#' @param relative Plot relative total biomass. Ref line indicates which reference point to use
#'
#' @return Plot total biomass from a stock assessment model as found in a NOAA
#' stock assessment report. Units of total biomass can either be manually added
#' or will be extracted from the provided file if possible. In later releases, model will not
#' @export
#'
plot_total_biomass <- function(dat,
model = "standard",
show_warnings = FALSE,
units = NULL,
scaled = FALSE,
scale_amount = 1000,
ref_line = c("target", "MSY", "msy", "unfished"),
end_year = NULL,
relative = FALSE
){

if(length(ref_line)>1){
ref_line = "target"
} else {
ref_line <- match.arg(ref_line, several.ok = FALSE)
}

# check units
if(!is.null(units)){
bu <- units
} else {
bu <- "metric tons"
}

if(model == "standard"){
output <- read.csv(dat)
totb <- output |>
dplyr::filter(label == "biomass",
module_name == "DERIVED_QUANTITIES" | module_name == "t.series") |> # SS3 and BAM target module names
dplyr::mutate(estimate = as.numeric(estimate),
year = as.numeric(year))
if (is.null(end_year)){
endyr <- max(totb$year)
}
# Select value for reference line and label
# update the target option later
if (any(grepl("target", output$label))) {
ref_line_val <- as.numeric(output[grep("(?=.*biomass)(?=.*target)", output$label, perl = TRUE),]$estimate)
ref_line_label <- "target"
if (scaled) {
ann_add <- ggplot2::annotate("text", x = endyr + 0.05, y=ref_line_val, label = bquote(B[target])) # this might need to change
} else {
ann_add <- ggplot2::annotate("text", x = endyr + 0.05, y=ref_line_val/1000, label = bquote(B[target]))
}
} else if (ref_line == "MSY" | ref_line == "msy") {
ref_line_val <- as.numeric(output[grep("(^biomass_msy)", output$label, perl = TRUE),]$estimate)
ref_line_label <- "MSY"
if (scaled) {
ann_add <- ggplot2::annotate("text", x = endyr + 0.05, y=ref_line_val, label = bquote(B[ref_line]))
} else {
ann_add <- ggplot2::annotate("text", x = endyr + 0.05, y=ref_line_val/1000, label = bquote(B[MSY]))
}
} else if (ref_line == "unfished") {
ref_line_val <- as.numeric(output[grep("(^biomass_unfished)", output$label, perl = TRUE),]$estimate)
ref_line_label <- "unfished"
if (scaled) {
ann_add <- ggplot2::annotate("text", x = endyr + 0.05, y=ref_line_val, label = bquote(B[unfished]))
} else {
ann_add <- ggplot2::annotate("text", x = endyr + 0.05, y=ref_line_val/1000, label = bquote(B[unfished]))
}
}
# Choose number of breaks for x-axis
x_n_breaks <- round(length(subset(totb, year<=endyr)$year)/10)
if (x_n_breaks <= 5) {
x_n_breaks <- round(length(subset(totb, year<=endyr)$year)/5)
}
if (relative) {
# plot relative TOTB
plt <- ggplot2::ggplot(data = subset(totb, year<=endyr)) +
ggplot2::geom_line(ggplot2::aes(x = year, y = estimate/ref_line_val), linewidth = 1) +
# ggplot2::geom_ribbon(ggplot2::aes(x = year, ymin = (value/ref_line_val - stddev/ref_line_val), ymax = (value/ref_line_val + stddev/ref_line_val)), colour = "grey", alpha = 0.3) +
ggplot2::geom_hline(yintercept = ref_line_val/ref_line_val, linetype = 2) +
ggplot2::labs(x = "Year",
y = paste("Biomass (", bu, ")", sep = "")) +
ggplot2::scale_x_continuous(n.breaks = x_n_breaks,
guide = ggplot2::guide_axis(minor.ticks = TRUE))
} else {
if(scaled){
plt <- ggplot2::ggplot(data = totb) +
ggplot2::geom_line(ggplot2::aes(x = year, y = estimate), linewidth = 1) +
ggplot2::geom_hline(yintercept = ref_line_val, linetype = 2) +
ggplot2::labs(x = "Year",
y = paste("Biomass (", bu, ")", sep = "")) +
ggplot2::scale_x_continuous(n.breaks = x_n_breaks,
guide = ggplot2::guide_axis(minor.ticks = TRUE))
plt <- plt + ann_add
} else {
plt <- ggplot2::ggplot(data = totb) +
ggplot2::geom_line(ggplot2::aes(x = year, y = estimate/1000), linewidth = 1) +
# ggplot2::geom_ribbon(ggplot2::aes(x = year, ymin = (value/1000 - stddev/1000), ymax = (value/1000 + stddev/1000)), colour = "grey", alpha = 0.3) +
ggplot2::geom_hline(yintercept = ref_line_val/1000, linetype = 2) +
ggplot2::labs(x = "Year",
y = paste("Biomass (", bu, ")", sep = "")) +
ggplot2::scale_x_continuous(n.breaks = x_n_breaks,
guide = ggplot2::guide_axis(minor.ticks = TRUE))
plt <- plt + ann_add
}
}
plt_fin <- add_theme(plt)
}
return(plt_fin)
}
6 changes: 1 addition & 5 deletions R/plot_recruitment.R
Original file line number Diff line number Diff line change
@@ -1,10 +1,6 @@
#' Plot Recruitment
#'
#' @param dat Stock assessment output file containing estimates parameters and
#' other associated metrics.
#' @param model Acryonym of the stock assessment model that produced the output file.
#' Currently, only stock synthesis (SS3) and Beaufort Assessment Model (BAM)
#' are options for plotting spawning biomass.
#' @inheritParams plot_recruitment
#' @param params Print/export the parameters of the stock recruitment function?
#' @param params_only Only export the stock recruitment function or both the parameters and the plot(s)?
#' @param units If units are not available in the output file, in metric tons,
Expand Down
43 changes: 43 additions & 0 deletions man/plot_total_biomass.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

Loading