基于无人机激光雷达的人工红松直径分布模型

Diameter Distribution Model of Korean Pine Plantation by Using UAV-LiDAR

摘要

林分直径分布是反映林分结构最重要的指标之一，有助于估算森林蓄积、材种出材量、评价林分稳定性。无人机激光雷达（UAV-LiDAR）技术通过提供高精度的三维空间信息，为大范围、连续性的林分直径分布预测提供了新的技术手段。以UAV-LiDAR特征因子和实测直径分布为数据源，利用3种参数预测方法（参数预测法、百分位参数回收法和矩估计参数回收法）构建Weibull直径分布模型，通过决定系数（R²）、均方根误差（R_MSE）、平均偏差（M_AE）、误差EI（E_IP和E_IR）等指标对模型进行拟合和检验，确定最优的人工红松直径分布模型。基于UAV-LiDAR的Weibull参数预测精度R²在0.57~0.71，3种参数预测方法建立的直径分布模型均能较好地预测林分直径分布，其中参数预测法表现最佳（E_IR=57.384，E_IP=0.319）。UAV-LiDAR技术能够有效弥补传统调查中数据获取效率低和覆盖范围小等问题，为实现可靠的林分直径分布预测提供技术支撑。

Abstract

Stand diameter distribution is one of the most important indicators reflecting stand structure， aiding in estimating forest biomass， timber yield， and evaluating stand stability. Unmanned aerial vehicle （UAV）-LiDAR technology provides high-precision three-dimensional spatial information， offering a new technical approach for large-scale and continuous predictions of stand diameter distributions. This study utilized UAV-LiDAR derived features and measured diameter distributions as data sources to develop Weibull diameter distribution models using three parameter prediction methods： the parameter prediction method， percentile-based parameter recovery method， and moment-based estimation parameter recovery method. The model’s performance was evaluated through multiple metrics including， R²， R_MSE， M_AE， E_IP， and E_IR， ultimately identifying the optimal Weibull model for Korean pine planatation diameter distributions. The precision of Weibull parameter predictions based on UAV-LiDAR data showed an R² ranging from 0.51 to 0.71. The diameter distribution models established by the three parameter prediction methods can predict the diameter distribution of the stand， and the parameter prediction method performed best （E_IR=57.384， E_IP=0.319）. UAV-LiDAR technology effectively addresses the limitations of conventional surveys， such as low data acquisition efficiency and restricted spatial coverage， providing a robust technical foundation for reliable predictions of stand diameter distributions.

Graphical abstract

关键词

无人机激光雷达 / 红松人工林 / Weibull函数 / 参数预测方法 / 直径分布模型

Key words

UAV-LiDAR / Korean pine plantation / Weibull function / parameter estimation method / diameter distribution model

0 引言

林分直径分布是指林分内各径阶树木的胸径分布密度特征，其反映了林分中不同大小树木的数量分配情况，是表征人工林林分结构的重要属性之一^［1-8］。同时，林分结构的动态变化能够有效反映林木在生长过程中受到的环境因子（如气候、土壤等）和经营措施（如间伐、抚育等）的综合影响。这种动态变化不仅体现了林木个体的生长响应，也反映了林分整体的发育过程，为森林经营决策提供了重要的科学依据^［9-10］。红松作为东北林区乡土优势树种之一，广泛分布于长白山脉以及小兴安岭等区域，通过建立该树种的直径分布模型可以揭示红松人工林直径结构的发展规律，为开展科学的经营管理提供依据^［11］。在传统人工林经营中，森林资源的动态监测过程主要通过长期固定样地调查，虽然该方法可直接获取相对准确的林分调查因子，但往往耗时费力，固定样地调查数据中往往存在同木异号、错测漏测等^［12］问题。而无人机激光雷达技术^［13］具有穿透植被冠层直接提供高精度三维信息的能力，在评估林分结构中展示出了巨大的优势。基于激光雷达的林分结构研究中，关于直径分布模型的构建上主要采用2种方法^［14］：单木检测法（individual tree delineation，ITD）和基于面方法（area-based approach，ABA）。对于ITD方法而言，建立直径分布模型过程是通过单木分割算法进行单木探测，利用所提取的单木参数构建直径分布模型。在复杂林分中单木分割算法的探测精度受限，主要问题在于单木信息提取错误，从而影响ITD方法的估测精度^［15］。相比之下，ABA方法是通过对与样地大小一致的每个栅格单元生成一系列林分特征因子来建立直径分布模型。ABA方法在直径分布模型中的预测精度通常高于ITD方法^［16］。因此，ITD方法通常作为ABA方法的辅助手段用于提升精度，弥补其在冠层下树木探测精度的不足^［17］。

