Optimal control of susceptible mature pest concerning disease-induced pest-natural enemy system with cost-effectiveness

Kunwer Singh Mathur; Bhagwan Kumar

doi:10.1515/cmb-2023-0120

Open Access Published by De Gruyter Open Access March 9, 2024

Optimal control of susceptible mature pest concerning disease-induced pest-natural enemy system with cost-effectiveness

Kunwer Singh Mathur and Bhagwan Kumar

From the journal Computational and Mathematical Biophysics

https://doi.org/10.1515/cmb-2023-0120

Abstract

This article addresses the pressing issue of pest outbreaks in India, which poses significant challenges for farmers and ecologists. A novel system is proposed for effective control that leverages natural enemies. Here, the pests are classified into juveniles and mature individuals, further categorized as susceptible or infected. The study introduces harvesting, incorporating external efforts and natural phenomena, in a pest-epidemic prey–predator system featuring a prey-stage structure. The model reveals three equilibria: trivial, boundary (indicating the absence of natural enemies), and interior equilibria. Notably, the trivial equilibrium is consistently unstable. As demonstrated by stability analysis, the survival or extinction of natural enemies hinges on control variables, including the harvesting rate, disease transmission rate, and natural death rate. Local stability is assessed using the Routh–Hurwitz criterion, while global stability is explored through the Lyapunov method. Furthermore, optimal control theory and Pontryagin’s maximum principle are applied for model optimization, unveiling crucial optimality conditions and determining the optimal harvesting rate for susceptible mature prey. Numerical computations validate theoretical insights, offering valuable guidance for formulating policies that optimize the control of susceptible adult pests within a disease-induced pest-natural enemy system, ensuring sustained cost-effectiveness.

Keywords: stage-structured pest; pest-natural enemy system; harvesting; optimal control; cost-effectiveness; eco-epidemic model

MSC 2010: 92-10; 92Bxx; 34Dxx

1 Introduction

The insect pests are unwanted living organisms that damage many crops, essential for human survival. In addition, a significant part of the economy of most developing countries, like India, depends on crop production. Several economic and ecological issues have arisen concerning pest outbreaks in a few decades. Thus, cost-effective pest control is needed to preserve crop production and enhance the economy of any developing country. The natural enemies are treated as the farmer’s friends as they can consume a large portion of pests naturally. Hence, the interaction of pests and natural enemies is reflected as prey–predator interaction. The mathematical study of pest-natural enemy prey–predator systems in any ecological environment is appealing to researchers due to its significance and universal acceptance. The history of the prey–predator system is not very old; biophysicist Lotka developed the first model in his study of “Element of physical biology” in 1925 [31], and the same model was also proposed by the prominent mathematician Vito Volterra in the same year with the help of a young zoologist Umberto d’Ancona, Volterra’s future son-in-law [17]. However, these two researchers were utterly unaware of each other’s work. Later on, the model was known as the Lotka–Volterra prey–predator model. This model globally attracted researchers due to its capacity to reveal the dynamics of two species, prey and predator, depending upon fundamental biological situations. In a number of research articles, several prey–predator models have been developed and studied with different factors based on some realistic conditions [9–11,14,15,21–25,35,40,43–45,47–49]. In last one decade, several attractive prey–predator models were proposed and studied. Panja and Mondal [39] proposed a three-species prey–predator model in a fishery system with Holling type II functional response. In this article, they studied the dynamic behaviors of phytoplankton, zooplankton, and fish and also determined the behavior of the system near the positive equilibrium point. A prey refuge model with the influence of anti-predator behavior due to the fear of predators with a Holling type II prey–predator model was investigated by Zhang et al. [53] and obtained the complex dynamics of the system along with Hopf bifurcation and limit cycle. The model is analyzed analytically and numerically to show the impact of the fear effect, which suggests that at the positive equilibrium, the density of predators not only reduced due to the fear effect but also stabilized the system. To protect a fixed quantity of prey, Das and Samanta studied a prey–predator model with the effect of additional food for predators and prey refuge in the study by Das and Samanta [11].

A pest-natural enemy epidemic model is a mathematical or ecological model used to study the dynamics of interactions between pests and their natural enemies within an ecosystem. Such models are often employed in agricultural and environmental research to understand how natural predators and parasites can help control pest populations. Researchers in the field of integrated pest management (IPM) have extensively studied various two-compartment models to enhance the understanding of these ecological dynamics [29,41,46,51]. Moreover, numerous pest species undergo a life cycle primarily divided into two stages: immature and mature, characterized by significant behavioral distinctions. For instance, a juvenile predator requires time to develop the ability to catch and hunt directly. Therefore, employing a mathematical model with a stage structure proves more effective in accurately representing biological scenarios than an ordinary model. Recently, several researchers have focused on extending models to incorporate diverse life stages within prey and predator populations, encompassing both immature and mature stages [6,14,23,28,30,36,37,52]. For instance, Ma et al. [33] studied and proposed a basic stage-structured single-species model with impulsive control and determined the dynamics of the system. However, the proposed single-species model has also been extended into a prey–predator model with a stage structure. Jatav and Dhar proposed a stage-structured prey–predator model with maturation delay for prey and gestation delay for predator [23]. A stage-structured prey–predator model with maturation and gestation delay is also extended by Dubey and Kumar in [14], with the assumption that the predator can consume immature prey in Holling type I manner, while the mature prey in Crowley–Maryin-type functional response manner. Thus, the stage-structured model has shown great interest in the prey and predator system and has worked to make it more realistic. The impacts of stage structure in the predator-prey ecosystem are essential topics. Furthermore, it is sensible to evaluate the effect of harvesting when studying a prey–predator model. The harvesting of populations is commonly practiced in fishery, forestry, and wildlife management [27]. Harvesting prey affects predator populations indirectly since it reduces the food population available in the area [2]. Prey harvesting also affects the stabilization or destabilization of prey–predator models. Several studies have dealt with the effect of harvesting on population dynamics, such as in previous studies [3,15,50]. In the study by Chakraborty et al. [7], an optimal control study was carried out for a prey–predator model with stage structure for prey and harvesting of mature prey and predator.

It is recognized that every species in the environment is susceptible to various infectious diseases, and these ailments influence the population size. The effect of infectious diseases on the prey-natural enemy system is significant and cannot be overlooked. In a prey–predator system, a prevalent scenario involves the transmission of parasites from one prey to another, leading to the division of the prey population into two distinct compartments: susceptible prey and infectious prey [8,19,20,42]. Consider aphids and beetles as an illustration. Aphids can become infectious due to parasitoids, which categorize them into two major forms: susceptible and infectious. In this context, beetles prefer consuming only susceptible aphids and avoiding attacking infectious aphids. The parasitoids disrupt the shape of aphids, leading to the beetles losing interest in preying on infectious aphids [18,26,38]. Moreover, aphids can be categorized into two significant life stages, determined by their behavior and reproductive capacity: immature and mature. The extensive scope of this classification underscores the relevance of the subject. Extending the conventional stage-structured pest-natural enemy model in alignment with this classification not only justifies the consideration’s applicability to reality but also presents a new perspective for a more realistic examination.

This article is organized as follows: section 2 introduces a mathematical model, while section 3 presents initial findings related to boundedness and the existence of equilibrium. The local stability of equilibrium is explored in section 4, followed by an analysis of global stability in section 5. Section 6 addresses the optimal control problem. This article concludes with numerical simulations in section 7 and final remarks in section 8.

