Influence of larval transport and temperature on recruitment dynamics of North Sea cod (Gadus morhua) across spatial scales of observation

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited. © 2020 The Authors. Fisheries Oceanography published by John Wiley & Sons Ltd Romagnoni and Kvile contributed equally to the article. 1Department of Biosciences, Centre for Ecological and Evolutionary Synthesis (CEES), University of Oslo, Oslo, Norway 2Department of Biology, Woods Hole Oceanographic Institution, Woods Hole, MA, USA 3Norwegian Meteorological Institute, Bergen, Norway 4Norwegian Institute for Water Research (NIVA), Oslo, Norway 5Department of Natural Sciences, Centre for Coastal Research (CCR), University of Agder, Kristiansand, Norway

variability in stock-recruitment models can help to explain parts of the observed recruitment variability. In fact, it has been suggested that ELS transport could be one of the drivers behind the unclear relationship between spawning stock size and recruitment (Huwer et al., 2016). While some studies have included proxies of larval transport in stock-recruitment relationships (Baumann et al., 2006;Zimmermann, Claireaux, & Enberg, 2019), few have included direct estimates of larval transport (but see Hidalgo et al., 2019).
Moreover, transport can influence connectivity among populations (e.g. through interannual variability in oceanographic current patterns, Huwer et al., 2016;Kvile, Romagnoni, Dagestad, Langangen, & Kristiansen, 2018) and recruitment dynamics across large geographic scales (Cadrin, Goethel, Morse, Fay, & Kerr, 2019;Henriksen et al., 2018;Hinrichsen, Von Dewitz, & Dierking, 2018) and thereby population management (Fogarty & Botsford, 2007;Hidalgo et al., 2019;Ramesh, Rising, & Oremus, 2019). Critically, the spatial scale of observation can affect the stock-recruitment relationship, providing contrasting results across scales (Chang, Chen, Halteman, & Wilson, 2016). The importance of environmental drivers for recruitment can also differ across subunits within a stock (Brosset et al., 2018). We therefore expect the importance of larval transport variability for recruitment to differ across spatial scales (i.e. basin vs. sub-basin) and between individual populations. In this study, we explore alternative approaches to explicitly include larval transport in stock-recruitment functions, and quantitatively assess the effect of transport on recruitment across spatial scales. We use annual estimates of larval retention and population connectivity, obtained through a coupled physical-biological model of larval drift, focusing on North Sea cod (Gadus morhua) as a case study.
We initially include retention anomaly as an environmental covariate in traditional parametric stock-recruitment model formulations and compare its effect to alternative covariates, namely sea surface temperature (SST) and the North Atlantic Oscillation index (NAO). In addition, we propose a novel approach for inclusion of the effect of larval transport in stock-recruitment models by weighting spawning stock biomass (SSB) according to yearly retention and advection rates, providing a measure of "effective biomass." Two alternative approaches are proposed to account for effective biomass: including only retention in the spawning area of origin (retention-only SSB, rSSB), and including retention and inflow of larvae from other areas (net drift SSB, ndSSB), effectively accounting for connectivity. Additionally, we quantify temporal patterns in population connectivity and their relationship with SST and NAO.

| Ocean circulation in the North Sea
Ocean circulation in the semi-enclosed North Sea basin is influenced by topography and inflow of North Atlantic water, separating the basin into a shallow southern and a deeper northern area.
The northern area is influenced by inflow of saline Atlantic water flowing along the western slope of the Norwegian Trench. This current transports the planktonic copepod, Calanus finmarchicus, an important food source for larval cod (Nicolas, Rochette, Llope, & Licandro, 2014) and other species, into the region. The current flows along the Norwegian Trench and into Skagerrak, where it enters the "Skagerrak loop." It follows a counterclockwise trajectory along the Skagerrak coast, and after mixing with the less saline Norwegian coastal current, flows north-westward along the eastern slope of the Norwegian Trench and into the Norwegian Sea (Huserbråten, Moland, & Albretsen, 2018). The southern North Sea is dominated by continental freshwater run-off and tidal patterns, which in combination with wind and wave turbulence and shallow topography result in permanent mixing. The intermediate saline current from the English Channel and the coastal, low saline Jutland Current flow along the continental coast and into the Skagerrak, entering the "Skagerrak loop" (Sundby, Kristiansen, Nash, & Johannessen, 2017).

