Analysis and Optimal Control of Fractional-Order Transmission of a Respiratory Epidemic Model.

The World Health Organization is yet to realise the global aim of achieving future-free and eliminating the transmission of respiratory diseases such as H1N1, SARS and Ebola since the recent reemergence of Ebola in the Democratic Republic of Congo. In this paper, a Caputo fractional-order derivative is applied to a system of non-integer order differential equation to model the transmission dynamics of respiratory diseases. The nonnegative solutions of the system are obtained by using the Generalized Mean Value Theorem. The next generation matrix approach is used to obtain the basic reproduction number R 0 . We discuss the stability of the disease-free equilibrium when R 0 < 1 , and the necessary conditions for the stability of the endemic equilibrium when R 0 > 1 . A sensitivity analysis shows that R 0 is most sensitive to the probability of the disease transmission rate. The results from the numerical simulations of optimal control strategies disclose that the utmost way of controlling or probably eradicating the transmission of respiratory diseases should be quarantining the exposed individuals, monitoring and treating infected people for a substantial period. © Springer Nature India Private Limited 2019.

Introduction

The respiratory syncytial virus, influenza virus and parainfluenza virus are some viruses that cause respiratory diseases. Influenza viruses cause more or less respiratory diseases [1]. Influenza viruses are group into Type A and Type B. The viruses are often transmitted from people and cause seasonal influenza epidemics each year. The evolving prevalence of infectious diseases is increasing every year. It occurs through the respiratory tract. The spread of the disease is rapid and widespread. Since 1980, the disease has appeared just like bird flu in Asia. The outbreak of the 2014 Ebola virus in West Africa, severe acute respiratory infectious (SARS) disease and the Middle East respiratory syndrome (MERS) caused severe impact on the health system. Most respiratory diseases have no vaccines, such as SARS and Ebola. These diseases spread quickly and may be re-infected. After the flu season, there are several different types of influenza (Types A and B) and subtypes (Type A) that circulate and cause disease. For example, the bird flu virus has several subtypes such as H5N1, H7N3, H7N7, H7N9, and H9N2, which can infect humans [1].

Epidemiological mathematical models have proven to be a valuable tool for understanding, analyzing influenza virus infection dynamics, disseminating and recommending control strategies. Although [2-6] has completed a great deal of work on dynamic modeling of influenza, it is limited to ordinary differential equations. However, currently, it has been found that the use of fractional differential equations to model many different fields of phenomena has been very successful [7-18]. For instance, in mathematical epidemiology, Ebola virus epidemic has been modeled with fractional-order differential equations by [19]. They used the SEIR fractional-order model to analyze data published by the WHO to provide projections of outbreaks in three countries in West Africa. The model considered fits precisely with the real data. Their findings revealed that the outbreak will last about two years, with an estimated 9 million infected people. Although individuals heredity play an important role in mortality rate of the disease, data analysis shows that the predicted death toll is very high.

Goufo et al. [20] took into account the stability analysis of a non-linear spread Ebola hemorrhagic fever epidemic model. They used conventional time derivatives to express models that contain new parameters that happen to be fractional. They proved that the model is well defined, poses non-negative solutions and also established the conditions for boundedness. The Routh–Hurwitz criterion was used to show the existence and stability analysis of Ebola virus model equilibrium states and showed that they strongly depend on non-linear propagation. Also, they provide the conditions for the persistence of the Ebola virus in the system. In addition, numerical simulations of non-linear propagation are provided and the results obtained are significant for combating and preventing Ebola hemorrhagic fever, which has so far caused the deaths of hundreds of families and continues to infect many people in West Africa.

