西藏高海拔地区气象数据缺失条件下的ET0计算研究

汤鹏程，徐 冰，高占义，高晓瑜

（1.中国水利水电科学研究院 流域水循环模拟与调控国家重点实验室，北京 100038；2.中国水利水电科学研究院 牧区水利科学研究所，内蒙古 呼和浩特 010020；3.中国农业大学，北京 100083；4.内蒙古农业大学，内蒙古 呼和浩特 010018）

1 研究背景

作物腾发量（ET）是水文循环过程中很重要的一部分，ET是水量平衡过程的关键环节，同时也是表面能量平衡过程的重要组成，对于决定灌水需求、进行区域水平衡研究和降雨-径流、生态模型的实施起着关键的作用。计算作物腾发量的第一步就是计算参照作物腾发量（ET0），所以不精确的ET0计算将会导致不合理的模型率定和不可行的地下水补给估算［1］。虽然已经有很多估算ET0的方法，但是各种方法的实施一般都局限于特殊的地理与气象条件［2-3］。很多不同气候区的结果都证实了FAO56 Penman-Monteith模型（以下简称PM）的可行性，所以联合国粮农组织将PM法作为一种在全球适用的计算ET0的标准方法。

广泛应用PM模型的主要障碍是复杂的计算过程和大量的气象数据获取，包括太阳辐射、风速、湿度、温度等［4］，尤其是在气象数据缺乏地区，很多情况下都无法得到连续且全面的数据。对于高海拔、地形陡峭等环境比较恶劣的地区，全套气象站的安装和维护较为昂贵与复杂，因此，探索精度高、计算简便的适合高海拔地区的ET0计算方法有很强的实际意义。数据缺失情况下，ET0计算的经验公式基本分为基于温度的、基于蒸发的、基于物质转换的和混合型的，其中基于温度的简化模型具有较为突出的简便性与计算数据易获得性，因而被广泛采用［5］。

Hargreaves和Samani于1985年采用最高、最低温度和外辐射数据计算得到太阳辐射，然后提出了基于温度的“Hargreaves-Samani（HS）模型”［6］。Almorox等在2015年采用4 362个气象站的数据，应用11种基于温度的ET0计算公式对ET0进行了估算，发现在干旱半干旱地区HS方法是最精准的简便计算方法［5］。Er-Raki等在2010年的研究也表明HS方法在美洲半干旱地区是基于温度的简化模型中最精确的方法［7］。国内一些学者对国际范围内常用的ET0经验模型在我国典型地区的适用性进行了分析：胡庆芳等基于全国105个气象站，在月时间尺度上评价了HS模型在我国不同气候区适用性［8］；除此之外，王声锋等（2008）［9］、王新华等（2006）［10］分别在我国内陆半干旱地区、西北干旱地区验证了HS模型的适用性。由于HS方法的简便性，对于ET0的估算，HS方法很受欢迎。然而Jensen等的研究发现在干旱地区HS会低估ET0，在湿润地区HS会高估ET0［11］。Martinez和Thepadia在2009年用佛罗里达72个站点的数据计算，发现HS方法过高估算了ET0［12］。Yoder等于2005年对美国东南湿润地区Cumberland高原的日ET0和周ET0进行了估算，发现相对于PM方法，HS方法计算的ET0偏大，HS方法更适合长时间尺度ET0的估算［13］。总体来说，国内外对HS等经验模型的适用性评价较多但修正研究相对较少，修正模型在国内外的推广和应用更加欠缺。

