孫英軍,唐為昊,王 成 ,李英德

(1.浙江省水文管理中心,浙江 杭州 310009;2.浙江工業(yè)大學(xué) 機(jī)械工程學(xué)院,浙江 杭州 310023)

河流重要水文站點(diǎn)的水位預(yù)測在水資源調(diào)控方面起著重要作用。水文預(yù)測模型在結(jié)構(gòu)上可大體分為確定性模型、黑箱模型和概念性模型[1]。確定性模型以實(shí)際物理過程為基礎(chǔ),用物理理論描述水文系統(tǒng)。以SWAM模型[2]和IHDM模型[3]為代表的分布式流域水文模型是應(yīng)用較為廣泛的確定性模型。黑箱模型是一種由水文數(shù)據(jù)驅(qū)動的預(yù)測模型,將水文系統(tǒng)內(nèi)部過程歸結(jié)為統(tǒng)計關(guān)系,以非線性函數(shù)擬合模型的輸入輸出關(guān)系,常用于缺乏具體產(chǎn)匯流背景資料的水文預(yù)測。單位過程線[4]、支持向量機(jī)[5]和BP神經(jīng)網(wǎng)絡(luò)[6]都屬于黑箱模型。概念性模型介于確定性模型和黑箱模型之間,是一種具有一定物理意義但在某些方面設(shè)置參數(shù)對水文物理過程進(jìn)行簡化的預(yù)測模型。常見的概念性模型有ARNO模型[7]和趙人俊等[8]提出的新安江模型。確定性模型和概念性模型都需要所研究流域具體的水文物理過程資料,例如應(yīng)用新安江模型需要完成蒸發(fā)能力折算系數(shù)等14項(xiàng)參數(shù)的率定[9-10]。雖然一些學(xué)者將搜索算法應(yīng)用于參數(shù)率定,并取得了一定成果,但是就我國廣泛分布的水系而言,考察流域的水文物理過程資料是一項(xiàng)成本較高的工作,且一些流域的水文測量數(shù)據(jù)無法滿足構(gòu)建確定或概念性模型的輸入數(shù)據(jù)需求,因此無法構(gòu)建相關(guān)預(yù)測模型。

黑箱模型作為一種統(tǒng)計方法,預(yù)測過程不依賴于流域的物理過程基礎(chǔ),是一種相對通用的水文預(yù)測模型。以往的黑箱模型常采用一些較為簡單的統(tǒng)計學(xué)習(xí)方法,例如以BP神經(jīng)網(wǎng)絡(luò)為代表的機(jī)器學(xué)習(xí)模型。這些方法將水文數(shù)據(jù)中隱含的物理過程用統(tǒng)計關(guān)系表達(dá),雖然在學(xué)習(xí)數(shù)據(jù)內(nèi)擬合效果良好,但是在外延預(yù)測過程中存在模型泛化能力不足等問題,使得預(yù)測效果大幅下降。近年來在機(jī)器學(xué)習(xí)基礎(chǔ)上發(fā)展而來的深度學(xué)習(xí)一定程度上解決了過擬合等問題,提升了泛化性能,在計算機(jī)視覺、自然語言理解和復(fù)雜過程預(yù)測等諸多場景展現(xiàn)了其強(qiáng)大的性能。水文預(yù)測領(lǐng)域亦有大量采用深度學(xué)習(xí)方法的相關(guān)研究,例如衣學(xué)軍等[11]將小波非線性自回歸網(wǎng)絡(luò)應(yīng)用于渭河流域水位站點(diǎn)的徑流量預(yù)測;Zhao等[12]為了解決高泥沙含量的河道水位不平衡預(yù)測問題,提出了一種混合機(jī)器學(xué)習(xí)框架;Xiang等[13]將Seq2Seq模型應(yīng)用于降雨量-徑流量關(guān)系回歸問題中;紀(jì)國良等[14]將循環(huán)神經(jīng)網(wǎng)絡(luò)應(yīng)用于水庫水位預(yù)測,預(yù)測效果優(yōu)于水動力學(xué)模型;倪漢杰等[15]將小波分析與長短期記憶網(wǎng)絡(luò)相結(jié)合,設(shè)計了一種航道水位智能預(yù)測模型。上述研究提出的黑箱模型均包含由Hochreiter等[16]設(shè)計的長短期記憶網(wǎng)絡(luò)(LSTM)。LSTM獨(dú)特的“門”結(jié)構(gòu)使其能自適應(yīng)地選擇數(shù)據(jù)中需要被“遺忘”以及“記憶”的部分,選擇性保留數(shù)據(jù)的循環(huán)結(jié)構(gòu)使其成為處理時間序列問題的首選深度學(xué)習(xí)模型之一。雖然大量研究驗(yàn)證了LSTM提取時序特征的能力,但其受限于循環(huán)單元結(jié)構(gòu)而難以處理空間特征(在各個循環(huán)步輸入的協(xié)變量)。例如水位預(yù)測任務(wù)中LSTM可很好地利用水位變化的時序特征,卻無法有效利用由降雨量組成的空間特征。為解決上述問題,設(shè)計了一種基于LSTM和一維卷積神經(jīng)網(wǎng)絡(luò)(CNN)的混合神經(jīng)網(wǎng)絡(luò)——卷積-序列到序列(Convolutional neural networks-Sequence to sequence,CNN-Seq2seq)黑箱水文預(yù)測模型。在流域物理過程資料未知的背景下,將設(shè)計的混合神經(jīng)網(wǎng)絡(luò)應(yīng)用于西苕溪流域水位預(yù)測任務(wù)中。模型在測試數(shù)據(jù)集上的預(yù)測準(zhǔn)確性相較于其他黑箱預(yù)測模型更優(yōu),證明了所設(shè)計的混合神經(jīng)網(wǎng)絡(luò)能有效改進(jìn)LSTM在提取空間特征上存在的缺陷,提高了河道水位預(yù)測的準(zhǔn)確性。

