杭艳红,苏 欢,于滋洋,刘焕军,※,官海翔,孔繁昌
(1. 东北农业大学公共管理与法学院,哈尔滨 150030;2. 中国科学院东北地理与农业生态研究所,长春 130012)
水稻作为重要的粮食作物,其生产关系到中国的粮食安全,而叶面积指数作为水稻生长的一项重要参数,能够提供水稻生长过程的动态信息,与地上部生物量的积累和产量形成密切相关[1],是评价作物长势及指导田间管理的重要指标。目前叶面积指数(Leaf Area Index,LAI)地面测量方法包括直接测量法和间接测量法。直接测量法精度相对精确,但费时费力,且离散的地面测量点不利于大规模的监测农作物长势[1],测量过程中还会对作物造成一定破坏;间接测量法是采用测量参数或光学仪器进行测量,常用的有LAI-2000植被冠层分析仪、Sunscan密度计等[2]。卫星遥感具有大面积同步观测的特点,被广泛应用于LAI估算中,但受分辨率和重访周期的限制,相比卫星遥感,无人机遥感不仅具有更高的时间和空间分辨率,而且能够快速地获得田块尺度作物的空间变异信息,能够准确反映地块内部的长势差异。
近年来,国内外学者针对植被指数估算LAI展开了大量研究,在机理上已经被证明植被指数能够较好地估算作物LAI[3]。根据绿度归一化植被指数(Green Normalized Difference Vegetation Index,GNDVI)和LAI的线性关系建立模型[4],但植被指数在LAI较大时会出现饱和现象[5],在估算LAI时存在一定局限性。近些年基于无人机影像,结合光谱特征和纹理特征估算作物生长状况的研究逐渐增多,研究发现将植被指数与纹理特征融合,可以有效提高生物量[6]、叶绿素含量[7]和氮营养指数[8]的估算精度,国内外结合光谱信息和纹理特征估算生物量[9-10]的研究有很多,但在估算作物LAI方面研究较少,且以往的研究多是直接将众多纹理特征值参与到模型中,缺少对纹理特征值的优化。同时,有研究发现在植被指数的基础上引入作物覆盖度也可用于估算LAI[11],但以往的研究多是在传统基于植被指数的LAI估算方法的基础上引进纹理特征或作物覆盖度等指标,缺少对光谱特征、纹理指数和作物覆盖度等指标共同估算LAI能力的研究。
综合以上研究存在的问题,该研究旨在利用具有高空间分辨率的无人机多光谱影像,通过对水稻冠层纹理特征进行组合运算来建立新的纹理指数,以期提高纹理特征值与水稻LAI的相关性,结合光谱特征、纹理指数和作物覆盖度,构建多指标结合的水稻LAI估算模型,以期提高水稻LAI估算精度。
水稻田间试验于2018年7-8月在黑龙江省五常市试验站(44°44′N,127°12′E)进行(图1)。研究区属温带大陆性季风气候,平均海拔约为450 m,年平均降雨量在500~800 mm,年平均气温3~4 ℃。采用随机区组设计进行了稻田试验。施肥处理采用移栽后施一次基肥,分蘖期和抽穗期之前进行两次追肥的施肥方式。基肥、分蘖肥和穗肥的施肥比例为5:2:3。施肥品种为磷酸氢二铵、尿素和硫酸钾,磷酸氢二铵的含P2O5量为46%,含N量为17%,尿素的含N量为46%,硫酸钾的含K2O量为50%(以上含量均为质量分数)。试验包括5种氮肥处理,每种施氮处理有3次重复试验,共30个样区,每个样区面积为50 m2,包括五优稻4号和松粳6号两个水稻品种。
1.2.1 地面数据获取
该研究于2018年水稻拔节期(7月22日)和抽穗期(8月8日)2个关键生育期分别采集了水稻LAI 数据,共有60个样本。LAI-2000测量的结果为总叶面积指数(Plant Area Index,PAI),实际上是包含LAI和木质面积指数(Woody Area Index,WAI),对于树林冠层,WAI的影响较大,不可用PAI直接替代LAI[12];而水稻属于水平均匀植被,WAI的影响较小,LAI-2000在连续、均质的冠层测量时具有较好的结果[13]。因此采用美国LI-COR公司生产的LAI-2000植物冠层分析仪测量LAI,测量时间为上午08:00-10:00。该研究在每个样区设置一个植被覆盖较均一且具有区域代表性的地点作为采样区域,每个生长阶段在不同的采样区进行LAI测定,每个采样区选择3处测量并记录LAI值,取3次测量平均值作为该样区的LAI值。每个采样区先测量冠层上方的1个天空光,然后测量冠层下方的4个测试目标值,样方的LAI值由 LAI-2000自带的软件计算。为了避免太阳光线直射带来的测量误差值,数据采集选择上午10:00之前直射光较弱的时段,测量时操作员背对阳光,并且在天气晴朗时使用遮盖帽进行遮光处理。
1.2.2 多光谱数据获取与处理
该研究采用大疆M600 Pro六旋翼高性能无人机搭载MicaSense RedEdgeTM3多光谱相机。大疆M600 Pro六旋翼高性能无人机最大飞行承载质量为6 kg,可承受最大风速8 m/s,可在-10~40℃的环境中工作,续航时间为25~35 min。MicaSense RedEdgeTM3多光谱相机拥有蓝、绿、红、红边、近红外5个波段(表1),配有光强传感器、全球定位系统(Global Positioning System, GPS)模块和白板[14]。MicaSense RedEdgeTM3是美国MicaSense公司生成的一款专业的多光谱相机,能够同时捕获5个不连续的光谱带,光谱采集范围为400~900 nm,拥有5个通道。MicaSense RedEdgeTM3在执行任务时,能够达到每秒一次捕捉。飞行过程中,光强传感器可用于校正太阳光线变化对影像造成的影响;GPS模块同时记录每张影像的位置信息;白板具有固定的反射率信息,可用来对影像反射率进行校正。
表1 多光谱相机波段 Table 1 Multispectral camera bands
数据采集时间为2018年7月22日和2018年8月8日,天气晴朗无云,风速小于3级,适合无人机飞行,选择的无人机遥感平台对地观测时间为8:00-10:00。为确保影像的完整性和准确性,地面站设置航向重叠为80%,旁向重叠为75%,主航线角度90°,飞行高度设定为110 m,图像空间分辨率为7 cm。
为去除姿态角度异常和由大气传输等因素导致的成像问题,对采集到的无人机影像先进行筛选;对筛选后的影像进行数据预处理,将影像导入Pix4D mapper软件,调节处理参数,采用自动生成的连接点与无人机内置的定位定向系统(Position and Orientation System, POS)数据进行空三运算,生成三维点云数据;最后对影像进行仿射变换,得到该研究区的无人机正射多光谱影像[15]。
1.3.1 植被指数计算
植被指数(Vegetation Index,VI)能够简单、有效地度量地表植被状况,植被指数是反映植被在可见光、近红外波段反射与土壤背景之间差异的指标,各植被指数能够在一定条件下用来定量表明植被的生长状况。归一化差值植被指数(Normalized Difference Vegetation Index,NDVI)对绿色植被表现敏感,常被用于LAI估测,但在植被覆盖度较高时呈现饱和状态,对植被检测灵敏度下降;GNDVI将绿波段代替NDVI中的红波段,对叶绿素a的敏感度更高;差值植被指数(Difference Vegetation Index,DVI)、修正三角植被指数(Modified Triangular Vegetation Index 2,MTVI2)和优化土壤调节植被指数(Optimized Soil-adjusted Vegetation Index,OSAVI)能够降低土壤和环境背景等影响因子的影响;红边叶绿素指数(Red-edge Chlorophyll Index,CIRE)与植被叶片的叶绿素浓度具有高度相关性,能够反应植被生长状况。该研究选取6个与LAI相关的植被指数对水稻LAI进行估算(表2)。
表2 植被指数及公式 Table 2 Vegetation indexes and formulas
1.3.2 纹理特征提取
在几种纹理算法中,常用灰度共生矩阵(Gray-Level Co-occurrence Matrix,GLCM)测试纹理分析。在进行辐射校正后,使用ENVI 5.3软件计算8类基于GLCM的纹理特征值,包括均值(Mean,mean)、方差(Variance,var)、协同性(Homogeneity,hom)、对比度(Contrast,con)、相异性(Dissimilarity,dis)、信息熵(Entropy,ent)、二阶矩(Second Moment,sm)和相关性(Correlation,corr)[22],共40个纹理特征值。无人机影像的分辨率为7 cm,大多数窗口大小涉及土壤背景和水稻植株,因此在该研究中,纹理分析选取最小的3×3窗口进行分析。该研究区的水稻田垄向接近90°,因此,在纹理分析时方向选取90°。
该研究所用3种纹理指数(Texture Index,TI)定义分别为:归一化差值纹理指数(Normalized Difference Texture Index,NDTI)[23],差值纹理指数(Difference Texture Index,DTI),比值纹理指数(Ratio Texture Index,RTI)。构造了5个波段(475、560、668、717和840 nm)与8类基于GLCM的纹理特征值相结合的所有可能的两种纹理测量组合,以探索它们的估算能力。利用Matlab软件实现对纹理指数NDTI,DTI,RTI的运算,并计算出不同纹理指数与LAI之间的相关性。
NDTI,DTI,RTI定义分别如下:
式中T1和T2为随机任意波段的纹理特征值。
1.3.3 作物覆盖度提取
作物覆盖度(Crop Coverage,CC)是指植被在地面的垂直投影面积占总面积的百分比。由于作物覆盖度与LAI之间具有显著的相关性,本文将提取作物覆盖度并将其输入水稻LAI估算模型中。本文采用像元二分法的思路,使用ArcGIS输出无人机影像的灰度图,再对灰度图应用OTSU算法,寻找令植被与非植被类间方差最大的灰度阈值,将非植被区域灰度值置零,提取出植被区域的像元,进而计算出该影像的作物覆盖度[24]。
该研究中,选取拔节孕穗期和抽穗期试验的数据进行模型构建,总样本量60,采用K-折交叉验证法作为验证方法。选取3种方法建立模型,分别为一元线性模型、多元逐步回归模型、人工神经网络模型。采用决定系数(Determination coefficient,R2)、调整后决定系数(Adjust determination coefficient,R2adj)和均方根误差(Root Mean Square Error,RMSE)评价模型的精度。一般来说,R2和R2adj越高,表明模型越稳定;RMSE越低,表明模型精度越高。
1.4.1 K-折交叉验证
本文采取K-折交叉验证的方式建模并评估模型精度。K-折交叉验证是应用最广泛的泛化误差估计方法之一,对于检验模型具有普适性,适用于样本集不大的情况[7],将数据集均分为K-份,每次训练取其中一份作为验证集,剩余部分作为测试集,重复K次,完成一次交叉验证[25]。K-折交叉验证保证了数据集中的每个数据都参与了建模与验证,能够让模型更好地学习到数据中的特征;同时K-折交叉验证建模的精度是K次训练精度的均值,能更好地估计模型的泛化误差。该研究采取5-折交叉验证,即训练集样本个数为48,测试集样本数为12,训练5次。
1.4.2 线性模型
植被指数法是光学遥感估算植被LAI常用的经验方法。分别选用以上6种植被指数、纹理指数和作物覆盖度与LAI进行相关性分析,并建立一元线性模型估算LAI。选择植被指数、纹理指数和作物覆盖度等因子作为变量,构建多元线性逐步回归模型估算LAI。
1.4.3 人工神经网络模型
基于误差反向传播(Back Propagation,BP)的人工神经网络是应用最广泛的神经网络,它的基本思想是利用梯度下降法,逐层求出目标函数对各神经元权值的偏导数,作为权值更新的依据,使模型学习达到期望的性能[26]。本文使用的神经网络主要包含三部分,输入层、隐藏层和输出层,其中输入层不参与运算,隐藏层和输出层为全连接层。本文以均方误差为目标函数,学习率为0.01建立了一个包含5层隐藏层的人工神经网络模型。选择植被指数、纹理指数和作物覆盖度等因子作为变量,输入到人工神经网络模型估算LAI。
植被指数是根据地物光谱反射率的差异进行运算,以此突出图像中植被的特征。对选取的6个植被指数和LAI的相关性进行分析。结果显示,植被指数与LAI之间的相关性均达到极显著水平(P<0.01),其中与LAI相关性由大到小依次是:OSAVI、MTVI2、DVI、GNDVI、NDVI、CIRE,相关系数依次是0.763、0.761、0.731、0.718、0.667、0.606。对1.3.3中提取的作物覆盖度与实测水稻LAI进行相关性分析,结果显示显著相关且相关系数为0.785。
单一纹理特征与LAI的相关性分析结果(表3)显示多数纹理特征与LAI的相关性并不高,只有少部分纹理特征显示出与LAI存在较高的相关性,其中LAI与近红外波段的纹理特征均值的相关性最高,相关系数为0.731。
表3 单一纹理特征和LAI的相关性分析 Table 3 Correlation analysis between single texture feature and leaf area index
由于纹理特征值与LAI的相关性较低,因此构建了由不同纹理特征值组成的指数NDTI,DTI,RTI以提高纹理特征值对LAI的估测能力。结果如图2显示,通过对纹理特征进行组合运算,在整体上明显提高了纹理特征值与LAI的相关性,3种指数的均值组合与LAI具有较高的相关性。优选出每种组合运算纹理指数中与LAI 相关性最高的纹理指数,按相关性由高到低依次是近红外波段和蓝波段均值的差值DTI(mean5,mean1)、近红外波段均值与绿波段信息熵的比值RTI(mean5,ent2)、近红外波段均值与绿波段信息熵的归一化值NDTI(mean5,ent2),相关系数依次是0.830、0.798、0.791。结果证明对纹理特征进行组合运算,能够明显提高纹理特征值与LAI的相关性,近红外波段均值与蓝波段均值的差值较近红外波段均值提高了13.54%。
根据多光谱影像提取的输入量,建立单一输入量与试验区实测LAI的一元线性模型,由结果(表4)可知,6种遥感植被指数与LAI建立的一元线性模型中,GNDVI模型R2最高,达0.603,R2adj为0.563,RMSE为0.541;优选纹理指数与LAI建立的一元线性模型中,DTI(mean5,mean1)模型R2最高,达0.668,R2adj为0.635,RMSE为0.447;作物覆盖度与LAI建立的一元线性模型R2为0.633,R2adj为0.596,RMSE为0.516。对比不同单一输入量与试验区实测LAI的一元线性模型结果发现,精度由高到低依次为DTI(mean5,mean1)、作物覆盖度、GNDVI。
本研究先将3类输入量分别两两类组合放入多元逐步回归模型中构建LAI估算模型,再将3类输入量全部放入多元逐步回归模型中构建LAI估算模型,对比不同输入量组合的模型精度并进行评价。由结果(表4和表5)可知,结合植被指数和纹理指数的多元逐步回归模型(R2=0.728,R2adj=0.668,RMSE=0.421)明显优于单一植被指数模型(R2=0.603,R2adj=0.563,RMSE=0.541),结合植被指数、纹理指数和作物覆盖度共同估算LAI模型(R2=0.866,R2adj=0.816,RMSE=0.308)精度最高,优于所有单一输入量和任意两类输入量组合。因此,结合植被指数、纹理指数和作物覆盖度共同估算LAI可以被认为是一种能够有效改善LAI估算精度的方法。
表4 LAI的一元线性估算模型与精度评价 Table 4 Unary linear estimation model and precision evaluation of LAI
为了对多元逐步回归模型筛选出的变量进行进一步分析,揭示各输入因子对模型的影响,本文引入“贡献度”这一指标。此处的贡献度为偏相关系数的绝对值与模型中个因子的偏相关系数绝对值总和之比[7]。如表5所示,基于多指标结合的LAI估算模型中,作物覆盖度和纹理指数DTI(mean5,mean1)对模型的贡献度均很大,植被指数GNDVI贡献度相对较小。作物覆盖度、DTI(mean5,mean1)、GNDVI对模型的贡献度分别为40.413 %、35.507%、24.08%。
表5 LAI的多元线性逐步回归估算模型与精度评价 Table 5 Multiple linear stepwise regression estimation model and accuracy evaluation of LAI
为了排除模型对文章结论的影响,证明本文结论的普适性,本文使用人工神经网络对结论进行进一步验证。本文将3类输入量依次放入人工神经网络模型中构建LAI估算模型,对比模型精度。结果如表6所示,将植被指数与纹理指数结合构建LAI模型,发现模型精度(R2=0.76,R2adj=0.707,RMSE=0.365)较单一植被指数构建的LAI模型(R2=0.654,R2adj=0.619,RMSE=0.431)精度得到明显提高;将植被指数、纹理指数和作物覆盖度共同构建LAI模型,发现模型(R2=0.836,R2adj=0.775,RMSE=0.286)精度最高。基于人工神经网络模型估算LAI得到的结论与基于多元逐步回归模型估算LAI得到的结论一致,再次验证结合植被指数、纹理指数和作物覆盖度共同估算LAI,是一种能够有效改善LAI估算精度的方法。
表6 LAI的人工神经网络预测模型与精度评价 Table 6 Artificial neural network forecast model and accuracy evaluation of LAI
在以往的研究中,多是利用光谱信息与地面实测数据来对LAI进行估算,但仅仅利用植被指数进行LAI估算存在一定的局限性。因此,从无人机影像中提取纹理特征和作物覆盖度来提高LAI估算精度具有重要意义。
纹理是物体表面的内在特性,不依赖于颜色和亮度而变化,能够抑制同谱异物和同物异谱现象的发生[27]。大多数纹理特征值与LAI的相关性较弱,该研究对纹理特征值进行组合运算得到新的纹理指数,有效提高了纹理特征估算LAI的性能。纹理归一化能够降低土壤背景、太阳角和传感器视角影响[28];纹理比值可以放大地物波谱特征间的细微差别,降低影像中地形、阴影等带来的影响,突出地物特征;纹理差值能够去除影像中的相同背景[27],遥感图像中由阴影产生的影响在纹理差值图像中表现得不明显,因此纹理信息与光谱特征的结合在一定程度上削弱了阴影对变化检测产生的影响。
该研究利用多元逐步回归模型,将植被指数、纹理指数和作物覆盖度结合起来建立水稻LAI估算模型。结合植被指数和纹理指数的多元逐步回归模型(R2=0.728,R2adj=0.668,RMSE=0.421)明显优于单一植被指数模型(R2=0.603,R2adj=0.563,RMSE=0.541),结合植被指数、纹理指数和作物覆盖度共同估算LAI(R2=0.866,R2adj=0.816,RMSE=0.308)精度最高,优于所有单一输入量和任意两类输入量组合。这主要是因为多指标结合模型综合了植光谱特征、纹理特征和作物覆盖度对LAI估算的共同贡献。其中GNDVI对叶绿素a较为敏感,绿叶的叶绿素在光照条件下发生光合作用产生植物干物质积累,使叶面积增大[29];水稻的覆盖度随着LAI的增加而增加,覆盖度可以直接反映从俯视图中提取的LAI值[30];纹理指数DTI(mean5,mean1)提供了作物的空间特征,一定程度上弥补了光谱信息的不足,其中纹理测量均值包含移动窗口中的目标和背景的平均值,可以平滑图像,使背景的干扰最小化,又因为绿色植被在蓝波段对光强吸收,在近红外波段由于冠层结构的漫反射导致近红外处反射率较大[31],近红外波段与蓝波段的差值能增强植被对光吸收和反射的差异,因此能够更好反映绿色植被冠层结构。光谱特征、纹理特征和作物覆盖度的结合为LAI估算模型提供了更多的信息。该研究尚存在部分局限性,考虑到影响水稻LAI的因素不止光谱特征、纹理特征和作物覆盖度,还需要综合考虑原始波段、叶绿素、土壤背景以及影像质量等因素的影响。
本文基于无人机多光谱影像,对光谱特征、纹理特征和作物覆盖度估算水稻叶面积指数(Leaf Area Index,LAI)的能力进行了分析,主要得到以下结论:
1)对纹理特征进行组合运算,能够明显提高纹理特征与LAI的相关性,近红外波段均值与蓝波段均值的差值较近红外波段均值提高了13.54%,且归一化差值纹理指数、差值纹理指数、比值纹理指数的均值组合与LAI具有较高的相关性,其中近红外波段和蓝波段均值的差值具有最高的相关性,相关系数为0.830。
2)基于一元线性模型构建的水稻LAI估算模型的精度由高到低为近红外波段和蓝波段均值的差值(决定系数R2=0.668,调整后决定系数R2adj=0.635,均方根误差RMSE=0.447)、作物覆盖度(R2=0.633,R2adj=0.596,RMSE=0.516)、绿度归一化植被指数(R2=0.603,R2adj=0.563,RMSE=0.541)。
3)对比单一植被指数、纹理指数、作物覆盖度,以及任意两类输入量的结合,发现多特征结合估算水稻LAI的能力最强(R2=0.866,R2adj=0.816,RMSE=0.308),人工神经网络模型结果对本文结论进行进一步验证,表明结合光谱特征、纹理指数、作物覆盖度共同估算LAI具有较好的精度,为监测水稻生长状况提供一种可行的方法。