一類隨機SIR傳染病模型的非標準數值離散化

譚偉 劉茂省

DOI：10.16783/j.cnki.nwnuz.2024.01.006

收稿日期：20221108;修改稿收到日期：20230522

基金項目：國家自然科學基金資助項目（12071445，12001501）

作者簡介：譚偉（1996—），男，河北張家口人，碩士研究生.主要研究方向為傳染病動力學.

Email：1848231722@qq.com

*通信聯系人，男，教授，博士.主要研究方向為傳染病動力系統及復雜網絡.Email：liumaoxing@126.com

摘要：研究一個具有飽和發生率和疫苗接種率的隨機離散SIR傳染病模型的穩定性.首先，引入一個具有飽和發生率和疫苗接種率的確定性SIR模型，考慮到隨機噪聲對疾病傳播有很大影響，應用非標準有限差分（NSFD）方法將模型離散化，最終得到一個隨機離散的SIR傳染病模型.這種離散方法是對系統的右側進行局部離散，得出離散模型，然后系統左側用廣義前向差分法對一階導數進行逼近，并且要選取恰當的分母函數.其次，應用Lyapunov函數和矩陣方法給出了系統平衡解穩定的充分條件，并且得到了非線性差分方程概率穩定的充分條件和線性差分方程漸近均方穩定的充分條件.最后，通過數值仿真對理論分析結論進行驗證.

關鍵詞：飽和發生率；隨機離散模型；非標準有限差分方法；Lyapunov函數；漸近均方穩定

中圖分類號：O 175.7??? 文獻標志碼：A??? 文章編號：1001-988Ⅹ（2024）01-0030-09

Nonstandard numerical discretization of

a stochastic SIR epidemic model

TAN Wei1，LIU Mao-xing1，2

（1.College of Mathematics，North University of China，Taiyuan 030051，Shanxi，China;

（2.College of Science，Beijing University of Civil Engineering and Architecture，Beijing 102616，China）

Abstract：This paper studies the stability of a stochastic discrete SIR epidemic model with saturated incidence and vaccination rate.A deterministic SIR model with saturated incidence and vaccination rate is introduced.Considering the great influence of stochastic noise on disease transmission，the model is discretized by nonstandard finite difference（NSFD） method，and finally a stochastic discrete SIR epidemic model is obtained.This discretization method is to locally discretize the right side of the system to obtain the discrete model，and then use the generalized forward difference method to approximate the first derivative on the left side of the system，and select the appropriate denominator function.The sufficient conditions for the stability of the equilibrium solution of the system are given by using Lyapunov function method and matrix method，and the sufficient conditions for the probabilistic stability of nonlinear difference equations and the sufficient conditions for the asymptotic mean square stability of linear difference equations are also proposed.Finally，the conclusion is verified by numerical simulation.

Key words：saturation incidence;stochastic discrete model;nonstandard finite difference method;Lyapunov function;asymptotic mean square stability

0? 引言

近幾十年來，數學模型被認為是理解傳染病傳播機制和評估其控制策略的關鍵工具［1-3］.在這方面，Yusuf等［4］引入并分析了一個具有疫苗接種和治療的SIR流行模型，其中假設發病率是雙線性的.但經過長期的驗證發現，雙線性發生率對于較大的種群規模來說并不符合實際，為了應對易感個體產生的心理或抑制效應，Capasso等［5］提出了飽和發病率βSI/（1+αI），其中β是疾病的傳染率，α是衡量由擁擠效應或社會行為變化引起的抑制效應的飽和因子.當前，飽和發病率已廣泛應用于確定性傳染病模型［6-9］和隨機傳染病模型［10-12］.本文研究具有飽和發生率和疫苗接種率的SIR確定性傳染病模型［13］：

dSdt=Λ-βSI1+αI-（μ+ν）S，

dIdt=βSI1+αI-（μ+γ+δ）I，

dRdt=νS+γI-μR.（1）

其中S≡S（t）為t時刻易感人群數量，I≡I（t）為t時刻已感染人群數量，R≡R（t）為t時刻染病后免疫人群的數量，且S（0）≥0，I（0）≥0，R（0）≥0；Λ為從外部進入種群的凈遷移率，μ為人口的自然死亡率，ν為接種疫苗的人口比率，γ為患病種群恢復的概率，δ為因為感染疾病而死亡的人口概率，且所有參數都是正的.

