土壤有机质的预测方法和设备与流程

本发明涉及计算机地理信息技术领域，更具体地，涉及一种土壤有机质的预测方法以及预测土壤有机质的设备。

背景技术：

土壤有机质作为农业生态系统的重要“指示器”而备受关注，也是环境碳循环和农业管理中空间异质性成因的主要自然和人为驱动力，同时影响着土壤生产力中作物系统的响应过程。为了克服土壤有机质的区域离散性对“空间面源”表达的制约，较高分辨率的土壤养分空间制图已成为其非点源空间变异估计的迫切需求，精确的有机质含量空间评估不仅可以避免不合理采样造成的土壤生态可持续被破坏，也进一步深化丰富着土壤科学的空间变异理论。

传统的土壤有机质空间预测方法已有很多，地统计学与克里格插值因其较好的地学解释能力而被广泛应用，如普通克里格模型的应用。普通克里格模型结构简练易于操作，但受研究区采样点密度、土壤属性空间自相关性及非线性等干扰较大，对土壤有机质空间分布的预测精度较低且泛化性不高。

近年来，混合地统计模型(线性或非线性算法与地统计方法的整合模型)被引入到土壤属性预测领域，例如应用广泛的回归克里格法与人工神经网络克里格法。然而，回归克里格法依然存在对土壤有机质与关联变量间非线性关系解析力弱而导致预测精度不高的问题，人工神经网络法虽然有效地提高了拟合效果，但需要复杂的模型参数调节，并且当制图分辨率提高时所需的运算时间随之增加，导致模型执行效率较低。

技术实现要素：

本发明要解决的技术问题在于提供一种新的土壤有机质的预测方法，能够提高土壤有机质的预测准确性。

根据本发明的一方面，提供了一种土壤有机质的预测方法，包括：

处理器从存储器获取目标区域中被测位置的与土壤有机质相关联的辅助变量数据，其中所述存储器存储有目标区域中各个位置的辅助变量数据；

该处理器将被测位置的辅助变量数据作为预先建立的极限学习机elm模型的输入值计算得到所述被测位置的土壤有机质预测值。

在可选实施例中，辅助变量包括以下变量的任意组合：遥感因子、地形因子、气候因子和土壤属性因子，其中，所述辅助变量至少包括所述遥感因子中的影像波段5和归一化植被指数、所述地形因子中的高程、所述气候因子中的月均温度和月均降水量，以及所述土壤属性因子中的土壤ph、土壤全氮、土壤有效磷和土壤速效钾。

在可选实施例中，所述遥感因子包括影像的波段2至波段7和归一化植被指数，所述地形因子包括高程、坡度、坡向和地形起伏度，所述气候因子包括月均降水量和月均温度，以及所述土壤属性因子包括土壤ph、土壤全氮、土壤有效磷和土壤速效钾。

在可选实施例中，所述预测方法还包括模型建立装置建立极限学习机elm模型，其包括：获取所述目标区域中多个采样点的土壤有机质实测值；根据所述多个采样点的土壤有机质实测值获取所述多个采样点的与土壤有机质相关联的辅助变量数据；根据所述多个采样点的与土壤有机质相关联的辅助变量数据并基于高维空间到低维空间的映射求取所述辅助变量数据组成的矩阵的特征值以及所述特征值的方差，当特征值的方差贡献率累加和大于预设定值时，确定出多个主成分；以及，将这多个主成分的数据作为输入值构建极限学习机elm模型。

在可选实施例中，所述预测方法还包括：模型建立装置建立极限学习机普通克里格elmok模型，其包括：在构建极限学习机elm模型之后，将多个采样点的与土壤有机质相关联的辅助变量数据作为elm模型的输入值计算得到所述多个采样点的土壤有机质预测值；计算所述多个采样点中各个采样点的土壤有机质实测值与土壤有机质预测值之间的残差值；根据各个采样点的残差值按照以下公式计算各个采样点的残差普通克里格插值：

其中，r(xi)是采样点xi的残差值，是残差普通克里格插值，λi是调节权重，n是采样点的数量；以及，通过多个采样点的残差普通克里格插值确定出所述目标区域中各个位置的残差普通克里格插值并存储到存储器。

进一步地，该方法还包括处理器从存储器获取所述被测位置的对应残差普通克里格插值，并且将所述被测位置的土壤有机质预测值与对应的残差普通克里格插值叠加求和得到所述被测位置的土壤有机质预测修正值。

