GEDI 不同算法组数据反演森林最大冠层高度和生物量精度比较

韩明辉，邢艳秋，李国元，黄佳鹏，蔡龙涛

（1.东北林业大学 森林作业与环境研究中心，黑龙江 哈尔滨 150040；2.自然资源部国土卫星遥感应用中心，北京 100048；3.辽宁工程技术大学 测绘与地理科学学院，辽宁 阜新 123000；4.桂林电子科技大学 电气工程与自动化学院，广西 桂林 541000）

森林生态系统是陆地生态系统的主体，是陆地碳循环最重要的碳库，在维护全球碳平衡和缓解全球气候变暖中发挥着巨大作用[1-4]。传统的森林参数测量方法需耗费大量的人力、财力，且效率较低。利用遥感技术对森林垂直分布和空间分布进行定量估计，对于提高森林生物量（Aboveground biomass，AGB）的估算精度，更好地了解森林生态系统功能、森林碳储量和全球碳循环具有重要意义[5-6]。激光雷达（Light detection and ranging，LiDAR）角分辨能力、距离分辨能力极高，抗干扰能力较强，具有很高的穿透性[7]，因而在探测森林垂直结构方面优于其他遥感技术，能够较好地获取森林内部结构信息和参数，获取精度较高的数据[8-9]。

星载激光雷达搭载在卫星上，轨道高度高、观测范围广，已被广泛运用于反演森林冠层高度等结构参数，目前常用的星载激光雷达包括线性体制的地球科学激光高度计系统（Geoscience laser altimeter system，GLAS）、全球生态系统动力学调查系统（Global ecosystem dynamics investigation，GEDI）和高分七号（GF-7），以及离散光子计数的先进地形激光高度计系统（Advanced topographic laser altimeter system，ATLAS）系统。于颖等[10]详细介绍了GLAS 波形数据的处理方法，建立了树高和生物量反演模型，最终结果表明利用激光雷达数据反演方法较传统方法精度明显提高。Adam等[11]在德国中部两个研究区评估了冠层覆盖、地形坡度、波束类型对GEDI 高度指标的影响，结果显示除陡坡地区外，在大多数条件下，GEDI 地面高程和冠层高度估算精度都较高。Iván等[12]在西班牙巴达霍斯省评估了5 个地中海生态系统中GEDI 冠层高度的精度，并建立了生物量反演模型，结果显示在稀疏的栎属为主的森林中相关性更强，生物量模型RMSE 范围为14.13～32.16 t/hm2，rRMSE 范围为38.17%～84.74%，这表明GEDI 数据在高度多层森林具有更高的不确定性。目前，虽然GEDI 数据已有一些研究，并取得了较好的结果，但对于GEDI 不同算法组数据的研究较少，而且GEDI L4A 数据自2021 年8 月份发布以来，并没有研究分析其生物量数据的准确性。

因此，为探究GEDI 不同算法组数据反演森林结构参数的能力，本研究拟将GEDI 的L2A 高程产品和L4A 生物量产品作为研究对象，分别提取L2A 中6 组数据和L4A 中7 组数据，利用机载验证数据，开展GEDI 不同算法组数据反演森林最大冠层高度和森林生物量的精度评价研究。

1 研究材料与方法

1.1 研究区概况

研究区为美国缅因州（Maine）的佩诺布斯科特实验森林（Penobscot experimental forest，PEF），其经纬度范围为（44.82°～44.87°N，68.58°～68.64°W），如图1 所示，左下角的红点是PEF 在缅因州的位置，黑色圆点是GEDI 光斑点在研究区的位置。研究区地面高程为9～91 m，地面坡度为0°～15°，总体比较平缓。研究区是大陆性湿润气候，夏季凉爽，冬季长而多雪，植被主要是针叶林，有云杉Picea asperata、冷杉Abies fabri、铁杉Tsuga chinensis以及一些松树Pinus和雪松Cedrus deodara[13]，冠层高度的范围约为3～30 m。

图1 Penobscot experimental forest 研究区及GEDI 光斑位置示意图Fig.1 Schematic diagram of the Penobscot experimental forest and the location of the GEDI laser points

1.2 研究数据

1.2.1 GEDI 数据

