典型喀斯特山地城市不透水面時空格局

母浩江， 崔文剛， 魯 霞

(貴州師范大學地理與環境科學學院， 貴州貴陽550025)

近年來，隨著我國城市化進程不斷加快，大量農村人口向城市集聚，在促進城市經濟發展的同時也導致城市用地壓力日益增大，越來越多的森林、 耕地等自然地表被不透水地表所替代。不透水地表又被稱為不透水面(impervious surface，IS)，是指常見的道路、 廣場、 建筑等水不能通過其下滲到土壤中的城市人工景觀[1-2]。IS的大量增加會改變地表下墊面性質[3]，使地表徑流和地下水的活動受到干擾，致使城市內澇、 地下暗河斷流等現象產生[4-5]，進而給生態環境帶來負面影響。及時、 準確地掌握城市IS的分布情況和時空演變對于城市發展規劃及生態環境保護具有重要意義。

早期IS信息的獲取主要依靠專業人員進行實地采集， 這樣的獲取方式成本高， 周期長， 時效性差， 使得IS信息的采集和應用十分困難[6]。 進入20世紀以來， 隨著遙感技術的飛速發展， 遙感數據獲取變得便利， 國內外大量學者開始將遙感技術引入IS的研究。 隨著研究的深入， 多種基于遙感技術定量獲取IS信息的方法相繼被提出， 較典型的方法可概括為光譜混合分析(spectral mixture analysis，SMA)法[7-8]、 決策樹模型[9-10]、 回歸模型[11]、 指數法[12-13]和分類法[14]5種， 其中線性光譜混合分析(linear spectral mixture analysis，LSMA)法可以有效解決中低分辨率遙感數據的混合像元問題， 被廣泛用于提取IS信息的研究[15]。 Wu等[7]使用LSMA法成功提取了美國大都市區哥倫布的IS信息； 唐菲等[8]基于陸地衛星7(Landsat 7)增強型專題制圖儀(enhanced thematic mapper plus, ETM+)影像， 利用LSMA法提取了北京等6個不同城市的IS信息， 效果良好； 李苗等[16]通過LSMA法提取區域多個時間節點的IS信息， 分析IS信息的時空變化規律， 進一步探討城市化進程的范圍和規模。 理論方法的日趨成熟也推動著IS的研究向不同的區域環境與應用領域擴展， 但是通過查閱文獻發現， 相關研究主要集中于經濟較為發達的沿海城市[17]和平原地區[18]， 對于城市化進程起步稍晚且地形條件更復雜的喀斯特地區， IS的研究極少。

貴陽市是貴州省省會， 同時也是世界上喀斯特地貌發育最典型的地區之一[19]， 境內喀斯特地貌占行政區域面積的73.3%， 地形起伏較大， 景觀破碎化嚴重， 生態環境脆弱且具有明顯的空間異質性。 近年來， 隨著經濟實力和城鎮化水平不斷提高， 貴陽已成為我國西部陸海貿易的新通道、 西南地區的中心城市之一。 在城市化狀況快速發展的同時， IS急速增加所造成的環境問題也引起人們的重視， 然而目前尚未發現關于該地區IS信息及其空間變化的研究。

本文中以貴陽市主城區為研究區，采用基于植被覆蓋度提取法和LSMA法分別提取研究區的IS信息，利用混淆矩陣對提取結果進行精度驗證，選出最優方法對研究區1996—2017年各時間節點的IS數據進行提取，依據結果定量分析21 a間研究區IS的時空分布和演化規律。

1 研究區概況與數據源

1.1 研究區概況

貴陽市位于我國的西南部(東經106°07′—107°17′， 北緯26°11′—27°22′)，地處長江和珠江兩大水系流域的分水嶺地段，四面環山，地勢西高東低，境內丘陵起伏，屬典型的喀斯特地形地貌，氣候類型為亞熱帶季風氣候[20]。本文中選取該市主城區為研究區，包括南明區、 云巖區、 觀山湖區、 白云區、 花溪區、 烏當區6個行政區，其中云巖區、 南明區大部分區域位于盆地，烏當區、 花溪區、 白云區和觀山湖區大部分區域以山地為主。

1.2 數據源與預處理

