Computer Systems Science & Engineering DOI:10.32604/csse.2022.019655 | [picture]
Article |

Deterministic and Stochastic Fractional Order Model for Lesser Date Moth

Moustafa El-shahed¹,* and Asmaa M. Al-Dububan²

1 Department of Mathematics, Unaizah College of Sciences and Arts, Qassim University, P.O. Box 3771, Unaizah, 51911, Saudi Arabia
2 Department of Mathematics, Buraydah College of Sciences and Arts, Qassim University, P.O. Box 1162, Buraydah, 51431, Saudi Arabia
Corresponding Author: Moustafa El-shahed. Email: elshahedm@yahoo.com
Received: 21 April 2021; Accepted: 23 May 2021

Abstract: In this paper, a deterministic and stochastic fractional order model for lesser date moth (LDM) using mating disruption and natural enemies is proposed and analysed. The interaction between LDM larvae, fertilized LDM female, unfertilized LDM female, LDM male and the natural enemy is investigated. In order to clarify the characteristics of the proposed deterministic fractional order model, the analysis of existence, uniqueness, non-negativity and boundedness of the solutions of the proposed fractional-order model are examined. In addition, some sufficient conditions are obtained to ensure the local and global stability of equilibrium points. The occurrence of local bifurcation near the equilibrium points is investigated with the help of Sotomayor’s theorem. Numerical simulations are conducted to illustrate the properties of the proposed fractional order model with respect to the intrinsic growth rate of the LDM larvae, natural enemy’s mortality rate, predation rate, sex pheromone trap parameter, fractional order and environmental noise. The impact of mating disruption on lesser date moth is demonstrated. Also, a numerical approximation method is developed for the proposed stochastic fractional-order model.

Keywords: Lesser date moth; stochastic; stability; natural enemies; Sotomayor’s theorem; mating disruption

1  Introduction

The palm tree is considered one of the oldest major and basic crops in Southwest Asia, North Africa and many other places of the world. Palm trees are affected by many agricultural pests that cause significant losses to the palm trees and their fruits, as well as affecting the age and growth of the palm tree if it is left without control [1]. The LDM is one of the most dangerous pests of young palm trees and immature palm fruits. The damage is mainly due to the way in which the LDM larvae feed, as soon as they leave the egg until the date of their entry the virgin dwelling develops feeding and causing tunnels in the affected part of the palm in all directions and depths without any early signs showing the infection [2–5]. One of the most promising strategies for controlling LDM is the use of a mating disruption using the sex pheromone traps [5]. Natural enemies’ use to stop pest infestation has long been recommended as a clean and environmentally friendly way to protect crops. The main natural enemies that are used in agricultural pest control are larval predators. Goniozus swirskiana can be considered as one of the most important natural enemies that can attack the LDM [6,7]. In the real world, plant diseases models are always affected by the environmental noise. Thus, the stochastic models may be a more appropriate way of modelling agricultural pests in many circumstances [8]. Recently, fractional calculus has been applied to describe different mathematical models, and it has been shown to be more accurate in some cases compared to the classical models [9]. The main objective of this paper is to propose and analyse a deterministic and stochastic fractional order LDM model with mating disruption and sex pheromone traps taking into consideration the effect of the natural enemy on LDM.

The paper is arranged as follows: In Section 2, the deterministic mathematical model is described as well as the existence, uniqueness, non-negativity, and boundedness of the solutions of LDM system are verified. The local and global stability of equilibrium points of the LDM system is analyzed in Section 3. Using Sotomayor’s theorem, the local bifurcation conditions are derived in Section 4. In Section 5, we extend the deterministic fractional-order LDM model to the stochastic case and a numerical approximation method developed for the proposed stochastic fractional-order case. Some numerical examples are given in Section 6 to illustrate the theoretical findings. Finally, the discussion and conclusion are given in Section 7.

2  Dynamic of the Deterministic Fractional Order Model

Following [10–14], the model of lesser date moth with mating disruption and sex pheromone

trap is describing by the following system

cαDαx=rx(1−xk)−α1x−βxza+x−μ1x,cαDαy=ϵα1x−α2y+δf−μ2y,cαDαf=α2y−δf−μ3f,cαDαm=(1−ϵ)α1x−γpmy+p−μ4m,cαDαz=eβxza+x−μ5z, (1)

where CDα is the Caputo fractional derivative of order α and 0<α<1 . The detailed explanation of the variables and parameters for system (1) are listed in Tab. 1.

2.1 Existence and Uniqueness

In this section, we investigate the existence and uniqueness of the solutions of the fractional order system (1) in the region Ω×(0,T] where