1 研究流域與數(shù)據(jù)

1.1 研究流域

研究預(yù)測對象為西苕溪流域內(nèi)河道水位。西苕溪流域如圖1所示。該流域位于浙江省北部安吉縣境內(nèi),東經(jīng)119°14～119°45,北緯30°22～30°45,是太湖的主要支流之一。流域面積約2 200 km2,平均海拔高度266 m。作為該流域內(nèi)的主要河流,西苕溪自西南山區(qū)流向東北部匯入太湖。該地區(qū)屬于亞熱帶季風(fēng)氣候,年平均溫度12.2～15.6 ℃,年平均降水量約1 500 mm,超過70%的降雨發(fā)生在豐水期(4—10月)。為了控制洪水同時兼具水資源調(diào)配等功能,分別于1958年和1972年在西苕溪上游修建了老石坎和賦石兩座大型水庫。

圖1 西苕溪流域地圖Fig.1 Map of the Xitiaoxi watershed

1.2 水文數(shù)據(jù)

西苕溪流域內(nèi)建有30個雨量測量站點(diǎn)、4個河道水位測量站點(diǎn)和2個水庫水位測量站點(diǎn)。采集了上述水文測量站點(diǎn)2014—2019年逐5 min的數(shù)據(jù)記錄,賦石、老石坎兩座水庫在此期間逐1 h的泄洪記錄。對雨量數(shù)據(jù)按照小時累計求和,水位數(shù)據(jù)逐小時求平均,共獲得46 451 h的水文數(shù)據(jù)記錄。按照水利部門相關(guān)預(yù)測需求,預(yù)測模型對“楊家埠”“梅溪”“橫塘村”“港口”4個水位測量站點(diǎn)的水位高度進(jìn)行預(yù)見期為12 h逐小時的預(yù)測。

2 主要方法

2.1 序列到序列結(jié)構(gòu)

長短期記憶網(wǎng)絡(luò)(Long short-term memory,LSTM)是一種為了解決長時依賴問題而改進(jìn)的循環(huán)神經(jīng)網(wǎng)絡(luò),通過在循環(huán)單元中加入幫助神經(jīng)網(wǎng)絡(luò)選擇“遺忘”不重要信息的遺忘門、增強(qiáng)對重要信息“記憶”的記憶門和使神經(jīng)網(wǎng)絡(luò)能綜合考慮輸入信息并生成輸出張量的輸出門。長短期記憶網(wǎng)絡(luò)結(jié)構(gòu)如圖2所示。網(wǎng)絡(luò)結(jié)構(gòu)中的xt表示在時間點(diǎn)t輸入網(wǎng)絡(luò)中的特征向量;ht表示在時間點(diǎn)t的網(wǎng)絡(luò)計算輸出;Ct被稱為循環(huán)單元的隱含狀態(tài),其中包含了在時間點(diǎn)t的網(wǎng)絡(luò)長時記憶信息。每個時間點(diǎn)中,LSTM的循環(huán)單元都需要更新6個參數(shù),具體計算過程為

圖2 長短期記憶網(wǎng)絡(luò)循環(huán)單元結(jié)構(gòu)Fig.2 Basic long short-term memory networkrecursive unit structure

ft=σ(Wf·[ht-1,xt]+bf)

(1)

it=σ(Wi·[ht-1,xt]+bi)

(2)

C′t=tanh(WC·[ht-1,xt]+bC)

(3)

Ct=ft×Ct-1+it×C′t

(4)

ot=σ(Wo·[ht-1,xt]+bo)

(5)

ht=ot×tanh(Ct)

(6)

式中:σ(·),tanh(·)分別表示Sigmoid和tanh激活函數(shù),其計算式為

(7)

(8)

