An SEIR model with modified saturated incidence rate and Holling type II treatment function

Shilpa Umdekar; Praveen Kumar Sharma; Shivram Sharma

doi:10.1515/cmb-2022-0146

Article Open Access

An SEIR model with modified saturated incidence rate and Holling type II treatment function

Shilpa Umdekar , Praveen Kumar Sharma and Shivram Sharma

Published/Copyright: May 8, 2023

Published by

Become an author with De Gruyter Brill

Submit Manuscript Author Information

From the journal Computational and Mathematical Biophysics Volume 11 Issue 1

Abstract

In this article, the behavior of an susceptible exposed infected recovered (SEIR) epidemic model with nonlinear incidence rate and Holling type II treatment function is presented and analyzed. Reproduction number of the model is calculated. Equilibrium points are determined. Disease-free equilibrium exists when R0 is below 1. Behavior of disease-free equilibrium is examined at R0 = 1. Endemic equilibrium exists when R0 crosses 1. Stability of both equilibrium points is investigated locally and globally. Simulation is provided to support the result.

Keywords: SEIR, basic reproduction number; Holling type II function; Lyapunov function

MSC 2010: 34D23; 93A30; 93D20

1 Introduction

Infectious diseases have been a part of human life for a long time already. Evidence tells us that epidemics often end up causing mass deaths. It was after the increase in healthcare, for a certain period of time, the health burden diminished of infectious diseases. However, in recent years, it has emerged that the challenge still exists, especially, in our rapidly changing world since every nation has limited resources to treat the infected. Emerging diseases pose a continuing threat, for example, human immunodeficiency virus in the twentieth century, acute encephalitis syndrome, malaria, cholera, and more recently COVID-19-coronavirus, causing mortality that has proven the necessity of having optimal resources to control an epidemic. Various mathematical models for infectious diseases proposed by many authors (see [1–9,16,18]). To figure out this problem, many treatment functions have been proposed by various researchers [13,14].

Zhang and Xianning [17] introduced the saturated treatment function for the better analysis of real system through the epidemic model. This function is widely known as Holling type II treatment function,

h ( I ) = aI 1 + bI ,

where a > 0 is a cure rate and b ≥ 0 measures the magnitude of the consequence of the infected person being held for treatment. The specialty of this function is that it is defined as continually differentiable and characterizes the situation of limited medical resources.

Whenever we talk about epidemic or pandemic, the incidence rate of the disease in a population is first to be discussed. Previously, bilinear incidence rate β SI was being used more often to assess the new infections per unit of time in a population of susceptible. But while pondering about incidence rate, some significant factors must always be considered, like the inhibitory effects caused by the sudden spread of disease in a population. Pathak et al. [11] proposed such saturated incidence rate. This modified incidence rate is defined using sociological and psychological parameters α 1 and α 2 , respectively, given by β I 1 + α 1 S + α 2 I . The introduction of these parameters by authors is a modified saturation effect.

In this article, we propose an susceptible exposed infected recovered (SEIR) model with modified saturated incidence rate and Holling type II treatment function. Then, for the model, we calculate R0, find out the equilibrium point, and discuss the local and global stability of disease free equilibrium (DFE) and endemic equilibrium (EE). We check the stability of equilibrium at R0 = 1 . In the next section, we provide simulation to assist theoretical results. In the last section, we discuss the aspects of the proposed model.

2 Proposed mathematical model

The epidemic model we propose and study in this work is an SEIR model with saturated incidence rate and Holling type II treatment function using nonlinear ordinary differential equations (Figure 1).

ⅆ S ⅆ t = A − β SI 1 + α 1 s + α 2 I − μ S ,

ⅆ E ⅆ t = β SI 1 + α 1 S + α 2 I − ( γ + μ ) E ,

ⅆ I ⅆ t = γ E − ( δ + μ ) I − aI 1 + bI ,

(2.1) ⅆ R ⅆ t = δ I − μ R + aI 1 + bI ,

where S ( t ) + E ( t ) + I ( t ) + R ( t ) = N ( t ) .

Figure 1

Transfer diagram of the presented model.

