Fishing, predation, and temperature drive herring decline in a large marine ecosystem

Abstract Since 1960, landings of Atlantic herring have been the greatest of any marine species in Canada, surpassing Atlantic cod and accounting for 24% of the total seafood harvested in Atlantic Canada. The Scotian Shelf‐Bay of Fundy herring fisheries (NAFO Division 4VWX) is among Canada's oldest and drives this productivity, accounting for up to 75% of the total herring catch in some years. The stocks’ productivity and overall health have declined since 1965. Despite management measures to promote recovery implemented since 2003, biomass remains low and is declining. The factors that drive the productivity of 4VWX herring are primarily unresolved, likely impeding the effectiveness of management actions on this stock. We evaluated potential drivers of herring variability by analyzing 52 time‐series that describe the temporal and spatial evolution of the 4VWX herring population and the physical, ecological, and anthropogenic factors that could affect them using structural equation models. Variation in herring biomass was best accounted for by the exploitation rate's negative effect and the geographic distribution of fishing and recruitment. Thermal phenology and temperature adversely and egg predation positively impacted the early life stage mortality rate and, ultimately, adult biomass. These findings are broadly relevant to fisheries management, but particularly for 4VWX herring, where the current management approach does not consider their early life stage dynamics or assess them within the ecosystem or climate change contexts.


| INTRODUC TI ON
The heavy exploitation combined with the socio-economic and ecological importance of forage fish globally has raised concerns about their conservation status (FAO, 2014). Hypotheses to explain forage fish population fluctuations and collapses have been diverse: heavy exploitation (Essington et al., 2015), high, early life stage mortality (Fassler et al., 2011), intrinsic density-dependent population regulation (Myers & Barrowman, 1996), predation (Kotterba et al., 2017;Richardson et al., 2011), and environmental variability (Brosset et al., 2018;Trochta et al., 2020). Traditional management approaches focus heavily on the role of exploitation in regulating recruitment and adult biomass and largely overlook the importance of ecosystem or climate factors. However, some forage fish stock collapses have no apparent explanation or have occurred naturally (Peck et al., 2014). The inability of exploitation alone to explain forage fish population collapses (McClatchie et al., 2017;Peck et al., 2014) or delayed recoveries (Essington et al., 2015) suggests that additional factors are critical in accounting for population variability.
Atlantic herring (Clupea harengus), one of Canada's oldest and largest fisheries, has been experiencing a long-term decline over the broad Scotian Shelf and Bay of Fundy region (NAFO Division 4VWX; Figure 1a,b). It accounted for 19% of Canada's total seafood landings between 1990 and 2018, the largest of any species ( Figure 1a). However, the productivity and population health of  Trochta et al. (2020), the decline in 4VWX herring stock status is exceptional: it has been in a collapsed state for 26 years , far longer than the average collapse duration of 11 years. The causes of the long-term decline in herring state ( Figure 1d) and the failure of the stock to respond to reduced exploitation are unknown, impairing management and conservation efforts. The cultural, economic, and ecological importance of herring within this region provides a strong incentive to understand the factors which underlie its declining status.
The 4VWX herring stock undergoes a complex seasonal cycle of spawning, overwintering, and feeding that involves separate geographic domains and differential mixing with neighboring populations. Mark and recapture studies indicate that adult herring (>2-3 years old) spawn near German Bank/Lurcher Shoals in August-November, overwinter ~700 km to the northeast in Chedabucto Bay in January-March, and feed on the Scotian Shelf and Bay of Fundy in April-July (Stobo & Fowler, 2009).
Upon hatching and before metamorphosis, larval herring remain within a well-defined, tidally-mixed "larval retention area" (LRA) off SWNS (Iles & Sinclair, 1982;Stephenson et al., 2015). Juvenile The 4VWX herring fishery targets spawning and nonspawning adult herring and juveniles (sardines) using multiple fishing gears, primarily purse seines but also weirs and gillnets. From 1985 to 2014, 69% of the effort was directed mainly to harvest roe during the spawning season (Singh et al., 2016). Unlike some roe fisheries, such as Pacific herring, which harvest roe after it is shed from females, adult herring are captured and then processed for their roe. The fishery is managed as four components, with the dominant Southwest Nova Scotia/Bay of Fundy (SWNS/BoF) component accounting for up to 93% of all landings in some years (DFO, 2015). Of these, the SWNS spawning component has directly accounted for a far greater fraction (44%) of the total herring landed annually in Division 4VWX  relative to the upper BoF (10%). Accordingly, due to its greater productivity and availability of datasets, our evaluation and discussion of 4VWX herring focused on the SWNS spawning herring complex.
Here, we evaluated how the physical, ecological, and anthropogenic conditions individually and synergistically affect the population dynamics of 4VWX herring across different life stages. To achieve this goal, we assembled time-series of the herring stock's intrinsic population characteristics across different life stages and coincident series of climate, weather, chemistry, biology, exploitation, and ecological dynamics (predation, prey availability, competition) that could plausibly influence population productivity. We used multivariate statistical approaches, including multi-model inference and structural equation models, to explore hypotheses for how these factors interact and affect herring population variability. We interpret our findings in the context of forage fish dynamics and suggest ways to improve forage fish management, particularly 4VWX herring.

