落叶松人工林三维绿量测算

王璐瑶，贾炜玮，孙毓蔓，胡 涛，郭昊天

（东北林业大学 a.林学院；b.森林生态系统可持续经营教育部重点实验室，黑龙江 哈尔滨 150040）

三维绿量又称三维绿色生物量（Tridimensional green biomass，TGB），简称绿量，是指生长中绿色植物的茎、叶所占据的空间体积[1]。三维绿量能够反映不同绿地垂直空间的结构差异，是准确分析生态环境效益的一个指标，因此准确估算三维绿量对于对城市森林健康发展具有重要的作用[2]。目前，传统的三维绿量测量方法大多采用光学遥感数据结合“冠形公式法”[3-8]，该方法将树冠冠形看作近似几何体，求几何体体积即视为三维绿量。该方法存在误差主要有以下2 方面原因：首先，光学遥感数据精度并不高；其次，树冠内部存在无效树冠部分，用树冠描述实心规则几何体并不准确[8]。

激光雷达具有密度大、精度高、非破坏性等特点，因此常被用于森林资源调查和植被三维结构参数推算。通过对地基激光雷达（Terrestrial laser scanning，TLS）关于森林测绘应用的研究发现，利用地基激光雷达提取树冠参数可有效解决光学遥感数据精度不高的问题[9-11]。刘鲁霞等[12]以地基激光扫描的白皮松为研究对象，采用体元法和凸包算法对单木因子及树冠体积进行估算，实现了利用地基激光雷达提取数据建立树冠体积模型的研究。韦雪花等[13]利用TLS 数据分别采用“体元模拟法”和传统树冠体积计算方法对树冠体积进行估算并进行了误差对比，结果证明“体元模拟法”计算树冠体积结果更接近真实三维绿量值。

为此，本研究以黑龙江省桦南县孟家岗林场落叶松人工林Larixgmelinii为对象，通过TLS 扫描LiDAR 360 软件处理可得到单木因子参数；采用“体元模拟法”通过Matlab 编程实现4 块样地（208 株）内单木三维绿量计算，与“冠形公式法”计算的三维绿量进行精度分析；以地基提取树冠因子与Matlab 计算的三维绿量构建单木三维绿量模型；根据构建的单木三维绿量模型计算42 块样地中3 366 株样木的每木三维绿量，通过样地中每木三维绿量累加可得到样地的三维绿量值，结合样地面积大小即可转化为林分三维绿量，再通过林分三维绿量与林分因子建立林分三维绿量模型。通过以上方法，三维绿量的估算就不仅仅局限于单木，可拓展为区域落叶松林三维绿量的测算，为今后研究区域三维绿量生态效益价值分析做铺垫。

1 研究区概况与数据源

1.1 研究区概况

研究区域位于黑龙江省佳木斯市桦南县孟家岗林场（图1），其地理坐标为130°32＇～130°52＇E，46°20＇～46°30＇N，海拔高度170～575 m，平均海拔250 m。林场地处温带湿润地区，夏季短暂湿热，全年气温变化范围在-34.7～35.6 ℃，年平均气温2.7 ℃，年降水量为550 mm。

图1 研究区域概况Fig.1 Overview of the study area

1.2 数据采集

1.2.1 外业数据采集

根据不同立地条件选择42 块落叶松林样地用于林分三维绿量模型构建，选择其中4 块样地（208株样木）用于单木三维绿量模型构建，样地基本情况如表1 所示。

表1 样地因子统计表Table 1 Statistics of the plot factors

2018—2020 年在黑龙江省孟家岗落叶松林场采集每木检尺数据，样地面积为20 m×30 m、20 m×20 m。根据每木检尺数据计算林分因子，包括平均胸径（Dg）、平均树高（H）、每公顷株树（n）、蓄积量（M）、林分密度指数（SDI）等。根据不同树冠大小选择其中的4 块落叶松样地，测定样地内每株单木的胸径、第一活枝高、树高、冠长（CL）和冠幅（CW）等单木因子及相对坐标。

1.2.2 激光雷达数据获取和处理

2018 年6 月，利用LiBackpack 背包激光雷达扫描系统对4 块样地进行数据采集。2019 年10 月，利用Trimble TX8 系统（Trimble，Sunnyvale，CA，USA）对4 块落叶松样地进行地基激光雷达数据采集。将采集配准处理后的点云数据导入LiDAR 360 软件进行预处理。

2 研究方法

2.1 单木因子提取

