基于監(jiān)測及Kriging方法的京津冀地區(qū)大氣污染物暴露分布研究

王占山， 李志剛， 錢 巖， 李曉倩， 郭 辰， 薛凱兵， 王孜曄， 朱曉晶， 魏永杰

中國環(huán)境科學研究院， 環(huán)境基準與風險評估國家重點實驗室， 北京 100012

隨著我國工業(yè)化和城市化進程的加快，城市生態(tài)環(huán)境日益惡化，成為制約經(jīng)濟發(fā)展、危害人體健康的重要因素[1-3]. 關于大氣污染物對人體暴露風險的評估成為當前環(huán)境領域研究的熱點和重點[4]. 空氣污染暴露風險評估依賴于連續(xù)一致的濃度空間分布，以及盡量精細化到每個個體的暴露水平估算，但受限于目前我國儀器設備、資金和人力等因素，很難實現(xiàn)采用每個個體佩戴暴露裝置的措施，同時，空氣質(zhì)量監(jiān)測點位的數(shù)量和密度只能觀測有限空間范圍內(nèi)的污染物濃度，不足以支撐科學全面的污染精細化暴露評估[5]. 以往在進行人群暴露水平研究時，常用以下兩種方式進行：①假設所有污染物在空間內(nèi)均具有相同濃度，用某一個站點的數(shù)據(jù)代表一個城市的暴露水平；②用全市多個站點的平均數(shù)據(jù)代表所關注人群的暴露水平. 這兩種方法均存在較大的誤差和不確定性.

為解決上述問題，國內(nèi)外學者通常使用數(shù)值模擬、空間插值、遙感反演、LUR (Land Use Regression)模型等方法對有限固定站點的污染物監(jiān)測濃度進行空間化[6]. 在數(shù)值模擬方面，Buonocore等[7]使用CMAQ模型模擬了美國電廠排放的PM2.5對人體健康的影響；Arunachalam等[8]使用CMAQ模型模擬了美國機場飛機起飛和降落時排放的PM2.5對人體健康的影響；CHEN等[9]使用CMAQ模型模擬了華北地區(qū)ρ(PM2.5)的空間分布和區(qū)域傳輸特征. 在空間插值方面，Pearce等[10]使用Kriging(克里金插值)方法研究了Cusco地區(qū)ρ(PM2.5)的空間分布特征；Sampson等[11]使用Kriging方法預測了美國地區(qū)ρ(PM2.5)；ZOU等[12]從多角度對比了基于Kriging模型和基于土地利用回歸(LUR)模型的ρ(PM2.5)空間分布特征. 在遙感反演方面，WANG等[13]基于衛(wèi)星遙感獲得了Montreal地區(qū)高分辨率的ρ(PM2.5)空間分布；Kloog等[14]基于衛(wèi)星遙感AOD數(shù)據(jù)反演了新西蘭地區(qū)ρ(PM2.5)并評價其暴露水平；HE等[15]基于衛(wèi)星遙感數(shù)據(jù)獲得了中國東部地區(qū)長時間系列的ρ(PM2.5). 在LUR模型方面，Jeanette等[16]基于LUR模型預測了美國地區(qū)基于面源和道路源的ρ(PM2.5)空間分布；Ho等[17]基于LUR模型預測了韓國地區(qū)顆粒物及其化學組分的空間分布；Moore等[18]使用LUR模型預測了洛杉磯地區(qū)的ρ(PM2.5)及空間分布.

上述方法均存在一定的不足，如數(shù)值模擬方法需要依賴高精度的大氣污染物排放清單和高性能運算資源，衛(wèi)星遙感方法難以解決有云遮擋的數(shù)據(jù)缺失問題，LUR方法對輸入基礎數(shù)據(jù)的要求較高，包括高精度的氣象數(shù)據(jù)、人口數(shù)據(jù)、污染源排放清單數(shù)據(jù)、地形數(shù)據(jù)等資料，數(shù)據(jù)的可獲取性較差. 因此，該研究擬以保定市作為京津冀地區(qū)試點城市，進行空氣污染采樣監(jiān)測，開展站點代表性和空間預測模型適用性的評估，旨在獲取科學有效且可操作性較強的適合中國地區(qū)的污染物濃度空間預測模型，為開展精細化大氣污染物人體健康暴露評估提供有力支撐.

1 研究方法

1.1 監(jiān)測點位及監(jiān)測方法

