城市热岛效应是指城市中的温度明显高于周边环境温度的现象^[1]，已成为现代城市气候的主要特征之一。影响热岛效应的因素一般有人为热量的排放、下垫面性质以及城市结构等^[2]。随着城市化进程的加快，大面积的自然地表被热容量更小、吸热更快的人造表面所替代，导致城市热岛现象的强度和影响范围越来越大，对城市的生态环境和居民生活质量造成了不良的影响。因此，研究热岛问题对城市的科学规划、生态城市的建设有着重要的参考意义。常用的研究方法有地面气象站观测^[3]、空间模拟^[4]和热红外遥感技术。相对于前两种依靠气象数据的技术，卫星热红外遥感有着能获取大范围的地表温度信息的优势，因此成为分析城市热岛的空间分布的主要手段^[5]。

改革开放以来，广州经历了快速的城市扩张^[6]，城市热岛问题日益凸显。已有学者利用Landsat卫星数据对广州区域进行热岛研究，但多数研究重点关注于广州热岛效应与土地利用类型之间的关联性^[7-8]、单一时相的热岛空间格局^[9]，而多时相、长时间序列的遥感分析还比较少见。因此，本文基于1995年~2019年的Landsat卫星数据，利用温度反演算法对广州地区24年来的热岛空间分布变化及其原因进行分析，为减缓热岛效应，合理规划城市发展和生态环境建设提供参考。

1 研究区及数据源概况 1.1 研究区概况

广州市是广东省省会，由荔湾、越秀、海珠、天河、白云、黄埔、番禺、花都、南沙、从化、增城11个行政区组成，总面积为7 434.4 km²。该城市具有典型的亚热带季风气候，年平均气温21.8 ℃~22.8 ℃，最高月平均气温发生在7、8月份，为28 ℃~30.3 ℃^[10]。改革开放以来，由于经济的快速增长，人口及城市规模迅速扩大。至2019年末，广州市常住人口已达1530.59万人，城镇化率达86.46%。本文将研究区划分成6个行政单元，研究区域如图1所示。

1.2 数据来源

本次研究使用的卫星影像为美国地质调查局(United States Geological Survey, USGS)提供的6期Landsat-8 OLI/TIRS 及9期Landsat-5 TM影像，每期均包含行列号为122/044和122/043的两景影像，其热红外波段空间分辨率均重采样至30 m^[11]。选择的影像为L1TP级，经过了地形校正，均为云量较少的影像，具体日期信息如表1所示。下载的影像再经过辐射校正^[12]、大气校正^[13]等预处理操作，消除大气造成的辐射误差。最后通过裁剪和镶嵌操作生成研究区的图像。

另外，本次研究使用的气象数据来自美国国家环境预报中心(National Centers for Environmental Prediction, NCEP)提供的全球地面观测数据以及欧洲中期天气预报中心(European Centre for Medium-Range Weather Forecasts, ECMWF)的全球陆地表面观测反馈记录报告，广州国家基本气象站编号为59287，观测的时间间隔为3 h。

2 基本原理与方法 2.1 土地利用分类

本文采用支持向量机(Support vector machines, SVM)的监督分类算法^[14]来获取土地利用分类图，土地利用类型包括：耕地、草地、林地、建筑用地、水体、未利用地。采用分层随机抽样的方法，来抽取验证样本进行精度评定。最终所有分类图的总体分类精度均在85%以上。

2.2 地表温度反演

根据不同的数据源，目前地表温度的反演算法主要包括大气校正法、单窗算法^[15]、劈窗算法^[16]、单通道法^[17]等。本研究中Landsat-8地表温度反演采用大气校正法^[18]，由于2000年前缺少大气参数，Landsat-5采用广泛应用的单窗算法。

2.2.1 大气校正法