在可选实施例中，目标区域面积在60平方千米到1000平方千米的范围时，采样点的数量范围在20个到600个。

根据本发明的另一方面，还提供了一种预测土壤有机质的设备，包括：存储器，用于存储有目标区域中各个位置与土壤有机质相关联的辅助变量数据；以及，处理器，用于从所述存储器获取所述目标区域中被测位置的辅助变量数据；以及，将所述辅助变量数据作为预先建立的极限学习机elm模型的输入值计算得到所述被测位置的土壤有机质预测值。

在可选实施例中，辅助变量包括以下变量的任意组合：遥感因子、地形因子、气候因子和土壤属性因子，其中，所述辅助变量至少包括所述遥感因子中的影像波段5和归一化植被指数、所述地形因子中的高程、所述气候因子中的月均温度和月均降水量，以及所述土壤属性因子中的土壤ph、土壤碱解氮、土壤有效磷和土壤速效钾。

在可选实施例中，该设备中的存储器，还用于存储有所述目标区域中多个采样点的土壤有机质实测值。设备还包括模型建立装置，用于构建极限学习机elm模型，其包括：从存储器获取所述目标区域中多个采样点的土壤有机质实测值；根据多个采样点的土壤有机质实测值获取所述多个采样点的与土壤有机质相关联的辅助变量数据；根据多个采样点的与土壤有机质相关联的辅助变量数据并基于高维空间到低维空间的映射求取所述辅助变量数据组成的矩阵的特征值以及所述特征值的方差，当特征值的方差贡献率累加和大于预设定值时，确定出多个主成分；以及，将这多个主成分的数据作为输入值构建极限学习机elm模型。

在可选实施例中，模型建立装置还用于建立极限学习机普通克里格elmok模型，其包括：在构建极限学习机elm模型之后，将所述多个采样点的与土壤有机质相关联的辅助变量数据所述elm模型的输入值计算得到所述多个采样点的土壤有机质预测值；计算多个采样点中各个采样点的土壤有机质实测值与土壤有机质预测值之间的残差值；根据各个采样点的残差值按照以下公式计算各个采样点的残差普通克里格插值

其中，r(xi)是采样点xi的残差值，是残差普通克里格插值，λi是调节权重，n是采样点的数量；

通过所述多个采样点的残差普通克里格插值确定出所述目标区域中各个位置的残差普通克里格插值并存储到所述存储器。

在可选实施例中，处理器，还用于从存储器获取被测位置的对应残差普通克里格插值，并且将被测位置的土壤有机质预测值与对应的残差普通克里格插值叠加求和得到被测位置的土壤有机质预测修正值。

通过本发明实施例的土壤有机质预测方法，可以克服现有非线性算法中非线性拟合欠佳、模型参数复杂、模型性能低效的问题，提高预测精度。

附图说明

图1是根据本发明一实施例的土壤有机质的预测方法的流程图。

图2是根据本发明另一实施例的土壤有机质的预测方法的流程图。

图3是根据本发明实施例的elm模型的逻辑结构原理图。

图4是根据本发明实施例的实现土壤有机质预测的设备的示意框图。

具体实施方式

下列公开提供了许多不同的实施例或例子，以实现本发明的不同特征。下面描述了部件和方法步骤的具体例子以简化本公开。当然，这些仅是例子并不旨在限制。应该意识到，可以将这里所公开的内容应用在多种广泛形式中，并且这里所公开的任何特定结构和/或功能仅仅是描述性的。基于这里所教导的内容，本领域技术人员应该意识到，可以独立于任何其他方面实现这里所公开的方面，并且可以以多种方式对这些方面中的两个或多个进行组合。

土壤有机质空间预测辅助变量的选择往往影响到模型拟合效果，更对模型结果驱动力的反馈解释能力至关重要。目前国内外选择的土壤有机质空间制图辅助变量多为环境类因子，如气候、地形等十分便于获取^[5～7]。但由于环境因子在区域尺度的空间变化有限且存在较弱的非线性解释能力，在预测土壤有机质含量时难以达到较高的精度。近年来，遥感波段数据也被用来作为土壤有机质预测的辅助因子并取得一定的预测效果，但是遥感因子难以较好地解释有机质与全球变化及碳汇的关系机理^[8]。因此，在预测模型中同时选取多光谱遥感要素、环境要素以及与土壤有机质密切相关的土壤属性要素(如全氮、有效磷、速效钾等)，共同作为土壤有机质空间预测的多源辅助变量十分必要。但在土壤有机质空间分布预测中，同时考虑环境、遥感和土壤属性要素作为辅助变量的应用很少。

