Modelling Bacteria-Inspired Dynamics with Networks of Interacting Chemicals

One approach to understanding how life-like properties emerge involves building synthetic cellular systems that mimic certain dynamical features of living cells such as bacteria. Here, we developed a model of a reaction network in a cellular system inspired by the ability of bacteria to form a biofilm in response to increasing cell density. Our aim was to determine the role of chemical feedback in the dynamics. The feedback was applied through the enzymatic rate dependence on pH, as pH is an important parameter that controls the rates of processes in cells. We found that a switch in pH can be used to drive base-catalyzed gelation or precipitation of a substance in the external solution. A critical density of cells was required for gelation that was essentially independent of the pH-driven feedback. However, the cell pH reached a higher maximum as a result of the appearance of pH oscillations with feedback. Thus, we conclude that while feedback may not play a vital role in some density-dependent behavior in cellular systems, it nevertheless can be exploited to activate internally regulated cell processes at low cell densities.


Introduction
Bacteria are arguably one of the most prolific forms of life on earth, capable of surviving the harshest of environmental conditions. The processes that allow bacteria to maintain homeostasis, adapt to changing environmental conditions, and switch between vastly different states in response to external stimuli are often driven by internal feedback mechanisms. Feedback, when a process is regulated by its output, plays a vital role in the functioning of living organisms [1].
Feedback is also utilized to drive a sharp switch from single cell to multi-cellular behavior in a phenomenon known as quorum sensing [2]. Bacteria communicate by release and amplification of signaling species resulting in population-level responses such as the formation of biofilm above a threshold number or density [3]. Biofilms consist of a sticky matrix of polysaccharides and other biomolecules that protect the bacteria against the external environment. The development of antimicrobial resistance or tolerance typically involves the formation of a biofilm driven by a quorum sensing mechanism.
The design of biologically-relevant reaction networks that can reproduce dynamical cell behaviors may provide an insight into the role of feedback in such phenomena [4,5]. Progress in systems and synthetic biology has resulted in oscillations and quorum sensing driven by autocatalytic networks in genetically modified organisms [6]. The systems approach to chemistry is also undergoing rapid development [7] but although feedback based on DNA [8], peptides [9] and other biomolecules [10,11] has been implemented, autocatalytic communication mechanisms have been mainly investigated in regard to inorganic redox processes [12,13]. There are relatively few robust examples of chemical feedback in which cell-to-cell communication can be explored.
Taking inspiration from bacteria, we developed a model to explore the role of feedback in a cellular system when the feedback is coupled to other process, such as gelation or precipitation of a substance in the external solution. The model consisted of an enzyme-catalyzed reaction network based on urease, a virulence factor produced by bacteria [14] and associated with biomineralization, the precipitation of inorganic salts by living organisms [15]. It has been demonstrated that the urease reaction can display feedback as a result of the enzyme rate dependence on pH [16,17]. The underlying chemistry, although biologically relevant, is not implicated in quorum sensing in bacteria. However, feedback through pH is of interest as pH plays an important role in the rates of many processes and oscillations in pH have been observed in living systems as well as implicated in the origin of life [18,19].
Here, we used the model to explore the influence of increasing density of cells on the dynamics and determined the importance of feedback on the resultant behavior. Our aim was to build a simple model that may shed light on the role of feedback in natural processes but also to aid in the design of synthetic cells that may be exploited in biotechnological applications. There have been great advances in the development of synthetic cells based on, for example, gene networks or enzymatic reactions in liposomes [20,21], but few examples exploiting feedback. Hence our model is based on an experimentally-realizable reaction network that might be implemented in a synthetic cell, rather than a realistic mimic of biological processes in bacteria.

Methods
The enzyme-catalyzed reaction network we considered here was inspired by the enzyme urease described in detail elsewhere [16,22], where S is the substrate and EH (active form of the enzyme) is the active form of the polyprotic enzyme which binds to the substrate to give product, N: The protonation equilibria constants of the enzyme in (1a) are K e1 and K e2 . At low and high pH, the dominant forms of the enzyme are EH 2 + and E − respectively, which are considered inactive.
Therefore, the rate of the enzyme-catalyzed reaction (1a + 1b) is given by a modified Michaelis-Menten expression, taking into account the fact that the activity depends on pH: where v max is the maximum enzyme rate and K M is the Michaelis constant. The constants for the enzyme-catalyzed reaction were v max = 0.0001 M/s; K e1 = 5 × 10 −6 M; K e2 = 2 × 10 −9 M and K M = 3 × 10 −3 M. The product, N, of the enzyme reaction is a weak base, and reacts with acid to form NH + : where the net rate is R2 and k 2 = 1 × 10 −3 s −1 , k 2r = 1 × 10 7 M −1 s −1 and pK a = 10. The addition of a large organic molecule, B, that can buffer changes in pH was also considered: with net rate R3 and k 3 = 1, k 3r = 1 × 10 7 and pKa = 7.
The key processes in the enzyme reaction network are illustrated in Figure 1a where the blue dashed box shows the enzyme protonation equilibria. The dependence of the enzyme rate on acid concentration gives rise to the typical bell-shaped curve expected for enzyme-catalyzed reactions, with a maximum rate for this reaction at pH = 7. Combining (2) and (3) with the water dissociation equilibrium: yields the product highlighted in the purple box and the buffer equilibrium in the red dashed box. Feedback arises through the pH: starting from low pH, the production of OH − leads to an increase in pH and hence an increase in enzyme rate as the reaction proceeds. The acid-dependent term in R1 was removed in simulations where no feedback was included, recovering the classic Michaelis-Menten expression relevant for buffered systems. We do not include + and − charges in the figure and the rest of the model development for simplicity.

