key: cord-0115991-d706gp56
authors: Andersson, J.; Ghersheen, S.; Kozlov, V.; Tkachev, V.; Wennergren, U.
title: Effect of density dependence on coinfection dynamics
date: 2020-08-23
journal: nan
DOI: nan
sha: 5bc11d5435c9745d85b0982dfa82e96638f44b26
doc_id: 115991
cord_uid: d706gp56

In this paper we develop an SIR model for coinfection. We discuss how the underlying dynamics depends on the carrying capacity $K$: from a simple dynamics to a more complicated. This can help in understanding of appearance of more complicated dynamics, for example, chaos etc. The density dependent population growth is also considered. It is presented that pathogens can invade in population and their invasion depends on the carrying capacity $K$ which shows that the progression of disease in population depends on carrying capacity. Our approach is based on a bifurcation analysis which allows to generalize considerably the previous Lotka-Volterra type models.

Two or more pathogens circulating in the same population of hosts can interact in various ways. One disease can, for instance, impart cross-immunity to the other, meaning that an individual infected with the first disease becomes partially or fully immune to infection with the second [7, 18] . One disease can mediate the progression of other disease in population.

Therefore it is important to understand the dynamics of coexistent pathogens. In epidemiology the interaction of strains of the same pathogen, such as influenza or interacting diseases such as HIV/AIDS and hepatitis is very common and involves many complexities. The central problem in studying such systems is the explosive growth in the number of state variables of the system with the linear increase in the number of strains or pathogens [13] . Mostly these strains or pathogens are interacting in a way which has limited the further understanding and dynamics of such systems in terms of limited analytical progress. In this regard, it is a challenge to understand the dynamics and evolution of pathogens in population. The complexity of multiple strain models allows a great variability in modelling strategies. Recently, attention has focused on understanding the mechanisms that lead to coexistence, competitive exclusion and co-evolution of pathogen strains in infectious diseases which is important from the management of disease prospective.

Several studies exist on the coinfection with specific diseases. There are also studies [19, 16, 17, 7, 14] which have addressed this issue in general. In [6] , a mathematical model has been studied and showed that strains with differing degree of infectivity extinct, except for those that have higher the basic reproduction number. Allen et al in [1] showed coexistence only occur when the basic reproduction number is large enough for persistence of strains. They numerically illustrate the existence of globally stable coexistence equilibrium point. In an other study, Allen et al [2] , studied an SI model of coinfection with application on hanta virus. They assumed a logistic growth with carrying capacity and horizontal transmission of both viruses and yet only vertical transmission of virus 2. The condition of coexistence of two strain is described.

In [4] , a SIR model with vertical and horizontal transmission and a different population dynamics with limited immunity is considered. It is shown that the competitive exclusion can occur which is independent of basic reproduction number but a threshold. The existence and stability of endemic equilibrium is also shown. Since coinfection involves many complexities, many studies are only restricted to numerical simulation to understand the dynamics.

Nevertheless, mathematical modelling is one of the effective tool to understand the dynamics of biological system. But the major challenge is to balance between the practicality and mathematical solvability of the model. The cost of realisticity in mathematical modelling is the diminution of mathematical machinery.

The way to deal with this challenge is to divide the model into different sub models. Difference between the models is due to different biological assumptions. There are two major advantages in that case. First is the understanding of the system completely under certain assumptions. It can help to apply it to some real-life situations, since the controlling strategies for a diseases sometimes moves the original system to more simple system. In those cases the complete information about such system is needed to deal with that type of unexpected situation from management prospective. The second is, by relaxing assumptions, one can understand the role of each new parameter and its effects on the dynamics of epidemic.

One of the important characteristics, to understand the coinfection dynamics, is transmission mechanism. In paper [12] we have developed an SIR model to understand the dynamics of coinfection. Limited transmission is considered and the competitive exclusion principle is observed. The transition dynamics is also observed when the equilibrium points exist in the form of branches for each set of parameters. The compete dynamics of the system for all set of parameters is described by using linear complementarity problem. It appeared that there always exist an equilibrium point which is globally stable. It is showed that the dynamics of the system changes when carrying capacity changes. There are certain assumptions on the transmission of coinfection in that model. It is assumed that the coinfection can only occur as a result of contact between coinfected class and susceptible class, coinfected class and single infected classes. Interaction between two single infected classes is not considered. Also the simultaneous transmission of two pathogens from coinfected individual to susceptible individual is assumed.

