羅燦坤,劉 昊,黃 鑫,邵 壯
(國網(wǎng)湖南綜合能源服務(wù)有限公司,湖南 長沙 410001)
河川徑流作為重要的水文要素之一,由于受到降水、氣溫和人類活動等諸多各種因素的影響,顯示出復(fù)雜、隨機(jī)、多維等特征[1]。人類活動以及全球氣候變暖等各種因素導(dǎo)致了水資源時空分布失衡態(tài)勢,從而增加了非線性、非平穩(wěn)性徑流預(yù)測難度[2]。同時,徑流預(yù)報不準(zhǔn)確也會為梯級水電站優(yōu)化調(diào)度帶來困難[2-3]。準(zhǔn)確的徑流預(yù)測對于解決旱澇災(zāi)害,梯級水電站優(yōu)化調(diào)度,保證水利設(shè)施正常運(yùn)行及提升其經(jīng)濟(jì)效益都具有重要意義[4]。因此,國內(nèi)外水文工作者致力于研究徑流特點,找尋預(yù)測結(jié)果更加精確的方法。目前,比較常見的徑流預(yù)測方法有很多,例如數(shù)理統(tǒng)計法、物理成因法等[5-6]。隨著計算科學(xué)技術(shù)的迅速發(fā)展以及關(guān)于徑流預(yù)報方法的深入研究,許多基于現(xiàn)代智能方法和數(shù)值天氣預(yù)報的綜合預(yù)報模型逐漸被提出,主要包括模糊分析[7]、灰色系統(tǒng)理論[8]、混沌理論[9]、小波分析理論[10]、人工神經(jīng)網(wǎng)絡(luò)[11]等,這些方法的應(yīng)用可以有效地提高徑流預(yù)測結(jié)果的可靠性和精度。
目前,國內(nèi)外的徑流預(yù)測大都偏向于通過建立新的預(yù)測模型來提高精度。但是,基于徑流時間序列的復(fù)雜性分析,單一的模型難以充分反映徑流形成的復(fù)雜性和影響因素的多樣性,導(dǎo)致對于整個徑流序列的擬合性較差,預(yù)測可靠性不高。由此,綜合采用多種預(yù)測方法的徑流預(yù)報方式近年來受到水文研究者的廣泛關(guān)注。綜合近年來水文工作者的研究成果分析,一般的混合模型多是采用各種分解方法對徑流數(shù)據(jù)進(jìn)行處理,而后代入預(yù)測模型進(jìn)行徑流預(yù)測。從數(shù)據(jù)結(jié)果分析,預(yù)測精度和穩(wěn)定性得到了較好的提升。近年來,對于徑流數(shù)據(jù)的時頻分析因其良好的效果而被研究者逐漸運(yùn)用于水文分析之中[12-13]?;跁r頻分析,將原始徑流序列分解,處理結(jié)果為多個不同頻率的子序列,通過這樣的方式可以使得徑流序列趨于平穩(wěn)化,而后將各子序列分別代入模型,這樣的方式對提高預(yù)測的精度效果顯著[14]。張敬平等[15]提出了經(jīng)驗?zāi)B(tài)分解與徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)混合徑流預(yù)測模型,模型預(yù)測結(jié)果具有較高的精度;Giulia等[16]建立EMD-ANN(Empirical Mode Decomposition-Artificial Neural Network)分解-集成模型,較好地減少了徑流預(yù)測誤差;劉艷等[17]通過瑪納斯河的局部信息運(yùn)用集合經(jīng)驗?zāi)B(tài)分解(EEMD,Ensemble Empirical Mode Decomposition)預(yù)處理,而后代入整合移動平均自回歸(ARIMA,Auto Regressive Integrated Moving Average Model)模型,驗證了模型預(yù)測的準(zhǔn)確性;趙力學(xué)等[18]運(yùn)用變分模態(tài)分解(VMD,Variational Mode Decomposition)和BP神經(jīng)網(wǎng)絡(luò)(Back Propagation Neural Network)組合模型來水位流量的非線性預(yù)測,取得較好的結(jié)果。上述研究表明,基于時頻分析的分解技術(shù)對于提升徑流預(yù)測精度效果明顯。但是,現(xiàn)有的處理方法大多是對徑流序列進(jìn)行處理后,采用相同預(yù)測模型對子序列進(jìn)行預(yù)測,沒能很好地考慮子序列的差別性和獨特性。
因此,提出了一種新的組合預(yù)測方法,基于不同時序分析方法,采用VMD方法對徑流數(shù)據(jù)進(jìn)行分解,利用LSTM方法對于長時間序列具有較好的泛化性和ARMA模型對于平穩(wěn)序列具有良好的預(yù)測效果等特點對不同頻率的模態(tài)進(jìn)行預(yù)測,通過VMD對徑流數(shù)據(jù)的處理結(jié)果分析,所得低頻的模態(tài)很好地繼承原始數(shù)據(jù)的時間特性,選用LSTM神經(jīng)網(wǎng)絡(luò)進(jìn)行處理,而高頻的模態(tài)通過了平穩(wěn)性檢驗,采用ARMA模型進(jìn)行預(yù)測,最后將不同模態(tài)的預(yù)測結(jié)果疊加,得到徑流預(yù)測數(shù)據(jù)。
VMD算法作為一種新的信號處理方法,是在經(jīng)驗?zāi)B(tài)分解(EMD,Empirical Mode Decomposition)的基礎(chǔ)上改進(jìn)所得。作為一種新的非遞歸、自適應(yīng)、準(zhǔn)正交信號分解方法,VMD算法通過參數(shù)預(yù)設(shè),將時間序列數(shù)據(jù)分解為指定個數(shù)的固有模態(tài)分量(IMF,Intrinsic Mode Function)。VMD算法將數(shù)據(jù)從時域轉(zhuǎn)化到頻域進(jìn)行分解,不僅可以很好地捕捉時間序列數(shù)據(jù)的非線性特征,還能避免變量信息重疊,其分解過程具有很強(qiáng)的魯棒性。
VMD算法將徑流數(shù)據(jù)的分解問題轉(zhuǎn)化為變分問題,見式(1):
(1)
式中f(t)——徑流數(shù)據(jù);δ(t)——沖激函數(shù);{u}={u1(t),u2(t),…,uk(t)}——分解后得到的K個模態(tài)分量;{ω}={ω1,ω2,…,ωK}——各模態(tài)分量所對應(yīng)的中心頻率;*——卷積運(yùn)算符。
然后,引入二次懲罰因子ɑ和拉格朗日算子λ(t),用于求解上述變分問題,得到擴(kuò)展的拉格朗日表達(dá)式,見式(2);并采用交替方向乘子法(Alternating Direction Method of Multipliers,ADMM)迭代搜索,迭代后的uk、ωk以及λ見式(3):
L=({uk},{ωk},λ)=
(2)
(3)
式中ω——頻率;γ——噪聲容忍度。
(4)
此時,徑流數(shù)據(jù)被分解為K個不同頻率的時間子序列。
作為特殊的循環(huán)神經(jīng)網(wǎng)絡(luò)(Recurrent Neural Network,RNN),LSTM(圖1)與傳統(tǒng)的時間序列算法相比,可以避免長依賴問題,可以解決時間序列中較長時間間隔和延遲的時間序列問題。
圖1 LSTM神經(jīng)元內(nèi)部結(jié)構(gòu)
其具體計算方法見式(5):
ft=σ(Wf·[ht-1,xt]+bf)
(5)
式中ft——t時刻的遺忘門;σ——sigmoid 函數(shù);bf——遺忘門的偏置項;Wf——權(quán)重矩陣;[ht-1,xt]+bf——2個向量拼接為1個新的向量。
當(dāng)前的輸入xt保存到單元狀態(tài)Ct的量由輸入門決定,見式(6)、(7):
it=σ(Wi·[ht-1,xt]+bi)
(6)
(7)
當(dāng)前時刻的單元狀態(tài)Ct的計算見式(8):
(8)
輸出門的計算見式(9)、(10):
σt=σ(W0·[ht-1,xt]+bo)
(9)
ht=Ot·tanh(Ct)
(10)
式中Ot——輸出;tanh——激活函數(shù)。
本文所搭建的網(wǎng)絡(luò)中,時間窗口步長設(shè)置為1 d,不僅調(diào)用了LSTM神經(jīng)網(wǎng)絡(luò)用來處理低頻子序列,還設(shè)置了多個用全連接層搭建的隱藏層和relu激活函數(shù),用來提取特征解決非線性問題,因此更加適合處理時間序列數(shù)據(jù)。
ARMA模型是基于線性自回歸模型(Autoregressive Model,AR)和滑動平均模型(Moving Average Model,MA)建立起的平穩(wěn)時間序列處理模型。ARMA模型可以描述為:
yn=φ1yn-1+φ2yn-2+…+φpyn-p+εn-θ1εn-1-…-θqεn-q
(11)
式中yn——徑流數(shù)據(jù)的預(yù)測值;{yn-1+yn-2+…+yn-p}——前p個時刻的徑流數(shù);{εn,n=0,1,2,…}——白噪聲序列;{φ1,…,φp}——AR模型的系數(shù);{θ1,…,θp}——MA模型的系數(shù),要確定p和q的取值,需要進(jìn)行多次嘗試與檢驗。
利用VMD對徑流進(jìn)行分解,可以得到若干個復(fù)雜度低的子模態(tài)。對所得模態(tài)進(jìn)行ADF單位根檢驗,用以區(qū)別平穩(wěn)與非平穩(wěn)模態(tài)。以往基于時序分解的組合模型未考慮分解子序列的差異性,因此,本文考量了分解模態(tài)的特征,選取LSTM與ARMA模型,將非平穩(wěn)及平穩(wěn)子模態(tài)分別代入進(jìn)行預(yù)測,最后將所得結(jié)果疊加即得到原始序列預(yù)測結(jié)果。組合模型建模流程見圖2。
圖2 模型流程
為衡量預(yù)測模型的精準(zhǔn)度和可靠性,本文選取了較為常用的2種評價指標(biāo):平均絕對誤差(Mean Absolute Error ,MAE )和均方根誤差(Root Mean Square Error ,RMSE)。
(12)
(13)
耒水屬于長江流域的湘江水系,干流全長439 km,流域面積11 905 km2。東江水文站設(shè)立在資興市東江鎮(zhèn)耒水左岸,是國家重要基本水文站,控制集水面積4 659 km2,站類為一類精度站。
采用東江水文站2020年實測徑流數(shù)據(jù),樣本的采樣周期為1 h,選取日期為2020年1月1日至2020年12月31日。原始徑流數(shù)據(jù)見圖3(為了更好地展示數(shù)據(jù)趨勢,只取前240個點)。
圖3 前10 d數(shù)據(jù)展示
從圖中分析,該徑流數(shù)據(jù)具有比較強(qiáng)的非平穩(wěn)性和非線性,采用一種時間序列預(yù)測模型(如LSTM、ARMA模型)難以保證預(yù)測結(jié)果的精度及可靠性。綜合考慮原始數(shù)據(jù)特點,本文利用VMD算法對于信號在頻域剖分的優(yōu)勢進(jìn)行原始徑流序列的處理。需要注意的是,模態(tài)數(shù)(K)的取值會影響VMD分解的效果:當(dāng)K的取值較小,徑流數(shù)據(jù)中的一些重要信息會丟失,從而影響后續(xù)的預(yù)測精度;當(dāng)K的取值較大時,鄰近的模態(tài)分量的中心頻率距離變小,會導(dǎo)致模態(tài)重疊或產(chǎn)生噪聲。因此,K的取值可以通過觀察來確定,即對相鄰模態(tài)分量進(jìn)行相關(guān)性分析或者觀察不同模態(tài)下的中心頻率的分布情況。
此外,為了確定分解后的模態(tài)是否是平穩(wěn)序列,可以對模態(tài)分量進(jìn)行單位根(ADF)檢驗,即確定時間序列中的單位根是否存在。若不存在單位根,則數(shù)據(jù)為平穩(wěn)的時間序列;若存在單位根,則數(shù)據(jù)為非平穩(wěn)的。VMD參數(shù)設(shè)置見表1。參數(shù)α、τ和ε的設(shè)置參考了文獻(xiàn)[19],K的選擇主要依據(jù)分解后各模態(tài)的頻率分布情況以及預(yù)測結(jié)果(例如若K=4,可以發(fā)現(xiàn)模態(tài)3、4的頻率分布幾乎沒有差別,而且預(yù)測效果比K=3時的預(yù)測結(jié)果差)。
表1 VMD算法參數(shù)設(shè)置
利用VMD算法將徑流數(shù)據(jù)分解出3個IMF,分別命名為IMF1、IMF2、IMF3。然后分別對3個模態(tài)進(jìn)行平穩(wěn)性檢驗。原假設(shè):檢驗的時間序列具有單位根,函數(shù)的返回值為檢驗統(tǒng)計量、p值和臨界值在1%、5%、10%置信區(qū)間。IMF1的p值為0.09,大于顯著性水平(α=0.05),因此可以認(rèn)為IMF1中的時間序列沒有通過ADF檢驗,是非平穩(wěn)的時間序列;然后對于IMF2、IMF3進(jìn)行平穩(wěn)性檢驗,p值同為0.001,小于顯著性水平(α=0.05),時間序列沒有單位根(特征根在單位圓之內(nèi)),因此可以認(rèn)為IMF2和IMF3中的時間序列通過了ADF檢驗,是平穩(wěn)的時間序列。
此外,分解后的3個模態(tài)中,IMF1頻率較低,繼承了原始徑流數(shù)據(jù)的變化趨勢,基本上剔除了隨機(jī)干擾噪聲的影響,具有明顯的時間序列特征,因此可以采用LSTM神經(jīng)網(wǎng)絡(luò)預(yù)測;而IMF2和IMF3數(shù)據(jù)幅值較小,頻率較大,呈現(xiàn)出一種類似于高斯白噪聲的數(shù)據(jù)類型,屬于平穩(wěn)的時間序列,可以利用ARMA模型預(yù)測,見圖4—6。因此,需要采用不同的預(yù)測算法對數(shù)據(jù)進(jìn)行處理,然后再對每個分量預(yù)測結(jié)果進(jìn)行重構(gòu),則可以得到更加準(zhǔn)確的預(yù)測效果。
圖4 VMD分解出的第一個分量
圖5 VMD分解出的第二個分量
圖6 VMD分解出的第三個分量
為了更好體現(xiàn)徑流受歷史徑流的影響,本文采用預(yù)測日前5天的數(shù)據(jù)作為訓(xùn)練樣本。其次,為了滿足水電站制定調(diào)度計劃的需求,對預(yù)測日當(dāng)?shù)臄?shù)據(jù)進(jìn)行預(yù)測,即將前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)練集和測試集,即訓(xùn)練集的特征樣本為6 149×120的矩陣,標(biāo)簽樣本6 149×24的矩陣;測試集的特征樣本為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)行時間序列分析。LSTM模型參數(shù):隱藏節(jié)點24,激活函數(shù)為relu,迭代次數(shù)為50。在建立LSTM神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)后,為了防止過擬合,在全連接層后面加入Dropout層,從而提高泛化能力。優(yōu)化器選用Adam,損失函數(shù)選用mae。當(dāng)訓(xùn)練集和測試集的訓(xùn)練誤差逐漸減少且趨于穩(wěn)定時,可以認(rèn)為LSTM神經(jīng)網(wǎng)絡(luò)擬合完畢。
利用LSTM算法對IMF1進(jìn)行預(yù)測,損失函數(shù)隨迭代次數(shù)的變化情況見圖7??梢钥闯?隨著迭代次數(shù)的增加,準(zhǔn)確率逐漸升高,損失函數(shù)逐漸減小。
圖7 LSTM損失函數(shù)隨迭代次數(shù)變化曲線
ARMA模型的輸入數(shù)據(jù)為IMF2和IMF3的前3×24個點的值預(yù)測接下來24個點的值。經(jīng)過多次嘗試和檢驗,本案例的ARMA模型的參數(shù)為:p=1,q=0。
為了得到徑流數(shù)據(jù)的預(yù)測結(jié)果,將IMF1、IMF2與IMF3各自的預(yù)測結(jié)果進(jìn)行疊加。為了驗證本模型對徑流數(shù)據(jù)預(yù)測的可靠性和精準(zhǔn)度,本文引入另外幾種常見的徑流預(yù)測模型進(jìn)行比對,包括LSTM神經(jīng)網(wǎng)絡(luò)、ARMA模型、隨機(jī)森林(Random Forest,RF)和支持向量機(jī)(Support Vector Machine,SVM),采用MAE和RMSE進(jìn)行定量分析,見表2、圖8。
表2 預(yù)測評價指標(biāo)對比 單位:m3/s
a)實際-組合模型
由表2可知,與其他4種算法的預(yù)測結(jié)果相比,本文提出的VMD-LSTM-ARMA組合模型的MAE與RMSE最小,徑流預(yù)測結(jié)果更為精準(zhǔn)。此外,為了直觀體現(xiàn)模型預(yù)測效果及差異,本文將原始徑流數(shù)據(jù)(取其中400個點)與VMD-LSTM-ARMA組合模型及其他4種對比模型預(yù)測值可視化,結(jié)果見圖8,4種對比模型的預(yù)測值與原始數(shù)據(jù)的擬合度較低,與本文提出的組合模型相差較大,說明VMD-LSTM-ARMA模型具有較高的預(yù)測精度,是一個有效的模型。
從徑流時間序列的非線性、復(fù)雜性特征角度考慮,采用時頻分析方法,基于VMD分解,充分考慮徑流數(shù)據(jù)分解所得模態(tài)的特點,考量LSTM對非平穩(wěn)低頻模態(tài)和ARMA模型對平穩(wěn)高頻模態(tài)的優(yōu)勢,建立起了VMD-LSTM-ARMA模型,以東江水文站逐小時流量數(shù)據(jù)為例進(jìn)行預(yù)測,并與LSTM、BP、隨機(jī)森林和SVM等常見單一模型進(jìn)行結(jié)果比對,性能分析,得出如下結(jié)論:①采用VMD分解方法可以有效提取徑流時間序列中不同頻率的固有信息,通過分解,較好地降低了原始數(shù)據(jù)的非平穩(wěn)性,降低了預(yù)測難度;②針對分解所得模態(tài)的差異性,分別選取LSTM與ARMA模型進(jìn)行預(yù)測,對于提升預(yù)測精度效果明顯,通過與比對模型比較,與直接將原始徑流序列代入模型預(yù)測相比,基于分解的組合模型預(yù)測效果更理想;③本文提出的組合模型將VMD、LSTM與ARMA模型的優(yōu)點有機(jī)結(jié)合起來,在預(yù)測時間序列的表現(xiàn)要優(yōu)于單一模型,具有較高的精度和可靠性。