基于堆叠集成模型对中国重庆市森林地上生物量的估算

《Frontiers in Plant Science》：Estimation of forest above-ground biomass based on stacked ensemble model in Chongqing, China

重庆地区森林地上生物量（AGB）估算研究，提出小胸径（DBH<5cm）树木估算方法，结合Sentinel-2影像和DEM数据，运用随机森林、XGBoost、CatBoost和SVM模型及堆叠集成模型，分析特征重要性及空间分布。包含小胸径树木时，堆叠模型R2达0.68，RMSE为25.45 Mg·ha?1，显著优于其他模型；排除小胸径树木时精度下降明显。光谱波段（B3、B5、B12）和地形指数（坡度、曲率）是主要预测因子，空间分布显示高山地区AGB较高，西部及东北部较低。研究证实小胸径树木不可忽视，堆叠模型提升精度，为西南山区森林生态管理提供数据支撑。

　　森林地上生物量（AGB）的准确区域估算对于有效的森林管理和陆地碳循环研究至关重要。然而，将多种机器学习模型（MLs）集成应用于山地森林的AGB估算仍存在局限。本研究提出了一种实用方法，利用国家森林清查（NFI）数据对DBH小于5厘米的幼树进行估算，并结合Sentinel-2遥感影像和数字高程模型（DEM）数据，采用随机森林（RF）、XgBoost、CatBoost和支持向量机（SVM）等个体模型以及堆叠集成方法，对重庆市的森林AGB进行了估算，分别考虑了是否包含DBH小于5厘米的幼树两种情景。结果表明，DBH估算方法具有高精度（R2=0.93，RMSE=1.46厘米），在交叉验证（CV）和外部验证（EV）中，堆叠集成模型表现出最佳性能（R2=0.65，RMSE=24.34 Mg·ha?1；R2=0.68，RMSE=25.45 Mg·ha?1），生成了10米分辨率的AGB地图，其空间模式具有一致性，适用于山地城市地形。该研究推动了中国西南地区山地森林AGB估算的进展，并为森林生态与管理提供了重要参考。

### 1. 引言

森林在生态系统服务中发挥着重要作用，不仅提供可再生资源和能源，还维持生物多样性、保护水源和防止土壤侵蚀。此外，森林在地球碳循环中也扮演着关键角色，植物光合作用占陆地生态系统碳储存的大约80%（Liu et al., 2021）。作为影响植物光能利用、周转和呼吸的重要因子，森林地上生物量（AGB）是评估森林成熟度和碳汇能力的关键指标（Hao et al., 2020；Shen et al., 2020）。

森林AGB估算方法通常分为实地测量、基于遥感建模和基于过程的模拟模型三类（Thornton et al., 2002；Zhu et al., 2020；Guo et al., 2023）。实地测量提供了直接的数据，是研究森林AGB的重要资源，但其成本高、耗时长，且在空间上难以进行精确建模。基于过程的模型包含详细的生态过程，能够模拟生物量的分配动态，但通常需要大量变量和特定的实地校准，因此在实际应用中面临挑战。由于传统森林清查方法和过程模型的局限性，基于遥感的估算方法在过去几十年中被广泛采用（Goetz et al., 2009；Zhu et al., 2020）。遥感技术取得了显著进展，许多机构启动了多传感器卫星项目（Qian et al., 2021）。准确的森林AGB估算依赖于三个因素：实地数据、影像获取和模型选择（Puliti et al., 2020；Wai et al., 2022；Feng et al., 2025）。

### 2. 材料与方法

#### 2.1 研究区域

重庆市位于长江上游地区，是中国西南部的重要经济中心（见图1）。其地理范围介于东经105°28′至110°19′，北纬28°16′至32°20′之间，总面积约为82,400平方公里。重庆市以其崎岖的山地地形而闻名，主要由丘陵和山地构成。该地区处于湿润亚热带气候区，具有典型的大陆性季风气候。年平均气温在17至18.8°C之间，年平均降水量为1000至1400毫米，年日照时数为1000至1400小时。重庆市的森林类型主要包括常绿阔叶林、次生和暖温针叶林、竹林以及常绿阔叶灌木林。主要树种包括马尾松（Pinus massoniana）、杉木（Cunninghamia lanceolata）、侧柏（Cypress）和壳斗科（Quercus）等。

