基于LSTM-GCN的PM2.5濃度預測模型

馬俊文，嚴京海，孫瑞雯，劉保獻

北京市生態環境監測中心，大氣顆粒物監測技術北京市重點實驗室，北京 100048

在環境監測部門開展空氣質量監測預報分析業務中，目前采用的主流空氣預報分析模型仍然主要是空氣數值預報模型。從研究大氣流體運動學和應用大氣內部流體物理學的規律關系出發，例如利用大氣流體動力學、熱力學等，用傳統數理分析方法建立相應的具體數學大氣物理規律模型，然后可以借助大型高性能計算機的計算能力，建立空氣污染物濃度的動態分布運輸和擴散數值預測模型。例如：XU等[1]基于尺度空氣質量模式系統(CMAQ)提出了CMAQMOS模型，很好糾正了CMAQ中平均污染物排放清單所導致的系統性預測誤差；CHUANG等[2]應用結合化學的天氣研究與預測模型(WRF-Chem-MADRID)對美國東南部地區空氣質量進行預測,結果顯示在臭氧和可吸入顆粒物預測中表現良好。

近年來，研究人員嘗試基于機器學習的深度學習算法預測PM2.5濃度。例如，在分析時間特征維度，王瑋等[3]開展了GRU神經網絡在可吸入顆粒物預測中的應用研究，白盛楠等[4]開展了LSTM神經網絡在可吸入顆粒物濃度預測中的應用研究。然而，由于空氣質量在時間維度具有很強的季節性，因此，研究人員開始將季節性分解算法和深度學習算法相結合，開展了PM2.5濃度預測模型的研究。例如，劉銘等[5]提出EMD分解算法和LSTM算法的組合算法EMD-LSTM，開展了PM2.5濃度預測模型研究；劉炳春等[6]提出小波分解變換算法Wavelet和LSTM算法的組合算法Wavelet-LSTM，開展了空氣污染物濃度預測模型研究；王曉飛等[7]提出Prophet分解方法和LSTM的組合算法Prophet-LSTM，開展了PM2.5濃度預測模型研究。

但是，空氣質量監測數據除了具有時間特征外，還具有空間特征，隨著污染物傳輸過程以及氣象因素的影響，各監測站的監測數據具有很強的空間傳輸演化特征。因此，研究人員開始關注基于時空特征的PM2.5預測模型。例如，宋飛揚等[8]利用KNN算法提取目標站點空間因素，提出KNN算法和LSTM算法的組合算法KNN-LSTM，開展了PM2.5濃度預測模型研究；劉旭林等[9]將地理空間網格化，利用CNN提取空間特征，提出CNN和Seq2seq的組合算法CNN-Seq2seq，開展了未來1 h PM2.5濃度的預測研究。

現有研究中時間特征提取基本是使用LSTM模型實現，空間特征提取基本是將非歐氏距離空間的監測站分布數據通過地理空間網格化后轉換到歐氏距離空間，然后利用面向歐氏距離空間的二維卷積神經網絡或KNN算法提取空間特征。2017年，Thomas Kpif提出了專門面向圖提取圖結構特征的圖卷積神經網絡(GCN)，可以直接用于處理非歐氏距離空間監測站地理分布數據。黃偉建等[10]利用圖卷積神經網絡改進門控循環單元網絡的輸入部分，提出面向空氣質量的時空混合預測模型。但是，該時空混合預測模型并未分別提取時間特征和空間特征。因此，本文基于LSTM算法構建時間特征提取組件，基于GCN算法構建空間特征提取組件，然后將時間特征和空間特征進行聯合，提出預測未來1 h PM2.5濃度的LSTM-GCN組合模型。

1 LSTM-GCN組合模型構建

1.1 LSTM模型

