Modelling and control of a fractional-order epidemic model with fear effect.

In this paper, we formulate and study a new fractional-order SIS epidemic model with fear effect of an infectious disease and treatment control. The existence and uniqueness, nonnegativity and finiteness of the system solutions for the proposed model have been analysed. All equilibria of the model system are found, and their local and also global stability analyses are examined. Conditions for fractional backward and fractional Hopf bifurcation are also analysed. We study how the disease control parameter, level of fear and fractional order play a role in the stability of equilibria and Hopf bifurcation. Further, we have established our analytical results through several numerical simulations. © The Joint Center on Global Change and Earth System Science of the University of Maryland and Beijing Normal University 2020.

Introduction

Infectious diseases have become one of the most threatening issues in todays lifestyle. The infectious diseases including chickenpox, measles, cholera, tuberculosis, influenza, SARS, COVID-19, etc., have massive impact than the other types of noninfectious diseases as these types of diseases can be transmitted from one individual to another, some are spread out by bites from insects or animals, and some diseases are acquired by consuming flagitious water or food or being exposed to organisms in the environment. Hence, an infectious disease may spread in a huge region throughout the globe within a very short time period. Due to the improvement of lifestyle of common people and massive enhancement in transportation and globalization, an infectious diseases become pandemic in a less amount of time compared to earlier days (for example Spanish flu in twentieth century had taken a long period of time to become a pandemic compared to the pandemic due to COVID-19 in ongoing time period). Thus, not only to study the dynamics of an infectious disease, but also to control or determining the procedure to control an infectious disease, researchers from various fields have engaged themselves.

Several mathematical tools induce a great affect in regulating many infectious diseases. Mathematical modelling is one of the popular and commonly used mathematical tools that can be utilized efficaciously to monitor several infectious diseases. The first effective research work on studying an infectious diseases using mathematical model was probably first done by Kermack and Mckendric (1927). There are several research works where the authors have been studied mathematical models on controlling infectious diseases (see Zhou et al. 2014; Jana et al. 2016, 2017a). In his book, Murray (2002) has analysed theoretical works on some simple SI,  SIS,  SIR epidemic models. Present days mathematicians formulate more complex epidemic model systems which are almost identical to the real-world problems. In mathematical studies of an epidemic model, proper control strategies are significant tools in monitoring and controlling the infectious diseases. In this regard, two most effective and commonly used control tools used in controlling diseases are vaccination and treatment (Kar and Jana 2013a, b). However, isolation (see Jana et al. 2017b), insecticide (to control vector-borne disease)(see Kar and Jana 2013b), etc., are also used in recent days.

In the present work, we consider fear factor due to the infectious disease which is an another important parameter for epidemic model. People generally got scared and make a significant distance to prevent the infectious disease. Hence, fear induced by the infectious disease compelled susceptible population to isolate which actually decreases the birth rate of population and also survival of adults is affected consequently. Due to the SARS outbreak in HongKong (SARS started on November 2002; peaked on March 2003; eliminated on June, 2003), the birth rate dramatically had fallen from 8.742 (2002) to 8.436 (2003) and again it increased to 8.558 in 2004 (see Worldbank (2018)). In our work, we formulate a new fractional-order SIS epidemic model with fear effect (Clinchy et al. 2013; Wang et al. 2016) and treatment control to eradicate the disease.