LSTM雖然能解決長時依賴問題,但是仍然存在幾個缺陷:1) 隨著預(yù)測時段時間步的增加,預(yù)測誤差大幅上升;2) 無法考慮降雨預(yù)測數(shù)據(jù)等預(yù)測時段的協(xié)變量信息;3) 訓(xùn)練過程收斂速度較慢。Cho等[17]提出了一種能有效解決上述問題的網(wǎng)絡(luò)結(jié)構(gòu),該結(jié)構(gòu)被稱為序列到序列(Seq2seq)或編碼器解碼器(Encoder-Decoder)。Seq2seq結(jié)構(gòu)包含兩組獨(dú)立的LSTM循環(huán)單元:一組循環(huán)單元將歷史時段的時序特征壓縮為隱含特征(如圖3中的LSTM循環(huán)單元1);另一組循環(huán)單元根據(jù)壓縮后的隱含特征和預(yù)測時段的協(xié)變量特征(如果存在)循環(huán)輸出最終預(yù)測結(jié)果(如圖3中的LSTM循環(huán)單元2)。Salinas等[18]在Seq2seq的基礎(chǔ)上設(shè)計了DeepAR概率預(yù)測模型,該模型在訓(xùn)練過程中加入“teacher forcing”策略,使訓(xùn)練過程的收斂速度有所加快。包含“teacher forcing”策略的Seq2seq模型訓(xùn)練過程和預(yù)測過程如圖3所示。

圖3 包含“teacher forcing”策略的Seq2seq模型訓(xùn)練和預(yù)測過程Fig.3 Seq2seq model training and prediction process with teacher forcing strategy

2.2 一維卷積神經(jīng)網(wǎng)絡(luò)

長短期記憶網(wǎng)絡(luò)雖然在時序數(shù)據(jù)特征提取上展現(xiàn)了強(qiáng)大的性能,但是在空間多特征學(xué)習(xí)上受到了內(nèi)部結(jié)構(gòu)的制約。卷積神經(jīng)網(wǎng)絡(luò)(Convolutional neural networks,CNN)因其空間特征提取能力和相對較少的模型參數(shù),被廣泛應(yīng)用于計算機(jī)視覺領(lǐng)域。根據(jù)卷積核的不同可將卷積神經(jīng)網(wǎng)絡(luò)分為一維、二維和三維卷積,其中一維卷積網(wǎng)絡(luò)(1D-CNN)常被用于提取時序數(shù)據(jù)的空間特征。

一維卷積網(wǎng)絡(luò)采用數(shù)個不同卷積核以一定步長在輸入數(shù)據(jù)的時間維上滑動與輸入數(shù)據(jù)相乘來提取數(shù)據(jù)的隱含特征。相較于全連接網(wǎng)絡(luò),一維卷積網(wǎng)絡(luò)的卷積計算方式具有更少的參數(shù),并能更好地提取出局部特征。卷積運(yùn)算過程為

x′=relu(Wk*x+bk)

(9)

式中:*表示卷積運(yùn)算;Wk表示卷積核權(quán)重矩陣;bk表示卷積偏置項(xiàng);relu(x)=max{x,0}表示非線性激活函數(shù)。

2.3 CNN-Seq2seq模型

將一維卷積網(wǎng)絡(luò)和Seq2seq模型相融合,提出了一種時間序列預(yù)測模型——卷積-序列到序列模型(Convolutional neural networks-Sequence to sequence,CNN-Seq2seq)。結(jié)合水位預(yù)測期長、特征多的特點(diǎn),將長短期記憶網(wǎng)絡(luò)和卷積神經(jīng)網(wǎng)絡(luò)的優(yōu)勢相結(jié)合,由一維卷積網(wǎng)絡(luò)提取以多站點(diǎn)的降雨量數(shù)據(jù)為代表的空間特征,再由Seq2seq模型提取水位數(shù)據(jù)變化情況的時序特征,最后由全連接網(wǎng)絡(luò)根據(jù)所提取到的隱含特征計算得出所需預(yù)測信息。

由于隨機(jī)因素的存在,水位預(yù)測總是存在預(yù)測誤差,因此相較于具體水位高度的預(yù)測,預(yù)測未來水位高度的置信區(qū)間更具現(xiàn)實(shí)意義,相關(guān)水利決策可以根據(jù)預(yù)測水位置信區(qū)間的上下界作出相應(yīng)調(diào)整。按照預(yù)測需求,所設(shè)計的模型輸出未來12 h內(nèi)逐小時河道水位的10%,50%和90%分位數(shù)。預(yù)測模型結(jié)構(gòu)如圖4所示。

圖4 CNN-Seq2seq模型結(jié)構(gòu)Fig.4 Model structure of CNN-Seq2seq

2.4 損失函數(shù)與評價指標(biāo)