| Cod populations in the North Sea
Although managed as one stock (ICES, 2018c), North Sea cod comprises a number of spatially segregated units, with limited overlap and varying degree of connectivity (Heath et al., 2014;Neat et al., 2014).
The main units are the Viking and the South populations. The latter is often separated into a South proper (centred around the Dogger Bank) and a Northwest unit; these two subpopulations are genetically homogenous but show contrasting demographic trends and limited adult connectivity, so their relationship is as yet unclear (Neat et al., 2014).
In this study, we considered alternative scenarios with three populations (Viking, South and Northwest), two populations (Viking and South including Northwest) and a single population (aggregating Viking, South and Northwest), the latter roughly corresponding to the current management unit. The populations' spatial extent (Figure 1) was based on ICES (2015). We calculated larval connectivity between the populations and assessed drift into the Skagerrak (which is excluded from our populations) and the Norwegian Sea, Scottish sea and English Channel (hereafter called "outside"; Figure 1). Particles leaving the study area (i.e. entering the "outside" area) were considered lost.

| Early-Life Stage (ELS) dispersal model
To quantify larval retention and connectivity between populations, we used a coupled physical-biological model (hereafter, ELS dispersal model) for the time period 1971-2014 and included the model output in statistical stock-recruitment models for the same years. The individual-based model (IBM) simulates development and transport of cod eggs and larvae based on earlier studies of larval cod (Kristiansen, Lough, Werner, Broughton, & Buckley, 2009;Kristiansen, Stock, Drinkwater, & Curchitser, 2014;Kristiansen, Vikebø, Sundby, Huse, & Fiksen, 2009). The IBM is integrated as a module in the open source Lagrangian particle tracking framework OpenDrift (github. com/opendrift; Dagestad, Röhrs, Breivik, & Ådlandsvik, 2018;Kvile et al., 2018), and the code for the cod eggs and larvae module is available on github.com/trond kr/KINO-ROMS/tree/maste r/Romag noni-2019-OpenD rift. To simulate transport with ocean currents and temperature-dependent development, the IBM was coupled offline to a reanalysis of the regional ocean circulation model ROMS (Shchepetkin & McWilliams, 2005) configured for ocean regions covering the Nordic Seas (including the North Sea) and parts of the Arctic Ocean, with 4 km horizontal resolution, 32 vertical layers and output stored daily (Lien, Gusdal, Albretsen, & Melsom, 2013).
For downloading options, see http://thred ds.met.no/thred ds/ nansen_daily.html. Further details on the characteristics and limitations of the ELS dispersal model are available in Kvile et al. (2018).
Due to long-term and interannual variation in the relative importance of spawning grounds (González-Irusta & Wright, 2016;Sundby et al., 2017) and the uncertainty in spawning ground locations early in the time series, we released particles representing cod eggs uniformly within the three populations' spatial extent ( Figure 1).
Although this could reduce the precision of connectivity estimates in some years, we considered this approach as more conservative when modelling larval transport over a long time period including years with unknown spawning ground distribution. To obtain uniform spatial distribution (0.12-0.14 eggs/km 2 ), we set the number of eggs released based on the sizes of the population areas: ~32,400 in the South (~270,000 km 2 ), ~22,950 in the Northwest (~170,000 km 2 ) and ~27,000 in the Viking area (~200,000 km 2 ), for a total of ~91,500 eggs.
We defined the timing of egg release using prior knowledge of the population spawning periods (Brander, 1994(Brander, , 2005Fox et al., 2008): between December 15th and April 15th for the South population, between January 1st and May 1st for the Northwest population and between February 1st and May 15th for the Viking population ( Figure 2a). The number of eggs released per day followed a Gaussian distribution, N(μ, σ 2 ), where μ = 1 and σ = 0.25, scaled to the length of the spawning season and the total number of particles defined per population area and with peaks that approximately matched the spawning peak described by Brander (1994). Setting a broader spawning season than observed in recent years accounts for uncertainty in the spawning season early in the time series. For example, spawning was allowed to start in December for the South population to account for the fact that the Southern Bight component, which spawns earlier than the German Bight and Central-west (Brander, 1994), was more abundant in the past. Eggs were released in equal numbers at 10 m depth intervals between 0 and 50 m (i.e. for a given population, an equal number of eggs was uniformly released at 0, 10, 20, 30, 40 and 50 m). After release, eggs and larvae were advected horizontally at fixed depths using an Euler interpolation scheme without horizontal diffusion and a 1-hr time step. The Euler scheme differed minimally compared to a more computationally costly Runge-Kutta scheme (Kvile et al., 2018). We used different drift depths to represent vertical movement within the depth range typically available in the North Sea, based on the finding that incorporating a more computationally costly vertical movement behaviour had limited effect on connectivity and retention of cod ELS at settlement in the North Sea (Kvile et al., 2018). Development time of planktonic eggs (d, days) was a function of the ambient sea water temperature (T, °C) according to the ocean model reanalyses, parameterised based on observations for cod eggs (Langangen, Stige, Yaragina, Vikebø, et al., 2014, based on data in Ellertsen, Fossum, Solemdal, Sundby, & Tilseth, 1987, Figure 2b): After completing the egg stage, the simulated individuals hatch into cod larvae. The simulated cod larvae grew with a growth rate (GR, percentage of larval weight/day) depending on larval weight (W, mg) and ambient temperature (T), as estimated experimentally for Atlantic cod larvae (Folkvord, 2005) (Figure 2c): F I G U R E 2 Cod egg and larvae IBM functions. (a) Number of eggs released per Julian day (counting from January 1st) per spawning ground; (b) egg development time (d) as a function of temperature; (c) growth rate of larvae day -1 (GR, contours) as a function of larval weight and temperature; (d) larval length (L) as a function of weight; and (e) mortality rate (m, day -1 ) for eggs (fixed at 0.2) and larvae as a function of larval weight Larvae were assumed to feed ad libitum, and their initial weight was set at 0.08 mg. Larval length (L, mm) was a function of weight (Folkvord, 2005) (Figure 2d): We assumed that cod larvae had no directional horizontal (swimming) movement. During the simulation, eggs were subject to a fixed daily mortality rate (m) of 0.2, which is within the range of mean values estimated in studies of cod eggs (0.1-0.32, Rijnsdorp & Jaworski, 1990; see Table 2 in Langangen, Stige, Yaragina, Vikebø, et al., 2014). For larvae, we set the mortality rate to decrease with weight ( Figure 2e) as parameterised for North Sea cod larvae in Akimova, Hufnagl, Kreus, and Peck (2016), based on the size-spectrum theory (Peterson & Wroblewski, 1984): The survival probability of each individual was updated throughout the simulation according to the mortality rate (i.e. individuals were not removed from the simulation), following a super-individual approach (Scheffer, Baveco, DeAngelis, Rose, & van Nes, 1995).
Larvae settled when reaching a length >49 mm (Bastrikin, Gallego, Millar, Priede, & Jones, 2014). Only larvae settling within known nursery areas for North Sea cod (based on Heath et al., 2014; see Figure 1) were considered to successfully settle and survive; larvae that reached settlement length outside nursery areas were considered dead (hereafter "not settling"). Larvae not reaching settlement length by the end of the simulation (set to 15th August for South and 29th September for Northwest and Viking) were considered dead (amounting to <1% of larvae, not included in the analysis).
The juvenile stage was not simulated since cod adopt a demersal lifestyle upon reaching settlement length.
For each population, we estimated the proportion of larvae (a) retained in a nursery area for the given population of origin; (b) drifting into the nursery area of another population; (c) drifting out of the study area (to the Skagerrak or "outside") and (d) reaching settlement size within any population area, but not within a nursery area ("not settling"). Annual values  for these metrics were included in the stock-recruitment analysis (see below). To test the robustness of the results of the stock-recruitment analysis to key assumptions in the larval dispersal model, we performed additional simulations where the mortality rate (for eggs and larvae) was adjusted by ±20% and separate simulations where settlement size was adjusted by ±20%. We ran these additional simulations for 1990 and 2010, two years with different climatic conditions (high and low NAO phase, respectively; Figure 3) and contrasting results of larval dispersal. Parameters included in the ELS model are summarised in Table 1.

| Observational data
We calculated population-specific estimates of SSB and recruit-  Figure 3). Although more accurate abundance estimates could be obtained by using standardised indices instead of raw data, these are only available from 1983. We instead used raw data to include a longer time series, spanning years when SSB was higher, to provide robustness to the stock-recruitment estimates.
We generated annual data of age-specific abundance to match the ICES statistical sub-rectangle using Catch-Per-Unit-Effort (CPUE) adjusted for swept area and gear catchability, assuming that the sample is representative of fish abundance. The swept area is a function of standardised tow length and net width during towing, which in turn is a function of tow depth (ICES, 2012). Average depth per sub-rectangle (NGDC, 1995) was used as a proxy for tow depth. Catchability coefficient for the survey gear by age (Fraser, Greenstreet, & Piet, 2007) was multiplied by the ratio of the swept area to the whole ICES rectangle area. SSB was calculated as abundance per age multiplied by maturity and weight per age and year, assumed homogeneous between populations lacking population-specific data (ICES, 2017). The recruitment index was estimated by back-calculating abundance of age 1 from averaged age 2 and 3 abundances scaled by the age-and year-specific natural mortality (ICES, 2015). As climate variables, we included sea

