中国农用地背景表层土壤有机碳对气候变化的土地用途依赖性敏感性

《CATENA》：Land use dependent sensitivity of background topsoil organic carbon to climate variability across China's agricultural lands

　　

• 背景耕地表层土壤有机碳（SOC）未呈现显著累积趋势。
• 背景表层土壤有机碳均值保持稳定，但存在结构性变化。
• 背景表层土壤有机碳的气候敏感性沿土地利用类型梯度而有所差异。

**1. 引言**

土壤有机碳（SOC）是陆地生态系统中最大的有机碳库（Lal, 2004），通过调节陆地与大气之间的碳、能量和水通量，在调节全球碳循环和减缓气候变化方面发挥着关键作用（Chaney et al., 2019; Crow et al., 2012）。此外，它是碳汇潜力的基石（Lal, 2004），也是维持土壤肥力和农业生产力可持续性的核心驱动力（Baldock, 2007; Chen et al., 2022）。

耕地作为与人类生产和生计紧密相关的生态系统（Zhang et al., 2024），既受气候变化影响，也受到持续的人为干预塑造（Ritchie et al., 2020; Pittelkow et al., 2015）。一方面，气候变化通过增强微生物活性影响土壤碳释放，形成潜在的正反馈，同时促进植物光合作用和植被碳输入（Wang et al., 2020; Liu et al., 2022b）。这些过程之间的平衡决定了生态系统碳汇或碳源的方向和幅度。另一方面，旨在提高产量、增强碳固存和减少温室气体排放的管理措施（如优化氮肥施用以减少N?O排放）也深刻影响着耕地SOC模式（Wang et al., 2020）。因此，阐明耕地SOC如何应对气候变化，不仅是中国实现碳达峰和碳中和目标的关键，也是保障国家粮食安全所必需的。

已有大量工作对耕地SOC进行了评估。例如，保护性耕作和秸秆还田已被证明能显著增强表层土壤SOC的累积，主要由活性碳组分驱动（Kang et al., 2026; Tao et al., 2026），而配对重采样研究揭示，尽管有机质输入增加导致表层SOC累积，但增温却导致深层SOC损失（Zhou et al., 2025b）。在国家尺度上，机器学习评估表明耕地SOC储量总体呈上升趋势（0.012 Pg C·yr?1）（Ai et al., 2025）。然而，这些评估将人为管理效应与气候变化对耕地SOC模式的影响混为一谈，而非将气候信号本身分离出来。集约化的人工管理投入掩盖了SOC对气候变异和气候驱动植被动态的内在响应；因此，目前的国家尺度评估未能将纯粹的气候驱动基线碳动态从人为管理信号中分离出来。

这一挑战在中国尤为突出，显著的温度纬度梯度与剧烈的降水经度变化（Li et al., 2016）塑造了SOC分布的高度空间异质性（Sun et al., 2026b; Li et al., 2024），使其成为研究全球变化背景下碳循环动态的天然实验室（Fu et al., 2015），但同时也使得在国家尺度上分离气候与管理效应格外困难。这一难度因三大农业用地类型的性质差异而进一步加剧（FAO, 2004）。水田受周期性淹水条件和厌氧分解影响（Zhu et al., 2023），旱地作物则在整个生长季保持好氧状态，这意味着相同的气候驱动力可能通过不同途径转化为碳累积或损失。相比之下，草地在直接农业干预（如集约化耕作和定期施肥）方面少于耕地，因此更直接地反映气候驱动力对SOC的影响，使其对气候波动天然更为敏感（Liu et al., 2021）。这些根本性差异意味着任何气候驱动的基线必须在不同土地利用类型之间分别解析，而非假设为统一一致。2000年至2024年为此类评估提供了一个关键窗口，涵盖了当代全球气候变暖最快的阶段（IPCC, 2023; Foster and Rahmstorf, 2026）。此外，自2000年以来遥感技术的进步提供了持续、高质量的环境协变量数据集（Zhang, 2017），使得对三类农业用地SOC空间显式模式的更可靠重构成为可能。

背景SOC含量定义为在最小管理投入条件下维持的内在碳含量，为解析这一挑战提供了独特基准。这类对照处理消除了施肥、秸秆还田和有机质改良等管理实践所引入的年际波动和方向性趋势，使得SOC动态主要由气候变异和气候相关的植被生产力所主导，使其成为量化气候敏感性的理想基线。尽管长期田间试验中的对照处理（CK）已被广泛用于量化特定管理实践的增量效应，但鲜有研究利用此类背景数据系统刻画不同土地利用类型下SOC的内在气候敏感性在国家尺度上的表现。草地和耕地CK样地均排除了刻意的碳增量投入，尽管其基线管理制度并非严格等同，但共同为评估气候敏感性提供了土地利用特异性基线。将水田、旱地作物和草地整合到统一的分析框架（农业用地；FAO, 2004）中，可以对不同土地利用类型下内在碳库如何响应共同的大气气候驱动力进行系统评估。

基于这一框架，我们认识到在最小管理条件下，SOC的年际变异源于两条耦合途径：气候对分解过程的直接调节以及由气候相关植被生产力所介导的间接碳输入变化。因此，本研究量化的SOC敏感性代表整合气候敏感性，涵盖直接气候效应和植被介导反馈。我们汇编了一个涵盖这三类农业用地的全国尺度背景表层（0–20 cm）SOC含量观测合成数据集。在本研究中，所有分析和预测模型严格针对该表层（0–20 cm）背景条件（为简洁起见，后文简称SOC）。以此统一基准，我们旨在回答两个基本科学问题：(1) 在气候变化背景下，中国主要农业用地的SOC如何动态演化？(2) SOC对气候变异的内在敏感性在水田、旱地作物和草地之间如何差异？

**2. 材料与方法**

为将气候驱动的SOC动态与混淆性人为干预分离，我们实施了基于对照处理数据筛选、土地利用稳定性掩模和静态模型参数化的三层解耦框架（文本S1）。

**2.1. 数据收集**

**2.1.1. 土壤有机碳含量数据**

