基于GBDT 的望奎县农田土壤有机碳主控因子研究

祝元丽,冯向阳,闫庆武,吴子豪*(.中国矿业大学公共管理学院,江苏 徐州 226；2.中国矿业大学矿区土地利用与生态安全研究中心,江苏 徐州 226；3.中国矿业大学中国资源型城市转型发展与乡村振兴研究中心,江苏 徐州 226)

东北黑土区是世界四大黑土区之一,具有优越的农业生产条件[1-3].然而由于高强度的人类利用和土壤侵蚀等原因,黑土地肥力逐年下降,表层土壤有机碳(SOC)含量降低等问题出现[4].因此,正确揭示黑土区农田土壤有机碳空间分异特征,并探究其分异的主控因子,阐明土壤有机碳与土壤理化性质、气象因子和农田管理措施变化等影响因素的关系,对于农田土壤有机碳的空间预测和碳库的精确估算具有重要意义.

近年来,很多学者利用多源数据对农田土壤有机碳的空间分异特征及其影响因素进行了深入地分析[5-9],研究指出,农业有机碳的变化与气温、降水、地形和土壤理化性质以及农艺管理措施,例如:轮作、秸秆还田、耕作制度以及化肥和动物粪肥的施用密切相关[10-12].土壤养分含量通过影响植被生长影响植被凋落物的输入量,从而影响SOC.Yu 等[13]、Li 等[14]和Yuan 等[15]的研究分别表明土壤氮、磷、钾元素含量与SOC 高度正相关.但是以上多数研究的数据来源于野外采样和国家土壤普查结果,受成本限制而难以扩展到大范围,其时效性也难以保障,这也导致了土壤氮、磷、钾等土壤养分数据在数字土壤制图研究中较少得到使用.近年来,随着开源土壤数据的逐渐丰富,例如世界土壤数据库(HWSD)、世界土壤信息数据(SoilGrids)以及中国高分辨率国家土壤信息格网基本属性数据集(2010～2018)等开源土壤数据网站,提供了我国土壤养分和土壤质地的空间分布信息.将以上开源多源土壤养分因子纳入解释变量指标体系,为更准确地探究SOC 空间分异及主控因子提供了新的契机.目前机器学习算法,包括决策树(DT)、人工神经网络(ANN)、支持向量机(SVM)、随机森林(RF)、Cubist 等模型,被广泛应用于土壤属性空间分异的影响因子探究中[16-19],例如,袁玉琦等[20]使用随机森林模型对福建闽东南地区SOC 的影响因子进行了相对重要性排序,进而确定该区域SOC 主控因子是水解性氮因素.这些机器学习算法具有宽松的建模假设,并且充分考虑了土壤属性与环境因子间复杂的非线性关系,其模型的拟合优度高于传统线性回归模型,可以更好确定土壤属性空间分异的主要影响因子.然而,现有研究大多止步于环境因子的相对重要性分析,鲜有研究对土壤属性与环境变量的非线性关系进行可视化,难以揭示环境变量对SOC 的影响机理.

望奎县位于典型黑土区,处于东北黄金玉米带,耕地占全县面积的75%.长期高强度耕作引发了严重的水土流失,侵蚀驱动的土壤物质迁移造成SOC的空间迁移与再分布,使SOC呈高度空间异质性.本研究以望奎县为典型范例,利用多源开源气候、土壤以及地形等数据,使用梯度提升决策树模型(GBDT)探究该区域SOC 与环境变量的非线性关系,确定SOC 主控因子,揭示SOC 空间分异机理.

1 材料与方法

1.1 研究区概况

望奎县位于黑龙江省绥化市的中部(46°32′07″N～47°08′24″N,126°10′23″E～126°59′00″E)(图1),辖区面积约为2314km2,地处于松嫩平原腹地.该县东部为丘陵漫岗区、中部为漫川漫岗区、西部为低洼平原区,地势东高西低,海拔在46～256m 之间.望奎县内南部有呼兰河流经、西部有通肯河环绕,水力资源丰富.该县属于中温带大陆性半湿润季风气候区,年均降水量约为650mm,年均气温为3.9℃.土壤类型以黑土为主,是国家商品粮食生产基地县和全国产量大县.耕地面积为1724km2,占全县面积的75%.随着望奎县黑土地的不断开发利用,带动了当地经济的快速发展,但也相对带来了土壤退化等一系列问题,对当地土地资源带来负面影响[21].

