三峽水庫中長期徑流預(yù)測(cè)及不確定性分析研究

黃華平，酈于杰，王 棟，靳高陽

（1.中水珠江規(guī)劃勘測(cè)設(shè)計(jì)有限公司，廣州 510610；2.浙江省水利水電勘測(cè)設(shè)計(jì)院，杭州 310002；3.長江水利委員會(huì)水文局，武漢 430010）

0 引 言

中長期徑流預(yù)報(bào)是指基于水文現(xiàn)象的變化規(guī)律，以已知信息作為輸入條件，采用成因分析及數(shù)學(xué)建模等手段對(duì)預(yù)見期超過3 天的徑流過程進(jìn)行定量或定性的預(yù)測(cè)［1］。作為一種非工程措施，中長期徑流預(yù)報(bào)在水旱災(zāi)害防治、水庫調(diào)度運(yùn)用、水資源優(yōu)化管理及生態(tài)環(huán)境保護(hù)等方面均發(fā)揮了顯著作用［2］。

目前，關(guān)于中長期徑流預(yù)報(bào)已有較多研究，依據(jù)原理的不同，可劃分為以下兩類—基于水文循環(huán)過程驅(qū)動(dòng)和基于相關(guān)因子數(shù)據(jù)驅(qū)動(dòng)。前者一般是指將氣象要素預(yù)測(cè)數(shù)據(jù)作為輸入項(xiàng)，驅(qū)動(dòng)水文模型來對(duì)徑流過程進(jìn)行預(yù)測(cè)。如劉甜等［3］以CFS降水?dāng)?shù)據(jù)為輸入項(xiàng)，驅(qū)動(dòng)SWAT 模型對(duì)丹江口水庫入庫流量進(jìn)行了相關(guān)預(yù)測(cè)研究；Singh 等［4］采用ESP 方法獲取未來氣象數(shù)據(jù)，將其導(dǎo)入Topnet 模型，預(yù)測(cè)了新西蘭南島四個(gè)典型流域的徑流過程。這一類方法具有較為堅(jiān)實(shí)的理論基礎(chǔ)，但預(yù)測(cè)過程中需要搜集大量數(shù)據(jù)來驅(qū)動(dòng)模型運(yùn)轉(zhuǎn)，且降水預(yù)測(cè)的不確定性在水文模擬過程中將進(jìn)一步放大，導(dǎo)致預(yù)測(cè)結(jié)果往往差強(qiáng)人意，故而該類方法相關(guān)應(yīng)用較少。而基于相關(guān)因子數(shù)據(jù)驅(qū)動(dòng)是指在分析前期徑流或氣候要素與未來徑流過程相關(guān)性基礎(chǔ)上，建立統(tǒng)計(jì)模型來實(shí)現(xiàn)徑流中長期預(yù)測(cè)。該類研究在早期一般采用灰度理論、小波分析、方差分析等數(shù)理統(tǒng)計(jì)方法［5-7］。隨著人工智能與數(shù)據(jù)挖掘技術(shù)的快速發(fā)展，大量機(jī)器學(xué)習(xí)算法被應(yīng)用于中長期徑流預(yù)測(cè)研究中，包括神經(jīng)網(wǎng)絡(luò)模型、隨機(jī)森林算法、支持向量機(jī)及深度信念網(wǎng)絡(luò)等［8-11］。但上述研究多數(shù)僅提供了徑流的確定性預(yù)測(cè)結(jié)果，缺乏預(yù)測(cè)不確定性方面的考慮。

本文將LGB算法與HUP模型相耦合，應(yīng)用于三峽水庫逐月流量中長期預(yù)測(cè)中。研究首先采用相關(guān)系數(shù)檢驗(yàn)及逐步回歸法篩選出合理的預(yù)測(cè)因子；其次，以預(yù)測(cè)因子為輸入項(xiàng)，驅(qū)動(dòng)LGB 算法對(duì)逐月徑流過程進(jìn)行預(yù)測(cè)；最終，采用HUP 模型對(duì)預(yù)測(cè)結(jié)果的不確定性進(jìn)行相應(yīng)分析。

