International Journal of Mathematics and Mathematical Sciences
Volume 2017 (2017), Article ID 8273430, 9 pages
https://doi.org/10.1155/2017/8273430

Research Article

Stability Analysis of a Fractional Order Modified Leslie-Gower Model with Additive Allee Effect

Agus Suryanto, Isnani Darti, and Syaiful Anam

Department of Mathematics, Brawijaya University, Jl. Veteran, Malang 65145, Indonesia

Correspondence should be addressed to Agus Suryanto; suryanto@ub.ac.id

Received 1 March 2017; Accepted 3 May 2017; Published 28 May 2017

Academic Editor: Shyam L. Kalla

Copyright © 2017 Agus Suryanto et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

We analyze the dynamics of a fractional order modified Leslie-Gower model with Beddington-DeAngelis functional response and additive Allee effect by means of local stability. In this respect, all possible equilibria and their existence conditions are determined and their stability properties are established. We also construct nonstandard numerical schemes based on Grünwald-Letnikov approximation. The constructed scheme is explicit and maintains the positivity of solutions. Using this scheme, we perform some numerical simulations to illustrate the dynamical behavior of the model. It is noticed that the nonstandard Grünwald-Letnikov scheme preserves the dynamical properties of the continuous model, while the classical scheme may fail to maintain those dynamical properties.

1. Introduction

The dynamical interaction of predator and prey is one of important subjects in ecological science. In recent years, one of the most important species interactions is predator-prey model [1]. One of well-known mathematical models which describe the dynamics of prey-predator interaction is the modified Leslie-Gower model proposed by Aziz-Alaoui and Okiye [2]. In this model, the growth rate of predator is in the form of logistics-type where its carrying capacity is proportional to the prey number and environment protection for predator. One of important parameters describing the prey-predator interaction is the functional response which describes the predator’s rate of prey consumption per capita. Aziz-Alaoui and Okiye [2] and Yu [3] have considered a modified Leslie-Gower model with Holling type II functional response, while Yu [4] considered the same model but with Beddington-DeAngelis functional response. In the normalized variables, the modified Leslie-Gower equation with Beddington-DeAngelis functional response can be written as

where and denote population densities of prey and predator at time , respectively. The parameters , , , , , and are positive constants.

One of other factors that influence the interaction of predator and prey is Allee effect, referring to a decrease in per capita fertility rate at low population densities. Allee effect may occur under several mechanisms, such as difficulties in finding mates when population density is low or social dysfunction at small population sizes. When such a mechanism operates, the per capita fertility rate of the species increases with density; that is, positive interaction among species occurs [5–8]. Recently Indrajaya et al. [9] investigate a modified Leslie-Gower equation with Beddington-DeAngelis functional response and additive (both weak and strong) Allee effect on prey

with initial conditions and . The criteria of the Allee effect are as follows [7, 8]:

(i)If , then the Allee effect is weak.(ii)If , then the Allee effect is strong.

It is shown in system (2) that the growth rates of both prey and predator depend only on the current state. In many situations, the growth rate is also dependent on the history of variable or its memory. With the rapid development of fractional calculus, fractional differential equations have been implemented in various fields including biological system. This is due the fact that fractional differential equations are naturally related to the real life phenomena with memory which exists in most of biological system [10–16]. To describe such memory effect, we first recall the definition of fractional integral operator as well as fractional differential operator.

Definition 1 (see [17]). The Riemann-Liouville -order fractional integral operator of any function , is defined by

where is the Euler Gamma function.

Definition 2 (see [17]). Let be an integer which satisfies . The Riemann-Liouville -order fractional derivative of function is defined as

where is the common -order derivative.

The Riemann-Liouville fractional derivative is historically the first concept of fractional derivative and theoretically well established. However, in the case of Riemann-Liouville fractional differential equation, the initial value is usually given in the form of fractional derivative, which is not practical. Consequently, one applies the Caputo fractional derivative which is defined as follows.

Definition 3 (see [17]). The Caputo fractional differential operator of order , with , , is defined by

For simplicity, the Caputo fractional derivative of function of order is denoted by .

From Definition 3, we see the -order fractional derivative at time t is not defined locally; it relies on the total effects of the commonly used m-order integer derivative on the interval . So it can be used to describe the variation of a system in which the instantaneous change rate depends on the past state, which is called the ‘‘memory effect’’ in a visualized manner [18].

