綜合多源不確定性的洪水概率預(yù)報(bào)試驗(yàn)

谷洪欽，劉開(kāi)磊，劉玉環(huán)，馬亞楠

(1. 國(guó)核電力規(guī)劃設(shè)計(jì)研究院有限公司，北京 100095；2. 淮河水利委員會(huì)水文局(信息中心)，安徽 蚌埠 233000；3. 河海大學(xué)水文水資源學(xué)院，江蘇 南京 210098)

在源遠(yuǎn)流長(zhǎng)的水文水資源學(xué)科發(fā)展歷史中，人們?cè)诜炙礄C(jī)制、坡面產(chǎn)匯流機(jī)理方面已經(jīng)形成了較多的經(jīng)驗(yàn)與知識(shí)儲(chǔ)備。然而，洪水預(yù)報(bào)的精準(zhǔn)程度受降雨輸入、土壤水文初始條件、參數(shù)概化能力、模型結(jié)構(gòu)等來(lái)源不確定性的影響顯著，導(dǎo)致傳統(tǒng)水文預(yù)報(bào)技術(shù)不能準(zhǔn)確反映流域產(chǎn)匯流物理機(jī)制，洪水預(yù)報(bào)結(jié)果往往明顯偏離實(shí)際，給實(shí)時(shí)洪水預(yù)報(bào)以及防洪減災(zāi)工作帶來(lái)極大的限制。

自19世紀(jì)40年代起，水文水資源學(xué)科便已出現(xiàn)有關(guān)誤差或水文不確定性對(duì)洪水預(yù)報(bào)的影響研究，尤其近年來(lái)針對(duì)來(lái)自降雨輸入、模型參數(shù)及模型結(jié)構(gòu)等來(lái)源不確定性影響的研究，取得了較多代表性成果。其中，對(duì)降雨、參數(shù)不確定性影響的研究多偏重于對(duì)誤差概率分布特征定量刻畫(huà)；對(duì)模型結(jié)構(gòu)不確定性的研究涉及對(duì)模型結(jié)構(gòu)本身的改進(jìn)以及最優(yōu)預(yù)報(bào)模型篩選策略、組合預(yù)報(bào)等方面。葉愛(ài)中等[1]于2007年提出由日降雨資料時(shí)間降尺度生成逐時(shí)降雨資料的方法，認(rèn)為所提出方法在保證逐日降水量總量一致的前提下，能在更細(xì)致的時(shí)間分辨率上提供可靠的降雨成果，對(duì)于降雨輸入誤差控制具有良好的借鑒價(jià)值；2016年梁忠民等[2]提出考慮降雨輸入不確定性的洪水概率預(yù)報(bào)方法，基于抽站法原理反推降雨輸入的概率分布，結(jié)合Monte-Carlo抽樣法驅(qū)動(dòng)確定性水文模型，實(shí)現(xiàn)對(duì)預(yù)報(bào)流量概率分布的估計(jì)；梁忠民等[3]利用貝葉斯理論對(duì)TOPMODEL模型的參數(shù)及預(yù)報(bào)不確定性問(wèn)題進(jìn)行探討，確認(rèn)了隨機(jī)抽樣方法對(duì)描述參數(shù)概率分布特征的有效性；賀新月等[4]針對(duì)水文模型的預(yù)報(bào)殘差序列的數(shù)值特征進(jìn)行分析，據(jù)此改進(jìn)傳統(tǒng)預(yù)報(bào)模型結(jié)構(gòu)，提高了融雪徑流的預(yù)報(bào)精度。

然而，現(xiàn)有研究往往將各種來(lái)源的不確定性割裂開(kāi)來(lái)，分別考慮對(duì)洪水預(yù)報(bào)結(jié)構(gòu)的影響，具有很大的理論與應(yīng)用局限性。實(shí)際上，降雨觀測(cè)與預(yù)報(bào)誤差、預(yù)報(bào)模型選擇、參數(shù)優(yōu)化等各方面的不確定性往往伴隨著一次洪水預(yù)報(bào)的整個(gè)過(guò)程，正是由于水文不確定性的客觀存在，導(dǎo)致人們不能夠完全精細(xì)掌握水文物理過(guò)程。本研究綜合考慮降雨輸入、參數(shù)、模型結(jié)構(gòu)3類典型的不確定性來(lái)源，在定量描述各要素不確定性程度的基礎(chǔ)上，融合各類型不確定性來(lái)源，實(shí)現(xiàn)洪水概率預(yù)報(bào)。

