## New Research In

### Physical Sciences

### Social Sciences

#### Featured Portals

#### Articles by Topic

### Biological Sciences

#### Featured Portals

#### Articles by Topic

- Agricultural Sciences
- Anthropology
- Applied Biological Sciences
- Biochemistry
- Biophysics and Computational Biology
- Cell Biology
- Developmental Biology
- Ecology
- Environmental Sciences
- Evolution
- Genetics
- Immunology and Inflammation
- Medical Sciences
- Microbiology
- Neuroscience
- Pharmacology
- Physiology
- Plant Biology
- Population Biology
- Psychological and Cognitive Sciences
- Sustainability Science
- Systems Biology

# Role of metabolic spatiotemporal dynamics in regulating biofilm colony expansion

Edited by William Bialek, Princeton University, Princeton, NJ, and approved March 12, 2018 (received for review May 5, 2017)

## Significance

Oscillatory dynamics commonly arises in a variety of multicellular biological systems. Bacterial colonies exploit such oscillations to control the interplay between growth and resource availability. We model a recent experiment that observes the oscillations of the growth rate during the expansion of a bacterial colony of *Bacillus subtilis*, elucidating the origin of the oscillations in terms of the spatiotemporal dynamics of the metabolic interactions between cells within the biofilm. Further, different bacterial cell populations are required in this model for the oscillatory behavior to arise, therefore demonstrating they are necessary for the colony survival. This approach provides a platform to model a large class of biological phenomena involving the formation of large aggregates of cells and/or a heterogeneous cell population.

## Abstract

Cell fate determination is typically regulated by biological networks, yet increasing evidences suggest that cell−cell communication and environmental stresses play crucial roles in the behavior of a cell population. A recent microfluidic experiment showed that the metabolic codependence of two cell populations generates a collective oscillatory dynamic during the expansion of a *Bacillus subtilis* biofilm. We develop a modeling framework for the spatiotemporal dynamics of the associated metabolic circuit for cells in a colony. We elucidate the role of metabolite diffusion and the need of two distinct cell populations to observe oscillations. Uniquely, this description captures the onset and thereafter stable oscillatory dynamics during expansion and predicts the existence of damping oscillations under various environmental conditions. This modeling scheme provides insights to understand how cells integrate the information from external signaling and cell−cell communication to determine the optimal survival strategy and/or maximize cell fitness in a multicellular system.

- biofilm expansion
- phenotypic differentiation
- metabolic codependence
- oscillations
- reaction−diffusion system

Great efforts have been spent in recent years toward understanding the role of gene regulatory networks in modulating the gene expression dynamics of cells (1⇓⇓–4). In these studies, the dynamic behavior of gene expression for a cell can be modeled by a set of coupled chemical rate equations that represent the combinatorial interactions among genes. Cell−cell communication and cellular microenvironment can, however, substantially affect this single-cell picture, such as in the cases of bacterial chemotaxis (5, 6), Notch signaling in normal development and tumorigenesis (7, 8), and the interplay among cell populations in synthetic systems (9).

The oscillatory expansion of bacterial colonies recently observed in a biofilm system exemplifies how intercellular communication plays a central role (10). Liu et al. (10) placed a colony of the Gram-positive organism *Bacillus subtilis* in a growth chamber that confines biofilm expansion within a 2D space, while a microfluidic device constantly provides glutamate and washes away metabolic products diffusing outside of the biofilm. Interestingly, the biofilm initially undergoes rapid and steady expansion until it reaches a certain size, when it exhibits oscillatory growth rate with periodic halting. Detailed analysis revealed that the spatial regulation of the biochemical pathways for cell metabolism is critical for the formation of this specific dynamics. Certainly, bacterial colonies must develop efficient communication mechanisms to react effectively to various environmental stresses and optimize the overall colony fitness.

