李榮庭，段 鵬，胡 瑞，范蓮靜

(云南民族大學 數(shù)學與計算機科學學院，云南 昆明 650000)

根據(jù)以往的經驗，我國冬季流感流行時間一般從12月份開始一直持續(xù)到來年的3—5月.經中國疾病預防控制中心病毒預防控制所首席專家武桂珍的介紹由于流感和新冠肺炎均可經呼吸道飛沫傳播及接觸傳播[1-4]，其均受氣候以及人際接觸等的影響較大：秋冬季節(jié)的低溫天氣會使呼吸道排入外界的病毒存活期延長，且在寒冷刺激下，人的抵抗力相對較弱，病毒趁虛而入，易于形成傳播；人員聚集，則是疫情擴大的另一必要條件[5].進入秋冬季之后，全國更大范圍的復工和學生開學，會形成廣泛的、不同程度的人員聚集，且秋冬季節(jié)里人們室內活動相對頻繁，室內通風情況不及夏季.加之目前全球范圍內許多其他國家的新冠肺炎疫情仍較嚴峻，尚未得到有效控制等，以上這些因素均可能導致我國在今年冬季面臨新冠肺炎疫情的反彈，其反彈時間也將很可能與冬季流感流行時間重疊[6].倘若我們能在流感發(fā)現(xiàn)早期及時預測疾病的增長趨勢，很有可能我們就可以采取相應的措施，及時止損.同時權威精準的預測也可以使人們盡早的把預防流感重視起來，從而進一步減少流感蔓延的規(guī)模，從而也從側面對新冠疫情的防控起到了積極作用.本課題雖然是以甲型H1N1流感為案例進行預測分析，但搭建構造的模型，其實亦適用與其他類似的傳染病，例如肺結核、新冠肺炎等.對于后續(xù)的相關研究，起到了一定的借鑒作用.從另一個角度出發(fā)，也可以用此模型對未來可能出現(xiàn)的潛在的傳染病威脅進行排查.

目前現(xiàn)代醫(yī)學技術取得重大發(fā)展，但諸多傳染性疾病仍是人類社會向前發(fā)展的重要阻力之一.查閱全世界范圍內已有的流感相關數(shù)據(jù)進行分析研究，多數(shù)預測模型應用多元線性回歸、LASSO回歸以及Ridge回歸模型結合相關檢索詞數(shù)據(jù)進行建模分析，探討回歸模型與流感疫情預測的相關性與可行性[7].常見的傳染病預測模型還有灰色模型、Markov模型、人工神經網(wǎng)絡模型、通徑分析模型等[8].但這都需要大量數(shù)據(jù)作為支持才能得到理想的結果，而針對突如其來的大規(guī)模傳染病，其預測數(shù)據(jù)量小、預測結果亟待使用的特點對預測工作的準確度和速度提出了雙重要求.基于上述要求，實驗選取少樣本數(shù)據(jù)集，利用數(shù)據(jù)變化的自回歸過程建立模型，使得搭建的模型具有較高的可信度及一定的實用性.

1 模型介紹

1.1 ARIMA模型

時間序列是利用統(tǒng)計學的基本原理,通過對數(shù)據(jù)的采集選用模型以近似估計,利用模型分析揭示數(shù)據(jù)的內在特性[9],達到推測發(fā)展趨勢規(guī)律的目的[10].自回歸差分移動平均模型(ARIMA)是常用的時間序列之一.該模型主要分為三個過程，平穩(wěn)過程(I)、自回歸過程(AR)和移動平均過程(MA).該模型反映了時間序列過去與現(xiàn)在，未來與現(xiàn)在之間的關系，適用于短期預測，需要的數(shù)據(jù)量小、參數(shù)少，這與少樣本數(shù)據(jù)集訓練有很高的鍥合度.基于上述特點，ARIMA模型在眾多領域的應用都有著不錯的變現(xiàn)，尤其是在與時間相關的預測方面[11].

1.2 SVM模型

