结合Sentinel影像与特征优选的山地城市不透水面提取

陈鑫亚，杨昆，王加胜

（1.云南师范大学信息学院；2.云南师范大学地理学部；3.西部资源环境地理信息技术教育部工程研究中心，云南昆明 650500）

0 引言

山区城镇化最显著的特征是大量不透水面取代了以植被为主的自然地表。不透水面是指道路、停车场等人工地表，不透水面能够阻止地表水下渗，从而破坏了地表能量平衡和生物多样性，容易导致洪涝灾害、热岛效应等问题。城镇化对山区生态环境的影响尤为显著，因此快速、准确地提取山区不透水面，可为优化山区发展规划、加强环境保护与治理提供决策支持。

遥感影像因其具有成本低、覆盖范围广等优点成为快速不透水面制图的重要数据来源。由于山区地物交错分布，低分辨率影像受空间分辨率和地物异质性的影响，混合像元占比较高。分辨率高的影像中的混合像元较少，Sentinel-2影像是目前免费的分辨率最高的光学影像，常被用于不同地物的遥感提取。受山区云雨天气的影响，光学影像的云量较高。同时，不透水面与裸土等地物容易出现光谱混淆的情况。因此，仅使用光学影像无法准确提取不透水面。有研究表明，融合Sentinel-2等光学影像和Sentinel-1合成孔径雷达（Synthetic Aperture Radar，SAR）影像等多源数据可以提高分类精度。

多源数据融合容易导致特征冗余，特征优选能够降低特征维度，提高分类性能。JM距离、随机森林（RF）和ReliefF算法是目前最为常用的3种优选方法。例如，文献［11］筛选出JM值处于一定范围的特征用于不同分类算法，以提取江苏夏收作物；文献［12］结合JM距离和RF算法对三七种植区进行遥感监测；文献［13］基于RF的袋外数据误差评估特征重要性，将特征从28个降至13个，以提取湖北洪湖湿地信息；文献［14］采用RF的平均不纯度减少衡量特征重要性后，对中南半岛进行了土地覆盖分类；文献［15］结合ReliefF算法和信息熵改进分离阈值算法进行特征优选，对华北平原地区作物进行分类；文献［16］使用基于ReliefF与粒子群优化算法的混合特征选择法用于土地分类，总体精度为85.88%。

目前，许多研究倾向于使用一种特征优选方法进行特征降维，并未对比不同优选方法在具体应用中的优劣和有效性。此外，大多数相关研究主要针对平原地区，特征优选方法在山区不透水面提取中的应用能力有待验证。为此，本文基于谷歌地球云计算平台（Google Earth Engine，GEE），以典型山地城市昆明市为例，比较了3种常用特征优选方法在不透水面提取中的有效性，并基于最优特征子集提取了昆明市2020年的不透水面。

1 研究区与数据源

1.1 研究区概况

昆明市（见图1）位于云南省中部，属低纬亚热带高原季风气候，全年干湿分明。地势由北向南呈阶梯状降低，大部分地区为海拔在1 500～2 800m之间的山区或半山区。云南省地处中国西南边境，与南亚、东南亚多国邻接，作为云南省的政治经济文化中心，昆明市城镇化水平发展迅速。

Fig.1 Overview of the study area图1 研究区概况

1.2 数据与预处理

（1）Sentinel-2 MSI数据。Sentinel-2携带多光谱成像仪（MultiSpectral Instrument，MSI），最高分辨率为10m。使用GEE上2020年Level-1C级别的影像共360景，数据已经过辐射校正和几何校正。筛选出云量低于10%的影像，并用具有云掩膜信息的QA60波段去云。为减少植被物候因素对分类的影响，将影像按归一化差异植被指数的最大值进行合成，并按研究区范围剪裁影像。

（2）Sentinel-1 SAR数据。Sentinel-1提供的合成孔径雷达数据最高分辨率为10m。为减少山区阴影和叠掩的影响，收集了2020年双轨道、双极化（VV+VH）和IW工作模式的影像共323景。影像都进行了热噪声去除、辐射校准和地形校正。将所有影像按均值进行合成，并按研究区范围进行剪裁。

