朱凌一，王 勇

（1.南京信息工程大學(xué)長望學(xué)院，江蘇 南京 210044；2.南京信息工程大學(xué)應(yīng)用氣象學(xué)院，江蘇 南京 210044）

目前研究所使用的夜間燈光數(shù)據(jù)主要有DMSP-OLS 與NPP-VIIRS數(shù)據(jù)，后者擁有更高的空間分辨率和時(shí)間分辨率，且已經(jīng)經(jīng)過輻射定標(biāo)。在對(duì)城市化的監(jiān)測(cè)應(yīng)用中，NPP-VIIRS數(shù)據(jù)也被證實(shí)其能力要優(yōu)于DMSP-OLS 數(shù)據(jù)[1]。隨著最新的NPP-VIIRS 夜光數(shù)據(jù)遙感影像發(fā)布，基于NPP-VIIRS數(shù)據(jù)提取城市建成區(qū)的方法體系還需完善[2]，目前急需基于NPP-VIIRS 數(shù)據(jù)提取城市群建成區(qū)的簡(jiǎn)單高效普適方法。因此，本文將基于2020年NPP-VIIRS夜光遙感數(shù)據(jù)，嘗試構(gòu)建基于差值平方函數(shù)的城市群建成區(qū)提取模型，并對(duì)蘇南地區(qū)建成區(qū)進(jìn)行特征分析，從面積、混淆矩陣和景觀指數(shù)等方面將所建模型與閾值法進(jìn)行對(duì)比分析，以期對(duì)蘇南城市群國土空間規(guī)劃制定、生態(tài)環(huán)境保護(hù)提供技術(shù)支持和決策參考。

1 研究區(qū)概況

本文選取蘇南地區(qū)城市群（南京市、蘇州市、鎮(zhèn)江市、無錫市、常州市）作為研究區(qū)域，該地區(qū)位于長江三角洲中心區(qū)域，東鄰上海，西接安徽，南靠浙江，北依長江；區(qū)內(nèi)低山、丘陵、平原、江河、湖泊縱橫交錯(cuò)，東部和中部較為低平，擁有廣袤的太湖平原，水網(wǎng)密集，西部南京、鎮(zhèn)江一帶丘陵起伏；土地總面積2.79 萬km2，占江蘇省土地總面積27.17%，其中平原占蘇南土地總面積50.45%，山地面積占28.4%，水域面積占21.15%。氣候呈現(xiàn)四季分明、冬冷夏熱，降水量在1 000 mm 以上，梅雨期內(nèi)暴雨頻發(fā)，強(qiáng)降水集中。蘇南地區(qū)是中國經(jīng)濟(jì)最發(fā)達(dá)的區(qū)域之一，常住人口3 378.09 萬人（2019年），城鎮(zhèn)化率超過70%，所有縣（市）都進(jìn)入全國綜合實(shí)力百強(qiáng)縣行列，近年來蘇南城市群建成區(qū)面積不斷擴(kuò)大。

2 數(shù)據(jù)來源與方法

2.1 數(shù) 據(jù)

本文采用美國國家大氣和海洋管理局國家地理數(shù)據(jù)中心提供的2020年的NPP-VIIRS夜光遙感影像，空間分辨率約500 m，能夠更精確、豐富地反映地表人類活動(dòng)的空間信息。VIIRS 數(shù)據(jù)產(chǎn)品有“VCM”和“VCMSL”2種類型，為保證結(jié)果的準(zhǔn)確性，本文選取“VCMSL”類型的VIIRS 影像，對(duì)其進(jìn)行投影轉(zhuǎn)換、重采樣、裁剪等預(yù)處理，轉(zhuǎn)換至分辨率為500 m 的阿爾伯斯等積圓錐投影坐標(biāo)系，并裁剪蘇南地區(qū)2020年月平均NPP-VIIRS 夜光遙感影像。此外，由于VIIRS數(shù)據(jù)產(chǎn)品并未去除陸地反射月光等背景噪聲，本文還去除VIIRS影像的背景值。