使用以下关键词组合在Web of Science（https://www.webofscience.com/）、中国知网（CNKI，https://www.cnki.net/）、Google Scholar（https://scholar.google.com/）和万方数据库（https://www.wanfangdata.com.cn/）进行了系统文献检索：（"soil organic carbon" OR "SOC" OR "SOM" OR "soil organic matter" OR "carbon sequestration" OR "carbon stock" OR "carbon density" OR "soil nutrient" OR "soil quality"）AND（"cropland" OR "dryland" OR "paddy soil" OR "grassland" OR "agricultural land" OR "rice" OR "upland"）。检索范围覆盖至2025年6月的出版物，无语言限制。为分离背景SOC含量，我们应用了严格的选择标准：仅纳入明确定义了CK的田间试验；排除了实验室孵化试验和模型模拟结果。对于耕地，CK定义为未施肥、未秸秆还田和未施有机改良剂的样地。对于草地，CK定义为未放牧围栏或在环境放牧压力下且无灌溉、施肥或播种等额外管理干预的样地。此外，数据收集仅限于表层土壤（0–20 cm），并识别并移除了报告相同试验站点和采样时期的重复发表记录以避免重复计算。当可从表格或文本中直接获取数据时，数据直接从中提取。在数据仅以图形格式呈现的情况下，使用WebPlotDigitizer（https://automeris.io/WebPlotDigitizer/）提取数值，该工具是高精度图表数据恢复的可靠工具。对于仅报告土壤有机质（SOM）的出版物，使用以下公式将数值转换为SOC（Pribyl, 2010）：

SOC = SOM / 1.724

为确保数据集的可比性，所有SOC值均标准化为g C·kg?1。同时，记录每个研究站点的精确地理坐标（经纬度）和采样年份，以便利与环境协变量进行时空匹配及后续空间分析。共从1432篇同行评议出版物中提取了3963条原始数据记录，覆盖1815个不同站点。表S1提供了100篇代表性来源出版物的完整文献信息（作者、发表年份、观测期、碳类型、DOI及语言）。该数据集包含2326条旱地作物记录、1005条水田记录和632条草地记录。观测覆盖广阔的空间范围（77°16′–133°31′ E，18°34′–50°15′ N；图1），并具有一致的1983至2024年时间序列（数据S1）。各研究站点之间，SOC在不同土地利用类型间存在显著差异（P < 0.001），草地平均SOC最高（21.68 g C·kg?1），显著高于水田（15.63 g C·kg?1）和旱地作物（9.77 g C·kg?1）（图S1）。

*下载：高分辨率图片（562KB）*
*下载：全尺寸图片*

图1. 研究站点的空间分布与气候代表性。(a) 中国水田、旱地作物和草地点位的地理分布。土地利用类型分布图来源于资源环境科学与数据中心；(b) 采样站点投影到基于年平均温度和年平均降水的Whittaker生物群落分类方案上。

**2.1.2. 环境协变量**

为考虑SOC变化的复杂驱动因素，分析中整合了一组综合的环境协变量，包括气候、地形、土壤和植被因素以及土地利用类型。这些变量的详细描述和数据来源总结于表1。

表1. 环境协变量列表。

| 变量类型 | 变量名称（缩写） | 分辨率 | 数据属性 | 数据来源 |
|---------|----------------|-------|---------|---------|
| 气候因子 | 年平均温度（Tmean） | 1 km | 动态 | 国家青藏高原科学数据中心 |
| 气候因子 | 年最高温度（Tmax） | 1 km | 动态 | 国家青藏高原科学数据中心 |
| 气候因子 | 年最低温度（Tmin） | 1 km | 动态 | 国家青藏高原科学数据中心 |
| 气候因子 | 年降水量（Pre） | 1 km | 动态 | 国家青藏高原科学数据中心 |
| 气候因子 | 干旱指数（AI） | 1 km | 静态 | 全球干旱指数和潜在蒸散发气候数据库（v3） |
| 地形因子 | 海拔（Altitude） | 1 km | 静态 | 资源环境科学与数据中心 |
| 地形因子 | 坡度（Slope） | 1 km | 静态 | 1 |
| 地形因子 | 坡向（Aspect） | 1 km | 静态 | 1 |
| 土壤因子 | 容重（BD） | 1 km | 静态 | 国家青藏高原科学数据中心 |
| 土壤因子 | 黏粒含量（Clay content） | 1 km | 静态 | 国家青藏高原科学数据中心 |
| 土壤因子 | 粉粒含量（Silt content） | 1 km | 静态 | 国家青藏高原科学数据中心 |
| 土壤因子 | 砂粒含量（Sand content） | 1 km | 静态 | 国家青藏高原科学数据中心 |
| 土壤因子 | pH值（pH） | 1 km | 静态 | 国家青藏高原科学数据中心 |
| 植被因子 | 年NDVI最大值（NDVI max） | 1 km（2000–2024） | 动态 | MOD13A3 v061 |
| 植被因子 | 年NDVI最大值（NDVI max） | 0.05°（1982–1999） | 动态 | 国家地球观测数据中心 |
| 植被因子 | 年NDVI均值（NDVI mean） | 1 km（2000–2024） | 动态 | MOD13A3 v061 |
| 植被因子 | 年NDVI均值（NDVI mean） | 0.05°（1982–1999） | 动态 | 国家地球观测数据中心 |
| 土地利用类型 | 土地利用分类（Land type） | 1 km | 静态 | 资源环境科学与数据中心 |

注：1. 该变量基于DEM（数字高程模型）计算得出。

气候变量的选取旨在直接量化SOC动态的热量和水分约束，这是本研究的核心目标。这些变量包括年平均温度（Tmean）、年度温度极端值（年最高和最低温度，Tmax和Tmin）、年降水量（Pre）（Peng, 2019; Peng, 2020a, Peng, 2020b, Peng, 2020c）和干旱指数（AI）（Zomer et al., 2022）。高分辨率月度温度和降水记录被聚合至年尺度，以匹配每条SOC观测的特定采样年份。具体而言，Tmean计算为每年各月值的算术平均数，Tmax和Tmin确定为该年内月记录中的相应极值，Pre则导出为月度累计总和。年AI值同样从月度记录中平均得到以对齐采样时期。

地形特征（包括海拔、坡度和坡向）也纳入环境协变量集，以考虑景观尺度对SOC再分布的控制作用。海拔（Altitude）直接从DEM（数字高程模型）中提取，而坡度和坡向则基于DEM通过表面几何分析推导得出。

我们从高分辨率中国土壤信息栅格数据（2010–2018）（Liu and Zhang, 2021）中整合了土壤因子。与其他全球土壤产品相比，该数据集为中国的区域提供了增强的细节和精度（Liu et al., 2022a）。我们选择了表层土壤pH、容重（BD）和土壤质地（黏粒、粉粒和砂粒含量）作为土壤因子。这些土壤变量代表了土壤环境的基本理化状态，为SOC空间分布建模提供了必要的约束条件。

