半城市化地區不透水面遙感提取優化方法研究

崔 雯 卓，林 晨，陳 曦，熊 俊 峰，侯 渲，馮 雪 嬌，吳 紫 靜,5

(1.中國科學院流域地理學重點實驗室，中國科學院南京地理與湖泊研究所，江蘇 南京 210008；2.中國科學院大學，北京 100049；3.中山大學地理科學與規劃學院，廣東 廣州 510275；4.自然資源部海岸帶開發與保護重點實驗室，江蘇 南京 210023；5.南京工業大學測繪科學與技術學院，江蘇 南京 211816)

0 引言

中國已進入快速城市化階段，大量具有城市化基礎的半城市化地區逐漸轉為城市地區[1-4]。不透水面(Impervious Surface Area,ISA)是典型的人工地表，目前已成為研究城市化進程及其發展強度最直觀的指標[5，6]。與城市化地區相似，半城市化地區主要由不透水面覆蓋[7-9]，但其建筑密度和功能與城市化地區存在明顯差距。半城市化地區作為城市發展的后備用地，肩負著承接城市職能和支持城市發展的重任[10]，因此，快速獲取半城市化地區的不透水面分布信息對城市規劃和管理具有重要意義。

遙感技術因具有探測范圍大、獲取資料速度快、實時性強等特點，被廣泛應用于不透水面提取[11]。目前，不透水面提取方法分為兩大類：1)基于輸入圖像的分類方法[12]，包含基于像素或面向對象的分類方法[13-16]，如光譜混合分析法(Spectral Mixture Analysis,SMA[17])、影像分類法[18,19]、基于人工神經網絡(Artificial Neural Network,ANN)[20]的分類方法等；2)通過閾值直接分割圖像的方法(指數法)，即利用不透水面的光譜影像特征，通過構建遙感光譜指數提取不透水面信息[11,21]，常采用影像的原始波段構建不透水面指數，通過對不透水面與其他地類之間有明顯區別的反射波段進行運算，以此抑制其他地類信息并增強不透水面信息，該方法無需預設訓練樣本，可操作性更強且自動化程度更高。目前，應用最廣泛的不透水面提取指數是歸一化建筑指數(Normalized Difference Build-up Index,NDBI)[22]，該指數初次應用于提取南京市不透水面且總體精度高達92.6%。由于城市地表具有復雜異質性，原始波段的組合無法精確地剔除不透水面以外的地物，一些學者采用其他指數組合或利用分量方法重構不透水面指數。例如：Sun等[23]提出組合建筑區指數(Combination Build-up Index,CBI)，通過組合主成分分析的第一分量和歸一化差異水體指數(Normalized Difference Water Index,NDWI)[24]提取不透水面信息；Xu[25]根據改進型歸一化差異水體指數(Modified Normalized Difference Water Index,MNDWI)、土壤調節植被指數(Soil-Adjusted Vegetation Index，SAVI)[26]和歸一化建筑指數(NDBI)提出遙感建筑用地指數(Index-based Built-up Index,IBI)，并在福州地區得到有效應用；Bhatti等[12]結合歸一化差異水體指數(NDWI)、歸一化植被指數(Normalized Difference Vegetation Index,NDVI)[27]和改進型歸一化差異水體指數(MNDWI)[28]構建了建筑指數(Built-up Area Extraction Method,BAEM)。

綜上，現有研究針對不透水面提取構建了較為成熟的指數，并在京津冀[29-31]、長江三角洲[32,33]等高度城市化地區得到廣泛應用，卻極少應用于半城市化地區。半城市化地區建筑密度較低，不透水面周圍存在大量耕地與裸土，不僅需要突出不透水面光譜信息，還需要抑制其他地物的光譜信息，而現有不透水面提取指數難以適用于半城市化地區，因此，本研究提出一種適用于半城市化地區的不透水面提取方法。

1 研究區概況與數據來源

1.1 研究區概況

目前我國城市正處于轉型期，江蘇省提出“1+3”功能區戰略助推高質量發展，江淮生態經濟區(圖1a)是全省經濟發展相對落后的區域，該區域城鎮化水平差異較大，城鎮體系結構不完整，城鎮綜合承載力低，亟須協同推進區域新型城鎮化建設，是城市轉型需要重視的區域。本文以區縣級行政單元作為半城市化地區，選取江淮生態經濟區的淮陰區和盱眙縣為研究區。1)盱眙縣位于江蘇省中西部(圖1b)，地處淮河中下游、洪澤湖南岸，作為丘陵區的代表，耕地多為水田，以坑塘養殖產業為主，區域經濟發展狀況位居江蘇省中等水平，根據第七次全國人口普查數據，盱眙縣城鎮化率達61.66%，近5年提高近10.46%，處于逐漸向城市轉型時期。2)淮陰區位于江蘇省北部平原的中心(圖1c)，耕地以旱地為主，處在江淮生態經濟區、淮河生態經濟帶、大運河文化帶戰略交匯處，又位于長三角一體化發展、長江經濟帶和南京都市圈的影響范圍內，目前經濟發展位居江蘇省中等偏高水平，根據第七次全國人口普查數據，淮陰區城鎮化率達57.61%，近5年提高了12.95%。為更形象地體現半城市化地區在遙感影像上的特點，選取A、B、C、D 4個區域(圖1d)進行展示，其中A和C地區除不透水面外還存在大量裸土，B和D地區建筑密度較低，周圍被大量耕地包圍。

