運(yùn)用空間計(jì)量模型偵測廣西結(jié)核病發(fā)病的氣象影響因素*

林定文， 崔哲哲，趙錦明， 區(qū) 進(jìn)

(廣西壯族自治區(qū)疾病預(yù)防控制中心結(jié)核病防制所，南寧 530028)

結(jié)核病是由攜帶結(jié)核分枝桿菌的微滴核進(jìn)入人體引起感染部位病變的慢性傳染病。結(jié)核病的流行已成為全世界，特別是包括中國在內(nèi)的發(fā)展中國家面臨的嚴(yán)重公共衛(wèi)生問題。2018年全球結(jié)核病報(bào)告指出，雖然估計(jì)發(fā)病率在緩慢下降，但中國是估計(jì)發(fā)病數(shù)僅次于印度的結(jié)核病高負(fù)擔(dān)國家之一[1]。在全國省級(jí)水平，廣西壯族自治區(qū)是結(jié)核病報(bào)告發(fā)病熱點(diǎn)地區(qū)，但其傳播和發(fā)病的特征和影響因素仍然不明?？臻g流行病學(xué)分析，可以有效地掌握疾病聚集的范圍和趨勢，加之空間計(jì)量統(tǒng)計(jì)對(duì)大氣等影響因素的探查，將對(duì)衛(wèi)生決策起到關(guān)鍵作用。本研究通過空間自相關(guān)統(tǒng)計(jì)量和空間回歸技術(shù)，對(duì)2016廣西壯族自治區(qū)活動(dòng)性肺結(jié)核患者空間特征和氣象影響因素進(jìn)行分析，為防控決策的制訂和衛(wèi)生資源的分配提供量化依據(jù)，現(xiàn)報(bào)道如下。

1 資料與方法

1.1一般資料 本研究是基于現(xiàn)住址定位的廣西壯族自治區(qū)全人群景觀流行病學(xué)研究。擬用2016年的分析結(jié)果闡述廣西壯族自治區(qū)傳染性結(jié)核病病例在該時(shí)期的空間分布和聚類情況，并運(yùn)用優(yōu)化空間回歸模型偵測影響結(jié)核病發(fā)病的氣象因素。人口數(shù)據(jù)來源于2017年廣西壯族自治區(qū)統(tǒng)計(jì)年鑒。PM2.5、PM10、SO2、NO2和O3數(shù)據(jù)來源于廣西壯族自治區(qū)環(huán)境保護(hù)局?jǐn)?shù)據(jù)庫。年降雨天數(shù)、氣溫和相對(duì)濕度數(shù)據(jù)來源于廣西壯族自治區(qū)氣象局。

1.2方法 研究地區(qū)廣西壯族自治區(qū)位于中國南疆，總?cè)丝? 613萬，占地236 700平方公里，下轄14個(gè)市和所轄的112個(gè)縣(區(qū))。從國家法定傳染病報(bào)告系統(tǒng)下載2016年所有確診的活動(dòng)性結(jié)核病患者資料，共計(jì)76 151例。根據(jù)每例患者的現(xiàn)住址的地區(qū)編碼進(jìn)行歸類統(tǒng)計(jì)，最小空間統(tǒng)計(jì)單元為市，即統(tǒng)計(jì)每個(gè)市的登記報(bào)告數(shù)。然后將每個(gè)單元的病例按照月和年再次歸類，形成2016年活動(dòng)性結(jié)核病數(shù)據(jù)庫。

從“Global Administrative area database” (GADM Inc,California,US) 中獲取中國地圖(Shape 文件)，從中剪取生成廣西Shape文件。將病例數(shù)據(jù)、環(huán)境數(shù)據(jù)與地圖文件根據(jù)地區(qū)編碼進(jìn)行關(guān)聯(lián)，生成結(jié)核病及氣象信息空間數(shù)據(jù)庫。

1.3統(tǒng)計(jì)學(xué)處理

