西藏地表溫度時空演變特征及影響因子

伍健恒， 孫彩歌， 樊風雷，2

（1.華南師范大學地理科學學院，廣東 廣州 510631；2.西藏大學高原地表遙感監測聯合實驗室，西藏 拉薩 850000）

0 引言

由于人類活動以及自然氣候波動的加劇，1880—2012年全球平均地表溫度（land surface temperature，LST）升高了0.85℃，1951—2012年全球平均LST升溫速率為0.12℃·（10a）-1，約為1880年以來的2倍，大量的觀測事實和實例研究證明，全球氣候呈現暖化的趨勢［1-2］。LST包括了地表水體表面溫度和陸地表面溫度［3］，它在地-氣系統間的能量交換和水熱平衡過程中扮演著重要角色，對水循環、植被生長、氣候變化產生深刻影響，是地表過程分析與模擬的關鍵參數［4］。

熱紅外遙感是獲取大區域尺度地表熱力信息的重要手段，能夠基于熱輻射平衡方程反演出高精度的LST。針對不同特性的熱紅外遙感數據，多種LST反演算法相繼被提出，主要包括單通道算法、分裂窗算法以及溫度發射率分離算法等［5-7］。MODIS與ASTER具有多個熱紅外波段，學者們分別利用分裂窗算法和溫度發射率分離算法生產了MODIS LST和ASTER LST地表溫度產品。隨后，學者們利用不同的驗證數據，如紅外溫度計數據、同步探空數據、氣象站氣溫數據以及同步高分辨率遙感影像開展了大量LST產品適用性評價工作［8］。已有研究指出，MODIS LST可以獲得區域尺度上溫度場變化，適用于熱力差異分析，而ASTER LST則能更精細地反映溫度變化的細節特征，適合于小區域的對比研究和精度驗證［9］。

由于MODIS LST覆蓋范圍廣、適用性強，學者們開展了以MODIS LST為數據基礎的相關研究，得到豐碩的成果。喬治等［10］利用MODIS LST從數量、形狀以及結構維度揭示北京地表熱力景觀時空分異；饒勝等［11］基于MODIS數據探究了珠三角地區的熱島分布，并且發現MODIS LST能較好反映城市化進程中區域熱環境的變化；趙冰等［12］以MODIS LST為數據源，系統分析了近15年中國LST的時空分布和變化格局，發現中國LST整體呈現微弱上升的趨勢，且暖化趨勢顯著不均；也有學者嘗試利用MODIS LST和氣象站點數據構建LST與近地氣溫關系模型，模型結果能夠較好適用于其研究區域［13-14］。以上研究說明基于MODIS LST的大區域尺度熱環境監測切實可行，但是具有長時間跨度、高時間分辨率，以中國西部生態脆弱地區為對象的LST時空演變研究仍然鮮見。

青藏高原是國內外學者研究的熱點地區，以往對其LST的時空分析主要依賴地面站點展開。胡軍等［15］利用西藏10個地面站點的逐月平均地溫資料，采用趨勢分析的方法研究了1971—2005年西藏近35年的年、季平均LST變化趨勢、氣候突變以及異常年份；卓嘎等［16］利用西藏地區38個觀測站點的0 cm地溫資料，采用經驗正交函數方法探討1980—2004年西藏年平均及不同季節地溫異常時空結構特征。隨著熱紅外遙感技術的發展，MODIS、GSM-5和AATSR等遙感影像成為了青藏高原LST研究的主要數據來源［17-19］，其中，MODIS憑借其相對較高的空間分辨率成為最為關鍵的數據源。楊成松等［20］基于MODIS LST數據通過時空信息融合的方法對高原地區缺失像元進行重構，并建立以年為周期的余弦函數模型來描述青藏高原LST在各年份的季節波動；Xu等［21］考慮到云層污染造成MODIS LST數據缺失的問題，采用時間序列諧波分析的方法進行LST時序影像的重構，以此為基礎探究了青藏高原LST的年際變化、峰值周期以及季節性波動。

