基于深度學習方法的PM2.5精細化時空估算模型

耿 冰,孫義博,曾巧林,商豪律,劉霄宇,單菁菁* (.中國社會科學院生態(tài)文明研究所,北京 0070；.中國環(huán)境科學研究院生態(tài)研究所,北京 000；.重慶郵電大學計算機科學與技術學院,重慶 0006；.中國科學院空天信息研究院數(shù)字地球重點實驗室,北京 0009；.中國信息通信研究院,北京 009)

為了應對日益嚴重的空氣污染問題,我國大規(guī)模建立地面PM2.5監(jiān)測站點,對重污染天氣進行監(jiān)測及預警[1-3].眾多學者也在此基礎上開發(fā)了 PM2.5濃度估算模型[4-8].但是,地面 PM2.5監(jiān)測站點僅能提供空間上“點”尺度的觀測,有限的空間范圍不足以代表PM2.5在空間上的異質(zhì)性, 而PM2.5的空間信息對于研究空氣污染與經(jīng)濟、地理及人口之間的關系至關重要.

近年來,隨著衛(wèi)星遙感技術的發(fā)展,使得區(qū)域尺度的污染信息獲取成為可能,利用衛(wèi)星反演的大氣氣溶膠光學厚度(AOD)估算大空間尺度的地表PM2.5濃度也已被廣泛采用[9-13].

目前,已有多種衛(wèi)星遙感數(shù)據(jù)反演的AOD產(chǎn)品被用于估算地表PM2.5濃度的時空分布[14-17].

與此同時,多種類型的統(tǒng)計模型被提出并應用于PM2.5濃度的估算中[18-22],這些模型均以AOD作為主要指示因子,結合氣象觀測及其他類型的統(tǒng)計參數(shù)估算地面 PM2.5的時空分布.例如,早期的研究采用一元線性回歸模型僅采用AOD作為指示因子來估算PM2.5濃度[23];或更為復雜的采用多元或廣義線性回歸模型,考慮更多的地表及氣象參數(shù)用以提高 PM2.5的估算精度[24-25].但是在真實的環(huán)境中,PM2.5濃度的分布是一個與多種因素有關的非線性過程,在時間和空間上存在著強烈的可變性,因此學者們開發(fā)了更為復雜的模型對PM2.5濃度與AOD關系的時空變異性進行描述,例如地理(和時間)加權回歸模型[10]、混合效應模型[26]以及廣義加權混合模型等[27].然而,本質(zhì)上這些統(tǒng)計模型仍然是線性的,模型內(nèi)部簡化了PM2.5與AOD及其他指示因子之間的復雜關系,使PM2.5濃度估算結果仍然存在較大的不確定性.隨著計算機技術的發(fā)展,機器學習(包括深度學習)方法以其強大的非線性建模能力越來越多地被用于 PM2.5濃度的估算中[28],例如支持向量回歸模型[29]、隨機森林模型[30]、人工神經(jīng)網(wǎng)絡模型[24]、貝葉斯方法[31]、廣義回歸神經(jīng)網(wǎng)絡模型[32]以及深度信念網(wǎng)絡[33]等,這些模型在對PM2.5濃度的估算方面均表現(xiàn)出比傳統(tǒng)統(tǒng)計模型更好的性能.在指示因子的選擇方面,這些機器學習模型除了采用AOD和常規(guī)的氣象觀測參數(shù)之外,還使用了包括相鄰時間和空間上觀測的 PM2.5信息、土地利用信息、植被指數(shù)信息、NO2濃度信息、人口密度、海拔高程[26,33]以及路網(wǎng)密度信息等,這些信息或多或少與PM2.5濃度分布相關.考慮的影響因子越多,越能夠提高PM2.5估算精度.但是,過多的人工設計的特征不僅耗時耗力,而且過于復雜的特征選擇也不利于模型的工程化實施.此外,目前的模型大多僅對日平均 PM2.5濃度進行估算,且空間分辨率相對粗糙(大于 3km).盡管此類模型可以有效降低目標函數(shù)的復雜性,但卻忽略了PM2.5濃度每小時的時空變異性.針對以上問題,為了有效地開展 PM2.5濃度精細化時空尺度(即每小時和1km的時空分辨率)估算,需要一種非線性表達能力更強并且容易實現(xiàn)工程化的模型.

