楊波, 于振華, 蔡遠利

(1.西安交通大學自動化科學與工程學院, 710049, 西安; 2.西安科技大學計算機科學與技術學院, 710054, 西安)

新型冠狀病毒(SARS-CoV-2)[1-2]出現(xiàn)以來,由于早期對其缺乏足夠的認識,在多個地區(qū)出現(xiàn)了擴散和大流行,對全世界產(chǎn)生了巨大的影響。這種病毒所引起的傳染性肺炎疫情被世界衛(wèi)生組織命名為“Coronavirus disease 2019”(COVID-19)[3]。為了遏制疫情的擴散,中國政府采取了一系列高效嚴格的防控措施[4],已經(jīng)完全控制住了疫情,但是在世界其他地區(qū),COVID-19疫情依然在繼續(xù),感染人數(shù)和死亡人數(shù)不斷增加,嚴重影響了人們的正常生活和工作,給各國社會經(jīng)濟造成巨大損失,給人們生命健康造成重大威脅。據(jù)世界衛(wèi)生組織的報告,截止2021年04月07日,全球已有131 487 572人感染COVID-19,因病死亡人數(shù)達2 857 702[5]。特別是2020年12月14日,英國報告了COVID-19的一種新的變異毒株,其傳染性可能更強[6],現(xiàn)已在世界范圍內(nèi)傳播。不僅如此,最近已發(fā)現(xiàn)多種變異毒株,給疫情防控帶來更大挑戰(zhàn)。目前,疫情仍在世界范圍內(nèi)廣泛傳播,針對COVID-19還沒有特效藥物,而疫苗供應量有限,不能滿足目前的需求,因此高效嚴格的非藥物防控措施仍是減緩疫情擴散的最有效途徑,特別是欠發(fā)達國家和地區(qū)。各國政府和人民依然在共同努力,加強防控,希望早日控制住COVID-19的進一步傳播,讓世界恢復昔日的繁榮。

數(shù)學模型是定量分析各種傳染病傳播動態(tài)(包括網(wǎng)絡病毒的傳播,例如文獻[7-9])的重要工具,它能模擬和預測疫情的發(fā)展趨勢,為制定相關的防控措施提供一定的參考。根據(jù)COVID-19疫情的特點,已建立了很多相關的數(shù)學傳播模型,例如文獻[10-21]等。但在這些模型中,文獻[12-14,17-20]中沒有考慮無癥狀感染者,而文獻[11,15]中沒有考慮無癥狀感染者可以發(fā)展為有癥狀的情況,文獻[16]中沒有考慮潛伏期。本文主要考慮政府的嚴格防控措施對疫情的影響,建立相應的非線性動力學傳播模型,模擬COVID-19的傳播動力學。模型除了考慮嚴格的防控措施外,還考慮了無癥狀感染者以及無癥狀感染者可以發(fā)展為有癥狀等實際情況。求出了模型的基本再生數(shù);應用武漢市的真實報告數(shù)據(jù),對模型參數(shù)和狀態(tài)初值進行了估計,用所建模型分別模擬了COVID-19在武漢市和湖北省除武漢市外地區(qū)的傳播動態(tài),并與實際報告數(shù)據(jù)對比。對參數(shù)做了靈敏度分析,并進行了仿真驗證,為疫情的有效防控提供一定參考。

1 模型的建立

對于COVID-19,還有很多未解之謎有待解開,但是以目前的研究和數(shù)據(jù)來看,COVID-19的傳播機理是基本清楚的。人們感染COVID-19后,會有一段時間的潛伏期,之后會發(fā)展為無癥狀或有癥狀的感染者。由于COVID-19的傳染力較強,故需要對患者和疑似患者進行隔離治療,對密切接觸者實施隔離觀察。除此之外,根據(jù)疫情的發(fā)展情況,政府部門還出臺了一些措施,包括交通管制,限制旅行、聚集,保持社交距離,商場停業(yè)、學校停課、工廠停產(chǎn)等,以減少易感人群感染COVID-19的機率。根據(jù)COVID-19的傳播特點和部分嚴格的疫情防控措施,將人群分為7組:易感人群(S)、潛伏期人群(E)、隔離觀察人群(Q)、無癥狀感染人群(A)、有癥狀感染人群(I)、感染后住院隔離治療人群(J)、康復人群(R)。