為了驗(yàn)證建成區(qū)提取效果，同時(shí)考慮到Google Earth平臺(tái)數(shù)據(jù)更新速度快，保證了數(shù)據(jù)的時(shí)效性，適合對(duì)成片的城市建成區(qū)進(jìn)行目視解譯[3]。本文采用Google Earth 軟件所提 供的2020 年12 月Landsat8 衛(wèi)星城市歷史影像，在研究區(qū)內(nèi)進(jìn)行隨機(jī)采樣，同時(shí)對(duì)研究區(qū)隨機(jī)采樣點(diǎn)進(jìn)行目視解譯判斷，并對(duì)本文所用方法的分類精度進(jìn)行驗(yàn)證。

本文以GlobeLand30-2020土地利用類型數(shù)據(jù)為基準(zhǔn)，從面積、混淆矩陣和景觀指數(shù)等方面對(duì)比不同方法的建成區(qū)精度。GlobeLand30-2020數(shù)據(jù)是我國研制的30 m分辨率全球地表覆蓋數(shù)據(jù)集，自然資源部對(duì)GlobeLand30數(shù)據(jù)集進(jìn)行的最新更新，包含了更加豐富詳盡的全球地表覆蓋空間分布信息，能更好地刻畫大多數(shù)人類土地利用活動(dòng)及其所形成的景觀格局，此數(shù)據(jù)是目前可免費(fèi)獲得的總體分類精度最高的地表覆蓋數(shù)據(jù)。

2.2 數(shù)據(jù)預(yù)處理

由于NPP-VIIRS傳感器較為敏感，城市建成區(qū)內(nèi)部燈光輻射值增速較快[4]，在分析過程中會(huì)掩蓋城郊過渡區(qū)的燈光變化較為平緩的特征。為了避免計(jì)算過程中像元DN值極端效應(yīng)帶來的影響，本文對(duì)原始NPP-VIIRS數(shù)據(jù)（圖2a）進(jìn)行篩選處理，結(jié)合先驗(yàn)知識(shí)與前人研究資料[5]，斷定DN值大于20處為城市建成區(qū)（DN值≥20 的像元賦值為20），DN值小于1 處為非城市建成區(qū)（DN值≤1 的像元賦值為1），以此消除原始數(shù)據(jù)中的城市內(nèi)部燈光突變和農(nóng)村極小值噪聲干擾，結(jié)果見圖2b。為利用差值平方函數(shù)來進(jìn)一步精確計(jì)算城市建成區(qū)與非建成區(qū)分界線提供數(shù)據(jù)基礎(chǔ)。

2.3 方 法

2.3.1 差值平方函數(shù)

統(tǒng)計(jì)學(xué)中，方差是用來描述數(shù)據(jù)離散程度的數(shù)學(xué)函數(shù)，常常用來衡量數(shù)值之間的差異程度。本文基于方差的思路，使用相鄰像元DN值之間的差值平方，來確定相鄰像元DN值之間的差異，尋找這種差異變化的突變點(diǎn)，并將其作為城市建成區(qū)與非城市建成區(qū)的邊界像元，從統(tǒng)計(jì)學(xué)意義上可以被視為是劃分建成區(qū)與非建成區(qū)的標(biāo)準(zhǔn)，即定義為差值平方函數(shù)。差值平方函數(shù)公式為：

式中，DNi為第i個(gè)像元的燈光指數(shù)；Δ 為統(tǒng)計(jì)學(xué)意義上相鄰2個(gè)像元的夜光指數(shù)的差值平方值。

2.3.2 驗(yàn)證方法

Google Earth影像空間比較法[6]在建成區(qū)提取精度分析方面具有操作簡(jiǎn)單和時(shí)效性強(qiáng)等優(yōu)勢(shì)。因此，本文選取Google Earth 平臺(tái)于2020 年12 月成像的Landsat8 遙感影像作為研究區(qū)建成區(qū)和非建成區(qū)的驗(yàn)證數(shù)據(jù)，依據(jù)城市建成區(qū)在影像中的表現(xiàn)特征，通過目視解譯的方式對(duì)建成區(qū)進(jìn)行識(shí)別，同時(shí)ArcGIS軟件在研究區(qū)內(nèi)隨機(jī)產(chǎn)生1 500個(gè)采樣點(diǎn)，根據(jù)建成區(qū)和非建成區(qū)對(duì)采樣點(diǎn)進(jìn)行分類，隨后與本文所用方法的分類結(jié)果進(jìn)行對(duì)比分析，再進(jìn)行精度評(píng)價(jià)。其中，分類精度的計(jì)算公式為：