The biofilm dynamics can be discerned by examining the glutamate/glutamine synthesis pathway in bacterial cells (Fig. 1). Since glutamate is the only available nitrogen source, the proliferation is supported by glutamine, which is synthesized from both glutamate and ammonium (11). Without ammonium in the growth media, cells need to produce ammonium from glutamate catalyzed by the glutamate dehydrogenase (GDH) enzyme (12). To explain the oscillatory dynamics of biofilm expansion, Liu et al. (10) proposed a model of metabolic codependence. According to their experimental observations, bacterial cells have the following two spatially distinct metabolic phenotypes. Cells at the internal region of the colony (interior) are responsible for ammonium synthesis by means of glutamate diffusing from periphery. Cells at the external part of the biofilm (periphery) do not produce ammonium by their own due to a futile cycle: if ammonium was synthesized in the periphery, it would rapidly diffuse to the growth media, making the process especially inefficient and thus unlikely (13, 14).

Therefore, peripheral cells exploit the ammonium produced in the interior and the glutamate diffused from the media to synthesize glutamine and support cell growth. The proliferation of peripheral cells consumes substantial amounts of glutamate and ammonium, leaving interior cells without substrate for ammonium synthesis. Since the interior stops ammonium production, peripheral cells do not have sufficient nutrients for growth. When growth halts, the external glutamate from the media becomes again accessible to the internal cells. Thus, the production of ammonium is restarted, leading to another round of growth.

To test this proposed mechanism, Liu et al. (10) established a quantitative model by considering two discrete spatial components: the internal cell population and the external cell population. The dynamics of the biofilm were described by a set of ordinary differential equations governing the time evolution of the size of the cell populations and the levels of the metabolites within these two regions. This design can capture the essential traits of the biofilm dynamics and reproduce the experimentally observed oscillations in the growth rate of the colony. It also correctly predicted the dynamics of biofilm growth under various conditions of growth media.

Much remains to be understood, however, regarding how these metabolic reactions are regulated by the spatial organization of the organisms inside the community. Without a direct characterization of the diffusion of the nutrient molecules in the biofilm, the original model is not capable of explaining why oscillations are detected only when the size of the colony exceeds a certain threshold and, more generally, how the metabolic interactions are regulated in different regions of the biofilm. These questions are hard to address without a more comprehensive framework that describes, at different spatial locations, both the temporal dynamics of the metabolic activities and the proliferation of the cells.

We introduced an improved scheme to study both temporal and spatial dynamics of a biofilm. Compared with the description proposed by Liu et al. (10), this model considers the diffusion of small metabolites by incorporating the spatial dynamics of the bacterial colony. The biofilm is described as a continuous 2D system where the more external layer represents bacteria with the peripheral phenotype and the more internal layer represents bacteria with the internal phenotype. These improvements allow us to investigate the spatial organization of cells with different phenotypes, and directly explore the repercussions of the glutamate/glutamine synthesis pathway on the biofilm development. The proposed model can explain the recurrent oscillatory cycles of the growth rate in terms of the space-dependent interplay between the internal and peripheral phenotypes and reproduce the observed growth dynamics in the presence of altered conditions of the growth media. Moreover, the occurrence of oscillations is insensitive to the radius of the biofilm; instead, it is sensitive to the width of the peripheral layer. These findings suggest that the onset of the oscillation of the biofilm expansion rate is due to the switch of the bacterial cells from interior to peripheral phenotype, and it is specifically triggered when the peripheral layer increases its width to a certain level. Finally, various types of growth dynamics, including dampened, stable, dissipating dynamics, are revealed by varying the ratio between interior and peripheral cells and modifying the biofilm’s external conditions.

## Results

### A Spatial−Temporal Model of Metabolic Interplay in the Biofilm.

In the original experiment by Liu et al. (10), a special microfluidic device was designed to limit biofilm expansion to two dimensions. With a growth media containing nutrients such as glutamate, *B. subtilis* grows and eventually forms a colony that is approximately centrosymmetric as no shear was applied to the system. Therefore, we described the biofilm dynamics in a 2D space, but the current framework can be naturally generalized to a 3D system. Due to the symmetry, we further simplified the model with a polar coordinate reference frame centered in the biofilm and considered the spatial dynamics along the radial coordinate only, hence reducing the biofilm to a 1D system.

