Computational Solutions of a Delay-Driven Stochastic Model for Conjunctivitis Spread

1  Introduction

Conjunctivitis, often called pink eye, is a highly contagious disease. Infection of the conjunctiva is generally a result of viruses, bacteria, and allergens. Some viruses can cause a precise conjunctivitis known as acute hemorrhagic conjunctivitis (AHC). Bacterial conjunctivitis causes the speedy onset of conjunctiva irritation, puffiness of the eyelid, and a muggy release. Usually, signs advance mainly in a single eye and then may extend to the next eye within 2 to 5 days. Three main viruses have been studied and built as the agents responsible for acute hemorrhagic conjunctivitis (AHC): coxsackievirus A24 variation (CA24v), enterovirus 70, and adenovirus 11. Acute hemorrhagic conjunctivitis (AHC) can only occur in a human population and is transmitted through human interaction with an infected specific object, like a towel used by a diseased person. The discomfort of the eyes, sensitivity to light, sore eyes, tears, and pus discharge are some symptoms. Antibiotic eye drops, sanitized passivity, quarantine, and permitting the infection to progress in 2–3 weeks are operative ways to break the extent of conjunctivitis infection [1]. Recently, a significant outbreak occurred in different cities in Pakistan. It began in Karachi and spread to major cities. Almost one hundred thousand cases were reported in Punjab. This massive outbreak of pink eye is a recent development in Pakistan’s history [2].

In 2024, Faisal et al. [3] examined the dynamics of conjunctivitis adenovirus using a unique tool such as bifurcation. At the end, to support the dynamical results, they simulated their numerical result very comprehensively. Gabriela et al. [4] conducted a study in 2024 using a basic Susceptible-Infectious-Susceptible (SIS) model with introducing uncertainty through transmission parameters. They attained conditions of extinction and existence using threshold constraints. They considered real-world problems and proposed suitable uncertain parameters. In 2024, Essak [5] considered three types of uncertainty in their study Lévy noise, Gaussian white noise, and telegraph noise to investigate their influence on the dynamics of disease. Albani et al. [6] used the fact that model parameters evolve with time and showed the randomness of sudden change. They proposed the jump-diffusion stochastic process. They presented their dynamical analysis and compared it with the variation and jumps. They used real data for their simulation and claimed their finding are accurate. In 2024, Ali et al. [7] considered a stochastic model for influenza dynamics and also used accurate parameter values for their comprehensive understanding of influenza dynamics. They used actual data for their simulation that supports their dynamic results. In 2014, Unyong et al. [8] suggested and investigated a mathematical system for the dynamical stability of conjunctivitis with nonlinear incidence terms considered. Finally, numerical results support the analytic results. In 2012, Sangsawang et al. [9] investigated the system for the spread of Hemorrhagic Conjunctivitis. Numerical simulations of the system demonstrate that the disease-free equilibrium and endemic are stable and unique. In 2006, Chowell et al. [10] modeled the diffusion dynamics of AHC in an eruption of pink eye in the host population in Mexico. Therefore, the effects of sustenance present a public health strategy and the efficient message of the eruption, and community health evidence press statements that initiate beings on evading infection and encourage them to pursue a clinical diagnosis. They studied epidemic models with well-known computational methods based on stochastic differential equations with artificial decay terms in [11–13]. They studied well-known stochastic models such as the stochastic within-host Chikungunya (CHIKV) virus model with saturated incidence rate, Extinction, and the stationary distribution of the stochastic hepatitis B virus model, and a Zika virus as a mosquito-borne transmitted disease with environmental fluctuations as presented in [14–16]. The approach aligns with recent works such as [17], where stochastic noise is interpreted as environmental randomness impacting disease spread. Foundational contributions such as those in [18] provide detailed bifurcation analysis of delay epidemic models, highlighting the critical role of time-lag in disease dynamics. Similarly, reviews in [19,20] summarize analytical techniques for deterministic and stochastic delay systems, particularly in biological contexts. Our model structure and analytical proofs draw on these theoretical foundations, especially in establishing positivity, stability, and global existence. The authors studied the Split-Step θ-Method for Stochastic Pantograph Differential Equations, focusing on convergence and mean-square stability analysis, in [21].

