張書源, 程全國, 邢紅彬
(沈陽大學 a. 環(huán)境學院, b. 大學生創(chuàng)新創(chuàng)業(yè)指導中心, 遼寧 沈陽 110044)
大氣污染是影響人類健康的主要環(huán)境危害因素之一。在公認的大氣污染物中,可吸入顆粒物PM10、PM2.5與人群健康效應終點的流行病學聯(lián)系最為密切[1-2]。近年來,對可吸入顆粒物的定量健康危害評價已成為國際社會關(guān)注的熱點之一[3]。
大量的研究結(jié)果表明,空氣可吸入顆粒物污染與居民每日死亡總?cè)藬?shù)、呼吸系統(tǒng)與心血管等疾病的超額死亡數(shù)顯著相關(guān)[4-5],即使在空氣顆粒物污染濃度低于標準的情況下,隨著顆粒物濃度的升高,呼吸系統(tǒng)和心腦血管疾病的發(fā)病率和死亡率亦呈現(xiàn)上升趨勢[6-8]。
數(shù)據(jù)挖掘是指從大量的數(shù)據(jù)中通過算法與技術(shù)搜索隱藏于其中的信息的過程。數(shù)據(jù)挖掘中經(jīng)常用到的算法和技術(shù)有統(tǒng)計分析、決策樹、貝葉斯網(wǎng)絡、粗糙集、人工神經(jīng)網(wǎng)絡等。隨著生態(tài)文明建設和生態(tài)環(huán)境保護進入了數(shù)據(jù)驅(qū)動的新時代,數(shù)據(jù)挖掘技術(shù)在生態(tài)環(huán)境保護領(lǐng)域中應用的廣度和深度都不斷加強。
目前國內(nèi)外研究人員運用流行病學和毒理學的研究方法探索大氣污染對死亡率的影響,并建立暴露反應關(guān)系模型[9-10],對于大氣污染與各死因類型的相關(guān)分析研究相對較少。本文以遼寧省本溪市2019年的死因數(shù)據(jù)、大氣污染和氣象監(jiān)測數(shù)據(jù)作為數(shù)據(jù)來源,通過數(shù)據(jù)挖掘技術(shù)分析PM2.5污染與各死因類型的顯著性與相關(guān)性,并以此為基礎(chǔ)建立PM2.5污染與每日死亡人數(shù)的暴露-反應關(guān)系模型,為促進數(shù)據(jù)挖掘技術(shù)在生態(tài)環(huán)境保護領(lǐng)域中的應用、制定有效的環(huán)境法規(guī)、保護易感人群提供科學依據(jù)。
分析所用數(shù)據(jù)包括數(shù)2019年1月1日—12月31日本溪市的死因數(shù)據(jù)、大氣污染和氣象數(shù)據(jù)。死因數(shù)據(jù)來源于本溪市疾控中心,死因按照國際疾病分類ICD-10進行分類;大氣污染數(shù)據(jù)來自中國空氣質(zhì)量在線監(jiān)測平臺(https:∥www.aqistudy.cn/historydata/)的監(jiān)測數(shù)據(jù),包括PM2.5、PM10、SO2、CO、NO2、O3的每日平均質(zhì)量濃度(μg·m-3);氣象數(shù)據(jù)來自國家氣象科學數(shù)據(jù)中心(http:∥data.cma.cn/)本溪站的監(jiān)測數(shù)據(jù),包括每日平均溫度(℃)、相對濕度(%)、平均氣壓(kPa)、平均風速(m·s-1)。
數(shù)據(jù)處理按照數(shù)據(jù)挖掘技術(shù)中統(tǒng)計分析的一般步驟,即數(shù)據(jù)清洗、缺失值填充、z-score數(shù)據(jù)標準化方法的順序?qū)υ紨?shù)據(jù)進行處理。由于死因數(shù)據(jù)中的ICD-10死因分類類別繁雜,為了避免特征維度的擴張與方便統(tǒng)計分析軟件識別各類死因,根據(jù)ICD-10分類標準對本溪市2019年1月1日—12月31日的死因進行分類,并將死因進行數(shù)值化編碼,如表1所示。
表1 死因分類數(shù)值化編碼Table 1 Numerical coding for the classification of causes of death
大氣污染流行病學研究中常用時間序列的廣義相加泊松(Poisson)回歸模型評價大氣污染短期暴露對人群的健康影響,由于對于總?cè)丝趤碚f,每日死亡人數(shù)是小概率事件,作為一種時間序列模型,其實際分布近似泊松分布。時間序列的廣義相加泊松回歸模型是對傳統(tǒng)廣義對數(shù)線性模型的進一步拓展,模型中除擬合普通的線性項(如大氣污染物濃度)外,還可將一些與因變量之間存在復雜非線性關(guān)系的變量(如時間序列資料中的長期趨勢、季節(jié)和其他一些與時間長期變異有關(guān)的混雜因子)以不同函數(shù)加和的形式擬合模型[11]。廣義相加Poisson回歸模型的一般形式為
lnE(Y)=βX。
式中:E(Y)為日死亡人數(shù)預測值;X為大氣污染物濃度及氣象因素數(shù)據(jù)等一些其他變量構(gòu)成的矩陣,包括截距列;β為模型自變量的系數(shù)向量。
通過PM2.5污染對各死因的差異顯著性分析,判斷PM2.5污染對不同類型的死因是否具有顯著性影響,在此基礎(chǔ)上將死亡人數(shù)的自然對數(shù)與大氣污染物指標、氣象指標進行相關(guān)性分析,并考慮滯后日效應,選擇自然對數(shù)與PM2.5和其他指標存在顯著關(guān)聯(lián)的變量進行線性回歸分析,從而建立PM2.5污染與居民每日死亡人數(shù)的暴露-反應關(guān)系。
本溪市2019年1月1日—12月31日的PM2.5日均質(zhì)量濃度在各死因分類(編碼1~20)中的分布情況如表2所示。
表2 PM2.5日均質(zhì)量濃度在各死因分類中的分布情況Table 2 Distribution of PM2.5 daily average mass concentration in the classification of death causes
從表2中可以看出,2019年本溪市居民死因中的腫瘤,內(nèi)分泌、營養(yǎng)和代謝疾病,循環(huán)系統(tǒng)疾病,呼吸系統(tǒng)疾病及消化系統(tǒng)疾病為主要死因,因循環(huán)系統(tǒng)疾病導致死亡的人數(shù)最多,因皮膚和皮下組織疾病導致死亡的人數(shù)最少,PM2.5質(zhì)量濃度在死因分類各組中的均值在35~50 μg·m-3之間波動。
顯著性P>0.05表示差異性不顯著;0.01≤P≤0.05表示差異性顯著;P<0.01表示差異性極顯著。本文在研究本溪市PM2.5污染對每日死亡人數(shù)的影響時,首先進行基于死因分組的大氣污染物指標均值差異顯著性分析,判斷PM2.5等大氣污染物對不同死因是否具有差異顯著性。
利用SPSS軟件,將數(shù)值編碼后的死因分類(編碼1~20)作為因子,將大氣污染物數(shù)據(jù)中PM2.5、PM10、SO2、CO、NO2、O3的每日平均質(zhì)量濃度作為因變量進行差異顯著性分析,結(jié)果如表3所示。從表3可以看出,PM2.5、PM10、SO2、CO、NO2污染的顯著性大于0.05,O3污染的顯著性小于0.05且接近0.01。分析結(jié)果說明只有O3污染對不同類型的死因具有顯著影響,PM2.5及其他大氣污染物對不同類型的死因影響不顯著,可以認為PM2.5對不同類型的死因影響相同。因此本文在分析本溪市PM2.5污染對每日死亡人數(shù)的影響時,應對全部死因類型的死亡人數(shù)進行分析。
表3 大氣污染物在各死因組間的均值差異顯著性分析Table 3 Significant analysis of the difference in the mean value of air pollutants among the death cause groups
差異顯著性分析除了對研究的總體均值的差異進行顯著性檢驗以外,還需要對2個獨立樣本所屬總體的總體方差的差異進行顯著性檢驗,統(tǒng)計學上稱為方差齊性檢驗。
利用SPSS軟件,基于死因分組的大氣污染物方差齊性檢驗結(jié)果如表4所示,從表4中可以看出,除O3外,其他大氣污染物的顯著性均大于0.05,可以認為除O3外的其他大氣污染物在各組死因中的方差相等,O3在各組死因中的方差差別較大,說明PM2.5在各組死因中數(shù)據(jù)波動一致,進一步證明PM2.5對不同類型的死因影響不顯著。
表4 大氣污染物在各死因組間的方差齊性檢驗Table 4 Test for the homogeneity of variance of air pollutants among the death cause groups
相關(guān)性分析是指對2個或多個具備相關(guān)性的變量元素進行分析,從而衡量2個變量因素的相關(guān)密切程度。相關(guān)系數(shù)是研究變量之間線性相關(guān)程度的量,較為常用的是皮爾遜相關(guān)系數(shù),系數(shù)為正值表示變量間為正相關(guān)關(guān)系,負值則表示為負相關(guān)關(guān)系,系數(shù)絕對值越大,表示關(guān)系緊密程度越高,本文相關(guān)性分析結(jié)果中的相關(guān)系數(shù)均為皮爾遜相關(guān)系數(shù)。
為探究大氣污染物指標、氣象指標與每日全因死亡人數(shù)(D)的相關(guān)性,利用SPSS軟件對本溪市2019年1月1日—12月31日的大氣污染物指標中的PM2.5、PM10、SO2、CO、NO2、O3每日平均質(zhì)量濃度以及氣象指標中的平均溫度、相對濕度、平均氣壓、平均風速,與每日全因死亡人數(shù)作為因子,進行雙變量相關(guān)性分析,并考慮滯后日效應。滯留日是在時間序列資料分析中常用的概念,其含義是用今日的健康指標與前期若干日的氣象或大氣污染質(zhì)量濃度值進行分析,從而研究前幾日的氣象情況或大氣污染對以后健康問題的影響。本文使用滯后0(即當天)~7d的大氣污染物數(shù)據(jù)、氣象數(shù)據(jù)與每日全因死亡人數(shù)進行相關(guān)性分析,分析結(jié)果如表5所示。
表5 滯后0~7 d日全因死亡人數(shù)(D)與各指標相關(guān)性分析Table 5 Analysis of the correlation between the number of daily all-cause deaths (D) and various indicators after a lag of 0~7 days
表5中,相關(guān)系數(shù)右上角的星號(*)代表顯著性水平,即假設變量間具有線性相關(guān)的可能性,有星號則表示可能有關(guān)系,否則表示無關(guān)系。通過表5可知,滯后0、3、6、7 d,日全因死亡人數(shù)與PM2.5的日平均質(zhì)量濃度顯著性均小于0.01,表明假設結(jié)果極顯著,具有相關(guān)性,日全因死亡人數(shù)與PM2.5污染具有明顯的統(tǒng)計學意義。同時,PM10、SO2、CO、NO2日平均質(zhì)量濃度及平均溫度、平均氣壓也與日全因死亡人數(shù)具有統(tǒng)計學意義。
為進一步探究日全因死亡人數(shù)與PM2.5污染的內(nèi)在關(guān)聯(lián),需要繼續(xù)進行偏相關(guān)分析。偏相關(guān)分析是指當2個變量同時與第3個變量相關(guān)時,將第3個變量的影響剔除,只分析另外2個變量之間相關(guān)程度的過程。通過偏相關(guān)分析,控制氣象因素(平均溫度、相對濕度、平均氣壓、平均風速)的影響,只分析日全因死亡人數(shù)與大氣污染物指標的內(nèi)在相關(guān)性,利用SPSS軟件分別對當天、滯后3、6、7 d的數(shù)據(jù)進行偏相關(guān)分析,分析結(jié)果如表6。通過表6可以看出,PM2.5污染的相關(guān)性由大到小為:當天、滯后6 d、滯后3 d、滯后7 d,顯著性由大到小為:滯后7 d、滯后3 d、滯后6 d、當天。當天PM2.5污染的顯著性相對其他組最小,且小于0.01,相關(guān)系數(shù)相對其他組最大。說明當天的PM2.5數(shù)據(jù)與日全因死亡人數(shù)具有最大線性相關(guān)的可能性,并且相關(guān)程度也最高。因此本文將采用死亡日當天的大氣污染物數(shù)據(jù)、氣象數(shù)據(jù)與日死亡人數(shù)進行回歸分析。
表6 日全因死亡人數(shù)(D)與各指標偏相關(guān)分析Table 6 Analysis of partial correlation between daily all-cause deaths (D) and various indicators
利用SPSS軟件將死因數(shù)據(jù)中日死亡人數(shù)的自然對數(shù)lnD的計算結(jié)果與大氣污染物指標及氣象指標進行相關(guān)分析,結(jié)果見表7。通過表7可以看出,lnD與大氣污染物指標中的PM2.5、PM10、SO2、NO2的日平均質(zhì)量濃度以及氣象指標中的平均溫度、平均氣壓P0的顯著性小于0.01,具有顯著的線性相關(guān)可能性與明顯的統(tǒng)計學意義;lnD與CO日平均質(zhì)量濃度、平均濕度的顯著性大于0.01,線性相關(guān)可能性軟弱;lnD與O3日平均質(zhì)量濃度、平均風速無線性相關(guān)性。因此,PM2.5、PM10、SO2、NO2日平均質(zhì)量濃度、日平均溫度、日平均氣壓P0均可作為預測日死亡人數(shù)模型的主要變量。
表7 ln D與各指標相關(guān)分析Table 7 Correlation analysis of ln D and each index
主成分分析是一種將多個指標化為少數(shù)幾個不相關(guān)的綜合指標(即主成分)的統(tǒng)計分析方法,該方法適于分析多指標的大量數(shù)據(jù)間的關(guān)系和趨勢,也是數(shù)學上用來降維的一種方法。主成分載荷是主成分分析中原始指標與主成分之間的相關(guān)系數(shù),主成分載荷的正、負決定原始指標對主成分影響的方向是正面還是負面的,主成分載荷的絕對值大小表示原始指標對主成分影響的強度,絕對值越大表示原始指標對主成分越具有代表性。
為使模型的自變量降維,從而達到簡化模型的作用,利用SPSS軟件對PM2.5、PM10、SO2、NO2的日平均質(zhì)量濃度、日平均溫度、日平均氣壓P06個評價指標進行主成分分析,提取出了2個不相關(guān)的綜合指標作為主成分進行分析,主成分載荷矩陣如表8所示。從表8中可以看出,在第1主成分中的最大載荷為ρPM2.5,為0.934,在第2成分中的最大載荷為P0,為0.955,說明第1主成分中ρPM2.5對該主成分最具有代表性,第2主成分中P0對該主成分最具有代表性。因此,為使模型自變量降維,選取PM2.5的日平均質(zhì)量濃度ρPM2.5和日平均氣壓P0與lnD進行回歸分析。
表8 評價指標主成分載荷矩陣Table 8 Evaluation index principal component load matrix
利用SPSS軟件的回歸分析功能,將本溪市2019年1月1日—12月31日的死亡數(shù)據(jù)中日死亡人數(shù)的自然對數(shù)lnD作為因變量,將評價指標中PM2.5日平均質(zhì)量濃度ρPM2.5、日平均氣壓P0作為自變量進行回歸分析,擬合模型系數(shù)及95%置信區(qū)間如表9所示。擬合模型調(diào)整后的R2為0.992,大于0.8,說明擬合模型能夠解釋因變量99.2%的變異,擬合精度高。表9顯示評價指標ρPM2.5回歸系數(shù)為0.001 1(0.000 3~0.001 9),評價指標P0的回歸系數(shù)為0.027 8(0.027 4~0.028 2),2個回歸系數(shù)估計值的95%置信限均分別在0的同側(cè),且2個回歸系數(shù)P值均小于0.05,表明2個回歸系數(shù)均具有顯著性與明顯的統(tǒng)計學意義。根據(jù)分析結(jié)果,以泊松回歸模型為基本模型,擬合本溪市PM2.5污染與死亡人數(shù)的暴露-反應關(guān)系模型表達式為
表9 系數(shù)及95%置信區(qū)間結(jié)果Table 9 Coefficient and 95% confidence interval results
lnD=0.001 1ρPM2.5+0.027 8P0。
本文研究發(fā)現(xiàn),本溪市PM2.5污染對各類死因的影響不顯著,O3對各類死因的影響顯著;PM2.5日平均質(zhì)量濃度、日平均氣壓與當日居民全因死亡人數(shù)的自然對數(shù)之間存在顯著關(guān)聯(lián)。當PM2.5日平均質(zhì)量濃度上升10 μg·m-3時,帶入擬合的日死亡人數(shù)泊松回歸模型中計算可得本溪市日死亡人數(shù)上升1.1%,而我國其他城市的同類研究結(jié)果在0.36%~0.85%之間[12-14]。研究結(jié)果表明,本溪市居民死亡率受PM2.5污染的影響程度相比國內(nèi)其他城市更為嚴重。