大量研究利用普通線性回歸（ordinary least square，OLS）或者地理加權回歸（geographically weighted regression，GWR）對LST的影響因子進行探析，表明LST受多種因素的綜合影響，與地形、地表覆蓋、濕度、夜間燈光以及風環境等因素密切相關［22-26］。文路軍等［27］利用逐步回歸分析構建了LST與DEM、坡度、坡向以及其他地表參數的回歸方程，探究了向陽面與背陽面下各地表參數對LST影響程度的差異；王佳等［28］基于Landsat土地覆蓋數據與MODIS LST數據采用OLS與GWR分析了林地、草地、濕地、耕地以及人工表面與LST的關系。然而，由于LST及其影響因子存在明顯的空間異質性，傳統的方法難以識別影響因子的空間尺度差異，捕捉的大量噪聲會造成較大的估計偏誤。

綜上所述，MODIS LST是大區域尺度熱環境演變研究的主要數據源。盡管已有研究表明MODIS LST在青藏高原有較好的適用性［29］，但是高原地區遙感影像普遍面臨云遮擋問題，年平均云量超過50%［30］，這大大降低了遙感數據的可用性。有限的數據和衛星傳感器本身的不確定性導致了LST的趨勢分析結果存在差異［31］，難以對高原地區進行大范圍長時間跨度的有效監測。此外，對LST影響因子的研究主要基于OLS與GWR，由于忽略了影響因子的空間作用尺度會造成回歸結果的不穩健。針對前者，近年來發展的谷歌云計算平臺（Google Earth Engine，GEE）通過云檢測以及影像疊加提取均值的方式可以獲取無云或者少云的遙感影像，解決了云層遮擋、數據缺失的問題［32］；針對后者，多尺度地理加權回歸（multi-scale geographically weighted regression，MGWR）通過允許帶寬具有特異性從而量化了影響因子作用尺度，解決了OLS和GWR無法識別影響因子作用尺度的問題，使其回歸結果更為可靠［33］。因此，本文基于GEE平臺采用均值合成的方法生成無云覆蓋的長時間序列MODIS LST遙感數據，利用趨勢分析、熱力空間分析以及熱力重心遷移等方法量化西藏LST年際變化特征及空間格局演變，并且選取高程（digital elevation model，DEM）、歸一化植被指數（normalized difference vegetation index，NDVI）、裸土指數（bare soil index，BI）、垂直不透水面指數（perpendicular impervious surface index，PISI）和纓帽變換的濕度分量（WET）作為LST影響因子，結合MGWR分析LST影響因子的作用效力與作用尺度。

1 研究區與數據

1.1 研究區概況

西藏地處青藏高原，地理范圍為26°50′～36°53′N、78°25′～99°06′E，東西橫跨2 000 km，南北相距1 000 km，平均海拔在4 000 m以上，三維地帶性分異顯著，由西北至東南分別為高原寒帶、高原亞寒帶、高原溫帶、亞熱帶和熱帶。西藏年日照時數1 443.5～3 574.3 h，各地年日照時數差異較大，年平均氣溫-2.4～12.1℃，氣溫由東南向西北遞減。西藏復雜且獨特的自然條件造成了LST的區域差異，為了細化分析西藏的熱環境特征，本文按照王小丹等［34］對西藏生態區的劃分，將研究區劃分為了7個生態區（圖1），各生態區體現不同的氣候和地貌特征。

圖1 西藏生態分區Fig.1 Ecological zoning in Tibet

1.2 數據來源與處理

本文借助GEE平臺獲得了551景的MOD11A2地表溫度產品，該產品是MOD11A1逐日地表溫度產品8天數據的算術平均。MOD11A1是基于分裂窗 算 法 利 用31（10.78～11.28 μm）和32（11.77～12.27 μm）熱紅外波段計算得到，已有研究表明，在地形平坦，地物均一且晴空的條件下，MODIS LST產品誤差低于1℃［35-36］。本文選取MOD11A2中的日間地表溫度數據利用QC質量控制波段識別云層并對數據進行云掩膜，最后通過影像疊加以及均值合成的方法［式（1）］生成覆蓋全區、無云遮擋的2000—2020年3—9月日間LST均值數據。此外，本文利用27景的MOD09A1影像分別計算了2020年3—9月均值NDVI［37］、WET［38］、PISI［39］和BI［40］，并且獲取集成于GEE平臺的30 m空間分辨率SRTM數字高程數據（digital elevation model，DEM）。為統一地表參數與LST的空間分辨率，將所有數據的空間分辨率重采樣為1 000 m。

