面雨量計算不確定性對新安江模型徑流模擬的影響研究

李 莎，張 濤，張行南,3,4，方園皓，程中陽

(1.河海大學水文水資源學院，南京 210098；2.長江水利委員會水文局，武漢 430010；3.河海大學水安全與水科學協同創新中心，南京 210098； 4.河海大學水資源高效利用與工程安全國家工程研究中心，南京 210098；5.水利部珠江水利委員會技術咨詢中心，廣州 510611)

0 引 言

流域水文模型是以產匯流過程為核心對水循環現象的描述，在防洪減災、水資源管理、水環境和水生態保護等方面發揮了重要作用。降雨是水文循環要素中最活躍的因子，與徑流的產生與匯聚有關，在水文模擬中發揮著控制水量平衡作用，準確描述降雨的空間分布是模型模擬成功的關鍵[1]。面雨量能客觀反映一個流域的降水情況，在分析、預報水情變化時面雨量應用非常廣泛，在水文模擬中是一個重要參數[2]。要精確獲取流域降雨空間分布特征需要建立密度極高的雨量站網，目前站網的資料水平遠達不到要求。因此需要從插值方法、站網密度、站點分布均勻度、降雨空間分布等方面研究面雨量計算不確定性及其對水文模擬的影響。石朋[3]從理論和應用方面比較、分析了距離平方倒數法、普通克里金方法、引入高程信息的協克里金方法插值面雨量的優缺點。井立陽等[4]以西江小流域為例，分析比較了自然流域分塊、相關分析、等雨量線3種降雨空間插值方法對新安江模型模擬精度的影響。張雪松等[5]從雨量站密度、雨量站點分布、降雨空間分布3方面研究了其對分布式水文模型SWAT模型徑流模擬的影響。薛豐等[6]研究了不同雨量站點密度下4種常用的降雨空間插值方法對SWAT模型徑流模擬結果的影響，指出在研究流域現有站網密度條件下距離平方倒數法和泰森多邊形法的徑流模擬結果較好，但在雨量站數量較少的流域應使用最鄰近法進行面雨量的空間插值。流域概念性水文模型被廣泛應用于水文預報、防洪減災、水資源開發利用等方面。新安江模型是我國自主研制的典型概念性水文模型，主要應用于濕潤和半濕潤地區，模型的模擬精度取決于輸入變量和參數率定[7]。肖楠等[8]在以相對誤差百分值為標準同倍比改變降雨輸入數據的基礎上，基于新安江模型研究山丘區小流域洪峰模擬值的誤差隨降雨輸入偏差的變化特點。國內學者對降雨空間插值方法與徑流模擬精度關系的研究大都在現有最高站網密度下進行，基于不同雨量站網密度條件研究降雨空間插值方法與徑流模擬精度關系的較少。

因此本文采用降雨空間插值方法中較為常用的泰森多邊形法和反距離權重法進行面雨量的空間插值計算，對比不同雨量站網密度條件下不同插值方法的面雨量計算結果，在統一模型參數前提下分析多種降雨輸入情境下新安江模型的模擬精度，研究站網密度對徑流模擬的影響以及不同站網密度下降雨空間插值方法與徑流模擬精度的關系，為涪江流域建立降雨空間描述方案提供參考。最后本文對比了面雨量計算相對偏差與徑流量模擬相對偏差，研究2者之間的關系。

1 研究區概況

天仙寺流域起源于涪江最大的一條支流梓潼江。梓潼江又名梓江，發源于龍門山東南麓江油縣藏王寨棋盤山鷹咀崖一帶，河長達340 km，多年平均徑流總量約4.7 億m3，僅為通口河徑流總量的1/7。天仙寺流域地理位置在30°59′N～32°12′N、104°48′E～105°40′E，位于涪江流域東部，流域面積4 976 km2，流域地形起伏較小，以小起伏地形和丘陵為主，流域東側基本為小起伏地形，西側基本為丘陵，北部和中部小部分地區被平原覆蓋，流域平均高程1 316 m。流域屬于亞熱帶濕潤性季風氣候區，雨量豐沛但時空差異大。天仙寺流域地理位置及雨量站點分布情況見圖1，流域內現有12個雨量站，站網密度為414.67 km2/站，雨量站點整體上分布較均勻但計算單元1內有3個雨量站的位置靠得很近，因此在隨機抽站時設定不能同時抽取這3個雨量站中的任意2個雨量站。

