董章功,宋 波,孟友新

(青島科技大學(xué),山東 青島 266061)

0 引 言

2019年12月以來(lái)，湖北陸續(xù)在華南海鮮市場(chǎng)發(fā)現(xiàn)不明冠狀病毒感染病例。2020年2月11日，世界衛(wèi)生組織宣布，將新型冠狀病毒肺炎命名為COVID-19[1]。據(jù)研究表明，該病毒具有人傳人的現(xiàn)象，因此，國(guó)家衛(wèi)健委發(fā)布1號(hào)通告，將該新冠肺炎納入乙類傳染病，但采取甲類傳染病預(yù)防、控制措施[2-4]。

截至2020年3月5日，全國(guó)累計(jì)確診人數(shù)80710人，其中湖北累計(jì)確診人數(shù)67592人，占全國(guó)的80%以上。面對(duì)這場(chǎng)疫情，國(guó)家科學(xué)的相關(guān)舉措顯得尤為重要。因此，對(duì)新冠肺炎感染人數(shù)的預(yù)測(cè)，有利于國(guó)家決策的制定。隨著疫情的發(fā)展，國(guó)家相關(guān)舉措的施行、人員管制、社區(qū)防疫站等措施以及新冠病毒的變異，使得原有的新冠肺炎預(yù)測(cè)數(shù)學(xué)模型存在差異。因此，本文構(gòu)建SEIR-ARIMA模型。

1 數(shù)據(jù)來(lái)源

自2020年1月20日以來(lái)，各省衛(wèi)健委官網(wǎng)每天公布感染人數(shù)、疑似人數(shù)、死亡人數(shù)等數(shù)據(jù)。各應(yīng)用平臺(tái)通過(guò)衛(wèi)健委官網(wǎng)對(duì)全國(guó)各省新冠肺炎數(shù)據(jù)進(jìn)行收集統(tǒng)計(jì)，建立新冠肺炎專題疫情地圖(http://ncov.dxy.cn/ncovh5/view/pneumonia)。本文主要是利用疫情地圖為建模提供數(shù)據(jù)支持。

2 模型與評(píng)估方法

本文所提出的SEIR-ARIMA混合模型，由傳統(tǒng)的傳染病模型SEIR和差分整合移動(dòng)平均自回歸模型ARIMA這2種模型結(jié)合組成。傳統(tǒng)的傳染病模型雖然能夠預(yù)測(cè)新冠疫情的發(fā)展趨勢(shì)，但是受外界因素的影響很大，無(wú)法根據(jù)外界因素的改變進(jìn)行模型外生變量(由模型以外的因素所產(chǎn)生，非模型本身產(chǎn)生的變量)的自我調(diào)節(jié)。因此，引入差分整合移動(dòng)平均自回歸模型ARIMA，ARIMA模型是由內(nèi)生變量(模型本身產(chǎn)生的變量)確定的，不受外生變量的影響。因此，ARIMA模型可以很好地解決SEIR的弊端。利用2種模型的預(yù)測(cè)結(jié)果，通過(guò)線性回歸能夠很好地表現(xiàn)出2種模型的預(yù)測(cè)結(jié)果與真實(shí)值之間的線性關(guān)系，通過(guò)不斷學(xué)習(xí)修改，從而提高模型預(yù)測(cè)的準(zhǔn)確性。

2.1 SEIR模型

新型冠狀肺炎是一種傳染性疾病，常見(jiàn)的傳染性模型按照傳染病類型分為SI、SIR、SEIR等模型，按照傳播機(jī)理又分為基于常微分方程、偏微分方程、網(wǎng)絡(luò)動(dòng)力學(xué)的不同類型[5]。據(jù)國(guó)家衛(wèi)健委與國(guó)家中醫(yī)藥管理局公布《新型冠狀病毒感染的肺炎診療方案(試行第四版)》文件表明：新型冠狀病毒與2003年的SARS病毒同源性達(dá)85%，都屬于冠狀病毒這個(gè)大家族[6]。有關(guān)SARS傳播動(dòng)力學(xué)研究多數(shù)采用的是SIR或SEIR模型[7]。但是據(jù)報(bào)道研究表明，新冠肺炎存在潛伏者，因此采用具有潛伏者的SEIR模型更合適。

2.1.1 SEIR模型介紹

SEIR模型是一種國(guó)際通用的傳統(tǒng)傳染病模型。該模型是研究傳染病的傳播速度、感染人數(shù)以及發(fā)展趨勢(shì)等問(wèn)題的一種非線性動(dòng)力學(xué)模型，對(duì)傳染病有效控制和預(yù)防具有指導(dǎo)的作用[8]。該模型引入潛伏人群的概念，將人群分為4類，分別是易感人群、潛伏人群、感染人群、康復(fù)人群[9]。其中潛伏人群指接觸過(guò)感染者，但暫無(wú)能力傳染給其他人的人群[10]。SEIR模型通過(guò)這4類人群之間的關(guān)聯(lián)關(guān)系建立微分方程，是一種典型的傳染病模型。

2.1.2 SEIR模型建立

如圖1所示，SEIR模型首先分別用符號(hào)S、E、I、R來(lái)表示易感人群、潛伏人群、感染人群、康復(fù)人群這4種人群的數(shù)量[11]。假設(shè)一個(gè)區(qū)域內(nèi)總?cè)藬?shù)為N，即N=S+E+I+R，每天感染者接觸的人數(shù)為P，健康人比例為S/N，其中易感人群按照比例A轉(zhuǎn)化為潛伏人群，潛伏人群按照比例B轉(zhuǎn)化為感染人群，感染人群按照比例C轉(zhuǎn)化為康復(fù)人群[12]。在此模型上可以建立如下微分方程：

