羅燦坤,劉 昊,黃 鑫,邵 壯

(國(guó)網(wǎng)湖南綜合能源服務(wù)有限公司,湖南 長(zhǎng)沙 410001)

河川徑流作為重要的水文要素之一,由于受到降水、氣溫和人類活動(dòng)等諸多各種因素的影響,顯示出復(fù)雜、隨機(jī)、多維等特征[1]。人類活動(dòng)以及全球氣候變暖等各種因素導(dǎo)致了水資源時(shí)空分布失衡態(tài)勢(shì),從而增加了非線性、非平穩(wěn)性徑流預(yù)測(cè)難度[2]。同時(shí),徑流預(yù)報(bào)不準(zhǔn)確也會(huì)為梯級(jí)水電站優(yōu)化調(diào)度帶來(lái)困難[2-3]。準(zhǔn)確的徑流預(yù)測(cè)對(duì)于解決旱澇災(zāi)害,梯級(jí)水電站優(yōu)化調(diào)度,保證水利設(shè)施正常運(yùn)行及提升其經(jīng)濟(jì)效益都具有重要意義[4]。因此,國(guó)內(nèi)外水文工作者致力于研究徑流特點(diǎn),找尋預(yù)測(cè)結(jié)果更加精確的方法。目前,比較常見(jiàn)的徑流預(yù)測(cè)方法有很多,例如數(shù)理統(tǒng)計(jì)法、物理成因法等[5-6]。隨著計(jì)算科學(xué)技術(shù)的迅速發(fā)展以及關(guān)于徑流預(yù)報(bào)方法的深入研究,許多基于現(xiàn)代智能方法和數(shù)值天氣預(yù)報(bào)的綜合預(yù)報(bào)模型逐漸被提出,主要包括模糊分析[7]、灰色系統(tǒng)理論[8]、混沌理論[9]、小波分析理論[10]、人工神經(jīng)網(wǎng)絡(luò)[11]等,這些方法的應(yīng)用可以有效地提高徑流預(yù)測(cè)結(jié)果的可靠性和精度。

目前,國(guó)內(nèi)外的徑流預(yù)測(cè)大都偏向于通過(guò)建立新的預(yù)測(cè)模型來(lái)提高精度。但是,基于徑流時(shí)間序列的復(fù)雜性分析,單一的模型難以充分反映徑流形成的復(fù)雜性和影響因素的多樣性,導(dǎo)致對(duì)于整個(gè)徑流序列的擬合性較差,預(yù)測(cè)可靠性不高。由此,綜合采用多種預(yù)測(cè)方法的徑流預(yù)報(bào)方式近年來(lái)受到水文研究者的廣泛關(guān)注。綜合近年來(lái)水文工作者的研究成果分析,一般的混合模型多是采用各種分解方法對(duì)徑流數(shù)據(jù)進(jìn)行處理,而后代入預(yù)測(cè)模型進(jìn)行徑流預(yù)測(cè)。從數(shù)據(jù)結(jié)果分析,預(yù)測(cè)精度和穩(wěn)定性得到了較好的提升。近年來(lái),對(duì)于徑流數(shù)據(jù)的時(shí)頻分析因其良好的效果而被研究者逐漸運(yùn)用于水文分析之中[12-13]?；跁r(shí)頻分析,將原始徑流序列分解,處理結(jié)果為多個(gè)不同頻率的子序列,通過(guò)這樣的方式可以使得徑流序列趨于平穩(wěn)化,而后將各子序列分別代入模型,這樣的方式對(duì)提高預(yù)測(cè)的精度效果顯著[14]。張敬平等[15]提出了經(jīng)驗(yàn)?zāi)B(tài)分解與徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)混合徑流預(yù)測(cè)模型,模型預(yù)測(cè)結(jié)果具有較高的精度;Giulia等[16]建立EMD-ANN(Empirical Mode Decomposition-Artificial Neural Network)分解-集成模型,較好地減少了徑流預(yù)測(cè)誤差;劉艷等[17]通過(guò)瑪納斯河的局部信息運(yùn)用集合經(jīng)驗(yàn)?zāi)B(tài)分解(EEMD,Ensemble Empirical Mode Decomposition)預(yù)處理,而后代入整合移動(dòng)平均自回歸(ARIMA,Auto Regressive Integrated Moving Average Model)模型,驗(yàn)證了模型預(yù)測(cè)的準(zhǔn)確性;趙力學(xué)等[18]運(yùn)用變分模態(tài)分解(VMD,Variational Mode Decomposition)和BP神經(jīng)網(wǎng)絡(luò)(Back Propagation Neural Network)組合模型來(lái)水位流量的非線性預(yù)測(cè),取得較好的結(jié)果。上述研究表明,基于時(shí)頻分析的分解技術(shù)對(duì)于提升徑流預(yù)測(cè)精度效果明顯。但是,現(xiàn)有的處理方法大多是對(duì)徑流序列進(jìn)行處理后,采用相同預(yù)測(cè)模型對(duì)子序列進(jìn)行預(yù)測(cè),沒(méi)能很好地考慮子序列的差別性和獨(dú)特性。