#### 2.2 数据来源与处理

本研究的实地数据收集于2017年4月至10月，依据国家森林清查（NFI）技术协议进行。采用机械采样方法设置固定样地，间距为4公里×4公里。每个样地为25.82米×25.82米的正方形，面积约为0.0667公顷。测量了所有DBH大于5厘米的树木以及样地内优势树种的高度，并记录了树种、数量、标尺类型、DBH、单株体积和年龄等信息。NFI数据还包含样地四个角点的坐标以及角点间的水平距离。经过数据清洗后，实际可用的样地数量为623个。2021年起，NFI数据收集改为年度测量方案，每年调查所有永久样地的五分之一，五年内完成一次完整的清查周期。因此，2017年测量的样地会在2021、2022或2023年被重新调查，从而可用于估算2017年时DBH小于5厘米的树木的DBH值。

Sentinel-2 MSI级别的1C产品从Copernicus数据空间生态系统下载，时间为2017年4月至10月。Sentinel-2的1C数据包括13个光谱波段，空间分辨率为10米、20米和60米。通过sen2cor（版本2.8.0）进行大气校正，以获得2A级产品。使用屏幕分类层（SCL）产品对云层进行掩膜处理。将20米和60米的影像进行镶嵌和重采样，最终得到10米空间分辨率的12个光谱波段（B10因大气校正被移除）。

DEM数据从欧洲航天局（ESA）下载，空间分辨率为30米，并根据Sentinel-2的瓦片网格进行裁剪和重投影到WGS84 UTM区48N和49N。在本研究中，选择了40个基于遥感的特征变量，通过递归特征消除与交叉验证（RFECV）方法筛选出20个关键变量，这些变量在四个机器学习模型（RF、XgBoost、CatBoost和SVM）中具有共同性。对这些20个遥感特征变量，计算了Pearson相关系数，以定量评估每个特征变量与总森林AGB之间的线性相关性。随后，使用这20个特征变量对每个模型进行训练，并通过100次重复的5折交叉验证进行特征重要性分析。

#### 2.3 方法

##### 2.3.1 DBH小于5厘米树木的估算

在本研究中，基于2017年NFI数据中所有DBH达到5厘米的树木，以及2021-2023年NFI数据中DBH达到5厘米的树木，估算了2017年时DBH小于5厘米的树木的DBH值（见图2）。首先，根据树木ID计算了80%保留树木的DBH生长率，然后按样地、树种和DBH等级（包括5-10厘米、10-15厘米、15-20厘米和大于20厘米）进行分组，估算出每组的平均生长率。剩余20%保留树木的估算DBH值与实际测量值进行比较，以验证方法的性能。结果显示，该方法具有良好的结果，可以用于进一步的分析（见图3）。最终，将该方法应用于所有保留树木，计算其生长率，并从2021-2023年NFI数据中减去这些生长率，以得到2017年时DBH小于5厘米树木的DBH值。

##### 2.3.2 样地级AGB计算

样地级AGB计算采用Zeng和Tang的方法（Zeng and Tang, 2011）。Zeng和Tang（2011）将他们的方法与全球不同研究进行了比较，表明理论参数值7/3能够描述不同树种AGB与DBH之间的平均全量化关系。AGB与DBH之间的系数可以通过将0.3乘以木材密度ρ来获得。计算公式为：AGB = 0.3 × ρ × D^(7/3)。样地级AGB是样地内所有树木AGB的总和。不同树种的ρ值如表2和补充表S1所示。大多数树种的ρ值参考了Zeng（2018）的研究，而其他未被涵盖的树种则采用《中国林业行业标准》LY/T 3256-2021中“主要树种（群）木材密度测试”提供的数据。

##### 2.3.3 特征选择与变量重要性