圖1 SEIR模型

(1)

(2)

(3)

(4)

2.2 ARIMA模型框架

ARIMA模型是一種時(shí)間序列預(yù)測(cè)方法，主要捕獲歷史數(shù)據(jù)的波動(dòng)，描述當(dāng)前值與歷史值之間的關(guān)系[13-14]，適用于含時(shí)間序列的新冠肺炎預(yù)測(cè)，而且其自身模型是借助內(nèi)生變量構(gòu)建，不受外生變量影響，可很好地解決外界因素的干擾。建模步驟如圖2所示。

圖2 ARIMA建模步驟圖

2.2.1 平穩(wěn)性檢測(cè)

平穩(wěn)時(shí)間序列是ARIMA模型預(yù)測(cè)的重要前提，通過(guò)時(shí)間序列上的平均值反映樣本數(shù)據(jù)的波動(dòng)情況，通過(guò)標(biāo)準(zhǔn)差反映樣本數(shù)據(jù)的離散程度，從而做出平穩(wěn)性假設(shè)。然后通過(guò)單位根檢驗(yàn)(ADF檢驗(yàn))在數(shù)據(jù)層面直觀地驗(yàn)證假設(shè)。

平均值公式：

(5)

標(biāo)準(zhǔn)差公式：

(6)

其中，xi(i=1,2,…,n)表示t時(shí)刻感染人數(shù)。

2.2.2 數(shù)據(jù)平穩(wěn)性處理

如果該時(shí)間序列數(shù)據(jù)存在離散性，無(wú)法順著現(xiàn)有的形態(tài)延續(xù)下去[15-16]，可能會(huì)增加新冠肺炎預(yù)測(cè)值與真實(shí)值之間的誤差，這也是對(duì)新冠肺炎相關(guān)數(shù)據(jù)進(jìn)行平穩(wěn)性檢驗(yàn)與處理的重要原因之一。數(shù)據(jù)離散可以通過(guò)采用差分的方法，使離散的時(shí)間序列數(shù)據(jù)平穩(wěn)化[17]。

差分運(yùn)算過(guò)程：

一階差分：

Δyt=yt+1-yt

(7)

其中，yt+1、yt分別表示時(shí)間序列在t+1與t的值。

二階差分：

Δ(Δyt)=Δ(yt+1-yt)=Δyt+1-Δyt

