Global dynamics of a two-strain flu model with a single vaccination and general incidence rate

In this his paper, we studied the global dynamics of a two-strain flu model with a single-strain vaccine and general incidence rate. Four equilibrium points were obtained and the global dynamics of the model are completely determined via suitable lyapunov functions. We illustrate our results by some numerical simulations.


Introduction
Seasonal influenza is an acute respiratory infection caused by influenza viruses which circulate in all parts of the world. Worldwide, these annual epidemics are estimated to result in about 3 to 5 million cases of severe illness, and about 290 000 to 650 000 respiratory deaths [1]. This infection can have an endemic, epidemic or pandemic behavior.
There were, three major flu pandemics during the 20th century, the so called Spanish flu in 1918 had, been the most devastating pandemic. It has been estimated that the Spanish flu claimed around 40-50 million deaths (as much as 3 % of the total population), and it also infected 20-40% of the whole population. In 1957-1958, the Asian flu or bird flu pandemic caused more than two million deaths [2]. Unlike the Spanish flu, this time the infection-causing virus was detected earlier due to the advancement of science and technology. A vaccine was made available but with limited supply. After a decade (in 1968), a flu pandemic that originated again from Hong Kong hit mankind. That flu pandemic also claimed one million lives. Beside these three major ones, there are some other flu pandemics spreading among nations on smaller scales. For instance, the 2009 H1N1 swine flu is one of the more publicized pandemics that attracted the attention of all scientists and health professionals in the world and made them very much concerned. The pandemic, however, did not result in great casualties like before. As of July 2010, only about 18,000 related deaths had been reported [2]. There are many methods of preventing the spread of infectious disease, one of them is vaccination. Vaccination is the administration of agent-specific, but relatively harmless, antigenic components that in vaccinated individuals can induce protective immunity against the corresponding infectious agent [3].
Influenza causes serious public-health problems around the world, therefore, we need to understand transmission mechanism and control strategies. Mathematical models also provided insight into severity of past influenza epidemics. Some models were used to investigate the three most devastating historical pandemics of influenza in the 20th century [4][5][6]. There are a lot of pathogens with several circulating strains. The presence of them is mostly due to incorrect treatment.
Rahman and Zou [2] proposed a two-strain model with a single vaccination, namely.
Global dynamics of a two-strain flu model with a single vaccination and general incidence rate
• 1 µ is the average time of life expectancy.
• r is the rate of vaccination with strain 1.
• k is the transmission coefficient of vaccinated individuals to strain 2.
• β1 is the transmission coefficient of susceptible individuals to strain 1.
• β2 is the transmission coefficient of susceptible individuals to strain 2.
is the average infection period of strain 1.
is the average infection period of strain 2.
• v1 is the infection-induced death rate of strain 1.
• v2 is the infection-induced death rate of strain 2.
The incidence rate of a disease measures how fast the disease is spreading and it plays an important role in the research of epidemiology. Rahman and Zou [2] used the bilinear incidence rate βSI. However, there are more realistic incidence rates than the bilinear incidence rate, For instance, Capasso and his co-workers observed in the seventies [7] that the incidence rate may increase more slowly as I increases, so they proposed a saturated incidence rate βIS 1+ζI . Baba and Hincal [8] studied an epidemic model consisting of three strains of influenza (I1 , I2, and I3) where we have vaccine for strain1 (V1) only, and force of infection βSI 1+ζS for strain 2. Baba et al. [9] studied an studied an epidemic model consisting of two strains of influenza (I1 and I2) where force of infection βSI 2 1+ζI 2 2 for strain 2.
We propose to study model (1) modifying the force of infection in the compartments I1 and I2, by extending the incidence function to a more general form: Which is based on the incidence rate studied in [10].
Thus, the resulting model is given by the following system: Whose state space is R 5 We make the following hypotheses on Fi, i = 1, 2.:

H3) lim
Ii exist and is positive for all S > 0.
The first of this hypotheses is a basic requirement for any biologically feasible incidence rate, since the disease cannot spread when the number of susceptible or infected individuals is zero. As for (H2), the condition ∂f i ∂S (S, Ii) > 0 ensures the monotonicity of fi(S, Ii) on S, is non-increasing with respect to Ii . In the case when fi monotonically increases with respect to both variables and is concave with respect to Ii , the hypothesis (H2) naturally holds. Concave incidence functions have been used to represent the saturation effectin the transmission rate when the number of infectives is very high and exposure to thedisease is virtually certain.
(H3) is needed only to ensure that the basic reproduction number is well defined. Some examples of incidence functions studied in the literature that satisfy (H1)-(H3) are as follows: (C1) F(S,I)=βSI [2].

