生态学报  2016, Vol. 36 Issue (13): 4109-4121

文章信息

王新云, 郭艺歌, 何杰
WANG Xinyun, GUO Yige, HE Jie.
基于HJ1B和ALOS/PALSAR数据的森林地上生物量遥感估算
Estimation of forest above-ground biomass based on HJ1B and ALOS/PALSAR remote sensing data
生态学报[J]. 2016, 36(13): 4109-4121
Acta Ecologica Sinica[J]. 2016, 36(13): 4109-4121
http://dx.doi.org/10.5846/stxb201410212062

文章历史

收稿日期: 2014-10-21
网络出版日期: 2015-10-29
基于HJ1B和ALOS/PALSAR数据的森林地上生物量遥感估算
王新云1,2, 郭艺歌1,2, 何杰3     
1. 宁夏大学 西北土地退化与生态恢复省部共建国家重点实验室培育基地, 银川 750021;
2. 宁夏大学 西北退化生态系统恢复与重建教育部重点实验室, 银川 750021;
3. 宁夏大学 资源与环境学院, 银川 750021
摘要: 森林地上生物量的精确估算能够减小碳储量估算的不确定性。为了探寻一种有效地提高森林生物量估算精度的方法,探讨了基于遥感物理模型和经验统计模型估算山地森林地上生物量的方法。首先,基于Li-Strahler几何光学模型和多元前向模式(MFM)进行模型模拟,结合查找表算法(LUT)从多光谱图像HJ1B估算贺兰山研究区的森林地上生物量。其次,采用统计方法建立了2种回归模型:(1)多光谱图像HJ1B进行混合像元分解(SMA),并与雷达图像ALOS/PALSAR进行图像融合建立生物量回归模型;(2)雷达图像ALOS/PALSAR后向散射系数和实测生物量建立了生物量回归模型。用实测数据对3种算法估算结果进行精度验证。研究结果表明:采用几何光学模型和MFM算法估算的森林地上生物量精度最好(决定系数R2=0.61,均方根误差RMSE=8.33 t/hm2P <0.001),其估算地上生物量与实测值一致性较好,估算生物量精度略优于SMA估算的精度(R2=0.60,RMSE=9.417 t/hm2);ALOS/PALSAR多元回归估算的精度最差(R2=0.39,RMSE=14.89 t/hm2)。由此可见,采用几何光学模型和混合像元分解SMA适合估算森林地上生物量,利用这2种方法进行森林地上生物量遥感监测研究具有一定的应用潜力。
关键词: 森林     地上生物量     环境卫星     ALOS/PALSAR     多元前向模式(MFM)     混合像元分解(SMA)    
Estimation of forest above-ground biomass based on HJ1B and ALOS/PALSAR remote sensing data
WANG Xinyun1,2, GUO Yige1,2, HE Jie3     
1. State Key Laboratory Breeding Base of Land Degradation and Ecological Restoration of Northwest China, Ningxia University, Yinchuan 750021, China;
2. Key Laboratory for Restoration and Reconstruction of Degraded Ecosystem in Northwestern China of Ministry of Education, Ningxia University, Yinchuan 750021, China;
3. School of Resources and Environment, Ningxia University, Yinchuan 750021, China
Abstract: Forest is one of the most widely distributed terrestrial ecosystems on earth. Global-scale biomass estimation has become a research hotspot. It is important to accurately estimate the spatial distribution of forest above-ground biomass (AGB) because its carbon budget forms part of the global carbon cycle and ecosystem productivity. Remotely sensed data have been widely used to quantitatively obtain the biophysical characteristics of vegetation canopy structure. The use of optical and microwave remote sensing in combination with field measurements can provide an effective method to improve the estimation of forest biomass over large regions. In order to improve the accuracy of estimating forest above-ground biomass from remotely sensed data, the methods for obtaining AGB data using a physically-based canopy reflectance model inversion approach and two other empirical statistical regression methods were introduced in this paper. A geometric-optical canopy reflectance model was run in multiple-forward mode (MFM) using multispectral HJ1B imagery to derive forest biomass at the Helan Mountain Nature Reserve region in the northwest of China. Structural parameters of the forest inventory were carried out in 50 separate 30 m by 30 m randomly distributed plots, and the data was used for either model development or validation. The two other empirical-statistical models were also established to estimate the biomass in the area. A multiple stepwise regression model was developed to estimate the forest above-ground biomass by integrating the field measurements of 30 sample plots with ALOS/PALSAR Synthetic Aperture Radar (SAR) backscatter remotely sensed data. The pre-processing of the HJ1B scenes included radiometric calibration, atmospheric correction, and georeferencing. Radiometric data were converted from radiance to reflectance. Additionally, spectral mixture analysis (SMA) was applied to decompose a mixture of spectral components of HJ1B into vegetation, soil, and shade fractions. The vegetation fraction image was fused with PALSAR data using the discrete wavelet transform (DWT) method. As a comparison, a regression model was also created by integrating field measurements with the fused image. Error levels for the three models and the field-measured data were analyzed. MFM predictions of AGB from HJ1B imagery were compared with the results from the SMA and PALSAR multiple stepwise regression models. Simultaneously, the estimation biomass using the three methods was evaluated for 20 field validation sites. The result shows that a good fit can be found between the AGB estimated by geometric-optical canopy reflectance model and the field-measured biomass with a R2 (coefficient of determination) and RMSE (root mean-square error) of 0.61 and 8.33 t/hm2, respectively. MFM provided the lowest error for all validation plots and its estimated accuracy is a little better than that of the SMA model (R2=0.60, RMSE=9.417 t/hm2). PALSAR multiple stepwise regression model has the worst estimation accuracy (R2=0.39, RMSE=14.89 t/hm2) and had a higher error. Consequently, it can conclude that geometric-optical canopy reflectance model and spectral mixture analysis (SMA) approach were considerably more suitable for estimating the forest biomass in mountainous terrain. Moreover, it demonstrates a good potential for monitoring the indicators of forest ecosystem by combined with the optical and polarimetric SAR remote sensing synergistic research.
Key words: forest     above-ground biomass (AGB)     HJ1B     ALOS/PALSAR     MFM     SMA    

