王亦斌，孫 濤，梁雪春,謝海洋

(1.南水北調(diào)東線江蘇水源有限責(zé)任公司，江蘇 南京 210019； 2.南京工業(yè)大學(xué)電氣工程與控制科學(xué)學(xué)院，江蘇 南京 211816)

準(zhǔn)確的水位和水量等水文時(shí)間序列預(yù)測(cè)是水資源管理的重要依據(jù)，對(duì)于調(diào)水監(jiān)測(cè)有著重要的作用。傳統(tǒng)的水文時(shí)間序列預(yù)測(cè)方法有回歸分析、時(shí)間序列分析等，還有模糊識(shí)別法、灰色模型、小波神經(jīng)網(wǎng)絡(luò)分析[1-5]。隨著機(jī)器學(xué)習(xí)方法的興起，很多學(xué)者使用機(jī)器學(xué)習(xí)方法進(jìn)行建模,相對(duì)于傳統(tǒng)的統(tǒng)計(jì)預(yù)測(cè)模型，機(jī)器學(xué)習(xí)方法具有更好的適應(yīng)性。面對(duì)數(shù)據(jù)量化不確定性問(wèn)題，Porter等[6]提出一種人工神經(jīng)網(wǎng)絡(luò)(artificial neural networks,ANN)方法來(lái)預(yù)測(cè)地下水水位；針對(duì)在臺(tái)風(fēng)期間受潮汐影響水位波動(dòng)較大問(wèn)題，Wei[7]提出了小波-支持向量機(jī)算法來(lái)預(yù)測(cè)每小時(shí)的水位數(shù)據(jù)，為臺(tái)風(fēng)襲擊期間的水位預(yù)測(cè)問(wèn)題提供了實(shí)用的解決方案；Moosavi等[8]提出了小波-自適應(yīng)神經(jīng)模糊推理系統(tǒng)(wavelet-adaptive neuro-fuzzy inference system，wavelet-ANIS)方法，比較了ANN、Wavelet-ANN、ANIS和Wavelet-ANIS在地下水位預(yù)報(bào)上的表現(xiàn)；Zhong等[9]分析了長(zhǎng)江潮汐水位預(yù)測(cè)的復(fù)雜性，在獲取了長(zhǎng)江下游(安慶、蕪湖、南京)的每日水位的歷史數(shù)據(jù)后，提出使用人工智能-Kalman濾波方法進(jìn)行水位預(yù)報(bào)；Wang等[10]分析南水北調(diào)中線陶岔渠首的水位數(shù)據(jù)情況，提出了一種改進(jìn)Cuckoo-搜索梯度增強(qiáng)樹算法(CCS-gradient boosting decision tree)進(jìn)行水位預(yù)報(bào)。近年來(lái)，隨著深度學(xué)習(xí)技術(shù)的發(fā)展，深度神經(jīng)網(wǎng)絡(luò)也逐漸被應(yīng)用于時(shí)間序列預(yù)測(cè)分析中，尤其是循環(huán)神經(jīng)網(wǎng)絡(luò)(recurrent neural networks，RNN)在時(shí)序數(shù)據(jù)處理上表現(xiàn)出很強(qiáng)的適應(yīng)性。Chang等[11]使用循環(huán)神經(jīng)網(wǎng)絡(luò)作城市水位預(yù)測(cè)和洪水預(yù)報(bào)。張軒等[12]為提高秦淮河流域東山站水位預(yù)報(bào)的精度，基于BP神經(jīng)網(wǎng)絡(luò)算法建立經(jīng)驗(yàn)預(yù)報(bào)模型，分別根據(jù)降雨歷時(shí)、起漲水位兩種模式對(duì)水位漲幅進(jìn)行預(yù)報(bào)。常見的RNN等算法存在梯度消失或者梯度爆炸的問(wèn)題，而長(zhǎng)短期記憶神經(jīng)網(wǎng)絡(luò)(long short-term memory networks，LSTM) 是一種引入了遺忘門、輸入門和輸出門3種門控單元和記憶細(xì)胞的特殊形式的RNN，這種特殊的結(jié)構(gòu)可以緩解梯度消失或者梯度爆炸問(wèn)題[13]。 LSTM通過(guò)門控結(jié)構(gòu)調(diào)整隱藏狀態(tài)中信息的流動(dòng)，輸入門可以識(shí)別需要添加到長(zhǎng)期項(xiàng)的狀態(tài)，遺忘門可以控制需要消除的長(zhǎng)期項(xiàng)狀態(tài)，可以使得模型能學(xué)習(xí)何時(shí)忘記以前的隱藏狀態(tài)，何時(shí)根據(jù)新的信息更新隱藏狀態(tài)，從而能夠成功地捕獲時(shí)間序列中的長(zhǎng)期模式，因此，被廣泛應(yīng)用于交通、金融和經(jīng)濟(jì)序列預(yù)測(cè)中[14-16]。LSTM模型在水文時(shí)間序列預(yù)報(bào)和模擬方面也受到了關(guān)注。為了有效提高流域內(nèi)站點(diǎn)徑流量預(yù)測(cè)的精度,馮銳[17]運(yùn)用LSTM技術(shù)建立九龍江流域站點(diǎn)徑流序列預(yù)測(cè)模型和數(shù)據(jù)模擬模型。劉亞新等[18]綜合使用梯度下降與Broyden-Fletcher-Goldfarb-Shanno算法來(lái)訓(xùn)練LSTM模型，并將模型用于葛洲壩水電站水位預(yù)測(cè)。以入湖流量和長(zhǎng)江干流流量作為輸入條件，郭燕等[19]采用LSTM模型預(yù)測(cè)鄱陽(yáng)湖的水位變化。閆佰忠等[20]基于地下水位、蒸發(fā)量、降水量、氣溫、氣壓、相對(duì)濕度、日照時(shí)長(zhǎng)、開采量等多個(gè)指標(biāo)構(gòu)建多變量的LSTM地下水位預(yù)測(cè)模型。通常水文時(shí)間序列是非平穩(wěn)且非線性的，為提高預(yù)測(cè)精度，陳旭等[21]提出了基于經(jīng)驗(yàn)?zāi)Ｊ椒纸?empirical mode decomposition，EMD)的AR預(yù)測(cè)模型用于水文站的年徑流量預(yù)測(cè)?；谏鲜鲅芯?，采用EMD和LSTM相結(jié)合的算法(EMD-LSTM)進(jìn)行河流水量水位預(yù)測(cè),首先采用信號(hào)處理中常用的經(jīng)驗(yàn)?zāi)Ｊ椒纸夥椒▽⑺臅r(shí)間序列數(shù)據(jù)分解為不同時(shí)間尺度的固有模態(tài)函數(shù)，然后對(duì)EMD分解的每個(gè)特征序列使用LSTM模型進(jìn)行預(yù)測(cè)，最后疊加每個(gè)特征序列預(yù)測(cè)值得到預(yù)測(cè)結(jié)果。

