## Abstract

The logistic growth model is one of the most frequently used formalizations of density dependence affecting population growth, persistence and evolution. Ecological and evolutionary theory and applications to understand population change over time often include this model. However, the assumptions and limitations of this popular model are often not well appreciated.

Here, we briefly review past use of the logistic growth model and highlight limitations by deriving population growth models from underlying consumer-resource dynamics. We show that the logistic equation likely is not applicable to many biological systems. Rather, density-regulation functions are usually non-linear and may exhibit convex or both concave and convex curvatures depending on the biology of resources and consumers. In simple cases, the dynamics can be fully described by the continuous-time Beverton-Holt model. More complex consumer dynamics show similarities to a Maynard Smith-Slatkin model.

Importantly, we show how population-level parameters, such as intrinsic rates of increase and equilibrium population densities are not independent, as often assumed. Rather, they are functions of the same underlying parameters. The commonly assumed positive relationship between equilibrium population density and competitive ability is typically invalid. As a solution, we propose simple and general relationships between intrinsic rates of increase and equilibrium population densities that capture the essence of different consumer-resource systems.

Relating population level models to underlying mechanisms allows us to discuss applications to evolutionary outcomes and how these models depend on environmental conditions, like temperature via metabolic scaling. Finally, we use time-series from microbial food chains to fit population growth models and validate theoretical predictions.

Our results show that density-regulation functions need to be chosen carefully as their shapes will depend on the study system’s biology. Importantly, we provide a mechanistic understanding of relationships between model parameters, which has implications for theory and for formulating biologically sound and empirically testable predictions.

## Introduction

Population regulation and density dependence of population growth are at the core of fundamental but also controversial research in ecology (see e.g., Turchin, 1999; Henle et al., 2004; Sibly et al., 2005; Herrando-Prez et al., 2012; Krebs, 2015; Saether et al., 2016). Density dependence of population growth is often captured by the logistic growth model (Verhulst, 1838) and its more complex extensions, such as the θ-logistic model (Gilpin and Ayala, 1973). Despite its widespread use, it is important to recall that the logistic model is an abstract description of population dynamics (Herrando-Prez et al., 2012). This level of abstraction makes the interpretation of parameters challenging and may lead to paradoxical behaviours (Ginzburg, 1992; Gabriel et al., 2005; Mallet, 2012). Such issues are especially apparent in the often used *r – K* formulation with *r*_{0} being the intrinsic rate of increase, *K* the carrying capacity and *N* the population density:

Additional challenges arise when these parameters, especially the carrying capacity (*K*), are interpreted in an evolutionary context (e.g., “*K*-selection” MacArthur, 1962). For instance, Luckinbill (1979) set out to test *r – K* selection theory using selection experiments in protist microcosms. Contrary to the expectation, he reported that r-selection actually led to higher carrying capacities, compared to the expected decrease in equilibrium population densities. Similar empirical evidence in the context of range expansions was reported by Fronhofer and Altermatt (2015) who showed that the interpretation of *K* as a parameter under selection and positively linked to competitive ability may be misleading. Recently, Reding-Roman et al. (2017) found positive *r – K* relationships in microbial systems counter to their initial hypothesis which led the authors to postulate ‘trade-ups’ and ‘uberbugs’ while discussing the relevance of these findings for cancer (see also Aktipis et al., 2013) and antibiotic resistance research (paper highlighted by Reznick and King, 2017). Although these and related issues have been discussed in detail by Matessi and Gatto (1984), Reznick et al. (2002), Rueffler et al. (2006) and Mallet (2012), to name but a few, current empirical work continues to expect negative *r – K* relationships (e.g., Fronhofer and Altermatt, 2015; Reding-Roman et al., 2017) and some theory continues to use “*K*” as an evolving trait (e.g., Lande et al., 2009; Burton et al., 2010; Fleischer et al., 2018, but see Engen and Sæther 2017).

In order to resolve some of the issues associated with the logistic growth model as described by Eq. 1, Mallet (2012), for instance, has promoted the use of Verhulst’s original r – *α* formulation of logistic growth (Verhulst, 1838, see also Kostitzin 1937). In comparison to the popular *r – K* formulation (Eq. 1), Verhulst’s model uses biologically interpretable parameters (Joshi et al., 2001; Ross, 2009; Mallet, 2012).

Namely, it includes *r*_{0} as the intrinsic rate of increase and *α* as the intraspecific competition coefficient:

From Eq. 2 it follows that the population density at equilibrium is . Similarly to Eq. 1, density dependence is assumed to act linearly, with *r*_{0} being the intercept (that is, the population growth rate when population density (*N*) is vanishingly small) and *α* representing the slope of population growth rate over population density.

Other authors have acknowledged the dynamic relationship between populations and their resources that causes density dependence by resorting to using more mechanistic consumer-resource models. For instance Matessi and Gatto (1984) show how resource dynamics and especially consumer traits have to be taken into account in order to understand density-dependent selection (see also Fronhofer and Altermatt, 2015). Such consumer-resource models provide a framework that can be used in an eco-evolutionary context (for a detailed discussion see McPeek, 2017) because model parameters linked to resource use (search efficiency, handling time) are related to real, individual-level traits that can be subject to evolutionary change (Rueffler et al., 2006; Govaert et al., 2019). Importantly, bottom-up population regulation due to renewing, depletable resources, as assumed in such consumer-resource models, is the most likely case according to Begon et al. (2006).

The disadvantage of these more mechanistic consumer-resource models is an increased complexity and number of parameters. Importantly, the quality and quantity of empirical data is often not sufficient for confronting such models to data. In an attempt to simplify, some studies have explored under what conditions population level growth models (e.g., the logistic, Eq. 1 or 2) can be used to describe the underlying consumer-resource dynamics. For instance, the consumer-resource dynamics underlying the logistic growth model have already been described by MacArthur (1970). A few years later, deriving the *r – K* logistic from the underlying resource dynamics, Schoener (1973) noticed that *r* and *K* share numerous parameters, implying that growth rates and resulting equilibrium densities may be linked through resource use traits (but see Getz, 1993). Similarly, Matessi and Gatto (1984) showed that selection for increased competitive ability (“K-selection”) does in fact not maximize equilibrium densities, but rather minimizes death rates and maximizes foraging rates and assimilation efficiencies. More recently, Abrams (2009b) compared the *θ*-logistic model (Gilpin and Ayala, 1973), an extension of the logistic model that allows for non-linear density dependence, to underlying dynamics in order to find plausible *θ* values (for an extension of this work to include multiple resources species see Abrams, 2009c). He showed that density dependence is non-linear if one assumes a Holling type II functional response for the consumer (see also Abrams, 2002), but that this non-linearity is different from the *θ* logistic model.

He finally extended his work to include a type III functional response and non-linear numeric responses. Reynolds and Brassil (2013) discussed, extended and reinterpreted these findings. Recently, O’Dwyer (2018) used the same approach to discuss the general applicability of Lotka-Volterra equations. Similar work deriving discrete-time population growth models has been conducted, for example, by Geritz and Kisdi (2004) and Brännström and Sumpter (2005). In parallel to this line of research on intra-specific density dependence, Lotka-Volterra models of inter-specific competition, and specifically the inter-specific competition coefficients (*α _{i,j}*), have been linked to resource utilization through the concept of limiting similarity for example by MacArthur and Levins (1967), Schoener (1974) and Abrams (1975).

