基于Landsat 8 OLI的滇西北高山亚高山针叶林地上生物量遥感估测模型构建及光饱和点确定

吴 勇，张海燕，卢腾飞，胥 辉，欧光龙

(西南林业大学 西南地区生物多样性保育国家林业局重点实验室，云南 昆明 650224)

森林生物量是反映森林生态系统变化规律的基础数据，在维持全球碳氧平衡、调节生态平衡等方面起着极其重要作用，这使得森林生物量估测方法的可靠性和精度变得极其重要[1-2].传统的生物量估测方法成本高、工作量大、周期长、代表性差以及不能实时反映实地情况等，因此众多学者开始使用遥感技术对森林生物量进行估测和模型模拟[3-6].遥感估测具有实时、动态、大面积同步监测和信息丰富等特点[7-9]，被众多学者广泛应用于不同地区、不同林分的生物量遥感估测中.近年来在使用Landsat 8 OLI 影像估测时生物量光饱和点的问题普遍存在.由于森林生态系统自身及分布区地形变化，遥感生物量估测模型、遥感数据源及其空间分辨率等都会使得森林生物量饱和值产生误差，估算误差范围多在5%～30%[10-11].生物量饱和现象是导致估测不确定性的重要因素，不同树种组成、年龄、空间和垂直冠层结构组成的植被类型反射率值不同，也会影响生物量光饱和点[12-13]，这使得如何提高生物量估测精度[14]成为研究的热点话题.

近年来，诸多学者对如何解决森林生物量数据饱和带来的问题进行了研究，付元元等[15]利用高光谱数据，将波段深度分析与最小二乘回归结合，结果能较好地克服生物量较大时存在的饱和问题，提高遥感估测生物量的精度.欧光龙等[16]、卢腾飞等[17]通过构建不同模型对森林生物量遥感进行估测以探索减少饱和现象造成的估测精度问题.赵盼盼等[18]、Zhao 等[19]基于Landsat TM 影像与地面样本数据，利用分层理论和曲线拟合确定浙江省南部区域的6 类树种的生物量估测所存在的饱和点的大致范围，以及通过分层法降低了数据饱和对浙江中西部主要植被生物量的影响，研究表明分层方法能有效提高估测精度.

滇西北地区属中国横断山脉，森林植被呈现明显的垂直分布，生长着大量云杉、冷杉、高山松、落叶松等.高山松（Pinus densata）为中国西部高山地区的特有树种，其主要分布在北纬27°～36°、东经92°～104°之间；云杉（Piceaspp.）和冷杉（Abiesspp.）均为高海拔地区特有树种，高山松、云杉、冷杉在滇西北广泛分布[20].目前对相关优势树种森林生物量光学遥感估测饱和点的研究及报道极少.

基于此，本文将以滇西北高山亚高山针叶林为研究对象，以森林资源二类调查数据和同期Landsat 8 OLI 遥感影像为数据源，利用分层理论构建遥感生物量估测模型，探索典型森林AGB 估算中的光学遥感数据源的饱和点阈值确定方法，并反演研究区高山松林和云冷杉林森林地上生物量.本研究可以实现对滇西北高山亚高山针叶林生物量的精确估测，减小森林生物量遥感估测中普遍存在的低值高估和高值低估问题的影响，提高生物量遥感估测精度，为滇西北森林生态效益、碳汇效益以及 科学经营与规划提供科学依据.

1 材料与方法

1.1 研究区概况研究区位于云南省西北部，在行政区划上涵盖丽江市、迪庆州、怒江州3 个州市，涵盖面积4.2 万km2.地理位置为北纬25°30′～29°30′，东经98°7′～101°30′，海拔高度在738～6 743 m，海拔绝对高差悬殊，造成气候显著差别，垂直气候变化十分明显，拥有南亚亚热带至寒带等多种气候类型（图1）.滇西北地区森林资源极其丰富，高山松、云杉、冷杉等作为滇西北地区优势树种在滇西北地区广泛分布，森林面积占全云南省森林面积的18.7%，森林生物量占全云南省森林生物量的24.2%，森 林覆盖率高达71.8%[20].

图1 研究区位置及优势树种分布图Fig.1 Location of the study area and the dominant species

1.2 数据收集与处理

1.2.1 森林资源二类调查数据收集与处理

（1）小班单位面积AGB 计算 本研究依据2016 年研究区的森林资源二类调查数据，并选择优势树种为高山松和云冷杉的林分（图1），参考胥辉等[20]采用生物量转换因子法计算小班单位面积AGB[21]（表1）.其计算公式为：

