基于遙感的遷西縣土地利用景觀格局變化及脆弱性分析

曹應舉 邢英梅 曹 琳 王 娜 郝雪景

(河北水文工程地質勘察院 石家莊 050000)

景觀是由多個生態系統構成的異質性地域或不同土地覆被方式的鑲嵌體[1]；景觀格局是大小和形狀各異的不同景觀要素在空間上的分布與組合，是自然、生態過程和人為因素相互作用的結果[2-3]。景觀格局變化能夠體現空間上的斑塊和景觀變化，與區域生態環境變化密切相關[4]。因此，對景觀格局動態變化的研究成為當今學者研究區域生態環境變化的必要前提。多年來，不少學者通過景觀格局的數量和空間變化分析景觀的演變特征，通過結合多種生態敏感性因子構建景觀脆弱性指數來評價景觀生態的好壞[5]。例如，曾光[6]等借助RS和GIS技術，對山西省1975～2014年間四個時間段的多源遙感數據采用景觀指數的方法，進行濕地景觀格局變化特征分析。結果表明，40年來山西省濕地面積減少、退化速率加快，濕地景觀的破碎度指數先增后減，山西省濕地景觀格局的變化歸根于人類活動的影響。張榮天[7]等借助Landsat TM/ETM+數據，以鎮江市為研究案例地，運用CA-Markov模型、動態度模型、景觀指數以及PCA方法，對1994～2012年研究區景觀格局演變特征、驅動機制以及變化趨勢進行研究，通過定量化分析了解了區域土地利用景觀格局演變的基本特征。董雅雯[8]等基于克里金插值法，通過構建的景觀脆弱度指數對三亞市景觀脆弱度進行評價，結果以不同等級景觀脆弱區的分布特征和數量大小反映出三亞市整體受外界干擾的適應程度大小。李慧穎[9]等以吉林省長吉示范區為例，利用1976～2015年多源遙感數據的景觀分類結果，對研究區森林景觀的變化特征以及驅動因素進行分析，結果表明人類活動是導致研究區森林景觀破碎化程度加劇的主要原因。

總結前人的研究成果，發現景觀格局變化的研究大多是對人口、經濟和氣候等外在變化驅動力的分析，通過景觀格局變化進行景觀脆弱性評價的研究還略顯不足?；诖?，本研究以遷西縣為例，基于土地利用動態度模型、景觀格局指數分析方法，對2009年、2013年以及2017年遙感影像景觀分類結果進行空間疊置對比分析，構建景觀脆弱性指數對三期數據的景觀脆弱性進行綜合評價，旨在進一步分析研究區景觀變化的動態特征以及發展方向。

1 材料與研究方法

1.1 研究區概況

遷西縣位于河北省東北部，地處東經118°6′～118°37′，北緯39°57′～40°27′，東西橫距39 km，南北縱距51 km，四季分明，屬于暖溫帶大陸性半濕潤的季風氣候。遷西縣是唐山市山區較多的縣，地勢南北高中間低，海拔70～831.3 m，境內有灤河、灑河、長河、清河等7條河流分布其中，土地利用類型多樣，林地覆蓋率達68 %，蘊含豐富的水資源、林地資源和礦產資源，其地貌呈典型的低山丘陵景觀，因其獨特的地形地貌特征被評選為國家級生態示范縣。該縣屬環渤海經濟圈，毗鄰北京、天津、秦皇島、承德等市，優越的地理位置和便捷的交通使得遷西縣經濟得以快速發展(圖1)。

圖1 研究區位置及高程

1.2 研究數據

研究所采用的數據來源于USGS美國地質調查局網站的Landsat TM/OLI影像，時段分別為2009年、2013年和2017年，獲取月份分別為8月、7月和6月，云量覆蓋率均小于1%且研究范圍內無云量覆蓋，影像投影分度帶為UTM Zone 50N，基準面為WGS-84。其他輔助數據包括遷西縣矢量邊界確定研究區精確范圍，Google Earth高分辨率歷史影像用于目視解譯選取高精度驗證樣本，ASTER GDEM 30 m分辨率高程數據，統計數據主要來源于2010～2016年《唐山統計年鑒》，用于對研究區土地利用分類結果進行對照分析。