深度學(xué)習(xí)模型的迭代優(yōu)化依賴于預(yù)測損失所提供的梯度信息,合理選擇損失函數(shù)可以提升預(yù)測準(zhǔn)確性,增強(qiáng)模型泛化能力。預(yù)測損失通常是無量綱的,因此需要更加符合直覺的評價指標(biāo)以評價模型的實(shí)際預(yù)測效果。區(qū)間預(yù)測包含兩種思路:概率回歸假定被預(yù)測對象服從某種概率分布,構(gòu)建的預(yù)測模型對該分布的參數(shù)進(jìn)行預(yù)測,并由預(yù)測得到的分布參數(shù)計算得出希望獲得的概率分布分位數(shù),該方法一般使用極大似然估計作為損失函數(shù);分位數(shù)回歸不需要假定被預(yù)測對象的概率分布,直接對有限個分位數(shù)進(jìn)行預(yù)測,這種方法使用分位數(shù)損失作為損失函數(shù)。水文預(yù)測具有高度不確定性和不平穩(wěn)性,難以假定所服從的概率分布,因此選擇分位數(shù)回歸進(jìn)行預(yù)測,分位數(shù)損失函數(shù)為

(10)

式中:Ω表示包含M個訓(xùn)練樣本的集合;Q表示輸出分位數(shù)的集合(按照實(shí)際需求,輸出Q={0.1,0.5,0.9});yt,j表示樣本yt在預(yù)見期j的水位高度;n表示預(yù)見期長度;(x)+=max(0,x)。

在進(jìn)行模型評價時,使用標(biāo)準(zhǔn)化分位數(shù)損失(q-risk)對模型在驗(yàn)證集和測試集上的性能進(jìn)行評估[23-24],q-risk代表對分位數(shù)“q”的計算結(jié)果,例如0.1-risk表示對10%分位數(shù)的標(biāo)準(zhǔn)化分位數(shù)損失。具體計算式為

(11)

使用區(qū)間覆蓋率(IC)驗(yàn)證預(yù)測置信區(qū)間的合理性,置信區(qū)間由所預(yù)測得到的0.1和0.9分位數(shù)決定。區(qū)間覆蓋率越接近于兩分位數(shù)數(shù)值之差越優(yōu)。區(qū)間覆蓋率偏大則說明預(yù)測置信區(qū)間過寬,覆蓋率偏小則說明預(yù)測置信區(qū)間過窄。區(qū)間覆蓋率計算式為

(12)

使用納什效率系數(shù)(NSE)、洪峰預(yù)報誤差(QL)和洪峰出現(xiàn)時間誤差(QT)對洪水過程預(yù)測效果進(jìn)行輔助評價,納什效率系數(shù)計算式為

(13)

洪峰水位誤差計算式為

(14)

洪峰出現(xiàn)時間誤差計算式為

(15)

洪峰水位誤差反映了洪水預(yù)測峰值和實(shí)測峰值的差距,洪峰出現(xiàn)時間誤差反映了洪水預(yù)測峰值出現(xiàn)時刻和實(shí)測峰值出現(xiàn)時刻的誤差。QL和QT在數(shù)值上越接近于0則模型對洪峰水位大小和出現(xiàn)時間預(yù)測越準(zhǔn)確。

3 對比實(shí)驗(yàn)

3.1 數(shù)據(jù)預(yù)處理

研究流域的水文數(shù)據(jù)記錄包含大量異常值和空缺值,要使用這些數(shù)據(jù)進(jìn)行模型訓(xùn)練需要剔除異常值并對空缺值進(jìn)行填補(bǔ),無法填補(bǔ)的樣本則棄而不用。剔除異常數(shù)據(jù),對水位數(shù)據(jù)按小時求均值、對雨量數(shù)據(jù)按小時求均值,若水位數(shù)據(jù)被求均值的區(qū)間包含但不全為空缺值,則對剩余數(shù)據(jù)求均值,考慮到短時強(qiáng)降雨等情況的存在,雨量數(shù)據(jù)若被求和的區(qū)間包含空缺值,則將求和結(jié)果置為空缺。一些水位站點(diǎn)在使用時存在基準(zhǔn)調(diào)整的過程,調(diào)整前后的測量數(shù)據(jù)存在整體偏差,選擇最長時間未進(jìn)行基準(zhǔn)調(diào)整的數(shù)據(jù)作為基準(zhǔn),其他時間段的數(shù)據(jù)按照此基準(zhǔn)進(jìn)行數(shù)值調(diào)整。對空缺間隔小于5 h的水位數(shù)據(jù)使用線性插值補(bǔ)全。

水位、雨量和放水量數(shù)據(jù)相差較大,對數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理可避免模型過度考慮某些數(shù)值上較大的特征,標(biāo)準(zhǔn)化過程為

(16)

按照70%,15%和15%的比例將輸入數(shù)據(jù)劃分為訓(xùn)練集、驗(yàn)證集和測試集,則xtrain=[h1:0.7n,f1:0.7n,y1:0.7n],xvalidate=[h0.7n:0.85n,f0.7n:0.85n,y0.7n:0.85n],xtest=[h0.85n:n-12,f0.85n:n-12,y0.85n:n-12]。xtrain表示訓(xùn)練集用于訓(xùn)練構(gòu)建的預(yù)測模型;xvalidate表示驗(yàn)證集用于驗(yàn)證訓(xùn)練好的模型的泛化能力并根據(jù)驗(yàn)證結(jié)果進(jìn)行模型調(diào)優(yōu);xtest表示測試集用于對比不同模型的泛化能力。