In this paper we reconsider system (2). By assuming that the growth rates of both prey and predator at time do not only depend instantaneously on the current state but also depend on the past state, we replace the first order derivatives in system (2) with the fractional order Caputo type derivatives:

where . Hence we have a system of fractional differential equation. In the following we discuss the dynamical properties of system (6). To study the stability of equilibrium points, we apply the following stability theorem.

Theorem 4 (see [19]). Consider the following autonomous nonlinear fractional order system:

The equilibrium points of the above system are solutions to the equation . An equilibrium point is locally asymptotically stable if all eigenvalues of the Jacobian matrix evaluated at equilibrium satisfy .

2. Equilibria and Their Stability

Based on Theorem 4, equilibria of model (6) can be determined by solving the following system:

System (8) has been solved by Indrajaya et al. [9], and the obtained equilibria are as follows:

(1)The equilibria of system (6) for the weak Allee effect case are

(a) trivial equilibrium , that is, the extinction of both prey and predator point,(b) two axial equilibria, that is, the prey extinction point and the predator extinction point where

(c) positive or coexistence equilibrium , where are all possible real positive solutions of cubic equation

 where

(2)The equilibria of system (6) for the strong Allee effect case are

(a) trivial equilibrium , that is, the extinction of both prey and predator point,(b) three axial equilibria that are the prey extinction point and two predator extinction points: and , where

(c) positive or coexistence equilibrium point , where are also all possible real positive solutions of cubic equation (10).

Using transformation , (10) can be reduced to

where and . Implementing Cardan’s method as performed by Cai et al. [7], we obtain the existence of the positive equilibria as follows.

Lemma 5 (existence of positive equilibria). Let be the interior equilibrium of model (6) for both weak and strong Allee effects where is a real positive solution of (10). Then the following statements hold:

(a) If , then (13) has a single positive root . As a result, model (6) has a unique positive equilibrium point, that is, , with .(b) Suppose that and , then

(b1)If , then (13) has a positive root of multiplicity two. Thus, model (6) has a unique positive equilibrium point, that is, .(b2)If , then (13) has two positive roots and . Thus, model (6) has two positive equilibrium points, namely, and , with .

(c) If and , then (13) has a unique positive root . As a result, model (6) has a unique positive equilibrium point, that is, with

Moreover, algebraic computations show that if (13) has two positive roots, then they are

If (13) has a positive root, then it must be

To check the local stability of each equilibrium point, we linearize system (6) around the equilibrium and verify all eigenvalues of the Jacobian matrix evaluated at the equilibrium. The stability properties of trivial and axial equilibrium points for the case of weak and strong Allee effect are, respectively, stated in Theorems 6 and 7.

Theorem 6. Stability of trivial and axial equilibrium for weak Allee effect :

(i)the trivial equilibrium and the axial equilibrium are always unstable;(ii)the axial equilibrium is asymptotically stable if .

Proof. (i) The Jacobian matrix at is , and the eigenvalues are and . It is clear that . Hence is unstable. The Jacobian matrix at is where one of its eigenvalues is and therefore is unstable because , .
(ii) The Jacobian matrix at is . The eigenvalues of are and . Thus and whenever . This proves part (ii).

Using the same argument as in the proof of Theorem 6, we obtain the following stability properties of equilibria for the case of strong Allee effect.

Theorem 7. Stability of trivial and axial equilibrium for strong Allee effect :

(i)the trivial equilibrium ; the axial equilibria: and are always unstable;(ii)the axial equilibrium is always asymptotically stable.

The stability properties of positive (coexistence) equilibrium for the case of weak and strong Allee effect are stated in Theorems 8 and 9.

Theorem 8. Stability of coexistence equilibrium for weak Allee effect :
suppose is the Jacobian matrix at coexistence equilibrium . Equilibrium is asymptotically stable if one of the following mutually exclusive conditions holds:

(i); and .(ii), and

Proof. The characteristics equation of is given by

(i)If ; ; and then ; hence and the result follows.(ii)If is an eigenvalue of and , then is also an eigenvalue. Using , we have that . Therefore the stability of follows.