很多ET0计算经验公式都只适用于特定的气候和区域［14-15］。本研究旨在基于HS温度法的基础上，找到一种适合于西藏高海拔极端环境地区（海拔2000 m以上地区）气象数据缺失条件下的ET0简便、精准的计算方法。海拔因子是在ET0计算中最容易获得参数，其不需要连续观测，计算时没必要针对不同时间尺度进行基础数据整理。同时海拔因子也是最容易被忽略的参数，在PM推荐标准计算公式中，海拔因子与γ（为湿度计常数）、Ra（天顶辐射，地球大气层顶部水平面吸收的太阳辐射）、Rn（太阳净辐射，地球表面吸收的能量）的计算均有直接函数关系。目前已有简化模型研究中多针对温度、大气相对湿度、风速、日照时数、降雨等气象因子与ET0的相关关系开展，往往忽略了ET0计算的空间变异性，因此直接导致简化模型在不同地区间应用推广的难度加大，往往区域不同经验模型里的很多参数就要重新校正。例如Allen指出HS模型中的温度系数取值在海拔超过1 500 m的地区并不完全合理；Annandale等认为HS模型中的温度系数、温度指数均应考虑不同地区的大气压而进行修正，然而大气压强的改变与海拔成明显负相关关系［16］。因此本文考虑海拔因素，对HS模型进行修正提高了原模型对不同区域空间变化的响应能力。

2 材料与方法

2.1 研究区概况与数据资料 西藏位于中国西南部，素有“世界屋脊”之称，平均海拔4 000 m，面积123万km2，低氧低压（不足海平面的2/3）、日照长（多在3 000 h以上）、辐射强（年太阳辐射6 000～8 000 MJ/m2）是其主要气候特点，由于西藏高海拔地区近地层冷热交换频繁，导致该地区空气温、湿度变化大，干湿季分明。本文根据西藏地理地貌与气候特点，充分考虑西藏全区海拔高程变化（整体趋势西高东低；本研究不考虑西藏西北及东南海拔高程6 000 m以上的喜马拉雅高山区），针对西藏主要农牧业生产区（海拔高程5 000 m以内），在西藏全区选取9个代表性站点开展研究（图1），9个典型站点海拔跨度2 000～5 000 m，其地理位置信息见表1。气象资料均来自国家气象信息中心，数据经过严格控制，质量较好。本文以1981—1990年逐日气象资料（n=32 868，无插补延长）进行先验研究，以1991—2000年逐日（n=32 877，无插补延长）、逐月（n=1 080）气象资料进行模型检验。

2.2 参考作物腾发量计算

2.2.1 FAO56 Penman-Monteith公式 国际粮农组织认为PM公式所有的计算程序都能够通过可得到的气象资料和时间尺度的计算公式得以标准化，相关计算成果可被项目管理者、项目咨询者、灌溉工程师、水文学者、农艺学者、气象学者作为参考标准使用。PM公式考虑了空气动力学项（ET0aero）

表1 西藏全区9个典型站点地理位置信息

图1 西藏全区9个典型站点分布

与辐射项（ET0rad），常被作为适合于大部分气候地区的标准方法使用。

式中：ET0为参照作物腾发量，mm·d-1；Rn为作物冠层表面净辐射，MJ·m-2·d-1；G为土壤热通量，MJ·m-2·d-1；T为平均温度，℃；u2为高度2 m处风速，m/s；es为饱和水汽压，kPa；ea为实际水汽压，kPa；Δ为饱和水汽压与温度曲线的斜率，kPa·℃-1；γ为湿度计常数，kPa·℃-1。

2.2.2 Hargreaves公式 HS方程是Hargreaves基于8年内实测蒸渗仪实测试验推导出的仅利用温度数据来计算参照作物腾发量（ET0）的方法，HS方程作为一种“温度法”巧妙地利用大气温度差（Tmax-Tmin，表明天顶辐射有多少能达到地球表面的一个指标）计算太阳辐射，国内外众多研究成果表明其作为气象资料缺失情况下估算ET0的方法，可以给出全球较为有效合理的ET0，但其仍保留着经验系数，如下式中的温度系数（0.0023）、温度常数（-17.8）与温度指数（0.5）。HS方程中的经验参数均基于多年实测试验由大量数据回归分析得到的，因此在每一个新的地区都应对HS方程中参数进行修正后使用。

式中：Tmax、Tmin、Tmean分别表示日最高、最低、平均温度；Ra为天顶辐射，mm/d。可根据地理位置与时间直接查表或计算。

2.3 高海拔地区ET0计算主要影响因素分析