注意到系統（1）的第一和第二個方程中不包含R，系統（1）的第三個式子不影響其動力學行為，因此，我們只需研究以下系統：

dSdt=Λ-βSI1+αI-（μ+ν）S，

dIdt=βSI1+αI-（μ+γ+δ）I.（2）

通過簡單計算可知，系統（2）中存在一個無病平衡點E0和一個地方病平衡點E+，并且無論在什么條件下，無病平衡點

E0=（S0，I0）=Λμ+ν，0

都存在.當R0>1時地方病平衡點E+=（S+，I+）存在，且

S+=αΛ+（μ+γ+δ）α（μ+ν）+β，

I+=（μ+ν）（R0-1）α（μ+ν）+β，

基本再生數

R0=Λβ（μ+ν）（μ+γ+β）.

由于白噪聲對疾病的傳播有很大的影響，所以與確定性模型相比，隨機模型更加符合實際.本文在模型（2）的基礎上加入隨機噪聲，得到

dSdt=Λ-βSI1+αI-（μ+ν）S+

σ1（S-S*）dB1（t），

dIdt=βSI1+αI-（μ+γ+δ）I+

σ2（I-I*）dB2（t），（3）

其中Bi（t）（i=1，2）為白噪聲，也就是常說的布朗運動，σ2i（i=1，2）為白噪聲Bi（t）的強度［14］.假設隨機擾動強度與狀態S（t）和I（t）偏離平衡點的偏差成正比，即當系統狀態偏離平衡點的偏差增大時，隨機擾動強度也隨之增加［15］.

傳染病模型通常是一個非線性微分方程組，因此很難得到其解析解，通常需要數值近似來得到可靠的結果.然而，許多標準的數值方法，如Euler方法、Runge-Kutta方法可能無法保持連續模型的動力學特性，例如收斂到錯誤的平衡點或數值不穩定［16-18］.為克服這些方法的缺點，Mickens［19］提出了非標準有限差分（NSFD）方法.

如果以下兩個條件之一被滿足，那么一階微分方程組的這種數值格式即為NSFD格式：

（i）用廣義前向差分近似系統中的一階導數

dundt≈u（n+1）-u（n）ψ，

其中un=u（tn），分母函數ψ=ψ（h）滿足ψ（h）=h+O（h2），h為時間步長.

（ii）非線性項是非局部近似的，例如

u2（tn）≈unun+1.

然而，注意到NSFD格式對一個微分方程來說并不是唯一的，所以為了確保NSFD格式與原始連續模型具有相同的動力學性態，我們需要對NSFD格式進行動力學特性的分析.

假定條件（i）成立，我們根據系統（3）提出以下NSFD格式：

S（n+1）-S（n）ψ=

Λ-βS（n）I（n）1+αI（n）-（μ+ν）S（n）+

σ11ψ（S（n）-S*）ω1（n+1），

I（n+1）-I（n）ψ=

βS（n）I（n）1+αI（n）-（μ+γ+δ）I（n）+

σ21ψ（I（n）-I*）ω2（n+1），（4）

其中分母函數ψ是與h相關的函數，初始條件為S（0）=φ1（0），I（0）=φ2（0）.

系統（4）是對系統（3）進行廣義前向差分和右側局部近似得到的，通過一些簡單的代數計算，我們可以很容易地驗證，具有NSFD格式的模型（4）與初始連續模型（2）具有完全相同的平衡點.數學期望我們通常由E來表示，ωi（n+1）（i=1，2）（n∈Z）是一個相互獨立且服從標準正態分布N（0，1）的隨機序列，并且對于ωi（n+1）（i=1，2，n∈Z），滿足

E［ωi（n）］=0， E［ω2i（n）］=1，

E［ωi（n）］≠E［ωj（n）］， i≠j.（5）

通過直接計算，容易驗證分母函數［20］

ψ（h）=1-e-μhμ.（6）

對離散NSFD系統（4）可重新排列，得到其顯式形式：

S（n+1）=S（n）+

Λ-βS（n）I（n）1+αI（n）-（μ+ν）S（n）ψ+