GEDI 于2018 年12 月5 日成功发射，质量为519 kg，表1 中记录了GEDI 的参数名称及参数大小。GEDI 由3 个激光器组成，其中2 个是全功率激光器，通过改变激光束指向，各自产生2 个光束。另外1 个激光器是覆盖激光器，被分成2 个光束，然后共产生4 个光束[14-15]，最终3 个激光器产生8 个地面光束轨迹（分别为BEAM0000、BEAM0001、BEAM0010、BEAM0011、BEAM0101、BEAM0110、BEAM1000、BEAM1011），GEDI 的光斑直径约为25 m，沿轨的光斑间距为60 m，相邻地面轨迹的间距为600 m，如图2 所示[16]。

图2 GEDI 波束采样模式Fig.2 GEDI beam sampling mode

表1 GEDI 参数[16]Table 1 GEDI parameters

GEDI 分为4 级产品，第一级别为原始波形数据（L1A）和定位波形数据（L1B），第二级别为冠层高度和剖面指标（L2A 和L2B），第三级别为网格冠层指标（L3），第四级别为网格生物量指标（L4）。本研究使用GEDI 的L2A 高程产品和L4A 生物量产品，L2A 级地理定位高程和高度度量产品由土地过程分布式主动档案中心（Land processes distributed active archive center，LPDAAC）分发，目的是提供波形解释，并且从每个L1B 接收到的波形中提取产品，包括地面高程、冠层高度和相对高度指标。L4A 数据产品来自橡树岭国家实验室分布式活动档案中心（Oak ridge national laboratory distributed active archive center，ORNL DAAC），包含单个GEDI 足印水平的地上生物量密度和相关预测不确定性。L2A 高程产品和L4A 生物量以HDF5 格式提供，空间分辨率为25 m，数据均可在https://search.earthdata.nasa.gov/search?q=GEDI 查询下载。

低能量的地面反射或较高的背景噪声等不利情况会使生产的GEDI L2A 高程数据存在误差，为了减小不利情况的影响，对GEDI L1B 的波形数据使用高斯平滑滤波器进行平滑处理，如公式（1）所示。

对L1B 的波形设定不同的噪声阈值、信号阈值、起始阈值和结束阈值，首先根据噪声阈值对原始波形的噪声进行拟合去噪，然后对拟合后的波形设置信号阈值检测是否存在误差。对于波形的起始信号和结束信号分别设置阈值，以便能够检测并改善虚假信号，其中包括设置阈值时在冠层顶部或地面以下检测到噪声，就会造成冠层高度过大或者地面高程过低[17]。为了减小上述情况的影响，并考虑白天夜晚的数据差异，设置6 组不同的算法阈值。根据不同的算法结果，生产了6组不同高程数据，存储在geolocation 子组中，具体的噪声、信号、波形起始信号和结束信号的阈值如表2 所示，σ是背景噪声的标准差，a1～a6表示6 个算法组，其中a1 是默认算法组，能满足大多数情况，可能遇到的其他情况包括：

表2 L2A 算法设置组[17]Table 2 L2A algorithm setting groups

1）GEDI 获取的最低高程低于真实地面，则需要用更高信号的结束阈值算法；

2）GEDI 的最低高程位于真实地面，但地面能量较弱，因而需要选择更低信号的结束阈值算法；

3）GEDI 最高回波低于顶部冠层，则需要更低信号的开始阈值算法；

4）GEDI 最高回波高于顶部冠层，具有更高信号的开始阈值算法可能更合适[18]。

1.2.2 G-LiHT 数据

美国宇航局戈达德的激光雷达，高光谱和热机载成像仪（NASA Goddard’s lidar,hyperspectral and thermal airborne imager，G-LiHT）是便携式的机载成像系统，以极高的空间分辨率同时测量植被结构、叶面光谱和表面温度[19]。G-LiHT 数据从0～4 级，本研究主要使用了由Riegl VQ-480 Scanning Lidar 生成的L3 级冠层高度模型（Canopy height model，CHM），其分辨率是1 m，数据可以通过https://glihtdata.gsfc.nasa.gov/ 获取。森林生物量数据集同样是G-LiHT 搜集，通过使用LiDAR 数据以及来自PEF 内604 个森林调查清单和分析图进行建模，在样地位置包含多个时间索引测量值来估计AGB 的变化，从而调整现场测量值和遥感数据之间的时间偏差，从而得出生物量[20]。该数据可以通过https://daac.ornl.gov/cgi-bin/dataset_lister.pl?p=33 查询下载。

1.3 研究方法