图1是根据本发明一实施例的土壤有机质的预测方法的流程图。如图1所示，该预测方法包括：

s101，模型建立装置从存储器获取目标区域中多个采样点的土壤有机质实测值。

s102，模型建立装置根据多个采样点的土壤有机质实测值从存储器获取多个采样点的与土壤有机质相关联的辅助变量数据。

在本发明实施例中，辅助变量包括以下变量的任意组合：遥感因子、地形因子、气候因子和土壤属性因子。优选地，辅助变量至少包括遥感因子中的影像波段5和归一化植被指数、地形因子中的高程、气候因子中的月均温度和月均降水量，以及土壤属性因子中的土壤ph、土壤全氮、土壤有效磷和土壤速效钾。

存储器存储有目标区域中各采样点的土壤有机质实测值，以及目标区域中各个位置的与土壤有机质相关联的辅助变量数据。

s103，模型建立装置根据多个采样点的与土壤有机质相关联的辅助变量数据确定出建立elm模型用的多个主成分。

具体地，模型建立装置将辅助变量数据基于高维空间到低维空间的映射求取辅助变量数据组成的矩阵的特征值以及特征值的方差，当特征值的方差贡献率累加和大于预设定值时，确定出建立elm模型用的多个主成分。在可选实施例中，与土壤有机质关联的辅助变量可以是5-25个，基于5-25个辅助变量可以确定出2-10个主成分。

s104，模型建立装置将多个主成分的数据作为输入值构建极限学习机elm模型。

s105，处理器从存储器获取目标区域中被测位置的与土壤有机质相关联的辅助变量数据。

s106，处理器将被测位置的辅助变量数据作为elm模型的输入值计算得到被测位置的土壤有机质预测值。

在本发明其他实施例中，还可以是在构建极限学习机elm模型之后，将多个采样点的与土壤有机质相关联的辅助变量数据作为elm模型的输入值计算得到所述多个采样点的土壤有机质预测值。然后，计算多个采样点中各个采样点的土壤有机质实测值与土壤有机质预测值之间的残差值，以及根据各个采样点的残差值按照以下公式计算各个采样点的残差普通克里格插值

其中，r(xi)是采样点xi的残差值，是残差普通克里格插值，λi是调节权重，n是采样点的数量。这样，可以通过多个采样点的残差普通克里格插值确定出目标区域中各个位置的残差普通克里格插值并存储到存储器中。

图2是根据本发明另一实施例的土壤有机质的预测方法的流程图。如图2所示，该预测方法包括：

s201，处理器从存储器获取目标区域中被测位置的与土壤有机质相关联的辅助变量数据。

在本发明实施例中，辅助变量包括以下变量的任意组合：遥感因子、地形因子、气候因子和土壤属性因子。优选地，辅助变量至少包括遥感因子中的影像波段5和归一化植被指数、地形因子中的高程、气候因子中的月均温度和月均降水量，以及土壤属性因子中的土壤ph、土壤全氮、土壤有效磷和土壤速效钾。

s202，处理器将辅助变量数据作为预先建立的极限学习机elm模型的输入值计算得到被测位置的土壤有机质预测值。

s203，处理器从存储器获取被测位置的对应残差普通克里格插值。

s204，处理器将被测位置的土壤有机质预测值与该被测位置的对应残差普通克里格插值叠加求和得到被测位置的土壤有机质预测修正值。

下面结合图1至图3进一步详细描述本公开的土壤有机质预测方法的实施例。

1、确定预测土壤有机质所需的辅助变量数据

1.1测得目标区域的土壤有机质含量

根据本发明实施例，在面积在大约450km²的目标区域内采集多个土壤表层(0～20cm)样本，同时记录gps地理坐标、土壤类型、种植年限、作物种类和地形等相关环境信息，然后采用油浴加热重铬酸钾氧化-容量法测得土壤有机质含量。

1.2确定与目标区域的土壤有机质预测相关联的辅助变量

