Disease Control in Age Structure Population

Abstract

We combine the Leslie model and its derivatives with the classical compartmental SIRS models to build a model of transmission of infected diseases, in a population of hosts, whether opened or closed systems. We calculate the basic reproductive rate R₀. Under certain conditions, when 𝑅₀<1, there is a disease-free equilibrium that is locally asymptotically stable. In contrast, when 𝑅₀>1, this equilibrium is unstable. Then, through an example, we show how we can define public health strategies to tackle an endemic. Finally we carry a global sensitivity analysis based on this basic reproduction rate to exhibit the most influential parameters of our model that are applied to influenza.

1. Introduction

During a certain period it was believed to have vanquished a large number of infectious diseases, until they redo their appearance, sometimes more dramatically (e.g., foot and mouth disease in the United Kingdom in 2001, highly pathogenic avian influenza pathogen in Europe in 2006) [1]. The same observation was made in Africa where, in the 70s, it was believed to have neutralized human African trypanosomiasis, long fought by the WHO until it revives later and often in the same historic foci [2]. This shows how much the prevention and control of transmission of infection diseases in the community requires constant vigilance. Many scientific disciplines have addressed this problem. Among these, the contribution of mathematical models has been very helpful. Specifically, mention may be made of introduction of the concept of 𝑅₀ and its consequences. Generally, 𝑅₀, the basic reproductive ratio, defined as the expected number of secondary infections that occur when one infective is introduced into a completely susceptible host population, is used to characterize the nature of the disease [3–9]. So, it is a threshold used to evaluate the conditions for starting an epidemic in a given area and for a specific disease. In general, when this threshold is less than 1, the disease will disappear. By cons, when it is greater than 1, the disease spreads permanently. Knowledge of the 𝑅₀ was used to evaluate the critical vaccination coverage, that is the minimum proportion of the population that must be immunized in order for the infection to die out in the population (𝑝_𝑐=1−1/𝑅₀). As an example, when 𝑅₀=3/2, the vaccination of the third of population is enough to stop the spread of the disease [10].

This threshold is even more used in most epidemiological models through compartmental representations [11–13]. In such a representation, compartments are discretized, each corresponding to all individuals in the population with specific epidemiological status (susceptible, latent, infected, removed, etc.). When the dynamics depend on time only, the balance of individuals who leave a compartment and return in different other compartments can be used in a deterministic way to model the phenomenon by a systems of ordinary differential equations (ODE) associated [12]. The spread of epidemics can also be studied in space or another parameter using a system of partial differentiation. Indeed, one can use partial differential equations (PDE) to include forms of population structure other than clinical status. These systems are studied analytically or numerically to determine the conditions under which the disease spreads in the population.

However, the proposed methods have limitations because very often they do not take into account the dynamics of the disease itself or the structure of the components of the study environment, such as the age structure. Indeed the structure may be very important for understanding the disease and its control; ignoring age-structured contacts is likely to result in misinterpretation of epidemiological data and potentially costly policy missteps [14–19]. For example, in a human population, the reproduction and the survival rates depend on the structure that we have chosen: a 20-year-old woman is more likely to procreate or to survive than a 60-year-old woman. For this reason, some mathematical tools, such as Leslie matrix (in discrete time) [20, 21] and its derivatives (in continuous time) have emerged [22]. These demographic models are still widely used in population dynamics and ecology because they are well suited to the methods of counting or census of individuals and their use involves relatively simple mathematical tools. An illustration of such a model is shown in Figure 1.

Figure 1 

Age groups with birth.

𝐴𝐺_𝑖 can be considered as individuals who are 𝑖 years old, assuming that nobody can survive after a maximum of 𝑛 years. But in general, when people's behaviors are substantially the same depending on the problem studied, several age groups will be grouped within the same single class, taking care to have disjoint classes. In this way, we have less age groups and the problem is more simplified. In this context, 𝑓_𝑖 is the rate of sexual maturation from age group 𝐴𝐺_𝑖 to age group 𝐴𝐺_𝑖+1 and 𝑏_𝑖 is the reproduction rate, assuming that there is no possible reproduction before a certain age.

In this paper, we are interested in the modeling of influenza. It is a contagious disease that is caused by the influenza virus. It attacks the respiratory tract in humans (nose, throat, and lungs). The influenza virus is transmitted in most cases by droplets through the coughing and sneezing of infected persons, but it can be transmitted as well by direct contact. Once recovered, an individual acquires temporary term protection and becomes susceptible again. For these reasons, we will use the SIRS compartmental structure to model the transmission. Moreover, since influenza can be highly contagious, particularly among persons without preexisting antibodies against influenza, such as young children during the interpandemic phase influenza and anyone during a pandemic, we will use a discretized age-class structure.

