论文标题：Estimation of aboveground biomass of Robinia pseudoacacia forest in the Yellow River Delta based on UAV and Backpack LiDAR point clouds

ABSTRACT

人工林是陆地碳汇的重要来源。黄河三角洲刺槐林是我国最大的人工生态防护林。然而，自1990年代以来，超过一半的森林出现了不同程度的枯死甚至死亡。及时准确地估算森林地上生物量（AGB）是研究森林碳循环的基础。光探测和测距 (LiDAR) 已被证明是森林生物量估算最有效的方法之一。然而，由于生长季节阔叶林冠层的形状不规则且重叠，因此很难分割单个树木并根据机载 LiDAR 数据估计树木生物量。在这项研究中，提出了一种基于无人机-光探测和测距（UAVLiDAR）与背包-LiDAR 相结合的新方法来解决刺槐林中的个体树木检测问题。所提出的方法主要包括以下步骤：（i）在地块级别，UAV-LiDAR 数据中的树木通过从 Backpack-LiDAR 数据中通过个体树分割（ITS）方法获得的种子点进行检测； (ii) 单棵树的高度和胸径 (DBH) 将分别从无人机和背包 LiDAR 数据中提取； (iii) 个体树的 AGB 将通过异速生长方程计算，并累积地块级的森林 AGB； (iv) 以地块级森林 AGB 为因变量，从 UAV-LiDAR 点云数据中提取的各种指标作为自变量，利用多元线性回归（MLR）和随机森林 (RF) 模型。结果表明：（1）Backpack-LiDAR提取的种子点可以显着提高个体树检测的整体精度（F = 0.99），从而提高森林AGB估计精度； (2) 与MLR模型相比，RF模型的估计精度更高(p < 0.05)； (3) MLR 和 RF 模型选择的 LiDAR 强度信息和激光穿透率 (LP) 在估计健康森林 AGB 中发挥了重要作用。

1. Introduction

人工造林被认为是通过吸收 CO2 和减缓气候变暖来增加碳封存的最生态有效的方法之一（Piao 等人，2009 年）。中国广泛开展了森林保护、植树造林等生态工程。中国人工林面积占全球人工林面积的73%，已成为增加中国陆地固碳的重要手段(Zhou et al, 2012)。黄河三角洲 (YRD) 拥有中国面积最大的人工刺槐林。由于土壤肥力低，土壤盐渍化普遍，长三角地区没有天然林。刺槐林具有适应性强、耐旱和一定耐盐性的特点，从20世纪70年代开始大量种植。然而，近60%的刺槐林出现了不同程度的枯死甚至死亡（Wang等，2015a）。及时准确估算刺槐林地上生物量（AGB）可为评估长三角地区人工林碳汇提供科学依据。

通过大空间尺度的野外调查测量森林生物量费时费力，难以推广（Hermosilla 等，2014；Van Leeuwen 和 Nieuwenhuis，2010）。光探测和测距 (LiDAR) 等遥感技术已证明其具有提供森林冠层结构的三维详细特征的潜力（Lefsky 等人，2002 年；Næsset 和 Gobakken，2008 年）。 LiDAR数据可以提取亚米级垂直精度的森林高度和水平分布信息，这在森林AGB估计中具有显着优势（e.g., Cao et al, 2016; Hudak et al, 2012）。

个体树分割 (ITS) 在 LiDAR 数据的森林 AGB 估计中具有重要意义。一旦树木被准确分割，树高和树冠直径等树木结构参数就可以直接高精度地提取出来(Solberg et al, 2006)。然后可以使用异速生长方程准确估算森林 AGB（Basuki 等，2009；Wang，2006）。 LiDAR数据衍生的ITS要么直接使用点云数据，要么使用从点云数据衍生的冠层高度模型（CHM）（Liu et al, 2019; Zhen et al, 2016）。基于 CHM 的分割仅使用激光点云的第一个回波，对于林下树木的检测是不切实际的。因此，基于点的分割近年来得到了迅速发展。例如，Reitberger 等人 (2009) 使用随机样本一致性算法来检测单个树干，准确率接近 70%。 Li等(2012)根据针叶树冠间距大于林下植被的特点，以树冠顶部激光点云为种子点，采用自上而下的区域生长法分离单株。 , 总准确率达到 90%。但是，当此方法应用于阔叶林时，准确性会降低。由于阔叶林冠层形状不规则，枝叶交叉重叠，点云数据难以从冠顶确定种子点。 Lu et al (2014)利用激光雷达在落叶季获取的点云数据的强度信息，提取树干与点云之间的拓扑关系，对树干进行分割，总准确率达到90%。但这种方法只能在落叶季节应用。这是因为在生长期，叶片强度值也很大，干扰了树干的提取过程。