ABA方法通常采用概率密度函数（probability density function，PDF）构建分布模型，该方法先提取样地尺度的激光雷达特征参数来预测林分株数密度及其PDF参数，最后结合实测数据建立PDF模型。近年来，许多学者针对不同树种尝试了多种概率密度函数来构建直径分布模型，其中包括Weibull分布、Johnson′s sb分布和Beta等^［18-21］建立相应的直径分布模型。萨其拉^［22］研究表明，使用2参数的Weibull函数建立的红松人工林直径分布模型相较于其他概率密度函数具有更高的预测精度。目前，Weibull参数的估计方法，主要分为参数预测法（parameter prediction method，PPM）和参数回收法（parameter recovery method，PRM）^［23］。其中PPM通常采用极大似然估计法（maximum likelihood estimate，MLE）估计直径分布的参数，随后基于林分信息采用最小二乘法（ordinary least squares，OLS）、最大似然估计回归（maximum likelihood estimation based regression，MLER）等预测模型的系数^［24］。PRM则使用一种间接的方法，首先估计与分布相关的其他参数，然后使用这些参数来预测Weibull参数。PRM可进一步分为基于百分位数的参数回收法（percentile-based parameter recovery method，PPRM）和矩估计参数回收法（moment-based parameter recovery method，MPRM）^［25］。PPRM在以往的研究中通常使用直径的24百分位数和93百分位数建立直径分布模型^{［24，26］}；MPRM主要通过林分平均胸径和算术平均胸径的转换关系预测Weibull参数^［27］。例如，Cao等^［26］利用机载激光雷达数据，比较了百分位参数回收法、矩估计参数回收法和参数预测法所构建的亚热带森林的材积分布预测效果，其中PPRM（E_IR=37.75）的效果最好。Gobakken等^［28］基于机载激光雷达数据分别构建PPRM和PPM的林分直径和断面积分布模型，结果表明PPRM（偏差：-4.8%~2.7%）比PPM方法（偏差：-4.7%~6.6%）更适用于估测林分直径和断面积的分布。Thomas等^［29］在加拿大安大略省中部北方森林中基于PPM方法建立了单峰双参数Weibull模型，该方法在特定林型中具有较高精度（R_{\mathrm{a}}^{\mathrm{2}}=0.65~0.88）。在构建Weibull直径分布模型时，通常采用全子集法对激光雷达特征变量进行筛选，以决定系数和方差膨胀因子等指标来确定最佳特征变量的组合，基于上述特征因子来预测Weibull分布参数和林分的株数密度^{［22，30］}。

因此，本研究旨在利用无人机激光雷达（UAV-LiDAR）数据，分别采用3种Weibull参数预测方法（参数预测法、百分位数回收法和矩估计参数回收法）建立人工红松直径分布模型，并与实测数据进行比较从而得到最优的直径分布模型。为基于UAV-LiDAR数据构建人工红松林直径分布模型提供科学依据和技术参考，有助于提升人工红松林的经营管理水平。

1 材料与方法

1.1 研究区概况

本研究使用的红松人工林数据来源于中国东北地区的5个代表性林场，分别是黑龙江省佳木斯市孟家岗林场（130°32′—130°52′E，46°20′—45°30′N）、尚志市尚志林场（127°39′—129°35′E，45°0′—45°23′N）、尚志市苇河林场（128°10′—128°47′E，44°01′—44°52′N）、牡丹江市东京城林业局（128°07′—130°02′E，43°30′—44°18′N）以及吉林省临江市闹枝林场（125°17′—127°51′E，41°20′—42°28′N）。其中，孟家岗林场位于完达山西麓余脉，东京城、尚志、苇河林场位于张广才岭西坡，闹枝林场则位于长白山主峰西麓。气候整体属于温带或中温带大陆性季风气候，年均气温2.3~3 ℃，年降水量为550~800 mm。

1.2 数据采集

1.2.1 实测数据收集

样地数据采集于2023年6—9月，分别在黑龙江省佳木斯市孟家岗林场、牡丹江市东京城林业局、尚志市尚志林场、尚志市苇河林场和吉林省临江市闹枝林场设置共220块样地，样地详细统计信息见表1。

1.2.2 无人机激光雷达数据收集

本研究于2023年7—9月期间，在5个林场开展激光雷达点云数据采集工作，以大疆DJI M600 Pro八旋翼无人机为飞行平台，搭载RIEGL VUX-1UAV激光雷达扫描仪。数据采集时无人机飞行高度为180 m，飞行速度为10 m/s，飞行时采用交叉航带飞行，旁向重叠率为50%，激光雷达扫描视场角为330°，脉冲频率为380 kHz，最大回波数5次。

1.2.3 点云数据处理

首先，对原始点云数据进行航带拼接，随后进行点云去噪，然后利用改进渐进加密三角网滤波算法^［31］（improved progressive TIN densification，IPTD）进行点云滤波，将点云分为地面点和非地面点，利用地面点进行高度归一化处理，消除地形起伏对植被高度的影响。利用实时动态载波相位差分技术（real time kinematic，RTK）获取的样地坐标从而确定样地边界，并提取每个样地的样地尺度激光雷达特征变量，共计101个，见表2。