圖1 天仙寺流域概況Fig.1 The map of Tianxiansi basin

2 研究方法

2.1 面雨量計算方法

2.1.1 泰森多邊形法

泰森多邊形法(Thiessen Polygon Method，THI)屬于比較簡單的面積加權平均法，只使用距離最近的單個的雨量站點來描述泰森多邊形區域的面平均雨量[3]。該方法由荷蘭氣象學家泰森(A.H. Thiessen)提出，旨在用離散點的性質描述區域的性質。具體做法為將雨量站點兩兩相連并作連線的中垂線，各中垂線相交形成若干個多邊形，這個多邊形被稱為泰森多邊形。在泰森多邊形中僅包含一個雨量站，而泰森多邊形區域的面平均雨量就用這個雨量站的實測值表示。該法僅考慮了面積加權沒有考慮到高程對降雨的影響，適用于地形起伏變化不大的流域。

2.1.2 反距離權重法

反距離權重法(Inverse Distance Weight Method，IDW)是地統計分析中的一種局部插值方法，也稱為距離反比加權法，是一種加權移動平均方法，以待估點與范圍內的雨量站點之間的距離倒數為權重，是一種確定性的內插方法，在采樣點分布均勻且未聚類時效果最佳[9]。反距離權重方法通用公式表示為：

(1)

式中：v0為未知點的雨量內插值；vi為第i(i=1,2,…，n)個雨量站點的雨量；di為雨量站點與未知點之間的距離；k為距離的冪，它顯著影響插值結果，通常取2。

文中將研究區域離散成0.05°×0.05°的網格，在此基礎上按上述2種插值方法進行面雨量計算，插值結果包括全流域面平均雨量和計算單元面平均雨量。

2.2 面雨量計算分析指標

一般認為雨量站密度越大其降雨量資料越能真實的反映流域降雨的空間分布情況[10]，因此本研究以現有12個雨量站點插值計算得到的面平均雨量作為面雨量“近似真值”。選用平均絕對偏差和均方偏差度量面雨量計算偏差。平均絕對偏差(Mean Absolute Error，MAE)反映了降雨理論值偏離近似真值的大小，并給出了可能的誤差范圍。均方偏差(Root Mean Square Error，RMSE)反映了理論值與近似真值系列的差異程度。兩者的計算公式如下：

(2)

(3)

式中：pc為雨量理論值；p0為雨量近似真值；n為降雨系列個數。

2.3 徑流模擬結果評價指標

本文采用相對誤差(RE)、均方根誤差(RMSE)、納什效率系數(Nash-Sutcliffe，NSE)作為新安江模型模擬精度的評價指標。其中相對誤差用于反映徑流總量模擬值與徑流總量實測值之間的差異程度；均方根誤差用于反映徑流模擬系列偏離徑流實測系列的平均程度；納什效率系數在水文預報中也稱為確定性系數，用于反映徑流模擬過程與徑流實測過程的吻合度，其值越接近于1，說明模擬過程越精確。各評價指標計算公式如下：

(4)

(5)

(6)

2.3 新安江模型構建及參數率定

新安江模型是分單元計算的，根據流域下墊面的水文、地形、地貌條件將流域劃分為若干個計算單元。在每個計算單元內進行產流計算和坡面、河網匯流計算，得到計算單元出口的流量過程，再將每個計算單元出口的流量過程沿河道采用馬斯京根法演算到流域出口，并在流域出口疊加得到整個流域的徑流過程。新安江模型采用3層蒸散發模式計算流域蒸散發，嚴格控制整個模型計算過程的水量平衡；以張力水蓄水容量曲線反映流域土壤缺水量空間分布的不均勻性，采用蓄滿產流理論計算流域產流量；模型采用自由水線性水庫進行三水源劃分，將徑流成分劃分為地面徑流、壤中流和地下徑流；河網匯流一般采用單位線法或滯時演算法，一般在大流域中常忽略河網匯流這一階段。

新安江模型輸入數據為實測降雨和蒸發皿蒸發數據，研究采用的數據均為日時間尺度數據。本研究將流域劃分成5個計算單元(計算單元劃分情況見圖1)，選取2007-2010年作為模型率定期，2011-2012年作為驗證期，以天仙寺水文站的日徑流系列作為實測資料，構建天仙寺流域徑流模擬模型。