根据本发明实施例，选取landsat7etm+遥感影像波段数据作为遥感辅助因子。为示例说明的目的，本实施例选取了2009年10月份的landsat7etm+遥感影像波段数据作为遥感辅助因子，其中遥感影像数据来自于美国宇航局(nationalaeronauticsandspaceagency，nasa)，分辨率为30m(其中第八波段为15m)。将etm+波段数据在envi5.2中经辐射校正和大气校正后，选取6个波段(波段2～波段7)，并计算出归一化植被指数(normalizeddifferencevegetationindex，ndvi)，共同作为土壤有机质预测的遥感辅助变量。可选地，也可以选取landsat8oli遥感影像波段数据作为遥感辅助因子。

根据本发明实施例，地形因子可以选取中科院地理空间数据云提供的先进星载热发射和反射辐射仪全球数字高程模型astergdemv2作为高程因子，分辨率为30m。在arcgis10.3中用3d分析工具分别计算出坡度、坡向和地形起伏度共同作为土壤有机质预测的地形辅助因子。

气候因子同样可以选取中科院地理空间数据云提供的modlt1m中国1km地表温度月合成产品，并计算出月均温度和月均降水量，在arcgis10.3中通过重采样转换成与dem一致的空间分辨率。

土壤属性因子包括土壤ph、土壤全氮、土壤有效磷和土壤速效钾。其中，土壤ph通过套筒式甘汞电极测得，土壤全氮采用半微量开氏法测定，土壤有效磷和速效钾分别通过0.5mol·l^-1nahco3和1mol·l^-1nh4oac的火焰光度法测定。将以上点因子数据分别进行反比距离加权插值法(inversedistanceweighted，idw)得到空间分布数据并转换分辨率为90m。

2、建立elm模型

2.1描述性统计分析：

对面积在大约450km²的目标区域内采样点的土壤有机质及17个辅助变量进行描述性统计，在以正负3倍标准差为基准剔除异常值后，共统计例如522个采样点值。发明人在研究中发现如果一定面积区域内的采样点数量太少，则所构建的elm模型的预测精度不好，即土壤有机质的预测值不精确；如果采样点数量太多，则构建elm模型的计算量明显增加，然而预测精度并不会提高很多，因此选择合适的采样点数量可以使得elm模型的预测精度较好同时减少不必要的计算量和计算复杂度。在一优选实施例中，面积为100km²的目标区域，采样点的数量在20-100个；在另一优选实施例中，面积为500km²的目标区域，采样点的数量在100-600个；在又一优选实施例中，面积为1000km²的目标区域，采样点的数量在600-1000个。

在本发明实施例中，如表1所示，通过522个采样点的土壤有机质值得到目标区域的土壤有机质的含量介于6.00～51.00gkg^-1之间，均值为24.74gkg^-1，经过spss22.0软件的柯尔莫诺夫-斯米尔诺夫(kolmogorov-smirnov，k-s)检验显示检验结果p值(即样本符合正态分布的概率大小，默认当p>0.05时，原样本的检验结果符合正态分布)为0.322，符合正态分布。如表1所示，522个采样点的17个辅助变量中，坡度、坡向、地形起伏度、土壤全氮、土壤速效钾和土壤有效磷的变异系数差异较大，波段2到波段7、ndvi、高程、月均温度、土壤ph和月均降水量的变异系数较低。总的来说，地形因子和土壤属性因子的变异程度较高，遥感辅助因子次之，气候因子的变异程度最低，分析其原因可能与研究区内土地覆盖类型多样、土壤酸度整体较高等因素有关。

表1采样点土壤有机质及辅助变量的描述性统计结果

在本发明实施例中，对目标区域所有采样点的土壤有机质实测值以及目标区域中各个位置(包括采样点)的土壤有机质关联的辅助变量数据进行描述性统计以用来判断辅助变量和土壤有机质这些数据的可用性及相关性。具体而言，如表1所示，可通过最大值、最小值、方差、变异系数等统计指标，发现并剔除异常值。通过初步分析辅助变量变异程度较大的原因并进行异常值剔除，可以避免变量间的数值异常对后续主成分分析造成干扰。

