地表温度与土地利用类型间的空间尺度依赖性——以成都为例

李喆，陈圣宾，陈芝阳

成都理工大学生态环境学院，四川 成都 610059

城市热岛效应（urban heat island effect，UHI）即城区气温明显高于周边郊区气温的现象（Schatz et al.，2016），热岛效应的产生不仅会对自然环境产生影响，而且还会严重影响城市居民的生活水平（Huang et al.，2012）。研究发现，不同土地利用类型由于自身热容量的不同，对热岛效应的影响也是不同的。

以往关于地表温度（Land Surface Temperature，LST）与土地利用类型间的相关性研究多数集中在建筑用地、绿地与湿地之间。研究认为，建筑区的环境温度最高，绿地和湿地环境温度较低（Mori et al.，2009；Coseo et al.，2014；Wu et al.，2014）。城市湿地作为城市的重要景观组成，在调节周边小气候方面具有显著效果，多数表现为在冬季和夜间增温，在白天和夏季降温（傅抱璞，1997）。不同城市湿地的降温效果存在差异，这与水体布局（聂冲等，2019）、形状水深（Deilami et al.，2017）、季节（Theeuwes et al.，2013）等多方面因素有关，是描述热岛效应的有效指标之一（Wu et al.，2019）。绿地经由植物的蒸腾作用对周围环境进行降温，植物的郁闭度（Xiao et al.，2018）、绿量（李英汉等，2011）、叶面积指数（Kong et al.，2017）等都会影响绿地的降温效果，不同植被类型组合对降温效应也不同（刘娇妹等，2008）。在遥感应用方面，研究者多采用NDVI指数探究植被与地表温度的关系，如Gallo et al.（1993）利用遥感影像对某地区的地表温度和归一化植被指数的关系进行了研究，发现两者呈负相关关系，之后这一结论也被广大学者认同（Rani et al.，2018）。NDBI作为归一化建筑指数，其值不随季节变化发生改变，且与地表温度存在明显的线性关系，是研究热岛效应的有效指标（Ogashawara et al.，2012）。Shahfahad et al.（2020）以印度四大城市为研究对象，发现地表温度在建筑密集区均偏高；Zhao et al.（2021）通过研究西安市气温异质性，发现NDBI是影响建成区与郊区温度差异的最重要驱动因子；Koko et al.（2021）通过研究尼日利亚阿布贾市在 1990—2019年间地表温度与NDBI的关系，发现NDBI与地表温度在不同时期均呈正相关关系。以上结果说明，若想缓解城市热岛效应，需要增加湿地与绿地面积，减少建筑用地面积，然而在有限的城市空间内大力发展湿地与绿地的空间面积是不现实的，如何经过合理布局发挥湿地与绿地的降温功效成为了新的研究热点。韩贵锋等（2011）以重庆市为例，研究了地表温度与NDVI的空间尺度特征，发现两者在空间上有较强的异质性；江颖慧等（2018）通过研究地表温度与NDVI的空间关联性的尺度效应，发现两者的相关尺度在300 m内均显著，并且呈现出明显的季节差异；李斌等（2017）对比了NDVI与NDWI在不同尺度下与LST的空间相关程度，发现随着尺度的增加，NDWI与LST的相关性增强，且相对NDVI，NDWI更适合对温度进行量化分析。

综上所述，土地利用类型对城市热环境具有不同的影响，有的可以缓解城市热环境，有的则会加剧城市热环境。在土地利用类型与地表温度的研究中，用遥感指数描述土地利用类型也是有效的量化指标。目前关于城市热岛与土地利用类型的空间尺度研究多是以植被为主，忽略了水体与建筑用地两大土地利用类型。因此，通过GIS手段提取研究区地表温度与NDBI、NDVI、NDWI值，运用半变异函数识别地表温度与3种遥感指数的空间相关性分析尺度，并结合空间相关性分析手段探究地表温度与各土地利用类型间空间尺度依赖度，在尺度方面认识和把握地表温度与土地利用类型简单空间关系及变化规律，有助于未来城市合理规划蓝绿空间，对有效推动城市生态文明建设提供科学参考。