森林是地球上分布范围最广的陆地生态系统类型之一,区域乃至全球尺度的森林地上生物量遥感估算一直是生态学研究的重要内容。森林作为陆地生态系统最大的碳库和土壤碳库的最重要输入源,因其巨大的碳储量和固碳能力而备受关注[1-2]。森林地上生物量的研究有助于对区域乃至全球的碳循环深入地理解和认识[3]。森林地上生物量不仅是森林固碳能力的重要标志,也是森林生态系统服务功能和生态资产评价的重要参数。利用遥感技术可以获取区域乃至全球尺度的陆表参数。光学、微波遥感数据和野外观测数据集成应用是提高区域尺度森林地上生物量估算精度的一种有效方法[4-5]。Basuki等分别采用混合像元分解、离散小波变换和Brovey变换将Landsat-7 ETM+和ALOS/PALSAR进行图像融合来估算热带森林地上生物量,结果表明采用光学和微波图像融合的方法能够提高生物量估算的精度[4]。Treuhaft等认为基于森林结构剖面的生物量估算比微波或光学遥感的方法更精确,通过融合高光谱和雷达数据所测得的叶面积指数和叶面积密度来提高生物量估算的精度[6]

由于地形对卫星接收的反射信号的影响,经验统计模型估算的森林地上生物量精度受到很大影响。利用植被指数和生物量建立的统计关系受土壤等下垫面背景的光谱特征影响比较明显[7]。统计模型只能用于特定的研究区域和数据源,同时植被冠层的反射特性也受到外界条件的影响,利用植被的光谱特征和植被结构参数建立的统计关系不适合在区域尺度和多传感器之间应用。为了克服下垫面背景的空间异质性对经验统计关系的影响,利用遥感物理模型估算森林生态参数是一种可行的方法[8]。遥感物理模型通过输入林分几何结构、地形特征、纯像元的光谱信息来模拟图像像元的光谱特性,输出林分的生态参数。Soenen等人采用几何光学模型和查找表方法从遥感图像SPOT-5估算了加拿大卡南纳斯基斯市落基山的森林参数,研究表明几何光学模型适合对复杂地形条件下的森林植被进行遥感监测研究[9]

贺兰山是我国西部温带草原与荒漠的分界线和连接青藏高原、蒙古高原及华北植物区系的枢纽[10-11]。国内研究人员对贺兰山的森林生长状况、群落特征、生物多样性在海拔梯度上的变化规律等进行了深入研究[12-15]。穆天民采用生态学统计的方法对贺兰山地区的青海云杉林地上生物量进行了研究,结果表明乔木层的生物量达到127.26 t/hm2,青海云杉林的地上生物量不亚于东北地区的云冷杉和落叶松林[16]。近年来很少见到利用遥感技术研究该地区森林分布状况、水文、生物量、碳循环、气候变化和生态资产评价等内容的相关报道。由于贺兰山受不合理土地利用和气候变化的影响,生态系统极其脆弱。因此,利用遥感技术对贺兰山的森林生态系统的生态参数和服务功能进行深入系统的研究尤为必要。

本文选择位于干旱半干旱地区的贺兰山自然保护区的乔木林(主要为针叶林、阔叶林和混交林)为研究对象,建立基于光学遥感HJ1B和微波遥感ALOS/PALSAR等多源遥感数据的森林地上生物量反演模型,并进行森林地上生物量制图,为进一步对该区域的森林生态系统服务功能和生态资产进行综合评价提供理论依据。

1 研究区与数据处理 1.1 研究区概况

本文的研究区选择位于宁夏和内蒙古自治区交界的贺兰山自然保护区境内,中心坐标为38°43′47.9″N,105°55′7.4″E(图 1)。贺兰山地处内陆,属于典型的大陆性气候,冬季严寒,夏季干燥炎热,降水量稀少,年平均气温8.2—8.6℃,年平均降水量250—438 mm,年平均蒸发量1600—1800 mm之间[15]。主要的乔木树种有青海云杉(Picea crassifolia)、油松(Pinus tabulaeformis)、山杨(Populus davidiana Dode)、灰榆(Ulmus glaucescens Franch)和杜松(Juniperus rigida)。灌木树种主要有小叶金露梅、蒙古扁桃、蒙古绣线菊、华西银露梅和忍冬。

图 1 研究区地理位置及样地在图像(HJ1B)上的分布 Fig. 1 Geographical location of the study areas
1.2 遥感图像数据的获取与预处理

