ESM 1 Description of model with ODD protocol

We follow the ODD (Overview, Design concepts, Details) protocol for describing individual-based and agent-based models (Grimm et al. 2006; 2010). The model we use here is an extension of the model of population dynamics of marble trout by Vincenzi et al. (2008a), accounting for variability in individual growth rates among individuals living in the same population.

Purpose

The purpose of the model is to explore how the distribution of individual body growth rates in a fish population may evolve in three different environments, namely constant with stochastic birth and death processes (CON), variable (i.e., with environmental stochasticity, VAR) and with occurrence of catastrophic events (i.e., floods) causing sudden and massive mortality (CAT). The fish population that we study has density-dependent survival and body growth and trade-off between growth and mortality (all described below).

Entities, state variables and scales

We model the population dynamics of female marble trout living in Zakojska (Slovenia). Individuals are characterized by the state variables: unique identifier, somatic growth parameter , age, body length, and sexual maturity. We refer to individuals in their first year of life (from zero to one winter) as juveniles, while referring to fish ≥ age-1 as immature and, after sexual maturity, adults (Fig. 1in the main text). We do not explicitly model metabolism (individual-oriented model, sensu Grimm and Railsback 2005).

The environment is a stretch of a stream, hosting a population of marble trout. We simulate evolution in three different environmental scenarios, constant (CON) and variable regimes (VAR), and the variable regime where catastrophic events in the form of flooding occur (CAT). Individual survival depends on the state of the environment. We assume the surface area of the stream to be constant over time and equal to 900 m2. There is no spatial variation within the stretch. We assume flood events not to modify stream environment. The population dynamics are simulated for a variable duration, as explained in Simulation scenarios.

Process overview and scheduling

The time step in the model is one year. We compute population size and the distribution of individual growth rates in the fall, before spawning. We consider populations with fewer than 5 females > age-0 extinct, and exclude such replicates from further analysis. At each time step, Bernoulli trials define if an individual survive or dies. Eggs are spawned in the fall and age-0 emerge in the following year in late spring-early summer. We assume flood events to impact the population after the spawning period.

Design concepts

Basic principles: The theoretical studies on effects of density-dependence and environmental stochasticity on life-histories and population dynamics usually assume infinite population sizes, so that effects of randomness can be ignored. However, since we study a small population, this randomness may be important. We therefore need an individual-oriented approach. Further, most theories and applications on environmental stochasticity focus on non-catastrophic events. Our study should therefore be of interest both for management of small fish populations in streams and for the general understanding of effects of environmental variation on population dynamics and life-histories.

Adaptation:The parameter is the single element in the strategy vector (sensuHuse et al. 2002) of the individual. It is heritable, with a small probability of mutation. It influences mortality risk and body growth, and thus age at sexual maturation and the reproductive output of the organism.

Due to natural selection and random drift, the distribution of  changes over time, with ecological consequences for population size, age-structure and distribution of phenotypes associated with .

Emergence, Learning, and Collectives:Population dynamics emerge from the life cycle of the individuals, which is represented by empirical rules describingsurvival at different life stages as probabilities. The body growth of fish is predicted by the Submodel “Body growth” from population density and the individual growth trait of the individual.

There is no individual learning in the model. The population size and the age-structure determine (together with the current state of the stream stretch in VAR and CAT simulations) the density-dependent aspects of growth and survival of the individual. Since we do not model behavior, these are passive processes.

Objectives, Prediction and Sensing: Therefore the model does not contain elements of prediction and sensing. However, since the distribution of  in future generationsis the result of reproduction in past generations, the implicit objective for the individual is to achieve the maximum reproductive success. Both survival and fecundity are functions of.

Interaction: We implicitly model three types of interactions: juvenile mortality increases with increasing juvenile densities (DU), mortality of fish ≥ age-1 increases with increasing densities of fish ≥ age-1 (DA), body growth decreases with increasing densities of fish ≥ age1 experienced during the first year of life. Densities (ind. m-2) are computed at each time-step by dividing twice the number of individuals (a 1:1 sex ratio) by stream surface area.

Stochasticity: In CON, we include only demographic stochasticity in the model. In VAR, we model environmental stochasticity by randomly drawing parameters values at each time step of the simulation from the uniform distribution [0.95 P,1.05 P], where P denotes the value of a parameter. We chose the uniform distribution for simplicity and the interval [0.95 P,1.05 P] to clearly separate the effect of environmental variability from the effect of catastrophic events (but see Sensitivity analysis). In CAT, we model catastrophes (i.e., flood events) by randomly drawing at each time step a rainfall value for October R (in mm) from the lognormal distribution fitted to rainfall data recorded in the study area from 1961 to 2007 (= 300 mm, sd= 328, both parameters are on the arithmetic scale) (Vincenzi et al. 2008a). Parameters values are reported in Table 1 of the main text.

