洪水預(yù)報實時校正技術(shù)比較及應(yīng)用研究

周 夢，陳 華，郭富強，許崇育

(武漢大學(xué) 水資源與水電工程科學(xué)國家重點實驗室，武漢 430072)

0 引 言

洪水預(yù)報是根據(jù)前期和現(xiàn)時的水文、氣象等信息，揭示和預(yù)測洪水的發(fā)生及其變化過程的應(yīng)用科學(xué)技術(shù)。進行洪水預(yù)報的水文模型主要是采用觀測到的歷史水文資料和數(shù)據(jù)，確定模型參數(shù)的合理值，然后對未來的水文數(shù)據(jù)進行預(yù)測從而進行洪水預(yù)報。當(dāng)系統(tǒng)的參數(shù)隨時間變化時，普通的參數(shù)估計將無法得到確定的最優(yōu)值，實時校正就是追蹤系統(tǒng)參數(shù)的變動過程尋找合適的參數(shù)來進行系統(tǒng)模擬[1]。

在實時洪水預(yù)報系統(tǒng)中，每次預(yù)報做出之前，實時校正模型根據(jù)當(dāng)時的實測信息(新息)對預(yù)報模型的參數(shù)、狀態(tài)變量、輸入向量或預(yù)報值進行某種校正，使其更符合客觀實際。現(xiàn)在模型洪水預(yù)報實時校正涉及兩方面的研究內(nèi)容：一是建立洪水實時預(yù)報校正模型；二是采用有效的實時校正技術(shù)[2]。實時洪水預(yù)報校正技術(shù)有典型的水文模型流量預(yù)報實時校正算法[1]，其在多站點的大流域和流量陡漲陡落的流域局限性比較大；誤差自回歸校正算法[3]關(guān)鍵在于其回歸系數(shù)的確定，確定之后無需改動模型，計算量較少，但是其使用效果跟預(yù)報的殘差與前期殘差序列的一致性;衰減記憶最小二乘算法[4]其重點在于實時洪水遺忘因子的確定，但是洪水過程隨機性較大，可以針對不同洪水選取合適的遺忘因子；卡爾曼濾波算法[5]可以有效濾掉實測及預(yù)報噪聲，關(guān)鍵在于模型噪聲和量測噪聲方差矩陣的確定；基于貝葉斯理論的自適應(yīng)實時校正模型[6]其能夠充分利用先驗信息以及實測信息得到預(yù)報的后驗分布，關(guān)鍵在于誤差分布的估計以及模型參數(shù)的確定。這些方法都能根據(jù)實測信息及時更新預(yù)報數(shù)據(jù)，來提高預(yù)報精度。

現(xiàn)有研究中對于某一流域不同方法的適用性比較較少，并且少有比較各種方法對于不同洪水評價指標(biāo)的表現(xiàn)的文章，所以本文在建溪流域進行實例研究，在新安江水文模型模擬的基礎(chǔ)上比較了誤差自回歸實時校正模型、卡爾曼濾波算法以及基于貝葉斯理論的實時校正方法的適用性。

1 研究區(qū)域與數(shù)據(jù)

建溪流域位于福建省的北部，流域面積14 787 km2，流域平均降水量1 600～2 600 mm，是福建省鋒面雨高值區(qū)。如圖1所示，建溪七里街水文站以上有崇陽溪、南浦溪和松溪三大支流。本文主要研究七里街及其上游的流域，涉及的水文站點有麻沙站，崇安站，管潭站，松溪站，建陽站，水吉站，新廠站以及七里街站。

本研究用到的數(shù)據(jù)為2005-2014年的小時觀測數(shù)據(jù)，取2005-2013年為率定期，2014年作為檢驗期，時間步長取1 h，對場次洪水進行計算。校正主要是對檢驗期新安江模型的模擬數(shù)據(jù)進行校正，選取的場次洪水共有10場。

圖1 建溪流域圖Fig.1 Jianxi river basin

2 模型及方法介紹

2.1 新安江模型

