基于HEC-HMS模型的蘭江流域徑流預測

唐中楠， 楊國麗， 李 軍， 劉鵬霄

(河北建筑工程學院 河北省水質工程與水資源綜合利用重點實驗室, 河北 張家口 075000)

政府間氣候變化專門委員會(IPCC)的第五次評估報告中指出，從1880年以來全球平均表面溫度上升0.85 ℃，21世紀初已成為有記錄以來最炎熱的時期[1]。隨著溫室氣體排放的增加，全球范圍內輻射強迫度和溫度持續(xù)上升，溫度的上升導致更高的蒸散速率，這極大影響了區(qū)域水文過程和水文事件(即洪水和干旱)發(fā)生頻率[2]，因此研究水文水資源對氣候變化的響應具有重要現(xiàn)實意義。大氣環(huán)流模式(general circulation models, GCMs)是目前廣泛應用于大尺度氣候研究的方法[3]。然而GCMs輸出受到低分辨率和缺乏區(qū)域氣候數(shù)據(jù)限制，降尺度方法常用于彌補GCMs對區(qū)域氣候預測的不足[4]。SDSM (statistical down scaling model)統(tǒng)計降尺度模型耦合了天氣發(fā)生器和多元線性回歸，在氣候變化情景研究中得到廣泛應用[5]。

氣候變化下的水文過程響應逐漸成為水文水資源研究領域的熱點[6]，水文模型是模擬流域水循環(huán)過程及氣候和人類活動變化對水資源響應，揭示徑流變化規(guī)律的重要工具和方法[7]。其中，與SDSM模型對接應用較為廣泛的水文模型有SWAT模型[8]、VIC模型等[9]。HEC-HMS(hydrologic engineering center and hydrologic modeling system)模型是由美國陸軍工程師兵團研發(fā)的半分布式水文模型，因其相對簡單高效的結構和明確的物理機制而被國內外學者廣泛應用[10]。Abdessamed等[11]利用HEC-HMS模型驗證阿爾及利亞艾因賽夫拉市防洪混凝土墻防洪效果及其對洪峰流量的削減效果，發(fā)現(xiàn)市區(qū)擋土墻的存在顯著減少了洪水區(qū)面積。Mahmood等[12]通過對該模型量化分析，發(fā)現(xiàn)了過去50 a余非洲第四大淡水湖乍得湖入庫流量迅速下降，其中人類活動因素占總比例的66%，氣候因素僅為34%。袁玉等[13]探究了秦淮河流域洪水特征對景觀格局的響應，發(fā)現(xiàn)景觀豐富度較高的區(qū)域洪水危害性較小。田競等[14]基于修正Morris法優(yōu)化了模型參數(shù)敏感性分析過程，并驗證了模型在官山河流域的適用性。基于上述分析，國內外學者對HEC-HMS模型的研究主要集中于山洪預報預警和雨洪對下墊面變化的響應機制研究。然而，目前將HEC-HMS模型與SDSM模型耦合并應用于未來氣候變化下的徑流響應的研究還較少，加深耦合模型的應用理解對未來徑流變化預測有現(xiàn)實意義。

蘭江流域受降雨和地形等因素影響，是浙江省內防洪抗汛的重點流域。歷史上有記錄的該流域發(fā)生的特大洪水超過56次。20世紀以來受到氣候變化等因素影響，洪水發(fā)展到平均每4 a一次，最近一次特大洪水出現(xiàn)在2017年6月下旬，蘭溪站監(jiān)測水位達到31.86 m[15]。因此，本文針對蘭江流域構建HEC-HMS模型，并將其與SDSM模型結合，依托CanESM2氣候模式研究未來不同RCPs情景下的流域徑流變化趨勢，探究在全球變暖氣候變化背景下徑流的響應過程，以期為區(qū)域防洪工程和水資源管理提供重要參考和一定的理論支持。

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

1.1 研究區(qū)概況

