TWO MATHEMATICAL MODELS FOR H 1 N 1 INFLUENZA WITH ANTIVIRAL TREATMENT

In this paper, we establish two mathematical models to study H1N1 influenza transmission dynamics. One model is for the case of concurrent treatment, in which we assume that untreated individuals are detected at random and moved to the treatment compartment at any time of their infected phase, and the other model deals with the case of early diagnosis, in which we assume that with some probability σ ∈ [0,1], individuals are diagnosed at the moment of infection and immediately moved to the treatment compartment. Both models are analyzed including the derivation of the basic and control reproduction numbers, the proof of global stability of disease-free equilibrium points, and demonstrating how the acquired reproduction number can be used to explain the adverse effects associated with antiviral treatment. This effect is also explained using a quantity termed the total control reproduction number. We also compare the differences between the two models in evaluating outcomes of influenza. Numerical simulations are conducted to verify the theoretical analysis results.


Introduction
Influenza (also known as flu) is a respiratory disease caused by certain RNA viruses of the Orthomyxoviridae family [1].RNA viruses lack means of correcting errors introduced during replication and influenza has multiple hosts with which mutant and recombinant genotypes may have selective advantage.Influenza viruses have coexisted with humans for centuries and have historically been a cause of excessive morbidity and mortality [2].
H1N1 influenza is a type A flu virus of swine origin that was first detected in April, 2009.
The virus infected people and spread from person-to-person, sparking a growing outbreak of flu in the United States [3].An increasing number of cases are being reported internationally as well.It is well-known that novel influenza A (H1N1) virus spreads in the same way as regular seasonal influenza viruses spread, mainly through coughs and sneezes of people who are infected.At the time of the outbreak in 2009, pandemic H1N1 was a new virus, to which most people did not have immunity, and illness may have been more widespread as a result [4].
There is a need for effective treatment and control strategies in the event of pandemics.The government of China took a series of emergency measures, such as setting up work mechanism to prevent H1N1 influenza from spreading; adopting strict inspection and quarantine measures in the ports of entry and exit; strengthening the management of close contacts; monitoring and reporting infected individuals; and constantly adjusting and improving the diagnosis and treatment of cases.
One of such strategies is treatment of affected patients with antiviral medication.Currently, the drugs for H1N1 influenza are the same as those used for treating seasonal influenza.Two types of antiviral drugs currently available are neuraminidase inhibitors and M2 channel blockers [5].Neuraminidase inhibitors such as oseltamivir, zanamivir and tamiflu play an important role in inhibiting the virus replication and dissemination.M2 channel blockers such as amantadine and rimantadine can effectively inhibit replication of the virus in infected cells, block the spread of the virus infection.Despite the effectiveness of these drugs in reducing influenzarelated morbidity and mortality, the emergence of drug resistance poses a critical limitation on their application.Incidence of viral resistance to M2 channel blockers has been associated with an increasing rate in seasonal and H1N1 influenza, possibly through widespread or indiscriminate use of the drugs [6].Neuraminidase inhibitors are less prone to be selected for resistant mutations [7], and therefore offer a better option for pandemic preparedness.However, the emergence of oseltamivir resistance in 2006 has raised concerns about our preparation for an influenza pandemic [8].Recently, studies show that amantadine and rimantadine do not work among individuals who have suffered drug abuse for a long time.But, oseltamivir, zanamivir and tamiflu can effectively alleviate symptoms in some individuals.This is particularly important for preventing pandemics caused by the emergence of resistant viral mutants.Some individuals who become resistant to single drug may continue to seek treatment with other drugs to alleviate pain.After individuals become resistant to multiple drugs, antiviral treatment will no longer be effective.
In this paper, we formulate two models to study the effect of drug treatment on the prevalence of drug resistance.One considers treatment while ill, and the other considers diagnosis at the moment of infection and immediate treatment with a given probability.We compare the differences between the two models, including their connection to parameters representing treatment strategies and evaluation of the impact of drug treatment on the prevalence of resistance.We examine the role of the acquired reproduction number, in the contribution of treatment failure to the adverse effect of drug treatment.
The paper is organized as follows.Sections 2 and 3 present Model I and Model II, respectively.Analyses of the two models are also included in these sections, including the derivation of the basic and control reproduction numbers, which are shown to determine the stability of equilibria.Section 4 is devoted to numerical simulations, and Section 5 compares the differences between two models.Finally, the findings and conclusions are summarized in Section 6.

