白 迪,周武略,趙繼軍
(青島大學復(fù)雜性科學研究所,山東 青島 266071)
瘧疾(malaria)是一種由于按蚊叮咬或輸入帶瘧原蟲者的血液而感染瘧原蟲所引起的蟲媒傳染病[1],寄生于人體的瘧原蟲共有4種,即間日瘧原蟲,三日瘧原蟲,惡性瘧原蟲和卵形瘧原蟲。其中間日瘧原蟲傳播范圍最廣,惡性瘧原蟲致死率最高[2]。自20世紀50年代開始中國對瘧疾逐步開展了有計劃的全面防治,到20世紀末大部分地區(qū)的瘧疾流行已得到控制[3]。雖然一些地區(qū)感染瘧疾的病例維持在了較低的水平,甚至部分省份自2019年開始陸續(xù)通過了省級消除瘧疾評估,但是這些地區(qū)的瘧疾傳播媒介按蚊仍普遍存在,瘧疾傳播的條件尚未消除[4],仍有再次流行的風險。江蘇省泗洪縣在20世紀90年代瘧疾發(fā)病數(shù)降低至歷史最低水平后,于2000年再次暴發(fā)了疫情,這次暴發(fā)與從鄰近省份輸入的瘧疾病例有關(guān)[5]。由于瘧疾傳播媒介的存在,輸入病例易成為傳染源,引起繼發(fā)傳播。此外對非洲高地瘧疾的研究中發(fā)現(xiàn),瘧疾的傳染率受降水量周期性變化的影響[6]。當降水量的年周期降低時,瘧疾的傳播程度逐漸增大。當滿足一些特定條件的情況下,例如足夠多的輸入病例、未來某段時間氣象因素的變化都有可能導(dǎo)致瘧疾的重新暴發(fā)。
對瘧疾在中國流行的研究中發(fā)現(xiàn),瘧疾的報告發(fā)病數(shù)與溫度、降水量等因素有關(guān)。在氣溫接近26℃與27℃時瘧疾的發(fā)病數(shù)相對較高[7-8];瘧疾發(fā)病數(shù)的時間序列與滯后0~3個月的降水量呈現(xiàn)強正相關(guān),與滯后4~9個月的降水量呈負相關(guān)[9]。這些研究是對瘧疾發(fā)病數(shù)與氣候因素之間的相關(guān)性研究。中國本地出現(xiàn)的瘧疾類型主要是間日瘧,其傳播媒介為中華按蚊[10]。相對于非洲高地關(guān)于岡比亞按蚊傳播惡性瘧傳染率的研究,中國通過中華按蚊傳播間日瘧傳染率的研究極其有限。傳染率是疾病傳播特性中的一個重要參數(shù),決定了疾病流行的程度(即發(fā)病數(shù)的高低)和流行的周期。本文以安徽省2004~2009年流行的間日瘧為例,分析間日瘧的傳播動態(tài),包括傳染率和傳染率季節(jié)性。
安徽省在地理位置上處于中國東南部,為亞熱帶和暖溫帶過渡氣候,四季分明,適宜間日瘧傳播媒介中華按蚊的繁殖,曾多次暴發(fā)間日瘧疫情[11]。2005~2015年間,安徽省間日瘧累計發(fā)病數(shù)達59 752例,居全國首位,占全國累計間日瘧發(fā)病數(shù)的51.39%[12]。本文將在研究安徽省間日瘧傳播動態(tài)的基礎(chǔ)上,分析氣候因素包括溫度、降水量對通過中華按蚊傳播的間日瘧傳染率的影響。本文的研究結(jié)果為進一步預(yù)測間日瘧在一個地區(qū)消除后重新暴發(fā)的可能性以及影響暴發(fā)的條件等方面的研究提供指導(dǎo)。
本文從中國疾病預(yù)防控制中心(www.phsciencedata.cn)獲得安徽省2004~2009年各月間日瘧的發(fā)病數(shù)。從安徽省統(tǒng)計局(http://tjj.ah.gov.cn/tjjweb/web/index.jsp)獲得2004~2009年安徽省人口出生率、死亡率。從中國氣象數(shù)據(jù)網(wǎng)(http://data.cma.cn/user/toLogin.html)獲得安徽省各地級市氣象站點2004~2009年溫度、降水量的數(shù)據(jù)。
建立Susceptible Infected Recovered-Susceptible Infected(SIR-SI)模型估算間日瘧傳染率矩陣并分析傳染率的季節(jié)性,進一步通過多元線性回歸模型分析影響間日瘧傳染率的因素。
圖1 SIR-SI模型示意圖
為研究安徽省間日瘧的傳播特性,我們建立病毒通過蚊子叮咬在人與蚊子之間傳播間日瘧的SIR-SI倉室模型。倉室模型是由微分方程構(gòu)成的確定性建模方法。對于蚊媒傳播的疾病來說,倉室模型將人群和蚊群分為幾種不同類型,每種類型作為一個“倉室”存在。本文中的模型包括兩個部分,SIR部分描述人類易感者(Sh)、感染者(Ih)、恢復(fù)者(Rh)各倉室中人群數(shù)量的動態(tài),SI部分描述蚊子易感(Sm)、感染(Im)各倉室中蚊子數(shù)量的動態(tài)。由于蚊子的生命周期較短,在這里我們只考慮人類從感染到恢復(fù),因此蚊子沒有恢復(fù)者倉室[13-14]。模型的傳播機制如圖1所示。
(1)
其中Xh,Yh,Zh和Nh分別表示t時刻人類易感者,感染者,恢復(fù)者和所有人類的數(shù)量,Xm,Ym和Nm分別表示t時刻易感蚊子,感染蚊子和所有蚊子的數(shù)量。r是蚊子叮咬率,Thm(≤1)是受感染的蚊子叮咬易感者后易感者被感染的概率,Tmh(≤1)是易感蚊子叮咬感染者后蚊子被感染的概率。γh是人的恢復(fù)率;νh是人的總出生率,νm是蚊子的總出生率;μh是人的死亡率,μm是蚊子的死亡率。參數(shù)值的設(shè)置見表1。
表1 參數(shù)取值
對于蟲媒傳播的傳染病(如瘧疾、登革熱),傳染率季節(jié)性起到了非常重要的作用。在這里我們考慮將模型中的季節(jié)性表示為蚊子叮咬率r,r是關(guān)于時間t的函數(shù)[13]:
r(t)=r0(1-r1cos(ωt+φ))
(2)
其中,r0是蚊子的平均叮咬率,ω是驅(qū)動頻率,φ是相位,r1是蚊子季節(jié)性叮咬率。
模型(1)中的傳染率矩陣為
(3)
通過最小化根據(jù)模型所得的數(shù)據(jù)與實際數(shù)據(jù)差值的平方和,即最小二乘法來得到最佳參數(shù)r0和r1的值,具體公式為
(4)
其中,DIi(t)為t時刻安徽省間日瘧的實際數(shù)據(jù),PIi(t)為t時刻根據(jù)模型仿真得到的數(shù)據(jù),T為仿真時長。利用R語言進行數(shù)據(jù)整理并建立模型,獲取使Sr最小的參數(shù)r0和r1,即是模型所估計的最佳參數(shù)值。
為了更準確地了解間日瘧傳播的動態(tài)范圍,我們計算基本再生數(shù)R0?;驹偕鷶?shù)R0是一項重要的流行病學指標,預(yù)示著疾病傳播的潛在危險程度和控制的難易程度[19]。對于間日瘧來說,R0是指一只被感染的蚊子會產(chǎn)生二次感染的數(shù)量?;驹偕鷶?shù)R0=1就是間日瘧流行與不流行的臨界水平,當R0>1時,即1個病例至少傳播產(chǎn)生1個繼發(fā)病例時,間日瘧才可以繼續(xù)傳播下去,否則傳播將趨向終止[20]。根據(jù)式(1)的均衡狀態(tài),可以推導(dǎo)出間日瘧的基本再生數(shù)為
(5)
本文選取安徽省所有地級市氣象站點的月氣溫和降水量兩個因素的數(shù)值并分別求平均,與模型估計出的蚊子叮咬率進行回歸分析,由于溫度和降水量具有不同的量綱單位,為了消除指標之間的量綱影響,我們對兩者進行數(shù)據(jù)標準化處理,保證溫度和降水量處于同一數(shù)量級,適合進行綜合對比評價。多元線性回歸函數(shù)為
r(t)=aTemperature(t)+bRain(t)+c+ε
(6)
其中,蚊子叮咬率r(t)為響應(yīng)變量,a、b為回歸系數(shù),c為回歸常數(shù),溫度Temperature(t)和降水量Rain(t)為回歸因子,ε是均值為零且方差不可觀測的隨機變量,稱為誤差項。
應(yīng)用最小二乘法對SIR-SI模型的參數(shù)進行估計,得出安徽省間日瘧的蚊子平均叮咬率r0為4.27e-06,季節(jié)性叮咬率r1為42.4%,相位φ為0.000 3。將蚊子叮咬r(t)帶入到傳染率矩陣中,即:
(7)
觀察到傳染率在每年的6、7月最高,12、1月最低。
對比模型仿真結(jié)果和實際數(shù)據(jù),我們發(fā)現(xiàn)模型結(jié)果與安徽省2004~2009年間日瘧發(fā)病趨勢基本符合,且具有每年7~10月發(fā)病季節(jié)高峰的特點(見圖2),但模型仿真結(jié)果與實際數(shù)據(jù)仍存在一定偏差。
將r0代入到式(5)中,計算出基本再生數(shù)R0的值為1.22。R0>1說明間日瘧在一定程度上流行,并觀察到發(fā)病數(shù)具有年周期性的特征。
多元線性回歸模型的結(jié)果表明溫度和降水量對間日瘧的蚊子叮咬率均有顯著影響,決定系數(shù)R2=0.73。其中溫度(℃)的回歸系數(shù)為2.609e-06(p=7.23e-11),降水量(mm)的回歸系數(shù)為9.951e-07(p=0.004),溫度的變化對于蚊子叮咬率的影響更顯著。此外由于溫度和降水量的方差膨脹因子VIF小于5,因此不存在顯著的共線性。
模型仿真結(jié)果與實際發(fā)病數(shù)的比較,-o-為實際值,——黑色為預(yù)測值。圖2 SIR-SI模型仿真結(jié)果
圖3 安徽省間日瘧預(yù)測感染者數(shù)(a)與預(yù)測易感者數(shù)(b)
圖4 安徽省間日瘧年報告發(fā)病數(shù)
假設(shè)2009年以后對間日瘧沒有進一步的控制措施,根據(jù)估計出的傳染率值應(yīng)用SIR-SI模型對發(fā)病數(shù)趨勢進行預(yù)測,觀察到在2018年后間日瘧感染者人數(shù)維持在較低的狀態(tài),在2034年后開始增加并于2048年達到峰值后又逐漸下降。預(yù)測的感染者人數(shù)不僅呈現(xiàn)出每年一個周期的特點,還存在多年一次大周期的特點(見圖3a)。在年感染者人數(shù)處于低谷期間,易感者人數(shù)則不斷累計,直至達到一定數(shù)值,感染者人數(shù)的均值和變化幅度又開始逐漸增加(見圖3b)。20世紀末的安徽省間日瘧年發(fā)病數(shù)呈現(xiàn)過再次上升的特點,在1993~1999年發(fā)病數(shù)較低,自2000年開始,間日瘧年發(fā)病數(shù)又不斷增加(見圖4)。
本文以安徽省間日瘧為例,估算蚊子叮咬率、間日瘧的傳染率以及它們的季節(jié)性,并進一步分析傳染率發(fā)生季節(jié)性變化的影響因素。采用人類與蚊子相耦合的SIR-SI模型對間日瘧的傳播動態(tài)建模,估算隨時間變化的蚊子叮咬率以及傳染率。并進一步建立多元線性回歸模型探究溫度和降水量對蚊子叮咬率和間日瘧傳染率的影響。根據(jù)對安徽省間日瘧傳染率的分析結(jié)果,我們發(fā)現(xiàn)間日瘧的傳染率在一年中存在大幅度變化(變化幅度為42.4%),傳染率在在6、7月最高,在12、1月最低。傳染病的傳染率季節(jié)性形成了發(fā)病數(shù)的周期性[13],因而間日瘧發(fā)病數(shù)的年周期性波動是由較高的傳染率的季節(jié)性導(dǎo)致的。一些兒童傳染病的傳染率峰值早于發(fā)病數(shù)峰值若干星期或若干月[21]。本文研究發(fā)現(xiàn),通過中華按蚊傳播的間日瘧的傳染率峰值也有類似的現(xiàn)象。間日瘧的蚊子叮咬率和傳染率又受溫度和降水量變化的影響,其中溫度對間日瘧傳染率的影響更大。有研究表明,溫度對中華按蚊的影響包括多個方面,例如繁殖、生存、活動行為等[22-23]。此外溫度對中華按蚊的吸血行為也存在顯著影響,在中華按蚊傳播間日瘧的溫度范圍內(nèi),隨著溫度的升高,雌蚊的叮咬頻率逐漸增加[24]。本文通過具有年季節(jié)性變化的叮咬率這個參數(shù)考慮了中華按蚊個體的叮咬頻率和中華按蚊數(shù)量變化對整體叮咬率的影響。
本文的研究結(jié)果顯示,間日瘧易感者人數(shù)存在多年周期性的變化趨勢,從而形成間日瘧感染者數(shù)多年的大周期。雖然間日瘧的感染者數(shù)能夠在多年維持較低的水平,但在易感者人數(shù)不斷累積并達到一定數(shù)量時,感染者數(shù)又會逐年增加,這能夠解釋90年代末安徽省間日瘧在保持若干年較低發(fā)病數(shù)后又逐年大幅增加的現(xiàn)象。值得注意的是90年代間日瘧開始保持較低的年發(fā)病數(shù)(<2 000),至2000年發(fā)病數(shù)開始快速增加,保持較低發(fā)病數(shù)的時間長度為6年。本文的仿真結(jié)果顯示,2016年后間日瘧被感染人數(shù)也呈現(xiàn)類似的結(jié)果,但不同點是仿真結(jié)果在2016年后有更長的時間(約20年)處于較低發(fā)病數(shù)的狀態(tài)??紤]可能的原因是2016年后對間日瘧的控制措施優(yōu)于90年代,降低了傳染率,這使得間日瘧的傳播速度相較于90年代更慢,因此2016年后間日瘧保持較低年發(fā)病數(shù)的時間更長。
本文中模型估計出的發(fā)病數(shù)與實際報告數(shù)據(jù)有一定偏差,產(chǎn)生的原因可能來自以下三個方面:1)最小二乘法對未知參數(shù)值的估計具有局限性,且仿真結(jié)果受初始值的影響較大,因此估計出的參數(shù)值并不一定是最優(yōu),未來我們考慮應(yīng)用馬爾科夫蒙特卡洛方法(Markov chain Monte Carlo,MCMC)進行估計,得到更為準確的參數(shù)值;2)間日瘧的傳播動態(tài)實際上更接近于隨機過程。本文所使用的是確定性模型,該模型能夠在長時間內(nèi)抓住間日瘧傳播動態(tài)的主要特征與總體趨勢,但對于具體年份上缺乏較好的描述。3)為了便于估計,我們沒有考慮各年間日瘧報告率的差異,報告率的差異會影響實際的發(fā)病數(shù)的多少。
安徽省作為中國間日瘧疫情不穩(wěn)定的省份[11],自90年代以來通過系列滅蚊措施,有效控制了間日瘧疫情,但是2000年后間日瘧發(fā)病數(shù)逐年上升并于2006年達到最高。2007年8月起安徽省在沿淮北地區(qū)的間日瘧重點防治縣實施預(yù)防服藥并采取多項綜合防治措施保證了間日瘧發(fā)病數(shù)穩(wěn)步下降[25]。盡管通過人工干預(yù)取得了消除間日瘧的階段性勝利,但是由于傳播媒介中華按蚊的存在以及間日瘧可在人體內(nèi)長期休眠等原因[26],本地病例反復(fù)出現(xiàn)的風險依舊存在,因此我們要保持對間日瘧的監(jiān)測,在通過省級消除瘧疾評估后仍要嚴防輸入再傳播,進一步防止本地疫情的反彈。此外據(jù)《第三次氣候變化國家評估報告》顯示,中國氣候變暖速率高于全球平均值。未來80年中國陸地區(qū)域增溫幅度將達到1.3~5.0℃[27]。由于溫度升高顯著增加間日瘧的傳染率,輸入性病例的增加以及累計的高易感者人數(shù)在增高的傳染率情況下是否會引起間日瘧的重新暴發(fā)是我們未來將要研究的內(nèi)容。