基于SWAT模型和降水隨機模擬的徑流預測

張金萍，王宇昊

（1.鄭州大學水利科學與工程學院，鄭州450001；2.鄭州大學黃河生態保護與區域協調發展研究院，鄭州450001）

SWAT 模型基于物理過程建立，且能夠充分反映流域內降水和下墊面要素空間變化對水文過程的影響，是研究流域水文循環的有效工具，是探究水文規律和機理的有效手段［1，2］。SWAT 模型在國內外應用廣泛，包括對水質、徑流、泥沙等方面的研究，其中最為常見的就是徑流模擬［3］。例如王博威等依據四期土地利用數據構建SWAT 模型，從而分析土地利用變化對流域徑流及其空間分布的影響［4］；廖亞一等將蒙特卡洛隨機采樣引入SWAT 模型，探討氣象數據輸入不確定性對徑流模擬結果的影響［5］；這些研究主要集中于通過構建SWAT 模型探究土地利用變化、人類活動、氣候變化等因素對徑流的影響［6，7］，但降水是氣象數據中對徑流模擬結果影響最大的因素，降水的總量和年內分布均會對徑流過程造成較大的影響［8］，而對不同量級和年內分布降水情景相應徑流過程的研究較少。

降水模擬包括降水概率的模擬和降水量的模擬，目前國內外運用比較多的降水概率模型主要有一階Markov 鏈、高階Markov鏈、間隔長度和半經驗分布法4種，而分布模型廣泛用于降水量模擬［9］。其中耦合Markov 鏈-Gamma 分布的隨機模擬模型既能利用Markov鏈將降水日和非降水日相聯系，也能在考慮降水年內分布不均勻的前提下利用降水的概率分布較好的模擬降水日的降水量［10］。例如王斌等應用模擬的降水序列，結合缺水率指標劃分農業旱情等級［11］；秦道清等將比例調整算法引入隨機模擬模型，生成日降水模擬序列［12］。但這些研究并沒有考慮到不同降水量級下降水年內分布規律和模擬參數的差異性，因此根據年降水量級的不同分別構建隨機模擬模型能有效提高模擬的代表性和準確性。

本文的創新點在于將降水隨機模擬和SWAT 模型相結合，基于不同降水量級下降水年內分布規律和模擬參數的差異性，通過構建不同量級降水下的隨機模擬模型，生成未來年降水量為800 mm（豐水）、550 mm（平水）、300 mm（枯水）相對應的日降水過程，并將其輸入SWAT模型預測徑流，從而進一步探究降水和徑流的關系，具體流程見圖1。

1 數據與方法

1.1 研究區域概況及數據來源

沙河屬于海河流域大清河水系，發源于山西省靈丘縣，流經河北省阜平境內，最后注入王快水庫。王快水庫以上流域控制面積為3 770 km2，主河道長166 km，占沙河流域6 420 km2的59%。該流域多為山區，地表植被覆蓋較差，水土流失較為嚴重。王快水庫以上流域地處中溫帶半濕潤半干旱大陸季風氣候區，多年平均氣溫12.5 ℃，最高氣溫41.7 ℃，多發生在7月份，最低氣溫-18.9 ℃，多發生在1月份［13］。

本文氣象數據選用CMADS 數據集［14，15］在研究區域內及周邊的15 個站點和阜平站的氣象資料（最高、最低氣溫，平均風速，日降水量，相對濕度，日照時數），徑流資料選用阜平水文站2008-2015年實測徑流資料，研究區域地理位置和站點分布見圖2。

圖2 研究區域和站點分布Fig.2 Study area and distribution of stations

1.2 SWAT模型

SWAT 是一種基于物理機制的半分布式水文模型，該模型以日為時間步長，能夠在不同的時間尺度（日、月和年）上進行連續的長時段的模擬，通常設計用于利用流域和流域尺度的輸入數據對連續的長期徑流、泥沙和農藥產量進行建模和模擬［16］。SWAT模型將流域水文過程模擬分為水循環陸面過程和水循環水面過程，其中前者綜合了氣候、冠層截留、蒸騰蒸發、入滲、植被、地下水等產流和坡面匯流的影響因素，后者綜合了子河道傳輸過程、主河道匯流、洼地滯蓄過程等河網匯流的影響因素。因此，相比集總式模型，SWAT模型對于流域下墊面非均質的描述更符合實際情況［17］，SWAT 模型的具體運行流程見圖3。

圖3 SWAT模型運行流程Fig.3 Operating process of SWAT model

1.3 耦合Markov鏈和Gamma分布的降水模擬

一階Markov 鏈模型將自然日分為降水日（P≥0.1 mm）和非降水日（P＜0.1 mm），模型的兩個基本參數是非降水日轉移到降水日的概率P（W/D）和降水日轉移到降水日的概率P（W/W）。研究表明當月的兩個狀態轉移概率值均與同月降水日出現頻率q密切相關且存在線性關系，且轉移概率可用下式估算［18］：

若連續型隨機變量X服從Gamma 分布，則X的密度函數為：

