胡書琦，吳文淵，陳喜芬，苗 瑩，許林霞，胡程浩

(1. 杭州師范大學(xué)遙感與地球科學(xué)研究院，浙江 杭州 311121； 2. 浙江省城市濕地與區(qū)域變化研究重點(diǎn)實(shí)驗(yàn)室，浙江 杭州 311121)

0 引言

在過(guò)去的幾十年中，城市擴(kuò)張以前所未有的速度發(fā)生，全球城市地區(qū)的增長(zhǎng)速度約為城市人口的兩倍[1-2].城市擴(kuò)張是土地利用/覆蓋變化的一種主要形式[3]，影響著本地和全球的生態(tài)系統(tǒng)[4]，不僅帶來(lái)了社會(huì)經(jīng)濟(jì)問(wèn)題的變化[5]，而且是重要的地理空間變化過(guò)程[6].因此對(duì)城市化的地理空間過(guò)程的研究尤為重要.

遙感(remote sensing，RS)提供了空間一致的數(shù)據(jù)集，可以覆蓋具有高空間細(xì)節(jié)和高時(shí)間頻率的大區(qū)域，因其在揭示城市土地?cái)U(kuò)張和描述時(shí)空格局方面的有效性而受到重視[7-10]，并與地理信息系統(tǒng)(geographic information system，GIS)相結(jié)合，廣泛用于城市擴(kuò)張研究[11],得到了許多高質(zhì)量的結(jié)果[12-16].

然而，由于20世紀(jì)80年代的衛(wèi)星遙感數(shù)據(jù)分辨率不高，使得城市的歷史變遷分析和格局變化研究等工作受到影像數(shù)據(jù)限制.隨著美國(guó)CORONA Keyhole衛(wèi)星數(shù)據(jù)的開始應(yīng)用，這一缺陷得以彌補(bǔ).CORONA Keyhole是美國(guó)偵察衛(wèi)星，從20世紀(jì)60年代開始提供高分辨率的全色數(shù)據(jù)[17]，這些數(shù)據(jù)在1995年被解密.歷史的具有良好空間分辨率的Keyhole影像采集以及廣闊的空間覆蓋范圍使它們廣泛應(yīng)用于各種科學(xué)研究中.

但在城市擴(kuò)張的研究中，由于獲取的5 m高空間分辨率的Keyhole衛(wèi)星數(shù)據(jù)僅有一個(gè)全色波段，不能在高分辨率的同時(shí)兼有豐富的光譜信息，所以無(wú)法發(fā)揮很大價(jià)值.這就需要通過(guò)多源遙感數(shù)據(jù)融合的方法來(lái)提高融合影像的空間分辨率，提升圖像的對(duì)比度和清晰度，增強(qiáng)影像的解譯和動(dòng)態(tài)監(jiān)測(cè)能力，從而更好地發(fā)揮鎖眼歷史衛(wèi)星的使用價(jià)值.這一方法可以將城市變遷的研究推進(jìn)至此前無(wú)法企及的20世紀(jì)60年代，并為區(qū)域土地利用與城市變遷提供更好的研究基礎(chǔ)，使其能夠?yàn)槌鞘邪l(fā)展與規(guī)劃提供更具價(jià)值的參考依據(jù).

本研究使用Keyhole衛(wèi)星數(shù)據(jù)和Landsat衛(wèi)星數(shù)據(jù)融合的方法，在影像融合與協(xié)同的基礎(chǔ)上對(duì)研究區(qū)進(jìn)行監(jiān)督分類及城市變遷與土地利用變化的分析，并利用GIS空間分析方法研究城市分布格局的變化，揭示古城山西大同的城市發(fā)展特點(diǎn).研究的內(nèi)容主要為：1)研究Keyhole全色影像(5 m空間分辨率)與Landsat多光譜數(shù)據(jù)(30 m空間分辨率)的融合；2)結(jié)合Spot遙感衛(wèi)星及高分一號(hào)衛(wèi)星探索過(guò)去40年中城市擴(kuò)張的時(shí)空特征；3)分析1980—2019年由于城市化形成的山西大同城區(qū)土地利用和土地覆蓋變化情況.

