不同目標函數對水文模型模擬效果影響研究

張如強，高紅凱，劉俊國,，田 勇

(1.北京林業大學自然保護區學院，北京 100083; 2.亞利桑那州立大學，亞利桑那州坦佩 AZ85282； 3.南方科技大學環境科學與工程學院，廣東 深圳 518055)

水文模型作為模擬水文過程和認識水文規律的有效工具，已廣泛應用于科學研究領域和解決實際問題，如洪水分析與預報、水資源管理等[1,2]。水文模型是對現實世界的簡化和抽象，模型中包含難以直接測量的參數，所以需要對模型參數進行率定，間接得到這些參數值[3,4]。另外，無論選用哪種水文模型，都需要對模型的模擬效果進行評價[5]。模型的率定和模擬效果的評價通常是對比模擬流量和觀測流量，使目標函數值最大或最小[6,7]。

模型的模擬效果與目標函數的選取直接相關[8]。納什系數(E)是模型率定過程中應用最為廣泛的目標函數，也是評價模擬效果的重要數值指標[7,9,10]。但是，納什系數存在不足之處：計算公式中觀測流量與模擬流量之差采用平方的形式，這導致洪峰的影響較基流大，即納什系數對洪峰過于敏感[7]。因此產生了納什系數的變形形式，主要包括以自然數為底的對數形式(Eln)、開方形式(Esqrt)、絕對值形式(Eabs)和相對值形式(Erel)[6,7,11,12]。盡管這些變形形式已經被學者應用于模型率定和評價[5,6,13]，但這些目標函數的選取對模擬效果影響的研究相對較少，選擇哪種變形形式更為合理這一問題尚未得到圓滿解答。

本研究擬通過將納什系數、各個變形形式設定為目標函數，以納什系數、相關系數、平均絕對誤差和平均值之差4個數值指標綜合分析模型模擬效果，探索不同目標函數的選取對模擬效果的影響。

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

黑河發源于祁連山，流經甘肅、青海、內蒙古自治區，最后流入居延海，干流全長821 km，是我國西北地區第二大內陸河。以鶯落峽和正義峽為界，黑河流域分為上、中、下游。本文的研究區為黑河上游，位于祁連山西北，面積約1.000 9 萬km2，海拔高度介于1 700～4 900 m；研究區覆被包括森林(20%)、草地(52%)、裸地(19%)、濕地或洪泛區(8%)和永久性冰川雪地(0.8%)。研究區降雨集中，超過80%大降雨分布在5-9月；多年平均降雨和潛在蒸散發為430和520 mm/a。

數據來源于寒區旱區科學數據中心。研究使用數據為日降雨、日平均氣溫和逐日的徑流量。其中，氣象數據為托勒站、祁連站、野牛溝站、張掖站四個氣象站1959-1979年和1990-2006年日降雨和日平均溫度，徑流數據為鶯落峽水文站相應時間觀測流量。研究區位置、覆被以及氣象站水文站位置見圖1。

圖1 研究區位置、黑河上游覆被及氣象、水文站點分布Fig.1 Location and Land cover of Upper Heihe River basin and distribution of the meteorological andhydrologic stations in and around upper Heihe

2 研究方法

2.1 模型簡介

地形與地質條件、土壤類型、氣候、地表覆被等密切相關，而這些正是影響水文過程的主要因素。地形綜合反映了流域下墊面信息，在不同尺度上影響了降雨-徑流過程[14]?；谶@一假設Gao等[15]開發了半分布式水文模型FLEX-Topo，并應用于黑河上游，取得了滿意的模擬效果[15]。因此，本文研究中采用FLEX-Topo模型[15]。依據地形因素將流域分為4種水文景觀類型。不同水文景觀類型的產匯流特點不同，因此，模型中采用了不同的產匯流機制，使模型具有了半分布式的特點。依據降雨形成徑流的過程，模型包含了5個水箱：積雪水箱(Sw)、植被冠層截留水箱(Si)、非飽和土壤層水箱(Su)、快速退水水箱(Sf)和慢速退水水箱(Ss)。整個模型通過5個水箱模擬降雨經過截留(I)、蒸散發(Ei&Ea)、優先流(RfB)、壤中流(Qf)、地下徑流(Qs)等最終流入河網的過程。模型結構見圖2。FLEX-Topo模型包含24個參數，參數含義及取值范圍見表1。

