Integro-Differential Age-Structured System for Influenza Transmission Dynamics

Influenza infection shows a wide range of severity and it is well known that a significant proportion of individuals is asymptomatic or experience mild infections. Also, It is widely accepted that influenza transmission dynamics depends on age distributions. An integro-partial differential system is considered for influenza transmission dynamics, which includes the standard Susceptible-Infected-Recovered (SIR) classes with a quarantine (Q) class and an asymptomatic class (A). In this work, we extend the previous model to an integro-partial differential model by including age-structure. We establish the existence of an endemic steady state distribution and its explicit expression. Then, an analytic expression for the basic reproduction number is obtained. Furthermore, we prove the local and global stability of the disease-free equilibrium. Some numerical simulations of the basic reproduction number have been carried out using age-dependent influenza parameter values. This study can provide effective interventions and implementing age-dependent countermeasures.


Introduction
Mathematical modeling has been a useful tool for studying the transmission dynamics and control of communicable human diseases. One of the most remarkable contributions in mathematical epidemiology is well-known as the Susceptible-Infectious-Recovered or SIR model introduced by Kermack and McKendrick [1,2]. Various infectious diseases such as measles, tuberculosis, rubella, chicken pox and influenza have been analyzed by SIR-models [3,4,5,6,7]. The inclusion of different epidemiological classes when modeling is in need of other factors, such as partial or permanent immunity after recovery, the introduction of control measures and vaccination programs. Recently, the inclusion of the quarantine-class has been a new innovation in modeling communicable diseases in order to gain a better understanding of multiple-outbreaks and emergent diseases like SARS or pandemic influenza [8,9,10]. The new model has been represented by the SIQR epidemiological model, where recurrent outbreaks are supported by this type of model [8,11,12].
Influenza is a common human disease that has several subtypes due to minor changes in its structure, called mutations (nucleotide substitutions in the HA molecule). Among all the influenza pandemics in history, the pandemic influenza A/H1N1 in 1918 is the most notorious [13]. There are three major different competing and coexisting strains of Influenza which are classified by the subtypes: A/H1N1, A/H2N2 and A/H3N2. Furthermore, research has proven that exposure and even preventive vaccination of a particular subtype does not provide any type of protection or cross-immunity against strains of a different subtype [14]. As a result of the accumulation of mutations within the influenza HA molecule, new strains of influenza A emerged within each subtype. Major genetic shifts are necessary to contribute to the formation of new subtypes. When the pathogen is changed, available treatments or vaccines may not be enough to prevent or control an epidemic outbreak, and consequently, quarantine or social distancing are effective preventive measures to control influenza [12,15,16]. The role of the quarantine-class and the asymptomatic-class has been investigated [17]. In that study, the existence of equilibrium, and the conditions for the existence of damped periodic solutions have been established.
Over the past two decades, age-structured models have been developed for different infectious diseases. Hethcote employed a demographic model with age groups under the assumption of proportional mixing between/within age-groups [18,19]. The previous study considered the age-structured SIQR model without the asymtomatic class [12]. In this manuscript, an integrodifferential system of partial differential equations is introduced to formulate a continuous age-structure SAIQR model. This age-structure model is based on the influenza dynamics established for the standard model SAIQR [17]. The existence of an endemic steady state age distribution has been assumed in order to find the basic reproduction number R 0 . We prove the local and global stability of the disease-free equilibrium when R 0 ≤ 1 by use of a Lyapunov functional. The population in consideration is divided into five epidemiological classes: S(a, t), susceptible; A(a, t), asymptomatic; I(a, t), infectious; Q(a, t), isolated (quarantined) and R(a, t), recovered individuals. The transferred proportion of individuals from the S(a, t)-class to the I(a, t)-class is denoted by p while the proportion transferred from S(a, t)-class to the A(a, t)-class is given by 1 − p. In this model, it is assumed that there is proportional mixing, which means where c(a) is the age-specific per-capita contact rate and N (a, t) is the total population of age a at time t with S(a, t) + A(a, t) + I(a, t) + Q(a, t) + R(a, t) = N (a, t).
This indicates that the population is asymptotically constant for any values of a and t. In fact, the previous relationship contributes to finding the steady state and the basic reproduction expressions for the SAIQR age-structure model. Contact rates in the heterogeneous population play an important role in epidemiological modeling for human diseases. For example, the contact rates are higher among elementary school children and much lower for elderly people, and other variables. Some assumptions and notations are introduced in order to formulate the model. The per-capita isolation quarantine rate θ(a) and µ(a) is the age-specific per-capita natural death. It was assumed that individuals from the A-class are infectious, however, they have a reduced per-capita infection rate β(a)σ, σ ∈ [0, 1]; The per-capita recovery rates A-class, I-class, and Q-class are γ 1 (a), γ 2 (a), and γ 3 (a), respectively. The measure of the infectiousness of the population at time t is the sum of the infection's age distribution I(a, t) + σA(a, t). The function λ(a, t) = β(a)λ(t) is called the force of infection. F(a) denotes the fertility per person of age a. Even though isolated (quarantined) individuals are assumed to have negligible contact with members of the overall population for a large population N , N (t) − Q(t) is a good approximation for N (t); then using the above definitions and assumptions leads to the following partial integro-differential equations age-structure system for the SAIQR influenza model: the boundary conditions at age a = 0 are System (1) with its boundary conditions corresponds to the initial boundary value problem for the age-structure SAIQR influenza model.
The following equation is obtained, by adding all the equations from System (1) with initial conditions and boundary condition where η(t) represents the births at time t and N 0 (0) = η(0). Equation (2) is known as the McKendrick-Von Foster equation and can be solved by using the method of characteristics. The solution can be expressed in the form Integro-Differential Age-Structured System for Influenza Transmission Dynamics where Equation (5) is a convolution integral equation with the following kernel Equation (5) is known as the renewal equation and is equivalent to the following expression Due to the complexity of System (1), it was assumed that the population age distribution approaches a steady state, and specific analytical expressions for the steady state, as well as the basic reproduction number, have been obtained.