南水北調(diào)工程對(duì)改善水環(huán)境、促進(jìn)區(qū)域經(jīng)濟(jì)和社會(huì)的可持續(xù)發(fā)展具有重要的戰(zhàn)略意義。在我國(guó)南水北調(diào)工程的重要河段中，準(zhǔn)確預(yù)測(cè)河流水位水量，對(duì)水庫(kù)調(diào)度、防洪、發(fā)電、灌溉等具有重要意義。本文擬使用EMD和LSTM模型對(duì)水位、水速和水量進(jìn)行預(yù)測(cè)。以南水北調(diào)工程某河流2017年11月到 2018年6月中每隔1 h瞬時(shí)流量、流速和水深監(jiān)測(cè)數(shù)據(jù)為研究對(duì)象，采用中值濾波對(duì)數(shù)據(jù)進(jìn)行預(yù)處理，基于EMD-LSTM進(jìn)行時(shí)間序列建模。

1 研究方法

1.1 EMD

EMD是Huang等[22]提出的自適應(yīng)信號(hào)時(shí)頻處理方法，它不需要預(yù)先設(shè)置任何基函數(shù)，根據(jù)數(shù)據(jù)本身的時(shí)間尺度特征將信號(hào)分解成有限數(shù)目的固有模態(tài)函數(shù)。對(duì)時(shí)間序列{Xt}的EMD分解過(guò)程如下:

步驟1：求出時(shí)間序列{Xt}的極大值和極小值，由三次樣條插值函數(shù)擬合所有極值點(diǎn)得到上包絡(luò)線{Ut}和下包絡(luò)線{Lt}。

步驟2：計(jì)算上包絡(luò)線和下包絡(luò)線的均值，記為

(1)

步驟3：令

(2)

步驟5：令殘差

(3)

(4)

式中：Rt為EMD分解的殘差。

1.2 LSTM