分别对于GEDI 的L2A 和L4A 数据进行不同的处理，通过L2A 产品可以根据高程数据反演森林最大冠层高度，L4A 产品提供森林生物量数据。具体研究方法流程如图3～4 所示，详细的研究过程包括L2A 获取相关字段、提取冠层高度、机载数据验证精度、L4A 获取相关字段、提取生物量和精度验证。

图3 GEDI L2A 反演最大冠层高度流程图Fig.3 Flow chart of GEDI L2A inversion of maximum canopy height

1.3.1 GEDI L2A 数据预处理

本研究利用GEDI L2A 高程数据，L2A 产品包括8 个波束和1 个元数据（METADATA），元数据集中记录了获取时间、目的、文件名、创建时间等综合信息。每个波束集中包含156 个字段，包括地理位置、陆地覆盖和数据处理等信息。

本研究主要提取了10 个字段，首先根据shot_number、lon_lowestmode和lat_lowestmode用来确定每个足印点的唯一索引及经纬度，elev_lowestmode和elev_highestreturn 对应着最高高程和最低高程；mean_sea_surface 是平均海平面，在进行高程数据处理时，需要在原有数据基础上根据这个参数进行修正。相对高度（Relative height，RH）定义为小于给定高度的接收激光波形强度百分比，其中给定高度是相对于波形中最低高度计算的[21]，如图5所示。最后quality_flag、degrade_flag 和sensitivity 是选择高质量足印的筛选指标，利用这3 个参数进行数据的筛选，quality_flag 有0 和1 两个值，1 表示质量良好，需要剔除0 的激光点；degrade_flag 同样有0 和1，含义不同，1 表示卫星处在降轨中，数据不准确，因而要保留0 的光斑点；sensitivity 表示的是冠层覆盖度，其取值范围为0～1，越接近1 表示质量越高，在森林区将其设置为sensitivity ≥0.95。这10 个字段具体的名称及其描述如表3 所示。

图4 GEDI L4A 反演生物量流程图Fig.4 Flow chart of GEDI L4A inversion of biomass

图5 L2A 相对高度示意图[16]Fig.5 Schematic diagram of L2A relative height

表3 GEDI L2A 字段介绍[21]Table 3 Introduction to GEDI L2A fields

经过筛选之后的高质量激光点，elev_lowestmode表示最低高程可作为地面高程数据，elev_highestreturn 表示最高高程，在林区内即可作为冠层顶部高度。最高高程减去最低高程可得到利用GEDI 高程数据推算的森林冠层高度，同时需要根据mean_sea_surface 进行大地水准面的校正，得到最大冠层高度。

对于GEDI 的RH 参数，选择精度最高的算法组，直接从该组中提取RH100，即最大冠层高度。为了分析该研究区最合适的相对高度，提取不同的RH，判断恰当的RH 反演森林最大冠层高度。

1.3.2 GEDI L4A 数据预处理

GEDI L4A 数据产品是森林生物量产品，主要包括森林生物量预测子集、地理定位子集、陆地覆盖子集和其他在根目录下的数据字段。L4A 算法使用L2A 相对高度指标作为参数线性模型的输入来预测AGB。对应于L2A 的6 组高程数据，L4A同样提供了6 组生物量数据，但相较于L2A，L4A中多了一个a10 组，a10 的设置是因为使用a5 检测到的最低信号精度较低，所以需要更高的信号来计算相对高度，从而获取更高精度的生物量[18]。

本研究共提取了11 个字段，根据shot_number、lon_lowestmode 和lat_lowestmode进行光斑的定位，然后根据l2_quality_flag、l4_quality_flag、algorithm_run_flag、degrade_flag、sensitivity进行光斑的筛选，degrade_flag、sensitivity 与L2A 的含义相同，algorithm_run_flag 表示算法运行标志，如果此标志设置为 1，则运行 L4A 算法。此标志选择具有足够波形保真度以进行生物量估计的数据。当algorithm_run_flag=1 且sensitivity ＞0.9时，l2_quality_flag=1；l4_quality_flag=1表示l2_quality_flag=1 且sensitivity ＞0.95。筛选后，将标识地上生物量的字段agbd、agbd_pi_lower、agbd_pi_upper 提取出来，然后进行生物量的反演，所提取的字段名称及其含义记录在表4 中。

表4 GEDI L4A 字段介绍[22]Table 4 Introduction to GEDI L4A fields

1.4 精度评价标准