Case 1: Cell in a Reservoir
In Case 1, we considered a model of a single cell containing enzyme and buffer molecules that cannot cross the cell wall ( Figure 1a). The other small molecules can be transferred between the cell and the reservoir: where A i,r is the reservoir concentration of the ith species. The concentrations of S r and H r were taken as constant and N was considered as dilute in the reservoir such that concentrations N r and NH r were effectively zero. The transfer rate coefficient k i was given by: where P i is the permeability coefficient (m s −1 ) that depends on the nature of the cell membrane and the species and r j is the radius of the jth cell. The values of permeability coefficients vary widely. Here, in line with earlier work considering transport of species from synthetic lipid vesicles [22], we took k s = 0.0014 s −1 ; k H = k OH = 0.008 s −1 and k N = 0.004 s −1 . The rate equation for H and OH can be combined [22] and the concentration of NH determined from mass balance resulting in a five variable model for case 1:

Case 2: Group of Identical Cells in a Fixed Volume of Solution
In the second case, we considered a scenario where a group of cells had secreted a substance, C, which accumulated in a droplet around the cells, where [C o ] is the concentration of C in the droplet. The cells of volume V j were contained in a fixed volume, V T , of outer solution (the droplet) in which base-catalyzed gelation or precipitation of C o can take place ( Figure 1b): where [A i,o ] is the concentration of the ith species in the outer solution, g([A i,o ]) contains reaction terms, k r is the transfer rate from the reservoir and A i,r is the concentration of the ith species in the reservoir. For n identical cells the outer solution concentration is given by: where d = nV j /V T is a cell volume fraction: it gives the change in concentration of the species upon transfer from the cells of total volume nV j to the large outer volume V T . The rate equations for the outer solution concentrations were given by: where R2 o is the net rate of reaction 2 with the concentrations of N o and NH o in the outer solution. We ignore the loss of N o and NH o to the reservoir as these were found to have negligible effect on the dynamics. Equations (15)(16)(17)(18)(19)(20) were coupled with rate equations for the species S, H, N, NH, B and BH in the cells, resulting in a twelve variable model for Case 2.
All equations were solved using XPP-AUT with the integration method of cvode [23]. Bifurcation diagrams were produced using the AUTO program in XPP. The ".ode" files are available in the supplementary information.

Results
In the enzyme-catalyzed reaction considered here, feedback was imposed through the dependence of the enzyme activity on the pH. The role of the feedback in the dynamics displayed in a single and multi-cell system connected to an external solution is explored below.

Case 1: Cell in a Reservoir
In a cell with initial substrate [S] = 0.3 mM, no buffer, pH = 3.7 and no mass transfer between the cell and reservoir we observed a lag of approximately 600 s before the appearance of basic product N (Figure 2a). The pH switched from 4 to 9 and was accompanied by an acceleration in the rate of consumption of substrate S. In Figure 2b, no feedback was present i.e., the terms in R1 that depended on H were not included, and the initial concentrations were the same as in Figure 2a. The reaction was faster and there was the sudden appearance of product N after a time lag of 40 s accompanied by a pH switch from low to high. There was no acceleration in the rate of loss of substrate S when the pH rapidly increased. In the latter two examples the reaction proceeded until all of the initial substrate was consumed. With the addition of mass transfer between the cell and the external reservoir, oscillations in substrate, product and pH were observed (Figure 3a). The transport of OH and N played a small role in the dynamics: setting the values of k OH or k N to zero resulted in only a slight change in the period and amplitude of the oscillations. However, oscillations occurred for a narrow range of values of the ratio of k H /k S (5 to 6.67) in agreement with earlier work [22]. The influence of the addition of the buffer (BH + B) on the oscillations is shown in Figure 3b. With an initial [BH] of 0.03 mM, a reduction of the pH maximum was observed, from pH 8 to 7. When the BH concentration was increased above the substrate concentration of 0.3 mM, the pH oscillations were suppressed. The oscillations in the concentration of the buffer components in time is shown in Figure 3c for the case where initial [BH] = 0.12 mM.
A bifurcation diagram is plotted in Figure 4a

