基于遥感的南四湖菹草群落时空演变特征及其原因分析

焦银合,于泉洲,2,刘恩峰,梁春玲,田晓飞,2,张保华,2

(1.聊城大学 环境与规划学院,山东 聊城 252059； 2.聊城大学 东平湖湿地研究所,山东 聊城 252059； 3.山东师范大学 地理与环境学院,济南 250014；4.商丘师范学院测绘与规划学院,河南 商丘 476000)

湖泊湿地具有独特的生态功能和复杂的生态系统结构,对维护区域生态平衡和支撑区域经济发展具有重要作用[1]。水生植被对草型湖泊的水质改善及生态系统稳定都起着不可忽视的调节作用。菹草(PotamogetoncrispusL.)是一种常见的多年生沉水单子叶植物。它多生长于池塘、湖泊湿地以及河流中,每年4—5月迅速繁殖露出水面,6月中上旬腐烂衰退。湖泊中含有适量的菹草对水质具有一定的净化作用[2],但大量繁殖不仅会引起水体酸碱性变化,更会对水质和湖泊生态环境造成破坏。因此,识别和监测菹草的时空变化特征对于湖泊生态系统的保护和管理具有重要意义。

传统的基于实地采样调查监测菹草的方法往往费时费力。近年来,随着遥感技术的发展,使得湖泊湿地植被种群识别监测的效率大大提高了。当前,国内外许多学者对于湖泊沉水植物的监测识别多是基于高分或高光谱数据,如,李俊生等[3]基于高光谱数据提取分析沉水植物的光谱特征差异,并对太湖的菹草进行了识别；刘克等[4]基于湿地植被的光谱特征,利用OMIS成像光谱数据提取了9个最佳子波段,很好地区分了浮水、湿生、沉水和挺水植被；袁琳等[5]借助高分辨率Quickbird反射率与实测的苦草盖度数据,经过回归建模反演了苦草盖度的分布特征；Matta等[6]利用高光谱MIVIS航拍图像,根据不同植物光谱差异绘制了加尔达湖内的水生植被分布图；Ozge等[7]利用Quickbird卫星数据识别了土耳其安纳托利亚地区的大型浅水湖Mogan中的水生植物；Underwood等[8]通过对HyMap影像的光谱混合分析将巴西水草和水葫芦进行分类提取,精度达到93%和73%。

当前,基于遥感数据光谱信息的沉水植物识别(特别是菹草的识别),往往是基于单一时相的数据进行的。若结合多时相中分辨率成像光谱仪(MODIS)影像的归一化植被指数(NDVI),建立典型水生植被的时序特征曲线,则能更好地指导Landsat数据对菹草的精确识别。目前,基于MODIS-NDVI时序曲线的地物识别研究主要应用于林业、农业[9-10]领域,在浅水湖泊沉水植被识别方面的研究还较少。譬如,Brozen等[11]借助MODIS-NDVI对美国第二大淡水湖奥基乔比湖中的入侵水生植被进行识别监测；张圳等[12]利用MODIS-NDVI数据,分析了太湖2000—2015年太湖水生植被的季节变化和年际变化。而借助NDVI时序对菹草的遥感识别仅见于梁莉莉等[13]对东平湖的研究,其研究方法也有待进一步完善和明确。

本研究利用MOD13Q1产品数据提取湿地不同植被群落的NDVI时序特征,并将其作为物候标志,分析菹草与其他植被的物候差异,然后利用这种差异并结合多期不同时相的Landsat 影像,提取了南四湖近30年来菹草群落的时空格局特征,并分析了原因。本研究在传统光谱识别方法的基础上,融合了沉水植被的物候信息,并借鉴了决策树方法,实现沉水植被的遥感提取。一方面为菹草这一特殊沉水植被的遥感识别提供了一种行之有效的方法；另一方面,也通过揭示菹草群落在南四湖的时空变化特征及驱动因素,为防控其生态风险、保障南水北调东线工程水质安全提供一定的科学支撑。

1 研究区概况和数据来源

1.1 研究区概况

南四湖位于山东省西南部,地理坐标为34°20′～35°20′N,116°30′～117°28′E,是华北地区最大的淡水湖泊湿地。该湖主要由上级湖(由南阳湖、昭阳湖和独山湖组成)和下级湖(微山湖)构成,湖区由西北向东南狭长延伸,南北长约120km,东西宽5～25km,湖泊北高南低,全湖面积约1 266km2。南四湖影像如图1所示。