深度學習方法[34]作為當前最先進的機器學習技術之一,以其優(yōu)異的非線性表達能力在許多領域都取得了超過傳統(tǒng)機器學習方法的顯著成果.目前,已有研究人員采用深度學習方法來估算 PM2.5濃度的時空分布[33,35-36],但是相關的模型規(guī)模仍然相對較小,并且很大程度上依舊依賴于人工特征選擇,并沒有充分利用深度學習方法通過更深更寬的網(wǎng)絡結構來表達高度復雜目標函數(shù)的優(yōu)點.因此,本文以北京市 2017年的 PM2.5觀測數(shù)據(jù)為基礎,提出了一種典型的深度學習模型(PM2.5-DNN),僅采用衛(wèi)星遙感反演的AOD數(shù)據(jù)以及常規(guī)的氣象觀測要素(例如氣溫、地表溫度、風速、風向、相對濕度、壓強以及能見度)來估算PM2.5濃度的時空分布.

1 數(shù)據(jù)與模型

1.1 研究區(qū)與數(shù)據(jù)來源

1.1.1 研究區(qū)概況 北京市位于華北平原,中心經(jīng)緯度為 116.41°E,39.92°N,總面積約 16410.54km2,2019年常住人口約2153.6萬人.其地處暖溫帶半濕潤半干旱季風區(qū),氣候四季分明,夏季炎熱多雨,冬季寒冷干旱,春秋短促.本文采用了北京市2017年全年觀測的每小時PM2.5濃度數(shù)據(jù)和氣象觀測數(shù)據(jù).研究區(qū)域地理空間范圍以及 PM2.5站點和氣象觀測站點的分布見圖1所示.

圖1 研究區(qū)及PM2.5監(jiān)測站和氣象觀測站點分布Fig.1 Distribution of study area and PM2.5 monitoring stations and meteorological observation stations

1.1.2 衛(wèi)星AOD數(shù)據(jù) 本文衛(wèi)星遙感數(shù)據(jù)采用葵花-8氣象衛(wèi)星數(shù)據(jù).葵花-8屬于第3代地球靜止氣象衛(wèi)星,其觀測范圍為東西 80°E～160°W,南北 60°N～60°S,距離地面高度 35800km,星下點位于 140.7°E[17].衛(wèi)星搭載的主要傳感器為AHI(高像素葵花成像儀),可見光最高分辨率為 0.5km,紅外-近紅外最高空間分辨率為 1km,最高時間分辨率為 10min,是目前全球最先進的氣象觀測傳感器之一.

葵花-8號衛(wèi)星 AOD數(shù)據(jù)產(chǎn)品的反演采用了Yang等提出的新暗目標算法(New-DT),該算法主要利用葵花-8衛(wèi)星的可見光和近紅外波段數(shù)據(jù)反演空間分辨率為1km的每小時AOD數(shù)據(jù)產(chǎn)品.

1.1.3 環(huán)境及氣象觀測數(shù)據(jù) 2017年PM2.5監(jiān)測站點數(shù)據(jù)來自國家環(huán)境監(jiān)測中心網(wǎng)站,共采用了北京市12個PM2.5監(jiān)測站點.2017年北京市氣象觀測數(shù)據(jù)來自于中國氣象數(shù)據(jù)中心網(wǎng)站,本文所用到的主要氣象觀測數(shù)據(jù)包括氣溫(A_temp)、地表溫度(S_temp)、風速(wing_S)、風速方向(wind_D)、相對濕度(RH)、地表壓強(SP)和能見度(VIS).為了獲得與PM2.5濃度相關性最佳的氣象觀測數(shù)據(jù),本文選用了在 PM2.5監(jiān)測站半徑為 5km區(qū)域內(nèi)的氣象觀測站,最后共選擇了8個氣象觀測站點.PM2.5監(jiān)測站點與氣象觀測站點的分布如圖1所示.