式中：LSTmean為日間LST均值；LSTi為第i張影像的日間LST值；n為影像張數。如果某一景影像某區域因云層覆蓋被掩膜，則該景的該區域像元不參與均值運算。

在沒有安裝紅外溫度測量設備和輻射測量儀器的區域，或無法得到長時間和大范圍觀測數據時，可以利用氣象站的氣溫數據對LST精度進行相對評價。本文選用2000—2020年西藏地區22個地面氣象觀測站點的3—9月平均氣溫數據對MODIS LST進行適用性評價，數據來源于中國氣象局國家氣象信息中心和NOAA數據搜索平臺（部分站點于2018年后缺失數據）。

2 研究方法

2.1 地表溫度趨勢分析

本文利用基于像元尺度的時間序列分析法探究研究區域LST的變化特征，采用趨勢線性斜率分析2000—2020年LST的演變趨勢，計算公式為

式中：Slopeij為第i行第j列像元的趨勢線斜率；n為時間序列長度；為第t年第i行第j列像元的LST。Slopeij＞0表示LST呈現上升趨勢，反之則為下降趨勢。

2.2 熱力空間分析

在熱環境研究中融入景觀生態研究的方法能夠為深入認識熱環境空間格局及其演變特征提供新的研究思路［41］。本文分別選取聚集度指數與景觀形狀指數對西藏熱環境空間格局進行定量描述。兩種景觀格局指數計算公式為

式中：AI為聚集度指數；gi為斑塊類型i的同類相鄰的像元數；maxgi為斑塊類型i的同類相鄰的最大鄰接數。AI越小表明聚集度越低，破碎度越高；AI越大表明聚集度越高，破碎度越低。

式中：LSI為景觀形狀指數；E為所有斑塊的邊界總長度；A為景觀總面積。景觀形狀指數用來衡量景觀形狀的復雜程度，其值越小表明景觀結構越緊湊、簡單，其值越大表明景觀形狀越復雜、扁長。

2.3 熱力重心遷移

熱力重心的轉移表征熱力斑塊的空間變化，為了深入分析LST的時空分布特征，采用熱力重心來描述LST的空間遷移過程。熱力重心坐標可以表達為

式中：X與Y分別為LST重心的橫坐標與縱坐標；wi為第i個柵格的LST；xi和yi分別為第i個柵格幾何中心的橫坐標與縱坐標。

根據力矩物理學原理，在獲得熱力重心基礎上計算各生態區對重心遷移的貢獻率C，以此量化生態區對熱力重心遷移的貢獻［42］。首先假設西藏整體熱力重心遷移軌跡為OA，遷移角度為∠NOA；其次將某一生態區的熱力重心遷移軌跡平行投影至O為原點的位置，投影后的軌跡記為OB，遷移角度記為∠NOB；最后過點B作OA的垂線交于點D，OD為OB對OA的貢獻率（圖2）。遷移貢獻率可以表示為

圖2 熱力重心遷移貢獻率計算Fig.2 Calculation of contribution rate of thermal centroid migration

2.4 多尺度地理加權回歸

MGWR允許各個變量的帶寬具有特異性，進而識別不同變量的影響尺度，能夠得到更為可靠的估計結果。MGWR可以表達為

式中：yi為第i個點的因變量；k為變量個數；β為第j個變量的回歸系數；bwj為第j個變量回歸系數所對應的帶寬；（ui，vi）為第i個點的坐標；xij為在第i個點的第j個變量；εi為第i個點的殘差。MGWR為廣義加性模型，使用后退擬合算法對各個平滑項進行擬合，初始化系數估計是由GWR回歸得到，并采用經典的殘差平方和變化比例作為收斂準則。MGWR模型的核函數為二次核函數，帶寬設定準則為校正的赤池信息準則AICc。

3 地表溫度時空演變特征及影響因子

3.1 地表溫度產品適用性評價

