梁 浩,黃生志,孟二浩,黃 強
(西安理工大學 西北旱區(qū)生態(tài)水利國家重點實驗室,陜西 西安 710048)
河川徑流作為水循環(huán)的關鍵環(huán)節(jié),對其進行準確預報十分重要[1]。當前水資源系統(tǒng)受到諸多因素的影響,徑流序列表現(xiàn)出多時間尺度變化特性和非平穩(wěn)特性,這給徑流序列的預測帶來了難度。尤其是在極端天氣頻發(fā)的背景下,徑流序列呈現(xiàn)更大的波動性[2]。當前,徑流預測模型已有不少,從不同的模型中優(yōu)選可靠的徑流預報模型,獲得更高的徑流預測精度,是水文預報工作者迫切開展的工作。
為了獲得未來幾天、幾個月甚至幾年的徑流預測值,水文預報人員已經提出了多種方法。這些水文預報模型方法可以大致分為過程驅動和數(shù)據(jù)驅動兩大塊[3]。過程驅動的模型以基于物理過程的水文模型為主,如王中跟等[4]通過SWAT模型利用地理信息系統(tǒng)和遙感的數(shù)據(jù)模擬流域上復雜的水文過程。數(shù)據(jù)驅動則主要通過針對水文站、氣象站等觀測的數(shù)據(jù),結合統(tǒng)計方法進行徑流模擬預測,如運用廣泛的多元線性回歸(MLR)[5-6]、人工神經網絡模型(ANN)[7-8]和支持向量機模型(SVM)[9-10]等諸多基于人工智能的模型[11]。水文模型要涉及較為復雜的水文過程,需要大量準確的水文過程資料,需要考慮不同時空尺度下水文系統(tǒng)規(guī)律的差異性,存在參數(shù)不確定性以及水文模型普適性較差等問題[12]。數(shù)據(jù)驅動模型不需要考慮水文過程的物理機制,對時間序列進行數(shù)學分析,在解決數(shù)據(jù)預測問題以及時間序列分析等方面有獨特的優(yōu)勢[13-14]。模型對不同數(shù)據(jù)有其適用性,本文首先優(yōu)選多元線性回歸、人工神經網絡和支持向量機3種數(shù)據(jù)驅動模型中對渭河流域徑流序列適用性最強的單一模型。
徑流過程本身具有一定的周期性、趨勢性以及隨機性特征,給徑流預測帶來難度,所以徑流的多時間尺度分析和認識對提高徑流預測精度十分必要。分解技術可以將徑流序列分解為一系列相對穩(wěn)定的子序列。近年來研究者基于分解技術進行的徑流預測研究已經有很多,對徑流進行分解重構,最終實現(xiàn)徑流預測的方法不斷被應用,有經驗模態(tài)分解法(EMD)、集合經驗模態(tài)分解法(EEMD)和小波分解(WD)等方法。張洪波等[15]探究了對EMD 分解預測后的分量采用不同的重構方式及對于徑流預報精度的影響。Tan 等[16]和杜懿等[17]采用EEMD 對徑流進行分解,更好地得到徑流的有用信息,提高預測精度。Mohammad 等[18]將EEMD 應用于韓國和伊朗河流的模擬預測,也驗證了其良好性能。周婷等[19]探討了基于小波分解的SVM 徑流模擬預測,徑流預測精度得到了提高。以上研究表明了基于分解技術的徑流研究對徑流預測精度的提高是有效的。然而,大多數(shù)研究者聚焦于采用單一分解方法的混合模型進行徑流預測研究,缺乏對基于多種分解方法的混合模型進行綜合對比分析及其適應性研究。
以徑流時間序列為研究對象的數(shù)據(jù)驅動型徑流預測模型,一般情況下模型的輸入項為降雨等氣象因子,同時考慮歷史的徑流序列。以輸入因子為依據(jù)的預測模型可以大致分為3類:基于氣象因素的徑流預測模型;基于歷史徑流序列以及氣象因素的徑流預測模型;僅考慮歷史徑流序列的徑流預測模型。白云天等[20]在對尼爾基水庫的入庫徑流預報時選取的預報因子為前期徑流以及降水。變化環(huán)境下水文序列發(fā)生了改變,給徑流預測研究帶來了困難。為得到較高的徑流預測精度,豐富輸入因子也是關鍵途徑之一。孟二浩等[21]在徑流預測時豐富了預報因子,加入北大西洋年代際震蕩、北極震蕩等因子對徑流預測精度提高起到了作用。楊開斌等[22]探究了在汛期時ENSO 事件對于徑流預測的影響。相關的研究表明大氣環(huán)流異常因子對徑流具有一定影響。值得一提的是,當前的研究主要集中于大氣環(huán)流異常因子對于徑流整體序列的影響,而沒有考慮分解技術下大氣環(huán)流異常因子對徑流子序列的影響。
基于此,本文對比分析了基于EMD、EEMD和WD方法構建的多種混合模型進行徑流預報的精度和效果,并融合大氣環(huán)流異常因子的信息,以期揭示不同混合模型對河川徑流的預測能力,為變化環(huán)境下徑流預報精度的提高提供新思路。
2.1 多元線性回歸多元線性回歸(MLR)是用于擬合多個自變量和因變量之間的線性關系,通過訓練確定得出模型需要的參數(shù),從而確定具體的MLR 方程,之后便可通過此方程來預測因變量[23]。其一般形式為:
式中:Y為因變量;βi為回歸系數(shù),i=0,1,…,k;Xi為自變量,i=1,2,…,k;μ為隨機誤差項,?;貧w系數(shù)β的求解采用最小二乘法,使誤差項的平方和最小。
2.2 人工神經網絡人工神經網絡(ANN)是一種通過模仿動物大腦神經的行為特征,從而進行大規(guī)模分布式并行信息處理的數(shù)學模型。ANN 由3部分組成:輸入項、包含激勵函數(shù)的隱含層和輸出項。其依據(jù)的主要思想是:通過對輸入輸出數(shù)據(jù)的訓練,建立能夠以任意精度逼近最佳輸出的非線性函數(shù)的模型[24]。
ANN以神經元作為基本處理單元。簡單的ANN是多輸入一輸出的,其表達式為:
式中:Ii為隱含層的值;yi為輸出值;xj(j=1,2,…,n) 為不同的輸入信號;θi為閾值;wji為從神經元j到神經元i的權值;f(·)為傳遞函數(shù)。
ANN 通過不斷地調整權值和閾值的大小,從而形成一種由輸入到隱含層的固有模式,然后通過傳遞函數(shù)由隱含層到輸出層。合理選擇隱含層數(shù)和神經元個數(shù)對ANN十分重要。
2.3 支持向量機支持向量機(SVM)通過引入核函數(shù),將原低維空間中非線性不可分的數(shù)據(jù)映射到高維特征空間。將兩類數(shù)據(jù)分類問題轉換為在高維特征空間中求最優(yōu)的線性分類面,即找到的分類平面與最靠近的數(shù)據(jù)點之間距離最大[25]。SVM 進行模擬預測主要思路:通過訓練期對輸入輸出數(shù)據(jù)的訓練,建立結構風險最小化的模型,然后給定輸入數(shù)據(jù)就可以得到相應的輸出值。通用的SVM 表達式如下:
式中:xi為輸入值;yi為與xi對應的輸出值;sgn(·)為符號函數(shù);為核函數(shù);為調控SVM確定最優(yōu)分類平面的參數(shù)。
核函數(shù)的選取是SVM 應用的關鍵之處。常用的核函數(shù)有:線性核函數(shù)、多項式核函數(shù)以及徑向基核函數(shù)等。核函數(shù)計算公式如下:
式中?函數(shù)將數(shù)據(jù)從低維空間映射到高維特征空間。
SVM 用于回歸分析,能夠對訓練數(shù)據(jù)進行學習,并實現(xiàn)其余數(shù)據(jù)的模擬預測,很好的解決了非線性擬合問題。對于樣本數(shù)據(jù),xi、yi分別為輸入變量和預測變量,n為數(shù)據(jù)量[26]。其回歸函數(shù)為:
式中:φ(x) 為輸入空間到高維空間的非線性映射;ω為權值向量;b為偏差。
SVM 基于結構風險最小化原理,而不是經驗風險最小化,提高了模型的泛化能力,求得風險最小的f(x):
式中:ε為損失系數(shù);ξi、為松弛變量;C為懲罰因子。
3.1 經驗模態(tài)分解經驗模態(tài)分解(EMD)簡言之就是將一段時間序列依據(jù)指定的條件分解為多個本征模態(tài)分量和一個趨勢項,使得每個本征模態(tài)分量的極值點和過零點數(shù)相等或者最多相差1,且相對于局部零值點是對稱的[27]。其具體步驟如下:
(1)輸入需要處理的原始時間序列y(t) ;
(2)找出序列中所有的極大值點和極小值點,根據(jù)三次樣條函數(shù)分別擬合出上包絡線yp(t)和下包絡線并求出上下包絡線的均值序列(t) 。計算公式如下:
(4)計算剩余項z1(t),以z1(t)作為原始序列,重復第(2)和第(3)步。如果最后的剩余項zi(t)滿足單調函數(shù)或至多具有一個局部極值點,此時zi(t)=q(t)為趨勢項,完成整個分解。否則繼續(xù)進行第(1)—(3)步。
式中r'i(t)為滿足本征模態(tài)分量條件的序列。
EMD分解原始序列y(t)得到n個本征模態(tài)分量和趨勢項??杀硎緸椋?/p>
式中:ri(t)為原始時間序列的本征模態(tài)分量,i=1,2,…,n;q(t)為原始時間序列的趨勢項。
3.2 集合經驗模態(tài)分解集合經驗模態(tài)分解(EEMD)與EMD 的目的是相同的,將原有的時間序列分解得到若干不同尺度的本征模態(tài)分量IMF 以及一個趨勢項。Huang 等[28]發(fā)現(xiàn)EMD 分解得到的某些本征模態(tài)分量IMF 包含了不止一個尺度下的信號,有些信號出現(xiàn)在不同的IMF 中。通過加入白噪音改進EMD 得到的EEMD 解決了模態(tài)混疊的問題,能夠更精確的分解得到不同尺度的本征模態(tài)分量和趨勢項。EEMD分解時間序列的步驟如下:
(1)輸入需要處理的原始時間序列y(t);
(2)產生m組白噪聲序列nj(t) (j=1,2,…,m),要求m組白噪聲的均值為零。分別加入原始序列中,序列變?yōu)閦j(t):
(3)對zj(t)進行EMD分解,得到n個本征模態(tài)分量rij(t)以及趨勢項qj(t):
(4)由第(3)步可以得到m組n個本征模態(tài)分量以及n個趨勢項,求各項均值為序列分解的最終分量;
可以看出,EEMD同樣將時間序列y(t)分為了若干的本征模態(tài)分量和趨勢項。公式如下:
式中:r'i(t)為原始時間序列的若干本征模態(tài)分量;q'(t)為原始時間序列的趨勢項。
3.3 小波分解小波分析(WD)是一種對時間序列在時域和頻域上進行逐級分析的方法[29]。小波是均值為0,長度有限的波形。滿足如下條件的函數(shù)為小波母函數(shù):
其中:ψ(ω)為小波母函數(shù);ψ(t)為在頻率ω處的Fourier變換;ψ*(ω)為ψ(ω)的共軛函數(shù)。
小波變換是在時間和頻率兩域對信號的局部特征進行表征。在R上平方可積函數(shù)構成的函數(shù)空間L2(R),f(t)∈L2(R)。則關于母小波函數(shù)ψ的連續(xù)小波變換公式為:
對于離散的信號,有離散小波變換,離散連續(xù)小波變換中的伸縮因子a和平移因子b得到,通常取。公式如下:
1989年,塔式算法融入小波變換中,實現(xiàn)信號的分解重構,即著名的Mallat 算法。分解和重構過程表達式為:
式中:cj為原始時間序列;h、g分別為低通、高通濾波器;cj-1、dj-1分別為分解尺度j-1下的低頻信號與高頻信號;分別為h、g的共軛反轉;的卷積;為對卷積的二次下抽樣;Ucj-1為對cj-1的二元上抽樣;U、D為二元上、下抽樣算子。
分解、重構后得到與原時間系列同長度的子序列,便于對時間序列進行預測。WD 方法應用于時間序列預測的步驟如下:(1)選擇合理的小波函數(shù),并選擇合適的時間尺度(或小波分解水平);(2)通過濾波器對時間序列進行分解,時間序列分解得到一個逼近信號序列和若干個細節(jié)信號序列;(3)利用小波重構方法可以得到與原序列長度相同子序列,最后進行時間序列模擬預測等工作。
結合以上模型方法,本文通過多站優(yōu)選MLR、ANN、SVM 模型,選取最優(yōu)的單一模型并基于此模型結合多種分解技術(EMD、EEMD、WD)構建混合預測模型,對比分析不同站點采用混合模型的徑流預測精度高低;探究變化環(huán)境下豐富輸入因子(大氣環(huán)流異常因子)對于混合模型徑流預測精度的影響。主要技術路線如圖1所示。
4.1 研究區(qū)域概況以渭河流域為研究區(qū)域,選取北道、林家村、華縣3 個水文站作為研究對象。
圖1 技術路線
其中,北道位于渭河上游,林家村位于渭河中游,華縣站位于渭河下游。北道位于甘肅省內,地處秦嶺北麓,集水面積1846 km2。渭河支流牛頭河匯集湯浴河、后川河等諸多小河于天水市北道區(qū)注入渭河。林家村水文站位于陜西省寶雞市,其集水面積30 661 km2,多年平均含沙量63.7 kg/m3,是渭河干流的控制站,屬于國家重點水文站。華縣站位于陜西省渭南市,集水面積106 498 km2,多年平均含沙量49.3 kg/m3,是國家重點水文站,同時也是渭河流域干流的控制站。渭河流域示意圖如圖2。
圖2 渭河流域
采用的數(shù)據(jù)為北道、林家村、華縣3 個水文站1960年1月至2010年12月徑流,3 個水文站所在區(qū)域內相應時間的降雨、蒸發(fā)以及大氣環(huán)流異常因子(北極濤動(Arctic Oscillation,AO)、太平洋十年濤動(Pacific Decadal Oscillation,PDO)、Nino3.4、太陽黑子)的月尺度數(shù)據(jù)。其中,以徑流作為預測目標,以降雨、蒸發(fā)和前期的徑流作為模型的輸入項,將大氣環(huán)流異常因子作為輸入項增加項。
本文收集的4 種具有代表性的大氣環(huán)流異常因子月尺度數(shù)據(jù),其中,AO 數(shù)據(jù)從美國國家海洋和大氣管理(National Oceanic and Atmospheric Administrator,NOAA)國家氣象數(shù)據(jù)中心獲取,網址為(http://www.ncdc.noaa.gov/teleconnections/ao.php);PDO 數(shù)據(jù)從東京氣候中心獲取,網址為(http://ds.data/jma.go.jp/tcc/products/elnino/decadal/pdo.html);Nino3.4 從NOAA 的地球系統(tǒng)研究實驗室獲得,網址為(http://www.esrl.noaa.gov/psd/data/correlation/nina34.data);太陽黑子數(shù)據(jù)從NOAA的國家地球物理數(shù)據(jù)中心獲得,網址為(https://www.esrl.noaa.gov/psd/gcos_wgsp/Timeseries/SUNSPOT/)。
由于降雨蒸發(fā)等因素對徑流的影響存在滯后效應,不同因子對徑流的影響時間存在差異,故本文基于相關系數(shù)對模型輸入項進行優(yōu)選。通過對于不同滯時下,各因子與目標因子的相關性計算,分析選取合理的因子,這對徑流的預測是至關重要的。由于在徑流預測中,考慮了滯時為一個月的徑流作為預測下一月徑流的輸入項,所以徑流預測的預見期為1個月。
4.2 徑流序列突變分析以及時間尺度特性分析采用Mann-Kendall 變異診斷法對各研究站點徑流序列進行趨勢分析以及變異點診斷,結果如圖3所示。
由圖3(a)北道站的UF曲線可知:自1970年9月起北道站徑流量有明顯的減少趨勢,且在1986年8月以后這種減少趨勢超過了0.05的顯著性水平,意味著北道站徑流量下降趨勢顯著。同時可以看到北道站UF與UB曲線交點位于1988年1—4月,不在信度線之間,即徑流量沒有發(fā)生顯著的突變。由圖3(b)林家村站的UF曲線可知:自1963年11月起林家村站徑流量有減少趨勢,在1989年11月后減少趨勢超過顯著性水平0.05的臨界線,并且自1992年10月開始趨勢出現(xiàn)更明顯的下降。林家村站UF與UB曲線交點位于1989年1月不在信度線之間,即徑流量于此月發(fā)生突變,但交點已很接近信度線。由圖3(c)華縣站的UF曲線可知:自1963年11月起華縣站徑流量有減少趨勢,在1980年7月至1984年11月之間有增加趨勢,之后徑流量為持續(xù)減少趨勢且減少趨勢超過顯著性水平0.05的信度線。華縣站UF與UB曲線交點位于1971年11月且交點位于信度線之間,即此月徑流量發(fā)生由多到少的突變。
較長時期的河川徑流變化受多種因素的影響,是氣候變化和人類活動因素(如修建大量水利工程等)綜合作用產生的結果[30]。渭河流域上的3個研究站點月徑流序列均有趨勢性變化,且華縣站月徑流序列發(fā)生了變異。渭河流域有大量灌區(qū)和水庫水電站,主要灌區(qū)集中在關中平原,關中九大灌區(qū)依渭河橫貫關中平原[31]。灌區(qū)灌溉用水量波動不大,且大多灌區(qū)建成時間較早。渭河流域的大中型水庫修建時間為1960—1980年代[32]。理論上認為模型可以在訓練期捕捉并學習到灌溉工程用水以及水庫運行等人類活動信息。選取此3個研究站點,以期探究不同混合模型對變化環(huán)境下徑流預測的效果及適用性,并且可以更真實的認識大氣環(huán)流異常對于分解后徑流子序列的影響。
圖3 各站徑流突變點檢驗UF、UB過程線
考慮滯后效應以及徑流變化可能對模型效果的影響,故劃分本文的訓練期為1960年1月至2000年12月,驗證期為2001年1月至2010年12月。在模型的訓練期,模型將捕捉并學習到氣候變化和人類活動因素的部分信息,并相應的體現(xiàn)在模型的參數(shù)變化上。在訓練期實現(xiàn)對于模型參數(shù)的確定,驗證期為預測驗證樣本進行徑流預測[33]。在徑流模型評價體系建立中,選取了納什系數(shù)和均方根誤差為衡量徑流預測精度的評價指標。
4.3 選用MLR、ANN、SVM 模型進行徑流預測本文基于相關系數(shù)選取合理的輸入因子,分別計算北道、林家村、華縣3個水文站徑流、降雨、蒸發(fā)與徑流的相關系數(shù)分析,結果如圖4。
圖4 基于相關系數(shù)的輸入因子優(yōu)選
基于相關系數(shù)優(yōu)選輸入因子(降雨、蒸發(fā)、前期的徑流),相關系數(shù)大于0.3即滿足95%置信度檢驗,滿足置信度檢驗的序列作為預選輸入項[25]。可以看出,華縣站徑流自滯時為1 時,相關系數(shù)大于0.3;降雨與徑流的滯時為0、1、2、7、8個月的相關系數(shù)絕對值均大于0.3;蒸發(fā)與徑流的滯時為1、2、3、7、8、9 個月的相關系數(shù)絕對值大于0.3。徑流預報模型輸入項需要高信息量,強相關性的因子。一些因子雖然與徑流相關性好,但它們的作用是重復的。這樣的因子并不能同其他因子配合共同對徑流起作用,結果可能使預測精度降低,故需要剔除[34]。本文通過窮舉法[35]獲得最佳的輸入項。如華縣站的最佳輸入因子有:滯時為1 的徑流,滯時為0、1、2 個月的降雨以及滯時為1、2、3、7 個月的蒸發(fā)。同理優(yōu)選得到北道站、林家村站的最佳輸入因子,然后以輸入因子、目標因子建立單一模型(MLR、ANN、SVM)。本文選用ANN 中的BP 三層神經網絡,模型訓練函數(shù)采用trainlm[36]。對于隱含層神經元個數(shù),本文在訓練過程預先選定一個較小的值,然后依次增加神經元數(shù),取驗證數(shù)據(jù)誤差達到最小時的隱層神經元數(shù)為最終隱層單元數(shù)。由于徑向基核函數(shù)具有穩(wěn)定性和去記憶性,故SVM 模型構建中,對比選用了徑向基核函數(shù),不敏感損失系數(shù)ε取0.001,采用網格尋優(yōu)方法確定懲罰因子C和核函數(shù)參數(shù)g[37]。在預測驗證期,將相應滯時下的前期徑流、降雨、蒸發(fā)作為輸入項,輸入已訓練好的模型,可得到預測徑流。基于不同單一模型的徑流預測結果,如表1所示。
通過表1 可以得出,對于不同站點,采用MLR、ANN、SVM 模型預測徑流精度有較大差距。林家村站體現(xiàn)尤為明顯,采用MLR 模型驗證期其納什系數(shù)僅為0.150,均方根誤差為0.733;采用ANN模型,驗證期其納什系數(shù)為0.332,均方根誤差為0.642。據(jù)調查可知,林家村站上游有大型灌區(qū)水庫,需要取水調水,人為因素對徑流的影響很強,使得林家村站徑流預測模型中,MLR 模型的泛化能力弱,相比采用描述線性關系的MLR 模型,采用模擬非線性關系的ANN 模型在處理徑流及其影響因素之間復雜模糊的關系中更具有優(yōu)勢,泛化能力更強,預測時能得到更高的精度。采用SVM 模型預測林家村站徑流,驗證期其納什系數(shù)為0.603,相較ANN 模擬結果提高了81.6%;相應的均方根誤差為0.495,降低了22.9%,其模擬效果優(yōu)于ANN的模擬效果,說明SVM相比ANN模型在預測徑流過程中更具有優(yōu)勢。從模型原理探究原因,ANN 模型通過訓練后建立了基于經驗風險最小化的網絡模型,存在著局部極小值以及具有不穩(wěn)定性等缺點,而SVM 模型引用核函數(shù),并且以結構風險最小化為原則,最終得到確定的唯一解,可以彌補ANN 的一些不足,所以模擬預測效果優(yōu)于ANN。北道站以及華縣站的徑流預測結果分析亦可得到相同的結論。
表1 基于MLR、ANN和SVM模型的徑流預測
4.4 基于混合模型的徑流預測本文對北道、華縣和林家村站的徑流分別采用EMD、EEMD 和WD技術處理。EMD分解方法自適應的將北道站、華縣站和林家村站分解得到的子序列個數(shù)分別為7、6、7。子序列分別包含各站徑流的高頻部分、低頻部分以及趨勢項。EEMD分解方法將北道站、華縣站、林家村站分解后,得到子序列數(shù)分別為7、8、8。由于徑流序列是離散的,所以采用離散小波分解徑流。對于WD方法,小波函數(shù)以及分解層數(shù)的選擇極為關鍵。通過對比優(yōu)選,本文選取的小波函數(shù)為Daubechies小波系中的db3[38],小波分解層數(shù)為5層,子序列有1個逼近信號(a5)和5個細節(jié)信號(d1、d2、d3、d4、d5)。
耦合不同分解方法與SVM 建立混合預測模型,即根據(jù)不同分解方法將徑流序列分解為多個子序列;將各子序列分別建立預測模型進行單獨模擬,最后將預測的各子序列值進行累加,從而得到最終的預測徑流序列。以降雨、蒸發(fā)以及前期徑流為輸入因子,考慮各因子與各徑流子序列的滯后性。對于各徑流子序列進行模擬預測時,采用窮舉法選取有效滯時下的影響因子作為模型的輸入因子。將預測精度最高的徑流子序列累加合成完整預測徑流序列。同樣,對于混合模型的預測徑流精度通過納什系數(shù)和均方根誤差來評價。基于混合模型的徑流預測結果見表2和圖5。
通過表2基于混合模型的徑流預測結果,可以看出不同研究站點的徑流序列經過不同的方法分解后,子序列通過SVM 預測合成后的徑流精度有明顯差別。EEMD-SVM 分解預測的精度要普遍高于EMD-SVM 分解預測的精度。北道、華縣、林家村徑流預測精度提高幅度分別為21.6%、1.9%、4.8%。EEMD 方法改進了EMD 方法所存在的模態(tài)混雜問題,在時間序列分解上更有優(yōu)勢。WD-SVM分解預測的精度要高于EEMD-SVM 分解預測精度。但是不同站點區(qū)別很大,EEMD-SVM 模型在華縣站徑流預測中泛化能力較弱,采用WD-SVM模型進行分解預報相較EEMD-SVM 模型徑流精度提高幅度為31.0%。而北道站,WD-SVM 與EEMD-SVM 徑流預測精度相差很小。姜浩等[39]的研究發(fā)現(xiàn):依據(jù)通量矩陣可得WD 方法在分解能力上優(yōu)于EMD。由以上結果可知,不同徑流序列對于分解方法有著不同的適應性,但一般來說預測精度高低符合如下表達式:WD-SVM>EEMD-SVM>EMD-SVM。
表2 基于混合模型的徑流預測
施能[40]的研究表明大氣環(huán)流是造成我國冬季氣候變化的重要原因之一。本文研究區(qū)域為渭河流域,海陸熱力的差異會給當?shù)貧夂蛟斐奢^大影響,同時可能會對當?shù)貜搅鳟a生影響。所以本文的輸入因子在降雨、蒸發(fā)、前期徑流的基礎上,考慮大氣環(huán)流異常因子(AO、PDO、Nino3.4、太陽黑子)作為徑流預測的輸入因子,并對增加大氣環(huán)流異常因子前后徑流預測精度進行了對比分析。原始徑流采用WD 方法分解得到子序列,通過相關系數(shù)計算大氣環(huán)流異常因子相對于各分解序列的滯時,優(yōu)選大氣環(huán)流異常因子作為輸入項,采用SVM 模型進行預測。表3 為基于WD-SVM 混合模型的大氣環(huán)流異常因子篩選結果。
如表3 所示,對于研究站點大氣環(huán)流異常因子的篩選結果可知:4 個大氣環(huán)流異常因子中,Nino3.4 對渭河徑流影響明顯,而AO、PDO 和太陽黑子的影響不明顯;北道、林家村和華縣站依次由西向東分布,季風對中國東部和中國西部的影響強度不同,且季風也會對Nino3.4 產生影響,故Nino3.4 對于不同水文站徑流的影響會不同[41];對北道站經小波分解后的徑流子序列d3 進行預測時,增加滯時為3、4、5、8、9、10 的Nino3.4 作為模型輸入項;累加分解預測后的徑流子序列,徑流精度得到了提高。針對小波分解后的各徑流子序列通過支持向量機進行預測時,支持向量機中懲罰因子C和核函數(shù)參數(shù)g結果如表4所示。
圖5 基于分解的混合模型及融合大氣環(huán)流異常因子的WD-SVM模型的徑流預測
最終徑流預測結果如表5 所示。驗證期混合模型以及融合大氣環(huán)流異常因子的WD-SVM 模型預測徑流與實測徑流對比情況見圖5。
由表5 基于WD-SVM 模型豐富輸入因子后的徑流預測,可知相比輸入因子只有降雨、徑流以及前期徑流,增加大氣環(huán)流異常因子對3站的徑流預測精度均有提高。北道站徑流預測中納什系數(shù)提高了0.25%,均方根誤差降低了0.78%;林家村站納什系數(shù)提高了0.23%,均方根誤差降低了1.0%;華縣站納什系數(shù)提高了0.7%,均方根誤差降低了16%。其次,對各研究站點徑流極大值點(10 個左右)的預測精度進行了分析。WD-SVM 混合模型融入大氣環(huán)流異常因子后,北道站其均方根誤差減小0.9%;林家村站其均方根誤差減小1.3%;華縣站其均方根誤差減小4.7%。結果印證基于分解技術下大氣環(huán)流異常因子對于徑流子序列有影響,豐富輸入因子可以提高徑流預測精度,尤其是對于徑流極大值點的預測精度提高更為明顯。
表3 基于WD-SVM混合模型的大氣環(huán)流異常因子篩選結果
表4 基于小波分解的徑流子序列支持向量機模型參數(shù)值
表5 基于WD-SVM混合模型融合大氣環(huán)流異常因子前后的徑流預測
綜上所述,豐富輸入因子后的WD-SVM 模型對研究站點的徑流預測精度均高于80%。林家村站和華縣站訓練期以及驗證期的徑流預測精度不低于北道站(水利工程影響較?。?。修建的大量水利工程并未對徑流預測造成較大影響。其他學者也做了不少有關渭河流域徑流模擬預測的研究,左德鵬等[42]基于SWAT 模型和SUFI-2 算法的研究結果中,林家村站和華縣站月徑流模擬結果為:林家村站校準期和驗證期納什系數(shù)為0.23 和0.75、華縣站校準期和驗證期納什系數(shù)為0.72 和0.78;Zhao 等[43]基于分布式水文模型DTVGM 得到咸陽站(咸陽站位于渭河干流,林家村站和華縣站之間)月模擬徑流和觀測徑流納什系數(shù)在0.6~0.7 之間;黎云云等[44]基于SWAT 模型的渭河流域徑流模擬中,林家村站校準期和驗證期納什系數(shù)為0.58 和0.80、華縣站校準期和驗證期納什系數(shù)為0.77 和0.82。而本文基于WD-SVM 混合模型徑流預測結果中,林家村站校準期和驗證期納什系數(shù)0.880 和0.868、華縣站校準期和驗證期納什系數(shù)為0.890 和0.855。通過單一預測模型和混合預測模型預測結果、以及與其他學者相關研究結果對比可知,基于分解技術來提高徑流預測精度的有效性以及變化環(huán)境下融合大氣環(huán)流異常因子會進一步提高徑流預測精度。
(1)本文采用MLR、ANN 和SVM 模型對渭河流域北道、林家村、華縣站徑流進行預測,結果表明,依次采用MLR、ANN、SVM 模型進行徑流預測,預測精度不斷提高。表明渭河流域上徑流與其影響因素之間復雜的非線性關系占較大比重,在對徑流預測過程中基于結構風險最小化的SVM 要優(yōu)于基于經驗風險最小化的ANN。
(2)混合模型的徑流預測研究中,3 種混合模型的徑流預測精度由高到低排序為:WD-SVM>EEMD-SVM>EMD-SVM。
(3)WD-SVM 混合模型中,預報因子加入大氣環(huán)流異常因子,各站點的模型預測精度均得到了提高,尤其是對于徑流極大值點的提升較為明顯。最終,北道、林家村、華縣站驗證期徑流的納什系數(shù)分別為0.818、0.868和0.855。
(4)在對渭河流域不同站點徑流預測研究中,相較于直接對徑流序列預測,基于分解的混合模型預測精度更高。采用徑流分解技術后可更好地呈現(xiàn)出徑流的周期、趨勢等本質特征,對提高徑流預測精度起促進作用。在接下來的研究工作中,將有3個突破點:①大氣環(huán)流異常因子對徑流影響機理的深入研究,實現(xiàn)作為預報因子的有效辨別和預處理;②徑流預測的模型進一步優(yōu)化,提高其學習、泛化能力;③分解方法的優(yōu)化,對徑流序列進行更為本質的剖析,提高徑流的預測精度。