支持向量機(support vector machines, SVM)是一種二分類模型，它的基本模型是定義在特征空間上的間隔最大的線性分類器，間隔最大使它有別于感知機；SVM還包括核技巧，這使它成為實質上的非線性分類器.SVM的的學習策略就是間隔最大化，可形式化為一個求解凸二次規(guī)劃的問題，也等價于正則化的合頁損失函數(shù)的最小化問題.SVM的的學習算法就是求解凸二次規(guī)劃的最優(yōu)化算法[12].相比于傳統(tǒng)的優(yōu)化算法，例如人工蜂群算法、蟻群算法、粒子群算法，SVM的優(yōu)化對象不僅可以是固定參數(shù)也可以是其他數(shù)值.ARIMA模型里的關鍵參數(shù)都是通過模型自身的差分值進行自調整的，無需再優(yōu)化，對于下文中動態(tài)殘差的優(yōu)化，SVM算法較為合適.

2 基于ARIMA-SVM模型的預測算法

2.1 ARIMA模型搭建

本文通過將2009年美國甲型流感新增病例數(shù)據(jù)與傳統(tǒng)ARIMA模型結合[13]，并在此基礎上增加模型評估環(huán)節(jié)以確定階數(shù)，繪制修改后的流程圖如下.

圖1 ARIMA模型構建流程示意圖

傳統(tǒng)ARIMA模型在階數(shù)選擇方面，一般采用置信度區(qū)間能夠完全覆蓋的數(shù)值，但是往往沒有容忍度的參數(shù)選擇會使模型錯過最佳參數(shù)，故本文加入模型評估，為確定階數(shù)增加一些彈性，使得模型在可接受范圍內忽略一些誤差從而找到最佳參數(shù).

2.1.1 平穩(wěn)數(shù)據(jù)

流感的傳播具有聚集性和爆發(fā)性，所以針對數(shù)據(jù)的區(qū)域特點和時間特點，考慮到模型對數(shù)據(jù)平穩(wěn)性的要求[14]，將美國2009年甲型H1N1新增病例數(shù)數(shù)據(jù)利用差分法進行平穩(wěn)，計算特征點數(shù)據(jù)方差，找到合適的差分階數(shù)，下圖2為差分比較圖.

通過圖2內容可以看出，一階差分后數(shù)據(jù)相對穩(wěn)定.又根據(jù)差分后的過濾掉空值缺失數(shù)據(jù)繪制效果展示圖.由圖3可知差分后的平穩(wěn)性要比原數(shù)據(jù)好很多，標準差反映出差分后的數(shù)據(jù)的離散程度在可接受的范圍[15],故確定ARIMA模型的參數(shù)d=1.

圖2 一、二階差分效果展示圖 圖3 差分數(shù)據(jù)效果展示圖

2.1.2 自相關圖和偏自相關圖的分析

本文利用自相關函數(shù)(ACF)和偏自相關函數(shù)(PACF)結合赤池準則(AIC)和貝葉斯信息準則(BIC)進行模型評估確定相應的階數(shù)p和q，下圖4中的陰影表示置信區(qū)間，可以看出不同階數(shù)自相關性的變化情況，從而選出p值和q值.p值代表預測模型中采用的時序數(shù)據(jù)本身的滯后數(shù)(lags),也叫做AR(auto-regressive)項.q值代表預測模型中采用的預測誤差的滯后數(shù)(lags)，也叫做MA(moving average)項[9].

針對p值和q值得選取，本模型提出利用Matlab軟件對自相關函數(shù)圖像和偏自相關函數(shù)圖像進行拉普拉斯(Laplace)濾波處理抽離置信區(qū)間和階數(shù)取值特征點，然后用cftool工具箱中的傅里葉函數(shù)和高斯函數(shù)分別進行擬合得到置信區(qū)間邊界函數(shù)s(x)和階數(shù)取值函數(shù)F(x).(由于置信度區(qū)間關于x軸大致對稱，所以本文選取x軸上方的部分，抽象出置信度區(qū)間邊界函數(shù)s(x)).

圖4 自相關函數(shù)和偏自相關函數(shù)選取p、q值

(1)