本研究采用了递归特征消除与交叉验证（RFECV）方法，这是一种稳健的特征选择方法，通过递归移除最不重要的特征，并在每一步中使用交叉验证评估模型性能，从而筛选出40个遥感指标中的20个关键变量，这些变量在四个机器学习模型（RF、XgBoost、CatBoost和SVM）中具有共同性（见图4）。对于这20个遥感指标，计算了Pearson相关系数，以定量评估每个特征变量与总森林AGB之间的线性相关性。随后，使用这20个特征变量对每个模型进行训练，并通过100次重复的5折交叉验证进行特征重要性分析。

##### 2.3.4 交叉验证与模型参数调优

本研究采用网格搜索结合交叉验证的方法进行超参数优化。网格搜索是一种穷举搜索方法，通过遍历所有可能的超参数组合寻找最优的超参数。不同模型的超参数范围见补充表S2。为了提高模型的可靠性和泛化能力，在超参数调优过程中采用了5折交叉验证策略。原始训练数据集被随机划分为五个相等的子集。在每次交叉验证迭代中，四个子集用于模型训练，剩余一个子集用于性能评估。此过程重复五次，每个子集恰好作为一次验证集。通过优化过程确定的最优超参数组合如图5和表3所示。

##### 2.3.5 堆叠集成模型的开发

本研究采用Ridge回归（RR）作为堆叠集成模型的元模型。在堆叠过程中，使用基于样本的交叉验证（CV）和外部验证（EV）进行性能评估。本研究采用5折交叉验证来测试模型的稳健性，所有站点样本被随机划分为五个子集。每次交叉验证中，模型在四个子集上进行训练，并在剩余子集上进行测试。相比之下，EV实验评估模型的泛化能力（即真实预测性能）在未参与任何模型训练过程的独立数据集上进行。在每次5折交叉验证迭代中，四个基础模型在相同训练集上并行训练，并在测试集上生成预测结果。这些基础模型的预测结果作为新特征，与目标变量一起训练元模型（RR）。随后，外部测试值输入基础模型进行训练，结果被输入到训练好的堆叠集成模型中进行预测。最终的预测结果经过验证（见图6）。

##### 2.3.6 模型评估与不确定性分析

模型性能通过决定系数（R2）、平均绝对误差（MAE）和均方根误差（RMSE）以及均方根误差百分比（RMSE%）进行评估。计算公式如方程2至5所示。其中，y和y?分别表示实际值和预测值，y?是实际值的平均值，n是训练数据集的数量。模型的不确定性由交叉验证中的预测值计算得出，计算公式如方程6所示。其中，n为5，p?是单折交叉验证的预测值，p?是5折交叉验证预测值的平均值。研究框架如图7所示。

### 3. 结果

#### 3.1 DBH小于5厘米树木的估算

2017年时，DBH小于5厘米的树木占样地数据中所有树木的近23%（14697株/60786株）（见图8A）。这些树木的估算DBH分布与NFI技术协议一致，大部分值集中在5厘米左右（见图8B）。分布的平均值和中位数分别为4.983厘米和4.533厘米，四分位数范围为3.646厘米至5.687厘米。一些超出10厘米的异常值可能是由于相邻NFI测量之间的边界定位误差或数据采集操作误差所致。

#### 3.2 样地级AGB的空间分布

本研究中使用的623个NFI样地的空间分布如图9所示。样地在重庆市分布相对均匀，以确保空间覆盖的代表性。通过全量化方程计算的样地级AGB值范围为5.15至329.81 Mg·ha?1，平均值为66.63 Mg·ha?1。这些样地级AGB值的大部分集中在33.38至90.48 Mg·ha?1之间（见表4）。排除DBH小于5厘米的树木后，样地级AGB的中位数和平均值相对较小，分别为55.08 Mg·ha?1和64.48 Mg·ha?1。整体范围从0.62 Mg·ha?1到329.81 Mg·ha?1。重庆市西部的样地级AGB值主要小于50 Mg·ha?1，而超过150 Mg·ha?1的AGB值主要集中在东北至西南方向的山地范围内（见图9）。考虑DBH小于5厘米的树木后，超过500个样地的AGB值变化在10%以内，其余样地的变化超过10%（见图10）。

#### 3.3 特征重要性与相关性

