无人机多光谱和LiDAR的红树林精细识别与生物量估算
2022-07-06吴培强任广波张程飞王浩刘善伟马毅
吴培强,任广波,张程飞,2,王浩,刘善伟,马毅
1.自然资源部第一海洋研究所,青岛 266061;
2.山东科技大学测绘与空间信息学院,青岛 266590;
3.中国石油大学(华东)海洋与空间信息学院,青岛 266580
1 引言
红树林是生长在热带和亚热带海岸潮间带、受到海水周期性浸淹的木本植物群落,是全球4大湿地生态系统中最具特色的一个湿地生态系统(林鹏,1997),具有巨大的生态价值和社会价值,不仅在御风消浪、固岸护堤和维持生物多样性等方面所发挥重要作用,也是珍稀濒危水禽重要栖息地,鱼、虾、蟹、贝类生长繁殖场所。特别是2004年印度尼西亚海啸的发生,让世界看到红树林在防灾减灾中的重要作用。另外,红树林生态系统具有强大的固碳能力,其碳汇能力远高于沿海滩涂盐沼,是地球碳循环至关重要的组成部分,以面积不足陆地总面积的0.1%,却埋置了大于10%的陆源可溶性有机碳,是重要的蓝色碳汇,在减缓全球气候变暖过程中发挥着重要作用(Duarte 等,2005)。碳储量高和生态服务价值大,使得红树林生态系统成为“碳中和”和生态系统服务付费策略研究的重点对象之一。
随着毁林填海造田、围海养殖、城镇建设等粗放型的人类开发活动的影响,全球红树林面积在过去的50年减少了30%—50%,在全球气候变化背景下,按此减少速度,红树林及其特有的生态功能在未来100年可能会从地球上消失,其快速减少的趋势必然对海洋中陆源可溶性有机碳产生深远影响(Duke 等,2007)。中国红树林在20世纪50—90年代初的40 多年间锐减了68.7%(廖宝文和张乔民,2014)。近年来,随着红树林人工栽种、保护等生态修复的开展,扭转了中国红树林面积急剧减少趋势,但红树林面积尚未得到完全的恢复,生境退化、生物多样性降低、外来物种入侵严重等问题还依然存在。因此,开展红树林生态系统监测和生物量的精确评估,对于红树林保护管理、生态修复具有重要的意义。
由于红树林生长在潮间带淤泥滩,潮水周期性淹没,现场调查和采样费时费力。遥感技术具有探测范围大、时效性强、受地面条件限制少等优势,已广泛的应用于红树林监测。特别是近年来,随着高分辨率多光谱卫星和无人机数据的渐多(Quickbird、IKONOS、WorldView-2 和WorldView-3等),使得红树林类型精准识别和生理参数高精度估算成为可能,(Wang 等,2004;Jia 等,2014;Manna等,2014;朱远辉等,2014;Zhu等,2015;唐焕丽等,2015;Pham 和Brabyn,2017;闻馨等,2020;徐逸等,2021),但由于红树类型间具有相似的光谱反射特性,仍然会发生错分和漏分现象,从而影响其生物量估算精度。激光雷达可以穿透植被冠层,测定地表和林冠的高度信息以及树冠面积、树木密度等信息,这些参数在红树林生物量估算起着重要的作用(Simard 等,2008;Olagoke等,2016;Feliciano等,2014,2017;Hickey等,2018;Fatoyinbo 等,2018;Navarro 等,2020;Salum 等,2020,2021);然而激光雷达由于缺少光谱信息在红树林种类识别方面不理想(Lu 等,2016)。近年来,大多数研究发现主动激光雷达数据获取的冠层高度、被动光学数据的光谱和纹理信息相融合可以有效克服单数据源的不足(Chadwick,2011;Li 等,2019;Wang 等,2019;Qiu 等,2019;Hu 等,2020;Zhu 等,2020),从而提高红树林种间类型识别和生物量估算精度。但是少有开展基于红树林单木识别的生物量估算研究(Li 等,2012;Yin 和Wang,2019;Silván-Cárdenas等,2020),且没有开展进一步的红树种类识别,这会影响红树林生物量高精度估算。
本文选择广西茅尾海红树林保护区为研究区,综合利用无人机多光谱遥感数据和机载激光雷达点云数据,结合地面实测的不同类型红树生长参数,开展基于无人机载主被动遥感手段相结合的红树林单木识别方法研究,并估算了红树林地上生物量。本文研究的技术方法和结果可为中国红树林生态系统的保护和修复,以及红树林碳储量的估算提供技术支撑。
2 数据与方法
2.1 研究区
茅尾海 (21°46′N—21°54′N, 108°28′E—108°37′E)地处广西壮族自治区钦州湾北部,钦江、茅岭江的径流带来了大量陆地泥沙,造就了大面积的淤泥质滩涂。该区域为南亚热带季风气候,潮汐为不规则全日潮,平均潮差为2.51 m,温暖湿润的气候以及丰富的滩涂资源为红树林生长提供了良好的条件。茅尾海红树林保护区成立于2005年,是以红树林生态系统保护为主的兼顾越冬鸟类栖息地保护的自然保护区,保护区内共有红树植物11 科16 种,其优势物种包括桐花树、秋茄、无瓣海桑和老鼠簕,其中还有白骨壤、卤蕨等非该区域优势种的红树类型。自2000年左右,保护区内较早开展了无瓣海桑的引种,随后多次大规模进行了秋茄、桐花树等本土红树类型的人工栽种,截止2020年,保护区内红树林面积相比2000年增加了约120%(王浩等,2020)。但由于广西最大的沙井大蚝(近江牡蛎)养殖区紧靠保护区,极大干扰了红树林的保护和修复。
图1 研究区示意图Fig.1 Location of study area
2.2 数据
2.2.1 遥感数据
2020年11月29日,在低潮位期间同步获取了无人机多光谱遥感和LiDAR 点云数据。无人机飞行数据获取过程中,航向重叠度设定为75%,旁向重叠度设定为50%。无人机多光谱数据获取使用的是大疆精灵4 多光谱版,该设备集成了1 个可见光镜头和5 个多光谱镜头(蓝光,绿光,红光,红边和近红外),分别负责可见光成像及多光谱成像。对于获取得到的多光谱数据利用ContextCapture软件进行图像拼接,对于拼接后的影像,基于ENVI 5.3 软件对其进行辐射校正等处理。无人机所搭载的LiDAR 为蜂鸟Genius 激光雷达系统,其360°的扫描视场和320 kHz 的高扫描频率,可以极大提高地物测绘效率。无人机激光雷达数据经过航带匹配、噪声点去除等预处理,得到较为完整的点云数据。无人机遥感数据的详细参数见表1,无人机激光雷达点云数据和剖面数据见图2。
表1 无人机遥感数据参数Table 1 The parameters of UAV data
图2 无人机激光雷达数据和点云剖面Fig.2 UAV-LiDAR and illustration of mangrove point cloud profile
2.2.2 现场调查数据
分别于2020年10月和2020年11月对广西茅尾海红树林保护区开展了现场样本采集,共获取了28 个站点的经纬度位置、类型等解译信息和植株高度、胸径、基径、冠幅等红树生长参数,包含4 个优势物种(秋茄、桐花、无瓣海桑和老鼠簕)以及与红树林混生的盐沼植物茳芏。对于桐花树、秋茄和无瓣海桑等3种乔木,现场调查样方设置为5 m×5 m,在样方内逐株测量株高、基径(离地面3—10 cm 处,避开有隆起的部位)、胸径(大于1.3 m的植株,测量第一个分叉处下方)、冠幅和分支数等。对于灌木老鼠簕,样方设置为1 m×1 m,测量其平均高度和株数,选取具有代表性的两株带回实验室烘干称量,作为地上生物量测算的现场样本依据。对于盐沼草本植物茳芏,样方设置为0.5 m×0.5 m,测量其平均高度和株数。现场调查站位分布见图3所示,红树生长参数统计见表2。
表2 茅尾海红树林生长参数Table 2 Growth parameters of mangrove in Maowei Sea
图3 现场调查站位图Fig.3 Site survey station map
2.3 方法
2.3.1 红树种间分类方法
红树林生长在淤泥质潮间带,环境复杂,且受潮汐影响频繁,现场调查困难,因此,对于开展红树林种间类型遥感监测,样本少是面临的普遍问题。针对此问题,本文选择使用建立在统计学习理论基础上的支持向量机SVM(Support Vector Machine)分类方法,可以有效解决分类过程中小样本问题,同时能通过核方法解决线性不可分问题,是目前遥感图像分类中应用最为广泛的方法之一(唐焕丽等,2015;甄佳宁等,2019)。本文选用ENVI 5.3软件中的SVM 模块开展红树林分类,核函数选用径向基核函数,惩罚参数为100。
2.3.2 红树单木分层分区距离判别聚类分割方法
对于LiDAR 点云数据,基于LiDAR 360 软件中的面向点云的距离判断聚类算法来提取单木的空间信息和冠幅数据。经现场调查发现,由于研究区内红树生长时间和环境的不同,导致红树生长状况也不尽相同,从而表现为树木的高度和冠幅存在较大差异,假如整景点云数据作为一个整体利用统一的阈值开展红树单木树高和冠幅的反演,易造成过分割和弱分割现象,从而影响红树的生物量估算精度,因此本文基于树高和类型将研究区分为若干区域,对于不同的区域选择合适的分割阈值提取红树单木结构信息,这样可以有效的避免由于同一分割尺度带来的误差。
2.3.3 精度验证
遥感影像分类精度验证旨在确定分类过程的准确程度,本文使用混淆矩阵法对红树林分类精度进行评价。具体指标包含:总体精度OA(Overall Accuracy)、制图精度PA(Producer Accuracy)、用户精度UA(User Accuracy)和Kappa 系数(贾明明,2014)。具体公式如下:
式中,Khat为Kappa系数,xii为第i行i列的样本数,xi+为分类结果中第i类总和,x+i为地表真实数据中第i类的总和,n为类别的数量,N为样本总数。
2.3.4 红树林地上生物量估算模型
本文通过查阅前人构建的红树生长参数—生物量异速生长方程(多是通过砍伐采样方式研究得到),计算红树地上生物量,将其作为真实的样本数据。由于调研的文献中所研究的红树林类型、树龄、地理位置和生长环境等均有差别,其所构建的异速生长方程也不尽相同,为此本文根据研究区内红树林类型、树龄、生长环境因素等特点,选取与研究区同等纬度、邻近区域和相同生长形态的红树异速生长方程,计算得到站位中每棵树的地上生物量。所选用红树林异速生长方程如表3所示。
表3 不同红树类型异速生长方程及茳芏回归模型Table 3 Algorithms of allometric growth equations for different mangrove species and Cyperus malaccensis of regression model
3 结果与分析
3.1 红树林种间类型分类
现场调查发现,无人机飞行区域主要分布有7 种地物类型,包括4 种红树:秋茄、桐花树、无瓣海桑和老鼠簕,1种本土盐沼植被茳芏,还包括裸露的潮滩和水域。本文基于现场调查站点数据和无人机低空飞行获取的照片,选取了316个独立的单一红树类型区域作为分类模型中训练和验证样本,这些样本均匀分布于整个研究区,包括红树不同的生长状况,其中训练样本和验证样本的比例为2∶1。利用LiDAR 360 点云处理软件,计算得到研究区的高度特征,提取了冠层高度模型(canopy height model)、高程标准偏差(standard deviation of elevation)、高程方差(variance of eleva⁃tion)、冠层起伏率(canopy relief ratio)、高程平均绝对偏差(absolute average deviation of elevation)和高程L2(elevation L2)(Li 等,2019)等6 个高度特征参量。将高度特征与多光谱特征结合用于红树林类型识别,对分类结果展开分析。
本文以无人机多光谱和LiDAR 点云数据为数据源,利用支持向量机分类方法开展了红树林种间类型识别研究,从总体精度、用户精度和生产者精度等方面进行了比较分析(表4),分类结果图见图6。其中,基于无人机多光谱和LiDAR 提取的高度特征相结合的分类精度最高,总体精度可达90.69%,较单数据源都有所提高;基于无人机多光谱数据的分类结果精度要高于基于LiDAR 点云的结果。可能是研究区内红树林冠层密度较高,能够穿透冠层的点云数量较少,形成的回波强度较为一致,红树林种间类型高度和冠幅存在较大的重叠区,仅利用LiDAR 点云数据难以做到红树林类型的高精度识别。从表4 的分类结果可以看出,秋茄的分类效果最好,其制图精度均在90%以上。桐花树识别精度大都在90%以上,但其中利用LiDAR 的高度特征识别桐花树时,用户精度最低,仅为51.16%,被错分为秋茄和老鼠簕较为明显。老鼠簕为灌丛结构,与外围的稀疏、矮小且高潮时受潮汐浸没的桐花树在纹理和光谱上差异较小,仅利用多光谱影像易造成错分现象,将多光谱和LiDAR 数据结合时,增加了老鼠簕高度特征,从分类精度提高可以看出多特征结合能解决这一问题。对于无瓣海桑,现场调查时,发现研究区内仅有4株,训练样本的缺少导致其分类精度最低。对于盐沼植被茳芏,其生长密度最大且高度均一,且无人机拍摄时,茳芏正处于物候期,植株干枯,从而基于高度特征分类效果要好过多光谱影像的光谱特征。潮滩无植被覆盖,地形起伏小,且高度最低,由此仅利用高程特征便可以很好对其区分,分类精度达到93.83%。由于水体对可见光和激光雷达的吸收性和衰减效应较其他地物要高,因此水域的分类精度都较高。
表4 基于不同遥感数据特征的红树林分类精度Table 4 Classification accuracy statistical table base on different characteristics of remote sensing data
本文还计算了多光谱波段和激光雷达高度特征的贡献率(图4)。可以看出,分类过程中,高度特征的分类精度要比多光谱波段高。在所有的高度特征中,冠层高度模型有很好的区分红树种类的能力,仅靠一个特征就达到了72.85%的分类精度,相对于总体精度的84.3%,其他高度特征发挥的作用要小一些。除此之外,高程L2、高程方差、高程标准偏差也都表现出比较重要的作用,这些特征的分类精度大都比多光谱波段要略高。在多光谱波段中,红波段贡献率最高,为55.78%。虽然各多光谱波段分类精度较低,但对于特定目标如水域、滩涂、茳芏等地物各个波段都有其优势。这些优势特征的融合所发挥的作用会出现超常效果,从而多光谱融合分类精度要高于激光雷达高度特征融合分类精度。
图4 不同特征参量的分类精度Fig.4 Classification accuracy based on different characteristic parameter
3.2 红树单木结构分割
基于本文提出的红树单木分层分区距离判别聚类分割方法,计算得到红树单木冠幅和树高(图6)。结合现场调查站点数据,共实测秋茄和桐花树分别为84棵和93棵,单木分割结果分别为72棵和56棵,精度分别为86.71%和60.21%。通过对比分析红树树高和冠幅等结构特征发现(图7),激光雷达估算红树高度的能力要远超冠幅。其中秋茄实测树高在1.8—4.4 m,中值为3.35 m,均值为3.36 m,而反演值在1.58—4.02 m,中值为2.99 m,均值为3.07 m,误差为0.4 m以内,比实测值略低。桐花树树高实测值在1.73—3.6 m,中值和均值分别为2.2 m 和2.24 m;反演值在1.94—3.2 m,中值和均值分别为2.38 m和2.47 m,精度要高于秋茄。从上述结果可以看出,激光雷达在估算树木高度方面具有很高的精度。但对于红树冠幅,反演值要比实测值偏大,秋茄和桐花树分割误差分别为21%和43%;通过对比发现,秋茄冠层边缘的老鼠簕和茳芏被错分为秋茄冠幅的一部分;桐花树由于生长茂密,树木冠幅交错,冠幅高低起伏差异较小,导致单木分割时出现弱分割,冠幅识别误差变大。可见对于复杂的红树林,LiDAR 点云单木分割还存在一定的局限性。
图6 基于无人机遥感影像的单木分割Fig.6 Single tree segmentation results base on UAV-multispectral
为提高研究区红树林地上生物量估算准确性,本文对单木分割结果进行了后处理,挑选出异常值区域,通过目视经验判读的方法对分类结果进行修正。对于秋茄冠层边缘的老鼠簕和茳芏进行重新提取和赋值;对于桐花树冠幅异常区域,利用邻近现场调查采集的桐花树密度、高度和冠幅等数据,重新定义该区域的单木数量和冠幅大小。
3.3 红树林地上生物量估算模型
对于研究区内植被地上生物量,考虑其生长特征和现场数据,本文分为以下情况开展估算。(1)对于秋茄、桐花树和无瓣海桑这3种乔木或灌木植被,基于现场采集的红树生长参数,结合已有研究中构建的不同类型红树异速生长模型,计算站位区域内各类型红树的单木地上生物量作为实测值,基于点云数据得到的红树单木结构特征,构建红树生长参数与生物量回归模型。(2)对于老鼠簕,多为单株灌木植被,构建归一化植被指数与现场采集生物量之间的关系模型来估算其地上生物量。(3)由于调查时间段内,茳芏正处于物候期,植株干枯,无法构建光谱与生物量之间的关系模型,因此采用单位面积生物量平均值乘以分布面积计算其地上生物量。
本文基于提取的秋茄、桐花树和无瓣海桑等3种红树类型的冠幅和树高等生长参数,结合站点的单木生物量数据,构建了包含树高H(Height)、冠幅C(Canopy)的单变量和双变量回归模型。由于研究区内无瓣海桑数量较少,现场调查时仅采集到其中两株的生长参数,因此未构建回归模型。
图5 不同遥感数据特征分类结果图Fig.5 Classification result of different characteristics of remote sensing data
本文构建了秋茄和桐花树地上生物量与树高和冠幅两种生长参数之间的线性回归模型,并使用决定系数R²来评价模型的精度。从表5 可以看出,双变量回归模型中,桐花树和秋茄的回归模型决定系数R²分别达到了0.832 和0.678。两种红树地上生物量与冠幅呈较高的正相关性,单变量模型中决定系数R2分别为0.824 和0.547;但两种红树的生物量与树木高度相关性却很低,决定系数R2只有0.019 和0.055。现场调查发现,研究区内桐花树和秋茄的植株高度与树龄的相关性小,大部分桐花树和秋茄的树高分别在2.3 m 和3 m 左右(图7),但其树龄却存在较大的差异,红树的基径、胸径和冠幅却随树龄的增加而变大。这也是表3中秋茄和桐花树的异速生长模型并未引入树高作为参数的原因。对于老鼠簕,构建了归一化植被和现场生物量之间的模型,决定系数R2为0.77。
图7 桐花树和秋茄生长参数实测值和反演值箱线图和散点图Fig.7 Box plots and scatter plots of field-measured and predicted values of growth parameters between Kandelia candel and Aegiceras corniculatum
表5 秋茄、桐花树和老鼠簕地上生物量估算模型Table 5 Estimation model of aboveground biomass of Kandelia candel,Aegiceras corniculata and Acanthus ilicifolius L
3.4 红树林空间分布特征分析
研究区内红树林多为人工种植,经过20 多年的生长,其空间分布具有较为明显的规律。研究区以桐花树和秋茄混生种植为主;但随着潮汐和
种间竞争的影响,红树林类型在垂直于堤坝呈现条带式分布:秋茄—老鼠簕/茳芏—桐花树。通过对比研究区地面高程发现,秋茄生长区域的平均高程比桐花树要高出大约15 cm。研究区内桐花树面积最大,为8.9 hm2,分布最为广泛,由于其耐盐和抗潮汐能力强,多分布于远离堤坝的潮沟两侧区域,树高较低且生长茂密。秋茄面积次之,为4.69 hm2,主要分布于受潮水影响较弱的邻近堤坝区域和被桐花树。茳芏环绕区域。无瓣海桑分布数量较少,整个研究区域内仅有4株,分布于潮水连通处,可能为种子自然扩散而来。老鼠簕和茳芏生长于秋茄和桐花树之间,其中老鼠簕更靠秋茄一侧多零散分布于秋茄树木之间,茳芏紧邻桐花树成片生长。
3.5 红树林地上生物量空间特征分析
根据图8 的研究区红树林地上生物量分布图,其空间分布规律与地物类型分布相关,研究区内地上生物量最小的区域是桐花树,从无人机现场照片解译发现,该区域为稀疏桐花树幼苗区,地上生物量为0.34 kg/m2,最大值为无瓣海桑的8.6 kg/m2。整体上,大部分地上生物量分布区为1.24—3.6 kg/m2。从各地物类型地上生物量分布来看,呈现无瓣海桑>桐花>茳芏>秋茄>老鼠簕。无瓣海桑作为高大乔木,在茅尾海红树林保护区内最大高度达20 m,且生长速度快,其他植被无法对其形成生存威胁,因此树木高大,枝干茂密,单位面积生物量最高。桐花树丛生,生长茂盛,地上生物量分布区间主要集中在1.76—3.73 kg/m2,均值为2.85 kg/m2。秋茄生长较为缓慢,其地上生物量分布区间为1.24—3.09 kg/m2,均值为2.32 kg/m2。老鼠簕地上生物量位于0.56—2.64kg/m2,其生物量均值为1.68 kg/m2;茳芏生物量位于2.00—3.73 kg/m2。
图8 红树林地上生物量分布图Fig.8 Spatial distribution of mangrove aboveground biomass
3.6 红树林监测精度的影响因素分析
在本文研究区无人机多光谱遥感图像的拍摄和拼接过程中,由于红树林树木之间的相似性和枝叶随风摇摆造成的同名点确定困难,导致了最终拼接而成的覆盖整个调查区域的遥感图像存在局部拉伸、畸变等情况。这种畸变和成像缺点在一定程度上影响了多光谱的分类效果和精度。同时由于遥感图像的空间分辨率高,对红树阴影区域的清晰成像,更进一步地影响了分类效果。将LiDAR 数据获取的高度特征与多光谱特征结合后,在很大程度上有效缓解了上述多光谱遥感图像的成像缺陷带来的分类精度不高的问题,因此结合无人机多光谱和LiDAR 点云数据在红树林精细识别方面具有较好的适用性。
由于飞行路线和倾角的原因,其激光雷达传感器正下方的点云密度要高于瞬时视场角两侧的点云密度,由此提取的树木高度特征也会随之变化,从而会影响红树种类的识别和高度的估算(Yin 和Wang,2019)。同时,现场调查时发现,稀疏的秋茄树下生长有较为茂密的老鼠簕和秋茄幼苗等植被,造成部分点云无法到达地面,影响了LiDAR 点云提取红树林木高度的精度,进而会降低红树林类型的分类精度和地上生物量估算准确性。再者,林下植被作为红树林生态系统的重要组成部分,能够体现整个生态系统的生物多样性、稳定性和健康状况,其地上生物量和种类丰富度是人工林的自然更新和演替的一个重要指标,在下一步研究中,将引入地面激光测量雷达,开展林下植被的研究。
图9 红树林林下结构Fig.9 Understory structure of mangrove
4 结论
本文选择广西茅尾海红树林保护区为研究区,以无人机多光谱和激光雷达影像为数据源,利用支持向量机分类方法对研究区红树林种类进行分类;在此基础上,结合实地测量的红树生长参数,提出了红树林分层分区距离判别聚类算法,实现了红树单木结构特征的准确提取,估算了研究区红树林地上生物量。得出以下结论:
(1)利用无人机多光谱和激光雷达数据在红树林种间类型识别是可行的,并且融合两种数据源的分类结果要比单一数据源或单一特征分类结果精度更高,从分类精度来看,融合无人机多光谱影像和激光雷达高程特征的分类精度最高,总体精度可达90.69%。Kappa系数为0.8819。
(2)本文提出的分层分区距离判别聚类算法能够较好的提取红树单木结构特征,秋茄和桐花树单木分割精度分别为86.71%和60.21%;红树树高的提取精度要高于冠幅,秋茄和桐花树结构特征中树木高度中误差0.36 m和0.18 m,而秋茄和桐花树冠幅分割误差分别为21%和43%。基于红树单木高度、冠幅构建了不同红树种类的地上生物量回归模型,双变量融合模型精度要高于单变量模型。其中,桐花树生物量估算模型精度最高,决定系数R²为0.832。
(3)通过构建的红树林地上生物量估算模型计算了研究区红树林地上生物量,发现无瓣海桑地上生物量最高,老鼠簕最低。其中桐花树和秋茄的地上生物量主要分布区间分别位于1.76—3.73 kg/m2和1.24—3.09 kg/m2,老鼠簕地上生物量位于0.56—2.64 kg/m2,茳芏地上生物量均值2.53 kg/m2。