Mathematical models based on ordinary differential equations depict the interactions between the population classes (e.g. susceptible-infected), and it is a classical approach in theoretical epidemiology. In recent days, researchers are concerned in developing mathematical model by the fractional-order differential equations because it is an important apparatus for the study of the memory and also some hereditary properties of several biological components and also it has a very close relation to the fractal theory. These are the main advantages of the fractional-order derivatives which are not included in the models based on the ordinary-order derivatives. Therefore, fractional-order derivatives are more naturalistic than the ordinary derivatives. The models constructed by the fractional-order derivatives have been widely applied in different field of research after the some famous books and research works on fractional-order differential equations (Diethelm and Braunschweig 2003; Hilfer 2000; Kilbas et al. 2006; Miller and Ross 1993; Petras 2011; Podlubny 1999; Sabatier et al. 2007; Sengupta et al. 2020; Yadav et al. 2020; Karthikeyan and Arul 2020). However, a fractional-order derivative may be defined in several ways. The most popular and commonly useful definitions of fractional-order derivatives are in the sense of Riemann–Liouville, Grünwald–Letnikov and Caputo definitions (Petras 2011; Podlubny 1999). Since, in the definition due to Caputo, the initial conditions can be expressed in a similar fashion as the integer-order differentiation, it is the most frequently used fractional-order derivative in mathematical modelling. There are very few theories to analyse the dynamical behaviour of the mathematical models with fractional-order derivatives (Delavari et al. 2012; Deshpande et al. 2017; Guo 2014; Li et al. 2016; Liang and Wu 2015; Li et al. 2010). In recent days, the linearization, Lyapunov method and Lyapunov direct method (due to Lyapunov function Garrappa 2010) have established to study the local and global stability of equilibrium point(s) of a fractional-order system.

Our research article studies the modelling and control of an epidemic problem using some advanced mathematical tools. In particular, this article focuses on approaches that integrate ecology and environment considerations and promote multidisciplinary solutions for tackling environmental problems and complex socio-economic problems. Present global scenarios show that human epidemiology has a significant impact on both the ecology and the environment. As this research article studies the modelling and control of an epidemic problem, the theoretical analysis of epidemiology is relevant to both ecology and environment (see Venturino et al. 2016; Maji et al. 2019; Khatua et al. 2020, etc.). Moreover, the modelling approach using fractional-order differential equations gives an additional novelty and dimension of the article.

The rest of the article is prepared as follows: in Sect. 2, we formulate an SIS-type fractional-order epidemic system and then discuss some preliminaries on fractional-order derivatives in Sect. 3. In Sect. 4, we describe the dynamics of the proposed model, including equilibria, their existence, bifurcation and both the local and the global behaviours. In Sect. 5, we present some computer simulation works to validate our theoretical observations. In Sect. 6, we provide some of the central findings of the paper and the last section is devoted to draw the conclusions of the paper.

Model formulation

In this section, we construct a new fractional-order SIS (susceptible-infected-susceptible)-type epidemic model with treatment control and fear effect due to the infectious disease. We separate the total populations into two time-dependable mutually exclusive classes, viz. susceptible S(t) (vulnerable for the infectious disease) and infected I(t) (already have pathogen of the disease and infectious for the population). Assume that birth rate of the susceptible population be r, natural death rate of both susceptible and infected population be a and d be the death rate due to the disease. Also, consider b as the death rate of susceptible population due to intraspecific competition (Ghirlanda et al. 2009). Note that b is smaller as compared to a, but it is not negligible for a large population. Further, let us take fear function as which reduces the growth rate of susceptible population. Again, we consider that the disease transmission function is g(I) and the disease control function is h(u). Consider m as the natural recovery rate of the infected population from the disease. Finally, assume that the recovery from infectious disease is not permanent and hence all the recovered individuals will be reinfected, i.e. they are susceptible to the infectious disease. We formulate our proposed model in the following form by the above conditions:Here is the fractional-order derivative (-order) in the sense of Caputo (Petras 2011; Podlubny 1999). In our present model, we consider fear function as with the level of fear as and the disease transmission function as , where is disease transmission rate and is half-saturation incidence rate. Further, we assume that the disease control function is . Here u is the treatment control. So our model becomesSystem (2.2) for different fractional differentiations is known as incommensurate fractional-order system. It is presented by following waySystem (2.2) is known as commensurate fractional-order system. Note that all the parameters considered in the formulations of the model are nonnegative.

Some preliminaries

In this section, we recall some essential definitions, effective lemmas for both types of commensurate and incommensurate fractional-order systems.

Definition 3.1

(Petras 2011) The Caputo-type fractional-order derivative of order for the function may be defined and denoted by:where is a space of n times continuously differentiable functions on , and be the Gamma function with (the set of all positive integers) such that . Particularly, if , the definition becomes:

Lemma 3.2