3.2 超參數(shù)選擇

要獲得性能優(yōu)異的預(yù)測模型,構(gòu)建理想的模型結(jié)構(gòu)與選擇合適的模型超參數(shù)同樣重要。對于構(gòu)建的深度學(xué)習(xí)算法A,其根本目的在于找到一個函數(shù)f使得給定損失函數(shù)的損失L(x;f)最小,x采樣自預(yù)測數(shù)據(jù)的真實(shí)分布Gx。A的學(xué)習(xí)過程可以看作將真實(shí)分布Gx的有限采樣xtrain映射為預(yù)測模型f。對于參數(shù)化的算法A,總是存在影響預(yù)測模型性能的超參數(shù)組Λ,算法A選用超參數(shù)λ∈Λ在數(shù)據(jù)集xtrain上訓(xùn)練得到預(yù)測模型f的過程可以表示為f=Aλ(xtrain)。例如對于一個全連接神經(jīng)網(wǎng)絡(luò),需要選擇調(diào)整的超參數(shù)包括網(wǎng)絡(luò)層數(shù)、每層網(wǎng)絡(luò)神經(jīng)網(wǎng)絡(luò)的神經(jīng)元數(shù)量、神經(jīng)元選擇的激活函數(shù)以及是否使用參數(shù)正則項(xiàng)、正則項(xiàng)的權(quán)重大小等。超參數(shù)優(yōu)化問題可以表示為

(17)

對式(17)的優(yōu)化是十分困難的:一方面深度學(xué)習(xí)模型的模型訓(xùn)練時間較長,意味著每次計算L(x;Aλ(xtrain))都需要較大的開銷,因此在搜索空間較大時不可能完成對所有超參數(shù)組合空間Λ的遍歷搜索;另一方面L(x;Aλ(xtrain))不存在導(dǎo)數(shù)信息,因此無法使用梯度下降等方法進(jìn)行求解。Bergstra等[25]通過大量實(shí)驗(yàn)證明了隨機(jī)搜索方法在相同時間開銷中對超參數(shù)優(yōu)化的搜索效率高于傳統(tǒng)的網(wǎng)格搜索或人工調(diào)整搜索。隨機(jī)搜索即在Λ中隨機(jī)選擇λ并評估L,在給定時間內(nèi)重復(fù)此過程并記錄最小的L和對應(yīng)的超參數(shù)組合。

采用隨機(jī)搜索方法對模型超參數(shù)進(jìn)行尋優(yōu),以q-risk作為損失函數(shù)計算L(xvalidate;Aλ(xtrain))評估超參數(shù)組合。構(gòu)建的CNN-Seq2seq模型可選的超參數(shù)包括:1) 歷史數(shù)據(jù)長度k(12,24,48,60,72);2) Dropout比例(0.1,0.2,0.3,0.4,0.5,0.6);3) Batchsize大小(32,64,128,256);4) 優(yōu)化器學(xué)習(xí)率(0.000 1,0.001,0.01);5) 最大梯度裁剪(0.01,0.1,1.0,100.0);6) LSTM網(wǎng)絡(luò)神經(jīng)元數(shù)量(16,32,64,128);7) LSTM網(wǎng)絡(luò)層數(shù)(1,3,6);8) CNN網(wǎng)絡(luò)卷積核數(shù)(8,16,32,64);9) CNN網(wǎng)絡(luò)層數(shù)(4,8,16,32);10) 訓(xùn)練迭代數(shù)(10,25,50,75,100,150,200)。

3.3 對比模型