式中，R為Google Earth與本文所用方法的提取結(jié)果相一致的點(diǎn)數(shù)；S為Google Earth 提取結(jié)果不一致的點(diǎn)數(shù)。通過對(duì)P值的分析，得出采樣點(diǎn)分類的正確率。

2.3.3 混淆矩陣

在遙感影像處理中，混淆矩陣可以衡量遙感影像地物分類結(jié)果與地表真實(shí)地物結(jié)果之間的差異，這些差異可以用來描述遙感影像的分類精度。

本文使用建成區(qū)與非建成區(qū)分類結(jié)果作為實(shí)際劃分建立建成區(qū)混淆矩陣，建成區(qū)混淆矩陣可以解讀為：有x11個(gè)像元被正確分類為建成區(qū)，x22個(gè)像元被正確分類為非建成區(qū)；有x12個(gè)建成區(qū)像元被錯(cuò)分為非建成區(qū)像元，有x21個(gè)非建成區(qū)像元被錯(cuò)分為建成區(qū)像元。

基于以上結(jié)果，本文將進(jìn)一步利用混淆矩陣，計(jì)算總體精度（OA）與Kappa 系數(shù)等分類精度指標(biāo)。其中總體精度（OA）表示分類正確的像元點(diǎn)占總像元點(diǎn)的百分比，Kappa 系數(shù)常被視為衡量模型分類結(jié)果與實(shí)際分類結(jié)果是否一致的指標(biāo)，值越大關(guān)系越好。具體公式如下：

式中，N為總樣本點(diǎn)數(shù)，本文取值為1 500；r為總行數(shù)，本文取值為2；xii為主對(duì)角線上的值；xi+與x+j分別表示第i行和第j列的樣本個(gè)數(shù)之和。本文將通過混淆矩陣計(jì)算總體精度與Kappa 系數(shù)，實(shí)現(xiàn)對(duì)提取建成區(qū)的精度評(píng)價(jià)。

2.3.4 景觀指數(shù)

景觀指標(biāo)是分析地形空間模式的重要指標(biāo)，也是驗(yàn)證城市區(qū)域制圖精度的有效方法。由于混淆矩陣精度驗(yàn)證法無法展示提取的城市建成區(qū)與驗(yàn)證數(shù)據(jù)提取的城市建成區(qū)在空間和組織上的對(duì)比性，引入景觀指標(biāo)指數(shù)對(duì)比能夠使精度驗(yàn)證更完善、全面、準(zhǔn)確[7]。

常用的景觀指數(shù)有景觀形狀指數(shù)（LSI）、周長面積比（PARA）和聚集指數(shù)（CONTIG）。其中LSI定義為景觀格局中斑塊的形狀指數(shù)，通過計(jì)算區(qū)域內(nèi)某斑塊形狀與相同面積的圓或正方形之間的偏離程度來測(cè)量形狀的破碎程度，因此當(dāng)建成區(qū)形態(tài)趨于破碎，LSI指數(shù)也隨之提高；PARA是以斑塊周長于面積的比值直接衡量形狀復(fù)雜程度的簡(jiǎn)單指標(biāo)，當(dāng)PARA 指數(shù)增加時(shí)，城市建成區(qū)邊界也更間斷、離散；CONTIG指數(shù)反映地物中不同斑塊類型的非隨機(jī)性或聚集程度[8]，其值越大表明建成區(qū)越連續(xù)、聚集。本文將GlobeLand30-2020 中的人造地表數(shù)據(jù)作為驗(yàn)證樣本，選取景觀形狀指數(shù)LSI、周長面積比PARA和聚集指數(shù)CONTIG等3個(gè)指數(shù)進(jìn)行景觀指數(shù)分析。

