中文 
English 
Spatial Heterogeneity Analysis of Forest Biomass Based on Spaceborne LiDAR ICESat–2/ATLAS and Geostatistics
-
摘要:
以ICESat–2/ATLAS数据为数据源，结合54块实测样地，构建机器学习模型并对光斑足迹的地上生物量进行预测，采用Moran's I和半变异函数对反演的森林AGB空间自相关和异质性进行研究。结果表明：梯度提升回归树（GBRT）模型具有较好的预测精度（R2=0.90，RMSE=11.08 t/hm2）；香格里拉市森林生物量的最佳拟合半变异函数模型为指数模型（C0=0.12，C0 + C=0.87，A0=
10200 m）；与普通克里格相比，序贯高斯条件模拟得到的AGB空间分布图具有较好的一致性（r=0.59**，d=0.70）。AGB的空间分异能够被地形因子解释，在解释力方面，海拔最大，坡向次之，坡度最小；基于星载激光雷达ICESat–2/ATLAS数据的森林AGB反演精度较高（Pp=81.43%），为地统计分析提供了可靠的数据源。因此，基于星载激光雷达与地统计学相结合的方法，能较好地实现森林AGB的空间异质性分析。Abstract:Using ICESat–2/ATLAS data as data source, combined with 54 measured plots, a machine learning model was built and the AGB of the spot footprint was predicted. Moran's I and semi-variogram function were used to study the spatial autocorrelation and heterogeneity of inverse forest AGB. The results showed that the Gradient Boost Regression Tree(GBRT) model had a great prediction accuracy(R2=0.90, RMSE=11.08 t/hm2). The best-fitting semi-variogram function model of forest biomass was exponential model in Shangri-La(C0=0.12, C0 + C=0.87, A0=
10200 m). Compared with ordinary Kriging, the spatial distribution of AGB obtained by the Sequential Gaussian Conditional Simulation had better consistency(r=0.59**, d=0.70). The spatial differentiation of AGB could be explained by topographic factors. In terms of explanatory power, elevation was the largest, slope was the second, and slope was the least. The inversion accuracy of forest AGB based on spaceborne LiDAR ICESat–2/ATLAS data was high(Pp=81.43%), which provided a reliable data source for geostatistical analysis. Therefore, the method based on spaceborne LiDAR and geostatistics can greatly analyze the spatial heterogeneity of forest AGB. -
森林地上生物量（AGB）是森林生态系统生产力的重要指标，也是碳循环的重要参数[1]，探究森林AGB的空间分异特征，对于研究森林生态系统碳循环具有重要意义[2-3]。森林生物量由地上生物量和地下生物量组成，但地下生物量的调查难以开展，所以大多研究都是围绕森林地上生物量这一对象开展[4]。现阶段，森林生物量估测常用方法包括样地调查法和基于遥感的方法[5]。与传统人工调查相比，遥感技术具有明显优势并成为在不同空间尺度下估算森林生物量的主要手段[5]。然而，光学遥感易受到天气和饱和特性的影响[6−8]。微波遥感易受到地形影响和存在饱和问题[9]。星载激光雷达数据具有高时空分辨率、低成本、大监测面积等特点，在大尺度AGB反演中优势明显[10−12]。Ice Cloud, and land Elevation Satelite–2（ICESat–2）所搭载的先进地形激光测高系统（ATLAS）采用了多波束、微脉冲、光子计数激光雷达技术[13−14]。其数据已被众多学者用于森林生物量的估测研究[15−19]。
空间异质性是90年代生态学研究的一个极为重要的理论问题[20]，作为地理现象的共同属性，它是指各种地理空间属性在一定地理区域内的不均匀分布[21−22]。林木生长与立地条件存在着密切的关系，环境条件差异影响着林木的生长及其分布[23−24]。森林生物量的遥感定量过程经常忽视森林生物量的空间关系。传统方法不能及时、准确地获取森林生物量作为空间属性，难以有效反映其空间异质性和空间自相关[25]。
香格里拉位于“三江并流”核心区，是全球公认的生物多样性热点地区。探索该区域的森林生物量空间分布特征具有重要意义。本研究以ICESat–2/ATLAS为主要数据源，结合54个实测样地数据，采用机器学习算法建模并预测森林生物量，使用地统计学对足迹生物量进行空间分析，以期为大尺度的森林生物量空间变异分析提供参考。
1. 研究区概况
研究区为云南省迪庆藏族自治州香格里拉市，位于北纬26°52′～28°52′，东经99°20′～100°19′。全市属山地寒温带季风气候类型，干湿季节明显，而四季变化不显著，夏秋季节降雨量较大，而冬春季节则相对干燥。区域地处高海拔（平均海拔约
3459 m）低纬度地带，海拔相对高差大，其地形主要由北向南延伸的山脉和北部较高、南部较低的地势构成。特殊的地理位置造成了明显的土壤垂直分布特征，海拔由低到高依次出现燥红土、红壤、棕壤、暗棕壤、棕色针叶林土、高山草甸土、亚高山草甸土和高山寒漠土等土壤类型。森林植被方面，该市的森林覆盖率为74.99%，主要的优势树种有云杉（Picea asperata）、高山栎（Quercus semecarpifolia）、冷杉（Abies fabri）、高山松（Pinus densata）等[26]。2. 材料与方法
2.1 数据来源
2.1.1 样地数据
54块基于ATLAS足迹的圆形样地（r=8.5 m, S=226.9 m2）于2021年11月在研究区建立。使用差分定位设备（Seeker SR3 Pro）确保样地中心点坐标与足迹中心点坐标匹配，匹配误差限制在采样误差（0.02 m）以内。生物量模型中，选取了相近树种或其他地区的生物量模型替代难以收集的部分树种的生物量公式，整理的研究区主要树种的生物量公式可见表1[27]。在所设立的样圆范围内记录生物量等主要信息，描述性统计可见表2。
表 1 研究区主要树种的立木生物量模型Table 1. Stand biomass model of main tree species in the study area树种 生物量模型 云南松（Pinus yunnanensis） MA= 0.070231 D2.10392 H0.41120 高山松（Pinus densata） MA= 0.0730 D2.3560 H0.1090 云杉（Picea asperata） MA= 0.09152 D2.2106 H0.25663 (D≥5)冷杉（Abies fabri） MA= 0.06127 D2.05753 H0.50839 (D≥5)华山松（Pinus armandii） MA= 0.009512 (D2H)1.138665 其他阔叶树 MA= 0.07806 D2.06321 H0.57393 (D≥5)注：MA代表立木生物量（kg）；D代表平均胸径（cm）；H表示平均树高（m）。其他阔叶树包括杨树、枫树等。 表 2 实测生物量的描述性统计Table 2. Descriptive statistics of measured biomassAGB统计项 值 数量 54 最小值 4.92 最大值 126.24 平均值 59.66 标准差 37.33 偏度 0.07 峰度 −0.99 2.1.2 ICESat–2/ATLAS数据的介绍及处理
ICESat–2数据产品目前分为4级共22种，以ATL00～22命名，本研究以ATL03、ATL08为主要数据源。ATL03作为生成其他产品的基础数据，记录光子事件发生的时间、椭球高度和经纬度等[28]。ATL08是对ATL03数据进行去噪分类后生成的。ATL08产品提供沿轨道方向的地形高度和森林冠层高度及其相关参数信息。本研究为使ATLAS数据在整个研究区内的均匀覆盖，使用了2020年6月1日后1 a时间段内的ATL03、ATL08数据产品。光子点云去噪方面，本研究采用的去噪算法是由基于密度差的空间聚类噪声算法（DDBSCAN）[29]和基于K−近邻的去噪算法（KNNB）[30]组成。除去噪外，在森林区域的研究中，需将光子点云分为地面光子点和冠层光子顶部点。本研究为克服不规则三角网络加密（PTD）法在复杂地形下的光子错分情况，使用改进的PTD方法，并将地面点定为离初始 TIN 最远点的下方高程最低点[31]。
研究区林地范围内的光斑采样足迹共计
58260 个。为了使足迹分布均匀、随机，对58260 个林地光斑足迹进行系统抽样（抽样间隔：50），最终抽样所得的1165 个星载激光雷达采样足迹被作为空间分析的预设样本。2.1.3 DEM 数据
数字高程模型（DEM）的空间分辨率为12.5 m，该数据是由ALOS（Advanced Land Observing Satellite-1）卫星的 PALSAR 传感器获取。DEM数据被重采样到17 m空间分辨率，并提取研究区的坡向、坡度和高程。为便于地形因子参与分析，需按照研究区地形特征和林业上的划分标准对地形因子进行赋值。海拔赋值如下：1对应低于
2200 m的海拔范围，2对应2200 ～2600 m的海拔范围，3对应2600 ～3000 m的海拔范围，4对应3000 ～3400 m的海拔范围，5对应3400 ～3800 m的海拔范围，6对应3800 ～4200 m的海拔范围，7对应4200 ～4600 m的海拔范围，8对应大于4600 m的海拔范围。坡度赋值如下：1对应平坡（ < 5°），2对应缓坡（5°～15°），3对应斜坡（15°～25°），4对应陡坡（25°～35°），5对应急坡（35°～45°），6对应险坡（ > 45°）。坡向赋值如下：无坡向赋值为1，北坡为2，东北坡为3，东坡为4，东南坡为5，南坡为6，西南坡为7，西坡为8，西北坡为9。2.1.4 研究区森林资源规划设计调查数据
森林资源规划设计调查的目的是查清森林资源现状，其成果是对森林资源进行科学经营管理的重要依据[32]。本研究中，该数据被用于提取ATLAS的林地光斑足迹和得到森林生物量空间分布的掩膜结果。
2.2 机器学习算法
随机森林（RF）是由Breiman[33]提出的，作为取一种较新的集成学习技术，其基本思想是将弱分类器组合成强分类器，弱分类器即使获得了错误的预测，其他弱分类器也可以纠正错误。RF开发之初作为分类和回归树（CART）的扩展，它可以生成聚合预测器来提高预测模型的性能[34]。
梯度提升回归树（GBRT）是一种集成学习方法，通过顺序聚合一组弱CART回归树子模型来构建强学习模型[35−36]。GBRT的关键概念是每个新的回归树子模型都是在残差缩减的梯度方向上建立的，以减少与先前模型的残差[37]。
构建的机器学习模型采用十折交叉验证进行的精度评价[38]，该方法的原理是：将数据集分成10份，每次轮流从中选1份作为测试集，剩下9份作为训练集，重复10次，并将10次验证结果平均后评价模型。均方根误差（RMSE）、决定系数（R2）、总体估测精度（Pp）被用于评估回归模型的拟合效果，相关计算方法可参考宋涵玥等[27]的研究。
2.3 地统计分析方法
2.3.1 空间自相关分析
Moran's I是空间自相关的常用指标。本研究用Moran's I进行空间自相关分析，Moran's I的取值区间为 [−1, 1]。当Moran's I的取值小于0，代表空间对象存在空间负相关，森林生物量呈现破碎化的空间分布；当Moran's I等于0或接近于0，空间自相关不存在；当Moran's I的取值大于0，说明存在空间自相关，空间对象呈现聚集分布。Moran's I计算方法如下：
$$ I=\frac{n\underset{i=1}{\displaystyle\sum ^{n}}\underset{j\ne 1}{\displaystyle\sum ^{n}}{W}_{ij}\left({X}_{i}- \overline{X}\right)\left({X}_{j}- \overline{X}\right)}{\underset{i=1}{\displaystyle\sum ^{n}}\underset{j\ne 1}{\displaystyle\sum ^{n}}{W}_{ij}\underset{i=1}{\displaystyle\sum ^{n}}{\left({X}_{i}- \overline{X}\right)}^{2}} $$ (4) 式中：n为观测值个数；
$ {X}_{i} $ 、$ {X}_{j} $ 是变量在空间位置i和j处的值；$ \overline{X} $ 为变量$ X $ 的平均值；$ {W}_{ij} $ 为空间权重函数。局部空间自相关能够阐明研究对象与其相邻单元之间存在的空间相关性，其计算方法参考杨晴青等[39]的研究。
2.3.2 空间异质性分析
表 3 半变异函数参数Table 3. Semi-variance function parameters参数 符号 描述 块金值 C0 反映了区域化变量内部随机性的可能程度 基台值 C0 + C 区域化变量总体特征的体现，用来衡量空间异质性程度 变程 A0 表示研究变量之间存在空间自相关性的平均最大距离 块金比 C0/(C0 + C) 反映随机部分引起的空间异质性占系统总变异的比例 半变异函数计算公式如下：
$$ r\left(h\right)=\dfrac{1}{2N\left(h\right)}\underset{i=1}{\displaystyle\sum ^{N\left(h\right)}[}{Z\left({x}_{i}\right)-Z\left({x}_{i} + h\right)]}^{2} $$ (5) 式中：
$ r\left(h\right) $ 代表生物量半变异函数；$ N\left(h\right) $ 为某一方向上间距为h的样本点总对数；$ Z\left({x}_{i}\right) $ 为空间位置$ {x}_{i} $ 处的生物量实测值；$ Z\left({x}_{i} + h\right) $ 为与点xi偏离距离为h处的生物量值。2.4 基于地统计学的空间插值
克里格（Kriging）是一种估计观测样点间内插值的地统计学方法，它基于区域化变量理论；当获得了变量的半变异函数的最佳拟合模型后，就可以利用样点观测值对研究区域上未取样点的区域化变量值进行估测[41]，进而可得到森林生物量的空间分布格局。
序贯高斯条件模拟（SGCS）是一种基于区域化变量理论和空间自相关生成研究变量的空间显式估计方法，基于已知数据构造高斯函数，将区域化随机变量
$ Z\left(x\right) $ 的每个值视为高斯函数$ F\left(x\right) $ 的随机实现。在每个模拟位置xm，使用已知数据$ Z\left({x}_{i}\right) $ （$ i $ =1, 2, ···, n）和先前模拟的值$ Z\left({x}_{j}\right) $ （j=1, 2, ···, m−1）作为条件累积条件概率密度函数推导$ F\left(x\right) $ [42]，然后利用累积条件概率分布通过序贯随机模拟生成空间预测。与Kriging相比，SGCS弱化了由Kriging产生的平滑效应，避免由单次拟合造成的不确定性[43]，其主要实现步骤可参考蒋威等[44]的研究。采用皮尔逊相关系数和一致性指数2个指标评价插值结果精度。皮尔逊相关系数和一致性指数的计算公式如下：
$$ r=\frac{\underset{i=1}{\displaystyle\sum ^{n}}\left[Z\left({x}_{i}\right)- \overline{M}\right]·\left[{Z}'\left({x}_{i}\right)-{ \overline{M}}'\right]}{\underset{i=1}{\displaystyle\sum ^{n}}{\left[Z\left({x}_{i}\right)- \overline{M}\right]}^{2}\sqrt{\underset{i=1}{\displaystyle\sum ^{n}}{\left[{Z}'\left({x}_{i}\right)-{ \overline{M}}'\right]}^{2}}} $$ (6) $$ d=1-\frac{{\left[{Z}'\left({x}_{i}\right)-Z\left({x}_{i}\right)\right]}^{2}}{\underset{i=1}{\displaystyle\sum ^{n}}{\left[\left|{Z}'\left({x}_{i}\right)- \overline{M}\right| + \left|Z\left({x}_{i}\right)- \overline{M}\right|\right]}^{2}} $$ (7) 式中：
$ n $ 为验证样本数，$ Z\left({x}_{i}\right) $ 和$ {Z}'\left({x}_{i}\right) $ 分别是位于$ {x}_{i} $ 处的森林生物量的预测值与实测值，$ \overline{M} $ 和$ \overline{M}'$ 分别表示实测值和预测值的平均值。2.5 基于地理探测器的森林生物量空间分异影响分析
地理探测器作为探测空间分异性，以及揭示其背后驱动力的统计学方法[45]，在应用时没有过多的假设条件，克服了统计方法处理变量时的局限，被广泛应用在多种因素的影响机理研究当中[46−47]。地理探测器的因子探测可以识别影响因子，交互探测可以解释影响因子对因变量的交互作用，是研究复杂地理因素驱动作用机理的有效工具[48-49]。
3. 结果与分析
3.1 机器学习模型
构建预测模型前，采用随机森林特征重要性评价法对ATLAS参数计算特征重要性。贡献度前10的参数被作为自变量（自变量及其详细信息见表4），实测森林生物量作为因变量建立机器学习模型。
表 4 ATLAS参数及说明Table 4. ATLAS parameters and their details特征项目 贡献度/% 描述 landsat_perc 19.99 陆地卫星冠层百分比 latitude 12.26 纬度 n_toc_photons 11.83 冠层顶部光子的数量 photon_rate_can 6.04 冠层光子比 segment_id_beg 5.75 第一光子在陆地部分的地理位置段号 h_dif_canopy 5.46 冠层高和中位数冠层高之差 asr 4.50 表面反射率 toc_roughness 3.52 在段内被归类为冠层顶部的所有光子相对高度的标准偏差 beam_azimuth 3.33 参考光子在本地ENU坐标系中的单位指向矢量的方位 h_min_canopy 3.29 段内个体相对冠层高度的平均值 使用默认参数设置的随机森林算法构建森林生物量估测模型，应用GBRT算法构建森林生物量估测模型。由图1a可知，随机森林模型R2=0.86、RMSE=15.52 t/hm2、Pp=73.99%。由图1b可知，GBRT模型R2=0.90、RMSE=11.08 t/hm2、Pp=81.43%。基于实测值与模型预测值绘制的散点图（图1）显示了两种机器学习模型预测值与实测值的比较。使用1∶1线和拟合线来直观地描述估计值和测量值之间的关系，模型的主要误差与高估森林生物量低值和低估森林生物量高值有关。
3.2 光斑足迹的森林地上生物量
预测精度较高的GBRT模型被用于预测研究区内
1165 个光斑内的森林生物量（图2），1165 个光斑对应的地形因子描述统计见表5。森林AGB的分布特征呈现正态分布，满足半变异函数结构分析的要求。表 5 光斑生物量及其地形因子的描述统计Table 5. Descriptive statistics of biomass and topographic factors of light spots项目 AGB/
(t∙hm−2)海拔/m 坡度/(°) 坡向/(°) 样本个数 1165 1165 1165 1165 最大值 119.89 4646 57.34 357.21 最小值 16.73 1821 2.43 0 平均值 60.20 3472 26.72 184.14 标准差 15.65 564 11.08 102.04 变异系数 0.26 0.16 0.42 0.55 分布类型 正态分布 正态分布 正态分布 平方根正态 3.3 森林AGB的空间自相关性
基于森林AGB预测值的空间自相关分析显示：Z得分=5.79和P < 0.01，说明Moran's I通过检验，置信度为99%。研究区森林AGB的Moran's I为正值（Moran's I=0.31），说明研究区域内森林AGB具有空间正相关性，属于空间集聚分布。
局部Moran's I散点图（图3）显示了每个森林生物量反演位置的空间自相关性。由图3可知，研究区中部和北部的生物量表现出显著的“高–高”聚类，而研究区西部和东部则表现出“高–低”型的分布特征。此外，“低–高”型分布主要集中在研究区的中部，“低–低”型的聚集分布模式主要集中在研究区的东部。
![]()
图 3 2021年香格里拉市森林生物量的局部Moran's I散点图Figure 3. Local Moran's I scatter plot of forest biomass in Shangri-La City in 20213.4 森林AGB的空间异质性
半变异函数主要参数可见表6。根据残差平方和（RSS）最小，决定系数最大的原则，研究区森林生物量的空间变异最适合用指数模型来描述。由表6可知，最优拟合的C0/(C0 + C)小于25%，表明研究区森林生物量异质性变化主要是由空间自相关引起。人类活动等外界产生的随机因素对研究区森林AGB空间变异的影响较小，当地的森林管理与保护较优。
表 6 研究区森林地上生物量半变异函数主要参数表Table 6. Main parameters of semi-variance function of forest aboveground biomass in the study area模型 C0 C0 + C A0 C0/(C0 + C)/% C/(C0 + C)/% R2 RSS 线性模型 0.79 0.98 10276 80.61 19.39 0.71 1.38 × 10−6 球状模型 0.41 0.87 7500 47.13 52.87 0.49 2.41 × 10−6 指数模型 0.12 0.87 10200 13.79 86.21 0.62 1.84 × 10−6 高斯模型 0.13 0.87 6408 14.94 85.06 0.51 2.37 × 10−6 根据表6中的相应参数可以获得指数模型表达式，指数模型表达式如下：
$$ \gamma \left(h\right)=\left\{\begin{array}{ccc}0& h=0& \\ 0.001\;2 + 0.007\;5\left(1-{\mathrm{e}}^{-\frac{h}{10\;200}}\right)& h > 0& \end{array}\right. $$ (8) 3.5 森林生物量的空间连续制图
基于ICESat–2/ATLAS光斑足迹的AGB反演结果，使用克里格和序贯高斯条件模拟（模拟次数：50次）两种方法绘制研究区森林生物量的空间分布图（图4）。
1165 个反演的生物量足迹被随机拆分为两部分（80%绘制空间分布，20%检验插值精度）。由图4可知，基于两种方法得到的空间分布基本一致。森林地上生物量高值区位于香格里拉北部、中部和南部；低值区主要分布在香格里拉西南部和东南部。由表7可知，相较于Kriging法，序贯高斯条件模拟法得到的空间分布图与验证集的相关系数相对提高5.36%，一致性系数相对提高2.94%，具有较强的一致性（r=0.59**，d=0.70）。选择精度较高的图4b得到研究区森林生物量，其估计值为8.24 × 107 t。![]()
图 4 研究区森林生物量空间分布图Figure 4. Spatial distribution map of forest biomass in the study area表 7 森林地上生物量插值方法的精度评价Table 7. Accuracy evaluation of interpolation methods for forest aboveground biomass插值方法 r d Kriging 0.56** 0.68 序贯高斯条件模拟 0.59** 0.70 注：**表示相关性极显著（P < 0.01）。 3.6 基于地理探测器的地形因子对森林生物量空间分异影响
基于地理探测器，对研究区森林生物量空间分异的地形因子进行分析。因子探测器的结果可见表8。因子探测器的q值由大到小的顺序：海拔（0.25） > 坡向（0.03） > 坡度（0.02）。可见，相较于其他地形因子，海拔对森林生物量空间分异的解释力最强。
表 8 地形因子解释力探测结果Table 8. Detection results of explanatory power of terrain factors统计值 海拔 坡度 坡向 q 0.25 0.02 0.03 p < 0.01 0.03 0.01 交互探测器可以评估地形因子之间共同作用是否会增加或减少对森林生物量的解释力[50]。由表9可知，森林生物量的空间分异是地形因子综合作用的结果。高程∩坡向两因子交互时解释力最大；其次是高程∩坡度，即本研究所采用的地形因子组合会增强因子对森林生物量间分异的解释力。此外，高程∩坡度、高程∩坡向、坡度∩坡向在影响森林生物量方面的作用是非线性增强的。
表 9 地形因子相互作用的探测结果Table 9. Detection results of terrain factor interactionC=A∩B A + B 比较 解释 海拔∩坡度=0.308 海拔(0.251) + 坡度(0.022)=0.273 C > A + B 海拔↑坡度 海拔∩坡向=0.309 海拔(0.251) + 坡向(0.030)=0.281 C > A + B 海拔↑坡向 坡度∩坡向=0.119 坡度(0.022) + 坡向(0.030)=0.052 C > A + B 坡度↑坡向 注：A↑B表示A和B非线性增强。 4. 结论与讨论
4.1 讨论
基于GBRT的机器学习算法构建的光斑足迹内的森林生物量估测模型精度较高，实现了使用较低成本实现较为密集的空间属性覆盖。然而，模型都存在高值低估的问题。基于不同森林类型的建模可以减少模型对训练样本的依赖并提高建模效果[51−52]。由于缺乏样本地块数据，无法区分森林类型或生物量水平以进行建模。将来，通过收集足够的样地数据，可以减少建模过程中的不确定性。此外，自然地理、生物气候等对估计森林参数至关重要[53−54]。因此，在后续研究中，应将遥感数据与森林生理过程模型耦合，提高模型的泛化性和准确性。
本研究中，研究区森林生物量的空间异质性主要是由结构性因素及水热条件等自然因素作用的结果。分析地形因子对森林参数空间异质性的影响可以更好地理解气候–林木之间的相互作用的机制，这是当前研究中经常被忽视的。然而，本研究缺乏地形与气候因子的关系以及两者对生物量空间分布的共同影响，未来需要进一步分析。
在已有的研究中，王继雄等[55]利用2016年森林资源规划设计调查数据，采用香格里拉不同类型优势树种的生物量模型估算结果为7.35 × 107 t。由于本研究结果来自2016年森林资源规划设计调查，故将其作为本研究估算精度的参考标准。本研究的森林生物量的估计值为8.24 × 107 t，与参考标准相同数量级，估计的绝对精度为77.6%。然而，本研究为了提升足迹分布的随机性，舍弃了大量光斑信息，插值的准确度受到一定影响。为了获取更多并保持随机性的空间观测足迹，今后的研究可引进另一个星载激光雷达GEDI数据。
4.2 结论
本研究基于星载激光雷达ICESat–2/ATLAS与地统计学相结合的方法，能够很好地反映云南省香格里拉市森林生物量的空间分布格局以及空间变异规律。研究表明：采用梯度提升树算法建立的森林生物量估测精度较高，模型决定系数为0.91，均方根误差为11.08 t/hm2；研究区域内森林生物量具有空间正相关性，属于空间集聚分布；森林AGB的空间变异特征最适合用指数模型描述（R2=0.62，RSS=1.84 × 10−6）；基于光斑足迹生物量的序贯高斯条件模拟能够较好地反映研究区森林生物量的空间分布（r=0.59**，d=0.70）。此外，地形因子会对森林AGB的空间分异造成影响，在解释力方面，海拔 > 坡向 > 坡度。
-
![]()
图 3 2021年香格里拉市森林生物量的局部Moran's I散点图
Figure 3. Local Moran's I scatter plot of forest biomass in Shangri-La City in 2021
![]()
图 4 研究区森林生物量空间分布图
Figure 4. Spatial distribution map of forest biomass in the study area
表 1 研究区主要树种的立木生物量模型
Table 1 Stand biomass model of main tree species in the study area
树种 生物量模型 云南松（Pinus yunnanensis） MA= 0.070231 D2.10392 H0.41120 高山松（Pinus densata） MA= 0.0730 D2.3560 H0.1090 云杉（Picea asperata） MA= 0.09152 D2.2106 H0.25663 (D≥5)冷杉（Abies fabri） MA= 0.06127 D2.05753 H0.50839 (D≥5)华山松（Pinus armandii） MA= 0.009512 (D2H)1.138665 其他阔叶树 MA= 0.07806 D2.06321 H0.57393 (D≥5)注：MA代表立木生物量（kg）；D代表平均胸径（cm）；H表示平均树高（m）。其他阔叶树包括杨树、枫树等。 
下载: 导出CSV
表 2 实测生物量的描述性统计
Table 2 Descriptive statistics of measured biomass
AGB统计项 值 数量 54 最小值 4.92 最大值 126.24 平均值 59.66 标准差 37.33 偏度 0.07 峰度 −0.99
下载: 导出CSV
表 3 半变异函数参数
Table 3 Semi-variance function parameters
参数 符号 描述 块金值 C0 反映了区域化变量内部随机性的可能程度 基台值 C0 + C 区域化变量总体特征的体现，用来衡量空间异质性程度 变程 A0 表示研究变量之间存在空间自相关性的平均最大距离 块金比 C0/(C0 + C) 反映随机部分引起的空间异质性占系统总变异的比例
下载: 导出CSV
表 4 ATLAS参数及说明
Table 4 ATLAS parameters and their details
特征项目 贡献度/% 描述 landsat_perc 19.99 陆地卫星冠层百分比 latitude 12.26 纬度 n_toc_photons 11.83 冠层顶部光子的数量 photon_rate_can 6.04 冠层光子比 segment_id_beg 5.75 第一光子在陆地部分的地理位置段号 h_dif_canopy 5.46 冠层高和中位数冠层高之差 asr 4.50 表面反射率 toc_roughness 3.52 在段内被归类为冠层顶部的所有光子相对高度的标准偏差 beam_azimuth 3.33 参考光子在本地ENU坐标系中的单位指向矢量的方位 h_min_canopy 3.29 段内个体相对冠层高度的平均值
下载: 导出CSV
表 5 光斑生物量及其地形因子的描述统计
Table 5 Descriptive statistics of biomass and topographic factors of light spots
项目 AGB/
(t∙hm−2)海拔/m 坡度/(°) 坡向/(°) 样本个数 1165 1165 1165 1165 最大值 119.89 4646 57.34 357.21 最小值 16.73 1821 2.43 0 平均值 60.20 3472 26.72 184.14 标准差 15.65 564 11.08 102.04 变异系数 0.26 0.16 0.42 0.55 分布类型 正态分布 正态分布 正态分布 平方根正态
下载: 导出CSV
表 6 研究区森林地上生物量半变异函数主要参数表
Table 6 Main parameters of semi-variance function of forest aboveground biomass in the study area
模型 C0 C0 + C A0 C0/(C0 + C)/% C/(C0 + C)/% R2 RSS 线性模型 0.79 0.98 10276 80.61 19.39 0.71 1.38 × 10−6 球状模型 0.41 0.87 7500 47.13 52.87 0.49 2.41 × 10−6 指数模型 0.12 0.87 10200 13.79 86.21 0.62 1.84 × 10−6 高斯模型 0.13 0.87 6408 14.94 85.06 0.51 2.37 × 10−6
下载: 导出CSV
表 7 森林地上生物量插值方法的精度评价
Table 7 Accuracy evaluation of interpolation methods for forest aboveground biomass
插值方法 r d Kriging 0.56** 0.68 序贯高斯条件模拟 0.59** 0.70 注：**表示相关性极显著（P < 0.01）。
下载: 导出CSV
表 8 地形因子解释力探测结果
Table 8 Detection results of explanatory power of terrain factors
统计值 海拔 坡度 坡向 q 0.25 0.02 0.03 p < 0.01 0.03 0.01
下载: 导出CSV
表 9 地形因子相互作用的探测结果
Table 9 Detection results of terrain factor interaction
C=A∩B A + B 比较 解释 海拔∩坡度=0.308 海拔(0.251) + 坡度(0.022)=0.273 C > A + B 海拔↑坡度 海拔∩坡向=0.309 海拔(0.251) + 坡向(0.030)=0.281 C > A + B 海拔↑坡向 坡度∩坡向=0.119 坡度(0.022) + 坡向(0.030)=0.052 C > A + B 坡度↑坡向 注：A↑B表示A和B非线性增强。
下载: 导出CSV
-
[1] 廖易, 张加龙, 鲍瑞, 等. 引入地形因子的高山松地上生物量动态估测 [J]. 生态学杂志, 2023, 42(5): 1243−1252. [2] 王书贤, 张加龙, 鲍瑞, 等. 香格里拉高山松地上生物量及遥感因子的时空变化特征研究 [J]. 西北林学院学报, 2023, 38(5): 43−48. DOI: 10.3969/j.issn.1001-7461.2023.05.06 [3] 郝晴, 黄昌. 森林地上生物量遥感估算研究综述 [J]. 植物生态学报, 2023, 47(10): 1356−1374. DOI: 10.17521/cjpe.2023.0008 [4] Lu D S. The potential and challenge of remote sensing-based biomass estimation [J]. International Journal of Remote Sensing, 2006, 27(7): 1297−1328. DOI: 10.1080/01431160500486732
[5] 毛学刚, 王静文, 范文义. 基于遥感与地统计的森林生物量时空变异分析 [J]. 北京林业大学学报, 2016, 38(2): 10−19. [6] Peduzzi A, Wynne R H, Fox T R, et al. Estimating leaf area index in intensively managed pine plantations using airborne laser scanner data [J]. Forest Ecology and Management, 2012, 270: 54−65. DOI: 10.1016/j.foreco.2011.12.048
[7] Wang L, Xu X, Yu Y, et al. SAR-to-optical image translation using supervised cycle-consistent adversarial networks [J]. IEEE Access, 2019, 7: 129136−129149. DOI: 10.1109/ACCESS.2019.2939649
[8] Chopping M, Moisen G G, Su L H, et al. Large area mapping of southwestern forest crown cover, canopy height, and biomass using the NASA Multiangle Imaging Spectro-Radiometer [J]. Remote Sensing of Environment, 2008, 112(5): 2051−2063. DOI: 10.1016/j.rse.2007.07.024
[9] Vatandaşlar C, Abdikan S. Carbon stock estimation by dual-polarized synthetic aperture radar (SAR) and forest inventory data in a Mediterranean forest landscape [J]. Journal of Forestry Research, 2022, 33(3): 827−838. DOI: 10.1007/s11676-021-01363-3
[10] Disney M. Terrestrial LiDAR: a three-dimensional revolution in how we look at trees [J]. New Phytologist, 2019, 222(4): 1736−1741. DOI: 10.1111/nph.15517
[11] Pitkänen T P, Raumonen P, Kangas A. Measuring stem diameters with TLS in boreal forests by complementary fitting procedure [J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2019, 147: 294−306. DOI: 10.1016/j.isprsjprs.2018.11.027
[12] Wang Y, Ni W J, Sun G Q, et al. Slope-adaptive waveform metrics of large footprint lidar for estimation of forest aboveground biomass [J]. Remote Sensing of Environment, 2019, 224: 386−400. DOI: 10.1016/j.rse.2019.02.017
[13] 朱笑笑, 王成, 习晓环, 等. ICESat–2星载光子计数激光雷达数据处理与应用研究进展 [J]. 红外与激光工程, 2020, 49(11): 76−85. [14] Magruder L A, Brunt K M. Performance analysis of airborne photon- counting lidar data in preparation for the ICESat–2 mission [J]. IEEE Transactions on Geoscience and Remote Sensing, 2018, 56(5): 2911−2918. DOI: 10.1109/TGRS.2017.2786659
[15] Song H Y, Xi L, Shu Q T, et al. Estimate forest aboveground biomass of mountain by ICESat–2/ATLAS data interacting cokriging [J]. Forests, 2022, 14(1): 13. DOI: 10.3390/f14010013
[16] Narine L L, Popescu S C, Malambo L. Using ICESat–2 to estimate and map forest aboveground biomass: a first example [J]. Remote Sensing, 2020, 12(11): 1824. DOI: 10.3390/rs12111824
[17] Guerra-Hernández J, Narine L L, Pascual A, et al. Aboveground biomass mapping by integrating ICESat–2, SENTINEL-1, SENTINEL-2, ALOS2/PALSAR2, and topographic information in Mediterranean forests [J]. GIScience & Remote Sensing, 2022, 59(1): 1509−1533.
[18] Narine L L, Popescu S C, Malambo L. Synergy of ICESat–2 and landsat for mapping forest aboveground biomass with deep learning [J]. Remote Sensing, 2019, 11(12): 1503. DOI: 10.3390/rs11121503
[19] Rodda S R, Nidamanuri R R, Fararoda R, et al. Evaluation of Height Metrics and Above-Ground Biomass Density from GEDI and ICESat–2 Over Indian Tropical Dry Forests using Airborne LiDAR Data[J]. Journal of the Indian Society of Remote Sensing, 2023: 1-16.
[20] Kareiva P. Special feature: space: the final frontier for ecological theory [J]. Ecology, 1994, 75(1): 1. DOI: 10.2307/1939376
[21] Fischer M M, Getis A. Handbook of applied spatial analysis: software tools, methods and applications [M]. Berlin Heidelberg: Springer-Verlag, 2010.
[22] Wang J F, Zhang T L, Fu B J. A measure of spatial stratified heterogeneity [J]. Ecological Indicators, 2016, 67: 250−256. DOI: 10.1016/j.ecolind.2016.02.052
[23] 孔宁宁, 曾辉, 李书娟. 四川卧龙自然保护区植被的地形分异格局研究 [J]. 北京大学学报(自然科学版), 2002(4): 543−549. [24] Hu Z W. Impacts of topography on the spatial pattern of the age of forest community [J]. Chinese Journal of Plant Ecology, 2007, 31(5): 814−824. DOI: 10.17521/cjpe.2007.0103
[25] 陈先刚, 张一平, 张小全, 等. 过去50年中国竹林碳储量变化 [J]. 生态学报, 2008, 28(11): 5218−5227. DOI: 10.3321/j.issn:1000-0933.2008.11.003 [26] Xu D F, Zhang J L, Bao R, et al. Temporal and spatial variation of aboveground biomass of Pinus densata and its drivers in shangri-La, CHINA [J]. International Journal of Environmental Research and Public Health, 2021, 19(1): 400. DOI: 10.3390/ijerph19010400
[27] 宋涵玥, 舒清态, 席磊, 等. 基于星载ICESat–2/ATLAS数据的森林地上生物量估测 [J]. 农业工程学报, 2022, 38(10): 191−199. DOI: 10.11975/j.issn.1002-6819.2022.10.023 [28] 俄相颖, 戴光耀, 吴松华. ICESat–2 ATL03数据预处理及校正方法[J]. 红外与激光工程, 2021, 50(6): 20211032. [29] Zhang J S, Kerekes J. An adaptive density-based model for extracting surface returns from photon-counting laser altimeter data [J]. IEEE Geoscience and Remote Sensing Letters, 2015, 12(4): 726−730. DOI: 10.1109/LGRS.2014.2360367
[30] Xia S B, Wang C, Xi X H, et al. Point cloud filtering and tree height estimation using airborne experiment data of ICESat–2 [J]. National Remote Sensing Bulletin, 2014, 18(6): 1199−1207. DOI: 10.11834/jrs.20144029
[31] Nie S, Wang C, Dong P L, et al. A revised progressive TIN densification for filtering airborne LiDAR data [J]. Measurement, 2017, 104: 70−77. DOI: 10.1016/j.measurement.2017.03.007
[32] 杨毅, 周汝良, 李靖, 等. 基于GIS与多扩展接口的林地资源数据库的研究与设计 [J]. 林业资源管理, 2003(4): 56−59, 69. DOI: 10.3969/j.issn.1002-6622.2003.04.015 [33] Breiman L. Random forests [J]. Machine Learning, 2001, 45(1): 5−32. DOI: 10.1023/A:1010933404324
[34] Breiman L. Bagging predictors [J]. Machine Learning, 1996, 24(2): 123−140.
[35] Opitz D, Maclin R. Popular ensemble methods: an empirical study [J]. Journal of Artificial Intelligence Research, 1999, 11: 169−198. DOI: 10.1613/jair.614
[36] Friedman J H. Greedy function approximation: a gradient boosting machine [J]. The Annals of Statistics, 2001, 29(5): 1189−1232. DOI: 10.1214/aos/1013203450
[37] Liu B, Gao L, Li B A, et al. Nonparametric machine learning for mapping forest cover and exploring influential factors [J]. Landscape Ecology, 2020, 35(7): 1683−1699. DOI: 10.1007/s10980-020-01046-0
[38] 范永东. 模型选择中的交叉验证方法综述[D]. 太原: 山西大学, 2013. [39] 杨晴青, 刘倩, 尹莎, 等. 秦巴山区乡村交通环境脆弱性及影响因素: 以陕西省洛南县为例 [J]. 地理学报, 2019, 74(6): 1236−1251. DOI: 10.11821/dlxb201906012 [40] 王政权, 王庆成, 李哈滨. 红松老龄林主要树种的空间异质性特征与比较的定量研究 [J]. 植物生态学报, 2000, 24(6): 718−723. DOI: 10.3321/j.issn:1005-264X.2000.06.013 [41] Bivand R S, Pebesma E, Gómez-Rubio V. Applied Spatial Data Analysis with R[M]. New York: Springer New York, 2013.
[42] 赵彦锋, 化全县, 陈杰. Kriging插值和序贯高斯条件模拟的原理比较及在土壤空间变异研究中的案例分析 [J]. 土壤学报, 2011, 48(4): 856−862. [43] Luo S L, Xu L, Yu J G, et al. Sampling estimation and optimization of typical forest biomass based on sequential Gaussian conditional simulation [J]. Forests, 2023, 14(9): 1792. DOI: 10.3390/f14091792
[44] 蒋威, 郜允兵, 刘玉, 等. 基于耕地图斑数据的土壤有机质估测方法研究: 以大兴区南部耕地为例 [J]. 中国农业大学学报, 2017, 22(11): 75−82. DOI: 10.11841/j.issn.1007-4333.2017.11.08 [45] 王劲峰, 徐成东. 地理探测器: 原理与展望 [J]. 地理学报, 2017, 72(1): 116−134. DOI: 10.11821/dlxb201701010 [46] 刘彦随, 杨忍. 中国县域城镇化的空间特征与形成机理 [J]. 地理学报, 2012, 67(8): 1011−1020. DOI: 10.11821/xb201208001 [47] 湛东升, 张文忠, 余建辉, 等. 基于地理探测器的北京市居民宜居满意度影响机理 [J]. 地理科学进展, 2015, 34(8): 966−975. [48] Wang J F, Hu Y. Environmental health risk detection with GeogDetector [J]. Environmental Modelling & Software, 2012, 33: 114−115.
[49] 吕晨, 蓝修婷, 孙威. 地理探测器方法下北京市人口空间格局变化与自然因素的关系研究 [J]. 自然资源学报, 2017, 32(8): 1385−1397. [50] 张舒瑾, 余珮珩, 白少云, 等. 面向国土空间规划的流域景观时空分异特征及驱动因子研究 [J]. 生态经济, 2020, 36(10): 219−227. [51] Xi Z L, Xu H D, Xing Y Q, et al. Forest canopy height mapping by synergizing ICESat–2, sentinel-1, sentinel-2 and topographic information based on machine learning methods [J]. Remote Sensing, 2022, 14(2): 364. DOI: 10.3390/rs14020364
[52] Xing Y Q, de Gier A, Zhang J J, et al. An improved method for estimating forest canopy height using ICESat-GLAS full waveform data over sloping terrain: a case study in Changbai Mountains, China [J]. International Journal of Applied Earth Observation and Geoinformation, 2010, 12(5): 385−392. DOI: 10.1016/j.jag.2010.04.010
[53] Su Y J, Guo Q, Xue B L, et al. Spatial distribution of forest aboveground biomass in China: estimation through combination of spaceborne lidar, optical imagery, and forest inventory data [J]. Remote Sensing of Environment, 2015, 173: 187−199.
[54] Fayad I, Baghdadi N, Alcarde Alvares C, et al. Terrain slope effect on forest height and wood volume estimation from GEDI data [J]. Remote Sensing, 2021, 13(11): 2136. DOI: 10.3390/rs13112136
[55] 王继雄, 唐文静. 基于森林资源规划设计调查数据的乔木林地上生物量及碳储量估算分析: 以香格里拉市为例 [J]. 绿色科技, 2021, 23(8): 14−16.
计量
- 文章访问数: 116
- HTML全文浏览量: 13
- PDF下载量: 45
