一類具有信息干預(yù)的SEPAIRM模型的動(dòng)力學(xué)分析與最優(yōu)控制

黃懷玉,向會(huì)立

(湖北民族大學(xué) 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,湖北 恩施 445000)

傳染病一直是人類生命安全和生態(tài)安全的巨大威脅,如何預(yù)防和控制傳染病的傳播一直是人們重要研究方向之一。建立貼近現(xiàn)實(shí)的數(shù)學(xué)模型對(duì)于更有效地預(yù)防和控制傳染病具有重大的意義。如在現(xiàn)實(shí)生活中許多疾病都存在無(wú)癥狀感染[1-2],而忽略了無(wú)癥狀傳播的影響將造成不必要的損失。Duong等[1]對(duì)登革熱病毒傳播的研究表明無(wú)癥狀感染者對(duì)蚊子的傳染性明顯高于有臨床癥狀的患者,研究結(jié)果從根本上改變了當(dāng)時(shí)對(duì)于登革熱流行病學(xué)的認(rèn)識(shí)和控制的模式。Hao等[3]通過(guò)病毒樣本模擬了COVID-19的傳播,推斷具有大量的無(wú)癥狀和癥狀前傳播,建立一類具有癥狀前感染和無(wú)癥狀感染的傳染病模型。在對(duì)傳染病傳播行為有一定的基礎(chǔ)之后,就需要采取一定的方法預(yù)防或者控制疾病的傳播,而最優(yōu)控制是重要的方法之一。Abbasi[2]提出了一個(gè)無(wú)癥狀感染區(qū)和隔離區(qū)的傳染病模型的最優(yōu)控制問(wèn)題研究并應(yīng)用到了新冠肺炎的研究。Saha等[4]提出了一個(gè)新的新冠模型的最優(yōu)控制問(wèn)題,即維持社會(huì)分化和疫苗接種程序隨時(shí)間的變化,結(jié)果表明,只有接種疫苗才能顯著降低整體感染人群。值得注意的是,媒體報(bào)道一直是預(yù)防和控制傳染病的重要措施[5-8]。Rai[7]提出了一個(gè)數(shù)學(xué)模型來(lái)評(píng)估社交媒體廣告對(duì)印度抗擊病毒大流行的影響,他們的結(jié)果表明媒體覆蓋是一種很重要的手段,為了徹底消滅疾病,應(yīng)該提高認(rèn)識(shí)水平,努力改變易感人群的行為。Buonomo 與Marca[8]提出了一種包括所有年齡段有效疫苗接種的行為流行病模型,疫苗接種依賴于信息。他們的結(jié)果表明地方性平衡點(diǎn)在R0>1與R0<1時(shí)都可能發(fā)生Hopf分叉。現(xiàn)實(shí)中疫苗是預(yù)防傳染病傳播最主要的措施之一,而加強(qiáng)人群意識(shí),促使易感人群盡快接種疫苗也與媒體覆蓋的強(qiáng)度有關(guān)。本文以受媒體信息影響的疫苗接種作為一種預(yù)防措施,通過(guò)加強(qiáng)媒體覆蓋的強(qiáng)度,進(jìn)而提高接種率,從而達(dá)到控制疾病傳播的目的。本文在文獻(xiàn)[3]的基礎(chǔ)上考慮了媒體覆蓋的因素,建立了如下模型:

(1)

其中:S, E, P, A, I, R, M分別是易感者,暴露者,癥狀前感染者,無(wú)癥狀感染者,癥狀感染者,移除者,媒體信息水平。Λ是招募率,w是無(wú)癥狀感染者和癥狀前感染者的傳播率,b是癥狀傳播者的傳播率,μ1是移除率,在某一時(shí)刻,暴露者轉(zhuǎn)移到癥狀前感染者的比率是de,癥狀前感染者到無(wú)癥狀感染者的比率是(1-r)dp, 癥狀前感染者到癥狀感染者的比率是rdp,無(wú)癥狀感染者,癥狀感染者到移除者的比率分別是da,di。δ代表媒體對(duì)于感染總數(shù)的反應(yīng)率,μ2代表媒體傳播中的信息損失率。參數(shù)Λ,w,b,r,de,dp,da,di,μ1,μ2,δ都是大于0的常數(shù)。