(Odibat and Shawagfeh 2007) Let us assume that and both the functions f(t) and its fractional derivative be elements of the metric space C[a, b] . If then the function f(t) is a monotone increasing, and the function is monotone decreasing if for all .

Lemma 3.3

(Li et al. 2016) Let us consider that be continuous function and satisfies the following:Then we have the inequality:for all . Here is the Mittag–Leffler function of one parameter.

Lemma 3.4

(Petras 2011) Assume the following fractional-order system of the order where and (set of all n × n real matrices). The system is said to be stable asymptotically iff satisfy for each of the eigenvalues of the matrix M and the system is said to be stable only iff for each of the eigenvalues of matrix M with the eigenvalues satisfying the critical condition must have geometric multiplicity one.

Lemma 3.5

(Petras 2011) Let us assume the fractional-order system of the order where . A stationary point of the system is called to be locally asymptotically stable iff for all eigenvalues of Jacobian matrix calculated at the corresponding stationary point.

The fractional-order system

Here, first we assume commensurate fractional-order system (2.2) by considering in the incommensurate system (2.3).

Equilibria

Now we analyse the existence criteria of nonnegative equilibria of constructed system (2.2). The system possesses the equilibrium points which are trivial equilibrium (TE) , disease-free equilibrium (DFE) and three endemic equilibria (EE) , and depending upon the positive roots , and of the following cubic equation.where ,where is known as the basic reproduction number which can be defined as the number of new infections made by an infected individual to the susceptible population during its infection period. Also, , and .

,

and

,

We rewrite the cubic equation (4.1) in following standard form of strum’s functionwhere , and

The strum’s functions of equation (4.2) may be represented asWe use the famous Strum’s theorem to verify the existence of number of distinct positive roots of equation (4.2) and hence of equation (4.1). It can be easily verified that implies equation (4.2) possesses three distinct real roots. The condition implies that . See that if , then . Now . Here the number of changes of signs of the successive strum’s functions is zero. Again . Here the number of changes of signs of the strum’s functions depends on the sign of B. If , then the number of changes sign of the successive strum’s functions is one, and hence equation (4.2) possesses exactly one positive root. If , then the number of changes sign of the successive strum’s functions is two, and hence equation (4.2) possesses exactly two positive roots. Since for any choice of parameters, and together imply that . Now implies that and implies that . Again together imply that . Now implies that . With this above discussion, we state the next theorem:

Theorem 4.1

Proposed model system (2.2) possesses at least one EE if .

Proposed model system (2.2) possesses exactly two EE points if and . Otherwise, system (2.2) possesses unique EE if , , and .

System (2.2) possesses no EE point if and .

From above Theorem 4.1, we can detect that system (2.2) possesses two endemic equilibrium points if , the basic reproduction number is less than unity, and this confirms that the system passes through a backward bifurcation at . Hence is not a sufficient condition to eradicate the disease from our system (Fig. 1).

Fig. 1

Curve for backward bifurcation of model (2.2)

The backward bifurcation diagram shows that when , there exist two endemic equilibrium points and one of which is asymptotically stable and other is unstable.

Theorem 4.2

(El-Saka 2015) Fractional-order system (2.2) undergoes backward bifurcation at iff and .

Here we choose as backward bifurcation parameter. If the condition (ii) of Theorem 4.1 is true (the conditions for the existence of two EE points), then system (2.2) possesses backward bifurcation at . This implies that system (2.2) goes through two EE points in the interval and notice that one EE is asymptotically stable and other is unstable. Here and can be calculated from .

Now the condition , i.e. implies thatAnd the condition implies thatA backward bifurcation of system (2.2) occurs at if and only if (4.3) and (4.4) are satisfied. If (4.4) holds, article by El-Saka (2015) implies that system (2.2) possesses a backward bifurcation at and also the system has two EE for the parameter in an interval. With the help of the critical value can be calculated from and which implies

Note: A similar condition of Theorem 4.2 can be found if we consider as the function of disease control parameter, i.e. .

Existence and uniqueness of solutions

In this stage, we state and prove the following lemma with the help of Li et al. (2010).

Lemma 4.3

