Ecotone formation through ecological niche construction: the role of biodiversity and species interactions

Rapid changes in species composition, also known as ecotones, can result from various causes including rapid changes in environmental conditions, or physiological thresholds. The possibility that ecotones arise from ecological niche construction by ecosystem engineers has received little attention. In this study, we investigate how the diversity of ecosystem engineers, and their interactions, can give rise to ecotones. We build a spatially explicit dynamical model that couples a multispecies community and its abiotic environment. We use numerical simulations and analytical techniques to determine the biotic and abiotic conditions under which ecotone emergence is expected to occur, and the role of biodiversity therein. We show that the diversity of ecosystem engineers can lead to indirect interactions through the modification of their shared environment. These interactions, which can be either competitive or mutualistic, can lead to the emergence of discrete communities in space, separated by sharp ecotones where a high species turnover is observed. Considering biodiversity is thus critical when studying the influence of species–environment interactions on the emergence of ecotones. This is especially true for the wide range of species that have small to moderate effects on their environment. Our work highlights new mechanisms by which biodiversity loss could cause significant changes in spatial community patterns in changing environments.

their individual goals and needs, behave and respond as a collective. Multiple mechanisms can plausibly lead to such large-scale self-organization in ecological communities.
A first possibility is a community of purpose, when all species have shared goals. This brings to mind the picture of a superorganism [1], a network of interactions carefully arranged toward self-perpetuation. But this picture has long been contested in ecology [2].
A second possibility is a community of means -a public good or marketplace tying together many individuals with different interests. When we can identify a large-scale ecosystem function, it is often tied to some shared resource -water, energy, basic elements. The existence of such a "common currency", through which a species can affect all others in a coherent fashion, is a widespread source of collective behavior.
This second possibility has already been largely explored in ecological theory, where resource competition plays a central role [3], and can indeed lead to collective organization [6]. Public goods models often assume that all agents benefit from accumulating some common resource; these works thus focus on setting up a tension between individual and collective means to achieve that profit.
Yet, the environment is no mere stockpile of resources: it is not only used, but also constructed and transformed. The environmental state -say, the concentration of various nutrients -differentially advantages some species over others, and the changes induced by a species need not be beneficial to itself. This can prompt many complex dynamics, such as a succession of different community stages, each benefitting from the outputs of the previous one.

A1.2 Model
We have S species, each with abundance N i , growth rate r i and interference competition α ij Their carrying capacities are given by their environmental niche with a maximum of k max i when the environmental variable equals the species' optimum E = C i , and a tolerance (niche width) of T i . For numerical stability and clarity of interpretation, we can decide of an extinction threshold θ (e.g. a single individual) such that we treat smaller carrying capacities as being zero Finally, each species affects the environment at a rate m i (the species' "engineering ability"), pushing it toward some value ε i , while the environment tends to return to its baseline value B with rate µ, Perfect engineer species will have ε i = C i and always draw the environment toward their own optimum. By contrast, imperfect engineers may create an environment which is suboptimal for themselves, ε i = C i , for instance by depleting resources that they need, or accumulating harmful byproducts.
To understand the long-term consequences of these dynamics, we will first study the equilibrium

A2 Perfect coexisting engineers
Throughout this section, we assume perfect engineer species (ε i = C i ) without interference competition, α ij = 0. Then, all species (with nonzero carrying capacity) can coexist, and at equilibrium A2.1 Potential landscape and equilibria

A2.1.1 Slow environment
If the dynamics of species abundances is much faster than that of the environment r i m i , µ, we expect that species quickly reach their carrying capacity for a given environment value, N i = K i (E), and hence the dynamics of the system is given by We now show that is in fact a gradient descent dynamics, similar to where U (x) is a potential function, with the dynamics always going toward the closest minimum of U (x). Indeed, notice that Thus, where the potential takes the form: We see it has two components: a parabolic well (E − B) 2 which has a single minimum at B, and a sum of Gaussian wells created by each of the engineer species. If m i are large enough, U (E) can have local minima corresponding to these engineered wells, and if m i µ, these wells are deeper than the parabola, so engineered states are more stable than the natural state.
The effective strength of a species' long term action on the environment is thus meaning that a species can be an important ecosystem engineer either through large engineering ability m i , large maximum abundance k max i , or wide niche T i . At a given patch, as E(t) changes over time, U (E(t)) will decrease until it reaches the bottom of the local basin: See Fig. 4 for a map of the potential in a simulation and how it can be used to predict the final environmental state.