為了驗(yàn)證構(gòu)建的CNN-Seq2seq深度學(xué)習(xí)預(yù)測模型在水位區(qū)間預(yù)測任務(wù)上相較于其他黑箱模型的優(yōu)勢,選擇4種常用的預(yù)測模型、1種預(yù)測基準(zhǔn)(Baseline)以及僅使用Seq2seq結(jié)構(gòu)的預(yù)測模型進(jìn)行比較。1) ARIMA:經(jīng)典的統(tǒng)計預(yù)測方法,被廣泛應(yīng)用于自回歸預(yù)測任務(wù),采用Hyndman-Khandakar算法求解最優(yōu)模型參數(shù);2) NN:全連接神經(jīng)網(wǎng)絡(luò)憑借其簡單的網(wǎng)絡(luò)結(jié)構(gòu)和靈活的輸入輸出關(guān)系被應(yīng)用于多種預(yù)測任務(wù),使用隨機(jī)搜索進(jìn)行隱含層層數(shù)和隱含神經(jīng)元數(shù)量尋優(yōu);3) LSTM:僅使用長短期記憶網(wǎng)絡(luò)提取時間序列的隱含特征,由全連接網(wǎng)絡(luò)根據(jù)提取出的特征輸出預(yù)測結(jié)果[31];4) DeepAR:采用Seq2seq結(jié)構(gòu)、“teacher forcing”訓(xùn)練策略的概率預(yù)測方法,是應(yīng)用較為廣泛的深度學(xué)習(xí)預(yù)測模型之一;5) Seq2seq:在構(gòu)建的CNN-Seq2seq模型中去除一維卷積結(jié)構(gòu),用于驗(yàn)證一維卷積提取輸入數(shù)據(jù)特征的能力,LSTM,DeepAR和Seq2seq模型與構(gòu)建的CNN-Seq2seq模型均具有類似的網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu),因此上述4種模型在隱含層層數(shù)和神經(jīng)元數(shù)量上保持一致,由CNN-Seq2seq模型的超參數(shù)尋優(yōu)結(jié)果確定;6) Baseline:不采用任何預(yù)測模型,以預(yù)測發(fā)生時刻的水位高度作為預(yù)測結(jié)果。Baseline的預(yù)測性能是預(yù)測模型性能的下限,上述模型的預(yù)測效果不應(yīng)差于Baseline。所選用的對比模型具有相同的輸入輸出關(guān)系,均以歷史m小時內(nèi)的水位高度為輸入,輸出未來n小時內(nèi)的水位預(yù)測結(jié)果。然而這些模型在能否考慮協(xié)變量信息、時序信息以及進(jìn)行不確定性預(yù)測(分位數(shù)或概率預(yù)測)上有所區(qū)別(表1)。圖5反映了CNN-Seq2seq模型和Baseline對“梅溪”站點(diǎn)2019年6月20日19時至次日7時水位的預(yù)測效果,圖5中短虛線為Baseline,長虛線和陰影部分為CNN-Seq2seq模型的預(yù)測結(jié)果。預(yù)測事件發(fā)生在圖5中的第48 h。

表1 對比模型所考慮的信息Table 1 Information considered in the comparison model

圖5 Baseline預(yù)測示意圖Fig.5 Prediction baseline schematic

結(jié)合所選擇的模型,整體實(shí)驗(yàn)過程包括4步:1) 整理水文數(shù)據(jù)記錄,對數(shù)據(jù)進(jìn)行清洗和填補(bǔ);2) 將水文數(shù)據(jù)以滑動窗口的方式構(gòu)造為數(shù)據(jù)樣本,并將數(shù)據(jù)樣本劃分為訓(xùn)練集、驗(yàn)證集和測試集;3) 使用訓(xùn)練集完成各個模型的訓(xùn)練,在驗(yàn)證集上調(diào)優(yōu)各個模型的超參數(shù),考慮到深度學(xué)習(xí)訓(xùn)練過程的隨機(jī)性,每個超參數(shù)組合對應(yīng)的模型重復(fù)訓(xùn)練5次,取其在驗(yàn)證集上的平均損失作為評價指標(biāo),完成50次搜索后選擇在驗(yàn)證集上損失最小的模型作為最終訓(xùn)練結(jié)果;4) 根據(jù)評價指標(biāo),對比各個模型在測試集上的預(yù)測效果。

4 實(shí)驗(yàn)結(jié)果與分析

能否準(zhǔn)確預(yù)測洪水過程是考察水位預(yù)測結(jié)果可靠性的關(guān)鍵因素。相較于非洪水過程,洪水水文過程具有不確定性高、非平穩(wěn)和無周期性等特點(diǎn),因此對其進(jìn)行準(zhǔn)確預(yù)測是相對困難的。為了檢驗(yàn)?zāi)Ｐ退活A(yù)測結(jié)果的可靠性,在測試集上分別對含有洪水和非洪水的整體水文過程和僅含洪水過程的數(shù)據(jù)樣本進(jìn)行水位預(yù)測準(zhǔn)確性評估。

對整體水文過程預(yù)測結(jié)果的評估可以反映模型區(qū)分洪水和非洪水水文特征的能力,驗(yàn)證模型是否會將非洪水過程預(yù)測為洪水過程。選擇標(biāo)準(zhǔn)化分位數(shù)損失作為整體水文過程預(yù)測結(jié)果的評估指標(biāo)。與模型所輸出的10%,50%和90%水位分位數(shù)預(yù)測保持一致,評估指標(biāo)包括0.1-risk,0.5-risk和0.9-risk。0.1-risk和0.9-risk兩項(xiàng)指標(biāo)檢驗(yàn)了模型所預(yù)測的區(qū)間上界和下界的準(zhǔn)確性,在這兩項(xiàng)指標(biāo)上較優(yōu)的模型能給出更合理的80%水位置信區(qū)間預(yù)測,0.5-risk反映了模型在水位真值預(yù)測上的相對偏差。