p值和q值的確定原則是，圖像中置信區(qū)間所涵蓋的最后一個超出其范圍的滯后數(shù)的取值.但在原有的模型中p值和q值的選取是沒有容忍度的，即只有完全符合置信區(qū)間要求才會被選取.本文在此基礎上提出利用容忍度函數(shù)，如公式(1)，篩選符合條件的最小的p值和q值以保證不會使模型陷入局部最優(yōu)而導致忽略全局最優(yōu).其中xi是待定p或q.根據(jù)上述圖像和公式本確定的p值為3，q值為2.

2.1.3 AR模型與MA模型

根據(jù)之前確定的差分階數(shù)還有p值和q值，建立如下數(shù)學模型：

(2)

(3)

(4)

其中，yt是當前值，μ是常數(shù)項，p是階數(shù)，γi是自相關系數(shù)，εt是誤差，θi能夠消除波動的最佳參數(shù).從實際情況出發(fā)由于自回歸模型本身的限制，在數(shù)據(jù)選取過程中既要保證穩(wěn)定性，還要保證自相關性，如果自相關系數(shù)φi小于0.5，則不宜采用[13].

2.1.4 模型殘差檢驗

基于以上模型，進行模型殘差檢驗，觀測ARIMA模型方差是否為常數(shù)的正態(tài)分布，以確保模型的可行性和合理性[16].并通過對其準確的預測采取下一步相應的預防措施.在選擇了合適的參數(shù)之后，對ARIMA模型的殘差進行初步預測，如圖5.

2.2 SVM模型殘差優(yōu)化

2.2.1 SVM優(yōu)化分析

根據(jù)上文中給出的相關圖像，不難看出預測數(shù)據(jù)與實際數(shù)據(jù)還存在這一定的誤差.故針對預測中存在的誤差本實驗欲通過支持向量機(Support Vector Machine,SVM)進行優(yōu)化[17].根據(jù)上文，本實驗選取原模型殘差作為對象進行優(yōu)化，并綜合兩個模型的實驗結果得出混合預測值作為最終預測值[11]，圖6為優(yōu)化過程流程圖.對于圖6混合預測值的詳細解釋參照下文殘差數(shù)據(jù)導入部分.

圖5 殘差預測效果圖展示 圖6 加入SVM模型優(yōu)化后的ARIMA模型流程示意圖

具體優(yōu)化的操作步驟主要分為3個部分，第1個部分根據(jù)之前ARIMA模型預測數(shù)據(jù)與實際數(shù)據(jù)做差得到新的誤差數(shù)據(jù)并利用numpy第3方庫進行矩陣化為誤差數(shù)據(jù)貼上時序標簽，如公式(5)；第2部分將得到的新的數(shù)據(jù)再做滑動處理擴增數(shù)據(jù)集，滑動次數(shù)按照參照之前ARIMA模型得出的，這樣就得到了新的矩陣，共計63行5列如圖7所示；第3部分將實際數(shù)據(jù)的誤差作為標簽利用SVM模型進行有監(jiān)督的學習，具體操作流程參下.

Eerror=yreal-yarima

(5)

2.2.2 導入SVR模塊

由于當前環(huán)節(jié)需要處理的數(shù)據(jù)是預測值與實際值的誤差，那么自變量就可以定做是出現(xiàn)新增病例之后的每一天，與此同時因變量也自然而然成了誤差本身.在此基礎上就不難推斷在進行優(yōu)化之前需要導入SVR模塊[15]，SVR()就是SVM算法來做回歸用的方法(即輸入標簽是連續(xù)值的時候要用的方法)，通過以下語句來確定SVR的模式(選取比較重要的幾個參數(shù)進行測試).

1) kernel：核函數(shù)的類型

一般常用的有rbf函數(shù)，linear函數(shù)，poly函數(shù)，如圖7所示，將數(shù)據(jù)導入函數(shù)擬合圖像發(fā)現(xiàn)rbf函數(shù)圖像與實際數(shù)據(jù)擬合度最高，故最佳核函數(shù)為rbf函數(shù)[18].

圖7 誤差數(shù)據(jù)最終形式 圖8 核函數(shù)與數(shù)據(jù)貼合度參考圖

2)C:懲罰因子