We implemented the glutamate metabolic pathway via a reaction−diffusion system that extends the previous modeling of Liu et al. (10) and describes the spatiotemporal dynamics of ammonium, glutamate, GDH, and housekeeping proteins (Eqs. **1**–**4**). Ammonium and glutamate can diffuse through the biofilm and its boundaries, while GDH and housekeeping proteins are not diffusible, because they represent the specific levels inside bacterial cells. Specifically, ammonium diffuses faster than glutamate (15, 16). Nutrient availability in the growth media was modeled via boundary conditions at the biofilm boundary. The cell metabolic phenotype is determined by a sigmoid signal that depends on the radial distance of the cell from the biofilm center (Eq. **5** and *SI Appendix*, Fig. S1). Cells within the inner layer, the interior, can activate the GDH enzyme and therefore catalyze ammonium production when there is sufficient glutamate. Conversely, cells in the external layer, the periphery, cannot synthesize the nutrient, due to their futile cycle of ammonia.

### Metabolic Oscillations in the Biofilm.

Although our main objective was to investigate the growth process of the biofilm, we studied a simplified “nongrowth model” first to understand the governing principles of the glutamate metabolic pathway. This reduced model considered the dynamics of the chemical species and increased nutrients’ consumptions for fast-growing cells, but disregarded biofilm expansion.

First, we discovered that the nongrowth model can exhibit different oscillatory dynamics by varying the sharpness of the transition between internal and peripheral phenotypes and the diffusion coefficients of the nutrients (*SI Appendix*, Fig. S2). Specifically, tuning the diffusion coefficient of ammonia and glutamate to measured values (15, 16) led to metabolic oscillations spread through the whole biofilm (Fig. 2 *A* and *B*). Further, such oscillations presented the typical traits of nutrients exchange and management reported by Liu et al. (10). The concentration of ammonium was higher in the interior, where the metabolite is produced, while the glutamate level was larger in the periphery due to the high external concentration (Fig. 2 *C* and *D*). Housekeeping proteins synthesized from glutamate and ammonium were mostly observed in the periphery (*SI Appendix*, Fig. S3).

Notably, a smooth interior−periphery phenotypic transition suppressed oscillations, strongly suggesting the need for two distinct cell phenotypes. Also, oscillations were insensitive to biofilm size, but required an intermediate size of the peripheral layer (*SI Appendix*, Fig. S4). Further, the system was robust against local variation of biochemical interaction parameters due to the architecture of the metabolic circuit (17) (*SI Appendix*, Figs. S4 and S5). Lastly, the oscillation period increased with the peripheral layer’s width, while it decreased with the glutamate diffusion coefficient (*SI Appendix*, Fig. S5), hence indicating that how fast glutamate diffuses through the periphery determines the timescale of metabolic oscillation.

### The Onset of Biofilm’s Oscillations.

Next, we developed a theoretical framework for the growth of the biofilm (Eqs. **6**–**8**). In the proposed scheme, the production of new biomass correlates to the level of housekeeping proteins. Specifically, we hypothesized that a basal level of proteins is employed in ordinary cellular activities and it is thus not available for the growth. The biofilm growth was coupled with the metabolism of glutamate by solving the reaction−diffusion equations of the glutamate biochemical network on the expanding biofilm area. To simulate the creation of new cells, the growth algorithm first computed the newly produced biomass and then updated the radius of the biofilm accordingly (see *Methods* for details).

Furthermore, we assumed that the peripheral layer grows proportionally to the radius of the biofilm when the colony is small but saturates to a constant value when the colony is large (Fig. 3*A* and Eq. **9**). This hypothesis ensures that there are enough peripheral cells to effectively sustain the growth at all stages of biofilm development and prevents the oscillation quenching observed for an excessively large periphery.

The model correctly reproduced the onset of oscillations and replicated the stable oscillations observed afterward (Fig. 3*B*). In the small biofilm, the glutamate diffuses through the thin peripheral layer and guarantees an adequate nutrient supply to interior cells. Under these circumstances, the growth rate does not oscillate. Once the critical colony diameter is exceeded, the trajectory of the growth rate presents a continuous transition to oscillations whose amplitude evolves from zero to a nearly steady value while the diameter of the biofilm—and thus the size of the peripheral layer—increases. In this regime, the time required by glutamate to diffuse through the peripheral layer is too long to directly supply the interior, as this nutrient is also employed in the periphery to synthesize housekeeping proteins. The possible starvation of interior cells is resolved through periodic decreases of the growth rate, the oscillation minima, that allow some glutamate to diffuse to the interior before being consumed by peripheral cells. Once the internal level of ammonium has been restored, the biofilm can start a new stage of fast growth, which in turn consumes the interior ammonium and sets the conditions for a new cycle. As the growth continues, the trajectory approaches a limit cycle and the growth rate oscillates with a nearly constant amplitude and a period around 2.5 h (Fig. 3*C*), in agreement with experiments (10).

