林业科学  2018, Vol. 54 Issue (2): 68-80   PDF    
DOI: 10.11707/j.1001-7488.20180208
0

文章信息

李晓红, 陈尔学, 李增元, 李世明
Li Xiaohong, Chen Erxue, Li Zengyuan, Li Shiming
综合应用多源遥感数据的面向对象土地覆盖分类方法
Object Based Land Cover Classification Method Integrating Multi-Source Remote Sensing Data
林业科学, 2018, 54(2): 68-80.
Scientia Silvae Sinicae, 2018, 54(2): 68-80.
DOI: 10.11707/j.1001-7488.20180208

文章历史

收稿日期:2016-02-22
修回日期:2016-06-28

作者相关文章

李晓红
陈尔学
李增元
李世明

综合应用多源遥感数据的面向对象土地覆盖分类方法
李晓红, 陈尔学, 李增元, 李世明    
中国林业科学研究院资源信息研究所 北京 100091
摘要:【目的】针对国家森林资源宏观监测业务对区域森林资源空间分布信息的迫切需求,发展一种基于国家森林资源连续清查固定样地数据,可充分发挥GF-1宽幅多光谱数据、MODIS遥感数据相应空间和时间分辨率优势的面向对象土地覆盖分类方法,以提高林地和森林资源的监测精度和自动化程度。【方法】以黑龙江省小兴安岭某林区为研究区,主要数据源包括GF-1宽幅多光谱数据、MODIS NDVI(250 m,8天合成)时间序列遥感数据、国家森林资源连续清查固定样地数据以及少量外业实地调查数据等。首先,基于GF-1宽幅多光谱数据进行多尺度影像分割,将研究区划分为许多区域性的分割对象;然后,以分割对象为分析单元,分别提取GF-1宽幅多光谱遥感影像的光谱特征、纹理特征、形状特征等以及MODIS NDVI时间序列遥感数据的时序特征,并采用随机森林算法进行特征选择;最后,利用训练样本建立基于分类回归树分类器完成面向对象的土地覆盖分类方法研究,分别比较单一GF-1 16 m宽幅多光谱数据、单一MODIS NDVI时间序列遥感数据以及综合多源数据的分类结果,并基于混淆矩阵对分类结果进行分析。【结果】精度检验和分析结果表明,面向对象的综合多源遥感数据分类方法总体分类精度达89.46%,Kappa系数为0.874,明显优于仅基于GF-1宽幅多光谱数据、MODIS NDVI时间序列遥感数据的分类方法。【结论】综合应用多源遥感数据的面向对象土地覆盖分类方法适用于综合GF-1与GF-4数据的土地覆盖类型分别制图,可有效提高主要土地覆盖类型的分类精度。针对国家森林资源连续清查的业务需求和特点,本文所发展的方法在分类对象生成、特征提取、特征选择、分类器训练和精度检验等关键环节均进行了优化设计,有利于提高森林资源连续清查业务中主要林地类型遥感分类制图的自动化、标准化程度。
关键词:多源数据    GF-1宽幅多光谱数据    MODIS NDVI遥感数据    随机森林    面向对象    土地覆盖分类    
Object Based Land Cover Classification Method Integrating Multi-Source Remote Sensing Data
Li Xiaohong, Chen Erxue , Li Zengyuan, Li Shiming    
Research Institute of Forest Resource Information Techniques, CAF Beijing 100091
Abstract: 【Objective】In order to meet the urgent needs of national forest inventory (NFI) for monitoring national forest resources on a macro-scale, an object-oriented regional land cover classification method was developed in this paper by using permanent forest plot data of NFI and integrating the corresponding spatial and temporal resolution advantages of GF-1 WFV multi-spectral data and MODIS remote sensing data.【Method】The test site is located in the central of the Xiaoxing'an mountain region in Heilongjiang Province. The GF-1 WFV multispectral data, time series MODIS NDVI product of 8 days synthetic (250 m spatial resolution), the permanent forest plot data collected by the NFI and some field survey data are employed as the key data sources. After image segmentation, spectral, texture and shape features from GF-1 WFV multi-spectral data and NDVI features from times series MODIS NDVI data are extracted for each object. Based on these features, the random forests algorithm is adopt to select the best features automatically and then the classification and regression tree is used to administer the supervised classification. The permanent forest plot data are used to build and validate the decision tree classifier. This method has been investigated with the data collected for the study site.【Result】Experimental results show that the overall accuracy and Kappa coefficient of the developed method combing multi-sources data can reach 89.46% and 0.874 respectively, with significant improvement compared with that using either GF-1 WFV multi-spectral data or MODIS NDVI time series data alone.【Conclusion】The land cover classification method developed in this study is appropriate for integrating GF-1 data and GF-4 data for land cover classification mapping and can improve accuracy of land cover classification effectively.Compared with existed research work, this paper has developed one land cover classification method of much more practical application value with some special features such as focusing on the operational application needs and the characteristics of NFI, optimizing the key procedures such as classification objects generation, feature extraction, feature selection, classifier training, accuracy validation, and so on.The land cover classification method developed is useful for improving the automation and standardization of the technique flow to produce forest land cover maps implemented by NFI.
Key words: multi-sources data    GF-1 WFV multi-spectral data    MODIS NDVI remote sensing data    random forest    objected-oriented    land cover classification    

