黄河流域植被物候遥感提取及空间分布

罗婉琳 韩 琳＊ 吴满玉 杨佳音 侯亚坤

（郑州师范学院地理与旅游学院，河南 郑州 450000）

植被物候是指植物在一年的生长中，随着气候的季节性变化植物自身发生萌芽、抽枝、展叶、开花、结果及落叶、休眠等规律性变化的现象，与地貌特征及地理空间分布密切相关，随降水、温度、光照等气候条件的变化而变化[1]。植被物候对气候变化具有敏感性，对于认识自然季节现象变化的规律，农作物的生产、全球气候变化、生态学的应用等方面具有重要的应用价值。

黄河流域地处我国干旱，半干旱和半湿润地区，生态环境脆弱，在全球变暖和人类活动共同影响下，流域气候及水文过程发生了显著变化，流域整体气候的暖干化和人类用水的不断增加使得黄河流域的水文气象变化不断加剧，该流域的生态问题变得越来越受人关注，了解黄河流域植被生长变化空间格局，对黄河流域生态环境保护具有重要意义。

近年来，遥感技术发展迅速，传感器能够提供多时相、多分辨率相对准确和连续的遥感数据，从而为实现大区域尺度、长时间序列的植被生长动态监测提供了可能。基于NDVI 时间序列具有较好的季节节律能够综合反映植物季相变化特征，结合地理信息系统技术进行空间格局及其动态变化分析，成为量化分析区域植被物候动态变化及其对气候变化响应研究的重要手段[2]。本文通过探讨植被物候遥感提取方法，对黄河流域植被生长开始期、植被生长结束期和植被生长长度等三个植被物候参数进行了提取，并对近20a 来流域植被物候空间分布进行了分析，为黄河流域生态环境治理提供信息参考。

1 研究区域与数据源

1.1 黄河流域概况。黄河流域位于96°～119°E，32°～42°N，东西长约1900 km，南北宽约1100 km，流域面积约为79.46×104km2。黄河流域幅员辽阔，地貌差异大，该流域发源于青海巴颜喀拉山，从西到东横跨青藏高原、内蒙古高原、黄土高原和黄淮海平原4个地貌单元。流域地势西高东低，西部河源地区平均海拔在4000m以上，由一系列高山组成；中部地区海拔在1000～2000m 之间，为黄土地貌，水土流失严重；东部海拔不超过100m，主要由黄河冲积平原形成。黄河流域属于大陆性气候，东南部属半湿润气候，中部属半干早气候，西北部属干旱气候。植被从高山草甸、灌木林到耕地林网，类型丰富。

1.2 NDVI 数据来源与处理。NDVI 序列影像采用2001～2019年中分辨率成像光谱仪（MODIS）陆地标准产品（MOD13Q1）中的归一化植被指数（NDVI）数据产品，进行物候分析。产品从NASA 官网获取（https://www.nasa.gov/）。MOD13Q1 中的NDVI 数据产品的空间分辨率为250m，是经过16d 最大值合成的反射率数据。运用MRT软件对MODIS 数据进行数据拼接、重投影和波段信息提取等，然后利用Python 在cmd 命令中调取相关Arcgis 工具进行数据批量的数据裁剪、投影转换、格式转换，最终输出成Albers 投影的NDVI数据。

2 植被物候提取方法

基于时间序列植被指数提取植被物候信息一般需要经过两个步骤[3]：一是NDVI 指数的平滑与重构，以消除数据获取时受云雾、大气等因素影响而产生的噪声；二是对平滑后的NDVI 序列数据曲线提取物候参数。

2.1 NDVI 植被植数。NDVI（Normal）是遥感图像研究植被特征的重要指标，植被信息主要通过绿色植物叶子和植被冠层的光谱特性及其差异、变化来反映。但对于复杂的植被遥感，仅用个别波段或多个单波段数据对比分析提取植被信息存在局限性，因此往往选用不同波段之间的运算来产生对植物长势、生物量等有一定知识意义的数值来分析植被变化。归一化植被指数（NDVI）是常采用的指数之一，其计算方法见公式（1）。

式中：Rnir指近红外波段地物的反射率，Rnir指可见光波段地物的反射率。

