基于Sentinel-2A影像干旱区棉花叶片SPAD数字制图

唐普恩,丁建丽,3,*,葛翔宇,张振华

1 新疆大学资源与环境科学学院, 乌鲁木齐 8300462 新疆大学绿洲生态教育部重点实验室, 乌鲁木齐 8300463 新疆大学智慧城市与环境建模自治区普通高校重点实验室, 乌鲁木齐 830046

新疆作为中国最大的棉花种植基地,其棉花产量能够左右国际棉价,在我国棉花产业中占有举足轻重的作用。而棉花作为新疆的主要经济作物,对新疆农业经济和社会发展也至关重要[1]。植被生长过程中,叶绿素是参与光合作用的重要色素,其既体现着植物的生理状态,又决定着光合作用的强弱,是生态系统活力的重要体现[2]。叶绿素浓度直接影响植物代谢活动、生长周期、产量形成,是评价植物长势的重要指标[3- 5]。叶绿素含量的高低直接或间接的反映自然环境的胁迫作用,特别是干旱区独有的土壤高盐背景下盐胁迫会直接影响植物的生理生长指标和离子吸收,抑制植被光合吸收速率。植被叶片的叶绿素相对含量(SPAD)可直接用以衡量其实际叶绿素含量的高低,植物叶片的SPAD值易于与高光谱数据实现准确对应[6]。因此,快速易行的掌握叶绿素含量空间状况对干旱区制定科学的水肥一体化管理、潜在灾害胁迫 监测具有重大现实意义。

传统实验室测量叶绿素的方法复杂、成本高、效率低,且不具备实时、大面积监测能力,遥感技术的实施应用恰逢其时,利用遥感技术可以快速、大尺度、高效的捕获农作物信息,很好的弥补了传统手段的不足[7- 8]。国内外研究学者利用高光谱技术对植被叶绿素含量进行监测方面取得了很多成果。Lu等认为由红色和近红外反射率组合计算的可见光与植被的生物物理和生物化学特性有很好的相关性[9]。郭超凡等建立叶面尺度下不同包络线去除衍生转换光谱也叶绿素含量估算模型[10]。田明璐等使用多远逐步回归与线性回归构建多种光谱参数建立与叶绿素相对含量的反演模型中,发现前者精度优于后者[11]。冯伟等研究了小麦叶片色素密度与冠层高光谱参数的定量关系,认为群体叶片色素密度的敏感波段主要分布在可见光区, 而红边区域导数光谱表现更显著[12]。李成等[13]发现红边和近红外波段构建的植被指数的苹果树冠层叶绿素含量估算模型更准确,更稳定。已被证明光谱与叶绿素含量间存在机理,特别是在红边波段,由于作物在红光波段的强吸收以及在近红外波段的强反射,在680—760 nm的光谱范围内,红光波段反射率下降,近红外波段陡然上升,形成“吸收谷”,是绿色作物光谱最为明显的特征。红边波段在实践中光谱指数通常被视为遥感光谱数据与制备生理化学参数检测桥梁。这些光谱指数的设计,将不良影响降低的同时具备较好简洁性、鲁棒性和准确性。然而,由于植物叶片结构和观测角度等因素的影响,基于经验统计的方法可能表征较差。另外,基于物理的辐射传输机理模型根据叶片和土壤特性,冠层结构和阳光传感器的几何形状模拟冠层内的光传播,并提供对光与植被之间相互作用的有价值的物理理解,也被广泛用于植被重要的生理参数的反演(如叶绿素、氮含量、叶面积指数(LAI)、类胡萝卜素等)[14- 16]。李振海等比较了仅使用SAIL+PROSPECT(PROSAIL)冠层反射率模型与农学先验知识(APK)关联的结果,经分析证实了PROSAIL模型反演生物物理变量的操作潜力[17]。但其容易受地域、植被类型、生长季变化等因素影响,致使其用于遥感像元尺度反演时存有较大误差,且模型普适性较差[18- 19],但植被中的其他生理化学参数与叶绿素间具有很强的协同性和相关性,如叶片是植被光合作用的主要器官,叶面积大小和作用时间都与健康状态有关,而叶绿素又是植被生长状况和光合作用的重要指标。类胡萝卜素与也叶绿素在光谱反射中存在遮蔽影响。

