-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathbasa_report.Rmd
183 lines (137 loc) · 15.5 KB
/
basa_report.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
---
title: "2023 Bayesian Age-structure Stock Assessment (BASA) Results for Prince William Sound (PWS) herring"
author:
- Joshua A. Zahner
- Trevor A. Branch
date: "`r format(Sys.time(), '%B %d, %Y')`"
geometry: margin=1in
fontsize: 11pt
output:
pdf_document:
keep_tex: yes
fig_caption: yes
template: default
---
# Executive Summary
The median expected spawning biomass of PWS herring in 2024 was estimated at approximately 27,050 metric tons, above the minimum threshold required for the opening of existing herring fisheries (19,958 metric tons). Taking into account uncertainity in this estimate, there is an approximately 11% probability that the true biomass of the PWS herring population is below this lower cutoff. The 2024 biomass estimate is the highest since 1993, and continues the recent trend of biomass growth observed since 2018. The estimated age-composition for 2023 shows a modest proportion of age-7 individuals, continuing to support the estimate of a strong 2016 cohort, as well as relatively large proportions of age-3 and age-4 fish, supporting a large 2020 and 2021 year-class. Based on recent trends in both the survey data and model estimates of biomass, the PWS herring stock appears to be recovering towards its mid-2000s level. While still a long way from the biomass levels sustained prior to the 1993 population crash, this is a welcome sign after nearly a decade of further biomass decline.
# Background
Before 2014, the Alaska Department of Fish and Game (ADF&G) ran an Excel-based age structured assessment (ASA) model to forecast PWS herring biomass for input into harvest control rule. The harvest control rule has a minimum biomass threshold at 19,958 metric tons, which is equivalent to 25% of the unfished biomass under equilibrium determined from simulations (Funk and Rowell 1995). When forecasted biomass is between 19,958 and 38,555 metric tons (22,000-42,500 short tons), the control rules scales the annual harvest rate from 0-20% (Botz et al. 2011). These reference points were last revised by the Alaska Board of Fisheries in 1994.
Since 2014, the ASA has been expanded to include a Bayesian formulation (BASA) that inherently weights the input data sources based on statistical probability distributions, and estimates uncertainty through the sampling of Bayesian posteriors (Muradian et al. 2017). Muradian et al. (2017) first demonstrated BASA as a more robust model to the previous ASA. Since then, BASA has been used in various studies to evaluate which historical input data were the most informative given the trade-off between information gain and cost (Muradian et al. 2019) and which ecological factors most likely regulate herring recruitment and natural mortality (Trochta and Branch 2021).
Various updates have been made to the equations of BASA since Muradian et al. (2017) to reflect changes in the understanding of herring biology and life history, as well as incorporate new data. These changes included, but were not limited to:
* Fitting hydroacoustic data to mature biomass instead of total age 3+ biomass (Trochta and Branch 2021)
* Estimating maturity (proportions mature at ages 3 and 4) over the entire modeling time period instead of two time periods split at 1997
* Incorporating aerial age-1 school count data to inform recruitment estimates before they enter the spawning population at age 3
* Integrating seroprevalence and disease data to inform natural mortality estimates (Trochta et al. _in review_)
Furthermore, a more efficient Markov Chain Monte Carlo (MCMC) algorithm called the No-U-turn Sampler (NUTS) has been introduced to rapidly sample posteriors of the parameters in BASA (Monnahan et al. 2019).
At present, BASA is primarily used to estimate spawning biomass and recruitment up to the most recent year with data. Persistent low levels of biomass and recruitment since the early 1990s continue to preclude consideration of reopening of fisheries under the current harvest strategy, and thus forecasts are not conducted. BASA has also been used as a research tool to investigate hypotheses and evaluate alternative models. In this report, we present the most recent fits and estimates from BASA for 2021 and summarize any modifications and alternative models explored.
# 2023 BASA Summary
To run the 2023 BASA model, the key software and versions used include:
* AD Model Builder v. 13.1
* R v. 4.3.2
* adnuts v. 1.1.2
The no-U-turn sampler (NUTS) was used within ADMB to sample the posterior distributions of BASA parameters and derived quantities. The 'adnuts' package and its dependencies were used to run NUTS and diagnostic checking from within R. Four NUTS chains were ran in total with the default arguments already supplied to 'sample_nuts()' (e.g. warmup=700, iter=2000), except for a higher target acceptance rate (adapt_delta=0.9) and using the inverse Hessian as the mass matrix (metric='mle'). Diagnostics supported convergence in all four chains (zero divergences and all R-hat convergence values < 1.05) and had sufficient sample size (estimated Bulk Effective Sample Size > 500 from merged chains). The total duration for running BASA was 1.8 minutes.
Results are shown from the BASA model fits to data up to and including 2023 (Figs. 1-3). The inner 95th percentiles of the posterior predictive distributions of the ongoing biomass survey data (Mile-days milt and PWSSC acoustic biomass) from BASA encompass all observations (Fig. 1). Fits of the discontinued data (egg deposition and ADF&G acoustic biomass) also fit well the historical time series.
```{r, out.width="85%", include=TRUE, fig.align="center", echo=FALSE}
knitr::include_graphics(here::here("figures/survey_fits.pdf"))
```
**Fig. 1. Estimated survey biomass from Bayesian age structured assessment (shading showing 50% and 95% posterior predictive intervals in dark and light gray, respectively) compared to indices of biomass in the population (points and lines showing observation CV).**
Posterior predictions of the juvenile aerial survey index (age-1 schools) bounded all observations, albeit with large uncertainty. BASA largely overestimated the 2017 index which was the largest in the available record, although the relative scale of this cohort (2016 age-0) agrees with the large proportions of age-3s, -4s, -5s, and -6s observed in 2019, 2020, 2021, and 2022 (Fig. 2). The consistent overestimation of schools since 2017 may be due to bias from a subjective standardization used to calculate this index; schools were numerated by four descriptive categories (small, medium, large, and extra large) and the largest three categories were converted to and summed as small school equivalents to calculate the index. Furthermore, the numbers of medium, large, and extra large schools in 2017 each represented the historical maxima in their respective categories, while the number of small schools was the third largest. Further investigation into the accuracy of this standardization is needed. The large aerial survey index in 2022 is also anomalous, in that two sequential years have never been observed to have near equally large numbers of age-1 schools. The exact reasons for this discrepency remain unknown at this time.
```{r, out.width="85%", include=TRUE, fig.align="center", echo=FALSE}
knitr::include_graphics(here::here("figures/age_compositions.pdf"))
```
**Fig. 2. Estimated age structure from the Bayesian age structured stock assessment (points = median, lines = 95% posterior predictive intervals) compared to the age composition data from catches and surveys (bars). Each color follows a single cohort as it ages through the fishery. Data are available only for ages-3 and above. The quantity 'J' in the upper-right corner of each panel is the Shannon-Weiner evenness index score of the age composition for each year (Shannon and Weaver, 1949).**
Posterior predictive intervals for the age composition data mostly show good fits, except for the age-3 classes in 1987 and 1998 (Fig. 2) as well as the age-4 and -7 class in 2023. The model misfit in 2023 could be due to a reduced sample size from the ASL survey being used. Estimates of the age-3 and age-4 year classes support strong 2019 and 2020 cohort. The median spawning biomass estimate in 2023 was approximately 28,700 metric tons which is above ADF&G's lower cut-off for fishing (Table 1). Additionally, uncertainty in this estimate indicates there was a 5% probability that 2023 spawning biomass was below this lower cutoff (Table 1).
The spawning biomass estimate for 2024 is approximately 27,050 metric tons (95% CI: 16,790 - 46,190), corresponding to a ~11% probability that spawning biomass is below the lower regulatory cutoff definde by ADF&G (Figure 3).
```{r, out.width="100%", include=TRUE, fig.align="center", echo=FALSE}
knitr::include_graphics(here::here("figures/management_outputs.pdf"))
```
**Fig. 3. Bayesian age structured assessment estimates of numbers of age-3 recruitment in millions, spawning biomass with 95% credibility intervals (light gray shading), total exploitation rate, and poterior probability density of pre-fishery biomass.**$\\$
\newpage
# 2023 BASA Diagnostics
A retrospective analysis was conducted to evaluate the effect that adding new data has on prior predictions of stock biomass. The Mohn's rho statistic, commonly used to assess the magnitude and direction of retrosepctive biad, was $\rho=-0.088$ in 2023, indicating a small negative bias in biomass estimates in previous years. This value of Mohn's rho is not considered indicative of substantial model misspecification, as it falls within the accepted interval of (-0.22, 0.30; Hurtado-Ferro et al. 2015).
```{r, out.width="100%", include=TRUE, fig.align="center", echo=FALSE}
knitr::include_graphics(here::here("figures/retrospective.pdf"))
```
**Fig. 4. 5-year retrospectives for the 2023 BASA model. Peels show no pervasive retrospective bias.**$\\$
\newpage
**Table 1.** Time series of posterior percentile (PCTL) estimates of the population from BASA. Pre-fishery spawning biomass (SB) is mature biomass of age-3+ fish. Catch (metric tons) includes all fisheries. Exploitation fraction is the catch divided by total age-3+ biomass. P(SB<20K) is the proportion of posterior SB samples that are less than 20,000 tons.
```{r, echo=FALSE,message=FALSE}
library(flextable,quietly = TRUE)
library(officer, quietly = TRUE)
#library(knitr,quietly = TRUE)
library(tidyverse,quietly = TRUE)
library(magrittr,quietly = TRUE)
set_flextable_defaults(fonts_ignore=TRUE)
ssb_est <- readr::read_csv(here::here("data_outputs/outputs-for-management.csv"))
names(ssb_est) <- paste0("col_",1:ncol(ssb_est))
tab_head <- tibble(col_keys=names(ssb_est),
line2=c("",
rep('Age-3 recruits\n(millions)',3),
rep("Spawning Biomass \n(1000 metric tons)",3),
'',
rep('Exploitation\nfraction',3),
''),
line3=c("Year",
'50th','2.5th','97.5th',
'50th','2.5th','97.5th',
'Catch',
'50th','2.5th','97.5th',
'P(SB<20K)'))
pgwid <- 8.5
ft_obj <- ssb_est %>%
regulartable() %>%
set_header_df(mapping = tab_head,key="col_keys") %>%
merge_h(part = "header", i=1) %>%
merge_v(part = "header", j=c(1,7,12)) %>%
colformat_double(j=c(1, 8), big.mark = "", digits = 0) %>%
colformat_double(j=2:7, big.mark = "", digits = 1) %>%
#theme_booktabs(bold_header = TRUE) %>%
bold(part="header") %>%
align(align="right", part = "all") %>%
valign(j=c(1,7,12),valign="top", part = "header") %>%
fontsize(size=9, part="all") %>%
hline_bottom(j=1:12,part = 'header', border = fp_border(width = 1.5)) %>%
hline_top(j=1:12,part = 'header', border = fp_border(width = 1.5)) #%>%
# set_formatter_type(fmt_double = paste0("%.0",decimal.places[j],"f")) %>%
# display(formatters = list(mpg ~ sprintf("%.01f", mpg) )) %>%
#autofit()
ft_obj %>%
width(width = c(0.3, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.8, 0.5, 0.5, 0.5, 1)) %>%
line_spacing(space=1, part="body") %>%
padding(padding=0,part="body")
# library(kableExtra,quietly = TRUE)
# ssb_est %>%
# kbl(caption = "Time series of posterior percentile (PCTL) estimates of the population from BASA. Pre-fishery spawning biomass (SB) is mature biomass of age-3+ fish. Catch includes all fisheries. Exploitation fraction is the catch divided by total age-3+ biomass. P(SB<20K) is the proportion of posterior SB samples that are less than 20000 tons.",
# col.names =c('Year',
# "2.5th PCTL",'Median','97.5th PCTL',
# '2.5th PCTL','Median','97.5th PCTL',
# 'Catch (thousand metric tons)',
# '2.5th PCTL','Median','97.5th PCTL',
# 'P(SB<20K)'),
# booktabs = T) %>%
# #kable_classic(full_width=T) %>%
# kable_styling(full_width=T) %>%
# add_header_above(c(" "=1,
# 'Age-3 recruits (millions)'=3,
# 'SB (thousand metric tons)'=3,
# ' '=1,
# 'Exploitation fraction'=3,
# ' '=1),
# bold =TRUE, line=TRUE) %>%
# row_spec(0,bold=TRUE,align='c') %>%
# #column_spec(c(1:7,9:11),width=paste0((pg_wid-2)/15,"in")) %>%
# column_spec(8,width=paste0((pg_wid-2)/15 * 3,"in")) #%>%
# #column_spec(12,width=paste0((pg_wid-2)/15 * 2,"in"))
```
\newpage
# References
Botz, J., G. Hollowell, J. Bell, R. Brenner, and S. D. Moffitt. 2010. 2009 Prince William Sound Area Finfish Management Report. Alaska Department of Fish and Game, Anchorage, Alaska.
Funk, F., and K. A. Rowell. 1995. Population model suggests new threshold for managing Alaska`s Togiak Fishery for Pacific herring in Bristol Bay. Alaska Fishery Research Bulletin 2:125-136.
Hurtado-Ferro, F., C.S. Szuwalski, J.L. Valero, S.C. Anderson, C.J. Cunningham, K.F. Johnson, R. Licandeo, C.R. McGilliard, C.C. Monnahan, M.L. Muradian, K. Ono, K.A. Vert-Pre, A.R. Whitten, and A.E. Punt. 2015. Looking in the rear-view mirror: bias and retrospective patterns in integrated, age-structured stock assessment models. ICES Jounal of Marine Science. 72(1): 99–110.
Monnahan, C. C., T. A. Branch, J. T. Thorson, I. J. Stewart, and C. S. Szuwalski. 2019. Overcoming long Bayesian run times in integrated fisheries stock assessments. ICES Journal of Marine Science 76:1477-1488.
Muradian, M. L., T. A. Branch, S. D. Moffitt, and P.-J. F. Hulson. 2017. Bayesian stock assessment of Pacific herring in Prince William Sound, Alaska. PLOS ONE 12:e0172153.
Muradian, M. L., T. A. Branch, and A. E. Punt. 2019. A framework for assessing which sampling programmes provide the best trade-off between accuracy and cost of data in stock assessments. ICES Journal of Marine Science 76:2102-2113.
Shannon, C. E., and Weaver, W. 1949. A Mathematical Theory of Communication, University of Illinois Press, Urbana.
Trochta, J. T., T. A. Branch, A. O. Shelton, and D. E. Hay. 2020. The highs and lows of herring: A meta-analysis of patterns and factors in herring collapse and recovery. Fish and Fisheries 21:639-662.
Trochta, J. T., and T. A. Branch. 2021. Applying Bayesian model selection to determine ecological covariates for recruitment and natural mortality in stock assessment. ICES Journal of Marine Science 78:2875-2894.
Trochta, J. T., M. Groner, P. Hershberger, and T. A. Branch. in review. A better way to account for disease in fisheries stock assessment: the powerful potential of seroprevalence data. ICES Journal of Marine Science.