Spatial dynamics of synthetic microbial mutualists and their parasites

A major force contributing to the emergence of novelty in nature is the presence of cooperative interactions, where two or more components of a system act in synergy, sometimes leading to higher-order, emergent phenomena. Within molecular evolution, the so called hypercycle defines the simplest model of an autocatalytic cycle, providing major theoretical insights on the evolution of cooperation in the early biosphere. These closed cooperative loops have also inspired our understanding of how catalytic loops appear in ecological systems. In both cases, hypercycle and ecological cooperative loops, the role played by space seems to be crucial for their stability and resilience against parasites. However, it is difficult to test these ideas in natural ecosystems, where time and spatial scales introduce considerable limitations. Here, we use engineered bacteria as a model system to a variety of environmental scenarios identifying trends that transcend the specific model system, such an enhanced genetic diversity in environments requiring mutualistic interactions. Interestingly, we show that improved environments can slow down mutualistic range expansions as a result of genetic drift effects preceding local resource depletion. Moreover, we show that a parasitic strain is excluded from the population during range expansions (which acknowledges a classical prediction). Nevertheless, environmental deterioration can reshape population interactions, this same strain becoming part of a three-species mutualistic web in scenarios in which the two-strain mutualism becomes non functional. The evolutionary and ecological implications for the design of synthetic ecosystems are outlined.


Introduction
The evolution of complexity is largely grounded in the emergence of new forms of cooperation capable of holding together higher-order entities from simpler ones. Cooperative interactions have played a great role in the so-called major transitions in evolution [1]. Cooperation pervades the rise of molecular systems capable of overcoming mutation thresholds, multicellular assemblies incorporating division of labour or the appearance of insect societies. Each of these structures incorporates new properties that cannot be observed at the level of its component parts. Despite the burden involved in sustaining the new, larger entity, the advantage of staying together can overcome, under some circumstances, the cost of the association.
Cooperation can be achieved in particular by means of closed catalytic loops. Mutualistic interactions pervade ecological communities at many different scales, from bacterial communities to microbiomes and large-scale ecosystems [2]. The presence of these reciprocal relations was already outlined by Charles Darwin in one of his memorable studies on the ecology of earthworms [3,4] and summarised by the diagram of fig 1a. Here earthworms improve soil porosity and organic content that helps plants to grow, which results in more organic matter and mechanisms of soil preservation (which favours the earthworm population). This is a simple, two-component (n = 2) diagram, but ecosystems are characterised by the presence of multiple feedback loops and thus interactions might be more complex, like the three-member (n = 3) loop shown in (fig 1b). Here vegetation is grazed by animals, whose activity enhances the survival of invertebrates, which in turn improve soil quality thus favouring plant growth. Because of their ecological and evolutionary relevance, cooperative interactions have also been a major topic in synthetic biology [5][6][7][8][9][10]. The possibility of engineering de novo cooperative interactions is of relevance for several reasons. On one hand, engineered mutualisms could be used to build desirable (even optimal) functionalities that require the presence of a tight metabolic dependence [11,16]. Moreover, the possibility of designing mutualistic interactions and even symbiotic pairs [13][14][15][16][17] provides a unique opportunity for exploring the emergence of cooperation in evolution under a 'synthetic" perspective [18].
Mutualistic interactions are also required to sustain stable communities, particularly when harsh conditions are present. An example (fig 1c) is provided by drylands [19] and in particular the interactions between the so-called biological soil crust (BSC) and vascular plants [20]. The BSC defines in itself a complex ecosystem enclosed within a Here we indicate in (a) the mutual support between vegetation (grasses) and earthworms and in (b) a more complex hypercycle composed by vegetation, cattle and earth worms (and other invertebrates). In (c) the image shows a small area within a semiarid ecosystem including a plant surrounded by biological soil crust. Formal models of these types of interactions are described by hypercycles. In (d) we display the basic logical scheme of interactions for a two-member hypercycle. In (e) we show an extended model where a parasitic species (colour circle) takes advantage of one of the species but gives no mutual feedback. These parasites can easily destroy the hypercycle, but this effect is reduced or suppressed under the presence of oscillations and spatial diffusion when spiral waves get formed (f). Here different colours indicate different molecular species in a n = 8 member hypercycle (adapted from Attolini and Stadler 2011). few centimetres of the topsoil, largely controlling the energy and matter flow through the soil surface, helping vegetation thrive under semiarid conditions. The soil microbiome plays a major role in sustaining plant diversity and its dynamics, with the latter often completely dependent on their microbial symbionts [21]. Since these ecosystems might experience sudden declines due to climate change [22,23] understanding their dynamics is crucial to predicting their future. In this context, it has been suggested that engineering new synthetic mutualistic loops in endangered ecosystems could help prevent catastrophic shifts [24,25].
Understanding cooperation, its rise and fall and how can it overcome competitive interactions is an important problem. A great insight has been obtained from both field and theoretical studies [2]. An elegant description of this class of cooperative loops is the hypercycle, first suggested within the context of prebiotic evolution [26-28, 30, 31]. Here a simple catalytic system is defined (as in figs 1a-b) forming a closed graph where the replication of each component is catalysed by a previous one in the loop, while it also catalyses the replication of the next. The simplest case is the one shown in fig 1d for a two-member syste [26,29]. If we indicate by Φ 1 and Φ 2 their population sizes, a pair of coupled equations allows us to represent the hypercycle model as follows: where K stands for the carrying capacity of the system, δ is the degradation/death rate of both species and the replication rates of the cross-catalytic loop are indicated by α ij . As defined, we can see that no proliferation of any of the two partners will occur in the absence of the other, as a consequence of the second-order kinetics that requires the product of the two concentrations. The hypercycle can outcompete other non-cooperative species [26,28] but a major drawback is that it can also be easily threatened by a parasite (figure 1e) capable of destabilising the whole system [32]. Interestingly, mathematical and computer models indicate that this problem can be limited by the presence of diffusion in a spatial domain [33][34][35][36]. Hypercycles displaying spatial structures ( fig. 1f) are obtained from n > 4 loops capable of exhibiting oscillations. In a nutshell, the spatial structure imposes a limitation to the spread of the parasite, and it can even go extinct if the inaccessibility of its target species, combined with its death rate, makes it non-viable. This suggestion has received considerable attention [37]. Although these models are reasonable for molecular systems, real populations of microbial mutualists spread in space as expanding fronts that impose particular constraints which can be analysed using front propagation theory [38,39]. It is in this context that previous theories on hypercycle propagation and parasitic interactions can be tested using synthetic ecologies.
In this paper, we address this problem by using an experimental design where populations of synthetic microbes are forced to cooperate as they expand on a two-dimensional substrate. In this context, we take into account the population spreading that leads to complex spatial structures, partially due to the cooperative loop but also to the physical impact of cellular shapes. Both populations had some potential for (Malthusian) growth in the absence of the mutualistic partner, provided that the necessary metabolites are present. Such response allows testing the conditions under which the hypercycle overcomes the effects of competition derived from Malthusian populations exploiting similar resources. As will be shown below, the contact surface between the engineered strains increases with the strength of mutualistic interactions. When the mutualistic interactions are neutralised, the synthetic strains display the same segregative dynamics described by competing invaders [52,53]. Moreover, a synthetic parasite was also designed to test the capacity of the spatial synthetic hypercycle to prevent it from spreading. Additionally, our parasitic strain has been engineered to degrade ampicillin from the medium, in order to explore the boundaries between parasitic and cooperative interactions,