1.3 研究方法

1.3.1 Weibull概率密度函数

本研究使用2参数Weibull概率密度函数作为直径分布模型。其公式为

f\left(x\right)=\left(\frac{c}{b}\right){\left(\frac{x}{b}\right)}^{c-\mathrm{1}}\mathrm{e}\mathrm{x}\mathrm{p}\left[-{\left(\frac{x}{b}\right)}^c\right],x\ge \mathrm{0},b,c\ge \mathrm{0}。（1）

式中：x为林分中的胸径；b、c分别是直径分布的尺度参数和形状参数。

1.3.2 Weibull参数拟合方法

本研究采用参数预测法（parameter prediction method，PPM）、百分位参数回收法（percentile-based parameter recovery method，PPRM）和矩估计参数回收法（moment-based parameter recovery method，MPRM）对基于UAV-LiDAR数据的直径分布模型进行参数预测。

1）参数预测法

参数预测法是一种直接求解Weibull分布参数的方法，首先基于样地实测数据通过最大似然估计法直接计算出Weibull分布参数b和c，然后构建激光雷达特征变量与参数b、c的模型。最大似然估计法公式为

\begin{array}{l}\mathrm{l}\mathrm{n}L=\left(c-\mathrm{1}\right){\displaystyle \sum _{i=\mathrm{1}}^n}{p}_i\mathrm{l}\mathrm{n}\left(d_i\right)-P\left(c-\mathrm{1}\right)\mathrm{l}\mathrm{n}\left(b\right)+\\ P\mathrm{l}\mathrm{n}\left(\frac{c}{b}\right)-{\displaystyle \sum _{i=\mathrm{1}}^n}{p}_i{\left(\frac{d_i}{b}\right)}^c\end{array}。（2）

式中：L为似然函数；p_i是第i棵树的频率；P为总频率；d_i为第i棵树的胸径。

2）矩估计参数回收法。

矩估计参数回收法是利用Weibull概率密度函数的参数与林分平均直径的关系建立联立方程，具体估测直径的一阶矩和二阶矩，然后通过估测的矩来求解尺度参数b和形状参数c。在基于激光雷达数据的矩估计参数回收法中，首先通过样地尺度激光雷达特征因子对一阶矩和二阶矩进行预估，进而间接计算得到Weibull模型参数。矩估计参数回收法通常使用算术平均直径和平均直径作为一阶和二阶矩。其公式为

D^{\mathrm{2}}-d^{\mathrm{2}}=\frac{d^{\mathrm{2}}}{{\Gamma }^{\mathrm{2}}\left(\mathrm{1}+\frac{\mathrm{1}}{c}\right)}\left(\Gamma \left(\mathrm{1}+\frac{\mathrm{2}}{c}\right)-{\Gamma }^{\mathrm{2}}\left(\mathrm{1}+\frac{\mathrm{1}}{c}\right)\right)。（3）

b=\frac{d}{\Gamma \left(\mathrm{1}+\frac{\mathrm{1}}{c}\right)}。（4）

式中：d和D分别是样地内林木算术平均直径和平均直径；Γ表示gamma函数。

3）百分位法参数回收法

基于百分位法的参数回收法是利用Weibull分布的已知百分位数与样本数据的百分位数进行匹配来估计参数。林木直径的百分位数p（按升序排列）可以通过胸径小于Q_p 的林木直径所占的百分比计算得到，Weibull的累积分布函数，公式为

F\left(Q_p\right)=\mathrm{1}-\mathrm{e}\mathrm{x}\mathrm{p}\left(-{\left(\frac{Q_p}{b}\right)}^c\right)。（5）

式中：Q_p 为整个林分直径按百分位划分的直径，将上述表达式进行转换得到Q_p 表达式。其公式为

Q_p=b{\left(-ln\left(\mathrm{1}-F\left(Q_p\right)\right)\right)}^{\frac{\mathrm{1}}{c}}。（6）

由表达式可知，Weibull模型参数（Q₂₄和Q₉₃）可以通过基于2个百分位点的联立方程求解。本研究使用24百分位和93百分位（Q₂₄和Q₉₃）来回收Weibull模型的2个参数。首先采用UAV-LiDAR提取的样地尺度特征因子作为自变量对Q₂₄和Q₉₃进行预估，进而间接计算Weibull模型参数b和c。其公式为