機器學習領域的深度學習網絡算法主要有3類：第一類是循環神經網絡，適用于處理序列數據；第二類是卷積神經網絡，適用于處理圖形圖像數據；第三類是全連接神經網絡，適用于解決分類問題。其中，LSTM網絡是循環神經網絡的主要變體之一,它在充分繼承了循環神經網絡遞歸運算優點的同時,能夠很好地解決目前循環神經網絡在應用深度機器學習算法過程中,產生的梯度遞歸消失和梯度遞歸爆炸問題。

LSTM結構的核心思想是在循環神經網絡基礎上，引入一個叫單元狀態的連接，最后的狀態不再簡單地存儲，而是通過LSTM神經網絡的訓練機制來選擇更新狀態。LSTM 的訓練機制是通過遺忘門、輸入門和輸出門來控制信息的刪減或增加[11]。

1.2 GCN模型

GCN是一種基于圖數據的深度學習方法，其核心思想是學習一個映射函數，可以使圖中的節點能夠聚合自身和鄰居節點的特征，生成節點的新表示。卷積神經網絡可以很好地處理歐氏空間數據以提取空間特征，對于非歐氏空間里的圖數據，可以使用GCN提取空間特征。

(1)

第二步是提取特征矩陣X的空間信息,具體換算方式見式(2)。

(2)

式中：σ(·)為激活函數；W(i)為第i層權值矩陣；H(i)為第i層的激活值；H(0)=X。

1.3 LSTM-GCN組合模型

設計LSTM-GCN組合模型結構如圖1所示。

圖1 LSTM-GCN模型結構Fig.1 LSTM-GCN model structure

圖1中，LSTM-GCN模型包含時間特征組件和空間特征組件。時間特征組件的輸入是多條時序數據，針對每條時序數據應用LSTM模型，然后將其輸出合并于特征聯合層，由多維輸入一維化層輸出所提取的時間特征；空間特征組件的輸入是鄰接矩陣A和特征矩陣X，然后經過GCN模型、全連接層后，由多維輸入一維化層輸出所提取的空間特征。最后，將時間特征和空間特征進行特征聯合后，經全連接層輸出預測值。

選擇均方根誤差(RMSE)和平均絕對誤差(MAE)為模型預測效果評價指標，具體換算見式(3)和式(4)。

(3)

(4)

2 時空趨勢分析和數據預處理

以北京市35個空氣質量監測站2018—2020年PM2.5小時監測數據、氣象監測數據、站點分布經緯度數據為樣本，首先分析PM2.5時空變化趨勢，然后構建LSTM-GCN模型輸入數據和預測數據。

2.1 時空趨勢分析

Prophet時序模型是Facebook公司于2017年發布的，它可以分解出時序數據的趨勢，且彌補了傳統時序模型對時序數據過于局限、需要填充缺失值、模型不靈活等不足。因此，基于該模型分解35個監測站PM2.5濃度數據的趨勢分量，然后按年計算趨勢分量的均值，繪制用于表征時空趨勢的雷達圖，如圖2所示。

圖2 時空趨勢圖Fig.2 Temporal and spatial trends

從圖2可知，2018—2020年，在時間趨勢方面，除了第15#站點的趨勢分量年均值在2018年較低以外，其余34個站點的趨勢分量年均值均隨時間逐年下降；在空間趨勢方面，站點之間的趨勢分量年均值相對變化趨勢各年基本一致。其中，第15#站點在2018年進行了停站改造，影響了數據完整性，導致趨勢分量年均值偏低。

2.2 時間特征組件輸入數據和預測數據構建

假設基于t時刻之前的連續24 h PM2.5濃度數據，預測t時刻PM2.5濃度，則設計構建時間特征組件輸入數據以及對應的模型預測數據的程序流程分6個步驟。

1)對PM2.5小時數據按時間順序排序，數據有35列，分別對應35個監測站。

2)剔除極大和極小異常值，按列執行歸一化，即對每個監測站的PM2.5小時數據分別進行歸一化。

3)定義存儲時間特征組件輸入數據的變量M，存儲LSTM-GCN模型預測數據的變量Y。

4)建立滑動窗口，窗口大小為25，并置于數據頂端。

