粳稻冠层叶绿素含量PSO-ELM高光谱遥感反演估算

于丰华，冯 帅，赵依然，王定康，邢思敏，许童羽

(1 沈阳农业大学 信息与电气工程学院，辽宁 沈阳 110866； 2 辽宁省农业信息化工程技术中心，辽宁 沈阳 110866)

粳稻是水稻的一种，在我国东北地区有着较大的种植面积。粳稻生育期较长，蛋白质含量较高，食味好[1]。叶绿素含量是表征粳稻生长状态的一个重要性状指标，无损、快速、大面积反演粳稻叶绿素含量一直都是国内外精准农业学者的重要研究方向[2-3]。无人机低空遥感平台是近几年快速发展的一种近地遥感手段，其平台的灵活性高、获取数据丰富、覆盖范围广，相较于其他遥感手段具有一定的优势，因此利用无人机高光谱遥感方式在近地尺度对粳稻冠层叶绿素含量进行反演估算对于辅助田间精准施肥、施药具有重要的现实意义。

目前，国内外在利用高光谱技术估算作物叶绿素含量已经取得了一定的研究成果。尼加提·卡斯木等[4]研究了820～940 nm的最大反射率及反射率综合所建立的偏最小二乘回归法(Patial leastsquares regression，PLSR)估算春小麦叶绿素含量的效果较为理想。毛博慧等[5]研究了冬小麦冠层325～1 075 nm范围的反射光谱，采取遗传算法对光谱特征参数寻优并结合相关分析结果，利用最小二乘－支持向量回归建立冬小麦叶绿素含量预测模型，为后续施肥决策提供支持。孙红等[6]按照叶片垂直分布位置采集马铃薯叶片的成像高光谱数据，选取随机蛙跳优化偏最小二乘(Random frog-partial least squares，RF-PLS)模型计算马铃薯叶片每个像素点的叶绿素含量，实现了不同叶位马铃薯叶片叶绿素含量无损检测以及分布可视化表达。裴浩杰等[7]利用无人机低空遥感平台，通过构建综合长势指标(Comprehensive growth index, CGI)建立小麦关键生育期叶片叶绿素含量反演模型，模型决定系数均大于0.7，可为小麦叶绿素含量长势监测提供参考。

现有的作物叶绿素含量反演研究中，多是通过构建植被指数建立统计回归模型来反演叶绿素含量，这对于特定地区的特定品种，会有比较理想的反演效果，但在模型的普适性等方面还存在一定的不足[8-10]。利用高光谱分析技术检测叶绿素含量方面，前人的工作主要集中在2个方面：建立各种植被指数，利用多元线性或非线性回归方法建立指数与叶绿素含量之间的反演模型；或者先对获取的水稻冠层高光谱数据的全部波段进行PCA、PLS等方法建模[11]。其中，第1种植被指数方法物理意义清晰，模型构建简单，但目前所构建的光谱指数数目繁多，单一植被指数不能很好地表征整个高光谱信息；第2种方法虽然考虑了全部高光谱信息，但每次反演均需要全部波段作为模型输入，计算量较大，模型应用复杂[12]。

本研究针对东北粳稻冠层叶绿素含量反演开展研究。采用无人机高光谱遥感平台，获取粳稻关键生育期冠层高光谱影像数据，采用连续投影算法(Successive projection algorithm，SPA)选择粳稻叶片叶绿素含量响应的敏感波段，建立基于粒子群优化的极限学习机反演模型，以期实现快速、精准、无损诊断粳稻叶绿素含量，提高反演精度。

1 材料与方法

1.1 试验区域与试验设计

为覆盖粳稻关键生育期如分蘖期、拔节期、抽穗期、灌浆期等，试验于2016—2017年的6—9月，每周周三在沈阳农业大学辽中水稻实验站(41°47′N、122°71′E)进行，水稻品种为‘沈稻529’。为避免阴雨多云天气对遥感数据采集造成误差，数据采集过程选择天气晴好的日期，如遇云量超过20%或不利于遥感数据采集天气则顺延采集。粳稻试验小区设计为4个氮肥梯度处理，分别为 CK、N1、N2、N3，每种处理 3 次重复，共设立96个面积相同的长方形试验小区。其中，CK为对照组，即不施用氮肥；N1为当地标准氮肥施用水平，氮肥施用量为45.0 kg/hm2；N2为低氮施肥水平，施用量为22.5 kg/hm2；N3为高氮施肥水平，施用量为67.5 kg/hm2。磷肥和钾肥的施用按照当地标准施用量进行，其中，磷肥标准施用量为51.75 kg/hm2，钾肥标准施用量为18 kg/hm2。图1为本研究试验区域与试验小区设计。