Model formation and its numerical simulations are the classical origins of studying the transmission dynamics of infectious diseases because they allow us to understand the complex transmission dynamics of any contagious disease. The trajectories of the deterministic model are idealistic, but the stochastic model is closer to reality. This outcome matters for policymakers and concerned individuals. So, this study took the deterministic model and reformulated it into a stochastic delayed conjunctivitis model, then studied its stochastic behavior to understand the dynamics of the spread of disease, and an exponential delay was introduced to understand the effect of delay on the spread of disease to get more optimal results to control the disease. This work makes several novel contributions to the modeling and control of conjunctivitis epidemics:

•   Model Innovation: The stochastic delayed compartment model for conjunctivitis is configured in a new manner such that it incorporates artificial delay effects to reflect realistic intervention strategies such as quarantine, delayed exposure, or behavioral response to symptoms.

•   Theoretical Advancement: The model’s dynamical properties, including positivity, boundedness, existence, and uniqueness, were studied rigorously. Also, the authors studied local and global stability conditions under biologically meaningful assumptions.

•   Numerical Methodology: A completely new and original stochastic nonstandard finite difference (NSFD) method has been developed and analyzed. This method retains important dynamical properties, ensures numerical positivity and convergence, and beats standard methods under larger time steps.

•   Biological Insight: Artificial delays significantly reduce the number of infected individuals and increase the susceptible population, with use as a control strategy. These findings are compatible with recent public health approaches and have practical relevance to the model.

•   Literature Integration: This paper contributes to filling the gap in existing literature by making an extended union of stochastic dynamics, artificial delays, and conjunctivitis modeling under one united framework, which is justified by the most recent epidemiological events and data.

•   Public Health Interpretation of Artificial Delay: The exponential delay term e−μT incorporated in the model captures realistic intervention strategies employed in controlling disease outbreaks. Specifically, this form of delay accounts for the time lag between exposure and infectiousness, which can be influenced by public health policies such as quarantine protocols, awareness campaigns, and social distancing mandates. For instance, if an infected individual is identified and quarantined within a specific time frame, the probability of disease transmission significantly drops, mimicking the effect of exponential decay in contact or transmission rate. Thus, the model not only reflects the underlying biology but also provides a quantitative framework for evaluating the timing and efficacy of such interventions.

The mathematical framework of this study is deeply rooted in functional analysis. This study employs key results from functional analysis, including fixed point theory, operator norms, and stability criteria in infinite-dimensional spaces, to establish existence, uniqueness, and positivity of solutions. These foundational tools ensure the rigorous formulation and analysis of the stochastic delay differential equations (SDDEs) that govern the proposed model.

The following subdivisions comprise five sections: Section 1 briefly introduces the research disease and literature review. Section 2 presents the derivation of the delayed mathematical model and its dynamical properties, equilibrium points, and local-global behavior of the delayed system. Section 3 describes the stochastic delayed model and the theorem-related unique global positive solution. Section 4 introduces the selective numerical method, and presents the derived stochastic NSFD-related theorems. Lastly, Section 5 discusses numerical simulations among selective techniques and their convergence behavior around the conjunctivitis present state with different step sizes. In addition, the discussion and conclusion are presented in different sections.

2  Derivation of the Delayed Conjunctivitis Mathematical Model

This model is present in [1] in the deterministic form. This study investigates artificial delay effects and stochastic behaviors to understand the realistic transmission dynamics. Delay-based interventions are used to control the spread of disease. Because of the exponential growth of the disease population, the exponential delay is employed to prevent disease transmission similarly. In [1], the human host population N(t) is divided into five compartments Sc(t),Ec(t),Ic(t),IR(t),R(t) at any time t≥0, represent the individual susceptible, individual exposed, individual infected, number of irritants, and individual recovered, respectively. The saturated incidence in the system (1)–(5) with delayed effect represents the rate of change of susceptible individuals to infected individuals and artificial delay e−μT that can be any artificial tactics be used to control the exponential growth of the transmission of the disease with counter artificial exponential-based interventions for all T≥0 (see Fig. 1 and Table 1).

Figure 1: Flow dynamics of the Conjunctivitis disease in a population

The delayed conjunctivitis mathematical model is derived under the following assumptions.

•   Direct and indirect transmission paths are taken.

•   Individuals who recover again become susceptible to the disease.

•   Individuals who lost their vision are considered.

•   The host population is equally mixed.

•   The natural death rate and birth rate are the same.

•   The natural death rate and natural expiry rate of bacteria/viruses are approximately equal.

