基于MODIS数据的神农架大九湖泥炭藓沼泽植被指数变化研究

庞毓雯,黄雨馨,俞立鹏,问静怡,吴玉环,徐俊锋,*

1 杭州师范大学遥感与地球科学研究院,杭州 311121 2 浙江省湿地与区域变化研究重点实验室,杭州 311121 3 浙江安吉小鲵国家级自然保护区,湖州 313300 4 杭州师范大学生命与环境科学学院，杭州 310036

泥炭藓(Sphagnum)是泥炭沼泽中最重要的固碳植物[1],对维持沼泽的结构和功能发挥关键作用[2]。全球一半以上的沼泽炭由泥炭藓形成,泥炭藓的固碳量约占土壤总碳量的15%[3]。每年由泥炭藓及其残体固定的碳比全球陆生植物积累的碳还要多[4- 5],这使得泥炭藓沼泽成为陆地生态系统中最重要的碳库之一,在全球碳循环中占据重要地位。泥炭藓沼泽一般分布于北半球极地、副极地地区[6],其发育对温度和水分有较高要求[7]。近几十年来,由全球变暖导致的水热格局变动,极端干旱、暴雨事件频发[8],改变了泥炭藓沼泽的地表水文和温度状况,增加了泥炭藓沼泽生态系统的脆弱性,削减了其碳汇能力[9]。已有研究指出,全球气候变化对泥炭藓沼泽的植被群落分布、植被生长等带来巨大挑战[10-11]。泥炭藓沼泽植被生长状况动态变化监测是了解气候变化对泥炭藓沼泽影响,预测泥炭藓沼泽演化[12]、生态系统稳定的良好指示器[13-14]。

目前,植被指数已被广泛应用于陆地植被生长状况监测以及植被生产力估算、碳源/碳汇估算的研究中。其中,归一化植被指数(Normalized Difference Vegetation Index, NDVI)和增强型植被指数(Enhanced Vegetation Index, EVI)应用最为广泛。国内学者,对于植被生长状况的监测进行了大量研究,从研究尺度上看,由全国[15-16]到省市级[17-18];植被类型以森林[19]、林草地[20]为主,湿地类型的植被监测较少涉及;在研究方向上,集中讨论了植被覆盖、生长状况时空变化及驱动因子[15- 20]。虽然NDVI/EVI被普遍认为是最佳的植被生长状况监测指标,但对于泥炭藓这类缺乏根系的植被类型,这两种指数的灵敏性评价研究仍需增强。

近年来,大量研究者开展了较多泥炭藓沼泽对气候变化的响应研究,泥炭藓已被广泛视为全球气候变化的生物放大器。Whinam和Copson等[21]通过历史资料分析了气候变化对泥炭藓生长的影响。研究指出在1999—2000年出现的干旱期破坏了麦夸里岛(Macquarie island)上泥炭藓的生长,以及未来的气候变化将影响泥炭藓的斑块大小和空间分布;Oke等[22]通过气候模型,模拟了2050年气候模式下,北美泥炭藓的生长状况,结论表明土壤水分缺乏和最干旱季的温度平衡将影响泥炭藓泥沼泽的分布,同时预测在不久的将来,泥炭藓沼泽适宜的气候面积可能会扩大。Loisel等[23]从已发表的文献中收集北方高纬泥炭沼泽泥炭藓生长的测量数据,以研究气候变化对泥炭藓生长的影响。研究发现由于全球变暖和生长季节的延长带来的光合有效辐射增加,在云量没有发生重大变化的情况下,将促进泥炭藓的生长。上述研究通过历史或预测资料分析了高纬地区泥炭藓沼泽在气候变化下的生长状况,得到了很好的结论,而不同区域泥炭藓对气候变化的响应存在差异,中纬地区也分布有大量泥炭藓沼泽,这些区域植被长势在气候变化格局下的时间周期规律尚未可知。