图1 本研究试验区域和试验小区设计Fig. 1 Design of experimental district and plot in this study

1.2 无人机高光谱影像数据获取

无人机高光谱平台采用深圳大疆创新公司的M600 PRO六旋翼无人机，高光谱成像仪选用四川双利合谱公司的GaiaSky-mini内置推扫式机载高光谱成像系统，高光谱的波段范围为400～1 000 nm，分辨率为3 nm，有效波段数为253个，高光谱成像系统如图2所示。

无人机高光谱遥感平台数据采集时间为每次试验的上午11:00—12:00，选择太阳光强较为稳定的时段，无人机飞行高度为150 m。运用ENVI5.3+IDL工具软件对获取的高光谱遥感影像进行小区高光谱数据提取，运用波谱角填图方法去除干扰地物光谱的影响，对每个小区的感兴趣区计算平均光谱，作为每个试验小区的高光谱信息。

1.3 粳稻田间数据获取与叶绿素含量测定

图2 无人机高光谱成像系统Fig. 2 UAV hyperspectral imaging system

粳稻田间冠层数据获取的关键生育期主要包括分蘖期、拔节期、孕穗期、抽穗期、开花期和灌浆期6个时期。粳稻田间取样在试验小区的中间长势均匀区域进行，每个小区采集4穴粳稻作为该小区的样本点，取冠层叶片样本50片放入密封袋，并迅速放入移动冰室进行冷藏保鲜，以确保后期试验数据的准确性。同时利用手持差分GPS设备进行采样点地理坐标位置采集，以便与无人机采集的高光谱影像数据一一对应[7,13-14]。

待粳稻样品送回实验室后，采用浸提法提取样品叶绿素含量。将丙酮、无水乙醇、蒸馏水3种溶液按9∶9∶2的体积比配置成混合溶液，选择粳稻的完全展开叶片，剪碎后混合均匀，称取0.4 g加入200 mL混合溶液；在实验室遮光环境下静置，直到样品完全发白，利用分光光度计进行比色，分别测定D663nm和D645nm，按公式(1)计算粳稻样品的叶绿素含量，ρ(叶绿素)[15-16]。

1.4 数据处理

本研究共进行了10次试验数据采集，共采集780组有效样本数据，选取其中的90%作为建模数据集(样本容量为702)，另外10%为验证数据集(样本容量为78)。本研究采样数据的统计特征如表1所示。

表1 试验小区粳稻叶绿素含量统计特征Table 1 Statistical characteristics of chlorophyll content in japonica rice in experimental plot

由表1可知，本研究的建模数据集和验证数据集2组数据除样本量不同之外，其他统计参量相差不大，变异系数大于40%，说明叶绿素含量数据的离散程度较大。

如图3所示，叶绿素吸收光的能力很强，有2个最强烈的吸收区，一个是波长为640～660 nm的红光部分，另一个是430～450 nm的蓝紫光部分。此外，在光谱的橙光、黄光和绿光部分有不明显的吸收带，所以本研究选择可见光(380～780 nm)之间的光谱信息作为高光谱区间[8,17]。

1.5 高光谱特征波段提取

本研究利用Savitzky-Golay卷积平滑算法[18]对试验小区获取的400～800 nm的反射率光谱数据进行平滑处理，降低噪声对叶绿素含量反演精度的影响，在平滑处理的基础上，利用双变量相关性分析方法对叶绿素含量与高光谱波段反射率之间分别进行相关性分析研究，同时选用皮尔森相关系数、距离相关系数、最大信息系数等3种方法进行单波段相关性分析。

为了能够有效的提取粳稻冠层叶绿素含量高光谱反演特征波段，本研究采用连续投影算法(SPA)提取叶绿素含量反演的特征波段。SPA是一种矢量空间共线性最小化的前向变量选择算法，SPA的流程由3个阶段构成[19]：

第1阶段，筛选出共线性最小的若干组备选波长变量子集。假设初始变量位置k(0)以及变量数N已经给出，该阶段具体步骤如下。

步骤1：第1次迭代之前(n=1)，将训练集光谱矩阵Xcal的第j列赋值给Xj，j=1，…，J；

步骤3：对所有的j，计算Xj在和Xk(n-1)正交的子空间上的投影，

式中：P为投影算子，

步骤 4：令xj=Pxj，j；

步骤5：令n=n+1。如果n＜N则返回步骤2结束：得到N个备选波长的位置：{k(n)；n=0，…，N－1}。