5.
A more thorough list can be found in [10].
This paper is organized as follows. In section 2, we describe the disease dynamics described by the model. In section 3, we calculate the basic reproduction number. In section 4, we establish the existence of equilibrium points. In section 5, we study the stability of the model. In section 6, provides some numeric simulations to illustrate our main theoretical results. The paper ends with a some remarks.  (2), the total population N = S + V1 + I1 + I2 + R satisfies: The comparison theorem then implies that lim and so are all components S(t), V1(t), I1(t), I2(t) and R(t).
Since the equation forṘ is actually decoupled from the rest in equation (2), we only need to consider dynamics of the following four-dimensional sub-system: Therefore, we have established the following.

Basic reproduction number
The basic reproduction number of infection of model (3), is a dimensionless quantity denoted by R0, and intuitively defined as the expected number of secondary infection cases caused by a single typical infective case during its entire period of infectivity in a wholly susceptible population. Then, referring to the method of [11]. Then , 0, 0 . The matrix F is non-negative and is responsible for new infections, while the V is invertible and is referred to as the transmission matrix for the model (3). It follows that, Where σi = ∂Fi(S0, 0) ∂Ii , for i = 1, 2. Thus, the basic reproduction number can be calculate as Where ρ(A) denotes the spectral radius of a matrix A. Let Global dynamics of a two-strain flu model with a single vaccination and general incidence rate 7.

Existence of equilibrium solutions
The four possible equilibrium points for the system (3) are: Disease-free equilibrium, single-strain (I1)-infection, single-strain (I2)-infection and endemic equilibrium. The system (3) has disease-free equilibrium E0 = Λ λ , rΛ µλ , 0, 0 for all parameter values. We will now prove the existence of the others equilibrium points. First we will show some lemmas.
Proof. By  Proof.
As Ω is a positively invariant set for model (4), it will be enough to show that if   if and only if R2 > 1. Also, if −α2rµ − α2µ 2 + kΛr < 0 then E2 is unique. While if −α2rµ − α2µ 2 + kΛr > 0 then the model (3) has at most one single-strain Proof. 1) If I2 = 0 and R1 > 1, we consider the system By (5) and (6) Global dynamics of a two-strain flu model with a single vaccination and general incidence rate
Note that S ≥ 0 if and only ifĪ1 ≤ Λ α 1 .Ī1 being determined by the positive roots of the equation.
Next, we shall show thatĪ1 is unique. From (6), it follows that Using (H2) and lemma 2, we have that −α 1 Which implies that G(Ī1) strictly decreases at any of the zero points of (7). Let us suppose that (7) has more than one positive root. Without loss of generality, we choose the one, denoted byĪ1 * , that is the nearest toĪ1. Because of the continuity of G(Ī1), we must have G ′ (Ī1 * ) ≥ 0, which results in a contraction with the strictly decreasing property of G(Ī1) at all the zero points.
Note thatS ≥ 0 if and only ifĨ2 ≤ Λ α 2 .Ĩ2 being determined by the positive roots of the equation.
See that Then Therefore H(Ī2) > 0 by I2 sufficiently small. Also then equation (11) has a positive root. Also if E2 exists then Note thatS < S 0 andṼ1 < V 0 1 . Then by lemma 2 and remark 1 Then R2 > 1. Next, we shall show thatĨ2 is unique if −α2rµ − α2µ 2 + kΛr < 0 and if −α2rµ − α2µ 2 + kΛr > 0 then the model (3) has at most one single-strain (I2)-infection in Furthermore, it can be found that Global dynamics of a two-strain flu model with a single vaccination and general incidence rate

