Determinism of Temporal Variability in Size at Maturation of Sardine Sardina pilchardus in the Bay of Biscay

Age and size at maturation appear as key parameters governing the dynamics of a population as they affect growth rate, fecundity, and survival. The expression of such life history traits is determined by genetic make-up and modulated by environmental factors mainly through phenotypic plasticity. Moreover, fishing, besides decreasing population size and changing demographic composition can alter allelic frequencies through fisheries-induced evolution by selecting for some particular traits. In the Bay of Biscay, a decreasing trend in both sardine body condition and size-at-age has recently been pointed out at the population level. The Probabilistic Maturation Reaction Norm (PMRN) approach was applied to help disentangle phenotypic plasticity and genetic changes. Based on the analysis of sardine spawning seasonality, PMRN was estimated by considering body condition as additional life-history state variable to predict the onset of maturation. The resulting PMRN was then used to investigate temporal trends in reaction norm midpoints to test whether changes in length at maturation can be explained by plastic and/or evolutionary adaptive change. Overall, our results emphasize for the first time that including sardine body condition as explanatory variable improves predictions of maturation probability. We found that better individual condition increases maturation probability. The assessment of temporal changes in length at maturation confirms the low plasticity in this trait for a species maturing mostly at age-1 and advocates for the use of a monthly time scale when investigating PMRNs for this species. Beside environmental variables included in this analysis (water temperature, chlorophyll-a, and population biomass) that only show a weak correlation with PMRN midpoints, our results reveal no evidence for recent fisheries-induced evolution in the sardine stock of the Bay of Biscay. They suggest that the short-term variability in length at maturation is strongly dependent upon individual growth which is likely driven by environmental factors. For sardine fisheries management, our study highlights the need to consider both the length-composition data and the seasonality within a stock assessment model. Finally, we discuss the fact that considering individual growth trajectories should improve our understanding of the relationship between environmental variability and changes in maturation for sardine.

Age and size at maturation appear as key parameters governing the dynamics of a population as they affect growth rate, fecundity, and survival. The expression of such life history traits is determined by genetic make-up and modulated by environmental factors mainly through phenotypic plasticity. Moreover, fishing, besides decreasing population size and changing demographic composition can alter allelic frequencies through fisheries-induced evolution by selecting for some particular traits. In the Bay of Biscay, a decreasing trend in both sardine body condition and size-at-age has recently been pointed out at the population level. The Probabilistic Maturation Reaction Norm (PMRN) approach was applied to help disentangle phenotypic plasticity and genetic changes. Based on the analysis of sardine spawning seasonality, PMRN was estimated by considering body condition as additional life-history state variable to predict the onset of maturation. The resulting PMRN was then used to investigate temporal trends in reaction norm midpoints to test whether changes in length at maturation can be explained by plastic and/or evolutionary adaptive change. Overall, our results emphasize for the first time that including sardine body condition as explanatory variable improves predictions of maturation probability. We found that better individual condition increases maturation probability. The assessment of temporal changes in length at maturation confirms the low plasticity in this trait for a species maturing mostly at age-1 and advocates for the use of a monthly time scale when investigating PMRNs for this species. Beside environmental variables included in this analysis (water temperature, chlorophyll-a, and population biomass) that only show a weak correlation with PMRN midpoints, our results reveal no evidence for recent fisheries-induced evolution in the sardine stock of the Bay of Biscay. They suggest that the short-term variability in length at maturation is strongly dependent upon individual growth which is likely driven by environmental factors. For sardine fisheries management, our study highlights the need to consider