在河北省保定市選擇了13個點位開展空氣污染監(jiān)測，點位布設以位于保定市區(qū)的河北大學為中心點，在東、南、西、北4個方位每隔3 km左右設置一個監(jiān)測點位(見圖1). 采樣時間為2018年1月3—18日，采樣期間空氣質(zhì)量包括4個等級：良(4 d)、輕度污染(3 d)、中度污染(4 d)和重度污染(4 d).

圖1 保定市監(jiān)測點位分布Fig.1 Locations of the monitoring sites in Baoding

每個點位在線監(jiān)測SO2、NO2、O3、CO、PM10和PM2.5等6項常規(guī)污染物監(jiān)測數(shù)據(jù)，并采集PM2.5濾膜和VOCs蘇瑪罐. 在線設備采用中國環(huán)境科學研究院自行研發(fā)生產(chǎn)的集成傳感器專利產(chǎn)品，主要部件為英國阿爾法傳感器. 顆粒物監(jiān)測原理為光散射方法，氣態(tài)污染物監(jiān)測原理為電化學方法. 經(jīng)過清華大學建筑環(huán)境檢測中心的測試，該傳感器對6項污染物的監(jiān)測結(jié)果與美國進口的TSI品牌傳感器監(jiān)測結(jié)果的相關系數(shù)達到0.92，相對誤差為8%，儀器間平行度為98%. 采集過程中設定時間分辨率1 h. PM2.5濾膜采樣使用SKC品牌采樣泵和采樣頭采集，使用Whatman品牌37 mm濾膜，每個點位采集石英膜和特氟龍膜各1張. VOCs使用ENTH品牌蘇瑪罐采集. 在線設備為24 h連續(xù)監(jiān)測，PM2.5濾膜每天采集23 h，VOCs為每天12：00采集瞬時樣品. 在實驗室分析中，PM2.5組分中的水溶性離子、OCEC和元素分析分別使用離子色譜儀(ICS-600, Dionex, USA)、熱光碳分析儀(RT-4, Sunset Lab, USA)和波長色散X射線熒光光譜儀(XRF-1800, 島津, 日本)，蘇瑪罐VOCs樣品使用氣相色譜質(zhì)譜聯(lián)用儀(GC-MS)(Nexis GC-2030, 島津, 日本)進行分析. 在實驗室分析的質(zhì)保質(zhì)控方面，各物種檢測規(guī)范必須嚴格遵守國家相關標準規(guī)范，至少設置10%的空白樣和平行樣，標準曲線相關系數(shù)大于0.999[19].

1.2 數(shù)據(jù)處理方法

對于SO2、NO2、O3、CO、PM10和PM2.5等6項常規(guī)污染物監(jiān)測數(shù)據(jù)，時間分辨率為小時，使用監(jiān)測時段內(nèi)小時數(shù)據(jù)的平均值獲得點位均值. 對于PM2.5組分監(jiān)測數(shù)據(jù)，使用23 h采樣時長的平均值作為日均值，使用監(jiān)測時段內(nèi)日均值的平均值獲得點位均值. 對于VOCs監(jiān)測數(shù)據(jù)，使用每天12：00采集的瞬時值作為日均值，使用監(jiān)測時段內(nèi)日均值的平均值獲得點位均值.

注： 虛線表示中心點濃度.圖2 不同點位φ(SO2)、φ(NO2)、φ(CO)、φ(O3)Fig.2 Concentrations of SO2, NO2, CO and O3 at different sites

2 結(jié)果與討論

將12個點位的各污染物監(jiān)測濃度與中心點濃度進行對比，該研究中認為各方向上污染物平均濃度與中心點濃度的絕對誤差的平均值與中心點濃度的偏差不超過20%時，說明中心點濃度可以代表此點位的濃度.

2.1 不同氣態(tài)污染物差異性比較

表1和圖2顯示了不同點位φ(SO2)、φ(NO2)、φ(CO)、φ(O3)以及與中心點的濃度差異. 對φ(SO2)來說，點位1(指各方向的第一個點位，下同)、點位2與中心點的濃度差均在20%以下，點位3與中心點的濃度偏差較大，在40%左右，特別是南3和西3兩個點位，其濃度明顯高于中心點. 對φ(NO2)來說，整體規(guī)律與φ(SO2)類似，點位1、點位2與中心點的濃度偏差均低于20%，點位3與中心點的濃度偏差高于20%，φ(NO2)較高的是北3和東3兩個點位. 對φ(CO)和φ(O3)來說，不同點位之間的濃度偏差沒有明顯差異和規(guī)律性，且均小于20%.

