劉旭林,趙文芳,唐 偉
1(北京市氣象探測中心,北京 100176)
2(北京城市氣象研究院,北京 100089)
3(北京市氣象信息中心,北京 100089)
4(中國氣象局發展研究中心,北京 100081)
E-mail:yoyozwf@sina.cn
PM2.5(fine particulate matter)是指一種懸浮于大氣中的空氣動力學直徑小于等于2.5μm的細顆粒物[1-3],它是構成霾的主要成分,其濃度越高,意味著霧霾污染越嚴重.霧霾天氣容易引起交通事故,引發多種呼吸道疾病,危害人類健康[4-7].近幾年,我國多地出現霧霾天氣,受到各級政府、部門及社會的廣泛關注.衡量霧霾污染程度的首要因子是PM2.5濃度,因此準確預測PM2.5濃度對大氣污染防御、空氣質量的監測和政府決策都具有重要意義.
目前PM2.5的預測主要有數值模式法和統計學預報法.數值模式法廣泛應用于氣象預報領域,主要基于空氣動力學理論和物理化學過程,使用數學方法建立大氣污染濃度的稀釋和擴散模型動態預測空氣質量和主要污染物的濃度變化.國內外常用的空氣質量數值模式包括美國環保局研發的多尺度空氣質量模式(The Community Multi-scale Air Quality,簡稱CMAQ)[8-11],美國NOAA預報系統實驗室研發的WRF-Chem(Weather Research Forecast-Chemical)[12,13],中國氣象局和中國氣象科學研究院研發的大氣化學模式(Chinese Unified Atmospheric Chemistry Environment,簡稱CUACE)[14]以及北京市氣象局研發的區域環境氣象數值預報系統(Beijing Regional Environmental Meteorology Pre-diction System,簡稱BREMPS)[15]等.這些模式通常依賴大量的計算,很多參數都是根據經驗估計的,很多條件也是假定為理想狀態,因此存在一定的局限性和不確定性.此外,由于這些模式對空氣質量數據和氣象數據的要求較高,相關數據往往很難獲取,所以數值模式方法在國內大多城市并不成熟[16].
統計學預報法是利用統計學方法建立模型開展天氣預報,常用的方法有多元線性回歸、隨機森林、貝葉斯、支持向量機(Support Vector Machine,簡稱SVM)、神經網絡等.基于這些算法建立的模型通常使用遺傳算法、粒子群算法和蟻群算法等優化算法獲取最優參數.已有大量學者利用空氣質量數據、氣象觀測數據和數值模式預報數據,應用一種或多種統計學方法建立預測模型,對PM2.5濃度和其他污染物濃度進行預報[17-23].劉杰等提出應用SVM和模糊粒化時間序列相結合預測PM2.5的方法;李龍等選擇綜合氣象指數、二氧化硫濃度、一氧化硫濃度、二氧化氮濃度和PM10濃度構成特征向量,利用特征向量和PM2.5濃度數據應用最小二乘支持向量機(Least Squares Support Vector Machine,簡稱LS-SVM)預測模型;戴李杰等聯合應用支持向量機(SVM)和粒子群優化(Particle Swarm Optimization,簡稱PSO)算法建立滾動預報模型,對PM2.5未來24小時濃度進行預報,同時對未來一天的晝、夜均值及日均值濃度進行預報;張恒德等將時間學列和卡爾曼濾波結合起來用在霧霾預報技術中.然而,除了氣象條件,污染物濃度還受排放量、交通條件、人口密度等因素的共同影響,使用單一統計方法很難建立準確度高的預報模型.
深度學習是人工智能領域中一種新穎的機器學習方法,可以對大量輸入數據的特征表示進行有效學習,為氣象時間序列的預測提供了新的研究思路和方法.深度學習的主要神經網絡模型主要有卷積神經網絡(Convolutional Neural Networks,CNN)、遞歸神經網絡(Recurrent Neural Network,RNN)、長短時記憶(Long Short Term Memory,LSTM)、對抗神經網絡等,已有一些學者利用這些模型開展氣象預測的研究.Vidushi Chaudhary 等提出了一個多層的LSTM模型預測未來空氣污染物的濃度[24];ZHAO等構建了LSTM-FC預測模型,使用歷史空氣質量數據、氣象數據和天氣預報數據預測未來48小時的PM2.5濃度[25].此外,還有一些學者將LSTM模型與特征空間相關性相結合應用于PM2.5濃度預報,如:Congcong Wen等提出了一種新的時空卷積長短期記憶神經網絡,將當前站和近鄰站點PM2.5濃度數據經過1維卷積運算后輸入到模型中[26];SOH等提出LSTM預測模型,使用當前站點和近鄰站的PM2.5濃度數據、氣象數據和地形數據,其中地形數據用以提取地形對空氣質量的影響[27].這些模型各有自的特點和適用場景,但都是根據輸入數據對未來某個時間的PM2.5濃度進行預測,而氣象預測和服務都需要根據輸入的歷史序列數據預測未來多個時間的PM2.5濃度.Seq2Seq是一種根據輸入序列生成輸出序列的模型,輸入序列的長度和輸出序列的長度可以不同,主要應用于機器翻譯和文本生成等應用.已有文獻利用Seq2Seq模型應用在語音識別、文本摘要、對話系統、圖像標題生成中,但是將Seq2Seq模型應用至PM2.5濃度預測的研究相對甚少.因此,本文選擇Seq2Seq進行建模,實現對未來多個時間段的PM2.5濃度預測.
綜上所述,基于PM2.5濃度預測的重要性和難度,考慮PM2.5濃度與氣象要素的時空相關性,本文將采用Seq2Seq建立PM2.5濃度預測模型,以歷史的空氣質量觀測數據、氣象觀測數據和PM2.5濃度數據作為輸入,先進行卷積操作提取出空間特征,再將結果輸入Seq2Seq模型,進而得到最終的預測結果.
為了研究北京地區霧霾天氣的季節變化特征,收集了北京地區2008年-2018年所有國家級氣象站的觀測數據,計算每個站逐年霧霾天氣出現的總次數,統計每個季節所有站的季平均霧霾天氣出現次數,結果如圖1所示.可以看出,近10年在北京地區霧霾天氣出現的平均次數從春節到冬季呈現先下降后上升的趨勢,冬季發生霧霾天氣次數最高,夏季最低.從2008年到2014年,冬季出現霧霾天氣現象的次數逐年遞增,2014年之后則逐年下降,從側面反映了近幾年空氣污染治理取得了良好效果.春季和秋季出現霧霾天氣的次數比較接近,逐年變化也不大.

