Introduction

Alternative tactics, utilized by many species at some point in their lives, describe discontinuous patterns of phenotypic variation, expressed by the same genotype, that result in divergent developmental programs (Roff 1996). For example, among adult dung beetles, individuals that have large horns behave as fighters, while those that do not become sneakers (Emlen 1997). In earwigs, small individuals have short forceps and large ones have long forceps, differences that originate due to the amount of proteins they acquired during development (Tomkins 1999). Salmonid fishes exhibit remarkable alternative life-history tactics. While some individuals migrate to the ocean to feed, grow to very large sizes, and then return to the natal stream to reproduce (migratory tactic), others complete their entire life cycle in freshwater streams (resident tactic) (Morita and Nagasawa 2010; Morita et al. 2014).

The coexistence of alternative tactics within a population is often explained as a status-dependent strategy (SDS; Gross 1996). If several tactics are available, and if the fitness of each tactic depends on the status (e.g., body size), the individual chooses the tactic with the highest expected fitness gain, given their status. Here, the phenotype of an individual switches based on whether the status is above or below a threshold. The value of the threshold should evolve to an evolutionarily stable strategy (Gross 1996).

A general trend observed throughout anadromous salmonids is that larger juveniles tend to remain in their natal streams (residents) (Piché et al. 2008). In a status-dependent strategy model for life-history tactics in salmonids (Morita et al. 2014), body size measurements, such as fork length and body mass, are considered to be the determinants for whether the fish stays in or leaves its natal stream (Dodson et al. 2013). Body size may depend on several factors, such as interactions among individuals, environmental stochasticity, and climate changes, and subsequently, these factors can influence the fraction of juveniles that employ each of the two tactics. For example, the growth of a juvenile can be suppressed by competition with other juveniles and adults in the stream, resulting in a decrease in the overall size and number of juveniles (Brown trout: Jenkins et al. 1999; Rainbow trout: Post et al. 1999). It can be predicted then that the number of juveniles developing to residents should decrease as the density increases.

Climate changes also affect migratory tactics in several different ways. Recently, global warming is considered to have affected environmental habitats, resulting in a shift of fish body sizes (Perry et al. 2005; Crozier et al. 2008; Daufresne et al. 2009; Satterthwaite et al. 2010). Naturally, we may imagine that a higher temperature might enhance the growth rate of juveniles due to the enhanced productivity of their environments (Crozier et al. 2008). However, Daufresne et al. (2009) demonstrated that an increase in the water temperature resulted in a decrease in the body size of aquatic species, due to their increased metabolic activity. The change in the growth rate caused the body sizes of the juveniles to change, which affected the life-history choices of the masu salmon, and subsequently, a change in the fraction of migrant fishes (versus resident fishes). Further, Finstad and Hein (2012) predicted that the fraction of migrants Arctic char that are important for the fishing industry will be reduced by increasing terrestrial productivity, under climate change scenarios.

In many models of life history choices, individuals are assumed to choose tactics based on their absolute status (such as body size), with a threshold value (Gross 1996). Tachiki and Koizumi (2016) proposed an alternative mechanism of life-history choices, called a relative assessment model, in which each juvenile makes a life-history choice based on their relative status within the juvenile population, estimated by the interactions among individuals. Tachiki and Koizumi (2016) showed that relative assessment stabilizes the fractions of the two tactics and that it is evolutionarily favored under a stochastic environment, in comparison with the absolute assessment model (Gross 1996).

In this paper, we consider the size of individuals affected by climate change, represented as a monotonic change of temperature, that might either increase or decrease their growth rate. We focus on salmonids, a taxon that exhibits alternative life-history tactics. We examine how their eco-evolutionary dynamics as well as the fraction of migrant tactics and the evolution of the threshold status responds to climate change. We first consider the case in which life-history choices are made based on absolute assessments, but we later examine the case with a relative assessment model. These two cases are very different and indicate the importance of understanding the underlying mechanisms of how individuals make life-history decisions in predicting the impact of climate change on wild life resources.

Life history choices affected by density-dependent growth

We begin with a mathematical model for the life-history choice of masu salmon (Oncorhynchus masou masou). Some males complete their entire life cycle in freshwater streams (residents), whereas others migrate to the ocean (migrants). In contrast, all females become migrants, as in several other salmonids (Dodson et al. 2013). Resident males can start reproduction in the year they hatch and continue to reproduce for multiple years (iteroparity). In contrast, migrant males mature in the ocean and return to their natal stream at the age of 2+ years, to reproduce once and then die (semelparity). Large juveniles tend to be residents (Piché et al. 2008), which has been regarded as the optimal tactic in the status-dependent strategy (SDS) model (Morita et al. 2014). The life-history choice of juveniles depends on whether their body size exceeds a threshold value, and the growth of juveniles is suppressed when there are many older resident males. An abundant number of residents in a year should reduce the growth of juveniles and result in fewer juveniles staying in the stream, leading to fewer resident males in the subsequent years. This causes a delayed feedback of fish density in the stream and may result in fluctuations of the juvenile body size and of the fractions of migrants and residents (Horita et al. 2018). In this section, we explain the population dynamic model for males and summarize the results that will be the basis of the analysis of the evolutionary outcomes.

The choice of life history by juveniles in absolute assessment model

Let x be the body size of a juvenile when the life-history choice is made. We assume that x has a normal distribution with mean \({\mu }_{t}\) and standard deviation \(\sigma\):

$${\widehat{p}}_{t}\left(x,{\mu }_{t}\right)=\frac{1}{\sqrt{2\pi }\sigma }\mathrm{exp}\left(-\frac{{\left(x-{\mu }_{t}\right)}^{2}}{2{\sigma }^{2}}\right)$$
(1)

The mean body size \({\mu }_{t},\) may vary between years, being affected by competition for food and space with other river residents, and by environmental conditions.

We assume that a juvenile becomes a resident if its body size x is larger than a threshold \({x}_{\mathrm{T}}\); it becomes a migrant if x is smaller than \({x}_{\mathrm{T}}\) (absolute assessment model). The probability to be a resident is