Although most parameter values are directly adopted from [1], several minor modifications were made to improve computational tractability and biological interpretability. Specifically, this study assumes that the birth and natural death rates are equal (0.05/day), a common simplification in epidemic models to maintain a stable total population in the absence of disease-related deaths. This assumption allows a clearer focus on the disease-induced dynamics without the confounding effects of demographic growth or decline.

The model of the conjunctivitis delay differential equations is as follows:

dSc(t)dt=bh−β1Sc(t)Ic(t−T)e−μT1+m0Ic(t−T)−β2Sc(t)IR(t)1+m1IR(t)−μSc(t)+σR(t)(1)

dEc(t)dt=β1Sc(t)Ic(t−T)e−μT1+m0Ic(t−T)+β2Sc(t)IR(t)1+m1IR(t)−(τ+μ)Ec(t)(2)

dIc(t)dt=τEc(t)−(γ+μ+ξ+δ)Ic(t)(3)

dIR(t)dt=ξIc(t)−ηIR(t)(4)

dR(t)dt=γIc(t)−(μ+σ)R(t)(5)

Given that Sc(0) ≥ 0, Ec(0) ≥ 0, Ic(0) ≥ 0, IR(0) ≥ 0 and R(0) ≥ 0, for all t ≥ 0 are the initial conditions, T<t where T is delay effect.

2.1 Essential Properties

Theorem 1: (Positivity) For any given initial condition (SC(0),Ec(0),Ic(0),IR(0),R(0))ϵR+5 then the solution of system (1)–(5) for all time t>0, (SC(t),Ec(t),Ic(t),IR(t),R(t))ϵR+5.

Proof: If SC(0)≥0,Ec(0)≥0,Ic(0)≥0,IR(0)≥0,R(0)≥0 and the norm ‖ϑ‖∞=SuptϵDϑ|ϑ|. Then, Eq. (1) yields:

dScdt=[bh−β1ScIC1+m0ICe−μT−β2ScIR1+m1IR−μSc+σR],∀ t≥0

dScdt≥−[β1IC1+m0ICe−μT+β2IR1+m1IR+μ]Sc,∀ t≥0

ln⁡Sc≥−[β1‖Ic‖∞1+m0‖Ic‖∞e−μT+β2‖IR‖∞1+m1‖IR‖∞+μ]t+c1,∀ t≥0

Sc≥S(0)e−[β1‖Ic‖∞1+m0‖Ic‖∞e−μT+β2‖IR‖∞1+m1‖IR‖∞+μ]t≥0,∀ t≥0,Sc(t)≥0,∀ t≥0(6)

Similarly, the following relationships from system (2)–(5) can be described:

Ec≥E(0)e−(τ+μ)t≥0,∀ t≥0(7)

Ic≥Ic(0)e−(γ+μ+ξ+δ)t≥0,∀ t≥0(8)

IR≥IR(0)e−ηt≥0,∀ t≥0(9)

R≥R(0)e−(μ+σ)t≥0,∀ t≥0(10)

Eqs. (6)–(10) are non-negative. Hence, the system (1)–(5) is positive ∀ t≥0. □

Theorem 2: (Boundedness) ∀t≥0, the trajectories of the system (1)–(5) with any given initial trajectories from the set Ω persist around the set:

Ω={(SC(t),Ec(t),Ic(t),IR(t),R(t))ϵR+5:Sc(t)+Ec(t)+Ic(t)+IR(t)+R(t)≤bhμ}if limt→∞Sup N(t)≤bhμ and μ=η.

Proof: After adding the system (1)–(5), the total population of the system is as follows:

dNdt=bh−μSc−μEc−μIc−δIc−ηIR−μR

dNdt≤bh−μSc−μEc−μIc−ηIR−μR

dNdt≤bh−μN,provided let μ=η

dNdt+μN≤bh

The general solution

N(t)≤bhμ+(N(0)−bhμ)e−μt(11)

limt→∞Sup N(t)≤bhμ

Now, any given initial trajectory N(0)≤bhμ, then limt→∞Sup N(t)≤bhμ. Hence, Ω is positively invariant. Therefore, system (1)–(5) is locally Lipschitz. □

Theorem 3: (Existence and Uniqueness) The unique solution of system (1)–(5) exists piece-wise if it satisfies linear growth and Lipschitz conditions and μ=η.

Proof: Define the norm ‖ϑ‖∞=SuptϵDϑ|ϑ|. If Li,L¯i and ωi be positive constants ∀i=1,2,3,4,5 for the system (12) for uniqueness and existence, Eqs. (13) and (14) are satisfied, respectively.