1.1.4 數(shù)據(jù)預處理 由于模型的構建需要將 PM2.5濃度數(shù)據(jù)、AOD數(shù)據(jù)以及氣象觀測數(shù)據(jù)一一對應,因此需要對所有數(shù)據(jù)進行預處理,使之形成時間和空間上一致的數(shù)據(jù)集,用于對所構建的PM2.5濃度估算模型進行訓練與驗證.其中 AOD數(shù)據(jù)提取 PM2.5監(jiān)測站點位置所在柵格的數(shù)值.對于氣象觀測數(shù)據(jù),由于 PM2.5監(jiān)測站點的分布與氣象觀測站點的分布不同,本文使用了PM2.5周圍5km范圍內(nèi)的氣象觀測站點觀測均值作為與之匹配的氣象觀測數(shù)據(jù).此外,由于較小的太陽高度角會導致大氣路徑變長,使得大氣漫反射的比例變大,從而導致衛(wèi)星遙感反演的AOD數(shù)據(jù)產(chǎn)品的精度變差.因此,為了有效避免太陽高度角對AOD數(shù)據(jù)反演精度的影響,本文僅使用北京時間9:00～16:00之間8h的葵花8號衛(wèi)星影像反演AOD數(shù)據(jù),并同時獲取與AOD數(shù)據(jù)時間上一致的 PM2.5濃度與氣象觀測數(shù)據(jù).最后,去掉數(shù)據(jù)中的缺失值以及由于降雨影響導致的無效數(shù)據(jù),并將數(shù)據(jù)進行歸一化處理(零均值,單位方差).經(jīng)過預處理后,共得到可用數(shù)據(jù)17059條.

1.2 模型構建方法

1.2.1 深度網(wǎng)絡模型(PM2.5-DNN) 深度學習屬于機器學習中的一類模型,它通過構建非常深的神經(jīng)網(wǎng)絡來學習輸入數(shù)據(jù)的多級表示特征,從而表達復雜及抽象的概念或模式[34].與傳統(tǒng)神經(jīng)網(wǎng)絡模型相比,深度學習模型通常由超過三層的隱含層構成.并且可以在沒有進行預訓練的情況下直接對深度網(wǎng)絡進行端到端的監(jiān)督訓練.

本文采用AOD數(shù)據(jù)和常規(guī)的氣象觀測參數(shù)作為輸入變量,通過深度學習模型直接構建輸入變量與 PM2.5濃度之間的高時空分辨率(每小時,1km)關系模型.模型的輸入?yún)?shù)為:葵花 8衛(wèi)星反演的AOD、空氣溫度(A_temp)、地表溫度(S_temp)、風速(wind_S)、風向(wind_D)、相對濕度(RH)、壓強(SP)、可見度(VIS)以及年積日(DOY),模型的輸出為PM2.5濃度的估算值.對模型隱含層神經(jīng)單元數(shù)據(jù)的確認采用啟發(fā)式搜索方法,將隱含層的神經(jīng)單元個數(shù)設為 10,并以 10為步長進行迭代的訓練和驗證,并統(tǒng)計驗證誤差,直至驗證誤差穩(wěn)定且不再降低,然后以相同的方式來確定多個隱含層的神經(jīng)元個數(shù),直至整個模型的驗證誤差不再降低為止.模型的訓練采用 ReLU作為激活函數(shù),輸出層采用線性函數(shù)作為激活函數(shù),最終確定的深度學習網(wǎng)絡模型結構如圖2所示,結構為9-300-300-100-20-1.

圖2 本文所采用的深度網(wǎng)絡模型結構Fig.2 The structure of the deep network model

模型的訓練和應用流程如圖 3所示,在訓練階段通過采用誤差反向傳播方法對輸入的訓練數(shù)據(jù)集進行學習,獲得能夠表征輸入數(shù)據(jù)時空特征的模型內(nèi)部參數(shù),并對模型的訓練結果進行驗證,在獲得可靠的驗證精度后即可應用訓練好的模型,通過輸入 AOD影像及柵格化后的氣象觀測數(shù)據(jù),對 PM2.5濃度的時空分布進行估算.

圖3 模型訓練和應用流程圖Fig.3 Model training and application flow chart

1.2.2 線性混合效應模型(LME) 線性混合效應模型(LME)是目前估算 PM2.5時空分布的主要方法之一,是既包含了固定效應參數(shù)又包含隨機效應的統(tǒng)計回歸模型.其中固定效應表示模擬因子對 PM2.5的多年平均影響狀態(tài),而隨機效應則用于解釋 PM2.5與AOD以及其他氣象因子之間的日變化關系,以隨機截距或者隨機系數(shù)的形式表示.線性混合效應模型可以表達為:

式中:PM2.5,it表示第i個監(jiān)測站點在時間t的 PM2.5濃度值;bo和bo,t分別代表模型的固定截距和隨機截距;b1～b9和b1,t～b8,t分別代表各自變量參數(shù)的固定效應斜率和各變量的隨機效應斜率;εit第i個監(jiān)測站點在時間t的隨機誤差項,b0,t、b1,t和b2,t為其參數(shù);Ψ為隨時間變化的隨機效應方差-協(xié)方差矩陣.

1.2.3 地理加權回歸模型(GWR) 地理加權回歸模型(GWR)區(qū)別于傳統(tǒng)回歸方法之處在于不同輔助變量的回歸系數(shù)不再是利用全局信息所獲得的常量進行估算,而是用鄰近觀測值進行局部加權回歸從而得到相關系數(shù),并考慮了數(shù)據(jù)的空間位置.其公式可表達為:

式中:PM2.5,it表示第i個監(jiān)測站在時間t的PM2.5濃度值,β0,i為模型的固定截距;β1,I～β9,i分別代表各自變量參數(shù)的回歸系數(shù),εi為第i個監(jiān)測站點的回歸殘差.

1.2.4 支持向量回歸模型(SVR) 支持向量回歸(SVR)模型是支持向量機模型在回歸問題上的應用模型,該模型會在盡可能擬合現(xiàn)有數(shù)據(jù)的情況下考慮模型的泛化誤差,從而盡量避免模型的過擬合,使得學習到的模型能夠在未知數(shù)據(jù)上具有良好的預測性能.SVR模型在PM2.5模擬方面已有了較多的應用,以下對其原理進行簡單敘述.

給定數(shù)據(jù)集{(xi,yi),i,…,m},x為輸入因子,y為輸出因子,SVR模型可以表示為:

式中:w是權重向量;φ(x)是將輸入數(shù)據(jù)從輸入空間映射到特征空間的核函數(shù);b是常數(shù)項.對模型的訓練是使風險函數(shù)最小,風險函數(shù)可表示為:

式中:等號右側第一項是對模型復雜程度的懲罰項;第二項是模型輸出值與真實值之間誤差的懲罰項;常數(shù) C是用來調(diào)整懲罰比例的系數(shù).本文使用徑向基函數(shù)(RBF)為核函數(shù),表示為:

在對SVR模型的訓練中,采用了格網(wǎng)搜索方法來確定模型超參數(shù)(σ,C),即設定σ和C的范圍在4到-4之間,以 0.8為間隔進行遍歷計算,找到模型驗證誤差最小的超參數(shù)組合即為最優(yōu)參數(shù).

1.2.5 隨機森林回歸模型(RFR) 隨機森林(RF)算法是通過集成學習的思想將多顆決策樹集成的一種分類與回歸算法.隨機森林引入集成學習思想和隨機子空間思想,通過實現(xiàn)樣本選取隨機性和特征選取隨機性,對樣本單獨構建決策子樹,結合集成學習思路將各決策子樹的結果按照一定規(guī)則匯總作為最后輸出.對于回歸模型而言,匯總規(guī)則為取平均值.

對一組由決策子樹{h(x,θt),t=1,2,…,T}構成的決策組合模型.其中θt為服從獨立同分布的隨機變量,x為自變量,T為決策子樹的個數(shù).回歸模型的估算結果為:

2 結果與討論

2.1 深度網(wǎng)絡模型PM2.5-DNN的訓練和驗證

將數(shù)據(jù)集隨機分隔為訓練數(shù)據(jù)集(80%)和測試數(shù)據(jù)集(20%)分別應用于模型的訓練階段和測試階段.對模型的訓練采用帶動量項的隨機梯度下降算法,并采用從均勻分布中采樣的方式對 PM2.5-DNN模型的初始權重進行初始化.此外,超參數(shù)的設置也會顯著的影響模型訓練結果,超參數(shù)主要包括學習速率,動量以及為防止過度擬合而引入的dropout和正則化因子[37].本文對超參數(shù)的設置采用了格網(wǎng)搜索的方法,即設置各超參數(shù)取值范圍,然后以一定的步長進行迭代訓練,找出獲得最佳驗證性能的數(shù)值作為后續(xù)模型訓練時采用的超參數(shù).本文最終確定的模型訓練超參數(shù)設置情況為:批量大小為 100,學習速率為0.04,動量為0.8,dropout為0.1,L2正則化因子為 0.0002,同時采用早停法來避免模型的過度擬合.對模型的訓練采用10折交叉驗證方法,即將訓練數(shù)據(jù)集隨機平均分為10等份,輪流將其中9份作為訓練數(shù)據(jù),1份作為驗證數(shù)據(jù)進行模型訓練,最后將10次訓練的驗證結果的均值作為模型的性能[37].

