楊瓊波,崔東文
(1.云南省水文水資源局紅河分局,云南 紅河 661199;2.云南省文山州水務(wù)局,云南 文山 663000)
主要從事水文水資源分析研究和水資源管理保護等工作;崔東文(通訊作者).
徑流時間序列預(yù)測是水文預(yù)報研究領(lǐng)域的熱點和難點。提高徑流時間序列預(yù)測精度對區(qū)域水資源開發(fā)利用、防洪抗旱規(guī)劃、水資源管理保護、水庫優(yōu)化調(diào)度等具有重要意義。由于徑流時間序列人類活動、受氣候變化、土地利用及植被覆蓋等多重因素的影響,常表現(xiàn)出多尺度、非線性、非平穩(wěn)性等特征,基于“分解—預(yù)測—重構(gòu)”模式的預(yù)測方法被廣泛應(yīng)用于徑流時間序列預(yù)測[1- 4],并取得了較好的預(yù)測效果。目前,常見的時間序列分解方法有經(jīng)驗?zāi)B(tài)分解(Empirical mode decomposition,EMD)、完全集合經(jīng)驗?zāi)B(tài)分解(Complete ensemble empirical mode decomposition,CEEMD)、小波分解(Wavelet decomposition,WD)、奇異譜分解(Singular spectrum analysis,SSA)、變分模態(tài)分解(Variational mode decomposition,VMD)等;預(yù)測模型有BP神經(jīng)網(wǎng)絡(luò)、廣義回歸神經(jīng)網(wǎng)絡(luò)(GRNN)、支持向量機(SVM)、相關(guān)向量機(RVM)、極限學(xué)習(xí)機(ELM)、LSTM神經(jīng)網(wǎng)絡(luò)、秩次集對分析法等。
WD是一種被廣泛應(yīng)用的時頻分析工具,主要通過伸縮平移運算對信號逐步進行多尺度細化,最終達到高頻處時間細分,低頻處頻率細分的目的,已在各領(lǐng)域得到廣泛應(yīng)用。然而,WD不足之處在于其僅對信號低頻部分詳細分解而忽略了信號的高頻部分。小波包分解(Wavelet packet decomposition,WPD)源于WD,它在分解信號低頻子集同時,對高頻子集繼續(xù)分解,且能夠根據(jù)分解需要自行設(shè)定分解層數(shù)和分解使用的小波函數(shù)[5-7],已在各種時序數(shù)據(jù)分解中得到廣泛應(yīng)用。人工水母搜索(Artificial Jellyfish Search,AJS)算法是Chou等人于2020年基于水母覓食行為而提出一種新型群體智能算法,已在工程結(jié)構(gòu)優(yōu)化中得到應(yīng)用[8]。數(shù)據(jù)分組處理方法(Group method of data handling,GMDH)是針對復(fù)雜系統(tǒng)自組織建模的一種前饋型多層神經(jīng)網(wǎng)絡(luò),網(wǎng)絡(luò)僅通過遞階分層取舍神經(jīng)元和調(diào)整閾值,實現(xiàn)對復(fù)雜非線性規(guī)律進行擬合,目前已在各種預(yù)測研究中得到應(yīng)用[9-11]。在實際應(yīng)用中,GMDH網(wǎng)絡(luò)權(quán)值、多項式系數(shù)、網(wǎng)絡(luò)層最大神經(jīng)元數(shù)等參數(shù)對GMDH網(wǎng)絡(luò)預(yù)測性能有著重要影響,目前主要有GMDH網(wǎng)絡(luò)權(quán)值優(yōu)化[12]、多項式系數(shù)優(yōu)化[13]兩種途徑,而鮮見于針對GMDH網(wǎng)絡(luò)層最大神經(jīng)元數(shù)、最大網(wǎng)絡(luò)層數(shù)和學(xué)習(xí)率等關(guān)鍵參數(shù)的優(yōu)化。
本文以云南省龍?zhí)墩?952年1月~2016年12月月徑流時間序列預(yù)測為例,基于WPD、AJS與GMDH網(wǎng)絡(luò)構(gòu)建WPD-AJS-GMDH徑流時間序列預(yù)測模型。內(nèi)容安排如下:①利用WPD及對比分解方法CEEMD、WD將月徑流時間序列數(shù)據(jù)分解為若干子序列分量;②介紹AJS算法,在不同維度條件下通過6個典型函數(shù)對AJS算法進行仿真驗證,采用AJS優(yōu)化GMDH網(wǎng)絡(luò)層最大神經(jīng)元個數(shù)、最大網(wǎng)絡(luò)層數(shù)、學(xué)習(xí)率、選擇壓力系數(shù),構(gòu)建WPD-AJS-GMDH模型,并構(gòu)建WPD-AJS-SVM、WPD-AJS-BP、WPD-GMDH、WPD-SVM、WPD-BP、CEEMD-AJS-GMDH、CEEMD-AJS-SVM、CEEMD-AJS-BP、CEEMD-GMDH、CEEMD-SVM、CEEMD-BP和WD-AJS-GMDH、WD-AJS-SVM、WD-AJS-BP、WD-GMDH、WD-SVM、WD-BP作對比分析模型;③利用所構(gòu)建的18種模型對實例月徑流時間序列進行檢驗。
WPD是在小波分解(WD)的基礎(chǔ)上改進和發(fā)展而來的一種有效分解方法,能同時對高頻、低頻部分實施信號分解。理論上,3層小波包分解能夠提取絕大部分信號中的有效信息,逼近任意非線性函數(shù),從而達到解決實際問題的目的。小波包分解算法公式[14-16]為
(1)
重構(gòu)算法為
(2)
1.2.1 AJS算法數(shù)學(xué)描述
AJS算法是基于水母覓食行為而提出一種新型元啟發(fā)式算法。該算法通過模擬水母洋流移動、群內(nèi)移動(主動移動和被動移動)和移動切換策略來達到求解優(yōu)化問題的目的。AJS算法基于3個理想化規(guī)則:①水母要么跟隨洋流移動,要么在群內(nèi)移動,“時間控制策略”控制著兩種移動之間的切換;②水母在海洋中尋找食物,他們更容易吸引到食物豐富的區(qū)域;③找到的食物數(shù)量取決于水母位置及其相應(yīng)的目標函數(shù)值[7]。
AJS算法數(shù)學(xué)描述簡述如下[7]:
(1)洋流移動。AJS通過計算海洋中每只水母到當(dāng)前最佳水母位置的平均矢量值確定洋流方向(趨勢)。 洋流移動數(shù)學(xué)描述
(3)
(4)
表1 AJS算法標準測試函數(shù)尋優(yōu)結(jié)果
式中,Xi(t+1)為第i只水母第t+1次迭代位置;Xi(t)為第i只水母第t次迭代位置;其他參數(shù)意義同上。
(2)水母群移動。在水母群中,水母分為被動移動(A型)和主動移動(B型)兩種類型。當(dāng)水母群剛形成時,大多數(shù)水母呈A型移動;隨著時間推移,大多數(shù)水母呈B型移動。A型移動表示水母在其自身位置周圍移動,其位置更新公式為
Xi(t+1)=Xi(t)+γ×rand(0,1)×(Ub-Lb)
(5)
式中,γ為運動系數(shù),本文取0.1;Ub、Lb分別為搜索空間上、下限值;其他參數(shù)意義同上。B型運動描述為:當(dāng)水母Xj所處區(qū)域食物量超過水母Xi所處區(qū)域食物量時,水母Xi向水母Xj移動;反知,水母Xi遠離水母Xj。其位置更新公式描述如下
(6)
式中,Xj(t)為第j只水母的第t次迭代位置;f(Xi)、f(Xj)分別為水母Xi、水母Xj的適應(yīng)度值;其他參數(shù)意義同上。
(3)時間控制切換策略。水母在水母群中的運動分為A型和B型,為更好地模擬水母在這兩種移動之間的切換,引入時間控制函數(shù)c(t)來描述這種情況,當(dāng)c(t)≥c0時,水母跟隨洋流移動;當(dāng)c(t) c(t)=|(1-t/Tmax)×[2×rand(0,1)-1]| (7) 式中,t為當(dāng)前迭代次數(shù);Tmax為算法最大迭代次數(shù);c0為時間控制常數(shù),本文取0.5。 1.2.2 AJS算法仿真驗證 在不同維度條件下,選取Sphere等6個測試函數(shù)對AJS算法進行仿真驗證(見表1)。AJS算法通過Matlab 2018a M語言實現(xiàn),利用函數(shù)20次尋優(yōu)平均值評估AJS算法的性能。設(shè)置AJS算法最大迭代次數(shù)Tmax=3 000,種群規(guī)模npop=100。 對于單峰函數(shù),AJS算法在不同維度條件下20次尋優(yōu)精度均在7.20×10-83以上,表現(xiàn)出較好的尋優(yōu)精度;對于多峰函數(shù),AJS算法在不同維度條件下20次尋優(yōu)精度均在1.22×10-12以上,具有較好全局尋優(yōu)能力??梢?,對于表1中的函數(shù),在不同維度條件下AJS算法均表現(xiàn)出較好尋優(yōu)效果。 GMDH網(wǎng)絡(luò)是針對復(fù)雜系統(tǒng)自組織建模的一種前饋型多層神經(jīng)網(wǎng)絡(luò),也稱多項式網(wǎng)絡(luò),目前已在電主軸熱位移[17]、經(jīng)濟短期預(yù)測[18]、碳價格預(yù)測[19]、潮汐預(yù)報[20]、用戶流失預(yù)測[21]等方面得到應(yīng)用。GMDH網(wǎng)絡(luò)由系統(tǒng)各輸入單元交叉組合產(chǎn)生活動神經(jīng)元,從已產(chǎn)生的一代神經(jīng)元中選擇若干與目標變量最接近的神經(jīng)元,通過強強結(jié)合再次產(chǎn)生新的神經(jīng)元,……。通過重復(fù)這一優(yōu)勝劣汰和逐漸進化的過程,直到獲得外準則值最低的最優(yōu)神經(jīng)元,即完成GMDH網(wǎng)絡(luò)模型訓(xùn)練過程[9-11]。 yi=f(xi1,xi2,…,xin),i=1,2,…,N (8) 式中,yi為實際輸出;xi1,xi2,…,xin為輸入;N為輸入變量的數(shù)量。 (9) (10) 在GMDH預(yù)測過程中,除網(wǎng)絡(luò)權(quán)值、多項式系數(shù)對GMDH預(yù)測效果產(chǎn)生直接影響外,GMDH網(wǎng)絡(luò)層最大神經(jīng)元數(shù)、最大網(wǎng)絡(luò)層數(shù)、學(xué)習(xí)率等關(guān)鍵參數(shù)對GMDH網(wǎng)絡(luò)預(yù)測性能同樣具有重要影響。為進一步提升GMDH網(wǎng)絡(luò)預(yù)測性能,拓展GMDH網(wǎng)絡(luò)優(yōu)化方向,本文利用AJS算法對GMDH網(wǎng)絡(luò)層最大神經(jīng)元數(shù)、最大網(wǎng)絡(luò)層數(shù)、學(xué)習(xí)率、選擇壓力系數(shù)進行優(yōu)化,以期進一步提高GMDH網(wǎng)絡(luò)預(yù)測性能。 WPD-AJS-GMDH模型預(yù)測實現(xiàn)如下: 步驟一,利用WPD將原月徑流時間序列分解為8子序列分量[3,0]~[3,7]。利用自相關(guān)函數(shù)法(AFM)確定各子序列分量的輸入、輸出向量,合理劃分訓(xùn)練、預(yù)測樣本[22],見圖1。 步驟二,構(gòu)建訓(xùn)練樣本均方誤差作為AJS-GMDH模型適應(yīng)度函數(shù) (11) 步驟三,設(shè)置AJS算法水母種群規(guī)模npop、最大迭代次數(shù)Tmax和算法終止條件。利用式(12)初始化水母空間位置Xi(i=1,2,…,npop),則 Xi+1=ηXi(1-Xi),0≤X0≤1 (12) 式中,Xi為第i只水母位置的logistic混沌值;X0用于生成水母的初始種群,X0∈(0,1),且X0?{0.00,0.25,0.50,0.75,1.00 };參數(shù)η設(shè)置為4.0。 步驟四,計算每只水母Xi適應(yīng)度值,找到并保存最佳水母適應(yīng)度fbest和最佳水母個體位置X*。 步驟五,利用式(7)計算時間控制函數(shù)c(t)。若c(t)≥0.5,利用式(4)更新水母位置;否則,利用式(5)、(6)更新水母位置。 步驟六,計算當(dāng)前每只水母Xi適應(yīng)度值;比較并保存當(dāng)前最佳水母適應(yīng)度值fbest和最佳水母個體位置X*。 步驟七,判斷算法是否達到終止條件,若是,則輸出最佳水母個體位置X*,算法結(jié)束。 步驟八,輸出最佳水母個體位置X*。即,GMDH網(wǎng)絡(luò)層最大神經(jīng)元數(shù)、最大網(wǎng)絡(luò)層數(shù)、學(xué)習(xí)率、選擇壓力系數(shù);利用AJS-GMDH模型對分量[3,0]~[3,7]進行預(yù)測,疊加重構(gòu)后即得到實例月徑流最終預(yù)測結(jié)果。 步驟九,利用平均絕對百分比誤差(MAPE)、平均絕對誤差(MAE)、納什系數(shù)(NSE)和均方根誤差(RMSE)對各模型有效性進行評估。即 (13) 圖1 實例月徑流預(yù)測流程 龍?zhí)墩疚挥谠颇鲜∥纳绞信手ㄦ?zhèn)龍?zhí)墩?,系紅河流域瀘江水系盤龍河干流控制站,控制徑流面積3 128 km2,為我國重要的水文站和中央報汛站。 2.1.1 WPD分解 基于demy小波包基,利用WPD將龍?zhí)墩?952年1月~2016年12月逐月徑流數(shù)據(jù)進行3層小波包分解,得到8個子序列分量數(shù)據(jù),見圖2。 從圖2來看,[3,0]~[3,3]為低頻空間數(shù)據(jù)波形,振幅較大、波長較短,大致反映月徑流時間序列的變化趨勢;[3,4]、[3,5]、[3,6]、[3,7]為高頻空間數(shù)據(jù)波形,從[3,4]分量到[3,7]分量,其振幅逐漸減小、波長逐漸變長,大致反映月徑流時間序列的波動情況[4]。 2.1.2 CEEMD分解 利用CEEMD方法對龍?zhí)墩緦嵗聫搅鲾?shù)據(jù)進行分解,分解結(jié)果見圖3。從圖3可以看出,高頻子序列RSE的頻率最高、幅度最大、波長最短;中頻子序列IMF4~IMF8和低頻子序列IMF1~IMF3的振幅逐漸減小,頻率逐漸降低,波長逐漸增大,周期性逐漸增強,符合經(jīng)驗?zāi)B(tài)分解規(guī)律[18]。 2.1.3 WD分解 利用db5小波將龍?zhí)墩緦嵗聫搅鲾?shù)據(jù)進行分解,得到1個低頻率分量rA5和5個高頻率分量rD1、rD2、rD3、rD4、rD5,見圖4。從圖4可以看出,高頻分量反映原始數(shù)據(jù)的變化特征,低頻分量反映原始數(shù)據(jù)的變化趨勢。 圖2 實例月徑流時間序列WPD分解3D效果 圖3 實例月徑流時間序列CEEMD分解3D效果 圖4 實例月徑流時間序列WD分解3D效果 本文采用自相關(guān)函數(shù)法(AFM)確定各分量的輸入、輸出向量,即通過計算各子序列分量的自相關(guān)系數(shù),在滯后數(shù)≥3情況下,選取自相關(guān)系數(shù)最大時所對應(yīng)的滯后數(shù)作為各子序列分量最優(yōu)輸入維數(shù),通過最優(yōu)輸入維數(shù)構(gòu)建預(yù)測模型輸入、輸出向量,見表2。本文利用實例前649~657個月徑流數(shù)據(jù)作為訓(xùn)練樣本,后120個月作為預(yù)測樣本[4]。 表2 各分量自相關(guān)系數(shù)、輸入維數(shù)及序列長度 從表2可以看出,經(jīng)WPD分解的各子序列分量自相關(guān)系數(shù)絕對值均在0.591以上;經(jīng)CEEMD分解的IMF1分量、經(jīng)WD分解的rD1分量自相關(guān)系數(shù)絕對值小于0.25,這在很大程度上制約分解效果和預(yù)測精度。 2.3.1 參數(shù)設(shè)置 (1)WPD-AJS-GMDH、WPD-GMDH、CEEMD-AJS-GMDH、CEEMD-GMDH、WD-AJS-GMDH、WD-GMDH模型:設(shè)置WPD-GMDH、CEEMD-GMDH、WD-GMDH模型網(wǎng)絡(luò)層最大神經(jīng)元數(shù)m=20;最大網(wǎng)絡(luò)層數(shù)r=5;選擇壓力系數(shù)α=0.6;學(xué)習(xí)率p=0.85;WPD-AJS-GMDH、CEEMD-AJS-GMDH、WD-AJS-GMDH模型AJS算法最大迭代次數(shù)Tmax=50,種群規(guī)模npop=30;網(wǎng)絡(luò)層最大神經(jīng)元數(shù)m搜索空間[5,50],最大網(wǎng)絡(luò)層數(shù)r搜索空間[2,10],選擇壓力系數(shù)α搜索空間[0.01,0.9],學(xué)習(xí)率p搜索空間[0.01,0.9],其他采用默認值。 (2)WPD-AJS-SVM、WPD-SVM、CEEMD-AJS-SVM、CEEMD-SVM、WD-AJS-SVM、WD-SVM模型:WPD-SVM、CEEMD-SVM、WD-SVM模型懲罰因子、核函數(shù)參數(shù)和不敏感損失系數(shù)通過人工試湊法確定;WPD-AJS-SVM、CEEMD-AJS-SVM、WD-AJS-SVM模型AJS算法最大迭代次數(shù)Tmax=50,種群規(guī)模npop=30;不敏感損失系數(shù)、懲罰因子、核函數(shù)參數(shù)搜索空間分別設(shè)置為[0.001,1]、[0.01,100]、[0.01,100],輸入數(shù)據(jù)采用[0,1]進行歸一化處理。 (3)WPD-AJS-BP、WPD-BP、CEEMD-AJS-BP、CEEMD-BP、WD-AJS-BP、WD-BP模型:設(shè)置WPD-BP、CEEMD-BP、WD-BP模型隱層傳遞函數(shù)、輸出層傳遞函數(shù)、訓(xùn)練函數(shù)分別為tansig、purelin、trainlm,設(shè)定期望誤差0.000 001,最大訓(xùn)練輪回1 000次,隱層節(jié)點數(shù)為輸入維數(shù)的2倍減1;WPD-AJS-BP、CEEMD-AJS-BP、WD-AJS-BP模型AJS算法最大迭代次數(shù)Tmax=50,種群規(guī)模npop=30,權(quán)值閾值搜索空間∈[-1,1],BP部分設(shè)置同WPD-BP、CEEMD-BP、WD-BP模型。輸入數(shù)據(jù)采用[0,1]進行歸一化處理。 2.3.2 實例預(yù)測結(jié)果分析 基于表2構(gòu)建各分量的輸入、輸出向量,并利用所建立的WPD-AJS-GMDH等18種模型對各分量進行預(yù)測,將預(yù)測結(jié)果疊加重構(gòu)后即得到實例月徑流最終預(yù)測結(jié)果。各模型采用上述MAPE、MAE、RMSE、NSE進行評估,結(jié)果見表3,預(yù)測效果見圖5、6。 由表3及圖5、6可以得出以下結(jié)論: (1)WPD-AJS-GMDH、WPD-AJS-SVM、WPD-AJS-BP模型的預(yù)測效果優(yōu)于WPD-GMDH、WPD-SVM、WPD-BP模型,遠優(yōu)于其他模型。其中,尤以WPD-AJS-GMDH模型預(yù)測的精度最高、效果最好,其平均絕對百分比誤差2.14%,平均絕對誤差0.25 m3/s,納什系數(shù)0.999 4,均方根誤差0.331 m3/s。將WPD-AJS-GMDH模型用于月徑流時間序列預(yù)測是可行和可靠的,模型具有更高的預(yù)測和更好的泛化能力。 (2)從不同分解方法預(yù)測結(jié)果對比來看,基于WPD分解的同種模型的預(yù)測效果要遠優(yōu)于基于CEEMD、WD分解的同種模型,表明WPD對月徑流時間序列的分解效果要優(yōu)于CEEMD、WD方法。 表3 實例月徑流時間序列預(yù)測結(jié)果對比 圖5 實例月徑流時間序列預(yù)測相對誤差3D 圖6 實例月徑流時間序列預(yù)測絕對誤差3D (3)基于AJS算法優(yōu)化的同種模型的預(yù)測效果要優(yōu)于未經(jīng)AJS算法優(yōu)化的同種模型,表明AJS算法能有效優(yōu)化GMDH網(wǎng)絡(luò)層最大神經(jīng)元數(shù)、最大網(wǎng)絡(luò)層數(shù)、學(xué)習(xí)率和選擇壓力系數(shù),以及SVM懲罰因子、核函數(shù)參數(shù)、不敏感損失系數(shù)和BP網(wǎng)絡(luò)權(quán)閾值,有效提高GMDH網(wǎng)絡(luò)、SVM、BP網(wǎng)絡(luò)的預(yù)測性能。 (4)從圖5、6可看出,WPD-AJS-GMDH模型預(yù)測的絕大多數(shù)樣本絕對百分比誤差、平均絕對誤差均在0值附近波動,具有較優(yōu)的預(yù)測效果;次之為WPD-AJS-SVM、WPD-AJS-BP模型;再次為WPD-GMDH、WPD-SVM、WPD-BP模型,其他12種模型則均較差。 本研究將WPD、AJS與GMDH網(wǎng)絡(luò)進行結(jié)合,構(gòu)建WPD-AJS-GMDH月徑流時間序列預(yù)測模型,利用云南省龍?zhí)墩?952年1月~2016年12月共780組月徑流時序數(shù)據(jù)進行實證分析,并建立基于CEEMD、WD分解的GMDH、SVM、BP模型做對比分析。結(jié)果如下: (1)在不同維度條件下,AJS算法均具有較好的尋優(yōu)效果,將AJS算法用于GMDH網(wǎng)絡(luò)層最大神經(jīng)元數(shù)、最大網(wǎng)絡(luò)層數(shù)、學(xué)習(xí)率和選擇壓力系數(shù)尋優(yōu)是可靠的。 (2)WPD-AJS-GMDH模型對實例后120個月月徑流時間序列預(yù)測的平均絕對百分比誤差為2.14%,平均絕對誤差為0.25 m3/s,納什系數(shù)0.999 4,均方根誤差0.331 m3/s,預(yù)測效果略優(yōu)于WPD-AJS-SVM、WPD-AJS-BP模型,優(yōu)于WPD-GMDH、WPD-SVM、WPD-BP模型,遠優(yōu)于其他模型。將WPD-AJS-GMDH模型用于月徑流時間序列預(yù)測是可行和可靠的,模型具有更高的預(yù)測和更好的泛化能力。 (3)實例驗證表明,WPD能將徑流時序數(shù)據(jù)分解為更具規(guī)律的子序列分量,有效弱化復(fù)雜環(huán)境對徑流時間序列的影響,降低預(yù)測復(fù)雜度,其分解效果明顯優(yōu)于CEEMD、WD方法。 (4)實例驗證表明,AJS算法能有效優(yōu)化GMDH網(wǎng)絡(luò)層最大神經(jīng)元數(shù)、最大網(wǎng)絡(luò)層數(shù)、學(xué)習(xí)率和選擇壓力系數(shù),不但克服了GMDH網(wǎng)絡(luò)人工選擇關(guān)鍵參數(shù)的不足,而且大大提高了GMDH網(wǎng)絡(luò)的預(yù)測精度和智能化水平,具有較好的實用價值。1.3 GMDH網(wǎng)絡(luò)
1.4 建模流程
2 實例應(yīng)用
2.1 時序數(shù)據(jù)分解
2.2 時間序列建模
2.3 參數(shù)設(shè)置及預(yù)測分析
3 結(jié) 論