1 研究區(qū)與數(shù)據(jù)

1.1 研究區(qū)概況

大同市城區(qū)研究區(qū)(113°11′—113°19′E，40°01′—40°08′N)隸屬于山西省大同市，位于大同盆地中心，是大同市的政治、經(jīng)濟(jì)和文化中心，總面積46 km2(圖1).城區(qū)西北略高，東南稍低，區(qū)內(nèi)無(wú)高山丘陵，基本為平原地帶，西北東三面環(huán)山；屬溫帶大陸性季風(fēng)氣候區(qū)，有嚴(yán)寒無(wú)酷暑，四季分明，多風(fēng)沙半干旱的氣候.

a：山西省大同的地理位置；b：研究區(qū)位置.

1.2 遙感影像數(shù)據(jù)

多時(shí)相衛(wèi)星圖像可有效監(jiān)測(cè)、測(cè)繪和評(píng)估城市擴(kuò)展與形態(tài)景觀特征之間相互作用的時(shí)空特征[18-19].本研究收集了1980、2000和2019年的衛(wèi)星遙感影像(表1)，以評(píng)估1980—2019年城市擴(kuò)張的時(shí)空特征.

表1 本研究中使用的數(shù)據(jù)Tab.1 The data used in this study

2 研究方法

2.1 圖像融合與協(xié)同

遙感圖像融合的目的在于集成或整合多源圖像中的優(yōu)勢(shì)信息和互補(bǔ)信息，利用優(yōu)勢(shì)互補(bǔ)的數(shù)據(jù)來(lái)提高圖像的信息可用度，同時(shí)也提高對(duì)研究對(duì)象解譯的可靠性.本研究首先使用地理信息系統(tǒng)通過(guò)精心選擇的地面控制點(diǎn)對(duì)圖像進(jìn)行幾何校正，將1980年的Keyhole全色影像和Landsat多光譜影像進(jìn)行精確的校正和配準(zhǔn)，使兩景影像中的地物精確重合；采用Brovey融合法、HSV(hue,saturation,value)融合法、Gram-Schmidt融合法和影像協(xié)同的方法進(jìn)行融合，選擇最適合Keyhole與Landsat TM影像的融合分辨率以及融合方法，作為后續(xù)目視解譯及監(jiān)督分類的輸入數(shù)據(jù).

采用數(shù)據(jù)配準(zhǔn)與數(shù)據(jù)分辨率重采樣等方法，將Landsat的7個(gè)多光譜影像波段采樣成5 m分辨率，并和Keyhole的1個(gè)5 m分辨率的全色波段生成Keyhole與Landsat的5 m空間分辨率協(xié)同影像.影像共有8個(gè)波段，綜合了Keyhole影像與Landsat影像的優(yōu)勢(shì)，用于后續(xù)的土地利用變化監(jiān)測(cè).

2.2 土地利用和土地覆蓋分類

根據(jù)土地利用分類標(biāo)準(zhǔn)及研究需要，將土地類型分為耕地、林地、水體、城鎮(zhèn)建設(shè)用地和草地5類.采用廣泛使用的監(jiān)督分類方法Maxi-mum Likelihood[20]來(lái)檢測(cè)土地覆蓋類型，并利用分類后處理工具對(duì)分類結(jié)果進(jìn)行精度評(píng)價(jià)，以土地分類專題數(shù)據(jù)為基礎(chǔ)，提取矢量數(shù)據(jù)，獲得研究區(qū)域城市擴(kuò)張系列數(shù)據(jù).

分別基于Keyhole全色影像、Landsat多光譜影像以及兩者的融合影像和協(xié)同影像進(jìn)行監(jiān)督分類與精度評(píng)價(jià)，得到1980年的城區(qū)土地利用分類圖，在此基礎(chǔ)上提取建設(shè)用地矢量數(shù)據(jù).然后將2000年的Spot數(shù)據(jù)與2019年的高分?jǐn)?shù)據(jù)進(jìn)行土地利用分類后提取建設(shè)用地矢量數(shù)據(jù).