新安江模型是由趙人俊教授領(lǐng)導(dǎo)的研究組在20世紀(jì)70年代提出的模擬徑流模型，在80年代由二水源模型發(fā)展成為了三水源模型，適用于濕潤地區(qū)與半濕潤地區(qū)。三水源新安江模型蒸散發(fā)計算采用三層蒸散發(fā)模型；產(chǎn)流計算采用蓄滿產(chǎn)流模型；用自由水蓄水庫結(jié)構(gòu)將總徑流劃分為地表徑流、壤中流和地下徑流三種；流域匯流計算采用線性水庫；河道匯流采用馬斯京根分段連續(xù)演算或滯后演算法[7]。本次率定模型所用的優(yōu)化算法是SCE-UA算法(Shuffled Complex Evolution),SCE-UA算法是一種全局優(yōu)化算法，適合新安江模型這類高維參數(shù)的全局優(yōu)化[8-10]。

新安江模型主要有15個參數(shù)需要率定[11],分段馬斯京根法對于不同的河段，參數(shù)k和x需要各自率定。現(xiàn)根據(jù)歷史蒸發(fā)、降雨以及徑流資料，可以用全局優(yōu)化方法率定參數(shù)，得到最優(yōu)的模擬模型。

2.2 AR算法

洪水預(yù)報誤差序列具有拖尾的特性，且為白噪聲，故洪水預(yù)報誤差序列為平穩(wěn)時間序列，可采用時間序列AR模型進行實時校正[12]。其主要是對模型殘差序列采用殘差自回歸估計式：

et+L=c1et+c2et-1+…+cpet-p+1+ξt+L

(1)

那么預(yù)報結(jié)果的校正式為：

(2)

其中：

et=Q(t)-QC(t)

2.3 Bayesian算法

(3)

假定天然徑流過程滿足自回歸滑動模型，先驗分布用自回歸滑動平均模型進行模擬；假定實測徑流和預(yù)報徑流服從線性關(guān)系，似然函數(shù)使用AR模型進行模擬。

2.4 卡爾曼濾波算法

卡爾曼濾波理論的基本思想即采用信號與噪聲的狀態(tài)空間模型，利用前一時刻的估計和觀測值來更新對狀態(tài)變量的估計，求現(xiàn)時刻的估計值[13]。其采用一個狀態(tài)方程和一個量測方程來完整的描述線性動態(tài)過程，前者描述系統(tǒng)狀態(tài)向量隨時間的變化規(guī)律，后者采用函數(shù)關(guān)系反映量測向量與系統(tǒng)狀態(tài)向量之間的依賴關(guān)系[14]。其遞推公式既可以得到濾波估計值，又可以得到誤差的方差陣，即可以完成自身的誤差分析[15]。其公式可以表示為：

狀態(tài)方程：

Xt+1=FtXt+BtUt+wt

(4)

量測方程：

Yt=HtXt+vt

(5)

式中：Xt為狀態(tài)向量；Yt為測量向量；Ut為系統(tǒng)的輸入；wt、vt分別為狀態(tài)噪聲和測量噪聲，且wt、vt為獨立一致的正態(tài)隨機變量；Bt、Ht為傳輸矩陣；Ft為系統(tǒng)從t到t+1時刻的n×n維狀態(tài)轉(zhuǎn)移矩陣[6]。

2.5 建溪流域洪水預(yù)報模型

本研究主要用到新安江模型及3種實時校正模型，如圖2所示，先根據(jù)歷史資料建立建溪流域新安江模型，確定匯流參數(shù)，以及區(qū)間的新安江模型參數(shù)。再用新安江模型進行檢驗期的流量模擬以及預(yù)報，預(yù)報時假設(shè)預(yù)報時刻以后流域降雨為零。預(yù)報值會分別經(jīng)過3種校正方法進行校正，校正前后的值會通過分段馬斯京根法匯流到七里街，區(qū)間用新安江模型進行預(yù)報，最后可以匯總得到七里街的流量預(yù)報。再次用相應(yīng)的方法進行校正。