Fractional diffusion mimics the human mobility network by simulating disease outbreaks [21]. Human mobility networks can smooth the spread of infectious diseases from the sidewalk to the flight route. The effort of control and removal depends on describing these networks in terms of individual connections and flux rates between the contact nodes. In some cases, transportation can be parametrized by a gravity-type model or approximated by diffuse random walks. Alternatively, they separated domestic commercial air traffic into a case study of the utility of non-diffusion heavy-tailed transport models. A new stochastic simulation of typical influenza-like infections was adopted, targeting the dense, highly connected America air travel network. They show that the mobility on the network can be defined mainly by a power law, which is consistent with the previous study. They observed that the global evolution of an outbreak on the network is precisely reproduced by a two-parameter space-fractional diffusion equation, which is determined by the air travel network. The dynamic changes of cytotoxic T-lymphocyte (CTL) responses in Ebola virus infection in vivo are described by a fractional-order model with time delay [22]. They introduced a time-delay during the CTL response period to represent the time required to simulate the immune system. Stability and Hopf bifurcation were obtained from the model through fractional Laplace transform. Moreover, the stability conditions show that the dynamics of the model can be improved by the fractional-order delay.

In [23] SIR epidemic model with non-integer order under the conditions of external noise is studied. The behavior of the system changes with the introduction of seasonality and noise force. It revealed that different non-integer order and parameters improve the dynamical behavior of the system. Multi-scale fuzzy entropy is used to study the complexity of stochastic models. They designed a hard limiter control system and simulation results which showed that with effective medical and health measures the proportion of infected people can be controlled to significantly small numbers.

Compared with integer order models, non-integer order models are more effective in modeling biological systems which have long-range and temporal memories [24]. In [25], Gonzlez-Parra et al. proposed a nonlinear order system to discuss the outbreaks of H1N1 influenza. The fractional model does not only rely on the present stage but also on all its history which makes fractional models more general than the integer-order models. They also determined that the nonlinear fractional epidemiological model can be well matched to provide numerical results that are in good agreement with the actual data for H1N1 influenza. The model suggested gives valuable facts for understanding, predicting and controlling the spread of different epidemics across the globe.

For these reasons, the fact that fractional (non-integer) order models of many phenomena have recently been shown to be more realistic than the integer order, and the fact that the fractional order also has long-range memories motivated this study. In addition, it is also obvious that modeling by fractional order derivative can capture any natural phenomena or a rich variety of dynamics observed in the system. Many of the past models focused only on the analysis of the influenza outbreaks in human populations. But, a model in which the order is an integer to describe the transmission of a respiratory epidemic has recently been proposed in [26]. However, modeling the transmission of a respiratory epidemic by fractional order model might provide a feasible alternative in controlling or probably eradicating the disease.

The present study introduces Caputo fractional-order into SEIR epidemic type of model proposed by [26]. The proposed model is analyzed without control measures to investigate its stability conditions. The sensitivity of the basic reproductive number is analyzed to determine the most sensitive parameters. The interpretation of the sensitivity led us to two control measures. The optimal control theory is then used to investigate the efficiency of incorporating the control measures, namely, quarantine of exposed population and monitoring and treating of the infected population.

The remaining part of the paper is structured as follows: We provide some basic and necessary definitions about the fractional calculus in section “Fractional Calculus”. The description of the model is discussed in section “The Model”. The nonnegative solution and the stability analysis of the fractional-order differential equations are discussed in Sections “Non-negative Solutions” and “Model Equilibrium States and Stability” respectively. The optimal control of the epidemic is discussed in section “Optimal Control”. A numerical solution of the fractional model using Atanackovic and Stankovic method is discussed in section “Numerical Method”. We support our theoretical analysis with numerical simulations in section “Numerical Simulation and Discussion”. Finally, the paper concludes in section “Conclusion”.

Fractional Calculus

In this section, we define fractional integral and fractional derivative of Riemann–Liouville and Caputo respectively which are applied in this work.

Definition 1

[27]. The order of fractional integral operator function defined byis called the Riemann–Liouville fractional integral.

Here and elsewhere is known as the Euler gamma function and is defined as

Definition 2

[12]. The derivative function of fractional order (where lie in the half-open interval (0,1]), defined byis called the Caputo fractional derivative, where a is the starting point.