In equation (2.1), S ( t ) , E ( t ) , I ( t ) , and R ( t ) denote the susceptible, exposed, infected, and recovered number of individuals at time t . Other parameters used in the system are infection rate β , recruitment rate A , natural death rate μ , natural recovery rate δ , progression rate γ , sociological parameter α 1 > 0 , psychological parameter α 2 > 0 , cure rate a > O , and magnitude of the consequence of delaying treatment for an infected individual by b ≥ 0 .

Since the first three equations are free from R , we can reduce system (2.1) into:

ⅆ S ⅆ t = A − β SI 1 + α 1 s + α 2 I − μ S ,

(2.2) ⅆ E ⅆ t = β SI 1 + α 1 S + α 2 I − ( γ + μ ) E ,

ⅆ I ⅆ t = γ E − ( δ + μ ) I − aI 1 + bI .

Lemma 2.1

The set Ω = ( S , E , I , R ) : S + E + I + R ≤ A μ , S ≥ 0 , E ≥ 0 , I ≥ 0 , R ≥ 0 is a positively invariant system of system (2.1).

From (2.2), we obtain:

ⅆ ⅆ t ( S + E + I + R ) = A − μ ( S + E + I + R ) − δ I − aI 1 + bI ≤ A − μ ( S + E + I + R ) ≤ A − μ N .

Solving the aforementioned equation and applying lim t → ∞ N ( t ) , we obtain N ( t ) < A μ .

Thus, feasible region invariant with system (2.2) is given by Lemma 2.1.

3 Main results

Reproduction number R0 is calculated using the method by van den Driessche and Watmough [12]. The DFE of system (2.2) is at D 0 = A μ , 0,0 . For calculating R0, we focus on infected compartments E and I only; thus, the system furthermore reduces to:

ⅆ E ⅆ t = β SI 1 + α 1 S + α 2 I − ( γ + μ ) E ,

(3.1) ⅆ I ⅆ t = γ E − ( δ + μ ) I − aI 1 + bI .

Let x = ( E , I ) T ,

F = β SI 1 + α 1 S + α 2 I 0 and V = ( γ + μ ) E − γ E + ( δ + μ ) I + aI 1 + bI .

Jacobian matrix of F and V at D 0 = A μ , 0,0 is given by:

F 1 = 0 β A μ + α 1 A 0 0 , V 1 = ( γ + μ ) 0 − γ ( δ + μ ) + a .

The next-generation matrix is given by:

F 1 V 1 − 1 = γ β A ( μ + α 1 A ) ( γ + μ ) ( δ + μ + a ) β A ( μ + α 1 A ) ( δ + μ + a ) 0 0 .

Reproduction number R0 is the dominant eigen value of F 1 V 1 − 1

(3.2) R 0 = γ β A ( μ + α 1 A ) ( γ + μ ) ( δ + μ + a ) .

3.1 Equilibriums of the system

Setting all rates to zero in system (2.2), that is, setting d S d t = 0 , d I d t = 0 , d E d t = 0 , we obtain:

A − β SI 1 + α 1 s + α 2 I − μ S = 0 ,

(3.3) β SI 1 + α 1 S + α 2 I − ( γ + μ ) E = 0 ,

γ E − ( δ + μ ) I − aI 1 + bI = 0 .

From (3.3), we obtain:

I = 0 or γ β S ( 1 + α 1 S + α 2 I ) ( γ + μ ) − ( δ + μ ) + a 1 + bI = 0 .

3.2 Disease-free equilibrium (DFE)

Taking I = 0 in the remaining equations of system (3.3), we obtain the values of S and E as follows:

S = A μ ,

E = 0 .

Thus, we determined the DFE:

(3.4) D 0 = A μ , 0 , 0 .

3.3 Endemic equilibrium (EE)

S * = ( γ + μ ) ( 1 + α 2 I * ) { ( 1 + b I * ) ( δ + μ ) + a } β γ ( 1 + b I * ) − ( γ + μ ) α 1 { ( 1 + b I * ) ( δ + μ ) + a } ,

E * = A − μ S * ( γ + μ ) .

Substituting the value of S * in the first equation of system (3.3) yields cubic equation of I * :

(3.5) ψ 1 I * 3 + ψ 2 I * 2 + ψ 3 I * + ψ 4 = 0 ,

where

ψ 1 = β γ b 2 ( γ + μ ) ( δ + μ ) − α 1 b 2 ( γ + μ ) 2 ( δ + μ ) 2 ,