Here, we expand on this work and use consumer-resource models to derive different forms of population growth models and to gain a better understanding of their parameters, that is, how these parameters may be interpreted in biological terms and how they are inter-related. We focus on consumers in food chains that are bottom-up regulated, as described previously. Importantly, our considerations are mechanistic as we derive consumer density-dependent population growth without assuming that the resources grow logistically in the first place to avoid circularity in the argument (see the Supplementary Material S1 for a discussion of Lakin and Van Den Driessche, 1977). We furthermore investigate relationships between parameters that may constrain evolutionary trajectories. We use these derivations to show how intrinsic rates of increase, competitive abilities and equilibrium population densities are non-independent and provide explicit relationships between those parameters. In addition, we highlight the potential of deriving climate driven population growth models by discussing examples of how to extend our results to include temperature-dependence based on the metabolic theory of ecology (Brown et al., 2004) and multiple interacting species. Finally, we confront our theoretically derived growth models with time series data from microbial model populations using Bayesian inference.

## Modelling populations using consumer-resource models

In order to derive population growth models and the density-regulation function capturing how per capita consumer population growth rates (*r*) change depending on population density (*N*), as well as to understand relationships between population-level parameters, we will use the following general consumerresource model in which *R* is the resource and *N* is the consumer population density:

In this consumer-resource model, the function *f* (*R*) captures the growth of the resources and the function *g*(*R*) captures the functional response of the consumer (Holling, 1959), that is, how much resources are harvested by consumers depending on resource density. Furthermore, the constant *e* is the assimilation coefficient which translates consumed resources into consumer offspring and the constant *d* is the consumer’s death rate.

Using a time-scale separation argument (for a critical discussion see O’Dwyer, 2018), that is, assuming that resources quickly equilibrate and solving Eq. 3a for *R*, we obtain the resource equilibrium density (piecewise defined, that is if there is no positive equilibrium). The per capita consumer dynamics (density-regulation function) then become:

In order to understand how the form of the density-regulation function (Eq. 4) depends on different consumer and resource characteristics (e.g., filter feeders, saturating feeders, respectively, abiotic or biotic resources) we explore multiple realizations of *f*(*R*) and *g*(*R*), including a chemostat model for resource growth (*f*(*R*)) as well as linear and type II (saturating) functional responses for the consumer (*g*(*R*)). For an overview of the model components used below see Table 1. For simplicity, our work does not consider the possibility of predator-dependent functional responses (for a discussion and overview see Abrams, 2014), although these can be included in principle, for instance using a Beddington-DeAngelis functional response (Geritz and Gyllenberg, 2012). Such a functional response may also be able to represent spatial variation and behavioural complexities (Cosner et al., 1999).

As a large part of the work introduced above, we stay in the realm of deterministic ordinary differential equations (ODEs), which have a long tradition in ecology and evolution. Of course, this implies that we are considering expected population sizes in continuous time systems and the equations may therefore be most appropriate for capturing the dynamics of biomass (Yodzis and Innes, 1992). Nevertheless, over finite time frames, trajectories of corresponding individual-based models are highly likely to remain close to tra jectories of the ODE whenever habitat size is sufficiently large (Kurtz, 1981). Furthermore, ODEs can be derived as a deterministic approximation of stochastic processes describing the births and deaths of individuals, where interactions emerge by the “collision” of two individuals in a well-mixed environment (see also Mobilia et al., 2007).

We start by exploring density-regulation functions that are appropriate for basal biotic consumers feeding on abiotic resources. Then, we use these density-regulation functions to describe biotic resources fed upon by higher tropic level consumers and derive density-regulation functions describing the dynamics of the latter consumers. In all cases we show how population level parameters, specifically intrinsic rates of increase (*r*_{0}) and equilibrium population densities , are impacted by and interrelated due to consumer traits such as parameters of the functional response (*g*(*R*)) as well as the assimilation coefficient (*e*) the consumer’s death rate (*d*).

## A simple case: non-saturating filter feeders

### Abiotic resources

For simplicity, we first assume that resources are abiotic, that is, resource population growth (*f*(*R*)) does not depend on resource population density, but rather on a fixed rate:

This resource model is often termed “chemostat model” with *ω* as the flow rate into and out of the system, and *R*_{0} as the resource concentration flowing into the system. The corresponding consumer dynamics are:

Chemostat models may follow different formulations and, for instance, assume no flow rate (see e.g., Abrams, 1988). This, however, does not impact our results qualitatively as the shape of the emerging density-regulation function as well as parameter relationships remain qualitatively unchanged (see Supplementary Material S2).

The amount of resources present at equilibrium can be obtained by setting Eq. 5 to zero, hence . We can now substitute into Eq. 6 to obtain the per capita growth rate of the consumer:

As Thieme (2003) notes, this result parallels the continuous-time version of the population growth model proposed by Beverton and Holt (1957) and derived by Schoener (1978) and Ruggieri and Schreiber (2005):
with *β* being the intraspecific competition coefficient in analogy to *α* in Verhulst’s model (Eq. 2). We will here always refer to the parameters scaling per capita population growth rate with population density as capturing intraspecific competitive ability.

In contrast to the logistic growth model (Eqs. 1 and 2) the density-regulation function described by Eq. 8 is not linear (see e.g., Pástor et al., 2016). Rather, it is convex, implying a decreasing strength of density regulation with increasing population density (Fig. 1; note that here and in the following a concave (convex) shape implies a (local) negative (positive) second derivative). Importantly, this implies that already for very simple consumer-resource systems the logistic growth model does not hold (but see Lakin and Van Den Driessche, 1977, and Supplementary Material S1 for a detailed discussion).

This derivation shows that non-saturating filter feeding consumer populations feeding on abiotic resources will generally have a convex per capita growth rate function (Eq. 7; Fig. 1) that is best described by a Beverton and Holt (1957) model (Eq. 8) and not by the logistic growth model given in Eq. 2. Following Jeschke et al. (2004), filter feeders that exhibit a type I functional response include branchiopods, some insect larvae, bryozoans, ascidians and molluscs, for example. Of course, the relevance of our derivation for any specific system depends on whether it fulfils relevant assumptions such that density is indeed regulated by food resources which may not be the case if densities are defined by spatial resources, for example. It is also important to note that we here assume a linear functional response for simplicity, and therefore no saturation of the consumer. Even if initially linear, functional responses of filter feeders will exhibit upper limits for high resource concentrations. Therefore, for clearly limited filter feeders, our considerations using type II functional responses below may be more appropriate.

We can also use Eq. 8 to study the relationship between population level parameters, such as the intrinsic rate of increase (*r*_{0}) and the equilibrium density of the consumer population. For nonsaturating filter feeders the intrinsic rate of increase is
and the equilibrium density is obtained as

Importantly, this shows that intrinsic rates of increase (*r*_{0}; Eq. 9), competitive abilities (; Eq. 8) and equilibrium population densities (; Eq. 10) depend on the same set of underlying parameters and are therefore not independent (see also Fig. 2). Competitive ability (*β*; Eq. 8) and intrinsic rates of increase (Eq. 9) are both linear functions of foraging rate (*a*). Note that the intrinsic rate of increase can nevertheless be independent of competitive ability, if the differences are driven by the assimilation efficiency (*e*). Counter to often made classical assumptions of a trade-off between growth rates and equilibrium densities, populations of consumers best characterized by this model will always exhibit a positive relationship between equilibrium densities and intrinsic rates of increase, regardless of whether this is due to a change in foraging rate (*a*), assimilation efficiency (*e*) or death rate (*d*; see Fig. 2). The increase will be linear for *e*, concave for *a* and convex for *d*. If density-independent mortality (*d*) is very small, the effect of foraging rate (*a*) on the equilibrium population density is negligible and only depends on the assimilation efficiency (*e*).