在干旱区复杂的下垫面和盐分胁迫的背景下,单独一种方法或一种属性的协变量是无法满足精准农业的需求。基于此,本文以囊括3个红边波段的Sentinel- 2A为本底,联合最适光谱指数和基于PROSAIL物理模型的生物协变量集成估算棉花叶片SPAD值。以期获取高精度、高分辨率区域叶绿素分布,服务以精准农业为抓手的一带一路核心区农业信息化。

1 数据与方法

1.1 研究区概况

研究区位于新疆维吾尔自治区昌吉回族自治州(包含玛纳斯县,沙湾县和石河子市),地处天山北麓中段,准格尔盆地中部古尔班通古特沙漠南缘,地理坐标范围44°20′—45°20′ N,85°20′—86°40′ E。气候、地貌、水文条件、植被具有明显的垂直地带性分布[20- 21],海拔由南向北递减,玛纳斯河流域主要形成“山地-绿洲-荒漠”系统流域,绿洲主要以片状或带状分散在山前冲积平原和冲洪积扇扇缘[22]。年平均气温4—7℃,极端最高气温43℃,极端最低气温-42℃,多年平均降雨量150—200 mm,多年平均蒸发量达1500—2100 mm,属于典型的温带大陆性干旱半干旱气候。

1.1.1土壤特征

研究区内的土壤类型十分丰富主要有灰漠土、棕钙土、风沙土、沼泽土、草甸土、盐土、灌耕土和灌淤土等。由于干旱少雨、植被稀疏,大多数土壤呈现出盐碱化、有机质少、沙性大和粘重板结等现象。但由于农业生产活动的发展,使得土壤逐渐熟化,有机质增加,土地承载力与环境容量提高,使得研究区成为了重要的农业生产基地。本研究主要以棉花为研究对象,其种植地土壤属性如表1所示。

表1 研究区内土壤属性

1.1.2植被类型

不同质地的土壤其理化性质差异性大,所以生长的植被也会有所不同。研究区内平原植被主要以人工种植作物棉花(Gossypiumhirssittum)、玉米(Corn)为主,天然植被主要是耐盐耐旱植物如红柳(Redwillow)、胡杨(Populuseuphratica)、猪毛菜(Salsola)、盐爪爪(Kalidiumfoliatum)等。

1.2 数据采集

1.2.1棉花叶片SPAD值

为估算模型构建和验证叶绿素精度,课题组在2019年7月进行野外调查,沿玛纳斯河流域根据地形地貌选取具有典型代表的104个样点并采样。利用日本美能达公司生产的SPAD快速叶绿素测定仪,在面积为10 m×10 m的样方内,从东、西、南、北方向的四个角和中间位置上对植株冠层的上、中、下采集20片健康、大小均匀的叶片测定SPAD值,并将五个点20个测定值的算术平均值作为对应样方内特定植被叶片的SPAD值。

1.2.2影像数据

