王燕芬，王旭峰，趙繼軍

(青島大學(xué)復(fù)雜性科學(xué)研究所，山東 青島，266071)

0 引言

圖1 2009年～2016年手足口病發(fā)病數(shù)時(shí)間序列圖Fig.1 Time series of the incidence of HFMD from 2009 to 2016

1 資料與方法

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

本文使用的手足口病監(jiān)測(cè)數(shù)據(jù)為2009年至2016年月報(bào)告發(fā)病數(shù)，來(lái)源于公共衛(wèi)生科學(xué)數(shù)據(jù)中心(www.phsciencedata.cn)，全國(guó)各省年人口總數(shù)及年出生率來(lái)源于中國(guó)統(tǒng)計(jì)局(www.stats.gov.cn)。氣候數(shù)據(jù)包括月平均氣溫、降雨量和相對(duì)濕度，來(lái)源于中國(guó)氣象信息中心(data.cma.cn)。

寒暑假時(shí)間的選取參照各省教育廳公布的歷年中小學(xué)開(kāi)學(xué)及放假時(shí)間通知，1月20日至2月20日為寒假時(shí)間，7月、8月為暑假。春運(yùn)共計(jì)40天，一般為春節(jié)前15天和節(jié)后25天，由于每年時(shí)間不同，根據(jù)最早開(kāi)始時(shí)間和最晚結(jié)束時(shí)間，選擇1月中旬至3月中旬作為春運(yùn)期間。

1.2 方法

本文首先建立時(shí)間序列SIR (Time Series Susceptible Infected Recovered，TSIR)模型估算各省手足口病月傳染率，模型中的參數(shù)由馬爾科夫蒙特卡洛方法(Markov chain Monte Carlo，MCMC)進(jìn)行估計(jì)，然后利用線性回歸模型分析月平均氣溫、降雨量、相對(duì)濕度、開(kāi)學(xué)放假和春運(yùn)對(duì)各省手足口病傳染率的影響。

(1)

其中，B(t)為出生人數(shù)，Yr(t)為報(bào)告發(fā)病人數(shù)，ρ為報(bào)告率，β(t)為傳染率，本文假設(shè)每年同一個(gè)月的傳染率相同。由于經(jīng)疾病監(jiān)測(cè)系統(tǒng)上報(bào)的發(fā)病人數(shù)僅占實(shí)際發(fā)病數(shù)的一定比例，在計(jì)算下一時(shí)刻易感者人數(shù)時(shí)候，需要除去實(shí)際感染者人數(shù)，所以需要計(jì)算相應(yīng)報(bào)告率來(lái)糾正發(fā)病人數(shù)。本文使用平滑樣條法對(duì)各地報(bào)告率進(jìn)行估算，模型中自由度選取經(jīng)驗(yàn)值2.5[23,27]。

實(shí)際生活中，疾病的發(fā)病人數(shù)通常服從泊松分布或負(fù)二項(xiàng)分布[19,25]，本文假設(shè)TSIR模型中感染者人數(shù)服從泊松分布：

Y(t)～Possion(λ)

(2)

其中，λ為感染力，且λ的期望值為：

E(λ)=β(t)Yα(t)X(t)

(3)

其中，α(0<α≤1)為人群混合程度，當(dāng)α<1，說(shuō)明人群是不完全混合的，即部分人之間的接觸高于他們與其他人群的接觸；當(dāng)α=1，表示人群是均勻混合。因?yàn)門SIR模型中步長(zhǎng)的選擇(兩周或是一個(gè)月)不影響最終結(jié)果[13]，而我們收集的數(shù)據(jù)為月發(fā)病數(shù)，所以采用時(shí)間步長(zhǎng)為1個(gè)月。利用TSIR模型，本文估計(jì)了各省、自治區(qū)月手足口病傳染率β(t)、人群混合程度α、易感者人數(shù)X(t)及報(bào)告率ρ共15個(gè)參數(shù)。

