周再文,何贞铭,蒋松谕,相龙伟,彭 李
(1. 长江大学地球科学学院,湖北 武汉 430100;2. 荆州市自然资源卫星应用技术中心,湖北 荆州 434000)
湖泊不仅是重要的国土资源,且具有巨大的环境调节功能和效益,在调节江河径流、防洪抗旱、发展灌溉、保护生物多样性等方面发挥着重要的作用[1]。湖泊对流域内环境变化及人类活动十分敏感,是区域气候敏感的指示器和调节器[2]。湖泊萎缩会导致生物多样性减少、水体污染、调洪蓄水能力减弱、水体富营养化加重等一系列生态环境问题[3]。对湖泊流域进行长时序的变化监测,在提高湖泊动态变化规律认识,揭示气候变化与人类活动对湖泊生态系统的影响具有重要的意义。
随着遥感技术的快速发展,遥感技术被广泛用于水体识别和监测的研究。基于遥感提取水体的方法主要有单波段阈值法、多波段谱间关系法、水体指数与阈值法、支持向量机、随机森林及深度学习等方法。其中,水体指数法被研究人员广泛采用[4-5]。文献[6]使用NDWI削弱土壤植被等非水体因素影响,提取大型开阔水域,但易将城区建筑用地错分为水体;文献[7]在NDWI的基础上,将近红外波段替换为短波红外波段,提出MNDWI并解决水体提取难以消除阴影的难题;文献[8]通过AWEI放大水体与背景地物区之间的差异实现对水体的准确提取;文献[9]基于半干旱地区水系与背景噪音反射特点,提出了增强型水体指数(EWI),有效解决了半干涸河道与背景噪声混淆的问题。
由于传统遥感方法数据获取难、数据量大、影像预处理过程耗时费力,对大尺度长时序的遥感影像使用较少。融合了谷歌服务器高性能的计算能力及集成了数据庞大的卫星影像的GEE(Google Earth Engine)平台的方法,为高时空分辨率、大尺度、长时序水体变化分析提供了一种新的分析工具[10]。
洞庭湖流域是长江流域重要的调蓄湖泊,在维系长江中下游流域江湖关系和生态平衡具有不可或缺的重要作用,也是长江经济带城镇化高速发展的经济活跃地区[11]。目前关于洞庭湖流域水体变化研究多以短时序的年际变化为主,缺少季度尺度的长期水体时空演变特征及其驱动因素研究[12]。因此,本文基于GEE平台的多源时序影像和高性能的计算能力,细化洞庭湖流域水体时空演变特征。
洞庭湖位于长江中游荆江南岸,湖南省东北部,不仅是我国长江流域的第2大淡水湖泊,也是我国仅存的两大自由通江湖泊之一。洞庭湖地处亚热带季风气候区,春、夏冷暖气流交替频繁,夏、秋晴热少雨,秋寒偏旱。多年平均气温为16.4℃~17℃,年平均降水量为1100~1400 mm,无霜期为258~275 d[13]。
洞庭湖区大致可分为东洞庭湖、南洞庭湖和西洞庭湖3部分。为更好研究流域内三大湖区长时序水体变化特点,综合流域内降水、气温、人口、土地利用等多方面因素,本文以3大湖区的外接矩形范围作为研究区(如图1所示)。
1.2.1 多源时序影像数据
基于GEE平台获取自1989年11月1日至2022年10月30日覆盖研究区且云量少于20%的多源遥感影像共874景,将其作为数据源(如图2所示),影像数据均已经过辐射校正、几何校正和去云预处理。为探究流域内水体季节性变化的特点,将历年影像数据分丰水期(5—10月)和枯水期(11—次年4月)两期进行水体提取[12]。部分年份丰、枯水期影像受云层遮挡影响,使用Landsat 7影像替代。Sentinel-2影像提取细小部分水体信息更为明显,整体提取效果优于Landsat 8影像,且空间分辨率更高,因此2016年以后的影像采用Sentinel-2数据[14]。
图2 多源时序影像信息
1.2.2 气象数据
从GEE平台上获取1989—2019年ERA5 Monthly Aggregates月均气温数据和1989—2022年CHIRPS Daily逐天降水数据[15],利用月均气温和逐天降水数据分别计算出研究区历年丰、枯水期的平均气温和降水。
1.2.3 人口密度和土地利用数据
为探究人类活动对洞庭湖流域水体面积的影响,收集了WorldPoP全球人口数据,以及文献[16]基于GEE的335 709景Landsat影像制作的1990—2020年30 m中国逐年土地覆盖数据集(annual Chinaland cover dataset, CLCD)。
基于多源遥感数据提取洞庭湖流域湖泊面积技术路线如图3所示,包括5个步骤:①影像预处理;②最优水体指数选取;③提取丰、枯水期水体;④计算水体面积;⑤水体变化分析。
图3 基于多源遥感数据计算水体面积技术路线
选取当前引用最多的归一化差分水体指数(NDWI)、改进的归一化差分水体指数(MNDWI)、自动水体提取指数(AWEIsh)及基于线性判别分析的水体指数(WI2015)4种水体指数提取研究区水体[17]。以10年为周期,选择1991、2001、2011和2021年遥感影像作为测试数据,分别计算4种水体指数,选择研究区内水体提取最优水体指数。2011年枯水期因云层遮蔽,影像缺失,使用2010年枯水期影像。最优水体指数的选取主要步骤为:①生成年度丰、枯水期合成指数数据。选择4种指数对应的波段编号,逐像素计算4种指数,取每个像元的指数中值,生成年度合成指数数据。②计算分割阈值和提取水体。采用大津算法[18](OTSU)提取水体和非水体的分割阈值,将水体指数大于阈值的栅格单元定义为水体,否则为非水体。使用StratifiedSample函数对分割后的结果进行分层抽样,提取训练样本,利用随机森林分类器对研究区进行监督分类,将研究区划分为水体和非水体。③精度评价。在研究区按丰、枯水期分别随机生成200个精度验证点(4年采取相同的样本点),按照每年分类结果赋值,对分类结果进行人工目视判读,构造混淆矩阵,再计算水体分类的总体精度(overall accuracy, OA)、F1和Kappa系数。④选取最优水体指数。根据试验结果,分别计算4年水体指数以提取OA、F1和Kappa系数的平均值(见表1),选择最优水体指数。
表1 水体指数精度统计
经计算,NDWI、MNDWI、AWEIsh和WI2015的OA均值在丰水期分别为0.971 4、0.977 5、0.967 5、0.956 3;在枯水期分别为0.907 7、0.932 5、0.914 3、0.907 5。F1均值在丰水期分别为0.971 8、0.976 8、0.967 6、0.956 4;在枯水期分别为0.878 5、0.885 1、0.868 3、0.861 2。Kappa系数均值在丰水期分别为0.944 8、0.954 8、0.934 9、0.912 5;在枯水期分别为0.825 0,0.837 5,0.806 1,0.793 7。MNDWI指数的3种评价指标在4年丰、枯水期试验中均表现良好,因此选择MNDWI提取洞庭湖流域水体。
2.2.1 水体面积萎缩比例
因1989、1990年历史影像缺失,选取1991年的丰、枯水期水体面积作为近30年研究区水体参考面积,以丰、枯水期季度为单位,计算洞庭湖流域水体面积萎缩比例[19],分析水体面积动态变化。计算公式为
(1)
式中,Ac为待评估年份的研究区水体面积;Ar为历史同期参考水体面积。
2.2.2 标准差椭圆方法
标准差椭圆(standard deviational ellipse, SDE)是空间统计方法中能够精确揭示要素空间分布多方面特征的方法。SDE通过椭圆重心、方位角、长短轴等参数对地理要素的空间分布及时空演化过程进行定量的描述[20]。
洞庭湖区季节性特征明显,其中,枯水期湖区面积较小,平均面积为1 524.40 km2;丰水期相对较大,平均面积为2 236.03 km2。湖区丰、枯水期水体面积总体表现为下降趋势。
枯水期内研究区在2002年湖面积最大,为2 184.42 km2,在1992年湖面积最小,为1 073.57 km2,分别比枯水期平均值高43%和低30%;丰水期内1996年湖面积最大,为3 459.70 km2,2022年研究区湖面积最小,为1 005.03 km2,分别比丰水期平均值高48%和低57%。湖体面积在丰、枯水期均表现出波动明显特征。1992—1999、2004—2018年枯水期水体面积持续减小,1992、1998年萎缩比例达0.36、0.31,面积较1991年同期减少570.98和489.50 km2;湖区面积在丰水期同样存在一定程度的缩减,平均减少面积93.27 km2,其中2022年丰水期萎缩比例达0.57,面积缩减1427.27 km2(如图4所示)。
图4 1989—2022年洞庭湖流域水体变化
在图4中,选取丰、枯水期ASR大于0.2的各5个年份,即枯水期:1992、1998、2007、2010、2013年;丰水期:2008、2011、2014、2015、2022年。利用GEE平台将这些年份提取的水体影像与1991年湖区季度水体影像叠加分析,突出显示典型时期湖区的空间形态变化特征。
图5为湖区典型丰、枯水期水体空间形态的时空差异(部分)。结果表明,东洞庭湖区西北、西南两条分叉型水槽,以及连接三大湖区的两条过渡型深水河槽为常年保存水体,洞庭湖流域的主体湖区为东洞庭湖。
图5 典型年份丰、枯水期水体空间形态变化(部分)
枯水期3大湖区的水体由四周向中心逐渐缩减,过渡河槽变窄,南洞庭湖区支离破碎,碎片化、板块化特征明显,呈现大面积河流浅滩交错区。在枯水期缩减最明显的是东洞庭湖区中心的三角洲浅滩、漉湖及南洞庭湖的横岭湖(如图5框内所示)。东洞庭湖区中心三角洲浅滩向四周扩张,湖区与过渡河槽交汇区缩减为狭长河道,漉湖面积缩减明显,与主湖区连接通道变窄。南洞庭湖的横岭湖大面积被河流浅滩覆盖(如图5框内所示)。丰水期水体面积缩减较严重的区域是漉湖、万子湖和大连湖,2022年丰水期东洞庭湖主湖区及漉湖、大连湖基本消失。
为更深刻研究洞庭湖流域水体空间形态演变规律,利用标准差椭圆分析流域内水体变化的空间全局特征。综合标准差椭圆参数(见表2—表3)及分布离散趋势信息(如图6所示)可以看出,枯水期湖区水体重心虽然在1992、1998及2007年向南移动,但大致方向依然是向西、北方向移动,长轴方向由东南-西北向转向东-西向,再转向东南-西北向最后又转向东-西向,趋势反映出枯水期南洞庭湖及东洞庭湖存在缩减,流域内东洞庭湖两条分叉型水槽成为主要湖区。椭圆面积和扁率变化均为先减小后增大最后减小,整体呈减小趋势,说明在枯水期湖区碎片化、板块化特征明显。在丰水期椭圆中心大致向西北方向移动,但1998—2008年,2015—2022年,中心向西南方向移动,长轴方向、椭圆面积和扁率均呈波动变化特征。
表2 典型年份枯水期水体标准差椭圆相关参数
图6 典型年份丰、枯水期水体空间变化趋势
3.3.1 气候变化
在自然条件下的水循环中,湖体面积受最直接的影响因素是大气降水和蒸发[21]。
大气降水作为湖泊水量的补给来源,直接影响湖泊水体的面积变化,降水量减少很容易导致湖泊水域面积减少甚至消失。气温则是决定蒸发的重要因素,温度越高蒸发越快,蒸发量的升高会导致湖泊水分减少[22]。
总体而言,湖区面积与丰、枯水期降水具有较高的一致变化趋势,呈正相关,相关系数分别为0.728和0.604,这验证了降水在一定程度上影响湖区面积的变化;湖区面积与丰、枯水期气温呈负相关,相关系数分别为-0.348和-0.408,这表明气温升高会加速湖泊表面水体蒸发,导致水体面积减少(如图7和图8所示)。在丰水期,2002、2012和2015年,湖区面积与降水和气温表现出同步变化,但趋势不明显,对比原始影像发现,2002年丰水期的Landsat 7影像因去云处理,导致一部分水体区域缺失,因此水体提取结果偏小。2012、2015年的5月月均降雨量分别为270.76和260.49 mm,高于历史同期降雨水平,同年丰水期其余月份降水量均值为122.13和116.06 mm,因5月的极端降水影响了丰水期的整体降雨量均值,导致面积与降水变化不一致。
图8 1989—2019年丰、枯水期水体面积与气温变化
在枯水期内,2015年湖区面积与降水变化不一致,分析枯水期降水数据为:2016年4月月均降水量251.89 mm,前5个月的月均降水量75.53 mm,4月的降水造成2015年枯水期降水均值与面积变化不同步异常。
3.3.2 人类活动影响
人口是社会发展中最活跃的因素,人口的空间密度在某种程度上反映了人类活动的强度,如农业生产、生态系统的改变等[23]。基于WorldPOP人口数据得到的湖区人口密度变化情况如图9所示,人口密度从362人/km2增长至394人/km2,表明近20年研究区内人口密度一直处于增长状态。
图9 2000—2020年洞庭湖流域人口密度变化
人类活动带来的最直接的影响是区域内土地利用变化,对研究区进行长时间的土地利用变化监测也能在一定程度上反映人类活动强度。以CLCD为参考,从中裁剪出研究区的土地利用情况,以10年为周期选取1990、2000、2010、2020年土地利用数据,进行土地利用转移矩阵的计算,提取湖区主要土地覆盖类型各时期的面积。
对比CLCD与本文提取的水体面积(如图10所示),发现两者总体趋势虽然均下降,但存在一些差异。这可能是因为CLCD以Landsat全系列影像作为数据源,同时本文基于历年丰、枯水期提取水体,并未按照年度变化提取。这也是本文方法在洞庭湖流域典型洪涝灾害年份(1996、1998、2002、2016、2020年)水体局部变化特征更明显的原因。
图10 不同方法下洞庭湖流域水体面积对比
在研究期内,水体主要转换为耕地,一部分转换为不透水面、裸地与草地。其中,湖区过度围垦占水体减少面积的95.3%、97.7%、95.1%,人工围垦、围湖造田成为湖区面积减少的主要影响因素(见表4)。对水体面积与耕地面积进行相关性分析,两者呈负相关,相关系数为-0.952。
表4 1990—2020年洞庭湖流域水体变化情况 km2
本文基于GEE平台对洞庭湖流域水体时空演变特征进行研究,并从气候变化和人类活动两个方面分析了水体变化的影响因素。结果表明:①1989—2022年,洞庭湖区丰、枯水期水体面积总体呈下降趋势,3大湖区均有不同程度的缩减,平均减少93.27和140.15 km2。湖区季节性特征明显,枯水期湖面支离破碎、板块化;丰水期3大湖区相互连通,水域面积相对较大。②降水和气温是影响水体面积变化的重要自然因素。其中,湖泊面积与季度降水量呈较强正相关,与季度气温呈负相关。③随着湖区社会经济发展,人口增长、围湖造田等人类活动导致的土地利用类型转移,也是流域水体面积减少的重要原因。