# A Mathematical Model for Nipah Virus Infection

Authors: Assefa Denekew Zewdie; Sunita Gakkhar
Journal: Journal of Applied Mathematics
 (2020)

---

## Abstract

It has been reported that unprotected contact with the dead bodies of infected individuals is a plausible way of Nipah virus transmission. An SIRD model is proposed in this paper to investigate the impact of unprotected contact with dead bodies of infected individuals before burial or cremation and their disposal rate on the dynamics of Nipah virus infection. The model is analyzed, and the reproduction number is computed. It is established that the disease-free state is globally asymptotically stable when the reproduction number is less than unity and unstable if it is greater than unity. By using the central manifold theory, we observe that the endemic equilibrium is locally stable near to unity. It is concluded that minimizing unsafe contact with the infected dead body and/or burial or cremation as fast as possible contributes positively. Further, the numerical simulations for the given choice of data and initial conditions illustrate that the endemic state is stable and the disease persists in the community when the reproduction number is greater than one.


---

## Body

## 1. Introduction

Nipah virus is first identified in the Malaysian village of Sungai Nipah during an outbreak of encephalitis and respiratory illness. Pig farmers got virus from pigs. There were more than 265 cases including 105 deaths from September 1998 to April 1999 [1–3]. Later on, it spreads to other south-east Asian regions. It is a zoonotic disease, transmitted from animals to humans caused by a Nipah virus which belongs to the Paramyxoviridae family. It causes a severe illness such as brain inflammations or respiratory infection in animal and human population followed by death. The disease is highly infectious and spreads through the community with infected animals; however, no prophylaxis or effective treatment is available. In the absence of effective drugs, its treatment is limited to a symptomatic treatment, hospitalization, and quarantine. Educating people and creating awareness is the only way to control the spread of infection due to nonavailability of any vaccine for Nipah virus [4].The natural reservoir of Nipah virus is fruit bats belonging to the Pteropus genus in Pteropodidae family [5]. Over 50 species of Pteropus bats are inhabiting the south and south-east Asian regions [4]. Many different pathways are observed for the Nipah virus transmission to human: from fruit bats to human, from fruit bats to animal then to human, and from human to human [2, 6]. Consumption of date palm sap is popular in countries of the south-east Asian region [4], and it is commonly collected in the cold season particularly in rural areas of Bangladesh [7]. Fruit bats also consume and contaminate the date palm and its sap with their secretions, e.g., saliva and urine. Human consumption of fruits or fruit products was one of the transmission route of Nipah from fruit bats to human in Bangladesh [6, 8]. The pigs will get infected by consuming dropped fruits that are partially eaten and infected by bats. These infected pigs transmit Nipah virus to human through fluid contact, while handling ill pigs or through consumption of pork [2]. In Malaysia and Singapore, most Nipah infection cases were observed in individuals having frequent contact with sick pigs [9, 10].Further, person to person transmission was observed during the Bangladesh Nipah outbreak of the year 2004 [11, 12]. In the outbreaks of Siliguri, India, the Nipah virus was not found in samples obtained from animals. However, samples of human fluid secretions found in health care settings were positive for Nipah virus. In particular, healthcare workers, family care givers, and hospital visitors who had exposure to secretions of Nipah virus patients with direct contact or through contaminated towels, bed sheets, etc., got infected [13, 14].Some infectious diseases, like Ebola virus disease [15], can transmit from person to person if there is unsafe direct contact with the dead body of infected individuals during funeral ceremonies or last rites. Unprotected contact of corpse of an infected individual before burial/cremation is another way of transmission of the Nipah virus from human to human [16, 17].Many researchers have studied the pathology and epidemiology of the Nipah virus disease, but very few models are available for it and presented as follows. Biswas studied the disease dynamics using SIR basic mathematical model [18]. He further investigated this model and studied the possible control and preventive strategies through optimal control [19, 20]. Optimal control is also carried out in another SIR model [21]. Mondal et al. have proposed an SEIR model to study the dynamics of the disease by incorporating two control parameters (number of quarantined individuals and enhanced personal hygiene) which have not yet used before in a Nipah dynamic model [22]. Shah et al. proposed an SEI model considering bat to human and human to human disease transmission [23]. The unprotected contact of the dead body of the Nipah virus-infected person is one way of transmitting disease from human to human which has not been considered so far. Durgesh et al. in [24] proposed an SVEIR model by considering bat-human interaction, and they analyzed that vaccination has its own role for controlling the disease to spread. Nita et al. also proposed an SEIHD epidemic model with bat-human interaction in [25] and incorporated control measures such as spray insecticides, buried bats, self-prevention, and hospitalization. They also analyzed the dynamics of the disease and optimal control. Agarwal and Singh in [26] proposed an SEI epidemic model for flying foxes and human. They consider a virus compartment incorporating the fractional order differential equation to analyze. Keeping this in view, an SIRD model is proposed and investigated by incorporating unsafe burial or unprotected contact of dead bodies of infected individuals capable of spreading Nipah virus and its disposal rate. The objective of this study is to investigate the role of safe burial or cremation of dead bodies in the control of Nipah virus. This paper is organized according to the following: Section 2 presents the formulation of the proposed model and description of the parameters. In Section 3, the model analysis including the basic reproduction number, equilibrium points, and their stabilities are presented. In Section 4, the numerical simulation is presented and the graphical results are discussed. Finally, the conclusion is presented in Section 5.


## 2. Model Formulation