1 多源不確定性的量化

本研究是對(duì)近年來(lái)淮河流域降雨誤差分析、參數(shù)規(guī)律、多模型集合預(yù)報(bào)領(lǐng)域最新研究成果的集成，在綜合考慮降雨輸入、參數(shù)、模型結(jié)構(gòu)3個(gè)方面不確定性的前提下實(shí)現(xiàn)洪水概率預(yù)報(bào)，并結(jié)合歷史洪水，尤其是2020年大洪水驗(yàn)證概率預(yù)報(bào)方法及預(yù)報(bào)結(jié)果的可靠性。

1.1 降雨輸入

TIGGE(THORPEX interactive grand global ensemble)等降雨預(yù)報(bào)產(chǎn)品在實(shí)時(shí)洪水預(yù)報(bào)應(yīng)用中的輸入誤差較為顯著，且不同的產(chǎn)品、不同預(yù)見(jiàn)期數(shù)據(jù)的輸入誤差較大。本研究擬采用GBM(generalized Bayesian model)法分析降雨預(yù)報(bào)的不確定性。在已知預(yù)報(bào)降水量x的前提下，降雨預(yù)報(bào)誤差的廣義概率密度函數(shù)如下：

(1)

式中：x——預(yù)報(bào)降水量；y——降水量預(yù)報(bào)誤差值；δ——均方根誤差值；βx,0——降雨真值為0的概率，可根據(jù)樣本估計(jì)；fy(y|Y>0,X=x)——X=x時(shí)的降水量預(yù)報(bào)誤差，可用截尾正態(tài)分布等統(tǒng)計(jì)描述。

Cai等[5-6]于2018年、2019年在大別山區(qū)進(jìn)行降雨誤差分析，將Weibull分布函數(shù)用于降雨的先驗(yàn)概率密度估計(jì)，采用截?cái)嘈驼龖B(tài)分布密度函數(shù)估計(jì)降雨預(yù)報(bào)誤差隨機(jī)變量相對(duì)于降雨預(yù)報(bào)值的條件概率密度函數(shù)。研究中以Weibull分布函數(shù)描述βx,0的特征參數(shù)，同時(shí)給出了降雨預(yù)報(bào)誤差的概率密度函數(shù)表達(dá)式：

(2)

式中：c——常數(shù)項(xiàng)；μ、σ——Weibull分布函數(shù)的特征參數(shù)；Φ——標(biāo)準(zhǔn)正態(tài)分布函數(shù)。

1.2 模型參數(shù)

參考陳竹青等[7]于2019年所提出的水文模型參數(shù)概率分布特征描述方法，利用s場(chǎng)洪水資料，分別以每一場(chǎng)洪水資料獨(dú)立的率定水文模型i的參數(shù)xi，得到參數(shù)xi的次優(yōu)解集，據(jù)此估計(jì)xi所服從的概率分布。在本研究中，采用Beta分布函數(shù)描述[7]水文模型參數(shù)的概率分布特征。

1.3 模型結(jié)構(gòu)

本研究借助于BMA(Bayesian model averaging)評(píng)估多個(gè)模型相對(duì)最優(yōu)的概率，用于表征不同模型結(jié)構(gòu)上的不確定性程度相對(duì)大小，模型相對(duì)最優(yōu)概率記為ωj(j=1，2，…，J)

(3)

式中：j——最優(yōu)模型序列號(hào)；ωj=p(M=j|Tobs)——在給定數(shù)據(jù)集Tobs的條件下，模型j為最優(yōu)的后驗(yàn)概率。ωj數(shù)值的估算需要基于歷史實(shí)測(cè)、預(yù)報(bào)流量序列，構(gòu)建高斯混合模型、利用期望最大化等算法求解，介紹具體求解方法的文獻(xiàn)[8-9]較多,此處不再贅述。

