嚴(yán)春華，張 玥，姜中清，梅 杰，高志發(fā)

(江蘇省水利勘測(cè)設(shè)計(jì)研究院有限公司，江蘇 揚(yáng)州 225000)

0 引 言

徑流量預(yù)測(cè)是洪水預(yù)報(bào)的關(guān)鍵環(huán)節(jié)，對(duì)于水庫(kù)大壩、堤防、水閘等水利工程的安全穩(wěn)定及正常運(yùn)行有著至關(guān)重要的作用，關(guān)系到流域人民生命財(cái)產(chǎn)安全[1]。隨著人工神經(jīng)網(wǎng)絡(luò)和深度學(xué)習(xí)的發(fā)展，基于數(shù)學(xué)模型的徑流量預(yù)測(cè)的研究成果也越來越多，如王軍等基于TCN-Attention模型對(duì)黃河徑流量進(jìn)行預(yù)測(cè)[2]，李代華構(gòu)建了基于PCA-SHO-SVM和PCA-SHO-BP模型進(jìn)行徑流預(yù)測(cè)[3]。這些組合模型可以較好地實(shí)現(xiàn)對(duì)徑流量的預(yù)測(cè)，但忽略了徑流量的時(shí)間序列的滯后性和誤差累積，從而導(dǎo)致誤差的增加和收斂速度的下降。因此，需要綜合考慮組合模型的特性，進(jìn)一步深度研究算法的優(yōu)缺點(diǎn)，進(jìn)而構(gòu)建高精度和穩(wěn)定性的徑流量預(yù)測(cè)數(shù)學(xué)模型。

為了解決徑流量時(shí)間序列的滯后性，本文引入偏自相關(guān)函數(shù)進(jìn)行研究[4]。偏自相關(guān)函數(shù)在處理時(shí)間序列數(shù)據(jù)時(shí)可以研究監(jiān)測(cè)序列和本身的k階滯后序列兩者的相關(guān)性，從而得出徑流量序列對(duì)某時(shí)刻的影響作用程度。在智能學(xué)習(xí)算法中，目前以神經(jīng)網(wǎng)絡(luò)模型(如CNN、RNN和LSTM等網(wǎng)絡(luò))的預(yù)測(cè)能耐較好，其缺點(diǎn)為運(yùn)算時(shí)間長(zhǎng)，如果訓(xùn)練和預(yù)測(cè)集過大的話，迭代速度會(huì)大大減慢、收斂時(shí)間長(zhǎng)。為了解決和克服常用神經(jīng)網(wǎng)絡(luò)算法的缺點(diǎn)，2014年Cho等[5]優(yōu)化了LSTM算法模型，構(gòu)建了GRU算法模型，借助設(shè)置隱藏層提高訓(xùn)練學(xué)習(xí)的速度。GRU相比于LSTM網(wǎng)絡(luò)優(yōu)化了網(wǎng)絡(luò)參數(shù)，計(jì)算收斂速度更快，并且兩個(gè)算法預(yù)測(cè)精度相差不大，所以本文選擇GRU網(wǎng)絡(luò)構(gòu)建徑流量預(yù)測(cè)模型。為了進(jìn)一步地提高徑流量預(yù)測(cè)模型的精度，需要解決預(yù)測(cè)過程中的誤差修正問題。由于徑流量的非穩(wěn)定性，任何模型模擬出的徑流量均具有一定偏差，并且誤差的積累是一種非線性問題，而1998年Huang等[6]提出了一種可以較為高效地解決非線性與非平穩(wěn)信號(hào)時(shí)頻分解難題的方法——經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)，因而采用EMD方法對(duì)殘差部分進(jìn)行預(yù)處理。該處理時(shí)間序列的方法被廣泛運(yùn)用到很多研究領(lǐng)域并有不錯(cuò)的結(jié)果精度，如測(cè)繪空間定位[7]、風(fēng)力預(yù)測(cè)[8]、軍事雷達(dá)[9]、工程爆破[10]研究等領(lǐng)域。因此，本文提出了以GRU為主要建模手段的徑流量預(yù)測(cè)方法，采用偏自相關(guān)方法確定時(shí)滯步長(zhǎng)的前處理和EMD進(jìn)行誤差修正的后處理，分別解決徑流量預(yù)測(cè)中存在的收斂速度慢、時(shí)序滯后性和誤差積累的問題。通過工程案例以及和其他預(yù)測(cè)模型進(jìn)行對(duì)比，驗(yàn)證了該徑流量預(yù)測(cè)模型的精度和可靠性。