In this paper we develop an SIR model for coinfection which is a relevant extension of model presented in [12] to understand the role of each new transmission parameter in the dynamics. Our aim here is to investigate how the dynamics changes due to a certain parameter, which in our case is the carrying capacity K, from a simple dynamics to a more complicated. This can help in understanding of appearance of more complicated dynamics for example chaos etc. Contrary to [12] , we could no more make use of the linear complementarity problem due to some additional term which appeared by relaxing the assumption of interaction between two single infected classes. We use the approach based on bifurcation analysis. The density dependent population growth is also considered. It is presented that pathogens can invade in population and their invasion depends on the carrying capacity K which shows that the progression of disease in population depends on carrying capacity.

The present model is displayed in Figure 1 . More precisely, we assume that the single infection cannot be transmitted by the contact with a coinfected person. According to Figure 1 , this process gives rise to the system of ODEs:

where we use the following notation:

• S represents the susceptible class,

• I 1 and I 2 are the infected classes from strain 1 and strain 2 respectively, • I 12 represents the co-infected class,

• R represents the recovered class.

Following [2, 6, 20] , we assume a limited population growth by making the per capita reproduction rate depend on the density of population. The recovery of each infected class is presented by the last equation in (1) . The fundamental parameters of the system are:

• r = b − d 0 is the intrinsic rate of natural increase, where b is the birthrate and d 0 is the death rate of S-class, • K is the carrying capacity (see also the next section), • ρ i is the recovery rate from each infected class (i = 1, 2, 3), where d 3 and d 4 correspond I 12 and R respectively, • d i is the death rate of each class, (i = 1, 2, 3, 4),

• α 1 , α 2 , α 3 are the rates of transmission of strain 1, strain 2 and both strains (in the case of coinfection), • γ i is the rate at which infected with one strain get infected with the other strain and move to a coinfected class (i = 1, 2),

• η i is the rate at which infected from one strain getting infection from a co-infected class (i = 1, 2);

Summing up all equations in (4) we have

where N = S + I 1 + I 2 + I 12 + R is the total population. We only need to consider the first four equations of (4) since R appears only in the last equation, hence it does not affect the disease dynamics. Rewrite the reduced system as

Furthermore, we only consider the case when the reproduction rate of the susceptible class is not less than their death rate, i.e.

Indeed, it is easy to see that the population will go extinct otherwise. The reduced system is considered under the natural initial conditions

Then it easily follows that any integral curve with (4) is well-defined and staying in the non negative cone for all t ≥ 0. Note also that since the variable R is not present in the first four equations, without loss of generality, we may consider only the first four equations of system (1).

It is convenient to introduce the notation

We shall always assume that the strains 1 and 2 are different, i.e. σ 1 = σ 2 . Then by change of the indices (if needed) we may assume that

Under this assumption, I 1 is the primary disease. Furthermore, let us assume for a moment that the susceptible class and only one infected class are nonzero. Let us suppose that only I i is nonzero zero. Then (3) reduces to

It is easy to see that there always exist two equilibrium points: the trivial equilibrium E 1 = (0, 0) and the disease-free equilibrium E 2 = (K, 0). Furthermore, if K > σ i then also exists (in the positive cone) the coexistence equilibrium E 3 = (σ i , r αi (1 − σi K )). Next, an elementary analysis shows that the following is true.

The trivial equilibrium state E 1 is always unstable. For any positive K = σ 1 there exists a unique locally stable equilibrium point E(K):

The reproduction number

can be used as a threshold. In other words, the transition from the disease-free equilibrium state to the disease equilibrium (the coexistence equilibrium point) occurs exactly when the reproduction number R 0 (I i ) of the corresponding infected class I i exceeds 1. We illustrate the transition by the diagram

The latter also clarifies the meaning of the parameter σ i . Namely, note that a more aggressive virus I has a greater value of R 0 (I). For a fixed value of the carrying capacity K this implies that a more aggressive virus I has a smaller value of σ (which, for example, means smaller recovery rate ρ or greater rate of transmission α).