表1 生物量转换因子法计算参数[21]Tab.1 The parameters using biomass conversion factor method

式中：B为小班单位面积地上生物量（t/hm2），V为小班单位面积蓄积量（m3/hm2），SVD 为基本木材密度（t/m3），BEF 为生物量转换因子（无量纲）.

（2）样本小班的确定 首先据二类调查数据提取出优势树种为高山松和云冷杉的小班斑块，优势树种在研究区分布如图1 所示，剔除小班平均树高、平均胸径、林木蓄积和林木株树为零的异常数据，筛出林木平均胸径大于等于5 cm 的小班.利用ArcGIS 的“子集要素”工具，等比例保留小班数据，再利用3 倍标准差法删除异常离群值，最终选出小班样本高山松1 094 块，云冷杉948 块，分别取80%小班进行建模，20%小班作为检验样本（表2）.

表2 样本小班单位面积地上生物量基本统计Tab.2 The basic information of the sampling sub-compartments

本研究选取坡向与龄组用不同分层方法（无分层、坡向分层、龄组分层以及坡向龄组结合）构建遥感生物量估测逐步回归模型，并统计计算了各分层 下的森林单位面积地上生物量值（图2）.

图2 分层数据统计Fig.2 Data analysis according to the different stratification

1.2.2 遥感数据收集与处理 本研究遥感数据（与二类调查数据同期）从地理空间数据云网站（http://www.gscloud.cn/）下载得到，其空间分辨率为30 m，平均云量都低于2%（表3）.对遥感影像数据进行地理配准、FLAASH 大气校正、地形校正、辐射定标等预处理，最后对研究区卫星影像进行图像 融合、镶嵌及裁剪（图1）.

表3 研究区Landsat 8 OLI 影像基本信息Tab.3 The basic information of Landsat 8 OLI images in study area

1.2.3 遥感因子提取与筛选 根据高山松林、云冷杉林小班面状矢量数据为单位，利用ArcGIS10.4软件的“分区统计”功能，将小班面状矢量数据和156 个遥感特征因子数据层相叠加，统计小班内各个遥感因子的平均值[22].光学遥感估测AGB 的关键在于确定自变量因子，一般采用统计检验的方法来确定影响因子[23].本研究提取了4 类遥感因子，分别为原始单波段因子、植被指数因子、K-T 和KL 因子、纹理特征.

（1）原始单波段[10]为b1、b2、b3、b4、b5、b6、b7；

（2）植被指数[22]为NDVI、ND43、ND67、ND563、DVI、SAVI、RVI、PVI、B、G、W、ARVI、MV15、MV17、SARV、II、MSAVI、VIS234、ALBEDO、b4/ALBEDO、SR、TVI、SAV12、NLI、MSR；

（3）K-T 和K-L 因 子[22]为KT1、KT2、KT3；PC1-A、PC2-A、PC3-A、PC1-B、PC2-B、PC3-B、PC1-P、PC2-P、PC3-P；（4）纹理特征[22]（提取分为3×3 和5×5 窗口，共计112 个遥感因子）：Mean、Variance、Homogeneity、Contrast、Dissimilarity、Entropy、Second Moment、C orrelation.

1.3 研究方法

1.3.1 相关性分析 本研究采用的是皮尔逊相关性分析，通过SPSS 软件对AGB 与遥感因子进行相关性分析，选择相关性显著（P≤0.05）的遥感因 子作为建模备选变量.

1.3.2 AGB 遥感估测的生物量饱和点确定 本文以生物量为自变量，以波段反射率为因变量，基于卢腾飞[17]、刘刚[9]等的研究发现，三次项模型拟合曲线方程的决定系数R2最高，故本研究采用三次项模型拟合曲线并求得其拐点的值，其拐点所对应的值即为云冷杉林与高山松林遥感估测AGB 饱和 值.

1.3.3 生物量估测模型构建 在森林生物量遥感估测模型中，线性回归模型是最常见的方法之一[7]，本文采用线性逐步回归进行模拟，其计算公式如下：

式中：Y为生物量，β0为常数项，β1,β2,···,βn为模型系数，x1,x2,···,xn为相关遥感因子，ε 为模型满足随机正态分布下的随机残差，n为自变量个数.

