青藏高原是维持中国乃至东亚地区生态系统的重要水塔,被称为世界的第三极,是全球变化响应敏感的地区(Solomon, 2007),也是研究全球变化的关键地区和热点(Freymueller, 2011; Wischnewski et al, 2011)。青藏高原上高山湖泊总面积超过44990 km2(姜加虎等, 2004; Ünal et al, 2011),约占中国湖泊总面积的49.5%。湖泊作为青藏高原最重要的地表水赋存形态,对气候的波动变化极为敏感, 是揭示全球气候变化与区域响应的重要信息载体(王苏民, 1998; Liu et al, 2009; 郭超等, 2014),并影响生态系统和人类活动的许多方面(Steinman, 2005)。

湖泊水量收支变化是青藏高原地表过程与环境变化的最重要组成部分之一,承载了大量的降水、蒸发、温度和湿度等气候变化的信息,反映着高原气候干湿的变化(鲁安新等, 2006),是认识流域水量响应气候变化尺度和强度的重要参照系(于革等, 2004)。青藏高原湖泊水量变化的动态监测,有助于理解青藏高原水文资料缺乏地区水循环过程,揭示高原水资源的动态变化及其对气候变化的响应特征,可为生态屏障的保护提供科学依据。然而,青藏高原地区气候环境条件恶劣、人迹罕至,常规的水文观测不足。遥感卫星技术的发展,为青藏高原地区湖泊变化监测提供了更多的数据支持,有效地弥补了高寒地区观测资料的缺乏。已有很多学者对青藏高原地区湖泊面积、水位的动态变化开展了研究(Wu et al, 2008; Zhang, Xie, Duan, et al, 2011; Zhang, Xie, Kang, et al, 2011; 李均力等,2011; Kropáček et al, 2012; Phan et al, 2012; 张鑫等,2014; 车向红等,2015; 张国庆, 2018)。就青海湖而言,当前基于遥感的研究也大多集中在湖泊面积或水位的变化方面(Zhang , Xie , Duan , et al, 2011; 李晓东等, 2012; Wang et al, 2013; 李燕等,2014; 骆成凤等, 2017; 祁苗苗等, 2018),而将卫星测高数据和卫星图像数据联合研究湖泊水量变化的工作还较少(张国庆等, 2013)。本文在前人研究的基础上,基于MODIS产品影像数据和卫星高度计数据构建青海湖的水位—面积关系,进而定量估算青海湖水量的变化,对湖泊流域的水量平衡模拟提供有效的验证,并初步探讨其与气候因子的关系,为该地区的水资源可持续发展及生态环境保护提供科学依据。

本文对青藏高原上面积最大的湖泊——青海湖(中国第一大咸水湖)进行研究。青海湖位于青藏高原东北隅(36°32'N~37°15'N, 99°36'E~100°47'E),祁连山东南部的内陆封闭盆地最洼处,地势从西北向东南倾斜。集水面积29661 km2,补给系数5.83(图1)。湖体呈椭圆形,入湖河流40余条,西、北部较多且水量较大(王苏民, 1998)。湖区处于中国东部季风区、西北干旱区和西南干旱区的交汇地带,干寒、少雨、多风、太阳辐射强、气温日差较大(中国科学院兰州分院, 1994),是典型的高寒半干旱气候。降水集中在夏季,5-9月的降水占全年的80%(Li, 2007)。湖水主要依赖地表径流和湖面降水补给(刘佳等, 2009)。

显示原图| 下载原图ZIP| 生成PPT

图1   青海湖地理位置及其水系特征

Fig.1   Geographical location of the Qinghai Lake and its drainage system

据《中国湖泊志》(王苏民, 1998)记载,1981年青海湖实测水位3193.92 m、面积4340 km2,与1961年测量结果相比,面积减少295 km2。岳辉等(2017)在其研究中也提到,近40年来青海湖面积整体呈现先缩小后趋于平缓而后增大的趋势。尤其是在2000-2010年间,湖泊水位面积呈逐年增加的趋势(李晓东等, 2012),湖泊边界变化区域集中在湖区东岸和西岸(刘瑞霞等, 2008; 刘宝康等, 2013)。但李林等(2011)预测了在2010-2020年间,青海湖水位总体仍可能以下降为主。

