Contrasting life-history responses to climate variability in eastern and western North Pacific sardine populations
All procedures accorded to administrative provision of animal welfare of the Fisheries Research Education Agency Japan. All statistical tests used in this study are two-sided.Otolith samplesFrom the western North Pacific, age-0 JP sardine were collected from samples taken during acoustic and sub-surface trawl surveys in the offshore Oyashio region conducted during 2006–2010 and 2014–2015. The surveys were conducted by Japan Fisheries Research and Education Agency every autumn since 2005 which aim to estimate the abundance of small pelagic species. The abundance of young-of-the-year sardine in the region in the season, approximately 10–15 cm in standard length (SL), is considered a proxy for the abundance of recruits of the Pacific stock and used to tune the cohort analysis in stock assessment4. As representatives of the young-of-the-year population in the region, 2–6 trawl stations each year that had relatively larger catch-per-unit-effort were selected (Supplementary Fig. 1), and 9–20 individuals were randomly selected from each station for otolith analyses (Supplementary Table 1). Age of fish was initially judged by SL (10–15 cm) and later confirmed by the counts of otolith daily increments.From the eastern North Pacific, archived otoliths of CA sardine captured in cruise surveys and in the pelagic fishery of the Southern California Bight during 1987, 1991–1998, and 2005–2007 were collected. Fish in the size range of 10–16 cm SL were regarded as age-1 individuals born in the previous year, following Takahashi and Checkley56. The number of individuals varied between year classes in the range of 4–20 (Supplementary Table 2).Otolith processing, microstructure and somatic growth analysisSagittal otoliths were cleaned to remove the attached tissue in freshwater and then air-dried. Otoliths of JP sardine were embedded in epoxy resin (Petropoxy 154, Burnham Petrographics LLC) on slide-glass, while those of CA were glued to slide-glass using enamel resin and then ground and polished with sandpaper to expose the core. For some otoliths of CA sardine, the polished surface was coated with additional resin to facilitate identification of the daily increment width. Using an otolith measurement system (RATOC System Engineering Co. Ltd.), the number and location of daily increments were examined along the axis in the postrostrum from the core. Although daily increments were clearly observed until the otolith edge for JP sardine, it was difficult to do this for CA sardine probably because they had experienced winter when otolith growth slowed down. Therefore, the rings were counted as far as possible for CA sardine, which typically resulted in more than 150 counts. The first daily increment was assumed to form after 3 days post hatch (dph) for JP and 8 dph for CA sardine following Takahashi et al.26 and Takahashi and Checkley56. The otolith radius at each age was calculated by adding all the increment widths up to that age. Standard lengths at each age were back-calculated assuming a linear relationship between otolith radius and standard length using the biological intercept method34 as follows:$${SL}_{n}=left({{SL}}_{{catch}}-{{SL}}_{{first}}right)times left({{OR}}_{n}-{{OR}}_{{first}}right)/left({OR}_{catch}-{{OR}}_{{first}}right)+{{SL}}_{{first}}$$
(1)
where SLn is the standard length at age n, SLcatch is the standard length at catch, SLfirst is the standard length at the age of first daily increment deposition fixed at 5.9 mm for JP sardine and 5.5 mm for CA sardine following the previous studies26,56, ORn is the otolith radius at age n, ORfirst is the otolith radius at the age of first daily increment deposition, and ORcatch is the otolith radius at catch. Based on rearing experiments of field collected eggs, Lasker57 showed the SL of CA sardine at 6–8 dph ranged between 3.8 to 6.5 mm, and Matsuoka and Mitani58 showed the total length at 2–4 dph ranged between 4.8 to 6.2 mm, corresponding to 4.7 to 6.1 mm in SL. To deal with these uncertainties regarding the size at the age of first daily increment deposition, we conducted Monte Carlo simulations (10,000 times) to estimate the uncertainties of back-calculated SL, assuming that the initial SLs fall between 3.8 to 6.5 mm for both sardines. Standard deviations of the temporal back-calculated SL at each age were presented as the uncertainty of each SLn estimation, which varied between 0.51 and 0.73 at the end of larval stage (JP: 45 dph, CA: 60 dph), between 0.34 and 0.64 at the end of early juvenile stage (JP: 75 dph, CA: 90 dph) and between 0.20 and 0.53 at the end of late juvenile stage (JP: 105 dph, CA: 120 dph). These values were significantly smaller than the variability of estimated SL among individuals assuming initial sizes of 5.9 and 5.5 mm for JP and CA sardine, respectively (standard deviations: 4.2, 8.1 and 8.3 in JP sardine and 5.5, 9.1 and 10.3 in CA sardine for the end of larval, early juvenile and late juvenile stages, respectively), suggesting that the back-calculated SL is robust to variations of initial size. Nevertheless, the biological intercept method assumes a constant linear relationship between fish and otolith size within individual59, which can vary depending on physiological or environmental conditions60,61. Therefore, to examine the relationships between temperature and growth, we used both otolith growth, which contains fewer assumptions, and back-calculated somatic growth as growth proxies. Since the use of the two proxies did not show remarkable differences in the relationships between temperature and growth (Supplementary Figs. 11, 12), we mainly used the back-calculated SL in the discussion, which has a more direct ecological implication.To more generally test whether growth trajectories are different between the western and eastern boundary current systems, otolith growth data of JP and CA sardines were compared with those of sardines in the east to south and west coasts of South Africa. The biological intercept method to back-calculate standard length could not be used in sardine from South Africa because the size at catch was large, some over 20 cm, and otolith radius and standard length were not linearly correlated for fish of this size. Therefore, the otolith radius and increment width were directly used as proxy for size and growth in this comparison, respectively. For visualisation (Fig. 2a), the means of year class mean otolith radii were estimated for JP and CA sardines. For CA sardine, otolith radii at ages were simply averaged within each year class. For JP sardine, to account for the variation in the number of individuals captured at the same station, otolith radii were first averaged within each station, and the station means were averaged within each year, weighted by catch-per-unit-effort. For South African sardine, data of otolith daily increment widths from hatch to 100 dph of 67 adults captured at six stations on the east to south coast ( >22oE), and 51 individuals captured at six stations on the west coast ( 0.05). Theoretically, the relationship between metabolism and temperature tends to show a linear trend after the metabolic rate is log-transformed79. Thus, we applied “identity (data without transformed)” and “log (data transformed)” links to evaluate if model shows a better linearity with data transformation. Based on AIC, however, the result showed Moto have a better linearity without data transformation (Supplementary Table 7). We, therefore, used “identity” links for the further model selection. Model selection base on AIC was performed for models including temperature, region (JP and CA sardines), life history stages (larvae, early juvenile and late juvenile) and interactions of these factors. The full model including all the interactions had the lowest AIC (Supplementary Table 7). As the diagnostic for the full model showed normality and homogeneity of residuals (Supplementary Fig. 9), we selected this model for interpretation. The CA sardine at the larval stage as the baseline, we found only JP sardine at early and late juvenile stages has relatively higher Moto values, and the temperature-dependent slope is significantly gentler in JP sardine at early and late juvenile stages (Supplementary Table 8).Next, the diversity of Moto across temperature range was assessed to estimate the optimal temperature in each stage. The relationship between the maximum metabolic rate and temperature is known to be parabolic, while that between the standard metabolic rate and temperature is logarithmic28,79. As the highest field metabolic rate would be constrained by maximum metabolic rate and the lowest field metabolic rate would be close to resting metabolic rate43, fish would have the most diverse metabolic performance at the optimal temperature with the widest aerobic scope. Thus, we modelled the highest and lowest Moto values in each 1 °C bin using a polynomial regression and a generalised linear model with Gaussian distribution and a log link for the 95th and 5th percentile values of each bin, respectively (Supplementary Fig. 10). The values of the bin that included less than four values were excluded from the regression analyses to reduce the uncertainty caused by under-sampled temperature bins. The gap between the two regression lines was considered as a proxy for the aerobic scope, and the temperature at which the gap reached the maximum was regarded as the optimal temperature.Statistical analyses for the relationships between temperature and growthTo understand how variation in ambient water temperature affects early life growth of sardines, we compared back-calculated standard length at around the end of the larval stage (hatch–35 mm; JP: 45 dph, CA: 60 dph), the end of the early juvenile stage (35–60 mm; JP: 75 dph, CA: 90 dph), and the end of the late juvenile stage (60–85 mm; JP: 105 dph, CA: 120 dph) and the mean seawater temperature from hatch to the ages. Median of each sampling batch were used as minimal data unit. Pearson’s r and p-values were first calculated for each comparison (Supplementary Table 9). As the relationship between mean temperature and standard length of JP at 75 dph seemed to be dome-shaped rather than linear, we introduced quadratic term of temperature and tested whether the term increased explanatory power using a linear model and stepwise model selection based on AIC. The model selection showed that the full model (Standard length ∼ Temperature2 + Temperature) was the best model, and the coefficients of the quadratic and linear terms were both significant (Supplementary Table 10). To account for these multiple tests, we corrected the p-values of the coefficients of the quadratic term in the linear model for JP sardine at 75 dph and of the Pearson’s r for the rest using the Benjamini-Hochberg procedure with α = 0.05, and selected the null hypotheses that could be rejected (Supplementary Table 9). To compare the temperature that allow maximisation of growth rate and optimal temperature derived from the analysis of Moto for each stage, median somatic growth rate and otolith increment width in each 1 °C bin was calculated together with its 3-window running mean (Supplementary Figs. 11, 12).Statistical analyses for the relationships between sea surface temperature and survival indexTo test whether habitat temperatures during the first 4 months after hatch affect the survival of sardines in the first year of life on a multidecadal scale, satellite-derived sea surface temperature (SST) since 1982 and survival of JP and CA sardines were compared. The log recruitment residuals from Ricker recruitment models (LNRR)13, representing early life survivals taking into account the effect of population density, were calculated based on the stock assessment data for JP and CA sardines as follows:$${LNR}{R}_{t}={ln}({R}_{t}/{S}_{t}) , – , (a+btimes {S}_{t})$$
(6)
where LNRRt is the LNRR at year t, Rt is the recruitment of year-class t, St is the spawning stock biomass in year t, and a and b are the coefficients of linear regression of ln(Rt/St) on St. Pearson’s r between the LNRR and the mean SST values from March to June for JP and from April to July for CA sardine was calculated for each grid points in the western and eastern boundaries of the North Pacific basin, derived from a SST product based on satellite and in situ observations80 (Global Ocean OSTIA Sea Surface Temperature and Sea Ice Reprocessed (https://resources.marine.copernicus.eu/product-detail/SST_GLO_SST_L4_REP_OBSERVATIONS_010_011/INFORMATION), accessed on 11th August and 28th October 2021). The correlations were generally negative and positive in the western and eastern regions, respectively (Supplementary Fig 13a, b). In particular, mean SST values in the area where eggs, larvae and juveniles of JP or CA sardines are mainly found in the months26,39,49,56,78,81,82 (dotted areas in Supplementary Fig 13a, b) were compared with LNRR values to test the relationship between habitat temperature and survival in the early life stages (Supplementary Fig 13c). It should be noted that the mean SST values were not significantly correlated with otolith-derived year-class mean temperatures of JP and CA sardines during the larval to late juvenile stages (JP: r = 0.01, p = 0.98, n = 7, CA: r = 0.29, p = 0.38, n = 11), likely due to the short periods analysed, patchy distribution and inter annual variation in larval and juvenile dispersal and migration patterns. Nevertheless, the regions included areas where SST showed weak to significant (p More