王 升，閆 妍，黃玉清1,**，胡寶清

(1.廣西百色國家農業(yè)科技園區(qū)管理委員會，廣西百色 533612；2.南寧師范大學,北部灣環(huán)境演變與資源利用教育部重點實驗室，廣西南寧 530001；3.中國科學院華南植物園， 廣東廣州 510650)

參考作物蒸散發(fā)(Reference Evapotranspiration，ET0)是水資源管理、灌溉制度設計以及水文模型的重要參數(shù)，反映了與下墊面(植被類型、地形及土壤水分等)無關的、在給定位置和時間的最大大氣蒸發(fā)能力[1,2]。聯(lián)合國糧食及農業(yè)組織(FAO)推薦采用基于空氣動力學及能量平衡原理建立的Penman-Monteith模型(FAO-56 P-M)作為ET0計算的標準方法，該模型還可用于校正其他ET0計算模型。然而FAO-56 P-M模型需要較完備的氣象數(shù)據(jù)(太陽輻射、氣溫、風速和相對濕度)，以及緯度和海拔高度，即使在發(fā)達國家，能同時測量這些氣象數(shù)據(jù)的氣象站點也很有限，而且FAO-56 P-M模型對數(shù)據(jù)質量有嚴格要求。如在全球范圍內，大多數(shù)氣象站均配備有氣溫傳感器[3]，但其他氣象變量，特別是太陽輻射，僅有少數(shù)氣象站點能夠監(jiān)測[3-5]，所以FAO-56 P-M模型的應用受到一定的限制。另外，F(xiàn)AO-56 P-M模型在計算過程中需要計算很多中間參數(shù)，計算量大且煩瑣[5,6]。因此，使用氣象因子較少、形式簡單且具有足夠精度的ET0計算模型(如Hargreaves、Priestley-Taylor、Ravazzani、Tralkovic、Makkink-Hansen、Oudin和Turc等經驗模型)受到較多關注。學者們對這些經驗模型在各地的適用性進行了檢驗及校正。如Gao等[4]評估了Hargreaves模型(需日平均氣溫、最高氣溫和最低氣溫)和McCloud模型(需平均氣溫)在貴州的適用性，并進行了校正，校正后的模型可替代FAO-56 P-M模型。Djaman等[7]分析了16個ET0計算模型在Senegal河谷的適用性，發(fā)現(xiàn)Hargreaves、修正的 Hargreaves、Ravazzani和Tralkovic模型高估了ET0，而Makkink-Hansen、 Oudin和Turc模型則低估了ET0，總體來看Valiantzas、Trabert、Romanenko、Schendel 和Mahringer模型精度較高，在該地區(qū)可代替FAO-56 P-M模型。隨著大數(shù)據(jù)時代的來臨以及計算能力的進步，不少機器學習方法被用于ET0計算研究[1,5,6,8]，然而其得到的為黑箱模型，沒有明確的算術表達式，不便于灌溉人員應用[9]，且通用性較弱[1]，也不能進行氣象因子敏感性分析[10]。

與以上純經驗模型方法不同，Valiantzas[11]利用全球2 000個氣象站點的氣象數(shù)據(jù)，基于數(shù)值模擬及數(shù)理統(tǒng)計方法，將FAO-56 P-M模型中一系列中間參數(shù)用氣象因子(輸入變量)的算術表達式表示，從而建立了不同氣象參數(shù)缺失下FAO-56 P-M的簡化模型(下文簡稱Valiantzas模型)。相比于FAO-56 P-M模型，Valiantzas模型更直觀、易用，計算更簡便，且便于進行氣象因子的敏感性分析，通用性更強，在一些地區(qū)發(fā)現(xiàn)其計算精度與FAO-56 P-M模型相當[12]。在Uganda濕潤氣候區(qū)的研究發(fā)現(xiàn)，缺少相對濕度數(shù)據(jù)的Valiantzas模型表現(xiàn)良好[13]，而缺少風速數(shù)據(jù)的Valiantzas模型由于使用了世界平均風速值替代，導致估算的ET0誤差較大。在伊朗的研究發(fā)現(xiàn)，Valiantzas模型最適用的氣象條件為太陽凈輻射大于24.2 MJ/(m2·d),氣溫為16-18℃，相對濕度為40%-50%，風速為1.50-2.50 m/s[14]。

