Aboveground biomass estimation of major coniferous forests in Qujing City based on Landsat images

云南省曲靖市位于24°21′00″~27°03′24″N，103°03′44″~104°32′36″E，土地总面积约28 900.0 km²，平均海拔约1 881.0 m。属亚热带季风气候，年平均气温为14.8 ℃，年均降水量为560.0~1 026.0 mm。主要植被类型为亚热带针叶林和常绿阔叶林，原生植被多已退化，现以天然次生植被及人工林为主。曲靖市森林资源较为丰富，2021年森林覆盖率达50.5%。优势树种主要包括云南松Pinus yunnanensis、杉木Cunninghamia lanceolata、华山松Pinus armandii、栎类Quercus spp. (如麻栎Quercus acutissima、青冈Quercus glauca等)及油茶Camellia oleifera等。

获取云南省曲靖市1992—2017年6期(5 a为1期)森林资源连续清查数据。云南省森林资源连续清查采用系统抽样方法^[15]，抽样间距为6 km × 8 km，方形样地面积为0.06 hm²。以曲靖市优势树种为研究对象，筛选出3类主要针叶林：云南松林、华山松林和针叶混交林，共计417个样地(其中云南松林264个、华山松林83个、针叶混交林70个)。以样地蓄积量为基础计算森林生物量^[16]，蓄积量由一类连续清查数据提供，依据样木胸径(DBH)、树高及现有生长量方程计算得到。通过蓄积量-生物量模型^[17]计算得到森林生物量为实测数据。计算公式为：B=aV+b。其中：B为森林生物量(t·hm⁻²)；V为森林蓄积量(m³·hm⁻²)；a、b为蓄积量与生物量之间的模型转换参数^[18]。

获取2007年云南省曲靖市森林资源二类调查数据^[19]，用于确定后续生物量制图范围。一类连续清查样地数据经几何校正后，重采样至30 m以匹配遥感影像。遥感数据来源于谷歌地球引擎(GEE)平台的Landsat TM/OLI影像(Collection 2 Tier 1)，时间范围为1992—2022年，每年1月1日至12月31日的影像，空间分辨率为30 m。该数据集已完成大气校正、几何校正和地形校正，经去云处理与时间序列插值，最终采用中值合成生成年度影像。地形数据来自地理空间数据云的DEM (空间分辨率为30 m)。

基于Landsat影像的光谱响应、冠层几何特征及灰度共生矩阵(GLCM)描述的林分空间异质性，提取与AGB显著相关的遥感因子。以Landsat 5 TM影像为例，包括①单波段因子^[20]：选取TM1~TM5和TM7共6个原始波段作为单波段因子，TM1区分土壤与植被，TM2区分植被类型，TM3 (红光)反映叶绿素吸收与冠层结构，TM4 (近红外光)对叶片结构和水分敏感，是生物量估算基础，TM5 (短波红外光)反映干物质积累并辅助判别岩土，TM7 (短波红外光)用于识别岩石与植被胁迫；②波段组合因子：比值指数(如 TM42、TM53、TM73等)可有效削弱地形光照干扰，增强植被信号，显著提升森林生物量反演精度；③DEM地形因子：坡度、坡向和高程，可调节光合有效辐射(PAR)并作为温度与水分的综合代理，从而间接影响生物量的空间分布。④纹理特征：对6个波段分别在3×3、5×5和7×7等3种窗口下提取GLCM的8类指标^[16]，包括均值、方差、均一性、对比度、相关性、相异性、角二矩和熵，纹理指标如均值、对比度和熵可量化林分内部树冠大小、郁闭度及空间复杂性，且已被证实与AGB高度相关，最终提取144个纹理变量；⑤植被指数涵盖多种典型指标^[21]：如归一化植被指数(NDVI)、增强型植被指数(EVI)、土壤调节植被指数(SAVI)、垂直植被指数(PVI)、土壤调整比值植被指数(SARV)等，其中NDVI和EVI与叶面积指数(LAI)及光合有效辐射吸收比例(fAPAR)密切相关，广泛用于AGB反演。⑥主成分变换因子：对原始6个波段进行K-L变换，保留信息量最大的前3个主成分(PCA1、PCA2、PCA3)，在降维的同时保留与生物量变化相关的综合光谱信息，从而提高模型的稳健性。变量命名规则：b1~b7表示Landsat TM影像的单波段因子；纹理特征命名规则中，数字前缀表示窗口大小，字母及后缀表示所用波段和纹理指标类型，例如R3Hom2表示3×3窗口下第2波段的均一性(homogeneity, Hom)纹理因子。

