基于稀疏型机载激光雷达数据的风景林参数估测

魏金龙，李明阳，赵邑晨，李 超，李盈昌

(南京林业大学 南方现代林业协同创新中心，江苏 南京 210037)

森林是陆地生态系统的主体，它具有净化空气、改善局部小气候、保持水土、观赏游憩和碳汇等功能，保护和经营管理具有一定景观质量和游憩功能的森林是当前森林质量提升工程的一个研究重点。随着森林康养游憩产业的发展，森林公园、风景区内风景林的参数提取成为森林资源调查与生态功能监测的重要内容。森林参数能有效反映森林资源的空间分布和结构特性[1]，是森林资源规划设计调查的重要内容。传统森林参数测定需要外业实地测量，不仅费时费力，而且只能获得小范围林区的信息，不利于研究大范围的森林参数。遥感作为一门新兴的应用科学，已经在许多行业发挥了重要的作用，在林业领域也不例外。激光雷达(light detection and ranging,LiDAR)作为主动遥感技术可以有效穿透森林，在获取森林垂直结构参数方面有着其他光学遥感无法比拟的优势[2]。获取高精度的森林结构参数，对于了解森林结构特征具有重要作用，机载激光雷达技术(airborn laser scanning,ALS)的出现使得提取高精度森林参数成为可能，同时它还能够提供森林冠层上表面的水平分布和垂直结构信息等[3]。庞勇等[4]从机载激光雷达数据归一化的植被点上提取四分位数处的高度，并结合实测树高估测林分树高，结果表明，使用分位数法能很好的估测林分树高，精度达到90.59%；曹林等[5]利用机载激光雷达全波形数据估测亚热带森林林分平均高，结果显示估测精度决定系数(R2)为0.85，均方根误差(RMSE)是实测平均高的7.26%；高婷等[6]基于机载激光雷达数据估测林分水平上森林参数，林分平均高估测精度最高(R2=0.800 7)，平均胸径估测精度仅次于平均树高；鲁林等[7]从机载激光雷达数据中提取点云统计特征和高度变量使用随机森林算法估测林分平均高，总体精度达到93.17%；J.Jensen等[8]从机载激光雷达数据提取高度参数估测林分平均胸径，R2=0.66，RMSE为6.31 cm；焦义涛等[9]基于机载LiDAR高程归一化后的植被点云数据建立林分平均树高线性回归模型，模型估测平均树高精度最高为99.81%；C.A.Silva等[10]采用主成分分析法选择LiDAR特征变量估算了热带桉树人工林林分蓄积，结果表明99%高度百分位数和变异系数对桉树人工林蓄积量的估测精度R2=0.84。

机载LiDAR点云密度决定了林木冠层结构描述的详细程度，高密度点云(点云密度>3个/m2)，通过一定拟合算法可以估测单木尺度森林参数[11-12]。稀疏型点云(点云密度<3个/m2)，受到点云数量的限制，无法精细刻画单木完整的树冠表面结构，多利用LiDAR数据高度和密度特征变量估测较大尺度森林参数[13-15]。E.Gonzalez-Olabarria[16]等利用机载LiDAR技术估测点云密度为8个/m2和0.5个/m2的林分平均树高、优势树高、林分蓄积量和生物量，结果表明点云密度稀疏程度对林分参数反演精度影响较小。S.G.Tesfamichael等[17]利用高密度点云数据(8个/m2)分别生成点云密度为0.25、1、2、3、4、5、6个/m2的数据，并使用多元线性回归模型(multiple linear regression)估测林分平均树高和优势高，结果显示，不同点云密度下林分平均树高的R2=0.93～0.94，林分优势高的R2均为0.95。P.Treitz等[18]利用点云密度分别为3.2、1.6个/m2和0.5个/m2的机载LiDAR数据分别提取了高度和密度特征变量，结合多元线性回归模型估测林分平均树高、优势树高和总地上生物量等，结果表明，林分平均树高估测精度R2=0.75～0.95，林分优势树高R2=0.74～0.98，总地上生物量R2=0.23～0.93，并且点云密度为0.5个/m2的机载LiDAR数据足以进行绘图和林分水平森林参数的估测。A.L.Montealegre等[13]利用点云密度为1个/m2的机载LiDAR数据结合多元线性回归模型对地中海松(Pinushalepensis)林分平均树高、平均胸径和单位蓄积量等进行了估测，结果显示，估测精度R2都在0.80以上，低密度机载LiDAR点云数据精确估测地中海松林分特征参数是有效的。