### Biotic resources

Since the population dynamics of filter feeding consumers of a first trophic level can be described by the Beverton-Holt model (Eq. 8), we can investigate population dynamics of the next tropic level. We start by considering a non-saturating filter feeding resource where *f*(*R*) can be described by Eq. 8 and a non-saturating filter feeding consumer where *g*(*R*) is linear. Without repeating the previously described derivation, we can use the principle of the inheritance of the curvature described by Abrams (2009b), stating that the curvature of the density-regulation function of a consumer with a linear functional response is identical to the curvature of the density-regulation function of its resource. Hence, the population dynamics of any non-saturating filter feeding consumers in a food chain with a basal abiotic resource will follow the Beverton-Holt model (Eq. 8).

## A more complex case: saturating feeders

### Abiotic resources

Up to now we have assumed a linear functional response for *g*(*R*) which has been shown to be likely only applicable to some non-saturating filter feeding organisms (Jeschke et al., 2004). While keeping resources abiotic (*f*(*R*) as a chemostat model; Eq. 5), we next explore the form of the density-regulation function, as well as how intrinsic rates of increase, competitive abilities and equilibrium population densities are linked, if the consumer follows a saturating (type II) functional response (, where *h* is the half-saturation constant).

Assuming resource equilibrium, we have to solve a quadratic equation (see Supplementary Material S3 for details). Substituting the resource equilibrium into the consumer equation results in a densityregulation function that is only slightly more complex than the Beverton-Holt model (Eq. 8; Fig. 1): with as the competitive ability, as the intrinsic rate of increase and as the equilibrium density.

While the density-regulation function may now exhibit both concave and convex sections (Fig. 3, Eq. 11 and see Supplementary Material S4) it only depends on 4 parameters. The most relevant parameter driving the extent of the concave portion of the density-regulation function is the parameter *δ* and therefore the half-saturation constant (*h* in Eq. S14 and *δ* in Eq. 11; Fig. 3C). The density-regulation function is concave at low densities only if *δ* > 0 (which is equivalent to *R*_{0} > h; see Supplementary Material S4). The smaller the half-saturation constant, the more the concave part of the density-regulation function will approach a threshold-like shape. Ultimately the foraging rate becomes independent of consumer population density (*δ* → 2) which leads to exponential growth before density regulation kicks in at high densities (see Eq. S16). By contrast, the larger the half-saturation constant is, the more the type II functional response will approach a linear shape (*δ* → – 2). The latter case approaches the Beverton-Holt model (Eq. S17; see Supplementary Material S3 for details). As discussed in Supplementary Material S5, Eq. 11 behaves in part similarly to existing density-regulation functions, specifically the continuous-time version of the Maynard Smith and Slatkin (1973) model. However, these models do not have fewer parameters and therefore do not represent sensible alternatives.

As in the simplest case of non-saturating filter feeders, intrinsic rates of increase (*r*_{0}; Eq. S12), competitive abilities (*β*; Eq. S15) and equilibrium population densities (; Eq. 12) depend on the same set of underlying parameters and are therefore not independent. All three parameters increase with the maximal foraging rate (*a*) and decrease with the half-saturation constant (*h*) and death rate (*d*). While competitive ability does not depend on the assimilation efficiency (*e*), equilibrium population densities and population growth rates increase with increasing assimilation efficiency. Therefore, as for non-saturating filter feeders, relationships between equilibrium population density and intrinsic rates of increase (*r*_{0}) will always be positive although the exact shape of the function will depend on which of the underlying consumer parameters is changing (Fig. 4).

### Biotic resources

Consumers can also exhibit saturating functional responses (*g*(*R*)), while feeding on biotic resources. In this case Eqs. S5a and S5b have to be adjusted as shown in the Supplementary Material S6 (Eqs. S27a–S27b) which yields an even more complex density-regulation function (Fig. 5).

While some simplifications are possible, the resulting density-regulation function remains unwieldy with its 7 parameters (see Eq. S37). Clearly, if such complex dynamics need to be analysed in detail, it might be more appropriate to directly rely on the consumer-resource model. However, Eq. S37 is a density-regulation function that, interestingly, behaves similarly to the well known, continuous-time version of the Maynard Smith and Slatkin (1973) model (see Supplementary Material S7 and Fig. S5):

In this model, the equilibrium density is . Via its shape exponent *γ* (Eq. 13) the Maynard Smith-Slatkin model is flexible enough (Bellows, 1981) to reproduce both the convex and concave parts of the density-regulation function (see Fig. S5 and S6) qualitatively. It is however not able to capture aspects like the asymmetry well (see Fig. S5).

Up to here, for abiotic and biotic resources and non-saturating filter feeding consumers, the equilibrium density was always a monotonically increasing function of the intrinsic rate of increase. In contrast to these simpler consumer-resource systems, our analyses show that for saturating consumers feeding on biotic resources this relationship can take many forms: unimodal, monotonically decreasing or increasing and and *r*_{0} can even be independent of each other (Fig. 6). Changes in the foraging rate (*a*) may lead to positive or negative relationships which are globally unimodal and concave (Fig. 6 A). Interestingly, changes in the assimilation efficiency (*e*) may not impact while increasing *r*_{0} (Fig. 6 A). Decreasing death rates (*d*) lead to a monotonically increasing relationship between and *r*_{0} (Fig. 6 B). Clearly, the unimodal pattern is centrally impacted by the half-saturation constant (*h*; Fig. 6 C): very low values of *h* can lead to negative relationships between and *r*_{0}. Such a negative relationship is otherwise only possible in the Lakin and Van Den Driessche (1977) model (see Supplementary Material S1).

## Confronting population growth models with data

Above, we show that many consumers that are bottom-up regulated will exhibit non-linear density regulation (Figs. 1, 3 and 5). For capturing dynamics at the consumer level, we offer alternatives to the logistic growth formulation: In the simplest case of non-saturating filter feeders and abiotic resources the dynamics of the consumer follow exactly the Beverton-Holt model (Eq. 8). For more complex cases involving saturating feeders (non-linear functional responses) or biotic resources Eq. 11 holds or the Maynard Smith-Slatkin model (Eq. 13) can be used as a heuristic description.

We will now consider how well these models perform when confronted with empirical data of a real biological system. We used populations of the freshwater ciliate model organism *Tetrahymena thermophila* as a consumer that feeds on the bacterium *Serratia marcescens*. Starting at low population densities, we recorded population growth trajectories over the course of two weeks for 7 different genotypes, each replicated 6 times. We fit the logistic, the Beverton-Holt model, Eq. 11, and the Maynard Smith-Slatkin model to these dynamics using a Bayesian approach (trajectory matching; i.e., assuming pure observation error; see Supplementary Material S8 and Rosenbaum et al. (2019) for details) to avoid the pitfalls of likelihood ridges (Clark et al., 2010; Delean et al., 2012) and compared fits using WAIC as an information criterion (McElreath, 2016).

Given our theoretical considerations, as well as the fact that the bacterial resources were replenished regularly in the microcosms (effectively mimicking abiotic resource flows), we expect either the Beverton-Holt model (BH) or Eq. 11 to best describe the dynamics and therefore produce the best fit to the time-series data compared to other candidate models of density-dependent dynamics. Indeed, over the 42 growth curves the logistic model had, on average, the worst fit (mean WAIC over all fits = 5.68; for individual results see Fig. S7 and Tab. S1), followed by the Beverton-Holt model (mean WAIC = 0.64) and the Maynard Smith-Slatkin model (mean WAIC = 0.22) while Eq. 11 performed best on average (mean WAIC = −0.16; see Fig. 7 for an example and Fig. S7 for all fits). In only 3 out of the 42 time-series the logistic was found to fit best (Fig. S7). The remaining growth curves were found to follow most often Eq. 11 (26 out of 42), the Beverton-Holt model (7 out of 42) or the Maynard Smith-Slatkin model (6 out of 42). This finding is in good accordance with work showing that *Tetrahymena may* follow a type II functional response (Fronhofer and Altermatt, 2015).