(8)

二階差分是在一階差分的基礎(chǔ)上進(jìn)行運(yùn)算，其中，Δyt+1、Δyt分別表示時(shí)間序列y在t+1與t的一階差分值。

2.2.3 ARIMA運(yùn)算過(guò)程

(9)

其中,φ表示AR的系數(shù)，θ表示MA的系數(shù),p為自回歸項(xiàng)，q為移動(dòng)平均項(xiàng)數(shù)[18],e為誤差。p、q暫時(shí)未知，需要通過(guò)自相關(guān)函數(shù)ACF(Autocorrelation Function)和偏自相關(guān)函數(shù)PACF(Partial Autocorrelation Function)確定。

1)自相關(guān)函數(shù)ACF，不僅受其他變量的影響，而且還與其他變量存在相關(guān)的關(guān)系，所以ACF里含有其他變量影響因子[19]。

(10)

2)偏自相關(guān)函數(shù)PACF，排除其他變量的相關(guān)影響下，k時(shí)間內(nèi)yt+k和yt之間的關(guān)系。

(11)

2.3 SEIR-ARIMA模型構(gòu)建

SEIR模型可以較好地預(yù)測(cè)疫情的發(fā)展趨勢(shì)，但是對(duì)于初始的傳染概率等外在因素也只是估算的，在國(guó)家采取疫情防控的措施下，傳染概率也會(huì)隨之降低，因此SEIR無(wú)法根據(jù)歷史的波動(dòng)數(shù)據(jù)，進(jìn)行模型的自我調(diào)節(jié)。ARIMA模型雖然是基于時(shí)間序列預(yù)測(cè)的模型之一，但是具有由內(nèi)生變量所決定不受外界因素干預(yù)的優(yōu)點(diǎn)，很好地彌補(bǔ)了SEIR模型的弊端。因此，為了更加準(zhǔn)確地預(yù)測(cè)疫情的發(fā)展趨勢(shì)，本文通過(guò)對(duì)SEIR模型預(yù)測(cè)結(jié)果和ARIMA模型預(yù)測(cè)結(jié)果進(jìn)行線性回歸建模，形成了SEIR-ARIMA混合模型。SEIR-ARIMA混合模型的優(yōu)勢(shì)是既學(xué)習(xí)到了一般傳染病的發(fā)展趨勢(shì)，又避免后期因管控措施、病毒變異等影響，造成外生變量的改變，而進(jìn)行的分階段劃分分析，可以根據(jù)歷史數(shù)據(jù)的發(fā)展趨勢(shì)，統(tǒng)籌地預(yù)測(cè)接下來(lái)疫情發(fā)展的優(yōu)點(diǎn)。SEIR-ARIMA模型是通過(guò)SEIR模型和ARIMA模型預(yù)測(cè)結(jié)果線性回歸建模，而線性回歸是研究因變量與自變量之間關(guān)系[20]，在線性回歸中，自變量可以是連續(xù)的，也可以是離散的[21]。通過(guò)一元線性方程或者多元線性方程，找到SEIR模型預(yù)測(cè)值與ARIMA模型預(yù)測(cè)值之間的內(nèi)部關(guān)系，形成最佳的線性擬合直線[22]。本文認(rèn)為SEIR模型和ARIMA模型存在線性相關(guān)性，使用SEIR-ARIMA模型進(jìn)行預(yù)測(cè)步驟如下：

1)SEIR模型預(yù)測(cè)。

通過(guò)對(duì)2020年1月20日至2020年2月10日疫情數(shù)據(jù)進(jìn)行分析，結(jié)合相關(guān)研究以及參數(shù)的優(yōu)化，得到其預(yù)測(cè)值f。

2)ARIMA模型預(yù)測(cè)。

將疫情歷史數(shù)據(jù)作為輸入值，檢驗(yàn)時(shí)間序列數(shù)據(jù)是否平穩(wěn)化，采用ADF驗(yàn)證假設(shè)。若不平穩(wěn)就通過(guò)一元差分、二元差分等方法，直到數(shù)據(jù)平穩(wěn)化[23]。通過(guò)計(jì)算ACF和PACF得到p、q的值，建立ARIMA模型。訓(xùn)練樣本集，優(yōu)化模型，得到預(yù)測(cè)值e。