将进行描述性统计分析后得到的采样点的土壤有机质实测值和土壤有机质关联的辅助变量数据存储到存储器中。在本发明实施例中，在确定出目标区域内采样点的土壤有机质关联的辅助变量数据后，可以进一步确定出其他非采样位置点的土壤有机质关联的辅助变量数据，例如，将研究区划分为90m×90m的格网(研究区面积450平方千米可划分为约55556个网格)，将17个辅助变量数据(如遥感因子、高程等面状图像数据)在arcgis软件中统一成90m大小的格网，就得到除去522个采样点之外(约55034个网格)的非采样位置点格网对应的辅助变量数据。选择90m左右×90m左右的格网可以尽可能将每个采样点落到相应的网格中，不会出现一个网格里有两个采样点的情况，此外还可以提高计算效率，这是因为格网越小，计算量将会成指数型增长。

2.2主成分分析降维处理：

在本发明实施例中，经过描述性统计的异常值剔除，采用主成分分析法(principalcomponentanalysis，pca)能减少变量间多重共线性的现象而且还能使模型结构简单化。jolliffeit.principalcomponentanalysis[m].springer-verlag,1986.此外，根据所述辅助变量数据和方差贡献率确定出八个主成分，可以较好地保留原始辅助变量的信息。采用主成分分析法是由于发明人在研究中发现：如果存在变量间高度线性相关，就会导致提取的主成分过度集中与那几个线性相关的变量，干扰对原始变量的综合提取，提取出来的主成分总贡献率较低也不能很好地反映原有变量信息。

在本发明实施例中，对17个多源辅助变量进行主成分分析，首先将17个辅助变量原始数据按行排列(一个辅助变量数据为一行)，组成新矩阵x，并对矩阵x进行标准化使其均值变为零；然后，求解矩阵的协方差，就得到k个特征向量，将每个特征向量对应的特征值大小由大到小排列，将k个特征向量按行组成矩阵h；最后，根据公式y＝hx计算，得到降维后数据y，并得到相应的方差贡献率。通常取方差总贡献率(每个特征值的方差贡献率累积之和即总贡献率)大于90％的主成分，能较好地保留多源辅助变量的原始信息。本发明中，将17个辅助变量在spss22.0软件中自动完成主成分分析过程，如表2所示，根据方差贡献率逐个主成分进行累加，取方差总贡献率大于90％得到8个主成分，对原始变量的解释程度达到91.46％。其中，第一主成分pc1解释了总变量的30.86％，波段2、波段3、波段4、波段5和波段7的相关性最高，而ndvi则在第五主成分，和月均降水量、月均温度具有较高的相关性，解释了总变量的6.41％；其他6个主成分(pc2＝18.72％、pc3＝11.79％、pc4＝8.97％、pc6＝5.52％、pc7＝5.09％和pc8＝4.10％)主要反映出土壤属性、气候、地形等环境辅助变量的密切关联性。所得的8个主成分可用于作为输入值建立elm模型。

表2主成分分析结果

3、elm模型建模

基于上述处理后的多个主成分据来构建elm模型以对目标区域的土壤有机质进行空间预测。

极限学习机(extremelearningmachine，elm)算法是由新加坡南洋理工大学黄广斌教授于2004年提出的一种由输入层、隐含层和输出层组成的单隐层前馈神经网络(single-hidden-layerfeedforwardnetworks，slfns)算法。huanggb,zhuqy,siewck.extremelearningmachine:anewlearningschemeoffeedforwardneuralnetworks[j].proc.int.jointconf.neuralnetw,2004,2:985-990.本公开发明人在研究中发现：elm算法在本质上属于神经网络算法体系，但其在算法信息传播过程中相比于传统的人工神经网络有着较大差异，首先在算法执行过程中它可以随机设定输入层与隐含层之间的权值和阈值，而且对网络的输入权值以及隐元的偏置无需反复调整迭代，最后产生唯一的最优解，避免了人工神经网络算法容易陷入局部最优解的情况，因此具有超快的学习速度和强泛化性。

elm的广义slfns输出函数，即公式(1)所示：

式中，x是模型的输入数据，fl(x)是模型的输出，β＝[β1,…,βl]^t是l个节点到m个输出节点(m≥1)之间的权向量集合，h(x)＝[h1(x),…,hl(x)]是elm隐含层输出行向量的非线性特征映射，hi(x)是第i个隐含层节点的输出函数，它包含常见的激活函数。图3示出了elm模型的逻辑结构示意图，其中隐含层节点的输出函数可以有很多种，不同的输出函数用于不同的隐含层神经元。elm算法最终推导出最优权值求解范式，即公式(2)所示：

式中，β*是elm模型隐含层与输出层之间的最优权值，是moore-penrose广义逆矩阵h，即公式(3)所示：