| Stock-recruitment models
We used the Cushing parametric stock-recruitment model formula- Similarly to the SSB-based models, models with ndSSB and rSSB were fitted with or without climate variables as covariates ( Table 2).
The models thus took the form (Akimova, Núñez-Riboni, et al., 2016): Where recruitment R was calculated as a function of the generic S (either SSB, rSSB or ndSSB). This was extended for inclusion of climate variable E as: E was any climate variable (SST, NAO or retention anomaly, RA).
RA was calculated as the annual deviation from the mean larval retention over the whole time series for a given population, determined from the ELS dispersal model.
The linear forms of the models (see Appendix S1) were tested for residuals assumptions, and outliers were removed from the analysis. Commonly used model comparison methods such as AIC and likelihood-based approaches could not be used since models with SSB, ndSSB or rSSB included different data in the predictor variable.
Models were therefore compared through their absolute fit to data using adjusted R 2 , with significance threshold set at .05. Adjusted R 2 allows highlighting the combinations of predictor and covariates with highest explanatory power, that is those that improve the model more than expected by chance, with penalisation of additional parameters.

| RE SULTS
We compared the performance of stock-recruitment models for North Sea cod across the three population levels and model formulations (SSB, rSSB, ndSSB, Table 3). At the three populations scale, models including SST as a covariate (models 2, 6 and 9) have In the single population case, ndSSB is not calculated.

