Global properties of virus dynamics with B-cell impairment

Ahmed M. Elaiw; Safiya F. Alshehaiween; Aatef D. Hobiny

doi:10.1515/math-2019-0113

Article Open Access

Global properties of virus dynamics with B-cell impairment

Ahmed M. Elaiw , Safiya F. Alshehaiween and Aatef D. Hobiny

Published/Copyright: December 13, 2019

Published by

Become an author with De Gruyter Brill

Submit Manuscript Author Information Explore this Subject

$Open Mathematics$

From the journal Open Mathematics Volume 17 Issue 1

Abstract

In this paper we construct a class of virus dynamics models with impairment of B-cell functions. Two forms of the incidence rate have been considered, saturated and general. The well-posedness of the models is justified. The models admit two equilibria which are determined by the basic reproduction number R₀. The global stability of each equilibrium is proven by utilizing Lyapunov function and LaSalle’s invariance principle. The theoretical results are illustrated by numerical simulations.

Keywords: virus dynamics; global stability; Lyapunov function; B-cell impairment

MSC 2010: 34D20; 34D23; 37N25; 92B05

1 Introduction

The study of within-host virus dynamics using mathematical modeling has been an interesting topic to research in the last decades. A proper model could provide insights of a better understanding of the virus dynamics and clinical treatments used to fight against it. In an infection process, the interaction between viruses and cells can be seen as an ecological system within the infected host. A wide of mathematical models focused on exploring the interaction between three basic compartments, uninfected cells (U), infected cells producing viruses (I) and viruses (P). A basic model of virus dynamics was originally developed by Nowak and Bangham [1] which has become highly used by experimentalists and theorists (see e.g., Nowak and May [2]). The model presented in [1] is given by:

U˙=ϱ−γU−ωUP, (1)

I˙=ωUP−βI, (2)

P˙=ϰI−ξP, (3)

where U, I and P are the concentrations of uninfected cells, infected cells and viruses, respectively. The parameters ϱ, γ, ω, β, ϰ and ξ are positive. The full description of the model was given in [1]. A huge number of papers have been published as extension of the basic model (see, e.g., [3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22]).

The immune response plays a critical role in controlling the virus spreading. The specificity and memory in adaptive immune responses are the responsibility of lymphocytes. B cells and T cells are the two main types of lymphocytes. The function of T cells is to recognize and kill infected cells, while the function of B cells is to produce antibodies which bind to virus particles and mark it as a foreign structure for elimination by other cells of the immune system. Antibody alone can neutralize, and thus protect against, viruses [23]. The virus dynamics model with B cell immune response was presented by Murase et al. [24] as

U˙=ϱ−γU−ωUP, (4)

I˙=ωUP−βI, (5)

P˙=ϰI−ξP−ρPC, (6)

C˙=εPC−μC, (7)

where C is the concentration of B cells. Many extended models are developed with B cell immune response (see, e.g., [25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36]).

In certain circumstances, some viruses can suppress immune response or even destroy it especially when the load of viruses is too high. Models with T cell immune impairment were studied several times (see, e.g., [37, 38, 39, 40]). In addition, there are factors affect B cell function and cause the impairment of B cell [41, 42, 43]. These factors include the following; malnutrition, tumors, cytotoxic drugs, irradiation, aging, trauma, some diseases (e.g., diabetes) and immunosuppression by microbes, e.g., malaria, measles virus but especially HIV [23]. In a very recent work, Miao et al. [44] have proposed a virus dynamics model which includes: humoral impairment, time delay, reaction-diffusion, and logistic growth of the target cells. Due to the complexity of the model presented in [44], the global stability analysis of the model’s equilibria did not studied. Studying the global stability of equilibria for virus dynamics models will give us a detailed information and enhances our understanding about the virus dynamics. Therefore, many mathematician have paid great efforts to study global stability of systems in virology (see, e.g., [7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19] and [45, 46, 47, 48, 49, 50, 51, 52, 53, 54]) and epidemiology (see, e.g., [55, 56, 57]).

In [44], the incidence rate of infection is given by bilinear. In reality, the bilinear incidence may not accurate to characterize the virus dynamics during different stages of infection especially when the concentration of the viruses is high [8]. Therefore, in the present paper, we propose viral infection model with B-cell impairment and with two nonlinear forms of the incidence rate, saturation and general. We show that the solutions of the model are nonnegative and bounded. The global stability of the equilibria is established by constructing Lyapunov functions and applying LaSalle’s invariance principle.

2 Model with saturation

In this section we propose a virus dynamics model including B-cell impairment and saturated incidence as:

U˙=ϱ−γU−ωUP1+αP, (8)

I˙=ωUP1+αP−βI, (9)

P˙=ϰI−ξP−ρPC, (10)

C˙=εP−μC−ϑCP, (11)

where, ϑ PC is the B-cell impairment rate and α ≥ 0 is a saturation constant.

2.1 Basic properties

We define the compact set

Ω=U,I,P,C∈R≥04:0≤U≤s1,0≤I≤s1,0≤P≤s2,0≤C≤s3 (12)

where s_i > 0, i = 1, 2, 3.

Proposition 1

The set Ω is positively invariant for model (8)-(11).

Proof

We have

U˙∣U=0=ϱ>0,I˙∣I=0=ωUP1+αP≥0, when U,P≥0,P˙∣P=0=ϰI≥0, when I≥0,C˙∣C=0=εP≥0, when P≥0.

Thus R≥04 is positively invariant for model (8)-(11). Let F(t)=U+I+β2ϰP+βξ4εϰC, then