1.3.1基本情況描述 對(duì)年度報(bào)告的活動(dòng)性結(jié)核病患者數(shù)據(jù)進(jìn)行總體描述并進(jìn)行趨勢χ2檢驗(yàn)；對(duì)各市的結(jié)核病報(bào)告發(fā)病率及大氣數(shù)據(jù)進(jìn)行中位數(shù)和均數(shù)描述。根據(jù)不同地區(qū)的結(jié)核病報(bào)告發(fā)病率繪制發(fā)病率專題圖，了解高發(fā)和低發(fā)地區(qū)分布，為空間自相關(guān)分析和空間掃描做準(zhǔn)備。

1.3.2空間自相關(guān)分析 本研究應(yīng)用統(tǒng)計(jì)學(xué)家帕特里克·阿爾弗雷德·皮爾斯·莫蘭創(chuàng)建的“Moran′s I”統(tǒng)計(jì)理論進(jìn)行空間相關(guān)性分析。該分析分為全局和局部自相關(guān)分析。

全局自相關(guān)分析通過方差歸一化，使研究區(qū)域內(nèi)總體相關(guān)系數(shù)(或稱莫蘭指數(shù))取值在-1～1。莫蘭指數(shù)大于0表示存在空間正相關(guān)性，指數(shù)值越大，空間相關(guān)性越顯著，也就是隨著相鄰的單元在空間上越聚集，相關(guān)性越明顯；莫蘭指數(shù)小于0表示存在空間負(fù)相關(guān)性，指數(shù)值越小，空間差異顯著，也就是相鄰的單元在空間上越離散，反而相關(guān)性變得顯著了；當(dāng)莫蘭指數(shù)等于0，空間呈隨機(jī)性，相鄰研究單位之間沒有相關(guān)性[2]。本研究將用總體莫蘭指數(shù)和象限圖進(jìn)行結(jié)果描述。象限圖將各相鄰坐標(biāo)點(diǎn)之間的相關(guān)系數(shù)用“high-high,高-高”“l(fā)ow-low,低-低”“high-low,高-低”和“l(fā)ow-high,低-高”在4個(gè)象限呈現(xiàn)。全局統(tǒng)計(jì)量公式如下：

式中W是空間權(quán)重系數(shù)，當(dāng)兩個(gè)相鄰地區(qū)(i和j)傳染性結(jié)核病均為高發(fā)(聚集)，W=1；若均為低發(fā)，W=-1。n為總研究單元數(shù)。若有顯著差異，W值則為0。S0為總體權(quán)重系數(shù)。

局部自相關(guān)分析是對(duì)每個(gè)研究個(gè)體的自相關(guān)分析及分解呈現(xiàn)。其在全局分析的基礎(chǔ)上，對(duì)局部地區(qū)的相鄰關(guān)系計(jì)算統(tǒng)計(jì)量，并根據(jù)不同的顯著性在LISA類聚地圖中呈現(xiàn)結(jié)果。一般認(rèn)為具有“高-高”關(guān)系的區(qū)域是顯著性高發(fā)類聚區(qū)(“熱點(diǎn)”地區(qū))，而具有“低-低”關(guān)系的區(qū)域是顯著性低發(fā)類聚區(qū)(“冷點(diǎn)”地區(qū))，“高-低”和“低-高”代表與周圍明顯差異的離群區(qū)域。本研究用“GeoDa 1.8.12” (Luc Anselin,University of ILLinois,Urbana-Champaign,US) 空間分析軟件進(jìn)行空間自相關(guān)分析。Anselin′s局部公式如下：

全局和局部自相關(guān)均用“經(jīng)驗(yàn)Bayes”標(biāo)化，以解決比率方差不穩(wěn)定的問題。