圖1 20個國家級氣象站逐年霾日數的季節變化趨勢
對所有氣象站的PM2.5濃度數據按春、夏、秋、冬季分類,然后在每個季節內再按照觀測時間0到23點進行分類排序,計算所有站點各季節0點到23點逐小時PM2.5濃度均值,結果如圖2所示.可以看出,春季和冬季,PM2.5濃度從凌晨到半夜呈現先下降再上升趨勢,波峰集中出現在0點和23點,波谷出現在上午9時左右.冬季PM2.5濃度均值在各時間段都保持在65μg/m3,最高突破90μg/m3,春季PM2.5濃度均值在各時間起伏較大,傍晚最低值不到60μg/m3,而0點峰值達到86μg/m3.夏秋兩季,PM2.5濃度隨時間的變化均不大,夏季峰值出現在上午9時左右,秋季峰值出現在0點和22點,濃度先下降,在上午8時左右達到一個峰值,隨后下降傍晚最低,再呈現上升趨勢到22時到達最大.

圖2 不同季節的PM2.5濃度日內逐小時變化趨勢
對所有氣象站選出與之最臨近的12個站點,按距離由近及遠排序,逐一計算每個站點和最臨近12個站點之間的PM2.5濃度相關性,結果如圖3(a)所示,其中,橫坐標為最臨近的12個站點,縱坐標為PM2.5濃度相關性.可以看出,所有站點的PM2.5濃度相關性隨距離增加呈下降趨勢,大部分的站點與最臨近6個站點的PM2.5濃度相關性大于0.5,也有的站點與距離最近站點的PM2.5濃度相關性小于0.5.計算站點PM2.5濃度與不同氣象要素的相關性,結果如圖3(b)所示.可以看出,溫度、最高溫度、平均相對濕度、能見度與PM2.5濃度呈現正相關;2米風、10米風、極大風速、小時降水量與PM2.5濃度呈負相關;PM10的濃度與PM2.5濃度相關性較強,而O3,SO2濃度與PM2.5濃度幾乎不相關.