Considering system (1)–(5) yields:

{Sc′=F1(Sc,Ec,Ic,IR,R,t)Ec′=F2(Sc,Ec,Ic,IR,R,t)Ic′=F3(Sc,Ec,Ic,IR,R,t)IR′=F4(Sc,Ec,Ic,IR,R,t)R′=F5(Sc,Ec,Ic,IR,R,t),t≥0(12)

∀i=1,2,3,4,5 then

|Fi(S,t)|2<Li(|Si|2+1).(13)

|Fi(S1,t)−Fi(S2,t)|<L¯i|S1−S2|.(14)

The uniqueness of Eq. (12) can be proved by satisfying Eq. (13) as follows:

|F1(Sc,Ec,Ic,IR,R,t)|2=|bh−β1ScIC1+m0ICe−μT−β2ScIR1+m1IR−μSc+σR|2

|F1(Sc,Ec,Ic,IR,R,t)|2≤(|bh|2+|β1|2‖Sc‖∞‖Ic‖∞e−2μT(1+m0‖Ic‖∞)2+|β2|2‖Sc‖∞‖IR‖∞(1+m1‖IR‖∞)2+|μ|2‖Sc‖∞+|σ|2‖R‖∞)

|F1(Sc,Ec,Ic,IR,R,t)|2≤(|β1|2‖Sc‖∞‖Ic‖∞e−2μT(1+m0‖Ic‖∞)2+|β2|2‖Sc‖∞‖IR‖∞(1+m1‖IR‖∞)2+|μ|2‖Sc‖∞)(1+|bh|2+|σ|2‖R‖∞(|β1|2‖Sc‖∞‖Ic‖∞e−2μT(1+m0‖Ic‖∞)2+|β2|2‖Sc‖∞‖IR‖∞(1+m1‖IR‖∞)2+|μ|2‖Sc‖∞))

The condition |bh|2+|σ|2‖R‖∞(|β1|2‖Sc‖∞‖Ic‖∞e−2μT(1+m0‖Ic‖∞)2+|β2|2‖Sc‖∞‖IR‖∞(1+m1‖IR‖∞)2+|μ|2‖Sc‖∞)<1,

|F1(Sc,Ec,Ic,IR,R,t)|2<L1(1+|Sc|2)(15)

Similar,

|F2(Sc,Ec,Ic,IR,R,t)|2=|β1ScIce−μT1+m0Ic+β2ScIR1+m1IR−(τ+μ)Ec|2

|F2(Sc,Ec,Ic,IR,R,t)|2≤(|(τ+μ)|2Ec∞)(1+|β1|2‖Sc‖∞‖Ic‖∞e−2μT(1+m0‖Ic‖∞)2+|β2|2‖Sc‖∞‖IR‖∞(1+m1‖IR‖∞)2(|(τ+μ)|2‖Ec‖∞))

The condition |β1|2‖Sc‖∞‖Ic‖∞e−2μT(1+m0‖Ic‖∞)2+|β2|2‖Sc‖∞‖IR‖∞(1+m1‖IR‖∞)2(|(τ+μ)|2‖Ec‖∞)<1,

|F2(Sc,Ec,Ic,IR,R,t)|2<L2(1+|Ec|2)(16)

|F3(Sc,Ec,Ic,IR,R,t)|2=|τEc−(γ+μ+ξ+δ)IC|2

|F3(Sc,Ec,Ic,IR,R,t)|2≤(|(γ+μ+ξ+δ)|2‖Ic‖∞)(1+|τ|2‖Ec‖∞|(γ+μ+ξ+δ)|2‖Ic‖∞)

The condition |τ|2‖Ec‖∞|(γ+μ+ξ+δ)|2‖Ic‖∞<1,

|F3(Sc,Ec,Ic,IR,R,t)|2<L3(1+|Ic|2)(17)

|F4(Sc,Ec,Ic,IR,R,t)|2=|ξIc−ηIR|2

|F4(Sc,Ec,Ic,IR,R,t)|2≤(|ξ|2‖Ic‖∞+|η|2‖IR‖∞)

|F4(Sc,Ec,Ic,IR,R,t)|2≤|η|2IR∞(1+|ξ|2‖Ic‖∞|η|2‖IR‖∞)

The condition |ξ|2‖Ic‖∞|η|2‖IR‖∞<1,

|F4(Sc,Ec,Ic,IR,R,t)|2<L4(1+|IR|2)(18)

|F5(Sc,Ec,Ic,IR,R,t)|2=|γIC−(μ+σ)R|2

|F5(Sc,Ec,Ic,IR,R,t)|2≤(|γ|2‖Ic‖∞+|(μ+σ)|2‖R‖∞)

|F5(Sc,Ec,Ic,IR,R,t)|2≤|(μ+σ)|2R∞(1+|γ|2‖Ic‖∞|(μ+σ)|2‖R‖∞)

The condition |γ|2‖Ic‖∞|(μ+σ)|2‖R‖∞<1,

|F5(Sc,Ec,Ic,IR,R,t)|2<L5(1+|R|2)(19)

max{|bh|2+|σ|2‖R‖∞(|β1|2‖Sc‖∞‖Ic‖∞e−2μT(1+m0‖Ic‖∞)2+|β2|2‖Sc‖∞‖IR‖∞(1+m1‖IR‖∞)2+|μ|2‖Sc‖∞),|β1|2‖Sc‖∞‖Ic‖∞e−2μT(1+m0‖Ic‖∞)2+|β2|2‖Sc‖∞‖IR‖∞(1+m1‖IR‖∞)2(|(τ+μ)|2‖Ec‖∞),|τ|2‖Ec‖∞|(γ+μ+ξ+δ)|2‖Ic‖∞,|ξ|2‖Ic‖∞|η|2‖IR‖∞,|γ|2‖Ic‖∞|(μ+σ)|2‖R‖∞}<1

Now, Lipschitz’s condition is verified for the existence as follows:

|F1(Sc1,t)−F1(Sc2,t)|=|(bh−β1Sc1IC1+m0ICe−μT−β2Sc1IR1+m1IR−μSc1+σR)−(bh−β1Sc2IC1+m0ICe−μT−β2Sc2IR1+m1IR−μSc2+σR)|

|F1(Sc1,t)−F1(Sc2,t)|≤(β1‖Ic‖∞e−μT1+m0‖Ic‖∞+β2‖IR‖∞1+m1‖IR‖∞+μ+ω1)|Sc1−Sc2|

|F1(Sc1,t)−F1(Sc2t)|≤L¯1|Sc1−Sc2|(20)

|F2(Ec1,t)−F2(Ec2,t)|=|(β1ScIce−μT1+m0Ic+β2ScIR1+m1IR−(τ+μ)Ec1)−(β1ScIce−μT1+m0Ic+β2ScIR1+m1IR−(τ+μ)Ec2)|

|F2(Ec1,t)−F2(Ec2,t)|≤|(τ+μ)+ω2||(Ec1−Ec2)|

|F2(Ec1,t)−F2(Ec2,t)|≤L¯2|(Ec1−Ec2)|(21)

|F3(Ic1,t)−F3(Ic2,t)|=|(τEc−(γ+μ+ξ+δ)Ic1)−(τEc−(γ+μ+ξ+δ)Ic2)|

|F3(Ic1,t)−F3(Ic2,t)|≤|(γ+μ+ξ+δ)+ω3||(Ic1−Ic2)|

|F3(Ic1,t)−F3(Ic2,t)|≤L¯3|(Ic1−Ic2)|(22)

|F4(IR1,t)−F4(IR2,t)|=|(ξIc−ηIR1)−(ξIc−ηIR2)|

|F4(IR1,t)−F4(IR2,t)|≤|η+ω4||(IR1−IR2)|

|F4(IR1,t)−F4(IR2,t)|≤L¯4|(IR1−IR2)|(23)

|F5(R1,t)−F5(R2,t)|=|(γIC−(μ+σ)R1)−(γIC−(μ+σ)R2)|

|F5(R1,t)−F5(R2,t)|≤|(μ+σ)+ω5||(R1−R2)|

|F5(R1,t)−F5(R2,t)|≤L¯5|(R1−R2)|(24)

Eqs. (15)–(19) indicate that systems (1)–(5) have unique solutions and Eqs. (20)–(24) have a sure existence. □

2.2 Model Equilibrium

To obtain the condition where conjunctivitis dies out from the system (1)–(5) or becomes endemic.

Conjunctivitis free equilibrium (CFE)=C1=(Sc1,Ec1,Ic1,IR1,R1)=(bhμ,0,0,0,0)(25)

Conjunctivitis present equilibrium (CPE)=C2=(Sc∗,Ec∗,Ic∗,IR∗,R∗)(26)

Sc∗=(bh+σR∗)(1+m0Ic∗)(1+m1IR∗)β1Ic∗e−μT(1+m1IR∗)+β2IR∗(1+m0Ic∗)+μ(1+m0Ic∗)(1+m1IR∗)

Ec∗=β1Sc∗Ic∗e−μT(1+m1IR∗)+β2Sc∗IR∗(1+m0Ic∗)(1+m0Ic∗)(1+m1IR∗)(τ+μ)

Ic∗=R0μ(τ+μ)Ec∗β1bhe−μT,IR∗=ξIc∗η,R∗=γIc∗(μ+σ)

2.3 Reproduction Number

To determine the system’s reproductive number [22] of (1)–(5) using the following generation matrix method.

Theorem 4: The reproductive number of the system (1)–(5) is given by

R0=τβ1bhe−μTμ(τ+μ)(γ+μ+ξ+δ)(27)

Proof: By the abovementioned method, three concerned compartments are considered, Ec,Ic, and R, under the process hypothesis and overlooked both remaining.

Now,

[Ec′Ic′R′]=[0β1bhe−μTμ0000000][EcIcR]−[τ+μ00−τγ+μ+ξ+δ00−γμ+σ][EcIcR]

here, Fc=[0β1bhe−μTμ0000000]

and Vc=[τ+μ00−τγ+μ+ξ+δ00−γμ+σ],

where Fc= Transmission matrix at C1, Vc= Transition matrix at C1.

Hence, by the largest eigenvalue of FcVc−1,

R0=ρ(FcVc−1)=τβ1bhe−μTμ(τ+μ)(γ+μ+ξ+δ). □

2.4 Stability Analysis

Now, the asymptotic behaviour of convergence for the system (1)–(5) is investigated locally and globally around defined equilibrium states with the constraints of the reproductive number.

Theorem 5: (Local stability at C1) The system (1)–(5) around the Conjunctivitis free equilibrium point is locally asymptotically stable if the reproductive number is less than one.

Proof: To prove this result, the Jacobian matrix of the right-side system (1)–(5) at C1=(bhμ,0,0,0,0) is given by:

J(C1)=[−μ0−β1bhe−μTμ−β2bhμσ0−(τ+μ)β1bhe−μTμβ2bhμ00τ−(γ+μ+ξ+δ)0000ξ−η000γ0−(μ+σ)]

For eigenvalue, |J−λI|=0,

|−μ−λ0−β1bhe−μTμ−β2bhμσ0−(τ+μ)−λβ1bhe−μTμβ2bhμ00τ−(γ+μ+ξ+δ)−λ0000ξ−η−λ000γ0−(μ+σ)−λ|=0(28)

By expanding (28), λ1=−μ and λ2=−(μ+σ). The characteristic equation evaluates the remaining eigenvalues of (28):

λ3+A1λ2+A2λ+A3=0(29)

where A1=−(k1+k2+k3), A2=k1k3+k2k3+k1k2(1−R0), A3=−k1k2k3(1−R0)−τξk1 and k1=−(τ+μ), k2=−(γ+μ+ξ+δ), k3=−η, k4=β1bhe−μTμ, k5=β2bhμ. Since A1A2>A3 and all Ai>0 for i=1,2,3 iff R0<1. Then, Routh-Hurwitz criteria for the 3rd-order polynomial indicate that all the eigenvalues of (29) have negative real parts. Hence, the system (1)–(5) at the Conjunctivitis-free equilibrium point C1 is locally asymptotically stable if R0<1. □

Theorem 6: (Local stability at C2) The system (1)–(5) around the Conjunctivitis present equilibrium point is locally asymptotically stable if the reproductive number exceeds one.

Proof: To prove this result, the Jacobian matrix of system (1)–(5) at C2=(Sc∗,Ec∗,Ic∗,IR∗,R∗) is given by:

J(Sc∗,Ec∗,Ic∗,IR∗,R∗)=[−β1Ic∗e−μT1+m0Ic∗−β2IR∗1+m1IR∗−μ0−β1Sc∗e−μT(1+m0Ic∗)2−β2Sc∗(1+m1IR∗)2σβ1Sc∗e−μT(1+m0Ic∗)2+β2Sc∗(1+m1IR∗)2−(τ+μ)β1Sc∗e−μT(1+m0Ic∗)2β2Sc∗(1+m1IR∗)200τ−(γ+μ+ξ+δ)0000ξ−η000γ0−(μ+σ)](30)