2.3 城市擴(kuò)張情況分析

城市擴(kuò)張速度是指某一時(shí)期內(nèi)城市的擴(kuò)張面積，在該時(shí)段內(nèi)的年均增長(zhǎng)量和絕對(duì)擴(kuò)張能直觀反映城市擴(kuò)張的變化程度，因此目前被廣泛應(yīng)用于城市擴(kuò)張?zhí)卣鞣治鯷21-23]，計(jì)算公式如下：

式中：V表示建設(shè)用地?cái)U(kuò)張速度(km2·a-1)；Sa表示研究時(shí)段初始年城鎮(zhèn)建設(shè)用地面積(km2)；Sb表示研究時(shí)段末年城鎮(zhèn)建設(shè)用地面積(km2)；T表示研究時(shí)間(a).

城市擴(kuò)張強(qiáng)度表示某一時(shí)期內(nèi)研究區(qū)城鎮(zhèn)建設(shè)用地的面積增長(zhǎng)比率即相對(duì)擴(kuò)張，用以表征單位時(shí)間內(nèi)城市擴(kuò)張的強(qiáng)弱和快慢，公式如下：

式中：R表示城市擴(kuò)張強(qiáng)度(%)；Sa表示研究期前建設(shè)用地面積(km2)；Sb表示研究期末建設(shè)用地面積(km2)，T表示研究時(shí)間(a).城市擴(kuò)張強(qiáng)度可劃分為5個(gè)等級(jí)，即：R>1.92，高速擴(kuò)張;1.05

城市空間重心轉(zhuǎn)移中的“重心”是一個(gè)物理概念，旨在表示某個(gè)區(qū)域在某些空間上的平衡點(diǎn)，通常用于指示城市演變、土地利用類型變化[25].通過(guò)不同時(shí)期城市重心遷移情況的對(duì)比，可分析城市空間變化規(guī)律、城市擴(kuò)張空間變化的方向及走勢(shì).重心模型是描述地理對(duì)象空間分布的一個(gè)重要指標(biāo)，能夠清晰客觀地反映區(qū)域地理對(duì)象在時(shí)間和空間上的變化軌跡[26].城市重心坐標(biāo)的計(jì)算公式如下：

式中：Xt、Yt分別表示第t年城市重心的經(jīng)緯度坐標(biāo)；Sti表示第t年第i個(gè)矢量圖斑的幾何面積(km2)；Xi、Yi分別表示第i個(gè)矢量圖斑的幾何中心坐標(biāo).

城市空間重心轉(zhuǎn)移距離計(jì)算公式如下：

式中：x、y分別表示某一時(shí)期地理單元重心的橫、縱坐標(biāo)；Dab表示從第a年到第b年城市空間重心轉(zhuǎn)移距離(km).

本研究中，將3個(gè)年份的建設(shè)用地矢量數(shù)據(jù)進(jìn)行轉(zhuǎn)移矩陣的計(jì)算，得到不同類型用地的轉(zhuǎn)變趨勢(shì)及面積百分比；然后，將建設(shè)用地?cái)?shù)據(jù)疊加并分別計(jì)算其城市重心、城市擴(kuò)張速度和城市擴(kuò)張強(qiáng)度，得到城市擴(kuò)張時(shí)間及空間上的分布規(guī)律，分析山西大同城區(qū)城市建設(shè)用地及重心的動(dòng)態(tài)變化.

3 結(jié)果

3.1 影像融合結(jié)果

比較Keyhole全色影像、Landsat多光譜影像以及兩者的融合影像結(jié)果(圖2)發(fā)現(xiàn)，Keyhole影像的紋理信息最為清晰，能夠較好地分辨道路和建設(shè)用地邊界；Landsat多光譜影像顏色信息豐富，可以區(qū)分河流與耕地；兩者的融合影像既有良好的紋理信息，又能夠通過(guò)顏色更好地辨別地物.