Results
Our model system to study mutualistic interactions is composed of the pair of bacterial engineered strains shown in Fig. 2a. The Istrain (depicted in yellow) cannot produce the isoleucine (iso) amino acid but overproduces and leaks leucine (leu), while L -(in blue) cannot produce leu but overproduces and leaks iso [6]. Therefore, the strains are able to engage in a cross-feeding mutualism that permits growth in coculture, in a minimal medium lacking both amino acids where neither Inor Lcan grow in monoculture (obligate mutualism scenario in Fig. 2b). However, both Iand Lare able to grow in monoculture when this same medium is supplemented with 10 −4 M of both iso and leu (competition scenario in Fig. 2b, see also S1 Fig and S2 Fig). We use a pair of engineered bacterial strains (yellow depicts Icells and blue stands for L -) that engage in mutualistic interactions by cross-feeding amino acids. b) Both strains are able to grow cooperatively in liquid cultures lacking both amino acids, but monocultures exhibit no growth in this conditions (Obl. Mutualism). When amino acids are supplemented at 10 −4 M (Competition), monocultures grow to comparable levels while the Lstrain overcomes its partner in cocultures. c) Time series (dots stand for average, shades for standard deviation) showing the coupled growth of the (obligate) hypercycle in liquid medium lacking both iso and leu. d) Invasion speed of the mutualistic strains according to a minimal reaction-diffusion model. The gray area indicates the domain where the cooperative interaction favours hypercyclic growth over Malthusian competition (parameter values extracted from observed growth for Iand Lin liquid cultures, see Supp. Info.). e) Front shape for the two strains for different self-reproduction rates (which models the effect of supplemented amino acids in the medium). The top panel shows the obligate (µ i = 0) hypercycle case: the coupled populations propagate as two travelling waves that approximately share the location of their fronts' edges. In the medium panel (µ i = µ Ci /2), the two species display interactions at the critical intersection that separate mutualism from competition: both strains travel at similar speeds, but the front edge of Iremains slightly behind one of Ldue to its smaller growth rate in the presence of amino acids. In the lower panel (µ i = µ Ci ), the faster replicator Lwins the competition by conquering the available space long before I -, which is progressively let behind until it is excluded from the population range expansion process. f) Observed front speed for cocultures spreading on agar surfaces, exhibiting a slowing-down in facultative mutualism scenarios that are not captured by the RD model.