2 Model development and proposition

In this section, we explore an epidemiological mathematical framework that incorporates stage structure in pest populations. The model encompasses two pivotal stages: juvenile and mature pests, further classifying the adult pest population into susceptible and infected mature pests. Additionally, external efforts contribute to harvesting the mature susceptible pest population.

Let x j ( t ) denote the immature pest population at time t , x m s ( t ) the susceptible mature pest population, x m I ( t ) the infectious mature pest population, and y ( t ) the population of natural enemies at time t . The assumption is made that only the susceptible mature pest contributes to growth, as infection suppresses the reproductive ability in the mature infectious pest population. The growth of the immature pest x j ( t ) is modeled with a rate r proportional to the mature susceptible pest population x m s ( t ) .

Natural enemies exhibit predation on mature susceptible pests with a predation rate β following a Holling type I manner. However, it is specified that they do not consume infected mature and juvenile pests. Notably, recent observations suggest predators may attack prey in a Holling type I manner. At the same time, the energy conversion due to this predation is represented through a Holling type II functional response, incorporating handling time and saturation effects [13,34].

The term λ β x m s y 1 + r 1 h 1 x m s captures the growth of natural enemies, considering energy conversion with a rate λ , a half-saturation constant 1 r 1 , and handling time h 1 . Furthermore, it is assumed that immature pests mature at a constant rate α , contributing to the growth of both susceptible and infectious mature pest populations.

The rate of change of susceptible mature pests increases with the rate α b , while the rate of change of infectious mature pests increases with the rate α ( 1 − b ) . Additionally, susceptible pest population becomes infectious with the transmission rate γ and dies out with disease mortality rate γ . Moreover, it is assumed that the susceptible population undergoes natural or manual harvesting at a constant proportion h of the existing mature susceptible pest population, recognizing the potential threat this population poses to crops.

Based on these assumptions, our mathematical model consists of the following differential equations:

(1) d x j d t = r x m s − α x j − d 1 x j , d x m s d t = α b x j − β x m s y − θ x m s x m I − d 2 x m s − h x m s , d x m I d t = α ( 1 − b ) x j + θ x m s x m I − ( d 2 + γ ) x m I , d y d t = λ β x m s y 1 + r 1 h 1 x m s − d 3 y ,

where d 1 , d 2 , and d 3 represent, respectively, the death rate of immature pests, mature pests, and natural enemies. All parameters are assumed positive, and system (1) has non-negative initial conditions.

To determine the stability analysis, the Routh–Hurwitz stability criterion is restated in the following theorem:

Theorem 2.1

Consider the fourth-degree polynomial with real constant coefficients:

(2) P ( λ ) = λ 4 + a 1 λ 3 + a 2 λ 2 + a 3 λ + a 4 .

All roots of polynomial (2) are negative or have a negative real part if the determinant of the Hurwitz matrix H 4 is positive, i.e., Det ( H 4 ) > 0 , where 4-Hurwitz matrix is given as:

H 4 = a 1 1 0 0 a 3 a 2 a 1 1 0 a 4 a 3 a 2 0 0 0 a 4 .

3 Analysis of the system

Due to the natural process, no species can be negative; hence, in the next, we will prove the positivity of the solution of system (1).

3.1 Positivity for Model (1)

From system (1), we obtain that

d x j d t x j = 0 = r x m s ,

d x m s d t x m s = 0 = α b x j ,

d x m I d t x m I = 0 = α ( 1 − b ) x j ,

d y d t y = 0 = 0 .

Here, all the rates are non-negative on bounding planes of the non-negative cone of R 4 . Thus, if we begin in the interior of this cone, we shall always remain in this cone as the direction of a vector field is inward on all the bounding planes. Hence, all the solutions of the model system (1) are non-negative.

Theorem 3.1

All solutions of system (1) that initiate in R 4 are bounded and enter into a region Γ defined by

Γ = ( x j , x m s , x m I , y ) ∈ R 4 : 0 ≤ x j ( t ) , x m s ( t ) , x m I ( t ) ≤ ( r − d 2 ) 2 4 r 1 d 2 h 1 μ , 0 ≤ y ( t ) ≤ λ ( r − d 2 ) 2 4 r 1 d 2 h 1 μ .

Proof

Assume ( x j ( t ) , x m s ( t ) , x m I ( t ) , y ( t ) ) be any solution of system (1).

Define V ( t ) = x j ( t ) + x m s ( t ) + x m I ( t ) + 1 λ y ( t ) with μ = min ( d 1 , h , d 2 + γ , d 3 λ ) .

The time derivative along a solution of (1) is obtained as:

d N d t = r x m s − α x j − d 1 x j + α b x j − β x m s y − θ x m s x m I − d 2 x m s − h x m s + α ( 1 − b ) x j + θ x m s x m I − ( d 2 + γ ) x m I + 1 λ λ β x m s y 1 + r 1 h 1 x m s − d 3 y ≤ r x m s − μ N − d 2 x m s − d 2 r 1 h 1 x m s 2 .

It follows that

d N d t + μ N ≤ r x m s − d 2 x m s − d 2 r 1 h 1 x m s 2 ≤ M 0 ,

where

M 0 = ( r − d 2 ) 2 4 r 1 d 2 h 1 .

Thus, for the initial value N ( 0 ) = N 0 , the comparison theorem [5] gives that limsup t → ∞ N ( t ) ≤ ( r − d 2 ) 2 4 r 1 d 2 h 1 μ . Hence, Γ is positively invariant and bounded.□

3.2 Existence and uniqueness of solution for Model (1)

Clearly, the interaction functions of system (1) are continuous and have continuous partial derivatives on the following positive four-dimensional space

Λ = { ( x j , x m s , x m I , y ) ∈ R 4 : x j ( 0 ) ≥ 0 , x m s ( 0 ) ≥ 0 , x m I ( 0 ) ≥ 0 , y ( 0 ) ≥ 0 } .

So, these functions are Lipschitzian on R 4 , and hence, the solution of system (1) exists and is unique.

3.3 Existence of equilibrium

Pest natural enemy epidemic models often reveal the existence of dynamic equilibria where the populations of pests and their natural enemies stabilize over time. A balance between pest reproduction and predation, parasitism, or disease transmission characterizes this equilibrium. In this section, we will determine all the equilibrium points of system (1). Clearly system (1) has three equilibrium points as follows:

The trivial equilibrium point P 0 ( x j , x m s , x m I , y ) = (0, 0, 0, 0) always exists.
Natural enemy-free equilibrium point P 1 ( x j ¯ , x m s ¯ , x m I ¯ , 0 ) exists, provided α r b > ( d 2 + h ) ( α + d 1 ) , where
x j ¯ = r ( d 2 + γ ) ( α r b − ( d 2 + h ) ( α + d 1 ) ) ( r α ( 1 − b ) + α r b − ( d 2 + h ) ( α + d 1 ) ) θ ( α + d 1 ) ,

x m s ¯ = ( d 2 + γ ) ( α r b − ( d 2 + h ) ( α + d 1 ) ) ( r α ( 1 − b ) + α r b − ( d 2 + h ) ( α + d 1 ) ) θ ,