研究中使用的图像数据包括HJ1B-CCD1多光谱图像和ALOS/PALSAR雷达图像数据。光学图像数据HJ1B—CCD1获取时间为2010年7月20日,卫星轨道17/68,图像分辨率为30 m。PALSAR是ALOS卫星上携带的一个L波段的主动式合成孔径雷达传感器,工作波段为L波段。原始图像为level 1.1双极化数据,极化方式为HH、HV极化,其中HH表示水平极化(horizontal polarization),HV表示交叉极化(cross polarization)。入射角为38.722°,数据类型为单视复数据(SLC),像元大小为9.368 m × 3.164 m,成像时间分别为2010年6月30日和2010年8月15日。

HJ1B-CCD1图像的预处理包括几何精校正、辐射定标和大气校正。利用HJ1B数据处理工具将灰度值(DN)值转为辐亮度值,利用野外采集的研究区GPS控制点在ENVI4.8中对影像进行几何精校正和大气校正,将辐亮度图像校正为反射率图像。

ALOS/PALSAR图像的预处理采用GAMMA软件进行,利用分辨率为90 m的SRTM(shuttle radar topography mission)将ALOS/PALSAR校正为像元大小为25m分辨率的图像。为了校正地形对雷达信号强度的影响,采用Zhou[17]提出的地形校正算法计算雷达成像斜距平面与地面的夹角,再进行几何平面的地形辐射纠正得到雷达地形辐射校正后的图像。最后将HJ1B和ALOS/PALSAR进行了图像配准,选择了20个控制点,配准精度在0.2个像元以内,均方根误差RMSE为0.1。

1.3 遥感图像融合分类

为了得到较高精度的分类图像,首先对配准后的HJ1B图像(图 2a)和ALOS/PALSAR SAR图像进行了图像融合处理,然后再进行图像分类。采用离散小波变换将预处理后的2010年6月30日HH极化SAR图像(图 2b)、HJ1B反射率图像、主成分变换图像和归一化植被指数(NDVI)图像进行了图像融合处理。

图 2 HJ1B反射率图像(a);6月30日ALOS/PALSAR HH极化图像(b);分类图(c) Fig. 2 HJ1B reflectance image (a); ALOS/PALSAR HH polarization image acquired on June 30 (b) and classification mapping (c)

决策树分类法已经广泛应用于遥感影像信息提取和土地利用/土地覆被分类研究中[18]。采用决策树分类将融合后的图像进行分类处理,分类结果如图 2c所示。图像共分为10类,分别为建筑、农田、针叶林、灌木林、阔叶林、混交林、沙地、草地、水体和道路。分类的总体精度为85.60%,Kappa系数为0.8269。

1.4 野外数据采集与处理

2012年7月13日—7月30日在研究区开展了野外实地考察和森林参数测量。在研究区内布设样地,选择地形坡度小、植被分布相对比较均匀,布设了50个有代表性的30 m×30 m标准样地(针叶林29个,平均树龄105a,郁闭度40%—80%;阔叶林2个,平均树龄38a,郁闭度30%—60%;混交林19个,平均树龄73a,郁闭度40%—70%),用手持GPS记录每个样地的四角坐标,选择胸径在5 cm以上的林木实行每木检尺,测量胸径和树高,随机选择几株林木测量枝下高,东西和南北两个方向的冠幅。在室内对测量数据进行处理。

1.5 样方生物量计算

森林生物量模型的建立采用Kittredge提出的CAR(Constant Allometric Ratio)算法[19],其数学表达式为:Y=aXb。穆天民[16]、刘斌[20]、梁咏亮[21]等人分别对贺兰山的青海云杉、天然油松林、灰榆疏林等的生物量模型进行了研究,建立了以胸径(DBH)、树高(H)、胸径平方乘树高(D2H)等为自变量的生物量幂指数回归模型,根据这些数学模型计算得到单株树的实测生物量。样地内所有的林木生物量累加得到各样地生物量总和。生物量总和除以样地的面积得到林分地上生物量(t/hm2),计算结果见表 1

表1 样地生物量 Table 1 Statistical field measured AGB of every sample plot
样方号 Quadrat number优势树种 Dominant Species株数 Tree number平均树高/m Average height平均胸径/cm Average DBH林分生物量 Stand AGB/(t/hm2)
0.1油松839.0716.9767.68
02青海云杉/山杨1309.6812.53103.74
03青海云杉/山杨1368.7410.6960.43
04油松1478.0713.7984.91
05油松/山杨1627.519.7945.65
06青海云杉1238.3714.5194.05
07青海云杉/山杨1209.2913.0988.15
08青海云杉/山杨1658.229.6964.66
09油松1047.6715.8654.49
10油松888.8714.1469.61
11油松1488.3515.81102.98
12油松/青海云杉919.0214.6972.66
13油松1107.0912.6964.27
14青海云杉/山杨1648.4910.5189.54
15油松837.2013.7539.13
16油松907.1913.6039.73
17青海云杉916.7412.5646.20
18青海云杉1018.0114.9877.21
19青海云杉9010.3414.7186.63
20青海云杉1578.5512.1095.91
21青海云杉/油松1008.3114.6967.97
22青海云杉/山杨1769.5410.65102.00
23青海云杉/山杨2088.949.5792.30
24山杨1897.9210.1676.02
25青海云杉2057.3811.03104.39
26青海云杉/山杨1747.9513.05119.88
27青海云杉/油松1058.5514.1568.26
28油松1048.9514.2674.98
29青海云杉/油松1208.3512.1870.58
30青海云杉/山杨1989.5012.20153.31
31青海云杉/山杨1237.5411.8454.64
32青海云杉/山杨1448.1310.5256.70
33青海云杉12410.349.1547.43
34青海云杉827.8312.3039.68
35青海云杉/山杨1438.639.8753.51
36青海云杉/山杨1197.7015.2078.12
37青海云杉949.2313.4658.91
38青海云杉8112.3410.3042.25
39青海云杉1019.5012.5057.88
40油松/山杨10411.5012.3067.14
41青海云杉10811.609.2345.62
42青海云杉1169.1214.5080.56
43青海云杉1208.6014.1076.47
44山杨1344.8016.771.07
45青海云杉/山杨9610.3513.565.85
46青海云杉12312.709.1054.44
47青海云杉/山杨8211.1012.1450.55
48青海云杉/山杨977.4015.6064.27
49青海云杉20313.29.4097.09
50青海云杉/山杨13514.3010.2077.50
2 反演方法 2.1 GOMS-MFM反演算法