Ω={(x,y,f,m,z)∈R+5:max(|x|,|y|,|f|,|m|,|z|)≤φ},

Theorem 1. For each X0=(x0,y0,f0,m0,z0)∈Ω , there exists a unique solution X(t)∈Ω of the fractional order system (1), which is defined for all t≥0 .

Proof. Define a mapping F(X)=(F1(X),F2(X),F3(X),F4(X),F5(X)) , in which

F1(X)=rx(1−xk)−α1x−βxza+x−μ1x,F2(X)=ϵα1x−α2y+δf−μ2y,F3(X)=α2y−δf−μ3f,F4(X)=(1−ϵ)α1x−γpmy+p−μ4m,F5(X)=eβxza+x−μ5z. (2)

For any X,X¯∈Ω, it follows from (2) that

where

Hence, F(X) satisfies the Lipschitz condition with respect to X . According to Cresson et al. [15], as F(X) locally Lipschitz. Then there exists unique local solution to the fractional order system (1).

2.2 Non-Negativity and Boundedness

The following results show the non-negativity of the solutions of the fractional order system (1). According to Cresson et al. [15], a model of the form dXdt=F(X) satisfies the positivity property if and only if for all i=1,….,5 , Fi(X)≥0 for all X∈R+5 such that Xi=0 . Thus, the solution of the integer-order model (1), with nonnegative initial conditions remains nonnegative. Also, the solution satisfies the Lipschitz condition, as stated in Theorem 1. By Theorem 5 and Theorem 6 in Cresson et al. [15], the solution of the fractional-order model (1) also satisfies the non-negativity. The boundedness of the solutions of model (1) are given in the following theorem.

Theorem 2. All the solutions of the fractional-order LDM Model (1) starting in R+5 are uniformly bounded.

Proof. The approach of [16,17] is utilized. Let (x(t),y(t),f(t),m(t),z(t)) to be any solution of the system (1) with non-negative initial conditions. Let M(t)=x(t)+y(t)+f(t)+m(t)+z(t), then

cαDαM(t)≤rx(1−xk)−μ1x−μ2y−μ3f−μ4m−μ5z−γpmy+p≤−rk(x−k2)2+rk4−νM,

where, ν<min{μ1,μ2,μ3,μ4,μ5} , thus, cαDαM(t)+νM≤rk4 . In accordance with Lemma 9 in Choi et al. [18], it follows that, 0≤M(t)≤M(0)Eα(−νtα)+rk4νtαEα,α+1(−νtα),

where Eα is the Mittag-Leffler function. According to Lemma 5 and Corollary 6 in Choi et al. [18], it follows

0≤M(t)≤rk4ν,ast→∞.

Hence all the solutions of fractional-order LDM model (1) that start in R+5 are uniformly bounded in the region, H={(x,y,f,m,z)∈R+5:M(t)≤rk4ν+ξ,for\ anyξ>0} .

3  Equilibria and Stability

The LDM model (1) has the following three equilibrium points:

1) E0=(0,0,0,0,0) , which always exists.

2) The free natural enemy equilibrium point E1=(x1,y1,f1,m1,0) , where

x1=k(1−1R0),y1=ϵ(δ+μ3)(1−1R0)δμ2+μ3(α2+μ2),f1=ϵα2(1−1R0)δμ2+μ3(α2+μ2), m1=(1−ϵ)α1k(y1+p)γp+μ4(y1+p)(1−1R0).

The free natural enemy equilibrium point exists positively if R0>1, where, R0=rα1+μ1 is the basic offspring number obtained by using the next generation method [19].

3) The coexistence equilibrium point E2=(x2,y2,f2,m2,z2) , where

x2=aμ5β−μ5,y2=ϵα1(δ+μ3)x2δμ2+μ3(α2+μ2),f2=α1α2ϵx2δμ2+μ3(α2+μ2), m2=(1−ϵ)α1x2(y2+p)γp+μ4(y2+p),z2=ea(α1+μ1)(eβ−μ5)(R0−1−arμ5k(α1+μ1)(eβ−μ5)).

4) The coexistence equilibrium point E2 exists if eβ>μ5 and R0>1+raμ5k(α1+μ1)(eβ−μ5).

The locally and globally asymptotically stable of equilibrium points of LDM system (1) are now investigated. The stability analysis of the equilibrium point E0=(0,0,0,0,0) is not considered because in this case all the population will go to extinction.

The stability of free natural enemy equilibrium point E1=(x1,y1,f1,m1,0) is investigated as follows.

Theorem 3. If R0<1+raμ5k(α1+μ1)(eβ−μ5), then the equilibrium point E1 is locally asymptotically stable.

