# Analysis of a Mathematical Model of Emerging Infectious Disease Leading to Amphibian Decline

Authors: Muhammad Dur-e-Ahmad; Mudassar Imran; Adnan Khan
Journal: Abstract and Applied Analysis
 (2014)

---

## Abstract

We formulate a three-dimensional deterministic model of amphibian larvae population to investigate the cause of extinction due to the infectious disease. The larvae population of the model is subdivided into two classes, exposed and unexposed, depending on their vulnerability to disease. Reproduction ratioℛ0 has been calculated and we have shown that if ℛ0<1, the whole population will be extinct. For the case of ℛ0>1, we discussed different scenarios under which an infected population can survive or be eliminated using stability and persistence analysis. Finally, we also used Hopf bifurcation analysis to study the stability of periodic solutions.


---

## Body

## 1. Introduction

Worldwide catastrophic declines in the amphibian population are perhaps one of the most pressing and discussed problems among the ecologists during the last two decades. Although many of these are attributable due to the habitat loss (see, e.g., [1, 2]), the majority have remained enigmatic till today. Many hypotheses responsible for this decline have been documented in the literature, such as adverse weather patterns [3, 4], acid precipitation [5], environmental pollution [3], increased ultraviolet (UV-B) radiation [6], introduction of predators or competitors [7], infectious disease [3, 8, 9], or a combination of these.Recently, infectious diseases have become one of the emergent factors behind rapid amphibian decline which often results in extinction of species, for example, the recent extinction of the golden toad in Costa Rica and some species of gastric-brooding frogs in Australia. Investigations reveal that the main cause behind the mass death of these species is two infectious diseases, chytridiomycosis in the rain forest of Australia and Central America and some parts of North America and iridoviral infections in United Kingdom, United States, and Canada [3, 10, 11].The larval stage of the amphibian is considered to be the most vulnerable towards the spread of infectious disease. Most of the larva population of tropical amphibian species remain alive for 12 to 18 months but some temperate species may also survive as long as 3 years before metamorphosing. Recently it has been observed, both in Australia and Central America, that larval amphibians infected with chytridiomycosis may exhibit disfigurement of their keratinized mouth parts and demonstrate a significant reduction in their growth and development which may eventually cause the complete extinction sometimes [3]. Further, in the case of reduced amphibian population, this infection causes prolonging of the existence ofBatrachochytrium [12] and implicates the life-cycle stage as a reservoir host for the pathogen. This kind of larval infection enhancing pathogen-mediated host population extinction has also been reported for invertebrates [13]. Therefore in this work, we keep our focus on the larval stage and investigate the dynamics of spread of disease among them.This work has been motivated by the recent work of [9], in which only susceptible and infected classes are considered. Since, as discussed earlier, larvae are the key players in the case of density-dependent disease incidence which eventually leads to the host extinction, therefore, it is imperative to keep this component intact while studying the disease based amphibian decline. Since the main focus of the model is the disease based larvae extinction, therefore, in order to look at the complete scenario, we assume that, at the early stage of their lives, the larvae are surrounded by a safe environment which is free of disease, like in a closed shell or some protected holes. For our model, we denote this class by L. At the later stage once these unexposed larvae enter into free environment they become vulnerable to disease and we call them susceptibles (denoted by S). These susceptibles can be infected with the disease and enter into the infective class, denoted by I. Further, these infected larvae can also recover from the disease and return to the susceptible class. The model is formulated in a way that the number of new infected larvae depends on the present number of susceptible and infective larvae in the sense that the more larvae (susceptible or infective) we have, the more chance of spreading the infections among the healthy or susceptible larvae will there be; that is, the disease can be transferred from one to another.We organize this paper as follows. First we describe the mathematical model along with its underlying assumptions. Then we introduce many threshold quantities, for example, the reproduction ratios and the critical host density for the disease establishments. We also derive the number of endemic equilibriums and discuss their stabilities. Next we also elaborate our findings using numerical simulations. Finally, we give the conclusion and discuss the potential limitations and impact of our results.


## 2. Formulation of Mathematical Model

In this mathematical model we assume that the disease only affects the larval population. So we divide the larvaeL into two categories, the susceptibles larvae S and the infected larvae I. Although this disease is transferable from one to another, the recovery is also possible and on recovering from the disease, an infected larva will reenter into the class of susceptible larvae. Further, for our model we also assume that only susceptibles can contribute to the reproduction process. Thus once the disease is spread, the whole population will go towards extinctions not only due to the illness but also due to the lack of reproduction process.The model is as follows:(1)L′=βS-π(L)L-μL,S′=π(L)L-νS+ρI-σSI,I′=σSI-νI-αI-ρI.
Here L represents the unexposed class of the larval population which is not yet entered into the susceptible class S. Here we assume that the transition from stage L to stage S requires a certain minimum size. This assumption is valid for many species, for example,Daphnia magna, a water flea, where the individuals typically have a length of 0.8 mm at birth and a length of 2.5 mm once they enter into the exposed class [14]. Further, in the case of some amphibians, the body size at metamorphosis is quite flexible (see, e.g., [15]); therefore, a certain minimum size still seems to be required in these cases. Thus based on these evidences, it is reasonable to assume that the scarcity of resources will prolong the length of the larval stage L. So, if there are more larvae in stage L and less resources, then it will take longer for a single larva in stage L to complete the transition into stage S. Therefore, in our model, we assume that π(L)L is nondecreasing in L in order to represent a strong negative feedback from the number of larval population L. It is clear from the model that this class will enter into the exposed susceptible class S with the rate μ which further can become infected with a rate ν. Since, in the case of closed population, the disease can transfer from one individual to the other, therefore by using the law of mass action, we express this by σSI, where σ is disease transmission rate. Finally the last term in the third equation represents the recovery from disease with the rate ρ.The detailed description of parameters of our model is given in Table1.Table 1
Description of parameters in model (1).
Variables
Description
β
Per capita birth rate of susceptibles.
μ
Per capita (natural) mortality rate of the larvae.
ν
Per capita (natural) mortality rate of the susceptible and infected larvae.
α
Per capita mortality rate of the infected larvae due to disease.
ρ
Recovery rate from disease.
π
(
L
)
Per capita transition rate. This is the rate at which the susceptible larva matures.
σ
Rate at which the disease is transferred from infected to susceptible class, when the contact between susceptible and infected individuals occurs.Assumption 1.
All the parameters are positive exceptα, which can also be zero. π(L) is a strictly decreasing nonnegative function of L such that π(L)→0 as L→∞.