1 徑流量預(yù)測(cè)模型構(gòu)建的基本理論及實(shí)現(xiàn)流程

1.1 徑流量時(shí)間序列相關(guān)性確定方法

為了尋取徑流量時(shí)間序列的滯后特性，采用偏自相關(guān)分析方法，并通過分析結(jié)果分析輸入量與預(yù)測(cè)值之間的關(guān)系[11]。偏自相關(guān)函數(shù)主要研究yt和yt+k之間的相關(guān)性，而排除了它們之間元素yt+1，yt+2，…，yt+k-1對(duì)yt的影響，其表達(dá)式如下[12]

(1)

1.2 GRU時(shí)序預(yù)測(cè)模型基本原理

GRU模型保留了LSTM模型的基本網(wǎng)絡(luò)構(gòu)造，同時(shí)增加隱藏層，模型的基本架構(gòu)由細(xì)胞狀態(tài)、更新門和重置門組成，該算法模型架構(gòu)如圖1所示。

圖1 GRU網(wǎng)絡(luò)模型的架構(gòu)

GRU算法的隱藏層的輸入由上一個(gè)時(shí)間節(jié)點(diǎn)的隱藏層狀態(tài)輸出h(t-1)和當(dāng)前時(shí)間序列實(shí)際輸入x(t)組成。模型內(nèi)部的更新門z(t)和重置門r(t)是優(yōu)化時(shí)序預(yù)測(cè)模型的內(nèi)部結(jié)果和減少參數(shù)數(shù)量的關(guān)鍵組成部分[13]。

更新門的作用是對(duì)當(dāng)前時(shí)間節(jié)點(diǎn)輸入值的梯度判斷，與LSTM算法輸入門的功能基本一致，能夠極大程度克服梯度消失的難題，其門結(jié)構(gòu)為

z(t)=σ(Wz·h(t-1)+Uzx(t)+bz)

(2)

式中，σ為激活函數(shù)，一般是sigmoid函數(shù)；Wz、Uz分別為更新門結(jié)構(gòu)中的權(quán)重；bz為更新門結(jié)構(gòu)中的偏置。

重置門的作用是判斷前一時(shí)間節(jié)點(diǎn)隱藏層狀態(tài)的h(t-1)對(duì)當(dāng)前時(shí)間節(jié)點(diǎn)的x(t)的影響程度。如果程度不大，從而打開重置門，其門結(jié)構(gòu)為

r(t)=σ(Wr·h(t-1)+Urx(t)+br)

(3)

式中，σ為激活函數(shù)，一般是sigmoid函數(shù)；Wr、Ur分別為重置門結(jié)構(gòu)中的權(quán)重；br為重置門結(jié)構(gòu)中的偏置。

確定更新隱藏層中的細(xì)胞狀態(tài)。細(xì)胞狀態(tài)的更新由更新門與重置門的輸出確定，其中隱藏層狀態(tài)中最右側(cè)的h～(t)為

h～(t)=tanh(W(r(t)°h(t-1))+Ux(t)+b)

(4)

式中，tanh為該部分的激活函數(shù)；W、U分別為更新h～(t)結(jié)構(gòu)中的權(quán)重；b為更新h～(t)結(jié)構(gòu)中的偏置。

當(dāng)前時(shí)間節(jié)點(diǎn)隱藏層的輸出h(t)，結(jié)合更新門的輸出更新細(xì)胞狀態(tài)得到當(dāng)前時(shí)間節(jié)點(diǎn)隱藏層的輸出h(t)

h(t)=(1-z(t))°h(t-1)+z(t)°h～(t)

(5)

計(jì)算預(yù)測(cè)值的輸出是對(duì)當(dāng)前時(shí)間節(jié)點(diǎn)的預(yù)測(cè)值輸出，則輸出門的輸出為

y(t)=σ(Vh(t)+c)

(6)

式中，V、c分別為權(quán)重參數(shù)及偏置。

1.3 EMD基本理論

EMD是一種自適應(yīng)分解法，把具有時(shí)序特性的信號(hào)分解成不同頻率特征的本征模態(tài)函數(shù)(Intrinsic Mode Function，IMF)和一個(gè)長(zhǎng)期趨勢(shì)項(xiàng)。EMD法根據(jù)原始信號(hào)的上包絡(luò)線和下包絡(luò)線的均值求解瞬時(shí)均值，進(jìn)而實(shí)現(xiàn)對(duì)原始信號(hào)的特征分解[8]。給定原始數(shù)據(jù)信號(hào)xt，EMD法的實(shí)現(xiàn)基本分為5個(gè)步驟。