廣西喀斯特地貌發(fā)育，因此，盡管該地區(qū)降雨充沛，但是由于特殊的地質背景(碳酸鹽巖分布范圍廣)，水文過程迅速[3]，地表水缺乏。準確估算ET0是水資源合理分配與管理的基礎。廣西面臨著氣象站點稀疏及氣象因子不完備等問題，且ET0在廣西喀斯特與非喀斯特地區(qū)具有不同的年際變化特征和影響因子[15]。本研究通過檢驗Valiantzas模型在廣西喀斯特與非喀斯特地區(qū)的適用性，并對相關模型進行校正，旨在為該地區(qū)農業(yè)水資源高效利用提供科學依據(jù)。

1 材料與方法

1.1 研究區(qū)概況與數(shù)據(jù)來源

廣西位于華南地區(qū)(104°26′-112°04′ E，20°54′-26°23′ N)，地處中國地勢第二階梯中的云貴高原東南邊緣，兩廣丘陵西部，地勢西北高、東南低，呈現(xiàn)西北向東南傾斜趨勢。陸地面積2.376×105km2，碳酸鹽巖面積約占廣西陸地總面積的33.9%。廣西屬亞熱帶和熱帶季風氣候，年平均氣溫為17.5-23.5℃，年平均降水量為841.2-3 387.5 mm。