因此,提出了一種新的組合預(yù)測(cè)方法,基于不同時(shí)序分析方法,采用VMD方法對(duì)徑流數(shù)據(jù)進(jìn)行分解,利用LSTM方法對(duì)于長(zhǎng)時(shí)間序列具有較好的泛化性和ARMA模型對(duì)于平穩(wěn)序列具有良好的預(yù)測(cè)效果等特點(diǎn)對(duì)不同頻率的模態(tài)進(jìn)行預(yù)測(cè),通過(guò)VMD對(duì)徑流數(shù)據(jù)的處理結(jié)果分析,所得低頻的模態(tài)很好地繼承原始數(shù)據(jù)的時(shí)間特性,選用LSTM神經(jīng)網(wǎng)絡(luò)進(jìn)行處理,而高頻的模態(tài)通過(guò)了平穩(wěn)性檢驗(yàn),采用ARMA模型進(jìn)行預(yù)測(cè),最后將不同模態(tài)的預(yù)測(cè)結(jié)果疊加,得到徑流預(yù)測(cè)數(shù)據(jù)。

1 算法原理

1.1 利用VMD算法將徑流數(shù)據(jù)分解

VMD算法作為一種新的信號(hào)處理方法,是在經(jīng)驗(yàn)?zāi)B(tài)分解(EMD,Empirical Mode Decomposition)的基礎(chǔ)上改進(jìn)所得。作為一種新的非遞歸、自適應(yīng)、準(zhǔn)正交信號(hào)分解方法,VMD算法通過(guò)參數(shù)預(yù)設(shè),將時(shí)間序列數(shù)據(jù)分解為指定個(gè)數(shù)的固有模態(tài)分量(IMF,Intrinsic Mode Function)。VMD算法將數(shù)據(jù)從時(shí)域轉(zhuǎn)化到頻域進(jìn)行分解,不僅可以很好地捕捉時(shí)間序列數(shù)據(jù)的非線性特征,還能避免變量信息重疊,其分解過(guò)程具有很強(qiáng)的魯棒性。

VMD算法將徑流數(shù)據(jù)的分解問(wèn)題轉(zhuǎn)化為變分問(wèn)題,見(jiàn)式(1):

(1)

式中f(t)——徑流數(shù)據(jù);δ(t)——沖激函數(shù);{u}={u1(t),u2(t),…,uk(t)}——分解后得到的K個(gè)模態(tài)分量;{ω}={ω1,ω2,…,ωK}——各模態(tài)分量所對(duì)應(yīng)的中心頻率;*——卷積運(yùn)算符。

然后,引入二次懲罰因子ɑ和拉格朗日算子λ(t),用于求解上述變分問(wèn)題,得到擴(kuò)展的拉格朗日表達(dá)式,見(jiàn)式(2);并采用交替方向乘子法(Alternating Direction Method of Multipliers,ADMM)迭代搜索,迭代后的uk、ωk以及λ見(jiàn)式(3):