It is natural to assume that the reproduction number of coinfection must be less than that of virus 1 and 2 respectively [15] . This makes it natural to assume the following hypotheses:

2.3 Some important notation

We shall also assume that (3) satisfies the following non-degenerate condition

This condition has a natural biological explanation: the virus strains 1 and 2 have different (co)infections rates. Let us also define

We also have

and

This implies an inequality which will be useful in the further analysis:

We shall further make use of the following relations:

On the other hand, one has

Remark 1 The parameters η * i can be thought of as the normalized co-infection rates. They play a distinguished role in the analysis of the thresholds given below.

The concise meaning of the parameter K becomes clear if we consider the limit case of (3) when the virus infection is absent, i.e. I 1 = I 2 = I 12 = 0. Then (1) reduces to the system

where the first equation (18) is the famous logistic (Verhulst) equation, r is the intrinsic rate of natural increase and K is the carrying capacity of the system. The carrying capacity K is one of the most fundamental parameters in population dynamics and it usually expresses the upper limit on the size of hypothetical populations, thereby enhancing mathematical stability. In basic ecology one defines carrying capacity as the equilibrium population size. Indeed, coming back to (3), we can see that K coincides with the healthy population size for the disease-free equilibrium. Mathematically this means that for any positive initial data, the corresponding solution of (18) converges to K as t → ∞. Furthermore, the equilibrium state G 2 := (K, 0, 0, 0) is the only possible equilibrium point of (3) with all I i = 0.

Equilibrium points of (3) are determined by the system

It is an elementary to see (see also Proposition 3 below for more explicit representations) that except for the trivial equilibrium point

and the disease-free equilibrium

there exist only 6 possible equilibrium points:

• three semi-trivial equilibria G 3 , G 4 , G 5 with only one nonzero infected class, i.e. I i = 0 for some i; • two coinfected semi-trivial equilibria G 6 , G 7 with I 12 = 0 but

Our main result extends the results obtained in [12] on the case of arbitrary values of γ i . More precisely, we will prove that we have the following possible scenarios for developing of an equilibrium point as a continuous function of increasing carrying capacity K:

Theorem 1 Let us assume that

Then there is exactly one locally stable nonnegative equilibrium point. Furthermore, changing the carrying capacity K, the type of this locally stable equilibrium point may be exactly one of the following alternative cases:

We consider the remained case

in the forthcoming paper [3] . This requires a delicate bifurcation analysis with application of methods similar to the principle of the exchange of stability developed in [8] ; see also [9] and [5] for recent applications in population analysis. We will show that in the remained cases one has the following two transition diagrams:

Furthermore, G 8 may loose stability for large K and small γ i in the latter case.

Remark 3 In particular, the above result implies that there are only three possible 'final destination' equilibrium states, namely G 3 , G 5 and G 8 . The latter may be interpreted as follows: taking the carrying capacity K as a random variable in (0, ∞), one concludes that a random choice of K corresponds to 'large enough' values of K, thus implying that, in average, the system (1) will have either of the three of the possible finishing scenarios G 3 , G 5 and G 8 . This is interesting from the epidemiological point of view: this means that if R 0 > 1, the most expected equilibrium states (in the aforementioned meaning) is either of the following: the state with the presence of only the first strain G 3 , purely coinfected case G 5 , or the state with the presence of all possible strains.

First we discuss some general results and equilibrium point analysis for (1).

In this section we discuss only stable equilibrium points with nonnegative coordinates. We denote Y = (S, I 1 , I 2 , I 12 ).

In what follows, by an equilibrium point we always mean an equilibrium Y of (3) point with nonnegative coordinates, Y = (S, I 1 , I 2 , I 12 ) ≥ 0.

In the next sections we identify all equilibria of the system (3) and determine their local stability properties. First, let us remark some useful relations which hold for any nonnegative equilibrium point of (3).

Lemma 1 Let Y = (S, I 1 , I 2 , I 12 ) = (0, 0, 0, 0) be a nontrivial equilibrium point of (3) with nonnegative coordinates. Then

and the right equality holds if and only if I 1 = I 2 = I 12 = 0, i.e. precisely when