a：Brovey融合影像；b：HSV融合影像；c：Gram-Schmidt融合影像；d：協(xié)同影像；e：多光譜影像；f：全色影像.

本文數(shù)據(jù)融合之前經(jīng)過(guò)精確配準(zhǔn)，使用Band7、Band4、Band3多光譜波段進(jìn)行的彩色合成接近自然色彩，既能夠較為真實(shí)地反映地物原本顏色，又在一定程度上突出顯示了植被、水體與城市建設(shè)用地，其中植被為亮綠色，水體為深藍(lán)色，城鎮(zhèn)建設(shè)用地為紅色，地物類型豐富.

依據(jù)圖2結(jié)果，從整體上看，融合影像提高了多光譜圖像的分辨率并且顏色、色調(diào)與原影像相似，說(shuō)明融合結(jié)果較好地保留了圖像的光譜特征.從紋理細(xì)節(jié)上看，HSV變換融合后的影像較原始多光譜圖像紋理內(nèi)容更加豐富，清晰度和空間細(xì)節(jié)表達(dá)能力有所提高，尤其是在原始全色影像上無(wú)法識(shí)別的田地與水體等地物的目視效果增強(qiáng).

圖像對(duì)比后還可以發(fā)現(xiàn)，1980年大同市城市東面的植被在鎖眼全色圖像中的呈現(xiàn)與湖泊幾近相同，難以用目視解譯的方法將兩者加以區(qū)分.而通過(guò)融合處理后，水體和耕地分別擁有了地物的自然顏色，對(duì)比結(jié)果見圖3、圖4.

a：全色影像中的水體；b：多光譜影像中的水體；c：融合影像中的水體.

a：全色影像中的耕地；b：多光譜影像中的耕地；c：融合影像中的耕地.

城區(qū)部分的融合影像(圖5)顯示，Keyhole數(shù)據(jù)可以清晰地勾勒出建筑物的形態(tài)邊界，而Landsat多光譜數(shù)據(jù)已經(jīng)失去可辨識(shí)性，除了大致輪廓及一些顏色信息外無(wú)法為目視解譯提供有效支持.融合后的影像則可以看清建筑及道路邊界，為后續(xù)的分類及城市擴(kuò)展研究提供基礎(chǔ).

a：全色影像中的城區(qū)；b：多光譜影像中的城區(qū)；c：融合影像中的城區(qū).

分別對(duì)Keyhole、Landsat、3種融合影像以及協(xié)同影像進(jìn)行監(jiān)督分類，選用同一批訓(xùn)練樣本并對(duì)其分類結(jié)果進(jìn)行精度評(píng)價(jià)，結(jié)果見表2和圖6.依據(jù)分類結(jié)果，協(xié)同影像的精度最高，HSV融合影像次之，兩者的準(zhǔn)確性均高于90%.Gram-Schmidt融合影像與Landsat影像分類結(jié)果的精度相似，準(zhǔn)確性達(dá)85%以上.Brovey融合影像分類精度略低，而Keyhole分類結(jié)果較差，準(zhǔn)確性僅為30.81%.從上述目視解譯及分類結(jié)果評(píng)價(jià)中不難看出HSV融合影像與協(xié)同影像的優(yōu)越性，融合與協(xié)同影像使得城市變遷研究中的監(jiān)督分類效果更好.

表2 融合前后影像分類精度評(píng)價(jià)Tab.2 Accuracy evaluation of image classification before and after fusion

a：Keyhole影像的分類結(jié)果；b：Landsat影像的分類結(jié)果；c：Brovey融合影像的分類結(jié)果；d：Gram-Schmidt融合影像的分類結(jié)果；e：HSV融合影像的分類結(jié)果；f：協(xié)同影像的分類結(jié)果.

