李 幔, 馬元婧
1(中國科學院 沈陽計算技術研究所, 沈陽 110168)
2(中國科學院大學, 北京 100049)
隨著我國經濟水平的提高, 交通運輸業的飛速發展以及工業現代化進程的不斷加快, 空氣質量水平作為可持續發展的重要組成部分, 也越來越受到人們的關注. 為及時精準地了解空氣質量問題時空分布狀況,提高動態反應空氣質量的能力, 空氣質量監測應運而生. 根據多年來的實驗表明, 環境監測在采樣過程中產生的誤差是總體誤差的重要組成部分[1]. 因此, 保證采樣過程中監測點位選擇的準確度就成為了空氣質量監測與表征空氣污染程度的關鍵[2]. 而監測點位的優化準確度依賴于空氣質量監測微子站(以下簡稱為微子站)的監測數據的有效完整性, 但在微子站監測過程中, 常常由于各種突發因素導致一段時間內的數據缺失. 因此, 在根據微子站監測數據進行點位優化前, 進行缺失數據處理成為了不可或缺的一部分.
近年來, 隨著空氣質量監測點位優化問題研究的不斷深入, 已有眾多發達國家與地區對于空氣質量監測點位優選提出了一系列的優化方法, 如相關性分析法、聚類分析法以及多目標優化等數學與統計學方法[3].我國對于空氣質量監測點位優化的研究也在逐步發展,如徐明德等人提出的凝聚層次物元法[4]、杜增榮采用的凝聚層次聚類法[5]、姜林等人提出的BP 分析法[6]應用于空氣質量監測點位優化. 這些方法各有特色, 可以很好地代表空氣質量的準確性與代表性. 與此同時,隨著對缺失數據處理問題的關注度不斷提高, 大多被提出的數據填充方案大概分為以下兩種填充思想, 第一種是對缺失數據鄰近數據進行特征分析, 利用均值, 眾數和聚類的方法表現數據特征進行填補的算法. 如李董等人提出的基于SMOTE 和KNN 的數據缺失填充算法[7]、袁兆祥等人提出的基于DBSCAN 聚類的數據補充算法[8]、楊挺等人提出的基于FSOM 神經網絡的數據補充方案[9], 另一種缺失數據的填補方式是利用神經網絡對數據整體特征進行分析與補充, 例如使用RNN來對缺失的數據進行填補操作. 如柯昊等人提出的BP 神經網絡補差法[10]、宋維等人提出的基于LSTM的活立木莖干水分缺失數據填補方法[11], 對于數據缺失問題都有很好的表現. 因此, 本文首先選取在空氣質量、水質、地質預測領域應用廣泛, 原理簡單, 計算量小的凝聚層次聚類法進行空氣質量監測點位優化. 鑒于凝聚層次聚類法易受到初始聚類中心以及缺失數據的影響, 單獨使用凝聚層次聚類法的聚類效果不理想.為了能夠在空氣質量監測點位監測數據缺失的情況下,擁有更好的聚類效果, 本文首先選取適用于具有時序特征缺失數據補充且能夠反映空氣質量監測數據雙向特征的BiLSTM 神經網絡對缺失數據進行填補, 然后應用凝聚層次聚類法對修復后的數據進行聚類分析,實現對空氣質量監測點位的優化, 獲取更加真實、精準的數據, 用以分析空氣污染成因. 同時, 為環境質量監測有關部門制定點位優化標準提供有利的技術支撐.
針對空氣質量監測點位優化的最終目的, 構建的基于BiLSTM 改進聚類的空氣質量監測點位優化模型, 如圖1 所示. 首先, 利用BiLSTM 數據補充模型進行數據補充解決監測數據缺失問題, 對于BiLSTM 中的正向LSTM 網絡和逆向LSTM 網絡具備相同的結構, LSTML部分是通過學習缺失數據之前時序的空氣質量數據的數據特征規律來產生缺失數據的候選補充結果; LSTMR是通過學習缺失數據之后時序的空氣質量數據進行同樣的操作, 產生另一候選補充結果. BiLSTM神經網絡通過使用缺失監測數據的前后m組數據對當前缺失的n組數據進行補充, 公式如下:

圖1 基于BiLSTM 改進聚類模型

其中,F為基于BiLSTM 的數據補充模型, [Xk···Xk+n?1]表示缺失的n組數據.
最后使用一個Softmax layer 來綜合正向LSTM和逆向LSTM 的候選結果并產生最終的缺失數據補充結果. 得到最終的輸出結果傳遞給層次聚類部分, 輸出聚類結果.
LSTM 神經網絡是一種循環神經網絡(RNN)的特殊變體. 在標準的RNN 中, 重復部分的模塊表現為只有單個tanh 層的簡單結構, 如圖2 所示. LSTM 的基本結構與基礎ANN 保持一致, 但是重復部分模塊包含具有4 個網絡層的特殊結構, 并且利用一種新的方式與其他部分實現交互[12]. 因此, 解決了RNN 的長期依賴以及梯度消失問題[13].