系統(tǒng)(1)初值條件為:

S(0)>0,E(0)>0,P(0)>0,A(0)>0,I(0)>0,R(0)>0,M(0)>0。

1 解的非負(fù)性與有界性

本節(jié)我們將驗(yàn)證系統(tǒng)(1)解的非負(fù)性和有界性。

定理1 系統(tǒng)(1)的解是非負(fù)的。

證明由系統(tǒng)(1)的第一個(gè)式子可得

(2)

其中

即

(3)

對(duì)(3)兩邊從0到t積分可得

(4)

由于S(0)>0, 所以S(t)非負(fù)。同理可以證明E(t),P(t),A(t),I(t),R(t),M(t)都是非負(fù)的。

定理2 系統(tǒng) (1) 的解有界。

證明由系統(tǒng) (1) 第一個(gè)和第二個(gè)式子可得

(5)

故

(6)

由定理1可得,S(t),E(t)非負(fù), 所以我們得出S(t),E(t)有界,同理我們得出P(t),A(t),I(t),R(t),M(t)有界。

2 基本再生數(shù)與穩(wěn)定點(diǎn)

接下來(lái)為了研究無(wú)病平衡點(diǎn)和地方性平衡點(diǎn)的穩(wěn)定性,本節(jié)我們將首先證明無(wú)病平衡點(diǎn)和地方性平衡點(diǎn)的存在性并給出基本再生數(shù)的表達(dá)式。

2.1 無(wú)病平衡點(diǎn)

引理1 系統(tǒng)(1)存在唯一的無(wú)病平衡點(diǎn)

證明令系統(tǒng)(1)右端為0且I=0,可得

2.2 基本再生數(shù)

基本再生數(shù)是指在易感者人群中,一個(gè)感染者所可以傳染的人數(shù),它是判斷疾病是否持續(xù)存在的關(guān)鍵閾值。根據(jù)文獻(xiàn)[9],基本再生數(shù)R0計(jì)算公式如下:

(7)

引理2 當(dāng)R0>1時(shí),系統(tǒng)(1)存在唯一的地方性平衡點(diǎn)X*。

證明令系統(tǒng)(1)右端為0,I≠0,計(jì)算可得

(8)

(9)

將這些公式帶入系統(tǒng)(1)中可得

令

(10)

顯然G(P)關(guān)于P是一個(gè)增函數(shù),而當(dāng)

R0>1時(shí),

(11)

則由介值定理可得存在P*>0, 使得

G(P*)=0,

(12)

所以,系統(tǒng)(1)存在唯一的地方性均衡點(diǎn)X*=(S*,E*,P*,A*,I*,R*, M*), 其中

3 穩(wěn)定性分析

本節(jié)我們將通過(guò)基本再生數(shù)R0研究無(wú)病平衡點(diǎn)和地方性平衡點(diǎn)的漸進(jìn)穩(wěn)定性。

3.1 無(wú)病平衡點(diǎn)的穩(wěn)定性

定理3 當(dāng)R0<1時(shí),無(wú)病平衡點(diǎn)X0是局部漸進(jìn)穩(wěn)定的。

證明非線性系統(tǒng)(1)在無(wú)病平衡點(diǎn)X0處的雅可比矩陣為

(13)

其中

g=wS0,d=bS0。

故J(X0)對(duì)應(yīng)的特征方程為

(α+μ1)2(α+μ2)T(α)=0,

(14)

其中

T(α)=(α+de+μ1)(α+dp+μ1)(α+da+μ1)(α+di+μ1)

gde(α+di+μ1)(α+da+μ1)