Proof Let S = 0. Then we have from the second equation of (20) that (η 1 I 12 + γ 1 I 2 +µ 1 )I 1 = 0, where the nonnegativity assumption gives η 1 I 12 +γ 1 I 2 +µ 1 ≥ µ 1 > 0, hence I 1 = 0. For the same reason, I 2 = 0, thus the last equation in (20) yields µ 3 I 12 = 0, hence I 12 = 0 too. This proves that Y = (0, 0, 0, 0), hence implying the left inequality in (22). Now assume that Y = (S, I 1 , I 2 , I 12 ) = (0, 0, 0, 0) is an equilibrium point. Since S = 0, we have from the first equation of (20) that

In particular, the nonnegativity of the left hand side in the latter identity implies that K − S ≥ 0, i.e. proving the right inequality in (22). On the other hand, summing up equations in (20) we obtain

Assuming that S = K and dividing (25) by (24) we get S = µ 1 I 1 + µ 2 I 2 + µ 3 I 12 α 1 I 1 + α 2 I 2 + α 3 I 12 which readily yields (23).

Corollary 1 For any equilibrium point Y = (0, 0, 0, 0) and Y = G 2 there holds K ≥ σ 1 .

Notice that for G 2 , all I i = 0, otherwise we have Corollary 2 If an equilibrium point Y is distinct from G 2 := (K, 0, 0, 0) then (24) implies the following a priori bound on the I-coordinates:

where r = r is the intrinsic rate of natural increase. In other words, any equilibrium point distinct from G 2 lies inside a block with sides depending only on the fundamental constants.

The trivial equilibrium point G 1 = (0, 0, 0, 0) is the equilibrium of no disease or susceptible and the standard (local asymptotic) stability treatment shows that this point is always unstable. The first nontrivial equilibrium point G 2 is the disease-free equilibrium, i.e G 2 = (K, 0, 0, 0) and it always exist (for any admissible values of the fundamental parameters). The argument of [12] is also applicable in the present case because the stability analysis for G 2 does not involve γ i , so it is literally equivalent to that given in [12] . Repeating this argument (see section 8 in [12] ) readily yields the following criterium.

The following three conditions are equivalent:

(a) the disease-free equilibrium point G 2 is locally stable; (b) the disease-free equilibrium point G 2 is globally (asymptotically) stable; (c) 0 < K < σ 1 .

The latter proposition is completely consistent with the dichotomy of the R 0 -number (the reproduction number, sometimes called basic reproductive ratio). Recall that in epidemiology, the basic reproduction number of an infection can be thought of as the number of cases one case generates on average over the course of its infectious period, in an otherwise uninfected population. In our case, using the formal definition (see for example [10] ), one has

using the fact that the first strain is the most inclined to spread. In this notation, R 0 < 1 corresponds exactly to the scenario when the infection will die out in the long run (i.e. the only asymptotically stable equilibrium state is the disease-free equilibrium point G 2 ), while R 0 > 1 means the infection will be able to spread in a population. Therefore, in what follows, we shall focus on the nontrivial case R 0 > 1 with different scenario admitting the equilibrium states with some of I 1 , I 2 , I 12 nonzero.

Coming back to (20) , note that the Bezout theorem yields (in generic setting) that a quadratic system with four equations and four independent variables has 2 4 = 16 distinct solutions (counting the identically zero solution (0, 0, 0, 0)). In fact, in our case we have only one-half of the relevant (the Bezout number) solutions. More precisely, we have Proposition 3 Except for the trivial equilibrium G 1 = (0, 0, 0, 0) and the disease-free equilibrium G 2 = (K, 0, 0, 0) there exist only the following equilibrium states:

where

and there may exist at most two distinct points of type G 8 .

Proof Let Y = (S, I 1 , I 2 , I 12 ) = G 1 , G 2 be an equilibrium point. Then by Lemma 1 S > 0 and by the assumption some of coordinates I 1 , I 2 , I 12 must be distinct from zero. First assume that I 12 = 0. Then the last equation in (20) implies I 1 I 2 = 0. By the made assumption this implies that exactly one of I 1 and I 2 is nonzero while another vanishes. This yields G 3 and G 4 in (27) and (39), respectively. Now, let I 12 = 0 but I 1 I 2 = 0. Then the last equation in (20) implies α 3 S + η 1 I 1 + η 2 I 2 − µ 3 = 0. An elementary analysis reveals exactly three possible points G 5 , G 6 and G 7 in (29)-(31). Finally, consider the case when all coordinates of Y are distinct from zero. Since Y is distinct from G 1 and G 2 , it must satisfy (24), (25). Also, since I 1 , I 2 = 0, we obtain from the second and the third equations (20) the following system:

Rewriting these four equations in the matrix form as follows

we conclude that (I 1 , I 2 , I 12 , 1) T is a 0-eigenvector of the matrix in the left hand side of (33), thus, the first coordinate S satisfies the determinant equation

and

In particular, it follows that P (S) is a quadratic polynomial in S, therefore there may be at most two distinct inner points of type G 8 . The condition P (S) = 0 is sufficient if γ 1 , γ 2 < ∆α α3 .

It follows from Proposition 3 that all the boundary (edge) points are uniquely determined and can be expressed very explicitly. The existence and uniqueness of coexistence (inner) points of type G 8 is more involved (in contrast with the Lotka-Volterra caseγ = 0) and depends on the value ofγ.

We study the existence and the local stability of inner points by a bifurcation approach in the forthcoming paper [3] . Notice also that in the particular case γ i = 0, the characteristic polynomial (34) becomes a linear function expressed explicitly by

where we used the notation in (13) . This considerably simplifies the analysis, see [12] .

(i) For each G j , j = 1, 2, 3, 5, there exists ε > 0 (depending on the fundamental parameters α i , µ i , η i and γ i ) such that G j − G 8 ≥ ε. (ii) Let G 4 is given by (28) and δ := α 1 S * − γ 1 I * 2 − µ 1 > 0 (or equivalently γ * < K/(K − σ 2 )). Then there exists ε(δ) > 0 such that G 4 − G 8 ≥ ε(δ). (iii) Let G 6 is given by (30) and δ := α 2 S * − η 2 I 12 − γ 2 I * 1 − µ 2 = 0. Then there exists ε(δ) > 0 such that G 6 − G 8 ≥ ε(δ). (iv) Let G 7 is given by (32) and δ :=

Proof (i) We prove the assertion for j = 5 since other cases are considered in a similar way. The second and the third equations in (20) near the point G 5 have the form

where ǫ = G 5 − G 8 . By the assumption (7), one of numbers α 1 K − µ 1 , α 2 K − µ 2 does not vanish and so the corresponding coefficient in (35) does not vanish also for small ǫ, which implies (i) for G 5 . Proofs of (ii)-(iv) use the same argument.

It turns out that the most natural way to study equilibrium points is to consider their dependence on the carrying capacity K (or, equivalently, on the modified carrying capacity K). We know by Proposition 2 that the diseasefree equilibrium point G 2 is the only stable equilibrium point for 0 ≤ K < σ 1 .

In this section we consider each equilibrium state separately and study their local stability for K ≥ σ 1 . We study first the local stability of each point individually and in the next sections consider the dependence on K.

Our main goal is to describe all possible continuous scenarios of how the locally stable equilibrium states of (3) depends on K provided that all other fundamental parameters α i , µ i , b, γ i remain fixed. To this end, we introduce the following concept.

Definition 1 By an equilibrium branch we understand any continuous in K ≥ 0 family of equilibrium points of (3) which are locally stable for all but finitely many threshold values of K.

We need to distinguish the threshold values of K in the above definition because, formally, the local stability (i.e. that the real parts of all the system characteristic roots are negative) fails when an equilibrium point change its type. On the other hand, a branch may be stable in the Lyapunoff sense even for the threshold values of K. Indeed, the latter holds at least for γ = 0, see [12] .

4 The equilibrium state G 3 : Proof of (i) Note that the next three boundary equilibriums G 3 , G 4 and G 5 have a priori non-zero coordinates; furthermore, the S-coordinate is a constant (independent on K). The first of these is the equilibrium point G 3 with the presence of only the first strain. Its explicit expression with the nonnegativity condition are given by (27) . Remark that when K = σ 1 , the globally stable point G 2 bifurcates into G 3 :

Using (27), we find the corresponding Jacobian matrix evaluated at G 3 = (σ 1 , I * 1 , 0, 0):