根據(jù)實際防控措施,做如下假設:①所有未感染COVID-19的人都是易感的;②E類人是被感染,沒有癥狀,但不具有傳染性,E類人群主要是感染者的密切接觸者(暴露者),一旦被跟蹤到,將被醫(yī)學隔離觀察,即成為Q類;③Q類人一旦確診將進入醫(yī)院隔離治療,即成為J類;④A類人群是被感染了,沒有癥狀,但具有傳染性,這是區(qū)別于E類人群的,且A類人可發(fā)展為有癥狀的I類;⑤A和I類人群一旦被確診,將進入醫(yī)院隔離治療,即成為J類,他們也可以自愈,直接成為R類人群;⑥重癥患者都就醫(yī)治療,故因病死亡者只包括在醫(yī)院治療無效死亡的人;⑦防護措施得當,操作規(guī)范、嚴格,因此Q和J類人群不會再傳染其他人;⑧感染者康復后,即R類人群,具有對COVID-19的免疫,即短時間內(nèi)不會被再次傳染。

圖1 人群演化示意圖Fig.1 Evolution of population

不同組群之間的具體轉換關系如圖1所示,參數(shù)說明見表1。建立如下非線性動力學模型

(1)

式中:M表示平均每天新增人口;S(t)、E(t)、Q(t)、A(t)、I(t)、J(t)和R(t)表示t時刻各人群數(shù);N(t)=S(t)+E(t)+Q(t)+A(t)+I(t)+J(t)+R(t)為t時刻總人口數(shù)。其初始條件為

(S(0),E(0),Q(0),A(0),I(0),J(0),R(0))≥0

(2)

狀態(tài)空間為

(3)

命題1設初值條件為

(S(0),E(0),Q(0),A(0),I(0),J(0),R(0))>0

(4)

則系統(tǒng)(1)的解滿足:?t>0,有

(S(t),E(t),Q(t),A(t),I(t),J(t),R(t))>0

(5)

證明:用反證法證明。若不然,設存在t1>0,使得S(t)、E(t)、Q(t)、A(t)、I(t)、J(t)、R(t)中至少有一個是非正的。由連續(xù)性可知,存在0

(b)若S(t0)>0,E(t0)=0,則

(6)

(c)若S(t0)>0,E(t0)>0,而Q(t0)=0,則

(7)

類似(a)的分析,這會導致與t0的取法矛盾,故而Q(t0)>0。同理可得A(t0)=0、I(t0)=0、J(t0)=0、R(t0)=0都會產(chǎn)生矛盾,故S(t)、E(t)、Q(t)、A(t)、I(t)、J(t)、R(t)都是正的。

命題2系統(tǒng)(1)是有界的。

證明:將系統(tǒng)(1)的所有方程加起來可得

(8)

上述不等式兩邊同乘以eυt,有

(9)

則有

(10)

(11)

上述命題表明