b=\mathrm{e}\mathrm{x}\mathrm{p}\left[\frac{\mathrm{l}\mathrm{n}\left(Q_{\mathrm{93}}\right)ln\left[-\mathrm{l}\mathrm{n}\left(\mathrm{1}-\frac{\mathrm{24}}{\mathrm{100}}\right)\right]-\mathrm{l}\mathrm{n}\left(Q_{\mathrm{24}}\right)\mathrm{l}\mathrm{n}\left[-\mathrm{l}\mathrm{n}\left(\mathrm{1}-\frac{\mathrm{93}}{\mathrm{100}}\right)\right]}{\mathrm{l}\mathrm{n}\left[-\mathrm{l}\mathrm{n}\left(\mathrm{1}-\frac{\mathrm{24}}{\mathrm{100}}\right)\right]-\mathrm{l}\mathrm{n}\left[-\mathrm{l}\mathrm{n}\left(\mathrm{1}-\frac{\mathrm{93}}{\mathrm{100}}\right)\right]}\right]。（7）

c=\frac{\mathrm{l}\mathrm{n}\left[-\mathrm{l}\mathrm{n}\left(\mathrm{1}-\frac{\mathrm{93}}{\mathrm{100}}\right)\right]}{\mathrm{l}\mathrm{n}\left(Q_{\mathrm{93}}\right)-\mathrm{l}\mathrm{n}b}。（8）

式中：Q₂₄和Q₉₃分别为直径排列的24百分位点和93百分位点。计算每个样地的Q₂₄和Q₉₃与激光雷达特征变量构建回归模型。

1.3.3 预估模型的构建与评价

建立直径分布模型时，首先基于激光雷达特征因子对参数预测法的参数b和c、百分位参数回收法的Q₂₄和Q₉₃以及矩估计参数回收法的d和D进行预估，用于计算Weibull参数。由于Weibull函数仅表示径阶上的概率分布，还需要株数密度（N）将概率分布转换为数量分布。因此需要建立［b，c，d，D，Q₂₄，Q₉₃，N］7个预估模型。首先提取与因变量相关性高的激光雷达特征因子之后使用全子集法以决定系数最大、变量相关性较高且方差膨胀因子（variance inflation factor，VIF）<10为标准确定基础模型。模型形式为

\mathrm{l}\mathrm{n}Y=k_{\mathrm{1}}{x}_{\mathrm{1}}+k_{\mathrm{2}}{x}_{\mathrm{2}}+k_{\mathrm{3}}{x}_{\mathrm{3}}+\dots +k_n{x}_n+\epsilon。（9）

式中：Y为因变量组合［b，c，d，D，Q₂₄，Q₉₃，N］中的任意变量，使用lnY是为了使自变量与因变量之间形成线性关系；k_{\mathrm{1}}、k_{\mathrm{2}}、k_{\mathrm{3}}、…k_n为自变量的系数；x_{\mathrm{1}}、x_{\mathrm{2}}、x_{\mathrm{3}}…x_n为筛选的激光雷达特征变量；\epsilon为截距项。

1）模型评价指标

对于模型的拟合优度评价，使用决定系数（R²）、均方根误差（R_MSE）、平均偏差（M_AE）来进行评价，各指标计算公式为

R^{\mathrm{2}}=\mathrm{1}-\frac{{\displaystyle \sum _{i=\mathrm{1}}^n}(y_i-{\widehat{y}}_i{)}^{\mathrm{2}}}{{\displaystyle \sum _{i=\mathrm{1}}^n}(y_i-{\overline{y}}_i{)}^{\mathrm{2}}}。（10）

M_{\mathrm{A}\mathrm{E}}=\frac{\mathrm{1}}{n}{\displaystyle \sum _{i=\mathrm{1}}^n}\left|y_i-\widehat{y_i}\right|。（11）

R_{\mathrm{M}\mathrm{S}\mathrm{E}}=\sqrt[]{\frac{\mathrm{1}}{n-\mathrm{1}}{\displaystyle \sum _{i=\mathrm{1}}^n} (y_i-{\widehat{y}}_i{)}^{\mathrm{2}}}。（12）

式中：n为观测值个数；y_i、{\overline{y}}_i和{\widehat{y}}_i分别为模型的观测值、观测值的平均值以及预测值。

2）模型的检验指标

对于基于UAV-LiDAR数据建立的预估模型，本研究采用十折交叉法对模型的预测效果进行检验，即将样本划分为10等份，轮流将一份数据作为测试集，依次轮流将其中一份数据作为测试集，其余9份作为训练集进行模型训练。通过比较十折交叉检验的决定系数与基于全部数据集的模型决定系数，对模型的拟合效果进行综合评估。

对于直径分布的检验，采用Packalén等^［32］提出的误差指数（E_{\mathrm{I}\mathrm{P}}）和Reynolds等^［33］提出的误差指数（E_{\mathrm{I}\mathrm{R}}）。其公式为

E_{\mathrm{I}\mathrm{P}}={\displaystyle \sum _{d=\mathrm{1}}^{n_d}}\mathrm{0.5}\left|\frac{f_d}{N}-\frac{\widehat{f_d}}{\widehat{N}}\right|。（13）

