Abstract
Fitness is typically represented in heavily simplified terms in evolutionary genetics, often using constant selection coefficients, to make it easier to infer or predict how type frequencies change over time. This excludes fundamental ecological factors such as dynamic population size or density-dependence from the most genetically-realistic treatments of evolution, a problem that inspired MacArthur’s influential but problematic r/K theory. Following the spirit of r/K-selection as a general-purpose theory of density-dependent selection, we develop a new model of density-dependent selection by generalizing the fixed-density classic lottery model of territorial acquisition to accommodate arbitrary population densities. We show that, with density dependence, co-existence is possible in the lottery model in a stable environment. Inspired by natural Drosophila populations, we consider co-existence under strong, seasonally-fluctuating selection coupled to large cycles in population density, and show that co-existence (stable polymorphism) is promoted via a combination of the classic storage effect and density-regulated population growth. We also show that the only significant bias introduced by selection at the different environmental extremes of Grime’s triangle is the relative importance of competitive ability at different densities, confirming an important role for phenotypic constraints in shaping “primary strategies”.
“…the concept of fitness is probably too complex to allow of a useful mathematical development. Since it enters fundamentally into many population genetics considerations, it is remarkable how little attention has been paid to it.” — Warren J. Ewens, Mathematical Population Genetics I, 2004, pp. 276
Introduction
Evolutionary models differ greatly in their treatment of fitness. In models of genetic evolution, genotypes are typically assigned constant (or occasionally frequency-dependent) selection coefficients describing the change in their relative frequencies over time. This simplified treatment of selection facilitates explicit time-dependent treatment of genotype frequencies, and can be justified over sufficiently short time intervals (Ewens, 2004, p. 276). The emphasis of population genetics is to infer past selection, migration and demographic change given a sample of nucleotide sequences, or to predict how allele frequencies change over time based on their relative fitness effects together with population structure, genetic drift and linkage. The resulting picture of evolution excludes basic elements of the ecological underpinnings of selection, including density dependence, and how selection affects population size. This complicates the inference of past selection, because demographic changes can look genealogically very similar to selective frequency changes (Barton, 1998).
By contrast, models of phenotypic trait evolution use absolute fitness functions to describe how some traits of interest affect survival and reproduction in particular ecological scenarios (Metz et al., 1992; Diekmann et al., 2004). These fitness functions can be quite problem-specific and often only account for a few traits at a time. The emphasis here is on the conditions for invasion from low frequencies and co-existence, rather than frequency or abundance trajectories over time. For instance, adaptive dynamics uses “invasion fitness” to explore the consequences of eco-evolutionary feedbacks (Diekmann et al., 2004).
It is challenging to generalize beyond particular traits or ecological scenarios to model fundamentally different forms of selection. Perhaps this is not surprising given that fitness is such a complex quantity, dependent on all of a phenotype’s functional traits (Violle et al., 2007) and its environment. A detailed, trait-based, predictive model of fitness would be enormously complicated and situation-specific. It is therefore easy to doubt the feasibility of a simplified, general mathematical treatment of fitness (Ewens, 2004, p. 276). For example, MacArthur’s famous r/K scheme (MacArthur, 1962; MacArthur and Wilson, 1967) is now almost exclusively known as a framework for understanding life-history traits, and judged on its failure in that role (Pianka, 1970; Stearns, 1977; Boyce, 1984; Reznick et al., 2002). However, the r/K scheme’s original purpose was to extend the existing population-genetic treatment of selection to account for population density (MacArthur, 1962). Few attempts have been made to develop it further along these lines.
However, empirical trait classification studies have suggested the existence of a few “primary strategies”, reflecting broadly distinct responses to selection (Winemiller et al., 2015). Grime famously considered (Grime, 1974, 1977, 1988; Westoby, 1998) two broad determinants of population density: stress (persistent hardship e.g. due to resource scarcity or unfavorable temperatures) and disturbance (intermittent destruction of vegetation e.g. due to trampling, herbivory, pathogens, extreme weather or fire). The extremes of these two factors define three primary strategies denoted by C/S/R respectively (Fig. 1): competitors “C” excel in low stress, low disturbance environments; stress tolerators “S” excel in high stress, low disturbance environments; and ruderals “R” excel in low stress, high disturbance environments. Population persistence is not possible in high-stress, high-disturbance environments. Grime showed that measures of C, S and R across a wide range of plant species are anti-correlated, so that strong C-strategists are weak S and R strategists, and so on, creating a triangular C/S/R ternary plot (Grime, 1974). Similar schemes were proposed for insects (Southwood, 1977), fishes (Winemiller and Rose, 1992), and zooplankton (Allan, 1976). More recently, modern hierarchical clustering techniques have also revealed, with greater statistical rigor, distinct trait clusters in corals analogous to Grime’s primary strategies (Darling et al., 2012). These empirical findings suggest that functional traits contribute to fitness predominantly via a few key factors such as intrinsic reproductive rate or stress-tolerance.
Here we explore the interplay between some “key factors” of fitness in a simplified, territorial model of growth, dispersal and competition. This broadly follows the original spirit of MacArthur’s r/K scheme. More specifically, our aim is to begin merging some ecological realism into population genetics’ time-dependent, genetically-focused view of evolution. We revisit the classic lottery model of Chesson and Warner (1981), which has two features that make it well suited for this role, but one critical flaw that we rectify here.
The first feature is that the lottery representation of competition is particularly concise. Mature individuals (“adults”) each require their own territory, whereas newborn individuals (“propagules”) disperse to, and subsequently compete for, territories made available by the death of adults. Territorial contest among propagules leaves a single victorious adult per territory, the victor chosen at random from the propagules present, with probabilities weighted by a coefficient for each type representing competitive ability, akin to a lottery (Sale, 1977). By comparison, coefficients for the pairwise effects of types on each other (e.g. the α coefficients in the generalized Lotka-Volterra equations and the associated concept of “α-selection”; Gill 1974; Case and Gilpin 1974; Joshi et al. 2001), or explicit resource consumption (Tilman, 1982), are much more complicated. The second feature is the close connection between the lottery model and one of the foundational models of population genetics, the Wright-Fisher model of genetic drift, which we discuss further below.
The critical flaw of the classic lottery model is that it breaks down at low densities (few propagules dispersing to each territory), precluding density-dependent behaviour. Our first task is to analytically extend the classic lottery model to correctly account for low density behavior (sections “Model” and “Mean field approximation”).
Using our extended lottery model, we then revisit Grime’s C/S/R scheme, and evaluate the extent to which selection favors fecundity, competitive ability or lower mortality under Grime’s environments (section “Primary strategies and Grime’s triangle”). This represents a simple “sanity check” on the invasion behavior of our extended lottery model under different extremes.
We then explore some time-dependent behavior of our extended lottery model. Taking an example inspired by recent studies of rapid, seasonal evolution in Drosophila (Bergland et al., 2014), we discuss how environmental fluctuations might stabilize polymorphisms when population density is cyclical.
Model
We assume that reproductively mature individuals (“adults”) each require their own territory to survive and reproduce (Fig. 2). All territories are identical, and the total number of territories is T. Time t advances in discrete iterations, each representing the time from birth to reproductive maturity. In iteration t, the number of adults of the i’th genotype is ni(t), the total number of adults is N (t) = ∑i ni(t), and the number of unoccupied territories is U (t) = T − N (t).
We assume that the ni and T are large enough that stochastic fluctuations in the ni (“drift”) can be ignored. We derive deterministic equations for the expected change in the ni over time, leaving the evaluation of drift for future work. This is an excellent approximation when the ni are all large. However, we also do not evaluate the initial stochastic behaviour of adaptive mutant lineages while they are at low abundance. When considering new mutations, we therefore restrict our attention to begin with the earliest (lowest ni) deterministic behavior of mutant lineages (the transition to deterministic growth occurs at an abundance ni of order equal to their inverse expected absolute growth rate; Uecker and Hermisson 2011).
Each iteration, adults produce new offspring (“propagules”), mi of which disperse to unoccupied territories. We assume that adults cannot be ousted from their territories, so that mi only includes propagules landing on unoccupied territories. Propagules disperse at random over the unoccupied territories, regardless of distance from their parents, and independently of each other. There is no interaction between propagules (e.g. avoidance of territories crowded with propagules). Loss of propagules during dispersal is subsumed into mi.
In general, mi will increase with ni, and will depend on population density N. For example, if bi is the number of successfully dispersing propagules produced per genotype i adult, then the loss of propagules due to dispersal to occupied territories implies mi = bi(1 − N/T)ni, akin to Levins’ competition-colonization model (Levins and Culver, 1971; Tilman, 1994). In section “Cyclical birth and death rates” we evaluate Eq. (4) numerically using this functional form for mi, with bi assumed to be constant.
In the sections “Invasion of rare genotypes and coexistence” and “Primary strategies and Grime’s triangle”, we assume the simpler form mi = bini, with constant bi, meaning that all propagules land on unoccupied territories (a form of directed dispersal). This simplifies the mathematics without affecting the results of those sections, which only depend on the low-frequency invasion behavior of Eq. (4). Note that due to our assumption of uniform dispersal, the parameter bi can be thought of as a measure of “colonization ability”, which combines fecundity and dispersal ability (Levins and Culver, 1971; Tilman, 1994; Bolker and Pacala, 1999).
The number of individuals of the i’th genotype landing in any particular territory is denoted xi. We assume that xi follows a Poisson distribution , where li = mi/U is the mean territorial propagule density. This is approximation becomes exact when the ni are large enough that drift in ni can be ignored (Appendix A).
When multiple propagules land on the same territory, the victor is determined by lottery competition: genotype i wins a territory with probability cixi/∑j cjxj, where ci is a constant representing relative competitive ability (Fig. 2).
In the classic lottery model (Chesson and Warner, 1981), unoccupied territories are assumed to be saturated with propagules from every genotype li ≫1. From the law of large numbers, the composition of propagules in each territory will then not deviate appreciably from the mean composition l1, l2, …, lG (G is the number of genotypes present), and so the probability that genotype i wins any particular unoccupied territory is approximately cili/∑j cjlj. Let Δ+ni denote the number of territories won by genotype i. Then Δ+n1, Δ+n2, …, Δ+nG follow a multinomial distribution with U trials and success probabilities , respectively. Genotype i is expected to win cili/∑j cjlj of the U available territories, and deviations from this expected outcome are small (since T is large by assumption), giving where is the mean propagule competitive ability for a randomly selected propagule, L = M/U is the total propagule density and M =∑j mj is the total number of propagules.
There is a close connection between the classic lottery model and the Wright-Fisher model of genetic drift (Svardal et al., 2015). In the Wright-Fisher model, genotype abundances are sampled each generation from a multinomial distribution with success probabilities wini/∑j wjnj, where w is relative fitness and the ni are genotype abundances in the preceding generation. Population size N remains constant. This is mathematically equivalent to the classic lottery model with non-overlapping generations (di = 1 for all i) and wi = bici. Thus, the classic lottery model allows us to replace the abstract Wright-Fisher relative fitnesses wi with more ecologically-grounded fecundity, competitive ability and mortality parameters bi, ci and di, respectively. Since birth and death rates affect absolute abundances, this allows us to evaluate selection at different densities (after appropriate extensions are made), in an otherwise very similar model to the canonical Wright-Fisher. We therefore expect that drift in realized values of ni in our extended lottery model should be similar to that in the Wright-Fisher model, but we leave this for future work.
In our extension of the classic lottery model, we do not restrict ourselves to high propagule densities. Eq. (1) is nonsensical if even a single type has low propagule),density (li ≪ 1): genotype i can win at most mi territories, yet Eq. (1) demands cili/∑j cjlj of the U unoccupied territories, for any value of U. Intuitively, the cause of this discrepancy is that individuals are discrete. Genotypes with few propagules depend on the outcome of contests in territories where they have at least one propagule present, not some small fraction of a propagule as would be implied by small li in the classic lottery model. In other words, deviations from the mean propagule composition l1, l2, …, lG are important at low density.
We expect that a fraction p1(x1) … pG(xG) of the U unoccupied territories will have the propagule composition x1, …, xG. Genotype i is expected to win cixi/∑j cjxj of these. Ignoring fluctuations about these two expectations (due to our no-drift, large T, large ni approximation), genotype i’s territorial acquisition is given by in our extended lottery model, where the sum only includes territories with at least one propagule present. Note that unlike the classic lottery model, not all unoccupied territories are claimed each iteration, since under Poisson dispersal a fraction e−L remain unoccupied.
For the majority of this manuscript we assume that mortality only occurs in adults (Fig 2; setting aside the juvenile deaths implicit in territorial contest), and at a constant, genotypespecific per-capita rate 0 ≤ di ≤ 1, so that the overall change in genotype abundances is
This seems reasonable in the absence of disturbances; when we come to consider the effects of disturbances (Section “Primary strategies and Grime’s triangle”), we will incorporate disturbance-induced mortality in competing juveniles.
Results
Mean Field Approximation
Eq. (2) involves an expectation over the time-dependent dispersal distributions pi, and is thus too complicated to give intuition about the dynamics of density-dependent lottery competition. We now evaluate this expectation using a “mean field” approximation.
Similarly to the high-li approximation of classic lottery model, we replace the xi with appropriate mean values, although we cannot simply replace xi with li. For a genotype with low propagule density li ≪ 1, we have xi = 1 in the territories where its propagules land, and so its growth comes entirely from territories which deviate appreciably from li. To account for this, we separate Eq. (2) into xi = 1 and xi > 1 parts. Our more general mean field approximation only requires that there are no large discrepancies in competitive ability (i.e. we do not have ci/cj ≫ 1 for any two genotypes). We obtain (details in Appendix B) where and
To supplement our analytical mean field derivation, we did numerical simulations of our exact our density-dependent lottery model, and verified that Eq. (4) is a good approximation (Appendix B). Thus, Eq. (4) describes how type abundances change over time in a lottery model where population density can itself vary with time.
Comparing Eq. (4) to Eq. (1), the classic lottery per-propagule success rate has been replaced by three separate terms. The first, e−L, accounts for propagules which land alone on unoccupied territories; these territories are won without contest. The second, represents competitive victories when the i genotype is a rare invader in a high density population: from Eq. (5), Ri → 0 when the i genotype is abundant (li ≫ 1), or other genotypes are collectively rare (L − li ≪ 1). The third term, , represents competitive victories when the i genotype is abundant: Ai → 0 if li ≪ 1. The relative importance of these three terms varies with both the overall propagule density L and the relative propagule frequencies mi/M. If li ≫ 1 for all genotypes, we recover the classic lottery model (only the term remains, and Ai → 1/L).
Primary strategies and Grime’s triangle
Here we describe how Grime’s “disturbed”, “stressful” and “ideal” environments can be captured mathematically in our model, and evaluate the role of selection in shaping the corresponding “primary strategies”. To proceed, we map these verbally-defined environments to quantitative parameter regimes in our model, summarized in Fig. 3.
Ideal environments are characterized by the near-absence of stress and disturbance, so that survival is easy (d ≪ 1) and growth and reproduction can be rapid (b ≫ d). The propagule density is determined by the magnitude of b/d, and hence L ≫ 1 (Appendix C).
Disturbed environments are characterized by short bursts of high extrinsic mortality d caused by physical destruction. We assume that disturbances are equally damaging to adults and juveniles, so that only (1 − di)Δ+ni (rather than Δ+ni) territories are secured by type i each iteration. In the heavily disturbed extreme of Grime’s triangle, many such bursts (which could each be spatially localized) will occur each iteration (i.e. over the time required for a propagule to reach reproductive maturity), which we approximate as a high constant mortality rate d ≈ 1. Population persistence requires b > d/(1 − d), which implies b ≫ 1 (severe disturbance comes with a high fecundity/dispersal requirement). We assume this requirement can just be met, so that b ≈ d/(1 − d) and L ≪ 1 (Appendix C).
Stressful environments are more ambiguous, and have been the subject of an extensive debate in the plant ecology literature (the “Grime-Tilman” debate; Aerts 1999 and references therein). Severe stress inhibits growth and reproduction, so that b ≪ 1. In Grime’s view, this means that the rate at which propagules successfully develop to adulthood cannot appreciably exceed the mortality rate (b/d ≈ 1); it is just possible to cope with the environmental stress (Grime, 1974, 1977). In our model, this implies L ≪ 1.
The alternative view is that stressed environments are highly competitive. In spite of inhibited growth and reproduction, density is actually high relative to carrying capacity of the stressful environment (Taylor et al., 1990). For example, if consumable resources are scarce, we expect intense resource competition (Davis et al., 1998). Thus, even though b ≪ 1, we still have b/d ≫ 1 and L ≪ 1.
The expected direction of evolution is determined by three factors. First, there may be pleiotropy between b, c, and d, which could be due to physiological trade-offs between them. Second, adaptive mutations in one trait may occur more frequently than in others (Yampolsky and Stoltzfus, 2001). Third, selection favors some trait changes over others: b, c and d will tend to evolve in the direction of greatest fitness increase because the probability that an adaptive mutant survives its initial low abundance phase and begins to grow deterministically rather than be lost to drift is proportional to its expected absolute growth rate (Haldane, 1927; Uecker and Hermisson, 2011).
Accounting for the first two factors is beyond the scope of this article. We only consider how selection shapes the direction of evolution for b, c and d, assuming that mutations affect them independently and with comparable distributions of beneficial fitness effects. Specifically, we consider the effect of mutations that change b, c, or d by a fixed proportion bj = (1 + ϵ)bi, cj = (1 + ϵ)ci or dj = (1 − ϵ)di, where j is is a novel mutant (Appendix C). Mutations with magnitude E are assumed to occur with the same frequency for all traits.
Not surprisingly, mutations which increase c are effectively neutral when L ≪ 1 and adaptive when L ≫ 1. However, selection in our density-dependent lottery model does not favor improvement in b versus d. The same is true for proportional c-mutations cj = (1 + ϵ)ci when L ≫ 1 in ideal and stressed (HD) environments.
In disturbed environments, where Δ+ni is replaced by (1 − di)Δ+ni, the argument for equal selection on b and d is more subtle and contingent. The bulk of the high mortality d ≈ 1 is environmental and unavoidable. This implies that E mutations in the remaining genetic part of mortality have a reduced impact on total d, suggesting the primacy of b under disturbance. However, this effect is easily matched by the fact that di ≈ 1 and that disturbance mortality affects both juvenile recruitment and adult mortality. In Appendix C, we argue that b and d mutations should contrive to have comparable effects on fitness in the disturbance case (assuming that d can evolve at all) because the fitness benefit of improving d declines as d improves.
Coexistence in constant and cyclical environments
In the previous section we only considered how b, c and d should respond to selection in Grime’s environmental extremes, based on invasion fitness. Here we further explore the low frequency behavior of Eq. (4) to determine which types can coexist in a constant environment, and then consider the full time-dependent behaviour of Eq. (4) in a cyclical environment.
In a constant environment, stable coexistence is possible in our extended lottery model. A b-specialist i and c-specialist j (bi > bj, cj > ci) can co-exist because then propagule density L is frequency-dependent, and so is the importance of competitive ability (Appendix D). This is a version of the classic competition-colonization trade-off (Tilman, 1994; Levins and Culver, 1971); the competitor (c-specialist) leaves many territories unoccupied (low L) due to its poor colonization ability (low b), which the colonizer (b-specialist) can then exploit. A similar situation holds for coexistence between high-c and low-d specialists; a “competition-longevity” trade-off (Tilman, 1994). These forms of co-existence require density dependence (being mediated by L), and are not present in the classic lottery model. Coexistence is not possible between b- and d-specialists in a constant environment (Appendix D).
Now suppose that birth and death rates vary periodically with amplitude sufficent to cause large changes in population density. This example is inspired by natural Drosophila populations, which expand rapidly in the warmer months when fruit is abundant, but largely die off in the colder months. Within this seasonal population density cycle, hundreds of polymorphisms also cycle in frequency (Bergland et al., 2014). Some of these polymorphisms may be adaptive and potentially millions of years old, suggesting stable coexistence (Bergland et al., 2014; Messer et al., 2016). Selection on allele frequencies thus occurs on the same time scale as population demography, a situation vastly more complicated than classical sweeps in demographically stable populations (Messer et al., 2016).
The classical population genetic treatment of fluctuating selection suggests that environmental fluctuations do not promote coexistence. Allele frequencies are successively multiplied by relative fitness values for each environmental iteration, and so two alleles favored in different enviroments can only stably coexist if the product of fitnesses for one type exactly equals the product for the other (Dempster, 1955). Thus, stable coexistence still requires frequency-dependent selection or heterozygote advantage (as is required in a constant environment).
This classical argument overlooks two general effects that promote coexistence in fluctuating environments. The first is the storage effect, which occurs when part of the population is protected from selection (due to overlapping generations in the lottery model; Chesson and Warner 1981). The second is the bounded population size effect of Yi and Dean (2013), which occurs when each environmental cycle involves growth from low to high density, with the time spent growing each cycle dependent on the fitness of the types present.
Fig. 4a-c shows the behavior of Eq. (4) for an example where b and d cycle between zero and positive values (“summers” with rapid growth and no mortality, and “winters” with mortality and no growth). Both the storage effect (adults are sheltered from selection during the summer growth phase) and the bounded density effect (expansion to high density occurs every cycle) are operating. Two types are present, a b-specialist, which is better at rapidly growing in the summer (higher b), and a d-specialist which is better at surviving the winter (lower d). Neither type has an advantage over a full environmental cycle, and they stably coexist. This is due to a combination of the storage and bounded density effects (recall that stable coexistence between b and d specialists was not possible in a constant environment).
The classic lottery model (Eq. 1) fails to give co-existence for these parameters because expansion to carrying capacity occurs immediately at the start of the summer (Fig. 4d-f). As a result, coexistence requires that the winter survivor’s b must be about 5 times smaller than required when we properly account for the growth in the abundance of each type using Eq. (4) (keeping the other parameters the same; Fig. 4g-i). Previous models of the promotion of genetic variation via the storage effect (Ellner and Hairston Jr, 1994) similarly assume that the total number of offspring per iteration is constant, and would produce a similar error.
Discussion
It is interesting to compare the predictions of the extended lottery model with earlier approaches, such as the r/K scheme, where r = b − d is the maximal, low-density growth rate (Pianka, 1972). Confusingly, the term “K-selection” sometimes refers generally to selection at high density (Pianka, 1972), encompassing both selection for higher saturation density (MacArthur and Wilson, 1967) and competitive ability (Gill, 1974). Contrary to predictions of an r/K trade-off, empirical studies have shown that maximal growth rate at low density and the high density at which saturation occurs (measured by abundance) are positively correlated, both between species/strains (Luckinbill, 1979; Kuno, 1991; Hendriks et al., 2005; Fitzsimmons et al., 2010), and as a result of experimental evolution (Luckinbill, 1978, 1979). From the perspective of our model, this positive correlation is not surprising since the saturation density, which is determined by a balance between births and deaths, increases with b.
There is support for a negative relationship between competitive success at high density and maximal growth rate (Luckinbill, 1979), consistent with a tradeoff between r and the competitive aspect of K. This could be driven by a tradeoff between individual size and reproductive rate. To avoid confusion with other forms of ”K-selection”, selection for competitive ability has been called “α-selection” after the competition coefficients in the Lotka-Volterra equation (Gill, 1974; Case and Gilpin, 1974; Joshi et al., 2001). However, competitive success as measured by α (i.e. the per-capita effect of one genotype on another genotype’s growth rate) is only partly determined by individual competitive ability — in the presence of age-structured competition and territoriality, it also includes the ability of each genotype to produce contestants i.e. b in our model. Our c is strictly competitive ability only — as such, changes in c do not directly affect population density (the total number of territories occupied per iteration is Δ+N = U (1 − e−L), which does not depend directly on the ci). The clean separation of a strictly-relative c parameter is particularly useful from an evolutionary genetics perspective, essentially embedding a zero-sum relative fitness trait within a non-zero-sum fitness model. This could have interesting applications for modeling the impacts of intra-specific competition on species extinction, for example due to clonal interference (Gerrish and Lenski, 1998; Desai and Fisher, 2007) between c-strategists on the one hand, and b- and d-strategists on the other.
K-selection in the narrow logistic sense of selection for a greater environmental carrying capacity for given r, sometimes referred to as “efficiency” (MacArthur and Wilson, 1967), could be represented in our model by smaller individual territorial requirements. To a first approximation, two co-occurring genotypes which differ by a small amount in their territorial requirements only should have the same fitness, since the costs or benefits of a change in the amount of unocupied territory is shared equally among genotypes via the propagule density per territory L. The situation is more complicated when the differences in territorial requirements become large enough that territorial contests can occur on different scales for different genotypes. We leave these complications for future work.
Our realization of Grime’s triangle (Fig. 1) differs from approaches which identify primary strategies as trait combinations that can stably co-exist (Bolker and Pacala, 1999), referring instead to the direction and rate of ongoing adaptive trait evolution under different regimes of stress and disturbance, which is closer in spirit to Grime’s arguments (Grime, 1974, 1977). Moreover, our formulation is mathematical, in contrast to Grime’s original verbal and descriptive approach, which is a recognized hindrance to the evaluation or broader application of the C/S/R scheme (e.g. Tilman 2007). However, section “Primary strategies and Grime’s triangle” suggests that, apart from the obvious irrelevance of territorial contests at low density, selection alone has no preference for fecundity, competitive ability or longevity. One or both of the factors we set aside — pleiotropy/trade-offs and mutation bias — are needed to get true “primary strategy” trait differentiation.
Nevertheless, it is interesting to note that ruderals, which are typically thought of as high fecundity dispersers (b-specialists), may also be strongly d-selected, which while unintuitive, is consistent with our findings. An effective way to reduce d in the face of unavoidable physical destruction is to shorten the time to reproductive maturity — short life cycles are a characteristically ruderal trait. Moreover, a recent hierarchical cluster analysis of coral traits did find a distinct “ruderal” cluster, but high fecundity was not its distinguishing feature. Rather, ruderals used brood-(as opposed to broadcast-) spawning, which could plausibly be a mechanism for improving propagule survivorship in disturbed environments (Darling et al., 2012).
One potential limitation of our model as a general-purpose model of density-dependent selection is its restriction to interference competition between juveniles for durable resources (lottery recruitment to adulthood), analogous to the ubiquitous assumption of viability selection in population genetics (Ewens, 2004, p. 45). In some respects this is the complement of consumable resource competition models, which restrict their attention to indirect exploitation competition, typically without age structure (Tilman, 1982). In the particular case that consumable resources are spatially localized (e.g. due to restricted movement through soils), resource competition and territorial acquisition effectively coincide, and in principle resource competition could be represented by a competitive ability c (or conversely, c should be derivable from resource competition). The situation is more complicated if the resources are well-mixed, since, in general, resource levels then need to be explicitly tracked. It seems plausible that explicit resource tracking may not be necessary when the focus is on the evolution of similar genotypes that use identical resources rather than the stable co-existence of widely differing species with different resource preferences (Ram et al., 2016). We are not aware of any attempts to delineate conditions under which explicit resource tracking is unnecessary even if it is assumed that community structure is ultimately determined by competition for consumable resources. More work is needed connecting resource competition models to the density-dependent selection literature, since most of the former has to date been focused on narrower issues of the role of competition at low resource availability and in the absence of direct interactions between organisms at the same trophic level (Aerts, 1999; Davis et al., 1998; Tilman, 2007).
While our model can be applied to species rather than genotypes (e.g. ecological invasions), our focus is genotype evolution i.e. the change in allele frequencies over time. Our assumption that there are no large c discrepancies (section “Mean field approximation”) amounts to a restriction on the amount of genetic variation in c in the population. Since beneficial mutation effect sizes will typically not be much larger than a few percent, large c discrepancies can only arise if the mutation rate is extremely large, and so the assumption will not be violated in most cases. However, this restriction could become important when looking at species interactions rather than genotype evolution.
In the introduction we mentioned the recurring difficulties with confounding selection and demography in population genetic inference. It seems that Eq. (4) or something similar (and hopefully more analytically tractable) is unavoidable for the analysis of time-course genetic data because, fundamentally, selective births and deaths affect both abundances and frequencies, not one or the other in isolation. Moreover, some aspects of allele frequency change are intrinsically density-dependent. In the classic lottery model, which as we have seen is essentially the Wright-Fisher model with overlapping generations, bi and ci are equivalent in the sense that the number of territorial victories only depends on the product bici (see “Model”). This is no longer the case in our extension, where b and c specialists can co-exist. This “colonization-competition trade-off” is well known in the co-existence literature (Tilman, 1994). It and similar forms of “spatial co-existence” in stable environments have previously been modeled either with Levin’s qualitative representation of competition (Levins and Culver, 1971; Tilman, 1994), as opposed to the quantitative c of lottery competition, or with a more sophisticated treatment of space (non-uniform dispersal; Shmida and Ellner 1984; Bolker and Pacala 1999). In cyclical environments, polymorphisms can be stabilized by the bounded density effect, which is completely lost if there is an exclusive focus on allele frequencies (Yi and Dean, 2013). We leave the details of how our model might be applied to inference problems, including the crucial issue of its genetic drift predictions (providing a null model for neutral sites), for future work.
Acknowledgments
We thank Peter Chesson and Joachim Hermisson for many constructive comments on this manuscript. This work was financially supported by the National Science Foundation (DEB-1348262).
Appendix
Appendix A: Poisson approximation
For simplicity of presentation, we have assumed a Poisson distribution for the xi as our model of dispersal. Strictly speaking, the total number of i propagules ∑xi (summed over unoccupied territories) is then no longer a constant mi, but fluctuates between generations for a given mean mi, which is more biologically realistic. Nevertheless, since we do not consider the random fluctuations in type abundances here, and for ease of comparison with the classic lottery model, we ignore the fluctuations in mi. Instead we focus, on Poisson fluctuations in propagule composition in each territory.
In the exact model of random dispersal, the counts of a genotype’s propagules across unnocupied territories follows a multinomial distribution with dimension U, total number of trials equal to mi, and equal probabilities 1/U for a propagule to land in a given territory.
Thus, the xi in different territories are not independent random variables. However, for sufficiently large U and mi, this multinomial distribution for the xi across territories is closely approximated by a product of independent Poisson distributions for each territory, each with rate parameter li (Arenbaev, 1977, Theorem 1). Since we are ignoring finite population size effects, we effectively have T → ∞, in which case U can be only be small enough to violate the Poisson approximation if there is vanishing population turnover, and then the dispersal distribution is irrelevant anyway. Likewise, in ignoring stochastic finite population size for the ni, we have effectively already assumed that mi is large enough to justify the Poisson approximation (the error scales as ; Arenbaev 1977).
Appendix B: Derivation of growth equation
We separate the right hand side of Eq. (2) into three components Δ+ni = Δuni +Δrni +Δani which vary in relative magnitude depending on the propagule densities li. Following the notation in the main text, the Poisson distributions for the xi (or some subset of the xi) will be denoted p, and we use P as a general shorthand for the probability of particular outcomes.
Growth without competition
The first component, Δuni, accounts for territories where only one focal propagule is present xi = 1 and xj = 0 for j ≠ i (u stands for “uncontested”). The proportion of territories where this occurs is lie−L, and so
Competition when rare
The second component, Δrni, accounts for territories where a single focal propagule is present along with at least one non-focal propagule (r stands for “rare”) i.e. xi = 1 and Xi ≥ 1 where Xi = ∑j≠i xj is the number of nonfocal propagules. The number of territories where this occurs is U pi(1)P (Xi ≥ 1) = binie−li (1 − e−(L−li)). Thus where denotes the expectation with respect to , and is the probability distribution of nonfocal propagule abundances xj after dispersal, in those territories where exactly one focal propagule, and at least one non-focal propagule, landed.
Our “mean field” approximation is to replace xj with its mean in the last term in Eq. (8),
Below we justify this replacement by arguing that the standard deviation (with respect to ), is much smaller than
We first calculate . Let X =∑j xj denote the total number of propagules in a territory and xi = (x1, …, xi−1, xi+1 …, xG) denote the vector of non-focal abundances, so that p(xi) = p1(x1) … pi−1(xi−1)pi+1(xi+1) … pG(xG). Then, can be written as and so
The inner sum over xi is the mean number of propagules of a given nonfocal type j that will be found in a territory which received X − 1 nonfocal propagules in total, which is equal to . Thus, where the last line follows from .
The exact analysis of the fluctuations in ∑j ≠ i cjxj is complicated because the xj are not independent with respect to . These fluctuations are part of the “drift” in type abundances which we leave for future work. Here we use the following approximation to give some insight into the magnitude of these fluctuations and also the nature of the correlations between the xj. We replace with , defined as the xi Poisson dispersal probabilities conditional on Xi ≥ 1 (which are independent). The distinction between with will be discussed further below. The approximation gives , and where C = 1 − e−(L−li) and j ≠ k.
The exact distribution assumes that exactly one of the propagules present in a given site after dispersal belongs to the focal type, whereas assumes that there is a focal propagule present before non-focal dispersal commences. As a result, predicts that the mean propagule density is greater than L (in sites with only one focal propagule is present) when the focal type is rare and the propagule density is high. This is erroneous, because the mean number of propagules in every site is L by definition. Specifically, if L − li ≈ L ≫ 1, then the mean propagule density predicted by is approximately L + 1. The discrepancy causes rare invaders to have an intrinsic rarity disadvantage (territorial contests under are more intense than they should be). In contrast, Eq. (12) correctly predicts that there are on average nonfocal propagules because accounts for potentially large negative covariances between the xj “after dispersal”. By neglecting the latter covariences, overestimates the fluctuations in ∑j ≠ i cjxj; thus gives an upper bound on the fluctuations. The discrepancy between and will be largest when L is of order 1 or smaller, because then the propagule assumed to already be present under is comparable to, or greater than, the entire propgaule density.
Decomposing the variance in ∑j ≠ i cjxj, and using the fact that and the first term in Eq. (13) are negative because C < 1, we obtain an upper bound on the relative fluctuations in ∑j ≠ i cjxj,
Suppose that the cj are all of similar magnitude (their ratios are of order one). Then Eq. (16) is ≪ 1 for the case when L − li ≪ 1 (due to the factor of C1/2), and also for the case when at least some of the nonfocal propagule densities are large lj ≫ 1 (since it is then of order ). The worst case scenario occurs when L − li is of order one. Then Eq. (16) gives a relative error of approximately 50%, which from our earlier discussion we know to be a substantial overestimate when L is of order 1. Our numerical results (Fig. 5) confirm that the relative errors are indeed small.
However, the relative fluctuations in ∑j≠i cjxj can be large if some of the cj are much larger than the others. Specifically, in the presence of a rare, extremely strong competitor (cjlj ≫ cj′lj′ for all other nonfocal genotypes j′, and lj ≪ 1), then the RHS of Eq. (16) can be large and we cannot make the replacement Eq. (9).
Substituting Eqs. (9) and (12) into Eq. (8), we obtain where Ri is defined in Eq. (5).
Competition when abundant
The final contribution, Δani, accounts for territories where two or more focal propagules are present (a stands for “abundant”). Similarly to Eq. (8), we have where is the probability distribution of both focal and nonfocal propagaule abundances after dispersal in those territories where at least two focal propagules landed.
Again, we argue that the relative fluctuations in ∑cjxj are much smaller than 1 (with respect to ), so that,
Following a similar procedure as for Δrni, where the vector of propagule abundances is denoted x, the mean focal genotype abundance is,
For nonfocal genotypes j ≠ i, we have
To calculate the relative fluctuations in ∑j≠i cjxj, we use a similar approximation as for Δrni: is approximated by , defined as the x dispersal probabilities in a territory conditional on xi > 2 (that is, treating the xj as indepenent). All covariances between nonfocal genotypes are now zero, so that , where for j ≠ i, and where D = 1 − (1 + li)e−li, and
Similarly to Eq. (16), the RHS of Eq. (23) is ≪ 1 for the case that L ≪ 1 (due to a factor of D1/2), and also for the case when at least some of the propagule densities (focal or nonfocal) are large — provided that ci and the cj are all of similar magnitude. Again, the worst case scenario occurs when li and L − li are of order 1, in which case Eq. (23) is around 35%, which is again where the approximation produces the biggest overestimate of the fluctuations in x. Similarly to Eq. (16), the RHS of (23) will not be ≪1 in the presence of a rare, extremely strong competitor.
Combining Eqs. (18) and (19), we obtain where Ai is defined in Eq. (6).
Comparison with simulations
Fig. 5 shows that Eq. (4) and its components closely approximate our density-dependent lottery model over a wide range of propagule densities (the latter is evaluated by direct simulations of uniform random dispersal and lottery competition). Two genotypes are present, one of which is at low frequency. The growth of the low-frequency genotype relies crucially on the low-density competition term , and also to a lesser extent on the high density competition term if l1 is large enough (Fig. 5b). On the other hand, is negligible for the high-frequency genotype, which depends instead on high density territorial victories (Fig. 5d). Fig. 3 also shows the breakdown of the classic lottery model at low propagule densities.
Appendix C: Mutant invasion under different environments
In this Appendix we evaluate the invasion of novel mutants in a population with a single resident type i (such that N = ni), which is in equilibrium.
For a single type i in equilibrium (Δni = 0), we have Ri = 0, , Ai = (1 − (1 + L)e−L)/L, and Eq. (4) becomes (Alternatively, Eq. (25) can be deduced directly from Eq. (2)). This implies L ≈ bi/di when bi/di ≫ 1 and L ≪ 1 when bi/di ≈ 1. Now suppose that a novel mutant j, which is initially rare, appears in the population. Then Aj/Rj ≪ 0, lj ≈ 0 and , and so, from Eq. (4), the mutant lineage’s fitness is where since lj ≪ 1.
We now consider mutant invasion under the different environments from section “Primary strategies and Grime’s triangle”. In our representation of Grime’s stressful environment, b ≪ 1, d ≈ b and L ≪ 1, and so Eq. (26) becomes Δnj/nj ≈ bj − dj. Mutations which improve b or d by the same fraction E, such that bj = (1 + ϵ)bi or dj = (1 − ϵ)di, yield identical fitness benefits.
In our representation of Grime’s ideal environment, d ≪ 1, bi/di ≫ 1 and L ≫ 1, and so Eq. (26) becomes where the last approximation follows from the fact that cj will be of similar magnitude to ci (ignoring mutations with very large effect sizes; incidentally, Eq. (1) yields the same expression by assuming L → ∞ even though this assumption is unrealistic for a rare mutant). Similarly to Grime’s stressful environment, improving b, c, or d by a fraction E yields identical fitness benefits, even though they may have very different absolute magnitudes, because the ancestral resident type’s density-regulated birth rate bici/L exactly balances its mortality di in equilibrium. Biases in the direction of evolution of b, c or d must again follow from biases in mutation availability or trade-offs/pleiotropy. The same argument applies for the high-L interpretation of the stress regime.
To represent Grime’s disturbed environment, Δ+ni is replaced with (1 − di)Δ+ni, and a single type reaches equilibrium when L/(1 − e−L) = (1 − di)bi/di (population persistence now requires (1 − di)bi > di rather than bi > di). We assumed di ≈ 1 and (1 − di)bi ≈ di, so that L ≪ 1. A mutant’s invasion fitness is then given by
Another distinct feature of disturbed environments is that the bulk of mortality is environmental and unavoidable. That is, di = die + dig, where the environmental part die does not evolve, and the genetic part dig is affected proportionally by mutations djg = (1 − E)dig.
The effects of improving b and d are now no longer identical for given values of die and dig. A b-mutation bj = (1 + E)bi gives Δnj/nj = E(1 − di)bi ≈ Edi ≈ E, whereas a d-mutation djg = (1 − E)dig gives Δnj/nj = Edig/(1 − di). Since overall mortality is high (di ≈ 1), dg can be small compared to de and yet still make a greater contribution to invasion fitness than b. The case where dig is so small that dig ≪ 1 − di (i.e. dig is an order of magnitude smaller than 1 − di which itself an order of magnitude smaller than die) effectively represents a hard constraint that d cannot evolve. We do not deny this possibility, but this brings us back to issues of mutational bias and constraint (the first two factors controlling the direction of evolution in section “Primary strategies and Grime’s triangle”), which we do not address in this manuscript.
Thus, assuming that d can evolve in this sense, dig is either comparable to or greater than 1 − di. Yet dig cannot be appreciably greater than 1 − di for long, because selection will then strongly favor reductions in dg over increases in b, thereby bringing dig closer to 1 − di. Thus, if b and d are subject to long-term evolution subject to external degradation at comparable rates (Bertram et al., 2017), then we expect selection on b and d to be of comparable strength in disturbed environments as well.
Appendix D: Coexistence in a stable environment
To determine whether coexistence is possible in a constant environment, we check for “mutual invasion”, that is, we check that j will invade an i-dominated population, but i will also invade a j-dominated population.
We consider the case of coexistence between a b-specialist i and a c-specialist j (bi > bj, cj > ci and di = dj). Suppose that bi is so large that L ≫ 1 when i is dominant, and bj is so small that L ≪ 1 when j is dominant. Then, when j is dominant, we have Δni/ni = bi − di = bi − dj = bi − bj > 0. When i is dominant, Eq. (27) applies, where Eq. (25) implies dj = di = bi(1 − e−L)/L ≈ bi/L, and so Therefore, coexistence occurs if cj/ci is sufficiently large. The analogous argument for d- and c-specialists (di < dj with L ≫ 1 when i dominates, L ≪ 1 when j dominates, and bi = bj) gives , which again implies coexistence if cj/ci is sufficiently large.
For b-and d-specialists (ci = cj), we have Δnj/nj ≈ bjdi/bi − dj when i dominates and Δni/ni ≈ bidj/bj − di when j dominates. Thus, either i or j grows when rare, but not both, and stable coexistence is not possible in a constant environment.