L=({uk},{ωk},λ)=

(2)

(3)

式中ω——頻率;γ——噪聲容忍度。

(4)

此時(shí),徑流數(shù)據(jù)被分解為K個(gè)不同頻率的時(shí)間子序列。

1.2 利用LSTM預(yù)測(cè)徑流數(shù)據(jù)的低頻子序列

作為特殊的循環(huán)神經(jīng)網(wǎng)絡(luò)(Recurrent Neural Network,RNN),LSTM(圖1)與傳統(tǒng)的時(shí)間序列算法相比,可以避免長(zhǎng)依賴問(wèn)題,可以解決時(shí)間序列中較長(zhǎng)時(shí)間間隔和延遲的時(shí)間序列問(wèn)題。

圖1 LSTM神經(jīng)元內(nèi)部結(jié)構(gòu)

其具體計(jì)算方法見(jiàn)式(5):

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

(5)

式中ft——t時(shí)刻的遺忘門(mén);σ——sigmoid 函數(shù);bf——遺忘門(mén)的偏置項(xiàng);Wf——權(quán)重矩陣;[ht-1,xt]+bf——2個(gè)向量拼接為1個(gè)新的向量。

當(dāng)前的輸入xt保存到單元狀態(tài)Ct的量由輸入門(mén)決定,見(jiàn)式(6)、(7):

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

(6)

(7)

當(dāng)前時(shí)刻的單元狀態(tài)Ct的計(jì)算見(jiàn)式(8):

(8)

輸出門(mén)的計(jì)算見(jiàn)式(9)、(10):

σt=σ(W0·[ht-1,xt]+bo)

(9)

ht=Ot·tanh(Ct)

(10)

式中Ot——輸出;tanh——激活函數(shù)。

本文所搭建的網(wǎng)絡(luò)中,時(shí)間窗口步長(zhǎng)設(shè)置為1 d,不僅調(diào)用了LSTM神經(jīng)網(wǎng)絡(luò)用來(lái)處理低頻子序列,還設(shè)置了多個(gè)用全連接層搭建的隱藏層和relu激活函數(shù),用來(lái)提取特征解決非線性問(wèn)題,因此更加適合處理時(shí)間序列數(shù)據(jù)。

1.3 ARMA模型

ARMA模型是基于線性自回歸模型(Autoregressive Model,AR)和滑動(dòng)平均模型(Moving Average Model,MA)建立起的平穩(wěn)時(shí)間序列處理模型。ARMA模型可以描述為:

yn=φ1yn-1+φ2yn-2+…+φpyn-p+εn-θ1εn-1-…-θqεn-q

(11)

式中yn——徑流數(shù)據(jù)的預(yù)測(cè)值;{yn-1+yn-2+…+yn-p}——前p個(gè)時(shí)刻的徑流數(shù);{εn,n=0,1,2,…}——白噪聲序列;{φ1,…,φp}——AR模型的系數(shù);{θ1,…,θp}——MA模型的系數(shù),要確定p和q的取值,需要進(jìn)行多次嘗試與檢驗(yàn)。

1.4 VMD-LSTM-ARMA模型構(gòu)建

利用VMD對(duì)徑流進(jìn)行分解,可以得到若干個(gè)復(fù)雜度低的子模態(tài)。對(duì)所得模態(tài)進(jìn)行ADF單位根檢驗(yàn),用以區(qū)別平穩(wěn)與非平穩(wěn)模態(tài)。以往基于時(shí)序分解的組合模型未考慮分解子序列的差異性,因此,本文考量了分解模態(tài)的特征,選取LSTM與ARMA模型,將非平穩(wěn)及平穩(wěn)子模態(tài)分別代入進(jìn)行預(yù)測(cè),最后將所得結(jié)果疊加即得到原始序列預(yù)測(cè)結(jié)果。組合模型建模流程見(jiàn)圖2。