大气校正法又称辐射传导方程法，该方法主要使用大气模型来模拟大气对地表热辐射的影响，然后把这部分大气的影响从卫星传感器观测到的热辐射总量中减去来求解地表温度。该方法表达式为^[19]

$ {I}_{\rm{sensor}}=[\varepsilon B({T}_{\rm{s}})+(1-\varepsilon ){I}^{\downarrow }]\tau +{I}^{\uparrow } $

(1)

式中， $ {I}_{\rm{sensor}} $ 为卫星高度上传感器测得的热辐射强度 ( ${\rm{W}} \cdot {\rm{m}}^{-2} \cdot {\rm{sr}}^{-1} \cdot {\text{μ}}{\rm{m}}^{-1}$ )； $ \varepsilon $ 为地表比辐射率； $ B\left({T}_{{\rm{s}}}\right) $ 是由Plank定律推导得到的黑体热辐射强度，其中 $ {T}_{{\rm{s}}} $ 为地表温度 (K)； $ \tau $ 为大气透射率， $ {I}^{\uparrow } $ 、 $ {I}^{\downarrow } $ 分别为大气的上行和下行热辐射强度( ${\rm{W}} \cdot {\rm{m}}^{-2} \cdot {\rm{sr}}^{-1} \cdot {\text{μ}}{\rm{m}}^{-1}$ )。 $ \tau $ 、 $ {I}^{\uparrow } $ 、 $ {I}^{\downarrow } $ 可以通过NASA的大气校正参数计算器( https://atmcorr.gsfc.nasa.gov)获取。

在求得地表比辐射率(见式(7)~式(10))后，便能得到 $ B\left({T}_{{\rm{s}}}\right) $ ，再通过式(2)求解地表温度^[20]

$ {T}_{{\rm{s}}}=\frac{{K}_{2}}{{\rm{In}}\left[1+\dfrac{{K}_{1}}{B\left({T}_{\rm{s}}\right)}\right]} $

(2)

式中，对于Landsat-8波段10， ${K}_{1}\!\!\!=\!\!\!774.89 ({\rm{W}} \!\!\cdot\!\!{\rm{m}}^{-2} \!\!\cdot\! \! {{\rm{s}}{\rm{r}}}^{-1} \!\!\cdot\!\! {{\text{μ}} }{\rm{m}}^{-1}$ )， ${K}_{2}=1\;321.08$ K。

2.2.2 单窗算法

根据覃志豪提出的单窗算法，地表真实温度反演公式为^[21]

$ \begin{split} {T}_{\rm{s}}= & \{{a}_{6}(1-{C}_{6}-{D}_{6})+ [{b}_{6}(1-{C}_{6}-{D}_{6})+\\&{C}_{6}+{D}_{6}]{T}_{6}-{D}_{6}{T}_{\rm{a}}\}/{C}_{6} \end{split} $

(3)

$ {C}_{6}={\tau }_{6}{\varepsilon }_{6} $

(4)

$ {D}_{6}=\left(1-{\tau }_{6}\right)\left[1+{\tau }_{6}\left(1-{\varepsilon }_{6}\right)\right] $

(5)

式中， ${a}_{6}=-60.326\;3$ ， ${b}_{6}=0.434\;36$ ， $ {T}_{\rm{s}} $ 为地表温度(K)， $ {T}_{6} $ 为卫星高度传感器探测的像元亮度温度(K)； $ {T}_{\rm{a}} $ 为大气平均作用温度， $ {\tau }_{6} $ 为大气透射率， $ {\varepsilon }_{6} $ 为地表比辐射率。 $ {T}_{\rm{a}} $ 的估算，根据热带平均大气经验公式计算

$ {T}_{\rm{a}}=17.976\;9+0.917\;15{T}_{0} $

(6)

式(6)中， $ {T}_{0} $ 为气象站的实测气温数据(K)。当大气水分含量在0.4~3.0 $ {\rm{g}}/{{\rm{c}}{\rm{m}}}^{2} $ 变动区间时， $ {\tau }_{6} $ 的估计方程如表2所示^[21]。

根据已有的研究表明，大气水分含量同地面露点温度存在相关性。地面露点温度可由气象站实测数据获得，再根据经验公式估计出大气水分含量^[22]。根据土地分类图可分成自然表面、水域和城镇表面，用以下公式计算地表比辐射率^[23]

$ {\varepsilon }_{\rm{water}}=0.995 $

(7)

$ {\varepsilon _{{\rm{nature}}}} = 0.962\;5 + 0.061\;4{P_{\rm{v}}} - 0.046\;1{P_{\rm{v}}}^2 $

(8)

$ {\varepsilon _{{\rm{built}} - {\rm{up}}}} = 0.958\;9 + 0.086{P_{\rm{v}}} - 0.067\;1{P_{\rm{v}}}^2 $

(9)

式中， $ {\varepsilon }_{\rm{water}} $ 、 $ {\varepsilon _{{\rm{nature}}}} $ 、 ${\varepsilon _{{\rm{built-up}}}}$ 分别代表水域、自然表面、城镇表面的比辐射率。通过公式(10)求植被覆盖度 $ {P}_{\rm{v}} $

$ {P}_{\rm{v}}=\frac{{\rm{NDVI}}-{{\rm{NDVI}}}_{\rm{s}}}{{{\rm{NDVI}}}_{\rm{v}}-{{\rm{NDVI}}}_{\rm{s}}} $

(10)

式中， ${{\rm{NDVI}}}_{\rm{v}}$ 、 ${{\rm{NDVI}}}_{\rm{s}}$ 分别为植被和裸土的归一化植被指数值，取经验值 ${{\rm{NDVI}}}_{\rm{v}}=0.7$ 和 ${{\rm{NDVI}}}_{\rm{s}}=0.05$ 。当某个像元的NDVI大于0.7时， ${P}_{\rm{v}}$ 取值为1；某个像元的NDVI小于0.05时， $ {P}_{\rm{v}} $ 取值为0。

2.3 热岛等级划分

不同时相影像的天气条件都存在差异，反演的温度值不能直接进行比较，但温度的空间分布是相对稳定的状态。采用划分强度等级的方法，既能将不同时相的结果统一到一个标准，也能更直观描述研究区地表温度差异情况。本文采用均值−标准差的方法将地表温度划分成5个等级，如表3所示。其中，低温区和较低温区统称冷岛区，较高温区和高温区统称热岛区。

3 结果 3.1 温度反演结果验证

本文所用Landsat影像获取时间均在GMT(Greenwich Mean Time)1:50 ~3:00之间，而气象数据观测时间为GMT3：00。由于时间相近，可近似将气象站所在位置的地表温度反演值与气象站观测值比较，用以检验温度反演结果,如图2所示。结果表明，温度反演值均比气象观测值大，这是由于广州在卫星过境时间处于地表升温阶段，且气温值是在距地1.5 m高处观测所得，地表温度比气温升温更快，温度更高。检验差值最大有3.8 ℃，大部分差值小于2 ℃，均方根误差为1.88，地气温差属于正常范围，温度反演结果可靠。

3.2 热岛动态变化分析

图3显示了1995年~2019年间共15期温度反演图的热岛和冷岛面积的变化，结果表明，广州热岛区域的面积在1995年~2001年间变化不大，呈相对稳定的状态；2001年~2014年进入快速增长阶段，整体呈增长趋势，由2001年的722 km²扩大至2014年的1176 km²，年均增长率为3.82%；自2014年起，热岛面积呈缓慢减少的趋势。而冷岛区域面积在2008年前处于波动状态，2008年后整体呈增长趋势，由2008年的733 km²增长到2019年的1139 km²。

图4显示了其中6年的地表温度等级分布图，可看出多年以来温度区的空间分布变化。1995年，城市主要热岛区域集中在珠江沿岸的港口区域以及海珠区西部、新塘镇；北部的从化、花都的热岛出现在裸地和因季节性变化而裸露的耕地(见图4(a))。1997年，沿珠江形成了荔湾−海珠−天河−黄埔的强热岛带。除此以外，各个市辖区的中心镇也开始形成较高温区为主的热岛中心：花都的新华街道、从化的街口街道、增城的荔城街道、番禺的市桥街道(见图4(b))。1997年到2004年间，主城区中心区域热岛效应加剧，南边番禺、南沙热岛面积增加，北边的花都新华镇和白云机场形成两大热岛中心(见图4(c))。2004年到2014年间，热岛区以沿珠江强热岛带为基础，往北沿主要交通道路呈辐射状蔓延至花都，往东向新塘镇、永和经济技术开发区扩张，往南向番禺和南沙扩散(见图4(d)和图(e))。2014年~2019年间，热岛区域无明显的大范围扩张(见图4(f))。

图5为广州市各区热岛占比统计图。可以看出，番禺、花都、主城区等区热岛面积都先呈上升的趋势，并于2014年后有所下降；从化和增城由于其地形和高植被覆盖率的生态环境，热岛占比维持在较低的水平；南沙则是持续增加的趋势。各区之间比较，主城区热岛面积最大，热岛占比最高可达37.9%，是热岛效应最严重的区域；番禺区仅次于主城区，在2014年可达31.2%；1997年~2014年间，主城区热岛面积增长速度最快，为12.7 $ {{\rm{k}}{\rm{m}}}^{2} $ /a，其次是花都区的10.3 $ {\;{\rm{k}}{\rm{m}}}^{2} $ /a。

3.3 土地利用类型与温度区分布的关系

通过计算不同温度水平下不同土地利用类型的面积比，可分析出不同土地利用类型对广州热岛效应的影响。表4计算了所有年份的土地利用类型的3种温度区的平均值，可以看出，热岛区占比最高的土地类型是建筑用地，有79.2%，其次是未利用地的8.29%。随着温度水平的降低，建筑用地所占的比例有所下降，林地的比例则开始增加。冷岛区占比最高是林地，高达94.18%，在比例上有着绝对的优势，其次是水域的4.69%。

结合谷歌地球和OpenStreetMap地图进行叠加分析，可以发现热岛区域在建筑用地中的分布特征：容易依附于主要的交通网络分布，如许广高速和广园快速路；较高温区广泛分布于人群密集的住宅区和商业圈，如员村和石牌等地；高温区基本以斑块状分布在工业区、产业聚集地和经济开发区，代表性区域为广州经济开发区、狮岭镇等。而冷岛区域则大面积聚集在森林茂密的广州北部，特别是拥有大量水库、湖泊和森林的从化。

4 讨论

从20世纪90年代初开始，广州进入了快速的城市化建设时期，围绕着中心城区进行向心聚集型的城市发展。2000年后，《广州市城市总体规划(2001~2010)》中提出城市空间发展的基本策略为：南拓、北优、东进、西联，确定南部、东部为城市发展的主要方向。城市发展模式转变为多中心式、组团式，并将旧城区的传统产业向黄埔至新塘一线集中迁移，形成密集的产业发展带。从图6显示的广州城市的扩张可知，除了中心城区的建筑密度增大外，南北方向和东方向的建筑密度均有明显的增加；外围离散的建筑均开始往中心城镇聚集，并沿交通网络向外延伸。这样的城市扩张的模式与政策建设的方向是一致的。与之对应的，图4显示的热岛空间分布格局呈现着相似的扩散模式。1995~2001年间，广州虽然处于快速的城市化建设进程中，但除了珠江沿岸的中心城区和番禺，其他地区的人口、建筑密度以及工业发展规模依然相对较小，因而还没有大范围的热岛面积增长，在图3中表现为变化不大的热岛曲线。

2001年以后，结合图5可以看出，花都、番禺热岛扩散速度最快，主城区也随着白云区的发展，热岛面积也有很大的增长。沿珠江强热岛带工业区形成的高温区有所减少，永和镇和新塘镇的热岛效应愈发强烈，各区中心城镇也因人口和产业的聚集形成局部热岛并通过交通网络向外扩散。

为分析热岛效应减弱原因，将1995年的城市建筑用地(旧城区)与2010年后新建城区的热岛面积做比较，如表5所示。可以看出2014年后，在建筑用地面积持续增长的情况下，新建城区的热岛面积仅有微弱的增长，而旧城区的热岛面积却大幅下降。这表明，城市建成区的热岛整体有得到减缓，且旧城区减少幅度明显。

在工业污染整治方面，此阶段广州市政府依旧重视“中调”的战略规划，继续推行“退二进三”的策略，将工业区外迁或将旧厂房改造^[24]。据《广州统计年鉴》显示，规模以上工业煤炭消费量从2007年的2.4千万吨减至2018年的1.3千万吨。广州多年以来不断压缩煤炭的使用量和进行清洁能源改造，使得越秀、海珠、荔湾、天河在2014年实现“无燃煤区”，减少了二氧化碳等温室气体的排放。建筑面积接近饱和的老城区，在减少工业区热源和温室气体减排的情况下，热岛效应能得到减缓。同时，工业的集聚化将城区分散的热源得到集中的控制，提高工业热源和能源的利用率，减少热量散失和释放，这是近年来新建城区热岛面积没有大幅度增长的原因之一。

另外，2010年颁布的《广州市城市绿地系统规划(2001~2020年)》和2012年颁布的《广州市城市总体规划(2011~2010)》中均强调了要改善城市生态系统空间结构，保护森林和建立绿道等。根据广州统计年鉴所示，建成区绿化覆盖率由2010年的40%增加到2018年的45%。而城市生态环境质量的改善也能有效减弱城市热岛效应，这也可能是造成2014年后热岛面积下降和冷岛面积持续扩大的原因。

5 结论

综上所述，根据广州市1995年~2019年的遥感影像反演出的温度分布图进行时空分析可以得出以下结论：

(1) 温度区的分布与土地利用类型相关，热岛区域主要分布在建筑用地，容易依附于主要交通线分布。其中较高温区多分布在人口和建筑密度大的住宅区和商业圈，而高温区主要出现在工业区、产业聚集地和经济开发区。冷岛多出现在水库、湖泊和森林区域。

(2) 广州20世纪90年代的向心聚集发展模式导致荔湾−海珠−天河−黄埔的沿珠江强热岛带的形成；2000年以后热岛向南北方向和东方向扩张，并随着市辖区的中心城镇的发展，形成以荔湾−海珠−天河−黄埔−新塘热岛带、白云向北的辐射状热岛带为主导，多个城镇热岛中心围绕的空间分布格局(见图4(f))。热岛变化的趋势与城市发展模式基本一致。

(3) 热岛面积在2000年后整体呈增长趋势。在2014年到达顶峰，之后以缓慢速度减少，热岛效应开始减弱。热岛减缓现象出现在旧城区，其变化原因可能与工业区外迁、燃煤消耗减少、绿化改善等原因有关；冷岛区域在2008年前无明显规律，之后呈增长趋势。

(4) 主城区、花都和番禺的热岛效应最严重，从化区最弱；1997年~2014年间，主城区热岛面积增长速度最快，其次是花都；南沙虽然热岛面积小，但依旧处于持续增长阶段。

致谢：感谢美国地质调查局(USGS) 免费提供Landsat-5/8卫星数据的使用。感谢美国国家环境预报中心(NCEP)提供的全球地面观测数据以及欧洲中期天气预报中心(ECMWF)的全球陆地表面观测反馈记录报告。