In: Advanced Math
Want a SEIRS model with stochastic transmission
project proposal under the following headings for a MSc student in
Applied Mathematical Modelling :
1)Topic
2) introduction
3) Background to the study
4)problem statement
5)aims
6)objectives
7)methodology
Abstract:
For an SEIRS epidemic model with stochastic perturbations on transmission from the susceptible class to the latent and infectious classes, we prove the existence of global positive solutions. For sufficiently small values of the perturbation parameter, we prove the almost surely exponential stability of the disease-free equilibrium whenever a certain invariant Rσ is below unity. Here Rσ<R, the latter being the basic reproduction number of the underlying deterministic model. Biologically, the main result has the following significance for a disease model that has an incubation phase of the pathogen: A small stochastic perturbation on the transmission rate from susceptible to infectious via the latent phase will enhance the stability of the disease-free state if both components of the perturbation are non-trivial; otherwise the stability will not be disturbed. Simulations illustrate the main stability theorem.
Introduction:
In recent years, a number of articles have been published on stochastic differential equation models of population dynamics of infectious diseases. In comparison with models described by ordinary differential equations (ode), the stochastic differential equation (sde) models provide of course a means of accommodating randomness in the model. Two themes of special interest in the modeling of population dynamics of a disease are the stability of equilibrium points and the optimal control of interventions such as vaccination, quarantine, public health education and others. For sde models, optimal control problems and solutions are presented in [1] of Cai and Luo, [2] of Ishikawa and in [3]. In the stochastic setting, stability of equilibria and the long term persistence or extinction of a disease in a population have been studied in most of the sde models in the literature. Such studies use different versions of stability. Stochastic perturbation has also been studied in multigroup models, such as in [4–6] for example. In many cases it has been proved that the introduction of stochastic perturbations into an ode epidemic model system can possibly render an unstable disease-free equilibrium of the ode system to become stable in the stochastic differential equation system. This phenomenon was highlighted in [7] by Chen et al., [8] by Gray et al. and in [9] for instance.
Since the basic models such as [10] by Li et al. on diseases of the SEIRS type, many variations have been presented in the literature, such as [11] of Melesse and Gumel. Starting with an ode model of SEIRS type, in this paper we study the effect of stochastic perturbations on the stability of the disease-free equilibrium of the system. The models in [7] and in [8] have perturbations of the transmission rate from the S-class to the I-class. The latter models do not include a latent infection compartment such as the E compartment in SEIR type models. The current paper is among the first studies of a disease model with a latent infection compartment, having a perturbation of the disease transmission. In the literature there are stochastic models such as in [4] and [5] having latent infection compartments, but with the stochastic perturbations not directly aimed at transmission. We prove the existence of solutions which are almost surely global and positive. We also study stability of the disease-free equilibrium. In particular, we introduce an invariant Rσ of the model that is related to the basic reproduction number R of the underlying deterministic model, with Rσ<R. With the given type of randomness in the system, we prove that there is a greater chance of the disease vanishing from the population. The main results are illustrated with simulations.
Preliminaries:
Notation-2.1
By Rn+ (resp. Rn++) we denote the set of points in Rn having only
non-negative (resp. strictly positive) coordinates.
We assume throughout the paper that we have a complete
probability space (Ω,F,P), equipped with a filtration, {Ft}t≥0,
that is right continuous and with F0 containing all the subsets
having measure zero. We consider a one-dimensional Wiener process
W(t) on this filtered probability space.
Consider the k-dimensional stochastic differential equation, for
some multi-dimensional Wiener process
B(t):dx(t)=f(x(t),t)dt+g(x(t),t)dB(t),with x(0)=x0.A solution is
denoted by x(t,x0).
Assume that f(0,t)=g(0,t)=0 for all t≥0, so that the origin point
is an equilibrium point of equation (2.1).
By L we denote the infinitesimal generator (see for instance [12])
associated with the function displayed in equation (2.1), defined
for a function V(t,x)∈C1,2(R+×Rk).
Definition-2.2
See [13]
The equilibrium x=0 of the system (2.1) is said to be almost surely
exponentially stable if for each initial value x0 in a given
subset, we have
lim supt→∞1tln∣∣x(t,x0)∣∣<0(a.s.).
The limit lim supt→∞1tlnx(t) is called the Lyapunov exponent of
x(⋅).
The following lemma was utilized in [9] and proved in [14]. For
completeness we include the simple proof.
Lemma-2.3
For k∈N, let X(t)=(X1(t),X2(t),…,Xk(t)) be a bounded Rk-valued
function and let (t0,n) be any increasing unbounded sequence of
positive real numbers. Then there is a family of sequences (tl,n)
such that, for each l∈{1,2,…,k}, (tl,n) is a subsequence of
(tl−1,n) and the sequence Xl(tl,n) converges to a chosen limit
point of the sequence Xl(tl−1,n).
Proof
Let b be an upper bound for the functions Xi(t). In the compact set
[0,b], we can choose a limit point in the closure of the set
{X1(t0,n)|n∈N} and select a convergent subsequence (t1,n) of (t0,n)
for which the limit is the chosen limit point. In the same way we
can start with the sequence (t1,n) and pick a subsequence (t2,n)
for which (X2(t2,n)) is convergent, etc. □
The function that we now introduce will be important in the
stability analysis. Consider any p∈[0,1] and let q=1−p. Now we
define the function
h:R++→R+by the rule x↦1x(p(1−x)+qx)2.
(2.2)
Let h∗:[0,1]→R+ be the function defined by
h∗(p)={4p(1−2p)(1−p)2for 0≤p≤13,for 13<p≤1.
Proposition-2.4:
Let p, h and h∗ be as above, and let h1=h|(0,1] be the (domain-)
restriction of h to (0,1]. Then h∗(p) is the absolute minimum of
h1.
Proof
If p=12, then h(x)=(4x)−1 and the result follows easily. Thus, for
the remainder of the proof we exclude the case p=12. Then we
observe that h tends to +∞ if x→0+ and also, h tends to +∞ as x→+∞.
Using calculus we find that h′(x) is continuous on R++ and has
exactly one root x0, which is x0=p⋅|q−p|−1. Therefore, the minimum
of h1 is h(x0) whenever x0≤1 and is h(1) otherwise. Further, x0≤1
if and only if 0≤p≤13. The rest of the proof follows readily. □
The model:
Melesse and Gumel [11] present a model for a disease of SEIRS type that may cause different stages of infectiousness in a patient. In a special case of the mentioned model, in this paper we study the effect of stochastic perturbations on the stability of the disease-free equilibrium. The population, which at any time t consists of N(t) individuals, is regarded as being divided into four compartments or classes. These are called the susceptible, exposed, infectious and removed classes. Their sizes, at any time t, are denoted by S(t), E(t), I(t) and R(t), respectively. The equations of motion of the system are assumed to be given by the system (3.1) of stochastic differential equations. If σ=0 then the system reduces to a system of ode, which can be called the underlying deterministic model or the underlying system of ode. For the system (3.1), the underlying system of ode coincides with a special case of the model in [11]. Inflow into the population is assumed to be all into the class of susceptibles, and it is at a rate μ0K. Additionally there is flow from the recovered class into the class of susceptibles at a rate α3R, due to loss of infection-acquired immunity. The mortality rates in the different classes are denoted by μi (i=0,1,2,3) and this allows for higher mortality rates in classes which have been affected by the disease, such as also in [15] of Beretta et al. Hence the condition (3.2) below. The symbol β denotes the effective contact rate. The parameters α1 and α2 determine the rates at which individuals in the population pass from classes E to I and (respectively) from I to R.
We further assume that W is a standard Brownian motion. The aim of the paper is to have stochastic perturbations on the transmission rate. We do this by introducing two complementary pairs of stochastic perturbation terms. The non-negative constants σ, p and q are such that σ determines the intensity of the perturbation, while p and q are the relative weights attached to the split parts of the perturbation. We assume that0≤p≤1andp+q=1(see Remark 3.3(b)).
The first pair of perturbation terms (−σpSEdW and +σpSEdW)
constitutes randomness in the transmission rate from the class S to
the class E. Let us explain the presence of the factor E (instead
of I) in this component of the perturbation. We note that, for any
equilibrium point P∗ of the underlying system of ode, there is a
proportionality,
I∗=α1α2+μ2E∗,
and this motivates the presence of the factor E in the first pair
of complementary perturbation terms. This form of the first pair of
terms is particularly significant since we are specifically
concerned with what happens near disease-free equilibrium. The
second pair of complementary perturbation terms can be understood
in view of the infection ultimately driving the susceptibles (via
the E class) into the I class. The shorter the average latent
period, the more relevant does the latter perturbation become. All
the parameters are non-negative or positive constants. So for
instance, if α3=0, then the model is said to be of the SEIR type,
but for α3≠0, the model is referred to as SEIRS.
The system of stochastic differential equations is as follows:
dS=[μ0(K−S)−βSI+α3R]dt−σS(pE+qI)dW,dE=[βSI−(α1+μ1)E]dt+σpSEdW,dI=[α1E−(α2+μ2)I]dt+σqSIdW,dR=[α2I−(α3+μ3)R]dt.
(3.1)
Throughout the paper we assume that
μ0≤min{μ1,μ2,μ3}.
(3.2)
The basic reproduction number of the underlying deterministic
model, see [11], is
R=α1βK(μ1+α1)(μ2+α2).
(3.3a)
The following invariant Rσ of the model (3.1) shall feature in the
main theorem on almost sure extinction of the I-class. In
describing Rσ we use the number h∗=h∗(p) from Section 2:
Rσ=α1βK(μ1+α1)(μ2+α2+12σ2K2h∗).
(3.3b)
We introduce the following set:
ΔK={x∈R4++|x1+x2+x3+x4≤K}.
Proposition-3.1:
Suppose that, for some T, there is a local solution
X(t)=(S(t),E(t),I(t),R(t))on t∈[0,T)
for the system, with X(t)∈R4+ for each t∈[0,T). If N(0)≤K, then
N(t)≤K for each t∈[0,T).
Proof
Given any such local solution X(t), then
d(N−K)dt=−μ0(N−K)−(μ1−μ0)E−(μ2−μ0)I−(μ3−μ0)R.
The condition X(t)∈R4++ together with μ0≤min{μ1,μ2,μ3} ensures
that
(μ1−μ0)E+(μ2−μ0)I+(μ3−μ0)R≥0.
Consequently,
d(N−K)dt+μ0(N−K)=−(μ1−μ0)E−(μ2−μ0)I−(μ3−μ0)R≤0.
Solution of the first order linear ordinary differential equation
reveals that if N(0)<K, then N(t)<K for all t∈[0,T). □
We now prove the existence of solutions which are almost surely
global and positive.
Theorem-3.2:
Given any initial value X0=(S0,E0,I0,R0)∈ΔK, then the system
(3.1) admits a unique solution X(t)=(S(t),E(t),I(t),R(t)) on t≥0,
and this solution remains in ΔK almost surely.
Proof
The coefficients of the system (3.1) are locally Lipschitz
continuous. By [13], Theorem 3.5, for the given initial value X0∈ΔK
there is a unique local solution X(t) over the interval t∈[0,τen),
where τen is the explosion time.
There is a number m0∈N which is sufficiently large to allow
S0,E0,I0,R0∈(1/m0,K). For each n∈N∩[m0,∞), let us write
Dn={t∈[0,τen):S(t)≤1n or E(t)≤1n or I(t)≤1n or R(t)≤1n}.
Then we define stopping times τn and τ∞ by taking τn to be the
infimum of Dn if Dn≠∅ and otherwise τn=∞. The set D∞ and the random
variable τ∞ are defined as
D∞={t∈[0,τen):S(t)≤0 or E(t)≤0 or I(t)≤0 or
R(t)≤0},τ∞=limn→∞τn=infD∞.
For each γ>0, let Ω(γ) be the subset of Ω defined thus:
Ω(γ)={∈Ω|τ∞(w)≤γ}.
We shall prove by contradiction that τen=∞ (a.s.). So let us assume
to the contrary that there exist T,C∈R with C>0, and with
T<τen such that P(Ω(T))=C.
Let us define the function V0(X), for X=(S,E,I,R), by the
formula:
V0(X)=lnKS+lnKE+lnKI+lnKR.
By Proposition 3.1, each of the four terms of V0(X(t)) are
non-negative for every t∈[0,τ∞).
We set up a contradiction by calculating upper and lower bounds on
expectations of V0. Firstly we calculate an upper bound. For every
u∈[0,τ∞∧T) we have
dV0(X(u))=−1S(u){[μ0(K−S(u))−βS(u)I(u)+α3R(u)]du−σS(u)(pE(u)+qI(u))dW(u)}−1E(u){[βS(u)I(u)−(α1+μ1)E(u)]du+σpS(u)E(u)dW(u)}−1I(u){[α1E(u)−(α2+μ2)I(u)]du+σqS(t)I(u)dW(u)}−1R(u){α2I(u)−(α3+μ3)R(u)}du+12{[σ(pE(u)+qI(u))]2+(σpS(u))2+(σqS(u))2}du.
We remove some negative terms and deduce the following inequality:
dV0(X(u))≤[βI(u)]du+σ(pE(u)+qI(u))dW(u)+(α1+μ1)du−σpS(u)dW(u)+(α2+μ2)du−σqS(u)dW(u)+(α3+μ3)du+12{(σ(pE(u)+qI(u)))2+(σpS(u))2+(σqS(u))2}du.
Now let
ρ=βK+(α1+μ1)+(α2+μ2)+(α3+μ3)+12(σ(p+q)K)2+12(σpK)2+12(σqK)2,
and for t∈[0,τ∞∧T], let M0(t) be
M0(t)=σ∫t0[pE(t)+qI(u)−pS(u)−qS(u)]dW(u).
Now we have the following inequality:
∫t0dV0(X(u))≤ρt+M0(t).
Therefore, for any k∈N∩[m0,∞) we have
V0(X(t∧τk))−V0(X(0))≤ρ(t∧τk)+M0(t∧τk)(a.s.).
The stochastic process M0(t) is a local martingale and therefore
for any m∈N∩[m0,∞) we have E[M0(t∧τm)]=M0(0)=0. Consequently,
E[V0(X(T∧τm))]≤ρ(T∧τm)+V0(X(0))≤ρT+V0(X(0)),
and we have the upper bound which we set out to find. We now search
for a lower bound for E[V0(X(T∧τm))]. Note that if w∈Ω(T) and we
evaluate V0(X(ζ)) for ζ=w(τm), then we get
V0(X(ζ))≥ln(mK).
We can deduce the lower bound:
E[V0(X(T∧τm))]≥Cln(mK).
These two bounds yield
Cln(mK)≤E[V0(X(T∧τm))]≤ρT+V0(X(0)).
We can choose a value of m sufficiently big, so that
Cln(mK)>ρT+V0(X(0)),
leading to a contradiction. Therefore we must have τ∞=∞ almost
surely. This completes the proof of Theorem 3.2. □
Remark 3.3
(a) In the remainder of this paper we assume that sample paths are
restricted to Ω0, which is defined as follows:
Ω0={ω∈Ω|(S(t,ω),E(t,ω),I(t,ω),R(t,ω))∈ΔK for all t≥0}.
(b) Let us briefly consider a slightly different form of the
stochastic perturbation. In the first equation of the model (3.1),
the (dS) equation, let us consider a perturbation of the form
(σ1E+σ2I)SdW, with σ1 and σ2 both non-negative and at least one of
them being non-zero. Now let σ=σ1+σ2. We set p=σ1/σ and q=σ2/σ.
Then p,q∈[0,1], p+q=1 and
(σ1E+σ2I)SdW=σS(pE+qI)dW.
The introduction of p and q simplifies the analysis when we get to
deal with the function h(⋅).
Stability theorems:
The concept of stability of a deterministic system of differential equations ramifies into different forms when dealing with stochastic differential equations. In this paper we shall focus on almost sure exponential stability, which is conceptually uncomplicated. We prove that when the basic reproduction number R of the underlying deterministic model is below unity, then the disease-free equilibrium is almost surely exponentially stable. We also prove stronger results on I(t) converging to zero, in terms of the analog Rσ of R.
Theorem 4.1
If R<1, then the disease-free equilibrium of the system (3.1) is
almost surely exponentially stable.
The proof of Theorem 4.1 will be presented following a discussion
which is relevant to all the stability results that we derive in
this paper.
Item 4.2
A construction and notation.
The following construction is crucial for the proofs of the
stability theorems. We fix a positive real number b and let a≥0 and
a1≥0. Let us write (S(t),E(t),I(t),R(t))=X(t). We define the
following stochastic processes:
z(X(t))=a(K−S(t))+bE(t)+I(t)+a1R(t),Q(t)=K−S(t)z(X(t)),Ez(t)=E(t)z(X(t)),Iz(t)=I(t)z(X(t))andRz(t)=R(t)z(X(t).
Note that for every t>0 we have
aQ(t)+bEz(t)+Iz(t)+a1Rz(t)=1.
(4.1)
Since (see Remark 3.3) we assume the sample paths to be in the
subset Ω0, it follows that z(X(t))>0 for all t>0.
LetV(X(t))=lnz(X(t)).
For every sample path w of the Wiener process W(t), there exists an
unbounded increasing sequence (τwn) of positive time values for
which
lim supt→∞LV(X(t,w))=limn→∞LV(X(τwn,w)).
Fix such a sequence. Then by Lemma 2.3 there exists a subsequence
(twn) for which the following limits
exist:q=limn→∞Q(X(twn,w)),f=limn→∞Ez(X(twn,w)),i=limn→∞Iz(X(twn,w)),r=limn→∞Rz(X(twn,w))ands=limn→∞S(twn,w).
(4.2a)
Let us write
Λ(w)=limn→∞L(V(X(twn,w))).
(4.2b)
Item 4.3
A useful inequality for LV(X(u)).
Using Itô’s formula we can express the stochastic process V(X(t))
as
V(X(t))=V(X(0))+∫t0LV(X(u))du+M1(t),
with M1(t) being the Itô integral
M1(t)=∫t0(σa)pS(u)E(u)+qS(u)I(u)z(X(u))+(σb)pS(u)E(u)z(X(u))+(σ)qS(u)I(u)z(X(u))dW(u).
Since
(σa)pS(u)E(u)+qS(u)I(u)z(X(u))+(σb)pS(u)E(u)z(X(u))+(σ)qS(u)I(u)z(X(u))≤σa(pK1b+qK)+σbpK1b+σqK,
and the latter is a (bounded) fixed number, it follows that
limt→∞1t∫t0[(σa)pS(u)E(u)+qS(u)I(u)z(X(u))+(σb)pS(u)E(u)z(X(u))+(σ)qS(u)I(u)z(X(u))]2du<∞.
Thus we may apply the strong law of large numbers for local
martingales, as from [16] for instance, and we deduce that
limt→∞1tM1(t)=0(a.s.).
Since also
limt→∞1tV(0)=0,
it follows that
lim supt→∞1tV(X(t))=lim supt→∞1t∫t0LV(X(u))du(a.s.).
(4.3)
Now we calculate LV(X(t)).
LV(X)=1z[−a(μ0(K−S)−βSI+α3R)+b(βSI−(α1+μ1)E)+(α1E−(α2+μ2)I)+a1(α2I−(α3+μ3)R)]−(σS)22z2(a2(pE+qI)2+2ab(pE+qI)pE+2a(pE+qI)qI)−(σS)22z2((pbE)2+2pqbEI+(qI)2).
For our further analysis it will suffice to have a suitable
function dominating LV. From the last equation we obtain
LV(X)≤1z[−a(μ0(K−S)−βSI)+b(βSI−(α1+μ1)E)+(α1E−(α2+μ2)I)+a1(α2I−(α3+μ3)R)]−(σS)22z2((pbE)2+2pqbEI+(qI)2).
Now we introduce notation from Item 4.2, to obtain the following inequality:
LV(X)≤−a(μ0Q−βSIz)+b(βSIz−(α1+μ1)Ez)+(α1Ez−(α2+μ2)Iz)+a1(α2Iz−(α3+μ3)Rz)−(σS)22(pbEz+qIz)2.
(4.4)
Remark 4.4
In the proofs of the stability theorems we shall need to prove
that, for paths w∈Ω0,
lim supt→∞1t∫t0LV(X(u,w))du<0.
To this end we note here that it suffices to prove that
lim supt→∞LV(X(t,w))<0.
Proof of Theorem 4.1
Since R<1 we have the inequality
cβK−(α2+μ2)<0,with c=α1μ1+α1.
Now we can choose a>0 sufficiently small such that
(c+a)βK−(α2+μ2)<0.
We can also choose positive real numbers a0 and a1 sufficiently
small such that
(c+a+a0α1+μ1)βK−(α2+μ2)+a1α2<0.
(4.5)
Now let b=c+a0(μ1+α1)−1, and with these values of a, b and a1 we
define z(X(t)) as in Item 4.2. It suffices to prove that Λ<0. We
modify the inequality (4.4) ignoring the last term, using S<K
and noting that b(α1+μ1)−α1=a0. Then we deduce the following
inequality:
LV(X)≤−aμ0Q+(a+b)βKIz−a0Ez−(α2+μ2)Iz+a1α2Iz−a1(α3+μ3)Rz.
(4.6)
We form limits as in equations (4.2a) and (4.2b), and after
rearranging the terms we obtain
The coefficients of q, f and r are negative, and by inequality
(4.5) the coefficient of i is also negative. Since aq+bf+i+a1r=1,
it follows that at least one of these limits {q,f,i,r} must be
non-zero. Hence Λ<0 and the proof is complete.
Theorem 4.1 implies that if a disease-free equilibrium is locally asymptotically stable with respect to the underlying ode-system, then it is almost surely exponentially stable with respect to the stochastic model, in particular, the perturbations do not disrupt the stability of the disease-free equilibrium.
Remark 4.5
In the sequel we shall use a special form of the inequality (4.4),
taking a=a1=0. Then z=bE+I and from (4.4) we obtain the
inequality.
(4.7)
We now present the main result of this paper, which proves that the
stochastic perturbation improves the stability of the disease-free
equilibrium for small values of the perturbation parameter.
Theorem 4.6
If the following conditions hold:
(1)
Rσ<1,
(2)
σ2≤cβKh∗ with c=α1μ1+α1,
then (E(t),I(t)) almost surely converges exponentially to 0.
Proof
Let us assume the conditions (1) and (2) of the theorem to hold. In
particular then, the condition (1) is equivalent to
cβK−(α2+μ2)−12σ2K2h∗<0.
Choose any positive number a3 to be sufficiently small such
that
(c+a3α1+μ1)βK−(α2+μ2)−12σ2K2h∗<0.
This can be written
bβK−(α2+μ2)−12σ2K2h∗<0,with b=c+a3α1+μ1.
(4.8)
Using this value of b in the inequality (4.7) while a=0=a1,
yields
LV(X)≤Izg(S)−a3Ez,where g(S)=bβS−(α2+μ2)−(σS)22h∗.
(4.9)
The quadratic function g(ζ) of equation (4.9) reaches an absolute
maximum when ζ=bβ(σ2h∗)−1. By assumption (2) it follows that
cβ(σ2h∗)−1≥K. Therefore also
bβσ2h∗>cβσ2h∗≥K.
Considering that 0<S≤K, it follows that g(S)≤g(K) and
therefore
LV(X)≤Izg(K)−a3Ez.
Therefore (see the notation in Item 4.2) we have
Λ≤ig(K)−a3f.
Now we observe that g(K) coincides with the left hand side of the
inequality (4.8). Since at least one of f and i must be non-zero,
it follows that Λ<0. This completes the proof of Theorem
4.6. □
In Section 5 we shall further reflect on Theorem 4.6. Also, Theorem 4.8 below combines very well with Theorem 4.6. However, while our main result Theorem 4.6 focused on small perturbations, let us briefly address the case of larger perturbations. The following theorem asserts that, for sufficiently large values of the perturbation parameter σ, the disease will eventually vanish from the population.
Theorem 4.7
The pair (E(t),I(t)) almost surely converges exponentially to 0
if
σ4>(cβ)2(α2+μ2)h2∗with c=α1μ1+α1.
Proof
Let a4 be sufficiently small to support the inequality:
σ4>[(c+a4μ1+α1)β]2(α2+μ2)h2∗,
and let b=c+a4(α1+μ1)−1. Now we revisit the construction presented
under Item 4.2, using the constant b as selected and a=0=a1. Then
similar to the proof of Theorem 4.6, we obtain an inequality,
LV(X)≤Izg(S)−a4Ez,where g(S)=bβS−(α2+μ2)−(σS)22h∗.
The inequality σ4>(bβ)2(α2+μ2)h2∗ is equivalent to
(bβσ2h∗)2−(α2+μ2)<0.
Therefore,
g(S)==<<−σ2h∗2(S2−2bβσ2h∗S)−(α2+μ2)−σ2h∗2(S−bβσ2h∗)2+(bβσ2h∗)2−(α2+μ2)−σ2h∗2(S−bβσ2h∗)20.
This implies that Λ≤ig(s)−a4f<0, and the proof is
complete. □
Theorem 4.8
If (E(t),I(t)) almost surely converges exponentially to 0, then
limt→∞S(t)=K(a.s.)andlimt→∞R(t)=0(a.s.).
Proof
This proof is by contradiction, so let us suppose, to the contrary,
that (on a subset Ωlim of Ω of positive measure) we have
limt→∞(K−S(t))+R(t)>0. Let z be as in Item 4.2, with a=b=a1=1.
Then, since limt→∞E(t)=0 while limt→∞(K−S(t))+R(t)>0, it follows
that f=0 on Ωlim. Similarly it follows that i=0 on Ωlim. Therefore
from the inequality (4.4) it follows that on Ωlim we have
Λ≤−μ0q−(α3+μ3)r.
Therefore, Λ<0. This implies that z converges to 0 and,
consequently, that limt→∞(K−S(t))+R(t)=0, which is a contradiction.
This completes the proof. □
Remark 4.9
(a) Theorem 4.6 is much more significant than Theorem 4.7 because
in disease modeling, in practice one is more interested in smaller
perturbations rather than the larger perturbations. Let us denote
the bounds on σ specified in Theorems 4.6 and 4.7 by θ1 and θ2
respectively. If θ1>θ2, then these theorems can be combined,
guaranteeing the disease-free equilibrium to be almost surely
exponentially stable irrespective of the magnitude of σ.
(b) Of course, Theorem 4.8 serves to extend Theorems 4.6 and 4.7.
Simulations:
Theorem 4.6 suggests that Rσ is an approximation for a threshold that decides stability in a way similar to the basic reproduction number. Simulations show that it is a rather useful approximation. For a non-negative stochastic process, almost sure convergence to 0 can be tested by computing the (approximation over finitely many paths, of the) mean of sample paths. If the mean of I is not asymptotically stable, then I is not almost surely exponentially stable. The simulations that were run produced trajectories of the mean of I which consistently appears to converge to a value which is smaller than, or at least not bigger than, in the deterministic case.
These simulations are obtained by considering an influenza infection of the type in [11] and [17]. The relevant parameter values for α1, α3, μ0, μ1 and μ3 are taken directly from [11] and other parameters values are derived. Our value for 1/α2 is obtained by taking the sum of the average times spent in the I(⋅)-compartments (of [11]). The value of μ2 is taken as 1.025μ0. The parameter β is not kept fixed in these simulations. Here we note that infections that are aerially transmitted will spread faster when people are in high density locations with poor ventilation. So for instance the same disease has a higher value for the effective contact rate, when considered in a concentration camp as compared to an ordinary small village or rural area.
Initial values in millions are S(0)=6, E(0)=2, I(0)=0.9 and R(0)=1.05.
Using the Euler-Maruyama scheme, we simulate the trajectories of one sample path of the I-class of the model,
the mean of I, taken over 1,000 sample paths and indicated as ‘I (1000ave)’, and
the I-class trajectory ‘I determ’ (broken line) of the underlying deterministic model,
with the parameter values as explained above. We take K=10 (in millions)and the values of the other parameters are all given, with one-day as time unit, in the accompanying Table 1
Parameter Value
μ0 (60×365)−1=4.566×10−5
μ1 4.566 × 10−5
μ2 1.025μ0=4.680×10−5
μ3 4.566 × 10−5
α1 1.9−1 = 0.5263
α2 0.2
α3 (83.33)−1=0.012
Table 1 Numerical values of the fixed parameters
Parameters such as p and σ are difficult to compute. We choose
p=0.333 for simulation. The parameters β and σ are varied in order
to obtain different values of Rσ.
In Figure 1 we show trajectories for the case β=0.0235 and
σ=0.05, for which we get Rσ=1.1429. In this case we cannot deduce
stability of the disease-free equilibrium from Theorem 4.6, since
condition (1) is not satisfied. In fact we observe the mean value
of I as seeming to converge to a positive value. The given I-path
also does not seem to converge to 0. We do note, however, that
after 400 days (a relatively long period) it gives a mean I value
(computed as 6,600), which is lower than in the deterministic case
(19,400).
Condition (1) of Theorem 4.6 is violated.
In Figure 2 we show a case in which we take β=0.021, σ=0.08, and
we calculate R=1.0497 and Rσ=0.9800≤1. However, condition (2) is
not satisfied. The I-trajectory shown does not seem to converge to
0, although the mean seems to converge to 0. This further
demonstrates the need for condition (2) in Theorem 4.6 (other than
condition (2) just being utilized in the proof). In this case the
deterministic model has a non-trivial equilibrium value I∗ for I.
For the stochastic case we observe that the mean seems to converge
to a value smaller than I∗.
Condition (2) of Theorem 4.6 is violated.
In Figure 3 we use β=0.0215, σ=0.05, and then we have R=1.0247 and
Rσ=0.99700. This time the parameters selection fulfills all the
conditions of Theorem 4.6, and indeed what we see appears to be in
line with the assertion of the theorem.
Conclusion:
we conclude that by constructed an SEIRS model, with stochastic perturbations which can be viewed as linked to the transmission rate out of the class of susceptibles. We proved that the system of stochastic differential equations permits solutions that are almost surely global and positive. The model permits a disease-free equilibrium which we showed to be almost surely exponentially stable whenever the basic reproduction number of the underlying deterministic model is below unity, and even slightly beyond under given conditions.
Biologically we observe, in particular, the following effect of a stochastic perturbation on the disease transmission in the case of a deterministic compartmental model which allows for a latently infectious class. Given a small stochastic perturbation on the transmission rate from susceptible to infectious via the latent phase, the stability of the disease-free state will be improved if both components of the perturbation are non-trivial. If any one of the components of the perturbation is zero, then the stability will not be disturbed. The simulations confirm the proven results and also provide further insights, such as about the behavior of the mean of the I-class trajectories.