Observation: For model analysis and testing, we record population-level variables at each time-step, e.g.distribution of individual growth rates in the population over time, population size over time, and age-distribution of the population.

Initialization

We start each simulation with 350 individuals randomly aged between 0 and 6, each individual having a randomly assigned growth parameter  from a uniform distribution bounded between 0 and 2 (see Sensitivity of model results).

Input

In CAT, at each time step we randomly draw an October rainfall value R from the empirical log-normal distribution of rainfall. If rainfall is greater than 400 mm, a flood event occurs and survival of fish is modified as described in the following section (Submodel Survival).

Submodels

Body growth: In marble trout, density conditions experienced during the first year of life have a lasting influence on the body growth of fish (Vincenzi et al. 2008c;see Harrison et al. 2011 for a review of carry-over effects and their effect on individual fitness).We model length-at-age as density-dependent von Bertalanffy growth. Thus, for a an individual of age a, length growth rate  and density of marble trout ≥ age-1 during its first year of life (density as underyearling, DU) and strength of density-dependence , length-at-age is:

(1)

We estimated model parameters by non-linear regression using data in Vincenzi et al. (2008c). Density at later years and maturation seem to have only a minor role in determining growth trajectories and thus for simplicity they were not considered in the model of body growth (Vincenzi et al. 2008c). The growth parameter  is a multiplier of growth rate k and equal to 1 for the population of marble trout (“base case”) upon which the model parameters were tuned (Fig. 2 in the main text). We assume that individuals living sympatrically may have different individual growth rates ranging from 0 to 2.

Reproduction: A tight relationship linkingfemale length to sexual maturity and egg production iscommonly observed in salmonids (e.g., Hendry and Stearns 2003).Vincenzi et al. (2008b) found from farm experiments that the number of eggs produced by marble trout E was related to body length L of female by the linear relationship:

E = + L (2)

where  and is length at sexual maturation. We estimated the form and parameters of Eq. 2 using data from experiments in a fish farm (Vincenzi et al. 2007). Based on Eq. (2), length at maturity is invariant to differences in growth rates,although the corresponding age at maturity can change substantially as a result of density-dependence in growth( Mangel and Abrahams 2001; Lorenzen 2005). A number of studies on salmonid species have reported a general correlation between rapid growth and early sexual maturation (e.g., Hutchings and Jones 1998; Utrilla and Lobon-Cervia 1999), with slowly growing fish delaying sexual maturation in order to reach the minimum size required for gonad development.

We assume that fish spawn annually from maturity until death. Approximately half of mature females reproduce successfully each year in the wild (Meldgaard et. 2007; AJ Crivelli, unpublished data) and we use Bernoulli trials (p= 0.5) to determine if eggs produced by the individual female are viable and thus included in the egg pool. Variation in growth patterns, either from differences in individual growth rates or population density, affects the reproductive value of individuals.

Transmission of growth parameter : In order to facilitate the interpretation of results, we assume that the newborn perfectly inherited the growth parameter  from the parent (see Sensitivity analysis), but we introduce “mutation” rates to allow for small variation in individual growth rates over the simulation time. Given p the growth parameter of the parent and o the growth parameter of the offspring:

(3)

where u is randomly draw from the uniform distribution over the interval [-0.2,0.2]. We performed simulations with reduced heritability of to assess the robustness of results (see Sensitivity analysis).

Survival: We assume that egg to age-0 survival SE is density-independent (Rubin et al., unpublished manuscript) and constant across , while juvenile (from age-0 to one winter) survival Sj is density-dependent:

(4)

where M is the annual mortality where densities approach zero, Dj is the density of juveniles and j describes the strength of density-dependence for juvenile survival. We estimated parameters of the model of first-year survival from the results of Vincenzi et al (2007).

The growth parameter determines the growth-mortality trade-off for the individual. The classic trade-offbetween growth and mortality arises in the case of foraging under predation risk, in which animals may grow at optimal, rather than maximal growth rates, due to predation risk (e.g., Arendt1997). However,several recent studies indicated that varioustrade-offs may limit the benefit of rapid growth (Arendt1997;Mangel and Stamps 2001; Post and Parkinson 2001; Munch and Conover2003; Biroet al. 2005; Stoks et al. 2005) andmany functionsor traits may be reduced or compromised or when individual growth rates are maximized (Roff 2002).For simplicity and in absence of empirically validated models for marble trout, we assume a linear effect of oninstantaneous mortality rate.

