何雪晴，韋煜明

(廣西師范大學(xué) 數(shù)學(xué)與統(tǒng)計學(xué)院，廣西 桂林 541006)

世界正面臨著新型冠狀病毒肺炎(COVID-19)帶來的嚴重威脅[1]，人類感染病毒后常見的體征有發(fā)燒、咳嗽、氣促和呼吸困難等，在較嚴重的病例中，感染可以導(dǎo)致肺炎、嚴重急性呼吸綜合征、腎衰竭，甚至死亡。這引起了許多學(xué)者的注意，他們通過數(shù)學(xué)建模對其進行研究[2-7]。例如，Wang等[5]建立了一個模型來估計武漢市的流行趨勢，對假設(shè)預(yù)防和控制措施足以或不足以控制疫情進行了討論。Hellewell等[6]構(gòu)建了一個傳播模型，發(fā)現(xiàn)在大多數(shù)情況下，高效地對接觸者進行追蹤和對染病者隔離，足以在三個月內(nèi)控制新的新冠肺炎疫情爆發(fā)。Mandal等[7]考慮了新冠肺炎的數(shù)學(xué)模型，其中的人類群體被細分為了5個關(guān)于時間t的類別，即易感類S(t)、暴露類E(t)、隔離類Q(t)、感染類I(t)和恢復(fù)類R(t)，假設(shè)當一個易受感染的人與一個暴露的人接觸時，新冠肺炎傳染病就會傳播。該模型由5個一階常微分方程組成，如下所示:

(1)

其中所有參數(shù)都是正數(shù)。A為易感人群的招募率，β為疾病的傳播率，ρ1(0<ρ1<1)為易感人群中保持適當預(yù)防措施的比例，ρ2(0<ρ2<1)為暴露人群中對疾病傳播采取適當預(yù)防措施的比例，α和b2分別為暴露類中被感染和隔離的比例，b1和c分別為隔離類轉(zhuǎn)移到易感類和感染類的比例，η和v分別為感染類和暴露類的恢復(fù)率，μ為自然死亡率，δ為因新冠肺炎誘發(fā)的死亡率。

但是，考慮到傳染病的傳播具有隨機性，生活中的環(huán)境白噪聲無所不在，所以確定性的COVID-19傳染病模型難以避免地會受到環(huán)境噪聲的影響，這使得模型中涉及的參數(shù)往往會隨著周圍環(huán)境的波動而隨機地在一些平均值附近波動，因此，有必要在確定性的COVID-19傳染病模型(1)中引入隨機波動。研究者通常用Gauss白噪聲來刻畫外部環(huán)境對系統(tǒng)的隨機擾動[8-9]，考慮到除了白噪聲之外，傳染病模型也可能受到彩色噪聲的干擾，導(dǎo)致系統(tǒng)從一種環(huán)境狀態(tài)切換到另一種環(huán)境狀態(tài)[10]，而且在一般情況下，環(huán)境狀態(tài)之間的切換是無記憶的，下一次切換的等待時間遵循指數(shù)分布[11]，所以這種切換可以用連續(xù)時間馬爾可夫鏈{r(t),t≥0}來描述，然后根據(jù)環(huán)境中的溫度、濕度、光照、降雨和風度等不同條件，將其分為N種不同狀態(tài)，其狀態(tài)空間記為={1,2,…,N}。因此，基于以上考慮，為了讓模型更符合實際情況，提出了以下具有Markov切換的隨機COVID-19傳染病模型

(2)

其中σi(i=1,2,3,4，5)用于刻畫白噪聲的強度，Bi(t)(i=1,2,3,4,5)是標準的Brown運動。設(shè)在有限狀態(tài)空間={1,2,…,N}中，Markov鏈 {r(t),t≥0}由轉(zhuǎn)移率矩陣Γ=(γij)N×N生成，即

dx(t)=b(x(t),r(t))dt+σ(x(t),r(t))dB(t),x(0)=x0,r(0)=r，

其中B(·)是d維Brown運動，r(·)是右連續(xù)的Markov 鏈，b(·,·):n×→n，σ(·,·):n×→n×d且滿足σ(x,k)σT(x,k)=(dij(x,k)),k∈。設(shè)V(·,k)是任意二次連續(xù)可微函數(shù)，線性算子L被定義為

由廣義的Ito公式[12]可知，如果x(t)∈n，那么有

dV(x(t),r(t))=LV(x(t),r(t))dt+Vx(x(t),r(t))σ(x(t),r(t))dB(t)。

引理1[13]設(shè)M={Mt}t≥0是一個實值的連續(xù)局部鞅，那么有