5/24
Malthusian growth can break-down the hypercycle during range expansions How the synthetic hypercycle would spread in different environments is a focal question we address here. We used different concentrations of supplemented amino acids in order to modify the interactions displayed by our synthetic pair Iand L -(as done in Refs. [43,44] for mutualistic yeast strains). Regarding the hypercycle minimal model, let us consider amino acid supplementation as a way to introduce Malthusian growth rates for the species in the system (see the observed Malthusian growth rates in S1 Fig for the competition scenario). Furthermore, to be able to explore growth on solid culture conditions, the simple model [Eq. set (1)], must be redefined to incorporate space. To this aim, let us consider the following Reaction-Diffusion (RD) approach [40,63] as a minimal model describing the spatiotemporal dynamics of the synthetic mutualistic replicators: where I and L stand for the population density of the Iand Lstrain respectively, t and r are the time and spatial coordinates (see Methods), D is the diffusion coefficient, µ i is the Malthusian growth rate of species i, α ij (≥ 0) is the growth rate of species i assisted by its mutualistic partner j, and k is the carrying capacity of the system. The above set of equations generalised the two-member hypercycle model by including, on the one hand, the spatial context (through the diffusion terms D∂ 2 /∂r 2 ) and, on the other, by considering both Malthusian (µ k ≥ 0) and mutualistic (α ij ≥ 0) growth terms. By considering the absence of either species in the set (2), we recover the one-species Fisher RD model [46,47] that leads to the well-known expression for the invasion speed: Moreover, the Fisher speed establishes the asymptotic invasion speed for our two-species system in(2) as µ i >> α ij (for i = I, L and i = j = I, L). In the case of two purely competing species (µ i > 0, and α ij = 0) we should expect the front to propagate at the speed of the faster competitor because this species will be more efficient at conquering the available space at the edge of the population front. In contrast, for the case of two purely mutualistic species (i.e., a pure hypercycle with µ i = 0, and α ij > 0), we derived the analytical solution for the invasion speed (see Methods): Our minimal model (2) thereby predicts two different invasion modes for our pair of mutualistic strains Iand L -. Indeed, in the competition scenario, the invasion speed (3) is governed by the growth rate at low population densities(which gives rise to a pulled front [REF]). In contrast, the carrying capacity k appearing in Eq. (4) is a hallmark of an invasion front governed by the growth dynamics at high population densities. This gives rise to a pushed front [48] : individuals at the edge of the front are pushed from the inside bulk where individuals reproduce at higher rates. Moreover, note that the invasion speed (4) is the same for the two mutualists Iand L -, consistent with their need for a mutualistic partner in order to grow and spread. Figure 2d shows how the transition between the two invasion modes takes place, according to the RD model. In the absence of Malthusian replication (µ i = 0), both strains spread at the same speed. As both µ I and µ L are increased towards their observed value (see S1 Fig) in the competition scenario, the front speed increases due to the corresponding enhancement in growth rates. However, once µ i induces stronger effects on the front than α ij , competition becomes important and the coupled advance of the two strains is replaced by two differentiated front speeds. At this point, further increasing the Malthusian growth rates µ i benefits the faster species (in this case, the Lstrain), while the second one is slowed down in a relatively abrupt way. This eventually leads the Istrain to be excluded from the front (which propagates at the Fisher's speed c L F as Malthusian growth rates approach the observed values in competition). The population density profiles in Fig 2e illustrate the change in the invasion front shape as the scenario transits from obligate mutualism to competition.
Self-reproduction rates of the auxotrophic strains can be experimentally tuned by supplementing the medium with different doses of the two amino acids (an analogous method was used in Ref [43] for yeast mutualists). Fig. 2f shows that the observed speeds for cocultures spreading on agar (See Methods) are not consistent with those observed in 2d, within the considered range of amino acid concentrations. The particularly low values of the front speed observed in the experimental transition, from the obligate mutualism to the competition scenario, revealed one of the limitations of the minimal RD model. According to the RD model, the minimum front speed predicted for the synthetic hypercycle should be the invasion speed of the obligate mutualists. In other words, even if one of the strains is slowed down because of competition, the edge of the front will keep travelling at the speed of the fastest strain (which should exceed the speed of the obligate hypercycle in order to overcome its partner species at the edge of the front). Thus, the decrease of the observed front speed as supplemented amino acids are increased indicates that other, more complex phenomena are driving the dynamics of the synthetic hypercycle. In particular, the physical embodiment of bacterial cells (not taken into account by the RD model) may affect their access to the extracellular amino acids, thus influencing the invasion speed.

Slowdown of hypercycle front speed under local resource depletion in moderately rich environments
In order to further study the role of spatial structure in range expansions from synthetic hypercycles and how it is influenced by the strength of mutualistic interactions, we seed the cross-feeding system on M63 plates with 1.2% agar and different concentrations of auxotrophic amino acids distributed as is depicted in figure 3a. Each point in the experimental setting represents a different value of α ij [43,44]. When no amino acids are supplemented into the medium, cells are only able to growth if mutualistic partners remain close enough, the population engages in an obligate mutualism, which leads to a self-organized distribution with a characteristic high intermixing of the two strains. The opposite scenario (iso and leu at 10 −4 M) in Fig. 3a reveals a remarkably different spatial structure. When the driving interaction is competition for space and resources, the invasion dynamics is governed by genetic drift [52], which leads to demixing of the population into wide (single-strain) patches. This experimental condition correlates with µ M i in the RD minimal model (4), the competition scenario. In these conditions, we would expect the slower replicator (lower µ i ) to go extinct at the edge of the advancing front (see the competition case in Fig 2e). However, in our experimental scenario (Fig.  3a, top-right panel), both strains are present at the edge of the front despite exhibiting significantly different growth rates µ I < µ L (see Supp Info). This result is consistent with the expected effects of genetic drift in population range expansions [53]. In between of the above two modes of invasion, we found the environmental conditions that allow a facultative mutualistic behaviour. Single-strain patches are wider than those observed in the absence of supplemented iso and leu, although genetic diversity is still preserved (the characteristic width of patches is preserved) as the front propagates, Fig. 3c. In other words, in the facultative scenario, the concentration of amino acids added to the media permit the strains to grow into wider patches (compared to those of obligate mutualists), but both strains still benefit from the cross-feeding. However, in the competition scenario, the high concentration of amino acids permits the two strains to spread at comparable speeds regardless of the presence of the mutualistic partner (see Supp Info). Once mutualism is suppressed, genetic drift becomes the governing mechanism at the edge of the front, leading to progressively wider (3c) single-strain patches as the front advances.
The scenarios in Fig 3a reveal a qualitatively identical interplay between mutualism and genetic drift in range expansions of yeast populations Ref. [43], which suggests that such feature is universal and independent of the biological organisms exhibiting the 8/24 mutualism. However, boundary domains between the yeast strains in Ref. [43] are significantly different than the ones in Fig 3a. The differences on the boundary domains could be explained by the different cell shape [50]. In order to develop an agent-based model for range expansions of the bacterial hypercycle that consider cell shape, we used the GRO package [51]. Our agent-based model takes into account the cross-feeding interactions between the two strains and captures the experimentally observed scenarios at the edge of the front as a function of the supplemented amino acid concentrations (see Fig. 3b). Moreover, the model also captures different boundary domains associated to cell shape (as we show in S3 Fig).
GRO simulations allow us to study whether local depletion of supplemented amino acids could effectively slow-down the invasion process similarly as observed in experimental conditions. Nutrients and amino acids are mainly consumed by cells at the edge of the front 3d, their depletion leaves a population of stagnant cells that effectively constitutes a fossil record of the invasion process [52]. In the obligate mutualism case, single-strain patches keep a characteristic width determined by the distance at which cells can sustain the cross-feeding mutualism (cells near the front can temporarily become stagnant when their location prevents an effective cross-feeding). This process shapes the spatial distribution of the population, leading to a relatively high fraction of active cells at the edge of the front (3d and e). However, in the case of facultative mutualism, the dynamics can be marked by episodes of opportunistic growth that exploits the available amino acids in the environment. During these periods, the dynamics are locally governed by genetic drift (single-strain sectors become wider). However, once the supplemented amino acids are locally depleted, a significant number of cells (remote to the boundary domains where cross-feeding is still effective) can become stagnant (arrow in Fig 3d). Figure 3e shows how the ratio of active cells is correlated with the invasion speed, suggesting that the dynamics in facultative mutualism scenarios can slow-down the invasion speed of the synthetic hypercycle.

Environmental deterioration can determine the survival of parasites during range expansions
Several processes (such as mutations or the arrival of foreign, invader species) may give rise to new organisms exploiting hypercycle feedbacks in a given ecosystem. The introduction of a new replicator organism that makes use of the limited resources in the medium will restrict the growth of the hypercycle, even more, if this new organism is a parasite (hereafter P cells), that takes advantage of the cross-feeding (Fig. 4e).
In order to experimentally study the ecological implications of such parasites, we used the synthetic parasitic strain P (see Methods) that exploits one of the cross-feeding amino acids (namely, iso). The coculture of those three organisms in well-mixed conditions, for both the obligate mutualism and the competition scenarios, give as result a restricted growth of Ior Lstrains, Fig. 4a shows (compare to Fig 2b). Moreover, for the competition scenario in Fig 4a, the P strain exhibits a relatively high Malthusian growth rate (see Supp Info) that leads it to overcome the growth of the hypercycle pair.
To test whether spatial structure can limit the parasitic exploitation of hypercycles, we coculture combinations of the three strains (I -, Land P ) on M63-agar plates. In the absence of supplemented amino acids, when Ior Lcells are lacking, no growth was observed. This means that P cells can be considered a hypercycle parasite, because they are unable to close an effective cross-feeding loop (see Fig. 2a) with either Ior Lcells. When the three strains are present (Fig. 4b), despite an initial success of the parasite at colonising available space (see fig. 4c, red line), the parasitic strain is progressively left behind as the range expansion takes place. This is because, in the spatial scenario, cell location determines a preferential access to the cross-feeding . Environmental conditions determine the fate of parasites during range expansions. a) Obligate mutualism scenario (absence of supplemented amino acid) leads the strain P to act as a parasite in well-mixed conditions, while competition is observed at 10 −4 M supplementation of iso and leu. b) Spatial structure leads the mutualists to conquer the edge of the population front, defeating the parasite P. Yellow arrows indicate regions where the parasite has been excluded from the population front (red arrow indicates one of the few regions in which the parasite still surf at the population wave). Note that the front curvature is enhanced at regions governed by the mutualists, a hallmark of an enhancement of the front speed at these regions. The grey rectangle indicates the magnified area on the right. c) Frequency of the P strain at the edge of the front for two different scenarios (0 and 100µM extracellular ampicillin). d) The P strain offers cross-protection to the mutualists when threatened by antibiotics, leading to the survival of the P strain at the edge of the front. e) Scheme of the complex mutualistic interaction (which involves cross-feeding and cross-protection) between the three species in the presence of antibiotics. Each species lacks a different ability needed to survive in the system, but the ensemble may be able to survive if able to develop the corresponding division of labour. f) Three-species spatial structure in a simulated heterogeneous environment with non-isotropic antibiotic concentration at t = 0. While the P strain is conserved in the areas where cross-protection is essential for the mutualistic ensemble, P cells are excluded from the front in areas where the antibiotic concentration does not reach the growth inhibition threshold. metabolites [8,57]. Thereby, the presence of a P patch increases the distance between Iand Land leads to restricted growth. This gives a significant advantage to mutualistic Iand Lneighbouring patches that engage in an efficient cross-feeding. Hence, spatial structure benefits the hypercycle species, eventually leading the hypercycle ensemble to overcome the parasite at the edge of the front (Fig 4b).
The ecological role of a species in a given community can be strongly dependent on its environment and transitions can occur between mutualism and parasitism as external conditions change [2,[54][55][56]58]. In our three-member microbial consortium, composed by I -Land P, we studied whether environmental deterioration can make this community to develop a more complex mutualistic network. In order to do this, the three-member microbial consortium was seeded on m63-agar plates containing a lethal concentration of ampicillin, for which P cells are resistant. The P cells are able to degrade extracellular Ampicillin (by secreting beta-lactamase). Now, two different mutualistic motives are present in this scheme (Fig 4e): (amino acids) cross-feeding and (antibiotic) cross-protection. Remarkably, the hypercycle trio was able to solve the complex environmental problem and develop the range expansion process on the corresponding agar layers. Figure 4d shows the observed spatial structure displayed by this new mutualistic ensemble in order to collectively invade the available space. In contrast to the previous parasitic case, the fraction of the P strain is approximately constant as the population front advances (see Fig 4c).
The definition of the three-member consortium as an agent-based model allows us to make some predictions on how the system would spread within heterogeneous environments and captures the main spatial dynamics features of the system (see Supp Info). Simulation in a heterogeneous environment, that presents an asymmetric spatial antibiotic distribution, allows us to see how the P strain remains present at the edge of the front in the top region of the colony, which is precisely where the population is exposed to higher doses of antibiotic. In contrast, in the lower region where the antibiotic dose is much lower, the P strain is excluded from the edge of the front (consistently with our previous results), (Fig. 4f).
This is an interesting result particularly within the context of bioengineering soils [24,25] by the rewiring of the ecological interactions within the biological soil crust (BSC). Here the vertical structure defines a heterogeneous set of conditions where different species and physicochemical spatial gradients are present. Both in the BSC and around the plant root system a complex microbiome exists. Soil engineering under a systems perspective is a promising domain to harness and restore different functionalities [59]. This approach could be complemented by designed microbiomes exploiting mutualistic ties following some of the basic findings reported here. Since different soil conditions might sustain different qualitative functional traits, the previous synthetic three-species ecosystem can inspire novel forms of improving soil communities and plant efficiency.

