The role of spatial distribution for growth and survival of juvenile cod Gadus morhua in the Barents Sea

C. The role of spatial distribution for growth and survival of juvenile cod morhua in the Barents accepted  July . While the importance of early life survival and growth variations for population dynamics is well documented, there is still a relatively limited understanding of how survival and growth is aﬀected by the species’ spatial distribution. Using Barents Sea spatial bottom survey data (– ), we study the spatiotemporal variability of the juvenile Atlantic cod ( Gadus morhua ) growth and survival. We used indices of the spatial distribution of juvenile cod at age- to study the role of distribution for the change in abundance and mean body size through their second winter of life (from age- to age-). Over the  years analysed, we found that the location where the age- cod are in the Barents Sea matters for their growth and survival. We found that year-classes growing up in the western Barents Sea have higher mortality but faster growth than year-classes distributed farther east. Our results indicate that the biotic and abiotic conditions encountered at the settlement location may inﬂuence the spatial survival and growth of age- cod and subsequently the population dynamics. Our results underscore the importance of distribution for survival and growth early in life and by providing this essential information has implications for stock


Introduction
Since individuals from the same year-class of a population are likely to experience different environmental conditions depending on their location, the spatial distribution of the year-class is important for its growth and survival (Kareiva et al., 1990;Ciannelli et al., 2007;Ciannelli et al., 2008). For many fish species, this is particularly true for young age classes with relatively low motility Ferreira et al., 2020). The distribution of fish juveniles may hence influence the population recruitment, i.e. the abundance at the age when the fish enter into the fisheries. To understand and eventually predict recruitment in these populations, it is necessary to investigate which factors influence distribution and in turn how distribution influences early life growth and survival. For fish species with pelagic early life stages, research has focused on the roles of spawning location and drift routes for growth and survival, e.g. in relation to temperature gradients (Opdal et al., 2008) and for closing the life cycle by reaching favourable nursery grounds (Agostini and Bakun, 2002;Sinclair and Power, 2015). Relatively less attention has been given to the role of the spatial distribution of the juveniles after settlement to a more demersal lifestyle.
The Northeast Arctic (NEA) stock of cod Gadus morhua is of high economic and ecologic importance in the Barents Sea. These cod mature at 7-11 years of age (Ottersen et al., 2006) and migrate to the coast of Norway to spawn in February to early May Langangen et al., 2019). The eggs, larvae, and early juveniles then drift to the Barents Sea, carried by the Norwegian Coastal Current and the North Atlantic Current (Tereshchenko, C International Council for the Exploration of the Sea 2021. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited.  J. M. Durant et al. 1980;Ellertsen et al., 1981;Mukhina et al., 2003;Hidalgo et al., 2012). In September-October, the juvenile cod (0-group) are settling . Their spatial extension is then at its maximum (Sundby, 1989;Ciannelli et al., 2007), going from the coasts of Spitsbergen and northern Norway to the central and eastern Barents Sea. The young 1 year of age cod (age-1) seems then to remain in the area where they settled during autumn as 0-group, i.e. at the end of their pelagic drift phase (Ottersen et al., 1998).
The distribution of the 0-group is depending on environmental conditions such as current strength and temperature  and possibly on numbers of offspring spawned that year. Water temperature in the Barents Sea is closely related to the amount of inflow of relatively warm Atlantic water (Ådlandsvik and Loeng, 1991). Enhanced inflow of warm water from the North Atlantic provides favourable conditions for young cod growth, and at the same time tends to transport the cod juveniles far eastwards in the Barents Sea (Ottersen and Loeng, 2000;Helle et al., 2002;Ciannelli et al., 2007). Hence, strong year classes tend to have a more easterly distribution in the Barents Sea than weak year classes (Shevelev et al., 1987;Ottersen et al., 1998;Ottersen et al., 2002). An expansion of abundant year classes into environmentally and geographically marginal areas in the eastern and northern Barents Sea may also be a result of density-dependent habitat selection, although the limited swimming capacity of the juveniles the first 1-2 years of life argues against this possibility . The distribution possibly also depends on the demographic structure of the spawning stock, as old and large females may spawn over a longer period of time, a larger area and have more viable offspring . Consistent with this hypothesis, NEA cod eggs and larvae have been found to be more widely distributed in years when the mean age and size of spawners is high (Stige et al., 2017;Endo et al., 2021).
Depending on where the cod juveniles are in the Barents Sea, they are confronted with different environmental conditions, e.g. in terms of temperature, that are influencing their growth and survival (Sundby, 2000;Ciannelli et al., 2007). For example, juveniles that are distributed far east in the Barents Sea tend to experience low temperatures and slow growth after settlement Ottersen et al., 2002;Ciannelli et al., 2007). An eastern distribution of juveniles in warm years with high larval and juvenile survival may explain a negative correlation between abundance and mean length of juvenile cod from age-1 onwards, although density-dependent effects could not be ruled out . Analyses of more recent data have confirmed a negative abundance-size association for NEA cod juveniles, but have also shown similar associations for three other fish stocks around the time of transition from pelagic to demersal habitats . Increased competition for space is likely to occur in such transition periods (Juanes, 2007), suggesting that density-dependent growth could be a common explanation for the associations. When it comes to spatial patterns in survival, cannibalism is probably an important driver. NEA cod juveniles in the eastern part of the distribution seem to have higher survival through the first winter (from age-0 to age-1) than juveniles farther west, most likely because of lower cannibalism rate Stige et al., 2019).
We here build on the earlier studies on NEA cod and extend these by including spatial bottom survey data from the last two decades. Our focus is on the transition between the two first surveyed age classes after settlement (age-1 and -2). Our aims were to assess the association between spatial distribution of age-1 cod juveniles and survival and growth to age-2 and to investigate whether age-1 distribution correlated with potentially influencing variables, i.e. cod population density, demographic structure, temperature and inflow. To do so we quantified how age-1 cod distribution metrics of location and dispersion, together with abundance and mean length at age-1, and temperature explained abundance and mean length at age-2. Abundance at age-2 is closely correlated with recruitment at age-3, whereas mean length at age-2 is likely to impact reproductive capacity later in life . We addressed the following hypotheses: High mean length at age-2 is associated with a south-western distribution of age-1 juveniles (H1a), high temperature in the Barents Sea (H1b), and low abundance at age-1 (H1c). High abundance at age-2 is associated with an eastern distribution of age-1 juveniles (H2a), and high temperature in the Barents Sea (H2b). An eastern distribution of age-1 juveniles is associated with high temperature in the Barents Sea (H3a), and strong inflow (H3b). A wide distribution of age-1 juveniles is associated with high population density (H4a), and high mean age and size in the spawning stock (H4b).