Let us assume the following fractional-order differential equation of the order of differentiation : be a function, if the Lipschitz condition is satisfied by f(t, x) with respect to x, then system (4.5) possesses unique solution on the interval .

Proof

We consider the region to prove the existence and uniqueness criterion of the solutions of system (2.2), where , and U and V are two finite positive real numbers. Let and be two points in and set the mapping by whereAny ,Therefore, the Lipschitz’s condition is satisfied by the function E(X) with respect to the variable . Lemma 4.3, implies that system (2.2) exhibits a unique solution in respect of the initial conditions .

Accordingly we can state the next theorem:

Theorem 4.4

For any initial point system (2.2) has a unique solution for any time .

Nonnegativity and boundedness of Solutions

The nonnegative and boundedness solutions are the only useful solutions in mathematical epidemiology. Let us consider , where is the set of all positive real numbers.

Theorem 4.5

Each solution of system (2.2) starting in is nonnegative and uniformly bounded.

Nonnegativity: First consider that be an initial solution of system (2.2). Then we establish that any solution is nonnegative. Assume T be a real number such that andThe first equation of (2.2) implies that . Lemma 3.2 implies that , and this contradicts to the assumption . Hence for any , we have . Similarly, we can also prove that for all .

Uniform boundedness: To prove the uniform boundedness of the solutions, assume the succeeding function . Then,Lemma 3.3 implies that as . Hence each solution of system (2.2) starting in the region is lying in the region .

Dynamical Behaviour

We now check the stability of the TE, , the DFE point and the unique EE point . The local stability analysis of model (2.2) may be studied by the linearization technique around the each equilibrium points. The Jacobian matrix J of model (2.2) at the point (S, I) is given by:At the TE, and DFE point the Jacobian matrix areHence, the eigenvalues of system (2.2) at TE, are and . See that if , otherwise , and , where . Hence by Lemma 3.4, system (2.2) is locally stable around the trivial equilibrium asymptotically iff , i.e. if natural death rate is higher than the growth rate.

,

,

Again the eigenvalues of system (2.2) at DFE, are and . See that as (existence criterion for DFE), otherwise , and if , i.e. , where . Hence by Lemma 3.4, the system is locally stable around the DFE asymptotically iff . Therefore we can state the following theorem relating the local stability of the TE and DFE points.

Theorem 4.6

The TE, of proposed model system (2.2) is locally asymptotically stable iff .

The DFE, of model system (2.2) is locally asymptotically stable iff .

Next our objective is to analyse the local stability behaviour of the unique EE point . Characteristic polynomial of Jacobian matrix J at the EE, point is given bywhere and . Let and be the eigenvalues. Now, will be negative ifAgain, we consider Equation (4.6) as function of disease control parameter, u. Then we havewhere and . Let and be the eigenvalues. Now, will be negative ifThe eigenvalues areThe different values of and are depending on the coefficients and . Thus we can state the next theorem regarding asymptotic local stability of the EE point.

If and , then both the eigenvalues , are negative. Therefore , and hence by Lemmas 3.4 and 3.5, the EE is asymptotically stable.

If and , then one of the eigenvalues or will be nonnegative. Therefore for or 2, where and hence by Lemmas 3.4 and 3.5, the EE is unstable.

If and , then both the eigenvalues and will be complex conjugate. Therefore, . Then the EE point will be stable asymptotically if . In this case we find an interval of differentiation for as where the EE is asymptotically stable.

If and , then . Then the EE, point will be stable asymptotically if . The interval of differentiation for is , where the EE, point is stable asymptotically.

If and , then . Therefore the EE is asymptotically stable.

Theorem 4.7

The asymptotic local stability conditions for EE point of system (2.2) are followed by

If and , then the EE is asymptotically stable.

If and , then the EE is unstable.

If and , then the EE is asymptotically stable.

If and , then the EE is asymptotically stable.

If and , then the EE is asymptotically stable.

