圖神經網絡驅動的流域洪水預報技術

馬森標

（福建中銳網絡股份有限公司研發中心，福建 福州，350108）

0 引言

流域洪水預報是防汛搶險和防洪系統調度運用的決策依據，傳統的洪水預報采用的機理模型[1]屬于半經驗半物理過程的預報模型，有一定的物理基礎，但模型參數的率定十分復雜和困難[2]，基于數據驅動的模型可根據海量的歷史水文數據挖掘其潛在變化規律，具有比較高的性能和價值。

由于數據驅動模型的優點，目前國內已有很多專家學者將數據驅動的人工智能模型引入到洪水預報中。劉冬英等[3]提出了改進的BP 網絡模型在洪水預報中的應用，取得較高精度的預報結果；金保明等[4]利用RBF 神經網絡預報山區流域洪水的方法，取得了不錯的預報精度效果；在閩江十里庵流量預報中應用BP 神經網絡，預報精度符合要求[5]；李曉麗等[6]提出了不確定支持向量機在洪水預報模型中的應用，獲得了較高的洪水預報精度；馬森標等[7]利用長短時記憶神經網絡（LSTM）附加Attention 機制進行水庫水位預報，取得了較好的性能。雖然以上算法在一定程度上提高了洪水的預報精確度，但是對于輸入的數據特征，尤其是關鍵性的降雨、水位數據特征挖掘明顯不足，關注了時序特征，卻忽略了其本身具有的空間分布特征。同時，降雨量對水位的影響受到流域的地形、土壤、植被等自然因素的作用，表現出了顯著的滯后性、持續性等特點。同樣，受地理位置和流域地形的作用，泄洪量對水位的影響，及上游水位對下游水位的影響同樣存在滯后性、持續性等特點。

圖卷積神經網絡（GCN）[8]是一種在圖結構數據上運行的基于深度學習的擴展方法，可以挖掘符合非歐幾里得空間域的流域數據結構，它能夠從底層數據中提取出網絡節點復雜的關聯性[9]。

基于上述出發點，本研究提出了一個圖指導的時空關聯預報模型（GSCPM，graph-guided spatiotemporal correlation prediction model），用以預報流域水位。具體而言，該模型首先通過多個長短期記憶網絡（LSTM）編碼每個監測點歷史屬性的時間關聯特征，隨后利用挖掘監測點間的地理空間依賴，最后通過一個全連接神經網絡進行回歸預報，可解決流域洪水預報中的時空關系建模和滯后影響問題。

1 模型概述

在本節中，我們將詳細介紹所提出的模型的實現細節，該模型的概述如圖1 所示。

圖1 GSCPM 模型

1.1 問題定義

流域水位預報問題的本質為多變量時間序列預報問題。設流域中有K 個監測點，第i 個監測點的可觀測數據被定義為，其中t 為當前時間，d為過去時間天數，x∈Rn為該觀測點包含的n 個變量，包括水位、降雨量和水庫泄洪量。給定所有監測點的可觀測數據，我們的目標是預報未來1 天的監測點水位

1.2 時間依賴編碼

在獲取K 個監測點的歷史數據后，我們先通過LSTM 編碼其時間依賴關系。由于每個輸入點的特征維度可能不盡相同，并且只使用一個LSTM 可能無法區分不同監測點數據的特異性，我們為每個監測點設置了一個LSTM 網絡，從而得到包含時間依賴的監測點嵌入：

其中，每個LSTM 的具體執行過程如下：

式中：Wf、Wi、WO、WC為可訓練權重矩陣；bf、bi、bO、bC為可訓練偏置矩陣；σ（·）和tanh（·）為非線性激活函數；ht-1為前一時刻的LSTM 單元的輸出；[ ]∵ 為向量連接操作。

1.3 空間依賴編碼

為編碼獲得時空融合特征，我們在時間依賴編碼后，使用GCN 對所得嵌入進行基于地理空間關系的特征聚合。

一個地理信息圖可由G=（V，E，C）描述。其中，V 表示圖的節點集合，其中節點vi∈V 表示第i 個監測點；E表示邊的集合，其中邊eij∈E 表示節點vi和vj中存在關系；C 為節點的特征集合，其中ci∈C 為第i 個節點的屬性值，我們將其設置為第i 個監測點的嵌入表示，即ci=ei。為使地理信息圖能表示實際監測點間的地理關系，我們提出了一個簡單有效的地理信息圖構造方法。具體而言，遍歷每個監測點，將其和流域河道中相鄰的監測點都進行連接，如果該條邊已經存在，則遍歷下一個監測點。圖2 展示了構圖的一個示例。

圖2 拓撲結構圖生成示例

經過上述構圖操作，流域監測點的地理信息被融入圖的拓撲結構之中。隨后，使用一個GCN 聚合節點屬性，從而將圖的拓撲信息融入節點嵌入中。一個圖卷積層的執行過程如下：