Proof. The first three eigenvalues of J(E1) are λ1=r−2rx1k−(α1+μ1)=(α1+μ1)(1−R0) , λ2=−(γpy1+p+μ4) and λ3=eβx1a+x1−μ5 . The other roots are determined by

λ2+(α2+δ+μ2+μ3)λ+δμ2+μ3(α2+μ2)=0. (3)

The roots of Eq. (3) have negative real parts. It can be observed that λ3<0, when eβx1a+x1<μ5 which equivalent R0<1+raμ5k(α1+μ1)(eβ−μ5). So, E1 is locally asymptotically stable if 1<R0<1+raμ5k(α1+μ1)(eβ−μ5).

Theorem 4. If R0<1+raμ5k(α1+μ1)(eβ−μ5) , then the natural enemy extinction equilibrium point E1 is globally asymptotically stable.

Proof. The following positive definite Lyapunov function is considered.

L1=eaa+x1(x−x1−x1ln⁡xx1)+z.

By calculating the time derivative of L1 along the solution of system (1), one obtains,

cαDαL1=ea(x−x1)a+x1[r(1−xk)−βza+x−(α1+μ1)]+[eβxa+x−μ5]z≤−ear(x−x1)2k(a+x1)+(eβx1a+x1−μ5)z.

Then cαDαL1≤0 , when eβx1a+x1<μ5 which equivalent R0<1+raμ5k(α1+μ1)(eβ−μ5). According to generalized Lyapunov–Lasalle’s invariance principle [20], the natural enemy extinction equilibrium point E1 is globally asymptotically stable when R0<1+raμ5k(α1+μ1)(eβ−μ5).

The stability of the coexistence equilibrium point E2=(x2,y2,f2,m2,z2) is investigated as follows.

The first eigenvalue of J(E2) is λ1=−γpp+y2−μ4 . The other eigenvalues are determined by

λ4+θ1λ3+θ2λ2+θ3λ+θ4=0, (4)

where

θ1=α2+δ+μ2+μ3+B11,θ2=α2μ3+eaβ2z2x2(a+x2)3+B11(α2+δ+μ2+μ3)+μ2(δ+μ3),θ3=α2(B11μ3(a+x2)3+eaβ2z2x2)+μ2(B11(a+x2)3(δ+μ3)+eaβ2x2z2)+eaβ2x2z2(δ+μ3)(a+x2)3,θ4=eaβ2x2z2(α2μ3+μ2(δ+μ3))(a+x2)3,

where B11=rx2k−βz2x2(a+x2)2 and B44=γpp+y2+μ4 . Then, the proposition proposed in Matouk [21] can be used to determine the stability conditions of the equilibrium point E2 . When B11>0 , then θi>0,i=1,2,3,4. Also, θ1θ2θ3−θ32−θ12θ4=(B11(a+x2)3Ω4+Ω2Ω3)(B11δμ2(a+x2)3+Ω1(α2+μ2)+eaβ2w2x2(δ+μ3))2(a+x2)9, where

Ω1=B11μ3(a+x2)3+eaβ2z2x2,Ω2=(a+x2)3−eaβ2z2x2,Ω3=α2μ3+μ2(δ+μ3)−1,Ω4=(α2+δ+μ2+μ3).

when Ω2Ω3>0 , then θ1θ2θ3−θ32−θ12θ4>0 , therefore all the eigenvalues of the Jacobian matrix J(E2) near the equilibrium point E2 have negative real parts. Thus, due to the Routh-Hurwitz criterion the equilibrium point E2 is locally asymptotically stable. The local stability of the coexistence equilibrium point E2 is given in the following theorem.

Theorem 5. If Ω2Ω3>0 and B11>0 , then the coexistence equilibrium point E2 is locally asymptotically stable.

The global stability of the coexistence equilibrium point is investigated in the following theorem.

Theorem 6. If βz2a(a+x2)<rk , then the coexistence equilibrium point E2 is globally asymptotically stable.

Proof. The following positive definite Lyapunov function is considered.

L2=H(x−x2−x2ln⁡xx2)+z−z2−z2ln⁡zz2.

By calculating the time derivative of L2 along the solution of system (1), one obtains,

cαDαL2=H(x−x2)[r(1−xk)−βza+x−(α1+μ1)]+(z−z2)[eβxa+x−μ5]≤H[βz2a(a+x2)−rk]+β(x−x2)(z−z2)(a+x)(a+x2)(ea−H(a+x2)).

Choosing H=eaa+x2, then cαDαL2<0 . According to generalized Lyapunov–Lasalle’s invariance principle [20], the coexistence equilibrium point E2 is globally asymptotically stable when βz2a(a+x2)<rk .

4  Bifurcation Analysis