t是训练数据的期望输出矩阵，即公式(4)所示：

由上述elm算法推导过程可知其根据连续概率分布随机产生隐含层和输入层之间的权值和阈值，而第二阶段的隐含层和输出层之间的权值转化为求解一个线性方程组，并且无需后续调整，因而相比于传统的人工神经网络其高学习效率和泛化性能更佳。

elm模型结构的这种单隐含层结构不同于传统的线性方程式表达，属于“黑箱”模型，模型预测的计算过程被隐含，本发明在matlab2013b软件中通过编程能够实现elm模型的计算过程。

在matlab2013b软件中通过编程，将表2中8个主成分作为elm模型的输入，通过公式(1)及图3“黑箱”过程的计算得到土壤有机质预测值。如表3所示，elm预测模型的最终结构为8-32-1(8个主成分输入即x1＝pc1,x2＝pc2,…,x8＝pc8，32个隐含层节点即l＝32，以及1个输出节点即土壤有机质som预测值)，激活函数为sigmoid核函数。

为了评估上述elm模型的预测精度和性能，选取平均误差(meanerror，me)即公式(5)，均方根误差(rootmeansquareerror，rmse)即公式(6)和相关系数(coefficientofdetermination，r²)即公式(7)作为评判依据。

式中，z(xi)是土壤有机质实测值，是elm模型的土壤有机质预测值，是土壤有机质实测值的平均值。其中，me和rmse越接近于零，r²越接近于1表明模型预测精度较高。

根据上述公式(5)、公式(6)和公式(7)计算elm模型的精度误差，得到其对土壤有机质预测的me和rmse分别为7.92gkg^-1和2.012gkg^-1，同时对目标含量值与辅助变量之间的非线性解释能力r²达到55.1％。作为对比，还可以根据上述公式(5)、公式(6)和公式(7)计算人工神经网络ann预测模型的精度误差，结果发现相比于ann预测模型，elm模型土壤有机质预测的me和rmse较低。

此外，由于在elm模型中初始权值随机赋值、无需权重调节及迭代等算法参数优势使得elm预测模型的运算时间大大快于传统的ann模型。

表3最佳elm土壤有机质预测模型参数及误差统计

4.不同辅助变量对土壤有机质预测的影响力分析

为了进一步探究17个辅助变量对土壤有机质预测的影响力，直接用原始的遥感类型因子、环境因子和土壤属性因子作为elm模型的输入，去预测土壤有机质含量。在matlab2013b软件中通过编程，将表1中17个原始辅助变量作为elm模型的输入，通过公式(1)及图3“黑箱”过程的计算得到土壤有机质预测的精度。经过逐步变量加入得到各因子对土壤有机质预测的精度大小序列，由表4所示，经elm模型的预测的17个辅助变量对土壤有机质的关系影响力从大到小依次为：土壤有效磷>土壤速效钾>月均温度>月均降水量>土壤全氮>高程>土壤ph>波段5>ndvi>波段4>坡向>坡度>地形起伏度>波段2>波段3>波段7>波段6。因此，17个辅助变量按对土壤有机质预测影响力较大的因子包括土壤有效磷、土壤速效钾、月均温度、月均降水量、土壤全氮、高程、土壤ph、波段5和ndvi等9个辅助变量，这些变量不同大小的驱动作用共同制约着土壤有机质的量效特征，在实际辅助因子选取过程中应对其加以着重考虑。

表4不同辅助变量对土壤有机质预测的精度统计结果

为进一步提高土壤有机质的预测精度，本发明实施例还提供了一种基于极限学习机普通克里格elmok模型的土壤有机质预测方法的实施例。

普通克里格插值是地统计学的一种空间局部插值方法，其求解结果无法通过明确的数学方程式进行表达，本发明实施例可以通过arcgis10.3软件的地统计分析模块完成对研究区土壤有机质的空间预测过程。

本发明实施例采用elmok模型预测土壤有机质含量方法包括以下三步：

首先，如上面所描述的方法构建土壤有机质与提取的主成分之间的elm模型。例如可以在matlab2013b软件中进行编程从而构建elm模型。然后基于elm模型将采样点的辅助变量数据作为elm模型的输入值计算得到采样点的土壤有机质预测值；

接着，根据公式(8)和公式(9)分别计算土壤有机质实测值与基于elm模型的土壤有机质预测值之间的残差值及残差普通克里格值；