2 研究區(qū)概況、數(shù)據(jù)及算法

2.1 研究區(qū)概況

史灌河是淮河南岸最大的支流，流經(jīng)河南省和安徽省，主要涉及河南省信陽(yáng)市、固始縣和安徽省六安市[10]。史灌河由上游的史河、灌河匯流而成，河流全長(zhǎng)261 km，流域面積6 889 km2。史灌河流域地形南高北低，南部最高峰太白峰海拔1 140 m，北部至淮河地面海拔一般為23 m左右。流域?qū)俦眮啛釒蚺瘻貛н^(guò)渡的季風(fēng)濕潤(rùn)區(qū)，多年平均降水量1 240 mm，多年平均水面蒸發(fā)量874 mm，多年平均水資源總量41.2億m3，其中地表水資源量36.1億m3。研究流域方位見(jiàn)圖1。

圖1 研究流域概況Fig.1 Conceptualization of experimental basin

2.2 數(shù)據(jù)介紹

史灌河流域蔣家集站2003年1月1日至2020年12月1日流量過(guò)程見(jiàn)圖2, 歷場(chǎng)洪水資料的時(shí)間步長(zhǎng)均為6 h。根據(jù)水文資料代表性的要求, 從現(xiàn)有資料中篩選出共16場(chǎng)大、中、小洪水；除2013年、2017年汛期流量、洪量均較小，且無(wú)典型洪水過(guò)程外，其余各年份均可篩選出一場(chǎng)洪水過(guò)程。具體的洪水場(chǎng)次信息見(jiàn)表1。其中，在2020年7月15—29日，史灌河流域產(chǎn)生特大暴雨，史灌河發(fā)生超保證水位洪水，下游支流白露河、淠河發(fā)生接近保證水位洪水。圖2中，各年度不同場(chǎng)次的洪水首尾相連繪制在同一張圖上，相鄰兩場(chǎng)洪水之間以藍(lán)色實(shí)線進(jìn)行分隔。

圖2 蔣家集站歷史洪水過(guò)程示意圖Fig.2 Schematic map of historical flood process in Jiangjiaji station

2.3 算法流程

綜合多源不確定性的洪水概率預(yù)報(bào)方法流程見(jiàn)圖3，其共分為6個(gè)步驟。

圖3 算法流程Fig.3 Flow chart of algorithm

步驟1：獲取面雨量預(yù)報(bào)誤差概率分布。一般要求流域內(nèi)降雨、蒸發(fā)、流量、水位數(shù)據(jù)資料條件較好，有至少10場(chǎng)典型洪水?dāng)?shù)據(jù)。

使用泰森多邊形法分別采用流域內(nèi)雨量站監(jiān)測(cè)數(shù)據(jù)和TIGGE氣象產(chǎn)品計(jì)算流域面雨量值[11-12]。其中，地面監(jiān)測(cè)相應(yīng)的面雨量作為真值，氣象產(chǎn)品預(yù)報(bào)降水量的面雨量數(shù)據(jù)為含誤差的序列，可參照式(1)(2)獲取預(yù)報(bào)誤差概率分布。

步驟2：構(gòu)建各模型參數(shù)的次優(yōu)解集，獲取水文模型參數(shù)的概率分布。分別選擇新安江模型[13]、BP神經(jīng)網(wǎng)絡(luò)模型作為預(yù)報(bào)模型，其中新安江模型以河網(wǎng)消退系數(shù)CS作為敏感參數(shù)率定獲得其次優(yōu)解集、最優(yōu)解；BP神經(jīng)網(wǎng)絡(luò)的參數(shù)即為其網(wǎng)絡(luò)結(jié)構(gòu)，直接利用全部場(chǎng)次洪水資料率定得到最優(yōu)的BP模型參數(shù)。參數(shù)優(yōu)選方面有眾多方法可供選擇，如SCE-UA、單純形法、客觀優(yōu)選法等。

表1 歷史洪水特征及預(yù)報(bào)精度指標(biāo)統(tǒng)計(jì)

