一种GSA-MLs相结合的小麦LAI估算方法

技术领域

本申请属于作物叶面积指数反演方法技术领域，尤其是涉及一种结合全局敏感性分析与机器学习的小麦LAI估算方法。

背景技术

叶面积指数(LAI)是冬小麦生长指标监测过程中重要的冠层结构参数，也是冬小麦估产模型中重要的输入参数。传统LAI测量方法依赖于人工定点观察或仪器测量，效率低下且耗费人力、物力、财力，难以应用于大面积的冬小麦生长监测。随着对地观测技术的发展，富含丰富光谱信息的区域乃至全球尺度的影像获取已经成为可能，为大范围冬小麦LAI的估算提供了数据基础。

众多卫星平台已经搭载多光谱谱段的传感器，提供可用于提取植被的长势信息的可见光-近红外波段。如Landsat 8卫星、MODIS卫星、Worldview-2卫星、RapidEye卫星等包含可见光-近红外波段或短波红外波段的影像。此外，相关研究表明，红边信息在植被卫星遥感监测过程中具有较大优势(参见DONG et al.Remote Sensing of Environment,2018,222,133-143)。如Worldview-2与RapidEye等商业卫星在可见光-近红外波段的基础上均增加了一个红边波段，结果证明在植被生长监测方面具有较好的效果(参见

在利用遥感影像估算植被LAI的过程中，常见的方法有经验或半经验模型、机器学习算法等。其中，经验或半经验模型往往是通过影像的波段反射率构建植被指数，增强植被信息，通过统计回归分析构建LAI的估算模型，实现LAI的预测；机器学习算法则是基于不同的模型框架，根据训练数据集进行学习，构建LAI估算模型，常见的有偏最小二乘法回归、高斯过程回归、支持向量回归、随机森林回归等。基于植被指数的经验或半经验模型通常具有较好的机理性、可操作性强，但易出现病态反演的问题。机器学习算法较强的数据解析能力，能够避免病态反演的问题且具有较高的精度，但模型的机理解释性尚待加强。

已有研究通过融合两类算法，在原始波段反射率的基础上，将多个植被指数作为机器学习算法的输入变量，改善植被长势参数的估算模型(参见LI,X C,etal..International Journal of Applied Earth Observation and Geoinformation,2016,44:104-112.)。但是过多输入变量的引入会增加计算负荷，降低大范围LAI估算效率，同时数据冗余的风险加大。因此，输入变量的筛选是必要且重要的。虽然部分机器学习算法自带变量筛选的功能，如PLSR、RF，但是该筛选功能仅聚焦于输入变量与目标变量之间的相关性。然而，在植被LAI的遥感估算方面，植被的光谱表现是一项综合指标，不仅仅受到LAI的影响，同时也会受到叶绿素、叶倾角、干物质等其他参数的影响。因此需要具有同时评估多个参数对光谱影响大小的方法来弥补机器学习的不足。

利用敏感性分析定量模型的各个输入参数对输出结果影响大小，显得尤为重要。常见的敏感性分析方法有两类，包括全局敏感性分析(GSA,Global SensitivityAnalysis)和局部敏感性分析(LSA,Local Sensitivity Analysis)；前者是对所有输入参数评估及分析对最终结果的贡献程度，常见的有扩展傅里叶振幅灵敏度试验(EFAST,Extended Fourier Amplitude Sensitivity Test)、Sobol法等；后者则在保持其他输入变量不变的情况下，逐个分析每个输入变量对输出结果影响的大小。

当前，已有学者尝试结合全局敏感性分析与植被辐射传输模型(PROSAIL模型等)开展植被长势参数对光谱变量影响大小的评估，尤其是多个长势参数之间的交互作用影响大小的评估，但是应用Sentinel-2卫星结合全局敏感性分析与植被辐射传输模型，并结合机器学习进一步应用于大范围冬小麦LAI的估算鲜有报道。

发明内容

为了提高LAI估算的精度和MLs的计算效率本发明提供了基于全局敏感性分析(global sensitivity analysis,GSA)与MLs(简称GSA-MLs)相结合的小麦LAI估算方法。

具体，本发明提供了一种估算小麦叶面积指数(LAI)的方法，包括以下步骤：

S1：基于PROSAIL模型模拟数据集，利用全局敏感性分析(GSA)量化不同光谱变量对植被生长参数的响应情况；

S2：基于不同光谱变量对植被生长参数响应的一阶灵敏度，利用S

S3：选择最优变量作为机器学习算法(MLs)的输入参数；

S4：用MLs整合筛选得到的光谱变量结合地面实测小麦LAI建立小麦LAI估算模型；

S5：通过计算决定系数

S6：通过最佳估算方法绘制了研究区LAI的分布图；

其中，S

S

S

S

优选地，上述方法中，S1所述光谱变量是采用Sentinel-2多光谱影像数据，通过获取Sentinel-2影像数据经辐射定标与大气校正，然后利用Sen2Res 1.0将Sentinel-2影像中20m空间分辨率的波段降尺度至10m得到。更进一步的，获取的Sentinel-2影像数据经过官方发布的Sen2Cor 2.5.5进行辐射定标与大气校正。随后利用Sen2Res 1.0将Sentinel-2影像中20m空间分辨率的波段降尺度至10m。

Sentinel-2包含两颗卫星，分别是2015年发射的Sentinel-2A与2017年发射的Sentinel-2B，两者均搭载多光谱传感器，能够获取具有13个不同空间分辨率波段的多光谱影像。

S1所述PROSAIL模型设置的参数包括：结构系数，叶绿素含量，类胡萝卜素含量、等效水厚度、褐色素的含量、平均叶倾角、叶面积指数、热点、土壤亮度指数、太阳天顶角、观察者天顶角以及方位；参数在确定范围后经马尔科夫链-蒙特卡罗方法(MCMC)进行采样。

进一步的，S1所述植被生长参数包括基于可见光-近红外构建的非红边植被指数，例如NDVI、RVI、DVI、GDVI、GRVI、GNDVI、EVI、SAVI、CIg等，加入红边波段后的红边植被指数，例如NDRE、RERVI、REDVI、REGDVI、REGRVI、REGNDVI、MEVI、SARE、CIred-edge等，以及基于短波红外波段构建的植被指数，例如NDWI、MSI等。

由于上述红边植被指数与短波红外相关的植被指数均只涉及一个红边波段或一个短波红外波段，而Sentinel-2影像具有三个红边波段与两个短波红外波段，因此，优选地，在计算红边植被指数与短波红外相关植被指数过程中，所述的红边植被指数分别用Sentinel-2的三个红边波段构建三次，而所有基于短波红外波段构建的植被指数分别用Sentinel-2的两个短波红外波段构建两次。

本发明所述的植被指数计算公式可以为本领域已知的计算公式，例如文献TUCKERC J.Remote Sensing and Environment,1979,8(2):127-150；SIMS,D A等.RemoteSensing of Environment,2002,81(2):337-354；MAIK,B等.Sensors,2019,19(24):5397；VINCINI M,FRAZZI E等.Precision Agriculture,2009；JORDAN C F等Ecology,1969,50(4):663-666；BUSCHMANN C等.International Journal of Remote Sensing,1993,14(4):711-722等公开的NDVI、NDRE、RVI、RERVI、DVI、REDVI、GDVI、REGDVI、GRVI等植被指数的计算公式，这些文献公开的内容通过引用结合到本发明中。

上述方法中，优选地，S1所述全局敏感性分析的方法包括EFAST、Sobol和Mirror，优选地为EFAST。

进一步的，S2所述最优变量是基于GSA所得不同光谱变量对植被生长参数响应的一阶灵敏度(First order sensitivity)结果(S

进一步的，S3所述机器学习算法包括采用偏最小二乘法(PLSR)、支持向量(SVM)和随机森林(RF)算法建立冬小麦LAI估算模型。

其中，偏最小二乘法(Partial Least Squares Regression，PLSR)综合了多元回归分析、典型相关分析和主成分分析的思想，能够在自变量存在严重多重相关性、样本点个数少于变量个数的条件下进行回归建模。

支持向量机(Support Vector Machine,SVM)是一种以非线性映射为理论基础的小样本机器学习方法。它避开了从归纳到演绎的传统过程，实现了从训练样本到预测样本的高效“转导推理”。径向基核函数被认为是具有良好分类功能，是最常用的核函数。本发明选用SVR_Epsilon模型、高斯径向基核函数(RBF)，通过调节拉格朗日乘上界，不敏感损失函数的参数、相对误差参数，实现模型最优解。SVR中的参数用交叉验证法获得，同时为防止“过学习”，对参数C进行适当调整。

随机森林算法(Random forests，RF)模型是建立在决策树基础上的一种集成学习方法，通过多次bootstrap抽样获得多个随机样本，并通过这些样本分别建立相对应的决策树，从而构成随机森林。

进一步的，S5所述决定系数

其中，n表示样本数；χ

其中，N表示样本数，LAI

进一步的，S6所述最佳估算方法为通过S

术语解释：

LAI：小麦叶面积指数；

MLs：机器学习算法

GSA：全局敏感性分析(global sensitivity analysis)

PLSR：偏最小二乘回归；

SVM：支持向量机；

RF：随机森林；

GSA-MLs方法：全局敏感性分析与机器学习算法相结合的方法。

本发明提出了结合全局敏感性分析与机器学习提高小麦LAI估算精度、效率的方法(GSA-ML)。首先进行光谱变量的筛选，在提高机器学习法估算LAI精度以及应用过程中计算机计算效率的同时，提高了机器学习应用过程中的机理性。综合对比LAI估算精度、计算机运行效率，结果显示GSA-RF具有最佳的表现。

本发明运用PROSAIL模型进行全局敏感性分析，不仅能够考虑目标参数对光谱变量的影响，还能将目标参数以外其他参数以及不同参数之间交互作用影响纳入考虑。既能获取单个参数对模型输出变量的敏感性，也能获取参数单独作用及其与其他参数之间的交互作用对模型输出变量的敏感性。

本发明提供的GSA-ML，利用全局敏感性分析综合评估现有的、敏感于LAI的光谱变量，通过四种不同的策略(SLAI、SLAI+SCab、SLAI-SInteraction和SLAI+SCab-SInteraction)筛选出对LAI敏感性高且对其他参数敏感性低的光谱变量，而后用机器学习法整合筛选得到的光谱变量来估算小麦LAI。试验证明，本发明选用的光谱变量大部分对LAI具有较好的敏感性。通过对比四种不同的光谱变量筛选策略以及不同机器学习法在LAI估算方面的表现，结果显示，相对于51个光谱变量的运行时间缩短了54.13％，有效地避免了数据冗余造成的精度影响与计算效率低的风险。

本发明提供的GSA-ML不仅适用于卫星平台，还适用于无人机与地面平台，尤其是普遍面临数据冗余风险的地面高光谱平台。

附图说明

下面结合附图和实施例对本申请的技术方案进一步说明。

图1是利用GSA-ML估算小麦LAI的方法的流程图；

图2是本发明实施例研究区及地面调查样点分布图；

图3是不同光谱变量全局敏感性分析结果图；

图4是在冬小麦LAI估算方面表现最佳的GSA-PLSR(a)、GSA-SVM(b)和GSA-RF(c)图。

图5是由GSA-RF生成的小麦LAI图。

具体实施方式

以下结合具体实施例对本发明的技术方案以及有益效果做进一步说明，应当理解的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。

实施例

1、研究区域

本实施例以研究区位于江苏省高邮与姜堰地区为例(如图2)。该区域的年均温度在15℃左右，年均降水量在1 032.3mm左右，土壤类型主要为壤土，研究区内冬小麦为扬麦158，在2017年10月下旬播种，于2018年6月成熟收获。试验于2018年4月19日开展，选择了46个冬小麦田样点作为地面调查的对象，该田块为稻麦轮作与秸秆还田，且是在低机械化的条件下播种，甚至人工播种。

2、Sentinel-2影像获取及预处理：

Sentinel-2包含两颗卫星，分别是2015年发射的Sentinel-2A与2017年发射的Sentinel-2B，两者均搭载多光谱传感器，能够获取具有13个不同空间分辨率波段的多光谱影像，波段信息如表1所示。

本实施例根据地面采样时间获取2018年4月19日的Sentinel-2卫星影像，获取的Sentinel-2影像数据经过官方发布的Sen2Cor 2.5.5进行辐射定标与大气校正。随后利用Sen2Res 1.0将Sentinel-2影像中20m空间分辨率的波段降尺度至10m。本实施例究尚未涉及用于大气成分监测且分辨率为60m的三个波段，包括B1、B9、B10。

表1 Sentinel-2波段信息

3、小麦LAI估算方法

在PROSAIL模型的基础上引入GSA，测试不同光谱变量对不同作物长势参数的响应情况，然后在GSA所得结果的基础上采用四种变量筛选的策略并将筛选得到的光谱变量应用于实测数据，作为机器学习估算小麦LAI过程中的输入变量，形成了用于小麦LAI估算方法，即GSA-ML，具体流程如图1所示。

(1)、PROSAIL模型模拟数据

为了定量不同长势参数对同一光谱变量影响大小并筛选出用于LAI估算的敏感变量，利用PROSAIL模型进行数据模拟。其中，模型输入参数的设置如表2所示。PROSAIL模型中的参数在确定范围后经马尔科夫链-蒙特卡罗方法(MCMC)进行采样，总计获取8 000个样本，用于后面EFAST的敏感性计算。

表2 PROSAIL模型的参数设置

(2)植被指数

用于植被LAI估算的植被指数，包括NDVI、RVI、DVI、GDVI、GRVI、GNDVI、EVI、SAVI、CIg等基于可见光-近红外构建的非红边植被指数、NDRE、RERVI、REDVI、REGDVI、REGRVI、REGNDVI、MEVI、SARE、CIred-edge等加入红边波段后的红边植被指数以及NDWI、MSI等基于短波红外波段构建的植被指数。在计算红边植被指数与短波红外相关植被指数过程中，所有的红边植被指数分别用Sentinel-2的三个红边波段构建三次，所有基于短波红外波段构建的植被指数分别用Sentinel-2的两个短波红外波段构建两次。

表3用于估算小麦LAI的植被指数

(3)全局敏感性分析筛选光谱变量策略

首先基于PROSAIL模型模拟数据，利用全局敏感性分析中的EFAST比较不同光谱变量对不同长势参数响应的敏感性分析。其中，马尔科夫链-蒙特卡罗方法被用于模型输入样本的生成，根据表2中各输入参数的范围与模型设置的采样参数N进行采样(本实施例中设置N＝1 000)，总计获取8 000个样本用于不同光谱变量对各参数响应的敏感性分析。然后，基于EFAST所得不同光谱变量对各参数响应的一阶灵敏度(First order sensitivity)结果(Sparamete)，利用以下四种不同策略进行机器学习输入变量的排序，分别为：

策略一：仅根据LAI的敏感性分析结果S

策略二：同时考虑LAI与Cab两者敏感性的总和，根据S

策略三：同时考虑高LAI敏感性与低交叉互作影响，根据S

策略四：为避免参数之间的交叉互作影响，同时考虑高LAI敏感性、高Cab敏感性与低交叉互作影响，根据S

根据以上四种变量筛选策略，分别筛选前10、20、30个变量用于机器学习估算小麦LAI，并对比不同策略估算小麦LAI的精度以及影像应用过程中计算机的计算速度。本实施例基于Python3.7环境下的Sklearn构建机器学习算法，计算机型号为联想ThinkStationP520,操作系统为Windows10，其GPU为NVIDIA Quadro P6000(32Gb)，CPU为Inter(R)Xeon(R)W-2125CPU(32Gb)。

(4)机器学习

偏最小二乘法(Partial Least Squares Regression，PLSR)综合了多元回归分析、典型相关分析和主成分分析的思想，能够在自变量存在严重多重相关性、样本点个数少于变量个数的条件下进行回归建模。基于PLSR建立的LAI光谱预测模型能有效地减少光谱维数，揭示最大LAI变化的主控因子，建立模型有更好的稳定性。

支持向量机(Support Vector Machine,SVM)是一种以非线性映射为理论基础的小样本机器学习方法。它避开了从归纳到演绎的传统过程，实现了从训练样本到预测样本的高效“转导推理”本申请选用SVR_Epsilon模型、高斯径向基核函数(RBF)，通过调节拉格朗日乘上界，不敏感损失函数的参数、相对误差参数，实现模型最优解。SVR中的参数用交叉验证法获得，同时为防止“过学习”，对参数C进行适当调整。

随机森林算法(Random forests，RF)模型是建立在决策树基础上的一种集成学习方法，通过多次bootstrap抽样获得多个随机样本，并通过这些样本分别建立相对应的决策树，从而构成随机森林。RF是通过随机采样点方式进行数据训练，因此能够有效避免过拟合，在对组分进行高效预测。本实施例应用Python中Scikitlearn库建立随机森林模型，回归树数量ntree为600，每棵树随机抽取特征变量mtry为6，其余参数选择默认设置。

(5)模型构建及精度评估

利用2017-2018年地面实测数据进行模型的构建，通过计算决定系数

其中，n表示样本数；χ

其中，N表示样本数，LAI

4、试验结果：

(1)不同光谱变量估算小麦LAI的表现

从相关性分析结果表4看，在各波段与LAI相关性中，红波段R、红边波段RE1、RE2以及近红外NIR、NNIR波段光谱反射率与LAI具有显著相关性，且五个波段的相关系数大都高于0.450，说明红边波段及近红外波段在冬小麦长势监测中具有很好的有效性。同时相关性最高的为近红外NIR波段，相关系数为0.643，RMSE为0.979。大部分植被指数与LAI具有很强的相关，其中相关系数均高于0.500的植被指数有10个，分别为SARE2、NDRE2、REDVI2、GDVI、CIre2、RERVI2、Cig、GRVI、GNDVI、SAVI、RVI、MTVI2，其中SARE2最佳，相关系数为0.707，其次为NDRE2，相关系数为0.681。由于RE-2相对于其他两个红边波段(包括RE-1，RE-3)位于红边区域的中心，大多数由RE-2构建的红边植被指数相对于无红边波段或由其他红边波段构建植被指数具有对LAI更高的敏感性。另外，NDRE相较于NDVI来说，可以克服在高LAI值下的饱和问题，同时对作物叶片水分和叶绿素含量表现出较弱的敏感性。总体来看，加入红边波段的植被指数与传统近红外植被指数相比，相关性表现出微弱的提升。

表4、单波段反射率及植被指数与LAI相关性

(2)不同光谱变量全局敏感性分析的结果

利用EFAST测试了PROSAIL模型中不同作物长势参数对Sentinel-2多光谱信息的影响情况，包括原始的波段反射率及其衍生出的多个植被指数，如图3所示。结果显示，大多数光谱变量都敏感于LAI，但是受到其他参数不同程度的影响，影响较大的有叶绿素、叶倾角以及参数之间的耦合作用影响。其中，平均叶倾角(Lidfa)虽然对红边波段(包括RE-1、RE-2、RE-3)及NIR具有明显的影响，但是通过植被指数的构建能够较大程度地降低这类影响；叶绿素(Cab)虽然对RE-2影响较小，但是植被指数的构建会放大其影响，包括3.1部分所得最佳表现的NDRE2与SARE2；而由SWIR波段及其构建的植被指数虽然不受叶倾角与叶绿素的影响，但是会遭受叶片含水量(Cw)的影响。综合以上情况，大部分植被指数在估算小麦LAI时会遭受其他参数影响，需要机器学习法整合多个变量的LAI敏感性来提高小麦LAI的估算精度。

(3)光谱变量排序结果

基于各长势参数对不同光谱变量影响大小的结果，即全局敏感性分析所得一阶灵敏度结果，利用四种变量排序的策略，包括SLAI、SLAI+SCab、SLAI-SInteraction和SLAI+SCab-SInteraction，进行光谱变量的排序(表5)，并分别筛选前10、20与30个光谱变量作为不同机器学习法(PLSR、SVM、RF)的输入变量，对比LAI估算的精度与机器学习法运算的效率。

结果显示，在不考虑Cab的(Strategy 1和Strategy 3)情况下，由于G、RE-1、RE-2波段构建的植被指数遭受Cab的影响具有较低的LAI敏感性未排列在前列，而大部分非红边植被指数由于遭受到的Cab影响较小排在前列(图4)；在同时考虑Cab与LAI敏感性(Strategy 2和Strategy 4)的情况下，部分由RE-2构建红边植被指数排在前列。值得注意的是，这两个策略排列与筛选得到的光谱变量组合在通过机器学习估算LAI过程中存在受植被叶绿素影响造成模型精度或普适性低的风险；除了参数各自敏感性影响外，参数之间的交互影响也是降低LAI估算模型精度低的一大原因。因此，本研究进一步考虑低参数互作影响(Strategy 3和Strategy 4)的情况。另外，在四种筛选策略中，由G波段构建的MTVI2表现较为稳定，这是因为它能够有效地降低Cab的影响并保留LAI的敏感性，然而在LAI估算过程中存在估算精度不如红边植被指数的表现(表4)。

表5不同策略光谱变量排序结果

(4)不同机器学习估算小麦LAI的变现

根据上述四种筛选策略的排序结果，选择前10、20、30个光谱变量作为机器学习法的输入变量，结合地面实测小麦LAI，分别采用偏最小二乘法(PLSR)、支持向量(SVM)和随机森林(RF)算法建立冬小麦LAI估算模型，结果如表6所示。

在不同GSA-ML对比方面，GSA-SVM，GSA-RF整体表现效果最佳，R

表6冬小麦LAI估算模型对比

表7 MLs法估算小麦LAI的比较

应用例：

基于上述对比结果，综合考虑LAI的估算精度与计算机运算效率，选用了变量数为10的三种GSA-ML进行进一步的对比。利用所有变量和基于全局敏感性分析获取的策略S

表8不同估算模型应用过程中的计算机运行时间(s)

进一步将LAI估算方面表现最佳的模型(RF)应用姜堰地区，基于4748×4231像素的Sentinel-2影像，绘制LAI的空间分布(图5)，期间机器学习的计算时长为618.920s。所得LAI的空间分布均处于合理值范围，且同实际生产具有一致性，说明本发明提供的方法在应用过程中具有可行性。

以上述依据本申请的理想实施例为启示，通过上述的说明内容，相关工作人员完全可以在不偏离本项申请技术思想的范围内，进行多样的变更以及修改。本项申请的技术性范围并不局限于说明书上的内容，必须要根据权利要求范围来确定其技术性范围。