选择过程中进行的投影操作次数为(N－1)(J－N/2)。

第2阶段，分别使用各子集中的变量建立多元线性回归(Multiple linear regression，MLR)模型，选出均方根误差(Root mean square error，RMSE)最小的子集。

第3阶段，对第2阶段选出的子集进行逐步回归建模，在尽量不损失预测准确度的前提下得到一个变量数目较少的集合，该集合中的波长变量即是所选有效波长[20]。

1.6 叶绿素含量PSO-ELM反演建模

本研究采用基于粒子群优化的极限学习机模型(Particle swarm optimization-extreme learning machine, PSO-ELM)反演粳稻冠层叶绿素含量。利用粒子群优化算法优化选择极限学习机的输入层权值和隐含层偏差，从而计算出输出权值矩阵。在本研究中考虑的PSO算法的参数主要包括种群规模(p)、惯性权重(w)、学习因子C1和C2(C1、C2分别调节微粒向个体最优或者群体最优方向飞行的最大步长，决定微粒个体经验和群体经验对微粒自身运行轨迹的影响，C1和C2的搭配不同，将会影响到PSO算法的性能、速度位置相关系数(m)。通过与传统的极限学习机反演算法比较，验证了基于粒子群优化的极限学习机在隐含层节点数目以及网络泛化性上的优势[21]。

本研究使用RMSE和模型决定系数(R2)作为无人机高光谱遥感反演粳稻冠层叶绿素含量精度的评价标准。

2 结果与分析

2.1 粳稻高光谱特征波段提取结果

利用Savitzky-Golay卷积平滑算法对试验小区获取的400～800 nm的反射率光谱数据进行平滑处理，处理后粳稻冠层高光谱曲线如图4所示。

图4 高光谱信息平滑降噪处理Fig. 4 Noise reduction of hyperspectral information smoothing

在平滑处理的基础上，利用双变量相关性分析方法对叶绿素含量与高光谱波段反射率之间分别进行相关性分析研究，本研究选用皮尔森相关系数、距离相关系数、最大信息系数(Maximum information coefficient, MIC)等3种方法进行单波段相关性分析(图5)。

图5 不同方法的单波段相关性分析性结果Fig. 5 Single band correlation analysis results of different methods

400～800 nm之间任意单一波段与粳稻冠层叶绿素含量相关性分析结果表明，基于MIC的相关性要低于基于距离相关性系数法。而利用皮尔森线性相关分析的结果在400～480 nm之间与距离相关系数的结果较为一致。而3种相关性分析方法在701～705 nm处的反射率与叶绿素含量的相关性要明显高于其他波段与叶绿素含量之间的相关性。其中702 nm波段反射率与叶绿素含量之间的相关性最高，3种方法的相关系数分别为0.66、0.64和0.53。由图5可知，400～800 nm波段反射率与叶绿素含量之间的相关系数介于0.50～0.65之间，在区间内不同波段之间也存在着很强的相关关系，难以提取能够有效建立叶绿素含量反演模型的高光谱特征波段[22]。

本研究采用SPA提取叶绿素含量反演的特征波段，结果如图6所示，综合分析来看，随着特征波段提取数目的增加，交叉验证根均方误差 (Root mean square error of cross validation, RMSECV)现出逐渐变小的走势，RMSECV值达到最低时，对应最优波段集合包含的波段数量，即为该时期最优子集波段数量。根据校正集的内部交叉验证RMSECV值确定最佳的光谱波段数为5个。如图7所示，最佳特征波段为 410、481、533、702和 798 nm。

图7 基于连续投影算法(SPA)的粳稻冠层高光谱特征波段选择Fig. 7 Selection of characteristic hyperspectral bands for japonica rice canopy based on successive projection algorithm (SPA)

2.2 叶绿素含量的反演模型

本研究通过SPA提取的粳稻冠层光谱特征波段为自变量，现场采样的粳稻冠层叶片叶绿素含量为因变量，利用ELM、PSO-ELM分别建立粳稻叶绿素含量估算模型，粒子群迭代次数为100。由于PSO算法的参数设置较多，本研究设计了PSO参数正交试验，结果如表2所示。将表2中各参数组合代入PSO-ELM算法(得到叶绿素含量预测值，并计算得出模型决定系数(R2)[22-24]。

表2 ELM和PSO-ELM建立粳稻叶绿素含量反演模型Table 2 Inversion models of chlorophyll contents in japonica rice established by ELM and PSO-ELM

表3为PSO-ELM参数正交试验分析结果，其中Wij为参数j的i水平对应的叶绿素含量反演模型的模型决定系数R2之和，以W11和W21为例，计算公式如下：

表3 PSO-ELM参数正交试验分析结果1)Table 3 Analysis results of orthogonal test of PSO-ELM parameters

其余PSO参量的计算方法同理可得。

由表 3 结果可知，W11、W21、W31、W41和W51之间的差别仅由PSO种群规模的变化而印象最终模型反演精度，与PSO算法的其他参数无关。可通过对比W11、W21、W31、W41和W51的大小来选取最优的PSO算法种群规模，其中R2最大值所对应的种群规模即为本研究所确定的PSO-ELM算法中的种群规模。因此，本研究PSO-ELM反演粳稻叶绿素含量模型的最优种群规模为70，同理可得最优惯性权重w为0.9～0.3线性递减，学习因子C1=2.80、C2=1.10，速度位置相关系数m为0.60。利用该组参数建立的PSO-ELM粳稻叶绿素含量反演模型，模型R2为0.887，RMSE为0.783，要明显优于ELM算法建立的叶绿素含量反演模型(R2=0.734、RMSE=1.824)，反演效果对比如图8所示。

图8 ELM和PSO-ELM模型的粳稻叶绿素含量反演模型检验结果Fig. 8 Test results of the inversion models for the chlorophyll contents in japonica rice using ELM and PSO-ELM

由反演结果可知，本研究采用的PSO优化的极限学习机粳稻叶绿素含量反演模型的预测效果，要明显优于单纯极限学习方法所建立的反演模型，展现了良好的预测能力。

3 讨论与结论

本研究尝试通过SPA方法提取无人机获取的粳稻冠层高光谱信息的特征波段，利用PSO优化的极限学习机算法建立粳稻叶绿素含量反演模型，以期能够实现对东北粳稻冠层叶绿素含量无人机遥感反演的快速监测和评估。传统的利用高光谱特征波段进行一元或者多元回归统计模型可能无法充分表达光谱信息与叶绿素含量之前的非线性数学关系，并在一定程度上制约了叶绿素含量反演精度[25-26]。基于此，本研究采用PSO优化的极限学习机方法，建立PSO-ELM叶绿素含量反演模型，选取了5个特征波段作为模型输入，建立粳稻冠层叶绿素含量反演模型。

所构建的PSO-ELM模型在反演精度上要高于单纯ELM反演精度，这可能是因为PSO算法优化了极限学习机的输入权值和隐层的阈值，避免了ELM在隐层节点数较少的情况下由于随机给定输入权值和隐层阈值带来的泛化能力差、标定精度不高等缺点。本研究所选择的特征波段与前人的结论在特征波段范围区间上基本是相同的，但具体的特征波段并不相同[27-32]。造成这样的原因主要是由于如下几个因素：1)由于试验品种不同造成的特征波段存在一定程度的差别。2)由于稻田内存在的例如水体、杂草、土壤等干扰物，在粳稻纯净高光谱信息提取过程中有可能会参杂上述光谱信息，也会造成特征波段提取的不同。