x m I ¯ = α r b − ( d 2 + h ) ( α + d 1 ) θ ( α + d 1 ) .
The interior equilibrium point P 2 ( x j * , x m s * , x m I * , y * ) exists, provided that
α b r λ β ( d 2 + γ ) + ( d 2 + h ) ( α + d 1 ) ( ( d 2 + γ ) d 3 r 1 h 1 + θ d 3 ) > α b r ( d 2 + γ ) d 3 r 1 h 1 + α b r θ d 3 + θ r α d 3 ( 1 − b ) ( α + d 1 ) + λ β ( d 2 + h ) ( α + d 1 ) ( d 2 + γ ) ,
where
x j * = r d 3 ( α + d 1 ) ( λ β − d 3 r 1 h 1 ) , x m s * = d 3 λ β − d 3 r 1 h 1 ,

x m I * = r α ( 1 − b ) d 3 λ β ( d 2 + γ ) − ( ( d 2 + γ ) d 3 r 1 h 1 + θ d 3 ) ,

y * = M β ( α + d 1 ) ( λ β ( d 2 + γ ) − ( ( d 2 + γ ) d 3 r 1 h 1 + θ d − 3 ) ) ,
with
M = α b r λ β ( d 2 + γ ) + ( d 2 + h ) ( α + d 1 ) − ( α b r ( d 2 + γ ) d 3 r 1 h 1 ( ( d 2 + γ ) d 3 r 1 h 1 + θ d 3 ) + α b r θ d 3 + θ r α d 3 ( 1 − b ) ( α + d 1 ) + λ β ( d 2 + h ) ( α + d 1 ) ( d 2 + γ ) ) .

4 Local stability analysis

Theorem 4.1

The trivial equilibrium point is locally asymptotically stable if ( α + d 1 ) ( d 2 + h ) > α r b .

Proof

The transition matrix for system (1) has been computed as:

(3) J = − ( α + d 1 ) r 0 0 α b − β y − θ x m I − d 2 − h − θ x m s 0 α ( 1 − b ) θ x m I θ x m s − ( d 2 + γ ) 0 0 λ β y ( 1 + r 1 h 1 x m s ) 2 0 λ β x m s 1 + r 1 h 1 x m s − d 3 .

The transition matrix (3) at trivial equilibrium P₀ is given as

(4) J [ P 0 ] = − ( α + d 1 ) r 0 0 α b − ( d 2 + h ) 0 0 α ( 1 − b ) 0 − ( d 2 + γ ) 0 0 0 0 − d 3 .

Clearly, equation (4) has two eigenvalues − d 3 , and − ( d 2 + γ ) , and the rest of the eigenvalues are given by a quadratic equation Λ 2 + A 1 Λ + A 2 = 0 , where

A 1 = α + d 1 + d 2 + h , A 2 = ( α + d 1 ) ( d 2 + h ) − α r b .

It reveals that A 1 > 0 , A 2 > 0 for ( α + d 1 ) ( d 2 + h ) > α r b . Hence, the Routh–Hurwitz criterion follows that the trivial equilibrium point is locally asymptotically stable for ( α + d 1 ) ( d 2 + h ) > α r b and unstable if ( α + d 1 ) ( d 2 + h ) < α r b .□

Theorem 4.2

The equilibrium point in the absence of natural enemies P 1 ( x j ¯ , x m s ¯ , x m I ¯ , 0 ) is locally asymptotically stable, provided that ( α + d 1 ) ( d 2 + h ) < α r and x m s ¯ < d 3 ( λ β − d 3 r 1 h 1 ) .

Proof

The transition matrix for system (1) at the equilibrium point in the absence of natural enemies is given as:

J [ P 1 ] = − ( α + d 1 ) r 0 0 α b − r α b ( α + d 1 ) − θ x m s − β x m s α ( 1 − b ) θ x m I − α ( 1 − b ) x j x m I 0 0 0 0 λ β x m s 1 + r 1 h 1 x m s − d 3 .

Clearly, the aforementioned matrix has a eigenvalue λ β x m s 1 + r 1 h 1 x m s ¯ − d 3 , which is negative, provided that x m s ¯ < d 3 ( λ β − d 3 r 1 h 1 ) . The rest of the eigenvalues can be calculated through a cubic equation:

Λ 3 + B 1 Λ 2 + B 2 Λ + B 3 = 0 ,

where

B 1 = α + d 1 + r α b α + d 1 + α ( 1 − b ) x j ¯ x m I ¯ ,

B 2 = α ( α + d 1 ) ( 1 − b ) + r α 2 b ( 1 − b ) α + d 1 x j ¯ x m I ¯ + θ 2 x m I ¯ x m s ¯ ,

B 3 = ( d 2 + γ ) ( α r b − ( d 2 + h ) ( α + d 1 ) ) 2 ( ( α + d 1 ) ( d 2 + h ) − r α ) + r α ( 1 − b ) ( d 2 + γ ) ( α r b − ( α + d 1 ) ( d 2 + h ) ) ( ( α + d 1 ) ( d 2 + h ) − r α ) .

It is easy to obtain that for ( α + d 1 ) ( d 2 + h ) < α r , B 1 > 0 , B 2 > 0 . Hence, the Routh–Hurwitz criterion follows that the natural enemy-free equilibrium point is locally asymptotically stable for ( α + d 1 ) ( d 2 + h ) < α r and x m s ¯ < d 3 ( λ β − d 3 r 1 h 1 ) .□

Theorem 4.3

The interior equilibrium point P 2 ( x j * , x m s * , x m I * , y * ) of system (1) is locally asymptotically stable, if it exists.

Proof

The transition matrix for system (1) at the interior equilibrium point is given as:

J [ P 2 ] = − ( α + d 1 ) r 0 0 α b − β y * − θ x m I * − d 2 − h − θ x m s * − β x m s * α ( 1 − b ) θ x m I * θ x m s * − d 2 − γ 0 0 λ β y * ( 1 + r 1 h 1 x m s * ) 2 0 λ β x m s * 1 + r 1 h 1 x m s * − d 3 .

Now, using the interior equilibrium equation, we obtain

J [ P 2 ] = − ( α + d 1 ) r 0 0 α b − α b x j * x m s * − θ x m s * − β x m s * α ( 1 − b ) θ x m I * − α ( 1 − b ) x j * x m I * 0 0 λ β y * ( 1 + r 1 h 1 x m s * ) 2 0 0 .

The characteristic equation of the matrix J [ P 2 ] is

(5) Λ 4 + C 1 Λ 3 + C 2 Λ 2 + C 3 Λ + C 4 = 0 ,

where

C 1 = α + d 1 + α x j * x m s * ,

C 2 = α r ( 1 − b ) + α 2 b ( 1 − b ) r ( α + d 1 ) + θ 2 x m s * x m I * + β 2 λ x m s * y * ( 1 + r 1 h 1 x m s * ) 2 ,

C 3 = ( α + d 1 ) θ 2 x m s * x m I * + ( α + d 1 ) λ β 2 x m s * y * ( 1 + r 1 h 1 x m s * ) 2 + α ( 1 − b ) λ β 2 x m s * y * ( 1 + r 1 h 1 x m s * ) 2 x m s * + r α θ ( 1 − b ) x m s * + r α θ ( 1 − b ) x m s * ,

C 4 = α ( 1 − b ) ( α + d 1 ) β 2 λ x j * y * ( 1 + r 1 h 1 x m s * ) 2 .

By the Routh–Hurwitz criterion, all the roots of the equation (5) have real parts less than zero if and only if C i > 0 , i = 1 , 2 , 3 , 4 and ( C 1 C 2 − C 3 ) C 3 > C 1 2 C 4 . Through direct calculations, it shows that C i > 0 , i = 1 , 2 , 3 , 4 whenever