為了精確反演地表輻射亮度值和反射率值，在ENVI5.1軟件中對遙感影像進行輻射定標和大氣校正，Landsat8-OLI影像的大氣校正結果與全色波段利用Gram-Schmidt變換方法進行圖像融合，再利用遷西縣矢量文件進行邊界裁剪，最終得到研究區范圍。后期為達到增強和提取有用信息的目的，對研究區影像進行了分段線性拉伸和主成分分析，主成分分析保留了最大主成分波段。

1.3 景觀分類及精度評價

為了系統分析遷西縣景觀格局動態特征，參照《土地利用現狀分類》(GB/T 21010-2017)，結合遷西縣景觀格局的實際情況，將研究區景觀類型分為：居民地、工礦用地、林地、水體、耕地、園地和其他用地七大景觀類型。研究以eCognition軟件為平臺，運用面向對象和決策樹相結合的分類方法，通過建立分類規則，區分開特征參數差異較大的地類，之后從易到難依次進行分類，根據制定的分類標準和分類規則，綜合考慮影像的光譜、紋理、形狀等特征，對遷西縣的3期遙感影像進行多尺度分割與分類，最后對分類結果進行后處理(圖2)。

圖2 2009～2017年遷西縣景觀分類結果

借助谷歌地球高分辨率歷史影像，隨機選取足夠數量且分布均勻的樣本點對2009～2017年分類結果逐一進行精度評價。精度驗證結果表明，遷西縣土地利用總體分類精度分別為91.3 %(2009年)、92.0 %(2013年)和90.6 %(2017年)，均達到最低允許判讀的精度，滿足研究需要。

1.4 研究方法

(1)土地利用動態度模型

土地利用動態度反映的是研究區在一定時間范圍內，土地利用類型變化的數量和速度[10]。應用土地利用動態度分析土地利用類型的動態變化，可以真實反映區域土地利用類型的變化快慢程度[11]。根據研究對象，可將土地利用動態度分為單一土地利用動態度和綜合土地利用動態度。

單一土地利用動態度表達的是在特定時期內某一土地利用類型變化的程度和速率[12]，其計算公式如式1所示。

(1)

式中：K為研究時段內某一土地利用類型的動態度，ua和ub分別為該研究時段初期和末期的土地利用類型面積，t為研究時間間隔，土地利用動態度K反映了研究區在t時間間隔內的土地利用變化速率。

綜合土地利用動態度表達的是土地利用變化的劇烈程度[13]，其計算公式如式2所示。

(2)

式中：LC為研究區土地利用類型的綜合土地利用動態度，Si為研究初期第i類土地利用類型的面積，ΔSi-j為研究時段內第i類土地利用類型轉為非i類土地利用類型面積的絕對值，t為研究時間間隔。

(2)景觀格局指數的選擇

景觀格局指數能夠高度濃縮景觀空間格局信息，并且能夠反映景觀結構組成和空間配置等方面的特征，是研究景觀格局時空變化的關鍵數據[14]。根據研究需要同時結合研究區特點，篩選出適合本研究的多種景觀指數。在斑塊類型水平上選取斑塊數量、最大斑塊指數、邊緣密度、斑塊密度、景觀形狀指數以及分散和并列指數；在景觀水平上選取蔓延度指數、景觀連接度指數、聚合指數和香農多樣性指數。以期對研究區的景觀破碎程度進行簡單度量，同時計算出豐度密度指數、景觀均勻度指數、面積加權平均分維數、斑塊類型面積4個重要的景觀指數，用以對研究區脆弱性進行時空變化分析。研究所用的指數通過景觀結構分析軟件Fragstats 4.2計算得到。

(3)景觀干擾度指數

景觀干擾度指數由景觀類型破碎度、分維數倒數和優勢度共同決定，其計算公式如式3所示。

Si=aFNi+bFDi+cDOi

(3)

