A mechanistic model of metabolic symbioses in microbes recapitulates experimental data and identifies a continuum of symbiotic interactions

Abstract

Microbial symbioses based on nutrient exchange and interdependence are ubiquitous in nature and biotechnologically promising; however, an in-depth mathematical description of their exact underlying dynamics from first principles is still missing. Hence, in this paper a novel mechanistic mathematical model of such a relationship in a continuous chemostat culture is derived. In contrast to preceding works on the topic, only parameters which can be directly measured and understood from biological first principles are used, allowing for a higher degree of mechanistic understanding of the underlying processes compared to previous approaches. The predictive power of the model is validated by demonstrating that it accurately recapitulates both the temporal dynamics as well as the final state of a previously published cross-feeding experiment. The model is then used to examine the influence of the biological traits of the involved organisms on the position and stability of the equilibrium states of the system using bifurcation analyses. It is additionally demonstrated how manipulating the external metabolite concentrations of the system can shift the species interaction on a continuous spectrum ranging from mutualism over commensalism to parasitism. This further reinforces the idea of a continuous spectrum of symbiotic interactions as opposed to static and discrete categories. Finally, the practical implications of the results for the biotechnological application of such microbial consortia are discussed.

Appendices

Appendix A: Parametrisation of the model for the validation with experimental data

In order to validate the derived model using previously published experimental data, the model was parametrised using parameter values provided in the experimental paper by Zhang and Reed (2014) or, when such parameter values have not been provided by the authors, estimated from literature values.

Specifically, the initial glucose concentration of the medium was set to $R(0) = 2 \; {\mathrm{{g/l}}}$ following the description of the culture medium given in the paper, while both initial metabolite concentrations were set to $M_1(0) = M_2(0) = 0$, as neither of the two exchanged metabolites (leucine or lysine, respectively) was supplemented into the growth medium. The washout rate was set to $\omega = 0$, as the coculture was not diluted over the course of the 75-h long experiment. Initial population densities were set to $N_1(0) = N_2(0) = 5 \times 10^6$ cells per millilitre assuming a 50:50 mixture of the two strains at the initial OD of 0.01 described by the authors.

The intrinsic growth rates of the two auxotrophy mutants were set to the values provided by the authors, i.e. $r_1 = 0.461 \; {\mathrm{{h}}}^{-1}, r_2 = 0.465 \; {\mathrm{{h}}}^{-1}$. The efflux rates of the two strains have been measured to be 0.027 mmol/(gDW $\times$ h) leucine, and 0.02 mmol/(gDW $\times$ h) lysine, respectively. Using a literature value of $3 \times 10^{-13} \; \mathrm{{g}}$ dry weight (DW) per E.coli cell (Neidhardt et al. 1990), and molar masses of 131.17 g/mol for leucine and 146.19 g/mol for lysine, efflux rates of $\epsilon _1 = 1.1 \; \mathrm{{fg/(cell}} \times {\mathrm{{h}}})$ and $\epsilon _2 = 0.9 \; \mathrm{{fg/(cell}} \times {\mathrm{{h}}})$ were estimated. Metabolite requirements were measured by the authors as 0.35 mmol/(gDW) lysine and 0.473 mmol/(gDW) leucine, respectively. Using the same conversion factors as before, these amount to $\gamma _1 = 15 \; {\mathrm{{fg/cell}}}$ lysine and $\gamma _2 = 19 \; {\mathrm{{fg/cell}}}$ leucine.

The Monod constants of the two strains for the respective two amino acids has not been measured by the authors, so it was estimated from Kerner et al. (2012), who measured the Monod constants of E. coli for the amino acids tyrosine and trypsine, to be approximately $K_1 = K_2 = 5 \; \upmu {\mathrm{{g/l}}}$. For the Monod constant of E. coli for glucose, the literature value of $L_1 = L_2 = 0.1 \; {\mathrm{{g/l}}}$ from Senn et al. (1994) was used. As a side note, this striking difference of several orders of magnitude between $K_1, K_2$ and $L_1, L_2$ indicating a very high affinity of E. coli for the exchanged amino acids does not come as a surprise, as bacteria have evolved highly specialised permeases which are able to take up even extremely small concentrations of externally present aromatic amino acids, most likely in order to avoid the significantly more energetically expensive de novo synthesis of these amino acids which is switched off via feedback inhibition, whenever the respective amino acid is sufficiently present externally (Ames 1964). Finally, from Fig. 2b of the paper, the requirement of glucose per one bacterial cell was estimated to be $\alpha _1 = \alpha _2 = 5000 \; {\mathrm{{fg/cell}}}$ by dividing the amount of utilised glucose by the final population density of the culture, assuming that bacterial death is neglectable for the short duration of the experiment of 75 h.

Appendix B: Proofs of Lemma 1 and 2

Proof

(Lemma 1—Invariance of the nonnegative orthant) Consider the behaviour of the five differential equations if the respective state variable equals zero. One gets:

$$\begin{aligned} \begin{aligned} \dot{N}_1(t) \arrowvert _{N_1(t)=0}&= 0 \\ \dot{N}_2(t) \arrowvert _{N_2(t)=0}&= 0 \\ \dot{R}(t) \arrowvert _{R(t)=0}&= \omega R_{{\mathrm{{in}}}}> 0 \\ \dot{M}_1(t) \arrowvert _{M_1(t)=0}&= \omega M_{1, {\mathrm{{in}}}} + \epsilon _1 N_1(t) - 0> 0\\ \dot{M}_2(t) \arrowvert _{M_2(t)=0}&= \omega M_{2, {\mathrm{{in}}}} + \epsilon _2 N_2(t) - 0 > 0\\ \end{aligned} \end{aligned}$$

(10)