蘭江流域位于錢塘江上游，浙江省中西部，水文上通常將馬金溪、衢江和蘭江統(tǒng)稱為蘭江。流域上游起源于安徽省休寧縣，干流全長約為303 km，流域境內總面積1.94×104km2。流域跨越東經118°40′—120°38′,北緯28°60′—29°30′之間，地形以丘陵和盆地為主，整體地勢呈現(xiàn)四周高、中間低的特點，氣候屬于東亞副熱帶季風區(qū)，年平均年降水約為1 200～1 700 mm，暴雨多集中在4—8月，多年平均溫度約在17～18 ℃間，流域徑流主要補給源是大氣降水，多年平均徑流量608 m3/s[16]。研究區(qū)域選取蘭江蘭溪水文站以上流域，區(qū)域地形和雨量站點分布見圖1。

圖1 蘭江流域雨量站及地形圖

1.2 數(shù)據(jù)來源及處理

1.2.1 地理空間數(shù)據(jù) 地理空間數(shù)據(jù)包括流域數(shù)字高程、土地利用和土壤分布數(shù)據(jù)。數(shù)字高程數(shù)據(jù)來自于地理空間數(shù)據(jù)云平臺(http:∥www.gscloud.cn)的30 m分辨率GDEMV2 30 M產品數(shù)據(jù)集。土地利用數(shù)據(jù)為中國科學院資源環(huán)境科學數(shù)據(jù)中心提供的1 km分辨率2018年土地利用現(xiàn)狀遙感監(jiān)測數(shù)據(jù)，通過蘭江流域矢量邊界提取(見圖2a)。土壤分布數(shù)據(jù)來自世界土壤數(shù)據(jù)庫HWSD (harmonized world soil database)數(shù)據(jù)集，邊界裁剪蘭江流域范圍見圖2b。

注：①LPk為黑色石灰薄層土；LVh1為簡育高活性淋溶土1；ACu為腐殖質低活強酸土；ACh1為簡育低活性強酸土1；RGc為石灰性疏松巖性土；CMo為鐵鋁性雛形土；CMx為艷色雛形土；RGd為不飽和疏松巖性土；ACh2為簡育低活性強酸土2；ALh為簡育高活性強酸土；ANh為簡育火山灰土；FLe為飽和沖積土；CMd為不飽和雛形土；ATc1為人為堆積土1；ATc2為人為堆積土2； ②ATc2比ATc1在土壤含沙量、壤土含量、土壤碳酸鈣降低，黏土含量上升；LVh2比LVh在土壤含沙量、壤土含量上升，黏土含量、土壤碳酸鈣下降；ACh2比ACh1在土壤含沙量上升，壤土含量、黏土含量和土壤碳酸鈣含量均下降。

1.2.2 水文氣象數(shù)據(jù) 構建HEC-HMS模型的水文氣象數(shù)據(jù)來自于歷年中華人民共和國水文年鑒第7卷第1冊，徑流數(shù)據(jù)選取2015—2018年的蘭江蘭溪站逐日平均流量表和洪水水文要素摘錄表，同時考慮到預測模型對下游蘭江水庫的安全性影響，故本著對水利工程最不利原則選取6場峰值流量大于3 000 m3/s的典型大規(guī)模洪水數(shù)據(jù)。降水量數(shù)據(jù)選取流域內15個雨量站點(表1)的降水量摘錄表。構建SDSM模型所需資料包括歷史氣象數(shù)據(jù)、NCEP(national centers for environment prediction)歷史逐日再分析數(shù)據(jù)和CanESM2大氣環(huán)流模式輸出數(shù)據(jù)。實測歷史氣象資料采用流域內東陽、衢州、開化、江山、金華和武義6個氣象站點1975—2005年逐日降水數(shù)據(jù)。NCEP歷史再分析數(shù)據(jù)選取1975—2005年日序列的26種大氣環(huán)流因子。CanESM2大氣環(huán)流模式數(shù)據(jù)來自加拿大氣候情景網，選取BOX_043 X_43 Y和BOX_044 X_43 Y網格下RCP2.6，RCP4.5和RCP8.53個排放情景的2030—2100逐日序列。

