劉彩花，景 浩，成一雄，牛 月

（太原理工大學(xué)水利科學(xué)與工程學(xué)院，山西 太原 030024）

國內(nèi)外統(tǒng)計(jì)資料顯示，地下滲漏引起大壩發(fā)生安全事故的概率為20%，失事原因?yàn)闈B流量增大、滲透坡降加大等[1]。我國1981年的調(diào)查結(jié)果表明，31.7%的大壩安全事故是由滲漏管涌引起的[2]。因此，對大壩進(jìn)行滲流監(jiān)測，通過建立恰當(dāng)模型來監(jiān)控大壩運(yùn)行狀態(tài)、預(yù)測大壩未來運(yùn)行性狀意義重大。

由于滲流序列體現(xiàn)出復(fù)雜的非線性、非平穩(wěn)特征，可能包含周期性、趨勢性及隨機(jī)波動(dòng)等成分，因此聯(lián)合運(yùn)用多種理論分析方法建立耦合模型取代單一模型預(yù)測成為研究趨勢。本文利用經(jīng)驗(yàn)?zāi)B(tài)分解法 （Empirical Mode Decomposition，EMD） 具有秉留原數(shù)據(jù)屬性的特點(diǎn)及時(shí)間序列ARI(Auto Regressive Integrated)模型具有預(yù)測精度較高與先進(jìn)統(tǒng)計(jì)性的優(yōu)勢，提出了EMD-ARI滲流耦合模型，將此模型應(yīng)用于汾河水庫左壩岸滲流量的預(yù)報(bào)中，并將預(yù)報(bào)結(jié)果與ARI單一時(shí)序模型的預(yù)報(bào)結(jié)果進(jìn)行比較和分析，以期為大壩安全運(yùn)行狀態(tài)監(jiān)控及模型構(gòu)建提供參考。

1 理論分析方法

1.1 EMD-經(jīng)驗(yàn)?zāi)B(tài)分解法

EMD[3]可將原始序列分解為若干階頻率不同的固有模態(tài)函數(shù) (Intrinsic Mode Function，IMF)和一階趨勢項(xiàng)，分解步驟為:

（1）找出原始序列x（t）全部極大值點(diǎn)，進(jìn)行3次樣條函數(shù)法插值，得各時(shí)刻上、下包絡(luò)線值xmax（t）、xmin（t） 及包絡(luò)線均值 m（t）

（2） 求類距平值函數(shù) h（t）

判斷h（t）是否為固有模態(tài)函數(shù)的依據(jù)為:①h（t）函數(shù)的極值點(diǎn)個(gè)數(shù)與跨零點(diǎn)個(gè)數(shù)相等或僅相差1；②m（t）恒等于零。滿足則h（t）為固有模態(tài)函數(shù)，否則將 h（t）作為原始序列重復(fù)步驟 （1）、（2），直至同時(shí)滿足①、②為止。得到的一階固有模態(tài)函數(shù)I1（t）及剩余值序列為

（3）以r1（t）為新序列，重復(fù)以上步驟，直至分解出所有階IMF分量及一階趨勢項(xiàng)r（t）。

1.2 求和自回歸模型ARI(p,d)

ARI[4,5]（p,d）模型是d次差分運(yùn)算與AR（p）模型的組合，表達(dá)式為:

式中，Xt為原序列；φ1,φ2,…,φp為自回歸系數(shù)；p 為自回歸模型階數(shù)；B為后移算子；B的次數(shù)表示后移期數(shù)；d為差分運(yùn)算次數(shù)；Φ(B)為自回歸系數(shù)多項(xiàng)式；at為零均值白噪聲序列。

2 汾河水庫E MD-ARI滲流模型開發(fā)

2.1 工程簡介

汾河水庫建于1958年，為大 (二)型水利樞紐，具有發(fā)電、防洪、灌溉等綜合效益。水庫于1990年建立了大壩滲流自動(dòng)觀測系統(tǒng)，監(jiān)測內(nèi)容包括左右壩基滲流、左壩岸滲流、古河床壩體繞滲、浸潤線等，其中，左壩岸滲流采用三角堰觀測，觀測頻率為每周一次，觀測設(shè)備完好。本文采用2010年7月2日～2012年11月17日汾河水庫左壩岸共125組滲流監(jiān)測數(shù)據(jù)進(jìn)行模型構(gòu)建；采用2012年11月24日～2012年12月29日汾河水庫左壩岸共6組滲流監(jiān)測數(shù)據(jù)檢驗(yàn)?zāi)Ｐ蛿M合結(jié)果與預(yù)報(bào)精度，部分滲流監(jiān)測數(shù)據(jù)如表1所示。

表1 汾河水庫左壩岸部分滲流量監(jiān)測數(shù)據(jù)L·s-1

