Abstract

We investigate the spatial dynamics of a predator-prey system with Allee effect. By using bifurcation analysis, the exact Turing domain is found in the parameters space. Furthermore, we obtain the amplitude equations and determine the stability of different patterns. In Turing space, it is found that predator-prey systems with Allee effect have rich dynamics. Our results indicate that predator mortality plays an important role in the pattern formation of populations. More specifically, as predator mortality rate increases, coexistence of spotted and stripe patterns, stripe patterns, spotted patterns, and spiral wave emerge successively. The results enrich the finding in the spatial predator-prey systems well.

1. Introduction

The Allee effect, named after the ecologist Warder Clyde Allee, has been recognized as an important phenomenon of positive density dependence in low-density population [15]. Allee effect can occur whenever fitness of an individual in a small or sparse population decreases as the population size or density also declines [6, 7]. Since the outstanding work of Allee [1], the Allee effect has been regarded as one of the central and highly important issues in the population and community ecology. And its critical importance has widely been realized in the conservation biology that Allee effect is most likely to increase the extinction risk of low-density populations. As a result, studies on Allee effect have received more and more attention from both mathematicians and ecologists.

Long time series of the density of both prey and predator is needed, so it is difficult to analyse their dynamics. As a result, it may provide useful information by constructing mathematical models to investigate the dynamical behaviors of predator-prey systems. There have been a large group of papers on predator-prey systems with Allee effect [813]. However, these previous works did not take into account the effect of space.

There are also some works done on spatial predator-prey systems with Allee effect [1416]. Petrovskii et al. found that the deterministic system with Allee effect can induce patch invasion [14]. Morozov et al. found that the temporal population oscillations can exhibit chaotic dynamics even when the distribution of the species in the space was regular [15]. Moreover, they found that the chaos accompanied with patch invasion even though the environments were heterogeneous [16]. However, their results were obtained by choosing particular initial conditions. Then, it is natural to ask what kind of patterns can be obtained in predator-prey systems with Allee effect by using other initial conditions. To understand that mechanism well, we will investigate a predator-prey system with Allee effect.

Because of the insightful work of many scientists over recent years, we can make research on pattern selection by using the standard multiple scale analysis [17, 18], in which the control parameters and the derivatives are expanded in terms of a small enough parameter. In the neighborhood of the bifurcation points (Hopf and Turing bifurcation points), the critical amplitudes follow the normal forms, and thus their general forms can be obtained from the methods of symmetry-breaking bifurcations.

The paper is organized as follows. In Section 2, we present a predator-prey system with Allee effect and give Turing region in parameters space. In Section 3, by using multiple scale analysis, we obtain amplitude equations. In Section 4, we show the spatial patterns by a series of numerical simulations. Finally, conclusions and discussions are presented in Section 5.

2. A Predator-Prey System with Allee Effect

We consider the following model of two-dimensional spatiotemporal system [1416, 19]:

where and are densities of prey and predator, respectively, at time and position . The function represents the intrinsic prey growth, represents predation term,is the food utilization coefficient, and are diffusion coefficients, and describes predator mortality.

It is assumed that the predation term is a bilinear form of prey and predator density and predator mortality is a nonlinear function of predator density. As a result, we choose and [20].

When the prey population obeys Allee dynamics, its growth rate can be parameterized as follows [14, 15, 21]: whereis the prey-carrying capacity,is the maximum per capita growth rate, andquantifies the intensity of the Allee effect. If , is a strong Allee effect; if , is a weak Allee effect; if , the Allee effect is absent.

In order to minimize the number of parameters involved in the model system, it is extremely useful to write the system in a nondimensionalized form. Although there is no unique method of doing this, it is often a good idea to relate the variables to some key relevant parameters. Introducing dimensionless variables we obtain the following equations:

where First of all, we need to investigate the dynamics of nonspatial model of systems (4a) and (4b)

Systems (6a) and (6b) have three boundary equilibrium named , , and and two interior equilibriums named and , where

where.

From a biological point of view, we are concerned with the dynamics ofand. The Jacobian matrix corresponding to the equilibrium point is that where

Diffusion-driven instability requires the stable, homogeneous steady state is driven unstable by the interaction of the dynamics and diffusion of the species; and therefore

