不同分辨率无人机数码影像的马铃薯地上生物量估算研究

2021-05-11冯海宽杨福芹杨贵军

光谱学与光谱分析 2021年5期

刘杨，冯海宽，孙乾，杨福芹，杨贵军

1.农业部农业遥感机理与定量遥感重点实验室，北京农业信息技术研究中心，北京 100097 2.山东科技大学测绘科学与工程学院，山东青岛 266590 3.国家农业信息化工程技术研究中心，北京 100097 4.北京市农业物联网工程技术研究中心，北京 100097 5.河南工程学院土木工程学院，河南郑州 451191

引言

地上生物量(above ground biomass，AGB)是作物生长发育过程中进行光合作用的产物，能够反映作物长势状况和营养状态，与作物产量有着密切的联系[1]。因此，快速精准地监测AGB信息对指导农田精准管理和挖掘生产潜力具有重要的意义[2]。传统的AGB测定主要以田间取样调查为主，虽能够达到一定的估算精度，但破坏性强、非动态性和花费成本高，且局限于小面积[3-4]。近年来，将遥感技术用于反演作物的生化参数，检测生态环境信息等方面得到广泛关注[5]。相较于其他遥感技术，无人机遥感具有运营成本较低、操作更为简单、可云下获取高分辨率影像等优势，更适合田块尺度观测，是现代化农业向定量精细化迈入的重要技术手段[6-7]。

目前，基于无人机遥感平台搭载的多种传感器(如高光谱和激光雷达)相继用于作物AGB获取，但因传感器价格昂贵，制约了无人机技术在作物AGB监测中的应用。数码相机传感器因价格低廉、空间分辨率较高且获取数据后续处理较为简单等特点，作为一种低成本获取作物AGB的方式，已经得到国内外许多学者的关注。崔日鲜等[2]基于无人机数码相机获取冬小麦影像，通过构建植被指数利用BP神经网络和多元线性回归估算冬小麦AGB，显示神经网络AGB估算模型优于多元线性回归模型。陶惠林等[8]基于无人机数码影像，采用逐步回归、偏最小二乘和随机森林3种方法以植被指数、植被指数加实测株高和植被指数加从作物表面模型(crop surface model，CSM)中获得的估测株高3种变量估算冬小麦AGB，结果表明每种方法中加入估测株高能够明显提高模型的精度。Zhang等[9]基于无人机数码影像，采用偏最小二乘、BP神经网络、支持向量机和随机森林4种方法以5种光学植被指数构建冬小麦AGB估算模型，结果表明偏最小二乘法相比其他方法构建的模型效果更优。以上成果表明，无人机数码影像能够很好的对作物AGB进行有效监测，但存在一定的不足，仅利用数码影像的光谱信息计算植被指数进行参数估算，模型的精度容易出现饱和现象。为了解决这一问题，众多学者通过引入作物株高和纹理特征分别结合光谱特征来估算作物参数，达到了理想的效果。相较于作物株高而言，纹理信息更具多面性，能够反映图像的灰度属性和空间位置关系，将其与光谱信息结合能够更好地减弱仅利用植被指数反演作物参数出现的低估现象，提高估算模型的适用性[10]。当前，利用数码影像以二者融合为变量估算作物AGB的研究较少，主要以高光谱和多光谱为主。刘畅等[11]利用无人机高光谱影像，比较分析光谱特征、纹理特征和二者结合对冬小麦AGB估算精度的影响，结果表明，光谱特征和纹理特征相结合能够解决AGB估算出现的过早饱和问题，而且模型估算效果优于单一光谱特征或单一纹理特征。Zheng等[12]基于无人机多光谱影像，以光谱信息和光谱信息结合纹理特征2种类型变量对水稻AGB进行估算，结果表明以后者为输入变量能够更好的估测AGB。Liu等[13]基于无人机多光谱影像，采用偏最小二乘和随机森林2种方法以光谱特征、光谱特征+纹理特征估算油菜AGB，经比较可知2种方法以融合指标构建的模型效果较优。纹理特征的提取十分依靠于获得的影像分辨率，以上成果并未考虑到不同飞行高度下所获得的影像特征对小麦、油菜和水稻AGB估算精度的影响。因此，探讨不同分辨率的无人机数码影像对作物AGB监测结果的影响，具有科学指导意义。