## 3. Existence, Positivity, and Boundedness of the Solution

Theorem 2.
Assume that Assumption1 is satisfied; then for nonnegative initial data, there are unique nonnegative solutions of the system which are defined for all nonnegative times. Further all the solutions to the system are uniformly eventually bounded.Proof.
LetF=(f1,f2,f3) and X=(x1,x2,x3)=(L,S,I), where
(2)f1=βS-π(L)L-μL,f2=π(L)L-νS+ρI-σSI,f3=σSI-νI-αI-ρI.
Notice that all the partial derivatives ∂fi/∂xj,i,j=1,2,3, are continuous. Further
(3)L=0⟹f1=βS≥0,S=0⟹f2=π(L)L+ρI≥0,I=0⟹f3=0.
So X(t)=(L(t),S(t),I(t))∈B3, where B=[0,∞) for all t≥t0≥0 for which it is defined and whenever X(t0)∈B3. Thus all the solutions of the given system are nonnegative by proposition A.1 of [14]. We conclude that there exists a unique solution to the above system with the values in ℝ+3 and it is defined in the interval [0,b), b∈(0,∞] (Theorem A.4 [14]) and if b<∞ then
(4)sup⁡0≤t<b{L(t)+S(t)+I(t)}=∞.
To prove the first part of the theorem, we need to show only that b=∞. Assume that b<∞; then from model (1) we get
(5)L′≤βS,S′≤π(L)L+ρI-σSI≤π(0)L+ρI-σSI,I′≤σSI.
Let η=max⁡{π(0),β,ρ}; then we get
(6)L′+S′+I′≤η(L+S+I)⟹L(t)+S(t)+I(t)≤[L(0)+S(0)+I(0)]eηt⟹sup⁡0≤t<bL(t)+S(t)+I(t)≤sup⁡0≤t<b[L(0)+S(0)+I(0)]eηt.
Since b<∞, therefore
(7)sup⁡0≤t<bL(t)+S(t)+I(t)≤C<∞.
This is a contradiction, so b=∞.
Now we show that the solutions are uniformly eventually bounded. Let us defineV=L+ξ(S+I), where ξ is to be determined later. Then by using the system (1), we have
(8)V′=βS-π(L)L-μL+ξ(π(L)L-νS-νI-αI)⟹V′=(β-ξν)S+[(ξ-1)π(L)-μ]L-ξ(ν+α)I.
By choosing ξ=(β/ν)+1, we get
(9)V′=-νS+(βνπ(L)-μ)L-ξ(ν+α)I.
Since we know that π(L)→0 as L→∞, therefore for each ϵ>0 there exist some L♯>0 such that, for all L≥L♯, π(L)<ϵ. So we can divide the proof into two cases.
Case
1  (L≥L♯). From (9) we get
(10)V′≤-νS+(βνϵ-μ)L≤-νS-μ2L≤-μ2L♯.
So if we take L♯>2/μ, then we get V′<-1.
Case
2  (L<L♯). Using ξ=(β/ν)+1, again from (9), we have
(11)V′≤-νS+βνπ(L)L-μL-ξ(ν+α)I≤-νS+βνπ(0)L♯-ξ(ν+α)I≤-ν(1ξ(V-L)-I)+βνπ(0)L♯-ξ(ν+α)I≤-νξV+νξL+νI+βνπ(0)L♯-(βν+1)(ν+α)I.
Simplification and neglecting the negative terms (involving I), yields
(12)⟹V′≤-νξV+(νξ+βνπ(0))L♯.
We find some V♯ such that V′<-1 whenever V≥V♯.
Therefore in either caseV′<-1. This implies that liminf⁡t→∞V(t)≤V♯. Otherwise since V′(t)≤-1 for all sufficiently large time and V would become negative in finite time, contradiction. Further since V≥L+S+I and we have shown that our system is weakly dissipative, so it is dissipative by Proposition 3.18 of [14].


## 4. Existence and Stability of Steady States

