唐佳佳 董 婧 杨永均 许木桑 雷少刚 华 夏

(1.中国矿业大学环境与测绘学院,江苏 徐州 221116;2.矿山生态修复教育部工程研究中心,江苏 徐州 221116;3.山东省煤田地质局采煤塌陷地与采空区治理工程研究中心,山东 济宁 272100)

生态修复是实现矿山固碳增汇和碳中和目标的重要途径[1-4]。生态修复后,植被覆盖度得到提高,水土保持、水源涵养等多项生态系统服务得到改善[5-6],同时,植被的恢复可以显著增加生态系统碳固存[7]。植被碳储量不仅是研究生态系统与大气间碳循环的基本参数,也是反映群落生态系统结构和功能特征的重要指标,准确评估碳储量是估算生态修复固碳效应的关键,对生态修复物种选取、生态修复后生态系统结构优化具有重要意义[8]。

目前,矿区植被碳储量估算得到了不少学者的关注,现有研究主要是通过实地调查直接获得碳储量数据,或利用基于土地覆被变化的因子法、基于遥感的CASA模型间接获得。在实地调查法方面,朴世龙等[9]根据森林清查数据,对植被碳库及其变化做了大尺度的研究,提高了生态系统的植物碳储量的估算精度。原野[10]运用土壤和植物采样分析方法,探讨了平朔露天矿复垦生态系统的碳存储机理与效应。然而,传统实地调查方法对于森林生态系统碳储量的估算结果在尺度转化上存在较大误差,存在空间局限性[11]。利用遥感影像数据,则能够实现大范围植被碳储量的快速评估。李帆等[12]基于Landsat影像,对伊敏矿区土地利用覆盖信息进行分类提取表述土地利用变化,进而监测了碳储量变化;HOU等[13]、吴国伟等[14]通过RS和GIS技术,采用改进的CASA模型对煤矿区植被碳储量进行了测算研究;张璐[15]基于遥感影像并结合InVEST模型对锡林浩特市大型煤电基地土地利用、覆盖时空动态变化进行了分析,定量评估了生态系统碳储量;刘英等[16]利用InVEST模型评估了新疆红沙泉矿区的碳储量变化,探究了采矿扰动下土地利用变化对固碳功能的影响。然而,现有的矿区植被地上碳储量遥感监测研究大多采用MODIS、Landsat和Sentinel等中低分辨率数据,监测结果分辨率较低[17-18]。矿区存在地形复杂、生态修复斑块小、植被配置空间异质性强等特征[19],中低分辨率遥感数据难以实现矿山生态修复区植物群落地上碳储量的精细反演与动态监测。

无人机机载激光雷达和高光谱遥感具备分辨率高、数据获取受时间与空间限制小等优势,逐渐得到了广泛应用,为矿区植被碳储量监测和评估提供了新途径[20]。为此,本研究针对矿山生态修复区,基于无人机遥感手段,提取植物群落遥感结构参数与光谱特征,建立植被地上碳储量反演模型,评估无人机机载激光雷达和高光谱反演矿区植被地上碳储量的能力。

1 研究区与数据

1.1 研究区概况

研究区位于大柳塔采煤沉陷区水土保持示范园,地处陕北黄土高原北侧和毛乌素沙漠东南缘(图1)。该区南北长2 km,东西宽约300 m,面积约0.6 km2。区内具有鲜明的严寒干燥、风沙频繁的季节性气候特征。地势北高南低,中间高而东西低。研究区于2005年开始进行生态修复,主要包括水土保持与植被重建措施,采用乔、灌、草结合的治理模式,经过生态修复,形成了以杨树、樟子松、沙棘、油蒿为主的植物群落。

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

1.2 数据采集与预处理

1.2.1 无人机遥感数据

