## Abstract

Critical slowing down (CSD) reflects the decline in resilience of equilibria near a bifurcation and may reveal early warning signals (EWS) of ecological phase transitions. We studied CSD in the recruitment dynamics of 120 stocks of three Pacific salmon (*Oncorhynchus* spp.) species in relation to critical transitions in fishery models. Pink salmon (*Oncorhynchus gorbuscha*) exhibited increased variability and autocorrelation in populations that had a growth parameter, *r*, close to zero, consistent with EWS of extinction. However, models and data for sockeye salmon (*Oncorhynchus nerka*) indicate that portfolio effects from heterogeneity in age-at-maturity may obscure EWS. Chum salmon (*Oncorhynchus keta*) show intermediate results. The data do not reveal EWS of Ricker-type bifurcations that cause oscillations and chaos at high *r*. These results not only provide empirical support for CSD in some ecological systems, but also indicate that portfolio effects of age structure may conceal EWS of some critical transitions.

## 1. Introduction

Critical transitions in ecosystems include catastrophic change to a new stable state [1,2], destabilization of populations into cyclical or chaotic dynamics [3,4], and extinction [5]. All are examples of bifurcations [6], where a small change in conditions induces large changes in the steady state and dynamical behaviours of a complex system. Such critical transitions may be foreshadowed by early warning signals (EWS) arising from the interaction of stochastic and nonlinear processes [5,7,8]. Typically, signs of an approaching transition are due to critical slowing down (CSD): the resilience of equilibria (defined as the rate of recovery from perturbation) declines near bifurcations, which allows noise to increase the magnitude, shape and autocorrelation (AC) of system variability [7].

However, some critical transitions do not express EWS [9,10], and the intrinsic heterogeneity of ecological systems may play a role. In population dynamics, portfolio effects describe how heterogeneity in habitat, ecotypes and life history may dampen the influence of environmental stochasticity on population fluctuations [11–13]. In some cases, however, life-history heterogeneity, via age– or stage–structure, may also amplify population fluctuations owing to resonance arising from interactions between random perturbations and population structure or nonlinear processes [14–16]. It is therefore likely that heterogeneity in ecological systems modulates EWS of critical transitions by dampening or amplifying the effects of stochastic perturbation.

Simple population models can exhibit a rich array of bifurcations as the population growth parameter (*r*) changes [3]. For the common transcritical bifurcation where *r* crosses zero, CSD has been shown to provide EWS of extinction in experimental populations [5]. However, some models, such as the Ricker equation [17], *N*_{t + 1} = *N _{t}*exp(

*r*−

*bN*), where

_{t}*r*is the population growth parameter,

*b*controls density-dependent mortality and

*N*is abundance [17–19], also have destabilizing bifurcations as

*r*increases [3]. Such instability underlies a leading explanation for how fishing magnifies fluctuations in marine fish populations [20], owing to interaction between nonlinearity and stochasticity [21–23]. Collectively, these properties suggest that CSD may occur in the vicinity of multiple, but not necessarily all, critical transitions in the parameter space of a complex system. Particularly, we predict that variability and AC in populations with Ricker-type density-dependence would be a U-shaped function of

*r*. Because the Beverton–Holt model,

*N*

_{t + 1}= exp(

*r*)[

*N*/(1 +

_{t}*bN*)], does not exhibit this destabilizing bifurcation, variability and AC in populations with density dependence of this type are predicted to decline monotonically as

_{t}*r*increases.

Ecological studies that examine CSD mostly have been theoretical analyses of simulated dynamics obtained from a stochastic dynamical system under continuous change in a forcing parameter [8,9,24,25]. Data from laboratory populations [5,26,27] and a whole-lake experiment [28] support the theory of EWS, but empirical characterizations of CSD indicators across a parameter range from wild populations are lacking (but see [29]). Such an investigation is data-intensive, because it requires replicated timeseries of ecological dynamics in which to study CSD indicators via at least two different scales of contrast: (i) variation in proximity to bifurcations among replicated timeseries of stationary systems; and (ii) as non-stationary systems transit (or approach) bifurcations within their individual timeseries. While the latter has been empirically investigated in the context of a marine ecosystem regime shift in the Baltic Sea [29], studies of the dynamics of replicated natural ecological systems are particularly rare. Our study addresses this empirical gap using Pacific salmon (*Oncorhynchus* spp.) as a model system.