第1步，求解原始信號(hào)xt的上包絡(luò)線u(t)與下包絡(luò)線l(t)。計(jì)算xt中所有的局部極值點(diǎn)，再采用三次樣條插值法對(duì)極值點(diǎn)插值處理，從而得到原始信號(hào)的上下包絡(luò)線。

第2步，計(jì)算瞬時(shí)均值。把第1步得到的上包絡(luò)線u(t)和下包絡(luò)線l(t)取均值，即為瞬時(shí)均值m(t)，

(7)

第3步，把原始信號(hào)xt減掉瞬時(shí)均值m(t)，可得新序列f1(t)，

f1(t)=xt-m(t)

(8)

第4步，得到IMF分量一定要有兩個(gè)條件：其一為序列f1(t)的局部極值點(diǎn)與過原點(diǎn)的個(gè)數(shù)的差值必須是0或1；其二為任何時(shí)刻，瞬時(shí)均值都為0?；跇?gòu)成IMF的條件對(duì)序列f1(t)進(jìn)行判別，如果達(dá)到IMF特征條件，序列f1(t)為IMF分量；如果序列f1(t)達(dá)不到條件，把f1(t)作為原始信號(hào)序列重復(fù)步驟1到3，一直分解得到的序列滿足IMF的特征條件。而且設(shè)置c1(t)=f1(t)，得到信號(hào)xt和c1(t)的差值，從而一階殘差量為r1(t)，再把殘差量作為原始信號(hào)重復(fù)步驟1到3。

第5步，當(dāng)有第1個(gè)IMF分量f1(t)，使c1(t)=f1(t)，把原始信號(hào)xt減掉c1(t)可以得到1個(gè)殘差量r1(t)，為

r1(t)=xt-c1(t)

(9)

上述5步完成了EMD法的首次信號(hào)分解。進(jìn)一步地，把殘差量r1(t)重新重復(fù)上述步驟從而可以分解得到若干個(gè)IMF分量。當(dāng)?shù)趎個(gè)殘差量rn(t)只有單個(gè)極值點(diǎn)或者為單調(diào)序列時(shí)分解停止，EMD分解流程結(jié)束。

原始信號(hào)xt經(jīng)過分解可以得到若干IMF分量為

(10)

式中，fi(t)為分解的第i個(gè)IMF分量；rn(t)為第n次分解的殘差量；n為原始信號(hào)的分解次數(shù)，同時(shí)是分解得到的IMF分量數(shù)。

1.4 徑流量預(yù)測(cè)模型實(shí)現(xiàn)步驟

綜上理論，構(gòu)建基于GRU的徑流量預(yù)測(cè)模型，該模型的基本實(shí)現(xiàn)流程如下：

第1步，確定徑流量序列的訓(xùn)練學(xué)習(xí)集和預(yù)測(cè)集。

第2步，通過偏自相關(guān)性方法確定徑流量序列的時(shí)滯步長(zhǎng)。

第3步，基于GRU算法進(jìn)行學(xué)習(xí)訓(xùn)練，輸出徑流量預(yù)測(cè)值。

第4步，提取殘差并基于EMD分解后利用GRU對(duì)各分量進(jìn)行殘差修正。

第5步，并與其他常用算法對(duì)比，驗(yàn)證本文模型預(yù)測(cè)結(jié)果的精度和模型的可靠性。

2 案例分析

2.1 徑流量監(jiān)測(cè)數(shù)據(jù)提取

將某河流的徑流量采集數(shù)據(jù)作為本章方法的建模對(duì)象，研究時(shí)段為2015年1月1日至2019年12月31日，共計(jì)5 a，該時(shí)段的徑流讀取頻次為每日一次，其徑流量原始數(shù)據(jù)如圖2所示。在進(jìn)行建模分析時(shí)，首先需要對(duì)徑流量數(shù)據(jù)進(jìn)行數(shù)據(jù)集的劃分，本文以2015年1月1日～2018年12月31日4年的日徑流量數(shù)據(jù)作為訓(xùn)練樣本集，2019年1月1日～2019年12月31日1年的日徑流量數(shù)據(jù)作為測(cè)試樣本，以分析提出的徑流量預(yù)測(cè)模型的有效性。