The population is divided into three groups: susceptible (𝑆), infected (𝐼), and removed (𝑅). Each group is made of two subgroups: those who can procreate and those who cannot yet procreate. We assume that the host population can vary but birth and mortality rates are independent of age within the age classes. For this reason, we use ODEs in place of PDEs, according to [23]. We further assume that inside the same age group, the population is perfectly homogeneous in terms of transmission. This transmission can vary across age group. We compute 𝑅₀ using the proposed van den Driessche method [7]. We show how to use 𝑅₀ to define some public health strategies. Finally, taking into account the very high number of parameters needed to calculate the basic reproductive rate, we do a global sensitivity analysis that highlights the most important parameters.

2. The Model

2.1. Description

Our population is subdivided into three disjoint clinical statuses. Each clinical status is then subdivided into two other classes, each class corresponding to one age group. The first age group is people who can not procreate while the second is those who may procreate. The number of hosts in each disjoint class or clinical status can vary with time. A schematic illustration of such the model is described in Figure 2.

Figure 2 

A 2 SIRS age groups model: Dynamics due to infection is represented by continuous lines, while that related to demography is shown in dotted lines. 𝑆𝑖, 𝐼𝑖, and 𝑅𝑖 correspond to susceptible, infected and removed states, respectively. Here, we have 𝐵𝑆=(1−𝜋)𝑏′𝐼2+𝑏𝑆2+𝑏′′𝑅2 and 𝐵𝐼=𝜋𝑏′𝐼2.

At any time an individual can be in one of the states Susceptible, infected, or removed. Individuals are said to be susceptible if they are not currently infected but can become infected [13]. They are infected if they can transmit the disease to others. As in the earlier models [11, 15], we use the structure of “Who Acquired Infection From Whom” WAIFW square matrix 𝛽 in which 𝛽_𝑖𝑗 is the infection rate of an individual of age class 𝑖 by another of the age class 𝑗. 𝛾 is used to represent the recovery rate. So, the infection state lasts an average period of period 1/𝛾 units of time. When an individual recovered from his illness, he enters a state of immunity and is said removed. This immunity lasts a certain period for the hosts who becomes susceptible again and the cycle repeats. The number of hosts is denoted by 𝑆_𝑖, 𝐼_𝑖, and 𝑅_𝑖 representing the number of susceptible, infected, and removed hosts within a clinical status of age group 𝑖. 𝑖=1 and 𝑖=2 are the values assigned to the first and second age groups, respectively. Taking into account the fact that the infection status of an individual may influence the birth rate, according to [24], we assume that birth rate varies with the epidemiological status: 𝑏 for susceptible hosts, 𝑏^′ for those who are infected, and 𝑏^′′ for removed people. Similarly, a transition rate is used to describe the transition from age group 1 to age group 2: for sake of generality, we consider 𝑓 for susceptible, 𝑓^′ for infected, and 𝑓^′′ for removed hosts; indeed, in animal species, worms, and snail trematode parasites of crustaceans barnacles infect both gonads and can therefore affect the rate of maturation and the birth rate [4, 5]. In addition, mortality can be applied differently depending on the status of infection and age group membership: 𝜇₁ (resp., 𝜇₂), 𝜇^′₁ (resp., 𝜇^′₂), and 𝜇₁^′′ (resp., 𝜇₂^′′) for susceptible, infected and removed hosts in age group 1 (resp. 2), respectively. Finally, we assume that it enters or leaves a constant number of susceptible hosts in the compartments 𝑆₁ and 𝑆₂ that we denote 𝑀₁ and 𝑀₂, respectively. It can be considered as a number representing a balance between input and output of the system. This model leads to the system of ordinary differential equations as follows: 𝑆̇₁=−𝛽₁₁𝑆₁𝐼₁−𝛽₁₂𝑆₁𝐼₂+𝑏𝑆₂+𝛿𝑅₁+𝑏^′′𝑅₂+(1−𝜋)𝑏^′𝐼₂−𝑓+𝜇₁𝑆₁+𝑀₁,𝑆̇₂=−𝛽₂₁𝑆₂𝐼₁−𝛽₂₂𝑆₂𝐼₂+𝑓𝑆₁+𝛿𝑅₂−𝜇₂𝑆₂+𝑀₂,𝐼̇₁=𝛽₁₁𝑆₁𝐼₁+𝛽₁₂𝑆₁𝐼₂+𝜋𝑏^′𝐼₂−𝛾+𝑓^′+𝜇^′₁𝐼₁,𝐼̇₂=𝛽₂₁𝑆₂𝐼₁+𝛽₂₂𝑆₂𝐼₂+𝑓^′𝐼₁−𝛾+𝜇^′₂𝐼₂,𝑅̇₁=𝛾𝐼₁−𝛿+𝜇₁^′′+𝑓^′′𝑅₁,𝑅̇₂=𝛾𝐼₂−𝛿+𝜇₂^′′𝑅₂+𝑓^′′𝑅₁, (1) where (𝑆₁(0),𝑆₂(0),𝐼₁(0),𝐼₂(0))≥0.

