高桂玲，呂錫宏，孫中興

上海市松江區(qū)疾病預防控制中心，上海 201620

新型冠狀病毒肺炎（簡稱新冠肺炎）為新發(fā)急性呼吸道傳染病，是全球性的重大公共衛(wèi)生事件，新型冠狀病毒傳播速度快，致病力強[1]。為做好新冠肺炎疫情防控工作，上海市于2020年1月24日啟動重大突發(fā)公共衛(wèi)生事件一級響應。疫情防控期間，全社會動員共參與，廣大居民養(yǎng)成了外出佩戴口罩，居家勤通風、勤洗手等良好個人衛(wèi)生習慣。全民健康素養(yǎng)的提高不僅有效減少了新冠病毒在國內(nèi)的傳播和流行，也對猩紅熱等其他呼吸道傳染病的發(fā)生與流行產(chǎn)生重大影響。有研究發(fā)現(xiàn)新冠肺炎疫情期間法定傳染病的報告發(fā)病率顯著下降[2-3]。本研究收集了2011—2020年上海市松江區(qū)猩紅熱月報告病例數(shù)據(jù)，使用2011—2019年猩紅熱月發(fā)病數(shù)據(jù)，構建季節(jié)性自回歸移動平均（seasonal autoregressive integrated moving average，SARIMA）模型，采用最佳模型預測2020年松江區(qū)猩紅熱疫情發(fā)生情況，與新冠肺炎疫情期間猩紅熱的實際發(fā)病情況進行比較分析，了解新冠肺炎疫情對松江區(qū)猩紅熱的流行所產(chǎn)生的影響，為探索猩紅熱防控政策提供依據(jù)。

1 材料與方法

1.1 資料來源 在中國疾病預防控制信息系統(tǒng)的“傳染病監(jiān)測”模塊中，收集發(fā)病日期為2011年1月1日—2020年12月31日，且現(xiàn)住址為上海市松江區(qū)的猩紅熱病例數(shù)據(jù)。以月報告病例數(shù)為基本單位構建猩紅熱SARIMA乘積季節(jié)預測模型。

1.2 方法

1.2.1 SARIMA模型 SARIMA模型主要用于擬合具有平穩(wěn)性或可以被轉換為平穩(wěn)序列的時間序列[4]。SARIMA模型能綜合考慮數(shù)據(jù)序列的趨勢變化、周期變化及隨機干擾，其公式為SARIMA（p,d,q）（P,D,Q）s，其中p、d、q分別表示時間序列的自回歸階數(shù)、差分階數(shù)和移動平均階數(shù)，P、D、Q分別表示季節(jié)自回歸階數(shù)、季節(jié)差分階數(shù)、季節(jié)移動平均階數(shù)，s表示季節(jié)步長（周期）。

1.2.2 構建SARIMA模型步驟 ①時間數(shù)據(jù)序列平穩(wěn)性檢驗：采用單位根檢驗方法（augmented dickey fuller，ADF）對原始數(shù)據(jù)序列做平穩(wěn)性檢驗，若檢驗結果P＜0.05，則序列平穩(wěn)，若序列不平穩(wěn)，可采用對數(shù)轉換、差分、季節(jié)性差分實現(xiàn)時間序列的平穩(wěn)性。②參數(shù)選擇：根據(jù)原始數(shù)據(jù)序列或差分后序列的自相關系數(shù)圖（ACF）和偏相關系數(shù)圖（PACF）初步識別模型參數(shù)。R語言程序包計算赤池信息量準則（Akaike information criterion,AIC）、貝葉斯信息準則（Bayesian information criterion,BIC），選擇AIC、BIC值最小的模型作為最優(yōu)模型，用平均絕對百分比誤差（mean absolute percentage error,MAPE）評價預測精度，同時用最大似然法估計模型參數(shù)，并對參數(shù)進行統(tǒng)計學檢驗[4]。③模型檢驗與診斷：采用R語言程序包中Ljung-Box對模型擬合后的殘差序列進行白噪聲檢驗。檢驗結果P＞0.05，殘差為白噪聲序列，表示信息提取充分，模型建立有效。④預測與評價：使用R語言擬合松江區(qū)2011年1月—2019年6月的猩紅熱月報告病例數(shù)，構建最佳SARIMA模型，預測2019年7—12月猩紅熱月病例數(shù)，并與實際月報告病例數(shù)比較，利用相對誤差對預測效果進行評價。

1.2.3 影響評估 利用最佳SARIMA模型預測2020年1—12月猩紅熱發(fā)病數(shù)據(jù)，與2020年實際發(fā)生情況比較，評估新冠肺炎疫情對松江區(qū)猩紅熱流行的影響。

1.3 統(tǒng)計學分析 使用Excel 2010整理猩紅熱疫情數(shù)據(jù)和制圖，利用R語言4.0.4軟件中的stata、forecast、tseries等程序包進行SARIMA模型構建和統(tǒng)計學分析。檢驗水準α=0.05。

2 結果