图1 南四湖2018年5月Landsat 8 OLI标准假彩色影像Fig.1 Landsat 8 OLI standard false color image of Nansi Lake in May 2018

南四湖属河迹洼地型湖泊,湖底浅平,平均水深1.5m[14]。南四湖水质肥沃,水草繁茂,分布着大量的芦苇荡和荷田[15],是典型的浅水草型湖泊。近年来,每到春季,南四湖菹草便爆发性增长,泛滥湖区,尤其自微山县2012年实施了压缩湖区网箱养殖政策之后,敞水面的释放更为菹草泛滥提供了良好的水域环境,使得菹草在每年5月初都成为南四湖沉水植物区的唯一优势种。菹草在6月中旬迅速腐烂导致水质污染严重,对湖泊湿地水环境和南水北调东线输水水质保护造成较大的压力[13]。

1.2 数据来源

本研究以南四湖地区MODIS及Landsat系列数据为基础。其中MODIS数据采用250m分辨率MOD13Q1(MODIS/Terra Vegetation Indices 16-Day L3 Global 250m SIN Grid)陆地产品,其时间分辨率为16d。由于MOD13Q1时间分辨率高,能满足菹草生长期内的监测要求,故可用年内23景MODIS数据建立湖区典型地物的物候特征曲线[16]。所用的Landsat数据包括Landsat 5 TM,Landsat 7 ETM+和Landsat 8 OLI的30m分辨率数据(表1),用于菹草群落的遥感提取。其中对于部分ETM+数据条带缺失问题,利用ENVI软件中的条带修复补丁进行修复。该过程主要包括缺失值修复以及修复后影像头文件信息的修改与添加,前者是基于单文件三角插值的条带填充,后者则是为了便于后续影像的辐射定标与大气校正。

表1 研究区Landsat数据介绍Tab.1 Landsat imagery used in this study

2 研究方法

2.1 总体处理流程

南四湖菹草在每年的4—5月开始快速生长并逐渐露出水面,4月底至5月中下旬大片泛滥并成为湖区沉水植物区的唯一优势种,进入6月逐渐衰亡腐烂,直至6月中旬消亡[17]。菹草这种特殊的物候特征可以被MODIS的时序归一化植被指数(NDVI)精确地记录。NDVI是植被生长状况和覆盖度变化的良好指标,其值被限定在[-1,1]之间[18-19]。本研究将2018年23景MODIS数据的NDVI波段进行叠加得到湖区的NDVI时序融合图像,并基于Landsat影像上典型地物的图像特征和文献获取的野外调查点[20]提取其NDVI时序曲线,分析物候特征的差异性并用于Landsat数据的菹草群落遥感识别[21],同时选取典型林地、水体、耕地、建设用地为训练样本,对研究区典型地物进行最大似然法监督分类,分类总体精度为99.47%,Kappa系数为0.98。利用分类得到的耕地掩膜去除动态监测分类结果中的耕地,最后得到湖区菹草群落的分布信息。整体流程如图2所示。

2.2 典型地物NDVI时序曲线提取

利用MCTK(Modis Conversion Toolkit)工具对MOD13Q1数据重投影为WGS-84坐标,将23景MODIS数据的NDVI波段叠加融合,得到时序融合影像。由于采集MODIS数据的过程中常常会由于诸如气溶胶和水汽等因素的影响,使得所提取的NDVI时序数据含有噪声,对后续的处理精度造成损失,故对图像进行S-G滤波(Savitzky-Golay Filter),改善曲线的整体平滑性以重构NDVI曲线。S-G滤波是一种在时域内基于局域多项式最小二乘法拟合的滤波方法[22],这种滤波的特点是可以同时去除噪声并保留原始信号的形状和宽度,效果如图3所示。

图2 菹草群落的遥感识别流程Fig.2 Flowchart of Potamogeton crispus L. remote sensing

基于前期调查数据和南四湖Landsat8标准假彩色影像特征[23],从MODISNDVI融合数据中提取典型地物的时序曲线。由于在5月前后未发现其他沉水植被大面积出露水面,故本研究提取出耕地、菹草、芦苇(PhragmitesCommunisTrin.)、荷田、水体5类典型地物的NDVI时序曲线,每类各5条曲线。

图3 S-G滤波处理前后NDVI曲线对比Fig.3 Comparison of NDVI curves before and after S-G filtering

在提取NDVI时序曲线时,对各类地物端元的选取尽量避免地物边缘混合像元的影响。最终NDVI时序曲线如图4所示。