## Implications and extensions

Our derivations and the empirical support for Eq. 11 and the Beverton-Holt model clearly show the value of considering alternatives to the logistic growth model when modelling population dynamics with density-dependent rates. Even more importantly, our work highlights the underlying consumer-resource parameters being responsible for changes of and relationships between population level growth parameters (Figs. 2, 4 and 6). In the following section we discuss evolutionary consequences of our findings, and showcase extensions of our models to include temperature-dependence or multiple resource and consumer species.

### Evolutionary consequences

Besides being crucial to ecology, density dependence and resulting density-dependent selection represent a central link between ecology and evolution (Travis et al., 2013), and is essential for the occurrence of eco-evolutionary feedbacks (Govaert et al., 2019). The shape of density dependence is also crucial for understanding the consequences of adaptive prey evolution (Abrams, 2009a), which is a central topic in eco-evolutionary dynamics research (Yoshida et al., 2003; Hiltunen et al., 2014).

From an evolutionary point of view, Matessi and Gatto (1984) discuss that density-dependent selection (“*K*-selection”) should minimize equilibrium resource availability rather than maximizing (note parallels to *R** theory, Tilman 1980 and MacArthur’s minimum principle, MacArthur 1969, Gatto 1990 and Ghedini et al. 2018). As a consequence, Matessi and Gatto (1984) show that, at the consumer level, should be minimized by evolution, as this allows to reduce resource availability due to high values of *e* and/or *a* and low values of *d*. Therefore, density-dependent selection can be predicted to increase *r*_{0} (see e.g., Eq. 9) and *β*, either leading to an increase or a decrease in equilibrium densities, depending on resource and consumer behaviour (see Fig. 8). Most importantly, a negative relation does not result from a trade-off between population growth rates and competitive ability, as both parameters are always positively related. The decrease in equilibrium densities is rather a consequence of increased competition. This is an important distinction to classical *r – K* selection assumptions.

Interestingly, Getz (1993) derives *r-K* models from underlying consumer-resource models and shows that *r*_{0} and *K* may be independent if one considers a parameter capturing the maximal growth rate. This parameter only acts on *r*_{0} and may lead to a difference between *r* and *K* selected populations. Some parallels can be found in our work: for both non-saturating filter and saturating feeders, we find that for vanishingly small density-independent mortalities (*d* → 0), the equilibrium population density is independent ofthe assimilationcoefficient (*e*). Bycontrast, ifwe assume abiotic resources and vanishingly small density-independent mortalities, the equilibrium population density becomes independent of the maximal foraging rate (*a*). Only in these cases evolution in the respective parameters (*a* or *e*) can lead to a change in growth rate (*r*_{0}) without affecting the equilibrium population density .

In a classical *r-K* selection study, Luckinbill (1979) investigated the consequences of adaptation to a low-density environment. Using protists as model organisms, he showed that, in contrast to *r-K* selection theory, *r*-selection did not only lead to an increase in *r*_{0} but simultaneously to an increase in equilibrium densities . This result can be explained by our work if the protists exhibit a linear functional response (Fig. 8A–C) or feed on abiotic resources (chemostat; Fig. 8D–F). Adaptation to high-density environments may be mainly driven by changes in foraging rates (*a*) as has been prominently investigated in *Drosophila* (Muelleret al., 1991). In line with the prediction of Matessi and Gatto (1984) that selection in high-density environments should minimize , Joshi and Mueller (1988) and Mueller (1990) report that selection in high-density environments increases feeding rate (bite size). Furthermore, Joshi and Mueller (1996) find a trade-off between foraging rate (*a*) and assimilation efficiency (*e*) in *Drosophila*. Similarly, Palkovacs et al. (2011) report that Trinidadian guppies from low-predation environments which have experienced high population densities show adaptations towards increased resource consumption rates.

Most recently, Abrams (2019) has clarified that, in contrast to widely held beliefs, especially in the eco-evolutionary dynamics literature (Hendry, 2017), that adaptation in the consumer usually leads to increases in (equilibrium) population density, maladaptive population density declines can be predicted to occur regularly. In our work, adaptive population density decline should occur especially in consumerresource systems characterized by biotic resources and saturating consumers (unimodal relationships; Fig. 6). Of course it should be noted that such declines may also happen in the other systems considered here: in Figs. 2, 4 and 6 we have only considered that maximally two consumer parameter may change at once, this of course must to be the case. Additionally, the picture may be even more complicated by trade-offs or relationships between the underlying consumer parameters.

Finally, on a larger geographical scale, a widely held assumption is that organisms should be most abundant and exhibit the highest densities in their range core (Burton et al., 2010). However, the generality of this pattern has recently come under scrutiny: in an experimental evolution study, Fronhofer and Altermatt (2015) showed the evolution of lower equilibrium densities in range core populations (see also Fronhofer et al., 2017). More generally, Dallas et al. (2017) could not find increased densities in range cores across 1400 species, including vertebrates and trees. Our current work provides a potential explanation for the lack of this pattern, and even for an inverse pattern: if species feed on biotic resources, exhibit non-linear functional responses, and evolution maximizes foraging rate (*a*) rather than assimilation efficiency (*e*), low equilibrium densities will be the result (Fig. S6).

### Temperature-dependence of population dynamics

Besides having evolutionary consequences, our work has implications for instance in the context of global change research (for a concrete example see Supplementary Material S9). Linking consumer-resource parameters and population level entities explicitly is at the centre of efforts to understand how populations behave under changing climatic conditions, specifically changing temperatures (Uszko et al., 2017). Equilibrium density relates population dynamics to community and ecosystem processes, because it provides a mechanistic link between individual performance, intraspecific interactions, and the total biomass and productivity of a species. Yet, as Gilbert et al. (2014) state, understanding the temperature-dependence of the equilibrium density has long remained a knowledge gap (see also Bernhardt et al., 2018). Important progress has been made recently by Uszko et al. (2017) who for the first time include a mechanistic formulation of prey equilibrium density as a function of temperature in a predator-prey model.

### Multi-species models

While all our considerations have up to now been focused on one consumer species, natural systems rather consist of communities of multiple consumers and resources. Based on our derivations, one can expand our work to multiple species (for some examples, see Supplementary Material S10). For example, Abrams (2009c) presents an extension of Abrams (2009b) to include multiple resource species (see also O’Dwyer, 2018). Extensions to include more complex trophic and non-trophic interactions via additional mortality terms are also possible. Of course, if such complexities are relevant to a biological system, it may be useful to model these interactions explicitly using an appropriate foodweb model instead of attempting a simplification by focusing on the consumer level. This would also resolve intrinsic weaknesses of our approach due to the issues associated with the time-scale separation argument (O’Dwyer, 2018). However, our models remain relevant for example for describing the lowest level represented in such a food web. The shape of resource density dependence has for instance been shown to be crucial when understanding the consequences of adaptive prey evolution (Abrams, 2009a).

### Discussion and conclusion