1 研究区域与数据源

成都市地处亚热带季风气候区，热量充足，雨量丰富，四季分明，雨热同期。除西北边缘部分山地以外，成都市大部分地区表现出的气候特点是：夏无酷暑，冬少冰雪，气候温和，夏长冬短，无霜期长，风速小，湿度大，云雾多，日照少。年平均温度为 15.7—17.7 ℃，年总降水量为 798.3—1541.0 mm，年平均日照时数为685.5—1002.9 h。

以成都市第二绕城高速以内的区域作为研究区（图1），成都市近几年的发展极为迅猛，城市化进程不断推进，成都都市圈的概念提出后，政府更会加大对该地区的发展建设，随之会产生更严重的热岛效应。因此，探究成都热环境的空间分异与成因，对改善成都热环境，提升居民生活水平，建设生态宜居现代化新天府、提高城市可持续发展的综合承载力具有重要意义。

采用的卫星遥感数据均来源美国 NASA的Landsat TM、Landsat ETM+和Landsat OLI-TIRS遥感影像，从美国地质勘探局网站（http://www.usgs.gov/）上获取。条代号为129/39，分辨率为30 m。在2010—2021年间选出冬夏两季共6期影像，冬季为2010年12月22日、2014年12月19日和2021年2月25日，夏季为2016年9月11日、2018年6月5日和2019年8月11日。由于研究区晴朗天气极少，因此不能在时间尺度上选择各年同期数据。

图1 研究区概况Figure 1 Overview of the study area

2 研究方法

2.1 地表温度反演

采用辐射传输方程法反演研究区的地表温度。首先计算卫星传感器所观测到的热辐射总量，然后以此估计出大气对地表热辐射的影响，再从卫星传感器所观测到的热辐射总量中减去这部分大气影响，从而得到地表热辐射强度，最后把这一热辐射强度转化为相应的地表温度（游绚等，2009）。

卫星传感器接收到的热红外辐射亮度值Lλ的表达式如下：

式中：

ε——地表比辐射率；

TS——地表真实温度，K；

B(TS)——黑体热辐射亮度；

τ——大气在热红外波段的透过率；

L↑——大气上行辐射亮度；

L↓——大气下行辐射亮度。

温度为T的黑体在热红外波段的辐射亮度为：

Ts可以用普朗克公式的函数获取：

对于 TM 卫星传感器，K1=607.76 W·m−2·µm−1·sr−1，K2=1260.56 K；对于 ETM+卫星传感器K1=666.09 W·m−2·µm−1·sr−1，K2=1282.71K；对于 TIRS 卫星传感器，K1=774.89 W·m−2·µm−1·sr−1，K2=1321.08 K。

利用中ENVI中的Bandmath工具计算式（2）和式（3）。大气剖面参数可在 NASA提供的网站（http://atmcorr.gsfc.nasa.gov/），输入成相时间和中心经纬度进行获取，查询结果见表1。

表1 影像基本信息Table 1 Basic image information

地表比辐射率通过NDVI阈值法计算：

其中：

Pv——植被覆盖度，通过以下公式计算：

其中：

NDVI——归一化植被指数；

NDVIsoil——完全裸露土壤或非植被覆盖区域的NDVI值；

NDVIveg——纯植被像元的 NDVI值。取经验值NDVIveg=0.70和NDVIsoil=0.05。

2.2 地表覆盖参数计算

归一化植被指数 NDVI的计算（杨春华等，2013）：

式中：

NIR——近红外波段的地表反射率；

Red——红光波段的地表反射率。

归一化水体指数 NDWI的计算（邓开元等，2021）：

式中：

Green——绿光波段的地表反射率。