圖2 徑流量原始變化數(shù)據(jù)

2.2 徑流量時(shí)間序列相關(guān)性分析

由理論部分可知，對(duì)河流徑流量建立GRU模型進(jìn)行徑流時(shí)序特性分析之前，首先應(yīng)當(dāng)確定時(shí)間序列的時(shí)滯步長(zhǎng)，這對(duì)于正確建立歷史徑流量與當(dāng)前徑流量之間的關(guān)系具有重要意義。若是確定的時(shí)滯步長(zhǎng)過小，則不足以有效建立徑流量的時(shí)序關(guān)系；若是時(shí)滯步長(zhǎng)過大，則將增加過多無用信息，GRU模型在學(xué)習(xí)徑流量時(shí)序關(guān)系時(shí)的效果也將變差。因此，本文按照理論部分提出的偏自相關(guān)分析方法進(jìn)行時(shí)滯步長(zhǎng)確定，若將徑流量的當(dāng)前時(shí)刻確定為t，本案例中分析了t-1，t-2，t-3，…，t-12共計(jì)12個(gè)時(shí)間點(diǎn)處徑流量與當(dāng)前時(shí)刻徑流量之間的偏自相關(guān)性，如圖3所示。

圖3 相關(guān)性分析

由圖3可以看出，當(dāng)徑流量時(shí)間序列的時(shí)滯步長(zhǎng)小于8時(shí)，歷史徑流數(shù)據(jù)與當(dāng)前徑流量之間的相關(guān)性是逐漸下降的；當(dāng)時(shí)滯步長(zhǎng)的取值大于8時(shí)，歷史徑流數(shù)據(jù)與當(dāng)前徑流量之間的相關(guān)性趨于平穩(wěn)且處于較低相關(guān)程度。因此，在后面的時(shí)序建模中，將時(shí)滯步長(zhǎng)的取值確定為8。

2.3 GRU建模分析和殘差提取

通過2.2節(jié)設(shè)定的徑流量時(shí)間序列時(shí)滯步長(zhǎng)，可以對(duì)原始徑流量時(shí)間序列進(jìn)行GRU建模。在此過程中，以2015年1月1日～2018年12月31日4年的日徑流量數(shù)據(jù)作為訓(xùn)練樣本集，2019年1月1日～2019年12月31日1年的日徑流量數(shù)據(jù)作為測(cè)試樣本，即訓(xùn)練集樣本數(shù)和測(cè)試集樣本數(shù)的比例為4∶1。首先通過訓(xùn)練集樣本訓(xùn)練GRU模型，然后將訓(xùn)練完成的模型應(yīng)用在測(cè)試集上，則徑流量GRU模型在測(cè)試集上的分析結(jié)果如圖4所示。

圖4 基于GRU模型的徑流量預(yù)測(cè)結(jié)果

通過圖4可以看出，該方法能夠大體模擬徑流量的變化趨勢(shì)，然而由于徑流量波動(dòng)性較大，導(dǎo)致分析結(jié)果不能夠準(zhǔn)確估計(jì)其波動(dòng)性，因此，本文提出了殘差修正的方法對(duì)訓(xùn)練集部分?jǐn)M合殘差進(jìn)行進(jìn)一步分析，訓(xùn)練集部分?jǐn)M合結(jié)果如圖5所示，由此得到的徑流量殘差如圖6所示。

圖5 基于GRU模型的訓(xùn)練集徑流量擬合結(jié)果

2.4 徑流量預(yù)測(cè)結(jié)果的殘差修正

由理論部分可知，本小節(jié)主要針對(duì)GRU分析結(jié)果殘差進(jìn)行修正，目的是得到更為準(zhǔn)確的徑流量預(yù)測(cè)結(jié)果。首先采用EMD方法對(duì)上圖中的徑流量殘差進(jìn)行分解，從而降低徑流量殘差序列的復(fù)雜性以提升建模分析結(jié)果的精度，EMD分析結(jié)果如圖7所示。

圖7 基于EMD的徑流量殘差分解結(jié)果