图4 南四湖典型地物NDVI时序曲线Fig.4 NDVI time-series curves of typical objects in Nansi Lake

2.3 基于NDVI时序的菹草群落遥感提取方法

2.3.1湿地典型地物NDVI曲线特征分析

对滤波重构后的NDVI曲线进行对比发现,菹草物候与其他植被明显不同,具体分析如下:

1) 菹草的NDVI在第100d之后开始迅速增加,约120d左右达到峰值,此时处于菹草生长的旺盛期；而后菹草迅速衰亡,至175d左右NDVI降至低谷；峰谷两个时期的NDVI变化值均高于0.4。菹草的这种特殊物候特征不仅与湖区其他植被有显著差异,而且经调查发现其与湖区的狐尾藻(MyriophyllumverticillatumL.)、金鱼藻(CeratophyllumdemersumL.)等沉水植物群落差异也较为明显,后者在6—9月开花结果,生长繁殖期在菹草衰亡之后[24],175d以后菹草已经消亡,在影像中原菹草对应像元的区域,其他沉水植被重新占据原菹草区域,故NDVI曲线于第200d左右出现另一峰值。

2) 耕地和菹草的NDVI时序曲线除年初值和年末值以外,趋势十分相似,尤其是4—7月,由于人为收割、耕地,NDVI的下降值约为0.5,这与菹草的NDVI峰谷差值相似,故在后续菹草提取过程中,需要将耕地信息去除以防止与菹草混淆。本研究用5月份Landsat数据进行监督分类后得到的耕地类型对研究区的耕地进行掩膜去除。

3) 湿生植被芦苇从第100d(4月底)开始,NDVI逐渐上升,达到0.6左右并在整个生长期内较稳定,直至第300d(10月底)枯萎,NDVI降至0.2,其在菹草峰谷期的NDVI差值略高于0.1。

4) 荷田NDVI在菹草的旺盛期显著低于菹草,但菹草衰亡之后,荷田的NDVI值迅速攀升至0.6以上,两个时期差值高于0.3。

5) 敞水体的NDVI与其他地物差异性最为显著,全年均低于0.1,大部分时期低于0,且呈波动状态,故其在菹草的峰谷期NDVI值的变化落在正负两个区间内,这与其他地物不同。

基于菹草生长的NDVI峰、谷值信息及与其他植被类型的区别,一方面,可以指导我们对Landsat影像时相的选取,即每年两景,一景为4月下旬至5月上旬,一景为6月中下旬；另一方面,可以采用图像直接比较法进行手工阈值设定,并结合监督分类法以提取菹草群落。

2.3.2基于Landsat动态监测方法提取菹草

根据对南四湖典型地物NDVI时序曲线的分析结果,每个研究年份采用对应菹草生长旺盛期和衰亡期的两景Landsat数据,并进行波段计算,得到两个时期的NDVI图像,对其进行遥感动态监测。所选取的动态监测方法为图像直接比较法[25],该方法主要的技术流程为：先对两个时相的NDVI数据求差,差值结果为ΔNDVI=NDVI衰亡-NDVI旺盛,然后以与NDVI曲线分析结果相匹配的阈值系列为标准,对ΔNDVI数据进行分类。其中,菹草的区分区间为小于-0.4；荷田的区分区间为大于0.3,芦苇的提取区间则位于区间(0,0.1]。将变化值介于水体与菹草耕地之间地部分界定为非典型地物。借鉴决策树方法设定整体阈值思路,其路线如图5所示。

值得注意的是,由于耕地和菹草的NDVI趋势与差值特征相似,故最终得到的变化检测的分类结果是菹草和耕地的混合类型,需要借助Landsat监督分类得到的耕地掩膜进行耕地去除。

图5 图像直接比较法手工阈值设定Fig.5 Manual threshold setting based on imagery comparison

3 结果与讨论

3.1 南四湖菹草 芦苇和荷田的优势度变化

对湖区典型水生植物菹草、芦苇和荷田进行优势度指数(辛普森多样性指数)计算。某种群优势度的计算公式为:

D=Nmax/N

(1)

式中：Nmax为该种群的个体数量,N为所有物种个体总数[26]。一般优势度大于0.1的为绝对优势种,大于0.01的为主要优势种。