C表征對離群點的重視程度，C大小與重視程度成正比，與對誤差分類的懲罰成正比，與margin成反比.當C趨近無窮的時候，表示不允許分類誤差的存在，容易過擬合；當C趨于0時，表示我們不再關注分類是否正確，容易欠擬合.

3)gamma

gamma是rbf函數(shù)的核系數(shù)且gamma的值必須大于0，其作用是保證模型不會過擬合.

綜上所述SVR模塊的幾個重要參數(shù)可以由上述規(guī)則基本確定：

SVR(C=10, cache_size=200, coef0=0.0, degree=2, epsilon=0.1, gamma=‘scale’, kernel=‘rbf’, max_iter=-1, shrinking=True, tol=0.001, verbose=False)

2.2.3 殘差值數(shù)據(jù)導入

首先要將每日的預測值，轉化為矩陣形式，方便將數(shù)據(jù)導入SVM.由于目前本實驗無法確定優(yōu)化后的效果優(yōu)秀與否，故要同時考慮原始ARIMA模型的預測值yarima與ARIMA-SVM模型的預測值yarima-svm，如公式4.然后將總預測值減去實際值，得到殘差值E，即公式5.在此基礎上將殘差值導入SVR模塊，為了篩選出最合適本數(shù)據(jù)集的C和gamma，還要引入一個新的變量[12]，記作最小殘差Emin.一方面初始Emin是為了拿來和E進行比較，縮小檢索C和gamma的范圍，保證篩選范圍不會太大，導致程序崩潰.另一方面，并將循環(huán)中出現(xiàn)的最小E賦值給Emin，循環(huán)最后找到Emin對應的C和gamma，記作Cbest和gammabest流程圖如圖8.

ysum=yarima+yarima-svm，

(6)

(7)

圖9 SVM參數(shù)調整流程圖

3 實驗與分析

實驗數(shù)據(jù)采用WHO官網(wǎng)公開的美國2009年甲型H1N1流感數(shù)據(jù)集，實驗選擇配置為AMD Ryzen 5中央處理器(主頻 3.4 GHz)，內存為8G的電腦為硬件平臺，在 Windows 10，64位操作系統(tǒng)下利用anaconda(python3.7.2)進行實驗.將上述過程中得到的最佳最佳參數(shù)分別導入到ARIMA和SVM模型中，便可以得到較為準確的預測結果了.

3.1 數(shù)據(jù)集

本文使用的數(shù)據(jù)來自世界衛(wèi)生組織WHO(World Health Organization)公開下載數(shù)據(jù)集.對數(shù)據(jù)進行了如下處理：①考慮到地方防疫措施和疫苗推廣對疫情的影響，本實驗選取2009年疫情初期4月24日到7月16日共計 83 d，美國甲型H1N1新增病例數(shù)信息(其中有 4 d 空值)；②數(shù)據(jù)清洗，取空值前后2天的平均值填充;③將數(shù)據(jù)集按4∶1比例分割為訓練集與驗證集.

3.2 基于ARIMA和SVM模型的預測

3.2.1 ARIMA模型預測結果

將2009年4月24日以來到7月6日之間美國甲型H1N1流感新增病例數(shù)數(shù)據(jù)導入到已經初步搭建好的時間序列模型，并利用python軟件繪制出模型差分值圖像，如圖9所示.并將差分后的數(shù)據(jù)導入傳統(tǒng)ARIMA模型預測結果如圖10.

圖10 ARIMA模型數(shù)據(jù)的差分值 圖11 ARIMA模型預測效果示意圖

很明顯預測結果存在較大偏差，其中原因有可能為以下2點：

1) 在偏離處當?shù)叵嚓P部門采取了有效的防控措施，使得實際新增病例數(shù)遠遠低于預測值.

2) 模型本身仍存在一定缺陷，需要進行進一步的優(yōu)化.

基于上述原因，本文提出了利用SVM模型對殘差值進行改進的方案.

3.2.2 ARIMA-SVM模型優(yōu)化結果

根據(jù)上文中SVM對殘差值的優(yōu)化，將新的預測模型預測值、舊的預測模型預測值、真實的新增病例利用python中的畫圖工具繪制甲流新增病例變化趨勢折線圖12用于進行相互比對.

3.2.3 模型結果評估