α b r λ β ( d 2 + γ ) + ( d 2 + h ) ( α + d 1 ) ( ( d 2 + γ ) d 3 r 1 h 1 + θ d 3 ) > α b r ( d 2 + γ ) d 3 r 1 h 1 + α b r θ d 3 + θ r α d 3 ( 1 − b ) ( α + d 1 ) + λ β ( d 2 + h ) ( α + d 1 ) ( d 2 + γ ) .

Therefore, the interior equilibrium point of system (1) is always locally asymptotically stable if it exists.□

5 Global stability analysis

Theorem 5.1

The trivial equilibrium point is globally asymptotically stable, provided that r d 2 + h < λ < 1 .

Proof

Consider the following function:

V 0 ( x j , x m s , x m I , y ) = x j + λ x m s + λ x m I + y .

Clearly, V 0 : R + 4 → R is a C 1 positive definite function. Then, by differentiating V 0 for time t and some algebraic manipulation, we obtain

d V 0 d t = d x j d t + λ d I m s d t + λ d I m I d t + d y d t = r x m s − ( α + d 1 − λ α ) x j − λ β x m s y − ( λ d 2 + λ h ) x m s − λ ( d 2 + γ ) x m I + λ β x m s y 1 + r 1 h 1 x m s − d 3 y , ≤ r x m s − ( α ( 1 − λ ) + d 1 ) x j − λ β x m s y − λ d 2 x m s − λ h x m s − λ ( d 2 + γ ) x m I + λ β x m s y − d 3 y , = − ( α ( 1 − λ ) + d 1 ) x j − ( λ d 2 + λ h − r ) x m s − λ ( d 2 + γ ) x m I − d 3 .

It is evident that d V 0 d t < 0 for r d 2 + h < λ < 1 . Hence, the trivial equilibrium point is globally asymptotically stable for r d 2 + h < λ < 1 .□

Theorem 5.2

The interior equilibrium point of system (1) is globally asymptotically stable if the following conditions hold:

(6) r > max { d 1 , β + d 2 + h − α b } , α ( 1 − b ) > d 2 + γ , λ > 1 , a n d β λ > d 3 .

Proof

Consider the following function:

V 3 ( x j , x m s , x m I , y ) = x j − x j * − x j * ln x j x j * + x m s − x m s * − x m s * ln x m s x m s * + x m I − x m I * − x m I * ln x m I x m I * + y − y * − y * ln y y * .

Here, it can be easily obtained that V 3 : R + 4 → R is a C 1 positive definite function. Then, by differentiating V 3 for time t and some algebraic manipulation, we obtain

d V 3 d t = ( x j − x j * ) x j ( r x m s − α x j − d 1 x j ) + ( x m s − x m s * ) x m s ( α b x j − β x m s y − θ x m s x m I − d 2 x m s − h x m s ) + ( x m I − x m I * ) x m I ( α ( 1 − b ) x j + θ x m s x m I − ( d 2 + γ ) x m I ) + ( y − y * ) y λ β x m s y 1 + r 1 h 1 x m s − d 3 y , = ( x j − x j * ) x j ( r x m s − α x j − d 1 x j − r x m s * + α x j * + d 1 x j * ) + ( x m s − x m s * ) x m s ( α b x j − β x m s y − θ x m s x m I − d 2 x m s − h x m s − α b x j * + β x m s * y * + θ x m s * x m I * + d 2 x m s * + h x m s * ) + ( x m I − x m I * ) x m I ( α ( 1 − b ) x j + θ x m s x m I − ( d 2 + γ ) x m I − α ( 1 − b ) x j * − θ x m s * x m I * + ( d 2 + γ ) x m I * ) + ( y − y * ) y λ β x m s y 1 + r 1 h 1 x m s − d 3 y − λ β x m s * y * 1 + r 1 h 1 x m s * + d 3 y * , = ( x j − x j * ) x j ( r ( x m s − x m s * ) − ( α + d 1 ) ( x j − x j * ) ) + ( x m s − x m s * ) x m s ( α b ( x j − x j * ) − β ( x m s y − x m s * y * ) − θ ( x m s x m I − x m s * x m I * ) − ( d 2 + h ) ( x m s − x m s * ) ) + ( x m I − x m I * ) x m I ( α ( 1 − b ) ( x j − x j * ) + θ ( x m s x m I − x m s * x m I * ) − ( d 2 + γ ) ( x m I − x m I * ) ) + ( y − y * ) y λ β ( x m s y 1 + r 1 h 1 x m s − x m s * y * 1 + r 1 h 1 x m s * ) − d 3 ( y − y * ) , ≤ − ( r − d 1 ) ( x j − x j * ) 2 − ( r + α b − β − d 2 − h ) ( x m s − x m s * ) 2 − ( α ( 1 − b ) − d 2 − γ ) ( x m I − x m I * ) 2 − ( λ β − d 3 ) ( y − y * ) 2 − β ( λ − 1 ) ( x m s y − x m s * y * ) 2 .

It is clear that d V 3 d t < 0 , if r > max { d 1 , β + d 2 + h − α b } , α ( 1 − b ) > d 2 + γ , λ > 1 , and β λ > d 3 . Hence, the interior equilibrium point of system (1) is globally asymptotically stable under condition (6).□

6 Optimal control problem

In this section, we develop an optimal control problem to show the impact of optimal harvesting of pest population using some pest-controlling agents on system (1). The aim of applying these controls to the model is to keep the number of susceptible mature pest populations at a minimum level under the low cost of implementation related to pest-controlling agents. Hence, this method will help lessen the susceptible mature pest population so that the cost incurred is minimal. Let us consider the following cost functional:

(7) J [ u ( t ) ] = ∫ 0 T [ ω 1 x m s + ω 2 u 2 ( t ) ] d t ,

where u ( t ) is the control variable, a Lebesgue measurable function for implementing pest-controlling agents considered such that 0 ≤ u ≤ 1 on a finite interval [ 0 , T ] , and T represents the termination time of the control strategy. We denote ω 1 x m s as the cost related to the susceptible mature population, encompassing expenses incurred by controlling agents. The cost associated with harvesting is expressed as ω 2 u 2 ( t ) , covering various expenditures such as additional harvesting costs. Here, additional harvesting cost refers to the extra expenses incurred in harvesting a susceptible mature pest population of interest as part of a control strategy. Harvesting, in this context, likely involves activities such as capturing, collecting, or otherwise managing the susceptible mature pest population, and the associated costs could include expenses for labor, equipment, materials, or any other resources required beyond the baseline costs. Including the quadratic nonlinearity u 2 ( t ) in the harvesting, cost is explored for its relevance and similarity to previous works [1,12]. The analysis focuses on a control problem involving discussed policies and aims to minimize the associated cost function.

The cost functional (7) will be minimized subjected to the following system of equations by considering the harvesting intervention in the second equation of Model (1) through pest controlling agent:

(8) d x j d t = r x m s − α x j − d 1 x j , d x m s d t = α b x j − β x m s y − θ x m s x m I − d 2 x m s − u ( t ) x m s , d x m I d t = α ( 1 − b ) x j + θ x m s x m I − ( d 2 + γ ) x m I , d y d t = λ β x m s y 1 + r 1 h 1 x m s − d 3 y ,

with initial conditions:

x j ( 0 ) ≥ 0 , x m s ( 0 ) ≥ 0 , x m I ( 0 ) ≥ 0 , y ( 0 ) ≥ 0 .