In this section we will discuss the existence and stability of the steady states of our model (1). We will have three steady states, the trivial steady states (0,0,0), the disease-free steady state of the form (L^,S^,0), and the interior steady state of the form (L*,S*,I*). The trivial steady state always exists. To find the disease-free steady state, we proceed as follows. From the I′ equation of our main model (1) we have
(13)(σS-ν-α-ρ)I=0.
From (13) we can say that either I=0 or S=(ν+α+ρ)/σ. Let us suppose the disease-free case(L^,S^,0) that is, I^=0; then by first and second equations of model (1) we have
(14)βS^-π(L^)L^-μL^=0π(L^)L^-νS^=0.
Now by adding these two equations we get
(15)S^=μ(β-ν)L^.
To find the value of L^, we substitute the value of S^ from (15) and I^=0 in L′ equation of model (1). We have
(16)βμβ-νL^-π(L^)L^-μL^=0[βμβ-ν-π(L^)-μ]L^=0⟹L^=0orβμβ-ν-μ=π(L^).
Notice that L^=0 will give us trivial state (0,0,0). In order to find the nontrivial state, we will consider the case when L^≠0. In this case we will have
(17)βμβ-ν-μ=π(L^)π(L^)=μνβ-ν.
Equivalently
(18)L^=π-1(μνβ-ν).
Thus we have two disease-free steady states, that is, (0,0,0) and (L^,(μ/(β-ν))L^,0), where L^>0 is given by (17). Observe that there is a threshold condition
(19)βμβ-ν-μ<π(0).
We can rewrite the condition (β-ν)π(0)<νμ as (β/ν)(π(0)/(π(0)+μ))<1. We refer this quantity as reproduction number and it is defined as
(20)ℛ0=βνπ(0)π(0)+μ<1.
Now we find the interior equilibrium (L*,S*,I*) with L*,S*,I*≠0. From the I′ equation of model (1), we get S*=(ν+α+ρ)/σ.Adding the first two equations of our original model (1) gives
(21)βS*-μL*-νS*+ρI*-σS*I*=0⟹(β-ν)S*+(ρ-σS*)I*-μL*=0.
Now by using the values of S*, we get
(22)(β-ν)(ν+α+ρ)σ+(ρ-ν-α-ρ)I-μL=0⟹(β-ν)(ν+α+ρ)+σ(-ν-α)I-σμL=0.
Therefore we have
(23)I*=σμL*-(β-ν)(ν+α+ρ)σ(-ν-α).
Now we need to find the value of L*. By substituting the values of S* and I* in the first equation of model (1), we will have
(24)βσ(ν+α+ρ)-(π(L*)+μ)L*=0(π(L*)+μ)L*=βσ(ν+α+ρ)=βS*.Now we need to solve (24) numerically (depending on the structure of the function π(L)) to find the value of L*. Therefore, in this case, the interior equilibrium is given by
(25)(L*,ν+α+ρσ,σμL*-(β-ν)(ν+α+ρ)σ(-ν-α)).
The L* and S* are nonzero. From the S′ equation of model (1) we get
(26)π(L*)L*-νS*+I*(ρ-σS*)=0.
If we solve this for I*, we get
(27)I*=π(L*)L*-νS*ν+α.
This implies that I*>0 if
(28)π(L*)L*>νS*.
This gives the condition for the existence of infected steady state. Observe that we can write (24) as
(29)σS*ν+α+ρ=1.
Since 1/(ν+α+ρ) is the average time a larva spends in the infected stage and σS* is the rate at which one average infected larva infects a susceptible larva, if they are at density S*, therefore we can express the left side of (29) as the replacement ratio of the infectious disease from the infected larvae to the susceptible larvae when they have density S*. That is,
(30)ℛpar(S*)=σS*ν+α+ρ.
The replacement ratio at the susceptible larval density is given by
(31)ℛ0par(S^)=σS^ν+α+ρ.From the above results we have the following theorem for the existence of the steady states.Theorem 3.
(a) The trivial steady state(0,0,0) exists all the time.
(b) The disease-free state(L^,S^,0) exists if and only if ℛ0>1.
(c) The infected state(L*,S*,I*) exists if and only if ℛ0>1 and π(L*)L*>νS*.
(d) Letℛpar(L^)>1. Then there is at least one steady state with the infected larvae.
(e) Letℛpar(L^)<1; then there is no or two steady states with the infective larvae.Proof.
Proof of part (b) and part (c) follows directly from the discussion above this theorem.
(d) The conditionℛpar(L^)>1 can be written as
(32)(π(L^)L^+μL^)>βν+α+ρσ.
The intermediate value theorem implies that there is some L~1∈(0,L^) which satisfies (24).
(e) The conditionℛpar(L^)<1 can be written as
(33)(π(L^)L+^μL^)<βν+α+ρσ.
If π(L)L+μL is strictly increasing function for L≥0, then there is no steady state in [0,L^]. Now assume that this function is not increasing and that there is at least one solution Lˇ1>0 of (24). We can choose Lˇ1∈(0,L1]. If Lˇ1=L1 then we have two solutions counting multiplicities. They are steady state with the infectives if Lˇ1>L^; otherwise that steady state does not exist with infectives. If Lˇ1<L1, then we can choose some L-∈(Lˇ1,L1). This implies
(34)π(L-1)L-1+μL-1>π(Lˇ1)Lˇ1+μLˇ1=β(ν+α+ρσ)π(L-1)L-1+μL-1>β(ν+α+ρσ).
Therefore by the intermediate value theorem, there exists a solution Lˇ2∈(Lˇ1,L-) of (20). By the previous theorem, Lˇ1 and Lˇ2 are components of equilibria with positive infective population.Theorem 4.
LetL*>0 be a solution of (24); then L* satisfy (29) if and only if L*<L^ with L^=π-1(μν/(β-ν)).Proof.
LetL*<L^. Since π(L) is strictly decreasing function, this is equivalent to
(35)π(L*)>π(L^).
Since π(L^)=νμ/(β-ν), so we have the equivalent expression
(36)π(L*)>νμβ-ν.
Multiplying both sides with L* and rearranging the terms imply
(37)(β-ν)π(L*)L*>νμL*.
At infected steady state (L*,S*,I*) the L′ equation of model (1) can be written as
(38)μL*=βS*-π(L*)L*.
So we get
(39)(β-ν)π(L*)L*>ν(βS*-π(L*)L*)π(L*)L*>νS*.
Similarly if L*≥L^ then (31) will not be satisfied.It is also clear that (24) has one or three solutions depending on the choice of the parameters ν, α, ρ, σ, and β.Now we will discuss the stability of the trivial steady state(0,0,0) and the persistence of the population. The persistence of a population is defined as follows.Definition 5.
The populationL(t),S(t),I(t) is called uniformly strongly persistent if there exists ϵ>0, independent of initial data, such that for all solutions of model (1) satisfying L(0)>0, S(0)>0, and I(0)>0 we have L(t),S(t),I(t)>ϵ for all sufficiently large t. It is robust uniform persistence, if model (1) depends continuously on a parameter η∈ℛ. If η0 is a parameter and uniform persistence holds with ϕ(t;,S,L,I) replaced by ϕ(t;,S,L,I,η), where ϕ(t,x) is semiflow, for all η in the neighborhood of  η0, then we say that the system (1) is robust uniformly persistent.Theorem 6.
(a) The trivial state is locally asymptotically stable if and only ifℛ0<1.
(b) Assume thatℛ0>1 and 𝒜^=π′(L^)L^+π(L^)+μ>0 (or equivalently (β/ν)π(0)+ν>-π′(L^)L^). Then the steady state (L^,S^,0) is locally asymptotically stable if and only if  S^<(ν+α+ρ)/σ or S^<S*.
(c) Assume thatℛ0>1 and 𝒜=π′(L*)L*+π(L*)+μ>0 (or equivalently (β/ν)π(0)+ν>-π′(L*)L*). Then the endemic equilibrium is locally asymptotically stable.
(d) Ifℛ0>1, then the disease is robustly uniformly persistent: that is, for every parameter vector η0 for the system (1), there exists a neighborhood 𝒩 of η0 and ϵ>0 such that
(40)liminf⁡t→∞min⁡(Sη(t),Iη(t))>ϵ,η∈𝒩,
for all solutions (Lη(t),Sη(t),Iη(t)) of model (1) corresponding to η with (Sη(0)+Iη(0))>0.Proof.
The Jacobian matrix associated with our model (1) is given by
(41)J=[-π′(L)L-π(L)-μβ0π′(L)L+π(L)-ν-σIρ-σS0σIσS-ν-α-ρ].
(a) This Jacobian matrix at steady state(0,0,0) is
(42)J1=[-π(0)-μβ0π(0)-νρ00-ν-α-ρ].
The three eigenvalues areλ1=-(ν+α+ρ), and the eigenvalues of the matrix
(43)A1=[-π(0)-μβπ(0)-ν].
The trace and determinant are given by
(44)trace(A1)=-(ν+μ+π(0))<0det⁡(A1)=νπ(0)+νμ-βπ(0)=νμ-(β-ν)π(0).
So the steady state (0,0,0) is locally asymptotically stable if and only if (β-ν)π(0)<νμ.
(b) The Jacobian matrix at(L^,S^,0) is
(45)J=[-π′(L^)L^-π(L^)-μβ0π′(L^)L^+π(L^)-νρ-σS^00σS^-ν-α-ρ].
One of the eigenvalues is λ1=σS^-ν-α-ρ, which is negative if S^<(ν+α+ρ)/σ. The other two eigenvalues come from
(46)J1=[-π′(L^)L^-π(L^)-μβπ′(L^)L^+π(L^)-ν].
Note that
(47)𝒯=trace(J1)=-A^-ν𝒟=det⁡(J1)=(σS^-ν-α-ρ)[-A^(ν+β)+βμ],
where A^=π′(L^)L+π(L^)+μ.
Notice that𝒟 is positive if A^>0.
(c) The endemic equilibrium is given by(L*,(ν+α+ρ)/σ,(σμL*-(β-ν)(ν+α+ρ))/σ(-ν-α)), where L* is the solution of (π(L)+μ)L=(β/σ)(ν+α+ρ).
The Jacobian matrix at this steady state is given by(48)𝒥*=[-π′(L*)L*-π(L*)-μβ0π′(L*)L*+π(L*)-ν-σI*ρ-σS*0σI*σS*-ν-α-ρ].
Since σS*-ν-α-ρ=0 and ρ-σS*=-(ν+α), so we get
(49)𝒥*=[-𝒜β0𝒜-μ-ν-σI*-(ν+α)0σI*0],
where 𝒜=π′(L*)L*+π(L*)+μ.
Now we have(50)𝒯=trace(𝒥*)=-(𝒜+ν+σI*)<0𝒟=det⁡(𝒥*)=-σI*𝒜(ν+α)<0.
Further we have
(51)𝒥*11=σI*(ν+α)𝒥*22=0𝒥*33=𝒜(ν+σI*)-β(𝒜-μ),
where 𝒥kk is the determinant of the matrix of size 2 resulting from 𝒥* by removing the kth row and the kth column:
(52)⟹∑k=13‍𝒥kk=σI*(ν+α)+𝒜(ν+σI*)-β(𝒜-μ)𝒯∑k=13‍𝒥*kk=-(𝒜+σI*)×[σI*(ν+α)+𝒜(ν+σI*)-β(𝒜-μ)]=-𝒜σI*(ν+α)-𝒜2(ν+σI*)+𝒜β(𝒜-μ)-(ν+σI*)σI*(ν+α)-(ν+σI*)𝒜(ν+σI*)+(ν+σI*)β(𝒜-μ).
Finally
(53)𝒟-𝒯∑k=13‍𝒥*kk=𝒜2(ν+σI*)-𝒜β(𝒜-μ)+(ν+σI*)σI*(ν+α)+(ν+σI*)𝒜(ν+σI*)-(ν+σI*)β(𝒜-μ).
Simplifying and rearranging the terms we get
(54)=𝒜2(ν+σI*)+(ν+σI*)σI*(ν+α)+𝒜(ν+σI*)2-β(𝒜-μ)[𝒜+ν+σI*]≥-β(𝒜-μ)[𝒜+ν+σI*].
Since 0<𝒜=(π(L)L)L=L*′+μ, so this implies μ-𝒜=-(π(L)L)L=L*′>0 by the assumption of the structure of π(L)L (i.e., the case when we have three interior equilibria and this is the most right one). Thus the endemic equilibrium is locally asymptotically stable.
(d) Letη be a fixed parameter and X={x∈ℛ3S=I=0}. Define M=V∩X, where V is the same region defined in Theorem 2 in which all the solutions are uniformly eventually bounded. Notice that both X and M are positively invariant and E0=(0,0,0)∈M that attracts all the solutions in X. Considering E0 as a periodic orbit (say, e.g., of period T=1), then our results will follow by applying the theorems discussed in [16] (see, in particular, Corollary (4.7), Proposition (4.1), and Theorem (3.2) of [16]). Now if we denote the last two equations of model (1) in the form of a 2×2 matrix A(x), then it is clear that A(E0) has a positive spectral radius provided ℛ0>1. So condition (9) of Corollary 4.7 in [17] is satisfied. Also, A(E0) is irreducible which, using Theorem 1.1 in [18], implies that etA(E0) has all positive entries for all t>0. This implies that condition (1) of the abovementioned corollary is also satisfied.Theorem 7.
L
(
t
)
,
S
(
t
), and I(t) converge to zero as t→∞ if and only if ℛ0<1.Proof.
We will use “La Salle’s Invariance Principle” to prove the theorem.
Case
1  ((β-ν)<0).Consider the Lyapunov function V=L+S+I. Since L,S,I≥0, so V≥0 is a continuously differentiable function with V(0)=0. Now consider
(55)V*=L˙+S˙+I˙=βS-π(L)L-μL+π(L)L-νS+ρI-σIS+σSI-νI-αI-ρI=(β-ν)S-μL-(ν+α)I≤0.
Now define
(56)E={(L,S,I)∣(β-ν)S-μL-(ν+α)I=0}.
Further, as we have
(57)(β-ν)S-μL-(ν+α)I=0⟹S=μL+(ν+α)Iβ-ν.
Since L and I are nonnegative and also β-ν<0, so we have S≤0. But from the condition that S should be nonnegative we get S=0. So we have
(58)μL+(ν+α)I=0.
As μ and ν+α are positive, therefore we must have L=I=0.
So(0,0,0) is the only subset of E which is also one of the steady states and so it is the largest invariant subset of E. Since we have shown above that all the solutions of the above system are bounded, thus, by the above result, all the solutions to the given system converge to (0,0,0).
Case
2  ((β-ν)=0). Again consider V*=L+S+I. In this case we have
(59)V*=-μL-(ν+α)I≤0.
Define
(60)E={(L,S,I)∣-μL-(ν+α)I=0}.
Since we have
(61)-μL-(ν+α)I=0⟹L=-(ν+α)Iμ⟹L≤0.
Thus from the nonnegativity of L, we conclude that L=0 and so we get I=0 and E={(0,S,0),S≥0}since, at (0,S,0), we get from the first equation of our main model that S=0.
Thus(0,0,0) is the only invariant subset of E; therefore all the bounded solutions of the given system will converge to (0,0,0).
Case
3  ((β-ν)>0). Consider the Lyapunov function V=νL+β(S+I). Then
(62)V*=νL˙+β(S˙+I˙)=ν(βS-π(L)L-μ(L))+β(π(L)L-νS+ρI-σSI)+β(σSI-νI-αI-ρI)=-νπ(L)L-νμL+βπ(L)L-β(ν+α)I=(β-ν)π(L)L-νμL-β(ν+α)I≤(β-ν)π(L)L-νμL≤((β-ν)π(0)-νμ)L≤0=(β-ν)π(L)L-νμL-β(ν+α)I≤0.
Define
(63)E={(L,S,I)(β-ν)π(L)L-νμL-β(ν+α)I=0}.
Since
(64)(β-ν)π(L)L-νμL-β(ν+α)I=0,(β-ν)π(L)-νμ≤(β-ν)π(0)-νμ≤0,
so we haveI≤0. Therefore, from the nonnegativity of I, we conclude that I=0. So we will have either L=0 or (β-ν)π(L)-νμ=0.
Notice that ifL≠0, then L>0, so from the strictly decreasing property of π, π(L)<π(0), and therefore
(65)(β-ν)π(L)-νμ<(β-ν)π(0)-νμ≤0.
Thus (β-ν)π(L)-νμ<0 and the only possibility is L=0 and hence (0,S,0) is the only subset of E. Again from the first equation of our model (1), L′=βS, which implies L′=0 if and only if S=0. Hence (0,0,0) is the only invariant subset of E; therefore all the bounded solutions of the given system will converge to (0,0,0).


