李晶晶,邢 鹏,屈鸿钧,杜凤兰
(1.二十一世纪空间技术应用股份有限公司,北京 100096;2.二十一世纪(广州)空间技术应用有限公司,广东 广州 510640)
为积极促进节约集约用地,确保已供土地合理开发利用,需做好土地供后监管工作。其中,已供国有建设用地开发建设情况的状态监测是土地供后监管的一项重要内容[1]。遥感影像以覆盖面积大、获取周期短的优势被广泛应用于变化检测[2]。高分辨率遥感影像具有空间细节清晰、丰富的特点,适用于检测特定地物(如建筑物)的变化[3],为建设用地变化检测提供了良好的契机与条件。
针对基于遥感影像的建设用地变化检测方法,早期学者多采用传统机器学习算法,如李悦[4]等结合旧时期土地利用矢量数据,采用随机森林和模糊C均值聚类方法提取了建设用地变化信息;Gamanya R[5]等以Landsat TM/ETM多时相遥感影像为数据源,采用基于模糊逻辑的面向对象自动分类方法提取了津巴布韦哈拉雷地区的建筑物目标;裴艳艳[6]等以Landsat TM/OLI遥感影像序列为数据源,采用CART分类方法分析了建设用地的变化情况;刘红超[7]等以Landsat8-OLI影像为数据源,提出了一种基于土地覆盖类型特征自适应确定阈值的遥感影像变化检测方法。上述机器学习方法是单一尺度的面向对象方法,适用于目标地物尺度接近的情况,否则将存在过分割或不完全分割的情况[8],导致提取结果中大尺寸地物过于破碎或小尺寸地物漏提。近年来深度学习技术得到快速发展,逐步应用于遥感影像建筑物变化检测领域,如吴海平[9]等结合全国土地利用遥感监测工程的实际需求,尝试利用深度学习方法对新增建设用地进行了自动提取,深度学习技术在新增建设用地变化提取方面具有更高的适用性和实用性;顾炼[10]等基于高分辨率遥感影像,采用FlowS-Unet深度学习方法实时、准确地检测区域内新建与扩建的建筑物;张翠军[11]等改进了U-Net模型,将变化检测问题转化为像素级二分类问题,提高了遥感图像检测建筑物变化的精度。然而,深度学习方法需具有足够数量的提取目标标注样本,且难以充分利用专家先验知识、遥感成像机理、遥感影像附载的地理学知识等解译知识[12]。本文提取的目标为新增建设用地,相较于其他建筑类型,建设用地在不同开发阶段呈现不同的形状、光谱、纹理特征,尤其在高分辨率遥感影像上,建设用地类型多样,小建筑物多且破碎,而且可用的样本集不多,传统机器学习和深度学习方法不适用于该类新增建设用地的提取。
鉴于此,本文提出了一种多层级双尺度的新增建设用地自动提取方法。该方法根据建设用地在各开发阶段的不同特征和类型,建立了新增建设用地类型层级,并选择不同尺度进行影像分割,从而通过不同分割尺度提升新增建设用地的提取精度。本文以北京二号遥感影像为实验数据,通过比较本文方法与两种单尺度决策分类树方法的提取结果精度,验证本文方法的有效性。
研究区位于广东省中山市板芙镇,地理范围为22°23′~22°24′N、113°17′~113°18′E,区域内地物类型包括建筑、水面、草地、道路、裸地等。实验采用的影像数据由北京二号卫星分别于2020年4月和2021年4月拍摄,空间分辨率为0.8 m,包括Blue(B)、Green(G)、Red(R)、Near Infrared(NIR)4个波段,位数为16 bit,实验区域大小为2 480×3 180像元。该实验区域按照R、G、B波段顺序合成的影像如图1所示,数据已完成几何校正、辐射校正、匀光匀色等处理,满足遥感智能解译的要求,图1a为变化前影像,即前时相影像,图1b为变化后影像,即后时相影像。此外,以人工提取结果为基准进行精度评价,人工提取结果由生产人员参考本文构建的新增建设用地类型层级和指定的上图面积,采用GIS软件人工勾画生成,经质检合格。
图1 研究区遥感影像
本文提出的多层级双尺度方法是采用面向对象的直接比较方法,根据新增建设用地的特征构建决策分类树进行分层信息提取。其具体思路为:①分析新增建设用地的特征,建立新增建设用地的类型层级;②分别选择适合各类新增建设用地的分割尺度,并对两幅影像同时进行分割,生成对应的影像对象层;③结合各类新增建设用地的特点,构建多特征约束的提取规则,实现不同分割尺度下的新增建设用地提取;④以人工提取结果为基准数据,进行精度评价。具体提取流程如图2所示。
图2 新增建设用地提取流程图
已供建设用地开发状态一般可分为场地平整、基础施工、地上主体开建和主体封顶4个阶段。根据研究区两期北京二号遥感影像上新增建设用地的特征,结合已供建设用地开发状态的监测需求,同时便于新增建设用地自动提取规则的特征描述,本文采用多层级方式提取新增建设用地,将新增建设用地作为一级类型,将新增彩钢房、新增高亮构筑物或硬化地面、绿地变裸地3类作为二级类型,将新增彩钢房再细分为新增红色彩钢房和新增蓝色彩钢房,如表1所示。
表1 新增建设用地的类型层级
通过影像分割来获取用于特征分析的影像对象层,由于实验区内新增建设用地各类型特征存在较大差异,因此需要设置多个分割尺度来实现不同类型新增建设用地的准确提取。本文采用的多分辨率分割算法[13],是一种基于异质性最小原则的自下而上的区域合并算法,参数包括波段权重、异质性因子和分割尺度。波段权重表示各波段被利用的信息比例,考虑到实验影像4个波段对新增建设用地提取具有同等重要作用,将波段权重均设置为1;异质性因子包括光谱因子、紧致度因子和平滑度因子,定义了影像对象合并生长的准则,根据经验[14],将光谱因子设置为0.9、紧致度因子设置为0.5、平滑度因子设置为0.5;分割尺度表示影像分割允许的影像对象的最大异质性,决定分割对象是否停止合并。本文通过固定波段权重、异质性因子,动态调整分割尺度,目视分析来确定合适的分割尺度。
新增建设用地特征提取需在各类新增建设用地对应的影像对象层上分别选择合适的特征组合,构建基于多特征约束的分类决策树提取规则,实现目标地物的提取。通过分析目标地物的光谱/形状特征分布图和特征值发现:
1)新增红色彩钢房影像对象在后时相影像R波段、前后时相影像R波段差值上,与其他影像对象具有明显差异;新增蓝色彩钢房对象在后时相影像B波段、前后时相影像B波段差值上,与其他影像对象具有明显差异;此外实际生产中对新增建设用地图斑最小上图面积有一定要求。因此,本文采用后时相影像对象R波段特征(R2)、R波段特征差值ΔR和新增图斑面积Area作为新增红色彩钢房的特征组合,这里新增图斑面积是指同类型相邻的新增影像对象合并后的图斑面积,下文同理;采用后时相影像对象B波段特征(B2)、B波段特征差值ΔB和新增图斑面积Area作为新增蓝色彩钢房的特征组合。
式中,meanR2、meanG2、meanB2分别为后时相影像对象在R、G、B波段的均值;meanR1、meanG1、meanB1分别为前时相影像对象在R、G、B波段的均值。
2)新增高亮构筑物或硬化地面影像对象在后时相影像亮度、前后时相影像亮度差值上,与其他影像对象具有明显差异。因此,本文采用前后时相影像对象亮度差值ΔHSI_I、后时相影像对象地物亮度特征Brightness2以及新增图斑面积Area作为该类型的特征组合。
式中,HSI_I1、HSI_I2分别为前后时相影像对象在HSI色彩空间的亮度分量;meanR2、meanG2、meanB2、meanNIR2分别为后时相影像对象在R、G、B、NIR波段的均值。
3)绿地变裸地影像对象在前时相影像植被指数、后时相影像植被指数、后时相影像土壤指数、后时相标准差均值上,与其他影像对象具有明显差异。因此,本文采用前时相影像植被指数NDVI1、后时相影像植被指数NDVI2、后时相影像土壤特征NDSI2、后时相影像标准差均值SD2和新增图斑面积Area作为特征组合。
式中,meanR1、meanNIR1分别为前时相影像对象在R、NIR波段的均值;meanR2、meanB2、meanNIR2分别为后时相影像对象在R、B、NIR波段的均值;SD G2、SD B2、SD R2、SDNIR2分别为后时相影像对象在G、B、R、NIR波段的方差。
为进一步验证多层级双尺度提取方法的有效性,本文进行以下实验分析。
1)确定分割尺度。通过初步实验分析,将分割尺度的取值范围设为20~150,步长设为10,其他参数设为§2.2中的固定值,分别对实验影像进行多分辨率分割,并分析分割效果。合适的分割尺度具有的特点为分割对象大小与地物目标(新增建设用地)大小尽量接近,特定地物类型可用一个或几个对象来表达,相同类型对象的异质性差异性较小,地物的边界轮廓不能太模糊[15]。以此为评价标准,本文确定的分割尺度为:新增彩钢房的分割尺度为30,新增高亮构筑物或硬化地面、绿地变裸地的分割尺度为100。
2)提取新增建设用地。根据上述确定的分割尺度和各新增建设用地类型的特征组合,提取各类新增建设用地。提取规则如表2所示,其中面积单位为m2,面积大小为经验值,可根据实际需要调整。合适的特征阈值是保证提取结果精度的关键,本文以所选择的阈值能够检测出大部分目标类型为阈值确定的基本原则[16],采用样本极值法确定特征阈值,再利用目标类型的多个特征依次筛选影像图斑,完成新增建设用地的提取。
表2 新增建设用地提取规则
样本极值法确定阈值要求选取样本时,尽量选取靠近目标与非目标临界点的样本(每类选取有代表性的样本即可)。根据不同类型新增建设用地对应的特征,本文采用样本的最大值或最小值作为阈值,样本极值法的计算公式为:
式中,Q为确定的阈值;k1、k2、…、k n分别为样本1、样本2、…、样本n在某一特征上的值;max()为取最大值函数;min()为取最小值函数。
两组新增彩钢房与所在位置的2020年影像、2021年影像叠加自动提取结果、2021年影像叠加人工提取结果如图3所示,新增高亮构筑物或硬化地面以及绿地变裸地与所在位置的2020年影像、2021年影像叠加自动提取结果、2021年影像叠加人工提取结果如图4所示。结合人工提取结果可知,本文方法可有效发现新增建设用地。
图3 新增彩钢房提取结果
图4 新增高亮构筑物或硬化地面和绿地变裸地提取结果
以人工提取的新增建设用地结果为基准数据,本文采用目视对比方法对自动提取结果进行精度评价。结合已供国有建设用地开发建设情况的状态监测需求,本文的研究目标是采用自动提取方式辅助发现新增建设用地。因此,精度评价采用自动提取结果与人工提取结果是否相交的方法来分析发现变化的准确性。另外,为了评价本文方法的有效性,采用单尺度决策分类树方法开展对比实验,分割尺度分别为30和100,具体统计结果如表3、4所示,可以看出,若考虑小尺度地物的提取,采用小尺度分割,则提取结果过于破碎,且会将一些小图斑错分为新增高亮构筑物或绿地变裸地,如单尺度方法(分割尺度为30)提取结果的图斑个数为145,新增高亮构筑物或硬化地面、绿地变裸地的图斑数量大于另外两种方法,单尺度方法(分割尺度为30)的整体正确率为83.45%,新增高亮构筑物或硬化地面、绿地变裸地的正确率明显低于本文方法;若考虑大尺度地物的提取,采用大尺度分割,则会漏提小尺度地物,如单尺度方法(分割尺度为100)中尺度较小的新增红色彩钢房的遗漏率为16.67%,高于另外两种方法,由于实验数据中新增彩钢房数量较少,若数量增加,遗漏现象将更加明显。本文方法采用的多个分割尺度能较好地改善上述问题,根据建立的新增建设用地类型层级,选择适合各类型的分割尺度,提取结果能用尽量少的图斑表示目标地物,同时可保证较高的正确率和较低的遗漏率。
表3 新增建设用地提取结果统计
本文以北京二号遥感影像为数据源,根据建立的新增建设用地类型层级选择不同的分割尺度,采用不同的特征组合,实现了新增建设用地的自动提取。多层级双尺度的新增建设用地变化监测方法能在一定程度上解决单尺度方法分割不充分或过度分割导致的问题,可在保证高正确率和低遗漏率的情况下,实现不同类型新增建设用地的自动提取,在土地供后监管的建设用地变化发现中具有一定的应用价值。