2.2 EMD分解

基于MATLAB平臺(tái)對汾河水庫左壩岸滲流序列進(jìn)行EMD分解，分解結(jié)果如圖1所示。分解得4階IMF分量與趨勢項(xiàng)，經(jīng)Hilbert變換及倒數(shù)運(yùn)算后求得 4階 IMF分量的平均周期分別為 23.8、43.6、148、222.7天，表明滲流序列可能存在23.8、43.6、148、222.7天這4個(gè)波動(dòng)周期；1階趨勢項(xiàng)總體呈上升趨勢，表明汾河水庫左壩岸滲流量有隨時(shí)間推移增加的趨勢。

圖1 汾河水庫左壩岸滲流序列EMD分解結(jié)果

2.3 IMF分量及趨勢項(xiàng)擬合

應(yīng)用Fourier級(jí)數(shù)對4階IMF分量進(jìn)行周期擬合[6]，不考慮均值項(xiàng)的擬合結(jié)果為

應(yīng)用線性方程對趨勢項(xiàng)進(jìn)行擬合，擬合結(jié)果為

式中，t為滲流量序列時(shí)間序號(hào)(1，2，…，125)。

2.4 剩余值序列ARI模型建立

水庫滲流序列經(jīng)提取可能周期與趨勢項(xiàng)后得剩余值序列，對其進(jìn)行時(shí)序分析建立ARI模型的步驟[7]為:

（1）剩余值序列平穩(wěn)性檢驗(yàn)。ADF單位根檢驗(yàn)結(jié)果表明剩余值序列顯著非平穩(wěn)，一階差分后平穩(wěn)。

（2）差分序列的純隨機(jī)性檢驗(yàn)。采用LB檢驗(yàn)統(tǒng)計(jì)量對差分后的剩余值序列進(jìn)行純隨機(jī)性檢驗(yàn),檢驗(yàn)結(jié)果為LB統(tǒng)計(jì)量的P（Pr>ChiSq）值均小于0.05，說明該序列為非白噪聲序列，滿足時(shí)序建模條件。

（3）差分序列擬合AR模型。模型階數(shù)由BIC準(zhǔn)則確定，BIC值最小時(shí)相應(yīng)模型階數(shù)最優(yōu)，其表達(dá)式如下

式中，BIC為信息準(zhǔn)則函數(shù)；N為序列長度；p為模型階數(shù)；C為常數(shù)；為模型殘差方差的極大似然估計(jì)值。

根據(jù)BIC準(zhǔn)則確定最優(yōu)模型為AR（2），模型參數(shù)經(jīng)檢驗(yàn)極顯著，還原差分后得ARI（2，1）模型表達(dá)式為

（2） 模型殘差檢驗(yàn)。 ARI（2，1）模型殘差的檢驗(yàn)結(jié)果為LB檢驗(yàn)統(tǒng)計(jì)量的P（Pr>ChiSq）值均大于0.05，表明ARI(2，1)模型適合滲流序列。

2.5 汾河水庫EMD-ARI模型建立

汾河水庫EMD-ARI模型的結(jié)構(gòu)為

式中，yt為汾河水庫左壩岸滲流量序列；aj、bj為相應(yīng)各階IMF分量的傅氏系數(shù)；ωj為相應(yīng)各階IMF分量的角頻率；ci為趨勢項(xiàng)各項(xiàng)系數(shù)；ut為汾河水庫左壩岸滲流量剩余值序列；φk為時(shí)間序列ARI模型參數(shù)；ut-k為t-k時(shí)刻的剩余值序列；εt為模型殘差項(xiàng)。

基于SAS軟件平臺(tái)運(yùn)用Gauss-Newton[9]（高斯-牛頓）法對全部參數(shù)進(jìn)行重構(gòu)，所得EMD-ARI模型參數(shù)重構(gòu)的最小二乘估計(jì)值如表2、3、4所示。

表2 IMF分量傅氏系數(shù)

表3 趨勢項(xiàng)系數(shù)

表4 ARI模型參數(shù)

汾河水庫EMD-ARI模型的回歸方程為

3 汾河水庫ARI滲流模型開發(fā)

水庫左壩岸滲流序列只進(jìn)行時(shí)序分析，經(jīng)平穩(wěn)性與隨機(jī)性檢驗(yàn)后確定最優(yōu)模型為ARI(2,1)，對模型參數(shù)進(jìn)行條件最小二乘檢驗(yàn)，檢驗(yàn)結(jié)果如表5所示。

表5 條件最小二乘結(jié)果

由表5可知，t檢驗(yàn)統(tǒng)計(jì)量的P(Pr>|t|)值全部小于檢驗(yàn)水平0.01，模型參數(shù)均極顯著，將一階差分還原后得ARI(2,1)模型的表達(dá)式為