INTRODUCTION
Understanding mechanisms that regulate a stock's reproductive potential is of fundamental importance for marine fisheries management as it corresponds to the ability of a fish population to produce a viable offspring, which in turn conditions the future recruitment (Trippel, 1999). Reproductive potential is itself influenced by stock productivity, which has declined for a wide range of commercially exploited fish species over the past decades (Christensen et al., 2003;Myers and Worm, 2003). Concomitant to those declines, several studies emphasized drastic changes in life history traits and, in particular, a decrease in both age and size at which fish start to reproduce [e.g., reviewed by Rijnsdorp (1993), Trippel (1995), and Jorgensen et al. (2007)]. Such adaptive changes may result from evolutionary and/or plastic responses, the relative importance of which are difficult to distinguish in exploited fish populations (Rijnsdorp, 1993;Law, 2000). However, since lifetime reproductive success and therefore stock productivity are closely linked to both age and size at maturation (Bernardo, 1993), it becomes essential to identify which underlying processes can induce changes in the maturity of those populations.
Maturation is a costly and complex physiological process that constitutes a key event in fish life history. Its schedule is determined by the individual genetic pool, which reflects selective pressures to which the population has been subjected to, and by environmental conditions that the individual experiences during maturation or at an earlier stage when the onset of the maturation process is triggered (abbreviated maturation "decision" in the following text) (Wright, 2007). Although it has long been described as dependent on size thresholds, a large corpus of knowledge underlines that environmental variations may also play a critical role through the modification of the individual's growth and energetic status at particular times of the life cycle (Thorpe et al., 1998;Wright, 2007;Tobin and Wright, 2011). Following initial growth, ageand size-at-maturation are determined by both the ability of fish to store energy and how they will allocate it between growth and other functions. Considering such dependence implies that external factors play an important role in the process of maturation either through phenotypic plasticity or evolutionary changes.
Among the potential drivers of these changes, several studies highlighted the role of water temperature (Grift et al., 2003;Kraak, 2007) and food availability (Trippel, 1995;Law, 2000) which may affect the probability of maturing either directly [e.g., Tobin and Wright (2011) for temperature], or indirectly by increasing growth rates (Sinclair et al., 2002) and body condition (Morgan, 2004;Grift et al., 2007). The social structure and the size composition of a fish population have also been pointed out as factors that may influence the probability of maturing at a given age and size (Diaz Pauli and Heino, 2013). Superimposed on the effects of environment, fishing is commonly known as a factor also affecting fish maturation. Fishing induces changes in the demographic structure of exploited populations (Law, 2000) and may finally lead to both phenotypically plastic (e.g., faster growth rates resulting from lower density-dependent effects) and genetic changes. Moreover, beside decreasing population sizes (Law, 2000) that may alter intra-and inter-specific competition, fishing may also induce direct or indirect (e.g., through habitat modification) changes in food availability (Rijnsdorp and Van Leeuwen, 1996) and cause evolutionary change by selecting for genotypes less affected by fishing (Law, 2000). Although the observed trends toward earlier maturation in fish population have been described as the result of a single factor or the combination of multiple factors (Marshall and McAdam, 2007), it remains unclear to what extent these changes are due to phenotypic plasticity on the one hand or evolutionary change on the other. Due to their different response time-scales, ranging from changes within a single generation (phenotypic plasticity) to those taking several generations before becoming detectable (evolutionary responses), disentangling genetic and plastic influences on maturation is of primary importance for the sustainability of marine resources.
The age and length at which an individual matures are variable and because of this variation at the individual level, maturation must usually be considered probabilistically: at a particular combination of age and size, some fish may mature while others do not. In this context, Probabilistic Maturation Reaction Norms (PMRNs), introduced by Heino et al. (2002b), are commonly used to describe the probability of maturation as a function of age and size and possibly other life history traits such as body condition (Grift et al., 2003(Grift et al., , 2007Barot et al., 2004;Mollet et al., 2007;Marty et al., 2014). Thus, using this method and providing some assumptions, one can estimate maturation probability independently of the influence of survival and growth .
The European sardine, Sardina pilchardus, is a small pelagic fish widely distributed along the Northeast Atlantic shelf, from the English Channel to Mauritania, and in the Mediterranean Sea (Parrish et al., 1989). This species, which has an observed lifespan of up to 14 years, is a multiple egg batch spawner and is characterized by a fast growth and an early maturation. Several studies emphasized an extended spawning activity for sardine in the Bay of Biscay, from autumn to spring with two spawning peaks occurring in April and in October (Coombs et al., 2006;Stratoudakis et al., 2007), and individuals mostly maturing during their first 2 years of life and in particular as 1 year old.
Over the last two decades, changes in life history traits of sardine have been brought to light in several regions of European waters (Silva et al., 2006(Silva et al., , 2013Doray et al., 2018;Saraux et al., 2019). In particular, Doray et al. (2018) highlighted a decreasing trend in both body length and weight for the Bay of Biscay sardine stock without emphasizing any correlations either with local changes in environmental conditions or with fishing intensity. Analyzing this decrease in details, Véron et al. (2020) recently noticed a concomitant decline in body condition at the population level. Moreover, recent analyses also suggest a decreasing trend in the proportion of mature individuals at age one which could be related to changes in sardine growth and/or maturation (Véron, 2020). Based on those observations, it becomes essential to identify if changes in the maturation of sardine have occurred over the last decades and which could be the underlying processes in the Bay of Biscay.
Despite a recent increase in fishing mortality rate (Ages 2-5) toward a value stabilized above both natural mortality (0.44 yr −1 ) and F MSY [0.51 yr −1 ; ICES (2019)], fishing mortality in this stock has remained historically relatively low (ICES, 2019). Consequently, based on the life history characteristics of this species and the relatively low exploitation rate for this stock over the study period, we do not expect recent fisheries-induced evolutionary changes in maturation.
Since the maturation process is strongly linked to spawning seasonality and because the duration of a population's spawning season plays a critical role in reproductive success and can be negatively impacted by the age truncation effects of fishing (Anderson et al., 2008;Wright and Trippel, 2009), we first investigated if temporal changes in sardine spawning seasonality occurred within the population over the period 2003-2018. With this first step we identified months that permitted accurate investigation of the sardine maturation process in the Bay of Biscay. Then, we used the PMRN approach to analyze the sardine maturation schedules in the Bay of Biscay and to disentangle genetic and plastic changes in these schedules. We focused on two aspects of the maturation process: we first investigated the relative importance of body length and condition in sardine maturation, and then we further examined the effects of environmental variables that could influence sardine maturation through growth-independent phenotypic plasticity. Subsequent to this analysis, the unexplained part of variability in sardine maturation process would potentially be attributed to genetic variability since this factor was not accounted for.

Biological Data Collection
We used biological data collected both by research vessels (pelagic and demersal trawlers) and commercial fleets ( et al., 2014), respectively] are yearly conducted by Ifremer aboard R/V Thalassa since 2000 and 1997, respectively. Since these surveys entirely cover the Bay of Biscay, the samples are assumed to be representative of the whole spawning fraction of the sardine population in this area. The design of the PelGas survey, especially dedicated to monitor the pelagic ecosystem on the continental shelf of the Bay of Biscay in springtime, makes it the main source of data used in this study. This dataset was complemented with fish market samples collected in the fishing harbors for sardine mainly located in the northern part of the Bay of Biscay (southwest Brittany). Finally, the dataset used in this study includes morphometric and physiological characteristics of sardine throughout the year from more than 24,000 individuals collected over the period 2003-2018 ( Table 1). Independently of the season considered, the distribution of these morphometric characteristics was assumed representative of the whole population.
For each fish, standard measurements were taken including total body length, weight, sex, and maturity stage. Length and weight were rounded down to the nearest 0.5 cm and the lower gram, respectively. Sex and maturity status were determined following ICES guidelines (ICES, 2008) through macroscopic inspection of the gonads. The maturity scale was based on a six-stage key defined as follows: (1) immature, (2) developing, (3) pre-spawning, (4) spawning, (5) partial post-spawning, and (6) post-spawning (ICES, 2008). Following the macroscopic inspection of gonads and the extraction of otoliths, age determination was carried out visually by binocular microscopy at the "Laboratory of Technology and Fisheries Biology" of Ifremer (Lorient, France) using the number of winter rings and considering January 1st as the nominal birthday.
Because of a limited number of samples in our dataset, we decided to pool males and females together in this analysis in order to reduce uncertainty in model estimates. This is not expected to create temporal changes in PMRN since sex ratio in the samples was relatively stable over the study period (Table 1). Moreover, we only considered individuals of length between 7.0 and 21 cm to reduce the number of length classes with only mature or immature individuals and therefore ensure model convergence and increase the goodness of fit around p = 0.5. In order to have enough individuals for each cohort, analyses were carried out for cohorts from 2002 to 2017. Since macroscopic observation of gonads does not allow distinguishing upcoming spawners from other individuals when the "decision" to mature takes place, maturity status was therefore used to separate individuals within the population. Individuals that had reached maturity stage 2 were considered mature while those in stages 3-5 were assumed to be spawning.

Consideration of Body Condition
The effect of body condition on sardine maturation was investigated using the relative condition index (Le Cren, 1951). This index, identified as one of the best morphometric indicators of fish condition in marine ecosystems (ICES, 2016), corresponds to the ratio between the measured weight of an individual and its predicted weight from a length-weight relationship. To account for variations in the weight of mature individuals within the spawning season, which might result from both feeding cycle and spawning activity, observed weights of mature individuals of different maturity stages were standardized as if they had all been collected when in maturity stage 2 (gonad developing stage), accounting for both endogenous characteristics of individuals (length, sex, and maturity stage) and year. The resulting corrected weight (termed W corr ) was obtained by adding to the observed weight the difference between observed weight and predicted weight for an individual with the same length but with maturity stage 2 in the same year. The Le Cren condition index was then computed as follow: whereŴ i,s is the predicted body weight of an individual i, of a given length L i,s obtained from a sex-specific length-weight relationship for sex s. A detailed description regarding the computation of the sex specific length-weight relationship used in this study is provided in Véron et al. (2020).

Application of the Probabilistic Maturation Reaction Norm
While maturity ogives describe the fraction of mature individuals within a population at a given age and/or size in a given year, maturation ogives refer to the probability that an immature individual, which has survived and grown until a certain age and size, sexually matures during a given time interval (Heino et al., 2002b). Maturity ogives therefore describe a population "state" without distinguishing between first-time spawners and repeat spawners whereas maturation ogives are more related to the process itself (Heino and Dieckmann, 2008). In this study, maturation propensity has been investigated through the PMRN approach which assumes that mature and immature individuals have the same growth and survival rates within age class and cohort. For sardine, newly mature individuals (first time spawners) and repeat spawners cannot be distinguished through visual examination of the gonads. Age and maturity status are therefore the only information that can be used to separate those individuals within a population. Taking into account that sardine maturation takes place during the first 2 years, with almost all individuals maturing within their first reproductive season (i.e., as age 0 or young age 1) and given that these individuals constitute the bulk of sardine spawning biomass in the Bay of Biscay (ICES, 2019), maturity ogives appear relevant to deduce sardine maturation probabilities. For age 0, maturation ogives can be approximated by maturity ogives. This assertion can be extended to age 1, as only a minority of sardine mature at age 0 and the contribution of repeat spawners at age 1 can be considered negligible.
Contrary to most studies investigating fish maturation process, the reaction norm approach applied to sardine here only considers length and condition as explanatory variables. The procedure for statistical analysis of sardine maturation probabilities in the Bay of Biscay involves the five following steps: (A) data selection and estimation of spawning ogives, (B) estimation of maturity ogives, (C) computation of the reaction norm midpoints and estimation of confidence limits, (D) analysis of differences between the two spawning seasons, and (E) testing the significance of trends in reaction norm midpoints. Each step of this procedure is described in detail below.

(A) Spawning ogives and data selection
Because data collected just prior to or during the spawning period are needed to analyze fish maturation, the sampling period used to compute sardine maturation probabilities was based on the analysis of sardine spawning seasonality. The probability of spawning of Atlantic sardine was modeled using generalized additive models (GAMs, as implemented in the R-package mgcv, 1.8-17, Wood, 2011). GAMs with a binomial error distribution and logit link function were fitted to the probability of spawning as an anisotropic bivariate function of fish length and month: where Sp is the probability of spawning, L corresponds to the length class and m is the month. Because Véron et al. (2020)
Since schedules for the start of vitellogenesis can differ among individuals of a given population and may depend upon lengthclass, differences in spawning seasonality within the population were assessed using these models to predict the probability of spawning for three length-classes throughout the year. These length-classes were assumed to correspond to young, medium and old individuals within the population. The outcomes of this analysis highlighted the existence of two spawning peaks within the year (end of April and mid October) which were used as a basis to analyze sardine maturation probabilities.

(B) Maturity ogives
Statistical analyses were conducted separately for each spawning season previously identified. Since maturity status can be reduced to a binary response variable (mature or maturing vs. immature), the fraction of mature fish (O) was estimated by means of logistic regression (using generalized linear models assuming a binomial error structure), with the proportion of mature individuals as the dependent variable and cohort and length (Eq. 3) or cohort, length and condition (Eq. 4) as independent variables. The linear predictor was linked to the fraction of mature individuals (i) using a logit link function [logit where length (L) and condition (Kn) are treated as continuous variables while cohort (c) is a factor. For both models, the α values correspond to the parameter estimates. In order to comprehensively assess the effects of cohort, length, and body condition on sardine maturation, performances of both models in explaining the fraction of mature individuals were compared by computing both a likelihood ratio-test and the Akaike Information Criterion.

(C) Reaction norm midpoint and confidence limits
Equation (5) was used to estimate L p 50 , defined as the midpoint of the Probabilistic Maturation Reaction Norm. This particular point, which corresponds to the length at which an immature fish has a 50% probability of sexually maturing, was used to summarize the reaction norm. For each cohort, L p 50 estimates were obtained by solving the equation of model (4) for length with a probability O equal to 50%. An increase in L p 50 will indicate an increase in length at first maturity. The reaction norm width, here defined as the interval between the lengths that correspond to maturation probabilities of 25% and 75% (L p 25 and L p 75 , respectively) was also considered to investigate the strength of the link between body length and maturation "decision." The narrower this width is the stronger the influence of length on maturation probabilities. Because the estimation of sardine maturation probabilities relies on several successive steps, confidence intervals for L p 50 were estimated by bootstrapping data (Barot et al., 2004). Individuals in the original data set were resampled 2000 times with replacement, with stratification by cohort, to generate distributions of L p 50 . Confidence intervals were then derived using the 2.5-97.5 percentiles of the bootstrap distribution.

(D) Trend analysis and differences between spawning seasons
Temporal trends in maturation probability were assessed with a linear regression of the PMRN midpoints (L p 50 ) against cohort (c) as a continuous variable: where ε is a normally distributed error term. The estimated midpoints were weighted with the inverse of the variance of each midpoint in order to lower the influence of imprecise estimates on the regression. Variance estimates were obtained from the bootstrap output. Furthermore, differences in reaction norm midpoints between spawning seasons were investigated using the following simple linear model, also weighted by the inverse of the variance of bootstrapped estimates of L p 50 , where s corresponds to the spawning season (spring or autumn), treated as a categorical variable and ε is a normally distributed error term. The comparison of the envelope widths between spawning seasons was also considered to test for potential differences in the strength of the influence of body length on sardine maturation between spawning seasons.

(E) Environmental impacts on reaction norm midpoints
When a linear temporal trend in reaction norm midpoints was found, growth-independent plasticity of sardine maturation was investigated considering the effect of environmental factors that may potentially affect sardine maturation. Sea surface temperature [SST ( • C), averaged over the first 30 m] and chlorophyll-a [Chl-a (mg.m −3 ), used as an indicator of food abundance] were extracted from the Copernicus Marine Service 1 and came from the Atlantic Iberian Biscay Irish (IBI) model and satellite data, respectively. These variables were averaged by intervals of 3 months over the continental shelf of the Bay of Biscay (Bathymetry < 330 m). Population abundances were considered as indicators of intra-and inter-specific competition that may affect sardine maturation. Long-term trends in sardine maturation process were investigated considering explanatory variables showing a significant temporal trend.
Since anchovy and sardine are commonly found together in the Bay of Biscay and feed on a similar wide variety of plankton organisms, abundance estimates for these two species obtained during the spring PelGas survey were considered in the analysis. The effect of abundance on maturation reaction norm midpoints was considered using two time series describing the total abundance of both species for the current year and the previous years.
Using separate linear models for each spawning season, L p 50,s,c were regressed against environmental factors e k,s, , abundance a s,c−d , and cohort c as continuous variables, where k identifies the environmental factor considered, is the time lag in months from the set of lags considered M (3 or 6 months before the peak of the given spawning season) and c − d is the year when abundance is estimated, with a time lag d (0 or 1 year). Similarly to previous models, ε s,c corresponds to the normally distributed error term. In order to account for nongrowth-mediated phenotypic plasticity in sardine maturation, before investigating possible evolutionary trends when analyzing long term trends in L p 50 , the cohort was only included after the selection of all other significant explanatory variables. Then, for each explanatory variable showing a significant effect on reaction norm midpoints over the study period, the absolute magnitude of change in L p 50 attributable to this effect was calculated as the product of the linear rate of change over time of this explanatory variable, the regression coefficient of L p 50 according to this variable and the number of years in the time interval.
Finally, since short-term fluctuations in both environmental conditions and population abundance may induce short-term variations in maturation propensity, fluctuations around the trend in L p 50 were investigated using Eq. (7). For each spawning season, we regressed the L p 50 [or its residuals from Eq. (5) if a long-term trend has been found] against explanatory variables for which the potential long-term trend was also removed (using residuals from regression of each explanatory variables against time).
For both analyses (long-term trend and short-term fluctuations), a stepwise procedure was implemented to remove non-significant variables with a threshold set at p-value = 0.05 (ANOVA; Zuur, 2009). The estimated reaction norm midpoints were weighted with the inverse of the variance obtained from the bootstrap procedure.

Spawning Modeling
Whatever the time period considered, our models described sardine spawning seasonality well as shown by the explained deviance ranging from 43.9 to 59.6% ( Table 2).
Because we wanted to investigate sardine maturation probabilities and since sardine mostly matures at 1 year old (represented here by the 16 cm length-class), the extent of spawning seasons has been delimited using a minimum threshold of 25% of individuals being spawning in the lower size-class. Based on this threshold, sardine spawning activity appears to be composed of two spawning seasons within the year with a first one extending from January to June and a second one occurring between September and December (Figure 1). Independently of the size-class, sardine exhibited a maximum spawning activity with a first peak in April-May and the second one in October-November. The magnitude of the spawning activity appeared size-dependent with small individuals showing both a shorter spawning season and a lower percentage of spawning individuals when compared with larger ones, in particular in autumn and early winter (Figure 1). Even if patterns in spawning seasonality showed no major difference across periods, our results suggested a slight shift toward latter and longer spring spawning season at the end of the study period (P3: 2012-2018). However, due to the small number of individuals collected in summer as reflected by the wide confidence interval (Figure 1), this predicted increase in the duration of the spring spawning season should be considered cautiously. Overall, our results do not show strong temporal shifts either in the peak or in the spawning activity pattern of sardine over the study period. Based on these outcomes, we characterized sardine spawning seasonality in the Bay of Biscay as being composed of two spawning seasons defined as follows: the first one during winter/springtime (from January to June, referred to as spring spawning season) and a second one occurring in autumn (from September to December). These seasons were used to derive maturity ogives at two different times of the year.

Maturity Ogives
For both seasons, models fitted the data well as indicated by the percentage of explained deviance. The fraction of mature individuals was better explained in springtime than in autumn with models accounting for around 60% and 35% of the deviance, respectively (Table 3). Cohort, size, and condition, as well as the 2-way interactions, significantly affect the probability of being mature when incorporated within the same model. As expected, whatever the season and the maturation model considered (length-or length and condition-based), body length was the best variable explaining sardine maturation (explained deviance ranging from 24 to 55.7%). Moreover, behind its strong dependence upon body length, sardine maturation process fluctuates among cohorts (explaining 5% and 7% of the total FIGURE 1 | Estimated probability of spawning for three length classes [16 cm (green), 20 cm (blue), and 23 cm (purple)] corresponding to the mean length at ages 1, 3, and 6 over the study period. The first three panels correspond to prediction considering the three identified periods in Véron et al. (2020Véron et al. ( ) (P1: 2003Véron et al. ( -2006P2: 2007-2011and P3: 2012 while the last one represents prediction considering the whole study period. deviance in spring and autumn, respectively) with a relationship between maturation and length that appeared cohort-dependent for both seasons, as indicated by the significant interaction between length and cohort.
Models in which sardine body condition was included performed significantly better than the model without it only in springtime, as shown by lower AIC (Akaike Information Criterion, significance based in chi-square test for likelihoodratio test statistic, p < 0.001). However, including body condition in addition to length resulted in a limited increase in predictive power (0.9% of the explained deviance in springtime). The significant interaction between condition and cohort reveals that the relationship between maturation and body condition is cohort-dependent for both spawning seasons and accounts for 0.5% and 0.76% of the deviance in spring and autumn, respectively.
Although only a small proportion of the deviance was accounted for by body condition, its positive effect on the probability to mature, in addition to size, remains statistically significant. Results are presented for one selected cohort in springtime only (Figure 2). In that case, there was a clear increase in the estimated probability of becoming sexually mature with an increase in body condition. Obviously, this increase remains strongly dependent on body length with an earlier increase for large individuals when compared with smaller ones (Figure 2A). While small individuals (13 cm length-class) do not exhibit drastic changes in the proportion of mature individuals under a "threshold" value of body condition (around 1.4), an increase of 0.4 units in body condition leads to a significant increase in the proportion of mature individuals in the 16 cm length-class (from around 25 to 75%, for condition values ranging from 0.7 to 1.1). When this proportion is estimated as a function of length for a small range of body condition (Figure 2B; 0.7, 0.9, and 1.1) corresponding to the median values of body condition observed during the three identified periods (P3, P2, and P1; respectively) in Véron et al. (2020), our results show a small increase in the proportion of mature individuals with increasing body condition which corresponds to a small decrease in length at maturation with an increase in condition.

"Spring" and "Autumn" Maturation Reaction Norms
Variability in maturation process between spawning season was investigated through the comparison of PMRN midpoints from models considering body condition. Our results show significant differences between seasons with higher average reaction norm midpoints in autumn than in spring indicating that at the same length, spring-spawning sardines have a higher maturation probability than autumn-spawning ones (ANOVA, F-value = 33.2, p = 2.44e-6; Figure 3A). On average, spring spawning sardines reached sexual maturity at 14.3 cm while in autumn L p 50 was 15.6 cm.
Envelope widths were significantly different between autumn and spring (two-sided paired Student's t-test: p = 0.003). While the average distance between L p 25 and L p 75 was around 4 cm for the former, it was less than 2 cm in springtime ( Figure 3B). Moreover, contrary to the spring reaction norm width which appears relatively stable over time, our analyses emphasize a larger variability in this envelope for autumn-spawners (onesided Fisher's F-test: p = 5.6 10 −9 ; Figure 3B). Altogether, these results highlight a stronger link between maturation "decision"   and body length in spring than in autumn underlying a more size-dependent probability to mature in spring than in autumn when sardine have lower and more variable tendencies to mature at a given size. While regression analyses of trends in maturation reaction norms (Eq. 5) do not show evidence for temporal trends in autumn L p 50 , a significant decreasing trend was found in spring reaction norms (midpoint: p < 0.01, width: p < 0.02; Figure 4A). The magnitude of the changes over the study period was quite substantial: from 2002 to 2018, L p 50 shifted from more than 15 cm to less than 14 cm.

Reaction Norm Midpoints and Environmental Effect
Among all explanatory variables tested in this analysis (temperature, chlorophyll-a, and fish abundance), significant trends over the study period were only found in both series of fish abundance (Figure 5; p < 0.036 for both series). These variables were consequently selected as potential factors affecting long-term trends in sardine spring maturation. No significant trend was found in series of temperature or chlorophyll-a, which were therefore only considered to investigate short term trends in both seasonal maturation reaction norms.
In contrast with autumn L p 50 which do not display any temporal trend, the long term decrease in spring L p 50 was significantly related to the total population abundance (anchovy and sardine) estimated the year prior to the observation of the maturity status (Table 4 and Figure 4A). Moreover, the inclusion of this variable explained enough of the variation in L p 50 to make the subsequently added cohort effect insignificant. From cohort 2002 to 2017, the rising population abundance experienced the year before maturation accounts for 15% of the decrease in L p 50 . Despite being non-significant, the cohort effect accounted for a reduction of 1.44 cm in spring sardine L p 50 corresponding to 69% of the downward trend in maturation reaction norm midpoints.
In contrast to spring, short term fluctuations in autumn L p 50 could be partly explained by short term variations in both the sea surface temperature and the chlorophyll-a 6 months before the spawning peak, and by total population abundance the year preceding the computation of the reaction norm (Table 4 and Figure 4B). The probability to mature at a given size in autumn increased significantly (p < 0.03) with an increase of both Chl-a and SST 6 months before the spawning peak and an increase in spring abundance the previous year.

DISCUSSION
The present study reveals no evidence that the recently observed decrease in length at maturity of sardine in the Bay of Biscay can be attributed to fisheries-induced evolution.
Instead, our results suggest that this decrease can be, at least, partly ascribed to an increase in the combined biomass of sardine and anchovy the year before. To our knowledge, this paper is the first to analyze sardine maturation processes in the Bay of Biscay by considering both spawning schedule and fish body condition. Overall, we show a significant increase in the estimated probability of becoming sexually mature with an increase in body condition. Our results serve to better understand the origin of the marked changes in the proportion of mature individuals at age 1 over the past decade.

Sardine Reproductive Period
Our results are consistent with the reproductive pattern described in the literature for the sardine population in the Bay of Biscay Lacroix, 1971, 1977;Coombs et al., 2006;Gatti, 2016). They confirm the extended reproductive period of sardine in the North East Atlantic and emphasize a synchronicity within the population with two spawning peaks occurring in spring and autumn (April-May and October-November, respectively). Taking into account both the key role of energy reserves in the reproductive process (Lloret et al., 2013) and the significant decline in sardine body condition over the last decades , the duration of the reproductive period was expected to be negatively impacted at the end of the study period. However, our results do not exhibit temporal shifts either in the duration of the spawning season or in the peak, which suggests that the investment toward the reproduction function has been kept constant over the study period.
The probability of spawning was both length-and seasondependent. In particular, smaller individuals exhibited lower probabilities of spawning than larger ones (whatever the season) and a shorter reproductive period. Thus, 1-year-old sardines maturing in spring appear to be smaller than those  Results are only presented for seasonal PMRN midpoints showing either long-term trend or short-term fluctuations. For each coefficient, F-statistic value (F) and associated p-value are given. Regarding the long-term trend analysis, the change in L p 50 , L p 50 (cm), that is attributable to each effect is also provided. Significant variables are shown in bold. maturing in autumn (Figure 3), the average difference in body length being 1.3 cm.

Impact of Body Condition on Maturation Process
Our analysis confirms that good condition has a significant positive effect on sardine maturation in the Bay of Biscay. These results are in line with the conclusions of various studies involving other fish species (Rowe et al., 1991;Marteinsdottir and Begg, 2002;Bromley, 2003;Morgan, 2004;Grift et al., 2007). Since energy allocation involves trade-offs between growth, reproduction and maintenance (Stearns and Koella, 1986), individuals exhibiting a higher body condition may devote more energy to reproduction making them able to make the "decision" to mature at a smaller size.
The influence of body condition on maturation was more substantial for sardine spawning in spring than those spawning in autumn. This is illustrated by a wider maturation envelope in autumn. Contrary to the autumn period when sardines have benefited from both spring and summer periods to store energy reserves (McBride et al., 2015), sardines emerge at the end of the overwintering period with, on average, a very low body condition (Gatti et al., 2018;Véron et al., 2020). The maturation envelope represents most of the combinations of size and condition at which sardine maturation can occur (Heino and Dieckmann, 2008) and the autumn maturation envelope could represent the strong variability in both growth rates and ability of individuals to store reserves. Alternatively, the size of the autumnal maturation envelope may be a consequence of the mixing of newly mature individuals (including both ages 0 and 1) and those that have already spawn. Assuming that the ability of an individual to enter sexual maturity relies on an energetic threshold value (Thorpe et al., 1998;Morita and Fukuwaka, 2006;Diaz Pauli and Heino, 2013), the amount of available energy in the beginning of the spring reproductive period will greatly influence the maturation of individuals. Such hypothesis is in agreement with several studies on salmonids which have pointed out how the accumulation of energy stores may act as a triggering signal for determining the onset of puberty (Rowe et al., 1991;Hutchings and Jones, 1998;Metcalfe, 1998).
The effect of body condition on maturity ogives appears furthermore dependent upon cohort and, to a large extent, on years considered over the study period. Indeed, when compared with the beginning of the study period (2002)(2003)(2004)(2005)(2006), when individuals exhibited a relatively good body condition , the lower body condition observed in recent years may have a stronger effect on sardine maturation propensity. This assumption is also supported by the observed difference between seasons in the significance of the cohort-dependent effect of body condition on maturation in Eq. (4). Indeed, the stronger significance of this interaction term in spring could mostly be explained by the fact that sardine enter the overwintering period with a lower body condition and thus a smaller amount of energy reserves than at the beginning of the time series , potentially resulting in a higher dependence of maturation upon individual condition.
Although we found that body condition has an effect on maturation, our results also show that this effect is lower than the effect of body length, which accounts for a much higher proportion of the explained deviance of maturation probability. Such outcome is in agreement with other studies showing similar results for other species (Marteinsdottir and Begg, 2002;Morgan, 2004;Grift et al., 2007). Predicted differences in the Lp 50 were relatively small over the observed range of body condition (0.7-1.1). This may indicate that despite the lower energy reserves observed in recent years the sardines have maintained their average size at maturation. Moreover, these results are in good agreement with a general prediction from the life-history theory that if fish are facing high mortality, selection will not only favor earlier reproduction but also higher reproductive effort at age and therefore length, at the expense of body growth and/or survival (Heino and Kaitala, 1999).

Differences in Maturation Lengths Between Reproductive Periods
Our results suggest that spring-spawning sardines had a significantly greater probability of maturing than autumn-spawning sardines of equivalent length (Figure 3). There are several reasons that may explain such differences. First, as confirmed in this study, sardine exhibit two spawning peaks throughout the year that occur before and after the main growth and energy storage season for sardines in the Bay of Biscay (spring and summer, McBride et al., 2015). Therefore and since sardine maturation takes place during the first 2 years of life, with most individuals maturing in their first year of life, individuals born in springtime and becoming mature in spring are more likely to be smaller when they arrive at their first reproductive period (age-1 fish) than those which will become mature in the following autumn (age-1.5 fish), after the productive planktonic season. Second, several studies emphasized the ability of fish to skip spawning if they have not reached both size and condition threshold values required at the time to make the "decision" to mature (Wright, 2007;Lowerre-Barbieri et al., 2011;Diaz Pauli and Heino, 2013). Therefore, taking into account that both body condition and growth have strongly declined over the study period , we hypothesized here that some small individuals skip their first spawning opportunity to increase their energy reserves in favor of both growth and postponed reproduction. As the "decision" to mature is made long before the maturation (Wright, 2007), and since it has been emphasized that the decrease in sardine body condition is mainly supported by both summer and autumn periods , such hypothesis could be relevant if the first reproductive period is in spring. This may partly explain the higher L p 50 highlighted by our results for the autumn spawning season.

Temporal Trends in Reaction Norm Midpoints
While changes in the proportion of young sardine mature individuals from cohort 2002 to 2017 have been detected (not shown), our analysis does not conclude in drastic changes in the PMRNs over time. We highlighted a quantitative decrease in spring L p 50 and no evidence for long term changes in autumn reaction norm midpoints has been found. As it has been underlined by several studies using PMRNs to disentangle the effect of phenotypic plasticity from the effect of genetic change (Heino et al., 2002a;Grift et al., 2003;Engelhard and Heino, 2004;Barot et al., 2005;Olsen et al., 2005), the cohort effect highlighted in this study could support the existence of evolutionary changes in spring maturation process. However, after including in the analysis a set of environmental variables capable of generating growth-independent plasticity in maturation, the effect of cohort on the spring reaction norm midpoints was no longer significant. The variability of environmental conditions affecting growth during the period between the date of the maturation "decision" and the date at which the maturity status is assessed appears as the most likely candidate to explain a major part of the interannual variability in L p 50 . The deviation from the assumption that negligible somatic growth occurs between these two dates is here too strong and the methodology carried out here does not allow to capture all growth-related variability in the maturation probability. This result has two implications. First, although genetic selection generated by fishing is commonly known as a potential factor contributing to changes in fish maturation (Trippel, 1995;Law, 2000;Wright, 2007), it suggests that temporal modifications in sardine maturity ogives do not result from such process over the period considered here. This is expected for a fish population that has not been subject to intense harvesting over the past decades. Second, it strengthens the need to consider environmental factors as additional explanatory variables that can drive plastic variation in PMRNs (Kraak, 2007;Marshall and McAdam, 2007), and in particular those impacting fish juvenile stages.

Phenotypic Plasticity of Sardine Maturation Process
Several studies emphasized that the onset of maturation may be affected by environmental factors such as temperature (Grift et al., 2003;Yoneda and Wright, 2005;Tobin and Wright, 2011), food availability or trends in abundance which may alter social structures (Kraak, 2007;Diaz Pauli and Heino, 2013). Moreover, these factors have been suspected to affect PMRN midpoints through their impacts during juvenile stages or more months/years prior to the maturation process (Mollet et al., 2007). Potential plastic responses of sardine maturation to environmental fluctuations were assessed at the population level by correlating PMRN midpoints with environmental variables (temperature, chlorophyll-a and total abundance of anchovy and sardine).
Analyzing long-term trends in spring reaction norm midpoints, our results emphasized a negative correlation between L p 50 and the spring abundance experienced the year prior to reaching sexual maturity. This indicates that abundance may have an effect on sardine maturation, driving individuals to mature earlier at a smaller size if population densities experienced during the critical "decision" phase are high. This counter-intuitive outcome is most likely another consequence of the deviation from the major assumption of negligible somatic growth between the time the "decision" is made and the observation of maturation, made in the PMRN approach. An increase in population abundance during the time interval between the "decision" to mature and maturation may have led to a decrease in growth rates, after the "decision" to mature, due to higher competition for food. Such process is in line with potential density-dependence effects that have already been suggested to occur within nursery grounds in the Bay of Biscay (Doray et al., 2018;Véron et al., 2020). Moreover, this result strengthens the idea that probabilistic maturation reaction norm approach for this case study might not completely account for growth-related plasticity (Morita and Fukuwaka, 2006;Dieckmann and Heino, 2007;Wright, 2007;Heino and Dieckmann, 2008) and in particular that occurring during the juvenile stages (Mollet et al., 2007;Diaz Pauli and Heino, 2013), when growth is very variable.
At shorter time scales, abundance seems to have a similar effect on the fluctuations in reaction norm midpoints. Our results emphasize that population densities experienced during springtime have a positive effect on the probability to mature in the following autumn reproductive period. Moreover, they also emphasized that short term fluctuations in autumn reaction norm midpoints could partly be explained by short term variations in water temperature and chlorophyll-a. Temperature (Charnov and Gillooly, 2004) and food availability are known to impact early life history stages, as they may affect age at maturity through their effect on juvenile growth rates. However, and as already highlighted by several studies (Grift et al., 2003;Dhillon and Fox, 2004;Mollet et al., 2007;Tobin and Wright, 2011), our results suggest potential effects of these environmental factors on maturation process at relatively short term. Since temperature has been demonstrated as a potential factor accelerating developmental rates other than through growth (Fuiman et al., 1998), the onset of maturation might therefore be affected by temperature in early life. However, like several studies which have examined the relationship between reaction norm midpoints and temperature (Grift et al., 2003;Kraak, 2007;Mollet et al., 2007;Morita et al., 2009), we cannot conclude whether temperature is having a direct effect on sardine maturation or is a proxy for other factors.

Toward New Perspectives to Improve Our Approach
Several of our findings should be interpreted with caution owing to the assumptions underlying the approach and the availability of data as potential explanatory factors.
First, a corpus of literature emphasized that the "decision" to mature relies on specific physiological threshold traits during the period of maturation (Thorpe et al., 1998;Wright, 2007;Tobin and Wright, 2011). Among them, fish energetic status before and during the critical phase of the developmental "decision to spawn" seems to play a critical role. Here, we used the relative condition index to account for the effect of sardine body condition on maturation process. Even if our results are in good agreement with those from several studies showing a positive effect of body condition on the probability to mature (Grift et al., 2007;Mollet et al., 2007;Diaz Pauli and Heino, 2013), they also emphasize that the decrease in sardine body condition over the study period cannot explain the whole temporal trend in maturation. Here, we highlight two main reasons for this. First, the timing of observations has been suggested as a potential source of misunderstanding of maturation process using the PMRN approach (Wright, 2007;Diaz Pauli and Heino, 2013). Since the developmental "decision to spawn" takes place long before spawning actually happens, the state of individuals when the initial maturation "decision" is made does not match the observed state in field-data. In particular, Wright (2007) highlighted that changes in energy status around the time of maturation "decision" cannot be accounted for in PMRN approach. In this context, Diaz Pauli and Heino (2013) notably suggested that the initiation of maturation stage for guppies (Poecilia reticulata) was closer to the maturation "decision" than to its completion. The authors argue that differences between initiation and completion of maturation might result from differential allocation of resources into growth and reproduction at the different stage and therefore conclude that PMRN defined around the initiation stage could better represent maturation schedule than that considering completion stage. Therefore, and since it has been suggested that maturation is controlled by successive inhibition through lipid-regulated switches during the critical period (Thorpe, 2007), we suggest here the use of those biochemical indices to identify the beginning of maturation in order to help understanding of maturation process. A second explanation of our results relies on the properties of the body condition index itself. In this study, data used to compute PMRN were selected during the identified spawning seasons. Since gonad development increases body weight at the time of spawning, it results (due to calculation process) in a higher body condition index. Such effect could therefore contribute to the positive effect of condition on the probability of maturation. This supports the difficulty to interpret the role of condition in the maturation process because with this in mind, a higher condition index can therefore be the consequence of maturation rather than the cause.
Second, our results confirm the existence of two spawning peaks for sardine in the Bay of Biscay. Moreover, they highlight significant differences in the reaction norm midpoints between the two reproductive seasons, with higher and more variable L p 50 in autumn than in spring. Those differences may be explained by the fact that sardines are multiple batch spawners and simple visual examination of the gonads or histological analyses do not allow distinguishing between first time spawners and individuals that have already spawned. Therefore, we could hypothesize that the higher PMRN midpoints may result from a mix between newly mature individuals and those that have already spawned in spring, grown during summer and skipped spawning on the second opportunity of the year, just before the autumn reproductive period, due to deficient diet and poor nutritional condition (Lowerre- Barbieri et al., 2011;Rideout and Tomkiewicz, 2011). With this result in mind, we suggest that the back-calculation of monthly growth trajectories from scales and otoliths of sardine may help to first set individuals on their own growth trajectory according to their birth season and second determine at which time they made the "decision" to mature through the identification of the shift in energy allocation between growth and reproduction.
This kind of back-calculation, that has already been used for some species (Engelhard et al., 2003;Baulier and Heino, 2008), may also help to better account for phenotypic plasticity in growth trajectories that can lead to the same age-size combination and which therefore constitutes a limit of the PMRN approach (Morita and Fukuwaka, 2006;Heino and Dieckmann, 2008;Diaz Pauli and Heino, 2013). Here we highlighted that consequence of such plasticity may be more important for species maturing mostly in their first year of life (old age 0 and young age 1). Indeed, as shown by Morita and Fukuwaka (2006), our results strengthen the idea that growth condition just before the initiation of maturation, here individual growth trajectories between the larval stage and age one, considerably influences maturation probabilities.
Finally, even if we considered a set of environmental factors that can influence sardine maturation process through growthindependent effects, the question remains as to what extent both the timing and the spatial windows were adequate with our analysis. This issue is strongly linked with the previous limit of our approach. In particular, Wright (2007) emphasized that temperature is likely to have the greatest direct influence during the "decision" period and therefore suggested that it should be examined when an individual initiates its maturation. In our case, we therefore recommend collecting individuals on nursery grounds some months before the spawning season in order to be able to understand the effect of environmental variables on sardine maturation process. Furthermore, due to the lack of data, our analysis was not spatially explicit. However, recent analyses of sardine population dynamics suggest different growth trajectories between the north and the south of the Bay of Biscay. Such a spatial pattern in growth may induce spatial patterns in the maturation process for sardine and could therefore involves strong repercussions on the dynamics of this stock and consequently its management if not accounted for.

Toward a Better Understanding of Stock Population Dynamics
This study clearly demonstrates that the strong phenotypic plasticity in sardine maturation process is tightly linked to variability in growth rates. Such variability originates from multiple sources among which environmental conditions play a prominent role and can consequently have a strong impact on both spawning stock biomass and population productivity (Dutil and Brander, 2003). Indeed, considering the current warming of sea surface temperature, climate change may lead to an increase in growth rate during the juvenile phase that may result in a decrease in the length at first maturation which in turn might impact adults growth since available energy has been routed toward reproduction at smaller size (Heino et al., 2002b). Moreover, even though a shift toward earlier maturation will positively affect spawning stock biomass, its consequences on the oncoming recruitment will not necessarily be as positive since fecundity and egg viability are positively linked with maternal size and condition (Trippel and Neil, 2011).
With all of this in mind, our results support the need to develop modeling tools allowing for the consideration of both temporal variations in biological processes and their relationships with biotic and abiotic drivers to better understand stock population dynamics. In the case of sardine, our study leads us to suggest the use of a length-based stock assessment model as a first step to account for both the variability in sardine maturation process and the spawning seasonality which may lead to several recruitments throughout the year and therefore strongly impact its population dynamics. Furthermore and more generally, considering that climate variation will impact fish biological processes such as maturation, a better understanding of its influence as well as its incorporation within fish stock assessment models will clearly help evaluate management options for sustainable fisheries under climate change (Hollowed et al., 2009). A future key step would therefore be to identify and quantify more precisely the environmental drivers of the maturation process in order to ultimately increase the predictive power of stock assessment models.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.