基于深度學習的溫度觀測數據長時間缺失值插補方法①

鄭欣彤,邊婷婷,張德強,賀 偉

1(中國科學院 地理科學與資源研究所 資源與環境信息系統國家重點實驗室,北京 100101)

2(中國科學院大學 資源與環境學院,北京 100049)

3(北京聯合大學 管理學院,北京 100101)

4(中國科學院 華南植物園鼎湖山森林生態系統定位研究站,廣州 510650)

1 前言

溫度是農業、生態系統研究非常重要的觀測量,農業作物生長的模擬、農業氣象災害監測、生態系統模擬中溫度是必不可少的輸入[1,2].隨著農業、生態模擬的精細化,要求溫度數據具有更高的精度,如農業氣象災害干熱風監測、林塊生態系統碳排放的模擬等[3,4],高精度的溫度觀測量是必不可少的輸入參數.溫度觀測數據一般通過野外氣象觀測站獲取,由于設備故障、惡劣環境或是認為操作失誤等原因,小氣象觀測難免會出現缺失[5],缺失數據插補或補全,是溫度觀測數據運用前必不可少的預處理工作.

本文針對中國廣州一個森林生態站長時間間隔溫度觀測數據缺失進行插補方法研究.由于該森林生態站處在雷雨區,小氣象站夏季容易因惡劣天氣損壞,容易造成較長時間的數據缺失.論文選擇了該生態站同時有自動觀測氣象站,作為比對觀測,該小氣象場還有人工溫度觀測設施.自動觀測氣象數據輸出的觀測產品時間頻率為30 分鐘,每天有48 條觀測記錄數據;人工觀測分早、中、晚每天3 次,產生3 條記錄.本文研究的實際應用問題:如何運用不同數據插補方法,通過每天低頻的人工溫度觀測數據,獲取完整的高精度半小時頻率溫度觀測數據.

數據插補是眾多學科數據分析前必不可少的預處理工作.目前以數據插值、統計分析和時間序列分析等為基礎,發展出了多種數據插補方法[6,7],但對高精度數據插補研究還處初步階段[5].高精度數據插補的要求是從已觀測數據中學習到數據的規律或缺值模式,從而實現對未觀測數據的準確估計.深度學習是機器學習領域一個新的研究方向,是人工智能領域的一項顛覆性技術創新,除了帶來圖像、語音和自然語言處理領域的突破,也成功應用到了眾多學科領域[8,9].深度學習旨在獲得樣本數據的內在規律和表示[10],和數據插補的需求非常契合.

深度學習技術已經在交通、醫療、傳感器網絡等多個領域的數據插補中獲得了成功應用[11-13],并發展出了GRU、LSTM、GAN 等不同結構的數據插補深度學習神經網絡[14].GRU(gate recurrent unit)和LSTM(long short-term memory)都是循環神經網絡的不同形式,可以解決RNN 網絡學習過程中的梯度消失或爆炸問題[15].這兩種結構神經網絡在數據插補應用中,不但可以從已觀測數據之中學習規律,也可從數據缺失值模式中進一步學習,提高數據插補精度[16,17].GAN(generative adversarial networks)網絡用于學習多變量時間序列的總體分布,從而對觀測數據中的缺失值進行插補[18].

在眾多深度學習時間序列數據插補模型中,一種基于序列-序列(Seq2Seq)的深度學習模型在多個標準樣本集數據插補都有很好的表現[19],該結構采用了雙向循環LSTM 網絡,在隨后的實際應用也進一步驗證了該結構適用于時間序列數據缺失插補問題[20].另外,基于Encoder-Decoder 結構的深度學習神經網絡在數據插補方面也獲得成功應用[21].這些不同結構的深度學習模型為本文研究提供了重要參考.

本文運用低頻人工溫度觀測數據,來插補高頻次機器觀測數據的長時間觀測值缺失問題.時間序列數據缺失值插補雖然已很豐富,但針對這一特定應用場景的數據插補方法研究文獻還是較少[22].下面是一個具體的生態臺站的小氣象觀測數據,該小氣象站同時具有半小時自動溫度觀測和每天3 次人工溫度觀測,但半小時自動溫度觀測數據存在較長時間的缺失值.為了實現半小時溫度觀測數據的高精度插補,文中詳細給出了一個編碼-解碼結構的序列-序列深度學習溫度插補模型的構建過程和數學公式,并將其與其他插補方法進行了插補精度對比分析.

2 研究數據介紹

本文研究采用了我國廣州鼎湖山森林生態系統國家野外科學觀測研究站的氣象溫度觀測數據.該生態系統觀測站開展有溫度觀測對比試驗,同時開展人工觀測和氣象機器自動觀測活動,有較長時間的溫度觀測數據記錄,表1 是用于本文研究的溫度人工觀測數據和自動機器觀測數據情況.