圖2 FLEX-Topo 模型結構[15]Fig.2 The structure of the FLEX-Topo model

參數名稱取值范圍含義FDD[1，8] 度日模型參數，表征當溫度高于Tt時，每升高1°每日融化形成的水量Tt[-2．5，2，5] 當溫度高于該參數時，產生截留Pt[5，35] 當有效降雨大于該參數時，發生霍頓超滲產流SumaxB[5，500] 裸地景觀土壤最大持水能力PmaxB[0．1，10] 裸地景觀的最大深層滲漏能力DB[0，1] 裸地景觀中將徑流分成流入快速退水水箱和慢速退水水箱的參數TlagB[0，5] 裸地景觀中徑流進入快速退水水箱的滯后時間KfB[2，50] 裸地景觀中存儲于快速退水水箱中水流入河網的時間ImaxFH[1，10] 林坡景觀最大截留量SumaxFH[100，1000] 林坡景觀土壤最大持水量βFH[0．1，5] 表征景觀的空間異質性，用于計算林坡的徑流系數D[0，1] 林坡和草坡景觀中將徑流分成流入快速退水水箱和慢速退水水箱的參數TlagFH[0，5] 林坡景觀徑流進入快速退水水箱的滯后時間KfH[1，20] 林坡和草坡景觀中快速退水水箱中的水流入河網的時間

續表1 FLEX-Topo模型參數及參數取值范圍

模型的輸入數據為DEM、日降雨、日平均溫度，輸出為逐日的蒸散發(Ei&Ea)、徑流量(Qf&Qs)、滲流量(RfB&f)。DEM主要用于生成HAND(Height Above the Nearest Drainage)[15-18]和坡度。黑河上游依據這兩個因素分為4種水文景觀類型，即裸地(bare soil/rock)、林坡(forest hillslope)、草坡(grass hillslope)、濕地或洪泛區(wetland/flood plain)，其中，裸地的產流機制主要為深層滲漏補給地下水(DP)和超滲產流(HOF)；林坡和草坡的產流機制為壤中快速徑流(SSF)，主要是優先流；濕地或洪泛區的產流機制為蓄滿坡面產流(SOF)。降雨和溫度是模型的驅動數據，模型采用Hamon方程[19]計算潛在蒸散發。黑河流域上游面積大、海拔高，而研究區周邊僅有4個氣象站，所以模型采用泰森多邊形法計算平均降雨量，并同時考慮海拔梯度對溫度和降雨的影響。徑流量是模型率定和驗證的主要數據。

2.2 模型率定算法

FLEX-Topo模型采用的參數率定算法為MOSCEM-UA(Multi-Objective Shuffled Complex Evolution Metropolis-University of Arizona)[20]。Vrugt等[20]開發了該算法(算法流程圖見圖3)，并將其用于參數數量、復雜程度不同的3個模型中，證明該算法的有效性和高效性。該算法已被國內外學者廣泛使用[21-25]。MOSCEM-UA算法需要設置3個參數：迭代次數、復合型個數以及隨機樣本數。根據Vrugt等[20]的建議本研究設定的3個參數值分別為100 000、24、2 304。

圖3 MOSCEM-UA算法流程圖[20]Fig.3 Flowchart of the MOSCEM-UA algorithm

2.3 研究方法

納什系數是模型率定過程中應用最為廣泛的目標函數，鑒于納什系數的不足之處，產生了變形形式，主要包括以自然數為底的對數形式(Eln)、開平方形式(Esqrt)、絕對值形式(Eabs)和相對值形式(Erel)4種。4種變形的計算公式見表2。