ψ 2 = 2 b β γ ( γ + μ ) ( δ + μ ) + ab β γ ( γ + μ ) − 2 α 1 b ( γ + μ ) 2 ( δ + μ ) 2 − 2 ab α 1 ( γ + μ ) 2 ( δ + μ ) + α 2 b μ ( γ + μ ) ( δ + μ ) − A β b 2 γ 2 + A α 1 γ ( γ + μ ) ( δ + μ ) b 2 ,

ψ 3 = β γ ( γ + μ ) ( δ + μ ) + a β γ ( γ + μ ) − α 1 ( γ + μ ) 2 ( δ + μ ) 2 − 2 a α 1 ( γ + μ ) 2 ( δ + μ ) − a 2 α 1 ( γ + μ ) 2 ( δ + μ ) + α 2 μ ( γ + μ ) ( δ + μ ) + α 2 μ ( γ + μ ) − 2 Ab β + 2 Ab α 1 γ ( γ + μ ) ( δ + μ ) + Aab α 1 γ ( γ + μ ) ,

ψ 4 = a + μ ( γ + μ ) − A β γ 2 + A α 1 γ ( γ + u ) ( δ + μ ) + Aa α 1 γ ( γ + μ ) .

Now, unique positive real root of equation (3.5) exists [15] if:

i. ψ 1 > 0 , ψ 2 > 0 , ψ 3 > 0 and ψ 4 < 0 ,

ii. ψ 1 > 0 , ψ 2 > 0 , ψ 3 < 0 and ψ 4 < 0 ,

iii. ψ 1 > 0 , ψ 2 < 0 , ψ 3 < 0 and ψ 4 < 0 ,

where ψ 1 > 0 . We can estimate S * and E * after we obtain the value of I * . Thus, there exists a unique EE if the aforementioned inequalities are satisfied under the following conditions.

Let us take m = ( γ + μ ) , n = ( δ + μ ) ,

ψ 2 > 0 , when 2 b β γ mn + ab β γ m + α 2 b μ mn + A α 1 γ mn b 2 > 2 α 1 b m 2 n 2 + 2 ab α 1 m 2 n + A β b 2 γ 2 < 0 , when 2 b β γ mn + ab β γ m + α 2 b μ mn + A α 1 γ mn b 2 < 2 α 1 b m 2 n 2 + 2 ab α 1 m 2 n + A β b 2 γ 2 ,

ψ 3 > 0 , when β γ mn + a β γ m + α 2 μ mn + α 2 μ m + 2 Ab α 1 γ mn + Aab α 1 γ m > α 1 m 2 n 2 + 2 a α 1 m 2 n + a 2 α 1 m 2 n + 2 Ab β < 0 , when β γ mn + a β γ m + α 2 μ mn + α 2 μ m + 2 Ab α 1 γ mn + Aab α 1 γ m < α 1 m 2 n 2 + 2 a α 1 m 2 n + a 2 α 1 m 2 n + 2 Ab β ,

ψ 4 > 0 , when a + μ m + A α 1 γ m ( n + a ) > A β γ 2 < 0 , when a + μ m + A α 1 γ m ( n + a ) < A β γ 2 .

Thus, EE D * = ( S * , E * , I * ) exists when the aforementioned conditions are satisfied.

3.4 Stability analysis of equilibria

Theorem 3.1

For R0 < 1 , DFE D 0 = A μ , 0 , 0 is locally asymptotically stable and unstable for R0 > 1 .

For this, we construct a Jacobian matrix of system (2.2) at DFE as follows:

(3.6) J ( X 0 ) = − μ 0 − β A μ + α 1 A 0 − ( γ + μ ) β A μ + α 1 A 0 γ − ( δ + μ ) − a .

Thus, | J ( X 0 ) − λ I | is given by:

− μ − λ 0 − β A μ + α 1 A 0 − ( γ + μ ) − λ β A μ + α 1 A 0 γ − ( μ + δ ) − a − λ = 0 .

That implies

(3.7) ( λ + μ ) { λ 2 + ( 2 μ + δ + γ + a ) λ + ( γ + μ ) ( μ + δ + a ) − ( 1 − R 0 ) } = 0 .

First eigen value of matrix is μ < 0 , and other two eigen values are zeros of the equation:

λ 2 + p λ + q = 0 ,

where p = 2 μ + γ + δ + a > 0 ,

q = ( μ + r ) ( μ + δ + a ) ( 1 − R 0 ) .

If R0 ≤ 1 , then q > 0 , so both zeros of equation are negative. So, all the eigen values of the Jacobian matrix at DFE are negative. Therefore, by the Routh–Hurwitz criterion, Theorem 3.1 is proven.

Now, we analyze the global stability of DFE. For this, define:

(3.8) R 0 a = γ β A ( μ + α 1 A ) ( γ + μ ) ( δ + μ + ( a / 1 + b ( A μ ) ) ) .

Theorem 3.2

DFE, D 0 , is globally asymptotically stable if R 0 a < 1 , implying R0 < 1 .

Proof

From the first equation of system (2.1), we have d S d t ≤ A − μ S . Solution of d x d t = A − μ x is a maximal solution of S ( t ) and x → A / μ as t → ∞ . By the comparison theorem, S ( t ) ≤ A / μ and from set Ω in Lemma 2.1, we obtain I ( t ) ≤ A / μ .

Let us define the Lyapunov function:

L = γ E + ( γ + μ ) I .

Therefore,

d L d t = γ β S 1 + α 1 S + α 2 I − ( γ + μ ) δ + μ + a 1 + bI I

≤ γ β S 1 + α 1 S + α 2 I − ( γ + μ ) δ + μ + a 1 + bI I

≤ γ β A μ + α 1 A − ( γ + μ ) δ + μ + a 1 + b ( A / μ ) I ≤ 0 ,

and d L d t = 0 iff I = 0 .

{ ( S , E , I ) ∈ Ω , d L d t = 0 } is the singleton set. Thus, by the Lasalle–Lyapunov theorem, DFE is asymptotically stable globally.□

Theorem 3.3

DFE, D 0 , is unstable at R0 = 1 ; thus, a positive equilibrium is found when R0 crosses 1 .

To check the stability of DFE, we make use of center manifold theory [12]. To do this, let us suppose:

x 1 = S , x 2 = E , and x 3 = I . Define X = ( x 1 , x 2 , x 3 ) T in such a way that system (2.2) can be written as ⅆ X ⅆ t = F ( x ) , where F = ( f 1 , f 2 , f 3 ) . Therefore,

ⅆ x 1 ⅆ t = A − β x 1 x 3 1 + α 1 x 1 + α 2 x 3 − μ x 1 = f 1 ,

(3.9) ⅆ x 2 ⅆ t = β x 1 x 3 1 + α 1 x 1 + α 2 x 3 − ( γ + μ ) x 2 = f 2 ,

(3.10) ⅆ x 3 ⅆ t = − γ x 2 − ( δ + μ ) x 3 − a x 3 1 + b x 3 = f 3 .

Now, at R0 = 1 , β = β * = ( μ + α 1 A ) ( γ + μ ) ( δ + μ + a ) γ A .

Let J * be the Jacobian matrix at R0 = 1 and β = β * . Then,

J * = − μ 0 − β * A μ + α 1 A 0 − ( γ + μ ) β * A μ + α 1 A 0 γ − ( δ + μ ) − a . Let w = [ w 1 , w 2 , w 3 ] and u = [ u 1 , u 2 , u 3 ] T be the left and right eigen vectors of J * correlating with zero eigen value. Then,

[ w 1 , w 2 , w 3 ] = 0 , γ ( γ + μ ) , 1 and [ u 1 , u 2 , u 3 ] T = − β * A μ + α 1 A , − β * A ( γ + μ ) ( μ + α 1 A ) , 1 .

Now, to obtain the bifurcation coefficients a 1 and a 2 , the nonzero partial derivatives of system (3.9) calculated at D 0 are as follows:

∂ 2 f 1 ∂ x 1 ∂ x 3 = ∂ 2 f 1 ∂ x 3 ∂ x 1 = ∂ 2 f 1 ∂ x 1 2 = − μ 2 β * ( μ + α 1 A ) 2 , ∂ 2 f 1 ∂ x 3 2 = 2 α 2 μ β * A ( u + α 1 A ) 2 , ∂ 2 f 2 ∂ x 3 2 = − 2 α 2 μ β * A ( u + α 1 A ) 2 , ∂ 2 f 2 ∂ x 1 ∂ x 3 = ∂ 2 f 2 ∂ x 3 ∂ x 1 = μ 2 β * ( μ + α 1 A ) 2 ,

