彭金莲，唐朝胜，蒋菊生，赵春梅

(1. 海南大学 计算机与网络空间安全学院， 海口 570228；2. 海南省农垦科学院， 海口 570206；3.中国热带农业科学院 橡胶研究所 海口 571101)

我国耕地地力评价主要方法有模糊数学法、层次分析法、指数和法等，我国第2次土壤普查多采用指数和法评定耕地地力等级。近年来，陈桂芬等[1]应用K-means聚类方法、Johnson粗糙集属性约简算法与C4.5决策树算法优化地力评价。蔡丽霞等[2]应用C4.5决策树、K-means和DBSCAN聚类算法，对土壤养分数据进行对比分析，依靠大量土壤肥力状况以预测未来土壤肥力的变化趋势。闫一凡[3]认为评价指标筛选目前尚无国际通用的标准，中国已建立了包括7类共64项指标的公用指标体系;当前主流的地力评价方法(如中国农业部推荐的特尔斐—层次分析法等)仍存在人为随意性等不足，推荐将不同评价方法(特别是特尔斐—层次分析法和分类与回归树模型法的组合)组合起来用于实际工作，可取长补短，进一步提高结果的可靠性。建议建立起统一的地力等级定义，引入分类与回归树、人工智能、机器学习等新方法，并将其与现有评价方法(如特尔斐法、多元回归等)集成，充分挖掘数据、实现优势互补，认为是未来耕地地力评价研究的重点和发展趋势[4-7]。海南农垦测土配方施肥项目从2006年开始试点，历时10年，实施范围为海南农垦原有的92个农场(分10个片区)，面积大约2万hm2，依据测土配方调查数据进行了橡胶园耕地地力评价。采用了专家打分法和特尔菲法，10个片区分别建立了耕地地力评价指标体系，10个片区采用的耕地地力评价指标不尽相同，耕地地力评价指标的选择与重要度由专家打分确定，人为主观影响明显，无法对海南垦区胶园耕地地力进行统一评价。笔者对海南垦区胶园大数据进行挖掘，旨在筛选出对海南胶园耕地类型有影响的地力指标，以便客观评价海南垦区胶园的耕地地力。

1 研究对象与技术方法

1.1 研究对象本研究通过整合海南农垦测土配方施肥10个片区(项目场)的属性数据库、空间数据库和测土配方施肥数据库，创建全垦区统一的耕地地力评价数据库；对海南垦区测土配方施肥原始数据进行数据探索，通过数据分析，探索土壤内在的分类规律，筛选出与垦区耕地地力相关的主要因素，以确定海南垦区的耕地地力指标；使用4种聚类模型，利用轮廓系统进行适用性评估，最终确定一种聚类算法作为耕地地力分类模型来划分胶园耕地类型，从而实现对海南橡胶园耕地地力进行统一的综合评价。

1.2 技术方法使用IBM公司开发的企业级大型数据挖掘工具SPSS Modeler进行数据挖掘，技术路线参照跨行业数据挖掘标准流程 CRISP-DM 模型设计而成，此数据挖掘标准流程由问题理解、数据理解、数据准备、建立模型、模型评估和成果发布等6个阶段构成。

1.2.1 分类型指标的数据预处理通过对原始数据的审核，发现能够作为地力评价潜在的分类指标为:地貌类型、地形部位、坡向、成土母质和质地等。但由于原始数据时间跨度大、涉及人员多，许多数据没有严格按照操作规程填写，导致各指标的类别数超出规程规定，产生了地貌类型71类，地形部位95类，坡向71类，成土母质41类，质地64类，所以数据质量较差，因此,必须按照农业部测土配方施肥野外调查操作规程要求调整。其中，地形部位原来只有其他部位和山腰，已经没有数据挖掘价值。调整前后各指标及其类别值见表1。调整后，地貌类型调整为3种类别，分别是山地、丘陵、平地。坡向调整为9种类别，分别是北、东北、东、东南、南、西南、西、西北和平地等。成土母质调整为8种类别，分别是安山岩、变质岩、花岗岩、浅海沉积物、砂页岩、石灰岩、玄武岩和云母片岩风化物。质地调整为6种类别，分别是粘土、粘壤土、中壤土、砂壤土、砂土和其他等。分类型地力指标的缺失值处理依据不同情况进行了不同的处理。地形部位只有4个片区农场有调查数据，缺失值占74%，已经没有数据挖掘价值。地貌类型的缺失值为11.44%，只有1个片区农场没有调查此类数据，通过C&RT模型进行缺失值插补。质地指标60%是其他类别，已经没有分类价值。