总结已有研究，国内外对高密度机载LiDAR点云数据研究的理论技术逐渐成熟并不断深入，研究的方法不断改进和优化，稀疏型机载LiDAR点云数据，由于点云密度低，难以对单木尺度的森林参数进行估测，在森林资源监测领域的潜在应用价值没有得到充分发挥。针对此问题，本研究充分考虑稀疏型机载LiDAR数据所提取的高度和密度等特征因子的特殊性，进行小班尺度森林参数估测的精度评价研究。本研究以南京紫金山国家森林公园为研究区，以2007年机载激光雷达Optech ALTM LiDAR数据、2007年森林资源规划设计调查数据为主要数据源，主要探讨：1)稀疏型机载LiDAR数据在小班尺度森林参数提取的应用；2)影响估测小班尺度森林参数的主要机载LiDAR数据特征变量；3)对比分析参数和非参数模型在小班尺度森林参数估测中的精度。本研究利用稀疏型机载LiDAR数据和森林资源规划设计调查数据探讨稀疏型机载激光雷达在林分尺度上的森林参数估测，可以为稀疏型机载LiDAR数据在林分尺度森林参数遥感估测中的应用提供科学参考。

1 材料与方法

1.1 研究区概况

南京紫金山国家森林公园位于南京市玄武区，地理坐标为32°01′57″-32°06′15″N，118°48′-118°53′04″E，是全国生态文明示范地、首批国家5A级旅游景区、中国首个被评为“国家级森林公园”的大型城市森林公园。紫金山国家森林公园辖区总面积3 008.8 hm2，森林覆盖率达到70.2%，森林气氛浓郁，林相整齐，生物多样性保存完好。植被类型可分为阔叶林、针叶林、针阔混交林、竹林四大类，现有种子植物1138种，其中木本植物499种、草本植物569种、藤本植物70种，蕨类植物25科36属78种，列入特、稀、危植物名录的50余种，古树名木924株。研究区常绿阔叶树有青冈(Cyclobalanopsisglauca)、石楠(Photiniaserratifolia)、冬青(Ilexchinensis)等；落叶树木有黄连木(Pistaciachinensis)、枫香(Liquidambarformosana)等，具有重要的观赏与科研价值。紫金山国家森林公园地处中纬度地区，属亚热带季风气候，区内降水丰富，四季分明，气候温和。研究区的地带性土壤为黄棕壤，在山北坡下部存在少量的紫色土壤。研究区年平均日照1 628.8 h，平均气温19.6℃；年平均降水量1 530.1 mm，无霜期322 d[19](图1)。

图1 研究区地理位置Fig.1 Geographic location map of the research area

1.2 数据来源

本研究所采用的数据为2007年紫金山国家森林公园森林资源规划设计调查632个有林地小班数据和2007年紫金山国家森林公园机载激光雷达点云数据。2007年紫金山国家森林公园森林资源规划设计调查数据，每个矢量化小班属性表包括海拔、坡度、坡位、林种、树种、平均树高、平均胸径、郁闭度、单位蓄积量等50多个属性因子。

2007年9月20日的机载激光雷达Optech ALTM LiDAR数据，根据航线设计软件中所确定的飞行高度为900 m、扫描角度为±25°、扫描频率(scanning frequency,SF)为70 s-1和激光脉冲发射频率(pulse reflection frequency,PRF)为100 kHz，平均地面光斑距离2 m，点云密度为0.25个点/m2。沿航线逐条采集地表三维坐标信息，再通过特定算法处理生成高精度的地表模型DTM。其设备参数见表2。