由圖7可以看出，徑流量殘差經(jīng)過EMD處理后，原本具有很強(qiáng)不規(guī)則特性的殘差變成了具有不同變化頻率的殘差分量，這對(duì)于殘差建模分析提供了基礎(chǔ)。與直接對(duì)徑流量包含殘差的建模分析相比可知，通過對(duì)具有較低頻且規(guī)律較為明顯的分量進(jìn)行建模，更容易學(xué)習(xí)其包含的徑流量信息，殘差修正的效果更明顯。通過對(duì)徑流量各殘差(IMF1～I(xiàn)MF11)分別建立GRU模型，并對(duì)各GRU模型的預(yù)測(cè)結(jié)果進(jìn)行綜合得到最終的徑流量預(yù)測(cè)結(jié)果，如圖8所示。

圖8 基于殘差修正的徑流量預(yù)測(cè)結(jié)果

對(duì)比圖4和圖8徑流量測(cè)試階段的預(yù)測(cè)結(jié)果，很顯然在經(jīng)過殘差修正后，徑流量的預(yù)測(cè)結(jié)果更加接近徑流量真實(shí)值，分析結(jié)果的波動(dòng)性較之殘差修正之前也更接近真實(shí)值，驗(yàn)證了本文提出的徑流量預(yù)測(cè)方法的有效性。

2.5 模型對(duì)比分析

上文已分析得到基于偏自相關(guān)分析和殘差修正的徑流量預(yù)測(cè)模型的預(yù)測(cè)結(jié)果，為說明該方法在徑流預(yù)測(cè)領(lǐng)域的優(yōu)越性和可靠性，選取ARIMA和LSTM這2種常用的時(shí)間關(guān)系分析模型進(jìn)行對(duì)比說明。其中在建立上述2個(gè)模型時(shí)，徑流數(shù)據(jù)集的劃分同上述案例。3種模型的預(yù)測(cè)結(jié)果對(duì)比如圖9所示。

圖9 3種方法徑流量預(yù)測(cè)結(jié)果

由圖9可以看出，ARIMA的預(yù)測(cè)結(jié)果較為平滑，且與徑流量的真實(shí)值相差甚遠(yuǎn)；而具有時(shí)序?qū)W習(xí)能力的LSTM方法相對(duì)來說具有一定模擬原始徑流量波動(dòng)程度的能力，相比之下文中方法預(yù)測(cè)結(jié)果與真實(shí)的徑流量變化過程一致，徑流預(yù)測(cè)精度高。上述3種方法分析得到的徑流量預(yù)測(cè)結(jié)果殘差對(duì)比圖，如圖10所示。由圖10可以看出，ARIMA的預(yù)測(cè)誤差明顯是最大的，甚至在某些時(shí)間內(nèi)其誤差高達(dá)3 000 m3/s，說明該模型并不適合應(yīng)用于徑流量預(yù)測(cè)問題；LSTM相對(duì)誤差較小，而文中提出的方法徑流量預(yù)測(cè)結(jié)果誤差最小。3種方法徑流量預(yù)測(cè)結(jié)果評(píng)價(jià)指標(biāo)如表1所示。從表1可以看出，文中方法分析得到的徑流量預(yù)測(cè)結(jié)果對(duì)應(yīng)的RMSE值和MAE值是3種方法中最小的，其R2值也明顯優(yōu)于另外2種方法。以上關(guān)于3種方法的對(duì)比分析說明了本文所提出方法的優(yōu)越性，這為洪水預(yù)報(bào)提供了高可靠度的技術(shù)支撐。

表1 3種方法徑流量預(yù)測(cè)結(jié)果評(píng)價(jià)指標(biāo)

圖10 徑流量預(yù)測(cè)結(jié)果殘差對(duì)比

3 結(jié) 論

本文通過偏自相關(guān)方法的前處理和EMD誤差的后處理，構(gòu)建了基于GRU模型的徑流量預(yù)測(cè)模型，其結(jié)論如下：

(1)基于GRU模型的徑流量預(yù)測(cè)結(jié)果更接近于徑流量真實(shí)值，計(jì)算結(jié)果的波動(dòng)性更接近真實(shí)值，預(yù)測(cè)誤差更小，說明該模型的工程實(shí)用性高。

(2)通過與ARIMA和LSTM預(yù)測(cè)模型的對(duì)比表明該模型的優(yōu)越性，驗(yàn)證了該徑流量預(yù)測(cè)模型的有效性和可靠度，對(duì)洪水預(yù)報(bào)提供了有力的技術(shù)支撐。

(3)由于本文采用的是固定時(shí)滯步長(zhǎng)，而河流在不同時(shí)間下的步長(zhǎng)是變化的，可以進(jìn)一步研究動(dòng)態(tài)滯后步長(zhǎng)對(duì)預(yù)測(cè)值的影響。