Long-term temporal characteristics of vegetation clumping index remote sensing products at mid/high latitudes in the Northern Hemisphere – a case study
-
摘要:
植被聚集指数(clumping index,CI( Ω ))是表征植被冠层聚集程度的重要结构参数,由于其定量化研究起步较晚,导致对CI季相变化特征的研究不充分,结论争议较大.为此,本文基于长时间序列的MODIS CI产品,从北半球中高纬度植被物候特征敏感区,在13个国际地圈-生物圈计划(IGBP)类型中,优选了84个高质量的代表性像元,开展典型像元CI季相变化特征的案例研究.以归一化植被指数(normalized difference vegetation index,NDVI)作为对比,提出改进的动态阈值法,结合离散Fourier变换方法,分别提取不同地类的生长季开始时间(start of season,SOS)及生长季结束时间(end of season,EOS),最终建立北半球中高纬度各地类生长季与休眠期的经验Ω.结果表明:CI具有较为明显的物候变化规律及季节变化特征,甚至能够识别出耕地的一年两熟迹象,但相对于NDVI相对稳定的季相变化特征,大部分地类的CI表现出较大的变化和不确定性,其中,SOS和EOS多分别在第100和第300天左右变化,生长季则多维持在200 d左右;提取物候特征参数的最佳阈值随提取时期、地物类别的变化而变化,其中提取SOS和EOS的最佳阈值多集中在40%~80%和80%~90%;经验Ω呈现出针叶林的聚集效应最强,耕地的聚集效应最弱的特征.本研究对于揭示不同地类CI季相特征及相关应用研究提供了有用的证据和参考.
Abstract:Clumping index (CI) is a vital structural parameter characterizing non-random spatial distribution of vegetation canopy leaves. Quantitative research on CI started relatively late compared to other vegetation structural parameters. This is especially true for the study of seasonal variations in CI, resulting in significant controversies. In this work, we conducted a case study on seasonal characteristics of typical CI pixels in mid/high latitudes vegetation phenology-sensitive regions in the Northern Hemisphere, using recently developed long-term MODIS CI products developed by us. We selected 84 high-quality representative pixels from 13 land cover types designated by the International Geosphere-Biosphere Program (IGBP). We applied an enhanced dynamic threshold method in conjunction with discrete Fourier transform to extract start of season (SOS) and end of season (EOS) for various land cover types, using normalized difference vegetation index (NDVI) as reference. We established empirical Ω values for both growing season and dormant period across various landcover types. CI exhibits distinct phenological and seasonal variations, to even detect biannual cropping signs in cropland. In comparison with the relatively stable seasonal variations of NDVI, CI exhibits somewhat large variability and uncertainty in most land cover types. SOS and EOS vary within approximately 100 and 300 days, respectively; the growing season typically spans around 200 days. The optimal threshold values for extracting phenological parameters vary dependent on extraction period and landcover type. The thresholds for SOS vary within the range from 40%-80%, but for EOS, from 80%-90%. Empirical Ω values show the strongest clumping effect in coniferous forests and the weakest in croplands.
-
Keywords:
- clumping index /
- time series /
- vegetation phenology /
- MODIS
-
0 引言
天然植被冠层的叶片并不是随机分布的,而是在树群、树冠以及叶簇等不同层次上呈现不同程度的聚集[1].植被聚集指数(clumping index,CI(Ω))作为表征植被冠层叶片在空间分布非随机程度的一个重要结构参数[2−3],量化了植被冠层叶片相对于随机分布的偏离程度[3].当叶片呈现随机分布时,Ω=1;当叶片呈现规则分布时,Ω>1;当叶片呈现聚集分布时,Ω<1,自然界中的植被叶片绝大部分都呈现聚集分布,其值随聚集效应的增强而减小.Ω被Chen等[4]定义为有效叶面积指数(effective leaf area index,LAIe)与真实叶面积指数(leaf area index,LAI)之比(即,Ω=LAIe/LAI),因而可用于提高叶面积指数反演的精度[5].在以下文本中,CI和Ω可以互换使用,分别表示更具描述性和算术意义的含义。CI和LAI是区分冠层光照和阴影叶片的基础,因而CI对于提高陆地碳循环模型的精度有着重要意义[6].例如,CI可用于提高全球总初级生产力(GPP)的估算精度.Chen等[7]研究表明,当直接使用LAI而忽略聚集效应时,全球GPP被高估12%,当直接使用LAIe且不考虑聚集效应时,全球GPP被低估9%.叶片的聚集效应同样影响着植被冠层内光能、降水的截获与分配,进而影响植被冠层的光合作用[8−10]及蒸腾作用[11].Lacaze [10]研究表明,当考虑叶片的聚集效应时,其黑云杉冠层的日光合作用的估测结果同不考虑聚集效应时的结果相差20%.总之,在区域和全球范围内,CI对光能传播、地表水文、光合作用、植被生产力等生态过程的研究具有重要作用.
目前获取CI的方法主要分为光学仪器实地测量、利用激光雷达技术提取、基于多角度星载数据进行反演.实地测量一般选用冠层结构和辐射分析仪(tracing radiation and architecture of canopies,TRAC)、多波段植被摄像仪(multiband vegetation imager,MVI)以及数字半球摄影(digital hemispheric photography,DHP)技术.但均需要较为严格的天气条件,其中TRAC需在晴空无云的天气状况下进行[12],MVI的光照条件需为阴天[13],DHP需在漫射光下获取影像并设置最佳曝光时间[14].Chen等[12]提出基于间隙率和间隙大小反演CI;Walter等[15]基于DHP技术,利用空间分段系数反演CI.激光雷达技术受天气影响较小,可快速获取植被三维结构,目前已被广泛应用于CI的反演中.Thomas等[16]基于机载激光雷达技术,利用地面实测CI和机载点云数据建立了一种新的计算CI的经验模型.此外,目前现有的CI产品主要是基于多角度星载数据进行估算,例如POLDER 数据[17−18]、MODIS 数据[19−20]以及多角度成像光谱仪(MISR)数据[21]等,其中Chen等[17]基于4尺度几何光学模型,模拟得到了不同情况下的归一化热点和暗点差异指数(NDHD)同CI的经验关系,随后该模型被广泛应用于CI产品的生产中.Chen等[17]利用POLDER 数据的红光和近红光波段生产了 6 km 分辨率的全球 8 个月的 CI 产品;He等[19]基于MODIS数据生产了空间分辨率为500 m的2006年全球CI产品,但每年仅给出了1个代表全年CI中位数的值,难以用于CI时相变化特征的分析;Jiao等[22]以NDHD-CI经验回归关系式为基础,基于AFX及BRDF原型库发展了备用算法对CI反演结果的异常值进行约束,生产了空间分辨率为500 m的每月和8 天(d)的全球MODIS CI产品,同时提高了该产品的时间及空间分辨率,为相关陆表及生态模型的应用、CI时空变化特性的研究提供了有益参考及数据支撑.
之前,由于缺乏对CI时空变化特性的理解,CI被认为是固定不变的[23−25],使得大多数CI反演算法及其产品应用尚未充分考虑其季相变化,对此,国内外部分学者开始关注CI的季相特征并展开了研究.部分研究针对不同植被类型的CI季相特征展开了研究:Ryu等[26]研究表明在温带落叶阔叶林中,CI及LAI均呈现出明显的季节效应;Sprintsin等[27]及Noe等[28]研究表明,即便是常绿针叶林,由于物候特性的驱动,其聚集程度也会随季节变化;Fang等[29]基于DHP法测定作物的CI,表明其具有明显的季节变化规律,CI值随作物的生长而减小,在衰老后期增加;Wei等[30]在全球尺度分析了CI的时空特性,研究表明落叶阔叶林和混交林等较为均一的作物,其季相特征强于其它植被类型.此外,部分学者还针对CI在区域间的季相变化展开了研究:He等[19]基于生产的全球CI产品,评估了某站点CI的季节变化,结果表明并无明显的季节波动现象;Pisek等[21]基于MISR数据反演得到区域CI产品,并首次对CI进行平滑,从而验证了该CI产品在季节变化时间序列里的适用性;朱高龙[31]利用MODIS BRDF产品反演得到了中国区域每8 d的500 m分辨率的CI产品,揭示了中国区域的CI具有明显的季节变化特征.Hu等[32]探究了三江平原CI时空变化特征及其驱动因素,表明三江平原的CI年际变化规律不明显,但季节变化规律显著;Yin等[33]通过对比基于不同版本MCD43产品反演得来的全球CI产品的时空变化,发现2类CI数据在全球不同纬度区域内的空间及季节变化相似,月CI数据更适于表征地表CI的季节变化.但前人针对CI的相关研究大多局限在某一区域平均CI的年内、年际变化的波动规律上,从提取物候特征参数方面定量分析CI的季相变化并建立其经验模型的相关研究还未见报道.
本文基于作者所在研究团队最近研发生产的2001—2020年全球500 m月合成CI产品,依据IGBP分类体系,在北半球中高纬度(23.5°~60.0°N)各地类(包括一年两熟耕地)中分别选取6个高质量典型像元进行CI季相变化特征的案例研究,既可以得到植被物候相对敏感的中高纬度区域中有代表性的CI季相变化规律,又可以暂时不考虑全球尺度上近20 a低质量像元CI反演结果的不确定性.研究中以较为稳定的NDVI作为对比参考,基于离散Fourier变换分别重建其时间序列,并使用本文改进的动态阈值法分别反演其关键物候参数,包括生长季开始时间(start of season,SOS)及生长季结束时间(end of season,EOS),利用MODIS物候产品MCD12Q2对反演结果进行间接验证,从定量和定性2个角度深入研究CI作为植被结构参数所表现的季相特征,并基于提取物候参数的最佳阈值建立北半球中高纬度各地类生长季与休眠期的经验Ω,以期更好地满足生态、气象及其他地表过程建模的需求.
1 数据与方法
1.1 数据收集
1)MODIS CI数据.本研究选用作者所在的研究团队研发生产的2001—2020年MODIS全球500 m月合成CI产品,由MODIS BRDF和反照率产品(MCD43A1和MCD43A2)及MODIS地表覆盖类产品(MCD12Q1)通过主算法和备用算法计算得出,并将冰雪(IGBP15)、无植被覆盖(IGBP16)、水体(IGBP17)等没有植被覆盖的区域进行掩膜处理[22].数据有效范围为3300~10000,尺度因子为0.0001,有效值为0.33~1.00,该数据集来源于http://www.geodata.cn.
2)NDVI数据.选用MOD13A1产品提供的NDVI数据集,该产品包含2个波段,分别为经过大气校正的双向表面反射率计算出来的NDVI及EVI数据,提供每16 d为周期、将全年划分为23个时段的500 m空间分辨率的3级正弦投影产品,其曲线可以反映植被的光谱及物候特征[34−35], 使用最大值合成法合成月尺度NDVI数据.
3)地表覆盖类型数据.选用美国NASA提供的2001—2020年的MCD12Q1数据集,该数据每年1期,覆盖全球,空间分辨率为500 m,时间分辨率为1 a,包括5种不同的土地覆盖分类体系[36],选用国际地圈-生物圈计划(IGBP)分类体系,对其进行镶嵌、波段提取等操作处理.
此外,本研究还选用地球系统科学数据共享平台提供的全球地表覆盖数据GlobalLand30作为本文选取典型像元的辅佐验证数据,用于进一步减小CI所在地类的不确定性,该数据空间分辨率为30 m,使用包含耕地、林地、草地等10个一级类型在内的分类系统,目前已发表2000、2010、2020版.
4)MODIS物候数据.选用美国NASA提供的2001—2020年MODIS土地覆盖动态产品(MCD12Q2),是MODIS全球陆地表面动态产品(MLCD),该产品每年1期,空间分辨率为500 m,采用增强型植被指数(EVI),利用Logistic函数法反演植被物候期[37].
1.2 研究方法
本研究以植被类型较为丰富、植被物候相对敏感的北半球中高纬度带作为研究区域,旨在通过在典型地物类中优选代表性像元,并提取关键物候参数,深入探讨研究区内CI作为植被结构参数所表现的季相特征,研究框架如图1所示.
1)典型像元选取规则.为深入探究CI的季相变化特征,本文以MCD12Q1数据集为主,辅以GlobalLand30数据集,对研究区内不同地类各选取6个典型像元,其中北美及东亚地区各地类均选取3个典型像元,共选取84个样本,分布在23.5°~60.0°N,包括基于MCD12Q2物候数据区分的一年一熟及一年两熟耕地,以进一步探究CI捕捉一年两熟植被生长状况的能力.通过对所选典型像元的样本量进行假设检验,发现在95%的置信水平下,所选研究区内的样本量能够满足本研究的需求.典型像元选取规则为:①2001—2020年MCD12Q1地表覆盖类型未变化;②2000、2010、2020这3年间的GlobalLand30同类型像元比例均占到MCD12Q1像元的70%及以上,从而确保本文选取典型像元的代表性,尽量避免因空间分辨率较大而引起的混合像元异质性导致季相特征不明显的问题.同时,为尽量减少噪声,降低近20 a长时间周期所可能带来CI反演结果的不确定性,统计典型像元周围3×3像元的CI及NDVI均值并对其进行研究分析,图2为典型像元的分布情况.
2)离散傅里叶变换.傅里叶变换将时间序列数据看作一系列不同频率的正弦和余弦函数的拟合.由于植被的生长在时间上呈现出周期性变化规律,因而傅里叶变换被广泛应用于物候信息提取、植被指数时间序列重构以及农作物种植信息提取等领域[38−41].本文研究数据较为离散且数据有限,受限于反演过程,CI时间序列噪声较多,而离散傅里叶变换作为有限长时间序列的算法在理论上十分重要,并且可以更好地去除噪声的影响,因而选用离散傅里叶变换对NDVI及CI的时间序列进行快速有效的拟合,其一般模型为
$$ F(n) = \frac{1}{N}{\sum\limits_{n = 0}^{N - 1} {g(t){\mathrm{e}}} ^{ - j\tfrac{{2{\text{π}} }}{N}nt}}, $$ 式中:$ F(n) $表示傅里叶变换;$ g(t) $为时间序列,$ t $为时域中的序列号;$ N $为序列长度;$ n $为频率域的序列号;$ n/N $为对应正弦函数的频率.$ F(n) $为一个复数,其中$ F(0) $等于$ g(t) $的平均值;$ {F_{{\mathrm{amplitude}}}} = {(F_{{\mathrm{real}}}^2 + F_{{\mathrm{imaginary}}}^2)^{1/2}} $为对应分量正弦函数的幅度(amplitude);${F_{{\mathrm{phase}}}} $ = $ \arctan (F_{{\mathrm{imaginary}}}^2/F_{{\mathrm{real}}}^2) $为对应分量正弦函数的相位(phase).
3)改进动态阈值法.动态阈值法采用动态比值的形式,与时间序列的季节变化幅度紧密相关,能够有效降低物候反演方法不确定性的影响[42],即当CI时间序列曲线下降或上升至当年CI振幅一定百分比的时刻,即为生长季开始或结束时间,NDVI则与之相反.但由于CI作为一种植被结构参数,其反演精度受到热点反射率[19−20]和混合像元异质性[17]等因素的影响,即使在经过平滑处理后,其时间序列曲线较NDVI仍存在部分噪声.故本研究基于CI时间序列曲线的特点,在原始动态阈值法的基础上加以改进(图3).对于CI波动幅度小于同侧振幅0.1的噪声,若按照原始动态阈值法,返青期应为SOS1,而改进的动态阈值法在检测到噪声后,继续提取到SOS2,并将生长季开始时间定义为SOS1与SOS2的平均值,从而达到进一步平滑噪声的目的.此外,动态阈值法的提出者Jonsson等[43]建议植被生长季开始日期阈值为20%左右,但阈值的选取依据植被指数、地表覆盖类型以及地理位置的不同而有所不同,故本文根据CI时相波动较大的特点,放宽域值范围,使用改进的动态阈值法,分别基于10%~90%的阈值提取物候参数,辅以MCD12Q2产品进行验证,确定CI提取物候参数时的最佳阈值,以期通过物候特征参数更好地反映CI的季相特征,并确定各地类的经验Ω.对于一年两熟的耕地种植区域,其CI及NDVI时间序列曲线应当含有2个波谷或波峰,同样使用10%~90%的阈值分别尝试提取2个生长季的物候特征参数.
4)经验Ω的建立.CI时间序列曲线除异常波动较多外,还有部分像元的部分年份无周期性变化规律,导致物候特征参数提取失败.为尽量排除噪声干扰,探究基于物候特征参数的CI季相特征,本文首先基于改进的动态阈值法进行筛选,定义若有超过30%的年份无法提取其物候特征参数,则该像元被定义为无明显季节变化规律的像元,予以剔除.并使用精度评价指标均方根误差(RMSE)及均值偏差(Bias)对提取结果进行精度评估.基于RMSE得到提取CI物候特征参数的最佳阈值,对2001—2020年各地类典型像元的CI时间序列进行划分,规定当其下降至SOS最佳阈值对应的Ω时,生长季开始,休眠期结束;当其上升至EOS最佳阈值对应的Ω时,生长季结束,休眠期开始.最后,分别对2001—2020年的CI时间序列进行处理,依据最佳阈值划定每年的生长季与休眠期范围,分别对范围内的Ω取平均,得到各典型像元20 a内每年生长季与休眠期的平均Ω,并分别求取其变异系数(coefficient of variation,CV(CV)),其含义为序列的标准差除以算术平均值.可将CV划分为5个等级:CV≤0.05为低变异;0.05<CV≤0.10为较低变异;0.10<CV≤0.15为中等变异;0.15<CV≤0.20为较高变异;CV>0.20为高变异.CV值越小,表明数据分布越集中,具有更好的稳定性[44].根据CV值的大小,决定是否对20 a内每年生长季与休眠期的平均Ω求取多年平均,进而建立不同地类的经验Ω.
2 结果与分析
2.1 研究区内典型像元CI时间序列曲线谷值统计
图4为2001—2020年间CI及NDVI时间序列在经过Fourier重建后的谷值及峰值数目及其对应月份,对于不超过总振幅20%的噪声,认为是由非植被生长规律导致的波动,其谷值或峰值的数目不予统计.由图4-a、b可得,大部分典型像元的CI谷值数目在20上下浮动,表明CI在20 a中存在20个物候变化周期;少数为0或过高,且大部分地类为CL2像元的谷值数在40上下浮动,表明一年两熟耕地的CI在20 a中存在40个物候变化周期,符合一年两熟作物的物候变化规律.此外,其谷值对应月份大多分布在4—9月,但较为分散,一年两熟的耕地像元第一生长季谷值月份位于2—5月份,第二生长季谷值月份位于7—9月份.其中北美地区典型像元的谷值月份分布迹象优于东亚地区,这可能是由于所选东亚地区作物种植种类在500 m像元尺度上较为分散破碎所导致的.由图4-c、d可得,绝大部分一年两熟耕地像元的NDVI峰值数目为40,除地类为OSh的部分像元外,其余像元的峰值数目均为20,且峰值月份基本位于6—7月,相较于CI时间序列,其分布更为集中、有序.由此可得,通过对典型像元的CI及NDVI时间序列曲线进行提取分析,表明CI时间序列曲线并不十分稳定,虽能捕捉其周期性变化规律,并且能够识别出耕地的一年两熟迹象,但同NDVI相比,其反映植被确切生长规律的优势并不明显.除此之外,统计可得CI时间序列平均谷值月份为6.14个月,NDVI时间序列平均峰值月份为6.87个月,说明代表聚集程度的CI最低值月份较代表绿度的NDVI最大值月份提前了0.73个月,这可能与植被生长期间,叶片中代表绿度的叶绿素浓度逐渐增大,而植被冠层叶片则更倾向于提前达到聚集状态有关[45−46].
2.2 基于CI及NDVI的植被物候特征参数提取
图5和图6分别为基于CI及NDVI提取的典型像元SOS及EOS同MODIS物候产品MCD12Q2作对比验证的精度评价结果,精度评价指标为RMSE和Bias,各地类的最小误差用圆圈标出,其中Bias最小误差的高估现象用红圈标出,低估现象用蓝圈标出.由图5-a、b可得:除PWe及CL1外,基于CI最佳阈值提取SOS的RMSE均<40,其最佳阈值多集中在40%及80%;对于EOS而言,除Wsa、PWe、CL1以及CVM外,其最佳阈值对应的RMSE均<35,最佳阈值多集中在80%~90%.由图5-c、d可得:基于CI最佳阈值提取的SOS,高估与低估现象相当,除CL2-1外,其余地类典型像元的Bias绝对值均<15,在此,CL2-1、CL2-2分别代表一年两熟耕地的第一及第二生长季,见图5图例;对于EOS而言,多呈现为低估现象,说明基于CI提取得到的EOS多呈现日期提前的现象,除CL1及CL2-1外,其余地类典型像元的Bias绝对值均<15.由图6-a、b可得:除Osh、PWe外,基于NDVI最佳阈值提取SOS的RMSE均<22,最佳阈值多集中在40%~70%;除OSh、GL、PWe及CL2-1外,基于NDVI最佳阈值提取EOS的RMSE均<28,最佳阈值多集中在10%~40%.此外,由图6-c、d可得,除ENF、CVM外,基于NDVI最佳阈值提取SOS的Bias绝对值均<3.6,高估与低估现象相当;EOS的Bias绝对值稍大,除CL2-1与CL2-2外,其余均<9,其中多呈现为高估现象,说明基于NDVI最佳阈值提取得到的EOS多呈现日期延后的现象.总体而言,基于改进的动态阈值法提取NDVI物候特征参数提取精度要优于CI,这是由于CI在反演过程中受到多种因素的影响,并不十分稳定,但二者误差均在可接受范围内,说明CI具有同NDVI类似的较为明显的物候变化规律及季节变化特征.此外,基于CI及NDVI提取SOS、EOS时的最佳阈值呈现较大差异,说明提取物候参数的最佳阈值不仅与提取时期有关,且随植被指数、地物类别的变化而变化.
2.3 最佳阈值提取的SOS及EOS时间分布特征
图7为不同地类各选取一个典型像元,其对应2001—2020年MCD12Q2产品的SOS、EOS,以及基于RMSE得到的最佳阈值提取CI和NDVI的SOS、EOS时间分布特征,可以发现三者相差不大,但基于CI得到的多年物候年均值的标准差较大,再次表明CI作为一种植被结构参数并不十分稳定,其物候提取日期年际差异较大.此外,DNF、DBF、Sav、CL1以及CL2-2的SOS>100 d,且EOS较其它地类较早,生长季(EOS-SOS)较短,均<180 d.CL2-2的SOS开始的最晚,CL2-1的EOS结束得最早,二者生长季时间相当,在100 d上下浮动.其余各地类的多年物候年平均值则相差不大,SOS多<100 d,EOS在300 d上下浮动,生长季则多维持在200 d左右.综上所述,CI时间序列曲线虽然能够反映较为明显的物候变化规律,但SOS及EOS年际间的差异相对较大,基于CI提取植被物候特征参数的相关研究还有待进一步深化.
2.4 基于最佳阈值提取不同地类生长季与休眠期的经验Ω
表1为提取不同地类CI物候特征参数的最佳阈值及其经验Ω,计算可得,各地类C V值均<0.06,其中有超过93%的典型像元为低变异状态,其余均为较低变异状态,说明20 a时间里,生长季与休眠期的Ω波动程度均较低,数据均较为稳定,故对其取平均,得到不同地类生长季与休眠期的经验Ω及其标准差.由表1可得,各地类的生长季经验Ω均小于休眠期经验Ω,说明在北半球中高纬地区,各地类的经验Ω均呈现出生长季时植被的聚集指数低、聚集效应强,休眠期时植被的聚集指数高、聚集效应差的规律.此外,由于针叶林的枝叶较其它植被类型更具有聚集作用,因此针叶林的经验Ω相对较低,均≤0.60,其中ENF的经验Ω最低,生长季与休眠期的经验Ω分别为0.48与0.54,DNF次之,分别为0.55与0.59;CL2-2、CL2-1、CL1以及CVM的经验Ω较高,均>0.75,且基本一致,生长季与休眠期的经验Ω分别为0.76~0.77和0.80~0.81;常绿林的生长季与休眠期的经验Ω也呈现出较为显著的季相变化,这主要是由于常绿林在冬半季会发生掉叶、掉针等现象,这与前人的研究结果相似[31].除OSh的生长季与休眠期的标准差分别为0.094与0.063外,其余地类标准差均<0.05,说明OSh的CI时间序列并不十分稳定,这可能与其聚集程度较低,分布较为稀疏,其反演过程受到限制有关.
表 1 提取不同地类CI物候特征参数最佳阈值及其经验ΩIGBP 地物类型 最佳阈值/% 经验Ω SOS EOS 生长季 休眠期 1 ENF 80 80 0.48±0.011 0.54±0.014 2 EBF 40 90 0.57±0.036 0.68±0.028 3 DNF 40 80 0.55±0.026 0.59±0.025 4 DBF 40 90 0.69±0.016 0.75±0.013 5 MF 50 90 0.63±0.016 0.68±0.012 6 CSh 40 90 0.68±0.028 0.77±0.049 7 OSh 20 80 0.74±0.094 0.81±0.063 8 Wsa 80 90 0.65±0.033 0.72±0.029 9 Sav 30 70 0.63±0.034 0.75±0.023 10 GL 50 60 0.72±0.013 0.81±0.017 11 PWe 80 40 0.74±0.030 0.78±0.027 12 CL1 20 90 0.76±0.023 0.80±0.014 12 CL2-1 80 90 0.77±0.035 0.81±0.020 12 CL2-2 40 90 0.77±0.028 0.80±0.028 14 CVM 50 90 0.76±0.032 0.80±0.024 3 结论
本文基于MODIS CI及MOD13A1等数据,使用离散Fourier变换分别重建CI及NDVI时间序列,并使用本文改进的动态阈值法分别反演不同地类的关键物候参数,将结果同MODIS物候产品MCD12Q2对比,并基于最佳阈值建立北半球中高纬度各地类生长季与休眠期的经验Ω,主要结论如下:
1) CI基本具有同NDVI类似的较为明显的物候变化规律及季节变化特征,其物候特征参数的最佳阈值随提取时期、地物类别等因素的变化而变化.其中:提取SOS的最佳阈值多集中在40%~80%,除PWe及CL1外,RMSE均<40,除CL2-1外,其余地类典型像元的Bias绝对值均<15;提取EOS的最佳阈值多集中在80%~90%,除Wsa、PWe、CL1以及CVM外,RMSE均<35,Bias多呈现为低估现象,除CL1及CL2-1外,其余地类典型像元的Bias绝对值均<15.
2) 除一年两熟耕地的CI谷值个数在40上下浮动外,其余地类的谷值个数多在20上下浮动,且谷值个月对应月份大致分布在4—9月,其平均谷值月份为6.14个月,较NDVI平均峰值月份提前了0.73个月,说明CI时间序列具有一定周期性,且能够识别出耕地的一年两熟迹象,但相较于NDVI而言较为分散,且植被叶片的聚集程度较绿度而言率先达到极值.此外,基于最佳阈值提取CI的SOS与EOS同NDVI及MCD12Q2产品相比,三者相差不大,但基于CI得到的多年物候年均值的标准差较大.其中:DNF、DBF、Sav、CL1及CL2-2的SOS>100 d,且EOS较早,生长季较短;其余各地类则相差不大,SOS多<100 d,EOS在300 d上下浮动,生长季则多维持在200 d左右.
3) 在北半球中高纬地区,各地类的生长季经验Ω均小于休眠期经验Ω,均呈现出生长季时植被的聚集效应强,休眠期时植被的聚集效应差的规律.针叶林的聚集效应最强,其中ENF的Ω最低,DNF次之;CL2-2 、CL2-1 、CL1及CVM的Ω较高,均>0.75,且基本一致,生长季与休眠期的经验Ω分别为0.76~0.77、0.80~0.81.
由于CI的定量化研究起步较晚,导致CI的季相变化特征研究结论争议较大,针对此类问题,本研究着重解决了目前对于CI存在争议的2个问题:CI是否具有物候变化特征;CI的物候变化具有什么样的特征.研究表明,受限于反演过程,CI捕捉植被确切生长规律的优势虽不如NDVI,但依旧具有较为明显的物候变化规律及季相特征。此外,部分学者[47−48]倾向于假设BRDF在时间上是缓慢变化的,甚至在其应用中假设其全年恒定;而本文所使用的MODIS CI产品为基于NDHD-CI经验关系生产的,由于NDHD使用热点及暗点方向的反射率计算产生,而热点和暗点是BRDF在主平面上的2个典型特征,因此本研究也间接证明了BRDF形状具有季节变化规律.同时,与对30~45°N内Ω均值进行初步探究的结果类似,在今后的研究中,将继续完善相关工作,并将结合本文研究成果生产MODIS CI季相产品,以满足相关地表过程建模的需求.同时,本研究以MCD12Q2产品为基准,确定了基于CI提取不同地类、不同物候特征参数的最佳阈值,并以此为基准建立了北半球中高纬地区不同地类生长季与休眠期的经验Ω,同前人研究成果及实际情况基本相符,进而证实了本文研究方法的可行性,为今后相关研究与应用提供了有益参考.
本文也存在一些不足:由于CI时间序列曲线噪声较多,本研究在综合考量不同时间分辨率的CI产品及参考相关文献[33]后,决定选用月CI产品,在一定程度上避免了噪声的影响,但同样掩盖了植被物候变化的细微信息,进而影响提取精度.因此,今后将进一步开展CI时序数据噪声消减方面的研究,以期针对更高时间分辨率的CI产品展开更精确的相关研究.此外,通过地面物候观测数据进行验证具有准确和客观的特点,作为直接验证方法,是结论可靠性的重要组成部分,但由于受到多种原因的限制,地面物候观测数据无法推广到大尺度区域,目前很难进行长时间序列连续的物候数据的获取.且无法得到与本研究经过地类限制等条件筛选得到的典型像元相匹配的地面物候观测数据,因此搜集地面物候观测数据进行精度验证也存在一定困难.但依然期望在未来可以将基于CI得到的结果同地面实测物候数据进行对比,并对全球CI产品全像元分区域进行时相分析,以期更充分地了解和验证基于CI的大尺度物候观测结果,更好地满足相关研究的需求.
-
表 1 提取不同地类CI物候特征参数最佳阈值及其经验Ω
IGBP 地物类型 最佳阈值/% 经验Ω SOS EOS 生长季 休眠期 1 ENF 80 80 0.48±0.011 0.54±0.014 2 EBF 40 90 0.57±0.036 0.68±0.028 3 DNF 40 80 0.55±0.026 0.59±0.025 4 DBF 40 90 0.69±0.016 0.75±0.013 5 MF 50 90 0.63±0.016 0.68±0.012 6 CSh 40 90 0.68±0.028 0.77±0.049 7 OSh 20 80 0.74±0.094 0.81±0.063 8 Wsa 80 90 0.65±0.033 0.72±0.029 9 Sav 30 70 0.63±0.034 0.75±0.023 10 GL 50 60 0.72±0.013 0.81±0.017 11 PWe 80 40 0.74±0.030 0.78±0.027 12 CL1 20 90 0.76±0.023 0.80±0.014 12 CL2-1 80 90 0.77±0.035 0.81±0.020 12 CL2-2 40 90 0.77±0.028 0.80±0.028 14 CVM 50 90 0.76±0.032 0.80±0.024 -
[1] CHEN J M. Optically-based methods for measuring seasonal variation of leaf area index in boreal conifer stands[J]. Agricultural and Forest Meteorology,1996,80(2/3/4):135
[2] CHEN J M,BLACK T A. Foliage area and architecture of plant canopies from sunfleck size distributions[J]. Agricultural and Forest Meteorology,1992,60(3/4):249
[3] NILSON T. A theoretical analysis of the frequency of gaps in plant stands[J]. Agricultural Meteorology,1971,8:25
[4] CHEN J M,BLACK T A. Defining leaf area index for non‐flat leaves[J]. Plant Cell and Environment,1992,15:421 doi: 10.1111/j.1365-3040.1992.tb00992.x
[5] 麻庆苗,李静,刘强,等. 混合像元聚集指数研究及尺度分析[J]. 遥感学报,2012,16(5):895 [6] CHEN J M,LIU J,LEBLANC S G,et al. Multi-angular optical remote sensing for assessing vegetation structure and carbon absorption[J]. Remote Sensing of Environment,2003,84:516 doi: 10.1016/S0034-4257(02)00150-5
[7] CHEN J M,MO G,PISEK J,et al. Effects of foliage clumping on the estimation of global terrestrial gross primary productivity[J]. Global Biogeochemical Cycles,2012,26(1):1
[8] PALMROTH S,HARI P. Evaluation of the importance of acclimation of needle structure,photosynthesis,and respiration to available photosynthetically active radiation in a Scots pine canopy[J]. Canadian Journal of Forest Research,2001,31(7):1235 doi: 10.1139/x01-051
[9] CHEN J M,LIU J,CIHLAR J,et al. Daily canopy photosynthesis model through temporal and spatial scaling for remote sensing applications[J]. Ecological Modelling,1999,124(2/3):99
[10] LACAZE R. Retrieval of vegetation clumping index using hot spot signatures measured by POLDER instrument[J]. Remote Sensing of Environment,2002,79(1):84 doi: 10.1016/S0034-4257(01)00241-3
[11] LIU J,CHEN J M,CIHLAR J. Mapping evapotranspiration based on remote sensing:an application to Canada’s landmass[J]. Water Resources Research,2003,39(7):1189
[12] CHEN J M,CIHLAR J. Plant canopy gap-size analysis theory for improving optical measurements of leaf-area index[J]. Applied Optics,1995,34(27):6211 doi: 10.1364/AO.34.006211
[13] KUCHARIK C J,NORMAN J M,MURDOCK L M,et al. Characterizing canopy nonrandomness with a multiband vegetation imager (MVI)[J]. Journal of Geophysical Research:Atmospheres,1997,102(D24):29455 doi: 10.1029/97JD01175
[14] LI J C,FAN W J,LIU Y,et al. Estimating savanna clumping index using hemispherical photographs integrated with high resolution remote sensing images[J]. Remote Sensing,2017,9(1):52 doi: 10.3390/rs9010052
[15] WALTER J M N,FOURNIER R A,SOUDANI K,et al. Integrating clumping effects in forest canopy structure:an assessment through hemispherical photographs[J]. Canadian Journal of Remote Sensing,2003,29(3):388 doi: 10.5589/m03-011
[16] THOMAS V,NOLAND T,TREITZ P,et al. Leaf area and clumping indices for a boreal mixed-wood forest:lidar,hyperspectral,and Landsat models[J]. International Journal of Remote Sensing,2011,32(23):8271 doi: 10.1080/01431161.2010.533211
[17] CHEN J M,MENGES C H,LEBLANC S G. Global mapping of foliage clumping index using multi-angular satellite data[J]. Remote Sensing of Environment,2005,97(4):447 doi: 10.1016/j.rse.2005.05.003
[18] PISEK J,CHEN J M,LACAZE R,et al. Expanding global mapping of the foliage clumping index with multi-angular POLDER three measurements:evaluation and topographic compensation[J]. ISPRS Journal of Photogrammetry and Remote Sensing,2010,65(4):341 doi: 10.1016/j.isprsjprs.2010.03.002
[19] HE L,CHEN J M,PISEK J,et al. Global clumping index map derived from the MODIS BRDF product[J]. Remote Sensing of Environment,2012,119:118 doi: 10.1016/j.rse.2011.12.008
[20] ZHU G L,JU W M,CHEN J M,et al. Foliage clumping index over China’s landmass retrieved from the MODIS BRDF parameters product[J]. IEEE Transactions on Geoscience and Remote Sensing,2012,50(6):2122 doi: 10.1109/TGRS.2011.2172213
[21] PISEK J,RYU Y,SPRINTSIN M,et al. Retrieving vegetation clumping index from Multi-angle Imaging SpectroRadiometer (MISR) data at 275 m resolution[J]. Remote Sensing of Environment,2013,138:126 doi: 10.1016/j.rse.2013.07.014
[22] JIAO Z T,DONG Y D,SCHAAF C B,et al. An algorithm for the retrieval of the clumping index (CI) from the MODIS BRDF product using an adjusted version of the kernel-driven BRDF model[J]. Remote Sensing of Environment,2018,209:594 doi: 10.1016/j.rse.2018.02.041
[23] BALDOCCHI D D,WILSON K B,GU L H. How the environment,canopy structure and canopy physiological functioning influence carbon,water and energy fluxes of a temperate broad-leaved deciduous forest:an assessment with the biophysical model CANOAK[J]. Tree Physiology,2002,22(15/16):1065
[24] HOUBORG R,ANDERSON M C,NORMAN J M,et al. Intercomparison of a ‘bottom-up’ and ‘top-down’ modeling paradigm for estimating carbon and energy fluxes over a variety of vegetative regimes across the U. S.[J]. Agricultural and Forest Meteorology,2009,149(11):1875 doi: 10.1016/j.agrformet.2009.06.014
[25] SAMPSON D A,JANSSENS I A,CEULEMANS R. Under-story contributions to stand level GPP using the process model SECRETS[J]. Agricultural and Forest Meteorology,2006,139(1/2):94
[26] RYU Y,NILSON T,KOBAYASHI H,et al. On the correct estimation of effective leaf area index:does it reveal information on clumping effects?[J]. Agricultural and Forest Meteorology,2010,150(3):463 doi: 10.1016/j.agrformet.2010.01.009
[27] SPRINTSIN M,COHEN S,MASEYK K,et al. Long term and seasonal courses of leaf area index in a semi-arid forest plantation[J]. Agricultural and Forest Meteorology,2011,151(5):565 doi: 10.1016/j.agrformet.2011.01.001
[28] NOE S M,KIMMEL V,HÜVE K,et al. Ecosystem-scale biosphere-atmosphere interactions of a hemiboreal mixed forest stand at Järvselja,Estonia[J]. Forest Ecology and Management,2011,262(2):71 doi: 10.1016/j.foreco.2010.09.013
[29] FANG H L,LI W J,WEI S S,et al. Seasonal variation of leaf area index (LAI) over paddy rice fields in NE China:Intercomparison of destructive sampling,LAI-2200,digital hemispherical photography (DHP),and AccuPAR methods[J]. Agricultural and Forest Meteorology,2014,198:126
[30] WEI S S,FANG H L,SCHAAF C B,et al. Global 500 m clumping index product derived from MODIS BRDF data (2001-2017)[J]. Remote Sensing of Environment,2019,232:111296 doi: 10.1016/j.rse.2019.111296
[31] 朱高龙. 2000−2013年中国植被叶片聚集度系数时空变化特征[J]. 科学通报,2016,61(14):1595 [32] HU K H,ZHANG Z,FANG H L,et al. Spatio-temporal characteristics and driving factors of the foliage clumping index in the Sanjiang Plain from 2001 to 2015[J]. Remote Sensing,2021,13(14):2797 doi: 10.3390/rs13142797
[33] YIN S Y,JIAO Z T,DONG Y D,et al. Evaluation of the consistency of the vegetation clumping index retrieved from updated MODIS BRDF data[J]. Remote Sensing,2022,14(16):3997 doi: 10.3390/rs14163997
[34] HUETE A,DIDAN K,MIURA T,et al. Overview of the radiometric and biophysical performance of the MODIS vegetation indices[J]. Remote Sensing of Environment,2002,83(1):195
[35] ZHANG X Y,FRIEDL M,SCHAAF C,et al. Monitoring vegetation phenology using MODIS[J]. Remote Sensing of Environment,2003,84:471 doi: 10.1016/S0034-4257(02)00135-9
[36] 夏文韬,王莺,冯琦胜,等. 甘南地区MODIS土地覆盖产品精度评价[J]. 草业科学,2010,27(9):11 doi: 10.3969/j.issn.1001-0629.2010.09.003 [37] 夏传福,李静,柳钦火. 植被物候遥感监测研究进展[J]. 遥感学报,2013,17(1):1 [38] GEERKEN R,ZAITCHIK B,EVANS J P. Classifying rangeland vegetation type and coverage from NDVI time series using Fourier Filtered Cycle Similarity[J]. International Journal of Remote Sensing,2005,26(24):5535 doi: 10.1080/01431160500300297
[39] 郑玉坤,庄大方. 多时相AVHRR数据的傅立叶分析[J]. 中国科学院研究生院学报,2003,20(1):62 [40] SELLERS P J ,Randall D A,Collatz G J,et al. A revised land surface parameterization(SiB2) for atmospheric GCMs. Part I:model formulation[J]. Journal of climate,1996,9(4):676
[41] 马超,郭增长,张晓克. 一种NDVI时间序列离散Fourier重建方法[J]. 测绘科学技术学报,2011,28(5):347 doi: 10.3969/j.issn.1673-6338.2011.05.009 [42] 王钊,彭艳,权文婷,等. 秦岭森林物候时空分布特征及对水热条件的响应[J]. 干旱区地理,2019,42(5):1048 [43] JONSSON P,EKLUNDH L. Seasonality extraction by function fitting to time-series of satellite sensor data[J]. IEEE Transactions on Geoscience and Remote Sensing,2002,40(8):1824 doi: 10.1109/TGRS.2002.802519
[44] 刘洋洋,章钊颖,同琳静,等. 中国草地净初级生产力时空格局及其影响因素[J]. 生态学杂志,2020,39(2):349 [45] CROFT H,CHEN J M,ZHANG Y. Temporal disparity in leaf chlorophyll content and leaf area index across a growing season in a temperate deciduous forest[J]. International Journal of Applied Earth Observation and Geoinformation,2014,33(1):312
[46] LI J,LU X H,JU W,et al. Seasonal changes of leaf chlorophyll content as a proxy of photosynthetic capacity in winter wheat and paddy rice[J]. Ecological Indicators,2022,140:109018 doi: 10.1016/j.ecolind.2022.109018
[47] CLAVERIE M,VERMOTE E,FRANCH B,et al. Evaluation of medium spatial resolution BRDF-adjustment techniques using multi-angular SPOT4 (Take5) acquisitions[J]. Remote Sensing,2015,7(9):12057 doi: 10.3390/rs70912057
[48] VERMOTE E,JUSTICE C O,BREON F M. Towards a generalized approach for correction of the BRDF effect in MODIS directional reflectances[J]. IEEE Transactions on Geoscience and Remote Sensing,2009,47(3):898 doi: 10.1109/TGRS.2008.2005977