We further tracked the growth rate in different points of the biofilm and showed that peripheral cells give the larger contribution to biofilm expansion (Fig. 3 *D*–*F*).

The growth trajectory was robust against local variation of the diffusion constants and the parameters of the growth model (*SI Appendix*, Fig. S6 *A*–*D*). Additionally, stress factors such as environmental pollution, pH variation, or waste accumulation could modulate cellular growth (18, 19). We modeled such factors as (*i*) decrease of the cellular growth rate constant or (*ii*) increase of housekeeping proteins threshold for cellular growth. Both approaches resulted in a slower growth but with robust oscillatory behavior (*SI Appendix*, Fig. S6 *E* and *F*).

### Oscillations Depend on the Peripheral Layer.

Different growth dynamics were observed (Fig. 4) by changing the maximum width of the peripheral layer (the parameter *A*). This parameter controls the fraction of cells with a peripheral phenotype in the biofilm. If the periphery did not reach the transition width for oscillation onset (vertical dashed line of Fig. 4*A*), the growth rate did not exhibit oscillations (Fig. 4 *A* and *B*). For intermediate values of *A* and *C*). As the parameter *A* and *D*). At even larger values, the oscillations could completely stop after some oscillatory cycles (Fig. 4 *A* and *E*).

According to the experimental observations of Liu et al. (10), the oscillations remain stable for a long period (more than a day), suggesting that an intermediate value of the parameter yields stable oscillations such as in Fig. 4*C*. This regime is characterized by a competition between the maximal growth and the effective protection of the interior cells, and supports the interesting interpretation of the oscillations as a trade-off between colony growth and sheltering of the interior population proposed by Liu et al. (10).

### Strategies to Control Biofilm Growth via Metabolic Interplay.

Ammonium administration in the growth media quenched oscillations as the periodical dearth of ammonium in the periphery was balanced by environmental supply (Fig. 5*A*). This finding agrees well with experimental observations (10). Interestingly, simulating biofilm growth for different ammonium levels at the biofilm boundary (*B*). Further, a low value of *C*). Similarly, GDH induction in the periphery quenched the oscillations via removing phenotypic diversification, as observed experimentally (10) (*SI Appendix*, Fig. S7).

A higher concentration of glutamate in the growth media (*A*) for periphery expansion is a linear function of glutamate (Eq. **10**). We further hypothesized that the deeper glutamate penetration inside the biofilm could allow a larger peripheral cell population, and therefore assumed that the peripheral layer could expand more (Eq. **11**). Under these assumptions, we observed faster biofilm growth and larger oscillation amplitude upon glutamate increase on a twofold change around the original concentration of 30 mM (15 mM to 60 mM) (Fig. 5*D*). The increased nutrient availability shifted the oscillation onset to a larger biofilm size (Fig. 5*E*), in good agreement with experiments (20). Intriguingly, the model predicted a maximal oscillation period at intermediate *F*). This finding can be interpreted as a balance of the following processes: a larger *SI Appendix*, Fig. S8).

## Discussion

We introduced a mathematical description to investigate the glutamate metabolic pathway in cells and explain the oscillatory growth of bacterial colonies of *B. subtilis* recently observed by Liu et al. (10). Improving the previous modeling scheme, this work regarded the biofilm as a continuum domain and investigated the spatiotemporal regulation of the metabolic interactions with a set of coupled reaction−diffusion equations. We considered a biofilm composed of two distinct cell phenotypes. The interior population accounts for ammonium production in the presence of active GDH enzymes, while the peripheral population utilizes the ammonium produced in the interior and the external glutamate to sustain colony growth. The phenotypic transition was based on the futile cycle of ammonia loss that makes ammonium homeostasis costly for peripheral cells. The choice of a space-dependent framework was advantageous as it allowed us to explicitly describe the diffusion of metabolites in the colony. Secondly, the metabolic interactions could be evaluated in different parts of the system, allowing us to precisely address the phenotypic distribution and/or heterogeneity of the biofilm. Therefore, this modeling scheme provides a perspective on how external signaling and cell−cell communication contribute to the dynamics of intracellular biological networks.

