曲延斌，王振锡*，吕金城，马琪瑶，郝康迪，葛 瑶

(1.新疆农业大学 林学与园艺学院，新疆 乌鲁木齐 830052；2.新疆教育厅 干旱区林业生态与产业技术重点实验室，新疆 乌鲁木齐 830052)

森林作为陆地生态系统中最大的碳库，在应对气候变化和实现碳中和目标过程中发挥着不可替代的作用[1]。森林植物吸收大气中的CO2并将其固定在植被或土壤中，从而减少CO2在大气中的浓度，具有显著的固碳能力。森林蓄积量是指一定森林面积上存在着的林木树干部分的总材积[2]，既可以反映一个地区的森林资源规模及水平，也能够衡量森林的健康状况、丰富程度和固碳能力，是森林资源调查中最为重要的指标之一。在传统森林资源调查过程中，森林蓄积量测算主要依赖于抽样技术和大量人工调查，具有工作量大、周期长、实效性差的缺点，快速、准确获取大面积森林蓄积量一直都是森林资源经营管理领域亟待解决的技术难题。

近年来，遥感技术凭借大面积同步观测和周期短的优势，广泛用于森林资源调查、规划和监测等，极大地提高了工作效率。尤其是机载激光雷达技术，突破了传统光学遥感的局限性，在获取地形和林分结构因子等信息方面具有不可比拟的优越性，为准确、高效获取森林蓄积量提供了新的技术途径[3-5]。许多学者利用机载激光雷达数据开展了森林结构参数信息提取研究[6-8]，不少学者在森林蓄积量估测方面进行了有益的探索[9-13]。王臣立等[9]利用LANDSAT-TM植被指数组合，结合森林清查数据，构建了热带人工林蓄积量的多元回归与神经网络模型；刘琼阁等[11]利用偏最小二乘法开展了基于TM影像的密云县森林蓄积量遥感估测，精度达到了90.1%；廖顺宝等[13]基于全国森林资源清查数据，采用不同森林类型面积和蓄积量之间的线性回归关系反映了公里网格尺度上的森林蓄积量空间变化。森林蓄积量的大小不仅与树高、林分密度等林分结构因子直接相关，而且地形因素对其也有着非常重要的影响。杨凤萍等[14]研究了油松(Pinustabuliformis)和华山松(P.armandii)蓄积量的垂直梯度变化规律，认为海拔梯度对林分蓄积量具有一定的影响。刘贵峰等[15]验证了海拔差异对天山云杉林蓄积测量有着显著影响。刘俊等[16]认为不同立地质量对森林蓄积量遥感估测有较明显的影响，且分立地质量类型建模精度要明显优于统一建模精度。R.Nelsonetal[17]利用脉冲激光系统获取森林树冠高度，以预测地面森林生物量与木材体积，模型预测的平均总蓄积量在地面值的2.6%以内。在天山云杉方面以机载激光雷达数据为数据源开展的森林蓄积量反演大多基于“树高-胸径”或“冠幅-胸径”关系模型，采用材积表法构建林分蓄积量反演模型[17-21]。这些研究，一方面没有综合考虑林分高度、林分密度以及地形条件等多变量因素与林分蓄积量的紧密关系，另一方面所构建的森林蓄积量模型大部分属于传统的基于回归关系的经验模型，而包含某种生物生长机理的生物生长模型较少。邓立斌等[22]以林分地位指数、林分年龄和林分密度为自变量构建模型对杉木的蓄积量进行预估，所建立模型精度较高。Schumacher模型是基于可变密度的收获预估模型，模型中充分考虑了对林分蓄积量具有决定作用的影响因子，除树高外，也包含年龄、林分密度、立地质量等因子，使林分蓄积量的反演精度得到较大的提高[23-24]。冉啟香等[25]，以 Richards和Schumacher模型为基础模型，采用不同密度指标，拟合油松林分的断面积与蓄积量，预估效果较好，且Schumacher模型的预估精度要更高于Richards模型，充分考虑到影响森林蓄积量的林分因子，大大提高了森林蓄积量以及收获量预估的准确度。

天山云杉(Piceaschrenkianavar.tianschanica)作为新疆天山山地最主要的地带性森林植被，其林分占新疆天然林乔木林地面积的44.9%，是构成天山乃至新疆森林生态系统的物质主体，在涵养水源、调节气候、应对气候变化、实现碳中和目标和维持山区生态系统平衡等方面发挥着重要的主导作用，为了更好地开展森林生态效益评价、建立生态效益补偿机制以及监测森林碳汇等工作，大面积天然林森林蓄积量的测定无疑成为区域生态综合监测的重要内容之一。由此，探索和研究利用现代信息技术进行天山云杉林蓄积量快速、准确测定的技术方法显得尤为迫切。