Populations of Pacific salmon (*Oncorhynchus* spp.) are known to vary among each other in their average population growth rates and therefore provide an opportunity to investigate if indicators of CSD are present among Pacific salmon populations. To theoretically characterize CSD in salmon populations, we first studied the influence of *r* on coefficient of variation (CV) and temporal AC in the Ricker and Beverton–Holt fisheries models. Both models contain a transcritical bifurcation at *r* = 0, but differ at high values of *r* where the Ricker model has unstable dynamics, but the dynamics of the Beverton–Holt model remain stable.

To briefly review these bifurcations, the positive equilibria for Beverton–Holt and Ricker models are *N** = (*e ^{r}* − 1)/

*b*and

*N**=

*r*/

*b*, respectively, and both models also have extinction equilibria at

*N**= 0. Local stability analysis reveals the stability of the equilibria based on their corresponding eigenvalues, , where

*f*(

*N*) is the recruitment function (i.e.

*N*

_{t + 1}=

*f*(

*N*)), and which measures how a small perturbation near an equilibrium, δ

*=*

_{t}*N*subsequently changes according to

_{t}– N**δ*

_{t + 1}=

*λδ*[30]. Equilibria are therefore stable if and unstable if [18]. The eigenvalues for the positive equilibria are for the Beverton–Holt and for the Ricker model, and both models have an eigenvalue at the extinction equilibrium that is . Thus, in both models, when

_{t}*r*crosses zero from above, the positive equilibrium is initially stable, and the extinction equilibrium is unstable, the two equilibria collide at

*r*= 0 where they exchange stability, so that the extinction equilibrium becomes stable (i.e. the definition of a transcritical bifurcation [6]). In addition, as

*r*approaches zero, the magnitude of

*λ*becomes larger, and therefore the rate of recovery from the perturbation modelled as

*δ*

_{t + 1}=

*λδ*becomes slower, and hence there is CSD. The models differ, however, when

_{t}*r*becomes large—the positive equilibrium in the Beverton–Holt model remains stable for all positive

*r*owing to its eigenvalue , but the Ricker model becomes unstable at

*r*= 2, where its eigenvalue becomes

*λ*= 1 −

*r*= −1, and therefore a flip (period-doubling) bifurcation occurs [6] where the equilibrium becomes unstable, and there is instead a period-two oscillation. CSD also precedes the flip bifurcation, because the magnitude of

*λ*becomes larger as

*r*approaches two, and hence rate of recovery is weaker, and in this case, is mixed with an oscillation, because

*λ*is negative. The Ricker model exhibits further period-doubling bifurcations and chaos as

*r*is increased further [3].

However, compared with the simple population models above, some salmon species exhibit heterogeneity in age-at-maturity [31], and such age– or stage–structure may have important consequences for fish population dynamics via nonlinear dynamics or portfolio effects [11,12,32,33]. To analyse this biocomplexity of salmon life histories in relation to CSD indicators, we first analyse simple Ricker and Beverton–Holt population models for the existence of CSD near bifurcations. Then, we study how stage-structured elaborations of Ricker and Beverton–Holt models may affect CSD indicators. Finally, we then examine historical timeseries of population abundance of Pacific salmon to characterize CSD indicators over an among-population gradient in average population growth rates, *r*. Results of these analyses revealed evidence of EWS of extinction in pink salmon (*Oncorhynchus gorbuscha*) and chum salmon (*Oncorhynchus keta*), but not sockeye salmon (*Oncorhynchus nerka*), where EWS may be obscured by heterogeneity in age-at-maturity.

## 2. Methods

### (a) Critical slowing down in Beverton–Holt and Ricker models

We studied CSD in the stochastic Ricker
2.1and stochastic Beverton–Holt
2.2models via computer simulation. In these models, *r* is the population growth parameter, *b* controls density-dependent mortality, *N* is abundance and *ɛ _{t}* is a normally distributed random variable with mean of zero and standard deviation,

*σ*, that represents environmental stochasticity [17–19,34]. For demonstration of concept, we set a carrying capacity of

*k*=

*r/b*= 1 000 000 and used

*σ*= 0.9, and simulated the model over a range of values of

*r*. The qualitative patterns are similar for other choices of