圖1 研究區概況Fig.1 Overview of the study area

1.2 數據及預處理

研究區存在大量較窄道路，選擇空間分辨率較高的遙感數據可以盡量減少混分情況。本文根據研究區實際，選用云量小于5%、空間分辨率為10 m的Sentinel-2A數據，數據來源于哥白尼數據開放訪問中心(Copernicus Open Access Hub)，獲取時間均為植被茂密、坑塘未排水的2021年7月30日。研究采用的波段為可見光波段(Band 2、Band 3、Band 4，10 m空間分辨率)、近紅外波段(Band 8，10 m空間分辨率)和短波紅外波段(Band 11、Band 12，20 m空間分辨率)。利用歐空局(ESA)發布的Sen2Cor工具箱對Sentinel-2A 1C級產品進行輻射定標、大氣校正，最后將數據重采樣至10 m。

2 半城市化地區不透水面提取方法

傳統BAEM指數[12]內部的約束條件并不符合半城市化地區情況，應用于半城市化地區時，不透水面提取效果較差。因此，本文從以下兩方面進行優化：1)考慮半城市化地區存在大量裸土，需要對現有的指數增加土壤調節植被指數(SAVI)進行優化，以區分裸土的影響；2)傳統方法采用雙窗口變步長閾值搜尋算法(DFPS)[26]提取閾值，考慮到半城市化地區存在混合像元的問題，對樣本選取要求方面進行閾值法改進。

2.1 傳統BAEM指數優化

BAEM指數將地物劃分為水體、植被與不透水面，利用各地物最具代表性的指數(歸一化建筑指數(NDBI)、歸一化植被指數(NDVI)和改進的歸一化差異水體指數(MNDWI))進行差值運算(式(1)-式(4))，以有效剔除不透水面以外的地物。然而，將傳統BAEM指數方法應用于半城市化地區時，易導致不透水面提取結果存在大量裸土與不透水面混分。

BAEM=NDBI-NDVI-MNDWI

(1)

(2)

(3)

(4)

式中：ρNIR、ρRed、ρGreen、ρSWIR分別為近紅外(Band 8)、紅光(Band 4)、綠光(Band 3)、短波紅外波段(Band 11)的地表反射率。

本研究將地表覆蓋物分為不透水面、水體、植被和裸土4類典型地物，由其光譜曲線(圖2)可以看出，不透水面與裸土的光譜特征相似，從綠光波段(Band 3)到短波紅外波段(Band 11)呈遞增趨勢，到Band 12反射率降低，植被呈明顯的紅邊效應，水體整體反射率較低。因半城市化地區存在大量裸土，需要在傳統的BAEM指數中增加SAVI指數(式(5))進行優化，以區分裸土的影響，其中，土壤調節因子L取0.5，最終得到改良后的BAEM指數(MBAEM)(式(6))。

圖2 典型地物反射率光譜曲線Fig.2 Reflective spectrum curves of typical surface features

(5)

MBAEM=NDBI-SAVI-NDVI-MNDWI

(6)

為驗證加入SAVI能否區分裸土與其他典型地物，選取盱眙縣和淮陰區具有代表性的典型像元進行分析。在每個區域中選擇較密集的建成區a、被大量耕地包圍的房屋b、道路c、水體d、植被e及裸土f，分別展示各個像元的光譜曲線以及4種指數(SAVI、NDVI、NDBI、MNDWI) (圖3-圖5)。由圖3b和圖4b可知，a、b、c 3類不透水面具有相似的光譜曲線特征，但其光譜曲線與水體d、植被e差別較大，與裸土f光譜曲線趨勢一致，因此將其像元的4種指數進行對比。如圖5所示，a、b、c 3類不透水面的4種指數變化趨勢基本一致；水體d的MNDWI指數與其他指數區分明顯；植被e的NDVI指數與其他指數區分明顯；裸土f的光譜曲線與不透水面相似，但在指數方面存在一定差異，因此可以利用4種指數的組合抑制植被、水體與裸土的信息，突出不透水面信息。