为了评价GEDI 反演森林最大冠层高度和生物量的精度，研究选择决定系数（R2）、均方根误差（RMSE）、绝对平均误差（MAE）作为评价指标，3 个指标的计算公式如（2）～（4）所示。其中，RMSE 和MAE 的范围都是大于0，但求出的结果往往是RMSE 略大于MAE，这是因为RMSE 是先对误差进行平方的累加后再开方，其实是放大了与较大误差之间的差距。而MAE 反映的就是真实误差。计算结果越接近于0，其意义越大，说明估测精度越高。R2取值0～1 之间，表示两种数据的线性关系，越接近于1说明两种数据相关性越高，模型的精度越高。

式中：n是样本数量，yi是第i个真实值，是第i个预测值，表示真实值的平均值。

2 结果与分析

2.1 GEDI 反演最大冠层高度

如图1 所示，研究区的光斑点共104 个，为探究不同算法组的GEDI L2A 数据反演森林最大冠层高度的精度，将6 组数据与机载数据进行验证，统计结果如表5 所示，GEDI L2A 与G-LiHT 反演最大冠层高度的散点图如图6 所示。

图6 不同算法组GEDI 数据与G-LiHT 反演森林最大冠层高度散点图Fig.6 Scatterplots of the max canopy height inversion of GEDI data of different algorithm groups and G-LiHT

表5 GEDI L2A 反演最大冠层高度精度统计Table 5 Accuracy statistics of the maximum canopy height retrieved by GEDI L2A

试验结果表明，6 个算法组数据R2的均值为0.81，RMSE 和MAE 的均值分别为2.54 和1.63 m。精度高于均值的分别是a1、a2、a3 和a4，a1 和a4的R2的精度分别是0.96 和0.97，RMSE 分别是1.15 和0.87 m，MAE 为0.50 和0.31 m，其中，a4组是精度最高的。a2 和a3 的精度也高于均值，a2的精度是R2=0.87，RMSE=1.71 m，MAE=0.84 m，a3 的精度为R2=0.85，RMSE=2.04 m，MAE=1.02 m。a5 和a6 的精度低于均值，精度最低的是a5 组，R2=0.44，RMSE=6.39 m，MAE=5.33 m，a6的精度分别为R2=0.77，RMSE=3.11 m，MAE=1.79 m。因此，6 组数据中，可以选择a4 组作为反演森林最大冠层高度的算法组。

为了研究相对高度反演冠层高度的精度，选择a4 组作为研究数据，提取了RH100、RH95、RH90、RH85 此4 组数据，其中RH100 即最大高程与最低高程之差，RH100 的精度即R2=0.97，RMSE=0.87 m，MAE=0.31 m。RH100、RH95、RH90、RH85 反演最大冠层高度的R2逐渐降低，RMSE 逐渐增加，MAE 呈现下降的趋势，但RH95 的MAE 最小，因此假设反演冠层高度精度最高的相对高度是在RH100～RH95。为了进一步确定精度最高的相对高度字段，提取RH99、RH98、RH97、RH96，最终相对高度反演森林最大冠层高度的精度估计结果如表6 所示。

表6 相对高度反演最大冠层高度精度统计Table 6 Accuracy statistics of the maximum canopy height retrieved from relative height

RH100～RH95 的结果有很小的差距，综合考虑R2、RMSE 和MAE 的结果，反演冠层高度的精度差异较小，可以作为GEDI 反演冠层高度的相对高度指标。在此区间内，RH99 的反演精度最高，R2=0.97，RMSE=0.85 m，MAE=0.25 m。从RH99～RH95，精度逐渐降低。R2从0.97 降低到0.96，RMSE 和MAE 分别从0.85 和0.25 m 增加到1.07 和0.30 m。因此，在该研究区可以判断RH99 可能是反演冠层高度最合适的RH 指标。

2.2 GEDI 反演森林生物量

为探究不同算法组的GEDI L4A 数据反演森林最大冠层高度的精度，将7 组数据与机载数据进行验证，统计结果如表7 所示，GEDI L2A与G-LiHT 反演生物量的散点图如图7 所示。试验结果表明，a1～a6 这6 个算法组样本数量的平均值（Mean1）是82，R2的均值为0.89，RMSE 和MAE 的均值分别为14.95 和8.09 t/hm2。除a1 组的精度低于均值，其精度分别为R2=0.80，RMSE=16.07 t/hm2，MAE=9.40 t/hm2。其他5 组数据的精度都高于均值，其中，a2 的精度是最高的，分别是R2=0.92，RMSE=10.80 t/hm2，MAE=4.99 t/hm2，同时a2 也是采集光斑点最多的算法组。