表1 蘭江流域雨量站點地理位置

1.3 研究方法

1.3.1 HEC-HMS水文模型 HEC-HMS模型包括流域模塊、控制模塊、氣象模塊和時間序列管理模塊4個部分[17]，流域模塊的降雨—徑流計算中又分為產流模塊、匯流模塊、基流模塊和河道匯流模塊，通過模塊式操作可以采用不同的計算方案模擬流域內的水文過程。

本研究通過HEC-GeoHMS 10.2對流域DEM數(shù)據(jù)進行處理，提取流域水系特征和地形參數(shù)并將流域劃分為18個子流域，利用泰森多邊形法計算每個子流域的雨量站權重，在此基礎上生成HEC-HMS工程文件(圖3)。產流模塊采用SCS-CN徑流曲線數(shù)法，率定參數(shù)包括降雨初損、CN值和不透水率；匯流模塊采用Snyder單位線法，參數(shù)包括流域滯時和峰值系數(shù)；基流模塊采用消退基流法，參數(shù)有初始基流、衰減常數(shù)和衰減閾值；河道匯流模塊采用馬斯京跟法，參數(shù)包括蓄量常數(shù)K，流量比重X和河段數(shù)。考慮到初始參數(shù)值輸入模型的模擬較差，故參數(shù)率定過程中采用手動試錯法結合模型內置的Nelder-Mead法及峰值加權均方根目標函數(shù)確定最優(yōu)參數(shù)[18]。驗證模型選用洪峰流量相對誤差(REQ)，洪量相對誤差(REW)，峰現(xiàn)時差(ΔT)，Nash效率系數(shù)(NSE)和相關性系數(shù)(R2)5個指標評價模型模擬結果，洪峰流量相對誤差REQ和洪量相對誤差REW誤差范圍應在20%以內，峰現(xiàn)時差ΔT誤差范圍應在3 h以內，Nash效率系數(shù)(NSE)應在0.6以上。REQ，REW和ΔT絕對值越低，NSE和R2越接近1代表模擬效果越好[19]。

根據(jù)上述方法構建蘭江流域HEC-HMS模型，洪水場次模擬選取4場暴雨洪水對模型參數(shù)進行率定，2場暴雨洪水對模型模擬效果進行驗證；逐日降雨徑流選取2015—2016年日徑流數(shù)據(jù)對模型參數(shù)進行率定，2017—2018年日徑流數(shù)據(jù)對模型模擬效果驗證。

圖3 蘭江流域劃分及模型結構

1.3.2 SDSM統(tǒng)計降尺度模型 SDSM模型是耦合多元線性回歸分析和隨機天氣發(fā)生器的統(tǒng)計降尺度模型，能夠將低分辨率的GCM大氣環(huán)流模式輸出數(shù)據(jù)轉換為站點或小尺度氣候要素的日序列，其核心為通過多元回歸分析方法建立預報因子與預報量之間的統(tǒng)計關系[20]。模型主要使用步驟包括：篩選預報因子、數(shù)據(jù)轉化及質量控制、率定驗證模型、天氣發(fā)生器和情景發(fā)生器，其中預報因子的選擇依據(jù)為預報因子與預報量間有強相關性和明確的物理過程、所選因子能夠被GCM準確模擬[21]。根據(jù)建立的統(tǒng)計關系率定及驗證模型后即可憑借于GCM輸出數(shù)據(jù)模擬預測未來日序列的氣象要素變化過程，其基本原理[22]如下所示：

(1)

(2)

式中：ωi為第i天降水概率;Pij為第i天第j個預報因子;α，β和γ為模型參數(shù);Ri為降水量;Ti為溫度變量;ei為誤差。