式中：景觀類型破碎度FNi，分維數倒數FDi，優勢度DOi，a,b,c分別為各個指數對景觀干擾度指數的貢獻量大小，本研究三個參數分別賦值為0.5、0.3和0.2。

(4)景觀敏感性指數

景觀敏感性指數是通過景觀干擾度指數和景觀類型易損度指數的乘積相疊加，用來表征景觀受外界干擾敏感程度的指標[16]，景觀敏感性指數的計算公式如式4所示。

(4)

式中：LSI為景觀敏感性指數，i為景觀類型，n為景觀類型數目，Si為景觀干擾度指數，Vi為景觀類型易損度指數。參考他人研究成果[17]，將各景觀類型按易損程度由低到高分別賦值：居民地為1，工礦用地為2，林地為3，園地為4，耕地為5，水體為6，其他用地為7，最后進行標準化處理得到每個景觀類型的易損度指數。

(5)景觀適應性指數

景觀脆弱性不僅與景觀敏感性有關，還與景觀的適應和恢復能力有關。一般情況下，景觀多樣性越大，結構分布越均勻，其景觀適應性越強。本文借鑒孫才志[18]等的研究，引入景觀適應性指數用來描述景觀的抗外界干擾及適應能力，景觀適應性指數LAI計算公式如式5所示。

LAI=PRD×SHDI×SHEI

(5)

式中：PRD為斑塊豐度密度指數，SHDI為香農多樣性指數，SHEI為香農均勻性指數。

(6)景觀格局脆弱性評價指數

為了直觀反映研究區2009～2017年三個時期的景觀格局空間變化對生態環境產生的負面影響，利用Fragstats 4.2軟件計算出的景觀指數構建景觀脆弱性指數。

系統的脆弱性通過系統面對外界干擾的敏感程度與適應能力來體現，根據引入景觀敏感性指數和景觀適應性指數構建可以評價景觀生態安全的景觀脆弱性指數LVI[15]，其計算公式如式6所示。

LVI=LSI×(1-LAI)

(6)

式中：LSI表示景觀敏感性指數，LAI表示景觀適應性指數。

2 結果與分析

2.1 景觀格局變化特征分析

基于獲得的土地利用分類結果，采用GIS空間分析功能統計各時期土地利用類型的變化面積，根據公式(1)和公式(2)，可以得出2009～2013年、2013～2017年遷西縣單一土地利用動態度和綜合土地利用動態度(表1)。

表1 2009～2017年遷西縣土地利用結構及其變化

從遷西縣土地利用動態變化表可以看出：2009～2013年間，遷西縣土地利用空間變化面積最大的地類是工礦用地，變化總量為2489.37 hm2，其次是林地，變化總量為-2375.29 hm2，其他用地面積變化最小，變化總量為-44.96 hm2。在單一動態度方面，工礦用地動態度最大，其數值為14.44%，其次分別是園地、水體、居民地、耕地和其他用地，其數值分別為10.50%、3.92%、3.91%、-2.36%和-1.35%，由于林地面積比重最大，使得其面積變化量較高、單一土地利用動態度最小，林地的單一土地利用動態度為-0.59%。2013～2017年間，面積變化最大的地類是水體，變化總量為-2033.35 hm2，其次是其他用地，變化總量為1916.06 hm2，居民地面積變化最小，變化總量為-254.52 hm2。在單一動態度方面，其他用地動態度最大，其數值為61.02%，其次分別是水體、園地、工礦用地、耕地和居民地，其數值分別為-8.77%、-3.91%、-3.16%、1.05%和-0.4%，林地的單一土地利用動態度最小為0.21%。以上數量變化表明在2009～2013年期間工礦用地增加速度最快，園地減少速度最快，工礦用地和居民地面積的大幅增加，主要由林地、耕地和園地轉化而來，其轉化原因與經濟增速發展和人口數量的增加有關。在2013年到2017年期間其他用地增加速度最快，水體減少速度最快。