-g(1-r)dedp(α+di+μ1)

-drdedp(α+da+μ1),

顯然α1=α2=-μ1,α7=-μ2是特征方程(14)的特征根且都具有負(fù)實(shí)部。下面我們將證明剩余的根具有負(fù)實(shí)部。不妨假設(shè)α=x+iy(x,y∈R)是特征方程(14)的一個(gè)特征根且x>0,則

=R0,

這與R0<1矛盾,故特征方程(14)的每個(gè)根都具有負(fù)實(shí)部。因此當(dāng)R0<1時(shí),無(wú)病平衡點(diǎn)X0是局部漸進(jìn)穩(wěn)定的。

3.2 地方性平衡點(diǎn)的穩(wěn)定性

定理4當(dāng)R0>1時(shí),若滿足

則地方性平衡點(diǎn)X*具有局部漸進(jìn)穩(wěn)定性,其中取

n=

a=(da+μ1),p=(dp+μ1),e=(de+μ1),i=(di+μ1),m=(l+μ1),

A1=m+p+a+i+μ2+e,

A2=l(p+a+i+μ2+e)+p(a+i+μ2+e)+a(i+μ2+e)+i(μ2+e)+μ2e,

A3=eμ2(i+a)+ai(e+μ2)+p(μ2e+ie+iμ2+ae+aμ2+ai)+l(eμ2+ie+iμ2+pa)+

l(ae+aμ2+ai+pμ2+pa)+

g*de(l+μ2+μ1+i+a)+deq*(l+μ2+μ1+i+p)+deq*(1-r)dp+rdpdeq*+ndeδ,

A4=eiμ2(a+p)+paμ2(e+i)+μ2l(ei+ea+ai)+pl(μ2e+μ2i+ia+aμ2+ae)+

deg*(l(i+μ2+a)+μ2(i+a)+ia)+qdedpdp(1-r)(l+μ2+i)-q*del(μ2+a+i)-

nδdel+ndeδ(dp-2rdp+i+a+l)+r2ded*dp(l+μ2+a),

A5=μ2ei(ap+la+lp)+lpa(eμ2+ie+iμ2)-nδ(i+a)del+

deg*(lμ2(i+a)+il(a+μ2)+al(μ2+i))+q*dpde(1-r)(l(μ2+i)+μ2i)+

rdpded*(l(μ2+a)+μ2a)-g*del(μ2(i+a)+ia)-nrplδd*de+

ndeδ[dp(1-r)(l+i)+rdp(a+l)+ai+il+al]-n(1-r)dpdelδ,

A6=aielpμ2+iaμ2delg*+il(1-r)dpdeg*μ2+raμ2lded*dp-ialg*deμ2-

nralδdpde-nilδ(1-r)dpde-niaδlde+ndeδ[(1-r)dpil-rdpam+ial],

證明非線性系統(tǒng)(1)在地方性平衡點(diǎn)X*的雅可比矩陣為

(15)

則J(X*)對(duì)應(yīng)的特征方程為

(α+μ1)(α6+A1α5+A2α4+A3α3+A4α2+A5α+A6)=0,

(16)

顯然α1=-μ1是特征方程(16)的一個(gè)特征根且具有負(fù)實(shí)部。經(jīng)過(guò)計(jì)算可得A1>0,A2>0,A3>0,A4>0,A5>0,A6>0。由Routh-Hurwitz準(zhǔn)則知,如果R0>1,且滿足

A1A2>A3,

則地方性平衡點(diǎn)X*是局部漸進(jìn)穩(wěn)定的。

4 最優(yōu)控制

為了有效控制疾病的傳播,本文考慮易感人群受媒體信息依賴的疫苗接種,以及通過(guò)治療來(lái)降低傳播者人數(shù)來(lái)控制病毒傳播,在系統(tǒng)(1)中引入兩個(gè)控制變量,u(t)代表媒體信息報(bào)道的強(qiáng)度,v(t)代表對(duì)感染人群的治療強(qiáng)度,建立了如下的最優(yōu)控制系統(tǒng):