Steady State Solution
The population age distribution approaches a steady state, and the population size approaches the exponential form e kt , where the parameter k is a population parameter to be evaluated (k value is the exponential growth or decay). In order to find an analytic expression for the solution of this equation, it was assumed that the population has reached an equilibrium age distribution and that N (a, t) = e kt U (a).
Substituting the partial derivatives of N (a, t) with respect to t and a, the solution for Equation (2) is A system of partial integro-differential equations is introduced to model the dynamics of the disease, a new set of re-scaled variables is defined as with the new variables, System (1) and its boundary conditions, become The following Lemma describes a nontrivial solution for the System (7) and its boundary conditions , this result will be used in the mathematical analysis to find an analytical expression for R 0 .

Proof
The standard theory of ordinary differential equations has been used to find the steady state solutions; if the time derivatives are set equal to zero, the first equation of System (7) can be written as the following separable differential equation where Λ(a) = a 0 λ(u)du; and the initial condition s 0 = 1, under the same assumptions, the second equation of System (7) becomes which corresponds to a linear differential equation, and its solution can be found using the integrating factor µ(a) = exp γ 1 (a)da . The solution for the previous equation is given bŷ with Γ 1 (a) = which is a linear differential equation with the integrating factor µ(a) = exp (γ 2 (a) + θ(a))da , and its solution is given by where Θ(a) = a 0 θ(u)du and Γ 2 (a) = a 0 γ 2 (u)du. Next, the fourth equation from System (7) becomes with the integrating factor µ(a) = exp γ 3 (a)da and the solution is given by q * (a) = e −Γ3(a) a 0 θ(y)e Γ3(y)−Γ2(y)−Θ(y) × p y 0 λ(x)e Γ2(x)+Θ(x)−Λ(x) dx dy, where Γ 3 (a) = a 0 γ 3 (u)du.