圖2 標準RNN 結構圖
LSTM 的網絡結構如圖3 所示, 其中包含輸入門、遺忘門、輸出門3 個邏輯單元, 輸入門控制記憶單元中當前輸入的狀態; 遺忘門對前一個記憶單元處理結果進行篩選保留; 輸出門控制記憶單元的輸出狀態[14].

圖3 LSTM 神經網絡結構圖
LSTM 的初始步驟在于判斷細胞狀態中信息的去留. 這一判斷由遺忘門層進行.
遺忘門的狀態更新公式為:

它接收ht–1和xt, 對于細胞狀態Ct–1中的每部分的輸出值都限制在0 到1 之間. 接受程度跟隨0 到1 數值大小遞增, 0 表示完全不接受, 1 表示完全接受.
下一步的目的在于判斷細胞狀態中需要更新的信息并建立新信息. 通過輸入門層決定需要被更新的數據. 然后, 通過一個 tanh形網絡層創建一個新的備選值向量Ct, 可以用來添加到細胞狀態.
輸入門狀態更新公式:

通過下述公式對細胞狀態的更新:

最后, 確定輸出值. 首先, 確定細胞狀態中的輸出部分. 然后, 我們把細胞狀態輸入tanh 與s 形網絡層的輸出值相乘, 實現輸出.
輸出門狀態更新公式如下:

對于LSTM 來說, 只能學習數據對之前數據的依賴性, 無法充分考慮對之后數據的依賴性[15]. 在面對數據缺失問題時, 需要結合數據缺失前后數據來共同決定當前缺失數據補充. 圖3 為BiLSTM 模型圖.
BiLSTM 可以通過雙向計算, 不僅能考慮到缺失數據之前的監測數據記錄, 還能捕捉到后續監測數據.將前向LSTM 模型和后向LSTM 模型構成BiLSTM模型來學習雙向監測數據信息.
本文采用聚類中的適用于少量特征、少量點位的凝聚層次聚類法(以下簡稱為層次聚類法)對沒有預先處理的微子站點位監測數據進行分類[16–19]. 如圖4 所示, 每個樣本先自成一類, 然后按距離準則逐步合并,減少類數. 首先, 將微子站分為n類(每個微子站為一類), 計算類間距離, 將類間距最為接近的兩類合并為一類, 循環計算n–1 類之間的間距, 不斷合并各個微子站, 直到所有類都完成歸類, 聚類結束[20]. 聚類分析根據類間距離的計算方法不同, 層次聚類法也不同. 本文采用最為常用的最小距離法計算類間距離.

圖4 BiLSTM 神經網絡結構圖
最小距離公式:

其中,D(Xm,Xn)為M類中的某個樣本Xm和N類中的某個樣本Xn之間的歐式距離.
歐式距離公式:

層次聚類法示意圖如圖5 所示.

圖5 層次聚類示意圖
為驗證本文提出的基于BiLSTM 改進的聚類算法的效果, 本實驗數據集來源為沈陽市渾南區18 個微子站2020 年1 月–2021 年8 月的SO2、NO2、O3、PM10、PM2.5、CO 六種大氣污染物的日均監測數據, 表1 為其中某月各站點6 種污染物日均監測數據.

表1 主要污染物濃度監測日均值 (μg/m3)
實驗中選取的18 個監測站點位置如圖6 所示.

圖6 監測站點位置
部分原始數據如圖7 所示, 在圖7 中橫軸代表的是部分月份各站點的CO 日均值原始數據組數, 縱軸代表的是CO 原始數據的濃度值. 從圖上可以看出部分原始數據存在缺失的情況. 因此, 在使用此數據集進行模型訓練和驗證之前需要對原始數據進行數據預處理, 補充缺失數據.

圖7 CO 原始數據
同時, 數據的標準化處理是數據處理過程中的至關重要的一環, 將原始觀測數據矩陣中的每類元素, 按照某種特定的運算把它變為一個新值, 使其更具有可比性. 因此, 在聚類分析之前, 對數據進行預處理, 將原始數據映射到(0, 1)之間, 即均勻化原始數據公式如下:

本文實驗所采用的硬件配置為: Intel(R) Core(TM)i5-1035G1 CPU@1.0 GHz 的處理器, 內存DDR4 (16 GB);軟件環境為: 編輯器PyCharm 2019.3.2, Python 3.7,PyTorch 1.5.
實驗過程如圖5 所示, 首先, 對數據進行標準化、歸一化處理. 由于原始監測數據存在缺失, 因此, 對序列缺失處, 補上日期索引并特殊值標記該時間步, 之后,Masking 層接收經過預處理后的序列, 為后續BiLSTM層過濾掉標記為缺失的時間步, 為保證監測數據季節完整性以及客觀地展示本文方法的補全效果, 將2020年1–12 月的數據作為訓練集, 2021 年1–4 月的數據作為驗證集, 輸入上述的BiLSTM 網絡模型中, 開始訓練, 設置誤差平方和E=0.001; 學習率 η=0.01. 經過589 次迭代訓練網絡達到要求. 最后, 將2021 年5–8 月的數據隨機選取300 個監測數據, 刪除其真實值作為缺失數據, 并將處理后數據作為測試樣本輸入已訓練好的網絡中進行估計, 得到缺失數據的估計值. 圖8 為部分缺失數據估計值與真實值的誤差曲線圖, 可以看出BiLSTM 神經網絡模型結果與實測值變化趨勢基本一致, 由此可得出該神經網絡模型訓練結果良好.