Similarly we have Theorem 9 for the stability of coexistence equilibrium for strong Allee effect case.

Theorem 9. Stability of coexistence equilibrium for strong Allee effect :
suppose is the Jacobian matrix evaluated at the coexistence equilibrium . Equilibrium is asymptotically stable if one of the following mutually exclusive conditions holds:

(i); ; and

(ii), , and

Based on the above theorems it can be seen that the stability properties of both trivial and axial equilibrium points are not dependent on (order of fractional derivative). But may influence significantly the stability of coexistence equilibrium point. Coexistence point can be asymptotically stable although the eigenvalue of Jacobian matrix has positive real part, provided that conditions of Theorem 8(ii) or Theorem 9(ii) are met. This is in contrast to the coexistence equilibrium point of the integer-order model (2) where coexistence point is asymptotically stable only if all real parts of the eigenvalues of the Jacobian matrix are negative [9].

3. Numerical Simulations

To solve system (6), we implement a nonstandard Grünwald-Letnikov scheme which is a combination of the Grünwald-Letnikov approximation [20] and the nonstandard finite difference (NSFD) method [21, 22]. According to [20], the explicit (or implicit) Grünwald-Letnikov (GL) approximation for a fractional differential equation with initial value

is given by

where ; ; and Here, represents the time step of numerical integration. The Grünwald-Letnikov approximation is proceeding iteratively but the sum in the scheme becomes longer and longer which represents the memory effects. Scherer et al. [20] have shown that the coefficient is positive and satisfies for . Observe that in the standard Grünwald-Letnikov approximation (21), the right hand side of (20) is approximated locally. We implement a nonstandard method which is adopted from the NSFD method [23]. A numerical scheme for an initial value problem

is called a NSFD method if at least one of the following conditions is satisfied [21, 22]:

(i) The left hand side is approximated by the generalization of forward difference scheme

 The nonnegative denominator function has to satisfy .(ii) The approximation of is nonlocal.

By implementing the Grünwald-Letnikov approximation for the fractional derivative and the nonlocal approximation for the right hand side of system (6) we get the following scheme:

Observe that scheme (24) is explicit and hence it is simple and easy to be implemented. Besides that, the nonstandard Grünwald-Letnikov scheme (24) also maintains the positivity solutions.

To verify our stability analysis as well as the effectiveness our numerical scheme, we perform some numerical simulations. First we use hypothetic values of parameters , , , , , , , , and . Model (6) with these parameters has four equilibrium points: , , , and . According to Theorem 7, only axial equilibrium is stable for any order of fractional derivative , where . This stability behavior is confirmed by our numerical solutions; see Figures 1 and 2. This shows that strong Allee effect may lead to an extinction of prey population. It is shown in Figure 2 that our numerical solutions for are convergent to the axial equilibrium , indicating that equilibrium is stable asymptotically for any order of fractional derivative. Detail observation shows that as the order of fractional derivative increases the convergence of solution is faster and the solution of system (6) closes to the integer-order model (2).

(a)

(b)

Figure 1: Phase portrait of system (6). The values of parameters are , , , , , , , , and .

(a)

(b)

Figure 2: Solutions of system (6) with initial value (1, 0.2) and parameter values: , , , , , , , , and , for various values of . All numerical solutions are convergent to axial equilibrium (0, 0.5). Solution of system (6) with larger order of fractional derivative has faster convergence compared to that with smaller .

Next, we set and and keep the rest of parameters as in Figure 1. This weak Allee case has a trivial equilibrium ; two axial equilibria: and ; and two interior points: and . Theorems 7 and 9 state that axial equilibrium is locally stable for and interior point is locally stable if , and other equilibria are always unstable. Such behavior is in accordance with our numerical results depicted in Figures 3(a) and 3(b). It is clearly seen in Figure 3(a) that produces bistable dynamic where depending on the initial values, solutions may be convergent to the extinction of prey point or to interior point. In other words, the solution of system (6) is highly sensitive to the initial conditions. An initially relatively small prey will converge to the prey extinction point. On the other hand, if the prey is initially relatively large then prey and predator will coexist. If we increase the order of derivative such that , the axial equilibrium is still locally stable but the interior point becomes unstable; see Figure 3(b). In latter case, there exists a stable limit cycle which shows that both prey and predator are fluctuating around the interior point. However, the appearance of limit cycle may be suppressed by increasing the coefficient of predator interference. For example we plot in Figure 4(a) the numerical solution using the same parameters as in Figure 3(b) but with . We see that the interior point is now stable while the axial equilibrium point is unstable. It can be said that relatively large predator interference can stabilize the interior point. On the other hand, strong Allee effect can destabilize or even remove the interior equilibrium and can cause the extinction of prey population. For example, we show numerical simulation using parameters the same as in Figure 3(b) except . This simulation shows that all initial values converge to the axial equilibrium which shows the extinction of prey population but predator species can still survive in the habitat because there is an enough environmental protection; see Figure 4(b).