The total populationNt is divided into three mutually disjoint compartments, namely, susceptible St, infected It, and recovered Rt such that N=S+I+R. The susceptible are those individuals who are not infected, but they will become infected when they come in contact with infectious individuals. Infected are those individuals who are infectious and may infect others. Recovered are individuals who have recovered from infectiousness either with treatment or by their own. Moreover, the dead bodies of Nipah-infected individuals can also transmit the disease. This model incorporates deceased body compartment Dt; it represents the number of unburied dead bodies of infected individuals. The natural death of susceptible and recovered is excluded from D compartment. The constant parameters Λ and μ are the natural birth and death rates, respectively, while θ is the rate of loss of temporary immunity acquired by recovered individuals. Considering the parameter β1 to be the effective unprotected contact rate of susceptible individuals to get infection from dead bodies of Nipah-infected individuals and β2 to be the effective contact rate of susceptible individuals to get infection from infected individuals, then the force of infection representing the effective transmission rate λ is given by
(1)λ=κβ1D+β2IN.Parameterκ represents a constant fraction of unsafe or unprotected handling of dead bodies leading to the spread of Nipah infection. Further, there is no possibility that susceptible individuals can get infection from the dead body of infected individuals when κ=0. The susceptible class is increased by Λ and θ rates and decreased by λ and μ rates. The infected compartment is increased by λ and decreased by γ, μ, and δ rates due to developing immunity, natural death, and disease-induced death, respectively. The recovered compartment is increased by γ rate and decreased by μ and θ rates. The deceased body compartment is increased by μ and δ rates released from infected class and buried with α rate. From those that are dead, κ rates of them can transfer the disease before burial due to unprotected contact. This model is described by the flow diagram given in (Figure 1) and nonlinear ordinary differential equations given in (2). The description and values of parameters used in the model are given in (Table 1).
(2)dSdt=Λ−λS−μS+θR≡f1S,I,R,D,t,dIdt=λS−μ+δ+γI≡f2S,I,R,D,t,dRdt=γI−μ+θR≡f3S,I,R,D,t,dDdt=μ+δI−αD≡f4S,I,R,D,t.Figure 1
Flow chart of the Nipah model.Table 1
Description and values of parameters.
ParameterDescriptionParameter valuesSourceΛRecruitment rate to susceptible class6102 day−1[22]β1Rate of infection from dead bodies of infected individuals0.65 day−1Assumedβ2Rate of infection from infected individuals0.75 day−1[22]θRate of loss of immunity0.85 day−1AssumedγRate of recovered individuals from infected class0.09 day−1[22]αRate of disposition (burial/cremation) of dead bodiesVariableAssumedμNatural death rate0.000038642 day−1[22]δRate of disease-induced death0.76 day−1[22]κFraction of dead bodies that are not handled safely[0, 1]AssumedThe following nonnegative initial conditions are associated with the dynamical system given in the model (2):
(3)S0=S0>0,I0=I0≥0,R0=R0≥0,D0=D0≥0.


## 3. Model Analysis

It can be easily show that the functionsfi,i=1,2,3,4 are sufficiently smooth and satisfy the Lipchitz condition in R+4. Accordingly, the IVP in the system (2) admits a unique solution in R+4.


### 3.1. Invariant Region and Positivity

Lemma 1.
The system in (2) with initial conditions given in (3) has a positive invariant solution in the regionΩ=S,I,R,D∈R+4:N=S+I+R≤Λ/μ,D≤Λμ+δ/αμ.Proof. All existing solutions starting from nonnegative initial conditions remain nonnegative for all timet≥0. It can be shown by contradiction [27]. Let there exists a time t1>0 such that
(4)St1=0,St≥0,It≥0,Rt≥0,Dt≥0,0≤t≤t1.There exists t2>0 such that
(5)It2=0,St≥0,It≥0,Rt≥0,Dt≥0,0≤t≤t2.There exists a t3>0 such that
(6)Rt3=0,St≥0,It≥0,Rt≥0,Dt≥0,0≤t≤t3.There exists t4>0 such that
(7)Dt4=0,St≥0,It≥0,Rt≥0,Dt≥0,0≤t≤t4.From the first equation of the system in model (2) and from case (4), we have
(8)S′t1=Λ−λSt1−μSt1+θRt1=Λ+θRt1>0.This meansSt<0, which contradicts the fact that St is initially nonnegative; it implies St is positive. Similarly, one can show that the solution of all state variables is positive.The total population of individuals is given byN=S+I+R and adding equations of the system in model (2), we have
(9)dNtdt=Λ−μNt−δIt≤Λ−μNt.It follows that:(10)Nt≤Λμ−Λμ−N0e−μt,where N0=N0 is the initial value of the total population and as t⟶∞,Nt≤Λ/μ and
(11)dDtdt=μ+δIt−αDtButNt=St+It+Rt≤Λ/μ⇒It≤Λ/μ as t≥0.It follows that(12)dDtdt=μ+δIt−αDt⇒Dt≤Λμ+δαμ−Λμ+δαμ−D0e−αt,and as t⟶∞,Dt⟶Λμ+δ/αμ.Hence, the regionΩ is positively invariant for the system (2) and it attracts all solutions of the equations in the system.


### 3.2. Reproduction Number

In the absence of the disease, the system given in model (2) has an equilibrium point E0=Λ/μ,0,0,0. The infectious class is X=I,DT, and the rate of appearance of new infection in each infectious class is denoted by FX and given by
(13)Fx=κβ1D+β2ISN,0T.The rate of other transition between infectious classes is denoted byVX and given by
(14)VX=γ+μ+δI,−μ+δI+αDT,where
(15)JF=∂F1∂I∂F1∂D∂F2∂I∂F2∂D=β2SN−κβ1D+β2ISN2κβ1SN−κβ1D+β2ISN200,JV=∂V1∂I∂V1∂D∂V2∂I∂V2∂D=γ+μ+δ0−μ+δα.At disease-free equilibrium pointE0, the matrices JF and JV will give us
(16)JFE0=F=β2κβ100,⇒JVE0=V=γ+μ+δ0−μ+δα,V−1=1γ+μ+δ0μ+δαγ+μ+δ1α,⇒FV−1=κδ+μβ1+αβ2αγ+μ+δκβ1α00.Now the spectral radiusρFV−1 of the matrix FV−1 is given by
(17)ρFV−1=κδ+μβ1+αβ2αγ+μ+δ.Applying next-generation matrix approach [28], the reproduction number RD is computed as
(18)RD=κδ+μβ1+αβ2αγ+μ+δ=R0+κδ+μβ1αγ+μ+δ,where R0 is the basic reproduction number of the SIR model when κ=0; that is, all infected dead bodies are safely buried. Further, it may be noted that RD increases as κ increases and it has an inverse relationship with α; that means the quicker the dead bodies are buried or cremated, the smaller will be the reproduction number.