13.
If −α2rµ − α2µ 2 + kΛr < 0, then H ′ (Ĩ2) < 0 which implies that H(Ĩ2) strictly decreases at any of the zero points of (11). Let us suppose that (11) has more than one positive root. Without loss of generality, we choose the one, denoted bỹ I2 * , that is the nearest toĨ2. Because of the continuity of H(Ĩ2), we must have If −α2rµ − α2µ 2 + kΛr > 0, Let us suppose that (11) has more than one positive root in Without loss of generality, we choose the one, denoted byĨ2 * , that is the nearest toĨ2. Note that H ′ (Ĩ2 * ) < 0 and Proof. Similar argument to the proof of Theorem 1 proof that the model (3) Then H ′ (Ĩ2) < 0 which implies that H(Ĩ2) strictly decreases at any of the zero points of (11). Let us suppose that (11) has more than one positive root. Without loss of generality, we choose the one, denoted byĨ2 * , that is the nearest toĨ2. Because of the  The model (3) can have endemic infection equilibrium E3 = (S * , V * 1 , I * 1 , I * 2 ). To find E3, we consider the system By (13), (14) and (15) Substituting in (12).
.Ī2 being determined by the positive roots of the equation.
I * 1 being determined by the positive roots of the equation.

Stability of equilibrium
In this section we will study the local and global stability of the equilibrium points.
Global dynamics of a two-strain flu model with a single vaccination and general incidence rate 15.
Proof. The Jacobian matrix of the model, we get as follows: Then Eq. (16) at the disease-free equilibrium E0 is Thus the eigenvalues of the above Eq. (17) are From (18), if R0 < 1, then λ3, λ4 < 0 and we obtain that the disease-free equilibrium E 0 of Model (3) is locally asymptotically stable. If R0 > 1, then the disease-free equilibrium loses its stability.
Proof. Then Eq. (16) at the equilibrium E1 is Where The last equality regarding A33 is that equation (6) implies that f1(S,Ī1) − α1 = 0. The corresponding characteristic polynomial is Then has an eigenvalue is A44 and the remaining ones satisfy Global dynamics of a two-strain flu model with a single vaccination and general incidence rate
The corresponding characteristic polynomial is Then (20) has an eigenvalue equal to B33 and the remaining ones satisfy Note that Then b1, b0 > 0. Also Global dynamics of a two-strain flu model with a single vaccination and general incidence rate
Proof. The result follows from an application of Theorem 4.6 in [12], with X1 = int(R 4 + ) and X2 = bd(R 4 + ) this choice is in accordance by virtue of Lemma 1 there exists a compact set Ω in which all solution of system (3) initiated in R 4 + ultimately enter and remain forever after. The compactness condition C4.2 is easily verified for this set Ω1. Denoting the omega limit set of the solution x(t, x0) of system (3) starting in x0 ∈ R 4 + by w(x0). Note that w(x0) is bounded (Lemma 1), we need to determine the following set: From the system equations (3) it follows that all solutions starting in bd(R 4 + ) but not on the I1 axis or I2 axis leave bd(R 4 + ) and that the axes I1 and I2 are invariant sets, which implies that Y2 = (S, V1, I1, I2) ∈ bd(R 4 + )|I1 = 0 or I2 = 0 . Therefore Ω2 = {E0, E1, E2}, then 3 i=1 {Ei} is a covering of Ω2, which is isolated (since Ei (i = 1, 2, 3) is a saddle point) and acyclic. It will be enough to show that Ei (i=1,2,3) is a weak repeller for X1.
Further, it is proved in [13] uniform persistence implies the existence of an interior equilibrium point. Therefore, we have established the following.
Global dynamics of a two-strain flu model with a single vaccination and general incidence rate 21.
Then E3 is locally asymptotically stable.