1.3.3空間回歸模型分析 結(jié)核病是傳染性疾病，缺乏空間獨(dú)立性，從而出現(xiàn)空間維度的播散，即空間計(jì)量統(tǒng)計(jì)學(xué)中提及的溢出效應(yīng)。用經(jīng)典線性回歸檢驗(yàn)影響因素時(shí)往往由于忽略了空間依賴性和空間差異效應(yīng)而出現(xiàn)失真。空間計(jì)量模型的基本思路是運(yùn)用空間距離權(quán)重矩陣，對(duì)普通的最小二乘法估計(jì)進(jìn)行修正。對(duì)于具備了空間屬性的變量，距離較近的觀察值比距離較遠(yuǎn)的觀察值更容易表現(xiàn)顯著的關(guān)系。本研究主要用于比較的空間計(jì)量模型有空間遲滯模型[3]和空間誤差模型[4]。其基本公式如下：

y=ρW1y+βx+μ

μ=λW2μ+ε

ε～N(0,δ2I)

式中的β表示自變量的空間回歸系數(shù);μ是空間誤差項(xiàng);ε是白噪聲；W1是結(jié)核病報(bào)告發(fā)病率的空間權(quán)重矩陣，W2是殘差空間權(quán)重矩陣；ρ為空間遲滯系數(shù)，取值0～1，越接近1，比鄰地區(qū)因變量值越相似；λ為空間誤差系數(shù)，取值0～1，越接近1，比鄰地區(qū)自變量取值越相似。因此，當(dāng)ρ=0，λ=0時(shí)，模型為普通的最小二乘法估計(jì)模型，模型中沒有空間特征；當(dāng)ρ≠0，λ=0時(shí)，模型為空間遲滯模型(或空間自回歸模型)，研究區(qū)域的因變量不但與本地區(qū)的自變量有關(guān)，還與比鄰區(qū)域因變量有關(guān)；當(dāng)ρ=0，λ≠0時(shí)，模型為空間誤差模型，研究區(qū)域的因變量與比鄰地區(qū)的因變量取值相互獨(dú)立，但比鄰地區(qū)同一種自變量殘差項(xiàng)不滿足獨(dú)立性，具有空間自相關(guān)性。

2 結(jié) 果

2.1廣西結(jié)核病及氣象一般情況 2016年廣西各醫(yī)療衛(wèi)生機(jī)構(gòu)共報(bào)告活動(dòng)性結(jié)核病患者47 701例。用2010-2016年度報(bào)告病例數(shù)進(jìn)行趨勢檢驗(yàn)，結(jié)果提示呈逐漸下降趨勢(χ2=159.76,P<0.01)。各市平均報(bào)告發(fā)病率為(109.33±29.46)/10萬，提示離散度較大。納入分析的大氣顆粒物中O3含量較高，其余均處于較低水平。年降雨天數(shù)，氣溫和相對(duì)濕度均符合濕熱地區(qū)氣候特征，見表1。

表1 廣西壯族自治區(qū)結(jié)核病報(bào)告發(fā)病和氣象基本情況

2.2結(jié)核病發(fā)病的空間分布和自回歸分析 2016年廣西壯族自治區(qū)以市為單位的結(jié)核病報(bào)告發(fā)病率呈現(xiàn)明顯的空間異質(zhì)性。西北部和西南部高于東部，見圖1。通過莫蘭空間自相關(guān)分析可知，空間發(fā)病分布具有正向自相關(guān)性(莫蘭指數(shù)為0.262)，提示實(shí)施影響因素分析時(shí)，考慮空間自相關(guān)效應(yīng)的必要性。比率經(jīng)過經(jīng)驗(yàn)貝葉斯標(biāo)化的局部聚類性分析，東部地區(qū)呈現(xiàn)顯著的低發(fā)聚類，中部地區(qū)呈現(xiàn)與周邊結(jié)核病疫情比較有顯著低發(fā)的離群態(tài)勢，未偵測到顯著的高發(fā)聚類，見圖2。