3)SEIR-ARIMA模型構(gòu)建。

SEIR-ARIMA模型構(gòu)建以上述所得的2種模型的預(yù)測(cè)值f、e為輸入值，真實(shí)數(shù)據(jù)為目標(biāo)值，采用線性回歸建模，計(jì)算SEIR-ARIMA混合模型的預(yù)測(cè)值。

(12)

求解系數(shù)k的過(guò)程：

(13)

即預(yù)測(cè)值與真實(shí)值的倍數(shù)，求解偏移數(shù)據(jù)b，

(14)

2.4 評(píng)估預(yù)測(cè)方法

本文采用4個(gè)指標(biāo)評(píng)估預(yù)測(cè)數(shù)據(jù)的精確度，分別是平均誤差ME(反映預(yù)測(cè)值與真實(shí)值之間的誤差)、平均絕對(duì)誤差MAD(體現(xiàn)真實(shí)值與預(yù)測(cè)值的真實(shí)誤差)、平均絕對(duì)百分比誤差MAPE(反映誤差大小的相對(duì)值)、決定系數(shù)R2(反映模型的可靠程度[24-25])。

(15)

(16)

(17)

(18)

3 模型的預(yù)測(cè)與評(píng)估

3.1 SEIR-ARIMA模型構(gòu)建

3.1.1 SEIR模型構(gòu)建

SEIR模型中未知變量為S、E、I、R、A、B、C、N。其中，假設(shè)易感人群為全部人群易感，湖北人口統(tǒng)計(jì)局統(tǒng)計(jì)湖北人口總數(shù)據(jù)為5700萬(wàn)，因此，S=N=57000000,初始潛伏者為0，感染者為1，康復(fù)者為0，因此，E=0、I=1、R=0。其次，假設(shè)新冠肺炎的傳染概率與2003年SARS相同為0.6，因此A=0.6。據(jù)相關(guān)研究報(bào)道，潛伏人群轉(zhuǎn)換為感染人群至少需要2 d～14 d，故潛伏者轉(zhuǎn)化為感染者概率為B=1/((2+14)/2)=0.125?？祻?fù)率為治愈人數(shù)/累計(jì)患病人數(shù)，根據(jù)2020年1月20日至2020年2月10日數(shù)據(jù)可以得到其平均值為C=0.143。根據(jù)以上數(shù)據(jù)建立相對(duì)應(yīng)的SEIR模型。

3.1.2 ARIMA模型構(gòu)建

1)判定數(shù)據(jù)是否平穩(wěn)化。

通過(guò)平均值和標(biāo)準(zhǔn)差進(jìn)行平穩(wěn)化判定。當(dāng)數(shù)據(jù)平穩(wěn)，在數(shù)據(jù)序列上，其平均值不隨時(shí)間有較大變化[26-27]。如圖3所示，其平均值隨時(shí)間呈現(xiàn)遞增的趨勢(shì)。因此，假設(shè)其不具有平穩(wěn)性。

圖3 平穩(wěn)性檢驗(yàn)圖

為了進(jìn)一步驗(yàn)證數(shù)據(jù)不平穩(wěn)，采取ADF檢驗(yàn)。ADF檢驗(yàn)又稱為單位根檢驗(yàn)，如果存在單位根，則序列平穩(wěn),否則，序列不平穩(wěn)[28]。所以，ADF檢驗(yàn)的假設(shè)就是存在單位根，如果得到的顯著性檢驗(yàn)統(tǒng)計(jì)量小于3個(gè)置信度(10%，5%，1%)，則對(duì)應(yīng)有(90%，95%，99%)的把握來(lái)拒絕原假設(shè)[29]。