本步驟推薦采用SCE-UA算法進(jìn)行參數(shù)率定，優(yōu)化目標(biāo)函數(shù)采用確定性系數(shù)指標(biāo)。SCE-UA算法、確定性系數(shù)指標(biāo)NSE的詳細(xì)內(nèi)容可參考文獻(xiàn)[14-15]。當(dāng)所優(yōu)化的參數(shù)值能夠使目標(biāo)函數(shù)確定性系數(shù)指標(biāo)最大時(shí)，認(rèn)為參數(shù)優(yōu)化完畢，當(dāng)前參數(shù)值即為所求。

步驟3：隨機(jī)生成面雨量、水文模型參數(shù)，驅(qū)動(dòng)各水文模型產(chǎn)生L1組初始預(yù)報(bào)流量過(guò)程，一般取L1大于50。確定各模型參數(shù)最優(yōu)解以及各模型為相對(duì)最優(yōu)的先驗(yàn)概率。將所有場(chǎng)次洪水資料劃分為率定期和驗(yàn)證期，以分別用于各水文模型的率定與驗(yàn)證，得到各水文模型參數(shù)應(yīng)用于所有場(chǎng)次的洪水預(yù)報(bào)時(shí)的綜合最優(yōu)參數(shù)；將已率定完畢的綜合最優(yōu)參數(shù)分別代入各水文模型，并利用BMA算法求解各模型為相對(duì)最優(yōu)的概率ω1、ω2、…、ωI。

步驟1～3分別用于獲取模型輸入、參數(shù)、結(jié)構(gòu)3種不同來(lái)源的概率分布特征，是對(duì)3種不確定性來(lái)源的不確定性程度的先驗(yàn)估計(jì)。綜合各模型在歷史洪水中表現(xiàn)的相對(duì)優(yōu)劣程度，估計(jì)模型為相對(duì)最優(yōu)的概率，并以BMA算法的結(jié)構(gòu)參數(shù)ω去表征概率值，體現(xiàn)了對(duì)模型結(jié)構(gòu)(模型選擇)不確定性的考量。

3 概率預(yù)報(bào)算法驗(yàn)證

將上述所選16場(chǎng)洪水中2003—2010年作為率定期，其8場(chǎng)洪水用于TIGGE降雨誤差分布特征參數(shù)估計(jì)，以及新安江模型、BP神經(jīng)網(wǎng)絡(luò)模型參數(shù)率定及兩模型相對(duì)最優(yōu)的權(quán)重的率定；2011—2020年作為驗(yàn)證期，其8場(chǎng)洪水用于驗(yàn)證降雨誤差分布特征參數(shù)準(zhǔn)確性，新安江模型、BP神經(jīng)網(wǎng)絡(luò)模型參數(shù)率定結(jié)果及相對(duì)權(quán)重值的合理性，以及檢驗(yàn)本研究所構(gòu)建洪水概率預(yù)報(bào)算法性能。其中，TIGGE產(chǎn)品中的中國(guó)氣象局CMA預(yù)報(bào)模式2008—2020年數(shù)據(jù)，時(shí)間分辨率6 h，空間分辨率0.5°×0.5°；對(duì)2008年之前以及2014年4月TIGGE降雨預(yù)報(bào)數(shù)據(jù)缺失的年份，以實(shí)測(cè)降水量替代。降雨預(yù)報(bào)誤差概率分布函數(shù)的特征參數(shù)沿用Cai等[6]在本流域前期的研究成果，即Weibull分布函數(shù)的參數(shù)值，分別為μ=11.537 3,σ2=295.635 3；率定得到新安江模型、BP神經(jīng)網(wǎng)絡(luò)模型相對(duì)最優(yōu)的概率分別為0.637和0.363。