Case 2. Group of Cells in a Fixed Volume of Solution
In case 2, n identical cells were coupled to an outer solution (droplet) of fixed volume V T in which the pH can change and a base-catalyzed gelation or precipitation of species C took place. The effect of increasing the value of d = nV j /V T on the dynamics was explored. Assuming the cells have constant volume, V j , then an increase in d corresponded to an increase in the number density of cells.
In Figure 5a, with d = 0.01 and [S r ] = 0.3 mM and pH r = 4.1, the internal cell pH rose initially but fell to a low pH as a result of transport of acid from the surrounding solution. The outer solution pH did not rise above pH 4 and the extent of conversion (p Cm ) of C to gelled form, C m , remained low. When the number density of cells was increased (d = 0.07), the cell pH oscillated between pH 4.5 and 6.5 (Figure 5b). However, the outer solution pH o remained low and no significant gelation was observed. If the density of cells was further increased (d = 0.3), then the cell pH rapidly jumped to pH 9 and the outer solution pH o rose after a lag time (Figure 5c). The increase in pH o catalyzed the gelation process and there was complete conversion of the soluble species C o to insoluble gel C m . In Figure 6a, the extent of conversion, p Cm , at 5000 s was plotted as a function of increasing number density of cells. There was a critical value of d above which gelation was observed in the outer solution. Below d crit , the value of C m increased at such a low rate that no net gelation was obtained. The value of d crit is plotted as a function of [S o ] with and without feedback in Figure 6b, i.e., when the rate of the enzyme reaction R1 does not depend on acid. It is clear that feedback did not influence the pH switch and gelation in the outer solution. The feedback nevertheless played an important role in the internal cell dynamics. As d was increased, oscillations appeared resulting in a higher pH maximum with feedback than without (Figure 6c). The cell pH is plotted as a function of d in Figure 6d, where the black line shows the upper and lower values of pH in the oscillatory regime. For the pH to reach values above 6, the critical value of d required with feedback is almost half of that without feedback (Figure 6d dashed lines).

