The Fractional SIRC Model and Influenza A

This paper deals with the fractional-order SIRC model associated with the evolution of influenza A disease in human population. Qualitative dynamics of the model is determined by the basic reproduction number, R0. We give a detailed analysis for the asymptotic stability of disease-free and positive fixed points. Nonstandard finite difference methods have been used to solve and simulate the system of differential equations.


Introduction
Influenza is transmitted by a virus that can be of three different types, namely A, B, and C 1 .Among these, the virus A is epidemiologically the most important one for human beings, because it can recombine its genes with those of strains circulating in animal populations such as birds, swine, horses, and so forth 2, 3 .Over the last two decades, a number of epidemic models for predicting the spread of influenza through human population have been proposed based on either the classical susceptible-infected-removed SIR model developed by Kermack  The notion of fractional calculus was anticipated by Leibniz, one of the founders of standard calculus, in a letter written in 1695.In recent decades, the fractional calculus and fractional differential equations have attracted much attention and increasing interest due to their potential applications in science and engineering see 8, 9 .
In this paper, we consider the fractional order SIRC model associated with the evolution of influenza A disease in human population.Qualitative dynamics of the model is determined by the basic reproduction number, R 0 .We give a detailed analysis for the asymptotic stability of disease-free and positive fixed points.Numerical simulations are presented to verify the obtained results.

Model Derivation
There are many definitions of fractional derivatives 8, 9 .Perhaps the best known is the Riemann-Liouvile definition.The Riemann-Liouville derivative of order α is defined as where Γ is the gamma function and n is an integer.An alternative definition was introduced by Caputo as follows, which is a sort of regularization of the Riemann-Liouville derivative: The most common definition is the Caputo definition, since it is widely used in real applications.The initial conditions for the fractional order differential equations with the Caputo's derivative are in the same form as for the integer-order differential equations.The Grunwald-Letnikov GL definition is given as This formula can be reduced to where h is the time step and ω α j are the Grunwald-Letnikov coefficients defined as ω α j 1 − 1 α /j ω α j−1 , j 0, 1, 2, . . ., and ω α 0 h −α .The model presented in 5 for the spread influenza disease in the human population classifies the population in four groups or classes: S t is the proportion of susceptibles at time t individuals that do not have the virus , I t is the proportion of infected at time t individuals that have the virus and can infect , R t is the proportion of recovered at time t individuals recovered from the virus and have a temporary immunity , and C t is the proportion of cross-immune individuals at time t.One of the main assumptions of this model is that the per capita birth rate is a constant μ > 0 and the birth rate is the same as death rate.Using the above assumptions Casagrandi et al. 5 introduced the following system: where the parameter β is the contact rate for the influenza disease also called the rate of transmission for susceptible to infected, γ −1 is the cross-immune period, θ −1 is the infectious period, δ −1 is the total immune period and σ is the fraction of the exposed cross-immune individuals who are recruited in a unit time into the infective subpopulation 5, 7 .
Recently great considerations have been made to the models of FDEs in different area of researches.The most essential property of these models is their nonlocal property which does not exist in the integer order differential operators.We mean by this property that the next state of a model depends not only upon its current state but also upon all of its historical states.Now we introduce fractional order into the ODE model by Casagrandi et al. 5 .The new system is described by the following set of fractional order differential equations: where D α t is the Caputo fractional derivative.Because model 2.6 monitors the dynamics of human populations, all the parameters are assumed to be nonnegative.Furthermore, it can be shown that all state variables of the model are nonnegative for all time t ≥ 0 see, for instance, 7, 10 .Proof.The fractional derivative of the total population, obtained by adding all the equations of model 2.6 , is given by The solution to 2.7 is given by N t N 0 E α,1 −μt α μt α E α,α 1 −μt α , where E α,β is the Mittag-Leffler function.Considering the fact that the Mittag-Leffler function has an asymptotic behavior 9, 11 , One can observe that N t converges to 1 when t → ∞.Therefore, all solutions of the model with initial conditions in Ω remain in Ω for all t > 0. Thus, region Ω is positively invariant with respect to model 2.6 .
In the following, we will study the dynamics of system 2.6 .

Equilibrium Points and Stability
To evaluate the equilibrium points let Then E 0 1, 0, 0, 0 .By 2.6 , a positive equilibrium and I 1 is the positive root of g I A 1 I 2 A 2 I A 3 1 − R 0 , where R 0 β/ μ θ and

3.3
The Jacobian matrix J E 0 for the system given in 2.6 evaluated at the disease-free equilibrium is as follows: Proof.The disease-free equilibrium is locally asymptotically stable if all the eigenvalues, λ i , i 1, 2, 3, 4 of the Jacobian matrix J E 0 satisfy the following condition 12, 13 : The eigenvalues of the Jacobian matrix We now discuss the asymptotic stability of the endemic positive equilibrium of the system given by 2.6 .The Jacobian matrix J E 1 evaluated at the endemic equilibrium is given as: The characteristic equation of J E 1 is where 3.9 Following 12, 14-16 , we have the proposition.i If the discriminant of Φ x , D Φ is positive and Routh-Hurwitz are satisfied, that is, D Φ > 0,a 1 > 0, a 3 > 0, a 1 a 2 > a 3 , then E 1 is locally asymptotically stable.

Numerical Methods and Simulations
Since most of the fractional-order differential equations do not have exact analytic solutions, so approximation and numerical techniques must be used.Several analytical and numerical methods have been proposed to solve the fractional-order differential equations.For numerical solutions of the system 2.6 one can use the nonstandard finite difference method NFDM .The nonstandard finite difference schemes were introduced by Mickens in the 1980s as a powerful numerical method that preserves significant properties of exact solutions of the involved differential equation 17 .The concept of the nonstandard finite difference method is discussed in 18 .Applying this method, the system 2.6 can be discretized as follows 7 :

4.1
Doing some algebraic manipulation to 4.1 yields the following relations: The approximate solutions S t , I t , R t , and C t are displayed in Figures 1, 2, 3, 4, and 5, respectively.In each figure three different values of α are considered.When α 1, system 2.6 is the classical integer-order system 2.5 .Figure 1 indicates behavior of the approximate solutions for system 2.6 obtained for the values of α 0.9.In Figure 2, the variation of S t versus time t is shown for different values of α 1, 0.9, 0.8 by fixing other parameters.It is revealed that increase in α increases with the proportion of susceptible while behavior is reverse after certain value of time.Figure 3 and McKendrick 4 .Casagrandi et al. 5 have introduced SIRC model by adding a new compartment C, which can be called cross-immune compartment, to the SIR model.This cross-immune compartment C describes an intermediate state between the fully susceptible S and the fully protected R one.They have studied the dynamical behaviors of this model numerically 6 .Jódar et al. 7 developed two nonstandard finite difference schemes to obtain numerical solutions of a influenza A disease model presented by Casagrandi et al. 5 .Very recently Samanta 6 considered a nonautonomous SIRC epidemic model for Influenza A with varying total population size and distributed time delay.

Lemma 2 . 1 .
The closed set Ω { S, I, R, C ∈ R 4 : S I R C 1} is positively invariant with respect to model 2.6 .
depicts I t versus time t, for Figures3, 4, and 5showing the similar variations of I t , R t , and C t with various values of α.