In this section the local bifurcations near the equilibrium points of LDM model (1) are inves-tigated with the help of Sotomayor’s theorem [22]. The Hopf bifurcation theorem given in Liu [23] is also presented to discuss the bifurcation analysis of the underlying system. One can compute

D2F(x,φ)(U,U)=((2aβz(a+x)3−2rk)ζ12−2aβζ1ζ5(a+x)200−2γpmζ22(p+y)3+2γpζ2ζ4(p+y)2eaβζ1ζ5(a+x)2−2eaβzζ12(a+x)3), (5)

where φ is any bifurcation parameter and U=(ζ1,ζ2,ζ3,ζ4,ζ5)T is any eigenvector.

Theorem 7. The LDM system (1) undergoes a transcritical bifurcation with respect to the bifurcation parameter r around E0 if r=r∗=α1+μ1 .

Proof. Let V1=(ν1,ν2,ν3,ν4,ν5)T be the eigenvector corresponding to the zero eigenvalue of the matrix (E0) r∗=α1+μ1 , hence J(E0)V1=0 , gives

V1=(ν1ϵα1(δ+μ3)ν1α2μ3+μ2(δ+μ3)ϵα1α2ν1α2μ3+μ2(δ+μ3)(1−ϵ)α1ν1γ+μ40),

where ν1 is any non zero real number. Similarly, suppose V2=(τ1,τ2,τ3,τ4,τ5)T be the eigenvector corresponding to the zero eigenvalue of the matrices J(E0) , thus J(E0)TV2=0 gives V2=(1,0,0,0,0)T . Consider, ∂F∂r=Fr(X,r)=(x(1−xk),0,0,0,0)T , thus, V2TFr(E0,r∗)=0. Therefore, according to Sotomayor’s theorem for local bifurcation, the LDM model (1) has no saddle-node bifurcation near E0 at r∗=α1+μ1 . One can note that r=α1+μ1 is equivalent to R0=1 . Now,

DFr(E0,r∗)=(1000000000000000000000000),

then V2TDFr(E0,r∗)V1=ν1≠0. Using (5), one obtains

V2TD2F(X,r)(V1,V1)=−2rkν12−2βaν1ν5≠0.

Thus, according to Sotomayor’s theorem, the LDM system (1) has a transcritical bifurcation at r∗=α1+μ1 as the parameter r passes through the value r∗ , thus the proof is complete.

Theorem 8. The LDM system (1) undergoes a transcritical bifurcation with respect to the bifurcation parameter μ5 around E1=(x1,y1,f1,m1,0) if μ5=μ5∗=eβx1a+x1 .

Proof. The Jacobian matrix of the LDM system (1) at the free enemy equilibrium point E1 with μ5=μ5∗ has zero eigenvalue takes the form

J(E1,μ5∗)=(A11000−βx1a+x1α1ϵ−A22δ000α2−A3300(1−ϵ)α1γpm1(p+y1)20−A44000000),

where A11=r−2rx1k−(α1+μ1),A22=(α2+μ2),A33=(δ+μ3) and A44=γpy1+p+μ4. J(E1)V3=0 , gives the eigenvector corresponding to the zero eigenvalue of the matrix J(E1) , hence

V3=(μ5ν5eA11ϵA33α1μ5ν5eA11(A22A33−δα2)ϵα1α2μ5ν5eA11(A22A33−δα2)α1(((1−ϵ)A22(y1+p)2+ϵγpm1)A33−(1−ϵ)δ(y1+p)2α2)μ5ν5e(y1+p)2A11A44(A22A33−δα2)ν5),

where ν5 is any non zero real number. Similarly, J(E1)TV4=0 gives the eigenvector corresponding to the zero eigenvalue of the matrix J(E1)TV4=0 , hence V4=(0,0,0,0,1)T. Consider, ∂F∂μ5=Fμ5(X,μ5)=(0,0,0,0,−z)T , thus, V4TFμ5(E1,μ5∗)=0. Therefore, according to Sotomayor’s theorem for local bifurcation, the LDM model (1) has no saddle-node bifurcation near E1 at μ5=μ5∗ . Now,

DFμ5(E1,μ5∗)=(000000000000000000000000−1),

then, V4TDFμ5(E1,μ5∗)V3=−ν5≠0. Using (5), one obtains

V4TD2F(X,μ5)(V3,V3)=eaβν1ν5(a+x1)2≠0.

Thus, according to Sotomayor’s theorem, the LDM system (1) has a transcritical bifurcation at μ5∗=eβx1a+x1 as the parameter μ5 passes through the value μ5∗ , thus the proof is complete.

In this part, we shall show that as the coexistence equilibrium loses stability, periodic solutions can bifurcate from the positive equilibrium. We first give the following lemma.

