Abstract
The complexity of the dynamical behaviors of interaction between prey and its predator is studied. The prey and predator relationship involves the age structure and intraspecific competition on predators and the nonlinear harvesting of prey following the Michaelis–Menten type term. Some biological validities are shown for the constructed model such as the existence and uniqueness as well as the non-negativity and boundedness of solutions. Three equilibrium points, namely the origin, axial, and interior points, are found including their global dynamics by employing the Lyapunov function along with the generalized Lassale invariant principle. The changes in dynamical behaviors driven by the harvesting and the memory effect are exhibited, including transcritical, saddle-node, backward, and Hopf bifurcations. The appearance of these interesting phenomena is strengthened by giving numerical simulations consisting of bifurcation diagrams, phase portraits, and their time series.
1. Introduction
Since Lotka and Volterra introduced the classical predator–prey model, theoretical studies of predation without age structure have attracted the attention of many authors, for example Deng et al. [1], Huang et al. [2], Tahara et al. [3], and Zeng et al. [4]. However, in nature, many species of plants and animals could have life histories that can simply be partitioned into two age stages: immature and mature stages. In each stage, individuals of species have identical biological characteristics, such as the ability to reproduce, motile, ingest food, and survive [5]. In particular, there are amphibians, insects, birds, and mammals with life cycles that can last from only several days or weeks to more than a century. For this reason, some researchers have developed the predator–prey model by incorporating age structure either in prey or/and predator population with other factors that also influence the dynamics of the predator–prey model, mainly restricted to the classical integer-order, stochastic, or delay equations [6–13].
In Wang and Chen [14] considered the predator–prey model with age structure for the predator population using time delays. If we ignore the effect of time delay, the model can be written as follows: Here x(t), y(t), and z(t) represent the population densities of prey, immature predator, and mature predator, at time t, respectively. Model (Equation 1) assumes that the prey grows logistically with r as the intrinsic growth rate, K is the carrying capacity; m is the linear Holling type I functional response, n is the conversion rate with which captured prey are converted to new immature predator, β is the maturity rate of the predator, δ1 and δ2 are the death rate of the immature and mature predators, respectively. It is also assumed that only the mature predator can feed the prey through the term mxz. If we do not consider the age structure of the predator population, then model (Equation 1) is reduced to the classical Lotka–Volterra model for which the positive equilibrium or the boundary equilibrium of this model is globally asymptotically stable. This means that the model has no periodic solution. On the other hand, Wang and Chen [14] prove that in the model (Equation 1), there exists an orbitally asymptotically stable periodic solution around the interior equilibrium point which suggests that the age structure can cause periodic oscillation of populations.
From the point of view of human needs, harvesting of populations generally occurs in wildlife, forestry, and fisheries management. When harvesting is integrated into the predator–prey model, there are three types of harvesting, namely constant harvesting [15], linear harvesting [16], and non-linear harvesting [17]. In this article, we assume that the predator is not a commercial species and there is intraspecific competition among immature predators. Therefore, the predator–prey model with age structure and intraspecific competition in predator (Equation 1), where the prey population is subject to Michaelis–Menten type harvesting, is given by An example of prey-predator interactions whose biological phenomena are described in the model (Equation 2) can be found in the African wild dog with its prey impala. The African wild dogs are a social structure that lives in packs. For 3–4 weeks, young African wild dogs were in the den with their mother. All adult members of African wild dogs are care for their young ones and provide food for them. The hunting members of the pack will return to the den where they regurgitate meat for the nursing female and young. In some cases, young ones fail to survive because the hunting member does not bring back sufficient food for the young, which leads to intraspecific competition in immature predator [18]. On the other hand, the prey, impala, even though there are no major threats to their survival, poaching has become significantly contributed to the decline in its number [19].
Note that the growth rates of the prey, immature, and mature predator populations in the model (Equation 2) depend only on the local state as the left-hand side is the integer-order derivative. On the other hand, most biological systems have properties where the current state is affected by all of the past states or it is called the memory effect. Therefore, modeling with memory effects can be done by analyzing the system using fractional-order calculus [20, 21]. The operators of the fractional-order derivative have non-local properties to make them more suitable for dynamical systems that have memory influences on their state variables.
After Riemann and Liouville generalized the concept of integer-order calculus to the fractional-order calculus over two decades ago, the studies about the predator–prey models with fractional-order differential equation have gained much attention, for example, Rahmi et al. [22], Owolabi [23], Barman et al. [24], Ghanbari and Djilali [25], Yousef et al. [26], Ghosh et al. [27], and Panigoro et al. [28]. The fractional-order derivatives are defined as an integration that provides the ability to store the whole memory over time, and hence, it could give an exact description of different ecological phenomena. For this reason, the new structure for the model (Equation 2) is given in the following form: The existence and local stability of all equilibrium points of the model (Equation 3) are discussed in Panigoro et al. [29]. However, to the best of our knowledge, the global dynamics and bifurcation analysis of the model (Equation 3), to this day, have not been investigated. Here, is the standard Caputo derivative for a continuous function f(x) ∈ ℂ([0, +∞), ℝ), which is defined as follows: where Γ(∗) is the gamma function, t ≥ 0, and 0 < α ≤ 1 is known as the order of the fractional derivative.
Based on the above discussion, we have organized our work in several sections: In Section 3, we develop the existence and uniqueness solution of the model (Equation 3). To check the biologically well-posedness of the model, we establish the non-negativity and boundedness of solutions of the model in Section 3. In Section 4, we derive some sufficient conditions to ensure the global asymptotic stability of each locally asymptotically stable equilibrium point by applying the Lyapunov functions. We then analyze the existing conditions of transcritical, saddle-node, backward, and Hopf bifurcations in Section 5. Some numerical simulations of our obtained results are carried out in Section 6. Finally, the conclusions are given in Section 7.
2. Existence and uniqueness
In this section, we will show that the model (Equation 3) has a unique solution. A similar manner given by Mahata et al. [30] is adopted. We first take the Riemann–Liouville integral (Definition 1 in Yavuz and Sene [31]) on both sides of Equation (3) to achieve the following Volterra-type integral equations. Now, we will show that the kernels Fi(x, y, z), i = 1, 2, 3, satisfy the Lipschitz condition. For ∥.∥ is the Euclidean norm, we suppose that ∥x(t)∥ ≤ a1, , ∥y(t)∥ ≤ a3, ∥(t)∥ ≤ a4, ∥z(t)∥ ≤ a5, and are bounded functions. For x, , y, , z, and , we have where , g2 = β + δ1 + ω(a3 + a4), and g3 = δ2. Therefore, we conclude that Fi(x, y, z), i = 1, 2, 3, satisfy the Lipschitz condition. Furthermore, it is clear that if 0 ≤ gi < 1, then Fi(x, y, z) are contractions for i = 1, 2, 3. Therefore, the following theorem is obtained.
Theorem 1. The kernel Fi(x, y, z), i = 1, 2, 3 satisfy the Lipschitz condition and contractions if 0 < gi < 1, i = 1, 2, 3.
Next, Equation (5) can be written as follows: Which can be written by the following recursive formula with initial conditions x0(t) = x(0), y0(t) = y(0), and z0(t) = z(0). Therefore, we have where , , and . Now, we evaluate the norm of Equation (6). We achieve From Theorem 1, we have that the kernel satisfy the Lipschitz condition and hence Equation (7) becomes The last inequality gives Finally, the existence of a solution is given by the following theorem.
Theorem 2. The solution of model (Equation 3) has a solution under the condition if we have t1 such that .
Proof. We assume that x(t), y(t), and z(t) are bounded and their kernel Fi, i = 1, 2, 3 satisfy the Lipschitz condition. According to Equation (8), we obtain which represent the existence and continuity of the system. To show that the solution of the model (Equation 3) can be set up from the functions in Equation (9), we assume where Qin, i = 1, 2, 3 are the remaining terms. Furthermore, the given terms would be demonstrated to hold Qi∞ → 0 ∀i = 1, 2, 3. Denoting that By applying a recursive pattern, we acquire At the point t1, we have By considering the results of Theorem 1 and applying n → ∞ to both sides, we have ∥Qi∞∥ → 0 ∀i = 1, 2, 3.
In the end, we will show that the solution is unique for each initial value by utilizing the contradiction approach. Suppose that there exists another solution of the model (Equation 3), namely x1(t), y1(t), and z1(t). Then, we have Applying the norm on both sides, we achieve By considering Theorem 1, we obtain Therefore, the following equations are obtained. As a result, we achieve ∥x(t) − x1(t)∥ = 0, ∥y(t) − y1(t)∥ = 0, and ∥z(t) − z1(t)∥ = 0, which impact x(t) = x1(t), y(t) = y1(t), and z(t) = z1(t). Then, we finally give the following theorem.
Theorem 3. The Caputo fractional-order predator–prey model (Equation 3) has a unique solution.
3. Non-negativity and boundedness
In this section, we will show that for any initial condition is in where The solution not only exists and is unique but also bounded and always in as t → ∞. Therefore, we have the following two theorems.
Theorem 4. If the initial condition in then both population densities of prey and predator given by model (Equation 3) remain in .
Proof. To prove this non-negativity condition, we apply reductio ad absurdum (contradiction method), which is also used in Barman et al. [24] and Maji [32]. We assume that there exists such that From the first equation in Equation (3) along with Equation (17), we have According to Lemma 3.1 in Barman et al. [33], we get which contradicts with Equation (17) where . This means that x(t) ≥ 0 for all t ∈ [0, ∞]. Similarly, we can show that y(t) ≥ 0 and z(t) ≥ 0 for all t ∈ [0, ∞]. In conclusion, we have the non-negative solution for model (Equation 3) when the initial values are in .
Theorem 5. The solution of model (Equation 3) is always bounded in for the initial condition in .
Proof. Since we work the population model, it is natural that the population must be bounded due to the limitation of their biological resources, which is also known as environmental carrying capacity. Thus, the boundedness of the solution of the model (Equation 3) is also important to learn and prove. From Theorem 4, we can define a positive function as follows: For any γ > 0, the following fractional-order differential equation holds. By choosing γ < min{δ1, δ2}, we obtain According to Lemma 3 in Panigoro et al. [34], we apply the comparison principle and obtain For t → ∞, we achieve , which means all solutions of model (Equation 3) with non-negative initial conditions are confined to where
4. Global dynamics
In this section, the global dynamics of model (Equation 3) are investigated. Note that all biological equilibrium points, their existence conditions, and their local stability are successfully described in Panigoro et al. [29], which can be rewritten by the following theorem.
Theorem 6.
The origin point always exists. It is locally asymptotically stable (LAS) if .
The axial point where is the positive root of , which has
an equilibrium point if .
a pair of equilibrium points if .
Moreover, it is LAS if and .
The interior point exists, if ai, i = 2, 3 in Panigoro et al. [29] satisfies the following statements.
An equilibrium point in if a3 < 0.
Two equilibrium points in if a2 < 0 and a3 > 0.
The LAS condition of can be seen in Theorem 4 in Panigoro et al. [29].
Note that all equilibrium points may attain local asymptotic stability with several biological conditions. Now, we will identify the biological properties to obtain globally asymptotically stable (GAS) for each equilibrium point. The analytical results are given by the following three theorems.
Theorem 7. The origin point is GAS if .
Proof. We define the positive definite Lyapunov function as follows: By applying Lemma 3.1 in Vargas-De-León [35], we compute the α−order derivative of along the solution of the model (Equation 3) as follows: From Equation (21), we have x ≤ σ and hence Therefore, for all , if . We also find that , if (x, y, z) = (0, 0, 0). This conveys that is the only invariant set on which . Obeying Lemma 4.6 in Huo et al. [20], obviously becomes the biological condition of to reach GAS.
Theorem 8. The axial point is GAS if .
Proof. We construct a positive definite Lyapunov function based on the Volterra equation as follows: The α-order derivative of along the solution of the model (Equation 3) given by Lemma 3.1 in Vargas-De-León [35] is given by Since , we have for all . It is also clear that if . This confirms that is the only invariant set on which . Therefore, is GAS due to Lemma 4.6 in Huo et al. [20]. This confirms the justifiability of Theorem 8.
Theorem 9. Let and . The interior point is GAS in ΩX.
Proof. Consider a positive definite Lyapunov function as follows: By applying Lemma 3.1 in Vargas-De-León [35] and Lemma 1 in Aguila-Camacho et al. [36], we obtain the α-order derivative of as follows: Applying Equation (21), we have Since , we achieve Thus, for all , when . We also confirm that if (x, y, z) = (x*, y*, z*) and hence is the only invariant set on which . Based on Lemma 4.6 in Huo et al. [20], the interior point is GAS in ΩX. This ends the proof.
5. Bifurcation analysis
In this section, we will study the occurrence of several phenomena namely transcritical, saddle-node, backward, and Hopf bifurcations. Two parameters are chosen, namely the harvesting rate (h) and the order of the derivative (α), as the memory index. For the analytical purpose, we define the following parameter. Next, the following theorem is given for describing the occurrence of transcritical bifurcation driven by the harvesting rate (h) as the bifurcation parameter.
Theorem 10. The origin point and the axial point exchange their stability via transcritical bifurcation when h passes through .
Proof. Since the axial consists of two equilibrium points, we focus on the axial point nearest to the origin point. When , the axial point merge with the origin point where the eigenvalues of the Jacobian matrix are: λ1 = 0, λ2 = (β + δ1), and λ3 = −δ2. We obtain arg(λ2,3) = π > απ/2 while arg(λ1) = απ/2. This means is non-hyperbolic. When , by applying Theorems 2 and 3 in Panigoro et al. [29], becomes LAS while the nearest becomes a saddle point. For , The origin becomes unstable and the nearest becomes unstable. This condition shows the existence of transcritical bifurcation, where h becomes the bifurcation parameter while is the bifurcation point.
Now, the existence of saddle-node bifurcation on axial will be proven by still regarding the harvesting rate (h) as the bifurcation parameter. As a result, the following theorem is proposed.
Theorem 11. Suppose that . The axial point undergoes saddle-node bifurcation when h passes through the bifurcation point .
Proof. According to Theorem 1 in Panigoro et al. [29], the axial point does not exist when . When , a unique equilibrium point occurs in axial where its Jacobian matrix has three eigenvalues: λ1 = 0 and . Since |arg(λ1)| = απ/2, this axial point is non-hyperbolic. When h < h*, two axial points occurs given by and , where and . It is easy to validate that both and are in and have different stability. As a consequence, all the given circumstances express the occurrence of saddle-node bifurcation.
Based on Theorems 10 and 11, we obtain more global bifurcation namely backward bifurcation given by the following lemma.
Lemma 1. The model (Equation 3) undergoes backward bifurcation driven by the harvesting rate (h).
Proof. From previous theorems, the axial point exists and is LAS, while is unstable when . When , becomes LAS, still exists and is LAS, and unstable occurs. The bistability condition is held for this interval of h which means that the convergence of the solution is very sensitive to the initial condition. Finally, those two axial points merge when and disappear when . This completes the proof.
Finally, we will show that the memory index in this case, that is, the order of the derivative (α), affects the dynamical behaviors of the model (Equation 3) indicated by the appearance of Hopf bifurcation around the interior point .
Theorem 12. Suppose the characteristic equation of the Jacobian matrix evaluated at can be written as , which has a pair of complex conjugate eigenvalues λ1,2 = ζ1 ± iζ2, where ζ1 > 0 and one real negative eigenvalue (λ3 < 0). The model (Equation 3) undergoes a Hopf bifurcation when the order of the fractional derivative α crosses out the critical value
Proof. From the earlier assumptions, we have . Therefore, the solution of is only when If we check the transversal condition: which is not equal to 0, we can assure that the sign of m(α) changes when the bifurcation parameter α passing by α*. It means that the equilibrium point is stable when α ∈ (0, α*) and is unstable for α* < α < 1.
6. Numerical simulations
In this section, we explore the dynamical behaviors of the model (Equation 3) numerically to support analytical findings, especially the occurrence of bifurcation. The predictor–corrector scheme given by Diethelm et al. [37] is employed. All of the parameters used in these simulations are assumptions matched with the biological conditions given by the previous analysis results. This decision was taken because this work does not specifically address an ecological case involving a particular organism.
To show the occurrence of several bifurcations driven by the harvesting rate (h), we first set the parameter as follows: By varying the harvesting rate in the interval 0 ≤ h ≤ 0.24, the bifurcation diagram is portrayed as in Figure 1. We have three types of dynamic behaviors around the axial point. When , we have unstable origin point and LAS . The origin point losses its stability via transcritical bifurcation when h crosses and the unstable axial point occurs simultaneously. These dynamics are maintained for interval . On the other hand, the stable branch of axial point is preserved for . The LAS point and unstable point of merge into the non-hyperbolic point when . The axial point finally disappeared when h passes through while the sign of does not change. Thus, we have saddle-node bifurcation on axial with as the bifurcation point. If we observe from a more global point of view, these interesting phenomena represent the existence of backward bifurcation marked by the occurrence of bistability condition. To show these dynamical behaviors, we choose three values of harvesting rate in each interval: h = 0.02, 0.12, and 0.18 and portray them as phase portraits and time series (see Figures 2–4). The interesting phenomenon called bistability is portrayed in Figure 3. Two equilibrium points LAS simultaneously impact the sensitivity of the convergence of the solution to the selection of the initial value. The two closest initial values are set which converge to the different equilibrium points. One of them convergent to the origin point and the other solution convergent to the axial point. This means, two conditions may arrive, namely the extinction of all populations and the only prey existence point. Several references show that the bistability condition occurs as the consequence of saddle-node bifurcation, see Adhikary et al. [38] and several references therein.
Figure 1
Figure 2
Figure 3
Figure 4
From the biological point of view, these numerical bifurcations describe the possibility for the prey to extinct or survive due to the change in the harvesting rate while the predator both mature and immature is always extinct (Figure 1). Three feasible conditions may happen. First, for any sufficiently small harvesting rate (), the prey population may maintain its existence in this ecosystem (Figure 2). Second, if the harvesting rate increases (), two possible conditions may occur namely the extinction of prey or the viability of prey. These circumstances depend on the initial condition, where if the initial condition is quite close to the origin point, the prey will be extinct, and for the initial condition is far enough from the origin point, the prey can survive extinction (Figure 3). Third, if the harvesting rate is again increased (), the population of prey will finally extinct and there are no population again in this ecosystem (Figure 4).
The next circumstance occurs in the interior point of the model (Equation 3), which demonstrates the influence of the order of the derivative as the memory index on the dynamical behaviors around the interior point. We set the parameter as follows: To identify the dynamical behaviors, we vary the values of α in the interval 0.76 ≤ α ≤ 1. As a result, we obtain the bifurcation diagram given in Figure 5. For α < α* ≈ 0.86, the interior point is LAS. To show this condition, we give the phase portraits by selecting α = 0.81 and α = 0.84 as given in Figures 6A, B. Nearby solution oscillates and convergent to . When α crosses α* ≈ 0.86, losses its stability via Hopf bifurcation which is indicated by the occurrence of periodic signal namely limit-cycle. The nearby solution stays away from and convergent to the limit-cycle. The evolution of the limit-cycle given in Figure 5 also shows that the diameter of the limit-cycle increases when alpha increases. We portray the phase portraits in Figures 6C, D to show the dynamics of solutions around for α = 0.87 and α = 0.9. It is shown that the densities of all populations are oscillated and finally converge to the limit cycle. The physical interpretations of Hopf bifurcation driven by α are closely related to the influence of the memory on the change of behaviors around the interior point. The stronger the influence of memory, the higher the ability of prey and predators to maintain their existence (α < α*). For less memory effect which is indicated by α > α*, all populations lose the ability to stabilize their number of individuals. The population density fluctuates according to seasonal patterns which indicates by the appearance of a limit cycle (Figures 6C, D), and the peak of each population also increases for less memory effect (Figure 7). Although the density for each population cannot tend to a constant number, in this case, the memory effect cannot cause the extinction of every population.
Figure 5
Figure 6
Figure 7
7. Conclusion
The dynamics of a predator–prey model incorporating four biological conditions, namely age structure, intraspecific competition, Michaelis–Menten type harvesting, and memory effect, have been studied. All biological validities have been presented such as the existence, uniqueness, non-negativity, and boundedness of the solution. The dynamics of the model have been explored by showing the global stability conditions for three equilibrium points namely the origin, the axial, and the interior points. We have shown that three bifurcations phenomena driven by the harvesting rate occur around the axial point namely transcritical, saddle-node, and backward bifurcations. The bistability condition exists as the impact of saddle-node bifurcation which states that the existence of prey depends on the initial condition. A bifurcation driven by the memory effect also has been shown around the interior point which is called Hopf bifurcation. All the bifurcation phenomena having an impact on the densities of the population not only may reduce their densities but also threaten the existence of several populations.
Statements
Data availability statement
The original contributions presented in the study are included in the article/supplementary material, further inquiries can be directed to the corresponding author/s.
Author contributions
All authors listed have made a substantial, direct, and intellectual contribution to the work and approved it for publication.
Funding
This research was funded by LPPM-UNG via PNBP-Universitas Negeri Gorontalo according to DIPA-UNG No. 023.17.2.677521/2021, under contract No. B/176/UN47.DI/PT.01.03/2022.
Acknowledgments
The authors intended to express our gratitude and appreciation to the editors, reviewers, and all those who have supported us in improving this article.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
References
1.
DengHChenFZhuZLiZ. Dynamic behaviors of Lotka-Volterra predator-prey model incorporating predator cannibalism. Adv Diff Equat. (2019) 2019:359. 10.1186/s13662-019-2289-8
2.
HuangQLinYZhongQMaFZhangY. The impact of microplastic particles on population dynamics of predator and prey: implication of the Lotka-Volterra model. Sci Rep. (2020) 10:1–10. 10.1038/s41598-020-61414-3
3.
TaharaTGavinaMKAKawanoTTubayJMRabajanteJFItoHet al. Asymptotic stability of a modified Lotka-Volterra model with small immigrations. Sci Rep. (2018) 8:1–7. 10.1038/s41598-018-25436-2
4.
ZengXLiuLXieW. Existence and uniqueness of the positive steady state solution for a Lotka-Volterra predator-prey model with a crowding term. Acta Math Sci. (2020) 40:1961–80. 10.1007/s10473-020-0622-7
5.
GuiZGeW. The effect of harvesting on a predator-prey system with stage structure. Ecol Modell. (2005) 187:329–40. 10.1016/j.ecolmodel.2005.01.052
6.
MagnússonKG. Destabilizing effect of cannibalism on a structured predator-prey system. Math Biosci. (1999) 155:61–75. 10.1016/S0025-5564(98)10051-2
7.
LiuCZhangQHuangJ. Stability analysis of a harvested prey-predator model with stage structure and maturation delay. Math Problems Eng. (2013) 2013:329592. 10.1155/2013/329592
8.
DubeyBAgarwalSKumarA. Optimal harvesting policy of a prey-predator model with Crowley-Martin-type functional response and stage structure in the predator. Nonlinear Anal Model Control. (2018) 23:493–514. 10.15388/NA.2018.4.3
9.
XiaoZLiZZhuZChenF. Hopf bifurcation and stability in a Beddington-DeAngelis predator-prey model with stage structure for predator and time delay incorporating prey refuge. Open Math. (2019) 17:141–59. 10.1515/math-2019-0014
10.
ZhangXLiuZ. Periodic oscillations in age-structured ratio-dependent predator-prey model with Michaelis-Menten type functional response. Physica D. (2019) 389:51–63. 10.1016/j.physd.2018.10.002
11.
ZhangFChenYLiJ. Dynamical analysis of a stage-structured predator-prey model with cannibalism. Math Biosci. (2019) 307:33–41. 10.1016/j.mbs.2018.11.004
12.
LuWXiaY. Periodic solution of a stage-structured predator-prey model with Crowley-Martin type functional response. AIMS Math. (2022) 7:8162–75. 10.3934/math.2022454
13.
LiJLiuXWeiC. The impact of role reversal on the dynamics of predator-prey model with stage structure. Appl Math Model. (2022) 104:339–57. 10.1016/j.apm.2021.11.029
14.
WangWChenL. A predator-prey system with stage-structure for predator. Comput Math Appl. (1997) 33:83–91. 10.1016/S0898-1221(97)00056-4
15.
ChowCHotiMLiCLanK. Local stability analysis on Lotka-Volterra predator-prey models with prey refuge and harvesting. Math Methods Appl Sci. (2018) 41:7711–32. 10.1002/mma.5234
16.
XiaoMCaoJ. Hopf bifurcation and non-hyperbolic equilibrium in a ratio-dependent predator-prey model with linear harvesting rate: analysis and computation. Math Comput Model. (2009) 50:360–79. 10.1016/j.mcm.2009.04.018
17.
ZhangZUpadhyayRKDattaJ. Bifurcation analysis of a modified Leslie-Gower model with Holling type-IV functional response and nonlinear prey harvesting. Adv Diff Equat. (2018) 2018:127. 10.1186/s13662-018-1581-3
18.
MulheisenMAllenCParrCS. Lycaon Pictus. (On-line), Animal Diversity Web. (2002). Available online at: https://animaldiversity.org/accounts/Lycaon_pictus/ (accessed on November 13, 2022)
19.
IUCN SSC Antelope Specialist Group. Aepyceros melampus. The IUCN Red List of Threatened Species 2016: eT550A50180828. (2016). Available online at: https://www.iucnredlist.org/species/550/50180828 (accessed on November 13, 2022)
20.
HuoJZhaoHZhuL. The effect of vaccines on backward bifurcation in a fractional order HIV model. Nonlinear Anal. Real World Appl. (2015) 26:289–305. 10.1016/j.nonrwa.2015.05.014
21.
MoustafaMMohdMHIsmailAIAbdullahFA. Stage structure and refuge effects in the dynamical analysis of a fractional order Rosenzweig-MacArthur prey-predator model. Progr Fract Different Appl. (2019) 5:49–64. 10.18576/pfda/050106
22.
RahmiEDartiISuryantoATrisilowatiPanigoroHS. Stability analysis of a fractional-order leslie-gower model with allee effect in predator. J Phys. (2021) 1821:012051. 10.1088/1742-6596/1821/1/012051
23.
OwolabiKM. Dynamical behaviour of fractional-order predator-prey system of Holling-type. Discrete Contin Dyn Syst S. (2020) 13:823–34. 10.3934/dcdss.2020047
24.
BarmanDRoyJAlamS. Modelling hiding behaviour in a predator-prey system by both integer order and fractional order derivatives. Ecol Inform. (2022) 67:101483. 10.1016/j.ecoinf.2021.101483
25.
GhanbariBDjilaliS. Mathematical and numerical analysis of a three-species predator-prey model with herd behavior and time fractional-order derivative. Math Methods Appl Sci. (2020) 43:1736–52. 10.1002/mma.5999
26.
YousefFBYousefAMajiC. Effects of fear in a fractional-order predator-prey system with predator density-dependent prey mortality. Chaos Solitons Fractals. (2021) 145:110711. 10.1016/j.chaos.2021.110711
27.
GhoshUPalSBanerjeeM. Memory effect on Bazykin's prey-predator model: stability and bifurcation analysis. Chaos Solitons Fractals. (2021) 143:110531. 10.1016/j.chaos.2020.110531
28.
PanigoroHSRahmiEAchmadNMahmudSLResmawanRNuhaAR. A discrete-time fractional-order rosenzweig-macarthur predator-prey model involving prey refuge. Commun Math Biol Neurosci. (2021) 2021:1–19. 10.28919/cmbn/6586
29.
PanigoroHSResmawanRSidikATRWalangadiNIsmailAHusunaC. A fractional-order predator-prey model with age structure on predator and nonlinear harvesting on prey. Jambura J Math. (2022) 4:355–66. 10.34312/jjom.v4i2.15220
30.
MahataAPaulSMukherjeeSDasMRoyB. Dynamics of caputo fractional order SEIRV epidemic model with optimal control and stability analysis. Int J Appl Comput Math. (2022) 8:1–25. 10.1007/s40819-021-01224-x
31.
YavuzMSeneN. Stability analysis and numerical computation of the fractional predator-prey model with the harvesting rate. Fractal Fract. (2020) 4:35. 10.3390/fractalfract4030035
32.
MajiC. Impact of fear effect in a fractional-order predator-prey system incorporating constant prey refuge. Nonlinear Dyn. (2022) 107:1329–42. 10.1007/s11071-021-07031-9
33.
BarmanDRoyJAlrabaiahHPanjaPMondalSPAlamS. Impact of predator incited fear and prey refuge in a fractional order prey predator model. Chaos Solitons Fractals. (2021) 142:110420. 10.1016/j.chaos.2020.110420
34.
PanigoroHSSuryantoAKusumahwinahyuWMDartiI. Dynamics of a fractional-order predator-prey model with infectious diseases in prey. Commun Biomath Sci. (2019) 2:105. 10.5614/cbms.2019.2.2.4
35.
Vargas-De-LeónC. Volterra-type Lyapunov functions for fractional-order epidemic systems. Commun Nonlinear Sci Num Simulat. (2015) 24:75–85. 10.1016/j.cnsns.2014.12.013
36.
Aguila-CamachoNDuarte-MermoudMAGallegosJA. Lyapunov functions for fractional order systems. Commun Nonlinear Sci Num Simulat. (2014) 19:2951–7. 10.1016/j.cnsns.2014.01.022
37.
DiethelmKFordNJFreedAD. A predictor-corrector approach for the numerical solution of fractional differential equations. Nonlinear Dyn. (2002) 29:3–22. 10.1023/A:1016592219341
38.
AdhikaryPDMukherjeeSGhoshB. Bifurcations and hydra effects in Bazykin's predator-prey model. Theor Populat Biol. (2021) 140:44–53. 10.1016/j.tpb.2021.05.002
Summary
Keywords
bifurcation, age structure, intraspecific competition, harvesting, memory effect
Citation
Panigoro HS, Rahmi E and Resmawan R (2023) Bifurcation analysis of a predator–prey model involving age structure, intraspecific competition, Michaelis–Menten type harvesting, and memory effect. Front. Appl. Math. Stat. 8:1077831. doi: 10.3389/fams.2022.1077831
Received
23 October 2022
Accepted
25 November 2022
Published
06 January 2023
Volume
8 - 2022
Edited by
Bapan Ghosh, Indian Institute of Technology Indore, India
Reviewed by
Huda Abdul Satar, University of Baghdad, Iraq; Prabir Panja, Haldia Institute of Technology, India
Updates
Copyright
© 2023 Panigoro, Rahmi and Resmawan.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Hasan S. Panigoro ✉ hspanigoro@ung.ac.id
This article was submitted to Mathematical Biology, a section of the journal Frontiers in Applied Mathematics and Statistics
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.