經(jīng)過超參數(shù)優(yōu)化后的模型在測試集上對整體水文過程的預(yù)測結(jié)果如表2所示。表2中粗體為各模型間準(zhǔn)確性最高的預(yù)測結(jié)果;“平均”列表示該模型在4個水位預(yù)測任務(wù)中的平均表現(xiàn);“p”列表示前一列平均值減去最優(yōu)結(jié)果的差與最優(yōu)結(jié)果的比值,反映了與最優(yōu)結(jié)果的相對差距。LSTM,NN和Baseline只能完成河道水位的點(diǎn)預(yù)測(輸出50%分位數(shù)),因此在0.1-risk和0.9-risk兩項(xiàng)指標(biāo)上無數(shù)據(jù)。

表2 標(biāo)準(zhǔn)化分位數(shù)損失(整體水文過程)Table 2 Standardized quantile loss (Overall hydrological process)

整體水文過程的預(yù)測準(zhǔn)確性評估結(jié)果表明:構(gòu)建的CNN-Seq2seq預(yù)測模型除了在“港口”站點(diǎn)90%分位數(shù)的預(yù)測準(zhǔn)確性上略差于Seq2seq模型外(劣化1.71%),其余站點(diǎn)的準(zhǔn)確性均優(yōu)于其他模型(在0.1-risk上優(yōu)于其他模型7%～46%,在0.5-risk上優(yōu)于其他模型8%～49%,在0.9-risk上優(yōu)于其他模型9%～28%)。

采用區(qū)間覆蓋率評估水位預(yù)測置信區(qū)間的預(yù)測準(zhǔn)確性。預(yù)測模型輸出水位預(yù)測的0.1和0.9分位數(shù),因此區(qū)間覆蓋率應(yīng)接近于0.80。模型在測試數(shù)據(jù)集上的區(qū)間覆蓋率統(tǒng)計結(jié)果如表3所示。結(jié)果顯示:ARIMA和Seq2seq模型的預(yù)測置信區(qū)間過窄,而DeepAR模型的置信區(qū)間過寬,CNN-Seq2seq模型的區(qū)間覆蓋率最接近于要求的80%置信區(qū)間。

表3 區(qū)間覆蓋率(整體水文過程)Table 3 Interval coverage (Overall hydrologic process)

為了評估洪水過程的水位預(yù)測準(zhǔn)確性,需在測試數(shù)據(jù)集中劃分洪水過程和非洪水過程。測試數(shù)據(jù)集包含了2018年10月—2019年8月的河道水位數(shù)據(jù),其中包括6次洪水過程,結(jié)果如圖6所示,黑框中為各個站點(diǎn)標(biāo)記的洪水過程。

圖6 測試數(shù)據(jù)集中的洪水時間段Fig.6 Flood period in the test dataset

洪水過程所在時間段分別為2019年2月12日18時—2019年3月10日0時,編號“2019-02-12”;2019年5月26日6時—2019年5月30日12時,編號“2019-05-26”;2019年6月19日18時—2019年6月27日18時,編號“2019-06-19”;2019年6月30日0時—2019年7月7日0時,編號“2019-06-30”;2019年7月9日0時—2019年7月21日0時,編號“2019-07-09”。開始于2019年8月10日的洪水過程因缺少水位、雨量等數(shù)據(jù),不參與評估。

NN,LSTM和Baseline無法進(jìn)行區(qū)間預(yù)測,表2的評估結(jié)果中CNN-Seq2seq,Seq2seq和DeepAR 3種模型的預(yù)測準(zhǔn)確性相近,而ARIMA模型的準(zhǔn)確性較差。為了兼顧水位的區(qū)間預(yù)測功能和水位點(diǎn)預(yù)測具有較高的準(zhǔn)確性,排除在整體水文過程中預(yù)測性能不佳的NN,LSTM,ARIMA和Baseline。采用納什效率系數(shù)、洪峰預(yù)報誤差和洪峰出現(xiàn)時間誤差3項(xiàng)指標(biāo),對CNN-Seq2seq,Seq2seq和DeepAR模型在上述洪水過程預(yù)見期為6 h和12 h的預(yù)測準(zhǔn)確性進(jìn)行評估,結(jié)果如表4～6所示。

表4 納什效率系數(shù)(洪水時間段,NSE)Table 4 Nash efficiency coefficient (flood period, NSE)

表5 實(shí)測洪峰出現(xiàn)時間誤差(QT)Table 5 The actual measurement of the flood peak appeared time error (QT)

表6 實(shí)測洪峰水位高度誤差(QL)Table 6 The actual measurement of the flood peak level error (QL)