FLEX-Topo模型及MOSCEM-UA算法的可靠性已經得到證明。本研究采用該模型和算法。首先，使用相同的參數范圍、迭代次數(100 000)、復合型個數(24)以及隨機樣本數(2 304)，設定不同的目標函數(納什系數及變形形式)進行參數率定(流程見圖3)，獲得5套參數。然后，將每一套參數代入模型進行驗證，獲得5套模擬結果。最后對5套模擬結果進行分析、比較。本研究中模型參數率定期設為1959-1978年，驗證期設為1990-2006年，兩階段分別取1959和1990年為模型預熱期。

表2 納什系數及變形形式計算公式[20]Tab.2 Mathematical formulation of Nash and Sutcliffe efficiency and its modified forms

衡量模擬效果通常有兩種途徑，一是比較直觀的圖形，二是采用數值指標[5,6]。常用的圖形包括水文過程圖和流量頻率分布曲線，從圖形中可以直觀地查看模擬結果。其中，水文過程圖強調流量隨時間的變化，而流量頻率分布曲線側重于對流量的高估與低估問題；按照有無量綱，數值指標可以分為相對指標和絕對指標兩大類，其中，相對指標中應用最為廣泛的是納什系數，絕對指標包括均方根誤差、平均絕對誤差等。Willmott[26]指出判斷模擬效果時納什系數這一相對指標不宜單獨使用。Legates等[7]通過研究表明，平均絕對誤差略優于均方根誤差，因此，根據Legates等[7]和Anderton等[27]的建議，在評價模型效果時選取納什系數(E)、相關系數(r)、平均絕對誤差(MAE)和平均值之差(DM)4個數值指標來衡量模擬效果。4個數值指標的計算公式見表3。

3 結果與討論

3.1 模擬結果

根據模擬結果繪制水文過程圖和流量頻率分布曲線。不同目標函數驗證期的模型水文過程圖和流量頻率分布曲線見圖4和圖5。圖4和圖5中，Qo為觀測流量，Qx為以x為目標函數的模擬流量。

由圖4可以看出，以納什系數及變形形式為目標函數的模擬結果中，模擬流量基本可以反映觀測流量隨時間的變化，但吻合程度不同：較好的是E和Esqrt，Eln和Eabs次之，Erel最差。流量頻率分布曲線反映模型對漲水段、洪峰、退水段和基流等的模擬效果。由圖5可以看出，對于洪峰流量，E出現高估，而Eln、Esqrt和Eabs較E略好，Erel則低估了洪峰流量；對于基流，各個目標函數的模擬效果相似，均低估了基流流量；漲水段和退水段除Erel最差外，E、Eln、Esqrt和Eabs表現差別不大，均小于觀測值。

表3 3個數值指標的計算公式Tab.3 Mathematical formulation of the three numerical criterion

圖4 驗證期鶯落峽水文過程圖Fig.4 Hydrograph of validation period of Yingluoxia

圖5 驗證期鶯落峽流量頻率分布曲線Fig.5 Flow duration curve (FDC) of validation period of Yingluoxia

雖然水文過程圖和流量頻率分布曲線可以直觀地看出模擬結果，但整體模擬效果難以客觀地直接衡量。因此，本文計算了5套模擬結果的4個數值指標，計算結果見表4。

從表4中可以看出：①4種變形形式的模擬結果中，相關系數一列的值相差不大；其余三列除Erel最差外，Eln、Esqrt和Eabs相差很小，Esqrt在各個數值指標中均為最優；②E與Esqrt相比，相關系數有所提高，這表明Esqrt較E更能反映流量的動態變化，而其余指標差別很??；③5個模擬結果中，Erel的各數值指標值均為最差；④5個模擬結果中，平均值之差一列均為正數，表明模擬結果均低估了總體水量。

表4 不同目標函數的數值指標Tab.4 Numerical criterionof different objective functions

3.2 結果分析與討論

Legates等[7]指出，納什系數的計算公式采用平方形式，導致洪峰流量的影響過大，即納什系數對洪峰流量過于敏感、忽視基流，這是欠妥的。為了解決這一問題，產生了對數形式、開方形式、絕對值形式、相對值形式。