表1 研究区森林参数统计Table 1 The statistics of forest parameters of the study area

1.3 数据的处理和提取特征参数

利用激光雷达数据进行森林参数反演，需要对机载激光雷达数据进行预处理。首先去除机载LiDAR原始点云数据的噪音点，使用LiDAR 360软件对原始数据去噪，所采用的算法将对每一个点搜索指定邻域点个数(为10)的相邻点，计算点到相邻点的距离平均值，计算这些平均值距离的中值和标准差，如果这个点的平均值距离大于最大距离(最大距离=中值+5倍标准差)，则认为是噪点，将被去掉；然后采用改进的渐进加密三角网滤波算法去除非地面点，将分类后的地面点数据使用不规则三角网(TIN)插值方法生成数字高程模型(DEM)；最后用DEM的高程值对植被回波点的高度进行归一化处理，即去除地形的高程，使植被点的高度值为相对于地面的高度值。

高度百分位数和密度特征变量能够很好地反演大多数森林结构参数[13-14,20]。因此，对紫金山机载LiDAR数据按小班单位提取高度变量，包括累积高度百分位数(AIH，15个)、峰度、中位数绝对偏差的中位数、最大值、最小值、平均值、中位数、二次幂平均、三次幂平均和高度百分位数(15个)以及10个密度变量(表3)。

表2 激光雷达设备参数Table 2 Parameters for airborne laser radar equipment

表3 LiDAR特征变量汇总Table 3 Summary of LiDAR characteristic variables

1.4 模型对比及精度评价

采用参数回归方法(多元线性回归(multiple linear regression))和支持向量机(support vector machine,SVM)、随机森林(random forest)2种非参数回归方法对紫金山国家森林公园林分内平均胸径、平均树高及单位蓄积量进行估测。

多元线性回归根据已知数据以多元线性回归分析法为理论基础，以求解求逆紧凑变换法为算法技巧，通过双检验法建立最优回归方程[21]。支持向量机算法根据样本数据在模型复杂性和学习能力之间寻求最佳值，以获得最好的结果。支持向量机模型中核函数的选择对估测精度的影响较大，常用的有线性核函数(line)、多项式核函数(polynomial)、径向基核函数(radius basis function,RBF)和多层感知机核函数(sigmoid)[22]。随机森林算法(random forest)是Breiman提出的一种基于分类与回归决策树的组合算法，随机森林利用多棵树对样本进行训练并预测，用Bagging方法生成具有独立同分布的训练样本集进行训练，生成的树的参数随机向量也是独立同分布的；随机森林中决策树的数量(ntree)和随机特征的数量(mtry)是2个非常重要的参数[23-24]。

本研究基于紫金山国家森林公园森林资源规划设计调查数据，采用10折交叉验证的方法评价估测模型精度。10折交叉验证将2007年紫金山国家森林公园森林资源规划设计调查小班数据随机分成10份，轮流将其中9份作为训练集，1份作为测试集进行估测。通过计算决定系数(R2)、均方根误差(RMSE)和相对均方根误差(rRMSE)3个常见指标10次结果的均值对估测模型的精度进行评价。3个指标的公式如下：

1)决定系数(R2)：R2值越大，则估测值与实测值之间的相关性越强。计算公式为：

(1)

2)均方根误差(RMSE)：RMSE越小，则表明模型预测的效果越好。计算公式为：

(2)

3)相对均方根误差(rRMSE)：rRMSE越小，表明模型预测的效果越好。计算公式为：

(3)

1.5 基于最优模型的森林参数反演

通过10折交叉验证选择每个森林参数估测精度最高的模型，利用R语言软件分别进行平均胸径、平均树高和单位蓄积量的反演制图，并与实测数据对比分析。

