基于模糊c均值聚类法的玉米农田管理分区研究

陈世超 杜太生 王素芬

(中国农业大学中国农业水问题研究中心， 北京 100083)

0 引言

“精准农业”概念从提出至今已逐渐成熟。同时大面积农田管理分区技术的应用也为农业生产过程中的精准调控提供了理论基础[1]。国内外众多学者基于农田土壤、作物的空间变异规律划分管理分区，进行分布式管理[2-4]，根据土壤要素与作物的空间变异性，对农田进行管理区域划分，以达到“精准管理”的目的。应用该方法可以明确各分区之间的土壤、作物等要素的差异，因地制宜，可显著提高水肥利用效率和管理精度，减少资源浪费[1,3]。

大面积农田内土壤类型、养分以及地形属性都存在空间变异性，这些因素会影响产量，从而造成产量分布不均匀[2]。国内外学者对农田进行分区管理的研究大都集中于土壤养分聚类分析。张泽等[3]和温鹏飞等[4]对土壤养分指标和作物指标进行统计分析与相关分析，利用模糊c均值聚类法进行聚类，基于GIS和RS技术最终实现管理分区；陈彦等[5]和BEHERA等[6]使用FCM进行管理分区，将分区结果作为变量施肥和养分精确管理的作业单元；GUASTAFERRO等[7]利用迭代自组织数据分析法(Iterative self-organizing data analysis)、FCM和非参数密度算法(Non-parametric density algorithm)对研究区域进行管理区划分，3种算法分区后子区域内各要素空间变异性减小，达到了分区后各因素均质化的目的，同样文献[8-10]使用模糊c均值聚类法得出了相似的结论；对于流域尺度的农田管理分区，周浩等[11]和杨建宇等[12]分别使用FCM和TOPSIS算法，划分了作物的适宜种植区和农田管理区。这些研究通过分析土壤养分含量以及金属阳离子交换量等多要素进行分区，阐述了分区后的均质性，但并未在时间尺度上考虑生育期内不同分区的土壤含水率变化和作物生长规律。同时，使用单一要素进行管理区划分的研究也有很多成果，表观电导率(ECa)作为主要的评价指标可以反映多种土壤阳离子的含量，可以作为管理区划分的依据，其分区结果也较为合理[13-15]，但单一要素反映产量空间变异性的信息少于多因素聚类分析法所得到的结果。

目前，农田分区后对不同区域土壤含水率(SWC)变化、作物生长规律以及产量的研究还鲜有报道，基于此，本文筛选对产量起决定作用的主控因子，基于主成分分析提取可反映各因子变异性的主成分，用模糊c均值聚类方法得到分区结果，分析分区后各因子的空间变异规律以及各分区的作物形态指标、土壤含水率在时间尺度上的差异，以期为农田分区管理提供科学依据和理论指导。

1 材料与方法

1.1 研究区概况

试验于2016年4—10月在甘肃省武威市凉州区黄羊河农场七队2号地(102°55′E，37°49′N)进行，研究区面积为6.84 hm2(120 m×570 m)。该区属典型的大陆性温带干旱荒漠气候，昼夜温差大，冬冷夏热、干燥、多风。全年日照时数3 000 h左右，无霜期150 d左右，年均气温8.8℃，大于0℃的年积温为3 550℃以上。多年平均降水量为163.2 mm，多年平均水面蒸发量为2 019.9 mm。灌溉水源为地下水，地下水埋深约40 m。

研究区供试作物为制种玉米，于2016年4月20日种植，9月5日收获。在研究区中划分30 m×30 m的均匀网格，另加设15个加密点，共91个取样观测点，用于取土测定土壤各项指标和作物产量，取样分布点如图1所示。灌溉施肥方式均采用膜下滴灌，使用“一膜两管四行”的铺设方式，滴灌带间距约75 cm，株距25 cm，植株与滴灌带的距离为15 cm，滴头间距30 cm，滴头流量3.0 L/h，滴灌带长度约为60 m。

图1 研究区取样点分布示意图Fig.1 Spatial distribution map of sampling points and elevation in study area

1.2 主要观测内容与方法

1.2.1土壤质地

播种前，在各取样点取土，风干后过2 mm筛，使用MaterSizer2000型激光颗粒分析仪(Malvern Instruments Ltd.，英国)测定土壤砂粒、黏粒、粉粒含量[16]，测量深度为0～20 cm、20～40 cm、40～60 cm。

1.2.2土壤含水率