表7 GEDI L4A 反演森林生物量精度统计Table 7 Accuracy statistics of the forest biomass retrieved by GEDI L4A

图7 不同算法组GEDI 数据与G-LiHT 反演生物量散点图Fig.7 Scatterplots of the biomass inversion of GEDI data of different algorithm groups and G-LiHT

续图7Continuation of Fig.7

根据表7 和图7 可知，a5 算法组反演最大冠层高度精度最低，且远低于均值，这可能是因为最低检测信号是噪声。为此，在进行生物量反演时，选择更高的检测信号阈值的数据来代替原来的a5数据，从而生成了额外的a10 组，a10 组的样本数量有91，相较于a5 组的47 个，多了44 个光斑点，样本数量大幅提高，反演森林生物量的精度为R2=0.91，RMSE=15.21 t/hm2，MAE=7.61 t/hm2，同样高于a5 组的精度。

取a1～a4，a6，a10 此6 组数据计算均值，记为Mean2，样本数量和反演精度都有所提升，其中，光斑平均数量多了7 个，R2保持不变，RMSE、MAE 分别提高了0.90、1.85 t/hm2，这表明a10 算法组在获取最终有效光斑数量和反演生物量的精度上，是优于a5 算法组的，a10 组可以有效改善a5 组数量少、精度低的不足，能够有效提高GEDI 反演森林生物量的精度。在此6 组数据中，精度最高的仍是a2 算法组所提取的。所以，7 组数据中，a2 是光斑点最多、精度最高的算法组，因此在数据验证或者生产时可以选择a2 组作为反演森林生物量的参考数据。

3 结论与讨论

为探究GEDI 数据反演森林最大冠层高度和生物量的精度，以美国的佩诺布斯科特实验森林作为研究区，以G-LiHT 的产品CHM 产品和衍生的生物量数据作为验证数据，分析了GEDI L2A 的6个算法组高程数据反演最大冠层高度的精度。同时针对L4A 生物量产品7 个算法组数据计算了估算森林生物量的精度，得出以下结论：

1）在GEDI L2A 的6 个算法设置组中，反演精度最低的是a5 组，精度最高的是a4 组，其精度是R2=0.97，RMSE=0.87 m，MAE=0.31 m。利用GEDI L2A 的相对高度指标反演冠层高度，可以得出RH99 是该研究区反演精度最高的，精度为R2=0.97，RMSE=0.85 m，MAE=0.25 m。相较于RH100，RH99 的RMSE 和MAE 分别提高了0.02和0.06 m。

2）GEDI L4A 数据集中设置有7 个算法组，多设置的第7 组a10 是用来替代a5 组数据的，改善了a5 数据量小、精度低的不足。相较于a5组，a10 组的光斑数量多了44 个，R2、RMSE、MAE分别提高了0.01、5.50 t/hm2、5.10 t/hm2。7 个算法组中，精度最高的是a2 组，R2=0.92，RMSE=10.80 t/hm2，MAE=4.99 t/hm2。

利用星载激光雷达反演森林最大冠层高度和生物量会受到低能量的地面反射或较高的背景噪声等不利情况的影响。为了减少这些不利情况的影响，针对不同的冠层和地面的场景，GEDI 设定了不同的噪声阈值、信号阈值、起始阈值和结束阈值共6 种组合供用户选择。在冠层高度方面，除a5 算法组数据外，其余5 组数据的精度明显优于Potapov等[23]得出的结果（RMSE=6.60 m;MAE=4.45 m,R2=0.62）。Liu等[24]推荐使用夜间的GEDI L2A 数据，因为RMSE 从日间5.02 m 降到3.56 m，但其精度也弱于除a5 组数据以外其余5 组数据的冠层高度反演精度。在生物量方面，本研究计算出的RMSE 结果中，a2 组和a3 组的精度优于Dorado-Roda等[12]得出的生物量区间14.13～32.16 t/hm2，其余5 组数据精度位于该生物量区间内，因此在选择L4A 数据时可以优先选择a2 和a3 算法组进行试验研究。由于本文所选研究区地形坡度较小，比较平缓，不能很好地反映坡度对GEDI 数据的影响。同时，森林类型和植被覆盖率也是影响反演精度的重要因素。因此，下一步的研究重点是针对不同森林类型、不同地表坡度和不同植被覆盖率对于GEDI 数据反演森林最大冠层高度和生物量的影响。