Notice that, J 3 has a block structure. The left upper 2 × 2-block is obviously stable. Therefore J 3 is stable if and only if the right lower block is so. By virtue

or, equivalently, using the expression I * 1 = r Kα1 (K − σ 1 ) and (9) we obtain

After some obvious manipulations we arrive at

Notice that the point G 3 remains nonnegative and locally stable for any K > σ 1 provided η * 1 ≤ 1. This provides us with the first (simplest) example of a branch. More precisely, we have Corollary 3 (Branch (i)) Let η * 1 ≤ 1. Then (a) for 0 < K < σ 1 the point G 2 is locally (in fact, globally) stable; (b) for K = σ 1 the point G 2 coincides with G 3 ; (c) for K > σ 1 the point G 3 is locally stable.

We display this schematically as

The latter corollary implies (i) in Theorem 1.

Corollary 3 completely describes all possible scenarios for 0 ≤ K < ∞ when η * 1 ≤ 1. In what follows, we shall always assume that η * 1 > 1. Then Proposition 4 tells us that G 3 remains locally stable for any σ 1 < K < K 1 . If we want to find a continuous equilibrium branch, we need to check which of the remained candidates G 4 , G 5 , G 6 , G 7 , G 8 becomes equal to G 3 for the right critical value K = K 1 .

An easy inspection shows that for a generic choice of the fundamental parameters there is only one possible candidate, namely G 6 . Thus, to construct the only possible scenario for a continuous equilibrium branch is when G 3 bifurcates into G 6 . In the next section we give stability analysis of G 4 and G 5 , and then continue with G 6 and construction of equilibrium branches.

The equilibrium point G 4 expresses the presence of only the second strain, see (28). It is nonnegative if and only if

Note that if G 4 is nonnegative then by virtue of (39) and (7), G 3 is nonnegative too. The Jacobian matrix evaluated at G 4 is

Note that, interchanging rows and columns of the matrix (40) only change the sign of the determinant of this matrix. Therefore, after an obvious rearrangement, the eigenvalues of J 4 solves the following equation:

(41) Again, one easily verifies that the left upper 2 × 2-block is stable, while the stability of the right down (lower-diagonal) block is equivalent to the negativity of the diagonal elements, i.e. to the inequalities

Thus the stability of G 4 is equivalent to the inequalities

where γ * := γ1 A3 . In summary, we have Proposition 5 The equilibrium point G 4 is stable and nonnegative iff • K 3 < K < σ2γ * γ * −1 when γ * > 1 and η * 2 > 1, or • K > σ2γ * γ * −1 when γ * > 1 and η * 2 < 1.

In this paper, we are primarily interested in the case of 'small' values of γ i . On the other hand, the latter proposition shows that G 4 may be stable only if γ 1 > A 3 , therefore this equilibrium is not stable for small values of γ 1 and will be eliminated from the subsequent analysis.

The equilibrium point G 4 is locally unstable if 0 ≤ γ * 1 < 1.

An equilibrium point in the presence of coinfection is given by (29).

The equilibrium point G 5 is stable and nonnegative iff

Furthermore, if the point G 5 is nonnegative and locally stable for a certain K 0 > 0 then it will be so for any K ≥ K 0 (provided that other parameters are fixed).

Proof By (29), I * 12 = r Kα3 (K − σ 3 ), hence the positivity of I * 12 is equivalent to

Next, the Jacobian matrix evaluated at G 5 is

The matrix has a block structure where the block −r σ3 K −α 3 σ 3 α 3 I * 12 0 is obviously stable, therefore the stability of J 5 is equivalent to the negativity of two diagonal elements:

First notice that stability of G 5 implies immediately that I * 12 > 0. Also, taking into account that I * 12 = r Kα3 (K − σ 3 ), the stability of G 5 is equivalent to the inequalities

In summary, we have (29). Finally, the last statement of the proposition follows immediately from the increasing (with respect to K) character of the second inequality in (43).

We emphasize that the stability of the equilibrium states G 2 , G 3 , G 4 and G 5 does not involve the interference constants γ 1 , γ 2 .

