(河海大学 水文水资源与水利工程科学国家重点实验室,江苏 南京 210098)
潜在蒸散发(Potential Evapotranspiration,记为ETp)是指充分供水条件下的区域蒸散发能力。该变量是影响区域水热平衡的重要因素,是水文过程中的重要环节,不仅可以指示一个区域的干湿状况,也对分析地区水文平衡过程、气候变化等具有重要意义。
不同研究者对于ETp的定义存在差异[1]。世界气象组织(WMO)定义ETp为特定气候条件下植物完全覆盖且充分供水时地表的蒸散发量[2],但对植物类型未做定义;联合国粮农组织(FAO)则以平均高度0.12 m,反照率为0.23,冠层阻力为70 s/m的植被作为“参考作物”,并将该植被在地表充分供水时的蒸散量称为“参考蒸散发”[3]。从这两个定义来看,FAO提出的“参考蒸散发”是将WMO对“潜在蒸散发”定义中的植被条件做了进一步的补充,可以认为“参考蒸散发”是一个特定条件下的“潜在蒸散发”。而关于ETp的估算,目前没有一个可以直接获取广阔自然水面或地表ETp的方法,通常是利用间接的手段进行估算[2],例如以蒸发皿蒸发量作为ETp观测值使用,或使用气象站点的观测资料结合各种蒸散公式(Penman-Moneith、Priestley-Taylor和Hargreaves-Samani等)来估算ETp[1]。而在全球50多种ETp的估算方法和模型中,修正的FAO56 Penman-Moneith(P-M)方程计算结果被认为是最可靠的[4],因此基于该方法计算的ETp结果常被作为标准值来评估其他方法的ETp估算结果。但无论是蒸发皿观测值还是公式的计算结果,当区域上站点分布不均且数量较少时,会使区域ETp的估算结果缺乏代表性。近几十年来,由于卫星遥感技术的发展 ,遥感蒸散估算方法开始使用并逐步发展完善[5]。
现有的蒸散发产品(包括实际蒸散发ETa和潜在蒸散发ETp)可以大致分为模型模拟和遥感数据反演两大类[6],前者的代表性产品有NASA基于全球气象数据以及多个陆面模型模拟得到的全球陆面数据同化系统(GLDAS)[7];而后者的代表性产品有由英国布里斯托大学地理科学学院的Miralles博士发布的GLEAM产品[8]、由美国航空航天局(NASA)发布的基于MODIS卫星遥感数据制作的全球陆地蒸散产品MOD16[9]和美国地质调查局(USGS)发布的基于业务化简化地表能量模型的SSEBop/ET产品[10]等。这些产品数据集都同时提供了ETa与ETp数据。目前,这些数据集中的ETa被广泛应用于全球的蒸散发时空变化及干湿变化规律研究当中,并在很多地区进行了评估,而对潜在蒸散发产品评估很少。例如赵燊等[11]在分析山东省蒸散及潜在蒸散时空变化的过程中,仅从时空相关性的角度分析了MOD16 ETp产品的精度,发现其与蒸发皿观测数据的空间相关系数为0.872,时间相关系数高达0.913,但对于该产品数据质量的时空特性缺乏系统性认识。
ETp作为许多水文和生态模型的重要输入参数,它的准确估算是展开相关研究的基础,并且,有些产品(例如GLEAM及SSEBop/ET)的ETa数据是基于ETp进一步计算得到的。因此,有必要对不同ETp产品在各地的运用和分析进行质量评估,为相关应用研究提供数据可靠性的支撑。本文利用2个版本的GLDAS气象驱动数据,通过FAO56 Penman-Monteith公式计算得到云南和贵州2000~2014年的逐月ETp,之后利用研究区内站点观测蒸发皿蒸发量数据,对基于Penman-Monteith公式的ETp计算结果以及GLDAS-2.1和MODIS MOD16的月尺度ETp产品,从时间和空间分布特征上进行数据质量的评估,为西南地区水资源管理和干旱监测等提供数据支撑。
贵州和云南是典型的喀斯特地貌地区,总面积约57万km2。该地区受到热带和亚热带季风的影响,气候类型多样并具有明显的区域差异。气温年较差较小,日较差较大,降水时空分布不均,干湿季分明。
由图1数字高程(DEM)数据(来源:http://srtm.csi.cgiar.org/SELECTION/inputCoord.asp)可看出整个云贵地区地形高差达到6 000 m。贵州地区的地势主要表现为西高东低,云南地区的地势则呈现出北高南低的特点。
1.2.1MODIS MOD16数据
由美国航空航天局(NASA)发布的基于MODIS卫星遥感数据制作的全球陆地蒸散量产品MOD16[12](http://files.ntsg.umt.edu/data/NTSG_Products/MOD16/)主要包含了陆面蒸散发数据、潜在蒸散发数据、潜热通量等数据。目前NASA已发布了版本6的MOD16数据,但因该版本的数据在云贵地区缺失值较多,故本文使用MOD16版本5的数据,由蒙大拿大学密苏拉分校地球动态数值模拟研究组(Numerical Terra dynamic Simulation Group,NTSG)提供。NTSG将空间分辨率为1 km、时间分辨率8 d的MOD16数据处理后,得到月和年尺度的0.05°×0.05°,0.5°×0.5°数据。本文使用的是2000~2014年0.05°×0.05°的月尺度潜在蒸散发量(记为MOD16_ETp)数据。
图1 云贵地区高程及气象站点分布Fig.1 Elevation and meteorological stations in Yunnan-Guizhou region
1.2.2全球陆面数据同化系统数据
由美国NASA戈达德空间飞行中心(GSFC)和美国海洋和大气局(NOAA)国家环境预报中心(NCEP)联合开发的全球陆面同化系统GLDAS(Global Land Data Assimilation System)用再分析气象数据离线驱动多个陆面模型(包括Noah模型、VIC模型等),得到全球范围的地表状态变量[13]。GLDAS-2产品(http://disc.sci.gsfc.nasa.gov/hydrology/data-holdings)包含了GLDAS-2.0和GLDAS-2.1两种数据集,其中,GLDAS-2.0由普林斯顿大学制作的全球气象数据,驱动Noah或CLM模型生成1948~2014年的1°×1°以及0.25°×0.25°的数据集[14];GLDAS-2.1产品则是基于NOAA/GDAS大气数据、全球降水气候学计划(GPCP)的降水数据以及由美国空军气象局农业气象模拟系统(AGRMET)生成的辐射数据,共同驱动Noah模型得到的2000年至今1°×1°和0.25°×0.25°的数据集。GLDAS-2.0和GLDAS-2.1产品都是按3h时间尺度生成各气象要素,并以此为基础合成月尺度数据[15]。本文所用数据为GLDAS 2.1的月尺度潜在蒸散发数据,空间分辨率为0.25°×0.25°,因潜在蒸散发数据的单位为W/m2,需利用潜热通量λ=2.45 MJ/kg将其换算,根据各月的日数将其换算为月潜在蒸散发量(记为GLDAS_ETp)。此外,由于GLDAS-2.0和GLDAS-2.1的气象驱动数据不同,且GLDAS-2.0只提供了实际蒸散发数据,因此本文结合GLDAS-2.0和GLDAS-2.1的3 h尺度产品中的长波辐射、短波辐射、气压、气温、风速和比湿等数据,运用FAO推荐的P-M公式求得逐日ETp(分别用PM_ETp_2.0和PM_ETp_2.1表示)后,再累加得到逐月PM_ETp_2.0和PM_ETp_2.1。
1.2.3蒸发皿数据
中国气象数据网地面气候资料日值数据集V3.0(http://data.cma.cn/site/index.html)提供了国家基准气象站1951~2017年逐日数据,由于大型蒸发皿数据在每年1~3月、11~12月存在缺失,故本文选取了云贵地区范围内42个站点(见图1)2000~2014年的逐日小型蒸发皿蒸发量,通过累加得到各站点的月蒸发皿蒸散发数据,并以该蒸发量作为参考值(记为Epan),用于检验其他潜在蒸散发数据的精度。
1.3.1潜在蒸散发的计算
采用P-M公式[3],利用GLDAS-2.0、GLDAS-2.1驱动数据中的长波辐射、短波辐射、气压、气温、风速和比湿等数据计算ETp,其表达式为
(1)
式中,Δ表示饱和水汽压-温度曲线的斜率,kPa·℃-1;γ表示干湿表常数,kPa·℃-1;Rn表示参考作物表面的净辐射,MJ·m-2·d-1;G表示土壤热通量密度,MJ·m-2·d-1,因为日均G相对Rn数值很小,一般默认G=0;T表示2 m高度处的气温,℃;u2表示2 m高度处的风速,m·s-1;es和ea分别表示饱和水汽压和实际水汽压,故es-ea表示饱和水汽压差,kPa。
1.3.2评估方法
采用相关系数(r)、平均偏差(MBE)、绝对误差(MRAE)和均方根误差(RMSE),基于站点分别从时间和空间分布特征上对MOD16、GLDAS 2.1以及基于P-M公式计算所得的ETp数据进行评估。各指标的计算公式如下:
(2)
(3)
(4)
(5)
式中,xi表示待评估ETp数据(包括MOD16_ETp、GLDAS_ETp、PM_ETp_2.0和PM_ETp_2.1),yi表示气象站实测值Epan,n为统计样本数。
42个站点实测Epan均值和待评估的4套ETp数据在云贵地区的逐月均值变化情况如图2所示。由图2可知,Epan与4套ETp数据的逐月变化基本一致,但GLDAS_ETp的数据明显高于其他ETp值,这种高估情况在3~5月最明显;PM_ETp_2.1与MOD16_ETp的一致性较好,两者在大部分时段约有±10 mm的偏差,在其余时段下的变化情况相似;PM_ETp_2.0比实测值存在低估的情况,在2002,2005,2009,2012等年份的冬季(11,12月),PM_ETp_2.0与Epan偏差可达约-20 mm,而其余时段内与Epan的匹配程度较好。
图2 2000~2014年云贵地区月潜在蒸散发变化情况Fig.2 Variation of monthly ETp during 2000~2014 in Yunnan-Guizhou region
综合对比发现,云南和贵州地区ETp数据大小排序为GLDAS_ETp,MOD16_ETp,PM_ETp_2.1,Epan,PM_ETp_2.0的特征,说明GLDAS 2.1和MOD16的ETp产品存在高估现象,PM_ETp_2.0最接近于实测Epan数据。进一步利用相关系数和均方根误差等指标对4套数据进行评估,结果列于表1,可以看出,4套数据与Epan都有较好的相关性,均能较好地反映地区ETp的时程变化特征,其中PM_ETp_2.0与Epan相关性最好,MBE及MRAE也最小,而MOD16和GLDAS-2.1的ETp产品则与Epan有较大的偏差,存在严重高估。
4套数据4个误差指标的逐月变化情况见图3,结
图3 不同ETp数据的r、RMSE、MBE、MRAE季节变化Fig.3 Seasonal variation of r , RMSE, MBE and MRAE of different ETp in comparison with Epan
果表明各套数据与Epan的相关性具有明显的季节变化特性。总体看,各月的PM_ETp_2.0与实测值的相关性不如其他数据,而2.1节中该数据与逐月Epan相关性最好,这说明PM_ETp_2.0能够较好地拟合出Epan的季节性波动情况,但在反映逐月Epan的年际波动情况上不如其他3套数据。4套数据相关系数的逐月变化显示,一年中,2,3月和8,9月的r均大于0.7,而在6月和10~12月,4套ETp数据的相关系数均出现了整体的下降。
表1 云贵地区4套ETp评估结果Tab.1 Evaluation results of four ETp data sets in Yunnan-Guizhou region
RMSE、MBE和MRAE的逐月变化情况较为相似,一致反映出GLDAS_ETp全年正偏,且误差的季节波动性较大,即在春季(3~4月)有大幅的正偏;MOD16_ETp和PM_ETp_2.1也全年正偏,夏半年(4~9月)的高估程度最严重;PM_ETp_2.0的偏差最小且季节波动最小,6~8月微有正偏,1,2月及9~12月则略有负偏。
从4套ETp数据与站点实测Epan的相关系数空间分布来看(图4),在贵州以及云南中、东部地区,各套数据与Epan均有较好的相关性,相关系数均大于0.8。而在海拔较高的云南西北部的站点,各套数据与实测值之间的相关性很低,大部分都低于0.4,这说明在海拔较高的云南地区,ETp数据与实测值的一致性较差。总体上相关关系表现出由东向西北递减的特征。
图4 不同ETp数据与Epan的相关系数r分布Fig.4 Correlation coefficient between Epan and ETp
由RMSE和MBE的空间分布可知(图5~6),PM_ETp_2.0总体在整个区域都存在负偏,负偏程度在云南的北部和南部站点较为明显,但总体精度较好,平均MBE为-6.1 mm/月;而GLDAS_ETp数据则整体高于Epan,其正偏程度由西部向东部降低。图6反映出位于云南西南部和贵州西部的部分站点正偏严重,MBE均大于50 mm/月。另一方面MOD16_ETp和PM_ETp_2.1的评估结果显示,两套数据在贵州地区比PM_ETp_2.0和GLDAS_ETp数据存在更大程度的高估,但总体上MOD16_ETp和PM_ETp_2.1比实际的Epan偏高了约48~49 mm/月,精度优于GLDAS_ETp。
MRAE的空间分布与RMSE、MBE的空间分布相似(图7),即绝大部分地区的MOD16_ETp、PM_ETp_2.1和GLDAS_ETp存在较大程度偏差,3套数据的MRAE均大于50%,在贵州北部的一些站点MOD16_ETp和PM_ETp_2.1的MRAE甚至达到了101%~130%;而全区域PM_ETp_2.0的MRAE则低于50%,呈现出均匀的空间分布状况。
综合4套数据的评估结合来看,PM_ETp_2.0与Epan偏差最小,数据质量最优,MOD16_ETp、PM_ETp_2.1及GLDAS_ETp在整个区域都存在高估的情况,尤其是GLDAS_ETp在春季严重高估,总体精度最差。
图5 不同ETp数据与Epan的均方根误差(RMSE)分布Fig.5 RMSE between Epan and ETp
图6 不同ETp数据与Epan的平均偏差(MBE)分布Fig.6 MBE between Epan and ETp
根据本文结果来看,基于GLDAS-2.0数据利用P-M方法计算得到的PM_ETp_2.0最为合理,而MOD16_ETp、GLDAS_ETp和PM_ETp_2.1产品计算存在高估的现象。下面从ETp计算方法及驱动数据两方面分析不同产品的数据质量差异原因。
从方法上看,MOD16_ETp、PM_ETp_2.0和PM_ETp_2.1均基于P-M方法,它们的逐月变化和误差的季节特征一致;而GLDAS_ETp是Noah模型模拟结果,存在较严重的高估问题。考虑到Noah模型以地表参数和气象驱动数据为输入数据,之后根据输入数据和一系列的参数化方案来计算地表温度、湿度、动量以及地表能量通量[16],因此用于计算ETa和ETp的一些参数是模型计算所得而非实际观测值,同时由于模型输入数据的误差,一些地区基于Noah模型所得的ETa会存在数值偏大的问题[17-18],同样也会造成ETp的高估问题。
图7 不同ETp数据与Epan的相对绝对误差(MRAE)分布Fig.7 MRAE between Epan and ETp
从驱动数据上看,ETp计算的驱动数据包括两类,一类是辐射数据,另一类是常规气象数据。而影响ETp数值大小的辐射驱动项主要是净辐射(Rn),该数据需由潜热通量减去显热通量和土壤热通量得到,故理论上净辐射应大于潜热通量(λE,水的汽化系数λ=2.45 MJ/kg)。GLDAS-2.0的Rn以及MOD16_ETp、PM_ETp_2.0所对应的潜热通量在2000~2014年的逐月变化过程见图8。图中显示,PM_ETp_2.0所对应的λE小于Rn,而MOD16_ETp的λE则明显高于Rn,这说明在使用MOD16和GLDAS-2.0数据计算ETp的过程中,前者所用净辐射数据明显大于后者,这是造成MOD16_ETp产品高估的主要原因。
而PM_ETp_2.0和PM_ETp_2.1存在较大差异的原因在于两个版本GLDAS的气象驱动数据存在偏差。如1.2.2节所述,GLDAS-2.0产品的全部气象驱动数据均为普林斯顿大学经过偏差校正的气象数据,而GLDAS-2.1的气象驱动产品则是不同来源的大气、降水、辐射数据。以2001年为例,对比GLDAS-2.0和GLDAS-2.1的气象驱动数据发现,GLDAS-2.0和GLDAS-2.1的气温、风速和比湿数据差别很小(相对误差约为0.3%,0.02%和0.4%)。两数据集的差异主要体现在,GLDAS-2.1与GLDAS-2.0相比短波辐射偏小(相对误差约-21.15%),长波辐射偏大(相对误差约21.86%),进而导致净辐射值偏小;气压偏大(相对误差约为0.14%),导致公式(1)的饱和水汽压差以及干湿表常数等参数发生变化,而这些差异最终造成了PM_ETp_2.1的计算结果偏大。相比而言,GLDAS-2.0的气象驱动数据精度较好,更适合于反映云贵地区ETp的变化。
图8 净辐射与潜热通量的逐月变化过程Fig.8 Monthly variation of net radiation and latent heat flux
本文基于2000~2014年云南和贵州地区42个站点的实测蒸发皿蒸发量数据Epan,对MODIS MOD16 与GLDAS-2.1潜在蒸散发数据产品(MOD16_ETp、GLDAS_ETp),以及结合P-M公式和GLDAS两个版本气象数据计算所得的潜在蒸散发数据(PM_ETp_2.0和PM_ETp_2.1),通过相关性和偏差等指标的变化,从空间分布、季节变化评估了各数据在云贵地区的质量,结果如下。
(1) 从全区面平均值来看,4套数据在贵州、云南中部和东部地区都与Epan有较好的相关性,逐月序列的总体相关性均大于0.9。4套数据的质量存在一定的季节性波动;GLDAS_ETp的季节波动性较大,总体与Epan的均方根误差约53.7 mm/月,在春季(3~5月)会有超过60%的大幅正偏;MOD16_ETp和PM_ETp_2.1全年正偏,且这种正偏情况在4~9月最明显,此时段内的MRAE均大于50%;PM_ETp_2.0的偏差最小且季节稳定性更好,仅在夏季(6~8月)有较小幅度正偏,而在1~2月以及9~12月略有负偏。
(2) 从空间分布来看,MOD16_ETp、GLDAS_ETp、PM_ETp_2.0和PM_ETp_2.1数据在贵州、云南中部和东部地区,与实测Epan都有较好的相关性,相关系数达0.6以上;而在海拔较高的云南西北部,4套ETp数据与Epan的一致性较差,相关系数仅为0.2~0.4。RMSE和MBE的空间分布表明,PM_ETp_2.0总体在整个区域存在较小程度的负偏,平均MBE为-6.1 mm/月;GLDAS_ETp则在云南西南部和贵州西部的部分站点严重正偏,MBE超过50 mm/月;而MOD16_ETp和PM_ETp_2.1虽然在贵州地区有较大程度正偏,但总体正偏量约为48~49 mm/月,精度优于GLDAS_ETp。
(3) 综合各类误差指标结果,PM_ETp_2.0与Epan最接近,数据质量最优,MOD16_ETp、PM_ETp_2.1及GLDAS_ETp在整个区域都存在高估的情况,尤其是GLDAS_ETp在春季严重高估,总体精度最差。