$$X\left(x|{x}_{\mathrm{T}}\right)=\left\{\begin{array}{c}1\ x\ge {x}_{\mathrm{T}}\\ 0\ x<{x}_{\mathrm{T}}\end{array}\right..$$
(2)

Let \({D}_{t}\) be the density of stream residents in year t. It is obtained by integrating the size distribution of individuals employing the resident tactic, multiplied by \({N}_{t}\), the number of juveniles in year t. We have

$${D}_{t}=\sum_{\tau =1}^{\infty }{{S}_{\mathrm{R}}}^{\tau }{N}_{t-\tau }{\int }_{-\infty }^{\infty }\left.{X\left(x|{x}_{\mathrm{T}}\right)\widehat{p}}_{t}(x,{\mu }_{t-\tau }\right)dx.$$
(3)

The size distribution of residents in year t is given as follows:

$${n}_{\mathrm{R}, t}\left(x|{x}_{\mathrm{T}}\right)=\sum_{\tau =1}^{\infty }\left.{{S}_{\mathrm{R}}}^{\tau }X\left(x-\tau \Delta l|{x}_{\mathrm{T}}\right){N}_{t-\tau }{\widehat{p}}_{t-\tau }(x-\tau \Delta l,{\mu }_{t-\tau }\right).$$
(4)

where \({S}_{\mathrm{R}}\) and \(\Delta l\) are annual survivorship and size increments of the residents, respectively.

Density-dependent growth and ecological dynamics

In the stream, when the density of residents is higher, the growth of juveniles becomes slower, due to competition for resources. Hence, the mean body size of juveniles is a decreasing function of the density of the residents (e.g., Brown trout: Jenkins et al. 1999; Rainbow trout: Post et al. 1999). Horita et al. (2018) considered the case when the mean body size of juveniles \({\mu }_{t}\) decreases with \({D}_{t}\), that the density of the resident males is as follows:

$${\mu }_{t}={\mu }_{0}+{\mu }_{1}\mathrm{erf}\left[k\left(\theta -{D}_{t}\right)\right],$$
(5)

where parameter \({\mu }_{0}\) determines the intermediate value of the body size \({\mu }_{t}\). Hereafter, we call it the intermediate size. We may regard \({\mu }_{0}\) as an index of juvenile growth rate. In all of the analyses in this paper, we assume that \({\mu }_{0}\) is larger (or smaller) if the juvenile growth rate is faster (or slower), due to climate changes. \({\mu }_{1}\) determines a range of the function \({\mu }_{t}\). \(\theta\) is the population density when the mean body size \({\mu }_{t}\) is equal to \({\mu }_{0}\). \({\mu }_{t}\) converges to \({\mu }_{0}+{\mu }_{1}\) when the density is very low, and to \({\mu }_{0}-{\mu }_{1}\) when it is very high. \(k\) determines the slope of the function, and we call k, the “strength of ecological feedback.” The strength of ecological feedback (k) affects the stability of the fraction of juveniles who become migrants. When k is small, the fraction of migrant males is stable. As \(k\) increases, the fraction of migrants becomes unstable, and fluctuates over a period of 2 years. Further increases in k make the amplitude of the fluctuation larger (Horita et al. 2018).

Since \({\mu }_{t}\) in Eq. (5) is a function of the density of male residents, the size distribution of juveniles in year t is a function of \({D}_{t}\), expressed as \({{p}_{t}\left(x,{D}_{t}\right)=p}_{t}\left(x,{\mu }_{t}\left({D}_{t}\right)\right)\). The residents in subsequent years comprise the resident juveniles of the current year who remain in the stream, plus the surviving residents in the current year. We define the population dynamics of \({D}_{t}\) as

$${D}_{t+1}={S}_{\mathrm{R}}N{\int }_{-\infty }^{\infty }X\left(x|{x}_{\mathrm{T}}\right)\left.{p}_{t}(x,{D}_{t}\right)dx+{S}_{\mathrm{R}}{D}_{t},$$
(6)

where the first term on the right-hand side indicates the contribution of juveniles from the current year who remain in the stream and \({S}_{\mathrm{R}}{D}_{t}\) represents that by the residents in the current year. Equation (6) indicates that the density of residents in year t + 1 can be calculated from the density in year t.

The population dynamics Eq. (6) may converge either to a stable equilibrium or to a 2-year cycle with a large amplitude (Horita et al. 2018). This can be explained in an intuitive manner as follows: when the threshold value \({x}_{\mathrm{T}}\) is clearly higher than the mean size of the juveniles, most juveniles migrate to the ocean and will be back to the stream later for reproduction. Then, only a small fraction of juveniles will remain in the stream and become residents in the following years, causing reduced growth suppression on juveniles. This creates a delayed negative feedback in the density of residents and may produce 2-year fluctuations of the mean body size of juveniles, the fraction of migrants, and the density of residents. In contrast, if the threshold body size \({x}_{\mathrm{T}}\) is rather small, compared with the mean body size of juveniles, a large fraction of juveniles remain in the stream and become residents in the following year.

Monotonic changes of the environment

We have considered the effect of the slow and steady change of the environment, which causes monotonic increases or decreases in the intermediate size, \({\mu }_{0}\). In the numerical analyses, we start with stable environmental conditions, in which the intermediate size \({\mu }_{0}\) is fixed for a very long time. Then, climate changes start to modify the intermediate size \({\mu }_{0}\), as follows:

$${\mu }_{0,t+1}={\mu }_{0,t}+mt.$$
(7)

The mean body size was stopped when it reached a particular value. We assume that the mean size is a linear function of time for simplicity (Prentice et al. 1993). The parameter m indicates the rate of change in \({\mu }_{0}\). m is positive if the environmental change enhances the juvenile growth, but negative if it reduces it. Gregory et al. (2017) showed that the body size of Atlantic salmon increased more than 10 (mm) in 20 years (0.5 [mm/year]). According to these literatures, we adopted the plausible value as the parameter of body size change (\(0.01 \le m\le 0.6\) [mm/year]). We investigated climate change for 100 years as in many other studies (e.g. Finstad and Hein 2012; Lehodey et al. 2013, 2015). We classified the dynamical behavior of the system in terms of ecology and evolution. In experimental studies, when the water temperature increases, the growth rate first increases and then starts to decrease due to metabolic disorder (Satterthwaite et al. 2010). In our simulation, we consider a monotonical change to simplify the model.

The fraction of migrants decreased with an increase in the intermediate body size \({\mu }_{0}\) (Fig. 1a). The change of stability in the ecological dynamics may occur both from equilibrium to fluctuation and from fluctuation to equilibrium. Hereafter, we call it “regime shift.” In Fig. 1b, the fraction of migrants was stable initially, but it became destabilized as \({\mu }_{0}\) changed. In Fig. 1c, the fraction of migrants fluctuated initially within the 2-year period, but it stabilized as \({\mu }_{0}\) changed. When the intermediate size decreased (Fig. 1d, e, and f), the fraction of migrants increased (Fig. 1d), and the ecological dynamics may change from an equilibrium to a 2-year cycle (Fig. 1e and f).

Fig. 1
figure 1

Fraction of male juveniles to be migrants when the threshold is fixed. Climate changes occurred between \({T}_{1}\) and \({T}_{2}\). The signs “E” and “C” indicate that the ecological dynamics are equilibrium and cycle over a 2-year period, respectively. In a, b, and c, intermediate size \({\mu }_{0}\) increased. In the other cases, intermediate size decreased. (a) The fraction of migrants decreased with an increase in \({\mu }_{0}\). Parameters: feedback strength \(k=0.1\), the rate of climate change \(m=0.6\), survivorship of residents \({S}_{\mathrm{R}}=0.05\), and that of migrants \({S}_{\mathrm{M}}=0.5\times {10}^{-4}\). (b) The stability of the fraction of migrants that changed from equilibrium to cycle. Parameters: \(k=0.5\), \(m=0.6\), \({S}_{\mathrm{R}}=0.02\), and \({S}_{\mathrm{M}}=7.0\times {10}^{-4}\). (c) The stability of the fraction of migrants that changed from cycle to equilibrium. Parameters: \(k=0.5\), \(m=0.6\), \({S}_{\mathrm{R}}=0.02\), and \({S}_{\mathrm{M}}=0.5\times {10}^{-4}\). (d) The fraction of migrants that increased with a decrease in \({\mu }_{0}\). Parameters: \(k=0.1\), \(m=-0.4\), \({S}_{\mathrm{R}}=0.05\), and \({S}_{\mathrm{M}}=0.5\times {10}^{-4}\). (e) The stability of the fraction of migrants that changed from equilibrium to cycle. Parameters: \({\mu }_{0,0}=160.0\), \(k=0.5\), \(m=-0.6\), \({S}_{\mathrm{R}}=0.05\), and \({S}_{\mathrm{M}}=0.5\times {10}^{-4}\). (f) The stability of the fraction of migrants that changed from cycle to equilibrium. Parameters: \(k=0.5\), \(m=-0.4\), \({S}_{\mathrm{R}}=0.02\), and \({S}_{\mathrm{M}}=0.5\times {10}^{-4}\). Other parameters were \({\mu }_{0,0}=100.0, {\mu }_{1}=50.0, m=0.01, \theta =0.4, \space\)\(\sigma =13.0,\Delta l=40.0, N=100.0, =0.1\), and \({c}_{1}=5.0\)

Evolutionary dynamics

To explore the evolution of the threshold trait, \({x}_{\mathrm{T}}\), which is expressed only in males, we consider the competition among males for mating opportunities. We adopted the fitness function formalized by Tachiki and Koizumi (2016), in which the fitness of a male is defined as a lifetime siring success for each tactic. Notably, population growth is independent of the threshold trait and population size does not change, because all females are migrants. We hence assumed that N, the number of new born individuals, is a constant. The evolutionary change in the threshold trait might affect the fraction of males who exhibit alternative life-history tactics.

Evolutionary change of the threshold trait

The spawning group of masu salmon consists of a female and several males including both migrants and residents (Watanabe et al. 2008; Esteve et al. 2011). Migrant males attack other males, whereas resident males adopt a sneaking tactic (Koseki and Maekawa 2000; Watanabe et al. 2008). Further, sneaker males attack each other, resulting in a dominance hierarchy among sneaker males (Koseki and Maekawa 2000; Esteve et al. 2011). Due to this interference competition, the competitive ability of resident males with body size x becomes a decreasing function of the density of resident males larger than x. We here assume that the competitive ability of a resident male with size x is given as follows:

$${C}_{t}\left(x|{x}_{\mathrm{T}}\right)={c}_{0}\mathrm{exp}\left[-{c}_{1}{\int }_{x}^{\infty }{n}_{\mathrm{R},t}\left(y|{x}_{\mathrm{T}}\right)dy\right],$$
(8)

where \({c}_{0}\) is the relative competitive ability of the strongest resident relative to a migrant and \({c}_{1}\) is a constant that mediates the strength of the interference of the competition among the residents. \({n}_{\mathrm{R},t}\left(y\right)\) is the size distribution of the residents (Eq. 4). Equation (8) could be interpreted as the product of \({c}_{0}\) and the probability that an individual with size x does not encounter larger individuals in the same spawning group. The number of spawning females and eggs produced does not change between years, because all females are migrants and the sex ratio and survivorship are constant. The migrants return to their natal stream to reproduce once after several years in the ocean for growth. Miyakoshi and Saitoh (2011) reported that the adult size and ocean survival of masu salmon are independent of their size at smoltification.

We assume here that the competitive ability of a migrant can be expressed as the product of constant survivorship, \({S}_{\mathrm{M}}\), and the size-independent competitive ability of the migrant, which is normalized to 1. We also consider that the situation is one where there are a few rare mutants in a population dominated by the wild type. Let the trait value of mutant \({x}_{\mathrm{T}}^{\prime}\) be slightly different from that of wild type \({x}_{\mathrm{T}}\). The fitness of a mutant (i.e., invasion fitness) is given as follows (Tachiki and Koizumi 2016):

$${W}_{t}\left({x}_{\mathrm{T}}^{\prime}|{x}_{\mathrm{T}}\right)={f}_{t}\left({x}_{\mathrm{T}}^{\prime}\right)/{f}_{t}\left({x}_{\mathrm{T}}\right),$$
(9)

where

$${f}_{t}\left({x}_{\mathrm{T}}^{\prime}\right)={S}_{\mathrm{M}}{\int }_{-\infty }^{\infty }\left(1-X\left(x|{x}_{\mathrm{T}}^{\prime}\right)\right)p\left(x,{D}_{t-2}\right)dx\ +{\int }_{-\infty }^{\infty }X\left(x|{x}_{\mathrm{T}}^{\prime}\right)\sum_{\tau =0}^{\infty }{{S}_{\mathrm{R}}}^{\tau }C\left(x+\tau \Delta l|{x}_{\mathrm{T}}\right)p\left(x,{D}_{t-\tau }\right)dx.$$
(10)

The first term on the right-hand side of Eq. (10) is the reproductive success of migrants. The survivorship, \({S}_{\mathrm{M}}\), is multiplied by the probability that a juvenile migrates to the ocean. Notably, the mating success of migrants when they successfully return to their natal stream is normalized as unity. The second term of the right-hand side of Eq. (10) is the contribution of residents. It is the product of survivorship during \(\tau\) years and the competitive ability summed over the years and averaged with the probability distribution.

Here, we adopt the evolutionary dynamics as follows:

$${{x}_{\mathrm{T}}}_{t+1}={{x}_{\mathrm{T}}}_{t}+\Lambda {\left.\frac{\partial }{\partial {x}_{\mathrm{T}}^{\prime}}\mathrm{ln}W\left({x}_{\mathrm{T}}^{\prime}|{x}_{\mathrm{T}}\right)\right|}_{{x}_{\mathrm{T}}^{\prime}={x}_{\mathrm{T}}}$$
(11)

where \(\Lambda\) is equal to the additive genetic variance (Iwasa et al. 1991) and mediates the rate of the evolution. Concerning the threshold size of maturation in salmonid fishes, there is no appropriate dataset estimating the additive genetic variance. However, there is a literature that estimates the genetic variance of threshold for smoltification in rainbow trout (Phillis et al. 2016). While there is substantial difference between threshold of smoltification and that of early maturation, we tested \(\Lambda\) for a wide range including the value estimated by Phillis et al. (2016) and examined how \(\Lambda\) modified the eco-evolutionary consequence.

Evolution of threshold size in a constant environment

Figure 2 indicates the selection gradient of the threshold trait \({\mathrm{x}}_{\mathrm{T}}\), i.e., the direction of evolutionary change. Here, all of the parameters are constant (no climate change). The horizontal axis is survivorship, \({S}_{\mathrm{R}}\) in Fig. 2. The evolutionary attractors may be single values, or two different values, depending on \({S}_\mathrm{R}\) (Fig. 2a). The system exhibits saddle-node bifurcation in which bistability of evolutionary outcomes occurs for a particular range of the parameter \({S}_{\mathrm{R}}\) (see Horita et al. 2018). These two evolutionary attractors differ in the stability of ecological dynamics (i.e., time series of the fraction of migrants). When the threshold trait is at the evolutionary attractor of a larger \({x}_{\mathrm{T}}\), the ecological dynamics exhibit a fluctuation with 2-year period. The evolutionary attractor of a smaller \({x}_{\mathrm{T}}\) corresponds to the ecological dynamics of a stable equilibrium (Horita et al. 2018).

Fig. 2
figure 2

Evolutionary consequences of threshold trait. The horizontal axis is the survivorship of resident \({S}_{\mathrm{R}}\). The vertical axis is the threshold trait. The signs “E” and “C” indicate the same as those in Fig. 1. Upper (lower) arrows indicate that the trait values increase (decrease) in evolution. Solid curves are the set of convergence stable points (attractor). The dotted curves are the set of evolutionarily unstable points (repeller). The ROB (range of bistability) indicates the range for which the evolutionary bistability occurs. (a) In the case of absolute assessment, for \({S}_{\mathrm{R}}<0.026\) or \(0.052<{S}_{\mathrm{R}}\), the evolutionary dynamics of \({x}_{\mathrm{T}}\) converged into a single stable point. For \(0.026\le {S}_{\mathrm{R}}\le 0.052\), the evolutionary dynamics of \({x}_{\mathrm{T}}\) converged into one of the two stable values, depending on the initial values. Open and solid triangles correspond to the parameter sets of Fig. 3a and b, respectively. Other parameters were \(k=0.5, {\mu }_{0}=100.0, {\mu }_{1}=50.0, \theta =0.4, \sigma =13.0,\)\(\space{S}_{\mathrm{M}}=2.5\times {10}^{-4},\Delta l=40.0, N=100.0, {c}_{0}=0.1\), and \({c}_{1}=5.0\). (b) In the case of the relative assessment, the evolutionary dynamics of \({z}_{\mathrm{T}}\) always converged to a stable equilibrium. The number of Bayesian update was \(n=100.0\). Other parameters were the same as those in Fig. 2a.

The threshold size evolves to the evolutionarily stable value. When the survivorship of resident males in the stream, \({S}_{\mathrm{R}}\), is high, the threshold \({x}_{\mathrm{T}}\) evolves to be small and the fraction of resident males evolves to be large, resulting in stable population dynamics (denoted as “E” in Fig. 2a). In contrast when the survivorship of resident males, \({S}_{\mathrm{R}}\) is small, \({x}_{\mathrm{T}}\) evolves to be large, and the fraction of migrant males fluctuates wildly between high and low values over a 2-year period (denoted as “C” in Fig. 2a). Interestingly, when the survivorship of residents is an intermediate value, the system is evolutionarily bistable. Depending on the initial value, either a large threshold size, \({x}_{\mathrm{T}}\) with a fluctuating fraction of migrants/residents evolves, or a low threshold size \({x}_{\mathrm{T}}\) with a stable and large fraction of migrants (i.e. a small fraction of residents) evolves.

Figure 3 also indicates the evolutionary outcomes, \({x}_{\mathrm{T}}\) in a constant environment in which the horizontal axis is the intermediate size of juveniles \({\mu }_{0}\). As \({\mu }_{0}\) increases, the threshold value \({x}_{\mathrm{T}}\) that is convergence stable increased accordingly (Fig. 3a and b).

Fig. 3
figure 3

Evolutionary consequences of threshold trait. The horizontal axis is the intermediate size \({\mu }_{0}\). The signs “E” and “C” indicate the same as those in Fig. 1. Arrows indicate the same as those in Fig. 2. (a) In the case of absolute assessment with only one evolutionary equilibrium (in Fig. 2a, it is donated by the open triangle), the convergence stable point increased with an increase in \({\mu }_{0}\). Survivorship of resident was \({S}_{\mathrm{R}}=0.02\). Survivorship of migrant was \({\mathrm{S}}_{\mathrm{M}}=0.5\times {10}^{-4}\). (b) In the case of the absolute assessment with three evolutionary equilibria (in Fig. 2a, it is donated by the solid triangle). Two attractors and one repeller increased with an increase in \({\mu }_{0}\). Parameters were \({S}_{\mathrm{R}}=0.05\) and \({S}_{\mathrm{M}}=0.5\times {10}^{-4}\). (c) In the case of relative assessment, the convergence stable threshold value did not correlate with \({\mu }_{0}\). Parameters were \({S}_{\mathrm{R}}=0.08\) and \({S}_{\mathrm{M}}=2.5\times {10}^{-4}\). The number of Bayesian update was \(n=100.0\). Other parameters were \(k=0.5, {\mu }_{1}=50.0, \theta =0.4, \sigma =13.0,\Delta l=40.0, N=100.0, {c}_{0}=0.1\), and \({c}_{1}=5.0\)

Eco-evolutionary dynamics in changing environments

In this section, we investigate the impact of slow and steady environmental changes over a period of years on the eco-evolutionary dynamics of the life-history choices of masu salmon. We ran numerical analyses of the eco-evolutionary dynamics under the conditions in which the intermediate size \({\mu }_{0}\) was changed according to Eq. (7).

Evolutionary trajectories of threshold \({x}_{\mathrm{T}}\) with climate change

Figure 4 illustrates the evolutionary trajectories of the system, in which \({\mu }_{0}\) changes due to environmental changes and these drive the evolutionary changes in \({x}_{\mathrm{T}}\). Here, the horizontal axis is for \({\mu }_{0}\), and the vertical axis is for \({\mathrm{x}}_{\mathrm{T}}\). The parameters were selected as those that correspond to Fig. 3b. For each value of \({\mu }_{0}\) (horizontal axis), there are two evolutionary stable values of \({x}_{\mathrm{T}}\) (evolutionary attractors), and one evolutionary unstable value (evolutionary repeller) located between them (Fig. 3b). Only the difference between Fig. 4a and b is the rate of change in \({\mu }_{0}\) (designated by m).

Fig. 4
figure 4

Evolutionary trajectories. The horizontal and vertical axis represent the intermediate size \({\mu }_{0}\) and threshold trait \({x}_{\mathrm{T}}\), respectively. Black arrows indicate the evolutionary trajectories of the system. The solid gray line indicates the convergence stable threshold, and the dotted gray line is the repeller. In the dashed shading region, which is donated by sign “C,” ecological dynamics fluctuated over a 2-year cycle, and in the white region, denoted by the sign “E,” ecological dynamics were at an equilibrium. (a) When the rate of the evolution was faster than the rate of climate change, the evolution could keep up with climate change. The rate of climate change was \(m =0.01\). (b) When the rate of the evolution was slower than the rate of climate change, the threshold trait converged to the lower attractor by crossing the repeller. The rate of climate change was  \(m=0.1\). Other parameters were \(\Lambda =100.0, k=0.5, {\mu }_{\mathrm{0,0}}=100.0, {\mu }_{1}=50.0, \theta =0.4, \sigma =13.0;\)\(\space {S}_{\mathrm{R}}=0.05, {S}_{\mathrm{M}}=0.5\times {10}^{-4},\Delta l=40.0, N=100.0, {c}_{0}=0.1;\)  and \({c}_{1}=5.0\)

If the environment does not change, \({\mu }_{0}\) remains at a particular value, and the evolution of the system should converge to either one of the two attractors, as shown in Fig. 3b. However, climate change causes an increase in \({\mu }_{0}\), making the system deviate from the attractor. As a result, the trajectory of the threshold trait (\({x}_{\mathrm{T}}\)), illustrated by the black arrows, always stayed a little below the attractors (gray lines), instead of exactly on them (Fig. 4). When climate change ends and \({\mu }_{0}\) becomes fixed to its final value, the evolution can make straight forward movements to the attractor, which is indicated as a vertical trajectory in Fig. 4a and b.

If the rate of change in \({\mu }_{0}\) is sufficiently small in comparison with the rate of evolution (Fig. 4a), the evolutionary trajectory during environmental change has a slope similar to the three lines indicating two attractors and a repeller. This implies that if the system starts from near the attractor, corresponding to stable equilibrium population dynamics, then the evolutionary trajectory kept close to the attractor. In contrast if the system started near the attractor corresponding to a 2-year cycle, then the system keeps close to the attractor with the 2-year cycle. Hence, we can conclude that the population dynamics remain unchanged during and after climate change, if the rate of the evolution is sufficiently fast compared with climate changes.

Figure 4b shows the case where the rate of the evolution is slower than that of climate change. Here, climate changes make the system move from the attractor corresponding to the fluctuating population dynamics to the attractor corresponding to the stable equilibrium dynamics. Hence, the population dynamics exhibit a transition of their dynamical behavior.

Regime shift in the eco-evolutionary dynamics caused by climate change

Figure 5 shows the eco-evolutionary dynamics under the same parameter set as that of Fig. 3a, in which there was only one evolutionary equilibrium which was convergence stable. When the rate of the evolution \(\Lambda\) was too small to follow the change in the intermediate size \({\mu }_{0}\) (i.e., the rate of climate change was too fast to keep up with), the evolution of the threshold trait could not follow climate change. Hence, the difference between the trait value adopted in the population and the evolutionary convergence stable strategy, under given \({\mu }_{0}\), increased while climate change occurred (\({T}_{1}<t<{T}_{2}\); in Fig. 5a). This results in a large fluctuation of the fraction of migrants among males (Fig. 5b). Even after climate change stopped (\(t>{T}_{2}\)), the evolution of the threshold still continued, and gradually converged to an equilibrium. Taking 4 or 5 times longer transient time than that of climate change, the threshold trait \({x}_{\mathrm{T}}\) finally reached the convergence stable strategy. During this long transient (Hastings et al. 2018), the stability of the ecological dynamics changed once the frequency of migrants converged to a particular value, and then fluctuated again and finally returned to stable 2-year cycle with the same amplitude with that before climate change (Fig. 5b). Such a complex pattern of dynamical behavior was realized because the combination of trait value \({x}_{\mathrm{T}}\) and the intermediate size \({\mu }_{0}\) was in the parameter region for a 2-year cycle most of the time, but in the one for stable equilibrium for some period of time (Fig. S1 in online resource 1). The fraction of migrants recovered to the same level as before the environmental change (Fig. 5a). Depending on the parameter set, there were various patterns of ecological dynamics generated with a change in intermediate size \({\mu }_{0}\). See Appendix A and Fig. S2 (Online Resource 1) for a detailed explanation.

Fig. 5
figure 5

Eco-evolutionary dynamics with climate change in the case of single equilibrium. The signs “\({T}_{1}\),” “\({T}_{2}\),” “E,” and “C” indicate the same as those in Fig. 1. (a) The solid line is the trajectory of threshold trait \({x}_{\mathrm{T}}\), and the dotted line is the intermediate size \({\mu }_{0}\). The threshold value increased with evolution with an increase in intermediate size. The interval of values of intermediate size and thresholds had been changed before and after climate changes. (b) The stability of the time series of the fraction of migrants was temporarily changed, such as fluctuated to equilibrium with the effects of climate change, and then the stability returned to the initial level by threshold evolution. Parameters were \(k=0.5, m=0.6, {\mu }_{\mathrm{0,0}}=100.0, {\mu }_{1}=50.0, \theta =0.4, \sigma =13.0, {\Lambda =10.0,}\)\({S}_{\mathrm{R}}=0.02, {S}_{\mathrm{M}}=0.5\times {10}^{-4},\Delta l=40.0, N=100.0, {c}_{0}=0.1\), and \({c}_{1}=5.0\)

Figure 6 shows the eco-evolutionary dynamics in the case of Fig. 3b, in which there are three evolutionary equilibria, one is the repeller (dotted line in Fig. 3b) and others are convergence stable (solid lines in Fig. 3b). Before climate change started (\(t<{T}_{1}\)), the threshold trait \({x}_{\mathrm{T}}\) was equal to an evolutionary attractor (Fig. 6a, which corresponds to the upper attractor in Fig. 4b,). In this state, the fraction of migrants oscillates over a 2-year period (Fig. 6b). When intermediate size \({\mu }_{0}\) increased with time (\({T}_{1}<t<{T}_{2}\)), the threshold trait once increased slightly but then suddenly decreased (Fig. 6a). This occurred when the threshold trait crossed the repeller (tipping point [Hastings et al. 2018]), and it attracted another evolutionary equilibrium (Fig. 4b). At the same time, the fraction of migrants was stabilized (Fig. 6b). Even after the environmental changes ceased (\(t>{T}_{2}\)), the evolution of the threshold trait \({x}_{\mathrm{T}}\) continued and finally stopped at the evolutionary equilibrium after long transient. In this case, the ecological dynamics were qualitatively different before and after the environmental changes, caused by the transition of the evolutionary attractor, from the upper one to the bottom one in Fig. 4a. Hereafter, we call it “evolutionary regime shift.” For a more detailed explanation see Appendix A and Fig. S3 (Online Resource 1).

Fig. 6
figure 6

Eco-evolutionary dynamics with climate change and the case of bistability. The signs “\({T}_{1}\),” “\({T}_{2},\)” “E,” and “C” indicate the same as those in Fig. 1. (a) The solid line and dotted line indicate the same as those in Fig. 5a. By increasing the value of the intermediate size, the regime shift of the threshold occurred and it converged to the other attractor. The interval of values of intermediate size and threshold had been changed before and after climate changes. (b) The stability of the timeseries of the fraction of migrants changed from unstable to stable. The fraction of migrants is finally fixed at a specific value, which is different from the initial situation. Parameters were \(k=0.5, m=0.6, {\mu }_{\mathrm{0,0}}=100.0, {\mu }_{1}=50.0, \theta =0.4, \sigma =13.0, {\Lambda =50.0,}\)\({S}_{\mathrm{R}}=0.05, {S}_{\mathrm{M}}=0.5\times {10}^{-4},\Delta l=40.0, N=100.0,{c}_{0}=0.1\), and \({c}_{1}=5.0\)

This evolutionary regime shift was caused by the balance between the rate of change in the intermediate size \(m\), and that of evolution \(\Lambda\). We examined whether the evolutionary regime shift occurred or not, by running the eco-evolutionary dynamics. When the rate of the evolution was large enough, the evolutionary regime shift did not occur (Fig. 7). If the evolution was too slow to follow climate change, the trait value deviated from the attractor converged formerly and then reached the repeller. When the trait value crossed the repeller, it was attracted by another evolutionary attractor.

Fig. 7
figure 7

The conditions under which the evolutionary regime shift occurs. The horizontal axis and vertical axis indicate the rate of threshold evolution \(\Lambda\), and rate of climate change \(m\), respectively. The black area indicates where evolutionary regime shifts (threshold converges to the other attractor) have occurred. The white area indicates where evolutionary regime shifts have not occurred. When the rate of the evolution was large, the evolutionary regime shift of the threshold by climate change hardly occurs. Parameters were \(k=0.5, {\mu }_{\mathrm{0,0}}=100.0, {\mu }_{1}=50.0, \theta =0.4, \sigma =13.0, \space\)\({S}_{\mathrm{R}}=0.05, {S}_{\mathrm{M}}=0.5\times {10}^{-4},\Delta l=40.0, N=100.0, {c}_{0}=0.1\), and \({c}_{1}=5.0\)

Relative assessment of the rules

Animals learn their relative status or social rank in a population from experiences by interacting with other individuals and adjust their behavior accordingly (Fricke and Fricke 1977; Warner and Robertson 1978; Kuwamura and Nakashima 1998; Munday et al. 2006). In the relative assessment models (Tachiki and Koizumi 2016), each individual is assumed to estimate its relative status in the population z by Bayesian inference. Let z be the difference of an individual’s own status to the population median status (\(z=0\) if its status is equal to the median of the population). We assume that individuals have a prior distribution for z that follows a normal distribution, with mean 0 and standard deviation \(\sigma\). It obtains a difference of its own status (i.e., body size) and that for a randomly chosen opponent, and then calculates posterior distribution. Each individual updates posterior distribution n times (Bayesian update). Each juvenile chooses the resident tactic if the estimated status is higher than a threshold value \({z}_{\mathrm{T}}\).

Let \(\Phi \left(x|{z}_{\mathrm{T}}\right)\) be the probability that an individual with true absolute status x, estimates its own relative status z, to be larger than the threshold \({z}_{\mathrm{T}}\). According to the derivation in Appendix B (Online Resource 1), \(\Phi\) becomes

$$\Phi \left(x|{z}_{\mathrm{T}}\right)=\frac{1}{2}\left(1+\mathrm{erf}\left[\frac{n\left(x-\mu \right)-\left(n+2\right){z}_{\mathrm{T}}}{2\sigma \sqrt{n+2}}\right]\right).$$
(12)

In the case of the relative assessment, we obtain the same statistics from Eqs. (3) and (4) by replacing \({x}_{\mathrm{T}}\) and \(X\left(x|{x}_{\mathrm{T}}\right)\) with \({z}_{\mathrm{T}}\) and \(\Phi \left(x|{z}_{\mathrm{T}}\right)\), respectively.

Equation (12) is a sigmoidal function of true size x, and the number of Bayesian update n mediates the slope of the function. In the limit of very large n, \(\Phi\) becomes a step function in which \(\Phi\) is 1 if \(x>{z}_{\mathrm{T}}+\mu\), and otherwise \(\Phi\) equals 0.

We have examined the relative assessment models by applying all the analyses adopted to absolute assessment models in previous sections. First, in relative assessment models, the change in the status of all the juveniles in the population simultaneously does not impact on the fraction of migrants, as illustrated in Fig. 8a. In addition, the ecological dynamics are always stable. Second, we can compute the evolution of threshold traits in the relative assessment models, with the same scheme as for the case of absolute assessment in the previous sections, but the variable to evolve is \({z}_{\mathrm{T}}\), instead of \({x}_{\mathrm{T}}\); the results were as follows:

Fig. 8
figure 8

Ecologically and eco-evolutionary dynamics with climate change in the case of the relative assessment. The signs “\({T}_{1}\),” “\({T}_{2}\),” and “E” indicate the same as those in Fig. 1. (a) When the threshold trait is fixed. The frequency of migration was not changed by climate change. (b) The dotted line is the intermediate value \({\mu }_{0}\), and the solid line is the threshold trait, \({z}_{\mathrm{T}}\). When the threshold trait evolves, the value of \({z}_{\mathrm{T}}\) hardly changed with an increase in intermediate size. (c) The fraction of migrants did not change with an increase in intermediate size in the case of b. Parameters were \(k=0.5, {\mu }_{\mathrm{0,0}}=100.0, {\mu }_{1}=50.0, m=0.6, \theta =0.4, \sigma =13.0, {\Lambda =50.0,}\) \({S}_{\mathrm{R}}=0.08, {S}_{\mathrm{M}}=2.5\times {10}^{-4},\Delta l=40.0, N=100.0, n=100.0, {c}_{0}=0.1,\)  and \({c}_{1}=5.0\)

[1] There was only one stable attractor (Fig. 2b). There was no parameter set that showed bistability of evolutionary attractors. The evolutionary consequence of the threshold trait \({z}_{\mathrm{T}}\), that was convergence stable, decreased with an increase in survivorship of resident \({S}_{\mathrm{R}}\). Ecological dynamics were always stable in this case.

[2] The threshold value \({z}_{\mathrm{T}}\) at the stable attractor was independent of the intermediate size \({\mu }_{0}\) (Fig. 3c).

[3] In eco-evolutionary dynamics, the threshold value \({z}_{\mathrm{T}}\) hardly changed with an increase in intermediate size \({\mu }_{0}\) (Fig. 8b). Slight changes in trait values were observed during the environmental changes (\({T}_{1}<t<{T}_{2}\)), and they were caused by the changes in the population structure of the residents. The changes in the intermediate size did not influence the ecological dynamics (Fig. 8c). These outcomes of the responses of eco-evolutionary dynamics to climate change (formalized as Eq. (12)) are greatly different from those of the absolute assessment models. This suggests that the responses to climate changes critically depend on the manner in which life-history choice is performed by the organism.

Discussion

Population respond to climate changes both ecologically and/or evolutionarily (Bradshaw and Holzapfel 2006; Charmantier et al. 2008; Chevin et al. 2010). Traits change as evolution operates on fitness differences in the environment and also by phenotypic plasticity depending on the conditions of the environment. These responses may enable the population to survive under environmental changes (Ashander et al. 2019). While some recent studies in chaos and dynamics in eco-evolutionary models under changing environment have focused on the phase transition (Gilpin and Feldman 2017) and predictability (Rego-Costa et al. 2018), the relationship between the rate of climate change and that of evolution was underexplored. We here explored the possibility that the irreversible regime shift was caused depending on the balance of these two.

Eco-evolutionary responses to change may result in long transient, or even irreversible changes to population dynamics

In this study, we have investigated the effects of climate change on population dynamics. In cases where individuals made decisions based on their own status, such as body size (absolute assessment), the fraction of migrants was changed because of climate change (Fig. 1). As the intermediate size \({\mu }_{0}\) increased, the number of juveniles becoming larger surpassed the threshold, resulting in a greater fraction of juveniles employing resident tactics. Climate change may also affect the qualitative behavior of the population dynamics (Fig. 1b and c).

These changes in the ecological dynamics could be reverted to those before climate change occurred, when we consider the evolution of threshold traits (Fig. 5), given that the rate of the evolution is fast enough when compared with the environmental change, and the threshold trait that was convergence stable, linearly increased with increasing intermediate size, \({\mu }_{0}\) (Fig. 3a). During climate changes, the threshold did not evolve linearly, and hence, the difference between the evolutionary attractor and the threshold traits adopted in the population increased (Fig. 5a). This caused considerable change in the ecological dynamics, in terms of both the fraction of migrants and the stability of the system (Fig. 5b). However, if there was only one evolutionary equilibrium, the evolutionary dynamics finally converged to the equilibrium (Fig. 5a), and this resulted in the recovery of the ecological dynamics, to the same qualitative behavior as that prior to climate changes (Fig. 5b).

On the other hand, when there were multiple equilibria in the threshold trait (Fig. 3b), increasing the intermediate size \({\mu }_{0}\) might cause the transition of trait values from one evolutionary equilibrium to another (Fig. 6a). This resulted in irreversible changes in the dynamic behavior of the system (Fig. 6b). This evolutionary regime shift occurred when the rate of climate change was sufficiently large in comparison with the rate of the evolution (Figs. 4a and 7). In order to demonstrate potential behavior of the model, in this paper we adopted arbitrary values as the rate of the evolution and climate change in our simulation, and based on the analyses, we concluded that the balance between these parameters is important in a factor realizing the evolutionary regime shift (Fig. 7). Phillis et al. (2016) evaluated the rate of the evolution of threshold trait of smoltification in the system of rainbow trout. Quantitative study of measuring the rate of evolution is important because it will allow us to know whether or not the evolutionary regime shift likely to occur under natural situations. These results highlight the importance of reducing the rate of climate change. Bradshaw and Holzaphel (2006) claimed that organisms respond evolutionarily to climate changes and they are able to adapt to the environment; however, the survival of organisms will be in danger if we do not reduce the rate of climate change. Thomas et al. (2004) also called attention to the importance of reducing the rate of climate change to avoid a crisis of extinction for organisms. In this study, we showed that the mismatch of a threshold trait in a given environment, caused by climate change, might consequently cause a regime shift of the ecological dynamics (Fig. 1). This would result in phenotypically different population dynamics. Since we cannot modify the rate of the evolution for natural populations without significant interventions, reducing the rate of climate change may be one potential policy to prevent irreversible changes to the natural populations.

We examined the distance between the upper attractor and repeller when evolutionary bistability occurred (Fig. S4 in Online Resource 1) and found that it increased with a decrease in survivorship of residents \({S}_{\mathrm{R}}\) and increased with an increase in the survivorship of migrants \({S}_{\mathrm{M}}\). In other words, when the survivorship of residents was small and that of the migrants was large, the evolutionary regime shift would be unlikely.

Relative versus absolute assessment rules

In the cases where individuals make the decision based on the estimated social ranks in the population (relative assessment), climate change had no effect on the evolutionary dynamics (Fig. 8a). The threshold trait hardly changed when we considered the evolution of the threshold with the occurrence of climate change (Fig. 8b), and there were no situations in which evolutionary bistability occurred. The ecological dynamics of the relative assessment model were less severely affected by climate change than those of the absolute assessment model. Relative assessment can be regarded as a form of phenotypic plasticity. In general, phenotypic plasticity is known to help the maintenance of a population in the face of climate change (Chevin et al. 2010). Charmantier et al. (2008) showed that the phenotypic plasticity might allow the population to cope with rapid climate change. Hence, it is important to know whether the salmonids make the decision by relative assessment or absolute assessment. In natural situations, organisms may adopt a mixture of absolute and relative assessments. To predict the outcomes of climate change in relation to eco-evolutionary dynamics, it is important to know how relative elements influence life-history decision making, and this should be the focus of more research in the future.

In this study, we assumed that climate change affects the body size only, because the size directly affects life-history choices. There are, however, documentations that climate changes affect organisms from multiple aspects simultaneously (Walther et al. 2002), such as body size (Crozier et al. 2008; Daufresne et al. 2009), habitat (Perry et al. 2005; Warren et al. 2001), phenology (Przybylo et al. 2000), survivorship (Ozgul et al. 2010; Satterthwaite et al. 2010), and species interactions (Tylianakis et al. 2008). These will change the eco-evolutionary dynamics, and it is necessary to consider which changes are significant when we study the effects of climate change on particular organisms. Furthermore, we here focused only on the life-history of males, because most female exhibits migrant tactic in masu salmon. There are however some species in which females also show alternative life-history tactics (e.g., brown trout: Jonsson et al. 2001; masu salmon: Morita & Nagasawa 2010). Since female directly affect the population dynamics, considering the decision making in both sexes enables us to explore the dynamical change and population persistence of whole population as consequences of environmental change.

Conclusion

The present study demonstrated the importance of understanding the mechanisms of life history choices in predicting the impacts of climate change. The way in which life history choice mechanisms respond to change also is important. Bradshaw and Holzapfel (2008) conjectured that both phenotypic plasticity and genetic evolution mitigate the harmful effects of climate change on organisms. However, the present study showed that when the rate of climate change is too fast, traits may shift to different evolutionary attractors, and these shifts may be irreversible. As a result, the ecological dynamics can become qualitatively different from the initial ones.