In order to meet fairly rigorous models, taking into account the assumption that birth can be contagious, a parent in 𝐼₂ can give birth to a susceptible newborn. This latter returns to 𝐼₁ and should instantly join the compartment 𝑆₁ that is (1−𝑏^′)𝐼₂. Other births with proportion 𝜋 been considered as infected for the newborn lead to the quantity 𝜋𝑏^′𝐼₂. Births from removed people are considered noncontagious. The other quantities are given following the principles of the traditional SIRS Model [12, 13].

Thus, 𝛽₁₁𝑆₁𝐼₁+𝛽₁₂𝑆₁𝐼₂ (resp., 𝛽₂₁𝑆₂𝐼₁+𝛽₂₂𝑆₂𝐼₂) corresponds to the susceptible of age group 1 (resp., age group 2) who become infected and (𝑓+𝜇₁)𝑆₁ (resp., (𝑓+𝜇₂)𝑆₂) to the output flows from this compartment. Similarly, (𝛾+𝑓^′+𝜇^′₁)𝐼₁ and (𝛾+𝜇^′₂)𝐼₂ correspond to output variations from 𝐼₁ and 𝐼₂, respectively.

2.2. The Analytical Expression of 𝑅₀

Our goal in this subsection is to give an analytical expression for 𝑅₀. The method of calculating the 𝑅₀ we use is that definition of van den Driessche and Watmough [7], which asks to reorder the system of ODEs (1) in order to show the first compartments of infected before the others. We have 𝑚=2 infected compartments: 𝐼₁ and 𝐼₂. A disease-free equilibrium (DFE) leads to the following system of equations: −𝑓+𝜇₁𝑆₁+𝑏𝑆₂=−𝑀₁,𝑓𝑆₁−𝜇₂𝑆₂=−𝑀₂. (2) The determinant of this system is given by (3) and we can note that if 𝐷 equals zero, there is an infinite number of DFE, which is irrelevant to us. For this reason, we assume throughout the rest of this paper that the parameters are defined so that 𝐷 is not equal to zero 𝐷=𝜇₂𝑓+𝜇₁−𝑏𝑓. (3)

So, a DFE has the form 𝑥₀=(0,0,𝑆^∗₀₁,𝑆^∗₀₂)^𝑡, where 𝑆^∗₀₁=𝜇₂𝑀₁+𝑏𝑀₂𝜇₂𝑓+𝜇₁,𝑆−𝑏𝑓_∗02=𝑓+𝜇₁𝑀₂+𝑓𝑀₁𝜇₂𝑓+𝜇₁.−𝑏𝑓 (4) A condition for this DFE to have a reasonable epidemiological interpretation is that sign(𝜇₂𝑀₁+𝑏𝑀₂)=sign(𝜇₂(𝑓+𝜇₁)−𝑏𝑓)=sign((𝑓+𝜇₁)𝑀₂+𝑓𝑀₁). Assuming it is the case, the Jacobian of the DFE is 𝐽(𝑥₀)=𝐹−𝑉 where 𝐹 is the part due to new infections while 𝑉 is the part coming from other transfers. We have 𝛽𝐹=¹¹𝑆^∗₀₁𝛽₁₂𝑆^∗₀₁+𝜋𝑏^′𝛽₂₁𝑆^∗₀₂𝛽₂₂𝑆^∗₀₂;𝑉=𝛾+𝑓^′+𝜇^′₁0−𝑓^′𝛾+𝜇^′₂. (5) If we set det𝑉=(𝛾+𝑓^′+𝜇^′₁)(𝛾+𝜇^′₂), because 𝛾>0, 𝛽>0 and the other epidemiological parameters are nonnegative, it follows that 𝑉⁻¹=1det𝑉𝛾+𝜇^′₂0𝑓^′𝛾+𝑓^′+𝜇^′₁. (6) For sake of simplicity, we set 𝑢=𝛾+𝜇^′₂𝑣=𝛾+𝑓^′+𝜇^′₁. (7) Then it follows that 𝐹𝑉⁻¹=1𝛽det𝑉¹¹𝑢𝑆^∗₀₁+𝛽₁₂𝑓^′𝑆^∗₀₁+𝜋𝑓^′𝑏^′𝛽₁₂𝑣𝑆^∗₀₁+𝜋𝑣𝑏^′𝛽₂₁𝑢𝑆^∗₀₂+𝛽₂₂𝑓^′𝑆^∗₀₂𝛽₂₂𝑣𝑆^∗₀₂=1det𝑉𝑀, (8) where 𝛽𝑀=¹¹𝑢𝑆^∗₀₁+𝛽₁₂𝑓^′𝑆^∗₀₁+𝜋𝑓^′𝑏^′𝛽₁₂𝑣𝑆^∗₀₁+𝜋𝑣𝑏^′𝛽₂₁𝑢𝑆^∗₀₂+𝛽₂₂𝑓^′𝑆^∗₀₂𝛽₂₂𝑣𝑆^∗₀₂. (9)