Discussion
Feedback plays a vital role in the functioning of living systems [1]. One of the important applications of feedback in bacteria is the release and detection of a small diffusible molecule, the autoinducer, in a phenomenon known as quorum sensing [24]. Quorum sensing is used by bacteria to activate population-level responses such as the formation of biofilm and a transition from single cell to multi-cellular behavior [3]. In biofilms, cells are bound together in a matrix of polysaccharides and other biomolecules that protect them from the external environment. Here, inspired by biofilm formation, we have developed a reaction network in a cellular system coupling feedback with other chemical processes such as gelation in the external solution. We chose to examine feedback through pH because pH is an important parameter that controls the rates of many processes. Our system is also based on experimentally-realizable processes that could be implemented in a synthetic cell.
In our simple model, the enzyme reaction took place in a cell coupled to an external solution containing substrate, S, and acid. Feedback arose through the characteristic bell-shaped rate-pH curve of the enzyme-catalyzed reaction. The reaction produced a weak base and the rate of the enzyme reaction depended on the pH, with a maximum at pH 7. In a cell with no external supply of species, we obtained a switch from low pH to high pH after a time lag accompanied by the sudden appearance of product, N, and rapid consumption of substrate. This is an example of pH-driven feedback with the basic product of the enzyme reaction acting as autocatalyst.
If the feedback was removed by making the enzyme rate independent of pH then we observed no rate acceleration in the loss of substrate. Nevertheless, the basic product, N, suddenly appeared after a time lag and there was a switch to high pH. In this acid-base reaction system, the initial acid reacted with the base and kept the pH low until sufficient base was produced to overcome the acid. This is similar to the substrate-depletive chemical clock reactions driven by inhibitor removal in the absence of feedback i.e., acid acted as an inhibitor suppressing formation of OH − [25].
When the mass transport of species from an external reservoir of constant concentrations was included, we obtained oscillations in pH. The differential transport of acid and substrate was essential for the observation of oscillations, i.e., the transport rate constant k H > k S . This is in agreement with earlier work on pH oscillations in some enzyme-catalyzed reactions [26][27][28]. Oscillations typically require positive feedback coupled with delayed negative feedback. The net rate of the transport of acid into the cell increased with increasing difference between the inner and outer pH, i.e., as the pH increased the rate of transport of acid into the cell increased. Therefore, the acid provided the necessary delayed negative feedback that terminated the autocatalytic, base-producing process.
There are few experimental examples of transport-driven instabilities in chemical systems, where the differential rates of transport of species cause oscillations in space or time [29]. The concept was first explored in Turing's pioneering modelling work on the basis of morphogenesis, in which spatial patterns were observed in a ring of coupled cells when the inhibitor transport was greater than that of the substrate in an autocatalytic process [30]. Differential transport may arise naturally in acid-base systems as acid, in the form H 3 O+, diffuses typically much quicker than other small ions in solution: the diffusion coefficient of acid can be up to nine times greater than small ions. In the model presented here, the chemical species were produced in a cell and crossed the cell wall with rates that were related to the permeability coefficients and radius of the cell. Permeability coefficients vary widely and depend on the nature of the membrane components, membrane thickness as well as the species involved [31]. In general, charged species do not cross cell membranes with appreciable rates. However, in a synthetic cell system constructed with phospholipids, fast acid transport could be facilitated by the incorporation of ion channels in the cell membrane [32].
Living cells are filled with buffering agents, so it is important to understand the influence that buffers might have on pH-driven feedback. We found that oscillations were, not surprisingly, completely suppressed for values of the buffer concentration greater than the substrate concentration. However, for lower concentrations the buffer primarily reduced the maximum of the oscillations. Thus, buffers could be used to control the amplitude of oscillations with a pH maximum less than 7.
We mapped the behavior of the single cell using bifurcation diagrams of pH as a function of substrate concentration in the reservoir, [S r ]. The transition from a low pH state to oscillations was obtained with feedback whereas only a single steady state pH was observed when no feedback was present. Nevertheless, we obtained a sharp switch in pH as [S r ] was increased in both cases. Thus the system may be considered ultrasensitive to small changes in substrate, even in the absence of feedback [1]. The acid acting as an inhibitor played an important role in generating this response.
In the model with multiple cells, we demonstrated that the increase in pH could be used to drive a base-catalyzed gelation or precipitation process in the external solution. This process was inspired by the quorum sensing of bacteria that leads to the formation of biofilm containing polysaccharides or the precipitation of inorganic minerals such as calcium carbonate, which occur in natural systems as the pH of the solution increases. We found that above a threshold density of cells, d = 0.3, the pH of the external solution increased, and gelation/precipitation occurred. A value of d = 0.01 corresponds to a number density of 2.7 million cells per ml of solution, assuming cells of radius 1 µm; comparable with bacterial cell cultures. The critical density of cells for a switch to high pH and gel formation was not influenced by the presence of feedback, rather the inhibition by acid played the major role in this switch-like behavior.
We note that the behavior here was obtained for relatively low transfer rates in which a decoupling between the internal and external system was observed: there was a difference in the timescales of the pH dynamics in the cell and the pH dynamics in the outer solution. No large amplitude oscillations in pH in the external environment were obtained for any values of the parameters explored. The transfer coefficients ranged from 0.0014-0.008 s −1 ; with a cell of the order 1 µm, this gives permeability coefficients from 10 −9 -10 −8 m s −1 , on the low side for transport of neutral molecules in synthetic vesicle systems. It may be that for larger transport coefficients the feedback plays a greater role in the switch for gelation. However, we demonstrated that the feedback played an important role in the internal cell dynamics. The emergence of oscillations above a critical density of cells resulted in a much higher cell pH than without feedback. Thus, we conclude that the presence of feedback could be used to drive an internal switch that activates other pH-regulated processes within the cell at low substrate or cell densities.

Conclusions
The design of reaction networks that mimic the feedback-driven behavior of biological systems remains a formidable challenge. There has been increasing interest in obtaining quorum sensing type behavior in synthetic systems [5,12,[33][34][35][36][37]. The underlying processes are vastly different however all examples exploit seemingly similar principles of an autocatalytic signaling process and communication between entities via a common surround. Here, taking inspiration from bacteria, we have proposed a model to determine the role of feedback in certain density-dependent behavior. The model involved an enzyme catalyzed reaction in cells which influenced the pH and the dynamics of an external pH-regulated process: the base-catalyzed formation of a gel or precipitate.
We found that a sharp switch in pH was obtained in the multi-cellular system above a critical number of cells and gelation/precipitation did not require the presence feedback. However, we only observed pH oscillations with feedback and the oscillations resulted in a higher pH maximum in the cells. We conclude that feedback may be used to activate internal cellular processes at low substrate concentrations and low cell densities, even when it plays no role in the external dynamics. Although the chemistry of the underlying feedback processes was not necessarily relevant to quorum sensing, similar dynamical responses may arise in cells. This model system was also based on processes that might be used for synthetic quorum sensing analogues in biotechnology applications [38,39].