從總體來看：2009～2017年間，遷西縣各土地利用類型之間變化較為頻繁，綜合土地利用動態度為1.64%。其中，耕地、水體、林地、園地的面積呈減少趨勢，且減少量依次上升，工礦用地、其他用地、居民地的面積呈增加趨勢，且增加量依次上升。從變化數量看，居民地的大幅增加和園地的大幅減少最為顯著，8年來居民地面積凈增1911.68 hm2，年變化量為238.96 hm2，園地凈減少1611.94 hm2，年變化量為201.49 hm2；從變化速度來看，其他用地的變化速度最快，平均每年增加28.18%，其次是園地，平均每年減少6.39%，工礦用地平均每年增加4.73%，水體平均每年減少3.11%，居民地平均每年增加1.72%，耕地平均每年減少0.7%，林地由于占總體面積比重較大，是變化速度最小的用地類型，平均每年減少0.19%。就整個研究區而言，2009～2013年土地利用行為相對穩定有序，2013～2017 年土地利用類型間流轉程度較為劇烈。

2.2 景觀指數分析

利用ArcGIS將土地利用分類結果轉化為柵格圖，為了保留原有數據的細節信息，柵格數據的像元大小設置為原始影像數據的分辨率大小，利用Fragstats 4.2軟件進行景觀指數的計算，最終繪制斑塊類型水平下各景觀指數的雷達圖(圖3)。

圖3 斑塊類型指數雷達變化圖

(1)林地

在研究時段內，林地的最大斑塊指數和邊緣密度最大，并且在2017年達到峰值，而林地的分散和并列指數變化不大，說明林地景觀一直是遷西縣的優勢景觀。同時，林地與其他斑塊類型的鄰接程度較高，斑塊數量持續增加說明林地景觀的破碎化程度也在增加。

(2)耕地

耕地的斑塊數量、邊緣密度、斑塊密度和景觀形狀指數一直上升，而分散和并列指數一直下降，表明耕地圖斑受地形起伏的影響形狀更加復雜，與其它景觀類型的嵌合度大，分散性強，團聚性差，破碎化程度逐年增加，相比較耕地與林地，雖然耕地的占有面積較少，但是由于整體破碎度高，空間結構更復雜，所以景觀形狀指數高于林地。

(3)水體

水體的分散和并列指數最大并且保持上升趨勢，表明水體景觀的破碎化程度較低，人類活動對水體景觀的影響不明顯，水體景觀維持在一個較為穩定的狀態，其他類型的斑塊受水體分布廣的影響，與水體斑塊鄰接程度很高，其分散和并列指數值也普遍較高。另外，2009年到2013年研究區內斑塊密度和邊緣密度都有所增加，主要是因為研究所用的2013年遙感影像獲取時間正值雨量充沛的7月份，溝渠、水庫蓄水量增大，淹沒了周圍區域，使得該景觀斑塊面積增大。從2013年分類圖上可以看出水體形狀變得狹長曲折，且比較復雜，這使得水體的景觀形狀指數增大；而2013年到2017年研究區內斑塊密度和邊緣密度都有所降低，形狀指數轉而減小，說明水體自身形成了較大的斑塊，水體的形狀趨于規則和簡單。

(4)居民地和工礦用地

居民地和工礦用地的斑塊數量、邊緣密度、斑塊密度和景觀形狀指數先上升后下降總體趨于上升，但分散和并列指數先下降后上升總體趨于下降，說明人類活動的干擾造成不同斑塊間的相互鑲嵌。居民地和工礦用地的破碎化程度一直在增加、斑塊形狀越來越復雜。但從2013年開始，居民地和工礦用地的破碎化程度有降低的趨勢，斑塊形狀趨向規則、簡單，說明城鎮化發展中經過合理的規劃引導和空間資源的有序配置，居民地和工礦用地景觀格局朝著向好的方向發展。

(5)指數分析

景觀水平指數主要識別整個景觀鑲嵌體的空間特征，為了了解遷西縣的整體景觀空間構型，計算2009、2013和2017年的景觀水平指數如表2所示。

表2 不同時期的景觀水平指數