2.1 猩紅熱流行概況 2011—2020年上海市松江區(qū)累計報告猩紅熱4 438例，年均報告444例，年均報告發(fā)病率23.27/10萬。其中2015年報告病例最多，為674例，年報告發(fā)病率35.64/10萬，分別是2020年報告病例數(shù)（94例）和年發(fā)病率（4.32/10萬）的7.17倍和8.25倍。2011—2020年月平均報告猩紅熱37例，其中2017年5月報告病例最多（140例，占當年報告病例的22.73%），報告病例最少的是2020年3月、4月和8月，均為零報告。2011—2019年松江區(qū)猩紅熱流行呈明顯季節(jié)性分布，每年有兩個發(fā)病高峰，4—6月為春季高峰，11月—次年1月為冬季高峰，兩個高峰累計報告病例數(shù)3 519例，占2011—2019年累計報告病例的79.29%。見圖1。

2.2 構建SARIMA模型

2.2.1 數(shù)據(jù)平穩(wěn)化處理 2011年1月—2019年6月松江區(qū)猩紅熱每年發(fā)病數(shù)基本持平，月報告病例數(shù)有明顯的季節(jié)性和周期性，見圖1。對原始數(shù)據(jù)進行ADF檢驗發(fā)現(xiàn)，原始時間序列數(shù)據(jù)不平穩(wěn)，需要對數(shù)據(jù)進行平穩(wěn)化處理?？紤]將原始數(shù)據(jù)進行自然對數(shù)變換，同時結合原始序列趨勢性和季節(jié)性，見圖2。在對數(shù)轉換的基礎上進行一階差分和一階12步季節(jié)差分，對差分處理后的數(shù)據(jù)序列進行影響因素分解，顯示松江區(qū)猩紅熱月報告病例數(shù)變化具有明確的趨勢性和季節(jié)性，見圖3。處理后的數(shù)據(jù)序列經(jīng)ADF檢驗（t=-8.96，P=0.01），具有統(tǒng)計學意義，可以認為差分后的數(shù)據(jù)序列平穩(wěn)。

圖1 2011—2020年上海市松江區(qū)猩紅熱月份報告病例情況

圖2 上海市松江區(qū)猩紅熱月報告病例對數(shù)時間序列圖

圖3 上海市松江區(qū)猩紅熱月報告病例數(shù)時間序列因素分解

2.2.2 模型識別與檢驗 通過對猩紅熱原始序列數(shù)據(jù)平穩(wěn)化過程得出，s、d、D值分別為12、1、1，建立初步模型SARIMA（p，1，q）（P，1，Q）12。繪制平穩(wěn)化后數(shù)據(jù)的自相關圖和偏相關圖，見圖4。觀察ACF圖和PACF圖，移動平均數(shù)參數(shù)q取值0、1、2，自回歸參數(shù)p取值1或2，由此建立待選模型。利用AIC、BIC、MAPE原則，結合模型參數(shù)的假設檢驗結果，綜合判定選出最優(yōu)模型SARIMA（1,1,0）（2,1,1）12，參數(shù)見表1。模型殘差的白噪聲Ljung-Box檢驗統(tǒng)計量的P值為0.96（P＞0.05），提示殘差序列為白噪聲，表明數(shù)據(jù)序列信息提示充分，模型建立有效。SARIMA（1,1,0）（2,1,1）12模型的殘差自相關和偏自相關圖，可見其值均在可信區(qū)間內(nèi)，也提示該擬合模型的殘差為白噪聲，見圖5。

圖4 上海市松江區(qū)猩紅熱平穩(wěn)化序列的自相關圖和偏自相關圖

圖5 SARIMA（1,1,0）（2,1,1）12模型的殘差自相關圖和偏自相關圖

表1 SARIMA（1,1,0）（2,1,1）12模型參數(shù)

2.2.3 模型擬合和預測 利用SARIMA（1,1,0）（2,1,1）12模型預測2019年7—12月松江區(qū)猩紅熱報告病例數(shù)。2019年猩紅熱實際發(fā)生值與預測值吻合度高，發(fā)病趨勢基本一致，8月、9月逐漸降低，11月、12月達到高峰。7—12月模型預測整體相對誤差為6%，其中8月預測偏差相對偏高，其他月份預測值均比較接近真實值，實際發(fā)生值均在95%可信區(qū)間，能夠真實地反映2019年7—12月松江區(qū)猩紅熱趨勢變化，見表2、圖6。

圖6 SARIMA（1,1,0）（2,1,1）12模型擬合與預測

表2 2019年7—12月上海市松江區(qū)猩紅熱月報告病例實際值與預測值比較

2.3 影響評估

2.3.1 2020年猩紅熱流行特征 2020年松江區(qū)累計報告病例94例，較2019年（515例）同期下降81.75%；年報告發(fā)病率4.56/10萬，較2019年（25.80/10萬）同期下降82.33%。月分布中2020年1月報告病例57例，占全年報告病例的60.64%；12月報告21例，占22.34%；3月、4月、8月均零報告。