對洪水過程的預(yù)測準(zhǔn)確性評估顯示:CNN-Seq2seq模型在絕大部分洪水的評估指標(biāo)上優(yōu)于其他模型,在未達(dá)最優(yōu)的指標(biāo)上與最優(yōu)指標(biāo)仍保持較小的差距。CNN-Seq2seq模型除“楊家埠”站點(diǎn)“2019-05-26”號洪水外的所有洪水過程均能做到在預(yù)見期內(nèi)對洪峰到來進(jìn)行準(zhǔn)確預(yù)測。3個模型關(guān)于“楊家埠”站點(diǎn)“2019-05-26”號洪水的預(yù)測結(jié)果均出現(xiàn)了明顯下降。CNN-Seq2seq模型預(yù)見期6,12 h的水位預(yù)測納什效率系數(shù)分別下降至0.6,0.4,峰現(xiàn)時間誤差分別達(dá)到了5,11 h。經(jīng)過對數(shù)據(jù)集的檢查,發(fā)現(xiàn)“楊家埠”站點(diǎn)洪水中的洪峰相較于其他站點(diǎn)在時間上有明顯偏移,推測該站點(diǎn)所記錄數(shù)據(jù)存在異常,致使水位變化情況無法被準(zhǔn)確預(yù)測。

對上述評估結(jié)果中各水文站點(diǎn)的所有洪水過程求平均,結(jié)果如表7所示。

表7 洪水過程平均評估指標(biāo)Table 7 Flood process average assessment indexes

從洪水過程預(yù)測結(jié)果的平均指標(biāo)來看,構(gòu)建的CNN-Seq2seq模型是所有模型中預(yù)測精度最高的。CNN-Seq2seq模型在所有指標(biāo)上均優(yōu)于其他兩個模型,表明構(gòu)建的CNN-Seq2seq模型在洪水預(yù)測上具有更好的預(yù)測精度。

選擇測試數(shù)據(jù)集中最后一場數(shù)據(jù)記錄完整的洪水過程繪制預(yù)測效果圖。CNN-Seq2seq,Seq2seq和DeepAR模型對編號“2019-07-09”的洪水以預(yù)見期6,12 h進(jìn)行預(yù)測的效果如圖7,8所示。

圖7 “2019-07-09”號洪水6 h預(yù)見期預(yù)測效果Fig.7 Flood prediction effect of “2019-07-09” with 6-hour forecasting horizon

圖8 “2019-07-09”號洪水12 h預(yù)見期預(yù)測效果Fig.8 Flood prediction effect of “2019-07-09” with 12-hour forecasting horizon

為了更清晰地反映各個模型在“2019-07-09”號洪水過程中水位高度點(diǎn)預(yù)測效果的差異,僅繪制CNN-Seq2seq模型的區(qū)間預(yù)測結(jié)果。CNN-Seq2seq和Seq2seq在該洪水過程中的水位高度點(diǎn)預(yù)測準(zhǔn)確性上較為接近,而DeepAR的預(yù)測誤差較大。例如在“楊家埠”站點(diǎn)12 h預(yù)見期第96 h的預(yù)測結(jié)果中,CNN-Seq2seq和Seq2seq模型的絕對預(yù)測誤差接近于0 m,而DeepAR模型的絕對預(yù)測誤差達(dá)到了0.4 m。

水文過程總是存在不確定性,因此無法對水位變化情況作出完全準(zhǔn)確的預(yù)測。河道水位區(qū)間預(yù)測是對水位預(yù)測誤差的補(bǔ)充,在水文狀態(tài)不明確時(例如洪峰前后),模型傾向于輸出范圍更大的水位預(yù)測置信區(qū)間。雖然在一些情況下水位的預(yù)測值和真值會出現(xiàn)較大誤差,但通過預(yù)測的置信區(qū)間仍能得到關(guān)于水位變化情況的估計。在防洪減災(zāi)等應(yīng)用上可將預(yù)測置信區(qū)間的上界考慮為水位預(yù)測結(jié)果,以確保更為穩(wěn)妥的決策方案。

5 結(jié) 論

為了實(shí)現(xiàn)缺乏水文流域物理過程資料的短時河道水位預(yù)測,設(shè)計了一種基于一維卷積網(wǎng)絡(luò)和長短期記憶網(wǎng)絡(luò)的河道水位實(shí)時預(yù)測模型CNN-Seq2seq。將歷史水文數(shù)據(jù)記錄作為預(yù)測模型的訓(xùn)練樣本,采用“teacher forcing”策略完成對模型的訓(xùn)練,使模型能預(yù)測流域下游4個河道水位測量站點(diǎn)未來12 h內(nèi)逐小時的10%,50%和90%河道水位分位數(shù)。對比ARIMA,DeepAR,LSTM,NN和Seq2seq模型的預(yù)測準(zhǔn)確性,實(shí)驗(yàn)結(jié)果表明:CNN-Seq2seq不僅整體水文過程的預(yù)測精度優(yōu)于其他模型,而且對洪水水位的預(yù)測同樣具有較高精度,并優(yōu)于其他模型。實(shí)驗(yàn)驗(yàn)證了CNN-Seq2seq模型在脫離產(chǎn)匯流機(jī)制的情況下仍能對短時河道水位變化過程作出較高精度的預(yù)測,經(jīng)過合理設(shè)計的黑箱水文預(yù)測模型兼具通用性和預(yù)測可靠性。