σ1ψ（S（n）-S*）ω1（n+1），

I（n+1）=I（n）+

βS（n）I（n）1+αI（n）-（μ+γ+δ）I（n）ψ+

σ2ψ（I（n）-I*）ω2（n+1）.（7）

近幾十年來，具有噪聲擾動的確定性模型的隨機流行病模型［21-23］已經被廣泛研究，并且越來越多的學者將差分方程模型［24，25］運用到傳染病學的建模中.然而，多數學者都是將隨機模型和差分格式分別建模來研究的［21］，但在傳染病學中對差分方程所描述的隨機離散模型的穩定性的研究很少，所以本文主要對隨機離散的傳染病模型進行研究.

1? 預備知識

對系統（7）的地方病平衡點進行平移變換，令u（n）=S（n）-S+，m（n）=I（n）-I+可得

u（n+1）=u（n）+

Λ-β（u（n）+S+）（m（n）+I+）1+α（m（n）+I+）-

（u+v）（u（n）+S+）ψ+

σ1ψu（n）ω1（n+1），??????? （8）

m（n+1）=m（n）+

β（u（n）+S+）（m（n）+I+）1+α（m（n）+I+）-

（u+γ+δ）（m（n）+I+）ψ+

σ2ψm（n）ω2（n+1）.

系統（8）的零解等價于系統（7）的正解E+=（S+，I+）.

在平衡點E+=（S+，I+）處對系統（7）進行線性化，得到線性近似系統

x（n+1）=1-βI+ψ1+αI+-（μ+v）ψx（n）-

βS+ψ（1+αI+）2y（n）+

σ1ψx（n）ω1（n+1），

y（n+1）=βI+ψ1+αI+x（n）+????? （9）

1+βS+ψ（1+αI+）2-

（μ+γ+δ）ψy（n）+

σ2ψy（n）ω2（n+1）.

類似地進行平移變換.對無病平衡點，令u（n）=S（n）-Λ/（μ+v），v（n）=I（n），則

u（n+1）=u（n）+

-βu（n）+Λμ+vm（n）1+αm（n）-

（u+v）u（n）ψ+

σ1ψu（n）ω1（n+1），??????? （10）

m（n+1）=m（n）+

βu（n）+Λμ+vm（n）1+αm（n）-

（u+γ+δ）m（n）ψ+σ2ψm（n）ω2（n+1）.

在無病平衡點處對系統（7）線性化，得到相應的線性近似系統為

x（n+1）=（1-（μ+v）ψ）x（n）-Λβψμy（n）+

σ1ψx（n）ω1（n+1），

y（n+1）=1+Λβψμ-（μ+γ+δ）ψy（n）+

σ2ψy（n）ω2（n+1）.（11）

令ρ（n）=（u（n），m（n））T，z（n）=（x（n），y（n））T，φ（n）=（φ1（n），φ2（n））T，T代表轉置.

下面引入了文獻［26］的一些定義和定理來研究系統的動態行為.

定義1（文獻［26］定義7.1）? 對初始函數S（0）=φ1（0），I（0）=φ2（0），如果對ε1>0和ε2>0，存在一個δ>0使得系統（8）或（10）的解ρ（n）=ρ（n，φ）滿足不等式Psupn∈Zρ（n）>ε1<ε2和P{φ<δ}=1，那么稱系統（8）或（10）是依概率穩定的.

定義2（文獻［26］定義1.2）? 如果對每個ε>0，存在δ（ε）>0，使得E［z（n）2］<ε，n∈Z，并且E［φ（n）2］<δ對任意初始條件都被滿足，那么線性系統（9）或（11）的零解稱為均方穩定的；如果系統（9）或（11）是均方穩定的且其解滿足limn→∞E［z（n）2］=0，則稱它們是漸近均方穩定的.

對非負函數Vi=V（i，z（0），z（1），…，z（i）），i∈Z，定義

ΔVi=V（i+1，z（0），z（1），…，z（i+1））-

V（i，z（0），z（1），…，z（i））.

定理A（文獻［26］定理1.1）? 對于線性系統（9）或（11），若存在一個滿足條件E［V（0，φ）］≤c1φ2和E［ΔVi］≤-c2E［z（i）2］（i∈Z）的非負函數Vi=V（i，z（0），z（1），…，z（i）），其中c1和c2是正常數，那么系統（9）或（11）的零解是漸近均方穩定的.