3 結(jié)果分析

3.1 城市群建成區(qū)提取模型構(gòu)建

3.1.1 模型技術(shù)方案

為了尋找一種操作流程少，計(jì)算量不大且簡(jiǎn)單高效普適性強(qiáng)的城市群建成區(qū)提取方法，本文基于夜光遙感影像中的紋理特征，結(jié)合差值平方函數(shù)，挖掘夜光遙感影像的城市地物分布特征，實(shí)現(xiàn)建成區(qū)和非建成區(qū)的分類，本文將其定義為差值平方法。具體技術(shù)路線（圖1）為構(gòu)建基于差值平方函數(shù)的蘇南城市群建成區(qū)提取模型。首先使用NPP-VIIRS夜光遙感影像數(shù)據(jù)，開展幾何校正、去除噪聲、提取研究區(qū)、篩選燈光亮度范圍等預(yù)處理，結(jié)合差值平方函數(shù)逐行逐列計(jì)算，提取了建成區(qū)范圍，然后與Google earth影像目視解譯的采樣點(diǎn)進(jìn)行對(duì)比，實(shí)現(xiàn)空間精度驗(yàn)證。

本文所建的蘇南城市群建成區(qū)提取模型的核心算法，就是采用差值平方函數(shù)檢測(cè)邊界像元。在遙感影像中，當(dāng)像元間波動(dòng)性達(dá)到一點(diǎn)程度，即2個(gè)像元DN值之間的差值平方達(dá)到局部極大值，可認(rèn)為用來計(jì)算差值平方的2個(gè)像元分別屬于兩類不同的地物類型，而像元的交界處就是兩類地物的邊界。圖1 給出經(jīng)過差值平方計(jì)算后如何劃分建成區(qū)與非建成區(qū)的示例，以一行為例，圖1 的上部分為一行像元DN值，下部分為一行相鄰兩像元DN值的差值平方值。當(dāng)DN值由較小值過度到較大值時(shí)，影像上具體表現(xiàn)為夜光遙感指數(shù)由低變高，即由非城市建成區(qū)過渡到城市建成區(qū)，而差值平方在DN值為10.877和16.433之間時(shí)達(dá)到最大值（30.87），故兩點(diǎn)之間為城市建成區(qū)與非城市建成區(qū)的分界點(diǎn)，以此方法實(shí)現(xiàn)對(duì)夜光遙感影像的逐行逐列計(jì)算，提取研究區(qū)內(nèi)的各個(gè)建成區(qū)邊界點(diǎn)，再采用泰森多邊形建立各個(gè)邊界點(diǎn)DN值的代表范圍，進(jìn)而采用柵格計(jì)算器實(shí)現(xiàn)對(duì)建成區(qū)范圍的提取，即將鄰近分界點(diǎn)的夜光指數(shù)柵格DN值≥分界點(diǎn)DN值，此柵格定為城市建成區(qū)，否則定為非城市建成區(qū)。

3.1.2 精度驗(yàn)證結(jié)果

本文基于蘇南地區(qū)2020年月平均NPP-VIIRS夜光遙感影像，采用上文所構(gòu)建的建成區(qū)提取模型，提取了蘇南城市群建成區(qū)范圍。采用ArcGIS軟件隨機(jī)產(chǎn)生1 500 采樣點(diǎn)和2020 年12 月Landsat8 遙感影像作為驗(yàn)證樣本，根據(jù)目視解譯對(duì)采樣點(diǎn)進(jìn)行建成區(qū)和非建成區(qū)分類，再采用混淆矩陣與差值平方法分類結(jié)果進(jìn)行對(duì)比，總體精度為93.6%，表明基于差值平方法提取的建成區(qū)與蘇南地區(qū)實(shí)際建成區(qū)分布結(jié)果具有高度的一致性。

3.2 蘇南城市群建成區(qū)特征分析