2 结果与分析

2.1 基于最优模型特征变量重要性分析

用R语言软件对3个森林参数特征变量的相对重要性和节点纯度进行分析，挑选对3个森林参数影响最大的10个变量对森林参数进行建模(表4)。由表4可知，高度百分位数和累计高度百分位数是影响小班平均胸径、平均树高和平均蓄积量的主要特征参数。percentile80%、percentile90%、percentile99%对3个森林参数均有较大影响。影响单位蓄积量的主要特征变量均影响平均树高和平均胸径，原因在于紫金山国家森林公园森林资源规划设计调查数据中单位蓄积量是由平均实验形数法计算出来，林分蓄积与平均胸径和平均树高高度相关。需要指出的是，3个森林参数的主要特征变量虽有相同的特征变量，但仍存在不同的特征变量，如除percentile80%、percentile90%、percentile99%外，影响平均树高的有average，sq，cube等，影响平均胸径的有AIH95%、madmedian等，影响单位蓄积量的有max、AIH1%、AIH99%等。

表4 影响3个森林参数估测的主要特征变量Table 4 Main characteristic variables affecting estimation of three forest parameters

2.2 基于3种模型反演森林参数精度评价

通过挑选出来的对3个森林参数影响最大的10个机载LiDAR特征变量，进行小班平均树高、平均胸径、单位蓄积量遥感估测。在R语言软件中不断调整模型参数以获取最优精度，最后得出：3个森林参数的随机森林模型均在决策树的数量为40、随机特征的数量为16时能达到最好的预测效果；在支持向量机模型中，3个森林参数估测时核函数均选择径向基核函数(radius basis function,RBF)、gamma取值为0.05时能达到最好的预测效果；通过逐步回归筛选变量，构建多元线性回归模型，在变量共线性诊断的基础上得到3个森林参数的多元线性回归模型(表5)。

采用10折交叉验证方法对3个回归模型进行精度验证(表6、图2)。从表6及图2可知，3种回归模型估测林分平均树高、平均胸径和单位蓄积量时，均是随机森林模型精度最高，支持向量机模型次之，多元线性回归模型最低。多元线性回归模型估测精度最低的原因是研究区海拔相对高差大，坡度、坡位变化复杂，同时地处中纬度地区，属亚热带季风气候，降水丰富，树高、胸径和蓄积量与激光雷达提取的特征变量不是简单的线性关系。非参数回归方法以严谨的分析和操作的灵活性，对3个森林参数的估测精度高于参数回归方法。支持向量机的关键

表5 3个森林参数多元线性回归模型Table 5 Stepwise regression models for three forest parameters

在于找到一个最优超平面对样本进行抽样，本研究选用径向基核函数以非线性关系对样本抽样，精度高于多元线性回归模型，但支持向量机对于大样本数据估测精度不是很好。随机森林不仅对样本进行抽样，还对特征变量进行抽样，生成多个树模型，再进行组合预测，估测精度高于支持向量机模型。因此，非参数回归模型精度高于参数回归模型，其中随机森林模型精度最高。

2.3 基于最优模型3个森林参数空间分布特征分析

通过上述3种估测模型对于林分平均高、平均胸径和单位蓄积量的反演精度评价对比，得出随机森林算法反演该3个森林参数精度最高(R2均>0.86)。因此，选择随机森林算法为最优估测模型以30 m×30 m的分辨率反演紫金山国家森林公园3个森林参数，并制作紫金山国家森林公园3个森林参数空间分布(图3)。

表6 3个森林参数不同模型估测精度比较Table 6 Precision comparison of three forest parameters using different estimation models

注：a：多元线性回归模型估测值与实测值散点图；b：支持向量机模型估测值与实测值散点图；c：随机森林模型估测值与实测值散点图；a1,b1,c1：平均树高；a2,b2,c2：平均胸径；a3,b3,c3：单位蓄积量。