式中：D和A分別為圖G 的度矩陣和鄰接矩陣；W和b 分別是可訓練的權重和偏置矩陣；l 和L 分別是GCN 執行的當前層數和最大層數。Zl為第l 層輸出節點嵌入矩陣，同時作為下一層的輸入。第一層GCN 的輸入Z0為節點特征集C 中每個元素拼接而成的矩陣，GCN 的最終輸出ZL被按行拆分K 個向量，代表每個監測點包含時空信息的嵌入表達。

1.4 回歸預測

為進行最終的水位預報，K 個監測點的嵌入向量被按列拼接為新的向量Z。隨后，采用一個全連接網絡進行回歸預測：

式中：WM和bM為可訓練的權重和偏置矩陣；y 為最終預報的水位值集合，包含未來1 天的K 個監測點的水位值。最后，采用均方誤差作為損失函數，用以反向傳播訓練：

式中：Y 為實際的水位值集合。

2 溪源流域水位預報研究

2.1 數據基本情況

本文選取福州溪源流域為研究對象，該流域共包含溪源水庫及下游的侯官和融僑3 個監測點，其中溪源水庫監測站點有泄洪量、降雨量數據，侯官和融僑站點有水位數據。通過章節1.3 所提出的構圖方法，對3 個監測點進行構圖，所得流域的時空示意圖[9]如圖3 所示。

圖3 局部流域時空示意圖

2.2 特征工程

通過機理和數學方法創造的特征，可以考慮到降雨量、泄洪量對水位的影響及上游水位對下游水位的影響都受到流域的地形、土壤、植被諸多自然因素的作用，構造的特征將融合滯后性、持續性特點，根據這些特征訓練的模型具有更好的性能，這也是數據驅動模型優于機理模型的重要原因之一。

2.2.1 滯后特征構造

目前，在機器學習領域，降雨量、泄洪量和上游水位對流域水位的滯后效應，很難通過計算得到，通過構造相應的滯后特征，來模擬這些滯后影響，簡化計算方法，提高建模效率。

一階和二階滯后特征的構造過程如圖4 所示。假設原序列為｛xt-d…，xt-1，xt｝，其m 階段滯后則為其向左平移m 個時間點的序列，表示為｛xt-d-m…，xt-1-m，xt-m｝。在此基礎上，設Rt(i)為t 時刻滯后（i-1）個時間單位的降雨量序列；Ft(i)為t 時刻滯后（i-1）個時間單位的泄洪量序列；Wut(i)為t 時刻滯后（i-1）個時間單位的流域上游水位序列；Wdt(i)為t 時刻滯后（i-1）個時間單位的流域下游水位序列，其中i≥2。

圖4 滯后因素構造過程示例

2.2.2 降雨量滯后特征與水位特征的相關性分析

流域監測點在某一時刻的水位跟之前的降雨量密切相關，相關程度可以通過相關性的強弱來間接表示（見圖5 圖6）。

圖5 雨量曲線圖

圖6 m 階降雨量滯后特征與上、下游水位相關系數曲線圖

從表1 的相關系數可以看出，降雨量對水位的影響具有滯后效應，降雨量對水位的影響是逐漸加強的，隨后逐漸減退，與蓄滿產流的自然規律具有一致性。

表1 m 階降雨量滯后特征與上、下游水位相關系數表

2.2.3 泄洪量滯后特征與水位特征的相關性分析

泄洪量對水位的影響巨大，在大壩某一時刻段的泄洪量，將持續影響流域監測點的水位，這個影響可以通過泄洪量和水位的相關性進行衡量，影響程度可以通過相關性的強弱來間接表示（見圖7 圖8）。

圖7 泄洪量曲線圖

圖8 m 階泄洪量滯后特征與上、下游水位相關系數曲線圖

從表2 的相關系數可以看出，泄洪量對上、下游水位的影響具有持續性，泄洪量對上、下游水位的持續影響是漸漸變弱的，但都具有強相關性，符合水庫泄洪河道水位陡漲的自然客觀規律。

表2 m 階泄洪量滯后特征與上、下游水位相關系數表

2.2.4 上游水位滯后特征與下游水位的相關性分析

在某一時刻流域下游監測點的水位跟該時刻前的上游水位密切相關，相關程度可以通過相關性的強弱來間接表示（見圖9 圖10）。

圖9 上、下游水位曲線圖

圖10 m 階上游水位滯后特征與下游水位相關系數曲線圖

從表3 的相關系數可以看出，下游水位跟上游水位具有較強的相關性，隨著時間的推移，相關性逐漸衰減，符合某一時刻的洪水波在傳播過程中發生衰減演進的客觀規律。

表3 m 階上游水位滯后特征與下游水位相關系數表

2.3 數據集構建及數據預處理