播种前，在各取样点取土使用干燥法测定，测量深度为0～20 cm、20～40 cm、40～60 cm；生育期内，使用土壤水分廓线仪(Diviner 2000型，Sentek Pty Ltd.，澳大利亚)测定，测量深度为100 cm，每20 cm一层，并用取土干燥法进行校准，测定频率为7～10 d，灌水前后和降雨后加测；同时安装土壤水分传感器(EC-5型)进行连续监测，探头埋深分别为20、40、60 cm(图1)。

1.2.3土壤速效氮含量

播种前，在各取样点取土，测定速效氮含量，测量深度为0～20 cm、20～40 cm、40～60 cm，使用氯化钾浸提法(2 mol/L KCl溶液)用AutoAnalyzer-Ⅲ型流动分析仪(SEAL Analytical Gmbh, 德国)测定[17-18]。

1.2.4土壤电导率

播种前，在各取样点取土，测定土壤电导率，测量深度为0～20 cm、20～40 cm、40～60 cm，采用1∶5的土水比法用SevenCompact S230型电导率仪(Mettler Toledo, 美国)测定[19]。

1.2.5作物生理指标和产量

玉米株高和叶面积指数(LAI)：在各取样点为中心划定2 m×2 m矩形范围内选择3株玉米，于关键生育期(苗期、拔节期、抽穗期、灌浆期、成熟期)用米尺测定玉米株高，用Sunscan冠层分析仪测量玉米LAI并以手工测量法进行校正。

产量：作物成熟后在各取样点为中心划定2 m×2 m的矩形范围内选择10株玉米，人工脱粒后在85℃条件下干燥至恒定质量并称量，换算成单位面积产量作为该取样点所代表网格的产量。

1.2.6取样点位置与地形属性

使用全球定位系统 GPS(Trimble Recon，美国)对取样点进行空间定位，取样点分布与高程图由ARCGIS 10.2生成，同时使用栅格数据计算、提取各取样点坡度、坡向数据。

1.3 主成分分析

影响作物产量的各因子之间存在相关性，使得各因子对产量的影响作用存在信息重叠。为了真实地反映各因子对产量的影响，需要将重复的因子进行压缩，达到“降维”的目的。主成分分析核心思想是降维，通过正交变换将一组可能存在相关性的变量转换为一组线性不相关的变量，转换后的这组变量即为主成分，主成分之间互相独立且包含原始变量中尽可能多的信息，依据各主成分贡献率判断其对因变量的影响程度[18]。

1.4 模糊聚类分析

在模糊聚类中，一个样本可以同时归属于不同的子集，其对不同子集归属程度通过隶属度来表示，隶属度越大，表示归属于某个子集的可能性就越大，反之则越小[20-22]。本文中，使用管理分区软件(Management zone analyst, MZA, version 1.0.1, 美国)对测定的各因子进行聚类分析[17]。MZA软件的核心是模糊c均值聚类算法，使用模糊性能指数(FPI)和归一化分类熵(NCE)作为评判聚类效果的评判标准，当FPI和NCE同时最小时的分区数即为最优[20,23]。

2 结果与分析

2.1 各因子统计与相关性分析

表1、2为各取样点的土壤质地和地形指标的描述性统计结果。结果表明，0～20 cm的SWC随土壤深度增加而增加，这是由于休耕期的表面土壤水分蒸发导致的；速效氮(AN)含量和电导率(EC1∶5)在0～40 cm大于40～60 cm，说明AN以及阴阳离子在上层土壤残留多。根据WARRICK等[24]对变异级别的分类，各因子均呈现中等变异性。SWC、AN和EC1∶5的变异系数均随土壤深度的增加而减小，这是由于降雨等外界条件对深层土壤影响较小，受自身因素影响较大，从而变异性变小。

表1 研究区土壤质地、土壤养分指标统计结果Tab.1 Descriptive statistics results of soil properties and nutrient concentration of study area

表2 研究区农田地形指标统计分析Tab.2 Descriptive statistics analysis of topography of study area

由表3可见，本研究区内，玉米产量与土壤粉粒含量、SWC和AN均存在极显著正相关关系，与土壤砂粒含量和EC1∶5均存在极显著负相关关系，与高程存在显著负相关关系。说明这些因素对玉米产量会产生较大影响，是玉米产量的主控因子，其结果与其他学者的研究相似[7,21]。各主控因子空间分布规律与产量分布相似，产量与土壤粉粒含量、SWC和AN空间分布趋势一致，而与砂粒含量和EC1∶5呈现相反的分布规律，这与相关性分析结果一致(图2)。