采用相關性系數(shù)(R2)、均方根誤差(RMSE)、平均預測誤差(MPE)以及相對預測誤差(RPF)作為模型性能的評價指標.

模型訓練結構如圖 4所示,在模型的訓練階段,模型的總體性能分別為:R2=0.95,RMSE=10.6μg/m3,MPE=1.64μg/m3,RPE=25.47%;在模型的測試階段,模型的總體性為:R2=0.88,RMSE=18.78μg/m3,MPE=0.73μg/m3,RPE=42.98%.通過對比訓練階段和測試階段的 R2可以看,出模型在訓練階段存在一定的過擬合現(xiàn)象,導致了所構建的PM2.5-DNN模型在測試階段對高濃度范圍PM2.5的低估和低濃度范圍PM2.5的高估.但同時測試結果的相關性系數(shù)和誤差水平仍然表現(xiàn)的非常穩(wěn)健,表明本文所構建的 PM2.5-DNN模型對PM2.5濃度的估算仍然非常有效.

圖4 PM2.5-DNN模型的總體性能評價Fig.4 Performance evaluation of PM2.5-DNN model

2.2 深度網(wǎng)絡模型與其他模型的對比

從表 1可見,機器學習方法(包括 SVR,RFR和PM2.5-DNN模型)表現(xiàn)出了比傳統(tǒng)方法(LME和GWR)更好的估算性能.PM2.5-DNN與SVR和RFR模型相比表現(xiàn)出了更加優(yōu)異的估算性能,這主要得益于PM2.5-DNN模型的深層結構,能夠從大量數(shù)據(jù)中直接學習高度復雜的函數(shù)關系.同時,本結果也表明,深度網(wǎng)絡模型在PM2.5濃度的估算方面具有極高的應用潛力,可以直接通過端到端的方式進行訓練,即直接使用容易獲得的觀測因子進行建模便能夠獲得最佳的估算性能.

表1 不同PM2.5估算模型之間的結果對比Table 1 Results comparison of different PM2.5 estimation models

2.3 北京市PM2.5濃度的精細化時空分布估算

2.3.1 每小時PM2.5濃度估算分析 為獲得能夠覆蓋整個北京市區(qū)域的氣象觀測插值數(shù)據(jù),本文采用了覆蓋整個京津冀區(qū)域內(nèi)的氣象觀測站點并采用基于反距離加權方法(IDW)進行插值獲得北京市氣象柵格數(shù)據(jù).然后將訓練好的PM2.5-DNN模型應用于衛(wèi)星AOD柵格數(shù)據(jù)及插值生成的氣象柵格數(shù)據(jù)生成2017年北京市每小時PM2.5濃度柵格圖.圖5中,最大PM2.5濃度分別為:89.92、82.37、91.42、89.89、89.73、83.27、80.6和91.62μg/m3.結合與地面PM2.5監(jiān)測站點的數(shù)據(jù)相對比, PM2.5-DNN模型估算的PM2.5濃度時空分布與地面實測結果均吻合較好,表明本文所構建的PM2.5-DNN模型具有較好的時空估算性能.

圖5 每小時PM2.5濃度模擬結果與監(jiān)測站點的對比Fig.5 Hourly PM2.5 concentration estimated by PM2.5-DNN model compared with the monitoring site

2.3.2 各月份PM2.5濃度分布估算分析 從各月份PM2.5濃度分布來看,北京市2017年1～12月份最大PM2.5濃度分別為:51.56、39.05、34.05、34、57.48、29.23、18.49、15.09、32.67、27.41、46.77 和39.47μg/m3.在PM2.5濃度時空分布中可知,冬季節(jié)污染最為嚴重,夏季空氣質(zhì)量最優(yōu).研究結果也清晰地展示了2017年5月發(fā)生的北京市PM2.5重度污染過程(圖6).