In the functional J , the positive weight constants ω 1 and ω 2 represent the proportionality constants for transforming the integral into currency expanded. The first sum in the cost functional represents the incurred cost due to susceptible mature pest populations, and the last term is the cost due to pest controlling agents. Our objective is to find an optimal control u * ( t ) such that

(9) J ( u * ) = min u ∈ Z J ( u ) ,

where the control set is defined as:

Z = { u ( t ) : 0 ≤ u ( t ) ≤ u max ; 0 ≤ t ≤ T , u ( t ) is Lebesgue measurable } .

6.1 Existence of optimal control

In this subsection, we guarantee the existence of optimal control u * in Z , which minimizes the cost-functional by using the result given in the study by Fleming and Rishel [16].

Theorem 6.1

There exists an optimal control u * ∈ Z such that J ( u * ) = min [ J ( u ) ] corresponding to the control systems (7) and (8).

Proof

The boundedness of the solution of system (8) asserts the existence of a solution to control the system using results [32]. Therefore, the set of controls and corresponding state variables is nonempty, and the control set is closed and convex by definition. Hence, the solution of system (7) is bounded above by a linear function in state and control. The integrand in the cost functional

ω 1 x m s ( t ) + ω 2 u 2 ( t )

is convex on control set Z , and also, there exists q 1 , q 2 > 0 and m > 1 such that,

ω 1 x m s ( t ) + ω 2 u 2 ( t ) ≥ q 1 + q 2 ∣ u 2 ( t ) ∣ m ,

where q 1 depends upon the upper bound of x m s ( t ) , and q 2 with q 2 = ω 2 .

Thus, the existence of an optimal control is established.□

6.2 Characterization of optimal control function

Pontryagin’s maximum principle [4] can be used for the differential systems of adjoint variables and characterization of optimal control. We define the Hamiltonian as:

(10) H ( x j , x m s , x m I , y , u , λ ) = L ( x j , x m s , x m I , y , u ) + λ 1 x j ˙ + λ 2 x m s ˙ + λ 3 x m I ˙ + λ 4 y ˙ = ω 1 x m s ( t ) + ω 2 u ( t ) 2 + λ 1 ( r x m s − α x j − d 1 x j ) + λ 2 ( α b x j − β x m s y − θ x m s x m I − d 2 x m s − u ( t ) x m s ) + λ 3 ( α ( 1 − b ) x j + θ x m s x m I − ( d 2 + γ ) x m I ) + λ 4 λ β x m s y 1 + r 1 h 1 x m s − d 3 y ,

where λ = ( λ 1 , λ 2 , λ 3 , λ 4 ) is known as an adjoint variable. Here, the cost-functional is minimized subjected to state variables of the model and the optimal control u * . Thus, it reveals that there exists an adjoint variable λ i satisfying the following canonical equations:

(11) d λ 1 d t = ( α + d 1 ) λ 1 − α b λ 2 − α ( 1 − b ) λ 3 , d λ 2 d t = − ω 1 − r λ 1 + ( β y + θ x m I + d 2 + u ( t ) ) λ 2 − θ x m I λ 3 − λ β y ( 1 + r 1 h 1 x m s ) 2 λ 4 , d λ 3 d t = θ x m s λ 2 − ( θ x m s − d 2 − γ ) λ 3 , d λ 4 d t = β x m s λ 2 − λ β x m s 1 + r 1 h 1 x m s − d 3 λ 4 ,

with transversality condition λ i ( T ) = 0 , i = 1 to 4 . Here, the Hamiltonian is minimized with respect to u at the optimal value u * . Now, from the optimality condition, we have ∂ H ∂ u = 0 at u = u * , which provides that

H = ω 2 u ( t ) 2 − u ( t ) x m s λ 2 + terms without u(t) .

Differentiation of H with respect to u gives that

∂ H ∂ u = 2 ω 2 u * − x m s * λ 2 .

Clearly, u * = λ 2 x m s * 2 ω 2 . Now, from the aforementioned findings along with the characteristics of control set Z , we have

u * = 0 , if λ 2 x m s * 2 ω 2 < 0 , λ 2 x m s * 2 ω 2 , if 0 ≤ λ 2 x m s * 2 ω 2 ≤ 1 , 1 , if λ 2 x m s * 2 ω 2 > 1 .

6.3 Optimal system

In this section, we collectively write the corresponding optimal system using the abovementioned optimal control. The optimal system with minimized Hamiltonian H * at ( x j * , x m s * , x m I * , y * , u * , λ i ) is given as:

(12) d x j d t = r x m s * − α x j * − d 1 x j * , d x m s d t = α b x j * − β x m s * y * − θ x m s * x m I * − d 2 x m s * − u * x m s * , d x m I d t = α ( 1 − b ) x j * + θ x m s * x m I * − ( d 2 + γ ) x m I * , d y d t = λ β x m s * y * 1 + r 1 h 1 x m s * − d 3 y * ,

with non-negative initial conditions. The corresponding adjoint system is

(13) d λ 1 d t = ( α + d 1 ) λ 1 − α b λ 2 − α ( 1 − b ) λ 3 , d λ 2 d t = − ω 1 − r λ 1 + ( β y * + θ x m I * + d 2 + u * ) λ 2 − θ x m I * λ 3 − λ β y * ( 1 + r 1 h 1 x m s * ) 2 λ 4 , d λ 3 d t = θ x m s * λ 2 − ( θ x m s * − d 2 − γ ) λ 3 , d λ 4 d t = β x m s * λ 2 − λ β x m s * 1 + r 1 h 1 x m s * − d 3 λ 4 ,

with transversality conditions λ 1 ( T ) = 0, λ 2 ( T ) = 0, λ 3 ( T ) = 0, and λ 4 ( T ) = 0 and u * as similar.

7 Numerical simulations and discussions

A pest natural enemy epidemic model is a mathematical or conceptual framework used to study the dynamics of a pest population and its interaction with its natural enemies (predators, parasites, pathogens, etc.). In this section, to confirm the feasibility of our analysis regarding existence and stability conditions, we will provide some numerical computations with the help of MATLAB using a set of parametric values.

7.1 Stability analysis

Let us assume that r = 2 , α = 0.5 , b = 0.5 , θ = 0.3 , h = 0.2 , γ = 0.02 , λ = 1.01 , h 1 = 0.4 , r 1 = 0.32 , and d 3 = 0.1 , with β = 0.12 , d 1 = 0.35 , and d 2 = 0.39 . In our system, as described by equations (1), it is observed that an equilibrium point exists in the absence of natural enemies. This equilibrium point is demonstrated in Figure 1, and further analysis reveals that it is locally asymptotically stable. The presence of this stable equilibrium underscores the dynamics of the system when natural enemies are not considered, providing valuable insights into the ecological behavior of the pest population. Once again, utilizing the same parameter values, specifically β = 0.5 , d 1 = 0.02 , and d 2 = 0.1 , it is observed that the system described by equations (1) exhibits the existence of an interior equilibrium point. This equilibrium point is demonstrated in Figure 2, and a detailed analysis reveals that it is locally asymptotically stable. This finding provides valuable insights into the system’s dynamics under specific parameter configurations, offering a nuanced understanding of the ecological interactions in the presence of natural enemies.

Figure 1

Local stability of natural enemy-free equilibrium point.

Figure 2

Local stability of interior equilibrium point.