2021年7月28日,利用大疆M600无人机采集了激光雷达和高光谱数据。其中,激光雷达数据采集共规划航线3条,飞行高度为90 m,飞行速度为5 m/s,水平视场角为360°,垂直视场角大于20°,平均点云密度为130个/m2。激光雷达数据处理采用Li-DAR 360和python软件。高光谱数据获取过程中,为保证区域全覆盖,设置高光谱传感器8条航线,飞行高度为140 m,空间分辨率为0.19 m,旁向重叠率为57%,光谱范围为400～1 000 nm,共112个波段。数据处理过程中,辐射定标、大气校正、几何校正以及影像波段合成等流程均借助ENVI软件实现。

1.2.2 实地调查数据

实地调查数据与无人机遥感数据同步进行采集,共设置样地60个,然后取2 m×2 m草本样方、5 m×5 m灌木样方、10 m×10 m 乔木样方,控制每个样方间距大于50 m,共设置了382个样方。利用GPS标记每个样方中心点坐标和高程,调查以乔、灌为优势种的植物群落的植物类型、高度、胸径、冠幅,钻取乔、灌植株树干、树枝、树皮、树叶4个部位的样品。采集草本群落的地上部分,去除植物脱落物和其他杂物后装袋标记。在实验室内将采集的乔、灌、草样品先置于105 ℃下进行杀青30 min,再在温度65 ℃下烘干48 h至恒重,用电子天平迅速称量,得到乔、灌、草群落碳储量测定样品与草本群落地上生物量样地实测数据。

将烘干后的样品通过粉碎机磨碎后再进行研钵研磨,经100目网筛筛选后用微量电子天平称取10 mg样品装入塑封袋后测定样品碳含量。参照《森林土壤有机质的测定及碳氮比的计算》(GB 7857-87)测定植被样品碳含量,取实测碳含量均值作为进一步估算碳储量的参数。利用乔、灌木树种异速生长方程估算生物量,再乘以含碳系数计算单株乔、灌木地上碳储量,结果见表1。利用草本群落地上生物量实测数据乘以含碳系数来计算草本群落地上碳储量。

表1 研究区主要植物生物量估算的异速生长方程Table 1 Allometric equation for main vegetation species in the study area

2 研究方法

2.1 特征变量提取

将激光雷达点云数据分为地面点、建筑物、低矮植被、中间植被和高层植被5类,根据点云数据中的首次回波点采用克里金插值法生成数字表面模型,基于地面点建立不规则三角网生成数字高程模型,最后将数字表面模型与数字高程模型作差得到冠层高度模型。同时,将激光雷达中的高层植被点、中间植被点和低矮植被点划分为乔、灌、草3种植物群落类型。

从激光雷达植被点云中提取了与植物群落高度相关的变量40个;从高光谱影像中提取了16个植被指数以及均值等26个纹理特征,作为植被地上碳储量估测模型的参数,分别构建乔、灌、草植物群落地上碳储量反演模型(表2)。在纹理特征提取过程中,为避免单一的像元值产生随机性误差,对图像进行3×3窗口的均值滤波处理,由最小噪声分离变换后高光谱数据的前两个波段参与后续纹理计算,分别命名为B1和B2。基于一阶及二阶概率统计法提取植被纹理特征,纹理特征因子命名规则定义为波段—纹理指数—概率统计阶数[27]。

表2 LiDAR与高光谱特征变量Table 2 Features of LiDAR and hyperspectral

2.2 碳储量反演方法

2.2.1 碳储量估算模型

本研究构建了LiDAR特征、高光谱特征、LiDAR+高光谱特征3个特征变量组合。基于变量投影重要性进行特征变量筛选,从而筛选出与植物群落地上碳储量相关性高、敏感性强的参数,选择重要性大于0.8的变量构建逐步多元回归模型,进而计算不同特征变量组合反演植被地上碳储量的精度。最后构建研究区乔、灌、草3种群落的碳储量反演模型,实现研究区植被地上碳储量制图。

植被地上碳储量CS估算模型采用Schumacher-Hall公式[28]构建:

式中,LIi为第i个激光雷达结构特征变量;SIj为第j个高光谱特征变量;c、ai和bj为回归系数。

为了便于建模,将式(1)模型进行自然对数变换后得到:

变量投影重要性(VIP)用来判断单个自变量在解释因变量的重要性。若自变量对因变量的重要性越大,则VIP值越大。有学者建议用0.8作为临界值来区分重要变量与不重要变量[29],其计算公式为