式中，yt為汾河水庫左壩岸滲流量序列；yt－1為t-1時(shí)刻滲流量序列；yt－2為 t-2時(shí)刻滲流量序列；yt－3為t-3時(shí)刻滲流量序列；εt為模型殘差項(xiàng)。

經(jīng)參數(shù)檢驗(yàn)后對ARI(2,1)模型進(jìn)行殘差檢驗(yàn)，檢驗(yàn)結(jié)果如表6所示。

從表6可以看出，ARI（2，1）模型殘差的LB檢驗(yàn)統(tǒng)計(jì)量的P（Pr>ChiSq）值均大于0.05，殘差序列為白噪聲，表明滲流序列可建立ARI（2,1)模型。

表6 ARI(2,1)模型殘差自相關(guān)檢驗(yàn)

表7 滲流模型預(yù)報(bào)誤差分析

圖2 EMD-ARI滲流模型擬合

圖3 ARI滲流模型擬合

4 汾河水庫滲流模型預(yù)報(bào)結(jié)果的對比分析

（1）EMD-ARI模型與ARI模型的擬合結(jié)果分別如圖2、3所示。從圖2和圖3可以看出，兩個(gè)模型的預(yù)報(bào)值與實(shí)測值的吻合度均較高、擬合效果均較好。

（2）EMD-ARI模型與ARI模型的預(yù)報(bào)結(jié)果如表7所示。由表7可知，①EMD-ARI模型其預(yù)報(bào)值與實(shí)測值的相對誤差E總體偏??；②EMD-ARI模型第6步預(yù)測的相對誤差E為7.56%，小于ARI模型的11.65%；④E<10時(shí)EMD-ARI模型的預(yù)報(bào)合格率為100%，高于ARI模型的83%；④E<15，E<20時(shí)兩模型的預(yù)報(bào)合格率均為100%，滿足相對誤差E小于20%的要求。

5 結(jié) 論

（1）汾河水庫左壩岸滲流序列經(jīng)EMD分解，得4階IMF分量與1階趨勢分量，4階IMF分量表征滲流序列可能具有的4個(gè)周期成分，1階趨勢分量表明滲流量隨時(shí)間推移增加。

（2）EMD-ARI模型與ARI模型的相對誤差呈先增后減再突變增趨勢，即相對誤差表現(xiàn)為跳躍特征，且第6步預(yù)測的相對誤差的突變增表明隨預(yù)測步長的增大模型預(yù)測精度降低。

（3）本文1步預(yù)測的時(shí)間間隔為7天，ARI模型第6步預(yù)測的相對誤差為11.65%，誤差較大，因此適用于進(jìn)行5步預(yù)測，即預(yù)測汾河水庫未來35天以內(nèi)的左壩岸滲流量；EMD-ARI模型第6步預(yù)測的相對誤差為7.56%，預(yù)測精度較高，可對汾河水庫左壩岸滲流量進(jìn)行6步預(yù)測。

（4）EMD-ARI滲流耦合模型的預(yù)報(bào)精度高于ARI滲流單一模型的預(yù)報(bào)精度，該法切實(shí)可行，可用于滲流預(yù)測。

［1］陳文燕,朱林,王文韜.大壩安全監(jiān)測的現(xiàn)狀與發(fā)展趨勢［J］.電力環(huán)境保護(hù),2009,25（6）:38-42.

［2］毛昶熙，段祥寶，李祖貽.滲流數(shù)值計(jì)算分析與程序應(yīng)用［M］.南京:河海大學(xué)出版社，1999.

［3］趙雪花，安莉莉，袁旭琦.基于HHT和R/S分析的黃河上游年徑流序列演變模式分析［J］.水電能源科學(xué)，2013，31（7）: 9-12.

［4］王燕.應(yīng)用時(shí)間序列分析［M］.北京:中國人民大學(xué)出版社，2005.

［5］王振龍.應(yīng)用時(shí)間序列分析（第二版）［M］.北京:中國統(tǒng)計(jì)出版社，2010.

［6］BOWERMAN B L,O’CONNELL R T.Forecasting and Time Seriers［M］.3rd ed.Brooks/Cole.Duxbury，United States America，1993.

［7］肖枝洪，郭明月.時(shí)間序列分析與SAS應(yīng)用 ［M］.武漢:武漢大學(xué)出版社，2009.

［8］楊金芳，翟永杰，王東風(fēng)，等.基于支持向量回歸的時(shí)間序列預(yù)測[J].中國電機(jī)工程學(xué)報(bào)，2005，25（17）: 110-114.

［9］高惠璇，等.SAS系統(tǒng) SAS/STAT軟件使用手冊 ［M］.北京:中國統(tǒng)計(jì)出版社，1997.