Data
Data of spatially explicit NEA cod Gadus morhua abundance at age-1 between 1994 and 2018 from the Barents Sea joint Norwegian/Russian winter survey (Mehl et al., 2019) were extracted from NMDC.no (access date 12.06.2020). Age-1 cod here refers to cod of approximately 10 months of age that were sampled in the winter survey between January and March. This winter survey is a collaboration between the Norwegian Institute of Marine Research (IMR) and the Polar branch of the Russian Federal Institute of Fisheries and Oceanography (PB VNIRO formerly PINRO) and provides both bottom trawl and acoustic estimates by age group for several demersal fish, among others the NEA cod. The survey area is divided into 23 strata, with since 2014 three more northern strata (Mehl et al., 2019). In the following analysis, we did not use these three additional strata since they were not sampled during the whole period. However, the inclusion of the new strata in 2014 was linked to a likely northward expansion of the distribution, which is not fully captured in our analyses. Here, we used abundance indices (numbers in millions) per strata per year from the bottom trawl survey estimated by StoX software (Johnsen et al., 2019).
We calculated the location (longitude, latitude) of the geographical centre of gravity for each of the strata i using the functions poly_center in the pracma library in R 3.6.1 (R Core Team, 2019) and correcting for earth flattening at the poles and difference of distance with latitude.
We calculated the weighted mean longitude and latitude centre of gravity (Lon and Lat), with the abundance of each stratum as weight. Using the weighted variance of the centre of gravity in longitude and latitude (var.lon and var.lat) using the SDMTools library it was possible to calculate an index of dispersion (Dis) each year by the formula: Average length in cm of age-1 cod (Length age1 ) and age-2 cod (Length age2 ) in the winter survey were taken from the 2019 ICES Arctic Fisheries Working Group (AFWG) report (Table A5 in ICES,  2019).
Abundance indices in millions for age-1 cod (Abundance age1 ) and age-2 cod (Abundance age2 ) in the winter survey were taken from the 2019 ICES AFWG report (Table A3 in ICES, 2019).
The Fugløya-Bear Island section (70 • 30 N and 20 • 00 E to 74 • 15 N and 19 • 10 E; see Figure 1) is in the western entrance of the Barents Sea. Moorings are deployed in this section to monitor the northward flow of warm Atlantic Water that reaches the Arctic Ocean through the Barents Sea. The volume flux in Sverdrups at Fugløya-Bear Island (Flux) was uploaded from NMDC.no (16.10.2020) (Ingvaldsen, 2019). Note that the data cover 20 years only (1997 to 2016).
The Kola meridian transect (section along the 33 • 30'E between 69 • 30'N and 77 • 00'N; see Figure 1), intersecting the Murman Current in the south-central Barents Sea, covers the inflow of Atlantic and Coastal water masses from the Norwegian Sea to the south-eastern Barents Sea. We used vertically averaged sea temperature data collected at five stations on the Kola section at the standard depth: 1, 10, 20, 30, 50, 75, 100, 150, and 200 m) (Bochkov, 1982;Tereshchenko, 1996). This sea temperature index is representative of the Atlantic water masses in the south-central Barents Sea (Ingvaldsen et al., 2003) and correlates positively with strong year-classes of NEA cod (reviewed by Ottersen et al., 2014).
Spawning stock biomass in thousand tonnes (SSB) was taken come from the 2019 ICES AFWG report (Table 3.18 in ICES, 2019). Mean age of the spawning stock biomass in years (mAge, Ottersen et al., 2006) was estimated using information from the 2019 ICES AFWG report (Tables 3.9, 3.11 and 3.16 in ICES, 2019).
Data used in the analysis can found in Supplementary Table S1.

Statistical modelling
To address hypotheses H1 and H2, we related the temporal variability in length of the 2 years old cod (age-2: Length age2 ) and of logtransformed abundance (ln(Abundance age2 )) to different explanatory variables using generalized additive models (GAM), via the mgcv library (Wood, 2006). All models started with a similar formulation: where dept variable is either Length age2 for investigating H1 or ln(Abundance age2 ) for investigating H2. ln(Abundance age1 ) is the abundance indices for age-1 cod log transformed. Length is the yearly average length of specific age class. ST is the yearly average temperature for the Kola section. Lon is the yearly average longitude of the centre of gravity of abundance of age-1 cod and Lat is the yearly average latitude of the centre of gravity of abundance of age-1 cod. f 1 , f 2 , and f 3 are smoothing spline functions and f 4 is a two-dimensional tensor product smooth. All functions f 1 , f 2 , f 3 and the two basis functions of f 4 are thin plate regression splines with extra shrinkage. ε is a normally distributed error term.
The mgcv GAM procedure automatically selects the degree of smoothing based on the generalized cross-validation (GCV) score. GCV is a proxy for the model's predictive performance. However, to avoid spurious and ecologically implausible relationships, we constrained the model to contain at maximum quadratic relationships, i.e. we set the maximum degrees of freedom to 2 for each smooth term (i.e. k = 3 in the GAM formulation).
We entered all candidate explanatory variables in the GAM model and conducted model selection using Akaike information criterion corrected for low sample size (AICc, Hurvich and Tsai, 1989). We retained the model using the set of explanatory variables providing the highest decrease in AICc (lowest AICc, see Table 1). The models were visualized (see Figure 2) using the vis.gam function with too.far = 0.2 (excluding area too far from the predictors).
To verify that the models selected by AICc were the optimal models and not overfitted, we used a cross-validation technique. We calculated the model root-mean-square error (RMSE) by training our model on a subset of 80% of the data, picked up randomly without replication, and testing it on the 20% remaining. We repeated this routine 500 times and average the RMSE (Table 1). A decreased RMSE with increased model complexity indicates a model that is not overfitted.
The quantification of an individual explanatory variable's contribution to a multiple regression model explaining the variation of the stock biomass was studied using the package relaimpo (Grömping, 2006) using the "proportional marginal variance decomposition" (pmvd).
To avoid statistical issues with collinearity between explanatory variables, we calculated variance inflation factors (VIFs) for all combinations. Since the VIFs <3, we do not consider collinearity to be a problem in this analysis. This conclusion was confirmed by an exploration of the selected models' concurvity using the function concurvity from the package mgcv. We did not detect any significant autocorrelation in the residuals (using autocorrelation function acf) in any selected model.
To explore what may affect the distribution of age-1 cod in the Barents Sea (hypotheses H3 and H4), we conducted a Pearson's product-moment correlation analysis on Lon, Lat and Dis against The deviance (Dev) of the best model selected by AICc is indicated in bold. The lowest RMSE is indicated in bold. RMSE = test root-meansquare error calculated by cross-validation (see methods). Note that the selected model by AICc displayed a low RMSE in comparison to the other length models. This is not the case for the abundance model, where a simpler model [using only f 1 (ln(Abundance age1 )] as an explanatory variable) gave a lower mean RMSE than the AICc selected model (note that a comparison of medians gave an opposite result, Supplementary Figure S). However, compared with a simple model (using only f 4 (Lon, Lat)), the selected model had a lower RMSE.
ST and ln(Abundance age1 ) the same year, and ln(SSB), mAge, and Flux the year before.

Results
The best models selected by AICc for the two dependant variables (length at age-2 for hypothesis H1, abundance at age-2 for H2) both keep the average distribution of age-1 cod as predictor term, and this term explains a substantial part of the deviance observed ( Figure 2 and Table 2). Model 1 (Length at age-2) explains 68.5% of the deviance observed (Table 2). Length age-1 explains 31.9% of the deviance, ST explains 12.2% while the remaining 24.4% are explained by the average distribution of age-1 cod (see Table 2 to see these values relative to the total deviance explained). We found that on average, the longer the cod are at age-1 the year before, the longer they are at age-2 (p < 0.01) with a slope of 1.09 ± 0.28 (SE). Temperature shows a positive and significant linear relationship with length (p < 0.05) with a slope of 1.51 ± 0.52 (SE). We found a significant relationship (p < 0.05) with the average location of age-1 cod. The more southwest the cod are at age-1, the higher is the length at age-2 ( Figure 2a). The difference in expected mean length between the most south-western and the most north-eastern distribution observed is about 3 cm (20.7 cm compared to 17.5 cm, Figure 2a).
Model 2 (Abundance at age-2) explains 86.8% of the deviance observed. The two terms, ln(Abundance age1 ) and average distribution of age-1 cod, explain respectively 62.4 and 24.4% of the deviance. We found that on average, the more abundant the cod are at age-1, the more abundant they are at age-2 (Figure 2b). This significant relationship (p < 0.001) is linear with a slope of 0.60 ± 0.07 (SE, Supplementary Figure S3). We found a significant relationship (p < 0.01) with the average location of age-1 cod (Figure 2b). The more southeast the cod are at age-1, the higher is the abundance at age-2. Abundance is predicted to be more than twice as high (ca. 1 on ln-scale) with a south-eastern than a north-western distribution (Figure 2b). The model also suggests a decline in abundance if the distribution is too far southeast, but this effect is caused by a single year with an extreme south-eastern distribution (2003). An alternative model (delta AICc < 2, Table 1), adding to the twodimensional location term the smooth effect of the dispersion in-dex, shows that the more age-1 cod are dispersed in the Barents Sea (higher value of Dis), the more abundant they are at age-2. However, this relationship is non-significant (p = 0.30) with a very flat slope of ca 0.19 ± 0.11. Note that the relationship becomes significant if the smooth effect of the dispersion index replaces the two-dimensional location term (Table 1).
Knowing that the winter survey only started to sample the northern Barents Sea in 2014, the centre of gravity location used in the analysis may be biased (too much south) as cod may have been present in this northern area earlier than 2014 ( Supplementary Figure S1). To check this potential bias, we have used the same model formulations (i.e. model 1 and model 2) on data taking into account the 3 omitted strata (see method) and found a similar direction of effect for the different explanatory variables.
The correlation analysis addressing hypotheses H3 and H4, summarized in Table 3, indicates significant positive relationships between the latitude (Lat) and the log-transformed spawning stock biomass (SSB) the previous year (Table 3). The location of the distribution in longitude, latitude or dispersion seem otherwise not to be affected by the variables tested, i.e. mAge, Flux, ST, ln(Abundance age1 ), and ln(SSB).

Effect of spatial distribution and other factors on length (H)
By accounting for cod length at age-1 in model 1, the effects of the other covariates can be interpreted in terms of growth rate and/or size-dependent survival (e.g. Stige et al. 2019). Our analysis indicates that the location where the age-1 NEA cod are in the Barents Sea matters for their growth and survival from age-1 to age-2 ( Figure 3). Similarly, the locations of the pelagic larvae of NEA cod have large effects on their growth potential earlier in life (Opdal et al., 2008), with a northeastern distribution associated with low ambient temperatures and slow growth. As proposed by hypothesis H1a, our results show that the more northeast in the Barents Sea the age-1 cod are settling the smaller they are at age-2 (Table 2 , Figure  2a). A plausible interpretation of this finding is that a north-eastern distribution of the cod juveniles is associated with lower ambient temperatures, which combined with physiological temperature lim-Figure 2. I Effect of age- cod distribution in the Barents Sea on two population dynamic variables: length and abundance at age-. The figure shows the predicted age- cod mean length (Panel a) and mean abundance (Panel b) depending on the location of the centre of gravity at age- and mean values of other predictor variables ( Table  ). The blue dots are the data (for panel b the size is proportional to the dispersion variable value). The colour scheme indicates the strength of the effect (from red to yellow: lower to stronger).
itations on growth may lead to lower mean size. In addition, spatial differences in prey availability may play a role. Our model also supports that an increase of sea temperature leads to a faster growth as suggested by H1b (Supplementary Figure S2). The analysis hence captures how cod growth is influenced by interannual differences in Barents Sea temperature as well as by temperature differences caused by shifts in cod distribution. However, we did not detect the effect of the abundance at age-1 on the length at age-2 that was suggested by H1c, thus providing no support for density-dependent growth.
Length at age is also potentially influenced by size-dependent survival linked to temperature. A theoretical study has shown that increasing temperature may increase predation risk on cod larvae in poor food conditions (Fouzai et al., 2015). Since mortality rates of juvenile fish are generally lower for large individuals than for small individuals (Sogard, 1997), high predation rates on small individuals in warm years can reinforce the positive effect of temperature on growth. However, juvenile survival was not correlated with temperature in our study, suggesting that the temperature-size association was mediated by temperature effects on growth rather than size-dependent survival. Further, size-dependent survival seemed not to play an important role for the association between the location of the centre of gravity and mean length. If this association was to be driven by size-dependent survival, the higher mortality suggested by model 1 in the southwestern Barents Sea would have to disproportionally affect large individuals, which is contrary to expectations from previous studies (Sogard, 1997).
In our data, most of the distribution centres of gravity are located in the mixed-water mass (0 • C < bottom temperature < 2 • C (Frainer et al., 2017)) that is moving northward with climate warming (Fossheim et al., 2015). This northward movement of water masses is associated to a fish community change (Fossheim et al., 2015) that may affect the food availability for young cod. The first 2 years of life is the period of diet diversification for the young cod. The diet shifts from euphausiids and other invertebrate food items such as copepods to a piscivorous diet (Holt et al., 2019). The diet of two years and older cod (>20 cm length) is then dominated by fish (> 60%, Holt et al., 2019). Knowing that temperature is affecting distribution of species and the community composition in the Barents Sea (Fossheim et al., 2015), it can also affect the growth of cod by changing the availability of their favourite prey during growth. For instance, the distribution of euphausiids changes with temperature in the Barents Sea (Zhukova et al., 2009). It remains to be known how fast the cod population, or different age groups of it, move northward compared to their prey (for instance euphausiids) and temperature conditions. Indeed, the distribution of juvenile cod moves along with the other fish species (cod was a key species in the boreal fish group in the analysis by Fossheim et al., 2015), but it was shown an age specific-distribution response to temperature with larger northerly distribution shifts for the older age classes. However, note that even in a warming situation, when area boundaries extend further north and east, most of the centres of gravity are located in the central Barents Sea (Figure 1).

Effect of spatial distribution and other factors on abundance (H)
By accounting for age-1 abundance in model 2, the effects of the other covariates can be interpreted in terms of mortality rate (see, e.g. Stige et al. 2019). Our results indicate that the more east the settlement occurs, the higher is the survival (as shown by the abundance at age-2) which conforms to H2a (Figure 3). Our results on survival to age-2 are similar to what was shown for the survival of age-0 cod to age-1  but based on data from 1980-2004 (i.e. with only 10 years of overlap with the current study). Ciannelli et al. (2007) explained their result on age-0 cod survival with the overlap with predators, e.g. adult cod, which are more abundant in the western part of the Barents Sea. Cannibalism of juveniles by older fish has been shown to be common in NEA cod (Yaragina et al., 2009;Holt et al., 2019) and to affect the population dynamics (Ohlberger et al., 2014) by affecting juvenile survival in some years. Cannibalism could then have a role in explaining the change in age-2 abundance (Dingsør et al., 2007). One can suggest that juveniles in the area with high abundance of age 3-6 cod (the most cannibalistic age-groups, Holt et al., 2019) would show less survival to age-2. Unfortunately, spatial cannibalism data (i.e. mortality of ages 0-2 cod induced by cannibalism) are very uncertain making such hypothesis difficult to test. Distribution changes can also impact survival through other mechanisms than through cannibalism. Opdal et al. (Opdal et al., 2008;Opdal et al., 2011) suggest that the recruitment success of NEA cod is more linked to shifts in spawning grounds than to interannual environmental variability, because larvae from southern spawning grounds experience higher temperature, faster growth and higher survival through the early life stages than larvae from northern spawning grounds. However, when taking into account spatial patterns in natural mortality rates of the larvae, offspring from central and northern spawning grounds appear to have similar survival probability (Langangen et al., 2014;Langangen et al., 2016). Current patterns during the pelagic drift may further impact both distribution and survival through the first year of life, as indicated by the increased survival of NEA cod juveniles to age-1 in years with an eastern distribution of the juveniles and of Atlantic water masses in the Barents Sea . This association, could be caused by a better suitability of prey associated with warmer water, for instance by increased availability of prey in warm water masses, such as Atlantic zooplankton species advected from productive areas farther south (Kvile et al., 2016). Our findings add to these by showing that the distribution is likely to impact survival also later in juvenile life.
It is to note that the interannual changes in distribution, i.e. in the location of the centre of gravity, are interlinked with the changes in the dispersion around the centre of gravity ( Supplementary Figure S1 and Table S1). When the distribution is wide the centre of gravity is located more east in the Barents Sea (p < 0.01, Supplementary Figure S1); in the coldest part of the Barents Sea. That is, when the distribution is wide, a higher proportion of individuals occupy colder areas to the east and north while when the distribution is concentrated, a higher proportion of the fish are located in the relatively warm southwestern Barents Sea (Supplementary Figure S5).
It seems that there is an opposite trend from west to east in the Barents Sea that describes a high growth-low survival to low growth-high survival continuum. Year-classes that grow up concentrated in the western Barents Sea have relatively high mortality but fast growth to age-2, whereas year-classes that are distributed in the south-eastern Barents Sea have lower mortality but also lower growth (see relationship between distribution and location in Supplementary Figure S1). Note that an alternative model to model 2 indicates a positive effect of a wider dispersion on survival (Table  1 and Figure 2b). While large mean body size prior to this age seems to have large positive consequences for mean survival of yearclasses , survival appears to be independent of mean length at age-1 and temperature as suggested by H2b (Table  1). Also survival from age-2 to recruitment at age-3 appears to be independent of mean length at age-2 . The reduced body size of year-classes that are distributed over large areas of the Barents Sea is therefore not expected to impair recruitment.

Correlations between environmental factors and population variables with distribution (H-H)
Several factors can cause a change in the distribution of age-1 cod in the Barents Sea. We suggested that the geographical location could be affected by physical factors (H3) and the dispersion by population variables (H4). None of the proposed hypotheses was confirmed. Instead, we found that, correlating the north-south distribution (latitude of the centre of gravity) to the population abundance, the bigger the spawning stock biomass is, the more north is the distribution of the age-1 cod (Table 3). This was complemented by a non-significant trend of northern distribution of age-1 cod with higher sea temperature. We found no other significant correlation explaining the distribution of age-1 cod, thus not supporting the possible effect of sea temperature (H3a), the strength of the water influx (H3b) or the population structure (H4b). Regarding the water influx, we note that the data available did not cover the whole studied period, with the missing years corresponding to the lowest temperatures recorded and relatively southern distributions (Supplementary Table S1). Using the three omitted strata in the correlation analysis (see Material and methods), we found that the dispersion index becomes significantly positively correlated with sea temperature, spawning stock biomass and mean age of the spawning stock (respectively H3a, H4a, and H4b), showing (i) that the Table 2. Models linking population dynamic variables to the change of age- distribution.
years 2014-2018 have a strong influence on our models and (ii) that the information of these three strata is important. Indeed, a geographically wide offspring distribution in years with high mean age and weight in the spawning stock has been shown for NEA cod eggs (Stige et al., 2017) and larvae (Endo et al., 2021). Furthermore, a wide distribution of cod larvae was associated with a high liver condition index of the spawners and high spawning stock biomass (Endo et al., 2021). We note that temperatures in the studied period mostly ranged from near-average to warm (Supplementary Table S1), and that low contrast in climate conditions may have contributed to the lack of significant effects. A factor that can affect the distribution of age-1 cod is the spawning location along the Norwegian coast. Opdal et al. showed that a northern spawning leads to a wider and more eastern drift in the Barents Sea (Opdal et al., 2011). According to our results, the juvenile distribution resulting from a northern spawning is associated with low growth but high survival of the juveniles. Spawning location may hence have contrasting effects on offspring survival at different life stages, with egg and larval survival being highest for southern spawning locations (Opdal et al., 2008;Opdal et al., 2011;but see Langangen et al., 2014;Langangen et al., 2016) and juvenile survival being highest for northern spawning locations (this study). Langangen et al. (2019) showed that between 2000 and 2016, the mean spawning latitude was significantly positively correlated to the Kola temperature and that the bigger fish were spawning in the northernmost spawning grounds. Our results underline that the consequences of such distribution shifts for recruitment depend on an interplay between processes at different pre-recruit life-stages.

Implication of our findings
While the central challenge in fisheries ecology still generally remains focussed on temporal variability of fish populations, more of this focus is shifted toward a better understanding of how species distribute and survive over space (Ciannelli et al., 2008). The current study adds to the later. Density-dependence in growth, reduction of growth rate when abundance is high because of intra-specific competition, has been widely reported as a key process in the regulation of many fish populations (Lorenzen and Enberg, 2002). For example, Stige et al.  showed that high abundance of NEA cod at age-0 and age-1 are associated with poor growth to age-1 and age-2, respectively. These authors found similar associations for three other fish populations. In our study, we did not find such relationship, the abundance of age-1 NEA cod was not explaining the length at age-2 (Table 2). Andersen et al. (2017) suggested that density-dependence can occur at different life stages depending on spatial distribution and habitat size. Marine fish species typically have complex life cycles composed of multiple life stages with different degrees of habitat preferences. However, marine fish are more constrained in space during the early phase of their life cycle than otherwise. When such constraint occurs the question becomes whether these habitats saturate (MacCall, 1990) potentially triggering a compensatory density dependence mechanism. In our study we found that one of the proximate factors affecting the changes of length at age-2 were linked to their distribution at age-1; the more the age-1 cod are concentrated in the southwestern Barents Sea, the bigger they are becoming at age-2. Our findings suggest that spatial aspects may play a role in how the density dependence operates. These results are in line with MacCall's basin model of habitat suitability, which considers that individuals first fill habitats of high quality (the ideal-free distribution of Fretwell and Lucas (Fretwell, 1969)) of a heterogeneous "basin" environment until some carrying level is reached (MacCall, 1990). In such a case, at high abundance, more individuals should occupy marginal habitats where growth conditions are poor (for instance cold eastern and northern parts of the Barents Sea). However, we found no such relationship between age-1 abundance and east-west distribution or crowding. Instead, we found a positive relationship of the spawning stock biomass with the north-south distribution of age-1 in the Barents Sea (Table 3). This finding suggests that the age-1 distribution may be a result of density-dependent habitat selection by spawners rather than by offspring, which is consistent with juvenile cod having limited capacity for self-propelled movement over long distances. A northern offspring distribution at high spawning stock biomass could be linked to the concurrent cod stock increase (Kjesbu et al., 2014) with a strong northward expansion in the recent years (Johansen et al., 2013). Further work is needed to investigate the causal mechanisms behind the distribution changes of the adults, and whether densitydependent spatial patterns in growth and survival of the offspring play a role. In summary, we hypothesize that at high spawning stock biomass, cod tend to use the northernmost spawning grounds more than at low spawning stock biomass, causing an expansion of the offspring into the colder, northern parts of the Barents Sea, resulting in low growth rates and low mean size of juveniles. At the same time, the distribution changes are associated with changes in survival, with the net effect depending on spatial patterns in growth and survival through the entire early-life period until the juveniles are large enough to determine their own distribution. It is to note that the recent increase of spawning stock biomass is coinciding with the period of warming in the Barents Sea (Kjesbu et al., 2014). Our results complement previous studies investigating the importance of distribution for pre-recruitment survival (e.g. for cod eggs, (Stige et al., 2017), larvae, (Endo et al., 2021), and age-0 juveniles, ). From age-2 to age-3 (recruitment age), the survival of NEA cod juveniles appears to be relatively insensitive to body size . Hence, poor growth from age-1 to age-2 probably has limited effects on recruitment, and mainly influence population dynamics through effects on reproductive potential at a later age.
The management of fisheries requires understanding the regulation mechanisms of fish abundance. Our results indicate that the spatial distribution of spawners may have complex impacts on offspring growth and survival, implying that the spatial distribution of fisheries may have unforeseen consequences for stock productivity ( Figure 3). Our results also suggest that stock assessments can become more accurate if spatial heterogeneity in natural mortality values is accounted for. For instance, the mortality of an age group could be calculated as a weighted mean of regional mortality estimates weighted by the spatial distribution.

Supplementary Data
Supplementary material is available at the ICESJMS online version of the manuscript.

Data availability
The data underlying this article are available in the article and in its online supplementary material.

Funding
This work was funded by the Norwegian Research Council through the project Spaceshift (no. 280468).