3个森林参数空间分布趋势一致，平均胸径大的地方平均树高和单位蓄积量也大，这符合林分生长规律。较高的平均树高、平均胸径和单位蓄积量主要分布在紫金山国家森林公园北坡和南坡的中上部。北坡交通不便、距离居民点较远、旅游景点较少，人为干扰活动较少；南坡中上部海拔较高，立地条件较好，因此这2处森林生长良好，林分平均胸径、平均树高、单位面积蓄积量较高。较低的平均树高、平均胸径和单位蓄积量主要分布在紫金山国家森林公园的山脊、南部边缘和东北角。山脊由于海拔较高，土层瘠薄、立地条件差，林木生长不良；南部边缘地带和东北角建筑物、道路、草地、水体等非林地面积比例大，人为干扰活动多。

注：a：平均树高分布；b：平均胸径分布；c：单位蓄积量分布。

3 结论与讨论

为了探究稀疏型机载激光雷达数据在森林平均树高、平均胸径及单位蓄积量估测方面的应用，本研究以江苏省南京市紫金山国家森林公园为研究对象，通过对稀疏型激光雷达点云数据进行处理，提取了高度和密度两种共48个特征变量，并结合紫金山国家森林公园森林资源规划设计调查实测数据，建立小班平均树高、平均胸径和单位蓄积量反演模型用于估测3个森林参数。结果表明，非参数回归方法优于多元线性回归，随机森林反演精度最高，随机森林估测3个森林参数相关系数R2都在0.86以上，其中林分平均树高估测精度最高(R2=0.89)，平均胸径和单位蓄积量估测精度R2分别为0.87，0.86；高度百分位数和累计高度百分位数是影响小班平均胸径、平均树高和平均蓄积量的主要特征参数。利用随机森林模型制作的紫金山国家森林公园3个森林参数分布图与林分实际状况符合。因此，在本研究中，应用稀疏型机载激光雷达数据结合森林资源调查数据估测林分平均树高、平均胸径以及单位蓄积量具有较高的精度。

相对于其他遥感技术，激光雷达是遥感技术领域的一场革命。由于激光雷达数据可直接获取地物的三维空间地理坐标信息，对地形和植被空间结构的探测能力和抗干扰能力强，因此在估测森林参数方面具有独特的优势。Lefsky等[25]将Landsat TM数据、高空间分辨率的机载LiDAR数据、机载高光谱的AVIRIS数据和机载激光雷达数据进行比较，结果表明机载激光雷达数据反演植被结构参数效果最好。高密度机载LiDAR点云数据(>3个/m2)可以精确地对单木尺度森林参数进行估测，但数据获取成本较高，数据处理相对复杂，很难应用于大尺度森林参数估测[26]。为了降低获取数据成本及数据处理的复杂度，并扩大机载LiDAR技术的应用范围，通常会获取一定精度的低密度机载LiDAR点云数据，但是由于单木上点云数量的限制，无法精细刻画完整的树冠表面结构，因此稀疏点云数据多利用其高度和密度特征变量进行较大尺度的森林参数估测[14]。本研究采取以小班尺度提取稀疏型机载LiDAR数据的高度和密度等特征因子建模的方法，取得了可靠的估测结果，为稀疏型机载LiDAR数据在森林资源规划设计调查中的推广应用提供科学参考。

目前的研究中，虽然机载激光雷达在森林参数提取和估测中有着诸多优势，但不可避免的仍有较多的不足与局限。直至今日，使用机载激光雷达获取数据的成本依然很高，难以做到覆盖大片区域，对于很多研究项目来说不能作为主要数据获取手段。而且，机载LiDAR技术提取森林参数反演的精度与点云密度和森林结构有着密不可分的关系，大部分情况下较低的点云密度很难做到对单木的信息提取。另外，由于LiDAR技术起步较晚，处于发展阶段，历史数据缺乏，导致很难对森林不同时期的监测数据做出对比与分析。