LSTM是一種引入了遺忘門、輸入門和輸出門3種門控單元和記憶細(xì)胞(可看作一種特殊的隱藏狀態(tài))的特殊形式的循環(huán)神經(jīng)網(wǎng)絡(luò)，這種精心設(shè)計(jì)的結(jié)構(gòu)緩解了RNN的梯度消失問(wèn)題, LSTM神經(jīng)元結(jié)構(gòu)如圖1所示，其中，Xt為當(dāng)前時(shí)刻輸入;Ht-1,Ht分別為前一和當(dāng)前時(shí)刻的隱藏狀態(tài);It、Ft、Ot分別為當(dāng)前時(shí)刻的輸入門、遺忘門和輸出門；Gt為當(dāng)前時(shí)刻的候選記憶細(xì)胞；ct-1、ct分別為前一和當(dāng)前時(shí)刻的記憶細(xì)胞；Yt為輸出。關(guān)于LSTM的更多細(xì)節(jié)可參考文獻(xiàn)[15]。

圖1 LSTM細(xì)胞單元結(jié)構(gòu)

設(shè)隱藏單元數(shù)為h, 在t時(shí)刻輸入Xt∈Rn×d(其中n為樣本數(shù)，d為輸入長(zhǎng)度)，前一時(shí)刻的隱藏狀態(tài)Ht-1∈Rn×h(可視為短期狀態(tài)),則在t時(shí)刻輸入門It∈Rn×h，遺忘門Ft∈Rn×h和輸出門Ot∈Rn×h表達(dá)式為

It=σ(XtWXI+Ht-1WHI+bI)

(5)

Ft=σ(XtWXF+Ht-1WHF+bF)

(6)

Ot=σ(XtWXO+Ht-1WHO+bO)

(7)

式中：It、Ft、Ot都由sigmoid激活函數(shù)控制, 它們的輸出值都在0～1之間。記憶細(xì)胞ct(可視為長(zhǎng)期狀態(tài)) 可由式(9)計(jì)算：

Gt=tanh(XtWXG+Ht-1WHG+bG)

(8)

ct=Ft?ct-1+It?Gt

(9)

式中：WXI、WXF、WXO、WXG∈Rd×h為4個(gè)全連接層關(guān)于它們的輸入向量的權(quán)重矩陣；WHI、WHF、WHO、WHG∈Rh×h為4個(gè)全連接層關(guān)于它們的短期狀態(tài)Ht-1的權(quán)重矩陣；bI、bF、bO、bG∈R1×h為4個(gè)全連接層的線性偏倚；σ(·)和tanh(·)為sigmoid函數(shù)和雙曲正切函數(shù)；?表示向量間的逐元素相乘。隱藏狀態(tài)和輸出的表達(dá)式為

Yt=Ht=Ot?tanh(ct)

(10)

1.3 經(jīng)驗(yàn)?zāi)Ｊ椒纸夥椒ê烷L(zhǎng)短期記憶網(wǎng)絡(luò)

2 試驗(yàn)步驟

2.1 數(shù)據(jù)

以南水北調(diào)工程某河流2017年11月3日到2018年6月30日每隔1 h的中水文數(shù)據(jù)為試驗(yàn)對(duì)象，所記錄的水文數(shù)據(jù)包括水深、流速和瞬時(shí)流量，共有5 727條記錄，水深、流速和瞬時(shí)流量的時(shí)間序列曲線見圖2。

圖2 瞬時(shí)流量、流速和水深數(shù)據(jù)序列

2.2 數(shù)據(jù)預(yù)處理和EMD分解

水流變化會(huì)受到天氣、測(cè)量等隨機(jī)因素的影響，原始數(shù)據(jù)包含一些噪聲數(shù)據(jù)，噪聲數(shù)據(jù)達(dá)到一定程度會(huì)影響預(yù)報(bào)的準(zhǔn)確性，為此采用中值濾波來(lái)消除噪聲數(shù)據(jù)。記{x1,x2,…,xT}為觀測(cè)時(shí)間序列, 其中xt為t時(shí)刻觀測(cè)值,令

(11)

2.3 LSTM模型構(gòu)建

為記號(hào)方便起見，記{y1,y2,…,yT}為L(zhǎng)STM模型需要處理的時(shí)間序列,輸入變量的時(shí)間窗口為K,即選擇當(dāng)前時(shí)刻的前K個(gè)時(shí)間觀測(cè)值作為輸入變量來(lái)預(yù)測(cè)，LSTM模型的輸入為

(12)

圖3 預(yù)測(cè)窗口長(zhǎng)度為12 h的LSTM模型對(duì)瞬時(shí)流量的EMD分量預(yù)測(cè)結(jié)果

輸出為

(13)

式中：L為預(yù)測(cè)窗口長(zhǎng)度, 也就是說(shuō)用前面的K個(gè)時(shí)刻的觀測(cè)值來(lái)預(yù)測(cè)接下來(lái)L個(gè)時(shí)刻的數(shù)據(jù)。使用均方根誤差ER:

(14)

作為損失函數(shù)來(lái)確定隱藏狀態(tài)的維數(shù)和式(5)～(10)的權(quán)重矩陣，采用均方根誤差ER和平均相對(duì)誤差EM來(lái)評(píng)價(jià)模型的預(yù)測(cè)精度，EM定義為

(15)

2.4 預(yù)測(cè)結(jié)果

把以南水北調(diào)工程某河流2017年11月3日到2018年6月30日每隔1 h瞬時(shí)流量、流速和水深3個(gè)時(shí)間序列分別按照3∶1作為訓(xùn)練集數(shù)據(jù)和測(cè)試集數(shù)據(jù)，即每個(gè)序列的前4 259個(gè)數(shù)據(jù)作為訓(xùn)練集，后1 432個(gè)數(shù)據(jù)作為測(cè)試集。

首先對(duì)中值濾波預(yù)處理后的瞬時(shí)流量數(shù)據(jù)進(jìn)行EMD分解為本征模函數(shù)部分和殘差部分，對(duì)EMD分解的各個(gè)本征模函數(shù)和殘差部分分別訓(xùn)練LSTM方法進(jìn)行預(yù)報(bào)，并將預(yù)報(bào)結(jié)果進(jìn)行疊加構(gòu)建瞬時(shí)流量的預(yù)測(cè)結(jié)果。LSTM實(shí)驗(yàn)在PyTorch環(huán)境下進(jìn)行，在數(shù)據(jù)預(yù)處理階段中值濾波的窗口長(zhǎng)度設(shè)置為2，LSTM模型訓(xùn)練中以均方根誤差作為損失函數(shù)，使用的梯度下降優(yōu)化器為Adam優(yōu)化器，輸入變量的時(shí)間窗口K=24 h，即每次預(yù)測(cè)需要前24 h的觀測(cè)數(shù)據(jù)，epoch(指所有的訓(xùn)練樣本輸入網(wǎng)絡(luò)中完成一次前向和反向傳播的過(guò)程)設(shè)為200，學(xué)習(xí)率為0.000 1，隱藏單元數(shù)為64。圖3為預(yù)測(cè)窗口長(zhǎng)度L=12 h的LSTM模型對(duì)瞬時(shí)流量的EMD分量(I(1),I(2),…,I(7)和殘差Rt)預(yù)測(cè)結(jié)果。

圖4 不同預(yù)測(cè)窗口長(zhǎng)度兩種方法對(duì)瞬時(shí)流量的預(yù)測(cè)結(jié)果

表1 兩種模型下瞬時(shí)流量預(yù)測(cè)誤差比較

窗口長(zhǎng)度為12 h關(guān)于瞬時(shí)流量LSTM預(yù)測(cè)結(jié)果和EMD-LSTM預(yù)測(cè)的結(jié)果比較。改變預(yù)測(cè)時(shí)間窗口長(zhǎng)度為6 h，對(duì)瞬時(shí)流量進(jìn)行預(yù)測(cè)，兩種模型的預(yù)測(cè)效果見圖4(c)(d)。由圖4可知，LSTM方法和EMD-LSTM方法都能取得較好的預(yù)測(cè)效果，特別當(dāng)預(yù)測(cè)時(shí)間較短時(shí)，預(yù)測(cè)效果更優(yōu)。為了更為顯著比較兩種方法優(yōu)劣，截取2018年6月11日0點(diǎn)到6月25日24點(diǎn)共360 h數(shù)據(jù)(該時(shí)段內(nèi)的數(shù)據(jù)變換形式較多，因此更能衡量預(yù)測(cè)算法優(yōu)劣。)在預(yù)測(cè)窗口長(zhǎng)度為12 h和 6 h兩種情況下的LSTM預(yù)測(cè)方法和EMD-LSTM預(yù)測(cè)方法的結(jié)果進(jìn)行比較，如圖5所示。由圖5可知，對(duì)瞬時(shí)流量而言，EMD-LSTM預(yù)測(cè)方法優(yōu)于LSTM方法，特別對(duì)幅度變化較大時(shí)段，EMD-LSTM方法優(yōu)勢(shì)更為顯著。同時(shí)，短時(shí)段的預(yù)測(cè)效果優(yōu)于長(zhǎng)時(shí)段，原因在于預(yù)測(cè)窗口長(zhǎng)度太長(zhǎng)，數(shù)據(jù)信號(hào)會(huì)產(chǎn)生新的變化，而這些變化很難從歷史數(shù)據(jù)中學(xué)習(xí)。

