三江源区蒸散量的时空分布特征

2021-07-14甘海洪金晓媚张绪财朱晓倩

现代地质 2021年3期

甘海洪，金晓媚，张绪财，朱晓倩

(中国地质大学(北京)水资源与环境学院，北京 100083)

0 引言

蒸散由地表水分蒸发和植物蒸腾两部分组成，是大气-植被-土壤中能量交换的主要途径[1]。三江源地区位于我国青海省南部，其面积约为35.2万km2，是中国海拔最高的天然湿地和面积最大的自然保护区，同时也是世界高海拔地区生物多样性最集中的自然保护区。然而三江源地区自然环境十分恶劣，生态系统脆弱敏感；更由于近年来气候变暖，三江源地区出现冰川消融、局部草地退化等生态环境问题。蒸散作用是干旱-半干旱地区地下水的主要排泄方式，因此准确计算地下水蒸散量不但是地下水资源量评价的关键，还对生态环境保护具有重要意义[2]。

目前已有部分学者[3-5]对三江源区蒸散、植被等做过一些研究，但针对三江源地区蒸散量变化趋势及其影响因素的研究还不足，少有的研究也只是利用了传统方法估算区域蒸散量，缺乏长时间序列蒸散量的时间与空间变化研究，无论从时间尺度还是空间尺度均不能准确地反映三江源区蒸散量的分布特征。随着遥感技术的发展，推出了许多地表蒸散量的估算模型，扩展了传统方法在时空尺度上的应用范围，减小了传统方法以点代面引发的误差。1983年，Seguin和Itier[6]尝试从卫星热红外资料中获得地面辐射温度与空气温度的差值，进而估算大尺度区域蒸散量。21世纪初Su[7-8]提出地表能量平衡系统(SEBS)，其目的是将卫星对地观测的可见光、近红外和热红外波段资料与实测气象或大气模式输出数据进行结合，从而更连贯地估算大气湍流通量，进而得到地表相对蒸散量。2008年，金晓媚等[9]利用SEBS模型对张掖盆地的区域蒸散量进行了估算，并用水均衡法进行验证，结果吻合较好；2014年，郭仁宏等[10-12]估算了柴达木盆地的区域蒸散量，并分析了其影响因素；2017年，杨笑天[13]利用Landsat-8数据和SEBS模型估算了格尔木河流域的区域蒸散量；2018年，张晓玉[14]利用SEBS模型估算了艾比湖流域的区域蒸散量，并研究了其时空变化规律。

本文将采用连续序列的MODIS数据(Moderate-resolution Imaging Spectroradiometer)和GLDAS (Global Land Data Assimilation System)数据，利用SEBS模型计算三江源地区2001—2017年17年间的区域蒸散量，并分析蒸散量的时空分布特征、变化趋势以及影响因素，以期对三江源地区水资源合理开发利用提供参考。

1 研究区概况

三江源地区位于青藏高原中东部、青海省南部，平均海拔4 300 m左右(图1)，其地理位置为北纬31°39′—36°12′、东经89°45′—102°23′，是青藏高原腹地，是我国三大主要河流长江、黄河和澜沧江的源头汇水区。研究区内地势西高东低，西部的长江源区平均海拔最高；山地是主要地貌，区内山脉分布众多，地形较为复杂。

图1 研究区地理位置图Fig.1 Geographic map of the study area

三江源地区表现为典型的高原大陆性气候。冷热两季分明，干湿两季交替；年际温差较小，昼夜温差较大。冷季时，由于受到青藏高原冷气流影响，温度低，降雨少，较为干燥；热季时，受到印度洋暖季风影响，降雨比较充沛，具有明显雨热同期特征。但是，三江源地区的降水量在区域上明显分布不均，其东部较为湿润，降雨较多，多年年平均降水量为672.6 mm；西部地区为干旱-半干旱地区，降雨较少，多年年平均降雨量为375.84 mm。全区气温年际变化不大，多年月均气温变化范围为-10.84～10.85 ℃；昼夜温差大，平均昼夜温差可达20 ℃以上[15]。