Boelman等[24]对以泥炭藓为主的苔原植被群落开展了15年的温室和施肥控制试验,并进行了2001—2003年为期3年的NDVI与生物量的关系研究。结论表明,NDVI是监测北极苔原植被群落地上植物生物量年际变化的有效的工具,同时指出,在利用NDVI估算与碳固存相关的生态参数时,必须了解感兴趣区存在哪些植被群落,以及以何种比例存在,强调在进行植被群落生长状况监测中,指标的选择具有特异性。May等[25]在研究短期地表水文变化对泥炭藓及苔藓群落NDVI、红光反射率、近红外反射率以及总初级生产力(Gross primary productivity, GPP)的影响时指出,NDVI和GPP与水文变化特征不匹配,用NDVI估算北方植被群落生产力可能存在问题,并强调了对苔藓群落盖度、组成和水分含量进行量化研究的必要性。Acunha等[26]评估了加拿大西部高度扰动泥炭藓沼泽的植被物候参数对沼泽恢复的响应。研究对比了MODIS的NDVI和蒸散量(Evapotranspiration, ET)在监测泥炭藓沼泽植被覆盖动态变化中的能力,指出NDVI比ET更有效地反映扰动事件(如沼泽火灾)的影响。这些研究都选用了NDVI指标进行泥炭藓沼泽植被生长状况监测,而针对不同的研究问题和研究区,NDVI的响应程度和灵敏度存在差异。

本文拟利用2001—2017年的MOD13Q1遥感植被指数(NDVI和EVI)产品,对我国罕有的中纬度亚高山泥炭沼泽——神农架大九湖泥炭藓沼泽进行植被指数时序分析,拟说明气候变化格局下,中纬泥炭藓沼泽植被的变化规律。本文研究内容包括：(1)Logistic模型在泥炭藓沼泽植被指数时间序列订正中的效果分析;(2)对比NDVI与EVI的泥炭藓沼泽植被生长状况监测能力;(3)得到泥炭藓沼泽植被季节生长规律及近18年的年际变化趋势。本文对比了两种最广泛的植被指数在表征不同时间周期上泥炭藓沼泽植被生长状况的能力,这将帮助我们改进从遥感植被指数数据集中反演其他重要生态变量的精度。另外,本文将作为泥炭藓沼泽植被变化研究的重要补充,能够进一步说明气候格局对不同纬度带泥炭藓沼泽植被生长状况的影响,为进一步分析气候变化对泥炭藓沼泽生态系统功能的发挥、空间分布的演变提供数据和理论支撑。

1 研究区

我国泥炭沼泽主要分布于大、小兴安岭地区等北方高寒地区[14,27],亚热带少量分布,如安徽徽州天湖山[28]、福建天宝岩国家级自然保护区[29]、浙江安吉龙王山[30]等。鄂西亚高山地带有较大面积的泥炭沼泽,且保存了中纬度地区少见的泥炭藓沼泽[31]。湖北神农架林区大九湖湿地在地球上已存在3万年以上,湿地核心区的藓类沼泽维持了最原始的样貌,泥炭藓盖度接近100%[32]。近几十年来,在人类活动影响下,大九湖泥炭藓沼泽受到了较为严重的干扰,目前在政府及相关部门的政策支持下,已开展大量沼泽恢复工作。杜耘等[33]针对神农架林区大九湖湿地退化问题进行了实地调查,指出20世纪70年代以来大规模的人为水文环境改造是湿地生态退化的主要原因。罗涛等[34]评估了大九湖湿地的植物群落分布现状,指出退耕还草、还泽等政策产生了积极效果,改善了湿地生态系统群落多样性。赵素婷等[35]利用中巴资源卫星与实地调查结合的方式,得到鄂西高山泥炭藓沼泽的分布、面积、物种资源状况。研究中,已建立相应的实地观测站点,黄咸雨等[32]指出在泥炭藓沼泽监测中,遥感能为地面观测站点提供有力补充。已有研究大多采用实地调查的方式勘察泥炭藓沼泽植被资源现状,研究区域受局限、时间周期不规律,不能够及时、全面的反映大九湖泥炭藓沼泽植被资源现状和发展趋势,而遥感技术的应用则能提供快速、历史回溯的资料,将极大提升研究质量和效率。