本文基于GEE平臺通過影像疊加提取均值的方式獲取無云干擾的西藏2000—2020年LST影像（圖3）。理想情況下，同一地點的LST與氣溫維持著穩定的溫度差距，如果在長時間序列中，LST與氣象站點氣溫數據始終保持相對固定的差距，兩者差距的波動較小，可以認為LST相對準確。為了對MODIS LST進行產品適用性檢驗，借助西藏地區22個地面氣象站點3—9月氣溫均值數據，計算得到地面氣象站點地-氣溫度差值（圖4），發現改則、當雄、拉孜、洛隆、左貢等站點地-氣溫度差值波動較大，阿里、普蘭、聶拉木、江孜、察隅等站點地-氣溫度差值波動較小。利用標準差的方法量化地-氣溫度差值的波動程度。標準差STD越小，越趨近于0，表明地-氣溫度差值波動程度越小，MODIS LST精度越高。表1顯示，各站點STD范圍在0.64～3.83℃，STD平均值為1.93℃，表明平均地-氣溫度差值變化幅度較小，MODIS LST具有適用性。

表1 西藏不同氣象站地-氣溫度差值的多年平均值Table 1 Multi-year average value of ground-air temperature difference in different meteorological stations of Tibet

圖3 2000—2020年西藏LST分布Fig.3 Distribution of LST in Tibet from 2000 to 2020

圖4 西藏不同氣象站地-氣溫度差值的年際變化Fig.4 Interannual variations of ground-air temperature difference in different meteorological stations of Tibet

3.2 地表溫度時間變化分析

由LST年際變化（圖5）可以看出，2000—2020年西藏LST均值緩慢上升，從2000年的18.72℃增高至2020年的20.28℃，增幅1.56℃，年均增幅0.09℃，整體上西藏LST呈現微弱上升的趨勢。

圖5 2000—2020年西藏LST年平均值的年際變化Fig.5 Interannual variation of annual mean LST in Tibet from 2000 to 2020

從全區角度看，2000—2020年西藏LST年際變化規律為波動上升—持續下降—迅速上升—大幅下降—小幅上升。其中，2000—2009年LST由18.72℃波動上升至21.77℃，年均增長0.34℃；2009—2012年LST出現持續下降的現象，年均下降0.91℃；2012—2015年LST迅速回升，由19.05℃上升至21.29℃，增幅達到2.24℃，年均增長0.75℃；2015—2017年LST急劇降低，由21.29℃大幅下降至18.53℃，年均下降1.38℃；2017—2020年LST呈現出小幅上升的趨勢。

從生態區角度看，各個生態區LST年際變化分異顯著。I區LST年際變化規律為波動上升—迅速下降—迅速上升—波動下降。其中，2000—2006年LST由16.13℃波動上升至18.29℃，年均增長0.36℃；2006—2010年呈現迅速下降的態勢，LST以每年0.77℃的速度下降至15.22℃；2010—2013年LST急劇上升，由15.22℃上升至18.48℃，升幅高達3.26℃；2013—2020年LST呈現出波動下降的特征。II區LST年際變化規律為，以2009年為轉折點，2009年以前表現為波動上升，2009年以后表現為波動下降，特別的是2009—2010年出現了急劇下降，由2009年的15.64℃下降至2010年的12.46℃。III區LST年際變化規律為波動上升—波動下降—迅速上升—波動下降。其中，2000—2007年LST由16.43℃波動上升至19.91℃；2007—2014年LST波動下降，年均下降0.54℃；2014—2016年LST大幅上升，升溫幅度高達3.12℃，年均增溫1.56℃；2016—2020年LST表現為波動下降，由2016年的19.25℃降至2020年的16.42℃。IV區與V區LST年際變化規律相似，2009年與2015年出現峰值，呈現出雙波峰形態，IV區兩個LST峰值分別為23.34℃與22.15℃，V區分別為24.06℃與24.78℃。VI、VII區均表現為震蕩變化，其中VI區LST表現為震蕩上升，VII區LST圍繞其多年均值24.71℃上下波動。

