A Fractional-Order Infectivity and Recovery SIR Model

The introduction of fractional-order derivatives to epidemiological compartment models, such as SIR models, has attracted much attention. When this introduction is done in an ad hoc manner, it is difficult to reconcile parameters in the resulting fractional-order equations with the dynamics of individuals. This issue is circumvented by deriving fractional-order models from an underlying stochastic process. Here, we derive a fractional-order infectivity and recovery Susceptible Infectious Recovered (SIR) model from the stochastic process of a continuous-time random walk (CTRW) that incorporates a time-since-infection dependence on both the infectivity and the recovery of the population. By considering a power-law dependence in the infectivity and recovery, fractional-order derivatives appear in the generalised master equations that govern the evolution of the SIR populations. Under the appropriate limits, this fractional-order infectivity and recovery model reduces to both the standard SIR model and the fractional recovery SIR model.


Introduction
The classic SIR epidemiological model was originally introduced by Kermack and McKendrick in 1927 [1].This ordinary differential equation (ODE) system models the spread of an epidemic through a population.The SIR model was generalised in the following decade by Kermack and McKendrick to allow for age dependencies in disease transmission [2,3].Since Kermack and McKendrick, the SIR model has become widely used for modelling a range of diseases and has been extended to allow for re-infection, latent infections and the interaction of species [4,5].More recently, there has been increased interest in the extension of SIR models through the incorporation of fractional derivatives [6][7][8][9][10][11][12].
In an SIR model, the population is split into three compartments, those susceptible (S) to the disease, those infected (I), and those recovered (R) from the disease [1].In a stochastic process view of an SIR model, individuals begin in the susceptible compartment and transition into the infected compartment with some probability after coming into contact with an infected individual.Infected individuals then transition probabilistically to the recovered compartment.For the standard SIR model, these probabilities are related to exponential waiting-time densities.The model can be constructed as a directed continuous-time random walk (CTRW) through the SIR compartments [10,11,13].
Some disease processes are dependent on both the current state of the system and its history [14].The classic SIR model cannot accomodate this; however, the age-structured models of Kermack and McKendrick [2,3] can model such diseases.Fractional time derivatives, which include integrals over the function's history, can also be used to incorporate the system's history [6][7][8][9][10][11].The generalisation of an integer-order derivative to a fractional derivative is not unique, with typical examples being the Caputo derivative and Riemann-Liouville derivative [15].Typically, fractional derivatives have been incorporated into compartment models by simply replacing the integer-order derivatives with Caputo derivatives [6].Whilst such models may be able to be fitted to data, the underlying assumptions, such as the positivity of certain parameters and dimensional agreement, can be violated [13,16,17].
Fractional derivatives have recently been included in compartment models in a physically consistent way by modelling the dynamics of transitions between compartments as a stochastic process, a CTRW [10,11,13].The parameters in the resulting fractional-order equations are well-defined and consistent with the dynamics of the individuals in the population.Moreover, these fractional-order equations are guaranteed to be dimensionally consistent.
In this paper, we extend the stochastic process derivation to allow for fractional-order infectivity and fractional-order recovery.In Section 2, we derive the governing master equations of an SIR model from a stochastic process with general history-dependent infectivity and recovery.In Section 3, we consider particular forms of the infectivity and recovery such that the governing equations will contain fractional derivatives.In Section 4, we consider the limits under which the fractional-order infectivity and recovery SIR model reduce back to the classic and fractional recovery SIR models.In Section 5, we derive the steady states of the fractional-order infectivity and recovery SIR model.