圖5 不同預(yù)測(cè)窗口長(zhǎng)度瞬時(shí)流量預(yù)測(cè)結(jié)果比較

進(jìn)一步對(duì)流速和水深分別使用LSTM模型和EMD-LSTM模型進(jìn)行預(yù)測(cè)，預(yù)測(cè)效果見圖6～9，預(yù)測(cè)誤差見表2和表3。從對(duì)瞬時(shí)流量、流速和水深這3種指標(biāo)的預(yù)測(cè)結(jié)果可以看出，兩種模型都可以很好地預(yù)測(cè)水文數(shù)據(jù)的變化，但EMD-LSTM模型的預(yù)測(cè)精度更高，特別對(duì)波動(dòng)較大時(shí)段，EMD-LSTM方法優(yōu)勢(shì)更為顯著，這是由于LSTM方法適合于預(yù)測(cè)相對(duì)比較規(guī)律的信號(hào)，高頻噪聲的干擾會(huì)極大地降低LSTM方法的預(yù)測(cè)精度，使用EMD-LSTM可以把信號(hào)中的高頻部分分離出來(lái)預(yù)測(cè)，極大地減少了高頻噪聲的干擾。另外，兩種模型下，預(yù)測(cè)窗口長(zhǎng)度越小，均方根誤差和平均相對(duì)誤差越小，預(yù)測(cè)越準(zhǔn)確，曲線擬合越好，即短時(shí)段的預(yù)測(cè)效果優(yōu)于長(zhǎng)時(shí)段，原因在于預(yù)測(cè)窗口長(zhǎng)度太長(zhǎng)，受氣候、環(huán)境和人類活動(dòng)影響，水文數(shù)據(jù)信號(hào)會(huì)產(chǎn)生新的變化，而這些變化很難從歷史數(shù)據(jù)中學(xué)習(xí)。

圖6 不同預(yù)測(cè)窗口長(zhǎng)度兩種方法對(duì)流速的預(yù)測(cè)結(jié)果

圖7 不同預(yù)測(cè)窗口長(zhǎng)度流速預(yù)測(cè)結(jié)果比較

圖8 不同預(yù)測(cè)窗口長(zhǎng)度兩種方法對(duì)水深的預(yù)測(cè)結(jié)果

圖9 不同預(yù)測(cè)窗口長(zhǎng)度水深預(yù)測(cè)結(jié)果比較

表2 兩種模型流速預(yù)測(cè)誤差比較

表3 兩種模型水深預(yù)測(cè)誤差比較

3 結(jié) 論

本文建立了基于EMD-LSTM神經(jīng)網(wǎng)絡(luò)的預(yù)測(cè)模型，并分別對(duì)水深、流速和瞬時(shí)流量序列預(yù)報(bào)，對(duì)每個(gè)序列預(yù)報(bào)的主要步驟是：①使用中值濾波消除噪聲數(shù)據(jù)；②對(duì)預(yù)處理后的數(shù)據(jù)進(jìn)行EMD分解為IMF部分和殘差部分；③對(duì)步驟②中的每一部分使用LSTM方法進(jìn)行預(yù)報(bào)，并將預(yù)報(bào)結(jié)果進(jìn)行疊加。利用LSTM和EMD-LSTM兩種模型分別對(duì)3個(gè)水文指標(biāo)數(shù)據(jù)序列在不同預(yù)測(cè)窗口下進(jìn)行預(yù)測(cè)，結(jié)果表明：①兩種模型下，預(yù)測(cè)窗口長(zhǎng)度越小，均方根誤差和平均相對(duì)誤差越小，預(yù)測(cè)越準(zhǔn)確，曲線擬合越好，即短時(shí)段的預(yù)測(cè)效果優(yōu)于長(zhǎng)時(shí)段；②兩種模型都可以很好地預(yù)測(cè)水文數(shù)據(jù)的變化,但EMD-LSTM模型的預(yù)測(cè)精度更高，特別對(duì)波動(dòng)較大時(shí)段，EMD-LSTM方法優(yōu)勢(shì)更為顯著。EMD-LSTM方法更適合于含有高頻噪聲干擾的水文時(shí)間序列。

本文提供的方法可以很好地對(duì)水文數(shù)據(jù)進(jìn)行短期預(yù)報(bào)，為水位預(yù)判和水資源的實(shí)時(shí)調(diào)度提供決策依據(jù)和技術(shù)支撐。本研究數(shù)據(jù)不足1 a，模型未考慮季節(jié)因素，后面隨著數(shù)據(jù)的積累，可以進(jìn)一步選取季節(jié)特征進(jìn)行研究。