烏魯木齊城市地表溫度時空變化與貢獻度研究

柴晉鵬，陳學剛, ，阿里木江·卡斯木, ，高鵬文

1. 新疆師范大學地理科學與旅游學院，新疆 烏魯木齊 830054；2. 新疆師范大學絲綢之路經濟帶城市化發展研究中心，新疆 烏魯木齊 830054

城市是人類社會發展的產物，是人類生產生活的主要場所。伴隨著城市化進程，城市局部氣候發生變化，首先表現為城市與郊區氣溫的差異，形成了城市熱島（宋艷玲等，2003），其直接表現為城市地表溫度的差異。土地利用是人類對自然環境的改變和干預的主要表現，其直接表現為土地的覆蓋狀況。城市地表溫度與城市土地利用/覆蓋狀況密切相關。隨著遙感反演地表溫度技術的日益成熟，使獲取較為準確的地表溫度成為可能，為城市熱島研究提供了技術支持。研究城市地表溫度時空演變規律其在空間上的集聚效應，不僅能夠深入理解土地利用/覆被變化下城市熱環境的空間特征和動態變化，而且對大氣污染、市政建設及土地合理規劃和利用具有重要意義（張心怡等，2005）。

近年來，國內外學者對城市熱島與土地利用/覆被的研究主要集中在兩個方面。一方面，探索其動態變化與空間分布特征，利用長時間多幅遙感影像，對城市熱島與土地利用的變化情況進行分析，揭示二者之間的關系（Ahmed，2018；歷華等，2008；蘇振華等，2015）。另一方面探索對其影響因素進行定量研究，（1）植被指數與地表溫度的關系研究（Melkamu，2019；唐曦等，2008；李成范等，2009；張曉娟等，2018）；（2）城市用地擴展與地表溫度的關系研究（Qian et al.，2006；Amiri et al.，2009；徐涵秋等，2004；貢璐，2007；黃偉嬌，2018）；（3）景觀指數與地表溫度的關系研究（游絢等，2009；徐雙等，2015；王雪等，2017）。然而，在西北干旱區綠洲城市特殊的地理背景下，基于空間尺度模式對城市區域熱島集聚程度并定量化評價各土地利用/覆蓋類型的對地表溫度的貢獻度指數進行研究的案例甚少。

快速城市化的背景下，作為新疆首府的烏魯木齊市經濟快速發展，隨著人口的不斷增長與建設用地的擴張，地表溫度也發生及其顯著的變化。選定烏魯木齊市規劃的中心城區為研究對象，研究用空間尺度模式，揭示研究區地表溫度的時空格局變化及其特征并定量評估土地利用/覆蓋類型對 LST的貢獻度，以期為干旱區綠洲城市建設與可持續發展提出建議，為有關部門決策提供參考。

1 研究區概況與數據來源

1.1 研究區概況

烏魯木齊地處亞歐大陸腹地，我國西北，新疆中部，北鄰準噶爾盆地，南接天山山脈，地處86°37′33″—88°58′24″E，42°45′32″—45°00′00″N，屬溫帶大陸性氣候，全年 7、8月最熱，平均氣溫在25.7 ℃，1月最冷，平均氣溫在-15.2 ℃。該區域地勢起伏懸殊，山地面積廣大。南部、東北部高，中部、北部低。市區三面環山，北部為開闊平原。全市轄7區1縣：天山區、沙依巴克區、高新技術開發區（新市區）、水磨溝區、經濟技術開發區（頭屯河區）、達坂城區、米東區、烏魯木齊縣，總面積 1.42×104km2，中心城區約為 1455 km2。2018 年，常住人口約 350.58萬人，人均生產總值達 88418元（https://baike.baidu.com/item/烏魯木齊/121593）。根據烏魯木齊市城市總體規劃（2014—2020年），選取中心城區進行研究，如圖1所示。

圖1 研究區位置示意圖及Landsat 8影像（2018-08-10）Fig. 1 Location of study area and Landsat 8 imege (2018-08-10)

1.2 數據來源及其預處理