根據(jù)基于差值平方函數(shù)的NPP-VIIRS夜光遙感影像建成區(qū)提取模型，實(shí)現(xiàn)了蘇南城市群建成區(qū)的提取，結(jié)果如圖2所示。圖2a為差值平方函數(shù)所得邊界點(diǎn)，已經(jīng)凸顯了城市建成區(qū)與非城市建成區(qū)之間的過渡帶，可以大致反映出城市內(nèi)部的建成區(qū)空間分布特征。

在所求得邊界點(diǎn)的基礎(chǔ)上，通過ArcGIS軟件中泰森多邊形實(shí)現(xiàn)各個(gè)邊界點(diǎn)DN值的代表范圍，進(jìn)而采用柵格計(jì)算器實(shí)現(xiàn)蘇南城市群建成區(qū)的提取，最終分類結(jié)果如圖2b所示。蘇南城市群中城市建成區(qū)總體呈現(xiàn)蘇錫常趨于一體、南京獨(dú)占一角的分布態(tài)勢(shì)，從數(shù)據(jù)上來看，南京、鎮(zhèn)江、常州、無錫、蘇州的城市像元所占比例分別為22.64%、19.45%、19.94%、28.36%與26.84%；按建成區(qū)面積指標(biāo)排序?yàn)樘K州（2 144.63 km2）>南京（1 703.12 km2）>無錫（1 305.53 km2）>常州（869.39 km2）>鎮(zhèn)江（747.80 km2）。差異的產(chǎn)生與各個(gè)城市不同的自然、社會(huì)經(jīng)濟(jì)條件及不同的城市規(guī)劃方案等有著密切的關(guān)系。與實(shí)際的建成區(qū)相對(duì)比發(fā)現(xiàn)，提取結(jié)果與實(shí)際建成區(qū)像元占比、建成區(qū)實(shí)際位置、城市實(shí)際土地利用類型較為一致，表明差值平方法在能對(duì)城市實(shí)際建成區(qū)的比例與城市空間范圍起到很好的解釋作用。

3.3 與閾值法對(duì)比

目前，閾值法是在城市建成區(qū)提取中最常見的方法，其利用二分法計(jì)算影像分類閾值。當(dāng)在某一閾值下的城市建成區(qū)面積最接近于統(tǒng)計(jì)數(shù)據(jù)的數(shù)值時(shí)即可認(rèn)為此閾值為分割閾值[7,10]，通過確定分割閾值進(jìn)行分類，實(shí)現(xiàn)城市建成區(qū)范圍的劃分與確定。為進(jìn)一步驗(yàn)證差值平方法是否有效提高了城市建成區(qū)的提取精度，本文以GlobeLand30-2020土地利用類型數(shù)據(jù)為基準(zhǔn)，對(duì)差值平方法與閾值法提取的建成區(qū)范圍進(jìn)行比較，并從面積、混淆矩陣和景觀指數(shù)等方面定量分析了建成區(qū)提取精度。

1）面積精度分析。本文統(tǒng)計(jì)了研究區(qū)內(nèi)各市的建成區(qū)面積，結(jié)果見表1。以GlobeLand30-2020 中的人造地表數(shù)據(jù)的面積作為驗(yàn)證數(shù)據(jù)，采用差值平方法的城市建成區(qū)面積準(zhǔn)確率平均值為97.71%，而高于采用率閾值法的92.44%，表明在夜間燈光影像，差值平方法提取的建成區(qū)精度更高，提取面積更接近真實(shí)城市建成區(qū)。差值平方法不僅可以保證城市建成區(qū)的整體形態(tài)，還可以剔除城市建成區(qū)中的非城市建成區(qū)影響，效果較為理想。

表1 基于差值平方法和閾值法的城市建成區(qū)面積精度評(píng)價(jià)/km2