E_{\mathrm{I}\mathrm{R}}={\displaystyle \sum _{d=\mathrm{1}}^{n_d}}\mathrm{100}\left|\frac{f_d-\widehat{f_d}}{N}\right|。（14）

式中：d为第i块样地的径阶；n_d 为这块样地的所有径阶；f_{\mathrm{d}}和\widehat{f_{\mathrm{d}}}分别为d径阶的株数密度的观测值和预测值；N和\widehat{N}分别是样地的株数密度的观测值和预测值。2个误差指数越小代表拟合效果越好，E_{\mathrm{I}\mathrm{P}}的值是0~1。

2 结果与分析

2.1 PPM、PPRM、MPRM的参数及株数密度预测模型的拟合结果

首先，对上述因变量进行对数转换，以提高模型线性化程度并降低异方差；然后，使用皮尔森相关系数筛选出相关性较高的激光雷达特征变量，用全子集法以决定系数最大且方差膨胀因子（VIF）<10为标准确定模型自变量，以此构建预估模型，见表3。7个模型都选择了强度变量和高度变量，Model-c还额外选择了密度变量（D₂）和叶面积指数（L_AI）。在高度相关变量中较低的高度分位数和累计高度百分位在各模型的变量中出现频率较高，其代表着林分中低层的树木高度，对于冠下红松具有一定的解释效果；而强度变量中的强度中位数绝对偏差中位数和变异系数出现频率较高。

通过R²、R_MSE、M_AE评估模型拟合效果，并采用十折交叉检验模型预测效果，见表4。7个模型中，除去N参数模型的R²为0.57之外，其他6组参数模型（b、c、d、D、Q₂₄、Q₉₃）的R²在0.64~0.71，表明对于7个因变量的预测效果整体较好。对于R_MSE和M_AE指标，Q₉₃的R_MSE和M_AE分别为0.12和0.10，是所有模型中最小的，说明该模型在预测Q₉₃误差较小；N和b的R_MSE为0.34和0.32，M_AE为0.26和0.25，说明2个模型的误差相较于其他模型较大，但仍在合理范围内。通过对7个模型进行十折交叉检验，评估了模型预测效果。具体而言，十折交叉检验的3个指标的10次平均值与使用全部数据拟合模型的指标之间差值的绝对值并不超过0.05，表明7个模型不存在明显的过拟合或欠拟合现象。

在建立预估模型时通常使用高度变量作为主要特征来构建预测模型，因为高度与胸径在生物学上具有重要意义，而激光雷达强度信息与胸径并无直接的关系^［34］。因此，本研究为验证强度变量对模型的贡献作用，在去除强度变量的情况下，使用相同的方法建立预估模型，模型形式以及精度见表5。结果表明，不含强度变量的预估模型精度降低（R²=0.19~0.41），与加入强度变量构建的模型相比，检验指标降低0.3~0.4。因此，在构建相关预估模型时，综合考虑高度变量和强度变量的作用至关重要，引入强度信息可以显著提升参数估测精度。

为进一步评估Weibull参数以及株数密度的预测效果，生成7个模型的标准化残差图，如图1所示。由图1可以看出，残差各模型残差分布均匀，无明显异方差现象。所建立的回归模型符合独立、正态的基本假设。

2.2 PPM、PPRM、MPRM建立的Weibull直径分布模型拟合效果

根据建立的预估模型，分别采用3种不同方法对直径分布进行了预测，使用E_IR和E_IP对各方法的预测精度进行了评价，见表6。由表6可知，对于整组样地数据，使用PPM所建立的直径分布模型对于直径分布的预测效果最好（E_IR=57.384，E_IP=0.319），而PPRM（式中记为P_PRM）和MPRM（式中记为M_PRM）对于林分的直径分布预测效果与PPM相差不大（P_PRM：E_IR=57.361，E_IP=0.324；M_PRM：E_IR=58.063，E_IP=0.324）。

为进一步分析3种方法对每个龄组中每块样地的直径分布预测效果，绘制了检验指标箱线图，如图2所示。对于E_IR指标来看（图2（a）），PPM在3个龄组中箱体较窄，说明对于每块样地的预测效果和实际分布差别较小，整体预测效果较好。对于E_IP指标来看（图2（b）），PPM在3个龄组中相较于PPRM、MPRM离散点更少，进一步表明PPM对于直径分布预测的精确性。