The last equation of System (7) becomes
with the solution In the next section, the number of secondary infections for the Age-Structure SAIQR influenza model will be evaluated, tedious mathematical calculations using multiple integrals were performed, the basic reproduction number was found in terms of integrals, and numerical simulations will be included in a future section.

The Basic Reproduction Number R 0 and Stability
We use the analytical expressions for the steady state solution from Lemma 3.1 under the proportionally mixing assumption in order to find an analytical expression for the basic reproduction number R 0 . This means that the contacts of a person of age a are distributed depending on the activity levels of the individuals with other ages [18], which implies λ(a) = τ c(a), at the disease free equilibrium c = 0 The following Lemma shows the analytical expression for the basic reproduction number R 0 depending on the parameters of the model. Lemma 4.1 The endemic integro-differential SAIQR age-structure model, represented by System (7) and its boundary conditions, has the following basic reproduction number R 0 with G 1 (a ) =

Proof
The basic reproduction number R 0 can be found by substituting the endemic steady solution corresponding to i * (a) andâ * (a) into the expression for λ(a, t) in System (7) where Since λ(a) = τ c(a), then the following relation follows under the assumption that λ(a) = τ c(a). The basic reproduction number R 0 can be found by setting c = 0 in the right side of Equation (22); Next, a Lyapunov function was determined to show the global stability of the disease-free steady state when R 0 ≤ 1. Since the two infectious classes in this model are I and A, and the initial system includes integrals, it seems natural to define a Lyapunov function as a linear combination of two integrals depending on I and A. The fact that the feasible set for the integro-differential System (7) consists of nonnegative fractions allowed us to find the appropriate linear combination between ∞ 0 i(a, t)da and ∞ 0â (a, t)da. Lemma 4.2 The steady state solution s * (a),â * (a), i * (a), q * (a), r * (a) corresponding to System (7) with boundary conditions, is globally stable if R 0 ≤ 1, for the Lyapunov functional where l 1 (a) and l 2 (a) are positive and bounded functions. If R 0 ≤ 1, Proof Assume that the Lyapunov function is a linear combination of the two integrals ∞ 0 i(a, t)da and ∞ 0â (a, t)da as follows: In order to determine l 1 and l 2 , we take the derivative of V with respect to t; we consider the following conditions for l 1 (a) and l 2 (a) in order to find their explicit expressions in terms of the parameters for the system: if l 1 (a) − (γ 2 + θ)l 1 (a) = 0, then the solution of this separable differential equation is given by Similarly, if l 2 (a) − γ 1 l 2 (a) = 0, then l 2 (a) = exp Γ 1 (a). Thus, if R 0 ≤ 1, thenV ≤ 0. Solutions for System (7) move toward the level set of V if they do not start on the set wherė V = 0. The feasible region withâ = 0 and i = 0 has the boundaryV = 0, also, on this boundary, dâ(a, t)/dt = 1 s and 20 Integro-Differential Age-Structured System for Influenza Transmission Dynamics di(a, t)/dt = 2 s, soâ and i move off this boundary unless s = 0. Ifâ = i = 0, then there are not asymptomatic or infectious individuals, so everybody will be susceptible. If there is a finite maximum ageâ such that there is a compact set that attracts all solutions, then all the paths in the feasible region approach the disease-free steady state because it is the only positively invariant subset of the set withV = 0. For the case R 0 > 1, given a certain age a, if s is close to 1, andâ > 0, i > 0, thenV > 0 for the points sufficiently close to the disease-free steady state, which implies that this disease-free steady state is unstable.