研究基于不同地物的NDVI阈值区间内的像元数代替物种数量以计算该指数。首先对2018年第80d到第180d的湖区不同植被NDVI统计数据进行拟合,得到该阶段内菹草、芦苇、荷田的优势度季节变化(图6)。结果显示,100d后,菹草逐渐成为整个湖区仅次于芦苇的绝对优势种,但由于芦苇和荷田为挺水植物,菹草为沉水植物,故这个时期菹草为沉水植物区的唯一优势种群,140d以后荷花取代菹草。再对1988—2018年第120d左右的湖区植被NDVI数据拟合,得到湖区优势度年际变化(图7)。可以看出,湖区中荷田的物种优势度经历了大幅下降又大幅上升的变化过程：芦苇一直以来都是湖区挺水植物的绝对优势种；菹草从1988年的优势度低于0.1逐年上升至优势度大于0.1。2000年以后,每年4—5月份菹草都是湖区沉水植物的绝对优势种,优势度指数维持在0.2以上。

3.2 南四湖菹草的时空演变特征

基于NDVI时序特征和Landsat系列数据,对1988—2018年间8个年份的菹草群落进行遥感提取,得到南四湖菹草群落时空分布图(图8),并对提取的结果统计菹草的面积和占比,具体情况如表2所示。

图6 2018年湖区主要植被优势度季节变化Fig.6 Main vegetation dominance change in 2018

图7 1988—2018年菹草旺盛期湖区植被优势度年际变化Fig.7 Vegetation dominance change in the prosperous period of Potamogeton crispus L.from 1988 to 2018

表2 1988—2018年南四湖菹草面积变化
Tab.2 Area change ofPotamogetoncrispusL. in NansiLake from 1988 to 2018

年份菹草面积/km2占湖区面积比/%占湿地面积比/%198825.223.611.98199632.874.712.58200040.325.773.16200673.3910.515.752010102.2614.648.022013188.7627.0314.802017193.7427.7415.192018171.4624.5513.44

图8 1988—2018南四湖菹草群落的时空演变Fig.8 Spatial and temporal distribution of Potamogeton crispus L.in Nansi Lake since 1988

1988—2006年间,南四湖中菹草呈局部聚集状零散分布于近湖岸区域,且前期独山湖和昭阳湖几乎没有分布。2006年前菹草面积由25.22km2增加至73.39km2,之后菹草逐渐蔓延至中部的独山湖和昭阳湖,分布区域主要在昭阳湖东部及独山湖,微山湖近岸区域也始有环状分布的大面积菹草群落。尤其到2010年后,湖区菹草开始大范围增加,由73.39km2增至171.46km2,面积占湖区百分比增加了14.04%。其中独山湖和昭阳湖于2013年菹草泛滥最为严重,菹草广布于除湖心区域之外的大面积水域,微山湖菹草于2017年后泛滥程度最高,由湖岸向湖心蔓延,其中东部湖区从微山岛周围至南部大面积水域,约1/2湖面被菹草占据,西部的养殖区也有广泛的分布。此时,相比于北部的南阳湖、独山湖和昭阳湖,微山湖的菹草泛滥面积更大,几乎由湖岸区蔓延至湖心。2013年之后菹草泛滥的情况在已有文献(如,马迎丽等[20])中已得到调查验证。2015年,研究者在33个南四湖调查样点中,有26个样点采集到菹草,且在湖区开阔水域的菹草生物量很高,这表明菹草在南四湖已大面积分布,有力地支撑了本研究的提取结果。

将提取的8个年份的菹草面积进行变化趋势二阶线性拟合[27],建立的曲线如图9所示,拟合曲线的R2=0.91。根据拟合效果可知,在当前人为管理和环境条件下,南四湖菹草群落仍会保持较大的分布规模。

图9 二阶线性拟合的南四湖菹草群落面积变化趋势Fig.9 Trend curves of Potamogeton crispus L. Area-Change-Line based on second-order linear fitting

3.3 南四湖菹草爆发的驱动因素

经过实地调研及查阅相关资料得知,2010年之后菹草的泛滥主要归因于南四湖湖水的富营养化[28],其次微山县政府开展的退渔还湖政策[34]以及湖区水位的变化也对菹草的生长产生了一定影响。

1) 水体富营养化。进入21世纪以来,南四湖已经成为典型的富营养型湖泊[29-33],水体富营养化指的是水体中N,P等营养盐含量过多而引起的水质污染现象。然而,菹草适宜在低浓度(总氮TN2.0～0.5mg/L,总磷TP0.20～0.05 mg/L)、中浓度(总氮TN2.0±0.3mg/L,总磷TP0.20±0.03 mg/L)内生长[35]。相关研究表明,在中低浓度下氮磷含量高更有利于菹草植株的生长[36]。在低氮低磷条件下,菹草叶片光合色素含量、植株可溶性糖含量显著提高；在高氮高磷条件下则相反。在低氮高磷和高氮低磷两种营养条件下,菹草繁殖体形成数增多[37]。根据相关文献整理的南四湖富营养化及水质数据如表3所示,由整理数据可知,虽然南四湖水体总磷(TP)、总氮(TN)含量呈现上升趋势,但仍大致保持在低中浓度范围内,有利于菹草的生长扩张。