Existence criteria of Hopf Bifurcation We can rewrite our proposed model system (2.2) as:where and , is a function from to with with parameter . For EE, point of system (4.10), suppose and are complex conjugate eigenvalues of the Jacobian matrix J of model system (4.10) at (which exists if and ). System (4.10) goes through a Hopf bifurcation at critical value [by equation (4.7)] such that following singularity condition[(i)] and transversality condition [(ii)] are satisfied (Deshpande et al. 2017).Note: a similar condition of Hopf bifurcation can also be found if we consider system (2.2) as function of the disease control parameter, u.

Global stability of the equilibria

In this portion, first we only state some lemmas in order to study the global asymptotic stability criterion of DFE and the EE points.

Lemma 4.8

(Li et al. 2010) Let us consider a function which is continuously differentiable. For any and all time , the following result holds:

Lemma 4.9

(Jingjing et al. 2015) Let be a closed and bounded subset in and be a continuous function which satisfies the condition: . Then for any solution of  the model system (4.5) starting from a point always remains in . Also assume that is the largest invariant subset of the set . Then any solution x(t) of the model (4.5) starting in converges to as . In peculiar case if , then as .

Now we introduce and establish the theorems corresponding to the global asymptotic stability criterion for DFF and EE points.

Theorem 4.10

If the basic reproduction number of system (2.2) is strictly less than unity, i.e. , then DFE point is asymptotically stable globally.

Using Lemma 4.9, we prove this theorem. We formulate the following positive definite Lyapunov function:Then with the help of  the order, the constructed positive definite function V and the solution of system (2.2), we haveSo by Lemma 4.9, any solution starting in converges to the largest invariant set . Hence . Moreover, if , then system (2.2) reduces to single eqationThe solution of (4.11) isAs , we have . Therefore the equilibrium point of limit set (4.11) is asymptotically stable globally, and hence DFE point of system (2.2), is asymptotically stable globally. Hence the proof.

Theorem 4.11

The EE point, of system (2.2) is asymptotically stable globally if it is asymptotically stable locally.

We construct the following Lyapunov function in order to prove this theorem:where . It can be verified easily that the Lyapunov function L(S, I) is positive definite for all . Since EE, is a fixed point of constructed system (2.2), we can calculate the following:Thus by Lemma 4.8, it follows thatThus if . Now see that for any choice of parameters. Hence the proof.

Incommensurate fractional-order model

Here, we assume incommensurate fractional-order system (2.3), for and where , with and are relatively primes for . Then the article by El-Saka et al. (2019), the fixed point E(S, I) of model (2.3) is locally asymptotically stable iff for all the eigenvalues of the following matrixwhere and J is the Jacobian matrix calculated at that fixed point E(S, I).

Numerical simulations

To perform the numerical simulations of our model, we use modified predictor–corrector method (Diethelm and Braunschweig 2003), and we mainly emphasize the effects of the parameters: level of fear , treatment control u, the fractional-order of system (2.2) and the fractional orders and of model system (2.3). First, we choose the parameters of our proposed model system (2.2) as and we see that backward bifurcation occurs at . Figure 1 depicts that if , then there exist two EE points and one of them is stable and another is unstable which means that only the condition is not sufficient to eliminate the disease from system. Again with the help of the parameters set , we solve the system (2.2). Here as , Theorem 4.6 confirms that the TE, is asymptotically stable. Figure 2 verifies this behaviour.

Fig. 2

Dynamical behaviour around TE

Now, if we choose the parameters set as , then as and , again the Theorem 4.6 supports that the DFE, is asymptotically stable. Figure 3 depicts this behaviour.

Fig. 3

Dynamical behaviour around DFE

Again, we consider the parameter set: as . Figure 4 shows that the EE(0.707, 0.205) is asymptotically stable. Also for the different parameter set: with fractional order , Fig. 5 shows that the EE(0.508, 0.034) is asymptotically stable (green curves). But if the increases to the value , then we will see limit cycles are formed around EE for this particular order of differentiation and hence it proves the existence of Hopf bifurcation. If the further increases, then limit cycles will be bigger. Figure 5 depicts these behaviours (red and blue curves).

Fig. 4

Dynamical behaviour around EE

Fig. 5

Limit cycle for fractional order