TSIR模型中的參數(shù)由馬爾科夫蒙特卡洛方法(Markov chain Monte Carlo，MCMC)進(jìn)行估計(jì)[28]，馬爾科夫鏈用Gibbs采樣法產(chǎn)生，各地區(qū)初始易感者人數(shù)X(0)、人群混合程度α和傳染率β(t)的先驗(yàn)參數(shù)都服從均勻分布，初始易感者人數(shù)X(0)取值范圍為總?cè)丝诘?%～40%，人群混合程度α取值范圍為[0.5, 0.999]。為了保證各省相應(yīng)的基本再生數(shù)在1至100范圍內(nèi)，其傳染率初值的取值范圍不同。在初始迭代1 000次后，馬爾科夫鏈迭代1 000萬(wàn)次，再對(duì)最后5 000個(gè)有效參數(shù)進(jìn)行采樣，為了避免采樣結(jié)果自相關(guān)，采樣間隔取50。由于手足口病疫苗于2016年才投入使用，且對(duì)發(fā)病數(shù)影響并不明顯，所以TSIR模型中可以不考慮免疫措施的影響。

最后，在獲得各省傳染率的基礎(chǔ)上，估算各省手足口病基本再生數(shù)，利用線性回歸模型分析月平均氣溫、降雨量、相對(duì)濕度、開(kāi)學(xué)放假和春運(yùn)期間對(duì)各省傳染率的影響。氣候因素為月數(shù)據(jù)變量，開(kāi)學(xué)放假(“學(xué)期”、“寒假”或“暑假”)和春運(yùn)(“春運(yùn)期間”或“平常時(shí)期”)均為分類變量。由于估計(jì)的傳染率為月數(shù)值，所以本文選用2月作為寒假時(shí)間，7月、8月作為暑假時(shí)間。每年春運(yùn)時(shí)間并不固定，考慮到人群接觸后疾病的傳播存在一定的時(shí)間滯后，本文選取有兩個(gè)星期滯后的2月和3月作為線性回歸模型中的春運(yùn)期間。我們首先只選取氣候因素進(jìn)行線性回歸，為了避免各因素之間相關(guān)性導(dǎo)致的模型中多重共線性問(wèn)題，對(duì)變量進(jìn)行篩選，將方差膨脹因子(方差膨脹因子越大，共線性越嚴(yán)重)大于5的因素依次去掉后，再找出對(duì)手足口病具有顯著影響的氣候因素。其次，應(yīng)用包括具有顯著影響的氣候因素、開(kāi)學(xué)放假和春運(yùn)期間在內(nèi)的線性回歸模型，分析各因素對(duì)手足口病傳染率的影響。

2 結(jié)論

2.1 MCMC-TSIR模型采樣結(jié)果

利用馬爾科夫蒙特卡洛方法估算TSIR模型中的參數(shù)，對(duì)最后5 000個(gè)有效參數(shù)進(jìn)行采樣(圖2)，傳染率β(t)、人群混合度α和易感者人數(shù)X(t)收斂在一定區(qū)間內(nèi)。

圖2 以湖南省為例，MCMC-TSIR模型最后5 000次迭代采樣結(jié)果Fig.2 Taking Hunan Province as an example, the results of the last 5 000 iterations of the MCMC-TSIR model

2.2 MCMC-TSIR模型估算結(jié)果及模型驗(yàn)證

中國(guó)各省手足口病傳染率的峰值有明顯的季節(jié)性(見(jiàn)圖3)，且在2月都有較大的增幅，這可能與春運(yùn)期間大規(guī)模人口流動(dòng)引起的人群接觸率快速增加有關(guān)，因?yàn)闆](méi)有發(fā)現(xiàn)氣候因素在2月有所變化。在暑假期間，多數(shù)省的傳染率處于較低水平。為了更好進(jìn)行描述，根據(jù)傳染率峰值發(fā)生時(shí)間，可以將全國(guó)各省歸類到為4個(gè)區(qū)域(見(jiàn)表1)：

圖3 全國(guó)各省手足口病傳染率峰值分布及各區(qū)域傳染率偏移值Fig.3 Peak distribution of HFMD transmission rate and the deviance of β(t) in each region

地區(qū)包含省份東南地區(qū)海南、廣東、廣西、浙江、江蘇、上海、安徽、云南、四川、重慶、貴州、山西、陜西、北京、天津、山東、江西、湖北、湖南、福建、河北、河南北部地區(qū)內(nèi)蒙古、吉林、遼寧和沈陽(yáng)西北地區(qū)新疆、青海、甘肅和寧夏西藏西藏