2.3.1 标准PM方程输入因子比较 PM标准方程计算过程中输入的连续数据系列包括日最高温度、日最低温度、平均温度、平均风速、平均相对湿度、日照时数等7个主要气象影响因子；以及非连续数据系列包括经度、纬度与海拔高程在内的3个地理位置信息。计算过程中需要参数较多，且需要气象数据连续，计算公式也较为复杂。

主成分分析是在损失较少有效信息的基础上将多指标转化为少数几个综合指标的多元统计分析方法，通常把转化生成的综合指标称为主成分，其中主成分都是原始变量的线性组合，且各个主成分之间互不相关［16］。本文基于西藏高海拔地区9个气象站点地理位置信息以及20年逐日气象资料，利用SPSS软件进行主成分分析及统计，通过该软件系统自动将原始数据标准化处理，消除各指标量纲与数量级的影响。结果表明：

（1）总方差解释：基于特征值不小于1的原则，由表2可知，第一、第二、第三主成分的累积方差贡献率可达77.29%，可知前3个主成分包含了原始变量的大部分信息，可以在一定程度上代替9个原始指标的影响关系。

（2）成分矩阵：由成分矩阵分析可知（表3），3个主成分与原始指标间的相关程度。经观察，第一主成分与日最低温度、日平均温度、日最高温度、海拔4个初始变量在同组之间具有显著的相关性。第二、第三主成分仅分别与日照时数、纬度具有较好的相关性。第一主成分的方差贡献率达到46.375%（表2），即该成分包含了接近一半的信息。

表2 总方差解释

表3 成分矩阵分析

综上所述日最低温度、日平均温度、日最高温度和海拔4个指标包含了原始变量多数信息且4项因子之间隶属于同一组（即第一主成分），是较为重要的输入指标。同时考虑到海拔因素对ET0的重要影响，且海拔因子不需要连续的观测与修正，易于获取，本研究在HS方法中引入海拔因子，构建适用于高海拔地区ET0计算的HS-E改进模型。

2.3.2 海拔及温度因素对ET0计算的影响 在图2中，横坐标表示海拔（海拔跨度介于2 000～5 000 m之间），纵坐标表示年累计ET0：

（1）4 000 m以下的地区，随着年平均日气温的逐渐降低（见表1），年累计ET0逐渐升高；4 000 m以上地区，随着年平均日气温的逐渐降低（见表1），年累计ET0呈现逐渐降低的趋势。年累积ET0最大值（平均约1 280 mm/a）会出现在泽当或拉孜地区（海拔4 000m左右），最小值（平均约950 mm/a）会出现在察隅（本研究所选典型站点中海拔最低点，2 000 m左右）与安多地区（本研究所选典型站点中海拔最高点，5 000 m左右）。

（2）当雄地区（海拔4 200 m）ET0及年平均日气温均出现反常，其年累积ET0相比较定日（海拔4 300 m）、拉孜（海拔4 000 m）较低，但由于受地形地貌影响，导致其多年内日平均温度（1.81℃）低于定日（3.12℃），拉孜（6.54℃），这直接导致当雄地区年累计蒸散发量较小。

综合2.3.2节高海拔地区ET0计算及2.3.1节主成分分析结果可知，ET0同时受温度与海拔因素影响明显，针对高海拔地区ET0计算，采用温度简化算法的同时引入海拔因素，完全可取且方向正确。

图2 高海拔地区年累积ET0与海拔因子响应关系

3 Hargreaves-Elevation改进模型

3.1 数学模型的建立 通过主成分分析，得知海拔因子是一个重要的输入指标，其获取相对容易且不需要连续观测，鉴于此在HS公式中引入海拔因子，改进后的数学模型为：

式中：f（H）为待求的海拔函数；a为温度常数，修正模型中为避免ET0HSE计算出现负数，a取近50年来海拔2 000 m以上地区的Tmean最小值；本研究充分考虑HS模型的构成及其参数来源的基础上对其进一步改进，其它符号同前。

3.2 温度常数确定与海拔函数推求

3.2.1 温度常数确定 FAO推荐当Ra单位为mm/d时，温度常数a=-17.8，然而在高海拔地区，尤其是海拔3 500 m以上地区1月、12月两个较为寒冷的月份，日Tmean经常出现低于-17.8℃的情况，直接导致由HS模型计算的ET0出现负值，该结果违背客观自然规律（实际蒸散发不可能为负），计算结果不具有参考性。