Numerical Simulations
We present numerical simulations to explore the effect of some age-dependent parameter values on the basic reproduction number given in (19). The effects of the inclusion of the Asymptomatic class is the main interest, hence, the most relevant parameters related with the infectious asymptomatic A-class are σ and p, where σ is the reduction constant in the transmission rate for the A-class and p is the proportion of symptomatic infected individuals. Another important parameter is the quarantine parameter (θ). Age-specific influenza parameters are taken from the previous work, which investigated the 2009 H1N1 agedependent influenza dynamics in Mexico [20]. The population is divided into six age-groups, and age-dependent parameters, including the transmission rates and the recovery rates (β i , γ 1i and γ 2i ), are varied. These parameter values are defined in Table  1. Age-specific recovery rate from A i (t)-class to R i (t)-class γ 2i Age-specific recovery rate from I i (t)-class to R i (t)-class γ 3i Age-specific recovery rate from Q i (t)-class to R i (t)-class σ Infectivity level constant of A i (t)-class p Proportion of individuals from S i (t)-class to I i (t)-class θ Proportion of quarantine from I i (t)-class For our numerical simulations, we select three different age-specific configurations, as shown in Figure 1. These age groups are based on the age-specific contact information for six age groups (0 − 5y, 6 − 12y, 13 − 19y, 20 − 39y, 40 − 59y, >= 60y) [20]. Scenario 1 assumes that the transmission rates β and the recovery rates γ1 and γ2 are constant for all the age-groups. In Scenario 2, it is assumed that the transmission rate β remains constant, but the recovery rates γ1 and γ2 are higher for the younger populations (groups 1, 2 and 3). Scenario 3 represents the higher transmission rates for the younger groups (groups 1, 2 and 3) with the same recovery rates as Scenario 2. Figure 2 illustrates the impact of the three parameter values p, σ, θ on the basic reproduction number. The results are shown as a function of transmission rates for each scenario. Regardless of different age-specific configurations, the basic reproduction number increases as the transmission rate increases. As expected, the basic reproduction number becomes smaller under Scenario 2 (higher recovery rates in the younger groups), while the basic reproduction number becomes largest under Scenario 3 (due to higher transmission rates in the younger groups). Overall, the basic reproduction number gets larger (smaller) as p or σ (θ) gets larger. Clearly, the basic reproduction number gets smaller as the level of quarantine increases (see the red triangles in the right most panels with θ = 0.3). Also, it is worth mentioning that even implementing a small level of quarantine reduces the basic reproduction number significantly (see the green circles using θ = 0.15 in the right most panels). Figure 3 shows the impact of the two parameter values p and σ on the basic reproduction number. Here, we keep the level of quarantine as θ = 0.15. As seen in Figure 2, the overall basic reproduction number increases as p and σ increase under all three scenarios. Note that the basic reproduction number can reach below one when p < 0.25 and σ < 0.1. Again, the basic reproduction gets the largest (smallest) under Scenario 3 (Scenario 2).

Conclusions
In this manuscript, we have developed an age-structure model for influenza that includes asymptomatic individuals. The inclusion of asymptomatic individuals as a new infectious class into the age-structure model can give valuable insights about the real dynamics of influenza. Also, it is pertinent to assume that the dynamics of influenza is dependent on age structures. We have obtained an analytic expression for the basic reproduction number (depending on the parameters of the model and integrals) as well as the endemic steady-state equilibrium.
It is shown that the force of infection λ becomes positive, which implies the existence of the endemic steady state solution. Moreover, we have proved that if 0 < 1, the disease-free steady state equilibrium is globally asymptotically stable by the use of a specific Lyapunov functional. One important contribution of this work is the explicit mathematical expressions for the basic reproduction number and the steady state solution for the integro-partial differential equations. Lastly, numerical simulations for the basic reproduction number are illustrated using age-dependent influenza parameters under different scenarios. Figure 1. Three configurations of age-specific transmission rates and infectious periods: Scenario 1 with constant rates in all age groups (top), Scenario 2 with higher recovery rates in the younger groups (middle), and Scenario 3 with higher transmission rates and higher recovery rates in younger groups (bottom). Six age groups are divided as 1 = 0 − 5y, 2 = 6 − 12y, 3 = 13 − 19y, 4 = 20 − 39y, 5 = 40 − 59y, and 6 = 60y or more.