2）混淆矩陣精度評(píng)價(jià)。在定量評(píng)價(jià)指標(biāo)方面，面積誤差和混淆矩陣常用于衡量圖像提取和分類精度。采用混淆矩陣進(jìn)行精度評(píng)價(jià)，計(jì)算出總體精度和Kappa 系數(shù)2個(gè)指標(biāo)，2 種方法的對(duì)比結(jié)果見表2。差值平方法與閾值法的平均總體精度為84.36%與84.40%，2 種方法差別很小，而它們的Kappa 系數(shù)分別為56.23%與55.38%，這表明總體上差值平方法提取效果優(yōu)于閾值法。

表2 基于差值平方法和閾值法的城市建成區(qū)混淆矩陣精度評(píng)價(jià)/%

3）景觀指數(shù)精度評(píng)價(jià)。根據(jù)差值平方法和閾值法的建成區(qū)提取結(jié)果，采用Fragstats 軟件進(jìn)行景觀指數(shù)分析?？紤]到GlobeLand30-2020數(shù)據(jù)空間分辨率較高（30 m）導(dǎo)致不同用地類型邊緣破碎的情況，本文將其進(jìn)行重采樣為500 m×500 m 的空間分辨率，使其與差值平方法和閾值法的提取結(jié)果柵格一致，再進(jìn)行景觀指數(shù)分析（表3）。從表3中的結(jié)果來看，差值平方法和閾值法的LSI 指數(shù)相近，與GlobeLand30-2020相距較大；差值平方法與GlobeLand30-2020土地利用類型數(shù)據(jù)的PARA 指數(shù)接近，而CONTIG指數(shù)完全一致，總體優(yōu)于閾值法，因此差值平方法提取的城市內(nèi)部建成區(qū)與GlobeLand30-2020土地利用類型數(shù)據(jù)中人造地表數(shù)據(jù)呈現(xiàn)更高的景觀一致性，更接近實(shí)際情況。

表3 基于差值平方法和閾值法的城市建成區(qū)景觀指數(shù)精度評(píng)價(jià)

4 結(jié)論

本文利用蘇南地區(qū)2020年月平均NPP-VIIRS夜光遙感影像，建立了基于差值平方函數(shù)的城市群建成區(qū)提取模型，實(shí)現(xiàn)了蘇南地區(qū)城市群建成區(qū)范圍的提取，結(jié)合GlobeLand30-2020土地利用分類數(shù)據(jù)，與閾值法進(jìn)行對(duì)比分析，具體結(jié)論如下：

1）本文對(duì)NPP-VIIRS 夜光遙感影像數(shù)據(jù)開展幾何校正、去除噪聲、提取研究區(qū)、篩選燈光亮度范圍等預(yù)處理，構(gòu)建了基于差值平方函數(shù)的建成區(qū)提取模型。并采用ArcGIS 隨機(jī)產(chǎn)生1 500個(gè)采樣點(diǎn)和Landsat8 遙感影像作為驗(yàn)證樣本，得出差值平方法提取的建成區(qū)與蘇南地區(qū)實(shí)際建成區(qū)分布結(jié)果具有高度的一致性。

2）基于差值平方法和NPP-VIIRS 夜光遙感影像所構(gòu)的建成區(qū)提取模型，實(shí)現(xiàn)了蘇南城市群建成區(qū)的提取。南京、鎮(zhèn)江、常州、無錫、蘇州的建成區(qū)像元比例分別為22.64%、19.45%、19.94%、28.36%與26.84%，按建成區(qū)面積指標(biāo)排序?yàn)樘K州>南京>無錫>常州>鎮(zhèn)江。與實(shí)際的建成區(qū)相對(duì)比發(fā)現(xiàn)，采用差值平方法能對(duì)城市實(shí)際建成區(qū)的占比與城市空間范圍起到很好的解釋作用。

3）從面積、混淆矩陣和景觀指數(shù)等方面對(duì)差值平方法與閾值法提取的建成區(qū)進(jìn)行精度分析，差值平方法提取的城市建成區(qū)面積準(zhǔn)確率高于率閾值法，差值平方法與閾值法的平均總體精度相近，而差值平方法的Kappa 系數(shù)高于閾值法，差值平方法的景觀指數(shù)更接近實(shí)際情況，進(jìn)一步證明差值平方法能夠有效提高城市建成區(qū)的提取精度。