表1 各點位污染物濃度與中心點濃度的差異

2.2 不同點位ρ(PM10)、ρ(PM2.5)及其主要組分差異性比較

表2和圖3顯示了各點位ρ(PM10)、ρ(PM2.5)與中心點的濃度差異，二者變化規(guī)律較為一致，即點位1、點位2與中心點的濃度偏差均在20%以下，點位3與中心點的濃度偏差均高于20%. 特別是ρ(PM10)，點位 3的濃度明顯高于中心點，ρ(PM2.5)空間分布更為均勻.

表2 各點位顆粒物濃度與中心點濃度的差異

注：虛線表示中心點濃度.圖3 不同點位ρ(PM10)和ρ(PM2.5)Fig.3 Concentrations of PM10 and PM2.5 at different sites

選取PM2.5組分中占比較高的SO42-、NO3-和NH4+[20]開展研究. 表3和圖4顯示了各組分濃度與中心點的濃度差異.ρ(SO42-)和ρ(NH4+)的變化規(guī)律較為一致，即點位1、點位2與中心點的濃度偏差均低于20%，點位3的濃度偏差高于20%. 對于ρ(NO3-)來說，各點位濃度分布較為均勻，濃度偏差均低于20%. 綜上，對于ρ(SO42-)和ρ(NH4+)而言，監(jiān)測點位的代表性為6 km左右；對于ρ(NO3-)而言，監(jiān)測點位的代表范圍更大.

綜合來看，對于顆粒物及其大部分組分而言，監(jiān)測點位的代表性為5～6 km.

表3 各點位顆粒物組分濃度與中心點濃度的差異

注：虛線表示中心點濃度.圖4 不同點位PM2.5中ρ(SO42-)、ρ(NO3-)和ρ(NH4+)Fig.4 Concentrations of SO42-, NO3- and NH4+ in PM2.5 at different sites

2.3 不同點位φ(VOCs)差異性比較

圖5顯示了各點位φ(VOCs)與中心點的濃度差異. 點位1、點位2、點位3與中心點的濃度偏差分別為17.5%、5.6%、10.0%，表明不同點位φ(VOCs)有所差異，但整體來看，各點位φ(VOCs)與中心點的平均濃度偏差維持在20%以內(nèi).

注： 虛線表示中心點濃度. 圖5 不同點位φ(VOCs)Fig.5 Concentrations of VOCs at different sites

2.4 5 km以外的污染物濃度預測方法對比

該研究得出，對于PM2.5等主要污染物來說，監(jiān)測點位具有5 km左右的代表性. 根據(jù)生態(tài)環(huán)境部發(fā)布的《環(huán)境空氣質(zhì)量監(jiān)測規(guī)范》(試行)[21]中的建議，空氣質(zhì)量監(jiān)測點位的空間代表性一般為500～4 000 m，與筆者所得結(jié)論較為一致. 由于目前我國的空氣質(zhì)量監(jiān)測網(wǎng)絡尚比較落后，并沒有實現(xiàn)每隔5 km都設有一個監(jiān)測點位[22]，因此在5 km以外沒有監(jiān)測點位的區(qū)域，需要尋求其他方法來獲得空氣污染物的濃度. 該研究中選取了國內(nèi)外研究中常用的Radial Basis Functions (徑向基函數(shù)，RBF)、Local Polynomial Interpolation(局部多項式插值，LPI)、Inverse Distance Weighting(反距離權重插值，IDW)、Kriging(克里根插值)、Kernel Smoothing(內(nèi)核平滑插值，KS)和Diffusion Kernel(內(nèi)核擴散插值，DK)等6種空間分析方法[23-26]，在保定市開展污染物濃度插值和驗證，選取準確率最高的方法. 使用13個點位中的中心點、點位1和點位3進行插值，使用點位2進行驗證，即提取空間插值的點位2的污染物濃度與監(jiān)測實況進行對比. 該研究以ρ(PM2.5)為例進行分析.

2.4.1RBF方法

表4和圖6顯示了基于RBF方法的ρ(PM2.5)預測值和監(jiān)測值的差異，可以看出，該模型預測的濃度高值區(qū)在研究區(qū)域的東南側(cè)，低值區(qū)集中在研究區(qū)域的中心位置，研究區(qū)域西側(cè)預測的ρ(PM2.5)處于中等水平. 提取4個點位濃度進行比對發(fā)現(xiàn)，模型對各點位ρ(PM2.5)的預測值偏高，其中南2點位的絕對偏差最小，為7.9%，西2點位的絕對偏差最大，為26.1%，4個點位的平均絕對偏差為17.9%.