*σ*in the range of 0.4–1.2 that is commonly seen in the analyses of salmon population dynamics. Other choices of

*σ*do not change the predicted relationships between CSD indicators and

*r.*For each simulation, we allowed an initialization period of 450 iterations (this appeared sufficient based on visual inspection) and then analysed the final 50 iterations. CSD indicators were calculated from the final 50 iterations in the same manner as they were calculated from the salmon recruitment data by calculating the CV and lag-1 AC from ln(

*N*), the log-transformed abundance. For each value of

_{t}*r*(ranging from 0.01 to 3.0), this was repeated 1000 times to generate a distribution of CV and AC estimates from which we estimated the mean and 95% confidence intervals on CV and AC for each value of

*r.*The Ricker model was also used to evaluate the effects of fishing and pre-spawn mortality on the relationship between CV and

*r*(electronic supplementary material, figure S1).

### (b) Critical slowing down in age-structured models

To evaluate the influence that heterogeneity in age-at-maturity, common in chum and sockeye salmon, has on CSD indicators, we studied the dynamics of the following matrix population model of female salmon with census time at spawning in autumn:
2.3which written in matrix form is
2.4where **B** is the transition matrix, *F* is the abundance of age-1 offspring in a river (fry), that survive to the next year with probability *θ* to become age-2 offspring that have left the river (*J*; post-smoults). After their first winter at sea, the female smoults either remain at sea as subadults (*U*_{1}) with probability *γ*_{1} or return to their river to spawn as 1-sea winter adults (*A*_{1}) with probability *μ*_{1}. Of the *S*_{1} females that spend their second winter at sea, they either remain at sea as subadults (*U*_{2}) with probability *γ*_{2} or return to spawn as 2-sea winter adults (*A*_{2}) with probability *μ*_{2}. All remaining females that spend three winters at sea then return to spawn (*A*_{3}) with probability *μ*_{3}. All females die after spawning. The fecundity function we used was for Ricker-type dynamics and for Beverton–Holt-type dynamics to represent density-dependent mortality from egg-laying to first year of life. In these functions, *f _{i}* is the production of age-1 offspring by spawning salmon in class

*A*

_{1},

*A*

_{2}or

*A*

_{3}. The parameter

*ɛ*is a normally distributed random variable with mean of zero and standard deviation,

_{t}*σ*, that represents environmental stochasticity. Note that

*ɛ*is drawn randomly once for each year and applied to fecundity for each spawning age group, so fecundity among age groups is stochastic but covaries among ages.

_{t}We simulated the age-structured model (equation (2.3)) with just environmental stochasticity entering through *ɛ _{t}* and also through a full-stochastic model where each parameter was randomly selected in each year according to uniform distributions (electronic supplementary material, table S1). In each simulation, we generated a timeseries of 500 recruit observations from which we extracted the last 50 data for further analysis. Recruitment that corresponds to a particular brood year

*t*is

*R*=

_{t}*A*

_{1,t + 3}+

*A*

_{2,t + 4}+

*A*

_{3,t + 5}. From the last 50 recruit data, we calculated the CV and lag-1 as well as lag-4 AC of the simulated ln(

*R*) recruitment timeseries (where recruitment was log-transformed as in the empirical analysis). Lag-1 AC is expected to have the highest power, but lag-4 was also studied to correspond to the 4 year age at maturity of chum and sockeye salmon. Note that lag-1 for pink salmon is 2 years owing to their data structure. Simulations were repeated 1000 times each over a range of values of

_{t}*f*to compare CSD indicators across a range of the population growth parameter. The population growth parameter for the stage-structured model (which relates spawners and the recruits that they produce) is

*r*= ln(

*R*

_{0}), where

*R*

_{0}is the dominant eigenvalue of the next-generation matrix

**F**(

**I**−

**T**)

^{−1},

**F**contains the linearized fecundity elements of

**B**, and

**T**contains the transition probabilities of

**B**so that

**F**+

**T**=

**B**when

**N**→ 0 [35,36]. The distribution of age-at-maturity was approximately 20%

*A*

_{1}, 60%

*A*

_{2}and 20%

*A*

_{3}.

### (c) Salmon data