欧洲委员会和欧空总署局联合倡议研发的全球环境与安全监测计划的第二颗卫星,Sentinel- 2A、B卫星分别于2015年6月23日、2017年3月7日发射成功,双星在轨运行的重返周期为5d,最佳分辨率可达10 m[23],其遥感影像拥有13个分辨率不同的光谱波段(10、20、60 m),包括可见光、红边、近红外、水汽、卷云以及短波红外波段(表2)。根据采样时间及云量(< 5%),择取2019年7月26日、28日两景影像(https：//scihub.copernicus. eu/dhus/#/home)。首先将下载的两景影像经Sen2Cor模块对影像实时辐射定标和大气校正,再利用SNAP软件使用三次内插卷积(Cubic convolution interpolation)的方法对大气校正后的波段进行空间分辨率为10 m的重采样,之后转入ENVI 5.1提取本研究所用9波段(B2、B3、B4、B5、B6、B7、B8a、B11、B12),然后将其在ArcGIS 10.0使用 (Gram-Schmidt Pan Sharpening)融合方法将研究区内两景影像拼接、融合,并按照研究区范围使用取子区(subset)功能将影像裁剪。

表2 Sentinel- 2A影像波段光谱参数

1.3 最适光谱指数

本研究选择Sentinel- 2A影像(B2、B3、B4、B5、B6、B7、B8a、B11、B12)9个光谱波段,其中包含3个红边波段,在农作物植株叶绿素含量估算中具有很大的应用潜力。因此,本研究在借鉴现有植被指数构建方法的基础上,引入红边波段替代可见光波段或者近红外波段,利用红边波段与近红外波段或者红波段的光谱反射特征表达不同叶绿素含量的差异。构建NDVI、DVI和RVI光谱指数,以采样点的棉花叶片SPAD值来反应Sentinel- 2A影像像元的数据,计算其与棉花叶片SPAD值相关性,以期找到评估叶绿素含量的最适光谱指数,如表3所示。

表3 植被指数

1.4 PROSAIL辐射传输模型

PROSAIL模型是将PROSPECT与SAIL模型的集成,已成功应用于不同的生态系统(如森林和作物)常用于植被生化物理参数的反演[26-29]。将PROSPECT模型的输出的叶片反射率和透射率,作为SAIL模型的输入参数,再结合土壤反射率、叶面积指数、叶倾角分布、太阳天顶角,观测天顶角等相关参数,即可得到植被冠层的反射率和透射率。遥感影像通过大气校正,也可以获取地表植被冠层反射率,将遥感影像与植被参数叶绿素通过物理过程联系起来,生成一个包含所有可能的反射率的数据集,完成了模型模拟过程[30]。

本研究所用波段的光谱反射率在490—870 nm之间,在该波段范围内LAI、叶绿素以及干物质含量均为敏感参数,将其确定为可变参数,类胡萝卜素在模型中表现为不敏感,设为固定值,反演参数主要通过野外实测、文献查阅、经验值及LOPEX94数据库获得[31],表4为参数设置。

表4 PROSAIL模型参数设置

1.5 棉花叶片SPAD值的模型构建

在构建棉花叶片SPAD值的反演模型,是通过R- 3.6.2中sample函数k-s划分70%(n=72)的建模集和30%(n=32)的验证集并用set.seed函数固定选中的数据集。

随机森林算法(Random forest,RF)作为一种较新的集成学习算法,是由Breiman[32]提出的基于决策树分类器融合算法,算法可以综合利用多种特征变量建模并评价,可以看成由多个弱预测器(决策树)集成的强预测器。在随机森林模型中[31],通过选择原始数据集和训练数据集(即校准数据集)的随机样本,使用确定性算法构建每棵树,每一棵树都依赖一个随机向量,通过对数据集的列变量和行变量观测进行随机化,生成多个分类树,最终将分类树结果进行综合。本文研究使用R- 3.6.2中随机森林工具包对棉花叶绿素SPAD值进行预测分类,并产生重要变量性估计,研究模型的重要性分类指标由平均下降精度参数(Mean Decrease Accuracy)提供[33]。

本文将野外测得的植被叶绿素SPAD值作为模型因变量,将最适光谱指数(DVI、NDVI、RVI),Sentinel- 2A中9个波段(B2、B3、B4、B5、B6、B7、B8a、B11、B12),基于物理模型生物协变量(叶面积LAI、冠层含水量LAI-Cw、冠层叶绿素LAI-Cab、植被覆盖度fcover、光合有效辐射fapar)共计17个特征参数分别作为自变量并构建3种方案(表5),在R- 3.6.2中利用随机森林按照设计方案对采样地的棉花建立植被叶片SPAD的估算模型。

表5 建模方案

1.6 模型精度评价

本研究为了评估上述3种预测模型的性能,选取了以下几种常用的模型验证标准。分别为决定系数(Determination coefficients,R2)、均方根误差(Root mean squared error, RMSE)、一致性相关系数(Concordance correlation coefficient,ρc),研究思路如图1所示。

图1 研究流程图Fig.1 Research flow process

图2 样本统计分析图Fig.2 Sample statistical analysis chart

2 结果与分析

由图2可知棉花叶片SPAD值全集均值为59.02,变异系数为0.092,变异系数小,说明数据离散程度小较为集中。如图2所示,棉花叶片SPAD值全集、建模集、验证集在统计特征中保持相似的统计形态,均值分别为：59.02、58.66、59.38,SD分别为5.43、5.25、5.43,由此可见样本划分合理。

2.1 棉花叶片SPAD值与高光谱植被指数的相关性

图3为棉花叶片SPAD值与Sentinel- 2A波段(B2、B3、B4、B5、B6、B7、B8a、B11、B12)两个随机光谱波段的2D指数(DVI、NDVI、RVI)关系。右侧的彩条表示r大小,深红、深蓝代表SPAD值与波段组合正、负相关性较高。结果表明,DVI、NDVI、RVI与棉花叶片SPAD值高度相关性的波段主要集中在红边(Red edge)波段(B5、B6)和近红边(Narrow NIR)波段B8a处。RVI与棉花叶片SPAD值相关性最高其为0.767,P**=0.195,其表达式B8a/B5。NDVI其r为0.720,DVI其r为0.692。说明红边波段对棉花叶绿素有较高的响应,3种指数在相应波段下对棉花叶片SPAD值均具有较高相关性。

图3 棉花叶片SPAD值与Sentinel- 2A(B2、B3、B4、B5、B6、B7、B8a、B11、B12)中任意两波段构建的DVI、NDVI、RVI的相关系数Fig.3 Correlation coefficient between SPAD value of cotton leaves and DVI、NDVI、RVI constructed in any two bands of sentinel- 2A (B2、B3、B4、B5、B6、B7、B8a、B11、B12) B2：蓝光 Blue band；B3：绿光 Green band；B4：红光 Red band；B5：植被红边Vegetation Red Edge band；B6：植被红边Vegetation Red Edge band；植被红边B7：Vegetation Red Edge band；B8a：近红外(窄)Narrow NIR band；B11：短波红外SWIR band；B12：短波红外SWIR band

2.2 特征变量重要性分析

本研究通过随机森林进行建模,并对参与不同方案建模的特征变量做出重要性分析,通过图4可以发现在三种方案中,红边波段的贡献率都大于6%,由红边波段参与构建的最适光谱指数在重要性占比中也较高,这也说明了红边波段与叶绿素响应度较强；此外三方案中RVI与LAI-Cab的贡献始终最大其占比均超过了10%,在方案3中重要性分别10.22%和10.56%,由此可见植被指数与物理模型对于估算叶绿素含量的贡献两者均较为重要,皆不可忽视。

图4 特征变量重要性Fig.4 The importance of characteristic variablesM1:方案1; M2:方案2; M3:方案3;LAI：叶面积指数 Leaf area index；LAI-Cw：冠层含水量 Canopy water content；LAI-Cab冠层叶绿素 Canopy chlorophyll；fcover：植被覆盖度 Vegetation coverage；fapar：光合有效辐射 Photosynthetically active radiation；

2.3 棉花叶片SPAD值的估算模型评价结果

根据表6可知在所有棉花叶片SPAD估算模型中,由验证集效果可知,方案3 >方案1 >方案2,其中方案3 (R2= 0.8258,RMSE=3.5660,ρc= 0.7205)的效果甚佳。比单一的使用最适光谱指数或物理模型生物协变量构建的模型精度提高0.05—0.07左右,这说明将Sentinel- 2A波段与植被指数、物理模型集成进行优势互补可以提高棉花叶片SPAD值估算精度。同时一致性ρc也为三个方案中最高的。

表6 棉花叶片SPAD预测值与实测值分布

图5 棉花叶片SPAD预测值与实测值相关性及精度验证Fig.5 Correlation and accuracy verification of SPAD predicted value and measured value of cotton leaves

图6 叶片SPAD遥感监测图Fig.6 Remote sensing monitoring map of blade SPAD

3 讨论

作物叶片的SPAD值反映了作物叶绿素含量的相对高低,通过了解叶绿素含量及其变化,对作物的营养诊断、长势监测及产量估算具有重要意义[34]。本研究以Sentinel- 2A为数据源,将PROSALL模型与最适光谱指数进行结合对棉花叶片SPAD值进行估算研究。

本文将3种方案进行对比后,选择方案3绘制图6叶片SPAD遥感监测图,由图6可知研究区内叶绿素的空间分布大致为中部高四周低,其主要受水源、土壤盐分、土壤有机质及灌溉方式影响,研究区的西北角及东北角属于沙漠地区,土壤水分少、盐分高、有机质含量低对棉花生长有一定抑制作用。研究区内的某些零星区域叶绿素偏高,是由于种植其他作物对其估算具有一定不准确性。同时通过表1可以发现土壤属性存在较大差异会对所生长植被的生理参数产生较大影响。因此,今后研究要对棉花叶绿素进行空间进行分析,除了要掌握主要水源分布,还应将土壤性质(盐分、水分、有机质等)作为重要的参考因素作为土壤特征变量参与建模,同时也要了解作物的种植分布特点。

Sentinel- 2A具有丰富的光谱波段信息,选取特征波段构建光谱指数在不同程度上抵消了由环境变化而引起的噪声,可以用来定量估算作物的理化参数[35]。因此,本研究利用Matlab,选取最优光谱波段构建3个2D最适光谱指数(DVI、NDVI、RVI),由于选取的最优波段均为红边波段与叶绿素具有较强的响应,且采样时棉花正处于花铃期,此时棉花生长旺盛,叶片呈深绿色,对光谱的反射率较高。这与洪帅[36]等研究的结果一致,棉花的生育时期冠层叶片的叶绿素含量与红边指数的相关性较高。但仅构建3种光谱指数无法充分利用高光谱数据,从而在一定程度上限制估算模型的精度[10],充分利用高光谱数据构建多种光谱指数之后进行最优筛选以提高模型精度值得今后研究。

在对3种模型方案不同的特征变量进行重要性分析,易知红边波段与基于红边波段构建的光谱指数对估算模型有着关键作用。生物协变量中LAI-Cab贡献最大,其在物理模型中就是表示叶绿素,同时也说明了物理模型对叶绿素估算起到一定作用。在以最适光谱指数和物理模型生物协变量构建的方案1、方案2中,RVI与LAI-Cab重要性占比为13.90%、15.24%,但在方案3中两者占比虽依旧为最大但却有所下降,可能是由于在特征变量增多的原因,未对变量进行筛选,往后的研究上应在优选特征变量的基础上进行建模[37]。另外,PROSAIL模型反演过程中受到自身参数不确定性及敏感性以及植被生理参数之间的相互作用可能带来误差影响,本研究由于采样条件有限,故敏感参数设置仅为LAI、Cab,没有考虑其他参数[38],从而会影响模型精度,今后研究可多获取参数实测值来解决。

本研究构建的三种模型方案,都是基于随机森林回归算法估算棉花叶片SPAD值的,出现了低值高估、高值低估的现象,同时给出不同变量对模型影响评分,在建模时进行参数优化,可以很好地学习特征变量与SPAD值关系,提高模型的预测精度。这与前人所做出的结果一致[39]。但构建的模型要想快速的监测棉花叶绿素含量,还需要获取多年数据进行改进。其次,本文研究区以新疆人工绿洲为主,它是由绿洲农田和绿洲聚落生态系统构成的人工生态系统主要受人类控制,在人类精心维护下,绿洲人工生态系统内部生态环境趋于优化,生态稳定性较强,而生态环境的影响通常以作物物候周期为尺度表现在作物生理化学参数上。基于此,本文尝试通过对棉花物候周期中的某一时期的SPAD高精度制图可能,在后续研究中将棉花的而整个物候周期作为研究的重点,从而了解生态环境对于棉花生长的影响,通过整个物候周期生理参数的变化来分析生态环境的质量发展变化,以期实现对干旱区作物生长状况的实时监测,为精准农业提供新的技术途径。

4 结论

本文以Sentinel- 2A为高光谱数据源,将Sentinel- 2A光谱波段、植被指数和PROSAIL模型三者进行组合构建棉花叶片SPAD值的估算模型。利用Sentinel- 2A光谱波段在构建植被指数时可以发现建立的2D指数(DVI、NDVI、RVI)与棉花叶片SPAD值相关系数最好的波段在红边波段(B5、B6)和近红边波段(B8a)处,说明红边波段与叶绿素含量响应明显,基于红边波段构建的植被指数对估算模型起到关键作用。在构建的3个方案模型中,利用实测数据对模型进行验证,结果表明：将Sentinel- 2A光谱波段、植被指数和PROSAIL模型三者利用随机森林算法构建的模型表现出了最佳SPAD预测能力,实现了棉花叶片SPAD值的遥感监测。综上所述,单一的以经验模型植被指数或物理模型构建的反演方法对植被叶绿素含量具有一定的预测能力,但将两者在随机森林中集成的模型可以提高干旱区和半干旱区棉花叶片SPAD值的预测精度,可以推广至其他植物的叶绿素监测,为合理水肥配置,精准种植,农作物生长过程监测提供数据支撑,促进精准农业蓬勃发展。