表1 溫度觀測數據集信息表

由于鼎湖山生態站位于中國南方的山區,自動觀測設備容易受雷雨季節影響而產生較長時間觀測記錄的缺失.圖1 是某一機器自動氣象觀測數據的數據缺失情況分布圖,從圖中可見在2020年7月有一次超過2 個月的溫度觀測數據缺失.

圖1 半小時溫度觀測數據缺失值分布圖

運用每天3 次的人工觀測數據對2020年7月超過2 個月的機器溫度觀測數據缺失進行插補,是本文方法研究的應用背景.由于人工觀測數據和機器觀測數據之間很容易建立線性相關關系,本文數據插補問題研究核心是如何運用低頻的人工溫度觀測數據對高頻的機器自動溫度觀測數據進行插補.

3 溫度觀測數據深度學習插補

3.1 基本定義

氣溫具有天的周期性,很自然的將長時間序列溫度觀測數據按照天為單位進行劃分,變為每天48 個觀測值的分段序列.為了更好的把研究集中到長時間間隔缺失值的插補,對時間序列中偶爾或短時間數據缺先采用上文的卡爾曼平滑的方法進行插補.這樣以天為分段的溫度時間序列包含兩種,即沒有缺失值的每天分段,記為,和只包含早中晚3 次觀測值的每天分段,記為,這樣以天為分段的時間分段序列可以表示為:

序列(1)表示長度為n天溫度時間序列中,缺失值窗口寬度為m天.半小時溫度觀測序列(1)為長度為48n,存在缺失值48m的半小時溫度觀測數據序列,以天分段后的表達形式.

為了表示時間序列(1)中缺失值的位置,對長度為L(48n)的半小時采樣溫度時間序列,構建相應長度為L的掩碼時間序列,其中:

現在以天為單位,對長度為L的半小時掩碼序列進行分段,沒有缺失值的掩碼每天分段記為和只包含早中晚3 次觀測值的掩碼每天分段,記為這樣就可以建立與式(1)對應的以天為分段的掩碼序列:

3.2 滾動窗口采樣

采用滾動窗的方法,基于以天為分段的時間序列為深度學習模型訓練構建樣本集.對長度為m(天)的缺失值進行插補,需構建樣本滾動窗口的長度大于m,并且在m的兩端各保留長度為s(天)的觀測數據,這樣滾動窗口長度w為m+2×s天.訓練樣本為適應序列-序列(Seq2Seq)的訓練方法來構建,對長度為w的訓練輸入樣本中溫度觀測序列為:

可通過訓練形成如下的時間序列結果輸出:

訓練樣本需要在沒有缺失值的溫度觀測序列基礎上構建,樣本中觀測值缺失的模式同實際情況一致,即每天只有早中晚3 次觀測值.表2 是訓練樣本中存在缺失值的某一天溫度數據及其對應的掩碼示例.

表2 樣本序列中缺失值窗口內某一天的數據示例

3.3 深度學習模型的設計

典型的基于Seq2Seq 的時間序列數據插補深度學習模型有SSIM,BRTS-I 等[19,21].本文吸收了這些模型的優點,將Seq2Seq 和Encoder-Decoder 深度學習架構結合起來,所設計的深度學習模型結構如下文所述.上面輸入序列(1)被記為x={x1,x2,…,xn},輸出序列(4)被記為y={y1,y2,…,yn},掩碼序列(2)被記為m={m1,m2,…,mn}.

(1)編碼

從圖2 中可見,深度學習結構中的編碼部分的基本結構是LSTM-I,該結構與BRTS-I 結構中的RTSI 結構相似,其中的循環神經網絡單元直接采用了長短期記憶單元;另外,本文溫度觀測缺失值部分,每天48 個半小時溫度值,只有3 個觀測值,有效值比較稀疏,所有沒有采用RTS-I 中的缺失值時間間隔的變量和相應的訓練公式.下面定義中LSTM 被簡化為一個簡單算子的形式,將LSTM-I 單元過程數學描述為:

式(5)將上一個LSTM 單元的隱狀態ht-1轉化為估計向量其中Wx、bx為模型參數;式(6)通過運用掩碼向量mt,把輸入向量xt中的缺失值替換為估計向量對應的值;式(7)通過LSTM 網絡單元把和隱狀態ht-1產生預測狀態ht;式(8)是LSTM-I 單元的估計誤差,為缺失值位置上觀測值與估計值絕對差的累計量.