2.3空間影響因素和模型篩選 將2016年P(guān)M2.5、PM10、SO2、NO2、O3、年降雨天數(shù)、氣溫和相對(duì)濕度數(shù)據(jù)納入無空間權(quán)重和空間效應(yīng)的最小二乘法模型和具有空間權(quán)重的空間遲滯模型、空間誤差模型，按照P<0.05的條件進(jìn)行結(jié)核病報(bào)告發(fā)病率的潛在影響因素篩選。結(jié)果顯示，最小二乘法模型未篩選出有統(tǒng)計(jì)學(xué)意義的影響因素；空間遲滯回歸模型篩選出O3、年降雨天數(shù)和氣溫3個(gè)有統(tǒng)計(jì)學(xué)意義的變量；空間誤差模型篩選出年降雨天數(shù)和氣溫2個(gè)有統(tǒng)計(jì)學(xué)意義的變量，見表2。表3顯示了3種模型的擬合度測試結(jié)果。空間誤差模型的R2和對(duì)數(shù)釋然比(Log likelihood)最高，赤池信息準(zhǔn)則(AIC)統(tǒng)計(jì)值和回歸標(biāo)準(zhǔn)誤(S.E of regression)最低，入選最佳模型。由回歸系數(shù)可知，該模型中降雨天數(shù)及氣溫與結(jié)核病報(bào)告發(fā)病率呈負(fù)相關(guān)關(guān)系，其中氣溫的回歸效應(yīng)最大。

圖1 2016年廣西壯族自治區(qū)各市結(jié)核病報(bào)告發(fā)病率空間分布

圖2 2016年廣西壯族自治區(qū)結(jié)核病報(bào)告發(fā)病率的莫蘭局部聚類和象限圖(比率經(jīng)驗(yàn)Bayes標(biāo)化)

表2 3種模型的潛在影響因素篩選結(jié)果

表3 3種模型的擬合度比較

3 討 論

有研究數(shù)據(jù)表明，廣西壯族自治區(qū)的活動(dòng)性結(jié)核病報(bào)告發(fā)病近年來一直處于時(shí)空高發(fā)聚類區(qū)域[5]。廣西壯族自治區(qū)不同地區(qū)間的報(bào)告發(fā)病率差別迥異[6]。作為結(jié)核病高負(fù)擔(dān)省份，廣西壯族自治區(qū)盡管投入了大量的防控資源，疫情下降仍然緩慢，每年報(bào)告的活動(dòng)性結(jié)核病例依然維持在4萬左右[7]，空間熱點(diǎn)地區(qū)的結(jié)核病疫情仍然持續(xù)高于其他地區(qū)，而冷點(diǎn)地區(qū)依據(jù)處于低發(fā)態(tài)勢[6]。作為主要通過空氣飛沫傳播的疾病，結(jié)核病具備一定的空間格局，特別是臨近區(qū)域的聚集發(fā)病，結(jié)核病的空間聚集性疫情時(shí)有報(bào)道。如果不及時(shí)控制，很有可能導(dǎo)致更大范圍的蔓延，促成突發(fā)公共衛(wèi)生事件等嚴(yán)重后果[8-9]。本研究通過空間描述和自相關(guān)分析證明了廣西壯族自治區(qū)的結(jié)核病報(bào)告發(fā)病同樣存在空間異質(zhì)性(莫蘭指數(shù)大于0)，東部地區(qū)出現(xiàn)了顯著的低發(fā)聚集區(qū)。對(duì)于廣西壯族自治區(qū)的這一空間布局，筆者需要探究其因，用以指導(dǎo)防控。防治結(jié)核病的關(guān)鍵措施除了針對(duì)個(gè)體進(jìn)行早期的診斷、治療和管理，宏觀層面的防控實(shí)施也是很有必要的。社會(huì)及自然環(huán)境在結(jié)核病傳播和發(fā)病驅(qū)動(dòng)力方面起著舉足輕重的作用。但此前針對(duì)廣西壯族自治區(qū)的環(huán)境影響因素研究尚屬空白。因此，本研究主要以空間自相關(guān)為出發(fā)點(diǎn)，結(jié)合回歸理論對(duì)廣西壯族自治區(qū)結(jié)核病的聚類及其氣象影響因素進(jìn)行分析，以期獲得大氣環(huán)境層面的宏觀驅(qū)動(dòng)因子，為結(jié)核病精準(zhǔn)防控提供量化依據(jù)。