為了考察新模型預測效果，本實驗引入平均絕對百分誤差(MAPE)和均方根誤差(RMSE)2個指標來測定模型預測的精度.其中MAPE表示預測結果每天的預測結果與當天實際情況的偏離指數(shù)，RMSE表示總體預測結果與實際情況的偏離指數(shù).故MAPE和 RMSE越小，預測精度越高.公式與對比見表1：

(8)

(9)

圖12 各模型預測結果對比圖

表1 新舊模型殘差評估表

通過上面的表格和對比圖，不難發(fā)現(xiàn)改進后的模型精準度大大增加了，疫情病發(fā)的 50 d 內，預測數(shù)據(jù)與實際數(shù)據(jù)基本重合，準確率非常高，50 d 之后預測與實際情況開始出現(xiàn)偏差，但走勢相對吻合.

查閱相關資料2009年6月11日(即疫情出現(xiàn) 50 d 以后)，WHO將甲型H1N1警戒級別從5級提高到最高級6級.這是WHO 40年來第1次把傳染病警戒級別升至最高級別.世衛(wèi)生組織確認全球75個國家和地區(qū)，共確診27737例患者，死亡141例.至此美國人民危機意識已經防疫意識有了本質上的覺醒，自覺在家隔離，相關部門加強衛(wèi)生管理.猜測這可能是導致疫情增長狀況出現(xiàn)拐點與預測數(shù)據(jù)產生偏差的本質原因.正如王菊[19]博士所說，獨立、透明、有效和及時的信息流對于國際社會控制不斷出現(xiàn)的疾病至關重要.隨著甲型HIN1流感的細節(jié)開始被披露，人們應對農業(yè)尤其對墨西哥的豬養(yǎng)殖業(yè)進行一個更加全面的評估，這對于為未來吸取經驗教訓至關重要.當然，人們目前的注意力集中在國際社會公共衛(wèi)生反應，為疫情作準備意味著為意外作準備，在不確定的情況下作出快速和有效的反應.正如應對禽流感的經驗教訓所表明的那樣，這可能需要采取更多的措施，而不只是從上至下的“積極和有力”的技術應對措施.本文的實證分析為上述衛(wèi)生結論提供有力證明.

4 結語

本研究基于2009年在美國和墨西哥等地大面積爆發(fā)的甲型H1N1流感新增病例數(shù)據(jù)，展開關于傳染性肺炎的研究，通過ARIMA -SVM模型預測得知甲型H1N1流感新增病例人數(shù)明顯高于實際值，這可能與當年當?shù)馗鲊t(yī)療工作者展開的防治工作有關，本研究結果也為甲型H1N1流感防治措施的有效性提供了間接證據(jù),為之后的傳染病預防工作提供參考，從而達到及時判斷、抉擇和控制的目的.

根據(jù)4月24日到7月16日數(shù)據(jù)本實驗針對傳統(tǒng)ARIMA模型進行改進：①提出利用圖像擬合函數(shù)得到了最佳p值和q值；②利用ARIMA預測值與實際值之間的殘差序列作為樣本集，導入支持向量機(SVM)預測殘差，建立ARIMA-SVM混合模型；③在SVM優(yōu)化過程中加入?yún)?shù)自循環(huán)得到最佳懲罰因子和支持向量數(shù)；④通過比較新的評估參數(shù)MAPE和RMSE，發(fā)現(xiàn)支持向量機有效擬合了原始序列的非線性部分，ARIMA-SVM混合模型的預測效果優(yōu)于單一ARIMA模型.

針對本模型可以做的進一步研究工作：①優(yōu)化ARIMA模型的數(shù)據(jù)平穩(wěn)過程；②本模型應對噪聲的自適應性較差，如果在預測的過程中數(shù)據(jù)被其它人為因素影響會導致預測誤差變大，故在去噪方面可以進一步改進；③ARIMA模型研究的對象只與時間有關，但在現(xiàn)實社會中，流行性傳染病的傳播受醫(yī)療水平、文化、經濟等多種因素的影響，在本模型的基礎上可以加入人口演化，醫(yī)療隔離等仿真模型，使得模型精度進一步提高.