圖3 盱眙縣典型地物的像元光譜曲線Fig.3 Pixel spectrum curves of typical surface features in Xuyi County

圖4 淮陰區典型地物的像元光譜曲線Fig.4 Pixel spectrum curves of typical surface features in Huaiyin District

圖5 盱眙縣和淮陰區典型地物像元的指數特征Fig.5 Pixel exponential characteristics of typical surface features in Xuyi County and Huaiyin District

2.2 閾值法改進

Chen等提出的雙窗口變步長閾值搜尋算法(DFPS)是在土地利用/覆被變化監測中提取變化像元的半自動閾值搜尋算法[34]。首先選取典型訓練樣區的樣本，其中選取的樣本最好為“島狀”，即典型變化區被非典型變化像元包圍，使典型區與非典型區存在明顯的差異性；其次設定搜尋范圍和步長，根據給定的步長在合理的搜索范圍內尋找獲得典型樣區變化檢測精度最高的閾值。搜索過程需要進行反復嘗試，每次循環對應的范圍和步長會依次縮小，保證搜索的變化閾值最優。該算法存在的問題是訓練樣本需要被其他地類訓練樣本包圍，樣本之間相互依賴，需要形成非獨立的雙窗口。但在半城市化地區，訓練樣本往往復雜交錯，仍存在混合像元問題，即使10 m空間分辨率也不能保證樣本完全獨立，因此本研究采用大津法(OTSU算法)[35]提取閾值，使圖像分割不受樣本影響，計算相對簡單。

OTSU算法[35]按圖像的灰度特性分為背景和前景兩部分，自動計算背景與前景的最大類間方差，從而將圖像進行二值化，其最大優勢在于不需要預設樣本，自動獲取最佳閾值，計算簡單。本研究以MBAEM指數灰度圖像為例對該算法進行展開介紹。理論上MBAEM圖像的灰度值范圍為-3～3，其中前景(不透水面)的灰度值大于背景(植被、水體以及裸土)的灰度值,利用OTSU算法分離不透水面和其他地類的計算公式見式(7)。

(7)

式中：σ2為前景和背景圖像的方差；μ1和μ0分別為前景和背景的平均值；μ為整幅影像像元平均值；w1和w0分別為前景和背景在整幅影像中的占比。

根據式(7)提取的閾值threshold，對MBAEM指數灰度圖像進行分割，得到半城市化地區不透水面影像(Impervious Surface,IS)(式(8))。本研究中MBAEM指數提取不透水面在盱眙縣閾值為-0.48，淮陰區閾值為-0.45，BAEM指數提取不透水面在盱眙縣和淮陰區的閾值分別為-0.28和-0.33。

(8)

式中：(i，j)表示MBAEM指數影像的像元行列號。

2.3 精度檢驗

為客觀驗證不透水面提取結果的準確性，本研究在Sentinel-2A原始影像上，采用目視解譯方法和分層隨機抽樣方法選取驗證樣本[36]，并在Google Earth上選取與成像時間接近的影像對驗證點進行比對和微調。最終，在淮陰區選取不透水面與非不透水面樣本數量分別為3 964個和4 272個，在盱眙縣選取不透水面與非不透水面數量分別為3 261個和3 796個。

3 不透水面提取結果與評價

3.1 優化方法與其他方法提取結果對比

3.1.1 優化方法與傳統方法提取結果對比 基于Sentinel-2A數據，采用BAEM指數和MBAEM指數提取盱眙縣和淮陰區的不透水面，分別用兩種閾值法提取閾值(圖6、圖7)。從空間分布上看，盱眙縣的不透水面主要集中在中部建成區，其他村鎮向四周呈現放射狀分布，村鎮之間有道路連接；淮陰區不透水面也主要位于中部建成區，周圍分布村鎮，道路縱橫交錯。從提取結果整體上看，BAEM方法得到的不透水面提取結果分布更多，因此，本研究進一步采用總體精度(OA)、Kappa系數、制圖精度 (PA)和用戶精度 (UA)進行精度檢驗(表1)，以客觀、具體展示4種方法的不透水面提取結果。由表1可知：盱眙縣和淮陰區利用MBAEM指數和OTSU算法得到的結果均為最佳，其總體精度達89.50%和89.73%，Kappa系數為0.7907和0.7950；相比BAEM指數，MBAEM指數提取的總體精度提高約10%，Kappa系數提高約15%。綜上，MBAEM指數提取結果比BAEM指數提取結果精度更高，OTSU算法的精度比DFPS算法的精度更高，因此本研究基于指數和閾值兩方面進行優化的精度均得到提升。

圖6 盱眙縣不透水面提取結果對比Fig.6 Comparison of extraction results of impervious surface in Xuyi County