The data on the population dynamics of Pacific salmon consist of timeseries of spawner–recruit data for 120 stocks: *n* = 43 for pink salmon (*O. gorbuscha*); *n* = 40 for chum salmon (*O. keta*) and *n* = 37 for sockeye salmon (*O. nerka*). These data give the abundance of spawners in each year for each stock and the corresponding abundance of recruits, which is the abundance of adult offspring generated by the spawners and which returned to spawn or were caught in fisheries. The average lengths of the timeseries are 28 years for pink salmon, 27 years for chum salmon and 39 years for sockeye salmon, all distributed among regions of Washington, British Columbia and Alaska. Each pink salmon stock was further divided into odd and even year runs, and analysed separately (pink salmon have a 2-year life cycle, which creates odd and even year populations that are reproductively isolated in each river) [37,38]. A full description of the data as well as the full spawner–recruit dataset itself are available online in reference [38], where among-population variation in average *r* and within-timeseries non-stationarity in *r* (associated with a climatic regime shift [39]) is quantified for each population.

### (d) Empirical analysis

As a preliminary analysis, we compared the magnitude of variation in *r* at among-population scales relative to non-stationarity in *r* within individual timeseries to evaluate if there was a relatively large amount of variation among populations in *r* to investigate variation in CSD indicators. This was done using the tabulated results of Kalman-filtered estimates of within-timeseries variation in *r* published by Dorner *et al.* [38] from which we calculated and compared the average standard deviation in *r* within timeseries to the standard deviation of average *r* among populations.

For each stock, we estimated the CV as well as the lag-1 and lag-4 temporal AC of log-recruitment directly from the time-series data. To estimate the population growth parameter of each stock, we fitted a Ricker model for spawner–recruit data for each population in the database using the log-transformed version of the Ricker model:
2.5where *R _{i,t}* is the number of recruits in stock

*i*produced by spawners,

*S*, in year

_{i,t}*t*. The population growth parameter of population

*i*is measured as

*r*, density-dependent mortality within each stock is determined by

_{i}*b*, leaving normally distributed residuals

_{i}*ɛ*, which represent environmental stochasticity, with mean of zero and variance,

_{j,t}*σ*

_{i}^{2}, to be estimated. The use of the Ricker model here does not imply that the data generating processes were in fact governed by the Ricker equation, only that the population growth parameter can be reasonably estimated from the logarithm of recruits per spawner [40].

Estimates of CV, AC, *σ* and *r* for each salmon stock then formed the basis of our analysis of CSD indicators for each species. For each species, we looked for a declining or U-shaped relationship between CV and *r* as well as between AC and *r* using generalized additive models (GAMs) with Gaussian error and identity link function. First, to place CSD indicators and *r* on the same scale, we log-transformed the recruitment data prior to calculating the CV and AC values. To account for conventional (non-CSD) effects of environmental variability that might confound the analysis, we compared GAM models containing only *σ* as an explanatory variable versus GAM models that contained *σ* and *r* as explanatory variables. The models were fit using maximum-likelihood in the mgcv package in R (v. 2.11.1, www.R-project.org) [41] and compared using Akaike's information criterion (AIC) [42]. For AC, we followed the same procedure. To account for uncertainty in the estimates of *r*, we weighted the data in the GAM models inversely relative to the standard errors of the *r* estimates. We report results from weighted and unweighted analyses. To control for potential non-stationarity within the timeseries associated with a 1976–1977 shift in climate, we repeated the analysis with the data constrained to brood years after 1977. Timeseries prior to the regime shift were frequently too short for estimating *r*, and so analysis of timeseries prior to the regime shift was not pursued. To summarize, the steps involved in the analysis are as follows:

— for each salmon stock, estimate the CV and temporal AC of the log-recruitment timeseries;

— for each salmon stock, estimate the population growth parameter (

*r*) and environmental stochasticity (*σ*) by fitting a Ricker model to the spawner–recruit data pairs;— assemble the CV, AC,

*r*and*σ*estimates from each stock into a dataset, one for each species of salmon—pink, chum and sockeye; and— fit and compare GAMs for the relationship between the CSD indicator and

*σ*versus the CSD indicator and*σ*and*r*using maximum-likelihood and AIC.

## 3. Results