以Landsat衛星影像數據作為反演地表溫度的基礎數據，數據影像來源于地理空間數據云（GS Cloud）和美國地質勘探局（USGS）。選取研究區的高溫月份6—8月，且天氣晴朗，云量少于10%，成像較好的影像，獲取的影像衛星過境均均在格林尼治時間04:00—05:00時之間（北京時間12:00—13:00）（表1）。借助ENVI 5.3對影像進行幾何校正、輻射定標、大氣校正、圖像鑲嵌、圖像裁剪等預處理。文中驗證地表溫度及計算反演地表溫度所需的大氣參數所用的氣溫、濕度等氣象數據均來自中國氣象數據網（http://data.cma.cn）。采用中國科學院發布的中國土地利用現狀遙感監測數據庫 30 m×30 m數據的作為土地利用/覆蓋數據，該數據是基于Landsat系列衛星影像經人工解譯而成，并通過了質量檢驗（http://www.resdc.cn/）。

表1 影像數據來源列表Table 1 List of image data sources

2 研究方法

通過對1990—2018年的4期影像數據預處理之后，采用地表溫度反演、重心轉移模型、空間聚類分析以及貢獻度指數計算等方法，定量研究烏魯木齊市城市化進程中地表溫度的時空格局變化及其在空間上的集聚效應。

2.1 地表溫度反演

利用Landsat衛星影像進行地表溫度反演已有諸多國內學者進行研究（JImenez-Munoz et al.，2014；Wang et al.，2016；覃志豪等，2001；胡德勇等，2015；宋挺等，2015；蔣大林等，2015），其主要包括輻射傳輸方程法、單窗算法，劈窗算法等。諸多學者也對多種算法的準確性進行了對比研究，結果表明針對Landsat數據單窗算法反演精度，可操作性均優于其算法（李斌等，2017；金點點等，2018），且Landsat OLI/TIRS第11波段具有較大的不確定行（徐涵秋，2016），因此主要采用覃志豪和胡德勇等（2015）提出的單窗算法對 Landsat TM和Landsat OLI/TIRS的第6和第10波段進行地表溫度反演。具體步驟如下：

2.1.1 大氣頂部（TOA）光譜輻射率

利用輻射亮度校正因子，對熱紅外波段進行輻射定標，將DN值轉換為大氣頂部（TOA）光譜輻射率，即：

式中，Lλ為波λ的大氣頂部光譜輻射率（W·m-2·sr-1·μm-1）；ML為增益參數，Qcal為像元 DN值（即灰度值），AL為偏移參數，其中，ML和AL參數可直接從影像元數據文件MTL.txt中獲取。

2.1.2 大氣頂部（TOA）亮度溫度

利用元數據文件中的熱常數值，根據普朗克定律將輻射亮度轉化為亮度溫度，即：

式中，Tb大氣頂部（TOA）亮度溫度（K）；L為熱紅外波段的亮度值；K1、K2為常數。對于Landsat TM，K1=607.76（W·m-2·sr-1·μm-1），K2=1260.56 K；TIRS 傳感器 Band 10，K1=774.89（W·m-2·sr-1·μm-1），K2=1321.08 K。

2.1.3 地表比輻射率的計算

利用經驗公式，采用一種較為簡單的混合像元構成比例估算方法來確定地表構成比例，從而進行地表比輻射率的計算（覃志豪等，2004a）。對于城市區域，我們簡單的將其分為水體、自然表面和建筑表面3種，因此針對混合像元尺度上的地表比輻射率通過下式來估算：

式中，ε為混合像元的地表比輻射率，Pv為植被覆蓋率，Rv為植被的溫度比率，Rm為建筑表面的溫度比率，εv為植被地表比輻射率，針對TM影像，εv=0.986，εm=0.972；對于OLI/TIRS第10波段數據，εv=0.98672，εm=0.96767，dε為輻射校正項，由經驗公式得；

利用覃志豪等（2004b）的研究結果，植被、建筑表面的溫度比輻射率為：

式中，植被覆蓋度Pv可由歸一化植被指數NDVI算得：

式中，采用覃志豪的經驗值 NDVIV=0.70，NDVIS=0.05來進行估計。

2.1.4 大氣透射率的確定