基于LiDAR数据的区域AGB估计通常是通过在LiDAR样本数据中提取的特征变量与森林清查属性之间建立模型来获得的。从 LiDAR 点云数据中提取的高度和密度变量已被证明与森林生物量密切相关（Hall 等人，2005 年；Næsset 和 Gobakken，2008 年）。然而，仅使用森林高度和密度信息的 LiDAR 指标不足以描述整体冠层和垂直方向的异质性，因此 Zhang 等（2017）增加了冠层高度分布（Weishampel 等，2007；Zhao 等，2009） )和枝叶剖面(Lovell et al, 2003)，描述了冠层剖面的特征变量，他们使用Weibull函数拟合参数来推断阔叶林AGB，结果很好(R2 = 0.66, RMSE = 26.67 Mg/ha). 虽然激光点云的第一个回波提供了上层冠层结构的信息，但最后一个回波分布描述了植被层中的最大穿透激光信号，可以区分不同退化程度的森林类型，提高生物量的估计准确性（Ioki 等，2014）。由于单光谱 LiDAR 传感器通常使用 1064 nm 近红外波长，这非常适合识别植物反射率的变化，因此健康的树木通常具有很强的反向散射（Lorenzen 和 Jensen，1988）。 Yoga 等人 (2017) 使用从 LiDAR 中提取的遥感图像和强度信息来消除随机森林模型分类识别的死树，并提高了森林存量的估计精度。因此，从LiDAR数据中提取的高度度量、密度度量、剖面特征度量、最后回波透射率和强度信息都被证明是森林AGB估计的有用特征变量。

基于LiDAR数据的区域AGB估计通常是通过在LiDAR样本数据中提取的特征变量与森林清查属性之间建立模型来获得的。从 LiDAR 点云数据中提取的高度和密度变量已被证明与森林生物量密切相关（Hall 等人，2005 年；Næsset 和 Gobakken，2008 年）。然而，仅使用森林高度和密度信息的 LiDAR 指标不足以描述整体冠层和垂直方向的异质性，因此 Zhang 等（2017）增加了冠层高度分布（Weishampel 等，2007；Zhao 等，2009） )和枝叶剖面(Lovell et al, 2003)，描述了冠层剖面的特征变量，他们使用Weibull函数拟合参数来推断阔叶林AGB，结果很好(R2 = 0.66, RMSE = 26.67 Mg/ha). 虽然激光点云的第一个回波提供了上层冠层结构的信息，但最后一个回波分布描述了植被层中的最大穿透激光信号，可以区分不同退化程度的森林类型，提高生物量的估计准确性（Ioki 等，2014）。由于单光谱 LiDAR 传感器通常使用 1064 nm 近红外波长，这非常适合识别植物反射率的变化，因此健康的树木通常具有很强的反向散射（Lorenzen 和 Jensen，1988）。 Yoga 等人 (2017) 使用从 LiDAR 中提取的遥感图像和强度信息来消除随机森林模型分类识别的死树，并提高了森林存量的估计精度。因此，从LiDAR数据中提取的高度度量、密度度量、剖面特征度量、最后回波透射率和强度信息都被证明是森林AGB估计的有用特征变量。

2. Study area and data sets

2.1. Study area

长三角位于中国山东省东营市黄河入海口（图 1a）。属暖温带大陆性季风气候，年平均气温11.7～12.6℃，年平均降水量530～630毫米。土壤盐分对树木的正常生长有负面影响（Zhang，2013）。刺槐是主要的适宜造林树种，自 20 世纪 70 年代以来已在该研究区广泛种植，形成了中国最大的人工林（Wang 等，2018）。长三角地区共有 4 个林区（图 1b），总面积 27.94 平方公里（Wang et al, 2015a）。在本研究中，孤岛林被选为我们的研究区域（图 1c）。

2.2. Datasets