研究采用自動參數率定和手動參數率定相結合的方法進行參數率定，自動率參基于PEST軟件，采用SCE-UA算法進行參數的校準和驗證，在參數自動率定過程中選用的目標函數為納什效率系數(Nash-Sutcliffe，NSE)、相對誤差(RE)、均方根誤差(RMSE)。

3 結果與分析

3.1 面雨量計算分析

考慮到計算單元1內有3個雨量站點分布緊密，設定隨機抽站時不能同時抽取到這3個雨量站中的任意2個雨量站，因此從現有雨量站點中抽取35%、50%、70%、85%的雨量站(化整處理后雨量站個數分別為4、6、8、10)，分別有378、462、117、3個總體。除抽取10個雨量站時考慮總體之外，其余情況均隨機抽取60次(經檢驗抽取60次時全流域和計算單元的面雨量計算MAE、RMSE累積平均值變化趨于穩定(見圖2，站點數量為6、8的計算結果與圖2相似)，此時抽取的樣本可以代表總體)，共183種站網分布情況。對這183種站網分布情境重復利用泰森多邊形法和反距離權重法進行面雨量插值計算。

在計算各插值方法下降雨插值結果的MAE和RMSE時，以現有12個雨量站點在該插值方法下插值得到的面平均雨量作為面雨量的“近似真值”。統計各站網分布情景基于泰森多邊形法和反距離權重法的MAE和RMSE，并計算各站點數目下MAE與RMSE的平均值，作平均值絕對偏差均值、均方根偏差均值與站點數目關系圖，見圖3。

由圖3可知：站點數量變化、插值方法的差異均對降雨數據插值結果有很大的影響。不同的站點數目降雨插值結果不同，計算單元和全流域面雨量插值偏差隨站點數目變化趨勢相同，均隨站點數目增多而減小，即隨著插值站點數目增多，面平均雨量插值結果的“量”和“過程”均更加接近面雨量“近似真值”的“量”和“過程”。這與實際應用是一致的，站點數量越多雨量站的分布越密集，其雨量資料越能真實描述流域的降雨空間分布。由圖3還能看出，在天仙寺流域，各站點數目下反距離權重法對降雨數據插值結果都優于泰森多邊形法，反距離權重法的這種優勢在站點數量較少時較為明顯；隨著站點數目增多，反距離權重法的優勢變得不明顯，其插值出的面雨量結果與泰森多邊形插值的面雨量結果區別不大。這是因為天仙寺流域為中小起伏山地和丘陵地形，參與插值的雨量站數目較少時，泰森多邊形內的地形起伏較大，多邊形內降雨空間分布不均勻性較大，由一個雨量站點的雨量表示多邊形面雨量將會產生較大的偏差。

圖2 抽站次數與統計偏差累計平均值關系(4站)Fig.2 The relationship between number of extracting rain-gauges and cumulative mean of statistical bias(4 gauges)

圖3 站點個數與偏差統計值關系Fig.3 The relationship between number of rainfall stations and statistical bias

3.2 徑流模擬結果分析

將上述368種面雨量計算結果作為模型的輸入進行水文模擬，為消除模型參數的影響，每種方案均采用由反距離權重法插值出的“面雨量近似真值”率定所得的參數進行模擬。

3.2.1 站網密度對徑流模擬的影響

統計各站點數目下年均評價指標的均值(其中年均相對偏差的均值為年均相對偏差絕對值的均值)，并計算各評價指標隨站點數目變化的幅度，變化幅度為評價指標變化范圍與平均值的比值，具體情況見表1。

由表1可知，2種插值方法下模擬結果的各評價指標隨站點個數變化趨勢一致，模型模擬精度隨著站點數目的增加而提高，納什效率系數隨站點數目增多而增大，均方根誤差隨站點數目增多而減小，相對誤差隨站點數目增多而減小但略有波動。流域雨量站站網密度加密其降雨資料更能表現流域內降雨的空間分布，插值出的面雨量也將更接近實際面雨量，徑流模擬結果也會更加精確。表1中10個站點的徑流量誤差低于12個站點的徑流量誤差，這是因為計算單元1的南北地形差異大，區域內地形變化較大，而12個站點時多出的2個雨量站聚集在計算單元1地勢偏低的東南角上，使得12個雨量站時面雨量計算偏差變大。說明雨量站不能聚集布設，尤其在地形起伏大較大的區域要遵從均勻分布的原則布設。另外由表1可知，相對誤差對站點數目的變化最為敏感，均方根誤差和納什效率系數對站點數目變化的敏感性較小，這是因為均方根誤差和納什效率系數受模型參數的影響較大，模型參數的調試將由雨量站網密度變化引起的偏差弱化了許多。