圖2 模型流程

2 模型評(píng)價(jià)

為衡量預(yù)測(cè)模型的精準(zhǔn)度和可靠性,本文選取了較為常用的2種評(píng)價(jià)指標(biāo):平均絕對(duì)誤差(Mean Absolute Error ,MAE )和均方根誤差(Root Mean Square Error ,RMSE)。

(12)

(13)

3 算例分析

耒水屬于長(zhǎng)江流域的湘江水系,干流全長(zhǎng)439 km,流域面積11 905 km2。東江水文站設(shè)立在資興市東江鎮(zhèn)耒水左岸,是國(guó)家重要基本水文站,控制集水面積4 659 km2,站類為一類精度站。

采用東江水文站2020年實(shí)測(cè)徑流數(shù)據(jù),樣本的采樣周期為1 h,選取日期為2020年1月1日至2020年12月31日。原始徑流數(shù)據(jù)見(jiàn)圖3(為了更好地展示數(shù)據(jù)趨勢(shì),只取前240個(gè)點(diǎn))。

圖3 前10 d數(shù)據(jù)展示

從圖中分析,該徑流數(shù)據(jù)具有比較強(qiáng)的非平穩(wěn)性和非線性,采用一種時(shí)間序列預(yù)測(cè)模型(如LSTM、ARMA模型)難以保證預(yù)測(cè)結(jié)果的精度及可靠性。綜合考慮原始數(shù)據(jù)特點(diǎn),本文利用VMD算法對(duì)于信號(hào)在頻域剖分的優(yōu)勢(shì)進(jìn)行原始徑流序列的處理。需要注意的是,模態(tài)數(shù)(K)的取值會(huì)影響VMD分解的效果:當(dāng)K的取值較小,徑流數(shù)據(jù)中的一些重要信息會(huì)丟失,從而影響后續(xù)的預(yù)測(cè)精度;當(dāng)K的取值較大時(shí),鄰近的模態(tài)分量的中心頻率距離變小,會(huì)導(dǎo)致模態(tài)重疊或產(chǎn)生噪聲。因此,K的取值可以通過(guò)觀察來(lái)確定,即對(duì)相鄰模態(tài)分量進(jìn)行相關(guān)性分析或者觀察不同模態(tài)下的中心頻率的分布情況。

此外,為了確定分解后的模態(tài)是否是平穩(wěn)序列,可以對(duì)模態(tài)分量進(jìn)行單位根(ADF)檢驗(yàn),即確定時(shí)間序列中的單位根是否存在。若不存在單位根,則數(shù)據(jù)為平穩(wěn)的時(shí)間序列;若存在單位根,則數(shù)據(jù)為非平穩(wěn)的。VMD參數(shù)設(shè)置見(jiàn)表1。參數(shù)α、τ和ε的設(shè)置參考了文獻(xiàn)[19],K的選擇主要依據(jù)分解后各模態(tài)的頻率分布情況以及預(yù)測(cè)結(jié)果(例如若K=4,可以發(fā)現(xiàn)模態(tài)3、4的頻率分布幾乎沒(méi)有差別,而且預(yù)測(cè)效果比K=3時(shí)的預(yù)測(cè)結(jié)果差)。

表1 VMD算法參數(shù)設(shè)置

利用VMD算法將徑流數(shù)據(jù)分解出3個(gè)IMF,分別命名為IMF1、IMF2、IMF3。然后分別對(duì)3個(gè)模態(tài)進(jìn)行平穩(wěn)性檢驗(yàn)。原假設(shè):檢驗(yàn)的時(shí)間序列具有單位根,函數(shù)的返回值為檢驗(yàn)統(tǒng)計(jì)量、p值和臨界值在1%、5%、10%置信區(qū)間。IMF1的p值為0.09,大于顯著性水平(α=0.05),因此可以認(rèn)為IMF1中的時(shí)間序列沒(méi)有通過(guò)ADF檢驗(yàn),是非平穩(wěn)的時(shí)間序列;然后對(duì)于IMF2、IMF3進(jìn)行平穩(wěn)性檢驗(yàn),p值同為0.001,小于顯著性水平(α=0.05),時(shí)間序列沒(méi)有單位根(特征根在單位圓之內(nèi)),因此可以認(rèn)為IMF2和IMF3中的時(shí)間序列通過(guò)了ADF檢驗(yàn),是平穩(wěn)的時(shí)間序列。