圖3 PM2.5濃度相關性分析
從以上的數據分析可以看出,北京霧霾天氣多發生在冬季和春季,PM2.5濃度小時峰值最容易出現在這兩個季節的凌晨和上午時間段.由于觀測站點分布的不均勻,PM2.5濃度不僅與站點之間距離存在相關性,還與不同氣象要素之間存在相關性,因此在建立預測模型時候必須充分考慮PM2.5濃度的時間相關性和空間相關性.

Seq2Seq屬于encoder-decoder結構的一種,基本思想就是利用兩個RNN,一個RNN作為encoder,另一個RNN作為decoder,實現從一個序列到另外一個序列的轉換.encoder負責將輸入序列壓縮成指定長度的向量,這個過程稱為編碼,decoder負責將encoder生成的固定向量再轉化成輸出序列,這個過程稱為解碼.
PM2.5濃度預測模型包含獲取空間特征的CNN與Seq2Seq模型,使用的數據包括北京市氣象局10個觀測站逐小時的空氣污染物濃度觀測資料、逐小時氣象要素觀測資料和北京地區網格化三維氣象要素客觀分析資料.本文選擇常用于序列模型的一維卷積來獲取與站點高度相關的短期空間特征.以站點為中心,從整個三維氣象要素客觀分析網格中取出10*10的子網格作為輸入,使用一維CNN分別獲取溫度、相對濕度、風要素的空間特征.Seq2Seq模型包含兩個RNN,一個RNN作為encoder,另一個RNN作為decoder.站點在某個時刻的觀測數據、最臨近6個站的PM2.5濃度數據經過預處理后與該站在這個時刻的空間特征作為當前時刻的特征向量輸入到Encoder結構,經過運算得到輸出和隱含狀態;Encoder的輸出和隱含狀態與站點的目標預測值作為當前時刻的Decoder的輸入,經過RNN運算得到預測結果.整個模型的結構如圖4所示.

圖4 PM2.5濃度預測模型架構
Chung等的研究成果指出,門控循環單元(Gated Recurrent Neural,簡稱GRU)在較小的數據集上比LSTM表現更突出,能更好解決RNN在長序列訓練中爆炸或消失梯度的問題[28],因此,本文選擇GRU作為Seq2Seq模型的RNN.另外,根據Geman等的研究,改變隱藏節點的數量可能會減少過度擬合和增加模型泛化.一個可能的經驗法則是采取隱藏節點的數量為輸入層和輸出層維度總和的2/3[29,30].按照這個法則,如果使用過去72小時數據預測未來24小時PM2.5濃度,隱藏節點數量應該設置為(72+24)*2/3=64,這個值也是大量文獻和研究中推薦使用的.與單層GRU相比,堆疊GRU可以增加模型的學習能力,但是網絡參數也會隨之增加,這對模型泛化能力和訓練時間都有直接影響.考慮到樣本數總量,本文使用2層GRU結構,每層有64個隱含節點.此外,為了防止過擬合,在Decoder結構中使用比例為0.15的Dropout.
模型優化算法選用Adam算法.為了評估模型的性能,使用均方根誤差(Root Mean Square Error,RMSE)、平均絕對誤差(Mean Absolute Error,MAE)和和平均絕對百分比錯誤(MAPE)作為評價指標.RMSE和MAE用于評估絕對誤差,而MAPE用于測量相對誤差.RMSE和MAE反映了預測的極值效應和誤差范圍值,MAPE反映了平均值的特異性預測值.
PM2.5預測模型算法的具體步驟如下:
1)提取站點的空間特征.以站點為中心,將原始的氣象要素客觀分析資料處理為10×10的子網格,按溫度、相對濕度、風要素網格分別存儲.使用一維卷積操作這些子網格作提取不同要素的空間特征,作為站點在該觀測時刻的空間特征.
2)對站點的空氣污染物濃度和氣象要素觀測數據進行預處理.首先,計算每類數據的均值和方差,將均值在縮小而偏差在增大的數據從樣本集中剔除.其次,缺測值使用均值替代,再進行歸一化處理.最后,將它們與同一觀測時刻的空間特征進行組合,作為站點在該觀測時刻的特征向量.
3)數據融合與格式化,設置滑動時間窗口,生成用于訓練和測試模型的時間序列數據,完成原始數據到張量的轉變.
4)訓練CNN-Seq2Seq模型,得到最優超參數和訓練步長.
5)依據RMSE、MAE和MAPE指標評估Seq2Seq模型的預測效果.
本文使用的實驗數據來源于北京市氣象局,分為3類,如表1所示.
表1 輸入LSTM模型的站點特征向量
Table 1 Input feature vectors of LSTM

