Abstract
We investigate the spatiotemporal dynamics of a bacterial colony model. Based on the stability analysis, we derive the conditions for Hopf and Turing bifurcations. Furthermore, we present novel numerical evidence of time evolution of patterns controlled by parameters in the model and find that the model dynamics exhibit a diffusion controlled formation growth to spots, holes and stripes pattern replication, which show that the bacterial colony model is useful in revealing the spatial predation dynamics in the real world.
1. Introduction
Spatial patterns which are formed by some kinds of bacterial colonies present an interesting structure during their growth conditions. In particular, colonies of bacterium bacillus subtilis can present a rich variety of structures [1–13]. The nature of the pattern exhibited depends on the particular bacterial species used and the environmental conditions imposed. Ohgiwari et al. [11] have shown that for a nutrient-poor solid agar, the bacterium colonies exhibit fractal morphogenesis similar to diffusion-limited aggregation (DLA). For softer agar medium, the colonies tend to show a dense-branching morphology (DBM) [7]. If both the nutrient concentration and the agar’s softness further increase, simple circular colonies grow almost homogeneously in space [14].
There are many mathematical models for explaining each characteristic colony pattern. Kawasaki et al. [7] have developed a reaction-diffusion model and have shown the patterns by using the computer simulations. Since in Kawasaki et al.’s model, all the nutrients must be consumed; L. Braverman and E. Braverman [4] have introduced a model of prey-predator type with Holling-II functional response under the situation of a renewable nutrient. In the present paper, motivated by the work of L. Braverman and E. Braverman, we consider the model with the consumption term of nutrient in a Holling III functional response.
Let us denote by and the nutrient concentration and the density of the bacterial cells at point , respectively. We consider the following system: where is the intrinsic nutrient growth rate, is the carrying capacity of the environment for the nutrient (prey), is the bacteria (predator) mortality rate, , , and are parameters of the Holling Type III functional response, and is the nutrient diffusion coefficient. Following [4, 7], we assumed that the diffusion coefficient is proportional to both nutrient and bacteria densities
Here we try to model the situation of a renewable nutrient. Then the system involves two reaction-diffusion equations of a predator-prey type with a Holling Type III functional response. Diffusive predator-prey systems were extensively studied; we mention here the recent papers [15–19], the monograph [20], and the references therein. In the present paper, it is to investigate the spatial pattern formation of system (1) which means the convergence of solutions to some stable spatially-in-homogeneous pattern as time tends to infinity. And in natural science, the pattern formation can reveal the evolution process of the species; it is, perhaps, the most challenging in modern ecology, biology, chemistry, and many other fields of science [21–38]. Thus, our basic concern is to find, if any, a spatially inhomogeneous equilibrium and periodic solutions that are stable in a certain sense. From the pioneer work by Turing [12], it is widely known that a reaction-diffusion system exhibits Turing instability if the homogenous steady state is stable to small perturbations in the absence of diffusion but unstable to small spatial perturbations when diffusion is present which implies the existence of spatially in-homogenous solutions. From the Hopf bifurcation analysis and the phrase transition theory developed by Ma and Wang [39–42], it is shown that the periodic solutions exist [43].
The paper is organized as follows. In Section 2, we give the analysis of the model and mathematical setup. In Section 3, we analyze the spatial model, we derive the conditions of the Turing bifurcation and Hopf bifurcation, and we give the existence of periodic solution. We give some computer simulations to illustrate the emergence of pattern formation in Section 4. Finally, some conclusions are given.
2. Modeling Analysis and Mathematical Setup
To obtain the dimensionless form of the system (1), we introduce the following:
Omitting the primes, we obtain the following nondimensional form of (1): with , .
Model (4) is to be analyzed under the following nonzero initial conditions: and Neumann boundary conditions: In the above, and denote the size of the system in square domain and is the outward unit normal vector of the boundary . The main reason for choosing such boundary conditions is that we are interested in the self-organization of the pattern and the Neumann conditions imply no external input [22].
It is known that only nonnegative solutions of (4) have biological significance. System (4) has two spatially homogeneous stationary solutions:(1)the bacteria-free equilibrium which implies that the nutrient is at the carrying capacity level;(2)coexistence equilibrium which represents a uniform distribution of bacteria, where and with and .
To consider the pattern formation of (4) from we make the translation Then, (4) are rewritten as where and , , , and are terms of high order.
Define two Hilbert spaces Then is dense and compact inclusion. where for .
Furthermore, denote that with where Here and are terms of high order.
Then are a family of parameterized bounded operators continuously depending on the parameter such that .
Then (9) can be written in the following operator form:
3. Bifurcation Analysis
Unless otherwise specified, in this section, we require that always exist; that is, and .
Consider the following eigenvalue problem of system (9): with the Neumann boundary condition (6).
Let and be the th eigenvalue and eigenvector of the Laplacian with Neumann boundary condition and with , .
Denote by the matrix given by Thus, all eigenvalues of (18) satisfy where is the eigenvector of corresponding to and is expressed as with Hence, the eigenvector of (18) corresponding to is where is as in (19).
3.1. Hopf Bifurcation Analysis
It is clear that with if and only if
Thus, we introduce one critical number where such that attains its minimum values. Consider
Theorem 1. Let be the number given in (26) such that (27) is satisfied. Then and are a pair of first complex eigenvalues of (18) near , and
3.2. Periodic Solution from Hopf Bifurcation
By Theorem 1, problem (4) undergoes a dynamic transition to a periodic solution from . To determine the types of transition we introduced a parameter as follows: where
Theorem 2. Let be the number given by (29), then the problem undergoes a transition to periodic solutions at , and the following assertions hold true.(1)When , the transition is continuous and the system bifurcates to a periodic solution on which is an attractor.(2)When , the transition is jump and the system bifurcates to a periodic solution on which is a repeller.
Proof. We will verify this theorem by using Theorem A.3 in [44]. The eigenvalues at in are given by . The eigenvectors and corresponding to satisfy
It is easy to see that
The conjugate eigenvectors and satisfy
It is easy to check that
It is known that functions and are given by
Because the first eigenvector space of (18) with (6) is invariant for the equations (4) with (6), the center manifold function vanishes; that is,
Therefore, we derive from (32) to (35) that
where
From the focus values in [39, 40, 43], we have that
is the same as in (29). Hence, by Theorem A.3 in [3] the system bifurcates from to a periodic solution; thus the proof is complete.
Remark 3. As an example, let , , , , and , then from (29), we compute . From Theorem 2, we can conclude that the transition is jump and the system bifurcates to a periodic solution on which is a repeller (see Figure 1).
3.3. Turing Bifurcation Analysis
In this subsection, we will state the Turing instability for the positive equilibrium of model (1). Mathematically speaking, the positive equilibrium is Turing instability, which was emphasized by Turing in his pioneering work in 1952 [12]. The Turing bifurcation occurs when and the wave-number satisfies Hence, Turing instability occurs when the condition either or is violated.
Since the positive equilibrium is stable without diffusion means that and hold, then is always true. Therefore, for the emergency of the diffusion-driven instability in model (1), it is needed for some . A necessary condition is otherwise for all since and . And we notice that achieves its minimum at the critical value where Summarizing the above calculation, we conclude.
Theorem 4. If then the positive equilibrium of model (1) is Turing unstable.
In Figure 2, based on the results of Theorem 4, we show the dispersal relation of with . The blue and red curves represent Hopf and Turing bifurcation curves, respectively. They separate the parametric space into three domains. The outside domain of the Hopf bifurcation curve is stable and the inside domain of the Turing bifurcation curve is unstable. Hence, among these domains, only the domain(III) satisfies the conditions of Theorem 4 and we call domain(III) as Turing space, where the Turing instability occurs and the Turing patterns may be undergone.
4. Pattern Formation
In this section, we perform extensive numerical simulations of the spatially extended model (4) in 2-dimensional spaces, and the qualitative results are shown here. Our numerical simulations employ the nonzero initial (5) and the zero-flux boundary conditions (6) with a system size of , with discretized through and , with . Other parameters are fixed as The numerical integration of (4) was performed by fourth-order Runge-Kutta scheme integration [45], with a time step of , and by using the standard five-point approximation for the 2D Laplacian with the zero-flux boundary conditions [46, 47]. More precisely, the concentrations at the moment at the mesh position are given by where is defined in (17).
Initially, the entire system is placed in the steady state , and the propagation velocity of the initial perturbation is thus on the order of space units per time unit. And the system is then integrated for 200 000 time steps, and the last images are saved. After the initial period during which the perturbation spreads, either the system goes into a time-dependent state or to an essentially steady state (time independent).
In the numerical simulations, different types of dynamics are observed and it is found that the distributions of predator and prey are always of the same type. Consequently, we can restrict our analysis of pattern formation to one distribution. In this section, we show the distribution of prey , for instance. And the parameters are located in the Turing space (cf., Figure 2), the region where Turing instability occurs. We have taken some snapshots with red (blue) corresponding to the high (low) value of prey .
Figure 3 shows the process of pattern formation for model (4) with and . In this case, the pattern takes a long time to settle down, starting with a homogeneous state (cf., Figure 3(a)), and the random perturbation leads to the formation of stripes and spots (cf., Figures 3(b) and 3(c)) and ends with stripes only (cf., Figure 3(d)), which is time independent.
(a)
(b)
(c)
(d)
In Figure 4, we show two spots-stripes patterns obtained with model (4) at 100 000 iterations; that is, . These two patterns are similar to each other. With , in this case, the equilibrium is and the spots-stripes pattern is relatively high (cf., Figure 4(a)), while with , the equilibrium is , at low prey densities (c.f., Figure 4(b)).
(a)
(b)
In Figure 5, we show the interesting and similartime-independent patterns which obtained by model (4) at 200 000 iterations. They consist of blue/red spots on a red/blue background. We refer to them as spots (cf., Figure 5(a)) and holes (cf., Figure 5(b)), respectively. In Figure 5(a), with , , the hot spots are isolated zones with high prey densities. In this case, the predators are in low density obviously. While with , , holes are isolated zones with low prey density (Figure 5(b)). In this case, the predators are in high density. From Figure 5(b), one can see that the predators apparently almost occupy the whole spatial domain.
(a)
(b)
5. Concluding and Remarks
In this paper, pattern formation of a spatial model for the growth of bacterial colonies with the two-dimensional space is investigated. Based on both mathematical analysis and numerical simulations, we have found that its spatial pattern includes periodic solutions from Hopf bifurcation and the spotted and striped patterns from Turing bifurcation.
It should be noticed that, if considered in a somewhat broader ecological perspective, our results have an intuitively clear meaning; there has been a growing understanding in the past regarding the dynamics of the system’s parameter. From this standpoint, it seems interesting to know that the dynamics vary when the parameter moves across the diagram. From our analysis, the parameters and play an important role in pattern formation. Our results show that the pattern formation formed by the bacterial colonies model represents rich spatial dynamics which will be useful for studying the dynamic complexity of bacterial ecosystems.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgements
This research was supported by the Natural Science Foundation of China (61373005 and 11371386), the Doctoral Program of Higher Education of China (20115134110001), the Scientific Research Foundation of CQUT(2012ZD37), and the Fund Project of Zhejiang Provincial Education Department (Y201223449).