此外,分解后的3個(gè)模態(tài)中,IMF1頻率較低,繼承了原始徑流數(shù)據(jù)的變化趨勢(shì),基本上剔除了隨機(jī)干擾噪聲的影響,具有明顯的時(shí)間序列特征,因此可以采用LSTM神經(jīng)網(wǎng)絡(luò)預(yù)測(cè);而IMF2和IMF3數(shù)據(jù)幅值較小,頻率較大,呈現(xiàn)出一種類似于高斯白噪聲的數(shù)據(jù)類型,屬于平穩(wěn)的時(shí)間序列,可以利用ARMA模型預(yù)測(cè),見(jiàn)圖4—6。因此,需要采用不同的預(yù)測(cè)算法對(duì)數(shù)據(jù)進(jìn)行處理,然后再對(duì)每個(gè)分量預(yù)測(cè)結(jié)果進(jìn)行重構(gòu),則可以得到更加準(zhǔn)確的預(yù)測(cè)效果。

圖4 VMD分解出的第一個(gè)分量

圖5 VMD分解出的第二個(gè)分量

圖6 VMD分解出的第三個(gè)分量

3.1 采用LSTM神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)IMF1

為了更好體現(xiàn)徑流受歷史徑流的影響,本文采用預(yù)測(cè)日前5天的數(shù)據(jù)作為訓(xùn)練樣本。其次,為了滿足水電站制定調(diào)度計(jì)劃的需求,對(duì)預(yù)測(cè)日當(dāng)?shù)臄?shù)據(jù)進(jìn)行預(yù)測(cè),即將前5天的數(shù)據(jù)作為特征,后1天的數(shù)據(jù)作為標(biāo)簽。將IMF1劃分訓(xùn)練集和標(biāo)簽,生成的訓(xùn)練樣本為8 784×120 的矩陣,生成的標(biāo)簽樣本為8 784×24 的矩陣。按照7∶3的比例劃分訓(xùn)練集和測(cè)試集,即訓(xùn)練集的特征樣本為6 149×120的矩陣,標(biāo)簽樣本6 149×24的矩陣;測(cè)試集的特征樣本為2 635×120的矩陣,標(biāo)簽樣本為2 635×24的矩陣。系統(tǒng)環(huán)境為Windows10-64位系統(tǒng),Python版本3.8.5,IDE為PyCharm,利用 TensorFlow的前端 Keras搭建 LSTM 神經(jīng)網(wǎng)絡(luò)進(jìn)行時(shí)間序列分析。LSTM模型參數(shù):隱藏節(jié)點(diǎn)24,激活函數(shù)為relu,迭代次數(shù)為50。在建立LSTM神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)后,為了防止過(guò)擬合,在全連接層后面加入Dropout層,從而提高泛化能力。優(yōu)化器選用Adam,損失函數(shù)選用mae。當(dāng)訓(xùn)練集和測(cè)試集的訓(xùn)練誤差逐漸減少且趨于穩(wěn)定時(shí),可以認(rèn)為L(zhǎng)STM神經(jīng)網(wǎng)絡(luò)擬合完畢。

利用LSTM算法對(duì)IMF1進(jìn)行預(yù)測(cè),損失函數(shù)隨迭代次數(shù)的變化情況見(jiàn)圖7?？梢钥闯?隨著迭代次數(shù)的增加,準(zhǔn)確率逐漸升高,損失函數(shù)逐漸減小。

圖7 LSTM損失函數(shù)隨迭代次數(shù)變化曲線