圖2 中神經網絡編碼部分由雙向的LSTM-I 神經網絡構成:一個是從時間序列的開始到結束讀取輸入,產生前向隱狀態向量序列另一個是從時間序列的結束到開始反向讀取輸入,產生后向隱狀態序列將前向和后向隱狀態序列拼接到一起,構成編碼層的編碼輸出其中向量hi為:

圖2 溫度值缺失值插補神經網絡結構圖

雙向編碼LSTM-I 編碼網絡誤差包括正向和逆向估計誤差兩部分.

(2)解碼

解碼層接收編碼輸出序列h,產生插補后的時間序列結果序列y.本文通過優選,采用了2 種解碼結構:一種是全連接層和一層線性層組合,如圖3所示;另一種是LSTM 和一層線性層的組合,如圖4所示.

圖3 全連接層和線性層組合的解碼詳細結構

圖4 LSTM 層和線性層組合的解碼結構

①基于兩層線性層組合解碼過程數學描述如下:

式(10)中的Dropout通過對解碼輸入隨機丟棄部分神經元,能夠起到預防過擬合的作用;式(11)為全連接層,g為激活函數,全連接層產生輸出狀態序列s={s1,s2,…,sn};由于溫度值是連續值,式(12)為最上層即線性變換層,輸出插補結果序列y;式(13)是解碼的插補結果誤差,為缺失值位置上觀測值與插補值絕對差的累計量結果.

②LSTM 層和線性層組合的解碼結構

LSTM 層和線性層組合的解碼結構如圖4所示.解碼層接收編碼層的輸出序列h,產生插補后的時間序列結果序列y.

LSTM 層和線性層組合的解碼過程數學描述為:

如式(14),解碼層底部是一個標準的LSTM 網絡,該網絡綜合編碼輸出序列h,產生包含更豐富信息輸出狀態序列s={s1,s2,…,sn};如式(15),由于溫度值是連續值,解碼層頂部采用了線性全連接層,輸出插補結果序列y.同式(13),式(16)是解碼層的插補結果誤差.

上述兩種解碼機制的數據插補深度學習網絡的誤差構成是相同的,神經網絡的誤差包括3 部分,即:

4 插值效果評估方法

本文采用多個指標評價不同數據插補方法的性能,評價指標的數值在測試樣本集上計算.包括均方根誤差(RMSE),平均絕對誤差(MAE),平均相對誤差(MRE)和皮爾遜相關系數(PCC),定義如下:

上面指標計算中,xi為所構造測試本中缺失值位置的實際觀測值,yi為缺失值位置對應的插補結果值.PCC評價指標計算如式(21),其中xˉ為樣本中缺失值位置實際觀測值的總體平均,yˉ為缺失值位置插補結果的總體平均.

5 方法結果

作為對比,本文選用了簡單的總體平均插補方法(Mean);基于時間序列分解的卡爾曼插補方法(Kalmanstruct);基于深度學習的BRTS-I 時間序列插補方法;本文設計的兩種編碼-解碼結構的雙向LSTM 網絡插補方法,解碼層為全連接(BiLSTM-FC-I)和解碼層為LSTM(BiLSTM-LSTM-I).

上述方法中,總體平均插補方法可以在整個數據集上計算獲取,而插補方法BRTS-I、BiLSTM-FC-I、BiLSTM-LSTM-I、Kalman-struct 均需要先把整個數據集分為訓練集和測試集,然后在相同的訓練集上進行訓練,在同一測試集上進行精度分析.基于深度學習的BRTS-I、BiLSTM-FC-I、BiLSTM-LSTM-I 雖然模型結構有所不同,但其關鍵參數LSTM 狀態的維度均相同,均取值為108;另外模型的訓練參數也相同,采用了相同的mini batch 參數和優化方法,優化方法均選用Adam,初始學習率取值為0.001;訓練終止策略均為連續10 輪測試精度均不超過訓練過程中測試精度的最優值.深度模型實現是以PyTorch 深度學習框架為基礎,以2020年07月缺失值窗口左側的觀測數據構建訓練集,右側的觀測數據構建測試集.深度學習插補方法構建了兩種訓練樣本,一種缺失值時間窗口設定為30 天,另一種缺失值時間窗口設定為60 天,這兩種訓練樣本缺失值窗口前后連續觀測值均設定為14 天.為了區別這兩種訓練樣本,下面將缺失值時間窗口長度作為了深度學習插補方法的后綴,將各種插補方法的結果總結如表3.

表3 時間序列插補方法結果統計表

從表3 中可見,深度學習方法要明顯優于簡單的總體平均方法.深度學習方法之間的精度也存在較大的差別,圖5 是卡爾曼插值方法,以及各種深度學習方法插值RMSE精度的比較圖.

圖5 不同插值方法RMSE 精度的對比圖

