陳高偉，常向榮，陳俊英

(西南交通大學(xué)材料科學(xué)與工程學(xué)院，四川 成都 610031)

1 引 言

2019年12月被發(fā)現(xiàn)的新型冠狀病毒肺炎( COVID-19)具有極強的空氣傳播能力，在人與人之間具有極高的傳染率。COVID-19病原體可與人體粘膜細胞上的血管緊張素轉(zhuǎn)化酶Ⅱ(angiotensin converting enzyme-Ⅱ, ACE-Ⅱ)特異性結(jié)合，引發(fā)呼吸系統(tǒng)病變[1]，世界衛(wèi)生組織于2020年1月30日將COVID-19疫情宣布為國際公共突發(fā)衛(wèi)生事件。基于疫情數(shù)據(jù)，使用數(shù)學(xué)模型可以有效地挖掘COVID-19的各種流行病學(xué)參數(shù)，對醫(yī)學(xué)研究與疾病防控具有重要意義[2]。本研究對傳統(tǒng)倉室SIR模型進行改進，采用Gauss函數(shù)模擬日傳染率、日隔離率的變化情況，根據(jù)2020年2月13日至2020年3月2日國家衛(wèi)生健康委員會發(fā)布的湖北地區(qū)實際疫情數(shù)據(jù)，反推出COVID-19疫情傳染、防控情況的演變，并依此推測了COVID-19的基本再生數(shù)與有效再生數(shù)。

2 模型的建立

本研究基于倉室SIR模型來模擬COVID-19的傳播情況。以λ表示日傳染率(平均每日每名感染者所能感染的易感者數(shù)量)，ρ表示日隔離率(平均每日可得到有效控制、移出感染體系的感染者數(shù)量)。模型建立在以下假設(shè)基礎(chǔ)之上：

(1)不考慮湖北省人口的自然增長與死亡，且數(shù)據(jù)始于2020年春節(jié)以后，疫情導(dǎo)致返工時間延遲，故不考慮人口遷移。

(2)患者一經(jīng)確診，立刻接受嚴格隔離(入院治療或居家隔離)，不再具備傳染性。

(3)感染者自被感染到發(fā)病所需時間為4 d(即COVID-19潛伏期中位數(shù)[3])。

為貼合實際情況，本研究對傳統(tǒng)SIR模型提出如下改進：

(1)重新解讀R、ρ的實際含義。根據(jù)假設(shè)(2)，確診者可以視為已退出傳染系統(tǒng)。因此，中國疾控中心發(fā)布的確診人數(shù)對應(yīng)于本模型的移出者(removed，R)類人群。同時，參數(shù)ρ的含義也等同于診出率。根據(jù)假設(shè)(3)，結(jié)合COVID-19在潛伏期同樣具備傳染性的事實，每日確診人數(shù)應(yīng)對應(yīng)于4 d前的實際感染人數(shù)。

(2)采用隨時間緩慢變化的Gauss函數(shù)代替常數(shù)以模擬λ、ρ的演變情況。

設(shè)s、i、r分別為三類人群所占的人口比例(湖北省人口取59 270 000)，據(jù)SIR模型，可得微分方程組：

(1)

其中，λ、ρ以類Gauss函數(shù)描述：

λ(t)=ae-t2+b;μ(t)=1-ce-dt2

(2)

式中a、b、c、d為待定系數(shù)，需依據(jù)實際疫情數(shù)據(jù)來確定。記ir,j為疫情第j天湖北省現(xiàn)存確診人數(shù)比，ic,j為根據(jù)改進的SIR模型計算所得的第j-4 d感染人數(shù)比,則可歸結(jié)為最小二乘問題：

(3)

式中k、n分別為有效數(shù)據(jù)的起、始點。在國家衛(wèi)健委發(fā)布的湖北省疫情數(shù)據(jù)中，2月13日發(fā)布數(shù)據(jù)(即2月12日0時至24時內(nèi))新增確診病例14 840例，增幅異常。事實上，為確保信息真實透明，武漢市于2月12日實施全面搜查，并將臨床初步診斷病例納入確診之列。因此，本研究取2月13日及以后的數(shù)據(jù)為有效數(shù)據(jù)。取1月24日為第1 d，則k=21，n=39。

3 模型的求解

3.1 初步結(jié)果

引入四階Runge-Kutta算法求解微分方程組，并使用MATLAB軟件內(nèi)置的fminsearch函數(shù)求解最優(yōu)化問題，取迭代初始值(a0,b0,c0,d0)=(0.1,0.5,0.7,0.001)，求解得待定系數(shù)(a,b,c,d)=(0.1663，0.9912，0.2188，0.0032)，累積最小二乘誤差為8.2947×10-10。數(shù)據(jù)擬合結(jié)果見圖1。

圖1 數(shù)據(jù)擬合結(jié)果(2月13日-3月2日))