TA B L E 2 Models used in analyses with their predictors and climate variables
The interannual variation in SSB, rSSB and ndSSB is largest for the Viking and Northwest units (Figure 4). The three indices show similar interannual patterns (but for example the peak around 1985 in Viking SSB Figure missing in rSSB and ndSSB), but can still result in different fit to data (Figure 5b), with, for example, higher fit using ndSSB than rSSB in Viking and Northwest populations (Table 3). Using traditional SSB with retention anomaly (model 4, Table 3), the effect of retention is captured at the single population scale but not at the two or three populations scale. In this case, the effect has a similar magnitude and effect as the inclusion of SST as covariate (Figure 5a,c).
The interannual variation in retention and connectivity is relatively low (Figure 6). Retention is higher in South (0.39 ± 0.08) and Northwest ( to the outside area is higher than the retention rate (0.39 ± 0.13).
The proportion of larvae remaining within the study area but not settling within a nursery area is high for all populations (0.46 ± 0.08, 0.48 ± 0.09 and 0.36 ± 0.08 for South, Northwest and Viking populations, respectively). Only drift from the South to the Viking population and to the Skagerrak and from South to Northwest populations significantly increase or decrease, respectively, in time ( Table 4).
The NAO and SST indices are significantly correlated with drift anomalies across population scales (Table 4)

| D ISCUSS I ON
In this study, we combined long-term observational data with modelled estimates of larval transport to quantitatively assess the effect of transport on recruitment across spatial scales of observation, and we propose a novel approach for measuring effective biomass contributing to recruitment. While the effect of transport on recruitment has previously been explored using coupled biological-oceanographic models (e.g. Daewel et al., 2015;Hinrichsen et al., 2016), direct inclusion of ELS dispersal model output in stockrecruitment models is less common (but see Hidalgo et al., 2019).
Some studies have used proxies for larval transport such as wind speed (e.g. Hare et al., 2015;Köster et al., 2003) Positive correlations are represented in light grey, negative correlation in dark grey.
White cells indicate non-significant values at the 0.05 level. modelled particles (Baumann et al., 2006). We instead incorporate estimates of the proportion of cod larvae retained within a population and the influx of larvae from neighbouring populations, that is a more direct proxy for the effect of larval transport (Hidalgo et al., 2019), and apply this approach to the North Sea cod.
Our results suggest that although larval drift appears to play a minor role in the recruitment dynamics of North Sea cod, the effect is comparable in magnitude to the well-established effect of SST on cod recruitment (Beaugrand & Kirby, 2010;Nicolas et al., 2014).

| Effective biomass
Estimating effective biomass is a novel approach to account for larval transport compared to using SSB with additive covariates. In traditional stock-recruitment models, a covariate allows higher (lower) asymptotic value, that is higher (lower) expected recruitment at a given SSB value, while maintaining the shape of the curve (Figure 5a, c). Incorporating retention anomaly as a covariate, the interpretation is that a positive anomaly (higher than usual retention) results in higher level of recruitment compared to the same level of SSB with a lower drift anomaly. Subbey et al. (2014)  Spawning stock biomass is a suboptimal variable for predicting recruitment, since it does not capture biological aspects such as age and size structure, sex ratio, total egg production, skipped spawning or interannual variability in fecundity or condition (Köster et al., 2003;Marshall et al., 2003;Marshall, Needle, Thorsen, Kjesbu, & Yaragina, 2006;Marteinsdottir & Begg, 2002;Mintevera et al., 2019). Marshall et al. (2006) and Köster et al. (2003) show that female-only spawner biomass and predicted potential egg production are better predictors of realised egg production than SSB in Northeast Arctic and Baltic cod stocks. Similarly, our study shows that effective SSB might be a better predictor at population scale for some populations, such as the Viking unit, charac-

| Effects of drift and climate variables on recruitment
We investigated the emergent relationships between climate variables (SST and NAO), connectivity and retention metrics and recruitment. High SST can influence recruitment through faster development and thus increased retention and survival to settlement (Heath, Kunzlik, Gallego, Holmes, & Wright, 2008). Additionally, both SST and NAO could be proxies for other phenomena acting at local scales, such as food availability (Capuzzo et al., 2018;Nicolas et al., 2014) and flow regimes (Henriksen et al., 2018). NAO and SST can furthermore be correlated to connectivity and retention (  (Nicolas et al., 2014), physiological constraints (Butzin & Pörtner, 2016;Nunez-Riboni, Taylor, Kempf, Pu, & Mathis, 2019) and predation by warm-water predators (Akimova, Hufnagl, et al., 2016). For the Viking population, our results suggest that SST could influence recruitment through drift both positively and negatively.
In fact, higher SST is associated with increased retention and inflow, but also increased outflow to Skagerrak. At present, the two effects seem to counterbalance each other: SST does not influence recruitment according to our model results. However, with increasing SST this equilibrium, which currently masks the underlying relationships, not rSSB, improves the stock-recruitment model fit. The two models involve the same variable (retention), but differing mechanisms, as described in the section "Effective biomass" above.
Overall, our results indicate that the key mechanisms affecting recruitment (summarised in Figure 7) include: SST in the South population through processes unrelated to larval transport, SST and transport through the same underlying phenomenon in the Northwest population, with inflow from the Viking population and retention being higher in low SST years, and inflow from other populations into the Viking population ( Figure 5b; Table 3).

| Drift patterns, retention and population connectivity
The retention and connectivity patterns estimated here broadly reflect known patterns for the area. The southern North Sea is characterised by a generally retentive system (Henriksen et al., 2018), while in the northern area there is a strong flow to the Skagerrak and the Norwegian Sea (Huserbråten et al., 2018). Consequently, the South and North populations are generally isolated, with limited connectivity (Heath et al., 2008). According to our results, connectivity between the Northwest and South units is higher, but declined from the 1970s to present, while connectivity between South and Viking units increased.
Drifting into a suitable nursery area, however, is not enough for granting survival to recruitment, as density dependence and predation after settlement might influence successful recruitment into the new populations (Akimova, Hufnagl, et al., 2016;Heath et al., 2014).
Some studies discriminate potential connectivity (estimated from modelled particle drift) from effective connectivity using genetic methods (e.g. Bode et al., 2019;Jahnke et al., 2017). In our study, effective connectivity is an emerging result of fitting stock-recruitment models to data after inclusion of drift anomaly. Our results highlight that effective connectivity only affects the Viking and Northwest populations.
Notably, we assess how larval drift influences recruitment, irrespective of whether individuals merge with the host population or return to the natal population after being accounted as recruits.
For example, our results indicate that larvae from the South unit enter the Viking area and survive until being accounted as recruitment of the Viking population (shown by higher fit with ndSSB than rSSB or SSB). However, the Viking and South units show genetic differences, generally considered incompatible with interbreeding between populations (Heath et al., 2014). We therefore speculate that juveniles from the South unit settle in the Viking area and survive until age 1, to then return to the population of origin. This mechanism, known as homing behaviour and site fidelity, is known for cod in the North Sea (Neat et al., 2014) and between the North Sea and Skagerrak Jonsson, Corell, André, Svedäng, & Moksnes, 2016), and is suggested for larvae drifting from the Norwegian Trench (within our Viking area) to the Norwegian Sea (Huserbråten et al., 2018).
Although drift between the North Sea and Skagerrak is well known , we show here for the first time, to our knowledge, that larval drift from the Viking and South units into the Skagerrak is potentially of the same order of magnitude, showing an increasing trend in time and positive correlation with SST (Table 4).
Although the effective contribution cannot be determined in this study, trends in larval influx from the North Sea might have implications for management and recovery of cod in offshore and coastal areas of the Skagerrak.
Our results are influenced by the assumptions and simplifications of the ELS dispersal model. However, in a previous study, Kvile et al. (2018) showed that the present model configuration yields comparable results to a more realistic but computationally costly alternative. Specifically, both the inclusion of vertical swimming behaviour and the use of a higher resolution ocean model that resolves tidal circulation had limited effects on larval drift patterns compared to interannual variations in ocean dynamics. Since our aim here was to quantify long-term interannual variation in population connectivity, we opted for a less computationally costly representation of vertical movement using fixed drift depths, and applied the coarser ocean model that was available for 44 years. Additionally, sensitivity analyses of the parameterisation of ELS mortality and settlement size, the latter related to temperature-dependent growth, confirmed the robustness of the results to these key parameters (Tables A2   and A3 in the Appendix S1). Finally, factors such as spatially explicit predation pressure and prey fields, variability in fecundity, juvenile mortality through predation and density dependence upon settlement might all affect recruitment dynamics, but are not accounted for in this study. These caveats need to be considered in the interpretation of results.

| Implications for management
Despite the relatively low prediction power and major assumptions (Subbey et al., 2014), stock-recruitment models are routinely applied in management for short-term advice (e.g. Punt, 2019), and there is increasing interest in including spatial structure in recruitment dynamics in stock assessment (Cadrin et al., 2019;Hidalgo et al., 2019;Punt, 2019). Although reliable ocean current forecasts are not available in advance, estimates of larval drift can be useful to inform short-term forecasts (Henriksen et al., 2018;Hidalgo et al., 2019). This effort is however constrained by the availability and rapid applicability of ocean models in the context of operational fisheries oceanography (Hidalgo et al., 2019).
We find relatively low fit to data in the stock-recruitment models for North Sea cod, and inclusion of indices for larval drift results in relatively small improvements. Considering the computational cost of running ELS dispersal models, one must therefore carefully consider the benefits of this approach for the specific case at hand. Regardless, our study highlights a novel approach for accounting for connectivity in stock-recruitment dynamics, with potential applications for fisheries assessment and management in stocks characterised by highly dynamic oceanographic conditions.
Adopting spawning output metrics that account for effective con- Future research should focus on how climate change can influence larval transport, survival of larvae drifting between units and homing behaviour. Understanding these aspects, and developing operational fisheries oceanography and its application to management, will improve our capacity to tailor management to the population structure in the context of a changing climate. search program in collaboration between The Norwegian Academy of Science and Letters, and Equinor. We thank an anonymous referee for valuable comments that substantially improved the article.

CO N FLI C T O F I NTE R E S T S
The authors declare to have no conflict of interest.

AUTH O R CO NTR I B UTI O N S
GR, KK, ØL and AME conceived the research idea; GR and KK ran data analyses and wrote the manuscript; KFD, TK, ØL, AME and NCS participated in discussions of the results and critically reviewed the manuscript.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are available from the corresponding author upon reasonable request. Spatially-resolved influence of temperature and salinity on stock and recruitment variability of commercially important fishes in the