### 3.3. Existence of Equilibrium Point

The system given in model (2) always admits a disease-free equilibrium E0=Λ/μ,0,0,0. Further, let the μ system have an endemic equilibrium point, denoted by E∗=S∗,I∗,R∗,D∗. Then, the force of infection from equation (1) at this equilibrium point is
(19)λ∗=κβ1D∗+β2I∗N∗,N∗=S∗+I∗+R∗.Substituting and simplifying yields(20)λ∗+aλ∗=0⇒λ∗=0,λ∗=−a,where
(21)a=θ+μαγ+μ+δ−κδ+μβ1+α+β2αμ+θ+γ=θ+μγ+μ+δμ+θ+γ1−RD.ForRD<1 equation (21) has only one nonnegative solution λ∗=0, which is a disease-free equilibrium. But it has two nonnegative solutions, one is the disease free and the other is endemic equilibrium for RD>1. This result gives us the endemic equilibrium point exists in the system of model (2) and the following lemma.Lemma 2.
ForRD>1, the system in model (2) has a unique endemic equilibrium point E∗ given by, E∗=S∗,I∗,R∗,D∗, where
(22)S∗=Λγ+μ+θδθ+μRD−1+μγ+θ+μRD,I∗=Λθ+μRD−1δθ+μRD−1+μγ+θ+μRD,R∗=ΛγRD−1δθ+μRD−1+μγ+θ+μRD,D∗=Λδ+μθ+μR0−1αδθ+μRD−1+μγ+θ+μRD.It is noted that the endemic point reduces to disease-free equilibrium point forRD=1.


### 3.4. Stability Analysis

#### 3.4.1. Local Stability of Equilibrium Points

Proposition 1.
The disease-free equilibrium point of the system given in model (2) is locally asymptotically stable if RD<1 and unstable ifRD>1.Proof. The Jacobian matrix of system (2) is
(23)JX=∂f1∂S∂f1∂I∂f1∂R∂f1∂D∂f2∂S∂f2∂I∂f2∂R∂f2∂D∂f3∂S∂f3∂I∂f3∂R∂f3∂D∂f4∂S∂f4∂I∂f4∂R∂f4∂D,where X=S,I,R,DT.At disease-free equilibrium point, the Jacobian matrix gives(24)JE0=−μ−β2θ−κβ10β2−γ+μ+δ0κβ10γ−θ+μ00δ+μ0−α.The Jacobian matrix (24) of the system in (2) at the disease-free equilibrium point gives the following polynomial characteristic equation:
(25)λ+µλ+θ+µλ2+Bλ+C=0,where
(26)B=α+γ+δ+µ−β2,C=αγ+δ+µ−β2−κδ+µβ1.Here, we haveλ1=−μ, λ2=−θ+μ, and other solutions are negative provided that B and C are positive, this is possible that when RD<1.
(27)B=α+γ+δ+µ−β2>γ+δ+µ−β2=κδ+μβ1+αβ2αRD−β2sinceRDin18=κδ+μβ1+αβ21−RDαRD>0,forRD<1,C=αγ+δ+μ−β2−κδ+μβ1=αγ+δ+μ1−RD>0,forRD<1.ForRD=1, the quadratic factor in equation (25) has one zero eigenvalue; it follows that the Jacobian matrix (24) has a simple zero eigenvalue and other eigenvalues are negative and the equilibrium point is nonhyperbolic at RD=1. Linearization does not show the stability of nonhyperbolic equilibrium points, so we will analyze it using the central manifold theory [29, 30] stated below.Theorem 1.
Consider the following general system of ordinary differential equations with a bifurcation parameterϕ,
(28)dxdt=fx,ϕ,f:ℝn×ℝ⟶ℝn,f∈ℂ2ℝn⟶ℝ.Without loss of generality, it is assumed that 0 is an equilibrium for system (28) for all values of the parameter ϕ; that is,f0,ϕ≡0for allϕ. Assume thatA1. A=Dxf0,0=∂fi/∂xj0,0is the linearization matrix of system ((28)) around the equilibrium 0 withϕevaluated at0. Zero is a simple eigenvalue ofA, and all other eigenvalues ofAhave negative real parts;A2. MatrixAhas a nonnegative right eigenvectorwand a left eigenvectorvcorresponding to the zero eigenvalue.Letfkbe thekthcomponent offand(29)a=∑k,i,j=1nvkwiwj∂2fk∂xi∂xj0,0,b=∑k,i,j=1nvkwi∂2fk∂xi∂xj0,0.The local dynamics of the ODE in ((28)) around 0 is totally determined byaandbgiven in ((29)).(i)
a>0,b>0. When ϕ<0 with ∣ϕ∣<<1, 0 is locally asymptotically stable, and there exists a positive unstable equilibrium; when 0<ϕ<<1, 0 is unstable and there exists a negative and locally asymptotically stable equilibrium.(ii)
a<0,b<0. When ϕ <0 with ∣ϕ∣<<1, 0 is unstable; when 0<ϕ<<1, 0 is locally asymptotically stable, and there exists a positive unstable equilibrium.(iii)
a>0,b<0. When ϕ<0 with ∣ϕ∣<<1, 0 is unstable, and there exists a locally asymptotically stable negative equilibrium; when 0<ϕ<<1, 0 is stable, and a positive unstable equilibrium appears.(iv)
a<0,b>0. When ϕ changes from negative to positive, 0 changes its stability from stable to unstable. Correspondingly, a negative unstable equilibrium becomes positive and locally asymptotically stable.Now settingS=x1, I=x2, R=x3, and D=x4, we can rewrite the system of model (2) as follows:
(30)dx1dt=Λ−κβ1x4+β2x2x1x1+x2+x3−μx1+θx3≔f1,dx2dt=κβ1x4+β2x2x1x1+x2+x3−γ+δ+μx2≔f2,dx3dt=γx2−θ+μx3:f3,dx4dt=δ+μx2−αx4≔f4.WithRD=1 corresponding to β1=β1∗ where
(31)β1∗=αγ+μ+δ−αβ2κδ+μ,which is assumed to be a bifurcation parameter. The disease-free equilibrium is
(32)X0=x10,x20,x30,x40=Λμ,0,0,0.The linearization matrix of system (30) around the disease-free equilibrium when β1=β1∗ is
(33)DXf=−μ−β2θ−κβ1∗0β2−γ+μ+δ0κβ1∗0γ−θ+μ00δ+μ0α.It is clear thatDxf has a simple zero eigenvalue and a right eigenvector corresponding to the zero eigenvalue is w=w1,w2,w3,w4, where
(34)w1=γμ+δ+μμ+θαμμ+θ,w2=α,w3=αγμ+θ,w4=μ+δ,and the left eigenvector associated with the zero eigenvalue satisfying w.v=1 is v=v1,v2,v3,v4.
(35)v1=0,v2=1α+γ+μ+δ−β2,v3=0,v4=1δ+μ.Based on the theoretical result given in Theorem1, we have to compute a and b at X0,β1∗ where
(36)a=∑k,i,j=14vkwiwj∂2fk∂xi∂xjX0,β1∗,b=∑k,i=14vkwi∂2fk∂xi∂β1X0,β1∗.Sincev1=v3=0, equations a and b in (36) becomes
(37)a=v2∑i,j=14wiwj∂2f2∂xi∂xjX0,β1∗+v4∑i,j=14wiwj∂2f4∂xi∂xjX0,β1∗,b=v2∑i=14wiwj∂2f2∂xi∂β1X0,β1∗+v4∑i=14wi∂2f4∂xi∂β1X0,β1∗.Substituting the eigenvectors and the computed partial derivatives of the system (30) at X0,β1∗ in the formula for a and b in (37), after some algebraic computation, we get
(38)a=−2μv2w2+w3Λw2β2+w4κβ1∗<0,b=v2w4κ>0.This indicates that atβ1=β1∗ for RD=1, the system exhibits a transcritical bifurcation (i.e., the disease-free equilibrium point changes its stability from locally stable for RD<1 to unstable for RD>1). Consequently, Theorem 1 is guaranteed for the following result.Theorem 2.
The system in (2) exhibits a forward bifurcation at β1=β1∗forRD=1 and a locally stable positive endemic equilibrium will appear whenever RD>1 but near to unity.