After calculating, the characteristic polynomial of the matrix 𝑀 is 𝑃_𝑥(𝑀)=𝑥²−tr(𝑀)𝑥+det(𝑀)=𝑥²−𝛽₁₁𝑢𝑆^∗₀₁+𝛽₁₂𝑓^′𝑆^∗₀₁+𝜋𝑓^′𝑏^′+𝛽₂₂𝑣𝑆^∗₀₂𝑥+det(𝑀), (10) where det(𝑀)=(𝛽₁₁𝛽₂₂−𝛽₁₂𝛽₂₁)𝑢𝑣𝑆^∗₀₁𝑆^∗₀₂−𝛽₂₁𝜋𝑢𝑣𝑏^′𝑆^∗₀₂. The discriminant for 𝑃_𝑥(𝑀)=0 is Δ=𝛽₁₁𝑢𝑆^∗₀₁+𝛽₁₂𝑓^′𝑆^∗₀₁+𝜋𝑓^′𝑏^′+𝛽₂₂𝑣𝑆^∗₀₂²−4𝛽₁₁𝛽₂₂−𝛽₁₂𝛽₂₁𝑢𝑣𝑆^∗₀₁𝑆^∗₀₂+4𝛽₂₁𝜋𝑢𝑣𝑏^′𝑆^∗₀₂. (11) It is easily shown that since all the parameters are all nonnegative, Δ≥0. Therefore, the spectral radius of 𝑀 is 𝜌(𝑀)=12𝛽₁₁𝑢𝑆^∗₀₁+𝛽₁₂𝑓^′𝑆^∗₀₁+𝜋𝑓^′𝑏^′+𝛽₂₂𝑣𝑆^∗₀₂+√Δ. (12) Finally, with 𝑢 an 𝑣 defined at (7), we obtain 𝑅₀=1𝛽2𝑢𝑣¹¹𝑢𝑆^∗₀₁+𝛽₁₂𝑓^′𝑆^∗₀₁+𝜋𝑓^′𝑏^′+𝛽₂₂𝑣𝑆^∗₀₂+√Δ. (13) From (7), we know that 1/𝑢 and 1/𝑣 are the average duration of infection of 𝐼₂ and 𝐼₁, respectively. Consequently, the first part of 𝑅₀ that is 1/2𝑢𝑣, is the geometric average of infection. The second part, 𝛽₁₁𝑢𝑆^∗₀₁+𝛽₁₂𝑓^′𝑆^∗₀₁+𝜋𝑓^′𝑏^′+𝛽₂₂𝑣𝑆^∗₀₂+√Δ can be considered as the total number of potential contacts. Seen this way, we can notice that there is an apparent similarity between the formulation of the 𝑅₀ obtained here and its standard formulation.

As shown in [7], if 𝑅₀<1 then the DFE is locally asymptotically stable, and if 𝑅₀>1 then the DFE is unstable.

3. Numerical Solutions