| Overview of approach
The herring stock has a long, detailed assessment history based on surveys designed to collect larval and adult herring and collation of landings data from each commercial fisheries sector (e.g. DFO, 2015). Various research initiatives have generated knowledge about distribution, growth, and maturity at different life stages. We used biologically relevant indicators from several data sources to define the spatial and temporal variability required to explore the drivers of herring variability across the different life stages. Before testing hypotheses, individual time-series were standardized to common spatial and temporal resolutions and measurement units for 1975-2005. The preprocessing steps were (1) defining the spatial, seasonal, and temporal domains, (2) data acquisition and calculation of the indicator time-series, (3) imputation, normalization, and standardization of time-series, and (4) statistical analyses.
2.1.1 | Defining the spatial, seasonal, and temporal domains The geographic boundaries used to define larval and adult herring habitat were identified by Boyce et al. (2019). The habitat domain of 4VWX larval herring (the larval retention area; LRA, see Figure 1b) was estimated using field observations of their abundance from the Bay of Fundy larval herring survey and particle transport simulations from WebDrogue (Hannah et al., 2001). The spatial domain of adult herring in 4VWX was previously defined by Boyce et al. (2019) using survey observations, tagging studies, and landings statistics and includes most of the 4VWX domain. Herring in this area could originate from other spawning complexes located on the offshore banks or in the United States and represent an unknown mixture of stock sub-components.
The biophysical data series that potentially govern 4VWX herring productivity at early life stages were restricted to the LRA and the July 15-October 31 spawning window. The biophysical data series used to explore adult herring dynamics were restricted to the adult spatial domain. We present available data for the 1970-2010 period; due to limitations of some of the series, our statistical modeling analyses concentrate on the 1975-2005 period. Much of the variation in herring population state over the past 50 years occurred during these periods (Boyce et al., 2019;DFO, 2018;Trochta et al., 2020), and these years have the most comprehensive data coverage.