屬于東南地區(qū)的各省(除山東省和北京市外)的手足口病傳染率峰值集中于2月～3月；屬于西北地區(qū)的省份傳染率峰值在4月；東南地區(qū)和西北地區(qū)各省的手足口病傳染率在秋季還有一個(gè)小高峰。北方地區(qū)傳染率峰值在5月；西藏傳染率峰值則在8月，這可能是因?yàn)槲鞑氐貐^(qū)手足口病的平均感染年齡在5歲左右，即學(xué)齡兒童更易感染疾病，因此學(xué)校開(kāi)學(xué)、放假導(dǎo)致的學(xué)生間接觸率的變化能夠影響其傳染率的季節(jié)性。

根據(jù)模型估計(jì)的結(jié)果(見(jiàn)圖4)：全國(guó)多數(shù)省份手足口病基本再生數(shù)在15至35之間，只有廣東省的手足口病基本再生數(shù)高于40，貴州省基本再生數(shù)則小于10；手足口病報(bào)告率在各地區(qū)也有較大差異，從整體上看，處于東南地區(qū)的省份手足口病報(bào)告率高于北部、西北地區(qū)，其中山西省、青海省、新疆和西藏自治區(qū)的報(bào)告率不足5%，廣西省和海南省的報(bào)告率卻高達(dá)35%，其余省份的報(bào)告率多處于10%至25%之間；各省手足口病易感者比例則集中在10%至35%之間。

圖4 MCMC-TSIR模型估算各省參數(shù)結(jié)果Fig.4 Results of estimated parameters in different province of MCMC-TSIR model

利用MCMC-TSIR模型所得參數(shù)估計(jì)發(fā)病數(shù)，對(duì)實(shí)際報(bào)告發(fā)病數(shù)和預(yù)測(cè)發(fā)病數(shù)取對(duì)數(shù)后進(jìn)行分析。模型預(yù)測(cè)的發(fā)病數(shù)和實(shí)際發(fā)病數(shù)之間的關(guān)系為圖5。

2.3 手足口病傳染率季節(jié)性影響因素

手足口病傳染率季節(jié)性在不同地區(qū)受不同因素影響(見(jiàn)表2)：東南地區(qū)各省份傳染率季節(jié)性僅與春運(yùn)有關(guān)；北部地區(qū)手足足病傳染率季節(jié)性則受平均氣溫(解釋貢獻(xiàn)率20.10%)、暑假(解釋貢獻(xiàn)率49.42%)和春運(yùn)(解釋貢獻(xiàn)率8.11%)的影響；西北地區(qū)傳染率季節(jié)性主要受到相對(duì)濕度(解釋貢獻(xiàn)率31.27%)和春運(yùn)(解釋貢獻(xiàn)率9.80%)的影響；西藏省手足口病傳染率季節(jié)性既不受氣候因素影響，也不受春運(yùn)的影響。

圖5 模型估計(jì)發(fā)病數(shù)及報(bào)告發(fā)病數(shù)關(guān)系Fig.5 The relationship between the estimated number of cases and the number of reported cases

影響因素估計(jì)值平均溫度相對(duì)濕度暑假春運(yùn)期間R2東南地區(qū)回歸系數(shù)---1.204e-060.642 2P值-< 2.2e-16北部地區(qū)回歸系數(shù) 7.565e-08--1.944e-069.748e-070.420 3P值4.93e-06-0.000 1630.075 522西北地區(qū)回歸系數(shù)--1.141e-07 -1.710e-060.414 9P值-0.010 85-0.068 26西藏回歸系數(shù)-----P值----

各地區(qū)手足口病均不受降雨量影響，因此表中去除了降雨量這一項(xiàng)。地區(qū)內(nèi)各省存在差異性，對(duì)各省的手足口病傳染率季節(jié)性影響因素進(jìn)行分析發(fā)現(xiàn)東南地區(qū)各省傳染率季節(jié)性都與春運(yùn)有顯著相關(guān)，其中廣東省傳染率季節(jié)性還與相對(duì)濕度有關(guān)，北京、天津傳染率季節(jié)性還與相對(duì)濕度和溫度相關(guān)；北部地區(qū)各省的傳染率季節(jié)性都與相對(duì)濕度有關(guān)，吉林、遼寧和黑龍江省傳染率季節(jié)性還受平均溫度影響；屬于西北地區(qū)的寧夏、甘肅手足口病傳染率季節(jié)性與相對(duì)濕度、溫度有關(guān)，青海省傳染率僅與春運(yùn)有關(guān)，新疆自治區(qū)傳染率季節(jié)性則與相對(duì)濕度和降雨量有關(guān)。

