Abstract
We investigate a spatial ratio-dependent predator-prey system with linear harvesting rate. By using linear stability and bifurcation analysis, we obtain the conditions for Hopf and Turing bifurcation in the spatial domain. In addition, we classify spatial pattern formations of the system by making use of numerical simulations. In fact, the numerical simulations unveil that the typical Turing patterns, such as spotted, spot-stripelike mixtures and stripelike patterns, can be observed even if the system has the linear harvesting rate. In order to analyze these patterns via the spatial frequency, the discrete Fourier transform is used. Moreover, we discuss that the linear harvesting system is more realistic than a predator-prey system with constant harvesting. Our results disclose that the spatially extended system with linear harvesting rate has more complex dynamic patterns in the two-dimensional space. It may help to understand the effects of harvesting on species in the real world.
1. Introduction
Since a classical Holling-Tanner type predator-prey system, so-called a ratio-dependent predator-prey system, has received great attention among mathematicians and mathematical biologists (see [1–4]), we will focus our concern on the following system of nonlinear coupled ordinary differential equations: where and stand for prey and predator, respectively. All parameters are positive constants. Especially, the parameters , , , and stand for the prey intrinsic growth rate, the carrying capacity of prey, the capturing rate, and the half-saturation constant, respectively. The predator grows logically with the growth rate and conversion efficiency . We adopt a ratio-dependent functional response in system (1) which roughly the per capita predator growth rate should be a function of the ratio of prey to predator abundant. One of the main reasons to consider such response function is that numerous laboratory experiments and observations showed that functional and numerical responses over typical ecological timescales ought to depend on the densities of both prey and predator, especially when predators must search for food (see [1, 5–8]).
Following Kuang and Beretta in [3] and Xiao et al. in [7, 8], with the next scaling we can arrive at the following equations containing dimensionless quantities by omitting bar notations;
Although such kinds of predator-prey systems play an important role in studying about dynamics of prey or predator populations, we cannot infer any information about spatial distributions of prey or predator from these ordinary differential systems. Thus, it is needed to investigate spatial systems. In this context, recently, many authors in [9–15] have investigated spatiotemporal pattern formations and bifurcation analysis of predator-prey systems with reaction diffusion. Thus, firstly, we will take into account the following a ratio dependent predator-prey system with reaction diffusion: where are the diffusion coefficients of prey and predator, respectively, and is the usual Laplacian operator in two-dimensional space .
In fact, Banerjee [10] and Wang et al. in [13] have studied the spatiotemporal complexity of a ratio-dependent predator-prey system similar to system (4).
From the point of view of human needs, the exploitation of biological resources and the harvest of population are commonly practiced in fishery, forestry, and wildlife management. Concerning the conservation for the long-term benefits of humanity, there is a wide range of interest in the use of bioeconomics modeling to gain insight in the scientific management of renewable resources like fisheries and forestries (see [7, 8]). Thus, it is necessary to investigate the effect of the harvesting of species in some models. In fact, many researchers have studied complex bifurcation phenomena and dynamics of the models with nonzero constant rate of harvesting of either species or both species simultaneously (see [7, 8, 16–18]). However, in general, human needs are not immutable for a long time. Actually, human needs could be changed according to the amount of biological resources. For instant, let us think about an outbreak of algal bloom caused by phytoplankton or zooplankton (see [12, 19]). In this case, harvesting the populations of phytoplankton or zooplankton is one of simple ways to control the algal bloom. It is an undoubted fact that a nonconstant harvesting method considering prey or predator populations is more efficient and economic than the constant harvesting method. For the reason, it is needed to take into account nonconstant harvesting rate. Such harvesting activity has an effect on both species, prey and predator, directly. Thus, with the ideas discussed above, in the paper, we consider the following ratio-dependent predator-prey system with reaction-diffusion and linear harvesting rates and : where are nonnegative constants.
However, since , the dynamics of system (5) are the same as the following system: where the boundary conditions are taken as on , which means that no external input is imposed from outside. Here, is the outward unit normal vector of the boundary .
Thus, from now on, we will focus on studying spatial dynamics of system (6) instead of system (5).
The main purpose of this paper is to investigate spatial dynamics of system (6). In Section 2, we investigate bifurcation phenomena and Hopf and Turing bifurcations by making use of linear stability theory and bifurcation analysis. In addition, we obtain the conditions of Hopf and Turing bifurcation in a spatial domain. In Section 3, using numerical simulations for system (6), we show that the spatial system (6) has various spatiotemporal patterns such as spotted and stripelike patterns. Also, we obtain the spatial frequency via discrete Fourier transforms of the spotted and stripelike patterns to analyze these patterns in terms of periodic components. Finally, we give conclusions and discussions in Section 4. Especially, we discuss that the linear harvesting rate is more realistic than the constant harvesting by comparing the above linear harvesting system with the constant harvesting system studied by Zhang et al. in [15].
2. Bifurcation Analysis
In the beginning, we will think about the spatially homogeneous system of the diffusive system (6) to look into spatial pattern formations of system (6). It is fundamental to mention that the nonspatial system of system (6) has at most three nonnegative equilibria, which correspond to spatially homogeneous equilibria of system (6), in the first quadrant , where and . Moreover, the point lies in the interior of the first quadrant if and .
It follows from [20] that the point is a saddle point and the point is either a saddle point when or a stable node when for the nonspatial system of system (6).
From the biological point of view, it is important to mainly take into account the dynamics of system (6) at the nontrivial stationary state , which corresponds to the coexistence of prey and predator. Throughout the paper, we assume that parameters in system (6) satisfy the conditions Under condition (8), the equilibria and are both saddle points.
In order to deal with the stability and bifurcation analysis of the spatiotemporal system (6), we linearize the dynamical system (6) around the spatially homogeneous fixed point for small space and time fluctuations. To accomplish this purpose, we assume that where and , and and are the wavenumber vector and frequency, respectively. Then we can obtain the corresponding characteristic equation , where and is the diffusion matrix and the matrix is given by where. From an elementary calculation, we can get Then the solutions of the characteristic equation yield the dispersion relation where
It is well known that reaction-diffusion systems have led to the characterization of two basic types of symmetry-breaking bifurcations, Hopf and Turing bifurcations, which are responsible for the emergence of spatiotemporal patterns (see [9–15, 21–27]).
The space-independent Hopf bifurcation breaks temporal symmetry of a system and gives rise to oscillations that are uniform in space and periodic in time (see [9, 13, 15, 20]). This bifurcation takes place only when the diffusion vanishes. Thus, mathematically speaking, the Hopf bifurcation occurs if the conditions hold. From simple calculation, it is easy to know that the critical value of Hopf bifurcation parameter is equal to From now on, for the existence of the critical value of Hopf bifurcation, it is assumed that the parameters are satisfied with the following inequality: Thus the nontrivial stationary state of system (6) is an unstable equilibrium if . Furthermore, as mentioned above, the bifurcation breaks the temporal symmetry of the system, which causes oscillations that give rise to uniform in space and periodic in time with the frequency and the corresponding wavelength is
On the other hand, Turing bifurcation is sometimes called Turing instability or diffusion-driven instability. It is well known that Turing bifurcation breaks spatial symmetry, leading to the formation of patterns that are stationary in time and oscillatory in space. Moreover, since this instability occurs only if the activator () diffuses more slowly than inhibitor (), the condition is assumed.
In fact, Turing bifurcation occurs when the positive steady state of system (6) is stable for the nonspatial ordinary differential system of system (6) and it is unstable for the diffusive system (6) (see [9, 10, 12, 13, 15, 21–23, 28]). The conditions for the steady state to be stable are given by Since , is always negative and hence the two roots of the characteristic equation cannot be positive at the same time, therefore Turing instability can happen only if because the stationary state is unstable to spatial perturbations for . To hold the condition , the coefficient of must be negative because of ; that is, However, condition (20) is necessary, not a sufficient condition. Thus we must require that the minimum value of is below zero. Actually, has the minimum value at which is called the wavenumber. Therefore the value of at must be negative; that is, In order to get the threshold of Turing bifurcation, let the inequality in (22) be the equality. Then we have the critical value of the bifurcation parameter for the Turing instability as follows:
At the Turing threshold , the spatial symmetry of the system is broken and the patterns are stationary in time and oscillatory in space with the wavelength
In order to discuss the bifurcations represented by the above formulae in the parameter space spanned by the parameters and , we will observe a bifurcation diagram of system (6) by letting , , , and in system (6). The bifurcation diagram shown in Figure 1 can be obtained using the linear stability analysis studied above. In this case, we consider the range to guarantee the existence of Turing instability. As seen in Figure 1, the bifurcation parametric space is separated into four distinct domains by the Hopf and the Turing bifurcation lines. In domain I located below two bifurcation lines, the steady state is the only stable solution of system (6). On the other hand, domains II and III are the regions of the pure Turing instability and the pure Hopf instability, respectively. Hopf and Turing instabilities take place in domain IV simultaneously.
3. Pattern Formation Analysis
In this section, we will investigate spatiotemporal pattern formations of the spatially extended system (6) via computational simulation in two-dimensional space. For this, let , , , and . In order to solve the system (6) numerically, we employ a finite difference scheme for the spatial derivatives, an explicit Euler method for the time integration with the time step , and the discrete domain of lattice sites with the space step . These numerical schemes for the diffusion equation (6) give stable solutions as long as the equalities are satisfied (see [29]).
It is well known that spatiotemporal dynamics of a diffusion-reaction system depend on the choice of initial conditions (see [12, 13, 28]). However, from the biological point of view, it is natural to take an initial condition as a small amplitude random perturbation around the steady state . We run the simulations until the numerical solutions either reach a stationary state or show a behavior that does not seem to change its characteristic anymore.
In the numerical simulation, we restrict our analysis of pattern formations to the distribution of prey since we have figured out that the distributions of predator and prey are always of the same type. Nevertheless, through the numerical simulation different types of spatiotemporal dynamics of the distribution of prey are observed.
Now, fixing the parameter value , the bifurcation thresholds, for the Hopf bifurcation, and for the Turing bifurcation can be obtained using (15) and (23), respectively (see Figure 1). In addition, due to Turing instability, three typical Turing patterns such as spotted, spot-stripe mixture, and stripelike patterns of prey () in system (6) can be observed accordingly as the harvesting rate varies increasingly. Snapshots that represent such phenomena are illustrated in Figures 2, 3, and 4 when taking the harvesting rate and , respectively.
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
To analyze these patterns by using the spatial frequency, we make use of the discrete Fourier transform. It is known that the discrete Fourier transform is one of useful mathematical tools to analyze a signal or an image by decomposing it into different periodic components. Especially, it has been widely used for spatial patterns (see [13, 26, 30]). Thus, we apply the spatial Fourier transform to the Turing patterns in Figures 2, 3, and 4 to get the spatial frequency as in Figure 5. The spatial frequency and direction of any component in the power spectrum are given in the length and direction, respectively, of a vector from the origin to the point on the circle (see [13, 26, 30]). We can infer from Figure 5 that Figures 2 and 3 have similar spatial frequency in the length of the space unit and presence of one mode with different wavelengths while Figure 4 has two modes with different wavelengths. Since the wavelength is defined by the inverse of the frequency, the small and large circles in Figure 5(c) have long and short wavelength, respectively.
If we select and in Domain III, which is the region of the pure Hopf instability, we can investigate influence of the Hopf instability phenomena on system (6). Of course, it is not easy to see the influence of the instability directly through snapshots as shown in Figure 6. However, an oscillatory pattern in time can be observed. From ecological point of view, this bifurcation gives rise to the coexistence of two species, prey and predator, in the long run. In order to investigate this phenomenon more clearly, it is needed to see a local phase portrait of system (6). For this, let us fix a spatial point ; then the local phase portrait of system (6) for this point can be obtained, exhibited in Figure 7. Therefore, we can assert from Figure 7 that the stable limit cycle occurs by virtue of the Hopf instability. Also, we can estimate numerically the frequency of the periodic oscillations in time and the corresponding wavelength . In fact, it is shown from (17) and (18) that the frequency of the periodic oscillations in time is and the corresponding wavelength is .
(a)
(b)
(c)
4. Conclusions and Discussions
In this study, we have investigated pattern formations of a ratio-dependent predator-prey system with linear harvesting rate. It has been shown from the bifurcation analysis that the system has Turing and Hopf bifurcations. For some fixed parameters, the numerical simulations illustrate that the spatial domain is divided into four regions by two bifurcation lines, Hopf and Turing bifurcation lines, as shown in Figure 1. In domain II of Figure 1, as the linear harvesting rate increases typical Turing patterns such as spotted, spot-stripelike mixtures, and stripelike patterns are emerging and in domain III, a stable limit cycle emerges due to the Hopf bifurcation. In addition, applying the spatial Fourier transform to such Turing patterns, we classify these patterns in terms of the spatial frequency.
As mentioned above, we think about spatiotemporal dynamics of a predator-prey system with linear harvesting rate while Zhang et al. in [15] have dealt with a predator-prey system with constant harvesting rate. Especially, they reported that some spatial patterns including isolated groups such as stripelike and patch-like could be shown when the system has constant harvesting of the predator (see system (26)) whereas typical Turing patterns could emerge when the system has constant harvesting of the prey (see system (25)). Also Fasani and Rinaldi in [31] showed that Turing instability could not arise from prey harvesting: In this paper, we also dealt with harvesting of the prey, but the harvesting rate depends linearly on the density of the prey. We illustrated that this system has also typical Turing patterns such as spotted, stripelike, or spot-stripelike mixtures patterns which are similar phenomena to those when Wang et al. in [13] probed into a ratio-dependent predator-prey system without any harvesting rate. Moreover, we can infer from Figures 2, 3, and 4 that the presence of linear harvesting rate for preys has dramatical effect on changing the spatial patterns from spotted patterns to the stripes pattern or coexistence of spotted pattern and stripes pattern of the prey. Biologically we have a conclusion that if the prey is harvested at a rate of its population, it has an opportunity to be centralized or assembled. From this, we could surmise that the prey has a tendency to take group defence against the harvest of it and the lower the harvesting rate is, the stronger the group defence of the prey is. Thus, we can deduce that the prey population is decreasing accordingly as the harvesting rate of the prey is increasing. On the contrary, in case of constant harvesting studied by Zhang et al. in [15], as the constant harvesting rate of the prey is increasing, the prey population is not decreasing but rather increasing (See Figure 8.) Such phenomena are out of range of general thinking that the harvesting species brings about reducing of its population.
(a)
(b)
(c)
Therefore if one can monitor population of prey and predator, it is clear that the harvesting control method of prey or predator considering their populations is more realistic and reasonable than the constant harvesting control method. Thus, this research will be useful for studying the dynamic complexity of ecosystems or physical system.
Conflict of Interests
The author declares that there is no conflict of interests regarding the publication of this paper.
Acknowledgment
The author greatly appreciates the efforts of the anonymous reviewers and their helpful comments and useful remarks, patently improving this work. This work was supported by research grants from the Catholic University of Daegu in 2012.