分层模型采取分别拟合不同坡向、不同龄组、坡向结合龄组3 种分层方式对生物量进行逐步回归分析.龄组分为幼龄林、中龄林、近熟林、成熟林和过熟林，由于云冷杉林幼龄林和中龄林样本数极少，将其就近归入近熟林；高山松林中过熟林样本数极少，将其归入成熟林.划分坡向时一般根据坡向所受日照时长和太阳辐射强度，将坡向划分为阳坡（南坡）、半阴坡（东坡、东北坡、西北坡）、阴坡（北坡）、半阳坡（西坡、东南坡、西南坡）[23]，本研究区处于高山峡谷地区，无无坡向小班.坡向和龄组结合分为幼龄林半阳坡/阳坡/半阴坡/阴坡、中龄林半阳坡/阳坡/半阴坡/阴坡、近熟林半阳坡/阳坡/半阴坡/阴坡、成熟林半阳坡/阳坡/半阴坡/阴坡、过熟林半阳坡/阳坡/半阴坡/阴坡.

以优选出的遥感因子作为解释变量，采用线性逐步回归法分别构建高山松林和云冷杉林不分层及3 种分层的森林生物量遥感估测模型，为消除多重共线性的影响，采用方差膨胀因子（VIF＜10）进行 变量筛选.

1.3.4 模型评价与检验 本研究采用调整决定系数（Radj2）和均方根误差（RMSE）评价拟合情况，选取相对均方根误差（rRMSE）、平均相对绝对误差值（MARE）和平均绝对误差（MAE）3 个指标作为模型的独立性检验指标.调整决定系数（Radj2）：

1.3.5 不同生物量段模型估测能力评价 本文对小班面积单位生物量进行分段，检验各模型对生物量高值和低值的估测能力，着力于探索改善数据饱和引起的生物量估测误差方法，采用的评价指标公式如下：

平均残差（ME）：

1.3.6 高山松林、云冷杉林生物量反演 依据森林资源二类调查数据获取高山松林以及云冷杉林的空间分布位置，基于Landsat 8 OLI 影像采用像元法计算各小班内云冷杉林及高山松林生物量值，再根据分层所构建的不同模型对高山松林、云冷杉林进行生物量反演，进而分析不同分层模型的反 演和预估表现.

2 结果与分析

2.1 生物量与遥感因子相关性分析本文共提取156 个遥感因子分别与高山松林和云冷杉林AGB进行皮尔逊相关性分析.从波段反射率与AGB 的相关性看，Band 4 波段与生物量有最强的相关性，高山松林、云冷杉林与Band 4 波段相关系数分别为0.616、0.672，因此本研究基于Band 4 波段进行分析和确定高山松林、云冷杉林的生物量饱和值.据相关性分析发现，有141 个遥感变量与高山松林生物量之间存在显著的相关性，包括107 个纹理特征因子，16 个植被指数因子，11 个缨帽特征和主成分变换特征，7 个原始单波段，其中135 个遥感变量在0.01 水平上显著相关.有146 个遥感变量与云冷杉林生物量之间存在显著的相关性，包括112 个纹理特征因子，16 个植被指数因子，11 个缨帽特征和主成分变换特征，7 个原始单波段，其中144 个遥感变量在0.01 水平上显著相关.由于遥感因子变量数量较多，限于篇幅原因，本研究不列出遥 感因子与生物量相关性系数表.

2.2 生物量光学遥感估测饱和值分析本研究基于Band 4 波段反射率与森林地上单位生物量拟合高山松林、云冷杉林生物量饱和曲线（图3），其拐点所对应的值即为该波段所对应的高山松林、云冷杉林生物量饱和值，分别为149.09 t/hm2和1 62.30 t/hm2.

图3 研究区针叶林地上部分生物量光饱和值拟合曲线Fig.3 Aboveground biomass saturation curves of the coniferous forests in the study area

2.3 森林生物量遥感估测模型构建及检验本研究采用线性逐步回归基于分层建立生物量回归模型，对遥感特征因子与不同植被类型生物量作相关性分析，以VIF 为检验指标进行多重共线性诊断，筛选出符合目标条件的遥感因子构建线性逐步回归模型.用Radj2评估预测值和实测值的符合程度，用RMSE 衡量预测值与实测值之间的偏差大小.

从回归模型的独立性检验看，高山松林与云冷杉林分层模型检验指标值总体优于无分层模型.其中坡向分层模型各检验指标值与无分层模型相差不大，总体预估水平相近.龄组分层整体估测精度优于不分层模型，高山松林各检验指标数值明显降低，拟合精度显著提高（表4），云冷杉林各检验指标值均 低于无分层，在一定程度上提高了拟合精度（表5）.

表4 高山松林线性逐步回归拟合结果Tab.4 Fitting results of Pinus densata forests using linear stepwise regression

表5 云冷杉林线性逐步回归拟合结果Tab.5 Fitting results of Picea -Abies forests using linear stepwise regression