在RF模型中，AGB预测主要依赖于光谱波段，植被指数则处于次要地位（见图11A）。值得注意的是，B12、B4、B3和B5的相对重要性分别为12.1%、9.5%、8.9%和7.7%，均与AGB具有显著的负相关性。由光谱波段导出的植被指数，如MVI、MSI和NDII，也表现出较高的重要性。相比之下，地形特征如曲率和高程的重要性相对较弱。

在XgBoost模型中，光谱波段在AGB预测中占据关键地位，B12、B5、B3和B4被识别为最重要的特征，分别贡献了19.9%、9.45%、8.6%和8.5%的预测能力（见图11B）。植被指数如MVI和MSI也表现出显著的重要性，分别贡献了5.9%和5.8%。地形特征如曲率和高程对AGB预测具有增强作用。此外，生物物理指数如CWC表现出中等重要性，对AGB估算具有正反馈效应。

在CatBoost模型中，光谱波段和植被指数仍然是最重要的特征（见图11C）。其中，B5（7.1%）和MVI（6.8%）表现出最高的相对重要性，但它们对AGB的反馈效应相互抵消。这些特征的贡献度随后由TPI主导，而其他地形特征如曲率和高程则表现出中等重要性。相比之下，B12、NDII和B3也被识别为高度重要的特征。

在SVM模型中，B11（16.6%）和NDVI（15.2%）是主要贡献因素，而植被指数如CIgreen（9.5%）和MCARI（9.4%）则处于次要地位。地形特征的重要性相对较低。值得注意的是，其他三个模型中表现出重要性的B2（1.8%）、B4（1.1%）和B5（0.9%）在SVM模型中表现出较低的重要性（见图11D）。

#### 3.4 AGB估算模型的性能

在包含DBH小于5厘米树木的情况下，CatBoost模型在所有个体模型中表现最佳，平均R2为0.64（四分位数范围：0.639–0.647），平均RMSE为25.15 Mg·ha?1（四分位数范围：24.99–25.28 Mg·ha?1）（见图12C）。相比之下，XgBoost表现出最低的准确性，平均R2为0.62（四分位数范围：0.616–0.633），平均RMSE为25.97 Mg·ha?1（四分位数范围：25.66–26.27 Mg·ha?1）（见图12B）。RF模型表现出中等性能，平均R2为0.64（四分位数范围：0.635–0.645），平均RMSE为25.26 Mg·ha?1（四分位数范围：25.11–25.39 Mg·ha?1）（见图12A）。SVM的性能略低于RF，平均R2为0.64（四分位数范围：0.634–0.643），平均RMSE为25.70 Mg·ha?1（四分位数范围：25.48–25.87 Mg·ha?1）（见图12D）。堆叠集成模型在预测准确性、泛化能力和稳健性方面有所提高，平均R2为0.65（四分位数范围：0.646–0.657），RMSE为24.38 Mg·ha?1（四分位数范围：24.22–24.56 Mg·ha?1）（见图12E）。

在排除DBH小于5厘米树木的情况下，堆叠集成模型在预测准确性方面优于个体模型，平均R2为0.65（四分位数范围：0.646–0.656），平均RMSE为25.58 Mg·ha?1（四分位数范围：25.38–25.78 Mg·ha?1）（见图13E）。在个体模型中，CatBoost的平均R2为0.64（四分位数范围：0.632–0.642），平均RMSE为26.41 Mg·ha?1（四分位数范围：26.23–26.59 Mg·ha?1）（见图13C）。RF模型表现出稍低的准确性，平均R2为0.63（四分位数范围：0.629–0.639），平均RMSE为26.53 Mg·ha?1（四分位数范围：26.32–26.68 Mg·ha?1）（见图13A）。XgBoost表现出最弱的性能，平均R2为0.61（四分位数范围：0.605–0.623），平均RMSE为27.04 Mg·ha?1（四分位数范围：26.71–27.36 Mg·ha?1）（见图13B）。SVM的平均R2和RMSE分别为0.62（四分位数范围：0.621–0.632）和26.80 Mg·ha?1（四分位数范围：26.66–26.93 Mg·ha?1）（见图13D）。

#### 3.5 AGB地图与不确定性分析