絕對值形式的納什系數改變了計算公式中采用的平方形式(見表2)。Legates等[7]建議使用絕對值形式的納什系數，因為平方形式將模擬流量與觀測流量的差值放大，而絕對值形式將賦予兩者差值相同的權重，因此減小了洪峰流量的影響。對數形式、開方形式和相對值形式則采用另一種方式：改變流量的絕對大小(見表2)。因為計算公式中采用這種處理，就降低了洪峰流量的影響。Pushpalatha等[5]、Krause等[6]、Houghton-Carr[11]、Oudin等[12]、Kim等[13]研究證明這些變形降低了納什系數對洪峰流量的敏感性。董磊華等[28]比較了對數均方誤和平方根均方誤為目標函數的模擬結果，文章指出：由對數函數曲線(y=lnx)和開方函數曲線(y=x1/2)可知，x值越小，y值遞增的梯度越大，同理，流量值越小模擬流量與實測流量的誤差放大越大，因此對數形式和開方形式削弱了洪峰流量的影響、增強基流流量作用。Krause等[6]提出了基于相對偏差的納什系數計算方式，即相對值形式。通過這種變形，觀測流量和模擬流量的差值用相對偏差來衡量，從而減小洪峰流量的影響、增強基流流量的影響。但是，從結果來看，這種變形形式的模擬效果在所有變形形式中最差。

上述研究表明變形形式彌補了納什系數對洪峰流量過度敏感這一不足。然而，本研究發現，變形形式對基流流量的模擬并不理想。造成這一結果的原因可能有兩個方面：一是觀測誤差的存在，包括降雨、溫度、徑流的觀測誤差；二是納什系數計算公式中參考值為平均流量。對于季節流量波動較小的流域，選擇平均值為參考值是合適的，但是對于流量波動大的流域則有待商榷。

4 結 語

本研究通過設定不同的目標函數，從納什系數、相關系數、平均絕對誤差和平均值之差4個指標評估模擬效果。對比各個目標函數驗證期的水文過程圖、流量頻率分布曲線和數值指標，可以看出以相對值形式為目標函數的模擬結果最差，而以納什系數和開方形式為目標函數的模擬結果最好；綜合比較以納什系數和開方形式為目標函數的模擬結果，兩者的數值指標差別很小，但前者會明顯高估洪峰流量，后者的相關系數更高，更能反映流量的動態變化。整體來講，當選擇納什系數及其變形形式為目標函數時，開方形式的模擬效果最好。

盡管如此，5個目標函數的模擬結果均低估了基流流量，同時，也低估了整體水量。造成這一結果的原因可能有兩個方面：一是觀測誤差；二是計算公式中參考值選擇為平均流量。鑒于觀測技術的限制以及觀測誤差存在的必然性，建議對納什系數公式中參考值的選擇開展進一步研究以提高模擬效果。

□

[1] 徐宗學,包為民,謝正輝,等. 水文模型[M].北京: 科學出版社, 2009.

[2] 劉金濤, 馮 杰, 張佳寶. 分布式水文模型在流域水資源開發利用中的應用研究進展[J]. 中國農村水利水電，2007,(2):142-144.

[3] Muleta M K. Model performance sensitivity to objective function during automated calibrations [J].Journal of Hydrologic Engineering, 2012,17(6):756-767.

[4] Vrugt J A, Gupta H V, Bouten W, et al. A shuffled complex evolution metropolis algorithm for optimization and uncertainty assessment of hydrologic model parameters [J]. Water Resources Research, 2003,39(8):1 201.

[5] Pushpalatha R, Perrin C, Moine N L, et al. A review of efficiency criteria suitable for evaluating low flow simulations [J]. Journal of Hydrology, 2012,420-421:171-182.

[6] Krause P, Boyle D P, B?se F. Comparison of different efficiency criteria for hydrological model assessment [J]. Advances in Geosciences, 2005,5(5):89-97.

[7] Legates D R, Gregory J, Mccabe JR. Evaluating the use of ‘goodness-of-fit’ measures in hydrologic and hydrologic model validation [J].Water Resources Research, 1999,35(1):233-241.