2.3.2 疫情影響 依據(jù)最佳模型SARIMA（1,1,0）（2,1,1）12預測2020年1—12月猩紅熱發(fā)生情況，預測年報告病例591例。2020年猩紅熱病例實際發(fā)生情況較模型預期發(fā)病數(shù)下降84.09%，除1月猩紅熱實際發(fā)生值與最佳模型預測值接近，其他月份實際發(fā)生值遠低于預測值，也不在預測值的95%可信區(qū)間。見表3、圖6。

表3 2020年上海市松江區(qū)猩紅熱發(fā)生與預測情況比較

3 討論

傳染病的預測對傳染病的早期識別、早期預警起重要作用，有利于相關部門及時實施有針對性的健康宣教和行政干預措施。目前用于傳染病預測研究的數(shù)理模型有微分方程模型、時間序列模型和多因素模型等[5-6]。SARIMA模型因充分考慮時間序列的趨勢變化、周期以及季節(jié)變化，通過反復識別修正，獲取最佳模擬模型，操作相對簡單，實用性強、精度高，已廣泛應用于肺結核、猩紅熱、手足口病等傳染病短期預測[6-11]，并獲得了良好的早期控制效果。模型構建要求時間序列長度在30個數(shù)據(jù)以上[4]，本研究以2011年1月—2019年6月的猩紅熱月報告病例數(shù)據(jù)，累計102期數(shù)據(jù)擬合SARIMA模型，序列相對較長，可以擬合構建模型。因數(shù)據(jù)序列不平穩(wěn)對原始數(shù)據(jù)進行平穩(wěn)化處理，最后擬合出最佳模型SARIMA（1,1,0）（2,1,1）12，對模型殘差檢驗，驗證殘差為白噪聲序列，表示本次擬合模型成功。利用最佳模型預測2019年7—12月猩紅熱發(fā)生情況，預測值和實際值曲線走勢基本一致，與預測值相比，實際發(fā)生病例數(shù)均在預測值的95%可信區(qū)間內(nèi)波動，表明本次所構建的模型預測效果良好，所以可將該模型外推2020年1—12月猩紅熱發(fā)病數(shù)據(jù)。

2011—2019年上海市松江區(qū)猩紅熱年均發(fā)病率為23.27/10萬，高于蘇州、廣州、濟寧和南京等地近年來年均發(fā)病水平[12-15]，也高于上海市猩紅熱平均發(fā)病水平[16]。這可能是因為松江區(qū)經(jīng)濟快速發(fā)展，人口導入速度快，學校等集體單位數(shù)量多，學生基數(shù)大，傳染病傳播概率高，容易造成人群間局部流行與暴發(fā)。松江區(qū)猩紅熱發(fā)病曲線呈現(xiàn)2個高峰，分別為4—6月和11月—次年1月，低谷為每年的2月、8月、9月。猩紅熱發(fā)病季節(jié)雙峰特征與南京、蘇州等附近城市的研究結果一致[10,13]。2月、8月發(fā)病較少可能因為猩紅熱易感人群為幼托兒童和小學生，而2月與8月正值寒暑假期，學校放假，學生聚集機會減少，猩紅熱的暴露風險低。

受新冠肺炎疫情影響，全國傳染病發(fā)生情況較往年同期均有變化。浙江省新冠肺炎疫情應急響應期間，省內(nèi)其他法定傳染病同期下降50%左右[2]。本研究利用SARIMA模型對2020年猩紅熱發(fā)生情況進行預測，結果顯示2020年松江區(qū)猩紅熱實際發(fā)生水平遠低于模型預測值。分析2020年猩紅熱月報告病例數(shù)據(jù)發(fā)現(xiàn)，1月猩紅熱病例實際報告57例，接近預測值68例，與往年同期流行水平也接近。分析原因可能是1月份松江區(qū)學校、幼兒園等集體單位尚未進入寒假，學生、幼托兒童等易感人群仍在校就讀，集體活動多，猩紅熱的暴露風險高，與往年猩紅熱流行水平基本持平。1月24日上海市啟動重大突發(fā)公共衛(wèi)生事件一級響應，在2—4月新冠肺炎防控期間，上海市取消聚集活動，居民外出減少，學校開學延遲，降低了猩紅熱等其他呼吸道傳染病在學生和幼兒間的傳播，極大地降低了猩紅熱的發(fā)病水平。秋季開學后，學校等集體單位繼續(xù)落實日常消毒、勤通風等防控措施，9—11月猩紅熱持續(xù)呈低水平流行，11月、12月猩紅熱發(fā)病數(shù)逐漸增加，但仍遠低于去年同期水平。

科學掌握傳染病的發(fā)生規(guī)律，不僅有利于傳染病的早期防控，更便于及時調整制定適宜的防疫措施。本研究通過分析新冠肺炎疫情對松江區(qū)猩紅熱發(fā)病的影響，為評估新冠肺炎疫情對猩紅熱發(fā)生與流行提供了科學依據(jù)，今后可繼續(xù)采用構建數(shù)理模型的方法開展新冠肺炎疫情對其他傳染病影響的相關研究。