It is found from direct calculations that is unstable and is stable. Denote .

Following the standard linear analysis of the reaction-diffusion equation [22], we consider a perturbation near the steady state: where , , and . Assume that where is the growth rate of perturbation in time , and represent the amplitudes, and and are the wave number of the solutions.

The characteristic equation of the systems (4a) and (4b) is where

As a result, we have characteristic polynomial: where .

The roots of (15) can be obtained by the following form:

When and , Hopf bifurcation will emerge. Then, we have that the critical value of Hopf bifurcation parameter- equals

When and , , Turing bifurcation will occur. Denote as the critical value ofas Turing instability occurs. Since the expression is complicated, we omit it here.

In Figure 1, we show the two critical lines in the parameter space spanned by and. The equilibria that can be found in the region, marked by (Turing space), are stable with respect to the homogeneous perturbations, but they lose their stability with respect to the perturbations of specific wave numbers . In this region, stationary patterns can be observed. To see the effect of parameter well, we plot in Figure 2 the dispersion relation corresponding to several values of while keeping the other parameters fixed. We see that the available Turing modes shift to higher wave numbers when decreases.

3. Spatial Dynamics of Systems (4a) and (4b)

In the following, we use multiple scale analysis to determine the amplitude equations when . Denote as the controlled parameters. When the controlled parameter is larger than the critical value of Turing point, the solutions of systems (4a) and (4b) can be expanded as with. and the conjugate are the amplitudes associated with the modes and .

Close to onset , one has that

Based on the center manifold near the Turing bifurcation point, it can be concluded that amplitudesatisfies

In order to obtain the amplitude equations, we first need to investigate the linearized form of systems (4a) and (4b) at the equilibrium point . By setting and , we have the following equations:

Close to onset, the solutions of systems (4a) and (4b) can be expanded as series form:

System (19) can be expanded as where is the eigenvector of the linearized operator.

From the standard multiple scale analysis, up to the third order in the perturbations, the spatiotemporal evolution of the amplitudes can be described as

Due to spatial translational symmetry, we have the following equation:

Comparing (25) with (26), one can find that the two equations hold only if . From the center manifold theory, we know that amplitude equations do not include the amplitude with unstable mode. As a result, we have the following equations: where and is a typical relaxation time.

In the following part, we will give the expressions of , , , and . Let Then systems (4a) and (4b) can be written as: where

We need to investigate the dynamical behavior whenis close to, and thus we expandas: where is a small enough parameter. We expand and as the series form of : Linear operator can be expanded as where

Let and is a dependent variable. For the derivation of time, we have that

The solutions of systems (4a) and (4b) have the following form: This expression implies that the bases of the solutions have nothing to do with time and the amplitude is a variable that changes slowly. As a result, one has the following equation:

Substituting the above equations into (29) and expanding (29) according to different orders of , we can obtain three equations as follows: where

We first consider the case of the first order of . Since is the linear operator of the system close to the onset, is the linear combination of the eigenvectors that corresponds to the eigenvalue zero. Since that we have that

As , we can obtain that by assuming .

Let then where and is the amplitude of the mode .

Now, we consider the case of the second order of. Note that

According to the Fredholm solubility condition, the vector function of the right hand of the above equation must be orthogonal with the zero eigenvectors of operator . And the zero eigenvectors of operator are

It can be found from the orthogonality condition that where and represent the coefficients corresponding to in and .

By investigating , one has

It can be obtained from the orthogonality condition that

By using the same methods, one has where

By solving the sets of the linear equations about , , , and , we obtain that where .

For the third order of, we have that where

Using the Fredholm solubility condition, we can obtain where

By using the same methods, we can obtain the other two equations. The amplitudecan be expanded as

As a result, we have

The other two equations can be obtained through the transformation of the subscript of. By calculations, we obtain the expressions of the coefficients of,,, andas follows: where and .

By using substitutions, we have where . In order to see the relationships between different parameters, we give the values of coefficients for different parameter sets in Table 1.

The dynamical systems (4a) and (4b) possess five kinds of solutions [23] as follows.(1) The stationary state , given by is stable for and unstable for .(2)Stripe patterns , given by are stable for, and unstable for .(3)Hexagon patterns are given by with or , and exist when The solution is stable only for and is always unstable.(4)The mixed states are given by with . They exist whenand are always unstable.