大氣透射率主要受大氣水含量的影響（覃志豪等，2003），提出了大氣透射率的估算方程，獲取影像的時間均為夏季，且為北京時間12:00時左右，氣溫較高，因而選取高氣溫的估算方程：

式中，ω為大氣水汽含量（g·cm-2），根據楊景梅等（2002）提出的算法利用地面濕度計算所得，τ為大氣透射率。

2.1.5 大氣平均作用溫度的計算

在天氣晴朗，沒有明顯的大氣垂直漩渦作用條件下，可由地面氣溫近似計算大氣平均作用溫度，研究區地處中緯度地帶，且影像獲取時間均為夏季，因而選取：

式中，T0為地面氣溫（K），Ta為大氣平均作用溫度。

2.1.6 地表溫度計算

式中，a和b為常量，a= -67.355351，b=0.458606；Ta為大氣平均作用溫度（K）；Tb為亮度溫度（K）；τ為大氣透射率；ε為地表比輻射率；C和D為大氣透射率τ和地表比輻射率τ計算所得。

2.1.7 相對地表溫度分級（RLST）

為減小獲取的影像具有時間上的不同的而帶來的差異，實現各年度地表溫度的可比性，采用Yu et al.（2019）的相對地表溫度分級（RLST）：

式中：RLSTi為第i年的相對地表溫度，LSTi為第i年的地表溫度，為研究區第i年的平均地表溫度，把 RLST劃分為 6個等級，分別為RLST＜0，0≤RLST＜2，2≤RLST＜4，4≤RLST＜6，6≤RLST＜8，RLST≥8。

2.2 重心轉移模型

空間重心是描述地理對象空間分布的一個重要指標，通過重心的轉移軌跡來反映溫度區在空間上的轉移特征（李文亮等，2010）。具體公式如下：

式中：Xt，Yt表示第t年某類型的重心坐標，n為第t年該類的斑塊數量，Cti為第t年該類斑塊第i個斑塊的面積，Xti，Yti分別表示第t年該類斑塊中第i個斑塊的幾何中心坐標。

2.3 空間聚類分析

為探討研究區地表溫度在空間上的集聚特征，采用平均最近鄰分析和多距離空間聚類分析對相對地表溫度進行分析。平均最近鄰分析和多距離空間聚類分析已被廣泛用于景觀格局研究（Jill，2006；高凱，2010）、物種種群空間分布研究（Sébastien et al.，1999）、區域經濟研究（王結臣等，2016）等方面。利用ArcGIS將RLST經行面點數據轉換，用以分析研究區熱環境的空間模式（康文敏等，2020）。

2.3.1 平均最鄰近分析

最近鄰距離是表示地理空間上的點要素相互之間鄰近程度的指標，能夠較為準確客觀的反映點數據的空間格局分布特征（余兆武等，2016）。平均最近鄰是基于空間距離，通過計算任意一個點到最近鄰點的距離來度量空間聚集性的一種空間分析方法。通過觀測均距離除以期望的平均距離計算得出平均最近鄰比率（NNR）來表征空間模式。計算公式如下：

式中：n為RLST點要素的數量，di為每個RLST點要素與其最鄰近點要素的距離，A為研究區面積。若 NNR=1，RLST點要素呈現完全隨機分布；若NNR＜1，RLST點要素呈現集聚分布；若NNR＞1，RLST點要素呈分散分布。

2.3.2 多距空間離聚類分析

在不同的空間尺度下，格局分布模式存在較大差異（王勁峰，2010；Wu，2004），因而引入 Ripley’sK函數分析多尺度下熱環境的空間分布格局（高凱等，2010）。Ripley’sK函數是點密度距離函數，可分析研究區內任意尺度下點狀地物的空間分布特征，是全域統計指標，按照給定的半徑，計算而生成關系函數（高超等，2015）。計算公式如下：

式中：A為研究區面積，n為RLST點要素的數量，Kij為權重，d為空間尺度大小。若L(d)＞0，RLST點要素在對應空間尺度d上呈集聚分布；若L(d)=0，RLST點要素在對應空間尺度d上呈完全隨機分布；L(d)＜0，RLST點要素在對應空間尺度d上呈分散分布。

2.4 貢獻度指數