The results presented hereafter correspond to the data resumed in Table 1 below. The demographic data are obtained from [25, 26]. Our time unit is the day. So, the infected (resp., removed) period lasts 21 (resp., 14) days; birthrate is considered the same for susceptible, infected, and removed hosts. We are interested here in a population where women are likely to deliver early, from about 15 years; for this reason, the age group I is made of people under 15 years old. The life expectancy is 75 years at birth, so the mortality at birth is set to 1/(75∗365) and is considered the same for all those in group age 1. Taking into consideration the fact that the disease can affect the life, the mortality for infected hosts is set to doubled that of susceptible. At each time step, the balance between those who leave the system and those who return is set to 10 for susceptible in age group 1 and 50 for those in age group 2. We have 𝑆¹₀=23∗80000, 𝑆²₀=77∗80000, 𝐼¹₀=1, 𝐼²₀=0 while 𝑅¹₀=0 and 𝑅²₀=0.

Table 1 

The parameters of the simulations.

3.1. Case 𝑅₀>1

With 𝛽₁₁=(4/3), 𝛽₁₂=2, 𝛽₂₁=4, and 𝛽₂₂=0.00000008, it follows that 𝑅₀=7.0708. Figure 3 presents the results we obtained. In all the simulations we have done with 𝑅₀>1, there was always an asymptotic equilibrium. We have noticed that even for only one susceptible host in the numerical simulations, the disease spreads out.

Figure 3 

Results of numerical solutions when 𝑅0>1.

3.2. Case 𝑅₀<1

In this case, We keep the same parameters as above, except each 𝛽_𝑖𝑗 is divided by 10. It follows that 𝑅₀=0.70708. Asymptotically, this case always leads to a disease-free equilibrium for the numerical model, as presented in Figure 4.

Figure 4 

Results of numerical solutions when 𝑅0<1.

3.3. Public Health Policy

The formula (13) can further be used to explore some public health interventions. For example, with the parameters above, when 𝑅₀>1, the disease spreads out as illustrated in Figure 3. If the demographic parameters are modified in an endemic situation so as to allow a positive balance of only 2 susceptible in age group 1 and 5 others in age group 2 instead of 10 and 50, respectively, 𝑅₀ will switch from 7.0708 to 0.81547 that is the disease will disappear after a while. This case is presented in Figure 5.

Figure 5 

Numerical simulations when applying health strategies.

4. Sensitivity Analysis

In this section, we carry a sensitivity analysis based on 𝑅₀ in order to determine the most important paremeters.

4.1. Theory

We consider a generic model with independent variables 𝑋_𝑖𝑌=𝑓𝑋₁,𝑋₂,…,𝑋_𝑛. (14)𝑉 is the variance of 𝑌, 𝑉_𝑖=𝑉(𝐸[𝑌∣𝑋_𝑖]) the resulting variance of 𝑌, taken over 𝑋_𝑖, and 𝑉_𝑖𝑗=𝑉(𝐸[𝑌∣𝑋_𝑖,𝑋_𝑗]) the average of 𝑉_𝑖 over all possible point 𝑋_𝑖.

Sensitivity analysis studies how uncertainty in the output of a model (numerical or otherwise) can be apportioned to different sources of uncertainty in the model input [27, 28]. The methods of sensitivity analysis can be grouped into three classes: screening methods [29], which consist of a qualitative analysis of the sensitivity of the output variable to the input variables; the methods of local analysis, which quantitatively assess the impact of a small variation around a set of input values data; methods of global sensitivity analysis [30], interested in the variability of the model output in its entire range of variation. In this paper we opt for methods of global sensitivity analysis not only because they are interested in both qualitative and quantitative aspects but they save the effort spent to effectively carry out a series of local sensitivity analyzes. Indeed, they determine the entries responsible for this variability by using sensitivity indices known as Sobol indices [31]. Three of them are widely used: the first-order, the second-order, and the total indices.

The first-order index generally denoted by 𝑆_𝑖 quantifies the sensitivity of the output 𝑌 to the input variable 𝑋_𝑖 that is the proportion of variance in 𝑌 due to the variable 𝑋_𝑖, while the second-order index 𝑆_𝑖𝑗 expresses the sensitivity of the variance of 𝑌 to the interaction of 2 variables 𝑋_𝑖 and 𝑋_𝑗. When the number 𝑝 of input variables is too large, the number of sensitivity indices explodes (≈2^𝑝−1). The estimation and interpretation of all these indices become almost impossible. In that case, it is convenient to use the total sensitivity index denoted by 𝑆_𝑇𝑖. The mathematical formulation of these indices is given in (15) where 𝐾_𝑖 represents the set of indices containing the index 𝑖: 𝑆_𝑖=𝑉_𝑖𝑉,𝑆_𝑖𝑗=𝑉_𝑖𝑗𝑉,𝑆_𝑇𝑖=_{𝑘∈𝐾𝑖}𝑆_𝑘. (15) The higher the index will be (close to 1), the more important the variable will be.