We use the Caputo derivative definition in this work. The initial conditions for Caputo defined fractional differential equation is the same as ordinary differential equations which is a core advantage.

The Model

Model Description

The model we studied in this work is proposed by [26]. The model considers four variables, namely, the population which is susceptible ( ), the population which is exposed (), the infected population (), and the recovered population (). According to [26], represents the rate of new susceptible people entering the population at time , represents the probability of disease transmission, represents the seroconversion rate, represents the recovery rate, represents the natural mortality rate, represents the disease induced death rate, and represents the rate at which the recovered return to the susceptible population (due to the loss of immunity). These assumptions lead to the following integer-order differential equation presented by [26], where is the effective risk population at time (or t).The system (4) in fractional order is given bywith initial conditionswhere is the Caputo fractional derivative of order . It is important to notice that when the fractional-order , the system (5) becomes the integer-order system (4).

Non-negative Solutions

Denote and let . We need the following Lemma in [28] to proof non-negative solutions of system (5).

Lemma 4.1

(Generalized Mean Value Theorem [28]). Let and    for   and thenwith

Remark 4.1

Assume and for . It follows from Lemma 4.1 that if then h(m) is non-decreasing , and if then h(m) is non-increasing .

Theorem 1

The solution of the fractional-order initial value problem given with (5)–(6) is unique and in .

Proof

Using Theorem 3.1 together with Remark 3.2 in [29], we can obtain the existence and uniqueness of the solution of (5)–(6) in . We prove that the domain for the model is positively invariant, sincearound non-negative domain or neighborhood on each hyperplane, the vector field points to .

Model Equilibrium States and Stability

For the equilibrium states of the system of fractional order model (5), letThen the equilibrium states are and whereThe diseases free equilibrium state , is where the infectives equal to zero and the endemic equilibrium state , is where the infectives is nonzero .

The basic reproduction number of the system (5) using the next generation matrix approach, given bywhere F is non-negative and V is a non-singular M-matrix.

Applying this method on the system (5), whereNow, finding F and V, at disease free equilibrium, , and using we haveTherefore, the basic reproduction number, , is the most dominant eigenvalue of , that is,The Jacobian matrix for the system (5) evaluated at is given by

Theorem 2

The disease free equilibrium state of system (5) is locally asymptotically stable if

The disease free equilibrium state is asymptotically stable locally given that the eigenvalues of satisfy the conditions [30, 31]:We can evaluate these eigenvalues by solving the following characteristic equationthis leads to the equationwhereThe characteristic equation gives the rootsObviously, , and if , then all the eigenvalues , satisfy the condition given by (9).

The basic reproduction number denoted by is value and is defined as the number of cases occurring in a population which is completely susceptible. Biologically, if is less than one, then the infection will disappear, but if it is more than one, the infection still exists.

For the discussion of the asymptotic stability of the persistence of the disease of system (5), we need the following definition and Lemma 5.1.

Definition 3

[32].The discriminant D(h) of a polynomialis defined by , where is the derivative of h. If ,R(h, f) is the determinant of the corresponding Sylvester matrix, the Sylvester matrix is formed by filling the matrix beginning with the upper left corner with the coefficients of h(y) and then shifting down one row and one column to the right side. The process is then repeated for the coefficients of f(y) .

Lemma 5.1

[32]. For the polynomial equation,the conditions displayed below make all the roots of (11) satisfy (9):

for , the condition for (11) is

for , the conditions for (11) are either Routh–Hurwitz conditions or , ;

for , if the discriminant of is positive, then Routh–Hurwitz conditions are the necessary and sufficient conditions for (11), i.e.

If ,  ,  , , then the condition (11) is satisfied. Also if , , , , then all roots of satisfies .

If , , , then condition (11) is satisfied for all .

For general m, is a necessary condition for condition (11) to be satisfied.

The Jacobian matrix calculated at the disease persistence equilibrium is given as:By using the characteristic equation , the linearized system above is in the formwhereLet D(H) be the discriminant of the polynomial H, then based on the definition 3, we obtain the discriminant of H as