Derivation
We incorporate both a fractional-order infectivity and recovery into an SIR model by deriving the master equations for a stochastic SIR model with age since infection dependences.We consider a generalised CTRW through three compartments, those susceptible (S) to the infection, those infectious (I) with the infection, and those recovered (R) from the infection.An individual is born into the S compartment.They wait a random amount of time in each compartment before moving to the next compartment.The individual may die in any compartment and be removed from consideration.Here, we derive the master equations for the time evolution of an ensemble of individuals undergoing these dynamics.
Considering an individual who has been infectious since time t , the probability this infectious individual will infect a particular susceptible person in the time interval t to t + δt is σ(t, t )δt + o(δt).The transmission rate per infected individual, σ(t, t ), is dependent on both how long the individual has been infectious, t − t , and the current time, t.If there are S(t) susceptible individuals at time t, then in the time interval t to t + δt, the expected number of new infections per infected individual will be σ(t, t )S(t)δt + o(δt).
The probability that an individual who is infected at time t is still infected at time t is given by the survival function Φ(t, t ).For an individual to become infected at time t, they must come into contact with an individual who has become infected already.The flux of individuals into the infected compartment I at time t is denoted by q + (I, t) and is therefore constructed recursively via Initial conditions are given as the number of individuals who are infected at time 0 and how long each individual has been infected.This is given by the function i(−t , 0) that represents the number of individuals that are still infected at time 0 who were originally infected at some earlier time t ; hence, Equation (1) can then be written as follows: We assume that the rate of infection, σ(t, t ), is a function of both the current time, t, and the time since infection, t − t .This accounts for both time-dependent extrinsic changes as well as the intrinsic change in the infectivity of the disease over its natural course.As such, we may write where ω(t) ≥ 0 is the extrinsic infectivity and ρ(t) ≥ 0 is the intrinsic infectivity.An individual may only leave the infected compartment by either dying or recovering from the disease.Assuming that these processes are independent, the survival function for remaining in the infectious compartment can be written as follows: Here, φ(t − t ) is the probability that an individual has not recovered and transitioned to the R compartment by time t given that they were infected at an earlier time t .Similarly θ(t, t ) is the probability that an individual has not died by time t given that they were infected at the earlier time t .We assume that the survival function of the death process takes the following form: (6) and hence, Individuals in the infected compartment at time t must have arrived in the compartment at some earlier time and not left the compartment.We can therefore express the number of individuals in the infectious compartment via the flux into the compartment and the survival function, to give the following: The function I 0 (t) gives the number of individuals who were infected at time 0 who are still infected at time t.In terms of the initial condition function, i(−t , 0), this can be written as follows: The master equations are derived by differentiating Equation ( 8).This yields where ψ(t) = − dφ(t) dt is the probability density function related to φ(t).Using Equations ( 3)-( 5), Equation ( 10) can be written as follows: A generalised master equation can be obtained by removing the dependence on q + (I, t) in the above equation.Using Equation (7), Equation ( 8) can be rewritten as As this is now in the form of a convolution, taking a Laplace transform from t to s, L{•}, then gives Again, using Equation ( 7), the first integral of Equation ( 11) can be rewritten using Laplace transforms as Making use of Equation ( 13), this becomes Here, we have defined the infectivity memory kernel as where L −1 {•} defines the inverse Laplace transform from s to t.Once again, using Equation ( 7), the third integral of Equation ( 11) can similarly be rewritten using Laplace transforms as Making use of Equation (13), this becomes L{ψ(t)}L q + (I, t) Here, we have defined the recovery memory kernel: Using Equations ( 15) and ( 18), Equation (11) becomes the master equation for the infectious compartment: This equation governs the time evolution of the number of individuals in the infectious compartment.All individuals who enter the infectious compartment must have previously been susceptible.Taking this into account, we may write the master equation for the susceptible compartment, with the addition of the vital dynamics, as where λ(t) ≥ 0 is the flux into the compartment due to births.The per capita death rate is assumed to be the same as that for the infectious compartment.Similarly, considering that individuals who enter the recovered compartment must have left the infectious compartment, we write the master equation for the recovered compartment as Together, Equations ( 21)-( 22) are the master equations for an SIR model with both time-sinceinfection-dependent infectivity and recovery.These equations are simplified by taking the initial conditions to be i(−t, 0) = i 0 δ(−t), where δ(t) is the Dirac delta function and i 0 is a constant.With these choices, we can write This leads to simplifications, and our full set of SIR master equations become Henceforth, we use the master equations with Dirac delta initial conditions for simplicity.

Fractional Infectivity and Recovery SIR
We incorporate fractional derivatives into both the infective and recovery terms by choosing ψ(t) to be power-law distributed and ρ(t) related to our choice of ψ(t).In particular, we take ψ(t) to be Mittag-Leffler distributed [18]: for 0 < α ≤ 1, where τ is a scaling parameter.This distribution has a power-law tail, such that ψ(t) ∼ t −1−α for large values of t.Here, E α,β (z) is the two-parameter Mittag-Leffler function, given by The corresponding survival function φ(t) is The Laplace transform of the recovery memory kernel, Equation (19), with Mittag-Leffer distributed ψ(t), is given by To express the inverse Laplace transform in the above equation, we consider the form of a Laplace transform of a fractional derivative [19]: for f (t) that is smooth around the origin, where is the definition of the Riemann-Liouville fractional derivative.With this, a convolution with the recovery memory kernel can be written as A fractional derivative can be incorporated into the infectivity if the infective memory kernel, Equation ( 16), has a Laplace transform similar to Equation (30).This is satisfied by taking ρ(t) of the following form: where φ(t) is defined in Equation ( 29).As we require ρ(t) ≥ 0, we must constrain α and β such that 0 < α ≤ β ≤ 1.This constraint is easily verifiable for β = α, as ρ(t) can be reduced to where ψ(t), as defined in Equation ( 27), and φ(t) are both positive functions.It is also possible to express Equation (34) using fractional derivatives as Using this form, it is clearer to see that the Laplace transform of the infectivity memory kernel becomes Using the relation between the Riemann-Liouville fractional derivative and its inverse Laplace transform, Equation (31), we are able to express the first integral of Equation (25) as Substituting Equations ( 33) and (38) into the master Equations ( 24)-( 26) yields the fractional-order infectivity and recovery SIR model:

Reduction to Classic and Fractional Recovery SIR Models
Both the classic SIR and fractional recovery model are special cases of our derived Equations (39)-(41).In this section, we consider the parameters required for the classic and fractional recovery SIR and how they relate to the generalised fractional model we have derived above.The classic SIR model can be obtained by taking constant functions for the birth-, death-and time-dependent infectivity rates, that is, λ(t) = λ, γ(t) = γ and ω(t) = ω, respectively, and by taking the limit as α, β → 1.Noting that the limit is lim we obtain the classic SIR equations: By considering the relationship between the fractional exponents α and β and Equations ( 29) and (36), we gain insight into the underlying stochastic process of the classic SIR model.For α = 1, the waiting-time function, Equation (29), reduces to an exponential function: Taking the limit β → 1 to the age of infection-dependent infectivity, Equation (36) becomes a constant: We note that this limit is independent of the form of φ(t).
In a similar fashion, we obtain the fractional recovery SIR model [10] by taking the limit β → 1 whilst leaving 0 < α ≤ 1. Making use of the limit in Equation (42) and the functional form of ρ(t) from Equation (47), we obtain While the fractional recovery SIR model can be obtained from the general fractional infectivity and fractional recovery SIR model, we are unable to obtain the fractional infectivity SIR model [11].The fractional infectivity SIR model requires α = 1 and 0 < β < 1; hence β < α violates our non-negativity conditions for ρ(t).A different form of ρ(t) was considered in [11].The form of the fractional infectivity in [11] could not readily be generalized to include a fractional recovery.Thus the model here with both fractional recovery and fractional infectivity provides an alternate form of fractional infectivity.The choice of which type of fractional infectivity model should be used could only be decided by comparisons with data.

Equilibrium State Analysis
The set of fractional-order infectivity and recovery SIR Equations (39)-(41) are a non-autonomous dynamical system because of both the history dependence of the fractional derivative and the time dependence of the parameters.To find the equilibrium states, we simplify the model by taking all time-dependent parameters to be constants, that is, λ(t) = λ, γ(t) = γ and ω(t) = ω.Hence the simplified master equations become The constant recovery rate allows us to write Equation ( 6) as For an equilibrium state (S * , I * , R * ) to exist, the following limits must exist: Taking the limit as t → ∞, Equations ( 51)-( 53) reduce to We use the result of [10] to evaluate the limit: These equations permit two distinct equilibrium states, a disease-free state: and an endemic state: For all valid values of the parameters, the disease-free state will give non-negative populations and hence be physically obtainable.The same is not true of the endemic state, which is only physically obtainable if λω (τγ) −α + 1 > τ β γ β+1 (67) In the case for which α = β = 1, the equilibrium states recover the equilibrium states of the standard SIR ODE model with vital dynamics.We expect that the endemic state will be asymptotically stable when it is physically obtainable in a similar manner to the endemic state for the fractional recovery SIR model [10].

Basic Reproduction Number
It is also possible to calculate the basic reproduction number for this model.This is defined as the expected number of individuals who will become infected from a single infectious individual in an otherwise uninfected population.This can be calculated from where N is the total equilibrium population.From Equation ( 9), with i(−t, 0) = δ(−t), we have It is then left to solve The fractional derivative of the Mittag-Leffler function is well known, and hence this can be simplified to This integral is now in a standard form and has the solution [20]: We can also rewrite the existence criterion for the endemic steady state, Equation (67), in terms of R 0 , by noting that the equilibrium population is N = λ γ , giving R 0 > 1 (73)

Summary and Discussion
In this work, we have derived a fractional-order infectivity and recovery model using a stochastic process.The fractional derivatives arise as a consequence of taking an age of infection-dependent infectivity and recovery to be power-law-distributed.In doing so, we have shown how to incorporate fractional derivatives into the model without violating the physicality of the parameters of the model.Under appropriate limits, we are able to simplify this generalised fractional model to the fractional recovery and classic SIR models.We have shown the conditions under which an endemic steady state exists.Whilst the fractional-order models here are well-defined, the fractional derivatives originate from power law assumptions, and these assumptions need to be tested by fitting to data.Some methods for fitting fractional-order equations to data have been proposed in [21], and this could be a useful starting point for further work.