陸 琪
(陜西省商洛水文局,陜西 商洛 726000)
中長期水文預(yù)報是根據(jù)前期水文氣象資料,用成因分析和數(shù)理統(tǒng)計方法, 對未來較長時期的水文情勢進行科學(xué)的預(yù)測。目前中長期水文預(yù)報的方法主要有數(shù)理統(tǒng)計方法、天氣學(xué)方法和宇宙- 地球物理分析方法,而各種方法又有多種不同的思路,但是對于一個具體流域的中長期預(yù)報問題需要通過分析、嘗試和檢驗,本文運用方差分析周期預(yù)報分析成果,對商洛市乾佑河青泥灣站年徑流預(yù)報進行實踐檢驗。
一個水文要素隨時間變化的過程盡管多種多樣,但總可以把它看成是有限個具有不同周期的周期波相互疊加而成的過程,假定這個水文要素的序列為y1,y2,…,yn,并認為非周期性變化的偶然因素影響不大,其數(shù)學(xué)模型為:
式中:y(t)為水文要數(shù)序列;y(i)為第i 個周期波序列;x(t)為誤差項,i 為周期波數(shù)。
如果由實測水文要數(shù)資料能夠分析識別出水文要數(shù)所含的周期,而且這些周期在預(yù)報區(qū)間內(nèi)比較穩(wěn)定,就可以將周期波外延,然后進行線性疊加進行預(yù)報。
需要說明的是,由于影響水文要素長期變化因素很多,也很復(fù)雜,這里所指的周期不像自然界中其它自然現(xiàn)象一樣具有嚴格物理意義的周期,而只是概率意義上的周期,即為某一水文現(xiàn)象出現(xiàn)之后,經(jīng)過一定的時間間隔,再次出現(xiàn)這種現(xiàn)象的可能性較大而已。
由上述基本原理可知,周期預(yù)報,首先需要解決的問題是判斷某一個水文要素序列是否存在周期,如果存在周期,周期為多長;其次是根據(jù)實測數(shù)據(jù)分析得到的周期,它的可靠性如何(周期檢驗)。
從歷史資料的樣本序列中識別周期時,將序列分成若干組,當(dāng)分組組數(shù)等于客觀存在的周期長度時,組內(nèi)各個數(shù)據(jù)的差異小,而組間各個數(shù)據(jù)的差異大;如果組間差異顯著大于組內(nèi)差異時,序列就存在周期,其長度就是組間差異最大而組內(nèi)差異最小的分組組數(shù),可以通過選擇不同的置信度α,用F 檢驗來判斷組內(nèi)差異比組間差異小的顯著程度進行判斷。
以乾佑河青泥灣站1967 年~2012 年具有46 年的年徑流系列數(shù)據(jù)來分析。在識別樣本系列是否有周期和計算周期波時,計算量較大,為了提高工作效率,利用Excel VBA 編程計算,算法如下:
①計算
②從y1開始每b 項當(dāng)一個試驗周期,當(dāng)n 為偶數(shù)時b=2,3,…,n/2;當(dāng)n 為奇數(shù)時b=2,3,…
③對每一個試驗周期,處于同一位相的y 編為一組,求組中要素之和平方的均值,這些均值之和為Zg。計算方差比Fb=根據(jù)選定的信度α 為0.05(0.1)、自由度f1(b-1),f2(n-b)查F 表得出F,當(dāng)Fb>F 時進行第④步。否則該序列無周期波,結(jié)束識別。
④以Fb最大的試驗周期為選定周期,同位相y 值的平均值即為周期波在該位相的值x1,x2,…,xb。
⑤延長x 序列至n 項,即為x1,x2,…,xn。
⑥原序列y1,y2,…,yn分別減去x 的延長序列,得出一新的序列,作為源序列的余波,仍記為y1,y2,…,yn。
⑦對余波序列y1,y2,…,yn重復(fù)①~⑥步驟,求得若干顯著周期波。
⑧將上述所得周期波全部延長至n+1 年并疊加,n+1 年的疊加值即為預(yù)報值。
通過分析計算,在信度α 為0.05 時有2 個明顯的周期,分別為22 年和9 年;當(dāng)信度α 為0.1 時有1 個周期,周期長度是13 年,見表1,各周期數(shù)據(jù)成果見表2,年徑流3 個周期波過程圖見圖1。
圖1 乾佑河青泥灣站年徑流周期過程圖
表1 周期識別統(tǒng)計表
表2 周期數(shù)據(jù)成果表 單位:億m3
按照預(yù)報模型,將3 個周期波線性疊加得到預(yù)報過程。模擬結(jié)果和實測對比圖見圖2,分析期模擬預(yù)報結(jié)果見表3。根據(jù)水文情報預(yù)報規(guī)范(GB/T 22482-2008)第7.5.4 條,“對中長期預(yù)報的精度評定規(guī)定如下:a)對于水位(流量)的特征值定量預(yù)報,按多年同期實測變幅的10%、其它要素按多年同期實測變幅的20%作為許可誤差。b)定性預(yù)報的評定分合格和不合格兩個等級;當(dāng)預(yù)報與實況在同一量級時為合格,否則為不合格?!?/p>
以本資料多年同期實測變幅(9.56-0.880) 的10%,即8.680×10%=0.868 億m3作為許可誤差,對預(yù)報數(shù)值進行評定,得到模擬分析期預(yù)報年徑流合率為=76.1%,年徑流實測和預(yù)報擬合過程線圖見圖2,合格率計算表3,成果檢驗表4。
圖2 青泥灣站年徑流方差分析周期外延疊加預(yù)報擬合圖
表3 青泥灣站周期波計算及外延疊加預(yù)報結(jié)果表(計算期)單位:億m3
表4 數(shù)據(jù)分析期(1967 年~2012 年)成果檢驗表
利用分析成果對近年(2013 年~2019 年)乾佑河年徑流預(yù)報,預(yù)報計算見表5,在7 年年徑流預(yù)報中,合格5 年,不合格2年,檢驗期合格率71.4%,見統(tǒng)計表6。
表5 (2013 年~2019 年)年徑流預(yù)報計算表 單位:億m3
表6 (2013 年~2019 年)年徑流預(yù)報結(jié)果檢驗統(tǒng)計表
采用本預(yù)測方案成果進行預(yù)測,得到模擬分析期(1967 年~2012 年)青泥灣水文站年徑流預(yù)報的合格率為76.1%,應(yīng)用檢驗期(2013 年~2019)預(yù)報合格率71.4%,總的(1967 年~2019年)合格率75.5%。通過檢驗,近年來周期較穩(wěn)定,效果較好,可以將周期波外延,然后進行線性疊加進行數(shù)值預(yù)報。
1)應(yīng)用周期外延疊加法進行預(yù)報時,實際上假定了分析得到的周期在未來一段時間內(nèi)是保持不變的。但是,水文要素的變化絕不會按照固定的周期循環(huán)往復(fù)。因此,分析得到的周期只能作為一段時間內(nèi)的預(yù)報依據(jù),不能無限外推。
2)應(yīng)用方差分析識別周期是通過F 檢驗來判斷的,而應(yīng)用F 檢驗需要滿足兩個條件:一是每組數(shù)據(jù)均服從正太分布;二是各組總體的均方差相等。水文要素一般來講很難滿足這些條件,特別是第二個條件,因此,所得到的分析結(jié)果是近似的。
3)進行周期識別時,要選擇適當(dāng)?shù)男哦取H绻哦圈?選得太高,可能識別不出周期;信度α 的標準取得太低,容易出現(xiàn)偽周期的現(xiàn)象。
4)一般來說,所識別的周期個數(shù)越多,疊加結(jié)果與歷史資料的擬合越好。但是,從預(yù)報的角度出發(fā),周期的穩(wěn)定性應(yīng)放在首位。因此不要過分追求歷史資料的擬合上,周期數(shù)不宜太多,以2~3 個為宜。
方差分析得到的周期用于預(yù)報時,關(guān)鍵是周期的穩(wěn)定性,目前還沒有判別周期穩(wěn)定的有效方法。一般是根據(jù)歷史資料,檢驗一下要用以預(yù)報周期的近幾年來的情況,如果這個周期近幾年比較穩(wěn)定,且能通過標準較高的F 檢驗,那么用它來預(yù)報時把握大一些;反之如果發(fā)現(xiàn)這個周期的顯著性有較大的波動且有明顯的下降趨勢,那么在預(yù)報時不宜應(yīng)用這種周期進行外推預(yù)報。