同時(shí)，對(duì)2000年與2019年影像進(jìn)行土地利用分類，得到了土地利用分類圖(圖7).2000年(圖7b)和2019年(圖7c)影像分類的準(zhǔn)確性分別為90.24%和95.05%，Kappa系數(shù)分別為0.86和0.88.為了做多期分類結(jié)果對(duì)比，將1980年分類精度最好的協(xié)同影像分類結(jié)果(圖7a)也放入圖中.

a：1980年土地利用圖；b：2000年土地利用圖；c：2019年土地利用圖.

3.2 城市的景觀變化

根據(jù)分類結(jié)果，完成了1980—2019年的土地利用/土地覆蓋轉(zhuǎn)換矩陣(表3).結(jié)合圖7比較相應(yīng)的結(jié)果，發(fā)現(xiàn)景觀變化較大，城市區(qū)域已有擴(kuò)展，而農(nóng)業(yè)用地部分則減少.表3結(jié)果同時(shí)揭示了土地利用/覆蓋轉(zhuǎn)換的關(guān)系，例如：城市土地顯示出較大的增長(zhǎng)趨勢(shì)，建設(shè)用地由1980年的60.60 km2變化為2019年的127.67 km2，增長(zhǎng)110.67%，主要是耕地和林地轉(zhuǎn)變?yōu)榫幼∮玫?；林地的面積也有所增長(zhǎng)，1980年持有50.86 km2，2019年增長(zhǎng)至68.89 km2，其中部分耕地轉(zhuǎn)換為林地與同一時(shí)期的退耕還林措施有關(guān)；水體的變化是各類用地變化中最微弱的，1980年和2019年占地面積分別為7.84和7.92 km2，其增長(zhǎng)率僅為1.02%；耕地與草地則呈現(xiàn)出較大幅度的負(fù)增長(zhǎng)，其中耕地面積由1980年的122.44 km2削減至2019年的72.45 km2，而草地從35.57 km2銳減至0.38 km2，變化率高達(dá)-98.93%.此外，城市建設(shè)用地的擴(kuò)張主要表現(xiàn)為東西方向的擴(kuò)展.在南北向上，城市南邊建設(shè)用地增加的面積比北邊的更多，尤其集中在城區(qū)的西南部，增加部分沿著城市河流與主干道路分布，由此可以看出交通驅(qū)動(dòng)因素對(duì)城市發(fā)展的重要性.

表3 1980—2019年土地利用/土地覆蓋變化轉(zhuǎn)換矩陣

3.3 城市擴(kuò)張情況

3.3.1 城市擴(kuò)張的時(shí)間屬性

表4 1980—2019年大同城區(qū)擴(kuò)張的數(shù)量特征Tab.4 Quantitative characteristics of urban expansion in Datong from 1980 to 2019

由表4可以看出，從1980年到2019年，大同市主城區(qū)城市建設(shè)用地面積總共增加67.07 km2，擴(kuò)張速度為1.92 km2·a-1，擴(kuò)張強(qiáng)度達(dá)3.16%(圖8).其中，建設(shè)用地面積從1980年到2000年增加了12.45 km2，擴(kuò)大了1.21倍，擴(kuò)張速度為0.78 km2·a-1；從2000年到2019年增加了54.62 km2，擴(kuò)大了1.73倍，擴(kuò)張速度為2.87 km2·a-1.第二階段的擴(kuò)張強(qiáng)度較第一階段有明顯增長(zhǎng)，兩個(gè)時(shí)期的城市發(fā)展呈現(xiàn)出“緩慢—加速”的擴(kuò)張?zhí)卣?大同正處于快速擴(kuò)張期，城市的規(guī)模不斷擴(kuò)大，城市發(fā)展的新格局正在進(jìn)一步建立，與之密切相關(guān)的是現(xiàn)行的經(jīng)濟(jì)政策和發(fā)展規(guī)劃.

a：1980年主城區(qū)建設(shè)用地；b：2000年主城區(qū)建設(shè)用地；c：2019年主城區(qū)建設(shè)用地.

3.3.2 城市擴(kuò)張強(qiáng)度的空間分布