归一化建筑指数 NDBI的计算（Mori et al.，2009）：

式中：

MIR——中红外波段的地表反射率。

2.3 半变异函数

半变异函数也称为半方差函数，是地统计学中用来描述局域化变量的空间连续变异的一个函数，该函数通过计算两点间空间属性的变异程度与两点间的距离关系来描述要素的分布。在满足二阶平稳或本征假设条件下，半变异函数的计算公式为（Mejia-Dominguez et al.，2012）：

式中：

γ(h)——半变异函数值；

Z(x)——各栅格中的LST或地表参数值；

h——两样本点的空间距离；

N(h) ——空间距离为h的样本点的数量；

Z(xi)、Z(xi+h)(i=1, 2, …,N(h)) ——研究区变量Z(x)在空间位置xi和xi+h上的数值。通过上述公式可以计算出任意给定距离h内样本的半变异值，绘出基于h-γ(h)的半变异方差图（图2）。半变异函数图有3个特征：基台值、变程和块金值，基台值为半变异函数随距离h增加而增加，直到达到平稳状态的值；变程是半变异函数从初始的块金值达到基台值所需要的距离；块金值反映的是最小抽样尺度以下变量的变异性及测量误差。

图2 半变异函数拟合曲线示例Figure 2 An example curve of the semi-variogram function

2.4 空间自相关分析

空间自相关是检验空间上一点的属性值是否与相邻点的属性值存在相关性的一种方法。有研究发现，遥感影像中相邻象元间存在很强的空间自相关性（柏延臣等，2004），但这种空间自相关性将随像元间距的增加而逐渐减弱，当像元间距增至某种程度时，这种空间相性将变得非常弱，以致不再具有相关性，不能提供有意义的遥感信息。Moran’sI作为常用的描述空间自相关性的指标，当 Moran’sI>0表示空间正相关性，其值越大，空间相关性越明显。Moran’sI<0表示空间负相关性，其值越小，空间差异越大，Moran’sI=0，空间关系呈随机性。本研究分别从单变量空间自相关性和双变量空间相关性两个角度出发，LST、NDBI、NDVI、NDWI各自的空间自相关性和 NDBI-LST、NDVI-LST、NDWI-LST的空间相关性。

单变量 Moran’sI的计算公式为（王劲峰，2006）：

式中：

I——单变量 Moran’sI值；

n——研究区内栅格总数；

xi和xj——相邻栅格i和栅格j的属性值，即LST或地表参数值；

——所有样本的平均值；

wij——空间单元和之间的空间连接矩阵。

双变量Moran’sI的计算公式为（Wartenberg，1985）：

式中：

Ixy——双变量 Moran’sI值；

xi——相邻研究区域i的x属性值；

yi——相邻区域j的y属性值；

和——样本中所有x和y的属性平均值；

wij——空间单元和之间的空间连接矩阵。

基于半变异函数确定了 LST与各地表覆盖参数的空间相关性尺度后，为了进一步研究 LST与NDBI、NDVI和NDWI之间空间关系，将研究区栅格数据进行重采样，从原始分辨率30 m的栅格数据转换为分辨率90、150、300、600、900、1200、2400 m的栅格数据，之后提取各栅格数据值在GeoDa软件中建立K=8的K最近邻空间权重矩阵进行空间相关分析，最终得到6期影像在7个分辨率下的LST与各地表参数的单变量Moran’sI值和LST与各地表参数两两之间的双变量Moran’sIxy值。

3 结果与分析

3.1 LST与各地表参数间的相关性尺度判定

通过对6期影像的反演，得到研究区LST与各地表参数分布图，以2018年6月5日的影像为例，从图3a可知，成都市中心城区的地表温度高于郊区的地表温度，说明成都市具有明显的热岛效应，NDBI的高值集中在绕城高速以内的中心城区；NDVI的高值分布在绕城高速以外的西侧、南侧和东侧的城郊区域，低值分布在中心城；NDWI的高值集中在绕城高速以内的区域，呈放射状分布。