Furthermore, to delve into the nonlinear stability behavior of the interior equilibrium within the specified parametric values, we consider four distinct initial conditions: Y 1 = (3.254, 2.2568, 4.25, 4.3658), Y 2 = (1.254, 2.256, 3.25, 3.3658), Y 3 = (1.254, 1.02568, 4.725, 4.93658), and Y 4 = (0.954, 0.82568, 2.925, 2.93658). The trajectories of the system, as described by equations (1), have been plotted for these initial conditions. The analysis of these solution trajectories reveals that the interior equilibrium is globally asymptotically stable. Refer to Figures 3 and 4 for a visual representation of the system’s solution trajectories, emphasizing the global stability of the interior equilibrium under diverse initial conditions.

$Figure 3 Phase plot of x j ‒ x m s ‒ x m I {x}_{j}‒{x}_{ms}‒{x}_{mI} showing the convergence to interior equilibrium point.$

Figure 3

Phase plot of x j ‒ x m s ‒ x m I showing the convergence to interior equilibrium point.

$Figure 4 Phase plot of x j ‒ x m s ‒ y {x}_{j}‒{x}_{ms}‒y showing the convergence to interior equilibrium point.$

Figure 4

Phase plot of x j ‒ x m s ‒ y showing the convergence to interior equilibrium point.

7.2 Sensitivity analysis

In this section, we conducted a sensitivity analysis of the model to understand the impact of various parameters on the compartments. Figures 5, 6, 7, 8 specifically examine the effect of the parameter θ on the compartments while maintaining the values of the other parameters constant. Notably, compartments y and x m I appear to be more significantly influenced by this parameter. Furthermore, Figures 9, 10, 11, 12 investigate the effect of parameter h on the compartments, with the other parameters held constant. In this case, compartment y emerges as more prominently affected by the parameter. These sensitivity analyses contribute valuable insights into how changes in specific parameters influence the dynamics of the model compartments, aiding a comprehensive understanding of the system’s behavior.

$Figure 5 Effect of parameters θ \theta on immature pest.$

Figure 5

Effect of parameters θ on immature pest.

$Figure 6 Effect of parameters θ \theta on susceptible mature pest.$

Figure 6

Effect of parameters θ on susceptible mature pest.

$Figure 7 Effect of parameters θ \theta on infectious mature pest.$

Figure 7

Effect of parameters θ on infectious mature pest.

$Figure 8 Effect of parameters θ \theta on natural enemies.$

Figure 8

Effect of parameters θ on natural enemies.

Figure 9

Effect of parameters h on immature pest.

Figure 10

Effect of parameters h on susceptible mature pest.

Figure 11

Effect of parameters h on infectious mature pest.

Figure 12

Effect of parameters h on natural enemies.

Here, it is observed that the natural enemy population plays a very crucial role in the susceptible mature pest population. The effects of the predation rate β , transmission rate θ , and harvesting rate h on the susceptible adult pest population are shown in Figures 13, 14 and 15. The results obtained from the figures indicate that an increase in predation rate leads to a reduction in the susceptible pest population. This finding aligns with ecological expectations, as higher predation rates contribute to more effective natural control of pest populations. Similarly, an elevated transmission rate is associated with a decreased susceptible pest population. This outcome underscores the importance of understanding the mechanisms behind disease transmission, as a higher transmission rate results in a more rapid decline in the pest population. Furthermore, an increase in the harvesting rate is found to have a similar effect on the susceptible pest population, leading to a reduction. Harvesting becomes more effective as a control measure when the rate of removal of pests from the ecosystem is increased. The maturation rate α also significantly affects susceptible mature pests. The effect of α on the susceptible adult pest population is shown in Figure 16. Here, it is observed that if the α is increased, the number of susceptible mature pest populations also increases.

$Figure 13 Effect of variation in β \beta on susceptible mature pest.$

Figure 13

Effect of variation in β on susceptible mature pest.

$Figure 14 Effect of variation in θ \theta on susceptible mature pest.$

Figure 14

Effect of variation in θ on susceptible mature pest.

Figure 15

Effect of variation in h on susceptible mature pest.

$Figure 16 Effect of variation in α \alpha on susceptible mature pest.$

Figure 16

Effect of variation in α on susceptible mature pest.

7.3 Optimal control problem

We employ the forward-backward sweep method to determine the optimal solution, solving the optimality system numerically with the specified parameter values. The process involves solving the system of state variables forward in time, followed by solving the system of adjoint variables backwards in time. The Runge–Kutta fourth-order iterative scheme is employed for these computations. In each iteration, the value of the control is updated, and this iterative process continues until the desired convergence is attained. This numerical approach allows us to efficiently navigate the complex dynamics of the system, providing insights into the optimal control strategy under the given set of parameters. We start with x j = 2.5096 , x m s = 5.5656 , x m I = 4.7488 , and y = 1.2330 . For the simulation of the optimal control system, let us assume that r = 2 , α = 0.5 , b = 0.5 , θ = 0.3 , h = 0.2 , γ = 0.02 , λ = 1.01 , h 1 = 0.4 , r 1 = 0.32 , d 3 = 0.1 , β = 0.5 , d 1 = 0.02 , and d 2 = 0.1 . In our optimization process, we select the weight factors linked to the susceptible mature pest population as follows: ω 1 = 1.23 , ω 2 = 1 , and u max = 1 . These weight factors play a crucial role in influencing the optimal control strategies. Subsequently, we conduct simulations for the populations of x j , x m s , x m I , and y under the influence of the time-dependent optimal control u ( t ) . Additionally, we compare these simulations with scenarios where no control strategies are implemented. The results of these simulations are depicted in Figure 17, providing a visual representation of the population dynamics under the influence of the time-dependent optimal control and the absence of control strategies.

Figure 17

Profiles of populations with applied optimal control.

8 Conclusion

Pest natural enemy epidemic models serve as indispensable tools in understanding the intricate dynamics between pests and their natural adversaries, providing a foundation for developing sustainable and environmentally friendly pest management strategies. These models are of paramount importance for ecologists, pest management professionals, and policymakers striving to devise effective measures for mitigating the impact of pests on agriculture and ecosystems. This study delves into the dynamic behavior of optimal harvesting within an eco-epidemic model featuring stage-structured dynamics. The model reveals three equilibrium points: the trivial point, the natural enemy-free point, and an interior equilibrium point. Local stability analyses are conducted for each equilibrium point, supplemented by a comprehensive examination of global stability at the trivial and interior equilibrium points. Utilizing Lyapunov’s stability theory, we derive sufficient conditions for the global stability behavior of the interior equilibrium, shedding light on the overarching system dynamics. Furthermore, an optimal control problem is explored using Pontryagin’s minimum principle to minimize the cost associated with harvesting interventions. This optimization approach contributes to more efficient and economically viable strategies for pest management.

In this study, it becomes evident that the population dynamics of natural enemies play a pivotal role in influencing the susceptible mature pest population. The intricate interplay between these components underscores the importance of considering ecological nuances for successfully implementing pest management strategies. Overall, this research contributes to the theoretical understanding of eco-epidemic models and offers practical insights for developing effective and cost-efficient pest control strategies in real-world scenarios.

ksmathur1709@gmail.com, bhagwaniisc@gmail.com

Acknowledgements

The authors are grateful to the anonymous reviewers for their valuable suggestions, which helped significantly improve the manuscript.