Discussion
Most experimental and theoretical studies concerning the dynamics of microbial populations are grounded in competition. However, cooperation is a crucial component of ecological dynamics on all scales, and much needed to truly understand the behaviour of a wide range of systems, from populations growing on biofilms to the gut microbiome. Moreover, it has been suggested that synthetic cooperation can help to design ecological circuits capable of preventing endangered ecosystems from collapsing [24,25].
Previous studies have analysed a family of models involving closed cooperative loops of cooperators. These systems are known as hypercycles, and because of their second-order kinetics, they are capable of hyperbolic growth. This kind of dynamics 11/24 allows the hypercycle to overcome the simple Malthusian replicators. Theoretical work shows that hypercycles can prevent their decay due to the presence of parasites by exploiting the constraints imposed by a spatially extended system. However, these models require some special properties concerning the nonlinear dynamics of hypercyclic sets, which are not feasible in realistic conditions.Instead, we have analysed the problem of hypercycle persistence and response to parasites by means of experimental setups where populations of engineered cooperators spread on a two-dimensional medium.
Our study reveals that, as predicted by theoretical models involving both linear (Malthusian) growth and hypercyclic cooperation, spatial dynamics makes a big difference when space is introduced as propagating fronts. This is favoured by both the microscopic impact of bacterial shapes (leading to characteristic fractal structures) and by the local correlations required to sustain cooperation, which favour a maximisation of contact domains between the two cell populations. Hypercyclic growth has been characterised using diverse sets of measures and the front speed mathematically derived from a diffusion model.
The experiments and models confirm the picture of spatial hypercycles as dynamical systems where the mutualistic tie forces the formations of complex structures that guarantee the propagation of the cooperative consortium. We have also studied the tradeoffs associated with Malthusian growth and the conditions pervading the breakdown of hypercyclic cooperation thus showing the presence of two phases: one associated with competitive interactions and a second phase associated with scarce resources promoting the cooperative feedback.
The second set of experiments and models are related to the impact of parasitic strains on the stability of the hypercycle. Parasites are easy to evolve under standard conditions of growth in cell cultures: when synthetic constructs have been added to microbial strains, the loss of one construct has an immediate impact on the metabolic burden, thus allowing the cell to replicate faster. However, the presence of a given dependency can create nontrivial interactions. We designed synthetic parasitic strains capable of exploiting a given amino acid while not completing the cooperation cycle. Such parasite (which has a small component of Malthusian growth) has been shown to overcome and kill the hypercycle under liquid conditions but becomes a much less harmful component under spatial constraints. It was recently shown that resource availability can modulate the interactions between microbial cross-feeding mutualists [44]. Our work is, as far as we know, the first experimental design of a synthetic ecological network showing how different contexts allow cooperation, competition or parasitism to succeed or even transition from one to the other in a spatially extended context. Further work should explore how these results translate into more realistic contexts, from the gut microbiome to soil ecosystems.