图3 研究区LST、NDBI、NDVI、NDWI空间分布Figure 3 Spatial distribution of the LST, NDBI, NDVI and NDWI in the study area

从图4可知，夏冬两个季节的LST和各地表参数的半变异方差图均呈现出先增大后平稳的趋势，并且变程均在300—600 m之间。在300 m尺度之前，基台值随距离增加较快，快速增加表明LST、NDBI、NDVI、NDWI在研究区内存在较强空间自相关性；300—600 m之间，基台值随距离增加的速率减小；600 m尺度后，基台值增幅明显减缓甚至不再增加，说明LST、NDBI、NDVI、NDWI均由空间异质性趋向于空间同质性。LST与各地表参数的半变异方差图表现出明显的季节差异，总体表现为夏季的半变异方差值大于冬季的半变异方差值，其中以LST、NDVI、NDWI表现最为明显。综上，研究区LST、NDBI、NDVI、NDWI的空间自相关分析尺度应小于600 m。

图4 研究区LST、NDBI、NDVI、NDWI的半变异函数图Figure 4 Semi-variograms of LST, NDBI, NDVI and NDWI in the study area

3.2 LST与各地表参数的单变量空间自相关性

从图5可知，研究区的各年份各尺度下LST的Moran’sI均为正值，说明 LST具有很强的空间自相关性，且夏季的LST的空间自相关性要大于冬季的LST的空间自相关性。随着尺度的增大，6条曲线均呈下降趋势，说明随着栅格的合并，栅格数据之间重组，重组后的栅格间的空间自相关性减弱。夏季的3条曲线在在150、300、600 m尺度范围内有3处明显转折，且在不同尺度下，夏季3年的LST空间自相关强度差异表现为2016年<2018年<2019年。冬季的3条曲线在150、300、600、900 m处有明显转折，2010年曲线除在300 m尺度范围内的空间自相关强度略低于2014年与2021年外，其余尺度下，3条曲线基本重合，说明2010—2021年间冬季LST空间自相关强度无明显变化。各年份各尺度下NDBI的Moran’sI值均为正值，并且随着尺度的增大，6条曲线均呈下降趋势，在300 m和900 m处有明显转折，1200 m后呈无规律变化。NDBI为归一化建筑指数，建筑物密度不应随季节发生变化，但如图所示，夏季的NDBI空间自相关强度大于冬季NDBI空间自相关强度，这可能与遥感提取NDBI的计算公式有关。有研究发现，由于夏季植被更为茂盛，裸地面积小，所有遥感影像提取城镇用地更为精确（杨智翔等，2010）。NDVI单变量空间自相关性在冬夏两季随分辨率增大的变化趋势与LST类似，且6条曲线在300 m尺度范围内出现一致的转折点，1200 m后变化无规律。各年份各尺度下的NDWI的Moran’sI值整体上呈下降趋向，6条曲线在300 m尺度范围内有明显转折，此后6条曲线无明显一致转折点，1200 m后6条曲线变化呈无规律状态。

图5 LST、NDBI、NDVI、NDWI的空间自相关性的尺度效应Figure 5 Scale effects of spatial autocorrelation of LST, NDBI, NDVI and NDWI in different years

3.3 LST与NDBI、NDVI、NDWI的双变量空间自相关分析

图6a中，冬夏两个季节的LST与NDBI的双变量Moran’sIxy值均为正值，说明冬夏两个季节的LST与NDBI呈现空间正相关关系，即建筑物密度越大，地表温度越高。有研究证明，建筑物是影像城市热岛效应的重要因素，建筑物的材料改变了地表热交换和大气动力学特征，更易吸收大量的辐射，使城市上空温度升高（周玄德等，2018)，这也解释了为何夏季的 Moran’sIxy值大于冬季的Moran’sIxy值。同样地，LST与 NDBI的双变量Moran’sIxy值在300 m尺度范围内逐渐下降，但在300 m尺度范围外，双变量Moran’sIxy值变化幅度变小，直至无规律。