字段名稱描述ID站號PM2.5PM2.5濃度PM10PM10濃度vis能見度tem平均氣溫tem_max最高氣溫tem_min最低氣溫win_2m2米風速win_10m10米風速win_max極大風速pre小時降水rh相對濕度rh_max最大相對濕度prs平均氣壓
1)10個站點的逐小時觀測資料,包括PM2.5 濃度、 PM10濃度、能見度和氣象要素,數據每小時更新一次,具體見表1;
2)北京地區的網格化三維氣象要素客觀分析資料,空間分辨率為1公里,每小時更新一次;主要包括溫度、風和相對濕度要素;
3)與每個站點最臨近K個站點的PM2.5濃度數據,每小時更新一次.
原始數據經過預處理后生成序列數據集,其中每個時間序列數據由45個特征向量組成,每個樣本則包含96個時間序列(輸入72,輸出24)數據,批量訓練一次輸入64個樣本.數據集的時間跨度為2016年-2018年,其中80%的數據用來訓練模型,20%的數據用來驗證模型,測試集的時間跨度為2018年10月-12月.
在參考大量文獻的基礎上,本文中對模型中的部分關鍵參數做出了合理的設定:隱藏層為2層,每層隱藏神經元個數為64,學習率為10-3,衰減率為0.995,參數初始化范圍設置為[-0.08,0.08],模型優化算法選用Adam算法.過大的訓練步長會引起模型對訓練數據過度擬合,也會消耗更多的時間,而過小的步長容易引起局部最優,因此,訓練步長參數需要通過測試找到最優值.圖5顯示了訓練集和驗證集數據的MAPE隨訓練步長的變化.訓練集和驗證集上,MAPE誤差都隨迭代次數增加逐步減少.迭代次數超過3000左右時,該模型似乎過度適應,不僅泛化能力沒有改善,而且還出現微弱波動.因此,設置迭代步長為3000.