∂ 2 f 3 ∂ x 3 2 = b 2 ,

a 1 = ∑ k , i , j = 1 3 w k u i u j ∂ 2 f k ∂ x i ∂ x j D 0 ,

a 1 = − 2 γ μ β * A ( γ + μ ) ( μ + α 1 A ) 2 β * ( μ + α 1 A ) + α 2 + b 2 < 0 .

(3.11) If a 1 satisfies b 2 < 2 γ μ β * A ( γ + μ ) ( μ + α 1 A ) 2 β * ( μ + α 1 A ) + α 2 ,

∂ 2 f 1 ∂ x 3 ∂ β * = − A ( μ + α 1 A ) , ∂ 2 f 2 ∂ x 3 ∂ β * = A ( μ + α 1 A ) ,

a 2 = ∑ k , i = 1 3 w k u i ∂ 2 f k ∂ x i ∂ β * D 0 ,

a 2 = γ A ( γ + μ ) ( μ + α 1 A ) = ( δ + μ + a ) β * > 0 .

Since a 1 < 0 when (3.11) is satisfied and a 2 > 0 , D 0 becomes unstable and a positive equilibrium is found when R0 crosses 1.

3.5 Endemic equilibrium (EE)

Theorem 3.4

EE, D * = ( S * , E * , I * ) , is locally asymptotically stable when R0 > 1

(3.12) iff μ + δ + a ( 1 + b I * ) 2 < β S * ( 1 + α 1 S * + α 2 I * ) 2 .

Jacobian matrix at EE:

J ( D * ) = − β I * ( 1 + α 2 I * ) ( 1 + α 1 S * + α 2 I * ) 2 − μ 0 − β S * ( 1 + α 1 S * ) ( 1 + α 1 S * + α 2 I * ) 2 β I * ( 1 + α 2 I * ) ( 1 + α 1 S * + α 2 I * ) 2 − ( γ + μ ) β S * ( 1 + α 1 S * ) ( 1 + α 1 S * + α 2 I * ) 2 0 γ − ( δ + μ ) − a ( 1 + b I * ) 2 .

Characteristic equation of J ( D * ) is given by:

κ 3 + s 1 κ 2 + s 2 κ + s 3 = 0 ,

where

s 1 = 3 μ + δ + γ + a ( 1 + b I * ) 2 + β I * ( 1 + α 2 I * ) ( 1 + α 1 S * + α 2 I * ) 2 ,

s 2 = 2 μ + δ + γ + a ( 1 + b I * ) 2 μ + β I * ( 1 + α 2 I * ) ( 1 + α 1 S * + α 2 I * ) 2 + μ μ + δ + a ( 1 + b I * ) 2 + γ μ + δ + a ( 1 + b I * ) 2 − β S * ( 1 + α 1 S * + α 2 I * ) 2 ,

s 3 = ( γ + μ ) μ + δ + a ( 1 + b I * ) 2 β I * ( 1 + α 2 I * ) ( 1 + α 1 S * + α 2 I * ) 2 + μ μ μ + δ + a ( 1 + b I * ) 2 + γ μ + δ + a ( 1 + b I * ) 2 − β S * ( 1 + α 1 S * + α 2 I * ) 2 .

Using the Routh–Hurwitz criterion, clearly, s 1 > 0 , also s 2 > 0 , s 3 > 0 , and ∆ = s 1 s 2 − s 3 > 0 when inequality (3.12) holds. Thus, D * = ( S * , E * , I * ) is locally asymptotically stable when R0 > 1 .

Theorem 3.5

EE, D * = ( S * , E * , I * ) , is globally asymptotically stable when R0 > 1 .

For system (2.2), Jacobian matrix is as follows:

J = − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − μ 0 − β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − ( γ + μ ) β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 0 γ − ( δ + μ ) − a ( 1 + bI ) 2 .

Second compound matrix is as follows:

J σ = − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − γ β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 γ − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − δ − a ( 1 + bI ) 2 0 0 β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − ( δ + γ + 2 μ ) − a ( 1 + bI ) 2 .