国家森林资源连续清查简称“一类调查”, 以固定样地连续调查监测进而估计总体为主要手段, 每5年进行1次。在5年调查期内, 每省确定1个固定年份开展固定样地调查, 而且必须在该年内完成。固定样地调查的第1个观测要素就是地类(地表覆盖/利用类型), 须按照《国家森林资源连续清查技术规定》确定的分类系统进行。根据所有样地的地类抽样调查结果, 采用成数估计法推算出调查总体(省)的各地类面积, 重点是统计林业用地下的2级地类(如乔木林地、灌木林地、疏林地等)及3级地类(如针叶林、阔叶林等)。遥感影像在一类调查中主要用于设立遥感判读样地, 为分层抽样设计提供辅助数据。虽然固定样地数据可为遥感影像分类提供地面实况数据, 有助于提高林业用地、有林地类型等遥感分类制图的自动化、标准化程度, 但目前将森林资源连续清查固定样地数据用于遥感影像分类的相关研究还比较缺乏。

森林扰动、林地更新或造林等是引起林地类型发生变化的主要因素, 所导致的变化通常发生在约0.5 hm2尺度, 因此适合省级全覆盖地类遥感分类制图的影像分辨率应具有中等空间分辨率(约30 m), 这也是目前国内外主要采用Landsat TM/ETM+、CBERS等多光谱影像辅助国家森林资源连续清查的原因。但由于多光谱遥感影像光谱、时间分辨率不高, 通常采用计算机分类方法进行分类的精度和自动化程度都比较低(章恒等, 2015)。高时间分辨率遥感数据能够提供反映各种地物随物候特征变化的信息, 有利于地表覆盖类型的分类(Wardlow et al., 2008), 自高时间分辨率遥感影像产生的NDVI时序数据已经成为开展大区域植被和土地覆盖分类的重要手段(Li, 1997)。李俊祥等(2005)利用时间跨度8个月的时间序列NOAA-AVHRR NDVI(空间分辨率1.1 km)最大值合成影像, 采用ISODATA迭代聚类算法对中国东部地区的五省一市植被进行分类, 结果分出了28种土地覆盖类型; 孙小添等(2013)基于MODIS 8天合成反射率数据(空间分辨率250 m)并结合陆地表面温度数据, 利用决策树方法对森林类型进行分类, 总体分类精度达87.80%。虽然利用高时间分辨率遥感数据可提高分类精度, 但由于空间分辨率较粗, 不太适合国家森林资源连续清查业务应用。

为了能够综合利用中等空间分辨率多光谱数据和粗空间-高时间分辨率遥感数据各自在空间和时间分辨率上的优势, 有学者提出将Landsat多光谱数据与MODIS反射率数据进行数据融合的算法(Gao et al., 2006), 进而提取与多光谱数据具有相同空间分辨率的NDVI时间序列数据用于地表覆盖分类方法研究(Thomas et al., 2009; 李玉东等, 2013), 但这类数据融合算法并没有考虑地形对反射率影响的校正问题。贾明明等(2014a)采用面向对象分析方法, 将我国环境星多光谱数据与MODIS时间序列数据结合, 对吉林省东部地区进行森林类型提取, 在提高运行效率的同时削弱了地形的影响。

我国2013年成功发射的高分一号(GF-1)遥感卫星, 其宽幅多光谱数据(空间分辨率16 m)幅宽可达200 km, 具有每年至少1次全国全覆盖观测的能力; 同时, 我国2015年底成功发射的高分四号(GF-4)凝视遥感卫星, 具有获取高时间分辨率NDVI时间序列数据(空间分辨率50 m)的能力。因此, 发展一种以国家森林资源连续清查固定样地数据为主要地面实况数据来源, 并可综合利用GF-1和GF-4各自空间和时间分辨率优势的分类方法, 对提高我国林地及森林资源的监测精度和自动化程度具有重要意义。在本文研究阶段, GF-4遥感卫星还未发射, 因此以MODIS作为替代数据源开展研究。本文通过对国内外的最新研究进展进行对比, 并考虑所采用遥感数据的特点, 选择了面向对象的总体分析方法(贾明明等, 2014b), 与已有研究相比, 本文特色在于针对国家森林资源连续清查的业务需求和特点, 对分类对象生成、特征提取、特征选择、分类器训练和精度检验等关键环节均进行了优化设计, 旨在发展一种具有较高实际应用价值的地表覆盖类型分类方法。

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

研究区为一景GF-1 16 m宽幅多光谱遥感影像数据所覆盖区域, 位于黑龙江省中部小兴安岭地区, 127°14′11″—130°24′59″E, 46°02′28″—48°18′54″N。研究区属北温带大陆性季风气候, 四季分明, 冬季严寒干燥, 夏季温热湿润; 年均温0 ℃左右, 1月均温-25 ℃, 7月均温22 ℃; 属于低山丘陵地形, 山势较为缓和; 土地覆盖/利用类型以林业用地为主, 非林业用地覆盖较少, 其中耕地/农田主要分布于洼地底部、山坡下部及西部平原地区, 主要作物为大豆(Glycine max)、水稻(Oryza sativa)和玉米(Zea mays); 林业用地中的有林地主要优势树种类型包括白桦(Betula platyphylla)、椴树(Tilia tuan)、落叶松(Larix gmelinii)、红松(Pinus koraiensis)和云杉(Picea asperata), 前3类为落叶树种, 后2类为常绿树种。

1.2 数据

研究所采用的多源遥感数据主要包括单个时相的GF-1宽幅多光谱数据和MODIS NDVI时间序列数据, 地面调查数据主要包括黑龙江省第八次国家森林资源连续清查固定样地数据及笔者于2015年9月获取的地类实况调查数据。

1.2.1 GF-1宽幅多光谱数据

