基于Fragstats4的景觀格局指數與地表溫度的相關性
——以石家莊市為例

楊鵬 高祺* 張艷品 劉思廷 齊曉華 曹春莉 程巳陽

(1 石家莊市氣象局，石家莊 050081; 2 中國氣象科學研究院災害天氣國家重點實驗室，北京 100081)

引言

隨著城市化的發展，土地利用/覆蓋和景觀格局都發生變化，水泥、瀝青等不透水地表組成的建成區持續擴張，水體、天然植被等構建的生態空間不斷萎縮，交通及路網用地密度提升，景觀格局完整性遭到破壞[1-5]。近年來，從景觀生態學的角度開展土地利用尺度效應的研究是一個逐漸興起的研究方向。土地利用景觀格局變化研究對優化城市布局和加快美麗城市建設至關重要，景觀的組成和結構特征是影響熱場分布的主要因素，組成特征對熱場的影響主要是由于各地類熱效應不同，結構特征主要影響斑塊間能量交換[6-8]。

國內外許多學者從不同方面、利用不同方法開展景觀格局指數與地表溫度的相關性研究。Callo等[9]最先利用NOAA/AVHRR數據證明NDVI與LST有負相關；Weng Qihao等[10]利用景觀格局指數研究不同土地利用景觀格局對熱環境的影響；Coseo等[11]通過對芝加哥城市熱島影響因素分析得出，不透水地表比率和數目可以解釋68%的氣溫變化；江學頂[12]采用遙感和MM5模式研究珠江三角洲城市群的城市熱力景觀及其日變化規律；馬勇剛等[13]以烏魯木齊為研究區，將景觀分類和景觀格局變化圖與地表溫度圖進行對比，研究景觀格局變化對城市熱島的影響；徐雙等[14]以Landsat7 ETM+數據和土地利用現狀數據為數據源，分析了2004年和2010年長沙市熱力景觀空間格局演變；吳未等[15]以無錫市2010年TM數據研究了快速城市化平原地區不同景觀指數隨粒度的變化特征及適宜尺度的選擇。

本文以石家莊市區作為研究區，選擇1987年、2004年和2019年Landsat遙感影像數據，以城市景觀格局為視角，深入研究景觀格局指數與地表溫度的相關性，可以掌握城市化發展對熱環境的影響，可以為城市的合理規劃、生態環境改善以及可持續發展提供理論依據和技術支撐。

1 研究區概況

本文以石家莊市區作為研究區(圖1)。石家莊市地處河北省中南部，太行山東麓，位于37°27′～38°47′N，113°30′～115°20′E之間，屬于溫帶半濕潤半干旱大陸性季風氣候， 太陽輻射的季節性變化顯著，四季分明。近年來，由于城市工業化進程加快，使得城市綠地越來越少，城市景觀結構發生改變，城市工業、交通運輸和居民向大氣中排放的人為熱和廢熱日益增加，城郊之間某些氣象要素已出現明顯差異， 城市熱島效應愈發顯著。通過研究景觀格局指數與地表溫度的相關關系，可以掌握城市景觀對城市環境的影響，為石家莊城市建設、城市規劃和生態平衡的改善提供一定的理論依據和技術支撐。

圖1 石家莊市區地理位置

2 資料和方法

2.1 Landsat遙感數據及處理

本文所需數據資料主要涉及了Landsat衛星遙感資料、氣象資料以及相關的地理信息資料等。遙感影像數據來源于國際科學數據服務平臺[16]，條帶號/行編號為124/34，數據格式為 .TIF，數據選取時間為1987年7月31日、2004年8月30日和2019年7月23日，其中1987年和2004年為Landsat5 TM數據，2019年為Landsat8 OLI_TIRS數據，具體遙感數據信息見表1所示；氣象資料來源于河北省氣象數據共享服務平臺；地理信息資料由河北省地理信息局提供，數據格式為.shp格式。

表1 Landsat遙感數據信息