1 計(jì)算方法

1.1 輕量梯度提升樹

輕量梯度提升樹（Light Gradient Boosting Machine，LGB）是微軟于2017年提出的一種以梯度提升決策樹（Gradient Boosting Decision Tree，GBDT）算法為基礎(chǔ)的輕量化串行增強(qiáng)算法［12］。這類算法的基本思想是將訓(xùn)練過程階梯化，每一輪訓(xùn)練以擬合前一輪計(jì)算結(jié)果的殘差值為目標(biāo)，直到殘差值足夠小或迭代次數(shù)達(dá)到最大值為止。具體模型結(jié)構(gòu)如式（1）所示，每輪僅訓(xùn)練一棵決策樹來擬合前期決策樹之和的殘差，然后將新的決策樹加入前期決策樹中構(gòu)成新模型。

與GBDT 算法相比，LGB 算法優(yōu)勢(shì)主要反映在以下4 個(gè)方面：

（1）采用直方圖算法對(duì)輸入數(shù)據(jù)進(jìn)行離散化，減小噪聲影響的同時(shí)降低運(yùn)行所需內(nèi)存，模型訓(xùn)練速度顯著提升；

（2）采用梯度單側(cè)采樣法，保留大梯度樣本和部分小梯度樣本，減小訓(xùn)練數(shù)據(jù)且不改變樣本分布，提高了模型計(jì)算效率；

（3）通過將互斥特征停留在不同直方圖中構(gòu)造特征束，對(duì)特征束進(jìn)行排序，從而達(dá)到合并容忍范圍內(nèi)互斥特征的目的，大幅降低遍歷的特征數(shù)量；

（4）采用了一種更為高效的按葉生長策略，在分裂前會(huì)遍歷所有葉節(jié)點(diǎn)，選擇具有最大增益的節(jié)點(diǎn)分裂，可降低分裂誤差，得到更高的模擬精度。

1.2 水文不確定性處理器

Krzysztofwicz R［13］于1999年提出了貝葉斯預(yù)報(bào)系統(tǒng)（Bayesian Forecasting System，BFS），用于解決水文預(yù)報(bào)中的不確定性問題。該系統(tǒng)包含降水不確定性處理器（PUP）和水文不確定性處理器（HUP），分別用于處理降水輸入產(chǎn)生的不確定性和除降水外其他所有不確定性。其中，HUP模塊基本原理如下所述：

記H0為已知實(shí)測(cè)流量過程，Hn與Sn分別為不同預(yù)見期n對(duì)應(yīng)的實(shí)際流量和預(yù)報(bào)流量。依據(jù)貝葉斯原理及實(shí)測(cè)流量過程，可推斷出預(yù)報(bào)值Sn=sn條件下，實(shí)際流量Hn=hn的后驗(yàn)密度函數(shù)，如式（2）。

式中：h0為已知實(shí)測(cè)流量過程；g(hn|h0)為hn的先驗(yàn)分布，與h0相關(guān)；fn(sn|hn，h0)為sn關(guān)于hn的似然函數(shù)，反映了水文模型的預(yù)報(bào)能力。

假設(shè)每個(gè)目標(biāo)點(diǎn)Tl(l=1,2,3,…,m×n)的面積均為Δm×Δn,如果目標(biāo)點(diǎn)被覆蓋,則目標(biāo)點(diǎn)的聯(lián)合感知概率為1,覆蓋面積為Δm×Δn,否則為0,所以目標(biāo)點(diǎn)Tl的覆蓋面積可表示為Il×Δm×Δn;同樣,整個(gè)區(qū)域T的總面積為AS=(m×n)(Δm×Δn)。節(jié)點(diǎn)部署后節(jié)點(diǎn)所覆蓋的面積占部署區(qū)域總面積的比值稱為節(jié)點(diǎn)覆蓋率ψ[16],計(jì)算如下:

依據(jù)實(shí)測(cè)流量的后驗(yàn)密度函數(shù)φn(hn|sn，h0)，結(jié)合觀測(cè)數(shù)據(jù)及任意時(shí)刻的預(yù)報(bào)流量，可以估計(jì)相應(yīng)時(shí)刻實(shí)際值的條件分布函數(shù)，從而對(duì)預(yù)報(bào)結(jié)果的不確定性分析進(jìn)行相應(yīng)分析。

1.3 精度評(píng)價(jià)指標(biāo)

3種預(yù)測(cè)精度指標(biāo)（納什效率系數(shù)、相關(guān)系數(shù)及平均絕對(duì)誤差）及三種不確定性區(qū)間指標(biāo)（覆蓋率、相對(duì)帶寬及平均偏移程度）被用于評(píng)價(jià)模型確定性預(yù)測(cè)結(jié)果及不確定性置信區(qū)間，計(jì)算公式下如下：

（1）納什效率系數(shù)（NSE）。

（2）相關(guān)系數(shù)（CC）。

（3）平均絕對(duì)誤差（MAPE）。

（4）覆蓋率（CR）。

（5）相對(duì)帶寬（RB）。

（6）平均偏移程度（RD）。

2 實(shí) 例

本文搜集了三峽水庫1965-2016年逐月流量數(shù)據(jù)，考慮到不同月份徑流特性與影響成因存在差異，研究對(duì)12個(gè)月份分別建立LGB 預(yù)測(cè)模型，并將GBDT 算法作為對(duì)比參照來驗(yàn)證前者精度，在確定性預(yù)報(bào)基礎(chǔ)上，采用HUP 模型對(duì)預(yù)測(cè)結(jié)果不確定性進(jìn)行了相關(guān)分析。

2.1 預(yù)報(bào)因子篩選

將國家氣象局提供的130 項(xiàng)氣候系統(tǒng)指數(shù)（88 項(xiàng)大氣環(huán)流指數(shù)、26項(xiàng)海溫指數(shù)和16項(xiàng)其他指數(shù)（http：//cmdp.ncc-cma.net/Monitoring/cn_index_130.php）與前期徑流過程作為預(yù)測(cè)因子，考慮其與當(dāng)月徑流過程間的遙相關(guān)性，將徑流發(fā)生前12個(gè)月氣候系統(tǒng)指數(shù)與前一個(gè)月徑流量納入預(yù)測(cè)因子初選范疇，采用相關(guān)系數(shù)顯著性檢驗(yàn)及逐步回歸法對(duì)上述因子進(jìn)行篩選，具體步驟如下：

（1）搜集各月徑流過程對(duì)應(yīng)的所有相關(guān)變量，計(jì)算各相關(guān)變量與當(dāng)月徑流過程間的相關(guān)系數(shù)；

（2）依據(jù)相關(guān)系數(shù)顯著性檢驗(yàn)表，設(shè)置置信度為0.05，篩選出相關(guān)系數(shù)顯著的變量作為備選因子；

（3）采用逐步回歸法對(duì)備選因子進(jìn)行分析，依據(jù)不同因子的方差貢獻(xiàn)率，篩選出10 個(gè)獨(dú)立性強(qiáng)，相關(guān)系數(shù)高的備選因子作為最終預(yù)測(cè)因子。

其中，1月及8月選定的預(yù)測(cè)因子如表1所示。

表1 典型月份對(duì)應(yīng)預(yù)測(cè)因子統(tǒng)計(jì)表Tab.1 Selected predictors for two typical months（January and August）

2.2 模型率定與驗(yàn)證

以1965-2001年為率定期，2002-2016年為驗(yàn)證期，將逐月預(yù)測(cè)因子導(dǎo)入模型中，對(duì)三峽水庫逐月徑流過程進(jìn)行模擬預(yù)測(cè)。模擬過程中，為防止出現(xiàn)過擬合問題，研究采用留一交叉驗(yàn)證法對(duì)模型參數(shù)進(jìn)行充分率定，并將實(shí)測(cè)系列與模擬系列間MAPE值最小作為率定目標(biāo)，預(yù)測(cè)結(jié)果如圖1、2及表2所示。