Our results, based on deriving density-regulation functions from underlying consumer-resource models, show that the parameters of population level growth models (e.g., *r*_{0} and *β*) and equilibrium population densities are not independent, but are all functions of consumer traits such as foraging rate and assimilation efficiency (Figs. 2, 4, 6 and 8). This is in good accordance with Matessi and Gatto (1984) and holds true regardless of the underlying consumer-resource dynamics. In contrast to previous work, we derive a priori expectations about the relationship between growth rates, competitive abilities and equilibrium densities that capture the essence of three different consumer-resource systems: non-saturating filter feeders, saturating feeders and consumers feeding on abiotic resources (Fig. 8). Contrary to the widely held belief that growth rates and equilibrium densities are negatively related, we show that various relationships must be expected, depending on the underlying consumer-resource dynamics (Fig. 8). In accordance with previous work (e.g., Thieme, 2003; Johst et al., 2008; Abrams, 2009b; Reynolds and Brassil, 2013) we show that (i) the logistic model (Eqs. 1 and 2) assuming linear density dependence may only be appropriate under very specific conditions, such as competition for spatial resources like nesting sites or territories (see Supplementary Material S1; see also O’Dwyer 2018), and that (ii) most ecological systems will rather follow concave or convex density-regulation functions (Fig. 8), because most systems are bottom-up regulated (Begon et al., 2006). As discussed in Abrams (2009b), the curvature of these density-regulation functions is different from what the θ-logistic model (Gilpin and Ayala, 1973) can achieve. We here show that non-saturating filter feeders, and generally organisms with linear functional responses (Jeschke et al., 2004), follow a convex density-regulation function that is exactly described by a continuous-time version of the Beverton-Holt model (Eq. 8; see also Thieme 2003 and Pástor et al. 2016, for example). For consumers feeding on abiotic resource and following a saturating functional response, we provide a mechanistically derived density-regulation function (Eq. 11). More complex consumer-resource systems can be approximately described by a continuous-time version of the Maynard Smith-Slatkin model (Eq. 13).

Based on our theoretical work, we predict non-saturating filter feeders to generally exhibit a positive relationship between growth rates, competitive abilities and equilibrium population densities (Fig. 2 and 8A–C). The form of this relationship will depend on whether higher growth rates are achieved by increasing foraging rates or by increasing assimilation efficiencies (Fig. 2 and 8C). For foraging strategies that are characterized by a type II, that is, saturating functional response (Fig. 8D–F) while keeping the resources abiotic, we show that the appropriate density-regulation function can be both concave and convex (Fig. 3 and 8E). Increasing foraging rates and assimilation coefficients will still increase both growth rates and equilibrium densities (Fig. 4 and 8F), while the effect of the half-saturation constant is the opposite. The relationship between growth rates and equilibrium densities may be different for more complex systems characterized by both biotic resources and non-linear functional responses (Fig. 8G–I). Specifically, changing foraging rates and half-saturation constants may lead to non-monotonic relationships between growth rates and equilibrium densities (Fig. 6). These considerations highlight that population growth and competitive ability are related, usually positively, and that equilibrium densities are a result of underlying consumer-resource dynamics that may change positively or negatively with population growth and competitive ability. This has important implications for evolutionary considerations as discussed above. Recent microbial work has started to explore the underlying genetic basis of *r – K* relationships (Wei and Zhang, 2019).

Previous work has investigated the shape of density-regulation functions empirically, by using times-series data from growth experiments and natural population dynamics. For example, Borlestean et al. (2015) investigated the curvature of density dependence in a *Chlamydomonas* chemostat. While the authors report to be surprised by the general convexity of the density-regulation function, these results are in perfect accordance with our predictions (Fig. 8B). In a comparative study, Sibly et al. (2005) used a large dataset to show that the relationship between growth rate and density is generally convex, that is, exhibits first a fast decrease and then slows down (*θ* < 1 in the *θ*-logistic model). This result corresponds to our scenario for non-saturating filter feeders and abiotic resources (Fig. 8B), which is rather surprising, given that the dataset included mammals, birds, fish and insects. However, as Clark et al. (2010) laid out in detail, the dominance of *θ* < 1 values may be due to technical difficulties in fitting the *θ*-logistic model. Nevertheless, this study, along with the findings of Eberhardt et al. (2008) who suggested that the *θ*-logistic usually outperforms the logistic in a number of vertebrate species, clearly highlight the general non-linearity of density-regulation functions. These examples show the value of our work, as it provides theoreticians and empiricists with theoretically grounded assumptions on the occurrence of specific shapes of density dependence and trait relationships. Besides the implications discussed above, the curvature of density dependence itself is highly relevant for a population’s response to stressors. Hodgson et al. (2017) demonstrated that when density dependence is concave the effect of stressors on focal populations is always amplified, while the response can be amplified or dampened for convex density dependencies. For a detailed discussion of the impact of the shape of density dependence in basic and applied ecology see Abrams (2009c).

Of course, it is important to keep in mind that our models are simplified representations of population dynamics and ignore, for instance, age or stage structure (Mueller, 1997). We do not include Allee effects (Allee, 1931; Courchamp et al., 2008), that is, reduced population growth at low densities. There are numerous mechanisms that can lead to Allee effects and we speculate that these may lead to different functional relationships between population growth and population density. At a descriptive level, demographic Allee effects can be included as shown in Kubisch et al. (2016) for discrete time systems, for example. Our considerations also do not include time-lags in density dependence (Ratikainen et al., 2008) which are relevant in both theoretical and applied contexts, for example for population stability. For reasons of space we also do not discuss discrete time models. Note that Turchin (2003) and Thieme (2003) treat this topic in detail and mechanistic derivations of the Ricker model have, for example, been used by Melbourne and Hastings (2009).

Finally, we would like to reiterate the point made by Mallet (2012): we cover a long-discussed topic in ecology. Multiple authors have noted difficulties with the *r – K* formulation of the logistic and potential non-linearities in density dependence. Advanced textbooks like Thieme (2003) and Pástor et al. (2016) have shown that for non-saturating filter feeders and abiotic resources the continuous-time version of the Beverton-Holt model can be derived (see also Abrams, 1977; Schoener, 1978). Mallet (2012) discusses how using the *r – α* formulation of the logistic, where appropriate, alleviates some of the problems encountered with the *r – K* formulation. While these considerations are highly relevant to both empiricists and theoreticians, they seem to remain largely ignored. We hope that the insights provided here, as well as our expansion of past work towards more complex consumer-resource interactions resulting in an understanding of functional relationships between population level parameters, will help to change how density regulation and relationships between equilibrium density and population growth are treated in ecology, evolution and beyond (Aktipis et al., 2013).

## Author contributions

E.A.F. conceived the study. All authors discussed the content and topic of the study. E.A.F., L.G. and S.J.S. analysed the mathematical models. E.A.F. gathered the empirical data. E.A.F. performed the statistical analyses including fitting of population growth models. E.A.F. wrote the manuscript and all authors commented on the draft.

## Data availability

Data will be deposited in Dryad and R-code will be made available via GitHub and a Zenodo DOI.

## Supplementary Material

### S1 Derivation of the logistic model by Lakin and Van Den Driessche (1977)