3.2 采用ARMA模型預(yù)測(cè)IMF2和IMF3

ARMA模型的輸入數(shù)據(jù)為IMF2和IMF3的前3×24個(gè)點(diǎn)的值預(yù)測(cè)接下來(lái)24個(gè)點(diǎn)的值。經(jīng)過(guò)多次嘗試和檢驗(yàn),本案例的ARMA模型的參數(shù)為:p=1,q=0。

3.3 預(yù)測(cè)結(jié)果對(duì)比分析

為了得到徑流數(shù)據(jù)的預(yù)測(cè)結(jié)果,將IMF1、IMF2與IMF3各自的預(yù)測(cè)結(jié)果進(jìn)行疊加。為了驗(yàn)證本模型對(duì)徑流數(shù)據(jù)預(yù)測(cè)的可靠性和精準(zhǔn)度,本文引入另外幾種常見(jiàn)的徑流預(yù)測(cè)模型進(jìn)行比對(duì),包括LSTM神經(jīng)網(wǎng)絡(luò)、ARMA模型、隨機(jī)森林(Random Forest,RF)和支持向量機(jī)(Support Vector Machine,SVM),采用MAE和RMSE進(jìn)行定量分析,見(jiàn)表2、圖8。

表2 預(yù)測(cè)評(píng)價(jià)指標(biāo)對(duì)比 單位:m3/s

a)實(shí)際-組合模型

由表2可知,與其他4種算法的預(yù)測(cè)結(jié)果相比,本文提出的VMD-LSTM-ARMA組合模型的MAE與RMSE最小,徑流預(yù)測(cè)結(jié)果更為精準(zhǔn)。此外,為了直觀體現(xiàn)模型預(yù)測(cè)效果及差異,本文將原始徑流數(shù)據(jù)(取其中400個(gè)點(diǎn))與VMD-LSTM-ARMA組合模型及其他4種對(duì)比模型預(yù)測(cè)值可視化,結(jié)果見(jiàn)圖8,4種對(duì)比模型的預(yù)測(cè)值與原始數(shù)據(jù)的擬合度較低,與本文提出的組合模型相差較大,說(shuō)明VMD-LSTM-ARMA模型具有較高的預(yù)測(cè)精度,是一個(gè)有效的模型。

4 結(jié)語(yǔ)

從徑流時(shí)間序列的非線性、復(fù)雜性特征角度考慮,采用時(shí)頻分析方法,基于VMD分解,充分考慮徑流數(shù)據(jù)分解所得模態(tài)的特點(diǎn),考量LSTM對(duì)非平穩(wěn)低頻模態(tài)和ARMA模型對(duì)平穩(wěn)高頻模態(tài)的優(yōu)勢(shì),建立起了VMD-LSTM-ARMA模型,以東江水文站逐小時(shí)流量數(shù)據(jù)為例進(jìn)行預(yù)測(cè),并與LSTM、BP、隨機(jī)森林和SVM等常見(jiàn)單一模型進(jìn)行結(jié)果比對(duì),性能分析,得出如下結(jié)論:①采用VMD分解方法可以有效提取徑流時(shí)間序列中不同頻率的固有信息,通過(guò)分解,較好地降低了原始數(shù)據(jù)的非平穩(wěn)性,降低了預(yù)測(cè)難度;②針對(duì)分解所得模態(tài)的差異性,分別選取LSTM與ARMA模型進(jìn)行預(yù)測(cè),對(duì)于提升預(yù)測(cè)精度效果明顯,通過(guò)與比對(duì)模型比較,與直接將原始徑流序列代入模型預(yù)測(cè)相比,基于分解的組合模型預(yù)測(cè)效果更理想;③本文提出的組合模型將VMD、LSTM與ARMA模型的優(yōu)點(diǎn)有機(jī)結(jié)合起來(lái),在預(yù)測(cè)時(shí)間序列的表現(xiàn)要優(yōu)于單一模型,具有較高的精度和可靠性。