本研究梳理1960—2015年内9个代表性站点（站点信息参考表1，海拔跨度2 000～5 000 m）最低日Tmean数据，发现1987年12月西藏高海拔地区接连出现极端低温气候，其中海拔4 414.9 m的改则地区于1987年12月26日出现近50余年内极端最低日Tmean=-36.6℃。综上所述，为避免新模型ET0计算出现负值错误（HS-E及HS模型中仅会在Tmean-a计算项中出现负数），本文将温度常数a为-36.6。

3.2.2 海拔函数推求 设：

以9个典型站点参照作物腾发量计算的月累计值（ET0PM/ET0X）为因变量，以海拔因子H为自变量，将ET0PM/ET0X的月累计值与实际海拔因子H带入计算公式，通过回归分析得到回归趋势方程如下（R2=0.26***，***代表在0.1%水平上显著；标准误差=1.58×10-4）：

综上所述，通过回归建立的HS-E改进模型计算公式为：

4 HS-E改进模型评价

4.1 模型验证 本研究中采用1981年到1990年10年的数据进行回归分析，得到一个考虑海拔因素的新方法HS-E模型，针对9个代表性站点1990—2000年不同时间尺度条件下ET0计算结果来分析HS-E模型在高海拔地区适用性。

验证HS-E模型适用性用到平均相对误差（MRE）、均方根误差（RMSE）和纳什系数（NSE）作为评价指标，各指标计算过程如下：

图3 HS-E改进模型回归分析（n=120×9）

式中：N表示误差对比点总数；Pi和Oi表示预测值和实测值（i=1、2、…、N）。

MRE是一个百分数，用来表征预测的精度性及误差范围；RMSE用来衡量观测值同真值之间的偏差，可以用来衡量一个数据集的离散程度，二者越接近于0，模型的预测质量越高。NES值越接近1，表示模型质量越好，模型可信度高；结果接近0，表示模拟结果接近观测值的平均值水平，即总体结果可信，但过程模拟误差大；结果远远小于0，则模型是不可信的。

4.2 逐日ET0对比分析 通过对不同海拔地区9个气象站点气象数据分析计算，分别应用PM、HS、HS-E模型计算得到逐日ET0（样本数N=3 653×9），将2种简化温度算法（HS、HS-E）计算结果与PM标准法计算结果进行RMSE、MRE、NSE分析。

新方法HS-E模型在仅基于温度数据、地理位置数据的基础上计算得到ET0，其结果较为接近标准值，具体如表4所示：

（1）HS-E模型NSE介于0.69～0.88之间，平均0.80，模型质量较好；HS模型NSE系数介于0.49～0.81之间，平均0.68（＜0.8），HS-E模型模拟结果优于HS模型，HS-E模型质量更高。

（2）HS-E模型RMSE介于0.41～0.65 mm/d，平均0.53 mm/d；HS模型RMSE介于0.52～0.94 mm/d，平均0.68 mm/d，经由HS-E模型的计算值同真值之间的偏差更小。

（3）HS-E模型MRE分析相对误差范围介于9.3%～22.91%之间，平均值（MRE绝对值平均）13.8%；HS-E模型MRE表现优于HS模型-39.96%～21.22%（绝对值平均19.3），新模型避免了负偏差的出现，即纠正了ET0出现负值的错误情况，预测精度更高，预测值更接近PM标准方程计算值。

（4）HS-E模型回归斜率介于1.05至1.24之间，平均1.14，略微呈现正偏差；对比HS模型（回归斜率0.63～1.21），HS-E模型表现更加稳定。

通过对比ET0逐日数据可知，HS-E修正式比HS模型在高海拔地区表现出更好的稳定性与适应性。

4.3 逐月ET0对比分析分别应用PM、HS、HS-E模型计算得到逐月ET0（样本数N=120×9）。将HS-E修正模型、HS模型分别与标准PM公式计算的逐月ET0值进行对比分析（具体结果如图4所示）。基于NSE模型质量分析、RMSE模型误差分析可知：