采用人机交互的方法从LiDAR360 软件中提取树高、胸径、冠幅和冠长，如图2 所示。保留第一活枝至梢头部分树冠点云，并导出为.xyz 格式便于后期导入Matlab 软件数据库，通过编程对落叶松的三维绿量进行计算。

图2 单木因子提取示意Fig.2 Schematic diagram of individual tree factor extraction

2.2 单木因子提取效果

为验证采用地基激光雷达提取的单木因子精度，将外业测量的单木因子与LiDAR 360 软件中提取的相应因子进行精度对比分析，计算公式如下：

式中：xi为第i棵落叶松地基提取树冠因子；Xi为第i棵落叶松实测因子；n为样地内提取的落叶松株树；p为精度。

2.3 单木三维绿量计算方法及误差分析

三维绿量的传统计算方法是采用“冠形公式法”，通过实地考察与查阅相关文献确定落叶松树冠形态与圆锥体相似，如图3 所示。因此，圆锥体体积即视为落叶松单木三维绿量[14]，传统三维绿量计算方法为：

图3 三维绿量传统计算方法示意Fig.3 Schematic diagram of the traditional tridimensional green biomass calculation method

式中：Vi为第i株落叶松三维绿量；xi为第i株落叶松的树冠冠幅，m；yi为第i株落叶松的树冠冠长，m；xi、yi均为地基激光雷达测量值。

2.3.1 三维绿量计算

“体元模拟法”是一种基于地基激光雷达数据，用固定大小的体元模拟不规则树冠形状的一种方法。根据韦雪花等[13]利用地基激光雷达数据计算树冠体积的方法，本研究采用“体元模拟法”通过Matlab 编程实现对落叶松三维绿量的计算。具体如下：首先，在Matlab 软件中导入.txt 格式的树冠点云数据坐标信息文件，通过编程将树冠激光点云数据沿垂直方向以间距k进行等距划分，得到若干个树冠切片；水平方向对每一层冠层进行边长为k的间距进行网格划分，树冠分割若干个k×k×k的小方格，计算每一层含有树冠点云的有效网格数，有效网格的体积之和为有效树冠体积[13,15-16]。在计算机可接受范围内，尽可能地缩小网格间距k，间距k越小，精度越高，体元模拟法如图4 所示。

图4 体元模拟法示意图Fig.4 Schematic diagram of the voxel simulation method

“体元模拟法”的体元大小直接影响落叶松三维绿量计算结果。因此，为寻找模拟落叶松树冠的最佳体元大小，从208 株落叶松中挑选出不同冠幅大小30 株落叶松，分析不同体元大小下的三维绿量值。k的取值分别为0.4、0.3、0.2、0.1、0.09、0.08、0.07、0.06 和0.05 m，进行多次运算，通过对比分析，选出适用于不同冠幅大小的树冠体积计算值趋于稳定时的最佳体元k值。

2.3.2 计算误差分析

根据刘常富等[5]和李乐[17]的研究结果可知传统“冠形公式法”测量三维绿量值明显大于真实三维绿量值[17]；通过刘红伟[18]的研究证实了对树冠进行无限切分组合的体积比“冠形公式法”所得体积更精确；通过韦雪花等[13]的研究证实了“体元模拟法”作为三维绿量真值的可行性。为验证基于TLS 扫描数据采用“体元模拟法”计算三维绿量准确程度，本研究以传统三维绿量估算方法“冠形公式法”为参照，采用2 种方法计算208 株落叶松的树冠体积，对比2 种方法测量体积的相对误差，计算公式见公式（3）。

式中：vi第i棵落叶松树冠“体元模拟法”计算所得体积；Vi为第i棵落叶松树冠传统方法计算所得体积；ΔEi为第i棵落叶松树冠体积计算相对误差。

2.4 林分三维绿量计算方法

林分三维绿量（TGBs）则根据已经建立好的单木三维绿量模型计算样地中的每木三维绿量，通过样地中每木三维绿量累加可得到样地的三维绿量值，结合样地面积大小即可转化为林分单位面积三维绿量即林分三维绿量[19]。见公式（4）。

式中：TGBsi为第i块样地的林分三维绿量（kg/hm2），TGBij为第i块样地中第j棵树的三维绿量值，n为样地株树，Si为第i块样地的样地面积，本研究样地面积大小为0.06、0.12 hm2。

2.5 三维绿量模型构建与检验

2.5.1 三维绿量模型选取

通过对208 株落叶松单木因子与三维绿量的散点趋势图可看出单木因子与三维绿量值均呈明显的幂函数关系，如图5 所示。因此，采用对数转换的线性回归模型构建单木三维绿量模型，见公式（5）[20]。

