Climat e-driv en c hang es in the timing of spa wning and the availability of walleye pollock ( Gadus chalcogrammus ) to assessment surv e ys in the Gulf of Alaska

Climate-driven changes in the timing of spawning or migration can affect the availability of fish to surveys designed to monitor their abundance, complicating ef for ts to assess stock status and sustainably manage fisheries. From 20 17 to 20 19, trends in biomass estimates from four surveys used to monitor Gulf of Alaska pollock diverged. These conflicting trends increased uncertainty in the stock assessment and occurred during a time of rapid environmental change. We h ypothesiz ed that changes in spawn timing affected availability of pollock to a winter survey that targets pre-spawning aggregations. To test this, we reconstructed relative spawn timing using two independent data sources: spring larval surveys and observations of spawning state in mature female pollock. We found that changes in spawn timing relative to survey timing explained a significant portion of recent and historical discrepancies between survey and model estimates of biomass. We then incorporated measures of spawn timing/survey timing mismatch as catchability covariates in an enhanced state-space stock assessment model. Including spawn timing-based catchability covariates significantly improved the model fit to survey data and provided a mechanistic explanation for recent survey discrepancies


Introduction
Changes in climate and ocean conditions are increasing the need for climate-informed approaches to fisheries management (Holsman et al. 2019, Karp et al. 2019 ).Shifts in species geographic distributions and productivity in response to climate have been well documented (Pinsky et al. 2013, Free et al. 2019 ), increasing risks to fishing communities dependent on fisheries resources (Rogers et al. 2019, Payne et al. 2021 ).Detecting and understanding the effects of climate on fish stocks is a key element of climate-informed fisheries management.Monitoring surveys are instrumental for this effort and can provide early warnings of shifts in fish stock abundance or distribution, as well as the ecosystem, food web, or life history information necessary to identify the mechanisms behind such shifts.For instance, monitoring surveys of adults and early life stages of Pacific cod ( Gadus macrocephalus ) in the Gulf of Alaska (GOA) detected a rapid decline in stock biomass and concurrent failed recruitment associated with a 2014-2016 marine heatwave, necessitating the closure of the federal directed fishery (Barbeaux et al. 2020 , Laurel andRogers 2020 ).
While surveys are important for assessing, detecting, and understanding the effects of climate on fish stocks, surveyderived quantities, such as estimates of biomass or abundance, may in turn be affected by climate-driven changes in fish stock distribution and migratory behavior (Karp et al. 2019, Bell et al. 2022 ).For instance, reduced sea ice in the Bering Sea resulted in a northward shift in the distribution of commercially important species, including walleye pollock ( G .chalcogram-mus ) and Pacific cod, reducing their availability to the longrunning bottom trawl survey in the region (Stevenson andLauth 2019 , O'Leary et al. 2022 ).Understanding and being able to account for climate-driven changes in availability of fish stocks to surveys is an important component of climateready fisheries management.
Shifts in phenology, or the seasonal timing of events, are a well-documented response to changes in climate (Parmesan and Yohe 2003 ).The timing of spawning of many fishes, in particular, has been shown to be sensitive to thermal conditions, often reflecting physiological or behavioral responses to temperature cues (e.g.North Sea mackerel, Jansen and Gislason 2011 ;Atlantic cod, McQueen and Marshall 2017 ;capelin, Carscadden et al. 1997 ).Beyond implications for recruitment and survival of offspring, changes in spawn timing can lead to practical challenges for estimating the abundance of fish stocks.In temperate to arctic systems, many fish species seasonally aggregate to spawn in contracted areas (Harden Jones 1981 ), which is convenient for both fisheries harvest (e.g. the Lofoten cod fishery of Norway) and biomass surveys, often conducted using acoustics (e.g.Kloser et al. 1996 ).However, annual changes in spawn timing can affect availability to surveys that are fixed in time and cover a limited geographic area.For instance, temperature-driven changes in the spawning migration timing of yellowfin sole ( Limanda aspera ) into nearshore areas of the Bering Sea affect their availability to bottom-trawl surveys monitoring their abundance (Nichol et al. 2019, Olmos et al. 2023 ).With ongoing climate change, changes in phenology may complicate efforts to assess stock status and sustainably manage fisheries.
Walleye pollock (hereafter pollock) are the target of an economically important fishery in the GOA, with an exvessel value of $34.7 million (2016-2020 average;Fissell et al. 2022 ).Pollock in the GOA are managed as a single stock, distinct from pollock stocks in the Bering Sea and Aleutian Islands.Catch limits for GOA pollock are determined by the North Pacific Fishery Management Council based on scientific advice derived from a peer-reviewed age-structured stock assessment model (Monnahan et al. 2021 ).The assessment model relies on fishery-independent estimates of biomass and age/length composition from four surveys, including a winter acoustic-trawl (AT) survey that targets pre-spawning pollock aggregated in Shelikof Strait, and three summer surveys over the broader GOA shelf.In recent years, biomass estimates from the four surveys have diverged, giving conflicting estimates of survey biomass and temporal trends.In particular, from 2017 to 2019, estimates of pollock biomass were record high in the winter Shelikof Strait AT survey, while the other surveys estimated near-record low biomass (Dorn et al. 2020 ).These conflicting trends increased uncertainty in the stock assessment, leading in part to the use of precautionary buffers when setting catch limits (Dorn and Zador 2020 ).
These divergent survey trends occurred during a time of rapid environmental change in the GOA, as one of the largest marine heatwaves on record altered thermal environments and ecosystem processes from 2014 to 2016 and again in 2019 (Rogers et al. 2021, Suryan et al. 2021 ).Simultaneously, the spawning stock was dominated by a single large year class spawned in 2012 (Dorn et al. 2020 ).During the 2018 winter Shelikof Strait AT survey, an unusually high proportion of females were in spawning or spent stages (Stienessen et al. 2019 ), suggesting that spawning was earlier than usual relative to the timing of the survey.These observations (divergent population biomass estimates and unusual spawn timing) led to an examination of the hypothesis that the availability of pollock to the winter Shelikof Strait AT survey may vary depending on the timing of spawning.Previous work found that GOA pollock spawning occurs earlier when temperatures are warmer and when spawners are older (Rogers and Dougherty 2019 ).
In this study, we investigate the hypothesis that changes in spawn timing relative to survey timing affect availability of pollock to the winter Shelikof Strait A T survey.W e test this hypothesis by reconstructing relative spawn timing from larval surveys and from observations of spawning state from the A T survey.W e then compare estimates of spawn timing/survey timing mismatch with residual errors from the age-structured stock assessment model to evaluate whether model lack-offit to survey biomass estimates is related to the timing of spawning relative to the timing of the survey .Finally , we test whether catchability covariates can improve the stock assessment model fit to the winter Shelikof AT indices by modeling survey catchability as a function of the timing covariates within the assessment model directly.

Walleye pollock in the GOA
Walleye pollock is a semi-pelagic gadid distributed throughout the continental shelves of the North Pacific Ocean.In the GOA, summer distributions of pollock are often widespread on feeding grounds across the shelf, whereas in winter, pollock aggregate in discrete locations and times to spawn.Shelikof Strait, near Kodiak Island, is a major spawning ground for GOA pollock ( Fig. 1 ); however, smaller spawning aggregations have been detected in other areas, such as the Shumagin Islands and near Prince William Sound (Wilson 1994, McCarthy et al. 2018 ).In the Shelikof Strait, pollock begin to aggregate in January or February prior to spawning (De Robertis et al. 2018 ).Spawning typically begins in March and can last into May, with the majority of eggs spawned in early to mid-April, although the timing and duration of spawning are variable from year to year (Rogers and Dougherty 2019 ).As multiple batch spawners, female pollock release 10-20 batches of eggs over a period of 3 weeks or more (Hinckley 1990 ).After spawning, pollock rapidly disperse from spawning grounds (Wilson 1994, De Robertis et al. 2018 ).Eggs hatch after ∼2 weeks and larvae rise to the upper 50 m of the water column, where they become entrained in coastal currents and are advected to the southwest.Pollock typically mature at ages 3 to 4, although 1-and 2-year-old fish are also found on the spawning grounds during the spawning season.

Winter acoustic-trawl survey and maturity assessment
The winter Shelikof Strait AT survey has been a key source of biomass estimates and biological information (age, length, and maturation) for the assessment of GOA pollock.The survey has been conducted annually since 1980; however, no survey was conducted in 1982, 1999, or 2011.The survey was timed to sample the aggregation of pre-spawning pollock in Shelikof Strait prior to the peak of spawning, with the hypothesis that pollock would disperse shortly after spawning and no longer be available to the survey (Nelson and Nunnallee 1985 ). Historically, the survey has typically occurred in the second half of March, although ship availability, weather delays, and proactive planning for potential earlier spawning have led to variations in survey timing, and recent surveys (2019)(2020)(2021) were conducted in the first half of March.Agestructured biomass and abundance estimates were derived from acoustic backscatter collected at a frequency of 38 kHz that is attributed to pollock, which is expanded between parallel transects and proportionally scaled to the acoustic target strength of the size of pollock caught in midwater trawls in the surveyed area (Stienessen et al. 2019 ).
On-board assessments of maturity have been used to estimate maturity-at-size.Male and female pollock were randomly sampled from hauls, measured for fork length (FL), and maturity examined macroscopically using a 5-or 8-stage key (Williams 2007 ).The proportion of females > 40 cm FL in a spawning or spent stage (as opposed to immature, developing, or pre-spawning stages) has been used as an indicator of survey timing relative to spawn timing (Stienessen et al. 2019 ), with an approximate target of < 10% in the spawning and spent stages combined.

Larval survey and spawn timing estimation
Early life stages of walleye pollock have been sampled since 1979, with spring ichthyoplankton surveys targeting larvae as they drift from presumed spawning grounds in Shelikof Strait.Rogers and Dougherty (2019) demonstrated that data on larval length, age, and catch per unit area from these surveys can be used to estimate pollock hatch date distributions.Timing of spawning can then be back-calculated using temperaturedependent egg development rates, which predict an incubation period of 12-16 days depending on thermal conditions.For this study, analyses in Rogers and Dougherty (2019) were extended to include two additional years of larval data (2017 and 2019), resulting in 34 years of estimated spawn timing distributions.These annual distributions were then characterized by the mean, median, and 20th percentile date of spawning.The 2019 ichthyoplankton survey was ∼10 days earlier than most historical surveys, which necessitated adjusting the fixed survey date window for samples included in the spawn timing analysis, from 17 May to 7 June in Rogers and Dougherty (2019) to 10 May to 7 June in this study.This had a negligible effect on estimates of mean and median spawn dates (estimated change of < 1 day for all years except 1989, when mean and median spawn timing were estimated to be 3 days earlier using the new temporal window).Since 2011, larval surveys have been carried out biennially (odd-numbered years only); thus, spawn timing could not be estimated in this way for recent even-numbered years.

Relative timing metrics
To test the hypothesis that mismatches in the timing of the winter Shelikof Strait AT survey relative to the timing of spawning affected the availability of pollock to the AT survey, we developed two approaches to estimate relative timing.
The first is based on direct observations of mature females in the winter Shelikof Strait AT survey to assess the proportion of fish that were in a spawning or spent stage.A higher proportion of fish in spawning or spent stages could indicate that the survey was more closely timed with peak spawning or possibly after peak spawning, while a lower proportion in spawning or spent stages could indicate the survey was early relative to the timing of spawning.The unweighted proportion of females > 40 cm FL in spawning and spent stages combined (as a proportion of all females > 40 cm) has been used as a metric of relative survey timing in the past (Stienessen et al. 2019 ).However, pollock often mature smaller than 40 cm and those individuals may also provide information on relative survey timing.We thus calculated spawning stage metrics using minimum size cutoffs of both 40 cm FL and 30 cm FL for comparison.
Because only a subset of females were sampled in each haul and because the total biomass associated with each haul varied, we calculated maturity proportions in three ways: (1) unweighted (all sampled females weighted equally; SP 30_fish, SP 40_fish ); (2) equal weighting by haul ( SP 30_haul, SP 40_haul ); and (3) abundance weighting (i.e.haul-level average maturity weighted by acoustically derived abundance estimates associated with that haul; SP 30_wt, SP 40_wt ).Proportions were logit transformed following Rogers and Dougherty (2019) to be on a continuous, unbounded scale prior to inclusion in models.In years prior to 2002, multiple survey passes of the study area were often conducted, in which case only maturation data from those trawls used to derive biomass and abundance estimates for the assessment model were included in this study.
The second approach to estimating survey timing relative to spawn timing uses the estimated spawning dates from the larval survey to calculate the relative temporal lag between the AT survey and the inferred timing of peak spawning.Survey reports and databases were used to identify the start, end, and mean dates of the AT survey, specifically identifying the dates for individual transect lines and trawls in Shelikof Strait and Sea Valley ( Fig. 1 ).Again, when multiple survey passes of the study area were conducted, only survey dates for those Downloaded from https://academic.oup.com/icesjms/advance-article/doi/10.1093/icesjms/fsae005/7591717 by NOAA Central Library user on 15 February 2024 passes used in the assessment model were included.The number of days between the mean date of the survey and the mean ( mismatc h mean ), median ( mismatc h med ), and 20th percentile ( mismatch 20 ) of spawning dates were calculated as metrics of relative survey timing.

Stock assessment model
The operational stock assessment model for GOA pollock is peer reviewed and fully described in Monnahan et al. (2021) , with relevant characteristics briefly summarized here.The custom age-structured ADMB model (Fournier andArchibald 1982 , Fournier et al. 2012 ) simultaneously fits to age compositions and indices of abundance from the winter AT survey in addition to data from the fishery and three other surveys mentioned previously in what is known as an integrated assessment model (Maunder and Punt 2013 ).While AT surveys started in 1980, only surveys since 1992 have been included in the stock assessment.The catchability for the winter AT survey is modeled as a time-varying process and follows a constrained random walk on the log scale, with annual deviations assumed to follow a N(0,σ q = 0.038) distribution.This formulation allows catchability to change smoothly over time in response to annual variation in the proportion of the stock spawning in Shelikof Strait, and was motivated by observations of smaller spawning aggregations in surrounding sites fluctuating relatively smoothly over time.Because there are many potential contributors to variation in catchability, this component is constrained by fixing the process error at a small value specifically to try to only account for the spatial availability component of variation in catchability; however, some portion of a change in availability due to spawn timing may also be absorbed by this random walk term.To more explicitly test the potential importance of spawn timing to variation in catchability, we therefore explored alternative formulations for catchability.To do so, it required porting the assessment model to a different platform for computational efficiency.
For this study, the stock assessment model was converted to the Woods Hole Assessment Model (WHAM) platform, which has the capability to estimate process errors by integrating out random effects via its backend Template Model Builder (Kristensen et al. 2016 , Stock andMiller 2021 ).WHAM has specific features for modeling covariate links to catchability (and other model quantities), and missing covariate values can be interpolated via a first-order autoregressive [AR(1)] time series smoother as part of the fitting process (see section 2.1.1.4 of Stock and Miller 2021 ).Our WHAM model is not identical to that used for management.It does not include an aging error matrix, has a different configuration for initial numbers at age, and has a different selectivity (age-based instead of descending logistic) for the Shelikof Strait AT survey.These differences are relatively minor and the models are similar enough for the purposes of this study.

Covariate assessment
To explore the potential effect of spawn timing/survey timing mismatch on survey catchability, we compared the timing metrics described above (see Relative timing metrics ) to residuals from the stock assessment model that measure the lack of fit of the model to survey-derived biomass estimates (ages 3 + ).These residuals are treated as a proxy for the relative availability of pollock to the winter Shelikof Strait AT survey.Specifically, we fit a model with constant catchability ( q ) to eliminate the possibility that the time-varying q was absorbing some spawn-timing-related variation in availability.Residuals were calculated as R y = log(biomass_survey_age3p y )log(biomass_predicted_age3p y ), with log-transformation used to reflect the log-normal assumptions in the stock assessment model.
We tested the ability of the timing covariates to explain the variance in the survey residuals using linear models with single covariates.Models were fit to a dataset that included only years with complete observations of all covariates (1992-1998, 2000-2010, 2013, 2015, 2017, and 2019) to facilitate model comparisons.A second set of models was fitted to the data preceding 2015 to test whether covariate relationships were present prior to the recent marine heatwave and dominance of the 2012 year class.Finally, models were refitted using all available data for each covariate, which in many cases was a longer dataset than available for the full model comparison.
In addition to the relative timing covariates described above, we tested temperature covariates alone to reflect a climate-informed but mechanistically incomplete model for changes in availability.Mean monthly sea-surface temperatures for January , February , and March ( SST Jan , SST Feb , and SST Mar ) from the NCEP Reanalysis (Kalnay et al. 1996 ) were extracted for the latitude/longitude box defined by 56.2 to 54.3 • N and 151.9 to 157.5 • W. Previously, March sea surface temperatures were found to be a significant predictor of spawn timing in pollock (Rogers and Dougherty 2019 ).Analyses were conducted using R 4.3.2(R Core Team 2023 ).

Including catchability covariates in the assessment model
We selected a small set of covariates to test in the WHAM stock assessment model following the analysis of survey residuals.We tested a mismatch covariate, a spawning state ( SP ) covariate, and a temperature ( SST ) covariate.We also simulated independent and identically distributed (iid) standard normal "noise" covariates, whose inclusion was used to test assessment behavior in the presence of no signal and compare it against the behavior of our real covariates, as has been done in previous studies (e.g.Deriso et al. 2008 ).
We evaluated four parameterizations of catchability within the assessment based on the full model, which model the log of catchability ( q ) as where μ is a global intercept, τ the annual random walk process with estimated process error σ q , and β the estimated fixed effect coefficient for the AR(1) smoothed environmental covariate E y .The four parameterizations were (1) constant catchability log( q y ) = μ, (2) a random walk log( q y ) = μ + τ y , (3) a covariate effect log( q y ) = μ + βE y , and (4) a random walk plus a covariate effect (Equation 1), where again, the random walk process error σ q is estimated.Formal inclusion and testing of environmental effects in stock assessments is a developing field, with many questions about best practices, including how to appropriately conduct model selection (Wilberg et al. 2010 , Stock andMiller 2021 ).Here, we hypothesized that the inclusion of a linear effect of timing covariates would have three important practical effects.First, the model fit to the abundance index would improve because the covariate explains the lack of fit caused by survey timing.Second, the total amount of variation in catchability explained would increase.Specifically, if the process error of the random-walk component ( σ q ) were estimated (as opposed to fixed), it could increase to account for multiple sources of variation, including survey timing.However, with an added covariate effect, the estimate of σ q would decrease because the timing portion of that variability would be explained by the covariate.Finally, metrics of model goodness of fit would improve, and metrics of predictive skill like AIC (Burnham and Anderson 2002 ) would be reduced by linking the covariates to catchability.We compared model fits, estimates of σ q (when included), and AIC values between these model versions to assess the ability of catchability covariates to improve the assessment model fit.Because AIC requires the data to be the same to compare models, we also included the covariate data in parameterizations ( 1) and ( 2) and fixed β = 0, so the covariate was smoothed but had no effect on catchability.

Acoustic-trawl survey timing
Since 1992, the AT survey of the Shelikof Strait has typically started around 16 March and ended around 24 March.However, the AT survey has started as early as 3 March and as late as 25 March ( Fig. 2 ).The mid-survey date was 20 March on average.Surveys in 2019, 2020, and 2021 were 1-2 weeks earlier than the long-term average, with mid-survey dates of 12, 11, and 7 March, respectively, which was in part an adaptive shift in response to perceptions of earlier spawn timing due to higher percentages of spawning/spent females, and in part a result of ship-time constraints.

Reconstructed spawn timing
As estimated from larval survey data, the spawn timing of pollock was found to vary from year to year, with an earliest mean (median) date of spawning on 2 April (28 March) and a latest mean date of spawning on 4 May (4 May) ( Fig. 2 ).Spawning during 2017 and 2019, which were not included in Rogers and Dougherty (2019) , was the earliest on record, peaking approximately a week earlier than in the next earliest year and 14-16 days earlier than the long-term average.Since 1992, estimated spawn timing has been more variable than AT survey timing ( Fig. 2 ; SD spawn_timing = 8.3 days, SD survey = 4.7 days).Data were not available for recent even-numbered years.

Spawn timing/survey timing mismatch
The AT survey has preceded the estimated timing of peak (mean or median) spawning in all years, consistent with its design as a pre-spawning survey ( Fig. 2 ).However, the number of days between the mid-survey date and the estimated median date of spawning varied from 14 to 44, with an average of 28 (SD: 9.1 days).Greater interannual variation was found in the number of days between the mid-survey date and the date when 20% of spawning was estimated to have occurred (range: 3 to 40, mean: 21, SD: 11.0).one another (Pearson's correlation coefficients ranged from r = 0.97 to r = 0.99).Results are reported here for SP 30_wt, but mirror those of the estimates using only females > 40 cm, as well as the other approaches to weighting fish-level observations.The proportion of females in spawning or spent stages ranged from < 0.001 in 2008 to 0.81 in 2018.In most years since 1992 (22 of 28), the proportion was under 0.1 ( Fig. 3 ).The surveys in 2017, 2018, and 2019 encountered a higher proportion of fish in spawning or spent stages than most other survey years, with proportions of 0.37, 0.81, and 0.21, respectively.The proportion of pollock in spawning or spent stages was highly correlated with estimates of spawn timing/survey timing mismatch (r = −0.82,P < 0.001 for SP 30_wt and mismatch med ) with more pollock in spawning or spent stages in years when the survey was estimated to be closer in timing to the date of peak spawning, suggesting both metrics captured changes in the timing of spawning relative to the survey ( Fig. 4 ).Note that this comparison includes all years with concurrent larval and winter AT surveys, including those that precede 1992.

Stock assessment residuals and covariate comparison
Survey residuals (R y ) from the assessment model fit to winter Shelikof Strait AT survey biomass estimates were large in magnitude and highly positive for 2017-2019, indicating an inability for the assessment model to fit to survey data in those years and a likely conflict with other data sources ( Fig. 5 a).Other years with relatively large residuals were 2008, 2009, and 2012 (negative residuals), while in earlier years residuals were smaller on average.All measures of relative timing, based on either reconstructed spawn dates from larval otoliths or the proportion of females in spawning and spent stages, explained significant variation in R y for the winter Shelikof Strait AT survey ( Table 1 , all P < 0.01).Comparing among common years, mismatch med had the strongest relationship, explaining 58% of the variation in survey residuals ( Fig. 5 b).The best maturity-based covariate, SP 30_wt , explained 46% of the variation ( Fig. 5 c).The estimates of relative timing using larval data ("mismatch" covariates) tended to explain more of the variation than the maturity/spawning state estimates ("SP" covariates).When examining only years prior to 2015, the rank order of covariates changed slightly, but all remained highly significant ( Table 1 ).Models with SST alone were also significant, particularly SST Mar ( Fig. 5 d), but explained less variation than the more explicit timing-related covariates.

Catchability covariates in assessment model
We fit stock assessment models that included either SP 30_wt , mismatch med, or SST Mar as a covariate on catchability for the winter AT survey.Due to missing observations, six data points were interpolated for mismatch med using an AR(1) smoother within the model.In models with a covariate and a random walk term, the two timing covariates had significant slopes; that is, 95% C.I. did not overlap 0: for SP 30_wt β 1 = 0.26 (0.12-0.40), and for mismatch med β 1 = -0.36(-0.49-0.23).The slope for SST Mar was positive but did overlap zero: SST Mar β 1 = 0.07 (-0.07-0.22).Models with covariates but no random walk had steeper slopes and all were significant ( Supplementary Fig. S1 ).The AIC values decreased by 8.3 and 12.3, respectively, with the inclusion of the timing covariates relative to models with only a random walk term, whereas the AIC increased by 1.0 for the assessment model with SST Mar ( Table 2 ).The estimated indices (and thus residuals) were similar for the versions of the WHAM model with a random walk term with and without each covariate, meaning that both a random walk and a random walk with a covariate were able to resolve the apparent conflict among the different surveys ( Fig. 6 ).However, the process variance for the random walk had to increase to fit the survey data in the absence of an informative covariate.There was a reduction in the process variance for the random walk of 54% and 75% for the models that included SP 30_wt and mismatch med covariates, and a reduction of 6% for the model with the SST Mar covariate ( Table 2 ).Models that only had catchability covariates and no random walk terms outperformed models with constant catchability in all cases ( Table 2 ) and provided an improved fit to the data even without the flexibility of a random walk ( Fig. 6 ).In comparison, models fit with either of the "noise" covariates had estimated slopes ( β 1 ) with CIs overlapping zero and higher AIC than models fit with a random walk term alone.The random walk sigma decreased by 0% and 6% with the inclusion of each noise covariate, respectively ( Supplementary Table S1 ).

Discussion
We found strong evidence that changes in the spawn timing of pollock relative to the timing of the winter AT Shelikof Strait survey are associated with the availability of adult (age 3 + ) pollock to this survey, and thus the proportion of the stock this survey encounters in a given year.While the survey has preceded larval-estimated peak spawning in all years studied, as it was designed to do, in years when the survey was closer in time to peak spawning (i.e. a later survey and/or earlier spawning), a higher proportion of the stock was estimated to have been surveyed.A change in spawn timing relative to sur-vey timing may at least partly account for the unusually high biomass estimates in 2017-2019, which were in contrast to the biomass estimates from other surveys.However, analysis of observations from years prior to 2015 showed these relationships as well, suggesting that it was not only the somewhat extreme 2017-2019 years that drove the estimated relationship between spawn/survey timing and availability.Accounting for spawn timing helps resolve an outstanding issue with divergent trends between different surveys of this stock, and has the potential to improve management advice if incorporated into the operational assessment.Historically, the conventional wisdom has been that the winter AT survey should occur in the second half of March to correspond with the highest densities of pre-spawning aggregations of pollock (Nelson andNunnallee 1985 , Wilson 1994 ).When multiple survey passes of the spawning grounds were conducted in 1981, 1983, 1984, and 1994, re-sults showed that biomass estimates decreased after ca.April 1, suggesting that post-spawning pollock were beginning to leave the spawning grounds by that time (Nelson and Nunnallee 1985, Nunnallee and Williamson 1989, Wilson 1994 ).More recently, De Robertis et al. ( 2018) reported concurring observations from echosounder moorings and repeated AT survey passes; taken together, these showed a rapid decrease in biomass in Shelikof Strait in April-May following high biomass from late February through the end of March.Based on these findings, it has been assumed that a late survey would underestimate the biomass.We found that availability to the survey apparently increased when the survey was closer in time to peak spawning (i.e. a later survey and/or earlier spawning, with a higher proportion of females spawning or spent), which would seem in conflict with these previous findings.One possible explanation is that the standard survey has historically never been "too late," and the hypothesized reduction in availability that would occur if the survey occurred post-spawning has not been observed during the time-frame of this study (1992 onwards).If it had, we would expect to detect a dome-shaped relationship between availability (survey residuals) and our relative timing metrics.The historical repeated survey passes do suggest an asymmetric risk, where the risk of being too late outweighs the risk of being too early, supporting the current survey timing despite evidence that an early survey also underestimates biomass (though possibly to a lesser degree than a late survey).
While evidence for changes in spawn timing seems robust, as it is based on two sources of data and reflects patterns observed in other related species (e.g.Atlantic cod), the exact mechanism linking spawn timing with survey availability may be more complicated than we have hypothesized.Changes in availability may reflect both seasonal migration timing, as well as changes in the overall proportion of the GOA stock that aggregates (and possibly spawns) in Shelikof Strait in a given year.This proportion has been suggested to vary through time (Ciannelli et al. 2007 ), although alternative spawning locations are not routinely sampled, so this has not been quantified.We currently lack an understanding of why interannual variation in the proportion of the stock aggregating in Shelikof would correlate with spawn timing (as estimated by the larval survey or the Shelikof maturation data), but such a result could be explained by shared environmental or demographic drivers of spawn timing and location choice.Additional field studies, such as extended spatial and temporal sampling of spawning aggregations or eggs, could help to resolve this uncertainty.
Climate-ready fisheries management requires reassessing our surveys and assessment methods to ensure they capture ongoing changes in species distributions, phenology, life history traits, and other characteristics necessary for monitoring stock dynamics.One possibility for adaptation is to adjust sur-  Results are reported for models using only the subset of years in common for all covariates ("All Common Years"), all common years prior to 2015 ("Prior to 2015"), and for all years that covariate values are available ("All Years").
vey timing or spatial extent to track species phenology or distribution.For instance, in the Bering Sea, a northern survey region is now routinely sampled following a northward distribution shift of many species during a recent warm period (O'Leary et al. 2020(O'Leary et al. , 2022 ) ).In systems with high interannual variability, pre-season forecasts of distribution or phenology could aid in adaptively planning surveys to track species (Payne et al. 2017 ).On the US West Coast, environmentally driven forecasts of the distribution of sardine (Zwolinski et al. 2011 ) and hake (Malick et al. 2020 ) could enable survey footprints to be dynamically adjusted within a year.Phenology metrics may also be forecasted, as in the case of Maine lobster (Mills et al. 2017 ), where the timing of the start of the high-volume fishing period is forecasted up to four months in advance based on temperature.Skillful forecasts alone are not sufficient, however, and need to be combined with operational flexibility in order to be able to adjust survey timing or locations from year to year.In the case of GOA pollock, spawn timing depends on winter temperatures and spawning stock age structure (Rogers and Dougherty 2019 ), both of which could likely be forecasted with some skill up to a few months in advance.However, optimizing survey timing to match the forecasted spawn timing is impractical given the constraints of planning and executing surveys in late winter in Alaska.
In the absence of dynamically adjusting survey timing or footprints, a second option for adaptation is to account for changes in availability by adjusting survey indices using analytical approaches (e.g.O' Leary et al. 2020, Webster et al. 2020 ) or by incorporating changing catchability directly within the stock assessment (Thorarinsson andJohannesson 1997 , Nichol et al. 2019 ).The latter is the approach that we demonstrate here by including catchability covariates directly in the assessment model.We evaluated multiple configurations for modeling time-varying catchability, and models with spawn timing-related covariates consistently outperformed those without covariates.The WHAM framework used here provides a flexible and powerful tool to incorporate such covariate effects in complex hierarchical models (Stock and Miller 2021 ) and thus has potential for future assessment needs.However, additional research on best practices for testing and incorporating environmental covariates is needed as climate-linked stock assessments are increasingly evaluated as tactical tools for providing management advice.
While metrics based on spawn dates estimated from larval surveys explained the most variance in stock assessment residuals and resulted in the greatest reduction in AIC when included in the assessment model directly, they are less practical in an applied framework.Larval surveys have not been conducted on an annual basis since before 2011, and funding constraints continue to result in reduced survey frequency.While WHAM can incorporate covariates with missing values, the lack of spawn timing estimates in recent even-numbered years limits the utility of this time series as a catchability covariate in the pollock stock assessment.The high correlation between estimated survey/spawn timing mismatch and maturationbased metrics, as well as the strong performance of regression models with maturation-based metrics as predictors, provide support for using on-board maturity assessments to develop catchability covariates.
Temperature is possibly the most common environmental covariate included in stock assessment models (Marshall et al. 2019 ), particularly as a catchability covariate, reflecting its fundamental role in biological processes and association with    with fixed σ q = 0.038, labeled "fixed RW." Estimates are shown for models with constant catchability ("constant q"), RW with estimated process error ("RW"), and combinations of single covariates with and without a RW term.
changes in species depth, distribution, and phenology, all of which can alter species' availability to surveys.For GOA pollock, warmer temperatures in March were previously found to be associated with earlier spawning, which we hypothesize should increase availability to the AT survey on average.Indeed, temperature did explain significant variation in survey residuals.However, when SST was included as a catchability covariate directly within the assessment model, the model was not a significant improvement over a model with a random walk term alone, and its inclusion only marginally reduced the estimated process error variance.Temperature alone provided an improved fit to survey biomass over a model with constant catchability, particularly prior to 2017, but could not explain the strong discrepancies during 2017-2019.This model failed to capture the changes in spawn timing due to population age structure, as well as the interannual variation in timing of the survey itself.In other words, this model was climate-informed but mechanistically naive and was outperformed by models that more directly incorporated the hypothesized mechanism controlling the availability of pollock to the survey (although note that other mechanisms may also contribute to the observed patterns, as discussed above).This highlights the importance of research focused on mechanisms so that we can move beyond a purely empirical approach (see also Mills et al. 2017 ), and thus a continued investment in process studies in order to adapt our monitoring and stock assessments to climate change.
Spawn timing estimates from this study and Rogers and Dougherty (2019) are not always consistent with previous reports for pollock in the Shelikof Strait.As mentioned above, repeat AT survey passes in 1981AT survey passes in , 1983AT survey passes in , 1984AT survey passes in , and 1994, and echosounder , and echosounder moorings in 2015 found that spawning activity peaked in late March to early April (Nelson and Nunnallee 1985, Wilson 1994, De Robertis et al. 2018 ), which is a few weeks prior to the larval-derived estimates of peak spawning (mid-April in most years).One possible explanation for this discrepancy is that the larvae sampled from the late May to June period are not fully representative of the Shelikof spawning output.Some limited samples from larval gear comparisons (Bailey et al. 1996 ) suggest that larvae > 15-18 mm Downloaded from https://academic.oup.com/icesjms/advance-article/doi/10.1093/icesjms/fsae005/7591717 by NOAA Central Library user on 15 February 2024 ( ∼50 days old) become less available to bongo nets used to sample larvae for this study.This could lead to an undersampling of the earliest-hatched larvae, which would bias the estimated distribution of spawning dates to be relatively late.A second possibility is that sampled larvae could also be derived from spawning activity outside of Shelikof Strait, which could occur after the timing of pre-spawning aggregations targeted by the AT survey.Studies of pollock egg distributions, which can be used as a proxy for spawning activity, have found that egg concentrations are highest in Shelikof Strait in late March-mid-April, but remain elevated across a broader region throughout April and into May (Dunn et al. 1984, Ciannelli et al. 2007 ).The time of maximum pollock biomass in Shelikof Strait, and thus availability to the AT survey, may therefore be earlier than the estimated dates of peak spawning derived from larval data.Despite these discrepancies, Rogers and Dougherty (2019) found that interannual variation in estimated mean spawn dates corresponded closely to relative timing estimated from maturation data, indicating that phenological shifts estimated from larvae captured in late May and June are reflective of interannual variation in spawn timing for the Shelikof spawning stock as a whole.
In general, many questions remain about the spawning behavior, migration timing, and population structure of GOA pollock.Larval-derived spawn dates are used here as a proxy for when pollock are in Shelikof Strait and thus available to the AT survey, but it is not known for how long individual spawners aggregate prior to spawning, how quickly they disperse afterwards (although see De Robertis et al. 2018 ), nor how this may change with temperature, spawning stock size, age composition, or other factors.Rogers and Dougherty (2019) found that larval-derived estimates of spawn timing were not only earlier but also more protracted in warmer years, which could further affect availability to the AT survey.Age-, cohort-, or sex-related differences in spawning behavior and timing could also play a role.Finally, Shelikof Strait is considered the most important, but not the only spawning location for GOA pollock, as noted previously, and changes in the relative use of different spawning locations through time could also affect survey availability.Further studies of pollock spawning behavior could improve our ability to interpret survey biomass estimates and monitor assess this stock, particularly as novel environmental conditions emerge with climate change.

A c kno wledgments
We gratefully acknowledge the many personnel who have led and contributed to the NOAA-AFSC winter acoustic-trawl and spring larval surveys over the last four decades.Annette Dougherty aged the larval pollock used in this study.Patrick Ressler and Ben Williams provided detailed feedback on an earlier version of this manuscript, and two anonymous reviewers contributed to the clarity of this manuscript.This paper is contribution EcoFOCI-1048 to NOAA's Ecosystems and Fisheries-Oceanography Coordinated Investigations Program.

Figure 1
Figure 1 Shelikof Strait region of the GOA indicating the typical winter Shelikof Strait acoustic-trawl (AT) survey track lines and larval survey sampling domain, as well as the locations of maturity samples taken during the AT survey.Dominant currents flow to the southwest along the shelf, resulting in the pollock larval distribution in May/June being displaced to the southwest from the core spawning grounds in Shelikof Strait.

Figure 2
Figure 2 Historical timing of the winter Shelikof Strait AT survey (orange) and estimated timing of spawning based on larval hatch dates (dark red with circles), showing the median and 20th percentile dates of spawning.The pollock stock assessment uses estimates of surv e y biomass from the AT surv e y beginning in 1992, indicated by the darker shading.

Figure 3
Figure 3The proportion of female pollock ( > 30 cm FL) by maturity stage, captured during the winter Shelikof Strait AT survey.Proportions were w eighted b y estimated biomass associated with tra wl-sampled specimens.T he proportion of females in spa wning and spent stages is used as an indicator of surv e y timing relative to spawn timing.

Figure 4
Figure 4The proportion of females > 30 cm FL in spawning or spent stages (logit transformed; SP 30_wt ) in relation to the estimated number of days between the winter Shelikof Strait AT survey and the median date of spawning as estimated from larval data ( Mismatch med ).

−Figure 5
Figure 5 Residuals (a) from the stock assessment model for the winter Shelikof Strait AT survey were used as a proxy for the relative availability of pollock (age 3 + ) to the surv e y.T he relationship betw een residuals and three catchability co v ariates is sho wn: (b) the mismatch in timing betw een the surv e y and estimated median spawn date; (c) the logit-transformed proportion of females spawning or spent; and (d) March SST.Estimated linear regressions with 95% confidence intervals are shown.

Figure 6
Figure6Observed (black points and vertical bars) and estimated (colored lines) winter Shelikof Strait AT survey biomass (ages 3 + ) for assessment models with different surv e y catchability formulations.The assessment model used for tactical advice currently uses a constrained random walk (RW) with fixed σ q = 0.038, labeled "fixed RW." Estimates are shown for models with constant catchability ("constant q"), RW with estimated process error ("RW"), and combinations of single covariates with and without a RW term.
Downloaded from https://academic.oup.com/icesjms/advance-article/doi/10.1093/icesjms/fsae005/7591717 by NOAA Central Library user on 15 February 2024 From 1983 to 2021, a total of 27 486 female pollock > 30 cm were visually assessed for spawning and maturity state during winter Shelikof Strait AT surveys ( n = 18 155 since 1992).The six analytical approaches to estimating the proportion of female pollock in spawning or spent stages resulted in time series of spawning status that were highly correlated with

Table 1 .
Downloaded from https://academic.oup.com/icesjms/advance-article/doi/10.1093/icesjms/fsae005/7591717 by NOAA Central Library user on 15 February 2024 Changes in the timing of spawning and the availability of walleye pollock to assessment surveys 9 Results of single-covariate regression models to explain the difference between measured pollock biomass in the winter Shelikof AT survey (age 3 Downloaded from https://academic.oup.com/icesjms/advance-article/doi/10.1093/icesjms/fsae005/7591717 by NOAA Central Library user on 15 February 2024

Table 2 .
Comparison of stock assessment models with different catchability formulations and covariates, including const ant catchabilit y (const ant q), an estimated random walk (RW), a linear covariate effect (covariate), or a RW plus a covariate (RW + covariate).
Table rows correspond to catchability formulations (1)-(4) in the main text.The AIC is the difference in AIC from the best supported model for each covariate tested, with a AIC of zero indicating the best supported model for that covariate.The % reduction in variance for the RW process that is achieved by adding the given covariate is shown.