Vincenzi et al. (2008a) did not find any evidence for size-selective mortality so that Sj is independent of size. Annualsurvival SIfor sexually immature trout is therefore:

(5)

where DA is the density of fish ≥ age-1 and A characterizes the strength of densitydependence for survival. We estimated parameters of Eq. 5 from mark-recapture data (Vincenzi et al. 2008a; S Vincenzi, unpublished manuscript). As is commonly observed in stream-dwelling salmonids, the strength of density-dependent survival is greater for juveniles than for adults (jA).

For sexually mature trout, we modify Eq. 5 to account for spawning mortality. :

(6)

where  characterizes the cost of reproduction. We obtained the value of from experiments performed in the fish farm (S Vincenzi and AJ Crivelli, unpublished data).

If  is the age of sexual maturity, survival from egg to age a in a constant environment is:

(7)

We model survival probability SF in case of flood event, common for juveniles, immature and adult fish, as follows:

(8)

where MF(R) is the mortality induced by the flood event, depending on the rainfall R. In case of R>400 mm, we assume mortality to be independent of growth parameter density and age. We developed the relationship between flood event and mortality MF by using empirical data and expert knowledge, with October rainfall recorded in the study area used as a proxy for flood event. For R> 400 mm:

(9)

were R is October rainfall (mm). Rainfall greater than 1000 mm (for whichMF(R) = Mmax) has a recurrence time of approximately 20 years.

Simulation scenarios

Number of replicates: We ran the model for 500replicates for each of the three scenarios CON, VAR and CAT.

End of simulations: If tdenotes the distribution of growth parameters at time twe stop the simulation when is quasi-stable (qs), characterized by the following. For two time steps 25 years apart (t, t+25):

  • both minimum and maximum value of (t, t+25) coincide;
  • the p-value of theKolmogorov-Smirnov (KS) test for differences between t andt+25is ≥ 0.7. The p-value was arbitrarily chosen after preliminary exploration of model simulations (but see Sensitivity of model results).

Metrics:Since changes in the distribution of a trait can occur as a result ofdirectional, stabilizing or disruptive selection (e.g., Bradshaw and Holzapfel 2006), drift, plastic response to environmental change (e.g., density-dependent body growth, Via et al. 1995) changing selection pressure (e.g., Lande 1976; Gienapp et al. 2008) and the distribution of traits in newborns, we report a variety of metrics in order to assess the population-level response in growth parameter . Preliminary explorations of replicates guided the choice of appropriate and informative metrics (see Fig.3 in the main text).

  1. We fit both a bi-modal and uni-modal Gaussian mixture model to qsand use AIC to assess which of the two models better explained the data. A Uni-modal distribution of qs may be seen as directional or stabilizing selection depending on the location of the peak, while the bi-modal distribution of qs may be seen as disruptive selection.
  2. We estimate the location of peak (or peaks, depending on the selection of uni- or bi-modality) in qs, as identified by the Gaussian mixture model.
  3. We record max(qs) and min(qs)
  4. We compute the range of growth parameters in the population at tqs, as max(qs) – min(qs).
  5. Finally, we record the time to reach a quasi-stable distribution tqs.

Sensitivity of model results: We assessed the robustness of results to alternative assumptions regarding the shape of the initial distribution of  (uniform, beta and lognormal) and the number of individuals (ranging from 300 to 500) in the first year of the simulation. We performed simulations with different p-values of the KS test ranging from 0.5 to 0.8 to test if alternative p-values substantially altered the time tqs to reach a quasi-stable distribution. We widened the interval of uniform distribution to [0.9P, 1.1P] of parameters estimates and tested whether it influenced the simulation results of scenario VAR. We tested the sensitivity of model results to variation in the strength of density-dependent body growth () and flood events () (Table 1 in the main text). We used logistic regression for the categorical variable qs (bi-modal  and uni-modal = 0, see Fig. 3 in the main text). We used ordinary least-square regression when analyzing the sensitivity of range of growth rates at tqs, max(qs) and min(qs). In all regressions, we started with value of and  and their interaction predictors and then proceeded with stepwise selection. Finally, we carried out simulations with reduced heritability of  (Electronic Supplementary Material 3)

References

Arendt JA (1997) Adaptive intrinsic growth rates: an integration across taxa j. Q Rev Biol72:149-177.