表1 評價指標隨站點個數變化情況Tab.1 Evaluating indicator changing with number of rainfall stations

3.2.2 降雨插值方法對徑流模擬的影響

統計2種插值方法的不同站點個數下的評價指標，并作評價指標與站點個數相關圖，如圖4所示。

圖4 2種插值方法下評價指標隨站點個數變化Fig.4 Map of evaluating indicator changing with the number of rainfall stations under two interpolation methods

從圖4可知:站點數目少于10個即站網密度大于497.6 km2/站時，反距離權重法在NSE、RE、RMSE方面的表現均優于泰森多邊形法，當站點繼續增多，泰森多邊形法的徑流模擬結果略優于反距離權重法。這與泰森多邊形法和反距離權重法在面平均雨量插值精度上表現一致。由此可見，在天仙寺流域雨量站密度較大的情況下，用反距離權重法插值的面平均雨量進行模擬具有一定的優勢；而站點數目較多時，泰森多邊形法略具優勢。

3.2.3 面雨量計算相對偏差與徑流量模擬相對偏差對比

以面雨量“近似真值”輸入得到的年均徑流量作為標準值，將上述183種站網分布情境下的年均徑流量與標準值對比，得到年均徑流量模擬相對偏差。統計站點數目為4、6、8、10時分別基于反距離權重法和泰森多邊形法的年均面雨量計算相對偏差與年均徑流量模擬相對偏差，各方案下年均徑流量模擬相對偏差與年均面雨量計算相對偏差的關系類似，因此這里只展示站點數目為6基于反距離權重法的年均徑流量相對偏差與年均面雨量計算相對偏差關系圖，見圖5(a)。為進一步分析面雨量計算相對偏差與徑流量相對偏差的關系，按面雨量計算相對偏差絕對值遞增作面雨量計算相對偏差絕對值與徑流量相對偏差絕對值的關系圖，見圖5(b)。

圖5 各站網分布下降雨偏差與徑流模擬誤差Fig.5 Rainfall deviation and runoff simulation error under rainfall station network

由圖5(a)可知，60種站網分布情境中大部分站網下的面雨量計算相對偏差低于徑流量模擬相對偏差，小部分站網的面雨量計算相對偏差高于徑流量模擬相對偏差，說明徑流量模擬相對偏差不一定低于面雨量計算相對偏差，這與站網分布有一定的關系。由圖5(a)還能看出不同的站網分布之間面雨量計算相對偏差相差不大，最大相差在10%以內，但徑流量模擬相對偏差相差較大，最大偏差與最小偏差相差約30%。這有可能是因為本研究所采用的模型參數沒有真實反映研究流域的產匯流規律。有的站網分布下面雨量計算相對偏差較大但徑流量模擬相對偏差可以很小，說明合適的站網分布可以很大程度上降低徑流量模擬相對偏差。由圖5(b)可知，隨著面雨量計算相對偏差的遞增，徑流量模擬相對偏差呈震蕩變化，2者的變化趨勢并不相同。

4 結 論

本文基于新安江模型，考慮站網密度的影響，利用泰森多邊形法和反距離權重法進行面雨量的計算，分別開展了站網密度、插值方法對降雨插值精度和徑流模擬精度影響的研究以及面雨量計算相對偏差與徑流量模擬相對偏差關系的研究，得出如下結論。

(1)插值站點的數量不同，面雨量計算結果不同，面雨量計算偏差差異較大，插值雨量站點越多，降雨插值結果的精度越高；空間插值方法對降雨空間插值結果也有較大影響，不同的空間插值方法面雨量計算偏差相差較大，在該研究流域反距離權重法的面雨量計算結果在雨量站網密度較低時優于泰森多邊形法。

(2)在徑流模擬方面，2種插值方法的徑流模擬精度均隨著站點數目的增加而提高。當站點數量較少時反距離權重法在天仙寺流域的模擬精度優于泰森多邊形法的模擬精度，站點數目較多時泰森多邊形法更優。

(3)徑流量模擬相對偏差不一定低于面雨量計算相對偏差，這主要與流域地形、模型參數有關，但站網的合理布設有利于降低徑流量模擬相對偏差。