Choosing

K = 1 0 0 0 E / I 0 0 0 E / I implies K − 1 = 1 0 0 0 I / E 0 0 0 I / E .

Thus, K l = 0 0 0 0 E ′ I − E I ′ I 2 0 0 0 E ′ I − E I ′ I 2 .

Therefore, K l K − 1 = 0 0 0 0 E ′ E − I ′ I 0 0 0 E ′ E − I ′ I .

And KJ σ K − 1 = − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − γ β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 I E β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 I E γ E I − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − δ − a ( 1 + bI ) 2 0 0 β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − ( δ + γ + 2 μ ) − a ( 1 + bI ) 2 .

Now,

B = K l K − 1 + KJ σ K − 1 implies

B = − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − γ β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 I E β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 I E γ E I E ′ E − I ′ I − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − δ − a ( 1 + bI ) 2 0 0 β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 E ′ E − I ′ I − ( δ + γ + 2 μ ) − a ( 1 + bI ) 2 .

In block form:

B = B 11 B 12 B 21 B 22 ,

where

B 11 = − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − γ , B 12 = β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 I E β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 I E , B 21 = γ E I 0 ,

B 22 = E ′ E − I ′ I − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − δ − a ( 1 + bI ) 2 0 β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 E ′ E − I ′ I − ( δ + γ + 2 μ ) − a ( 1 + bI ) 2 .

If vectors of R 3 are given by ( p , q , r ) . Then, norm in R 3 can be selected as | p , q , r | = max ⁡ { | p | , | q | , |r| } . Let F ( B ) ≤ sup { j 1 , j 2 } , where F denotes Lozinski measure [10] and j 1 = F ( B 11 ) + | B 12 | , j 2 = F ( B 22 ) + | B 21 | .

| B 12 | = β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 I E , | B 21 | = γ E I ,

F ( B 11 ) = − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − γ ,

F ( B 22 ) = E ′ E − I ′ I − ( δ + 2 μ ) − a ( 1 + bI ) 2 ,

j 1 = − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − 2 μ − γ + β S ( 1 + α 1 S ) ( 1 + α 1 S + α 2 I ) 2 I E ,

j 2 = γ E I + E ′ E − I ′ I − ( δ + 2 μ ) − a ( 1 + bI ) 2 .

Second and third equations of system (2.2) can be rewritten as:

(3.13) E ′ E + γ + μ = β SI ( 1 + α 1 s + α 2 I ) E ,

(3.14) I ′ I + δ + μ + a 1 + bI = γ E I .

Substitution of (3.13) and (3.14), respectively, in j 1 and j 2 results in:

j 1 = E ′ E − μ − β I ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 − β S ( 1 + α 2 I ) ( 1 + α 1 S + α 2 I ) 2 I E ≤ E ′ E − μ ,

j 2 = E ′ E − μ + a ( 1 + bI ) 2 ≤ E ′ E − μ ,

iff a ( 1 + bI ) 2 < μ .

Thus, F ( B ) ≤ sup { j 1 , j 2 } ≤ E ′ E − μ .

And so,

y 2 ̅ = 1 t ∫ 0 t F ( B ) ds ≤ ln E ( t ) E ( t ′ ) + 1 t ∫ 0 t ′ F ( B ) d s − μ ,

implies, y 2 ̅ = lim t → ∞ sup sup 1 t ∫ 0 t F ( B ) d s < − μ < 0

implies y 2 ̅ < 0 .

Thus, D * = ( S * , E * , I * ) is globally asymptotically stable.

4 Example

For understanding the model better, we provide an example through the simulation of data. First, we consider the scenario when R0 < 1 . For this, we take A = 7 , μ = 0.2 , β = 0.5 , δ = 0.12 , γ = 0.03 , α 1 = 0.3 , α 2 = 0.3 , a = 1 , and b = 1 when R0 = 0.1504 . The initial values of all the four compartments in population are S ( 0 ) = 10 , E ( 0 ) = 5 , I ( 0 ) = 5 , and R ( 0 ) = 5 . Then, we obtain D 0 = ( S 0 , I 0 , E 0 ) = ( 35 , 0,0 ) . The following results can be seen in Figure 2. In this case, it can be clearly seen that DFE D 0 is approached and sustains as per Theorems 3.1 and 3.2. A whole lot of population is susceptible in less time when an infectious disease emerges but due to lower reproduction number, the average number of exposed and infected population remains low and approaches DFE.