考慮以下系統：

x（n+1）=a11x（n）+a12y（n）+

σ1ψx（n）ω1（n+1），

y（n+1）=a21x（n）+a22y（n）+

σ2ψx（n）ω2（n+1），（12）

其中σ1，σ2是隨機項中的常數，ψ表示分母函數，ωi（n+1）（i=1，2）是相互獨立的隨機序列且滿足條件（5）.可以看出，系統（12）為系統（9）和（11）的一般形式，因而，要探究系統（9）和（11）解的穩定性，只需研究系統（12）解的穩定性.

令

A=a11a12a21a22，

D=d11d12d12d22，

P=p11p12p12p22，（13）

其中D和P是對稱矩陣.對于實對稱矩陣P和Q，如果P-Q是一個正定矩陣，那么P>Q.

定理B（文獻［26］定理5.1）? 假設存在正定矩陣P，使得形為（13）的半正定矩陣D滿足矩陣方程ATDA-D=-P，同時P要滿足條件

P>d11σ21ψ00d22σ22ψ，（14）

那么系統（12）的零解是漸近均方穩定的.

證明? 由（13）式，系統（12）可表示為

z（n+1）=（A+θ（ω（n+1）））z（n），

其中

z（n）=x（n）y（n），

θ（ω（n+1））=

σ1ψω1（n+1）00σ2ψω2（n+1）.

考慮Lyapunov函數V（n）=zT（n）Dz（n），則

ΔV（n）=V（n+1）-V（n）=

zT（n+1）Dz（n+1）-zT（n）Dz（n）.

然后計算ΔV的期望，得到

E［ΔV（n）］=E［zT（n+1）Dz（n+1）-

zT（n）Dz（n）］=

E［zT（n）］·［（A+θ（ω（n+1）））T×

D（A+θ（ω（n+1）））-D］z（n）=

E［zT（n）］［ATDA-D+θT（ω（n+1））×

Dθ（ω（n+1））］z（n）=

EzT（n）［-P+θT（ω（n+1））×

Dθ（ω（n+1））］z（n）.

由關系式（5），可以求得

E［θT（ωω（+1）Dθ）Dθ（+1）］=

Eσ1ψw1（n+1）00σ2ψw2（n+1）×

d11d12d12d22×

σ1ψw1（n+1）00σ2ψw2（n+1）=

d11σ21ψ00d22σ22ψ，

然后根據（14）式，可以得到

E［ΔV（n）］=EzT（n）-P+

d11σ21ψ00d22σ22ψz（n）≤

-cE［z（n）2］，

其中c是正數.依據定理A可知系統（12）的零解是漸近均方穩定的.? 】

2? 地方病平衡點的穩定性

因為系統（9）是系統（12）的特殊形式，可以推導出系統（9）的系數矩陣

A=A11A12A21A22，（15）

其中

A11=1-βI+ψ1+αI+-（μ+v）ψ，

A12=-βS+ψ（1+αI+）2，

A21=βI+ψ1+αI+，

A22=1+βS+ψ（1+αI+）2-（μ+γ+δ）ψ.

從矩陣方程ATDA-D=-P可解出

-p11=（a211-1）d11+2a11a21d12+a221d22，

-p12=a11a12d11+（a21a12+a11a22-1）d12+

a21a22d22，

-p22=a212d11+2a12a22d12+（a222-1）d22.（16）

將矩陣A的元素帶入（16）式可得

p11=-1-βI+ψ1+αI+-（μ+v）ψ2-1d11-

2βI+h1-βI+ψ1+αI+-（μ+v）ψ1+αI+ψd12-

βI+ψ1+αI+2d22，

p12=1-βI+ψ1+αI+-（μ+v）ψ×

βS+ψ（1+αI+）2d11+

βI+ψ1+αI+·βS+ψ（1+αI+）2d12-

-1-βI+ψ1+αI+-（μ+v）ψ×

1-βS+ψ（1+αI+）2-（μ+γ+δ）ψ-1d12-

βI+ψ1+αI+1+βS+ψ（1+αI+）2-

（μ+γ+δ）ψd22，