圖8 監測數據補充值與真實值對比
為了更好地驗證BiLSTM 神經網絡模型的監測數據插補方法的效果, 本文將提出的數據缺失補充模型與以往在空氣質量數據缺失情況下常用的均值補插法以及使用單向LSTM 神經網絡進行對比. 為了評價模型描述實驗數據的精確度與反映補充值誤差的實際情況, 本文采用計算簡便的均方根誤差(RMSE) 與具有更強魯棒性的絕對平均誤差(MAE)來共同分析補充結果和真實測量值之間的偏差. 表2 是3 種模型關于部分站點數據補充結果的評價指標的對比, 一般來說當預測值與真實值之間的偏差越小, 也就是說MAE和RMSE 的值越小, 模型的補充效果越好.

表2 不同模型下缺失數據填補實驗結果對比
由以上分析可以得出本文提出的模型在空氣質量數據缺失的補充效果表現更優. 進而在此基礎上, 基于補充數據進行聚類分析, 從18 個采樣點中可選出1、2、3、5、6、8、9、10、13、14 共10 個優化點, 如圖9 所示.

圖9 填補數據層次聚類結果圖
為驗證算法選取的優化點位準確度, 對比單獨使用原始數據進行層次聚類法選取的優化結果, 如圖10所示, 發現與單獨使用聚類分析法選擇的1、2、4、3、13、16、8、5、6、10、14、17 等12 個點位基本保持一致.

圖10 原始數據層次聚類結果圖
本文通過相關性檢驗, 推斷點位優化前后點位的可靠性, 通過選取沈陽市渾南區2021 年4–8 月6 種主要污染物濃度日均實測值進行相關性檢驗, 判斷優選后的微子站能否客觀地反映空氣質量.
相關系數公式:

其中,Xi,Yi指優化前后污染物濃度值.
在給定α=0.05 時, 相關性表如表3.

表3 各污染物相關性系數表
根據表1 所示的相關性檢驗, 表明原點位與優化后點位的各污染物濃度密切相關, 無明顯差異性.
為進一步驗證優化前后監測數據的一致性, 選取部分監測數據采用方差齊性檢驗-F 檢驗法, 對優化前后數據的顯著性差異進行檢驗. 在顯著性差異不明顯的前提下, 繼續進行t檢驗, 確定兩組數據是否具有一致性. 結果如表4 所示, 可以得出優化前后各污染物濃度之間具有一致性. 由此表明優化點位具有較好的代表性, 點位選擇較為準確, 表明優化后的點位可以替代原點位.

表4 一致性檢驗結果
聚類屬于無監督學習中的一種, 原始數據不具有標簽, 需要一些定量的指標來進行評估算法的好壞[19,21].根據是否提供樣本的標簽信息, 相關的指標可以分為外部方法與內部方法, 內部方法指的是不需要數據的標簽, 僅僅從聚類效果本身出發, 而制定的一些指標[16].因為點位優化屬于無標簽數據, 故選擇內部方法的2 種評價指標:CH指數(Calinski-Harabaz index,CH)、戴維森堡丁指數(Davies-Bouldin index,DBI)作為評價標準, 對基于補充數據后數據集與基于原始數據集的聚類效果進行評價.
CH指數, 綜合考慮了簇間距離和簇內距離, 計算公式如下:

CH的數值大小與簇內距離成負相關, 簇間距離越大, 聚類效果越好.
戴維森堡丁指數, 公式如下:

其中, avg(C)表示聚類簇的緊密程度, 公式如下:

計算該聚類簇內樣本點的距離,d表示不同聚類簇中心點之間的距離, 公式如下:

聚類簇間距離越遠, 聚類簇內距離越近,DBI指數越小, 聚類算法性能越好.
BiLSTM 神經網絡填補監測數據與原始數據兩種狀態下, 進行層次聚類后, 聚類效果評價指標對比如表5所示, 可以看出相比較于對原始數據進行聚類操作, 本文所采用的基于BiLSTM 改進的聚類算法聚類效果有所提升.

表5 評價指標對比
根據BiLSTM 神經網絡補充數據的改進聚類算法從18 個監測點位中選出了10 個監測點位, 與在原始數據上使用層次聚類法選出的12 個點位基本保持一致, 對比單獨使用聚類算法選出的12 個監測點位. 此模型選擇的點位更少, 并且聚類效果更好, 選擇的點位代表性與準確性均得到提升, 為之后監測點位布設與優化提供一個更加高效準確的模型. 由表3 的計算結果可以推斷本文提出的方法選擇的點位是準確且具有代表性的.
本文針對空氣質量監測點位優化的最終目的, 通過層次聚類法對BiLSTM 神經網絡補充后的微子站監測數據進行聚類處理, 在選擇更加少量點位的基礎上,提高了點位選擇準確度, 提出了一種基于BiLSTM 神經網絡改進聚類的點位優化算法.