圖2 建溪流域洪水預(yù)報模型結(jié)構(gòu)圖Fig.2 Technology roadmap

3 結(jié)果與討論

3.1 模型模擬結(jié)果

根據(jù)子流域區(qū)間降雨蒸發(fā)數(shù)據(jù)以及上游流量匯流，新安江模型可以對子流域流量數(shù)據(jù)進行模擬，然后采用分段馬斯京根法進行七里街匯流計算，區(qū)間繼續(xù)采用新安江模型進行模擬，得到的預(yù)報結(jié)果如表1所示。

表1 率定期及檢驗期新安江模型模擬結(jié)果Tab.1 Simulated results for calibration and validation period

表1列出了新安江模型模擬的結(jié)果，率定期為2005-2013年，檢驗期為14年。可以看到，率定期的洪水平均納西效率系數(shù)都達到了0.9以上，子流域的模擬情況普遍比整個大流域的模擬情況要好。檢驗期最小的納西效率系數(shù)都不到0.9，子流域的平均納西效率系數(shù)還是保持在0.9以上。對于各場洪水模擬的結(jié)果來看，80%的子流域洪水模擬納西效率系數(shù)都高于0.9，對于個別洪水，有的子流域降雨很少，洪峰不明顯，所以模擬效果不夠好，但是對整個流域的誤差影響較小。

3.2 各校正方法校正結(jié)果比較分析

子流域的預(yù)測及校正完成后，用校正后的預(yù)測數(shù)據(jù)，采用分段馬斯京根法進行七里街匯流計算，區(qū)間繼續(xù)采用新安江模型進行模擬預(yù)報，得到的預(yù)報結(jié)果如表2所示。

由表2可以看到，3種校正方法都相較于新安江模型模擬值明顯提高了精度，納西效率系數(shù)能增加0.15左右，水量平衡誤差也大大降低。隨著預(yù)見期的增加，預(yù)報精度降低，就研究洪水的平均水平來說，對于前幾個小時的預(yù)見期，貝葉斯方法的納西效率系數(shù)和水量平衡誤差這兩個指數(shù)都優(yōu)于其余兩種方法。在預(yù)見期增加后，AR模型校正的水量相對誤差比較小，在預(yù)見期大于9之后，卡爾曼濾波的納西效率系數(shù)較高，而貝葉斯方法的相對誤差較小。綜合來看，貝葉斯方法在納西效率系數(shù)和水量誤差方面在各場洪水的總體水平上上表現(xiàn)較好。

為了具體分析各種方法不同預(yù)見期的表現(xiàn)，現(xiàn)在對于每場洪水的納西效率系數(shù)進行統(tǒng)計，找出每場洪水最高和最低的納西效率系數(shù)如表3所示。

各種方法在不同預(yù)見期時納西效率系數(shù)最高和最低的頻率如圖3所示，由柱狀圖的分布情況可以明顯發(fā)現(xiàn)，貝葉斯方法對短預(yù)見期表現(xiàn)最好，模擬出納西效率系數(shù)比另外兩種方法好的頻率很高，但是隨著預(yù)見期的增加，貝葉斯方法的優(yōu)勢呈降低趨勢，在預(yù)見期超過6 h時已呈劣勢。相反地，卡爾曼濾波校正方法對短預(yù)見期的校正效果不如貝葉斯方法和AR方法，但是在預(yù)見期擴大時，開始顯現(xiàn)出明顯的優(yōu)勢。這說明卡爾曼濾波校正技術(shù)比貝葉斯方法更適合長預(yù)見期的洪水。而AR方法很穩(wěn)定，頻數(shù)分布比較均勻，對各預(yù)見期表現(xiàn)沒有明顯變化趨勢。綜合比較，在預(yù)見期較短時，可以選擇貝葉斯方法，而預(yù)見期變長之后，可以選擇卡爾曼濾波方法。