不同年度的土地利用格局存在較大的空間異質性，其對地表溫度的影響也呈現出較大的空間異質性（孫宗耀等，2018）。通過貢獻度指數（Contribution Index，CI）來定量衡量土地利用/覆蓋對熱環境的貢獻程度：

式中：i為i類土地利用類型，CIi為i類土地利用類型對熱環境的貢獻度，為i類土地利用類型的平均地表溫度，為研究區平均地表溫度，Si為i類土地利用類型的面積，S為研究區總面積。若CIi≥0，則i類土地利用類型對地表溫度的升高呈正相關；若 CIi≤0，則類土地利用類型對地表溫度的升高呈負相關。

3 結果與分析

3.1 地表溫度的反演與分級

為保證地表溫度反演的準確性，采用2018年8月10日的Landsat 8 OLI/TIRS影像反演的地表溫度，與研究區內氣象站點的大氣溫度做相關分析。因衛星過境時間為格林尼治標準時間2018年8月10日04:55，北京時間2018年8月10日12:55，因而獲取2018年8月10日13時研究區內17個氣象站點的氣溫數據。通過實地調研，對 17個站點的位置進行校正（圖2）。通過ArcGIS軟件以站點位置為中心構建100 m緩沖區，將緩沖區內的平均溫度與站點溫度做相關性分析，發現其在0.05水平呈相關性顯著正相關，皮爾遜相關性為0.55，因而表明通過單窗算法反演的 LST在研究區內具有一定的可靠性，采用相對地表溫度對反演結果進行分級，結果如圖3所示。從圖中可以看出，研究區內相對地表溫度等級發生了較為顯著的變化。

圖2 地表溫度驗證氣象站點位置圖Fig. 2 Location map of weather station for verification of land surface temperature

圖3 1990—2018年相對地表溫度等級分布圖Fig. 3 Rlative land surface temperature zones from 1990 to 2018

從4期溫度等級分區做轉移矩陣發現（表2），近28年來，研究區內RLST在＜0、0—2、＞8等級上面積增加，2—4、4—6、6—8等級上面積減少。其中，面積增長最大的是RLST＜0的等級，凈增長63.99 km2，總變化率10.51%，面積減少最多的是6—8等級，凈減少72.61 km2，總變化率為-40.76%。究其原因是由于人類活動的影響，城市的擴展，人類活動改變原本的生態環境，城市綠地和不透水面不斷增加，使生態環境向好發展，適應人類生存的需要。從整體的空間分布上看，溫度區分布有從中間向外圍遞增的現象，RLST＜0的區域位于研究區西北部，下墊面多為耕地、水域與城市中心建成區，6—8、＞8的區域多位于東北，西南部，下墊面多為稀疏草地；城市中心建成區RLST＜0呈現向北和西北方向逐年擴大，下墊面多為建設用地，低溫面積逐年擴大。在干旱區綠洲城市特殊的氣候環境下，形成了不同于內地城市的“冷島”現象。

表2 1990—2020年相對地表溫度分區面積統計Table 2 Area statistics of rlative land surface temperature from 1990 to 2018

3.2 相對地表溫度等級重心轉移

利用重心轉移來分析各溫度分區在空間上的變化（圖4）。0—2、2—4、4—6區域的重心均表現為先向東南，再向西，再向東北轉移的趨勢，RLST＞8區域的重心先向西南，再向東北轉移，RLST＜0的區域溫度重心一直北移，空間上呈現等級越高越靠東的現象，RLST＜0分布在中間，RLST＞8區域的重心轉移距離最長，0—2區域的重心轉移距離最短。這與其處于干旱區綠洲的環境下，溫度區分布有從中間向外圍遞增的現象以及城市發展方向、城市建成區面積的擴大關系密切。

圖4 1990—2018年地表溫度分級重心轉移圖Fig. 4 Land surface temperature type center of gravity shift chart from 1990 to 2018

3.3 地表溫度空間集聚分析

基于RLST等級分布圖（圖3），將RLST＞2的區域定義為區域熱島（Qiu et al.，2017；Santamouris et al.，2018）。利用ArcGIS軟件提取RLST＞2的區域，提取RLST＞2的區域，計算NNR（表3）。由表3可知，1990—2018年的NNR均小于1（P＜0.01），說明研究區內的LST在空間上呈現出很強的集聚效應，其原因在于烏魯木齊特殊的地理環境，城市周圍稀疏草地占地面積廣，植被覆蓋度低。