其中：

在式(8)中，r(xi)是采样点xi的残差值；z(xi)是采样点xi土壤有机质实测值；是基于elm模型的采样点xi的土壤有机质预测值。

在式(9)中，r(xi)是采样点xi的残差值，是采样点xi的残差普通克里格插值，λi是调节权重，n是以采样点xi为圆心半径为例如2.04km内的圆中存在的已知采样点的数量。具体而言，通常假设以土壤有机质预测值未知的采样点x(这个采样点的土壤有机质实测值已知)为圆心设定一个半径h为例如2.04km(这个半径根据研究区实际情况可以在arcgis10.3软件的地统计分析模块中确定)形成一个圆，在这个圆的范围内会有n个土壤有机质预测值已知的采样点，用这n个已知采样点和未知采样点x可以建立一个一阶或多阶的线性方程。通过未知采样点x的土壤有机质实测值，就可以求出一组调节权重λ，进而用这些权重去预测目标区域内其他未知采样点的土壤有机质含量。权重的调节优化、函数模型选取等过程本领域普通技术人员能够通过arcgis10.3软件的地统计分析模块自动完成。

然后，根据公式(10)将基于elm模型的土壤有机质预测值和残差普通克里格结果求和最后得到采样点的土壤有机质修正预测值：

在式(10)中，是基于elm混合地统计模型的采样点xi的土壤有机质修正预测值；是采样点xi的残差普通克里格插值；是基于elm模型的采样点xi的土壤有机质预测值。

在本发明实施例中，还可以根据公式(5)、(6)和(7)评估本发明实施例提出的极限学习机普通克里格(elmok)模型、以及作为对比方法的普通克里格(ok)模型、多元线性回归普通克里格(mlrok)模型和人工神经网络普通克里格(annok)模型的土壤有机质预测精度和性能。

由表5可以得出，elmok模型相比于传统的普通克里格以及其他混合地统计模型具有较低的rmse(1.801gkg^-2)和较高的r²(0.662)，表明其具有更加精确的预测水平。从4种模型的土壤有机质预测值与真实值范围比较来看，普通克里格法预测结果与真实值范围差异最大，局部区域有机质含量分布的临界跳跃性较大，其它3种模型整体的预测效果与实测值较为接近，而elmok模型的土壤有机质预测结果范围与实测值最为接近，整体的预测效果较佳，对土壤有机质含量具有更深的非线性解析能力。

表5不同模型下土壤有机质预测的精度统计结果

在本发明实施例中，对于522个采样点，可以将其分成380个训练样本和142个测试样本，其中380个训练样本的土壤有机质实测值以及土壤有机质关联的辅助变量数据用于构建elm模型和elmok模型，另外142个测试样本的土壤有机质实测值以及土壤有机质关联的辅助变量数据用来验证elm模型或elmok模型的预测精度。

现有的土壤有机质空间预测方法主要包括诸如多元线性回归模型等的经典统计学模型、诸如简单克里格法、普通克里格法等的地统计学模型以及诸如回归克里格模型、神经网络克里格模型等混合地统计学模型三类。其中，回归克里格模型、神经网络克里格模型这种混合地统计学模型，通过非线性机器学习算法的高拟合精度结合克里格插值的高效性，是目前土壤有机质预测领域较为先进的技术方案。

然而，传统的地统计学方法对土壤有机质空间分布的预测精度有限，如普通克里格模型虽然结构简练易于操作，但难以解决土壤有机质与辅助变量间复杂非线性关系引起的低拟合度问题。此外，常用的混合地统计学方法如回归克里格法依然存在对土壤有机质与关联变量间非线性关系解析力弱而导致预测精度不高的问题，人工神经网络克里格法虽然有效地提高了拟合效果，但需要复杂的参数调节，并且当制图分辨率提高时所需的运算时间随之增加，导致模型执行效率不高。

本发明实施例提出的采用极限学习机普通克里格组合模型，并结合环境因子、遥感因子及土壤属性等多源辅助变量，挖掘区域尺度土壤有机质空间分布的新预测方法，其相比于现有技术具有以下优点及有益效果：