表4 HS模型与HS-E模型在高海拔地区逐日ET0分析比较

（1）在所有不同海拔的9个不同站点间，综合对比NSE（0.84，月尺度＞0.80，日尺度）、RMSE（11.90 mm/月优于0.53 mm/d）、MRE（12.50＜13.80）不同时间尺度条件下的误差分析结果，可知新模型HS-E逐月ET0计算结果更优于逐日ET0计算结果，因此随着计算时段增长HS-E模型误差更小，且对比HS模型优势更加明显。

（2）HS-E模型NSE介于0.57～0.98之间，平均0.84，接近于1；HS模型NSE系数介于-0.86～0.90之间，平均0.54（＜0.84），HS-E模型模拟结果优于HS模型，HS-E模型质量更高。

（3）HS-E模型RMSE介于6.55～21.92 mm/月，平均11.90 mm/月；HS模型RMSE介于11.53～45.43 mm/月，平均20.00 mm/月；经由HS-E模型的逐月ET0计算值同真值之间的离散程度更小。

（4）HS-E模型MRE分析相对误差范围介于-24.48%（安多）至24.92%（那曲）之间，平均值（MRE绝对值平均）为12.50%；HS模型-53.11%（安多）至19.73%（察隅），平均值（MRE绝对值平均）21.40%，新模型误差范围较小，计算结果更接近实际情况。

图4 HS模型与HS-E模型在高海拔地区逐月ET0分析比较

5 讨论与结论

5.1 讨论 现有研究中对HS模型的改进存在一定不足，已有成果多为提高HS模型计算精度而引入大气平均相对湿度、日照时数等连续系列气象参数作为模型新增输入项［17］，ET0计算精度提高的同时导致模型计算对基础气象数据需求度更高以及计算过程的复杂化，忽略了温度法简便、少参的本质。

温度计算法HS模型巧妙地利用大气温度差（Tmax-Tmin）作为表达天顶辐射有多少能达到地球表面的一个指标，然而由于地理位置的特殊性导致高海拔地区强辐射的同时温度并不是很高，准确地说是由于高海拔地区大气稀薄，对地面长波辐射的吸收较少，导致热量大量散失，大气温度较低；同时海拔高的地方云层较少，白天云层吸收地面长波辐射较少，夜间云层对地面的逆辐射作用被削弱导致保温作用较差，这些原因均削弱了温度差（Tmax-Tmin）对辐射的敏感性，因此引入海拔因子对HS模型进行改进，可显著提高模型计算精度。

本文研究区域主要集中于西藏高海拔地区，建立的HS-E方程主要适用范围也同样具有地区性。但本文提出的HS-E修正模型通过引入海拔因素对现有的温度法HS模型进行修正，以期为今后在更大范围的高寒区域完善和拓展“ET0简便温度计算法”提供一种思路与方法。为此本文提出HS-E修正模型：在其中增设经验系数a、b、c。在西藏海拔2 000 m以上的大部分地区经验系数a为-8×10-6；b为0.07；c为5。但在在西藏西部喜马拉雅山区（温度骤降，常出现水汽冻结的现象）以及其他高海拔地带该系数需要回归修正后采用。

5.2 结论

（1）基于西藏高海拔地区9个气象站点地理位置信息以及20年逐日气象资料，通过主成分分析可知，在海拔2 000～5 000 m地区，日最低温度、日平均温度、日最高温度和海拔4个指标包含了PM方程原始变量多数信息，且隶属于同一组（即第一主成分），是较为重要的输入指标。鉴于海拔因子不需要连续观测与修正，易于获取，将其引入HS方法进行高海拔地区ET0计算较适宜。

（3）HS-E修正模型改进了HS模型中的温度常数（-17.8），经过梳理西藏典型地区1960年-2015年56年内最低日Tmean数据，将温度场数修订为-36.6，避免了原HS模型在高海拔地区ET0计算出现负值的情况出现，提升ET0计算结果的实用性与精度。