（3）DEM数据。使用GEE上用于完成航天飞机雷达地形测绘任务的数字高程模型（Digital Eevation Model，DEM），数据分辨率为30m，将其按研究区范围进行裁剪。使用GEE的ee.Algorithms.Terrain（input）函数计算出海拔、坡度、坡向作为地形特征。

（4）样本数据。使用Google Earth进行线上样本采集。首先通过分层随机采样生成2 000个随机均匀分布的样本点，再基于Google Earth上2020年的影像，通过目视解译将样本标记为透水面和不透水面两类。透水面包括植被、水体和裸土，分别得到透水面和不透水面样本1 748与252个。

（5）土地覆盖分类数据。目前可免费获取的基于其他方法的土地覆盖数据有限，因此收集了2020年30m分辨率的全球地表覆盖数据集（GLC30）、2018年30m的全球人工不透水面数据集（GAIA）和2017年10m的全球土地覆盖观测与监测数据集（FROM-GLC10）用于结果精度验证，并将产品重分类为透水面和不透水面两类。

2 研究方法

本文基于GEE平台，首先获取并对遥感数据进行预处理，完成光谱、光谱纹理、地形及SAR纹理特征构建；然后，分别基于JM距离、RF模型、ReliefF算法优选出3组特征子集，评价不同特征子集的分类精度；最后，基于最佳的特征优选方法与RF分类算法提取昆明市2020年的不透水面分布信息。

2.1 特征构建

（1）光谱特征。通过Sentinel-2影像提取的光谱特征包括蓝波段（Blue）、绿波段（Green）、红波段（Red）、近红外波段（NIR）、短波红外（SWIR1、SWIR2）以及归一化差异建筑指数（NDBI）、改进的归一化差异水体指数（MNDWI）、归一化差异植被指数（NDVI）和裸土指数（BSI）。

（2）光谱纹理特征。纹理特征能够反映地物大小形态、疏密程度和规整程度。鉴于不透水面在短波红外波段中具有较高反射率，基于灰度共生矩阵对Sentinel-2影像的B12波段计算纹理特征。使用GEE函数glcmTexture（）得到常用的5个纹理特征：角二阶矩（B12_ASM）、对比度（B12_Contrast）、逆差矩（B12_IDM）、熵（B12_ENT）和方差（B12_VAR）。

（3）地形特征。地形因素通常限制着山区城市的分布与扩张。因此，采用由DEM数据计算得到的海拔（ELEVATION）、坡度（SLOPE）和坡向（ASPECT）3个地形特征分量参与特征集构建。

（4）SAR纹理特征。相关研究表明，SAR的纹理特征相比光谱纹理特征具有更丰富的纹理信息。本文通过glcmTexture（）函数计算得到VV和VH极化方式的常用纹理特征：VV_asm、VV_contrast、VV_idm、VV_ent、VV_var、VH_asm、VH_contrast、VH_idm、VH_ent、VH_var。

2.2 特征优选方法

特征选择即将高维空间的特征通过映射或变换方法转换到低维空间，剔除冗余和不相关特征。特征选择的原则是在不显著降低分类精度、不影响类分布的前提下优选出尽可能小且稳定的特征子集。本文构建了4类特征：光谱特征10个，光谱纹理特征5个，地形特征3个，SAR纹理特征10个，共28个原始特征。若不进行优选，冗余特征将会影响分类性能。

2.2.1 基于JM距离的特征优选

JM距离（Jeffries-Matusita distance）是反映类间可分性的重要指标。在特征正态分布的前提下，计算类间样本的JM距离可衡量不同类别基于某一特征的可分离性。对于类别

i

和类别

j

，其JM值的计算方法见式（1）和式（2）。

2.2.2 基于RF模型的特征优选

随机森林（Random Forest，RF）是集成多棵决策树的学习算法，每棵决策树通过随机选取样本和特征独立生成。决策树节点在分裂时，由最优分裂特征得到的子节点样本在理想状态下属于同一类。利用基尼指数（Gini index）衡量样本集合的不纯度，Gini index越小，样本的错分概率越小。基于RF模型的特征优选通过对比特征在决策树上的平均贡献率来评估特征重要性，通常使用Gini index或袋外数据错误率作为衡量标准。本文基于Gini index进行评估，并将所有特征贡献率的均值作为阈值，优选出贡献率大于均值的特征用于分类。