Choose the next parameter set: with fractional-order . Figure 6 shows that the EE(0.819, 0.039) is asymptotically stable (blue curves). But if the decreases to the value , then we will see limit cycles are formed around EE for this particular and hence it proves the existence of Hopf bifurcation. Figure 6 shows this behaviour (red curves).

Fig. 6

Limit cycle for level of fear

Again, assume the parameter set: with fractional-order , Fig. 7 shows that the EE(0.703, 0.043) is asymptotically stable (blue curves). But if the u decreases to the value , then we will see limit cycles are formed around EE for this particular and hence it confirms that the system undergoes a Hopf bifurcation. Figure 7 shows these phenomena (red curves).

Fig. 7

Limit cycle for disease control parameter u

Next, we choose the parameter set: with , and using this, we solve system (2.3). Since for these parameters the DFE(0.17,0) is asymptotically stable. Figure 8 (left) confirms this behaviour. Finally, using the parameter set: and , we solve system (2.3). Since for these parameters the EE(0.7,0.13) is asymptotically stable. Figure 8 (right) verifies this.

Fig. 8

Dynamical behaviour of system (2.3) for different

Discussions

In this present paper, we have proposed and analysed a new fractional-order SIS-type epidemic model with fear effect of the infectious disease and also the effect of treatment control. The fractional-order system is more naturalistic, and it brings on more affects as compared to the ordinary system or commonly known as integer-order system. The features of most of the real-world problems (or systems) assure that their model possibly more exactly formed with the fractional-order system. A higher-order system can also be modelled to a lower-order system by using fractional calculus. In case of control theory, the systems which are modelled by fractional calculus produce much better affect than the ordinary or integer-order system of differential equations.

In our proposed model systems, i.e. commensurate system (2.2) and incommensurate system (2.3), we deduced the criterion for the existence and uniqueness of various equilibrium points of the systems in several fractional-order cases. Then we examined the conditions for nonnegativity of solutions and also uniform boundedness of that solutions. Furthermore, we established sufficient conditions to check that DFE, and the EE, of system (2.2) are asymptotically stable globally by developing suitable Lyapunov functions. Further, the theoretical results are asserting by some numerical simulations. We restrict our numerical works at the EE point of proposed systems (2.2) and (2.3).

We have studied that , the basic reproduction number of our model acts a significant purpose in local and also global stability of DFE point (Fig. 3). Also, treatment control u and the parameter of the level control have key role to local and also the global asymptotic stability of EE point (Fig. 4). We have demonstrated the criterion for the existence of Hopf bifurcation for commensurate system (2.2), against the parameter of the level of fear and also noticed that similar types of conditions can be found for the existence of Hopf bifurcation against the disease control parameter u. We have also presented the Hopf bifurcation against the parameters and u in Figs. 5, 6 and 7, respectively. In case of incommensurate system (2.3), Fig. 8 depicts the dynamical behaviour for the different fractional orders and .

In this numerical simulation, we have applied the predictor–corrector (Predict, multi-term (Evaluate, Correct), Evaluate) method Diethelm et al. (2002), Diethelm and Braunschweig (2003) and Garrappa (2010) with the help of MATLAB software. We use the implicit fractional linear multistep methods (FLMMs) as solver function for our proposed systems (2.2) and (2.3) with the fractional differential equations, and applying this method, we have described the phase portrait. Notice that is modified iteration formula of the method, PECE (Predict, Evaluate, Correct, Evaluate) and Adams–Moulton algorithm.

Conclusions

We conclude that our proposed model yields efficient results and the model may serve as a tool for a wide range of utilizations in the mathematical epidemiology. In our future works, we develop and study several fractional-order epidemic models with the help of various numerical methods like Taylor series approximations method, homotopy perturbation method, Diethelma’s method, Adomian decomposition method, etc. Also in future, we stretch SIS-type epidemic model to SIR, SEIR, SEIRS model with treatment and vaccination controls in the presence of fear of the disease and apply fractional optimal control analysis of the disease control parameters which helps us to optimize the cost used in the controlling of the disease. Moreover, this mathematical model may be applied to examine the nature and also dynamics of almost every infectious disease.