$Figure 2 Behavior of S , E , I , and R S,E,I,{\rm{and}}{R} when R0 < 1 1 .$

Figure 2

Behavior of S , E , I , and R when R0 < 1 .

For the case, when R0 > 1 , A = 7 , μ = 0 . 2 , β = 1 , δ = 0 . 12 , γ = 0 . 3 , α 1 = 0 . 3 , α 2 = 0 . 3 , a = 0 . 3 , b = 0 . 4 results in R0 = 2 . 9453 with initial values S ( 0 ) = 10 , E ( 0 ) = 5 , I ( 0 ) = 5 , R ( 0 ) = 5 , which results in EE D * = ( S * , E * , I * , R * ) = ( 3.14 , 12.74 , 10.07 , 9.05 ) when inequality (3.12) is satisfied. Susceptible population decreases over time as a result of increased exposed and infected population, which also achieves a constant rate as they approach EE, which can be seen in Figure 3.

$Figure 3 Behavior of S , E , I , and R S,E,I,{\rm{and}}{R} when R0 < 1 \lt 1 .$

Figure 3

Behavior of S , E , I , and R when R0 < 1 .

In this article, we have called out Holling type II function as a better strategy to read out treatment of infected population since there is a significant difference. In the same time frame, while using the function, the infected number of individuals is less as compared to not using the function (Figure 4). As cure rate increases, we see a significant drop in infected population. Thus, in order to eradicate a disease, higher cure rate is absolutely essential (Figure 5). The longer it takes to treat the infected, the higher the number of infected; as a result of slow recovery, infectious disease spreads more widely (Figure 6).

Figure 4

Effect of the introduced treatment on infected.

Figure 5

Effect of cure rate.

Figure 6

Effect of delay in treatment.

We also discuss the difference in results with Holling type II function and nonlinear incidence rate, which is used in this article for all the compartments with the classic SEIR model in Figure 7. We can observe that with the introduction of the aforementioned functions when R0 > 1 , there is a significant increase in susceptible and recovered population, whereas there is a noticeable drop in exposed and infected. Thus, taking psychological and sociological effect into consideration helps us to curtail the infection by reducing the exposed population. Also, timely treatment reduces the infected individuals and joins the healthy recovered population.

$Figure 7 Comparison of the introduced model with classic SEIR {\rm{SEIR}} model.$

Figure 7

Comparison of the introduced model with classic SEIR model.

5 Discussion

In this work, we investigated the SEIR model with saturated incidence rate and Holling type II treatment rate. With the aid of the expression, we derived in (3.8), we measured the R0 value to determine the probability of epidemic and checked the behavior of population in the same time frame using graphs with the help of MATLAB. Note that the reproduction number R0 relies only on social behavior parameter α 1 related to susceptible individual, not on psychological impact α 2 of infected in population. This indicates that the reduction of social awareness in infected and rise in protective measure in susceptible are likely to have an effect on transmission of disease in the population. In order to grasp the psychological effect of the model, we make some further estimates by checking the existence and stability of EE. R0 also relies on the cure rate a of infected individuals, implying higher the value of a , lesser the R0 and lesser the amount of infected. So, cure rate should be higher to control epidemic. Sooner the infected people receive appropriate medical help, the quicker they heal and lot quicker we can control the epidemic. In this model, R0 does not depend on the magnitude of the consequence of the infected individual being delayed for treatment b ; to understand this impact, we did further assessment through EE and found that more delay in treatment of the infected will lead to slow recovery and they remain infected for longer time period causing the breakout of disease. Thus, to control an epidemic, authorities should be well equipped with proper resources to avoid the delay.

The model can be used to describe the behavior of infectious diseases with latency period, such as COVID-19, influenza, and smallpox.

Funding information: No funding.
Conflict of interest: Authors state that no conflict of interest.
Ethical approval: The conducted research is not related to either human or animal use.

References

[1] Brauer, F., & Castillo-Chavez, C. (2012). Mathematical models in population biology and epidemiology (Vol. 40). New York: Springer.10.1007/978-1-4614-1686-9Search in Google Scholar

