1988―2019年中国农田表土有机碳库密度变化及其主要影响因素

韩天富，都江雪，曲潇琳，马常宝，王慧颖，黄 晶，柳开楼，刘立生，谢建华，吴远帆，张会民*

（1 中国农业科学院农业资源与农业区划研究所/耕地保护国家工程研究中心，北京 100081；2 中国农业科学院农业资源与农业区划研究所祁阳农田生态系统国家野外科学观测研究站，湖南祁阳 426182；3 农业农村部耕地质量监测保护中心，北京 100125；4 江西省红壤及种质资源研究所，江西南昌 331717；5 湖南省土壤肥料工作站，湖南长沙 410005）

近200个国家在第26届气候变化大会上初步达成了《格拉斯哥气候公约》，要求各国政府基于最佳的科学知识制定有效的气候行动和决策，以实现21世纪中叶前后达到碳中和的目标[1]。土壤有机碳(SOC)是全球碳循环关键的碳库之一，对农业生产、土壤肥力和地表与大气温室气体通量具有重要的调控作用[2]，而土壤中碳的固定对全球气候变化和粮食安全至关重要[3]。单位面积内土壤保持的有机碳量称为土壤(有机)碳密度。深入剖析SOC密度时空分布特征及其变化的影响因素，对提升SOC密度、提高土壤质量、缓解气候变化以及保护生态环境具有重要意义，尤其是受人为因素干扰较大和含量相对较高的表层SOC密度变化是当前研究的热点[4–5]。

近几十年来，围绕农田土壤碳密度变化国内外展开了大量的研究，包括基于土壤普查数据的全国尺度估计[6–7]，单个省份或地理区[8–9]，不同农作区[10]或观测点位[11]的变化分析。从土壤背景SOC水平看，我国是SOC密度低且分布不均匀，普遍低于欧美等发达国家[6]。SOC密度的空间分布主要受气候、植被以及人类活动的影响，且在地带性分布基础上人类活动加剧了SOC密度的变化幅度[7]。于严严等[12]提取了发表文献中数据并与全国二次土壤普查数据作对比，判明2000年左右中国农田土壤SOC平均密度略有增加，具体是华北的潮土和褐土区增幅较大，而东北区呈降低趋势，西北区没有明显变化。Zhou等[9]的研究表明，碳密度地理变异的主要驱动力是土地利用类型，特别是华北和东北地区，而气候变化因素较为普遍，只是贡献程度不同。根据Yu等[13]研究，中国半干旱地区SOC变化分别与温度、降水量和pH等因素相关。此外，由于我国农田土壤类型多、分布广、土地利用方式不一、施肥等管理措施差异较大[14]，导致目前缺乏全国尺度的数据探究土壤有机碳密度的时空变化。

遥感和GIS技术的信息开发为利用大尺度数据估测农田土壤有机碳密度提供了技术支撑。但要精确估计投入碳密度变化，尚需区域代表性的实测数据支撑[15–16]和土壤容重与SOC含量间的关联关系[17]。缺乏长期连续监测下的土壤有机碳动态变化数据，也是目前我国农田有机碳密度估计存在巨大不确定性的因素。所幸，中国农业农村部在1988―2019的30年间，在我国典型农田土壤类型区陆续布设了一批肥料和肥力监测试验点，主要监测农民常规管理下作物产量、肥料投入、气候、表层(0―20 cm)土壤理化性质等指标，为进一步分析农田土壤有机碳库密度时空变化以及关键调控因素提供了条件。考虑到土地利用类型对土壤有机碳库的影响较大[10]，本研究以旱地、水田和水旱轮作3种利用方式为切入点，分别从时间和空间尺度，结合差异性分析和多元分析，深入探明中国农田表层土壤有机碳库的变化特征及调控因素，为我国耕地地力提升和粮食安全提供理论支撑，同时为农田土壤固碳减排和应对气候变化的策略提供参考。

1 材料与方法

1.1 数据来源