(17)

其中:p0代表接種疫苗后基本接種率的正常數(shù)。

允許控制集為:

Uad={u,v|0≤u≤1,0≤v≤1,u,v是可測(cè)函數(shù)}

在對(duì)傳染病的預(yù)防控制中,人們希望盡可能降低因?yàn)橐呙缃臃N,感染治療與通過(guò)媒體覆蓋而產(chǎn)生的經(jīng)濟(jì)與時(shí)間成本。為此,定義最優(yōu)目標(biāo)泛函為:

(18)

式中:k1,k2,n1,n2是相應(yīng)的權(quán)重系數(shù)。

4.1 最優(yōu)控制的存在性

在這一節(jié)中,我們將首先證明最優(yōu)控制的存在唯一性。

定理5 對(duì)系統(tǒng)(17),存在唯一的最優(yōu)控制對(duì)(u*,v*)和所對(duì)應(yīng)的(S*,E*,P*,A*,I*,R*,M*),使得

(19)

證明為了證明最優(yōu)控制的存在唯一性,應(yīng)用文獻(xiàn)[10]的結(jié)果。主要分為以下幾步:

1)允許狀態(tài)集非空;

2)允許控制集Uad是凸且閉集;

3)系統(tǒng)(17)的右側(cè)被狀態(tài)和控制變量的線性函數(shù)所限定;

4)被積函數(shù)

是凸函數(shù)。

5)存在常數(shù)η1>0.η2>0,η3>0 使得被積函數(shù)L滿足

現(xiàn)在來(lái)逐個(gè)驗(yàn)證以上5項(xiàng):

1)由定理2方法可得,系統(tǒng)(17)的右端有界, 進(jìn)一步由系統(tǒng)(17)右端是連續(xù)的,所以根據(jù)文獻(xiàn)[11]得,系統(tǒng)(17)至少存在一個(gè)解,因此允許狀態(tài)集非空。

2)任意U1,U2∈Uad及λ∈[0,1],其中U1=(u1,v1),U2=(u2,v2),取U3=λU1+(1-λ)U2,則

U3=(u3,v3)=(λu1+(1-λ)u2,λv1+(1-λ)v2)

由0≤ui≤1且0≤vi≤1(i=1,2),故u3∈[0,1],v3∈[0,1],故U3∈Uad,即U是凸集。另一方面,若Un=(un,vn)∈Uad(n∈N+),且

由0≤un≤1,0≤vn≤1,可得0≤u≤1,0≤v≤1且u,v是可測(cè)函數(shù),故Uad是閉集。

3)由0≤u≤1, 0≤v≤1,從而通過(guò)系統(tǒng)(17),我們可知

(20)

其中:K1,K2的值由系統(tǒng)(17)決定, 因此(3)成立。

4)顯然被積函數(shù)是凸集,在此便不過(guò)多證明。

5)取

則

4.2 最優(yōu)控制的特征

在這一部分,我們將研究通過(guò)Pontryagin最大值原理研究最優(yōu)控制及其協(xié)態(tài)變量的表達(dá)式。為此,我們定義系統(tǒng)(17)的拉格朗日算子為

(21)

及所對(duì)應(yīng)的哈密頓函數(shù)

(22)

定理6 在系統(tǒng)(17)中,對(duì)最優(yōu)控制對(duì)(u*,v*)與狀態(tài)變量(S*,E*,P*,A*,I*,R*,M*)存在相應(yīng)的協(xié)態(tài)變量(λ1,λ2,λ3,λ4,λ5,λ6,λ7) 且

(23)

所對(duì)應(yīng)的橫截條件為

λ1(T)=0,λ2(T)=0,λ3(T)=0,

λ4(T)=0,λ5(T)=0,λ6(T)=0,λ7(T)=0。

此外最優(yōu)控制變量u*,v*的表達(dá)式為:

(24)

其中,

證明將(17)代入(22)式可得:

λ4((1-r)dpP-daA-μ1A)+

λ7(δ(P+A+I)-μ2M)。

對(duì)任意t∈[0,T], 根據(jù)Pontryagin最大值原理,那么存在協(xié)態(tài)變量λ1,λ2,λ3,λ4,λ5,λ6,λ7,滿足

故得出式(23)。

根據(jù)最優(yōu)條件,當(dāng)u=u*,v=v*時(shí),

經(jīng)過(guò)計(jì)算得出

(24)

根據(jù)允許控制集的定義, 有

(25)

5 數(shù)值擬合

5.1 敏感性分析

敏感性分析揭示了每個(gè)參數(shù)對(duì)疾病傳播的重要性,這些信息對(duì)研究傳染病模型的非線性動(dòng)力學(xué)行為至關(guān)重要,對(duì)如何控制和預(yù)防疾病傳播,如何采取有效而又經(jīng)濟(jì)的控制措施都具有重要的意義。敏感性分析主要分為兩種,一種是局部敏感性分析,一種是全局敏感性分析。局部敏感性分析主要通過(guò)每次改變一個(gè)輸入?yún)?shù),同時(shí)將其他參數(shù)保持在中心值來(lái)檢查輸出的局部影響,而全局敏感性(見(jiàn)文獻(xiàn)[12])分析則是通過(guò)大范圍擾動(dòng)模型輸入?yún)?shù),從而量化模型輸入對(duì)輸出的影響。本文主要通過(guò)全局敏感性分析來(lái)評(píng)價(jià)參數(shù)對(duì)基本再生數(shù)的相對(duì)重要性。

利用超立方采樣(Latin hypercube sampling, LHS)計(jì)算各參數(shù),并用這些參數(shù)的局部秩相關(guān)系數(shù)進(jìn)行全局敏感性分析。選取參數(shù)Λ=1,μ1=0.05,b=0.40,w=0.55,r=0.23,de=0.35,dp=0.43,da=0.35,di=0.35,可以得出不同R0的局部秩相關(guān)系數(shù)值,如圖1所示。

由圖1可知,參數(shù)μ1,da,dp對(duì)基本再生數(shù)R0有消極影響,de,r,dp,S0,b,w對(duì)基本再生數(shù)R0有積極影響。值得注意的有以下幾個(gè)方面:第一,疾病傳播率依然是影響最大的因素,其中包括無(wú)癥狀感染者的傳播率與癥狀傳播者的傳播率,無(wú)癥狀感染的傳播率影響稍大這也符合實(shí)際,這提醒人們應(yīng)該提高自身安全意識(shí),避免無(wú)癥狀感染。第二,易感者的人數(shù)對(duì)于基本再生數(shù)的影響也比較大,這代表在不同人口基數(shù)的城市,病毒傳播的速率可能大有不同。第三,考慮到疾病傳播函數(shù)與易感者人數(shù)和傳播率有關(guān),通過(guò)媒體覆蓋的影響改變易感者的行為從而延緩病毒的傳播非常重要。

5.2 穩(wěn)定性分析

在式(7)中,選取參數(shù)Λ=1,μ1=0.01,μ2=0.05,b=0.007,w=0.55,r=0.23,de=0.35,dp=0.43,da=0.35,di=0.35, 此時(shí)R0=0.02,S0=100,圖2和圖3驗(yàn)證了定理3,無(wú)病平衡點(diǎn)是漸進(jìn)穩(wěn)定的。

圖2 R0=0.02時(shí),S(t) 隨時(shí)間變化

圖3 R0=0.02時(shí),E(t),P(t),A(t),I(t) 隨時(shí)間變化