p22=-βS+ψ（1+αI+）22d11+2βS+ψ（1+αI+）2×

1+βS+ψ（1+αI+）2-（μ+γ+δ）ψd12-

1+βS+ψ（1+αI+）2-（μ+γ+δ）ψ2-1d22.（17）

再由（14）式可導出

p11-d11σ21ψ>0， p11p22-p212>0，

（p11-d11σ21ψ）（p22-d22σ22ψ）-p212>0.（18）

根據以上的分析我們得出以下推論.

推論1? 若存在一個正定矩陣P，使得條件（18）對于矩陣方程（17）的解（d11，d12，d22）被滿足，并且矩陣D是半正定的，那么系統（9）的零解是漸近均方穩定的.因此，系統（8）的零解是依概率穩定的，這等同于系統（7）的正平衡點Ee是依概率穩定的.

當R0>1時，系統（7）存在地方病平衡點E+=（S+，I+）.對于數值仿真，假設初值

S（0）=5， I（0）=2.（19）

其它參數設置為

Λ=4，α=0.01，β=0.5，μ=0.5，

γ=0.1， ν=0.2，δ=0.1，h=0.1.（20）

通過這些參數能夠求出分母函數為ψ=0.0976，基本再生數R0=4.0816>1，從而地方病平衡點存在，且E+=（1.4597，4.2629）.通過對（15）求解得到矩陣

A=0.7322-0.06550.19950.9972；

假設P為單位矩陣，即

P=1001，

求解（17）式得到半正定矩陣

D=7.86907.38447.384412.3981.

由推論1可知，若σ1，σ2滿足條件

σ1<1d11ψ=1.1410，

σ2<1d22ψ=0.9091，（21）

那么，系統（9）的零解是漸近均方穩定的.因此，系統（8）的零解是依概率穩定的，這等同于系統（7）的地方病平衡點E+是依概率穩定的.

下面借助數值模擬證明之前的分析.（19）和（20）式給出了初值和一些參數值.圖1顯示了噪聲擾動強度σ1=σ2=0的情況，系統（9）和線性系統（7）的解曲線分別由圖1（a）和圖1（b）表示.從中可以看出，線性系統（9）的解曲線最終收斂到0，系統（7）的解曲線最終收斂于地方病平衡點E+=

（1.4597，4.2629）.因此，系統（9）的零解是漸近均方穩定的，系統（7）的地方病平衡點E+是依概率穩定的.

對具有隨機擾動的系統（7）的解進行模擬.（19）和（20）式已給出初值和其它的參數值.假設σ1，σ2的值滿足條件（21），我們取σ1=σ2=0.8，得到具有噪聲擾動的系統（7）的解軌跡，見圖2.

從圖2可以看出，隨機噪聲并不影響曲線的最終走向，系統（7）的解軌跡最終還是收斂到地方病平衡點E+，這就表明系統（7）在正平衡點處是依概率穩定的.

圖2? σ1=σ2=0.8時系統（7）的解軌跡

Fig 2Solution trajectories of system （7） when σ1=σ2=0.8

3? 無病平衡點的穩定性

由系統（11）可以寫出它的系數矩陣

A=1-（μ+v）ψ-Λβψμ

01+Λβψμ-（μ+γ+δ）ψ.（22）

同理，由方程（16）解得

p11=-（（1-（μ+ν）ψ）2-1）d11，

p12=（1-（μ+ν）ψ）·Λβψμd11-

（1-（μ+ν）ψ）1+Λβψμ-

（μ+γ+δ）ψ-1d12，

p22=-Λβψμ2d11+

2Λβψμ1+Λβψμ-（μ+γ+δ）ψd12-

1+Λβψμ-（μ+γ+δ）ψ2-1d22.（23）

再由（18）式，我們同樣能推導出以下結論.

推論2? 對于（22）式中的矩陣A，如果存在一個正定矩陣P，使條件（18）被方程（23）的解（d11，d12，d22）所滿足，且矩陣D為半正定的，那么稱系統（11）的零解是漸近均方穩定的.因此，系統（10）的零解是依概率穩定的，這等同于是依概率穩定的系統（7）的邊界平衡點E0.

當R0<1時，系統（7）不存在地方病平衡點，只有一個無病平衡點.取參數值