Theoretical invasion speed of a 2-species hypercycle
Our theoretical RD model for the two-species hypercycle considers that the dynamics of the system is governed by diffusion and population growth as: For convenience, we rewrite this set of equations in terms of dimensionless variables I * = I/k, L * = L/k, t * = α IL kt and r * = (α IL k/D) 1/2 r, and dimensionless parameters α * = α LI k/α IL . Then, the new set reads: Let us assume that there exist travelling wave-shaped solutions of the previous equations of the form: with s > 0, b > 0, a > 0, and z = r − ct (where c is the speed of the travelling wave, i.e. the front speed of the hypercyclic population). Using (7) can be rewritten as: Developing the derivatives U I and U I , Eq. (10) reads: where η = (1 + ae bz ). Neglecting the trivial solution (ε I = 0) for Eq. (12), and reorganising terms according to powers of e bz , we obtain the characteristic equation for the front speed c: Solutions for the travelling wave have to be valid ∀z, and thus each line in Eq. (13) gives an independent expression that must necessarily vanish. Analysing the terms in the last line in Eq. (13) leads to the necessary condition s < 2. This leads to s = 1 because we only consider solutions with s > 0. Then, considering s = 1, we develop the conditions given by the different powers of e bz in Eq. (13), which leads to: and b = c.
Combining Eqs. (14)- (16) leads to: With an analogous procedure to the one performed above for Eq. (10), analysis of Eq. (11) leads to: Combining Eqs. (17) and (18) we obtain the expressions for the species abundances in the travelling front: Replacing terms from Eq. (19) into Eq. (18), we obtain the analytical solution for the front speed in dimensionless variables: Finally, recovering dimension variables, the speed of the front reads: The agent based model Our approach to the study of hypercycles reveals the importance of considering cells as embodied entities, both as interacting elements on a microscopic scale and as spatially extended populations. Moreover, cells need to incorporate the molecular circuits associated to the specific regulatory mechanisms along with chemical reactions, spatial diffusion and molecular signalling. To this goal, we used the specification language gro [51] as the platform for individual-based simulation of growing populations. Our model integrates the main physical features of bacterial shape and growth [51], as well as the cross-feeding and cross-protection interaction between I -Land P strains. We used a very simple approach that considers a few step (Heavyside) functions to emulate cell behaviour. A list of the considered cell behaviour features follows: 1. Sensing: at each time step, each cell senses the extracellular concentration of three kinds of molecules: amino acids (Icells sense iso, while Land P cells sense leu), food (this category embraces any other nutrients that cells may need to grow), and antibiotic (i.e., ampicillin).
2. Growth: cells grow (increase their cell volume) and divide at the realistic speed proposed in Ref. [51], provided that: (a) food concentration exceeds a given threshold value g f .
(b) the corresponding amino acid (according to cell strain) exceeds a given threshold value g am . (c) antibiotic concentration is below a given inhibitory threshold g at .
Accordingly, cell growth is arrested whenever any of the above conditions are violated.
3. Cells absorb extracellular food and release amino acid (or β-lactamase) at constant rates, provided that extracellular food exceeds g f . Specifically, Icells release leu, Lcells release iso, and P cells release the betalactamase enzime (that degrades the antibiotic) to the extracellular medium. Provided that growth conditions are satisfied, cells will also absorb the amino acid they need.
The corresponding logical loop experienced by a given Lcell at each time step is illustrated in Fig 5. Iand P cell dynamics follow analogous logical schemes.
Furthermore, in order to consider a fitter parasitic strain that evades the cost of the mutualism in antibiotic-free scenarios, we consider the growth rate of P cells to be higher (by a 10% difference) than that of Iand Lcells. As shown above, the hypercycle was able to escape the parasite despite such faster growth rate.
Admittedly, actual cell dynamics is far more complex than this Heavyside representation. However, our goal for the agent-based model was to use a minimal set of assumptions, in order to provide an easy understanding of the key features governing the system dynamics. Remarkably, the Heavyside-based cell behaviour is enough to capture the essential dynamics, as discussed in the Results section. The source code and additional details on specific values for metabolic rates and concentration threshold values can be found in the Supp. Info.