把蔣家集站實(shí)測(cè)流量以及5%～95%概率預(yù)報(bào)結(jié)果的置信區(qū)間同時(shí)繪制在圖2上，用于展示2003—2020年期間的概率預(yù)報(bào)結(jié)果。橫坐標(biāo)上的紅色點(diǎn)用于分隔率定期與驗(yàn)證期。從圖2概率預(yù)報(bào)試驗(yàn)的流量過(guò)程可以看到，率定期、驗(yàn)證期洪水概率預(yù)報(bào)結(jié)果與實(shí)測(cè)結(jié)果的數(shù)量級(jí)、變化趨勢(shì)基本保持一致；且相對(duì)于率定期，概率預(yù)報(bào)結(jié)果的置信區(qū)間寬度在驗(yàn)證期并未有嚴(yán)重的發(fā)散現(xiàn)象。從圖2還可以看到，對(duì)于2003年、2016年、2020年汛期的預(yù)報(bào)洪水過(guò)程，其置信區(qū)間相對(duì)其他場(chǎng)次洪水的置信區(qū)間略有增寬，這在一定程度上表明概率預(yù)報(bào)對(duì)于較大量級(jí)洪水的不確定性程度會(huì)略有增大。上述研究結(jié)果證明，概率預(yù)報(bào)結(jié)果較為穩(wěn)定，且研究中所采用的降雨誤差分布特征參數(shù)、新安江模型及BP神經(jīng)網(wǎng)絡(luò)模型參數(shù)、模型相對(duì)最優(yōu)概率等參數(shù)相對(duì)合理。進(jìn)一步以表格形式統(tǒng)計(jì)概率預(yù)報(bào)結(jié)果的覆蓋率指標(biāo)CR(cover rate)及評(píng)估均值預(yù)報(bào)精度的確定性系數(shù)NSE(Nash-Sutcliffe efficiency coefficient)指標(biāo),見(jiàn)表1。

從表1可以看到，2003—2020年，2020年、2003年兩場(chǎng)洪水的洪峰流量最大，分別為4 350 m3/s和3 380 m3/s；相對(duì)地，2006年、2011年的洪峰流量?jī)H分別達(dá)到336 m3/s和308 m3/s。從新安江模型、BP神經(jīng)網(wǎng)絡(luò)模型的預(yù)報(bào)精度來(lái)看，模型預(yù)報(bào)結(jié)果的NSE指標(biāo)在驗(yàn)證期、率定期的均值分別為0.88、0.86和0.79、0.77，預(yù)報(bào)精度略有降低，降低幅度不大，證明了所采用降雨預(yù)報(bào)、模型參數(shù)的合理性。從概率預(yù)報(bào)結(jié)果的精度評(píng)價(jià)指標(biāo)分析，概率預(yù)報(bào)結(jié)果在2011年、2019年兩場(chǎng)較小量級(jí)的洪水中，90%置信區(qū)間預(yù)報(bào)結(jié)果的覆蓋率指標(biāo)CR僅分別達(dá)到67%、73%，同時(shí)也可以發(fā)現(xiàn)新安江模型、BP神經(jīng)網(wǎng)絡(luò)模型兩類單一模型的預(yù)報(bào)精度也均低于0.81，確定性模型、概率預(yù)報(bào)精度均較低，對(duì)于量級(jí)較小的洪水水文模型或概率預(yù)報(bào)模型的預(yù)報(bào)精度相對(duì)會(huì)有降低。但是從總體上來(lái)看，CR的均值達(dá)到了93%，即從2003—2020年的多場(chǎng)洪水過(guò)程綜合分析，概率預(yù)報(bào)結(jié)果的90%區(qū)間能夠覆蓋93%以上的實(shí)測(cè)洪水過(guò)程，所提供概率預(yù)報(bào)的可靠性程度相對(duì)較高，可以有效地避免偏大誤報(bào)或偏小漏報(bào)的可能性，提高洪水預(yù)報(bào)信息的參考價(jià)值。