First, to describe the among-population and within-population variation in *r*, we compared among-population variation in average *r*-values with non-stationarity in *r* within individual timeseries from a Kalman filter analysis of these data published in Dorner *et al.* [38]. Among-stock variation in the average growth parameter was larger than the average of non-stationary variation in *r* within timeseries (measured as the average among stocks in their within timeseries standard deviation of the growth parameter, ): pink salmon = 0.57 versus = 0.24; chum salmon = 0.43 versus = 0.35; sockeye salmon = 0.44 versus = 0.31. Thus, the key ingredient we need to study CSD indicators—variation in *r* that spans parameter space near and far from bifurcation—is present in the timeseries.

In Ricker and Beverton–Holt models with no age structure, variability and temporal AC increased with proximity to the transcritical bifurcation at *r* = 0 (figure 1). Further, instability in the Ricker model at the flip bifurcation at *r* = 2 [3,6] also induced CSD (figure 1), giving rise in simulations to the hypothesized U-shaped relationship between variance and *r*. The effect of fishing, or other pre-spawn mortality, which occur in the salmon populations we studied, can affect the location of the minimum and curvature of the overall U-shaped pattern (electronic supplementary material, figure S1). Reddening the noise term (i.e. increasing AC) did not change the overall predicted relationships between CSD indicators and *r*, but it did cause the U-shaped relationship between CV and *r* to shift rightward and, unsurprisingly, also increased the overall level of AC (electronic supplementary material, figure S2).

Adding stage-structure to the population model to account for variation in age-at-maturity—common to chum and sockeye salmon [31]—can reduce EWS near *r* = 0 (where *r* = ln(*R*_{0}) and *R*_{0} is the dominant eigenvalue of the linearized next-generation matrix; see Methods; figure 2) although CSD owing to instability of Ricker-type dynamics at higher *r* persists*.* The effect of heterogeneity in age-at-maturity occurred in both the fully stochastic model (where all parameters were stochastic) as well as the model where stochasticity entered only through fecundity and fecundity among age classes was 100% correlated (electronic supplementary material, figure S3). Reddening the noise in the stage-structured model had the effect of slightly recovering patterns of CSD indicators and *r* observed in the unstructured model (electronic supplementary material, figures S4 and S5).

For the empirical analyses, the most reliable empirical results are those based on post-1977 data (after the regime shift) with data weighted by uncertainty in *r* estimates (column D of table 1), because these minimize influences of non-stationarity within timeseries and account for uncertainty. For the bifurcation at *r* = 0, pink salmon, which do not exhibit variable age-at-maturity, provided the most data near the bifurcation and revealed a negative relationship in among-population variation between CV and *r* and between AC and *r*, as predicted (table 1 and figure 3; electronic supplementary material, figure S6). However, sockeye salmon, which exhibit life-history heterogeneity via variable age-at-maturity [31], did not reveal evidence of CSD near *r* = 0 at the among-population level consistent with predictions from the stage-structured models (figure 3). Chum salmon, which also exhibit variation in age-at-maturity did, however, reveal evidence for a negative relationship between CV and *r* at low *r*-values as well as between lag-4 AC and *r* at low *r-*values (table 1 and figure 3). There was no evidence for a positive relationship between CV and *r* nor AC and *r* for any species (table 1 and figure 3), as would be expected owing to Ricker-type instability at high values of *r*.

## 4. Discussion

Our theoretical predictions of EWS of ecological phase transitions for salmon population dynamics indicate that heterogeneity in age-at-maturity may reduce EWS of extinction for both Ricker and Beverton–Holt-type dynamics. The mechanism is that, within one recruitment year, the effects of environmental stochasticity are averaged over several age classes. This averaging effectively dampens interannual variability caused by environmental stochasticity akin to how a diversified set of assets dampens the overall fluctuations of an investment portfolio. This so-called portfolio effect is relatively new to population dynamics [12] and is similar also to the statistical inevitability of diversity–stability relationships in community ecology [43]. Such portfolio effects also dampen the variance in recruitment during Ricker-type destabilizing bifurcations when the population growth parameter, *r*, is large and increasing, although the qualitative nature of increasing CV as *r* is increased through the bifurcation is preserved for the age-structured Ricker model. Thus, there is a variety of qualitative forms that CSD indicators in salmon population dynamics may theoretically take in relation to *r* that is mediated by the heterogeneity that is inherent to the species life history.