本研究采用加拿大氣候中心研發(fā)的CanESM2模式的未來輸出降水數(shù)據(jù)進行研究。林朝暉等[23]對17個GCMs模式在中國東部歷史降水模擬能力進行了分析，發(fā)現(xiàn)只有4種GCMs模式能夠同時還原小雨至暴雨雨量的年代增加，CanESM2模式是其中之一。Xuan等[24]分析了18個GCMs模式在浙江省區(qū)域氣候變化模擬能力，發(fā)現(xiàn)CanESM2模式是溫度模擬、降水模擬和風速模擬的最佳模式之一，溫度模擬絕對偏差小于1 ℃，降水偏差小于5%。因此，CanESM2氣候模式可以認為適用于浙江區(qū)域。

2 結果與分析

2.1 HEC-HMS模型結果

2.1.1 洪水模擬結果 洪水模擬結果見表2。由表2可以看出6場洪水的洪峰相對誤差和洪量相對誤差均小于20%，其中洪峰流量相對誤差絕對值均值為6.82%，洪量相對誤差絕對值均值為7.55%；峰現(xiàn)時差除20170411號洪水外均小于2 h；Nash效率系數(shù)為0.78～0.92，相關性系數(shù)均高于0.8。

從模擬與實測流量對比結果(圖4)中看出，模型模擬洪水過程與實測過程趨勢基本一致。率定期20170614，20170625，20180305這3場洪水Nash效率均接近于0.9，洪量誤差均接近或小于5 %，峰現(xiàn)表現(xiàn)趨于雨停峰現(xiàn)，模擬效果較好。洪號20180430，20170411洪水洪量相對誤差稍大，具體表現(xiàn)為局部模擬過程和實測洪水過程有所差距，其原因可能為實測流量數(shù)據(jù)為摘錄洪水數(shù)據(jù)線性插值處理后的結果，坦化后實際流量過程線導致局部模擬過程誤差增大；同時通過泰森多邊形法計算流域面降水量不能全部還原降水的空間和時間過程，從而進一步增大洪水過程線和峰現(xiàn)時間誤差。總體上看，率定期和驗證期共6場洪水的相對誤差、效率系數(shù)等大體符合要求，表明HEC-HMS模型可以用于蘭江流域洪水模擬。

表2 HEC-HMS模型洪水模擬結果

圖4 蘭江流域率定期(a-d)與驗證期(e-f)洪水過程模擬結果

2.1.2 日徑流模擬結果 表3為HEC-HMS模型日徑流模擬結果評價，率定期逐日徑流相關系數(shù)和效率系數(shù)均為0.79，相對誤差為6.11%；驗證期相關系數(shù)為0.74，效率系數(shù)為0.73，相對誤差Re為14.02%，總體上率定期模擬效果優(yōu)于驗證期，模擬值稍大于實測值。圖5為率定期和驗證期模擬逐日徑流和實測徑流過程對比，可以看出過程線擬合效果較好。模型對于極端徑流模擬值偏小，實測最高流量日為2017年6月25日，日平均徑流量達到12 900 m3/s，模型模擬結果僅為6 653.50 m3/s；2015年6月19日，實測日平均徑流量達到8 570 m3/s，模擬結果為4 907.50 m3/s；2015年6月30日，實測日平均徑流量達到6 290 m3/s，模擬結果為4 829.30 m3/s。綜合效率系數(shù)、相關系數(shù)和相對誤差來看，HEC-HMS模型在蘭江流域有較好的適用性，可以用于該地區(qū)的徑流模擬。

表3 模型模擬結果評價

圖5 蘭江流域驗證期和實測期徑流模擬過程

2.2 SDSM模型結果

2.2.1 SDSM模型率定與驗證 選取1975—1995年為SDSM模型率定期，1996—2005年為模型驗證期。采用平均解釋方差(E)、標準誤差(SE)和相對誤差評價模擬結果。考慮到HEC-HMS模型僅需輸入降水序列，故本研究中僅對6個氣象站點的日降水量進行模擬驗證，結果見表4。各站點標準誤差為0.35～0.40 mm，解釋方差為0.34%～0.47%，說明所選的預報因子能夠解釋降水量超過34%～47%的誤差。本研究中對降水量的模擬結果與郝麗娜等[25]對河西走廊降水模擬的解釋方差0.09%～0.64%,初祁等[26]對太湖流域降水模擬的解釋方差0.11%～0.26%和婁偉等[27]對涇河流域降水模擬中解釋方差為0.35%～0.45%等結果相似，同時年降水量統(tǒng)計結果表明模擬與實測值基本相近，因此認為構建的SDSM模型可適用于蘭江流域未來徑流模擬。