從景觀水平上來看，研究時段內蔓延度指數呈遞減趨勢，在2017年達到最小，其值為65.083，說明研究區中林地景觀斑塊類型連接性減弱，研究區內不同斑塊類型的團聚程度在下降，空間分布逐漸趨向分散，研究區景觀的異質性增強；景觀連接度指數均在99.4%以上，并呈現出逐年遞增的趨勢，說明整個研究區是以少數大斑塊為主體構建起來的，其他不同斑塊之間的連接度較為緊密，研究區內景觀破碎度有所改善；2009年到2013年聚合指數和Shannon多樣性指數呈增長趨勢，2013年以后呈降低趨勢，說明2013年前遷西縣景觀類型趨于豐富多元化，各種景觀類型的發展變化，使得景觀中存在較多的小斑塊，景觀破碎化程度增加；2013年后景觀格局連通性向好的趨勢發展，景觀中不同斑塊類型聚集程度有所減小，各種景觀類型之間逐漸形成較好的連接性，且景觀格局破碎度緩和，朝著均勻化和多樣化方向發展。

2.3 景觀脆弱性分析

利用網格分析法對研究區景觀脆弱性進行評價，綜合考慮研究區的大小和采樣區的工作量，將研究區分割為4 km*4 km的單元格網，按照從上到下從左至右的順序依次編號，生成編號為1-119的119個景觀類型圖文件。同時，每個方格網中的景觀類型圖在ArcGIS中批量轉換為柵格文件，在Fragstats 4.2中導入前期在Excel中制作好的batch文件批量計算景觀指數，計算每個格網的脆弱性指數并將其作為每個格網中心點的脆弱性值，以每個格網中心點的ID標識碼和格網編號為連接字段建立屬性連接關系，利用統計學中的克里金空間插值法繪制景觀格局脆弱性分布圖，最后按照相等間隔將研究區劃分為5個等級分區：低脆弱區、較低脆弱區、中等脆弱區、較高脆弱區和高脆弱區。該尺度下的景觀脆弱性分布如圖4所示，統計研究區5個等級分區所占面積比例見表3。

圖4 4 km*4 km采樣區景觀格局脆弱性分布圖

表3 景觀脆弱分區面積比例

結合圖4和表3可以看出：2009～2017年研究區景觀脆弱性整體呈現出以較低脆弱區、中等脆弱區和較高脆弱區分布較多，而低脆弱區和高脆弱區分布較少的空間布局，且整體脆弱性程度從西北到東南大致呈現降低趨勢。

從景觀脆弱性不同等級來看，研究區2009年低脆弱區和較低脆弱區都不具優勢，占整個研究區的面積比例分別為5.72%、17.86%，主要分布在研究區東部及南部邊緣區域，該區域零散分布有林地、耕地、居民地和工礦用地，景觀適應性提高，敏感性和脆弱性較低；中等脆弱區占研究區面積比例為21.71%，主要呈半包圍狀分布于研究區中部地區，該區域林地、居民地和工礦用地相互嵌套，耕地景觀較為聚集，景觀適應性處于一般水平，多樣性不明顯，導致該區域有一定的脆弱性；較高脆弱區和高脆弱區占研究區面積比例分別為37.34%、17.37%，主要分布在研究區北部、西北部及中部地區，以林地為主，耕地為輔，其景觀類型相對單一，該區域也是灤河流經的區域，灘涂等其他用地面積變化較大，導致景觀脆弱性值較大。

與2009年相比，研究區2013年低脆弱區和較低脆弱區面積比例都有所增加，空間分布從研究區東部及東南部邊緣區域向中部及南部擴展，占據2009年部分中等脆弱區和較高脆弱區；中等脆弱區面積比例較2009年有所增加，其變化走向朝西擠壓并向東北方向延伸，占據了一些景觀較高脆弱區和高脆弱區；處于較高脆弱區和高脆弱區的面積比例大幅減少，較高脆弱區空間位置整體向西部靠攏，高脆弱區變動劇烈，由北部及西北部區域遷移至西南區域，在東南部區域也出現小范圍的較高脆弱區和高脆弱區。