本研究在天山云杉林实地调查的基础上，结合机载激光雷达数据，综合考虑林分年龄、林分密度和立地质量等因子，采用Schumacher收获预估模型建立多变量天山云杉林蓄积量反演模型，旨在获得便捷、快速提取森林蓄积信息的技术方法，为探索山地天然林精准监测与评价提供技术途径。

1 材料与方法

1.1 研究区概况

研究区位于天山山脉中段的新疆农业大学实习林场(86°46′-86°57′E，43°16′-43°26′N)，该地域地处头屯河上游，乌鲁木齐市西南方位110 km(图1)，东西长达15 km，南北长达17.5 km。林场属于温带大陆性气候，平均海拔2 200 m，年平均降水量达600 mm，平均温度3 ℃，年日照大于1 300 h，无霜期约为140 d。森林主要分布在北坡、东北坡和西北坡，以天山云杉纯林为主，混杂有少部分天山桦(BetuIatianschanica)、山杨(Populusdavidiana)、天山花楸(Sorbustianschanica)等树种。

1.2 数据来源

1.2.1 影像数据来源 激光雷达影像数据拍摄于2020年10月15日，采用华测导航AS-1300HL多平台激光雷达系统(表1)，该雷达可以有效穿透植被，获取植被的垂直结构信息，更好地进行点云分类与森林参数提取。选择天气晴朗的中午时间起飞，在最大程度上消除树木阴影的影响。

表1 激光雷达主要参数

1.2.2 外业数据获取 2020年6月，按照森林资源二类调查技术规程，以无人机航拍图像为底图，结合实际地形、林分结构特征等布设样地。采用系统抽样法，在林场范围内以东西30 m和南北50 m为间隔确定样方中心，样方规格为30 m×30 m，通过遥感判读，保留100个天山云杉纯林样方。样方围测采用罗盘法，相对闭合差小于1/200，在样方边界涂刷红油漆进行标识。对样方内胸径大于5 cm，树高大于3 m的天山云杉全部每木检尺，并采用天山云杉二元材积表法计算样地蓄积量。采用华测i90惯导RTK，基于载波相位动态实时差分(Real-time kinematic)方法实时获得厘米级的实测样地中心点和4个角桩的坐标，同时记录样地的坡度、坡向、海拔等基本信息(表2)。

表2 样方每木检尺数据

2 特征参数提取与结果

利用GIS软件提取与天山云杉林蓄积量相关的林木树高、林木每公顷株数、样地坡度、坡向、海拔这些GIS因子建立蓄积量模型，提高对蓄积量的估测精度。

2.1 特征参数的提取

2.1.1 树高与每公顷株数的提取 对点云数据进行预处理，通过设定的高度阈值进行去航带重叠点与噪声点处理，基于Kraus滤波算法[24]结合中值滤波分离出地面点，再收集地面点进行地面高程的重建，插值生成数字高程模型(digital elevation model,DEM)。收集去噪后的第1次回波点插值生成的数字表面模型(digital surface model,DSM)，利用DSM减去DEM获得研究区林分的冠层高度模型(CHM，canopy height model)。冠层高度模型构建原理见图1。

图1 冠层高度模型构建原理

树高的提取是在冠层高度模型(CHM)的基础上对其进行局部最大值提取，利用ArcGIS软件空间分析中的焦点统计功能，对邻域内象元的最大值赋予焦点，利用栅格计算工具将焦点统计所得结果与CHM相减，差值为0的点为树木的最高点，因天山云杉为针叶树种，其树木最高点即为树冠顶点，运行重分类命令，去掉非冠层顶点，并将冠层顶点由栅格数据转换为矢量数据。

对获得的矢量树冠顶点与DEM、DSM和CHM的高程数据进行对比，结合树冠顶点密度情况，对伪树冠顶点进行剔除。并且根据研究区域实际树高情况，对低于2 m的树冠顶点也进行剔除，剩余的树冠顶点位置则为林木顶点，在ArcGIS软件中将剩余树冠顶点与CHM图层叠加，获取林木树高，并根据树冠顶点数，得到样地每公顷株数。提取结果见图2。

图2 树冠顶点高度提取