2.2.3 基于ReliefF算法的特征优选

ReliefF特征优化算法能够解决分类及回归问题，算法效率高且不限制数据类型，为较常用的特征选择方法之一。ReliefF算法每次从训练样本集中随机取出一个样本

x

，并从与

x

同类的样本集中找出

x

的k个近邻样本，从与

x

不同类的样本集中找出k个近邻样本，然后更新每个特征权重，重复m次。特征权重计算方法见式（3）。

式中，

Diff

（）表示样本在特征

f

上的距离，

H

（

x

）、

M

（

x

）分别为与

x

同类和非同类的最邻近样本，

p

（

x

）表示类别

x

的概率，

m

为迭代次数，

k

为最近邻样本个数。特征权重越大，表示该特征对样本的区分效果越好。本文将所有特征权重均值作为阈值选择最优特征子集。

2.3 分类方法

本文运用RF分类算法提取不透水面。RF是对多棵决策树投票结果进行集成的分类器，采用bootstrap抽样策略选取训练样本用于生成决策树，具有分类精度较高、不易出现过拟合等优点，较少受噪声和特征选择的影响，适用于多类别与多特征分类。决策树数量和特征选取个数是RF建立的关键参数，通过以50为间隔调整参数的方法，选取分类精度最高的决策树数量。通过多次实验，将决策树数量设置为200，特征个数设置为特征维数的平方根。

2.4 精度评价方法

通过构建混淆矩阵可对基于不同优选方法的分类结果进行精度评价，具体评价指标有生产者精度（Producer′s Accuracy，PA）、用户精度（User′s Accuracy，UA）、总体精度（Overall Accuracy，OA）和Kappa系数。其中，OA和Kappa系数可反映总体分类精度，PA、UA用来衡量不透水面的错分和漏分误差。运用交叉验证方法，随机选取约7/10的样本用于模型训练，约3/10的样本作为测试样本用于精度验证。对每个特征子集进行10次精度验证，并对10次验证结果求平均值。

3 结果与分析

3.1 特征优选结果

3.1.1 基于JM距离的特征优选结果

计算出单个特征的JM距离以及同类特征JM距离的均值，并将同类特征中JM距离低于均值的特征剔除，完成同类特征优选（见表1），再将优化后的不同类特征进行组合，计算JM距离（见表2）。

Table1 Features of the same category and their JM distance after optimization表1 优选后的同类特征及其JM距离

Table 2 JM distance of combinations of features belong to different categories表2 不同类特征组合的JM距离

由表1可知，对同类特征进行优选后，共减少了13个特征。表2表明将不同类特征进行组合能增大JM距离，纹理特征的贡献率较高。其中，S+T+ST和S+T+ST+D的JM距离最大。由于S+T+ST特征数量更少，因此选择Blue、Green、Red、SWIR2、NDVI、BSI、B12_asm、B12_contrast、B12_var、VV_asm、VV_var、VH_contrast、VH_var共13个特征用于构建特征子集。

3.1.2 基于RF模型的特征优选结果

基于RF模型的特征重要性如图2所示，计算出特征重要性均值为0.035 7。其中，B12_idm、B12_asm、B12_ent、B12_contrast、VV_contrast、BSI、NDVI、Blue的贡献率超过3%；Red、NDBI、SLOPE、VV_idm、VV_ent、VV_asm和VH_contrast的贡献率在2.5～3%之间，表明纹理特征的贡献率最高，且VV纹理特征贡献率高于VH纹理特征。贡献率高于均值的特征包括Blue、NDVI、BSI、B12_asm、B12_contrast、B12_ent、B12_idm、VV_contrast共8个特征，用于构建特征子集。

Fig.2 Feature contribution rate based on RF model图2 基于RF模型的特征贡献率

3.1.3 基于ReliefF算法的特征优选结果

基于ReliefF算法的特征权重如图3所示，特征权重的均值为0.017。NDVI、NDBI、BSI、B12_ent、B12_idm的权重超过3%，Blue、VV_ent、VV_idm、ASPECT、SLOPE的权重处于2～3%之间，其余特征的权重低于2%。因此，权重高于均值的特征包括Blue、NDBI、NDVI、BSI、B12_ent、B12_idm、VV_ent、VV_idm、SLOPE、ASPECT共10个特征，用于构建特征子集。

Fig.3 Feature weight based on ReliefF algorithm图3 基于Relief F算法的特征权重

3.2 最优特征选择方法判定

基于JM距离、RF模型、ReliefF算法构建的特征子集分别用于训练RF分类器，并提取不透水面，将提取结果分别命名为JM_2020、RF_2020和ReliefF_2020。基于不同特征子集的分类精度和分类模型训练时间见表3。为进一步展现3种方法在不透水面提取上的差异，选取昆明市3个典型区域，对比JM_2020、RF_2020和ReliefF_2020的分类精度，见图4。

Table 3 Evaluation index of different feature selection methods表3 不同特征优选方法评价指标

Fig.4 Extraction results of impervious surface based on different feature optimization methods图4 基于不同特征优选方法的不透水面提取结果

由表3 可知，相比原始特征，特征优选剔除了大量冗余或无关特征，提高了分类效率和精度。其Kappa 值和PA值有所提高，模型训练时间骤减。基于ReliefF 算法的优选特征数量较少，相比JM 距离和RF 模型，Kappa 值分别提高了2%和3%，PA 值分别提高了7%和4%，模型训练时间分别减少了84s和16s。

分析图4 可知，（a）区域为居民区，RF_2020 和JM_2020 在居民楼密集区都存在不透水面高估现象。在（b）、（c）区域中，JM_2020的不透水面低估情况相比其他两种方法更严重，漏分了大量不透水面。根据Sentinel-2 的地面真实影像，与RF_2020 和JM_2020 相比，ReliefF_2020 更接近不透水面的真实分布，提取效果明显优于其他两种方法。因此，本研究将ReliefF 作为最优的特征选择方法，ReliefF_2020为最优的不透水面提取结果。

3.3 最优提取结果与其他土地覆盖数据比较

为进一步验证ReliefF_2020 的精度，选取昆明市内具有不同城镇化程度的3 个典型区域，将ReliefF_2020 与基于其他方法获取的土地覆盖数据精度进行比较，见图5。

区域（a）中FROM_GLC10和GAIA 都存在低估现象，区域（b）中FROM_GLC10 和ReliefF_2020 精度较高，区域（c）中的大棚作物种植区（Sentinel-2 中亮白色区域）易与不透水面混淆，FROM_GLC10 和GAIA 中的混淆现象严重。GLC30 具有最低精度，ReliefF_2020 与Sentinel-2 目视解译效果基本一致。因此，ReliefF_2020 的精度仅次于FROM_GLC10，但高于GAIA 和GLC30。

Fig.5 Accuracy comparison of land cover data图5 土地覆盖数据精度比较

4 结语

本研究融合多源遥感数据，对比了3 种特征优选方法（JM 距离、RF 模型和ReliefF 算法）在山地城市不透水面提取中的有效性，得到以下结论：

（1）融合光学卫星和雷达卫星数据能有效提高山区不透水面的提取精度。在所有特征中，光谱纹理特征与SAR 纹理特征的重要性较高。

（2）特征优选方法能减少大量冗余特征，减少模型训练时间，提高分类效率和精度。基于ReliefF 算法优选出的特征较少，模型训练时间最少，不透水面提取结果具有最高的Kappa和PA 值，因此成为最优的特征优选方法。

（3）最终提取结果中存在不透水面的高估与低估现象，主要与特征优选时阈值的不确定性、研究区的特殊性、样本的准确性，以及山区的光谱混淆、混合像元等问题有关。

以特征权重均值作为阈值具有一定的不确定性。未来研究将根据特征类型与数量，对比不同阈值的优选效果或寻找自动化的阈值确定方法。同时，将样本的线上选取与线下检验相结合以提高样本准确度。