2.2.1. UAV-LiDAR data

2017年6月，实施GreenValley（GreenValley, International, USA, 2019）激光雷达系统，采集激光雷达数据。以八旋翼无人机为平台，飞行高度为120 m，速度为4.8 m·s−1，飞行半径为2 km。无人机平台上安装了 Velodyne Puck VLP-16 双回波激光扫描仪、IMU（Novate）和双频 GPS（Novatel）。 LiDAR 系统配置为发射近红外波段的激光脉冲，扫描角度为距最低点 ± 30°；激光发散角为0.5 mrad；光斑直径约50mm，一个激光脉冲最多有4个回波，平均点密度70m-2，测距精度10mm。采用WGS84坐标系和UTM投影。一般来说，强度归一化可以提高森林属性的估计精度，但这种提高非常小，因此我们没有进行强度归一化(You et al, 2017)。表 1 显示了激光返回密度和强度的汇总。

图 1 研究地点及样地分布。 (a) 东营市黄河三角洲的位置。 (b) 使用 2013 年 6 月 11 日获取的 Landsat 8 OLI 图像作为背景的四个刺槐林的分布。 (c ) 以 2013 年 6 月 9 日获取的 IKONOS 图像为背景，位于孤岛森林三个不同健康水平的采样地块（有关此图图例中颜色参考的解释，读者可以参考网络版本文）。

表格1 研究地点的激光返回密度和强度摘要。

2.2.2. Backpack-LiDAR data

GreenValley Backpack LiDAR 系统由激光雷达扫描仪 (Velodyne Puck VLP-16)、定位系统 (POS) 和手持触摸板组成，在本研究中得到应用。在孤岛森林中，Backpack LiDAR 系统扫描了八个 30 m×30 m 的样地（图 1），具有三个健康水平（即健康、中度枯死、重度枯死）（Wang et al, 2015b）。我们设计了一条“S”形条带路径，并在样地中心放置了一个已知坐标的人工标记（图 2）。 Backpack LiDAR激光扫描距离为100米；扫描频率为300,000 pts·s-1；激光波长为903nm；平均点密度为7135 m-2；水平视角360°，垂直视角15°。数据也采用WGS84坐标系和UTM投影。

2.2.3. Field data

正如我们之前的研究 (Wang et al, 2015a)，在每个 30 m × 30 m的地块中，五个10 m × 10 m的子地块被部署在四个角，每个地块的中心有一个。在每个子图中，选择一棵标准树。 2017年6月，天宝GEOXT6000 GPS定位仪记录了各样地和子样地中一棵标准树的地理坐标。分别使用激光测高仪和卷尺测量了40棵标准树的树高和胸径。

图 2. Backpack-LiDAR 点云在一个样地（30 m × 30 m）内采集的轨迹。

3. Methods

图 3 显示了估计 AGB 的工作流程概览，主要分为两部分：（1）树木分割和观察到的 AGB 计算，以及（2）基于 UAV-LiDAR 的变量提取和估计模型比较。我们首先对无人机和背包激光雷达点云进行匹配和归一化。其次，使用比较最短路径算法（CSP）（Tao et al，2015）对Backpack-LiDAR点云进行分割以获得种子点和单树DBH。使用这些种子点，检测个体树木并从无人机激光雷达点云中获得树高，然后计算并累积在八个样地观测到的 AGB。第二部分提取LiDAR指标后，比较了MLR和RF模型对森林AGB的预测结果。下面给出了有关此方法的更多详细信息。

图 3. 使用 LiDAR 数据估算森林地上生物量的工作流程概述。

3.1. Matching between Backpack-LiDAR and UAV-LiDAR

由于派生的 Backpack-LiDAR 点云对尺度方差和平面或垂直偏差敏感（图 4（a）），我们首先使用改进的渐进三角不规则网络致密化过滤算法生成数字地形模型（DTM）（Zhao et al, 2016)，然后我们使用反距离加权算法对剩余的地面点进行插值。通过减去地面高度对点云高度进行归一化。归一化后，Backpack-LiDAR 和 UAV-LiDAR 数据的地面点位于同一平面上（图 4（b））。然后我们使用 Polewski 等人（2019）建议的方法来执行两个数据匹配。至少在 Backpack-LiDAR（图 4（c））和 UAVLiDAR（图 4（d））的范围内手动选择了至少三对同名点。最后对匹配后的数据再次归一化，消除点云匹配中z值差异带来的影响。表2显示了两类点云数据的配准精度。