图1 基于GBDT 的土壤有机碳模型建立流程Fig.1 Soil organic carbon modeling process based on GBDT

1.2 梯度提升决策树模型

梯度提升决策树(GBDT)是基于CART 算法的增强集成学习模型[22].它利用弱分类器(CART 树)构建集成模型,通过反复迭代训练多个决策树模型,并且将其进行叠加,达到提高预测准确性的目标.在每次迭代过程中,在上一次迭代的基础上,计算损失函数(式1),求得伪残差(式2),并求得该迭代,进而构建一个新的决策树,然后通过梯度下降的方式,根据决策树的权重(式3),将所有生成的决策树进行加权融合(式4).GBDT 模型将决策树与集成思想进行了有效的结合,提高了模型的预测精度.GBDT 算法中,通常包括以下步骤:

GBDT 具有较高的准确性,可以在建模过程中识别并纠正误差.但是GBDT 对异常值比较敏感,在多次迭代中,GBDT 模型会尝试着去拟合异常值,导致过拟合的风险增加.为防止过拟合现象发生并提高预测性能,本研究导入数据进行超参数调优,以便获取参数最优解.本研究初始设置了5000 个弱分类器,6 个节点,采用十倍交叉验证导入数据进行拟合(图1).

运用一元线性回归模型与GBDT 模型构建了SOC 与解释变量间定量关系,并对比了一元线性回归模型与GBDT 模型用于探究SOC 与解释变量间定量关系的解释程度,建模与拟合过程均在 R studio 软件中完成,使用决定系数R2来评估模型的拟合效果.最后,通过比较GBDT 模型中的单因素与双因素的R2来探究因子间的交互效应,若R2(X1∩X2) = R2(X1) + R2(X2),则因素X1与因素X2独立,若R2(X1∩X2)

1.3 SOC 和环境变量获取