Accordingly, the system will never be able to leave ${\mathbb {R}}_{\ge 0}^5$. $\square$

Proof

(Lemma 2—Boundedness) By inspection of the third differential equation, one can easily see the upper bound $\forall t: R(t) < R_{{\mathrm{{in}}}}$.

Now introduce the quantity $C(t) {:}{=}\alpha _1 N_1(t) + \alpha _2 N_2(t) + R(t)$, describing the complete concentration of resource in the culture vessel in form of free resource or microbial individuals. With some arithmetics, one may deduce that

$$\begin{aligned} \dot{C}(t) = \omega \cdot (R_{{\mathrm{{in}}}} - C(t)), \end{aligned}$$

(11)

from which directly follows that $N_1(t), N_2(t)$ are bounded by the amount of resource influx below some values $N_1^{\max }, N_2^{\max } \in {\mathbb {R}}^+$.

From this, one may obtain the following upper bounds of $\dot{M}_1(t)$:

$$\begin{aligned} \forall t: \dot{M}_1(t) \le \omega M_{1, {\mathrm{{in}}}} + \epsilon _1 N_1^{\max } - \omega M_1(t) \end{aligned}$$

(12)

Thus, $M_1(t)$ is bounded and similarly, $M_2(t)$ is as well, which concludes the proof. $\square$

Appendix C: Jacobian of the system

The Jacobian of the system is given by

$$\begin{aligned} \mathbf{J } = \begin{bmatrix} \dfrac{M_2 R r_1}{(K_1+M_2)(L_1+R)} - \omega & 0 & \dfrac{L_1 M_2 N_1 r_1}{(K_1+M_2)(L_1+R)^2} & 0 & \dfrac{K_1 N_1 R r_1}{(K_1+M_2)^2 (L_1+R)} \\ 0 & \dfrac{M_1 R r_2}{(K_2+M_1)(L_2+R)} - \omega & \dfrac{L_2 M_1 N_2 r_2}{(K_2+M_1)(L_2+R)^2} & \dfrac{K_2 N_2 R r_2}{(K_2+M_1)^2 (L_2+R)} & 0 \\ - \dfrac{\alpha _1 M_2 R r_1}{(K_1 + M_2)(L_1+R)} & - \dfrac{\alpha _2 M_1 R r_2}{(K_2 + M_1)(L_2+R)} & - \dfrac{\alpha _1 L_1 M_2 N_1 r_1}{(K_1+M_2)(L_1+R)^2} - \dfrac{\alpha _2 L_2 M_1 N_2 r_2}{(K_2+M_1)(L_2+R)^2} - \omega & - \dfrac{\alpha _1 K_2 N_2 R r_2}{(K_2+M_1)^2 (L_2+R)} & - \dfrac{\alpha _2 K_1 N_1 R r_1}{(K_1+M_2)^2 (L_1+R)} \\ \epsilon _1 & - \dfrac{\gamma _2 M_1 R r_2}{(K_2 + M_1)(L_2+R)} & - \dfrac{\gamma _2 L_2 M_1 N_2 r_2}{(K_2+M_1)(L_2+R)^2} & - \dfrac{\gamma _2 K_2 N_2 R r_2}{(K_2+M_1)^2 (L_2+R)} - \omega & 0 \\ - \dfrac{\gamma _1 M_2 R r_1}{(K_1 + M_2)(L_1+R)} & \epsilon _2 & - \dfrac{\gamma _1 L_1 M_2 N_1 r_1}{(K_1+M_2)(L_1+R)^2} & 0 & - \dfrac{\gamma _1 K_1 N_1 R r_1}{(K_1+M_2)^2 (L_1+R)} - \omega \end{bmatrix}. \end{aligned}$$

(13)

As long as no metabolites are manually added to the culture, the Jacobian of the extinction state always reduces to

$$\begin{aligned} \mathbf{J } = \left[ \begin{matrix} - \omega & 0 & 0 & 0 & 0 \\ 0 & - \omega & 0 & 0 & 0 \\ 0 & 0 & - \omega & 0 & 0 \\ \epsilon _1 & 0 & 0 & - \omega & 0 \\ 0 & \epsilon _2 & 0 & 0 & - \omega \end{matrix} \right] . \end{aligned}$$

(14)

For the non-trivial case, where all system variables exceed zero, using Eq. (6) the Jacobian of the system simplifies to

$$\begin{aligned} \mathbf{J } = \begin{bmatrix} 0 & 0 & {\mathcal {L}}_1 & 0 & {\mathcal {K}}_1 \\ 0 & 0 & {\mathcal {L}}_2 & {\mathcal {K}}_2 & 0 \\ - \alpha _1 \omega & - \alpha _2 \omega & - \alpha _1 {\mathcal {L}}_1 - \alpha _2 {\mathcal {L}}_2 - \omega & - \alpha _1 {\mathcal {K}}_2& - \alpha _2 {\mathcal {K}}_1 \\ \epsilon _1 & - \gamma _2 \omega & - \gamma _2 {\mathcal {L}}_2 & - \gamma _2 {\mathcal {K}}_2 - \omega & 0 \\ - \gamma _1 \omega & \epsilon _2 & - \gamma _1 {\mathcal {L}}_1 & 0 & - \gamma _1 {\mathcal {K}}_1 - \omega \end{bmatrix}, \end{aligned}$$

(15)

where ${\mathcal {L}}_1 = \frac{L_1 N_1 \omega }{(L_1 + R) R}$, ${\mathcal {L}}_2 = \frac{L_2 N_2 \omega }{(L_2 + R) R}$, ${\mathcal {K}}_1 = \frac{K_1 N_1 \omega }{(K_1 + M_2) M_2}$, and ${\mathcal {K}}_2 = \frac{K_2 N_2 \omega }{(K_2 + M_1) M_1}$.