5)提取滑動窗口內的數據，前24條數據追加到變量M，第25條數據追加到變量Y。

6)如果窗口滑動到達數據底端，則終止程序，否則，滑動窗口向下一小時滑動，返回第5步。

時間特征組件輸入數據是35個站點t時刻之前的連續24 h PM2.5濃度數據，用Mk表示，則其數據格式如式(5)所示。

(5)

式中：k=1,2,…,35。

LSTM-GCN的預測數據是35個站點t時刻PM2.5濃度預測值,用Yt表示，則其數據格式如式(6)所示。

(6)

2.3 空間特征組件輸入數據構建

空間特征組件用于提取35個監測站在空間分布上的演化趨勢特征。其中，監測站在北京市的空間分布情況如表1所示。

表1 35個站點空間分布Table 1 Spatial distribution of 35 stations

空間特征組件的輸入數據有站點分布圖和站點特征數據。站點分布圖是以站點為頂點，以站點之間的距離倒數為邊權構建的加權無向圖，站點分布圖用鄰接矩陣表示，其數據格式如式(7)所示。

(7)

式中：yij為站點i和站點j之間地理距離的倒數，yii=0；i,j=1,2,…,n；n=35。

站點之間的地理距離可以由站點經緯度數據計算得到，假設站點P的經緯度為(P1，P2)，站點Q的經緯度為(Q1，Q2)，則計算站點P和站點Q之間的地理距離D，具體換算方式見式(8)。

D=R×arccos[sin(P2)×sin(Q2)+

cos(P2)×cos(Q2)×cos(P1-Q1)]×π/180

(8)

式中：R為赤道半徑，取值為6 371.004 km。

站點特征數據是t-1時刻每個站點所具有的特征，假設每個站點有m個特征，則形成的特征矩陣如式(9)所示。

(9)

式中：zij是站點i的第j個特征；n=35。

因為每個監測站點除了監測PM2.5外，還監測相關氣象數據，且氣象條件的變化對PM2.5有影響[12]，因此，使用PM2.5、風速、風向、溫度、濕度構建監測站點的特征數據。數據預處理方法：對PM2.5、風速、溫度、濕度4列數值型數據，剔除異常數據后分別執行歸一化；對風向數據，按照東、南、西、北、西北、東北、東南、西南8個方向，使用獨熱編碼(One-Hot Encoding)方法，編碼為由0和1組成的8位數據，每個數據位代表1個方向。例如，若8位數據中代表方向東的數據位為1，代表其他方向的數據位為0，則表示風向為東風。以上2類數據按列合并后，形成包含12列的特征數據。

按照時間維度依次提取35個站點同一時刻的特征數據，形成特征矩陣，則其階數為35×12階，行代表35個監測站點，列代表監測站點特征。

LSTM-GCN模型輸入和輸出矩陣的階數如圖3所示。

圖3 模型輸入/輸出矩陣階數Fig.3 Input/output matrix order of model

圖3中，LSTM-GCN模型輸入有3個矩陣，分別是24×35階的PM2.5時間序列矩陣、35×35階的站點分布圖鄰接矩陣、35×12階的監測站特征矩陣，輸出預測值是1×35階的矩陣。

3 仿真實驗

本文實驗基于TensorFlow開發，它是一個基于Python語言的較為成熟的開放源代碼的深度學習框架，應用較為廣泛，尤其在自然語言處理、推薦系統、圖形分類、音頻處理等領域應用較為豐富[13]。

模型化的開發工具可以選擇Keras來實現,它是一個可以在Tensorflow里被導入、能夠被用于開發深度學習模型的高級API,在極大程度簡化模型的基本構建設計過程和模型訓練方法的同時,還可以完全兼容在Tensorflow里已有的模型開發相關的絕大部分功能,具有易于模塊化、可彈性組合、用戶友好、易于進行擴展的四大優點[14]。此外，針對復雜模型，可以使用Keras函數式API進行定義，然后把模型看作層，使用傳遞張量的方法來調用模型。