#### 3.4.2. Global Stability of Disease-Free Equilibrium Point

Proposition 2.
The disease-free equilibrium point of the system given (2) is globally asymptotically stable if RD<1.Proof. Leta be an arbitrary chosen positive constant. Consider a nonnegative definite function on the invariant set Ω as
(39)V=I+aD,a≥0.Taking the derivative ofV gives V̇=İ+aḊ⇒V̇=κβ1D+β2IS/N−γ+δ+μI+αμ+δI+αDSinceSt≤Nt for all t≥0 in the domain Ω, it follows that
(40)V̇≤κβ1D+β2I−γ+µ+δI+aµ+δI−aαD=β2−γ+µ+δ+aµ+δI+κβ1−aαD.We can choosea such that κβ1−aα=0; that is, a=κβ1/α.
(41)⇒V̇≤β2−γ+μ+δ+κβ1αμ+δI=γ+µ+δI=γ+µ+δIκδ+μβ1+αβ2αγ+δ+μ−1=γ+μ+δIRD−1≤0,forRD≤1Since all parameters and variables in the model (2) are nonnegative, the derivative of the Lyapunov function is V̇≤0 if RD≤1 and V̇=0 if and only if I=0. It follows that from LaSalle’s invariance principle [31], the only invariant singleton set in Ω is the set E0.Hence, the disease-free equilibrium point is globally asymptotically stable forRD<1 in the invariant region Ω.


## 3.1. Invariant Region and Positivity

Lemma 1.
The system in (2) with initial conditions given in (3) has a positive invariant solution in the regionΩ=S,I,R,D∈R+4:N=S+I+R≤Λ/μ,D≤Λμ+δ/αμ.Proof. All existing solutions starting from nonnegative initial conditions remain nonnegative for all timet≥0. It can be shown by contradiction [27]. Let there exists a time t1>0 such that
(4)St1=0,St≥0,It≥0,Rt≥0,Dt≥0,0≤t≤t1.There exists t2>0 such that
(5)It2=0,St≥0,It≥0,Rt≥0,Dt≥0,0≤t≤t2.There exists a t3>0 such that
(6)Rt3=0,St≥0,It≥0,Rt≥0,Dt≥0,0≤t≤t3.There exists t4>0 such that
(7)Dt4=0,St≥0,It≥0,Rt≥0,Dt≥0,0≤t≤t4.From the first equation of the system in model (2) and from case (4), we have
(8)S′t1=Λ−λSt1−μSt1+θRt1=Λ+θRt1>0.This meansSt<0, which contradicts the fact that St is initially nonnegative; it implies St is positive. Similarly, one can show that the solution of all state variables is positive.The total population of individuals is given byN=S+I+R and adding equations of the system in model (2), we have
(9)dNtdt=Λ−μNt−δIt≤Λ−μNt.It follows that:(10)Nt≤Λμ−Λμ−N0e−μt,where N0=N0 is the initial value of the total population and as t⟶∞,Nt≤Λ/μ and
(11)dDtdt=μ+δIt−αDtButNt=St+It+Rt≤Λ/μ⇒It≤Λ/μ as t≥0.It follows that(12)dDtdt=μ+δIt−αDt⇒Dt≤Λμ+δαμ−Λμ+δαμ−D0e−αt,and as t⟶∞,Dt⟶Λμ+δ/αμ.Hence, the regionΩ is positively invariant for the system (2) and it attracts all solutions of the equations in the system.