本工作利用马铃薯块茎增长期不同飞行高度的无人机数码影像，结合地上部实测生物量数据，分析不同分辨率影像下获取的光谱特征和纹理特征对马铃薯植株地上部生物量估算结果的影响，为无人机遥感监测AGB时，对飞行高度和影像分辨率的选择提供技术参考。

1 实验部分

1.1 研究区域与实验设计

在北京市昌平区小汤山镇国家精准农业研究示范基地的马铃薯试验田开展试验，该处位于北纬40°10′34″，东经116°26′39″，平均海拔为36 m，气候类型为暖温带半湿润大陆性季风气候。试验为小区随机试验设计，采用了中薯5(P1)和中薯3(P2)这2个均为早熟的马铃薯品种，试验田共计48个小区，每个小区面积为6.5 m×5 m。试验区共设密度试验(N区)、氮素试验(S区)、钾肥试验(K区)3个试验区，每种试验3种重复，每种重复进行不同程度的密度、氮素和钾肥处理，其中密度设置了3种水平(T1：60 000株·hm-2、T2：72 000株·hm-2、T3：84 000株·hm-2)，氮素设置了4种水平(N0：0 kg尿素·hm-2、N1：244.65 kg尿素·hm-2、N2：489.15 kg尿素·hm-2(正常处理，15 kg纯氮)、N3：733.5 kg尿素·hm-2)，钾肥设置了3种水平(K0：0 kg钾肥·hm-2、K1：970.5 kg钾肥·hm-2(N区和S区均为K1处理)、K2：1 941 kg钾肥·hm-2)。为了更好的获取试验田的位置，在试验小区周围均匀布控11个GCP(k1—k11由0.3 m×0.3 m的木板和埋于地下的木桩组成,其上有黑白标志的塑料板，目的是准确确定木板的中心位置)并用差分GPS测定其三维空间位置，详细的试验设计见图1所示。

图1 马铃薯研究区域和试验设计Fig.1 The experimental area and design for potato study

1.2 地面数据获取

马铃薯AGB通过收获法获取，在每个小区中选取能够代表小区长势水平的3棵植株，进行取样并迅速带回实验室将其茎叶分离，随后用流水洗净，105 ℃杀青，80 ℃烘干48 h以上，直到质量恒定再进行称量。将植株茎和叶的干质量求和得到样本干质量，最终通过群体密度和样本干质量，得到每个小区的马铃薯地上生物量。计算公式如式(1)

AGB=[(茎干重+叶干重)/取样株数]×种植密度

(1)

1.3 无人机数码影像获取

采用大疆公司生产的精灵系列4A产品在马铃薯块茎增长期进行无人机遥感数据获取，无人机参数见表1所示。无人机飞行时间为中午12:00，此时太阳光照强度稳定。人为设置了10，20，30，40和50 m共5种无人机飞行高度，相应的地面分辨率约为0.43，0.85，1.30，1.83和2.05 cm。

表1 无人机部分参数Table 1 Partial parameters of UAV

1.4 影像预处理

通过利用俄罗斯Agisoft LLC公司研发的基于运动三维重建(structure from motion，SFM)技术的PhotoScan 扫描软件对获取的不同高度的数码影像进行拼接处理，得到马铃薯试验区域的数字正射影像(digital orthophoto map，DOM)。首先将带有经纬度、高度、偏航角、俯仰角和旋转角共6种POS信息的每张影像导入软件中，还原影像拍摄时的姿态和空间位置；将GCP添加到软件中，利用GCP的三维坐标信息优化照片对齐，生成飞行区域的稀疏点云；利用稀疏点云构建密集点云，进而生成飞行区域的空间格网并建立纹理信息；生成马铃薯试验田的DOM。

1.5 影像指数选取