Funding information: The work of Bhagwan Kumar is financially supported by the Council of Scientific andIndustrial Research, India (Grant No.: 09/0150(16285)/2023-EMR-I).
Author contributions: Both authors contributed equally to the work.
Conflict of interest: The authors have no conflict of interest to disclose.
Ethical approval: This research did not require ethical approval.
Data availability statement: Data sharing is not applicable to this article as no data sets were generated or analyzed during the current study.

References

[1] Amalia, R. D., Arif, D. K., Windarto, W., & Fatmawati, F. (2018). Optimal control of predator-prey mathematical model with infection and harvesting on prey. Journal of Physics: Conference Series, 974, 012050. IOP Publishing. 10.1088/1742-6596/974/1/012050Search in Google Scholar

[2] Ávila-Vales, E., Estrella-González, Á., & Esquivel, E. R. (2017). Bifurcations of a Leslie Gower predator prey model with Holling type III functional response and Michaelis-Menten prey harvesting. arXiv: http://arXiv.org/abs/arXiv:1711.08081. Search in Google Scholar

[3] Bellier, E., Sæther, B.-E., & Engen, S. (2021). Sustainable strategies for harvesting predators and prey in a fluctuating environment. Ecological Modelling, 440, 109350. 10.1016/j.ecolmodel.2020.109350Search in Google Scholar

[4] Pontryagin, L. S. (2018). Mathematical theory of optimal processes. London: Routledge. 10.1201/9780203749319Search in Google Scholar

[5] Birkhoff, G. & Rota, G. (1982). Ordinary differential equation, Boston, New York: Ginn. and Co.Search in Google Scholar

[6] Cao, Q., Chen, G., & Yang, W. (2023). The impact of fear effect on the dynamics of a delayed predator-prey model with stage structure. International Journal of Biomathematics, 16(8), 2250139. 10.1142/S179352452250139XSearch in Google Scholar

[7] Chakraborty, K., Chakraborty, M., & Kar, T. K. (2011). Optimal control of harvest and bifurcation of a prey–predator model with stage structure. Applied Mathematics and Computation, 217(21), 8778–8792. 10.1016/j.amc.2011.03.139Search in Google Scholar

[8] Charles, R., Makinde, O. D., & Kungaro, M. (2022). A review of the mathematical models for the impact of seasonal weather variation and infections on prey predator interactions in serengeti ecosystem. Open Journal of Ecology, 12(11), 718–732. 10.4236/oje.2022.1211041Search in Google Scholar

[9] Chattopadhyay, J. & Arino, O. (1999). A predator-prey model with disease in the prey. Nonlinear analysis, 36, 747–766. 10.1016/S0362-546X(98)00126-6Search in Google Scholar

[10] Comins, H. N. & Blatt, D. W. (1974). Prey-predator models in spatially heterogeneous environments. Journal of Theoretical Biology, 48(1), 75–83. 10.1016/0022-5193(74)90180-5Search in Google Scholar PubMed

[11] Das, A. & Samanta, G. (2020). A prey–predator model with refuge for prey and additional food for predator in a fluctuating environment. Physica A: Statistical Mechanics and its Applications, 538, 122844. 10.1016/j.physa.2019.122844Search in Google Scholar

[12] Das, H. & Shaikh, A. A. (2021). Dynamical response of an eco-epidemiological system with harvesting. Journal of Applied Mathematics and Computing, 65, 67–91. 10.1007/s12190-020-01379-8Search in Google Scholar

[13] Dawes, J. & Souza, M. (2013). A derivation of Holling’s type i, ii and iii functional responses in predator-prey systems. Journal of theoretical biology, 327, 11–22. 10.1016/j.jtbi.2013.02.017Search in Google Scholar PubMed

[14] Dubey, B. & Kumar, A. (2019). Dynamics of prey–predator model with stage structure in prey including maturation and gestation delays. Nonlinear Dynamics, 96(4), 2653–2679. 10.1007/s11071-019-04951-5Search in Google Scholar

[15] Firdiansyah, A. L. (2021). Effect of prey refuge and harvesting on dynamics of eco-epidemiological model with Holling type III. Jambura Journal of Mathematics, 3(1), 16–25. 10.34312/jjom.v3i1.7281Search in Google Scholar

[16] Fleming, W. H. & Rishel, R. W. (2012). Deterministic and stochastic optimal control (Vol. 1). New York: Springer-Verlag. Search in Google Scholar

[17] Freedman, H. I. (1980). Deterministic mathematical models in population ecology. New York: Marcel Dekker, Inc. Search in Google Scholar

[18] Gallé, R., Császár, P., Makra, T., Gallé-Szpisjak, N., Ladányi, Z., Torma, A., Ingle, K., & Szilassi, P. (2018). Small-scale agricultural landscapes promote spider and ground beetle densities by offering suitable overwintering sites. Landscape Ecology, 33, 1435–1446. 10.1007/s10980-018-0677-1Search in Google Scholar

[19] Haque, M. (2010). A predator-prey model with disease in the predator species only. Nonlinear Analysis: Real World Applications, 11(4), 2224–2236. 10.1016/j.nonrwa.2009.06.012Search in Google Scholar

[20] Hethcote, H. W., Wang, W., Han, L., & Ma, Z. (2004). A predator-prey model with infected prey. Theoretical Population Biology, 66(3), 259–268. 10.1016/j.tpb.2004.06.010Search in Google Scholar PubMed

[21] Huang, Y., Chen, F., & Zhong, L. (2006). Stability analysis of a prey–predator model with Holling type iii response function incorporating a prey refuge. Applied Mathematics and Computation, 182(1), 672–683. 10.1016/j.amc.2006.04.030Search in Google Scholar

[22] Jana, S. & Kar, T. (2013). A mathematical study of a prey–predator model in relevance to pest control. Nonlinear Dynamics, 74, 667–683. 10.1007/s11071-013-0996-3Search in Google Scholar

[23] Jatav, K. S. & Dhar, J. Global behavior and Hopf bifurcation of a stage-structured prey–predator model with maturation delay for prey and gestation delay for predator. Journal of Biological systems, 23(1), 57–77. 10.1142/S0218339015500047Search in Google Scholar

[24] Kar, T. K. (2005). Stability analysis of a prey–predator model incorporating a prey refuge. Communications in Nonlinear Science and Numerical Simulation, 10(6), 681–691. 10.1016/j.cnsns.2003.08.006Search in Google Scholar

[25] Kermack, W. O. & McKendrick, A. G. (1927). A contribution to the mathematical theory of epidemics. Proceedings of the royal society of London. Series A, Containing papers of a mathematical and physical character, 115(772), 700–721. 10.1098/rspa.1927.0118Search in Google Scholar

[26] Kim, T. N., Bukhman, Y. V., Jusino, M. A., Scully, E. D., Spiesman, B. J., & Gratton, C. (2022). Using high-throughput amplicon sequencing to determine diet of generalist lady beetles in agricultural landscapes. Biological Control, 170, 104920. 10.1016/j.biocontrol.2022.104920Search in Google Scholar

[27] Lee, J. & Baek, H. (2017). Dynamics of a Beddington-Deangelis type predator-prey system with constant rate harvesting. Electronic Journal of Qualitative Theory of Differential Equations, 2017(1), 1–20. 10.14232/ejqtde.2017.1.1Search in Google Scholar

[28] Li, F. & Li, H. (2012). Hopf bifurcation of a predator-prey model with time delay and stage structure for the prey. Mathematical and Computer Modelling, 55(3–4), 672–679. 10.1016/j.mcm.2011.08.041Search in Google Scholar