如圖4所示，其中Test Statistic表示檢驗(yàn)值，p-value表示檢驗(yàn)結(jié)果，critical value表示置信度,看到Test Statistic的值-1.364767要大于置信度10%，所以無(wú)法拒絕原假設(shè)，另外，p-value的值也很大。因此，時(shí)間序列不平穩(wěn)。

圖4 ACF檢驗(yàn)置信度圖

2)平穩(wěn)化處理。

如圖5所示，一階差分之后數(shù)據(jù)并未平穩(wěn)，通過(guò)二階差分后，數(shù)據(jù)基本平穩(wěn)。

圖5 差分對(duì)比圖

3)模式識(shí)別，求p、q值。

通過(guò)建立ACF和PACF得出p、q值。

如圖6所示，根據(jù)表1的取值規(guī)則可見(jiàn)，ACF的2階衰減趨于0，PACF在2階后處于置信區(qū)間，因此p=2,ACF在1階后處于置信區(qū)間，PACF在1階后衰減趨于0，因此q=1。

圖6 ACF和PACF置信區(qū)間圖

表1 p和q取值規(guī)則

3.1.3 SEIR-ARIMA模型構(gòu)建

將湖北2020年1月20日至2020年2月10日疫情數(shù)據(jù)作為訓(xùn)練集，輸入樣本訓(xùn)練ARIMA模型，通過(guò)訓(xùn)練好的ARIMA模型，預(yù)測(cè)90 d內(nèi)感染人數(shù)。SEIR模型通過(guò)輸入傳染概率、最初感染者、康復(fù)概率等，構(gòu)建SEIR模型，預(yù)測(cè)90 d內(nèi)數(shù)據(jù)變化。最后，構(gòu)建SEIR-ARIMA混合模型是通過(guò)ARIMA模型預(yù)測(cè)值和SEIR模型預(yù)測(cè)值作為輸入，真實(shí)值為目標(biāo)，構(gòu)建ARIMA模型和SEIR模型預(yù)測(cè)結(jié)果與真實(shí)數(shù)據(jù)之間的線性模型，完成SEIR-ARIMA混合模型構(gòu)建。SEIR-ARIMA混合模型通過(guò)不斷學(xué)習(xí)修改預(yù)測(cè)值與真實(shí)值之間的線性關(guān)系，形成最佳的線性擬合直線，從而提高模型預(yù)測(cè)的準(zhǔn)確性。

3.2 模型評(píng)估與對(duì)比

本文使用2020年1月20日至2020年2月10日湖北新冠疫情數(shù)據(jù)，對(duì)湖北疫情2020年1月20日至2020年4月20日進(jìn)行預(yù)測(cè)。從圖7可以看出，單一的ARIMA模型預(yù)測(cè)效果較差，與真實(shí)值差距較大，SEIR模型預(yù)測(cè)現(xiàn)有確診人數(shù)發(fā)展趨勢(shì)相對(duì)較好。但是，與總的確診人數(shù)偏差較大，預(yù)測(cè)效果也存在較大誤差。

圖7 湖北省新冠肺炎預(yù)測(cè)模型對(duì)比圖