图5 树冠因子与单木三维绿量的关系趋势图Fig.5 Trend diagram of the relationship between the canopy factor and tridimensional green biomass of individual trees

式中：xj是各单木因子；βj为待估参数。

2.5.2 模型拟合与检验

以208 株落叶松为研究对象构建单木三维绿量模型，对因变量和自变量作相关性分析。采用逐步回归方法进行多元三维绿量模型构建，以4个单木因子lnDBH、lnH、lnCW、lnCL为自变量，以lnTGB为因变量。采用逐步回归分析中的逐步筛选法筛选方程变量，根据赤池信息准则（AIC）依次选择变量，采用方差膨胀因子（Variance inflation factor，VIF）诊断法排除多重共线性可能，临界值设置为5，仅选择参数具有显著性且VIF值小于5 的变量作为有效变量并进行最优模型拟合。输出结果包括单因子最优模型、双因子最优模型、三因子最优模型。

以42 块落叶松样地为研究对象构建林分三维绿量模型，根据相关文献研究，以实测林分因子lnM、lnSDI、lnN、lnDg、lnHt为备选自变量[11,19,21]；以单木三维绿量累积结合样地面积转换所得林分三维绿量lnTGBs为因变量。同样对因变量和自变量作相关分析，采用逐步回归分析中的逐步筛选法筛选方程变量并构建林分三维绿量模型。

本研究将样本数据按3∶1 的比例随机抽取分别用于建模和检验，使用R 语言对模型进行线性回归拟合。选择确定系数（R2）、调整后的确定系数（Ra2）、均方根误差（RSME）和赤池信息准则（AIC）评价模型的拟合优度，其中模型的R2越大，RSME和AIC越小，模型的拟合效果越好；选择平均绝对误差（MAE）、平均绝对百分比误差（MAPE）和预估精度（P）对模型进行检验，MAE和MAPE越小，P%越大，模型预估效果越好。

3 结果与分析

3.1 单木因子提取结果

基于208 株地基激光雷达数据，提取胸径、冠幅、冠长和树高，结果见表2。其中，胸径提取精度整体较高，平均提取精度为97.14%，树高平均提取精度为90.96%，冠幅提取精度变化比较大，平均提取精度为89.25%，主要原因是样地内树冠之间存在遮挡，单木分割存在缺冠或过分分割树冠的情况。冠长平均提取精度为77.92%，主要原因是部分样地密度较大，单木数量较多，扫描时受天气的影响，扫描影像模糊不清，地基提取第一活枝高效果不佳，因此导致冠长提取精度差别较大。单木因子提取精度最低为77.92%，因此，在可允许的误差范围内，采用地基激光雷达获取样地树冠因子数据取代实测数据具有实际意义。

表2 单木因子提取结果Table 2 Extraction results of individual tree factors

3.2 体元模拟法体元大小分析

为寻找模拟落叶松树冠的最佳体元大小，本研究从208 株落叶松中挑选出30 株不同冠幅大小的落叶松，对比分析r与V之间的变化趋势，结果如图6 所示。发现V与r成幂函数关系，V随r的增大而减小，到达一定界限V值趋于稳定，且冠幅大小对指标r有影响，冠幅越大，三维绿量趋于稳定时的r值越大。

图6 体元模拟法计算树冠三维绿量大小随r 值变化趋势Fig.6 The variation trend of canopy tridimensional green biomass with r values calculated by the voxel simulation method

由指标公式k=K/r（K为树冠冠幅，k为体元大小，r为指标公式）可知，通过V值趋于稳定时r值的大小，即可推算出最佳体元k值的取值范围。因此，本研究以1 m 为区间对树冠进行分类，通过分析不同冠幅大小情况下三维绿量值变化规律确定最佳体元k值。从表3 可以看出，当冠幅范围在1～2 m 时，k取值小于0.1；当冠幅范围在2～3 m 时，k取值小于0.11，当冠幅范围在3～4 m 时，k取值小于0.12；当冠幅范围在4～5 m 时，k取值小于0.11，当冠幅范围在5～6 m 时，k取值小于0.10。通过取平均值，三维绿量计算值趋于稳定的体元大小值应小于等于0.108 m。同时，体元大小影响计算机计算三维绿量运行时间，由图7可以看出，体元大小与计算时长成反比，体元越小，计算时间越长，随着体元变大，计算时间逐渐减少。综上所述，基于对三维绿量计算的准确性和计算时间的考虑，在单株落叶松三维绿量计算时长约为600 s 时，取0.1 m 确定为最佳k值。