從NSE指標(biāo)分析，率定期概率預(yù)報(bào)結(jié)果的NSE指標(biāo)均值為0.88，驗(yàn)證期的NSE指標(biāo)均值為0.87，洪水概率預(yù)報(bào)的均值預(yù)報(bào)結(jié)果在驗(yàn)證期預(yù)報(bào)精度與率定期相比變化更小，概率預(yù)報(bào)模型的預(yù)報(bào)結(jié)果相對(duì)更為穩(wěn)定。分析預(yù)報(bào)NSE的最大最小值可以看到，概率預(yù)報(bào)結(jié)果的NSE值變幅在0.75～0.92，而新安江、BP神經(jīng)網(wǎng)絡(luò)模型的NSE值變幅在0.66～0.95、0.69～0.91，概率預(yù)報(bào)的均值預(yù)報(bào)結(jié)果精度變化區(qū)間更小、且變化區(qū)間的下限比兩個(gè)單一模型的下限值更高，進(jìn)一步證明了概率預(yù)報(bào)結(jié)果所提供均值預(yù)報(bào)的可靠性較高。

4 結(jié) 語(yǔ)

本文提出了綜合考慮降雨、模型參數(shù)、模型結(jié)構(gòu)3類不確定因素的洪水概率預(yù)報(bào)方法。通過(guò)對(duì)預(yù)報(bào)洪水過(guò)程的對(duì)照，以及對(duì)預(yù)報(bào)結(jié)果覆蓋率CR、確定性系數(shù)NSE指標(biāo)的分析，可以發(fā)現(xiàn)，概率預(yù)報(bào)結(jié)果一定程度上依賴于其所采用的水文模型精度，如當(dāng)水文模型在2011年、2019年這樣的小水年預(yù)報(bào)精度較低時(shí)，概率預(yù)報(bào)的CR、NSE指標(biāo)相應(yīng)也處于較低水平。

從CR指標(biāo)上來(lái)看，概率預(yù)報(bào)結(jié)果所提供的90%置信區(qū)間預(yù)報(bào)結(jié)果能夠覆蓋93%的實(shí)測(cè)洪水過(guò)程，實(shí)測(cè)洪水過(guò)程基本上都落在90%置信區(qū)間之內(nèi)，根據(jù)洪水概率預(yù)報(bào)結(jié)果做出的防汛決策能夠充分考慮洪水可能的量級(jí)，對(duì)于準(zhǔn)確評(píng)估洪水發(fā)展情勢(shì)具有重要參考價(jià)值。概率預(yù)報(bào)結(jié)果相對(duì)新安江模型、BP神經(jīng)網(wǎng)絡(luò)模型的確定性預(yù)報(bào)結(jié)果，能以可靠的概率分布描述洪水可能的發(fā)展趨勢(shì)，降低洪水預(yù)報(bào)調(diào)度的誤判風(fēng)險(xiǎn)。從NSE指標(biāo)來(lái)看，洪水概率預(yù)報(bào)所提供均值預(yù)報(bào)結(jié)果，雖然并不能夠在每一場(chǎng)洪水中都達(dá)到預(yù)報(bào)精度高于單一模型，但是其預(yù)報(bào)精度總是高于預(yù)報(bào)精度相對(duì)較低的水文模型，NSE指標(biāo)更接近于水文模型中預(yù)報(bào)精度相對(duì)較高者，所提供均值預(yù)報(bào)結(jié)果的可靠性程度也較高。

考慮到不同的流域之間產(chǎn)匯流規(guī)律的巨大差異性，以及降雨數(shù)據(jù)產(chǎn)品、水文模型在流域的適應(yīng)性不同，不同來(lái)源不確定性因素的顯著性程度在流域間也存在較大差異，綜合多源不確定性的洪水概率預(yù)報(bào)應(yīng)當(dāng)根據(jù)洪水預(yù)報(bào)輸入數(shù)據(jù)、流域特征等因素來(lái)定制研究。除此之外，如何科學(xué)描述不確定性在輸入、模型運(yùn)算到結(jié)果展示過(guò)程的傳遞與變化，描述不同來(lái)源不確定性因素之間相互加強(qiáng)或削弱作用，并量化評(píng)估不同來(lái)源的不確定性因素對(duì)最終預(yù)報(bào)結(jié)果不確定性的貢獻(xiàn)程度，仍有待進(jìn)一步研究，這也是通過(guò)本研究可以推進(jìn)洪水預(yù)報(bào)技術(shù)發(fā)展的有價(jià)值的研究方向。