表1 地力评价各分类指标调查数据调整前后类别数

1.2.2 连续型地力指标的预处理连续型地力指标的预处理相对分类型指标来说要复杂很多，包括离群值、极值、缺失值的数据预处理，也包括数据分布和数据分箱处理。连续型地力指标共有10个，分别是pH值、有机质、全氮、有效磷、速效钾、北纬、东经、海拔、常年平均降雨量和有效土层厚度。从描述性统计分析，查看均值、标准差、最大值和最小值。一般情况下，若标准差远远大于均值，可粗略判定数据存在异常值。采用3倍标准差上最大值替换离群值和丢弃极值的方法。 为了提高数据挖掘效率，基于对连续型地力指标进行离散化处理和适应某些模型只能处理分类型变量的要求，减少变量取值数以实现样本量的缩减。参考《海南岛土系概论》[8]，对海南农垦胶园土壤肥力的养分指标进行分组，如表2所示,分组给出的区间范围是大于左边的区间数和小于右边的区间数。在SPSS Modeler中对连续数据分组的方法有固定宽度、分位数、等级、平均值/方差和最优等5种方法，根据调查样本的数据分布，确定海拔高度和有效土层厚度采用样本量固定宽度分级，年平均降雨量、东经和北纬采用数值固定宽度分组，各立地指标分组标准如表3所示，分组给出的区间范围是大于左边的区间数和小于右边的区间数。

表2 土壤养分指标分组

表3 立地指标分组标准表

1.2.3 聚类模型算法聚类分析将一批样本数据，在没有先验知识的前提下，根据数据的诸多特征，按照其在特征上亲疏自动分组，使得组内个体的结构特征具有较大的相似性，组间个体特征相似性较低。聚类与分类不同，聚类分析的输入数据集是一组未标记的对象，也就是说此时输入的对象还没有被进行任何分类，聚类的目的是根据一定的规则合理地进行分组或聚类，并用显式或隐式的方法描述不同的类别。由于分析可以采用不同的算法，所以对于相同的数据集合可能有不同的划分。从聚类的原理角度来看，聚类的算法可以分为:划分聚类、层次聚类、基于密度的聚类和基于网格聚类等方法。本研究分别使用K-Means，Two-Step，Kohonen和TwoStep-AS等4种聚类算法，对海南农垦测土调查数据进行聚类分析，最终，根据轮廓系数(Silhouette Coefficient)的大小确定TwoStep-AS聚类模型为海南垦区耕地进行分类。

1.2.3.1 TwoStep-AS聚类模型算法介绍[9]SPSS Modeler中的TwoStep-AS聚类算法，它的优势至少表现在以下几个方面：1)可同时基于类别变量和连续变量进行聚类；2)可自动确定最终的分类个数；3)可处理大型数据集。

TwoStep-AS聚类基本步骤如下：

1)预聚类：采用“贯序”方式将样本粗略划分成L个子类。开始阶段视所有数据为一个大类。读入一个样本数据后，根据“亲疏程度”决定该样本应派生出一个新类，还是应合并到已有的某个子类中。这个过程反复进行，最终形成L个类。预聚类过程聚类数目不断增加。

2)聚类：在预聚类的基础上，再根据“亲疏程度”决定哪些子类可以合并，最终形成L类。聚类数目不断减少的过程，随着聚类的进行，类内部的差异性将不断增大。

在距离测量有“欧氏距离”和“对数似然”两种算法作为聚类变量相似度的测量形式；如果聚类变量均为数值型，采用欧氏距离，若为混合型采用对数似然距离。在聚类准则有贝叶斯(BIC)和Akaik信息准则(AIC)，作为聚类个数的判断依据。同其他统计方法一样，TwoStep-AS聚类算法也有严苛的适用条件，它要求模型中的变量独立，类别变量是多项式分布，连续变量符合正态分布。