1.3.1 SOC 数据的获取及处理 SOC 数据来源于国家地球系统科学数据中心(http://www.geodata.cn),首先对望奎县2012年的1000多个采样点进行整理,剔除异常值和错误点,之后采用ArcGIS 软件将余下的采样点进行矢量化,并根据关键字段与SOC 数据进行关联,使用Kriging 法进行插值;按分层随机抽样方法,剔除掉望奎县内非耕地,林地与草地区域,在望奎县境内随机生成2000个样本点,每个样本点相距至少1km,在面积较小的土壤类型(如:雏形土)与土地利用类型(如:林地和草地)中生成更多的点,以提高样本的代表性,使样点覆盖整个研究区(图2).

图2 研究区位置和SOC 空间分布Fig.2 Location and SOC spatial distribution of the study area

1.3.2 环境变量数据的获取及处理 基于Scorpan 理论模型探究SOC 的影响因子,选取了包括土壤属性(S)、气候因子(C)、有机体(O)、地形(R)、和空间位置(N)5 类25 个影响因子(表1),空间分布状况见图3.

表1 环境变量分类体系及来源Table 1 Classification system and source of environmental variables

2 结果与分析

2.1 不同农用地利用类型的土壤环境因子分析

由图4(a)可以看出,水田的平均SOC 含量为20.28%,旱地的平均SOC 含量为20.41%,林地的平均SOC含量为20.37%,草地的平均SOC含量为20.35%,出现旱地>林地>草地>水田的情况;而土壤pH 值在各个农用地利用类型中均呈碱性,且差异较小,对SOC 含量影响较小.除碱解氮含量外,旱地的SOC、有效磷、速效钾、黏粒含量以及pH 值均高于水田、林地以及草地,并且旱地的碱解氮含量仅低于林地.通过对比可知,SOC 含量直接与碱解氮、有效磷、速效钾含量有关,旱地的SOC 含量以及碱解氮、有效磷、速效钾含量高是施肥作用的效果,而水田为了防止烧苗环境污染,施肥较少且慢速释放肥料,在水里水肥共生,肥料缓慢分解,难于深入到土壤中,导致肥效降低.对比不同农用地利用类型的SOC、碱解氮、有效磷、速效钾、pH 值以及黏粒含量的差异性,可以看出农用地利用类型对SOC 差异表现出有效性,而图5 中GBDT 模型自变量的相对重要性可以得出土地利用这一指标在整个环境变量指标体系相对重要性排序中相对较低,表明农用地利用类型对于土壤养分含量影响的主要方式为施肥这一途径.

图4 不同农用地利用类型的土壤环境因子统计Fig.4 Statistical map of soil environmental factors of different land uses

图5 GBDT 模型中环境变量的相对重要性Fig.5 Relative importance of environmental variables in GBDT

2.2 模型预测效能

表2 结果表明,GBDT 模型相较于一元线性回归模型在所有的自变量上均具有更好的预测效果,R2提升了13.4%～3671.94%,尤其是地形粗糙指数、到道路的距离、增强型水体指数、归一化植被指数、坡度、裸地指数等自变量,模型的拟合效果都有明显的提升.在一元线性回归模型和GBDT 模型中,各变量的个重要性排序基本一致.以上结果表明,SOC 与选取的各环境变量的关系为非线性关系,GBDT 算法可以更灵活地拟合复杂的非线性模式,因此预测效果比ULR 更好.因此后续研究中仅展示GBDT 模型的结果.

表2 GBDT 模型与一元线性回归模型拟合效果对比Table 2 Comparison of fitting effect between GBDT model and ULR model

2.3 土壤有机碳的主控因子分析

GBDT 模型总体的R2为0.958,表明该模型能够很好地解释目标变量的变异,与实际观测值较为一致.

如图5 所示,利用GBDT 模型探究SOC 与解释变量间定量关系的解释程度时,各个变量对于SOC的预测能力存在很大差异.其中,碱解氮、有效磷、气温、土壤微生物多样性和速效钾是所有变量中影响SOC含量最为显著的几个因素.如果考虑变量的度量单位和量纲,可以看出:碱解氮以33.85%的影响权重位居最高,这表明当提高土壤的碱解氮含量时,通常呈正向关联的情况下,SOC含量往往也会随之显著上升;其次,有效磷以24.88%的权重位列次席.碱解氮和有效磷对于作物生长的影响很大程度上决定了土壤中营养元素的供给能力,这些变量与作物生长发育密切相关,可以直接或间接影响作物的产出效益.

年最低气温和土壤微生物多样性以18.01%和11.72%的权重位列第3 和第4 位,这表明气温的高低和微生物丰富度是控制SOC 分解的重要调节因素.其中气温对于土壤中微生物代谢和分解活动具有重要作用,温度升高时,微生物活动率也会提高,这会影响SOC的分解速率.速效钾的相对重要性为6.8%,虽然在SOC 预测模型中的相对重要性较低,但其是植物生长运转所必需的关键元素之一,影响植物吸收养分的能力和土壤中其它营养元素的释放率.碱解氮、有效磷、气温、土壤微生物多样性和速效钾的相对重要性加和达到95.27%,其余变量均小于1.5%.

2.4 土壤有机碳与环境变量的非线性关系分析

图6 中,主要环境协变量与SOC 含量之间的非线性关系以相对重要性的降序呈现.环境变量并不是在其全部范围内都影响SOC 含量,在以上7 个环境变量中,均观察到对SOC 的阈值效应,6 个环境变量具有双阈值.7 个环境变量与SOC 之间的关系可以分类3 类:第1 类是碱解氮和有效磷,两者均与SOC 含量成正相关.碱解氮在170～180mg/kg 之间迅速增加,而在其余区间是稳定的,有效磷和碱解氮的发展趋势类似,在有效磷含量<40mg/kg 的范围内,SOC 含量保持平稳,在40～44g/kg 范围内,SOC 含量快速增长,随后回归平稳.

图6 环境变量对SOC 的阈值效应Fig.6 Threshold effects of environmental variables on SOC

第2 类是气温和土壤微生物多样性,SOC 含量先随温度升高,含量逐渐上升,-28.3℃后后急剧下降,最后趋于平稳,而土壤微生物多样性也是如此;温度升高,微生物代谢速率增加,能够参与有机质分解的微生物数量也会增多,但当超过一定的温度门槛时,微生物生长和活动阵亡极高,SOC 的分解反而会减缓;微生物指标小于108 时,固碳作用占主导,微生物会通过代谢过程吸收和存储有机碳并且积累在自己的组织中,从而促进了 SOC 的增加,大于108 之后代谢作用占主导,将土壤中的有机物转换成无机物质从而使SOC 含量降低.

最后一类是速效钾、土壤pH 值和降水量,与SOC 含量呈负相关.随着速效钾的含量上升,SOC 含量先保持平稳,当速效钾的含量超过195mg/kg 之后,SOC 含量快速下降,最后回归平稳状态.降水量与土壤有机碳之间存在着非线性关系,随着降水量的增加,SOC 含量先是处于平稳阶段,当达到降水量610mm 时,SOC 含量迅速下降进入第二个平稳期,直至降雨量达到640mm 左右,SOC 含量再次下降再次进入稳定阶段.

2.5 环境变量交互效应对土壤有机碳的影响

GBDT 模型中前7 个重要变量对SOC 的双因素交互效应见图7.任何2 个变量都不是完全独立的,均具有相互作用的关系.任何两个环境变量的相互作用效应均为负值,而仅有降水量与速效钾、降水量与微生物、土壤pH 值与速效钾的交互效应为正.其中碱解氮与微生物解释了89.9%的有机碳变异,碱解氮和有效磷解释了87.5%的有机碳变异行为,而碱解氮与其他变量也均解释78.3%以上的有机碳变异.综上,碱解氮对于SOC 空间变异的解释程度最强.

图7 GBDT 模型中环境变量的双因素交互效应Fig.7 Interaction effects between environmental variables in GBDT

3 讨论

3.1 环境变量对土壤有机碳的影响机理

土壤碱解氮、有效磷和速效钾含量等土壤养分含量是土壤有机碳空间分布的主要影响因子,本研究中土壤氮磷钾均对土壤有机碳产生较大程度的影响,其中碱解氮和有效磷对土壤有机碳的影响为正效应,表明土壤氮和磷元素在一定程度上对土壤有机碳含量的增加产生积极效应.一方面,氮添加会影响土壤pH 值、养分有效性和土壤微生物多样性,进而影响土壤有机碳的分解.氮元素的输入对土壤有机碳的分解速率也会产生影响,因为研究表明氮元素的富集会导致土壤酸化,通过改变土壤的理化性质而抑制微生物的代谢活动,进而抑制有机碳的分解[33-35].人为活性氮的增加能够显著提升生态系统中土壤有机碳的固定,从而增加土壤有机碳的含量[35].另一方面,随着施氮量的增加耕层土壤有机碳含量显著提高[13,37],耕层土壤有机碳的输入主要是种植的作物生物量决定的[38],而大量氮肥的输入能够显著提高地上作物生物量和根生物量,从而显著增加土壤有机碳的含量[39],这与本研究的结果一致.本研究发现土壤氮元素在所有的环境变量中对土壤有机碳的影响最大,这可能与土壤氮元素从影响土壤有机碳的输入和输出两个方面来对总量产生影响.

本研究中,磷元素对于土壤有机碳的影响仅次于氮元素,呈明显的正相关关系,说明磷元素对土壤有机碳的促进作用可能主要与植物碳输入的增加有关[40-41].有研究表明,磷元素的增加也会对土壤有机碳的分解产生显著的负效应[42],但有研究表明,磷元素含量对土壤有机碳含量的影响并不显著[43],这证明了关于磷元素对于有机碳分解的微生物调控记住目前尚不明确.Li 等[14]发现土壤氮和磷含量的高低虽然对土壤有机碳的分解以及固存没有直接的影响,但却能通过碳氮磷三者的比值增加或减少微生物的生物量,从而达到影响微生物活性和微生物碳利用效率的效果,从而影响土壤碳动态.

本研究采用的环境变量数据以及土壤有机碳数据均为多源开源数据,仍得出与很多研究相似的结果,侧面证明了多源开源获取数据进行土壤有机碳空间分布的主控因子研究具有有效性,可能受到栅格数据分辨率的限制,部分影响因子的R2不高,这表明之后结合更精细的遥感影像数据获取环境变量有望进一步提升模型的精度.并且受数据可获取性的限制,本研究未考虑风向风速、土壤母质、施用农药、杀虫剂和有机肥等因子的影响,对模型的精度造成了影响.本研究得出土壤氮元素对于土壤有机碳含量的影响非常显著,证明控制土壤C/N 对于耕地土壤有机碳的可持续利用具有重要意义.

3.2 环境变量的阈值效应分析

本研究发现,环境协变量与土壤有机碳之间存在着非线性关系,这与前人研究一致[44-45],模型评估结果表明GBDT 模型为何优于线性回归模型,即GBDT 模型的估计结果更为可靠,而线性回归模型假设数据间的关系是线性的,因此以上结果证实了在实际情况下,数据间的关系常常是非线性的,而GBDT 算法可以灵活地拟合复杂的非线性模式.Lamichhane 等[46]和Siewert[47]的研究也证实了随机森林等非线性模型优于传统线性回归模型.因此,后续研究推荐使用机器学习算法来探究土壤属性和环境变量间复杂的非线性关系.

本研究发现碱解氮和有效磷这两种土壤养分与土壤有机碳含量之间的非线性关系存在阈值效应,即引起土壤有机碳含量发生变化的范围区间是相对固定的,更多量的氮磷元素输入并不会引起SOC 的增加.Brown 等[48]研究也发现,过量施氮对土壤有机碳没有影响,只是补充了在农业生产中流失的氮含量,当施氮量超过需求量时,多余的氮素会以N2O 和NO3的形式排出,此研究佐证了阈值效应的存在,但并没有进一步得到阈值点.以上研究结果对于当地科学施肥和精确农业有重要意义,如何确定当前我国连续耕作的农田土壤的肥料投入,以最低肥料成本达到最大程度增加土壤中的养分和有机碳含量的目的,避免过度施肥造成的资源浪费和农田生态污染问题.化学氮肥的大量输入对土壤碳循环影响显著,引起了广泛关注[49],同时不合理的施肥是影响土壤温室气体排放的重要因素.综上所述,不合理的化肥施用容易对大气环境造成了严重的危害,同时也造成土壤性质的恶化以及肥力的退化.

4 结论

4.1 相较于一元线性回归模型,GBDT 模型在所有的自变量上均具有更好的预测效果,可以更灵活地拟合SOC 和环境变量之间复杂的非线性模式,R2提升了13.4%～3 671.94%不等.

4.2 7个主要协变量可以解释98.011%的SOC浓度变化,其中碱解氮和有效磷浓度起主要作用.

4.3 7 个主导环境变量对SOC 含量均存在阈值效应,其中6 个环境变量存在双阈值.阈值效应结果表明,在碱解氮(170～180mg/kg)和有效磷(40～44mg/kg)的位置有机碳含量迅速升高,在土壤微生物多样性(108)、速效钾(195～215mg/kg)、降水量(610mm 和640mm)的位置土壤有机碳含量下降,确定了协变量对土壤有机碳含量的双向和三向交互作用.

4.4 在确定环境变量的有效范围、SOC 空间分布的潜在过程和提高农田SOC 变化的解释力时应考虑阈值和相互作用效应.

致谢：感谢国家科技基础条件平台-国家地球系统科学数据中心(http://www.geodata.cn)提供数据支撑.感谢国家青藏高原科学数据中心(http://data.tpdc.ac.cn)提供数据支撑.