圖5 MAPE隨步長變化趨勢
為了方便模型之間的對比,除了CNN-Seq2Seq,本文還設計了以下幾種模型:基于機器學習的SVR模型、Just-Seq2Seq模型和Local-Seq2Seq模型.SVR廣泛用于分類或者回歸目的,在小樣本、非線性、高維模式識別等問題的解決上表現出許多特有優勢,具有較好的泛化能.本文使用的SVR模型參考文獻[31]調整參數,以獲得最佳性能.Just-Seq2Seq和Local-Seq2Seq結構與CNN-Seq2Seq的Seq2Seq部分完全相同,使用同樣的超參數,輸入序列長度也都相同,只是序列數據中包含的特征向量不同.Just-Seq2Seq模型的輸入僅包含用于生成時間特征的特征向量,而local-Seq2Seq模型的輸入僅包含站點自身的觀測數據,和最臨近K個站的數據,這里K取值為6.
4.4.1 幾種模型預測的誤差分析
基于驗證集,選擇RMSE、MAE和MAPE對SVR、Just-Seq2Seq、local-Seq2Seq和CNN-Seq2Seq這4個模型的預測結果進行評價.表2列出了4個模型對10個站點預測得到的整體精度.可以看出,盡管三個深度學習模型各項指標較接近,CNN-Seq2Seq模型預測效果仍最佳,在兩個指標中均優于所有其他模型.相比僅使用本站觀測數據的Local-Seq2Seq模型,對于PM2.5濃度,CNN-Seq2Seq對MAPE和RMSE的改進分別至少提高15.38%和26.48%.圖6(a)~圖6(b)顯示了四個模型在10個站上的RMSE和MAE誤差分布情況.CNN-Seq2Seq在十個站上的RMSE在[15,24]范圍區間,MAE在[11,16]范圍區間;Just-Seq2Seq在十個站上的RMSE在[20,38]區間,MAE在[11,25]范圍區間,僅有1個站的RMSE小于25;而Local-Seq2Seq模型的RMSE和MAE分布與Just-Seq2Seq類似,峰值比Just-Seq2Seq略高;SVR模型的預測誤差明顯高于其他三個深度學習的模型.參考前人基于深度學習的PM2.5小時濃度預測研究,當RMSE≤25、MAE≤15時,可以認為模型的預測效果較為理想,由此可見,CNN-Seq2Seq模型的預測效果是4個模型中最好的.
表2 不同模型預測結果的誤差比較
Table 2 Error comparison among different prediction models

模型MAPERMSESVR55.1945.79Just-Seq2Seq32.2321.42Local-eq2Seq Seq2Seq35.2524.87CNN-Seq2Seq29.8317.55

圖6 不同模型在10個站點上的PM2.5濃度預測誤差對比
4.4.2 不同深度學習模型的預測結果比較分析
從SVR、Just-Seq2Seq、local-Seq2Seq和CNN-Seq2Seq4個模型在驗證數據集上的預測的誤差分布可看出,三個深度學習模型的誤差均比SVR模型小,因此,在測試集上,僅對三個深度模型的預測準確性進行對比.選擇2018年11月30日-12月3日霧霾天氣過程為例,CNN-Seq2Seq模型、Just-Seq2Seq和Local-Seq2Seq模型下的ID為10的站點PM2.5逐小時濃度預測曲線隨時間的變化如圖7所示.由圖7可知,三種模型的預測曲線與實際觀測曲線的趨勢都基本保持一致.CNN-Seq2Seq模型的預測曲線與實際觀測曲線最為接近,尤其是在濃度波峰附近,表明該模型能較好的對PM2.5小時濃度峰值進行預測,這對PM2.5濃度的短時預測有著十分積極的作用.Just-Seq2Seq模型的預報曲線在PM2.5濃度波谷部分擬合有波動,對未來48小時后PM2.5濃度的預測偏離較大,預測效果不理想.local-Seq2Seq模型的預報曲線起伏最大,對PM2.5濃度峰值的預測值比實際觀測值偏大較多,和實際觀測曲線偏離較大.相比之下,CNN-Seq2Seq模型預測效果最好.

圖7 不同模型預測逐小時PM2.5濃度的誤差對比
由此可見,對于PM2.5未來24小時內的逐時濃度預測而言,考慮時空特征模型的準確度最高,包含周邊鄰近站點數據的時間序列模型的準確度高于僅考慮站點自身的時序模型,這說明空間相關性對PM2.5濃度預測很重要.
本文在分析北京地區PM2.5濃度隨季節變化規律、逐日變化趨勢以及氣象要素與PM2.5濃度相關性基礎上,提出了一種使用深度學習模型進行預測的方法,利用PM2.5觀測實況數據和多種氣象實況資料,建立CNN-Seq2Seq預測模型對PM2.5未來24小時的逐時濃度進行預測.通過幾個模型預測結果和誤差的對比表明,CNN-Seq2Seq預測模型能有效的獲取時空特征,適合解決時空序列數據的短時預測問題.在后續的研究中,將接入氣象數值模式預報產品和氣象格點預報產品,通過多種模型的對比分析,進一步提高模型的效率和精度.