由于遙感影像在采集過程中不可避免地會受到傳感器自身光學性質和結構特點等內部因素以及太陽輻射、地形條件等外部因素的影響，為保證研究結果的可靠性，需要對遙感數據進行預處理，主要包括波段合成、輻射定標、大氣校正、幾何校正和邊界圖裁剪等[17]。遙感數據預處理操作過程均在ENVI 5.2下進行，經過處理后，1987年、2004年和2019年遙感影像圖見圖2所示。

圖2 1987年(a)、2004年(b)和2019年(c)石家莊衛星遙感影像

2.2 監督分類方法

監督分類方法又稱訓練場地法，是以建立統計識別函數為理論基礎，依據典型樣本訓練方法進行分類的技術[18]。即根據已知訓練區提供的樣本，通過選擇特征參數，求出特征參數作為決策規則，建立判別函數以對各待分類影像進行的圖像分類，是模式識別的一種方法。判別準則若滿足分類精度要求，則此準則成立；反之，需重新建立分類的決策規則，直至滿足分類精度要求為止。目前比較常見的監督分類方法主要有平行六面體法、最大似然法、最小距離法、馬氏距離法、神經網絡法、支持向量機法等。本文利用最大似然法的監督分類方法對遙感影像數據進行分類處理，根據研究需要和4種類型的影像圖的表征特點(表2)，將景觀類型分為綠地、水體、不透水地表和未利用地等四種類型，分類結果見圖3所示。通過精度計算得到1987年、2004年和2019年下墊面分類的總體精度分別為98.23%、98.79%和96.48%，Kappa系數分別為95.92%、97.43%和95.12%，總體精度和Kappa系數均在95%以上，表明分類精度很高，能夠滿足研究要求。

表2 4種下墊面類型的樣本影像

圖3 1987年(a)、2004年(b)和2019年(c)石家莊4類下墊面分類結果

2.3 地表溫度反演算法

本文針對Landsat5 TM數據，采用單窗算法反演地表溫度[19]；針對Landsat8 OLI數據，采用劈窗算法反演地表溫度[20]。具體方法及參數請參考文獻。其中， Landsat5 TM數據的熱紅外通道為第6波段，Landsat8 OLI數據是第10和第11波段。本文通過相同時間點的地面氣象觀測站點的地溫數據作為驗證數據源對反演精度進行檢驗， 發現單窗算法和劈窗算法的誤差絕對值分別在1.2 K和0.8 K之內，反演結果能夠滿足研究需要，數據可用。

2.4 景觀格局指數

景觀格局指數能夠比較精確地反映景觀的結構組成和空間配置信息，是一個簡單的定量指標。根據描述對象的等級，可分為斑塊水平指數、斑塊類型水平指數和景觀水平指數[21-22]。本文利用Fragstats 4.2軟件計算景觀格局指數，在類型水平上選取8種景觀格局指數，在景觀水平上選取9種景觀格局指數，具體描述見表3、表4所示。

表3 類型水平上選取的景觀指數定義

表4 景觀水平上選取的景觀指數定義

2.5 移動窗口

移動窗口(Moving window)是通過對窗口內所選中景觀指標進行統計計算，輸出對應于所選景觀指數柵格圖，移動窗口算法[23]可以從空間角度展示出景觀格局的動態變化過程。移動窗口算法的基本原理是將窗口(用戶指定形狀和大小)放置在每個柵格像元上，計算選定的度量并將值返回給柵格像元。因此，柵格像元周圍的每個窗口都被視為子景觀，并對每個單元重復此過程，移動窗口算法原理的示意,如圖4所示，窗口a、b、c計算結果分別寫入柵格1、2、3中，以次類推。

圖4 移動窗口算法原理(窗口a、b、c計算結果分別寫入柵格1、2、3中)

3 結果與分析

3.1 地表溫度反演

