馬思婕，黃珈銘，印英東，曾楚儀

（廣州大學(xué)數(shù)學(xué)與信息科學(xué)學(xué)院，廣東 廣州 510006）

0 引言

2019年末，新冠突襲無疑是人類史上的一大災(zāi)難，短短幾日便讓人們看到了它的“威力”。2020年3月11日，世衛(wèi)組織將新冠肺炎疫情認定為全球大流行。實際上，全球大流行并不是第一次出現(xiàn)，西班牙流感、H1N1豬流感、西非埃博拉疫情……這也側(cè)面告訴人們：在以往經(jīng)歷中吸取經(jīng)驗，尋求合適的預(yù)測模型進行科學(xué)防疫尤為重要。

在疫情早期便有許多學(xué)者開始了與COVID-19相關(guān)的研究，不過早期的數(shù)據(jù)、現(xiàn)象還不夠成熟，過早過快地下結(jié)論會造成較大偏差。本文在相對成熟的數(shù)據(jù)及現(xiàn)象下進行了分析預(yù)測，對疫情走勢進行有效預(yù)測及傳染力判斷，希望通過對COVID-19的研究，得到較成熟的模型，為日后類似突發(fā)情況提供幫助。

1 模型建立

1.1 人群流動分析

在考慮居家隔離、醫(yī)學(xué)隔離、入院治療等防疫政策和醫(yī)學(xué)干預(yù)的情況下，本文對各艙室人群的流動情況進行了分析（見圖1）。

圖1 人群流動變化

在模型建立的過程中，為了排除其他因素干擾，將不考慮非COVID-19所引起的人口變化，如自然死亡等。在將模型實踐到具體地區(qū)時，假設(shè)該地區(qū)在所研究時間段內(nèi)人口恒定，系統(tǒng)內(nèi)各類人口總和不發(fā)生改變。浦江創(chuàng)新論壇系列活動“特別對話：科技創(chuàng)新與全球健康共治”［1］中，復(fù)旦大學(xué)附屬華山醫(yī)院感染科主任張文宏表示：在現(xiàn)實中存在復(fù)陽、二次感染情況，但屬于極小概率事件，故在模型中將不考慮二次感染情況，即被感染者人群痊愈后歸入康復(fù)者類中，不再流動轉(zhuǎn)變。

1.2 各類人群關(guān)聯(lián)性的微分方程建立

根據(jù)1.1中各類人群的流動分析，可建立如下各類人群關(guān)聯(lián)性的微分方程：

該模型中引進了一些參數(shù)以表現(xiàn)群體的變化流動。式中：c代表一個暴露的具有感染力的患者的接觸率；β代表每次接觸的傳播概率；q代表接觸者被追蹤到隔離的概率；σ代表受感染個體轉(zhuǎn)向受感染類的轉(zhuǎn)化率；λ代表隔離的未感染接觸者被釋放到更廣泛的社區(qū)的速度；δI代表有癥狀的感染者轉(zhuǎn)變?yōu)樽≡夯颊叩谋嚷?；δq代表隔離易感人群轉(zhuǎn)變?yōu)樽≡夯颊叩谋嚷?；θ代表無癥狀感染者相對于有癥狀感染者傳播力的比值；γI代表未住院有癥狀感染者的康復(fù)率；γA代表無癥狀感染者的康復(fù)率；γH代表住院患者的康復(fù)率；α代表死亡率。

2 模型求解

本文分別以北京市總?cè)丝凇⑽錆h市總?cè)丝跒槌跏家赘腥巳鹤鳛檠芯繉ο?，通過北京市衛(wèi)生健康委員會網(wǎng)站［2］、湖北省衛(wèi)生健康委員會網(wǎng)站［3］公開發(fā)布的官方數(shù)據(jù)為原始數(shù)據(jù)，收集2020年1月23日—3月17日的疫情信息并對所獲得的數(shù)據(jù)進行處理與篩選。數(shù)據(jù)處理與篩選后出現(xiàn)數(shù)據(jù)缺失的情況，因此在得到正確的實驗數(shù)據(jù)前本文進行了數(shù)據(jù)的缺失度檢驗，發(fā)現(xiàn)數(shù)據(jù)缺失度為1.37%，缺失度較低。經(jīng)分析可知缺失部分為MAR隨機缺失，可直接選用模型方法對數(shù)據(jù)進行填補。對于處理后得到的理想數(shù)據(jù)進行去噪處理與匯總，得到兩城市較準確的COVID-19感染人數(shù)變化數(shù)據(jù)，為本文模型的訓(xùn)練預(yù)測擬合提供重要基礎(chǔ)。

2.1 參數(shù)估計相關(guān)計算

（1）由疫情防控指南了解到［4］，一般隔離周期為14 d，故可得到隔離的未感染接觸者被釋放到更廣泛社區(qū)的速度λ=1/14。

（2）由世衛(wèi)組織的相關(guān)研究［4］可知，COVID-19的平均潛伏期為7 d，故可得到受感染個體轉(zhuǎn)向受感染類的轉(zhuǎn)化率σ=1/7。

