楊祖榮,劉 昆,譚 新,吳克堅(jiān)
(1空軍軍醫(yī)大學(xué)流行病學(xué)教研室,西安 710032;2西安市鄠邑區(qū)疾病預(yù)防控制中心;3空軍軍醫(yī)大學(xué)數(shù)理教研室;*通訊作者,E-mail:wukejian1983@163.com)
腎綜合征出血熱(hemorrhagic fever with renal syndrome,HFRS)又稱(chēng)流行性出血熱,是由漢坦病毒引起,以高熱、出血和腎功能損害為主要臨床特征[1],是嚴(yán)重危害人類(lèi)健康的自然疫源性疾病和重點(diǎn)傳染病之一[2,3]。我國(guó)作為世界上HFRS發(fā)病人數(shù)最多的國(guó)家,擁有占世界報(bào)道病例數(shù)90%以上的累計(jì)病人數(shù),內(nèi)陸31個(gè)省、自治區(qū)、直轄市均有HFRS的病例報(bào)告[4],特別是在我國(guó)華北、東北的部分地區(qū)HFRS已成為當(dāng)?shù)匚:ψ顕?yán)重的傳染病之一[5]。陜西省是我國(guó)HFRS疫情高發(fā)區(qū),西安地區(qū)的鄠邑、周至、長(zhǎng)安等區(qū)(縣)是疫情的主要分布地,其中鄠邑區(qū)是HFRS疫情國(guó)家級(jí)監(jiān)測(cè)點(diǎn)[6]。鄠邑區(qū)HFRS疫情的空間分布相對(duì)集中、人群分布以青壯年農(nóng)民為主。自1994年對(duì)鄠邑區(qū)HFRS發(fā)病率較高的平原鄉(xiāng)鎮(zhèn)啟動(dòng)了HFRS疫苗接種計(jì)劃,接種對(duì)象以15-60歲的人群為主[7]。為探究所采取措施是否具有針對(duì)性,本文基于求和自回歸移動(dòng)平均(autoregressive integrated moving average,ARIMA)模型對(duì)鄠邑區(qū)1971-01~2012-12的HFRS疫情數(shù)據(jù)進(jìn)行分析。
本研究1971-01~2012-12的HFRS疫情月度數(shù)據(jù)來(lái)自西安市鄠邑區(qū)疾病預(yù)防控制中心。
1.2.1 ARIMA模型 ARIMA模型來(lái)自20世紀(jì)70年代美國(guó)Wisconsin-Madison大學(xué)的Box和Jenkins所著的TimeSeriesAnalysis:ForecastingandControl一書(shū),主要用于分析時(shí)間序列,并對(duì)未來(lái)趨勢(shì)進(jìn)行預(yù)測(cè)。在ARIMA(p,d,q)模型中,AR代表自回歸模型,I表示差分,MA代表滑動(dòng)平均模型;p表示自回歸階數(shù),d表示對(duì)含有長(zhǎng)期趨勢(shì)、季節(jié)變動(dòng)、循環(huán)變動(dòng)的非平穩(wěn)時(shí)間序列進(jìn)行差分處理的次數(shù),q表示滑動(dòng)平均的階數(shù)[8]。ARIMA(p,d,q)模型的基本形式如下[9]:
(1)
1.2.2 模型建立 利用鄠邑區(qū)1971-01~2012-12每月HFRS發(fā)病數(shù)據(jù)建立數(shù)據(jù)庫(kù),通過(guò)R語(yǔ)言“tseries”和“forecast”軟件包進(jìn)行數(shù)據(jù)處理與分析,并建立ARIMA模型。根據(jù)時(shí)間序列圖或自相關(guān)系數(shù)判斷時(shí)間序列的平穩(wěn)性,對(duì)非平穩(wěn)時(shí)間序列進(jìn)行差分處理平穩(wěn)化,并利用自相關(guān)函數(shù)(autocorrelation function, ACF)圖和偏自相關(guān)函數(shù)(partial autocorrelation function, PACF)圖擬定模型p、q的階數(shù),對(duì)模型和參數(shù)的顯著性進(jìn)行檢驗(yàn),判斷殘差序列是否為白噪聲序列;若不是白噪聲序列則需再次調(diào)整差分階數(shù)以充分提取樣本信息。逐個(gè)檢驗(yàn)估計(jì)的參數(shù)是否顯著非零,若參數(shù)顯著不為零,則認(rèn)為構(gòu)建的模型合理,可進(jìn)行下一步的預(yù)測(cè)[11]。
1.2.3 模型比較 將1971-01~2012-12年每月HFRS發(fā)病數(shù)據(jù)代入對(duì)應(yīng)的ARIMA模型進(jìn)行擬合及預(yù)測(cè),比較所建立模型的擬合優(yōu)度及預(yù)測(cè)殘差。
1971-2012年,陜西省西安市鄠邑區(qū)共報(bào)告人間HFRS患者12 702例,HFRS年發(fā)病率在1971-1990年間呈波動(dòng)增長(zhǎng)趨勢(shì)。從1994年開(kāi)始年發(fā)病率逐步降低,但近幾年有回升趨勢(shì)。其中,1984年年發(fā)病率出現(xiàn)峰值,高達(dá)298.65/10萬(wàn)(1 439例),2005年出現(xiàn)最低值9.53/10萬(wàn)(55例),見(jiàn)圖1。疾病發(fā)生存在明顯的季節(jié)性和雙高峰特征,主高峰期出現(xiàn)在10-12月份(9 586例),占全年總病例數(shù)的75.47%,次高峰期出現(xiàn)在7月份,3月份發(fā)病率最低(見(jiàn)圖2)。
圖1 1971-2012年鄠邑區(qū)HFRS年發(fā)病率變化趨勢(shì)Figure 1 The annual incidence of HFRS in the Huyi District from 1971 to 2012
2.2.1 數(shù)據(jù)平穩(wěn)化 根據(jù)疫苗接種情況將鄠邑區(qū)HRFS月度發(fā)病數(shù)據(jù)分成1971-01~1993-12、1971-01~2012-12和1994-01~2012-12三組,分別進(jìn)行差分平穩(wěn)化處理。
2.2.2 模型的參數(shù)估計(jì) 通過(guò)分析數(shù)據(jù)預(yù)處理后的自相關(guān)函數(shù)和偏自相關(guān)函數(shù)圖,發(fā)現(xiàn)基于數(shù)據(jù)1971-01~2009-12建立的模型(記為模型Ⅰ)自相關(guān)函數(shù)在p=1,2,3時(shí)較大,偏自相關(guān)函數(shù)在q=1,2時(shí)較大,經(jīng)計(jì)算確定函數(shù)模型為ARIMA(1,0,1)(0,1,0)12,根據(jù)參數(shù)估計(jì)的結(jié)果構(gòu)建模型I為:
圖2 1971-2012年鄠邑區(qū)HFRS月發(fā)病人數(shù)分布Figure 2 The monthly incidence of HFRS in the Huyi District from 1971 to 2012
(2)
其中xt代表HFRS病例發(fā)病人數(shù)。
基于1971-01~1990-12發(fā)病數(shù)據(jù)建立的模型,記為模型Ⅱ,經(jīng)計(jì)算確定函數(shù)模型為ARIMA(1,0,1)(0,1,0)12,根據(jù)參數(shù)估計(jì)的結(jié)果構(gòu)建模型Ⅱ?yàn)?
(3)
類(lèi)似地,基于1994-01~2009-12發(fā)病數(shù)據(jù)建立的模型,記為模型Ⅲ,確定函數(shù)模型為ARIMA(2,0,1)(0,1,0)12,具體是:
(4)
2.2.3 模型的擬合優(yōu)度檢驗(yàn) 我們將時(shí)間序列原始數(shù)據(jù)和相應(yīng)的擬合數(shù)據(jù)的差值定義為模型的殘差序列,用殘差序列判斷模型的優(yōu)劣:假如殘差序列的自相關(guān)系數(shù)和偏自相關(guān)系數(shù)均在95%可信限以?xún)?nèi),則殘差為白噪聲序列。經(jīng)計(jì)算,模型Ⅰ、模型Ⅱ和模型Ⅲ擬合殘差是一個(gè)圍繞0波動(dòng)的平穩(wěn)序列,大部分殘差值在2倍標(biāo)準(zhǔn)差范圍內(nèi),且Box-Pierce統(tǒng)計(jì)量的P值均大于0.05,表明殘差為白噪聲序列,模型擬合程度較好。
2.3.1 利用所建的模型Ⅰ對(duì)1971-01~2009-12鄠邑區(qū)HFRS月發(fā)病率進(jìn)行模擬,并用實(shí)際發(fā)病率擬合,然后預(yù)測(cè)2010-01~2012-12的月發(fā)病率,利用預(yù)測(cè)值的殘差(即真實(shí)值-預(yù)測(cè)值)進(jìn)行模型檢驗(yàn)。由圖3可知,1971-01~2009-12月發(fā)病率的實(shí)際值與擬合值存在一定差異,但都在合理范圍內(nèi)波動(dòng)。2010-01~2012-12的預(yù)測(cè)月發(fā)病率與實(shí)際月發(fā)病率在短期內(nèi)基本相符,但隨著時(shí)間的增加,模型的預(yù)測(cè)誤差有增大的趨勢(shì),因此模型Ⅰ能比較合理地預(yù)測(cè)短期HFRS發(fā)病趨勢(shì)。
圖3 鄠邑區(qū)HFRS月發(fā)病率模型Ⅰ預(yù)測(cè)圖Figure 3 Prediction of the monthly incidence of HFRS in Huyi District by model Ⅰ
2.3.2 根據(jù)建立的模型Ⅱ?qū)?971-01~1990-12鄠邑區(qū)HFRS發(fā)病率進(jìn)行模擬,用實(shí)際發(fā)病率擬合,并預(yù)測(cè)1991-01~1993-12的發(fā)病情況。同樣地,用預(yù)測(cè)值的殘差進(jìn)行模型評(píng)價(jià)。1971-01~1991-12鄠邑區(qū)HFRS發(fā)病率觀(guān)測(cè)值與擬合值基本相符,波動(dòng)差異較小且都在合理范圍內(nèi);HFRS發(fā)病率預(yù)測(cè)值曲線(xiàn)與實(shí)際值曲線(xiàn)吻合程度高,表明模型Ⅱ可以很好地模擬預(yù)測(cè)未接種疫苗時(shí)鄠邑區(qū)HFRS的發(fā)病情況(見(jiàn)圖4)。
圖4 鄠邑區(qū)HFRS月發(fā)病率模型Ⅱ預(yù)測(cè)圖Figure 4 Prediction of the monthly incidence of HFRS in Huyi District by model Ⅱ
2.3.3 利用所建立的模型Ⅲ擬合1994-01~2009-12的HFRS月發(fā)病率,并對(duì)2010-01~2012-12的發(fā)病趨勢(shì)進(jìn)行預(yù)測(cè),見(jiàn)圖5。1994-01~2009-12鄠邑區(qū)HFRS月發(fā)病率觀(guān)測(cè)值與擬合值在擬合前期(1994-1999年)有一定波動(dòng),產(chǎn)生了較小的差異,但都在合理范圍內(nèi);而在擬合后期(2000-2009年)隨著模型逐漸趨于穩(wěn)定,觀(guān)測(cè)值與擬合值基本相符;HFRS月發(fā)病率預(yù)測(cè)值曲線(xiàn)與實(shí)際曲線(xiàn)基本吻合,表明模型Ⅲ可以更好地預(yù)測(cè)該時(shí)期鄠邑區(qū)HFRS的發(fā)病趨勢(shì)。
圖5 鄠邑區(qū)HFRS月發(fā)病率模型Ⅲ預(yù)測(cè)圖Figure 5 Prediction of the monthly incidence of HFRS in Huyi District by model Ⅲ
通過(guò)對(duì)模型Ⅰ和模型Ⅱ的綜合對(duì)比分析,發(fā)現(xiàn)對(duì)于同一個(gè)ARIMA模型(1,0,1)(0,1,0)12而言,模型Ⅱ的AIC值(1 922.72)明顯小于模型Ⅰ的AIC值(3 452.76),由圖3、圖4也能發(fā)現(xiàn)模型Ⅱ擬合程度優(yōu)于模型Ⅰ,其原因可能是模型Ⅰ是基于整體時(shí)間序列構(gòu)建,未考慮疫苗的接種情況;而模型Ⅱ是基于疫苗接種前的時(shí)間序列而構(gòu)建,未包含接種疫苗這一干預(yù)措施可能對(duì)時(shí)間序列造成的影響。由圖6可知,模型Ⅱ預(yù)測(cè)殘差的波動(dòng)范圍小,與HFRS的實(shí)際月發(fā)病率近似吻合,可見(jiàn)該模型可以對(duì)鄠邑區(qū)接種疫苗前的HFRS疫情做到很好的擬合及預(yù)測(cè)。
對(duì)模型Ⅰ和模型Ⅲ進(jìn)行對(duì)比分析,發(fā)現(xiàn)接種HFRS疫苗后建立的模型Ⅲ,其AIC值、擬合程度明顯優(yōu)于模型Ⅰ(見(jiàn)表1)。通過(guò)比較模型Ⅰ和模型Ⅲ的預(yù)測(cè)發(fā)病率,并利用方差來(lái)描述預(yù)測(cè)殘差的離散程度,發(fā)現(xiàn)模型Ⅲ的預(yù)測(cè)發(fā)病率更符合實(shí)際發(fā)病率(見(jiàn)圖7,8)。
圖6 模型Ⅱ HFRS月發(fā)病率預(yù)測(cè)殘差圖Figure 6 Prediction residuals of the monthly incidence of HFRS with model Ⅱ
表1模型Ⅰ與模型Ⅲ的擬合信息
Table1FittinginformationbetweenmodelⅠandmodelⅢ
模型AIC值殘差方差預(yù)測(cè)方差 模型Ⅰ3452.7690.3022.96 模型Ⅲ948.937.346.37
圖7 模型Ⅰ HFRS月發(fā)病率預(yù)測(cè)殘差圖Figure 7 Prediction residuals of the monthly incidence of HFRS with model Ⅰ
圖8 模型Ⅲ HFRS月發(fā)病率預(yù)測(cè)殘差圖Figure 8 Prediction residuals of the monthly incidence of HFRS with model Ⅲ
本文利用陜西省西安市鄠邑區(qū)疾病預(yù)防控制中心提供的HFRS疫情數(shù)據(jù),對(duì)1971-2012年鄠邑區(qū)HFRS新發(fā)病例的時(shí)間流行病學(xué)特征進(jìn)行了描述性分析和模型分析。通過(guò)構(gòu)建不同時(shí)間段的模型探討了疫苗接種干預(yù)效果,分析了HFRS疫苗接種情況對(duì)HFRS疫情的防控作用。
根據(jù)鄠邑區(qū)報(bào)告HFRS新發(fā)病例及其趨勢(shì)圖,可以發(fā)現(xiàn)自1994年開(kāi)始接種疫苗,鄠邑區(qū)HFRS發(fā)病率明顯降低,疫情得到控制。然而近幾年HFRS發(fā)病率有回升趨勢(shì),我們分析可能與外來(lái)務(wù)工人員的快速增長(zhǎng)有關(guān)。隨著經(jīng)濟(jì)發(fā)展,關(guān)中地區(qū)近年來(lái)外來(lái)務(wù)工人員數(shù)量龐大,該群體HFRS疫苗接種意識(shí)不強(qiáng)、疫情防護(hù)知識(shí)匱乏。這些因素提示我們應(yīng)在6-7月份和10-12月的發(fā)病高峰前期開(kāi)展針對(duì)性防病知識(shí)宣傳教育,加大宣講力度,并擴(kuò)大HFRS疫苗接種范圍尤其是對(duì)外來(lái)務(wù)工人員進(jìn)行疫苗接種。
通過(guò)模型擬合及預(yù)測(cè)比較,發(fā)現(xiàn)模型Ⅰ對(duì)接種疫苗后的HFRS發(fā)病數(shù)據(jù)擬合程度變差,只有重新建模型才能更好地?cái)M合、預(yù)測(cè)接種疫苗后的發(fā)病情況。當(dāng)模型Ⅰ的自回歸階數(shù)p由1階升到2階后,模型再次較好地?cái)M合疫苗接種后的HFRS發(fā)病數(shù)據(jù),表明接種疫苗使HFRS疫情周期性分布的循環(huán)周期發(fā)生改變[12],因此可以認(rèn)為疫苗的接種改變了HFRS疫情的發(fā)病趨勢(shì),對(duì)HFRS發(fā)病率的降低可能存在顯著影響。
由于ARIMA模型可以綜合考慮序列的趨勢(shì)變化、周期變化及隨機(jī)干擾,借助模型參數(shù)進(jìn)行量化表達(dá),能較精確地反映時(shí)間序列中所包含的動(dòng)態(tài)依存關(guān)系,因而具有明顯的優(yōu)勢(shì)和特色[13],同時(shí)也適用于各種復(fù)雜的時(shí)間序列,是一種實(shí)用性強(qiáng)、精確度高的短期預(yù)測(cè)方法[14]。本研究通過(guò)考慮時(shí)間序列與HFRS疫情的相關(guān)性、分組建模討論了接種HFRS疫苗對(duì)疫情的影響,但建立的模型只適合進(jìn)行短期較為準(zhǔn)確的預(yù)測(cè),并且疫情的傳播是多種因素共同作用的結(jié)果,各因素之間存在著錯(cuò)綜復(fù)雜的聯(lián)系,很難運(yùn)用結(jié)構(gòu)式的因果模型進(jìn)行預(yù)測(cè)。因此下一步,我們將綜合考慮鄠邑區(qū)HFRS疫情發(fā)生的季節(jié)性、接種疫苗、溫度、濕度、降水、地域特征和人群特征等因素,不斷積累、添加新的監(jiān)測(cè)數(shù)據(jù),以便能獲得更加準(zhǔn)確的ARIMA模型,為鄠邑區(qū)HFSR疫情的防控提供理論指導(dǎo),同時(shí)本研究也清晰地表明接種HFRS疫苗對(duì)HFRS疫情的防控存在明顯的積極作用,但進(jìn)一步加強(qiáng)HFRS防治工作,開(kāi)展疫情監(jiān)測(cè)和深入研究工作仍是HFRS防控工作的主要任務(wù)。