一景成像时间为2013年9月30日的GF-1宽幅多光谱影像, 具有近红外(NIR)、红光(R)、绿光(G)和蓝光(B)4个波段, 空间分辨率16 m。首先对该多光谱影像进行辐射定标, 然后应用ENVI软件的大气校正模块(FLAASH)进行大气校正, 最后对其进行正射校正处理。正射校正处理所需控制点以一景已经做过正射校正处理的Landsat ETM+全色影像(空间分辨率15 m)为参考, 通过尺度不变特征变换(SIFT)算法自动寻找不变特征匹配点104个, 并在自动寻找的控制点分布较少区域人工选择控制点18个, 控制点的高程信息自研究区的ASTER DEM(已由其原始的30 m格网重采样为16 m)提取。以获取的控制点三维坐标信息为输入, 对GF-1宽幅多光谱影像自带的RPC(rational polynomial coefficients)几何校正模型进行参数优化解算, 优化后的正射模型东西向误差为0.52 m, 南北向误差为0.47 m。利用该模型和ASTER DEM对GF-1宽幅多光谱影像进行最近邻法重采样处理, 得到最终的多光谱正射校正影像。

1.2.2 MODIS NDVI时间序列数据

从NASA网站(http://ladsweb.nascom.nasa.gov)下载2013年8天合成的MODIS数据集产品MOD09Q1, 共46个时相, 空间分辨率250 m, 每个时相数据包括近红外和红光2个波段。对下载的数据进行如下处理:1)分别对每个时相数据, 利用MRT(MODIS Reprojection Toolbox)软件进行地图投影转换, 由原来的正弦曲线投影转换为通用横轴墨卡托投影(UTM WGS84); 2)分别对每个时相数据, 由近红外与红光波段反射率计算得到NDVI; 3)将处理好的NDVI数据按照时间顺序进行波段叠加, 生成一个包含46个波段的NDVI时间序列数据; 4)采用S-G(Savtzky-Golay)滤波算法(郭芬芬等, 2011)对NDVI时间序列数据进行滤波; 5)将46个波段的NDVI数据进行立方卷积重采样处理, 得到16 m×16 m像元大小的NDVI数据, 使NDVI时间序列数据与GF-1宽幅多光谱数据具有相同的像元大小, 以方便后续处理和分析。

图 1所示为一个落叶林NDVI时间序列曲线滤波前后的对比。滤波前曲线变化比较剧烈, 特别是在森林生长旺季, 这显然不符合平稳生长期森林NDVI随时间变化比较平缓的事实, 因此特别突然的NDVI变化可认为是MODIS 8天合成遥感数据本身缺陷引起的“噪声”。S-G滤波算法可降低“噪声”的影响, 使得NDVI随时间的变化特征体现了地物本身生物物理特性随物候的变化规律。

图 1 S-G滤波前后NDVI随时间变化曲线比较 Figure 1 Comparison of the NDVI changing curves with time before and after S-G filtering
1.2.3 国家森林资源连续清查固定样地数据

采用覆盖研究区域的黑龙江省第八次国家森林资源连续清查固定样地数据, 调查时间为2010年, 虽然与本研究所用遥感影像(2013年获取)相差2年, 但由于固定样地地类发生变化的可能性较小, 因此认为2010年调查的固定样地地类可代表 2013年的实际地类状况。图 2展示了本研究所使用样地在GF-1宽幅多光谱影像上的空间分布情况, 从北到南靠近中间的区域, 森林分布较为集中, 但样地分布却比较稀疏, 样地之间的间隔为8 km×8 km; 左下角、右下角和右上角部分区域, 林业用地与非林业用地镶嵌分布, 有林地、灌木林地等分布较为破碎, 但样地分布却比较紧密, 样地之间的间隔为8 km×4 km。固定样地总数为791个。

图 2 研究区国家森林资源连续清查固定样地数据空间分布 Figure 2 Spatial distribution of the permanent plot data collected by the national forest inventory
1.2.4 外业调查数据

为了更好地辅助影像分类、分类结果的分析及精度检验, 笔者于2015年9月10日—10月6日在研究区开展了外业调查工作。采用手持GPS记录采样点的地理坐标, 记录采样点及其周边地物覆盖实况并拍摄照片, 最终得到主要地表覆盖类型的野外调查样本50个。

2 研究方法 2.1 分类系统及分类类别

《国家森林资源连续清查技术规定》确定的分类系统包含2个一级地类:林地和非林地。林地包含8个二级地类:乔木林地、灌木林地、竹林地、疏林地、未成林造林地地、苗圃地、迹地和宜林地。非林地包括5个二级地类:耕地、牧草地、水域、未利用地和建设用地。该分类系统是为摸清国家林地及森林资源的宏观现状而确定的, 调查人员在调查现场容易确定样地所属的土地利用/覆盖类型, 但利用本文拟研究的遥感影像计算机分类技术很难对“土地利用类型”进行有效分类, 更适用于“土地覆盖类型”的分类。同时, 本文要分的类别应该具有较高的可分性, 即便是综合利用多源遥感数据, 也只能对有限的几个类别进行有效分类。基于以上考虑, 确定了本研究的目标分类类别, 包括林地和非林地2个一级地类, 其中林地包括乔木林地和灌木林地2个二级地类; 本试验区虽然只出现乔木林, 但为了突出对森林类型的分类, 将乔木林划分为常绿林和落叶林, 作为2个自定义的三级地类。非林地下确定了耕地、建设用地和水域3个二级地类, 不能分为以上6个地类的地物为其他地类。最终确定本研究的目标分类类别共有7个:常绿林、落叶林、灌木林地、耕地、建设用地、水域和其他。

2.2 分类总体技术方法

本研究确定的分类技术流程如图 3所示, 总体上是一个面向对象的监督分类方法, 主要包括基于GF-1宽幅多光谱影像进行影像分割、最佳分割尺度的确定、以对象为分类单元提取分类特征、采用随机森林算法进行特征选择、利用训练样本建立决策树分类器并对所有对象进行分类、采用精度评价样本对分类精度进行评价等步骤。

图 3 土地覆盖分类技术流程 Figure 3 Flow chart of the land cover classification technique
2.2.1 影像分割

影像分割在 eCognition 软件中实现。将含有4个波段的GF-1 16 m宽幅多光谱遥感影像、含有46个波段的MODIS NDVI时间序列遥感数据依次载入到eCognition软件中, 生成的待分类影像共50个波段, 其中波段1-4为GF-1 16 m宽幅多光谱遥感数据的第1波段(Band 1)至第4波段(Band 4), 波段5-50为MODIS NDVI时间序列遥感数据的第1时相(第1天)至第46时相(第361天)。

在多尺度分割中, 分割尺度的选取至关重要, 不同分割尺度将得到不同的分割结果。一般而言, 分割尺度越大, 所生成的对象层内多边形面积就越大且数目越小, 反之亦然(张俊等, 2011)。最优的分割尺度应该是所分割的结果能通过1个或者多个影像对象来表达该地物, 同时又不会出现过于破碎或边界模糊的现象。

影像分割的目的是为了得到合理的分析单元, 使其既能充分利用GF-1宽幅多光谱数据的高空间分辨率优势, 又能充分利用MODIS NDVI时间序列遥感数据的高时间分辨率优势。本研究中由于GF-1宽幅多光谱数据空间分辨率为16 m, 最小分辨单元面积为256 m2, 而MOD09Q1时间序列数据空间分辨率为250 m, 最小分辨单元面积为62 500 m2, 是GF-1宽幅多光谱数据的244倍, 差别很大, 为方便后续按对象提取MODIS NDVI特征, 应使分割得到的绝大部分对象都至少包含1个MODIS像元。因此, 在对GF-1宽幅多光谱遥感数居进行分割时, 本文选择了80、100和120共3种较大的分割尺度进行对比分析。最优尺度的选取和尺度参数定量选取实际上是一个不断对比试验的过程, 即通过选取不同的尺度及参数设置, 目测选择最为合适的分割尺度或与已知数据进行叠加, 分析分割结果对象边界与原始数据的吻合度来实现最优尺度的选取。

以GF-1宽幅多光谱遥感影像4个波段数据为输入, 进行多尺度分割, 4个波段权重都设为1。分割时, 参考经验值, 发现对于研究区及所采用的数据, 将形状因子、光滑度和紧致度分别设置为0.2、0.5和0.5, 得到的多边形能够较好显示地物的边界。对比试验发现, 当将形状因子、光滑度和紧致度作为固定的常量时, 分别设置为3种(80、100和120)不同分割尺度参数变量, 得到的分割效果对比如图 4所示。当分割尺度为100时, 对象边界与原始数据的吻合度最高, 产生的对象平均大小尽量小的同时又不过于破碎, 既可以满足多源数据综合的需求, 又能够考虑到自身地物特征。

图 4 不同分割尺度分割效果对比 Figure 4 Comparison of the segmentation results with different segmentation scale parameters (Left:80, Center:100, Right:120) (左:80, 中:100, 右:120)

研究中, 分割尺度参数的确定, 不仅基于传统的目视判读方式, 而且还采用了面积统计方法, 如图 5所示即为从尺度100的分割结果中随机抽选500个对象计算得到对象面积的统计直方图。图 5中竖线与横轴的交点对应的对象面积等于1个MODIS像元面积(250 m×250 m), 竖线右边表示面积大于1个MODIS像元的对象个数, 有96%的对象面积大于1个MODIS像元面积, 剩余的是一些比较细小、破碎的对象, 主要包括建筑用地和一些小的水域。

图 5 分割对象面积统计直方图(尺度:100) Figure 5 Statistical histogram of the segmentation objects area (Scale:100)
2.2.2 特征提取

以对象为单元分别从GF-1宽幅多光谱影像和MODIS NDVI时间序列数据中提取分类特征, 包括自GF-1宽幅多光谱影像提取的植被指数、光谱特征、统计特征、纹理特征、形状特征和位置特征共28个以及自MODIS NDVI时间序列数据提取的46个NDVI特征(每个时相1个)(表 1), 特征总数为74。后面将利用训练样本采用随机森林算法在74个特征中选择对分类最为有用的特征。

表 1 按对象提取的分类特征 Tab.1 Classification features extracted for each object
2.2.3 训练和精度检验样本的获取

分类器训练和精度检验样本均基于国家森林资源连续清查固定样地数据确定。将高空间分辨率影像分割后的矢量多边形图层与森林资源清查样地数据矢量图层(点)求交集, 得到对应森林资源清查样地点的多边形矢量图层。图 6左图为所有森林资源清查样地对应的多边形矢量图层在遥感影像上的叠加显示, 将这些多边形按照2:1比例分为训练样本和精度检验样本, 如图 6右图所示, 其中红色为训练样本多边形(527个)、黄色为检验样本多边形(264个)。图 7所示为局部放大的多边形样地示意图。图 8显示了各个地类的样本数(不分训练样本和精度检验样本), 显然各地类的样本数存在严重的不均衡现象, 如灌木林地及水体类别过少, 阔叶林类别过多。待分类样本数不均衡将会导致分类规则过度拟合或欠拟合问题, 降低分类精度。为解决此问题, 本文参考外业调查获取的地面实况数据, 对主要类别的样本数进行了调整, 最终得到样本共计838个, 按照2:1比例, 随机分为训练样本558个、精度检验样本280个, 样本数调整前后的对比如图 8所示。

图 6 按对象获取的分类训练与精度检验样本 Figure 6 Classification training and accuracy test samples extracted at classification object level (Red:training samples, Yellow:accuracy test samples) (红:训练样本, 黄:检验样本)
图 7 多边形样本数据 Figure 7 The polygon samples data
图 8 调整前后的各类别样本数 Figure 8 Sample number of all the classes before and after adjustment
2.2.4 特征选择

随机森林(random forest, RF)算法是Leo Breiman于2001年提出的一个组合分类算法, 其通过自助法(bootstrap)重采样技术, 从原始训练样本集Z中有放回地、重复地随机抽取K个样本生成新的训练样本集, 然后根据自助样本集生成K个分类树组成随机森林, 待测数据的分类结果按分类树投票所得分数确定(Breiman, 2001)。其实质是对决策树算法的一种改进, 与传统决策树算法相比, 有更强的泛化能力和更好的分类效果。由于RF算法可以对每个变量的重要性进行分析, 所以本研究尝试采用该算法进行特征选择。利用RF算法进行特征选择的基本原理是将待考察属性随机替换, 然后评价替换前后预测器的预测正确率变化, 从而对变量进行评分。依据评分, 对变量进行排序从而判断在预测中起到的作用大小, 并筛选出相对重要的变量。

特征选择结果如表 2所示。其中, IncMSE代表自变量出现在袋外数据时的模型均方误差(MSE)增量, IncNode Purity代表自变量出现在袋外数据时对模型树节点纯度的影响力, 二者值越大, 说明特征变量的贡献度越大。本研究中特征选择的前提条件是IncMSE大于15, 对于出现的相邻时相NDVI特征变量, 参考图 9所示的MODIS NDVI时间序列变化曲线, 若二者对各典型地类响应大致相同, 只选择其中之一。随机森林特征选择个数统计如表 3所示, 主要包括GF-1 NDVI、GF-1第1/3/4波段标准差、亮度值、最大差、GF-1 4个波段反射率均值及MODIS第1/14/29/35时相的NDVI。

表 2 随机森林特征选择结果 Tab.2 Feature extraction results based on random forest
图 9 各典型地类的MODIS NDVI时间序列变化 Figure 9 The characteristic of MODIS NDVI time series curves for the key land cover types
表 3 随机森林特征选择个数统计 Tab.3 The statistics of feature selected based on random forest

特征选择结果中, 贡献率最大的为NDVI特征值, NDVI特征值既包括GF-1 16 m宽幅多光谱遥感数据的NDVI特征, 又包括MODIS时间序列数据的NDVI特征, 其中, 筛选出MODIS时间序列的NDVI特征主要包含时间序列中的第1时相、第14时相、第29时相和第35时相。结合图 9所示的各典型地类NDVI时间序列变化曲线可以发现, 其第1/14/29/35时相可分别代表时间序列曲线的基准值、开始生长时间、峰值(生长期的中期时刻)和停止生长时间4个体现不同地物物候变化特性的参数。

GF-1 16 m宽幅多光谱遥感影像各个波段的标准差、亮度值、平均差分值和亮度值等特征, 均出现在特征选择结果中, 说明其在区分各典型地类的分类过程中同样发挥着重要作用。

通常情况下, 纹理特征的加入, 可为高空间分辨率遥感影像的目视解译带来额外显著增益; 但是本研究发现, 针对GF-1 16 m宽幅多光谱遥感影像的纹理特征, 随机森林特征选择后的结果贡献较小。其原因一方面可能与遥感影像数据本身有密切关系, GF-1宽幅多光谱遥感数据本身的空间分辨率为16 m, 属于较高空间分辨率遥感影像, 纹理特征相对高空间分辨率遥感影像不明显; 另一方面可能是研究区位于黑龙江省小兴安岭地区, 林区覆盖面积70%以上, 影像分割后, 各典型地类以对象为单元提取的纹理信息差异不显著, 无法有效识别各典型地物类别。对于GLCM、GLDV这2种纹理特征提取方法本身而言, 可能存在的问题包括: 2种方法本身与人类视觉模型脱节, 缺少全局信息的利用, 难以研究纹理尺度间像素的遗传或依赖关系且计算复杂度很高, 制约了其实际应用(刘丽等, 2009)。与纹理特征类似, GF-1宽幅多光谱数据的对象形状特征和位置特征也仅提供了有限的信息。

2.2.5 决策树训练及分类

决策树是一种直观的知识表示方法, 其以信息论为基础, 将复杂的决策形成过程抽象成易于理解和表达的规则或判断(潘琛等, 2008)。CART(classification and regression tree, 分类回归树)是目前被广泛采用的决策树分类算法, 其基本原理是通过对由测试变量和目标变量构成的训练数据集的循环分析, 形成二叉树形式的决策树结构(张晓娟等, 2010)。本研究以CART决策树为分类器, 利用分类训练样本(对象)地类属性和所选择的分类特征训练分类器, 建立分类决策树, 进而用于每个待分类对象的分类。

2.2.6 分类精度检验

以精度检验样本地类为对象地类真值, 对精度检验样本对象的分类正确性进行检验, 计算分类混淆矩阵, 得到各地类分类精度、总分类精度及Kappa系数。

3 结果与分析 3.1 各地类的NDVI时间变化规律分析

为了说明本研究所构建的MODIS NDVI时间序列遥感数据对典型地类的可分性, 以图 9为例分析典型地类对NDVI时间序列的响应关系。可以发现, 常绿林全年NDVI均值在0.70左右, 曲线较为平缓; 落叶林从第14时相(4月初)开始发芽生长, NDVI急剧增加, 第25/26时相(7月)NDVI达到顶峰, 枝叶繁茂, 第35时相(9月)缓慢下降, 开始落叶, 第39时相(10月)曲线平缓; 建筑用地全年NDVI均值在0.13左右, 值为正数, 曲线较为平缓; 水域全年NDVI均值在0上下浮动, 值有正有负; 耕地生长周期从第20时相(5月)开始, 峰值从第25时相(7月)持续到第30时相(8月), 此后急剧下降; 灌木林地与落叶林类似, 但是在非生长周期(第1-13、40-46时相)NDVI远小于落叶林, 均值在0.10左右。

3.2 分类结果及分析

采用决策树分类算法, 分别对单一GF-1宽幅多光谱数据、单一MODIS NDVI时间序列数据及多源遥感数据进行分类, 结果如图 10所示。

图 10 分类结果对比 Figure 10 Comparison of classification results

比较图 10中3种不同分类结果发现, 单一采用GF-1多光谱数据的分类结果(图 10a), 研究区左下角建筑用地与耕地混分严重, 左侧水域与灌木林地出现错分现象, 无法准确区分有林地; 单一使用MODIS NDVI时间序列数据的分类结果(图 10b), 由于可反映不同地类随物候变化特点的时间序列NDVI特征的加入, 地类错分、误分现象得到改善, 并能够区分出部分落叶林与常绿林, 然而, 由于MODIS数据空间分辨率较低, 对于面积较小或相对破碎的地类不能精确识别, 如在研究区中部的部分常绿林和建筑用地出现了错分现象; 图 10c为综合利用2种数据源的分类结果, 由于综合利用了2种数据源在空间和时间分辨率上的优势, 总体分类结果得到了很大改善。

图 11所示为分类结果局部放大示意图。第1排遥感影像对应区域的主要地类为耕地和建筑用地。图 11a为GF-1多光谱原始影像, 黑色圆圈内呈亮白色的地类为建筑用地, 而黑色方框内地类为耕地, 该区域主要包括2种耕地:呈红色的地块为未收割农田, 呈靛蓝色的地块为已收割农田; 图 11b为只利用GF-1多光谱影像的分类结果, 耕地与建筑用地错分严重, 尤其是已收割农田, 有很多被错分成建筑用地; 图 11c为只利用MODIS NDVI时间序列数据的分类结果, 大部分耕地可以被正确分类, 但部分未收割农田错分为有林地, 且对建筑用地的分类能力较差; 图 11d为综合利用2种遥感数据源的分类结果, 有效避免了单一数据源的分类错误, 可准确区分耕地与建筑用地。

图 11 局部区域分类结果对比 Figure 11 Comparison of classification result of some regions

图 11第2排遥感影像对应区域的主要地类为有林地(包括常绿林、落叶林2类)。图 11e为GF-1多光谱原始影像, 地形起伏变化明显, 植被生长受地形影响较大, 主山脉为东北—西南走向, 次山脉多为东西走向; 图 11f为只利用GF-1多光谱影像的分类结果, 未能有效区分常绿林与落叶林; 图 11g为只利用MODIS NDVI时间序列数据的分类结果, 对于面积较大的常绿林能够正确分类, 但很多面积较小或比较破碎的落叶林被误分为常绿林, 且林地细节信息保留较少; 图 11h为综合利用2种遥感数据源的分类结果, 改善了单一采用MODIS NDVI时间序列数据不能有效区分的面积较小或比较破碎常绿林的分类问题。

多源数据综合分类结果混淆矩阵如表 4所示, 不同分类结果的精度对比如表 5所示。综合应用多源遥感数据进行分类, 除灌木林地外, 各地类分类精度均有显著提高, 总分类精度达89.46%, 其中常绿林精度为91.34%, 落叶林精度为87.11%, 灌木林地精度为81.25%, 建筑用地精度为93.08%, 耕地精度为93.12%, Kappa系数为0.87, 总分类精度及Kappa系数远高于单一GF-1 16 m宽幅多光谱数据分类结果及单一MODIS NDVI时间序列数据分类结果。可见, 研究中综合应用多源数据的分类方法优于传统的基于单一遥感数据源方法, NDVI时间序列数据的加入能够有效提高分类精度。

表 4 多源数据分类结果混淆矩阵 Tab.4 Confusion matrix of multi-sources data classification result
表 5 不同分类结果的精度对比 Tab.5 Accuracy validation results of different using different remote sensing data

表 4中, 制图精度最高的为耕地(93.08%)。综合应用多源数据分类后, 总分类精度较单一GF-1多光谱数据提升了15%左右, Kappa系数提高了0.18左右, 总分类精度较单一MODIS NDVI时间序列数据提升了2.8%左右, Kappa系数提高了0.05左右。对于单一GF-1多光谱数据分类结果, 常绿林、水域及其他地类的分类结果较差, 也正是由于这3类地类的分类精度远低于综合多源数据分类方法的分类精度, 因此拉低了单一GF-1多光谱数据分类的总精度。

表 4中第3、4列方框内类别数量统计可以看出, 落叶林和灌木林地的混分程度高于其他地类。出现这种现象的原因:一方面可能是训练样本中灌木林地样地仍然偏少导致, 因为理论上, 当用于分类的特征不变时, 样本数增加有利于提高分类精度或分类模型的推广能力, 样本数相对少的类别, 由于代表性不够, 可能导致其精度降低; 另一方面, 研究区内灌木林地主要包括有兴安杜鹃(Rhododendron dauricum)、毛赤杨(Alnus sibirica)、筐柳(Salix linearistipularis)及少量的桦树等植被类型, 均属于落叶性植被, 光谱曲线和NDVI时间序列曲线与落叶林差异较小, 不易区分。同样可以用于解释表 5多源数据综合分类结果中灌木林地的分类精度与仅采用GF-1宽幅多光谱影像分类相比有小幅下降的原因。

表 4中第5、6列方框内类别数量统计还可以看出, 耕地和建筑用地同样存在一定程度的混分现象。主要原因包括:一是耕地和建筑用地的位置关系, 通常情况下, 二者紧密相连, 且均分布在较为平坦的区域, 即以建筑用地为中心, 周边一般存在大量耕地, 由此导致影像分割结果边界没有完全按照实际情况准确区分二者; 二是由于研究区所在区域以及GF-1 16 m宽幅多光谱数据成像时间、建筑用地和耕地(尤其是已收割耕地)光谱差异较小, 导致混分现象出现。

以上定性和定量分析结果表明, 本文所发展的综合中等空间分辨率GF-1宽幅多光谱影像和粗空间分辨率、高时间分辨率MODIS NDVI时间序列数据的面向对象分类方法可充分发挥2种遥感数据分别在空间、时间分辨率上的优势, 有效提高了主要土地覆盖类型的分类精度及总分类精度。

林业是我国高分辨率对地观测系统GF-1/6和GF-4卫星的主要用户之一, 本文发展的相关技术流程适用于综合GF-1/6与GF-4数据的土地覆盖类型分类制图, 有助于提高国家森林资源连续清查业务中林业用地、有林地等地类遥感分类制图的自动化、标准化程度。

4 讨论

随着高分系列卫星的成功发射及遥感影像数据的逐渐丰富, 综合应用多源数据影像分类方法将是遥感影像解译中重要的研究领域。虽然经过定性和定量评价分析, 验证了本文分类方法的有效性, 但还有很多值得进一步探索的问题:

1) 本文采用了一种基于中等空间分辨率多光谱影像分割“对象”综合应用多光谱、MODIS数据进行分类的技术框架。要综合利用这2种遥感数据源的优势, 还有其他技术框架, 如可将MODIS光谱反射率数据与GF-1 16 m宽幅多光谱数据相融合, 形成16 m的时间序列融合数据, 再以像元为分析单元进行监督分类。该技术框架与本文方法相比是否可得到更好的分类效果, 应在今后的工作中进一步分析评价。