(a)

(b)

Figure 3: Phase portrait of system (6). The values of parameters are , , , , , , , , and .

(a);

(b);

Figure 4: Phase portrait of system (6). The values of parameters are , , , , , and .

Finally, we compare our numerical results obtained by the NSGL scheme to those obtained by the standard GL scheme using parameters , , , , , , , and ; see Figure 5. We see that both NSGL and GL schemes using produce solutions where their difference cannot be observed in the scale of Figure 5. Using , the numerical solutions of both NSGL and GL schemes are in excellent agreement with those of both schemes using . If we take time step , both schemes have comparable solutions which are initially distorted from solutions with much smaller time step; see Figure 5(a). In Figure 5(b), we plot solutions using relatively large time step ( and ). Although the NSGL scheme has solutions which are quantitatively different from solution with , nevertheless those solutions still have the same behavior as before; that is, they are always positive and convergent to the correct equilibrium point. However, the GL scheme in this case gives unrealistic negative value for prey population. If the time step is further increased, the GL scheme will be unstable and leads to blowing up solutions.

(a)

(b)

Figure 5: Phase portrait of system (6) calculated using nonstandard (NSGL) and standard (GL) schemes. The values of parameters are , , , , , , , , and .

4. Conclusion

The dynamic of a fractional order modified Leslie-Gower model with Beddington-DeAngelis functional response and additive Allee effect has been analyzed. Our model has four types of equilibria that are the trivial (extinction of both prey and predator) equilibrium, two axial equilibria (the prey extinction point and the predator extinction point), and the interior (coexistence) point. The trivial and the predator extinction for both weak and strong Allee effects are always unstable. For the case of weak Allee effect, the prey extinction is conditionally stable while for that of strong Allee effect, the prey extinction is always stable. Our analysis also shows that the interior point for both weak and Allee effects is conditionally stable. The order of fractional derivative may influence the stability of interior point. Here, when the order is larger than critical order , then the interior point may be destabilized. These dynamical properties are confirmed by our NSGL schemes which shows the effectiveness of NSGL scheme. It is also shown that the NSGL scheme preserves the positivity of numerical solutions. Furthermore, our numerical results show that the NSGL scheme produces numerical solutions which satisfy the dynamical behavior of our model. However, the standard GL scheme may fail to preserve such properties; for example, it can produce nonrealistic negative solutions.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

This work was supported by the Directorate of Research and Community Service, The Directorate General of Strengthening Research and Development, and the Ministry of Research, Technology and Higher Education (Brawijaya University), Indonesia, Contract no. 137/SP2H/LT/DRPM/III/2016 dated March 10, 2016, and Contract no. 460.18/UN10.C10/PN/2017 dated April 18, 2017.

References

A. A. Berryman, “The origins and evolution of predator-prey theory,” Ecology, vol. 73, no. 5, pp. 1530–1535, 1992.

M. A. Aziz-Alaoui and M. D. Okiye, “Boundedness and global stability for a predator-prey model with modified leslie-gower and holling-type II schemes,” Applied Mathematics Letter, vol. 16, no. 7, pp. 1069–1075, 2003.

S. Yu, “Global asymptotic stability of a predator-prey model with modified Leslie-Gower and Holling-type II schemes,” Discrete Dynamics in Nature and Society, vol. 2012, Article ID 208167, 2012.

S. Yu, “Global stability of a modified Leslie-Gower model with Beddington-DeAngelis functional response,” Advances in Difference Equations, 2014:84, 14 pages, 2014.