Li-Strahler几何光学模型(geometric-optical mutual shadowing,GOMS)[22]主要用于遥感影像像元尺度的森林结构参数的反演。多元前向模式(multiple forward mode,MFM)是利用几何光学模型和查找表(LUT)反演生态结构参数的常用方法,它由遥感物理模型正向模拟计算建立参数查找表数据库[23-24]。利用MFM算法进行参数反演,需要从查找表里搜索和图像像元反射率比较接近的值作为模型的最优解。为了得到比较好的反演结果,需要建立比较大的查找表[25]。国外一些学者利用几何光学模型和MFM算法从多光谱图像数据反演了森林结构参数[8-9, 26-27]

GOMS模型的输入参数包括树的结构、地形和光谱信息(表 2)。树的结构参数包括:树密度(λ,trees/m2)、树冠短半轴(r,m)、树冠长半轴(b,m)、树冠中心高度(h,m)、树高分布(dh,m)。根据样地实测数据统计结果输入树的结构参数的最小值、最大值和增加值。地形信息包括坡度(°)和坡向(°),利用分辨率为30 m的ASTER GDEM数据提取研究区的地形坡度和坡向图,统计坡度和坡向的变化范围。模型的光谱信息包括:光照冠层、光照背景和阴影等纯像元组分端元的光谱信息。采用ENVI软件纯像元指数(pixel purity index,PPI)从HJ1B反射率图像提取3种森林类型(阔叶林、针叶林和混交林)的纯像元的组分端元光谱信息。根据GOMS模型的输入格式要求,将3种森林类型输入参数的变化范围和增加值在表 2中分别列出。

表2 GOMS-MFM建立查找表森林结构和地形输入文件 Table 2 GOMS-MFM forest structural and terrain inputs used in creating a set of Look-up table (LUT)
林型Forest type结构参数Structural parameter树密度λ density树冠短半轴rHorizontal crown radius树冠长半轴bVertical crown radius树冠中心高度h Height to crown center树高分布Height distribution坡度Slope坡向Spect查找表大小Size of parameter set
针叶林最小值0.050.50.545001382400
Coniferous forest最大值0.57.510.5202570315
增加值0.0512251045
阔叶林最小值0.050.50.54500860160
Broadleaved forest最大值0.56.510.5182070315
增加值0.0512251045
混交林最小值0.050.50.545001209600
Mixed forest最大值0.56.510.5202570315
增加值0.0512251045

Soenen等人以树冠表面积(SA)和测量的单木的生物量(B,t/株)建立了单木生物量幂指数生长方程来计算林分地上生物量[9]。为了利用遥感图像估算像元尺度的森林地上生物量,根据野外测量的3种森林类型(阔叶林、针叶林和混交林)结构参数数据,以实测树高H为自变量建立3种林型单木生物量和树高H之间的幂指数函数回归模型如(1) 式所示:

(1)

式中,B为单株树地上生物量(kg/株),H为树高(H,m),b0b1为回归方程的系数。由公式(1) 计算得到单株树的地上生物量。利用GOMS模型模拟计算建立参数查找表时,根据树密度参数λ和公式(1) 建立的树高和生物量回归模型可以得到像元尺度的森林地上生物量。

利用HJ1B图像和GOMS模型反演森林参数的过程包括两步:(1) 光谱模拟和林分地上生物量数据库建立。由公式(1) 建立的生物量模型计算得到单木的地上生物量B(B,kg/株),像元尺度的森林地上生物量AGB(t/hm2)为:AGB=λBλ为树密度。GOMS模型正向模拟输出4个波段反射率和生物量信息,然后建立参数查找表文件。(2) 模型解的查找。采用光谱距离的均方根误差(RMSE)优化技术从查找表中搜索图像中每个像元4个波段反射率和模型模拟的对应波段反射率最接近的值,作为模型的最优解。通过查找模型的最优解便得到每个像素的生物量。RMSE的计算公式如式(2) 所示[28-29]

(2)

式中,ρi表示图像的反射率,ρ^i表示模拟的反射率,nb为图像的波段个数。

2.2 线性光谱混合分解

线性光谱混合分解模型(linear spectral mixing model,LSMM)是混合像元分解(spectral mixture analysis,SMA)的常用方法,它假设图像中每个像元的反射率是由该像元中每种地物(纯像元)反射率以及它们的地物类型在该像元内所占像元面积的比例作为权重系数的线性组合[30-31]。它的数学表达式如式(3,4) :

(3)
(4)