## 5. The Disease Transmission Model Special Case

In this section we will assume thatπ(L)=ξ(1-ψL) with ξ, ψ>0. Also assume that π(L)=(1-L) for L≤1. π(L) is strictly decreasing as long as it is strictly positive. Let us define L•=ψL, S•=ψS, and I•=ψI. This implies
(66)L=L•ψ,S=S•ψ,I=I•ψ.
Now by using these substitutions, our main model will become
(67)L•′ψ=βS•ψ-ξ(1-L•)L•ψ-μψL•S•′ψ=ξ(1-L•)L•ψ-νψS•+ρψI•-σψ2S•I•I•′ψ=σψ2S•I•-(ν+α+ρ)I•ψ.
On simplifying these equations, we get
(68)L•′=βS•-ξ(1-L•)L•-μL•S•′=ξ(1-L•)L•-νS•+ρI•-σψS•I•I•′=σψS•I•-(ν+α+ρ)I•.
First, we will make this system dimensionless. Assume that the β, ξ, μ, ν, ρ, and α are measured in 1/time. L,S, and I are measured in the units of population size, say p. So ψ is measured in 1/p. Also σ is measured in 1/((p)(time)). Now dividing by ξ both sides our system becomes
(69)1ξL•′=β~S•-(1-L•)L•-μ~L•1ξS•′=(1-L•)L•-ν~S•+ρ~I•-σ~ψS•I•1ξI•′=σ~ψS•I•-(ν~+α~+ρ~)I•,
where β~=β/ξ and vice versa. Rewriting the above system after renaming the parameters, we will have
(70)L′=βS-(1-L)L-μLS′=(1-L)L-νS+ρI-σSII′=σSI-(ν+α+ρ)I.Here we measured time inξ-1 and renamed the parameters β~,μ~,ν~,ρ~, and α~ as β,μ,ν,ρ, and α, σ~/ψ as σ, and L•,S•, and I• as L,S, and I. The above system (70) is dimensionless. Note that for the system ℛ0=(β/ν)(1/(1+μ)).Again we will have three types of steady states, that is,(0,0,0),(L^,S^,0), and (L*,S*,I*).Theorem 8.
(a) The trivial state(0,0,0) always exists.
(b) The disease-free steady state of the form(L^,S^,0) exists and is unique if and only if ℛ0>1.
(c) The infected state(L*,S*,I*) exists if and only if ℛ0>1 and (1-L*)L*>νS*.
(d) Ifμ<1, then for the system (70) there exist one or more than one interior steady state.
Ifμ>1 then for the system (70) there exists only one interior steady state with the L component given by l1.Proof.
(b) Since we know that the “L” component of disease-free steady state comes from
(71)π(L^)=μνβ-ν,
this implies
(72)1-L^=μνβ-ν,
which gives
(73)L^=1-μνβ-ν,
and the “S” component is given by
(74)S^=μβ-νL^=μβ-ν(1-μνβ-ν).
This steady state will exist if and only if
(75)0<μν<β-ν.
(c) It directly follows from part (c) of Theorem3.
(d) From the third equation of the system (70), we have
(76)S*=ν+α+ρσ.
As calculated earlier, the “I” component is given by
(77)I*=(β-ν)S*-μL*ν+α=π(L*)L*-νS*ν+α.
It is clear that the feasibility condition for the existence of that steady state is
(78)π(L*)L*>νS*,
where L*<L^ as proved in Theorem 3. This is equivalent to
(79)(1-L*)L*>νS*.
Observe that this holds only if L*<1, since for L*≥1, we get νS*<0, which is not true. So we must have L*<1. In this case we get
(80)L*2-L*+νS*<0.
We can write the left-hand side of (80) as
(81)(L*-K1)(L*-K2)<0,
where
(82)K1=1-1-4νS*2K2=1+1-4νS*2.
Notice that K1<K2 and both of them will exist only if S*≤1/4ν. Further from (81) it is clear that L*-K1 and L*-K2 have the opposite signs. So we have the following two cases.
Case
1  (L*-K1>0 andL*-K2<0). This means L*>K1 and L*<K2 or equivalently K1<L*<K2.
Case
2   (L*-K1<0 and L*-K2>0). This means L*<K1 and L*>K2. Since K1<K2, so this case is not possible. So we can have only one choice for L*; that is, K1<L*<K2 (with K1 and K2 given in (82)) for the existence of I*>0. The “L” component of the endemic steady state is the solution of the equation
(83)π(L*)L*+μL*=βS*.
By using the value of π(L), we get
(84)(1-L)L+μL=βS*.
This implies that for L<1(85)L2-L(μ+1)+βS*=0,
which is a quadratic equation in L and the roots of this equation are given by
(86)l1=(1+μ)-(1+μ)2-4βS*2(87)l2=(1+μ)+(1+μ)2-4βS*2.
It is clear that 0<l1<l2. Also l1 will exist only if l1<L^=1-(μν/(β-ν))<1. That is, we have
(88)(1+μ)-(1+μ)2-4βS*<2.
Simplifying and rearranging the terms yield
(89)(μ-1)<(1+μ)2-4βS*.
We have the following two cases.
Case
i  (μ<1). For this case inequality (89) is automatically satisfied.
Case
ii   (μ≥1). This is equivalent to μ-1≥0. Then on squaring both sides of inequality (89), we get
(90)1+μ2-2μ<1+μ2+2μ-4βS*.
This implies that
(91)βS*<μ.
Now consider l2; this will also exist only if l2<1; from (89), we get
(92)(1+μ)+(1+μ)2-4βS*<2
which can be rewritten as
(93)(μ-1)<-(1+μ)2-4βS*.
From the above equation notice that l2 will exist only if μ<1. Taking square of both sides of this equation we get
(94)1+μ2-2μ>1+μ2+2μ-4βS*
which implies that μ<βS*.Theorem 9.
(i) The steady states of the form(L^,S^,0) are stable if and only if ℛ0>1 and μ>2L^-1.
(ii) The endemic statel1 is locally stable and l2 is unstable whenever they exist.Proof.
(a) The proof directly follows from Theorems3 (b) and 6 (b).
(b) The Jacobian matrix is given by(95)J=[-(1-2L)-μβ0(1-2L)-ν-σIρ-σS0σIσS-(ν+α+ρ)].
At the steady state (L*,S*,I*) we get
(96)J*=[-(1-2L*)-μβ0(1-2L*)-ν-σI*-(ν+α)0σI*0].
Assume that (L*)=(1-2L*)+μ, so we have
(97)J*=[-ϕ(L*)β0ϕ(L*)-μ-ν-σI*-(ν+α)0σI*0];
the trace and determinant are given by
(98)𝒯=-(ϕ(L*)+ν+σI*)𝒟=-σI*(ν+α)ϕ(L*).
Observe that
(99)ϕ(l1)=(1-2l1)+μ=1-[(1+μ)-(1+μ)2-4βS*]+μ.
This implies
(100)ϕ(l1)=(1+μ)2-4βS*>0.
Similarly
(101)ϕ(l2)=-(1+μ)2-4βS*<0.
Now we evaluate the determinant for l1 and l2. At L*=l1, 𝒯<0 and also 𝒟<0.
AtL*=l2, we have 𝒟>0, so at least one of the eigenvalues is positive and thus this equilibrium is unstable.
As calculated earlier(102)𝒟-𝒯∑k=13‍Jkk=ϕ2(L*)(ν+σI*)+(ν+σI*)σI*(ν+α)+ϕ(L*)(ν+σI*)2-β(ϕ(L*)-μ)×[ϕ(L*)+ν+σI*].
This implies
(103)𝒟-𝒯∑k=13‍Jkk>ϕ(L*)(ν+σI*)2-β(ϕ(L*)-μ)[ϕ(L*)+ν+σI*].
Since the only point of interest here is L*=l1 (because l2 is unstable) and ϕ(l1)>0, therefore we have
(104)𝒟-𝒯∑k=13‍Jkk>-β(ϕ(l1)-μ)[ϕ(l1)+ν+σI*].
Now for the local asymptotic stability, we want
(105)ϕ(l1)-μ<0⟹(1+μ)2-4βS*<μ⟹1+μ2+2μ-4βS*<μ2;
simplification yields
(106)1+2μ<4βS*⟹μ<4βS*-12.
It can be easily shown that we can find such a μ which satisfies (106). Further this is clear from (101) that
(107)2μ<4βS*-1<(1+μ)2-1⟹μ2>0.Special Case. Here we will discuss the special case when there is no disease, that is, when I=0. Also assume that π(L)=(1-L) for L≤1. Under this condition our system given by (70) will be reduced to a two-dimensional system and is given by
(108)L′=βS-(1-L)L-μLS′=(1-L)L-νS.
Observe that (0,0) is a steady state for this system. We can write the steady state form of this system into matrix form as
(109)[-(1-L)-μβ(1-L)-ν][LS]=0
for any nontrivial steady state, L*,S*>0; the determinant of this matrix should be zero; that is,
(110)ν(1-L)+μν-β(1-L)=0.
On rearranging the terms, we get
(111)(ν-β)(1-L)+μν=0.
This implies
(112)μν=(β-ν)(1-L)μ(1-L)=βν-1.This steady state(L*,S*) will exist if β>ν. Now the Jacobian matrix at (L*,S*) is given by
(113)J=[-(1-2L*)-μβ(1-2L*)-ν].
The trace and determinant are given by
(114)𝒯=-(1-2L*)-μ-ν=-[(1-2L*)+μ+ν]𝒟=ν(1-2L*)+μν-β(1-2L*)=-(β-ν)(1-2L*)+μν=-(β-ν)[(1-L*)-L*]+μν=-(β-ν)(1-L*)+μν+(β-ν)L*.
Notice that, from (111), (β-ν)(1-L*)-μν=0. So we get
(115)𝒟=(β-ν)L>0.
So this internal steady state is locally asymptotically stable if (1-2L*)+μ+ν>0, that is, when L*<(1/2)(1+μ+ν), and unstable if (1-2L*)+μ+ν<0, that is, when L*>(1/2)(1+μ+ν). Thus we have that L⋄=(1/2)(1+μ+ν) is the point of bifurcation.Theorem 10.
Assume that the positive parametersμ and ν can be chosen such that L*>(1/2)(1+μ+ν); then L* becomes the first coordinate of a unique nontrivial steady state that is unstable spiral point. Every solution that does not start at the nontrivial steady state or the origin converges towards a periodic orbit. In particular there exists an orbitally stable periodic orbit.Proof.
Since we have shown earlier that all the solutions are bounded and we have only one interior steady state for this limiting system which is also nondegenerate, that is, all the eigenvalues of the associated Jacobian matrix evaluated at this steady state point are nonzero, so by Theorem A-15 of [19], there exists a locally asymptotically stable periodic orbit.So far we have shown that, with proper choice of the values ofμ and ν, we can find L*>0 which becomes the first coordinate of an unstable spiral point or in particular there exists a periodic solution. Here we will use the approach of the theory of Hopf bifurcation and in particular the concepts of supercritical and subcritical bifurcation to discuss the stability of these periodic orbits. Actually we need to find a Jacobian matrix of the form
(116)[0-γγ0]
with γ>0.Again consider the system equation (86). Define V=S+(ν/β)L. This implies
(117)V′=(1-L)L-νS+νβ[βS-(1-L)L-μL]=(1-L)L-νS+νS-νβ(1-L)L-μνβL.
This implies
(118)V′=Lβ[(β-ν)(1-L)-νμ].
Similarly
(119)L′=β(V-νβL)-(1-L)L-μL=βV-νL-(1-L)L-μL.
This implies
(120)L′=βV-(ν+μ)L-(1-L)L.
Let u=(β/γ)V. So we have
(121)u′=βγ(Lβ[(β-ν)(1-L)-νμ])
or
(122)u′=Lγ[(β-ν)(1-L)-νμ]L′=γu-(ν+μ)L-(1-L)L.
Observe that the L component of the interior steady state L* of this system satisfies (β-ν)(1-L*)-νμ=0. The Jacobian matrix is given by
(123)J=[01γ[(β-ν)(1-L)-νμ]+Lγ[-(β-ν)]γ-[(ν+μ)+(1-2L)]].
Evaluating the Jacobian matrix at the interior steady state (u*,L*) we get
(124)J=[0-L*γ(β-ν)γ-[(ν+μ)+(1-2L*)]].
Since we also know that at L*=L⋄ (i.e., Hopf bifurcation point) we have (ν+μ)+(1-2L⋄)=0, so at this point the Jacobian matrix will become
(125)J=[0-L⋄γ(β-ν)γ0].
Now we need to choose γ in such a way that γ=(L⋄/γ)(β-ν) which gives
(126)γ=L⋄(β-ν)>0.
So with this choice of γ, we have Hopf bifurcation. Now assume that (F1,F2) is the vector field associated with (122); that is,
(127)F1=Lγ[(β-ν)(1-L)-νμ]F2=γu-(ν+μ)L-(1-L)L.
The stability of the bifurcating periodic orbit is determined by the sign of the following number:
(128)a=γ[Fuuu1+FuLL1+FuuL2+FLLL2]+FuL1(Fuu1+FLL1)-FuL2(Fuu2+FLL2)-Fuu1Fuu2+FLL1FLL2.
Now we calculate these partial derivatives:
(129)Fuuu1=0,FuLL1=0,FuuL2=0,FLLL2=0FuL1=0,FuL2=0,Fuu1=0FL1=1γ(β-ν)(1-L)-νμ+Lγ[-(β-ν)]⟹FLL1=-1γ(β-ν)-1γ(β-ν)=-2γ(β-ν)⟹FL2=-(ν+μ)-(1-2L)⟹FLL2=2.
So by substituting these values we get
(130)a=-4γ(β-ν).
Since a<0, so the bifurcating periodic orbits are asymptotically stable; that is, we have the case of supercritical bifurcation.