2.4 不同生物量分段残差检验经分段残差分析发现，2 类森林的4 种模型都存在低值高估和高值低估的情况（表6）.高山松林与云冷杉林线性逐步回归模型在低生物量段（0～30 t/hm2）时，高山松林坡向分层ME（19.724）、MRE（-1.141）与云冷杉林坡向分层ME（10.680）、MRE（-0.609）最低，与无分层的模型比较，ME 绝对值分别降低了4.938 和2.525，MRE 绝对值分别降低了0.234 和0.135.高山松林与云冷杉林龄组分层模型在生物量段（30～90 t/hm2）有较好的表现.坡向龄组结合对于高山松林与云冷杉林高生物量段（＞150 t/hm2）的估测误差值最小，高山松林与云冷杉林高生物量段（＞150 t/hm2）ME绝对值分别降低了10.291 和24.322，MRE 绝对值分别降低了0.05 和0.148，且云冷杉林在生物量段（30～150 t/hm2）的估测误差绝对值均为最小.由此可见，考虑龄组和坡向分层在一定程度上降低了AGB 遥感估测中低值高估和高 值低估问题的影响.

表6 分段残差分析结果Tab.6 Segmented residual analysis

2.5 森林生物量反演由图4、5 可知，滇西北高山松林与云冷杉林生物量主要集中在范围30～150 t/hm2和60～150 t/hm2，少部分处在30 t/hm2以下和150 t/hm2以上.图4 中（b）、（c）、（d）与（a）以及图5 中（b）、（c）、（d）与（a）分别比较后发现，各分层模型生物量反演高值数量明显增加，高山松林与云冷杉林坡向分层模型在低生物量段值（0～30 t/hm2）数量增加，高山松林坡向龄组结合分层模型低生物量值（0～30 t/hm2）和高生物量值（大于150 t/hm2）数量都有所增加，云冷杉林坡向龄组结合分层模型高生物量段（120～150 t/hm2）值数量显著增加.总体来说，分层模型相较于无分层模型在低生物量值（0～30 t/hm2）和高生物量值（大于150 t/hm2）区间具有较多的分布.这使得生物量估测范围增大，在一定程度上能够反映出分层估计可以减小生物量遥感估测中低值高估和高值低估的影响.

图4 高山松林不同建模方案生物量反演图Fig.4 Biomass inversion map of Pinus densata forests with different modeling schemes

图5 云冷杉林不同建模方案的生物量反演图Fig.5 Biomass inversion map of Picea– Abies forests with different modeling schemes

3 讨论与结论

3.1 讨论光学影像对于冠层结构复杂、AGB 值较大的森林，其反射率出现饱和现象往往导致对高生物量值的低估.此外，由于冠层密度较低，幼林的反射率值也可能受到灌木、草、裸地等林下植被的影响，导致低生物量值被高估.目前，数据饱和被认为是影响森林生物量估测精度的主要原因，大量的遥感估测生物量均出现了数据饱和现象[12,24-25].在不同区域、不同类型的森林植被中，因林分因子与地形因子的不同，其生物量饱和值也会有所变化，通常还会受到遥感数据本身空间和光谱辐射分辨率差异的影响[9,20].本文基于滇西北地区森林资源二类调查数据计算AGB 值与遥感影像所提取小班面状数据和遥感因子进行曲线拟合，求得高山松林、云冷杉林生物量饱和值为149.09 t/hm2和162.3 t/hm2；此外，通过无分层模型的分段残差检验也可以看出，高山松林在120～150 t/hm2以及＞150 t/hm2的ME值较高，而云冷杉林生物量在120～150 t/hm2区间的ME 值低于12.772 t/hm2，且MRE 为0.093，但在＞150 t/hm2区间段ME 值达到40 t/hm2以上，这在一定程度上也说明了云冷杉林的生物量光学遥感估测的光饱和值高于高山松林.Ou 等[26]通过LM、LMC、RF 和ANN 4 种模型对高山松林AGB 进行范围估算，并得到了大致的高值低估区间，其区间阈值略高于赵盼盼[10]等通过半变异函数计算的针叶林生物量饱和值.这可能是由于滇西北地区高山松林和云冷杉林及其分布环境的异质性高，从而使得高山松林、云冷杉林的森林生物量遥感估测的光饱和值较高，也符合Lu 等[12]和Zhao 等[19]森林异质性程度在一定程度上影响森林生物量遥感估测的饱和值的解释.