引理2[14]對幾乎所有的ω∈Ω，存在正常數(shù)t0=t0(ω)，使得當t≥t0時，有

由引理2可知，集合X是模型(2)解的一個吸引域，下面要討論的內(nèi)容僅在此區(qū)域內(nèi)研究。

1 全局正解的存在唯一性

定理1對任意給定的初始值(S(0),E(0),Q(0),I(0),R(0),r(0))∈X，模型(2)在t≥0時，存在唯一的全局解(S(t),E(t),Q(t),I(t),R(t))，并且該解以概率1位于Γ中。

證明由于模型(2)的系數(shù)滿足局部Lipschitz條件，所以對任意給定的初始值(S(0),E(0),Q(0),I(0),R(0),r(0))∈Χ，模型(2)存在唯一的局部解(S(t),E(t),Q(t),I(t),R(t),r(t))，t∈[0,τe)，其中τe表示爆炸時間。要證隨機模型(2)存在唯一的全局解，只需證τe=+∞，a.s.。

P{τk≤T}≥ε,k≥k1。

(3)

V(S(t),E(t),Q(t),I(t),R(t),r(t))=(S(t)-1-lnS(t))+(E(t)-1-lnE(t))+(Q(t)-1-lnQ(t))+

(I(t)-1-lnI(t))+(R(t)-1-lnR(t))。

因為對任意的x>0，有x-1-lnx≥0，所以V正定。由Ito公式可知

dV(S,E,Q,I,R)=LV(S,E,Q,I,R)dt+σ1r(t)(S-1)dB1(t)+σ2r(t)(E-1)dB2(t)+

σ3r(t)(Q-1)dB3(t)+σ4r(t)(I-1)dB4(t)+σ5r(t)(R-1)dB5(t)。

(4)

≤Ar(t)+βr(t)(1-ρ1r(t))(1-ρ2r(t))E+5μr(t)+b1r(t)+b2r(t)+αr(t)+vr(t)+cr(t)+ηr(t)+δr(t)+

這里的K1是一個正常數(shù)，所以由(4)式整理可得

dV(S,E,Q,I,R)≤K1dt+σ1r(t)(S-1)dB1(t)+σ2r(t)(E-1)dB2(t)+σ3r(t)(Q-1)dB3(t)+

σ4r(t)(I-1)dB4(t)+σ5r(t)(R-1)dB5(t)。

(5)

對(5)式從0到τk∧T積分并取期望得

(6)

V(S(τk,ω),E(τk,ω),Q(τk,ω),I(τk,ω),R(τk,ω))≥k-1-lnk，

則

V(S(0),E(0),Q(0),I(0),R(0))+K1T≥[IΩk(ω)V(S(τk,ω),E(τk,ω),Q(τk,ω),I(τk,ω),R(τk,ω))]

其中IΩk表示Ωk的示性函數(shù)。令k→+∞，則有+∞>V(S(0),E(0),Q(0),I(0),R(0))+K1T=+∞，矛盾，故τe=+∞a.s.，定理1得證。

2 COVID-19傳染病滅絕的充分條件

定理2設(shè)(S(t),E(t),Q(t),I(t),R(t))是初始值為

證明由模型(2)的第二個等式及Ito公式可知

對任意的k∈，根據(jù)引理2，得到

(7)

對(7)式兩邊從0到t積分并同時除以t，得到

(8)

由r(t)的遍歷性可知

(9)

(10)

3 數(shù)值模擬

這一部分將利用Milstein方法[16]以及來自印度的馬哈拉施特拉邦(印度西部邦)和德里的真實新冠肺炎數(shù)據(jù)，說明關(guān)于模型(2)的一個主要理論結(jié)果。

假設(shè)連續(xù)時間Markov鏈{r(t),t≥0}只有兩個環(huán)境狀態(tài)，即={1,2}，平穩(wěn)分布π=(0.5,0.5)，選取初始值為(S(0),E(0),Q(0),I(0),R(0),r(0))=(10,8,5,3,1,1)，結(jié)合表1中Mandal等[17]給出的參數(shù)值，對模型(2)進行兩種情況模擬。

表1 來自印度的各參數(shù)取值

1)如果固定環(huán)境狀態(tài)e∈，則模型(2)的閾值參數(shù)變?yōu)槟敲喘h(huán)境在狀態(tài)1的情況下，取β1=3.5，A1=0.12，σ1=0.15，σ2=0.25，σ3=0.01，σ4=0.02，σ5=0.02，則根據(jù)定理2可知疾病將會滅絕，如圖1所示。

4 結(jié)論及展望