式中,Ri,j是第i个纯像元、第j个波段的反射率,Ki是第i个纯像元的面积百分比,vj是第j个波段的误差。m为组分的个数,p为图像波段的个数。

2.3 混合像元分解SMA模型

利用地形校正模型和ASTER GDEM数据对配准后的HJ1B反射率图像进行地形校正,消除了地形变形对遥感图像反射率的影响,将坡地光谱信号校正为水平地面信号[32]。对地形校正后的图像HJ1B进行线性光谱混合分解,分解后得到植被组分图像veg(图 3a)。分别将植被组分图像(veg)和雷达图像ALOS/PALSAR的HH和HV极化的雷达后向散射系数图像进行离散小波融合。利用S-Plus软件的线性回归将样地实测生物量和不同的融合图像(图 3b)对应的图像像元值建立了回归模型,从建立的回归模型中选择精度最好的模型,作为最后的回归模型如式(5) :

图 3 植被组分图像(a)和离散小波融合图像(b) Fig. 3 Vegetation fraction image resulting from spectral mixture analysis (a) and fusion of vegetation fractions with PALSAR HH polarization (b)
(5)

式中,veg_hh表示HJ1B混合像元分解后的植被组分图像veg和2010年6月30日的PALSAR的HH极化融合后的图像。决定系数R2=0.61,均方根误差RMSE=17.11 t/hm2,F值为20.0752,P<0.01。

2.4 PALSAR多元逐步回归模型

采用S-Plus软件的多元逐步回归算法,将两个时相雷达图像HH和HV极化雷达后向散射系数和样地实测生物量进行多元逐步回归分析,建立的地上生物量模型如式(6) :

(6)

式中,σHV0630σHH0815分别表示6月30日获取图像的HV极化雷达后向散射系数和8月15日获取图像的HH极化雷达后向散射系数。决定系数R2=0.63,均方根误差RMSE=12.94 t/hm2,F值为21.62,P<0.01。

3 结果与分析 3.1 3种森林类型地上生物量模型

根据3种不同的森林类型,分别从实测数据中选择200株树,采用公式(1) 计算像元尺度地上生物量,拟合3种森林(阔叶林、针叶林和混交林)地上生物量模型的回归系数、R2、RMSE(kg/株)和F值如表 3所示(P<0.01) 。利用该生物量回归模型模拟计算得到不同森林类型的林分地上生物量(t/hm2)。

表3 三种森林类型生物量模型回归系数 Table 3 Tree level regression coefficients of three forest types used in the estimation tree biomass from field data
森林类型Forest typeb0b1R2RMSE(kg/株)F
阔叶林Broadleaved forest-1.5152.3170.510.39205.948
针叶林Coniferous forest-1.4562.4280.740.57567.706
混交林Mixed forest-2.8652.7890.580.60276.333
3.2 生物量制图

利用本文介绍的3种算法对贺兰山研究区的森林进行生物量制图,研究区的森林地上生物量分布制图分别如图 4所示。根据遥感图像分类图(图 2c),将分类图上非森林覆盖区域(包括建筑物、道路、草地、灌木)合并为一类,定义为非林地。图 4a为利用GOMS模型和MFM方法估算的生物量分布图,其中非林地约占27.57%,乔木林地面积占整个图像区域的72.43%,乔木林最大生物量为332.87 t/hm2。生物量在范围0—50、50—100、100—150、150—200、200—250、>250 t/hm2,分别占总图像覆盖区域面积的28.54%、19.49%、12.07%、6.17%、3.03%、3.12%,其中0—150 t/hm2约占60.11%,>150 t/hm2约占12.32%。图 4b为利用SMA模型估算的生物量分布图,其中非林地约占29.39%,乔木林最大生物量为241.93 t/hm2。生物量在范围0—50、50—100、100—150、150—200、200—250 t/hm2,分别占总图像覆盖区域面积的12.00%、29.79%、24.97%、3.70%、0.14%,其中0—150 t/hm2约占66.76%,>150 t/hm2约占3.84%。图 4c为利用PALSAR回归模型估算的生物量分布图,其中非林地约占34.25%,乔木林最大生物量为177.35t/hm2。生物量在范围0—50、50—100、100—150、150—200 t/hm2,分别占总图像覆盖区域面积的24.35%、32.50%、8.64%、0.25%,其中0—150 t/hm2约占65.49%,>150 t/hm2约占0.25%。

图 4 3种算法估算的森林生物量分布图 Fig. 4 Mapping of forest AGB distribution based on MFM method (a); Mapping of forest AGB distribution based on SMA method (b); Mapping of forest AGB distribution based on PALSAR regression(c)

由3种算法的统计结果可以看到:3种算法估算的生物量主要分布在范围0—150 t/hm2,不同方法估算的生物量结果之间存在一定的差异。采用混合像元分解和PALSAR雷达后向散射系数图像估算的生物量大于150 t/hm2时达到饱和,建立的回归模型无法正确反演,反演的结果存在较大的误差[6]。SMA方法将PALSAR和光学图像进行融合处理减少了生物量估算的误差。

3.3 估算误差分析