Stage–structure or multidimensional systems, in general, may also affect the existence of EWS in other ways. In particular, EWS owing to CSD near a catastrophic bifurcation in a simple predator–prey model is detectable only in the juvenile component of the prey population, but not in the adult [44]. Most fishery datasets only track abundance or catch of adults, and so in cases where the model in reference [44] is representative, the bifurcation may not be presaged, because the relevant variable has not be measured. Note that the mechanism in reference [44] is not a consequence of statistical averaging of asynchronous variability among population components (i.e. portfolio effects) but rather via a different mechanism: the orientation of eigenvectors and the value of their corresponding eigenvalues that characterize the dynamical system near the bifurcation point allow some variables to not express CSD [44]. Thus, our results provide a new mechanism that adds to EWS theory in that EWS in multidimensional systems are not universal, but rather dependent on the structure of the system.

Comparing CSD indicators across populations in relation to their average *r*-values revealed patterns of EWS of extinction that were consistent with theoretical predictions in two of the three species. For pink salmon, which have minimal variation in age-at-maturity, the data revealed evidence for CSD, both for CV and AC indicators. For chum salmon, we also observed evidence for EWS of extinction in CV as well as lag-4 AC, even though chum salmon have variability in age-at-maturity. By contrast, for sockeye salmon, which also have variable age-at-maturity, there was no evidence for EWS of extinction consistent with theoretical predictions. Our simulation study suggests that AC rises sooner than CV when *r* approaches zero, suggesting this may be a more reliable EWS of extinction. The empirical results on CSD were also stronger for AC than CV, suggesting that AC may be a more reliable EWS of extinction. Evidence that heterogeneity in age-at-maturity may diminish EWS of extinction mostly rests on our theoretical results; the power to detect effects of heterogeneity in age-at-maturity in data is limited by few replicates of species that vary in heterogeneity in age-at-maturity as well as variation among species in the quantity of data at low *r-*values.

The results did not reveal evidence for CSD at high values of *r* as was predicted by the Ricker model, but not by the Beverton–Holt model. However, this result is not strong evidence that instability of fish population dynamics at high *r*-values does not occur. For salmon population dynamics, potential overcompensation that would lead to Ricker-type dynamics could occur through overcrowding of spawning habitat in natural populations. Such dynamics may rarely occur in the data we studied, because fisheries harvests are likely to often reduce abundance of spawners, so that overcompensation occurs infrequently. This idea is consistent with our simulations where fishing makes CSD indicators flatten out at high *r-*values making CSD less detectible there, but also makes CSD at low *r-*values more detectable. More generally, instability of fish stocks owing to increased *r* is a main explanation for increased amplitude of fluctuations of exploited fish species [20,45].

The main objective of our study was to test for predicted directional relationships between CSD indicators and population growth—a qualitative test. Although this qualitative objective was achieved, we failed to observe clear quantitative correspondence between the observed and predicted values of CSD indicators. Why? One possibility is that other factors, such as variation in the local magnitude of environmental stochasticity may cause quantitative estimates of CSD indicators from the theoretical models to differ from observed estimates. Additionally, parameter values used in the theoretical models were chosen to reflect known variation in those parameters and to generate qualitatively testable predictions on the relationships between CSD indicators and population growth rates. The observed directional relationships between CSD indicators and population growth rates were indeed consistent with predictions, but not in all cases. Finally, because the empirical results are based on correlations from field data and not experimental manipulations, it is possible that there are other mechanisms that underlie the results.

We did not analyse the data for within-timeseries non-stationarity in *r* for EWS of the 1977 climate-associated regime shift [39]. We did not expect to see EWS of the regime shift in the salmon data, because the regime shift is an external environmental change that occurred abruptly relative to salmon population dynamics. Indeed, even if this regime shift caused *r* to approach or cross a bifurcation, then we would not expect to observe EWS of the transition, because the change would have been abrupt. Furthermore, according to theory one only expects to see EWS of the regime shift expressed via changing *r*-values in salmon-recruitment dynamics if the shift is a smooth transition over multiple spawning periods. These arguments highlight the need to better understand how the rate of environmental change and the associated speed of approach to bifurcation affect the expression of CSD.