Multiple authors have reported that one can derive the logistic growth model (Eqs. 1 or 2) for consumers from underlying consumer-resource dynamics (e.g. MacArthur, 1970; Thieme, 2003; Abrams, 2009; Reynolds and Brassil, 2013). Most of these derivations, however, assume that the underlying resource dynamics (Eq. 3a in the main text) follow a logistic growth model (Eqs. 1 or 2 in the main text) which, given the principle of inheritance of the curvature of density dependence reported by Abrams (2009), appears circular. Reynolds and Brassil (2013) justify their usage of the logistic for describing resource dynamics by referring to work by Lakin and Van Den Driessche (1977) who indeed derive the logistic from a chemostat model. In the Lakin and Van Den Driessche (1977) model resource dynamics (*R*) are described as follows:
with *ω* being the flow rate into and out of the system, and *R*_{0} the resource concentration flowing into the system. The consumer dynamics (*N*) are thought to follow
where *e* is the assimilation efficiency of the consumer (*N*) and a its foraging rate. As described in the main text (Eq. 4), one can assume that the resources are at equilibrium and solve Eq. S1 for *R* to obtain the equilibrium resource population density . Then the per capita consumer dynamics follow

This formulation is identical to the linear density-regulation function assumed in the Verhulst (1838) model (Eq. 2 in the main text) with setting *r*_{0} = *eaR*_{0} and (see also Fig. S1). While this result may seem satisfactory, it is worth noting that in Eq. S1 the consumer harvests a constant amount of resources (*aN*; linear functional response), regardless of resource availability. Biologically, this seems a rather unrealistic assumption.

Despite its artificiality, this model already hints at the non-independence of population growth rates (*r*_{0}) and competitive abilities (here *α*, the slope of the density-regulation function). As depicted in Fig. S1 population growth rates and competitive abilities will covary. The relationship is either linear or quadratic, depending on whether the assimilation efficiency (*e*; Fig. S1B) or the foraging rate (*a*; Fig. S1A) is the underlying driver (since r_{0} = *eaR*_{0} and ). It is also interesting to note that the equilibrium density, here defined as is not a function of the assimilation efficiency (*e*) and decreases with increasing foraging rate (*a*).

### S2 Consumer dynamics with a linear functional response and abiotic resources – an alternative chemostat model

As an alternative to Eq. 5 and following for instance Abrams (1988) we may assume the following relationship which does not include a flow rate: while keeping Eq. 6 for the consumer.

The equilibrium resource concentration is then given by . Substituting this into Eq. 6 and rearranging gives us exactly the Beverton-Holt model (Eq. 8) with and the competition coefficient . The shape of resulting relationships are summarized in Fig. S2.

### S3 Consumer dynamics with a type II functional response and abiotic resources

For a consumer (*N*) with a type II functional response feeding on an abiotic resource (*R*; chemostat model) we assume the following dynamics:

For any consumer density (*N*), the growth rate of the resource is a continuous, decreasing function of *R* and takes on the value *ωR*_{0} > 0 at *R* = 0 and becomes negative for sufficiently large *R*. Thus, there is a unique positive value at which for the given consumer density *N*. When the resource dynamics occurs at a much faster time scale than the consumer dynamics (see, e.g, Rinaldi and Muratori, 1992; Hek, 2009), this resource density corresponds to a quasi-steady state. To find this quasi-steady state, we have to solve the following equation:

Assuming that *h* + *R* = 0, we can multiply both sides with *h* + *R*. Rearranging terms results in the following quadratic equation:

This quadratic equation has two solutions:
with *D* = (*ωR*_{0} – *aN*)^{2} + 2*hω*(*ωR*_{0} + *aN*) + *h*^{2}*ω*^{2}. By uniqueness and positivity of , we have . Substituting into Eq. S5b results in

Theorem 1 of Schreiber (1996) implies that the consumer-resource system (S5a)–(S5b) always has a unique globally stable equilibrium – in particular, there are no consumer-resource oscillations. Consequently, the one-dimensional system (S8) provides a reasonable approximation of the consumer dynamics on the slower time scale for all parameter values and all initial conditions.

Equation (S8) is of the following form
with *X* = *ωR*_{0} – *aN* – *hω*. We can then use the binomial product (*a* + *b*)(*a* – *b*) = *a*^{2} – *b*^{2} to simplify this equation by multiplying the numerator and denominator with . This results in:
where . Substituting *X = ωR*_{0} – *aN – hω* and *D* = (*ωR*_{0} – *aN*)^{2} + 2*hω*(*ωR*0 + *aN*) + *h*^{2}*ω*^{2} into *X*^{2} – *D* gives

For the denominator we get . For the second term we already know that it equals to –4*hω*^{2}*R*_{0}. Thus the denominator becomes:

We thus simplify Eq. S8 to:

Now we can solve Eq. S10 to find the density *N* at equilibrium . We first need to assume that . Note that we can rewrite Eq. S10 to an expression of the form by rearranging Eq. S10 to

Now squaring both sides and solving to *N* results in:

From Eq. S10 we can calculate *r*_{0} by assuming *N* = 0:

From this we find that which we can substitute into Eq. S10. This then gives the density-regulation function for consumer dynamics with a type II functional response and abiotic resources:

We can simplify this equation by defining the competitive ability *β* as in Eq. 8 as :
with . If *δ* → 2, then the denominator of Eq. S15 simplifies to . Hence, when 1 – *βN* ≥ 0, this simplifies Eq. S15 to

This implies exponential growth as a result of *h* → 0, that is, foraging being resource independent. By contrast, if *δ* → –2, then the denominator of Eq. S15 simplifies to . This reduces Eq. S15 to
which results into the Beverton-Holt model (Eq. 8 in main text). This result is expected, because *δ* → –2 is a result of *h* → ∞ which leads to a linear functional response. Solving equation (S15) to equilibrium density then gives:

### S4 Curvature of the density-regulation function – type II functional response and abiotic resources

We will use the following general consumer-resource model
where *R* is the resource and *N* the consumer population density. The function *f*(*R*) captures the growth of the resources and the function *g*(*R*) captures the functional response of the consumer. In order to have a general idea of the form of the per-capita growth function of the consumer , we can take the first and second derivative with respect to *N* of *r*(*N*). This then equals:

Hence, if we have information on the sign of and , this will provide information on the sign of and . But as *g* is a function of *R* and *R* is a function of *N*, the first and second derivative of *g* to *N* can be rewritten as:

So we need to find expressions for and . To do so, note that at resource equilibrium , we find that

If we set this fraction to a function *v*(*R*), and implicitly differentiate to *N*, this yields:
as *R* is a function of *N*. From which it follows that
where . A second implicit differentiation yields:
which gives after rearranging terms

Eq. S23 determines the curvature of *R* as a function of *N*. Hence, Eq. S22 and Eq. S23 gives an expression for the first and second derivative of *R* in function of *N*, which we can use to substitute in the first and second derivative of *g* to *N*:

For abiotic resources with type II functional response, function *f* and *g* become:

In this case, v(R) and its first and second derivative equal:

From this it follows that and . Thus, the resource equilibrium density is a decreasing convex function of consumer population density, just as in the case of the linear functional response.

From this we can determine the curvature of the per-capita growth rate by substituting the information given in Eq. S25 into Eq. S24. We get:

From this Eq. S24 becomes

As *v*’ < 0, it follows that . Hence 2*v*’ + (*h* + *R*)*v*″ determines the sign of . Using information from Eq. S26, we can replace *v*’ and *v*″ in 2*v*’ + (*h* + *R*)*v*″ and obtain:

Now . Hence, *R* decreases with *N* from *R*_{0} at *N* = 0 until *N* is such that *R*^{3} = *h*^{2}*R*_{0}. At *N* = 0, *R* = *R*_{0}, and we find that . Note that the concavity in the density-regulation function stems from the concavity of the type II functional response of the consumer, not from the resource level itself, as the resource is a decreasing convex function of consumer population density.

### S5 Fitting known density-regulation models to consumer dynamics with a type II functional response and abiotic resources