5. Application

Formula (13) shows that the model output can be represented only by a simple function given by 𝑅₀. Therefore, the model output sensitivity is reduced to the sensitivity of this function relative to variations of 16 model parameters among 21. Indeed, 5 parameters 𝑓^′′, 𝑏^′′, 𝜇₁^′′, 𝜇₂^′′, and 𝛿 which are the output of compartments 𝑅₁ and 𝑅₂, are not involved in the formula 𝑅₀. This situation is analogous to the classical SIRS model, where the output parameters of the compartment 𝑅 which are not input to 𝑆, are not involved in the formulation of the basic reproductive ratio 𝑅₀. For this reason, we can write: 𝑅₀=𝑓(𝛽₁₁,𝛽₁₂,𝛽₂₁,𝛽₂₂,𝛾,𝑓,𝑓^′,𝑏,𝑏^′,𝜇₁,𝜇₂,𝜇^′₁,𝜇^′₂,𝜋,𝑀₁,𝑀₂). (16) To carry out the sensitivity analysis, we use the Monte-Carlo based numerical procedure for computing the full set of first-order and total-effects indices for a model of 16 factors. The method that is used here is due to Saltelli [32]. Each parameter has a standard uniform distribution 𝑈(0,1) except 𝑀₁ and 𝑀₂ that follow the uniform distribution 𝑈(−50,50). For the Saltelli model with 16 parameters, 100 000 perturbations have been generated.

The first-order sensitivity indices (𝑆_𝑖) for all 16 parameters are calculated numerically. It is obvious from Figures 6 and 7 that the single indices are very small as comparing to the total indices. Indeed, the sum of first-order effects is approximately 2.37×10⁻³, while the sum of the total indices is 4.55; as these two sums are both different from 1, this implies high importance of parameter interactions. Moreover, given that all the parameters have total indices greater than their first orders, we conclude that all of them are taking part in interactions.

Figure 6 

Single sensitivity indices for 𝑅0.

Figure 7 

Total sensitivity indices for 𝑅0.

The pictures given in Figure 8 show the values of 𝑅₀ across the 100 000 perturbations that have been generated. In these figures, an attempt at clustering shows that the infection spreads in approximately 99% of cases while only 1% corresponds to cases where it recedes. This may explain the fact that flu rife in an endemic scale. There are almost 42% that correspond to the prohibited cases in the model above, that is, the cases where 𝑅₀<0. In addition, the total indices give more information about model sensitivity and factors determining it. From Figure 7, we can distinguish three groups of parameters. As the maturation rate (𝑓), birth rate of susceptible (𝑏), and the mortality rates of susceptible individuals (𝜇₁ and 𝜇^′₁) retain the leading role with the total index of about 1, they are the most important parameters. The total indices for each of these parameters belongs to [0.95,1.03]. It is important to note that these four parameters are the only ones involved in formula (3) of the determinant at the DFE. Another group is made of the migration parameters (𝑀₁ and 𝑀₂) which are the second most influential parameters with total indices up to 0.75. However, all the 10 remaining parameters are less important but still influential, with total indices close less than 0.22 for each of them.

(a)

(b)

Figure 8 

The variability of 𝑅0.

6. Discussion

In this paper, we have coupled models that are well known and widely used not only in population dynamics but also in the transmission of infected diseases, namely, the Leslie model and its derivatives with the SIRS compartmental model. Applying our model to influenza, under certain assumptions, we have deduced an analytical expression of the basic reproductive rate 𝑅₀ for which the DFE is locally asymptotically stable when it is below 1, but unstable when it is greater than 1. From an example, we have shown how to define public health policies to influence or even eliminate a disease that spreads out. This last result has already been highlighted in [33] in the context of a vector-born transmission disease. In addition, we have made a global sensitivity analysis, allowing us to accurately determine the parameters most influencing the asymptotic dynamics represented by 𝑅₀. Finally, we were able to note that the numerical simulations confirm the analytical results we have exhibited in this paper. Another track will be to complete our model in order to extent it to more age groups and other compartmental models, such as SEIR or SEIRS [13, 34, 35]. Finally the effect of spatialization in such models can be investigated as it has been done in [33, 34].

Acknowledgment

Support from Institut de Recherche pour le Developpement and Institut de la Francophonie pour l'Informatique is gratefully acknowledged.