注意，与其他的聚类分析算法一样，所得到的模型可以在一定程度上依赖于训练数据的顺序。重新排序数据和重建模型可能会导致不同的最终群集模型。对记录顺序的解决方案的鲁棒性可以通过多次拟合一个模型来评估，每次使用不同的随机顺序记录。

评价模型质量的方法采用轮廓系数(Silhouette Coefficient)[10-14]，是聚类效果好坏的一种评价方式。最早由Leonard kaufman，Peter J.Rousseeuw 在1986年提出。它结合内聚度和分离度两种因素。可以用在相同原始数据的基础上用来评价不同算法，或者算法不同运行方式对聚类结果所产生的影响。

对数据中每一个数据样本与本结果簇中其他样本间的相似性以及该数据样本与其他结果簇中样本的相似性进行定量分析，并将定量结果进行某种形式的组合，得到聚类结果优劣的评价。

对于第i个样本，计算该样本到其所属簇中所有样本的平均距离，记为a(i)，对于第i个样本和不包含该样本的任意簇，计算该样本到给定簇的中所有样本的平均距离，所有簇的平均距离最小者记为b(i)。对于第i个样本，其轮廓系数的计算如下公式(1)所示。

(1)

式中：a(i) = average(i向量到所有它属于的簇中其他点的距离)；b(i) = min(i向量到所有非本身所在簇的点的平均距离)。

a(i)越小，说明样本i越应该被聚类到该簇。b(i)越大，说明样本i越不属于其他簇。

轮廓系数s(i)=+1时，表示 样本i与其他簇中的对象相异性较大, 聚类合理。

轮廓系数s(i)=0时，表示样本i分类不明显,在两个簇的边界上。

轮廓系数s(i)= -1时，表示样本i被分配到一个错误的簇中，轮廓系数的值在-1和+1之间变化。

所有样本的s(i)的均值称为聚类结果的轮廓系数，是该聚类是否合理、有效的一种度量方式。

1.2.3.2 评价指标的重要度算法 考察地力指标的重要性，可使用Pearson样本相关系数、似然比(Likelihood Ratio)卡方值和Cramer’V系数。SPSS Modeler计算Cramer’sV系数[9]，见公式2。

(2)

式中：χ2是Pearson卡方统计量，R和C分别表示列联表的行数和列数。Cramer’s V系数取值范围在0～1，越接近1，越说明输入变量与输出变量的相关性越强，输入变量越重要。重要度在0～ 0.2 范围为差，重要度在0.2～ 0.6范围为中， 重要度在0.6～1范围为良。

2 结果与分析

2.1 地力指标的筛选结果初次运算，输入地力特征参数包括成土母质、地貌类型、坡向、pH值、全氮、速效钾、有机质、年平均降水量、北纬、东经、有效土层厚度、海拔等在内的14个指标，特征重要度排名前10的指标如图1所示。年平均降雨量特征重要度为1，全氮特征重要度为0.67，有效土层厚度特征重要度为0.58，纬度特征重要度为0.58，成土母质特征重要度为0.51，有机质特征重要度为0.46，经度特征重要度为0.36,地貌类型特征重要度为0.21,坡向特征重要度为0.18，有效磷特征重要度为0.11。

图1 初次运算的地力指标的特征重要度

根据特征重要度的大小，选择特征重要度>0.2以上的地力指标作为模型聚类的特征值，经过反复迭代筛选，以轮廓系数最优为依据，最终确定地力指标为年平均降雨量、全氮、有效土层厚度、纬度和成土母质(图2)，年平均降雨量重要度为1，纬度的重要度为0.9，全氮的重要度为0.75，有效土层厚度的重要度为0.59，成土母质的重要度为0.57。从图2可以看出，影响海南胶园地力的5个指标中，土壤方面有3个，分别是全氮、成土母质、有效土层厚度，这3个因子一起综合在一定程度上代表了不同土壤的养分总量，因有机质含量与全氮有较大的相关性，选择特征重要度>0.2重新迭代时就隐含了；地理方面有1个指标，就是纬度，这是一个影响气温和热量的十分重要的因子；气象方面有一个指标，就是年均降雨量，影响的是土壤水分。可见，水、热、养分，作为热带树种橡胶树来说，是其生长发育和高产高效的前提条件，也就是构成地力的重要因子。