为了进一步对3种方法估算的生物量进行误差分析,根据野外实测生物量数据,将MFM、混合像元分解SMA和PALSAR多元逐步回归等3种方法估算的森林地上生物量进行分析。将模型估算地上生物量和实测生物量的差值作为模型估算误差。选择其余的20个实测样地(plots)生物量数据作为真值,计算它们的平均值和标准差(SD)。把实测值落在区间avg-SD和avg+SD内的值记为plots <1 SD,反之则记为plots > 1 SD。根据这个划分标准,这20个实测生物量值分成两组,其中plots <1 SD共有11个,plots > 1 SD共有9个数据。然后分别计算实测生物量、模型估算的生物量、模型估算生物量误差、模型估算生物量plots <1 SD误差(11个数据)和plots > 1 SD(9个数据)误差的平均值、标准差、最小值、最小值和变化范围,统计结果如表 4图 5所示。

表4 估算生物量与实测生物量误差分析 Table 4 Error analysis for estimators and measured AGB
误差分析Error analysis样方数Plots平均值Biomass avg/(t/hm2)标准差SD/(t/hm2)最小值Min/(t/hm2)最大值Max/(t/hm2)范围Range/(t/hm2)
样地实测生物量Field data2085.915.860.4119.959.4
MFM法估算样地生物量Estimated AGB from MFM2082.712.665.5108.342.8
MFM估算的总误差Overall error for all plots208.06.60.125.024.9
plots <1 SD误差Error for plots<1 SD117.15.90.122.722.5
plots >1 SD误差Error for plots >1 SD99.17.12.825.022.2
SMA法估算样地生物量Estimated AGB from SMA2086.214.252.7116.964.2
SMA估算的总误差Overall error for all plots208.25.91.223.322.1
plots <1 SD误差Error for plots<1 SD117.56.11.223.322.1
plots >1 SD误差Error for plots >1 SD99.05.22.020.618.7
PALSAR多元逐步回归估算样地生物量Estimated AGB from PALSAR stepwise regression2072.218.036.3109.773.5
样地总误差Overall error for all plots2015.113.51.239.738.5
plots <1 SD误差Error for plots<1 SD1117.613.91.739.738.1
plots >1 SD误差Error for plots >1 SD912.111.31.236.535.3
SD表示标准差indicated standard deviation;MFM表示多元前向模式算法multiple forward mode;SMA表示混合像元分解算法spectral mixture analysis,respectively

图 5 三种方法生物量估算误差对比 Fig. 5 Comparison of biomass estimation error for MFM,SMA and regression model

根据表 4图 5的误差统计结果分析可以看到:MFM、SMA和PALSAR多元逐步回归等3种方法的平均误差分别为8.0、8.2、15.1 t/hm2,MFM和SMA估算误差比较接近,MFM误差最小;SMA法估算生物量的最小值、最大值和变化范围比实测生物量略宽,而MFM法估算生物量变化范围比实测值变化范围窄,PALSAR回归模型估算生物量变化范围最宽;MFM法和SMA法估算生物量的平均值和实测结果比较接近,而PALSAR回归法估算生物量平均值比实测结果小;MFM法和SMA法估算生物量标准差比实测结果略低,而PALSAR回归模型估算生物量标准差要大于实测结果;11个数据误差(plots <1 SD),MFM法和SMA法的平均误差和标准差都接近总误差(20个数据),而PALSAR回归比总误差(20个数据)大;9个数据误差(plots >1 SD),与<plots 1 SD情况正好相反。

由以上的分析可以看出,采用经验统计和基于遥感物理模型的反演方法从遥感图像数据估算森林地上生物量,不同的反演算法估算的森林地上生物量的精度有一定的差异,原因主要有以下几点:(1) 森林下垫面背景对卫星图像信息有一定的影响。数学模型很难真实的描模拟森林的三维结构及下垫面背景,会造成一定的反演误差。(2) 地形对反演精度的影响。地形改变了太阳、地表和卫星传感器三者所构成的相对观测几何变化,从而造成方向性反射的变化[33],地形对遥感图像的质量和反演精度的影响比较大。(3) 野外实测数据误差的影响。由于地形和林下植被的影响,在野外很难保证实测样地大小为标准样地,同时卫星接收的图像数据受地形坡度和坡向的影响,野外样地的位置很难与图像完全匹配,这对反演的验证精度有一定的影响。

3.4 估算结果验证

为了检验参数反演过程的过拟合现象,有必要结合野外实测数据对模型进行验证。选择20个样地实测生物量和3种方法估算的生物量进行了对比,用实测生物量值验证模型,结果见图 6图 6分别验证了MFM和SMA法估算生物量结果,R2分别为0.61和0.60,实测值与估算值之间的RMSE分别为8.33、9.417 t/hm2,图中各点均分布在1∶1线附近。图 6还验证了PALSAR回归法估算生物量结果,R2=0.39,RMSE=14.89 t/hm2

图 6 生物量实测值与估算值对比值散点图 Fig. 6 Scatter plot on measured and MFM estimated AGB; Scatter plot on measured and SMA estimated AGB; and Scatter plot on measured and PALSAR regression model estimated AGB

3种模型估算的生物量和实测生物量的误差变化范围和模型的验证结果进一步表明:MFM和SMA两种方法估算的生物量结果比较接近实测值,而PALSAR多元逐步回归模型估算的生物量精度最差。由此可见,采用PALSAR多元回归模型估算的山地森林生物量与实测的生物量相差比较大,精度比较低,决定系数小于0.5,该方法估算的生物量精度明显低于其它两种方法估算的生物量。这表明采用几何光学模型和混合像元分解SMA更适合估算山地森林地上生物量。

3.5 讨论