神农架大九湖湿地公园,位于湖北省神农架林区西北部,地理坐标为31°24′—31°33′N,109°56′—110°11′E。大九湖湿地发育于高山盆地,总面积为1645 hm2,其中泥炭沼泽有779 hm2[36]。泥炭藓沼泽分布于大九湖一凹型盆地中,沼泽的主要植被类型为红穗苔草-泥炭藓群落(Com.Carexargyi-Sphagnumpalustre)、紫羊茅-泥炭藓群落(Ass.Festucarubra-Sphagnumpalustre)等[31,35]。

本文在大九湖湿地公园中选择泥炭藓沼泽实验斑块两处,分别记为Bog1和Bog2(图1)。其中,根据文献[31,34]确定泥炭藓沼泽斑块Bog1(31°28′51.19″N, 110°0′51.19″E),该处主要植物物种为泥炭藓、紫羊茅等;根据文献[35],由中巴资源卫星遥感提取的泥炭藓沼泽斑块为Bog2(31°30′45.413″N, 110°0′7.273″E),该处主要植物物种为泥炭藓、紫羊茅等。两处实验泥炭藓沼泽斑块植被类型均属紫羊茅-泥炭藓群落(Ass.Festucarubra-Sphagnumpalustre),草本层的平均盖度为45%,优势种紫羊茅,平均高度为0.6 m,地被层平均盖度为100%,其中泥炭藓平均盖度为95%,藓丘的平均高度为0.22 m。从研究尺度上来看,两处实验泥炭藓沼泽斑块分别代表实地勘察的小尺度和遥感水平中等尺度,两斑块同时分析既能够交叉验证,也能够说明遥感泥炭藓沼泽植被生长监测是否受到尺度效应影响。

2 材料和方法

2.1 遥感数据

本文所用遥感数据为美国国家航空航天局(NASA)提供的EOS MODIS(Terra)产品系列中的MOD13Q1,该数据集包括了NDVI和EVI,空间分辨率为250 m,时间分辨率为16 d。数据集经过辐射校正、几何校正、大气校正等处理有效提高了数据质量,广泛应用于区域植被变化监测中。