LST趨勢分析結果（圖6）顯示，2000—2020年西藏LST增溫趨勢呈現西北高、東南低的分布格局。從生態區角度看，I區和II區Slope分別為-0.30～0.29與-0.33～0.42，表明區域LST年際變化差異較小。兩區均值分別為0.01與0.03，LST表現為微弱上升的趨勢。III區和IV區Slope分別為-0.14～0.36與-0.59～0.87，反映出前者區域LST年際變化差異較小，后者年際變化差異較大。兩區均值均為0.06，說明III區和IV區具有相同的升溫態勢。V區 和VI區Slope分 別 為-0.88～0.36與-1.11～0.28，均表現出明顯的LST年際變化差異。兩區均值分別為0.08與0.04，前者表現出顯著的增溫趨勢，后者表現出一定程度的升溫傾向。VII區Slope取值范圍在-0.70～0.22之間，均值為-0.02，說明VII區LST年際變化差異較大，整體上表現出微弱的降溫態勢。

圖6 2000—2020年西藏LST變化趨勢Fig.6 Variation trend of LST in Tibet from 2000 to 2020

3.3 地表溫度空間變化分析

本文在獲取西藏LST基礎上，對LST進行了歸一化處理［式（9）］，然后以0.2為間隔進行等分［43］，將西藏熱力空間劃分了5個溫度分區，分別為低溫區（0≤NLST＜0.2）、次低溫區（0.2≤NLST＜0.4）、中溫區（0.4≤NLST＜0.6）、次高溫區（0.6≤NLST＜0.8）以及高溫區（0.8≤NLST≤1）。

式中：NLST為歸一化后的LST值；LSTi為第i個像元的LST值；LSTmin為LST的最小值；LSTmax為LST的最大值。

西藏LST每個年份具有相似的熱力空間分布（圖7）與熱力空間結構（圖8）：熱力空間分布表現出西北高、東南低的特征；熱力空間結構主要以中溫區及以上為主，約占研究區域面積的70%。低溫區主要集中在II區；次低溫區零散分布于II與III區；中溫區聚集在III區和V區的北部區域；次高溫區鑲嵌于各個生態區；高溫區廣泛分布于VII區、V區西部以及IV區南部。對比2000—2020年不同溫區面積占比變化可以發現：面積比例增加最為明顯的是高溫區，由18.79%上升至20.75%，增長幅度為1.96%；次高溫區面積占比增長僅次于高溫區，增長幅度為1.02%；低溫區面積占比略有增加，由13.41%上升至14.09%；次低溫區與中溫區面積占比呈現下降趨勢，分別減少了1.11%和2.55%。

圖7 2000—2020年西藏熱力空間分布Fig.7 Thermal spatial distribution in Tibet from 2000 to 2020

圖8 2000—2020年西藏熱力空間結構Fig.8 Thermal spatial structure in Tibet from 2000 to 2020

進一步計算景觀格局指數AI和LSI，用以描述LST的空間形態變化（圖9）。整體上看，AI由小到大排序依次為次低溫區、中溫區、次高溫區、低溫區、高溫區，說明高溫區連片分布，空間聚集程度最高，而次低溫區空間分布破碎，空間聚集程度最低。從變化趨勢上看，2000—2020年各個溫區AI均呈現上升的趨勢，說明各溫區分布趨于集中。其中，低溫區與次高溫區AI變化最為明顯，由2000年的88.46和82.24增加至2020年的89.59和83.49，AI分別增加了1.13和1.25，表現出了較強的空間聚合傾向。

圖9 2000—2020年聚集度指數與景觀形狀指數的年際變化Fig.9 Interannual variations of aggregation index（AI）and landscape shape index（LSI）from 2000 to 2020

LSI由小到大排序依次為低溫區、高溫區、次高溫區、次低溫區、中溫區，其中低溫區與高溫區LSI遠遠小于其余三個溫區，說明低溫區與高溫區形狀較為簡單、規則，而次高溫區、次低溫區、中溫區形狀較為復雜、扁長。從變化趨勢上看，2000—2020年高溫區LSI由51.17波動上升至51.73，漲幅僅為0.56，表明高溫區形狀變化穩定。低溫區、次低溫區、中溫區和次高溫區LSI分別下降了3.73、7.93、7.83以及5.63，均表現出波動下降的趨勢，其形狀在一定程度趨于簡單化與規則化。