本研究以沈阳农业大学辽中水稻实验站的‘沈稻529’为研究对象，利用2016—2017年的无人机高光谱粳稻冠层遥感数据，并结合响应的地面粳稻叶绿素含量的实地测量，采用连续投影法提取能够表征粳稻冠层高光谱信息的特征波段，在此基础上通过PSO-ELM算法建立粳稻冠层叶绿素含量反演模型，得到以下结论：1)通过SPA算法提取了能够表征400～800 nm范围的高光谱特征波段，本研究所选取的特征波段分别为410、481、533、702和798 nm。2)通过正交试验确定PSO算法优化ELM模型的参数最优值，其中最优种群规模为80，同理可得最优惯性权重为0.9～0.3线性递减，学习因子1和2分别为2.80和1.10，速度位置相关系数为0.60。3)PSO-ELM粳稻叶绿素含量反演模型的R2=0.887、RMSE=0.783，明显优于单纯极限学习机算法建立的粳稻叶绿素含量反演模型(R2=0.734、RMSE=1.824)。该模型展现了良好的预测能力，能够为粳稻叶绿素含量反演研究及养分诊断提供一定的数据支撑和模型参考。

本研究采用无人机高光谱遥感平台，存在一定的采集误差，同时由于无人机平台的限制地面样本数量还相对有限，并且所建立的叶绿素含量反演模型仅是针对试验粳稻品种建立的，对于其他品种叶绿素含量反演的适用性还需要进一步完善。因此在未来的研究中将增加试验品种，同时针对粳稻不同生育期分别建立叶绿素含量反演模型，提高模型精度和普适性。