表4 基于RBF方法的不同點位ρ(PM2.5)預測值與監(jiān)測值的差異

圖6 基于RBF方法的空間插值Fig.6 Spatial interpolation based on the RBF method

2.4.2LPI方法

表5和圖7顯示了基于LPI方法的ρ(PM2.5)預測值和監(jiān)測值的差異，可以看出，該模型預測濃度的低值區(qū)集中在研究區(qū)域的中心位置，隨著距離向外擴散預測濃度有所升高，高值區(qū)集中在研究區(qū)域的邊界位置. 提取4個點位濃度發(fā)現(xiàn)，模型對北、南、西點位的ρ(PM2.5)預測值偏高，對東側(cè)點位的預測值偏低，其中東2點位的絕對偏差最小，為7.8%，西2點位的絕對偏差最大，為13.7%，4個點位的平均絕對偏差為11.2%.

表5 基于LPI方法的不同點位ρ(PM2.5)預測值與監(jiān)測值的差異

圖7 基于LPI的空間插值Fig.7 Spatial interpolation based on the LPI method

2.4.3IDW方法

表6和圖8顯示了基于IDW方法的ρ(PM2.5)預測值和監(jiān)測值的差異，可以看出，該模型預測的濃度高值區(qū)集中在研究區(qū)域南側(cè)、東側(cè)和北側(cè)的3個點位周邊，低值區(qū)集中在研究區(qū)域的中心位置. 提取4個點位的濃度發(fā)現(xiàn)，模型對北、西、東點位的ρ(PM2.5)預測值偏高，對南側(cè)點位的預測值偏低，其中東2點位的絕對偏差最小，為9.9%，南2點位的絕對偏差最大，為18.4%，4個點位的平均絕對偏差為12.1%.

表6 基于IDW方法的不同點位ρ(PM2.5)預測值與監(jiān)測值的差異

圖8 基于IDW方法的空間插值Fig.8 Spatial interpolation based on the IDW method

2.4.4Kriging(克里根插值)

表7和圖9顯示了基于Kriging方法的ρ(PM2.5)預測值和監(jiān)測值的差異，可以看出，該模型預測的中心位置濃度最低，東側(cè)和南側(cè)濃度最高，其他區(qū)域濃度空間分布較為均勻. 提取4個點位的濃度發(fā)現(xiàn)，模型對北、西、東點位的ρ(PM2.5)預測值偏高，對南側(cè)點位的預測值偏低，其中北2點位的絕對偏差最小，為7.7%，西2點位的絕對偏差最大，為15.0%，4個點位的平均絕對偏差為9.3%.

表7 基于Kriging方法的不同點位 ρ(PM2.5)預測值與監(jiān)測值的差異

圖9 基于Kriging方法的空間插值Fig.9 Spatial interpolation based on the Kriging method

2.4.5KS方法

表8和圖10顯示了基于KS方法的ρ(PM2.5)預測值和監(jiān)測值的差異，可以看出，該模型預測濃度的高值區(qū)處于研究區(qū)域的南側(cè)、東側(cè)和北側(cè)，低值區(qū)位于中心位置偏西的區(qū)域. 提取4個點位的濃度發(fā)現(xiàn)，模型對各點位的ρ(PM2.5)預測值偏高，其中南2點位的絕對偏差最小，為4.3%，西2點位的絕對偏差最大，為13.7%，4個點位的平均絕對偏差為10.4%.

表8 基于KS方法的不同點位ρ(PM2.5)預測值與監(jiān)測值的差異

圖10 基于KS方法的空間插值Fig.10 Spatial interpolation based on the KS method

2.4.6DK方法

表9和圖11顯示了基于DK方法的ρ(PM2.5)預測值和監(jiān)測值的差異，可以看出，該模型預測的濃度空間分布呈現(xiàn)自西北部向東南部依次升高的趨勢，中間區(qū)域濃度分布較為均勻. 提取4個點位的濃度發(fā)現(xiàn)，模型對北、西、東點位的ρ(PM2.5)預測值偏高，對南側(cè)點位的預測值偏低，其中東2點位的絕對偏差最小，為5.8%，西2點位的絕對偏差最大，為21.7%，4個點位的平均絕對偏差為12.0%.