研究表明，分层模型预估效果整体上优于无分层模型.当生物量过低或过高时，这两种模型都存在较大的估测误差，分层模型相对于无分层模型降低了估测误差，增大了生物量估测范围，在一定程度上提高了生物量估测精度，减小了森林生物量饱和所造成的低值高估和高值低估的影响.坡向和龄组作为影响森林生物量分布的两个重要因素[27-29].本研究选用坡向与龄组因子，采用分层法（龄组分层、坡向分层、坡向与龄组结合分层）构建模型以期解决AGB 估测饱和问题.研究发现研究区坡向分层森林生物量均值差异较小（图2），但是由于滇西北地处高山峡谷，地形差异较大，不同坡向光谱反射率不同，使得光谱差异大，从而对森林生物量估测误差产生影响[3].研究表明坡向分层对低生物量段（0～30 t/hm2）有更好的预估效果.高山松林、云冷杉林平均残差（ME）绝对值分别降低了4.938和2.525，平均相对残差（MRE）绝对值分别降低了0.234 和0.135；龄组分层模型在高山松林生物量值30～150 t/hm2与云冷杉林生物量值60～150 t/hm2时预估效果较好；龄组分层后可以在一定程度上提高估测精度[26]，但不同龄组间的森林生物量均值差异较大，使得采用平均模型预估不同龄组森林生物量时产生较大误差[22,30].坡向龄组结合分层模型在高生物量段（120～150 t/hm2）对AGB 估测效果最好，高山松林、云冷杉林ME 绝对值分别降低了10.291 和24.322，平均相对残差（MRE）绝对值分别降低了0.05 和0.148.通过分段残差检验可以看出，分层估计的方法在一定程度上降低了遥感估测AGB 饱和值的影响，获得了较好的拟合和预估效果，该方法对相似地区及针叶林进行遥感生物量估测及降低森林生物量饱和值的影响提供参考和借鉴作用，为滇西北森林生态效益、碳汇效益以及科学经营与规划提供科学依据.此外，研究还发现考虑不同坡向分层，从阳坡到阴坡，高山松林模型决定系数变大，云冷杉林模型决定系数变小，其原因可能是与优势树种的生长特性有关.高山松是阳性树种而云杉和冷杉是阴性树种，树种的耐荫性会影响到其分布区差异[31-32]，这将有待进一步探索和研究.

当然，由于本研究采用的都是参数模型，分层法虽提高了拟合精度，但总体估测精度依然不高.参数之间通常会存在比较复杂的非线性关系，本研究所采用的模型不足以描述其非线性关系，且当样本量不足或森林生物量与选择变量之间的相关性较弱时，多元逐步回归模型不能解释森林生物量与变量之间的复杂关系，从而使得模型估测精度降低[8,33].在今后的研究中将进一步探索基于分层的非线性参数模型与非参数模型估测方法.

另外，由于缺乏足够多的样地信息，尤其是在高海拔、多山地的滇西北地区，不易开展外业，测量森林生物量因子不易获取[34].本文基于滇西北地区森林资源二类调查数据，二类调查数据的统计数据精度与各小班的调查允许误差直接相关[35].利用充分的样本信息以及二类调查数据具有分布范围广、数据量大等优点，采用二类调查数据作为森林生物量遥感估测的地面调查数据可以在一定程度上解决样本数量不足的问题.周律等[9]利用森林资源二类调查数据作为基本地面调查数据，构建遥感生物量估测模型，取得了较好的估计精度，这在一定程度上说明了森林资源二类调查数据作为地面调查数据的可用性.本研究采用二类调查数据作为地面调查数据，构建的估测模型，并进行森林生物量的反演，其反演结果较好地反映了研究区高山松林和云冷杉林的森林生物量分布，这也再次验证了二类调查的小班数据作为森林生物量遥感估测 地面调查数据源的可靠性.

3.2 结论本文通过皮尔逊相关系数比较筛选出Landsat 8 OLI Band 4 波段变量与AGB 进行曲线拟合，求得高山松林与云冷杉林遥感生物量估测饱和值分别为149.09 t/hm2和162.3 t/hm2.通过对研究区高山松林和云冷杉林地上单位面积生物量与提取遥感因子进行皮尔逊相关性及共线性分析，并考虑龄组和坡向分层，据优选出的遥感因子进行线性逐步回归分析，发现分层模型相较于无分层模型有着较好的估测结果.坡向分层对低生物量段（0～30 t/hm2）、龄组分层对中生物量段（30～90 t/hm2）和坡向龄组结合分层模型对高生物量段（＞150 t/hm2）有更好的预估效果，坡向龄组结合分层模型在整体上估测AGB 效果较好.总体来说，考虑龄组和坡向分层可在一定程度上降低AGB 遥感估测中低值高估和高值低估问题的影响.