圖7 淮陰區不透水面提取結果對比Fig.7 Comparison of extraction results of impervious surface in Huaiyin District

表1 不透水面精度驗證結果Table 1 Accuracy validation of impervious surface extraction results

3.1.2 MBAEM指數與其他常用指數提取結果對比 在研究區選取A、B、C、D 4個區域，對NDBI、IBI和CBI指數提取結果與MBAEM指數提取結果進行對比(圖8)，所有指數閾值均采用OTSU法確定。可以看出，MBAEM指數的提取結果更接近真實情況；IBI指數的提取結果分布相對較少且稀疏，存在大量漏提；CBI指數的提取結果更密集，但也有許多非不透水面被錯分為不透水面；在細節上，NDBI指數則漏分或錯分了一些半城市化地區中的典型地物(如裸地)。由盱眙縣和淮陰區的提取精度(表2)可知：MBAEM指數的總體精度、Kappa系數、制圖精度和用戶精度均最高。

圖8 MBAEM與3種指數不透水面提取結果對比Fig.8 Comparison of extraction results of impervious surfaces using MBAEM and three indexes

表2 MBAEM指數與其他3個指數提取結果精度指標對比Table 2 Accuracy comparison of impervious surface extraction results using MBAEM and three indexes

3.2 優化方法與傳統方法提取結果細節對比

通過對比發現，閾值上OTSU算法得到的結果優于DFPS算法的結果，因此選取OTSU算法的結果進行細節對比。將原始BAEM指數提取結果與優化后的MBAEM指數提取結果局部進行對比(圖9、圖10)，發現優化后的指數提取結果在水體和裸土方面效果更好。1)水體(圖9)。水體在真彩色影像呈藍色系線條，而瀝青混凝土等道路在真彩色影像上呈灰色條狀特征，因此二者呈現的特征相似。圖9A中，道路呈淺色線狀分布，河流呈深藍色，MBAEM指數可以準確提取道路，但BAEM指數將水體誤分為道路；圖9B與圖9C中，MBAEM指數均將水體與不透水面更好地區分開，而BAEM指數卻將水體誤分為不透水面。因此，MBAEM指數較好地解決了水體被誤分為不透水面的問題且能清晰提取出道路。2)裸土(圖10)。裸土在真彩色合成影像中呈現棕色系，大部分呈矩形塊狀分布，少量呈條狀零星分布，且顏色偏亮。傳統提取方法未考慮裸土與不透水面各波段均具有較高反射率的特點，因此未將二者區分開。如圖10A所示，MBAEM指數較好地區分開建設用地周圍未開發利用的裸土區域，而BAEM指數卻將裸土與不透水面混分；在圖10B和圖10C中，MBAEM指數能有效將耕地周圍未耕種的裸土與不透水面區分開。由此可見，通過增加土壤調節植被指數(SAVI)，能較準確地提取半城市化地區不透水面。

圖9 水體的局部對比Fig.9 Local contrast for water

圖10 裸土的局部對比Fig.10 Local contrast for bare soil

4 結論與展望

本文以盱眙縣和淮陰區作為半城市化地區的典型區，利用Sentinel-2A遙感數據，在原有BAEM指數的基礎上優化指數選擇和閾值方法，提出一種新的不透水面綜合指數MBAEM，該指數綜合考慮了不透水面與水體、裸土和植被之間的差異，并采用OTSU閾值法提取半城市化地區的不透水面。結論如下：1)利用MBAEM指數和OTSU閾值法提取效果最好，盱眙縣和淮陰區的總體精度分別為89.50%和89.73%，Kappa系數分別為0.7907和0.7950； 2) 相較于傳統BAEM指數提取方法，MBAEM指數較好解決了部分不透水面與水體、裸土混分的問題；相比其他提取不透水面指數，總體精度和Kappa系數均有不同程度的提升，且優化方法的提取結果在空間分布上更接近真實結果。

本研究提出的MBAEM指數可以有效識別遙感影像下具有復雜環境的半城市化地區的不透水面，部分解決了傳統指數裸土與不透水面混分的問題，且能挖掘不透水面的特征信息，進而較好地獲取半城市化地區不透水面。但本研究仍存在局限性：1)當選用不同的影像時，利用OTSU算法得到的每期影像閾值不同，存在不確定性，未來需綜合考慮多區域多時序，以得到一個相對具有代表性的閾值；2)該方法目前僅在小區域得到實證，在更大區域尺度的應用有待進一步研究。

感謝國家重點研發計劃(2021YFD1700600)項目為本研究提供資金支持；感謝國家科技資源共享服務平臺-國家地球系統科學數據中心湖泊-流域分中心(http://lake.geodata.cn)提供數據支持！