植被变量的纳入旨在表征植被生产力作为气候变异对SOC影响效应的潜在中介因子。植被生产力由年归一化差异植被指数（NDVI）的最大值（NDVI max）和均值（NDVI mean）表征。2000至2024年的数据通过Google Earth Engine平台从MODIS产品（MOD13A3 v061）中导出。我们实施了线性插值以填补时间缺口，并采用最近邻插值处理边界值以确保数据连续性。

对于2000年之前的记录，我们使用了一个与MODIS兼容的NDVI数据集（1982–2000年），该数据集提供了时间延伸的、与MODIS兼容的记录（Zhang et al., 2025）。该数据集采用优化的机器学习算法生成，以同步复杂的时空相关性并考虑多个传感器（AVHRR和MODIS）之间的异质性。所有月度记录最终被汇总为年度最大值（NDVI max）和均值（NDVI mean），分别代表峰值和平均植被活力。土地利用信息（土地类型；Xu et al., 2018）被用于将研究站点划分为三种主要管理类别：稻田、旱地和草地。为确保SOC估值的代表性，我们筛选了2000年至2020年间土地类型未发生变化的像元。仅保留这些稳定的像元以构建研究区域的最终基础底图（这些稳定区域在2020年共同代表了这些类别总土地面积的79%）。通过将分析限制在这些时空一致的网格单元中，我们有效消除了土地利用转变可能引入的偏差，并将研究聚焦于管理措施稳定的区域。

所有动态预测变量均与每个SOC观测值的具体采样年份和位置进行了时空匹配。对于没有确切年份或跨越多年窗口的记录，预测变量以其相应的时间均值表示。这一策略生成了一个稳健的、时空明确的、与3963个SOC观测值完全对齐的环境驱动因子数据库。为便于空间分析，所有协变量图层均被重投影并重采样为统一的研究区域1公里分辨率，形成后续机器学习建模的一致空间框架。

2.2. 模型模拟与分析

为准确模拟和分析SOC的时空分布，比较了四种候选算法：随机森林（RF）、分位数回归森林（QRF）及其残差克里金变体（RFK和QRFK），以建立SOC与环境协变量之间的复杂关系。

RF是一种稳健的数据驱动算法，通过对多棵树的预测结果取平均来生成预测值（Breiman, 2001a; Sun et al., 2026a）。QRF扩展了RF框架，描述了预测值的整个条件分布（Meinshausen and Ridgeway, 2006）；在本研究中，我们提取中位数（第50百分位数）作为最终输出，以增强预测稳定性。这些模型的一个显著优势是它们能够处理高维数据集，无需对数据分布做严格的参数假设，也无需排除存在多重共线性的变量（Breiman, 2001a; Breiman, 2001b; Meinshausen and Ridgeway, 2006）。因此，这些算法已被广泛应用于碳循环研究（Warner et al., 2019; Wang et al., 2023; Sun et al., 2023; Dharumarajan et al., 2024）。为进一步提高预测精度，将残差克里金与RF和QRF相结合，以同时考虑环境驱动因子和空间自相关（Wu et al., 2024; Kong et al., 2026）。

所有四种算法均在嵌套交叉验证（cv）框架下进行训练和评估。10折交叉验证用于评估整体模型性能，内部5折交叉验证用于根据均方根误差（RMSE）自动调整超参数。该方法确保参数优化与模型评估相互独立，有效减少过拟合。模型性能通过决定系数（R2）和均方根误差（RMSE）进行评估，计算公式如下：

R2 = ∑(i=1 to N)(X_model,i ? X?_obs)2 / ∑(i=1 to N)(X_obs,i ? X?_obs)2

RMSE = √(1/N × ∑(i=1 to N)(X_obs,i ? X_model,i)2)

其中，N为样本量，X_obs,i为第i个实际观测值，X_model,i为第i个预测值，X?_obs表示观测值的平均值。基于这些指标选择最优建模方案。

在最优超参数配置下，所有四种模型均表现出稳健的预测性能，在10折交叉验证下R2值范围为0.70至0.73，RMSE值为6.96至7.41 g C kg?1（图2）。RFK模型在通过内部5折交叉验证确定的最优超参数下表现最佳（R2 = 0.73，RMSE = 6.96 g C kg?1，ntree = 800，mtry = 4）（图2b），因此被选为最终建模方案，并使用完整数据集重新训练以进行空间外推。对于年度SOC预测（2000–2024年），最终RFK模型反复应用年度特定的气候和NDVI协变量。在该预测框架中，拟合的模型参数、静态环境变量和分类土地利用属性在这些稳定像元内随时间严格固定，有效将气候驱动的变化与土地利用转变或管理措施变化分离开来（补充文本S1）。这种固定关系固有地将原始SOC测量中存在的异质性吸收到模型参数中，而不是在任何特定年份中将其作为额外误差重新引入。因此，预测SOC的年际变化预计主要反映时空明确的输入变量的变化。

下载：下载高分辨率图像（472KB）下载：下载全尺寸图像

图2. 基于10折交叉验证，使用决定系数（R2）、均方根误差（RMSE）和P值比较观测与预测的背景表层土壤有机碳（SOC）含量（g C kg?1）。散点图展示了四种模型的性能：(a) 随机森林（RF），(b) 随机森林与残差克里金结合（RFK），(c) 分位数回归森林（QRF）和(d) 分位数回归森林与残差克里金结合（QRFK）。蓝色虚线代表1:1线，红色实线代表线性回归拟合线。颜色梯度表示观测点密度。（对于本图例中颜色参考的解释，读者请参阅本文的网络版本。）

最终最优模型的空间不确定性通过非参数自助法程序（10次迭代）进行评估（Sun et al., 2023; Yang et al., 2025; Sun et al., 2026b）。具体而言，完整数据集被随机重采样，每次模拟选取80%的观测值。每次重采样后，最优模型被重建并投影到空间环境协变量上。对于动态变量，采用多年平均值（2000–2024年）表示稳定的环境基线。为表征1公里分辨率下的预测稳定性，计算了10次模拟输出之间逐像元的变异系数（CV），作为不确定性的空间指标。由于自助法训练子集直接取自原始观测数据集，该不确定性估计固有地纳入了基础SOC测量中存在的异质性。CV计算公式如下：

CV = (σ / μ) × 100%