Bacterial strains
Both the Iand the Lstrains are from E. coli strain DH1 (National BioResource Project, National Institute of Genetics, Shizuoka, Japan) and were genetically modified to cross-feed as described in [6]. The I -(L -) strain carries the dsred.T3 (gfpuv5) gene that provides the corresponding fluorescence labelling.
Cloning for the P strain was carried out using the Biobrick assembly method and the parts: B0014, J23100, B0032 and E0020, from the Spring 2010 iGEM distribution 15/24 assembled into a low copy number plasmid pSB4A5. A complete description of the construction protocols can be found at [61,62].

Culture conditions
All regular cultures and amplifications were done at 37 ℃ in well-mixed media Lysogeny Broth (LB). Bacterial strains were cryopreserved in LB-glycerol 20% (v/v) at -80 ℃. 12h before each experiment, 5 colonies of each type of cell were selected and grown separately in LB plus 10 − 4 M of auxotrophic amino acid at 37 ℃. After 12h, a 100-fold dilution with fresh amino acid supplemented LB plus antibiotics (for auxotrophic cells) or 500-fold (for the parasite), with fresh LB plus cabenicillin (100 µg/ml), were done and grown until OD660∼0.4.

Range expansions on agar surfaces
Fresh cultures, with a OD660∼0.4, were washed twice by centrifugation and resuspension with M63 medium. The optical density of each culture was settled to OD660 nm = 0.15. Equitable volumes were mixed to generate cocultures with 2 or 3 different cells. Finally, 0.4 µL of the mono or coculture volume was seeded on the center of a m63 (1.2% agar). Cells were grown 4 days at 37 ℃ and humidity 90%. colonies were observed using a Leica TCS SP5 AOBS (inverted) confocal microscope.