Ω={(S(t),E(t),Q(t),A(t),I(t),J(t),R(t)):(S(t),E(t),Q(t),A(t),I(t),J(t),R(t))≥0;S(t)+E(t)+Q(t)+

(12)

為系統(tǒng)(1)的正不變集。

表1 模型(1)參數(shù)說明

2 基本再生數(shù)

模型的基本再生數(shù)R0是指一個感染者(病人)在平均感染期內(nèi)所傳染的人數(shù),是研究傳染病的一個重要概念,是判斷傳染病是否會消亡的閾值。當R0<1時表明一個感染者在平均感染期內(nèi)能傳染的人數(shù)小于1,意味著傳染病將逐步消亡;反之,當R0>1,傳染病將始終存在,不會消亡[22]。因此要消滅傳染病,需要通過分析R0,制定相應的防控措施,使其減小到1以下。

首先求出模型的無病平衡點。令模型(1)中7個方程全等于0,由A(t)=I(t)=0,可得E(t)=Q(t)=J(t)=0,從而得無病平衡點

(13)

應用下一代矩陣方法[23-24]求R0。設

(14)

(15)

其中

(16)

F和V在P0處的雅克比矩陣分別為

(17)

(18)

其中

(19)

(20)

(21)

則

(22)

從而可得基本再生數(shù)

(23)

分析R0的表達式,不難發(fā)現(xiàn)采取如下措施可以使其減小:①加大對密切接觸者(暴露者)檢疫隔離力度,則遺漏在易感人群中的感染者(A和I類)將變少,即σ1增大,此時δ2和δ3會減小;②加強無癥狀感染者的隔離和治療,即σ1和σ3增大;③醫(yī)療資源允許的情況下,讓有癥狀感染者盡可能入院治療,即γ1和γ2增大;④加強對易感人群的管控,如加強疫情防控宣傳,控制人員流動,保持一定的社交距離,關閉人群容易聚集的場所等,即u減小。

上述措施中①和④是最有效果的,可以明顯減小R0,這與實際相符,特別是措施④。

3 模型驗證

本節(jié)利用COVID-19在武漢市傳播的相關實際報告數(shù)據(jù)擬合模型參數(shù)并驗證模型。為了說明模型具有廣泛的適用性,將模型應用于湖北省除武漢市外的地區(qū)。

3.1 參數(shù)擬合

用模型(1)來模擬COVID-19在中國武漢市的傳播。利用武漢市的COVID-19相關數(shù)據(jù)對模型中的參數(shù)和狀態(tài)初值進行估計。初期由于大量人員被感染,致使醫(yī)療資源被嚴重擠兌,包括病床數(shù)量、醫(yī)護人員、醫(yī)療防護物資等嚴重缺乏,很多病人無法入院治療。輕癥患者被要求在家自我隔離,這導致家庭成員感染增加,無法阻斷COVID-19的傳播。2020年2月上旬開始,在中國政府的強有力干預下,迅速改造了大量現(xiàn)有醫(yī)院,新建了兩所大型臨時醫(yī)院,用于收治重癥病例;從其他省調(diào)集了大量醫(yī)護人員和醫(yī)療物資馳援武漢市;將會展中心等大型公共場所改造為簡易臨時醫(yī)院,即所謂的“方艙醫(yī)院”,集中收治輕癥患者;將酒店、學校等設置為醫(yī)學隔離觀察點,用于隔離觀察密切接觸者等。這些措施的目的是將感染者盡可能的都集中起來,進行醫(yī)學隔離治療或觀察,徹底切斷COVID-19的傳染源。

本文模型正是基于這些強有力的措施而建立的。前期由于疑似病人較多,而每天核酸檢測能力有限,有很多病例沒有得到及時確診,因此前期數(shù)據(jù)不夠準確?；谝陨显?選取武漢市2020年2月14日～3月4日共20天的數(shù)據(jù)對模型進行校準,即對模型中的參數(shù)和狀態(tài)初值進行擬合,使用的實際報告數(shù)據(jù)來源于中國國家衛(wèi)生健康委員會官網(wǎng)[25]。隨著醫(yī)護人員和醫(yī)療物資的補充,診療經(jīng)驗的積累,病死率逐漸減小,治愈率不斷提高,因此假設θ=θ(t),μ=μ(t)都是時間的函數(shù),這更加符合實際。衛(wèi)生部門實時報告數(shù)據(jù)主要有:確診人數(shù)(住院人數(shù))、死亡人數(shù)、出院人數(shù)等。首先根據(jù)實際報告數(shù)據(jù),對θ(t)和μ(t)進行擬合,如圖2所示。圖中紅色點是由實際報告數(shù)據(jù)計算得到的值,曲線是擬合的。分別用二次函數(shù)和負冪指函數(shù)擬合每天的治愈率和病亡率,得

θ(t)=0.000 099 363t2+0.001 939 74t+

0.006 688 72

(24)

μ(t)=0.007 743 668 47t-0.183 432 3-

0.003 384 146

(25)

(a)治愈率

(b)病亡率圖2 θ(t)和μ(t)的擬合(武漢市)Fig.2 Fitting of θ(t) and μ(t) (in Wuhan)

確診病例基本都是有癥狀的感染者,且一旦確診就會進入醫(yī)院隔離治療,即確診人數(shù)就是住院人數(shù),包括醫(yī)院和方艙醫(yī)院所有病人。利用現(xiàn)有確診病例數(shù)(當前在醫(yī)院隔離治療的病例數(shù))的實際數(shù)據(jù)與模型的解J(t)做最小二乘擬合[16],估計模型中參數(shù)值和狀態(tài)初值。由于政府實施的強力防疫措施和民眾的積極配合,取u=0.3。且此時人口流動幾乎沒有,20天內(nèi)的出生人口有限,故可以認為M=0。武漢市人口大約1 000萬,受疫情防控措施的影響以及對病毒的認識,絕大多數(shù)人都盡量不出門,不與外界接觸,故假設實際易感人數(shù)初值為10萬。以2020年2月14日為J(t)時刻,武漢市現(xiàn)有確診病例數(shù)擬合仿真如圖3所示,其中紅色點表示真實報告數(shù)據(jù),藍色曲線表示模型解J(t)。模型參數(shù)的估計值見表2,狀態(tài)初始值見表3。

表2 模型(1)參數(shù)值及其靈敏度指數(shù)

表3 模型(1)狀態(tài)初值(武漢市)

圖3 現(xiàn)有確診病例數(shù)(武漢市)Fig.3 Currently confirmed case number (in Wuhan)

用平均相對誤差來衡量擬合的效果。平均相對誤差定義如下

(26)

3.2 模擬COVID-19在武漢市的傳播

利用累計確診病例數(shù)、治愈出院病例數(shù)和死亡病例數(shù)的實際報告數(shù)據(jù)來驗證模型。根據(jù)模型(1),第t天累計確診病例數(shù)(住院病例數(shù))為

(27)

第t天累計治愈出院病例數(shù)Z(t)和病死數(shù)D(t)分別為

(28)

(29)

其中Jc(0)、Z(0)和D(0)分別為2020年2月14日實際報告的累計確診病例數(shù)、治愈出院病例數(shù)和死亡病例數(shù),見表3,武漢市COVID-19的實際與模擬傳播動態(tài)如圖4所示。與真實數(shù)據(jù)相比,各類病例數(shù)的平均相對誤差見表4。從圖4和表4中可以看出,模型(1)與實際數(shù)據(jù)基本吻合,能反映在政府強有力干預下COVID-19在武漢市的傳播動態(tài)。此時,基本再生數(shù)R0為0.085 66,已遠小于1了,可以說已基本控制住了COVID-19的傳播,符合實際情況。

(a)累計確診病例數(shù)

(b)累計治愈出院病例數(shù)

(c)累計死亡人數(shù)圖4 COVID-19的實際與模擬傳播動態(tài)(武漢市)Fig.4 Real and simulated spread dynamics of COVID-19 (in Wuhan)

表4 各類病例數(shù)的平均相對誤差

3.3 模擬COVID-19在湖北省(除武漢市)的傳播

模型(1)不僅適用于COVID-19在武漢市的傳播,也適用于其他對疫情嚴格防控的地區(qū)。本小節(jié)以湖北省除武漢市外的其他地區(qū)為研究對象,驗證模型(1)。由于疫情中心是武漢市,故假設湖北省其他地區(qū)的防控力度稍稍弱于武漢市,即取u=0.4。模型其他參數(shù)采用3.1節(jié)中的擬合參數(shù)(見表2),結合實際數(shù)據(jù),只需對部分狀態(tài)初值進行擬合。時間范圍同樣取2020年2月14日～3月4日內(nèi)的20天。類似地,首先根據(jù)實際數(shù)據(jù)的分布,分別用一次函數(shù)和常數(shù)擬合湖北省(除武漢市)每天的治愈率θ(t)和病亡率μ,擬合圖如圖5所示,結果如下

θ(t)=0.006 971t+0.012 384 6

(30)

μ=0.001 677 41

(31)

圖6 現(xiàn)有確診病例數(shù)(湖北省除武漢市外其他地區(qū))Fig.6 Currently confirmed case number (in Hubei except Wuhan)

同樣利用現(xiàn)有確診病例數(shù)的實際數(shù)據(jù)與模型的解J(t)做擬合,現(xiàn)有確診病例數(shù)(湖北省(除武漢市))如圖6所示,狀態(tài)初值見表5。類似地,計算得到現(xiàn)有確診病例擬合的平均相對誤差為EMRE=0.833%。

(a)治愈率

20天內(nèi)新增確診病例數(shù)、累計確診病例數(shù)Jc(t)、累計治愈出院病例數(shù)Z(t)和累計死亡病例數(shù)D(t)的傳播動態(tài)如圖7所示,平均相對誤差見表4。

(b)病亡率圖5 θ(t)和μ(t)的擬合(湖北省除武漢市外其他地區(qū))Fig.5 Fitting of θ(t) and μ(t) (in Hubei except Wuhan)

表5 模型(1)狀態(tài)初值(湖北省除武漢市外其他地區(qū))

(a)累計治愈出院病例數(shù)

(b)累計死亡病例數(shù)

(c)累計確診病例數(shù)圖7 COVID-19的實際與模擬傳播動態(tài)(湖北省除武漢市外其他地區(qū))Fig.7 Real and simulated spread dynamics of COVID-19 (in Hubei except Wuhan)

從圖7和表4可以看出,累計治愈出院病例Z(t)和累計死亡病例數(shù)D(t)與實際報告數(shù)據(jù)較吻合,而從圖7c可以看出,2月23日以前累計確診病例實際報告數(shù)據(jù)波動較大,這主要是因為前期疑似病例較多,而核酸檢測能力有限,為了盡快消化疑似病例,衛(wèi)生部門制定了臨床診斷標準來區(qū)分確診病例。后期核酸檢測能力充裕了,對這些臨床診斷病例再進行核酸檢測,非COVID-19感染病例就會核減掉?？傮w上說,模型(1)能真實反映COVID-19在湖北省除武漢市外的地區(qū)的傳播動態(tài)。

4 靈敏性分析

為了降低COVID-19的傳播,需要考慮采取哪些措施最有效,即要分析哪些參數(shù)對基本再生數(shù)R0的影響最大。這里采用文獻[26]中的靈敏度分析方法,計算基本再生數(shù)R0對各參數(shù)的靈敏性指數(shù)。

定義[26]:設函數(shù)Φ關于變量x可偏導,則Φ關于x的靈敏度指數(shù)為

(32)

為了充分說明參數(shù)的調(diào)整對OVID-19播動態(tài)的影響,根據(jù)上述靈敏性分析,調(diào)整部分參數(shù)進行仿真對比。以武漢市為例,假設防控措施可以更嚴格一些,即u可以更小,取u=0.2,隔離率δ1由0.767 18增加到0.867 18,其他參數(shù)不變,調(diào)參前后各狀態(tài)變化趨勢分別如圖3、圖4和圖8所示。從圖8可以看出,調(diào)整兩個參數(shù)的值后,未被發(fā)現(xiàn)的暴露者E(t)、無癥狀感染者A(t)和有癥狀感染者I(t)的數(shù)量比之前有所減少;而隔離觀察者Q(t)(圖8b)明顯變多了,會導致新增確診病例增加,從而累計確診病例Jc(t)、累計治愈出院病例Z(t)和現(xiàn)有確診病例J(t)均有所增加。這些都是希望看到的結果,人們希望所有感染者都能盡快確診,一旦確診了就會被隔離治療,避免了COVID-19的進一步傳播。這些仿真驗證了靈敏度分析的結果。

(a)暴露(潛伏期)人群E(t)

(b)隔離觀察人群Q(t)

(c)無癥狀感染人群A(t)

(d)有癥狀感染人群I(t)圖8 模型(1)狀態(tài)變化趨勢(武漢市)Fig.8 Trends of states of model I (in Wuhan)

5 結 論

距離COVID-19首次發(fā)現(xiàn)至今已有一年多時間了,然而COVID-19在全球范圍內(nèi)仍然處于大流行階段,而且還出現(xiàn)了傳染力更強的變異病毒,這給全球防疫帶來挑戰(zhàn)。本文建立數(shù)學模型模擬COVID-19在一些嚴格防疫措施下的傳播,具體工作和結論總結如下:

(1)利用實際報告的有關數(shù)據(jù)對模型參數(shù)進行了擬合,并用模型(1)分別對COVID-19在武漢市和湖北省除武漢市外的地區(qū)的傳播進行了模擬仿真,仿真結果基本與實際數(shù)據(jù)吻合,說明模型(1)是合理的,能較準確地反映COVID-19在嚴格防控措施下的傳播動態(tài);

(2)通過對模型參數(shù)的靈敏性分析,給出各參數(shù)對COVID-19傳播影響的排序,其中對COVID-19傳播影響最大的是u和δ1,對制定相應的防控措施提供了理論依據(jù),仿真也證實了靈敏性分析結果。

通過本文的建模和分析,希望能在制定防疫措施阻止COVID-19的傳播方面提供一些有益建議。