3.2. Tree segmentation of LiDAR data

3.2.1. Tree segmentation of Backpack-LiDAR

我们应用 GreenValley LiDAR360 商业软件 (LiDAR360, 2018) 对无人机和 Backpack 获取的数据进行预处理。基于密度的应用噪声空间聚类算法 (Wu et al, 2013) 用于分割单个树干。提取了垂直长度为 10 cm、高度为 1.3 m 的切片并将其用作 CSP 算法的输入（Tao 等人，2015），然后计算了 DBH（图 5）。去除激光雷达点云中的现场工作人员、倒下的树木和参考杆等噪声数据后，获取包含每棵树的 X、Y 坐标和 DBH 值的最终种子点文件，并用于从中分割单棵树UAV-LiDAR 点云。

3.2.2. Tree segmentation of UAV-LiDAR

去噪后，使用 3.1 节中相同的方法再次对 UAV-LiDAR 点云进行归一化。将点云分割 (PCS) 方法（Li 等人，2012）应用于分割树冠。然而，在本研究中，并未将树冠的局部最大点作为树冠提取的树尖（种子点），而是将Backpack-LiDAR提取的树干的X和Y坐标作为种子点。由于 Backpack-LiDAR 的数据质量高，我们将匹配和归一化的 Backpack-LiDAR 数据设置为树分割精度评估的参考。三个统计参数（Goutte 和 Gaussier，2005），分别是树木的检测率 r（“召回率”）、检测到的树木的检测准确率 p（“精度”）和整体准确率 F（F 分数） , 用于评估 ITS 算法的性能。三个参数定义如下：

其中，TP 是图中检测到的树木数量； FN 是单个树木分割遗漏的树木数量，FP 是图中错误检测到的树木数量。同时，分别用Backpack-LiDAR和UAV-LiDAR提取的子图中每棵标准树的胸径和高度值，以及相应的实地测量值进行了评估。

图 4.背包-激光雷达和无人机-激光雷达匹配示意图。(a)原始背包-激光雷达(黑色)和UAV-激光雷达(彩色)点云的叠加图，(b)归一化背包-激光雷达和UAV-激光雷达点云的叠加图，©和(d)分别在UAV-激光雷达和背包-激光雷达点云中标记相同树点云的编号。

表2 八个地块(H:健康地块，M:中度枯死地块，S:重度枯死地块)的配准精度。

图 5. 1.3 m 树干切片的演示结果和计算的 DBH。

图 6 单地块分割示意图。

3.3. Biomass calculation of Robinia pseudoacacia

将不同健康水平的8个30 m×30 m样地对应的LiDAR数据进行切割，每个样地切割成9个10 m×10 m的样地，共72个样地（图6）。在每个子样地统计UAV-LiDAR和Backpack-LiDAR提取的树高和胸径，利用国家林业局公布的刺槐林异速生长方程（见公式（4-6））计算计算树干 (WS)、树枝 (WB) 和叶 (WL) 的生物量，然后将所有子图的森林 AGB 总结在表 3 中。

其中 D 是胸径 (cm)，H 是树高 (m)。

表 3 具有三个不同健康水平的 72 个子地块的实地估计森林特征摘要。

表 4 LiDAR 指标摘要。

3.4. UAV-LiDAR metrics

根据之前的研究（Lim 等人，2003 年；Næsset，2002 年；Næsset 和 Gobakken，2008 年；Ioki 等人，2014 年；Thomas 等人，2006 年），本研究中从归一化 UAV-LiDAR 点云数据中提取的指标包括高度-相关指标、密度相关指标、从所有回波中提取的强度相关指标以及从最后一个回波中提取的 LP（表 4）。高度相关指标描述了与点云高度（H1、H5、H25、H50、H75、H95、H99）相关的高度百分位数，例如高度均值（Hmean）、高度变异系数（Hcv）和高度方差（Hvariance）等。与密度相关的指标描述了冠层回归密度的比例，即百分位数以上的点云数量与总点数的比值。强度相关指标与高度相关指标类似，将点的高度值替换为点的强度值进行计算。作为表征植被冠层结构的基本参数，叶面积指数（LAI）被定义为每单位表面积所有叶片表面积的一半（Chen和Black，1991）。由于刺槐树的枯死是从顶冠开始的，因此健康林和枯死林的 LAI 是不同的。 LAI 指标是根据比尔朗伯定律推导出来的（Richardson 等人，2009 年）：