由于三江源地区海拔较高，空气稀薄，太阳辐射强烈，年日照时间超过2 300 h[15]。三江源水系发达，植被分布状况不均，黄河源区东北部和澜沧江源区为山区，植被覆盖度最好；长江源区的中部地形平坦，分布有覆盖度较高的高原草甸；而长江源区的西部以及黄河源区的西部气候干旱，植被稀疏，有较大面积的荒漠分布。

2 数据来源及研究方法

2.1 数据来源

2.1.1 MODIS数据

本文选用的是2001—2017年间的MODIS数据，包含地表反射率(MOD09A1)、地表比辐射率与地表温度(MOD11A2)。MOD09A1和MOD11A2两种数据时间分辨率均为8 d，空间分辨率分别为500 m和1 000 m。为保证数据一致性，对两种数据进行重采样处理，得到统一空间分辨率500 m。本文选用MOD09A1和MOD11A2数据共3 128景，均源于美国NASA网站(https://reverb.echo.nasa.gov/reverb)。

2.1.2 GLDAS数据

GLDAS数据全称为全球陆面数据同化系统数据，选用的是2001—2017年GLDAS_NOAH025_M数据中比湿 (kg·kg-1)、风速 (m·s-1)、气温(K)和气压(Pa)、土壤(0～10 cm)湿度(kg·m-2)5个波段，空间分辨率为0.25°，时间分辨率为1个月，最后将空间分辨率调整为500 m。数据共204景，源于https://disc.gsfc.nasa.gov/datasets。

2.1.3 DEM数据

文中采用了空间分辨率为90 m的SRTM-DEM数据，包括三江源地区地表高程、坡度和坡向数据。为使研究数据一致，对DEM数据也进行重采样，将空间分辨率与GLDAS数据和MODIS数据统一。DEM数据源于中国科学院地理空间数据云网站(http://www.gscloud.cn)。

2.1.4 气象数据

为分析气象因素对蒸散量的影响，本文选取三江源地区内14个气象站的气温、降水等气象数据。数据源于中国气象数据网(http://data.cma.cn)的“中国地面气候资料月值数据集”。

2.2 研究方法

2.2.1 表面能量平衡系统(SEBS)

利用SEBS模型计算区域蒸散量主要分为4步：(1)对遥感图像进行处理，获得一系列地表物理参数(如反照率、比辐射率、植被覆盖度等)；(2)建立热传导粗糙度模型[16]；(3)运用总体相似理论[17]确定摩擦速度、显热通量和奥布霍夫稳定度；(4)利用地表能量平衡指数[18]计算蒸发比，从而得到区域日蒸散量。

在任一时刻，地表的净辐射在感热通量、潜热通量和地面热通量方面都处于平衡，地表能量平衡方程为:

Rn=G0+H+λE

(1)

式中:Rn为净辐射通量，W·m-2；G0为土壤热通量，W·m-2；H为感热通量，W·m-2；λ为水的汽化潜热，J·kg-1；E为实际水分蒸散总量，mm·d-1；λE为潜热通量，W·m-2。

蒸发比可利用最干和最湿两种极限状态下的感热通量特征得到:

(2)

式中:EF为相对蒸发比，是实际蒸散量与最大蒸散量之比；Hdry为干极限状态下的感热通量，W·m-2；Hwat为湿极限状态下的感热通量，W·m-2。

日蒸散量为：

(3)

2.2.2 Sen+Mann-Kendall趋势分析

Sen 趋势度分析是一种利用中位数函数研究长时间序列的方法，但因其本身不能实现序列趋势显著性判断，而 Mann-Kendall方法则无须样本数据服从特定的分布，也不受少数异常值的影响，因此引入该方法对序列趋势进行显著性检验，得到Sen+Mann-Kendall趋势分析方法[19-20]。

Sen趋势度β计算公式为:

(4)

趋势度β可判断时间序列X(x1，x2，x3，…)趋势的升降，当β> 0时，时间序列呈上升的趋势；β< 0 时，时间序列呈下降的趋势。

Mann-Kendal公式如下:

对于时间序列X(x1，x2，x3，…)，

(5)

其中:n为数据样本长度；sgn为符号函数。sgn定义如下:

(6)

当n<10时，直接使用统计量S进行双边趋势检验；当n≥10时，统计量大致服从正态分布，则使用检验统计量Z进行趋势检验。

Z值计算公式为:

(7)

同样，采用双边趋势检验，在给定显著性水平α下，当|Z|≤Z1-α/2时，接受原假设，即趋势不显著；当|Z|>Z1-α/2时，则拒绝原假设，即认为趋势显著。

3 结果分析

3.1 区域蒸散量分布特征

3.1.1 蒸散量时间分布特征

SEBS系统计算出的为日蒸散量:取计算区域内所有像元蒸散量值的平均值作为区域日蒸散量值；再对该月多幅图像日蒸散量值平均，平均日蒸散量值乘以相应月份的天数得到月蒸散量；最后将该年12个月蒸散量值相加即得到年蒸散量值。经过计算得到2001—2017年期间三江源地区逐年蒸散量，得到研究区蒸散量年际变化趋势图(图2)和研究区多年平均蒸散量年内变化趋势图(图3)。从图2和图3可以看出，长江源、黄河源、澜沧江源蒸散量的时间变化趋势是基本一致的。

图2 研究区蒸散量年际变化趋势Fig.2 Annual variation of actual ET in the study area

图3 研究区多年平均蒸散量年内变化趋势Fig.3 Variation of the annual average ET in the study area

2001—2017年，长江源、黄河源、澜沧江源年平均蒸散量分别为414.18 mm、424.67 mm、482.87 mm，整个三江源区年平均蒸散量为425.99 mm。17年间，三江源区年总蒸散量整体呈增长趋势；2001年至2009年蒸散量增加明显，增加速率约为23.4 mm/a；2009年至2017年蒸散量总体平稳，略有减小。三江源区年蒸散量最大值为552.55 mm(2009年)，最小值为275.03 mm(2002年)。长江源、黄河源、澜沧江源区年蒸散量最小值分别为238.04 mm(2002年)、287.98 mm(2002年)、275.03 mm(2002年)，最大值分别为523.78 mm(2009年)、613.61 mm(2009年)、597.76 mm(2014年)。从多年月平均蒸散量趋势图(图3)中可以看出，蒸散量主要集中在5月至9月，这5个月的蒸散量占全年80%以上；6月蒸散量最高，平均可达80 mm，1月和12月的蒸散量很小，接近0 mm。

3.1.2 蒸散量空间分布特征

根据2017年研究区年蒸散量空间分布图(图4)可知，研究区内湖泊水体的年蒸散量最高，可达900 mm以上；裸土地区蒸散量较小，年蒸散量基本小于100 mm。三江源三个源区多年年平均蒸散量从大到小依次为:澜沧江源(482.87 mm)>黄河源(424.67 mm)>长江源(414.18 mm)。澜沧江源蒸散量明显大于黄河源和长江源，黄河源和长江源蒸散量相差不大。研究区内黄河源区的东北部主要为山区，由于山区的蒸散量误差较大，因此本文将山区部分剔除。对于地势相对平坦的地区，长江源西部荒漠地带蒸散量较小，为30～100 mm，但湖泊水体的蒸散量大，年蒸散量最高可达1 100 mm左右；中部地区有高山草甸分布，年蒸散量为350～500 mm；东部地区植被较好，有高原草甸分布，年蒸散量较大，为500～850 mm。黄河源西部分布有扎陵湖、鄂陵湖，湖泊水体年蒸散量最高可达1 200 mm；西部其余地区年蒸散量为350～550 mm；南部地区分布有高山草原和草甸，年蒸散量为550～850 mm。澜沧江源有覆盖度较高的高山草甸分布，年蒸散量较大，为500～800 mm。

图4 2017年研究区年蒸散量空间分布图Fig.4 Spatial distribution of annual ET in the study area in 2017

3.2 蒸散量趋势分析

为分析研究区蒸散量的连续空间变化趋势及其变化的显著性，本文利用Sen+Mann-Kendall显著性检测法对研究区17 a(2001—2017年)的蒸散量叠加数据进行分析。首先，使用Sen 趋势度分析方法得到研究区蒸散量空间变化趋势图(图5)；其次，根据Sen 趋势度分析结果和研究区实际情况，取0≤β≤ 6时，蒸散量呈基本不变趋势；再根据趋势度β来判断时间序列趋势的升降，检验统计量Z来进行趋势检验；最后，利用变化趋势分类表(表1)将变化趋势分为5类，分别为显著减小、轻微减小、基本不变、轻微增长、显著增长，得到研究区蒸散量显著性变化分类图(图6)。

表1 变化趋势分类

图5 研究区蒸散量空间变化趋势图Fig.5 Spatial ET variation trend in the study area

由研究区蒸散量显著性变化分类图(图6)可知，在2001—2017年期间，三江源地区绝大部分地区蒸散量变化呈增长趋势，显著增长地区面积占比达62.62%，轻微增长地区面积占比为28.03%，两者面积占比之和超过90%；基本不变地区面积占比为8.99%，轻微减小地区面积占比为0.37%，显著减小地区面积占比远小于0.01%，可忽略不计。由此可知，三江源区绝大多数地区蒸散量都在增加。

图6 研究区蒸散量显著性变化分类图Fig.6 Classification of significant ET change in the study area

从研究区蒸散量变化趋势面积统计结果(表2)来看，显著增长趋势表现为澜沧江源(68.37%)>长江源(63.25%)>黄河源(58.73%)；轻微增长趋势表现为黄河源(32.69%)>长江源(28.21%)>澜沧江源(23.37%)；基本稳定趋势表现为黄河源(8.32%)>长江源(8.22%)>澜沧江源(7.63%)；轻微减小趋势表现为澜沧江源(0.63%)>长江源(0.31%)>黄河源(0.26%)。三个源区中，轻微减小地区面积占比均小于1%，显著减小地区面积占比均远小于0.01%。

表2 研究区蒸散量变化趋势面积统计

4 蒸散量影响因素分析

4.1 气象因素

本次研究选用五道梁、兴海、托托河、杂多、曲麻莱、玉树、班玛、清水河、玛沁、达日、河南、久治、囊谦、班玛共14 个气象站的地面月值数据。分别对2001—2017年长江源、黄河源、澜沧江源和三江源区内的多个气象站数值进行平均后分析。

4.1.1 降水

通过对气象站降水数据分析，得到研究区降水量年际变化趋势图，如图7所示；并对全区气象站降水数值进行平均求取三江源区2001—2017年期间的月均降水值，利用每月的降水值与对应月的蒸散量进行相关性分析，得到降水对蒸散量影响的分析图，如图8(a)所示。

图7 研究区降水量年际变化趋势Fig.7 Annual variation of precipitation in the study area

图8 降水(a)和气温(b)对蒸散量的影响Fig. 8 Effect of precipitation (a) and air temperature (b) on evapotranspiration

结果表明:三江源区降水量最大值为593.63 mm(2009年)，最小值为419.14 mm(2002年)，与三江源区年蒸散量最大值和最小值年份相对应。多年年平均降水量:澜沧江源(548 mm)>黄河源(545 mm)>长江源(468 mm)，与蒸散量空间分布特征一致。研究区蒸散量与降水的相关性很好，确定系数达0.89。蒸散量总体与降水量变化趋势一致，呈正相关关系，即蒸散量随降水量的增大而升高。

4.1.2 气温

通过对气象站气温数据分析，得到研究区年均气温和月均气温变化趋势图(图9)，并对全区气象站气温数值进行平均求取三江源区2001—2017年期间的月均气温值，利用每月的气温值与对应月的蒸散量进行相关性分析，得到气温对蒸散量影响的分析图，如图8(b)所示。

图9 研究区年均气温 (a)、月均气温 (b)变化趋势Fig.9 Annual (a) and monthly (b) variation of air temperature in the study area

结果表明:除2008年较异常之外，2001年至2009年三江源区年平均气温明显增加，2009年至2017年气温相对平稳，与蒸散量变化趋势基本一致。三江源三个源区多年月平均气温从大到小依次为:澜沧江源>黄河源>长江源，与蒸散量空间分布特征一致。长江源、黄河源、澜沧江源多年月平均气温最大值分别为9.9 ℃、11.3 ℃、13.2 ℃，最小值分别为-11.7 ℃、-11.3 ℃、-6.7 ℃。当气温低于零下10 ℃时，蒸散量基本上接近0 mm；当温度处于-10～-5 ℃之间时，蒸散量基本上小于10 mm。研究区蒸散量与气温具有很好的正相关关系，确定系数达0.80。因三江源区蒸散量在2001—2009年期间增加相对明显，而在2009—2017年期间相对平稳，这与气温的变化较为一致，故相对于降水，蒸散量对气温的响应更明显。

4.2 植被指数

通过归一化植被指数(NDVI)可以了解植物的生长状态、植被覆盖程度。NDVI的取值范围为-1～1，<0为水体。NDVI 值越大，说明植被长势越好，植被覆盖率越高；NDVI 越小，植被长势越差，植被覆盖率越低。本文将2001—2017年的连续序列共204个植被指数MODIS NDVI数据与相对应的月蒸散量进行相关性分析，得到植被指数对蒸散量影响分析图(图10)。结果表明:两者吻合较好，蒸散量波动趋势与植被指数NDVI变化趋势趋于一致；两者相关性较好，确定系数为0.62，呈正相关。

图10 植被指数对蒸散量的影响Fig.10 Effect of NDVI on evapotranspiration

4.3 用地类型

用地类型不同，蒸散量的差异较大。研究区主要用地类型较为简单，主要分为灌木丛、草地、裸土、水体4类。为了分析对比4种用地类型的蒸散量状况，求取每种用地类型年蒸散量平均值，对三江源地区2017年年蒸散量进行统计得到研究区不同用地类型蒸散量柱状图(图11)。结果表明:水体的蒸散量最大，约为860 mm；灌木丛地区蒸散量其次，约为610 mm；草地地区蒸散量约为415 mm；裸土地区蒸散量最小，约为65 mm。

图11 研究区不同用地类型蒸散量柱状图 Fig.11 Histogram of ET of different land-use types in the study area

4.4 土壤湿度

土壤湿度对蒸散量也具有一定影响。由于缺乏土壤湿度实测数据，本文采用GLDAS土壤(0～10 cm)湿度数据，分析土壤湿度对蒸散量的影响。将年内12组GLDAS土壤(0～10 cm)湿度数据进行平均，得到该年平均土壤湿度数据，绘制出研究区土壤湿度与蒸散量年际变化趋势图(图12)。对2001—2017年共204个GLDAS土壤(0～10 cm)湿度数据与计算得出的逐月蒸散量进行相关性分析，得到土壤湿度对蒸散量影响分析图(图13)。结果表明:2001—2017年，蒸散量的变化趋势与土壤湿度的变化趋势基本一致，两者具有一定的相关性，确定系数为0.60。