图2 最终确定的地力指标的特征重要度

2.2 TwoStep-AS模型聚类结果根据最终确定的地力指标(常年降雨量、纬度、全氮、有效土层厚度和成土母质)建立TwoStep-AS聚类模型。模型参数信息如下：最小常规聚类数2类，最大常规聚类数15类。自适应特征选择，信息标准采用贝叶斯信息标准 (BIC)，距离测量采用对数似然。最终模型常规聚类数为3个，即将海南垦区十大片区的耕地划分为3种类型。3种耕地类型的概要特征如表4所示，每种耕地类型最重要的特征在表中以粗体显示。第1种耕地类型最重要的特征是年平均降雨量在1 600～1 700 mm，全氮质量含量在0.5～0.65 g·kg-1,有效土层厚度在70～100 cm，成土母质为玄武岩。第2种耕地类型最重要的特征是年平均降雨量大于2 000 mm，全氮质量含量在0.05～0.10 g·kg-1,有效土层厚度在55～70 cm，成土母质为花岗岩。第3种耕地类型最重要的特征是年平均降雨量在1 700～1 800 mm，全氮质量含量在0.65～0.85 g·kg-1,有效土层厚度在70～100 cm，成土母质为花岗岩。

表4 3种耕地类型的聚类中心

注： 最重要的特征用粗体表示

Note: The most important characteristics are in bold.

2.3 模型评估在SPSS Modeler中，提出了2种度量聚类效果的指标，一个是优度，另一个是重要度。优度是一种测量聚类内聚性和分离性的指标。优度在-1～ 0.2范围为差，优度在0.2～0.5范围为中，优度在0.5～1范围为良。重要度是一种测量聚类内聚性的指标。重要度在0～0.2 为差，重要度在0.2～0.6范围为中， 重要度在0.6～1范围为良。轮廓系数是聚类效果好坏的一种评价方式，轮廓系数在-1～0.2范围为差，在0.2～0.5范围为中，在0.5～1范围为良。用TwoStep-AS聚类算法为基于海南农垦测土调查数据进行分类时，总体模型优度(平均轮廓系数)为0.38。各聚类的优度和重要度如表5所示。从表5的数据来看，各聚类类型的优度为中等，而重要度为良好，说明用TwoStep-AS聚类算法对海南农垦测土调查数据进行聚类时内聚性一般，而分离性很好。

表5 模型质量

3 结 论

1)本研究使用基于大数据的TwoStep-AS聚类模型对海南农垦测土调查数据进行聚类分析，将海南垦区十大片区划分为3种耕地类型。第1种耕地类型最重要的特征是年平均降雨量在1 600～1 700 mm，全氮质量含量在0.5～0.65 g·kg-1范围,有效土层厚度在70～100 cm范围，成土母质为玄武岩。第2种耕地类型最重要的特征是年平均降雨量大于2 000 mm，全氮质量含量g·kg-1在0.05～0.10 g·kg-1范围,有效土层厚度在55～70 cm，成土母质为花岗岩。第3种耕地类型最重要的特征是年平均降雨量在1 700～1 800 mm，全氮质量含量在0.65～0.85 g·kg-1范围,有效土层厚度在70～100 cm，成土母质为花岗岩。2)影响海南胶园耕地的最重要的地力指标分别是年平均降雨量、全氮、北纬、有效土层厚度、成土母质。3)基于数据挖掘的地类划分没有人为主观因素的干扰，与作物产量无关，以耕地本身特征内在的凝聚和分离为依据，排除了人为主观因素，这种分类将更有助于人们对耕地性质的认识，有利于土地的区划和综合利用。4)本研究结果为建立统一的海南垦区胶园地力评价指标体系提供了初步的方法和结果，并为下一步开展系统性的研究打下了基础，其结果的准确性和适用性还有待今后在深入研究和应用中得到检验。