分别采用Pearson相关分析和Spearman秩相关分析对遥感提取的菹草面积与TP、TN、综合富营养化程度和水质级别的关系进行定量研究,其中人为定义综合富营养化程度中的中度营养化为1、轻度富营养化为2、中度富营养化为3。

研究结果表明,南四湖菹草群落的扩张与南四湖水质级别不具备相关性(R<0.4,P>0.05)；但与水体TP、TN浓度和富营养化程度3个指标呈现显著正相关关系(R>0.8,P<0.05),说明南四湖水体TP、TN浓度和水体富营养化程度对菹草群落的扩张均有显著影响。菹草面积与相关变量的相关性结果如表4所示。

表3 南四湖富营养化程度及水质情况统计Tab.3 Statistics of eutrophication degree and water quality of Nansi Lake

表4 菹草面积与相关变量的相关性结果Tab.4 Relevant results of variables related to Potamogeton crispus L.

2) 退渔还湖。由于南四湖是南水北调东线工程最重要的调蓄湖和疏水廊道,为了保证南水北调输水水质,近年来湖区先后进行了渔业养殖箱和养殖网的压缩撤回,累计缩减率达到80%,缩减面积达129.4km2,累计释放的敞水面达到200km2[34]。该政策主要成效区域为南部的微山湖,结合图8可知,压缩养殖区后,大面积的敞水面空间可能使得菹草群落进一步在微山湖蔓延。

3) 水位与湿地演化。有文献表明,2000年以后南四湖水位不断升高,加之南水北调东线工程,导致水位进一步升高。水位变化使湖区植被生境条件改变,部分挺水植被如芦苇被浮水或沉水植被取代[38-39]。加之近几十年来受到人为干扰和水位波动影响,南四湖湿地在2014年加大治理措施之前,湿地呈现衰退演化趋势,期间南四湖湿地植被生境破碎、多样性降低,菹草等沉水植被随之扩张,一定程度上促进了湖区菹草的泛滥[40]。

3.4 不确定性分析

由于不同年份气候和水文环境的波动,菹草生长的旺盛期与衰亡期出现的时间可能会有一定的差异,对Landsat数据时相的选择造成一定的限制。另外,本研究所使用的方法适用于对某时期内较大浅水湖泊中泛滥成灾的优势水生植被菹草的识别,而对于一些小型湖泊或水生植被复杂多样的湖泊,并不适合利用250m分辨率的MODIS数据提取NDVI时序曲线,这些问题也是将来湿地植被遥感研究需要进一步细化解决的问题。

4 结论

1) 1988—2018年间,南四湖菹草群落呈现持续增长趋势,其时空演变主要分为两个阶段。第一阶段为2006年之前:菹草面积缓慢上升,该阶段菹草面积从25.22km2增加到73.39km2,占湖区面积的10.51%,总体涨势不突出；菹草群落主要分布于北部的南阳湖和南部的微山湖,分布范围相对较小,呈聚集小斑块状分布。第二阶段为2006年至今:菹草面积占比又增加了14.04%,2017年面积一度达到193.74km2的峰值(湖区面积占比达27.74%),2018年略降至171.46km2；该阶段菹草的分布开始呈现条带状蔓延,由南阳湖、微山湖向中部的独山湖和昭阳湖蔓延,且后期逐渐由湖边缓冲区域向湖心扩张。

2) 南四湖主要水生植物优势度表现出年际和季节的不同变化特征。年际尺度上,菹草从优势种变为绝对优势种,荷花的优势度先降后升,而芦苇一直是南四湖的绝对优势种。季节尺度上,菹草在每年的110—140d前后为湖区绝对优势种,随后腐烂,芦苇和荷花的优势度则较为稳定。

3) 南四湖菹草的泛滥受多方面因素影响,其中水体富营养化程度的增加显著促进了湖区菹草的扩张与泛滥,湖区管理政策可能加剧了菹草的扩张,水位的上升和南四湖湿地退化也在一定程度上促进了湖区菹草群落的扩张。