2.1.2 海拔、坡度、坡向的提取 使用ArcGIS10.4软件，将DEM数据导入，利用3D Analyst工具栅格表面执行等值线命令，在DEM的基础上生成等高线图(图3A)，执行坡向命令，生成坡向图(图3B)，执行坡度命令，以5°为间隔进行处理生成坡度图(图3C)，利用数据管理工具中的栅格属性获取等高线、坡度、坡向的平均值，为接下来的蓄积量的计算提供参数。

注:A等高线图像；B坡向图像；C坡度图像。

2.2 Schumacher收获预估模型

采用Schumacher(1939)收获模型构建天山云杉林蓄积量模型。Schumacher模型的蓄积量收获曲线为

M=α0e-α1/t

(1)

(2)

转换的一般形式为：

lnM=β0+β1t-1+β2f(SI)+β3f(SD)

(3)

式中:M表示单位面积上的林分收获量；t表示林分年龄；f(SI)表示地位指数SI的函数；f(SD)表示林分密度SD的函数；α0、α1、β0、β1、β2、β3为方程参数。

由于生态、气候等随机因素的影响，树高生长是一个随机过程，这个过程可以用林分树高生长的全体样本函数空间表示，并且受到立地质量的影响。但林分立地质量与林分密度不能当作随机因子来分析处理，随着林分年龄的增大，它对林分树高的影响愈发明显。对于天山云杉树种，其材积生产潜力与树高生长呈正相关，且受林分密度和间伐的影响不大，因此，本研究用天山云杉林分平均高来估计立地质量。

3 结果与分析

3.1 林分因子提取精度分析

将100个样地的树高、每公顷株数、坡度、坡向、海拔进行提取，为了验证基于高光谱遥感数据的样地因子提取精度，对其进行与实测值的对比检验。

基于激光雷达影像对树高和每公顷株数进行克里金插值法，对坡度、坡向、海拔进行裁剪、属性提取、样地因子属性平均数分析后，所获得的树高提取精度为89.64%，每公顷株数提取精度为85.13%，坡度提取精度为84.26%，坡向提取精度为82.63%，海拔提取精度为97.25%(表3)。其中提取每公顷株数时，由于天山云杉的冠幅形状与不同林木冠幅大小的差距，进行焦点统计时选择圆形，并将半径设置为1，以避免冠幅较小的云杉有遗漏，但会造成较大冠幅的云杉形成2个以上的伪顶点，导致每公顷株数提取数量高于实测数量，精度降低，除此之外，其他提取因子与实际测量所得数据的精度均在80%以上，因选取样地具有代表性且基数较大，有较高的适用性，所以使用遥感影像来对样地因子进行估测提取是适用于天山云杉的(图4)。

图4 立地因子提取精度对比折线

表3 立地因子提取值精度统计

3.2 基于激光雷达数据的Schumacher蓄积收获模型构建

将70个样地提取的林分因子数据代入Schumacher蓄积收获模型中，进行基于激光雷达影像数据的模型构建获得以下函数

lnM=5.027-99.217/t+0.047SI+0.001×A

(5)

SI=2.951-(-0.000 127)×B-0.004×C-0.000 174×D+0.179×t

(6)

t=5.6435×H-14.707

(7)

式中:M为样地蓄积量；t为林分平均年龄；A为林分每公顷株数；B为样地坡向；C为样地坡度；D为样地海拔；H为林分平均高。

样地提取因子代入模型后，所获蓄积量lnM的模型拟合度R2=0.8；立地质量SI函数拟合度R2=0.99；林分年龄t的函数R2=0.83。

3.3 模型检验

将剩余30个样地的林分提取因子数据代入模型中，求出样地的估测蓄积量，并与实地调查所得蓄积量进行比较分析，检验模型的适用性以及对天山云杉蓄积量反演的精度。得到样地估测蓄积量与实地测量蓄积量相比，其中精度最高为99.04%，最低为78.52%，平均精度为90.22%。

将样地实测获得的蓄积量作为自变量x，根据模型所估测获得的蓄积量作为因变量y，建立线性回归方程：y=0.964 1x+0.749 1，其中相关系数R=0.875(图5)。

图5 样地蓄积量实测值与预测值线性回归

将实测蓄积量与代入模型检测获得的估测蓄积量放入SPSS软件中，根据实测值与估测值建立残差分布图(图6)，并进行T检验，从残差图可以看出数据的残差值大多数在(-2，+2)残差区间内，说明回归模型能较好地符合原始数据，且进行配对T检验后，表4中T值为1.741，显著性概率P=0.092>0.05。因此，对天山云杉蓄积量的估测值与实测值差异并不显著，建立的蓄积量收获预估回归模型符合预期效果。