为直观展示不同方法对直径分布的预测效果，在每个龄组中根据E_IR和E_IP值选择了具有代表性的样地进行绘制预测图，如图3所示。在幼龄林阶段，3种方法均能较好地捕捉到直径分布的峰值位置以及各径阶株数的变化趋势。其中，PPM表现出最佳的预测效果。在中龄林阶段，3种方法都可以捕捉峰值的位置，其中PPM在峰值大小的预测更接近真实值，显示出更高的精度。然而，3种方法都不能很好地预测高径阶部分，存在一定的欠缺。在近熟林阶段，PPM和PPRM都可以很好地预测直径分布峰值以及位置，在对于径阶两端（即高径阶和低径阶部分）PPRM相较于PPM的预测效果存在一定的低估。总体而言，PPM在不同龄组中的预测效果表现均最佳，优于PPM和MPRM。

3 讨论

在传统的森林经营管理中，直径分布模型虽然为森林经营管理者提供了便捷的手段以了解森林生长状况，但传统林业数据调查通常会耗费大量的人力、物力和时间，而且面对一些复杂地形和大面积区域，实施难度较大。相比之下，基于激光雷达技术建立的直径分布模型能够实现对大面积、连续性森林资源的有效管理。本研究将地面调查数据和激光雷达数据相结合，建立直径分布模型，通过E_IR和E_IP值评估了Weibull模型的预测精度，从而验证了利用激光雷达进行大范围预测的可行性。基于全子集法对激光雷达特征因子的筛选结果表明，引入强度变量能有效提升预测精度（R²提高0.3~0.4），这一结果与Goerndt等^［34］的研究结论一致，主要原因在于强度变量与林分年龄存在一定相关性：不同年龄树木的叶片生化特性及其在林冠中的水平和垂直分布差异，导致其反射强度信息有所不同，这种差异通过年龄异质性间接影响株数密度与直径分布，从而导致在强度信息引入模型之后能够有效提升估测精度。而对于常用的高度变量来看，以往研究往往使用的是百分数在99处的高度变量^［15］，而本研究筛选结果则是百分位之1~25的高度变量，导致此结果的原因为红松幼龄林往往是冠下红松，林分中存在着其他的高大乔木从而拔高了林分整体的高度，较低的高度变量可以代表冠下红松的高度，从而对直径分布进行预测。

在对预估模型的残差图（图1）进行分析时发现，ln（N）图的左侧区域出现了一个明显的断层。通过比对实测数据发现，这些位于图左侧的样本主要来自吉林省临江市闹枝林场的样地数据。这可能由于该林场与其他地区（黑龙江省）的样地在地理位置上存在较大差异。通过对临江市林业局档案资料的调查发现，闹枝林场对红松人工林实施了较多的经营管理干预措施，这可能导致吉林省红松的直径分布特征与黑龙江省存在一定差异。残差分析表明，不同地区林场数据的区域不同的差异性可能对模型预测产生干扰。为了解决这一问题，建议引入混合效应模型或贝叶斯模型等^［35-37］高级的数学建模方法，从而提高Weibull参数的预测效果。

在对林分中的直径分布预测结果分析时发现，MPRM对林分直径分布预测效果在3个龄组中的预测效果较差，可能是研究数据中牡丹江市东京城林业局、尚志市苇河林场中红松人工林存在着其他杂木，导致使用激光雷达特征因子对该方法使用的林分因子（林木算术平均直径和平均直径）的预估效果较差。进一步分析发现，利用激光雷达数据对幼龄林尤其是冠下造林的红松幼龄林进行直径分布预测具有较大难度。主要是树冠下激光雷达点云的获取往往受到林分中高大乔木的遮挡导致点云数据精确性下降，从而降低对冠下幼树的探测能力，进而影响直径分布模型的预测精度。针对这一问题，可以通过对林分进行分层预测^［38］，根据林分不同林层建立直径分布模型，引入冠层结构^［39］，使用更加详细的林分冠层特征以及林分中激光雷达特征因子，提高Weibull参数的预测效果；融合地基和无人机激光雷达数据^［40-41］，划分冠层和树种有助于提高林分水平的特征变量精度。

4 结论

本研究基于激光雷达数据，采用了参数预测法、矩估计参数回收法和百分位参数回收法3种方法，建立了红松人工林的直径分布模型。研究结果表明：1）3种方法构建的直径分布模型都可以对红松人工林进行很好的预测，其中PPM建立的直径分布模型在对林分的直径分布预测具有较高的精度（E_IR=57.384，E_IP=0.319）；2）引入强度信息的参数模型其估测精度有较大提升，R²最高可提升0.4。

参考文献

[1] GOBAKKEN T， NÆSSET E.Estimation of diameter and basal area distributions in coniferous forest by means of airborne laser scanner data［J］.Scandinavian Journal of Forest Research，2004，19（6）：529-542.

[2] VALBUENA R， EERIKAINEN K， PACKALEN P，et al.Gini coefficient predictions from airborne lidar remote sensing display the effect of management intensity on forest structure［J］.Ecological Indicators，2016，60（1）：574-585.