所用数据利用美国地质勘探局(USGS)开发的数据平台AppEEARS(Application for Extracting and Exploring Analysis Ready Samples)(https://lpdaacsvc.cr.usgs.gov/appeears)获得。该平台提供了用户自定义研究点位(Points)或研究区(Area)的两种感兴趣类型的定制NASA卫星产品的时空快速分析及可视化服务。

MOD13Q1的NDVI和EVI产品为16 d数据,在数据生产过程中会受到云雪、雨水土壤背景、传感器自身等影响,出现异常值或缺失数据。因此必须进行植被指数订正才能真实反映植被的季节和年际变化规律。

2.2 遥感植被指数的时间序列订正

植被指数是遥感监测植被生长过程的主要媒介,植被指数时间曲线可以表征植被生长过程[37]。在实际应用过程中,卫星传感器在获取地表信息时,太阳和传感器角度、云层和天气状况,以及地表水体和冰雪覆盖等[38]都会对植被指数产品的数据质量产生影响,使得原始植被指数时间曲线表现为与植物生长过程不符的非平稳变化的锯齿状(图2虚线)。

植物的生长过程一般先上升,再逐渐下降,这种规律可以通过Logistic模型来模拟。Zhang等[39]提出分段逻辑斯蒂克回归分析模型(Stepwise logistic)来模拟植物的季节变化,目前该模型被NASA采纳作为监测地表植物物候的核心算法。模型为：

(1)

式中,y(t)为t时刻的NDVI/EVI值,a,b为拟合参数,d为NDVI/EVI初始背景值,c+d为最大NDVI/EVI值[39]。本文利用该模型思想对NDVI、EVI时间序列进行订正,进而由此提取泥炭藓沼泽植被的生长规律(图2)。

图2 泥炭藓沼泽的季节周期Logistic订正Fig.2 Logistic correction of seasonal period of Sphagnum palustre bog从图中可以提取植被的物候参数：生长季开始(季节性振荡的20%)、生长季最大值、生长季结束

2.3 趋势分析

为了分析泥炭藓沼泽植被生长状况的变化趋势,本文利用一元线性回归模型统计植被指数(NDVI和EVI)在2000—2017年的年际变化趋势。其中,分别提取植被指数年均值(代表逐年泥炭藓沼泽植被的平均生长水平,这是最常用的趋势分析指标)和植被指数年最大值(表征泥炭藓沼泽植被年内的最佳生长状态,是更灵活的植被生长监测指标)。由此,分别统计泥炭藓沼泽植被指数年均值和年最大值的年际变化趋势K—为一元线性回归模型的斜率,其计算公式分别为：

(2)

(3)

式中,Kavg代表植被指数年均值的年际变化趋势线的斜率,Kmax代表植被指数年最大值的年际变化趋势线的斜率;n为遥感产品累计的时间周期,本文时间周期为2000—2017年,即n=18;aVIi代表第i年植被指数年均值;mVIi代表第i年植被指数年最大值。若K>0,说明时间周期内泥炭藓沼泽植被生长状况呈增加趋势,反之,则呈下降趋势。随后,对线性回归的植被指数年际变化趋势进行检验,其中0.01

3 结果与分析

3.1 植被指数订正

3.1.1季节变化与订正

本文分别对NDVI和EVI多年逐月均值应用Logistic模型,得到结果如图2(实线)。季节周期的植被指数时间序列曲线能够反映泥炭藓沼泽植被一年的生长规律。订正后的生长曲线更为平缓,能够获得泥炭藓沼泽植被任意时间结点的植被指数值,从而估计沼泽植被各物候期或生长过程关键节点。图2中NDVI曲线在6月中旬出现一个谷,图2中EVI振荡出现在年初及年末,8月出现的异常高值,经过Logistic订正后这些噪声都得到去除。

3.1.2年际变化与订正

植被指数产品在原始处理时以年为单位,因此在相邻年份衔接上,数值振荡和差异明显。本文分别对两个泥炭藓沼泽斑块的NDVI及EVI进行2000—2017年的年际周期时间序列订正分析(图3),经Logistic订正后,相邻年份间指数值变化平缓,更符合真实规律。从整体上看,Logistic订正后的植被指数时间序列变化趋势平缓,无跳跃值,可反映泥炭藓沼泽植被长时间序列的变化趋势。

图3 泥炭藓沼泽的年际周期Logistic订正Fig.3 Logistic correction of annual period of Sphagnum palustre bog

3.2 泥炭藓沼泽植被生长遥感监测指标选取

陆地生态系统中植被的长时间序列变化研究多借助遥感卫星数据。用于植被生长监测最广泛的两种遥感植被指数为归一化植被指数(NDVI)和增强型植被指数(EVI)。国内外学者大量讨论了NDVI与EVI在不同植被群落、生态系统中的响应能力[40-42],研究结果表明NDVI在监测较高覆盖度植被时容易发生饱和,而EVI通过调节参数克服了易饱和的缺点,在植被覆盖度较高的区域监测精度优于NDVI。

泥炭藓是一种缺乏根系的苔藓植物,其叶片和冠层结构与其他针叶和阔叶维管植物不同,NDVI和EVI在泥炭藓植物中的适用性还需要进一步分析。本文比较了订正后的NDVI与EVI应用于泥炭藓沼泽植被生长状况监测中的适应性和灵敏度。

3.2.1NDVI与EVI季节对比分析

图4 泥炭藓沼泽NDVI与EVI季节生长周期 Fig.4 NDVI and EVI seasonal growth cycle of Sphagnum palustre bog

泥炭藓沼泽植被的季节生长周期贯穿4月中旬至10月中旬(图4)。基于订正后的NDVI和EVI多年月份均值的沼泽植被周期曲线(图4),可知：1)沼泽植被为一年生,NDVI与EVI生长曲线为单峰,且生长季时间跨度较长,集中在夏季,与雨水充沛促进泥炭藓生长有关;2)NDVI与EVI得到的生长季时间大致相等,即生长季开始时间和结束时间近似。由于NDVI在年初起始值较高,随着植被生长旺盛,NDVI值趋于饱和,对生长季最大值变化检测不敏感,此时EVI更适用;3)EVI拟合后的整体效果优于NDVI。在之后的泥炭藓群落物候期规律定量变化研究中,应使用EVI作为监测指标。