表3 1990—2018年平均最近鄰比率Table 3 Average near neighbor ratio from 1990 to 2018

進一步分析不同時期不同RLST等級下的NNR值（圖5），橫向看，2000、2018年研究區內的NNR值隨著RLST等級的遞增而減小，1990、2010年的 NNR值先增大后減少，其從整體上均表現出等級越高集聚性越強的現象；縱向看，在 2—4、4—6等級上，NNR值均呈現先減少后增大，2000年的NNR值低，集聚性最強，2018年的NNR值最高，聚集性略有降低。在6—8、＞8等級上，NNR值呈現先增大再減少的趨勢，1990年的NNR值最低，集聚性最強；2010年的NNR值最高，聚集性略有降低。在 6—8、＞8等級上 2010—2020年的NNR值有下降趨勢，且在2018年各階段上最低，為 2018年集聚程度最高的區域，因此，在未來城市的發展過程中應該重點關注此區域。

圖5 1990—2018年不同RLST等級的NNR值Fig. 5 NNR of different RLST grades from 1990 to 2018

利用 ArcGIS空間統計工具的多距離空間聚類分析功能，將設置空間尺度為 1 km，置信區間為90%，對LST進行L(d)函數計算（圖6）。將LST在空間上出現隨機分布即L(d)=0時的空間尺度定義為臨界值（高超等，2015）。由圖6可看出，1990年空間尺度小于臨界值18 km時，L(d)＞0，LST呈現集聚趨勢，大于臨界值18 km時，L(d)＜0，LST呈現分散趨勢；2000年空間尺度小于臨界值 15.1 km 時，L(d)＞0，LST呈現集聚趨勢，大于臨界值15.1 km時，L(d)＜0，LST呈現分散趨勢；2010年空間尺度小于臨界值13.7 km時，L(d)＞0，LST呈現集聚趨勢，大于臨界值13.8 km時，L(d)＜0，LST呈現分散趨勢；2018年空間尺度小于臨界值 13.5 km 時，L(d)＞0，LST呈現集聚趨勢，大于臨界值13.5 km 時，L(d)＜0，LST 呈現分散趨勢；1990—2018年空間尺度臨界值逐年下降，從15 km降低到13.5 km，L(d)值最大從 2.22降低到 1.12，表明研究區LST的空間集聚范圍縮小，集聚程度降低。這主要是由于城市建設過程中溫度較低的區域從城市中心向外圍蔓延的原因。近年來，烏魯木齊市在“南控，北擴，西延、東進”的中心城區空間發展策略指導下，建設用地發展向北發展（《烏魯木齊市城市總體規劃（2014—2020年）（2017年修訂）》），建設用地面積大幅度增長，為搭建適合人類生存的城市生態，其過程中對自然環境的改變較為強烈，LST的集聚程度不斷降低，對區域熱島具有一定的緩解作用，使城市人居環境不斷向好發展。

圖6 1990—2018年 RLST 的 Ripley’s L(d)函數Fig. 6 Ripley’s L(d) functions of RLST from 1990 to 2018

3.4 不同土地利用/覆蓋的貢獻度指數時空差異

獲取1990、2000、2010、2018年中國科學院土地利用/覆蓋分類，將4期土地利用/覆被圖做轉移矩陣（表4）發現，近 28年來，研究區內土地利用/覆被面積增加的為建設用地與水域，其余類型均在減少。面積增長最大的為建設用地，從 237.17—658.79 km2，凈增長 421.63 km2，總變化率達到177.78%；土地利用/覆被面積減少最大的耕地面積，從603.53 km2減少到410.36 km2，凈減少193.17 km2，總變化率為 32.01%；未利用地在 2010年已經不存在，說明城市的快速擴展優先占據了未利用土地。其中，在2000—2010年期間，耕地面積急劇較少，2010年之后，由于耕地紅線的劃定，烏魯木齊做出土地利用總體規劃的調整，遏制了耕地減少。

