基于Sentinel-2數據的不透水面識別方法與應用研究

徐 瑞，屈 冉，劉 曉 龍，馬 慧 榕，史 正 濤

(1.云南師范大學地理學部，云南 昆明 650500；2.生態(tài)環(huán)境部衛(wèi)星環(huán)境應用中心，北京 100094；3.內蒙古工業(yè)大學信息工程學院，內蒙古 呼和浩特 010051)

0 引言

不透水面[1](如瀝青道路、房屋、廣場、停車場、塑料大棚等)在擴張過程中占用了大量耕地、林地及水域，不僅導致生物多樣性銳減，影響人與自然和諧相處[2]，還影響城市水循環(huán)系統(tǒng)，改變地表徑流，增加城市內澇災害風險[3,4]，且加劇了城市“熱島”形成，導致城市生態(tài)系統(tǒng)惡化[5-7]。因此，快速、準確獲取不透水面的空間分布信息對于制定生態(tài)環(huán)境保護政策、合理規(guī)劃土地資源、構建新型生態(tài)城市以及海綿城市等具有重要意義。

遙感技術時效性高、覆蓋范圍廣，被廣泛應用于不透水面提取。面向對象法[8]、影像分類法[9]、光譜混合分析(Spectral Mixture Analysis，SMA)[10]、分類回歸樹(Classification and Regression Trees，CART)[11]、人工神經網絡(Artificial Neural Network，ANN)[12,13]、深度學習[14]及指數法[15,16]常用于不透水面識別。面向對象法基于訓練樣本對影像進行分割分類，從而提取不透水面，主要適用于高分辨率遙感影像[8];影像分類法(如最大似然法、平行管道法、支持向量機等)依據不透水面光譜信息利用傳統(tǒng)分類器進行自動分類，操作性強[17]，主要針對中低空間分辨率光學影像(如MODIS[18,19]、Landsat[20,21])，但存在大量混合像元，降低了提取精度;SMA將影像像元分為高、低反射率對象、土壤和植被4種端元，通過生成各端元的灰度圖提取不透水面，能解決混合像元導致分類精度低的問題[10]，適用于大規(guī)模城市(如京津冀、長三角等城市群)的不透水面提取[22-24];ANN適用于識別范圍較小的城市不透水面，且分類精度優(yōu)于傳統(tǒng)的統(tǒng)計分類方法[12,13]，但存在大量隱含層數和參數，求取最佳學習速率和動量系數耗時長[14];深度學習能輔助ANN處理參數，提取高層次的特征數據，常被應用于不透水面提取[25,26];CART算法通過訓練樣本形成的規(guī)則預測連續(xù)變量提取不透水面[11]。上述方法均建立在模型學習的基礎上，對訓練樣本依賴性較強，一直面臨模型訓練成本與精度權衡這一主要問題[27];而指數法利用不透水面的光譜響應特征，通過構建遙感光譜指數提取不透水面信息[28]，無需訓練樣本，操作簡便、計算效率高[15,16,28]。常用的遙感指數有歸一化建筑指數(Normalized Difference Built-up Index，NDBI)[29]、遙感建筑用地指數(Index-based Built-up Index，IBI)[30]、歸一化差值不透水面指數(Normalized Difference Impervious Surface Index，NDISI)[31]、生物物理綜合指數(Biophysical Composition Index， BCI)[32]、比值居民地指數(Ratio Resident-area Index，RRI)[33]、不透水面綜合指數(Combinational Build-up Index，CBI)[34]、改進型歸一化差值不透水表面指數(Modified Normalized Difference Impervious Surface Index，MNDISI)[15]等。

