diff --git a/.Rbuildignore b/.Rbuildignore index a32bacb..d5655a0 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -2,3 +2,7 @@ ^\.Rproj\.user$ ^doc$ ^Meta$ +analysis +archive +past_years +plots diff --git a/0_update_condition_data.Rmd b/0_update_condition_data.Rmd index 2a09674..67d9263 100644 --- a/0_update_condition_data.Rmd +++ b/0_update_condition_data.Rmd @@ -114,17 +114,17 @@ ESP_SETTINGS <- list(ESP_SPECIES = n_knots = c(400, 400, 400, 600, 750, 750), - fl_min = c(0, 421, - 0, 461, - 0, 461), - fl_max = c(420, 1e7, - 460, 1e7, - 460, 1e7)) + fl_min = c(0, 504, + 0, 581, + 0, 581), + fl_max = c(503, 1e7, + 580, 1e7, + 580, 1e7)) ) ``` -# 1. Update stratum-biomass weighted condition inidicators +# 1. Update stratum-biomass weighted condition indicators ```{r get_ai, message=TRUE, warning=TRUE} # Aleutian Islands @@ -277,9 +277,7 @@ Update sysdata.rda with raw data for the condition indicator and write raw data ```{r save_to_sysdata, message=FALSE, warning=FALSE} EBS_INDICATOR <- list( FULL_REGION = as.data.frame( - # dplyr::full_join( - ebs_sbw$full_sbw) %>% #, - # ebs_vast_df)) %>% + ebs_sbw$full_sbw) |> dplyr::filter(common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$EBS]), STRATUM = as.data.frame( dplyr::filter( @@ -291,9 +289,7 @@ EBS_INDICATOR <- list( NBS_INDICATOR <- list( FULL_REGION = as.data.frame( - # dplyr::full_join( - nbs_sbw$full_sbw) %>% #, - # nbs_vast_df) %>% + nbs_sbw$full_sbw) |> dplyr::filter( common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$NBS]), STRATUM = as.data.frame( @@ -306,9 +302,7 @@ NBS_INDICATOR <- list( GOA_INDICATOR <- list( FULL_REGION = as.data.frame( - # dplyr::full_join( - goa_sbw$full_sbw) %>% #, - # goa_vast_df)) %>% + goa_sbw$full_sbw) |> dplyr::filter(common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$GOA]), STRATUM = as.data.frame( dplyr::filter( @@ -320,10 +314,7 @@ GOA_INDICATOR <- list( AI_INDICATOR <- list( FULL_REGION = as.data.frame( - # dplyr::full_join( - ai_sbw$full_sbw) %>% #, - # ai_vast_df) - # ) %>% + ai_sbw$full_sbw) |> dplyr::filter(common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$AI]), STRATUM = as.data.frame( dplyr::filter( @@ -335,22 +326,18 @@ AI_INDICATOR <- list( PCOD_ESP <- list( FULL_REGION_EBS = as.data.frame( - # dplyr::full_join(ebs_sbw$full_sbw, ebs_vast_df) %>% dplyr::filter(ebs_sbw$full_sbw, common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name )), FULL_REGION_GOA = as.data.frame( - # dplyr::full_join(goa_sbw$full_sbw, goa_vast_df) %>% dplyr::filter(goa_sbw$full_sbw, common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name )), FULL_REGION_AI = as.data.frame( - # dplyr::full_join(ai_sbw$full_sbw, ai_vast_df) %>% dplyr::filter(ai_sbw$full_sbw, common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name )), FULL_REGION_NBS = as.data.frame( - # dplyr::full_join(nbs_sbw$full_sbw, nbs_vast_df) %>% dplyr::filter( nbs_sbw$full_sbw, common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name )), diff --git a/0_update_indicator.R b/0_update_indicator.R new file mode 100644 index 0000000..87fa4a3 --- /dev/null +++ b/0_update_indicator.R @@ -0,0 +1,253 @@ +library(akfishcondition) + +# Make ESR/ESP settings ---- + +ESR_SETTINGS <- list(ESR_SPECIES = data.frame(common_name = c( + "walleye pollock", "walleye pollock (100-250 mm)", "walleye pollock (>250 mm)", "Pacific cod", + "Pacific cod (juvenile)", "Pacific cod (adult)", "Atka mackerel", "arrowtooth flounder", + "flathead sole", "yellowfin sole", "northern rock sole", "southern rock sole", "Alaska plaice", + "Pacific ocean perch", "dusky rockfish", "northern rockfish", "Dover sole", "rex sole", "shortraker rockfish", "rougheye rockfish", "blackspotted rockfish", "sharpchin rockfish"), + species_code = c(21740, 21741, 21742, 21720, 21721, 21722, 21921, 10110, 10130, + 10210, 10261, 10262, 10285, 30060, 30152, 30420, 10180, 10200, 30576, 30051, 30052, 30560), + AI = c(FALSE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE, FALSE, FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE), + GOA = c(FALSE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, TRUE, TRUE, FALSE, TRUE, TRUE, FALSE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE), + EBS = c(FALSE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, TRUE, TRUE, TRUE, TRUE, FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE), + NBS = c(FALSE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE)), + VAST_SETTINGS = data.frame(species_code = c(30420, 30060, 10262, 21720, 21740, 10110, 30152, + 30420, 30060, 10262, 21720, 21740, 21921, 10110, 21741, 21742, + 10130, 10210, 10285, 21740, 21720, 10110, 10261, 21741, 21742, + 21740, 21720, 10210, 10285, 21741, 21742 + ), + region = c("GOA", "GOA", "GOA", "GOA", "GOA", "GOA", "GOA", + "AI", "AI", "AI", "AI", "AI", "AI", "AI", "AI", "AI", + "EBS", "EBS", "EBS", "EBS", "EBS", "EBS", "EBS", "EBS", "EBS", + "NBS", "NBS", "NBS", "NBS", "NBS", "NBS"), + ObsModel_1 = c(2, 2, 2, 2, 2, 2, 2, + 2, 2, 2, 2, 2, 2, 2, 2, 2, + 2, 2, 2, 2, 2, 2, 2, 2, 2, + 2, 2, 2, 2, 2, 2), + ObsModel_2 = c(1, 1, 1, 1, 1, 1, 1, + 1, 1, 1, 1, 1, 1, 1, 1, 1, + 1, 1, 1, 1, 1, 1, 1, 1, 1, + 1, 1, 1, 1, 1, 1), + ObsModel_3 = c(3, 3, 3, 3, 4, 4, 4, + 3, 3, 3, 3, 4, 3, 4, 4, 4, + 3, 3, 3, 4, 3, 4, 3, 4, 4, + 4, 4, 3, 4, 4, 4), + ObsModel_4 = c(3, 3, 3, 3, 4, 4, 4, + 3, 3, 3, 3, 4, 3, 4, 4, 4, + 3, 3, 3, 4, 3, 4, 3, 4, 4, + 4, 4, 3, 4, 4, 4), + fl_min = c(0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 100, 251, + 0, 0, 0, 0, 0, 0, 0, 100, 251, + 0, 0, 0, 0, 100, 251), + fl_max = c(1e7, 1e7, 1e7, 1e7, 1e7, 1e7, 1e7, + 1e7, 1e7, 1e7, 1e7, 1e7, 1e7, 1e7, 250, 1e7, + 1e7, 1e7, 1e7, 1e7, 1e7, 1e7, 1e7, 250, 1e7, + 1e7, 1e7, 1e7, 1e7, 250, 1e7), + n_knots = c(400, 400, 400, 400, 400, 400, 400, + 400, 600, 400, 400, 400, 400, 400, 400, 400, + 400, 600, 400, 400, 750, 400, 400, 400, 400, + 200, 200, 200, 200, 200, 200)) +) + +ESP_SETTINGS <- list(ESP_SPECIES = + data.frame( + common_name = c("Pacific cod (juvenile)", "Pacific cod (adult)"), + species_code = c(21721, 21722), + AI = c(TRUE, TRUE), + GOA = c(TRUE, TRUE), + EBS = c(TRUE, TRUE), + NBS = c(FALSE, FALSE)), + VAST_SETTINGS = data.frame(species_code = c(21721, 21722, + 21721, 21722, + 21721, 21722), + region = c("GOA", "GOA", + "AI", "AI", + "EBS", "EBS"), + ObsModel_1 = c(2, 2, + 2, 2, + 2, 2), + ObsModel_2 = c(1, 1, + 1, 1, + 1, 1), + ObsModel_3 = c(3, 3, + 4, 4, + 3, 3), + ObsModel_4 = c(3, 3, + 4, 4, + 3, 3), + n_knots = c(400, 400, + 400, 600, + 750, 750), + fl_min = c(0, 504, + 0, 581, + 0, 581), + fl_max = c(503, 1e7, + 580, 1e7, + 580, 1e7)) +) + +devtools::install() + +# Get data --- + +library(akfishcondition) + +channel <- akfishcondition:::get_connected(schema = "AFSC") + +akfishcondition::get_condition_data(channel = channel) + + +# Generate data visualizations --- + +akfishcondition:::make_data_summary(dat_csv = here::here("data", "nbs_all_species.csv"), region = "NBS") +akfishcondition:::make_data_summary(dat_csv = here::here("data", "ebs_all_species.csv"), region = "EBS") +akfishcondition:::make_data_summary(dat_csv = here::here("data", "ai_all_species.csv"), region = "AI") +akfishcondition:::make_data_summary(dat_csv = here::here("data", "goa_all_species.csv"), region = "GOA") + + +# Calculate residual condition indicator ---- + +ai_sbw <- akfishcondition:::run_sbw_condition_multimodel( + region = "AI", + covariates_to_use = c('sex', 'stratum'), + min_n = 10) + +goa_sbw <- akfishcondition:::run_sbw_condition_multimodel( + region = "GOA", + covariates_to_use = c('sex', 'stratum'), + min_n = 10) + +ebs_sbw <- akfishcondition:::run_sbw_condition_multimodel( + region = "EBS", + covariates_to_use = c('sex', 'stratum'), + min_n = 10) + +nbs_sbw <- akfishcondition:::run_sbw_condition_multimodel( + region = "NBS", + covariates_to_use = c('sex', 'stratum'), + min_n = 10) + + +# Add updated indicator to built-in data sets ---- + +EBS_INDICATOR <- list( + FULL_REGION = as.data.frame( + ebs_sbw$full_sbw) |> + dplyr::filter(common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$EBS]), + STRATUM = as.data.frame( + dplyr::filter( + ebs_sbw$stratum_sbw, + common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$EBS]) + ), + LAST_UPDATE = Sys.Date() +) + +NBS_INDICATOR <- list( + FULL_REGION = as.data.frame( + nbs_sbw$full_sbw) |> + dplyr::filter( + common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$NBS]), + STRATUM = as.data.frame( + dplyr::filter( + nbs_sbw$stratum_sbw, + common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$EBS]) + ), + LAST_UPDATE = Sys.Date() +) + +GOA_INDICATOR <- list( + FULL_REGION = as.data.frame( + goa_sbw$full_sbw) |> + dplyr::filter(common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$GOA]), + STRATUM = as.data.frame( + dplyr::filter( + goa_sbw$stratum_sbw, + common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$GOA]) + ), + LAST_UPDATE = Sys.Date() +) + +AI_INDICATOR <- list( + FULL_REGION = as.data.frame( + ai_sbw$full_sbw) |> + dplyr::filter(common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$AI]), + STRATUM = as.data.frame( + dplyr::filter( + ai_sbw$stratum_sbw, + common_name %in% ESR_SETTINGS$ESR_SPECIES$common_name[ESR_SETTINGS$ESR_SPECIES$AI]) + ), + LAST_UPDATE = Sys.Date() +) + +PCOD_ESP <- list( + FULL_REGION_EBS = as.data.frame( + dplyr::filter(ebs_sbw$full_sbw, + common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name + )), + FULL_REGION_GOA = as.data.frame( + dplyr::filter(goa_sbw$full_sbw, + common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name + )), + FULL_REGION_AI = as.data.frame( + dplyr::filter(ai_sbw$full_sbw, + common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name + )), + FULL_REGION_NBS = as.data.frame( + dplyr::filter( + nbs_sbw$full_sbw, common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name + )), + STRATUM_EBS = as.data.frame( + dplyr::filter( + ebs_sbw$stratum_sbw, common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name + )), + STRATUM_GOA = as.data.frame( + dplyr::filter( + goa_sbw$stratum_sbw, common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name + )), + STRATUM_AI = as.data.frame( + dplyr::filter( + ai_sbw$stratum_sbw, common_name %in% ESP_SETTINGS$ESP_SPECIES$common_name + )), + stratum_NBS = NA, + LAST_UPDATE = Sys.Date() +) + +usethis::use_data(EBS_INDICATOR, overwrite = TRUE) +usethis::use_data(AI_INDICATOR, overwrite = TRUE) +usethis::use_data(GOA_INDICATOR, overwrite = TRUE) +usethis::use_data(NBS_INDICATOR, overwrite = TRUE) +usethis::use_data(PCOD_ESP, overwrite = TRUE) + +save( + EBS_INDICATOR, + NBS_INDICATOR, + GOA_INDICATOR, + AI_INDICATOR, + PCOD_ESP, + ESR_SETTINGS, + ESP_SETTINGS, + file = here::here("R", "sysdata.rda") +) + +AI_raw <- ai_sbw$input_data +GOA_raw <- goa_sbw$input_data +EBS_raw <- ebs_sbw$input_data +NBS_raw <- nbs_sbw$input_data + +save(AI_raw, + GOA_raw, + EBS_raw, + NBS_raw, + file = here::here("inst", "extdata", "raw_lw_bio.rda")) + +# Check update ---- + +library(akfishcondition) + +print(akfishcondition::AI_INDICATOR$LAST_UPDATE) +print(akfishcondition::GOA_INDICATOR$LAST_UPDATE) +print(akfishcondition::EBS_INDICATOR$LAST_UPDATE) +print(akfishcondition::NBS_INDICATOR$LAST_UPDATE) \ No newline at end of file diff --git a/AI_GroundfishCondition_2022.Rmd b/AI_GroundfishCondition_2022.Rmd deleted file mode 100644 index e6f2d25..0000000 --- a/AI_GroundfishCondition_2022.Rmd +++ /dev/null @@ -1,244 +0,0 @@ ---- -title: Aleutian Islands Groundfish Condition -author: -- affiliation: RACE - description: Research Fisheries Biologist - email: cecilia.oleary@NOAA.gov - name: Cecilia O'Leary -output: word_document -fontsize: 12pt -bibliography: EBS_references.bib -csl: fish-and-fisheries.csl -addr: - l1: 7600 Sand Pointe Way, NE - l2: NMFS RACE Division, Groundfish Assessment Program - l3: Seattle, WA 98115 ---- - -```{r setup, include=FALSE} -library(akfishcondition) -library(knitr) - -region <- "AI" -# Unzip map files -unzip(system.file("extdata/2022_ESR.zip", package = "akfishcondition")) - -if(!dir.exists(here::here("plots", region))) { - dir.create(here::here("plots", region), recursive = TRUE) -} -``` - -Contributed by Cecilia O'Leary^1^ and Sean Rohan^1^ -^1^ Resource Assessment and Conservation Engineering Division, Groundfish Assessment Program, Alaska Fisheries Science Center, National Marine Fisheries Service, NOAA, Seattle, WA -**Contact**: cecilia.oleary@noaa.gov -**Last updated**: October 2022 - -**Description of Indicator**: Morphometric condition indicators based on length-weight relationships characterize variation in somatic growth and can be considered indicators of prey availability, growth, general health, and habitat condition [@Blackwell2000; @Froese2006]. This contribution presents two morphometric condition indicators based on length-weight relationships: a new relative condition indicator that is estimated using a spatiotemporal model and the historical indicator based on residuals of the length-weight relationship. - -The new model-based relative condition indicator (VAST relative condition) is the ratio of fish weight-at-length relative to the time series mean based on annual allometric intercepts, _a~year~_, in the length-weight equation (_W_ = _aL_^_b_^; _W_ is mass (g), _L_ is fork length (cm)), i.e., $condition = a_{year}/\overline{a}$. Relative condition greater than one indicates better condition (i.e., heavier per unit length) and relative condition less than one indicates poorer condition (i.e., lighter per unit length). - -The historical length-weight indicator based on residuals of the length-weight relationship represents how heavy a fish is per unit body length relative to the time series mean [@Brodeur2004]. Positive length-weight residuals indicate better condition (i.e., heavier per unit length) and negative residuals indicate poorer condition (i.e., lighter per unit length) [@Froese2006]. Fish condition calculated in this way reflects realized outcomes of intrinsic and extrinsic processes that affect fish growth which can have implications for biological productivity through direct effects on growth and indirect effects on demographic processes such as, reproduction, and mortality (e.g., @Rodgveller2019; @Barbeaux2020). - -The model-based relative condition indicator was estimated using a spatiotemporal model with spatial random effects, implemented in the software VAST v3.8.2 [@Thorson2019; @Gruss2020a]. Allometric intercepts, _a~year~_, are estimated as fixed effects using a multivariate generalized linear mixed model that jointly estimates spatial and temporal variation in a and catch per unit effort (CPUE; numbers of fish per area swept). Density-weighted average _a~year~_ is a product of population density, local a, and area. Spatial variation in _a~year~_ was represented using a Gaussian Markov random field. The model approximates _a~year~_ using a log-link function and linear predictors [@Gruss2020a]. Parameters are estimated by identifying the values that maximize the marginal log-likelihood. - - -```{r map, include = TRUE, out.width = "200%", fig.cap = "\\label{fig:figs}Figure 1. National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center (AFSC) Resource Assessment and Conservation Engineering Groundfish Assessment Program (RACE-GAP) Aleutian Islands summer bottom trawl survey area with International North Pacific Fisheries Commission (INPFC) statistical fishing strata delineated by the purple lines.", echo = FALSE} -knitr::include_graphics("MapAI.png") -``` - -The historical indicator was estimated from residuals of linear regression models based on a log-transformation of the exponential growth relationship for all years where data were available from 1984 to 2022. A unique slope (_b_) was estimated for each stratum to account for spatial-temporal variation in growth and bottom trawl survey sampling. Strata were delineated based on International North Pacific Fisheries Commission (INPFC) stratum boundaries for the Southern Bering Sea, Eastern Aleutian Islands, Central Aleutian Islands, and Western Aleutian Islands (Figure 1). Length-weight relationships for 100–250 mm fork length (1–2 year old) walleye pollock were established independent of the adult life history stages caught. Bias-corrected weights-at-length (log scale) were estimated from the model and subtracted from observed weights to compute individual residuals per fish. Length-weight residuals were averaged for each INPFC stratum and weighted in proportion to stratum biomass based on stratified area-swept expansion of summer bottom trawl survey CPUE. Average length-weight residuals were compared by stratum and year to evaluate spatial variation in fish condition. Combinations of stratum and year with <10 samples were used for length-weight relationships but excluded from indicator calculations. - -Both condition indicators were calculated from paired fork length and weight measurements of individual fishes that were collected during bottom trawl surveys of the Aleutian Islands (AI) which were conducted by the Alaska Fisheries Science Center’s Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP). Analyses focused on walleye pollock (_Gadus chalcogrammus_), Pacific cod (_Gadus macrocephalus_), arrowtooth flounder (_Atheresthes stomias_), southern rock sole (_Lepidopsetta bilineata_), Atka mackerel (_Pleurogrammus monopterygius_), northern rockfish (_Sebastes polyspinis_), and Pacific ocean perch (_Sebastes alutus_) collected in trawls with satisfactory performance at standard survey stations (Figure 1). - -**Methodological Changes**: The historical length-weight residual indicator (used in 2020 and 2021) and new VAST relative condition indicator [@Gruss2020a] are both presented this year to allow comparison between methods. Overall, trends were similar between historical and new indicators based on the strong correlation (_r_ > 0.87) between indicators for most species (Figs. 3–4). An exception to this stronger correlation between the two indices was southern rock sole (_r_ = 0.57), where there was a large difference between historical and new indices in 1998, a year when all length-weight samples (n = 10) were collected at a single station and southern rock sole are almost exclusively caught in the southeastern Bering Sea in Aleutian years, a restricted spatial distribution likely impacting the historical condition indices and leading to a more precise current index as it accounts for spatially and temporally unbalanced sampling. Mean estimates and confidence intervals for the new condition indicator are likely more reliable than the historical indicator because the new indicator affords more precise expansion of individual samples to the population and better accounts for spatially and temporally unbalanced sampling that is characteristic of historical bottom trawl survey data due to changes in sampling protocols (e.g., transition from historical sex-and-length stratified to current random sampling). - -**Status and Trends**: Body condition varied amongst survey years for all species considered (Figure 2). Prior to 2010, the AI bottom trawl survey was characterized by condition cycling between positive and negative values through the years. Condition of most species since 2012 has primarily been below the long term average or neutral. Exceptions occur for 100–250 mm walleye pollock in 2016 and Atka mackerel in 2012 where the residual body condition is neutral or slightly positive. Overall, walleye pollock have fluctuated around the mean and are near the average in 2022. Walleye pollock > 250 mm had above average condition in 2010 and declined from 2010 - 2016, but were near average condition in 2018 and 2022. Atka mackerel showed above average condition in 2010, declined to below average by 2014, but have been near average since 2016. Southern rock sole residual body condition is trending positive in the Aleutians since 2012. Southern rock sole are near the time series mean in 2022. Pacific cod, northern rockfish, arrowtooth flounder, and Pacific ocean perch were above or near average condition in 2010, but subsequently had declining conditions through 2018. These species were in better condition in 2022 than 2018, but are still below their time series means. Notably, in 2022, residual body condition remained at neutral or became slightly more positive than the condition values since 2010 for all species considered, although the body condition of all species besides southern rock sole remain below the long term average for both the historical and model-based index. - - -```{r figure 2 grid, include = TRUE, echo = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 2. Weighted length-weight residuals for seven groundfish species collected during AFSC/RACE GAP standard summer bottom trawl surveys of the Aleutian Islands, 1986-2022. Filled bars denote weighted length-weight residuals using this year's indicator calculation. Error bars denote standard errors, thin black lines are 2 standard errors and thick blue boxes are 1 standard error.", message = FALSE, warning = FALSE} -fig2 <- plot_anomaly_timeseries(x = AI_INDICATOR$FULL_REGION, - region = "AI", - fill_color = "#0085CA", - var_y_name = "vast_relative_condition", - var_y_se_name = "vast_relative_condition_se", - var_x_name = "year", - y_title = "VAST relative condition", # expression("VAST Condition "~((grams/cm)^beta)), - plot_type = "box", - set_intercept = 1, - format_for = "rmd") - -fig2_png <- plot_anomaly_timeseries(x = AI_INDICATOR$FULL_REGION, - region = "AI", - fill_color = "#0085CA", - var_y_name = "vast_relative_condition", - var_y_se_name = "vast_relative_condition_se", - var_x_name = "year", - y_title = "VAST relative condition", # expression("VAST Condition "~((grams/cm)^beta)), - plot_type = "box", - set_intercept = 1, - format_for = "png") - -print(fig2 + theme_condition_index()) -``` - -```{r figure 2 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} -png(here::here("plots", region, "AI_condition.png"),width=6,height=7,units="in",res=600) -print(fig2_png + theme_blue_strip()) -dev.off() -``` - -```{r figure 3 grid, include = TRUE, echo = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 3. Length-weight residuals for groundfish species and age 1–2 walleye pollock (100–250 mm) collected during AFSC/RACE GAP summer bottom trawl surveys of the Aleutian Islands from 1986-2022. Black triangles denote the historical lenght-weight residual condition indicator and blue circles indicate the VAST relative condition indicator. Reported r values are the results of the Pearson's correlation.", message = FALSE, warning = FALSE} -fig3 <- plot_two_timeseries(x_1 = akfishcondition::AI_INDICATOR$FULL_REGION, - x_2 = akfishcondition::AI_INDICATOR$FULL_REGION, - region = "AI", - series_name_1 = "Length-weight residual (Historical)", - series_name_2 = "VAST relative condition (New)", - var_y_name_1 = "mean_wt_resid", - var_x_name_1 = "year", - var_y_name_2 = "vast_relative_condition", - var_x_name_2 = "year", - y_title = "Condition anomaly (Z-score)", - scale_y = TRUE, - fill_colors = c("#001743", "#0085CA"), - format_for = "rmd") - -fig3_png <- plot_two_timeseries(x_1 = akfishcondition::AI_INDICATOR$FULL_REGION, - x_2 = akfishcondition::AI_INDICATOR$FULL_REGION, - region = "AI", - series_name_1 = "Length-weight residual (Historical)", - series_name_2 = "VAST relative condition (New)", - var_y_name_1 = "mean_wt_resid", - var_x_name_1 = "year", - var_y_name_2 = "vast_relative_condition", - var_x_name_2 = "year", - y_title = "Condition anomaly (Z-score)", - scale_y = TRUE, - fill_colors = c("#001743", "#0085CA"), - format_for = "png") - -print(fig3 + theme_condition_index() + - theme(legend.position = "bottom", - legend.title = element_blank())) -``` - -```{r figure 3 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} -png(here::here("plots", region, "AI_VAST_SBW_timeseries.png"),width=6,height=7,units="in",res=600) -print(fig3_png + - theme_blue_strip() + - theme(legend.position = "bottom", - legend.title = element_blank())) -dev.off() -``` - -```{r figure 4 set up, include=FALSE, fig.height=4, fig.width=4, message=FALSE, warning=FALSE} -fig4 <- akfishcondition::plot_xy_corr( - x_1 = akfishcondition::AI_INDICATOR$FULL_REGION, - x_2 = akfishcondition::AI_INDICATOR$FULL_REGION, - region = region, - series_name_1 = "Length-weight residual (Historical)", - series_name_2 = "VAST relative condition (New)", - var_y_name_1 = "mean_wt_resid", - var_y_se_name_1 = "se_wt_resid", - var_x_name_1 = "year", - var_y_name_2 = "vast_relative_condition", - var_y_se_name_2 = "vast_relative_condition_se", - var_x_name_2 = "year", - format_for = "rmd") - -fig4_png <- akfishcondition::plot_xy_corr( - x_1 = akfishcondition::AI_INDICATOR$FULL_REGION, - x_2 = akfishcondition::AI_INDICATOR$FULL_REGION, - region = region, - series_name_1 = "Length-weight residual (Historical)", - series_name_2 = "VAST relative condition (New)", - var_y_name_1 = "mean_wt_resid", - var_y_se_name_1 = "se_wt_resid", - var_x_name_1 = "year", - var_y_name_2 = "vast_relative_condition", - var_y_se_name_2 = "vast_relative_condition_se", - var_x_name_2 = "year", - format_for = "png") -``` - -```{r figure 4 grid, include = TRUE, echo = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 4. Length-weight residual condition based on length-weight residuals versus VAST condition for the Aleutian Islands. Points denote the mean, error bars denote two standard errors.", message = FALSE, warning = FALSE} -print(fig4 + - theme_condition_index()) -``` - -```{r figure 4 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} -png(here::here("plots", region, "AI_VAST_VS_SBW_scatterplot.png"), width = 6, height = 7, units = "in", res = 600) -print(fig4_png + theme_blue_strip()) -dev.off() -``` - -```{r sbw_legacy_plots, include=FALSE, message=FALSE, warning=FALSE} -# Stratum-biomass weighted residual time series anomaly -sbw_timeseries <- plot_anomaly_timeseries(x = AI_INDICATOR$FULL_REGION, - region = region, - fill_color = "#0085CA", - var_y_name = "mean_wt_resid", - var_y_se_name = "se_wt_resid", - var_x_name = "year", - y_title = "Length-weight residual (ln(g))", - plot_type = "box", - format_for = "png", - set_intercept = 0) - -png(here::here("plots", region, "AI_SBW_timeseries.png"), width = 6, height = 7, units = "in", res = 600) -print(sbw_timeseries + akfishcondition::theme_blue_strip()) -dev.off() - -# Stratum-biomass weighted residual stratum stacked bar plots -stacked_bar <- akfishcondition::plot_stratum_stacked_bar(x = AI_INDICATOR$STRATUM, - region = region, - var_x_name = "year", - var_y_name = "stratum_resid_mean", - y_title = "Length-weight residual (ln(g))", - var_group_name = "inpfc_stratum", - fill_title = "Stratum", - fill_palette = "BrBG") - -png(here::here("plots", region, "AI_SBW_stratum_stacked_bar.png"), width = 6, height = 7, units = "in", res = 600) -print(stacked_bar + akfishcondition::theme_blue_strip()) -dev.off() - -# AI single species stratum plots -akfishcondition::plot_species_stratum_bar(x = AI_INDICATOR$STRATUM, - region = region, - var_x_name = "year", - var_y_name = "stratum_resid_mean", - var_y_se_name = "stratum_resid_se", - y_title = "Length-weight residual (ln(g))", - var_group_name = "inpfc_stratum", - fill_title = "Stratum", - fill_palette = "BrBG", - write_plot = TRUE) - -# EBS single species anomaly plots -akfishcondition::plot_species_bar(x = AI_INDICATOR$FULL_REGION, - region = "AI", - var_x_name = "year", - var_y_name = "vast_relative_condition", - var_y_se_name = "vast_relative_condition_se", - y_title = "VAST relative condition", - fill_color = "#0085CA", - set_intercept = 1, - write_plot = TRUE) -``` - -**Factors causing observed trends**: Several factors could affect morphological condition, including water temperature. Since the Warm Blob in 2014 [@Bond2015; @Stabeno2019a], there has been a general trend of warming ocean temperatures in the survey area through 2022 that could affect fish growth conditions there. The influence of temperature on growth rates depends on the physiology of predator species, prey availability, and the adaptive capacity of predators to respond to environmental change through migration, changes in behavior, and acclimatization. Thus, the factors underpinning the negative or neutral condition remain unclear. - -Other factors that could affect morphological condition include survey timing, stomach fullness, fish movement patterns, sex, and environmental conditions [@Froese2006]. Changing ocean conditions along with normal patterns of movement can cause the proportion of the population resident in the sampling area during the annual bottom trawl survey to vary. The date that the first length-weight data are collected is generally in the beginning of June and the bottom trawl survey is conducted throughout the summer months moving from east to west so that spatial and temporal trends in fish growth over the season become confounded with survey progress. We can expect some fish to exhibit seasonal or ontogenetic movement patterns during the survey months. Effects of survey timing on body condition can be further compounded by seasonal fluctuations in reproductive condition with the buildup and depletion of energy stores [@Wuenschel2019]. Another consideration is that fish weights sampled at sea include gut content weights so variation in gut fullness may influence weight measurements. Since feeding conditions vary over space and time, prey consumption rates and the proportion of total body weight attributable to gut contents may also be an important factor influencing length-at-weight. - -Finally, although condition indicators characterizes spatial and temporal variation in morphometric condition of groundfish species in the Aleutian Islands, it does not inform the mechanisms or processes behind the observed patterns. - -**Implications**: Fish morphometric condition can be considered an indicator of ecosystem productivity with implications for fish survival, maturity, and reproduction. For example, in Prince William Sound, the pre-winter condition of herring may determine their overwinter survival [@Paul1999], differences in feeding conditions have been linked to differences in morphometric condition of pink salmon in Prince William Sound [@Boldt2004], variation in morphometric condition has been linked to variation in maturity of sablefish [@Rodgveller2019], and lower morphometric condition of Pacific cod was associated with higher mortality and lower growth rates during the 2014–2016 marine heat wave in the Gulf of Alaska [@Barbeaux2020]. The condition of Aleutian Islands groundfish may similarly contribute to survival and recruitment and provide insight into ecosystem productivity, fish survival, demographic status, and population health. The condition of Aleutian Islands groundfish may similarly contribute to survival and recruitment and provide insight into ecosystem productivity, fish survival, demographic status, and population health. - -Survivorship is likely affected by many factors not examined here. As future years are added to the time series, the relationship between length-weight residuals and subsequent survival will be examined further. It is important to consider that residual body condition for most species in these analyses was computed for all sizes and sexes combined. Requirements for growth and survivorship differ for different fish life stages and some species have sexually dimorphic growth patterns. It may be more informative to examine life-stage (e.g., early juvenile, subadult, and adult phases) and sex specific body condition in the future for more insight into individual health and survivorship [@Froese2006]. - -The trend toward lowered body condition for many Aleutian Islands species from 2010 to 2018 RACE/AFSC GAP bottom trawl surveys (i.e., increasingly negative length-weight residuals) is a potential cause for concern. However, the increase in body condition for all groundfish species in 2022 may portend a reversal of the trend. Recent downward trends in body condition could indicate poor overwinter survival or may reflect the influence of locally changing environmental conditions depressing fish growth, local production, or survivorship. Indications are that the Warm Blob [@Bond2015; @Stabeno2019a] has been followed by years with elevated water temperatures - (e.g., @Barbeaux2020) which may be related to changes in fish condition in the species examined. As we continue to add years of fish condition to the record and expand on our knowledge of the relationships between condition, growth, production, and survival, we hope to gain more insight into the overall health of fish populations in the Aleutian Islands. - -**Research priorities**: The new model-based condition indicator (VAST relative condition) will be further explored for biases and sensitivities to data, model structure, and parameterization. Specifically, the 100-250 mm walleye pollock VAST relative condition indicator model does not converge, and so further model structure and parameterization research is needed. Research is also being planned and implemented across multiple AFSC programs to explore standardization of statistical methods for calculating condition indicators, and to examine relationships among putatively similar indicators of fish condition (i.e., morphometric, bioenergetic, physiological). Finally, we plan to explore variation in condition indices between life history stages alongside density dependence and climate change impacts [@Bolin2021; @Oke2022]. - -## References diff --git a/AI_GroundfishCondition_2024.Rmd b/AI_GroundfishCondition_2024.Rmd new file mode 100644 index 0000000..cd387de --- /dev/null +++ b/AI_GroundfishCondition_2024.Rmd @@ -0,0 +1,123 @@ +--- +title: Aleutian Islands Groundfish Condition +author: +- affiliation: RACE + description: Research Fisheries Biologist + email: cecilia.oleary@NOAA.gov + name: Cecilia O'Leary +output: word_document +fontsize: 12pt +bibliography: EBS_references.bib +csl: fish-and-fisheries.csl +addr: + l1: 7600 Sand Point Way, NE + l2: NMFS RACE Division, Groundfish Assessment Program + l3: Seattle, WA 98115 +--- + +```{r setup, include=FALSE} +library(akfishcondition) +library(knitr) + +region <- "AI" +# Unzip map files +unzip(system.file("extdata/2022_ESR.zip", package = "akfishcondition")) + +if(!dir.exists(here::here("plots", region))) { + dir.create(here::here("plots", region), recursive = TRUE) +} +``` + +Contributed by Cecilia O'Leary^1^ and Sean Rohan^1^ +^1^ Resource Assessment and Conservation Engineering Division, Groundfish Assessment Program, Alaska Fisheries Science Center, National Marine Fisheries Service, NOAA, Seattle, WA +**Contact**: cecilia.oleary@noaa.gov +**Last updated**: October 2022 + +**Description of Indicator**: Morphometric condition indicators based on length-weight relationships characterize variation in somatic growth and can be considered indicators of prey availability, growth, general health, and habitat condition [@Blackwell2000; @Froese2006]. This contribution presents two morphometric condition indicators based on length-weight relationships: a new relative condition indicator that is estimated using a spatiotemporal model and the historical indicator based on residuals of the length-weight relationship. + +To calculate indicators, length-weight relationships were estimated from linear regression models based on a log-transformation of the exponential growth relationship, W = aLb, where W is weight (g) and L is fork length (mm) for all areas for the period 1991–2024. Unique unique intercepts and (a) and slopes (b) were estimated for each survey stratum, sex, and interaction between stratum and sex to account for sexual dimorphism and spatial-temporal variation in growth and bottom trawl survey sampling. Length-weight relationships for 100–250 mm fork length walleye pollock (corresponding with ages 1–2 years) were calculated separately from adult walleye pollock (> 250 mm). Residuals for individual fish were obtained by subtracting observed weights from bias-corrected weights-at-length that were estimated from regression models. Individual length-weight residuals were aggregated and averaged for each stratum, weighted based on the proportion to total biomass in each stratum from area-swept expansion of bottom-trawl survey catch per unit effort (CPUE; i.e., design-based stratum biomass estimates). Variation in fish condition was evaluated by comparing average length-weight residuals among years. To minimize the influence of unrepresentative samples on indicator calculations, combinations of species, stratum, and year with a sample size < 10 were used to fit length-weight regressions but were excluded from calculating length-weight residuals. Morphometric condition indicator time series, code for calculating the indicators, and figures showing results for individual species are available through the akfishcondition R package and GitHub repository (https://github.com/afsc-gap-products/akfishcondition). + + +```{r map, include = TRUE, echo = FALSE, out.width = "100%", fig.cap = "\\label{fig:figs}Figure 1. National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center (AFSC) Resource Assessment and Conservation Engineering Groundfish Assessment Program (RACE-GAP) Aleutian Islands summer bottom trawl survey area with International North Pacific Fisheries Commission (INPFC) statistical fishing strata delineated by the purple lines."} +knitr::include_graphics("MapAI.png") +``` + +**Methodological Changes**: In 2022, historical stratum-biomass weighted residuals condition indicators were presented alongside condition indicators that were calculated using the R package VAST following methods that were presented for select GOA species during the Spring Preview of Ecological and Economic Conditions in May 2020. The authors noted there were strong correlations between VAST and stratum biomass weighted condition indicators for most species (r = 0.79–0.98). The authors received the following feedback about the change from the BSAI Groundfish Plan Team meeting during their November 2022 meeting: +"The Team discussed the revised condition indices that now use a different, VAST-based condition index, but felt additional methodology regarding this transition was needed. The Team recommended a short presentation next September to the Team to review the methods and tradeoffs in approaches. The Team encouraged collaboration with the NMFS longline survey team to develop analogous VAST indices." +Based on feedback from the Plan Team, staff limitations, and the lack of a clear path to transition condition indicators for longline species to VAST, the 2024 condition indicator was calculated from stratum-biomass weighted residuals of length-weight regressions. + +**Status and Trends**: Body condition varied amongst survey years for all species considered (Figure 2). Prior to 2010, the AI bottom trawl survey was characterized by condition cycling between positive and negative values through the years. Condition of most species since 2012 has primarily been below the long term average or neutral. Exceptions occur for 100–250 mm walleye pollock in 2016 and Atka mackerel in 2012 where the residual body condition is neutral or slightly positive. Overall, walleye pollock have fluctuated around the mean and are near the average in 2022. Walleye pollock > 250 mm had above average condition in 2010 and declined from 2010 - 2016, but were near average condition in 2018 and 2022. Atka mackerel showed above average condition in 2010, declined to below average by 2014, but have been near average since 2016. Southern rock sole residual body condition is trending positive in the Aleutians since 2012. Southern rock sole are near the time series mean in 2022. Pacific cod, northern rockfish, arrowtooth flounder, and Pacific ocean perch were above or near average condition in 2010, but subsequently had declining conditions through 2018. These species were in better condition in 2022 than 2018, but are still below their time series means. Notably, in 2022, residual body condition remained at neutral or became slightly more positive than the condition values since 2010 for all species considered, although the body condition of all species besides southern rock sole remain below the long term average for both the historical and model-based index. + + +```{r figure 2 grid, include = TRUE, echo = FALSE, message = FALSE, warning = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 2. Weighted length-weight residuals for seven groundfish species collected during AFSC/RACE GAP standard summer bottom trawl surveys of the Aleutian Islands, 1991--2024. Filled bars denote weighted length-weight residuals using this year's indicator calculation. Error bars denote standard errors, thin black lines are 2 standard errors and thick blue boxes are 1 standard error."} +fig2 <- plot_anomaly_timeseries(x = AI_INDICATOR$FULL_REGION, + region = "AI", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "rmd") + +fig2_png <- plot_anomaly_timeseries(x = AI_INDICATOR$FULL_REGION, + region = "AI", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "png") + +print(fig2 + theme_condition_index()) +``` + + +```{r figure 2 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", region, "AI_condition.png"),width=6,height=7,units="in",res=600) +print(fig2_png + theme_blue_strip()) +dev.off() +``` + + +```{r figure 3 grid, include = TRUE, echo = FALSE, message = FALSE, warning = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 3. Residual body condition index for Aleutian Islands groundfish species collected on the National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP) standard summer bottom trawl survey (1991--2024) grouped by International North Pacific Fisheries Commission (INPFC) statistical sampling strata."} +fig3 <- plot_stratum_stacked_bar(x = AI_INDICATOR$STRATUM, + region = "AI", + fill_palette = "BrBG", + var_y_name = "stratum_resid_mean", + var_x_name = "year", + var_group_name = "inpfc_stratum", + fill_title = "Stratum", + y_title = "Length-weight residual (ln(g))") + +print(fig3 + + theme_condition_index() + + theme(legend.position = "bottom", + legend.text = element_text(size = 19), + legend.title = element_blank())) + +``` + +```{r figure 3 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", region, "AI_condition_by_stratum.png") ,width=6,height=7,units="in",res=300) +print(fig3 + theme_blue_strip() + theme(legend.text = element_text(size = 8))) +dev.off() +``` + +**Factors causing observed trends**: Several factors could affect morphological condition, including water temperature. Since the Warm Blob in 2014 [@Bond2015; @Stabeno2019a], there has been a general trend of warming ocean temperatures in the survey area through 2022 that could affect fish growth conditions there. The influence of temperature on growth rates depends on the physiology of predator species, prey availability, and the adaptive capacity of predators to respond to environmental change through migration, changes in behavior, and acclimatization. Thus, the factors underpinning the negative or neutral condition remain unclear. + +Other factors that could affect morphological condition include survey timing, stomach fullness, fish movement patterns, sex, and environmental conditions [@Froese2006]. Changing ocean conditions along with normal patterns of movement can cause the proportion of the population resident in the sampling area during the annual bottom trawl survey to vary. The date that the first length-weight data are collected is generally in the beginning of June and the bottom trawl survey is conducted throughout the summer months moving from east to west so that spatial and temporal trends in fish growth over the season become confounded with survey progress. We can expect some fish to exhibit seasonal or ontogenetic movement patterns during the survey months. Effects of survey timing on body condition can be further compounded by seasonal fluctuations in reproductive condition with the buildup and depletion of energy stores [@Wuenschel2019]. Another consideration is that fish weights sampled at sea include gut content weights so variation in gut fullness may influence weight measurements. Since feeding conditions vary over space and time, prey consumption rates and the proportion of total body weight attributable to gut contents may also be an important factor influencing length-at-weight. + +Finally, although condition indicators characterizes spatial and temporal variation in morphometric condition of groundfish species in the Aleutian Islands, it does not inform the mechanisms or processes behind the observed patterns. + +**Implications**: Fish morphometric condition can be considered an indicator of ecosystem productivity with implications for fish survival, maturity, and reproduction. For example, in Prince William Sound, the pre-winter condition of herring may determine their overwinter survival [@Paul1999], differences in feeding conditions have been linked to differences in morphometric condition of pink salmon in Prince William Sound [@Boldt2004], variation in morphometric condition has been linked to variation in maturity of sablefish [@Rodgveller2019], and lower morphometric condition of Pacific cod was associated with higher mortality and lower growth rates during the 2014–2016 marine heat wave in the Gulf of Alaska [@Barbeaux2020]. The condition of Aleutian Islands groundfish may similarly contribute to survival and recruitment and provide insight into ecosystem productivity, fish survival, demographic status, and population health. The condition of Aleutian Islands groundfish may similarly contribute to survival and recruitment and provide insight into ecosystem productivity, fish survival, demographic status, and population health. + +Survivorship is likely affected by many factors not examined here. As future years are added to the time series, the relationship between length-weight residuals and subsequent survival will be examined further. It is important to consider that residual body condition for most species in these analyses was computed for all sizes and sexes combined. Requirements for growth and survivorship differ for different fish life stages and some species have sexually dimorphic growth patterns. It may be more informative to examine life-stage (e.g., early juvenile, subadult, and adult phases) and sex specific body condition in the future for more insight into individual health and survivorship [@Froese2006]. + +The trend toward lowered body condition for many Aleutian Islands species from 2010 to 2018 RACE/AFSC GAP bottom trawl surveys (i.e., increasingly negative length-weight residuals) is a potential cause for concern. However, the increase in body condition for all groundfish species in 2022 may portend a reversal of the trend. Recent downward trends in body condition could indicate poor overwinter survival or may reflect the influence of locally changing environmental conditions depressing fish growth, local production, or survivorship. Indications are that the Warm Blob [@Bond2015; @Stabeno2019a] has been followed by years with elevated water temperatures + (e.g., @Barbeaux2020) which may be related to changes in fish condition in the species examined. As we continue to add years of fish condition to the record and expand on our knowledge of the relationships between condition, growth, production, and survival, we hope to gain more insight into the overall health of fish populations in the Aleutian Islands. + +**Research priorities**: The new model-based condition indicator (VAST relative condition) will be further explored for biases and sensitivities to data, model structure, and parameterization. Specifically, the 100-250 mm walleye pollock VAST relative condition indicator model does not converge, and so further model structure and parameterization research is needed. Research is also being planned and implemented across multiple AFSC programs to explore standardization of statistical methods for calculating condition indicators, and to examine relationships among putatively similar indicators of fish condition (i.e., morphometric, bioenergetic, physiological). Finally, we plan to explore variation in condition indices between life history stages alongside density dependence and climate change impacts [@Bolin2021; @Oke2022]. + +## References diff --git a/AI_GroundfishCondition_2024.docx b/AI_GroundfishCondition_2024.docx new file mode 100644 index 0000000..802760f Binary files /dev/null and b/AI_GroundfishCondition_2024.docx differ diff --git a/DESCRIPTION b/DESCRIPTION index ce5c472..a58931f 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Package: akfishcondition Type: Package Title: Groundfish morphometric condition indicator -Version: 4.0.2 +Version: 4.1.1 Authors@R: c(person("Sean", "Rohan", email = "sean.rohan@noaa.gov", role = c("aut", "cre")), person("Cecilia", "O'Leary", email = "cecilia.oleary@noaa.gov", role = c("aut")), person("Bianca", "Prohaska", email = "bianca.prohaska@noaa.gov", role = c("ctb")), @@ -22,15 +22,15 @@ Description: This package calculates groundfish condition based on length-weight License: GPL (>=2) Encoding: UTF-8 LazyData: false -RoxygenNote: 7.2.3 +RoxygenNote: 7.3.1 Imports: colorspace, viridis Suggests: knitr, rmarkdown, RODBC, getPass, - VAST (>= 3.8.2), - FishStatsUtils (>= 2.11.0), ggthemes, akgfmaps +Remotes: James-Thorson-NOAA/VAST (>= 3.8.2), + James-Thorson-NOAA/FishStatsUtils (>= 2.11.0) VignetteBuilder: knitr diff --git a/EBS_GroundfishCondition_2023.Rmd b/EBS_GroundfishCondition_2023.Rmd index 79b3ea4..e55b4a2 100644 --- a/EBS_GroundfishCondition_2023.Rmd +++ b/EBS_GroundfishCondition_2023.Rmd @@ -239,7 +239,7 @@ fig5_png <- akfishcondition::plot_stratum_stacked_bar(x = NBS_INDICATOR$STRATUM, fill_title = "Stratum", fill_palette = "PuBu") -print(fig5 + theme_condition_index()+ theme(strip.text = element_text(size = 20))) +print(fig5 + theme_condition_index() + theme(strip.text = element_text(size = 20))) ``` @@ -249,96 +249,6 @@ print(fig5_png + theme_blue_strip() + theme(legend.position = "right", legend.ti dev.off() ``` - - -```{r sbw_legacy_plots, include=FALSE, message=FALSE, warning=FALSE} -# # EBS stratum-biomass weighted residual time series anomaly -# ebs_sbw_timeseries <- plot_anomaly_timeseries(x = EBS_INDICATOR$FULL_REGION, -# region = region, -# fill_color = "#0085CA", -# var_y_name = "mean_wt_resid", -# var_y_se_name = "se_wt_resid", -# var_x_name = "year", -# y_title = "Length-weight residual (ln(g))", -# plot_type = "box", -# format_for = "png", -# set_intercept = 0) -# -# png(here::here("plots", region, "EBS_SBW_timeseries.png"), width = 6, height = 7, units = "in", res = 600) -# print(ebs_sbw_timeseries + akfishcondition::theme_blue_strip()) -# dev.off() -# -# # EBS stratum-biomass weighted residual stratum stacked bar plots -# ebs_stacked_bar <- akfishcondition::plot_stratum_stacked_bar(x = EBS_INDICATOR$STRATUM, -# region = region, -# var_x_name = "year", -# var_y_name = "stratum_resid_mean", -# y_title = "Length-weight residual (ln(g))", -# var_group_name = "stratum", -# fill_title = "Stratum", -# fill_palette = "BrBG") -# -# png(here::here("plots", region, "EBS_SBW_stratum_stacked_bar.png"), width = 6, height = 7, units = "in", res = 600) -# print(ebs_stacked_bar + theme_blue_strip() + theme(legend.position = "right", legend.title = element_text(size = 9))) -# dev.off() -# -# # EBS single species stratum plots -# akfishcondition::plot_species_stratum_bar(x = EBS_INDICATOR$STRATUM, -# region = "EBS", -# var_x_name = "year", -# var_y_name = "stratum_resid_mean", -# var_y_se_name = "stratum_resid_se", -# y_title = "Length-weight residual (ln(g))", -# var_group_name = "stratum", -# fill_title = "Stratum", -# fill_palette = "BrBG", -# write_plot = TRUE) -# -# # NBS stratum-biomass weighted residual time series anomaly -# nbs_sbw_timeseries <- plot_anomaly_timeseries(x = NBS_INDICATOR$FULL_REGION, -# region = region, -# fill_color = "#0085CA", -# var_y_name = "mean_wt_resid", -# var_y_se_name = "se_wt_resid", -# var_x_name = "year", -# y_title = "Length-weight residual (ln(g))", -# plot_type = "box", -# format_for = "png", -# set_intercept = 0) -# -# png(here::here("plots", region, "NBS_SBW_timeseries.png"), width = 6, height = 7, units = "in", res = 600) -# print(nbs_sbw_timeseries + akfishcondition::theme_blue_strip()) -# dev.off() -# -# # NBS stratum-biomass weighted residual stratum stacked bar plots -# nbs_stacked_bar <- akfishcondition::plot_stratum_stacked_bar(x = NBS_INDICATOR$STRATUM, -# region = region, -# var_x_name = "year", -# var_y_name = "stratum_resid_mean", -# y_title = "Length-weight residual (ln(g))", -# var_group_name = "stratum", -# fill_title = "Stratum", -# fill_palette = "PuBu") -# -# png(here::here("plots", region, "NBS_SBW_stratum_stacked_bar.png"), width = 6, height = 7, units = "in", res = 600) -# print(nbs_stacked_bar + theme_blue_strip() + theme(legend.position = "right", legend.title = element_text(size = 9))) -# dev.off() -# -# -# # NBS single species stratum plots -# akfishcondition::plot_species_stratum_bar(x = NBS_INDICATOR$STRATUM, -# region = "NBS", -# var_x_name = "year", -# var_y_name = "stratum_resid_mean", -# var_y_se_name = "stratum_resid_se", -# y_title = "Length-weight residual (ln(g))", -# var_group_name = "stratum", -# fill_title = "Stratum", -# fill_palette = "PuBu", -# write_plot = TRUE) - -``` - **Factors influencing observed trends**: Temperature appears to influence morphological condition of several species in the EBS and NBS, so near-average cold pool extent and water temperatures in 2023 likely played a role in the near-average condition (within 1 S.D. of the mean) for most species in the EBS and NBS. Historically, particularly cold years tend to correspond with negative condition, while particularly warm years tend to correspond to positive condition. For example, water temperatures were particularly cold during the 1999 Bering Sea survey, a year in which negative condition was observed for all species that data were available. In addition, spatiotemporal factor analyses suggest the morphometric condition of age-7 walleye pollock is strongly correlated with cold pool extent in the EBS [@Gruss2021]. In recent years, warm temperatures across the Bering Sea shelf, since the record low seasonal sea ice extent in 2017–2018 and historical cold pool area minimum in 2018 [@Stabeno2019a], may have influenced the positive trend in the condition of several species from 2016 to 2019. However, despite near-average temperature in 2023 large walleye pollock (>250 mm) condition in the EBS was the second lowest recorded over the time series. Although warmer temperatures may increase growth rates if there is adequate prey to offset temperature-dependent increases in metabolic demand, growth rates may also decline if prey resources are insufficient to offset temperature-dependent increases in metabolic demand. The influence of temperature on growth rates depends on the physiology of predator species, prey availability, and the adaptive capacity of predators to respond to environmental change through migration, changes in behavior, and acclimatization. For example, elevated temperatures during the 2014–2016 marine heatwave in the Gulf of Alaska led to lower growth rates of Pacific cod and lower condition because available prey resources did not make up for increased metabolic demand [@Barbeaux2020]. diff --git a/EBS_GroundfishCondition_2023.docx b/EBS_GroundfishCondition_2023.docx index 61e94c6..f212b74 100644 Binary files a/EBS_GroundfishCondition_2023.docx and b/EBS_GroundfishCondition_2023.docx differ diff --git a/EBS_references.bib b/EBS_references.bib index 930820d..add7daf 100644 --- a/EBS_references.bib +++ b/EBS_references.bib @@ -125,22 +125,6 @@ @article{Gruss2021 year = {2021} } -@article{Haberle2023, -abstract = {Individual performance defines population dynamics. Condition index – a ratio of weight and some function of length – has been louded as an indicator of individual performance and recommended as a tool in fisheries management and conservation. However, insufficient understanding of the correlation between individual‐level processes and population‐level responses hinders its adoption. To this end, we use composite modelling to link individual's condition, expressed through the condition index, to population‐level status. We start by modelling ontogeny of European pilchard ( Sardina pilchardus , Clupeidae) as a function of food and constant temperature using Dynamic Energy Budget theory. We then provide a framework to simultaneously track the individual‐ and population‐level statistics by incorporating the dynamic energy budget model into an individual‐based model. Lastly, we explore the effects of fishing pressure on the statistics in two constant and food‐limited environmental carrying capacity scenarios. Results show that, regardless of the species' environmental carrying capacity, individual condition index will increase with fishing mortality, that is, with reduction of stock size. Same patterns are observed for gilthead seabream ( Sparus aurata , Sparidae), a significantly different species. Condition index can, therefore, in food‐limited populations, be used to (i) estimate population size relative to carrying capacity and (ii) distinguish overfished from underfished populations. Our findings promote a practical way to operationally incorporate the condition index into fisheries management and marine conservation, thus providing additional use for the commonly collected biometric data. Some real‐world applications, however, may require additional research to account for other variables such as fluctuating environmental conditions and individual variability.}, -author = {Haberle, Ines and Bav{\v{c}}evi{\'{c}}, Lav and Klanjscek, Tin}, -doi = {10.1111/faf.12744}, -issn = {1467-2960}, -journal = {Fish and Fisheries}, -keywords = {dynamic energy budget,fisheries management,fishing mortality,individual-}, -month = {jul}, -number = {4}, -pages = {567--581}, -title = {{Fish condition as an indicator of stock status: Insights from condition index in a food‐limiting environment}}, -url = {https://onlinelibrary.wiley.com/doi/10.1111/faf.12744}, -volume = {24}, -year = {2023} -} - @techreport{Hurst2021, author = {Hurst, Thomas P. and O'Leary, Cecilia A. and Rohan, Sean K. and Siddon, Elizabeth C. and Thorson, James T. and Vollenweider, Johanna J.}, doi = {10.25923/p1yd-0793}, diff --git a/GOA_GroundfishCondition_2023.Rmd b/GOA_GroundfishCondition_2023.Rmd index 8a7b6e4..2ae1d1b 100644 --- a/GOA_GroundfishCondition_2023.Rmd +++ b/GOA_GroundfishCondition_2023.Rmd @@ -25,10 +25,6 @@ unzip(system.file("extdata/2022_ESR.zip", package = "akfishcondition")) dir.create(here::here("plots", region)) - -# Load data -# goa_ann_mean_resid_df <- GOA_INDICATOR$FULL_REGION -# goa_stratum_resids_df <- GOA_INDICATOR$STRATUM ``` Contributed by Cecilia O'Leary^1^ and Sean Rohan^1^ @@ -39,13 +35,13 @@ Contributed by Cecilia O'Leary^1^ and Sean Rohan^1^ **Description of Indicator**: Length-weight residuals represent how heavy a fish is per unit body length and are an indicator of somatic growth variability (Brodeur et al., 2004). Therefore, length-weight residuals can be considered indicators of prey availability, growth, general health, and habitat condition (Blackwell et al., 2000; Froese, 2006). Positive length-weight residuals indicate better condition (i.e., heavier per unit length) and negative residuals indicate poorer condition (i.e., lighter per unit length) (Froese, 2006). Fish condition calculated in this way reflects realized outcomes of intrinsic and extrinsic processes that affect fish growth which can have implications for biological productivity through direct effects on growth and indirect effects on demographic processes such as, reproduction, and mortality (e.g., Rodgveller (2019); Barbeaux et al. (2020)). -```{r map, include = TRUE, out.width = "200%", fig.cap = "\\label{fig:figs}Figure 1. National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP) Gulf of Alaska summer bottom trawl survey area with International North Pacific Fisheries Commission (INPFC) statistical fishing strata delineated by the red lines.", echo = FALSE} +```{r map, include = TRUE, out.width = "100%", fig.cap = "\\label{fig:figs}Figure 1. National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP) Gulf of Alaska summer bottom trawl survey area with International North Pacific Fisheries Commission (INPFC) statistical fishing strata delineated by the red lines.", echo = FALSE} knitr::include_graphics("MapGOA.png") ``` The groundfish morphometric condition indicator is calculated from paired fork lengths (mm) and weights (g) of individual fishes that were collected during bottom trawl survey of the Gulf of Alaska (GOA) which were conducted by the Alaska Fisheries Science Center biennial Resource Assessment and Conservation Engineering (AFSC/RACE) - Groundfish Assessment Program’s (GAP). Fish condition analyses were applied to walleye pollock (*Gadus chalcogrammus*), Pacific cod (*Gadus macrocephalus*), arrowtooth flounder (*Atheresthes stomias*), southern rock sole (*Lepidopsetta bilineata*), northern rockfish (*Sebastes polyspinis*), Pacific ocean perch (*Sebastes alutus*), dusky rockfish (*Sebastes variabilis*), shortraker rockfish (*Sebastes borealis*), rougheye rockfish (*Sebastes aleutianus*),sharpchin rockfish (*Sebastes zacentrus*),flathead sole (*Hippoglossoides elassodon*),dover sole (*Microstomus pacificus*), and rex sole(*Glyptocephalus zachirus*) collected in trawls with satisfactory performance at standard survey stations. Data were combined in the former International North Pacific Fisheries Commission (INPFC) strata; Shumagin, Chirikof, Kodiak, Yakutat and Southeast (Figure 1). -To calculate indicators, length-weight relationships were estimated from linear regression models based on a log-transformation of the exponential growth relationship, W = aLb, where W is weight (g) and L is fork length (mm) for all areas for the period 1984–2023. UniqueA unique intercepts and (a) and slopes (b) wereas estimated for each survey stratum, sex, and interaction between stratum and sex to account for sexual dimorphism and spatial-temporal variation in growth and bottom trawl survey sampling. Length-weight relationships for 100–250 mm fork length walleye pollock (corresponding with ages 1–2 years) were calculated separately from adult walleye pollock (> 250 mm). Residuals for individual fish were obtained by subtracting observed weights from bias-corrected weights-at-length that were estimated from regression models. Individual length-weight residuals were aggregated and averaged for each stratum, weighted based on the proportion to total biomass in each stratum from area-swept expansion of bottom-trawl survey catch per unit effort (CPUE; i.e., design-based stratum biomass estimates). Variation in fish condition was evaluated by comparing average length-weight residuals among years. To minimize the influence of unrepresentative samples on indicator calculations, combinations of species, stratum, and year with a sample size < 10 were used to fit length-weight regressions but were excluded from calculating length-weight residuals. Morphometric condition indicator time series, code for calculating the indicators, and figures showing results for individual species are available through the akfishcondition R package and GitHub repository (https://github.com/afsc-gap-products/akfishcondition). +To calculate indicators, length-weight relationships were estimated from linear regression models based on a log-transformation of the exponential growth relationship, W = aLb, where W is weight (g) and L is fork length (mm) for all areas for the period 1984–2023. Unique unique intercepts and (a) and slopes (b) were estimated for each survey stratum, sex, and interaction between stratum and sex to account for sexual dimorphism and spatial-temporal variation in growth and bottom trawl survey sampling. Length-weight relationships for 100–250 mm fork length walleye pollock (corresponding with ages 1–2 years) were calculated separately from adult walleye pollock (> 250 mm). Residuals for individual fish were obtained by subtracting observed weights from bias-corrected weights-at-length that were estimated from regression models. Individual length-weight residuals were aggregated and averaged for each stratum, weighted based on the proportion to total biomass in each stratum from area-swept expansion of bottom-trawl survey catch per unit effort (CPUE; i.e., design-based stratum biomass estimates). Variation in fish condition was evaluated by comparing average length-weight residuals among years. To minimize the influence of unrepresentative samples on indicator calculations, combinations of species, stratum, and year with a sample size < 10 were used to fit length-weight regressions but were excluded from calculating length-weight residuals. Morphometric condition indicator time series, code for calculating the indicators, and figures showing results for individual species are available through the akfishcondition R package and GitHub repository (https://github.com/afsc-gap-products/akfishcondition). **Methodological changes**: In 2022, historical stratum-biomass weighted residuals condition indicators were presented alongside condition indicators that were calculated using the R package VAST following methods that were presented for select GOA species during the Spring Preview of Ecological and Economic Conditions in May 2020. The authors noted there were strong correlations between VAST and stratum biomass weighted condition indicators for most species (r = 0.79–0.98). The authors received the following feedback about the change from the BSAI Groundfish Plan Team meeting during their November 2022 meeting: "The Team discussed the revised condition indices that now use a different, VAST-based condition index, but felt additional methodology regarding this transition was needed. The Team recommended a short presentation next September to the Team to review the methods and tradeoffs in approaches. The Team encouraged collaboration with the NMFS longline survey team to develop analogous VAST indices." @@ -147,12 +143,18 @@ fig3b <- plot_stratum_stacked_bar(x = dplyr::filter(GOA_INDICATOR$STRATUM, fill_title = "Stratum", y_title = "Length-weight residual (ln(g))") -print(fig3a + theme_condition_index()) +print(fig3a + theme_condition_index() + + theme_condition_index() + + theme(legend.position = "bottom", + legend.title = element_blank())) ``` ```{r figure 3b grid, include = TRUE, echo = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 3B. Residual body condition index for sixteen Gulf of Alaska groundfish species collected on the National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP) standard summer bottom trawl survey (1984--2023) grouped by International North Pacific Fisheries Commission (INPFC) statistical sampling strata.", message = FALSE, warning = FALSE} -print(fig3b + theme_condition_index()) +print(fig3b + theme_condition_index() + + theme_condition_index() + + theme(legend.position = "bottom", + legend.title = element_blank())) ``` ```{r figure 3 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} diff --git a/GOA_GroundfishCondition_2023.docx b/GOA_GroundfishCondition_2023.docx new file mode 100644 index 0000000..c478595 Binary files /dev/null and b/GOA_GroundfishCondition_2023.docx differ diff --git a/R/data.R b/R/data.R index 7b00566..d00b62f 100644 --- a/R/data.R +++ b/R/data.R @@ -1,6 +1,6 @@ -#' Aleutian Islands (1986-2022) +#' Aleutian Islands (1991-2024) #' -#' Morphometric condition indicators based on the allometric intercept estimated using VAST and residuals from a length-weight regression for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), southern rock sole, arrowtooth flounder, Atka mackerel, northern rockfish, and Pacific ocean perch in the Aleutian Islands. +#' Morphometric condition indicators based on residuals from length-weight regressions for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), southern rock sole, arrowtooth flounder, Atka mackerel, northern rockfish, and Pacific ocean perch in the Aleutian Islands. #' #' @format A list containing two data frames (indicator for the full region, indicator by stratum) and a character vector. #' \describe{ @@ -11,24 +11,17 @@ #' \item{common_name}{: Year} #' \item{mean_wt_resid}{: Mean residual for the full region} #' \item{se_wt_resid}{: Standard error of the indicator for the full region} -#' \item{vast_condition}{: Estimate of the allometric intercept,a, of the length-weight equation (W = aL^b), estimated using VAST.} -#' \item{vast_condition_se}{Standard error of a, estimated using VAST.} -#' \item{species_code}{: RACE/GAP species code} -#' \item{scaled_vast_condition}{: Anomaly (Z-score) of vast_condition} -#' \item{vast_relative_condition}{: Relative condition based on the allometric intercept estimated using VAST, i.e., allometric intercept divided by the time series mean allometric intercept.} -#' \item{vast_relative_condition_se}{: Standard error of vast_relative_condition.} #' } #' \item STRATUM (data frame): Residuals by stratum #' \itemize{ -#' \item{year}{: Year} #' \item{common_name}{: Species Common name} -#' \item{species_code}{: RACE/GAP species code} +#' \item{year}{: Year} #' \item{inpfc_stratum}{: INPFC Stratum} #' \item{stratum_resid_mean}{: Unweighted stratum mean residual} #' \item{n}{: Sample size for the stratum} #' \item{stratum_resid_se}{: Standard error of the stratum mean residual} #' \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} -#' \item{weighted_resid_mean}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} +#' \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} #' } #' \item LAST_UPDATE: Date when indicator and data were last updated. #' } @@ -37,31 +30,30 @@ #' @export "AI_INDICATOR" -#' Gulf of Alaska (1984-2021) +#' Gulf of Alaska (1990-2023) #' -#' Morphometric condition indicators based on residuals from a length-weight regression for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), southern rock sole, arrowtooth flounder, Atka mackerel, northern rockfish, and Pacific ocean perch in the Gulf of Alaska. +#' Morphometric condition indicators based on residuals from length-weight regressions for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), Pacific cod, northern rock sole, southern rock sole, arrowtooth flounder, Dover sole, rex sole, flathead sole, Atka mackerel, blackspotted rockfish, rougheye rockfish, northern rockfish, shortraker rockfish, shaprchin rockfish, and Pacific ocean perch in the Gulf of Alaska. #' #' @format A list containing two data frames (indicator for the full region, indicator by stratum) and a character vector. #' \describe{ #' \itemize{ #' \item FULL_REGION (data frame): Residuals for the full region #' \itemize{ -#' \item{YEAR}{: Year} -#' \item{COMMON_NAME}{: Year} +#' \item{year}{: Year} +#' \item{common_name}{: Year} #' \item{mean_wt_resid}{: Mean residual for the full region} #' \item{se_wt_resid}{: Standard error of the indicator for the full region} #' } #' \item STRATUM (data frame): Residuals by stratum #' \itemize{ -#' \item{YEAR}{: Year} -#' \item{COMMON_NAME}{: Species Common name} -#' \item{SPECIES_CODE}{: RACE/GAP species code} -#' \item{INPFC_STRATUM}{: INPFC Stratum} +#' \item{common_name}{: Species Common name} +#' \item{year}{: Year} +#' \item{inpfc_stratum}{: INPFC Stratum} #' \item{stratum_resid_mean}{: Unweighted stratum mean residual} #' \item{n}{: Sample size for the stratum} #' \item{stratum_resid_sd}{: Standard error of the stratum mean residual} #' \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} -#' \item{weighted_resid_mean}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} +#' \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} #' } #' \item LAST_UPDATE: Date when indicator and data were last updated. #' } @@ -72,7 +64,7 @@ #' Eastern Bering Sea continental shelf (1999-2022) #' -#' Morphometric condition indicators based on the allometric intercept estimated using VAST and residuals from a length-weight regression for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), northern rock sole, arrowtooth flounder, flathead sole, Alaska plaice, and yellowfin sole in the eastern Bering Sea. +#' Morphometric condition indicators based on length-weight regressions for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), northern rock sole, arrowtooth flounder, flathead sole, Alaska plaice, and yellowfin sole in the eastern Bering Sea. #' #' @format A list containing two data frames (indicator for the full region, indicator by stratum) and a character vector. #' \describe{ @@ -83,24 +75,18 @@ #' \item{common_name}{: Year} #' \item{mean_wt_resid}{: Mean residual for the full region} #' \item{se_wt_resid}{: Standard error of the indicator for the full region} -#' \item{vast_condition}{: Estimate of the allometric intercept,a, of the length-weight equation (W = aL^b), estimated using VAST.} -#' \item{vast_condition_se}{Standard error of a, estimated using VAST.} #' \item{species_code}{: RACE/GAP species code} -#' \item{scaled_vast_condition}{: Anomaly (Z-score) of vast_condition} -#' \item{vast_relative_condition}{: Relative condition based on the allometric intercept estimated using VAST, i.e., allometric intercept divided by the time series mean allometric intercept.} -#' \item{vast_relative_condition_se}{: Standard error of vast_relative_condition.} #' } #' \item STRATUM (data frame): Residuals by stratum #' \itemize{ -#' \item{year}{: Year} #' \item{common_name}{: Species Common name} -#' \item{species_code}{: RACE/GAP species code} -#' \item{inpfc_stratum}{: INPFC Stratum} +#' \item{year}{: Year} +#' \item{stratum}{: Survey Stratum} #' \item{stratum_resid_mean}{: Unweighted stratum mean residual} -#' \item{n}{: Sample size for the stratum} #' \item{stratum_resid_sd}{: Standard error of the stratum mean residual} +#' \item{n}{: Sample size for the stratum} #' \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} -#' \item{weighted_resid_mean}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} +#' \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} #' } #' \item LAST_UPDATE: Date when indicator and data were last updated. #' } @@ -111,24 +97,31 @@ #' Northern Bering Sea (2010-2022) #' -#' Morphometric condition indicators based on the allometric intercept estimated using VAST and residuals from a length-weight regression for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), Alaska plaice, and yellowfin sole in the northern Bering Sea. +#' Morphometric condition indicators based on length-weight regressions for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), Alaska plaice, and yellowfin sole in the northern Bering Sea. #' #' @format A list containing a data frames (indicator for the full region) and a character vector. #' \describe{ #' \itemize{ #' \item FULL_REGION (data frame): Residuals for the full region #' \itemize{ -#' \item{YEAR}{: Year} -#' \item{COMMON_NAME}{: Species common name} +#' \item{year}{: Year} +#' \item{common_name}{: Species common name} #' \item{mean_wt_resid}{: Mean residual for the full region} #' \item{se_wt_resid}{: Standard error of the indicator for the full region} -#' \item{vast_condition}{: Estimate of the allometric intercept,a, of the length-weight equation (W = aL^b), estimated using VAST.} -#' \item{vast_condition_se}{Standard error of a, estimated using VAST.} #' \item{species_code}{: RACE/GAP species code} -#' \item{scaled_vast_condition}{: Anomaly (Z-score) of vast_condition} -#' \item{vast_relative_condition}{: Relative condition based on the allometric intercept estimated using VAST, i.e., allometric intercept divided by the time series mean allometric intercept.} #' \item{vast_relative_condition_se}{: Standard error of vast_relative_condition.} #' } +#' \item STRATUM (data frame): Residuals by stratum +#' \itemize{ +#' \item{common_name}{: Species Common name} +#' \item{year}{: Year} +#' \item{stratum}{: Survey stratum} +#' \item{stratum_resid_mean}{: Unweighted stratum mean residual} +#' \item{stratum_resid_sd}{: Standard error of the stratum mean residual} +#' \item{n}{: Sample size for the stratum} +#' \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} +#' \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} +#' } #' \item LAST_UPDATE: Date when indicator and data were last updated. #' } #' } @@ -149,8 +142,6 @@ #' \item{common_name}{: Species common name} #' \item{mean_wt_resid}{: Mean residual for the full region} #' \item{se_wt_resid}{: Standard error of the indicator for the full region} -#' \item{vast_condition}{: Estimate of the allometric intercept,a, of the length-weight equation (W = aL^b), estimated using VAST.} -#' \item{vast_condition_se}{Standard error of a, estimated using VAST.} #' } #' \item STRATUM_* (data frame): Residuals by stratum #' \itemize{ @@ -162,7 +153,7 @@ #' \item{n}{: Sample size for the stratum} #' \item{stratum_resid_sd}{: Standard error of the stratum mean residual} #' \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} -#' \item{weighted_resid_mean}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} +#' \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} #' } #' \item LAST_UPDATE: Date when indicator and data were last updated. #' } diff --git a/R/make_plots.R b/R/make_plots.R deleted file mode 100644 index 4596b89..0000000 --- a/R/make_plots.R +++ /dev/null @@ -1,950 +0,0 @@ -#' Make condition time series plot -#' -#' @param x Input data.frame -#' @param region Region. AI, BS, or GOA. Character (1L). -#' @param fill_color Fill color to use for points. Character (1L). -#' @param var_y_name Name of the y variable in the data.frame. Character (1L). -#' @param var_y_se_name Name of the standard error for the y variable. Character (1L). -#' @param var_x_name Name of the x (time variable). Character (1L). -#' @param y_title Y-axis title. Character (1L). -#' @param format_for "rmd" or "png" -#' @param set_intercept Intercept to use for plots. Numeric (1L). 0 for residual, 1 for vAST relative condition, NULL for VAST a -#' @return A ggplot object of the time series -#' @examples -#' # EBS anomaly timeseries plot -#' -#' names(EBS_INDICATOR$FULL_REGION) -#' -#' plot_anomaly_timeseries(x = EBS_INDICATOR$FULL_REGION, -# 'region = "BS", -#' fill_color = "#0085CA", -#' var_y_name = "mean_wt_resid", -#' var_y_se_name = "se_wt_resid", -#' var_x_name = "year", -#' y_title = "Length-weight residual (ln(g))", -#' format_for = "png") -#' @export - -plot_anomaly_timeseries <- function(x, - region, - fill_color = "#0085CA", - var_y_name = "mean_wt_resid", - var_y_se_name = "se_wt_resid", - var_x_name = "year", - y_title = "Length-weight residual (ln(g))", - format_for = "rmd", - plot_type = "point", - set_intercept = NULL) { - - region <- toupper(region) - - point_rel_size <- c(5, 2.5)[match(tolower(format_for), c("rmd", "png"))] - - stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'BS', 'AI', or 'GOA'" = (region %in% c("BS", "AI", "GOA"))) - - x$display_name <- akfishcondition::set_plot_order(x$common_name, - region = region) - - names(x)[which(names(x) == var_x_name)] <- "var_x" - names(x)[which(names(x) == var_y_name)] <- "var_y" - names(x)[which(names(x) == var_y_se_name)] <- "se_var_y" - - x <- dplyr::filter(x, !is.na(var_y)) - - # Anomalies - x_anomaly <- x %>% - dplyr::group_by(common_name, - display_name) %>% - dplyr::summarise(mean_var_y = mean(var_y), - sd_var_y = sd(var_y)) - - if(!is.null(set_intercept)) { - x_anomaly$mean_var_y <- set_intercept - } - - # Make plot - if(plot_type == "point") { - p1 <- ggplot() + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y), - linetype = 1, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y + sd_var_y), - linetype = 2, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y - sd_var_y), - linetype = 2, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y + 2*sd_var_y), - linetype = 3, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y - 2*sd_var_y), - linetype = 3, - color = "grey50") + - geom_linerange(data = x, - aes(x = var_x, - ymax = var_y + 2*se_var_y, - ymin = var_y - 2*se_var_y)) + - geom_point(data = x, - aes(x = var_x, - y = var_y), - fill = fill_color, - color = "black", - shape = 21, - size = rel(point_rel_size)) + - facet_wrap(~display_name, ncol = 2, scales = "free_y") + - scale_x_continuous(name = "Year", breaks = scales::pretty_breaks(n = 4)) + - scale_y_continuous(name = y_title) - } else if(plot_type == "bar") { - - if(is.null(set_intercept)) { - set_intercept <- 0 - } - - p1 <- ggplot() + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y), - linetype = 1, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y + sd_var_y), - linetype = 2, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y - sd_var_y), - linetype = 2, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y + 2*sd_var_y), - linetype = 3, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y - 2*sd_var_y), - linetype = 3, - color = "grey50") + - geom_linerange(data = x, - aes(x = var_x, - ymax = var_y + 2*se_var_y, - ymin = var_y - 2*se_var_y)) + - geom_bar(data = x, - aes(x = var_x, - y = var_y), - stat = "identity", - fill = fill_color, - color = "black", - width = 1) + - facet_wrap(~display_name, ncol = 2, scales = "free_y") + - scale_x_continuous(name = "Year", - breaks = scales::pretty_breaks(n = 4)) + - scale_y_continuous(name = y_title, - trans = scales::trans_new("shift", - transform = function(x) {x-set_intercept}, - inverse = function(x) {x+set_intercept})) - } else if(plot_type == "box") { - if(is.null(set_intercept)) { - set_intercept <- 0 - } - - p1 <- ggplot() + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y), - linetype = 1, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y + sd_var_y), - linetype = 2, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y - sd_var_y), - linetype = 2, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y + 2*sd_var_y), - linetype = 3, - color = "grey50") + - geom_hline(data = x_anomaly, - mapping = aes(yintercept = mean_var_y - 2*sd_var_y), - linetype = 3, - color = "grey50") + - geom_linerange(data = x, - mapping = aes(x = var_x, - ymax = var_y + 2*se_var_y, - ymin = var_y - 2*se_var_y), - color = "black") + - geom_rect(data = x, - mapping = aes(xmin = var_x+0.4, - xmax = var_x-0.4, - ymax = var_y + se_var_y, - ymin = var_y - se_var_y), - fill = fill_color) + - geom_segment(data = x, - aes(x = var_x+0.4, - xend = var_x-0.4, - y = var_y, - yend = var_y, - group = var_x), - color = "black") + - facet_wrap(~display_name, ncol = 2, scales = "free_y") + - scale_x_continuous(name = "Year", - breaks = scales::pretty_breaks(n = 4)) + - scale_y_continuous(name = y_title, - trans = scales::trans_new("shift", - transform = function(x) {x-set_intercept}, - inverse = function(x) {x+set_intercept})) - } - - return(p1) - -} - - - -#' Make condition time series stacked bar plot -#' -#' @param x Input data.frame -#' @param region Region. AI, BS, or GOA. Character (1L). -#' @param var_y_name Name of the y variable in the data.frame. Character (1L). -#' @param var_y_se_name Name of the standard error for the y variable. Character (1L). -#' @param var_x_name Name of the x (time variable). Character (1L). -#' @param var_group_name Name of the variable to use for grouping (i.e., stratum). Character (1L). -#' @param y_title Y-axis title. Character (1L). -#' @param fill_title Name of the fill variable to use for plotting. Character (1L). -#' @return A ggplot object of the time series -#' @examples -#' # EBS strata stacked bar plot -#' -#' names(EBS_INDICATOR$STRATUM) -#' -#' plot_stratum_stacked_bar(x = EBS_INDICATOR$STRATUM, -#' region = "BS", -#' fill_palette = "BrBG", -#' var_y_name = "stratum_resid_mean", -#' var_x_name = "year", -#' var_group_name = "stratum", -#' fill_title = "Stratum", -#' y_title = "Length-weight residual (ln(g))") -#' @export - -plot_stratum_stacked_bar <- function(x, - region, - fill_palette = "BrBG", - var_y_name = "mean_wt_resid", - var_x_name = "year", - var_group_name = "stratum", - fill_title = "Stratum", - y_title = "Length-weight residual (ln(g))") { - - region <- toupper(region) - - stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'BS', 'AI', or 'GOA'" = (region %in% c("BS", "AI", "GOA"))) - - x$display_name <- akfishcondition::set_plot_order(x$common_name, - region = region) - - names(x)[which(names(x) == var_x_name)] <- "var_x" - names(x)[which(names(x) == var_y_name)] <- "var_y" - names(x)[which(names(x) == var_group_name)] <- "var_group" - - x <- dplyr::filter(x, !is.na(var_y)) - - p1 <- ggplot(data = x, - aes(x = var_x, - y = var_y, - fill = set_stratum_order(trimws(var_group), - region = region))) + - geom_hline(yintercept = 0) + - geom_bar(stat = "identity", - color = "black", - position = "stack", - width = 1) + - facet_wrap(~display_name, ncol = 2, scales = "free_y") + - scale_x_continuous(name = "Year") + - scale_y_continuous(name = y_title) + - scale_fill_brewer(name = fill_title, palette = fill_palette) - - return(p1) - -} - - - -#' Make single species stratum bar plots -#' -#' Make condition time series bar plots with separate panels for each stratum. -#' -#' @param x Input data.frame -#' @param region Region. AI, BS, or GOA. Character (1L). -#' @param var_x_name Name of the x (time variable). Character (1L). -#' @param var_y_name Name of the y variable in the data.frame. Character (1L). -#' @param var_y_se_name Name of the standard error for the y variable. Character (1L). -#' @param y_title Y-axis title. Character (1L). -#' @param fill_title Name of the fill variable to use for plotting. Character (1L). -#' @param fill_pallete Character vector denoting which color pallette to use. Must be a valid name for RColorBrewer::brewer.pal(name = {fill_palette}) -#' @param var_group_name Name of the variable to use for grouping (i.e., stratum). Character (1L). -#' @param write_plot Should plots be written to the /plot/ directory? -#' @return A list of bar plots as ggplot objects -#' @examples -#' # EBS single species stratum bar plots -#' -#' names(EBS_INDICATOR$STRATUM) -#' -#' plot_species_stratum_bar( -#' x = EBS_INDICATOR$STRATUM, -#' region = "BS", -#' var_x_name = "year", -#' var_y_name = "stratum_resid_mean", -#' var_y_se_name = "stratum_resid_se", -#' y_title = "Length-weight residual (ln(g))", -#' var_group_name = "stratum", -#' fill_title = "Stratum", -#' fill_palette = "BrBG", -#' write_plot = FALSE) -#' @export - -plot_species_stratum_bar <- function(x, - region, - var_x_name, - var_y_name, - var_y_se_name, - y_title = "Length-weight residual (ln(g))", - var_group_name = "stratum", - fill_title = "Stratum", - fill_palette = "BrBG", - write_plot = TRUE) { - - region <- toupper(region) - - stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'EBS', 'NBS', 'AI', or 'GOA'" = (region %in% c("EBS", "NBS", "AI", "GOA"))) - - if(region %in% c("EBS", "NBS")) { - region_dir <- "BS" - region_order <- "BS" - } else { - region_dir <- region - region_order <- region - } - - names(x)[which(names(x) == var_x_name)] <- "var_x" - names(x)[which(names(x) == var_y_name)] <- "var_y" - names(x)[which(names(x) == var_group_name)] <- "var_group" - names(x)[which(names(x) == var_y_se_name)] <- "var_y_se" - x <- dplyr::filter(x, !is.na(var_y)) - - x$display_name <- akfishcondition::set_plot_order(x$common_name, region = region_order) - - unique_groups <- unique(x$display_name) - - out_list <- list() - for(ii in 1:length(unique_groups)) { - - if(is.na(unique_groups[ii])) { - next - } - - spp_dir <- gsub(pattern = ">", replacement = "gt", unique_groups[ii]) - - if(!dir.exists(here::here("plots", region_dir, spp_dir))) { - dir.create(here::here("plots", region_dir, spp_dir), recursive = TRUE) - } - - p1 <- - ggplot(data = x %>% - dplyr::filter(display_name == unique_groups[ii]), - aes(x = var_x, - y = var_y, - fill = set_stratum_order(trimws(var_group), region = region_order), - ymin = var_y - 2*var_y_se, - ymax = var_y + 2*var_y_se)) + - geom_hline(yintercept = 0) + - geom_linerange() + - geom_bar(stat = "identity", - color = "black", - position = "stack", - width = 1) + - facet_wrap(~set_stratum_order(trimws(var_group), - region = region_order), - ncol = 2, - scales = "free_y") + - ggtitle(unique_groups[ii]) + - scale_x_continuous(name = "Year", breaks = scales::pretty_breaks()) + - scale_y_continuous(name = y_title) + - scale_fill_brewer(name = fill_title, - palette = fill_palette, - drop = FALSE) + - theme(legend.position = "none", - title = element_text(hjust = 0.5)) - - - if(write_plot) { - - png(here::here("plots", region_dir, spp_dir, - paste0(region, "_STRATUM_CONDITION_", spp_dir, ".png")), - width = 6, height = 7, units = "in", res = 600) - print(p1 + theme_blue_strip() + - theme(legend.position = "none", - title = element_text(hjust = 0.5))) - dev.off() - } - - out_list[ii] <- p1 - - } - - - return(out_list) - -} - - - -#' Make single species bar plots -#' -#' Make condition time series bar plots for single species. -#' -#' @param x Input data.frame -#' @param region Region. AI, BS, or GOA. Character (1L). -#' @param var_x_name Name of the x (time variable). Character (1L). -#' @param var_y_name Name of the y variable in the data.frame. Character (1L). -#' @param var_y_se_name Name of the standard error for the y variable. Character (1L). -#' @param y_title Y-axis title. Character (1L). -#' @param fill_colors Fill color to use for bars. Character (1L). -#' @param set_intercept Intercept to use for plots. Numeric (1L). 0 for residual, 1 for vAST relative condition, NULL for VAST a -#' @param write_plot Should plots be written to the /plot/ directory? -#' @return A list of bar plots as ggplot objects -#' @examples -#' # EBS single species stratum bar plots -#' -#' names(EBS_INDICATOR$STRATUM) -#' -#' plot_species_stratum_bar( -#' x = EBS_INDICATOR$STRATUM, -#' region = "BS", -#' var_x_name = "year", -#' var_y_name = "stratum_resid_mean", -#' var_y_se_name = "stratum_resid_se", -#' y_title = "Length-weight residual (ln(g))", -#' write_plot = FALSE) -#' @export - -plot_species_bar <- function(x, - region, - var_x_name = "year", - var_y_name, - var_y_se_name, - fill_color = "#0085CA", - y_title = "Length-weight residual (ln(g))", - set_intercept = NULL, - write_plot = TRUE) { - - region <- toupper(region) - - stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'EBS', 'NBS', 'AI', or 'GOA'" = (region %in% c("EBS", "NBS", "AI", "GOA"))) - - if(region %in% c("EBS", "NBS")) { - region_dir <- "BS" - region_order <- "BS" - } else { - region_dir <- region - region_order <- region - } - - x$display_name <- akfishcondition::set_plot_order(x$common_name, - region = region_order) - - names(x)[which(names(x) == var_x_name)] <- "var_x" - names(x)[which(names(x) == var_y_name)] <- "var_y" - names(x)[which(names(x) == var_y_se_name)] <- "var_y_se" - x <- dplyr::filter(x, !is.na(var_y)) - - # Anomalies - x_anomaly <- x %>% - dplyr::group_by(common_name, - display_name) %>% - dplyr::summarise(mean_var_y = mean(var_y), - sd_var_y = sd(var_y)) - - if(!is.null(set_intercept)) { - x_anomaly$mean_var_y <- set_intercept - } - - unique_groups <- unique(x$display_name) - - - - out_list <- list() - for(ii in 1:length(unique_groups)) { - - if(is.na(unique_groups[ii])) { - next - } - - spp_dir <- gsub(pattern = ">", replacement = "gt", unique_groups[ii]) - - if(!dir.exists(here::here("plots", region_dir, spp_dir))) { - dir.create(here::here("plots", region_dir, spp_dir), recursive = TRUE) - } - - if(is.null(set_intercept)) { - set_intercept <- 0 - } - - sel_anomaly <- dplyr::filter(x_anomaly, display_name == unique_groups[ii]) - - p1 <- - ggplot() + - geom_hline(data = sel_anomaly, - mapping = aes(yintercept = mean_var_y), - linetype = 1, - color = "grey50") + - geom_hline(data = sel_anomaly, - mapping = aes(yintercept = mean_var_y + sd_var_y), - linetype = 2, - color = "grey50") + - geom_hline(data = sel_anomaly, - mapping = aes(yintercept = mean_var_y - sd_var_y), - linetype = 2, - color = "grey50") + - geom_hline(data = sel_anomaly, - mapping = aes(yintercept = mean_var_y + 2*sd_var_y), - linetype = 3, - color = "grey50") + - geom_hline(data = sel_anomaly, - mapping = aes(yintercept = mean_var_y - 2*sd_var_y), - linetype = 3, - color = "grey50") + - geom_linerange(data = x %>% - dplyr::filter(display_name == unique_groups[ii]), - mapping = aes(x = var_x, - ymax = var_y + 2*var_y_se, - ymin = var_y - 2*var_y_se), - color = "black") + - geom_rect(data = x %>% - dplyr::filter(display_name == unique_groups[ii]), - mapping = aes(xmin = var_x+0.4, - xmax = var_x-0.4, - ymax = var_y + var_y_se, - ymin = var_y - var_y_se), - fill = fill_color) + - geom_segment(data = x %>% - dplyr::filter(display_name == unique_groups[ii]), - aes(x = var_x+0.4, - xend = var_x-0.4, - y = var_y, - yend = var_y, - group = var_x), - color = "black") + - # geom_bar(stat = "identity", - # color = "black", - # fill = fill_color, - # width = 1) + - # geom_linerange() + - facet_grid(~display_name) + - scale_x_continuous(name = "Year", breaks = scales::pretty_breaks(n = 4)) + - scale_y_continuous(name = y_title, - trans = scales::trans_new("shift", - transform = function(x) {x-set_intercept}, - inverse = function(x) {x+set_intercept})) + - scale_fill_manual(name = ) - theme(legend.position = "none", - title = element_text(hjust = 0.5)) - - - if(write_plot) { - - png(here::here("plots", region_dir, spp_dir, - paste0(region, "_CONDITION_", spp_dir, ".png")), - width = 6, height = 4, units = "in", res = 600) - print(p1 + theme_blue_strip() + - theme(legend.position = "none", - title = element_text(hjust = 0.5))) - dev.off() - } - - out_list[ii] <- p1 - - } - - - return(out_list) - -} - - - -#' Two condition time series and correlation -#' -#' @param x_1 A data.frame for the first time series. -#' @param x_2 A data.frame for the second time series. -#' @param region Region. AI, BS, or GOA. Character (1L). -#' @param fill_colors Fill colors to use for points. Character (2L). -#' @param var_y_name_1 Name of the y variable in the data.frame. Character (1L). -#' @param var_x_name_1 Name of the x (time variable). Character (1L). -#' @param var_y_name_2 Name of the y variable in the data.frame. Character (1L). -#' @param var_x_name_2 Name of the x (time variable). Character (1L). -#' @param scale_y Logical. Should y variables be Z-score transformed? -#' @param year_break Optional. Year to split a times series for geom_line(); denotes a break in the time series (i.e., EBS and NBS in 2020) -#' @param y_title Y-axis title. Character (1L). -#' @param x_offset Offset value to add to an x variable (e.g. year) from one of the time series to improve interpretability. Numeric (1L). -#' @param fill_title Name of the fill variable to use for plotting. Character (1L). -#' @param fill_color Fill colors to use for points. -#' @param shapes Shapes to use for points. -#' @param format_for "rmd" or "png" -#' @return A ggplot object of the time series -#' @export - -plot_two_timeseries <- function(x_1, - x_2, - region, - series_name_1 = "Historical", - var_y_name_1 = "mean_wt_resid", - var_x_name_1 = "year", - series_name_2 = "VAST", - var_y_name_2 = "mean_wt_resid", - var_x_name_2 = "year", - var_group_name = "common_name", - y_title = "Condition (Z-score)", - fill_title = "Method", - scale_y = TRUE, - year_break = NULL, - x_offset = 0.25, - fill_colors = c("#001743", "#0085CA"), - shapes = c(24, 21), - format_for = "rmd") { - - region <- toupper(region) - - point_rel_size <- c(5, 2.5)[match(tolower(format_for), c("rmd", "png"))] - - stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'BS', 'AI', or 'GOA'" = (region %in% c("BS", "AI", "GOA"))) - - shared_cols <- c("var_x", "var_y", "var_y_se", "display_name", "common_name", "series") - - # Setup time series 1 - - names(x_1)[which(names(x_1) == var_group_name)] <- "var_group" - - if(var_group_name == "common_name") { - x_1$display_name <- akfishcondition::set_plot_order(x_1$var_group, - region = region) - } else { - x_1$display_name <- x_1$var_group - } - - names(x_1)[which(names(x_1) == var_x_name_1)] <- "var_x" - names(x_1)[which(names(x_1) == var_y_name_1)] <- "var_y" - x_1$series <- series_name_1 - x_1 <- x_1[, which(names(x_1) %in% shared_cols)] - x_1 <- dplyr::filter(x_1, !is.na(var_y)) - - # Setup time series 2 - names(x_2)[which(names(x_2) == var_group_name)] <- "var_group" - - if(var_group_name == "common_name") { - x_2$display_name <- akfishcondition::set_plot_order(x_2$var_group, - region = region) - } else { - x_2$display_name <- x_2$var_group - } - - names(x_2)[which(names(x_2) == var_x_name_2)] <- "var_x" - names(x_2)[which(names(x_2) == var_y_name_2)] <- "var_y" - x_2$series <- series_name_2 - x_2 <- x_2[, which(names(x_2) %in% shared_cols)] - x_2 <- dplyr::filter(x_2, !is.na(var_y)) - - # Scale y variables - if(scale_y) { - x_1 <- x_1 %>% - dplyr::group_by(common_name, display_name) %>% - dplyr::mutate(var_y = scale(var_y)[,1]) - x_2 <- x_2 %>% - dplyr::group_by(common_name, display_name) %>% - dplyr::mutate(var_y = scale(var_y)[,1]) - } - - - # Correlation between timeseries - corr_df <- dplyr::inner_join( - x_1 %>% - dplyr::ungroup() %>% - dplyr::rename(var_y_1 = var_y) %>% - dplyr::select(display_name, - var_x, - var_y_1), - x_2 %>% - dplyr::ungroup() %>% - dplyr::rename(var_y_2 = var_y) %>% - dplyr::select(display_name, - var_x, - var_y_2), - by = c("display_name", "var_x")) %>% - dplyr::group_by(display_name) %>% - dplyr::summarise( - r = round(cor(var_y_1, - var_y_2, - use = "complete.obs", - method = "pearson"), 2)) - - x_2$var_x <- x_2$var_x + x_offset - x_combined <- dplyr::bind_rows(x_1, x_2) - - if(!is.null(year_break)) { - x_combined$grp <- x_combined$var_x < year_break - } else { - x_combined$grp <- 1 - } - - max_x <- max(x_combined$var_x) - min_x <- min(x_combined$var_x) - lab_x <- max_x - (max_x - min_x)*0.12 - - min_y <- min(x_combined$var_y) - max_y <- max(x_combined$var_y) - lab_y <- min_y + (max_y - min_y)*0.07 - - # Make plot - p1 <- ggplot() + - geom_hline(yintercept = 0, - linetype = 1, - color = "grey50") + - geom_hline(yintercept = c(-1,1), - linetype = 2, - color = "grey50") + - geom_hline(yintercept = c(-2,2), - linetype = 3, - color = "grey50") + - geom_point(data = x_combined, - aes(x = var_x, - y = var_y, - fill = series, - shape = series), - color = "black", - size = rel(point_rel_size)) + - geom_text(data = corr_df, - aes(x = lab_x, - y = lab_y, - label = paste0("r = ", format(r, nsmall = 2))), - size = rel(point_rel_size*1.9)) + - scale_fill_manual(name = fill_title, - values = fill_colors) + - scale_shape_manual(name = fill_title, - values = shapes) + - facet_wrap(~display_name, ncol = 2) + - scale_x_continuous(name = "Year", breaks = scales::pretty_breaks(n = 4)) + - scale_y_continuous(name = y_title, limits = c(range(x_combined$var_y))) - - return(p1) - -} - - - -#' Plot of condition versus condition -#' -#' @param x_1 A data.frame for the first time series. -#' @param x_2 A data.frame for the second time series. -#' @param region Region. AI, BS, or GOA. Character (1L). -#' @param fill_colors Fill colors to use for points. Character (2L). -#' @param var_y_name_1 Name of the y variable in the data.frame. Character (1L). -#' @param var_y_se_name_1 Name of the standard error for the y variable. Character (1L). -#' @param var_x_name_1 Name of the x (time variable). Character (1L). -#' @param var_y_name_2 Name of the y variable in the data.frame. Character (1L). -#' @param var_y_se_name_2 Name of the standard error for the y variable. Character (1L). -#' @param var_x_name_2 Name of the x (time variable). Character (1L). -#' @param scale_y Logical. Should y variables be Z-score transformed? -#' @param y_title Y-axis title. Character (1L). -#' @param format_for "rmd" or "png" -#' @return A ggplot object of the time series -#' @export - -plot_xy_corr <- function(x_1, - x_2, - region, - series_name_1 = "Historical", - var_y_name_1 = "mean_wt_resid", - var_y_se_name_1 = "se_wt_resid", - var_x_name_1 = "year", - series_name_2 = "VAST", - var_y_name_2 = "mean_wt_resid", - var_y_se_name_2 = "se_wt_resid", - var_x_name_2 = "year", - format_for = "rmd") { - - region <- toupper(region) - - point_rel_size <- c(5, 2.5)[match(tolower(format_for), c("rmd", "png"))] - - stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'BS', 'AI', or 'GOA'" = (region %in% c("BS", "AI", "GOA"))) - - shared_cols <- c("var_x", "var_y", "var_y_se", "display_name", "common_name", "series") - - # Setup time series 1 - x_1$display_name <- akfishcondition::set_plot_order(x_1$common_name, - region = region) - names(x_1)[which(names(x_1) == var_x_name_1)] <- "var_x" - names(x_1)[which(names(x_1) == var_y_name_1)] <- "var_y" - names(x_1)[which(names(x_1) == var_y_se_name_1)] <- "var_y_se" - x_1$series <- series_name_1 - x_1 <- x_1[, which(names(x_1) %in% shared_cols)] - x_1 <- dplyr::filter(x_1, !is.na(var_y)) - - # Setup time series 2 - x_2$display_name <- akfishcondition::set_plot_order(x_2$common_name, - region = region) - names(x_2)[which(names(x_2) == var_x_name_2)] <- "var_x" - names(x_2)[which(names(x_2) == var_y_name_2)] <- "var_y" - names(x_2)[which(names(x_2) == var_y_se_name_2)] <- "var_y_se" - x_2$series <- series_name_2 - x_2 <- x_2[, which(names(x_2) %in% shared_cols)] - x_2 <- dplyr::filter(x_2, !is.na(var_y)) - - # Correlation between timeseries - x_combined <- dplyr::inner_join( - x_1 %>% - dplyr::ungroup() %>% - dplyr::rename(var_y_1 = var_y, - var_y_se_1 = var_y_se) %>% - dplyr::select(display_name, - var_x, - var_y_1, - var_y_se_1), - x_2 %>% - dplyr::ungroup() %>% - dplyr::rename(var_y_2 = var_y, - var_y_se_2 = var_y_se) %>% - dplyr::select(display_name, - var_x, - var_y_2, - var_y_se_2), - by = c("display_name", "var_x")) %>% - dplyr::ungroup() - - - corr_df <- x_combined %>% - dplyr::group_by(display_name) %>% - dplyr::summarise( - r = round(cor(var_y_1, - var_y_2, - use = "complete.obs", - method = "pearson"), 2)) - - corr_df <- corr_df %>% - dplyr::inner_join( - x_combined %>% - dplyr::group_by(display_name) %>% - dplyr::summarise(min_x = min(var_y_1-2*var_y_se_1), - max_x = max(var_y_1+2*var_y_se_1), - min_y = min(var_y_2-2*var_y_se_2), - max_y = max(var_y_2+2*var_y_se_2)) %>% - dplyr::mutate(lab_x = max_x - (max_x - min_x)*0.12, - lab_y = min_y + (max_y - min_y)*0.07), - by = "display_name") - - # Make plot - p1 <- ggplot() + - geom_linerange(data = x_combined, - aes(x = var_y_1, - ymin = var_y_2 - 2 * var_y_se_2, - ymax = var_y_2 + 2 * var_y_se_2), - color = "black") + - geom_linerange(data = x_combined, - aes(y = var_y_2, - xmin = var_y_1 - 2 * var_y_se_1, - xmax = var_y_1 + 2 * var_y_se_1), - color = "black") + - geom_point(data = x_combined, - aes(x = var_y_1, - y = var_y_2), - color = "black", - size = rel(point_rel_size)) + - geom_text(data = corr_df, - aes(x = lab_x, - y = lab_y, - label = paste0("r = ", format(r, nsmall = 2))), - size = rel(point_rel_size*1.9)) + - facet_wrap(~display_name, ncol = 2, scale = "free") + - scale_x_continuous(name = series_name_1, breaks = scales::pretty_breaks()) + - scale_y_continuous(name = series_name_2) - - return(p1) - -} - - - -#' Generate tables with b coefficient and VAST Settings -#' -#' Generate tables containing settings for VAST and the VAST estimate of the allometric slope of the length-weight (b in W = a*L^b) relationship. -#' -#' @param region One or more region as a character vector ("EBS", "NBS", "AI", "GOA") -#' @param write_table Should the table be written as a .csv file to /plots/[region]/[region]_VAST_summary_table.csv? -#' @export - -make_vast_table <- function(region, write_table = TRUE) { - - library(VAST) - - output_region_dir <- c("BS", "BS", "AI", "GOA")[match(region, c("EBS", "NBS", "AI", "GOA"))] - - sel_region <- region - - sel_spp <- dplyr::filter(akfishcondition::ESR_SETTINGS$VAST_SETTINGS, - region %in% sel_region ) - - out_df <- data.frame() - - for(ii in 1:length(sel_spp)) { - vast_file <- here::here("results", - sel_spp$region[ii], - sel_spp$species_code[ii], - paste0(sel_spp$species_code[ii], "_VASTfit.rds")) - - if(file.exists(vast_file)) { - message("make_vast_tables: Reading VASTfit.rds for species ", sel_spp$species_code[ii], " in ", sel_spp$region[ii], ".") - vast_fit <- readRDS(file = vast_file) - - if(!is.na(vast_fit$parameter_estimates$SD$par.fixed['log_sigmaPhi2_k'])) { - bb <- vast_fit$parameter_estimates$SD$par.fixed['lambda2_k'] - se <- exp(vast_fit$parameter_estimates$SD$par.fixed['log_sigmaPhi2_k']) - } else { - fixed_summary <- summary(vast_fit$parameter_estimates$SD, "fixed") - bb <- fixed_summary[which(row.names(fixed_summary) == "lambda2_k"), 1] - se <- fixed_summary[which(row.names(fixed_summary) == "lambda2_k"), 2] - } - - - out_df <- dplyr::bind_rows(out_df, - dplyr::bind_cols(sel_spp[ii,], - data.frame(b_est = bb, - b_se = se, - max_gradient = vast_fit$parameter_estimates['max_gradient']$max_gradient, - convergence_check = vast_fit$parameter_estimates['Convergence_check']$Convergence_check, - knots_match = vast_fit$spatial_list$n_x == sel_spp$n_knots[ii]))) - - } else { - warning("\nmake_vast_tables: VASTfit file for species ", sel_spp$species_code[ii], " in ", sel_spp$region[ii], " cannot be found in ", vast_file, "\n") - } - } - - if(write_table) { - if(length(region) > 1) { - out_loc <- here::here("plots", paste0(paste0(region, collapse = "_"), "_VAST_summary_table.csv")) - } else { - out_loc <- here::here("plots", output_region_dir, paste0(region, "_VAST_summary_table.csv")) - } - - message("make_vast_tables: Writing table to ", out_loc) - - write.csv(x = out_df, file = out_loc, row.names = FALSE) - - } - - row.names(out_df) <- 1:nrow(out_df) - - return(out_df) - -} \ No newline at end of file diff --git a/R/make_vast_table.R b/R/make_vast_table.R new file mode 100644 index 0000000..450d581 --- /dev/null +++ b/R/make_vast_table.R @@ -0,0 +1,72 @@ +#' Generate tables with b coefficient and VAST Settings +#' +#' Generate tables containing settings for VAST and the VAST estimate of the allometric slope of the length-weight (b in W = a*L^b) relationship. +#' +#' @param region One or more region as a character vector ("EBS", "NBS", "AI", "GOA") +#' @param write_table Should the table be written as a .csv file to /plots/{region{/{region}_VAST_summary_table.csv? +#' @export + +make_vast_table <- function(region, write_table = TRUE) { + + library(VAST) + + output_region_dir <- c("BS", "BS", "AI", "GOA")[match(region, c("EBS", "NBS", "AI", "GOA"))] + + sel_region <- region + + sel_spp <- dplyr::filter(akfishcondition::ESR_SETTINGS$VAST_SETTINGS, + region %in% sel_region ) + + out_df <- data.frame() + + for(ii in 1:length(sel_spp)) { + vast_file <- here::here("results", + sel_spp$region[ii], + sel_spp$species_code[ii], + paste0(sel_spp$species_code[ii], "_VASTfit.rds")) + + if(file.exists(vast_file)) { + message("make_vast_tables: Reading VASTfit.rds for species ", sel_spp$species_code[ii], " in ", sel_spp$region[ii], ".") + vast_fit <- readRDS(file = vast_file) + + if(!is.na(vast_fit$parameter_estimates$SD$par.fixed['log_sigmaPhi2_k'])) { + bb <- vast_fit$parameter_estimates$SD$par.fixed['lambda2_k'] + se <- exp(vast_fit$parameter_estimates$SD$par.fixed['log_sigmaPhi2_k']) + } else { + fixed_summary <- summary(vast_fit$parameter_estimates$SD, "fixed") + bb <- fixed_summary[which(row.names(fixed_summary) == "lambda2_k"), 1] + se <- fixed_summary[which(row.names(fixed_summary) == "lambda2_k"), 2] + } + + + out_df <- dplyr::bind_rows(out_df, + dplyr::bind_cols(sel_spp[ii,], + data.frame(b_est = bb, + b_se = se, + max_gradient = vast_fit$parameter_estimates['max_gradient']$max_gradient, + convergence_check = vast_fit$parameter_estimates['Convergence_check']$Convergence_check, + knots_match = vast_fit$spatial_list$n_x == sel_spp$n_knots[ii]))) + + } else { + warning("\nmake_vast_tables: VASTfit file for species ", sel_spp$species_code[ii], " in ", sel_spp$region[ii], " cannot be found in ", vast_file, "\n") + } + } + + if(write_table) { + if(length(region) > 1) { + out_loc <- here::here("plots", paste0(paste0(region, collapse = "_"), "_VAST_summary_table.csv")) + } else { + out_loc <- here::here("plots", output_region_dir, paste0(region, "_VAST_summary_table.csv")) + } + + message("make_vast_tables: Writing table to ", out_loc) + + write.csv(x = out_df, file = out_loc, row.names = FALSE) + + } + + row.names(out_df) <- 1:nrow(out_df) + + return(out_df) + +} \ No newline at end of file diff --git a/R/plot_anomaly_timeseries.R b/R/plot_anomaly_timeseries.R new file mode 100644 index 0000000..d6852bf --- /dev/null +++ b/R/plot_anomaly_timeseries.R @@ -0,0 +1,204 @@ +#' Make condition time series plot +#' +#' @param x Input data.frame +#' @param region Region. AI, BS, or GOA. Character (1L). +#' @param fill_color Fill color to use for points. Character (1L). +#' @param var_y_name Name of the y variable in the data.frame. Character (1L). +#' @param var_y_se_name Name of the standard error for the y variable. Character (1L). +#' @param var_x_name Name of the x (time variable). Character (1L). +#' @param y_title Y-axis title. Character (1L). +#' @param format_for "rmd" or "png" +#' @param set_intercept Intercept to use for plots. Numeric (1L). 0 for residual, 1 for VAST relative condition, NULL for VAST a +#' @return A ggplot object of the time series +#' @examples +#' \dontrun{ +#' # EBS anomaly timeseries plot +#' +#' names(EBS_INDICATOR$FULL_REGION) +#' +#' plot_anomaly_timeseries(x = EBS_INDICATOR$FULL_REGION, +# 'region = "BS", +#' fill_color = "#0085CA", +#' var_y_name = "mean_wt_resid", +#' var_y_se_name = "se_wt_resid", +#' var_x_name = "year", +#' y_title = "Length-weight residual (ln(g))", +#' format_for = "png") +#' } +#' @export + +plot_anomaly_timeseries <- function(x, + region, + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + format_for = "rmd", + plot_type = "point", + set_intercept = NULL) { + + region <- toupper(region) + + point_rel_size <- c(5, 2.5)[match(tolower(format_for), c("rmd", "png"))] + + stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'BS', 'AI', or 'GOA'" = (region %in% c("BS", "AI", "GOA"))) + + x$display_name <- akfishcondition::set_plot_order(x$common_name, + region = region) + + names(x)[which(names(x) == var_x_name)] <- "var_x" + names(x)[which(names(x) == var_y_name)] <- "var_y" + names(x)[which(names(x) == var_y_se_name)] <- "se_var_y" + + x <- dplyr::filter(x, !is.na(var_y)) + + # Anomalies + x_anomaly <- x |> + dplyr::group_by(common_name, + display_name) |> + dplyr::summarise(mean_var_y = mean(var_y), + sd_var_y = sd(var_y)) + + if(!is.null(set_intercept)) { + x_anomaly$mean_var_y <- set_intercept + } + + # Make plot + if(plot_type == "point") { + p1 <- ggplot() + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y), + linetype = 1, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y + sd_var_y), + linetype = 2, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y - sd_var_y), + linetype = 2, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y + 2*sd_var_y), + linetype = 3, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y - 2*sd_var_y), + linetype = 3, + color = "grey50") + + geom_linerange(data = x, + aes(x = var_x, + ymax = var_y + 2*se_var_y, + ymin = var_y - 2*se_var_y)) + + geom_point(data = x, + aes(x = var_x, + y = var_y), + fill = fill_color, + color = "black", + shape = 21, + size = rel(point_rel_size)) + + facet_wrap(~display_name, ncol = 2, scales = "free_y") + + scale_x_continuous(name = "Year", breaks = scales::pretty_breaks(n = 4)) + + scale_y_continuous(name = y_title) + } else if(plot_type == "bar") { + + if(is.null(set_intercept)) { + set_intercept <- 0 + } + + p1 <- ggplot() + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y), + linetype = 1, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y + sd_var_y), + linetype = 2, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y - sd_var_y), + linetype = 2, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y + 2*sd_var_y), + linetype = 3, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y - 2*sd_var_y), + linetype = 3, + color = "grey50") + + geom_linerange(data = x, + aes(x = var_x, + ymax = var_y + 2*se_var_y, + ymin = var_y - 2*se_var_y)) + + geom_bar(data = x, + aes(x = var_x, + y = var_y), + stat = "identity", + fill = fill_color, + color = "black", + width = 1) + + facet_wrap(~display_name, ncol = 2, scales = "free_y") + + scale_x_continuous(name = "Year", + breaks = scales::pretty_breaks(n = 4)) + + scale_y_continuous(name = y_title, + trans = scales::trans_new("shift", + transform = function(x) {x-set_intercept}, + inverse = function(x) {x+set_intercept})) + } else if(plot_type == "box") { + if(is.null(set_intercept)) { + set_intercept <- 0 + } + + p1 <- ggplot() + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y), + linetype = 1, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y + sd_var_y), + linetype = 2, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y - sd_var_y), + linetype = 2, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y + 2*sd_var_y), + linetype = 3, + color = "grey50") + + geom_hline(data = x_anomaly, + mapping = aes(yintercept = mean_var_y - 2*sd_var_y), + linetype = 3, + color = "grey50") + + geom_linerange(data = x, + mapping = aes(x = var_x, + ymax = var_y + 2*se_var_y, + ymin = var_y - 2*se_var_y), + color = "black") + + geom_rect(data = x, + mapping = aes(xmin = var_x+0.4, + xmax = var_x-0.4, + ymax = var_y + se_var_y, + ymin = var_y - se_var_y), + fill = fill_color) + + geom_segment(data = x, + aes(x = var_x+0.4, + xend = var_x-0.4, + y = var_y, + yend = var_y, + group = var_x), + color = "black") + + facet_wrap(~display_name, ncol = 2, scales = "free_y") + + scale_x_continuous(name = "Year", + breaks = scales::pretty_breaks(n = 4)) + + scale_y_continuous(name = y_title, + trans = scales::trans_new("shift", + transform = function(x) {x-set_intercept}, + inverse = function(x) {x+set_intercept})) + } + + return(p1) + +} \ No newline at end of file diff --git a/R/plot_species_bar.R b/R/plot_species_bar.R new file mode 100644 index 0000000..83b404d --- /dev/null +++ b/R/plot_species_bar.R @@ -0,0 +1,166 @@ +#' Make single species bar plots +#' +#' Make condition time series bar plots for single species. +#' +#' @param x Input data.frame +#' @param region Region. AI, BS, or GOA. Character (1L). +#' @param var_x_name Name of the x (time variable). Character (1L). +#' @param var_y_name Name of the y variable in the data.frame. Character (1L). +#' @param var_y_se_name Name of the standard error for the y variable. Character (1L). +#' @param y_title Y-axis title. Character (1L). +#' @param fill_colors Fill color to use for bars. Character (1L). +#' @param set_intercept Intercept to use for plots. Numeric (1L). 0 for residual, 1 for VAST relative condition, NULL for VAST a +#' @param write_plot Should plots be written to the /plot/ directory? +#' @return A list of bar plots as ggplot objects +#' @examples +#' # EBS single species stratum bar plots +#' +#' names(EBS_INDICATOR$STRATUM) +#' +#' plot_species_stratum_bar( +#' x = EBS_INDICATOR$STRATUM, +#' region = "BS", +#' var_x_name = "year", +#' var_y_name = "stratum_resid_mean", +#' var_y_se_name = "stratum_resid_se", +#' y_title = "Length-weight residual (ln(g))", +#' write_plot = FALSE) +#' @export + +plot_species_bar <- function(x, + region, + var_x_name = "year", + var_y_name, + var_y_se_name, + fill_color = "#0085CA", + y_title = "Length-weight residual (ln(g))", + set_intercept = NULL, + write_plot = TRUE) { + + region <- toupper(region) + + stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'EBS', 'NBS', 'AI', or 'GOA'" = (region %in% c("EBS", "NBS", "AI", "GOA"))) + + if(region %in% c("EBS", "NBS")) { + region_dir <- "BS" + region_order <- "BS" + } else { + region_dir <- region + region_order <- region + } + + x$display_name <- akfishcondition::set_plot_order(x$common_name, + region = region_order) + + names(x)[which(names(x) == var_x_name)] <- "var_x" + names(x)[which(names(x) == var_y_name)] <- "var_y" + names(x)[which(names(x) == var_y_se_name)] <- "var_y_se" + x <- dplyr::filter(x, !is.na(var_y)) + + # Anomalies + x_anomaly <- x |> + dplyr::group_by(common_name, + display_name) |> + dplyr::summarise(mean_var_y = mean(var_y), + sd_var_y = sd(var_y)) + + if(!is.null(set_intercept)) { + x_anomaly$mean_var_y <- set_intercept + } + + unique_groups <- unique(x$display_name) + + + + out_list <- list() + for(ii in 1:length(unique_groups)) { + + if(is.na(unique_groups[ii])) { + next + } + + spp_dir <- gsub(pattern = ">", replacement = "gt", unique_groups[ii]) + + if(!dir.exists(here::here("plots", region_dir, spp_dir))) { + dir.create(here::here("plots", region_dir, spp_dir), recursive = TRUE) + } + + if(is.null(set_intercept)) { + set_intercept <- 0 + } + + sel_anomaly <- dplyr::filter(x_anomaly, display_name == unique_groups[ii]) + + p1 <- + ggplot() + + geom_hline(data = sel_anomaly, + mapping = aes(yintercept = mean_var_y), + linetype = 1, + color = "grey50") + + geom_hline(data = sel_anomaly, + mapping = aes(yintercept = mean_var_y + sd_var_y), + linetype = 2, + color = "grey50") + + geom_hline(data = sel_anomaly, + mapping = aes(yintercept = mean_var_y - sd_var_y), + linetype = 2, + color = "grey50") + + geom_hline(data = sel_anomaly, + mapping = aes(yintercept = mean_var_y + 2*sd_var_y), + linetype = 3, + color = "grey50") + + geom_hline(data = sel_anomaly, + mapping = aes(yintercept = mean_var_y - 2*sd_var_y), + linetype = 3, + color = "grey50") + + geom_linerange(data = x |> + dplyr::filter(display_name == unique_groups[ii]), + mapping = aes(x = var_x, + ymax = var_y + 2*var_y_se, + ymin = var_y - 2*var_y_se), + color = "black") + + geom_rect(data = x |> + dplyr::filter(display_name == unique_groups[ii]), + mapping = aes(xmin = var_x+0.4, + xmax = var_x-0.4, + ymax = var_y + var_y_se, + ymin = var_y - var_y_se), + fill = fill_color) + + geom_segment(data = x |> + dplyr::filter(display_name == unique_groups[ii]), + aes(x = var_x+0.4, + xend = var_x-0.4, + y = var_y, + yend = var_y, + group = var_x), + color = "black") + + facet_grid(~display_name) + + scale_x_continuous(name = "Year", breaks = scales::pretty_breaks(n = 4)) + + scale_y_continuous(name = y_title, + trans = scales::trans_new("shift", + transform = function(x) {x-set_intercept}, + inverse = function(x) {x+set_intercept})) + + scale_fill_manual(name = ) + theme(legend.position = "none", + title = element_text(hjust = 0.5)) + + + if(write_plot) { + + png(here::here("plots", region_dir, spp_dir, + paste0(region, "_CONDITION_", spp_dir, ".png")), + width = 6, height = 4, units = "in", res = 600) + print(p1 + theme_blue_strip() + + theme(legend.position = "none", + title = element_text(hjust = 0.5))) + dev.off() + } + + out_list[ii] <- p1 + + } + + + return(out_list) + +} \ No newline at end of file diff --git a/R/plot_stratum_species_bar.R b/R/plot_stratum_species_bar.R new file mode 100644 index 0000000..00f20c3 --- /dev/null +++ b/R/plot_stratum_species_bar.R @@ -0,0 +1,126 @@ +#' Make single species stratum bar plots +#' +#' Make condition time series bar plots with separate panels for each stratum. +#' +#' @param x Input data.frame +#' @param region Region. AI, BS, or GOA. Character (1L). +#' @param var_x_name Name of the x (time variable). Character (1L). +#' @param var_y_name Name of the y variable in the data.frame. Character (1L). +#' @param var_y_se_name Name of the standard error for the y variable. Character (1L). +#' @param y_title Y-axis title. Character (1L). +#' @param fill_title Name of the fill variable to use for plotting. Character (1L). +#' @param fill_palette Character vector denoting which color palette to use. Must be a valid name for RColorBrewer::brewer.pal(name = {fill_palette}) +#' @param var_group_name Name of the variable to use for grouping (i.e., stratum). Character (1L). +#' @param write_plot Should plots be written to the /plot/ directory? +#' @return A list of bar plots as ggplot objects +#' @examples +#' # EBS single species stratum bar plots +#' +#' names(EBS_INDICATOR$STRATUM) +#' +#' plot_species_stratum_bar( +#' x = EBS_INDICATOR$STRATUM, +#' region = "BS", +#' var_x_name = "year", +#' var_y_name = "stratum_resid_mean", +#' var_y_se_name = "stratum_resid_se", +#' y_title = "Length-weight residual (ln(g))", +#' var_group_name = "stratum", +#' fill_title = "Stratum", +#' fill_palette = "BrBG", +#' write_plot = FALSE) +#' @export + +plot_species_stratum_bar <- function(x, + region, + var_x_name, + var_y_name, + var_y_se_name, + y_title = "Length-weight residual (ln(g))", + var_group_name = "stratum", + fill_title = "Stratum", + fill_palette = "BrBG", + write_plot = TRUE) { + + region <- toupper(region) + + stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'EBS', 'NBS', 'AI', or 'GOA'" = (region %in% c("EBS", "NBS", "AI", "GOA"))) + + if(region %in% c("EBS", "NBS")) { + region_dir <- "BS" + region_order <- "BS" + } else { + region_dir <- region + region_order <- region + } + + names(x)[which(names(x) == var_x_name)] <- "var_x" + names(x)[which(names(x) == var_y_name)] <- "var_y" + names(x)[which(names(x) == var_group_name)] <- "var_group" + names(x)[which(names(x) == var_y_se_name)] <- "var_y_se" + x <- dplyr::filter(x, !is.na(var_y)) + + x$display_name <- akfishcondition::set_plot_order(x$common_name, region = region_order) + + unique_groups <- unique(x$display_name) + + out_list <- list() + for(ii in 1:length(unique_groups)) { + + if(is.na(unique_groups[ii])) { + next + } + + spp_dir <- gsub(pattern = ">", replacement = "gt", unique_groups[ii]) + + if(!dir.exists(here::here("plots", region_dir, spp_dir))) { + dir.create(here::here("plots", region_dir, spp_dir), recursive = TRUE) + } + + p1 <- + ggplot(data = x |> + dplyr::filter(display_name == unique_groups[ii]), + aes(x = var_x, + y = var_y, + fill = set_stratum_order(trimws(var_group), region = region_order), + ymin = var_y - 2*var_y_se, + ymax = var_y + 2*var_y_se)) + + geom_hline(yintercept = 0) + + geom_linerange() + + geom_bar(stat = "identity", + color = "black", + position = "stack", + width = 1) + + facet_wrap(~set_stratum_order(trimws(var_group), + region = region_order), + ncol = 2, + scales = "free_y") + + ggtitle(unique_groups[ii]) + + scale_x_continuous(name = "Year", breaks = scales::pretty_breaks()) + + scale_y_continuous(name = y_title) + + scale_fill_brewer(name = fill_title, + palette = fill_palette, + drop = FALSE) + + theme(legend.position = "none", + title = element_text(hjust = 0.5)) + + + if(write_plot) { + + png(here::here("plots", region_dir, spp_dir, + paste0(region, "_STRATUM_CONDITION_", spp_dir, ".png")), + width = 6, height = 7, units = "in", res = 600) + print(p1 + theme_blue_strip() + + theme(legend.position = "none", + title = element_text(hjust = 0.5))) + dev.off() + } + + out_list[ii] <- p1 + + } + + + return(out_list) + +} \ No newline at end of file diff --git a/R/plot_stratum_stacked_bar.R b/R/plot_stratum_stacked_bar.R new file mode 100644 index 0000000..b086a14 --- /dev/null +++ b/R/plot_stratum_stacked_bar.R @@ -0,0 +1,68 @@ +#' Make condition time series stacked bar plot +#' +#' @param x Input data.frame +#' @param region Region. AI, BS, or GOA. Character (1L). +#' @param var_y_name Name of the y variable in the data.frame. Character (1L). +#' @param var_y_se_name Name of the standard error for the y variable. Character (1L). +#' @param var_x_name Name of the x (time variable). Character (1L). +#' @param var_group_name Name of the variable to use for grouping (i.e., stratum). Character (1L). +#' @param y_title Y-axis title. Character (1L). +#' @param fill_title Name of the fill variable to use for plotting. Character (1L). +#' @return A ggplot object of the time series +#' @examples +#' \dontrun{ +#' # EBS strata stacked bar plot +#' +#' names(EBS_INDICATOR$STRATUM) +#' +#' plot_stratum_stacked_bar(x = EBS_INDICATOR$STRATUM, +#' region = "BS", +#' fill_palette = "BrBG", +#' var_y_name = "stratum_resid_mean", +#' var_x_name = "year", +#' var_group_name = "stratum", +#' fill_title = "Stratum", +#' y_title = "Length-weight residual (ln(g))") +#' } +#' @export + +plot_stratum_stacked_bar <- function(x, + region, + fill_palette = "BrBG", + var_y_name = "mean_wt_resid", + var_x_name = "year", + var_group_name = "stratum", + fill_title = "Stratum", + y_title = "Length-weight residual (ln(g))") { + + region <- toupper(region) + + stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'BS', 'AI', or 'GOA'" = (region %in% c("BS", "AI", "GOA"))) + + x$display_name <- akfishcondition::set_plot_order(x$common_name, + region = region) + + names(x)[which(names(x) == var_x_name)] <- "var_x" + names(x)[which(names(x) == var_y_name)] <- "var_y" + names(x)[which(names(x) == var_group_name)] <- "var_group" + + x <- dplyr::filter(x, !is.na(var_y)) + + p1 <- ggplot(data = x, + aes(x = var_x, + y = var_y, + fill = set_stratum_order(trimws(var_group), + region = region))) + + geom_hline(yintercept = 0) + + geom_bar(stat = "identity", + color = "black", + position = "stack", + width = 1) + + facet_wrap(~display_name, ncol = 2, scales = "free_y") + + scale_x_continuous(name = "Year") + + scale_y_continuous(name = y_title) + + scale_fill_brewer(name = fill_title, palette = fill_palette) + + return(p1) + +} \ No newline at end of file diff --git a/R/plot_two_timeseries.R b/R/plot_two_timeseries.R new file mode 100644 index 0000000..9dafbc6 --- /dev/null +++ b/R/plot_two_timeseries.R @@ -0,0 +1,165 @@ +#' Two condition time series and correlation +#' +#' @param x_1 A data.frame for the first time series. +#' @param x_2 A data.frame for the second time series. +#' @param region Region. AI, BS, or GOA. Character (1L). +#' @param fill_colors Fill colors to use for points. Character (2L). +#' @param var_y_name_1 Name of the y variable in the data.frame. Character (1L). +#' @param var_x_name_1 Name of the x (time variable). Character (1L). +#' @param var_y_name_2 Name of the y variable in the data.frame. Character (1L). +#' @param var_x_name_2 Name of the x (time variable). Character (1L). +#' @param scale_y Logical. Should y variables be Z-score transformed? +#' @param year_break Optional. Year to split a times series for geom_line(); denotes a break in the time series (i.e., EBS and NBS in 2020) +#' @param y_title Y-axis title. Character (1L). +#' @param x_offset Offset value to add to an x variable (e.g. year) from one of the time series to improve interpretability. Numeric (1L). +#' @param fill_title Name of the fill variable to use for plotting. Character (1L). +#' @param fill_color Fill colors to use for points. +#' @param shapes Shapes to use for points. +#' @param format_for "rmd" or "png" +#' @return A ggplot object of the time series +#' @export + +plot_two_timeseries <- function(x_1, + x_2, + region, + series_name_1 = "Historical", + var_y_name_1 = "mean_wt_resid", + var_x_name_1 = "year", + series_name_2 = "VAST", + var_y_name_2 = "mean_wt_resid", + var_x_name_2 = "year", + var_group_name = "common_name", + y_title = "Condition (Z-score)", + fill_title = "Method", + scale_y = TRUE, + year_break = NULL, + x_offset = 0.25, + fill_colors = c("#001743", "#0085CA"), + shapes = c(24, 21), + format_for = "rmd") { + + region <- toupper(region) + + point_rel_size <- c(5, 2.5)[match(tolower(format_for), c("rmd", "png"))] + + stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'BS', 'AI', or 'GOA'" = (region %in% c("BS", "AI", "GOA"))) + + shared_cols <- c("var_x", "var_y", "var_y_se", "display_name", "common_name", "series") + + # Setup time series 1 + + names(x_1)[which(names(x_1) == var_group_name)] <- "var_group" + + if(var_group_name == "common_name") { + x_1$display_name <- akfishcondition::set_plot_order(x_1$var_group, + region = region) + } else { + x_1$display_name <- x_1$var_group + } + + names(x_1)[which(names(x_1) == var_x_name_1)] <- "var_x" + names(x_1)[which(names(x_1) == var_y_name_1)] <- "var_y" + x_1$series <- series_name_1 + x_1 <- x_1[, which(names(x_1) %in% shared_cols)] + x_1 <- dplyr::filter(x_1, !is.na(var_y)) + + # Setup time series 2 + names(x_2)[which(names(x_2) == var_group_name)] <- "var_group" + + if(var_group_name == "common_name") { + x_2$display_name <- akfishcondition::set_plot_order(x_2$var_group, + region = region) + } else { + x_2$display_name <- x_2$var_group + } + + names(x_2)[which(names(x_2) == var_x_name_2)] <- "var_x" + names(x_2)[which(names(x_2) == var_y_name_2)] <- "var_y" + x_2$series <- series_name_2 + x_2 <- x_2[, which(names(x_2) %in% shared_cols)] + x_2 <- dplyr::filter(x_2, !is.na(var_y)) + + # Scale y variables + if(scale_y) { + x_1 <- x_1 |> + dplyr::group_by(common_name, display_name) |> + dplyr::mutate(var_y = scale(var_y)[,1]) + x_2 <- x_2 |> + dplyr::group_by(common_name, display_name) |> + dplyr::mutate(var_y = scale(var_y)[,1]) + } + + + # Correlation between timeseries + corr_df <- dplyr::inner_join( + x_1 |> + dplyr::ungroup() |> + dplyr::rename(var_y_1 = var_y) |> + dplyr::select(display_name, + var_x, + var_y_1), + x_2 |> + dplyr::ungroup() |> + dplyr::rename(var_y_2 = var_y) |> + dplyr::select(display_name, + var_x, + var_y_2), + by = c("display_name", "var_x")) |> + dplyr::group_by(display_name) |> + dplyr::summarise( + r = round(cor(var_y_1, + var_y_2, + use = "complete.obs", + method = "pearson"), 2)) + + x_2$var_x <- x_2$var_x + x_offset + x_combined <- dplyr::bind_rows(x_1, x_2) + + if(!is.null(year_break)) { + x_combined$grp <- x_combined$var_x < year_break + } else { + x_combined$grp <- 1 + } + + max_x <- max(x_combined$var_x) + min_x <- min(x_combined$var_x) + lab_x <- max_x - (max_x - min_x)*0.12 + + min_y <- min(x_combined$var_y) + max_y <- max(x_combined$var_y) + lab_y <- min_y + (max_y - min_y)*0.07 + + # Make plot + p1 <- ggplot() + + geom_hline(yintercept = 0, + linetype = 1, + color = "grey50") + + geom_hline(yintercept = c(-1,1), + linetype = 2, + color = "grey50") + + geom_hline(yintercept = c(-2,2), + linetype = 3, + color = "grey50") + + geom_point(data = x_combined, + aes(x = var_x, + y = var_y, + fill = series, + shape = series), + color = "black", + size = rel(point_rel_size)) + + geom_text(data = corr_df, + aes(x = lab_x, + y = lab_y, + label = paste0("r = ", format(r, nsmall = 2))), + size = rel(point_rel_size*1.9)) + + scale_fill_manual(name = fill_title, + values = fill_colors) + + scale_shape_manual(name = fill_title, + values = shapes) + + facet_wrap(~display_name, ncol = 2) + + scale_x_continuous(name = "Year", breaks = scales::pretty_breaks(n = 4)) + + scale_y_continuous(name = y_title, limits = c(range(x_combined$var_y))) + + return(p1) + +} \ No newline at end of file diff --git a/R/plot_xy_corr.R b/R/plot_xy_corr.R new file mode 100644 index 0000000..d0fc842 --- /dev/null +++ b/R/plot_xy_corr.R @@ -0,0 +1,130 @@ +#' Plot of condition versus condition +#' +#' @param x_1 A data.frame for the first time series. +#' @param x_2 A data.frame for the second time series. +#' @param region Region. AI, BS, or GOA. Character (1L). +#' @param fill_colors Fill colors to use for points. Character (2L). +#' @param var_y_name_1 Name of the y variable in the data.frame. Character (1L). +#' @param var_y_se_name_1 Name of the standard error for the y variable. Character (1L). +#' @param var_x_name_1 Name of the x (time variable). Character (1L). +#' @param var_y_name_2 Name of the y variable in the data.frame. Character (1L). +#' @param var_y_se_name_2 Name of the standard error for the y variable. Character (1L). +#' @param var_x_name_2 Name of the x (time variable). Character (1L). +#' @param scale_y Logical. Should y variables be Z-score transformed? +#' @param y_title Y-axis title. Character (1L). +#' @param format_for "rmd" or "png" +#' @return A ggplot object of the time series +#' @export + +plot_xy_corr <- function(x_1, + x_2, + region, + series_name_1 = "Historical", + var_y_name_1 = "mean_wt_resid", + var_y_se_name_1 = "se_wt_resid", + var_x_name_1 = "year", + series_name_2 = "VAST", + var_y_name_2 = "mean_wt_resid", + var_y_se_name_2 = "se_wt_resid", + var_x_name_2 = "year", + format_for = "rmd") { + + region <- toupper(region) + + point_rel_size <- c(5, 2.5)[match(tolower(format_for), c("rmd", "png"))] + + stopifnot("Invalid region in plot_anomaly_timeseries. Must be 'BS', 'AI', or 'GOA'" = (region %in% c("BS", "AI", "GOA"))) + + shared_cols <- c("var_x", "var_y", "var_y_se", "display_name", "common_name", "series") + + # Setup time series 1 + x_1$display_name <- akfishcondition::set_plot_order(x_1$common_name, + region = region) + names(x_1)[which(names(x_1) == var_x_name_1)] <- "var_x" + names(x_1)[which(names(x_1) == var_y_name_1)] <- "var_y" + names(x_1)[which(names(x_1) == var_y_se_name_1)] <- "var_y_se" + x_1$series <- series_name_1 + x_1 <- x_1[, which(names(x_1) %in% shared_cols)] + x_1 <- dplyr::filter(x_1, !is.na(var_y)) + + # Setup time series 2 + x_2$display_name <- akfishcondition::set_plot_order(x_2$common_name, + region = region) + names(x_2)[which(names(x_2) == var_x_name_2)] <- "var_x" + names(x_2)[which(names(x_2) == var_y_name_2)] <- "var_y" + names(x_2)[which(names(x_2) == var_y_se_name_2)] <- "var_y_se" + x_2$series <- series_name_2 + x_2 <- x_2[, which(names(x_2) %in% shared_cols)] + x_2 <- dplyr::filter(x_2, !is.na(var_y)) + + # Correlation between timeseries + x_combined <- dplyr::inner_join( + x_1 |> + dplyr::ungroup() |> + dplyr::rename(var_y_1 = var_y, + var_y_se_1 = var_y_se) |> + dplyr::select(display_name, + var_x, + var_y_1, + var_y_se_1), + x_2 |> + dplyr::ungroup() |> + dplyr::rename(var_y_2 = var_y, + var_y_se_2 = var_y_se) |> + dplyr::select(display_name, + var_x, + var_y_2, + var_y_se_2), + by = c("display_name", "var_x")) |> + dplyr::ungroup() + + + corr_df <- x_combined |> + dplyr::group_by(display_name) |> + dplyr::summarise( + r = round(cor(var_y_1, + var_y_2, + use = "complete.obs", + method = "pearson"), 2)) + + corr_df <- corr_df |> + dplyr::inner_join( + x_combined |> + dplyr::group_by(display_name) |> + dplyr::summarise(min_x = min(var_y_1-2*var_y_se_1), + max_x = max(var_y_1+2*var_y_se_1), + min_y = min(var_y_2-2*var_y_se_2), + max_y = max(var_y_2+2*var_y_se_2)) |> + dplyr::mutate(lab_x = max_x - (max_x - min_x)*0.12, + lab_y = min_y + (max_y - min_y)*0.07), + by = "display_name") + + # Make plot + p1 <- ggplot() + + geom_linerange(data = x_combined, + aes(x = var_y_1, + ymin = var_y_2 - 2 * var_y_se_2, + ymax = var_y_2 + 2 * var_y_se_2), + color = "black") + + geom_linerange(data = x_combined, + aes(y = var_y_2, + xmin = var_y_1 - 2 * var_y_se_1, + xmax = var_y_1 + 2 * var_y_se_1), + color = "black") + + geom_point(data = x_combined, + aes(x = var_y_1, + y = var_y_2), + color = "black", + size = rel(point_rel_size)) + + geom_text(data = corr_df, + aes(x = lab_x, + y = lab_y, + label = paste0("r = ", format(r, nsmall = 2))), + size = rel(point_rel_size*1.9)) + + facet_wrap(~display_name, ncol = 2, scale = "free") + + scale_x_continuous(name = series_name_1, breaks = scales::pretty_breaks()) + + scale_y_continuous(name = series_name_2) + + return(p1) + +} \ No newline at end of file diff --git a/R/run_sbw_condition.R b/R/run_sbw_condition.R index 25a0dcd..c150e0a 100644 --- a/R/run_sbw_condition.R +++ b/R/run_sbw_condition.R @@ -22,7 +22,7 @@ run_sbw_condition <- function(region, stratum_col = NULL, biomass_col = NULL, co } if(is.null(biomass_col)) { - biomass_col <- c("area_biomass", "biomass", "biomass", "area_biomass")[region_index] + biomass_col <- "biomass" } biomass <- read.csv(file = here::here("data", paste0(region, "_stratum_biomass_all_species.csv"))) @@ -32,8 +32,8 @@ run_sbw_condition <- function(region, stratum_col = NULL, biomass_col = NULL, co lw <- read.csv(file = here::here("data", paste0(region, "_all_species.csv"))) if(region == "EBS") { - lw <- lw %>% - dplyr::filter(stratum < 70) %>% + lw <- lw |> + dplyr::filter(stratum < 70) |> dplyr::mutate(stratum = floor(stratum/10)*10) } @@ -41,27 +41,27 @@ run_sbw_condition <- function(region, stratum_col = NULL, biomass_col = NULL, co lw <- dplyr::mutate(lw, stratum = 999) } - lw <- lw %>% - dplyr::filter(!is.na(length_mm)) %>% + lw <- lw |> + dplyr::filter(!is.na(length_mm)) |> dplyr::mutate(ID = paste(cruise, vessel, haul, specimenid, sep = "_")) names(lw)[which(names(lw) == stratum_col)] <- "survey_stratum" - nsamp <- table(lw$survey_stratum, lw$common_name) %>% + nsamp <- table(lw$survey_stratum, lw$common_name) |> as.data.frame() |> dplyr::rename(n_raw = Freq) - dat <- lw %>% + dat <- lw |> dplyr::inner_join(biomass, - by = c("year", "species_code", "survey_stratum")) %>% + by = c("year", "species_code", "survey_stratum")) |> dplyr::select(-vessel, -cruise, -haul, -hauljoin) # Checks (stratum count match should be TRUE) - qc_df <- table(dat$survey_stratum, dat$common_name) %>% - as.data.frame() %>% - dplyr::rename(n_filtered = Freq) %>% - dplyr::inner_join(nsamp) %>% + qc_df <- table(dat$survey_stratum, dat$common_name) |> + as.data.frame() |> + dplyr::rename(n_filtered = Freq) |> + dplyr::inner_join(nsamp) |> dplyr::mutate(equal = n_filtered == n_raw) count_match <- sum(qc_df$equal) == nrow(qc_df) @@ -76,14 +76,14 @@ run_sbw_condition <- function(region, stratum_col = NULL, biomass_col = NULL, co message("Splitting Pacific cod and pollock for ESP and ESR length cutoffs.") # Create separate data.frames for adult and juvenile pollock and cod - pcod <- dat %>% dplyr::filter(species_code == 21720) + pcod <- dat |> dplyr::filter(species_code == 21720) pcod$species_code[pcod$length < cod_juv_cutoff_mm] <- 21721 pcod$species_code[pcod$length >= cod_juv_cutoff_mm] <- 21722 pcod$common_name[pcod$species_code == 21721] <- "Pacific cod (juvenile)" pcod$common_name[pcod$species_code == 21722] <- "Pacific cod (adult)" pollock <- dplyr::filter(dat, species_code == 21740) dat$species_code[dat$species_code == 21740 & dat$length_mm >= 100 & dat$length_mm <= 250] <- 21741 - dat$common_name[dat$species_code == 21741] <- "walleye pollock (100–250 mm)" + dat$common_name[dat$species_code == 21741] <- "walleye pollock (100-250 mm)" dat$species_code[dat$species_code == 21740] <- 21742 dat$common_name[dat$species_code == 21742] <- "walleye pollock (>250 mm)" @@ -134,22 +134,22 @@ run_sbw_condition <- function(region, stratum_col = NULL, biomass_col = NULL, co stratum_resids <- NULL - ann_mean_resid_df <- dat %>% - dplyr::group_by(common_name, year) %>% + ann_mean_resid_df <- dat |> + dplyr::group_by(common_name, year) |> dplyr::summarise(mean_resid = mean(resid, na.rm = TRUE), se = sd(resid, na.rm = TRUE)/n(), - n = n()) %>% + n = n()) |> dplyr::filter(n >= min_n) } else { # Estimate mean and std. err for each stratum, filter out strata with less than min_n samples - stratum_resids <- dat %>% - dplyr::group_by(common_name, species_code, year, survey_stratum, stratum_biomass) %>% + stratum_resids <- dat |> + dplyr::group_by(common_name, species_code, year, survey_stratum, stratum_biomass) |> dplyr::summarise(stratum_resid_mean = mean(resid_mean), stratum_resid_sd = sd(resid_mean), - n = n()) %>% - dplyr::filter(n >= min_n) %>% + n = n()) |> + dplyr::filter(n >= min_n) |> dplyr::mutate(stratum_resid_se = stratum_resid_sd/sqrt(n)) # Weight strata by biomass @@ -167,12 +167,12 @@ run_sbw_condition <- function(region, stratum_col = NULL, biomass_col = NULL, co } # Biomass-weighted residual and SE by year - ann_mean_resid_df <- stratum_resids %>% - dplyr::group_by(year, common_name) %>% + ann_mean_resid_df <- stratum_resids |> + dplyr::group_by(year, common_name) |> dplyr::summarise(mean_wt_resid = mean(weighted_resid_mean), se_wt_resid = mean(weighted_resid_se)) - stratum_resids <- stratum_resids %>% + stratum_resids <- stratum_resids |> dplyr::select(-stratum_biomass, - stratum_resid_sd) names(stratum_resids)[which(names(stratum_resids) == "survey_stratum")] <- stratum_col @@ -183,6 +183,8 @@ run_sbw_condition <- function(region, stratum_col = NULL, biomass_col = NULL, co names(lw)[which(names(lw) == "survey_stratum")] <- stratum_col names(biomass)[which(names(biomass) == "stratum_biomass")] <- biomass_col + stratum_resids <- stratum_resids |> dplyr::ungroup() |> dplyr::select(-species_code) + return(list(full_sbw = ann_mean_resid_df, stratum_sbw = stratum_resids, input_data = list(LW = lw, diff --git a/R/run_sbw_condition_multimodel.R b/R/run_sbw_condition_multimodel.R index d58fb5a..e56c946 100644 --- a/R/run_sbw_condition_multimodel.R +++ b/R/run_sbw_condition_multimodel.R @@ -46,7 +46,7 @@ run_sbw_condition_multimodel <- function(region, } if(is.null(biomass_col)) { - biomass_col <- c("area_biomass", "biomass", "biomass", "area_biomass")[region_index] + biomass_col <- "biomass" } biomass <- read.csv(file = here::here("data", paste0(region, "_stratum_biomass_all_species.csv"))) @@ -59,36 +59,32 @@ run_sbw_condition_multimodel <- function(region, lw$day_of_year <- lubridate::yday(lw$start_time) if(region == "EBS") { - lw <- lw %>% - dplyr::filter(stratum < 70) %>% + lw <- lw |> + dplyr::filter(stratum < 70) |> dplyr::mutate(stratum = floor(stratum/10)*10) } - # if(region == "NBS") { - # lw <- dplyr::mutate(lw, stratum = 999) - # } - - lw <- lw %>% - dplyr::filter(!is.na(length_mm)) %>% + lw <- lw |> + dplyr::filter(!is.na(length_mm)) |> dplyr::mutate(ID = paste(cruise, vessel, haul, specimenid, sep = "_")) names(lw)[which(names(lw) == stratum_col)] <- "survey_stratum" - nsamp <- table(lw$survey_stratum, lw$common_name) %>% + nsamp <- table(lw$survey_stratum, lw$common_name) |> as.data.frame() |> dplyr::rename(n_raw = Freq) - dat <- lw %>% + dat <- lw |> dplyr::inner_join(biomass, - by = c("year", "species_code", "survey_stratum")) %>% + by = c("year", "species_code", "survey_stratum")) |> dplyr::select(-vessel, -cruise, -haul, -hauljoin) # Checks (stratum count match should be TRUE) - qc_df <- table(dat$survey_stratum, dat$common_name) %>% - as.data.frame() %>% - dplyr::rename(n_filtered = Freq) %>% - dplyr::inner_join(nsamp) %>% + qc_df <- table(dat$survey_stratum, dat$common_name) |> + as.data.frame() |> + dplyr::rename(n_filtered = Freq) |> + dplyr::inner_join(nsamp) |> dplyr::mutate(equal = n_filtered == n_raw) count_match <- sum(qc_df$equal) == nrow(qc_df) @@ -103,14 +99,14 @@ run_sbw_condition_multimodel <- function(region, message("Splitting Pacific cod and pollock for ESP and ESR length cutoffs.") # Create separate data.frames for adult and juvenile pollock and cod - pcod <- dat %>% dplyr::filter(species_code == 21720) + pcod <- dat |> dplyr::filter(species_code == 21720) pcod$species_code[pcod$length < cod_juv_cutoff_mm] <- 21721 pcod$species_code[pcod$length >= cod_juv_cutoff_mm] <- 21722 pcod$common_name[pcod$species_code == 21721] <- "Pacific cod (juvenile)" pcod$common_name[pcod$species_code == 21722] <- "Pacific cod (adult)" pollock <- dplyr::filter(dat, species_code == 21740) dat$species_code[dat$species_code == 21740 & dat$length_mm >= 100 & dat$length_mm <= 250] <- 21741 - dat$common_name[dat$species_code == 21741] <- "walleye pollock (100–250 mm)" + dat$common_name[dat$species_code == 21741] <- "walleye pollock (100-250 mm)" dat$species_code[dat$species_code == 21740] <- 21742 dat$common_name[dat$species_code == 21742] <- "walleye pollock (>250 mm)" @@ -123,97 +119,67 @@ run_sbw_condition_multimodel <- function(region, # Calculate length weight residuals for(i in 1:length(spp_vec)) { - # if(region == "NBS") { - # dat$resid[dat$species_code == spp_vec[i]] <- akfishcondition::calc_lw_residuals_multimodel( - # len = dat$length_mm[dat$species_code == spp_vec[i]], - # wt = dat$weight_g[dat$species_code == spp_vec[i]], - # sex = null_flag(use = use_sex, var = dat$sex[dat$species_code == spp_vec[i]]), - # year = dat$year[dat$species_code == spp_vec[i]], - # day_of_year = null_flag(use = use_doy, var = dat$day_of_year[dat$species_code == spp_vec[i]]), - # stratum = null_flag(use = use_stratum, var = dat$survey_stratum[dat$species_code == spp_vec[i]]), - # make_diagnostics = TRUE, # Make diagnostics - # include_ci = FALSE, - # bias_correction = TRUE, # Bias correction turned on - # outlier_rm = TRUE, # Outlier removal turned off - # region = region, - # species_code = dat$species_code[dat$species_code == spp_vec[i]] - # ) - # - # } else { - # Separate slope for each stratum. Bias correction according to Brodziak, no outlier detection. - raw_resid_df <- akfishcondition::calc_lw_residuals_multimodel( - len = dat$length_mm[dat$species_code == spp_vec[i]], - wt = dat$weight_g[dat$species_code == spp_vec[i]], - sex = null_flag(use = use_sex, var = dat$sex[dat$species_code == spp_vec[i]]), - year = dat$year[dat$species_code == spp_vec[i]], - day_of_year = null_flag(use = use_doy, var = dat$day_of_year[dat$species_code == spp_vec[i]]), - stratum = null_flag(use = use_stratum, var = dat$survey_stratum[dat$species_code == spp_vec[i]]), - make_diagnostics = TRUE, # Make diagnostics - bias_correction = TRUE, # Bias correction turned on - outlier_rm = TRUE, # Outlier removal turned off - region = region, - species_code = dat$species_code[dat$species_code == spp_vec[i]] - ) - - dat$resid_mean[dat$species_code == spp_vec[i]] <- raw_resid_df$lw.res_mean - dat$resid_lwr[dat$species_code == spp_vec[i]] <- raw_resid_df$lw.res_lwr - dat$resid_upr[dat$species_code == spp_vec[i]] <- raw_resid_df$lw.res_upr - # } - } - - # if(region == "NBS") { - # - # stratum_resids <- NULL - # - # ann_mean_resid_df <- dat %>% - # dplyr::group_by(common_name, year) %>% - # dplyr::summarise(mean_resid = mean(resid, na.rm = TRUE), - # se = sd(resid, na.rm = TRUE)/n(), - # n = n()) %>% - # dplyr::filter(n >= min_n) - # - # } else { - - # Estimate mean and std. err for each stratum, filter out strata with less than min_n samples - stratum_resids <- dat %>% - dplyr::group_by(common_name, species_code, year, survey_stratum, stratum_biomass) %>% - dplyr::summarise(stratum_resid_mean = mean(resid_mean), - stratum_resid_sd = sd(resid_mean), - n = n()) %>% - dplyr::filter(n >= min_n) %>% - dplyr::mutate(stratum_resid_se = stratum_resid_sd/sqrt(n)) + # Separate slope for each stratum. Bias correction according to Brodziak, no outlier detection. + raw_resid_df <- akfishcondition::calc_lw_residuals_multimodel( + len = dat$length_mm[dat$species_code == spp_vec[i]], + wt = dat$weight_g[dat$species_code == spp_vec[i]], + sex = null_flag(use = use_sex, var = dat$sex[dat$species_code == spp_vec[i]]), + year = dat$year[dat$species_code == spp_vec[i]], + day_of_year = null_flag(use = use_doy, var = dat$day_of_year[dat$species_code == spp_vec[i]]), + stratum = null_flag(use = use_stratum, var = dat$survey_stratum[dat$species_code == spp_vec[i]]), + make_diagnostics = TRUE, # Make diagnostics + bias_correction = TRUE, # Bias correction turned on + outlier_rm = TRUE, # Outlier removal turned off + region = region, + species_code = dat$species_code[dat$species_code == spp_vec[i]] + ) - # Weight strata by biomass - for(i in 1:length(spp_vec)) { - stratum_resids$weighted_resid_mean[stratum_resids$species_code == spp_vec[i]] <- - weight_lw_residuals(residuals = stratum_resids$stratum_resid_mean[stratum_resids$species_code == spp_vec[i]], - year = stratum_resids$year[stratum_resids$species_code == spp_vec[i]], - stratum = stratum_resids$survey_stratum[stratum_resids$species_code == spp_vec[i]], - stratum_biomass = stratum_resids$stratum_biomass[stratum_resids$species_code == spp_vec[i]]) - stratum_resids$weighted_resid_se[stratum_resids$species_code == spp_vec[i]] <- - weight_lw_residuals(residuals = stratum_resids$stratum_resid_se[stratum_resids$species_code == spp_vec[i]], - year = stratum_resids$year[stratum_resids$species_code == spp_vec[i]], - stratum = stratum_resids$survey_stratum[stratum_resids$species_code == spp_vec[i]], - stratum_biomass = stratum_resids$stratum_biomass[stratum_resids$species_code == spp_vec[i]]) - } + dat$resid_mean[dat$species_code == spp_vec[i]] <- raw_resid_df$lw.res_mean + dat$resid_lwr[dat$species_code == spp_vec[i]] <- raw_resid_df$lw.res_lwr + dat$resid_upr[dat$species_code == spp_vec[i]] <- raw_resid_df$lw.res_upr - # Biomass-weighted residual and SE by year - ann_mean_resid_df <- stratum_resids %>% - dplyr::group_by(year, common_name) %>% - dplyr::summarise(mean_wt_resid = mean(weighted_resid_mean), - se_wt_resid = mean(weighted_resid_se)) - - stratum_resids <- stratum_resids %>% - dplyr::select(-stratum_biomass, - stratum_resid_sd) - - names(stratum_resids)[which(names(stratum_resids) == "survey_stratum")] <- stratum_col - - # } + } + + # Estimate mean and std. err for each stratum, filter out strata with less than min_n samples + stratum_resids <- dat |> + dplyr::group_by(common_name, year, species_code, survey_stratum, stratum_biomass) |> + dplyr::summarise(stratum_resid_mean = mean(resid_mean), + stratum_resid_sd = sd(resid_mean), + n = n()) |> + dplyr::filter(n >= min_n) |> + dplyr::mutate(stratum_resid_se = stratum_resid_sd/sqrt(n)) + + # Weight strata by biomass + for(i in 1:length(spp_vec)) { + stratum_resids$weighted_resid_mean[stratum_resids$species_code == spp_vec[i]] <- + weight_lw_residuals(residuals = stratum_resids$stratum_resid_mean[stratum_resids$species_code == spp_vec[i]], + year = stratum_resids$year[stratum_resids$species_code == spp_vec[i]], + stratum = stratum_resids$survey_stratum[stratum_resids$species_code == spp_vec[i]], + stratum_biomass = stratum_resids$stratum_biomass[stratum_resids$species_code == spp_vec[i]]) + stratum_resids$weighted_resid_se[stratum_resids$species_code == spp_vec[i]] <- + weight_lw_residuals(residuals = stratum_resids$stratum_resid_se[stratum_resids$species_code == spp_vec[i]], + year = stratum_resids$year[stratum_resids$species_code == spp_vec[i]], + stratum = stratum_resids$survey_stratum[stratum_resids$species_code == spp_vec[i]], + stratum_biomass = stratum_resids$stratum_biomass[stratum_resids$species_code == spp_vec[i]]) + } + + # Biomass-weighted residual and SE by year + ann_mean_resid_df <- stratum_resids |> + dplyr::group_by(year, common_name) |> + dplyr::summarise(mean_wt_resid = mean(weighted_resid_mean), + se_wt_resid = mean(weighted_resid_se)) + + stratum_resids <- stratum_resids |> + dplyr::select(-stratum_biomass, - stratum_resid_sd) + + names(stratum_resids)[which(names(stratum_resids) == "survey_stratum")] <- stratum_col names(biomass)[which(names(biomass) == "survey_stratum")] <- stratum_col names(lw)[which(names(lw) == "survey_stratum")] <- stratum_col names(biomass)[which(names(biomass) == "stratum_biomass")] <- biomass_col + stratum_resids <- stratum_resids |> dplyr::ungroup() |> dplyr::select(-species_code) + return(list(full_sbw = ann_mean_resid_df, stratum_sbw = stratum_resids, input_data = list(LW = lw, diff --git a/R/sysdata.rda b/R/sysdata.rda index 018d17c..6daf4c3 100644 Binary files a/R/sysdata.rda and b/R/sysdata.rda differ diff --git a/R/themes.R b/R/themes.R index db4714c..c3e7fac 100644 --- a/R/themes.R +++ b/R/themes.R @@ -13,7 +13,7 @@ theme_condition_index <- function() { legend.text = element_text(size = 24), legend.position = "right", panel.grid.minor = element_blank(), - strip.text = element_text(size = 24, + strip.text = element_text(size = 20, color = "white", face = "bold", margin = margin(0.5, 0, 0.5, 0, "mm")), @@ -78,7 +78,7 @@ set_plot_order <- function(common_name, region) { return(factor(common_name, levels = c("walleye pollock", "walleye pollock (>250 mm)", - "walleye pollock (100–250 mm)", + "walleye pollock (100-250 mm)", "Pacific cod", "Pacific cod (juvenile)", "Pacific cod (adult)", @@ -89,7 +89,7 @@ set_plot_order <- function(common_name, region) { "Pacific ocean perch"), labels = c("walleye pollock", "walleye pollock (>250 mm)", - "walleye pollock (100–250 mm)", + "walleye pollock (100-250 mm)", "Pacific cod", "Pacific cod (juvenile)", "Pacific cod (adult)", @@ -102,7 +102,7 @@ set_plot_order <- function(common_name, region) { return(factor(common_name, levels = c("walleye pollock", "walleye pollock (>250 mm)", - "walleye pollock (100–250 mm)", + "walleye pollock (100-250 mm)", "Pacific cod", "Pacific cod (juvenile)", "Pacific cod (adult)", diff --git a/R/utils.R b/R/utils.R index 33bc999..43c9645 100644 --- a/R/utils.R +++ b/R/utils.R @@ -87,7 +87,7 @@ get_condition_data <- function(channel = NULL) { names(dat_lw) <- casefold(names(dat_lw)) names(dat_cpue) <- casefold(names(dat_cpue)) names(dat_biomass) <- casefold(names(dat_biomass)) - combined_cpue_lw <- dplyr::bind_rows(dat_lw, dat_cpue) #, pollock_cpue, pollock_lw, cod_cpue, cod_lw) + combined_cpue_lw <- dplyr::bind_rows(dat_lw, dat_cpue) write.csv(combined_cpue_lw, paste0(getwd(),"/data/",i, "_all_species.csv"), row.names = FALSE) write.csv(dat_biomass, paste0(getwd(),"/data/",i, "_stratum_biomass_all_species.csv"), row.names = FALSE) @@ -182,39 +182,39 @@ select_species <- function(species_code, region, remove_outliers = TRUE, bonferr #' @noRd make_data_summary <- function(dat_csv, region) { - dat <- read.csv(file = dat_csv) %>% - dplyr::mutate(start_time = as.POSIXct(start_time)) %>% - dplyr::mutate(yday = lubridate::yday(start_time)) %>% + dat <- read.csv(file = dat_csv) |> + dplyr::mutate(start_time = as.POSIXct(start_time)) |> + dplyr::mutate(yday = lubridate::yday(start_time)) |> dplyr::arrange(year) - n_spp_by_year <- dat %>% - dplyr::filter(!is.na(length_mm)) %>% - dplyr::group_by(common_name, year) %>% + n_spp_by_year <- dat |> + dplyr::filter(!is.na(length_mm)) |> + dplyr::group_by(common_name, year) |> dplyr::summarise(n = n()) - out <- list(n_cpue_by_year = dat %>% - dplyr::filter(!is.na(cpue_kg_km2)) %>% - dplyr::group_by(common_name, year) %>% - dplyr::summarise(n = n()) %>% - tidyr::pivot_wider(names_from = c("common_name"), values_from = "n", values_fill = 0) %>% + out <- list(n_cpue_by_year = dat |> + dplyr::filter(!is.na(cpue_kg_km2)) |> + dplyr::group_by(common_name, year) |> + dplyr::summarise(n = n()) |> + tidyr::pivot_wider(names_from = c("common_name"), values_from = "n", values_fill = 0) |> data.frame(), - n_specimen_by_year = n_spp_by_year%>% - tidyr::pivot_wider(names_from = c("common_name"), values_from = "n", values_fill = 0) %>% - dplyr::arrange(year) %>% + n_specimen_by_year = n_spp_by_year|> + tidyr::pivot_wider(names_from = c("common_name"), values_from = "n", values_fill = 0) |> + dplyr::arrange(year) |> data.frame(), - n_specimen_by_vessel =dat %>% - dplyr::filter(!is.na(length_mm)) %>% - dplyr::group_by(vessel, cruise) %>% - dplyr::summarise(n = n()) %>% - tidyr::pivot_wider(names_from = c("vessel"), values_from = "n", values_fill = 0) %>% - dplyr::arrange(cruise) %>% + n_specimen_by_vessel =dat |> + dplyr::filter(!is.na(length_mm)) |> + dplyr::group_by(vessel, cruise) |> + dplyr::summarise(n = n()) |> + tidyr::pivot_wider(names_from = c("vessel"), values_from = "n", values_fill = 0) |> + dplyr::arrange(cruise) |> data.frame()) dir.create(here::here("output", region), recursive = TRUE) saveRDS(object = out, file = here::here("output", region, paste0(region, "_sample_tables.rds"))) - yday_df <- dat %>% + yday_df <- dat |> dplyr::filter(!is.na(length_mm)) yday_range <- range(yday_df$yday) @@ -235,7 +235,7 @@ make_data_summary <- function(dat_csv, region) { for(jj in 1:ceiling(length(yday_years)/6)) { png(file = here::here("output", region, paste0(region, "_ecdf_samples_by_year_", jj, ".png")), width = 7, height = 7, units = "in", res = 120) - print(ggplot(data = yday_df %>% + print(ggplot(data = yday_df |> dplyr::filter(year %in% yday_years[(1+6*(jj-1)):min(c((6+6*(jj-1)), length(yday_years)))]), aes(x = yday, color = common_name)) + stat_ecdf(size = rel(1.2)) + @@ -260,10 +260,10 @@ make_data_summary <- function(dat_csv, region) { map_layers <- akgfmaps::get_base_layers(select.region = tolower(map_region), set.crs = "EPSG:3338") - lw_dat.sub <- read.csv(file = here::here("data", paste0(region, "_all_species.csv"))) %>% - dplyr::filter(!is.na(length_mm)) %>% - dplyr::select(year, species_code, latitude, longitude, common_name) %>% - unique() %>% + lw_dat.sub <- read.csv(file = here::here("data", paste0(region, "_all_species.csv"))) |> + dplyr::filter(!is.na(length_mm)) |> + dplyr::select(year, species_code, latitude, longitude, common_name) |> + unique() |> akgfmaps::transform_data_frame_crs(coords = c("longitude", "latitude"), in.crs = "+proj=longlat", out.crs = "EPSG:3338") @@ -277,7 +277,7 @@ make_data_summary <- function(dat_csv, region) { geom_sf(data = map_layers$survey.area, fill = NA, size = 0.7) + - geom_point(data = lw_dat.sub %>% + geom_point(data = lw_dat.sub |> dplyr::filter(common_name == ii), aes(x = longitude, y = latitude), size = 0.4, @@ -336,7 +336,7 @@ knitr::kable(out$n_specimen_by_vessel, caption = 'Table 3. Number of length-weig add_common_name <- function(x) { spp_df <- data.frame(common_name = c( - "walleye pollock", "walleye pollock (100–250 mm)", "walleye pollock (>250 mm)", "Pacific cod", + "walleye pollock", "walleye pollock (100-250 mm)", "walleye pollock (>250 mm)", "Pacific cod", "Pacific cod (juvenile)", "Pacific cod (adult)", "Atka mackerel", "arrowtooth flounder", "flathead sole", "yellowfin sole", "northern rock sole", "southern rock sole", "Alaska plaice", "Pacific ocean perch", "dusky rockfish", "northern rockfish", "Dover sole", "rex sole", "shortraker rockfish", "rougheye rockfish", "blackspotted rockfish", "sharpchin rockfish"), diff --git a/R/weight_lw_residuals.R b/R/weight_lw_residuals.R index 76de44a..bf92bff 100644 --- a/R/weight_lw_residuals.R +++ b/R/weight_lw_residuals.R @@ -45,20 +45,20 @@ weight_lw_residuals <- function(residuals, year, stratum = NA, stratum_biomass = } } else { # Weight residuals by stratum with biomass expansion (2020 ESR) - biomass_df <- data.frame(stratum_biomass, stratum, year) %>% + biomass_df <- data.frame(stratum_biomass, stratum, year) |> unique() - biomass_proportion_df <- biomass_df %>% - dplyr::group_by(year) %>% - dplyr::summarise(year_biomass = sum(stratum_biomass)) %>% - dplyr::inner_join(biomass_df) %>% - dplyr::mutate(stratum_weight = stratum_biomass/year_biomass) %>% + biomass_proportion_df <- biomass_df |> + dplyr::group_by(year) |> + dplyr::summarise(year_biomass = sum(stratum_biomass)) |> + dplyr::inner_join(biomass_df) |> + dplyr::mutate(stratum_weight = stratum_biomass/year_biomass) |> dplyr::select(year, stratum, stratum_weight) residuals_df <- data.frame(residuals, year, - stratum) %>% - dplyr::inner_join(biomass_proportion_df) %>% + stratum) |> + dplyr::inner_join(biomass_proportion_df) |> dplyr::mutate(weighted_residuals = residuals * stratum_weight) wtlw.res <- residuals_df$weighted_residuals diff --git a/data/AI_INDICATOR.rda b/data/AI_INDICATOR.rda index ef14bd7..6a12ff0 100644 Binary files a/data/AI_INDICATOR.rda and b/data/AI_INDICATOR.rda differ diff --git a/data/EBS_INDICATOR.rda b/data/EBS_INDICATOR.rda index 3a38de0..26c6d00 100644 Binary files a/data/EBS_INDICATOR.rda and b/data/EBS_INDICATOR.rda differ diff --git a/data/GOA_INDICATOR.rda b/data/GOA_INDICATOR.rda index 270a86d..0f66579 100644 Binary files a/data/GOA_INDICATOR.rda and b/data/GOA_INDICATOR.rda differ diff --git a/data/NBS_INDICATOR.rda b/data/NBS_INDICATOR.rda index 4db5cd7..69fd9c4 100644 Binary files a/data/NBS_INDICATOR.rda and b/data/NBS_INDICATOR.rda differ diff --git a/data/PCOD_ESP.rda b/data/PCOD_ESP.rda index b3fc695..fa9dcec 100644 Binary files a/data/PCOD_ESP.rda and b/data/PCOD_ESP.rda differ diff --git a/inst/extdata/raw_lw_bio.rda b/inst/extdata/raw_lw_bio.rda index 80c9ee1..f417226 100644 Binary files a/inst/extdata/raw_lw_bio.rda and b/inst/extdata/raw_lw_bio.rda differ diff --git a/inst/sql/VAST_ai_cpue.sql b/inst/sql/VAST_ai_cpue.sql index 188f281..0260fbc 100644 --- a/inst/sql/VAST_ai_cpue.sql +++ b/inst/sql/VAST_ai_cpue.sql @@ -1,12 +1,28 @@ /* Query to retrieve CPUE in kilograms per square kilometer from all Aleutian Islands strata for VAST. */ -select a.hauljoin, a.year, c.start_latitude latitude, c.start_longitude longitude, c.start_time, - b.species_code, b.common_name, a.number_fish, a.effort effort_km2, a.wgtcpue cpue_kg_km2 - from ai.cpue a, racebase.species b, racebase.haul c, race_data.cruises f, race_data.surveys g - where a.species_code = b.species_code - and f.survey_id = g.survey_id - and g.survey_definition_id = 52 - and c.cruisejoin = f.racebase_cruisejoin - and a.hauljoin = c.hauljoin - and a.species_code in (21740,21741,21720,10262,10110,30060,30152,21921,30420) - and (c.cruise >= 198401 and c.cruise != 198901) +SELECT + a.hauljoin, + a.start_latitude latitude, + a.start_longitude longitude, + a.start_time, + a.stationid, + b.species_code, + floor(a.cruise/100) year, + b.cpue_kgkm2 cpue_kg_km2, + b.count number_fish, + b.area_swept_km2 effort_km2, + d.common_name +FROM + racebase.haul a, + gap_products.cpue b, + race_data.v_cruises c, + racebase.species d +WHERE + a.cruise >= 199100 + AND a.hauljoin = b.hauljoin + AND a.haul_type = 3 + AND a.cruisejoin = c.cruisejoin + AND c.survey_definition_id = 52 + AND b.species_code = d.species_code + AND b.species_code IN (21740,21741,21720,10262,10110,30060,30152,21921,30420) +ORDER BY hauljoin; diff --git a/inst/sql/VAST_ebs_cpue.sql b/inst/sql/VAST_ebs_cpue.sql index f3c2fa6..cff35ee 100644 --- a/inst/sql/VAST_ebs_cpue.sql +++ b/inst/sql/VAST_ebs_cpue.sql @@ -1,13 +1,28 @@ /* Query to retrieve CPUE in kilograms per square kilometer from all eastern Bering Sea continental shelf strata for VAST. */ -select a.hauljoin, a.start_latitude latitude, a.start_longitude longitude, a.start_time, a.stationid, b.species_code, floor(a.cruise/100) year, -b.cpue_kgkm2 cpue_kg_km2, b.count number_fish, b.area_swept_km2 effort_km2, d.common_name - from racebase.haul a, gap_products.cpue b, race_data.v_cruises c, racebase.species d - where a.cruise > 199900 - and a.hauljoin = b.hauljoin - and a.haul_type in (3,13) - and a.cruisejoin = c.cruisejoin - and c.survey_definition_id = 98 - and b.species_code = d.species_code - and b.species_code in (21740, 21720, 10110, 10210, 10130, 10261, 10285) - order by hauljoin +SELECT + a.hauljoin, + a.start_latitude latitude, + a.start_longitude longitude, + a.start_time, + a.stationid, + b.species_code, + floor(a.cruise/100) year, + b.cpue_kgkm2 cpue_kg_km2, + b.count number_fish, + b.area_swept_km2 effort_km2, + d.common_name +FROM + racebase.haul a, + gap_products.cpue b, + race_data.v_cruises c, + racebase.species d +WHERE + a.cruise > 199900 + AND a.hauljoin = b.hauljoin + AND a.haul_type IN (3,13) + AND a.cruisejoin = c.cruisejoin + AND c.survey_definition_id = 98 + AND b.species_code = d.species_code + AND b.species_code IN (21740, 21720, 10110, 10210, 10130, 10261, 10285) +ORDER BY hauljoin diff --git a/inst/sql/VAST_ebs_length_weight.sql b/inst/sql/VAST_ebs_length_weight.sql index 52661a3..f082e0f 100644 --- a/inst/sql/VAST_ebs_length_weight.sql +++ b/inst/sql/VAST_ebs_length_weight.sql @@ -1,15 +1,36 @@ /* Query to retrieve length and weight samples for all eastern Bering Sea continental shelf strata for VAST. */ -select a.hauljoin, a.vessel, a.cruise, a.haul, a.stratum, floor(a.cruise/100) year, a.start_latitude latitude, a.start_longitude longitude, a.start_time, - b.species_code, b.specimenid, c.common_name, b.length length_mm, b.weight weight_g, b.sex, b.age - from racebase.haul a, racebase.specimen b, racebase.species c, race_data.v_cruises d - where a.hauljoin = b.hauljoin - and a.cruise >= 199900 - and b.species_code in (21740, 21720, 10110, 10210, 10130, 10261, 10285) - and b.species_code = c.species_code - and a.cruisejoin = d.cruisejoin and d.survey_definition_id = 98 - and a.haul_type in (3,13) - and a.region = 'BS' - and length != 0 - and weight != 0 - order by hauljoin \ No newline at end of file +SELECT + a.hauljoin, + a.vessel, + a.cruise, + a.haul, + a.stratum, + floor(a.cruise/100) year, + a.start_latitude latitude, + a.start_longitude longitude, + a.start_time, + b.species_code, + b.specimenid, + c.common_name, + b.length length_mm, + b.weight weight_g, + b.sex, + b.age +FROM + racebase.haul a, + racebase.specimen b, + racebase.species c, + race_data.v_cruises d +WHERE + a.hauljoin = b.hauljoin + AND a.cruise >= 199900 + AND b.species_code IN (21740, 21720, 10110, 10210, 10130, 10261, 10285) + AND b.species_code = c.species_code + AND a.cruisejoin = d.cruisejoin + AND d.survey_definition_id = 98 + AND a.haul_type IN (3,13) + AND a.region = 'BS' + AND length != 0 + AND weight != 0 +ORDER BY hauljoin; \ No newline at end of file diff --git a/inst/sql/VAST_goa_cpue.sql b/inst/sql/VAST_goa_cpue.sql index d0210e8..93b6e9d 100644 --- a/inst/sql/VAST_goa_cpue.sql +++ b/inst/sql/VAST_goa_cpue.sql @@ -1,11 +1,28 @@ /* Query to retrieve CPUE in kilograms per square kilometer from all Gulf of Alaska strata for VAST. */ -select a.hauljoin, a.year, c.start_latitude latitude, c.start_longitude longitude, c.start_time, - b.species_code, common_name, a.wgtcpue cpue_kg_km2, a.number_fish, a.effort effort_km2 - from goa.cpue a, racebase.species b, racebase.haul c, race_data.cruises f, race_data.surveys g - where a.species_code = b.species_code - and a.hauljoin = c.hauljoin - and f.survey_id = g.survey_id - and g.survey_definition_id = 47 - and c.cruisejoin = f.racebase_cruisejoin - and a.species_code in (21740,21741,21720,30420,10262,10110,30060,30152) \ No newline at end of file +SELECT + a.hauljoin, + a.start_latitude latitude, + a.start_longitude longitude, + a.start_time, + a.stationid, + b.species_code, + floor(a.cruise/100) year, + b.cpue_kgkm2 cpue_kg_km2, + b.count number_fish, + b.area_swept_km2 effort_km2, + d.common_name +FROM + racebase.haul a, + gap_products.cpue b, + race_data.v_cruises c, + racebase.species d +WHERE + a.cruise >= 199000 + AND a.hauljoin = b.hauljoin + AND a.haul_type = 3 + AND a.cruisejoin = c.cruisejoin + AND c.survey_definition_id = 47 + AND b.species_code = d.species_code + AND b.species_code IN (21740,21720,30420,10262,10110,30060,30152,10261,10180,10200,10130,30576,30051,30052,30560) +ORDER BY hauljoin; \ No newline at end of file diff --git a/inst/sql/VAST_nbs_cpue.sql b/inst/sql/VAST_nbs_cpue.sql index bf321ae..c3ed02c 100644 --- a/inst/sql/VAST_nbs_cpue.sql +++ b/inst/sql/VAST_nbs_cpue.sql @@ -1,13 +1,28 @@ /* Query to retrieve CPUE in kilograms per square kilometer from all northern Bering Sea strata for VAST. */ -select a.hauljoin, a.start_latitude latitude, a.start_longitude longitude, a.start_time, a.stationid, b.species_code, floor(a.cruise/100) year, -b.cpue_kgkm2 cpue_kg_km2, b.count number_fish, b.area_swept_km2 effort_km2, d.common_name - from racebase.haul a, gap_products.cpue b, race_data.v_cruises c, racebase.species d - where a.cruise > 200900 - and a.hauljoin = b.hauljoin - and a.haul_type in (3,13) - and a.cruisejoin = c.cruisejoin - and c.survey_definition_id = 143 - and b.species_code = d.species_code - and b.species_code in (21740, 21720, 10210, 10285) - order by hauljoin +SELECT + a.hauljoin, + a.start_latitude latitude, + a.start_longitude longitude, + a.start_time, + a.stationid, + b.species_code, + floor(a.cruise/100) year, + b.cpue_kgkm2 cpue_kg_km2, + b.count number_fish, + b.area_swept_km2 effort_km2, + d.common_name +FROM + racebase.haul a, + gap_products.cpue b, + race_data.v_cruises c, + racebase.species d +WHERE + a.cruise >= 199100 + AND a.hauljoin = b.hauljoin + AND a.haul_type IN (3, 13) + AND a.cruisejoin = c.cruisejoin + AND c.survey_definition_id = 143 + AND b.species_code = d.species_code + AND b.species_code IN (21740, 21720, 10210, 10285) +ORDER BY hauljoin; \ No newline at end of file diff --git a/inst/sql/VAST_nbs_length_weight.sql b/inst/sql/VAST_nbs_length_weight.sql index d3ee4e2..c414989 100644 --- a/inst/sql/VAST_nbs_length_weight.sql +++ b/inst/sql/VAST_nbs_length_weight.sql @@ -1,14 +1,35 @@ /* Query to retrieve length and weight samples for all northern Bering Sea for VAST. */ -select a.hauljoin, a.vessel, a.cruise, a.haul, a.stratum, floor(a.cruise/100) year, a.start_latitude latitude, a.start_longitude longitude, a.start_time, - b.species_code, b.specimenid, c.common_name, b.length length_mm, b.weight weight_g, b.sex, b.age - from racebase.haul a, racebase.specimen b, racebase.species c, race_data.v_cruises d - where a.hauljoin = b.hauljoin - and a.cruise > 201000 - and b.species_code in (21740, 21720, 10210, 10285) - and b.species_code = c.species_code - and a.cruisejoin = d.cruisejoin and d.survey_definition_id = 143 - and a.region = 'BS' - and length != 0 - and weight != 0 - order by hauljoin \ No newline at end of file +SELECT + a.hauljoin, + a.vessel, + a.cruise, + a.haul, + a.stratum, + floor(a.cruise/100) year, + a.start_latitude latitude, + a.start_longitude longitude, + a.start_time, + b.species_code, + b.specimenid, + c.common_name, + b.length length_mm, + b.weight weight_g, + b.sex, + b.age +FROM + racebase.haul a, + racebase.specimen b, + racebase.species c, + race_data.v_cruises d +WHERE + a.hauljoin = b.hauljoin + AND a.cruise > 201000 + AND b.species_code IN (21740, 21720, 10210, 10285) + AND b.species_code = c.species_code + AND a.cruisejoin = d.cruisejoin + AND d.survey_definition_id = 143 + AND a.region = 'BS' + AND length != 0 + AND weight != 0 +ORDER BY hauljoin; \ No newline at end of file diff --git a/inst/sql/ai_biomass.sql b/inst/sql/ai_biomass.sql index cfcf944..373604f 100644 --- a/inst/sql/ai_biomass.sql +++ b/inst/sql/ai_biomass.sql @@ -1,9 +1,17 @@ /* Query to retrieve stratum biomass for Aleutian Islands INPFC strata */ -select a.year, a.species_code, a.area_biomass, a.biomass_var, b.inpfc_area inpfc_stratum -from ai.biomass_inpfc a, (select inpfc_area, summary_area, sum(area) inpfc_Stratum_area from goa.goa_strata where nvl(stratum,1) != 0 -and survey = 'AI' group by inpfc_area, summary_area) b where (a.year >= 1984 and a.year != 1989) -and a.species_code in -(21740,21741,21720,30420,10262,10110,30060,21921) -and a.summary_Area = b.summary_area -order by a.year, a.species_code, b.inpfc_Area \ No newline at end of file +SELECT + b.species_code, + b.year, + a.area_name inpfc_stratum, + b.biomass_mt biomass, + b.biomass_var +FROM + gap_products.biomass b, + gap_products.akfin_area a +WHERE + b.species_code in (21740, 21720, 30420, 10262, 10110, 30060, 21921) + AND b.survey_definition_id = 52 + AND b.area_id in (299, 799, 3499, 5699) + AND a.survey_definition_id = b.survey_definition_id + AND a.area_id = b.area_id; \ No newline at end of file diff --git a/inst/sql/ai_length_weight.sql b/inst/sql/ai_length_weight.sql index f62ad3c..0cf703a 100644 --- a/inst/sql/ai_length_weight.sql +++ b/inst/sql/ai_length_weight.sql @@ -1,14 +1,37 @@ /* Query to retrieve length-weight samples for Aleutian Islands bottom trawl surveys */ -select a.haul, a.vessel, a.cruise, floor(a.cruise/100) year, b.species_code, a.region, c.inpfc_area inpfc_stratum, a.start_latitude latitude, a.start_longitude longitude, d.specimenid, a.start_time, -e.common_name, d.sex, d.length length_mm, d.weight weight_g, d.age -from racebase.haul a, racebase.catch b, goa.goa_strata c, racebase.specimen d, racebase.species e, race_data.cruises f, race_data.surveys g -where a.region = 'AI' and (a.cruise >= 198401 and a.cruise != 198901) -and b.species_code in (21740,21741,21720,30420,10262,10110,30060,21921) -and a.hauljoin = b.hauljoin -and b.hauljoin = d.hauljoin -and f.survey_id = g.survey_id -and g.survey_definition_id = 52 -and a.cruisejoin = f.racebase_cruisejoin -and b.species_code = e.species_code and b.species_code = d.species_code and a.stratum = c.stratum -and a.region = c.survey and d.length != 0 and d.weight != 0 \ No newline at end of file +SELECT + h.haul, + c.vessel_id vessel, + c.cruise, + c.year, + s.species_code, + a.area_name inpfc_stratum, + h.latitude_dd_start latitude, + h.longitude_dd_start longitude, + h.date_time_start start_time, + s.specimen_id specimenid, + t.common_name, + s.sex, + s.length_mm, + s.weight_g, + s.age +FROM + gap_products.akfin_specimen s, + gap_products.akfin_area a, + gap_products.akfin_haul h, + gap_products.akfin_cruise c, + gap_products.akfin_taxonomic_groups t, + gap_products.akfin_stratum_groups sg +WHERE + s.species_code IN (21740,21741,21720,30420,10262,10110,30060,21921) + AND c.survey_definition_id = 52 + AND a.area_id IN (299, 799, 3499, 5699) + AND a.survey_definition_id = c.survey_definition_id + AND h.hauljoin = s.hauljoin + AND h.cruisejoin = c.cruisejoin + AND h.stratum = sg.stratum + AND sg.area_id = a.area_id + AND s.weight_g != 0 + AND s.length_mm != 0 + AND s.species_code = t.species_code; \ No newline at end of file diff --git a/inst/sql/ebs_biomass.sql b/inst/sql/ebs_biomass.sql index 54eb8f1..fd598c2 100644 --- a/inst/sql/ebs_biomass.sql +++ b/inst/sql/ebs_biomass.sql @@ -1,5 +1,12 @@ -SELECT species_code, year, area_id*10 stratum, biomass_mt biomass, biomass_var -FROM gap_products.biomass -WHERE (((species_code) in (21740, 21741, 21720, 10110, 10210, 10130, 10261, 10285))) and -survey_definition_id = 98 -and area_id between 1 and 6 \ No newline at end of file +SELECT + species_code, + year, + area_id*10 stratum, + biomass_mt biomass, + biomass_var +FROM + gap_products.biomass +WHERE + species_code IN (21740, 21741, 21720, 10110, 10210, 10130, 10261, 10285) + AND survey_definition_id = 98 + AND area_id BETWEEN 1 AND 6; \ No newline at end of file diff --git a/inst/sql/ebs_length_weight.sql b/inst/sql/ebs_length_weight.sql index 27ab7ef..5593a78 100644 --- a/inst/sql/ebs_length_weight.sql +++ b/inst/sql/ebs_length_weight.sql @@ -1,12 +1,33 @@ /* Query to retrieve length and weight samples for eastern Bering Sea continental shelf strata 10-62. */ -select a.haul, a.vessel, a.cruise, floor(a.cruise/100) year, b.species_code, a.region, a.start_time, -decode(a.stratum, 31, 30, 32, 30, 61, 60, 62, 60, 41, 40, 42, 40, 43, 40, a.stratum) stratum, -a.start_latitude latitude, a.start_longitude longitude, d.specimenid, -e.common_name, d.sex, d.length, d.weight -from racebase.haul a, racebase.catch b, racebase.specimen d , racebase.species e, race_data.v_cruises c -where a.region = 'BS' and a.cruise >= 199900 -and b.species_code in (21740, 21720, 10110, 10210, 10130, 10261, 10285) -and a.hauljoin = b.hauljoin and b.hauljoin = d.hauljoin -and b.species_code = e.species_code and b.species_code = d.species_code and d.length != 0 and d.weight != 0 -and a.cruisejoin = c.cruisejoin and c.survey_definition_id = 98 and a.stratum not in (82, 90) -and a.stratum between 10 and 62 + +SELECT + h.haul, + c.vessel_id vessel, + c.cruise, + c.year, + s.species_code, + h.stratum, + h.latitude_dd_start latitude, + h.longitude_dd_start longitude, + h.date_time_start start_time, + s.specimen_id specimenid, + t.common_name, + s.sex, + s.length_mm, + s.weight_g, + s.age +FROM + gap_products.akfin_specimen s, + gap_products.akfin_haul h, + gap_products.akfin_cruise c, + gap_products.akfin_taxonomic_groups t +WHERE + s.species_code IN (21740, 21720, 10110, 10210, 10130, 10261, 10285) + AND c.survey_definition_id = 98 + AND h.hauljoin = s.hauljoin + AND h.cruisejoin = c.cruisejoin + AND h.stratum BETWEEN 10 and 62 + AND s.weight_g != 0 + AND s.length_mm != 0 + AND s.species_code = t.species_code + AND c.year >= 1999; diff --git a/inst/sql/goa_biomass.sql b/inst/sql/goa_biomass.sql index 169a411..d3ddd17 100644 --- a/inst/sql/goa_biomass.sql +++ b/inst/sql/goa_biomass.sql @@ -1,9 +1,18 @@ /* Query to retrieve stratum biomass for Gulf of Alaska INPFC strata */ -select a.year, a.species_code, a.area_biomass, a.biomass_var, b.inpfc_area inpfc_stratum -from goa.biomass_inpfc a, (select inpfc_area, summary_area, sum(area) inpfc_Stratum_area from goa.goa_strata where nvl(stratum,1) != 0 -and survey = 'GOA' group by inpfc_area, summary_area) b where (a.year >= 1984 and a.year != 1989 and a.year != 1985) -and a.species_code in -(21740,21741,21720,30420,10262,10110,30060,30152,10261,10180,10200,10130,30576,30051,30052,30560) -and a.summary_Area = b.summary_area -order by a.year, a.species_code, b.inpfc_Area \ No newline at end of file +SELECT + b.species_code, + b.year, + b.area_id, + a.area_name inpfc_stratum, + b.biomass_mt biomass, + b.biomass_var +FROM + gap_products.biomass b, + gap_products.akfin_area a +WHERE + b.species_code in (21740,21741,21720,30420,10262,10110,30060,30152,10261,10180,10200,10130,30576,30051,30052,30560) + AND b.survey_definition_id = 47 + AND b.area_id in (919, 929, 939, 949, 959) + AND a.survey_definition_id = b.survey_definition_id + AND a.area_id = b.area_id; \ No newline at end of file diff --git a/inst/sql/goa_length_weight.sql b/inst/sql/goa_length_weight.sql index 08afab4..1bfdc55 100644 --- a/inst/sql/goa_length_weight.sql +++ b/inst/sql/goa_length_weight.sql @@ -1,13 +1,37 @@ /* Query to retrieve length-weight samples for Gulf of Alaska bottom trawl surveys */ -select a.haul, a.vessel, a.cruise, floor(a.cruise/100) year, b.species_code, a.region, c.inpfc_area inpfc_stratum, a.start_latitude latitude, a.start_longitude longitude, a.start_time, -d.specimenid, e.common_name, d.sex, d.length length_mm, d.weight weight_g, d.age -from racebase.haul a, racebase.catch b, goa.goa_strata c, racebase.specimen d, racebase.species e, race_data.cruises f, race_data.surveys g -where a.region = 'GOA' and (a.cruise >= 198401 and a.cruise != 198901) -and b.species_code in (21740,21741,21720,30420,10262,10110,30060,30152,10261,10180,10200,10130,30576,30051,30052,30560) -and a.hauljoin = b.hauljoin and b.hauljoin = d.hauljoin -and f.survey_id = g.survey_id -and g.survey_definition_id = 47 -and a.cruisejoin = f.racebase_cruisejoin -and b.species_code = e.species_code and b.species_code = d.species_code and a.stratum = c.stratum -and a.region = c.survey and d.length != 0 and d.weight != 0 \ No newline at end of file +SELECT + h.haul, + c.vessel_id vessel, + c.cruise, + c.year, + s.species_code, + a.area_name inpfc_stratum, + h.latitude_dd_start latitude, + h.longitude_dd_start longitude, + h.date_time_start start_time, + s.specimen_id specimenid, + t.common_name, + s.sex, + s.length_mm, + s.weight_g, + s.age +FROM + gap_products.akfin_specimen s, + gap_products.akfin_area a, + gap_products.akfin_haul h, + gap_products.akfin_cruise c, + gap_products.akfin_taxonomic_groups t, + gap_products.akfin_stratum_groups sg +WHERE + s.species_code IN (21740,21720,30420,10262,10110,30060,30152,10261,10180,10200,10130,30576,30051,30052,30560) + AND c.survey_definition_id = 47 + AND a.area_id IN (919, 929, 939, 949, 959) + AND a.survey_definition_id = c.survey_definition_id + AND h.hauljoin = s.hauljoin + AND h.cruisejoin = c.cruisejoin + AND h.stratum = sg.stratum + AND sg.area_id = a.area_id + AND s.weight_g != 0 + AND s.length_mm != 0 + AND s.species_code = t.species_code; \ No newline at end of file diff --git a/inst/sql/nbs_biomass.sql b/inst/sql/nbs_biomass.sql index ceabdf6..6943935 100644 --- a/inst/sql/nbs_biomass.sql +++ b/inst/sql/nbs_biomass.sql @@ -1,5 +1,12 @@ -SELECT species_code, year, area_id stratum, biomass_mt biomass, biomass_var -FROM gap_products.biomass -WHERE (((species_code) in (21740, 21741, 21720, 10210, 10285))) and -survey_definition_id = 143 -and area_id < 90 \ No newline at end of file +SELECT + species_code, + year, + area_id stratum, + biomass_mt biomass, + biomass_var +FROM + gap_products.biomass +WHERE + species_code IN (21740, 21741, 21720, 10210, 10285) + AND survey_definition_id = 143 + AND area_id < 90; \ No newline at end of file diff --git a/inst/sql/nbs_length_weight.sql b/inst/sql/nbs_length_weight.sql index 07f1d72..1426c71 100644 --- a/inst/sql/nbs_length_weight.sql +++ b/inst/sql/nbs_length_weight.sql @@ -1,13 +1,32 @@ /* Query to retrieve length and weight samples for northern Bering Sea strata. */ -select a.haul, a.vessel, a.cruise, floor(a.cruise/100) year, b.species_code, a.region, -decode(a.stratum, 31, 30, 32, 30, 61, 60, 62, 60, 41, 40, 42, 40, 43, 40, a.stratum) stratum, -a.start_latitude latitude, a.start_longitude longitude, a.start_time, d.specimenid, -e.common_name, d.sex, d.length length_mm, d.weight weight_g -from racebase.haul a, racebase.catch b, racebase.specimen d , racebase.species e, race_data.v_cruises c -where a.region = 'BS' and a.cruise >= 201000 -and b.species_code in (21740, 21720, 10210, 10285) -and a.hauljoin = b.hauljoin and b.hauljoin = d.hauljoin -and b.species_code = e.species_code and b.species_code = d.species_code and d.length != 0 and d.weight != 0 -and a.cruisejoin = c.cruisejoin and c.survey_definition_id = 143 and a.stratum not in (82, 90) -and a.stratum between 70 and 99 \ No newline at end of file +SELECT + h.haul, + c.vessel_id vessel, + c.cruise, + c.year, + s.species_code, + h.stratum, + h.latitude_dd_start latitude, + h.longitude_dd_start longitude, + h.date_time_start start_time, + s.specimen_id specimenid, + t.common_name, + s.sex, + s.length_mm, + s.weight_g, + s.age +FROM + gap_products.akfin_specimen s, + gap_products.akfin_haul h, + gap_products.akfin_cruise c, + gap_products.akfin_taxonomic_groups t +WHERE + s.species_code IN (21740, 21720, 10210, 10285) + AND c.survey_definition_id = 143 + AND h.hauljoin = s.hauljoin + AND h.cruisejoin = c.cruisejoin + AND h.stratum BETWEEN 70 and 81 + AND s.weight_g != 0 + AND s.length_mm != 0 + AND s.species_code = t.species_code; \ No newline at end of file diff --git a/man/AI_INDICATOR.Rd b/man/AI_INDICATOR.Rd index 1fd05a0..63f7dbc 100644 --- a/man/AI_INDICATOR.Rd +++ b/man/AI_INDICATOR.Rd @@ -3,7 +3,7 @@ \docType{data} \name{AI_INDICATOR} \alias{AI_INDICATOR} -\title{Aleutian Islands (1986-2022)} +\title{Aleutian Islands (1991-2024)} \format{ A list containing two data frames (indicator for the full region, indicator by stratum) and a character vector. \describe{ @@ -14,24 +14,17 @@ A list containing two data frames (indicator for the full region, indicator by s \item{common_name}{: Year} \item{mean_wt_resid}{: Mean residual for the full region} \item{se_wt_resid}{: Standard error of the indicator for the full region} - \item{vast_condition}{: Estimate of the allometric intercept,a, of the length-weight equation (W = aL^b), estimated using VAST.} - \item{vast_condition_se}{Standard error of a, estimated using VAST.} - \item{species_code}{: RACE/GAP species code} - \item{scaled_vast_condition}{: Anomaly (Z-score) of vast_condition} - \item{vast_relative_condition}{: Relative condition based on the allometric intercept estimated using VAST, i.e., allometric intercept divided by the time series mean allometric intercept.} - \item{vast_relative_condition_se}{: Standard error of vast_relative_condition.} } \item STRATUM (data frame): Residuals by stratum \itemize{ - \item{year}{: Year} \item{common_name}{: Species Common name} - \item{species_code}{: RACE/GAP species code} + \item{year}{: Year} \item{inpfc_stratum}{: INPFC Stratum} \item{stratum_resid_mean}{: Unweighted stratum mean residual} \item{n}{: Sample size for the stratum} \item{stratum_resid_se}{: Standard error of the stratum mean residual} \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} - \item{weighted_resid_mean}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} + \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} } \item LAST_UPDATE: Date when indicator and data were last updated. } @@ -44,6 +37,6 @@ A list containing two data frames (indicator for the full region, indicator by s AI_INDICATOR } \description{ -Morphometric condition indicators based on the allometric intercept estimated using VAST and residuals from a length-weight regression for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), southern rock sole, arrowtooth flounder, Atka mackerel, northern rockfish, and Pacific ocean perch in the Aleutian Islands. +Morphometric condition indicators based on residuals from length-weight regressions for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), southern rock sole, arrowtooth flounder, Atka mackerel, northern rockfish, and Pacific ocean perch in the Aleutian Islands. } \keyword{datasets} diff --git a/man/EBS_INDICATOR.Rd b/man/EBS_INDICATOR.Rd index 636e822..0ddc53f 100644 --- a/man/EBS_INDICATOR.Rd +++ b/man/EBS_INDICATOR.Rd @@ -14,24 +14,18 @@ A list containing two data frames (indicator for the full region, indicator by s \item{common_name}{: Year} \item{mean_wt_resid}{: Mean residual for the full region} \item{se_wt_resid}{: Standard error of the indicator for the full region} - \item{vast_condition}{: Estimate of the allometric intercept,a, of the length-weight equation (W = aL^b), estimated using VAST.} - \item{vast_condition_se}{Standard error of a, estimated using VAST.} \item{species_code}{: RACE/GAP species code} - \item{scaled_vast_condition}{: Anomaly (Z-score) of vast_condition} - \item{vast_relative_condition}{: Relative condition based on the allometric intercept estimated using VAST, i.e., allometric intercept divided by the time series mean allometric intercept.} - \item{vast_relative_condition_se}{: Standard error of vast_relative_condition.} } \item STRATUM (data frame): Residuals by stratum \itemize{ - \item{year}{: Year} \item{common_name}{: Species Common name} - \item{species_code}{: RACE/GAP species code} - \item{inpfc_stratum}{: INPFC Stratum} + \item{year}{: Year} + \item{stratum}{: Survey Stratum} \item{stratum_resid_mean}{: Unweighted stratum mean residual} - \item{n}{: Sample size for the stratum} \item{stratum_resid_sd}{: Standard error of the stratum mean residual} + \item{n}{: Sample size for the stratum} \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} - \item{weighted_resid_mean}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} + \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} } \item LAST_UPDATE: Date when indicator and data were last updated. } @@ -44,6 +38,6 @@ A list containing two data frames (indicator for the full region, indicator by s EBS_INDICATOR } \description{ -Morphometric condition indicators based on the allometric intercept estimated using VAST and residuals from a length-weight regression for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), northern rock sole, arrowtooth flounder, flathead sole, Alaska plaice, and yellowfin sole in the eastern Bering Sea. +Morphometric condition indicators based on length-weight regressions for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), northern rock sole, arrowtooth flounder, flathead sole, Alaska plaice, and yellowfin sole in the eastern Bering Sea. } \keyword{datasets} diff --git a/man/GOA_INDICATOR.Rd b/man/GOA_INDICATOR.Rd index 31e19da..75db315 100644 --- a/man/GOA_INDICATOR.Rd +++ b/man/GOA_INDICATOR.Rd @@ -3,29 +3,28 @@ \docType{data} \name{GOA_INDICATOR} \alias{GOA_INDICATOR} -\title{Gulf of Alaska (1984-2021)} +\title{Gulf of Alaska (1990-2023)} \format{ A list containing two data frames (indicator for the full region, indicator by stratum) and a character vector. \describe{ \itemize{ \item FULL_REGION (data frame): Residuals for the full region \itemize{ - \item{YEAR}{: Year} - \item{COMMON_NAME}{: Year} + \item{year}{: Year} + \item{common_name}{: Year} \item{mean_wt_resid}{: Mean residual for the full region} \item{se_wt_resid}{: Standard error of the indicator for the full region} } \item STRATUM (data frame): Residuals by stratum \itemize{ - \item{YEAR}{: Year} - \item{COMMON_NAME}{: Species Common name} - \item{SPECIES_CODE}{: RACE/GAP species code} - \item{INPFC_STRATUM}{: INPFC Stratum} + \item{common_name}{: Species Common name} + \item{year}{: Year} + \item{inpfc_stratum}{: INPFC Stratum} \item{stratum_resid_mean}{: Unweighted stratum mean residual} \item{n}{: Sample size for the stratum} \item{stratum_resid_sd}{: Standard error of the stratum mean residual} \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} - \item{weighted_resid_mean}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} + \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} } \item LAST_UPDATE: Date when indicator and data were last updated. } @@ -38,6 +37,6 @@ A list containing two data frames (indicator for the full region, indicator by s GOA_INDICATOR } \description{ -Morphometric condition indicators based on residuals from a length-weight regression for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), southern rock sole, arrowtooth flounder, Atka mackerel, northern rockfish, and Pacific ocean perch in the Gulf of Alaska. +Morphometric condition indicators based on residuals from length-weight regressions for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), Pacific cod, northern rock sole, southern rock sole, arrowtooth flounder, Dover sole, rex sole, flathead sole, Atka mackerel, blackspotted rockfish, rougheye rockfish, northern rockfish, shortraker rockfish, shaprchin rockfish, and Pacific ocean perch in the Gulf of Alaska. } \keyword{datasets} diff --git a/man/NBS_INDICATOR.Rd b/man/NBS_INDICATOR.Rd index 9545e8a..0ae54a7 100644 --- a/man/NBS_INDICATOR.Rd +++ b/man/NBS_INDICATOR.Rd @@ -10,17 +10,24 @@ A list containing a data frames (indicator for the full region) and a character \itemize{ \item FULL_REGION (data frame): Residuals for the full region \itemize{ - \item{YEAR}{: Year} - \item{COMMON_NAME}{: Species common name} + \item{year}{: Year} + \item{common_name}{: Species common name} \item{mean_wt_resid}{: Mean residual for the full region} \item{se_wt_resid}{: Standard error of the indicator for the full region} - \item{vast_condition}{: Estimate of the allometric intercept,a, of the length-weight equation (W = aL^b), estimated using VAST.} - \item{vast_condition_se}{Standard error of a, estimated using VAST.} \item{species_code}{: RACE/GAP species code} - \item{scaled_vast_condition}{: Anomaly (Z-score) of vast_condition} - \item{vast_relative_condition}{: Relative condition based on the allometric intercept estimated using VAST, i.e., allometric intercept divided by the time series mean allometric intercept.} \item{vast_relative_condition_se}{: Standard error of vast_relative_condition.} } + \item STRATUM (data frame): Residuals by stratum + \itemize{ + \item{common_name}{: Species Common name} + \item{year}{: Year} + \item{stratum}{: Survey stratum} + \item{stratum_resid_mean}{: Unweighted stratum mean residual} + \item{stratum_resid_sd}{: Standard error of the stratum mean residual} + \item{n}{: Sample size for the stratum} + \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} + \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} + } \item LAST_UPDATE: Date when indicator and data were last updated. } } @@ -32,6 +39,6 @@ A list containing a data frames (indicator for the full region) and a character NBS_INDICATOR } \description{ -Morphometric condition indicators based on the allometric intercept estimated using VAST and residuals from a length-weight regression for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), Alaska plaice, and yellowfin sole in the northern Bering Sea. +Morphometric condition indicators based on length-weight regressions for Pacific cod, walleye pollock (> 250 mm), walleye pollock (100-250 mm), Alaska plaice, and yellowfin sole in the northern Bering Sea. } \keyword{datasets} diff --git a/man/PCOD_ESP.Rd b/man/PCOD_ESP.Rd index 52cb4a1..5df35af 100644 --- a/man/PCOD_ESP.Rd +++ b/man/PCOD_ESP.Rd @@ -14,8 +14,6 @@ Eight data frames (indicator for the full region and indicator by stratum for th \item{common_name}{: Species common name} \item{mean_wt_resid}{: Mean residual for the full region} \item{se_wt_resid}{: Standard error of the indicator for the full region} - \item{vast_condition}{: Estimate of the allometric intercept,a, of the length-weight equation (W = aL^b), estimated using VAST.} - \item{vast_condition_se}{Standard error of a, estimated using VAST.} } \item STRATUM_* (data frame): Residuals by stratum \itemize{ @@ -27,7 +25,7 @@ Eight data frames (indicator for the full region and indicator by stratum for th \item{n}{: Sample size for the stratum} \item{stratum_resid_sd}{: Standard error of the stratum mean residual} \item{weighted_resid_mean}{: Stratum mean Residual weighted in proportion to stratum biomass} - \item{weighted_resid_mean}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} + \item{weighted_resid_se}{: Standard error of the stratum mean residual weighted in proportion to stratum biomass} } \item LAST_UPDATE: Date when indicator and data were last updated. } diff --git a/man/make_vast_table.Rd b/man/make_vast_table.Rd index a881235..a61bc63 100644 --- a/man/make_vast_table.Rd +++ b/man/make_vast_table.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/make_plots.R +% Please edit documentation in R/make_vast_table.R \name{make_vast_table} \alias{make_vast_table} \title{Generate tables with b coefficient and VAST Settings} @@ -8,8 +8,6 @@ make_vast_table(region, write_table = TRUE) } \arguments{ \item{region}{One or more region as a character vector ("EBS", "NBS", "AI", "GOA")} - -\item{write_table}{Should the table be written as a .csv file to /plots/[region]/[region]_VAST_summary_table.csv?} } \description{ Generate tables containing settings for VAST and the VAST estimate of the allometric slope of the length-weight (b in W = a*L^b) relationship. diff --git a/man/plot_anomaly_timeseries.Rd b/man/plot_anomaly_timeseries.Rd index 15d1844..044dd19 100644 --- a/man/plot_anomaly_timeseries.Rd +++ b/man/plot_anomaly_timeseries.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/make_plots.R +% Please edit documentation in R/plot_anomaly_timeseries.R \name{plot_anomaly_timeseries} \alias{plot_anomaly_timeseries} \title{Make condition time series plot} @@ -34,7 +34,7 @@ plot_anomaly_timeseries( \item{format_for}{"rmd" or "png"} -\item{set_intercept}{Intercept to use for plots. Numeric (1L). 0 for residual, 1 for vAST relative condition, NULL for VAST a} +\item{set_intercept}{Intercept to use for plots. Numeric (1L). 0 for residual, 1 for VAST relative condition, NULL for VAST a} } \value{ A ggplot object of the time series @@ -43,6 +43,7 @@ A ggplot object of the time series Make condition time series plot } \examples{ +\dontrun{ # EBS anomaly timeseries plot names(EBS_INDICATOR$FULL_REGION) @@ -55,3 +56,4 @@ var_x_name = "year", y_title = "Length-weight residual (ln(g))", format_for = "png") } +} diff --git a/man/plot_species_bar.Rd b/man/plot_species_bar.Rd index c723fbc..0200155 100644 --- a/man/plot_species_bar.Rd +++ b/man/plot_species_bar.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/make_plots.R +% Please edit documentation in R/plot_species_bar.R \name{plot_species_bar} \alias{plot_species_bar} \title{Make single species bar plots} @@ -29,7 +29,7 @@ plot_species_bar( \item{y_title}{Y-axis title. Character (1L).} -\item{set_intercept}{Intercept to use for plots. Numeric (1L). 0 for residual, 1 for vAST relative condition, NULL for VAST a} +\item{set_intercept}{Intercept to use for plots. Numeric (1L). 0 for residual, 1 for VAST relative condition, NULL for VAST a} \item{write_plot}{Should plots be written to the /plot/ directory?} diff --git a/man/plot_species_stratum_bar.Rd b/man/plot_species_stratum_bar.Rd index a55cea9..a7e4485 100644 --- a/man/plot_species_stratum_bar.Rd +++ b/man/plot_species_stratum_bar.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/make_plots.R +% Please edit documentation in R/plot_stratum_species_bar.R \name{plot_species_stratum_bar} \alias{plot_species_stratum_bar} \title{Make single species stratum bar plots} @@ -34,9 +34,9 @@ plot_species_stratum_bar( \item{fill_title}{Name of the fill variable to use for plotting. Character (1L).} -\item{write_plot}{Should plots be written to the /plot/ directory?} +\item{fill_palette}{Character vector denoting which color palette to use. Must be a valid name for RColorBrewer::brewer.pal(name = {fill_palette})} -\item{fill_pallete}{Character vector denoting which color pallette to use. Must be a valid name for RColorBrewer::brewer.pal(name = {fill_palette})} +\item{write_plot}{Should plots be written to the /plot/ directory?} } \value{ A list of bar plots as ggplot objects diff --git a/man/plot_stratum_stacked_bar.Rd b/man/plot_stratum_stacked_bar.Rd index f8a5a10..27551b9 100644 --- a/man/plot_stratum_stacked_bar.Rd +++ b/man/plot_stratum_stacked_bar.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/make_plots.R +% Please edit documentation in R/plot_stratum_stacked_bar.R \name{plot_stratum_stacked_bar} \alias{plot_stratum_stacked_bar} \title{Make condition time series stacked bar plot} @@ -39,6 +39,7 @@ A ggplot object of the time series Make condition time series stacked bar plot } \examples{ +\dontrun{ # EBS strata stacked bar plot names(EBS_INDICATOR$STRATUM) @@ -52,3 +53,4 @@ var_group_name = "stratum", fill_title = "Stratum", y_title = "Length-weight residual (ln(g))") } +} diff --git a/man/plot_two_timeseries.Rd b/man/plot_two_timeseries.Rd index 746f2c7..c31fc56 100644 --- a/man/plot_two_timeseries.Rd +++ b/man/plot_two_timeseries.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/make_plots.R +% Please edit documentation in R/plot_two_timeseries.R \name{plot_two_timeseries} \alias{plot_two_timeseries} \title{Two condition time series and correlation} diff --git a/man/plot_xy_corr.Rd b/man/plot_xy_corr.Rd index 92d8403..f4b1b95 100644 --- a/man/plot_xy_corr.Rd +++ b/man/plot_xy_corr.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/make_plots.R +% Please edit documentation in R/plot_xy_corr.R \name{plot_xy_corr} \alias{plot_xy_corr} \title{Plot of condition versus condition} diff --git a/past_years/2023/EBS_GroundfishCondition_2023.Rmd b/past_years/2023/EBS_GroundfishCondition_2023.Rmd new file mode 100644 index 0000000..79b3ea4 --- /dev/null +++ b/past_years/2023/EBS_GroundfishCondition_2023.Rmd @@ -0,0 +1,360 @@ +--- +title: Eastern and Northern Bering Sea Groundfish Condition +author: +- affiliation: RACE + description: Research Fisheries Biologist + email: Sean.Rohan@NOAA.gov + name: Sean Rohan +output: word_document +fontsize: 12pt +bibliography: CopyOfEBS_references.bib +csl: fish-and-fisheries.csl +addr: + l1: 7600 Sand Point Way NE + l2: NMFS RACE Division, Groundfish Assessment Program + l3: Seattle, WA 98115 +--- + +```{r setup, include=FALSE} +# Load packages +library(knitr) +library(akfishcondition) + +region <- "BS" + +make_map <- TRUE +pkg_version <- packageVersion("akfishcondition") + +# Unzip map files +unzip(system.file("extdata/2022_ESR.zip", package = "akfishcondition")) + +if(!dir.exists(here::here("plots", region))) { + dir.create(here::here("plots", region), recursive = TRUE) +} + +if(make_map){ + library(akgfmaps) + + ebs_layers <- akgfmaps::get_base_layers(select.region = "ebs", + set.crs = "EPSG:3338") + + + # Simple version + survey_strata <- ebs_layers$survey.strata %>% + dplyr::mutate(agg_stratum = Stratum) %>% + dplyr::mutate(agg_stratum = replace(agg_stratum, agg_stratum %in% c(31,32), 30), + agg_stratum = replace(agg_stratum, agg_stratum %in% c(41,42,43), 40), + agg_stratum = replace(agg_stratum, agg_stratum %in% c(61,62), 60), + SURVEY = replace(SURVEY, SURVEY == "EBS_SHELF", "EBS Shelf"), + SURVEY = replace(SURVEY, SURVEY == "NBS_SHELF", "NBS")) %>% + dplyr::group_by(agg_stratum, SURVEY) %>% + dplyr::summarise() + + plot_ebs_nbs_survey_stations <- ggplot() + + geom_sf(data = survey_strata, + aes(fill = SURVEY), + color = "black") + + geom_sf_text(data = sf::st_centroid(survey_strata), + aes(label = agg_stratum), + color = "black", + size = rel(6)) + + geom_sf(data = sf::st_centroid(ebs_layers$survey.grid %>% + dplyr::mutate(lab = "Station Location")), + aes(shape = lab)) + + ggplot2::geom_sf(data = ebs_layers$akland, + fill = "grey85", + color = "black") + + geom_sf(data = ebs_layers$graticule, + size = 0.2, + alpha = 0.5) + + ggplot2::coord_sf(xlim = ebs_layers$plot.boundary$x, + ylim = ebs_layers$plot.boundary$y) + + ggplot2::scale_x_continuous(name = "Longitude", + breaks = ebs_layers$lon.breaks) + + ggplot2::scale_y_continuous(name = "Latitude", + breaks = ebs_layers$lat.breaks) + + scale_shape_manual(values = "x") + + scale_fill_manual(values = c("#0085CA", "#54ADDB")) + + theme_bw() + + ggplot2::theme(axis.title = element_blank(), + axis.text = element_text(color = "black"), + axis.ticks = element_line(color = "black"), + panel.border = element_rect(color = "black", fill = NA), + panel.grid = element_blank(), + panel.background = element_rect(color = "black", fill = "white"), + legend.title = element_blank(), + legend.margin = margin(-12,0,0,0), + legend.position = c(0.2,0.1), + legend.background = element_blank()) + + ragg::agg_png(filename = here::here("plots", region, "EBS_NBS_survey_area.png"), width = 5, height = 5, units = "in", res = 600) + print(plot_ebs_nbs_survey_stations) + dev.off() +} + +``` + +Contributed by Bianca Prohaska^1^ and Sean Rohan^1^ + +^1^Resource Assessment and Conservation Engineering Division, Alaska Fisheries Science Center, National Marine Fisheries Service, NOAA +**Contact**: sean.rohan@noaa.gov +**Last updated**: October 2023 + +**Description of Indicator**: Length-weight residuals represent how heavy a fish is per unit body length and are an indicator of somatic growth variability [@Brodeur2004]. Therefore, length-weight residuals can be considered indicators of prey availability, growth, general health, and habitat condition [@Blackwell2000; @Froese2006]. Positive length-weight residuals indicate better condition (i.e., heavier per unit length) and negative residuals indicate poorer condition (i.e., lighter per unit length) [@Froese2006]. Fish condition calculated in this way reflects realized outcomes of intrinsic and extrinsic processes that affect fish growth which can have implications for biological productivity through direct effects on growth and indirect effects on demographic processes such as, reproduction and mortality (e.g., [@Rodgveller2019; @Barbeaux2020]). + +```{r map, include=TRUE,out.width="200%",fig.cap="\\label{fig:figs}Figure 1. AFSC/RACE GAP summer bottom trawl survey strata (10-90) and station locations (x) in the eastern Bering Sea and northern Bering Sea. ", echo=FALSE} + knitr::include_graphics(path = here::here("plots", region, "ebs_nbs_survey_area.png")) +``` + +The groundfish morphometric condition indicator is calculated from paired fork lengths (mm) and weights (g) of individual fishes that were collected during bottom trawl surveys of the eastern Bering Sea (EBS) shelf and northern Bering Sea (NBS), which were conducted by the Alaska Fisheries Science Center’s Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP). Fish condition analyses were applied to walleye pollock (_Gadus chalcogrammus_), Pacific cod (_Gadus macrocephalus_), arrowtooth flounder (_Atheresthes stomias_), yellowfin sole (_Limanda aspera_), flathead sole (_Hippoglossoides elassodon_), northern rock sole (_Lepidopsetta polyxystra_), and Alaska plaice (_Pleuronectes quadrituberculatus_) collected in bottom trawls at standard survey stations (Figure 1). For these analyses and results, survey strata 31 and 32 were combined as stratum 30; strata 41, 42, and 43 were combined as stratum 40; and strata 61 and 62 were combined as stratum 60. Northwest survey strata 82 and 90 were excluded from these analyses. + +To calculate indicators, length-weight relationships were estimated from linear regression models based on a log-transformation of the exponential growth relationship, _W_ = _aL_^_b_^, where _W_ is weight (g) and _L_ is fork length (mm) for all areas for the period 1997–2023 (EBS: 1997–2023, NBS: 2010, 2017, 2019, 2021-2023). Unique intercepts (*a*) and slopes (*b*) were estimated for each survey stratum, sex, and interaction between stratum and sex to account for sexual dimorphism and spatial-temporal variation in growth and bottom trawl survey sampling. Length-weight relationships for 100–250 mm fork length walleye pollock (corresponding with ages 1–2 years) were calculated separately from adult walleye pollock (> 250 mm). Residuals for individual fish were obtained by subtracting observed weights from bias-corrected weights-at-length that were estimated from regression models. Length-weight residuals from each stratum were aggregated and weighted proportionally to total biomass in each stratum from area-swept expansion of mean bottom-trawl survey catch per unit effort (CPUE; i.e., design-based stratum biomass estimates). Variation in fish condition was evaluated by comparing average length-weight residuals among years. To minimize the influence of unrepresentative samples on indicator calculations, combinations of species, stratum, and year with a sample size <10 were used to fit length-weight regressions, but were excluded from calculating length-weight residuals for both the EBS and NBS. Morphometric condition indicator time series, code for calculating the indicators, and figures showing results for individual species are available through the *akfishcondition* R package and GitHub repository (https://github.com/afsc-gap-products/akfishcondition). + +**Methodological Changes**: In Groundfish Morphometric Condition Indicator contributions to the 2022 Eastern Bering Sea and Aleutians Islands Ecosystem Status Reports, historical stratum-biomass weighted residuals condition indicators were presented alongside condition indicators that were calculated using the R package VAST following methods that were presented for select GOA species during the Spring Preview of Ecological and Economic Conditions (PEEC) in May 2020. The authors noted there were strong correlations between VAST and stratum-biomass weighted condition indicators for most EBS and NBS species (r = 0.79–0.98). The authors received the following feedback about the change from the BSAI Groundfish Plan Team meeting during their November 2022 meeting: + +*"The Team discussed the revised condition indices that now use a different, VAST-based condition index, but felt additional methodology regarding this transition was needed. The Team recommended a short presentation next September to the Team to review the methods and tradeoffs in approaches. The Team encouraged collaboration with the NMFS longline survey team to develop analogous VAST indices."* + +Based on feedback from the Plan Team, staff limitations, and the lack of a clear path to transition condition indicators for longline survey species to VAST, analyses supporting the transition to VAST were not conducted during 2023. Therefore, the 2023 condition indicator was calculated from statum-biomass weighted residuals of length-weight regressions. + +Stratum-biomass weighted residuals for NBS strata are presented for the first time in 2023. NBS length-weight samples were previously pooled across strata to calculate region-wide length-weight residuals because of the lack of samples from regular survey sampling prior to 2017. The authors have opted to present stratum-biomass weighted residuals for the NBS in 2023 because of the accumulation of regular length-weight samples in recent years. + +**Status and Trends**: Fish condition, indicated by length-weight residuals, has varied over time for all species examined in the EBS (Figure 2 & 3). In 2023 a downward trend in condition from 2022 was observed for all species in the EBS, with large walleye pollock (>250 mm), and arrowtooth flounder decreasing since 2019; however, all species were still within one standard deviation of the mean except for large walleye pollock (>250 mm) which was negative but within two standard deviations of the mean (Figure 2). Large walleye pollock (>250 mm) exhibited the second worst condition observed over the full time series with the lowest observed condition occurring in 1999 (Figure 2). In 2019, an upward trend in condition was observed for most species relative to 2017–2018 with positive weighted length-weight residuals relative to historical averages for large walleye pollock (>250 mm), northern rock sole, yellowfin sole, arrowtooth flounder, and Alaska plaice; however, in 2021 condition had a downward trend in most species examined. In 2022 in the EBS, conditions were near the historical mean, or positive for all species examined except for arrowtooth flounder and large walleye pollock (>250 mm). While their conditions were below average, the mean for both groups fell within one standard deviation of the historical mean (Figure 2). + +In the EBS in 2023, condition was negative for large walleye pollock (>250 mm), arrowtooth flounder, and flathead sole across most strata (Figure 3). In 2023, there was a divergence in small walleye pollock (100-250 mm) condition among strata with more positive condition observed on the inner shelf (Stratum 10), and more negative condition observed on the middle shelf (Stratum 30). + +```{r figure 2 grid, include=TRUE, echo=FALSE, fig.height=14,fig.width=12,fig.cap="\\label{fig:figs}Figure 2. Morphometric condition of groundfish species collected during AFSC/RACE GAP standard summer bottom trawl surveys of the eastern Bering Sea shelf (1999 to 2023) based on residuals of length-weight regressions. The dash in the blue boxes denote the mean for that year, the box denotes one standard error, and the lines on the boxes denote two standard errors. Lines on each plot represent the historical mean, dashed lines denote one standard deviation, and dotted lines denote two standard deviations.", message=FALSE, warning=FALSE} +# Set factor levels for plotting order +fig2 <- plot_anomaly_timeseries(x = EBS_INDICATOR$FULL_REGION, + region = "BS", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "rmd") + +fig2_png <- plot_anomaly_timeseries(x = EBS_INDICATOR$FULL_REGION, + region = "BS", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "png") + +print(fig2 + theme_condition_index()) +``` + +```{r figure 2 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", region, "EBS_condition.png"),width=6,height=7,units="in",res=600) +print(fig2_png + theme_blue_strip()) +dev.off() +``` + +```{r figure 3 grid, include=TRUE, echo=FALSE, fig.height=14,fig.width=12,fig.cap="\\label{fig:figs}Figure 3. Length-weight residuals by survey stratum (10-60) for seven eastern Bering Sea shelf groundfish species and age 1–2 walleye pollock (100–250 mm) sampled in the AFSC/RACE GAP standard summer bottom trawl survey, 1999-2023. Length-weight residuals are not weighted by stratum biomass.", message=FALSE, warning=FALSE} +# EBS stratum-biomass weighted residual stratum stacked bar plots +fig3 <- akfishcondition::plot_stratum_stacked_bar(x = EBS_INDICATOR$STRATUM, + region = region, + var_x_name = "year", + var_y_name = "stratum_resid_mean", + y_title = "Length-weight residual (ln(g))", + var_group_name = "stratum", + fill_title = "Stratum", + fill_palette = "BrBG") + + +fig3_png <- akfishcondition::plot_stratum_stacked_bar(x = EBS_INDICATOR$STRATUM, + region = region, + var_x_name = "year", + var_y_name = "stratum_resid_mean", + y_title = "Length-weight residual (ln(g))", + var_group_name = "stratum", + fill_title = "Stratum", + fill_palette = "BrBG") + +print(fig3 + theme_condition_index()+ theme(strip.text = element_text(size = 20))) +``` + +```{r figure 3 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", region, "EBS_SBW_stratum_stacked_bar.png"), width = 6, height = 7, units = "in", res = 600) +print(fig3_png + theme_blue_strip() + theme(legend.position = "right", legend.title = element_text(size = 9))) +dev.off() +``` + +In the NBS in 2023, positive condition was observed for large walleye pollock (>250 mm), which has been increasing since 2021. The remaining species exhibited near-average condition in the NBS in 2023, except for yellowfin sole which exhibited negative condition, and has been declining since 2019 (Figure 4). + +In 2023 large walleye pollock (>250 mm) condition was positive in all NBS strata, whereas condition was previously negative in all strata from 2021-2022 (Figure 5). Pacific cod, small walleye pollock (100-250 mm), Alaska plaice, and yellowfin sole condition have been consistently negative across all strata since 2021, with a notable exception in 2023 of positive condition for Pacific cod in the inner southern NBS shelf, and Alaska plaice in the northern inner NBS shelf and Norton Sound (Figure 5). + +```{r figure 4 grid, include=TRUE, echo=FALSE, fig.height=14,fig.width=12,fig.cap="\\label{fig:figs}Figure 4. Morphometric condition of groundfish species collected during AFSC/RACE GAP standard summer bottom trawl surveys of the northern Bering Sea shelf (2010, 2017, 2019 and 2021-2023) based on residuals of length-weight regressions. The dash in the blue boxes denote the mean for that year, the box denotes one standard error, and the lines on the boxes denote two standard errors. Lines on each plot represent the historical mean, dashed lines denote one standard deviation, and dotted lines denote two standard deviations.", message=FALSE, warning=FALSE} +# Set factor levels for plotting order +fig4 <- plot_anomaly_timeseries(x = NBS_INDICATOR$FULL_REGION, + region = "BS", + fill_color = "#54ADDB", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "rmd") + +fig4_png <- plot_anomaly_timeseries(x = NBS_INDICATOR$FULL_REGION, + region = "BS", + fill_color = "#54ADDB", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "png") + +print(fig4 + theme_condition_index()) + +``` + +```{r figure 4 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", region, "NBS_condition.png"), width=6, height=7, units="in", res=600) +print(fig4_png + theme_blue_strip()) +dev.off() +``` + +```{r figure 5 grid, include=TRUE, echo=FALSE, fig.height=14,fig.width=12,fig.cap="\\label{fig:figs}Figure 5. Length-weight residuals by survey stratum (70, 71 and 81) for four northern Bering Sea shelf groundfish species and age 1–2 walleye pollock (100–250 mm) sampled in the AFSC/RACE GAP standard summer bottom trawl survey during 2010, 2017, 2019 and 2021-2023. Length-weight residuals are not weighted by stratum biomass.", message=FALSE, warning=FALSE} +# NBS stratum-biomass weighted residual stratum stacked bar plots +fig5 <- akfishcondition::plot_stratum_stacked_bar(x = NBS_INDICATOR$STRATUM, + region = region, + var_x_name = "year", + var_y_name = "stratum_resid_mean", + y_title = "Length-weight residual (ln(g))", + var_group_name = "stratum", + fill_title = "Stratum", + fill_palette = "PuBu") + +fig5_png <- akfishcondition::plot_stratum_stacked_bar(x = NBS_INDICATOR$STRATUM, + region = region, + var_x_name = "year", + var_y_name = "stratum_resid_mean", + y_title = "Length-weight residual (ln(g))", + var_group_name = "stratum", + fill_title = "Stratum", + fill_palette = "PuBu") + +print(fig5 + theme_condition_index()+ theme(strip.text = element_text(size = 20))) + +``` + +```{r figure 5 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", region, "NBS_SBW_stratum_stacked_bar.png"), width = 6, height = 7, units = "in", res = 600) +print(fig5_png + theme_blue_strip() + theme(legend.position = "right", legend.title = element_text(size = 9))) +dev.off() +``` + + + +```{r sbw_legacy_plots, include=FALSE, message=FALSE, warning=FALSE} +# # EBS stratum-biomass weighted residual time series anomaly +# ebs_sbw_timeseries <- plot_anomaly_timeseries(x = EBS_INDICATOR$FULL_REGION, +# region = region, +# fill_color = "#0085CA", +# var_y_name = "mean_wt_resid", +# var_y_se_name = "se_wt_resid", +# var_x_name = "year", +# y_title = "Length-weight residual (ln(g))", +# plot_type = "box", +# format_for = "png", +# set_intercept = 0) +# +# png(here::here("plots", region, "EBS_SBW_timeseries.png"), width = 6, height = 7, units = "in", res = 600) +# print(ebs_sbw_timeseries + akfishcondition::theme_blue_strip()) +# dev.off() +# +# # EBS stratum-biomass weighted residual stratum stacked bar plots +# ebs_stacked_bar <- akfishcondition::plot_stratum_stacked_bar(x = EBS_INDICATOR$STRATUM, +# region = region, +# var_x_name = "year", +# var_y_name = "stratum_resid_mean", +# y_title = "Length-weight residual (ln(g))", +# var_group_name = "stratum", +# fill_title = "Stratum", +# fill_palette = "BrBG") +# +# png(here::here("plots", region, "EBS_SBW_stratum_stacked_bar.png"), width = 6, height = 7, units = "in", res = 600) +# print(ebs_stacked_bar + theme_blue_strip() + theme(legend.position = "right", legend.title = element_text(size = 9))) +# dev.off() +# +# # EBS single species stratum plots +# akfishcondition::plot_species_stratum_bar(x = EBS_INDICATOR$STRATUM, +# region = "EBS", +# var_x_name = "year", +# var_y_name = "stratum_resid_mean", +# var_y_se_name = "stratum_resid_se", +# y_title = "Length-weight residual (ln(g))", +# var_group_name = "stratum", +# fill_title = "Stratum", +# fill_palette = "BrBG", +# write_plot = TRUE) +# +# # NBS stratum-biomass weighted residual time series anomaly +# nbs_sbw_timeseries <- plot_anomaly_timeseries(x = NBS_INDICATOR$FULL_REGION, +# region = region, +# fill_color = "#0085CA", +# var_y_name = "mean_wt_resid", +# var_y_se_name = "se_wt_resid", +# var_x_name = "year", +# y_title = "Length-weight residual (ln(g))", +# plot_type = "box", +# format_for = "png", +# set_intercept = 0) +# +# png(here::here("plots", region, "NBS_SBW_timeseries.png"), width = 6, height = 7, units = "in", res = 600) +# print(nbs_sbw_timeseries + akfishcondition::theme_blue_strip()) +# dev.off() +# +# # NBS stratum-biomass weighted residual stratum stacked bar plots +# nbs_stacked_bar <- akfishcondition::plot_stratum_stacked_bar(x = NBS_INDICATOR$STRATUM, +# region = region, +# var_x_name = "year", +# var_y_name = "stratum_resid_mean", +# y_title = "Length-weight residual (ln(g))", +# var_group_name = "stratum", +# fill_title = "Stratum", +# fill_palette = "PuBu") +# +# png(here::here("plots", region, "NBS_SBW_stratum_stacked_bar.png"), width = 6, height = 7, units = "in", res = 600) +# print(nbs_stacked_bar + theme_blue_strip() + theme(legend.position = "right", legend.title = element_text(size = 9))) +# dev.off() +# +# +# # NBS single species stratum plots +# akfishcondition::plot_species_stratum_bar(x = NBS_INDICATOR$STRATUM, +# region = "NBS", +# var_x_name = "year", +# var_y_name = "stratum_resid_mean", +# var_y_se_name = "stratum_resid_se", +# y_title = "Length-weight residual (ln(g))", +# var_group_name = "stratum", +# fill_title = "Stratum", +# fill_palette = "PuBu", +# write_plot = TRUE) + +``` + +**Factors influencing observed trends**: Temperature appears to influence morphological condition of several species in the EBS and NBS, so near-average cold pool extent and water temperatures in 2023 likely played a role in the near-average condition (within 1 S.D. of the mean) for most species in the EBS and NBS. Historically, particularly cold years tend to correspond with negative condition, while particularly warm years tend to correspond to positive condition. For example, water temperatures were particularly cold during the 1999 Bering Sea survey, a year in which negative condition was observed for all species that data were available. In addition, spatiotemporal factor analyses suggest the morphometric condition of age-7 walleye pollock is strongly correlated with cold pool extent in the EBS [@Gruss2021]. In recent years, warm temperatures across the Bering Sea shelf, since the record low seasonal sea ice extent in 2017–2018 and historical cold pool area minimum in 2018 [@Stabeno2019a], may have influenced the positive trend in the condition of several species from 2016 to 2019. However, despite near-average temperature in 2023 large walleye pollock (>250 mm) condition in the EBS was the second lowest recorded over the time series. + +Although warmer temperatures may increase growth rates if there is adequate prey to offset temperature-dependent increases in metabolic demand, growth rates may also decline if prey resources are insufficient to offset temperature-dependent increases in metabolic demand. The influence of temperature on growth rates depends on the physiology of predator species, prey availability, and the adaptive capacity of predators to respond to environmental change through migration, changes in behavior, and acclimatization. For example, elevated temperatures during the 2014–2016 marine heatwave in the Gulf of Alaska led to lower growth rates of Pacific cod and lower condition because available prey resources did not make up for increased metabolic demand [@Barbeaux2020]. + +Other factors that could affect morphological condition include survey timing, stomach fullness, fish movement patterns, sex, and environmental conditions [@Froese2006]. The starting date of annual length-weight data collections has varied from late May to early June and ended in late July-early August in the EBS, and mid-August in the NBS. Although we account for some of this variation by using spatially-varying coefficients in the length-weight relationship, variation in condition could relate to variation in the timing of sample collection within survey strata. Survey timing can be further compounded by seasonal fluctuations in reproductive condition with the buildup and depletion of energy stores [@Wuenschel2019]. Another consideration is that fish weights sampled at sea include gut content weights, so variation in gut fullness may influence weight measurements. Since feeding conditions vary over space and time, prey consumption rates and the proportion of total body weight attributable to gut contents may be an important factor influencing the length-weight residuals. + +Finally, although the condition indicators characterize temporal variation in morphometric condition for important fish species in the EBS and NBS, they do not inform the mechanisms or processes behind the observed patterns. + +**Implications**: Fish morphometric condition can be considered an indicator of ecosystem productivity with implications for fish survival, maturity, and reproduction. For example, in Prince William Sound, the pre-winter condition of herring may determine their overwinter survival [@Paul1999], differences in feeding conditions have been linked to differences in morphometric condition of pink salmon in Prince William Sound [@Boldt2004], variation in morphometric condition has been linked to variation in maturity of sablefish [@Rodgveller2019], and lower morphometric condition of Pacific cod was associated with higher mortality and lower growth rates during the 2014–2016 marine heat wave in the Gulf of Alaska [@Barbeaux2020]. Condition can also be an indicator of stock status relative to carrying capacity because morphometric condition is expected to be high when the stock is at low abundance and low when the stock is at high abundance because of the effects of density-dependent competition [@Haberle2023]. Thus, the condition of EBS and NBS groundfishes may provide insight into ecosystem productivity as well as fish survival, demographic status, and population health. However, survivorship is likely affected by many factors not examined here. + +Another important consideration is that fish condition was computed for all sizes of fishes combined, except in the case of walleye pollock. Examining condition of early juvenile stage fishes not yet recruited to the fishery, or the condition of adult fishes separately, could provide greater insight into the value of length-weight residuals as an indicator of individual health or survivorship [@Froese2006], particularly since juvenile and adult walleye pollock exhibited opposite trends in condition in the EBS this year. + +The near-average condition for most species in 2023 may be related to the near historical average temperatures observed. However, trends in recent years such as prolonged warmer water temperatures following the marine heat wave of 2014-16 [@Bond2015] and reduced sea ice and cold pool areal extent in the eastern Bering Sea [@Stabeno2019a] may affect fish condition in ways that have not yet been determined. Additionally, periods of high fishing mortality that reduce population biomass are likely to increase body condition because of the compensatory alleviation of density-dependent competition [@Haberle2023]. As we continue to add years of length-weight data and expand our knowledge of relationships between condition, growth, production, survival, and the ecosystem, these data may increase our understanding of the health of fish populations in the EBS and NBS. + + +**Research priorities**: Research is being planned and implemented across multiple AFSC programs to explore standardization of statistical methods for calculating condition indicators, and to examine relationships among putatively similar indicators of fish condition (i.e., morphometric, bioenergetic, physiological). Research is underway to evaluate connections between morphometric condition indices, temperature, and density-dependent competition. + + +## References \ No newline at end of file diff --git a/past_years/2023/EBS_GroundfishCondition_2023.docx b/past_years/2023/EBS_GroundfishCondition_2023.docx new file mode 100644 index 0000000..61e94c6 Binary files /dev/null and b/past_years/2023/EBS_GroundfishCondition_2023.docx differ diff --git a/past_years/2023/ESP_Pacific_cod.Rmd b/past_years/2023/ESP_Pacific_cod.Rmd new file mode 100644 index 0000000..99f6d7f --- /dev/null +++ b/past_years/2023/ESP_Pacific_cod.Rmd @@ -0,0 +1,127 @@ +--- +title: ESP Pacific cod morphometric condition +author: +- affiliation: RACE + description: Research Fisheries Biologist + email: Sean.Rohan@NOAA.gov + name: Sean Rohan +output: word_document +fontsize: 12pt +bibliography: EBS_references.bib +csl: fish-and-fisheries.csl +addr: + l1: 7600 Sand Point Way NE + l2: NMFS RACE Division, Groundfish Assessment Program + l3: Seattle, WA 98115 +--- + +# Eastern Bering Sea Pacific cod + +## Adults (≥46 cm FL) + +*Status and Trends:* In 2022, the condition of adult Pacific cod in the EBS was neutral (within one standard deviation of the time series mean), which continues the trend of neutral morphometric condition since 2018. The neutral condition in recent years (2018–2022) represents an increase from the three prior years with negative to neutral condition from 2015–2017. This year, morphometric condition indicator was estimated using VAST [@Thorson2019; @Gruss2020a] instead of the stratum biomass weighted length-weight residual method from previous years. Trends in the historical and new indicator are similar based on the strong correlation between the historical and new indicator (r = 0.87), although there are notable differences for some years. The year with the lowest condition for the old indicator was 1999 (a ‘cold’ year with an early survey start), with an anomaly greater than two standard deviations from the mean. Based on the new VAST relative condition indicator, 1999 was a neutral year and the year with the lowest condition was 2012. Despite these differences, new indicator trends generally match the trend from the old indicator with condition increasing from 1999 to 2003, decreasing from 2003 to 2006, then fluctuating around neutral from 2007 to 2022, aside from negative (>1 standard deviation below the mean) years in 2009, 2012, 2015, and 2017. + +## Juveniles (<46 cm FL) + +*Status and Trends:* In 2022, the condition of juvenile Pacific cod in the EBS was neutral (within one standard deviation of the time series mean), which continues the trend of neutral morphometric condition since 2017. This year, the morphometric condition indicator was estimated using VAST [@Thorson2019; @Gruss2020a] instead of the stratum biomass weighted length-weight residual method from previous years. Trends in the historical and new indicator are similar based on the strong correlation between the historical and new indicator (r = 0.91), although there are noteable differences in some years. Specifically, 2017 was a negative year using the old indicator (~1 standard deviation below the mean) but a neutral year with the new indicator, negative years in 2009 and 2012 are still negative but the anomaly is larger, and the anomaly in 1999 decreased from 3.2 standard deviations below them mean to 1.8 standard deviations below the mean (a ‘cold’ year with an early survey start). Despite these differences, new indicator trends generally match the trend from the old indicator with condition increasing from 1999 to 2004, decreasing from 2005 to 2010, then fluctuating around neutral from 2011 to 2022, aside from a negative year in 2012 and positive year in 2016. + +## Influential Factors + +*Influential Factors:* Many factors contribute to variation in morphometric condition so it is unclear which specific factors contributed to neutral condition of adult Pacific cod in the EBS in 2022. Factors that may contribute to variation in morphometric condition include environmental conditions that affect prey quality and temperature-dependent metabolic rates, survey timing, stomach fullness of individual fish, fish migration patterns, and the distribution of samples within survey strata. Temperature is an important factor that can influence the morphometric condition of Pacific cod by influencing metabolic rates, prey availability, and prey quality. Historically in the eastern Bering Sea (EBS), ‘cold’ years (with a small cold pool) were associated with negative morphometric condition (e.g., 1999, 2012) and warm years (e.g., 2002-2005) were associated with positive morphometric condition. However, during recent (2018–2021) exceptionally warm recent years, the morphometric condition of Pacific cod was neutral for adult and juvenile Pacific cod and this trend continued into the average temperature year in 2022. Temperature can negatively affect growth rates if prey resources are insufficient to make up for increased metabolic demand. In GOA, elevated temperatures during the 2014–2016 marine heatwave in the Gulf of Alaska were associated with lower growth rates of Pacific cod and lower morphometric condition in 2015 (adults and juveniles combined), likely because diminished prey resources during the heatwave were insufficient to make up for increased metabolic demand [@Barbeaux2020]. Additional information about the groundfish morphometric condition indicator and factors that can influence estimates of morphometric condition are described in the EBS Groundfish Morphometric Condition contribution in the 2022 Eastern Bering Sea Ecosystem Status Report (Rohan et al., _In prep_). + + +```{r steup_data, include=FALSE} +library(akfishcondition) + +akfishcondition:::PCOD_ESP$FULL_REGION_EBS %>% + dplyr::filter(common_name == "Pacific cod (juvenile)") %>% + dplyr::mutate(indicator_name = "Summer_Pacific_Cod_Condition_Juvenile_EBS_Survey", + data_value = vast_relative_condition) %>% + dplyr::select(year, indicator_name, data_value) %>% + write.csv(file = here::here("output", "ESP_Summer_Pacific_Cod_Condition_Juvenile_EBS_Survey.csv"), + row.names = FALSE) + +akfishcondition:::PCOD_ESP$FULL_REGION_EBS %>% + dplyr::filter(common_name == "Pacific cod (adult)") %>% + dplyr::mutate(indicator_name = "Summer_Pacific_Cod_Condition_Adult_EBS_Survey", + data_value = vast_relative_condition) %>% + dplyr::select(year, indicator_name, data_value) %>% + write.csv(file = here::here("output", "ESP_Summer_Pacific_Cod_Condition_Adult_EBS_Survey.csv"), + row.names = FALSE) + +print(paste0("Last Update: ", akfishcondition:::EBS_INDICATOR$LAST_UPDATE)) +``` + +```{r ebs_anomaly, include=TRUE, echo=FALSE, fig.height=6,fig.width=12, fig.cap="\\label{fig:figs}Figure 1. VAST relative condition for adult and juvenile Pacific cod collected during AFSC/RACE GAP standard summer bottom trawl surveys of the eastern Bering Sea shelf, 1999 to 2022. The dash in the blue boxes denote the mean for that year, the blue box denotes one standard error, and the lines on the boxes denote two standard error. Horizontal lines on each plot represent the historical mean, dashed lines denote one standard deviation, and dotted lines denote two standard deviations.", message=FALSE, warning=FALSE} +fig1_rmd <- plot_anomaly_timeseries(x = akfishcondition::PCOD_ESP$FULL_REGION_EBS, + region = "BS", + fill_color = "#0085CA", + var_y_name = "vast_relative_condition", + var_y_se_name = "vast_relative_condition_se", + var_x_name = "year", + y_title = "VAST relative condition", + plot_type = "box", + set_intercept = 1, + format_for = "rmd") + +fig1_png <- plot_anomaly_timeseries(x = akfishcondition::PCOD_ESP$FULL_REGION_EBS, + region = "BS", + fill_color = "#0085CA", + var_y_name = "vast_relative_condition", + var_y_se_name = "vast_relative_condition_se", + var_x_name = "year", + y_title = "VAST relative condition", + plot_type = "box", + set_intercept = 1, + format_for = "png") + +print(fig1_rmd + theme_condition_index()) +``` + + +```{r ebs_anomaly_png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", "BS", "ESP_EBS_PCOD_condition.png"),width=6,height=3,units="in",res=600) +print(fig1_png + theme_blue_strip()) +dev.off() +``` + +```{r ebs_vast_sbw, include=TRUE, echo=FALSE, fig.height=6,fig.width=12,fig.cap="\\label{fig:figs}Figure 3. Time series of VAST relative condition and length-weight residual condition anomalies for adult and juvenile Pacific cod in the eastern Bering Sea. Triangles denote the length-weight residual, while points denote the VAST relative condition. Lines represent the historical mean, dashed lines denote one standard deviation, and dotted lines denote two standard deviations. The Pearson correlation coefficient (_r_) is shown at the bottom right of each panel.", message=FALSE, warning=FALSE} +fig2_rmd <- plot_two_timeseries(x_1 = akfishcondition::PCOD_ESP$FULL_REGION_EBS, + x_2 = akfishcondition::PCOD_ESP$FULL_REGION_EBS, + region = "BS", + series_name_1 = "Length-weight residual (Historical)", + series_name_2 = "VAST relative condition (New)", + var_y_name_1 = "mean_wt_resid", + var_x_name_1 = "year", + var_y_name_2 = "vast_relative_condition", + var_x_name_2 = "year", + y_title = "Condition anomaly (Z-score)", + scale_y = TRUE, + year_break = 2020, + fill_colors = c("#001743", "#0085CA"), + format_for = "rmd") + +fig2_png <- plot_two_timeseries(x_1 = akfishcondition::PCOD_ESP$FULL_REGION_EBS, + x_2 = akfishcondition::PCOD_ESP$FULL_REGION_EBS, + region = "BS", + series_name_1 = "Length-weight residual (Historical)", + series_name_2 = "VAST relative condition (New)", + var_y_name_1 = "mean_wt_resid", + var_x_name_1 = "year", + var_y_name_2 = "vast_relative_condition", + var_x_name_2 = "year", + y_title = "Condition anomaly (Z-score)", + scale_y = TRUE, + year_break = 2020, + fill_colors = c("#001743", "#0085CA"), + format_for = "png") +print(fig2_rmd + theme_condition_index() + theme(legend.position = "bottom")) +``` + +```{r ebs_vast_sbw_png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(file = here::here("plots", "BS", "ESP_PCOD_VAST_vs_SBW_time_series.png"), width = 6, height = 3, units = "in", res = 600) +print(fig2_png + theme_blue_strip() + theme(legend.position = "bottom")) +dev.off() +``` + +# References \ No newline at end of file diff --git a/ESP_Pacific_cod.docx b/past_years/2023/ESP_Pacific_cod.docx similarity index 100% rename from ESP_Pacific_cod.docx rename to past_years/2023/ESP_Pacific_cod.docx diff --git a/past_years/2023/ESP_Pacific_cod_EBS.Rmd b/past_years/2023/ESP_Pacific_cod_EBS.Rmd new file mode 100644 index 0000000..548fee0 --- /dev/null +++ b/past_years/2023/ESP_Pacific_cod_EBS.Rmd @@ -0,0 +1,106 @@ +--- +title: Morphometric condition of Eastern Bering Pacific cod for 2023 Ecosystem and Socioeconomic Profile +author: +- affiliation: RACE + description: Research Fisheries Biologist + email: Sean.Rohan@NOAA.gov + name: Sean Rohan +output: word_document +fontsize: 12pt +bibliography: EBS_references.bib +csl: fish-and-fisheries.csl +addr: + l1: 7600 Sand Point Way NE + l2: NMFS RACE Division, Groundfish Assessment Program + l3: Seattle, WA 98115 +--- + +# Eastern Bering Sea Pacific cod + +## Adults (≥46 cm FL) + +*Status and Trends:* In 2023, the condition of adult Pacific cod in the EBS was neutral (within one standard deviation of the time series mean), which continues the trend of neutral morphometric condition observed since 2018. + +## Juveniles (<46 cm FL) + +*Status and Trends:* In 2023, the condition of juvenile Pacific cod in the EBS was neutral (within one standard deviation of the time series mean), which continues the trend of neutral morphometric condition observed since 2017. + +## Influential Factors + +*Influential Factors:* Many factors contribute to variation in morphometric condition so it is unclear which specific factors contributed to neutral condition of adult Pacific cod in the EBS in 2023. Factors that may contribute to variation in morphometric condition include environmental conditions that affect prey quality and temperature-dependent metabolic rates, survey timing, stomach fullness of individual fish, fish migration patterns, and the distribution of samples within survey strata. Temperature is an important factor that can influence the morphometric condition of Pacific cod by influencing metabolic rates, prey availability, and prey quality. Historically in the eastern Bering Sea (EBS), ‘cold’ years (with a small cold pool) were associated with negative morphometric condition (e.g., 1999, 2012) and warm years (e.g., 2002-2005) were associated with positive morphometric condition. However, during exceptionally warm years from 2018–2021, the morphometric condition of Pacific cod was neutral for adult and juvenile Pacific cod and this trend continued into the average temperature years in 2022-2023. Temperature can negatively affect growth rates if prey resources are insufficient to make up for increased metabolic demand. Additional information about the groundfish morphometric condition indicator and factors that can influence estimates of morphometric condition are described in the EBS Groundfish Morphometric Condition contribution in the 2023 Eastern Bering Sea Ecosystem Status Report (Prohaska and Rohan, _In prep_). + +## Implications + +*Implications:* In the Gulf of Alaska, elevated temperatures during the 2014–2016 marine heatwave were associated with lower growth rates of Pacific cod and lower morphometric condition in 2015 (adults and juveniles combined), likely because diminished prey resources during the heatwave were insufficient to make up for increased metabolic demand [@Barbeaux2020]. + + +```{r setup_data, include=FALSE} +library(akfishcondition) + +akfishcondition:::PCOD_ESP$FULL_REGION_EBS %>% + dplyr::filter(common_name == "Pacific cod (juvenile)") %>% + dplyr::mutate(indicator_name = "Summer_Pacific_Cod_Condition_Juvenile_EBS_Survey", + data_value = mean_wt_resid) %>% + dplyr::select(year, indicator_name, data_value) %>% + write.csv(file = here::here("output", "ESP_Summer_Pacific_Cod_Condition_Juvenile_EBS_Survey.csv"), + row.names = FALSE) + +akfishcondition:::PCOD_ESP$FULL_REGION_EBS %>% + dplyr::filter(common_name == "Pacific cod (adult)") %>% + dplyr::mutate(indicator_name = "Summer_Pacific_Cod_Condition_Adult_EBS_Survey", + data_value = mean_wt_resid) %>% + dplyr::select(year, indicator_name, data_value) %>% + write.csv(file = here::here("output", "ESP_Summer_Pacific_Cod_Condition_Adult_EBS_Survey.csv"), + row.names = FALSE) + +print(paste0("Last Update: ", akfishcondition:::EBS_INDICATOR$LAST_UPDATE)) +``` + +```{r ebs_anomaly, include=TRUE, echo=FALSE, fig.height=6,fig.width=12, fig.cap="\\label{fig:figs}Figure 1. VAST relative condition for adult and juvenile Pacific cod collected during AFSC/RACE GAP standard summer bottom trawl surveys of the eastern Bering Sea shelf, 1999 to 2022. The dash in the blue boxes denote the mean for that year, the blue box denotes one standard error, and the lines on the boxes denote two standard error. Horizontal lines on each plot represent the historical mean, dashed lines denote one standard deviation, and dotted lines denote two standard deviations.", message=FALSE, warning=FALSE} +fig1_rmd <- plot_anomaly_timeseries(x = akfishcondition::PCOD_ESP$FULL_REGION_EBS, + region = "BS", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "rmd") + +fig1_png <- plot_anomaly_timeseries(x = akfishcondition::PCOD_ESP$FULL_REGION_EBS, + region = "BS", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "png") + +print(fig1_rmd + theme_condition_index()) +``` + + +```{r ebs_anomaly_png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", "BS", "ESP_EBS_PCOD_condition.png"),width=6,height=3,units="in",res=600) +print(fig1_png + theme_blue_strip()) +dev.off() +``` + + + +```{r make_txt, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} + con <- file(here::here("output", "ESP_EBS_PCOD_Adult.txt")) + + writeLines(text = paste(akfishcondition::PCOD_ESP$FULL_REGION_EBS$mean_wt_resid[akfishcondition::PCOD_ESP$FULL_REGION_EBS$common_name == "Pacific cod (adult)"], collapse = ", "), + con = con) + close(con) + + con <- file(here::here("output", "ESP_EBS_PCOD_Juvenile.txt")) + + writeLines(text = paste(akfishcondition::PCOD_ESP$FULL_REGION_EBS$mean_wt_resid[akfishcondition::PCOD_ESP$FULL_REGION_EBS$common_name == "Pacific cod (juvenile)"], collapse = ", "), + con = con) + close(con) +``` + +# References \ No newline at end of file diff --git a/past_years/2023/ESP_Pacific_cod_GOA.Rmd b/past_years/2023/ESP_Pacific_cod_GOA.Rmd new file mode 100644 index 0000000..300a08d --- /dev/null +++ b/past_years/2023/ESP_Pacific_cod_GOA.Rmd @@ -0,0 +1,105 @@ +--- +title: Morphometric condition of Gulf of Alaska Pacific cod for 2023 Ecosystem and Socioeconomic Profile +author: +- affiliation: RACE + description: Research Fisheries Biologist + email: Sean.Rohan@NOAA.gov + name: Sean Rohan +output: word_document +fontsize: 12pt +bibliography: EBS_references.bib +csl: fish-and-fisheries.csl +addr: + l1: 7600 Sand Point Way NE + l2: NMFS RACE Division, Groundfish Assessment Program + l3: Seattle, WA 98115 +--- + +# Gulf of Alaska Pacific cod + +## Adults (≥42 cm FL) + +*Status and Trends:* In 2023, the condition of adult Pacific cod in the GOA was neutral (within one standard deviation of the time series mean), which continues the trend of neutral morphometric condition observed since 2017. + +## Juveniles (<42 cm FL) + +*Status and Trends:* In 2023, the condition of juvenile Pacific cod in the EBS was neutral (within one standard deviation of the time series mean), which continues the trend of neutral morphometric condition observed since 2019. + +## Influential Factors + +*Influential Factors:* Many factors contribute to variation in morphometric condition so it is unclear which specific factors contributed to neutral condition of adult and juvenile Pacific cod in the GOA in 2023. Factors that may contribute to variation in morphometric condition include environmental conditions that affect prey quality and temperature-dependent metabolic rates, survey timing, stomach fullness of individual fish, fish migration patterns, and the distribution of samples within survey strata. + +## Implications + +In the Gulf of Alaska, elevated temperatures during the 2014–2016 marine heatwave were associated with lower growth rates of Pacific cod and lower morphometric condition in 2015 (adults and juveniles combined), likely because diminished prey resources during the heatwave were insufficient to make up for increased metabolic demand [@Barbeaux2020]. Additional information about the groundfish morphometric condition indicator and factors that can influence estimates of morphometric condition are described in the GOA Groundfish Morphometric Condition contribution in the 2023 Gulf of Alaska Ecosystem Status Report (O'Leary and Rohan, _In prep_). + + +```{r setup_data, include=FALSE} +library(akfishcondition) + +akfishcondition:::PCOD_ESP$FULL_REGION_GOA %>% + dplyr::filter(common_name == "Pacific cod (juvenile)") %>% + dplyr::mutate(indicator_name = "Summer_Pacific_Cod_Condition_Juvenile_GOA_Survey", + data_value = mean_wt_resid) %>% + dplyr::select(year, indicator_name, data_value) %>% + write.csv(file = here::here("output", "ESP_Summer_Pacific_Cod_Condition_Juvenile_GOA_Survey.csv"), + row.names = FALSE) + +akfishcondition:::PCOD_ESP$FULL_REGION_GOA %>% + dplyr::filter(common_name == "Pacific cod (adult)") %>% + dplyr::mutate(indicator_name = "Summer_Pacific_Cod_Condition_Adult_GOA_Survey", + data_value = mean_wt_resid) %>% + dplyr::select(year, indicator_name, data_value) %>% + write.csv(file = here::here("output", "ESP_Summer_Pacific_Cod_Condition_Adult_GOA_Survey.csv"), + row.names = FALSE) + +print(paste0("Last Update: ", akfishcondition:::EBS_INDICATOR$LAST_UPDATE)) +``` + +```{r ebs_anomaly, include=TRUE, echo=FALSE, fig.height=6,fig.width=12, fig.cap="\\label{fig:figs}Figure 1. Stratum-biomass weighted morphometric condition of adult and juvenile Pacific cod collected during AFSC/RACE GAP standard summer bottom trawl surveys of the eastern Bering Sea shelf, 1984 to 2023. The dash in the blue boxes denote the mean for that year, the blue box denotes one standard error, and the lines on the boxes denote two standard error. Horizontal lines on each plot represent the historical mean, dashed lines denote one standard deviation, and dotted lines denote two standard deviations.", message=FALSE, warning=FALSE} +fig1_rmd <- plot_anomaly_timeseries(x = akfishcondition::PCOD_ESP$FULL_REGION_GOA, + region = "GOA", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "rmd") + +fig1_png <- plot_anomaly_timeseries(x = akfishcondition::PCOD_ESP$FULL_REGION_GOA, + region = "GOA", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "png") + +print(fig1_rmd + theme_condition_index()) +``` + + +```{r ebs_anomaly_png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", "GOA", "ESP_GOA_PCOD_condition.png"),width=6,height=3,units="in",res=600) +print(fig1_png + theme_blue_strip()) +dev.off() +``` + + +```{r make_txt, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} + con <- file(here::here("output", "ESP_GOA_PCOD_Adult.txt")) + + writeLines(text = paste(akfishcondition::PCOD_ESP$FULL_REGION_GOA$mean_wt_resid[akfishcondition::PCOD_ESP$FULL_REGION_GOA$common_name == "Pacific cod (adult)"], collapse = ", "), + con = con) + close(con) + + con <- file(here::here("output", "ESP_GOA_PCOD_Juvenile.txt")) + + writeLines(text = paste(akfishcondition::PCOD_ESP$FULL_REGION_GOA$mean_wt_resid[akfishcondition::PCOD_ESP$FULL_REGION_GOA$common_name == "Pacific cod (juvenile)"], collapse = ", "), + con = con) + close(con) +``` + +# References \ No newline at end of file diff --git a/ESP_Pacific_cod_GOA.docx b/past_years/2023/ESP_Pacific_cod_GOA.docx similarity index 100% rename from ESP_Pacific_cod_GOA.docx rename to past_years/2023/ESP_Pacific_cod_GOA.docx diff --git a/past_years/2023/GOA_GroundfishCondition_2023.Rmd b/past_years/2023/GOA_GroundfishCondition_2023.Rmd new file mode 100644 index 0000000..8a7b6e4 --- /dev/null +++ b/past_years/2023/GOA_GroundfishCondition_2023.Rmd @@ -0,0 +1,186 @@ +--- +title: Gulf of Alaska Groundfish Condition +author: +- affiliation: RACE + description: Research Fisheries Biologist + email: cecilia.oleary@noaa.gov + name: Cecilia O'Leary +output: word_document +fontsize: 12pt +addr: + l1: 7600 Sand Point Way NE + l2: NMFS RACE Division, Groundfish Assessment Program + l3: Seattle, WA 98115 +--- + +```{r setup, include=FALSE} + +library(akfishcondition) +library(knitr) + +region <- "GOA" + +# Unzip map files +unzip(system.file("extdata/2022_ESR.zip", package = "akfishcondition")) + +dir.create(here::here("plots", region)) + + +# Load data +# goa_ann_mean_resid_df <- GOA_INDICATOR$FULL_REGION +# goa_stratum_resids_df <- GOA_INDICATOR$STRATUM +``` + +Contributed by Cecilia O'Leary^1^ and Sean Rohan^1^ + +^1^ Resource Assessment and Conservation Engineering Division, Groundfish Assessment Program, Alaska Fisheries Science Center, National Marine Fisheries Service, NOAA, Seattle, WA +**Contact**: +**Last updated**: September 2023 + +**Description of Indicator**: Length-weight residuals represent how heavy a fish is per unit body length and are an indicator of somatic growth variability (Brodeur et al., 2004). Therefore, length-weight residuals can be considered indicators of prey availability, growth, general health, and habitat condition (Blackwell et al., 2000; Froese, 2006). Positive length-weight residuals indicate better condition (i.e., heavier per unit length) and negative residuals indicate poorer condition (i.e., lighter per unit length) (Froese, 2006). Fish condition calculated in this way reflects realized outcomes of intrinsic and extrinsic processes that affect fish growth which can have implications for biological productivity through direct effects on growth and indirect effects on demographic processes such as, reproduction, and mortality (e.g., Rodgveller (2019); Barbeaux et al. (2020)). + +```{r map, include = TRUE, out.width = "200%", fig.cap = "\\label{fig:figs}Figure 1. National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP) Gulf of Alaska summer bottom trawl survey area with International North Pacific Fisheries Commission (INPFC) statistical fishing strata delineated by the red lines.", echo = FALSE} +knitr::include_graphics("MapGOA.png") +``` + +The groundfish morphometric condition indicator is calculated from paired fork lengths (mm) and weights (g) of individual fishes that were collected during bottom trawl survey of the Gulf of Alaska (GOA) which were conducted by the Alaska Fisheries Science Center biennial Resource Assessment and Conservation Engineering (AFSC/RACE) - Groundfish Assessment Program’s (GAP). Fish condition analyses were applied to walleye pollock (*Gadus chalcogrammus*), Pacific cod (*Gadus macrocephalus*), arrowtooth flounder (*Atheresthes stomias*), southern rock sole (*Lepidopsetta bilineata*), northern rockfish (*Sebastes polyspinis*), Pacific ocean perch (*Sebastes alutus*), dusky rockfish (*Sebastes variabilis*), shortraker rockfish (*Sebastes borealis*), rougheye rockfish (*Sebastes aleutianus*),sharpchin rockfish (*Sebastes zacentrus*),flathead sole (*Hippoglossoides elassodon*),dover sole (*Microstomus pacificus*), and rex sole(*Glyptocephalus zachirus*) collected in trawls with satisfactory performance at standard survey stations. Data were combined in the former International North Pacific Fisheries Commission (INPFC) strata; Shumagin, Chirikof, Kodiak, Yakutat and Southeast (Figure 1). + +To calculate indicators, length-weight relationships were estimated from linear regression models based on a log-transformation of the exponential growth relationship, W = aLb, where W is weight (g) and L is fork length (mm) for all areas for the period 1984–2023. UniqueA unique intercepts and (a) and slopes (b) wereas estimated for each survey stratum, sex, and interaction between stratum and sex to account for sexual dimorphism and spatial-temporal variation in growth and bottom trawl survey sampling. Length-weight relationships for 100–250 mm fork length walleye pollock (corresponding with ages 1–2 years) were calculated separately from adult walleye pollock (> 250 mm). Residuals for individual fish were obtained by subtracting observed weights from bias-corrected weights-at-length that were estimated from regression models. Individual length-weight residuals were aggregated and averaged for each stratum, weighted based on the proportion to total biomass in each stratum from area-swept expansion of bottom-trawl survey catch per unit effort (CPUE; i.e., design-based stratum biomass estimates). Variation in fish condition was evaluated by comparing average length-weight residuals among years. To minimize the influence of unrepresentative samples on indicator calculations, combinations of species, stratum, and year with a sample size < 10 were used to fit length-weight regressions but were excluded from calculating length-weight residuals. Morphometric condition indicator time series, code for calculating the indicators, and figures showing results for individual species are available through the akfishcondition R package and GitHub repository (https://github.com/afsc-gap-products/akfishcondition). + +**Methodological changes**: In 2022, historical stratum-biomass weighted residuals condition indicators were presented alongside condition indicators that were calculated using the R package VAST following methods that were presented for select GOA species during the Spring Preview of Ecological and Economic Conditions in May 2020. The authors noted there were strong correlations between VAST and stratum biomass weighted condition indicators for most species (r = 0.79–0.98). The authors received the following feedback about the change from the BSAI Groundfish Plan Team meeting during their November 2022 meeting: +"The Team discussed the revised condition indices that now use a different, VAST-based condition index, but felt additional methodology regarding this transition was needed. The Team recommended a short presentation next September to the Team to review the methods and tradeoffs in approaches. The Team encouraged collaboration with the NMFS longline survey team to develop analogous VAST indices." +Based on feedback from the Plan Team, staff limitations, and the lack of a clear path to transition condition indicators for longline species to VAST, analyses supporting the transition to VAST were not conducted during 2023. Therefore, the 2023 condition indicator was calculated from stratum-biomass weighted residuals of length-weight regressions. + +**Status and Trends**: Residual body condition varied among survey years for all species considered (Figure 2). Fish condition indicators for six of the seven species were below average in 2023, but with the same condition or reduction in magnitude for most species in 2023 relative to 2021. The exception to this was rougheye rockfish, which had an above average fish condition in 2023 and improved condition from 2021. Residual body condition for pollock, shortraker rockfish, flathead sole, dover sole, rex sole, and arrowtooth flounder remained relatively constant compared to 2021. Residual body condition for Pacific cod, dusky rockfish, and sharpchin rockfish declined in 2023. Southern rock sole residual body condition improved over the last four years, but then declined again in 2023 and the final two years remained a constant below average condition. Northern rock sole had the opposite trend, where it declinded in residual body condition over the last four years and then impvoed again in 2023, but remains below average condition. Residual body condition for northern rockfish, Pacific Ocean Perch, and blackspotted rockfish also improved, but are still below average. Rougheye rockfish residual body condition improved back to above average condition. + +```{r figure 2a grid, include = TRUE, echo = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 2A. Biomass-weighted residual body condition index across survey years (1984-2023) for sixteen Gulf of Alaska groundfish species collected on the National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP) standard summer bottom trawl survey. Filled bars denote weighted length-weight residuals, error bars denote two standard errors.", message = FALSE, warning = FALSE} + +group_1 <- (set_plot_order(common_name = unique(GOA_INDICATOR$FULL_REGION$common_name), region = "GOA") |> + factor() |> + levels() |> + as.character() |> + tolower())[1:8] +group_2 <- (set_plot_order(common_name = unique(GOA_INDICATOR$FULL_REGION$common_name), region = "GOA") |> + factor() |> + levels() |> + as.character() |> + tolower())[9:16] + +fig2a <- plot_anomaly_timeseries(x = dplyr::filter(GOA_INDICATOR$FULL_REGION, + tolower(common_name) %in% group_1), + region = "GOA", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "rmd") + +fig2a_png <- plot_anomaly_timeseries(x = dplyr::filter(GOA_INDICATOR$FULL_REGION, tolower(common_name) %in% group_1), + region = "GOA", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "png") + +print(fig2a + theme_condition_index()) +``` + +```{r figure 2b grid, include = TRUE, echo = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 2B. Biomass-weighted residual body condition index across survey years (1984-2023) for sixteen Gulf of Alaska groundfish species collected on the National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP) standard summer bottom trawl survey. Filled bars denote weighted length-weight residuals, error bars denote two standard errors.", message = FALSE, warning = FALSE} +fig2b <- plot_anomaly_timeseries(x = dplyr::filter(GOA_INDICATOR$FULL_REGION, tolower(common_name) %in% group_2), + region = "GOA", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "rmd") + +fig2b_png <- plot_anomaly_timeseries(x = dplyr::filter(GOA_INDICATOR$FULL_REGION, tolower(common_name) %in% group_2), + region = "GOA", + fill_color = "#0085CA", + var_y_name = "mean_wt_resid", + var_y_se_name = "se_wt_resid", + var_x_name = "year", + y_title = "Length-weight residual (ln(g))", + plot_type = "box", + format_for = "png") + +print(fig2b + theme_condition_index()) +``` + +```{r figure 2 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", region, "GOAbyyear_1.png"),width=6,height=7,units="in",res=300) +print(fig2a_png + theme_blue_strip()) +dev.off() + +png(here::here("plots", region, "GOAbyyear_2.png"),width=6,height=7,units="in",res=300) +print(fig2b_png + theme_blue_strip()) +dev.off() +``` + +The general patterns of above and below average residual body condition index across recent survey years for the GOA as described above were also apparent in the spatial condition indicators across INPFC strata (Figure 3). The relative contribution of stratum-specific residual body condition to the overall trends (indicated by the height of each colored bar segment) does not demonstrate a clear pattern. Although, for many species, the direction of residual body condition (positive or negative) was synchronous among strata within years. For example, residual body condition for pollock , Pacific Ocean Perch, and dusky rockfish in Southeast was positive while the majority of other locations for other fish trended negative. Exceptions include rougheye rockfish in Chrikof and Kodiak and rex sole in Kodiak. While Pacific cod residuals trended negative again, residual body condition in the Kodiak strata remained positive. Patterns of fish distribution are also apparent in the stratum condition indexes. For example, northern rockfish have primarily been collected from the Shumagin and Chirikof strata in recent surveys. + +```{r figure 3a grid, include = TRUE, echo = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 3A. Residual body condition index for sixteen Gulf of Alaska groundfish species collected on the National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP) standard summer bottom trawl survey (1984--2023) grouped by International North Pacific Fisheries Commission (INPFC) statistical sampling strata.", message = FALSE, warning = FALSE} +fig3a <- plot_stratum_stacked_bar(x = dplyr::filter(GOA_INDICATOR$STRATUM, + tolower(common_name) %in% group_1), + region = "GOA", + fill_palette = "BrBG", + var_y_name = "stratum_resid_mean", + var_x_name = "year", + var_group_name = "inpfc_stratum", + fill_title = "Stratum", + y_title = "Length-weight residual (ln(g))") + +fig3b <- plot_stratum_stacked_bar(x = dplyr::filter(GOA_INDICATOR$STRATUM, + tolower(common_name) %in% group_2), + region = "GOA", + fill_palette = "BrBG", + var_y_name = "stratum_resid_mean", + var_x_name = "year", + var_group_name = "inpfc_stratum", + fill_title = "Stratum", + y_title = "Length-weight residual (ln(g))") + +print(fig3a + theme_condition_index()) + +``` + +```{r figure 3b grid, include = TRUE, echo = FALSE, fig.height = 14, fig.width = 12, fig.cap = "\\label{fig:figs}Figure 3B. Residual body condition index for sixteen Gulf of Alaska groundfish species collected on the National Marine Fisheries Service (NMFS) Alaska Fisheries Science Center Resource Assessment and Conservation Engineering (AFSC/RACE) Groundfish Assessment Program (GAP) standard summer bottom trawl survey (1984--2023) grouped by International North Pacific Fisheries Commission (INPFC) statistical sampling strata.", message = FALSE, warning = FALSE} +print(fig3b + theme_condition_index()) +``` + +```{r figure 3 grid png, include=FALSE, echo=FALSE, message=FALSE, warning=FALSE} +png(here::here("plots", region, "GOA_condition_by_stratum_1.png") ,width=6,height=7,units="in",res=300) +print(fig3a + theme_blue_strip()) +dev.off() + +png(here::here("plots", region, "GOA_condition_by_stratum_2.png") ,width=6,height=7,units="in",res=300) +print(fig3b + theme_blue_strip()) +dev.off() + +plot_species_stratum_bar( + x = GOA_INDICATOR$STRATUM, + region = "GOA", + var_x_name = "year", + var_y_name = "stratum_resid_mean", + var_y_se_name = "stratum_resid_se", + y_title = "Length-weight residual (ln(g))", + var_group_name = "inpfc_stratum", + fill_title = "Stratum", + fill_palette = "BrBG", + write_plot = TRUE) +``` + +**Factors causing observed trends**: Factors that could affect residual fish body condition presented here include temperature, trawl survey timing, stomach fullness, movement in or out of the survey area, or variable somatic growth. Following an unprecedented warming event from 2014--2016 (Bond et al., 2015; Stabeno et al., 2019; Barbeaux et al., 2020), there has been a general trend of warming ocean temperatures in the survey area and sea surface temperature anomaly data continue to reflect temperatures above average historical conditions through 2023 (NOAA 2023); these warmer temperatures could be affecting fish growth conditions in this region. Changing ocean conditions along with normal patterns of movement can cause the proportion of the population resident in the sampling area during the annual bottom trawl survey to vary. Recorded changes attributed to the marine heatwave included species abundances, sizes, growth rates, weight/body condition, reproductive success, and species composition (Suryan et al., 2021). Warmer ocean temperatures can lead to lower energy (leaner) prey, increased metabolic needs of younger fish, and therefore slower growth for juveniles, as observed in Pacific cod (Barbeaux et al., 2020). Additionally, spatial and temporal trends in fish growth over the season become confounded with survey progress since the first length-weight data are generally collected in late May and the bottom trawl survey is conducted throughout the summer months moving from west to east. In addition, spatial variability in residual condition may also reflect local environmental features which can influence growth and prey availability in the areas surveyed (e.g., local differences in IPHC region average cross-shelf transport of heat via eddies reported this year; NOAA 2023). + +**Implications**: Variations in body condition likely have implications for fish survival. The condition of GOA groundfish may contribute to survival and recruitment. As future years are added to the time series, the relationship between length-weight residuals and subsequent survival will be examined further. It is important that residual body condition for most species in these analyses was computed for all sizes and sexes combined. Requirements for growth and survivorship differ for different fish life stages and some species have sexually dimorphic or even regional growth patterns. It may be more informative to examine life-stage (e.g., early juvenile, subadult, and adult phases) and sex-specific body condition in the future. + +The trend toward below average body condition for many GOA species over the last four to five RACE/AFSC GAP bottom trawl surveys is a potential cause for concern. It could indicate poor overwinter survival or may reflect the influence of locally changing environmental conditions depressing fish growth, local production, or survivorship. Indications are that the Warm Blob (Bond et al., 2015; Stabeno et al., 2019) has been followed by subsequent years with elevated water temperatures (e.g., Barbeaux et al., 2020; NOAA, 2021, 2023) which may be related to changes in fish condition in the species examined. It should be noted that for while many GOA species fish condition remained below average this year, many of the species (with the exception of southern rock sole, dusky rockfish, Pacific cod, and sharpchin rockfish) fish condition improved relative to 2021. As we continue to add years of fish condition to the record and expand on our knowledge of the relationships between condition, growth, production, and survival, we hope to gain more insight into the overall health of fish populations in the GOA. + +**Research priorities**: Research is being planned and implemented across multiple AFSC programs to explore standardization of statistical methods for calculating condition indicators, and to examine relationships among putatively similar indicators of fish condition (i.e., morphometric, bioenergetic, physiological). Finally, research is underway towe plan to evaluate connections betweenexplore variation in condition indices andbetween life history stages alongside effects if density dependentce and temperature dependent growth climate change effects (Bolin et al., 2021; Oke et al., 2022). \ No newline at end of file diff --git a/plots/AI/AI_SBW_stratum_stacked_bar.png b/plots/AI/AI_SBW_stratum_stacked_bar.png deleted file mode 100644 index 899280f..0000000 Binary files a/plots/AI/AI_SBW_stratum_stacked_bar.png and /dev/null differ diff --git a/plots/AI/AI_SBW_timeseries.png b/plots/AI/AI_SBW_timeseries.png deleted file mode 100644 index bb01e3b..0000000 Binary files a/plots/AI/AI_SBW_timeseries.png and /dev/null differ diff --git a/plots/AI/AI_VAST_SBW_timeseries.png b/plots/AI/AI_VAST_SBW_timeseries.png deleted file mode 100644 index 9cba728..0000000 Binary files a/plots/AI/AI_VAST_SBW_timeseries.png and /dev/null differ diff --git a/plots/AI/AI_VAST_VS_SBW_scatterplot.png b/plots/AI/AI_VAST_VS_SBW_scatterplot.png deleted file mode 100644 index ff8c421..0000000 Binary files a/plots/AI/AI_VAST_VS_SBW_scatterplot.png and /dev/null differ diff --git a/plots/AI/AI_condition.png b/plots/AI/AI_condition.png index 3b4e703..f98b01b 100644 Binary files a/plots/AI/AI_condition.png and b/plots/AI/AI_condition.png differ diff --git a/plots/AI/AI_condition_by_stratum.png b/plots/AI/AI_condition_by_stratum.png new file mode 100644 index 0000000..9692dee Binary files /dev/null and b/plots/AI/AI_condition_by_stratum.png differ diff --git a/plots/BS/EBS_NBS_survey_area.png b/plots/BS/EBS_NBS_survey_area.png index fffc8c4..e4e5b49 100644 Binary files a/plots/BS/EBS_NBS_survey_area.png and b/plots/BS/EBS_NBS_survey_area.png differ diff --git a/plots/BS/EBS_SBW_stratum_stacked_bar.png b/plots/BS/EBS_SBW_stratum_stacked_bar.png index 4bdb13f..f24dc0a 100644 Binary files a/plots/BS/EBS_SBW_stratum_stacked_bar.png and b/plots/BS/EBS_SBW_stratum_stacked_bar.png differ diff --git a/plots/BS/EBS_condition.png b/plots/BS/EBS_condition.png index 6d59871..efb52cd 100644 Binary files a/plots/BS/EBS_condition.png and b/plots/BS/EBS_condition.png differ diff --git a/plots/BS/ESP_EBS_PCOD_condition.png b/plots/BS/ESP_EBS_PCOD_condition.png new file mode 100644 index 0000000..78a783c Binary files /dev/null and b/plots/BS/ESP_EBS_PCOD_condition.png differ diff --git a/plots/BS/NBS_SBW_stratum_stacked_bar.png b/plots/BS/NBS_SBW_stratum_stacked_bar.png index f3410fd..d82e140 100644 Binary files a/plots/BS/NBS_SBW_stratum_stacked_bar.png and b/plots/BS/NBS_SBW_stratum_stacked_bar.png differ diff --git a/plots/BS/NBS_condition.png b/plots/BS/NBS_condition.png index 4aca714..570a427 100644 Binary files a/plots/BS/NBS_condition.png and b/plots/BS/NBS_condition.png differ diff --git a/plots/GOA/Dover sole/GOA_STRATUM_CONDITION_Dover sole.png b/plots/GOA/Dover sole/GOA_STRATUM_CONDITION_Dover sole.png index 80a3d66..d1fce15 100644 Binary files a/plots/GOA/Dover sole/GOA_STRATUM_CONDITION_Dover sole.png and b/plots/GOA/Dover sole/GOA_STRATUM_CONDITION_Dover sole.png differ diff --git a/plots/GOA/ESP_GOA_PCOD_condition.png b/plots/GOA/ESP_GOA_PCOD_condition.png new file mode 100644 index 0000000..fdf4aa8 Binary files /dev/null and b/plots/GOA/ESP_GOA_PCOD_condition.png differ diff --git a/plots/GOA/GOA_condition_by_stratum_1.png b/plots/GOA/GOA_condition_by_stratum_1.png index 4ff7f8b..bc2b71b 100644 Binary files a/plots/GOA/GOA_condition_by_stratum_1.png and b/plots/GOA/GOA_condition_by_stratum_1.png differ diff --git a/plots/GOA/GOA_condition_by_stratum_2.png b/plots/GOA/GOA_condition_by_stratum_2.png index 1289fe3..853c0c7 100644 Binary files a/plots/GOA/GOA_condition_by_stratum_2.png and b/plots/GOA/GOA_condition_by_stratum_2.png differ diff --git a/plots/GOA/GOAbyyear_1.png b/plots/GOA/GOAbyyear_1.png index 4341a35..9faa45c 100644 Binary files a/plots/GOA/GOAbyyear_1.png and b/plots/GOA/GOAbyyear_1.png differ diff --git a/plots/GOA/GOAbyyear_2.png b/plots/GOA/GOAbyyear_2.png index 5a1275e..26a3fc2 100644 Binary files a/plots/GOA/GOAbyyear_2.png and b/plots/GOA/GOAbyyear_2.png differ diff --git a/plots/GOA/Pacific cod/GOA_STRATUM_CONDITION_Pacific cod.png b/plots/GOA/Pacific cod/GOA_STRATUM_CONDITION_Pacific cod.png index c51210e..ce0fbbf 100644 Binary files a/plots/GOA/Pacific cod/GOA_STRATUM_CONDITION_Pacific cod.png and b/plots/GOA/Pacific cod/GOA_STRATUM_CONDITION_Pacific cod.png differ diff --git a/plots/GOA/Pacific ocean perch/GOA_STRATUM_CONDITION_Pacific ocean perch.png b/plots/GOA/Pacific ocean perch/GOA_STRATUM_CONDITION_Pacific ocean perch.png index 53d2626..499957c 100644 Binary files a/plots/GOA/Pacific ocean perch/GOA_STRATUM_CONDITION_Pacific ocean perch.png and b/plots/GOA/Pacific ocean perch/GOA_STRATUM_CONDITION_Pacific ocean perch.png differ diff --git a/plots/GOA/arrowtooth flounder/GOA_STRATUM_CONDITION_arrowtooth flounder.png b/plots/GOA/arrowtooth flounder/GOA_STRATUM_CONDITION_arrowtooth flounder.png index 7db66c7..275b08e 100644 Binary files a/plots/GOA/arrowtooth flounder/GOA_STRATUM_CONDITION_arrowtooth flounder.png and b/plots/GOA/arrowtooth flounder/GOA_STRATUM_CONDITION_arrowtooth flounder.png differ diff --git a/plots/GOA/blackspotted rockfish/GOA_STRATUM_CONDITION_blackspotted rockfish.png b/plots/GOA/blackspotted rockfish/GOA_STRATUM_CONDITION_blackspotted rockfish.png index 76c0eeb..21e9fe6 100644 Binary files a/plots/GOA/blackspotted rockfish/GOA_STRATUM_CONDITION_blackspotted rockfish.png and b/plots/GOA/blackspotted rockfish/GOA_STRATUM_CONDITION_blackspotted rockfish.png differ diff --git a/plots/GOA/dusky rockfish/GOA_STRATUM_CONDITION_dusky rockfish.png b/plots/GOA/dusky rockfish/GOA_STRATUM_CONDITION_dusky rockfish.png index 0bd12ab..eb23825 100644 Binary files a/plots/GOA/dusky rockfish/GOA_STRATUM_CONDITION_dusky rockfish.png and b/plots/GOA/dusky rockfish/GOA_STRATUM_CONDITION_dusky rockfish.png differ diff --git a/plots/GOA/flathead sole/GOA_STRATUM_CONDITION_flathead sole.png b/plots/GOA/flathead sole/GOA_STRATUM_CONDITION_flathead sole.png index 3e31fb7..cfcc0e7 100644 Binary files a/plots/GOA/flathead sole/GOA_STRATUM_CONDITION_flathead sole.png and b/plots/GOA/flathead sole/GOA_STRATUM_CONDITION_flathead sole.png differ diff --git a/plots/GOA/northern rock sole/GOA_STRATUM_CONDITION_northern rock sole.png b/plots/GOA/northern rock sole/GOA_STRATUM_CONDITION_northern rock sole.png index 8b2b067..0e0b16e 100644 Binary files a/plots/GOA/northern rock sole/GOA_STRATUM_CONDITION_northern rock sole.png and b/plots/GOA/northern rock sole/GOA_STRATUM_CONDITION_northern rock sole.png differ diff --git a/plots/GOA/northern rockfish/GOA_STRATUM_CONDITION_northern rockfish.png b/plots/GOA/northern rockfish/GOA_STRATUM_CONDITION_northern rockfish.png index 99a6609..e875eb4 100644 Binary files a/plots/GOA/northern rockfish/GOA_STRATUM_CONDITION_northern rockfish.png and b/plots/GOA/northern rockfish/GOA_STRATUM_CONDITION_northern rockfish.png differ diff --git a/plots/GOA/rex sole/GOA_STRATUM_CONDITION_rex sole.png b/plots/GOA/rex sole/GOA_STRATUM_CONDITION_rex sole.png index 9b21a6d..715589f 100644 Binary files a/plots/GOA/rex sole/GOA_STRATUM_CONDITION_rex sole.png and b/plots/GOA/rex sole/GOA_STRATUM_CONDITION_rex sole.png differ diff --git a/plots/GOA/rougheye rockfish/GOA_STRATUM_CONDITION_rougheye rockfish.png b/plots/GOA/rougheye rockfish/GOA_STRATUM_CONDITION_rougheye rockfish.png index dcd9a1a..2e80e15 100644 Binary files a/plots/GOA/rougheye rockfish/GOA_STRATUM_CONDITION_rougheye rockfish.png and b/plots/GOA/rougheye rockfish/GOA_STRATUM_CONDITION_rougheye rockfish.png differ diff --git a/plots/GOA/sharpchin rockfish/GOA_STRATUM_CONDITION_sharpchin rockfish.png b/plots/GOA/sharpchin rockfish/GOA_STRATUM_CONDITION_sharpchin rockfish.png index 737b80a..eb4b6f2 100644 Binary files a/plots/GOA/sharpchin rockfish/GOA_STRATUM_CONDITION_sharpchin rockfish.png and b/plots/GOA/sharpchin rockfish/GOA_STRATUM_CONDITION_sharpchin rockfish.png differ diff --git a/plots/GOA/shortraker rockfish/GOA_STRATUM_CONDITION_shortraker rockfish.png b/plots/GOA/shortraker rockfish/GOA_STRATUM_CONDITION_shortraker rockfish.png index 1d707e0..14de08d 100644 Binary files a/plots/GOA/shortraker rockfish/GOA_STRATUM_CONDITION_shortraker rockfish.png and b/plots/GOA/shortraker rockfish/GOA_STRATUM_CONDITION_shortraker rockfish.png differ diff --git a/plots/GOA/southern rock sole/GOA_STRATUM_CONDITION_southern rock sole.png b/plots/GOA/southern rock sole/GOA_STRATUM_CONDITION_southern rock sole.png index 31cbdf2..585b44d 100644 Binary files a/plots/GOA/southern rock sole/GOA_STRATUM_CONDITION_southern rock sole.png and b/plots/GOA/southern rock sole/GOA_STRATUM_CONDITION_southern rock sole.png differ diff --git a/plots/GOA/walleye pollock (100-250 mm)/GOA_STRATUM_CONDITION_walleye pollock (100-250 mm).png b/plots/GOA/walleye pollock (100-250 mm)/GOA_STRATUM_CONDITION_walleye pollock (100-250 mm).png index 8c3fad8..4d05634 100644 Binary files a/plots/GOA/walleye pollock (100-250 mm)/GOA_STRATUM_CONDITION_walleye pollock (100-250 mm).png and b/plots/GOA/walleye pollock (100-250 mm)/GOA_STRATUM_CONDITION_walleye pollock (100-250 mm).png differ diff --git a/plots/GOA/walleye pollock (gt250 mm)/GOA_STRATUM_CONDITION_walleye pollock (gt250 mm).png b/plots/GOA/walleye pollock (gt250 mm)/GOA_STRATUM_CONDITION_walleye pollock (gt250 mm).png index 1dc2e63..f4be837 100644 Binary files a/plots/GOA/walleye pollock (gt250 mm)/GOA_STRATUM_CONDITION_walleye pollock (gt250 mm).png and b/plots/GOA/walleye pollock (gt250 mm)/GOA_STRATUM_CONDITION_walleye pollock (gt250 mm).png differ