图6 不同年份下LST与土地利用类型的空间相关性的尺度效应Figure 6 The spatial correlation between LST and land use classification in different years

图6b中，冬夏两个季节的LST与NDVI的双变量Moran’sIxy值均为负值，说明冬夏两个季节的LST与 NDVI呈现空间负相关关系，并且夏季的Moran’sIxy值小于冬季的 Moran’sIxy值，即植被指数越高，周围地表温度越低。有研究证明，植被覆盖对城市热环境具有积极的影响，能缓解城市热环境问题（Maimaitiyiming et al.，2014）。同样地，NDVI与LST的双变量Moran’sIxy值在300 m尺度范围内逐渐上升，但在300 m尺度范围外，双变量Moran’sIxy值变化幅度变小，且冬季的变化幅度要小于夏季的变化幅度，这是因为夏季植被的叶面积指数较大，提高了植被的降温率（刘海轩等，2015）。

图6c中，冬夏两个季节的LST与NDWI的双变量 Moran’sIxy值均为正值，同样地，夏季的Moran’sIxy值大于冬季的 Moran’sIxy值且降幅明显大于冬季，在300 m空间尺度内，夏季的NDWI与LST的空间相关性显著，而冬季Moran’sIxy值随空间尺度的变化不明显。这与以往的大多数研究结果不同，以往的研究结果认为水体由于自身的比热容较大，因此会吸收周围的热量，从而降低周边环境的温度（刘洪杰等，2003）。这与本研究结果不同，通过观察图3d可发现，水体在三环以内分布相对密集，但并没有对热岛效应起到缓解作用，通过查阅资料发现，成都三环内的水体主要以河流的形式存在，主干河道包括府河、南河、沙河、清水河以及江安河水域，水域与陆地的温差主要由水面与陆地的反射率、比热容和水域的水量来决定（傅抱璞，1997）。就夏季而言，水体在正午时分吸收并储存了大量比陆地更多的热量，且成都地区夏季气候湿润，水陆蒸发差异不大，从而使水体对周边环境温度产生了增温效应，这与崔林林在 2018年研究成都热岛效应与下垫面关系时的结果相似（聂冲等，2019）。就冬季而言，水体在冬季正午更容易蓄积太阳辐射热和人类活动热来缓和与周边空气的温差，但这种效果微乎其微（牛少凤等，2008）。

4 结论与展望

以成都市为例，运用遥感和GIS方法与空间相关性软件对地表温度与3地表参数的空间相关性尺度进行了探讨，得出如下结论：

（1）LST、NDBI、NDVI、NDWI的单变量空间自相关性存在明显的尺度效应，300 m尺度范围内，自相关性尺度效应变化明显，300—600 m尺度范围内自相关性尺度效应减弱。夏季的LST单变量空间自相关性最强，NDVI与NDWI整体上表现出夏季自相关性尺度大于冬季。

（2）LST 与NDBI、LST与NDVI、LST 与NDWI的双变量空间相关性存在显著的尺度效应，且表现出明显季节差异，夏季的相关性强于冬季的相关性强度，夏季相关性随尺度变化幅度大于冬季相关性随尺度变化幅度。

（3）LST与NDBI、NDVI、NDWI均具有良好的空间自相关性，因此可以通过合理的绿地、湿地规划和建成区改造缓解城市热岛效应，改善城市热环境。

本研究只选取了成都市作为研究对象，今后可扩充其他地区的典型城市进行探讨，并且可以与当前城市规划模式相结合，进一步分析相关性尺度的机理。本研究从地学统计角度分析了LST和各地表参数间的空间相关性的尺度效应进行定性分析，揭示现象并总结规律，对于LST和各地表参数两两间相互作用的空间机理还有待进一步探究。