圖5給出了1987年、2004年和2019年地表溫度反演結果?？梢钥闯?，1987年、2004年和2019年研究區地表溫度空間分布均表現出中心高、周邊低的現象，城市熱島效應比較明顯。隨著石家莊城市建設的快速發展，市區熱島發展較快，熱島強度和面積不斷增加，東-西向發展明顯；1987年和2004年滹沱河河道及附近存在一條比較顯著的“高溫帶”，隨著滹沱河的綜合治理及蓄水工程的實施，“高溫帶”逐漸消失，現已成為石家莊市區的“生態冷源”。

圖5 反演地表溫度:(a)1987年，(b)2004年，(c)2019年

3.2 景觀粒度的選擇

景觀粒度是景觀中可識別的最小斑塊單元，有關研究發現，景觀格局指數會隨粒度的變化而發生變化[24-25]。本文將研究區重采樣為空間分辨率為30 m、60 m、90 m、120 m、150 m、180 m、210 m、240 m、

270 m和300 m的粒度大小，選取蔓延度指數(CONTAG)、景觀導度指數(COHESION)、香農多樣性指數(SHDI)和香農均度指數(SHEI)等景觀格局指數，研究不同粒度下景觀格局指數的變化。

圖6給出了1987年、2004年和2019年不同粒度下景觀格局指數的變化。從圖6可以看出，1987年、2004年和2019年不同年份的景觀格局指數大小不同，但變化趨勢基本一致；四種景觀格局指數均隨粒度的不同而發生變化，蔓延度指數和景觀導度指數變化較為明顯，香農多樣性指數和香農均度指數變化相對較小。研究表明，景觀指數隨粒度變化的第1尺度域是確定粒度大小，進行景觀格局分析的適宜取值范圍，根據研究區各景觀指數變化曲線的拐點和躍變區間確定最佳適宜粒度[26]，從圖6可以看出，在120 m、180 m、240 m和270 m處出現拐點，從圖6b景觀導度指數(COHESION)看，90 m～150 m范圍為第1尺度域，210 m～270 m范圍為第2尺度域，由此可以判定，90 m～150 m范圍內為最佳適宜粒度，由于熱紅外波段的空間分辨率為120 m，經綜合考慮，最終確定該研究區的適宜粒度為120 m。

圖6 不同空間粒度下景觀格局指數的變化:(a)蔓延度，(b)景觀導度，(c)香農多樣性，(d)香農均度

3.3 景觀格局指數計算

景觀格局指數能夠比較精確地反映景觀的結構組成和空間配置信息，是一個簡單的定量指標。根據描述對象的等級，可分為斑塊水平指數、斑塊類型水平指數和景觀水平指數。斑塊水平指數往往被作為計算其他景觀指數的基礎，而其本身對了解整個景觀的結構并沒有很大的解釋價值，通常采用的是斑塊類型水平指數和景觀水平指數。在獲得了各種類型水平指數和景觀水平指數的基礎上，對單個時相的景觀格局和不同時相的景觀格局變化進行分析，說明景觀的變化情況及變化趨勢。利用Fragstats 4.2軟件分別計算斑塊類型水平指數和景觀水平指數，結果見表5、表6所示。從表5可以看出，1987—2019年綠地的斑塊類型面積(CA)整體呈逐年減少趨勢，不透水地表和水體的斑塊類型面積(CA)在逐年增加。1987—2019年綠地的斑塊數量(NP)在逐年增加，綠地的破碎化指數(SPLIT)基本呈逐年增加趨勢，綠地的最大斑塊面積指數(LPI)在1987年、2004年和2019年分別為68.27%、12.40%和9.54%，呈逐漸減少趨勢，這都說明研究區范圍內綠地斑塊在不斷地裂化和減少。綠地的斑塊形狀指數(LSI)在1987年、2004年和2019年分別為19.69、22.15和36.37，說明綠地斑塊形狀不規則化程度加劇。綠地的聚集度指數(AI)在1987年、2004年和2019年分別為89.44%、86.51%和71.66%，說明研究區內綠地景觀斑塊的離散程度增大，綠地斑塊的空間優勢性大大降低。不透水地表的最大斑塊面積指數(LPI)在1987年、2004年和2019年分別為12.09%、29.58%和54.18%，呈逐年增加趨勢，不透水地表的破碎化指數(SPLIT)在1987年、2004年和2019年分別為68.15、11.41和3.40，呈逐年減少趨勢，說明研究區范圍內不透水地表斑塊在不斷的增加和聚合，形成優勢斑塊，說明城市化進程對城市景觀格局影響較大。