[3] MAGNUSSEN S， NÆSSET E， GOBAKKEN T.Prediction of tree-size distributions and inventory variables from cumulants of canopy height distributions［J］.Forestry，2013，86（5）：583-595.

[4] SAARINEN N， KANKARE V， VASTARANTA M，et al.Feasibility of Terrestrial laser scanning for collecting stem volume information from single trees［J］.ISPRS Journal of Photogrammetry & Remote Sensing，2017，123（6）：140-158.

[5] BOLLANDSAS O M， NÆSSET E.Estimating percentile-based diameter distributions in uneven-sized Norway spruce stands using airborne laser scanner data［J］.Scandinavian Journal of Forest Research，2007，22（1）：33-47.

[6] MALTAMO M， SUVANTO A， PACKALÉN P.Comparison of basal area and stem frequency diameter distribution modelling using airborne laser scanner data and calibration estimation［J］.Forest Ecology and Management，2007，247（1/3）：26-34.

[7] BERGSENG E， ØRKA H O， NÆSSET E，et al.Assessing forest inventory information obtained from different inventory approaches and remote sensing data sources［J］.Annals of Forest Science，2015，72（1）：33-45.

[8] COOMES D A， ALLEN R B.Mortality and tree-size distributions in natural mixed-age forests［J］.Journal of Ecology，2010，95（1）：27-40.

[9] LÄMÅS T， SÄNGSTUVALL L， ÖHMAN K，et al.The multi-faceted Swedish Heureka forest decision support system：Context，functionality，design，and 10 years experiences of its use［J］.Frontiers in forests and global change，2023，6（18）：1-19.

[10] BAŞKENT E Z， KELEŞ S.Developing alternative wood harvesting strategies with linear programming in preparing forest management plans［J］.Turkish Journal of Agriculture and Forestry，2006，30（1）：67-79.

[11] JIN X J， PUKKALA T， LI F R，et al.Developing growth models for tree plantations using inadequate data – a case for Korean pine in Northeast China［J］.Silva Fennica，2019，53（4）：10217-10234.

[12] SA Q L， JIN X J， PUKKALA T，et al.Developing Weibull-based diameter distributions for the major coniferous species in Heilongjiang Province，China［J］.Journal of Forestry Research，2023，34（6）：1803-1815.

[13] SUN Y S， JIN X J， PUKKALA T，et al.Two-level optimization approach to tree-level forest planning［J］.Forest Ecosystems，2022，9（1），1125-1147.

[14] PACKALEN P， PUKKALA T， PASCUAL A.Combining spatial and economic criteria in tree-level harvest planning［J］.Forest Ecosystems，2020，7（1）：1-13.

[15] HAO Y S， WIDAGDO F R A， LIU X，et al.A hierarchical region-merging algorithm for 3D segmentation of individual trees using UAV-LiDAR point clouds［J］.IEEE Transactions on Geoscience and Remote Sensing，2022，60：1-16.

[16] RÄTY J， PACKALEN P， KOTIVUORI E，et al.Fusing diameter distributions predicted by an area-based approach and individual-tree detection in coniferous-dominated forests［J］.Canadian Journal of Forest Research，2020，50（2）：113-125.

[17] SHIN J， TEMESGEN H.Generating tree-lists by fusing individual tree detection and nearest neighbor imputation using airborne LiDAR data［J］.Open Journal of Forestry，2018，8（4）：500-531.

[18] GORGOSO-VARELA J J， OGANA F N， IGE P O.A comparison between derivative and numerical optimization methods used for diameter distribution estimation［J］.Scandinavian Journal of Forest Research，2020，35（3/4）：156-164.

[19] GORGOSO-VARELA J J， PONCE R A， RODRÍGUEZ-PUERTA F.Modeling diameter distributions with six probability density functions in Pinus halepensis mill.plantations using low-density airborne laser scanning data in Aragón （Northeast Spain）［J］.Remote Sensing，2021，13（12）：2307-2324.

[20] 叶涛，臧颢.基于气候因子的川中柏木人工林直径分布预测模型［J］.西北林学院学报，2024，39（3）：59-66.

[21] YE T， ZANG Y.Prediction model of diameter distribution of Cupressus funebris plantations in central Sichuan based on climatic factors［J］.Journal of Northwest Forestry University，2024，39（3）：59-66.

[22] 熊子月，周春国.目标树抚育间伐对老山林场天然次生林直径分布的影响［J］.国土与自然资源研究，2024（1）：67-71.

[23] XIONG Z Y， ZHOU C G.Effects of goal tree thinning on diameter distribution of natural secondary forest in Laoshan forest farm［J］.Territory & Natural Resources Study，2024（1）：67-71.

[24] 萨其拉.黑龙江省典型人工针叶林直径分布模型［D］.哈尔滨：东北林业大学，2023.