3.2 殘差分析

殘差為擬合數(shù)據(jù)與實際數(shù)據(jù)之差。記t時刻的擬合數(shù)據(jù)為ic(t)，實際數(shù)據(jù)為ir(t)，定義標準化擬合殘差e*為：

(4)

σ為殘差標準差。根據(jù)擬合模型基本假設(shè)，殘差服從零均值的正態(tài)分布，故取樣本標準差作為σ的無偏估計：

(5)

式中n為殘差樣本容量。現(xiàn)對標準化殘差e*的期望是否為0進行檢驗。構(gòu)造t統(tǒng)計量：

(6)

取顯著性水平α=0.05，計算得：

(7)

因此，接受e*的期望為0假設(shè)。由此可知，本模型對數(shù)據(jù)的擬合程度相當高，能有效地反映疫情變化的過程。

4 分析與討論

4.1 λ、ρ的演變情況

依據(jù)數(shù)據(jù)擬合得出的待定系數(shù)a、b、c、d，作出λ、ρ的函數(shù)曲線，見圖2。

圖2 感染速率λ、診出率ρ函數(shù)曲線

由圖可知，傳染率λ從起初的1.1575迅速下降至0.9912并保持穩(wěn)定，本模型所選取的時間區(qū)間為Gauss函數(shù)的下降至末尾與后期穩(wěn)定部分。實際上，1月24日前，有關(guān)部門已采取一定措施遏制疫情發(fā)展。1月20日，武漢市加強進出人員的管控力度。因此，本模型得出λ的演變情況與實際相吻合。診出率ρ起初停留在較低的水平(78%)，且發(fā)展緩慢，主要由于臨床醫(yī)生缺乏診斷經(jīng)驗。隨著診斷技術(shù)與診斷標準的進步，由單純核酸檢測逐步發(fā)展為結(jié)合核酸檢測、胸部影像學(xué)檢測、臨床診斷為一體的檢測方案[3]，ρ加速增長且趨向于1。因此，本模型所得診出率的演變情況與實際吻合較好。

圖3 整體擬合曲線(1月24日～3月2日)

4.2 基本再生數(shù)與有效再生數(shù)的推測

基本再生數(shù)(Basic Reproduction Number)R0為一個感染者在無干預(yù)的情況下可以傳染的易感者人數(shù)的均值，是現(xiàn)代傳染病學(xué)的一大核心指標，用于刻畫傳染病的最大潛在感染能力，對疫情防控政策的制定與實施具有重要指導(dǎo)意義?；诜N群增長的生態(tài)模型，Wallinga等[4]提出了使用指數(shù)增長率與不同年齡的群體生育率計算R0。推廣至基于微分方程組的傳染病SIR模型中，結(jié)合Dreessche等[5]建立的利用SIR模型參數(shù)計算R0的公式，可得：

(8)

考慮到R0定義中強調(diào)感染的無干預(yù)性，λ在計算時選取最大值(1.1575)，ρ選取最小值(0.7812)，得R0=2.4817。

對于COVID-19的基本再生數(shù)，各學(xué)者已做出初步的研究與估計?；诒灸Ｐ偷耐茰y結(jié)果，數(shù)值水平與之接近，具有參考價值。現(xiàn)列舉部分R0估計值如下：

表1 COVID-19基本再生數(shù)參考Table 1 Reference for the basic reproduction number of COVID-19

由于R0是對無干預(yù)感染系統(tǒng)的估算，加之方法種類繁多、統(tǒng)計樣本偶然性較大，故另采用有效再生數(shù)(Effective Reproduction Number)Re來刻畫傳染病的傳播過程[11-12]。Re定義為：在考慮防控措施的情況下，一個感染者可傳染的易感人數(shù)的均值。結(jié)合本模型λ的定義，設(shè)T為疫情持續(xù)的總時間，可得：

(9)

代入擬合所得的λ數(shù)據(jù)，得Re=0.9912。

5 結(jié)論

本研究對新型冠狀病毒在湖北地區(qū)2020年2月13日至3月2日的疫情數(shù)據(jù)進行分析，利用數(shù)學(xué)模型計算出COVID-19的基本再生數(shù)R0與有效再生數(shù)Re。R0=2.4817，結(jié)果與當下各學(xué)者的結(jié)論基本吻合，體現(xiàn)了COVD-19在自然條件下的高度傳染性；Re=0.9912，表明在當下防控措施下，平均一名感染者每日仍可以感染約1名易感者。擬合殘差分析表明，本模型擬合程度相當高，證明了Gauss函數(shù)能更好擬合病毒的實際感染率，及隔離率隨時間的真實變化情況。在本研究建模過程中，未考慮疫情潛伏期的隨機分布性，無法準確地模擬感染者從感染至發(fā)病、確診的時間延遲。未來可引入時滯動力學(xué)參數(shù)對疫情做進一步分析。