|Home | About | Journals | Submit | Contact Us | Français|
Spatial patterns are a subfield of spatial ecology, and these patterns modify the temporal dynamics and stability properties of population densities at a range of spatial scales. Localized ecological interactions can generate striking large-scale spatial patterns in ecosystems through spatial self-organization. Possible mechanisms include oscillating consumer–resource interactions, localized disturbance–recovery processes, and scale-dependent feedback. However, in this paper, our main aim is to study the effect of tide on the pattern formation of a spatial plant-wrack model. We discuss the changes of the wavelength, wave speed, and the conditions of the spatial pattern formation, according to the dispersion relation formula. Both the mathematical analysis and numerical simulations reveal that the tide has great influence on the spatial pattern. More specifically, typical traveling spatial patterns can be obtained. Our obtained results are consistent with the previous observation that wracks exhibit traveling patterns, which is useful to help us better understand the dynamics of the real ecosystems.
The study of complexity in ecosystems has been recognized as one of the frontiers of ecology . Many ecosystems exhibit remarkable complexity in both structure and dynamics, which has repeatedly frustrated attempts to predict the response of ecosystems to human interference . Complexity within ecosystems, however, has also been recognized as the basis for the adaptive capacity of ecosystems to future changes on both regional and global scales [3–5]. A key property of complex systems is emergence, where properties of the ecosystem at large spatial scales result from feedback interactions between components occurring at smaller scales [5, 6]. These emergent properties result from the capacity of systems to self-organize [5, 7–31].
More recently, many ecologists have paid more and more attention to the experimental investigation of regular spatial patterning in Carex stricta. Carex stricta, the tussock sedge, is a species with widespread distribution in freshwater marshes of North America. Spatial dispersals of vegetation (through tillers) and wrack (resulting from dead plant leaves dropping to the soil surface and movement by the tides) are modeled using a diffusion approximation. The model, which describes the interaction of the plant and wrack, is as follows :
where P is plant the biomass, W is the wrack biomass, F(P) is a function describing the positive effect of plant biomass on its own growth, s is the specific rate of plant senescence, I(P, W) is a function describing the inhibiting effect of wrack on plant growth as a function of plant and the wrack biomass, b is the decay rate of wrack, and dP and dW are diffusion constants describing lateral movement of plants and wrack. Here, is the usual Laplacian operator in two-dimensional space.
However, another type of dynamics will emerge when the mobility of the tides exhibits a correlated motion towards a certain direction, which is the cyclic rising and falling of Earth’s ocean surface caused by the tidal forces of the Moon and the Sun acting on the oceans. Tides are important natural forces that cause changes in the depth of marine and estuarine water bodies and produce oscillating currents (see http://en.wikipedia.org/wiki/Tide for more detail), which may have great impact on the spatial patterns of plant and wrack population dynamics.
In the present paper, we want to see the effect of tidal current on the pattern formation of both the plant and wrack. To exclude the Turing pattern, we set dP=dW=1. Additionally, system (1) will be changed as:
Advection by tidal current is represented by the gradient operator multiplied by the advection constant δ . For this simple model, we assume that tide only has an effect on the wrack for the reason that the wrack is much lighter than the plants in weight.
In this model, we focus on the following specific type with an indirect facilitation of growth by the root mound by lowering of inhibition by the wrack, which leads to:
where K/(P+K) is added to the inhibition term, lowering inhibition as P increases, K is the level of plant biomass where inhibition is lowered by half, and a is an inhibition coefficient . As a result, kinetics of the system are
In order to provide guidelines on the dynamics of the full reaction–diffusion–advection system, it is important to consider the local dynamics of the system . It is, naturally, the dynamics in the biologically meaningful region P≥0, W≥0 that are of interest. By considering the nullclines f=0, g=0, and the intersection of these curves in phase space, linear stability analysis reveals that we have a stationary point at (0, 0) for the kinetics. The stability is determined by s, i.e., it is stable if s>1 and unstable if s<1 (which is a saddle). Furthermore, we want to find the positive steady state of the non-spatial model, (P *,W *), which corresponds to the coexistence of plant and wrack. By direct calculation, we find that, if s<1, there exists a unique positive steady state, and
where M=b2−2absK+2Kb2− 2sb2+a2s2K 2+2asbK 2+2as2Kb+K 2b2−2Kb2s+s2b2. Additionally, if s>1, there is no positive steady state.
From the biological point of view, we are interested in studying the stability behavior of the positive steady state (P *,W *). The Jacobian corresponding to this equilibrium point is
where the elements are the partial derivatives of the reaction kinetics evaluated at the stationary state (P*,W*). For the sake of brevity, we do not write down the expression for the elements.
It is easy to see that model (2) has the form of a generalized reaction–diffusion–advection model. To consider pattern formation for system (2), we should look for the dispersion relation, which gives the condition for pattern formation. Following , we will assume periodic (or zero-flux) boundary conditions and prescribed initial conditions. A symmetry breaking occurs when a homogeneous steady state solution of system (2) is linearly stable to perturbations in the absence of the diffusion and advection terms, but linearly unstable to small spatial perturbations in the presence of diffusion and advection. It is a simple matter to find the exact algebraic formulation of the dispersion relation when these equations are subject to certain boundary conditions and stability properties. Here, we consider an approach based on Laplace transform methods, which is convenient for finding symmetry breaking (Turing instability) conditions in the case when the systems include diffusion and advection. The Laplace transform methods are extensively used in the reaction–advection–diffusion equations including the system with fractional diffusion and fractional reaction[35–38].
We make the following substitution:
into the kinetic equations and assume . Then, in the linear approximation, we have
The initial conditions are assumed as
where the functions and decay rapidly for . Following the standard approach, let us now perform a Laplace transformation of the linearized equations over the two independent variables and t. For , we use the so-called two-sided version of the transformation. The relations for the forward and backward transforms are
where s and q are complex variables. Additionally, s is the Laplace transform variable and q is the Fourier transform variable. That is to say that q=ik, or (ik,il), corresponding to one- and two-dimensional space, and the wave numbers k and l are real numbers. In formula (10), for the backward transformation, the integration contour in the q plane is the imaginary axis. In the s plane, the contour is parallel to the imaginary axis and located to the right of all singularities of the integrand.
After this transformation, the kinetic equations read
where H(q) and Z(q) are the transforms of and . To reveal the presence of an instability and disclose its character, it is sufficient to consider one variable. The temporal growth of the perturbations can now be found by inverting the Laplace transforms, which follows directly after factorizing the denominator. In the presence of advection, the instability is either absolute or convective [39–45]. If the advection coefficient is small, the diffusion carries the perturbation upstream faster than the advection sweeps. As a result, the perturbation invades the whole space. This situation is referred to as an absolute instability. If the advection coefficient is large, the advection acts faster than the diffusions do. Therefore, the perturbation asymptotically decays towards zero at every fixed point. This is the case of the convective instability. The transition between these two cases happens at the critical advection coefficient that can be determined from the set of equations [46, 47, 37] (see also references therein):
Then, we obtain the linear stability of this state, which is described by the dispersion relation
For the sake of convenience, we first consider to the one-dimensional case. By setting q=ik, we have that
Straightforward manipulation of (17) yields
Under these conditions, the spatial mode q (in one- or two-dimensional space) is unstable, and thus, system (2) grows into a pattern, that is, Re(λ)>0. For the non-spatial model, the equilibria should be locally stable, i.e., a11+a22<0 and a11a22−a12a21>0. By direct calculation, we find that the condition for emergent pattern formation is
The dynamical behavior of the spatial model cannot be studied by using analytical methods or normal forms. Thus, we have to perform simulations by computer. In this section, we give the numerical simulations of the partial differential equations (PDEs) (2) in the one- and two-dimensional space. We have taken the approach of discretizing the PDEs in space, applying upwinding to the convective term and Euler integration of the finite difference equations, which means approximating the value of the concentration at the next time step based on the rate of change of the concentration at the previous time step. Here, the time stepsize Δt=0.05 and space stepsize (lattice constant) Δh=1. The scale of the space and time is averaged for the Euler method. We used a rectangular spatial grid consisting of 50×50 cells, with a unidirectional water flow in the y axis direction.
All our numerical simulations employ periodic (zero-flux) boundary conditions. As far as we know, there are two different initial conditions used for analysis of the spatial extended systems [48–51]. The first one, which is more general from the biological point of view, is random spatial distribution of the species. The other one is a special choice, which is made by taking the species community in a horizontal layer as decreasing gradually and the vertical distribution of species as homogeneous. In this paper, we choose the first method. The initial density distribution corresponds to random perturbations around the trivial stationary state (P *, W *) in model (2), with a variance significantly lower than the amplitude of the final patterns, which seems to be more general from the biological point of view.
In most cases, we used the explicit scheme. In order to avoid numerical artifacts, we checked the sensitivity of the results to the choice of the time and space steps, and their values have been chosen as sufficiently small. Also, some of the results were reproduced by means of using a more advanced alternate directions scheme. Both numerical schemes are standard; hence, we do not describe them here; details and particulars can be found in .
It must be mentioned that, since the system under study depends on a relatively large number of parameters (five for the system (2)), its detailed numerical investigation in the whole parameter space is virtually impossible. Instead, we choose one controlling parameter and consider the changes in the pattern formation subject to its variation. In this paper, we choose the advection coefficient δ as the controlling factor and keep all other parameters fixed. The choice of δ as the controlling factor is justified by the results of field observations, showing that intensity of advection affects the pattern formation. Also, this choice is consistent with our objective in this paper, i.e., finding the effect of the advection on the regular pattern formation in freshwater systems.
Figure 1 shows the evolution of the spatial pattern of both plant and wrack at t=385, 490, and 705, with small random perturbation of the stationary solution P * and W * of the spatially homogeneous systems when the parameters can ensure Re(λ)>0. In this case, one can see that there are travelling spatial patterns with δ=2.5. In Fig. 2, we show the typical spatial pattern formation of both the plant and wrack with δ=4, which is similar to Fig. 1. However, the wavelength in Fig. 2 is larger than that in Fig. 1. That is to say, the single parameter of the system (3), namely, δ, can lead to dramatic changes in the qualitative dynamics of the solutions. The parameter δ is the coefficient of the tidal current. From the biological point of view, how fast the tide diffuses can have a profound effect on the fate of the system.
In order to clearly see the evolution of the pattern, we show spatial plots in Fig. 3. It is well known that the basic arguments in spiral stability analysis can be carried out by reducing the system to one-dimensional space [53–56]. The method of generating space-time plots is to let Y be a constant from each pattern snapshot (here, Y=25, the corresponding to the center line), and to stack these lines in temporal order. As seen from Fig. 3, a typical traveling spatial pattern is shown, which is consistent with Fig. 2.
Our numerical analysis suggests that system (2) does not have any stationary patterns, with the homogeneous steady states being the only stationary solution. Instead, the patterns move at a constant speed, which means that the patterns are a periodic traveling wave solution. For the sake of convenience, we first investigate patterns in the one-dimensional space, which have the mathematical form P (X,t)=P (ξ) and W(X,t)=W(ξ), where ξ=X+ct and c is the migration speed. The results are also valid in two-dimensional space. Substituting these solution forms into (3) gives the ordinary differential equations
where ξ is now treated as a time- or space-like variable. This second-order ordinary differential equation predicts spatial periodic wave forms. Equation 25 can then be written as a system of four first-order ordinary differential equations, which are as follows:
It is easy to obtain the steady-state solutions to this system, which are of the form (P,u,W,v)=(P *,0,W *,0), where (P *,W *) is the spatially uniform steady-state solution of the nonspatial equation (2).
In the following, we track the family of traveling wave solutions of (26) by using the software package AUTO , which is a software tool partly designed for continuation and bifurcation problems in ordinary differential equations. We can use it to locate a period of that periodic traveling wave that varies as we gradually change the equation parameters. The first calculation is to locate the wave speed at which there is a Hopf bifurcation in (26), and the starting point is the homogeneous steady state (26).
We plot the typical bifurcation diagrams in Fig. 4, which shows that the changes of stability occur via Hopf bifurcations with respect to the speed c, from which arises a single branch that links the two Hopf points. When the speed is above the critical value of speed, the stability of the equilibrium will be lost and a periodic orbit is born
The bifurcation diagrams of solution amplitudes are represented in Fig. 5, where we use the standard definition of L2-norm
where L is the spatial period and Vk(ξ) are four variables of four first-order ordinary differential equations about (26). Figure 5 shows that the periodic orbits are formed by bifurcation with an intermediate range of speeds.
From the practical point of view, one may also be concerned with the effect of other parameters, such as s,a,K, and b. Using AUTO, it is possible to track the locus of the Hopf bifurcation lines, which we show in Fig. 6.
In this paper, we investigate spatial pattern formation in a plant-wrack model combined with tide effects. By obtaining the dispersion relation by numerical calculation, we discuss the effects of tide on the dispersion relation including the wavelength, and conditions of spatial pattern formation by using linear analysis. Furthermore, we perform nonlinear analysis by using the bifurcation package AUTO to confirm our numerical prediction resulting from the linear analysis.
In , they investigate the plant-wrack model only with the diffusion term, and obtain typical Turing patterns. However, when the former system is combined with the tidal current term, typical travelling patterns emerge. In the real world, wrack exhibits travelling patterns [58, 59]. That is to say, our obtained results are consistent with the real world.
Large-scale spatial patterns widely exist inecosystems. Despite abundant theoretical literature, studies revealing spatial self-organization in real ecosystems are very few [23–27]. Further research is necessary to understand and predict regular pattern formation in ecosystems .
This work is supported by the National Natural Science Foundation of China under Grant No. 60771026, Program for New Century Excellent Talents in University (NCET050271), the Special Scientific Research Foundation for the Subjects of Doctors in University (20060110005), Graduate Students’ Excellent Innovative Item of Shanxi Province No. 20081018, and the US National Science Foundation Bio-Complexity Program and the University of California Agricultural Experiment Station.
Open Access This article is distributed under the terms of the Creative Commons Attribution Noncommercial License which permits any noncommercial use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.
Gui-Quan Sun, Email: nc.moc.oohay@nusnauqg.
Li Li, Email: moc.361@311138ilil.
Zhen Jin, Email: ten.362@nhznij.
Bai-Lian Li, Email: email@example.com.