植被的生态结构参数不能通过卫星图像数据直接获取,需要利用经验统计或遥感物理模型反演的方法从卫星图像估算得到。利用Li-Strahler几何光学模型反演森林生物量不需要较多的地面验证数据,模型的输入参数不需要每棵树的详细信息。利用几何光学模型模拟建立查找表时,根据野外直接测量的数据,如树高(H)、胸径(DBH)、胸径平方乘树高(D2H)等数据与生物量建立回归模型,然后通过中间变量(胸径或树高)从遥感图像上估算森林地上生物量。

由于地形引起雷达成像几何关系的变化,进而导致雷达后向散射系数的变化。利用雷达图像PALSAR估算森林结构参数时,地形对估算精度影响较大。因此,采用雷达图像估算森林参数精度较低,反演的结果存在较大的误差。结合雷达图像和光学图像能够减少生物量估算的误差,改善地形对森林参数估算的精度的影响。采用离散小波变换将混合像元分解后的组分图像和雷达图像进行融合,然后建立生物量多元逐步回归模型,该模型比利用雷达后向散射系数建立的生物量模型估算精度更高,这是由于组分图像改善了地形对雷达后向散射系数的影响。因此,利用混合像元分解和离散小波融合可以有效的提高森林生物量的估算精度。

采用光谱距离的均方根误差(RMSE)最小从查找表中搜索模型的最优解,由于测量误差和模型模拟数据的问题,方程的解不唯一或无解,从而导致方程的病态问题[34-35],在实际应用中可以结合一些其它的植被指数如归一化植被指数NDVI或比值植被指数(SR)来查找方程解。采用查找表进行查找时,如果图像的覆盖范围比较大,同时建立的查找表数据库比较大时,在查找过程中查找速度会比较慢,需要对查找表结构和查找算法进行优化。

4 结论

本文探讨了利用GOMS几何光学模型和经验统计模型估算山地森林地上生物量的方法。利用几何光学模型模拟时将可直接测量的参数,如胸径(DBH)或树高(H)作为中间变量,建立树高和森林地上生物量回归模型,计算像元尺度的林分地上生物量,通过MFM方法从卫星图像HJ1B数据反演研究区的森林地上生物量,同时将反演结果和经验统计模型方法进行了对比。研究结果表明,MFM方法比经验统计模型估算森林地上生物量精度更高。采用几何光学模型反演的生物量精度最好,利用混合像元分解和图像融合可以改善地形对森林参数估算的精度的影响,利用PALSAR图像估算的森林生物量的精度最差。因此,几何光学模型和混合像元分解模型适合估算森林地上生物量。

利用几何光学模型进行前向模拟的输入参数不需要较多的先验信息和野外样地调查数据,该方法更适合于通过多时相多传感器的遥感图像反演区域尺度的植被生态结构参数。因此,本文介绍的方法进行森林地上生物量遥感监测研究具有一定的应用潜力。然而,本研究所用的野外观测样方的数量相对偏少,采用统计方法所建立的生物量模型的适用性需要进行一步验证。在今后的研究中,可以考虑通过增加野外观测的数量和小样本回归的方法,比如机器学习算法建立回归模型,并对模型的适用性进行进一步的验证。