[29] Li, W., Chen, Y., Huang, L., & Wang, J. (2022). Global dynamics of a Filippo predator-prey model with two thresholds for integrated pest management. Chaos, Solitons & Fractals, 157, 111881. 10.1016/j.chaos.2022.111881Search in Google Scholar

[30] Liu, Q., Jiang, D., Hayat, T., Alsaedi, A., & Ahmad, B. (2020). Dynamical behavior of a stochastic predator-prey model with stage structure for prey. Stochastic Analysis and Applications, 38(4), 647–667. 10.1080/07362994.2019.1710188Search in Google Scholar

[31] Lotka, A. J. (1925). Elements of physical biology. Philadelphia: Williams & Wilkins. Search in Google Scholar

[32] Lukes, D. L. (1982). Differential equations: Classical to controlled. University of Virginia, Charlottesville, Virginia: Elsevier. Search in Google Scholar

[33] Ma, Z., Yang, J., & Jiang, G. (2010). Impulsive control in a stage structure population model with birth pulses. Applied Mathematics and Computation, 217(7), 3453–3460. 10.1016/j.amc.2010.09.012Search in Google Scholar

[34] Mathur, K. S. (2016). A prey-dependent consumption two-prey one predator eco-epidemic model concerning biological and chemical controls at different pulses. Journal of the Franklin Institute, 353(15), 3897–3919. 10.1016/j.jfranklin.2016.07.012Search in Google Scholar

[35] May, R. M. (1972). Limit cycles in predator-prey communities. Science, 177(4052), 900–902. 10.1126/science.177.4052.900Search in Google Scholar PubMed

[36] Mbava, W., Mugisha, J., & Gonsalves, J. W. (2017). Prey, predator and super-predator model with disease in the super-predator. Applied Mathematics and Computation, 297, 92–114. 10.1016/j.amc.2016.10.034Search in Google Scholar

[37] Mortoja, S. G., Panja, P., & Mondal, S. K. (2018). Dynamics of a predator-prey model with stage-structure on both species and anti-predator behavior. Informatics in Medicine Unlocked, 10, 50–57. 10.1016/j.imu.2017.12.004Search in Google Scholar

[38] Mugala, T., Brichler, K., Clark, B., Powell, G. S., Taylor, S., & Crossley, M. S. (2023). Ground beetles suppress slugs in corn and soybean under conservation agriculture. Environmental Entomology, xx(xx), 1–9. 10.1093/ee/nvad047Search in Google Scholar PubMed

[39] Panja, P. & Mondal, S. K. (2015). Stability analysis of coexistence of three species prey–predator model. Nonlinear Dynamics, 81, 373–382. 10.1007/s11071-015-1997-1Search in Google Scholar

[40] Pearce, C. (1970). A new deterministic model for the interaction between predator and prey. Biometrics, 26(3), 387–392. 10.2307/2529095Search in Google Scholar

[41] Prakash, D. B. & Vamsi, D. D. Stochastic time-optimal control and sensitivity studies for additional food provided prey–predator systems involving Holling type-iv functional response. Frontiers in Applied Mathematics and Statistics, 9, 1122107. 10.3389/fams.2023.1122107Search in Google Scholar

[42] Saha, S. & Samanta, G. (2020). A prey–predator system with disease in prey and cooperative hunting strategy in predator. Journal of Physics A: Mathematical and Theoretical, 53(48), 485601. 10.1088/1751-8121/abbc7bSearch in Google Scholar

[43] Sen, M., Banerjee, M., & Morozov, A. (2012). Bifurcation analysis of a ratio-dependent prey–predator model with the allee effect. Ecological Complexity, 11, 12–27. 10.1016/j.ecocom.2012.01.002Search in Google Scholar

[44] Sih, A. (1987). Prey refuges and predator-prey stability. Theoretical Population Biology, 31(1), 1–12. 10.1016/0040-5809(87)90019-0Search in Google Scholar

[45] Smith, R. & Mead, R. (1974). Age structure and stability in models of prey–predator systems. Theoretical Population Biology, 6(3), 308–322. 10.1016/0040-5809(74)90014-8Search in Google Scholar PubMed

[46] Srinivasu, P., Vamsi, D., & Ananth, V. (2018). Additional food supplements as a tool for biological conservation of predator-prey systems involving type iii functional response: A qualitative and quantitative investigation. Journal of Theoretical Biology, 455, 303–318. 10.1016/j.jtbi.2018.07.019Search in Google Scholar PubMed

[47] Sugie, J., Kohno, R., & Miyazaki, R. (1997). On a predator-prey system of Holling type. Proceedings of the American Mathematical Society, 125(7), 2041–2050. 10.1090/S0002-9939-97-03901-4Search in Google Scholar

[48] Wang, X. & Zou, X. (2017). Modeling the fear effect in predator-prey interactions with adaptive avoidance of predators. Bulletin of Mathematical Biology, 79, 1325–1359. 10.1007/s11538-017-0287-0Search in Google Scholar PubMed

[49] Wangersky, P. J. and Cunningham, W. (1957). Time lag in prey–predator population models. Ecology, 38(1), 136–139. 10.2307/1932137Search in Google Scholar

[50] Wikan, A. & Kristensen, O. (2021). Compensatory and overcompensatory dynamics in prey–predator systems exposed to harvest. Journal of Applied Mathematics and Computing, 67(1), 455–479. 10.1007/s12190-020-01484-8Search in Google Scholar

[51] Yang, J., Tang, S., & Cheke, R. A. (2013). Global stability and sliding bifurcations of a non-smooth Gause predator-prey system. Applied Mathematics and Computation, 224, 9–20. 10.1016/j.amc.2013.08.024Search in Google Scholar

[52] Yu, X., Zhu, Z., Lai, L., & Chen, F. (2020). Stability and bifurcation analysis in a single-species stage structure system with Michaelis-Menten-type harvesting. Advances in Difference Equations, 2020(1), 1–18. 10.1186/s13662-020-02652-7Search in Google Scholar

[53] Zhang, H., Cai, Y., Fu, S., & Wang, W. (2019). Impact of the fear effect in a prey–predator model incorporating a prey refuge. Applied Mathematics and Computation, 356, 328–337. 10.1016/j.amc.2019.03.034Search in Google Scholar

Received: 2023-08-31

Revised: 2023-12-27

Accepted: 2023-12-30

Published Online: 2024-03-09

This work is licensed under the Creative Commons Attribution 4.0 International License.

Optimal control of susceptible mature pest concerning disease-induced pest-natural enemy system with cost-effectiveness

Abstract

1 Introduction

2 Model development and proposition

Theorem 2.1

3 Analysis of the system

3.1 Positivity for Model (1)

Theorem 3.1

Proof

3.2 Existence and uniqueness of solution for Model (1)

3.3 Existence of equilibrium

4 Local stability analysis

Theorem 4.1

Proof

Theorem 4.2

Proof

Theorem 4.3

Proof

5 Global stability analysis

Theorem 5.1

Proof

Theorem 5.2

Proof

6 Optimal control problem

6.1 Existence of optimal control

Theorem 6.1

Proof

6.2 Characterization of optimal control function

6.3 Optimal system

7 Numerical simulations and discussions

7.1 Stability analysis

7.2 Sensitivity analysis

7.3 Optimal control problem

8 Conclusion

Acknowledgements

References

Journal and Issue

Articles in the same Issue