本文在提取湖泊水体面积时使用MODIS地表反射率8天合成产品MOD09A1,空间分辨率为500 m,时间范围为2001-2016年,青海湖涉及两幅MODIS影像,故解译时需收集影像1468景(2001年影像较其他年份少2景影像)。另外,在提取的过程中还需要使用MODIS土地覆盖类型产品MCD12Q1,其空间分辨率为500 m,每年需1景数据。湖泊水位数据采用的是法国LEGOS (Laboratoire d' Etudes en Géophysique et Océanographie Spatiales)实验室卫星产品数据,该数据整合了Topex/Poseidon、ERS-1/ERS-2、ENVISAT、JASON-1及GFO一系列卫星雷达高度计数据。最终的动态水位数据,是以Topex/Poseidon数据为基础将不同卫星数据成果融合而成的(卢善龙等, 2010)。

气象数据来自中国科学院青藏高原研究所开发的《中国区域高时空分辨率地面气象要素驱动数据集》(CMFD)和《中国地面气候资料日值数据集(V3.0)》。CMFD是一种接近地表的气象和环境要素再分析数据集,空间分辨率为0.1°,时间分辨率是3 h。《中国地面气候资料日值数据集(V3.0)》包含了中国824个基准、基本气象站自1951年1月以来各站的气压、气温、降水量、蒸发量、相对湿度、风向风速、日照时数和0 cm地温要素的日值数据,本文中所使用的数据为刚察气象站((37°20' N, 100°08' E)的气温、蒸发、降水和风速数据。

基于上述MODIS和LEGOS遥感数据分别提取青海湖面积、水位长序列信息,然后通过概化湖泊形态构建青海湖水位—面积关系,进而估算湖泊水量的变化。

3.2.1 湖泊面积与水位信息的提取

水体面积提取的关键是水陆分割的阈值选择,阈值选择的主观性可影响水体面积提取的精度。针对这一问题,本文采用矢量边界辅助MODIS数据的水体自动提取方法。该方法将实现图像分割阈值的自动化选取,对不同的图像给予不同的阈值,并且给定谷值出现的阈值区间,在该区间内寻找直方图最小值即为图像分割的阈值。这种方法比单一阈值法的精度有了大幅提高,不需要人工干预即可获得适合的阈值。根据这一思路,基于IDL平台实现了相关算法,自动提取了青海湖不同时间点的边界信息用以估算湖泊水面面积。另外,为消除云覆盖的影响,提高提取精度,基于矢量边界,将初始边界向内腐蚀,将腐蚀后边界内的范围划分为水体,得到湖泊水体面积数据(图2)。LEGOS实验室卫星产品数据是完整的水位产品数据,其水位数据为一颗卫星或几颗卫星连续几次测量数据的平均值整合而成。因此,湖泊水位数据在提取出研究区的相关数据后只需对时间格式进行简单的转化。

显示原图| 下载原图ZIP| 生成PPT

图2   湖泊面积提取步骤

Fig.2   Flowchart of the lake area extraction process

3.2.2 基于水位—面积关系的湖泊水量平衡估算

湖泊水量平衡(水量净收支)变化最终将表现为湖泊水位和面积的相应变化。当水量净收支为正直时,水位上升,湖泊水面面积扩大;反之相反。基于此,当湖泊水位发生微小变化时,其蓄水量的微小变化可表示为：

dV=AdL(1)

即：

V=∫L1L2AdL(2)

式中：V为水位从 L1变化到 L2时,湖泊蓄水量的变化量,反映了这一过程湖泊的水量净收支。显然,对式(2)进行求解首先必须知道湖泊面积A与水位L的函数关系,即湖泊的水位—面积关系。

由于青海湖湖面近似呈椭圆形,假设湖体可概化为倒椭圆台(图3),则对任一横截面的椭圆面 积为：

A=παβ=kπα2(3)

显示原图| 下载原图ZIP| 生成PPT

图3   湖泊形态概化

Fig.3   Hypothetical lake morphology

式中： A代表面积,α和β分别为椭圆的长半轴和短半轴,k(=β/α)近似为湖泊的宽长比。当湖泊水位由 L1变化到 L2,则湖泊面积由 A1变化到 A2,且有：

A2=kπα22=kπ(α1+∆L∙ctgθ)2(4)

式中： θ为湖岸倾角, ∆L=L2-L1为水位变化量。由式(4)可以得出,在湖泊宽长比与湖岸倾角确定的情况下,湖泊面积的变化量与湖泊水位的变化量的平方成正比。基于此,为简化计算,进一步假设湖泊水位—面积关系可以用二次方程加以拟合,即：

A=aL2+bL+c(5)

将湖泊的水位—面积关系式(5)代入式(2),则湖泊蓄水量的变化量可表示为：

V=∫L1L2(aL2+bL+c)dL=13aL3+12bL2+cL|L2L1=13L2-L1A1+A2+aL1L2+12bL2+L1+c(6)

必须指出,在湖泊面积—水位关系的构建中,针对不同的湖泊形态、数据与研究目的,不同的研究可有不同的方法。比如,Song等(2013)和Zhu等(2014)在研究中均通过构建面积水位的线性关系来估算水量变化;而Duan等(2013)在估算米德湖(美国)及塔纳湖(埃塞俄比亚)水量变化时,分别构建了2个湖泊的面积水位二次拟合方程来估算湖泊水量变化。此外,如果忽略式(5)中的一次项和常数项,即湖泊水位—面积关系表示为 A=KL2(K为系数),则湖泊蓄水量变化的计算公式可简化为：

$V = \frac{1}{3}(L_{2} - L_{1})[A_{1} + A_{2} + \sqrt{A_1A_2}] \quad (7)$

式(7)即为Taube(2000)提出的计算公式,也是张国庆等(2018)采用的计算方法。根据推导过程,式(6)更符合湖泊的形态,估算湖泊水量变化精度也会更好,因此,本文采用式(6)来计算湖泊水量平衡的变化。

湖岸边水与陆地的混合像元影响着水体面积提取结果,因此不同分辨率的影像对于湖岸线的提取的效果不同。在分辨率越高的情况下,混合像元所占整体水体像元比就越小,提取的湖泊面积就约接近于实际值。本文将利用相同或相近日期的上述数据、算法提取的面积结果与分辨率为30 m的Landsat遥感数据提取结果进行比较,其中利用Landsat数据提取时应用的方法为目视解译,比较结果见表1。

表1   青海湖面积、水位数据精度验证

Tab.1   Verification of the accuracy of the Qinghai Lake area and water level data

新窗口打开

从表1中可以看出,本文与利用分辨率较高的Landsat数据提取的湖泊面积结果相比差距不大,RMSE=11.32 km2,仅为青海湖面积的0.27%。另外与相同时期提取结果(骆成凤等, 2012)相比,2016年8月份的面积相差34.68 km2,相对误差为0.79%。

卫星测高数据由于与实测水文站采用的基准面不一样,所以在湖泊高程的提取结果上存在一定差异。赵云等(2017)的研究认为,高度计获得的水位数据要减去约40 cm的垂直偏差(焦文海等, 2002; 郭海荣等, 2004)。表1将水位结果与文献(Zhu et al. 2014)中所用实测水文站水位数据进行了比较,其中最大的水位差为0.19 m,而最小的水位差仅为0.05 m。

基于MODIS和LEGOS等遥感信息,应用前述算法,得出青海湖2001-2016年面积与水位的变化(图4)。由图4a可以看出,2001-2014年青海湖水位整体呈现上升态势,变化速率为0.10 m/a; 同一时期,青海湖湖泊水面面积扩张了187.93 km2,变化速率为11.56 km2/a(图4b)。进一步分析可以发现,2004年为湖泊水位和面积变化的转折年。2001-2004年,湖泊年均水位由3193.77 m下降到3193.64 m,降幅达13 cm。2004-2014年,湖泊水位持续上升,升幅达1.28 m。相应地,湖泊水面面积由2001年的4135.63 km2缩小到2004年的4122.76 km2。2004-2016年,除个别年份有所波动外,湖泊年均面积整体持续扩张。本文的研究结果与既有的相关研究有较好的一致性。例如：骆成凤等(2017)的基于Landsat的遥感数据研究也发现2004年是青海湖面积变化的转折年;张国庆等(2013)基于ICESat卫星高度计数据的研究则发现,青海湖年均水位从2003-2009年上升了0.41 m,与本文同一时段内的水位变化量(0.44 m)基本相当。

显示原图| 下载原图ZIP| 生成PPT

图4   青海湖面积水位拟合情况

Fig.4   Relationship between water levels and surface areas of the Qinghai Lake

由于寒冷季节湖泊面积的提取会受到结冰情况的影响,因此在分析湖泊水位—面积关系时仅使用6-10月份的面积、水位数据。图4c给出了青海湖水位—面积的关系,其拟合方程为：

A=2.1267L2-13501.5L+21431919.1≈2.1267(L-3174.515)2,L>3193m(8)

式中： A代表的是湖泊面积月均值, L代表的是湖泊水位月均值。决定系数 R2=0.83( p 500 km), Yamzhog Yumco had seriously shrunk, and it has continued to shrink in recent years; Qinghai Lake had shrunk during the period, but some new researches indicated that it has been expanding since the year 2004; Siling Co, Nam Co, and Chibuzhang Co had expanded in the period. We divided the newborn lakes into six categories according to their forming reasons, including river expansion, wetland conversion, etc. The changes of natural conditions led to the death of four lakes, and human exploitation was the main reason for the death of Dalianhai Lake in Qinghai. We picked out three regions which were sensitive to the change of climate and ecological environment: Nagqu Region, Kekexili Region, and the source area of the Yellow River (SAYR). Lakes in both Nagqu and Kekexili have been expanded; meanwhile, most lakes in the SAYR have obviously been shrunk. These regional patterns of lake changes were highly related to variations of temperature, glacier, precipitation, and evaporation. Our investigation and analysis will provide references for researches related to lake changes in the QTP and the response to climate fluctuations. [44] Wang X W, Gong P, Zhao Y Y, et al.2013.

Water-level changes in China's large lakes determined from ICESat/GLAS data

https://doi.org/10.1016/j.rse.2013.01.005      URL      [本文引用: 1]      摘要

Terrestrial and aquatic responses to climate change and human impact on the southeastern Tibetan Plateau during the past two centuries

https://doi.org/10.1111/j.1365-2486.2011.02474.x      URL      [本文引用: 1]      摘要

The response of lake-glacier variations to climate change in Nam Co Catchment, central Tibetan Plateau, during 1970-2000

https://doi.org/10.1007/s11442-008-0177-3      URL      [本文引用: 1]      摘要

Water level variation of Lake Qinghai from satellite and in situ measurements under climate change

https://doi.org/10.1117/1.3601363      URL      [本文引用: 1]      摘要

Monitoring lake level changes on the Tibetan Plateau using ICESat altimetry data (2003-2009)

https://doi.org/10.1016/j.rse.2011.03.005      URL      [本文引用: 2]      摘要

Trend of salt lake changes in the background of global warming and tactics for adaptation to the changes

https://doi.org/10.1111/j.1755-6724.2004.tb00198.x      URL      摘要

Monitoring the fluctuation of Lake Qinghai using multi-source remote sensing data

https://doi.org/10.3390/rs61110457      URL      [本文引用: 3]      摘要