致谢: 中国林业科学研究院资源信息研究所谭炳香老师帮助写作,特此致谢。
参考文献
[1] 徐德福, 王让会, 李映雪, 张慧.陆地生态系统碳收支及其主要影响因素分析.中国农业气象,2009,30(4): 519–524.
[2] 刘双娜, 周涛, 舒阳, 戴铭, 魏林艳, 张鑫.基于遥感降尺度估算中国森林生物量的空间分布.生态学报,2012,32(8): 2320–2330.
[3] 何红艳, 郭志华, 肖文发.遥感在森林地上生物量估算中的应用.生态学杂志,2007,26(8): 1317–1322.
[4] Basuki T M, Skidmore A K, Hussin Y A, Duren I V.Estimating tropical forest biomass more accurately by integrating ALOS/PALSAR and Landsat-7 ETM+data.International Journal of Remote Sensing,2013,34(13): 4871–4888.
[5] Rauste Y.Multi-temporal JERS SAR data in boreal forest biomass mapping.Remote Sensing of Environment,2005,97(2): 263–275.
[6] Treuhaft R N, Asner G P, Law B E.Structure-based forest biomass from fusion of radar and hyperspectral observation.Geophysical Research Letters,2003,30(9): 1472–1475.
[7] Elvidge C D, Lyon R J P.Influence of rock-soil spectral variation on the assessment of green biomass.Remote Sensing of Environment,1985,17(3): 265–279.
[8] Scarth P, Phinn S.Determining forest structural attributes using an inverted geometric-optical model in mixed eucalypt forest, Southeast Queensland, Australia.Remote Sensing of Environment,2000,71(2): 141–157.
[9] Soenen S A, Peddle D R, Hall R J, Coburn C A, Hall F G.Estimating aboveground forest biomass from canopy reflectance model inversion in mountainous terrain.Remote Sensing of Environment,2010,114(7): 1325–1337.
[10] 梁存柱, 朱宗元, 王炜, 裴浩, 张韬, 王永利.贺兰山植物群落类型多样性及其空间分异.植被生态学报,2004,28(3): 361–368.
[11] 郑敬刚, 董东平, 赵登海, 何明珠, 李新荣.贺兰山西坡植被群落特征及其与环境因子的关系.生态学报,2008,28(9): 4559–4567.
[12] 俞益民, 赵登海, 梅曙光, 董仁才, 孙萍, 张源润, 周惠玉.贺兰山地区青海云杉生长与环境的关系.西北林学院学报,1999,14(1): 16–21.
[13] 朱源, 康慕谊, 刘全儒, 江源, 和克俭.贺兰山油松和青海云杉交错区的物种丰富度研究.植被资源与环境学报,2007,16(3): 1–6.
[14] 朱源, 康慕谊, 刘全儒, 苏云, 江源.贺兰山针叶林结构特征与种类组成的比较.地理研究,2007,26(2): 305–313.
[15] Jiang Y, Kang M Y, Zhu Y, Xu G C.Plant biodiversity patterns on Helan Mountain, China.Acta Oecologica,2007,32(2): 125–133.
[16] 穆天民.贺兰山区青海云杉森林群落的生物量.内蒙古林业科技,1982,1(1): 34–45.
[17] Zhou Z S, Lehmann E, Wu X, Caccetta P, McNeil S, Mitchell A, Milne A, Tapley I, Lowell K.Terrain slope correction and precise registration of SAR data for forest mapping and monitoring//International Symposium for Remote Sensing of the Environment.Sydney, Australia,2011 1–4.
[18] Mclver D K, Friedl M A.Usin prior probabilities in decision-tree classification of remotely sensed data.Remote Sensing of Environment,2002,81(2/3): 253–261.
[19] Kittredge J.Estimation of amount of foliage of trees and shrubs.Journal of Forest Research,1944,42(1): 905–912.
[20] 刘斌, 刘建军, 任军辉, 杜承星.贺兰山天然油松林单株生物量回归模型的研究.西北林学院学报,2010,25(6): 69–74.
[21] 梁永亮.贺兰山灰榆梳林单株生物量回归模型的研究.林业资源管理,2012(5): 98–104.
[22] Li X W, Strahler A H.Geometric-optical bidirectional reflectance modeling of the discrete crown vegetation canopy:Effect of crown shape and mutual shadowing.IEEE Transactions on Geoscience and Remote Sensing,1992,30(2): 276–292.
[23] Kimes D S, Knyazikhin Y, Privette J L, Abuelgasim A A, Gao F.Inversion methods for physically-based models.Remote Sensing Review,2000,18(2/4): 381–440.
[24] Peddle D R, Franklin S E, Johnson R L, Lavigne M B, Wulder M A.Structural change detection in a disturbed conifer forest using a geometric optical reflectance model in multiple-forward mode.IEEE Transactions on Geoscience and Remote Sensing,2003,41(1): 163–166.
[25] Darvishzadeh R, Skidmore A, Schlerf M, Atzberger C.Inversion of a radiative transfer model for estimating vegetation LAI and chorophyll in a heterogeneous grassland.Remote Sensing of Environment,2008,112(5): 2592–2604.
[26] Peddle D R, Johnson R L, Cihlar J, Latifovic R.Large area forest classification and biophysical parameter estimation using the 5-scale canopy reflectance model in Multiple-Forward-Mode.Remote Sensing of Environment,2004,89(2): 252–263.
[27] Peddle D R, Johnson R L, Cihlar J, Leblanc S G, Chen J M, Hall F G.Physically based inversion modeling for unsupervised cluster labeling, independent forest classification, and LAI estimation using MFM-5-Scale.Canadian Journal of Remote Sensing,2007,33(3): 214–225.
[28] Weiss M, Baret F, Myneni R B, Pragnère A, Knyazikhin Y.Investigation of a model inversion technique to estimation canopy biophysical variables from spectral and directional reflectance data.Agronomie,2000,20(1): 3–22.
[29] Combal B, Baret F, Weiss M, Trubuil A, Macé D, Pragnère A, Myneni R, Knyazikhin Y, Wang L.Retrieval of canopy biophysical variable from bidirectional reflectance:Using prior information to solve the ill-posed inverse problem.Remote Sensing of Environment,2003,84(1): 1–15.
[30] Jacob B G, Mwangangi J M, Mbogo C M, Novak R J.A taxonomy of unmixing algorithms using Li-Strahler geometric-optical model and other spectral endmember extraction techniques for decomposing a QuickBird visible and near infra-red pixel of an Anopheles arabiensis habitat.The Open Remote Sensing Journal,2011,4(1): 1–25.
[31] Tompkins S, Mustard J F, Pieters C M, Forsyth D W.Optimization of endmembers for spectral mixture analysis.Remote Sensing of Environment,1997,59(3): 472–489.
[32] Riaño D, Chuvieco E, Salas J, Aguado I.Assessment of different topographic corrections in landsat-TM data for mapping vegetation types.IEEE Transactions on Geoscience and Remote Sensing,2003,41(5): 1056–1061.
[33] 王新云, 郭艺歌, 过志峰, 覃文汉, 孙国清.地形因子对森林冠层BRF模拟的影响分析.遥感信息,2012,27(4): 82–85.
[34] Bacour C, Jacquemoud S, Leroy M, Hautecoeur O, Weiss M, Prévot L, Bruguier N, Chauki H.Reliability of the estimation of vegetation characteristics by inversion of three canopy reflectance models on airbone PolDER data.Agronomy,2002,22(6): 555–565.
[35] Meroni M, Colombo R, Panigada C.Inversion of a radiative transfer model with hyperspectral observation for LAI mapping in polar plantation.Remote Sensing of Environment,2004,92(2): 195–206.