为降低高维特征带来的计算负担与过拟合风险，提升AGB估测的精度与稳定性，研究在初步提取187个遥感特征变量基础上，采用二阶段特征选择策略进行优化。首先，通过皮尔逊(Pearson)相关性分析筛选与AGB显著相关(P＜0.05)的变量，剔除无效信息；随后结合RF的特征重要性，利用递归特征消除(RFE)算法优选10个关键变量用于建模^[22]。

RF算法是BREIMAN^[23]于2001年提出的一种集成学习方法，由多个决策树构成。该模型通过在训练过程中对样本和特征进行随机抽样，构建多个相互独立的决策树，并最终通过集成(如投票或平均)各棵树的预测结果，提升模型的整体预测精度和稳定性。作为一种非参数模型，RF具有出色的抗过拟合能力、良好的抗噪声性能和较强的泛化能力，同时支持并行计算，在计算效率和预测精度方面均表现优异，在回归任务中，RF通过多棵回归树的集成，有效提高对连续变量的建模能力^[24]。

XGBoost算法是CHEN等^[25]在传统梯度提升方法的基础上提出的一种高效的集成学习算法。在模型结构与计算机制上均进行优化的增强型梯度提升算法，具备显著的训练效率与预测性能优势，尤其适用于处理大规模、高维度的数据集。在传统梯度提升决策树(GBDT)框架的基础上引入正则化项，有效控制模型复杂度，从而缓解过拟合问题。同时，XGBoost采用基于二阶导数的信息进行损失函数的优化更新，进一步提升了模型的拟合精度与收敛速度。

GBRT算法与RF类似，是一种基于集成学习思想的回归建模方法，由FRIEDMAN^[26]首次提出，并在后续研究中不断优化与拓展^[27]。GBRT通过集成弱回归树并逐步拟合残差，基于损失函数梯度优化预测性能。相比随机森林，其更适用于处理非线性、高维和复杂特征关系，且对异常值具有鲁棒性。在遥感估测、生物量反演及碳储量建模等复杂生态问题中，GBRT具有优越的预测与泛化能力。

PLS算法是一种结合主成分分析(PCA)和多元线性回归思想的多变量建模方法。不仅保留主成分分析在降维处理中的优势，有效解决自变量之间多重共线性问题，还进一步考虑各主成分与因变量之间的相关性，在提取特征时同时最大化自变量与因变量之间的协方差，从而提高模型的预测能力^[28]。

CatBoost是一种基于梯度提升决策树的高效机器学习框架，针对类别特征做出优化，无需手动编码，自动高效处理类别变量，提升建模效率与预测精度。采用有序目标统计方法和先验项，缓解目标泄露引起的过拟合；引入对称树结构和特有损失函数优化，提高训练速度与稳定性。核心算法基于梯度提升框架，结合有序Boosting策略和类别特征编码，调整样本顺序敏感性。CatBoost凭借优异的泛化能力和计算性能，在处理复杂结构、类别变量多和数据分布不均的问题中表现突出。

将优选的特征变量带入模型的构建，采用随机抽样的方法，以其中80%的样地数据作为训练数据，其余20%的样地数据作为检验数据。除PLS外，其他模型均采用网格搜索结合五折交叉验证进行超参数调优，通过遍历参数空间，筛选最优参数，确保模型具备良好的鲁棒性和泛化能力。此外，为加快模型训练与参数优化过程，所有可并行算法均设置n_jobs=−1，充分调用多线程计算资源。依据各模型算法特点和树种数据特征设置合理的参数范围。GBRT模型主要调节基础学习器数量(n_estimators)；RF模型进一步考虑最大树深(max_depth)与节点最小分裂样本数(min_samples_split)；XGBoost模型引入学习率(learning_rate)与树结构复杂度控制；PLS模型采用固定成分数(n_components=5)，CatBoost模型则基于迭代次数(iterations)进行性能优化。

采用多种评价指标对模型拟合效果进行综合评估，包括调整决定系数(R²)、相对均方根误差(rRMSE)、均方根误差(RMSE)^[29]以及平均绝对误差(MAE)^[30−31]。其中R²用于衡量解释能力并校正自变量影响；RMSE反映预测误差大小；rRMSE为标准化指标，便于不同模型误差比较；MAE表示所有预测误差的绝对值的平均值，能直观反映模型预测值与真实值之间的平均偏离程度^[31]。