表2 不同預(yù)見期新安江模型與各校正方法的平均納西效率系數(shù)以及平均水量誤差Tab.2 Average NSE and RE of data forecasted by Xinanjiang model and three correction methods in different prediction period

表3 各校正方法校正之后的納西效率系數(shù)比較Tab.3 Comparison of NSE after correction of three correction methods

圖3 各種方法校正之后得到最高和最低納西效率系數(shù)的頻率比較Fig.3 Frequency comparison of the highest and lowest NSE after correction of three methods

圖4 各方法校正之后對應(yīng)不同洪峰流量的最高最低納西效率系數(shù)頻率比較Fig.4 Frequency comparison of maximum and minimum NSE corresponding to different peak flow after different methods correction[16] Maximum NSE Minimum NSE

除了考慮不同校正方法跟預(yù)見期的關(guān)系，還考慮到不同方法的校正效果可能跟洪峰流量有一定的關(guān)系。子流域以及整個流域各場洪水的洪峰流量以及6個預(yù)見期時各種校正方法校正效果如圖4所示，可見AR方法在6個預(yù)見期時表現(xiàn)較好，得到最優(yōu)納西效率系數(shù)的次數(shù)比較穩(wěn)定，貝葉斯方法對于低流量洪水表現(xiàn)較好，卡爾曼濾波比較適合高流量的洪水。

洪水預(yù)報中，對洪水洪峰的預(yù)測也十分重要，現(xiàn)在對不同預(yù)見期預(yù)測的洪水的洪峰誤差進行統(tǒng)計，得出比較圖如圖5所示。

圖5 各校正方法不同各預(yù)見期的洪峰誤差比較Fig.5 Comparison of flood peak errors at different forecast periods with different correction methods

由圖5可以看到，各個校正方法在前兩個預(yù)見期預(yù)報精度都比較高，洪峰誤差都在±3%以內(nèi)，在預(yù)見期延長時，誤差都在擴大，三種方法校正效果從圖中可以看到，卡爾曼濾波和貝葉斯方法，較于AR模型洪峰誤差較小，并且從分布區(qū)間和中位數(shù)來看，卡爾曼濾波表現(xiàn)更好，洪峰誤差較小。

圖6以20140725這場洪水為例，展示了不同預(yù)見期3種方法預(yù)報的結(jié)果，可以看到，隨著預(yù)見期的增加，預(yù)報值偏離實測值越來越遠。在3個預(yù)見期時，3種方法的校正曲線比較接近實測值，但是已經(jīng)不再光滑，在到12個預(yù)見期時，已經(jīng)偏離較遠了，這是因為12個預(yù)見期，在用新安江模型預(yù)報時，假設(shè)沒有降雨，12個預(yù)見期沒有降雨這在雨量方面的誤差過大，這個時候的預(yù)報精度已經(jīng)不能很好的保證了，只能做一個大致的參考。

4 結(jié) 語

建溪流域雨量豐沛，氣候濕潤，新安江模型本身就可以較好地模擬洪水過程。AR模型，卡爾曼濾波以及貝葉斯方法3種校正方法都能明顯提高流量模擬的精度，對于不同的評判指標(biāo)，3種校正方法校正效果有所不同。對于水量相對誤差，貝葉斯方法校正后效果最好，AR模型與卡爾曼濾波方法校正結(jié)果不相上下。AR模型和卡爾曼濾波方法校正之后的洪水過程納西效率系數(shù)保持很高，卡爾曼濾波方法對于長預(yù)見期或者流量較大的洪水表現(xiàn)得比較好，可以得到更好的納西效率系數(shù)，而貝葉斯方法對于流量較小的洪水表現(xiàn)得更好，AR模型對于各個量級的洪水都能有比較好的表現(xiàn)。綜上，3種方法都能在實際應(yīng)用中達到較為理想的效果，根據(jù)不同的預(yù)報需求，可以選擇相應(yīng)的方法進行校正。

圖6 不同預(yù)見期3種方法預(yù)測值與實測比較Fig.6 Comparison between predicted and measured values of three methods in different forecast periods