為了避免由影像數據源不同、 拍攝季節不同、 云霧遮蔽等原因造成的誤差， 保證研究數據的可靠性， 本文中選取研究區內無云遮蔽且拍攝月份相近的3期Landsat影像作為數據源， 影像數據包含4景陸地衛星5號(Landsat 5)專題成像儀(TM)影像和2景陸地衛星8號(Landsat 8)陸地成像儀(OLI)影像(見表1)。 影像數據來自中國科學院地理空間數據云平臺(http://www.gscloud.cn)。前期需要對影像數據進行預處理， 包括輻射校正、 大氣校正、 拼接以及按研究區范圍裁剪。

表1 研究區遙感數據

2 研究方法

2.1 IS信息的提取

2.1.1 基于植被覆蓋度提取法

基于植被覆蓋度提取IS數據的方法最早由Carlson等[21]提出。該方法根據IS與植被覆蓋度的負相關關系提取IS信息，計算方法相對簡單，并且能夠快速獲取大范圍連續的IS信息，已經成為現階段應用最多的IS信息提取方法。該方法的表達式[13-14，21]為

I=1-F，

(1)

(2)

(3)

式中：I為IS指數，%；F為植被覆蓋度，%；N為歸一化植被指數；Nsoil為N在裸土或無植被區的值；Nveg為N在純植被區的值；NIR為近紅外波段的反射率，%；R為紅光波段反射率，%。由于影像存在不可避免的噪聲，因此在一般情況下，Nveg和Nsoil取一定置信度范圍內的最大值與最小值。在參考前人經驗[22]并多次試驗后，本文中選取累積分數95%～5%作為置信范圍，得到N的最大值和最小值。

在利用植被覆蓋度提取IS信息前，還需對研究區內的水體信息進行掩膜。徐涵秋[23]提出的改進型歸一化差異水體指數(modified normalized difference water index，MNDWI)能夠清楚地分出陰影和水體，是現階段提取效果最好的水體提取法，因此本文中使用該指數對水體信息進行提取、 掩膜。該指數的表達式為

(4)

式中：M為歸一化水體指數；G為綠光波段的反射率，%。

2.1.2 LSMA法

LSMA法通過計算混合像元中各端元的組合比例，將混合像元的光譜分解為多種端元光譜的線性組合[8]。該方法假設一個像元的反射率由影像拍攝瞬間各地物的純凈反射率(端元)按一定比例組合而成，即

(5)

式中：Rx為波段x的光譜反射率，%；n為端元個數；fi為端元i在像元中所占的比例，%；Ri,x為端元i在波段x的反射率，%；ex為殘差。

理論上， 端元的數量在一定程度上決定了LSMA法的擬合程度， 但這并不意味著端元種類越多，模型效果越優。Wu等[7]通過實際研究發現，選取4種端元時IS信息的反演效果最優。本文中依據前人經驗[7]對影像進行端元采集時選擇由玻璃、混凝土為代表的高反射率地物，傳統的舊式瓦片屋頂、瀝青路面為代表的低反射率地物，以及植被和土壤4種端元類型。除此之外，研究[2]發現，水體在LSMA法中的表現為低反射率。為了避免水體對研究結果的影響，本文中在提取純凈端元之前使用式(4)中的修正歸一化水體指數，對研究區內的水體進行剔除，將最終得到的高反射率蓋度圖和低反射率蓋度圖相加得到研究區IS蓋度圖[24]。

2.2 精度驗證

對影像分類提取結果驗證的最好方式是隨機提取采樣點進行實地驗證；但是，實地驗證成本較高，對過往年份展開實地驗證也不切實際，因此，本文中先利用地理信息系統(GIS)工具隨機生成驗證點，參照驗證點的空間位置對基準影像(2007、 2017年基準影像采用同期分辨率為1.07 m的Google影像，1996年的Google影像因效果不佳而采用原始影像)進行目視解譯，將解譯結果與二值化后的IS結果進行對比，對比結果通過混淆矩陣來反映提取精度。

2.3 IS的空間分布與演變

2.3.1 標準差橢圓分析法

為了探究研究區IS的空間變化情況，嘗試利用標準差橢圓分析法進行分析。通過分析標準差橢圓，可以從橢圓中心、 分布范圍、 方向和形狀等多個角度分析要素的空間分布特征及空間演變過程[25]。在地學領域，主要用橢圓的長軸半徑長度與扁率反映空間格局總體要素的集中密度，用長軸偏角反映分布格局的主導方向。

2.3.2 IS的擴展指標

IS擴展指標能夠全面反映研究區IS擴展的情況。 本文中選擇IS密度、 IS擴展速度這2個指標來反映各區域IS擴展趨勢[26]。 IS密度是指某區域中IS的面積占該區域總面積的比例； IS擴展速度是指在一段特定時間段內， 以年為單位的IS面積增長情況。

IS密度計算公式為

(6)

式中：ρa為區域a的IS密度，%；Aa為區域a的IS的面積， km2；Aa,tot為區域a的總面積， km2。

IS擴展速度計算公式為

(7)

式中：v為區域a的IS擴展速度， km2/a；t為時間間隔，a；Apri和Aend分別為區域a研究初期和末期的IS面積， km2。

3 結果與分析

3.1 精度驗證結果分析

將2種IS信息提取方法獲取的研究區2017年的IS影像選擇合適閾值進行二值化處理。利用生成的驗證點，分別將2幅二值化IS影像與基準影像疊在一起進行驗證，結果見圖1。基于植被覆蓋度提取法提取的IS總體分類精度為0.89，Kappa系數為0.53，而LSMA法得到的總體分類精度為0.95，Kappa系數為0.74。此外，在目視判斷過程中發現基于植被覆蓋度提取IS信息的方法很容易將大片裸地誤分為IS，而用LSMA法時該現象很少出現。由此可知，采用LSMA法在喀斯特地區具有較好的可行性，并且提取效果優于基于植被覆蓋度提取法。

在后續研究中選擇LSMA法繼續提取研究區1996、 2007年的IS信息， 并進行精度驗證， 結果如表2所示。 由表可知， 本文中提取的IS信息均能滿足研究要求， 可用于分析研究區IS空間格局演變。

表2 不同方法對貴陽市主城區不透水面的提取精度驗證

3.2 IS時空分布

利用GIS技術對得到的1996、 2007、 2017年的IS二值化圖(見圖2)進行面積統計， 由圖2可知， 研究區1996、 2007、 2017年的IS面積分別為50.82、 126.56、 234.08 km2， 21 a間貴陽市IS面積共增加183.26 km2，增長了近5倍，整體呈現出高速擴張趨勢，其中1996—2007年IS面積增加了75.76 km2，年平均增長量為6.89 km2/a；2008—2017年IS面積增加量達107.49 km2，年平均增長量為10.75 km2/a，前10 a的平均增長量明顯高于后10 a的。在2008—2017年中， 貴陽市內、 外部交通快速發展， 內部以環城高速和中環為代表的各等級道路竣工通車， 外部貴遵附線、 花安高速以及貴廣、 滬昆高速鐵路建成通車， 大力推動了貴陽市城市的發展。 除此之外， 觀山湖區作為緩解舊城區人口稠密、 交通擁堵的發展新區， 是城市擴展的重點方向， 在大量文化商貿區、 居住區、 公共基礎設施建設的同時， 也導致了IS的增長。

將提取的IS信息按行政區進行分區統計， 結果見表3。 由表可看出， 1996—2007年IS面積增長量最大的是花溪區， 增長了25.72 km2，增長量最小的是云巖區，僅為8.59 km2。通過觀察各區IS面積增長率發現，觀山湖區的增長率最高，達到了610.67%， 其次是花溪區、 烏當區， 增長率分別為377.13%、 369.85%，而南明區的增長率最低，僅為45.91%。2008—2017年IS面積增長最多的是觀山湖區，其次是花溪區，增長最少的是南明區。對比1996—2007年的數據，各區IS面積增長率明顯降低且區域間降低幅度差異較大，其中觀山湖區增長率下降超400%，但其增長率仍居各行政區值之首，南明區21 a間IS面積增長率都最低，但其一直處于平穩增長狀態。

表3 1996—2007年貴陽市主城區各區域不透水面面積的變化情況

從貴陽市各行政區IS面積的總體變化趨勢可以看出： 1996—2007年觀山湖區、 花溪區、 烏當區處于城市化擴張爆發期， 大量的IS井噴式出現； 南明區處于城市穩定發展期， IS格局基本穩定， 增長相對緩慢； 云巖區、 白云區則屬于典型的中速發展期， IS平穩擴展。 2008—2017年觀山湖區、 烏當區仍處于高速發展階段， 強度明顯減弱； 花溪區與云巖區、 白云區一樣進入城市中速發展期， IS繼續平穩擴張。

為了更直觀地展現各階段研究區IS面積變化的空間分布情況， 利用GIS工具計算研究區不同階段之間IS面積差值， 由此得到不同階段IS面積變化的圖斑， 結果見圖3。 由圖可以看出， 1996年以前研究區IS分布較集中， 主要分布在以南明區、 云巖區為主的市中心； 1996—2007年IS面積變化主要分布在各行政區中心周圍， 同時城市中心的IS面積增加也很顯著； 2008—2017年， IS的增加使得之前許多破碎的斑塊連接起來， IS開始出現成片聚集現象， 還有很多新增IS以“飛地”形式出現在各行政區中。

3.3 IS格局演變分析

3.3.1 IS空間變化特征

通過分析研究區1996、 2007、 2017年IS的標準差橢圓，探討該區域IS的整體變化趨勢， 如圖4所示。由圖可以看出：1996—2007年，標準差橢圓面積明顯增大且整體向南移動，表明此階段IS擴展重心正在向南傾斜；主軸方向偏角從9.87°減小到2.30°， 意味著該區域IS的分布趨勢由東北—西南方向往正北—正南方向發展； 標準差橢圓的長半軸由19.11 km增大到22.58 km， 短半軸由10.46 km增大到13.51 km，表明IS的空間分布在主軸方向(南北方向)和輔軸方向(東西方向)都呈現較弱的離散化現象。

圖3 貴陽市主城區不同時期不透水面的變化

2008—2017年研究區IS標準差橢圓范圍逐漸縮小且向東北方向移動,擴展重心逐漸向西北方向移動； 主軸方向的偏角從2007年的2.30°增至17.02°， 該階段貴陽市主城區IS的分布方向轉變為東北—西南方向， 與1996年相同且趨勢更加顯著； 長半軸由22.58 km減至22.2 km， 短半軸由13.51 km增至14.67 km， 可知IS的空間分布在主軸方向(南北方向)的離散化趨勢有所緩解，在輔軸方向(東西方向)離散化程度進一步加強。

(a)1996年 (b)1996—2007年 (c)1996—2017年

縱觀研究區IS的變化情況， 與研究區地形密切相關。 貴陽市屬于典型的喀斯特山地城市， 地表景觀破碎化嚴重。 南明區、 云巖區是貴陽市的老城區， 大部分區域位于貴陽—中曹司盆地內， 四面環山， 受地形嚴重限制， 迫使其城市范圍不斷尋求向外擴展， 進而致使IS向其他行政區擴張， 但內部填充式擴張仍是其擴展的重要方式， 這種狀況一直保持到2000年。 2000年后， 大量的交通干線開始建設， 以埡口和河谷為走廊逐步向地勢平坦的花溪區、 觀山湖區延伸[27]。 2001—2006年間，花溪區成為貴陽市交通的主要擴展地區， 交通建設量占同期工程建設總量的45%。 從2006 年開始， 貴陽市城市建設發展的主方向轉移至觀山湖區， 其次是花溪區， 這2個方向擴展的比例達到同期工程建設總量的80%[27]， 這些變化情況都與研究區IS標準差橢圓的變化相呼應。

3.3.2 IS擴展趨勢

1996—2017年貴陽市主城區各區域IS的擴展指標如圖5所示。 從圖5(a)中可以看出， 1996、 2007、 2017年研究區的IS密度分別為2.14%、 5.32%、 9.84%， 22 a間增加了7.7%。各行政區的IS密度變化存在明顯差異，其中云巖區的IS密度在2個研究階段都是最高的，1996、 2007、 2017年的IS密度分別為12.23%、 21.62%、 36.33%，而烏當區在整個研究階段IS密度最低， 增長最緩慢。從圖5(b)中可以看出：1996—2007年研究區IS的擴展速度為6.89 km2/a，各行政區IS擴展極不均衡，其中花溪區的擴展最快，擴展速度為2.34 km2/a，是各行政區擴展速度平均值的2倍，烏當區的IS擴展最慢，擴展速度僅為0.67 km2/a；2008—2017年，研究區IS擴展速度整體上升至10.75 km2/a，花溪區雖然保持較高的擴展速度，但是擴展速度呈倒退趨勢，其余各區擴展速度明顯提升，特別是觀山湖區擴展速度上升最多，達到了2.65 km2/a，擴展速度最低的南明區也較前一階段提升了0.4 km2/a。

(a)IS密度

4 結論

1)通過對LSMA法和基于植被覆蓋度提取法得到的貴陽市1996—2017年IS信息進行精度驗證可知，LSMA法的提取結果明顯優于基于植被覆蓋度提取法的，說明在喀斯特地貌環境下使用LSMA法提取IS信息是適用的，并且效果更優。

2)1996—2017年貴陽市IS面積從50.82 km2增加到210.12 km2，該市的IS擴展可分為3個階段：前期(1996年)IS分布比較集中，主要分布在城市中心區域；中期(1996—2007年)主要呈現出以各行政區為發散中心向四周擴張，同時伴隨著城市中心區域的內部填充；后期(2008—2017年)IS的擴充得到進一步發展，擴展方式主要呈現出多中心發散的特點，并出現IS“飛地”，其中觀山湖區擴張尤為突出。

3)研究區IS面積的標準差橢圓偏角從1996年的9.87°增加到2017年的17.02°，表明研究區IS的分布趨勢表現為東北—西南方向，并且隨著城市的發展，該分布態勢得到進一步加強。

4)本文中將LSMA法應用于喀斯特地區IS信息的提取研究，與基于植被覆蓋度提取法提取的IS信息相比，LSMA法更適用于喀斯特地區。本文中還提取了研究區多個時間節點的IS信息，對IS信息的時空演變進行了進一步探討，研究結果為喀斯特地區IS研究和貴陽市后期城市發展規劃、 生態文明城市建設提供了參考。對提取到的IS信息開展相關應用研究將是后期研究的重點。