## 3.2. Reproduction Number

In the absence of the disease, the system given in model (2) has an equilibrium point E0=Λ/μ,0,0,0. The infectious class is X=I,DT, and the rate of appearance of new infection in each infectious class is denoted by FX and given by
(13)Fx=κβ1D+β2ISN,0T.The rate of other transition between infectious classes is denoted byVX and given by
(14)VX=γ+μ+δI,−μ+δI+αDT,where
(15)JF=∂F1∂I∂F1∂D∂F2∂I∂F2∂D=β2SN−κβ1D+β2ISN2κβ1SN−κβ1D+β2ISN200,JV=∂V1∂I∂V1∂D∂V2∂I∂V2∂D=γ+μ+δ0−μ+δα.At disease-free equilibrium pointE0, the matrices JF and JV will give us
(16)JFE0=F=β2κβ100,⇒JVE0=V=γ+μ+δ0−μ+δα,V−1=1γ+μ+δ0μ+δαγ+μ+δ1α,⇒FV−1=κδ+μβ1+αβ2αγ+μ+δκβ1α00.Now the spectral radiusρFV−1 of the matrix FV−1 is given by
(17)ρFV−1=κδ+μβ1+αβ2αγ+μ+δ.Applying next-generation matrix approach [28], the reproduction number RD is computed as
(18)RD=κδ+μβ1+αβ2αγ+μ+δ=R0+κδ+μβ1αγ+μ+δ,where R0 is the basic reproduction number of the SIR model when κ=0; that is, all infected dead bodies are safely buried. Further, it may be noted that RD increases as κ increases and it has an inverse relationship with α; that means the quicker the dead bodies are buried or cremated, the smaller will be the reproduction number.


## 3.3. Existence of Equilibrium Point

The system given in model (2) always admits a disease-free equilibrium E0=Λ/μ,0,0,0. Further, let the μ system have an endemic equilibrium point, denoted by E∗=S∗,I∗,R∗,D∗. Then, the force of infection from equation (1) at this equilibrium point is
(19)λ∗=κβ1D∗+β2I∗N∗,N∗=S∗+I∗+R∗.Substituting and simplifying yields(20)λ∗+aλ∗=0⇒λ∗=0,λ∗=−a,where
(21)a=θ+μαγ+μ+δ−κδ+μβ1+α+β2αμ+θ+γ=θ+μγ+μ+δμ+θ+γ1−RD.ForRD<1 equation (21) has only one nonnegative solution λ∗=0, which is a disease-free equilibrium. But it has two nonnegative solutions, one is the disease free and the other is endemic equilibrium for RD>1. This result gives us the endemic equilibrium point exists in the system of model (2) and the following lemma.Lemma 2.
ForRD>1, the system in model (2) has a unique endemic equilibrium point E∗ given by, E∗=S∗,I∗,R∗,D∗, where
(22)S∗=Λγ+μ+θδθ+μRD−1+μγ+θ+μRD,I∗=Λθ+μRD−1δθ+μRD−1+μγ+θ+μRD,R∗=ΛγRD−1δθ+μRD−1+μγ+θ+μRD,D∗=Λδ+μθ+μR0−1αδθ+μRD−1+μγ+θ+μRD.It is noted that the endemic point reduces to disease-free equilibrium point forRD=1.


## 3.4. Stability Analysis

### 3.4.1. Local Stability of Equilibrium Points