The spatiotemporal model can reproduce the observed growth dynamics of a *B. subtilis* biofilm. Specifically, improving the preexisting scheme, this model can predict the onset of the oscillations. When the size of the biofilm is small, the growth does not present oscillations, as diffusion ensures nutrient availability everywhere in the colony. As the biofilm expands, the dual role of glutamate, required both for ammonium synthesis in the interior and new biomass creation in the periphery, generates a conflict that requires an ad hoc management of the available resources. These competitive processes are resolved with oscillations of the growth rate, where the periodic halting of the expansion allows the glutamate to diffuse in the biofilm interior. The variation of the peripheral layer size led to a rich ensemble of growth dynamics, including stable growth, stable oscillations, and partial and complete damping. Currently available experiments suggest that oscillations can remain stable for a long period (10). It has been suggested that the metabolic interplay between interior and peripheral populations represents an optimal trade-off that can resolve the conflict between sheltering interior cells and promote growth of the colony at once (10, 21). This model’s predictions strengthen this hypothesis, since the stable oscillations correspond to an intermediate size of the peripheral layer which balances a maximal growth rate, achieved with a small periphery, and an effective protection of interior cells, accomplished with a large periphery. The validation of these predictions would further improve our understanding of the organization and survival strategies of bacterial communities.

Biofilm proliferation is certainly regulated by biophysical processes such as local variations of pH, accumulations of waste products, or pollution (18, 19) that we modeled as simple rescaling of the growth parameters. However, understanding how these processes impact the biofilm dynamics represents an exciting challenge at both the experimental and theoretical levels.

Seemingly counterintuitive, providing excessive nutrients suppresses the metabolic interplay of the colony: the supplement of ammonium allowed peripheral cells to grow autonomously, hence decoupling them from the interior. This “overfeeding” strategy results in strong nutrient consumption in the periphery and a consequent starvation of the interior. Although sporulation can be a response to growth slowdown and nutrient shortage (22, 23), such an energy-consuming process can be quite challenging in the dense biofilm environment if the level of nutrients diminishes too rapidly. As suggested by Liu et al. (10), this could thus represent an effective strategy to suppress the interior population of the colony. Investigating, both experimentally and theoretically, cell fate dynamics in response to overfeeding could play a major role in developing new therapies to fight biofilm infections. Further, a low exposure to ammonium conserved oscillations but could increase the oscillation onset compared with the ammonium-free biofilm.

Increasing the glutamate level in the growth media resulted in a larger growth rate. Nonetheless, metabolic oscillations were still observed, because ammonium was still self-produced inside the biofilm. Additionally, the oscillation onset increased with the external glutamate, in agreement with experiments (20), while the oscillation period was maximized at an intermediate exposure. Quantitatively validating these predictions would require tracking biofilm growth for different nutrients conditions, which, in turn, would offer further data to calibrate the mathematical model.

This work presented a theoretical formalism that can go beyond the expansion of bacterial biofilm and may be applied to a broader class of biological systems that involve the competition over limited resources, the development of large aggregates of cells, and/or a heterogeneous cell population. Cells within the biofilm integrate the information from cell−cell signaling and microcellular environment to design the optimal survival strategy and maximize cell fitness. In the considered resource-limited environment, the bacterial community resolves the shortage of nutrients through collective metabolic oscillations. Cell fate determination is influenced by intrinsic fluctuations in gene expression dynamics and extrinsic noise resulting from the surrounding microenvironment at a single-cell level (24⇓–26). Signaling mechanisms among neighboring cells represent a driving force, as well, in establishing phenotypic expression (21, 27), especially when external conditions force microorganisms to cooperate to maximize the overall fitness level (28, 29). The traditional population-based approach already provided important insights in the system biology modeling over the last decades (30, 31). This approach, however, needs to be expanded to include the spatial dynamics to capture the heterogeneity of the cell−environment and cell−cell interactions (32, 33). Here, the spatiotemporal modulation of the glutamate metabolism in *B. subtilis* was addressed within a multicellular system, thus taking a step further toward understanding the connection between cell−cell signaling, metabolism, and biofilm expansion.