(1)土壤有机质含量预测的极限学习机普通克里格(extremelearningmachineordinarykriging，elmok)模型对比普通克里格(ordinarykriging，ok)、多元线性回归普通克里格(multiplelinearregressionordinarykriging，mlrok)和人工神经网络普通克里格(artificialneuralnetworkordinarykriging，annok)模型具有更低的均方根误差(rmse＝1.801gkg^-2)及更高的相关系数(r²＝0.662)，对研究区土壤有机质含量与多源辅助变量之间的非线性关系具有较强的捕捉及解析能力，模型拟合精度较高且稳定性较好；

(2)极限学习机普通克里格模型相比人工神经网络普通克里格模型的参数调节简单、不易陷入局部最优解、模型结构简单且执行效率较高，十分便于土壤有机质预测及土壤肥力评估。

根据本发明实施例的城市更新改造的选址方法的设备或装置可以是计算机、或与计算机相关的实体等。该设备进一步可以包括处理器，配置成用于实施根据图1或图2示出的土壤有机质的预测方法中的各动作。如图4所示，根据本发明一个实施例的预测土壤有机质的设备40，包括：存储器42，用于存储有目标区域中各个位置与土壤有机质相关联的辅助变量数据；以及，处理器44，用于从存储器42获取目标区域中被测位置的辅助变量数据；以及，将该辅助变量数据作为预先建立的极限学习机elm模型的输入值计算得到被测位置的土壤有机质预测值。

在本发明实施例中，辅助变量包括以下变量的任意组合：遥感因子、地形因子、气候因子和土壤属性因子。辅助变量至少包括所述遥感因子中的影像波段5和归一化植被指数、地形因子中的高程、气候因子中的月均温度和月均降水量，以及土壤属性因子中的土壤ph、土壤碱解氮、土壤有效磷和土壤速效钾。

在可选实施例中，存储器42，还用于存储有所述目标区域中多个采样点的土壤有机质实测值。设备40还包括模型建立装置46，用于构建极限学习机elm模型。模型建立装置46包括从存储器42获取目标区域中多个采样点的土壤有机质实测值；根据多个采样点的土壤有机质实测值从存储器获取多个采样点的与土壤有机质相关联的辅助变量数据；以及根据多个采样点的与土壤有机质相关联的辅助变量数据并基于高维空间到低维空间的映射求取所述辅助变量数据组成的矩阵的特征值以及所述特征值的方差，当特征值的方差贡献率累加和大于预设定值时，确定出多个主成分。进一步地，模型建立装置46将多个主成分的数据作为输入值构建极限学习机elm模型。

在另一实施例中，模型建立装置46还用于建立极限学习机普通克里格elmok模型。具体而言，模型建立装置46在构建极限学习机elm模型之后，将多个采样点的与土壤有机质相关联的辅助变量数据作为elm模型的输入值计算得到多个采样点的土壤有机质预测值，然后计算多个采样点中各个采样点的土壤有机质实测值与土壤有机质预测值之间的残差值；以及，根据各个采样点的残差值按照以下公式计算各个采样点的残差普通克里格插值

其中，r(xi)是采样点xi的残差值，是残差普通克里格插值，λi是调节权重，n是采样点的数量。

最后，处理器44通过多个采样点的残差普通克里格插值确定出目标区域中各个位置的残差普通克里格插值并存储到存储器42。

在本发明实施例中，处理器44，还用于从存储器42获取被测位置的对应残差普通克里格插值，并且将被测位置的土壤有机质预测值与对应的残差普通克里格插值叠加求和得到被测位置的土壤有机质预测修正值。

通过本发明实施例的设备预测目标区域各个位置的土壤有机质可进一步参考本公开前面土壤有机质的预测方法中相对应步骤的描述，在此不再赘述。应当理解，可以以硬件、软件、固件、中间件、代码或其任何恰当组合来实现这里描述的实施例。对于硬件实现，处理器可以在一个或多个下列单元中实现：专用集成电路(asic)、数字信号处理器(dsp)、数字信号处理设备(dspd)、可编程逻辑器件(pld)、现场可编程门阵列(fpga)、处理器、控制器、微控制器、微处理器、设计用于实现这里所描述功能的其他电子单元或其组合。

上面概述了几个实施例的特征使得本领域技术人员可较好地理解本公开的方面。本领域技术人员应当理解他们可容易地使用本公开作为基础以设计或修改其他工艺和结构以实行相同目的和/或实现在此介绍的实施例的相同优点。本领域技术人员也应意识到这种等同构造没有脱离本公开的精神和范围内，并且他们在没有脱离本公开的精神和范围情况下可以做各种改变、代替和更改。