2) 面向对象分类方法的一大优势, 即基于分割对象可以提取比基于像元更丰富的分类特征, 由于众多分类特征服从的统计分布一般都未知, 因此应首选非参数的分类器进行分类, 本文的分类结果也表明CART是一种性能优良的分类器。但是, 未采用其他非参数分类器进行尝试。因此, 比较评价多种分类器或多分类器融合方法的分类效果为下一步的研究方向。

3) 目前GF-4卫星已经成功发射, 于2016年7月开始获取数据, 考虑到数据处理方法、影像成像时间等问题, 预计2017年底即可获取1年的完整GF-4 NDVI时间序列数据, 后续应将本文基于MODIS NDVI替代数据源研发的技术方法尽快应用于GF-4数据, 并根据GF-4的特点对分类技术细节进行适当调整, 形成一套适用于我国森林资源宏观监测、融合GF-1和GF-4多光谱数据优势的地表覆盖类型(包括林地类型、森林类型)遥感分类制图方法。

5 结论

本研究发展了一种综合应用多源遥感数据空间和时间分辨率优势、采用决策树的面向对象土地覆盖分类方法, 该方法具有以“对象”为分析单元提取多源遥感数据的光谱、空间、时相等特征信息, 采用随机森林算法进行特征选择, 充分利用国家森林资源连续清查固定样地记录的地类信息训练分类器和检验分类效果等特点。对所发展的综合分类方法进行定量精度评价和分类效果比较分析, 表明本文所发展的多源遥感数据面向对象土地覆盖分类方法比分别基于GF-1宽幅多光谱数据、MODIS NDVI时间序列数据的分类结果精度有大幅提高, 有助于提高国家森林资源连续清查业务中林业用地、有林地等地类遥感分类制图的自动化、标准化程度。主要研究工作及结论如下:

1) 研究采用面向对象思想, 首先要得到最优的分析单元, 在多尺度分割过程中, 影像分割参数的选取一直是个难点。研究通过目视对比现有矢量数据及分割后面积统计, 探索出了适合研究区的分割尺度, 使得到的分割结果所含信息量最为丰富, 且能够最大程度地包含GF-1宽幅多光谱数据信息以及MODIS NDVI的物候信息。

2) 区别于传统基于知识的人工勾绘感兴趣区域方法, 研究中训练与检验样本是基于国家森林资源清查样地数据产生的。将影像分割后结果与清查样地点叠加, 得到面状森林资源清查样地数据, 结合外业调查样点数据, 对面状森林资源清查样地数据逐一筛选, 最终得到进入分类器以及分类结果检验的样本数据。通过面状样地数据, 统计分析各典型地类的NDVI时间序列曲线, 发现各地类间不同的物候节点信息, 如常绿林与落叶林可通过第25、39时相的NDVI区分出来, 这2个时相是落叶林生长和衰落的节点, 常绿林全年均值在0.70左右; 灌木林地与落叶林可通过非生长周期(第1—13、40—46时相)来区分, 此时灌木林地NDVI远小于落叶林, 均值在0.10左右; 耕地与灌木林地可通过第14、25时相信息来区分, 灌木林地的生长期为第14时相, 远早于耕地的第25时相。