Lemma 9 . The characteristic Eq.(4) has a pair of purely imaginary roots and the remaining roots have negative real parts if and only if z2=r(a+x2)2kβ and θ1θ2θ3−θ32−θ12θ4=0.

Suppose (4) has two eigenvalues which have negative real parts and two complex conjugates eigenvalues (call them λ=m(φ)±in(φ) ) such that m(φ∗)=0,n(φ∗)>0, dmdφ|φ=φ∗≠0. Substituting λ=m(φ)±in(φ) into (4), and separating the real and imaginary, we get

m4+θ1m3+θ2m3+θ3m+θ4−(6m2+3θ1m+θ2)n2+n4=0, (6)

4m3+3θ1m2+2θ2m+θ3−(4m+θ1)n2=0, (7)

Substituting (6) into (7), differentiating with respect to φ and utilizing m(φ∗)=0 and n(φ∗)≠0, we have

dmdφ=[ddφ(θ1θ2θ3−θ32−θ12θ4)2θ1(4θ4−θ1θ3−θ22)]φ=φ∗≠0.

Theorem 10. For the coexistence equilibrium point E2 of the integer order LDM system (1), the system around E2 enters into the Hopf bifurcation when φ passes φ∗ if the coefficients θj(φ)(j=1,2,3,4) at φ=φ∗ satisfying the following condition:

1.     Φ(φ∗)=[θ1(φ)θ2(φ)θ3(φ)−θ32(φ)−θ12(φ)θ4(φ)]|φ=φ∗=0,

2.     (4θ4−θ1θ3−θ22)|φ=φ∗≠0 ,

3.     dΦ(φ)dφ|φ=φ∗≠0 .

According to Theorem 10, there exists a Hopf bifurcation in the LDM model (1) with α=1 , where the Hopf bifurcation is controlled by φ .

5  Stochastic Fractional Order Model

This section extends the deterministic fractional-order LDM model (1) to the following stochastic fractional-order model.

CαDαx=rx(1−xk)−α1x−βxza+x−μ1x+σ1xdW1dt,CαDαy=ϵα1x−α2y+δf−μ2y+σ2ydW2dt,CαDαf=α2y−δf−μ3f+σ3fdW3dt,CαDαm=(1−ϵ)α1x−γpmy+p−μ4m+σ4mdW4dt,CαDαz=eβxza+x−μ5z+σ5zdW5dt, (8)

where Wi(i=1,2,3,4,5) are independent standard Brownian motions with Wi(0)=0 and σi>0 denote the intensities of the white noise. The stochastic fractional-order LDM model (8) can be written in the general form:

CαDαX(t)=F(X)+g(X)dWdt, (9)

where F(x) is given in (2), g(x)=(σ1x,σ2y,σ3f,σ4m,σ5z) and dWdt=(dW1dt,dW2dt,dW3dt,dW4dt,dW5dt)T. Applying Riemann–Liouville integral to both sides of (9), one can obtain the following stochastic Volterra integral equation.

X(t)=X0+∫0tF(X)(t−s)α−1Γ(α)ds+∫0tg(X)(t−s)α−1dW(s)Γ(α)ds. (10)

According to Wang et al. [24,25], under some conditions on the coefficient functions, the global existence and uniqueness of solutions for the stochastic fractional-order system (8) can be investigated. Because Grunwald-Letnikov’s definition is the most straightforward from the point of view of numerical implementation, so we will use it to solve the LDM system of fractional order stochastic differential equations. Grunwald-Letnikov ( GLDα ) fractional derivative of order α defined by Aminikhah et al. [26,27]

GLαDαf(t)=Limh→0⁡h−α∑j=0[t−ah](−1)j(αj)f(t−jh), (11)

where [t−ah] means the integer part of t−ah . This formula (11) can be reduced to

GLαDαf(tn)≈h−α∑j=0nwjαf(tn−j), (12)

where h is the time step, tn=nh and wjα are the Grunwald-Letnikov coefficients satisfy the following recurrence relationship

w0α=1,wjα=(1−1+αj)wj−1α,j=1,2,3,…

If f(t) is continuous function and f′(t) is integrable function in the interval [0,T] , then the relation between Caputo and Grunwald-Letnikov fractional derivative takes the form [28,29]

CαDαf(t)=GLαDαf(t)−f(0)t−αΓ(1−α)=GLαDαf(t)−Limn→∞(−1)nhα(α−1n)f(0)≈1hα∑j=0nwjα(f(tn−j)−f(0)). (13)

Now, the fractional order stochastic LDM model (8) can be written as