表3 研究区地形指标、土壤质地指标、土壤含水率、速效氮含量、电导率和产量的相关性系数Tab.3 Correlation analysis between topography, soil properties, SWC, AN, EC1∶5 and maize yield in study area

注:*和** 分别表示显著性(P<0.05)和极显著性(P<0.01)水平。

图2 研究区高程、土壤粉粒、砂粒含量、土壤速效氮含量、土壤含水率和电导率的空间分布图Fig.2 Spatial distribution maps of elevation, silt, sand, AN, SWC, and EC1∶5 in study area

2.2 主成分分析

由相关性分析结果可以得出，土壤粉粒含量、砂粒含量、SWC、AN、EC1∶5和高程为产量的主控因子，对其标准化后进行主成分分析。由分析结果可以看出，第3主成分(PC3)特征值接近1，且前3个主成分累积贡献率为88.62%(>85%)，故选取前3个主成分(PC1、PC2和PC3)作为主成分，可以反映各因子空间变异性88.62%的原始信息(表4)。土壤粉粒含量和SWC对PC1的荷载分别为0.91和0.67，砂粒含量的荷载为-0.88；EC1∶5对PC2的荷载为0.75，而对PC1的荷载为-0.59；高程对PC3的荷载为0.65(表5)。对比主成分得分和产量的空间分布图(图3)可以看出，PC1得分分布与产量较为一致，在产量高的区域其得分也越高；PC2得分分布与产量的空间分布相反；PC3得分空间分布与高程较为一致。这与主成分中荷载大的因子分布相关：在PC1中，粉粒含量和SWC的正荷载占据主导地位，而在PC2中EC1∶5荷载大于其他因子，这与相关性分析结果一致。结果说明，PC1主要反映土壤质地与养分的信息，PC2则主要反映EC1∶5(盐碱度)的信息，PC3则主要反映研究区高程的信息。

表4 主成分特征值和贡献率Tab.4 Characteristic values and contribution rate of principal components

表5 研究区主成分因子载荷Tab.5 Factor loadings of principle components in study area

图3 研究区玉米产量、管理区和主成分得分空间分布图Fig.3 Spatial distribution maps of maize yield, management zones and principal components values in study area

2.3 基于模糊聚类的管理区划分

将前3个主成分作为模糊c均值聚类分析的输入变量[23]。依据模糊效果指数(FPI)和归一化分类墒(NCE)同时最小即为最优的原则确定最优分区数[13,25]。由图4可以看出，当分区数为3时满足FPI和NCE均为最小值，故将本研究区域分为3个管理区，分别命名为M1、M2和M3。

图4 不同分区数时的模糊效果指数(FPI)和归一化分类墒(NCE)Fig.4 Calculated fuzziness performance index (FPI) and normalized classification entropy (NCE) in study area

由分区图(图3)可见，分区情况与实测产量的分布有很大的相似性。对3个管理区各因子与产量进行方差分析，结果显示不同管理区间除高程外其他因子的差异均达极显著水平(P<0.01)(表6)。对各分区内产量的主要影响因子进行统计分析(表6)发现，M1管理区内土壤粉粒含量、SWC和AN在3个分区中最高，砂粒含量和EC1∶5最低，是玉米生长较为适宜的土壤环境，其产量为7 244 kg/hm2；反之，M3管理区内SWC和AN含量低，砂粒含量和EC1∶5高于其他两个分区，土壤环境对玉米生长有不利影响，其产量为5 502 kg/hm2；M2管理区各因子均值处于M1和M3之间。这与相关性分析和主成分分析的结果一致。同时，大部分因子的变异系数较分区前下降0.01%～38.10%，说明在分区内各因子的分布趋于均质化(表6)，达到了农田管理分区的目的。

2.4 生育期内不同分区土壤含水率与作物生长指标对比分析

根据分区结果，对不同分区内的土壤含水率和玉米株高、LAI进行时间尺度上的对比分析，可以看出在各因子之间差异最大的位于M1和M3区域内，SWC变化以及植株形态指标差异显著；同时3个管理区内株高与LAI由大到小顺序均为M1、M2、M3(图5a)。在拔节期之前，玉米株高和叶面积指数在3个管理区中差异不大，但在拔节开始后产生明显差异。造成作物形态指标存在明显差异的原因是土