Attarchi等[35-38]學者嘗試應用合成孔徑雷達(SAR)和夜間燈光等遙感影像提取不透水面。其中，SAR數據具有穿云透霧、不受天氣影響等優(yōu)點，但在復雜地形區(qū)提取不透水面仍面臨挑戰(zhàn)[39]；夜間燈光傳感器(如DMSP-OLS、VIIRS)基于其獨特的光學放大能力[40,41]，可以根據燈光亮度閾值提取不透水面，但在偏遠地區(qū)和欠發(fā)達地區(qū)仍存在不足[42]。與常用的陸地衛(wèi)星相比，Sentinel-2空間分辨率更高(Red、Green、Blue、NIR波段空間分辨率為10 m，SWIR波段空間分辨率為20 m)，重訪周期短(單星重訪周期為10 d，雙星重訪周期為5 d)。因此，本文基于Sentinel-2數據在不透水面提取方面的良好性能以及指數法的優(yōu)點[43,44]，提出一種結合多種指數和主成分分析的不透水面識別方法(Impervious Surface Identification Method，ISIM-S2)，以滿足快速、準確獲取局部或大范圍不透水面分布信息的需求。

1 研究區(qū)和數據

1.1 研究區(qū)概況

不透水面的分布與擴張受地形地貌、氣候(降水)的影響，與地理環(huán)境特征密切相關。為驗證ISIM-S2方法的普適性，選取中國內陸不同氣候區(qū)、地貌類型、城市環(huán)境的4個典型城市(北京、長春、廣州和昆明)作為研究區(qū)。北京位于華北平原北部(115°42′～117°24′E，39°24′～41°36′N)，地勢西北高、東南低，屬于暖溫帶半濕潤半干旱季風氣候，冬季嚴寒干燥，夏季高溫多雨；長春地處東北平原地區(qū)(124°18′～127°5′E，43°5′～45°15′N)，屬于溫帶大陸性半濕潤季風氣候，雨熱同期，四季分明，地勢起伏小(地表相對高差不超過50 m)，有利于城市基礎設施建設；廣州地處南部沿海地區(qū)(112°57′～114°3′E，22°26′～23°56′N)，中部以中低山和丘陵為主，東部以山地為主，南部以平原為主，屬于海洋性亞熱帶季風氣候，冬暖夏涼，雷雨較多；昆明地處云貴高原中部(102°10′～103°40′E，24°23′～26°22′N)，總體地勢北高南低，屬于亞熱帶高原山地季風氣候，全年溫差較小，干、濕季分明。

1.2 數據及預處理

1.2.1 實驗數據 1)Sentinel-2影像源自哥白尼開放訪問中心(https://scihub.copernicus.eu/)，選用無云或云量較少的影像，北京、長春、廣州影像獲取衛(wèi)星為Sentinel-2A，獲取時間分別為2017年7月7日、2017年6月28日、2017年12月31日，昆明影像獲取衛(wèi)星為Sentinel-2B，獲取時間為2018年1月7日。本文使用歐空局發(fā)布的Sen2Cor處理模塊(http://step.esa.int/main/third-party-plugins-2/sen2cor/)對Sentinel-2 1C級產品進行輻射定標和大氣校正，并采用最近鄰法將20 m空間分辨率的短波紅外(SWIR)波段數據重采樣至10 m。2)高程數據(空間分辨率為90 m)和氣候分區(qū)數據分別來源于https://srtm.csi.cgiar.org/srtmdata/和K?PPEN-GEIGER氣候分類世界地圖(http://koeppen-geiger.vu-wien.ac.at/)。

1.2.2 驗證數據 為評價不透水面提取結果的準確性，本文按隨機分層抽樣原則選擇驗證樣本。盡管當前有不透水面分類產品(如GAIA[45]、GLASS-GLC[46])和高分辨率影像(如Google Earth中的遙感影像)可作為驗證樣本，但裸土和城市不透水面的快速變化使上述產品的時效性不足，不能直接作為本文提取結果的驗證參考。為此，本文通過對Sentinel-2影像進行目視解譯選取驗證樣本，然后根據成像日期接近的Google Earth影像對驗證樣本進行調整，以確保不透水面驗證樣本的準確性。最終獲取北京、長春、廣州和昆明4個研究區(qū)的驗證樣本數分別為217個、198個、213個和214個。

2 不透水面提取方法

ISIM-S2方法主要包括不透水面靶區(qū)(含裸土和不透水面)提取和靶區(qū)內不透水面識別兩部分，具體流程如圖1所示。

圖1 ISIM-S2方法流程Fig.1 Flow chart of ISIM-S2 method

2.1 不透水面靶區(qū)提取

遙感指數可量化地表生物的物理特征，增強地類間的差異性，但由于不同地表覆蓋類型(如裸土和不透水面)存在混淆，導致單一遙感指數的有效性受到限制[47]。為最大限度地提高遙感指數性能，本文采用多個遙感指數逐步分離不同土地覆蓋類型。1)水體識別。改進的歸一化水體指數(Modified Normalized Difference Water Index，MNDWI)[48]不僅能有效提取城市區(qū)域內的水體，還能區(qū)分水體中的陰影以及水體中的懸浮沉積物等細微特征，因此，本文選擇MNDWI指數提取水體(式(1))，即將MNDWI值大于0的像元識別為水體。2)云/云陰影識別。采用Fmask算法[49]識別Sentinel-2 1C級影像中的云/云陰影，該算法下載自GitHub官網(https://github.com/GERSL/Fmask)，算法參數均為默認值。將上述已識別的水體、云/云陰影進行掩膜后，研究區(qū)內剩余植被、裸土和不透水面3種地類。3)不透水面靶區(qū)提取。NDBI指數主要根據不同地類的光譜特征差異增強不透水面信息[29]，雖然低覆蓋植被、裸土及不透水面之間存在異物同譜現象，導致NDBI指數將部分裸土與不透水面和低覆蓋植被混淆，但NDBI指數具有易操性[47]；而SAVI[50]能抑制裸土和背景噪聲，與NDVI相比，更適用于分離植被和裸土。因此，本文采用NDBI-SAVI[47]提取裸土和不透水面(式(2)、式(3))。

(1)

(2)

(3)

式中：ρGreen、ρSWIR、ρNIR、ρRed分別為綠光(Green)、短波紅外(SWIR)、近紅外(NIR)和紅光(Red)波段的地表反射率；L為土壤調節(jié)因子，本文取0.5。

NDBI-SAVI圖像的灰度值范圍為-2.5～2.5，其中不透水面靶區(qū)像元的灰度值大于植被像元的灰度值,為分離不透水面靶區(qū)和植被，應設置最佳閾值。OTSU[51]算法通過閾值將待處理圖像分為前景與背景，取最佳閾值時，分類效果較好。因此，本研究采用OTSU算法確定NDBI-SAVI灰度圖像的最佳閾值。

2.2 靶區(qū)內不透水面識別

主成分分析(PCA)是一種有效的特征提取工具，可在特征數量顯著減少的同時，保證每個主成分丟失的原始數據信息較少，且不同主成分包含的信息唯一。因此，本文利用PCA和OTSU閾值分割方法將裸土和不透水面分離。實驗表明：北京、長春采用PCA2，廣州、昆明采用PCA3，OTSU閾值均為0.01。PCA系數可能隨影像的變化而變化(如場景范圍、土地覆蓋類型、波段數、圖像采集日期等)，表明PCA衍生的特征不穩(wěn)定，主要是因為地表覆蓋物(如植被、水體)在某些特定波段或所有波段具有高反射率(或低反射率)特征，而ISIM-S2方法基于分層分類思想去除裸土和不透水面以外的地物，可降低PCA衍生特征的不穩(wěn)定性。

3 不透水面提取結果評價

本文采用BCI、CBI、RRI指數基于OTSU閾值分割法提取研究區(qū)不透水面并與ISIM-S2方法進行對比，采用總體精度(OA)、用戶精度(UA)、生產者精度(PA)和Kappa系數進行精度評估(圖2)，并對各研究區(qū)全局和局部區(qū)域不透水面空間分布進行比較分析(圖3)。

圖2 ISIM-S2與3種指數在不同研究區(qū)的分類精度對比Fig.2 Classification accuracies comparison of ISIM-S2 with three indexes in different study areas

由圖2可知，ISIM-S2方法在4個研究區(qū)的總體分類精度(OA)均大于90%，Kappa系數均大于0.80，且該方法的OA和Kappa系數均高于其他3種指數，提取不透水面的精度最高，顯著減少了裸土與不透水面的混淆。但ISIM-S2方法也存在如下不足：會將部分暗色不透水面(水泥、瀝青等)誤分為透水層(陰影、水體)；面積較小的不透水面邊界常與其他地物形成混合像元，從而易被誤分為透水層(圖3長春研究區(qū)中)；部分有色屋頂(以藍色為主)與植被在紅光波段光譜值差異較小[52]，易被誤分為透水層。

圖3 ISIM-S2與3種指數在不同研究區(qū)的不透水面提取結果對比Fig.3 Result comparison of impervious surface extraction using ISIM-S2 with three indexes in different study areas

由4個研究區(qū)內提取結果對比(圖3)可知：1)由于云與不透水面等部分高反射率地物的反射率相似，導致3種指數均會將云識別為不透水面，BCI和RRI指數的誤識別情況更嚴重，而ISIM-S2在提取靶區(qū)時先用Fmask算法將云/云陰影去除，降低了云對不透水面的影響，識別錯誤較少。2)3種指數均存在將水體誤識別為不透水面的情況，其中RRI指數的誤識別最嚴重，而ISIM-S2在靶區(qū)提取階段將水體去除，降低了誤分率。3)4種方法均存在將不透水面與裸土混淆的情況，其中，CBI將裸土誤識別為不透水面的概率最高，ISIM-S2針對靶區(qū)范圍的裸土和不透水面進行主成分分析，混淆程度最低。

綜上，相比BCI、CBI、RRI指數，ISIM-S2方法降低了云、水體對不透水面的影響，有效改善了裸土和不透水面提取混淆情況，且在4個研究區(qū)內總體識別精度均屬最高。

4 結論

不透水面與裸土的混淆一直是基于遙感數據識別不透水面的主要難題，為降低裸土對不透水面提取的影響，本文提出一種基于Sentinel-2數據的不透水面識別方法(ISIM-S2)。首先利用Fmask、MNDWI、NDBI和SAVI 指數剔除云/云陰影、水體和植被，得到不透水面靶區(qū)；然后采用主成分分析法對靶區(qū)內的裸土和不透水面進行分離，得到不透水面提取結果。將該方法應用于4個不同地理環(huán)境下的研究區(qū)(北京、長春、廣州和昆明)，并與BCI、CBI、RRI指數的提取結果進行對比，得出以下結論：1)ISIM-S2方法在4個研究區(qū)內的不透水面總體分類精度(OA)均大于90%，Kappa系數均大于0.80，不透水面提取精度較高。2)相比BCI、CBI、RRI指數，ISIM-S2預先去除云/云陰影和水體，有效抑制了高(低)反射率地物像元值變化導致PCA衍生特征的不穩(wěn)定性，不僅降低了云、水體對不透水面的影響，還能有效改善裸土與不透水面提取混淆問題，在4個研究區(qū)提取的總體精度(OA)和Kappa系數均高于BCI、CBI、RRI指數。

本文提出的ISIM-S2方法無需訓練樣本，操作簡單、高效且精度較高，可為不透水面識別提供一種新途徑。但由于不同地理環(huán)境下城市分布、不透水面像元占比和裸土反射特性存在差異，如何實現靶區(qū)內PCA裸土和不透水面區(qū)分閾值的自適應設置，后期需要深入研究。