[25] SA Q L.Diameter distribution models for the typical coniferous plantations in Heilongjiang Province［D］.Harbin：Northeast Forestry University，2023.

[26] LEI Y.Evaluation of three methods for estimating the Weibull distribution parameters of Chinese pine （Pinus tabulaeformis）［J］.Journal of Forest Science，2008，54（12）：566-571.

[27] CAO Q V.Predicting parameters of a Weibull function for modeling diameter distribution［J］.Forest Science，2004，50（5）：682-685.

[28] LIU C， ZHANG S Y， LEI Y，et al.Evaluation of three methods for predicting diameter distributions of black spruce （Picea mariana） plantations in central Canada［J］.Canadian Journal of Forest Research，2004，34（12）：2424-2432.

[29] CAO L， ZHANG Z N， YUN T，et al.Estimating tree volume distributions in subtropical forests using airborne LiDAR data［J］.Remote Sensing，2019，11（1）：97-124.

[30] SIIPILEHTO J.Modelling stand structure in young Scots pine dominated stands［J］.Forest Ecology and Management，2009，257（1）：223-232.

[31] GOBAKKEN T， NÆSSET E.Weibull and percentile models for lidar-based estimation of basal area distribution［J］.Scandinavian Journal of Forest Research，2005，20（6）：490-502.

[32] THOMAS V， OLIVER R D，LIM K，et al.LiDAR and Weibull modeling of diameter and basal area［J］.The Forestry Chronicle，2008，84（6）：866-875.

[33] 郝元朔.基于UAV-LiDAR的落叶松人工林生物量多尺度估算［D］.哈尔滨：东北林业大学，2022.

[34] HAO Y S.Mapping larch （Larix olgensis） plantation aboveground biomass using UAV-LiDAR data［D］.Harbin：Northeast Forestry University，2022.

[35] WANG M，IM J， ZHAO Y H，et al.Multi-platform LiDAR for non-destructive individual aboveground biomass estimation for Changbai larch （Larix olgensis Henry） using a hierarchical bayesian approach［J］.Remote Sensing，2022，14（17）：4361.

[36] PACKALÉN P， MALTAMO M.Estimation of species-specific diameter distributions using airborne laser scanning and aerial photographs［J］.Canadian Journal of Forest Research，2008，38（7）：1750-1760.

[37] REYNOLDS M R， BURK T E， HUANG W C.Goodness-of-fit tests and model selection procedures for diameter distribution models［J］.Forest Science，1988，34（2）：373-399.

[38] GOERNDT M E， MONLEON V J， TEMESGEN H.Relating forest attributes with area- and tree-based light detection and ranging metrics for western Oregon［J］.Western Journal of Applied Forestry，2010，25（3）：105-111.

[39] TIAN X L， SUN S C， MOLA-YUDEGO B，et al.Predicting individual tree growth using stand-level simulation，diameter distribution，and Bayesian calibration［J］.Annals of Forest Science，2020，77（2）：57-74.

[40] HAO Y S， WIDAGDO F R A， LIU X，et al.Estimation and calibration of stem diameter distribution using UAV laser scanning data：A case study for larch （Larix olgensis） forests in Northeast China［J］.Remote Sensing of Environment，2022，268：1127-1169.

[41] PEUHKURINEN J， TOKOLA T， PLEVAK K，et al.Predicting tree diameter distributions from airborne laser scanning，SPOT 5 satellite，and field sample data in the perm region，Russia［J］.Forests，2018，9（10）：639-654.

[42] 王金池，严铭海，黄清麟，等.基于分层的半天然马尾松阔叶混交林直径分布［J］.林业科学研究，2021，34（3）：72-80.

[43] WANG J C， YAN M H， HUANG Q L，et al.Diameter distribution of semi-natural mixed forest of Pinus massoniana and broadleaved trees based on stratification［J］.Forest Research，2021，34（3）：72-80.

[44] GALLAGHER-DUVAL X， VAN LIER O R， FOURNIER R A.Estimating stem diameter distributions with airborne laser scanning metrics and derived canopy surface texture metrics［J］.Forests，2023，14（2）：287-305.

[45] 张丛凯，于颖.ICESat-2/ATLAS 数据地面高程及植被冠层高度反演精度验证［J］.森林工程，2023，39（3）：1-11.

[46] ZHANG C K， YU Y.Accuracy verification of ground elevation and vegetation canopy height inversion from ICESat-2 / ATLAS data［J］.Forest Engineering，2023，39（3）：1-11.

[47] 黄旭，贾炜玮，王强，等.背包式激光雷达的落叶松单木因子提取［J］.森林工程，2019，35（4）：14-21.

[48] HUANG X， JIA W W， WANG Q，et al.Study on individual tree factor extraction of Larix olgensis in backpack lidar［J］.Forest Engineering，2019，35（4）：14-21.

基金资助

国家自然基金面上项目(32071758)