A2.1.2 Fast environment
We now consider the opposite limit, when the dynamics of the environment variable E are faster than the species'. Crucially, equilibria are independent from the relative timescale of environment and species dynamics. Therefore, equilibria must always be minima of the potential U (E), even for fast environment dynamics. However, in that case, there is no guarantee that these dynamics can be approximated by gradient descent, meaning that E(t) will not necessarily remain within the initial basin of attraction.
To move out of the initial basin, there must be some time during which the dynamics are climbing up the potential landscape, i.e. dU/dt > 0. Notice that hence dU/dt > 0 requires the second term in parentheses to be larger than the slope dU/dE of the potential. Given a potential barrier, the dynamics may climb it if: • there exist basins beyond the barrier (in the same direction as the slope, i.e. same sign of dU/dE and C i − E) that are created by species with N i > K i (E) • there exist basins in the other direction (opposite sign of dU/dE and C i − E) created by species with N i < K i (E) In other words, faster environmental dynamics will be able to climb out of shallow wells provided that there are overabundant species attracting them, or, more likely for barren initial conditions, if the species creating these wells remain at low abundance N i < K i (E) for sufficiently long. Thus, the dynamics are most likely to settle in a deep (and not necessarily wide) basin. By contrast, slow environment dynamics might favor the widest basin, which is more likely to contain the initial condition E(0) (see Discussion).

A2.1.3 Effective species interactions for fast environment
Let us set µ = 1 for simplicity (i.e. m i is measured in units of µ). Then the equilibrium environment value is If the environment quickly reaches this value for any species abundance N i , the dynamical equation for species becomes

A2.2 Effective facilitation and competition
For m i 1, we can do a series expansion and thus At equilibrium Once again, we do a Taylor expansion to get which is equivalent to the equilibirum of a Lotka-Volterra model A2.2 Effective facilitation and competition A2.2.1 When does an engineer species create an equilibrium?
Let us consider a single engineer species and ask when it can create an equilibrium with a value of E distinct from B. The equilibrium criterion is dU/dE = 0, i.e.
For simplicity, let us consider the case C 1 − B > T 1 . Then, we can approximate our question by asking whether the maximum of the right-hand term at some particular value E m is larger in absolute value than the left-hand term µ(E m − B) at that point. The maximum of the right-hand term is given by a zero of its derivative , and a sufficient condition for species 1 to create a new equilibrium is We see that it is harder to create a new equilibrium far away from the baseline environment value, i.e. when C 1 − B is large.