In order to compare the behaviour of Eq. S15 to existing population-regulation models, we take the same approach as Abrams (2009) and fitted these models to Eq. S15 using a least-squares approach. The respective ODEs were solved (function ‘ode’ of the ‘deSolve’ package in R version 3.4.4) and the model was fitted by minimizing the sum of squared residuals (Levenberg-Marquardt algorithm implemented in function ‘nls.lm’ of the ‘minpack.lm’ package in R). Specifically, we here consider the Beverton-Holt model (Eq. 8), the *θ*-logistic model (Gilpin and Ayala, 1973), the Hassel model (Hassell et al., 1976) and the Maynard Smith-Slatkin model (Maynard Smith and Slatkin, 1973) as possible candidates.

Because the Maynard Smith-Slatkin model approximates the dynamics of Eq. S15 sufficiently well, we can use this model to analyse numerically how changes in consumer-resource parameters impact the parameters used in the Maynard Smith-Slatkin model (Fig. S4). We recapture the results laid out in the main text and show that strong concavity (determined by the shape paramete*r γ*) mainly occurs at low values of the half-saturation constant (*h*) and specifically for low foraging rates (*a*).

### S6 Consumer dynamics with a type II functional response and biotic resources

For a consumer-resource system with biotic resources (Eq. 8) and a consumer exhibiting a type II functional response, we can write:

Consumer-resource dynamics of this type have been studied extensively (see, e.g., Cheng et al., 1982; Kuang and Freedman, 1988). Let *α* = *r*_{0},*R*+ *d _{R}*. The consumer and resource can coexist, in the sense of permanence (Schreiber, 2000), if and only if

*α*>

*d*>