Global stability of equilibria
In this section, we study the global properties of the equilibria. We use Lyapunov function to show the global stabilities. Such Lyapunov functions all take advantage of the properties of the function. g which is positive in R+ except at x = 1, where it vanishes.
Theorem 9. The DFE E0 is globally asymptotically stable if, Global dynamics of a two-strain flu model with a single vaccination and general incidence rate
Global dynamics of a two-strain flu model with a single vaccination and general incidence rate 25.
Note that ThenV ≤ 0. Furthermore, dV dt = 0 if and only if S =S and I1 =Ī1, which implies that S →S, I1 →Ī1 and I2 → 0 as t → ∞. By LaSalle's invariant principle, this implies that all solution in Ω1 approach the plane S =S, I1 =Ī1 and I2 = 0 as t → ∞. Also, All solution of (3) contained in such plane satisfyV1 = rS − µV1, which implies that V1 → rS µ =V1 as t → ∞, that is, all of these solution approach E1. Therefore we conclude that E1 is globally asymptotically stable in Ω1.
Global dynamics of a two-strain flu model with a single vaccination and general incidence rate
Now we will show that S(t) →S, V1(t) →Ṽ1 and I1(t) →Ĩ1 Consider the Lyapunov function Now, we need to showV < 0. Note that We concludeV < 0. Therefore E3 is globally asymptotically stable.
Global dynamics of a two-strain flu model with a single vaccination and general incidence rate 29.
By the relation of geometric and arithmetic means, we concludeV < 0. Therefore E3 is globally asymptotically stable.

Numerical simulations
In this section, we present some numerical simulations of the solutions for system (3) (3) indicates that E 3 is globally asymptotically stable. 9, we see that the disease-free equilibrium E0 is globally asymptotically stable. Numerical simulation illustrates our result (see Fig. 1).
• Example 6.4. In system (3) By theorem 12, we see that the E3 is globally asymptotically stable. Numerical simulation illustrates our result (see Fig. 5).

Concluding remarks
In this paper,we studied a system of ordinary differential equations to model the disease dynamics of two strains of influenza with only one vaccination for strain 1 being implemented, and general incidence rate for strain 1 and strain 2. We obtained four equilibrium points: • E0 disease free equilibrium, I1 and I2 are both zero.
We have investigated the topics of existence and non-existence of various equilibria and their stabilities. We also used next generation matrix method to obtain two threshold quantities R1 and R2, called the basic reproduction ratios for strain 1 and 2 respectively. It was shown that the global stability of each of the equilibrium points depends on the magnitude of these threshold quantities. More precisely, we have proved the following: • If R0 < 1 the disease free equilibrium E0 is globally asymptotically stable. If R0 > 1, then E0 is unstable.

35.
In order to discuss the meaning of our mathematical results, let us rewrite the two key indirect parameters R1 and R2 in terms of the direct model parameters as shown below: Also the derivative of R2 with respect to r is, Note that R1(r) is decreasing and R2(r) depends on ∂f 2 Λ µ ,0 ∂S . Now we will analyse some cases of incidence rate.  Note that for (C1) and (C3) R2(r) is increasing if βi < k, R2(r) is decreasing if βi > k and R2(r) is constant if βi = k. For (C2) R2(r) is increasing if βi ≤ k (ζ = 0). If βi > k R2(r) is increasing if ζ i kΛ β i −k − µ < r and decreasing if ζ i kΛ β i −k − µ > r. Also for if the force of infection of strain 1 is (C2), then R1 = β 1 α 1 (1+ζ 1 S 0 ) , note that R1 is decreasing in ζ1. If the force of infection of strain 2 is (C2), then R2 = β 2 Λ α 2 (λ+ζ 2 Λ) + krΛ α 2 µλ , note that R2 is decreasing in ζ2.
With the above information and the results in Section 5, the vaccination is always beneficial for controlling strain 1, its impact on strain 2 depends on the force of infection. If the forced of infection of strain 2 is (C2), the impact of vaccination depends of values of β2, k and ζ2. If ζ2 = 0; if β2 > k it plays a positive role, and if β2 < k, it has a negative impact in controlling strain 2. This is reasonable because larger k (than β2) means that vaccinated individuals are more likely to be infected by strain 2 than those who are not vaccinated, and thus, is helpful to strain 2. Smaller k (than β2) implies the opposite. If ζ2 = 0; if β2 > k, it plays a positive role and if β2 < k, not necessarily has a negative impact in controlling strain 2, because R2(ζ2) is decreasing, i.e., for ζ2 sufficiently large it can play a positive role. This is reasonable because larger k (than β2) means that vaccinated individuals are more likely to be infected by strain 2 than those who are not vaccinated, but if ζ2 is large it means that the population is taking precautions to avoid the infection of strain 2.
Finally, we remark that our model can be improved and generalized. For example, the model can be modified to contain two vaccinations, also we can consider the effect of time delay on vaccine-induced immunity and incorporate the diffusion of individuals. We leave these problems for future investigation.