3.2.2NDVI与EVI年际对比分析

统计2000—2017年2个泥炭藓沼泽斑块的年际NDVI和EVI变化规律,如图5所示：1)NDVI整体高于EVI,NDVI高值处于0.7—0.8之间,EVI高值位于0.5—0.6附近。两者总体趋势一致,在相邻时间周期内,表现出相似的变化特征;2)从长时间变化来看,相较于EVI,NDVI的年最大值变化较平缓,EVI的峰值在年际间的振荡更明显,说明EVI能更灵活监测泥炭藓沼泽植被年际最佳生长状况的变化;3)在曲线的谷值处,NDVI振荡强于EVI。总体而言,EVI在年际间泥炭藓沼泽植被动态监测中优于NDVI。

图5 泥炭藓沼泽斑块的NDVI和EVI年际变化Fig.5 NDVI and EVI annual change of Sphagnum palustre bog

3.3 泥炭藓沼泽植被时间序列变化分析

3.3.1泥炭藓沼泽植被季节规律分析

泥炭藓沼泽植被的季节EVI变化曲线(图4)可知：泥炭藓沼泽植被为一年生,生长季开始于4月中旬,结束于10中旬,生长季长度为6个月。一年中的1月至4月泥炭藓沼泽的EVI值处于较低水平,随着冰雪融化,雨季的到来,EVI开始缓慢上升,进入生长季,直至10月中旬结束生长季。泥炭藓对地表及地下水文变化敏感,从生长周期结点变化能够灵敏反映短期的水文格局变化,体现一年中极端降水和干旱事件。通过长时间的物候和沼泽地表水文状况联动分析,能够灵敏表征泥炭藓沼泽在气候变暖引起水热格局变化中的响应。

3.3.2泥炭藓沼泽植被生长年际变化分析

本文分别对2000—2017年大九湖两个泥炭藓沼泽斑块的EVI年均值及年最大值进行线性回归,结果如图6—7所示：1)无论是年均值还是年最大值泥炭藓沼泽的EVI都表现为平稳小幅增长趋势,年际变化率K均大于0。2)泥炭藓沼泽EVI年最大值的变化曲线的趋势程度强于年均值曲线,即年际变化率Kmax=0.0036>Kavg=0.0018,且对应R2更大(R2=0.33,P<0.05)。3)对比两个泥炭藓沼泽斑块,斑块2的趋势程度及相关性都优于斑块1。在年际变化中,两个斑块的变化规律大体一致,但在个别年份上存在差异,该现象在EVI最大值中表现明显。例如,斑块1中2007年EVI年最大值出现一个相邻年份间的极大值,而在斑块2中,该年份与相邻年份EVI值近似。该现象出现可能是由于泥炭藓沼泽植被生长对微地貌敏感或受到遥感监测尺度效应的影响。

图6 泥炭藓沼泽斑块1 EVI的变化趋势Fig.6 The trend of EVI in the Sphagnum palustre bog1

图7 泥炭藓沼泽斑块2 EVI的变化趋势Fig.7 The trend of EVI in the Sphagnum palustre bog2