整體上看，低溫區與高溫區空間聚集程度較高，形狀簡單、規則；次低溫區、中溫區以及次高溫區空間分布破碎，形狀復雜。其原因是高溫區主要位于人工硬化地表等城市不透水面和戈壁荒漠，低溫區主要位于高山積雪區域，兩個溫區地表覆蓋均呈現連片分布特點，聚合程度較高，形狀較為規則，其余溫區主要分布于河流、湖泊的邊緣或者鑲嵌分布于荒漠戈壁，形狀復雜且破碎。

3.4 熱力重心遷移分析

本文從熱力重心遷移的角度刻畫2000—2020年西藏LST時空演變過程（圖10），并對比分析各生態區對熱力重心遷移的貢獻程度（圖11）。

由圖10可知，各溫區重心的標準差橢圓較為扁平，表明低溫區至高溫區重心分布具有明顯的方向性，觀察各年份熱力重心位置可以發現各溫區重心遷移軌跡具有顯著差異。2000—2017年低溫區重心以每年3.98 km的速度向西北移動，2017—2020年低溫區重心呈現快速向東南移動的態勢，遷移速度達22.24 km·a-1；次低溫區重心遷移軌跡與低溫區相似，2000—2012年以每年9.43 km的速度向西北遷移，2012—2020年向東南移動，移動速度為20.79 km·a-1；2000—2005年中溫區重心以7.59 km·a-1的速度持續向西北移動，2005—2020年其遷移軌跡轉向東南方向，速度下降至5.33 km·a-1；2000—2016年次高溫區重心整體呈現向東南移動的趨勢，年均移動距離3.72 km，2016—2020年以每年17.54 km的速度持續向西北轉移；高溫區重心遷移軌跡與次高溫區類似，2000—2017年以每年9.61 km的速度向東南移動，2017—2020年高溫區重心向西北轉移，年均轉移距離63.84 km；全區熱力重心分布格局與遷移軌跡與高溫區高度相似。綜合以上分析可知，低溫區重心與高溫區重心遷移軌跡呈現出由相向而行到背向而行的轉變，結合前文發現LST空間分布具有東南低西北高的特征，說明了研究區東西部區域LST差距經歷了由縮小到擴大的過程。

圖10 2000—2020年西藏LST重心遷移Fig.10 LST centroid migration in Tibet from 2000 to 2020

逐年計算各個生態區的熱力重心遷移貢獻率，求取平均得到2000—2020年各生態區平均熱力重心遷移貢獻率（圖11）。不同生態區對熱力重心遷移的貢獻差異明顯：低溫區重心遷移主要受V區與VI區的影響，兩區貢獻率均在60%以上；次低溫區重心遷移主要受到V區正向影響以及I區負向影響；中溫區重心遷移幾乎由III區與V區的熱環境狀況控制，兩區貢獻率均高達150%；次高溫區重心遷移貢獻最大的兩個生態區是II區與VI區，貢獻率分別為54.77%和45.36%；高溫區重心遷移受到V區與VI區的影響較為明顯。西藏全區熱力重心貢獻率由低到高排序分別為VII區＜I區＜III區＜V區＜IV區＜VI區＜II區，其中II區、IV區、V區和VI區為正貢獻，I區、III區和VII區為負貢獻。總體來看，IV區、V區以及VI區對不同溫區重心遷移均產生較強的正向影響。I區、II區、III區以及VII區分別對特定溫區的重心遷移具有一定程度的負向影響。

圖11 2000—2020年西藏各生態區平均熱力重心遷移貢獻率Fig.11 Average contribution rate of thermal centroid migration in various ecological zones of Tibet from 2000 to 2020

3.5 地表溫度影響因子解析

以2020年數據為例，利用漁網工具在研究區內生成等間距的6 355個樣本點，在進行數據標準化基礎上利用MGWR對LST進行影響因子解析（圖12），揭示各個因子對LST的影響程度以及尺度差異（P＜0.05）。由表2可知，與OLS相比，MGWR與GWR具有更高的R2以及更低的AICc。與GWR相比，MGWR具有更低的AICc，說明考慮到影響因子尺度差異的MGWR優于OLS與GWR。