圖1 三峽水庫1965-2016年逐月流量過程對(duì)比Fig.1 The observed and simulated monthly runoff series from 1965 to 2016 for the Three Gorges Reservoir

表2 模型率定期及驗(yàn)證期精度評(píng)價(jià)Tab.2 Performance Indices for the calibration and validation periods

由表2 可知，率定期內(nèi)，GBDT 和LGB 模型關(guān)于三峽水庫逐月流量的模擬結(jié)果精度較為接近，其中MAPE值分別為16.1%和15.7%，CC值均為0.91，NSE值分別為0.84 和0.86；而驗(yàn)證期內(nèi)，除CC值較為接近外，其他兩種指標(biāo)差異顯著增大，其中MAPE值分別為25.8%和23.0%，NSE值分別為0.71 和0.75。整體來看，兩種模型均具有較好的預(yù)測(cè)精度，與GBDT 模型相比，LGB 模型在率定期雖然與前者精度差異不大，但在驗(yàn)證期精度明顯更高。

圖1 提供了三峽水庫實(shí)測(cè)流量與模擬流量過程對(duì)比，其結(jié)果表明：率定期內(nèi)，豐水年的模擬值偏小，尤其對(duì)于1997 和1998年兩場特大洪水而言，汛期模擬值與實(shí)測(cè)值間絕對(duì)偏差超過30%，而對(duì)于平水年和枯水年的模擬結(jié)果略偏大；驗(yàn)證期內(nèi)，兩種模型模擬結(jié)果整體大于實(shí)測(cè)系列，但偏差程度不大。圖2提供兩種模型模擬結(jié)果與實(shí)測(cè)系列間散點(diǎn)分布圖，依據(jù)其特征可以得出同樣的結(jié)論，率定期內(nèi)，兩種模型位于下半部分的散點(diǎn)均勻地分布在45°線兩側(cè)，而上半部分散點(diǎn)則略低于45°線，驗(yàn)證期下半部分散點(diǎn)分布特征與率定期類似，而上半部分散點(diǎn)則略高于45°線。對(duì)比兩種模型結(jié)果可發(fā)現(xiàn)，圖1 與圖2 均表明：與GBDT模型相比，LGB模型對(duì)于三峽水庫逐月流量過程的模擬性能更佳，特別是對(duì)汛期而言，其模擬結(jié)果與實(shí)測(cè)值的偏差程度更小，對(duì)應(yīng)散點(diǎn)也更為接近45°線。

圖2 三峽水庫實(shí)測(cè)與模擬月流量散點(diǎn)圖Fig.2 Scatter plots of observed and predicted monthly runoff for the Three Gorges Reservoir

2.3 不確定性分析

HUP 模型處理過程中，需先采用正態(tài)分位數(shù)變換將實(shí)測(cè)值與預(yù)測(cè)值分別轉(zhuǎn)換至正態(tài)空間中，然后依據(jù)貝葉斯理論與線性似然假設(shè)來推求實(shí)際流量的后驗(yàn)概率密度函數(shù)，從而對(duì)預(yù)報(bào)結(jié)果進(jìn)行不確定性分析。與上述預(yù)測(cè)過程類似，研究將1965-2001年作為率定期，2002-2016年作為驗(yàn)證期，來對(duì)HUP 模型后處理精度進(jìn)行相應(yīng)驗(yàn)證。