為了對(duì)SEIR-ARIMA混合模型進(jìn)行評(píng)估，本文引入了常用于新冠疫情預(yù)測(cè)的LSTM模型、Logistic模型。LSTM模型是一種長(zhǎng)短期記憶網(wǎng)絡(luò)模型，王瑞等[30]對(duì)新冠肺炎疫情研究表明，LSTM模型對(duì)國(guó)內(nèi)新冠疫情發(fā)展趨勢(shì)的預(yù)測(cè)相對(duì)較好。Logistic模型作為傳統(tǒng)的機(jī)器學(xué)習(xí)模型的一種，在國(guó)內(nèi)新冠肺炎預(yù)測(cè)方面也有較好的預(yù)測(cè)效果，而且，其具有訓(xùn)練高效、準(zhǔn)確率好的優(yōu)點(diǎn)，一般可以作為模型訓(xùn)練的基準(zhǔn)，用來(lái)衡量其他算法的性能。為了更全面地證明SEIR-ARIMA模型結(jié)合的優(yōu)勢(shì)，本文實(shí)現(xiàn)了SEIR-LSTM模型和SEIR-Logistic模型并對(duì)它們進(jìn)行比較。如圖7所示，LSTM和Logistic模型、SEIR-Logistic模型預(yù)測(cè)的趨勢(shì)與真實(shí)趨勢(shì)基本一致，但是預(yù)測(cè)人數(shù)與真實(shí)值差距較大。而SEIR-LSTM模型預(yù)測(cè)人數(shù)接近真實(shí)值，但是相比于SEIR-ARIMA模型，SEIR-ARIMA模型更加接近真實(shí)值，預(yù)測(cè)效果更加精確。綜上所述可以得出，SEIR-ARIMA混合模型預(yù)測(cè)相比于LSTM模型、Logistic模型、ARIMA模型、SEIR模型、SEIR-LSTM和SEIR-Logistic模型，預(yù)測(cè)效果更好。如表2所示，SEIR-ARIMA模型的4個(gè)評(píng)估指標(biāo)平均誤差(ME)、平均絕對(duì)誤差(MAD)、平均絕對(duì)百分比誤差(MAPE)、決定系數(shù)(R2)分別為-574、931、0.03、0.90。

表2 模型預(yù)測(cè)誤差對(duì)比

雖然SEIR-ARIMA模型具有對(duì)歷史數(shù)據(jù)的學(xué)習(xí)能力，但是為了驗(yàn)證在不同時(shí)間段、不同地點(diǎn)上獲得的數(shù)據(jù)仍然具有較好的預(yù)測(cè)效果，選取2021年于7月20日在南京爆發(fā)的疫情數(shù)據(jù)做預(yù)測(cè)，比較SEIR-ARIMA預(yù)測(cè)效果。

如圖8所示，可以看出SEIR-ARIMA預(yù)測(cè)效果較好，SEIR-LSTM次之。

圖8 南京新冠肺炎預(yù)測(cè)模型對(duì)比圖

綜上所有對(duì)比實(shí)驗(yàn)可以得出，SEIR-ARIMA預(yù)測(cè)效果較好，SEIR-LSTM次之，SEIR-Logistic略好于其他單一模型，其他單一模型在長(zhǎng)期的疫情發(fā)展趨勢(shì)預(yù)測(cè)是較好的，但是其預(yù)測(cè)值還是與真實(shí)值有較大差距。因此可以看出，SEIR-ARIMA模型不管與其他模型進(jìn)行比較，還是不同時(shí)間、不同地點(diǎn)的預(yù)測(cè)效果的比較，SEIR-ARIMA模型都具有更好的預(yù)測(cè)效果。

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

針對(duì)新冠疫情趨勢(shì)的預(yù)測(cè)，本文提出了一種SEIR-ARIMA混合模型預(yù)測(cè)的方法。通過(guò)傳統(tǒng)的傳染病模型SEIR建模，刻畫(huà)了傳染病一般的發(fā)展規(guī)律，同時(shí)結(jié)合了ARIMA模型，使整個(gè)混合模型在管控措施下，傳染概率等內(nèi)因改變的條件下，可以根據(jù)歷史因素，在模型內(nèi)部進(jìn)行調(diào)節(jié)。通過(guò)與Logistic模型、LSTM模型、單一的SEIR模型、ARIMA模型、SEIR-Logistic模型和SEIR-LSTM模型在不同時(shí)間、不同地點(diǎn)預(yù)測(cè)結(jié)果對(duì)比，可以看出，本文提出的SEIR-ARIMA混合模型有更好的預(yù)測(cè)效果，對(duì)國(guó)家未來(lái)預(yù)防其他類型傳染病具有很好的應(yīng)用價(jià)值。但是，本文提出的SEIR-ARIMA混合模型并未對(duì)康復(fù)者是否還會(huì)再感染等因素進(jìn)行實(shí)驗(yàn)。基于本模型具備的優(yōu)良特性，下一步可以通過(guò)研究康復(fù)者的感染情況，證明能否達(dá)到預(yù)期縮小誤差的效果。