Λ=4，μ=0.01，γ=0.9，β=0.01，

δ=0.2，ν=0.1，h=0.1.（24）

此時R0=0.167<1，所以只有一個無病平衡點存在，且E0=（20，0）.借助（6）式求得分母函數ψ=0.0995，給定初值為S（0）=10，I（0）=5.由以上參數能夠求出（22）式中的矩陣

A=0.9801-0.039800.9204.

仍取P為單位矩陣，借助（23）式求得半正定矩陣

D=25.3807-10.1124-10.112411.4208.

并由推論1可知，若σ1，σ2滿足條件：

σ1<1d11ψ=0.6293，

σ2<1d22ψ=0.9381，（25）

那么系統（11）的零解是漸近均方穩定的，因此，是依概率穩定的系統（10）的零解，這等同于是依概率穩定的系統（7）的無病平衡點.

根據（24）中的參數值和給定的初值，當隨機擾動強度為0時，我們得到系統（7）和（11）的解曲線，見圖3.從中可以看出，線性系統（11）和系統（7）的解曲線最終分別收斂于它們的平衡點.因此，通過圖像，我們也可以得到系統（11）的零解是漸近均方穩定的，而系統（7）的無病平衡點E0是依概率穩定的.

當隨機噪聲強度不為0，且σ1和σ2的值滿足條件（25）時，我們取σ1=σ2=0.5可以得到具有隨機噪聲擾動的系統（7）的解軌跡，見圖4.從圖4可以直觀地看到，白噪聲只在開始的一段時間內對疾病傳播具有較小的影響，但系統（7）的最終走向仍舊收斂到邊界平衡點E0=（20，0），這也表明系統（7）的邊界平衡點在噪聲擾動下依然是依概率穩定的.

when σ1=σ2=0

when σ1=σ2=0.5

當R0>1時，考慮系統（7）中地方病平衡點存在時無病平衡點的穩定性.仍?。?9）和（20）式中的初值和參數，此時求得（22）式中的矩陣

A=0.9317-0.390401.3221.

P仍為單位矩陣，推導出矩陣

D=7.511511.786911.786913.4013，

通過計算可知，矩陣D不是半正定的，所以根據推論1可知，系統（7）的無病平衡點不穩定.

當隨機噪聲強度為0時，圖5中代表染病者的曲線I趨向于無窮，所以此時系統（11）的零解是不穩定的.當加入一個較小的噪聲強度后，系統的解軌跡是非常不規則的（圖6）.因此，當一個正平衡點存在時，無論是否存在隨機擾動，系統（7）的邊界平衡點都是不穩定的.

4? 結論

本文討論了一個具有飽和發生率和疫苗接種率的連續SIR傳染病模型，鑒于疾病在傳播過程中會受到各種外部因素的影響，例如隨機白噪聲的影響，我們在系統中加入隨機擾動項.由于離散模型的動力學性態比連續模型表現的更加豐富和準確，我們還構造了保持連續模型基本性質的非標準有限差分格式，并選擇了一個合適的分母函數.隨后用Lyapunov函數方法探究了所構建的隨機差分模型的動力學性質，建立了模型在平衡點處穩定的充分條件.數值模擬表明，當所選參數和噪聲強度滿足給定條件時，適當的噪聲強度不會影響系統在平衡點處的穩定性；并且當正平衡點存在時，無論隨機擾動強度是多少，邊界平衡點的穩定性都不會被改變，在這種情況下，邊界平衡點是不穩定的.

參考文獻：

［1］? CASTILLOCHAVEZ C.Mathematical Models in Population Biology and Epidemiology［M］.New York：Springer，2012.

［2］? BRAUER F.Mathematical epidemiology：past，present，and future［J］.Infectious Disease Modelling，2017（2）：113.

［3］? REZAPOUR S，MOHAMMADI H，SAMEI M E.SEIR epidemic model for COVID-19 transmission by Caputo derivative of fractional order［J］.Advances in Difference Equations，2020，2020（1）：490.

［4］? YUSUF T T，BENYAH F.Optimal control of vaccination and treatment for an SIR epidemiological model［J］.World Journal of Modelling & Simulation，2012，8（3）：194.