[8] Gan T Y, Dlamini E M, Biftu G F. Effects of model complexity and structure, data quality, and objective functions on hydrologic modelling [J]. Journal of hydrology, 1997,192:81-103.

[9] W?hling T, Samaniego L, Kumar R. Evaluating multiple performance criteria to calibrate the distributed hydrological model of the upper Neckar catchment [J]. Environmental Earth Sciences, 2013,69(2):453-468.

[10] Schaefli B, Gupta H V. Do Nash value have value? [J]. Hydrological Processes, 2007,21:2 075-2 080.

[11] Houghton-carr H A. Assessment criteria for simple conceptual daily rainfall-runoff model [J]. Hydrological Sciences, 1999,44(2):237-261.

[12] Oudin L, Aneréassian V, Mathevet T, et al. Dynamic averaging of rainfall-runoff model simulations from complementary model parameterizations [J]. Water Resources Research, 2006,42(7):887-896.

[13] Kim H S, Lee S. Assessment of a seasonal calibration technique using multiple objectives in rainfall-runoff analysis [J].Hydrological Processes, 2014,28:2 159-2 173.

[14] Savenije H H G. HESS Opinions “Topography driven conceptual modelling (FLEX-Topo)”[J]. Hydrology and Earth System Sciences, 2010,14:2 681-2 692.

[15] Gao H, Hrachowitz M, Fenicia F, et al. Testing the realism of a topography-driven model (FLEX-Topo) in the nested catchments of the Upper Heihe, China[J].Hydrology and Earth System Sciences, 2014,18:1 895-1 915.

[16] Renn? C D, Nobre A D, Cuartas L A, et al. HAND, a new terrain descriptor using SRTM-DEM: Mapping terra-firme rainforest environment in Amazonia [J]. Remote Sensing of Environment, 2008,112:3 469-3 481.

[17] Nobre A D, Cuartas L A, Hodnett M, et al. Height Above the Nearest Drainage-a hydrologically relevant new terrain model [J]. Journal of Hydrology, 2011,404:13-29.

[18] Gharari S, Hrachowitz M, Fenicia F, et al. Hydrological landscape classification: investigating the performance of HAND based landscape classifications in a central European meso-scale catchment [J]. Hydrology and Earth System Sciences, 2011,15:3 275-3 291.

[19] Hamon W R. Estimating potential evapotranspiration [J]. Journal of the Hydraulics Division-Asce,1961,87:107-120.

[20] Vrugt J A, Gupta H V, Bastidas L A, et al. Effective and efficient algorithm for multiobjective optimization of hydrologic model [J]. Water Resources Research, 2003,39(8):1-19.

[21] Schoups G, Addams C L, Gorelick S M. Multi-objective calibration of a surface water-groundwater flow model in an irrigated agricultural region: Yaqui Valley, Sonora, Mexico [J]. Hydrology and Earth System Sciences, 2005,9:549-568.

[22] Tang Y, Reed P, Kollat J. Parallelization strategies for rapid and robust evolutionary multiobjective optimization in water resources applications [J]. Advances in Water Resources,2007,30(3):335-353.

[23] Parajka J, Merz R, Bl?schl G. Uncertainty and multiple objective calibration in regional water balance modelling: case study in 320 Austrian catchments [J]. Hydrological Processes, 2007,21(4):435-446.

[24] Fenicia F, Solomatine D P, Savenije H H G, et al. Soft combination of local models in a multi-objective framework [J]. Hydrology and Earth System Sciences, 2007,11:1 797-1 809.

[25] 王宇暉,雷曉輝,蔣云鐘.基于MOSCEM-UA算法的水文模型參數優化[J].人民長江，2010,41(19):95:99.

[26] Willmott C J. On the validation of models [J]. Physical Geography, 1981,2(2):184-194.

[27] Anderton S, Latron J, Gallart F. Sensitivity analysis and multi-response, multi-criteria evaluation of a physically based distributed model [J]. Hydrological Processes, 2002,16:333-353.

[28] 董磊華, 熊立華. 水文模型中不同目標函數的影響比較分析[J].水文，2009, 29(3):24-27.