2.2 NDVI 数据降噪与重构。NDVI 反映植被生长过程变化的时间序列曲线理论上是连续且平滑的，虽然经过16d 最大合成的MODIS 数据虽减弱了云的影响，但由于传感器性能、天气等随机因素的干扰，NDVI 数据仍然有一些噪声的存在，使NDVI 时序数据出现不规则情况，因此对数据进行平滑降噪处理是必不可少的一个工作。常用的方法有双Logistic 曲线拟合、非对称高斯函数拟合和Savitzky-Golay（S-G）滤波法。许多学者对这三种常用方法进行了研究，杨恒等[4]运用S-G 滤波的方法对江西省植被覆盖度时空变化进行了相关分析研究。丁潇等[5]对比分析了Savitsky-Golay 滤波法、Logistic 调合函数模型和非对称Gaussian 模型3 种方法对植被指数曲线拟合重构的效果，最后得出Savitsky-Golay 滤波法在与原始曲线近似度方面优于另外两种拟合方法。本研究利用TIMESAT3.3 软件，基于S-G 滤波方法对NDVI 数据进行降噪处理，从而对时间序列进行重构。滤波方法见公式（2）。

式中，Yj为拟合之后的序列数据，Yj+i为原始序列数据，Ci为滤波系数，N 为滑动窗口的大小（2m + 1）。本研究利用Savitsky-Golay滤波对NDVI 时间序列进行去噪，窗口大小为设置3，并采用NDVImean=0.05 作为排除异常阈值。图1 为某一像元基于S-G 滤波前后的NDVI 时间序列对比。

图1 时间序列重建效果图

2.3 基于序列NDVI 植被物候信息提取

基于遥感技术的植被物候信息提取是通过在时间序列NDVI形成的曲线，来定义植被物候参数。植被物候各参数在NDVI 时间序列曲线上含义见图2 所示。则是生长季期间NDVI 的积分;h+i 为大积分值;j, k 分别为NDVI 左右导数。

图2 植被物候参数[6]

当前，基于时间序列NDVI 植被物候参数提取方法有动态阈值法、滑动平均法、导数法、曲线拟合法、最大变化率、Logistic 函数拟合法等，其中阈值法是大尺度区域植被物候提取最为常用的方法[3]。本研究结合先验知识及区域特点，重点讨论动态阈值法植被物候信息提取。

动态阈值法是基于遥感图像的统计学特征而确定植被的物候变量的，即将NDVI 增长达到当年NDVI 振幅一定百分比的时刻定义为生长季的开始时间，而NDVI 降低到当年NDVI 振幅一定比例的时刻定义为生长季的结束时间。它能有效地识别不同土地覆盖类型的特点，得出有针对性的阈值。动态阈值法在进行大区域多年数据的物候信息提取分析时，有效的避免了在大区域不同植被类型植被指数差异导致合理阈值范围选取困难的问题，另外该方法使得不同年份NDVI 变异也不会相互影响，往往能够得到较准确的结果。

利用S-G 滤波后重建的长期序列数据，借助TIMESAT3.2 软件，利用动态阈值法提取2001 至2019 年植被生长季开始时间（start of season，SOS）、生长季结束时间（end of season, EOS）和生长季长度（Length of season, LOS）三个物候参数，对黄河流域近20a 来植被物候空间分布进行分析。根据本研究区的实测物候数据对比并与前人的研究结果分析，通过反复试验最后本文的植被SOS 和植被EOS 提取阈值为30%、40%。即将NDVI 曲线上升阶段，距离最小值为最大值与最小值间距离的30%的时间点定义为植被SOS，将NDVI 曲线下降阶段，距离最大值为最大值与最小值间距离的40%的时间点定义为植被EOS。

3 结果与分析

3.1 黄河流域不同年份植被SOS 空间分布。图3（a）～（d）分别是黄河流域2016～2019 近4 年植被SOS 空间分布图，表1 是各年份SOS 不同时间区间占整个流域植被覆盖区的面积占比。由图3 和表1 分析可见，黄河流域植被SOS 主要集中在第100-160 天，该区间平均约占整个流域植被覆盖区的71.8%，其中，2017 年的占比最高，约占75.2%；2018 年的占比最少，约占65.7%。流域植被SOS 在第120-140 天的占比最大，该区间植被主要分布在黄土高原中部地区。植被SOS 在小于第80 天的占比最小，主要集中在黄河下游及中游关中平原部分区域，该区域主要以农作物种植为主，2018 年该区间空间分布占比明显增大。流域SOS 在第80-100 天的空间分布区域在2018 年、2019 年较2016 年、2017 年发生了明显增大，而在第160-180 天的空间分布区域则明显减小。

图3 黄河流域不同年份SOS 空间分布

表1 黄河流域SOS 不同时间区间面积占比（%）