由圖3可知,無(wú)癥狀感染者和癥狀感染者都將經(jīng)歷一個(gè)短暫的上升再迅速下降,這個(gè)現(xiàn)象說(shuō)明即便R0<1,在一段時(shí)間內(nèi),感染人數(shù)依舊會(huì)增加,所以在疾病出現(xiàn),還是應(yīng)當(dāng)注意防護(hù)。其他參數(shù)保持不變,將b取為b=1.27,可得R0=4.04。

利用圖4和圖5,可以驗(yàn)證引理2和定理4,當(dāng)R0>1時(shí), 存在一個(gè)地方性平衡點(diǎn),而且這個(gè)地方性平衡點(diǎn)是漸進(jìn)穩(wěn)定的。

圖4 R0=4.04時(shí),S(t) 隨時(shí)間變化隨時(shí)間變化

圖5 R0=4.04時(shí),E(t),P(t),A(t),I(t)隨時(shí)間變化

由圖5可知,無(wú)癥狀感染者和癥狀感染者都將經(jīng)歷一個(gè)上升的階段后逐步穩(wěn)定,但是無(wú)癥狀感染者上升稍稍快一點(diǎn),這可能與無(wú)癥狀傳播的隱蔽性有關(guān)。

5.3 最優(yōu)控制

通過(guò)bvp4c 函數(shù)法進(jìn)行模擬,取Λ=1,μ1=0.01,μ2=0.05,b=1.27,w=0.55,r=0.23,de=0.35,dp=0.43,da=0.35,di=0.35, 得到了最優(yōu)控制函數(shù)圖(圖6和圖7)。

圖6 無(wú)癥狀感染者

圖7 癥狀感染者

由圖6和圖7可知,在加入最優(yōu)控制策略之后,無(wú)論是無(wú)癥狀感染者還是癥狀感染者人數(shù)都將發(fā)生顯著的下降,但是癥狀感染者的下降最為明顯,這說(shuō)明加強(qiáng)治療措施非常有效,因?yàn)樵诂F(xiàn)實(shí)生活中,癥狀感染者的治療與恢復(fù)非常重要。從另一角度考慮,治療措施主要針對(duì)的是癥狀感染者,那么對(duì)于預(yù)防控制無(wú)癥狀傳播,加強(qiáng)媒體覆蓋的強(qiáng)度這一措施就尤為重要。

由圖6可知,無(wú)癥狀感染者人數(shù)雖然最終也會(huì)得到很大程度的控制,但是依然會(huì)上升一段時(shí)間,這與現(xiàn)實(shí)生活中,無(wú)癥狀感染者的隱蔽性也是符合的。而由圖8與圖9可知,一方面間歇性的信息覆蓋比較合理,一方面在感染人數(shù)大幅度下降之后,應(yīng)該繼續(xù)保持一段時(shí)間的控制措施,以避免疫情再次上升。

圖8 信息強(qiáng)度(u)

圖9 治療強(qiáng)度

6 結(jié)束語(yǔ)

本文考慮了媒體傳播對(duì)于傳染病傳播的影響,建立了SEPAIRM 模型,推導(dǎo)了基本再生數(shù)R0,并且證明了當(dāng)R0<1時(shí), 系統(tǒng)的無(wú)病平衡點(diǎn)是局部漸進(jìn)穩(wěn)定的,當(dāng)R0>1時(shí),存在唯一的地方性平衡點(diǎn),并且得出了地方性平衡點(diǎn)局部漸進(jìn)穩(wěn)定的條件。在最優(yōu)控制研究中,通過(guò)Pontryagin最大值原理得出了最優(yōu)控制變量的表達(dá)式。在最后的數(shù)值模擬中,驗(yàn)證了無(wú)病平衡點(diǎn)和地方性平衡點(diǎn)的穩(wěn)定性,同時(shí)通過(guò)最優(yōu)控制函數(shù)圖也發(fā)現(xiàn),最優(yōu)控制策略對(duì)于降低感染人數(shù)是有效的,媒體覆蓋對(duì)于預(yù)防控制無(wú)癥狀傳播非常重要,而且間歇性的控制比較合適。