While we compared empirical patterns of CSD with theoretical predictions from Beverton–Holt and Ricker models (and their stage-structured elaborations), our analysis does not require that these models were the data generating processes of salmon population dynamics. The true data generating processes are probably far more complex. Rather, these models illustrate and contrast how CSD indicators may be qualitatively related to the population growth parameter for unstructured and stage-structured populations, and for compensatory and overcompensatory density dependence. Similarly, the parameter values used in simulations were not empirically derived, but rather chosen to give representative values for fecundity, smoult survival and distribution of age at maturity [31]. Nevertheless, these models are routinely used in contemporary fisheries science and population ecology [16,18,19,46]. Indeed, a key feature of the theory of CSD is that its predictions are qualitative—a full understanding of the nonlinear dynamical mechanisms is not required for the theory to be applied.

While these results indicate that a rise in CSD indicators may be useful as EWS, it is important to note that the practicality of this approach also depends on the pace of change in *r*, which must be slow compared with the timescale of recruitment*.* Furthermore, our results pertain to bifurcations away from a stable equilibrium and not among unstable equilibria, where dynamics are more complex and the utility of EWS is far from clear. For EWS theory to be applied practically, one needs to look within non-stationary timeseries using a moving-window analysis to detect temporal increases in CSD indicators. Such an approach is more data intensive than the salmon data can support. This is because many recruitment observations are needed to calculate CSD indicators, which takes time, and it is possible that fast environmental changes could cause bifurcations to occur without the accumulation of data necessary to reveal EWS. This is partly the basis of why we do not expect salmon population dynamics to reveal EWS of shifts in the Pacific decadal oscillation. As such, age structure is not the only mechanism by which EWS may be eroded; fast environmental change is another. Overall, the practical use of early warning systems will depend on the speed at which the critical transition is approached, the severity of the nonlinearity in the approach to bifurcation, the magnitude of intrinsic noise, and the precision of the data. However, analyses from laboratory populations indicated that CSD could be detected as early as eight generations prior to the critical transition in an experimentally deteriorating model ecosystem [5].

Our main empirical results revealed patterns of CV and AC among timeseries from replicated wild populations of pink and chum salmon (but not sockeye) that were consistent with predictions of CSD at low values of *r* in theoretical models. However, we did not directly test for changes in CSD indicators within individual timeseries during a bifurcation, and as such it departs from other studies that typically change a parameter, so that the system transits a bifurcation and then analyses data for EWS of the bifurcation [8,9,24,25]. By contrast, we empirically characterized parameter space of salmon population dynamics using variation among distinct populations in their average distance from locations in parameter space where bifurcations are predicted by theory. EWS theory, according to our simulation results, yields testable predictions at this among-population scale that were evident in two of the three species we studied.

EWS of critical transitions may be useful for anticipating changes in system state or dynamics, including extinction [5]. Of course, if one has adequate data and well-validated models, then the approach to bifurcations may be anticipated directly by trends in parameter estimates alone—for example through state–space modelling [38,47]. In most situations, however, sufficient models and data will not be available, so that any inferences about impending system change must be based on more limited information. Our results indicate that partial data (i.e. recruitment data alone as opposed to full spawner–recruit data) in the absence of detailed dynamical models may be helpful, particularly for anticipating extinction in species with simple life cycles. We have also shown, however, that species with intact heterogeneity owing to life-history variation may conceal EWS of extinction. Thus, extinction may be less predictable in ecological systems with intact heterogeneity than in systems that have been simplified by human activities such as exploitation or habitat modification. However, although intact heterogeneity may conceal EWS of extinction it may also reduce the probability of extinction owing to reduced variability at low *r*. Furthermore, populations in degraded habitats may be more variable owing to the erosion of portfolio effects [12] as well as CSD associated with a decline in *r.* In conclusion, our results may not only be the first observations of CSD across a parameter range from natural populations, but also caution that age–structure effects may conceal EWS of some critical transitions.

## Funding statement

This work was supported by funding from the Natural Sciences and Engineering Research Council of Canada (M.K.) and the Odum School of Ecology, the Sarah Moss Fellowship at the University of Georgia, and a grant from the Leverhulme Trust (J.M.D.).

## Acknowledgements

We thank Randall Peterman for contributing the Pacific salmon spawner–recruit dataset as well as RM and referees for helpful comments that improved the work.

- Received December 9, 2013.
- Accepted March 25, 2014.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.