其中，σ表示10次自助法迭代的标准差，μ为这些模拟的算术平均值。较高的CV值通常指示预测可靠性较低的区域。为评估预测不确定性的空间模式是否反映了样本数据覆盖的限制，将观测样本在平均气温（Tmean）和降水（Pre）空间中的分布与整个空间域进行了比较（MAT和MAP，基于2000至2024年多年平均温度和降水从研究区域随机选取的10万个像元表示）。整体技术流程图如图3所示。

下载：下载高分辨率图像（345KB）下载：下载全尺寸图像

图3. 研究工作流程示意图。

为进一步评估预测SOC的年际波动是否反映真实的气候信号而非模型噪声，采用所选RFK模型（相同超参数配置，ntree = 800，mtry = 4）进行了时间分割验证。模型使用1983至2018年的观测值进行训练（数据集的78%），并使用2019至2024年的独立样本进行验证（数据集的22%）。计算每种土地利用类型的平均误差（ME）和均方根误差（RMSE），以评估当模型外推到训练数据中未包含的时间段时可能存在的系统偏差。

ME计算公式如下：

ME = 1/N × ∑(i=1 to N)(X_model,i ? X_obs,i)

其中，N为样本量，X_obs,i为第i个实际观测值，X_model,i为第i个预测值。

我们采用SHAP（SHapley Additive exPlanations）值增强模型可解释性，并捕获因变量与环境协变量之间的复杂非线性关系（Lundberg and Lee, 2017; Huang et al., 2025）。公式如下：