图6 云杉蓄积量残差分布

表4 样地蓄积量估测值与实测值配对样本检验

4 结论与讨论

4.1 结论

1)采用激光雷达进行遥感数据采集，不仅提高工作效率，在进行单木分割时，使用基于分水岭算法及点云聚类等多种分割算法，实现对测区整体激光点云的单木分割，获得清晰的分割效果，在进行因子提取时有较高的精度，其中样地海拔提取平均精度达到97.25%，树高提取精度为89.64%。

2)在对蓄积量进行建模过程中，对立地质量进行参数化，使用坡度、坡向、海拔等已知因子构建立地质量函数，其拟合度R2达到0.99，说明对立地质量进行在参数化的方法是可行的。

3)基于Schumacher蓄积收获模型，以林分蓄积量的倒数lnM作为因变量，将林龄的倒数作为预测变量，构建出的模型为lnM=5.027-99.217/t+0.047SI+0.001SD，林分的蓄积量随着年龄t的增大而不断增大，与立地质量和林分密度也有着正相关关系。

4)建立的收获预估模型与以往Schumacher蓄积收获模型的研究有着共性，且包含的分模型精度也都较高，实际值与模型计算出估测值差异性不显著，总模型对天山云杉林蓄积量估测精度达到90%，为以后的研究提供了一定的参考价值，利用Schumacher蓄积收获模型能够对天山云杉的蓄积量估测提供科学的理论基础，对天山云杉的营林措施制定也有着较好的指导意义。

4.2 讨论

相比较以往的遥感技术，激光雷达数据具有更加便捷、准确的特征，本次试验工作时间选在天气良好，天空无云，12:00左右以消除云朵、山体或云杉本身阴影对影像提取时可能产生的影响因子。

在对因子进行提取时，误差的因素包括在中龄林或幼龄林林分中，由于树木顶端生长的枝条过于幼嫩，顶端较细，在激光雷达进行数据采集时，激光探测信号不能在达到林木顶点时成功反射回在信号采集器中，处理点云数据时，选取到林木尖端以下较粗的枝干作为树木的顶点，导致提取树高低于实测树高，在后期蓄积量计算过程中，由于小树的蓄积较小，且本试验中，对林木树高计算方面，采用群体水平绘制树高曲线图的方法获取样地整体的平均树高，因此单株小树树高误差对试验结果蓄积量的计算影响较小。坡向会由于一个样地中存在着山脊或者山谷的情况，导致坡向有较大的起伏，提取坡向时会产生误差值。

对蓄积量反演过程中，需要考虑到树木的生长依赖于林地的立地条件，与立地质量和林分密度等因素是密切相关的，因此，在构建模型过程中，要获取立地质量和林分密度等因子信息数据，建立合适的模型。在加入坡向、坡度、海拔等地形因素构建模型估测蓄积量方面，冯凯等[26]在基于landset8遥感影像以立地因子为自变量，采用BP神经网络与逐步线性回归为建模方法，构建蓄积量反演模型，精度达到78.5%，苏迪等[27]在基于无人机影像的基础上，对地形因子进行提取，使用偏最小二乘回归法构建蓄积量估测模型，估测精度达到88.43%。江希钿等[28]以林分年龄、地位指数SI、林分密度指数SDI为输入变量，从模型间的相容性出发，建立了杉木人工林可变密度的全林分模型。

本试验基于激光雷达数据，获得的模型效果与蒋伊尹等[29]对大兴安岭兴安落叶松的天然林分收获量预估结果规律一致，林分的蓄积量与该林分的立地质量和林分密度有着密切的关系。

将立地指数与林分密度因素引入到Schumacher方程中，建立出天山云杉的蓄积收获预估模型，回避了过往的以胸径树高反演蓄积量的模式，由于树木胸径在遥感影像上无法直接获取，只能通过树高或者冠幅等其他因子进行建模转换，二次转换导致误差环节增多，在激光雷达点云数据下，直接获取了林分的树高与株数，结合密度与立地质量对林分的蓄积量进行反演，减免了通过中间变量换算林木胸径，从而减小了误差存在的几率，对于林分的蓄积量估测精度远高于以往含有胸径因素参与其中的蓄积模型，更有益于构建具有普适性的林分调查蓄积反演模型。