1.1.1 NDVI数据 采用由西部环境与生态科学数据中心下载(http://westdc.westgis.ac.cn)的美国国家航空航天局戈达德航天中心(Goddard Space Flight Center, GSFC)全球监测与模型研究组制作的1986-2015年15 d最大化合成的8 km GIMMS NDVI 3g数据集(版本3)。该数据集已经过辐射校正、几何校正且消除了因传感器更换、卫星轨道漂移和太阳高度角变化对数据质量的影响^[28], 相关研究证明了NDVI 3g数据比NDVI g数据在中高纬度地区具有更高的质量^[29]。因此本研究利用的GIMMS NDVI 3g数据能较好地保证数据质量, 来分析我国草地活动状况且具有一定的可靠性。

需要说明的是:下载的GIMMS NDVI 3g数据在叠加中国草地数据的基础上, 每期均有小于-1的值。而理论上, NDVI的范围在-1~1。在冰雪覆盖或水域地区, NDVI常常为负值; 在戈壁沙漠等非植被盖地表, NDVI又不稳定, 国际上通常采用0.05或0.1的NDVI年均值作为阈值^{[30, 31]}来排除非植被因素的影响。考虑到高纬度、高海拔地区的高山亚高山草地和荒漠草地也是本次研究重要的草地类型, NDVI值年均较低。基于上述情况又参考梁爽^[32]对中国北方草地物候遥感监测时的研究, 本研究采用0.05的NDVI年均值作为阈值, 将年均NDVI小于0.05的像元排除在中国草地研究区之外, 投影定义为Albers投影, 以便进行下一步分析。

1.1.2 植被类型数据 本研究使用的草地类型数据来源于国家自然科学基金委员会中国西部环境与生态科(http://westdc.westgis.ac.cn)的GLC2000数据。是由中科院遥感所承担开发的基于SPOT4遥感数据的全球土地覆盖数据中国子集, 数据的空间分辨率为1 km, 共分为24种土地覆盖类型。其中6种草地类型(图1)包括高山亚高山草甸(alpine and sub_alpine meadow)、山地草原(slope grassland)、平地草原(plain grassland)、荒漠草原(desert grassland)、草甸(meadow)和高山亚高山平地草原(alpine and sub-alpine plain grassland)。此外, GLC2000产品在1:10万的尺度上与中国植被覆盖数据的一致性最高^[33]。定义为Albers投影, 利用栅格属性提取全国6种草地植被类型。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 中国草地类型分布Fig.1 Distribution map of grassland types in China

1.1.3 中国草地物候验证资料 本研究采用的野外物候观测资料来自国家生态系统观测研究网络数据共享网(http://cemdisl.cern.ac.cn/index. html, CNERN)发布的海北、奈曼、策勒、鄂尔多斯、阜康、临泽、沙坡头7个站点的物候观测数据, 1986-1988年《中国动植物物候观测年报》第10~11号中嫩江、伊春、牡丹江、虎林、石河子、乌鲁木齐等12个站点的数据, 1986-1993和2011年《保定动植物物候观测年报》实测数据以及来自文献资料的武威^[34]荒漠生态站点的草地物候观测资料, 草地物候相关论文的研究成果^{[35, 36, 37, 38, 39, 40]}。

首先将日期型的站点物候观测资料转换成数值型数据, 根据站点的经纬度坐标, 提取包括站点所在像元及邻近共9个像元的平均值, 作为该站点物候期的遥感提取结果; 然后根据侯学会^[36]的研究剔除记录时间与其他数据相差30 d以上的记录, 以增加验证数据的稳定; 再把开始展叶期和开始黄枯期定为草地的SOS和EOS; 最后提取站点地区的像元均值, 作为该站点物候期的遥感提取数据。

为了减少噪声对分析结果的影响, 本研究首先逐像元计算了中国草地植被的NDVI时间累积曲线, 然后对时间累积的NDVI曲线利用全局S-G滤波曲线拟合, 并用动态阈值法提取物候参数。并将该方法简记为iNDVI-SG方法, 其中iNDVI表示NDVI累积时间曲线。对基于上述NDVI年均临界值剔除后的数据进行S-G滤波时间序列重建, 具体步骤为:首先在曲线重构过程中, 采用Spike Method来剔除NDVI曲线的无效点; 然后在原始的NDVI数据上构建窗口, 进行加权平均, 并将得到的拟合值代替原始数据; 最终通过对原始NDVI时序数据窗口平滑实现对整条曲线的拟合。经反复试验最后设置窗口大小为3, 有理多项式次数为3。要说明的是, 由于Timesat软件中在n年中只会提取每个物候参数的n-1个特征时期, 为了克服这个问题, 本研究在研究时间的开始和结束时各添加半年的虚拟数据。滤波前后数据比较见图2。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 S-G滤波时间序列重建前后数据对照Fig.2 Comparison of data before and after reconstruction of S-G filter time series

动态阈值法的提出者Jö nsson等^[17]以及大多数研究者^{[41, 42]}, 建议植被生长季开始日期阈值为20%左右, 但是对草地植被生长季结束期提取时的阈值设定各有不同^{[21, 43, 44]}。根据现有实地记录的典型物候相关信息并经过多次试验的基础上, 本研究最后将动态阈值生长季开始时间和长季结束时间阈值分别设定为20%和50%。

基于一元线性回归分析计算3种草地物候参数变化率G_slope, 并用pearson相关系数R衡量定距变量间的线性关系。在分析草地物候变化幅度和显著性时:根据G_slope的大小判断物候变化程度, 当G_slope< 0时, 物候趋势变化为提前或缩短趋势, 当G_slope> 0时, 为推迟或延长趋势; 再由R分析物候和时间的相关性, 利用t检验对相关系数R的显著性进行分析, 令P< 0.1时, 为变化趋势显著。在对提取物候结果精度分析时, 采用pearson相关系数R分析遥感提取物候参数和站点实测物候参数的相关性。

为了准确评价草地植被物候遥感提取结果, 利用站点实测数据对草地植被物候遥感提取结果进行了精度检验(图3)。草地植被SOS的均方根误差(root-mean-square error, RMSE)为9.03 d, 而偏差(bias)为2.05。草地植被EOS的提取结果与站点实测数据相关系数为0.756, RMSE为9.44 d, bias为0.65。两者均通过99.9%的显著性检验, 为极显著相关关系。除个别观测数据外, SOS和EOS结果与观测数据误差大都控制在10 d以内。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 物候监测结果与站点实测数据的比较Fig.3 Comparison of phonological monitoring results with site measured data

对比本研究与相关研究同地区的草地物候提取结果(表1), 发现本研究的物候提取信息与文献研究同地区草地物候参数的平均值相差不大。虽与高亚敏等^[35]和侯学会^[36]对通辽草地研究和国志兴等^[37]对中国东北草地研究的EOS结果有较大差别, 但本研究提取的生长季末期均值介于两种研究结果中间, 且通辽与东北地区相邻, 因此结果处于合理的范围之内。考虑到使用的原始数据为15 d合成的, 以及与相关数据来源和研究方法的不同, 故这一误差仍在可接受范围之内。验证结果说明本研究应用的物候提取方法对中国草地植被进行物候的遥感监测结果是可靠的。

表1

Table 1

表1(Table 1)

表1 物候监测结果与其他研究结果的比较 Table 1 Comparison of phenological monitoring results with other studies 研究区域 Research area 文献 Document (SOS, d) 本研究 This research (SOS, d) 文献 Document (EOS, d) 本研究 This research (EOS, d) 数据 Data resources 研究时段 Study period (yr) 文献来源 References 通辽Tongliao city 105~130 100~120 235~265 260~280 物候观测数据Phenological observation data 1981-2016 [35] 中国北方Northern China 120~150 100~130 250~280 260~290 GIMMS NDVI& SPOT-VGT NDVI 1982-2012 [36] 中国东北 Northeast China 110~140 100~130 280~300 270~280 GIMMS NDVI 1982-2003 [37] 青海湖流域 Qinghai Lake watershed 110~150 120~140 240~300 260~280 MODIS EVI 2001-2012 [38] 青藏高原Tibetan Plateau 115~125 110~150 260~300 260~300 GIMMS NDVI 3g 1982-2013 [39] 黄河源区Source region of the Yellow River 126~146 110~140 277~290 260~290 MODIS NDVI 2000-2012 [40]

表1 物候监测结果与其他研究结果的比较 Table 1 Comparison of phenological monitoring results with other studies

1986-2015年中国天然草地SOS、EOS和LOS的年均空间分布(图4)。我国草地SOS主要集中在第100~140天, 即4月下旬到5月下旬(图4a)。北疆和大兴安岭以西的内蒙古地区草地SOS较早在第90~100天, SOS较晚的地区主要在青藏高原地区在第120天之后, 而我国东南部的草地, 呈明显的由南向北逐渐延后的趋势。从图4b可以看出, 中国草地的EOS较为集中在一年中第260~290天, 即从9月下旬到10月上旬。对于青藏高原地区, 冈底斯山和喜马拉雅山脉东段草地EOS出现的时间最晚, 在我国的东南部EOS大致呈自南向北递减趋势。图4c显示出我国草地生长季持续时间集中在130~170 d, 但存在较高的空间分异, 持续时间较短的地区主要分布在昆仑山脉、横断山脉和青海高原地区大致在130 d左右, 北疆大部分的草地和大兴安岭以西的内蒙古地区由于物候开始的时间比较早, 生长季持续时间也较长。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 中国草地物候多年均值的空间分布Fig.4 Spatial distribution of multi-year mean value of grassland phenology in China

不同类型草地的SOS、EOS和LOS也存在较大差异(图5)。对于SOS, 高山亚高山草甸和高山亚高山平地草原较晚, 时间大约为第120~130天, 所占面积分别约为58%和48%, 由于两种草地面积所占较大, 所以在全国范围内的SOS在第120~130天这一时段内贡献占较大比重。按物候开始时间早晚排序, 山地草原< 平地草地< 草地< 荒漠草原< 高山亚高山草甸< 高山亚高山平地草原; 就EOS而言, 各草地类型差异不明显, 按物候结束时间排列, 荒漠草原< 草地< 平地草地< 高山亚高山平地草原< 高山亚高山草甸< 山地草原; LOS受各类草地的SOS和EOS的影响也呈现不同差异, 分布于低海拔地区的山地草地和草地的草地植被生长季持续较长在第150天以后的分别约占90%和80%。按物候生长季长度排列, 荒漠草原< 高山亚高山平地草原< 高山亚高山草甸< 平地草地< 草地< 山地草原。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 6种草地类型的多年平均物候各时段所占像元的百分比 Ⅰ :高山亚高山草甸Alpine and sub_alpine meadow; Ⅱ :山地草原Slope grassland; Ⅲ :平地草原Plain grassland; Ⅳ :荒漠草原Desert grassland; Ⅴ :草甸Meadow; Ⅵ :高山亚高山平地草原Alpine and sub-alpine plain grassland.下同The same below.Fig.5 Percentage of pixels in 6 grassland types per year average phenological period

2.3.1 中国草地物候不同时段变化趋势的空间分布 参照前人研究结论^{[16, 26, 27]}, 为研究中国天然草地植被物候期变化规律, 选取1998年为分界, 对1998年前后草地SOS、EOS和LOS进行逐像元回归分析(图6)。近30年内, 中国大部分草地生长季SOS变化幅度为-0.3~0.3 d· yr^-1, 变化趋势不显著(P> 0.1); SOS变化趋势在1998年前后存在较大转折, 在转折之前, SOS平均提前0.5~1.0 d· yr^-1的草地所占面积较大, 特别是在内蒙古东北部、北疆和青藏高原大部分地区呈显著提前趋势; 转折之后SOS空间变化趋势格局发生转变, 除部分地区草地略有提前外整体呈不明显的推迟趋势。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 中国草地在不同时段物候变化趋势的空间分布Fig.6 Spatial distribution of phenological trends in different periods of grassland in China

从草地EOS变化趋势来看, 30年间全国草地EOS的变化特征除内蒙古高原和鄂尔多斯高原推迟趋势显著提前外, 普遍呈推迟和提前不明显的趋势; 1998年之前全国65%的草地EOS呈推迟趋势; 1998年之后出现转折除昆仑山脉EOS明显推迟外其他地区变化趋势不明显。

30年间青藏高原地区LOS主要为不明显的缩短趋势, 内蒙古高原、鄂尔多斯高原和北疆表现为增长趋势, 主要变化时间集中在0.5~1.0 d· yr^-1; 1986-1998年大部分地区草地LOS呈增长趋势, 为大于0.3 d· yr^-1的变化趋势; 1999-2015年全国草地LOS变化为不明显的提前趋势。

2.3.2 6种草地类型不同时段的物候变化趋势 不同草地类型由于其生境和生物学特性的不同, 物候变化趋势也呈现出一定的差异性。图7显示了3个时段6种草地类型的平均物候变化速率。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 各草地类型不同时段的物候变化趋势 A:1986-2015 SOS; B:1986-2015 EOS; C:1986-2015 LOS; D:1986-1998 SOS; E:1986-1998 EOS; F:1986-1998 LOS; G:1999-2015 SOS; H:1999-2015 EOS; I:1999-2015 LOS.红色加框表示趋势显著(P< 0.1)。Fig.7 Phenological change trends of different grassland types at different time periods red boxed indicates significant trend (P< 0.1).

在30年时间尺度上, SOS、EOS和LOS大多数类型草地物候特征未有明显变化。尽管如此, 高山亚高山草甸和草甸的物候指标趋势显著性水平达90%。特别地, 具有显著突变的草地类型如:高山亚高山草甸和高山亚高山平地草原SOS在1998年前呈现显著提前趋势, 分别为0.5、0.8 d· yr^-1, 之后呈显著推迟趋势; 高山亚高山草甸的EOS在1998年有0.1~0.2 d· yr^-1的推迟现象, 1998年后提前幅度达到0.1~0.2 d· yr^-1; 对于LOS高山亚高山草甸、荒漠草原和高山亚高山平地草原在1998年前的时段, 有大于0.3的增长趋势, 而1998年后则有大于0.3的提前趋势。该结果也印证了遥感物候观测中不同物种在全球变化下表现出不同程度的变化趋势^[45]。