2.2.2 未來情景預測 將2030—2100年分為3個時段，即2030s(2030—2049年)，2050s(2050—2079年)和2080s(2080—2100年)，以1975—2005年作為基準期對比。應用已構建的SDSM模型輸入CanESM2模擬下RCP2.6，RCP4.5和RCP8.53種排放情景數(shù)據(jù)得到6個氣象站點未來降水逐日序列，并根據(jù)泰森多邊形法計算流域面降水量，結果見表5。

表4 SDSM模型降水模擬結果統(tǒng)計值

RCP2.6，RCP4.5和RCP8.5濃度路徑下流域降雨量較于基準期分別增幅-0.82%，6.18%和18.17%。RCP2.6情景下未來降水量呈現(xiàn)輕微下降趨勢，平均減幅為1.97 mm/10 a，21世紀末期較于基準期下降1.41 %；RCP4.5情景下降水量增長至頂峰后減少，增幅為14.84 mm/10 a，2050s時期平均多年降水量達到最高值1 860.75 mm，較于基準期增幅10.76 %，21世紀末期回落至1 767.36 mm，較于基準期增幅5.20%；RCP8.5情景下降水量持續(xù)快速增長，增幅為43.60 mm/10 a，21世紀中前期已明顯高于基準期，21世紀末期平均多年降水量達到2 338.70 mm，較于基準期上升39.21 %，相比于其他情景變化更為劇烈。總體上看，蘭江流域未來降水量呈上升趨勢，降水量增幅程度隨著輻射強迫度的上升而增大。

表5 蘭江流域各情景下未來降水預測

2.3 未來徑流模擬

將SDSM模型降尺度生成的未來較好序列輸入到已校準的HEC-HMS模型中，得到未來RCPs 3種排放情景下的流域徑流變化過程(圖6)。結合表6可知未來不同排放情景各時段下平均多年徑流量，對比基準期(2015—2018年)未來徑流量均成上升趨勢，RCP2.6，RCP4.5和RCP8.5 3種情景下多年平均徑流相較于基準期分別增幅為17%，26.22%和41.93%。

未來各時段徑流量變化規(guī)律與降水變化過程基本一致，RCP2.6情景下徑流量基本趨于穩(wěn)定，較于基準期增長幅度15%～20%，預測期內平均每10 a減少1.6 m3/s；RCP4.5情景下徑流量增長至峰值后緩慢下降，較于基準期增幅20%～32%，預測期內平均每10 a上升4.65 m3/s；RCP8.5情景下徑流量持續(xù)上升，從2030s時段較于基準期增幅17%持續(xù)上升至2080s時段增幅71%，預測其內平均每10 a上升49.49 m3/s，同時RCP8.5情景下各年份徑流量起伏較大，變化過于劇烈的徑流過程代表旱澇事件頻率和強度將明顯提升。

圖6 蘭江流域3種排放情景下未來徑流變化

表6 蘭江流域未來情景下流域年平均流量

圖7所示為各時段平均月徑流模擬結果。RCP8.5情景下各時段徑流除6，10，11月外均出現(xiàn)不同程度的增幅，其中豐水期(4—7月)漲幅顯著高于枯水期(9—12月)，且隨著時間的推移，豐水期和枯水期的徑流差也不斷擴大，至21世紀末期，豐水期已占全年徑流71.93%。

RCP4.5情景和RCP2.6情景下徑流年內分配均趨于穩(wěn)定，相較于基準期1—5月徑流出現(xiàn)明顯上升，7—9月和12月徑流出現(xiàn)了小幅上升，而6，10，11月徑流則出現(xiàn)了不同程度的下降，至21世紀末期2種情景汛期徑流占比分別達到60.06%和61.63%，枯水期徑流量則有減小的態(tài)勢。總體來看預測結果表明10月和11月以外的月份平均徑流均呈明顯增加趨勢，其中2—5月平均徑流量增幅明顯，豐水期徑流增幅程度均大于枯水期，表明未來流域徑流年內分配可能趨于參差，需要更加注意豐水期防洪水利工程。