随后,将泥炭藓沼泽斑块1和2的年际EVI数据汇总进行全局时序分析,结果如图8所示：整体上,泥炭藓沼泽植被在18年来呈现稳定的上升趋势,EVI年最大值的变化趋势较EVI年最大值更强(Kmax=0.0058>Kavg=0.0038),且趋势相关性更好(R2=0.47,P<0.01)。

图8 泥炭藓沼泽斑块1和2的EVI变化趋势Fig.8 The trend of EVI in the Sphagnum palustre of bog1 and bog2

4 总结

本文利用2000—2017年MODIS NDVI和EVI数据,对我国独特的中纬度亚高山大九湖湿地中的两个泥炭藓沼泽斑块的植被生长状况进行了时间序列趋势分析,得到以下结论：

(1)Logistic模型能够有效消除原始植被指数时间序列中的噪声,在季节和年际曲线订正中表现出良好效果。中纬度高山地区年均云层覆盖率高,原始的植被指数时间周期产品有较大的噪声。无论是季节周期上夏季的异常低值,还是相邻年份间的数值跳动,logistic模型均能够很好订正这些误差,得到平滑的泥炭藓沼泽植被生长规律曲线及年际间的周期规律;

(2)泥炭藓沼泽植被生长状况时间序列变化动态监测中,EVI优于NDVI,表现在：1)季节生长周期监测中,EVI能够更精确的反映植物物候变化节点。NDVI在年初的初始值较高,随着植被生长的旺盛,值域变化空间较小,相比于EVI对植被生长季变化的感应能力较弱,无法准确捕捉泥炭藓沼泽植被的成熟期。由于,泥炭藓沼泽地表在生长季旺盛期泥炭藓藓丘覆盖度较高,约为90%,在这种状况下,NDVI检测植被生长变化的能力弱于EVI。随着植被的枯黄,NDVI下降,有了更大的值域响应范围,因此NDVI与EVI均能够反映植被的衰老期。2)在植被长势的年际周期分析中,NDVI在年际间的变化规律十分微弱,而EVI能够体现年份间的植被长势的强弱状况。由于NDVI在泥炭藓沼泽植被上易饱和的不足,不能够及时体现泥炭藓植被的生长变化,此时EVI更适用;

(3)经趋势分析可得,神农架大九湖泥炭藓沼泽植被表现为微弱增长,EVI年均值增长率达3.8‰(R2=0.46,P<0.01),泥炭藓沼泽植被生长稳定。本文的趋势分析指标分别为EVI年均值和年最大值,年最大值的趋势程度为5.8‰,强于年均值,表明泥炭藓沼泽植被的季节生长周期对气候变化的响应更为灵敏。泥炭藓是一种特殊的藓类,没有根系和维管组织对外界环境变化极其敏感,本文选取的指标(EVI年均值和年最大值)不够精细,可能未能及时捕捉到泥炭藓对气候变化的快速反应,因此在后续研究中应引入物候规律分析,包括植被的返青期、成熟期、衰老期及生长季长度的变化等。另外,本文基于MODIS产品的植被指数变化分析,时间周期为2000—2017年,在未来的研究中,可以叠加更多的分析数据(如AVHRR GIMMSS产品等),从更长的时间周期上,探究气候变化对泥炭藓沼泽植被生长的影响;

(4)此外,本文选取的两个泥炭藓沼泽斑块通过遥感的手段将传统的研究单元由点尺度提升到面尺度,由小区域跨越到大区域。两个斑块的趋势分析体现了较好的一致性,均表现了泥炭藓沼泽植被的微弱的增长趋势(Bog1年均值增长率1.8‰,Bog2年均值增长率5.8‰),但在趋势程度上也表现出了略微差异。这既展现了遥感在泥炭沼泽生态监测中的应用能力,也体现了在不同尺度分析上,遥感的尺度效应。因此在后续研究中,可进一步分析由尺度效应引起遥感监测效果差异,及遥感应用于泥炭藓沼泽生态监测的最佳研究尺度。

致谢：贾文晓博士和宋垚彬博士帮助植被群落表述,特此致谢。