## 6. Discussion

In this paper, we have introduced and analyzed a model of disease based amphibian decline. We found basic reproduction numberℛ0, which guarantees the extinction of disease population in the presence of disease as shown in Figure 1. We proved that only the disease-free steady state (0,0,0) is stable and the stability of the other two steady states is conditional and depends on some other environmental factors whenever ℛ0>1. We also have shown that not only is ℛ0 the sufficient condition for the survival of infected population but also for the case when ℛ0>1 disease may or may not persist depending on certain constraints (see Figures 1, 2, and 3). In our model we also have considered a small population size of L. The large population case will be relatively easy to handle in terms of the function π(L). In this case, since π(L) is decreasing, therefore for large population size π(L*)≈0 in (24) and this will take the form
(131)L*=βσμ(ν+α+ρ)=βμS*.Figure 1
This figure shows the time series of model (1) when ℛ0<1. It is clear that all the population will extinct in this case. The parameters used are β=0.002;μ=0.05;ν=0.02;ρ=0.1;σ=0.001;α=0.01.Figure 2
The case ofℛ0>1. It is clear that, under certain parameter values, only larvae and susceptible population will survive while the infected population will be extinct. All parameter values are the same except β=0.208.Figure 3
Phase portrait ofL versus S. The existence of a periodic solution for model (1) is evident. Parameters values are the same as of Figure 2.This is the case when we have large number of susceptibles.Further, in our model, we just consider the simple case when only the larvae in stageL are subject to intrastage competition. A future direction for this model is to consider a rather more complicated case when both secure larvae L and exposed larvae S will compete for the resources. In this case, the term π(L)L will be replaced by π(L,S)L, which may cause some oscillatory effects due to the density dependent development rate for both L and S.



---

*Source: 145398-2014-04-28.xml*