基于Sentinel-2A與NPP-VIIRS夜間燈光數據的城市建成區提取

劉智麗，張啟斌，岳德鵬，郝玉光，蘇 凱

(1.北京林業大學林學院，北京 100083；2.中國林業科學研究院沙漠林業實驗中心，巴彥淖爾 015200)

0 引言

近年來，中國的城市化水平大幅度提升。經國家統計局統計(http://www.stats.gov.cn/tjsj/)，2017年中國的城市化率已達58.52%，相比2010年增長了11.93%。然而快速城市化在提升人類生活水平的同時，也給城市帶來了一系列環境及社會問題，如生態惡化、熱島效應、環境污染、交通擁擠等[1]。城市建成區是一片具有建成環境的連續區域，包括城市用地(道路、建筑及工業設施)及其周邊的生態用地(公園綠地、行道樹、水體等)，作為衡量城市發展的一個重要指標，其范圍及空間分布不僅可以反映城市發展的規模及趨勢，如城市擴張水平、人口密度等[2]，還對城市發展規劃起著重要作用。因此，對城市建成區的提取得到越來越多的關注[3]。

隨著遙感技術的不斷發展，城市建成區提取的自動化程度與精度都在不斷提升[4-5]。夜間燈光數據(nighttime light data，NTL)能夠直接反映城市建成區的夜間亮度，分布連續[6]，但受城市邊緣亮度影響較大，因而提取準確度不高。因此，許多研究將NTL與光學遙感影像結合起來提取城市建成區，例如Zhang等[7]將NTL與歸一化植被指數(normalized difference vegetation index，NDVI)相結合，得出基于植被校正的城市夜間燈光指數(vegetation adjusted NTL urban index，VANUI)，有效豐富了城市邊緣信息[8]，但VANUI會放大水體的反射率，混淆城市邊緣地帶的水體與城市建成區，且空間分辨率不夠高，提取結果仍無法滿足更高精度的研究需求。

歐洲航天局發射的Sentinel-2A光學遙感衛星覆蓋了可見光、近紅外以及短波紅外范圍內的13個波段，空間分辨率最高可達10 m[9]，能夠準確識別城市建成區的邊界，可為城市邊界提取提供良好的數據支持。但目前將Sentinel-2A數據應用于城市建成區識別、提取的研究較少，因此，本文對VANUI指數加以改進，基于Sentinel-2A與NTL數據提出基于建筑校正的城市夜間燈光指數(building adjusted NTL urban index，BANUI)，以包頭市南部為研究區，進行城市建成區的提取。

1 研究區概況與數據源

1.1 研究區概況

包頭市位于內蒙古自治區中西部，總面積約為27 768 km2，平均海拔為2 000 m，中部地勢較高，由中部向南北兩側地勢逐漸變低，橫貫其中部的陰山山脈將包頭市分成中部山區、北部高原草場與南部平原地區3個部分，其中南部地區包括包頭市中心城區與土默特右旗。近年來，包頭市經濟發展迅速，南部平原地區城市規模不斷擴張，在2000—2016年間，包頭市城區面積增長了71.10%。但當地氣候干燥、地形起伏較大，生態環境脆弱，而城市的不斷擴張又加劇了對生態環境的破壞，由此引發了一系列環境問題。包頭市南部地區聚集了全市89.69%的人口，城市分布密集，具備進行城市建成區遙感監測的基礎，而中部和北部地區大多為高原山地和草原，人口稀疏、城市數量屈指可數，缺乏城市遙感監測的條件，因此將包頭市南部平原地區作為本次研究的研究區(圖1)。

圖1 研究區遙感影像Fig.1 Image of study area

1.2 數據源

本文所使用的數據有Sentinel-2A遙感影像數據、NPP-VIIRS NTL數據以及研究區2016年土地利用數據。遙感影像選用處理后的2016年6月份研究區Sentinel-2A Level-2A數據，來源于歐洲航天局網站(https://scihub.copernicus.eu/)，經Sen2Cor模型處理，包含12個波段。

目前使用較多的NTL數據有DMSP-OLS和NPP-VIIRS 2種，其中，NPP-VIIRS數據的空間、時間分辨率都比較高，并已經過輻射定標，更適宜應用于分析與研究[12]。因此，本文選用2016年6月份包頭NPP-VIIRS數據，下載自美國國家海洋大氣局網站(https://www.ngdc.noaa.gov/eog/viirs/download_dnb_composites.html)。數據預處理時將亮度值小于0的異常像元賦值為0。

研究區2016年土地利用調查數據來源于包頭市環境保護局，但其中缺少城市建成區的地面實況數據，因此基于此數據通過手動標注來獲取城市建成區的范圍，作為本次研究的驗證數據。

本文中使用的所有地理空間數據投影均為WGS1984 UTM投影。

2 研究方法

2.1 技術路線

本文提出BANUI指數用于研究區城市建成區的提取。首先，基于Sentinel-2A數據得到研究區的歸一化建筑指數(normalized difference build-up index，NDBI)，將其與NPP-VIIRS數據結合得到BANUI數據；然后，基于BANUI利用分水嶺分割算法提取包頭市城市建成區；最后，對提取結果進行對比分析和精度驗證。技術路線如圖2所示。

圖2 技術路線Fig.2 Technology roadmap

2.2 VANUI

為了將城市與農村更好地區分開來，增大兩者亮度值的差異，Zhang等[14]根據DMSP/OLS數據與NDVI提出VANUI指數，計算公式為

VANUI=(1-NDVI)NTL，

(1)

式中NDVI為歸一化植被指數。當VANUI接近0時，區域植被覆蓋較好；VANUI的值較大時，植被相對稀疏。但是計算NDVI時使用的MODIS數據空間分辨率僅為500 m，NTL數據中的細節信息無法加以利用；且在VANUI中，非植被覆蓋地區的亮度值會被過度放大，水體反射與其他異常亮度點也會對提取精度有所影響。

2.3 BANUI

2.3.1 NDBI的提取

由TM影像中提取的NDBI可以準確描述城市建筑用地信息，計算公式為

(2)

式中TM4和TM5分別為TM影像的近紅外和短波紅外波段的反射率。一般情況下，NDBI=0時，為耕地或林地；NDBI<0時，為水體；NDBI>0時，為建筑用地。

處理后的Sentinel-2A數據包含12個波段，不同波段的組合會導致城市建成區范圍有所差異，因此分析不同地類(建筑、植被、水體、未利用地)上的光譜特征并確定最佳波段組合具有重要意義。本文利用感興趣區(region of interest，ROI)工具選取研究區各地類的代表性像元，并利用ArcGIS10.2軟件進行分區統計，結果如圖3所示。

圖3 Sentinel-2A數據各波段的平均亮度值Fig.3 Mean brightness of each band in Sentinel-2A

參照式(2)中所使用波段的波長，本文從Sentinel-2A數據的第8、第11與第12波段中選出合適的波段進行組合。由圖3可知，在第11和12波段上，建筑用地的像元亮度值較高，水體和植被的亮度值遠小于建筑用地；在第8波段上，植被的亮度值達到極值，水體與建筑的亮度值相對較??；在第11、第8波段或第12、第8波段的亮度差值中，建筑為正數，植被為負數，水體為接近于0的負值。利用這3個波段的組合構建了2個NDBI，分別記作NDBI11和NDBI12，公式為

(3)

(4)

式中：S8,S11,S12分別為Sentinel-2A數據中第8、第11與第12波段的反射率。

為了直觀地對比NDBI11與NDBI12在不同地物上的表達差異，選取涵蓋了建筑、植被、水體與未利用地4種類型的矩形區域，并截取對應的NDBI圖像，結果如圖4所示。圖4顯示，NDBI11和NDBI12都能夠較好地反映建筑用地的范圍，建筑的亮度值明顯高于圖中其他2種地類，水體的亮度值在圖中并不高，植被的亮度值也明顯低于建筑，且植被覆蓋程度越高的地方亮度值越低。此外，2幅圖植被亮度差異不明顯，NDBI12中建筑的亮度略高，因此NDBI12中建筑與植被更容易區分；但NDBI12中水體的DN值略高于NDBI11的水體，這在一定程度上會干擾建筑用地的提取。為定量分析2種波段組合所獲取的NDBI的優劣，本文計算出各種地類的NDBI概率分布并繪制其正態分布曲線，結果如圖5所示。由圖5可知，2幅圖中建筑用地的NDBI值最高，水體次之，植被的NDBI值最低，其值遠小于0。2幅圖像中建筑用地與植被的差別都比較明顯，但在圖5(b)中，建筑與水體有較多重疊區域，此外建筑用地與未利用地在2幅圖中都難以區分。

(a)NDBI11(b)NDBI12

圖4 NDBI11和NDBI12的代表性區域

Fig.4RepresentativeareasofNDBI11andNDBI12

(a)NDBI11(b)NDBI12

圖5 不同地類在NDBI中的概率分布及其正態分布曲線

Fig.5ProbabilitydistributionandnormaldistributioncurvesofdifferentlandtypesinNDBI

基于上圖中正態分布曲線以及色譜分離分析理論可以計算建筑與其余各種地類之間的分離度R，目前該指數被廣泛應用于色譜柱中2個相鄰元素間分離程度的計算。一般當R<1時，2波峰有部分重疊；R=1時，其分離程度可達98%；R=1.5時，分離程度接近99.7%，此時2個相鄰元素完全分離；R=2時，2波峰從基線起就完全分離。通常情況下，利用色譜峰的保留時間以及基線上的峰寬可計算分離度R，本文則借助各地類正態曲線峰值所對應的NDBI值計算分離度，公式為

(5)

式中：BR1和BR2分別為2條正態曲線的峰值所在NDBI的值；W1和W2分別為2條正態曲線在基線上的峰寬。計算結果如表1所示。

表1 NDBI11與NDBI12中建筑與各地類之間的分離度RTab.1 R between buildings and the other land types in NDBI11 and NDBI12

由表1計算結果可知，建筑與植被的分離度均大于2，表明2種地類的正態曲線完全分離，在NDBI中差異明顯，但NDBI12中兩者的分離度更高；建筑與水體在NDBI11中的分離度高于NDBI12中的分離度，因此NDBI11中建筑與水體差別明顯；建筑與未利用地的分離度都比較低、不易區分。由此可以推斷，除了建筑與未利用地，其他地類信息在NDBI11中都得到了不同程度的抑制，而NDBI12對水體信息的抑制作用并不明顯，NDBI11與NTL數據的結合可以有效消除噪聲的影響、并增強建筑信息的表達。因此本文選擇Sentinel-2A數據中第8波段和第11波段來提取NDBI。

2.3.2 BANUI的提取

在VANUI中，水體的值相對偏離正常范圍，且空間分辨率較低。為解決這一問題，本文利用NDBI與NPP-VIIRS數據構建了BANUI指數，公式為

BANUI=(1+NDBI)VNTL，

(6)

式中VNTL為NPP-VIIRS數據的亮度值。經過重采樣后的NPP-VIIRS數據的空間分辨率為20 m，最終BANUI的分辨率同樣為20 m。為檢驗BANUI的提取效果，同樣利用Sentinel-2A與NPP-VIIRS數據計算得到VANUI，其空間分辨率為10 m。

選取E109.35°～110.85°范圍內一橫截面，該橫截面可以最大程度地經過包頭市南部城市建成區，信息豐富，基于此可檢驗BANUI中建筑信息的表達。橫截面上VANUI與BANUI的值如圖6所示。由圖6可知，兩者變化趨勢比較相似；橫截面經過包頭市一公園，公園內部植被覆蓋度較高，且因NDVI對植被更敏感，導致VANUI的值要比BANUI更低些。

圖6 橫截面上VANUI與BANUI的值Fig.6 Values of VANUI and BANUI on the cross section

表2數據為研究區內VANUI和BANUI的像元亮度值范圍、像元亮度均值、像元亮度值標準差，以及水體、植被、建筑的像元亮度均值。由于未利用地距離城區較遠，其面積只有153.47 km2(研究區總面積為4 965.78 km2)，對于城市建成區的提取影響甚微，因此未利用地的數據并未包含在其中。

表2 VANUI和BANUI的統計數據Tab.2 Statistical data of VANUI and BANUI

從表2可以看出，BANUI像元亮度值覆蓋范圍較廣、標準差更高，其變化幅度更大，所涵蓋的信息更多；水體在BANUI中的像元亮度均值較低，建筑在BANUI中的像元亮度均值較高，表明BANUI更好地凸顯了建筑信息。

圖7為包頭市城市建成區邊緣地帶，包含了水體、植被以及連續分布建筑等信息，因此將其作為檢驗BANUI區分建筑用地效果的典型區域。

(a)典型區域位置 (b)典型區域影像

圖7 典型區域

Fig.7Typicalarea

該區域水體、建筑、植被的概率分布以及正態分布情況如圖8所示。

(a)VANUI (b)BANUI

圖8 典型區域內不同地類在VANUI與BANUI中的概率分布及正態分布曲線

Fig.8ProbabilitydistributionandnormaldistributioncurvesofdifferentlandtypesofVANUIandBANUIinthetypicalarea

圖8顯示，BANUI中建筑與植被的像元亮度值要高于VANUI中建筑與植被的亮度值，水體的亮度值在兩幅圖像中相差不大，導致圖8(b)中建筑與水體距離更遠，二者更容易區分。

計算VANUI和BANUI中建筑與其他地類的分離度，結果顯示，建筑與植被的分離度都超過了1.5(1.76和1.52)；而BANUI中建筑與水體的分離度(0.63)明顯高于VANUI中二者的分離度(0.29)，表明BANUI相比VANUI能夠更好地將建筑與水體區分開來。由此可見，BANUI可以豐富NPP-VIIRS的信息、更好地抑制水體反射造成的影響，并增強建筑信息的表達，因此將其視作提取城市建成區的有效數據。

2.4 分水嶺分割算法

為檢驗BANUI在城市建成區提取應用中的優勢，本文利用分水嶺分割算法從NPP-VIIRS與VANUI數據中提取城市建成區。

分水嶺分割算法是一種基于區域進行影像分割的有效方法，可用于城市建成區范圍的提取[15]。BANUI數據的像元被視作地形表面，每個像元的亮度值為地形表面上相應像元的高程值，亮度值高的像元為高地，亮度值低的像元為低地，城市建成區的邊界為地形中集水盆地的邊界即分水嶺[16]。

主要步驟如下：①首先，對影像進行平滑處理，去除BANUI影像中椒鹽噪聲的影響，自適應加權中值濾波方法可在去除噪聲的同時保留更多信息，而高斯濾波、中值濾波等常用方法則無法滿足這一需求；②利用形態重構提高影像中同類物體的同質性，并降低不同物體之間的異質性，由此得到形態重構影像；③計算影像的形態梯度，通常情況下像元亮度值變化程度較高的地區，影像梯度較大；④借助形態重構影像求得局部極大值并標記為前景對象，通過歐氏距離倒數法標記背景影像；⑤利用分水嶺分割算法劃分標有前景與背景對象的待處理影像，提取出城市建成區。

2.5 對比分析與精度驗證

為了檢驗本文提取城市建成區所采取方法的有效性，選取包頭市2016年的土地利用調查數據作為對照，通過人機交互方式在其基礎上繪制出城市建成區的范圍作為驗證數據，此范圍具有高精度、連續性高等特點。

通過實驗結果與驗證數據之間逐像元的對比建立混淆矩陣，進而計算總體精度、Kappa系數、生產者精度和用戶精度等指標進行精度檢驗[17]。此外，需借助景觀指數來描述結果空間結構上的特征，包括總面積(total area，TA)、景觀形狀指數(landscape shape index，LSI)和聚集指數(aggregation index，AI)等。

LSI可衡量整個景觀的形狀與面積相同的正方形或圓形之間的偏差程度，以此來表示該景觀整體形態的復雜程度[18]，取值范圍為[1，+∞)。LSI值越大，該景觀形狀與正方形的偏差越大。公式為

(7)

式中：0.25為正方形校正常數；E為景觀中所有斑塊的邊緣總長度，100 m；A為景觀總面積，hm2。

AI表示景觀內部斑塊之間相臨的格網單元數目與在這些數中的最大值的比值，它反映了景觀的聚集程度與連續性[19]。公式為

(8)

式中gi為景觀內斑塊i相鄰的格網單元數。AI的取值范圍為(0,100%]，當景觀內部斑塊無限分散時，AI的值趨向于0；隨著斑塊的逐漸集聚，AI的值也隨之增加。

3 成果分析與檢驗

3.1 提取結果

利用分水嶺分割算法由不同數據中提取的城市建成區(圖9)。由圖9可知，利從BANUI中提取的城市建成區與驗證數據相似度較高，而從VANUI和NPP-VIIRS中提取的城市建成區相比于驗證數據范圍較大。由BANUI中提取的城市建成區包括城市用地，以及水體、公園等其他地類，空間分布連續、成分完整，與上文對城市建成區的定義相吻合。

(a)NTL (b)VANUI (c)BANUI

圖9 由不同數據提取的城市建成區

Fig.9Extractedurbanbuilt-upareasfromdifferentdata

3.2 對比分析與精度驗證

通過上文所述精度驗證方法得到不同數據的提取精度，如圖10所示。從圖10可以看出，由BANUI提取出的城市建成區總體精度最高，為93.61%，Kappa系數為0.793 4，城市建成區的生產者精度與用戶精度分別為85.34%和81.34%，總體上提取的精度較高；VANUI位列其次，NPP-VIIRS精度最低。

圖10 不同數據源提取的城市建成區的精度對比Fig.10 Accuracy comparison of urban built-up areas extracted by different data sources

表3為不同數據源提取的城市建成區景觀格局指數對比，其中，TA.Diff，LSI.Diff，AI.Diff分別表示TA，LSI，AI與驗證數據的誤差。由表4可知，BANUI提取結果的景觀指數與驗證數據十分接近，提取的城市建成區總面積比驗證數據多了4.91%，LSI和AI的差距也較小，提取效果最佳；而其余2個數據源所提取結果與驗證數據相差較大，面積分別多出了24.12%和16.55%，精度較低。

表3 不同數據源提取的城市建成區景觀格局指數對比Tab.3 Comparison of landscape pattern indexes of urban built-up areas extracted by different data sources

4 結論

1)基于Sentinel-2A與NPP-VIIRS數據，對傳統的VANUI進行改進，提出了用于提取城市建成區的BANUI指數。

2)通過對Sentinel-2A數據各個波段的特征進行分析，選擇Sentinel-2A影像的第8和第11波段計算研究區NDBI，從而構建BANUI指數

3)實驗結果表明，由BANUI所提取的研究區城市建成區的面積為971.25 km2，相比驗證數據僅相差4.91%，總體精度高達93.61%，Kappa系數為0.793 4，明顯優于由VANUI和NPP-VIIRS數據的提取結果。

4)BANUI在城市建成區的提取上具有結果完整、覆蓋范圍廣、提取精度高等優點，應用前景廣闊，可為進一步研究NTL數據對城市建成區的提取提供思路，也可用于城市發展的規劃與監測。