A2.2.2 Two species
To understand the long-term interactions between two species through their engineering capabilities, we can study the environmental variable's potential landscape U (E) and ask: does each species create a potential well (are there as many alternative stable states as there are species)?
The contribution of engineers to U (E) is a sum of Gaussian terms, where λ i is defined in (13). If we assume for now µ = 0, our question becomes: when is the sum of two Gaussians unimodal? Consider two species that are equal in every respect, save their optimum We find three cases (see Fig. 1) • Pure competition: if ∆C 2T , the two Gaussians are well separated and each species forbids the other from existing • Mixed facilitation and competition: if ∆C 2T (e.g. ∆ ∈ [2T, 8T ]), the two species allow each other to exist, and even facilitate each other to some extent (making the other's potential well deeper, and thus more likely to overcome environment inertia or other competitors). Still, the sum of Gaussians remains bimodal: there are two possible equilibria, each favoring one of the species. This is a form of moderate competition between facilitators.
• Coalescence: if ∆C < 2T , the sum of Gaussians becomes unimodal, with a peak halfway between the species optima. Now, the two species act together as a single, more influential species.
If species differ in other parameters, the same qualitative picture holds, although the weaker species (smaller λ) will need a larger ∆C to maintain its own distinct peak, rather than be absorbed in the stronger species'. There exists a general quantitative criterion for bimodality in a mixture of two arbitrarily different Gaussians [5], which we report in the next section, but it does not easily generalize to more species.

A2.2.3 Criterion for bimodality of a Gaussian mixture
According to [5], a weighted sum of two gaussians is bimodal if µ > µ 0 with and p ∈ [p 1 , p 2 ] where 1 with y 1 and y 2 the roots of the equation with 0 < y 1 < y 2 < µ. Otherwise, the sum is unimodal.

A2.2.4 Many species
For many species, no exact results exist but we can provide a scaling estimate of the average number of alternative states. Given S the number of species, with their optima distributed over interval [0, L], the average distance between their optima is If the optima are uniformly distributed, the number of optima within a certain interval follows a Poisson distribution, and thus, the probability that ∆C > 2T (avoiding coalescence) is the probability of having no optima within a span of 2T , i.e.
Hence, the typical number of clusters of coalesced species scales like This approximation will only hold up to ∆C ∼ T , i.e. S ∼ L/2 T , after which adding more species will typically not contribute more equilbria. Now recall that species i on its own can create an alternative equilibrium despite the natural environmental dynamics only if meaning that the potential well created by the species is deep enough to compensate the recovery of the environment, which gets faster as E moves away from B. This means that species can only contribute to an equilibrium if their optimum falls within a range L ≤ min L, 2 mkT µ .
Only the fraction L /L of species clusters with optima within that range can create new equilibria.

A3 Interference competition and imperfect engineers
Let us now consider the case of direct interference competition α ij = 0 and imperfect engineers ε i = C i . Now, some species can go extinct even while their carrying capacity is nonzero. Recall the equilibrium condition

Supporting Information
The abundance of surviving species N * i = 0 is given by or in vector form, where I is the identity matrix, α * is the matrix of interactions restricted to the S * surviving species, and K * (E) the vector of carrying capacities of surviving species. For convenience, define the matrix

A3.1 Equivalence
Taking once again the limit of fast species dynamics, we now have which we can rewrite as Thus, we see that direct competition appears equivalent, in terms of its equilibrium effect, to imperfect engineering with effective values of species engineering capability m i and target environment value ε i . An important consequence is that even perfect engineers (ε i = C i ) will behave like imperfect ones if they also interact directly. It is, however, crucial to note that the calculation above involves summing only on the S * species that survive the direct competitive interaction. While the matrix V ij depends only on interactions α ij and not on the environment, it does depend on who survives, which is controlled by the carrying capacities as well. Thus, direct competition cannot simply be replaced by imperfect engineering, except in the regime where all species coexist (i.e. for weak direct interactions). In particular, if direct competition allows for alternate equilibria (mutual exclusion), each will correspond to a different equation (42).
Perhaps counter-intuitively, it is the action of species i on others that appears in the effective parameters above. For instance, if direct interactions are weak, and we see that a species' effective engineering abilitym i decreases due to its competitive effect on others, α ji , while its effective target environmental valueε i moves away from the optima ε j of the species it affects.

A3.2 Skewed potential
From equation (42), we can construct the corresponding potential by noticing that The last term is the only one that differs significantly from the equation in the case of perfect engineers, (7). Thus, we can separate the resulting potential into two contributions: first, the usual potential for perfect engineers, obtained here with the S * surviving species and effective engineering ratesm i , and second, a correction ∆U (E) coming from the last term above. We write where is the contribution that is similar to the perfect engineer case. Since K i (E) is Gaussian, its integral is an error function, and the correction to the potential takes the form Error functions are sigmoidal and comprised between 0 and 1, so ∆U (E) will have the general shape of a "staircase", i.e. a sum of step-like functions going up or down, with a step height of This could create new potential wells, if two (or more) species push E in the direction of each other's optimum, giving rise to a new type of interaction: obligate facilitation (or mutual stabilization), where each species degrades its environment from its own perspective, but improves it from the perspective of the other.
Ifε i − C i has the same sign for many species (e.g. all species tend to degrade complex sugars into simpler ones, pushing the environment variable E in a constant direction), the effect will be to create a general slope in that direction, and thus, dynamics akin to succession.

A4.1 Ecotones on an environmental gradient
Let us assume an environmental (e.g. latitudinal or altitudinal) gradient, reflected in the fact that the baseline value of the environmental variable B = B(x) now depends on position x along the gradient.
If we choose a solution of (5) and follow it along the gradient, as we progressively change parameters such as B(x), we do not expect sharp ecotones (transition zones between communities with different species compositions and abundances). The only option for a singular transition is to have alternative stable states, with the transition occurring when one of the states loses its stability.
Whenever two or more attractors exist for the same patch x, see black lines on Fig. 4 top-right, there is potential for hysteresis.  which is also the initial state at each point. Increasingm/µ from 0.1 (a) to 10 (b) to 100 (c), we go from an environment that closely follows B, and thus a continuous turnover of species, to the existence of alternate stable states engineered by these species, and separated by sharp transitions. Finally, reducing the rate of species dynamics from r = 10 (c) to r = 0.1 (d) increases the difference between the pure gradient descent prediction (red line, right panels) and the observed environment value E at each position x. Possible equilibria are independent of r, but which equilibrium is reached does depend on it: as we explain in Fig. 4, high r entails gradient descent toward the closest equilibrium, while low r allows the dynamics to climb a potential wall, and tends to favor the deeper basins rather than the ones closest to the initial condition.

A4.2 Succession trajectories
In this model, including direct competition and imperfect engineering, one can imagine three different succession scenarios. The reference scenario is the usual competition-colonization tradeoff setting [7]. In that case, species at later stages are expected to have slower growth but stronger competitive ability, either through direct competition, engineering ability, or a combination of both. As a consequence, later stages will be longer, but transitions between stages will also be slower.
Another scenario is succession driven by the environmental variable slowly descending down the potential landscape created by the engineer species. In that case, there is no implication that later stages will be longer. If the landscape is shaped by perfect engineers, we can expect a rather smooth change of the environmental variable. If it is shaped by imperfect engineers, it is possible to have long transition periods where all species have low abundance, separating shorter periods where a set of imperfect engineers dominates. There may be stabilization at low abundances if two sets of engineers are pushing the environment toward each other.
The third scenario is perturbation-driven succession: jumps between equilibria, either due to random noise, or to directed perturbations (e.g. a gradual increase in the baseline environmental value B). It is only in this scenario that succession will generically exhibit discrete stages separated by sudden transitions. Under random perturbations, succession will proceed on average toward deeper wells, and thus later stages will be longer on average (deeper wells resist perturbations longer), but there may be reversions to earlier stages.

A5.1 Slow and fast environment
An intuitive aspect of the speed of environment change is its inertia upon removal of some engineer species. In one limit, the environmental state may remain the same for long times (for instance, peat created by Sphagnum mosses can remain for thousands of years [8]), long enough that the species could potentially recolonize at the same abundance before any significant change occurred. In the other limit, the environmental state may revert suddenly, even instantaneously when the engineering results from physical properties of the species themselves (e.g. shielding of light by the canopy).
As noted in Sec. A2, when the environment dynamics are slow, the environmental variable effectively follows a gradient descent. This means that the species (or group of species) creating the widest basin of attraction control the dynamics for a broad range of initial conditions. By contrast, when the environment dynamics are fast, they are drawn toward the optimum of the species with the largest carrying capacity and best engineering abilities, even if its niche is narrow. We thus predict a prevalence of generalists in slow environments, and specialists in fast environments.

A5.2 Multiple environmental variables
A single environmental variable E may not suffice to accurately represent the ways in which species interact through modifications of their surroundings. On the other hand, our modelling approach for ecosystem engineers is most relevant if the number of environmental variables is limited, and small compared to the number of species -otherwise, it may be simpler to directly model pairwise species interactions.
There is one important qualitative feature that distinguishes the outcomes of this model, and those of a model with multiple environmental "dimensions": with only one dimension E, there can be at most as many equilibria as there are engineer species. This stops being the case with more environmental variables. A corresponding mathematical result states that a mixture of Gaussian components can have more maxima than components in dimension d > 1 [4,9]. This means that some equilibria could not be assigned to, nor expected from, the action of any given species on its own.