依據(jù)HUP模型確定的實(shí)際流量后驗(yàn)概率密度分布，可以提取逐月徑流過程對(duì)應(yīng)的50%分位數(shù)流量值（Q50值）。表3和圖3 分別提供了Q50 值的精度指標(biāo)及其與實(shí)測(cè)過程的對(duì)比。兩者結(jié)果表明：與LGB 模型預(yù)報(bào)結(jié)果相比，Q50 值與實(shí)測(cè)系列間擬合程度更佳，各精度指標(biāo)在率定期及驗(yàn)證期也顯著優(yōu)于LGB 模型結(jié)果。為進(jìn)一步分析兩者精度差異的年內(nèi)變化特征，研究還計(jì)算了驗(yàn)證期內(nèi)實(shí)測(cè)值與模擬值間的相對(duì)誤差，如圖4 所示。由圖4，不難看出，兩者在枯水期（1-4月及12月）相對(duì)誤差普遍偏小，大體上處于-35%～35%之間，而在汛期，模擬值較實(shí)測(cè)值明顯偏大，在部分年份的8-10月內(nèi)，相對(duì)誤差甚至達(dá)到了90%，而對(duì)比圖4（a）與圖4（b），不難看出與LGB 模型預(yù)測(cè)結(jié)果相比，Q50 值在大部分月份的預(yù)測(cè)精度呈現(xiàn)一定程度提高，特別是對(duì)汛期而言，改善效果更為顯著。

圖3 三峽水庫1965-2016年逐月流量觀測(cè)系列與Q50系列對(duì)比圖Fig.3 The observed runoff and Q50 values series for the Three Gorges Reservoir

圖4 驗(yàn)證期逐月流量模擬結(jié)果相對(duì)誤差圖Fig.4 The relative error of predicted results for the validation period

表3 HUP模型率定期及驗(yàn)證期精度評(píng)價(jià)（Q50預(yù)報(bào)值）Tab.3 Performance Indices for the calibration and validation periods（Q50 predictions）

除Q50 預(yù)測(cè)值外，研究還提取了預(yù)測(cè)結(jié)果的90%不確定性區(qū)間，其中驗(yàn)證期不確定性區(qū)間如圖5 所示。為定量評(píng)估預(yù)測(cè)結(jié)果的不確定性，研究選取了覆蓋率（CR），相對(duì)帶寬（RB）和平均偏移（RD）幅度分別對(duì)率定期及驗(yàn)證期的不確定性區(qū)間進(jìn)行相應(yīng)分析，對(duì)應(yīng)指標(biāo)結(jié)果如表4 所示。由表4 及圖5 結(jié)果，不難看出提取的90%不確定性區(qū)間覆蓋了率定期與驗(yàn)證期絕大多數(shù)的實(shí)測(cè)點(diǎn)據(jù)，兩者覆蓋率分別達(dá)到0.86 和0.83，均較為接近0.9，說明提取的置信區(qū)間是較為可靠的。而相對(duì)帶寬和平均偏移幅度結(jié)果表明，與率定期相比，驗(yàn)證期的相對(duì)帶寬和平均偏移幅度更大，說明后者的預(yù)測(cè)精度要差于前者，且對(duì)應(yīng)不確定性也更大。

圖5 三峽水庫月流量不確定性分析結(jié)果Fig.5 Uncertainty analysis of the predicted monthly runoff sequences of the Three Gorges Reservoir

表4 不確定性區(qū)間精度評(píng)價(jià)Tab.4 Performance Indices for uncertainty confidence intervals

3 結(jié) 論

（1）預(yù)測(cè)精度指標(biāo)結(jié)果表明，GBDT 模型與LGB 模型對(duì)于三峽水庫逐月流量過程的預(yù)測(cè)結(jié)果在率定期和驗(yàn)證期內(nèi)均取得較好效果，但后者的精度要顯著優(yōu)于前者。

（2）采用HUP 模型對(duì)LGB 模型預(yù)測(cè)結(jié)果進(jìn)行后處理，提取的Q50 值系列對(duì)應(yīng)的MAPE、CC及NSE值均優(yōu)于LGB 模型預(yù)測(cè)結(jié)果，說明HUP模型具有一定程度的校正能力。

（3）90%不確定性區(qū)間在率定期及驗(yàn)證期的覆蓋率分別為0.86和0.83，較為接近0.9，且相對(duì)帶寬和平均偏移幅度較小，說明基于HUP模型提取的90%不確定性區(qū)間是較為合理可靠的?！?/p>