石河子绿洲地表温度时空变化与城镇发展协调性分析

张雪玲，阿里木江·卡斯木,2,3①，梁洪武

(1.新疆师范大学地理科学与旅游学院，新疆 乌鲁木齐 830054；2.新疆师范大学丝绸之路经济带城镇化发展研究中心，新疆 乌鲁木齐 830054；3.新疆干旱区湖泊环境与资源重点实验室，新疆 乌鲁木齐 830054)

地表温度(land surface temperature, LST)作为气候要素和变量，是探究区域热环境的重要指标，对于地表热量平衡和气-地相互作用具有推动作用[1]。绿洲是干旱-半干旱区的独特景观单元，由于中心绿洲区与外围荒漠区具有不同的下垫面类型[2]，绿洲地区呈现中心温度低、四周温度高的环状分布特点[3]，且绿洲相对于荒漠具有更加湿冷的“冷岛效应”特征[4]。热岛效应会形成区域持续性高温，造成能量的高消耗，给资源可持续利用和人居环境质量带来沉重负担；而冷岛效应在一定程度上能够遏制水体、植被过度的蒸发与蒸腾，这对于干旱地区资源环境大有裨益[5]。但是，随着城市化飞速发展，城镇范围极速扩张，使得绿洲的冷岛强度不断降低，LST的持续升高可能会引发不可抗力的自然灾害[6]。2006年颁布的《国家中长期科学和技术发展规划纲要(2006—2020年)》[7]提出，我国要切实推进全球环境公约，以加强对于极端气候变化的应对能力。“十四五”规划中也明确提出要促进生态文明建设发展，这意味着LST状况已成为改善城市生态环境问题的关键。

当前国内外学者关于地表热环境的研究主要包括以下几个方面。首先，地表下垫面覆被物性质对LST空间分布具有主导作用[8]。通常来说，土地利用强度越高、景观空间分布格局越复杂的地区，其热岛效应越强烈[9]。AYANLADE等[10]基于贡献指数和景观指数探究“汇”“源”景观对城市热岛的贡献率，其结果表明以水体和绿地等为主体的“蓝绿景观”是缓解城市热的重要手段。乔治等[11]采用景观类型比例(PLAND)、景观形状指数(LSI)和聚合度指数(AI)等景观指数揭示景观破碎化程度与斑块形状大小对地表热力景观的影响规律。其次，地表热环境是由多种因素综合构成。ZHU等[12]采用主成分分析，利用多个影响因子降维构建生态环境影响指数，以此来分析LST的协调性和贡献度。沈中健等[13]结合夜间灯光数据采用耦合态势模型和协调性模型探究不同因素与LST的耦合协调关系。此外，地表热环境研究还包括GDP发展水平[14]、人口密度[15]以及海拔和坡度[16]等自然要素对LST的影响。

地理探测器模型适宜于分析类型数据，能够用于独立探测因子的非线性关系，分析连续数据或比值数据时，则需要先进行一定离散化处理。相较于单一离散法，R软件GD包在建模前采用最优参数方法对数据进行离散，该方法处理后的数据仍保留原有特点，探测结果也更准确。以新疆典型绿洲城市——石河子市为研究对象，分析不同时相LST的空间聚集性和空间分异性，揭示不同土地利用/覆被类型以及高程、坡度和夜间灯光数据与LST的相关性，探究城镇发展与LST变化的协调程度，旨在为促进地区生态环境质量稳步提高，加快城镇化建设，建设和谐、友好的低碳型社会提供一定理论依据。

1 研究区概况与数据预处理

1.1 研究区概况

石河子市位于新疆维吾尔自治区天山北麓中段，地理坐标为北纬43°26′～45°20′，东经84°58′～86°24′，地处准噶尔盆地南部，北部、西部和南部均被沙湾市环绕，东部以玛纳斯河为界，与玛纳斯县毗邻相望[17]。石河子市地势南高北低，南部为山地区，故而形成大范围山前冲积平原，地形起伏较小，地势较为平缓。研究区为典型温带大陆性气候区，水源充足，自然资源和物种资源丰富，区内农业较发达，土地开发利用水平较高。