其中 ang 是平均扫描角度； GF 是间隙分数； k 是消光系数。在这项研究中，我们假设叶角的分布是球形的并且 k = 0.5 (Richardson et al, 2009)。 anglei 为无人机las文件中记录的第i个激光雷达点的扫描角度； nground 是地面点的数量，n 是 LiDAR 点的数量。通过使用 LiDAR360 软件从 UAV-LiDAR 中提取的总共 53 个指标（表 4）通过 MLR 和 RF 模型准备用于孤岛森林中的 AGB 预测。

3.5. Model development

不同的建模方法对结果质量有不同的影响（Straub 等人，2010 年）。在这项研究中，参数和非参数方法进行了比较。 MLR 和 RF 回归模型都用于估计每个子图的森林 AGB，并分别使用测量的森林 AGB 和 LiDAR 数据导出指标作为因变量和自变量运行。 MLR 是一种非常简单的参数方法，能够处理预测变量的依赖性或相关性，并且经常用于 AGB 估计（例如，Fassnacht 等人，2014 年；Morin 等人，2019 年）。以前的研究已经使用因变量和自变量的对数变换来提高模型的拟合能力（Næsset 等，2005）。由于我们的 LiDAR 指标中存在负数，因此未进行对数变换以避免信息丢失。在采用 MLR 模型之前，使用 WEKA 软件（Hall 等人，2009 年）来选择最佳的 LiDAR 衍生指标。 WEKA 中的 CfsSubsetEval 评估器评估每个属性（即本研究中的指标）的预测能力及其相互冗余，并倾向于选择与目标属性（即本研究中的因变量）高度相关但与每个属性相关性较低的属性其他（霍尔，1998 年）。我们采用从空属性集开始的前向搜索来过滤属性子集。在选择最优候选指标后，进行MLR以获得最优AGB模型。留一法交叉验证 (LOOCV) 方法（Bengio 和 Grandvalet，2004 年）使用所有其他样本估计这些样本中的每一个，用于评估校准 MLR 模型的性能。

RF 是一种基于决策树的无分布（非参数）分类算法，可以避免过度拟合问题，并且对异常值和噪声具有鲁棒性（Breiman，2001）。对于RF，需要设置两个参数ntree和mtry。 ntree 表示回归模型中运行的树的总数，mtry 表示树的每个节点上可以拆分的变量数 (Mutanga et al, 2012; Yu et al, 2011)。在本研究中，为了获得更好的ntree和mtry用于预测不同健康水平的森林AGB，我们根据误差分布和解释率对两个参数进行了优化，ntree取值为1000，mtry取值为三分之一的预测变量被收购。通过计算去除该变量后模型均方误差的增加，对每个指标的相对重要性进行排序。精度评估也通过使用 LOOCV 方法进行。

回归模型的准确性通过决定系数（R2）、均方根误差（RMSE）和相对均方根误差（rRMSE）来评估，表示如下：

其中 xi 是图 i 的测量值； x^i 是地块 i 的估计值； x̄i 是子图 i 的观测平均值； x̄ 是所有地块的观察平均值； n 是子图的数量。

6. Conclusions

本研究利用UAV-LiDAR和Backpack-LiDAR两种点云数据对刺槐林不同健康水平下的单株检测和AGB估计进行了研究。实验结果得出以下三个结论：

(1) 从 Backpack-LiDAR 数据中提取的树干位置可以用作使用 UAV-LiDAR 数据进行个体树分割 (ITS) 的种子点。该方法可有效提高生长季阔叶林树木分割精度（总精度F=0.99）。随着激光雷达技术的发展，激光探测器和承载平台变得更加轻量化和多样化。无人机与Backpack LiDAR数据的结合可以有效解决阔叶林中使用ITS方法分割不足或过度的问题。

(2) 反映植被健康状况的冠层覆盖度、LP 和强度相关指标在刺槐林的AGB 估计中具有重要作用。

(3) RF 模型预测森林 AGB 的准确性高于 MLR 模型。