Lemma 5.2

From condition (6) in Lemma 5.1 the positive equilibrium point is locally asymptotically stable, since the polynomial as given in (12) has a coefficients , , , and positives.

A Sensitivity Analysis of

Here, we investigate the response of to parameter changes and determine the effect of each parameter on and the potential for effective control and elimination of the disease. It is straightforward to calculate the partial derivatives of the value of using Eq. (8) with respect to the parameters and the recovery rate . With all other parameters held constant, the elasticity (or the variable’s normalized forward sensitivity index) approximates the fractional change in that results from a unit fractional change in parameter x, defined asThis index shows how sensitive is to changes of parameter x. Specifically, a positive (negative) index shows that an increase in the parameter value results in an increase (decrease) of [33].

The elasticities for the quantities of interest are

Figure 1 indicate that, is most sensitive to the disease transmission rate, followed by the recovery rate. The seroconversion rate and the natural death rate have the same sensitivity index. It can also be observed that, is least sensitive to , the disease induced death rate.

Fig. 1

The sensitivity analysis of the basic reproductive number

In detail, the sensitivity indexes for , , , , and , are found to be , 1, 0.0058, , and respectively, once all parameters are fixed at their baseline values (Fig. 1). Thus, for instance, if the rate of recovery were to increase or (decrease) by , then the value of would decrease or (increase) by . Likewise, a increase or (decrease) of the disease transmission rate would correspond to a decrease or (increase) of the , increase or (decrease) of seroconversion rate would decrease or (increase) the by , increase or (decrease) of the disease induced death rate would correspond to decrease or (increase) the value of by and increase or (decrease) of the natural death rate would correspond to a decrease or (increase) in the value of .

Therefore, the above interpretations recommend that control strategies that can efficiently decrease the probability of disease transmission , natural death rate , disease induced death rate , should be used to control the disease transmission effectively. Additionally, increase the rate of recovery will lead to a decrease in , thus all control strategies that can effectively help reduce the transmission of the respiratory diseases should be applied. The mathematical perspective for these strategies would be detailed in our subsequent model.

Optimal Control

In this section, we extend our model in Eq. (6) by introducing two time-dependent control measures, namely (quarantine of exposed population groups) and (monitoring and treatment of infected populations ). It is assumed that the exposed population is reduced by the factor as they are quarantine. Furthermore, the infected population is reduced by a factor of as they are monitored and treated by health professionals. The model system (6) becomeswith the given objective functionwhere E is the exposed population and I is the infected population. T is the final time and the coefficients are positive weights. Our aim is to minimize the exposed and infected population while minimizing the cost of control . Thus, we search for an optimal control , such thatwhere the control set isThe terms and represent the cost of reducing the exposed and infected population respectively, while is the cost of quarantine and also, is the cost of monitoring and treatment. The necessary conditions that an optimal control must satisfy come from the Pontryagin’s Minimum Principle [34-37]. This principle converts Eqs. (6) and (18) into a problem of point-wise minimizing a Hamiltonian M with respect to stated as follows:where , , and , adjoint variables or co-state variables [34-37]The transversality conditions are

On the interior of the control set, where , for we haveWe obtain

Theorem 3

The control parameters that minimizes over is given bywhere and are the adjoint variables satisfying (6) and the following transversality conditions: and

Numerical Method

Atanackovic and Stankovic [38] numerical method FDE is discussed in this section. This method indicates that the fractional derivative of a function with order may be stated aswherewith the following properties:By using K terms with the sums appearing in (20), we can approximate asThe above equation can be written aswhereWe setfor We can rewrite system (7) aswhereNow we can rewrite (23) and (25) as follows:with the initial conditionsThe system (27) with (28) can be solved numerically by using the Runge-Kutta fourth order method.

Numerical Simulation and Discussion

In this section, we present numerical simulations to confirm the theoretical results obtained in the preceding section.