式中,i为第i个自变量;k为自变量个数;ci为第i个相关自变量提取的主成分;r(y,ci)为因变量和主成分的相关系数,表示主成分对y的解释能力;wij为自变量在主成分上的权重。

2.2.2 精度评价方法

反演精度评价采用决定系数(R2)和均方根误差(RMSE)来检验,R2越接近1,反映模型精度越高;RMSE值越小,说明回归模型更准确。相关公式为

式中,Yi为实测植被地上碳储量;yi为预测植被地上碳储量;为估测植被地上碳储量均值;n为样本数量。

3 结果与分析

3.1 特征变量优选结果

将LiDAR点云和高光谱影像所提取的参数作为自变量,输入SIMCA 14.1中进行参数VIP值排序,筛选VIP＞ 0.8的特征参数用于后续构建乔、灌、草植物群落地上碳储量反演模型。其中,使用单一LiDAR特征组合构建乔、灌、草群落地上碳储量模型时,分别筛选出16、32、39个参数;使用单一高光谱特征组合构建植被地上碳储量模型时,分别筛选出28、28、23个参数;使用LiDAR+高光谱特征组合构建乔、灌、草群落地上碳储量模型时,分别筛选出51、46、61个参数。数据筛选结果表明:LiDAR和高光谱特征中的高度变量和绿色指数与地上碳储量的敏感度较高,这两种参数分别能体现植物群落垂直结构特征和光谱特征,能在反演植物群落地上碳储量中发挥关键性作用,对反演精度贡献较大。

3.2 反演模型构建结果

3.2.1 反演模型

3种特征变量组合下的植被碳储量估算模型及精度见表3。在LiDAR特征组合中,构建乔木群落地上碳储量估算模型时,HP20th、d9两个参数入选;构建灌木群落地上碳储量估算模型时,Hskew、Hsms等5个参数入选;构建草本群落地上碳储量估算模型时,d1、Hcv等5个参数入选。在高光谱特征组合中,构建乔木群落地上碳储量估算模型时,GI、MSAVI等9个参数入选;GI、mSR750等7个参数入选参与构建灌木群落地上碳储量估算模型;GI、NIR_G两个参数入选构建草本群落地上碳储量估算模型。在LiDAR+高光谱特征组合中,HP95th、HP99th等5个参数入选构建乔木群落地上碳储量估算模型,HP1st、HP99th等7个参数入选构建灌木群落地上碳储量估算模型,d1、HP90th等5个参数用于构建草本群落地上碳储量估算模型。

表3 不同特征组合构建的植被地上碳储量反演模型Table 3 Inversion models of aboveground carbon storage of vegetation based on different feature combinations

由表3可知:LiDAR+高光谱特征组合在植被碳储量估算中表现最佳,R2最高,乔、灌、草植物群落地上碳储量精度分别为0.87、0.73、0.72。而使用单一LiDAR特征时,构建的乔灌草群落地上植被碳储量反演模型精度分别为0.64、0.44、0.36;使用单一高光谱特征组合时,碳储量估算模型精度R2分别为0.42、0.43、0.41。因此,尽管LiDAR数据或高光谱数据可以单独估算植被地上碳储量,但植物群落碳储量不仅与从LiDAR数据中提取的群落结构特征有关,还与高光谱信息中反映的光合固碳特征有关,因而结合LiDAR和高光谱两种数据能显著提高植物群落碳储量的反演精度。

3.2.2 模型精度

将在3种特征组合下使用偏最小二乘法筛选得出的特征变量与地面实测样地数据,分别作为自变量和因变量输入线性多元逐步回归模型,结果如图2所示。

图2 3种特征组合下植物群落地上碳储量实测值与估测值Fig.2 Measured and estimated aboveground carbon stock of plant communities under 3 feature combinations