3) 由于采用多源数据, 分类特征的维数很大, 属于高维数据分类问题。但由于获取客观真实的训练和精度检验样本数据需要消耗大量的人力、物力, 因此样本的数量总是有限的。对于有限训练样本集, 类别分布参数的估计精度会随着特征数量的增加而降低, 因此特征选择尤为重要。区别于传统基于知识的人工特征选择方法, 本研究采用随机森林算法, 对提取出的74个特征进行特征选择, 按照自变量出现在袋外数据时的模型MSE增量对特征进行排序, 并保证其筛选出来的每个自变量出现在袋外数据时的模型MSE增量均在20以上。特征选择后的特征个数为14个, 既包含GF-1宽幅数据特征, 又包含MODIS NDVI时间序列数据特征。

4) 本文所发展的综合中等空间分辨率GF-1宽幅多光谱影像和粗空间分辨率、高时间分辨率MODIS NDVI时间序列数据的面向对象分类方法, 可充分发挥2种遥感数据源在空间和时间分辨率上的优势, 有效提高了主要土地覆盖类型的分类精度及总分类精度。

参考文献(References)
贾明明, 任春颖, 刘殿伟, 等. 2014a. 基于环境星与MODIS时序数据的面向对象森林植被分类[J]. 生态学报, 34(24): 7167-7174.
(Jia M M, Ren C Y, Liu D W, et al. 2014a. Object-oriented forest classification based on combination of HJ-1 CCD and MODIS-NDVI data[J]. Acta Ecologica Sinica, 34(24): 7167-7174. [in Chinese])
贾明明, 王宗明, 张柏, 等. 2014b. 综合环境星与MODIS数据的面向对象土地覆盖分类方法[J]. 武汉大学学报:信息科学版, 39(3): 305-310.
(Jia M M, Wang Z M, Zhang B, et al. 2014b. Land cover classification of compositing HJ-1 and MODIS data based on object-based method[J]. Geomatics and Information Science of Wuhan University, 39(3): 305-310. [in Chinese])
郭芬芬, 范建容, 边金虎, 等. 2011. 基于MODIS NDVI时间序列数据的藏北草地类型识别[J]. 遥感技术与应用, 26(6): 821-826.
(Guo F F, Fan J R, Bian J H, et al. 2011. Grassland types identification based on time-series MODIS-NDVI data in northern Tibet[J]. Remote Sensing Technology And Application, 26(6): 821-826. DOI:10.11873/j.issn.1004-0323.2011.6.821 [in Chinese])
李俊祥, 达良俊, 王玉洁, 等. 2005. 基于NOAA-AVHRR数据的中国东部地区植被遥感分类研究[J]. 植物生态学报, 29(3): 436-443.
(Li J D, Da L J, Wang Y J, et al. 2005. Vegetation classification of east China using multi-temple NOAA-AVHRR data[J]. Chinese Journal of Plant Ecology, 29(3): 436-443. DOI:10.17521/cjpe.2005.0058 [in Chinese])
李玉东, 黄永喜. 2013. 三种TM与MODIS数据融合方法在山区的适用性研究[J]. 测绘与空间地理信息, 36(9): 124-133.
(Li Y D, Huang Y X. 2013. Application research of three fusion data methods of TM and MODIS in mountain areas[J]. Geomatics & Spatial Information Technology, 36(9): 124-133. [in Chinese])
刘丽, 匡纲要. 2009. 图像纹理特征提取方法综述[J]. 中国图像图形学报, 14(4): 622-635.
(Liu L, Kuang G Y. 2009. Overview of image textural feature extraction methods[J]. Journal of Image and Graphics, 14(4): 622-635. [in Chinese])
孙小添, 邢艳秋, 李增元, 等. 2013. 基于MODIS影像的决策树森林类型分类研究[J]. 西北林学院学报, 28(6): 139-144.
(Sun X T, Xing Y Q, Li Z Y, et al. 2013. Forest type classification by decision tree based on MODIS images[J]. Journal of Northwest Forestry University, 28(6): 139-144. [in Chinese])
潘琛, 杜培军, 张海荣. 2008. 决策树分类算法及其在遥感图像处理中的应用[J]. 测绘科学, 33(1): 208-211.
(Pan C, Du P J, Zhang H R. 2008. Decision tree classification and its application in processing of remote sensing images[J]. Science of Surveying and Mapping, 33(1): 208-211. [in Chinese])
张俊, 朱国龙, 李妍. 2011. 面向对象高分辨率影像信息提取中的尺度效应及最优尺度研究[J]. 测绘科学, 36(2): 107-109.
(Zhang J, Zhu G L, Li Y. 2011. Scale effect and optimal scale in object-oriented information extraction of high spatial resolution remote sensing image[J]. Science of Surveying and Mapping, 36(2): 107-109. [in Chinese])
章恒, 王世新, 周艺, 等. 2015. 多源遥感影像红树林信息提取方法比较[J]. 湿地科学, 13(2): 145-151.
(Zhang H, Wang S X, Zhou Y, et al. 2015. Comparison of different methods of mangrove extraction from multi-source remote sensing images[J]. Wetland Science, 13(2): 145-151. [in Chinese])
张晓娟, 杨英健, 盖利亚, 等. 2010. 基于CART决策树与最大似然比法的植被分类方法研究[J]. 遥感信息, 2(7): 88-92.
(Zhang X J, Yang Y J, Gai L Y, et al. 2010. Research on vegetation classification method based on combined decision tree algorithm and maximum likelihood ratio[J]. Remote Sensing Information, 2(7): 88-92. [in Chinese])
Breiman L. 2001. Random forests[J]. Machine Learning, 45(1): 5-32. DOI:10.1023/A:1010933404324
Gao F, Masek J, Schwaller M, et al. 2006. On the blending of the Landsat and MODIS surface reflectance:predicting daily Landsat surface reflectance[J]. IEEE Transactions on Geosciences and Remote Sensing, 44(8): 2207-2218. DOI:10.1109/TGRS.2006.872081
Li X B. 1997. Application of NOAA-AVHRR data on the research of land cover change[J]. Earth Sci Frontiers, 4(1/2): 1-6.
Thomas H, Michael A W, Nicholas C C, et al. 2009. A new data fusion model for high spatial-and temporal-resolution mapping of forest disturbance based on Landsat and MODIS[J]. Remote Sensing of Environment, 113(8): 1613-1627. DOI:10.1016/j.rse.2009.03.007
Wardlow B D, Egbert S L. 2008. Large-area crop mapping using time-series MODIS 250 m NDVI data:an assessment for the U[J]. S. Central Great Plains. Remote Sensing of Environment, 112(3): 1096-1116.