（3）其他參數(shù)的估測方法。對于SEIR模型中的參數(shù)估計，本文采用基于分支過程的Monte-Carlo算法［5］，通過建立隨機并基于分支過程描述COVID-19的傳播過程（見圖2），從而對數(shù)據(jù)進行預(yù)測（見表1）。

表1 參數(shù)的計算估計取值

圖2 分支過程

根據(jù)分支過程理論，病毒傳播過程中人與人接觸產(chǎn)生的傳播是隨機獨立的。對這種大量隨機行為的描述符合計算機Monte-Carlo算法的應(yīng)用要求。

本文以收集到感染人數(shù)數(shù)據(jù)的前80%作為訓(xùn)練集擬合目標(biāo)函數(shù)，后20%作為參數(shù)估計集用于優(yōu)化函數(shù)中的參數(shù)取值。為了增強算法應(yīng)用的準確性，本文在執(zhí)行模型的過程中，不斷調(diào)整步長以求得較優(yōu)的參數(shù)取值。對實際數(shù)值和模擬數(shù)值用最小二乘法衡量其相似程度，不斷優(yōu)化參數(shù)取值。由多種參數(shù)的取值并兼顧疫情發(fā)展帶來的參數(shù)變化，進行疫情發(fā)展過程模擬。模擬結(jié)果與實際數(shù)據(jù)進行對比，發(fā)現(xiàn)模擬的貼合度達到91.49%，有理由認為算法擬合出的模擬情況從概率意義上表明了未來實際疫情的發(fā)展趨勢。

2.2 基本再生數(shù)R0的計算結(jié)果

由于疫情暴發(fā)中經(jīng)歷了醫(yī)療資源不足以及相關(guān)隔離政策頒布的過程，因此COVID-19傳播前期的R0與傳播后期的R0有較大區(qū)別［6］。根據(jù)現(xiàn)實情況可知，由于醫(yī)療資源供應(yīng)在2020年2月13日左右充足，因此本文以2020年2月13日為分界點，對兩城市的R0進行前期和后期計算，得到結(jié)果，如表2所示。

表2 R0分界計算結(jié)果

2.3 基于模型的感染者數(shù)據(jù)走勢預(yù)測

通過上文中相關(guān)參數(shù)的估計，對后續(xù)感染者的數(shù)據(jù)走勢進行了預(yù)測，預(yù)測情況，如圖3所示。

圖3 預(yù)計累計感染人數(shù)與實際情況比對

3 靈敏度分析

為研究本模型是否能較為靈敏地進行預(yù)測，且誤差是否會對模型結(jié)果造成顯著影響，本文通過改動在SEIR模型中對傳染影響較大的兩個參數(shù)c，q，以觀察最終擬合曲線的改變情況，結(jié)果如圖4所示。

圖4 基于參數(shù)c，q的模型靈敏度分析

改動參數(shù)c，即暴露的具有感染力的患者的接觸率發(fā)生了改變，這是對疫情下防疫政策的直觀體現(xiàn)。觀察圖4可以看到，模型預(yù)測的累計感染人數(shù)產(chǎn)生較大波動，但參數(shù)變化造成的預(yù)測值變化在合理范圍內(nèi)。表明本模型在進行預(yù)測時對于參數(shù)c的變動感應(yīng)較為靈敏，可以依據(jù)不同的防疫政策、措施進行合理預(yù)測，適用度較廣泛。

同時，本文在保持c取值為14.8的情況下，對參數(shù)q進行靈敏度分析。

改動參數(shù)q，即接觸者被追蹤到隔離的概率，這是對防疫體系成熟與否的直觀鑒定，相對成熟的防疫體系下，參數(shù)q自然相對較高，反之則較低。觀察圖4可以看到，改變參數(shù)q的值，參數(shù)變化造成的預(yù)測值變化在合理范圍內(nèi)。表明本模型在進行預(yù)測時對于參數(shù)q的變動感應(yīng)較靈敏，再次驗證本模型可以為不同的情境、防疫體系進行合理預(yù)測。

綜上可知，在不同病毒、不同的防疫防護政策下，模型反應(yīng)靈敏，具有時效性，可為不同類型傳染性病毒、不同地區(qū)及國家不同的防疫政策、措施提供模型基礎(chǔ)。

4 結(jié)語

針對COVID的傳染力研究，本文建立了改進的SEIR病毒模型，并且以醫(yī)療資源以及相關(guān)隔離政策頒布情況前后進行分析。以2020年2月13號作為分界日期，2月13號之前為前期，之后為后期，將模型實踐到北京市及武漢市，得到北京市前期基本再生數(shù)R0=2.135 0、北京市后期R0=1.973 0；武漢市前期R0=3.127 0、武漢市后期R0=2.429 0。同時，本文研究結(jié)果表明，北京市的疫情高潮會在2020年2月3日左右到來，武漢市的疫情高潮將會在2020年2月6日左右到來，該段時間病毒傳播的概率最高，要注意防控，嚴防嚴控才能有效降低病毒的傳染力。