表2 MGWR與GWR、OLS模型指標對比Table 2 Comparison of model index among MGWR，GWR and OLS

圖12 LST多尺度地理加權回歸系數空間分布Fig.12 Spatial distribution of LST MGWR coefficients

由表3可知，常數項帶寬為43，表示捕獲鄰近43個樣本用作模型估算，占樣本數量的0.68%，作用尺度最小。常數項系數在空間上呈現出西高東低的分布格局，系數范圍在-1.27～0.89，在不同生態區具有明顯的影響性質差異：在III、IV、V、VI、VII區表現出明顯的增溫效應，而在I、II區表現為較強的降溫作用。從系數均值絕對值上看，常數項的影響強度僅次于DEM。

表3 多尺度地理加權回歸系數統計Table 3 Statistical description of MGWR coefficient

BI、PISI和NDVI帶寬分別為394、293和290，占樣本數量的4.56%～6.20%，具有相近的作用尺度。BI與NDVI回歸系數均表現為東南高、西北低的分布特征，PISI回歸系數大致呈現南高北低的趨勢。BI、PISI和NDVI的 回 歸 系 數 分 別 為0.09～0.68、-0.27～0.46和0.02～0.56，對LST均存在正向的影響。從系數均值絕對值上看，其影響強度居中。

WET帶寬為878，占樣本數量的13.82%，具有最大的作用尺度，系數在空間上較為平穩。WET回歸系數大部分在V區顯著，系數由南向北逐漸減小。WET回歸系數為-0.26～0.14，對LST整體表現出微弱的負向影響。從系數均值絕對值上看，其影響強度最弱。

DEM帶寬為52，占樣本數量的0.82%，作用尺度相對較小，僅次于常數項。DEM回歸系數由東南向西北呈現減小趨勢。DEM回歸系數范圍為-1.94～0.27，對LST產生極為顯著的負向影響。從系數均值絕對值上看，DEM對LST的影響強度在所有因子中最大。

4 結論與討論

本文基于GEE平臺獲取長時間序列的MODIS數據，探究西藏LST的時空演變特征，并且利用帶有統計推斷的MGWR對LST影響因子的空間分異和空間尺度差異進行探討，旨在深入了解西藏熱環境特征以及LST的影響機制。通過研究分析得到如下結論：

（1）時間上，2000—2020年西藏LST年際變化規律為波動上升—持續下降—迅速上升—大幅下降—小幅上升。整體上LST呈現波動上升的趨勢，年均增長0.09℃。

（2）空間上，西藏LST空間分布呈現西北高東南低的分布格局，增溫趨勢由東南向西北逐漸增強。LST空間結構以中溫區及以上溫區為主，約占研究區域面積的70%。5個溫區中，低溫區與高溫區空間聚集程度較高，形狀簡單、規則；次低溫區、中溫區以及次高溫區空間分布破碎，形狀復雜。

（3）5個溫區的重心分布均具有明顯的方向性，并且各溫區重心遷移軌跡具有顯著差異。低溫區重心與高溫區重心遷移軌跡呈現出由相向而行到背向而行的轉變，反映出2000—2020年西藏東部與西部區域LST差距經歷了由縮小到擴大的變化過程。

（4）LST影響因子具有不同的作用尺度，其尺度由小到大分別為常數項、DEM、NDVI、PISI、BI、WET。DEM和WET負向影響LST、NDVI、PISI以及BI正向影響LST，而常數項在不同的生態區其影響性質存在差異。DEM是影響LST的主要因子，作用效力最強；常數項的影響強度僅次于DEM；NDVI、PISI以及BI對LST的影響強度居中；WET對LST的影響最弱。

