Processing math: 100%
Research article Special Issues

Mathematical model to investigate transmission dynamics of COVID-19 with vaccinated class

  • The susceptible, exposed, infected, quarantined and vaccinated (SEIQV) population is accounted for in a mathematical model of COVID-19. This model covers the therapy for diseased people as well as therapeutic measures like immunization for susceptible people to enable understanding of the dynamics of the disease's propagation. Each of the equilibrium points, i.e., disease-free and endemic, has been proven to be globally asymptotically stable under the assumption that R0 is smaller or larger than unity, respectively. Although vaccination coverage is high, the basic reproduction number depends on the vaccine's effectiveness in preventing disease when R0>0. The Jacobian matrix and the Routh-Hurwitz theorem are used to derive the aforementioned analysis techniques. The results are further examined numerically by using the standard second-order Runge-Kutta (RK2) method. In order to visualize the global dynamics of the aforementioned model, the proposed model is expanded to examine some piecewise fractional order derivatives. We may comprehend the crossover behavior in the suggested model's illness dynamics by using the relevant derivative. To numerical present the results, we use RK2 method.

    Citation: Mdi Begum Jeelani, Abeer S Alnahdi, Rahim Ud Din, Hussam Alrabaiah, Azeem Sultana. Mathematical model to investigate transmission dynamics of COVID-19 with vaccinated class[J]. AIMS Mathematics, 2023, 8(12): 29932-29955. doi: 10.3934/math.20231531

    Related Papers:

    [1] A. Q. Khan, Ibraheem M. Alsulami . Discrete Leslie's model with bifurcations and control. AIMS Mathematics, 2023, 8(10): 22483-22506. doi: 10.3934/math.20231146
    [2] Xiongxiong Du, Xiaoling Han, Ceyu Lei . Dynamics of a nonlinear discrete predator-prey system with fear effect. AIMS Mathematics, 2023, 8(10): 23953-23973. doi: 10.3934/math.20231221
    [3] Fatao Wang, Ruizhi Yang, Yining Xie, Jing Zhao . Hopf bifurcation in a delayed reaction diffusion predator-prey model with weak Allee effect on prey and fear effect on predator. AIMS Mathematics, 2023, 8(8): 17719-17743. doi: 10.3934/math.2023905
    [4] Abdul Qadeer Khan, Ayesha Yaqoob, Ateq Alsaadi . Neimark-Sacker bifurcation, chaos, and local stability of a discrete Hepatitis C virus model. AIMS Mathematics, 2024, 9(11): 31985-32013. doi: 10.3934/math.20241537
    [5] Xiaoming Su, Jiahui Wang, Adiya Bao . Stability analysis and chaos control in a discrete predator-prey system with Allee effect, fear effect, and refuge. AIMS Mathematics, 2024, 9(5): 13462-13491. doi: 10.3934/math.2024656
    [6] Ahmad Suleman, Rizwan Ahmed, Fehaid Salem Alshammari, Nehad Ali Shah . Dynamic complexity of a slow-fast predator-prey model with herd behavior. AIMS Mathematics, 2023, 8(10): 24446-24472. doi: 10.3934/math.20231247
    [7] Weili Kong, Yuanfu Shao . Bifurcations of a Leslie-Gower predator-prey model with fear, strong Allee effect and hunting cooperation. AIMS Mathematics, 2024, 9(11): 31607-31635. doi: 10.3934/math.20241520
    [8] Ali Al Khabyah, Rizwan Ahmed, Muhammad Saeed Akram, Shehraz Akhtar . Stability, bifurcation, and chaos control in a discrete predator-prey model with strong Allee effect. AIMS Mathematics, 2023, 8(4): 8060-8081. doi: 10.3934/math.2023408
    [9] Abdul Qadeer Khan, Ayesha Yaqoob, Ateq Alsaadi . Discrete Hepatitis C virus model with local dynamics, chaos and bifurcations. AIMS Mathematics, 2024, 9(10): 28643-28670. doi: 10.3934/math.20241390
    [10] Fethi Souna, Salih Djilali, Sultan Alyobi, Anwar Zeb, Nadia Gul, Suliman Alsaeed, Kottakkaran Sooppy Nisar . Spatiotemporal dynamics of a diffusive predator-prey system incorporating social behavior. AIMS Mathematics, 2023, 8(7): 15723-15748. doi: 10.3934/math.2023803
  • The susceptible, exposed, infected, quarantined and vaccinated (SEIQV) population is accounted for in a mathematical model of COVID-19. This model covers the therapy for diseased people as well as therapeutic measures like immunization for susceptible people to enable understanding of the dynamics of the disease's propagation. Each of the equilibrium points, i.e., disease-free and endemic, has been proven to be globally asymptotically stable under the assumption that R0 is smaller or larger than unity, respectively. Although vaccination coverage is high, the basic reproduction number depends on the vaccine's effectiveness in preventing disease when R0>0. The Jacobian matrix and the Routh-Hurwitz theorem are used to derive the aforementioned analysis techniques. The results are further examined numerically by using the standard second-order Runge-Kutta (RK2) method. In order to visualize the global dynamics of the aforementioned model, the proposed model is expanded to examine some piecewise fractional order derivatives. We may comprehend the crossover behavior in the suggested model's illness dynamics by using the relevant derivative. To numerical present the results, we use RK2 method.



    The well-known Leslie-Gower model with a Holling Ⅱ functional response is expressed as follows:

    {dxdt=rx(1xK)bxya+x,dydt=sy(1yx), (1.1)

    where the biological interpretation of each parameter is provided in Table 1. It is noteworthy that the Leslie-Gower model incorporates a response function indicating the predator's carrying capacity, which is proportional to the prey population, a feature absent in the Lotka-Volterra model. This distinction has been extensively investigated in many references, such as [1,2,3,4,5].

    Table 1.  Biological meaning of parameters in model (1.1).
    Parameter Interpretation
    x prey population density
    y predator population density
    r intrinsic growth rate of the prey
    K maximum prey carrying capacity of the environment
    b maximum predation rate when prey is abundant
    a prey population at which predator attack capability saturates
    s intrinsic growth rate of the predator

     | Show Table
    DownLoad: CSV

    Previous studies predominantly focused on direct predation effects. In 2016, Wang et al. [6] demonstrated through experiments that predator-induced fear (indirect effects) in prey leads to a reduction in prey birth rates. They introduced the fear factor F(k,y)=11+ky, where k represents the intensity of fear driving antipredator behavior. To further investigate the influence of fear on population dynamics, we enhance the system (1.1) by incorporating the fear factor F(k,y). The modified model is formulated as follows:

    {dxdt=rx1+ky(1xK)bxya+x,dydt=sy(1yx). (1.2)

    For additional significant findings on the impact of fear in predator-prey models, refer to [7,8,9].

    There has been increasing recognition that traditional integer-order differential equations may not sufficiently capture the complexity of biological systems [10,11,12,13,14]. In reality, the behaviors of most organisms in nature are influenced by their historical context. Fractional-order derivatives, which extend the concept of integer-order differentiation, offer a more flexible and accurate framework for modeling memory and hereditary properties in population dynamics. Therefore, the authors of this paper aim to apply the Caputo fractional derivative to (1.2), thereby extending it into a fractional model. The Caputo fractional derivative of a function u(t) of order α(0,1] is given in [15] as follows

    Dαtu(t)=1Γ(1α)t0u(τ)(tτ)αdτ,

    where Γ is the Gamma function. By replacing the integer-order derivative with the Caputo fractional derivative, the following model is obtained:

    {Dαtx=rx1+ky(1xK)bxya+x,Dαty=sy(1yx). (1.3)

    Specifically, when α=1, the fractional-order model (1.3) reduces to the integer-order model (1.2), demonstrating that the fractional-order model serves as a generalization of the integer-order model. Moreover, in fractional calculus, the rate of change at any given moment, expressed by the fractional-order derivative, depends on the population density over a specified time interval. This feature gives the fractional-order model (1.3) a distinct advantage in capturing memory effects within populations.

    When studying populations with nonoverlapping generations or small sizes, mathematical models are often expressed in discrete terms. Although some variables may change continuously in real-world scenarios, the recording of these changes typically happens at specific time intervals during data collection. Thus, employing discrete systems to examine the dynamic behavior of biological populations is highly practical and significant. In references [16,17], the authors employed the piecewise constant approximation method to discretize continuous fractional predator-prey models and explored their dynamical properties. Singh and Sharma [18] examine a discrete prey-predator model with Holling type Ⅱ functional response and prey refuge, identifying bifurcations and controlling chaos through state feedback, pole placement, and hybrid techniques. Berkal and Almatrafi [19] used the exponential piecewise constant argument to discretize continuous fractional activator-inhibitor system. Their analysis includes stability assessments, investigations of Neimark-Sacker and period-doubling bifurcations, and numerical simulations that validate the theoretical findings on the system's dynamics. For a more comprehensive exploration of discrete model studies, readers are advised to consult references [20,21,22,23,24,25] and the related literature cited therein. However, current literature lacks studies on fractional discrete-time predator-prey Leslie-Gower systems that incorporate the fear effect in prey population. Using the same method as in references [16,17] to discretize model (1.3), we can obtain the following discrete model:

    {xn+1=xn+hαΓ(α+1)[rxn1+kyn(1xnK)bxnyna+xn],yn+1=yn+hαΓ(α+1)[syn(1ynxn)]. (1.4)

    Here h>0 represents the time interval of production.

    Furthermore, the key contributions and findings of this study are summarized as follows:

    ● The Caputo fractional derivative of order (0,1] is utilized to incorporate the memory effect into the dynamical behavior of the proposed model.

    ● The existence and stability of fixed points are investigated.

    ● Conditions for the occurrence and direction of period-doubling bifurcation and Neimark–Sacker bifurcation at the positive fixed point are established.

    ● State feedback and hybrid control strategies are employed to manage bifurcations and chaotic behavior in the model.

    ● To validate the accuracy of our theoretical findings, numerical examples for the fractional-order discrete-time Leslie-Gower model with a fear factor are provided.

    The remainder of this paper is structured as follows: In Section 2, we investigate the existence and stability of the equilibrium points of model (1.4). Section 3 analytically demonstrates that model (1.4), under specific parametric conditions, undergoes period-doubling or Neimark-Sacker bifurcation. Section 4 explores the control of chaos toward an unstable equilibrium point using feedback control or hybrid control approaches. Section 5 presents a quantitative analysis of the dynamics of model (1.4) to validate our analytical findings. Finally, Section 6 offers brief conclusions.

    This section explores the existence of equilibrium points in model (1.4) and assesses their stability by evaluating the eigenvalues of the Jacobian matrix at these points. The following definition and lemma are introduced to assist in this stability analysis of equilibrium points.

    Definition 1. [26] Let λ1 and λ2 denote the two roots of the characteristic equation F(λ)=λ2+pλ+q=0 associated with the Jacobian matrix J(x,y). The equilibrium point (x,y) is termed

    (1) sink if |λ1|<1 and |λ2|<1, and the sink is locally asymptotically stable;

    (2) source if |λ1|>1 and |λ2|>1, and the source is locally unstable;

    (3) saddle if |λ1|>1 and |λ2|<1 (or |λ1|<1 and |λ2|>1);

    (4) non-hyperbolic if either |λ1|=1 or |λ2|=1.

    Lemma 1. [26] Let F(1)>0 in F(λ)=λ2Mλ+N, where λ1,λ2 are the two roots of F(λ)=0. Then, the following results hold true:

    (1) |λ1|<1 and |λ2|<1 if, and only if, F(1)>0 and N<1;

    (2) |λ1|<1 and |λ2|>1 (or |λ1|>1 and |λ2|<1) if, and only if, F(1)<0;

    (3) |λ1|>1 and |λ2|>1 if, and only if, F(1)>0 and N>1;

    (4) λ1=1 and |λ2|1 if, and only if, F(1)=0 and M0,2;

    (5) λ1 and λ2 are conjugate complex and |λ1,2|=1 if, and only if, M2<4N and N=1.

    It is evident that the equilibrium points of model (1.4) satisfy the following equations:

    {x=x+hαΓ(α+1)[rx1+ky(1xK)bxya+x],y=y+hαΓ(α+1)[sy(1yx)].

    This algebraic system is satisfied if x=K and y=0, indicating that the model (1.4) has a boundary equilibrium point E0(K,0) for all model parameters. To find the interior equilibrium point E1, we will solve the following system simultaneously:

    r1+ky(1xK)bya+x=0,s(1yx)=0.

    From the second equation, we find y=x. Substituting y=x into the first equation yields:

    (rK+bk)x2+(arK+br)xar=0,

    Given that all model parameters are positive, we obtain the following results through direct calculations.

    Theorem 1. The model (1.4) always has a boundary equilibrium point E0(K,0) and a positive equilibrium point E1=(x,y), where

    x=y=(arK+br)+(arK+br)2+4ar(rK+bk)2(rK+bk).

    The Jacobian matrix of model (1.4) evaluated at any point (x,y) is as follows:

    J(x,y)=[1+hαΓ(α+1)[r1+ky(12xK)aby(a+x)2]hαΓ(α+1)[rkx(1+ky)2(1xK)+bxa+x]hαΓ(α+1)sy2x21+hαΓ(α+1)(s2syx)]. (2.1)

    Theorem 2. The boundary equilibrium point E0(K,0) exhibits the following behaviors:

    (1) It is a source point if r>2Γ(α+1)hα;

    (2) It is a saddle point if 0<r<2Γ(α+1)hα;

    (3) It is a non-hyperbolic point if r=2Γ(α+1)hα.

    Proof. The Jacobian matrix (2.1) evaluated at the boundary equilibrium point E0(K,0) is given by:

    J(E0)=[1hαrΓ(α+1)hαΓ(α+1)bKa+K01+hαsΓ(α+1)].

    The eigenvalues of J(E0) are λ1=1hαrΓ(α+1) and λ2=1+hαsΓ(α+1). Clearly, |λ2|>1 and:

    |λ1|={<1if r>2Γ(α+1)hα,>1if 0<r<2Γ(α+1)hα,=1if r=2Γ(α+1)hα.

    By applying Definition 1, the proof is complete.

    Next, we analyze the local dynamics of model (1.4) at the positive equilibrium point E1(x,y). The Jacobian matrix (2.1) evaluated at E1(x,y) is expressed as:

    J(E1)=(1+Aa11Aa12As1As), (2.2)

    where

    A=hαΓ(α+1),a11=r1+kx(12xK)abx(a+x)2,a12=rkx(1+kx)2(1xK)bxa+x.

    The characteristic equation of J(E1) is given by:

    λ2Mλ+N=0, (2.3)

    where

    M=2+(a11s)A,N=1+(a11s)A(a11+a12)sA2.

    Let F(λ)=λ2Mλ+N, then:

    F(0)=N,F(1)=1+M+N,F(1)=1M+N.

    As a result, we establish the following theorem.

    Theorem 3. Let E1 be the unique positive equilibrium point of model (1.4). Then:

    (1) E1 is a sink point if |M|<1+N<2;

    (2) E1 is a source point if |M|<|1+N| and |N|>1;

    (3) E1 is a saddle point if M2>4N and |M|>|1+N|;

    (4) E1 is non-hyperbolic if |M|=|1+N| or N=1 and |M|<2.

    Proof. (1) According to Lemma 1, E1 is a sink point if and only if F(1)>0, F(1)>0, and N<1. This condition is satisfied when |M|<1+N<2. Similarly, the proofs for Theorem 3 (2)–(4) follow straightforwardly from the definitions and properties of M and N.

    In this section, we will analyze the period-doubling and Neimark-Sacker bifurcation behaviors of model (1.4) at the positive equilibrium point E1(x,y). The reason for not analyzing the boundary equilibrium point is that, at E0(K,0), the predators have become extinct, leaving only the prey.

    Assume M2>4N and

    s=s1:=4+2Aa112A+A2(a11+a12),

    and we can ascertain that the characteristic equation (2.3) satisfies

    F(1)=1+M+N=1+[2+(a11s1)A]+[1+(a11s1)A(a11+a12)s1A2]=4+2a11A[2A+(a11+a12)A2]s1=0.

    Therefore, the eigenvalues of the characteristic equation (2.3) are

    λ1=1,λ2=1+M=3+(a11s1)A.

    Let

    s2=a11+2A,s3=a11+4A.

    The condition |λ2|1 implies that s1s2,s3.

    Based on the above analysis, we can conclude that when the parameters vary within the set

    P.D={(h,α,r,k,K,b,a,s):M2>4N,s=s1andss2,s3}, (3.1)

    model (1.4) will undergo period-doubling bifurcation at E1(x,y).

    Select parameters (h,α,r,k,K,b,a,s)P.D and consider s as a small perturbation of s, i.e., s=ss1, where |s|1. After this perturbation, model (1.4) can be represented as follows:

    {xn+1=xn+hαΓ(α+1)[rxn1+kyn(1xnK)bxnyna+xn],yn+1=yn+hαΓ(α+1)[(s1+s)yn(1ynxn)]. (3.2)

    Let un=xnx and vn=yny, transforming the equilibrium point E1(x,y) into the origin O(0,0). The model (3.2) can thus be expressed as:

    {un+1=un+hαΓ(α+1)[r(un+x)1+k(vn+y)(1un+xK)b(un+x)(vn+y)a+(un+x)],vn+1=vn+hαΓ(α+1)[(s1+s)(vn+y)(1vn+yun+x)]. (3.3)

    The Taylor expansion of model (3.3) around (un,vn)=(0,0) yields the following form:

    [un+1vn+1]=[1+Aa11Aa12As11As1][unvn]+[f(un,vn,s)g(un,vn,s)], (3.4)

    where

    f(un,vn,s)=c13u2n+c14unvn+c15v2n+c16u3n+c17u2nvn+c18unv2n+c19v3n+O((|un|+|vn|+|s|)4),g(un,vn,s)=c23u2n+c24unvn+c25v2n+c26u3n+c27u2nvn+c28unv2n+c29v3n+d1uns+d2vns+d3u2ns+d4unvns+d5v2ns+O((|un|+|vn|+|s|)4),

    and

    c13=ArK(1+kx)+abAx(a+x)3,c14=Ark(1+kx)2(12xK)abA(a+x)2,c15=Ark2x(1+kx)3(1xK),c16=Aabx(a+x)4,c17=ArkK(1+kx)2+abA(a+x)3,c18=Ark2(1+kx)3(12xK),c19=Ark3x(1+kx)4(1xK),c23=c25=As1x,c24=2As1x,c26=c28=As1x2,c27=2As1x2,c29=0,d1=A,d2=A,d3=d5=Ax,d4=2Ax.

    Define

    T=[Aa12Aa122Aa11λ21Aa11],

    then

    T1=1Aa12(1+λ2)[λ21Aa11Aa122+Aa11Aa12].

    Using the transformation:

    [unvn]=T[˜un˜vn],

    model (3.4) is transformed into the following form:

    [˜un+1˜vn+1]=[100λ2][˜un˜vn]+[˜f(˜un,˜vn,s)˜g(˜un,˜vn,s)], (3.5)

    where

    ˜f(˜un,˜vn,s)=(λ21Aa11)c13Aa12c23Aa12(1+λ2)u2n+(λ21Aa11)c14Aa12c24Aa12(1+λ2)unvn+(λ21Aa11)c15Aa12c25Aa12(1+λ2)v2n+(λ21Aa11)c16Aa12c26Aa12(1+λ2)u3n+(λ21Aa11)c17Aa12c27Aa12(1+λ2)u2nvn+(λ21Aa11)c18Aa12c28Aa12(1+λ2)unv2n+(λ21Aa11)c19Aa12c29Aa12(1+λ2)v3nd11+λ2unsd21+λ2vnsd31+λ2u2nsd41+λ2unvnsd51+λ2v2ns+O((|un|+|vn|+|s|)4),
    ˜g(˜un,˜vn,s)=(2+Aa11)c13+Aa12c23Aa12(1+λ2)u2n+(2+Aa11)c14+Aa12c24Aa12(1+λ2)unvn+(2+Aa11)c15+Aa12c25Aa12(1+λ2)v2n+(2+Aa11)c16+Aa12c26Aa12(1+λ2)u3n+(2+Aa11)c17+Aa12c27Aa12(1+λ2)u2nvn+(2+Aa11)c18+Aa12c28Aa12(1+λ2)unv2n+(2+Aa11)c19+Aa12c29Aa12(1+λ2)v3n+d11+λ2uns+d21+λ2vns+d31+λ2u2ns+d41+λ2unvns+d51+λ2v2ns+O((|un|+|vn|+|s|)4),

    and

    un=Aa12˜un+Aa12˜vn,vn=(2+Aa11)˜un+(λ21Aa11)˜vn.

    Next, we apply the center manifold theorem [26] to analyze the dynamics around the equilibrium point (˜un,˜vn)=(0,0) at s=0. According to the theorem, the model (3.5) has a center manifold, which can be represented as:

    Wc(0,0,0)={(˜un,˜vn,s)R3+:˜vn=w(˜un,s),w(0,0)=0,Dw(0,0)=0}.

    Assume that

    w(˜un,s)=η1˜u2n+η2˜uns+η3(s)2+O((|˜un|+|s|)3).

    Then, the center manifold must satisfy

    w(˜un+˜f(˜un,w(˜un,s),s),s)λ2w(˜un,s)˜g(˜un,w(˜un,s),s)=0.

    By comparing the coefficients, it can be obtained that

    η1=Aa12[(2+Aa11)c13+Aa12c23](2+Aa11)[(2+Aa11)c14+Aa12c24]1λ22+(2+Aa11)2[(2+Aa11)c15+Aa12c25]Aa12(1λ22),η2=(2+Aa11)d2Aa12d1(1+λ2)2,η3=0.

    Thus, the model (3.5), when restricted to the center manifold Wc(0,0,0), is given by:

    ˜W:˜un+1=˜un+m1˜u2n+m2˜uns+m3˜u3n+m4˜u2ns+m5˜un(s)2+O((|˜un|+|s|)4) (3.6)

    where

    m1=Aa12[(λ21Aa11)c13Aa12c23]1+λ2(2+Aa11)[(λ21Aa11)c14Aa12c24]1+λ2+(2+Aa11)2[(λ21Aa11)c15Aa12c25]Aa12(1+λ2)m2=d1Aa121+λ2+d2(2+Aa11)1+λ2m3=2Aa12η1[(λ21Aa11)c13Aa12c23]1+λ2+Aa12η1(λ232Aa11)[(λ21Aa11)c14Aa12c24]1+λ22η1(2+Aa11)(λ21Aa11)[(λ21Aa11)c15Aa12c25]Aa12(1+λ2)+(Aa12)2[(λ21Aa11)c16Aa12c26]1+λ2Aa12(2+Aa11)[(λ21Aa11)c17Aa12c27]1+λ2+(2+Aa11)2[(λ21Aa11)c18Aa12c28]1+λ2(2+Aa11)3[(λ21Aa11)c19Aa12c29]Aa12(1+λ2)m4=2Aa12η2[(λ21Aa11)c13Aa12c23]1+λ2+η2(λ232Aa11)[(λ21Aa11)c14Aa12c24]1+λ22η2(2+Aa11)(λ21Aa11)[(λ21Aa11)c15Aa12c25]Aa12(1+λ2)Aa12η1d11+λ2(λ21Aa11)η1d21+λ2(Aa12)2d31+λ2+Aa12(2+Aa11)d41+λ2(2+Aa11)2d51+λ2m5=Aa12η2d11+λ2(λ21Aa11)η2d21+λ2.

    In order for Eq (3.6) to undergo period-doubling bifurcation, it is necessary that the following two quantities possess nonzero values.

    β1=(2˜W˜uns+12˜Ws2˜W˜u2n)|(0,0)=m2β2=[163˜W˜u3n+(122˜W˜u2n)2]|(0,0)=m3+m21.

    Summarize the above analysis into the following theorem.

    Theorem 4. If β1β20, then model (1.4) undergoes a period-doubling bifurcation at the positive equilibrium point E1(x,y) when parameters vary in a small neighborhood of P.D. Additionally, when β2>0 (respectively, β2<0), model (1.4) bifurcates from the equilibrium point E1(x,y) to a stable (respectively, unstable) 2-periodic orbit.

    Assume M2<4N and

    s=s4:=a111+A(a11+a12),

    and we can determine that the eigenvalues of the characteristic equation (2.3) satisfy

    |λ1,2|2=|M2±4NM22i|2=1+(a11s4)A(a11+a12)s4A2=1+a11A[1+(a11+a12)A]As4=1.

    Namely, Eq (2.3) has two complex conjugate roots with unit modulus. It is clear from the above discussion that

    d|λ1,2|ds|s=s4=A2[1+(a11+a12)A]0.

    In addition, it is crucial that when s=s4, λθ1,2(s4)1(θ=1,2,3,4), which is equivalent to M(s4)2,1,0,2. Since M2<4N, we deduce M(s4)2,2. Additionally, we necessitate that M(s4)0,1, which leads to

    s4s5:=a11+3A,s4s6:=a11+2A.

    Based on the preceding analysis, we conclude that when the parameters vary within a small neighborhood of the set

    N.S={(h,α,r,k,K,b,a,s):M2<4N,s=s4,ss5,ss6}, (3.7)

    model (1.4) undergoes a Neimark-Sacker bifurcation at E1(x,y).

    Next, assuming (h,α,r,k,K,b,a,s)N.S. and |s|1 represents a small perturbation of s4, model (1.4) can be described as follows:

    {xn+1=xn+hαΓ(α+1)[rxn1+kyn(1xnK)bxnyna+xn],yn+1=yn+hαΓ(α+1)[(s4+s)yn(1ynxn)]. (3.8)

    Let un=xnx and vn=yny, then we have

    [un+1vn+1]=[1+Aa11Aa12As41As4][unvn]+[F(un,vn)G(un,vn)] (3.9)

    where

    F(un,vn)=c13u2n+c14unvn+c15v2n+c16u3n+c17u2nvn+c18unv2n+c19v3n+O((|un|+|vn|)4),G(un,vn)=c23u2n+c24unvn+c25v2n+c26u3n+c27u2nvn+c28unv2n+c29v3n+O((|un|+|vn|)4),

    and c13, c14, c15, c16, c17, c18, c19, c23, c24, c25, c26, c27, c28, c29 are given in (2.3) by substituting s1 for s4+s.

    In order to obtain the normal form of model (3.9) at s=0, we use the following transformation:

    [unvn]=[ω1+Aa11ρ0As4][ˆunˆvn],

    with

    ρ=M2,ω=4NM22.

    Using this transformation, model (3.9) will transform as follows:

    [ˆun+1ˆvn+1]=[ρωωρ][ˆunˆvn]+[ˆF(ˆun,ˆvn)ˆG(ˆun,ˆvn)], (3.10)

    where

    ˆF(ˆun,ˆvn)=As4c13+(ρ1Aa11)c23As4ωu2n+As4c14+(ρ1Aa11)c24As4ωunvn+As4c15+(ρ1Aa11)c25As4ωv2n+As4c16+(ρ1Aa11)c26As4ωu3n+As4c17+(ρ1Aa11)c27As4ωu2nvn+As4c18+(ρ1Aa11)c28As4ωunv2n+As4c19+(ρ1Aa11)c29As4ωv3n+O((|un|+|vn|)4),ˆG(ˆun,ˆvn)=1As4(c23u2n+c24unvn+c25v2n+c26u3n+c27u2nvn+c28unv2n+c29v3n)+O((|un|+|vn|)4),

    and

    un=ωˆun+(1+Aa11ρ)ˆvn,vn=As4ˆvn.

    Next, a nonzero real number is defined as follows:

    Ω=Re[(12λ1)λ221λ1ξ11ξ12]12|ξ11|2|ξ21|2+Re(λ2ξ22) (3.11)

    where

    ξ11=14[2ˆFˆu2n+2ˆFˆv2n+i(2ˆGˆu2n+2ˆGˆv2n)]|s=0ξ12=18[2ˆFˆu2n2ˆFˆv2n+22ˆGˆunˆvn+i(2ˆGˆu2n2ˆGˆv2n22ˆFˆunˆvn)]|s=0ξ21=18[2ˆFˆu2n2ˆFˆv2n22ˆGˆunˆvn+i(2ˆGˆu2n2ˆGˆv2n+22ˆFˆunˆvn)]|s=0ξ22=116[3ˆFˆu3n+3ˆFˆunˆv2n+3ˆGˆu2nˆvn+3ˆGˆv3n+i(3ˆGˆu3n+3ˆGˆunˆv2n3ˆFˆu2nˆvn3ˆFˆv3n)]|s=0.

    Through some complicated calculations, we get

    2ˆFˆu2n|s=0=2ωAs4[c13As4+c23(ρ1Aa11)],2ˆFˆunˆvn|s=0=1As4[c14(As4)2As4(ρ1Aa11)(2c13c24)3c23(ρ1Aa11)2],2ˆFˆv2n|s=0=1ωAs4[2c15(As4)22(As4)2(ρ1Aa11)(c14c25)2As4(c24c13)(ρ1Aa11)2+2c23(ρ1Aa11)3],3ˆFˆu3n|s=0=6ω2As4[c16As4+c26(ρ1Aa11)],3ˆFˆu2nˆvn|s=0=ωAs4[2c17(As4)2As4(ρ1Aa11)(6c162c27)6c26(ρ1Aa11)2],3ˆFˆunˆv2n|s=0=1As4[2c18(As4)2+2(As4)2(c282c17)(ρ1Aa11)+6As4(c16c27)(ρ1Aa11)2+6c26(ρ1Aa11)3],3ˆFˆv3n|s=0=6(1+Aa11ρ)ωAs4[c18(As4)3+(As4)2(c28c17)(ρ1Aa11)+As4(c16c27)(ρ1Aa11)2+c26(ρ1Aa11)3],

    and

    2ˆGˆu2n|s=0=2c23ω2As4,2ˆGˆunˆvn|s=0=ωAs4[c24As42c23(ρ1Aa11)],2ˆGˆv2n|s=0=2As4[c25(As4)2As4c24(ρ1Aa11)+c23(ρ1Aa11)2],3ˆGˆu3n|s=0=6c26ω3As4,3ˆGˆu2nˆvn|s=0=ω2As4[2c27As46c26(ρ1Aa11)],3ˆGˆunˆv2n|s=0=ωAs4[2c28(As4)24c27As4(ρ1Aa11)+6c26(ρ1Aa11)2],3ˆGˆv3n|s=0=6(1+Aa11ρ)As4[c26(1+Aa11)2+c27As4(1+Aa11)+c28(As4)2(2(1+Aa11)c26+As4c27)ρ+c26ρ2].

    Based on the aforementioned calculations, we establish the theorem regarding the existence and direction of Neimark-Sacker bifurcation.

    Theorem 5. If (h,α,r,k,K,b,a,s)N.S. and Ω0, then model (1.4) undergoes a Neimark-Sacker bifurcation at the equilibrium point E1(x,y) when the parameter s varies in the vicinity of s4. Moreover, if Ω<0 (respectively, Ω>0), then an attracting (respectively, repelling) closed invariant curve bifurcates from the equilibrium point for s>0 (respectively, s<0).

    Chaos often has detrimental effects on biological systems, disrupting the ecological balance of populations and directly influencing long-term population growth projections. Implementing effective control policies not only safeguards the size of ecological populations but also establishes a strong foundation for the sustainable exploitation of ecological resources [27,28]. This section explores two control methods aimed at effectively managing the chaos generated by model (1.4).

    In this subsection, the state feedback control method [22] will be employed to regulate the chaos exhibited by model (1.4). To achieve this, we introduce the following controlled model.

    {xn+1=xn+hαΓ(α+1)[rxn1+kyn(1xnK)bxnyna+xn]+Sn,yn+1=yn+hαΓ(α+1)[syn(1ynxn)], (4.1)

    which corresponds to model (1.4). The feedback controlling force is defined as

    Sn=p1(xnx)p2(yny), (4.2)

    where (x,y) represents the positive equilibrium point of the model (1.4), and p1, p2 stand for the feedback gains. The Jacobian matrix of the controlled model (4.1) evaluated at the positive equilibrium point E1(x,y) is given by

    J1(x,y)=(1+Aa11p1Aa12p2As1As), (4.3)

    where the variables A, a11, and a12 are defined in Eq (2.2). The corresponding characteristic equation of the Jacobian matrix J1(x,y) is

    λ2(2+Aa11Asp1)λ+(1+Aa11p1)(1As)(Aa12p2)As=0. (4.4)

    Let λ1 and λ2 be the roots of the Eq (4.3), then

    λ1λ2=(1+Aa11p1)(1As)(Aa12p2)As. (4.5)

    The lines of marginal stability l1, l2, and l3 are derived by solving λ1λ2=1, λ1=1 and λ2=±1, respectively. These conditions ensure that |λ1,2|=1. Then, we derive the marginal stability lines as follows:

    l1:(1As)p1Asp2=A2s(a11+a12)+A(a11s), (4.6)
    l2:p1+p2=A(a11+a12), (4.7)
    l3:(2As)p1Asp2=A2s(a11+a12)+2A(a11s)+4. (4.8)

    Therefore, l1, l2, and l3 in the (p1,p2)-plane form a triangular region which leads to |λ1,2|<1.

    Theorem 6. If p1 and p2 lie within a triangular region bounded by the lines l1, l2, and l3, it can be concluded that the model (4.1) is stable.

    Next, we apply the hybrid control approach proposed by [23] to control chaos. The controlled model of (1.4) with the hybrid control approach is depicted below:

    {xn+1=ρxn+ρhαΓ(α+1)[rxn1+kyn(1xnK)bxnyna+xn]+(1ρ)xn,yn+1=ρyn+ρhαΓ(α+1)[syn(1ynxn)]+(1ρ)yn. (4.9)

    where 0<ρ<1, and the controlled strategy in (4.9) combines feedback control and parameter perturbation. By appropriately choosing the controlled parameter ρ, the chaotic behaviors of the equilibrium point (x,y) of the controlled model (4.9) can be accelerated (delayed) or even entirely eliminated. The Jacobian matrix of the controlled model (4.9), evaluated at the positive equilibrium point (x,y), is given by

    J2(x,y)=(1+Aρa11Aρa12Aρs1Aρs), (4.10)

    where the variables A, a11, and a12 are defined in Eq (2.2). Then, the positive equilibrium point (x,y) of the controlled model (4.9) is locally asymptotically stable if the roots of the characteristic polynomial of (4.10) lie within the open unit disk. According to the Jury condition, the equilibrium point of the model remains stable if, and only if, the following conditions are met:

    |2+Aρa11Aρs|<1+(1+Aρa11)(1Aρs)(Aρa12)(Aρs)<2.

    Theorem 7. If |2+Aρa11Aρs|<1+(1+Aρa11)(1Aρs)(Aρa12)(Aρs)<2 can hold, it can be concluded that the model (4.9) is stable.

    In this section, by exploring specific cases of model (1.4), we confirm the theoretical analysis above and discover new intriguing complex dynamic behaviors. Additionally, we validate the effectiveness of linear feedback techniques and hybrid control strategies for chaos control through numerical simulations.

    Take

    h=0.69,α=0.85,r=2.27,k=1.64,K=2.16,b=0.24,a=2.02. (5.1)

    After straightforward calculations, we determine the positive equilibrium point E1=(1.7498,1.7498) and the critical bifurcation value s1=2.8472. The Jacobian matrix of model (1.4) evaluated at E1 with s=s1 is given by

    J(E1)=(0.67330.14972.19651.1965),

    whose characteristic equation is

    F(λ)=λ2+0.5232λ0.4768=0. (5.2)

    The roots of (5.2) are λ1=1 and λ2=0.4768. Moreover, we have β1=0.79590 and β2=1.0080>0. According to Theorem 4, model (1.4) undergoes period-doubling bifurcation at E1 as s passes through s1. This behavior, verified through corresponding bifurcation diagrams shown in Figure 1, utilizes initial conditions (x0,y0)=(1.74,1.74) and varies s in the range [2.7,3.7].

    Figure 1.  Bifurcation diagrams of model (1.4) with parameter values as given in (5.1) and initial conditions (1.74,1.74).

    Take

    h=0.75,α=0.8,r=4.2,k=0.5,K=3.5,b=3.6,a=8.1. (5.3)

    After straightforward calculations, we determine the positive equilibrium point E1=(2.1746,2.1746) and the critical bifurcation value s4=1.1872. The Jacobian matrix of model (1.4) evaluated at E1 with s=s4 is given by

    J(E1)=(0.07120.98841.01260.0126),

    whose characteristic equation is

    F(λ)=λ20.0586λ+1=0. (5.4)

    The roots of (5.4) are λ1,2=0.02930±0.9996i. Moreover, we have d=0.39120 and Ω=1392.0464<0. According to Theorem 5, model (1.4) undergoes Neimark-Sacker bifurcation at E1 as s passes through s4. This behavior, verified through corresponding bifurcation diagrams shown in Figure 2, utilizes initial conditions (x0,y0)=(2.2,2.2) and varies s in the range [1.1,1.6].

    Figure 2.  Bifurcation diagrams of model (1.4) with parameter values as given in (5.3) and initial conditions (2.2,2.2).

    We chose the parameter values as follows:

    h=0.75,α=0.8,r=4.2,k=0.5,K=3.5,b=3.6,a=8.1,s=1.58. (5.5)

    By Figure 3, we can get that the variables xn and yn in the model (1.4) are in a chaotic state.

    Figure 3.  Plots for model (1.4) with parameter values as given in (5.5) and initial conditions (2.2,2.2).

    We use the same parameter values as in (5.5). In Figure 4, the triangular region defined by Theorem 6 bounds the parameters p1 and p2. Inside this region, the chaotic behavior produced by model (1.4) is effectively managed, resulting in asymptotic convergence toward the equilibrium point E1=(2.1746,2.1746).

    Figure 4.  Stability region for the controlled model (4.1).

    In this case, with feedback gains set to p1=0.57 and p2=0.38 and the controller activated at the 3000th iteration of the model, Figure 5 demonstrates the control effect. A chaotic trajectory is successfully stabilized at the equilibrium point E1=(2.1746,2.1746). This indicates that the feedback control approach is effective in mitigating bifurcation and chaos.

    Figure 5.  Plots for controlled model (4.1) with parameter values as given in (5.5), p1=0.57, p2=0.38, and initial conditions (2.2,2.2).

    Finally, using the same parameter values as in (5.5), the Jacobian matrix of the controlled model (4.9), evaluated at E1, is given by:

    (10.9288ρ0.988439ρ1.34765ρ11.34765ρ). (5.6)

    The characteristic polynomial of (5.6) is given by

    λ2(22.27645ρ)λ+2.58377ρ22.27645ρ+1=0. (5.7)

    The roots of (5.7) lie within the open unit disk if, and only if, 0<ρ<0.881058. Additionally, the plots for xn and yn of the controlled model (4.9) are shown in Figure 6 with ρ=0.87. From Figure 6, it is clear that the positive equilibrium point E1 is stable. Therefore, it can be concluded that employing the hybrid control approach is effective in mitigating bifurcation and chaos.

    Figure 6.  Plots for controlled model (4.9) with parameter values as given in (5.5), ρ=0.87, and initial conditions (2.2,2.2).

    In this study, we proposed a fractional Leslie-Gower model with a Holling type Ⅱ functional response and antipredator behavior. By employing the piecewise constant approximation method, we derived the discrete model (1.4) and analyzed its dynamical behavior, including the existence and stability of equilibrium points and the possibility of local bifurcations. The following conclusions can be drawn from our research:

    (1) The discrete model (1.4) has two equilibrium points: E0(K,0) and E1(x,y). The only positive coexistence equilibrium point, E1(x,y), reflects the coexistence of predators and prey.

    (2) Our theoretical analysis and numerical simulations of the positive equilibrium point E1(x,y) indicate that the model undergoes period-doubling bifurcation and Neimark-Sacker bifurcation under specific parameter conditions. Figures 1 and 2 illustrate how these bifurcations can lead to chaotic behavior at E1.

    (3) By applying state feedback control and hybrid control methods, we effectively managed the chaotic behavior generated by the model (1.4), as shown in Figures 36. These interventions mitigated the adverse effects of chaos and bifurcations, consequently enhancing ecosystem resilience.

    This study advances our understanding of the complex dynamics in ecological models influenced by fear effects and provides practical techniques for controlling chaotic behavior in such models. Future work could explore different discretization methods for the model, and new parameters could be chosen to study the influence of various ecological effects on population dynamics.

    Yao Shi: Formal analysis, validation, writing-original draft & editing, visualization, methodology; Zhenyu Wang: Conceptualization, investigation, writing-review & editing, software, supervision. All authors have read and approved the final version of the manuscript for publication.

    This work was supported by the Innovation Foundation of Hebei University of Engineering (Grant Nos. SJ2401002097), the National Natural Science Foundation of China (Grant Nos. 12401519), and the Central Guidance on Local Science and Technology Development Fund of Hebei Province (Grant Nos. 246Z1825G).

    The authors declare no conflict of interest.



    [1] Naming the coronavirus disease (COVID-19) and the virus that causes it, Available from: World Health Organization (WHO), 2020, https://www.who.int/emergencies/diseases/novel-coronavirus-2019/technical-guidance/naming-the-coronavirus-disease-(covid-2019)-and-the-virus-that-causes-it.
    [2] S. Zhao, S. S. Musa, Q. Lin, J. Ran, G. Yang, W. Wang, et al., Estimating the unreported number of novel coronavirus (2019-nCoV) cases in China in the first half of January 2020, a data-driven modelling analysis of the early outbreak, J. Clin. Med., 9 (2020), 388. https://doi.org/10.3390/jcm9020388 doi: 10.3390/jcm9020388
    [3] I. Nesteruk, Statistics based predictions of coronavirus 2019-nCoV spreading in mainland China, MedRxiv, 4 (2020), 1988–1989. https://doi.org/10.1101/2020.02.12.20021931 doi: 10.1101/2020.02.12.20021931
    [4] D. S. Hui, E. I. Azhar, T. A. Madani, F. Ntoumi, R. Kock, O. Dar, et al., The continuing 2019-nCoV epidemic threat of novel coronaviruses to global health–The latest 2019 novel coronavirus outbreak in Wuhan, China, Int. J. Infect. Dis., 91 (2020), 264–266. https://doi.org/10.1016/j.ijid.2020.01.009 doi: 10.1016/j.ijid.2020.01.009
    [5] S. Zhao, Q. Lin, J. Ran, S. S. Musa, G. Yang, W. Wang, et al., Preliminary estimation of the basic reproduction number of novel coronavirus (2019-nCoV) in China, Int. J. Infect. Dis., 92 (2020), 214–217. https://doi.org/10.1016/j.ijid.2020.01.050 doi: 10.1016/j.ijid.2020.01.050
    [6] K. Shah, R. Din, W. Deebani, P. Kumam, Z. Shah, On nonlinear classical and fractional order dynamical system addressing COVID-19, Results Phys., 24 (2021), 104069. https://doi.org/10.1016/j.rinp.2021.104069 doi: 10.1016/j.rinp.2021.104069
    [7] A. J. Lotka, Contribution to the theory of periodic reactions, J. Phys. Chem., 14 (1910), 271–274. https://doi.org/10.1021/j150111a004 doi: 10.1021/j150111a004
    [8] N. S. Goel, S. C. MAITRA, E. W. MONTROLL, On the Volterra and other nonlinear models of interacting populations, Rev. Mod. Phys., 43 (1971), 231. https://doi.org/10.1103/RevModPhys.43.231 doi: 10.1103/RevModPhys.43.231
    [9] P. Zhou, X. L. Yang, X. G. Wang, B. Hu, L. Zhang, W. Zhang, et al., A pneumonia outbreak associated with a new coronavirus of probable bat origin, Nature, 579 (2020), 270–273. https://doi.org/10.1038/s41586-020-2012-7 doi: 10.1038/s41586-020-2012-7
    [10] Q. Li, X. Guan, P. Wu, X. Wang, L. Zhou, Y. Tong, et al., Early transmission dynamics in Wuhan, China, of novel coronavirus infected pneumonia, N. Engl. J. Med., 382 (2020), 1199–1207. https://doi.org/10.1056/NEJMoa2001316 doi: 10.1056/NEJMoa2001316
    [11] I. I. Bogoch, A. Watts, A. Thomas-Bachli, C. Huber, M. U. G. Kraemer, K. Khan, Pneumonia of unknown aetiology in Wuhan, China: potential for international spread via commercial air travel, J. Travel Med., 27 (2020), taaa008. https://doi.org/10.1093/jtm/taaa008 doi: 10.1093/jtm/taaa008
    [12] C. Lu, H. Liu, D. Zhang, Dynamics and simulations of a second order stochastically perturbed SEIQV epidemic model with saturated incidence rate, Chaos Soliton. Fract., 152 (2021), 111312. https://doi.org/10.1016/j.chaos.2021.111312 doi: 10.1016/j.chaos.2021.111312
    [13] X. Liu, L. Yang, Stability analysis of a SEIQV epidemic model with saturated incidence rate, Nonlinear Anal. Real, 13 (2012), 2671–2679. https://doi.org/10.1016/j.nonrwa.2012.03.010 doi: 10.1016/j.nonrwa.2012.03.010
    [14] A. B. Gumel, S. Ruan, T. Day, J. Watmough, F. Brauer, P. van den Driessche, et al., Modelling strategies for controlling SARS out breaks, Proc. R. Soc. Lond. B, 271 (2004), 2223–2232. https://doi.org/10.1098/rspb.2004.2800 doi: 10.1098/rspb.2004.2800
    [15] A. Atangana, S. I. Araz, Modeling third waves of Covid-19 spread with piecewise differential and integral operators: Turkey, Spain and Czechia, Results Phys., 29 (2021), 104694. https://doi.org/10.1016/j.rinp.2021.104694 doi: 10.1016/j.rinp.2021.104694
    [16] Z. Zhang, A. Zeb, S. Hussain, E. Alzahrani, Dynamics of COVID-19 mathematical model with stochastic perturbation, Adv. Differ. Equ., 2020 (2020), 451. https://doi.org/10.1186/s13662-020-02909-1 doi: 10.1186/s13662-020-02909-1
    [17] A. Atangana, S. I. Araz, Mathematical model of COVID-19 spread in Turkey and South Africa: theory, methods, and applications, Adv. Differ. Equ., 2020 (2020), 659. https://doi.org/10.1186/s13662-020-03095-w doi: 10.1186/s13662-020-03095-w
    [18] N. H. Alharthi, M. B. Jeelani, A Fractional model of COVID-19 in the frame of environmental transformation with caputo fractional derivative, Adv. Appl. Stat., 88 (2023), 225–244. https://doi.org/10.17654/0972361723047 doi: 10.17654/0972361723047
    [19] M. B. Jeelani, Stability and computational analysis of COVID-19 using a higher order galerkin time discretization scheme, Adv. Appl. Stat., 86 (2023), 167–206. https://doi.org/10.17654/0972361723022 doi: 10.17654/0972361723022
    [20] C. A. B. Pearson, F. Bozzani, S. R. Procter, N. G. Davies, M. Huda, H. T. Jensen, et al., COVID-19 vaccination in Sindh Province, Pakistan: A modelling study of health impact and cost-effectiveness, PLoS Med., 18 (2021), e1003815. https://doi.org/10.1371/journal.pmed.1003815 doi: 10.1371/journal.pmed.1003815
    [21] R. P. Curiel, H. G. Ramírez, Vaccination strategies against COVID-19 and the diffusion of anti-vaccination views, Sci. Rep., 11 (2021), 6626. https://doi.org/10.1038/s41598-021-85555-1 doi: 10.1038/s41598-021-85555-1
    [22] J. T. Wu, K. Leung, G. M. Leung, Nowcasting and forecasting the potential domestic and international spread of the 2019-nCoV outbreak originating in Wuhan, China: a modelling study, The Lancet, 395 (2020), 689–697. https://doi.org/10.1016/S0140-6736(20)30260-9 doi: 10.1016/S0140-6736(20)30260-9
    [23] M. S. Arshad, D. Baleanu, M. B. Riaz, M. Abbas, A novel 2-stage fractional Runge-kutta method for a time-fractional logistic growth model, Discrete Dyn. Nat. Soc., 2000 (2000), 1020472. https://doi.org/10.1155/2020/1020472 doi: 10.1155/2020/1020472
    [24] T. Abdeljawad, Q. M. Al-Mdallal, F. Jarad, Fractional logistic models in the frame of fractional operators generated by conformable derivatives, Chaos Soliton. Fract., 119 (2019), 94–101. https://doi.org/10.1016/j.chaos.2018.12.015 doi: 10.1016/j.chaos.2018.12.015
    [25] O. A. Omar, R. A. Elbarkouky, H. M. Ahmed, Fractional stochastic modelling of COVID-19 under wide spread of vaccinations: Egyptian case study, Alex. Eng. J., 61 (2022), 8595–8609. https://doi.org/10.1016/j.aej.2022.02.002 doi: 10.1016/j.aej.2022.02.002
    [26] S. Saha, A. K. Saha, Modeling the dynamics of COVID-19 in the presence of Delta and Omicron variants with vaccination and non-pharmaceutical interventions, Heliyon, 9 (2023), e17900. https://doi.org/10.1016/j.heliyon.2023.e17900 doi: 10.1016/j.heliyon.2023.e17900
    [27] H. M. Ahmed, R. A. Elbarkouky, O. A. M. Omar, M. A. Ragusa, Models for COVID-19 daily confirmed cases in different countries, Mathematics, 9 (2021), 659. https://doi.org/10.3390/math9060659 doi: 10.3390/math9060659
    [28] F. Liu, K. Burrage, Novel techniques in parameter estimition for fractinal dynamical models arising from biological systems, Comput. Math. Appl., 62 (2011), 822–833. https://doi.org/10.1016/j.camwa.2011.03.002 doi: 10.1016/j.camwa.2011.03.002
    [29] M. T. Hoang, O. F. Egbelowo, Dynamics of a fractional-order hepatitis b epidemic model and its solutions by nonstandard numerical schemes, In: Mathematical modelling and analysis of infectious diseases, Cham: Springer, 2020,127–153. https://doi.org/10.1007/978-3-030-49896-2_5
    [30] A. Zeb, E. Alzahrani, V. S. Erturk, G. Zaman, Mathematical model for coronavirus disease 2019 (COVID-19) containing isolation class, BioMed Res. Int., 2020 (2020), 3452402. https://doi.org/10.1155/2020/3452402 doi: 10.1155/2020/3452402
    [31] K. Shah, A. Ali, S. Zeb, A. Khan, M. A. Alqudah, T. Abdeljawad, Study of fractional order dynamics of nonlinear mathematical model, Alex. Eng. J., 61 (2022), 11211–11224. https://doi.org/10.1016/j.aej.2022.04.039 doi: 10.1016/j.aej.2022.04.039
    [32] S. Boccaletti, W. Ditto, G. Mindlin, A. Atangana, Modeling and forecasting of epidemic spreading: The case of Covid-19 and beyond, Chaos Soliton. Fract., 135 (2020), 109794. https://doi.org/10.1016/j.chaos.2020.109794 doi: 10.1016/j.chaos.2020.109794
    [33] E. Atangana, A. Atangana, Facemasks simple but powerful weapons to protect against COVID-19 spread: Can they have sides effects, Results Phys., 19 (2020), 103425. https://doi.org/10.1016/j.rinp.2020.103425 doi: 10.1016/j.rinp.2020.103425
    [34] S. T. M. Thabet, M. S. Abdo, K. Shah, T. Abdeljawad, Study of transmission dynamics of COVID-19 mathematical model under ABC fractional order derivative, Results Phys., 19 (2020), 103507. https://doi.org/10.1016/j.rinp.2020.103507 doi: 10.1016/j.rinp.2020.103507
    [35] A. Al Elaiw, F. Hafeez, M. B. Jeelani, M. Awadalla, K. Abuasbeh, Existence and uniqueness results for mixed derivative involving fractional operators, AIMS Mathematics, 8 (2023), 7377–7393. https://doi.org/10.3934/math.2023371 doi: 10.3934/math.2023371
    [36] A. Moumen, R. Shafqat, A. Alsinai, H. Boulares, M. Cancan, M. B. Jeelani, Analysis of fractional stochastic evolution equations by using Hilfer derivative of finite approximate controllability, AIMS Mathematics, 8 (2023), 16094–16114. https://doi.org/10.3934/math.2023821 doi: 10.3934/math.2023821
    [37] J. T. Machado, V. Kiryakova, F. Mainardi, Recent history of fractional calculus, Commun. Nonlinear Sci., 16 (2011), 1140–1153. https://doi.org/10.1016/j.cnsns.2010.05.027 doi: 10.1016/j.cnsns.2010.05.027
    [38] F. C. Meral, T. J. Royston, R. Magin, Fractional calculus in viscoelasticity: an experimental study, Commun. Nonlinear Sci., 15 (2010), 939–945. https://doi.org/10.1016/j.cnsns.2009.05.004 doi: 10.1016/j.cnsns.2009.05.004
    [39] L. M. Richard, Fractional calculus in bioengineering, part 1, Critical Reviews in Biomedical Engineering, 32 (2004), 104. https://doi.org/10.1615/CritRevBiomedEng.v32.i1.10 doi: 10.1615/CritRevBiomedEng.v32.i1.10
    [40] M. Dalir, M. Bashour, Applications of fractional calculus, Appl. Math. Sci., 4 (2010), 1021–1032.
    [41] A. S. Alnahdi, M. B. Jeelani, H. A. Wahash, M. A. Abdulwasaa, A Detailed Mathematical Analysis of the Vaccination Model for COVID-19, Computer Modeling in Engineering Sciences, 135 (2022), 1315–1343. https://doi.org/10.32604/cmes.2022.023694 doi: 10.32604/cmes.2022.023694
    [42] K. Dehingia, M. B. Jeelani, A. Das, Artificial intelligence and machine learning: A smart science approach for cancer control, In: Advances in deep learning for medical image analysis, Boca Raton: CRC Press, 2022. https://doi.org/10.1201/9781003230540
    [43] M. B. Jeelani, A. S. Alnahdi, M. A. Almalahi, M. S. Abdo, H. A. Wahash, N. H. Alharthi, Qualitative analyses of fractional integro-differential equations with a variable order under the Mittag-Leffler power law, J. Funct. Space., 2022 (2022), 6387351. https://doi.org/10.1155/2022/6387351 doi: 10.1155/2022/6387351
    [44] R. L. Magin, Fractional calculus in bioengineering: A tool to model complex dynamics, In: Proceedings of the 13th International Carpathian Control Conference (ICCC), High Tatras, Slovakia, 2012,464–469. https://doi.org/10.1109/CarpathianCC.2012.6228688
    [45] Y. A. Rossikhin, M. V. Shitikova, Applications of fractional calculus to dynamic problems of linear and nonlinear hereditary mechanics of solids, 50 (1997), 15–67. https://doi.org/10.1115/1.3101682
    [46] A. Carpinteri, F. Mainardi, Fractals and fractional calculus in continuum mechanics, Vienna: Springer, 1997. https://doi.org/10.1007/978-3-7091-2664-6
    [47] L. M. Richard, Fractional calculus models of complex dynamics in biological tissues, Comput. Math. Appl., 59 (2010), 1586–1593. https://doi.org/10.1016/j.camwa.2009.08.039 doi: 10.1016/j.camwa.2009.08.039
    [48] M. Shimizu, W. Zhang, Fractional calculus approach to dynamic problems of viscoelastic materials, JSME International Journal Series C Mechanical Systems, Machine Elements and Manufacturing, 42 (1999), 825–837. https://doi.org/10.1299/jsmec.42.825 doi: 10.1299/jsmec.42.825
    [49] F. Mainardi, An historical perspective on fractional calculus in linear viscoelasticity, Fract. Calc. Appl. Anal., 15 (2012), 712–717. https://doi.org/10.2478/s13540-012-0048-6 doi: 10.2478/s13540-012-0048-6
    [50] Z. Dai, Y. Peng, H. A. Mansy, R. H. Sandler, T. J. Royston, A model of lung parenchyma stress relaxation using fractional viscoelasticity, Med. Eng. Phys., 37 (2015), 752–758. https://doi.org/10.1016/j.medengphy.2015.05.003 doi: 10.1016/j.medengphy.2015.05.003
    [51] M. A. Matlob, Y. Jamali, The concepts and applications of fractional order differential calculus in modeling of viscoelastic systems, Critical Reviews in Biomedical Engineering, 47 (2019), 249–276. https://doi.org/10.1615/CritRevBiomedEng.2018028368 doi: 10.1615/CritRevBiomedEng.2018028368
    [52] W. Grzesikiewicz, A. Wakulicz, A. Zbiciak, Non-linear problems of fractional calculus in modeling of mechanical systems, Int. J. Mech. Sci., 70 (2013), 90–98. https://doi.org/10.1016/j.ijmecsci.2013.02.007 doi: 10.1016/j.ijmecsci.2013.02.007
    [53] C. Celauro, C. Fecarotti, A. Pirrotta, A. C. Collop, Experimental validation of a fractional model for creep/recovery testing of asphalt mixtures, Constr. Build. Mater., 36 (2012), 458–466. https://doi.org/10.1016/j.conbuildmat.2012.04.028 doi: 10.1016/j.conbuildmat.2012.04.028
    [54] W. Adel, A. Elsonbaty, A. Aldurayhim, A. El-Mesady, Investigating the dynamics of a novel fractional-order monkeypox epidemic model with optimal control, Alex. Eng. J., 73 (2023), 519–542. https://doi.org/10.1016/j.aej.2023.04.051 doi: 10.1016/j.aej.2023.04.051
    [55] A. El-Mesady, A. Elsonbaty, W. Adel, On nonlinear dynamics of a fractional order monkeypox virus model, Chaos Soliton. Fract., 164 (2022), 112716. https://doi.org/10.1016/j.chaos.2022.112716 doi: 10.1016/j.chaos.2022.112716
    [56] N. Ahmed, A. Elsonbaty, A. Raza, M. Rafiq, W. Adel, Numerical simulation and stability analysis of a novel reaction-diffusion COVID-19 model, Nonlinear Dyn., 106 (2021), 1293–1310. https://doi.org/10.1007/s11071-021-06623-9 doi: 10.1007/s11071-021-06623-9
    [57] A. M. R. Elsonbaty, Z. Sabir, R. Ramaswamy, W. Adel, Dynamical analysis of a novel discrete fractional SITRS model for COVID-19, Fractals, 29 (2021), 2140035. https://doi.org/10.1142/S0218348X21400351 doi: 10.1142/S0218348X21400351
    [58] A. El-Mesady, A. Waleed Adel, A. A. Elsadany, A. Elsonbaty, Stability analysis and optimal control strategies of a fractional-order Monkeypox virus infection model, Phys. Scr., 98 (2023), 095256. https://doi.org/10.1088/1402-4896/acf16f doi: 10.1088/1402-4896/acf16f
    [59] M. M. Khalsaraei, An improvement on the positivity results for 2-stage explicit Runge-Kutta methods, J. Comput. Appl. Math., 235 (2010), 137–143. https://doi.org/10.1016/j.cam.2010.05.020 doi: 10.1016/j.cam.2010.05.020
    [60] Z. J. Fu, Z. C. Tang, H. T. Zhao, P. W. Li, T. Rabczuk, Numerical solutions of the coupled unsteady nonlinear convection-diffusion equations based on generalized finite difference method, Eur. Phys. J. Plus, 134 (2019), 272. https://doi.org/10.1140/epjp/i2019-12786-7 doi: 10.1140/epjp/i2019-12786-7
    [61] A. Atangana, S. I. Araz, New concept in calculus: piecewise differential and integral operators, Chaos Soliton. Fract., 145 (2021), 110638. https://doi.org/10.1016/j.chaos.2020.110638 doi: 10.1016/j.chaos.2020.110638
    [62] Current information about COVID-19 in Pakistan, 2021, Available from: https://www.worldometers.info.
    [63] Pakistan COVID-19 Corona tracker, 2021, Available from: https://www.coronatracker.com/country/pakistan/.
    [64] F. Chamchod, N. F. Britton, On the dynamics of a two-strain influenza model with isolation, Math. Model. Nat. Phenom., 7 (2012), 49–61. https://doi.org/10.1051/mmnp/20127305 doi: 10.1051/mmnp/20127305
    [65] Pakistan population, Available from: Worldometer, https://www.worldometers.info/world-population/pakistan-population/.
    [66] S. Ahmad, A. Ullah, Q. M. Al-Mdallal, H. Khan, K. Shah, A. Khan, Fractional order mathematical modeling of COVID-19 transmission, Chaos Soliton. Fract., 139 (2020), 110256. https://doi.org/10.1016/j.chaos.2020.110256 doi: 10.1016/j.chaos.2020.110256
    [67] Vaccines, Available from: UNICEF Pakistan, https://www.unicef.org/pakistan/topics/vaccines
    [68] S. Mwalili, M. Kimathi, V. Ojiambo, D. Gathungu, R. Mbogo, SEIR model for COVID-19 dynamics incorporating the environment and social distancing, BMC Res. Notes, 13 (2020), 352. https://doi.org/10.1186/s13104-020-05192-1 doi: 10.1186/s13104-020-05192-1
  • This article has been cited by:

    1. Yao Shi, Xiaozhen Liu, Zhenyu Wang, Bifurcation Analysis and Chaos Control of a Discrete Fractional-Order Modified Leslie–Gower Model with Nonlinear Harvesting Effects, 2024, 8, 2504-3110, 744, 10.3390/fractalfract8120744
    2. Rutaba Saleem, Ahmad Javid, Interaction solutions and qualitative analysis of the Painlevé integrable (2+1)-dimensional extended Kadomtsev–Petviashvili equation, 2025, 128, 11100168, 611, 10.1016/j.aej.2025.05.013
  • Reader Comments
  • © 2023 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(1581) PDF downloads(62) Cited by(0)

Figures and Tables

Figures(18)  /  Tables(2)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog