基于Sentinel-2A和Landsat8的城市不透水面的提取
Impervious surface extraction based on Sentinel-2A and Landsat8
通讯作者: 钟凯文(1972-),男,学士,研究员,研究方向为城市遥感。Email:zkw@gdas.ac.cn。
责任编辑: 张仙
收稿日期: 2020-07-16 修回日期: 2020-10-29 网络出版日期: 2021-06-15
基金资助: |
|
Received: 2020-07-16 Revised: 2020-10-29 Online: 2021-06-15
作者简介 About authors
赵 怡(1992-),女,博士研究生,研究方向为城市环境遥感。Email:
城市不透水面相关研究对城市的发展至关重要。为了提取城市不透水面盖度,通常采用线性光谱混合分析方法,在亚像元尺度上计算混合像元内的不透水面面积比例。由于端元光谱曲线存在误差,导致不透水面盖度提取精度较低,因而提出端元优化方案,通过Sentinel-2A影像选择比较纯净的端元,利用其光谱信息优化从Landsat8影像中选择的端元的光谱曲线,提高纯净像元光谱曲线精度。此外,结合解混结果优化方案,利用归一化植被指数(normalized differential vegetation index,NDVI)和干旱裸土指数(dry bare-soil index,DBSI),对解混结果进行优化。采用WorldView-2遥感影像进行样本验证,结果显示,该方法所提取不透水面盖度的精度比传统方法提高了20%,为端元选取和不透水面提取提供可靠的理论支持。
关键词:
The extraction of impervious surface (IS) is very important for the development of cities, and linear spectral mixture analysis is commonly adopted to calculate the fraction of IS in the mixed pixel to improve the extraction of the urban IS at the subpixel scale. Owing to errors in the spectra of pure pixels selected from remote sensing images, incorrect fractions of different land cover types often emerge after unmixing. In this paper, the modified endmember selection was proposed to improve the accuracy of the spectral information of endmembers. Sentinel-2A images were applied to selected endmembers to get the spectral, which was used to modify the spectral information of the endmembers from Landsat8. In addition, the optimization scheme of LSMA results in which the normalized differential vegetation index (NDVI) and dry bare-soil index (DBSI) thresholds are used to optimize the mixing results was applied to improve the accuracy of LSMA results. With the WorldView-2 remote sensing image for sample verification, the results showed that the accuracy of IS fraction extracted by the method in this paper was 20% higher than that of the traditional method, providing reliable theoretical support for endmember selection and IS extraction.
Keywords:
本文引用格式
赵怡, 许剑辉, 钟凯文, 王云鹏, 胡泓达, 吴萍昊.
ZHAO Yi, XU Jianhui, ZHONG Kaiwen, WANG Yunpeng, HU Hongda, WU Pinghao.
0 引言
在城市不透水面遥感提取中,线性光谱混合分解(linear spectral mixture analysis, LSMA)是常用的方法[14,15],该方法基于亚像元尺度和植被-不透水面-土壤(vegetation-impervious-surface, V-I-S)模型[16],假设遥感影像中混合像元的反射率是该像元中所有地物端元反射率的线性组合[17]。首先从混合像元内各类端元的纯净像元获取端元的光谱信息,然后通过光谱解混计算,可以得到混合像元中高反射率和低反射率不透水面的面积比例,得到不透水面的面积比例,即不透水面盖度。由于环境影响,以及在特定波段上不透水面的光谱特征与其他地物相似,导致不透水面的提取精度较低[18]。在LSMA解混过程中,端元选择是至关重要的一步[19],原因在于LSMA是基于某一端元光谱特性,判断混像元内该类端元的面积比例,因此,所选取端元的光谱信息是进行解混混合像元的基础。
端元的纯净像元光谱信息来源有两个渠道: 一是从标准地物光谱库中直接选取目标地物的光谱曲线,二是从实际目标影像中选择纯净像元,得到不同地物的光谱曲线。由于地物的光谱特征曲线会受到环境因素的影响,且不同地物之间也会产生相互作用,因此,标准地物光谱库中的光谱曲线并不能较好地表示目标地物,故而大部分研究选择第二种方式获取端元的光谱曲线[20]。传统的LSMA在端元中采用最小噪声分离变换(minimum noise fraction rotation,MNF)[21]将遥感影像的信息集中到前几个波段,达到去除噪声的目的; 然后通过计算纯净像元指数(purity pixel index, PPI)[22]筛选出纯净度较高的像元; 最后通过N维可视化选择纯净像元。此外,N-FINDAR[23]算法也常用于纯净端元的选择。
为了描述城市进程发展,长时间序列的遥感影像必不可少。在目前的遥感发展中,能满足长时间序列要求的影像大都是中等空间分辨率的。然而在中等空间分辨率的遥感影像(例如30 m空间分辨率的Landsat8 OLI影像)中通常存在大量的混合像元,纯净像元选择存在不确定性,其精度直接影响不透水面的提取精度。而在较高空间分辨率的遥感影像中(例如10 m空间分辨率Sentinel-2A MSI影像),可以选择更多的纯净像元。在同一时期,目标地物在特定波段上的反射率是一定的,这就意味着,高空间分辨率中端元的光谱特性与中等空间空间分辨率中相似。为了获取更为可靠的纯净像元,本研究提出结合Sentinel-2A MSI影像改进从Landsat8 OLI 选择端元的光谱曲线。
此外,由于地物的光谱曲线受外界影响,且端元选择存在误差,因此,经过LSMA解混的不透水面盖度仍然存在错分现象,即在土壤分量中可能含有少量的不透水面,而在不透水面分量中也存在土壤或植被盖度。为了解决以上问题,本研究提出解混结果优化方案,即在像元尺度上,基于阈值分割法,利用归一化植被指数(normalized differential vegetation index,NDVI)和干旱裸土指数(dry bare-soil index,DBSI)识别不透水面、植被和土壤,从而分离出错分的像元,提高不透水面盖度提取精度。
1 研究区概况与数据源
1.1 研究区概况
图1
图1
研究区Landsat8 B4(R),B3(G),B2(B)合成影像
Fig.1
Study area image combined with Landsat8 B4(R),B3(G),B2(B)
1.2 数据源
Sentinel-2A 卫星发射于2015年6月23日,是“全球环境与安全监测”计划的第二颗卫星。该卫星的重访周期为十天,幅宽为290 km, 携带的多光谱成像仪(multispectral instrument,MSI)有13个波段,从可见光和近红外到短波红外,具有不同的空间分辨率[24,25]。本研究选取2017年11月1日的Sentinel-2A MSI遥感影像的B2,B3,B4,B8,B11,B12作为研究目标,采用Sen2cor处理器进行大气校正。此外,通过哨兵系列卫星产品综合应用平台(Sentinel Application Platform,SNAP),将B11,B12的遥感影像进行重采样,使其空间分辨率为10 m。
Landsat 8卫星发射于2013年2月11日,是美国陆地卫星计划的第8颗卫星。该卫星的重访周期为16 d,携带的陆地成像仪(operational land imager, OLI)有9个波段,其中包括7个常用波段(B1—B7,空间分辨率为30 m),一个全色波段(B8,空间分辨率为15 m)和一个云检测波段(B9,空间分辨率为30 m)。本研究选取2017年10月23日的Landsat8 OLI遥感影像的B1—B7作为研究目标,通过辐射定标和大气校正,将影像的DN值转换为地表反射率。地表发射率产品是在美国地质调查局(United States Geological Survey,USGS)地球资源观测(Earth Resources Observation and Science,EROS)和科学中心科学处理机构(Center Science Processing Architecture,ESPA)(
本研究采用的Sentinel-2A MSI 和Landsat8 OLI数据参数表1所示。
表1 Landsat8 OLI 和Sentinel-2A MSI 数据介绍
Tab.1
Landsat8 OLI | Sentinel-2A MSI | |||||
---|---|---|---|---|---|---|
波段 | 空间分 辨率/m | 波长/μm | 波段 | 空间分 辨率/m | 波长/μm | |
B1 | 30 | 0.433~0.453 | ||||
B2 | 30 | 0.450~0.515 | B2 | 10 | 0.458~0.523 | |
B3 | 30 | 0.525~0.600 | B3 | 10 | 0.543~0.578 | |
B4 | 30 | 0.630~0.680 | B4 | 10 | 0.650~0.680 | |
B5 | 30 | 0.845~0.885 | B8 | 10 | 0.785~0.900 | |
B6 | 30 | 1.560~1.660 | B11 | 20 | 0.855~0.875 | |
B7 | 30 | 2.100~2.300 | B12 | 20 | 2.100~2.280 |
此外,采用2017年WorldView-2遥感影像(空间分辨率为1.8 m),通过对样本区域内的不透水面进行矢量化,计算不透水面面积,作为精度验证数据。
2 研究方法
本研究所提出的方法分为3个步骤: 首先,从Sentinel-2A和Landsat8 OLI 遥感影像中选取纯净像元,获取线性光谱解混所需的端元光谱曲线; 然后利用传统线性光谱解混方法对两种数据进行光谱解混; 最后,利用NDVI和DBSI两种光谱指数,对解混结果进行后处理,得到精度较高的城市不透水面比例。本文整体流程如图2所示。
图2
2.1 端元选择优化方案
端元选择是LSMA过程中的重要步骤。理论上,作为光谱解混的基础,所选取的端元应是遥感影像上的纯净像元(只包括植被,裸土,高反射率不透水面与低反射率不透水面4类地物中的一类),然后利用所选纯净像元的光谱信息,对混合像元进行光谱解混,计算得到混合像元内每一种地物的面积比例。因此,实验中所选端元的纯净度直接影响解混的精度。Sentinel-2A MSI影像的分辨率较高,像元大小为10 m×10 m,与Landsat8 OLI影像(像元尺寸为30 m×30 m)相比,其像元的纯净度较高,纯净像元的选取也比较容易。此外,在城市环境中,土壤分布比较破碎,很少有大面积的裸土呈现,因此,像元面积越小,裸土纯净像元就越多,也就是说在Sentinel-2A影像上裸土的纯净像元更多。
为了提高端元光谱的选去精度,本研究提出端元选取优化方案: 结合Sentinel-2A MSI影像优化从Landsat8选取的端元光谱信息。具体过程如下:
1)采用传统方法对Landsat8 OLI影像进行端元选取,将所选纯净像元的光谱信息以文件的形式输出。
2)利用选择感兴趣区域(region of interest,ROI)的方式对Sentinel-2A影像进行端元选取,选择植被、裸土、高反射率不透水面与低反射率不透水面4类地物的纯净像元,并统计其光谱信息,生成端元光谱曲线,并以文件的形式输出。
3)利用从Sentinel-2A MSI 影像中选取的的纯净像元在每个波段(B2,B3,B4,B8,B11,B12)上的反射率替换Landsat8 OLI影像中选取的纯净像元在相应波段(B2—B7)的反射率,生成一个新的端元光谱信息文件,作为线性光谱混合分解的基础(从表1可以看出,Sentinel-2A MSI影像的6个波段与Landsat8 OLI的前6个波段的波谱范围比较接近,可以视为相同的波段)。
2.2 线性光谱混合分解(LSMA)
式中: i = 1, 2…, M, M 为波段数; n为端元数目; Ri为混合像元在i 波段的反射率; fk为端元k的面积比例; ERi为i波段的计算残差; RMS为计算残差,以此评价解混结果的优劣[29]。
由于混合像元中多类端元的光谱之间会发生相互作用,此外,环境等因素也会影响端元的光谱信息。因此在解混过程中,会存在不透水面被错误估计的问题,例如在不透水面覆盖度低于20%的区域,采用LSMA解混混合像元时,不透水面比例会被严重高估[3]。
由于不透水面的低反射率部分与水体部分反射率都比较低,十分接近,为了减少解混误差,在数据预处理阶段,利用改进的归一化水体指数(modified normalized difference water index, MNDWI)[30]将水体去除,计算公式为:
式中: BGREEN为影像的绿色波段(Landsat8 B3),BSWIR1为影像的中红外波段(Landsat8 B6)。其中,MNDWI的阈值选取采用OTSU阈值分割算法[31]选取。
2.3 解混结果优化方案
在LSMA结果中,解混误差普遍存在,具体表现为: ①在土壤盖度中含有不透水面,即不透水面被错分为土壤; ②在不透水面盖度中(高反射率与低反射率盖度之和)含有植被和土壤,即少量植被和土壤被错分为不透水面。为了降低以上误差,提高城市不透水面的提取精度,本研究结合像元尺度,利用DBSI和NDVI,通过二者的阈值,提取解混过程中被错分的像元,将其归还到相应的盖度中去。其中,DBSI和NDVI的计算公式分别为:
式中: BNIR为影像的近红外波段(Landsat8 B5); BR为影像的红色波段(Landsat8 B4)。
优化过程分为3步,具体如下:
1)从土壤盖度分量中提取错分的不透水面盖度。利用DBSI阈值筛选土壤盖度中的不透水面像元,借助DBSI影像的直方图,通过实验确定0.1作为该步DBSI的阈值,即在土壤盖度中,DBSI值小于0.1的像元被认为是不透水面,将这部分像元的解混值认为是不透水面盖度,对应图2中的不透水面盖度(1);
2)将LSMA解混结果中的高反照率和低反照率不透水面盖度相加,得到LSMA解混遥感影像得到的不透水面盖度,然后将上一操作中的不透水面盖度(1)与之相加,得到不透水面盖度(2);
3)去除上一操作中的不透水面盖度(2)中含有的少量植被和土壤盖度,得到最终的不透水面盖度。这部分误差来源于LSMA本身的解混误差。理论上,NVDI值越接近1,说明目标像元越接近植被,负值代表水体,0代表土壤或者沙漠化土地[32]。为了保证识别的像元为纯净的植被和土壤,借助DBSI和NDVI影像的直方图, 确定DBSI和NDVI阈值分别为0.2和0.3,即DBSI<0.2且NDVI>0.3的像元被认为是植被,DBSI>0.2且NDVI<0.3的像元被认为是土壤。
在以往研究中,DBSI>0.26的像元被认为是裸土,NDVI值接近于1被认为是植被。在步骤1)中,从土壤盖度中提取出不透水面,为了尽量保留LSMA结果的土壤盖度中土壤的真实面积比例,DBSI值应稍小于理论值,因此DBSI值选择为0.1。而在步骤3)中,从不透水面中去除植被和土壤,为了尽量保留真实不透水面的面积比例,DBSI和NDVI阈值应稍大于理论值,因此NDVI阈值被确定为0.3; 然而在实际实验过程中,经过反复试验得知,DBSI阈值为0.2能较好地区分不透水面与其他地物。
2.4 精度验证
为了验证本文所提出的方法对提取城市不透水面的有效性,随机选取了171个样本区域,样本大小为480 m×480 m,并且采用WorldView-2 影像对样本区域进行矢量化,得到不透水面的样本真实面积比例。通过系统误差(system error,SE)、平均绝对误差(mean absolute error,MAE)、均方根误差(root mean square error,RMSE)和决定系数R2来衡量实验结果的精度。前3个指标的计算公式分别为:
式中: Xi为不透水面盖度的实验结果;
3 结果与分析
3.1 端元选取
图3
3.2 结果比较分析
图4
图5
表2 不同方法提取的城市不透水面盖度精度分析
Tab.2
研究方法 | SE | MAE | RMSE | R2 |
---|---|---|---|---|
传统LSMA | 0.195 | 0.240 | 0.298 | 0.735 |
结合端元优化方案的LSMA | 0.190 | 0.217 | 0.265 | 0.851 |
结合端元优化和解混结果优化方案的LSMA | -0.066 | 0.103 | 0.133 | 0.879 |
综上所述,本研究所提出的端元优化方案可以使LSMA更加准确地识别不透水面,降低了不透水面盖度被低估的问题(在不透水面覆盖度较高的区域),提高城市化程度较高的区域中不透水面盖度的提取精度。
但是图5(a)和(b)同时存在一个问题,即在不透水面实际盖度较低的区域(横轴接近于0的位置),观测值明显比真值大。也就是说,在不透水面覆盖度较低的区域,传统LSMA和结合端元光谱优化方案的LSMA会高估不透水面盖度。图5(c)显示了结合端元优化方案和解混结果优化方案的LSMA提取的水面盖度与真值的线性拟合结果,该图表明,在不透水面盖度较低的区域,观测值明显降低。如表2所示,结合端元优化方案和解混结果优化方案的LSMA所提取的不透水面盖度的R2比传统LSMA高20%,比仅结合端元优化方案的LSMA高3%,且SE为-0.066,说明整体系统误差较小,该方法所得到的不透水面盖度与实际情况接近。以上结果说明通过解混结果优化方案,大大改善了不透水面被高估的问题(在不透水面盖度较低的区域)。
图6
4 结论
1)本研究针对Landsat8 OLI遥感影像的端元选择问题提出了端元选择优化方案; 针对LSMA后处理问题,结合解混结果优化方案,以此提高不透水面的提取精度。
2)与Landsat系列影像相比,Sentinel-2A MSI影像的空间分辨率较高,可以快速选择相对比较纯净的像元作为LSMA的端元,然后将选取像元的光谱信息代替对应在波段上从Landsat8 OLI影像上选择的端元的光谱信息,优化了光谱解混所必需的端元光谱信息,从而达到了提高LSMA对混合像元解混精度的目的。
3)对于LSMA后的各地物盖度分量,结合解混结果优化方案,利用DBSI和NDVI两个指数作为识别植被和土壤的重要参数,可以比较准确地从土壤盖度分量中提取被错分的不透水面盖度分量,且剔除了不透水面盖度分量中被错分的植被和土壤盖度分量,整体提高了不透水面盖度的提取精度。
以上两种优化方案解决了传统LSMA在不透水面盖度较高的区域(例如繁华的商业区,聚集的居民地等区域)被低估的问题,以及在不透水面盖度较低的区域(例如大面积的森林区域等)被高估的问题。然而,本研究在DBSI和NDVI阈值选择过程中,大部分操作依赖统计实验,未来可以从两个指数的阈值入手,为不透水面盖度提取提供可靠的理论方法,从而达到提高不透水面盖度提取精度的目的。
参考文献
Global distribution and density of constructed impervious surfaces
[J].DOI:10.3390/s7091962 URL [本文引用: 1]
Estimating impervious surfaces by linear spectral mixture analysis under semi-constrained condition
[C]//
Improving urban impervious surface mapping by linear spectral mixture analysis and using spectral indices
[J].DOI:10.1080/07038992.2015.1112730 URL [本文引用: 1]
Measuring spatio-temporal dynamics of impervious surface in Guangzhou,China,from 1988 to 2015,using time-series Landsat imagery
[J].DOI:10.1016/j.scitotenv.2018.01.155 URL [本文引用: 1]
Methods to extract impervious surface areas from satellite images
[J].DOI:10.1080/17538947.2013.866173 URL [本文引用: 1]
LSMA结合NDBI提取广州市部分城区不透水面的方法
[J].
Extraction of urban impervious surface in Guangzhou by LSMA with NDBI
[J].
Analysis of impervious surface and its impact on urban heat environment using the normalized difference impervious surface index (NDISI)
[J].
Large-area assessment of impervious surface based on integrated analysis of single-date Landsat-7 images and geospatial vector data
[J].DOI:10.1016/j.rse.2009.03.012 URL [本文引用: 1]
Thermal sharpening of land surface temperature maps based on the impervious surface index with the TsHARP method to ASTER satellite data:A case study from the metropolitan Kuala Lumpur,Malaysia
[J].DOI:10.1016/j.measurement.2018.04.092 URL [本文引用: 1]
Impervious surface area extraction from IKONOS imagery using an object-based fuzzy method
[J].DOI:10.1080/10106049.2010.535616 URL [本文引用: 1]
基于空间自回归模型的广州市NDVI和NDBI与气温关系研究
[J].
Relationship of air temperature to NDVI and NDBI in Guangzhou City using spatial autoregressive model
[J].
Improving impervious surface estimation:an integrated method of classification and regression trees (CART) and linear spectral mixture analysis (LSMA) based on error analysis
[J].
Estimating impervious surface distribution by spectral mixture analysis
[J].DOI:10.1016/S0034-4257(02)00136-0 URL [本文引用: 1]
Linear spatial spectral mixture model
[J].
Exploring a V-I-S (vegetation-impervious surface-soil) model for urban ecosystem analysis through remote sensing:Comparative anatomy for cities?
[J].DOI:10.1080/01431169508954549 URL [本文引用: 1]
A sub-pixel analysis of urbanization effect on land surface temperature and its interplay with impervious surface and vegetation coverage in Indianapolis,United States
[J].
联合丰度信息与空谱特征的高光谱影像分类研究
[D].
Research on hypercpectral imagery classification by combing abundance information and spectral-spatial feature
[J].
Fractional vegetation cover estimation based on an improved selective endmember spectral mixture model
[J].DOI:10.1371/journal.pone.0124608 URL [本文引用: 1]
基于线性混合模型的端元提取方法综述
[J].
A review on endmember extraction algorithms based on the linear mixing model
[J].
Enhancement of high spectral resolution remote-sensing data by a noise-adjusted principal components transform
[J].
The application and evaluation of two methods based on LSMM model:A case study in Guangzhou
[J].
一种改进的快速N-FINDR端元提取算法
[J].
An improved fast N-FINDR endmember extraction algorithm
[J].
Sentinel-2:ESA’s optical high-resolution mission for GMES operational services
[J].DOI:10.1016/j.rse.2011.11.026 URL [本文引用: 1]
Characterization of Sentinel-2A and Landsat-8 top of atmosphere,surface,and nadir BRDF adjusted reflectance and NDVI differences
[J].DOI:10.1016/j.rse.2018.04.031 URL [本文引用: 1]
Coupling modified linear spectral mixture analysis and soil conservation service curve number (SCS-CN) models to simulate surface runoff:Application to the main urban area of Guangzhou,China
[J].DOI:10.3390/w8120550 URL [本文引用: 1]
地表不透水面信息遥感的主要方法分析
[J].
Remote sensing-based retrieval of ground impervious surfaces
[J].
Normalized spectral mixture analysis for monitoring urban composition using ETM+ imagery
[J].DOI:10.1016/j.rse.2004.08.003 URL [本文引用: 1]
A spectral mixture analysis and mapping of impervious surfaces in built-up land of Fuzhou City
[J].
A study on information extraction of water body with the modified normalized difference water index (MNDWI)
[J].
A threshold selection method from gray-level histograms
[J].
Data on identification of desertified regions in Anantapur district,Southern India by NDVI approach using remote sensing and GIS
[J].DOI:10.1016/j.dib.2020.105560 URL [本文引用: 1]
/
〈 |
|
〉 |