表9 基于DK方法的不同點位ρ(PM2.5)預測值與監(jiān)測值的差異

圖11 基于DK方法的空間插值Fig.11 Spatial interpolation based on the DK method

2.5 研究結(jié)果的適用性

該研究選取了京津冀地區(qū)污染較為嚴重的保定市作為研究對象，以河北大學校區(qū)為中心點位，每隔3 km在南、北、東、西方向布設采樣點，避免了顯著常年主導風向(東南風向)的影響. 試圖通過密集布點并同時采樣、分析的方法，考察大氣中主要污染物在與中心點位間不同距離下的濃度差異. 研究通過對各點位氣態(tài)污染物、顆粒物(PM2.5和PM10)以及PM2.5中主要組分的濃度差異性比較發(fā)現(xiàn)，對于SO2、NO2、顆粒物及其組分來說，距離不超過6 km點位的濃度與中心點濃度偏差均小于20%，在6 km以上時其偏差均大于20%，空氣質(zhì)量監(jiān)測點位的代表性在一般為5～6 km；對于CO、O3和VOCs來說，不同地區(qū)其濃度的空間分布更為均勻，各點位與中心點的濃度偏差均小于20%，監(jiān)測點位的代表性范圍更大. 該結(jié)果與基本科學認知相符：在冬季，SO2、CO、VOCs以及顆粒物與一次排放關聯(lián)性較大，因此在小范圍區(qū)域內(nèi)具有較好的一致性；而對于O3、NO2、PM2.5中的硝酸鹽，則因其主要來自于二次生成，在較大范圍內(nèi)具有一致性.

為了解決6 km內(nèi)無監(jiān)測站點情形下大氣污染物濃度水平的估算，該研究進行了6種空間分析方法對ρ(PM2.5)預測值的比較研究，結(jié)果表明，基于Kriging的空間分析方法對保定市ρ(PM2.5)的預測效果最好，預測值和監(jiān)測值的偏差小于10%，小于其與中心點位之間的偏差，可以用于估算監(jiān)測點位6 km以外沒有其他監(jiān)測站點區(qū)域的空氣污染物濃度. 考慮到一定的保險系數(shù)和實際應用中的使用習慣，該研究將6 km人為確定為5 km.

注：5 km以內(nèi)用監(jiān)測站點濃度代替，5 km以外使用Kriging方法插值.圖12 應用綜合方法在北京市進行大氣污染物濃度預測的結(jié)果Fig.12 The sketch map of spatial distribution of air pollutants in Beijing based on the comprehensive method

該研究結(jié)果表明，在進行某城市大氣污染物暴露水平研究時，如果人群聚集范圍在監(jiān)測點位周邊5 km 以內(nèi)的區(qū)域時，可直接使用該站點的監(jiān)測數(shù)據(jù)，監(jiān)測點位周邊5 km以外的區(qū)域可使用基于Kriging的空間分析方法獲得，圖12為以該方法在北京市進行應用的示意圖. 另外，該研究是基于冬季相對污染較為嚴重的北方城市的實測和模型模擬，從污染物排放和大氣化學反應的規(guī)律考慮，研究結(jié)果可能更適合于冬季的北方城市. 今后將在更多城市和不同季節(jié)進行相關研究.

3 結(jié)論

a) 通過在保定市的現(xiàn)場監(jiān)測發(fā)現(xiàn)，對于SO2、NO2、顆粒物及其主要組分，空氣質(zhì)量監(jiān)測點位的代表性一般為5 km左右；對于CO、O3和VOCs來說，不同地區(qū)其濃度的空間分布更為均勻，監(jiān)測點位的代表性范圍更大.

b) 基于Kriging的空間分析方法對ρ(PM2.5)的預測效果最好，預測值和監(jiān)測值偏差小于10%，可以用于獲取監(jiān)測點位5 km以外沒有站點區(qū)域的空氣污染物濃度.

c) 在進行某點位污染物濃度預測時，若該點位周邊5 km以內(nèi)有空氣質(zhì)量監(jiān)測點位，則可用該點位的監(jiān)測濃度代替；若5 km以內(nèi)沒有空氣質(zhì)量監(jiān)測點位，則可基于該市監(jiān)測點位的污染物濃度進行Kriging空間插值，從而獲得該點位的濃度.

d) 該研究結(jié)果基本可適用于所有城市的各季節(jié)，但更適合于北方的冬季.