Abstract
We introduce an Interaction and Tradeoffbased EcoEvolutionary Model (ITEEM), in which species are competing in a wellmixed system, and their evolution in interaction trait space is subject to a lifehistory tradeoff between replication rate and competitive ability. We demonstrate that the shape of the tradeoff has a fundamental impact on ecoevolutionary dynamics, as it imposes four phases of diversity, including a sharp phase transition. Despite its minimalism, ITEEM produces a remarkable range of patterns of ecoevolutionary dynamics that are observed in experimental and natural systems. Most notably we find selforganization towards structured communities with high and sustained diversity, in which competing species form interaction cycles similar to rockpaperscissors games.
https://doi.org/10.7554/eLife.36273.001eLife digest
A patch of rain forest, a coral reef, a pond, and the microbes in our guts are all examples of biological communities. More generally, a community is a group of organisms that live together at the same place and time. Many communities are composed of a large number of different species, and this diversity is maintained for long times.
Although diversity is a key feature of biological communities, the mechanisms that generate and maintain diversity are not well understood. Research had hinted at links between diversity and the tradeoffs that species are subject to. For instance, there is a tradeoff between competitiveness and reproduction: if there are limited resources in the environment a species may either produce many offspring that are not very competitive, or fewer, more competitive offspring.
Farahpour et al. have now simulated the development of communities of organisms that reproduce, compete, and die in a uniform environment. Crucially, these computational simulations introduced a tradeoff between competitive ability and reproduction.
The simulations show that the form of tradeoff has a fundamental impact on diversity: moderate tradeoffs favor diversity, whereas extreme tradeoffs suppress diversity. The simulations also revealed mechanisms that underlie how diversity is generated. In particular, cyclic relationships emerge where one species dominates another but is also dominated by a third, similar to the rockpaperscissors game.
Since Farahpour et al. used a barebone model with only a few essential features the results could apply to a larger class of communitylike systems whose evolution is driven by competition. This includes economic and social systems as well as biological communities.
https://doi.org/10.7554/eLife.36273.002Introduction
We observe an immense diversity in natural communities (Hutchinson, 1961; Tilman, 1982; Huston, 1994), but also in controlled experiments (Maharjan et al., 2006; Gresham et al., 2008; Kinnersley et al., 2009; Herron and Doebeli, 2013; Kvitek and Sherlock, 2013), where many species continuously compete, diversify and adapt via ecoevolutionary dynamics (Darwin, 1859; Cody and Diamond, 1975). However, the basic theoretical models (Volterra, 1928; Tilman, 1982) predict that both ecological and evolutionary dynamics tend to decrease the number of coexisting species by competitive exclusion or selection of the fittest. This apparent contradiction between observations and theory gives the stunning biodiversity in communities the air of a paradox (Hutchinson, 1961; Sommer and Worm, 2002) and hence has begotten a long, ongoing debate on the mechanisms underlying emergence and stability of diversity in communities of competitive organisms (Hutchinson, 1959; Huston, 1994; Chesson, 2000; Sommer and Worm, 2002; Doebeli and Ispolatov, 2010).
To identify candidate mechanisms that could resolve the problem of generation and maintenance of diversity, the basic theoretical ecological and evolutionary models have been extended by numerous features (Chesson, 2000; Chave et al., 2002), including spatial structure (Mitarai et al., 2012; Villa Martín et al., 2016; Vandermeer and Yitbarek, 2012), spatial and temporal heterogeneity (Caswell and Cohen, 1991; Fukami and Nakajima, 2011; Hanski and Mononen, 2011; Kremer and Klausmeier, 2013), tailored interaction network topologies (Melián et al., 2009; Mougi and Kondoh, 2012; Kärenlampi, 2014; Laird and Schamp, 2015; Coyte et al., 2015; Grilli et al., 2017), predefined niche width (Scheffer and van Nes, 2006; Doebeli, 1996), adjusted mutationselection rate (Johnson, 1999; Desai and Fisher, 2007), and lifehistory tradeoffs (Rees, 1993; Bonsall et al., 2004; de Mazancourt and Dieckmann, 2004; Gudelj et al., 2007; Ferenci, 2016; Posfai et al., 2017). However, it is still unclear which features are essential to explain biodiversity. For instance, diversity is also observed under stable and homogeneous conditions (Gresham et al., 2008; Kinnersley et al., 2009; Maharjan et al., 2012; Herron and Doebeli, 2013; Kvitek and Sherlock, 2013).
So far, models of ecoevolutionary dynamics have been developed in three major categories: models in genotype space, like population genetics (Ewens, 2012) and quasispecies models (Nowak, 2006); models in phenotype space, like adaptive dynamics (Doebeli, 2011) and webworld models (Drossel et al., 2001); and models in interaction space, like LotkaVolterra models (Coyte et al., 2015; Ginzburg et al., 1988) and evolving networks (Mathiesen et al., 2011; Allesina and Levine, 2011). Each of these categories has strengths and limitations and emphasizes particular aspects. However, in nature these aspects are entangled by ecoevolutionary feedbacks that link genotype, phenotype, and interaction levels (Post and Palkovacs, 2009; Schoener, 2011; Ferriere and Legendre, 2013; Weber et al., 2017). In a closed system of evolving organisms mutations, that is, evolutionary changes at the genetic level (Figure 1a), can cause phenotypic variations if they are mapped to novel phenotypic traits in phenotype space (Figure 1b)(Soyer, 2012). These variations have ecological impact only if they affect biotic or abiotic interactions of species (Figure 1c); otherwise they are ecologically neutral. The resulting adaptive variations in the interaction network change the species composition through population dynamics. Finally, frequencydependence occasionally selects strategies that adapt species to their new environment (Schoener, 2011; MoyaLaraño et al., 2014; Hendry, 2016; Weber et al., 2017).
Thus, we have a link from interactions to ecoevolutionary dynamics, suggesting that we do not need to follow all evolutionary changes at the genetic or phenotypic level if we are interested in macroecoevolutionary dynamics, but only those changes that affect interactions. In this picture, evolution can be considered as an exploration of interaction space, and modeling at this level can help us to study how complex competitive interaction networks evolve and shape diversity. This neglect of genetic and phenotypic details in interactionbased models (Ginzburg et al., 1988; Solé, 2002; Tokita and Yasutomi, 2003; Shtilerman et al., 2015) equals a coarsegraining of the ecoevolutionary system (Figure 1). This coarsegraining not only reduces complexity but it should also make the approach applicable to a broader class of biological systems.
Interactionbased evolutionary models have received some attention in the past (Ginzburg et al., 1988; Solé, 2002) but then were almost forgotten, despite remarkable results. We think that these works have pointed to a possible solution of a hard problem: The complexity of evolving ecosystems is immense, and it is therefore difficult to find a representation suitable for the development of a statistical mechanics that enables qualitative and quantitative analysis (Weber et al., 2017). Modeling at the level of interaction traits, rather than modeling of detailed descriptions of genotypes or phenotypes, coarsegrains these complex systems in a natural way so that this approach may be helpful for developing a biologically meaningful statistical mechanics.
The first ecoevolutionary interactionbased model was introduced by Ginzburg et al. (1988) based on Lotka–Volterra dynamics for competitive communities. Instead of adding species characterized by random coefficients, taken out of some arbitrary species pool, they made the assumption that a new mutant should be ecologically similar to its parent, which means that phenotypic variations that are not ecologically neutral generate mutants that interact with other species similar to their parents (Figure 1). Thus, speciation events were simulated as ecologically continuous mutations in the strength of competitive interactions. This model, although conceptually progressive, was not able to produce a large stable diversity, possibly because diversity requires components not included in this model. Therefore subsequent interactionbased models supplemented it with ad hoc features to specifically increase diversity, such as special types of mutations (Tokita and Yasutomi, 2003), addition of mutual interactions (Tokita and Yasutomi, 2003; Yoshida, 2003), enforcement of partially connected interaction graphs (Kärenlampi, 2014), or imposed parentoffspring niche separation (Shtilerman et al., 2015). While these models generated, as expected, higher diversity than the original Ginzburg model, they could not reproduce key characteristics of real systems, for example emergence of large and stable diversity, diversification to separate species and mass extinctions. Of course, the use of ad hoc features that deliberately increase diversity also cannot explain why diversity emerges.
An essential component missing in the previous interactionbased models had been a constraint on strategy adoption. In real systems such constraints prevent the emergence of Darwinian Demons, that is, species that develop in the absence of any restriction and act as a sink in the network of population flow. Among all investigated features responsible for diversity, mentioned above, lifehistory tradeoffs that regulate energy investment in different lifehistory strategies are fundamentally imposed by physical laws such as energy conservation or other thermodynamic constraints, and thus present in any natural system (Stearns, 1989; Gudelj et al., 2007; Del Giudice et al., 2015). These physical laws constrain evolutionary trajectories in trait space of evolving organisms and determine plausible evolutionary paths (Fraebel et al., 2017; Ng'oma et al., 2017), i.e. combinations of strategies adopted or abandoned over time. Roles of tradeoffs for emergence and stabilization of diversity have been investigated in previous ecoevolutionary studies (Posfai et al., 2017; Rees, 1993; Bonsall et al., 2004; de Mazancourt and Dieckmann, 2004; Ferenci, 2016; Gudelj et al., 2007) and experiments (Stearns, 1989; Kneitel and Chase, 2004; Agrawal et al., 2010; Maharjan et al., 2013; Ferenci, 2016). It has been shown, for example, that if metabolic tradeoffs are considered, even at equilibrium and in homogeneous environments, stable coexistence of species becomes possible (Gudelj et al., 2007; Beardmore et al., 2011; Maharjan and Ferenci, 2016).
Here, we introduce a new, minimalist model, the Interaction and Tradeoffbased EcoEvolutionary Model (ITEEM), with simple and intuitive ecoevolutionary dynamics at the interaction level that considers a lifehistory tradeoff between interaction traits and replication rate, that means, better competitors replicate less (Jakobsson and Eriksson, 2003; Bonsall et al., 2004). To our knowledge, ITEEM is the first model which joins these two elements, the interactionspace description with a lifehistory tradeoff, that we deem crucial for an understanding of ecoevolutionary dynamics. We use ITEEM to study development of communities of organisms that diversify from one ancestor by gradual changes in their interaction traits and compete under LotkaVolterra dynamics in wellmixed, closed system.
We show that ITEEM dynamics, without any ad hoc assumption, not only generates large and complex biodiversity over long times (Herron and Doebeli, 2013; Kvitek and Sherlock, 2013) but also closely resembles other observed ecoevolutionary dynamics, such as sympatric speciation (Tilmon, 2008; Bolnick and Fitzpatrick, 2007; Herron and Doebeli, 2013), emergence of two or more levels of differentiation similar to phylogenetic structures (Barraclough et al., 2003), occasional collapses of diversity and mass extinctions (Rankin and LópezSepulcre, 2005; Solé, 2002), and emergence of cycles in interaction networks that facilitate species diversification and coexistence (Buss and Jackson, 1979; Hibbing et al., 2010; Maynard et al., 2017). Interestingly, the model shows a unimodal (‘humpback’) course of diversity as function of tradeoff, with a critical tradeoff at which biodiversity undergoes a phase transition, a behavior observed in nature (Kassen et al., 2000; Smith, 2007; Vallina et al., 2014; Nathan et al., 2016). By changing the shape of tradeoff and comparing the results with a notradeoff model, we show that diversity is a natural outcome of competition if interacting species evolve under physical constraints that restrict energy allocation to different strategies. The natural emergence of diversity from a barebone ecoevolutionary model suggests that a unified treatment of ecology and evolution under physical constraints dissolves the apparent paradox of stable diversity.
Model
ITEEM is an individualbased model (Black and McKane, 2012; DeAngelis and Grimm, 2014) with simple intuitive updating rules for population and evolutionary dynamics. A simulated system in ITEEM has ${N}_{s}$ sites of undefined spatial arrangement (no neighborhood), each providing permanently a pool of resources that is sufficient for the metabolism of one organism. The community is wellmixed, which means that the probability for an encounter is the same for all pairs of individuals, and that the probability of an individual to enter a site (i.e. to access resources) is the same for all individuals and sites.
We start an ecoevolutionary simulation with individuals of a single strain occupying a fraction of the ${N}_{s}$ sites, and then carry out long simulations for millions of generations. Note that in the following, to facilitate discourse, we use the term strain for a group of individuals with identical traits, whereas the term species denotes a monophyletic cluster of strains with some intraspecific diversity (for a discussion on application of these terms in this study see Appendix 1, Species and strains). Over time $t$, measured in generations, the number of individuals, ${N}_{ind}\left(t\right)$, number of strains, ${N}_{st}\left(t\right)$, and number of species, ${N}_{sp}\left(t\right)$, change by ecological (birth, death, competition) and evolutionary dynamics (mutation, extinction, diversification).
Every generation or time step consists of ${N}_{s}$ sequential replication trials of randomly selected individuals, followed at the end by a single death step. In the death step all individuals that have reached their lifespan at that generation will vanish. Lifespans of individuals are drawn at their births from a Poisson distribution with overall fixed mean lifespan $\lambda $. This is equivalent to an identical per capita death rate for all strains. For comparison, simulations with no attributed lifespan ($\lambda =\infty $) were carried out, too; in this case the only cause of death is defeat in a competitive encounter.
At each replication trial, a randomly selected individual of a strain $\alpha $ can replicate with probability ${r}_{\alpha}$. Age of individuals plays no role in their reproduction and thus a newborn individual can be selected and replicate with the same probability as adult individuals. With a fixed probability $\mu $ the offspring mutates to a new strain ${\alpha}^{\prime}$. Then, the newborn individual is assigned to a randomly selected site. If the site is empty, the new individual will occupy it. If the site is already occupied, the new individual competes with the current holder in a lifeordeath struggle. In that case, the surviving individual is determined probabilistically by the ‘interaction’ ${I}_{\alpha \beta}$, defined for each pair of strains $\alpha $, $\beta $. ${I}_{\alpha \beta}$ is the survival probability of an $\alpha $ individual in a competitive encounter with a $\beta $ individual, with ${I}_{\alpha \beta}\in [0,1]$ and ${I}_{\alpha \beta}+{I}_{\beta \alpha}=1$ (Grilli et al., 2017). All interactions $I}_{\alpha \beta$ form an interaction matrix $I\left(t\right)$ that encodes the outcomes of all possible competitive encounters in this probabilistic sense. Row $\alpha $ of $I$ defines the ‘interaction trait’ ${\mathbf{T}}_{\alpha}=\left({I}_{\alpha 1},{I}_{\alpha 2},\dots ,{I}_{\alpha {N}_{st}\left(t\right)}\right)$ of strain $\alpha $, with ${N}_{st}\left(t\right)$ the number of strains at time $t$.
If strain $\alpha $ goes extinct, its interaction elements must be removed, i.e. the $\alpha $th row and column of $I$ are deleted. Conversely, if a mutation of $\alpha $ generates a new strain ${\alpha}^{\prime}$, its trait vector is obtained by adding a small random variation to the parent trait, that is $\mathbf{T}}_{{\alpha}^{\mathrm{\prime}}}={\mathbf{T}}_{\alpha}+\mathit{\eta$, where $\mathit{\eta}=\left({\eta}_{1},\cdots ,{\eta}_{{N}_{st}\left(t\right)}\right)$ is a vector of independent random variations, drawn from a zerocentered normal distribution of fixed width $m$. With this, $I$ grows by one row and column. The new elements of the matrix are:
where $\beta =1,\cdots ,{N}_{st}\left(t\right)$ and thus ${\alpha}^{\prime}$ inherits its interactions from $\alpha $, but with a small random modification. Evolutionary variations in ITEEM generate mutants that are ecologically similar to their parents. Such variations can represent any phenotypic variation that influences interactions of strains with their community and thus changes their relative competitive abilities (Thompson, 1998; Thorpe et al., 2011; Bergstrom and Kerr, 2015; Thompson, 1999). With Equation 1 we assume that all the interaction terms of the new mutant can change independently.
To implement tradeoff between competitive ability and fecundity, we introduce a relation between competitive ability $C$, defined as average interaction
and replication ${r}_{\alpha}$ (for fecundity). When ${N}_{st}=1$, competitive ability of that single strain is set to zero. To study the influence of tradeoff between competitive ability and replication, we systematically change its shape by varying a parameter $\delta $ $(0\le \delta <1)$ (Figure 2). For details of tradeoff function and its effect on trait distribution and relative fitness see Appendix 1, Tradeoff. Tradeoff functions can be concave ($\delta <0.5$), linear ($\delta =0.5$), or convex ($\delta >0.5$). The tradeoff function ties better competitive ability to lower fecundity and vice versa. The extreme case $\delta =0$ makes $r=1$ and thus independent of $C$, which means no tradeoff.
We compare ITEEM results to the corresponding results of a neutral model (Hubbell, 2001), where we have formally evolving trait vectors ${T}_{\alpha}$ but fixed and uniform replication probabilities and interactions. Accordingly, the neutral model has no tradeoff.
ITEEM belongs to the wellestablished class of generalized LotkaVolterra (GLV) models in the sense that the populationlevel approximation of the stochastic, individualbased ecological dynamics of ITEEM leads to the competitive LotkaVolterra equations (Appendix 1, Generalized Lotka–Volterra (GLV) equation). Thus the results of the model can be interpreted in the framework of competitive GLV equations that model competition for a renewable resource pool and summarize all types of competition (Gill, 1974; Maurer, 1984) in the elements of the interaction matrix $I$ (see above), i.e. these elements represent the resultant negative effect of all competitor populations on each other.
Our model also allows to study speciation in terms of network dynamics. The interaction matrix $I$ defines a complete dominance network between coexisting strains. In this network the nodes are strains ($\alpha ,\beta $), and the directed edges connecting them indicate direction and strength of dominance, i.e. sign and size of ${I}_{\alpha \beta}{I}_{\beta \alpha}$, respectively. Thus, the elements of the weighted adjacency matrix of this network are defined as either ${W}_{\alpha \beta}={I}_{\alpha \beta}{I}_{\beta \alpha}$, if $\alpha $ is the superior competitor in the pairwise encounter with $\beta $ (${I}_{\alpha \beta}>{I}_{\beta \alpha}$), or otherwise as ${W}_{\alpha \beta}=0$. With this definition all ${W}_{\alpha \beta}$ are in $[0,1]$. Accordingly, for the dominance network of species, we computed directed edges between any two species, $i$ and $j$, by averaging over edges between all pairs of strains belonging to these species, that is $W}_{ij}^{sp}={\overline{W}}_{\alpha \beta$ for all strains $\alpha $ and $\beta $ in the $i$th and $j$th species, respectively. The strength and direction of dominance edges indicate the effective flow of population between species.
As we consider a tradeoff between replication and competitive ability in the framework of GLV equations, we can distinguish between $r$ and $\alpha $selection (Gill, 1974; Kurihara et al., 1990; Masel, 2014). $r$selection selects for reproductive ability, which is beneficial in low density regimes, while $\alpha $selection selects for competitive ability and is effective at high density regimes under frequencydependent selection. $\alpha $selection, first introduced by Gill (Gill, 1974), can be realized by acquisition of any kind of ability or mechanism that increases the chance of an organism to take over resources, to prevent competitors from gaining resources (Gill, 1974), or helps the organism to tolerate stress or reduction of contested resource availability (Aarssen, 1984). $\alpha $selection is different from $K$selection; although both are effective at high density, the latter is limited to investments in efficient and parsimonious usage of resources (Masel, 2014).
The source code of the ITEEM model is freely available at GitHub (Farahpour, 2018; copy archived at https://github.com/elifesciencespublications/ITEEM).
Results
Generation of diversity
Our first question was whether ITEEM is able to generate and sustain diversity. Since we have a wellmixed system with initially only one strain, a positive answer implies sympatric diversification: the emergence of new species by evolutionary branching without geographic isolation or resource partitioning. In fact, we observe that during longtime ecoevolutionary trajectories in ITEEM new, distinct species emerge, and their coexistence establishes a sustained high diversity in the system (Figure 3a).
Remarkably, the emerging diversity has a clear hierarchical structure in the phylogeny tree and trait space: at the highest level we see that the phylogenetically separated strains (Figure 3a and Appendix 1, Species and strains) appear as wellseparated clusters in trait space (Figure 3b) similar to biological species. Within these clusters there are subclusters of individual strains (Barraclough et al., 2003). Both levels of diversity can be quantitatively identified as levels in the distribution of branch lengths in minimum spanning trees in trait space (Appendix 1, SMST and distribution of species and strains in trait space). This hierarchical diversity is reminiscent of the phylogenetic structures in biology (Barraclough et al., 2003).
Overall, the model shows evolutionary divergence from one ancestor to several species consisting of a total of hundreds of coexisting strains (Figure 3c). This evolutionary divergence in interaction space is the result of frequencydependent selection without any further assumption on the competition function, for example a Gaussian or unimodal competition kernel (Dieckmann and Doebeli, 1999; Doebeli and Ispolatov, 2010), or predefined niche width (Scheffer and van Nes, 2006). In the course of this diverging sympatric evolution, diversity measures typically increase and, depending on tradeoff parameter $\delta $, high diversity is sustained over hundreds of thousands of generations (Figure 3d, and Appendix 1, Diversity over time). This observation holds for several complementary measures of diversity, no matter whether they are based on abundance of strains or species, or on functional diversity, i.e. quantities that measure the spread of the population in trait space (Appendix 1, Functional diversity (FD), functional group and functional niche).
The observed pattern of divergence contradicts the longheld view of sequential fixation in asexual populations (Muller, 1932). Instead, we see frequently concurrent speciation with emergence of two or more species in quick succession (Figure 3a), in agreement with recent results from longterm bacterial and yeast cultures (Herron and Doebeli, 2013; Maddamsetti et al., 2015; Kvitek and Sherlock, 2013).
ITEEM systems selforganize toward structured communities: the interaction matrix of a diverse system obtained after many generations has a conspicuous block structure with groups of strains with similar interaction strategies (Figure 3e), and these groups being wellseparated from each other in trait space (Figure 3b) (Sander et al., 2015). This fact can be interpreted in terms of functional organization as the interaction trait in ITEEM directly determines the functions of strains and species in the community (Appendix 1, Functional diversity (FD), functional group and functional niche). This means that the block structure in Figure 3e corresponds to selforganized, wellseparated functional niches (Whittaker et al., 1973; Rosenfeld, 2002; Taillefumier et al., 2017), each occupied by a cluster of closely related strains. This niche differentiation among species, which facilitates their coexistence, is the result of frequencydependent selection among competing strategies. Within each functional niche the predominant dynamics, determining relative abundances of strains in the niche, is neutral. Speciation can occur when random genetic drift in a functional group generates sufficiently large differences between the strategies of strains in that group, and then selection forces imposed by biotic interactions reinforce this nascent diversification by driving strategies further apart.
We observe as characteristic of the dynamics of the dominance network $W$ (see Model) the appearance of strong edges as diversification increases trait distance (or dissimilarity) between species (Figure 3f) (Anderson and Jensen, 2005).
Emergence of intransitive cycles
Three or more directed edges in the dominance network can form cycles of strains in which each strain competes successfully against one cycle neighbor but loses against the other neighbor, a configuration corresponding to rockpaperscissors games (Szolnoki et al., 2014). Such intransitive dominance relations have been observed in nature (Buss and Jackson, 1979; Sinervo and Lively, 1996; Lankau and Strauss, 2007; Bergstrom and Kerr, 2015), and it has been shown that they stabilize a system driven by competitive interactions (Allesina and Levine, 2011; Mathiesen et al., 2011; Mitarai et al., 2012; Laird and Schamp, 2015; Maynard et al., 2017; Gallien et al., 2017). We find in ITEEM networks that the increase of diversity coincides with growth of mean strength of cycles (Figure 3d,g and Appendix 1, Intransitive dominance cycles). Note that these cycles emerge and selforganize in the evolving ITEEM networks without any presumption or constraint on network topology.
Formation of strong cycles could also hint at a mechanistic explanation for another phenomenon that we observe in long ITEEM simulations: Occasionally diversity collapses from medium levels abruptly to very low levels, usually followed by a recovery (Figure 3d). Remarkably, dynamics before these mass extinctions are clear exceptions of the generally strong correlation of diversity and average cycle strength. While the diversity immediately before mass extinctions is inconspicuous, these events are always preceded by exceptionally high average cycle strengths (Appendix 1, Collapses of diversity). Because of the rarity of mass extinctions in our simulations we currently have not sufficient data for a strong statement on this phenomenon, however, it is conceivable that the emergence of new species in a system with strong cycles likely leads to frustrations, i.e. the newcomers cannot be accommodated without inducing tensions in the network, and these tensions can destabilize the network and discharge in a collapse. The extinction of a species in a network with strong cycles will probably have a similar effect. This explanation of mass extinctions would be consistent with related works where collapses of diversity occur if maximization of competitive fitness (here: by the newcomer species) leads to a loss of absolute fitness (here: breakdown of the network) (Matsuda and Abrams, 1994; Masel, 2014). This is a special case of the tragedy of the commons (Hardin, 1968; Masel, 2014) that happens when competing organisms under frequencydependent selection exploit shared resources (Rankin and LópezSepulcre, 2005), as it is the case in ITEEM.
Impact of tradeoff and lifespan on diversity
The ecoevolutionary dynamics described above depends on lifespan and tradeoff between replication and competitive ability. This becomes clear if we study properties of dominance network and trait diversity. Figure 4a relates properties of the dominance network to the tradeoff parameter $\delta $, at fixed lifespan $\lambda $. Specifically, we plot two indicators of community structure against tradeoff parameter $\delta $, namely mean weight of dominance edges $\u27e8W\u27e9$, and mean strength of cycles $\rho $. Figure 4b summarizes the behavior of diversity as function of $\delta $ and lifespan $\lambda $. For this summary, we chose ten parameters that quantify different aspects of diversity, for example richness, evenness, functional diversity, and trait distribution, and then averaged over their normalized values to obtain an overall measure of diversity (color bar in the figure). The full set of parameters is detailed in Appendix 1, Diversity indexes and parameters of dynamics for different tradeoffs and lifespans. The resulting phase diagram gives us an overview of the community diversity for different tradeoff parameters $\delta $ and lifespans $\lambda $. The diagram shows a weak dependency of diversity on $\lambda $ and a strong impact of $\delta $, with four distinct phases (IIV) from low to high $\delta $ as described in the following.
Without tradeoff ($\delta =0$), strains do not have to sacrifice replication for better competitive abilities. Any resident community can be invaded by a new mutant with relatively higher $C$ that does not have to compensate with a lower $r$. These mutants resemble Darwinian Demons (Law, 1979), i.e. strains or species that can maximize all aspects of fitness (here $C$ and $r$) simultaneously and would exist under physically unconstrained evolution. Such Darwinian Demons can then be outcompeted by their own mutant offspring’s that have higher $C$ and the same $r$. Thus we have sequential predominance of such strategies with constantly changing traits and improving competitiveness, but no diverse network emerges. As we increase $\delta $ from this unrealistic extreme into phase I ($0<\delta \lesssim 0.2$) coexistence is facilitated. However, the small $\delta $ still favors investing in relatively higher competitive ability as a lowcost strategy to increase fitness. In this phase $\u27e8W\u27e9$ and $\rho $ (Figure 4a) slightly increase: biotic selection pressure exerted by interspecies interactions starts to generate diverse communities (left inset in Figure 4b, Appendix 1, Diversity indexes and parameters of dynamics for different tradeoffs and lifespans).
When $\delta $ increases further (phase II), tradeoff starts to force strains to choose between higher replication or better competitive abilities. Extremes of these quantities do not allow for viable species: sacrificing $r$ completely for maximum $C$ stalls population dynamics, whereas maximum $r$ leads to inferior $C$. Thus strains seek middle ground values in both $r$ and $C$. The nature of $C$ as mean of interactions (Equation 2) allows for many combinations of interaction traits with approximately the same mean. Thus, in a middle range of $r$ and $C$, many strategies with the same overall fitness are possible, which is a condition of diversity (Marks and Lechowicz, 2006). From this multitude of strategies, sets of trait combinations emerge in which strains with different combinations keep each other in check, for example by the competitive rockpaperscissorslike cycles between species described above. An equivalent interpretation is the emergence of diverse sets of nonoverlapping compartments or functional niches in trait space (Figure 3b,e). Diversity in this phase II is the highest and most stable (middle inset in Figure 4b, Appendix 1, Diversity indexes and parameters of dynamics for different tradeoffs and lifespans).
As $\delta $ approaches $0.7$, $\u27e8W\u27e9$ and $\rho $ plummet (Figure 4a) to interaction values comparable to the noise level $m$ (see Model), and a cycle strength typical for the neutral model (horizontal light green ribbon in Figure 4a), respectively. The sharp drop of $\u27e8W\u27e9$ and $\rho $ at $\delta \approx 0.7$ is reminiscent of a phase transition. As expected for a phase transition, the steepness increases with system size (Appendix 1, Size of the system). For $\delta \gtrsim 0.7$, weights of dominance edges never grow and no structures, for example cycles, emerge. Diversity remains low and close to that of a neutral system. The sharp transition at $\delta \approx 0.7$ which is visible in practically all diversity measures (between phases II and III in Figure 4b, see also Appendix 1, Diversity indexes and parameters of dynamics for different tradeoffs and lifespans) is a transition from a system dominated by biotic selection pressure to a neutral system. In high tradeoff phase III, a small relative change in $C$ produces a large relative change in $r$ (Appendix 1, Strength of tradeoff function). For instance, given a resident strain $R$ with $r$ and $C$, a closely related mutant $M$ increases the fitness by adopting a relatively high $r$ while paying a relatively small penalty in $C$ (see Appendix 1, Strength of tradeoff function for the relative impacts of the traits), and therefore will invade $R$. Thus, diversity in phase III will remain stable and low, and is characterized by a group of similar strains with no effective interaction and hence no diversification to distinct species (right inset in Figure 4b and Appendix 1, Diversity indexes and parameters of dynamics for different tradeoffs and lifespans). In this high tradeoff regime, lifespan comes into play: here, decreasing $\lambda $ can make lives too short for replication. These hostile conditions minimize diversity and favor extinction (phase IV).
Tradeoff, resource availability, and diversity
There is a wellknown but not well understood unimodal relationship (‘humpback curve’) between biomass productivity and diversity: diversity as function of productivity has a convex shape with a maximum at middle values of productivity (Smith, 2007; Vallina et al., 2014). This productivitydiversity relation has been reported at different scales in a widerange of natural communities, for example phytoplankton assemblages (Vallina et al., 2014), microbial (Kassen et al., 2000; HornerDevine et al., 2003; Smith, 2007), plant (Guo and Berry, 1998; Michalet et al., 2006), and animal communities (Bailey et al., 2004). This behavior is reminiscent of horizontal sections through the phase diagram in Figure 4b, though here the driving parameter is not productivity but tradeoff. However, we can make the following argument for a monotonic relation between productivity and tradeoff shape. First we note that biomass productivity is a function of available resources (Kassen et al., 2000): the larger the available resources, the higher the possible productivity. This allows us to argue in terms of available resources. For ecoevolutionary systems with scarce resources, species with high replication rates will have low competitive ability because for each individual of the numerous offspring there is little material or energy available to develop costly mechanisms that increase competitive ability. On the other hand, if a species under these resourcelimited conditions produces competitively constructed individuals it cannot produce many of them. This argument shows a correspondence between a resourcelimited condition and high $\delta $ for tradeoff between replication and competitive ability. At the opposite, rich end of the resource scale, evolving species are not confronted with hard choices between replication rate and competitive ability, which is equivalent to low $\delta $. Taken together, the tradeoff axis should roughly correspond to the inverted resource axis: high $\delta $ for poor resources (or low productivity) and low $\delta $ for rich resources (or high productivity); a detailed analytical derivation will be presented elsewhere. The fact that ITEEM produces this frequently observed humpback curve proposes tradeoff as underlying mechanism of this productivitydiversity relation.
Frequencydependent selection
Observation of ecoevolutionary trajectories as in Figure 3 suggested the hypothesis that speciation and extinction events in ITEEM simulations do not occur at a constant rate and independently of each other, but that one speciation or extinction makes a following speciation or extinction more likely. Such a frequencydependence occurs if emergence or extinction of one species creates the niche for emergence and invasion of another species, or causes its decline or extinction (Herron and Doebeli, 2013). Without frequencydependence such evolutionary events should be uncorrelated.
To test for frequencydependent selection we checked whether the probability distribution of interevent times (time intervals between consecutive speciation or extinction events) is compatible with a constant rate Poisson process, i.e. a purely random process, or whether such events are correlated (Appendix 1, Frequencydependent selection). We find that for long interevent times the decay of the distribution in ITEEM simulations is indistinguishable from that of a Poisson process. However, for shorter times there are significant deviations from a Poisson process for speciation and extinction events: at interevent times of around ${10}^{4}$ the probability decreases for a Poisson process but significantly increases in ITEEM simulations. Thus, the model shows frequencydependent selection with the emergence of new species increasing the probability for generation of further species, and the loss of a species making further losses more likely. This behavior of ITEEM is similar to microbial systems where new species open new niches for further species, or the loss of species causes the loss of dependent species (Herron and Doebeli, 2013; Maddamsetti et al., 2015).
The above analysis illustrates a further application of ITEEM simulations. Ecoevolutionary trajectories from ITEEM simulations can be used to develop analytical methods for the inference of competition based on observed diversification patterns. Such methods could be instrumental for understanding the reciprocal effects of competition and diversification.
Effect of mutation on diversity
Mutations are controlled in ITEEM by two parameters: mutation probability $\mu $, and width $m$ of trait variation. In simulations, diversity grew faster and to a higher level with increasing mutation probability ($\mu ={10}^{4},5\times {10}^{4},{10}^{3},5\times {10}^{3}$), but without changing the overall structure of the phase diagram (Appendix 1, Mutation probability). One interesting tendency is that for higher $\mu $, the lifespan becomes more important at the interface of regions III and IV (high tradeoffs), leading to an expansion of region III at the expense of the hostile region IV: long lifespans in combination with high mutation probability establish low but viable diversity at large $\delta $. The humpback curve of diversity over $\delta $ is observed for all mutation probabilities. Thus, the diversity in ITEEM is not a simple result of a mutationselection balance but tradeoff plays an important role in shaping diversity in trait space.
The width of trait variation, $m$, influences both the speed of evolutionary dynamics and the maximum variation inside species, i.e. clusters of strains. The smaller $m$ the slower the dynamics and the smaller the clusters. However extreme values of $m$ can completely suppress the diverging evolution: Very small variations are wiped out by rapid ecological dynamics, and very large variations disrupt selection forces by imposing big fluctuations.
Comparison of ITEEM with neutral model
The neutral model introduced in the Model section has no meaningful interaction traits, and consequently no meaningful competitive ability or tradeoff with fecundity. Instead, it evolves solely by random drift in trait space. Similarly to ITEEM, the neutral model generates clumpy structures of traits (Appendix 1, Neutral model), though here the clusters are much closer and thus the functional diversity is much lower. This can be demonstrated quantitatively by the size of the minimum spanning tree of populations in trait space that are much smaller for the neutral model than for ITEEM at moderate tradeoff (Appendix 1, Neutral model). The clumpy structures generated with the neutral model do not follow a stable trajectory of divergent evolution, and, hence, niche differentiation cannot be established. In a neutral model, without frequencydependent selection and tradeoff, stable structures and cycles cannot form in the community network, and consequently, diversity cannot grow effectively (Appendix 1, Neutral model). The comparison with the neutral model points to frequencydependent selection as a promoter of diversity in ITEEM. For high tradeoffs (region III in Figure 4b), diversity and number of strong cycles in ITEEM are comparable to the neutral model (Figure 4a).
Discussion
Phenotype traits and interaction traits
In established ecoevolutionary models, organisms are described in terms of one or a few phenotype traits. In contrast, the phenotype space of real systems is often very highdimensional; competitive species in their evolutionary arms race are not confined to few predefined phenotypes but rather explore new dimensions in that space (Maharjan et al., 2006; Maharjan et al., 2012; Zaman et al., 2014; Doebeli and Ispolatov, 2017). Coevolution systematically pushes species toward complex traits that facilitate diversification and coexistence (Zaman et al., 2014; Svardal et al., 2014), and evolutionary innovation frequently generates phenotypic dimensions that are completely novel in the system (Doebeli and Ispolatov, 2017). Complexity and multidimensionality of phenotype space have recently been the subject of several experimental and theoretical studies with different approaches that demonstrate that evolutionary dynamics and diversification in highdimensional phenotype trait space can produce more complex patterns in comparison to evolution in lowdimensional space (Doebeli and Ispolatov, 2010; Gilman et al., 2012; Svardal et al., 2014; Kraft et al., 2015; Doebeli and Ispolatov, 2017). For example, it has been shown that the conditions needed for frequencydependent selection to generate diversity are satisfied more easily in highdimensional phenotype spaces (Doebeli and Ispolatov, 2010). Moreover, the level at which diversity saturates in a system depends on its dimensionality, with higher dimensions allowing for more diversity (Doebeli and Ispolatov, 2017), and the probability of intransitive cycles in species competition networks grows rapidly with the number of phenotype traits. The conventional way to tackle this problem is to use models with a larger number of phenotype traits. However, this is not really a solution of the problem because this still confines evolution to the chosen fixed number of traits, and it also makes these models more complex and thus computationally less tractable. As will be discussed below, interactionbased models such as ITEEM offer a natural solution to this problem by mapping the system to an interaction trait space that can dynamically expand by the emergence of novel interaction traits as ecoevolutionary dynamics unfolds.
Ecoevolutionary dynamics in interaction trait space
Interactionbased ecoevolutionary models rely on the assumption that phenotypic evolution can be coarsegrained to the interaction level (Figure 1). This means that regardless of the details of phenotypic variations, we just study the resultant changes in the interaction network. In an ecoevolutionary system dominated by competition this is justified because phenotypic variations are relevant only when they change the interaction of organisms, directly or indirectly; otherwise they do not impact ecological dynamics. The interaction level is still sufficiently detailed to model macroevolutionary dynamics that are dominated by ecological interactions.
A transition from phenotype space to interaction space requires a mapping from the former to the latter, based on the rules that characterize the interaction of individuals with different phenotypic traits. As a concrete example, we might consider the competition kernel of adaptive dynamics models (Doebeli, 2011) that determines the competitive pressure of two individuals with specific traits. That formalism describes well how, after mapping phenotypic traits to the interaction space, ecological outcome eventually is determined by interactions between species. In Appendix 1, Phenotypeinteraction map, some properties of this mapping are discussed.
Interactionbased models
In the first interactionbased model by Ginzburg et al. (1988), emergence of a new mutant was counted as speciation, and it was shown that simulating speciation events as ecologically continuous mutations in the strength of competitive interactions resulted in stable communities. However the Ginzburg model produced stable coexistence of only a few similar interaction traits, without branching and diversification to distinct species. As outlined in the introduction, subsequent interactionbased models tried to solve this problem by supplementing the Ginzburg model with some ad hoc features. For example, Tokita and Yasutomi (2003) mixed mutualistic and competitive interactions, and showed that only local mutations, i.e. changes in one pairwise interaction rate, can produce stable diversity. Recently, Shtilerman et al. (2015) enforced diversification in purely competitive communities by imposing a large parentoffspring niche separation. To our knowledge, ITEEM is the first interactionbased model in which, despite its minimalism and without ad hoc features, diversity gradually emerges under frequencydependent selection by considering physical constraints of ecoevolutionary dynamics.
In all previous interactionbased models, ecoevolutionary dynamics has been divided into iterations over two successive steps: each first step of continuous population dynamics, implemented by integration of differential equations, was followed by a stochastic evolutionary process, namely speciation events and mutations, as a second step. However, in nature these two steps are not separated but intertwined in a single nonequilibrium process. Hence, the artificial separation necessitated the introduction of model components and parameters that do not correspond to biological phenomena and observables. In contrast, individualbased models like ITEEM operate with organisms as units, and efficiently simulate ecoevolutionary dynamics in a more natural and consistent way, with parameters that correspond to biological observables.
Tradeoff anchors ecoevolutionary dynamics in physical reality
Lifehistory tradeoffs, like the tradeoff between replication and competitive ability, now experimentally established as essential to living systems (Stearns, 1989; Agrawal et al., 2010; Masel, 2014), are inescapable constraints imposed by physical limitations in natural systems. Our results with ITEEM show that tradeoffs fundamentally impact ecoevolutionary dynamics, in agreement with other ecoevolutionary models with tradeoff (Huisman et al., 2001; Bonsall et al., 2004; de Mazancourt and Dieckmann, 2004; Beardmore et al., 2011). Remarkably, we observe with ITEEM sustained high diversity in a wellmixed homogeneous system. This is possible because moderate lifehistory tradeoffs force evolving species to adopt different strategies or, in other words, lead to the emergence of wellseparated functional niches in interaction space (Gudelj et al., 2007; Beardmore et al., 2011).
Given the accumulating experimental and theoretical evidence, the importance of tradeoff for diversity is becoming more and more clear. ITEEM provides an intuitive and generic conceptual framework with a minimum of specific assumptions or requirements. This makes the results transferable to different systems, for example biological, economical and social systems, wherever competition is the driving force of evolving communities. Put simply, ITEEM shows generally that in a barebone ecoevolutionary model withal standard population dynamics (birthdeathcompetition) and a basic evolutionary process (mutation), diverse set of strategies will emerge and coexist if physical constraints force species to manage their resource allocation.
Power and limitations of ITEEM
Despite its minimalism, ITEEM reproduces in a single framework several phenomena of ecoevolutionary dynamics that previously were addressed with a range of distinct models or not at all, namely sympatric and concurrent speciation with emergence of new niches in the community, mass extinctions and recovery, large and sustained functional diversity with hierarchical organization, spontaneous emergence of intransitive interactions and cycles, and a unimodal diversity distribution as function of tradeoff between replication and competition. The model allows detailed analysis of ecoevolutionary mechanisms and could guide experimental tests.
The current model has important limitations. For instance, the tradeoff formulation was chosen to reflect reasonable properties in a minimalist way. This should be revised or refined as more experimental data become available. Secondly, individual lifespans in this study came from a random distribution with an identical fixed mean. Hence we have no adaptation and evolutionarybased diversity in lifespan. This limits the applicability of the current model to communities of species that have similar lifespans, and that invest their main adaptation effort into growth or reproduction and competitive ability. Furthermore, our model assumes an undefined pool of steadily replenished shared resources in a wellmixed system. This was motivated by the goal of a minimalist model for competitive communities that could reveal mechanisms behind diversification and niche differentiation, without resource partitioning or geographic isolation. However, in nature, there will in general be few or several limiting resources and abiotic factors that have their own dynamics. For this scenario, which is better explained by a resourcecompetition model than by the GLV equation, it is possible to consider resources as additional rows and columns in the interaction matrix $I$ and in this way to include abiotic interactions as well as biotic ones.
In an interactionbased model like ITEEM the interaction terms of the mutants change gradually and independently (Equation 1). This assumption of random exploration of interaction space can be violated, for example, in simplified models with few fixed phenotypic traits. Further studies are necessary to investigate the general properties and restrictions of the map between phenotype and interaction space. In Appendix 1, Phenotypeinteraction map we briefly introduced and discussed some properties of this map.
Appendix 1
Species and strains
There is no universally accepted definition of species (Zachos, 2016), especially for asexual populations (Zachos, 2016; Birky Jr and Barraclough, 2009; Richards, 2013). In the present work, we follow RossellóMora and Amann (2001) and use the concept of phylophenetic species applicable to asexual populations. A phylophenetic species is defined as a monophyletic cluster of strains that show a high degree of overall similarity with respect to many independent characteristics (RossellóMora and Amann, 2001). We used genealogical trees to define species as group of strains that share a most recent common ancestor and are separated by longlasting gaps in the tree. Each of these clusters that is branching off from a point of divergence in the tree was counted as an individual species if it has existed for more than a certain number of generations ($7\times {10}^{4}$ in the manuscript), considering all of its subbranches. Changing the threshold in a range from $2\times {10}^{4}$ to $10}^{5$ had no profound impact on the results. Branches that lasted shorter than this threshold were counted as strains of their parents (Appendix 1—figure 1).
Distribution of strains in trait space and their diversification (see Figure 3b and c of the main text) shows that these monophyletic clusters are also the wellseparated clusters in functional space which together allow us to consider them as phylophenetic species.
See Birky Jr and Barraclough, (2009); Ereshefsky (2010); Richards (2013); Wilkins (2018) for a discussion on pros and cons of this definition.
Tradeoff
Tradeoff function
Linear, concave and convex tradeoffs among different lifehistory traits have been reported in many studies (Jessup and Bohannan, 2008; Maharjan et al., 2013; Saeki et al., 2014; Ferenci, 2016). It has been shown that different forms of tradeoff reproduce different diversity and coexistence patterns (Levins, 1968; Maharjan et al., 2013; Kasada et al., 2014; Ehrlich et al., 2017). The shape of tradeoff is determined by various factors like quantitative relationship between resource allocations in lifehistories (Saeki et al., 2014), physiological mechanisms (Bourg et al., 2017), and environment (Jessup and Bohannan, 2008).
In this study, to implement a tradeoff between reproduction ${r}_{\alpha}$ and competitive ability $C\left({T}_{\alpha}\right)$ with a variable form, we used a function with one shape parameter $s$:
where $\mathbf{T}}_{\alpha$ is the trait vector of strain $\alpha$. This tradeoff function maps competitive ability of species to reproduction probability in the range $[0,1]$. By changing the exponent $s$ between simulations, we can study the effect of tradeoff shape on ecoevolutionary dynamics. For a systematic scan of tradeoff shapes (Figure 2 in the main text), we formulated the shape parameter $s$ as
with tradeoff parameter $\delta $ covering $[0,1]$ in equidistant steps. Of course, other functional forms of the tradeoff are conceivable.
Tradeoff and explored interaction trait space
Tradeoffs between lifehistory traits are constraints imposed by fundamental resourceallocation principles. They confine evolutionary adaptations and innovations to a permissible subspace of all trait combinations. However the outcome of evolution – determined by the underlying mechanisms of the system and selection forces – is a subset of this permissible subspace, occupied by the selected coexisting organisms, in which all organisms should have more or less the same fitness to be able to coexist. Thus, in each community and at each time point, just a small part of this permissible space is usually occupied (Bourg et al., 2017).
In ITEEM, the competitive ability $C$ of a strain (Equation 2 of the main text) quantifies how successfully individuals of this strain compete against individuals of all coexisting strains in direct encounters. Thus, $C$ is a relative and densitydependent component of the fitness. In the course of evolutionary dynamics, $C$ never explored extreme regimes, which means that we never observed $C\approx 0$ (organisms that fail nearly in all encounters) or $C\approx 1$ (organisms that defeat nearly all the rivals). Instead, we saw a distribution around middle values, $C\approx 0.5$. In the low tradeoff regime, emergence of strategies with relatively high competitiveness, without a considerable cost in reproduction, drives the system to low diversity by outfighting the competitors. In this case, as $C$ is a relative, interactionbased measure (Equation 2 of the main text) extinction of species with low competitive ability pushes the distribution of $C$ again toward $0.5$. In the high tradeoff regime, even small gains in $C$ come with a severe penalty in $r$, that is, increasing competitive ability to high values is very unlikely; thus, strategies with a relatively high $r$ can prevail by a negligible decrement in their relative competitive ability, which again limits diversity of strategies and brings them back to $C\approx 0.5$. For moderate tradeoff values between the above limiting cases, $C$ also stabilizes around 0.5, as described in the main text (Results, Impact of tradeoff and lifespan on diversity).
Strength of tradeoff function
At first sight, the strength of a tradeoff function – how strongly changes in one trait influence the other trait – seems to be just a synonym for the slope (=first derivative) of the tradeoff function. For instance, consider two traits $x$ and $y$ that both contribute to the fitness and are related by tradeoff function $y=f\left(x\right)$ with first derivative ${f}^{\prime}\left(x\right)$. A small change $\Delta x$ in trait $x$ will cause a change $\Delta y\approx {f}^{\prime}\left(x\right)\Delta x$ in trait $y$, so that, obviously, for given $\mathrm{\Delta}x$ the slope ${f}^{\prime}\left(x\right)$ determines the change $\Delta y$. However, the effect of such a change will very much depend on the community context: the same change $\Delta x$ or $\Delta y$ may be relatively large or relatively small, depending on the actual values of $x$ and $y$. For example, a $\Delta y=0.1$ will change $y=0.1$ by 100%, but a larger $y=0.8$ by a mere 12.5%, and accordingly the change $\Delta y$ will have different effects on the fitness of the affected strain. Therefore, we define as tradeoff strength $\theta (x,y)$ the ratio ${\scriptscriptstyle \frac{\Delta y}{y}}/{\scriptscriptstyle \frac{\Delta x}{x}}$, or, in the limit of $\Delta x,\Delta y\to 0$ as
To demonstrate that the tradeoff strength $\theta $ is a more meaningful quantity than the slope to characterize the effect of the tradeoff, we discuss in the following a few characteristic cases.
Assume the same strain with a trait value $x$ in two different systems, 1 and 2, with different tradeoff functions ${f}_{1},{f}_{2}$, respectively. In the first system, $x$ may be mapped to a large value of trait ${y}_{1}={f}_{1}\left(x\right)$, while in the second it may be mapped to a small value of trait ${y}_{2}={f}_{2}\left(x\right)$. Emergence of a mutant with a small change $\Delta x$ in trait $x$ causes different variations in the two systems, namely $\mathrm{\Delta}{y}_{1}\approx {f}_{1}^{\prime}\left(x\right)\mathrm{\Delta}x$ among large traits ${{\displaystyle \overline{y}}}_{1}$, and $\mathrm{\Delta}{y}_{2}\approx {f}_{2}^{\prime}\left(x\right)\mathrm{\Delta}x$ among small traits ${{\displaystyle \overline{y}}}_{2}$. Thus the same $\Delta x$ impacts the two systems differently, even if the slopes ${{\displaystyle {f}^{\prime}}}_{i}\left(x\right)$ would be the same. The tradeoff strength ${\scriptscriptstyle \frac{\Delta y}{y}}/{\scriptscriptstyle \frac{\Delta x}{x}}$ captures this difference as it is smaller for system 1 and larger for system 2.
$\theta $ is also expressive if the two systems have different values of the first trait ${x}_{1}<{x}_{2}$ that are mapped to the same second trait $y$ and have the same slope of tradeoff function in the respective range (${f}_{1}^{\prime}\left({x}_{1}\right)={f}_{2}^{\prime}\left({x}_{2}\right)$). The same changes in the first trait $\Delta x$ have different impacts on relative fitness of the two systems, which is again captured by tradeoff strength ${\scriptscriptstyle \frac{\Delta y}{y}}/{\scriptscriptstyle \frac{\Delta x}{x}}$. Only in the special case ${x}_{1}={x}_{2}$ and ${y}_{1}={y}_{2}$, the slope $\frac{\mathrm{\Delta}y}{\mathrm{\Delta}x}$ is sufficient to compare the effects of tradeoffs on fitness and dynamics.
For $x=C$ and $y=r$ the derivative ${\scriptscriptstyle \frac{dy}{dx}}$ is the derivative of the tradeoff function in Figure 2 of the main text, explicitly formulated in Equation 3.
In ITEEM, as explained in the previous section, the competitive ability $C$ is typically in the middle range, i.e. organisms with low or high competitive ability are rare. In this middle range, tradeoff strength $\theta (C,r)={\scriptscriptstyle \frac{dr}{dC}}{\scriptscriptstyle \frac{C}{r}}$ increases with increasing tradeoff parameter $\delta $, as shown in Appendix 1—figure 2 below.
Generalized Lotka–Volterra (GLV) equation
As an individualbased model, ITEEM simulates systems consisting of distinct, interacting organisms, and thus can model nonequilibrium dynamics, demographic fluctuations, effects of diverse lifespans, and other features of real systems, as discussed in the main text. If these features were not of concern we could replace the ecological dynamics of ITEEM by the corresponding populationlevel model. In the following we show that such an abstraction of ecological interactions of ITEEM leads to the competitive generalized LotkaVolterra (GLV) equation.
We start from the main equation of population dynamics for our model:
In which ${x}_{\alpha}={\scriptscriptstyle \frac{{n}_{\alpha}}{{N}_{s}}}$ is the relative abundance or probability of finding strain $\alpha $ in the system (${N}_{s}$ is the number of sites in the system). The first term on the right side of Equation 6 is the growth of the population of strain $\alpha $ when it produces progeny that is able to find an empty space in the system. The second term shows the growth of population $\alpha $ when after reproduction its offspring is able to invade a site occupied by another individual, the third term is the decrease of population $\alpha $ due to invasion by offspring of other strains and the last term is the decrease of population $\alpha $ due to the intrinsic death rate because of the attributed life span ($d=1/\lambda $). We can rewrite the equation as follows:
In which we used $1{I}_{\alpha \beta}={I}_{\beta \alpha}$ (see Equation 1 in the main text). The last equation above, which we can rewrite compactly as
is the GLV equation. $g=rd$ is the effective population growth rate and $A$ is the community matrix; its elements ${A}_{\alpha \beta}=({r}_{\alpha}+{r}_{\beta}){I}_{\beta \alpha}$ are always negative in our system which shows that ITEEM strains and species are competing.
The close relationship of our individualbased ecological dynamics with the GLV equation shows that organisms are competing in the sense that they expand their populations at the expense of their competitors populations to secure resources and increase fitness. This similarity of GLV with ITEEM ecological dynamics also explains why ITEEM individualbased dynamics corresponds to a wellmixed system. A ‘site’ in the model is not a patch of space or a piece of a spatially structured resource – neighborhood has no meaning, as in the GLV model. Instead, a ‘site’ stands for a discrete portion of the steadily replenished resource pool that is equally accessible to all extant individuals, and sufficient for their respective metabolisms. Being wellmixed means that any individual meets any other individual and site with the same probability. A difference between the individualbased ITEEM and the populationlevel version in Equation 7 is that the former models encounters at the level of individuals whereas the latter maps encounters to interactions between populations.
As outlined above, ITEEM simulations can be interpreted in the framework of the competitive GLV equation. In evolving systems governed by this equation, fitness is determined by reproduction rate, carrying capacity and competitive abilities (Gill, 1974; Masel, 2014). In the present model, the carrying capacity is the same for all strains and species, and the fitness at the low density limit (corresponding to the initial phase of the simulations) is determined by replication $r$. At the high density limit typically simulated in the present work, the product of $r$ and $C$ determines the fitness (Masel, 2014).
Classical multidimensional scaling (CMDS)
Multidimensional scaling (MDS) algorithms take sets of points in $N$dimensional space and represent them in a lowerdimensional space (typically 2dimensional) so that the original distances in $N$dimensional space are preserved as well as possible. The lowerdimensional representation can then be easily visualized.
Classical MDS (CMDS) is a member of the family of MDS methods (Cox and Cox, 2000; Borg and Groenen, 2005; Wang, 2012). The algorithm takes as input an $N\times N$ distance matrix, where the distances could for example be dissimilarities between pairs of $N$ objects, and outputs a coordinate matrix that determines positions of the points in a lowerdimensional (often 2dimensional) space with the condition of minimizing the loss function that measures discrepancy between the algorithm’s output and the real distances. The quality of the lowerdimensional representation can be assessed from the eigenvalues of a factor analysis that shows the fraction of variation in the data explained by each dimension. CMDS is mathematically closely related to principal component analysis (PCA) (Wang, 2012).
In our analysis, objects are interaction traits of strains and our aim is to visualize the distribution of them in trait space. To this end we first calculate the Euclidean distances between all trait vectors
and then, by applying the cmdscale function of the R software (version 3.3.0) to that distance matrix, we project our trait space into 2dimensional plots. Thus, each point in the 2dimensional CMDS plot represents a trait vector, i.e. a strain. The larger the distance between two points, the more different the traits of corresponding strains. While CMDS is a handy tool for the visualization of evolutionary processes in trait space (Figure 3b and f of the main text), all quantitative analyses were performed in the original highdimensional space.
In the very early stages of evolution (Appendix 1—figure 3top) strains are very similar so that two dimensions are not sufficient to represent their dissimilarities accurately (relatively high eigenvalues beyond the 2nd eigenvalue in top right of Appendix 1—figure 3). But when evolutionary speciation’s and branching’s occur, the 2dimensional space is more appropriate, as the eigenvalues from the factor analysis show (Appendix 1—figure 3bottom and Appendix 1—figure 4).
SMST and distribution of species and strains in trait space
The minimum spanning tree (MST) has been used as a tool to characterize the distribution of species and strains in interaction trait space. The MST of a graph is a tree that connects all nodes of that graph so that the total edge weight is minimized. For $N$ nodes, the MST has $N1$ edges. Specifically, if we have $N$ strains in interaction trait space as nodes, we compute as MST a tree of $N1$ edges that links all strains with a minimum sum of edge lengths (= distances between strains in interaction trait space [Equation 8]). For a more familiar example think of nodes as $N$ cities on a map and the MST as a tree of $N1$ edges that links all cities, so that the sum of edge lengths (= distances between pairs of cities) is minimized.
The sum of edges of the MST (SMST) can be used as a quantitative measure that characterizes the distribution of points. The MST of an evenly distributed structure, with $N$ nodes, has $N1$ edges with more or less the same length and thus, the SMST scales with $N$. For a hierarchical structure, MST consists of long edges between clusters and short edges that connect the nodes inside each cluster; in this case, the SMST scales with the number of clusters. The SMST increases by divergent evolution. The left panel of Appendix 1—figure 5 shows an example MST in the trait space of our simulation, represented in 2D for the sake of visualization.
The right panel of Appendix 1—figure 5 plots the lengths of the lengthsorted edges of the MST versus their ranks for 500 simulation snapshots ($\delta =0.5$, $\lambda =300$, $\mu =0.001$, $m=0.02$). From this loglogscaled plot we see that there are two clearly different scales in the size of edges. ${R}_{1}$ is a representative value for the size of clusters (distance between strains within a typical species) and ${R}_{2}$ is a representative value for the scale of the trait space (typical distance between different species). $N$ in Appendix 1—figure 5 represents approximately the number of distinct clusters (species) in the system.
Diversity indexes and parameters of dynamics
A diversity index quantifies a certain aspect of diversity in a single number. Since diversity is itself complex, no single diversity index is sufficient to describe the diversity of a community. For example richness, i.e. number of species, has no information about the distribution of the population among species. Evenness or Shannon entropy takes into account this distribution but does not inform about the diversity of the trait of species, i.e. how diverse is a system with respect to the functionality of its species. Functional diversity indexes focus on this aspect but none of them exhaustively describes properties of trait space (Mouchet et al., 2010). For a comprehensive assessment of diversity and community dynamics, information about the density of species over resources, rate of extinction and emergence, and also details of community structure, for example interaction of species and topology of the network, should be considered, too.
Diversity over time
The next plots show how SMST, an index for functional diversity, changes over time. Note that evolutionary collapses (mass extinctions) occasionally occur (see Appendix 1, Collapses of diversity) with a probability that depends on the tradeoff parameter, lifespan and mutation probability. The plots in Appendix 1—figure 6 follow SMST (see Appendix 1, SMST and distribution of species and strains in trait space) over time for different tradeoffs $\delta $ but the same lifespan $\lambda $ (and mutation probability) in each plot.
Appendix 1—figure 7 follows SMST over time for different lifespans $\lambda $ but the same tradeoff $\delta $ in each plot. Comparison of Appendix 1—figure 7 and Appendix 1—figure 6 confirms that diversity and dynamics are strongly associated with tradeoff without a noticeable effect of lifespan (except for very short lifespans, upper left panel of Appendix 1—figure 6 and bottom panel of Appendix 1—figure 7).
Functional diversity (FD), functional group and functional niche
Univariate diversity indexes that are defined based on abundance of species, like richness and evenness, are routinely used to quantify diversity in biological communities. These indexes are most expressive if species are equal in their effect on their community and ecosystem functioning. However, in the last decades ecologists are increasingly realizing that without information on variety of functions in a community, diversity can not be correctly evaluated (Mouchet et al., 2010), and that traitbased measures that reflect the importance, essentiality, or redundancy of species may be more meaningful than abundancebased measures (Cadotte et al., 2011). Inspired by the concept of Hutchinsonian niche, functional diversity (FD) was introduced by Rosenfeld as distribution of species in functional space (Whittaker et al., 1973; Rosenfeld, 2002). The axes of this space represent the functional features of species (Mouchet et al., 2010) which are usually measurable characteristics (traits) that are indicators of organismal performance, and that are associated with species fitness and their ecological function (Violle et al., 2007; Májeková et al., 2016). In ITEEM, we operate with the interaction trait, which is already an optimal indicator of function of strains and species: the whole vector of interactions that determines the role or function of strain/species in the actual community (Hooper et al., 2002; Sander et al., 2015). Thus, the distribution of interaction traits in trait space determines variety of functions in the system.
Different measures of FD have been introduced, each quantifying and explaining one facet of trait distribution in trait or functional space, very similar to SMST (Appendix 1, SMST and distribution of species and strains in trait space). In the following we also use three other indexes: functional dispersion, Rao index and functional evenness (Mouchet et al., 2010; Mason et al., 2013).
Distinct, wellseparated clusters in functional (trait) space mean that species are diversified to different functional groups. A functional group is defined in ecology either as a set of species with similar effect on their environment, or as cluster in trait space (Hooper et al., 2002). In ITEEM, by using the framework of interactionbased models, these two definitions are interchangeable.
The positions of functional groups in functional space define their functional niches. The notion of functional niche was first introduced by Elton as the place of an animal in its community or its biotic environment (Elton, 1927). Then Clarke (Clarke, 1954) noted that the functional niche stresses the function of the species in the community, which is different from its physical niche, the latter determining its place in the habitat. A suitable definition of functional niche is the area occupied by a species in the $n$dimensional functional space (Clarke, 1954; Whittaker et al., 1973; Rosenfeld, 2002). This concept was subsequently differentiated by Odum who considered the habitat as the organism’s ‘address’ and the niche as its ‘profession’ (Odum, 1959). Following this picture and considering that there is no physical niche or habitat in our wellmixed model, we can say that in ITEEM, the position of trait vectors in functional space determines the profession (function) of species. In the ecoevolutionary dynamics of ITEEM, distinct functional groups with different professions/roles/functions emerge in a community of competitive organisms.
Diversity indexes and parameters of dynamics for different tradeoffs and lifespans
For the phase diagram in Figure 4 of the main text we have synthesized a descriptive dimensionless diversity parameter by averaging over normalized values of several diversity indexes, namely richness, Shannon entropy, standard deviation of replication $r$, maximum distance in trait space, standard deviation of interaction terms, sum of squared lengths of minimum spanning tree of trait space, functional diversity indexes (functional dispersion, Rao index and functional evenness, all three in two versions: with and without abundance), and strength of cycles. This phase diagram gives a good overview about different characteristics of communities with different $\delta $ and $\lambda $, but, of course, the averaging procedure leads to a loss of detailed information. Therefore, we report in Appendix 1—figure 8 some of the most important indexes of community state computed from our simulations for different tradeoffs and lifespans. Each index value is averaged over $5\times {10}^{6}$ generations. The four phases described in Figure 4 of the main text can be seen in nearly all the parameters. Functional diversity indexes (functional dispersion, functional evenness and Rao’s quadratic entropy) are calculated using the dbFD function in Rpackage FD, version 1.0–12.
Intransitive dominance cycles
Flow of energy/mass between strains in ITEEM community is determined by the dominance matrix ${W}_{\alpha \beta}$:
The corresponding network is a directed network with one directed edge between each pair of nodes (strains), pointing from the dominating to the dominated one, with weight between 0 and 1 according to (Equation 9). Three or more directed edges in the dominance network can form cycles of strains in which each strain competes successfully against one cycle neighbor but loses against the other neighbor, a configuration corresponding to the rockpaperscissors game. Even in a completely connected random dominance network, a randomly selected triplet of nodes forms a cycle with a probability of ${\scriptscriptstyle \frac{1}{4}}$. We are interested in characterizing evolved networks in ITEEM in comparison to random networks. Hence, we compare number, ${N}_{cyc}^{Network}$, and average strength, ${S}_{cyc}^{Network}$, of cycles of the evolving network at each time step with number, ${N}_{cyc}^{Random}$, and average strength, ${S}_{cyc}^{Random}$ of cycles of its equivalent shuffled random networks. For this purpose we
average over cycles: We select at random $3$ nodes of the network and check if they form a cycle of size $3$. If yes, the number of 3cycles, ${N}_{cyc}^{Network}$, of the network increases by one and the minimum weight among its edges (corresponding to the limiting edge in that cycle for energy/mass flow) is the strength of that cycle. This procedure is repeated many times ($>{10}^{5}$), and then we average over the strengths of all cycles to obtain average strength of cycles, ${S}_{cyc}^{Network}$.
build the equivalent random networks: We shuffle the edges of original network to obtain a random network. Then we apply the procedure described in step one on this network to measure the number of cycles and their average strengths. This step is done several times ($>10$) to sample different random networks, and by averaging over them we obtain ${N}_{cyc}^{Random}$ and ${S}_{cyc}^{Random}$.
normalize the values: Number and average strength of ITEEM network are normalized to the number and average strength of the corresponding random network, respectively: ${N}_{cyc}^{Network}={\scriptscriptstyle \frac{{N}_{cyc}^{Network}}{{N}_{cyc}^{Random}}}$ and ${S}_{cyc}^{Network}={\scriptscriptstyle \frac{{S}_{cyc}^{Network}}{{S}_{cyc}^{Random}}}$.
The results of these three steps are plotted in Figure 3g and Figure 4a of the main text.
Collapses of diversity
Collapses of diversity occur occasionally in ITEEM simulations. The probability of collapses depends on the tradeoff parameter $\delta $, attributed lifespan $\lambda $, and mutation probability $\mu $. In our analysis, a diversity collapse is defined as a sharp decrease in diversity, i.e. a sudden drop of the SMST, larger than half of the temporal average of the SMST in ${10}^{4}$ generations (Appendix 1—figure 9a). In this way we excluded small or gradual decreases in a diversity measure.
In order to find a qualitative explanation for diversity collapses in ITEEM, we examined the relation between diversity and average cycle strength (Appendix 1—figure 9b). During simulations, these two quantities are usually correlated, but this correlation is blurred by the stochastic nature of ecoevolutionary dynamics. Sometimes, diversity increases faster than cycle strength or vice versa. If we highlight the time steps before the sudden collapses, we see that they always lie in the right part of the cycle strength distribution, which means that cycle strengths are larger than expected at these time points.
Size of the system
We checked the effect of system size by comparing simulations with sizes ${N}_{s}=1\times {10}^{4},3\times {10}^{4},1\times {10}^{5},3\times {10}^{5}$, $\lambda =\infty $ and a set of tradeoff parameters ($0\le \delta <1$). One of the diversity indexes (SMST) is plotted in Appendix 1—figure 10 for different sizes. With increasing ${N}_{s}$ the final diversities in trait space increase, and the drop at $\delta \approx 0.7$ becomes sharper, which is typical of a phase transition.
Plotting diversity versus size of the system for the middle rage of tradeoff parameter $\delta $ in a loglog plot reveals a scaling relation with exponent around 2: $D\sim {S}^{1.97}$ (Appendix 1—figure 11).
Frequencydependent selection
Frequencydependent selection mediated by interaction of species could be a source of temporal correlation between ecoevolutionary ‘events’, for example speciation, invasion, and extinction of species. In the absence of such biotic selections, speciation and extinction of different species occur randomly with a constant rate without any autocorrelation in time. To examine if there is such a correlation we used the distribution of interevent times, that is, the distribution of intervals between occurrence of consecutive events. For a completely random (Poisson) process – which is the null hypothesis for correlated speciation and extinction – this distribution follows an exponential distribution. Deviation from the exponential distribution is a signature of correlation between events. Appendix 1—figure 12 compares the interevent distribution of ITEEM data with the best fit of a geometric distribution (discrete version of exponential distribution) to the data. The clear deviation from the Poisson process supports that speciation and extinctions are not just random events, but that after occurrence of an event, with a delay ($\approx 10000$ generations), the probability of observing another event is higher than in a random process.
Mutation probability
Appendix 1—figure 13 shows the behavior of a typical diversity index (SMST) as function of $\delta $ and $\lambda $ for different mutation probabilities $\mu $. The overall dependency on tradeoff and lifespan is the same for a wide range of $\mu $, but the value of diversity indexes depend on it: the smaller the $\mu $ the lower the diversity in community. $\mu $ also affects the rate of increase in diversity (Appendix 1—figure 14).
Neutral model
To compare the diversity generated by genetic drift (neutral model) with diversity generated under selection pressure induced by competition with moderate tradeoffs, we simulated neutral models in which all strains compete equally for resources, that means ${I}_{\alpha \beta}=0.5$ for all pairs of strains, but traits evolve by mutation as before. In the neutral model, reproduction probabilities should also be the same for all strains, hence we attributed the same replication in each simulation to all strains. We carried out simulations for $r=0.1,0.5,0.9$. The distribution of strains over trait space (Appendix 1—figure 15) shows that genetic drift is able to spread trait vectors in trait space and to produce a small cloud of strains, but the diversity generated in this way is much smaller than that of communities evolved under biotic selection pressure mediated by competition under moderate tradeoffs (compare scale to that of bottom left panel of Appendix 1—figure 3).
In order to show clearly the difference between the diversity produced in both models we also studied diversity measures and other parameters. The to panel in each column of Appendix 1—figure 16 follows changes of a typical diversity index (SMST) over time for a simulation of genetic drift for one lifespan (with three different reproduction probability), and compares these changes with those in a typical simulation with competitive selection pressure for a moderate tradeoff ($\delta =0.56$).The bottom panels illustrate the corresponding comparison for cycle formation. The relative strengths of cycles in these simulations have large fluctuations around a value less than one, without any stable pattern over time. This means that community dynamics in the neutral model is determined by fluctuations, as expected from a model dominated by random genetic drift.
Phenotypeinteraction map
A map between two spaces, for example interaction and phenotype space, can be constructed by the rules and laws that link them. The interaction of two individuals is a function of their phenotypic traits. Generally, this can be a complex relation with different functionality of different traits. When this function is known, any phenotypic variation can be mapped into the interaction space. To generally investigate this map, we borrow the term competition kernel from adaptive dynamics theory as a phenotypebased model (Doebeli, 2011). A competition kernel measures the competitive impact of two individuals from different strains with different traits, i.e. for two individuals from strains $\alpha $ and $\beta $, the competition kernel is a function $a\left({\mathbf{x}}_{\alpha},{\mathbf{x}}_{\beta}\right)$ where ${x}_{\alpha}=({x}_{1\alpha},{x}_{2\alpha},\dots )$ and ${x}_{\beta}=({x}_{1\beta},{x}_{2\beta},\dots )$ are the phenotypic trait vectors of the corresponding strains. Elements of these vectors can be any relevant phenotype like size, color, expression of a gene, etc. Considering that in our model the carrying capacity is fixed and equal for all traits, interaction terms of ITEEM are proportional to the competition kernel, ${I}_{\alpha \beta}\propto a({x}_{\alpha},{x}_{\beta})$.
The dimension of the phenotype space is equal to the number of traits with each axis representing a trait, and strains are distributed according to their traits over this space (Appendix 1—figure 17a). Interaction space, on the other hand, has one axis for each strain, and individuals are distributed based on their interactions with the strains that represent the axes (Appendix 1—figure 17b). The dimension of interaction space is dynamic because it increases with new strains and shrinks as strains go extinct.
A new, phenotypical mutant strain appears in phenotype space close to its parent (Appendix 1—figure 18a and c). Any phenotype variation that is not ecologically neutral and produces a new ‘interaction’ mutant, i.e. a strain with a novel interaction vector, adds a new dimension to the interaction space (Appendix 1—figure 18 b and d) while it is still close to its parent, as expected from their ecological similarity.
A relevant question about interaction basedmodels could be if the random Gaussian variations in the interaction traits are biologically meaningful or not. The key to answering this question is the competition kernel, that is, how phenotypes are translated to the interactions, and hence depends on its functionality. However, if we consider that a mutant should be ecologically similar to its parent, the assumption of random Gaussian variations appears justified; it is also in line with models at the phenotype and genotype levels. Evolution is exploring the trait space, which, depending on the model, could be genotype, phenotype or interaction space. This exploration can be random (neutral evolution), or directional (adaptive evolution). In evolutionary models in phenotype space, usually phenotypic variations are drawn from a Gaussian distribution around the parent’s trait vector. The fate of these mutants are then determined by either genetic drift, or the fitness landscape, or the community, i.e. when frequencydependent selection is the driving evolutionary force. If we use such Gaussian distributions to model phenotype traits, we can use the competition kernel to map them into interaction space and test if such variations produce random variations in the interaction traits.
Consider a system with $P$ phenotypic traits in which each strain $\alpha $ is described by its trait vector ${x}_{\alpha}=({x}_{1},{x}_{2},\dots ,{x}_{P})$. A mutant offspring ${\alpha}^{\prime}$ should have a trait vector that is a random variation of the trait vector of its parent $\alpha $, i.e. ${\mathbf{x}}_{{\alpha}^{\mathrm{\prime}}}={\mathbf{x}}_{\alpha}+\mathit{\nu}=\left({x}_{1}+{\nu}_{1},{x}_{2}+{\nu}_{2},\dots ,{x}_{n}+{\nu}_{P}\right)$, where elements ${\nu}_{i}$ ($i=1,\dots ,P$) of $\mathit{\nu}$ are drawn independently from a normal distribution. To map this system to the interaction space, we need the competition kernel. In a system with $N$ strains, the interaction trait vector of parent strain $\alpha $ is ${T}_{\alpha}=({I}_{\alpha 1},\dots ,{I}_{\alpha N})=\left(a\right({x}_{\alpha},{x}_{1}),\dots ,a({x}_{\alpha},{x}_{N}\left)\right)$ while the interaction trait vector of mutant ${\alpha}^{\prime}$ is ${\mathbf{T}}_{{\alpha}^{\mathrm{\prime}}}=\left({I}_{{\alpha}^{\mathrm{\prime}}1},\dots ,{I}_{{\alpha}^{\mathrm{\prime}}N}\right)=\left(a\left({\mathbf{x}}_{{\alpha}^{\mathrm{\prime}}},{\mathbf{x}}_{1}\right),\dots ,a\left({\mathbf{x}}_{{\alpha}^{\mathrm{\prime}}},{\mathbf{x}}_{N}\right)\right)=\left(a\left({\mathbf{x}}_{\alpha}+\mathit{\nu},{\mathbf{x}}_{1}\right),\dots ,a\left({\mathbf{x}}_{\alpha}+\mathit{\nu},{\mathbf{x}}_{N}\right)\right)$. As the phenotypic variations are small, we can approximate the $\beta $th element of ${T}_{{\alpha}^{\prime}}$ with a Taylor series expansion:
$g={\scriptscriptstyle \frac{d}{d{x}_{\alpha}}}a({x}_{\alpha},{x}_{\beta})$ and $H={\scriptscriptstyle \frac{d}{d{x}_{\alpha}}}g$ are the gradient vector and the Hessian matrix (first and second order derivatives in single variable case) of the competition kernel and ${\eta}_{\beta}$ is a random number with normal distribution (Equation 1 in the main text). For the last equation we have used the central limit theorem, which states that the sum of many independent random variables is approximated well by a normal distribution.
Despite the fact that the precise influence of random phenotypic variations on the interaction trait depends on the functionality of competition kernel, the above approximation (Equation 10) shows that those variations can be mapped to random variations in the interaction terms. This means that if a mutant emerges with random variation of its parent’s phenotype, its interactions with the extant strains are random variations of the interactions of the parent. It is important to mention that neither in phenotype nor in interaction space, random trait variation between parent and mutant offspring leads necessarily to random independent characters of parent and offspring. The latter is true only if evolution is governed by neutral drift. In adaptive evolution, the fate of a mutant is determined by the interaction of that mutant with the community or the environment.
One important aspect of the interaction level modeling is that the organisms are defined in this space by their interaction traits. Thus, variations that occur in different phenotype traits are coarsegrained into the interaction terms. Interaction space is not restricted by the dimensions of the phenotype trait and hence allows for evolutionary innovations that happen due to emergence of new phenotypic dimensions, for example a novel metabolic pathway activated by epigenetic changes. However, this coarsegraining neglects phenotypic variations that do not affect ecological interactions, and thus maps different phenotypic mutations that lead to the same effective ecological interactions to the same interaction term, for example if those phenotypic variations yield the same ecological dominance (Appendix 1—figure 19). This is similar to the genotypephenotype map if several genotypes are mapped to the same phenotype. The noninjectivity between the phenotype and interaction space is not an issue when the ecological outcome of the ecoevolutionary dynamics is studied.
Data availability
The source code of the model is freely available at https://github.com/BioinformaticsBiophysicsUDE/ITEEM; copy archived at https://github.com/elifesciencespublications/ITEEM).
References

BookTradeoffs and negative correlations in evolutionary ecologyIn: Bell M. A, Futuyama D. J, Eanes W. F, Levinton J. S, editors. Evolution Since Darwin: The First 150 Years. Sinauer Associates, Inc. pp. 243–268.

Network properties, species abundance and evolution in a model of evolutionary ecologyJournal of Theoretical Biology 232:551–558.https://doi.org/10.1016/j.jtbi.2004.03.029

Stochastic formulation of ecological models and their applicationsTrends in Ecology & Evolution 27:337–345.https://doi.org/10.1016/j.tree.2012.01.014

Sympatric speciation: models and empirical evidenceAnnual Review of Ecology, Evolution, and Systematics 38:459–487.https://doi.org/10.1146/annurev.ecolsys.38.091206.095804

BookModern Multidimensional Scaling: Theory and ApplicationsSpringer Science & Business Media.

Competitive networks: nontransitive competitive relationships in cryptic coral reef environmentsThe American Naturalist 113:223–234.https://doi.org/10.1086/283381

Beyond species: functional diversity and the maintenance of ecological processes and servicesJournal of Applied Ecology 48:1079–1087.https://doi.org/10.1111/j.13652664.2011.02048.x

BookCommunities in patchy environments: a model of disturbance, competition and heterogeneityIn: Jurek Kolasa S. T. P, editors. Ecological Heterogeneity. Springer Verlag.. pp. 97–122.

Comparing classical community models: theoretical consequences for patterns of diversityThe American Naturalist 159:1–23.https://doi.org/10.1086/324112

Mechanisms of maintenance of species diversityAnnual Review of Ecology and Systematics 31:343–366.https://doi.org/10.1146/annurev.ecolsys.31.1.343

TradeOff geometries and FrequencyDependent selectionThe American Naturalist 164:765–778.https://doi.org/10.1086/424762

Individualbased models in ecology after four decadesF1000Prime Reports 6:39.https://doi.org/10.12703/P639

The Handbook of Evolutionary PsychologyLife history theory and evolutionary psychology, The Handbook of Evolutionary Psychology, Wiley.

Diversity and coevolutionary dynamics in HighDimensional phenotype spacesThe American Naturalist 189:105–120.https://doi.org/10.1086/689891

A quantitative genetic competition model for sympatric speciationJournal of Evolutionary Biology 9:893–909.https://doi.org/10.1046/j.14209101.1996.9060893.x

The influence of predatorprey population dynamics on the longterm evolution of food web structureJournal of Theoretical Biology 208:91–107.https://doi.org/10.1006/jtbi.2000.2203

Microbiology and the species problemBiology & Philosophy 25:553–568.https://doi.org/10.1007/s1053901092119

SoftwareMuller Plot: Muller Plots from Population/Abundance/Frequency Dynamics Data, version 0.1.2R Package Version.

Tradeoff mechanisms shaping the diversity of BacteriaTrends in Microbiology 24:209–223.https://doi.org/10.1016/j.tim.2015.11.009

Ecoevolutionary feedbacks, adaptive dynamics and evolutionary rescue theoryPhilosophical Transactions of the Royal Society B: Biological Sciences 368:20120081.https://doi.org/10.1098/rstb.2012.0081

The effects of intransitive competition on coexistenceEcology Letters 20:791–800.https://doi.org/10.1111/ele.12775

Evolution of community structure: CompetitionJournal of Theoretical Biology 133:513–523.https://doi.org/10.1016/S00225193(88)803382

Constraints on microbial metabolism drive evolutionary diversification in homogeneous environmentsJournal of Evolutionary Biology 20:1882–1889.https://doi.org/10.1111/j.14209101.2007.01376.x

Bacterial competition: surviving and thriving in the microbial jungleNature Reviews Microbiology 8:15–25.https://doi.org/10.1038/nrmicro2259

Biodiversity and Ecosystem Functioning: Synthesis and Perspectives195–208, Species diversity, functional diversity and ecosystem functioning, Biodiversity and Ecosystem Functioning: Synthesis and Perspectives, Oxford University Press on Demand.

BookThe Unified Neutral Theory of Biodiversity and BiogeographyPrinceton University Press.

Homage to Santa Rosalia or why are there so many kinds of animals?The American Naturalist 93:145–159.https://doi.org/10.1086/282070

The approach to mutationselection balance in an infinite asexual population, and the evolution of mutation ratesProceedings of the Royal Society B: Biological Sciences 266:2389–2397.https://doi.org/10.1098/rspb.1999.0936

Symmetry of interactions rules in incompletely connected random replicator ecosystemsThe European Physical Journal E 37:56.https://doi.org/10.1140/epje/i2014140567

Evolutionary branching under asymmetric competitionJournal of Theoretical Biology 197:149–162.https://doi.org/10.1006/jtbi.1998.0864

Coexistence in a variable environment: Ecoevolutionary perspectivesJournal of Theoretical Biology 339:14–25.https://doi.org/10.1016/j.jtbi.2013.05.005

Competitive intransitivity, population interaction structure, and strategy coexistenceJournal of Theoretical Biology 365:149–158.https://doi.org/10.1016/j.jtbi.2014.10.010

Optimal Life Histories Under AgeSpecific PredationThe American Naturalist 114:399–417.https://doi.org/10.1086/283488

BookEvolution in Changing Environments: Some Theoretical ExplorationsPrinceton University Press.

The form of a tradeoff determines the response to competitionEcology Letters 16:1267–1276.https://doi.org/10.1111/ele.12159

Alternative designs and the evolution of functional diversityThe American Naturalist 167:55–66.https://doi.org/10.1086/498276

A guide for using functional diversity indices to reveal changes in assembly processes along ecological gradientsJournal of Vegetation Science 24:794–806.https://doi.org/10.1111/jvs.12013

Ecosystems with mutually exclusive interactions selforganize to a state of high diversityPhysical Review Letters 107:188101.https://doi.org/10.1103/PhysRevLett.107.188101

Runaway evolution to selfextinction under asymmetrical competitionEvolution 48:1764–1772.

Diversity begets diversity in competition for space. NatEcology and evolution 1:0156.

Emergence of diversity in a model ecosystemPhysical Review E 86:011929.https://doi.org/10.1103/PhysRevE.86.011929

BookAdvances in Ecological Research: EcoEvolutionary Dynamics, Volume 50Academic Press.

How to get the most bang for your buck: the evolution and physiology of nutritiondependent resource allocation strategiesProceedings of the Royal Society B: Biological Sciences 284:20170445.https://doi.org/10.1098/rspb.2017.0445

Metabolic TradeOffs promote diversity in a model ecosystemPhysical Review Letters 118:1–5.https://doi.org/10.1103/PhysRevLett.118.028103

Ecoevolutionary feedbacks in community and ecosystem ecology: interactions between the ecological theatre and the evolutionary playPhilosophical Transactions of the Royal Society B: Biological Sciences 364:1629–1640.https://doi.org/10.1098/rstb.2009.0012

The Species Problem  Ongoing IssuesThe species problem: A conceptual problem?, The Species Problem  Ongoing Issues.

The species concept for prokaryotesFEMS Microbiology Reviews 25:39–67.https://doi.org/10.1111/j.15746976.2001.tb00571.x

What can interaction webs tell us about species roles?PLOS Computational Biology 11:e1004330.https://doi.org/10.1371/journal.pcbi.1004330

Emergence of structured communities through evolutionary dynamicsJournal of Theoretical Biology 383:138–144.https://doi.org/10.1016/j.jtbi.2015.07.020

Microbial diversityproductivity relationships in aquatic ecosystemsFEMS Microbiology Ecology 62:181–186.https://doi.org/10.1111/j.15746941.2007.00381.x

Biological Evolution and Statistical Physics312–337, Modelling macroevolutionary patterns: An ecological perspective, Biological Evolution and Statistical Physics, Berlin Heidelberg, Springer.

BookCompetition and Coexistence, 161Springer Science & Business Media.https://doi.org/10.1007/9783642561665

Cyclic dominance in evolutionary games: a reviewJournal of the Royal Society Interface 11:20140735.https://doi.org/10.1098/rsif.2014.0735

Rapid evolution as an ecological processTrends in Ecology & Evolution 13:329–332.https://doi.org/10.1016/S01695347(98)013780

The evolution of species interactionsScience 284:2116–2118.https://doi.org/10.1126/science.284.5423.2116

Interactions among plants and evolutionJournal of Ecology 99:729–740.https://doi.org/10.1111/j.13652745.2011.01802.x

Resource competition and community structureMonographs in Population Biology 17:1–296.

BookSpecialization, Speciation, and Radiation: The Evolutionary Biology of Herbivorous InsectsUniversity of California Press.

Emergence of a complex and stable network in a model ecosystem with extinction and mutationTheoretical Population Biology 63:131–146.https://doi.org/10.1016/S00405809(02)000382

Selforganized spatial pattern determines biodiversity in spatial competitionJournal of Theoretical Biology 300:48–56.https://doi.org/10.1016/j.jtbi.2012.01.005

Ecoevolutionary model of rapid phenotypic diversification in SpeciesRich communitiesPLOS Computational Biology 12:e1005139.https://doi.org/10.1371/journal.pcbi.1005139

Variations and Fluctuations of the Number of Individuals in Animal Species living togetherICES Journal of Marine Science 3:3–51.https://doi.org/10.1093/icesjms/3.1.3

Geometric Structure of HighDimensional Data and Dimensionality Reduction115–129, Classical multidimensional scaling, Geometric Structure of HighDimensional Data and Dimensionality Reduction, Springer.

Evolution in a community context: on integrating ecological interactions and macroevolutionTrends in Ecology & Evolution 32:291–304.https://doi.org/10.1016/j.tree.2017.01.003

Evolutionary dynamics of species diversity in an interaction web systemEcological Modelling 163:131–143.https://doi.org/10.1016/S03043800(02)004179
Decision letter

Michael DoebeliReviewing Editor; University of British Columbia, Canada

Ian T BaldwinSenior Editor; Max Planck Institute for Chemical Ecology, Germany

Yaroslav IspolatovReviewer; Universidad de Santiago de Chile, Chile
In the interests of transparency, eLife includes the editorial decision letter and accompanying author responses. A lightly edited version of the letter sent to the authors after peer review is shown, indicating the most substantive concerns; minor comments are not usually included.
[Editors’ note: a previous version of this study was rejected after peer review, but the authors submitted for reconsideration. The first decision letter after peer review is shown below.]
Thank you for submitting your work entitled "Tradeoff shapes diversity in ecoevolutionary dynamics" for consideration by eLife. Your article has been reviewed by four peer reviewers, one of whom is a member of our Board of Reviewing Editors, and the evaluation has been overseen by a Reviewing Editor and a Senior Editor. The following individuals involved in review of your submission have agreed to reveal their identity: Yaroslav Ispolatov (Reviewer #3).
Our decision has been reached after consultation between the reviewers. Based on these discussions and the individual reviews below, we regret to inform you that we have decided to reject your paper. However, we would like to offer the possibility of resubmission of a substantially revised manuscript. Our detailed assessment is as follows.
The paper has been carefully read by four expert reviewers. All of the reviewers found your paper thought provoking and potentially interesting. However, three reviewers raised substantial concerns. Two reviewers question the biological relevance of the models considered in the paper, and they think the model is too abstract to be useful for biologists thinking about longterm evolution. Since eLife is a "generalist" journal, it is imperative that papers have a high relevance not only for a narrow theoretical audience, but also for a general audience, in this case primarily in ecology and evolution. The reviewers also raised concerns that the work was not put into the context of existing theoretical work. In particular, one reviewer was concerned that a very similar approach has been taken in a previous paper, which is cited in the manuscript, but not discussed at all. All these concerns would need to be addressed in a revision. Specifically, the authors need to convince a mainly biological audience that the work presented is interesting and relevant, and they need to be clear about what advances the current work represents over previous work. Consequently, a serious effort would be required for revising the paper, and there is no guarantee for a positive outcome. Should you choose to resubmit a revised version, you would also need to take into account all the other reviewer comments in a constructive manner.
If you feel you can make the necessary revisions, we will make every effort to return a new version of the manuscript to the same editors and reviewers. Please note that this would be treated as a new submission.
Reviewer #1:
This is an interesting paper presenting a "minimal" model for the evolution of diversity in communities governed by competitive interactions. The model is minimal in the sense that only the outcome of interactions is modelled, but not their mechanistic basis, which greatly simplifies the description.
I have two main concerns, and a number of more technical points.
Essential revisions:
1) The work presented appears to be very similar to the theory presented in Shtilerman et al., (2015). The authors need to clarify to what extent their work is novel compared to that earlier paper.
2) In fact, I think the logistic equations presented by the authors in the appendix are similar to the ones used in Shtilerman et al., 2015. I think the authors need to consider using these equations, rather than their individualbased simulations, to derive their results. It appears to be straightforward to implement the evolutionary dynamics based on the deterministic logistic equations, with one differential equation per species, and with evolution being implemented by adding new equations to the system (and deleting those ones in which the frequency falls below a threshold). This would yield a computationally more efficient model, whose salient features would nevertheless be essentially the same as in the individualbased model (perhaps barring the nearly neutral variation within species clusters). If the results are not the same, then this needs to be known, as it would probably point to an important role of stochasticity. I think this line of analysis is necessary to obtain a complete picture of the system studied by the authors, which is required for eLife.
Essential revisions:
– Subsection “Model”: the procedure for the individualbased model is not clear: the text states that "we try 𝑁𝑖𝑛𝑑(𝑡) (number of individuals at that time step) replications of randomly selected individuals. Each selected individual of a strain 𝛼 can replicate with rate r_α…" this is imprecise; what does is it mean that an individual that has already been selected for replication replicates with a rare r_\α? Isn't it rather the case that individuals get selected to produce an offspring with a probability that is proportional to r_α?
– Also, it is not clear how the discreteness of time is dealt with exactly. For example, the paper says that each offspring is assigned a randomly selected site, then competes against the occupant of the site and probabilistically takes over. What if one offspring takes over, but then a second offspring is chosen (by chance) to compete for the same site? Is that second offspring then competing against the previous occupant, or against the new occupant, i.e., the offspring that has already taken over from the previous occupant?
– And when do individuals die? Before or after offspring production?
– The authors mention the existence of rockpaperscissor interactions, but it is not at all clear what the importance of such interactions are for the overall dynamics of the system. Are such interactions very common, essential, just a quirk?
– Subsection “Impact of tradeoff and lifespan on diversity”: it is not clear what "average cycle strength" is. Which cycles were used to calculate the average of what?
– Subsection “Impact of tradeoff and lifespan on diversity”: shouldn't the average interaction strength, i.e. the average values of all I's, be ca. 0.5?
– I failed to understand Figure 2. Please explain better.
– Subsection “Impact of tradeoff and lifespan on diversity”: the Darwinian demon doesn't really make sense: with no diversification, there is only one strain present, and it is impossible to say what kind of competitive ability this strain has, because no other strain is present… so what are you actually trying to say?
– I don't find the section about the relationship between tradeoff and productivity very convincing, because it is too vague. I think the authors need to be more specific and detailed here.
– Section "Frequencydependent selection": I think a more mechanistic explanation for differences of speciation and extinction rates compared to Poisson processes is needed. What is the cause of these differences? When is this mechanism at play, and when not? In particular, why do speciation rates increase with the diversity in the system? Shouldn't this only happen at certain intermediate levels of diversity, because when the community reaches a highlevel diversity, speciation rates presumably decline?
– I would like to bring the authors' attention to the recent paper by Doebeli and Ispolatov, (2017), which also presents models for community assembly and evolution of diversity are presented, with some parallels to the work discussed here. These parallels could be taken up in the Discussion section.
Reviewer #2:
The authors present a model for studying diversification dynamics as an outcome of, as they claim, competition for a single resource. I have substantial concerns regarding the insight that could be gained from this work, mainly due to the type of model and the lack of a mechanistic interpretation.
Subsection “Model” and subsection “Frequentdependent selection”: Your trait space is Nspdimensional, where Nsp is the number of extant strains. In particular, the trait space grows explicitly with every strain added. Is it surprising that you get coexistence of Nsp strains in an Nspdimensional trait space? Further, is it surprising that "emergence of new species increases the probability for generation of further species"? It seems to me that this is quite expected, given your model. This also explains why you get these mass extinction events; there's a positive feedback loop between speciation rate and species number that is easy to recognize even without simulations.
Subsection “Generation of Diversity”: Please explain how you define "functional diversity" in your model at first mention. Appendix 6 was uninformative. The only information I found was in the caption of Figure 1, where you define functional diversity "in terms of the size of minimum spanning tree (SMST) in trait space". But your trait space is continuous ([0,1]Subsection “Generation of Diversity”: Please explain how you define "functional diversity" in your model at first mention. Appendix 6 was uninformative. The only information I found was in the caption of Figure 1, where you define functional diversity "in terms of the size of minimum spanning tree (SMST) in trait space". But your trait space is continuous ([0,1]Subsection “Generation of Diversity”: Please explain how you define "functional diversity" in your model at first mention. Appendix 6 was uninformative. The only information I found was in the caption of Figure 1, where you define functional diversity "in terms of the size of minimum spanning tree (SMST) in trait space". But your trait space is continuous ([0,1]Subsection “Generation of Diversity”: Please explain how you define "functional diversity" in your model at first mention. Appendix 6 was uninformative. The only information I found was in the caption of Figure 1, where you define functional diversity "in terms of the size of minimum spanning tree (SMST) in trait space". But your trait space is continuous ([0,1]^Nsp) so a priori there is no tree structure connecting species (apart from phylogeny). So how do you define a spanning tree, and how is this spanning tree not a function purely of the number of strains?
Subsection “Generation of Diversity”: Related to the previous point. What is a "functional niche" in your model? Since you did not discuss mechanisms underlying the interaction matrix, it is not clear what a function is.
Subsection “Generation of Diversity”: Saying your model has no geographical isolation nor resource partitioning sounds meaningless. You have not specified the underlying mechanisms for your "interaction trait", so it is hard to draw a comparison to other more mechanistic models (i.e. where the physiological/metabolic traits are explicit).
Subsection “Model”: Since replication is nonsexual, what is the purpose of defining "species" separately from "strain" in your model? In AppendixI you mention that you define "species" operationally based on divergence time and using some arbitrary cutoff threshold, so this sounds a bit analogous to operational taxonomic units (OTUs) in microbiology. However, you model mutation/selection/growth dynamics at the strain level. Please explain early on what additional insight one may get from counting "species".
Discussion section: You claim that modeling in terms of interaction traits (which in your case means in terms of outcome probabilities of local competitive exclusion) "coarsegrains these complex systems in a natural, biologically meaningful way.". But many of your interpretations and comparisons to other models or data require translating your Abstract "interaction traits" to functions or life histories. You did not discuss at all what real traits could possibly give rise to your interaction trait matrix. Your interaction matrix seems to loosely represent competitive interactions; but then it can only explain diversity within a single trophic level (e.g. in the case of animals) or a single metabolic niche (e.g. in the case of bacteria); yet, you keep referring to "functional diversity".
Discussion section: Related to the previous comment. You say that your formulation was "chosen to reflect reasonable properties" and that you "have assumed a single, limiting resource in a wellmixed system". However, you did not provide any plausible mechanism for how such an interaction matrix could arise from competition for a single resource in a wellmixed system. This is actually a big question: how can one obtain a Nspdimensional trait space through competition for a single resource pool?
Reviewer #3:
The manuscript "Tradeoff shapes diversity in ecoevolutionary dynamics" presents a quite original evolutionaryecological individualbased model based on a compromise between growth rate and competitiveness. The model is simple but exhibits quite rich evolutionary properties, which, while being not entirely unpredictable, are intriguing and provide useful insights and generalizations. This combination of simplicity of the rules and relevance of the dynamics usually distinguishes successful and longliving models from the rest and makes them understandable and appealing to a wide audience. Besides, I cannot see any potentially "hidden" flaws in the model and interpretation of results that could cast doubts on the main conclusions. Thus, in my opinion, the manuscript may be published in eLife after the following and, perhaps, other comments are addressed:
1) The complexity of algebra in (3) and subsequent definition of s is definitely unwarranted by an otherwise quite accessible level of math in the rest of the paper. Furthermore, it apparently confuses even the Authors when they classify the domains of weak and strong tradeoffs:
It looks like the definitions of hightradeoff (\δ ~ 1) and intermediate tradeoff (\δ ~ 1/2) are misleading. The strongest dependence of r on C appears to be when \δ=1/2, which also follows from Figure 1. The phrase "In hightradeoff phase III, any small change in C changes r drastically" is simply wrong for all but very small C. I would call both the \δ=1 and \δ=0 limits as small tradeoffs as they look perfectly symmetric, or, even better, choose another, more heuristic and transparent, form of parametrization of r vs. C.
2) Would it be possible to say anything about population dynamics of "species"? Is it cyclic?
3) At least qualitatively, what are mechanisms of mass extinction?
4) An explicit plot of diversity vs. system size (perhaps just for the δ optimal for diversity) and, ideally, an estimate of corresponding scaling, would be very revealing.
5) Similarly, how the level of diversity and the typical number of traits in a species depends on the mutation amplitude m?
6) For a general reader of eLife, the MDS algorithm needs to be explained and properly referenced. It plays a major role in interpreting the results, however, is presented only by the name of the corresponding function in R.
7) A qualitative explanation about the minimal spanning tree analysis would help as well.
Reviewer #4:
I've carefully read the paper "Tradeoff shapes diversity in ecoevolutionary dynamics" by Farnoush Farahpour and colleagues. In it, they use an individualbased ecoevolutionary model to understand the emergence of community structure based on evolving interactions. The model produces some interesting patterns, such as clustering in traitspace and the evolution of intransitive competitive loops.
There are a number of aspects of this paper that I liked: the focus on evolution of interactions, the emergence of intransitive loops, the dimension reduction applied to the vectorvalued traits, and the comparison with a neutral model variant. These are all creative contributions to the modeling literature.
While I was intrigued by the idea of modeling the evolution of interactions directly, it was hard for me to connect it to real ecological systems. The interactions between species are not determined by phenotypic traits of the organisms but evolve independently. It's based on an unstated assumption that species interactions are totally idiosyncratic and unpredictable. I feel that evolution in this model is too unconstrained, despite the tradeoff between competitive ability and reproduction, resulting in the prevalence of intransitive loops. As a complexsystems researcher, I'm fascinated, but as an evolutionary ecologist, I'm skeptical.
The authors need to do a better job putting their work in the context of the extensive literature on ecoevolutionary dynamics. The text had many statements that had me scratching my head. Examples:
1) Framing the problem in terms of resources, but in the model there are no actual resources.
2) The competitive exclusion principle only holds at equilibrium (Armstrong, Levins) and must count resources plus shared predators (Levin, 1970).
3) The "ecoevolutionary models" cited in the Introduction seem to be just ecological models.
4) In the Introduction, the "observed ecoevolutionary dynamics" cited that the model "closely resembles" aren't empirical patterns observed in real systems, but just results of other models.
5) Discussion of speciation overlooks that the species here are clonal, so what's hard about speciation?
Also, some of the sentences throughout were hard to understand the meaning of. E.g., "Evolutionary changes at the genetic level influence ecology if they cause phenotypic variations that affect biotic or abiotic interactions of species which in turn changes the species composition and occasionally forces species to evolve their strategies."
Some of the details of the model implementation weren't clear. For example:
How exactly do births happen (subsection “Model”)?
Is mu a mutation "rate" or probability of mutation during a replication event (subsection “Model”)?
Why is lifespan drawn from a Poisson distribution (Subsection “Model”) and how can that be infinite (Figure 2)?
If each individual stays in a site, is it really wellmixed?
Does mutation of one species' interaction coefficients end up changing another species' reproduction rate through the tradeoff (2)?
Could you not get at the same questions more efficiently using deterministic LotkaVolterra dynamics?
[Editors’ note: what now follows is the decision letter after the authors submitted for further consideration.]
Thank you for resubmitting your work entitled "Tradeoff shapes diversity in ecoevolutionary dynamics" for further consideration at eLife. Your revised article has been evaluated by Ian Baldwin (Senior editor), a Reviewing editor, and 4 reviewers (the same reviewers as for the initial submission, reports below).
The reviewers agree that substantial effort was put into the revision. However, there are a number of remaining issues that call for a further indepth revision. Reviewers 1, 2 and 3 (who was previously reviewer 4) still have substantial concerns about aspects of the paper. These concerns are all related to issues raised in the initial reports and will need to be addressed in a constructive manner if the paper is to be published in eLife.
Essential revisions:
1) Individualbased vs deterministic systems: both reviewer 1 and 4 (formerly reviewer 3) point out that it is not at all clear that these different approaches would yield different results. In other words, it is not clear that salient results reported in the paper depend crucially on the presence of noise, as appears to be the contention of the authors. This needs to be explored at least to some extent. This issue cannot be dealt with by simply "citing it away".
2) The assumption of independent evolution of the elements of the interaction matrix needs to be discussed in more detail and clarity, and in the context of biological realism. This point was raised by both reviewers 1 and 3 (formerly reviewer 4).
3) The "single resource" issue is related to point 2, and also needs a revised treatment.
4) I am sympathetic with reviewer 3's (formerly reviewer 4) concern about some of the references to previous work, particularly with regard to tradeoffs.
5) Please also address the remaining points, e.g. the definition of "species" raised by reviewer 2, as well as other minor points.
Reviewer #1:
The authors did a very good and thorough job in revising their paper. Almost all of my concerns have been addressed satisfactorily. The one remaining issue is that I don't buy the authors case for only using individualbased simulations (my original comment 2). The authors state in their rebuttal that "It has been extensively discussed in the literature that continuum approaches are unsuitable in cases of nonequilibrium dynamics,.…", and they somehow conclude from this that the deterministic LotkaVolterra description (which in my mind is the same as a "logistic" description) would not be appropriate for the problem at hand. But that is exactly the question: in some general sense, one would expect that in the limit of large population sizes, the individualbased model used by the authors would converge to some deterministic model, and my guess is that this model would be at least close to the "meanfield" LotkaVolterra model. It then becomes important to understand just which features of the individualbased models can be understood by studying the much simpler meanfield model. The authors present no arguments why there even are *any* features of their individualbased model that could not be observed in the deterministic model. They refer to nearneutrality, but that's exactly one of the features that was observed in similar deterministic LotkaVolterra models by Shtilerman et al., (2015). I am not convinced that the salient results reported in this paper could not also be obtained with deterministic models. Just claiming that some results cannot be obtained in that way on general grounds does not make it true in this particular case. The obvious advantage of using deterministic models would be that such models are much more tractable analytically (e.g. from a statistical physics point of view), and it is therefore important to know how far one can get using them. I think it would not be too onerous to at least do some tests using the deterministic models to either confirm or refute the claim that they can produce similar results as the individualbased models. The question is: does stochasticity really play a major role in producing the results reported in this paper? If so, then this would be important to know (but this point would need to be made based on more than just a vague statement that their "system falls into the category of those better modeled by individual based models"). If not, then the deterministic models should do a good job reproducing these results.
Reviewer #2:
The authors have partly addressed my concerns and those of the other reviewers. I do however still have two major concerns:
1) I agree that a lowdimensional phenotype space (e.g. pertinent to exploitation of/competition for a single resource) can give rise to an Nst x Nst interaction matrix that encodes the competitive interaction between strains. However, the crucial assumption of the authors is that each term in this matrix (well, half of the terms) can vary independently. How this could come about in reality is unclear to me.
In other words, if "P" is the underlying phenotype space (solely related to consumption of/competition for the common resource) and "I" is the space of possible interaction matrices, what could the mapping f:P>I possibly look like, such that f(P) is an Nst x Nst/2 dimensional manifold?
I strongly recommend:
(a) Avoid any comparison to "single resource" models or real systems.
(b) Acknowledge early on that an important assumption of the model is that the terms of the interaction matrix (well, half of them) can in principle vary independently (i.e. are not constrained explicitly due to genetics or ecology). Whether this assumption is met in reality is an open question.
(c) Clarify in the discussion that this paper does not address the important question of how such a highdimensional interaction trait space (i.e. with Nst x Nst/2 independent axes) might arise, or provide a plausible example.
2) The definition of "species" by the authors is still confusing and of questionable relevance. The authors define "species" operationally based on a cutoff threshold in phylogenetic distance. While this is common practice in microbial ecology (where such clusters are called Operational Taxonomic Units), few would claim that the emergence and disappearance of OTUs is comparable to "speciation" dynamics in sexually reproducing organisms.
What I also found confusing is that in their "response to reviewers", the authors explain that "species" are "wellseparated nontransient clusters in trait space". This does not align with the definition provided in their manuscript (Appendix 1), where species are defined as "clusters of strains separated by longlasting gaps in a phylogenetic tree". Are these definitions equivalent in your model?
While the emergence of clusters in trait space is indeed interesting, I would recommend not calling these clusters "species", since clusters in trait space need not always be monophyletic and could in principle also consist of distantly related strains that happen to have converged in trait space.
Reviewer #3:
This is reviewer 4 from the original submission again. This remains an interesting yet frustrating manuscript. The authors resisted many of the good suggestions from the other reviewers and myself in how they can place their manuscript in the broader context. In the end, it's the authors' manuscript, but I still think they could do a better job in the introduction and discussion to not confuse potential readers.
To me, the most interesting part of the manuscript is the idea that species interactions might be so high dimensional that it is best to focus on interaction traits that summarize many idiosyncratic phenotypes. This is described in the discussion but should also be highlighted more in the Introduction. The relationship between phenotypic traits and interaction traits should be clarified to better address comments 1 and 67 of reviewer 2. Maybe could be described as a rugged "phenotypeinteraction map", in analogy to the idea of "genotypephenotype maps"? By the way, this is a big assumption of the model, not an established empirical fact, but still an interesting basis for the theory.
Reviewer 2 and I were confused by statements about limiting resources and the competitive exclusion principle. In the revision, the authors still make statements like "GLV equations model competition over renewable resources" (Subsection “Model”), "we observe high diversity in a wellmixed homogenous system without violating the competitive exclusion principle" (subsection “Tradeoff anchors ecoevolutionary dynamics in physical reality”) and "we have assumed a single, limiting resource" (Subsection “Power and limitations of ITEEM”). Such statements will misdirect many readers into thinking about resource competition, R* rules, and the impossibility of coexistence of more species than resources. This is not appropriate, because in this model, the species interactions are direct (interference competition) and idiosyncratic. Allelopathy among plants or microbes would be a more relevant example than resource competition. The authors should remove all mentions about resources in the paper, because they will only confuse readers.
Concerning my previous comment 4, please don't portray other theoretical results as empirical support for your new model. Keep them clearly distinct.
It's a big stretch to say claim "lifehistory tradeoffs" are a missing ingredient in existing theory (Introduction). Almost all existing ecoevolutionary theory is built around tradeoffs. Models without tradeoffs are the exception, not the rule.
The references cited in subsection “Model”, do not represent current, mainstream ecological thinking.
Reviewer #4:
I think that the authors addressed the main points from my previous review, just a few minor issues remain.
I still disagree with the definition of “strong' and “weak' tradeoff limits, seeing, for example, plots in Figure 1 as completely symmetric. In my opinion, the strongest dependence of birthrate on competitiveness happens at the line in the middle of the plot, presumably for δ=0.5. I guess it's more a terminological discussion, however, I find the definition of the strong tradeoff adopted by the Authors rather confusing.
Caption to Figure 2, “Disc diameter scales with total abundance of species' Does it mean that it scales with the number of individuals in a species? Or the number of species in the system? What kind of scaling is it?
Subsection “Generation of diversity”, “Occasionally diversity collapses from medium levels abruptly to very low levels, usually followed a recovery”, should it read “by a recovery”?
Appendix 1 (Eq.8) Is γ the index of summation, running from one to N_st? If not, what is the index?
I don't know if I should dwell on that, yet I also strongly disagree with the Authors' reply to the second comment of the first reviewer, and especially, with the apparent misuse of the term “mean field'.
In short, I believe that both implementations of this process, the individualbased and continuouspopulations models should yield very similar phenomenology. Both those implementations are mean field in their nature as neither has any spatial correlations (in phenotypic or geographical space) or long temporal memory. The only difference between those is the presence of some stochasticity or noise in the individualbased model. If the Authors truly believe that such noise is the necessary source of the observed phenomenology, it should be clearly stated in the manuscript. Which, I believe, would have strongly depreciated the generality of conclusions. However, I don't think this is the case; on contrary, a continuous population version of this model would have enabled one to get “cleaner' results, speeding up the simulations and expanding the scaling range by including more species. The main distinct features that the model develops, such as temporal changes in the population, interaction cycles, speciations, etc., do not appear to be fluctuationdominated. A minor related comment, contrary to what is said int he appendix, the per capita death rate can be included into the continuous description (which is often also called the logistic model as all elements of matrix A are negative) by simply reducing the birth rate.
https://doi.org/10.7554/eLife.36273.029Author response
[Editors’ note: the author responses to the first round of peer review follow.]
Our decision has been reached after consultation between the reviewers. Based on these discussions and the individual reviews below, we regret to inform you that we have decided to reject your paper. However, we would like to offer the possibility of resubmission of a substantially revised manuscript. Our detailed assessment is as follows.
The paper has been carefully read by four expert reviewers. All of the reviewers found your paper thought provoking and potentially interesting. However, three reviewers raised substantial concerns. Two reviewers question the biological relevance of the models considered in the paper, and they think the model is too abstract to be useful for biologists thinking about longterm evolution. Since eLife is a "generalist" journal, it is imperative that papers have a high relevance not only for a narrow theoretical audience, but also for a general audience, in this case primarily in ecology and evolution. The reviewers also raised concerns that the work was not put into the context of existing theoretical work. In particular, one reviewer was concerned that a very similar approach has been taken in a previous paper, which is cited in the manuscript, but not discussed at all. All these concerns would need to be addressed in a revision. Specifically, the authors need to convince a mainly biological audience that the work presented is interesting and relevant, and they need to be clear about what advances the current work represents over previous work. Consequently, a serious effort would be required for revising the paper, and there is no guarantee for a positive outcome. Should you choose to resubmit a revised version, you would also need to take into account all the other reviewer comments in a constructive manner.
If you feel you can make the necessary revisions, we will make every effort to return a new version of the manuscript to the same editors and reviewers. Please note that this would be treated as a new submission.
Reviewer #1:
This is an interesting paper presenting a "minimal" model for the evolution of diversity in communities governed by competitive interactions. The model is minimal in the sense that only the outcome of interactions is modelled, but not their mechanistic basis, which greatly simplifies the description.
I have two main concerns, and a number of more technical points.
Major concerns:
1) The work presented appears to be very similar to the theory presented in Shtilerman et al., (2015). The authors need to clarify to what extent their work is novel compared to that earlier paper.
We have cited Shtilerman et al., (2015) in the previous and the current version of the manuscript as one of the most recent papers that use interactionbased modeling. Shtilerman et al., (2015) introduce their model as a modified version of the Ginzburg et al., (1988) model, a LotkaVolterra model with mutation events. The main flaw – and insight – of the Ginzburg et al., (1988) model had been its failure to generate the stable high diversity observed in actual competitive communities. Thus, it became clear that important aspects of reality had been missed. One way of fixing this problem is the introduction of ad hoc modifications likely to stabilize complex communities, which is what Tokita and Yasutomi (2003), Shtilerman et al., (2015) and others did. The main ad hoc modification in Shtilerman et al., (2015) is that they introduced explicit speciation events coupled to a large and sudden change of interaction between parent and o spring species, so that competition pressure between a parent species and its o spring species is much smaller than the pressure within each group. However, in real systems evolution happens gradually, and, whenever a mutation occurs, the resulting o spring species is most similar to its parent species, and experiences the highest competitive pressure from the latter.
In contrast, our Interaction and Tradeo based EcoEvolutionary Model (ITEEM) does not make these assumptions. In particular, it does not enforce speciation with an explicit mechanism, but novel species emerge out of the ecoevolutionary dynamics.
In the dynamics simulated in ITEEM, parent and o spring strain are in fact most similar, as in natural communities. The main conceptual advance that comes with ITEEM is the following. If we acknowledge the natural fact that ecoevolutionary dynamics is bound by physical constraints – in ITEEM acknowledged as tradeo between the abilities to compete and to replicate – major features of real competitive communities, such as stable high diversity, just emerge.
In a direct comparison with ITEEM, the model by Shtilerman et al., (2015) corresponds to the limiting case of zero tradeoff, though with ad hocfeatures to enforce diversity. In terms of evolution, the Shtilerman et al., (2015) model could be considered closer to cladogenesis or punctuated equilibrium, while ITEEM models phyletic gradualism in an interactionbased formalism.
Besides, we study further aspects not addressed by Shtilerman et al., (2015) in their model, such as the effects of lifespan, system size, and mutation rate on diversity. Moreover, we think that the investigation of the role of cycles in the manuscript brings us closer to understanding the mechanisms behind the dynamics of diversity, be it in the growth of stable complex communities or collapses.
In the new version of the manuscript we have expanded the fourth paragraph of theDiscussionwithadescriptionofthehistoryofinteractionbasedmodels, including a description of Shtilerman et al., (2015).
2) In fact, I think the logistic equations presented by the authors in the appendix are similar to the ones used in Shtilerman et al., 2015. I think the authors need to consider using these equations, rather than their individualbased simulations, to derive their results. It appears to be straightforward to implement the evolutionary dynamics based on the deterministic logistic equations, with one differential equation per species, and with evolution being implemented by adding new equations to the system (and deleting those ones in which the frequency falls below a threshold). This would yield a computationally more efficient model, whose salient features would nevertheless be essentially the same as in the individualbased model (perhaps barring the nearly neutral variation within species clusters). If the results are not the same, then this needs to be known, as it would probably point to an important role of stochasticity. I think this line of analysis is necessary to obtain a complete picture of the system studied by the authors, which is required for eLife.
As explained in the Model section and Appendix 3, and now expanded in the Introduction and the last paragraph of section Model, ITEEM is an interaction based model for LotkaVolterra dynamics of competitive communities. It is not a logistic model.
A standard way of simulating LotkaVolterra dynamics has been to treat populations as continuous quantities, and to integrate continuous LotkaVolterra equations over time. This avenue was chosen by Ginzburg et al., (1988) and all descendant models, including Shtilerman et al., (2015). A second standard approach is the use of individualbased models (IBMs) (Black and McKane, 2012; DeAngelis and Grimm,2014) that simulate the LotkaVolterra dynamics at the level of individuals. We have used the second approach but we show in the Supplementary material that, if we apply a meanfield approximation to our model and treat the population as continuous, we formally obtain the wellknown LotkaVolterra equations.
It is important to realize that although the two approaches are equivalent in the meanfield limit, but correspond to distinct models of reality. Clearly, communities are composed of individuals, and disregarding this fact by replacing them by a continuum will lead to artifacts. It has been extensively discussed in the literature that continuum approaches are unsuitable in cases of nonequilibrium dynamics, if demographic variations influence the dynamics, when size of the system or discreteness of individuals matter or when chaotic dynamics can cause unexpected behavior, and that then stochastic, individualbased models are better models of nature (McKane and Newman, 2004; DeAngelis and Mooij, 2005; Black and McKane, 2012; DeAngelis and Grimm, 2014).
Since our system falls into the category of those better modeled by individual based models, we chose such a model right from the beginning and do not see the necessity to reiterate the above studies. Specifically, the following points of our study required an individualbased model: the ecoevolutionary dynamics is of a nonequilibrium type; we wanted to study the effects of demographic variations of lifespans; we had to study the effect of systemsize, e.g. to test the observed phase transition.
To illustrate these points further, we refer to the important fact were the Reviewer herself/himself sees a difference between continuous and individualbased model, namely that the numerical solution of the continuum equations fail to produce neutral variations within species (“perhaps barring the nearly neutral variation within species clusters” in reviewer’s comment). However, as explained in the manuscript, this variation plays an important role in the first stage of speciation inside each cluster (see e.g. Vellend, (2006)). Similarly, stochasticity is known to impact extinction (Melbourne and Hastings, 2008).
To make it more clear that the choice of an individualbased model is wellfounded, we have provided the key arguments in the first paragraph in the Appendix 3.
Essential revisions:
1) Subsection “Model”: the procedure for the individualbased model is not clear: the text states that "we try 𝑁𝑖𝑛𝑑(𝑡) (number of individuals at that time step) replications of randomly selected individuals. Each selected individual of a strain 𝛼 can replicate with rate r_α…" this is imprecise; what does is it mean that an individual that has already been selected for replication replicates with a rare r_\α? Isn't it rather the case that individuals get selected to produce an offspring with a probability that is proportional to r_α?
The selected individuals do not reproduce deterministically but with a probability r_{↵}, as is customary in individualbased models. In the new version of the Model section we have clarified the description of the procedure by using the term “replication trial” (third and fourth paragraph of Model section quoted here for convenience; see esp. the respective first sentence):
“Every generation or time step consists of N^{ind}_{(t)} replication trials of randomly selected individuals, followed at the end by a single death step. In the death step all individuals that reached their lifespan at that generation will vanish. Lifespans of individuals are drawn at their births from a Poisson distribution with overall fixed mean lifespan. This is equivalent to an identical per capita death rate for all species. For comparison, simulations with no attributed lifespan (At each replication trial, a randomly selected individual of a strain= 1) were carried out, too. _{↵} can replicate with probability r _{↵}. With a fixed probability µ the o spring mutates to a new strain ↵0. Then, the newborn individual is assigned to a randomly selected site. If the site is empty, the new individual will occupy it. If the site is already occupied, the new individual competes with the current holder in a lifeordeath struggle. In that case, the surviving individual is determined probabilistically by the “interaction” I_{↵},defined for each pair of strains ↵,. I_{↵}is the survival probability of an ↵ individual in a competitive encounter with a. All interactions ↵ form an interaction matrix individual, with that encodes the outcomes of all I_{↵}2 [0,^{1]} and I_{↵}+ I_{↵}^{= 1}. _{I} _{I}_{(t)} possible competitive encounters.”
2) Also, it is not clear how the discreteness of time is dealt with exactly. For example, the paper says that each offspring is assigned a randomly selected site, then competes against the occupant of the site and probabilistically takes over. What if one offspring takes over, but then a second offspring is chosen (by chance) to compete for the same site? Is that second offspring then competing against the previous occupant, or against the new occupant, i.e., the offspring that has already taken over from the previous occupant?
This is now clarified in the fourth paragraph of the revised subsection “Model” (second quoted paragraph in previous response).
3) And when do individuals die? Before or after offspring production?
Individuals die after o spring production, see newly formulated third paragraph of subsection “Model” (first quoted paragraph in response to point 1).
4) The authors mention the existence of rockpaperscissor interactions, but it is not at all clear what the importance of such interactions are for the overall dynamics of the system. Are such interactions very common, essential, just a quirk?
Rockpaperscissor interactions are synonymous to cycles in the dominance network (see the response to the 6th technical point of the same reviewer). Such cycles are commonly emerging in competitive systems in ITEEM (not in the neutral model) and in nature as essential stabilizers of diversity. In the revised manuscript, the relationship between diversity and cycles is shown in Figure 2D, 2G and Figure 3A, and described in the corresponding text, especially subsection “Generation of diversity” and subsection “Impact of tradeoff and lifespan on diversity”. See also the new section in Appendix 8 “Collapses of Diversity” which describes the role of cycles in evolutionary collapses.
5) Subsection “Impact of tradeoff and lifespan on diversity”: it is not clear what "average cycle strength" is. Which cycles were used to calculate the average of what?
In the revised manuscript, we have provided a detailed description of the algorithm in the Supplementary material, Appendix 7. Shortly, we compute the average cycle strength in the dominance network as average over a large number of cycles between randomly selected triplets of nodes (=individuals).
6) Subsection “Impact of tradeoff and lifespan on diversity”: shouldn't the average interaction strength, i.e. the average values of all I's, be ca. 0.5?
To improve the presentation of the network results, we define in the revised manuscript the dominance networkand detail its computation (subsection “Generation of diversity” and Appendix 7). Briefly, a dominance edge points from an individual ↵ to an individual if ↵ dominates, i.e. if ↵ wins against in an encounter. The weight of the edge is given by the interactions as I^{↵}_{I}^{↵}. The mean dominance is now reported also in Figure 3a, instead of average interaction strength. In ITEEM system that starts from monomorphic population, deviation of the mean dominance from 0 and the corresponding mean of neutral model shows that selection generates species with significant pairwise dominance
7) I failed to understand Figure 2. Please explain better.
In the revision we have clarified the terminology (dominance network, see previous response), reworked the Figure (now Figure 3) and its caption, reformulated and extended the corresponding description in the main text (subsection “Impact of tradeoff and lifespan on diversity”) and Appendix 6.
8) Subsection “Impact of tradeoff and lifespan on diversity”: the Darwinian demon doesn't really make sense: with no diversification, there is only one strain present, and it is impossible to say what kind of competitive ability this strain has, because no other strain is present… so what are you actually trying to say?
We have reformulated the corresponding paragraph of the Results section (“Without tradeo (= 0).…”). To summarize the argument here: We follow in our definition of the Darwinian demon the notion in Law (1979), i.e. a hypothetical organism which can maximize all aspects of fitness simultaneously and would exist if the evolution of species was entirely unconstrained. In our model fitness is composed of reproductive and relative competitive fitness. In the zero tradeo limit, strains with high relative fitness emerge and quickly conquer all resources. If established, a mutant of that type could acquire yet another advantageous variation that would make it superior to its ancestor without any penalty in another trait. Thus, we see in phase I of Figure 3 a sequential dynamics of emergence, invasion, extinction with constantly changing traits. This makes the phase I different from the phase III, in which no diversification occurs and strains remain closely related variants of the same species.
9) I don't find the section about the relationship between tradeoff and productivity very convincing, because it is too vague. I think the authors need to be more specific and detailed here.
We have extended the description of this relationship in the revised manuscript.
10) Subsection "frequencydependent selection": I think a more mechanistic explanation for differences of speciation and extinction rates compared to Poisson processes is needed. What is the cause of these differences? When is this mechanism at play, and when not? In particular, why do speciation rates increase with the diversity in the system? Shouldn't this only happen at certain intermediate levels of diversity, because when the community reaches a high level diversity, speciation rates presumably decline?
In a system with balanced speciation and extinction rates a null hypothesis is that these events can occur randomly with a constant rate. This happens when ecological interactions and frequencydependent selection are not the underlying mechanisms for diversification, speciation and extinction. By comparison of distribution of time intervals between branching points with the Poisson distribution we reject this null hypothesis. In ITEEM when one species (as a group of similar traits) emerges in a branching event, it provides a new function to the community and this increases the probability for the generation of further species. The most plausible explanation for the observed deviation from Poisson process is frequencydependent selection which is mediated by interaction of species and causes autocorrelation between evolutionary events over time. In the new version of the manuscript we have expanded the subsection “Frequencydependence selection” and also modified Appendix 10 to clarify these points. Please also see the response to the first major concern of reviewer 2 about speciation rate.
Regarding the question about the increase in speciation rate: We did not explicitly report the speciation rate versus diversity and so did not claim that speciation rate increases with diversity. We agree with the reviewer that in the early stage of simulation, speciation rate increases (the positive feedback loop between number of species and speciation that is discussed in the first major concern of reviewer 2), and then will decrease when the number of species in the system saturates to a more or less stable value (Rabosky, 2013). In ITEEM simulations, this saturation depends on system size and mutation rate (Figure 2A and 2D, Appendix 6 (Figures 4 and 5) where diversity increases very fast and then saturates to a constant value). This saturation confirms the reviewer’s point about the decline in speciation rate. See also the answer to the first concern of reviewer 2.
– I would like to bring the authors' attention to the recent paper by Doebeli and Ispolatov, (2017), which also presents models for community assembly and evolution of diversity are presented, with some parallels to the work discussed here. These parallels could be taken up in the Discussion section.
We are grateful for this reference and have included it in our discussion.
Reviewer #2:
The authors present a model for studying diversification dynamics as an outcome of, as they claim, competition for a single resource. I have substantial concerns regarding the insight that could be gained from this work, mainly due to the type of model and the lack of a mechanistic interpretation.
1) Subsection “Model” and subsection “Frequentdependent selection”: Your trait space is Nspdimensional, where Nsp is the number of extant strains. In particular, the trait space grows explicitly with every strain added. Is it surprising that you get coexistence of Nsp strains in an Nspdimensional trait space? Further, is it surprising that "emergence of new species increases the probability for generation of further species"? It seems to me that this is quite expected, given your model. This also explains why you get these mass extinction events; there's a positive feedback loop between speciation rate and species number that is easy to recognize even without simulations.
This is a complex comment that we split up into three questions: (1) Is the coexistence of N_{sp} species in N_{sp} dimensional trait space surprising? (2) Is it surprising that emergence of new species increases probabilities for generation of further species? (3) Is there a positive feedback loop between speciation rate and species that explains mass extinctions?
The first question may have been arisen because of a misunderstanding of the term trait. We have to emphasize that in our model the term trait does not carry the usual meaning of a phenotypic trait (wikipedia: “A phenotypic trait is an obvious, observable, and measurable trait; it is the expression of genes in an observable way”), but that it is an interaction trait. Even in a system of n strains that are distinct with only one phenotypic trait, each strain has n interaction traits, quantifying its interactions with all strains in the system. In ITEEM we only consider interaction traits and are agnostic about the number of phenotypic traits. Thus, the underlying eco evolutionary system could well have only a few phenotypic traits, and we consider it surprising to see under these conditions the emergence of stable diversity (As it has been shown that a simple GLV model without further assumptions fails to produce diversity. See the fourth paragraph of Discussion section).
However, a crucial point is that even in a system with an N^{st} dimensional trait vector, diversity only emerges if we account for physicochemical constraints that appear in ITEEM as tradeo between replication rate and competitive ability. Depending on the tradeo strength, the diversity can be low or high, stable or unstable (see the phase diagram in Figure 3 of manuscript), although the trait vector has always more or less the same dimension N_{st}. As we emphasized in the manuscript (section Generation of diversity) diversity is not measured as the number of strainsbut as the number of species, i.e. coexisting clusters of strains that emerge via diversification.
We hope that the new version of the manuscript clarifies these points (see Materials and methods section and the first paragraph of Discussion section.)
The second question of this comment by reviewer 2 was about whether it is surprising that new species increase the probability for further species. This question might have been arisen by a confusion between strain and species, possibly because of the unfortunate symbol N_{sp} for the dimension of the interaction trait vector T_{↵}. In the new version we have replaced this symbol by N_{st}, clarifying that the dimension of T_{↵}is given by the number of strains (not the number of species). We hope that the more clear nomenclature, especially the use of N_{st} instead of N_{sp} introduced above should sharpen the distinction between species and strain that we make, in contrast to other works. Strains emerge under certain conditions (tradeo regimes) continuously as individuals successfully replicate with mutation, while speciation events are more rare and require a branching of the phylogenetic tree and diversification in trait space. Given our model it is indeed “quite expected” that the number of strains N_{st} initially increases in the almost empty system, and that we then have “a positive feedback loop between speciation (emergence of strains) and species (strain) number”. But this early stage comes quickly (after less than as the system saturates with individuals. The number of⇠ 2 ⇥ 10 strains^{5} generations) to an end is then more or less stable for long times. Likewise, the number of speciessaturates at values that depend on tradeoff, system size, and mutation rate (see Figure 2D and Appendix 61).
Our analysis focuses not on the initial stage, but on the stage after onset of saturation when we have a balance between speciation and extinction. In such a system, the null hypothesis is that speciation and extinction occur by random drift. However, in the manuscript we can quantitatively reject this null hypothesis. Instead, our analysis demonstrates that in ITEEM the emergence of new species (cluster with similar traits) by branching phenomena provides a new function to the community, which increases the probability for the generation of further species (Herron and Doebeli, 2013). We hope that the new subsection “Frequencydependent selection” and Appendix 10 provide a better explanation than the first version of the manuscript.
For the mass extinctions, the third point in the comment, similar arguments apply: The phenomenon is occasionally observed in the saturation phase, and it is not just a trivial consequence of the basic model but depends on tradeo strength and other parameters. To better explain this complex phenomenon we have improved the analysis and presentation of mass extinctions, especially in subsection “Generation of diversity”and the new Appendix 8.
Please also see the answer to the tenth technical point of reviewer 1 about speciation rate and the response to the third concern of reviewer 3 about mass extinction.
2) Subsection “Generation of Diversity”: Please explain how you define "functional diversity" in your model at first mention. Appendix 6 was uninformative. The only information I found was in the caption of Figure 1, where you define functional diversity "in terms of the size of minimum spanning tree (SMST) in trait space". But your trait space is continuous ([0,1]^Nsp) so a priori there is no tree structure connecting species (apart from phylogeny). So how do you define a spanning tree, and how is this spanning tree not a function purely of the number of strains?
We have added in the revised manuscript in caption of Figure 2d and in the text referencing this Figure (end of the subsection “Generation of diversity”) as a short explanation that the functional diversity (FD) quantifies the spread of the population in interaction trait space. On both occasions we refer the reader to Appendix 6. The revised and expanded Appendix 6 covers evaluation of FD including further explanations and references.
On the minimum spanning tree (MST): The Reviewer gives a good example for a tree structure in a continuous space induced by a distance structure on species, namely phylogenetic trees. In practice, phylogenetic trees are often inferred from estimated evolutionary distances (i.e. continuous quantities) between pairs of all considered taxa, assuming minimum evolution (see e.g. Saitou and Nei, (1987); Felsenstein, (1997)). In our manuscript, the MST is computed in an analogous way. Here we have between each pair of strains a distance in interaction trait space, and from the set of these distances we infer a minimum spanning tree (MST). The MST is a tree that links all strains and has a minimum sum of edges (sum of edges of MST = SMST). For the computation of the SMST see Appendix 5.
The SMST quantifies the spread of the total population of strains in interaction trait space, i.e. the space that matters for function in the community. To us it was therefore the most natural way of expressing functional diversity. As we detail in Appendix 63 there are a number of FD measures used in the literature that we have computed (Mouchet et al., 2010). However, all FD measures show the same strong dependency of tradeoff (Figure 6 of Appendix) and are not simply reflecting the number of strains, to answer the last question of the reviewer. For example in Figure 1 of the Appendix, we compare a set of strains evolved by neutral evolution (no difference in the interactions and no tradeo regime) with a set of strains at an intermediate tradeoff, both with the same number of strains. The former set is compact (low FD) while the latter is dispersed in interaction trait space (high FD, note the different scaling).
3) Subsection “Generation of Diversity”: Related to the previous point. What is a "functional niche" in your model? Since you did not discuss mechanisms underlying the interaction matrix, it is not clear what a function is.
In ITEEM the function of a strain is its interaction trait vector, describing its interaction with each member of the community. A functional niche is a cluster of strains in the space of interaction trait vectors. This definition is given in subsection “Generation of diversity” (paragraph ITEEM systems selforganize toward structured communities), and it is in line with the notion of functional niche developed in the literature (Elton, 1927; Clarke, 1954; Whittaker et al., 1973; Rosenfeld, 2002; Odum, 1959) as detailed in the new version of Appendix 62, subsection “Functional diversity, functional group and functional niche”.
4) Subsection “Generation of Diversity”: Saying your model has no geographical isolation nor resource partitioning sounds meaningless. You have not specified the underlying mechanisms for your "interaction trait", so it is hard to draw a comparison to other more mechanistic models (i.e. where the physiological/metabolic traits are explicit).
As we explain in the revised subsection “Model”, the system can be considered well mixed so that no individual is assigned to a specific spatial location and each individual accesses the same resources. In this respect ITEEM is not different from referenced interactionbased models, the LotkaVolterra competition equation, or the replicator equation.
5) Subsection “Model”: Since replication is nonsexual, what is the purpose of defining "species" separately from "strain" in your model? In AppendixI you mention that you define "species" operationally based on divergence time and using some arbitrary cutoff threshold, so this sounds a bit analogous to operational taxonomic units (OTUs) in microbiology. However, you model mutation/selection/growth dynamics at the strain level. Please explain early on what additional insight one may get from counting "species".
In ITEEM we observe reliably the emergence of clusters of strains of similar interaction traits, and we find that these clusters play important roles in ecoevolutionary dynamics. Thus, it seemed natural to us to introduce the term species to simplify the text and to offer a biological interpretation. We follow there the common usage of species in models of asexual populations as a term for separated clusters of very similar genotypes or phenotypes, depending on the model, see e.g. Dieckmann and Doebeli, (1999); Bonsall et al., (2004); Sevim and Rikvold, (2005); Scheffer and van Nes, (2006); Mathiesen et al., (2011); Ispolatov et al., (2016). See also the expanded Appendix 1.
6) Discussion section: You claim that modeling in terms of interaction traits (which in your case means in terms of outcome probabilities of local competitive exclusion) "coarsegrains these complex systems in a natural, biologically meaningful way.". But many of your interpretations and comparisons to other models or data require translating your Abstract "interaction traits" to functions or life histories. You did not discuss at all what real traits could possibly give rise to your interaction trait matrix. Your interaction matrix seems to loosely represent competitive interactions; but then it can only explain diversity within a single trophic level (e.g. in the case of animals) or a single metabolic niche (e.g. in the case of bacteria); yet, you keep referring to "functional diversity".
In our coarsegrained model the interaction trait of a strain summarizes all negative effects of the competing populations on the population of that strain, no matter whether these strains compete for the common resources by interference or by exploitation, and, in this sense, our interaction trait matrix is similar to the community matrix of the wellestablished Generalized LotkaVolterra equation (see penultimate paragraph of revised Model section). Specifically, the coarsegrained interaction trait summarizes over all phenotype traits, abilities and mechanisms that increase the chance of organism to take over resources, help the species to tolerate reduction in contested resource availability (Aarssen, 1984), or prevent competitors from gaining resources (Gill, 1974). The current version of Appendix 3 takes up this point again.
On the last point of the comment, the use of “functional diversity”: We agree with the reviewer that the presence of different trophic levels implies functional diversity. However, the reverse conclusion is not necessarily true if we understand function as the role in a community, not limited to the mode of resource consumption. We follow here this broader interpretation of function (see e.g. Petchey and Gaston, (2006)) that, for instance, considers interaction traits. If we accept this broader interpretation of function, functional diversity can be measured as spread in interaction trait space (see also responses to comments 2 and 3 of this reviewer).
7) Discussion section: Related to the previous comment. You say that your formulation was "chosen to reflect reasonable properties" and that you "have assumed a single, limiting resource in a wellmixed system". However, you did not provide any plausible mechanism for how such an interaction matrix could arise from competition for a single resource in a wellmixed system. This is actually a big question: how can one obtain a Nspdimensional trait space through competition for a single resource pool?
The comment of the reviewer implies that we are modeling a competition over the consumption of limited resources. However, the situation modeled here is that of a competition over one common, renewable resource pool (see section Model). There is a long history of models for this situation, for instance the LotkaVolterra equation (Volterra, 1928) and interactionbased ecoevolutionary models (Ginzburg et al., 1988); in our case, we are using the LotkaVolterra approach, which models both exploitative and interference competition. In these models, an interaction matrix quantifies how competitive organisms influence, by any means, N_{st}_{⇥}N_{st} each other’s populations in pairwise dominance relations. The final question by the reviewer is closely related to his/her first comment, and we refer the reviewer to our response there.
Reviewer #3:
The manuscript "Tradeoff shapes diversity in ecoevolutionary dynamics" presents a quite original evolutionaryecological individualbased model based on a compromise between growth rate and competitiveness. The model is simple but exhibits quite rich evolutionary properties, which, while being not entirely unpredictable, are intriguing and provide useful insights and generalizations. This combination of simplicity of the rules and relevance of the dynamics usually distinguishes successful and longliving models from the rest and makes them understandable and appealing to a wide audience. Besides, I cannot see any potentially "hidden" flaws in the model and interpretation of results that could cast doubts on the main conclusions. Thus, in my opinion, the manuscript may be published in eLife after the following and, perhaps, other comments are addressed:
1) The complexity of algebra in (3) and subsequent definition of s is definitely unwarranted by an otherwise quite accessible level of math in the rest of the paper. Furthermore, it apparently confuses even the Authors when they classify the domains of weak and strong tradeoffs:
It looks like the definitions of hightradeoff (\δ ~ 1) and intermediate tradeoff (\δ ~ 1/2) are misleading. The strongest dependence of r on C appears to be when \δ=1/2, which also follows from Figure 1. The phrase "In hightradeoff phase III, any small change in C changes r drastically" is simply wrong for all but very small C. I would call both the \δ=1 and \δ=0 limits as small tradeoffs as they look perfectly symmetric, or, even better, choose another, more heuristic and transparent, form of parametrization of r vs. C.
We are grateful to the reviewer for pointing out that the description of the tradeo needs clarification. To this end we have moved Eq (3) from section Model to a revised Appendix 21, and, instead, added to subsection Model the new Figure 1, which depicts actual shapes of the tradeo function in a more accessible way. We have also made the behavior of the tradeo function and the meaning of more transparent (subsection “Model”).
The new section Appendix 23 should now make it more clear why the tradeo
increases from = 0 to = 1. To summarize the main point here: The strength of the tradeo is not simply given by the slope of the ∆rr/∆CC=drdC.Cr tradeoff curve, but the computation of the actual strength of the tradeo requires consideration of the community context. Specifically, the tradeo strength can be approximated by, with the replication rate r and competitive ability C of the current community, slope drdC of the tradeo curve, and perturbations of reproduction rate r and competitive ability C due to emergence of a new mutant. This tradeo strength increases with tradeo parameter (Appendix 23).
2) Would it be possible to say anything about population dynamics of "species"? Is it cyclic?
We have looked for cyclic population dynamics by Fourier analysis but could not see clear signals. One problem may be our limited time resolution (10^{4} generations as sampling interval) so that we do not see short cycles, but the dynamics may also be chaotic because we have many interacting and evolving species.
3) At least qualitatively, what are mechanisms of mass extinction?
In the revised manuscript we sketch a mechanism possibly underlying mass extinctions in Results, subsection “Generation of diversity” (paragraph “Formation of strong cycles…”). Appendix 8 provides more details.
A problem for this part of the study was that we see mass extinctions only occasionally so that we do not have a strong statistics. Moreover, we have a limited time resolution of 10^{4} generations. We do not see conspicuous diversity values preceding mass extinctions. However, one remarkable pattern that occurred prior to all observed mass extinctions was high cycle strength (Appendix 8, Figure 7B). It could be that new species emerge or species go extinct, resulting in conflicting forces in the strong network that discharge in a collapse, similar to a crystal defect leading to geometrical frustration and destabilization of a highly ordered crystal lattice.
4) An explicit plot of diversity vs. system size (perhaps just for the δ optimal for diversity) and, ideally, an estimate of corresponding scaling, would be very revealing.
Appendix 9 provides a set of plots of diversity vs. system size for different tradeoffparameters. The new Figure 9 in Appendix 9 shows the scaling of diversity with size for middle values of.
5) Similarly, how the level of diversity and the typical number of traits in a species depends on the mutation amplitude m?
We have added a description of the effect of variation in m to the Results, subsection “Effect of mutation on diversity”. Briefly, the smaller m, the slower the dynamics and the smaller the extension of individual species (= clusters of strains) in trait space. Large m leading to big trait changes suppress the establishing of diverging evolution and selection
6) For a general reader of eLife, the MDS algorithm needs to be explained and properly referenced. It plays a major role in interpreting the results, however, is presented only by the name of the corresponding function in R.
Explanations of MDS in Appendix 4 and the caption of Figure 2 (previously Figure 1) of the main text have been expanded. However, MDS is only used for visualization. All quantitative analyses are based on evaluations of the original high dimensional traitspace.
7) A qualitative explanation about the minimal spanning tree analysis would help as well.
The explanations on the minimum spanning tree analysis in Appendix 5 have been expanded.
Reviewer #4:
I've carefully read the paper "Tradeoff shapes diversity in ecoevolutionary dynamics" by Farnoush Farahpour and colleagues. In it, they use an individualbased ecoevolutionary model to understand the emergence of community structure based on evolving interactions. The model produces some interesting patterns, such as clustering in traitspace and the evolution of intransitive competitive loops.
There are a number of aspects of this paper that I liked: the focus on evolution of interactions, the emergence of intransitive loops, the dimension reduction applied to the vectorvalued traits, and the comparison with a neutral model variant. These are all creative contributions to the modeling literature.
While I was intrigued by the idea of modeling the evolution of interactions directly, it was hard for me to connect it to real ecological systems. The interactions between species are not determined by phenotypic traits of the organisms but evolve independently. It's based on an unstated assumption that species interactions are totally idiosyncratic and unpredictable. I feel that evolution in this model is too unconstrained, despite the tradeoff between competitive ability and reproduction, resulting in the prevalence of intransitive loops. As a complexsystems researcher, I'm fascinated, but as an evolutionary ecologist, I'm skeptical.
The authors need to do a better job putting their work in the context of the extensive literature on ecoevolutionary dynamics. The text had many statements that had me scratching my head. Examples:
1) Framing the problem in terms of resources, but in the model there are no actual resources.
There is one durable resource or resource pool in the model for which individuals compete (see the first paragraph of Model section and the last paragraph of Appendix 3). In this respect ITEEM is very similar to established models in ecology and evolution that are based on LotkaVolterra equations (Gill, 1974; Masel, 2014) which model exploitative and interference competition for durable resources.
Note that these models are different from the consumerresource model, which models indirect competition for consumable resources, and from the logistic equation, which models exploitative competition for durable resources.
2) The competitive exclusion principle only holds at equilibrium (Armstrong, Levins) and must count resources plus shared predators (Levin, 1970).
We agree with the reviewer that the principle holds only at equilibrium although it is still referred to in nonequilibrium systems and even in the presence of disturbance, adaptation and evolution (Posfai et al., 2017; Germerodt et al., 2016; Kinnersley et al., 2009; Pfei er and Bonhoeffer, 2004). We had referenced the principle as a contribution to the historical context of the diversity debate. But since it might lead to confusion and because it is not the concern of the manuscript, we have removed the reference from the Introduction.
3) The "ecoevolutionary models" cited in the Introduction seem to be just ecological models.
In the revised manuscript we have rephrased the paragraph accordingly.
4) In the Introduction, the "observed ecoevolutionary dynamics" cited that the model "closely resembles" aren't empirical patterns observed in real systems, but just results of other models.
The references are mixed empirical (see for example Maynard et al., (2017); Kvitek and Sherlock, (2013); Bolnick and Fitzpatrick, (2007); Coyne, (2007); Herron and Doebeli, (2013)) and theoretical of welldescribed evolutionary and ecological phenomena such as sympatric speciation, emergence of two or more levels of differentiation similar to phylogenetic structures, large and complex biodiversity over long times, evolutionary collapses and extinctions and emergence of cycles.
5) Discussion of speciation overlooks that the species here are clonal, so what's hard about speciation?
Speciation would indeed be trivial if it was the generation of new (clonal) strains. However, speciesare not clonal (see e.g. subsection “Model” or Figure 2). The model continuously produces new strains, but the diversification and branching of species, i.e. wellseparated nontransient clusters in trait space, are rare and heavily depends on tradeoff.
6) Also, some of the sentences throughout were hard to understand the meaning of. E.g., "Evolutionary changes at the genetic level influence ecology if they cause phenotypic variations that affect biotic or abiotic interactions of species which in turn changes the species composition and occasionally forces species to evolve their strategies."
We have split this sentence up to make it more intelligible.
7) Some of the details of the model implementation weren't clear. For example:
a) How exactly do births happen (subsection “Model”)?
This has been clarified in the revised manuscript. See also response to 1st and 2nd technical point of reviewer 1.
b) Is mu a mutation "rate" or probability of mutation during a replication event (subsection “Model”)?
This has been clarified in the revised manuscript. See also response to 1st technical point of reviewer 1.
c) Why is lifespan drawn from a Poisson distribution (Subsection “Model”) and how can that be infinite (Figure 2)?
We used the Poisson distribution as a simple positive valued probability distribution that for large means approaches a normal distribution.
In the case of an infinite lifespan, no lifespan is attributed to individuals. In that case the only cause of death is defeat in a competitive encounter. The corresponding sentence in the subsection “Model” has been expanded for clarification.
If each individual stays in a site, is it really wellmixed?
Yes, because there is no neighborhood relation between sites. A site is just a portion of the durable resources, sufficient for the sustenance of an individual. The property of being wellmixed means that arbitrary encounters between new individuals and sites happen with equal probabilities. The subsection “Model” and Appendix 3 have been modified to clarify these points.
Does mutation of one species' interaction coefficients end up changing another species' reproduction rate through the tradeoff (2)?
Indirectly, yes. Adding a mutant strain introduces a new interaction trait and thus in general influences the competitive ability C. A change in the competitive ability will then change reproduction rate r through the tradeoff. This is qualitatively similar to reproductive plasticity observed in natural systems (Claridge and Franklin, 2002; Buffi et al.,, 2013; Goldstein et al.,, 2016). However, when ITEEM systems reach high diversity (i.e. complex interaction trait vectors), this effect is usually small.
Could you not get at the same questions more efficiently using deterministic LotkaVolterra dynamics?
Please see the response to the point 2 of reviewer 1.
[Editors' note: the author responses to the rereview follow.]
The reviewers agree that substantial effort was put into the revision. However, there are a number of remaining issues that call for a further indepth revision. Reviewers 1, 2 and 3 (who was previously reviewer 4) still have substantial concerns about aspects of the paper. These concerns are all related to issues raised in the initial reports and will need to be addressed in a constructive manner if the paper is to be published in eLife.
Among the issues to be addressed are:
1) Individualbased vs deterministic systems: both reviewer 1 and 4 (formerly reviewer 3) point out that it is not at all clear that these different approaches would yield different results. In other words, it is not clear that salient results reported in the paper depend crucially on the presence of noise, as appears to be the contention of the authors. This needs to be explored at least to some extent. This issue cannot be dealt with by simply "citing it away".
To avoid misunderstandings, we would like to state our notion of the terminology used in the referee report and in our response below. In the referee reports, four nonsynonymous terms have been used to denote alternative ecoevolutionary approaches that should or could be compared to ITEEM: (1) deterministic models, (2) continuum approaches or continuous models, (3) population models, and (4) mean field models.
The fourth term, meanfield models, is used with different meanings in the comments of reviewers 1 and 4. In the response to the comment 5 of the fourth reviewer we explain our understanding of the different usages of the term “meanfield” in statistical mechanics and in population biology. To avoid misunderstandings, we have replaced in the manuscript the term “meanfield” with “populationlevel”, a term that should be more clear for the target audience.
Thus, we are left with three nonsynonymous terms used by editor and reviewers. We hope that the reviewing team agrees with our understanding of the terms:
1) Deterministic model: A mathematical model that unfolds completely deterministically, in contrast to models that make use of randomness.
2) Continuous model (in time): A mathematical model with variables changing continuously in time without any abrupt transition between discrete states, in contrast to models with discrete time steps.
3) Population model: A mathematical model that studies dynamics of populations and their interactions. Such models often describe populations as continuous quantities, in contrast to models with discrete individuals and individual interactions.
The editor takes up comments by reviewers 1 and 4 on the comparison of ITEEM with deterministic/continuous/population/meanfield models. Since both reviewers solely refer to Shtilerman et al., (2015) as reference for this latter class of models, we assume that they propose a comparison of ITEEM with the model by Shtilerman et al., (2015), or equivalent models. Note that the model by Shtilerman et al., (2015) is a hybrid with stochastic and deterministic elements that are applied iteratively: each stochastic evolutionary “kick” is followed by an ecological relaxation of fixed length in which the deterministic LotkaVolterra population dynamics is integrated. This “kickrelax” model as a whole is neither deterministic, nor continuous in time, nor continuous in population.
One reason suggested by the reviewers for a comparison of such a kickrelax model to ITEEM is that the former would speed up the simulations, thus enabling to get cleaner results with larger systems. We have therefore developed a kickrelax model that is as close as possible to ITEEM and compared its output to the ITEEM output. We have further extended subsection “Interactionbased models” on this topic.
To not water down by technical details the main point of the manuscript, i.e. the effect of tradeoff on ecoevolutionary dynamics, we provide the following detailed information here in our response, but we have not included it in the manuscript.
Tests with a kickrelax code show that the method is much less efficient than ITEEM (though with a caveat, see below). At the same system size of 10^{5}, i.e. number of sites in ITEEM and carrying capacity in the kickrelax model, mutation rate 0.001 and lifespan 10^{5}, we achieve 15×10^{4} generations with ITEEM and 200 generations with the kickrelax code per CPU hour. The main reason is that the integration of a LotkaVolterra system of hundreds or thousands of strain equations is computationally far more demanding than the simple operations on an array of individuals in ITEEM. Because of these computational problems we could not collect as many generations for kickrelax as for ITEEM: with kickrelax we sampled at substantial computational cost about 6 × 10^{5} generations for each of 9 simulations (9 tradeoff parameters and one lifespan); with ITEEM we sampled 5 × 10^{6} generations for each of 540 simulations (3 runs for 18 tradeoff parameters and 10 lifespans). Our kickrelax code is an Rscript with the CPU time critical integration step performed by a fast FORTRAN routine. We have made the kickrelax and ITEEM source codes freely available (see the end of subsection “Model” of the manuscript).
The kickrelax model is also very sensitive to the length of the relax interval, i.e. the time interval between two successive kicks. To obtain results close to the nonequilibrium model of ITEEM, this time interval must be small which makes the simulation slow. Increasing this time interval 10fold abolishes the adaptive evolutionary dynamics and gives rise to communities with no speciation and diversity.
The caveat mentioned above is that a clean comparison in terms of efficiency requires that the results are the same. However, the results are qualitatively different. One important biological result of ITEEM is the phase diagram of diversity as function of tradeoff 𝛿, and the humpback shaped crosssection through that phase diagram (Figure 4 of the current version of the manuscript and Figure 8 of Appendix). The Figure below opposes this humpback curve, in some diversity indices, with the corresponding curve obtained with the kickrelax model.
The reason for the discrepancy lies in the interface of the kick and relax stages: We have to decide in the evolutionary kick stage what happens to which population, how and when we introduce a mutant, and when we consider a population extinct in the ecological relaxation stage and remove it from the pool for evolution. This interface is not motivated by biology but it is a consequence of the technical realization of the model and therefore necessarily artificial. In our kickrelax simulations we set the extinction threshold to 1. In the kick stage, the number of new mutants that each extant strain 𝛼 could produce during the relax interval is calculated from the LV equations as
⌊" close="⌋" separators="">μrαnα1∑βxβ+∑βrαxαIαβxβ,
where the symbols have their usual meaning and the outermost brackets indicate the “floor” function, i.e. we take the next smaller integer number below the value of the term enclosed by these brackets. Then these mutants are added to the community, each with an initial population of 1.
In the early stages of ecoevolutionary simulations, ITEEM and kickrelax model differ due to individuality of the agents. With the large system size and the negligible competitive dominance of a mutant to its parent, a mutant can invade very slowly. In kickrelax, it then can produce a new mutant of its own if its population reaches a large value (typically around 2000 in our simulations) to be able to compensate the small product of mutation rate and competitive ability. Taken together, this means that an adaptive character displacement takes a long time. However, during this time the parent produces numerous novel mutants. Thus, we have a slow and more or less even spread in trait space, effectively without branching. In ITEEM, on the other hand, an individual mutant has a nonnegligible chance (proportional to its reproduction rate and competitive ability) to double its population in one generation, and in the process, to produce a mutant of its own. Hence, from the early stages of the simulation onwards, adaptive forces are active.
Shtilerman et al., solve the above problem technically with an arbitrary assumption: they switch a certain percentage (e.g. 5%) of a mother strain population selected for mutation in one step to a daughter strain which is different from its parent by an arbitrary overlap reduction parameter ℎ introduced to enforce niche separation. The artificial character of the kickrelax model may also have forced Shtilerman et al., to make further assumptions, leading to more parameters with seemingly arbitrary values, such as initial population of a mutant (explained above), an extinction threshold 𝑛_{0} below which a species is removed forever, and a relax time interval 𝑇_{𝑠}between stochastic kicks.
Considering that we are trying to understand features of a complex nonequilibrium ecoevolutionary dynamics, it is not surprising that models like kickrelax that implement a number of technically motivated artificial assumptions concerning exactly this nonequilibrium dynamics will in general yield different features of this dynamics than a model that does not make these assumptions.
In contrast to the kickrelax model, the nature of ITEEM requires only a few elementary assumptions. It imposes no separation of ecological and evolutionary steps or time scales. In this conceptual sense it is more continuous than kickrelax models.
The new short paragraph added to the Discussion section explains briefly that why we adopt the individualbased approach.
2) The assumption of independent evolution of the elements of the interaction matrix needs to be discussed in more detail and clarity, and in the context of biological realism. This point was raised by both reviewers 1 and 3 (formerly reviewer 4).
To respond to the first comments of reviewers 2 and 3 (we guess “1 and 3” was a typo), we have added now introductory context, and a description and discussion of the mapping from phenotype to interaction. Moreover, we have added a substantial new Appendix13 (“Phenotypeinteraction map”) where we explain how a random variation in phenotype space yields a random variation in interaction space. However, we emphasize that despite the randomness of interaction term variations at the timepoint of their Introduction (in accordance with random variations of phenotypic traits), the evolutionary dynamics in ITEEM is adaptive, and the fate of mutants are determined by frequencydependent selection (in accordance with phenotypebased models like adaptive dynamics).
3) The "single resource" issue is related to point 2. and also needs a revised treatment.
We have taken precautions to not misdirect readers into thinking ITEEM was a consumerresource model. Please see the response to the reviewer's comment.
4) I am sympathetic with reviewer 3's (formerly reviewer 4) concern about some of the references to previous work, particularly with regard to tradeoffs.
We have exchanged or dropped references accordingly. Please see the response to the corresponding comments of the third reviewer.
5) Please also address the remaining points, e.g. the definition of "species" raised by reviewer 2, as well as other minor points.
Please find below our responses to all comments of the reviewers.
Reviewer #1:
The authors did a very good and thorough job in revising their paper. Almost all of my concerns have been addressed satisfactorily. The one remaining issue is that I don't buy the authors case for only using individualbased simulations (my original comment 2). The authors state in their rebuttal that "It has been extensively discussed in the literature that continuum approaches are unsuitable in cases of nonequilibrium dynamics,…", and they somehow conclude from this that the deterministic LotkaVolterra description (which in my mind is the same as a "logistic" description) would not be appropriate for the problem at hand. But that is exactly the question: in some general sense, one would expect that in the limit of large population sizes, the individualbased model used by the authors would converge to some deterministic model, and my guess is that this model would be at least close to the "meanfield" LotkaVolterra model. It then becomes important to understand just which features of the individualbased models can be understood by studying the much simpler meanfield model. The authors present no arguments why there even are *any* features of their individualbased model that could not be observed in the deterministic model. They refer to nearneutrality, but that's exactly one of the features that was observed in similar deterministic LotkaVolterra models by Shtilerman et al., (2015). I am not convinced that the salient results reported in this paper could not also be obtained with deterministic models. Just claiming that some results cannot be obtained in that way on general grounds does not make it true in this particular case. The obvious advantage of using deterministic models would be that such models are much more tractable analytically (e.g. from a statistical physics point of view), and it is therefore important to know how far one can get using them. I think it would not be too onerous to at least do some tests using the deterministic models to either confirm or refute the claim that they can produce similar results as the individualbased models. The question is: does stochasticity really play a major role in producing the results reported in this paper? If so, then this would be important to know (but this point would need to be made based on more than just a vague statement that their "system falls into the category of those better modeled by individual based models"). If not, then the deterministic models should do a good job reproducing these results.
As described in the response to the editor, the model of Shtilerman et al., (2015) is not a deterministic model but a kickrelax model with periodic stochastic evolutionary kicks followed by deterministic population relaxations. All previous interactionbased models including that of Shtilerman et al., (2015) are of this kickrelax type. We presume that the reasons for reviewer #2 the adoption of kickrelax models were historical: the scientists in the field were familiar with LotkaVolterra continuous models for ecological systems, and they transferred them therefore to ecoevolutionary problems, adding just periodic stochastic kicks for the evolutionary part. Not only do such kickrelax models artificially disentangle evolution and ecology and therefore are of questionable value for studying dynamics, but the artificial split necessitates a number of additional assumptions and parameters. This includes for instance the discrete extinction threshold (value not reported), and length of the relax interval (i.e. interval between two stochastic kicks) of 2000 generations. Another value that is difficult to justify is the low mutation rate of about 10^{−8} per replication. Diversification is imposed ad hoc by a discontinuous flipping of 5% of a parent population into a mutant population with artificial niche separation, described by another ad hoc parameter ℎ. The arguments given in the introduction of Shtilerman et al., (2015) explains the necessity of such discontinuous branching in the model to enforce diversity. It is noteworthy that most of these technically motivated parameters in kickrelax models do not correspond to observables. Conversely, in ITEEM the number of parameters is small, they correspond to observables, the model does not artificially split the dynamics into alternating evolutionary and ecological stages, and diversity is emerging in a simple, natural process as a function of a tradeoff that links ecoevolutionary dynamics to physical constraints.
As described in the response to the editor's first comment, we have nevertheless developed our own kickrelax model to match ITEEM as closely as possible. It is much less efficient than ITEEM and the results are qualitatively different (see response to editor's comment above). Given the artificial character of the kickrelax protocol, we find the latter not surprising.
The kickrelax model is also not suitable to answer the question about the role of stochasticity because in the kick step speciation events or mutations are introduced randomly. However, as described in the response to the editor, the discrepancy of some diversity indices between kickrelax model and ITEEM can be explained by the individuality of the agents and stochasticity of population dynamics that allow for faster adaptive character displacement in the early stage of the simulations.
Reviewer #2:
The authors have partly addressed my concerns and those of the other reviewers. I do however still have two major concerns:
1) I agree that a lowdimensional phenotype space (e.g. pertinent to exploitation of/competition for a single resource) can give rise to an Nst x Nst interaction matrix that encodes the competitive interaction between strains. However, the crucial assumption of the authors is that each term in this matrix (well, half of the terms) can vary independently. How this could come about in reality is unclear to me.
In other words, if "P" is the underlying phenotype space (solely related to consumption of/competition for the common resource) and "I" is the space of possible interaction matrices, what could the mapping f:P>I possibly look like, such that f(P) is an Nst x Nst/2 dimensional manifold?
I strongly recommend:
(a) Avoid any comparison to "single resource" models or real systems.
(b) Acknowledge early on that an important assumption of the model is that the terms of the interaction matrix (well, half of them) can in principle vary independently (i.e. are not constrained explicitly due to genetics or ecology). Whether this assumption is met in reality is an open question.
(c) Clarify in the discussion that this paper does not address the important question of how such a highdimensional interaction trait space (i.e. with Nst x Nst/2 independent axes) might arise, or provide a plausible example.
On recommendation (a): In the revision, we avoid the expressions “single” or “limiting” resource and it is clarified that ITEEM is not a consumerresource model with a single limiting resource (Discussion section). With ITEEM we provide a minimalist framework for ecoevolutionary dynamics that shows the qualitative behavior of evolutionary biological systems. On (b): Variations of the interaction matrix appear at random but only a few survive the frequencydependent selection that emerges in the system. This is akin to the principle observed by Luria and Delbrück for genotypes (Luria and Delbrück, 1943): mutations appear at random but only a few are selected. These random variations generate mutants that are ecologically similar to their parent, i.e. close to the parent trait in trait space. Please see the Appendix13on the random variations of interaction terms that are mapped from random variations in phenotypic traits. Please also see the response to the second point summarized by the reviewer. On (c): As described in the response to the second comment by the editor, we have added substantial explanations concerning the mapping of phenotype space to interaction space, including an extended example in Appendix13.
2) The definition of "species" by the authors is still confusing and of questionable relevance. The authors define "species" operationally based on a cutoff threshold in phylogenetic distance. While this is common practice in microbial ecology (where such clusters are called Operational Taxonomic Units), few would claim that the emergence and disappearance of OTUs is comparable to "speciation" dynamics in sexually reproducing organisms.
What I also found confusing is that in their "response to reviewers", the authors explain that "species" are "wellseparated nontransient clusters in trait space". This does not align with the definition provided in their manuscript (Appendix 1), where species are defined as "clusters of strains separated by longlasting gaps in a phylogenetic tree". Are these definitions equivalent in your model?
While the emergence of clusters in trait space is indeed interesting, I would recommend not calling these clusters "species", since clusters in trait space need not always be monophyletic and could in principle also consist of distantly related strains that happen to have converged in trait space.
We have revised the manuscript, especially Appendix1 (see also subsection Model), to clarify the term species. Shortly, we follow the widelyused phylophenetic species definition of Rosselló Mora, (2001), i.e. species are monophyletic clusters in trait space. This definition is biologically meaningful and applies also to asexual populations. To implement the species definition we identified longlasting monophyletic clusters of strains, i.e. strains that have the same parent (as inferred from the genealogical tree), and that have been separated from this parent in trait space for a long time (see Figures 3A and 3C of the current version and the new figure in Appendix1 for examples).
Reviewer #3:
This is reviewer 4 from the original submission again. This remains an interesting yet frustrating manuscript. The authors resisted many of the good suggestions from the other reviewers and myself in how they can place their manuscript in the broader context. In the end, it's the authors' manuscript, but I still think they could do a better job in the introduction and discussion to not confuse potential readers.
To me, the most interesting part of the manuscript is the idea that species interactions might be so high dimensional that it is best to focus on interaction traits that summarize many idiosyncratic phenotypes. This is described in the discussion but should also be highlighted more in the introduction. The relationship between phenotypic traits and interaction traits should be clarified to better address comments 1 and 67 of reviewer 2. Maybe could be described as a rugged "phenotypeinteraction map", in analogy to the idea of "genotypephenotype maps"? By the way, this is a big assumption of the model, not an established empirical fact, but still an interesting basis for the theory.
We have reorganized introduction and discussion to help the reader to see this work in the broader context. This includes also a new Figure 1 that should make more clear the idea of coarsegraining to the interactionlevel. Moreover, we have provided a substantial new Appendix13 on the genotypephenotype map. Please see also the response to the second comment of the editor.
Reviewer 2 and I were confused by statements about limiting resources and the competitive exclusion principle. In the revision, the authors still make statements like "GLV equations model competition over renewable resources" (Subsection “Model”), "we observe high diversity in a wellmixed homogenous system without violating the competitive exclusion principle" (subsection “Tradeoff anchors ecoevolutionary dynamics in physical reality”) and "we have assumed a single, limiting resource" (Subsection “Power and limitations of ITEEM”). Such statements will misdirect many readers into thinking about resource competition, R* rules, and the impossibility of coexistence of more species than resources. This is not appropriate, because in this model, the species interactions are direct (interference competition) and idiosyncratic. Allelopathy among plants or microbes would be a more relevant example than resource competition. The authors should remove all mentions about resources in the paper, because they will only confuse readers.
To avoid misdirection of readers we have in the revision stressed that the model is not a consumerresource model, we avoid the expression “limiting resource” (except for one appropriate occasion at the end of the Discussion section), and we describe in the subsection “Model” how the model deals with resources. We cannot avoid the term “resource” completely because tradeoff as a key concept of the manuscript is about options for resource allocation (competitive ability vs replication).
Concerning my previous comment 4, please don't portray other theoretical results as empirical support for your new model. Keep them clearly distinct.
We removed all the citations to theoretical models and replaced them with empirical studies, except for those that contain both empirical and theoretical results.
It's a big stretch to say claim "lifehistory tradeoffs" are a missing ingredient in existing theory (Introduction). Almost all existing ecoevolutionary theory is built around tradeoffs. Models without tradeoffs are the exception, not the rule.
We revised the Introduction accordingly, see paragraphs starting “To identify candidate mechanisms …” and “An essential component missing.
The references cited in subsection “Model” do not represent current, mainstream ecological thinking.
We agree with the reviewer, but we think that these publications are relevant in the current context and that they should be known more widely.
Reviewer #4:
I think that the authors addressed the main points from my previous review, just a few minor issues remain.
I still disagree with the definition of “strong' and “weak' tradeoff limits, seeing, for example, plots in Figure 1 as completely symmetric. In my opinion, the strongest dependence of birthrate on competitiveness happens at the line in the middle of the plot, presumably for δ=0.5. I guess it's more a terminological discussion, however, I find the definition of the strong tradeoff adopted by the Authors rather confusing.
To avoid confusion, we no longer use the term tradeoff strength in the main text. Appendix2 explains that the outcome of ecoevolutionary dynamics is not only determined by the slope of the tradeoff function (dependence of birthrate on competitiveness) but also by the community context.
Caption to Figure 2, “Disc diameter scales with total abundance of species' Does it mean that it scales with the number of individuals in a species? Or the number of species in the system? What kind of scaling is it?
The disc diameter 𝑑_{𝑖}of species 𝑖 indicates the relative abundance of the respective species:
di=N1ind∑αϵinα=∑αϵixα,
with 𝑛_{𝛼}the abundance of strain 𝛼 ∈ 𝑖, 𝑁_{𝑖𝑛𝑑}the total number of individuals in the system, and 𝑥_{𝛼}the fraction of individuals belonging to strain 𝛼.
This size of 𝑑_{𝑖}is proportional to the relative abundance of species in Muller diagram (Figure 3A in the current version). For a better visualization, all 𝑑_{𝑖}sizes are scaled with a same arbitrary factor to avoid disc overlaps or very small disc diameters. To clarify this, a short description has been added to the caption of Figure 3 (in the current version).
Subsection “Generation of diversity”, “Occasionally diversity collapses from medium levels abruptly to very low levels, usually followed a recovery”, should it read “by a recover”?
Corrected.
Appendix 1 (Eq.8) Is γ the index of summation, running from one to N_st? If not, what is the index?
Yes, and it is now corrected.
I don't know if I should dwell on that, yet I also strongly disagree with the Authors' reply to the second comment of the first reviewer, and especially, with the apparent misuse of the term “mean field'.
In short, I believe that both implementations of this process, the individualbased and continuouspopulations models should yield very similar phenomenology. Both those implementations are mean field in their nature as neither has any spatial correlations (in phenotypic or geographical space) or long temporal memory. The only difference between those is the presence of some stochasticity or noise in the individualbased model. If the Authors truly believe that such noise is the necessary source of the observed phenomenology, it should be clearly stated in the manuscript. Which, I believe, would have strongly depreciated the generality of conclusions. However, I don't think this is the case; on contrary, a continuous population version of this model would have enabled one to get “cleaner' results, speeding up the simulations and expanding the scaling range by including more species. The main distinct features that the model develops, such as temporal changes in the population, interaction cycles, speciations, etc., do not appear to be fluctuationdominated. A minor related comment, contrary to what is said int he appendix, the per capita death rate can be included into the continuous description (which is often also called the logistic model as all elements of matrix A are negative) by simply reducing the birth rate.
Regarding the usage of the term “meanfield” it is necessary to mention that we had used it in the spirit of statistical mechanics, which is different from the usage in some population biology papers where the term refers to nonspatial models (McKane and Newman, 2004). In statistical mechanics a “meanfield model” is a model that replaces all interactions between constituents (here: individuals) by an average interaction, i.e. an interaction between individuals and other components (Brézin, 2010). In this approach the cross correlations between individual degrees of freedom are neglected and the product of independent random variables can be approximated by the product of their means. Such a mean field model is also applicable to spatial models (McKane and Newman, 2004). Deterministic models at the populationlevel, like LotkaVolterra equations, logistic models and resourcecompetition models, are, in the statistical mechanics spirit, the meanfield approximation of their individualbased equivalents. To avoid confusion, we have replaced in the revised manuscript the term “meanfield” by “populationlevel”.
Regarding the minor comment about death rate, we agree with the reviewer. For a better comparison between the individualbased ecological dynamics and the presented populationlevel model we have introduced this term into Equation 6 in Appendix3.
Regarding the point about comparison of individualbased and continuouspopulations models, please see the responses to the first comment by the editor, and to the first comment by reviewer 1.
https://doi.org/10.7554/eLife.36273.030Article and author information
Author details
Funding
No external funding was received for this work.
Acknowledgements
We thank S MoghimiAraghi for helpful suggestions on the tradeoff function. We also thank the reviewers for their comments and insights, which helped us to improve the paper.
Senior Editor
 Ian T Baldwin, Max Planck Institute for Chemical Ecology, Germany
Reviewing Editor
 Michael Doebeli, University of British Columbia, Canada
Reviewer
 Yaroslav Ispolatov, Universidad de Santiago de Chile, Chile
Publication history
 Received: February 28, 2018
 Accepted: August 3, 2018
 Accepted Manuscript published: August 17, 2018 (version 1)
 Version of Record published: September 6, 2018 (version 2)
Copyright
© 2018, Farahpour et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Metrics

 3,206
 Page views

 383
 Downloads

 12
 Citations
Article citation count generated by polling the highest count across the following sources: Crossref, Scopus, PubMed Central.