圖6 PM2.5-DNN模型生成的2017年各月份PM2.5濃度分布Fig.6 PM2.5 concentration distribution per month in 2017 generated by PM2.5-DNN model

2.3.3 各行政區(qū)PM2.5濃度分布分析 北京市共16個行政區(qū),但僅有12個觀測站點,無法實現(xiàn)觀測站點的行政區(qū)全覆蓋,利用PM2.5-DNN模型可生成更加精細的 PM2.5濃度時空分布數(shù)據(jù).從計算結果可知(表 2),2017年北京市 PM2.5濃度最大值出現(xiàn)在房山和豐臺,其次為昌平、海淀、石景山,延慶PM2.5濃度最低.城六區(qū)(東城、西城、朝陽、海淀、豐臺、石景山)年平均 PM2.5濃度均高于其他地區(qū),生態(tài)涵養(yǎng)區(qū)(門頭溝、平谷、懷柔、密云、延慶、昌平和房山的山區(qū))PM2.5濃度較低.

表2 PM2.5-DNN模型生成的2017年各行政區(qū)PM2.5濃度分析Table 2 Analysis of PM2.5 concentration in each administrative district in 2017 generated by PM2.5-DNN model

同時,本文將所用12個觀測站點的年度平均值與其所在柵格單元的 PM2.5-DNN模型模擬結果的年度平均進行了對比,其中 R2=0.89,RMSE=20.04μg/m3,與模型的測試結果一致,進一步論證了所構建的PM2.5-DNN模型的可靠性(圖7).

圖7 各站點年平均觀測值與估算值對比Fig.7 Comparison of annual average observation values and estimated values of each station

2.4 討論

與 LME、GWR、SVR以及RFR相比,本文提出的基于深度網(wǎng)絡模型的 PM2.5-DNN模型僅采用容易獲得的觀測因子就可以獲得最佳的估算性能.結合衛(wèi)星遙感反演的AOD數(shù)據(jù),采用PM2.5-DNN模型實現(xiàn)了PM2.5濃度1km逐時的時空精細化模擬.

在模型應用方面,由于模型的建立主要是基于2017年北京市內(nèi)的12個PM2.5監(jiān)測站點、8個氣象觀測站點以及葵花-8號衛(wèi)星的AOD數(shù)據(jù)產(chǎn)品,考慮到其他年份或區(qū)域會出現(xiàn)不同的PM2.5時空模式,模型應用的最佳方式是進行內(nèi)插,因此本文僅對 2017年北京市PM2.5濃度的時空分布進行了模型估算.對于深度網(wǎng)絡模型而且,隨著模型的深度(隱含層數(shù))和廣度(每層神經(jīng)單元個數(shù))的增加以及模型構建時輸入數(shù)據(jù)量的增大,所訓練的深度網(wǎng)絡模型就越能更好地對復雜 PM2.5時空變化模式進行模擬.因此,采用更多年份及更多觀測站點的數(shù)據(jù)進行深度網(wǎng)絡模型的構建,實現(xiàn)全國區(qū)域多年份的PM2.5濃度時空精細化模擬是本文下一步的研究重點.此外,由于云層覆蓋造成的AOD數(shù)據(jù)在空間上的數(shù)據(jù)缺失,造成了PM2.5濃度時空分布數(shù)據(jù)在部分空間上的不連續(xù),這也是本文存在的不足之一.

3 結論

3.1 提出一種基于深度學習方法的地面PM2.5濃度時空估算模型 PM2.5-DNN,該模型僅需要常規(guī)的氣象觀測數(shù)據(jù)(包括氣溫,地表溫度、風速、風向、相對濕度、壓強以及能見度)結合衛(wèi)星遙感反演的AOD數(shù)據(jù),就可以對地表 PM2.5濃度進行高時空分辨率的估算.

3.2 深度學習模型在構建復雜關系模型中具有較強的性能,與線性混合效應模型、地理加權回歸模型、支持向量回歸模型、隨機森林回歸模型相比,PM2.5-DNN模型表現(xiàn)出更高精度的估算性能,其測試結果R2可以達到0.88.