本實驗除了實現LSTM-GCN模型外，還分離時間特征組件和空間特征組件，形成2個獨立的LSTM模型和GCN模型，作為評估LSTM-GCN模型的對比模型。

模型訓練和預測程序流程圖如圖4所示。

圖4 模型訓練和預測程序流程圖Fig.4 Flow chart of modeltraining and prediction

從圖4可見，樣本數據預處理后，按照時間順序排序，劃分前80%為訓練數據，后20%為測試數據，設定模型激活函數為ReLU[15]，損失函數為均方誤差(MSE)，優化方法為自適應運動估計算法(Adam)。

同時，將時空地理加權回歸模型(GTWR)作為對比模型，GTWR是一種通過構建時空非平衡關系的局部線性回歸模型，設置模型解釋變量為t時刻PM2.5濃度、風速、溫度、濕度，因變量為t+1時刻PM2.5濃度。

根據每個目標站點上實際值和4個目標模型的預測值,分別計算所對應的4個模型的評估指標RMSE和MAE，計算結果如圖5和圖6所示。

圖5 各監測站RMSE對比圖Fig.5 RMSE comparison of monitoring stations

從圖5和圖6可知，4個模型對35個監測站的預測精度具有較為一致的波動趨勢，且LSTM-GCN模型預測精度最高，GCN模型預測精度最差。

圖6 各監測站MAE對比圖Fig.6 MAE comparison of monitoring stations

計算35個站點RMSE和MAE的平均值，結果如表2所示。

表2 模型評估結果Table 2 Model evaluation results

從表2可見，LSTM-GCN模型評估結果優于LSTM模型、GCN模型、GTWR模型，相較于LSTM模型RMSE、MAE分別降低了11.68%、7.34%；相較于GCN模型RMSE、MAE分別降低了40.22%、36.37%；相較于GTWR模型RMSE、MAE分別降低了17.52%、23.69%。

以1#站點為例，抽取同一時間的真實值和4個模型的預測值，繪制曲線如圖7所示，以展示時間維度模型預測效果。

從圖7可知，在時間維度，4個模型的預測值基本隨著實際值波動，且趨勢較為一致，預測效果較好。

圖7 1#站PM2.5預測效果Fig.7 PM2.5 prediction effect of 1# station

隨機選擇某時刻，獲取35個監測站實際值和預測值，繪制曲線如圖8所示，以展示空間維度模型預測效果。

從圖8可知，在空間維度，4個模型預測值的波動幅度和演化趨勢同實際值較為一致，預測效果較好。

圖8 35個站PM2.5預測效果Fig.8 PM2.5 prediction effect of 35 stations

應用LSTM-GCN模型預測北京市35個監測站2021年1—7月的PM2.5數據，隨機抽取部分時間的預測值和實際值，繪制對比圖，其中6個站的對比圖如圖9所示。

圖9 6個監測站的LSTM-GCN預測值和實際值對比圖Fig.9 Comparison between LSTM-GCN predicted valueand actual value of six monitoring stations

從圖9可知，6個監測站的預測值很好地擬合了實際值，說明預測效果較好。

4 結論

本文提出的LSTM-GCN模型同時考慮空氣質量數據時間特征和空間特征，從模型預測和評估結果來看，能夠進行PM2.5數據的預測，且預測精度比LSTM模型、GCN模型、GTWR模型均有所提升。LSTM-GCN模型在2021年1—7月的PM2.5濃度預測效果表明，該模型具有一定的可用性。

但是，LSTM-GCN模型沒有分解PM2.5的季節性，后續可以再研究疊加一種適用的季節性分解算法，并增加LSTM和GCN的層數，進一步優化模型。

LSTM-GCN模型相比傳統數值模型具有數據需求少、運行速度快、計算資源需求少等優勢，有望在短時預測、數據審核、異常數據判定等業務中成為一種技術支撐手段。