xn=x0+hα(rxn−1(1−xn−1k)−(α1+μ1)xn−1−βxn−1zn−1a+xn−1+σ1xn−1hζ1n)−∑j=1nwjα(xn−j−x0)yn=y0+hα(ϵα1xn−(α2+μ2)yn−1+δfn−1+σ2yn−1hζ2n)−∑j=1nwjα(yn−j−y0)fn=f0+hα(α2yn−(δ+μ3)fn−1+σ3fn−1hζ3n)−∑j=1nwjα(fn−j−f0)mn=m0+hα((1−ϵ)α1xn−γpmn−1yn+p−μ4mn−1+σ4mn−1hζ4n)−∑j=1nwjα(mn−j−m0)zn=z0+hα(eβxnzn−1a+xn−μ5zn−1+σ5zn−1hζ5n)−∑j=1nwjα(zn−j−z0), (14)

where, σi and ζin represent real constants and a 5D Gaussian white noise processes, respectively, i=1,2,3,4,5. ζi satisfy the follows:

⟨ζj(t)⟩=0,(j=0,1,2,3,4,5)and⟨ζi(t1)ζj(t2)⟩=δijδ(t1−tj),

δij is Kronecker delta and δ(t1−tj) is the Dirac delta function.

6  Numerical Simulations

In this section, we simulate the fractional-order LDM system (1) and stochastic fractional-order LDM (8) by the following parameters:

r=1.8,ϵ=0.6,k=30,μ1=0.01,μ2=0.01,μ3=0.01,μ4=0.01,μ5=0.01, α1=0.5,α2=0.5,δ=0.5,a=15,γ=0.4,p=0.2,e=0.9,β=0.1.

To show the effect of the intrinsic growth rate of the LDM larvae, we draw the bifurcation diagram concerning r as a bifurcation parameter. It can be seen that a transcritical bifurcation occurs at r = 0.51 as shown in Fig. 1 and stated in Theorem 7. It can also be observed that when r > 0.51, the coexistence equilibrium point E2 = (1.6667, 25.2475, 24.7525, 25.36, 157.407r − 85) is locally asymptotically stable. According to Theorem 10, it can be seen that the supercritical Hopf bifurcation value localized at r = 1.31143. For r > 1.31143 the LDM system (1) undergoes limit cycle behaviour.

To show the effect of the natural enemy’s mortality rate around the coexistence equilibrium points, we draw the bifurcation diagram for µ5 as a bifurcation parameter. It can be seen that the Hopf bifurcation value localized at µ5 = 0.0284692 as shown in Fig. 2. It can also be observed that when µ5 > 0.0284692, the coexistence equilibrium point E2 is locally asymptotically stable. For µ5 < 0.3235, the system undergoes limit cycle behaviour. According to Theorem 7, for µ5 > 0.0642346 the natural enemy goes extinct from the system and a transcritical bifurcation occurs at µ5 = 0.0642346 and the equilibrium point E1 = (26.94, 408.101, 400.009, 528.446, 0) is locally asymptotically stable.

In order to show the effect of the predation rate, we draw the bifurcation diagram with respect to β as a bifurcation parameter. It can be seen that a transcritical bifurcation occurs at β = 0.0501342 as shown in Fig. 3. It can also be observed that when β < 0.0501342, the free enemy equilibrium point E1 = (22.35, 338.569, 331.93, 436.688, 0) is locally asymptotically stable as indicated in Theorem 8. For 0.0501342 < β < 0.152449, the coexistence equilibrium point E2 is locally asymptotically stable. It can be seen that the supercritical Hopf bifurcation value localized at β = 0.152449. For β > 0.152449 the LDM system (1) undergoes limit cycle behaviour.

From Fig. 4a, it can be seen that the sex pheromone trap parameter p is important in that it affects the population density of the LDM male. One can observe from Fig. 4b. that the population density of LDM male decrease with increasing p. We conclude that the dynamics of LDM can be controlled by sex pheromone trap parameters p.

Fig. 5 shows that the deterministic fractional-order remains stable for different values of fractional-order α though solutions reach to equilibrium point E2(7.5, 113.614, 111.386, 140.149, 189) more slowly for a smaller value of fractional-order α. It is important to notice that when α = 1 the fractional order model for lesser date moth (1) reduces to the classical integer-order model [10–14].