表4 1990—2018年土地利用/覆蓋類型面積變化統計Table 4 Area statistics of land use/cover types from 1990 to 2018

以各年度土地利用/覆蓋類型為單元，計算各類型在不同年度的熱貢獻指數（圖7）。各土地利用/覆蓋類型的熱貢獻度從大到小表現為草地，建設用地，未利用地，林地，水域，耕地，其中草地，建設用地，未利用地為正貢獻，林地，水域，耕地為負貢獻。耕地由于其面積占比大，且獲取影像時間多為農作物的生長季，大面積的農作物的光合作用吸收并遮擋太陽輻射，致使地表保持較低的溫度，成為CI值最低的地類；2018年耕地的CI絕對值最低是耕地面積減少而造成。較少的面積林地、水域是導致其 CI并不很低的原因，林地隨著面積的逐年較少，CI整體呈現增大趨勢，水域也是隨著其面積的變化而變化。草地 CI值最高是由于研究區內的草地面積廣，且多為地覆蓋度草地，地表裸露度較高。建設用地的CI值整體并不是很高，其在2000年達到最高 0.251，其后逐年下降，這與人類活動對城市的改造有密切的關系，1990—2000年之間只是追求城市發展速度的發展，而忽略了城市發展的生態環境，造成建設用地 CI升高的現象。未利用地CI值低與面積少有直接的關系。

圖7 1990—2018年各地類貢獻度指數Fig. 7 Contribution index of different land use/cover types from 1990 to 2018

4 結論與討論

4.1 結論

基于1990—2018年的4期Landsat衛星影像，通過地表溫度反演、重心轉移模型、空間聚類分析以及計算貢獻度指數等方法，定量研究烏魯木齊市城市化進程中地表溫度的空間格局變化及其在空間上的集聚效應。主要結論如下：

（1）1990—2018年間，研究區內相對地表溫度等級發生了較為顯著的變化。相對地表溫度等級在＜0、0—2、＞8等級上面積增加，2—4、4—6、6—8等級上面積減少。面積增長最大的是相對地表溫度等級＜0的等級，面積減少最多的是6—8等級，從空間分布上看，溫度區分布有從中間向外圍遞增的現象。

（2）1990—2018年間，研究區內相對地表溫度等級在 0—2、2—4、4—6區域的重心均表現為先向東南，再向西，再向東北轉移的趨勢；相對地表溫度等級＞8區域的重心先向西南，再向東北轉移；相對地表溫度等級＜0的區域溫度重心一直北移。空間上呈現等級越高越靠東的現象，相對地表溫度等級＜0分布在中間，相對地表溫度等級＞8區域的重心轉移距離最長，相對地表溫度等級在 0—2區域的重心轉移距離最短。

（3）研究區內的LST在空間上呈現出很強的集聚效應，從整體上均表現為等級越高集聚性越強的現象，在RLST≥6等級上NNR值有下降趨勢，在未來城市的發展過程中應該重點關注此區域；1990—2018年間，區域熱島的集聚程度在不斷降低，空間尺度臨界值逐年下降。

（4）1990—2018年間，研究區內土地利用/覆被面積增加的為建設用地與水域，其余類型均在減少，面積增長最大的為建設用地，面積減少最大的耕地面積。各土地利用/覆蓋類型的熱貢獻度從大到小表現為草地、建設用地、未利用地、林地、水域、耕地，其中草地、建設用地、未利用地為正貢獻，林地、水域、耕地為負貢獻。CI值大小與面積密切相關。

4.2 討論

根據上述研究，給出以下 3點建議：（1）合理配置林地面積，提倡在低植被覆蓋區植樹造林，發揮林地對熱環境的緩解作用同時較少低植被覆蓋區的面積；（2）嚴格落實耕地紅線，完善耕地保護法律法規體系，延緩耕地的減少速度；（3）在城市未來的發展過程中應該重點關注RLST≥6的區域。

文中僅從宏觀角度，分析土地利用/覆蓋與地表溫度的時空變化，對城市建成區內部土地利用類型、城市功能對城市地表溫度的影響并未考慮，此后的研究應該更多地從城市內部結構、建筑高度和建筑密度等城市空間形態參數方面探討地表溫度的響應。