本研究主要引用了普通的最小二乘法模型(不考慮空間依賴性)、空間遲滯回歸和空間誤差模型進(jìn)行氣象相關(guān)影響因素的篩選和模型擬合度比較。在近期關(guān)于結(jié)核病發(fā)病影響因素的空間回歸模型研究中不同地區(qū)和不同的測試變量均有可能選擇不同的優(yōu)化模型。有研究者在對(duì)2005-2014年影響北京市結(jié)核病發(fā)病的社會(huì)經(jīng)濟(jì)因素進(jìn)行空間模型擬合時(shí)，發(fā)現(xiàn)空間誤差模型是最優(yōu)模型[5]。也有學(xué)者對(duì)2004-2012年四川結(jié)核病空間影響因素調(diào)查發(fā)現(xiàn)，空間遲滯模型擬合度更好[10]。本研究通過實(shí)證研究發(fā)現(xiàn)，空間誤差模型的擬合度優(yōu)于最小二乘法和空間遲滯模型(R2=0.763，Log likelihood=-59.444；AIC=136.888，S.E of regression=14.334)。因此，將該模型作為最優(yōu)模型進(jìn)行因素篩選。

通過氣象潛在影響因素篩選，年降雨天數(shù)和氣溫2個(gè)變量有統(tǒng)計(jì)學(xué)意義(P<0.05)，與結(jié)核病報(bào)告發(fā)病率呈負(fù)相關(guān)。降雨量較少而且溫度較低的地區(qū)，結(jié)核病發(fā)病風(fēng)險(xiǎn)相對(duì)較高。這一結(jié)果符合結(jié)核分枝桿菌存活的自然優(yōu)勢環(huán)境，也與其他的研究結(jié)果一致[11]。結(jié)核分枝桿菌主要以微滴核的方式在自然環(huán)境中存活，一般在干痰中存活6～8個(gè)月，如處于3 ℃的低溫環(huán)境下，可存活1年。若黏附于塵埃上，可以保持傳染性8～10 d。降雨可在一定程度上清洗飄浮在空氣中的帶菌塵埃，從而降低傳染和發(fā)病的可能性。因此對(duì)于廣西壯族自治區(qū)北部溫度和降雨量相對(duì)較少的地區(qū)應(yīng)加強(qiáng)結(jié)核病的監(jiān)測和防控，在干旱季節(jié)做好人工降雨等工作。由大氣顆粒物監(jiān)測數(shù)據(jù)可知，廣西壯族自治區(qū)的PM2.5、PM10、SO2、NO2均處于較低水平，經(jīng)換算的空氣質(zhì)量指數(shù)屬于良好的范圍[12-13]，可正常戶外活動(dòng)，也未對(duì)結(jié)核病的發(fā)病構(gòu)成危險(xiǎn)。但O3在大氣中的含量超過了100 μg/m3，空間遲滯回歸提示了O3有影響結(jié)核病發(fā)病的可能。在今后應(yīng)針對(duì)這些顆粒物進(jìn)一步加強(qiáng)監(jiān)測和相關(guān)性研究。此外，也有研究提示結(jié)核分枝桿菌的基因多態(tài)性與所處的環(huán)境可能存在一定的交互作用，從而導(dǎo)致傳播的差異性[14-16]。但目前這一觀點(diǎn)仍然缺乏更多的實(shí)證研究。結(jié)核病在分子維度和宏觀環(huán)境維度的交互作用及驅(qū)動(dòng)因子研究將是下一步的重點(diǎn)探索方向。