Proposition 1.
The disease-free equilibrium point of the system given in model (2) is locally asymptotically stable if RD<1 and unstable ifRD>1.Proof. The Jacobian matrix of system (2) is
(23)JX=∂f1∂S∂f1∂I∂f1∂R∂f1∂D∂f2∂S∂f2∂I∂f2∂R∂f2∂D∂f3∂S∂f3∂I∂f3∂R∂f3∂D∂f4∂S∂f4∂I∂f4∂R∂f4∂D,where X=S,I,R,DT.At disease-free equilibrium point, the Jacobian matrix gives(24)JE0=−μ−β2θ−κβ10β2−γ+μ+δ0κβ10γ−θ+μ00δ+μ0−α.The Jacobian matrix (24) of the system in (2) at the disease-free equilibrium point gives the following polynomial characteristic equation:
(25)λ+µλ+θ+µλ2+Bλ+C=0,where
(26)B=α+γ+δ+µ−β2,C=αγ+δ+µ−β2−κδ+µβ1.Here, we haveλ1=−μ, λ2=−θ+μ, and other solutions are negative provided that B and C are positive, this is possible that when RD<1.
(27)B=α+γ+δ+µ−β2>γ+δ+µ−β2=κδ+μβ1+αβ2αRD−β2sinceRDin18=κδ+μβ1+αβ21−RDαRD>0,forRD<1,C=αγ+δ+μ−β2−κδ+μβ1=αγ+δ+μ1−RD>0,forRD<1.ForRD=1, the quadratic factor in equation (25) has one zero eigenvalue; it follows that the Jacobian matrix (24) has a simple zero eigenvalue and other eigenvalues are negative and the equilibrium point is nonhyperbolic at RD=1. Linearization does not show the stability of nonhyperbolic equilibrium points, so we will analyze it using the central manifold theory [29, 30] stated below.Theorem 1.
Consider the following general system of ordinary differential equations with a bifurcation parameterϕ,
(28)dxdt=fx,ϕ,f:ℝn×ℝ⟶ℝn,f∈ℂ2ℝn⟶ℝ.Without loss of generality, it is assumed that 0 is an equilibrium for system (28) for all values of the parameter ϕ; that is,f0,ϕ≡0for allϕ. Assume thatA1. A=Dxf0,0=∂fi/∂xj0,0is the linearization matrix of system ((28)) around the equilibrium 0 withϕevaluated at0. Zero is a simple eigenvalue ofA, and all other eigenvalues ofAhave negative real parts;A2. MatrixAhas a nonnegative right eigenvectorwand a left eigenvectorvcorresponding to the zero eigenvalue.Letfkbe thekthcomponent offand(29)a=∑k,i,j=1nvkwiwj∂2fk∂xi∂xj0,0,b=∑k,i,j=1nvkwi∂2fk∂xi∂xj0,0.The local dynamics of the ODE in ((28)) around 0 is totally determined byaandbgiven in ((29)).(i)
a>0,b>0. When ϕ<0 with ∣ϕ∣<<1, 0 is locally asymptotically stable, and there exists a positive unstable equilibrium; when 0<ϕ<<1, 0 is unstable and there exists a negative and locally asymptotically stable equilibrium.(ii)
a<0,b<0. When ϕ <0 with ∣ϕ∣<<1, 0 is unstable; when 0<ϕ<<1, 0 is locally asymptotically stable, and there exists a positive unstable equilibrium.(iii)
a>0,b<0. When ϕ<0 with ∣ϕ∣<<1, 0 is unstable, and there exists a locally asymptotically stable negative equilibrium; when 0<ϕ<<1, 0 is stable, and a positive unstable equilibrium appears.(iv)
a<0,b>0. When ϕ changes from negative to positive, 0 changes its stability from stable to unstable. Correspondingly, a negative unstable equilibrium becomes positive and locally asymptotically stable.Now settingS=x1, I=x2, R=x3, and D=x4, we can rewrite the system of model (2) as follows:
(30)dx1dt=Λ−κβ1x4+β2x2x1x1+x2+x3−μx1+θx3≔f1,dx2dt=κβ1x4+β2x2x1x1+x2+x3−γ+δ+μx2≔f2,dx3dt=γx2−θ+μx3:f3,dx4dt=δ+μx2−αx4≔f4.WithRD=1 corresponding to β1=β1∗ where
(31)β1∗=αγ+μ+δ−αβ2κδ+μ,which is assumed to be a bifurcation parameter. The disease-free equilibrium is
(32)X0=x10,x20,x30,x40=Λμ,0,0,0.The linearization matrix of system (30) around the disease-free equilibrium when β1=β1∗ is
(33)DXf=−μ−β2θ−κβ1∗0β2−γ+μ+δ0κβ1∗0γ−θ+μ00δ+μ0α.It is clear thatDxf has a simple zero eigenvalue and a right eigenvector corresponding to the zero eigenvalue is w=w1,w2,w3,w4, where
(34)w1=γμ+δ+μμ+θαμμ+θ,w2=α,w3=αγμ+θ,w4=μ+δ,and the left eigenvector associated with the zero eigenvalue satisfying w.v=1 is v=v1,v2,v3,v4.
(35)v1=0,v2=1α+γ+μ+δ−β2,v3=0,v4=1δ+μ.Based on the theoretical result given in Theorem1, we have to compute a and b at X0,β1∗ where
(36)a=∑k,i,j=14vkwiwj∂2fk∂xi∂xjX0,β1∗,b=∑k,i=14vkwi∂2fk∂xi∂β1X0,β1∗.Sincev1=v3=0, equations a and b in (36) becomes
(37)a=v2∑i,j=14wiwj∂2f2∂xi∂xjX0,β1∗+v4∑i,j=14wiwj∂2f4∂xi∂xjX0,β1∗,b=v2∑i=14wiwj∂2f2∂xi∂β1X0,β1∗+v4∑i=14wi∂2f4∂xi∂β1X0,β1∗.Substituting the eigenvectors and the computed partial derivatives of the system (30) at X0,β1∗ in the formula for a and b in (37), after some algebraic computation, we get
(38)a=−2μv2w2+w3Λw2β2+w4κβ1∗<0,b=v2w4κ>0.This indicates that atβ1=β1∗ for RD=1, the system exhibits a transcritical bifurcation (i.e., the disease-free equilibrium point changes its stability from locally stable for RD<1 to unstable for RD>1). Consequently, Theorem 1 is guaranteed for the following result.Theorem 2.
The system in (2) exhibits a forward bifurcation at β1=β1∗forRD=1 and a locally stable positive endemic equilibrium will appear whenever RD>1 but near to unity.


### 3.4.2. Global Stability of Disease-Free Equilibrium Point

Proposition 2.
The disease-free equilibrium point of the system given (2) is globally asymptotically stable if RD<1.Proof. Leta be an arbitrary chosen positive constant. Consider a nonnegative definite function on the invariant set Ω as
(39)V=I+aD,a≥0.Taking the derivative ofV gives V̇=İ+aḊ⇒V̇=κβ1D+β2IS/N−γ+δ+μI+αμ+δI+αDSinceSt≤Nt for all t≥0 in the domain Ω, it follows that
(40)V̇≤κβ1D+β2I−γ+µ+δI+aµ+δI−aαD=β2−γ+µ+δ+aµ+δI+κβ1−aαD.We can choosea such that κβ1−aα=0; that is, a=κβ1/α.
(41)⇒V̇≤β2−γ+μ+δ+κβ1αμ+δI=γ+µ+δI=γ+µ+δIκδ+μβ1+αβ2αγ+δ+μ−1=γ+μ+δIRD−1≤0,forRD≤1Since all parameters and variables in the model (2) are nonnegative, the derivative of the Lyapunov function is V̇≤0 if RD≤1 and V̇=0 if and only if I=0. It follows that from LaSalle’s invariance principle [31], the only invariant singleton set in Ω is the set E0.Hence, the disease-free equilibrium point is globally asymptotically stable forRD<1 in the invariant region Ω.


## 3.4.1. Local Stability of Equilibrium Points