Biro PA, Post JR, Abrahams MV (2005) Ontogeny of energy allocation reveals selective pressure promoting risk-taking behaviour in young fish cohorts. P R Soc B 272:1443-1448.

Bradshaw WE, Holzapfel CM (2006) Evolutionary response to rapid climate change. Science 312:1477-1488.

Gienapp P, Teplitsky C, Alho J, Mills J, Merilä J (2008) Climate change and evolution: disentangling environmental and genetic responses. Mol Ecol 17:167-178.

Grimm V, Railsback S (2005) Individual-based Modeling and Ecology. Princeton University Press, Princeton.

Grimm V, Berger U, Bastiansen F, Eliassen S, Ginot V, Giske J, Goss-Custard J, Grand T, Heinz SK, Huse G,Huth A, Jepsen JU, Jørgensen C, Mooij WM, Müller B, Pe’er G, Piou C, Railsback SF, Robbins AM, Robbins MM, Rossmanith E, Rüger N, Strand E, Souissi S, Stillmann R, Vabø R, Visser U & DeAngelisDL (2006) A standard protocol for describing individual-based and agent-based models. EcolModel 198:115-126.

Grimm V, Berger U, DeAngelis D, Polhill G, Giske J, Railsback S (2010) The ODD protocol: a review and first update. EcolModel 221:2760-2768.

Harrison Xa, Blount JD, Inger R, Norris DR, Bearhop S (2011) Carry-over effects as drivers of fitness differences in animals. J Anim Ecol. 80: 4-18.

Hendry AP, Stearns SC (2003) Evolution Illuminated: Salmon and Their Relatives. Oxford University Press.

Hutchings J, Jones M (1998) Life history variation and growth rate thresholds for maturity in Atlantic salmon, Salmo salar. Can J Fish Aquat Sci 55:22-47.

Lande R (1976) Natural selection and random genetic drift in phenotypic evolution. Evolution 30:314-334.

Lorenzen K (2005) Population dynamics and potential of fisheries stock enhancement: practical theory for assessment and policy analysis. Philos T Roy Soc B 360:171-89.

Mangel M, Abrahams M (2001) Age and longevity in fish, with consideration of the ferox trout. Experimental Gerontology 36:765-790

Mangel M, Stamps J (2001) Trade-offs between growth and mortality and the maintenance of individual variation in growth. Evol Ecol Res 3:583-593.

Meldgaard T, Crivelli AJ, Jesensek D, Poizat G, Rubin J, Berrebi P (2007) Hybridization mechanisms between the endangered marble trout (Salmo marmoratus) and the brown trout (Salmo trutta) as revealed by in-stream experiments. Biol Conserv 136:602-611.

Munch SB, Conover DO (2003) Rapid growth results in increased susceptibility to predation in Menidia menidia. Evolution 57:2119-2127.

Post JR, Parkinson Ea (2001) Energy Allocation Strategy in Young Fish: Allometry and Survival. Ecology 82:1040-1051.

Roff D (2002) Life History Evolution. Sinauer Associates, Sunderland, MA.

Stoks R, De Block M, McPeek MA (2005) Alternative growth and energy storage responses to mortality threats in damselflies. Ecol Lett 8:1307-1316.

Utrilla CG, Lobon-Cervia J(1999) Life-history patterns in asouthern population of Atlantic salmon. J. Fish Biol.55: 68–83.

Via S, Gomulkiewicz R, DeJong G, Scheiner SM, Schlichting CD,VanTiendern PH (1995) Adaptivephenotypic plasticity: consensus and controversy. Trends EcolEvol 10: 212–217.

Vincenzi S, Crivelli AJ, Jesensek D, Rubin J-F, De Leo GA (2007) Early survival of marble trout Salmo marmoratus: Evidence for density dependence? Ecol Freshw Fish 16:116-123.

Vincenzi S, Crivelli AJ, Jesensek D, De Leo GA (2008a) The role of density-dependent individual growth in the persistence of freshwater salmonid populations. Oecologia 156:523-534.

Vincenzi S, Crivelli AJ, Jesensek D, Rubin J, Poizat G, De Leo GA (2008b) Potential factors controlling the population viability of newly introduced endangered marble trout populations. Biol Conserv 141:198-210.

Vincenzi S, Crivelli AJ, Jesensek D, De Leo GA (2008c) Total population density during the first year of life as a major determinant of lifetime body-length trajectory in marble trout. EcolFreshw Fish 17:515-519.