10米分辨率的AGB地图揭示了包含和排除DBH小于5厘米树木时的明显差异。包含DBH小于5厘米树木时，RF、XgBoost、CatBoost、SVM和堆叠集成模型的平均AGB值分别为61.25 Mg·ha?1、59.15 Mg·ha?1、60.76 Mg·ha?1、61.48 Mg·ha?1和60.38 Mg·ha?1，对应的总AGB值分别为3.58×10? Mg、3.26×10? Mg、3.35×10? Mg、3.39×10? Mg和3.33×10? Mg（见图16A-E）。排除DBH小于5厘米树木时，RF、XgBoost、CatBoost、SVM和堆叠集成模型的平均AGB值分别降至59.29 Mg·ha?1、57.08 Mg·ha?1、58.61 Mg·ha?1、61.11 Mg·ha?1和58.70 Mg·ha?1，对应的总AGB值分别为3.27×10? Mg、3.15×10? Mg、3.23×10? Mg、3.37×10? Mg和3.24×10? Mg（见补充图S1A-E）。包含和排除DBH小于5厘米树木的AGB地图的空间分布模式相似：重庆市西部的AGB值较低，而东北至西南方向的山地区域AGB值较高。此外，包含DBH小于5厘米树木时，模型的不确定性较低，RF、XgBoost、CatBoost和SVM的平均不确定性分别为2.68 Mg·ha?1、5.44 Mg·ha?1、2.87 Mg·ha?1和4.30 Mg·ha?1。堆叠集成模型的不确定性进一步降低至3.04 Mg·ha?1（见图17A-E）。在排除DBH小于5厘米树木的情况下，所有模型的R2值显著降低，RMSE值显著增加，表明预测准确性下降。这些结果突显了DBH小于5厘米树木在AGB估算中的重要性。

### 4. 讨论

#### 4.1 DBH小于5厘米树木应纳入AGB估算

根据NFI技术协议，所有树木（包括经济树木）、竹子（包括非竹林中的竹子）和树状灌木（不包括灌木状树种）在永久样地中应测量其DBH是否达到5厘米。此前利用NFI数据进行森林AGB估算的研究通常未关注DBH小于5厘米树木的影响（Zhu et al., 2020；Qian et al., 2021；Li et al., 2022；Zhang et al., 2023），因为研究人员通常认为这些树木的生物量贡献可以忽略不计。然而，本研究通过比较两个连续的NFI数据集，并采用开发的DBH估算方法，估算出2017年时DBH小于5厘米树木的数量。结果表明，DBH小于5厘米的树木占总树木数量的23%，其中一些样地的AGB变化超过10%（见图9和图10）。尽管DBH小于5厘米树木的AGB变化较小，但这些树木主要由阔叶树种组成，如刺柏（Cupressus funebris）、青冈栎（Cyclobalanopsis glauca）、马尾松（Pinus massoniana）、杉木（Cunninghamia lanceolata）和壳斗科（Quercus）等。Yang等人（2012）和Wen等人（2015）的研究发现，在中国南方和西南地区，随着森林进入中年阶段，新的野生树种开始出现在马尾松、青冈栎等阔叶林中，并在成熟森林中存在一定程度的群落演替。根据本研究的2017年森林资源调查数据，中年林、近成熟林和成熟林共同占森林小班的显著比例（42%）。此外，排除DBH小于5厘米的树木显著降低了所有模型的预测准确性。因此，在利用NFI数据进行森林AGB遥感反演研究时，对DBH小于5厘米树木的校准是至关重要的。

#### 4.2 特征对AGB估算的重要性