Proposition 1.
The disease-free equilibrium point of the system given in model (2) is locally asymptotically stable if RD<1 and unstable ifRD>1.Proof. The Jacobian matrix of system (2) is
(23)JX=∂f1∂S∂f1∂I∂f1∂R∂f1∂D∂f2∂S∂f2∂I∂f2∂R∂f2∂D∂f3∂S∂f3∂I∂f3∂R∂f3∂D∂f4∂S∂f4∂I∂f4∂R∂f4∂D,where X=S,I,R,DT.At disease-free equilibrium point, the Jacobian matrix gives(24)JE0=−μ−β2θ−κβ10β2−γ+μ+δ0κβ10γ−θ+μ00δ+μ0−α.The Jacobian matrix (24) of the system in (2) at the disease-free equilibrium point gives the following polynomial characteristic equation:
(25)λ+µλ+θ+µλ2+Bλ+C=0,where
(26)B=α+γ+δ+µ−β2,C=αγ+δ+µ−β2−κδ+µβ1.Here, we haveλ1=−μ, λ2=−θ+μ, and other solutions are negative provided that B and C are positive, this is possible that when RD<1.
(27)B=α+γ+δ+µ−β2>γ+δ+µ−β2=κδ+μβ1+αβ2αRD−β2sinceRDin18=κδ+μβ1+αβ21−RDαRD>0,forRD<1,C=αγ+δ+μ−β2−κδ+μβ1=αγ+δ+μ1−RD>0,forRD<1.ForRD=1, the quadratic factor in equation (25) has one zero eigenvalue; it follows that the Jacobian matrix (24) has a simple zero eigenvalue and other eigenvalues are negative and the equilibrium point is nonhyperbolic at RD=1. Linearization does not show the stability of nonhyperbolic equilibrium points, so we will analyze it using the central manifold theory [29, 30] stated below.Theorem 1.
Consider the following general system of ordinary differential equations with a bifurcation parameterϕ,
(28)dxdt=fx,ϕ,f:ℝn×ℝ⟶ℝn,f∈ℂ2ℝn⟶ℝ.Without loss of generality, it is assumed that 0 is an equilibrium for system (28) for all values of the parameter ϕ; that is,f0,ϕ≡0for allϕ. Assume thatA1. A=Dxf0,0=∂fi/∂xj0,0is the linearization matrix of system ((28)) around the equilibrium 0 withϕevaluated at0. Zero is a simple eigenvalue ofA, and all other eigenvalues ofAhave negative real parts;A2. MatrixAhas a nonnegative right eigenvectorwand a left eigenvectorvcorresponding to the zero eigenvalue.Letfkbe thekthcomponent offand(29)a=∑k,i,j=1nvkwiwj∂2fk∂xi∂xj0,0,b=∑k,i,j=1nvkwi∂2fk∂xi∂xj0,0.The local dynamics of the ODE in ((28)) around 0 is totally determined byaandbgiven in ((29)).(i)
a>0,b>0. When ϕ<0 with ∣ϕ∣<<1, 0 is locally asymptotically stable, and there exists a positive unstable equilibrium; when 0<ϕ<<1, 0 is unstable and there exists a negative and locally asymptotically stable equilibrium.(ii)
a<0,b<0. When ϕ <0 with ∣ϕ∣<<1, 0 is unstable; when 0<ϕ<<1, 0 is locally asymptotically stable, and there exists a positive unstable equilibrium.(iii)
a>0,b<0. When ϕ<0 with ∣ϕ∣<<1, 0 is unstable, and there exists a locally asymptotically stable negative equilibrium; when 0<ϕ<<1, 0 is stable, and a positive unstable equilibrium appears.(iv)
a<0,b>0. When ϕ changes from negative to positive, 0 changes its stability from stable to unstable. Correspondingly, a negative unstable equilibrium becomes positive and locally asymptotically stable.Now settingS=x1, I=x2, R=x3, and D=x4, we can rewrite the system of model (2) as follows:
(30)dx1dt=Λ−κβ1x4+β2x2x1x1+x2+x3−μx1+θx3≔f1,dx2dt=κβ1x4+β2x2x1x1+x2+x3−γ+δ+μx2≔f2,dx3dt=γx2−θ+μx3:f3,dx4dt=δ+μx2−αx4≔f4.WithRD=1 corresponding to β1=β1∗ where
(31)β1∗=αγ+μ+δ−αβ2κδ+μ,which is assumed to be a bifurcation parameter. The disease-free equilibrium is
(32)X0=x10,x20,x30,x40=Λμ,0,0,0.The linearization matrix of system (30) around the disease-free equilibrium when β1=β1∗ is
(33)DXf=−μ−β2θ−κβ1∗0β2−γ+μ+δ0κβ1∗0γ−θ+μ00δ+μ0α.It is clear thatDxf has a simple zero eigenvalue and a right eigenvector corresponding to the zero eigenvalue is w=w1,w2,w3,w4, where
(34)w1=γμ+δ+μμ+θαμμ+θ,w2=α,w3=αγμ+θ,w4=μ+δ,and the left eigenvector associated with the zero eigenvalue satisfying w.v=1 is v=v1,v2,v3,v4.
(35)v1=0,v2=1α+γ+μ+δ−β2,v3=0,v4=1δ+μ.Based on the theoretical result given in Theorem1, we have to compute a and b at X0,β1∗ where
(36)a=∑k,i,j=14vkwiwj∂2fk∂xi∂xjX0,β1∗,b=∑k,i=14vkwi∂2fk∂xi∂β1X0,β1∗.Sincev1=v3=0, equations a and b in (36) becomes
(37)a=v2∑i,j=14wiwj∂2f2∂xi∂xjX0,β1∗+v4∑i,j=14wiwj∂2f4∂xi∂xjX0,β1∗,b=v2∑i=14wiwj∂2f2∂xi∂β1X0,β1∗+v4∑i=14wi∂2f4∂xi∂β1X0,β1∗.Substituting the eigenvectors and the computed partial derivatives of the system (30) at X0,β1∗ in the formula for a and b in (37), after some algebraic computation, we get
(38)a=−2μv2w2+w3Λw2β2+w4κβ1∗<0,b=v2w4κ>0.This indicates that atβ1=β1∗ for RD=1, the system exhibits a transcritical bifurcation (i.e., the disease-free equilibrium point changes its stability from locally stable for RD<1 to unstable for RD>1). Consequently, Theorem 1 is guaranteed for the following result.Theorem 2.
The system in (2) exhibits a forward bifurcation at β1=β1∗forRD=1 and a locally stable positive endemic equilibrium will appear whenever RD>1 but near to unity.