g(x') = φ? + ∑(i=1 to M) φ?

其中，g(x')是样本x'的预测SOC值，φ?表示基准值，代表整个数据集的平均模型输出，φ?表示Shapley值，量化第i个特征对预测的边际贡献。

值得注意的是，由于SHAP值需要单一预测模型，SHAP分析在最终RFK模型的RF组分上进行。RF是捕获环境协变量响应的机器学习组分，而残差克里金仅建模空间自相关的残差，因此RF组分的SHAP结果有效代表了RFK框架内协变量与SOC之间的关系。

2.3. 统计分析

Shapiro-Wilk检验证实，研究站点的SOC呈现显著的非正态分布（图S1a）；因此，采用非参数Kruskal-Wallis H检验评估不同土地利用类型之间的总体差异，随后使用Dunn事后检验进行两两比较（α = 0.05），并采用Bonferroni校正。所有报告的描述性统计（均值、中位数、四分位距（IQR）、标准差（SD）、标准误（SE）和偏度）及显著性字母均源自完整数据集。为可视化目的，密度图和箱线图的坐标轴限制在99百分位数，以最小化极端异常值的影响。为确保视觉清晰同时保持代表性，散点叠加图采用分层随机采样：对每个类别，随机选取最多10,000个像元；少于10,000个像元的类别保留所有可用像元。

采用Theil-Sen斜率估计量（Sen, 1968; Theil, 1992）计算2000至2024年SOC时间序列的变化率。作为一种基于所有可能成对斜率中位数的非参数方法，它对外部值具有很高的稳健性，在存在较大年际变异性时，比传统线性回归提供更可靠的趋势表征。随后，通过修正Mann-Kendall（MMK）检验确定观测趋势的显著性，该方法考虑了时间序列内的序列自相关，以避免显著性膨胀（Mann, 1945; Kendall and Gibbons, 1962; Hamed and Rao, 1998; Gu et al., 2025）。为进一步表征SOC的时间演变，Theil-Sen斜率估计量和MMK检验也应用于高SOC和低SOC组的面积比例时间序列。高SOC和低SOC含量区域使用预测多年（2000–2024年）平均SOC的中位数作为空间基准来定义：全国尺度为15.58 g C kg?1，稻田、旱地和草地分别为15.40、10.06和19.21 g C kg?1。

2.4. SOC对气候变化的敏感性

为量化SOC对气候变化的敏感性，我们使用了从2000–2024年时空SOC数据集中计算的变异系数（CV）（SOC CV）。在该框架中，较高的SOC CV值表示碳库对年际气候波动及其级联效应的更高敏感性。为便于空间比较和结构分析，基于数据分布的四分位数将计算的SOC CV值分类为四个不同的敏感性水平。具体而言，使用25、50和75百分位数作为切分点，将SOC CV值分为四个等级：高度稳定、稳定、中等波动和高度可变。利用Kruskal-Wallis H检验评估不同土地利用类型之间SOC敏感性的差异，随后使用Dunn事后检验进行两两比较（α = 0.05），并采用Bonferroni校正。

为归因SOC敏感性观测到的时间模式，我们量化了气候驱动因子的相对重要性及其植被介导的级联效应。对于全球和本地分析，使用了一组相同的预测变量，每个变量以其2000–2024年期间的CV表示，以捕获其年际波动的幅度。这些预测因子包括热因子（Tmean CV、Tmax CV和Tmin CV）、降水可利用性（Pre CV）和植被介导的气候效应（NDVI max CV、NDVI mean CV）。鉴于SOC站点没有接受管理投入，且土地利用类型在时间上保持不变，NDVI的年际波动主要由气候变化驱动。因此，植被CV被解释为气候变化通过植被通道产生的级联效应的综合替代指标——气候调控的碳输入路径。

采用RF回归结合自助重采样来量化变量相对重要性。对于全局分析，对每种土地利用类型进行30次独立自助迭代，每次随机有放回抽取5000个像元。使用平均平方误差的百分比增加（%IncMSE）量化变量重要性，原始值被归一化为相对重要性（总和为100%），以计算均值和标准差（SD）。

在局部动态分析中，沿年均降水量（MAP）梯度应用了双重嵌套移动窗口随机森林（RF）框架（窗口宽度=100 mm，步长=10 mm）。在每个窗口内，通过30次自助法迭代对3000个网格像素进行动态重抽样，以训练局部随机森林模型。所得的相对重要性和标准差（SD）值通过局部回归（LOESS）技术进行汇总和平滑处理，以追踪驱动因子沿MAP梯度的变化轨迹。值得注意的是，为确保变异系数（CV）的数学有效性，并防止统计异常——例如当摄氏度（°C）均值接近零或低于零时CV值趋于无穷大或失去意义——所有温度变量均被转换为开尔文（K）。所有统计分析和数据可视化均使用R软件（R Core Team 2025，版本4.5.2）在RStudio集成开发环境（Posit team 2026，版本2026.01.0）中完成。地理地图使用ArcMap 10.6（ESRI 2011，Environmental Systems Research Institute，美国加利福尼亚州Redlands市）生成。

3. 结果

3.1 模型性能与空间不确定性

所选的RFK模型取得了稳健的预测性能（图2）。模型中相对较高的均方根误差（RMSE）主要由异常大的标准差（SD，28.34 g C kg?1；图S1）以及草地地理采样密度不均匀（图1a）所致。按土地利用类型的验证进一步证实了模型的可靠性：水稻田的平均误差（ME）为?0.37 g C kg?1，旱地为0.19 g C kg?1，草地为?0.11 g C kg?1，表明在所有土地利用类型中均存在极小的系统性偏差（表S2）。所有区域合并后的平均误差为0.0035 g C kg?1（表S2），表明模型在全国尺度上并未系统性地高估或低估土壤有机碳（SOC）。

空间不确定性评估揭示了预测可靠性的显著地理异质性，变异系数（CV）值范围为0.43%至64.12%，平均值为6.99%；76.18%的像素值集中在0.43%至8.92%之间（图4a、b，S2）。最高的CV值集中在西北地区和青藏高原（图4a）。在各土地利用类型中，草地的平均CV最高（8.21%），其次是旱地（5.48%）和水稻田（3.99%）（图4b）。观测样本与完整空间域在气候空间中的对比表明，样本充分覆盖了主要气候梯度，但在高海拔寒冷干旱地区——年际均温（MAT）低于?5 °C且年均降水量（MAP）低于500 mm的区域——覆盖相对稀疏（图4c）。与旱地和水稻田相比，草地的样本-像素比最低（图4d）。

时间分割验证显示，所有土地利用类型的平均误差（ME）均较小（水稻田：1.73 g C kg?1；旱地：0.0086 g C kg?1；草地：?0.77 g C kg?1；表S3），表明当模型外推至训练数据中未包含的时间段时，不存在显著的系统性偏差。对于水稻田和旱地，时间分割的RMSE值（分别为5.61和5.63 g C kg?1；表S3）与或低于10折交叉验证的RMSE（6.96 g C kg?1，图2b），支持了模型对耕地预测的时间稳定性。草地较大的RMSE（26.08 g C kg?1；表S3）归因于其固有的高变异性（SD = 28.34 g C kg?1；表S2）和有限的晚期样本（n = 85；表S3）。总体而言，这些结果表明观测到的SOC年际波动包含真实的气候驱动成分，而非主要由模型噪声主导。

3.2 SOC含量的时空变化

通过RFK模型生成的年度空间预测图显示了2000至2024年间中国SOC高度一致且稳定的异质性分布格局（图5a）。在区域上，高SOC值显著且持续地分布于东北地区和青藏高原，这些地区分别以温带草原和高寒草甸生态系统为主，这归因于其低温条件和有限的微生物分解。相反，低SOC值集中在西北广大干旱半干旱地区和华北平原的部分区域，这些地区以集约化旱作农业为主，对应于9.77 g C kg?1的基线值（图S1a，表S2）。

在时间上，全国平均SOC在25年的研究期间内呈现显著下降趋势，从2000年的17.71 g C kg?1下降至2024年的17.04 g C kg?1，累计减少3.78%（P < 0.05），多年平均值为17.53 g C kg?1（图5b）。除绝对浓度的下降外，空间频率分布的趋势分析揭示了中国农业用地SOC库的显著结构性变化（图5c）。具体而言，高SOC区域（> 15.58 g C kg?1）的比例以0.04% yr?1的速率持续收缩（P < 0.001）。

按土地利用类别细分，三种土地利用类型在数量和时间结构方面表现出独特的变化特征（图6，表S2）。水稻田的多年平均SOC为15.61 g C kg?1，表现出同步退化模式，平均SOC（?0.0041 g C kg?1 yr?1，P < 0.001）和高SOC（> 15.40 g C kg?1）面积比例均以0.039% yr?1的速率显著下降（P < 0.05）。相比之下，旱地的多年平均SOC最低（11.05 g C kg?1），其平均SOC相对稳定，但表现出独特的结构优化，高SOC区域（> 10.06 g C kg?1）的比例以0.022% yr?1的速率扩大（P < 0.001）。值得注意的是，草地尽管保持了最高的多年平均值（20.96 g C kg?1）和最高的高/低SOC分组阈值（19.21 g C kg?1），但在三种土地利用类型中经历了最快的结构收缩，其高SOC面积以0.047% yr?1的速率缩小（P < 0.05）。高SOC面积收缩速率遵循土地利用类型的SOC等级顺序（图6，表S2）。这种差异凸显了不同土地利用类型间SOC动态的异质性，这很可能由基线碳储量对气候变率响应的差异所驱动。

进一步的时空区域分层揭示，SOC动态对气候的响应在不同生态系统类型间存在显著差异（图7）。在25年的研究期间内，旱地表现出全国范围内最高的SOC显著积累比例（35.63%，图7b），最显著的相对增加集中在北温带（74.55%）和南亚热带（70.58%）（图7f）。相反，草地经历了最高的相对空间收缩，其总面积的38.64%显示出显著的SOC下降（图7c），主要集中在中热带（99.28%）、南温带（51.46%）和暖亚热带（47.51%）的严重退化（图7g）。需要强调的是，这些子区域百分比反映的是各特定气候区内SOC趋势的内部比例，而非绝对地表面积。因此，在较小气候区观察到的高比例变化，本质上与该土地利用类型在这些特定区域的空间范围和基线足迹相关。

3.3 SOC对气候变率响应的分化

SOC敏感性空间格局通过SOC变异系数（CV）量化，在中国表现出显著的地域梯度（图8a）。空间异质性与400 mm降水量等值线高度一致，后者构成一条关键的生态界线。该线以东地区通常表现出较低的SOC敏感性，"高度稳定（0-2.84%）"和"稳定（2.85%–5.67%）"类别主导了景观格局，表明在气候变率下碳库相对稳定。相反，该线以西地区表现出明显升高的敏感性，"高度可变（> 10.40%）"斑块广泛分布于西北地区和高原地带（图8a）。

差异分析证实，土地利用类型显著调节了SOC敏感性（P < 0.001）。根据各土地利用类型的平均SOC CV，形成了明确的脆弱性等级：草地（9.02%）> 全国平均（7.85%）> 旱地（7.02%）> 水稻田（3.15%）（图8b）。草地小提琴图的分布范围更宽且峰值更高，进一步证明其气候响应范围更广，对气候变率的易感性高于耕地（图8b）。

在全国尺度上，SOC库主要由"稳定"（34.42%）和"中度波动"（28.22%）类别主导，两者合计占总面积的60%以上（图8c）。然而，分层分析揭示了三种土地利用类型在结构组成上的鲜明对比。水稻田表现出最高的气候缓冲能力，超过一半的面积被归类为"高度稳定"（53.27%），另有40.03%为"稳定"，而"高度可变"的比例微乎其微（0.69%）。

旱地处于中间位置，其中"稳定"类别占比最大（45.18%），"高度稳定"占20.42%、"中度波动"占21.46%。相比之下，草地表现出显著的结构脆弱性；"高度变化"（27.86%）与"中度波动"（35.12%）区域的合计比例超过60%，而"高度稳定"群体则被严重压缩至仅8.70%（图8c）。这一结构性差异进一步印证了以下结论：草地的土壤有机碳（SOC）对气候驱动退化异常敏感，而管理型水田生态系统在气候变异性下则表现出显著的稳定性。4. 讨论4.1. 预测性能与环境控制因素基于综合数据集，RFK模型展现出稳健的预测性能（R2 = 0.73，RMSE = 6.96 g C kg?1，图2b）。该预测精度与中国SOC制图研究中所报告的精度相当或更优（R2 = 0.55–0.72，Zhang等，2023；Xue等，2025；Zhou等，2026）。选取这三项研究作为基准，是因为它们代表了最近期的SOC制图工作，其模型架构与协变量结构均与本研究决定。本研究采用CK处理（对照处理）的SOC数据，以消除复杂管理措施对SOC的混杂效应，并整合残差克里金插值与随机森林方法以捕获空间自相关性，二者共同促进了模型性能相较于以往研究的提升。尽管模型对草地中的极端变异性进行了部分平滑处理（标准差：观测值28.34 g C kg?1 vs. 预测值11.47 g C kg?1；表S2），但这一效应并未损害均值预测及其趋势的可靠性。这从根本上取决于随机森林的自助法原理——自助法通过对大量噪声较大但近似无偏的树取平均来降低预测方差，而平均后树的期望仍与任何单棵树的期望相同（Hastie等，2009）。此外，基于Theil-Sen估计器的时间趋势检测本质上对这种方差压缩具有鲁棒性，因为该非参数方法评估的是单调方向性而非绝对量级（Sen，1968；Hamed和Rao，1998），从而有效过滤随机噪声以揭示一致的时间信号。基于最终RFK模型中随机森林（RF）组件的SHAP分析，用于解释SOC的空间分异，将年均温（Tmean）确定为SOC空间分布的主要驱动因素，这与将热量条件与土壤有机质分解-积累平衡相联系的理论框架相一致（Davidson和Janssens，2006；图9a）。除热力控制因素外，海拔（Altitude）作为次要驱动因素表现出显著的非线性阈值效应：在海拔约2000米以下，海拔对SOC的影响微乎其微；而在此阈值以上，其正向贡献显著增强并随海拔升高而加速（图9a、c）。高海拔地区观测到的极高SOC值与青藏高原高山草甸和冻土生态系统中已被充分记录的碳积累特征相一致，低温抑制了微生物分解过程（Zhong，2025），这进一步佐证了模型所识别的海拔-SOC耦合关系的生态学有效性。土壤理化因素（包括pH值、容重和黏粒含量）通过矿物-有机质络合与团聚体保护机制进一步调节了SOC的分布（Blanco-Canqui和Lal，2004）。SOC与其余环境变量之间的全面关系详见图S3。下载：下载高分辨率图片（444KB）下载：下载全尺寸图片图9. 基于SHAP值的随机森林（RF）模型可解释性分析。SHAP值代表某一特征值相对于所有样本平均预测值对模型预测的贡献度。（a）SHAP摘要图，展示各特征效应对背景表土SOC含量的全局重要性及分布。颜色梯度代表特征值；(b-f) SHAP依赖性图，显示SOC对各关键驱动因素的边际效应和非线性响应。不同颜色的点和实线表示不同土地利用类型的观测值及其由广义可加模型（GAM）拟合的相应平滑趋势，阴影区域代表95%置信区间。黑色虚线水平线表示SHAP值为零，高于此基线表示对SOC有正向贡献，低于此基线表示负向贡献（Tmean：年均温；Tmax：年最高温）。预测不确定性的空间分布与SOC对气候变异的敏感性分布具有显著相似性（图4a、8a），表明SOC本征敏感性较高的区域同时也是模型预测更具挑战性的区域。然而，驱动这两种空间模式的潜在机制存在根本性差异：基于多年平均环境协变量进行bootstrap重采样所得的较高预测不确定性，仅反映了观测覆盖不足所引发的不稳定性（图4c、d；Stell等，2021），与年际气候变异性无关。因此，未来研究应优先在这些观测不足的区域开展针对性采样以提高预测可靠性，并探索替代建模方法以进一步提升预测精度。相比之下，这些区域的高SOC敏感性既整合了稀疏采样导致的模型不稳定性，也包含了气候强迫驱动的真正年际SOC波动。时间序列分割验证支持了以下观点：年际SOC波动中包含真正的气候驱动成分；验证结果在所有土地利用类型中均显示出较小的平均误差（范围从?0.77至1.73 g C kg?1；表S3），表明在超出训练期进行外推时不存在显著的系统性偏差。然而，该验证并未量化气候强迫与残余模型不确定性对SOC敏感性的相对贡献，这仍是未来工作待解决的重要问题。4.2. 与既往研究的比较本研究汇编了1983年至2024年中国主要农业用地背景表土SOC含量的综合时空数据库。与以往整合了不同管理处理观测的全国尺度SOC评估研究相比（Zhou等，2025a；Ai等，2025；Wang等，2026），本研究提供了一种量化SOC对气候变化响应变异的框架，弥补了现有研究中人为信号与气候信号天然混杂的关键局限。尽管模拟的空间格局与既往SOC文献中报告的宏观地理分布总体吻合（Zhang等，2023；Song等，2025；Xue等，2025；Dong等，2026），但本研究识别出的时间轨迹却呈现出根本性的差异。具体而言，过去数十年间大多数常规评估均强调全国耕地SOC积累了广泛的上升趋势（Cheng等，2009；Han等，2022；Ding等，2023；Ai等，2025）。相反，我们的结果揭示耕地（水田和旱地）的SOC变化趋势相对稳定（图S4）。这一关键差异直接源于我们的目标数据库：通过严格采用长期对照处理，我们的框架成功地将背景SOC含量动态从混杂效应中剥离出来。在此，背景SOC含量的概念被明确定义为在无额外管理投入的情况下形成并维持的SOC基线水平；它代表了一种未经有意碳增强措施缓冲的生态状态。需要承认的是，"最小化管理"在耕地与草地系统中的具体表现形式存在差异：耕地CK样地保留了耕作和灌溉等基线人为干扰，而草地CK样地则主要受自然放牧或封育条件影响。因此，耕地和草地的背景SOC含量不应被解读为严格等价的生态状态，而应视为具有共同属性——即均排除了有意碳增强措施——的土地利用类型特定基线。在此背景下，我们的发现与包含管理措施的评估之间的对比尤具启发性。这种差异强烈表明，中国耕地表观的碳汇功能主要由有意的管理干预——包括化学施肥、秸秆还田和有机物料改良（Dong等，2021；Liao等，2025）——所维持，而非反映了气候驱动的内在积累过程。一旦去除这些管理带来的增益，背景碳库在耕地总体尺度上并未表现出显著的净变化。然而，这一总体信号不应被解读为气候驱动的SOC损失在中国所有耕地中普遍存在的证据。相反，旱地中SOC高值区域的显著扩张（图6e）表明，在特定条件下气候驱动的SOC积累是可以发生的。因此，尽管已报道的碳汇估计值可能高估了在持续管理投入缺失的情况下耕地SOC的内在韧性，但这一结论在应用于不同土地利用类型和区域时应保持适当的区分。包含管理措施的SOC轨迹与我们仅考虑气候的基线之间的比较，为定量分离管理和气候对SOC变化相对贡献提供了一条路径。然而，我们提醒，这两个组分的定量划分需要各研究之间采用协调一致的SOC指标和一致的土层深度信息，这仍是一个尚未解决的挑战。未来在统一框架下将背景SOC含量与包含管理措施的SOC产品相结合的研究，有望为将碳汇动态归因于气候还是管理驱动因素提供更为稳健的依据。4.3. SOC对气候变化的差异性响应尽管在SOC均值和SOC高值区域比例方面观测到了时间趋势，但2000至2024年间中国农业用地SOC的空间格局始终保持显著稳定性（图5a）。这种空间稳定性表明，在25年的研究期内，气候变异性并未改变SOC的空间分布格局。虽然气候变异性在这一稳定的空间框架内驱动了可检测的时间波动和结构性变化，但高SOC和低SOC区域的地理配置在整个研究期内持续存在。尽管三种土地利用类型的全国SOC均值均显著下降（图5b），但该总体信号代表的是异质性农业用地系统的一种涌现属性，而非所有土地利用类型在统一方向上的变化。全国均值为三种土地利用类型的面积加权平均值，其显著下降主要由水田的方向性损失信号驱动（?0.0041 g C kg?1 y?1，P < 0.001，图6a）。在耕地内部，水田和旱地表现出相反的反应——显著下降与稳定——这一对比可能反映了淹水与非淹水系统对气候强迫不同敏感性的差异。不同生态系统间的独特时间特征反映了生物气候异质性、空间抵消效应以及控制微生物分解与碳稳定的基本生物地球化学约束之间的复杂相互作用。水田主要分布于中国南方（图1a），气候变异性相对较低（Xu等，2017），这种较低的空间异质性可能限制了补偿性缓冲效应，导致全国范围出现空间上一致且可检测的SOC下降。此外，水田SOC的保存严重依赖于连续的厌氧条件以抑制氧化性微生物酶活性（Ma等，2026）。然而，区域增暖叠加年际降水波动的加剧，可能加剧蒸散发作用和干湿交替循环（Condon等，2020；Zhang等，2021）。表土的暂时性通气将分子氧引入根际环境，使微生物代谢转向效率更高的有氧呼吸过程，从而解除酶活性限制并加速活性有机质的氧化矿化（Keiluweit等，2017）。因此，空间缓冲效应的缺乏与增暖导致的通气性损失共同驱动了水田中SOC均值与高SOC覆盖面积的同步退化。相比之下，草地全国SOC均值的表观稳定掩盖了深层的内部空间缓冲效应：草地在横跨中国广阔的生态环境梯度上分布（图1a；Shen等，2016），敏感区域的碳损失被不太敏感区域的增益或相对稳定性所掩盖，这与复杂生态系统中观察到的统计平均效应相一致（Zhao等，2022）。关键在于，即使草地全国均值未出现可检测的变化，SOC高值像素的逐步缩减仍代表一种严峻的结构性退化。SOC高值草地主要分布于青藏高原和中国东北的寒冷高山草甸和温带草甸区域，这些环境以低温、充足的水分供应和高植被生产力为特征（Post等，1982；Wang等，2025）。在这些寒冷且碳含量高的环境中，气候变暖加速了微生物代谢动力学并促进了活性有机质的分解（Crowther等，2016），从而推动了这些高碳库的结构性快速退化。相反，旱地时间上的稳定性与SOC高值区域的局部扩张则源于基质限制与矿物-基质稳定化的协同效应。旱地在所有土地利用类型中SOC含量最低（图

（S1；表S2）中可能包含相对有限的易分解碳库，可供微生物分解利用。因此，气候驱动的下层SOC损失的幅度受到底物可获得性的严格约束（Fissore等，2013年）。综合来看，尽管全国平均水平保持稳定或得到缓冲，高SOC区域的渐进式收缩仍表明，陆地碳库中深刻的结构性重组可能在表观宏观平衡之下悄然发生。SOC平均轨迹与空间结构变化之间的这种脱耦现象，凸显了仅依赖基于平均值的指标来评估陆地碳库对气候变化脆弱性的局限性。

SOC敏感性沿400毫米等降水线（图8a）所呈现的观测一致性与该线在中国长期被认可的基本生态分界角色相符，它将半湿润生态系统与半干旱生态系统分隔开来（高X等，2020年；江X等，2026年）。许多生态过程，包括植被生产力、生态系统稳定性、土壤水分状况和生物地理分布，均在该水气候阈值两侧表现出急剧的转变（马X等，2016年；宋X等，2023年；江X等，2026年）。在全国尺度上，草地的SOC敏感性显著高于耕地（图8b），这与以往将草地确定为全球变化下高度气候敏感生态系统的研究结果一致（Liu等，2021年）。尤为关键的是，在控制背景气候条件后，我们的结果表明土地利用的SOC敏感性存在明显的纬度排序转移（表S4）。在温带至北亚热带区域，草地的SOC敏感性显著高于耕地。相反，从中亚热带到热带区域（以及高原气候区），旱地耕地表现出最高的SOC敏感性。这表明，在温暖潮湿的背景气候下，高强度的SOC耕作显著放大了对热力的敏感性。在所有主要气候区（除中亚热带区和高原气候区外），水田始终表现出最低的SOC敏感性（表S4）。总体而言，对比全国尺度分析与按气候分层的分析结果，突显了广泛尺度评估中固有的空间聚合偏差。这些发现表明，不同类型土地利用对气候变化的响应高度依赖于具体背景条件，并深受各气候区基线水热制度的控制。

在全国尺度上，温度波动在草地和旱地耕地中的SOC敏感性中占主导地位，而降水变率在水田中成为主要驱动因素（图S5a）。然而，沿年降水均值梯度，所有三种生态系统均表现出水文控制与热力控制之间的明显交替（图S5b）。在低降水条件下，降水变率在所有土地利用类型中对SOC敏感性的控制作用最强，表明水分限制构成干旱和半干旱环境中首要的生态约束（Noy-Meir，1973年；Ahlstr?m等，2015年）。这种局部水文主导性得到机制实验的强烈支持，实验证明在干旱和半干旱生态系统中，水分可获得性而非温度对土壤生物地球化学过程具有主导控制作用。具体而言，离散降水事件及其导致的水分变率决定了SOC的季节性矿化和微生物动态，驱动了显著高于升温影响的碳损失脉冲（Waring等，2021年）。随着年降水均值的增加，驱动因素并未简单地从水文主导线性过渡到热力主导。相反，降水变率与温度波动的相对重要性在整个气候梯度上呈现出非线性的振荡行为（图S5b）。值得注意的是，在高降水区域，降水变率仍然保持高度影响力，甚至出现反弹（图S5b）。在湿润环境中反复出现的水文主导性意味着，即使在年均值上水分并非限制性因素，短期降水波动仍可能关键性地调节SOC敏感性，这与Knapp等（2008年）提出的生态水文框架一致，该框架认为年降水内分配格局的改变可深刻重塑土壤水分动态和生态系统碳循环，而不必改变年均降水量。

滑动窗口分析表明，SOC敏感性的主要驱动因素并非空间静止的，而是沿降水梯度动态重组。这一发现表明，调控土壤碳稳定性的机制高度依赖于具体背景，无法仅用单一全国尺度平均值来充分刻画。需要指出的是，本研究未纳入滞后的年际变率指标。在某些生态系统中，这些滞后效应可能很重要，我们将在未来研究中加以探索，以进一步完善对气候-碳耦合机制的理解。

4.4. 局限性

尽管本研究提供了有益见解，但仍有若干局限性需要指出。首先，对照处理（CK）的操作定义不可避免地因土地利用类型而异，从而在数据集中引入了一定程度的概念异质性。耕地CK样地排除了施肥和有机改良措施，但保留了常规耕作和灌溉；而草地CK样地则代表未放牧或常规放牧条件下且无额外投入的状态。这种定义上的不对称意味着，不同类型土地利用间的背景SOC估算反映的是不同的基线扰动制度，因此跨系统直接定量比较应谨慎解读。未来采用跨土地利用类型标准化方案的研究将有助于解决这一定义模糊性。

其次，作为环境协变量的土壤理化性质源自静态数据集（2010–2018年），可能无法完全反映2000–2024年研究期间土壤性质的时间变化。使用时不变土壤数据可能给时空预测引入不确定性，特别是在环境快速变化的区域。未来研究应优先开发和纳入动态土壤性质数据集，以量化预测偏差并进一步提升模型保真度。

第三，观测样本在研究域内的空间分布本质上是不均匀的。尽管我们致力于汇编全国3963个观测点的数据集，但在SOC变率最高的区域（图8a），如中国西北高海拔冷干区和青藏高原（图1a），采样密度仍然相对稀疏。这种不均匀覆盖可能限制模型在未充分代表区域的性能，并限制这些区域趋势估计的普适性。在气候极端和生态敏感区域扩展观测网络将显著改善未来大尺度SOC评估的可靠性。

第四，所汇编的文献跨越四十年，且各引用研究中使用SOC测定方法各异，这可能引入无法在此规模文献综合中直接量化的研究间测量偏差。鉴于并非所有来源出版物均具备完整的方法学元数据，我们无法完全排除更细微的方法相关伪影。未来的综合分析将受益于系统提取方法学元数据并采用协调统一的方案，以实现对方法间偏差的显式定量校正。

5. 结论

本研究建立了包含3963个背景表层土壤有机碳含量观测值的综合全国数据库，最大限度地降低了集约化人为管理的影响，以更好地揭示低投入条件下气候相关的SOC动态。利用整合机器学习框架，我们重构了2000年至2024年背景SOC含量的时空轨迹并量化了其气候敏感性。我们揭示背景表层SOC含量动态从根本上取决于土地利用类型，仅依靠基于平均值的指标不足以刻画碳库脆弱性，突显了在土壤碳评估中纳入空间结构指标的必要性。此外，以往研究报告的中国耕地表观碳汇功能可能主要依赖于持续的管理干预而非内在的气候驱动积累，这意味着如果管理措施得不到维持，未来气候变化下耕地SOC的韧性可能被高估。最后，SOC敏感性的主要气候驱动因素被证明在空间上是非平稳的，沿降水梯度呈现出热力控制与水文控制之间的非线性交替。然而，由于观测覆盖稀疏，中国西北和青藏高原的预测不确定性较高，这些区域估计的SOC动态应谨慎解读。总体而言，本研究为准确评估陆地生态系统中自然碳-气候反馈提供了新视角。

贡献声明

何楠：原稿撰写、方法学、正式分析、数据整理。赵芬：审阅与编辑、监督、方法学、经费获取、概念化。胡国政：审阅与编辑、经费获取、数据整理、概念化。哈斯巴根·甘珠日华：验证、方法学、概念化。丁杰：验证、方法学、数据整理。黄黎明：可视化、软件、方法学。高庆柱：审阅与编辑、经费获取、概念化。

何楠|赵芬|胡国政|哈斯巴根·甘珠日华|丁杰|黄黎明|高庆柱