自1988年以来，农业农村部先后在我国粮食主产区有代表性田块上，陆续布设了一批肥力和肥效监测试验，起步探索阶段有46个点，始于1988年，之后逐步增加了监测点，包括1998年增加了121个，2004年增加了176个，2016年增加了955个，到2019年共计1298个持续监测点。监测点的地理位置、耕作制度、土壤类型、作物类型、分布面积、管理水平等均有较好的代表性，施肥、灌溉、除草等均按照当地常规方式进行[17]。监测内容主要为：作物产量、氮磷钾肥料每季用量，土壤理化指标(土壤有机质含量、有效磷含量、速效钾含量、pH和容重)。根据各监测点提供的肥料类型、用量和养分含量进行计算获得纯氮、磷和钾养分投入量。年均温和年降雨数据部分由监测点直接提供，部分通过国家气象科学数据中心(http://data.cma.cn/data/)下载。2000年之前，所有监测点均不进行秸秆还田，2000年之后，所有监测点陆续开展秸秆还田。

1.2 样品采集与测试

各监测点每年最后一季作物收获后，按“随机”、“等量”和“多点混合”的原则，采集5个0―20 cm土壤样品，手动除去肉眼可见的根茬及秸秆碎屑，混合均匀后风干研磨过1和0.15 mm筛，送至省级土壤测试中心，按照《土壤分析技术规范》[18]进行各指标的测定，其中土壤有机质含量采用重铬酸钾外加热法测定，土壤容重采用环刀法测定。

1.3 数据处理与分析

农田表层有机碳库即SOC密度(t/hm2)的计算公式为[19]：

式中，BD为土壤容重(g/cm3)；SOC含量为土壤有机碳含量(g/kg)；0.2为土壤深度(m)；10为转换系数。

为了探究农田表层土壤有机碳库在时间维度上的变化特征，我们基于监测时间的跨度，同时兼顾不同阶段部分点位的调整和增加，将监测数据每隔10年左右进行划分，共分为3个阶段，具体分为1988―1999 (鉴于前期数据较少，此阶段包含12年的数据)、2000―2009和2010―2019，这样划分在一定程度上可以降低因点位过少而导致的有机碳库的调查差异。为了探究表层土壤有机碳库在空间维度上的差异特征，我们结合监测点自然地理区域和粮食主产区将监测点划分为：东北地区(黑龙江、吉林、辽宁、内蒙古东北部共401个)、华北地区(北京、天津、河北、河南、山东和山西共276个)、西北地区(内蒙古西北部、甘肃、宁夏、陕西、青海和新疆共134个)、西南地区(重庆、四川、贵州、云南、广西和西藏共104个)、长江中游地区(湖北、湖南和江西共152个)、长江下游地区(安徽、江苏、上海和浙江共134个)、华南地区(广东、广西、福建和海南共97个)。

鉴于相同区域的气候相对一致，加之监测点在最初选定时考虑到了主要土壤类型和管理措施的代表性，进而弱化了由于新增点位和外源碳的输入量不同而造成较大的差异，能在一定程度上满足我们探究全国及各区域SOC密度随时间变化规律以及关键驱动因素。考虑到气候(年均降雨和年均温度)、氮磷钾养分年均投入量和土壤性质不仅能够直接或间接影响SOC的输入，而且可间接影响SOC的分解和转化[11]；加之提升回归树模型(BRT，一种基于机器学习的数据挖掘算法)能够分析不同类型的变量和变量之间的相互作用效应，适用于非线性关系，且能回答各自变量对同一因变量变化影响的相对重要性或贡献的百分比[20]，因此，利用BRT探究上述各因素对SOC密度变化的相对重要性[21]。

监测数据统一用Excel 2016整理并制表，运用Sigmaplot 14.0进行方程(单直线、双直线和二次曲线)拟合，借助SPSS 26.0进行差异显著性检验，利用R语言(4.0.5)中的gbm包进行BRT分析。

2 结果与分析

2.1 全国尺度上土壤有机碳密度统计分析

1988―2019年间全国农田表土层(0―20 cm，下同) SOC密度平均值和标准差分别为35.13和15.91 t/hm2，其中水田、旱地和水旱轮作下的平均值分别为 46.18、30.15 和 37.65 t/hm2(表1)，标准差分别为15.16、14.37和12.33 t/hm2，水田和水旱轮作较旱地分别平均高53.2%和24.9%。近3年(2017―2019)全国农田表土层SOC密度平均值和标准差分别为36.88和17.07 t/hm2，其中水田、旱地和水旱轮作下的平均值分别为45.15、36.06和38.01 t/hm2，标准差分别为19.08、16.47和15.19 t/hm2，水田和水旱轮作较旱地分别平均高25.2%和5.4%。与初始阶段(1988―1990)相比，近3年全国农田表土层SOC密度的增幅为1.7%，其中水田、旱地和水旱轮作下的增幅分别为10.3%、–3.0%和0.2%。

表1 土壤有机碳密度(t/hm2 )描述性统计分析Table 1 Descriptive statistics of soil organic carbon (SOC) density

2.2 全国尺度上土壤有机碳密度随时间变化特征

1988―2019年间全国农田表土层SOC密度在3个时间段呈先降低后显著增加趋势(图1)，旱地和水旱轮作下3个时间段变化趋势同全国尺度上一致，而在水田下无显著差异。方程拟合结果表明，整体上SOC密度在2003年之前随时间延长呈显著下降趋势，之后呈显著上升趋势，且上升速率低于下降速率(表2)；就水田而言，SOC密度在2000年之前随时间延长呈显著下降趋势，之后呈显著上升趋势，且上升速率低于下降速率；就旱地而言，SOC密度在1998年下降到最低点，之后呈增加趋势；就水旱轮作而言，SOC密度在2004年之前随时间延长呈显著下降趋势，之后呈显著上升趋势，且上升速率低于下降速率。

表2 年均土壤有机碳密度(y)与试验时间(x，年)的方程拟合Table 2 The fitted equations of yearly average SOC density (y) and duration (x, year)

图1 全国尺度下土壤有机碳密度随时间变化特征Fig.1 Characteristic variation in SOC density with monitoring time at the national level

2.3 各区域土壤有机碳密度随时间的变化特征

各区域SOC密度在3个时间段间的差异特征存在一定的异质性(图2)。就水田而言，东北地区Ⅲ阶段(2010―2019) SOC密度平均值显著高于Ⅰ阶段(1988―1999)，华南地区呈相反趋势，西南和长江中游地区Ⅱ阶段(2000―2009) SOC密度平均值低于Ⅰ和Ⅲ阶段，长江下游地区各时间段无显著差异。方程拟合结果表明东北和华南SOC密度分别随时间呈增加和降低趋势(表3)，西南和长江中游地区呈先降低后增加趋势。就旱地而言，除了西南和长江中游地区3个时间段无显著差异以外，其他各区域Ⅲ阶段SOC密度均显著高于Ⅰ阶段，且华北和华南地区Ⅱ阶段SOC密度均显著高于Ⅰ阶段和显著低于Ⅲ阶段。方程拟合结果表明东北、长江中游和长江下游地区SOC密度随时间呈先降低后增加趋势，华北、华南和西北地区呈增加趋势，而西南地区呈降低趋势。就水旱轮作而言，西南地区Ⅱ和Ⅲ阶段SOC密度均显著高于Ⅰ阶段，长江中游呈相反趋势，长江下游Ⅰ和Ⅱ阶段显著低于Ⅲ阶段。方程拟合结果表明长江中游和下游地区SOC密度随时间呈先降低后增加趋势。近3年的结果表明，长江中游和东北地区水田SOC密度较高；除了东北旱地SOC密度较高以外，其他区域旱地SOC密度均较低；西南地区水旱轮作SOC密度高于其他地区。

表3 各区域年均土壤有机碳密度(y)与监测年限(x)的方程拟合Table 3 The fitted equations of yearly average SOC density (y) and duration (x) in the main regions

图2 各区域不同时间段土壤有机碳密度随时间的变化特征Fig.2 Characteristic variation in SOC density across the monitoring time in each region

2.4 影响土壤有机碳密度变化的因素分析

由图3可知，SOC与全氮之间存在极显著的线性相关关系(P<0.001)，直线斜率表明旱地C/N (9.543)高于水田(8.706)和水旱轮作(8.656)。综合考虑碳氮较为强烈的耦合关系，以及在进行SOC密度计算时已经包含土壤容重因子，在进一步分析各因素对SOC密度的重要性时未包含全氮和容重的影响(图4和表4)。

图3 土壤有机碳和全氮含量之间的关系Fig.3 The relationship between SOC and total nitrogen contents

提升回归树分析结果(图4)表明，pH是水田SOC密度下降阶段差异的最重要的解释变量，且随pH增加SOC密度呈先增加后降低的趋势(P<0.01)(仅分析重要性排名前5的相应指标与有机碳密度之间的相关性)，其次是年均降雨和速效钾；速效钾是SOC密度上升阶段差异最重要的解释变量，且二者呈显著正相关关系(P<0.01)，其次是年均温度和全钾。年均温度是解释旱地SOC密度下降阶段差异的最重要的变量，且二者呈显著负相关关系(P<0.01)，其次是全磷和有效磷；年均降雨是SOC密度上升阶段差异的最重要的解释变量，且二者呈显著正相关关系 (P<0.05)，其次是有效磷和pH。磷肥是解释水旱轮作SOC密度下降阶段差异的最重要变量，且随磷肥增加SOC密度呈先增加后降低的趋势(P<0.01)，其次是氮肥和年均降雨；有效磷是SOC密度上升阶段差异最重要的解释变量，且二者呈显著正相关关系(P<0.05)，其次是年均降雨和年均温度。

图4 各因素对土壤有机碳密度的重要性分析Fig.4 Analysis of the relative importance of each index to SOC density

就不同区域而言，年均温度是东北、西南和华南地区水田SOC密度上升阶段差异的最重要的解释变量(表4)；钾肥是西南和长江中游SOC密度下降阶段最重要的解释变量；速效钾是长江中游SOC密度上升阶段差异最重要的解释变量。有效磷和年均降雨分别是东北旱地SOC密度下降和上升阶段差异最重要的解释变量；年均温度是西南地区SOC密度差异的最重要解释变量；年均降雨为华北和华南，氮肥是长江中游和下游SOC密度上升阶段差异最重要的解释变量，西北地区在前期和后期上升阶段最重要的解释变量分别是有效磷和年均温度。年均降雨和年均温度分别是西南和长江中游地区水旱轮作SOC密度转折点分别为2012和2001年之前差异最重要的解释变量，转折点之后均为有效磷；长江下游地区在转折点(2013年)之前和之后SOC密度差异最重要的解释变量分别是年均降雨和氮肥。

表4 各因素对不同区域土壤有机碳密度的重要性(%)分析Table 4 The relative importance (%) of each index to SOC density in the main regions

3 讨论

3.1 中国农田表层土壤有机碳密度现状

本研究通过对农业农村部的监测数据深入分析，全面了解了中国主要农田有机碳密度的时空变化特征。研究结果显示全国农田有机碳密度在近30年和近3年的平均值分别为35.13和36.88 t/hm2，而Tang等[22]通过采样调查分析了2011―2015年全国58个典型县域农田有机碳密度分布特征，发现表层0―20 cm有机碳密度的平均值(32.76 t/hm2)低于本研究结果，He等[23]通过分析1991―2012年IPNI(International Plant Nutrition Institute)项目的数据，得出的结果也低于本研究结果，主要原因可能是本研究在计算有机碳密度时忽略了土壤中的砾石占比；其次是本研究所选监测点位管理措施均相对较完善，土壤肥力水平均相对较高；此外，IPNI项目主要是根据经纬度网格化采样获得的数据，土壤肥力差异较大，最终致使本研究结果相对偏高；此外，部分点位的容重是根据经验公式计算获取的，在一定程度上影响有机碳密度的统计准确性。

随监测年限的延长，SOC密度整体上呈先降低后增加的趋势，这与国内部分研究结果存在微弱的差异。如于严严等[12]提取了发表文献中数据并与全国二次土壤普查数据作对比，判明2000年左右中国农田土壤SOC平均密度略有增加；邓祥征等[24]基于遥感反演的1988年与2000年土壤数据，发现年农田土壤有机碳库在1988―2000年呈增加趋势，而在2000―2012年呈降低趋势。主要原因是前人研究农田SOC密度时间变化规律主要以第二次土壤普查数据(1980年左右)为参照，而本研究的监测数据始于二次土壤普查10年后，期间随着作物产量的持续提升[25]，根茬和秸秆还田的有机物质也相应增加，在一定程度上提高了本研究的基础参照标准。此外，在1988―1999年，化肥用量显著增加，而有机物料投入的比例总体呈减小趋势[26]，进而导致土壤环境不断恶化，尤其是造成土壤酸化[27]，最终影响不同碳库之间的转换和土壤固碳能力[28]，图4结果也显示土壤pH在水田SOC密度下降过程中有重要的调控作用。在2000―2009年，随着高产品种的持续推广[29]、测土配方施肥技术[30]和秸秆还田进一步推进[31]，保证我国粮食产量整体呈逐年增加趋势[29]，在一定程度上增加了凋落物、根茬及根系分泌物等，进而为增加我国农田表层SOC密度提供了重要保障，这也是我国农田表层SOC密度在2000―2009年中后期整体上增加的主要因素。图4结果也显示土壤钾素在水田SOC密度上升过程中有重要的调控作用，这也进一步证明持续推进秸秆还田的重要性，我国水稻土整体钾含量相对较低且呈亏缺状态，随着氮、磷肥的大量施用和农民忽视钾肥的施用[32]，致使秸秆还田在缓解水田钾素缺乏[33]和保证水稻稳产和高产方面的作用尤为重要。

3.2 驱动各区域表层土壤有机碳密度变化的因素分析

SOC密度变化主要受其输入和输出的调控，输入主要来源于地上部凋落物、根茬及根系分泌物，而输出包括根系和微生物消耗及自身矿化作用[34]。这些过程均会受到作物耕作制度、施肥、灌溉、气候(温度和降雨)及其它土壤环境等因素影响[35]。大尺度SOC密度的时空分布主要受气候因素控制；而在相对小的尺度上，土地利用及管理措施等人为因素的作用强度往往超过气候因素[14]。如SOC密度在东北地区相对较高、西北地区最低(表3)，这与对应区域气候特点密切相关[36]。土地利用方式能够通过影响生态系统的分布与结构，进而影响碳循环和不同碳库之间的转换[37]。本研究结果显示，水田和水旱轮作下SOC密度比旱地平均高53.2%和24.9%，主要原因如下：与旱地相比，外源碳(如秸秆、残茬、凋落物等)投入到水田时能够更容易转化成SOC[38]；水田特殊的厌氧环境能够抑制微生物活性，进而阻碍土壤固有的有机质分解和促进SOC含量的稳定和提高[39]；此外，有机质中各组分的稳定性和物理化学保护机制的差异，导致其分解、矿化难易程度对环境变化的响应也存在较大不同[40]。

东北水田SOC密度呈逐年增加趋势，这与该地区水稻土类型和水稻生产特点密切相关，研究表明，在1980―2000年，黑龙江小麦种植区逐渐被水稻取代而水稻成为该区主要作物类型[41]，且种植水稻对土壤的碳汇效果高于小麦[42]；随着机械化、秸秆还田的推广和产量的持续增加[26,43]，在一定程度上促进了碳的开源和节流。此外，随着气候变暖，更加有利于水稻的生产和秸秆、残茬、凋落物等外源碳的输入[44–45]，有利于SOC密度的持续增加。而东北旱地在2008年之前SOC密度呈降低趋势，这与前人研究结果[46]一致，这可能与该区土壤过度垦殖而造成养分耗竭、水土流失加剧和外源养分的投入不足密切相关[47]，加之东北地区土壤初始SOC含量较高和气候快速变暖，进而促进碳的分解[48]和SOC密度逐年降低。而随着国家的不断重视和持续投入，如加大科研投入力度、奖励农民秸秆还田和兴修水利等，在一定程度抑制了东北地区SOC密度的下降[49]。

华北地区和华南地区旱地SOC密度分别呈逐年增加和降低趋势，且受年均降雨量影响较大(图4和表4)。尽管华北地区灌溉措施相对较完善，但是随着地下水位的不断下降致使作物生产对降雨的需求不断增加[50]，加之华北地区主要的种植制度是小麦玉米轮作，高度集约化的生产增加了对土壤的扰动和削弱对碳的保护，严重阻碍土壤有机碳的提升；而华南地区的旱地主要以红壤为主，肥力贫瘠，降雨充沛，过多的降雨将加速养分的淋失(尤其是钾素)、造成土壤酸化等，不利于作物的生长和有机碳的固存，且相对较高的温度在一定程度上促进了土壤中有机碳的分解[51]。研究表明，温度和降雨还能通过调控土壤pH进而间接影响(相关系数−0.35～−0.27)土壤有机碳的储存[52]。另外，华南地区随着双季或三季稻种植面积的不断缩减[53]，在一定程度上导致随作物残茬和秸秆输入[54]碳的量在减少。此外，华南地区的潜育化稻田分布较多，随着国家高标准农田建设的不断推进，使南方丘陵区一些长期淹水的冷浸田得到有效改善[55]，在一定程度上削弱水田的碳汇功能，进而导致SOC密度呈降低趋势。

磷是西南地区水田和水旱轮作SOC密度变化最重要的解释变量(图4和表4)，研究表明西南地区有效磷含量(<12 mg/kg)显著低于其他区域，难以保证作物正常的磷素需求[56]，应当重视磷肥的施用。西北地区SOC密度同样受土壤磷素影响较大，这与该区土壤中较高钙镁含量密切相关，大量的钙镁易固定土壤和肥料中磷[57]，阻碍磷素的有效性和作物产量的提升，而随着科学施肥(尤其是肥料类型)不断深入和秸秆还田持续推进[31]，促使该地区磷素的有效性呈逐渐增高趋势[56]。而长江中游地区由于风化程度高而导致土壤中铁铝氧化物含量较高[58]，易吸附土壤和肥料中的磷而降低其有效性，进而凸显磷肥投入对抑制SOC密度下降的重要作用。随着高产品种的大力推广，氮磷肥施用量逐渐增加和钾肥的施用易被忽视，加之该区域土壤钾素的本底值因主要含高岭石粘土矿物而较低，促使对外源钾素的需求不断增加[32]，导致缺钾成为该区域稻田产能提升的关键限制因素，进而造成土壤速效钾成为水田SOC密度上升阶段主要的驱动因素。然而土壤速效钾对水旱轮作下SOC密度变化的影响减弱，主要是因为水旱轮作下的干湿交替能促进土壤矿物钾的释放[59]，进而弱化外源钾的重要性。与水田不同，该地区旱地SOC密度下降阶段分别受土壤速效钾和pH的影响较大，该区域除了缺钾以外，旱地土壤pH显著低于其他区域[60]，严重阻碍作物正常生长，进而阻碍碳的输入和有机碳的提升。而氮肥对长江中下游地区旱地SOC密度下降阶段影响较大(表4)，因为长江中下游地区的旱地主要是一年两、三熟或者两年三熟，对养分的需求量大，足量氮肥能够在保证作物产量的同时增加作物残茬、根系和根际碳沉淀的输入[15,61]。此外，本研究缺少碳投入方面的详细数据，对研究SOC密度变化的影响因素有一定的影响，希望在以后的监测指标中尽量收集碳投入指标(如有机物料类型和投入量等)的数据，以便今后更准确的探究我国农田SOC密度变化特征。

4 结论

全国农田表层土壤有机碳平均密度为35.13 t/hm2，以水田>水旱轮作>旱地。总体而言，水田、旱地和水旱轮作土壤有机碳密度分别在2000、1998和2004年之前呈下降趋势，之后呈上升趋势。东北水田和西北、华北和华南旱地SOC密度呈逐年增加趋势，西南、长江中游水田和东北、长江中下游旱地呈先降低后增加趋势，华南水田SOC密度呈逐年降低趋势。

影响各区域农田表层SOC密度提升的主要因素：东北和西南地区的水田和旱地均为年均温；长江中、下游地区的水田分别为钾肥投入和土壤速效钾含量，旱地均为氮肥的投入；华南地区的水田和旱地分别为年均温和年均降雨；华北和西北地区的旱地分别为年均降雨和土壤有效磷含量；西南和长江中游水旱轮作区均为土壤有效磷含量，而长江下游为土壤速效钾含量。