3.2 黄河流域不同年份植被EOS 空间分布。图4（a）～（d）分别是黄河流域2016～2019 近4 年植被EOS 空间分布图，表2 是各年份EOS 不同时间区间占整个流域植被覆盖区的面积占比。由图4和表2 分析可见，黄河流域植被EOS 主要集中在第290 天之后，该区间平均约占整个流域植被覆盖区的83.9%，其中，2019 年的占比最高，约占87.9%；2018 年的占比最少，约占79.0%。流域植被EOS在第300-310 天的占比最大，该区间植被主要分布在黄土高原中南部以及流域东部的部分区域。植被EOS 在小于第270 天的占比最小，主要集中在黄河下游及内蒙古河套平原区域，该区域主要以农作物种植为主。

表2 黄河流域EOS 不同时间区间面积占比（%）

图4 黄河流域不同年份EOS 空间分布

3.3 黄河流域近20a 植被物候均值空间分布特征

3.3.1 植被SOS 均值空间分布特征。黄河流域2001～2019 年植被SOS 均值空间分布见图5 所示。黄河流域由于植被随地形以及水热条件的差异而不同，不同区域植被物候参数差异较大，植被SOS 呈现出由东南向西北逐渐推迟的空间分布格局，整个流域植被SOS 主要集中在第100-160 天，即3 月初到5 月上旬。其中，流域上游海拔较高的青藏高原区植被SOS 主要集中在第140-160天，该区植被生长期开始季最晚；流域中部的黄土高原区植被SOS大部分集中在第120-140 天；黄土高原南部及黄河下游大部分区域植被SOS 主要集中在第100 天以内，该区域植被生长季开始最早，该区域植被主要以人工种植植被为主。

图5 黄河流域2001-2019 年植被SOS 空间分布

3.3.2 植被EOS 均值空间分布特征。黄河流域2001～2019 年植被EOS 空间分布见图6 所示。流域植被EOS 空间分布与SOS 几乎相反，即EOS 呈现由西向东、由北向南逐渐推迟的空间分布格局。就整个流域而言，大部分地区自然植被EOS 集中在第270-310天左右，即10 月上旬至11 月上旬。其中，流域东部的青藏高原区植被EOS 主要集中在270-290 天，该区结束期开始最早；流域中部黄土高原区中部及北部自然植被EOS 主要集中在第290-310天，部分区域EOS 则在第270 天之前；黄土高原南部及下游区域自然植被EOS 主要集中在第300 天之后，该区域植被结束期最晚。

图6 黄河流域2001-2019 年植被EOS 空间分布

3.3.3 植被LOS 均值空间分布特征。黄河流域2001～2019 年植被LOS 空间分布见图7 所示。流域植被LOS 呈现由西向东、由北向南逐渐增大的空间分布格局。流域大部分区域植被LOS 在110-200 天左右。其中，流域东部的青藏高原区植被LOS 主要在130 天左右，该区自然植被生长季长度相对较短；流域中部黄土高原中部及北部地区植被LOS 主要集中在150-170 天左右；黄土高原南部及下游区域植被LOS 大于170 天，该区域自然植被生长季最长。

图7 黄河流域2001-2019 年植被LOS 空间分布

4 结论

本文结合黄河流域尺度特点，探讨了基于遥感技术的植被物候信息遥感提取方法。在此基础上，基于2001～2019 年MODIS NDVI 数据，利用动态阈值法提取了近20a 来黄河流域植被物候信息，对流域植被物候信息进行空间分布分析。结论如下：

4.1 黄河流域植被SOS 主要集中在第100-160 天，其中上游海拔较高的青藏高原区植被SOS 主要集中在第140-160 天，黄土高原中部及北部区域植被SOS 大部分集中在第120-140 天，黄土高原南部及黄河下游大部分区域植被SOS 主要集中在第100 天以内。流域植被SOS 呈由东南向西北逐渐推迟的空间分布格局。

4.2 黄河流域植被EOS 主要集中在第270-310 天，即10 月上旬至11 月上旬。其中流域东部的青藏高原区植被EOS 主要集中在270-290 天，流域中部黄土高原区自然植被EOS 主要集中在第290-310 天，黄土高原南部及下游区域自然植被EOS 主要集中在第300 天之后。流域植被EOS 的空间分布格局与SOS 的分布基本相反，即由西北向东南逐渐推迟。

4.3 黄河流域植被LOS 主要集中在110-200 天左右。其中，流域东部青藏高原区植被LOS 主要在130 天左右，流域中部黄土高原中部及北部地区植被LOS 主要集中在150-170 天左右；黄土高原南部及下游区域植被LOS 大于170 天。流域植被LOS 整体呈现由西向东、由北向南逐渐增大的空间分布格局。