2017年，研究區內處于低脆弱區面積比例較2013年大幅下降，由14.67%下降至3.10%，較低脆弱區面積比例有所增加，由23.34%上升到25.19%，兩者的空間位置從中部及南部區域逐步向東部邊緣地區擴展，出現在2013年東南部區域的較高脆弱區和高脆弱區完全消失，被低脆弱區和較低脆弱區占據。結合景觀分類圖可以發現，造成該片區域景觀脆弱性發生劇烈變化的原因是水體景觀和耕地景觀變化導致的生態環境失衡、脆弱性升高，經過綜合整治后的水體景觀和耕地景觀朝向好的方向發展，景觀適應性提高，從而脆弱性降低；研究區中等脆弱區的面積比例較2013年持續增加，由27.48%上升到29.69%，造成該變化的原因是，自2013年以來耕地景觀的合理規劃加上對土地資源的高效整合，造成中等脆弱區面積增大，空間分布呈半包圍趨勢朝西南方向擴展；相較2013年、2017年研究區處于較高脆弱區的面積比例減小，由30.22%下降至21.16%，高脆弱區面積比例大幅增加，由4.29%上升到20.86%，遷西縣中西部山巒環抱，灤河干流縱貫遷西縣境，兩山夾一河的地形條件在加快當地旅游業發展的同時也在影響著各類景觀要素的穩定性，該區域受人類活動的影響較大。因此，對外界的干擾較為敏感，較高脆弱區大多分布于該區域；同時結合景觀分類圖可以看出，研究區高脆弱區域分布有大量的其他用地景觀，景觀類型的頻繁更替導致該區域的景觀脆弱性達到最高。

3 結語

對遷西縣景觀變化及脆弱性分析應考慮以下方面：雖然描述景觀格局的指數有很多，但許多景觀指數之間存在很大的相關關系，因而在景觀格局研究中，選取信息不重復的景觀指數是值得探討的問題。自然因素和人為因素是遷西縣景觀格局脆弱性發生變化的主要原因，隨著一系列生態保護政策的實施，遷西縣景觀脆弱性會逐漸穩定。此外，本文所構建的景觀脆弱性評價指標中，低脆弱性并不表示生態環境的質量優良，僅代表了區域生態處于相對穩定的狀態，可恢復性較好。

遷西縣土地利用從變化數量來看，以居民地的大幅增加和園地的大幅減少最為顯著；從變化速度來看，其他用地的變化速度最快，然后依次為園地、工礦用地、水體、居民地、耕地、林地，林地是變化速度最小的用地類型。但大體上遷西縣耕地、水體、林地和園地的面積仍呈逐漸減少趨勢，工礦用地、其他用地和居民地的面積逐步增加。

通過斑塊水平和景觀水平對遷西縣景觀特征進行分析，遷西縣林地景觀和耕地景觀的斑塊數量一直處于上升趨勢，景觀邊界越來越復雜，破碎化程度逐漸增長；水體景觀受降雨量等的影響，形狀趨于簡單和規則，破碎化程度逐漸降低；居民地景觀和工礦用地景觀在人類活動的干預下，形狀呈現出由復雜向簡單過渡的趨勢，破碎化程度逐漸得到緩和；其它景觀破碎化程度變化不大。在景觀水平上，連接度指數和Shannon多樣性指數增加，蔓延度指數和聚合指數減少，表明遷西縣景觀趨于豐富多元化，景觀異質性減少，破碎化程度朝著均勻化和多樣化方向發展。

遷西縣景觀脆弱性呈現出由西北到東南逐步降低的趨勢，且脆弱性西部高于東部，北部高于南部，中西部分布有大量的其他用地景觀，景觀類型的頻繁更替導致高脆弱區面積有所擴張。但從整體來看，較低脆弱區和中等脆弱區面積比例上升，較高脆弱區所占面積比例下降，景觀脆弱性較嚴重的地區趨于分散，逐漸被脆弱性較低的區域所占據，遷西縣景觀脆弱性整體有所緩和，但局部區域仍需加強治理，保證該區域的可持續發展。