城市建設(shè)用地不同年份的中心點(diǎn)坐標(biāo)是通過(guò)空間分析的方法計(jì)算得出的，同時(shí)得到了不同時(shí)期城市重心的水平坐標(biāo)、垂直坐標(biāo)和移動(dòng)的距離.圖9顯示的是3個(gè)年份建設(shè)用地的矢量疊加，可以看到主城區(qū)即原古城墻所在區(qū)域建設(shè)用地沒有太大變化，密度與面積都維持在一定水平.擴(kuò)張強(qiáng)度最大的是主城區(qū)東側(cè)區(qū)域.東西差異的主要原因是東邊地形較為平坦，多耕地，而西部以山區(qū)林地為主，不利于建設(shè)用地的擴(kuò)張.1980年至2000年，大同市建設(shè)用地空間重心向西北方向移動(dòng)了2.90 km.2000年至2019年，建設(shè)用地空間重心則向東南移動(dòng)1.39 km.總體來(lái)看，大同市建設(shè)用地空間重心向西北移動(dòng)1.51 km，呈現(xiàn)西北緩慢延伸、東南跨越發(fā)展的擴(kuò)張?zhí)卣?

圖9 1980—2019年大同城區(qū)建設(shè)用地及重心動(dòng)態(tài)變化Fig.9 Dynamic changes of construction land and center of gravity in Datong from 1980 to 2019

4 結(jié)論與討論

本研究基于1980年的Keyhole衛(wèi)星和Landsat影像融合數(shù)據(jù)、2000年的Spot遙感影像及2019年的高分衛(wèi)星遙感數(shù)據(jù)，提取了大同市主城區(qū)建設(shè)用地信息，在計(jì)算擴(kuò)張速度、擴(kuò)張強(qiáng)度等基礎(chǔ)上，分析了研究區(qū)的時(shí)空擴(kuò)張?zhí)卣鳎瑸榇笸谐鞘袛U(kuò)張規(guī)律及驅(qū)動(dòng)力研究提供了一定的影像基礎(chǔ)和分析參考.結(jié)果表明:山西大同市城區(qū)40 年的變化較大，1980年的遙感影像顯示為四方形的古城，而2019年卻擴(kuò)張為周邊城市群.研究區(qū)土地利用類型以耕地、林地為主，建設(shè)用地持續(xù)增長(zhǎng).1980—2000年及2000—2019年兩個(gè)時(shí)間段內(nèi)的城市擴(kuò)張量不同，但擴(kuò)張過(guò)程明顯，總體呈現(xiàn)出階段性擴(kuò)張、區(qū)域性增長(zhǎng)的特點(diǎn).城市內(nèi)各區(qū)域擴(kuò)張模式不同，中心城區(qū)主要以填充式擴(kuò)張為主，城市化水平趨于成熟，城市擴(kuò)張緩慢.在周邊城區(qū)中，農(nóng)業(yè)用地大幅轉(zhuǎn)化為城市用地，是大同市城市擴(kuò)張的主要模式.城市形態(tài)在西北角較為穩(wěn)定，東南面則變化較大，城市邊界從曲折向整齊過(guò)渡.整個(gè)城區(qū)重心以0.04 km·a-1的速度向西北方向轉(zhuǎn)移.

另外，依據(jù)1980—2019年城市建設(shè)用地及重心的動(dòng)態(tài)變化，大同的城市擴(kuò)張與中國(guó)其他大城市相似，其發(fā)展大致沿城市中心的主要交通線呈放射狀趨勢(shì)[27]，且在城市化邊緣地區(qū)表現(xiàn)出區(qū)域擴(kuò)張的特點(diǎn)[8].

本研究結(jié)果也表明，1980年的Keyhole衛(wèi)星和Landsat影像融合數(shù)據(jù)在1986年Spot1發(fā)射之前提供了高空間分辨率的影像，在城市擴(kuò)張、濕地監(jiān)測(cè)等方面可以有較好的應(yīng)用.