表6 研究区内不同分区主控因子均值、变异系数统计与方差分析Tab.6 Mean value, coefficient of variation and variance analysis of master factors for different management zones in study area

壤各组分存在变异性，且3个管理区内土壤含水率差异明显：在M1管理区内土壤持水能力强，土壤水分向深层移动速度缓慢，在下一次灌水前可为作物持续提供水分，植株可持续利用的水量多，利于作物吸收水分和养分；在M3管理区内土壤砂粒含量高使土壤透水性强，土壤水分向深层移动速度快，造成作物主要根系层没有充足的水分与养分供给，从而出现水分和养分胁迫影响作物正常生长；M2管理区的土壤水分的变化趋势在M1和M3之间(图5b)。

图5 M1、M2和M3管理区内玉米株高、叶面积指数和土壤含水率变化对比曲线Fig.5 Changing curves of leaf area index (LAI), crop height and SWC in M1, M2 and M3 zones during growth period

3 讨论

农田管理分区是实现精准农业的重要环节，是农田分区和分布式管理的基础。大面积农田内土壤属性、地形属性的空间变异性，导致了作物生长和产量也有较强的空间变异特征[25-27]。李翔等[28]验证了多种传统分区方法的实际分区效果，发现利用空间连续性聚类算法(SCKM)所得分区结果的聚集度和破碎度均为最优，有利于田间实际操作；土壤质地对作物有着很大影响，土壤质地的差异以及各层土壤组分的不同对土壤水分的分布产生较大影响，同时土壤养分与地形因素也是不可忽略的因素[29-30]；土壤质地的差异对土壤水分、养分运移与分布的影响也造成了不同分区内作物产量的差异[31-32]。本研究中，高程作为玉米产量的主控因子，其影响比其他主控因子小，且在前两个主成分中荷载明显小于其他因子，与其他学者的研究不尽相同[33-36]，这主要是因为本研究区域内地势平坦，高程的空间变异性弱，导致高程、坡度和坡向在本研究中对产量的影响不显著。

对于大面积农田，在灌水、施肥等管理措施均一致的情况下，由于土壤质地的空间变异性造成土壤持水性和土壤入渗性能的空间变异，对土壤水分、养分的运移产生影响最终导致产量存在极大差异[37-38]。因此对本研究中得到的3个管理分区，可以采用不同的方式进行管理。对于M3区域，其土壤质地砂粒含量高，土壤持水能力弱，因此需要在生育期内提高灌水频率，以保证作物不受水分胁迫；同时土壤速效氮含量较低、电导率较高，因此在种植前应进行滴灌压盐，并在生育期内增施氮肥。M1区域土壤质地主要为壤土或粉质壤土，持水性能良好，因此可减少生育期内灌溉次数以节省劳动成本，且土壤速效氮含量较高可适当减少施氮量。M2区域内各项指标均处于M1和M3之间，在管理过程中可减少生育期内灌溉次数同时适量施肥以满足作物需求量。工程方面，3个分区面积均为1.8～2.7 hm2，便于集中分区管理和田间机械化作业，因此在进行滴灌工程设计时，可以对不同分区进行变径设计，以实现灌溉的精确管理，从而为分区制定不同的水肥一体化制度提供理论依据和方法支撑。同时磷、钾等主要营养元素对作物生长有关键作用，因此在进一步的研究中，应对各管理区内土壤磷和钾等元素在生育期内的变化进行分析，以得到更为科学的管理建议。

4 结论

(1)研究区域内土壤粉粒含量和高程均呈现弱变异性，土壤粘粒、砂粒含量呈现中等变异性，SWC、AN和EC1∶5呈现中等变异性；玉米产量与土壤粉粒含量、SWC和AN存在显著正相关关系，与土壤砂粒含量、高程和EC1∶5存在显著负相关关系；通过主成分分析可知，土壤粉粒含量、SWC对产量空间变异有主要的正贡献，砂粒含量和EC1∶5对产量有主要负贡献。

(2)根据3个主控成分得分的空间分布规律，利用模糊c均值聚类的方法，将农田分为3个管理区，管理区之间各因子有很大差异，管理区内大部分因子呈弱变异性；不同管理区内作物生长和产量有明显差异，主要是由于土壤质地的空间变异性使SWC变异性增强，造成了生长与产量的较大差异。

(3)可以根据分区结果指导农业的精准管理，为规模化农田滴灌工程的设计提供科学指导；同时依据分区内不同的的土壤属性与养分，制定不同的灌水施肥制度，以提高大面积农田管理的针对性。