S1 Fig
Malthusian and hypercycle growth rates for the synthetic strains Time series for the fluorescence of the Istrain, when cultured in M63 medium supplemented with 100 µM of both iso and leu. Coloured dots stand for the average values across 9 replicates (three technical replicates from each of three biological replicates), shaded area indicates standard deviation. The Malthusian growth rate µ I was obtained by linear regression (black solid line) to the data during the exponential growth regime (region delimited by the vertical dashed lines), as described in S1 Text. b) Malthusian growth rate for the Lstrain (growth conditions as in a)). c) Malthusian growth rate for the P strain (growth conditions as in a)). Hyperbolic growth rates α IL and α LI were obtained from the observed growth at low population densities (region between dashed lines), as described in S1 Text. The time series correspond to the growth of both Iand Lstrains in coculture, in M63 medium with no supplemented amino acids.

Growth rates in well-mixed conditions
In order to measure the Malthusian growth rates associated with competition scenarios, we cultured each of the three strains in M63 medium supplemented with 100 µM of both iso and leu amino acids. Fluorescence measures showed consistency with the expected exponential growth regime, that precedes growth saturation when the population reaches its carrying capacity, as shown in S1 Fig. Malthusian growth rates for the three species were obtained through linear regression of the observed growth data, according to the Malthusian growth model: where F stands for the fluorescence value, µ j is the Malthusian growth rate of species j, t stands for the time and β is a constant value for t = 0. Thus we obtained the values µ I = (9.1 ± 0.1) × 10 −2 hr −1 , µ L = (2.18 ± 0.02) × 10 −1 hr −1 and µ P = (3.75 ± 0.02) × 10 −1 hr −1 .
In order to obtain an estimate for the hyperbolic growth rates, we considered the well-mixed version of Eq. (2) in the Main text. This correspond to the following set of equations, which do not account for the diffusion process: where we have also assumed that we deal with the obligate mutualism scenario (µ I = µ L = 0).
Considering an approach for the growth at low population densities, it is easy to obtain a solution for I(t) and L(t). Thus, let us neglect the carrying capacity effects [last term on the right hand side of the set of Eqs. (23)], which leads us to: The above set (24) permits to write I(t) in terms of L(t) as: with δ = I 0 − α IL L 0 /α LI (where I 0 and L 0 stand for the population density numbers at the initial instant t 0 ). Replacing Eq. (25) into (24), and solving the corresponding differential equation we get the following solution for L(t): where with a = 1/δ and b = −aα IL /α 2 LI . The above Eqs. (25) and (27) constitute a set that we can adjust to the observed growth in well-mixed conditions, with α IL and α LI as the only adjustable parameters. By applying a least squares algorithm to the the fluorescence time series for the growth of the obligate mutualists, we obtained α IL = (4.4 ± 0.1) × 10 −2 hr −1 and α IL = (6.2 ± 0.1) × 10 −2 hr −1 , as S1 Fig shows.