B. Sahoo, “Role of additional food in eco-epidemiological system with disease in the prey,” Applied Mathematics and Computation, vol. 259, pp. 61–79, 2015.

P. J. Pal and P. K. Mandal, “Bifurcation analysis of a modified Leslie-Gower predator-prey model with Beddington-DeAngelis functional response and strong ALLee effect,” Mathematics and Computers in Simulation, vol. 97, pp. 123–146, 2014.

Y. Cai, C. Zhao, W. Wang, and J. Wang, “Dynamics of a Leslie-Gower predator-prey model with additive Allee effect,” Applied Mathematical Modelling. Simulation and Computation for Engineering and Environmental Systems, vol. 39, no. 7, pp. 2092–2106, 2015.

M.-H. Wang and M. Kot, “Speeds of invasion in a model with strong or weak Allee effects,” Mathematical Biosciences, vol. 171, no. 1, pp. 83–97, 2001.

D. Indrajaya, A. Suryanto, and A. R. Alghofari, “Dynamics of modified Leslie-Gower predator-prey model with Beddington-DeAngelis functional response and additive Allee effect,” International Journal of Ecology and Development, vol. 31, no. 3, pp. 60–71, 2016.

E. Ahmed, A. M. El-Sayed, and H. A. El-Saka, “Equilibrium points, stability and numerical solutions of fractional-order predator-prey and rabies models,” Journal of Mathematical Analysis and Applications, vol. 325, no. 1, pp. 542–553, 2007.

S. Das and P. K. Gupta, “A mathematical model on fractional Lotka-Volterra equations,” Journal of Theoretical Biology, vol. 277, pp. 1–6, 2011.

Y. Ding, Z. Wang, and H. Ye, “Optimal control of a fractional-order HIV-immune system with memory,” IEEE Transactions on Control Systems Technology, vol. 20, no. 3, pp. 763–769, 2012.

F. A. Rihan, “Numerical modeling of fractional-order biological systems,” Abstract and Applied Analysis, vol. 2013, Article ID 816803, 2013.

M. Javidi and N. Nyamoradi, “Dynamic analysis of a fractional order prey-predator interaction with harvesting,” Applied Mathematical Modelling. Simulation and Computation for Engineering and Environmental Systems, vol. 37, no. 20-21, pp. 8946–8956, 2013.

F. A. Rihan, D. Baleanu, S. Lakshmanan, and R. Rakkiyappan, “On fractional SIRC model with Salmonella bacterial infection,” Abstract and Applied Analysis, Article ID 136263, Art. ID 136263, 9 pages, 2014.

A. E. Matouk, A. A. Elsadany, E. Ahmed, and H. N. Agiza, “Dynamical behavior of fractional-order Hastings-Powell food chain model and its discretization,” Communications in Nonlinear Science and Numerical Simulation, vol. 27, no. 1-3, pp. 153–167, 2015.

I. Podlubny, Fractional Differential Equations, vol. 198 of Mathematics in Science and Engineering, Academic Press, San Diego, Calif, USA, 1999.

Y. S. Mishura, Stochastic Calculus for Fractional Brownian Motion and Related Processes, vol. 1929 of Lecture Notes in Mathematics, Springer, Berlin, Germany, 2008.

D. Matignon, “Stability results for fractional differential equations with applications to control processing,” in Computational Engineering in Systems Applications, pp. 963–968, Lille, France, 1996.

R. Scherer, S. L. Kalla, Y. Tang, and J. Huang, “The Grnwald-Letnikov method for fractional differential equations,” Computers and Mathematics with Applications, vol. 62, no. 3, pp. 902–917, 2011.

R. E. Mickens, “Nonstandard finite difference schemes,” in Applications of nonstandard finite difference schemes (Atlanta, {GA}, 1999), pp. 1–54, World Sci. Publ., River Edge, NJ, 2000.

R. E. Mickens, Application of Nonstandard Finite Difference Schemes, World Scientific Publishing Co Pte. Ltd, Singapore, 2000.

A. J. Arenas, G. González-Parra, and B. M. Chen-Charpentier, “Construction of nonstandard finite difference schemes for the SI and SIR epidemic models of fractional order,” Mathematics and Computers in Simulation, vol. 121, pp. 48–63, 2016.