## Methods

The reaction scheme presented in Fig. 1 is implemented via a reaction−diffusion system that describes the spatial and temporal dynamics of the concentrations of the four considered variables: ammonium (*A*), glutamate (*G*), GDH (*H*), and housekeeping proteins (*R*), as a function of the radial position

Ammonium is synthesized from glutamate, and the reaction is catalyzed by the enzyme GDH, so the production rate of ammonium is expressed as

where

The conditions at biofilm boundary (

The growth rate of biofilm area (

where

where

The thickness of the peripheral layer

It expands linearly for small

where

with glutamate threshold

The algorithm for the growth of the biofilm includes three steps. First, the newly produced biomass is computed by solving Eqs. **1**–**4** while the radius of the biofilm **6**–**8**. Lastly, the width of the peripheral layer **9**–**11**, and the interior radius **1**–**4** are solved. Therefore, the concentrations of chemicals are rescaled to the new biofilm area to ensure the same amount of total nutrients in the biofilm area before and after the growth step; *X*_{new} = *X*_{old} *S*_{old}^{2}/*S*_{new}^{2}, where *X* = *A*, *G*, *H*, *R*, and *S* is the radius of the biofilm. The code used for the numerical solution of the model is freely available on the github page of F.B. (https://github.com/federicobocci91).

The numerical values of the parameters and further details on the mathematical model are presented in *SI Appendix*, sections S1–S5.

## Acknowledgments

We have benefited from fruitful conversations with Gurol Suel and Jintao Liu. The work at the Center for Theoretical Biological Physics was sponsored by the National Science Foundation (Grants PHY-1427654 and MCB-1241332). M.L. is partially supported by the National Cancer Institute of the National Institute of Health under Award P30CA034196.

## Footnotes

- ↵
^{1}To whom correspondence should be addressed. Email: jonuchic{at}rice.edu.

Author contributions: F.B., Y.S., M.L., and J.N.O. designed research; F.B., Y.S., and M.L. performed research; F.B., Y.S., and M.L. analyzed data; and F.B., Y.S., M.L., and J.N.O. wrote the paper.

The authors declare no conflict of interest.

This article is a PNAS Direct Submission.

Data deposition: The code used for the numerical solution of the model is freely available at https://github.com/federicobocci91.

This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1706920115/-/DCSupplemental.

Published under the PNAS license.

## References

- ↵
- ↵
- ↵.
- Ben-Jacob E,
- Lu M,
- Schultz D,
- Onuchic JN

- ↵
- ↵.
- Ellison D, et al.

- ↵.
- Levine H,
- Kessler DA,
- Rappel W-J

- ↵
- ↵.
- Bocci F, et al.

- ↵.
- Chen Y,
- Kim JK,
- Hirning AJ,
- Josić K,
- Bennett MR

- ↵
- ↵
- ↵.
- Belitsky BR,
- Sonenshein AL

*Bacillus subtilis*glutamate dehydrogenase genes. J Bacteriol 180:6298–6305. - ↵
- ↵
- ↵.
- Stewart PS

- ↵
- ↵.
- Huang B, et al.

- ↵
- ↵.
- Asally M, et al.

- ↵.
- Liu J, et al.

- ↵
- ↵.
- Narula J, et al.

- ↵.
- Schultz D,
- Wolynes PG,
- Ben Jacob E,
- Onuchic JN

- ↵.
- Maamar H,
- Raj A,
- Dubnau D

*Bacillus subtilis*. Science 317:526–529. - ↵.
- Elowitz MB,
- Levine AJ,
- Siggia ED,
- Swain PS

- ↵.
- Narula J,
- Devi SN,
- Fujita M,
- Igoshin OA

*Bacillus subtilis*sporulation decision. Proc Natl Acad Sci USA 109:E3513–E3522. - ↵
- ↵
- ↵
- ↵.
- Kussell E,
- Leibler S

- ↵.
- Chuang JS,
- Rivoire O,
- Leibler S

- ↵
- ↵

## Citation Manager Formats

### More Articles of This Classification

### Biological Sciences

### Systems Biology

### Related Content

- No related articles found.

### Cited by...

- No citing articles found.