Team:ETH Zurich/Modeling/Reaction-diffusion

"What I cannot create I do not understand."
- Richard Feynmann

Reaction-diffusion Models

Introduction

While single-cell models are useful for correctly implementing and debugging chemical reaction models, they are not sufficient to fully understand the real-life functionality of our system. Since an essential part of our system is increasing the perceived concentrations of lactate and AHL through co-localization, it is necessary to model the concentrations the chemical species though a reaction-diffusion system.

Figure 1: "Doughnut" model for our system.

"Doughnut" model in COMSOL

Geometry: Abstraction of a population of discrete E. coli to a connected uniform space

The well in which our reaction takes place is a box of length 200μm. In the middle of the well, we have a circle representing our target cell. Surrounding this is a ring of appropriate radius (see below) representing the union of all E. coli cells bound to the target cell. Finally, the remaining space represents both the medium and free-floating E. coli cells.

Logistic growth of E. coli

If we let \(n_\text{bulk}:[0,tsim]\longrightarrow \mathbb N\) be the function giving the number of E. coli in the bulk at time \(t\), with an initial population of size \(n_0\), a carrying capacity of \(K\), and a growth rate of \(R\), then $$n_\text{bulk}(t) = \frac{n_0 K e^{Rt}}{L + n_0(e^{Rt}-1)}$$ so with an E. coli doubling time of \(t_2\) (which we assume is 30 minutes), we can solve for \(R\) $$R = \frac{\log2 + \log(K+n_0) - \log(K-2n_0)}{t_2}$$

Concentration correction for differing volumes

Since we are assuming a fixed number of bound E. coli to the target cell and since diffusion occurs almost instantly in our well, the concentrations in the doughnut will be accurate if we set its radius such that the area is the correct value. \begin{align*} A_\text{doughnut} &= \pi(r_\text{target} + r_\text{doughnut})^2 - \pi r_\text{target}^2 = n_\text{bound}\pi r_\textit{E. coli}^2 = n_\text{bound}A_\text{bound}\\ \Rightarrow r_\text{doughnut} &= \sqrt{n_\text{bound} r_\textit{E. coli}^2 + r_\text{target}^2} - r_\text{target} \end{align*} Unfortunately, the same principle does not apply for the bulk due to the logistic growth of the E. coli, so a more nuanced approach is necessary. Let \(n\) represent the number of chemical species and \(t_\text{sim}\) be our total simulation time. Let \(\mathbf X:[0,tsim]\longrightarrow \mathbb R^n\) be a function representing the molar concentrations of our chemical species over the simulation period. Then our system of non-linear ordinary differential equations (ODEs) can be represented by the following equation $$\left.\frac{d\mathbf X}{dt}\right|_{(\mathbf X(t),t)} = f(\mathbf X(t),t)$$ The units of \(X_i\) are \(\frac{\text{mol}}{L_\textit{E. coli}}\), we can define a new function \(\mathbf Y(t) := \nu(t)\mathbf X(t)\) representing the concentrations of the species within our simulated bulk, where $$\nu(t) := \frac{n_\text{bulk}A_\textit{E. coli}}{A_\text{bulk}}$$ Our original ODE system after this change of variables is then $$\left.\frac{d\frac{\mathbf Y}{\nu}}{dt}\right|_{(\mathbf Y(t),t)} = f\left(\frac{\mathbf Y(t)}{\nu(t)},t\right)$$ If we solve for \(\frac{d\mathbf Y}{dt}\), we get $$\left.\frac{d\mathbf Y}{dt}\right|_{(\mathbf Y(t),t)} = \nu(t)\left(f\left(\frac{\mathbf Y(t)}{\nu(t)},t\right) + \mathbf Y\frac{d}{dt}(\nu(t))^{-1}\right)$$

Diffusion and transport of chemical species

Under alkaline conditions, E. coli actively import lactate via a proton-motive symporter. Thus, a cross-membrane transport reaction had to be implemented. Since this is not possible directly in COMSOL, we had to model lactate in two states. Suppose our reference is the subspace of the interiors of the E. coli. We then defined the two states \(Lac_\text{int}\) and \(Lac_\text{ext}\), denoting intracellular and extracellular lactate, respectively. \(Lac_\text{ext}\) is produced by the target cell and can diffuse freely though the medium and all membranes. \(Lac_\text{int}\) is in equilibrium with \(Lac_\text{ext}\) with rate constants set to maintain a 20-fold difference of lactate concentration between interior and exterior. $$ Lac_\text{ext} \mathop{\mathop{\xrightarrow{\hspace{4em}}}^{\xleftarrow{\hspace{4em}}}}_{k_\text{int}}^{k_{\mathrm{ext}}} Lac_\text{int} \qquad \frac{k_\text{int}}{k_\text{ext}}\approx 20 $$ In addition, only the \(Lac_\text{int}\) state can react with the other chemical species in the E. coli.

AHL is able to freely diffusion in the medium and across membranes. All other chemical species are only able to diffuse intracellularly. The effective diffusion coefficients of AHL and Lactate through the E. coli membrane \(D_e\) were approximated by the method proposed by [Stewart 2003] as a fraction of their respective diffusion coefficients in water \(D_{aq}\) by the relation $$\frac{D_e}{D_{aq}}\approx 0.25$$

Four cases

To test whether our system acts as an AND gate on our two inputs (higher lactate production and co-localization signals), we combinatorially tested our system in environments with high vs. low lactate production and E. coli co-localization vs. dispersion.

Results

Parameter search on a simplified reaction-diffusion model in MATLAB