When the strength of environmental noise is very close to zero, the fractional-order stochastic LDM system (8) behaves like the deterministic model. Fig. 6a indicates the dynamical behavior of stochastic LDM (8) without noise (i.e., σi = 0), which gives the deterministic model (1). In Fig. 6a, the dynamical behavior of the LDM system was unstable, which corresponds with the results of Theorem 5. Fig. 6b shows that the LDM smooth oscillations’ dynamical behavior when the strength of the noise was low ( σi = 0.02). However, with an increase in the strength of noise, such as the medium-noise situation shown in Fig. 6c, the dynamical behavior of the LDM became more complex, and they tended to extinction. Fig. 6c represents the dynamical behavior of the model (1) when the noise strength was high ( σi = 0.05). The natural enemy can die out due to the white noise stochastic disturbance. By comparing Figs. 6a and 6b, one can realize that if the noise is not strong, the stochastic perturbation does not cause sharp changes in the LDM model (8). However, when the environmental noise σi is sufficiently large (see Fig. 6c), the noise can force the natural enemy to become extinct.

Figure 1: Bifurcation diagram of LDM model (1) with respect to r

Figure 2: Bifurcation diagram of LDM model (1) with respect to µ5

Figure 3: Bifurcation diagram of LDM model (1) with respect to β

Figure 4: (a) Plot of the LDM male versus time with p = 0.1, 0.3, 0.8, (b) Plot of the LDM male versus time with 0 ≤ p ≤ 1

Figure 5: Time series of the LDM female with α = 0.94, 0.97, 1

Figure 6: Time series of the fractional order stochastic model (10) with σ = 0, 0.02, 0.05 (a) σi = 0 (b) σi = 0.02 (c) σi = 0.05

7  Discussion and Conclusion

In this paper, we consider a deterministic and stochastic fractional order model for LDM using mating disruption and natural enemies. We obtain some sufficient conditions that ensure the local and global stability of equilibrium points. We conclude that sex pheromone trap parameters can control dynamics of LDM. The occurrence of local bifurcation near the equilibrium point is performed using Sotomayor’s theorem. Numerical simulations are performed to support and illustrate the theoretical findings. From the numerical results, one can realize that if the noise is not strong, the stochastic perturbation does not cause sharp changes in the dynamics of the LDM model. However, when the environmental noise is sufficiently large, the noise can force the population to become extinct.

Funding Statement: The authors gratefully acknowledge Qassim University, represented by the Deanship of Scientific Research, on the financial support under the number (cosao-bs-2019-2-2-I-5469) during the academic year 1440 AH/2019 AD.

Conflicts of Interest: The authors declare that they have no conflicts of interest to report regarding the present study.

References

 1.  L. I. El-Juhany, “Degradation of date palm trees and date production in Arab Countries: Causes and potential rehabilitation,” Australian Journal of Basic and Applied Sciences, vol. 4, no. 8, pp. 3998–4010, 2010.

 2.  J. Al-Khayri, “Date palm phoenix dactylifera micropropagation,” in Protocols for Micropropagation of Woody Trees and Fruits, Springer, pp. 509–526, 2007.

 3.  A. Levi-Zada, D. Fefer, L. Anshelevitch, A. Litovsky, M. Bengtsson et al., “Identification of the sex pheromone of the lesser date moth, Batrachedra amydraula, using sequential SPME auto-sampling,” Tetrahedron Letters, vol. 52, no. 35, pp. 4550–4553, 2011.

 4.  M. A. Al-Deeb and H. A. Al-Dhaheri, “Use of a pheromone-baited trap to monitor the population of the lesser date moth Batrachedra amydraula (Lepidoptera: Batrachedridae) in the UAE,” Journal of Entomology and Zoology Studies, vol. 5, no. 6, pp. 2572–2575, 2017.

 5.  M. Hoddle, A. H. Al-Abbad, H. El-Shafie, J. Faleiro, A. Sallam et al., “Assessing the impact of areawide pheromone trapping, pesticide applications, and eradication of infested date palms for Rhynchophorus ferrugineus (Coleoptera: Curculionidae) management in Al Ghowaybah, Saudi Arabia,” Crop Protection, vol. 53, pp. 152–160, 2013.

 6.  E. Sadeghi, V. Baniameri and A. Marouf, “Oviposition behaviour of Goniozus swirskiana (Hymenoptera: Bethylidae: Bethylinae) a parasitoid of Batrachedra amydraula Meyrick from the warmest desert of Iran,” World Applied Sciences Journal, vol. 20, no. 11, pp. 1493–1498, 2012.

 7.  M. Abbas, S. Al-Khatry, R. Shidi and N. A. Al-Ajmi, “Natural enemies of the lesser date moth, Batrachedra amydraula Meyrick (Lepidoptera: Batrachedridae) with special reference to its parasitoid Goniozus sp,” Egyptian Journal of Biological Pest Control, vol. 24, no. 2, pp. 293–296, 2014.

 8.  L. Addison, “Analysis of a predator-prey model: A deterministic and stochastic approach,” Journal of Biometrics and Biostatistics, vol. 8, no. 4, pp. 359–368, 2017.

 9.  B. Ahmad, S. K. Ntouyas, A. Alsaedi and A. F. Albideewi, “A study of a coupled system of Hadamard fractional differential equations with nonlocal coupled initial-multipoint conditions,” Advances in Difference Equations, vol. 2021, no. 1, pp. 1–16, 2021.