利用数码影像指数可以较好的反映作物的AGB信息，因此，通过DOM，利用ENVI软件提取每个试验小区影像的平均红波段值R、绿波段值G和蓝波段值B。并将R，G和B三个波段的DN值进行归一化处理，分别得到r，g和b。计算公式如下

r=R/(R+G+B)

(2)

g=G/(R+G+B)

(3)

b=B/(R+G+B)

(4)

根据文献中选取的能够很好反映AGB信息的15个数码影像指数，加上定义的3个，总共18个影像指数用来估算马铃薯地上部生物量，具体如表2所示。

表2 数码影像指数Table 2 Digital image indexes

1.6 纹理特征选取

纹理特征通过像素与空间邻域的灰度分布变化来表现，能够反映目标地物与周围环境的空间位置关系，是遥感影像中一种重要的图像信息。通过ENVI软件，利用灰度共生矩阵提取红、绿、蓝3个波段的0°，45°，90°和135°4个方向的纹理特征，不同方向的纹理特征分别求平均值，得到每个波段的8个纹理特征，分别为对比度(contrast，con)、相关性(correlation，cor)、差异性(dissimilarity，dis)、熵(entropy，ent)、同质性(homogenetity，hom)、二阶矩(second moment，sm)、方差(variance，var)和均值(mean)。红波段的纹理特征为con_R，cor_R，dis_R，ent_R，hom_R，sm_R，var_R和mean_R。绿波段的纹理特征为con_G，cor_G，dis_G，ent_G，hom_G，sm_G，var_G和mean_G。蓝波段的纹理特征为con_B，cor_B，dis_B，ent_B，hom_B，sm_B，var_B和mean_B。

1.7 数据分析方法

1.7.1 多重共线性分析

采用方差膨胀因子(variance inflation factor,VIF)评价输入变量之间存在的自相关性，当方差膨胀因子值较大时，输入变量之间存在高度自相关，会造成数据冗余和计算困难，因此有必要对数据进行降维。一般当VIF<10可视为输入变量之间不存在自相关性；当10≤VIF≤20时，认为输入变量之间存在一定的自相关性；当VIF>20时，认为输入变量之间自相关性较为严重。

1.7.2 主成分分析

主成分分析(principal components analysis，PCA)，又称K-L变换，本质上是对数据进行降维，即用很少的不相关变量来表达原始信息，使信息更简单化，是当前应用最为广泛的数据降维方法之一。现实中，为了能描述研究问题，往往选择多个自变量，变量之间却不完全独立，导致所表达的信息重叠，为了简化信息和提升运算效率，需要对数据进行降维处理。以光谱和纹理特征构建AGB估算模型前，首先对光谱信息和纹理特征进行K-L变换，依据特征值大于1的原则，分别选取前n个主成分来建立估算模型。

1.7.3 多元线性回归

多元线性回归(multiple linear regression，MLR)可以同时使用多个自变量来预测因变量，要求自变量与因变量之间要有较好的相关性。分别以光谱特征、纹理特征以及二者结合作为自变量，马铃薯AGB作为因变量来综合构建模型。

1.8 精度评价指标

挑选2/3样本数据(32个，重复1和重复2)作为建模集，1/3样本数据(16个，重复3)作为验证集以此来构建马铃薯AGB估算模型。将采用决定系数(coefficient of determination,R2)、均方根误差(root mean square error,RMSE)、标准均方根误差(normalized root mean square error,NRMSE)来评价模型的精度。R2越接近于1，RMSE和NRMSE越低，其估测模型的精度就越高。

2 结果与讨论

2.1 影像指数与AGB的相关性

将表2选取的影像指数与马铃薯AGB进行相关性分析，将相关系数绝对值按从大到小排列，选取不同飞行高度下前10个影像指数作为模型变量，其结果见表3。由表可知，不同飞行高度下，前10个影像指数与马铃薯AGB均达到0.01显著水平(p<0.01)，相关系数绝对值都在0.606以上。随着飞行高度的增加，VARI，GRVI，NDI，MGRVI，EXGR，g，r+b和EXG共8个影像指数与AGB的相关性均在前10以内，相关性顺序略有变动。当飞行高度为10和30 m时，影像指数EXR和RGBVI相关性都位于前10；当飞行高度为20和50 m时，影像指数RGBVI和g-b相关性也都位于前10；而当飞行高度为40 m时，影像指数g-b和EXR的相关性虽位于前10但都处于后两位。因此可以看出，不同飞行高度下，与AGB相关性较好的前10个影像指数类别差异小，不同分辨率的数码影像对光谱信息组成的影像指数与AGB的相关性排序差别不大。