_{R}, ea*d*, and . When this coexistence occurs there is a unique equilibrium supporting both species. The stability of this equilibrium depends on the slope ofthe non-trivial resource nullcline (given by

_{N}*N*= (

*h+R*)(

*α*/(1+

*β*)–

_{R}*d*)/

_{R}*a*) at the point at which it intersects the non-trivial consumer nullcline (given by

*R*=

*d*/(

_{N}h*ea–d*). Iftheslope isnegative, then Cheng et al. (1982) prove the consumer-resource equilibrium is globallystable. If the slope is positive, then Kuang and Freedman (1988) prove that the equilibrium is unstable and there is a unique limit cycle at which the consumer and resource ultimately coexist.

_{N}If we assume that the resource dynamics occur at a faster time scale then the consumer dynamics, then we can solve for potential quasi-steady states for the resource at any consumer density *N*. Unlike the case of the abiotic resource dynamics, there does not always exist a unique, positive resource quasi-steady state. To understand when such a state does or does not exist, we need to solve for *N* which gives the non-trivial resource nullcline

Whether or not *f* (*R*) for *R* ≥ 0 is hump-shaped or just decreasing in *R*, depends on *f*’(0) = (*α* – *d _{R}* –

*αβ*)/

_{R}h*a*is positive or negative. If

*f*‘(0) ≤ 0,

*f*(

*R*) decreases for all

*N*and the maximum value of

*f*(

*R*) is

*N*=

_{max}*f*(0) =

*h*(

*α*–

*d*)/

_{R}*a*. If

*f*’(0) > 0, then there is a critical value at which

*f*’(

*R*) = 0 in which case the maximum value of

_{c}*f*(

*R*) is

*N*=

_{max}*f*(

*R*).

_{c}If *N* > *N _{max}*, solving

*dR/dt*= 0 for the quasi-steady value of

*R*always yields

*R*= 0 as the only nonnegative solution. Hence, for these values of the consumer density, the consumer is driving the resource extinct. If

*N*<

*N*, there will be at least one positive solution to solving

_{max}*dR/dt*= 0 for

*R*. Solving this equation yields two roots of a quadratic equation: with . If

*f*’(0) ≤ 0 and

*N*<

*N*, then . Hence, in the fast-slow limit, the resource density quickly approaches the quasi-steady state . The same conclusion applies when

_{max}*f*’(0) > 0 and

*N*≤

*f*(0). When

*f*’(0) > 0 and

*f*(0) <

*N*<

*N*, both quasi-steady states are positive. However, in the fast-slow limit, the resource density quickly approaches the quasi-steady state whenever its initial density is and quickly approaches 0 otherwise. are stable at the fast-time scale, while is unstable at this time scale.

_{max}Putting all of these results together, we have the following approximation for the consumer dynamics:

However, this approximation is only valid as long as *N* < *N _{max}*. When the equilibrium of the full system (S27a)–(S27b) is globally stable, this constraint holds for all

*t*providing that

*N*(0) <

*N*. Alternatively, when the equilibrium of the full system (S27a)–(S27b) is unstable, solutions

_{max}*N*(

*t*) of (S29) with

*N*(0) <

*N*will in finite time increase to value

_{max}*N*at which time the approximation is no longer valid.

_{max}Using , we can simplify Eq. S29 using similar steps as in Supplementary Material S3. Setting *X* = – *β _{R}d_{R}h* –

*β*+

_{R}aN*r*, Eq. S29 can be rewritten as follows:

_{0,R}Multiplying the numerator and denominator with and simplifying gives:

In order to find the equilibrium density , we can solve Eq. S31 for *N*, giving the solutions: and
*r*_{0} can be obtained by setting *N* = 0 in Eq. S31, which gives:

Substituting *r*_{0} into by using Eq. S33 givesg

The previous equation can then be simplified by dividing by *β _{R}d_{R}h* +

*r*

_{0},

*R*: and by setting and we can rewrite the previous equation as follows: which equals

### S7 Fitting known density-regulation models to consumer dynamics with a type II functional response and biotic resources

For details of the fitting procedure see Supplementary Material S5.

Fitting the Maynard Smith-Slatkin model to a wide range of parameter combinations highlights that relations between intrinsic rates of increase (*r*_{0}) and competition coefficients do not change qualitatively, while relations with the equilibrium density do. Up to here, for abiotic or biotic resources and non-saturating filter feeding consumers, the equilibrium density was always a monotonically increasing function of foraging rates. We here find that this relationship is unimodal or monotonically decreasing for biotic resources and saturating consumers (Fig. S6). The abrupt changes from concave to convex of the density-regulation function (see Fig. 5) lead to high values of the shape parameter (*γ* 1) of the fitted Maynard Smith-Slatkin model (Fig. S6 K–O).

### S8 Fitting population growth models to data

#### Study organism, experimental procedure and data collection

We used 7 strains of the freshwater ciliate *Tetrahymena thermophila* (strains: SB3539, B2086.2, A*III, CU438.1, A*V, CU427.4, CU428.2) originally obtained from the Cornell Tetrahymena Stock Center. These strains were kept for maintenance and during the experiments in protist medium (1 g dried *Lactuca sativa* powder in 1.6 L Volvic water) with *Serratia marcescens* as a food resource.

Growth experiments were carried out in 20 mL vials (Sarstedt) at 20 °C. 6 replicated growth curves of each *Tetrahymena* strain were started with 250 μL of protist culture (4 d old, at equilibrium) in 20 mL bacterized medium (10% *Serratia marcescens* from a 3 d old culture that had reached equilibrium). The growth experiment was followed over the course of 2 weeks and 2 mL of medium was replaced in each microcosm three times per week with freshly bacterized new medium.

Data was collected using video recording and analysis. At regular intervals (twice per day for the first two days, subsequently once per day or once every 2–3 days) the microcosms were sampled and the samples were placed on a counting slide (height: 0.5 mm) under a stereomicroscope (Perfex Pro 10) at a 2-fold magnification. Using a microscope camera (Perfex SC38800) videos were recorded for a total duration of 10s at a rate of 15 frames per second imaging a volume of 31 μL. Videos were subsequently analysed using a customized version of the ‘bemovi’ package (Pennekamp et al., 2015) in R.

#### Model fitting

We used the ‘rstan’ package in R to solve the ODEs and fit the 4 potential deterministic growth models (logistic, Eq. 1; Beverton-Holt, Eq. 8; Eq. 11; Maynard Smith-Slatkin, Eq. 13). We use trajectory matching, that is, we assume pure observation error for simplicity. For a detailed description see Rosenbaum et al. (2019). As in such a trajectory matching fitting exercise the fit depends very strongly on the first population density value (*N*_{0}) we alleviate this problem to some degree by estimating N0 as an additional parameter. To facilitate fitting, both data and model parameters (except *δ* in Eq. 11) were log-transformed, and, besides fitting the respective parameters of the ODEs we also fitted the initial density as a free parameter. Priors were chosen based on published data on *Tetrahymena* (Fronhofer and Altermatt, 2015; Fronhofer et al., 2017; Altermatt and Fronhofer, 2018). See below for the rstan model code.

#### Stan code

### S9 Temperature scaling of equilibrium densities

#### Overview

In order to analyse how the equilibrium densities change with temperature, we used metabolic theory of ecology (Brown et al., 2004) that states that a given biological rate *I* will depend on temperature as follows:
with *E _{a}* as the activation energy,

*k*the Boltzmann constant and

*T*as the absolute temperature. While this relationship has been directly applied to population growth rates and competition coefficients, the effect of temperature on the equilibrium density has been less straight-forward to infer and justify theoretically (Gilbert et al., 2014; Bernhardt et al., 2018). Similarly to Uszko et al. (2017), our work provides a mechanistic way forward, at least in the case of abiotic resources. For non-saturating filter feeders, we can take the expression for equilibrium density in Eq. 10 and substitute Eq. S38 for all biological rates, that is, the feeding rate

*a*and the death rate

*d*, except for the efficiency

*e*, which has been shown to be temperature independent (Del Giorgio and Cole, 1998; López-Urrutia and Morán, 2007). For simplicity, we assume the parameters of the abiotic resource to be temperature-independent and that all activation energies

*E*are identical. Then

_{a}This result implies that, as temperature increases, equilibrium density decreases as shown in Fig. S8B for non-saturating filter feeding consumers feeding on abiotic resources (for details of the derivation see the Supplementary Material S9).

We can apply the same thoughts to the entire density-regulation function which then behaves as depicted in Fig. S8A. The simultaneous increase of population growth rate and competition coefficient with temperature decreases the equilibrium density and makes the density-regulation function steeper around the equilibrium density.

The negative scaling of equilibrium density with temperature holds more generally, for all cases of abiotic and biotic resources and shapes of functional responses investigated here (see Supplementary Material S9). Clearly, changing the temperature scaling law captured by Eq. S38 along the lines suggested by Uszko et al. (2017), that is, to be unimodal, because some rates may decline beyond certain critical temperature values, would change these results. While an in-depth analysis is beyond the scope of our current work, this extension highlights the potential of our framework. Of course, our derivations can be used to investigate the effects of other scaling laws, such as body size scaling, which we will not analyse any further here.

#### Detailed derivation

Substituting Eq. S38 for all biological rates, except the assimilation coefficient, e, which has been shown to be temperature independent (Del Giorgio and Cole, 1998; López-Urrutia and Morán, 2007), into Eq. 10 gives the expression of temperature dependent equilibrium population density of a non-saturating filter feeding consumer feeding on abiotic resource:

Assuming that all activation energies *E _{a}* are identical, we can simplify this to

More generally, in terms of density dependence, we obtain for the capita growth rate of the consumer:

Under the same assumptions that all activation energies are identical, we can simplify the previous equation to:

In analogy, temperature dependence of the equilibrium density can be investigated for type II functional responses and abiotic resources. In this case, we can use Eq. S11 and Eq. S38 to obtain:

Similarly, we can use Eq. S32 and Eq. S38 to obtain the temperature dependence of a consumer characterized by a type II functional response feeding on biotic resources. Without repeating the above derivation we obtain .

### S10 Multi-species extensions

Our consideration below should be understood as representing examples and potential that can be developed. In general, developments should consider a systems logical consistency as detailed in Kuang (2002). We here again do not include predator-dependence of the functional responses (see Abrams, 2014, for a discussion). Further developments could also follow ideas laid out by Bastolla et al. (2005).

#### Overview

Considering our simplest case, abiotic resources and a non-saturating filter feeding consumer, we can add *n* consumer species to the system (see Ruggieri and Schreiber, 2005, for an example with two consumers). After solving for resource equilibrium equilibrium, substitution and some reparametrization (see Supplementary Material S10), we obtain the consumer dynamics for species *i*:

Based on Eq. S45, coexistence is only possible if growth and death rates are identical, because all consumers feed on one resource and do not interact in any other way. In analogy to Matessi and Gatto (1984) the species that minimizes excludes the others (see also Hofbauer and Sigmund, 1998).

If there is any interspecific interaction that is not resource-mediated, an additional interaction term may be added in analogy to a death term as follows:

This additional density-dependent mortality term allows for coexistence of more than one species (McPeek, 2012). Possible mechanisms behind the additional interaction term include cannibalism at the intraspecific level and chemical compounds (that is, allelopathy) at both the intra- and interspecific levels.

Finally, for more than one (abiotic) resource and more than one (non-saturating filter feeding) consumer we obtain for species i after simplification:
where *b _{i,r}* =

*e*

_{i,r}a_{i,r}R_{0,r}is the resource specific birth term of species

*i*and is the competition coefficient. See below for details.

#### Multiple consumers

We here consider consumer-resource dynamics, where we have two consumer species feeding on one abiotic resource (similar to Eqs. 5 and 6 in the main text). In analogy to our previous considerations, we assume the following dynamics:

Solving the resource for equilibrium we get:

We can now substitute Eq. S49 into Eq. S48b and get:

Similarly, substituting equation (S49) into equation (S48c) gives:

Hence, for a system with *n* consumer species, uniquely indexed with *i* ∈ {1, …, *n*}, we get the consumer dynamics of species *i*:
with .

#### Multiple resources

There may also exist consumer-resource systems with multiple resources. We here consider a system with two abiotic resources, and one consumer species. This results in the following system of equations:

Solving this system to resource equilibrium gives for and . Substituting both resource equilibria into Eq. S52C gives:

Hence, for a system with multiple resources, uniquely indexed with *r* ∈ {1, …, *m*}, we would find the following consumer dynamics:

However, when determining *r*_{0} by setting *N* = 0, we get: . Note that because of the different denominators we cannot easily substitute this into Eq. S53.

#### Multiple resources and multiple consumers

Solving to resource equilibrium gives and Substituting these resource equilibria into Eq. S54c gives:

Similarly, substituting and into Eq. S54d gives:

From this it follows that for a system with *m* resources and *n* consumers, the consumer dynamics of a consumer species *i* is given by the following equation:
where . This is Eq. S47 found in the main text.

## Acknowledgements

We thank Peter Abrams, François Massol as well as three anonymous referees for valuable comments on the previous versions of the text. F.A. received funding from the Swiss National Science Foundation Grant No PP00P3_I79089 and the University of Zurich Research Priority Programme “URPP Global Change and Biodiversity”. S.J.S. is supported by US NSF DMS-1716803. This is publication ISEM-YYYY-XXX of the Institut des Sciences de l’Evolution – Montpellier.