F˙(t)=ϱ−γU−ωUP1+αP+ωUP1+αP−βI+β2I−βξ2ϰP−βρ2ϰPC+βξ4ϰP−βξμ4εϰC−βξϑ4εϰPC=ϱ−γU−β2I−βξ4ϰP−βρ2ϰ+βξϑ4εϰPC−βξμ4εϰC≤ϱ−γU−β2I−βξ4ϰP−βξμ4εϰC≤ϱ−σU+I+β2ϰP+βξ4εϰC=ϱ−σF(t),

where σ=minγ,β2,ξ2,μ. Then,

F(t)≤ϱσ+F(0)−ϱσe−σt.

Then, 0 ≤ F(t) ≤ s₁, if F(0) ≤ s₁ for t ≥ 0 wheres₁ = ϱσ . Hence, 0 ≤ U(t), I(t) ≤ s₁, 0 ≤ P(t) ≤ s₂ and 0 ≤ C(t) ≤ s₃ for all t ≥ 0 if U(0)+I(0)+β2ϰP(0)+βξ4εϰC(0)≤s1, where s2=2ϰs1β and s3=4εϰs1βξ. This guarantees that the solutions of the model are bounded.□

The basic infection reproduction number for model (8)-(11) is given by:

R0=ϱωϰβξγ.

Lemma 1

Consider model (8)-(11), we have

if R₀ ≤ 1, then the model has only one equilibrium point EP₀,
if R₀ > 1, then the model has two equilibria EP₀ and EP₁.

Proof

At any equilibrium EP(U, I, P, C) we have

ϱ−γU−ωUP1+αP=0, (13)

ωUP1+αP−βI=0, (14)

ϰI−ξP−ρPC=0, (15)

εP−μC−ϑCP=0. (16)

From equations (13)-(16) we get an infection-free equilibrium EP₀ = (U₀, 0, 0, 0), where U0=ϱγ and a unique endemic equilibrium EP₁ = (U₁, I₁, P₁, C₁), where

U1=ϱ1+αP1γ+αγP1+ωP1,I1=ϱωP1(1+αP1)βγ+αγP1+ωP1,P1=−b+b2−4ac2a,C1=εP1ϑP1+μ,

where

a=βξγαϑ+βξϑω+βγρεα+βρεω,b=βξγϑ+βξγμα+βρεγ+βξμω−ϱωϑϰ,c=βξγμ1−R0.

Then the equilibrium EP₁ exists when R₀ > 1.□

2.2 Global properties

Define a function G(u) = u – 1 – ln u. Clearly G(u) ≥ 0, for u > 0 and G(1) = 0. The global stability analysis of the two equilibria of model (8)-(11) will be established in the next theorems.

Theorem 1