1.2 数据源与预处理

选取2006年7月31日和2011年7月29日Landsat TM影像，以及2017年7月29日和2021年7月29日Landsat TRIS影像，空间分辨率为30 m，成像清晰且云量均小于1%，数据来源于美国地质勘探局(https:∥glovis.usgs.gov/)，并对数据进行辐射定标、大气校正和图像裁剪等预处理工作。夜间灯光数据用于表征人口分布状况，2006和2011年采用DMSP/OLS数据，空间分辨率为1 km；2017和2021年采用NPP/VIIRS数据，空间分辨率为500 m，数据来源于Google Earth Engine (GEE)地理计算云平台(https:∥earthengine.google.com)。数字高程模型(DEM)来源于地理空间数据云平台(https:∥www.gscloud.cn/)，空间分辨率为30 m，并基于DEM提取坡度数据。

2 研究方法

2.1 土地利用分类

采用遥感和Google Earth影像，使用目视解译法对训练样本和验证样本进行采样，样点总数量高于100个。根据研究区现状，将土地利用类型划分为草地、耕地、林地、建设用地、水域和未利用地6种1级类型。

随机森林法由BREIMAN[18]提出，是基于分类决策树的机器学习算法。分类决策树采用Bagging方法，基于选取的训练样本构建每一棵决策树，由全体决策树投票确定最终的分类结果[19]。随机森林法基于原始的全体训练样本M，采用重采样随机有放回抽样选取N个训练样本，再从全体训练样本中随机选取m个分类特征(m

2.2 地表温度反演

2.2.1辐射传输方程法

基于Landsat热红外波段(TM B6和TIRS B10波段)数据，运用辐射传输方程对LST进行反演[13]。

首先，基于归一化植被指数(NDVI，INDV)阈值法计算地表比辐射率(ε)，计算公式为

Pv=[(INDV-INDV，Soil)/(INDV，Veg-INDV，Soil)]，

(1)

ε=0.004Pv+0.986。

(2)

式(1)～(2)中，Pv为植被在整体混合像元中所占的比例，即植被覆盖度；INDV，Soil和INDV，Veg分别为没有植被(即裸土地区)和完全是植被地区的归一化植被指数值，取值分别为0.05和0.70。

其次，计算由Landsat传感器接收到的大气热红外辐射亮度值(Lλ)，计算过程如下：

Lλ=[εB(Ts)+(1-ε)L↓]τ+L↑。

(3)

式(3)中，ε为地表比辐射率；Ts为地表真实温度，℃；B(Ts)为黑体辐射亮度，W·m-2·μm-1·sr-1；τ为热红外波段大气透过率；L↓为大气向下辐射亮度，W·m-2·μm-1·sr-1；L↑为大气向上辐射亮度，W·m-2·μm-1·sr-1。τ、L↓和L↑数据均由NASA官网(http:∥atmcorr.gsfc.nasa.gov/)获取。黑体辐射亮度计算公式为

B(Ts)=[Lλ-L↑-τ(1-ε)L↓]/τε。

(4)

最后，用普朗克公式获取地表真实温度，计算公式为

Ts=K2/ln[K1/B(Ts)+1]-273。

(5)

式(5)中，对于TM B6波段，K1=607.76 W·m-2·μm-1·sr-1，K2=1 260.56 K；对于TIRS B10波段，K1=774.89 W·m-2·μm-1·sr-1，K2=1 321.08 K。

2.2.2地表温度等级划分

由于不同时相不同年份大气状况不同，因此对反演得到的LST进行归一化处理，并采用“均值-标准差”法将LST分为极高温、高温、中温、低温和极低温5个等级(表1)，并将低温和极低温区定义为冷岛区。“均值-标准差”法利用研究区温度平均值与标准差来划分温度梯度，这在一定程度上可以避免由于时相不同而带来的误差[21]。

Hi=(Ti-Tmin)/(Tmax-Tmin)。

(6)

式(6)中，Hi为第i个经归一化处理的地表温度；Ti为第i个像元温度值，℃；Tmax和Tmin分别为温度最大值和最小值，℃。

表1 地表温度等级划分

2.2.3冷岛比例指数

由于绿洲特殊的自然地理条件，形成独特的绿洲“冷岛效应”。根据徐涵秋等[22]提出的热岛比例指数(URI)，采用马勇刚等[23]构建的绿洲冷岛比例指数(OCRI，IOCR)表示热环境变化强度，计算研究区内冷岛强度。

(7)

式(7)中，m为归一化等级指数，将地表温度划分为5个级别，即m=5；n为绿洲冷岛范围等级数，将低温和极低温区定义为冷岛区，因此n=2；w为权重值，为第i级级值；p为温度等级i所对应LST的百分比；i为低温或极低温区等级值，即为4或5。

2.3 空间自相关分析

空间自相关分析分为全局空间自相关和局部空间自相关，用于解释LST空间聚集特征，其相关程度用Moran′sI来表征[24]。其中，全局空间自相关能够用来分析地块之间的聚集程度，计算公式为

(8)

局部空间自相关可以衡量单个空间地块的相关性，且基于Z-Score的LISA图能够清晰明确地表达局部空间变化特征，计算公式为

(9)

式(9)中，I局部为局部空间自相关指数。

2.4 遥感指数计算

不同遥感指数能够反映不同地物覆被特征，选取归一化植被指数(NDVI)用于表征植被生长状况，归一化建筑指数(NDBI)用于表征建成区密度及建筑用地面积比例，改进归一化水体指数(MNDWI)用于表征水域信息，裸土指数(SI)用于表征研究区自然裸土分布状况，归一化水汽指数(NDMI)用于表征植被含水量(表2)。以5种遥感指数作为影响LST的因子来描述地表土地利用类型与地表热环境之间的相关性[25]。

表2 遥感指数计算

2.5 地理探测器与协调性模型

2.5.1地理探测器

地理探测器是由王劲峰等[26]提出，用于探测空间分异程度，揭示自变量与因变量之间相关性的统计学模型。采用地理探测器中的因子探测和交互探测2个方面解释LST与各影响因子之间的关系。

因子探测用q值表达不同类型影响因子对LST空间分异性的解释程度，计算公式为

(10)

式(10)中，h为地表温度或影响因子的分层或分区；Nh和N分别为h层和全区样本单元数；σh2和σ2分别为h层和全区样本的离散方差。q取值为[0,1]，其值越大，表明影响因子对热环境的解释力越强；反之则越弱。

交互探测是识别两两影响因子共同作用时会加强还是减弱对地表温度的解释力。分别计算2个影响因子的q值〔q(X1)和q(X2)〕，再计算它们相叠加后的q值〔q(X1∩X2)〕，并将q(X1)、q(X2)与q(X1∩X2)进行比较。当q(X1∩X2)max〔q(X1),q(X2)〕时，交互探测结果为双因子增强；当q(X1∩X2)=q(X1)+q(X2)时，交互探测结果为独立；当q(X1∩X2)>q(X1)+q(X2)时，交互探测结果为非线性增强。

2.5.2协调性模型

协调性模型可以清晰反映2个因素之间变化发展的协调水平，并以此反映LST变化与城镇发展的协调性过程[13]，计算公式[27]为

(11)

式(11)中，O为协调性指数；i为LST和NDBI逐个对应像元；Mi和Ni分别为地表温度和城镇发展的年平均增长率。O取值范围为[0,1]，O越接近于1，协调程度越强，将0.9

表3 协调关系类型划分

3 结果与分析

3.1 土地利用变化特征分析

石河子市土地利用分类见图1。2006、2011、2017和2021年4期影像总体精度分别为98.67%、98.56%、93.58%和94.08%；Kappa系数分别为0.972、0.969、0.867和0.873，分类效果较好，可靠性较高。

图1 石河子市土地利用分类

由图1和表4可知，建设用地扩张速度激增，由2006年的83.42 km2增长至2021年的164.05 km2，增加17.63%。由于石河子市地处天山北坡经济带，城镇化发展迅速，逐渐形成以主城区为中心向四周扩张的发展趋势。石河子市林地面积逐年增长，除建设用地和林地外，其余类型面积均呈下降趋势。石河子市农业发达，耕地面积占比最高，但下降速度也最快，由2006年的272.33 km2下降至2021年的204.46 km2，耕地面积的减少主要是由于建成区扩张所造成。草地主要位于石河子市南郊，地形为山岭，多为极低覆盖度植被，面积变化幅度较小，呈不明显下降趋势。石河子市以东侧玛纳斯河为主要水源，玛纳斯河年均径流量逐年下降，故水域面积逐年降低。未利用地面积较小，至2021年，未利用地面积所占比例仅为0.19%。

3.2 地表温度分布特征分析

采用2021年同日MODIS地表温度产品MOD11A1数据对Landsat反演得到的LST数据进行对比验证，将数据空间分辨率统一为1 km，并随机提取200个点位数据进行拟合，结果表明2种数据平均温差为±0.94 ℃，最高温差为±5.67 ℃。如图2所示，2种温度数据拟合效果较好，精度较高。

表4 石河子市土地利用类型面积与比例

由图1和图3可知，与大部分城市高温中心分布在建设用地不同，石河子市高温中心位于南部草地。由于石河子市南部草地植被类型主要为荒漠草原，覆盖度与密度均极低，因此极高温和高温区主要分布在该地区，中温区大多分布在建设用地和低密度覆被的耕地，而低温和极低温区则主要分布在高密度覆被的耕地和水域。对2021年LST进行分区统计发现，面积较大的草地、耕地和建设用地平均LST分别为44、34和40 ℃，呈中间温度低、四周温度高的状况，这是典型绿洲冷岛效应特征。但是，随着建设用地面积增加，中温、高温区面积不断扩大，导致“冷岛”范围不断缩小。

图2 温度验证结果

图3 石河子市地表热环境空间分布

由表5～6可知，2006年极低温区面积占比最高，冷岛指数高达0.53，冷岛效应强烈。2011年，城镇开始快速扩张，极低温区面积大幅下降，并导致中、高温区面积增加，冷岛强度下降，冷岛指数为0.47。2017年，中、高温区面积持续上升，城镇扩张速度达到最快，冷岛强度下降速度也到达峰值，冷岛指数下降至0.39。2021年，中温区面积升为最高，这主要是由于土地利用类型由极低温的耕地转为中温的建设用地，虽然冷岛指数持续下降，但是冷岛强度衰减速度开始减缓。绿洲冷岛效应受到城市扩张速度与强度的剧烈影响，在发展中应重视“蓝绿空间”建设，以缓解“城镇热”的问题。

3.3 地表温度空间自相关特征分析

石河子市LST全局空间自相关Moran′sI值见表7，LST局部LISA图见图4。表7显示，4期全局自相关Moran′sI值分别为0.805、0.731、0.821和0.804，且P值均小于0.001，这表明LST具有显著空间正相关关系，且2017年相关性最高。

由图4可知，低-低聚集区主要分布在耕地，即低温区，区域随着耕地范围发生变化，由南北聚集转为西北聚集。高-高聚集区主要集中在南部草地和中部建设用地，即高温和极高温区。研究区中部高-高聚集区呈向南北方向转移的趋势，这与城镇发展方向一致；南部草地高-高聚集区多年来基本保持不变。高-低聚集和低-高聚集区分布较少，呈零星斑块状分布。耕地植被覆盖度较高，具有显著缓解温度的作用，区域LST较低，因而呈现低-低聚集效应；草地地区植被类型为极低覆盖植被，地面升温较快，因而呈现高-高聚集效应。

表5 地表温度等级面积与冷岛比例指数

表6 地表温度转移矩阵

表7 LST全局空间自相关Moran′s I统计

3.4 基于地理探测器的LST相关性研究

采用R 4.1.2软件地理探测器模型GD软件包探测9个影响因子对LST空间分异性的解释力大小，并对各影响因子进行交互探测。将各影像数据空间分辨率统一重采样至1 km，并进行离散化处理，采用gdm函数得到地理探测器探测结果。

由表8可知，2006—2021年，NDBI、NDMI和SI的q值最大，对LST的解释力最高，表明建筑物密集程度、植被含水量和裸土状况对LST具有较大影响。而坡度、夜间灯光和MNDWI的q值最低，对LST的解释力较弱，这是由于石河子市地势平坦，市区人口约为35万～40万，人口密度较低，且水域面积分布较小，因而对LST的影响较弱。不同因素对LST的解释力各不相同，从而构成了石河子绿洲独特的热环境。各影响因子q值由大到小依次为NDBI、NDMI、SI、NDVI、土地利用、DEM、MNDWI、坡度和夜间灯光。

图4 LST局部空间自相关LISA图

表8 因子探测q值

对不同影响因子进行交互探测(图5)，判断两两影响因子相互作用对LST的解释程度。由图5可知，石河子绿洲LST空间分异性不是由某一单因子控制，两两因素在交互作用下大多为双因子增强，表明两两因子交互对LST的影响远远大于单因素的影响。根据q值大小，4个年份起主导作用的交互因子主要为DEM∩NDMI和DEM∩NDBI，表明海拔与植被水汽和建筑的交互作用对LST的解释力最高。石河子市地处塔里木盆地，平均海拔为450 m，适合人口居住和农业发展，区内农作物以棉花、玉米为主，植被水汽含量高，两者因素结合对LST具有较大影响。此外，4个年份坡度∩夜间灯光和DEM∩夜间灯光均呈非线性增强；坡度、夜间灯光和DEM的q值均小于0.5，对LST的解释力较低，但在交互作用下对LST的影响力却大大加强。在石河子市生产发展过程中，自然因素与人为因素的结合对LST具有较大影响。

DEM为地理高程模型，NDVI为归一化植被指数，NDBI为归一化建筑指数，MNDWI为改进归一化水体指数，SI为裸土指数，NDMI为归一化水汽指数。*为非线性增强；无标记为双因子增强。

3.5 城镇发展与LST协调性分析

NDBI反映建筑物状态，在一定程度上可表示城镇发展，且地理探测器探测因子探测结果表明NDBI对LST的解释力最强。因此，采用协调性模型计算3个时期NDBI与LST的协调发展程度。对NDBI和LST数据进行归一化和格网化处理，将研究区划分为300 m×300 m的网格，提取每个网格单元内NDBI和LST像元平均值，计算两两因素之间年平均增长率的协调性(图6)。

协调型反映城镇发展与LST变化协调增长且速度较快，沥青、混凝土等不透水面扩张的同时会引发LST协同升高，因此协调型区域主要分布在城镇扩张方向；而在原有中心城区，城镇扩展和地表升温幅度比较有限，其协调性水平也相对较低。磨合型反映城镇发展与LST变化处于良好发展阶段，发展速度适中；拮抗型反映城镇发展与LST变化幅度较小，或不显著。从整体来看，3个时期协调型区域稳定在27%左右，变化起伏较小；磨合型区域占比最大，且呈先增后减趋势；拮抗型区域占比最小，且呈先减后增趋势(图6)。

2006—2011年，城镇开始向南北方向扩张，LST快速升高并超前于城镇发展，因此协调超前型区域占主导地位；在原有中心城区，冷岛效应强烈且发展空间有限，因而拮抗滞后型区域所占比例同样较高。2011—2017年，城镇开始向四周扩张，磨合超前型区域开始增长并占据主导地位，同时，由于城市冷岛效应衰减速度达到最快，拮抗型区域占比下降，而磨合滞后型区域占比上升。2017—2021年，城镇发展到达后期，城镇发展与LST均受到限制，LST变化滞后于城镇发展，磨合超前型区域占比达到最低，磨合滞后型区域占比有所上升，且冷岛效应衰减速度放缓，拮抗型区域面积有所上升。石河子市全区建设用地类型占有较大范围，LST受城镇发展影响较大，因此LST与城镇发展的协调性变化与土地利用变化的趋势基本一致。

图6 2006—2021年石河子市归一化建筑指数(NDBI)与地表温度(LST)的协调类型比例

4 讨论与结论

4.1 讨论

绿洲是干旱-半干旱地区空间和地域分异的结果，是一种独特的景观单元，也是维持干旱-半干旱地区人类栖居和繁荣发展的重要场所。由于荒漠与绿洲下垫面性质的差异，其地表热力分布具有非均衡性，中心绿洲低温和外围荒漠高温，造就绿洲地区具有典型的冷岛效应[28]。因此，研究绿洲冷岛效应对于评估绿洲生态环境质量具有积极意义[29]。但是，伴随着城市化进程快速发展，建筑物更加繁密紧凑，城镇扩展速率与日俱增[30]，使得研究区绿洲冷岛强度不断降低，这与高原等[17]的研究结果基本一致。2017年，研究区冷岛指数下降速度最快，2021年，冷岛指数下降速度开始大幅减缓。这是由于近年来，国家对生态环境保护愈发重视，但绿洲生态系统的稳步恢复仍需要较长时间的努力。

一般来说，地表温度空间分布由多方面因素共同决定[31]。地理探测器因子探测结果(表8)发现，研究区归一化建筑指数(NDBI)对LST的解释力最高。PENG等[32]基于最小二乘回归模型发现夏季引起LST变化的主导因子为NDBI，其贡献率也最高，而过渡季节影响力最高的则为NDVI。笔者研究结果与之基本一致。大量研究证明，植被对LST具有显著降温作用，但在干旱地区，植被类型大多为干草类和灌木类荒漠植被，植被呈斑块状分布在裸露地面，叶片面积较小而根茎发达，以此来适应严苛的生存条件。因此，研究区地面升温迅速，地表温度较高，形成特殊的干旱区热力景观。石河子市空间结构分布清晰，城镇发展对LST的影响极大，各时期城镇发展大多超前于LST变化，这与沈中健等[13]的研究结果基本一致。而在城市景观格局复杂多变的地区，NDBI对LST的影响会被削弱。

夏季植物种类丰富，植被覆盖度较高，对LST的影响较为典型，因此笔者仅对夏季LST进行研究，在后续研究中可以研究多季节不同因素与LST的相关性。笔者以多个遥感指数为变量探讨LST驱动因素，在未来研究中可以从景观生态学角度进行探讨，分析景观异质性与LST的关系，从更深层次对城镇发展与LST变化进行探讨。

4.2 结论

(1)随着城镇化水平不断提高，城镇地区不断向四周扩张，研究区建设用地由2006年的70.3 km2增长至2021年的164.73 km2，石河子市林地面积逐年增长，而耕地、草地和水域面积均呈下降趋势。土地利用类型对LST的分布影响较大，极高温、高温和中温区集中分布在建设用地、草地和部分低密度覆被的耕地，并呈逐年增加趋势；低温和极低温区分布在耕地和水域，面积逐年减少。2006、2011、2017和2021年研究区冷岛指数分别为0.53、0.47、0.39和0.35，冷岛强度在2017年下降速度最快，2021年下降速度开始减缓。

(2)2006、2011、2017和2021年LST全局自相关Moran′sI分别为0.805、0.731、0.821和0.804，LST具有显著空间正相关性。高-高聚集区主要集中在高温的建设用地和草地，低-低聚集区则集中在低温的耕地。

(3)各影响因子探测结果q值由大到小依次为NDBI、NDMI、SI、NDVI、土地利用、DEM、MNDWI、坡度和夜间灯光，NDBI对LST的解释力最强。交互探测结果表明，4个年份海拔与NDMI和NDBI的交互作用对LST的影响占主导地位。

(4)2006—2021年研究区NDBI与LST协调性的3种类型分布较为均衡，协调型区域变化比较稳定，磨合型和拮抗型区域分别呈先增加后减少和先减少后增加趋势。2006—2011年，城镇发展初期，LST变化超前于城镇发展，呈协调超前和磨合超前型；2011—2017年，冷岛强度减弱，中心城市LST协调性有所上升，磨合型范围扩大；2017—2021年，城镇发展后期，城镇扩散与LST发展具有局限性，LST变化滞后于城镇发展，拮抗型面积有所上升。