基于样地数据，采用Pearson相关性分析对187个遥感特征变量与实测森林地上生物量AGB之间的相关性进行检验，剔除P≥0.05的变量，初步筛选出显著相关因子。在此基础上，引入基于RFE进行进一步特征优化。结果显示(图1)：针叶混交林权重最高的变量为R7Hom2、b2、R7Sec2等；云南松林为aspect、slope、ND54等；华山松林为ND53、TM54、ARVI等，不同林分间特征变量组合存在差异。

Figure 1.  Correlation analysis of remote sensing characteristic variables of 3 coniferous forests

如图2所示：GBRT模型在针叶混交林、云南松林和华山松林3种林型的AGB估测中均表现良好(R²＞0.77)，其中华山松林拟合最佳(R²=0.82)，云南松林误差最小(RMSE=8.30 t·hm⁻²)；RF模型对针阔混交林、云南松林和华山松林等3种林型的预测精度略低于GBRT，R²分别为0.71、0.62、0.74，RMSE分别为9.64、11.57和11.27 t·hm⁻²；XGBoost模型中，针叶混交林R²=0.73，估测趋势较好，云南松林R²=0.61、RMSE=11.75，部分样点偏差较大，华山松林拟合略优但在高值区出现高低估交叉；PLS模型的3种林型R²普遍偏低，拟合能力较差，变异系数分别为0.36 (针叶混交林)、0.28 (云南松林)和0.50 (华山松林)，表明生物量差异明显；CatBoost模型在针叶混交林表现最佳(R²=0.90，RMSE=5.66)，拟合精度高，散点集中，优于其他模型，云南松林与华山松林R²分别为0.60和0.69，RMSE为11.88和12.37 t·hm⁻²。模型差异可能受样地变异性、遥感特征质量及样本代表性等因素影响。综上所述，模型适应性因树种而异，CatBoost在针叶混交林中表现最佳(R²=0.90，RMSE=5.66)，具备优异的非线性拟合能力；云南松林和华山松林则以GBRT效果最佳，稳定且解释力强。可见，模型选择需结合树种特性和数据结构，优选最优算法以提升生物量反演准确性和可靠性。综上，CatBoost适用于针叶混交林地上生物量估算，GBRT更适合云南松林与华山松林的遥感建模与全域逐像元反演制图。

Figure 2.  Scatter plots of the results of each model for biomass estimation

基于Landsat影像数据，采用CatBoost和GBRT模型对曲靖市主要针叶林地上生物量进行30 m×30 m分辨率的空间估算。如图3所示：将森林生物量分为低生物量区(1.78~20.00 t·hm⁻²)、中低生物量区(＞20.00~40.00 t·hm⁻²)、中等生物量区(＞40.00~60.00 t·hm⁻²)、中高生物量区(＞60.00~80.00 t·hm⁻²)及高生物量区(＞80.00 t·hm⁻²)。1992—2022年，曲靖市针叶林地上生物量呈持续增长态势，空间分布格局显著优化。1992年，区域生物量集中于20.00~40.00 t·hm⁻²及以下，高生物量区分布零散，反映当时森林质量较差，受林地退化和过度采伐等人为干扰影响明显。1997年，东南部部分区域生物量回升，40.00~60.00 t·hm⁻²区比例有所增加，表明中幼龄林逐步成熟，林分结构改善及天然林保护政策的实施促进碳储量提升。2002年以后，生物量增长加速，中高生物量区(＞60.00 t·hm⁻²)开始连片分布，低生物量区面积明显缩减，显示封山育林和人工造林等生态修复措施取得显著成效。2007年，高生物量区(＞80.00 t·hm⁻²)面积大幅扩展，集中于西部和南部高山区，生态系统生产力显著增强，该阶段变化得益于“退耕还林”和“生态公益林”等国家重点生态工程的持续推进。2012年，生物量分布趋于均衡，区域内部差异缩小，人工林进入成材期，天然林保护成果逐步显现。2017年，生物量整体继续提升，高生物量区范围进一步扩大，尤其在西部和南部山区呈现更为明显的连片分布，低生物量区进一步收缩，显示生态恢复效果持续深化，森林碳汇能力不断增强。至2022年，针叶林地上生物量达到历史最高水平，高生物量区在西部及南部山区形成大面积连片，区域碳储能力显著提升。相比1992年，生物量空间格局由“低量为主、斑块状高量”转变为“中高量为主、分布连续均衡”。

Figure 3.  Estimation results of aboveground biomass of major coniferous forests in Qujing City