Let R₀ < 1, then the infection-free equilibrium EP₀ of model (8-(11) is globally asymptotically stable.

Proof

Construct a Lyapunov function L₀(U, I, P, C) as

L0=U0GUU0+I+βϰP+βξεϰ1−R0C.

Calculating dL0dt as:

dL0dt=1−U0Uϱ−γU−ωUP1+αP+ωUP1+αP−βI+βϰϰI−ξP−ρPC+βξεϰ1−R0εP−μC−ϑCP=−γ1−U0UU−U0−βρϰPC−βξϑεϰ1−R0PC+ωU01+αP−βξϰ+βξϰ(1−R0)P−βξμεϰ1−R0C=−γ(U−U0)2U−βρϰ+βξϑεϰ1−R0PC−βξμεϰ1−R0C−βξαR0ϰ1+αPP2.

Since R₀ < 1, then for all U, P, C > 0 we have dL0dt ≤ 0. Moreover, dL0dt = 0 when U(t) = U₀ and P(t) = C(t) = 0. Let D0=(U,I,P,C):dL0dt=0 and M₀ be the largest invariant subset of D₀. The trajectory of model (8)-(11) tends to M₀ [58]. All the elements of M₀ satisfyU(t) = U₀ and P(t) = C(t) = 0. Then Eq. (10) we get

P˙(t)=0=ϰI(t), ⇒I(t)=0.

Hence, M₀ = {EP₀}. From LaSalle’s invariance principle, we derive that if R₀ < 1, then EP₀ is globally asymptotically stable.□

Theorem 2

Let R₀ > 1, then the endemic equilibrium EP₁ of model (8)-(11) is globally asymptotically stable.

Proof

Construct a Lyapunov function L₁(U, I, P, C) as

L1=U1GUU1+I1GII1+βϰP1GPP1+βρ2ϰε−ϑC1C−C12.

Note that from the equilibrium condition Eq. (16) that

ε−ϑC1=μC1P1>0.

Then dL1dt is given by:

dL1dt=1−U1Uϱ−γU−ωUP1+αP+1−I1IωUP1+αP−βI+βϰ1−P1PϰI−ξP−ρPC+βρϰε−ϑC1C−C1εP−μC−ϑCP=1−U1Uϱ−γU+ωU1P1+αP−ωUP1+αPI1I+βI1−βξϰP−βρϰPC−βP1PI+βξϰP1+βρϰP1C+βρϰε−ϑC1C−C1εP−μC−ϑCP.

From the equilibrium conditions, we have:

ϱ=γU1+ωU1P11+αP1,βI1=ωU1P11+αP1,ϰI1=ξP1+ρP1C1,εP1=μC1+ϑP1C1.

Utilizing the conditions of EP₁, we get

dL1dt=1−U1UγU1+ωU1P11+αP1−γU+ωU1P1+αP−ωUP1+αPI1I+ωU1P11+αP1−βξϰP−βρϰPC−βP1PI+βξϰP1+βρϰP1C+βρϰε−ϑC1C−C1εP−μC−ϑCP−εP1+μC1+ϑC1P1.

Simplifying the result, we obtain

dL1dt=−γU−U12U+ωU1P11+αP12−U1U+ωU1P1+αP−ωUP1+αPI1I−βξϰP−P1−βρϰP−P1C−ωU1P11+αP1P1IPI1+βρϰε−ϑC1C−C1εP−μC−ϑCP−εP1+μC1+ϑC1P1+ϑC1P−ϑC1P+βρϰP−P1C1−βρϰP−P1C1=−γU−U12U+ωU1P11+αP12−U1U−P1IPI1+ωU1P11+αP1P1+αP1P11+αP−ωU1P11+αP11+αP1UPI11+αPU1P1I−βϰP−P1ξ+ρC1−βρϰP−P1C−C1+βρεϰε−ϑC1C−C1P−P1−βρμϰε−ϑC1C−C12−βρϑϰε−ϑC1C−C12P−βρϑC1ϰε−ϑC1C−C1P−P1=−γU−U12U+ωU1P11+αP13−U1U−P1IPI1−1+αP1UPI11+αPU1P1I+ωU1P11+αP1P1+αP1P11+αP−ωU1P11+αP1PP1−βρϰP−P1C−C1+βρ(ε−ϑC1)ϰ(ε−ϑC1)(C−C1)(P−P1)−βρμ+ϑPϰε−ϑC1C−C12=−γU−U12U+ωU1P11+αP14−U1U−P1IPI1−1+αP1UPI11+αPU1P1I−1+αP1+αP1−βρμ+ϑPϰε−ϑC1C−C12+ωU1P11+αP1P1+αP1P11+αP−PP1−1+1+αP1+αP1=−γU−U12U+ωU1P11+αP14−U1U−P1IPI1−1+αP1UPI11+αPU1P1I−1+αP1+αP1−βρμ+ϑPϰε−ϑC1C−C12−αωU1P−P121+αP1+αP12.

Using geometrical mean (GM) and arithmetical mean (AM) inequality

AM≥GM, (17)

we get

4≤U1U+P1IPI1+1+αP1UPI11+αPU1P1I+1+αP1+αP1.

Thus for all U, I, P, C > 0 we have dL1dt ≤ 0. In addition dL1dt = 0 when U = U₁, I = I₁, P = P₁ and C = C₁. Let D1=W1(U,I,P,C):dL1dt=0 and M₁ be the largest invariant subset of D₁. Clearly M₁ = {EP₁}. Applying LaSalle’s invariance principle we obtain that if R₀ > 1, then EP₁ is globally asymptotically stable.□

3 Model with general incidence rate

In this section we propose a model with more general incidence rate function Θ(U, P) as:

U˙=ϱ−γU−Θ(U,P), (18)

I˙=Θ(U,P)−βI, (19)

P˙=ϰI−ξP−ρPC, (20)

C˙=εP−μC−ϑPC, (21)

We need the following Assumptions of the function Θ(U, P):

Θ(U, P) is continuously differentiable, Θ(U, P) > 0, and Θ(0, P) = Θ(U, 0) = 0 for all U > 0 and P > 0,
∂Θ(U,P)∂U>0, ∂Θ(U,P)∂P>0, and ∂Θ(U,0)∂P>0 for all U > 0 and P > 0,
ddU∂Θ(U,0)∂P>0 for all U > 0,
Θ(U,P)P is decreasing with respect to P for all P > 0.

One can show that the set Ω given by Eq. (12) is positively invariant for model (18)-(21).

Lemma 2

Assume that Assumptions (A1)-(A4) are satisfied, then there exists a threshold parameter R0G > 0 such that:

if R0G < 1, then the model has only one equilibrium point EP₀; and
if R0G > 1, then the model has two equilibria EP₀ and EP₁.

Proof

At any equilibrium EP(U, I, P, C) we have

ϱ−γU−Θ(U,P)=0, (22)

Θ(U,P)−βI=0, (23)

ϰI−ξP−ρPC=0, (24)

εP−μC−ϑCP=0. (25)

From Eq. (25), we have

C=εPμ+ϑP, (26)

and from Eq. (24), we get

I=ξPϰ+ρεϰP2μ+ϑP. (27)

Now from Eqs. (27) and (22)-(23), we obtain

U=ϱγ−βξγϰP+βρεγϰP2μ+ϑP. (28)

Let

Ψ(P)=ϱγ−βξγϰP−βρεγϰP2μ+ϑP.

Therefore, we can write U as U = Ψ(P). Note that Ψ(0) = ϱγ .

From Eqs. (27) and (22)-(23), we have

Θ(Ψ(P),P)−βξPϰ+ρεϰP2μ+ϑP=0. (29)

Observe that, P = 0 is a solution of Eq. (29). Then from Eqs. (26)-(28), we have U = U₀, I = 0, and C = 0. Then we get an infection-free equilibrium EP₀ = (U₀, 0, 0, 0).

Let

H(P)=Θ(Ψ(P),P)−βξPϰ+ρεϰP2μ+ϑP,

then H(0) = 0. Let P be such that Ψ(P) = 0, i.e.,

U0−βξγϰP¯−βρεγϰP¯2μ+ϑP¯=0,

which gives

(βρε+βξϑ)P¯2+(βξμ−γϰϑU0)P¯−γϰμU0=0. (30)

Thus, the positive solution of Eq. (30) is given by

P¯=(γϰϑU0−βξμ)+(βξμ−γϰϑU0)2+4γϰμU0(βρε+βξϑ)2(βρε+βξϑ).

We can see that

H(P¯)=Θ(0,P¯)−βξP¯ϰ+ρεϰP¯2μ+ϑP¯=−βξP¯ϰ+ρεϰP¯2μ+ϑP¯<0.

Moreover,

H′(P)=∂Θ(U,P)∂P+Ψ′(P)∂Θ(U,P)∂U−βξϰ−βρε2μP+ϑP2ϰ(μ+ϑP)2.

Assumption (A1) implies that ∂Θ(U0,0)∂U=0, then

H′(0)=∂Θ(U0,0)∂P−βξϰ=βξϰϰβξ∂Θ(U0,0)∂P−1.

Therefore, if ϰβξ∂Θ(U0,0)∂P>1, then H^′(0) > 0 and ∃P₁ ∈ (0, P) such that H(P₁) = 0. Let us define

R0G=ϰβξ∂Θ(U0,0)∂P,

which represents the basic reproduction number. Now, let

g(U)=ϱ−γU−Θ(U,P1)=0.

Then we have g(0) = ϱ > 0 and g(U₀) = –Θ(U₀, P₁) < 0. Assumption (A2) implies that g(U) is a strictly decreasing function of U, and then there exists a unique U₁ ∈ (0, U₀) such that g(U₁) = 0. Moreover, from Eqs. (26) and (27), we have

C1=εP1μ+ϑP1>0,I1=ξP1ϰ+ρεϰP12μ+ϑP1>0.

Therefore, an endemic equilibrium EP₁ = (U₁, I₁, P₁, C₁) exists if R0G > 1.□

3.1 Global stability of equilibria

The global stability analysis of the two equilibria of model (18)-(21) will be investigated in this section.

Theorem 3

Let R0G > 1, then the infection-free equilibrium EP₀ of model (18)-(21) is globally asymptotically stable.

Proof

Construct a Lyapunov function Z₀(U, I, P, C) as

Z0=U−U0−∫U0UlimP→0+Θ(U0,P)Θ(η,P)dη+I+βϰP+βξεϰ1−R0GC.

Calculating dZ0dt as:

dZ0dt=1−limP→0+Θ(U0,P)Θ(U,P)ϱ−γU−Θ(U,P)+Θ(U,P)−βI+βϰϰI−ξP−ρPC+βξεϰ1−R0GεP−μC−ϑPC=1−limP→0+Θ(U0,P)Θ(U,P)ϱ−γU+Θ(U,P)limP→0+Θ(U0,P)Θ(U,P)−βξϰP−βρϰPC+βξϰ1−R0GP−βξμεϰ1−R0GC−βξϑεϰ1−R0GPC=ϱ−γU1−limP→0+Θ(U0,P)Θ(U,P)+Θ(U,P)limP→0+Θ(U0,P)Θ(U,P)−βξR0GϰP−βξμεϰ1−R0GC−βρϰ+βξϑεϰ1−R0GPC=γU01−UU01−∂Θ(U0,0)/∂P∂Θ(U,0)/∂P+βξR0GϰϰΘ(U,P)βξR0GP∂Θ(U0,0)/∂P∂Θ(U,0)/∂P−1P−βξμεϰ1−R0GC−βρϰ+βξϑεϰ1−R0GPC.

From the Assumptions, we have the first term is less than or equal to zero. In addition,

Θ(U,P)P≤limP→0+Θ(U,P)P=∂Θ(U,0)∂P,

for all U > 0. Then

ϰΘ(U,P)βξR0GP∂Θ(U0,0)/∂P∂Θ(U,0)/∂P≤ϰβξR0G∂Θ(U0,0)∂P=1.

It implies that

dZ0dt≤γU01−UU01−∂Θ(U0,0)/∂P∂Θ(U,0)/∂P−βξμεϰ1−R0GC−βρϰ+βξϑεϰ1−R0GPC.

Therefore, if R0G < 1, then dZ0dt≤0 for all U, P, C > 0. Similar to the proof of Theorem 1, one can show that EP₀ is globally asymptotically stable.□

Theorem 4

Let R0G > 1, then the endemic equilibrium EP₁ of model (18)-(21) is globally asymptotically stable.

Proof

Define Z₁(U, I, P, C) as

Z1=U−U1−∫U1UΘ(U1,P1)Θ(η,P1)dη+I1GII1+βϰP1GPP1+βρ2ϰε−ϑC1C−C12.

Then dZ1dt can be calculated as:

dZ1dt=1−Θ(U1,P1)Θ(U,P1)ϱ−γU−Θ(U,P)+1−I1IΘ(U,P)−βI+βϰ1−P1PϰI−ξP−ρPC+βρϰε−ϑC1C−C1εP−μC−ϑPC.

Using the equilibrium condition, εP₁ – μC₁ – ϑP₁C₁ = 0, we have

dZ1dt=1−Θ(U1,P1)Θ(U,P1)ϱ−γU+Θ(U,P)Θ(U1,P1)Θ(U,P1)−I1IΘ(U,P)+βI1−βξϰP−βρϰPC−βP1PI+βξϰP1+βρϰP1C+βρϰε−ϑC1C−C1εP−μC−ϑPC−εP1+μC1+ϑP1C1+ϑPC1−ϑPC1+βρϰP−P1C1−βρϰP−P1C1=1−Θ(U1,P1)Θ(U,P1)ϱ−γU+Θ(U,P)Θ(U1,P1)Θ(U,P1)−I1IΘ(U,P)+βI1−βξϰP−βρϰPC−βP1PI+βξϰP1+βρϰP1C+βρεϰε−ϑC1C−C1P−P1−βρμϰε−ϑC1C−C12−βρϑPϰε−ϑC1C−C12−βρϑC1ϰε−ϑC1C−C1P−P1+βρϰP−P1C1−βρϰP−P1C1.

From the equilibrium conditions, we have

βI1=Θ(U1,P1),ϱ=γU1+βI1.

Applying these conditions, we obtain

dZ1dt=1−Θ(U1,P1)Θ(U,P1)γU1+βI1−γU+βI1Θ(U,P)Θ(U,P1)−I1IΘ(U,P)+βI1−βξϰP−P1−βρϰP−P1C−βP1PI+βρ(ε−ϑC1)ϰε−ϑC1C−C1P−P1−βρ(μ+ϑP)ϰε−ϑC1C−C12+βρϰP−P1C1−βρϰP−P1C1=1−Θ(U1,P1)Θ(U,P1)γU1−γU+βI12−Θ(U1,P1)Θ(U,P1)+βI1Θ(U,P)Θ(U,P1)−βI1I1Θ(U,P)IΘ(U1,P1)−βϰP−P1ξ+ρC1−βρϰP−P1C−C1−βI1P1IPI1+βρϰC−C1P−P1−βρ(μ+ϑP)ϰε−ϑC1C−C12=1−Θ(U1,P1)Θ(U,P1)γU1−γU+βI13−Θ(U1,P1)Θ(U,P1)−I1Θ(U,P)IΘ(U1,P1)−P1IPI1+βI1Θ(U,P)Θ(U,P1)−PP1−βρμ+ϑPϰε−ϑC1C−C12=1−Θ(U1,P1)Θ(U,P1)γU1−γU+βI14−Θ(U1,P1)Θ(U,P1)−I1Θ(U,P)IΘ(U1,P1)−P1IPI1−PΘ(U,P1)P1Θ(U,P)+βI1Θ(U,P)Θ(U,P1)−PP1−1+PΘ(U,P1)P1Θ(U,P)−βρμ+ϑPϰε−ϑC1C−C12=γU11−UU11−Θ(U1,P1)Θ(U,P1)+βI14−Θ(U1,P1)Θ(U,P1)−I1Θ(U,P)IΘ(U1,P1)−P1IPI1−PΘ(U,P1)P1Θ(U,P)+βI1Θ(U,P)Θ(U,P1)−PP11−Θ(U,P1)Θ(U,P)−βρμ+ϑPϰε−ϑC1C−C12.

From Assumptions (A2) and (A4) we have

1−UU11−Θ(U1,P1)Θ(U,P1)≤0,Θ(U,P)Θ(U,P1)−PP11−Θ(U,P1)Θ(U,P)≤0.

Therefore, using inequality (17) we get that for all U, I, P, C > 0 we have dZ1dt ≤ 0 and dZ1dt = 0 if and only if U = U₁, I = I₁, P = P₁ and C = C₁. Applying LaSalle’s invariance principle, we obtain that if R0G > 1, then EP₁ is globally asymptotically stable.□

4 Numerical simulations

We conduct numerical simulations for model (18)-(21) with specific incidence rate function

Θ(U,P)=ωUP1+α1P+α2U.

Then we get following model with Beddington-DeAngelis functional response:

U˙=ϱ−γU−ωUP1+α1P+α2U, (31)

I˙=ωUP1+α1P+α2U−βI, (32)

P˙=ϰI−ξP−ρPC, (33)

C˙=εP−μC−ϑPC, (34)

where ω is a positive parameter, while α₁ and α₂ are non-negative parameters. We note that if α₁ = α₂ = 0, then we obtain a model with bilinear incidence, if α₁ ≠ 0 and α₂ = 0, then we get saturated incidence which given in model (8)-(11), and if α₁ = 0 and α₂ ≠ 0, then we obtain Holling type-II. We can easily see that Θ(U, P) is continuously differentiable function. Moreover, Θ(U, P) satisfying the following conditions:

We have

∂Θ(U,P)∂U=ωP+α1ωP2(1+α1P+α2U)2,∂Θ(U,P)∂P=ωU+α2ωU2(1+α1P+α2U)2,

then Θ(U, P) is continuously differentiable. Moreover, Θ(U, P) > 0, and Θ(0, P) = Θ(U, 0) = 0 for all U > 0 and P > 0. Thus (A1) is satisfied.

Since ∂Θ(U,P)∂U>0,∂Θ(U,P)∂P>0, and ∂Θ(U,0)∂P=ωU1+α2U>0 for all U > 0, then (A2) is satisfied.

We have

ddU∂Θ(U,0)∂P=ω(1+α2U)2>0 for all U≥0,

then (A3) is satisfied.

Finally we have

∂∂PΘ(U,P)P=−α1ωU(1+α1P+α2U)2<0, for all P≥0,

then (A4) is also satisfied.

The basic reproduction number of model (31)-(34) is given by

R0G=ϰωU0βξ(1+α2U0).

In the numerical simulations we fix the values of parameters ϱ = 10, ϰ = ξ = 3, ρ = 0.1, μ = γ = 0.01, β = 0.3, ε = 0.2 and vary ω and ϑ

Case(1): Effect of ω on the stability of equilibria.

For this case, we take α₁ = α₂ = 0 and ϑ = 0.01. We choose three different initial conditions as:

IC1: U(0) = 700, I(0) = 5, P(0) = 5, C(0) = 0.5,

IC2: U(0) = 400, I(0) = 10, P(0) = 10, C(0) = 1,

IC3: U(0) = 300, I(0) = 20, P(0) = 15, C(0) = 1.5.

We consider two values of the parameter ω as:

ω = 0.0001, then we compute R0G = 0.3333 < 1. Figure 1 shows that, for all IC1-IC3, the solution of the model tends to EP₀ = (1000, 0, 0, 0). It means that, EP₀ is globally asymptotically stable.

$Figure 1 Solution trajectories of system (31)-(34) in case α1 = α2 = 0.$
Figure 1
Solution trajectories of system (31)-(34) in case α₁ = α₂ = 0.
ω = 0.001, then we compute R0G = 3.3333 > 1. Figure 1, shows that the solutions of the model converge to the equilibrium EP₁ = (482.9, 17.23, 10.7, 18.29) for all IC1-IC3. Then, EP₁ is globally asymptotically stable.

Case(2): Effect of the saturation infection on the virus dynamics.

In this case, we take α₂ = 0 and ϑ = 0.01. We choose ω = 0.001, and α₁ varied. Moreover we consider the initial condition IC2. Figure 2 shows that as α₁ is increased, the concentrations of the uninfected target cells is increased, while the the concentration of infected cells, virus particles and B cells are decreased. We note that the parameter α₁ has no effect on the stability of equilibria

$Figure 2 Solution trajectories of system (31)-(34) for different values of α1 when α2 = 0.$

Figure 2

Solution trajectories of system (31)-(34) for different values of α₁ when α₂ = 0.

Case(3): Effect of Holling type-II.

For this case, we take α₁ = 0, ω = 0.001, and ϑ = 0.01 then Θ(U, P) represents the Holling type-II. Let us choose the initial condition IC2. We suggest different values of α₂ to see its effect on the model as we can see in Figure 3. Moreover, we have the following cases:

$Figure 3 Solution trajectories of system (31)-(34) for different values of α2 when α1 = 0.$

Figure 3

Solution trajectories of system (31)-(34) for different values of α₂ when α₁ = 0.

EP₁ is globally asymptotically stable when 0 ≤ α₂ < 0.0023,
EP₀ is globally asymptotically stable when α₂ > 0.0023.

This means that α₂ can play the role of controller which can be designed to stabilize the system around the infection-free equilibrium EP₀.

Case(4): Effect of the B cell impairment parameter ϑ.

In this case, we take α₁ = 0.01 and α₂ = 0.002. We choose ω = 0.001, and ϑ varied. Moreover, we consider the following initial condition

IC4: U(0) = 700, I(0) = 10, P(0) = 10, C(0) = 10.

Figure 4 shows that as ϑ is increased, the concentrations of infected cells and virus particles are increased, while the concentration of uninfected cells is decreased. We note that the parameter ϑ has no effect on the stability of equilibria.

$Figure 4 Solution trajectories of system (31)-(34) for different values of ϑ.$

Figure 4

Solution trajectories of system (31)-(34) for different values of ϑ.

Acknowledgement

This project was funded by the Deanship of Scientific Research (DSR), King Abdulaziz University, Jeddah, Saudi Arabia under grant no. KEP-MSc-33-130-40. The authors, therefore, acknowledge with thanks DSR technical and financial support.

References

[1] Nowak M., Bangham C., Population dynamics of immune responses to persistent viruses, Science, 1996, 272, 74-79.10.1126/science.272.5258.74Search in Google Scholar PubMed

[2] Nowak M., May R., Virus dynamics: mathematical principles of immunology and virology, Oxford university press, 2000.Search in Google Scholar

[3] Perelson A., Nelson P., Mathematical analysis of HIV-1 dynamics in vivo, SIAM Rev., 1999, 41, 3-44.10.1137/S0036144598335107Search in Google Scholar

[4] Callaway D., Perelson A., HIV-1 infection and low steady state viral loads, B. Math. Biol., 2002, 64, 29-64.10.1006/bulm.2001.0266Search in Google Scholar PubMed

[5] Roy P., Chatterjee A., Greenhalgh D., Khan Q., Long term dynamics in a mathematical model of HIV-1 infection with delay in different variants of the basic drug therapy model, Nonlinear Anal-Real., 2013, 14, 1621-1633.10.1016/j.nonrwa.2012.10.021Search in Google Scholar

[6] Leenheer P., Smith H., Virus dynamics: A global analysis, SIAM J. Appl. Math., 2003, 63, 1313-1327.10.1137/S0036139902406905Search in Google Scholar

[7] Song X., Neumann A., Global stability and periodic solution of the viral dynamics, J. Math. Anal. Appl., 2007, 329, 281-297.10.1016/j.jmaa.2006.06.064Search in Google Scholar

[8] Huang G., Takeuchi Y., Ma W., Lyapunov functionals for delay differential equations model of viral infections, SIAM J. Appl. Math., 2010, 70, 2693-2708.10.1137/090780821Search in Google Scholar

[9] Georgescu P., Hsieh Y., Global stability for a virus dynamics model with nonlinear incidence of infection and removal, SIAM J. Appl. Math., 2006, 67, 337-353.10.1137/060654876Search in Google Scholar

[10] Li M., Shu H., Global Dynamics of an in-host viral model with intracellular delay, B. Math. Biol., 2010, 72, 1492-1505.10.1007/s11538-010-9503-xSearch in Google Scholar PubMed

[11] Li M., Shu H., Impact of intracellular delays and target-cell dynamics on in vivo viral infections, SIAM J. Appl. Math., 2010, 70, 434-2448.10.1137/090779322Search in Google Scholar

[12] Hattaf K., Yousfi N., Tridane A., Mathematical analysis of a virus dynamics model with general incidence rate and cure rate, Nonlinear Anal-Real., 2012, 13, 1866-1872.10.1016/j.nonrwa.2011.12.015Search in Google Scholar

[13] Korobeinikov A., Global properties of basic virus dynamics models, B. Math. Biol., 2004, 66, 879-883.10.1016/j.bulm.2004.02.001Search in Google Scholar PubMed

[14] Elaiw A., Global properties of a class of HIV models, Nonlinear Anal-Real., 2010, 11, 2253-2263.10.1016/j.nonrwa.2009.07.001Search in Google Scholar

[15] Elaiw A., Global properties of a class of virus infection models with multitarget cells, Nonlinear Dynam., 2012, 69, 423-435.10.1007/s11071-011-0275-0Search in Google Scholar

[16] Elaiw A., Almuallem N., Global properties of delayed-HIV dynamics models with differential drug efficacy in cocirculating target cells, Appl. Math. Comput., 2015, 265, 1067-1089.10.1016/j.amc.2015.06.011Search in Google Scholar

[17] Elaiw A., Almuallem N., Global dynamics of delay-distributed HIV infection models with differential drug efficacy in cocirculating target cells, Math. Method. Appl. Sci., 2016, 39, 4-31.10.1002/mma.3453Search in Google Scholar

[18] Elaiw A., Hassanien I., Azoz S., Global stability of HIV infection models with intracellular delays, J. Korean Math. Soc., 2012, 49, 779-794.10.4134/JKMS.2012.49.4.779Search in Google Scholar

[19] Elaiw A., Azoz S., Global properties of a class of HIV infection models with Beddington-DeAngelis functional response, Math. Method. Appl. Sci., 2013, 36, 383-394.10.1002/mma.2596Search in Google Scholar

[20] Elaiw A., Global dynamics of an HIV infection model with two classes of target cells and distributed delays, Discrete Dyn. Nat. Soc., 2012, Article ID 253703.10.1155/2012/253703Search in Google Scholar

[21] Gibelli L., Elaiw A., Alghamdi M., Althiabi A., Heterogeneous population dynamics of active particles: Progression, mutations, and selection dynamics, Math. Method. Appl. Sci., 2017, 27, 617-640.10.1142/S0218202517500117Search in Google Scholar

[22] Bellomo N., Tao Y., Stabilization in a chemotaxis model for virus infection, Discrete Cont. Dyn-S., 2020, 13, 105-117.10.3934/dcdss.2020006Search in Google Scholar

[23] Lydyard P., Whelan A., Fanger M., Instant notes in immunology, Taylor & Francis e-Library, 2005.10.4324/9780203488287Search in Google Scholar

[24] Murase A., Sasaki T., Kajiwara T., Stability analysis of pathogen-immune interaction dynamics, J. Math. Biol., 2005, 51, 247-267.10.1007/s00285-005-0321-ySearch in Google Scholar PubMed

[25] Wang T., Hu Z., Liao F., Stability and Hopf bifurcation for a virus infection model with delayed humoral immunity response, J. Math. Anal. Appl., 2014, 411, 63-74.10.1016/j.jmaa.2013.09.035Search in Google Scholar

[26] Elaiw A., AlShamrani N., Global stability of humoral immunity virus dynamics models with nonlinear infection rate and removal, Nonlinear Anal-Real., 2015, 26, 161-190.10.1016/j.nonrwa.2015.05.007Search in Google Scholar

[27] Wang T., Hu Z., Liao F., Ma W., Global stability analysis for delayed virus infection model with general incidence rate and humoral immunity, Math. Comput. Simulat., 2013, 89, 13-22.10.1016/j.matcom.2013.03.004Search in Google Scholar

[28] Elaiw A., AlShamrani N., Stability of a general delay-distributed virus dynamics model with multi-staged infected progression and immune response, Math. Method. Appl. Sci., 2017, 40, 699-719.10.1002/mma.4002Search in Google Scholar

[29] Elaiw A., Alshaikh M., Stability analysis of a general discrete-time pathogen infection model with humoral immunity, J. Differ. Equ. Appl., 2019, 10.1080/10236198.2019.1662411.Search in Google Scholar

[30] Elaiw A., AlShamrani N., Stability of an adaptive immunity pathogen dynamics model with latency and multiple delays, Math. Method. Appl. Sci., 2018, 36, 125-142.10.1002/mma.5182Search in Google Scholar

[31] Elaiw A., AlShamrani N., Stability of a general adaptive immunity virus dynamics model with multi-stages of infected cells and two routes of infection, Math. Method. Appl. Sci., 2019, 10.1002/mma.5923.Search in Google Scholar

[32] Elaiw A., Elnahary E., Analysis of general humoral immunity HIV dynamics model with HAART and distributed delays, Mathematics, 2019, 7, Article number: 157.10.3390/math7020157Search in Google Scholar

[33] Elaiw A., Almatrafi A., Hobiny A., Hattaf K., Global properties of a general latent pathogen dynamics model with delayed pathogenic and cellular infections, Discrete Dyn. Nat. Soc., 2019, Article ID 9585497.10.1155/2019/9585497Search in Google Scholar

[34] Lin J., Xu R., Tian X, Threshold dynamics of an HIV-1 virus model with both virus-to-cell and cell-to-cell transmissions, intracellular delay, and humoral immunity, Appl. Math. Comput., 2017, 315, 516-530.10.1016/j.amc.2017.08.004Search in Google Scholar

[35] Pan S., Chakrabarty S., Threshold dynamics of HCV model with cell-to-cell transmission and a non-cytolytic cure in the presence of humoral immunity, Commun. Nonlinear Sci., 2018, 61, 180-197.10.1016/j.cnsns.2018.02.010Search in Google Scholar

[36] Miao H., Teng Z., Kang C., Muhammadhaji A., Stability analysis of a virus infection model with humoral immunity response and two time delays, Math. Method. Appl. Sci., 2016, 39, 3434-3449.10.1002/mma.3790Search in Google Scholar

[37] Hu Z., Zhang J., Wang H., Ma W., Liao F., Dynamics analysis of a delayed viral infection model with logistic growth and immune impairment, Appl. Math. Model., 2014, 38, 524-534.10.1016/j.apm.2013.06.041Search in Google Scholar

[38] Regoes R., Wodarz D., Nowak M., Virus dynamics: the effect to target cell limitation and immune responses on virus evolution, J. Theor. Biol., 1998, 191, 451-462.10.1006/jtbi.1997.0617Search in Google Scholar PubMed

[39] Wang S., Song X., Ge Z., Dynamics analysis of a delayed viral infection model with immune impairment, Appl. Math. Model., 2011, 35, 4877-4885.10.1016/j.apm.2011.03.043Search in Google Scholar

[40] Elaiw A., Raezah A., Azoz S., Stability of delayed HIV dynamics models with two latent reservoirs and immune impairment, Adv. Differ. Equ-Ny., 2018, Article number: 414.10.1186/s13662-018-1869-3Search in Google Scholar

[41] Milito A., B Lymphocyte Dysfunctions in HIV Infection, Curr. HIV Res., 2004, 2, 11-21.10.2174/1570162043485068Search in Google Scholar PubMed

[42] Amu S., Ruffin N., Rethi B., Chiodi F., Impairment of B-cell functions during HIV-1 infection, AIDS, 2013, 27, 2323-2334.10.1097/QAD.0b013e328361a427Search in Google Scholar PubMed

[43] Chiodi F., Scarlatti G., Editorial: HIV-Induced damage of B cells and production of HIV neutralizing antibodies, Front. immunol., 2018, 9, Article number: 297.10.3389/fimmu.2018.00297Search in Google Scholar PubMed PubMed Central

[44] Miao H., Abdurahman X., Teng Z., Zhang L., Dynamical analysis of a delayed reaction-diffusion virus infection model with logistic growth and humoral immune impairment, Chaos Soliton. Fract., 2018, 110, 280-291.10.1016/j.chaos.2018.03.006Search in Google Scholar

[45] Wang W., Ma W., Feng Z., Dynamics of reaction-diffusion equations for modeling CD4⁺ T cells decline with general infection mechanism and distinct dispersal rates, Nonlinear Anal-Real., 2020, 51, Article number: 102976.10.1016/j.nonrwa.2019.102976Search in Google Scholar

[46] Wang W., Ma W., Feng Z., Complex dynamics of a time periodic nonlocal and time-delayed model of reaction-diffusion equations for modeling CD4⁺ T cells decline, J. Comput. Appl. Math., 2020, 367, Article number: 112430.10.1016/j.cam.2019.112430Search in Google Scholar

[47] Elaiw A., Abukwaik R., Alzahrani E., Global properties of a cell mediated immunity in HIV infection model with two classes of target cells and distributed delays, Int. J. Biomath., 2014, 7, Article ID 1450055.10.1142/S1793524514500557Search in Google Scholar

[48] Elaiw A., AlShamrani N., Global properties of nonlinear humoral immunity viral infection models, Int. J. Biomath., 2015, 8, Article ID 1550058.10.1142/S1793524515500588Search in Google Scholar

[49] Elaiw A., Alshehaiween S., Hobiny A., Global properties of delay-distributed HIV dynamics model including impairment of B-cell functions, Mathematics, 2019, 7, Article number: 837.10.3390/math7090837Search in Google Scholar

[50] Elaiw A., Alshaikh M., Stability of discrete-time HIV dynamics models with three categories of infected CD4⁺ T-cells, Adv. Differ. Equ-Ny., 2019, Article number: 407.10.1186/s13662-019-2338-3Search in Google Scholar

[51] Elaiw A., Raezah A., Stability of general virus dynamics models with both cellular and viral infections and delays, Math. Method. Appl. Sci., 2017, 40, 5863-5880.10.1002/mma.4436Search in Google Scholar

[52] Hobiny A., Elaiw A., Almatrafi A., Stability of delayed pathogen dynamics models with latency and two routes of infection, Adv. Differ. Equ-Ny., 2018, Article number: 276.10.1186/s13662-018-1720-xSearch in Google Scholar

[53] Elaiw A., Elnahary E., Raezah A., Effect of cellular reservoirs and delays on the global dynamics of HIV, Adv. Differ. Equ-Ny., 2018, Article number: 85.10.1186/s13662-018-1523-0Search in Google Scholar

[54] Elaiw A., AlAgha A., Global dynamics of reaction-diffusion oncolytic M1 virotherapy with immune response, Appl. Math. Comput., 2019 (accepted).10.1016/j.amc.2019.124758Search in Google Scholar

[55] Gao N., Song Y., Wang X., Liu J., Dynamics of a stochastic SIS epidemic model with nonlinear incidence rates, Adv. Differ. Equ-Ny., 2019, article number: 41.10.1186/s13662-019-1980-0Search in Google Scholar

[56] Song Y., Miao A., Zhang T., Wang X., Liu J., Extinction and persistence of a stochastic SIRS epidemic model with saturated incidence rate and transfer from infectious to susceptible, Adv. Differ. Equ-Ny., 2018, Article number: 293.10.1186/s13662-018-1759-8Search in Google Scholar

[57] Wang Y., Cao J., Alofi A., AL-Mazrooei A., Elaiw A., Revisiting node-based SIR models in complex networks with degree correlations, Physica A., 2015, 437, 75-88.10.1016/j.physa.2015.05.103Search in Google Scholar

[58] Hale J., Verduyn Lunel S., Introduction to functional differential equations, Springer Science & Business Media, 2013.10.1007/978-1-4612-4342-7Search in Google Scholar

Received: 2019-04-07

Accepted: 2019-10-04

Published Online: 2019-12-13

This work is licensed under the Creative Commons Attribution 4.0 Public License.

Articles in the same Issue

https://doi.org/10.1515/math-2019-0113

Keywords for this article

virus dynamics; global stability; Lyapunov function; B-cell impairment

Creative Commons

BY 4.0