表3 k 稳定值Table 3 k stable values

图7 体元大小与单木三维绿量计算时长关系Fig.7 Relationships between volume size and calculation time of individual wood tridimensional green biomass

3.3 体元模拟法计算结果分析

以208 株落叶松单木为研究对象，分别采用“体元模拟法”与“冠形公式法”计算三维绿量，以传统三维绿量估算方法为参照，对比分析“体元模拟法”三维绿量相对误差，对比分析结果见表4。传统“冠形公式法”与“体元模拟法”相对误差最小值为0.43%，最大值为78.80%，平均相对误差为37.25%。结果表明，传统三维绿量计算方法与“体元模拟法”计算三维绿量具有显著差异，原因主要有以下几点：

表4 体元模拟法与传统公式法测算结果对比分析Table 4 Comparative analysis of the calculation results of the voxel simulation method and traditional formula method

1）传统“冠形公式法”在测量冠幅时采用目测法，即通过观察树冠冠幅轮廓最外侧点的位置确定。目测本身存在一定误差，所测出的冠幅可能存在偏差，本研究采用的“体元模拟法”通过机器识别点云空间坐标位置，可有效避免人为因素的误差。

2）传统三维绿量计算方法是将落叶松树冠形状视为规则的圆锥体，而实际落叶松人工林树冠冠形并非完整的圆锥体。其次，阔叶林叶片大，枝叶密度高，落叶松叶片较小，且树冠内部存在较大空隙。按传统计算方法计算三维绿量包含的无效树冠部分较多，使得2 种方法计算三维绿量值相对误差也会偏大。

通过以上分析，可以得出采用“体元模拟法”比传统“冠形公式法”计算三维绿量更精确；同时也说明了传统方法技术粗糙，不能满足当今精准森林资源检测要求。因此，本研究采用“体元模拟法”计算树冠点云体积稳定值作为三维绿量的真值。

3.4 单木三维绿量模型拟合与检验

树冠单木因子与单木三维绿量的相关性分析结果如表5 所示。胸径、冠幅、冠长和树高与三维绿量的相关系数均大于0.700，说明4 个单木因子均与三维绿量存在较强的相关性。其中，胸径与三维绿量相关性最强，相关系数为0.885；其次，冠幅相关系数为0.876，树高相关系数为0.844，冠长的相关系数最低，为0.716。

表5 树冠参数与三维绿量相关性分析Table 5 Correlation analysis between canopy parameters and tridimensional green biomass

采用R 语言中逐步回归分析中的逐步筛选法对树冠参数lnDBH、lnH、lnCW、lnCL与lnTGB进行多元线性拟合，方程的预估参数和拟合指标见表6（未全部列出）。其中，单因子三维绿量模型的R2、Ra2均在0.78 以上，AIC值为-345.480；当增加到2 个变量时模型精度提高，确定系数R2、均在0.84以上，AIC值减小为-398.940；当增加到3个变量时模型确定系数R2、Ra2均在0.846以上，AIC值减小为-402.960；当增加到4 个变量时，模型AIC值增加，无法进入模型，因此排除4 因子模型。综上所述，模型的拟合效果为模型3 ＞模型2 ＞模型1，模型3 调整后的确定系数Ra2值最高为0.846，均方根误差RMSE值最小为0.250，AIC值最小为-402.960，因此模型3 拟合效果最好。

表6 单木三维绿量模型参数估计及拟合指标Table 6 Parameter estimation and fitting indexes of individual wood tridimensional green biomass model

通过模型的检验结果表7 可以看出，模型3的平均绝对误差（MAE）值为0.210，平均绝对百分比误差（MAPE）值为0.099，精度（P）值为97.314%。在3个模型中，模型3的MAE和MAPE最小，P值最高，图8 为单木模型检验残差分布，可以看出模型3 的残差值较均匀地分布于0 值上下两侧，分布区间为-0.75～0.75，因此模型3 的预估效果最好。

表7 单木三维绿量模型检验结果Table 7 Test results of individual wood tridimensional green biomass model

图8 单木三维绿量模型检验残差分布Fig.8 Residual distribution of individual wood tridimensional green biomass model

综上所述，模型3 为最优模型，落叶松单木三维绿量模型表示为：

lnTGB=1.374lnx2+0.282lnx3+0.911lnx4-2.299。

式中：TGB为单木三维绿量，x2为冠幅，x3为冠长，x4为树高。

3.5 林分三维绿量模型拟合与检验