［5］? CAPASSO V，SERIO G.A generalization of the Kermack-Mc Kendrick deterministic epidemic model［J］.Mathematical Biosciences，1978，42（1/2）：43.

［6］? XU R，MA Z，WANG Z.Global stability of a delayed SIRS epidemic model with saturation incidence and temporary immunity［J］.Computers & Mathematics with Applications，2010，59（9）：3211.

［7］? JANA S，NANDI S K，KAR T K.Complex dynamics of an SIR epidemic model with saturated incidence rate and treatment［J］.Acta Biotheoretica，2016，64：65.

［8］? MATHUR K S，NARAYAN P.Dynamics of an epidemic model with vaccination and saturated incidence rate［J］.Nature Public Health Emergency Collection，2018，4（5）：1.

［9］? LIU J.Bifurcation analysis for a delayed SEIR epidemic model with saturated incidence and saturated treatment function［J］.Journal of Biological Dynamics，2019，13（1）：461.

［10］? RAJASEKAR S P，PITCHAIMANI M，ZHU Q.Dynamic threshold probe of stochastic SIR model with saturated incidence rate and saturated treatment function［J］.Physica A：Statistical Mechanics and its Applications，2019，535：122300.

［11］? RAJASEKAR S P，PITCHAIMANI M，ZHU Q.Progressive dynamics of a stochastic epidemic model with logistic growth and saturated treatment［J］.Physica A：Statistical Mechanics and its Applications，2019，536：122649.

［12］? LI D，ZHAO Y，SONG S.Dynamic analysis of stochastic virus infection model with delay effect［J］.Physica A：Statistical Mechanics and its Applications，2019，528：121463.

［13］? SURYANTO A，DARTI I.On the nonstandard numerical discretization of SIR epidemic model with a saturated incidence rate and vaccination［J］.AIMS Mathematics，2020，6（1）：141.

［14］? MAO X.Stochastic Differential Equations and Applications［M］.New York：Academic Press，2006.

［15］? BERETTA E，KOLMANOVSKII V，SHAIKHET L.Stability of epidemic model with time delays influenced by stochastic perturbations［J］.Mathematics and Computers in Simulation，1998，45（3/4）：269.

［16］? SURYANTO A.A Dynamically consistent nonstandard numerical scheme for epidemic model with saturated incidence rate［J］.Int J Math Comput，2011，13：112.

［17］? SURYANTO A.Stability and bifurcation of a discrete SIS epidemic model with delay［C］//Proceedings of the 2nd International Conference on Basic Sciences，Malang，Indonesia，2012：1.

［18］? HU Z，TENG Z，JIANG H.Stability analysis in a class of discrete SIRS epidemic models［J］.Nonlinear Analysis：Real World Applications，2012，13（5）：2017.

［19］? MICKENS R E.Nonstandard Finite Difference Models of Differential Equations［M］.Singapore：World Scientific，1993.

［20］? MICKENS R E.Calculation of denominator functions for nonstandard finite difference schemes for differential equations satisfying a positivity condition［J］.Numerical Methods for Partial Differential Equations，2010，23（3）：672.

［21］? LIN Y，JIANG D.Long-time behaviour of a perturbed SIR model by white noise［J］.Discrete & Continuous Dynamical Systems（Series B），2013，18（7）：1873.

［22］? GRAY A，GREENHALGH D，HU L，et al.A stochastic differential equation SIS epidemic model［J］.SIAM Journal on Applied Mathematics，2011，71（3）：876.

［23］? DANG H N，TRAN K，DIEU N T，et al.Extinction and permanence in a stochastic SIRS model in regime-switching with general incidence rate［J］.Nonlinear Analysis：Hybrid Systems，2019，34：121.

［24］? SHAIKHET L.Stability of equilibrium states for a stochastically perturbed exponential type system of difference equations［J］.Journal of Computational & Applied Mathematics，2015，290：92.

［25］? PALMER P.Application of a discrete It formula to determine stability（instability） of the equilibrium of a scalar linear stochastic difference equation［J］.Computers & Mathematics with Applications，2012，64：2302.

［26］? SHAIKHET L.Lyapunov Functionals and Stability of Stochastic Difference Equations［M］.London：Springer Science & Business Media，2011.

（責任編輯? 馬宇鴻）