表3 不同高度下影像指数和AGB的相关性Table 3 Correlation between image indexes of different altitudes and AGB

2.2 纹理特征与AGB的相关性

将提取的不同飞行高度下的影像纹理特征与AGB进行相关性分析，将相关系数绝对值按从大到小排列，选取前8个纹理特征作为模型变量，其结果见表4。由表可知，不同飞行高度下，前8个纹理特征与AGB均达到0.01显著水平(p<0.01)，相关系数绝对值在0.366以上。随着飞行高度的增加，筛选的前8个纹理特征类别差异较大，同种纹理特征在不同飞行高度下与AGB的相关性排序波动明显。当飞行高度为50 m时，纹理特征与AGB相关系数相较于前4种高度比较小，可知纹理特征的提取非常依赖于数码影像的分辨率。

表4 不同高度下纹理特征和AGB的相关性Table 4 Correlation between texture features of different altitudes and AGB

2.3 不同飞行高度下光谱特征估算马铃薯AGB

在构建AGB估算模型前，首先对每种飞行高度下的前10个影像指数进行共线性诊断，判断输入变量之间的共线性。依据各个飞行高度下的共线性诊断VIF>20，得知输入变量之间存在严重共线性，因此有必要对前10个影像指数进行PCA处理。根据主成分特征值大于1的原则，确定不同飞行高度下的最佳主成分为1，1，1，1和2个。基于PCA-MLR方法，得到不同飞行高度下利用光谱特征估算马铃薯AGB的结果如表5所示。由表可知，在飞行高度10～50 m之间，随着影像分辨率降低，光谱信息估算马铃薯AGB的效果整体上呈现下降趋势变化，建模R2从0.69递减到0.65，建模RMSE从226.07 kg·hm-2递增到240.69 kg·hm-2，建模NRMSE从21.22%递增到25.49%；模型验证效果同样呈下降趋势，R2从0.73递减到0.70，RMSE从199.21 kg·hm-2递增到239.25 kg·hm-2，NRMSE从17.74%递增到18.89%。

表5 不同高度下影像光谱信息的AGB估算结果Table 5 Prediction results of AGB from image spectral information with different altitudes

2.4 不同飞行高度下纹理特征估算马铃薯AGB

依据各个飞行高度下的共线性诊断VIF>20，得知输入的纹理特征之间存在严重共线性，因此通过PCA降维处理，确定不同飞行高度下的最佳主成分为1，1，1，1和3个。基于PCA-MLR方法，得到不同飞行高度下利用纹理特征估算马铃薯AGB的结果如表6所示。由表可知，随着飞行高度的增加，纹理特征估算马铃薯AGB的效果逐渐变差，建模R2从0.75递减到0.68，RMSE从208.86 kg·hm-2递增到230.69 kg·hm-2，NRMSE从17.22%递增到21.44%；验证R2从0.83递减到0.73，RMSE从161.65 kg·hm-2递增到201.26 kg·hm-2，NRMSE从16.39%递增到18.99%。以上结果说明，飞行高度的增加会使影像分辨率降低，进而影响纹理特征的提取，造成马铃薯AGB的估算效果不断变差。所有模型中，当飞行高度为10 m及对应的影像分辨率为0.43 cm时，估算效果最好。

表6 不同高度下影像纹理特征对生物量的反演结果Table 6 Prediction results of biomass from image texture features with different altitudes

2.5 不同飞行高度下光谱与纹理特征结合估算马铃薯AGB