在所有模型中，光谱波段是AGB预测的主要特征，B3、B5和B12表现出一致的高重要性，这一发现与Wai等人（2022）的研究一致（见图11）。B5，红边波段，表现出检测关键植被生理参数（如叶绿素含量和冠层结构）的强大能力，这对于准确的AGB估算至关重要（Zhang et al., 2023）。相比之下，其他红边波段（B6、B7）在先前研究中表现出更大的重要性（Yang et al., 2012）。B3位于叶绿素反射峰（500–600纳米），对植被“绿色度”高度敏感，并有效反映叶绿素含量和光合作用活动等生理状态。同样，B2和B4也显示出在RF和XgBoost模型中的显著重要性（Wen et al., 2015）。短波红外（SWIR）波段B12在AGB估算中的有效性已被广泛证明。B12对植被含水量敏感，随着生物量和植被覆盖度的增加，植被或水体对光谱的吸收会减少反射，而在低生物量区域，反射度较高，这解释了B12与AGB之间的负相关性（Wai et al., 2022）。这些发现强调了SWIR、可见光和红边波段对生物量的高敏感性，突显了它们在生物量评估中的关键作用。

植被指数的重要性在不同模型中有所差异。MVI在所有三个模型中都表现出显著的重要性，主要反映了森林冠层特征，对生物量积累有重要贡献。MSI和NDII也显示出显著的重要性。MSI通常表现出对AGB的负反馈效应，这是由于植被对水分的依赖以及水分胁迫对生长和生理过程的影响。NDII和SR主要反映了叶片叶绿素和类胡萝卜素含量的变化，反映了光合作用活动，成为AGB遥感反演的关键参数（Richardson et al., 2002；Merzlyak et al., 2003；Main et al., 2011）。

地形特征，如TPI、曲率和高程，在CatBoost模型中表现出较强的预测能力。重庆市的山地地形具有多样化的地貌特征，影响阳光照射和水分保持，从而影响植被生长（Wai et al., 2022）。TPI和曲率有效捕捉了这些地表变化。高程和坡向也显示出显著的重要性，与先前研究一致（Chen et al., 2019b；Wang et al., 2021）。在山地地区，高程和坡向的变化影响湿度、温度和物种丰富度，从而影响植被生物量（Shen et al., 2014；Cong et al., 2019）。

#### 4.3 AGB估算模型的性能

在个体模型评估中，CatBoost表现出最佳性能，达到R2=0.66，随后是RF和XgBoost，而SVM的性能相对较低（见图12）。此前的研究表明，基于树的机器学习方法特别适合生态遥感研究（Belgiu and Dragut, 2016）。本研究中评估的RF、CatBoost和XgBoost模型均为基于决策树的集成方法，它们生成的AGB地图在空间上具有高度一致性。与RF和XgBoost相比，CatBoost是估算AGB的更好选择，因为它采用有序提升机制，降低了过拟合的风险，并通过特定顺序处理训练样本来减少噪声数据的影响。由于NFI数据通常包含内在的噪声和变量性，CatBoost在处理分类特征和提高泛化能力方面的优势使其特别适合在复杂和噪声数据集中准确估算AGB。基于个体模型，开发了一个使用RR作为元模型的堆叠集成模型。堆叠集成模型显著提高了预测准确性和泛化能力，达到R2=0.68，这明显高于排除DBH小于5厘米树木时的R2=0.59。这些发现为探索混合机器学习方法在森林AGB估算中的潜在应用奠定了基础。

### 5. 结论

本研究开发了一种基于NFI数据的DBH估算方法，用于估算DBH小于5厘米的树木，并结合遥感影像和地形数据，比较了包含和排除DBH小于5厘米树木的个体和堆叠集成模型的性能，最终生成了重庆市的10米分辨率森林AGB地图。所开发的DBH估算方法表现出高精度（R2=0.93，RMSE=1.46厘米）。考虑到DBH小于5厘米的树木占总树木数量的23%，其排除显著影响了AGB预测的准确性。因此，在利用NFI数据进行森林AGB遥感反演研究时，对DBH小于5厘米树木的校准是至关重要的。

在所有模型中，光谱波段是AGB预测的主要特征，而植被和地形指数在不同模型中的重要性表现出显著差异。因此，针对特定模型选择不同的特征变量有助于提高预测准确性。堆叠集成模型的性能优于个体模型。尽管四个个体模型在包含DBH小于5厘米树木时的R2值在0.64至0.66之间，但堆叠集成模型有效减少了模型间变异性，提高了预测准确性（R2=0.68），这明显高于排除DBH小于5厘米树木时的R2=0.59。这些发现为探索混合机器学习方法在森林AGB估算中的潜在应用奠定了基础。