表5 1987、2004、2019年類型水平上的景觀指數

景觀水平指數可以從全局的角度對研究區的景觀水平進行分析，通過表6可以看出，1987—2019年最大斑塊面積指數(LPI)呈現先減少、后增加的趨勢，破碎化指數(SPLIT)與之相反，呈現先增加、后減少的趨勢，說明隨著城市化進程的加快，綠地斑塊面積在逐漸減少和裂化，綠地的景觀優勢在不斷下降，不透水地表的面積在逐漸增加和聚合，不透水地表的景觀優勢在不斷加強，逐漸形成優勢景觀，這一結論合理解釋了研究區范圍內優勢景觀由綠地景觀逐漸向不透水地表景觀演變的過程。另外，1987—2019年聚集度(AI)逐漸減少，表明景觀中不同類型的斑塊離散程度在逐漸增大。

表6 1987、2004、2019年景觀水平上的景觀指數

3.4 移動窗口尺寸選取

景觀格局指數的分析僅以數值描述景觀格局，難以直接描述景觀格局的空間形態，利用移動窗口可以實現景觀格局指數的空間可視化，這就更加有利于統計分析景觀格局指數與對應地表溫度的關系。本文利用Fragstats 4.2的移動窗口功能得到景觀格局指數空間分布圖，為了更加科學地確定研究區移動窗口的網格尺度大小，本研究以2019年為例，統計分析了不同網格尺度下不同景觀類型占比(PLAND)與地表溫度的關系(圖7)。從圖7可以看到，綠地和不透水地表與地表溫度的相關性較高，在360×360～1080×1080(m)范圍內，基本呈穩定增加趨勢，到1320×1320(m)則開始下降；水體則呈逐漸下降趨勢，未利用地呈波動趨勢，在600×600(m)處相關性最高，但未利用地與地表溫度的相關性最差，這是與未利用地的不規則分布有關，在本文分析中將不再考慮未利用地的分析部分。經反復測試，600×600和360×360(m)相差不大，由于水體的跨度在1400 m左右，相關研究[23-25]表明一般選取其尺度的1/2作為窗口尺寸為宜，經綜合考慮，最終確定以600×600(m)的正方形移動窗口進行景觀格局分析。

圖7 不同網格尺度下4種下墊面類型占比與地表溫度的相關系數(絕對值)

3.5 不同景觀格局指數之間的相關關系

以2019年為例，利用移動窗口法對2019年的遙感影像圖進行計算，得到9種景觀格局指數柵格圖，并對其進行相關性分析，結果見表7。從表7可以看出，很多景觀格局指數之間呈現出較高的相關性，如香農多樣性指數(SHDI)和最大斑塊面積指數(LPI)的相關性為-0.910，SHDI值越高，表明景觀種類越多，斑塊破碎化程度越高，其最大斑塊面積指數LPI越??；斑塊破碎化指數(SPLIT)和最大斑塊面積指數(LPI)的相關性為-0.946，說明破碎化程度與最大斑塊呈顯著負相關，斑塊越破碎導致LPI指數越??；斑塊形狀指數(LSI)和聚集度(AI)的相關性為-0.897，表明LSI與AI呈顯著負相關，斑塊形狀指數增加使其不規則度增加，進而導致聚集度減少；斑塊形狀指數(LSI)和香農多樣性指數(SHDI)的相關性為0.891，LSI增加使得其不規則度增加，最終導致斑塊破碎，香農多樣性指數SHDI增加。根據表7并考慮不同景觀格局指數的生態學意義，本文選取景觀占比(PLAND)、最大斑塊指數(LPI)、斑塊破碎化指數(SPLIT)和聚集度指數(AI)進行分析。