## 3.4.2. Global Stability of Disease-Free Equilibrium Point

Proposition 2.
The disease-free equilibrium point of the system given (2) is globally asymptotically stable if RD<1.Proof. Leta be an arbitrary chosen positive constant. Consider a nonnegative definite function on the invariant set Ω as
(39)V=I+aD,a≥0.Taking the derivative ofV gives V̇=İ+aḊ⇒V̇=κβ1D+β2IS/N−γ+δ+μI+αμ+δI+αDSinceSt≤Nt for all t≥0 in the domain Ω, it follows that
(40)V̇≤κβ1D+β2I−γ+µ+δI+aµ+δI−aαD=β2−γ+µ+δ+aµ+δI+κβ1−aαD.We can choosea such that κβ1−aα=0; that is, a=κβ1/α.
(41)⇒V̇≤β2−γ+μ+δ+κβ1αμ+δI=γ+µ+δI=γ+µ+δIκδ+μβ1+αβ2αγ+δ+μ−1=γ+μ+δIRD−1≤0,forRD≤1Since all parameters and variables in the model (2) are nonnegative, the derivative of the Lyapunov function is V̇≤0 if RD≤1 and V̇=0 if and only if I=0. It follows that from LaSalle’s invariance principle [31], the only invariant singleton set in Ω is the set E0.Hence, the disease-free equilibrium point is globally asymptotically stable forRD<1 in the invariant region Ω.


## 4. Numerical Analysis

In this section, some numerical simulations of model (2) are carried out for the data shown in Table 1.The main objective is to understand and investigate the impact of unsafe handlingκ and disposition rate α of the dead bodies of infected individuals during Nipah outbreak. The proposed model is also compared with that of the SIR model where the dead body has no contribution towards transmission of infection. Using the reproduction number RD in (18) and solving RD=1, the critical value of κ is computed as 0.20249, by assuming parameters as given in Table 1 with α=1day−1. From RD versus κ curve, drawn in Figure 2, it is observed that RD<1 for 0≤κ<0.20249 and RD>1 for 0.20249<κ≤1. It means that the high ratio of unsafe burial practice increases the reproduction number, and this is leading to instability of disease-free state (Proposition 1) and existence of endemic state (Lemma 2).Figure 2
Graph of reproduction numberRD as a function of κ. Other parameter values are given in Table 1.Furthermore, if we assume all dead body contacts are unprotectedκ=1 and for our choice of data, we have got the critical value of α=4.9383day−1. This is the value where the reproduction number is close to unity. As shown in Figure 3, the reproduction number is inversely related with the rate of disposing dead bodies. For α>4.9383day−1, the reproduction number is less than one and the infection decreases through time and finally dies out from the community.Figure 3
Graph of reproduction numberRD as a function of α. Other parameter values are given in Table 1.Infection curves for various values ofκ with α=1day−1 are drawn in Figure 4. For different values of κ with RD>1, the infection increases with time up to a peak point then decreases, converging to endemic level. For instance, when κ=1, the infection is sharply increasing up to a peak value and then decreasing to stabilize to the endemic equilibrium level I∗=2107.5. A similar pattern is observed for κ=0.5 and 0.25, with I∗=2107 and I∗=2103.5, respectively. The peak is delayed and decreases with decreasing κ. It may be noted that there is a significant change in the peak values for different values of κ. The equilibrium level is also decreasing with decreasing κ although not at the same rate as for the peak value. Thus, the increase in safe burial decreases the reproduction number which slows down the spread of infection and the peak of infection as well as the equilibrium level is also lowered.Figure 4
The dynamics of infectious population forα=1 with different values of κ and other parameter values are given in Table 1.Furthermore, the remaining values ofκ we used in the graph make the value of RD<1. So, broken line curves in Figure 4 show the infection decreases and stabilizes at the disease-free state with time. When safe burial increases, the infection decreases faster and will die out from the community.In Figure5, it is seen the infection curves for different values of α and assuming all contacts with the dead body are unprotected κ. As the rate of disposal α increases for RD>1, the peak value of the infection decreases with time and finally comes to stabilize at endemic state. For RD<1, the infection decreases fast and stabilize at disease-free state.Figure 5
The dynamics of infectious population forκ=1 with different values of α and other parameter values are given in Table 1.


## 5. Conclusion

In most of human culture, there may be contact with the dead body during last rites due to bathing, expressing love or respect by hugging, touching heads and faces, etc. of the dead body particularly very close relatives. Researchers confirmed such contact of the dead body of Nipah-infected individuals is one way of transmitting disease from person to person. They recommend safe and protected contact for healthcare workers and family caregivers when handling the dead body. In this paper, we have proposed the SIRD model and presented the impact of unprotected contact with the dead bodies of Nipah-infected individuals and their disposal (burial or cremation) on the spread of Nipah infection. It is seen the numerical and analytical behavior of the dynamics that the high fraction of unprotected contact increases the reproduction number. It causes the disease to persist in the community and contributes some difficulty for the control of the outbreak. For the given choice of data, it is observed that the reproduction numberRD is greater than one for 0.20249<κ≤1. However, for 0≤κ<0.20249, RD is less than one. For κ=0, the handling of the human dead body is well protected and no possibility of transmitting disease through it. In this case, the model is reduced to SIR.Furthermore, the rate of disposal of dead bodies is inversely related with the reproduction number. In the absence of protected contactκ=1 and for the choice of data, RD is less than unity when α>4.9383 day−1 and greater than unity when 0<α<4.9383 day−1. If the handling of the human dead body to be buried is unprotected, it should be removed within a quarter of a day to minimize outbreak. This result shows that applying some effort to reduce unprotected contact and removing the dead body of infected individuals within a short period of time will make the reproduction number less than one and it helps in controlling the Nipah outbreak.



---

*Source: 6050834-2020-09-30.xml*