本研究表明：CatBoost模型在针叶混交林中估算性能最优(R²=0.90)，该结果与LUO等^[32]在吉林省基于RFE特征选择与CatBoost的AGB估测结论一致，CatBoost的精度均优于XGBoost和RF等估测模型。CatBoo st模型能够自适应学习高维非线性特征，并有效抑制过拟合，尤其适用于针叶混交林冠层结构复杂和物种组成异质性高的场景^[32]。在结构相对单一、冠层性状较一致的云南松林和华山松林中，GBRT模型表现出更稳健的估算能力，R²分别达0.80和0.82，此类林分中林木个体在高度、胸径和冠幅等表型特征差异较小，导致遥感信号在像元尺度上内部异质性低、统计分布平稳，有利于提升模型拟合效果。这一发现与ZHANG等^[33]在结构单一的高山松Pinus densata林中研究结果一致，其GBRT模型R²达0.94，显著优于其他算法，验证GBRT在结构一致性高的针叶林中具有更强的稳健性与泛化能力。进一步验证了“林型-算法”适配策略的有效性，针对特定林型结构特征选择适宜的机器学习算法，可显著提升区域森林生物量估算的精度与可靠性，为构建差异化生物量反演模型提供关键的方法学依据。

1992—2022年，曲靖市针叶林AGB总体持续增长，空间格局由零散低值逐步演变为中高值连片分布。低生物量区明显缩减，中高值斑块自2000年起快速扩展并趋于连通，与史川等^[34]基于2002—2017年多期遥感和样地数据的研究结果一致，滇中地区(含曲靖市)针叶林碳储量总体持续上升，处于显著的碳汇增强阶段。2012年后高生物量区进一步扩展并趋于连通，体现出持续的管护与经营增汇效应；但陡坡及立地条件较差区域，增幅有限，呈现明显的空间异质性。2009—2012年极端干旱事件曾导致中龄林生长停滞和高生物量区出现短期回落，揭示了森林恢复过程存在显著的空间异质性和气候敏感性^[35]。

本研究采用二阶段特征选择策略，借鉴相关性初筛结合RFE的集成特征选择策略^[32]，先按Pearson初筛将187个原始特征缩减至106个，显著缓解高维共线性问题^[33]；随后在此基础上引入RFE，筛选出与各林型显著相关的10个关键变量，在实现高变量压缩率(＞70%)与增强模型稳健性的同时，也表明该策略具有良好的跨区域、跨林型普适性与可迁移性。其次，时间序列遥感数据通过整合多时相信息，有效削弱单一时相数据中不确定性因素的干扰，提升生物量估算模型的准确度与鲁棒性^[36]。然而，生物量空间异质性显著，其估算结果仍存在一定不确定性，与岳彩荣^[37]关于香格里拉高山松林生物量估算的研究结论一致，生物量的空间分布属性随时间动态变化，这种动态特征成为不确定性的重要来源，在长时序研究中尤为突出。本研究存在若干局限性：所使用的Landsat TM/OLI影像在传感器差异、空间分辨率偏低以及云处理后时间插值等方面存在局限，对预测结果会产生一定影响^[17]；森林资源清查样地数据年代偏旧且样本有限，后续研究可引入更新、更丰富的数据源；当前研究仅估算针叶林乔木层的AGB，未涵盖其他植被类型及土壤碳库，导致区域碳储量结果偏低。这一偏差得到已有研究的支持，如张加龙等^[38]指出，忽略地下生物量与枯死木碳库可能导致碳储量被低估约10%~15%。此外，GBRT与CatBoost模型在森林碳储量估算中的研究较少^[39]，其可行性与适用性仍需进一步探讨。后续应综合考虑树种差异与地形异质性，整合GF-6、LiDAR等多源高分辨率遥感数据，及多源遥感特征变量对模型性能的影响，以期提升碳储量空间分布估算的准确性。

不同树种针叶林对遥感因子的响应存在明显差异，呈现特异性变量组合，其中，针叶混交林主要依赖高阶纹理和原始光谱，云南松林更敏感于地形因子和植被指数，华山松林则体现植被指数与纹理的复合响应，揭示其光谱结构和地形依赖的异质性。基于GBRT、RF、XGBoost、PLS和CatBoost等5种模型构建生物量估测，针叶混交林以CatBoost表现最佳，云南松林和华山松林则以GBRT模型效果优异，可见需结合树种特性选择相应合适的模型。1992—2022年，曲靖市针叶林生物量持续增长，空间分布显著优化，早期低值斑块广泛，碳汇能力弱，2002年起生态修复促使中高值区扩展，2022年生物量达峰，显示森林生态系统由退化向恢复转变，体现了政策与自然恢复共同提升碳汇能力。