By using the well known generalized Euler method (GEM) [39] with values in Table 1, we simulate system (5).

Table 1

Parameter values

Parameter	Discription	Value
b^σ(ξ)	Recruitment rate into the susceptible population at time ξ	0.061/day [26]
β	Probability of disease transmission	1.1/day [40]
ν^σ	Rate of seroconversion (from expose to infectious)	0.004107/day [26]
α^σ	Rate of recovery	0.7222/day [26]
κ^σ	Rate of loss of immunity	0.95 [26]
μNσ	Natural death rate	0.000024/day [26]
μDσ	Disease induced death rate	0.00000088/day [26]

For the parameter values in Table 1 above and by calculation, we obtained and the endemic equilibrium .

We obtained , and the simulations in Fig. 2 show that the endemic equilibrium is positive and locally asymptotically stable for , , and , satisfying condition 6 of Lemma 5.1. It can obviously be seen in Fig. 2 that, compared with the situation of order and , the trajectory of the system with order is nearer to the trajectory of the system with the order . Thus, the bigger of the trajectory difference, the distant from to 1.0.

Fig. 2

The figures show the trajectories for the state variables () of system (5) for different order values

It can be seen from Fig. 3 when system (5) is the classical integer-order system (4). The display of the trajectories indicates the behavior of the approximate solutions for system (5) obtained for the value of .

Fig. 3

The approximate solution of , against time for

The Effects of Optimal Control Strategies on the Exposed and Infected Populations

It can be seen from Fig. 4 that the optimal control and has a significant effect on the exposed and the infected populations respectively. Infection level is reduced rapidly but not eliminated. This suggests that monitoring and treatment strategies that can allow the immune response to rebuilding should also be well-thought-out (Figs. 5, 6).

Fig. 4

The simulations displaying the result of (quarantine of the exposed population groups) and (monitoring and treatment of the infected people) on a exposed populations, b infected populations

Fig. 5

Simulation displaying the optimal control profile and . The indicates control profile of quarantine of the exposed population’s whiles indicates control profile of monitoring and treating the infected populations

Fig. 6

This shows the phase portrait for fractional-order model with

Conclusion

We introduced Caputo fractional-order into a classical integer-order model proposed by [26] to model the transmission dynamics of respiratory disease. The nonnegative solution of the model is provided by using the generalized mean value theorem. We obtained the basic reproductive number , which perform as a threshold parameter in the disease control. We established and investigated the stability analysis of the fractional-order model with respect to the values of . The disease-free equilibrium is locally asymptotically stable if . For , using Lemma 5.1 and condition 6, we investigated the local stability of the positive endemic equilibrium state. Sensitivity analysis shows that is most sensitive to the disease transmission rate . This recommends that periodic monitoring by medical professionals and researchers should be done to control the transmission of the disease. Additionally, we investigated the optimal control problem by the application of the optimal control theory. We used the Pontryagin’s Minimum Principle to provide the necessary conditions needed for the existence of the optimal solution to the optimal control problem. Also, we applied the Atanackovic and Stankovic method to provide a numerical solution to system (5). Lastly, the theoretical results were verified by numerical simulations to measure the efficacy and impact of control on the transmission of the respiratory diseases. From the numerical simulation, the size of the infected population is significantly reduced under the controlled conditions. This proposes that if all two control measures (quarantine of exposed population groups) and (monitoring and treatment of infected populations) are employed for the same period of time and continue for a considerable period of time, the future of free from transmission of a respiratory disease could be achieved. In this manner, the fractional-order optimal control method can progress the value of the treatment (Fig. 7).

Fig. 7

This shows the phase portrait for fractional-order model with

The major advantages of our proposed fractional order model which cannot be exhibited by classical order model are:

It’s highly effective and efficient which help us to obtain better results.

It’s easy to implement.

It provides improved precision of the process model by offering more flexibility in model identification.

By modeling system by fractional order we can model a higher system by low order model.

It has the effect of memory, which is essential factor in many biological processes.