Analysis of the remaining three equilibrium points G 6 , G 7 and G 8 is more delicate and now also involves the coinfection constants γ 1 , γ 2 . Let us consider the boundary equilibrium point

see (30). First notice that the coordinates of G 6 are nonnegative if and only if the two conditions hold: K 1 > 0, what is equivalent to η * 1 > 1, and also σ 1 < S * < σ 3 .

We see that G 6 is nonnegative if and only if

(Note that the bilateral inequality is inconsistent with (43)). Now let us study the local stability of G 6 . Using (30), the Jacobian matrix for G 6 is found as

with S * , I * 1 , I * 12 given by (30). Using the block structure of J 6 , it follows that G 6 is locally stable if and only if

• and the matrix below is stable:

The stability ofJ is equivalent to the stability of the last matrix factor in (47). An easy application of the Routh-Hurwitz criteria [11] confirms thatJ is always stable. Hence, the stability of G 6 is equivalent to the condition (46). Using (30), we can rewrite it as follows:

see (8) . Let us defineŜ

We have by using (9)-(10)

Consider first the case ∆ α + α 3 γ 2 = 0. Then by (15) it follows that ∆ µ + γ 2 µ 3 > 0 therefore (48) holds automatically true in this case, and G 6 is locally stable.

Next consider the case ∆ α + α 3 γ 2 < 0. Then it follows from (48) that G 6 is stable whenever S * >Ŝ 1 . On the other hand, (50) implies in this caseŜ 1 < σ 1 , therefore using (23) we see that

whenever S * is nonnegative. Therefore in this case G 6 is locally stable whenever (45) are fulfilled. Note also that under the made assumption ∆ α +α 3 γ 2 < 0 one necessarily has η 2 A 1 > η 1 A 2 . Indeed, if η 2 A 1 ≤ η 1 A 2 then (16) implies ∆ α > 0, therefore ∆ α + α 3 γ 2 > 0, a contradiction. Finally, assume that

Then by (48) the point G 6 is locally stable if and only if S * <Ŝ 1 , i.e.

Under assumption (54), (50) impliesŜ 1 > σ 1 . On the other hand, we havê

On the other hand, in the latter case, the inequality η * 2 ≥ η * 1 by virtue of (16) that in fact ∆ α > 0, therefore (54) holds automatically true in this case. Combining (55) with the nonnegativity condition (45), and summarizing the above observations we arrive at Proposition 7 The equilibrium point G 6 is nonnegative stable iff η 1 > A 1 and the following conditions hold:

Now we are ready to describe the equilibrium branch for η * 1 > 1.

Corollary 5 Let η * 1 > 1. Then (a) for 0 < K < σ 1 the point G 2 is locally (in fact, globally) stable; (b) for K = σ 1 the point G 2 coincides with G 3 ; (c) for σ 1 < K < K 1 the point G 3 is locally stable; (d) for K = K 1 the point G 3 coincides with G 6 ; (e) for K 1 < K < Q σ1 K 1 the point G 6 is locally stable, where Q is defined by (57).

We display this schematically as

Proof The first three items are obtained by combining Proposition 4 with Proposition 2. Note that the upper bound in (c) here is smaller than that in (c) in Corollary 3. When K = K 1 = σ1η * 1 η * 1 −1 , it follows that the I 12 -coordinate of G 6 vanishes (see (30)), i.e. G 6 = G 3 , which proves (d). Next, Proposition 7 yields (e).

With Corollary 3 and Corollary 5 in hand, it is natural to ask: What happens with an equilibrium branch when η * 1 > 1 and K > K 1 ? So far, we see that any continuous equilibrium branch develops uniquely determined accordingly (58). But at G 6 the situation becomes more complicated: this point may a priori bifurcate into different points.

In this paper we only consider the particular case (ii), i.e. when 1 < η * 1 < η * 2 . This yields by (57) that Q = σ 3 , hence (56) implies that G 6 is locally stable for

The upper critical value K 2 substituted in (30) implies that I * 1 = 0, hence G 6 naturally bifurcates into G 5 . It is easy to see that the corresponding I *