S2 Fig
Cell concentration scales linearly to fluorescence for the three species a) Cell concentration in liquid cultures of the Istrain according to their fluorescence. The value of a indicates the slope (in ml −1 ) obtained by linear regression of the data points. b) In agreement with cell concentration, optical density also scales linearly to fluorescence for theIstrain. c) and d) show the same analysis as in a) and , but for the L -(while e) and f) correspond to analogous results for the P strain).

S1 Table
Relevant parameters in the agent-based model The table shows the main parameters of the agent-based model, as well as the main processes they affect. Unless stated otherwise in the text, the parameter values used in simulations correspond to those in the source code (S2 Text).

S2 Text
Agent based simulations source code Source code used to run our simulations in the GRO package [51].

S3 Fig
Cell shape influences mesoscopic boundary domains a) Fractal dimension for the boundaries between Iand Lpatches in the obligate mutualism scenario. Bars indicate average values, while vertical lines indicate standard deviation from three different simulations. b) A snapshot showing the patches of the Istrain (in white), when de division size parameter is set to 2.0, for a colony with approximately 1.6 × 10 4 individuals. c) A snapshot showing the patches of the Istrain (in white), when de division size parameter is set to 3.5, for a colony with approximately 1.6 × 10 4 individuals.

S4 Fig
Agent-based simulations capture the spatial dynamics of hypercycle range expansions a) Agent-based simulations show analogous scenarios to those observed in Fig 3a. Values on the vertical and horizontal axis indicate the parameter values for the initial extracellular concentration of amino acids (I0 and L0, respectively, see S1 Table). b) Patch width in simulated range expansions, for a different initial extracellular concentration of amino acids (initial nutrient concentration F 0 = 90). c) A biological replicate for each of the cases presented in Fig. 3c in the Main Text.

S5 Fig
Fraction of P strain in range expansions a) In silico, fraction of territory colonized by P cells in three-species population range expansions. Three different scenarios are shown: no ampicillin (Ampi0 = 0.0, see S1 Table), moderate ampicillin concentration (Ampi0 = 2.0), and high ampicillin concentration (Ampi0 = 4.0). b) Biological replicate for the two scenarios in Fig. 4c.

Front speed for one-species hypercycles
In the Main Text, we have derived an analytical solution for the front speed of two-species hypercycles. For completion, we here present the theoretical speed that would correspond to a one-species hypercycle. Analogous theoretical front speeds for one species mutualistic populations are also described in Refs. [40,63].
Let's consider an expanding one-species hypercyclic population (u), modelled under the following reaction-diffusion approach: We are interested on the invasion speed to the hypercycle. As a first approach, we will consider that the front is planar, and then the one-dimensional speed is a good approximation of the speed of the front in two dimensions. We rescale t and x in order to work with dimensionless variables: t * = rt, and x * = x( r D ) 1/2 . Then, Eq. (28) becomes: We look for front propagation solutions, so let us assume that there exist solutions to Eq. (29) with the propagating wave form: u(x * , t * ) = U (z) = 1 (1 + ae bz ) s , with b, s > 0 and z = x * − ct * . Using, u x = U z = U , and u t = −cU z = cU (where the subscripts denote the corresponding partial derivatives), we obtain the expressions for the following partial derivatives: u xx = s(s + 1)a 2 b 2 e 2bz (1 + ae bz ) −s−2 − sab 2 e bz (1 + ae bz ) −s−1 (32) u t = csabe bz (1 + ae bz ) −s−1

19/24
We rewrite Eq. (29) as: Rewriting Eq. (34), and reorganising terms as powers of e bz we obtain: e 2bz (s(s + 1)a 2 b 2 − sa 2 b 2 − csa 2 b) +e bz (−sab 2 − csab) The above Eq. (35) has to be equal to zero ∀z. Then ,the coefficients of e 0 , e b z, and e 2 bz must all be identical to zero. Taking into account that we look for travelling waves solution of the form (34) (with s > 0), it is easy to show that the only value that s can take is s = 1. Using this into Eq. (34) leads to: which leads to the dimensionless front speed: Recovering dimension variables from the front speed c we obtain: which leads to the dimensionless front speed: v = rD/2.