式中：α為形狀參數，決定了Gamma 分布的形狀；β為尺度參數，決定了Gamma 分布的分散度。一階Markov 鏈模型的兩個轉移概率值和Gamma分布的兩個參數合稱為降水模擬參數。

隨機模擬步驟如下：

（1）以1月為例，假設1月1日為非降水日，在確定模擬參數的基礎上，將計算機產生的［0，1］之間的隨機數和P（W/D），P（W/W）進行比較。

（2）當前一日為非降水日時，若隨機數大于當月的P（W/D），則當日仍為非降水日，反之則為降水日；前一日為降水日時，若隨機數小于當月的P（W/W），則當日仍為降水日，反之則為非降水日；

（3）當日為降水日時，將當月降水的Gamma 分布結合隨機數產生當日降水量。以此類推，即可產生若干年的逐日降水序列［19］。

2 SWAT模型的構建

2.1 基礎數據庫的構建

除氣象數據和徑流數據以外，SWAT 模型的構建還需要流域的地形高程數據、土地利用數據、土壤分類數據等資料。其中地形數字高程圖（DEM）來源于地理空間數據云，分辨率為30 m；土地利用數據選取2015年全國土地利用數據，分辨率為30 m；土壤數據來自HWSD 全球土壤數據庫1∶100 萬中國土壤數據集。此外，由于水庫位于流域出口位置且研究區域多為山區，人類活動對研究區域內水文循環過程的影響較弱，因此在SWAT模型中可忽略不計。

根據DEM 圖提取坡向、坡度數據，生成王快水庫以上流域河網，并根據設定的臨界集水面積將王快水庫以上流域劃分為14 個子流域（圖4），其中阜平水文站位于5 號子流域的出口；對流域內的土地利用/覆蓋數據進行重分類，共分為耕地、林地、草地、水域、城鄉工礦居民用地、未利用地等6個種類；以HWSD土壤數據集為基礎，結合SPAW軟件計算土壤層有效持水量、土壤濕密度等參數，完成土壤數據庫的制作并對土壤進行重分類，共分為26 個種類，其中土地利用分類和土壤分類見圖5。根據土地利用、土壤、坡向設置，共劃分水文響應單元（HRU）119個。

圖4 子流域劃分Fig.4 Watershed Subdivision of study area

圖5 土地利用分類和土壤分類Fig.5 The classification of land use and soil

2.2 參數敏感性分析

在模型初步構建完成后，應當在其運行結果的基礎上進行參數敏感性分析，從而確定該模型中對徑流模擬影響較大的參數，以此為依據求解使模型模擬效果最佳的參數組合。本文依托SWAT-CUP 軟件的SUFI2算法進行參數敏感性分析，參數的相對顯著性可由t-Stat檢驗值確定，其中t-Stat檢驗值的絕對值越大，參數敏感性越高［20］。選取影響徑流的28 個參數［21］，將模型輸出結果和阜平水文站的實測徑流數據結合，進行1 000 次迭代運算，挑選敏感性排名前17 的參數進行下一步的率定，敏感性分析見表1。

表1 參數敏感性分析Tab.1 Parameter sensitivity analysis

由表1 可得：各參數對徑流模擬結果均有不同程度的相關性，其中最為敏感的參數是濕潤條件下的初始SCS 徑流曲線數，該參數的大小代表著下墊面不透水性的強弱，即CN2 值越大，徑流量越大；第二敏感的參數為地下水延遲時間，反映的是基流快慢，對水文過程有重要影響；此外，土壤飽和容重、最大冠層截留量等參數對模型輸出結果的影響也較為顯著。

2.3 參數率定與驗證

將2008，2009 作為預熱年份，選用阜平水文站2010.1.1-2012.12.31 的實測流量對模型進行參數率定，本文選取決定系數（R2）、納什效率系數（NSE）、相對誤差（Re）3 個指標用于評價模型的模擬精度，R2和NSE越接近1，Re越接近0，模擬效果越好。一般來說，當模擬結果滿足R2＞0.6，NSE＞0.5，|Re|≤20%時視為模擬結果可信；當模擬結果滿足R2＞0.7，NSE＞0.65 時視為模擬結果良好［22］。計算公式分別如下：

式中：Qk為第k個實測流量；Sk為第k個模擬流量為實測流量的均值為模擬流量的均值。

在參數敏感性分析的基礎上，進行多次迭代，使模擬結果與實測結果較為擬合，此時月徑流R2為0.92，NSE為0.90，Re為-9.1%。選取2013.1.1-2015.12.31 的實測流量數據作為驗證資料，將參數組代入SWAT模型進行模擬驗證，模擬結果與實測徑流擬合較好，此時月徑流R2為0.85，NSE為0.81，Re為8.6%，精度達到模型的良好標準，率定期和驗證期的模擬結果如圖6。

圖6 阜平站實測月徑流與模擬徑流對比Fig.6 Comparison of measured monthly runoff and simulated runoff at Fuping Station

3 基于降水隨機模擬的徑流預測

3.1年際-年內降水統計分析