for G 5 and G 6 coincide when K = K 2 holds. This observation combined with Proposition 6 implies that in this case for any K > K 2 the point G 5 will be locally stable, hence we arrive at Corollary 6 (Branch (ii)) Let η * 2 ≥ η * 1 > 1 hold. Then (a) for 0 < K < σ 1 the point G 2 is locally (in fact, globally) stable; (b) for K = σ 1 the point G 2 coincides with G 3 ; (c) for σ 1 < K < K 1 the point G 3 is locally stable; (d) for K = K 1 the point G 3 coincides with G 6 ; (e) for K 1 < K < K 2 the point G 6 is locally stable; (f ) for K = K 2 the point G 6 coincides with G 5 ; (g) for K > K 2 the point G 5 is locally stable.

We display this schematically as

Thus, we are remained to study the case when

hold. Notice that in fact by virtue of (16) the latter inequality implies

We know by (e) in Corollary 5 that G 6 is locally stable for

Substituting the corresponding critical value K = K 0 such that

η * 1 − 1 in (30) reveals that the coordinates G 6 do not vanish, i.e. G 6 does not change its type. Instead it losts its local stability because the determinant of J 6 vanishes at this moment. To continue the equilibrium branch (58) beyond G 6 we need to find an appropriate candidate for a stable point. By the continuity argument (because G 6 keeps all coordinates nonzero for K = K 0 ), the only possible candidate for a continuous equilibrium branch is a point of type G 8 . Since we do not have any explicit expression of G 8 , the analysis in this case is more complicated and involves a certain bifurcation technique which we develop in a forthcoming paper [3] .

It is natural, from biological point of view, to relax the constancy condition on the transmission rates α i and assume that in general they may depend on the carrying capacity. Indeed, a larger carrying capacity corresponds to a larger size of the population (or the susceptible class), thus it must imply a slower spread of strains, i.e. α i = α i (K) must be a non-increasing function. One natural assumptions is the following relation:

This implies for the other fundamental constants

The main consequence of (62) is that the coordinates of a stable equilibrium point is no longer bounded and develop as K increases. For example, under assumption (62) one has from (23) merely

This, in particular implies that already the first bifurcation S 2 → S 3 is completely different. Indeed, it follows from Proposition 2 that G 2 becomes stable for all K > 0 provided s 1 ≥ 1. In the nontrivial case s 1 < 1, G 2 is never stable. In general, Proposition 4 and Corollary 6 instead imply Corollary 7 We have the following stability analysis:

(i) If s 1 ≥ 1 then G 2 is stable for all K > 0;

(ii) If s 1 < 1 and 0 < η * 1 ≤ 1 1−s1 then G 3 stable for all K > 0; Let now s 1 < 1, η * 2 > η * 1 > 1 1−s1 hold. Then (iii) if s 3 ≥ 1 or s 3 < 1 and η * 1 < 1 1−s3 then G 6 stable for all K > 0; (iv) if s 3 < 1 and η * 1 > 1 1−s3 then G 5 stable for all K > 0.

Thus, we have a complete description in the cases η * 1 ≤ 1 and η * 2 ≥ η * 1 > 1. The remained case η * 1 ≥ max{1, η * 2 } will be considered in [3] .

Competitive exclusion and coexistence for pathogens in an epidemic model with variable population size

The dynamics of two viral infections in a single host population with applications to hantavirus

Effect of density dependence on coinfection dynamics, the bifurcation analysis

Global analysis of multi-strains sis, sir and msir epidemic models

Introducing a population into a steady community: the critical case, the center manifold, and the direction of bifurcation

A competitive exclusion principle for pathogen virulence

On the relationship between evolution of virulence and host demography

The principle of exchange of stability

Stability and bifurcation analysis of Volterra functional equations in the light of suns and stars

On the definition and the computation of the basic reproduction ratio R 0 in models for infectious diseases in heterogeneous populations

The theory of matrices. Vols. 1, 2. Translated by

Dynamical behaviour of sir model with coinfection: the case of finite carrying capacity

Dynamics and selection of many-strain pathogens

Effects of metapopulation mobility and climate change in si-sir model for malaria disease

The role of coinfection in multidisease dynamics

Coinfection and the evolution of parasite virulence

Evolution of virulence: a unified framework for coinfection and superinfection

Threshold effects for two pathogens spreading on a network

Superinfection and the evolution of parasite virulence

Population size dependent incidence in models for diseases without immunity