值得注意的是，西藏LST年際變化極小值出現在2008年，極大值出現在2009年與2015年，分別對應拉尼娜年與厄爾尼諾年，極值出現的年份反映出高原LST對于極端天氣事件的高度敏感。2008年出現的拉尼娜事件導致高原地區出現了顯著的降水現象，2007—2008年22個氣象站點數據顯示其平均氣溫由10.53℃下跌至10.00℃，下跌幅度為0.53℃，顯著的降水與氣溫驟降造成西藏LST明顯下降；與之相反，2009年與2015年出現了較強的厄爾尼諾事件。2008—2009年22個氣象站點數據顯示其平均氣溫由10.00℃大幅上升至10.86℃，增幅高達0.86℃，氣溫的升高是2009年LST急速上升的主要原因。22個氣象站點數據顯示2014—2015年其平均氣溫由10.24℃微弱上升至10.30℃，上升幅度僅為0.06℃，氣溫并非2015年LST上升的主要原因。已有資料顯示厄爾尼諾年亞洲季風邊緣區夏季風呈現顯著衰退，2015年青藏高原中部出現了顯著干旱［44］，據此推測干旱條件壓迫植被生長環境，造成植被被裸巖裸土等對LST正向影響較高的地表覆蓋所替代，導致LST快速上升，同時這也可能是與歷年相比2015年高溫區占比最大以及具有最大AI的主要原因。另外，熱力重心遷移分析表明東西部區域LST差距經歷了由縮小到擴大的過程。在全球暖化背景下冰雪消融初期，西北部高山積雪融水增加，改善了荒漠區域植被生長環境，原本裸露的土壤、沙礫被植被覆蓋，高溫區重心向東南方向遷移。2006—2016年位于東南部的山南與林芝地區經濟發展迅猛，GDP年均增長12.07%，比其他地區高約0.62%，這是森林砍伐、農業開墾、產業發展帶來的結果［45］。同時，這也導致了城市不透水面的擴張，促使低溫區重心向西北移動，因此該時期東西部區域LST差距正在縮小。然而，隨著高原冰川退縮、積雪完全消融使得巖石裸露于地表。已有學者研究指出，青藏高原荒漠土地面積自東南向西北逐漸增加，荒漠化程度不斷加重［46］，阿里-昆侖區域更是中國五大荒漠化區域之一［47］，西北部的巖漠化加重使得高溫區重心快速向西北移動。而隨著生態文明建設的開展，封山育林、退耕還林使得東南部的林地區域得到有效保護，城市擴張得以有效規劃，城市不透水面無序蔓延得以遏制，低溫區重心逐步向東南遷移。因此，該時期東西部區域LST差距逐漸擴大。此外，MGWR結果顯示常數項在III、IV、V、VI、VII區表現出明顯的增溫效應，而在I、II區表現為較強的降溫作用，原因可能是常數項在一定程度上可以捕獲未納入到MGWR模型中的氣溫、降水、光照等影響LST的因子，一方面由于隨著全球暖化西藏西部地區凍土退化現象加劇，淺層凍土吸收的地表熱量減少導致地表增溫，另一方面是由于西藏東南部地區山高谷深，地形起伏產生的陰影遮蔽地表，并且夏季印度洋西南季風受到喜馬拉雅山脈的阻擋，水汽被迫集聚在西藏東南部形成了豐沛的降水，在陰影與降水的共同作用下造成LST下降。最后，本文發現BI、PISI和NDVI均對LST產生正向影響，可能是干燥的土壤和不透水面比熱容小，具有較快的增溫速率，而高原地區低矮茂密植被的阻礙作用減弱了強風的影響從而減少熱量的散失。

在生態分區基礎上科學識別LST時空演變特征，準確掌握LST影響因子及其機制，是有效監測高原地區熱環境狀況的關鍵，能夠為分區管理與保護高原生態環境，防范生態系統服務功能退化，維護區域生態安全提供科學依據。本文基于GEE生成無云干擾的MODIS LST遙感數據，實現了以往由于云層影響而缺乏的大范圍長時間跨度高原LST監測，并且結合MGWR量化了LST影響因子的作用效力和作用尺度，彌補了過往LST回歸分析中忽略空間異質性以及尺度差異的不足，得到了更為合理的回歸結果。但是本文在地表覆蓋遙感指數選取上存在兩個方面的不足，一方面，利用BI與PISI無法精確區分土壤與自然不透水面，造成一定程度的混淆現象；另一方面，由于現階段缺少能夠恰當表征冰川、凍土以及積雪的遙感指數，未能將一般意義上認為負向影響LST的地表覆蓋納入模型中，對研究結果產生了一定的不利影響。