由圖7 可得：1958-2017年間，阜平站年降水量在均值線附近上下波動，21世紀以來，年降水量的波動幅度逐漸趨于平緩。汛期降水量占全年降水量比重大多為80%～90%，且年降水量過大或過小時，汛期降水占比相比均值存在較大的偏差，即年降水量的量級和降水年內分布的規律存在一定的相關性。因此針對不同量級年降水量進行年內分布模擬更準確，根據多年降水量均值和均方差將年降水量劃分為3 個狀態，劃分標準見表2。

圖7 1958-2017年阜平站年降水量Fig.7 Annual precipitation at Fuping Station from 1958 to 2017

為合理預測未來不同量級年降水量下的徑流量，設置300 mm 為某一枯水年的年降水量，550 mm 為某一平水年的年降水量，800 mm 為某一豐水年的年降水量。根據狀態I內22年的逐日降水量進行300 mm年降水量的年內分布模擬，同理根據狀態II、狀態III 內逐日降水量分別進行550、800 mm年降水量的年內分布模擬。

表2年降水量狀態劃分標準Tab.2 Standard for dividing annual precipitation status

3.2 降水隨機模擬

根據阜平站1958-2017年的逐日降水數據，分別計算各狀態下各月降水日的出現頻率，結合式（1）、（2）估算降水狀態轉移概率，其中各月P（W/D）見圖8。

圖8 各月降水狀態轉移概率P（W/D）Fig.8 Probability of monthly precipitation state transition P（W/D）

從圖8 中可以看出P（W/D）總體上均呈現先增大后減小的趨勢，且在汛期（5-9月）較大，在非汛期（10-4月）較小。當年降水量為狀態III 時，P（W/D）最大值出現在8月；當年降水量為狀態I 或II 時，P（W/D）最大值出現在7月。以各狀態1月、7月、12月降水量為例，其Gamma 分布概率Q-Q 圖見圖9，同理可得不同狀態年降水量的各月降水基本符合Gamma分布，滿足構建耦合Markov 鏈-Gamma 分布隨機模擬模型的要求，并計算相應的分布參數見圖10。

圖9 各狀態1、7、12月降水量Gamma分布Q-Q圖Fig.9 Q-Q chart of Gamma distribution of precipitation in January，July and December under different states

圖10 不同狀態年降水量各月Gamma分布參數Fig.10 Gamma distribution parameters of annual precipitation in different states

以800 mm年降水量的年內分布模擬為例，耦合狀態III 中各月的降水轉移概率和Gamma 分布參數，根據生成的隨機數，連續進行100年的模擬，綜合考慮年降水總量接近800 mm、汛期降水量占比和各月降水天數與樣本均值相近等條件選取模擬年，對其進行適當修正后即可得800 mm年降水量在日尺度上的分配。同理，求得550 mm 和300 mm年降水量的年內分配見圖11。

圖11 不同量級年降水量年內分配Fig.11 Annual distribution of precipitation in different magnitudes

3.3 不同量級降水下的徑流預測

將模擬的日降水輸入到SWAT 模型中，結合已驗證的參數進行輸出，即可得出3種降水情景下的月徑流過程線見圖12。

由圖12 可得，年降水量為800、550、300 mm 的隨機降水模擬結果輸出的阜平水文站月徑流過程中，三者均呈現減小-增大-減小的變化趨勢，且1 至3月的徑流量基本相等，這是由于汛期降水量占了年降水量較大的比重，不同量級的年降水量在1 至3月份的降水分配相差不大；年降水量為800 mm 相對應的月徑流過程中，月徑流量在8月份達到最大，為17.81 m3/s，全年平均徑流量為5.97 m3/s；年降水量為550 mm 相對應的月徑流過程中，月徑流量在9月份達到最大，為7.01 m3/s，全年平均徑流量為2.79 m3/s；年降水量為300 mm 相對應的月徑流過程中，月徑流量在9月份達到最大，為4.08 m3/s，全年平均徑流量為1.51 m3/s。

圖12 不同量級年降水量下的徑流模擬Fig.12 Runoff simulation under annual precipitation of different magnitudes

4 結 論

基于2015年土地利用數據構建相應的SWAT 模型，可得CN2、GW_DELAY、SOL_BD等參數對徑流模擬結果的影響較為顯著；根據流域內阜平水文站2008-2015年的實測徑流數據對模型進行率定和驗證，率定期和驗證期的R2和NSE均高于0.8，|Re|均小于10%，表明該SWAT 模型對于研究區域徑流的模擬結果良好。

根據阜平站1958-2017年逐日降水數據可得，不同降水量級對應的年內分布規律存在一定差異。在此基礎上，分別構建了耦合Markov-Gamma 分布的降水隨機模擬模型，并生成了年降水量為800、550、300 mm 對應的日降水過程，輸入模型后得三種降水情景下的年徑流量分別為5.97、2.79、1.51 m3/s，最大月徑流量分別為17.81、7.01、4.08 m3/s。這一結論體現了不同量級降水對徑流過程的影響，可為后續降雨-徑流關系的研究提供理論基礎。此外，耦合Markov-Gamma 分布的降水隨機模擬模型沒有充分考慮特大暴雨等極端降水事件，從而在一定程度上影響徑流的模擬，因此對極端降水的發生過程進行研究并和降水模擬相結合是一個研究方向。