| Data acquisition and indicator time-series
We assembled data that included spatial and temporal variation in 4VWX herring across different life stages, the exploitation dynamics, and the biological and environmental factors that may influence herring. Data sources related to changes in the mean states, seasonal dynamics (timing, amplitude), and community composition of the plankton or larval assemblages that herring interact with, the type and intensity of fishing pressure, predation and competition were obtained or calculated from the stock assessments, peerreviewed publications, at-sea surveys, ships of opportunity, and remote sensing are listed in Table 1. From these data, 102 indicator time-series consisting of annual observations that could be averaged over a spatial domain (e.g., Lurcher Shoal SST) or from a single site (e.g., St. Andrew's SST), were initially developed. However, series missing >35% of their values (i.e., more than ten annual data points)  Since independent data sources were unavailable, three of the five intrinsic series (spawning stock biomass; SSB, recruitment; r, and recruitment rate) combine observations (e.g. industry and DFO acoustic surveys) with an assessment model. The SSB series was calculated from calibrated VPA and acoustic estimates using methods described in Boyce et al. (2019). The recruitment series (abundance at age 1) were derived from the VPA-based 4VWX stock assessments (Power et al., 2006;Singh et al., 2015). The recruitment rate was calculated using the approach of Platt et al. (2003) as the recruitment is standardized by SSB 3 years prior. Dahlke et al., 2020) and has been found to yield significantly less bias in subsequent analyses than case-wise deletion (Ellington et al., 2015). Statistical descriptions of MICE are readily available (Azur et al., 2011;Patrician, 2002;Schafer, 1999;Slade & Naylor, 2020 Figure S1).

| Statistical analyses
Understanding how the 52 indices interact to affect the temporal variability in herring productivity is complex, involving the consideration of lagged effects, collinear predictors, and numerous permu- Structural equation models were estimated as a network of interacting linear models within which variables can function as both predictors and responses and where relationships between unobservable (latent) processes of interest can be estimated. SEMs are valuable for distinguishing between processes that are of interest but cannot be directly measured or observed (latent constructs) from measurements that are useful but imperfect proxies for these processes (observed variables). Whereas traditional statistical models can automatically search for and evaluate correlative relationships between predictors and the response variable, by requiring a priori model specification and enabling multiple pathways between variables to be simultaneously tested, SEMs facilitate a more careful and deeper consideration of causation. SEMs are becoming increasingly common in ecology (e.g. Boyce et al., 2017) and can lead to a more rigorous causal inference network than can be achieved with traditional linear models or correlative approaches (Pearl, 2009). They also permit an accounting of time-lagged effects and multi-collinearity. Accordingly, SEMs can also yield different results than correlative statistical approaches (Pearl, 2009). We estimated the SEM models using the complete time-series of imputed data   insensitive to the use of imputed versus unimputed values. Figure   S2 illustrates the steps and workflow to calculate and analyze the indicator time-series.

| Temporal variability of indicators
The 52 standardized and imputed indicator time-series varied considerably between 1970 and 2010, featuring a broad division of variance between long and short periods (Figure 2). For example, we found that herring SSB has ~80% of its variance in periods >8 years
The geographic distribution of fishing was also a significant predictor of SSB, with a more geographically dispersed fishing footprint associated with higher SSB. Over the long term, the geographic foot- County accounted for 50% of the total 4VWX fishery and after that declined to their current level of 30%. Alternatively, from their minimum in ca. 1920, herring landings in Yarmouth County, primarily a roe fishery, accounted for 5% of the total 4VWX fishery, rose to about 35% in ca. 1985; after a slight decline, it currently accounts for ~55% of all landings (green in Figure 4c, Figure S6). Yarmouth County is adjacent to the 4VWX spawning locations, where herring aggregate in dense predictable concentrations and are subject to high mortality during spawning. The timing of this apparent shift from a juvenile fishery to one that primarily targets spawners for their role in the mid-1980s also coincides with the timeline of declining herring status (Figure 1d). Taken (Brosset et al., 2018;Johannessen, 1980;Kotterba et al., 2017;Richardson et al., 2011). While the environmental and biological factors did not strongly affect SSB, indices related to the temperature, timing of the seasonal temperature development (phenology), and egg predation by haddock significantly influenced both the recruitment rate and absolute recruitment, which then affected adults ( Figure 3).
The strong, consistent, and positive effects of haddock predation on the survivorship of early life stage herring (recruitment rate) were unexpected. We evaluated this relationship's robustness through sensitivity analyses and assessed numerous alternative hypotheses to explain it but could not identify any. However, several lines of evidence suggest the positive effect of haddock predation on herring may be valid. Egg mortality and predation (Johannessen, 1980), particularly by haddock (Bowman, 1922;Richardson et al., 2011;Toresen, 1991), has been emerging as an essential driver of herring population productivity in other ecosystems. Predator exclusion experiments in the Baltic Sea reported that predation accounted for 42% of herring egg mortality (Kotterba et al., 2017), while in a Norwegian fjord, Atlantic cod (Gadus morhua) alone were reported to consume 40-60% of spawned herring eggs (Johannessen, 1980).
The importance of egg predation in our analyses is also consistent with results from a study that reported a decoupling of herring SSB and larval density on the adjacent Georges Bank due to the interaction between larval predation by haddock and exploitation pressure (Richardson et al., 2011). The geographic distribution of haddock in July further reinforces the hypothesis that they are important predators of herring eggs. Using standardized survey observations (1970-2018), we found that the areas where haddock reach peak abundances (>95th percentile) overlap substantially with the spawning areas of 4VWX herring and with Georges Bank (Figure 5A), where significant egg predation by haddock has occurred (Richardson et al., 2011).
Contrary to the findings of Richardson et al. (2011) andKotterba et al. (2017), which reported adverse effects of egg predation on adult production, we observed a positive effect of egg predation on recruitment rate, which then propagated to recruitment and ultimately SSB lagged by three years (Figures 3 and 5). Both SEMs ( Figure S5) and univariate relationships ( Figure 5) suggested that the positive effect of egg predation on the recruitment rate (0.39, Figure 3) may operate by regulating the density of herring eggs (high predation →low larval density →greater larval size →higher recruitment rate). The relationships between haddock predation and larval density (r = −.49), larval length (r = .59), recruitment rate (r = .37), and SSB lagged by three years (r = .42) could affect recruitment rate via several pathways (Figure 5b-e). First, with all else equal, egg consumption by haddock leads to fewer eggs hatched into larvae (reduced larval density), and thus, to greater resources (planktonic prey) per individual larvae, contributing to larger larvae. Larger larvae are, in turn, better able to capture prey (Blaxter, 1962), avoid predation, and maintain position within the water column, thereby enhancing larval retention near the spawning locations where conditions are favorable (Frank, 1988;Iles & Sinclair, 1982;Stephenson et al., 2015), leading to increased survivorship. Second, herring eggs are spawned on a gravelly bottom substrate in dense aggregations that can lead to hypoxia-driven mortality; hypoxia can also reduce hatching success by 24%-80% and lead to reduced larval sizes (DePasquale et al., 2015). Egg consumption by haddock can reduce egg density, thus reducing the extent and severity of hypoxia experienced by herring eggs, potentially increasing the hatching success of the remaining eggs and the size of subsequent larvae. Lastly, haddock egg predation can impact the deposition pattern of herring eggs, significantly affecting hatching distribution and larval size at hatching (Munk & Tosenthal, 1983). These chains of interactions from eggs to larvae, juveniles, and adults under low and high egg predation are depicted quantitatively (Figure 5b-e, and Figure S5) and schematically (Figure 5f,g).
The effect of predation on herring eggs is likely nonlinear, depending on herring eggs' concentration and predation intensity.
For example, low levels of egg predation may yield positive effects, while at some higher levels of predation, adverse effects on recruitment rate and recruitment would occur. It is possible that over the core focal period of consideration in this study , the F I G U R E 3 Structural equation model effects between herring population variability across different life stages and environmental and anthropogenic processes that explain it. Observed variables are depicted as smaller gray circles and are used to understand the unobservable latent processes (larger colored circles). Arrows depict the directed model relationships. Negative effects between unobservable latent processes are shown in red and positive effects in blue; the strength of these effects is in units of variance change in the response per unit variance increase in the predictor. Model effects are displayed with an asterisk denoting statistical significance (*p < .05, **p < .01, ***p < .0001). Latent processes related to intrinsic herring population status are black circles, while those related to environmental, anthropogenic, and ecological factors are turquoise circles   (Darbyson et al., 2003) may be important. We also observed an interaction between egg predation and exploitation: the largest SSB levels were predicted at the lowest exploitation rates and highest egg predation. Despite the observed effects of egg predation by haddock, the direct predation on adult herring by dogfish (Squalus acanthias), pollock (Pollachius pollachius), silver hake (Merluccius bilinearis), white hake (Urophycis tenuis), and Atlantic cod (Gadus morhua) was weak ( Figure S4) and nonsignificant ( Figure 3). However, due to the complex life history and many ecosystem interactions of 4VWX herring and the incomplete time-series, the haddock egg predation index's importance may have arisen through an alternate pathway that we could not evaluate. Due to its novelty and potentially overarching impacts on herring SSB, the effect of haddock predation on early life stage herring merits further scrutiny and should be a priority for future study.
Most ocean temperature time-series suggested warming from the 1970s throughout the area occupied by herring, particularly during autumn across the LRA. These temperature changes were negatively related to the early life stages, that is, to recruitment rate and the number of recruits but did not have a significant direct relationship with SSB ( Figure 3). This finding is consistent with the reports that herring may have lower thermal tolerances when spawning and during their early life stages (larvae, juveniles) than when adults (Fassler et al., 2011;Payne et al., 2009). Due to their greater mobility and range, adult herring can also better avoid extreme warm temperatures than larvae. The more substantial temperature effects on early life stages of herring than on adults are consistent with studies that have reported narrower thermal ranges during early life stages (embryos and larvae) and reproduction (Dahlke et al., 2020;Portner & Farrell, 2008;Portner & Peck, 2010).
Bioenergetics models suggest that temperature may impact herring larvae at a critical period, during yolk-sac absorption and first feeding (Hufnagl & Peck, 2011); further, experiments suggest that such temperature increases could lead to a reduced size of newly hatched herring (Ware, 1975). Therefore, it is possible that the increasing temperatures are contributing to the reduced size of herring eggs and larvae with consequent effects on their fitness and survivorship. Indeed, warming trends over the study area support this hypothesis. Although annual average temperature trends have been moderate since 1970 (Figure 2), warming trends during the fall spawning period for herring when annual temperatures reach a maximum were stronger ( Figure S7a). For example, in the LRA, the average peak SST during the herring spawning season was 13.3°C in the 1980s but increased by 3.4-16.7°C in 2012 ( Figure S7b). Further, as these peak SST values are averages, observed temperatures are higher in some locations and days. Based on the reported optimal temperature for herring larvae in the eastern Atlantic derived from experiments and field studies (Moyano et al., 2020), it is likely that temperature conditions are becoming increasingly stressful to larvae and spawning adults in the LRA. The cardiac function of larval herring declines at temperatures above 16°C (Moyano et al., 2020), while growth rates have been reported to decline at temperatures above 17°C (Moyano et al., 2020) or 17.5°C (Fey, 2001). The maximum predicted temperature value for any single year and location ( Figure S7a Table S1). This disproportionate warming during the autumn may partly explain why the state of autumn spawning herring is declining (Boyce et al., 2019), whereas spring spawning herring in the nearby inner Bay of Fundy is not (DFO, 2015). If this warming trend continues, herring Proportion of 4VWX landings 19104VWX landings 19204VWX landings 19304VWX landings 19404VWX landings 19504VWX landings 19604VWX landings 19704VWX landings 19804VWX landings 1990  Charlotte county Yarmouth county larvae and recruitment will likely become increasingly affected, and adults may soon experience direct physiological stress from these autumn temperature increases. The seasonal warming trends may also induce a shift in the phenology of autumn spawning herring or plankton development, leading to a mismatch between larval herring and their food supply (Cushing, 1990;Platt et al., 2003). In addition to temperature, the size of herring eggs can also be affected by the size and fitness of spawners (Blaxter & Hempel, 1963;Óskarsson et al., 2019;Ware, 1975), which has steadily declined since the 1980s (Boyce et al., 2019).
In undertaking this analysis, we evaluated numerous factors that could plausibly underlie the observed variability in 4VWX herring stock dynamics. Despite this, due to limited data availability and some incomplete time-series, it was not possible to examine all variables that might reasonably play a significant role in this stock.
For example, the zooplankton series that we compiled contained too many missing values to be included. It is possible that changes in the availability, timing, or quality of zooplankton prey may have affected herring dynamics, particularly during the early life stages (Brosset et al., 2018). Further, we were unable to assess the effect of some crucial predators such as whales and seabirds on herring.
Notwithstanding these caveats, this analysis emphasizes the overarching importance of exploitation as a regulator of herring productivity. However, it also highlights the importance of consid-  (Baum & Fuller, 2016;. Adapting fisheries management to include climate and ecosystem factors is a high priority for fisheries agencies worldwide and an objective within Fisheries and Oceans Canada. Despite this, analyses suggest that few fisheries in Canada are currently including these factors. A recent review reported that only 11% of 729 fisheries assessments in Atlantic Canada and the Eastern Arctic published between 2000 and 2020 mentioned climate change . There is a growing suite of approaches and methods for furthering the consideration of climate change in fisheries management (e.g. Busch et al., 2016;Gattuso et al., 2015;Lawler et al., 2010;Ojea et al., 2017;Pinsky & Mantua, 2014). For example, climate vulnerability assessments (Greenan et al., 2019;Stortini et al., 2015) have been widely advocated as an approach for increasing understanding of species and fisheries vulnerability to climate changes and deploying climate adaption resources (Barange et al., 2018;Busch et al., 2016;Hare et al., 2016). Management strategy evaluations can optimize harvest rules to ensure that they are robust to future climate scenarios, population and ecosystem dynamics and other uncertainties. Dynamic management can set harvest rates based on near real-time forecasts or respond to rapidly changing conditions (Dunn et al., 2016). The US National Oceans and Atmospheric Administration employs ecosystem models that include multiple species interactions and environmental effects to address the impact of exploitation and climate changes on the dynamics of exploited species (Holsman et al., 2017) include temperature-dependent weightat-age functions and temperature-specific predation interactions.
Integration of these approaches into scientific advice would be ideal.
For example, the Alaska Eastern Bering Sea Integrated Ecosystem assessment program employs climate forecasts and projections developed by regional ocean modeling systems, food web and multispecies assessment models, and scientific surveys to support and inform fisheries decision-making in the North Pacific (NOAA, 2021).
In addition to incorporating climate considerations, our findings also emphasize the critical importance of adopting precautionary management principles. During the 1970s and 1990s, ecosystem conditions (e.g., ocean temperature, haddock predation) were shifting, the fishery's geographic distribution was contracting, the fishery's nature was changing, the number of spawning locations was reduced, and the productivity and population health of 4VWX herring declined. These conditions served to heighten uncertainty over the herring stock's status. Overall, our results agree with recent studies that emphasize the critical importance of anthropogenic, climate (Trochta et al., 2020), and ecosystem (Kotterba et al., 2017;Richardson et al., 2011) factors in determining the early life stage dynamics and emergent adult biomass of Atlantic herring. These findings suggest that a more comprehensive ecosystem approach must be considered.

ACK N OWLED G M ENTS
The authors are grateful to all data providers and Rabindra Singh

CO N FLI C T S O F I NTE R E S T
The authors declare no conflict of interest. Writing -review & editing (equal).

DATA AVA I L A B I L I T Y S TAT E M E N T
The data used in this study are available through the Dryad digital repository (doi: https://doi.org/10.5061/dryad.gtht7 6hnm).