4. Spatial Pattern of Systems (4a) and (4b)

In this section, we perform extensive numerical simulations of the spatially extended systems (4a) and (4b) in two-dimensional spaces. All our numerical simulations employ the zero-flux boundary conditions with a system size of 200 × 200. The space step is , and the time step is .

In Figure 3, we show the spatial pattern of prey population at different time. In the parameter set, , , and , we find that , which means that there is coexistence of spotted and stripe patterns. As shown in this figure, our theoretical results are consistent with the numerical results.

By setting , , and , one can obtain that . In Figure 4, we show the spatial pattern of prey population whenequals 0, 50, 100, 200, 500, and 1000. At the initial time, the prey population shows patched invasion. As time increases, stripe pattern appears and the structure does not change a lot. While keeping other parameters fixed and increasing, we find that stripe pattern will occupy the whole space. However, some stripe patterns connect with each other and cause the emergence of spotted patterns which are shown in Figure 5.

Figure 6 shows the evolution of the spatial pattern of prey population at , 100, 300, 500, 1000, and 2000 iterations, with small random perturbation of the stationary solution of the spatially homogeneous systems (4a) and (4b). The corresponding parameters values are , , and . By the amplitude equations, we can conclude that there are spotted patterns of prey population for this parameter set. In this case, one can see that for the systems (4a) and (4b), the random initial distribution leads to the formation of an irregular transient pattern in the domain. After these forms, it grows slightly and spotted patterns emerge. When the time is large enough, the spotted patterns prevail over the two-dimensional space. As time further increases, the pattern structures of the prey population do not undergo any further changes.

5. Conclusion and Discussion

Allee effect has been paid much attention due to its strong potential impact on population dynamics [24]. In this paper, we investigated the pattern dynamics of a spatial predator-prey systems with Allee effect. Based on the bifurcation analysis, exact Turing pattern region is obtained. By using amplitude equations, the Turing pattern selection of the predator-prey system is well presented. It is found that the predator-prey systems with Allee effect have rich spatial dynamics by performing a series of numerical simulations.

It should be noted that our results were obtained under the assumption that predation is modeled by the bilinear function of the prey and predator densities. However, this function has limitations to describe many realistic phenomena in the biology. By numerical simulations, we find that the system exhibits similar behaviour when the functional response is of other types, such as Holling-II and Holling-III forms.

To compare the spatial dynamics for different parameters, we give the spatial patterns of populationwhen the parameter values are out of the domain of Turing space. For this parameter set, systems (4a) and (4b) have Hopf bifurcation, and spiral waves occupy the whole domain instead of stationary patterns, which is shown in Figure 7. The stability of spiral wave can be done by using the spectrum theory analysis [25, 26]. In the further study, we will use the spectrum theory to show the stability of spiral wave.

In [15], they found that a spatial predator-prey model with Allee effect and linear death rate could increase the system’s complexity and enhance chaos in population dynamics. However, in this paper, we showed that a spatial population model with Allee effect and nonlinear death rate can induce stationary patterns, which is different from the previous results.

From a biological point of view, our results show that predator mortality plays an important role in the spatial invasions of populations. More specifically, low predator mortality will induce stationary patterns (cf. Figures 36), and high predator mortality corresponds to travelling patterns (cf. Figure 7). When the populations exhibit wave distribution in space, the dynamics of populations may be accompanied with chaotic properties [27, 28]. If the chaotic behavior occurs, it may lead to the extinction of the population, or the population may be out of control [29, 30]. In that case, we need to find out the best way to control the chaos or change the chaotic behavior.

Acknowledgments

This research was partially supported by the National Natural Science Foundation of China under Grant nos. 11301490, 11301491, 11331009, 11147015, 11171314, 11305043, and 11105024; Natural Science Foundation of Shan’xi Province Grant nos. 2012021002-1 and 2012011002-2, the Opening Foundation of Institute of Information Economy, Hangzhou Normal University, Grant no. PD12001003002003; and the specialized research fund for the doctoral program of higher education (preferential development) Grant no. 20121420130001.