10. R. Anguelov, C. Dufourd and Y. Dumont, “Mathematical model for pest-insect control using mating disruption and trapping,” Applied Mathematical Modelling, vol. 52, pp. 437–457, 2017.

11. J. P. Ntahomvukiye, A. Temgoua and S. Bowong, “Study of the population dynamics of Busseola fusca, maize pest,” Acta Biotheoretica, vol. 66, no. 4, pp. 379–397, 2018.

12. S. Xiang, Y. Pei and X. Liang, “Analysis and optimization based on a sex pheromone and pesticide pest model with gestation delay,” International Journal of Biomathematics, vol. 12, no. 5, pp. 1950054, 2019.

13. M. D. Tapi, L. Bagny-Beilhe and Y. Dumont, “Miridae control using sex-pheromone traps. Modeling, analysis and simulations,” Nonlinear Analysis: Real World Applications, vol. 54, pp. 103082, 2020.

14. M. El-Shahed and A. Al-Dubiban, “Mathematical modelling of lesser date moth using sex pheromone traps and natural enemies,” Mathematical Problems in Engineering, vol. 2021, pp. 1–14, 2021.

15. J. Cresson and A. Szafran´ska, “Discrete and continuous fractional persistence problems the positivity property and applications,” Communications in Nonlinear Science and Numerical Simulation, vol. 44, pp. 424–448, 2017.

16. H. Li, L. Zhang, C. Hu, Y. Jiang and Z. Teng, “Dynamical analysis of a fractional-order predator- prey model incorporating a prey refuge,” Journal of Applied Mathematics and Computing, vol. 54, pp. 435–449, 2016.

17. M. Sambath, P. Ramesh and K. Balachandran, “Asymptotic behavior of the fractional order three species prey-predator model,” International Journal of Nonlinear Sciences and Numerical Simulation, vol. 19, pp. 721–733, 2018.

18. S. K. Choi, B. Kang and N. Koo, “Stability for Caputo fractional differential systems,” Abstract and Applied Analysis, vol. 2014, pp. 1–6, 2014.

19. P. V. Driessche and J. Watmough, “Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission,” Mathematical Biosciences, vol. 180, pp. 29–48, 2002.

20. J. Huo, H. Zhao and L. Zhu, “The effect of vaccines on backward bifurcation in a fractional order HIV model,” Nonlinear Analysis: Real World Applications, vol. 26, pp. 289–305, 2015.

21. A. Matouk, “Stability conditions, hyperchaos and control in a novel fractional order hyperchaotic system,” Physics Letters A, vol. 373, no. 25, pp. 2166–2173, 2009.

22. L. Perko, Differential equations and dynamical systems. Vol. 7, Springer Science & Business Media, 2013.

23. W. M. Liu, “Criterion of Hopf bifurcations without using eigenvalues,” Journal of Mathematical Analysis and Applications, vol. 182, no. 1, pp. 250–256, 1994.

24. W. Wang, S. Cheng, Z. Guo and X. Yan, “A note on the continuity for Caputo fractional stochastic differential equations,” Chaos: An Interdisciplinary Journal of Nonlinear Science, vol. 30, no. 7, pp. 073106, 2020.

25. T. Doan, P. Huong, P. Kloeden and A. Vu, “Euler-Maruyama scheme for Caputo stochastic fractional differential equations,” Journal of Computational and Applied Mathematics, vol. 380, pp. 112989, 2020.

26. H. Aminikhah, A. H. R. Sheikhani, T. Houlari and H. Rezazadeh, “Numerical solution of the distributed-order fractional Bagley-Torvik equation,” IEEE/CAA Journal of Automatica Sinica, vol. 6, no. 3, pp. 760–765, 2017.

27. M. Hou, J. Yang, S. Shi and H. Liu, “Logical stochastic resonance in a nonlinear fractional-order system,” European Physical Journal Plus, vol. 135, no. 9, pp. 747, 2020.

28. I. Podlubny, Fractional Differential Equations. New York, NY, USA: Academic Press, 1999.

29. C. Huang and M. Stynes, “Error analysis of a finite element method with GMMP temporal discretisation for a time-fractional diffusion equation,” Computers & Mathematics with Applications, vol. 79, no. 9, pp. 2784–2794, 2020.

[picture] | This work is licensed under a Creative Commons Attribution 4.0 International License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.