本研究依據(jù)廣西地形地貌，選取10個氣象站點(5個為喀斯特站點，5個為非喀斯特站點[15])，共5年半(2010年1月-2015年6月)的逐日氣象數(shù)據(jù)(來源于中國氣象數(shù)據(jù)網http://data.cma.cn/)，包括日照時數(shù)n(h)、最高氣溫Tmax(℃)、最低氣溫Tmin(℃)、平均氣溫T(℃)、平均相對濕度RH、距地面2 m高處風速u2(m/s)。各站點位置及氣象因子平均值如表1所示。

1.2 ET0計算模型

采用表2中的模型對ET0進行模擬計算。其中，太陽凈輻射值由日照時數(shù)、地理緯度及日序數(shù)按經驗公式求得。Valiantzas-M1-Valiantzas-M5為不同氣象因子缺失下的Valiantzas簡化模型[11,12]，Turk[16]、Makkink[5]、Priestley-Taylor和Irmark-Allen為常用的經驗模型。

表1 氣象站點位置及氣象因子平均值(2010-2015年)Table 1 Location of meteorological stations and average value of meteorological factors (2010-2015)

表2中，Δ為飽和水汽壓-氣溫關系曲線在平均氣溫T處切線的斜率(kPa/℃)，γ為干濕計常數(shù)(kPa/℃)，λ為水的汽化潛熱(MJ/kg)，G為土壤熱通量[MJ/(m2·d)]，Tmax、Tmin和T分別為最高、最低和平均氣溫(℃)，TR=Tmax-Tmin(℃)，Tdew為露點溫度(℃)，α為常數(shù)0.23，RS為太陽輻射[MJ/(m2·d)]，Ra為大氣頂層輻射[MJ/(m2·d)]，Rn為太陽凈輻射[MJ/(m2·d)]，RH為平均相對濕度(%)，u2為距地面2 m高度處風速(m/s)，Z為氣象站海拔高度(m)，φ為緯度(rad)。Valiantzas-M4中，RH>50%時，Cu=0.076-0.0119(RH-50)0.2；RH≤50%時，Cu=0.076+0.0084(50-RH)0.2。Turk模型中，RH≤50%時，WRH取1；RH>50%時，WRH取0。

表2 參考作物蒸散量計算模型Table 2 Calculation model reference evapotranspiration

續(xù)表

Continued table

1.3 評價指標

本研究采用5種常用的統(tǒng)計指標：Kling-Gupta Efficiency系數(shù)(KGE)、Nash效率系數(shù)(NSE)、決定系數(shù)(R2)、均方根誤差(RMSE)以及平均偏差(PBIAS)，計算公式如下：

(1)

(2)

(3)

(4)

(5)

2 結果與分析

2.1 Valiantzas模型在喀斯特地區(qū)的適用性

表3給出了喀斯特5個站點在不同氣象因子缺失情況下，5個Valiantzas模型(含1個修正的模型)與FAO-56 P-M模型模擬ET0精度對比情況。Valiantzas-M1模型與FAO-56 P-M模型所需氣象因子相同，使用了所有的氣象數(shù)據(jù)，具有較高的KGE(0.735)、NSE (0.841)和R2(0.969)，以及最小的誤差(RMSE和PBIAS分別為0.557 mm/d和-40.785%)，因此其能更穩(wěn)健地模擬ET0。

Valiantzas-M2模型為原始缺失風速的模型(用2.00 m/s的世界平均風速來代替)，Valiantzas-M2-modify模型為基于廣西10個站點的平均風速(1.34 m/s)修正的Valiantzas-M2模型。由表3可見，基于Valiantzas-M2-modify得到的ET0精度得到了提高，R2變化較小，KGE和NSE分別提高了9.9%和16.6%，RMSE和PBIAS分別降低了27.4%和33.3%。表明廣西平均風速小于世界平均風速(2.00 m/s)，用廣西平均風速對Valiantzas-M2模型校正后能顯著提高ET0模擬效果。

Valiantzas-M3模型為缺失太陽輻射和風速的模型，其模擬精度最低。這可能是因為太陽輻射和風速是該區(qū)域ET0最主要的兩個影響因素[1]，同時缺失這兩個氣象因子導致ET0模擬出現(xiàn)較大誤差。

Valiantzas-M4模型也是缺失風速情況下的ET0計算模型，主要適用于濕潤地區(qū)[12]，其假定相對濕度較高時，空氣動力項對ET0的影響較小，反之亦然。據(jù)此引入依賴于相對濕度的調節(jié)系數(shù)waero，應用于風速缺失情況。從NSE、RMSE和PBIAS來看，Valiantzas-M4模型優(yōu)于Valiantzas-M2模型，NSE提高了1.9%，RMSE和PBIAS分別降低了1.4%和3.4%，表明在濕潤地區(qū)缺失風速情況下，Valiantzas-M4模型比Valiantzas-M2模型更適用。

Valiantzas-M5模型是相對濕度和風速缺失下的FAO-56 P-M簡化模型，通過建立相對濕度與氣溫的回歸關系并結合Valiantzas-M1模型得到，具有足夠高的精度。

表3 Valiantzas模型在喀斯特站點總體模擬精度Table 3 Overall simulation accuracy of Valiantzas models at the karst stations

由圖1可見，F(xiàn)AO-56 P-M模型與Valiantzas系列模型計算的ET0值在年內的變化趨勢較一致，均呈“幾”字形分布。各個模型計算的ET0值均在8月最大，分別為114.9 mm、136.8 mm、147.1 mm、139.2 mm、145.0 mm、146.6 mm和142.3 mm；在1月最小，分別為41.9 mm、43.7 mm、49.4 mm、44.2 mm、52.4 mm、49.2 mm和42.1 mm。Valiantzas系列模型均高估了各個月的ET0值，Valiantzas-M1、Valiantzas-M2、Valiantzas-M2-modify、Valiantzas-M3、Valiantzas-M4和Valiantzas-M5分別高估了16.8%、26.8%、17.9%、33.8%、25.9%和17.6%，可見Valiantzas-M1和Valiantzas-M5模型模擬的ET0值與FAO-56 P-M最接近。

圖1 2010/1-2015/6喀斯特地區(qū)FAO-56 P-M模型與Valiantzas模型模擬ET0的各月均值Fig.1 Monthly average value of ET0 stimulated by FAO-56 P-M model and Valiantzas models in karst area during Jan,2010 to Jun,2015

2.2 Valiantzas模型在非喀斯特地區(qū)的適用性

表4為不同氣象因子缺失下Valiantzas模型在非喀斯特地區(qū)的計算精度。Valiantzas-M1模型使用了所有的氣象數(shù)據(jù)，因此KGE (0.722)、NSE (0.796)和R2(0.964)均較高，而誤差(RMSE和PBIAS分別為0.641 mm/d和-50.264%)較低，因此Valiantzas-M1可靠性較高。

Valiantzas-M2、Valiantzas-M2-modify和Valiantzas-M4模型均為風速缺失情況下的ET0計算模型?；趶V西氣象站點風速修正后的Valiantzas-M2-modify模型能更好地模擬ET0，相比于Valiantzas-M2模型，其KGE和NSE分別提高了8.7%和17.4%，RMSE和PBIAS分別降低了24.9%和29.3%。Valiantzas-M4模型的模擬精度與Valiantzas-M2相當。綜合來看，缺失風速時，Valiantzas-M2-modify模型在非喀斯特地區(qū)較為適用。

同時缺失太陽輻射和風速、僅使用氣溫和相對濕度的Valiantzas-M3模擬精度較低，NSE最低(僅為0.411)且誤差較大(RMSE為0.988 mm/d，PBIAS為-82.747%)。

Valiantzas-M5模型僅使用氣溫和日照時數(shù)(相對濕度和風速缺失)，具有足夠的精度，NSE和R2分別達到0.753和0.953。

表4 Valiantzas模型在非喀斯特站點總體模擬精度Table 4 Overall simulation accuracy of Valiantzas models at the non-karst stations

由圖2可見，F(xiàn)AO-56 P-M模型與Valiantzas系列模型計算的ET0值在年內的變化趨勢較一致，均呈“幾”字形分布。各個模型計算的ET0值均在7月最大，分別為116.3 mm、142.2 mm、149.2 mm、142.2 mm、143.3 mm、150.0 mm和146.9 mm；在1月最小，分別為44.1 mm、48.5 mm、52.7 mm、47.8 mm、56.5 mm、52.0 mm和48.1 mm。Valiantzas系列模型均高估了各個月的ET0值，Valiantzas-M1、Valiantzas-M2、Valiantzas-M2-modify、Valiantzas-M3、Valiantzas-M4和Valiantzas-M5分別高估了19.5%、27.6%、19.5%、31.8%、26.9%和22.1%，可見Valiantzas-M1和Valiantzas-M2-modify模型模擬的ET0值與FAO-56 P-M模型最接近。

圖2 2010/1-2015/6非喀斯特地區(qū)FAO-56 P-M模型與Valiantzas模型模擬ET0的各月均值Fig.2 Monthly average value of ET0 stimulated by FAO-56 P-M model and Valiantzas models in non-karst area during Jan,2010 to Jun,2015

2.3 Valiantzas模型與其他經驗模型模擬精度比較

4個經驗模型在喀斯特站點的精度如表5所示。其中Turk模型與Valiantzas-M2、Valiantzas-M2-modify及Valiantzas-M4模型具有相同的輸入變量(Tmax，Tmin，n，RH)，其KGE和NSE均高于這3個Valiantzas簡化模型，且RMSE和PBIAS的絕對值均低于這3個Valiantzas簡化模型，說明在喀斯特地區(qū)缺失風速的情況下，Turk模型模擬ET0的可靠性高于相同輸入下的Valiantzas簡化模型。Makkink模型、Priestley-Taylor模型、Irmark-Allen與Valiantzas-M5模型具有相同的輸入變量(Tmax，Tmin，n)，可見Makkink精度最高，其次為Valiantzas-M5模型、Priestley-Taylor模型和Irmark-Allen模型，表明只有溫度和日照時數(shù)時，Makkink模型適用于喀斯特地區(qū)的ET0估算。

4個經驗模型在非喀斯特站點的精度如表6所示。其中Turk模型與Valiantzas-M2模型、Valiantzas-M2-modify模型、Valiantzas-M4模型具有相同的輸入變量(Tmax，Tmin，n，RH)，根據(jù)NSE、RMSE，Valiantzas-M2-modify模型精度最高，其次為Turk模型，說明在缺失風速的情況下，Valiantzas-M2-modify模型更適用于非喀斯特地區(qū)的ET0模擬。同理可見，Makkink模型、Priestley-Taylor模型、Irmark-Allen模型與Valiantzas-M5模型具有相同的輸入變量(Tmax，Tmin，n)，Makkink模型精度最高(具有最高的KGE、NSE和絕對值最低的RMSE、PBIAS)，其次為Valiantzas-M5模型、Priestley-Taylor模型和Irmark-Allen，表明只有溫度和日照時數(shù)時，Makkink模型是非喀斯特地區(qū)ET0的最佳估算模型。

表5 經驗模型在喀斯特站點總體模擬精度Table 5 Overall simulation accuracy of empirical models at the karst stations

表6 經驗模型在非喀斯特站點總體模擬精度Table 6 Overall simulation accuracy of empirical models at the non-karst stations

3 討論

桂西北以喀斯特峰叢洼地為基本景觀單元，該地貌水源涵養(yǎng)能力低、易旱、易澇，且“十里不同天”，氣象因子(降水量、溫度、濕度和光照等)時空變異大，這使得ET0經驗模型的使用受到限制。Valiantzas等[11,12]根據(jù)氣象因子之間存在一定的統(tǒng)計關系這一事實，對FAO-56 P-M模型進行了具有物理意義的簡化，因此與其他經驗模型相比具有較強的普適性。例如，從表3-表6可見，Valiantzas系列模型在喀斯特與非喀斯特地區(qū)表現(xiàn)較為穩(wěn)定，而Turk模型在喀斯特地區(qū)的表現(xiàn)優(yōu)于非喀斯特地區(qū)。這可能是由于Turk模型具有調節(jié)因子WRH(RH≤50%時，WRH取1；RH>50%時，WRH取0)，而RH在廣西喀斯特地區(qū)(桂西北)通常小于非喀斯特地區(qū)(桂東南)[1]，這使得該模型在喀斯特地區(qū)引入相對濕度的概率更大，從而在該地區(qū)具有更高的精度。

Wang等[3]統(tǒng)計了2010-2014年廣西24個氣象站點的多年平均氣象數(shù)據(jù)，結果顯示大多數(shù)站點(都安、鳳山、北海和潿洲島)的平均風速均小于世界平均風速(2.00 m/s)。本研究使用廣西喀斯特和非喀斯特10個站點的平均風速來代替風速缺失值，而非使用世界平均風速，發(fā)現(xiàn)能夠提高ET0的預測精度。Djaman等[18]也發(fā)現(xiàn)，在缺失風速數(shù)據(jù)的Sahelian地區(qū)，使用2.00 m/s的世界平均風速作為代替值時，Valiantzas模型高估了該地區(qū)的ET04%-35%，其RMSE為0.59-2.11 mm/d，誤差較大；使用該地區(qū)平均風速校正后，ET0模擬精度顯著提高，特別是對于較低的ET0值(ET0<6 mm/d)。下一步擬使用喀斯特和非喀斯特地區(qū)更多氣象站點和更長時間序列的風速值來進一步降低預測誤差。

基于5個模型評價指標，如喀斯特地區(qū)Valiantzas-M1模型、Valiantzas-M2模型和Valiantzas-M2-modify模型，R2分別為0.969，0.977和0.978，而RMSE值分別為0.557 mm/d、0.776 mm/d和0.563 mm/d，PBIAS分別為-40.785%、-65.119%和-43.460%；且Valiantzas-M2-modify模型和Valiantzas-M4模型的R2值(均為0.978)相同，但RMSE(分別為0.563 mm/d和0.765 mm/d)和PBIAS(分別為-43.460%和-62.921%)相差較大，可見R2對ET0模擬敏感性較差，不適合單獨用于ET0預測模型的評價，須結合NSE、KGE、RMSE和PBIAS等指標綜合評價。與喀斯特地區(qū)類似，在非喀斯特地區(qū)，Valiantzas-M2模型、Valiantzas-M2-modify模型和Valiantzas-M4模型的R2值均相同(0.980)，而NSE(0.691，0.811和0.703)、KGE(0.667，0.725和0.661)、RMSE(0.820 mm/d、0.616 mm/d和0.916 mm/d)和PBIAS(-70.675%、-49.939%和-68.980%)則相差較大，同樣說明R2值不適合單獨用于ET0預測模型的評價，這主要是由于R2值忽略了模型在低值區(qū)的表現(xiàn)[19]。因此，在評價ET0預測模型時，除了R2，還應結合其他評價指標，如NSE、KGE等。

從PBIAS值以及圖1、圖2可見，所有的Valiantzas模型(包括Valiantzas-M2的修正模型)均高估了ET0值。Valiantzas系列模型是基于全球2 000個氣象站點數(shù)據(jù)，并通過數(shù)值模擬方法簡化FAO-56 P-M模型得到的[11,12]，因此其ET0計算值與FAO-56 P-M模型的計算結果相比具有一定的誤差，表現(xiàn)為高估了ET0值。而在Senegal河流域的研究發(fā)現(xiàn)，Valiantzas系列模型低估了ET0值[20]，這可能是由于在不同氣候區(qū)，ET0的敏感性氣象要素不同。

4 結論

使用完備氣象要素的Valiantzas-M1模型在喀斯特和非喀斯特地區(qū)均表現(xiàn)最好，且其形式簡單，便于進行ET0對氣象因子的敏感性分析，可代替FAO-56 P-M模型。缺失風速時，在喀斯特地區(qū)，Turk模型較適用(RMSE為0.537 mm/d)；在非喀斯特地區(qū)，Valiantzas-M2-modify模型較適用(RMSE為0.616 mm/d)。同時缺失太陽輻射和風速時，Valiantzas-M3可用于ET0估算，在喀斯特和非喀斯特地區(qū)其RMSE分別為0.991 mm/d和0.988 mm/d。同時缺失相對濕度和風速時，Makkink模型最佳。此外，由于Valiantzas系列模型是基于全球氣象站點數(shù)據(jù)簡化FAO-56 P-M模型而來的，因此高估了廣西地區(qū)的ET0值，在應用時需進行校正。