從圖5 中比較BRTS-I、卡爾曼方法、BiLSTM-I 三種方法,本文設計的兩種BiLSTM-I 深度學習氣溫插補方法要優于其他兩種方法;基于LSTM 解碼的BiLSTM-I深度學習方法精度優于基于全連接解碼的BiLSTMI 深度學習方法;BRTS-I 深度學習時間序列插補方法精度最低.

卡爾曼平滑方法時間序列插補方法精度取決于狀態方程是否準確表達了時間序列,Kalman-struct 假定時間序列的趨勢性和季節性的成分可通過基本的線性方程進行擬合.深度學習方法,沒有對時間序列的表達形式進行任何假設,而是通過訓練數據集,自動學習時間序列的準確表達形式.從測試結果看,BiLSTM-I 方法比Kalman-struct 卡爾曼方法,以及BRTS-I 深度學習方法更有可能獲取了時間序列的準確表達,從而獲得了更高的數據插補精度.

本文設計的兩種解碼結構的BiLSTM-I 模型插補精度都高于BRTS-I 模型.BiLSTM-I 模型與BRTS-I模型的區別主要有兩點:首先從模型結構上,BiLSTM-I采用了Encoder-Decoder 結構,而BRTS-I 只相當于BiLSTM-I 模型的Encoder 部分,模型結構上BiLSTM-I有利于充分學習數據的潛在分布規律,從而可以提高數據插補精度;其次模型代價函數有區別,BiLSTM-I和BRTS-I 的代價函數均由3 部分構成[19],前兩部分是相同的,BRTS-I 模型代價函數的第3 部分是前向和后向LSTM 網絡對缺失估計值的差;BiLSTM-I 模型代價函數的第3 部分則為最后的估計值和真實觀測值的差,BiLSTM-I 模型誤差函數對插補結果的評價更直接,模型收斂誤差和插補精度直接對應,從而確保模型收斂時插補誤差也能達到最小.

兩種解碼結構的BiLSTM-I 模型,基于LSTM 解碼的模型插補精度優于全連接解碼結構模型插補精度,這主要是LSTM 解碼時,不但可以利用當前的編碼輸出信息,還可以利用之前的編碼輸出信息;而全連接編碼則只能利用當前的編碼信息,LSTM 解碼方法對編碼信息的利用更為充分.

圖5 中兩種解碼結構的BiLSTM-I 模型,缺值窗口分別為30 天和60 天的測試精度基本一致.深度學習方法應用中模型的泛化能力非常重要,本文問題體現在模型對不同寬度缺失值窗口插補精度是否一致.為了對這一點進行檢驗,我們用缺失值窗口為30 天的模型對缺失值為60 天時間溫度觀測時間序列進行插補,然后用缺失值窗口為60 天的模型對缺失值為30 天時間溫度觀測時間序列進行插補,表4 是這兩種情況下,兩種解碼結構的插補方法結果的精度統計表.

從表4 可見,無論是缺失值窗口為60 天的模型應用到缺失值情況為30 天,還是缺失值窗口為30 天模型應用到缺失值情況為60 天,兩種解碼結構的模型精度的各項指標都非常穩定,這表明這兩種解碼結構的BiLSTM-I 深度學習模型都對不同缺失值窗口有較好的泛化能力.

表4 BiLSTM-I 模型分別應用到缺失值為30、60 天的插值精度統計表

6 總結

本文運用不同深度學習數據插補方法,通過每天低頻的人工溫度觀測數據,獲取完整的高精度半小時頻率溫度觀測數據.本文采用序列-序列的時間序列插補方法,基于編碼-解碼結構的深度學習模型(BiLSTM-I),編碼層采用雙向LSTM-I 網絡,解碼層分別采用LSTM解碼結構與全連接解碼結構,設計了兩種解碼結構的深度學習數據插補模型.

試驗分析結果表明,本文設計的BiLSTM-I 深度學習氣溫插補方法要優于其他方法.基于LSTM 解碼結構的深度學習模型,缺失值時間為30 天的測試集,測試結果精度RMSE為0.47℃;缺失值時間為60 天的測試集,測試結果精度RMSE為0.49℃.基于全連接解碼結構的深度學習模型,缺失值時間為30 天的測試集,測試結果精度RMSE為0.76℃;缺失值時間為60 天的測試集,測試結果精度RMSE為0.70℃.

最后,文章還分析了BiLSTM-I 深度學習插補方法對不同時間溫度缺失長度的適應能力.分別用缺失值時間長度為30 天的訓練模型,對缺失值為60 天的測試集進行插補;以缺失值時間長度為60 天的訓練模型,對缺失值為30 天的測試集進行插補,結果表明兩種解碼結構的深度學習訓練模型對不同的溫度缺失時間長度具有泛化能力.