3 討論

手足口病傳染率在中國(guó)各省都有明顯的季節(jié)性，且在2月都有大幅度增加，這與春運(yùn)期間，整體人群接觸率變化有關(guān)。根據(jù)傳染率峰值的時(shí)間分布，可以將全國(guó)各省歸類為四個(gè)區(qū)域：東南地區(qū)、西北地區(qū)、北部地區(qū)及西藏。東南地區(qū)的多數(shù)省份傳染率高峰集中在2月和3月，西北地區(qū)，包括新疆、青海、甘肅和寧夏，手足口病傳染率高峰發(fā)生在4月，北方地區(qū)，包括內(nèi)蒙古、遼寧、吉林和沈陽(yáng)的傳染率在5月達(dá)到最大值，而西藏地區(qū)的傳染率在8月達(dá)到峰值。

根據(jù)MCMC-TSIR模型，全國(guó)多數(shù)省份手足口病基本再生數(shù)在15至35之間，這可能與人口密度、易感者比例等因素有關(guān)。根據(jù)估算的基本再生數(shù)，手足口病的免疫覆蓋率應(yīng)高于93.3%。從整體上看，多數(shù)省份的報(bào)告率在10%至25%之間，處于東南地區(qū)的省份手足口病報(bào)告率高于北部地區(qū)和西北地區(qū)，而西藏自治區(qū)報(bào)告率僅2.5%，這可能與當(dāng)?shù)蒯t(yī)療水平、公共衛(wèi)生健康監(jiān)測(cè)系統(tǒng)的完善程度有關(guān)；各省手足口病易感者比例則集中在10%至35%之間。

手足口病傳染率的季節(jié)性影響因素主要分為兩類：社會(huì)因素及氣候因素。東南地區(qū)的大多數(shù)省份的手足口病傳染率季節(jié)性僅受春運(yùn)期間的顯著影響。東南地區(qū)多數(shù)省份經(jīng)濟(jì)發(fā)達(dá)，流動(dòng)人口數(shù)目高于其他地區(qū)，春運(yùn)期間，人口遷移量更是急劇增加，使得各年齡組的接觸率遠(yuǎn)高于平常時(shí)間，封閉、擁擠的旅途環(huán)境使得疾病的傳播可能性大大增加。北方地區(qū)的傳染率季節(jié)性則與平均溫度和暑假相關(guān)，西北地區(qū)的傳染率季節(jié)性受到相對(duì)濕度和春運(yùn)的影響；這可能是因?yàn)槲鞅薄⒈辈康貐^(qū)的月氣候變化較為顯著，氣候變化會(huì)影響病毒在環(huán)境中的生存率，還會(huì)在一定程度上影響人的行為方式，間接影響人群接觸率變化。西藏地區(qū)傳染率不受氣候、學(xué)期和春運(yùn)的影響。在寒暑假期間，各省的手足口病傳染率也有所增加，尤其西藏地區(qū)，其手足口病傳染率峰值在8月。這可能是因?yàn)橹袊?guó)各省手足口病平均感染年齡不同，西藏地區(qū)手足口病的平均感染年齡在5歲左右，學(xué)齡兒童更易感染疾病[4]，因此學(xué)校開(kāi)學(xué)、放假導(dǎo)致的兒童間接觸率的變化對(duì)其傳染率的季節(jié)性的影響較其他省更加顯著。

為了更好了解手足口病傳染率的季節(jié)性及其影響因素，需要更進(jìn)一步研究，在未來(lái)研究中，我們將具體分析不同年齡組手足口病傳染率季節(jié)性及其影響因素，并將影響因素代入模型中，通過(guò)曲線下的面積(Area Under Curve，AUC)來(lái)進(jìn)行評(píng)估。