[2] Castillo-Chavez, C., & Song, B. (2004). Dynamical models of tuberculosis and their applications. Mathematical Biosciences and Engineering, 1, 361–404.10.3934/mbe.2004.1.361Search in Google Scholar PubMed

[3] Coddington, E., & Carlson, R. (1997). Linear ordinary differential equations. Philadelphia: Society for Industrial and Applied Mathematics.10.1137/1.9781611971439Search in Google Scholar

[4] Diekmann, O., Heesterbeek, J., & Metz , J. (1990). On the definition and the computation of the basic reproduction ratio R0 in models for infectious diseases. Journal of Mathematical Biology, 35, 503–522.10.1007/BF00178324Search in Google Scholar PubMed

[5] Dubey, B., Dubey, P., & Dubey, U. S. (2015). Dynamics of an SIR model with nonlinear incidence and treatment rate. Applications and Applied Mathematics, 10(2), 718–737.Search in Google Scholar

[6] Hethcote, H. (2000). The mathematics of infectious diseases. Society for Industrial and applied Mathematics, 42(4), 599–653.10.1137/S0036144500371907Search in Google Scholar

[7] Li, Y., & Muldowney, J. (1993). On Bendixson’s criterion. Journal of Differential Equations, 106(1), 27–39.10.1006/jdeq.1993.1097Search in Google Scholar

[8] Li, M. (2018). An introduction to mathematical modeling of infectious diseases. Cham, Switzerland: Springer International Publishing. doi: 10.1007/978-3-319-72122-4.Search in Google Scholar

[9] Li, M., & Muldowney, J. (1995). Global stability for the SEIR model in epidemiology. Mathematical Biosciences, 125(2), 155–164.10.1016/0025-5564(95)92756-5Search in Google Scholar PubMed

[10] Li, M., & Muldowney, J. (1996). A geometric approach to global-stability problems. SIAM Journal on Mathematical Analysis, 27(4), 1070–1083.10.1137/S0036141094266449Search in Google Scholar

[11] Pathak, S., Maiti, A., & Samanta, G. (2010). Rich dynamics of an SIR epidemic model. Nonlinear Analysis: Modelling and Control, 15, 71–81.10.15388/NA.2010.15.1.14365Search in Google Scholar

[12] van den Driessche, P., & Watmough, J. (2002). Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission. Mathematical Biosciences, 180, 29–48.10.1016/S0025-5564(02)00108-6Search in Google Scholar PubMed

[13] Wang, W. (2006). Backward bifurcation of an epidemic model with treatment. Mathematical Biosciences, 201, 58–71.10.1016/j.mbs.2005.12.022Search in Google Scholar PubMed

[14] Wang, W., & Ruan, S. (2004). Bifurcation in an epidemic model with constant removal rate of infectives. Journal of Mathematical Analysis and Applications, 291(2), 775–793.10.1016/j.jmaa.2003.11.043Search in Google Scholar

[15] Wang, X. (2004). A simple proof of Descartes’s rule of signs. The American Mathematical Monthly, 111(6), 525–526. doi: 10.2307/4145072.Search in Google Scholar

[16] Xiao, D., & Ruan, S. (2007). Global analysis of an epidemic model with nonmonotone incidence rate. Mathematical Biosciences, 208, 419–429.10.1016/j.mbs.2006.09.025Search in Google Scholar PubMed PubMed Central

[17] Zhang, X., & Xianning, L. (2008). Backward bifurcation of an epidemic model with saturated treatment function. Journal of Mathematical Analysis and Applications, 348, 433–443.10.1016/j.jmaa.2008.07.042Search in Google Scholar

[18] Zhang, J., Jia, J., & Song, X. (2014). Analysis of an SEIR Epidemic Model with Saturated Incidence and Saturated Treatment Function. The Scientific World Journal, 2014, 11. doi: 10.1155/2014/910421.Search in Google Scholar PubMed PubMed Central

Received: 2022-08-18

Revised: 2023-03-04

Accepted: 2023-03-06

Published Online: 2023-05-08

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

Articles in the same Issue

https://doi.org/10.1515/cmb-2022-0146

Keywords for this article

SEIR, basic reproduction number; Holling type II function; Lyapunov function

Creative Commons

BY 4.0