表7 2019年不同景觀格局指數之間的相關系數

3.6 景觀格局指數與地表溫度的相關性分析

本文選取景觀占比(PLAND)、最大斑塊指數(LPI)、斑塊破碎化指數(SPLIT)和聚集度指數(AI)，對2019年的遙感影像圖進行計算和分析，從景觀優勢度、破碎化以及聚集程度等角度，分析了研究區范圍內綠地、水體和不透水地表的景觀格局指數與平均地表溫度之間的相關性，結果見表8所示。

表8 4種類型景觀格局指數與平均地表溫度的相關系數

從表8可以看出，綠地、水體和不透水地表的景觀格局指數與平均地表溫度均表現出較好的相關關系，均為極顯著性相關，其中綠地和不透水地表的相關性最好，水體次之。斑塊百分比指數(PLAND)、最大斑塊指數(LPI)和聚集度(AI)與平均地表溫度表現出一致的極顯著相關關系(綠地和水體為負相關，不透水地表為正相關)，破碎化指數(SPLIT)則恰恰相反(綠地和水體為正相關，不透水地表為負相關)。綠地和水體占比與平均地表溫度相關系數分別為-0.844和-0.598，不透水地表占比與平均地表溫度相關系數為0.911，綠地和不透水占比與平均地表溫度的相關性高于水體，表明綠地為城市地表熱環境效應抑制的主要貢獻源，不透水地表是造成城市熱島的主要因素。平均地表溫度與斑塊百分比指數(PLAND)和最大斑塊指數(LPI)的相關系數明顯高于平均地表溫度與聚集度(AI)和破碎化指數(SPLIT)的相關系數，這一研究結果可以說明，一個優勢斑塊(景觀)對平均地表溫度的影響效果明顯大于幾個比較分散或破碎的斑塊(景觀)。

由此可以得出，一個較大規模的綠地(或水體)產生的冷島效應明顯強于幾個較小的分散綠地(或水體)；同理，規模較大或連片的不透水地表會比幾個較為分散的不透水地表產生更強的熱島效應。由此可見，合理優化城市景觀的空間布局是改善城市熱環境、減輕城市熱島效應最有效的途徑和手段。

4 結論

本文以石家莊市為例，基于1987年、2004年和2019年Landsat遙感影像數據，運用遙感、景觀生態學、地理信息系統和統計分析等理論方法，分析景觀格局指數與地表溫度的相關性，得出結論如下：

(1)1987—2019年綠地CA逐年減少，NP、LSI和SPLIT逐年增加，LPI和AI逐年減少，說明研究區范圍內綠地斑塊在不斷地裂化和減少，綠地斑塊形狀不規則化程度加劇，綠地景觀斑塊的離散程度增大，空間優勢性不斷降低；不透水地表和水體CA逐年增加，LPI逐年增加，SPLIT逐年減少，說明研究區范圍內不透水地表斑塊在不斷地增加和聚合，形成優勢斑塊，城市化進程對城市景觀格局影響較大。

(2)1987—2019年LPI呈現先減少、后增加的趨勢，SPLIT先增加、后減少，說明綠地斑塊面積在逐漸減少和裂化，綠地的景觀優勢在不斷下降，不透水地表的面積在逐漸增加和聚合，不透水地表的景觀優勢在不斷加強，解釋了優勢景觀由綠地景觀逐漸向不透水地表景觀演變的過程。

(3)LST與PLAND和LPI的相關系數明顯高于LST與AI和SPLIT的相關系數，說明一個優勢斑塊(景觀)對LST的影響效果明顯大于幾個比較分散或破碎的斑塊(景觀)，一個較大規模的綠地(或水體)產生的冷島效應明顯強于幾個較小的分散綠地(或水體)。