钟嫣然,李晓龙
(1.武汉大学 测绘学院,湖北 武汉430079;2. 重庆市地理信息中心,重庆 401120)
及时准确的城市扩张度量对于了解一个城市的经济发展情况、协调城乡之间的发展平衡是非常必要的[1]。遥感技术作为一种多时相、覆盖广、快速实时的信息获取手段,近年来在资源环境、城市扩张监测等方面得到了广泛应用[2]。迄今为止,已有许多学者基于不同分辨率的可见光与红外波段遥感影像、夜间灯光遥感数据、合成孔径雷达数据等,提出或使用了多种建成区即“不透水面”[3]的提取方法,本文主要讨论多波段影像数据所使用的方法。其中,较常见的算法是各种光谱指数,如归一化差异建筑指数(NDBI)[4]、建成区提取指数(BAEI)[5]、基于指数的建筑指数(IBI)[6]和生物物理成分指数(BCI)[7]。然而,由于地物间存在同物异谱、异物同谱以及混合像元的情况,使得这些指标在区分建成区和背景,尤其是裸土时,并没有获得较为理想的精度[8]。除此以外,也有基于纹理特征提取的算法,如Pesaresi M[9]等针对SPOT 卫星的2.5 m分辨率全色波段影像提出了适用于5 m 甚至更低分辨率影像的建成区存在指数PanTex;沈小乐[10]等利用采样轮廓波变换来描述图像的纹理特征,并引入视觉注意机制,最终得出该方法适用于高分辨率遥感影像的结论;NING X G[11]等以甚高分辨率遥感影像为对象,通过探测几何元素(如直角边和角)并计算其密度来提取建成区。
上述方法主要是通过提取影像特征来达到提取建成区的目的。近年来,随着机器学习的兴起,能自行训练分类模型的机器学习监督分类算法成为研究热点。陈亮[12]基于Landsat TM/ETM+影像,提出了基于最小二乘支持向量机模型(GA-LS-SVM)和遗传算法的不透水面提取方法,结果较理想;饶萍[13]等基于高空间分辨率影像,以4 种归一化指数、SVM 分类结果、河流缓冲区结果为决策树分类器的特征结点,提取了Landsat 8 影像中的建设用地,其精度高于其中任意单独特征的分类结果。然而,过量且不必要的数据可能反而会使分类精度下降[14]。对于建成区提取而言,光谱特征仍是许多研究中选择采用的特征,原因在于不透水面的反射特征与植被、水体等地物的差异较明显,其提取结果有很大的参考价值。此外,建成区是由若干典型尺寸范围的物体和投射阴影组成的,纹理特征相对明显[9],因此本文选取光谱特征和纹理特征对建成区提取进行分析。同时,为了比较单一特征和多特征提取的效果,本文将两种特征结合起来作为第三种提取方法,并根据总体精度和局部效果分析比较了这3 种方法的特点。
本文选取的研究区为广东省广州市南沙区的局部区域。随着近年来GDP 的增长,该区域内的建成区组成由单一的农村住宅群转化为新城区与旧农村房屋并存,区域内包括林地、农田、水域、建成区和裸土5 种地类,种类较齐全。
本文采用的实验数据为2017 年10 月的Landsat8-OLI 影像和同年12 月的Sentinel-2A 影像,成像时间较接近,区域内建成区基本没有较大变化。Landsat 8 卫星的OLI 传感器能获得空间分辨率为15 m 的全色波段以及空间分辨率为30 m 的其他短光谱段影像,因此对其他波段数据和全色波段进行融合可使其空间分辨率提高至15 m。Sentinel-2A 卫星可获得可见光、近红外和短波红外波段的影像,前两者空间分辨率为10 m,后者空间分辨率为20 m,在免费公开数据中属于中高分辨率影像[15]。本文采用的Landsat 8 和Sentinel-2A数据均已经过几何和大气校正、重采样、图像裁剪和灰度拉伸等预处理,得到的数据如图1 所示,其中Sentinel-2A 影像的空间分辨率为10 m,Landsat8 多波段影像融合后的空间分辨率为15 m,没有融合的为30 m。另外,本文的建成区参考结果为Sentinel-2A 影像的人工目视解译提取样本,最后得到的样本数据为2 000 个建成区矢量点和2 000 个非建成区矢量点,并将其按2∶3 的比例分别随机选取800 个点作为训练决策树模型的样本,另外1 200 个点参与分类后的精度评定。
图1 研究区影像图
本文首先计算影像的CBI 和PanTex,再把二者作为决策树分类法的节点进行提取,最后得到3 种方法的提取结果,如图2 所示。根据Valdiviezon J C[8]等的实验,对多种建筑指数计算裸土和建成区之间的光谱区分指数(SDI)可知,CBI 在这两类地物之间的SDI 值是最高的,因此本文采用CBI。另外,本文研究区为小面积住宅居多、建筑物既有集中又有分散的农村地域,且实验数据最高空间分辨率为10 m,建筑物个体在图上难以清晰展现,因此几何特征的提取方法采用PanTex。特征集合方法采用决策树分类法,因为样本点数量不多,采用树形算法可能比神经网络算法的效果要好。
图2 研究流程图
CBI 由主成分变换得到的第一主分量(PC1)、调整土壤亮度的植被指数(SAVI)和归一化差异水体指数(NDWI)3 个子指数组成[16],其中主成分变换是一种常用的数据降维、影像增强方法,利用正交变换原理得到原图像的主要信息(如亮度、绿度、湿度等),尤其是对应最大方差的第一主分量[17];SAVI 与经典的NDVI 相比,更能从裸土背景中提取出较为破碎的植被斑块[18],有利于削弱建成区内小范围植被的影响。其计算公式为:
式中,L为土壤调节系数,取值范围为0~1,对应植被覆盖度由低至高[19]。由图2 可知,研究区内有不少被建成区割碎的块状植被,因此这里L取0.5。
NDWI 的计算公式为:
在计算CBI 之前,需先用水体指数对上述各指数结果进行掩膜提取,裁掉水体部分。实验表明,改进的归一化差异水体指数(MNDWI)计算结果中的建成区基本为负值[20],其效果比NDWI 更理想,因此本文采用MNDWI 来提取水体部分。其计算公式为:
裁掉水体后,将3 个指数进行归一化,再进行CBI 计算。其计算公式为:
最后得到的CBI 计算结果介于-1~1 之间,其中正值部分为建成区,接近0 的部分为裸土,负值部分为植被[16]。
PanTex 基于一幅灰度图的灰度共生矩阵,根据模糊规则对多方向的纹理特征值进行融合,即进行旋转不变的各向异性纹理分析,再根据纹理特征值提取建成区。
对于Sentinel-2A 数据,本文选取可见光波段中各像元的最大灰度值作为该像元在新的灰度图中的值,某些研究会选取可见光3 个波段的平均值或对图像作HIS 转换等形式得到新的灰度图,但考虑了一些不必要的波段而导致实验结果不太理想[21];对于Landsat 8数据,可直接使用经过线性拉伸的全色波段。
计算GLCM 时包括滑动窗口大小、矢量位移(方向和距离)、图像灰度级[9]3 个关键参数。为了避免数据量过大,在预处理时将Sentinel-2A 影像从16 位降为8 位整型;根据图中建筑斑块的面积,本文选择的窗口大小为3×3;方向选取了常用的8 个方向,即0°、45°、90°和135°及其对角线方向;距离的选择取决于纹理的粗细,本文设定距离d=1,着眼于最邻近像元的灰度变化情况。
Pesaresi M[21]等通过对不同空间分辨率的图像进行实验发现,对比度(CON)是最能把建成区和非建成区区分开的特征量,计算公式为:
式中,Ng为图像灰度值的最大值;i、j为GLCM 中的行列号,即任意两个灰度值。
若要获得每个像元对比度的唯一值,需对所有方向的值进行模糊逻辑处理。实验证明,相比于最小值,平均值和原始各方向值都倾向于使农田、乡村道路等非建成区呈现较高的亮度值[9,21],因此本文采用所有方向值中的最小值作为最终的对比度值。
决策树分类法是机器学习中的一种基于空间知识挖掘、无需先验知识的监督分类方法,可对训练样本进行归纳学习,自行生成决策树或决策规则,即从“无序”中找出“有序”的规则,再对输入的数据进行分类。决策树可用树型结构表示,每个决策树包含一个根节点、若干个内部结点和叶节点,每个内部结点对应一个非类别属性或属性的集合。本文将计算得到的CBI 和PanTex 结果作为节点输入到决策树分类器中,再进行整个图像的分类,上述过程均利用Python 代码实现。
基于3 种空间分辨率数据,本文分别利用3 种方法提取的建成区结果如图3~5 所示,可以看出,双特征决策树方法提取的结果区域间连通性最强,其次为CBI,最后为PanTex。
图3 Sentinel-2A 数据(10 m)提取结果
基于Sentinel-2A 数据提取的建成区结果局部图如图6 所示,蓝色为PanTex 提取结果,紫色为CBI 提取结果。由图6 上方3 幅图可知,PanTex 提取结果与不规则分布的农村住宅群较为吻合,且它对建成区和裸土的区分度优于CBI 提取结果;CBI 作为一种光谱指数,其结果难免掺杂不少非建成区;PanTex 对道路等条状地物的敏感度不高,是由纹理不明显导致的。由图6 下方3 幅图可知,PanTex 提取结果里存在部分高亮建筑没有被识别为建成区,这是由于高亮建筑在滑动窗口内成分单一、亮度变化不明显所造成的;采用双特征决策树分类法可以改善PanTex 的这种现象。
图4 Landsat8-OLI 融合后数据(15 m)提取结果
图5 Landsat8-OLI 数据(30 m)提取结果
图6 基于Sentinel-2A 影像的提取结果局部图
基于3 种数据,3 种方法得到的9 个结果的总体精度如表1 所示。总体而言,双特征决策树分类法和PanTex 的提取结果均优于CBI,即按总体精度对3 种方法的排序为:双特征决策树分类法>PanTex>CBI。此外,图像空间分辨率越高,提取结果精度越高。
表1 提取结果的总体精度
本文以2017 年Landsat8-OLI 和Sentinel-2A 的影像数据为数据源,以广州市南沙区内局部区域为研究区,利用CBI、PanTex 和双特征决策树分类法3 种方法,基于相同时间内3 种空间分辨率的数据,得到了9 个提取结果。结果表明:①3 种方法提取结果的精度排序为双特征决策树分类法>PanTex>CBI;②在数据源和提取方法相同的情况下,PanTex 提取结果的面积均小于其他两种方法,斑块较破碎;③影像的空间分辨率越高,提取结果的总体精度也越高;④PanTex在识别占地面积较小的农村住宅以及区分裸土和建成区方面更优,而CBI 识别高亮度、面积较大建筑和道路的效果更好。
本文没有专门比较裸土与建成区CBI 和PanTex 的差别规律,在后续研究中可对二者作较详细的分析;两个指标都是通过重复试验来自定义阈值的,因此分类结果易受解译者的主观影响。此外,采用双特征决策树分类法得到的结果中仍包含错分区域,可能是由于非建成区样本在其他几类地物中分布不均匀,导致缺乏更全面的信息,使得建成区与其他地物区分开来。本文仅用CBI、PanTex 结果代表建成区的物理特征和几何特征是较片面的,在后续研究中应各自多选取几种指标来综合比较,同时增加其他类型的数据(如DEM、水文信息等)作为决策树分类法的节点,使提取效果更优。