Model I
The model is based on the transmission dynamics of H1N1 influenza.Assume that, for an individual under treatment with a single drug, if drug resistance has developed, the individual may continue to receive treatment with other drugs.For those individuals, multi-drug resistance may develop.Let f 1 and f 2 denote the treatment rate of a drug-sensitive individual and individuals who are resistant to a single drug, respectively.Assume that f i = p i k i (i = 1, 2), where k i denote the rate at which an infected patient seeks a doctor, which can happen at any time of the infection, and p i denote the proportion of patient who receive treatment.
Let N denote the number of the population, which is divided into five subclasses: susceptible S; infected but sensitive to drugs I 1s ; resistant to a single drug I 1r ; multi-drug resistant strains I 2r ; and recovered R. Based on the above assumptions, we establish the following model(called (2.1) Average life-span Fraction of treated individuals who recover Average immunity period due to infection Some of the main assumptions made in the formulation of the Model I are as follows: (i) Homogeneous mixing in the populations; (ii) Individuals in the I 1r class can still receive drug treatment but there is no antiviral drug to treat I 2r individuals, and assume that β 1r > β 2r ; (iii) Transmission rate from individuals in the I 1r class is lower than that from individuals in the I 1s class due to the antiviral treatment, and so β 1s > β 1r ; (iv) Disease-induced mortality is ignored.
Note that the total population size has the following properties: It is easy to show that the following biologically-feasible region of Model I is positively-invariant and attracting.

Global stability and reproduction numbers of Model I
We derive the reproduction numbers using the next generation matrix [14].The three infected variables are I 1s , I 1r , I 2r .Let F and V denote the matrices corresponding to the new infection terms and the transitions between stages, respectively.And for simplicity we introduce the following notations: Note that The control reproduction number is where, (i) The biological interpretations of these quantities R 1sc , R 1rc and R 1m are as follows.R 1sc and R 1rc represent the numbers of secondary sensitive cases produced by a drug-sensitive and single drug-resistant strain, respectively, during the period of infection in a susceptible population.
While R 1m explains the number of secondary multiple drug-resistant cases, during the period of infection in a susceptible population, which is the basic reproduction number for the strain with multi-drug resistance.From the formulation of R 1sc and R 1rc , we know treatment can reduce infection period.
(ii) In the absence of treatment, i.e., f i = 0(i = 1, 2), we can get the basic reproduction number where The basic reproduction number R r1 of drug-resistant strain can be written as R r1 = max{R 1r , R 1m }.
R 1s denotes the basic reproduction number of drug-sensitive strains.
Theorem 1.If R 1c < 1, then the disease-free equilibrium E 10 is locally asymptotically stable, and unstable if R 1c > 1.
Proof.The Jacobian of system (2.1) at E 10 is respectively, that means all eigenvalues of J| E 10 have negative real parts, hence E 10 is locally asymptotically stable.
Recall that total population size N(t) satisfies the equation dN dt = Λ − µN and N(t) → Λ µ as t → ∞.Using results from Castillo-Chavez and Thieme [15] and Mischaikow et al. [16], we can obtain analytical results by considering the following limiting system (2.2) is locally asymptotically stable.In order to verify the global stability of E 10 , we only need to prove that E 10 is a global attractor.
From the first equation in (2.2), it follows that By the comparison principle, we have S(t) ≤ Λ µ + (S(0) − Λ µ )e −(µ+ω)t .Without loss of generality, we can assume that S(t) ≤ Λ µ .Then it follows from the second and third equations of (2.2) that      Then by the comparison principle [15], it is easy to show that I 1s (t) → 0 and I 1r (t) → 0 as From the fourth equation of (2.2), we have Since I 1s (t) → 0, I 1r (t) → 0 as t → +∞ and R 1m < 1, we have I 2r (t) → 0 as t → +∞.Similarly, from the fifth equation of (2.2) it is easy to obtain that R(t) → 0 as t → +∞.Substitution of these into the first equation of(2) gives S(t) → Λ µ as t → +∞.This implies that the disease free equilibrium E 10 is a global attractor.Then we have the following theorem.
Theorem 2. If R 1c < 1, the disease-free equilibrium E 10 is globally asymptotically stable.

More on reproduction numbers of Model I
We define R 1sc and R 1rc as drug-sensitive and single drug-resistant strain control reproduction number, respectively.R 1m denotes the basic reproduction number of multi-drug resistant strains.
Then we rewrite the expression of R 1sc and R 1rc as follows Considering the partial derivative of R 1sc and R 1rc with respect to the proportion of patients who receive treatment, which yields Clearly, R 1sc and R 1rc are decreasing functions of p i (i = 1, 2), that mean the more infected patients to seek antiviral treatment, the easier disease and the epidemic size may be controlled.
At the same time, we notice that R 1m is independent of p i (i = 1, 2), that implies if drug resistance has developed, antiviral treatment will no longer be effective.But this can not capture the nonlinear relationship with p i (i = 1, 2), which suggests that more about reproduction number need to be considered.The paper [12] gives a new idea about acquired reproduction number, which can be used to explain the adverse effects associated with antiviral treatment.When a person is infected with H1N1 drug-sensitive strain, the number of secondary cases consists of five components given below: If the total population size N 0 is sufficiently large so that the number of infected cases is relatively small, then S(t) N(t) can be closely approximated by S 0 N 0 , and S 0 N 0 can be closely approximated 1.The component 1 represents the number of new sensitive cases if the person is untreated.
The component 2 represents the number of new sensitive cases if the person is treated and recover.The component 3 represents the number of cases who have not been cured and do not continue to seek treatment.The component 4 represents the cases of single drug resistant when the person continue to seek treatment and recover.The component 5 represents the cases who continue to seek treatment and finally develop multi-drug resistant strains.We denote the sum of components 3 , 4 and 5 by R AR , which refers to the acquired reproduction number and can be used to explain the adverse effects associated with antiviral treatment.So we get Define R TC as the total control reproduction number, which is equal to the sum of both R 1sc and R AR , then we have and let Thus, the dependence of R TC on p 1 is nonlinear, and there may be a critical value p 1c such that R TC decreases for 0 < p 1 < p 1c and increases for p 1c < p 1 < 1.In addition, and Then, the dependence of R TC on p 2 is nonlinear, there may be a critical value p 2c such that R TC increases for 0 < p 2 < p 2c and decreases for p 2c < p 2 < 1.
Thus, the dependence of R TC on p i (i = 1, 2) is nonlinear, there may be two critical values p *

Model II
(3.1) The meanings of parameters µ, γ i (i = 1, 2, 3) are the same as those in Model I.
Note that the total population size has the following properties: It is easy to show that the following biologically-feasible region of Model II is positively-invariant and attracting.

Reproduction numbers of Model II
We derive the reproduction numbers using the next generation matrix [9].The four infected variables are I 1U , I 1T , I 1r , I 2r , the matrices F and V denote the matrices corresponding to the new infection terms and the transitions between stages, respectively.
Note that The control reproduction number is where and Here, R 2m represents the number of secondary infection cases produced by multi-drug resistant strains during the period of infection in a susceptible population. Let The biological interpretations of these quantities R 2U and R 2T are as follows.R 2U and R 2T represent the numbers of secondary infection cases produced by an untreated and treated case, respectively, during the period of infection in a susceptible population.

Global stability and more on reproduction numbers of Model II
As discussed in Model I, we consider more about the reproduction number to capture the nonlinear relationship with σ and show that adverse effects associated with antiviral treatment.
In some cases, the level of infection may actually get higher when treatment rate is increased.
To prove this, we consider the partial derivative of R 2sc with respect to σ .
From 0 < θ 1 < 1 and γ 2 > γ 1 , we know that R 2T < R 2U , and hence, ∂ R 2sc ∂ σ < 0. Thus, R 2sc is a decreasing function of σ .Referring to the expression of R r2 , we know that R r2 does not depend on σ .Similarly, this can not capture the nonlinear relationship between R 2sc and σ , which suggests that more about reproduction number need to be considered [12].
When a person is infected with H1N1 influenza, the number of secondary cases consists of the following four components: If the total population size N 0 is sufficiently large so that the number of infected cases is relatively small, then S(t) N(t) can be closely approximated by S 0 N 0 and closely approximated 1.The component 1 represents the cases of new infected H1N1 influenza cases who are untreated.
The component 2 represents the cases of new infected H1N1 who are treated and do not develop drug resistance.The component 3 represents the cases of new infected H1N1 influenza with drug resistant due to antiviral treatment.The component 4 represents the cases of new infected H1N1 influenza with multi-drug resistant strains due to antiviral treatment.Clearly, the quantity R 2sc is the sum of the first two components.Therefore, R 2sc underestimates the number of secondary infections by a sensitive case.For ease of reference, we denote the sum of components 3 and 4 by R AR , which refers to the acquired reproduction number, and hence, Let R TC denote the sum of both R 2sc and R AR , it follows that Obviously, R TC is a function of σ .Notice that which is independent of σ .Thus, R TC is a linear function of σ .Furthermore, R TC is always a monotone function of σ .While, paper [12] finds that this relationship can not capture the nonlinear relationship with σ , that suggests that more generations of infections need to be considered.
For the purpose of presentation, denote R TC as the second generation in which new infections are produced.Let R TC denote the number of tertiary infected cases, from the paper of Qiu and Feng [12], we know the expression of R TC as follows the first and second items in the above expression represent the numbers of tertiary sensitive and resistant cases, respectively, produced in the third generation by a typical sensitive case.The and We can show that if we have Similarly, if we get Thus, the dependence of R TC on σ is nonlinear, there may be a critical value σ c such that R TC can provide a more accurate description on how treatment can negatively impact the disease dynamics.
We observe that Model I does not separate untreated and treated individuals, and there are more parameters than that in Model II, and we can not capture the nonlinear relationship of treatment until the second generation.However, Model II separates untreated and treated individuals immediately, and we can not capture the nonlinear relationship of treatment until the tertiary generation.There are some other differences between the two models, which will be given in the section 4. Now, for convenience we denote µ + γ i (i = 1, 2, 3, 4) by µ i (i = 1, 2, 3, 4).
Theorem 3. If R 2c < 1, the disease-free equilibrium E 20 is locally asymptotically stable, and Proof.The Jacobian of system (3.1) at E 20 is .
Obviously, we can calculate that four of the eigenvalues are −µ, −(µ + ω), µ 3 (R 2r − 1) < 0, µ 4 (R 2m − 1) < 0, and the other two eigenvalues are determined by the following quadratic equation From and we can conclude that both λ 1 and λ 1 have negative real parts, thus all eigenvalues of J| E 20 have negative real parts, hence E 20 is locally asymptotically stable.
The total population size N(t) satisfies the equation dN dt = Λ − µN and N(t) → Λ µ as t → ∞.Similar to the proof of Theorem 2.2, we can obtain the following limiting system If R 2c < 1, then R 2sc < 1 and R 2m < 1.From Theorem 3.1, the disease-free equilibrium E 20 is locally asymptotically stable.In the following, we only need to prove that E 20 is a global attractor.
It follows from the first equation of (4) that By the comparison principle, we have S(t) ≤ Λ µ + (S(0) − Λ µ )e −(µ+ω)t .Without loss of generality, we can assume that S(t) ≤ Λ µ .Then it follows from the second and third equations of (3.2) that      Based on the comparison principle [15], it is easy to show that I 1U (t) → 0 and I 1T (t) → 0 as From the fourth and fifth equation of (3.2), we get Since I 1U (t) → 0, I 1T (t) → 0 as t → +∞ and R 2m < 1 holds, we obtain I 1r (t) → 0, I 2r (t) → 0 as t → +∞.Similarly, from the sixth equation of (3.2), it is easy to get that R(t) → 0 as t → +∞.
Substitution of these into the first equation of (3.2) gives S(t) → Λ µ as t → +∞.This implies that the disease-free equilibrium E 20 is a global attractor.Thus, we obtain the following theorem.
Theorem 4. If R 2c < 1, the disease-free equilibrium E 20 is globally asymptotically stable.

Numerical simulations
In this section, we present some numerical simulation results, which confirm or extend the analytic results and illustrate the effect of new infections with sensitive strain and treatment ( f i (i = 1, 2)) on controlling the infection.
We consider the situation in which the population size has reached the steady-state Λ µ = 10 5 .Assume that the life span is 1 60 , then µ ≈ 0.00005.According to [13], the immunity obtained from infection about one year, so we choose ω = 1 365 ≈ 0.003.Estimating of the basic reproduction number for H1N1 influenza ranges from 1.75 to 3.96 in 2009-2010 year [14].In this paper, we set R 10 = 1.75.The baseline transmission coefficient (β 1s ) for drug-sensitive case can be calculated from the formula for R 10 , which gives β 1s = 0.3.During H1N1 influenza, several percentage of individuals who receive adamantanamine treatment may develop drugresistant and when drug-resistant strain develops, some individuals can not recover, they will continue to seek other drugs for treatment, so we choose reasonable c 1 and c 2 .We assume an average period of infection to be 6 days, so the recovery rate γ i (i = 1, 2, 3) = 0.1667 [12].The basic reproduction number of the drug-resistant strain R r1 can be either greater or smaller than that of the sensitive strain R 1s .Here, we consider the case R r1 R 1s < 1, and assume R r1 = 0.9 * R 1s = 1.5, those parameter values are summarized in Table 2.
Firstly, we only consider the Model I, and perform some simulations.When the treatment rate is greater than 0.12, the basic reproduction number is less than one, which implies the control measures can effectively prevent disease outbreaks.So we consider the case that the treatment rate is less than 0.12.Firstly, we consider f 2 = 0 and vary the treatment rate f 1 .So we consider the case that the treatment rate is less than 0.4.With the increase of treatment rate σ , infected sensitive cases decrease and the time of peak arriving can be delayed.However, both the number of infected single and multi-drug resistant cases increase, that can be explained for adverse effects on treatment.Furthermore, there may appear rebound phenomena among infected multi-drug resistant.

Comparison between Model I and Model II
We now discuss the difference between two models.

Discussion and Conclusion
In this paper, we study two mathematical models of H1N1 influenza transmission dynamics, that incorporate drug-sensitive and drug-resistant strains.The main purpose of this study is to examine the impact of antiviral treatment on the prevalence of drug resistance.A detailed stability analysis of the disease-free equilibria, and derivation of the basic and control reproduction number are presented.We also Demonstrate how the acquired reproduction number can be used to explain the adverse effects associated with antiviral treatment.This effect is interpreted using a quantity termed the total control reproduction number.
One of the interesting findings is that, despite the key role that control reproduction numbers play, they do not provide sufficient measures for examining the effect of antiviral treatment or reflect the adverse effects of treatment.To reveal the role of antiviral treatment in controlling H1N1 influenza, we need to consider the acquired reproduction number and the total control reproduction number.
The two models considered in this study have several different features, some of which allow for unique insight into the disease dynamics including the influence of drug resistance on the reproduction numbers.For example, in Model II, we derived the quantity R [3] TC for the control reproduction number, which has different properties than the conventional control reproduction number due to the contribution of acquired resistance R AR .We showed that, under certain conditions, there may exist a threshold value σ * such that R TC decreases with σ for σ < σ * but increases with σ for σ > σ * .This suggests that the number R TC may provide a better quantity than R 2sc for elevating the effect of treatment on the level of infection.
We develop two models with different assumptions about the timing of treatment.It is important to examine the critical assumptions and better understand their possible impact on model outcomes.It often happens that, when a model is formulated, certain assumptions are made without consideration of their consequences.In this paper, we use the exponential waiting time in disease stages.That is, the survival probability is described by a negative exponential function.Due to the memoryless property of exponential distribution, this leads to ODE models that are easy to be analyzed.But models with exponentially distributed infectious stage can result in misleading or incorrect evaluations of effectiveness.In the future, we will pay further attention on the impact of non-exponentially distributed disease stages to epidemic control.

2 )
denote the rate at which an infected patient seeks a doctor, that can happen at any time of the infection.Next, we study the model dealing with the case of early diagnosis, in which we assume that with some probability σ ∈ [0, 1] individuals are diagnosed at the moment of infection and immediately moved to the treatment compartment.As usual the population is divided into susceptible S, infected I and recovered R. Here, the infected are further subdivided into untreated individuals I 1U , individuals I 1T who receive treatment, individuals I 1r who are resistant to a single drug and multi-drug resistant strains I 2r .Based on the above assumptions, we develop the following model(called Model II) TC decreases for 0 < σ < σ c and increases for σ c < σ < 1.The key difference between R [3] TC and R TC in terms of their functional relationships with σ suggests that R [3]

Fig 2 :
Fig 2: The change of infected drug sensitive, single and multi-drug resistant individuals with different treatment rate f 1 .

Fig 3 :
Fig 3: The change of infected drug sensitive, single and multiple drug-resistant individuals with different treatment rate f 2 .

Fig 4 :
Fig 4: The change of I 1U + I 1T , I 1r and I 2r with treatment rate σ = 0.

Fig 5 :
Fig 5: The change of infected drug sensitive, single and multi-drug resistant individuals with different treatment rate σ .

Fig 7 : 1 From
Fig 7: The proportion of cumulative incidence in the total population with treatment rate σ = 0.1

Table 1
Definitions of frequently used symbols

Table 2
Parameter values for system(2.1)