将整合后的输入变量进行PCA处理，确定不同飞行高度下的最佳主成分为1，1，1，2和3个。基于PCA-MLR方法，得到不同飞行高度下利用光谱信息+纹理特征估算马铃薯AGB的结果如表7所示。由表可知，随着飞行高度的增加，光谱信息+纹理特征估算AGB效果逐渐变差，建模R2从0.82递减到0.70，RMSE从183.05 kg·hm-2递增到225.85 kg·hm-2，NRMSE从15.61%递增到19.64%；同样验证R2从0.85递减到0.75，RMSE从146.51 kg·hm-2递增到184.29 kg·hm-2，NRMSE从14.36%递增到18.11%。当无人机飞行高度为10 m及对应的地面分辨率为0.43 cm时，估算效果最好。此外，与表5和表6结果相比，同种高度下，通过融合指标构建的AGB估算模型精度最高，效果最佳。

表7 不同高度下影像光谱信息+纹理特征对AGB的估算结果Table 7 Prediction results of AGB from image spectral information and texture features with different altitudes

本研究构建AGB估算模型前，发现当飞行高度为50 m及对应影像分辨率为2.05 cm时，纹理特征与AGB的相关性较低，主要原因是新生叶片的长度和宽度小于2 cm，此分辨率下难以获取完整的影像纹理特征。以光谱信息、纹理特征以及二者结合分别估算马铃薯AGB的结果表明：同种高度下光谱信息+纹理特征估算效果最优，单一纹理特征次之，单一光谱信息最差。此研究结果与刘畅等[11]估算冬小麦AGB时，得出以光谱信息和纹理特征融合估算精度最高，单一纹理特征次之，单一光谱最差结论一致。与陈鹏等[15]估算马铃薯叶绿素含量时，得出结论略有不同，相同之处都是光谱信息+纹理特征融合指标估算精度高于单一光谱、纹理特征指标，不同在于，本工作得出光谱信息估算效果最差，而陈鹏等得出纹理特征估算效果最差。主要原因可能是二者估算的作物理化参数不同，其次是采用的传感器和影像分辨率不同。在探究不同飞行高度下数码影像对马铃薯AGB估算精度的影响时，本研究人为设置了共5种飞行高度并获取相应地影像分辨率，增加了3种变量对AGB估算能力的分析，结果发现同种飞行高度下光谱信息+纹理特征估算AGB的精度相较于单一光谱和纹理特征分别提高了18.84%，11.76%，10.45%，5.97%，7.69%和9.33%，4.13%，2.78%，2.89%，2.94%，表明仅利用数码影像的光谱信息估算马铃薯AGB，模型的精度容易出现饱和现象，而影像的纹理特征不受环境的影响，所得到的模型效果同种分辨率下优于光谱信息，综合考虑了光谱和纹理信息对AGB估算的贡献性，从而结合二者的优势，在一定程度上弥补了单一指标估算AGB能力的不足，提高了模型适用性。本研究还发现无人机数码影像分辨率在0.43～2.05 cm之间变化，光谱信息、纹理信息以及光谱+纹理信息会随影像分辨率提高，其估算AGB的精度逐渐变好，这与贾丹等[10]的研究结果略有不同。贾丹等研究结果表明，随着影像分辨率提高，光谱信息对氮素浓度反演结果影响不大，而基于纹理信息以及光谱+纹理信息的氮浓度反演精度随分辨率增强而提高，主要原因其一可能是实验设置的飞行高度不同，其二研究所用的传感器和作物对象不同，贾丹等研究利用多光谱估算冬小麦氮素浓度，其三本文通过设计氮素和密度梯度试验，会导致每个小区的植被覆盖度发生变化，植株叶片的表现形态和冠层空间结构会明显不同，这会导致作物纹理特征和冠层光谱有所差别。综上可知，对于同种传感器来讲，要想获得高分辨率影像，必须降低无人机飞行高度，但是这样会降低影像采集效率。因此，在对AGB进行监测时，需要根据现实需求并考虑模型精度来选择一定的飞行高度。本研究也存在一些缺陷，仅使用一年的马铃薯块茎增长期的地上部生物量数据进行探讨，得到的模型适用性还需要进一步验证，今后应采集更多年份、生育期及不同作物的样品数据，对影像分辨率如何影响AGB估算结果进行深入分析。