圖7 蘭江流域多年月平均徑流基準期與實測期比較

3 討論與結論

3.1 結 論

(1) HEC-HMS模型在蘭江流域有較好的適用性。通過構建適用于蘭江流域的HEC-HMS模型并率定校準模型，場次洪水模擬平均效率系數(shù)0.86，平均相關系數(shù)0.89，平均洪峰相對誤差6.82%，平均洪量相對誤差7.56%；日徑流量模擬平均效率系數(shù)0.76，平均相關系數(shù)0.77，模型可以較為準確地還原研究區(qū)水文過程。

(2) 未來蘭江流域徑流量呈上升趨勢，增幅程度隨輻射強迫度同步增大。RCP2.6情景下徑流趨于穩(wěn)定，流域徑流量每10 a減少1.6 m3/s；RCP4.5情景下徑流量增加至峰值后減少，每10 a上升4.65 m3/s；RCP.8.5情景下徑流量持續(xù)攀升，每10 a徑流量上升49.49 m3/s。RCP8.5情景下年徑流量變化劇烈，旱澇事件發(fā)生概率和強度均有所增加。

(3) 未來徑流量年內分配將更趨于極端，汛期徑流量占全年比例將有所上升。RCP2.6，RCP4.5和RCP8.5 3種情景下豐水期徑流量上升幅度均明顯高于枯水期，增幅程度與輻射強迫度增大趨勢一致。預測期汛期徑流占比總體高于基準期，至21世紀末期RCP2.6，RCP4.5和RCP8.53種濃度路徑下分別由基準期的46.59%上升到60.06%，61.63%和71.93%。

3.2 討 論

本研究所用的CanESM2模式是國際耦合模式比較計劃第5階段CMIP5發(fā)布的氣候模式之一，CMIP5與第3階段CMIP3相比在氣候模擬能力上已有顯著提高[28]，但是GCMs模式仍然不能完全模擬大氣要素之間的互相作用，在未來氣候變化的預測中有極大的不確定性，尤其是影響因素較多且過程復雜的降水事件[29]，不同氣候模式由于機制原理、物理氣候結構、模式分辨率和輻射強迫度情景設計的較大差異，對未來氣候的模擬能力也不盡相同，通過比較多個GCMs模式能夠盡量減少不確定性。Zhang等[30]通過HadCM3，ECHAM5和CCSM3這3個模式研究發(fā)現(xiàn)錢塘江流域未來降水量極有可能呈現(xiàn)上升趨勢，Xia等[31]通過GFDL等3個GCMs模式預測錢塘江流域豐水期降水量上升而枯水期下降，以上多種氣候模式研究和本研究一致發(fā)現(xiàn)錢塘江流域和其上游蘭江流域的降水量在未來有上升趨勢，因此本研究基于SDSM模型分析蘭江流域未來降水量將有所上升的結論具有可信度。

降尺度方法也是區(qū)域尺度的未來氣候變化研究中的不確定性因素之一。SDSM統(tǒng)計降尺度方法雖然可以減少空間分辨率尺度上帶來的誤差，但其預測因子的選擇存在主觀性，模型情景生成與基線一致這一結構導致模擬結果年際降水量變化較大[32]。

本研究通過分析未來降水變化并結合水文模型研究了蘭江流域徑流的變化趨勢，但僅通過單出口觀測流量進行模型的率定與驗證難以說明流域整體的徑流變化得到了準確的還原，進一步分析不同區(qū)域的徑流量演變趨勢將是未來研究的重點。同時，土地利用類型的改變會對流域徑流產生重大影響[33]，充分考慮流域未來土地利用、土壤分布和高程的改變對徑流的影響是值得進一步研究的科學問題。