基于章節2.2 中所展示的結論，我們添加相關度最高的滯后特征作為模型的輸入。添加特征后，在溪源水庫監測點包含泄洪量、泄洪量的1 階滯后特征、降雨量和降雨量的6 階滯后特征，共四個特征值；在上游流域監測點，共包含流域水位和流域水位的1 階滯后特征兩個特征值；在下游流域監測點，包含流域水位一個特征值。

由于水庫水文數據包含水位、泄洪量、降雨量等數據，數值差異大，需要進行量綱的統一，以此提高模型的收斂速度和計算精度，因此，本文對序列的每一個特征值進行最大最小值歸一化處理，即：

式中：x 為預處理序列的取值；xmin為序列的最小值；xmax為序列的最大值。

2.4 評價指標

本研究所采用的評價指標包括平均絕對誤差（MAE，Mean Absolute Error）、均方誤差（MSE，Mean Square Error）、平均絕對百分比誤差（MAPE，Mean Absolute Percentage Error）三種，其定義如下所示：

2.5 實現細節

本研究的所有實驗均用Pytorch（V1.7.0）[10]進行實現。為找到最優參數，我們使用網格搜索進行尋優。具體而言，學習率在{0.01，0.005，0.0025}區間進行搜索；正則化系數在{1e-5，1e-7，1e-9}區間搜索；歷史時間長度在{2，4，16}中搜索；嵌入維度和LSTM 單元數在{8，16，32，64，128}區間搜索。此外，每次訓練執行300 個epoch，并使用Adam[11]優化器優化模型參數。

2.6 實驗結果

2.6.1 對比分析

為驗證所提出的GSCPM 的有效性，我們將其與RNN、LSTM 和GRU 網絡進行對比，表4 展示了這4 種方法的對比結果?？梢钥闯?，GSCPM 在三個指標上均明顯優于其他三個方法。在所有方法中，RNN 表現得最差，LSTM 在解決RNN 梯度爆炸問題后得到了性能提升是明顯的，GRU 在簡化LSTM 后也帶來一定的性能提升。然而，這些方法都無法考慮不同監測點間的空間相關性，因而限制了性能。將GSCPM 與LSTM 對比，LSTM 在添加GCN 后帶來的性能提升是明顯的，說明了時空關聯預報對流域水位預報的重要性。

表4 GSCPM 及其對比方法的實驗結果

2.6.2 消融分析

為驗證空間依賴性和所構建的滯后特征的有效性，本研究設置了4 個子方法來與進行對比：

1）GSCPM（-rain）：將降雨量滯后特征刪除后的GSCPM。

2）GSCPM（-flood）：將泄洪量滯后特征刪除后的GSCPM。

3）GSCPM（-up）：將上游水位滯后特征刪除后的GSCPM。

4）GSCPM（-all）：刪除所有滯后特征后的GSCPM。

表5 展示了所有子方法的實驗結果。從表中可看出，GSCPM 展示了最佳的性能和添加滯后特征的有效性。將GSCPM（-rain）、GSCPM（-flood）、GSCPM（-up）和GSCPM（-all）進行對比，它們三者都得到了不同程度的性能提升，說明所提議的三種滯后特征都能單獨提升模型性能。此外，將GSCPM 與GSCPM（-rain）、GSCPM（-flood）、GSCPM（-up）進行對比，GSCPM 仍然具有性能提升，說明滯后特征的組合使用會帶來更大的性能增益。

表5 GSCPM 及其子方法的實驗結果

2.6.3 可視化分析

為進一步理解所提出的GSCPM 模型，我們將該模型及對比模型在測試集上的預報結果進行了可視化，其呈現如圖11～圖18 所示?？梢钥闯?，該模型在向前一步預報的情況下，基本擬合出了未來的變化趨勢，只在少數突變情況下造成較大誤差。值得注意的是，下游水位的真實值波動劇烈，然而GSCPM 仍能捕捉到大部分波動，形成較為精準的預報。相比于其他三個方法的可視化結果，GSCPM 能更好地擬合真實水位曲線，也能更好地預報出突變情況，這也表明為什么我們的方法在眾多評價指標上表現得最好。

圖11 GSCPM 上游預報結果可視化

圖12 GSCPM 下游預報結果可視化

圖13 RNN 上游預報結果可視化

圖14 RNN 下游預報結果可視化

圖15 LSTM 上游預報結果可視化

圖16 LSTM 下游預報結果可視化

圖17 GRU 上游預報結果可視化

圖18 GRU 下游預報結果可視化

3 結論

本研究提出了一種圖指導的時空關聯預報網絡模型，用以針對性地解決流域洪水預報中的時空關系建模問題和滯后影響問題。實踐證明，該模型可以有效降低對流域數據維度的要求，并挖掘流域洪水跟降雨量、泄洪量、上游水位和下游水位之間復雜的時空關聯性。本文提供的是流域2 個節點的洪水預報模型，接下來將研究通過圖神經網絡節點的擴展，通過一個圖神經網絡模型，對整個流域的洪水進行預報的可行性問題，以促進流域洪水預報的廣泛實踐。