由表3和图2可知:LiDAR特征组合中,乔木群落地上碳储量模型反演精度最高,为0.64,表明Li-DAR提供的高度百分位变量特征可以有效反映乔木群落的垂直结构信息,从而提高地上碳储量反演精度。在高光谱特征组合中,灌木群落地上碳储量模型反演精度最高,为0.43,说明高光谱影像能体现灌木群落丰富的光谱及纹理特征,对于提升地上碳储量反演精度的贡献较大。在LiDAR+高光谱特征组合中,融合激光雷达与高光谱数据提取高度百分位变量特征和绿色指数共同参与植被地上碳储量建模,乔、灌、草3种群落的植被地上碳储量模型反演精度分别提升了23%、29%、31%,反演精度的提升受益于LiDAR检测到的结构特征和高光谱捕获的植被光谱信息的融合互补。

以上反演结果表明:采用数据融合、特征消冗、参数优选技术可有效提高植物群落地上碳储量的反演精度。首先采用偏最小二乘算法对两种数据3种组合方式进行总体参数筛选,分别筛选出一定数量的回归参数;然后使用逐步多元回归算法对参数进行了优选,两次参数优选过程能去除与植物群落碳储量相关性低、敏感性差的参数[30]。综合反演结果可知:在两种数据的3种组合方式下,单一使用两种数据反演均有较大的局限性,采用LiDAR和高光谱两种数据相结合可发挥各自优势,能同时表征植物群落水平和垂直结构,体现植物群落的三维结构特征,进而提高植物群落碳储量的反演精度。

3.3 碳储量反演结果

本研究采用线性多元逐步回归算法,反演得到的研究区植物群落地上碳储量分布如图3所示。

图3 3种特征变量组合反演的植被碳储量Fig.3 Inversion results of vegetation carbon stock using three characteristic combination models

由图3可知:群落碳储量高值分布在研究区东北部的杨树群落与研究区中部沙棘群落,图3(a)中灌木群落碳储量较为均一,道路两侧有少量高值分布,乔、灌、草3种群落之间碳储量差异较大,研究区西部草本群落碳储量最低。图3(b)中,乔木与灌木群落碳储量较大,草本群落区域碳储量较低。研究区乔灌草群落碳储量分布表现出明显的差异性,碳储量高值分布在沙棘林及杨树林区域,各灌木群落反演结果均较好,但乔木群落有明显的低估现象。图3(c)中各植物群落碳储量的空间差异性明显,植被地上碳储量高值区主要分布在研究区北部、中部的乔木林附近,灌木群落地上碳储量较为均一,低值区主要分布在研究区主干道草本群落附近,可以有效地反映出不同群落间地上碳储量的差异。

综合分析可知:使用LiDAR与高光谱数据融合反演植物群落地上碳储量时,植被地上碳储量的空间分布趋势与研究区实际相符,有效地反映了复杂地形下植被地上碳储量的空间异质性。其中,研究区北部及中部区域多为人工种植的杨树或樟子松林,生长状况较好,树高与胸径均较大,因此地上碳储量较高;在研究区域的主干道附近种植有行道树,因而有少量碳储量高值点分布。

4 结 论

(1)评估了激光雷达和高光谱遥感预测植被地上碳储量的能力,并构建了一种联合LiDAR和高光谱数据的矿山生态修复区乔、灌、草3种植被类型地上碳储量的反演模型。

(2)采用LiDAR和高光谱两种数据相结合可发挥各自优势,进而提高植物群落碳储量反演精度,乔、灌、草群落地上碳储量R2分别达到0.87、0.73、0.72。LiDAR 和高光谱特征中的高度百分位变量和绿色指数分别与地上碳储量的敏感度最高,对提升地上碳储量反演精度贡献最大。碳储量的提高不仅与从LiDAR数据中提取的群落结构特征有关,还与高光谱信息中反映的群落光合固碳特征有关。

(3)总体来说,融合激光雷达与高光谱数据在估算矿区复杂地形和植被配置下的植被地上碳储量方面具有一定的优势,可为准确评估生态修复区植物生态系统碳储量提供技术支撑。今后可利用无人机载激光雷达与高光谱对矿区植被碳储量开展持续监测、时序变化分析、固碳效益评估等相关工作。