以42 块落叶松样地为研究数据，通过逐步回归分析中的逐步筛选法对林分因子lnN、lnSDI、lnM、lnDg、lnHt与lnTGBs进行多元线性拟合。模型拟合过程中发现林分因子密度指数（SDI），平均胸径（Dg），平均数高（Ht）存在共线性问题，采用逐步筛选法选择拟合效果最好的双因子林分三维绿量模型，模型拟合检验结果见表8。其中，确定系数（R2）为0.706，均方根误差（RMSE）为0.138 m3；平均绝对误差（MAE）值为0.120，平均绝对百分比误差（MAPE）值为0.013，精度（P）值达到99.142%。

表8 林分三维绿量值模型参数拟合结果Table 8 Parameter fitting results of the tridimensional green biomass model of stands

综上所述，落叶松林分三维绿量的模型形式为：

lnTGBs=4.749+0.372lnx1+0.386lnx2。

式中：TGBs为林分三维绿量，x1为每公顷株树，x2为每公顷蓄积量。

4 结论与讨论

4.1 结 论

本研究选取黑龙江省佳木斯市孟家岗林场的落叶松人工林为研究对象，获取了样地的地基激光雷达数据和每木检尺数据，提取单木因子参数，结合实测数据对比其提取精度；利用TLS 数据结合“体元模拟法”计算单木三维绿量，结合“冠形公式法”计算的单木三维绿量进行误差分析。通过单木三维绿量与提取单木因子构建单木三维绿量模型；结合42 块样地每木检尺数据通过累积转换得到林分三维绿量，结合林分因子构建林分三维绿量模型。研究结论如下：

1）单木参数提取精度。TLS 提取树冠参数的精度较高。其中，胸径平均提取精度为97.14%，树高平均提取精度为90.96%，冠幅平均提取精度为89.25%，冠长平均提取精度为77.92%，树冠因子最低精度达到77.92%，提取精度满足要求，可用于单木三维绿量模型构建。

2）最佳体元大小。采用“体元模拟法”计算单木三维绿量，三维绿量值随k值的减小而减小，到达临界值后趋于稳定。单木三维绿量计算值趋于稳定的平均时长为600 s，确定k=0.1 m 为最佳体元大小。

3）单木三维绿量模型。相关性分析结果表明，lnTGB与lnDBH、lnCW、lnH和lnCL存在强相关性，皮尔逊相关系数分别依次为0.885、0.876、0.844和0.716；通过单木三维绿量模型的拟合与检验，确定模型3 的拟合效果最好。模型3 调整后的确定系数Ra 2值最高为0.846，均方根误差RMSE值最小为0.250 m3，AIC值-402.960 最小；平均绝对误差（MAE）值最小为0.210，平均绝对百分比误差（MAPE）值最小为0.099，精度（P）值最高为97.314%，确定为最优单木三维绿量模型。

4）林分三维绿量模型。通过林分三维绿量模型的拟合与检验，确定双因子林分三维绿量模型的拟合效果最好。模型拟合检验指标如下：确定系数（R2）最大为0.706，均方根误差（RMSE）最小为0.138 m3；平均绝对误差（MAE）值最小，为0.120，平均绝对百分比误差（MAPE）值最小为0.013，精度（P）值达到99.142%。

4.2 讨 论

本研究采用“体元模拟法”[14-15]和Matlab 编程实现树冠点云体积的计算，该方法有效地剔除了树冠内部无植物茎叶的部分，实现了三维绿量计算方法优化。在三维绿量数据采集方面，采用实测每木检尺数据和高精度的TLS 数据，相比传统光学遥感数据更精确。本研究通过构建单木三维绿量模型及林分三维绿量模型，为估算区域三维绿量提供新尝试，对森林资源调查和森林资源动态检测具有重要意义。

本研究只采用了黑龙江省孟家岗地区落叶松人工林样地数据，所以建立的单木三维绿量模型与林分三维绿量模型适用范围有限。本研究只采用了实测数据与TLS 数据进行模型构建，利用TLS 数据研究单木模型的适用性更强，对于大尺度林分三维绿量模型的构建最好能够结合机载激光雷达数据。另外，在样地密度较大、样地单木数量较多时，树冠顶端识别可能出现树高测量不准确的情况。今后可利用TLS 数据结合机载激光雷达数据对大尺度区域进行数据采集，提高树冠参数的采集精度。本研究在构建林分三维绿量模型时样本量比较少，林分三维绿量最优模型的拟合参数也只有每公顷蓄积与每公顷株树，林分三维绿量模型还可能受其他因素的影响，在后续的研究中可以考虑扩大样本数量、丰富样本的林分类型以提高林分模型的拟合精度。