（4）HS-E模型在不同时间尺度条件下（日、月）计算精度均优于HS模型，新模型NSE分析分别达到0.80（日）与0.84（月），结果较优，模型可信度较高；同时综合对比NSE、RMSE、MRE在不同时间尺度条件下的分析结果（0.84＞0.80；11.90 mm/月优于0.53 mm/d；12.50＜13.80），相比于短时间段逐日ET0计算结果，长时间段逐月ET0计算结果更接近于PM标准公式，误差更小，因此计算时间尺度越大HS-E模型结果越优。

参 考 文 献：

［1］ ALLEN R G，PEREIRA L S，RAES D，et al.Crop evapotranspiration guidelines for computing crop water re⁃quirements［M］.FAO Irrigation and Drainage.Paper 56，Roma，Italy，1998.

［2］ IRMAK S，ALLEN R G，WHITTY E B.Daily grass and alfalfa-reference evapotranspiration estimates and alfalfato-grass evapotranspiration ratios in Florida［J］.Journal of Irrigation and Drainage Engineering，2003（129）：360-370.

［3］ 刘钰，PEREIRA L S.气象数据缺测条件下参照腾发量的计算方法［J］.水利学报，2001（3）：1-17.

［4］ VALIANTZAS J D.Simplified limited data Penman’s ET0 formulas adapted for humid locations［J］.Journal of Hydrology，2015（524）：701-707.

［5］ ALMOROX J，QUEJ V H，MARTI P.Global performance ranking of temperature-based approaches for evapo⁃transpiration estimation considering Koppen climate classes［J］.Journal of Hydrology，2015（528）：514-522.

［6］ HARGREAVES G H，ALLEN R G.History and evaluation of Hargreaves evapotranspiration equation［J］.Jour⁃nal of Irrigation and Drainage Engineering，2003，129（1）：53-63.

［7］ ER-RAKI S，CHEHBOUNI A，KHABBA S，et al.Assessment of reference evapotranspiration methods in semiarid regions：can weather forecast data be used as alternate of ground meteorological parameters［J］.Journal of Ar⁃id Environment，2010（74）：1587-1596.

［8］ 胡庆芳，杨大文，王银堂，等.Hargreaves公式的全局校正及适用性评价［J］.水科学进展，2011，22（2）：160-167.

［9］ 王声峰，段爱旺，张展羽.半干旱地区不同水文年Hargreaves和P-M公式的对比分析［J］.农业工程学报，2008，24（7）：29-33.

［10］ 王新华，郭美华，徐中民.分别利用Hargreaves和PM公式计算西北干旱区ET0的比较［J］.农业工程学报，2006，22（10）：21-25.

［11］ JENSEN M E，BURMAN R D，ALLEN R G.Evapotranspiration and irrigation water requirements［M］//ASCE manuals and reports on engineering practice NO.70，ASCE，NewYork，1990.

［12］ MARTINEZ C J，THEPADIA M.Estimating reference evapotranspiration with minimum data in Florida［J］.Jour⁃nal of Irrigation and Drainage Engineering，2009，136（7）：494-501.

［13］ YODER R E，ODHIAMBO LO，WRIGHT W C.Evaluation of methods for estimating daily reference crop evapo⁃transpiration at a site in the humid southeast united states［J］.Applied Engineering in Agriculture，2005，21（2）：197-202.

［14］ 李晨，崔宁博，魏新平，等.改进Hargreaves模型估算川中丘陵区参考作物蒸散量［J］.农业工程学报，2015，31（11）：129-134.

［15］ JABLOUM M，SAHLI A.Evaluation of FAO56 methodology for estimating reference evapotranspiration using lim⁃ited climatic data application to Tunisia［J］.Agricultural Water Management，2008（95）：707-715.

［16］ 杨永红，张展羽.改进Hargreaves方法计算拉萨参考作物蒸发蒸腾量［J］.水科学进展，2009，20（5）：614-618.

［17］ ANNANDALE J G，JOVANIC N Z，BENADE N，et al.Software for missing data error analysis of Penman-Mon⁃teith reference evapotranspiration［J］.Irrig Sci，2002，21（2）：57-67.