城市热环境风险预测及时空格局分析

乔 治, 孙宗耀,孙希华,徐新良,杨 俊

1 天津大学环境科学与工程学院,天津 300072 2 山东师范大学,地理与环境学院,济南 250014 3 中国科学院地理科学与资源研究所,资源与环境信息系统国家重点实验室,北京 100101 4 辽宁师范大学人居环境研究中心,大连 116029

城市热岛作为城市化进程的产物,是城区温度明显高于郊区的现象[1]。城市发展不断改变下垫面形态特征和土地利用/覆被,造成潜热通量减小和显热通量增加,城市热岛效应日益严重[2]。同时,为消除城市热岛效应产生负面效应所付出的经济代价也是极为巨大的。现阶段,热环境安全作为城市生态安全的重要组成部分,成为影响城市生态环境和可持续发展的重大风险问题[3]。保障城市热环境安全的前提是通过热环境时空格局和过程模拟与预测实现热环境风险的识别、评估和预控。但是,当前城市热环境风险研究没有引起足够重视,热环境风险分级标准的缺失也使得城市热环境安全防范和调控措施相对滞后[4]。

学者先后通过数理统计模型、能量平衡模型、数值模型、解析模型和物理模型等理论和技术方法模拟与预测城市热环境时空格局和过程[5]。数理统计模型通常利用统计规律建立城市热环境时间序列趋势方程,但并未真正揭示热环境空间格局变化过程[6- 9]；能量平衡模型更加关注城市与周围地区进行的热量传递、交换关系,对于城市热环境在微观尺度上的形成机理研究较为深入,但因所需参数过多导致宏观尺度热环境时空变化过程分析难度较大[10-11]；数值模型是以热力学和动力学为基础,模拟和预测不同时空尺度下城市热环境格局和过程。尽管数值模型弥补了统计模型和能量平衡模型的不足,但下垫面性质(如城市人为热源、城市建筑物)及其所处的地形特征会使模拟结果异常敏感[12- 15]；解析模型通过重新建立热环境模型参数、初始条件和其他输入信息以及模拟时间和结果之间的关系,可以直观揭示城市热环境及其影响因素的波动性质,但模型所设置的多种理想条件难免会与现实条件有一定偏离[16]；物理模型中最常见的风洞试验,将研究区场景按一定比例微缩,模拟气体和热量交换过程,更加真实地反映城市热岛变化以及不同参数变化对热岛强度的影响,但其研究成本高昂。当前城市热环境模拟预测的诸多方法无法平衡空间精度和模型复杂度的要求,简单的数理统计方法无法实现城市热环境空间格局的预测,而较复杂的能量平衡模型、解析模型和物理模型受制于硬件环境和人为经验等难以广泛应用。因此,如何实现城市热环境风险的快速精准预测是未来研究的重点。

当前,GIS和遥感技术的发展则有效推动了城市热环境时空格局和过程模拟与预测研究,在提高空间精度和简化模拟过程中具有显著意义。尤其通过土地利用/覆被及时空格局变化揭示城市地表能量吸收和传输机制,成为当前研究城市热环境时空异质性和格局过程预测的重要研究方法和手段。本文拟构建基于MARKOV-CA(Cellular Automata)的城市热环境风险预测模型,以北京市为研究区,根据土地利用和城市热环境等级的空间匹配性,计算热环境等级与土地利用及其空间变化的概率关系,并将其作为CA模型迭代因子层。从空间和时间两个维度上,建立基于土地利用类型的城市热环境时空演化预测方法,大大降低了热环境模拟的复杂过程,实现城市热环境快速精准预测。在此基础上,构建基于现状风险与预测风险的城市热环境风险划分标准,预测城市热环境风险时空格局,为城市总体规划和城市生态安全保护提供决策支持,最终实现城市热环境风险防范和调控。

1 研究区域与数据来源

1.1 研究区概况

北京市地处华北平原北部,面积1.68×105km2,地势西北高、东南低,属温带大陆性季风气候,夏季炎热多雨,冬季寒冷干燥,年平均气温为10—12℃,年平均降水量为600—700 mm。作为我国政治、文化和国际交往中心,近50年来,北京市基本建设投资增加近390倍,城市人口增加6.5倍,线性增温速率达到0.36℃/10a[17]。2012年北京市热岛足迹面积比2001年增加1.4倍[18],严重干扰其生态坏境可持续发展[19- 21]。北京市单核心圈层式发展模式极易影响其城市热环境,而其城市空间发展模式又代表着我国大部分城市发展共性。研究表明,当建设用地面积比例超过50%时,该区域会产生显著的热岛效应[18]。因此,从土地利用角度预测北京市城市热环境风险,可为我国城市生态安全防控提供理论和技术支撑。

1.2 数据来源

由于MODIS具有高时间分辨率优势,可保证区域地表温度年际和季节变化规律对比分析的时间一致性。本文使用2005年和2015年夏季(6—8月)MODIS地表温度8天合成产品(MOD11A2)表征城市热环境时空格局和过程,每年包括13期数据,覆盖研究区需要4景影像进行镶嵌处理。MODIS地表温度产品使用分裂窗算法,由MODIS的第31通道(10.780—11.280 μm)和第32通道(11.770—12.270 μm)地表比辐射率和亮温作为输入条件反演得到,空间分辨率为1 km。大量研究表明,MODIS分裂窗算法反演得到的地表温度达到了1K的精度[22]。原始MODIS地表温度产品为HDF格式,需使用MODIS Reprojection Tools软件进行几何纠正和重采样等预处理,数据存储为GEOTIFF格式,输出投影Geographic,基准面是WGS84。对预处理后的数据,需要进一步将其像元灰度值转为地表实际温度值。在计算季节地表平均温度时需进行云掩膜处理,对被云覆盖的像元不进行计算,从而去除温度的离群值,最终得到北京市2005和2015年夏季平均地表温度空间分布数据。

图1 北京市2015年土地利用现状图Fig.1 Land use type in Beijing in 2015

本文土地利用数据使用中国科学院地理科学与资源研究所生产的土地利用现状遥感监测数据产品(2005年、2010年和2015年),空间分辨率为1 km。土地利用类型划分为耕地、林地、草地、水域、城乡建设用地、未利用地等6种一级土地利用类型和25种二级土地利用类型,其中城乡建设用地类型包括城市建设用地、农村居民点和工矿用地等3种二级土地利用类型。通过验证分析,土地利用类型综合评价精度可达92.9%[23],该数据集已广泛应用于城市热环境时空格局分析和生态环境综合评估等领域(图1)。

2 研究方法

首先根据2005年和2010年城市热环境空间格局建立城市热环境分级空间转移矩阵,并进一步与2005年土地利用类型空间格局耦合,构建热环境分级变化与土地利用类型空间关系矩阵；其次,将所构建的城市热环境分级空间转移矩阵、热环境分级变化与土地利用类型空间关系矩阵分别作为MARKOV和CA规则,预测2015年热环境空间格局并通过2015年热环境空间格局现状进行验证,最终预测2020年北京城市热环境风险空间格局(图2)。

图2 总体技术路线图Fig.2 General technical routeMARKOV,马尔科夫；CA,元胞自动机 Cellular automata

2.1 地表温度等级划分

根据已有研究结果,北京市夏季城市热岛范围和强度最为显著[18]。因此,本文主要对北京市夏季城市热岛时空格局及热环境风险进行预测。为消除温度变化年际影响,首先对2005—2015年3期地表温度归一化处理并采用平均值标准差法进行定量分级[24](表1)。

(1)

式中,ti为归一化后的像元值；Ti为原始地表温度值；Tmax为研究区地表温度最大值；Tmin为研究区地表温度最小值。

表1 地表温度等级划分

tmean为所有像元平均值,tstd为所有像元标准差;LST,地表温度 land surface temperature

2.2 城市热环境空间格局预测

城市热环境空间格局预测需要明确不同热环境等级的像元数量及其空间位置。不同热环境等级的像元数量基于MARKOV过程模型定量计算[5]；同时,根据两期热环境等级变化及其与土地利用类型的空间关系建立CA空间过程规则[9, 25],从而预测目标年份热环境等级空间分布。

在MARKOV过程模型中,热环境不同等级可以看作MARKOV过程中的可能状态。通过计算前后两期不同热环境等级间转换的面积和概率,建立热环境等级间转移的数量规则,最终预测目标年份不同热环境等级的像元数量。

(2)

C=[C低温,C次低温,C中温,C次高温,C高温]

(3)

式中,C代表像元数量。

Ct+a=P×Ct

(4)

式中,t代表过程年份,a代表过程年份与初始年份的时间间隔。

目标年份城市热环境等级空间分布是基于CA空间模型进行预测。根据MARKOV预测得到的不同城市热环境等级像元数量,CA空间模型建立了不同热环境等级像元的空间转移规则。这种转移规则是基于前后两期城市热环境等级空间转移及其与初始年份土地利用类型之间的空间相关性所构建的概率矩阵,可得到不同土地利用类型的热力等级转移概率(表2)。

表2 土地利用类型热力等级转移概率

P高(耕|高)表示在所有初始年份和目标年份均为高温等级的像元,耕地类型所占的概率。其他依此类推

2.3 热环境风险等级划分标准

基于前后两期热环境等级建立城市热环境风险等级标准。风险等级划分为极高风险、高风险、临界风险、较安全和理想安全5级,不同风险等级具有不同的指示意义(表3)。

3 结果与分析

3.1 基于MARKOV过程的城市热环境分级时空规则

通过分析2005—2010年北京市热力等级时空格局特征(图3),2005年高温区域主要分布于东城区、西城区以及朝阳区西南部、大兴区西部、丰台区中部、石景山区南部和顺义区首都国际机场区域,高温斑块分布较为分散；朝阳区、海淀区、石景山区、丰台区、大兴区以及顺义区几乎均处于次高温区,通州北部、昌平区中东部、房山区东部和延庆区西部也有大片次高温区热力斑块分布；低温区域主要分布于怀柔区中部和密云区密云水库周边。2010年,北京市建成区城市热岛效应进一步加剧,丰台区东部和大兴区北部出现新的高温区域,并且原有的高温区域进一步扩大。

表3 城市热环境风险等级划分标准

初始年份热环境等级与目标年份热环境等级采用并集规则；同时满足两个风险等级的像元依据较高风险等级标准

分析2005—2010年北京市热环境分级时空转移特征 (表4),低温和中温热力等级像元数量减少,次低温、次高温和高温热环境等级像元数量增加,整体呈现升级趋势,其中高温等级增幅为85.81%,主要来源于次高温等级；次高温等级中有747个像元转向中温等级,551个像元转向高温等级,分别占次高温等级像元的12.16%和8.97%；中温等级有26.38%转向次高温等级,有58.19%维持中温等级；次低温等级中仅有161个和414个像元转变为低温等级和中温等级,分别占4.18%和10.75%；低温等级面积减幅为25.18%,其中有454个像元转向次低温等级。

表4 北京市2005—2010年城市热环境分级转移矩阵

3.2 基于CA过程的土地利用类型与城市热环境分级变化空间规则

耦合土地利用类型和城市热环境等级空间变化关系,分析北京市2005—2010年城市热环境分级变化与土地利用类型相关关系(表5)。林地在保持低温情况下比例为88.66%,在低温区域转向次低温区域比例高达94.27%；耕地主要分布于中温条件和次高温条件,在中温转向次高温区域和保持次高温比例都超过50%；城市建设用地主要分布在高温区域,其中保持高温情况下,城市建设用地比例为69.73%。根据土地利用类型与城市热环境分级变化空间规律,建立2015年北京市城市热环境分级空间分布概率(图4)。

表5 2005—2010年北京市城市热环境分级变化与土地利用类型空间关系

图4 2015年北京市城市热环境分级空间分布概率Fig.4 The probability of spatial distribution of urban thermal environment levels in Beijing in 2015

图5 2015年北京市城市热环境分级空间分布现状(a)和预测图(b)Fig.5 The status and prediction map of urban thermal environment in Beijing in 2015

3.3 城市热环境空间格局预测及验证

首先,通过2005—2010年城市热环境分级时空变化MARKOV过程计算2015年北京市热环境各等级像元数量；在此基础上,耦合土地利用类型与城市热环境分级变化CA过程所生成的2015年北京市热环境分级空间分布概率,预测2015年北京市热环境分级空间分布(图5)。北京市高温等级区域集中在以东西城区为中心的城市中心区域并有向南延伸的趋势,其中朝阳区中西南大部分地区、海淀区东南部、石景山区中南部、丰台区中东部、大兴区西北部及顺义区西南部热岛效应进一步加剧；平谷区、昌平区、顺义区次高温等级区域继续扩大；门头沟区西北部、延庆区西部、怀柔区中部和密云水库地区仍保持低温等级。对比北京市2015年热环境分级空间现状,热环境分级空间分布预测与真实结果具有较高一致性(Kappa系数为0.73)。

3.4 城市热环境风险时空格局变化特征

根据城市热环境风险等级划分标准,分析2015—2020年北京市城市热环境风险时空格局变化特征(图6和表6)。

2015年北京市热环境高风险区比例最高,达到36.16%,主要分布在通州区、大兴区和房山区东部、丰台区西部、朝阳区东北部、顺义区、平谷区中南部、密云区南部和怀柔区东南部地区；极高风险区比例为9.66%,集中于东西城区、朝阳区、丰台区、石景山区、海淀区东部和大兴区西北部；临界安全区和较安全区占比分别为26.64%和22.51%,主要分布在高风险区外围的房山区、门头沟区、昌平区、怀柔区、延庆区和密云区；热环境风险理想安全区域仅占全市面积的5.02%,分布在怀柔区中部、门头沟区和延庆区西部。

图6 2015—2020年北京市城市热环境风险空间分布图Fig.6 The spatial distribution of urban thermal environment risk in Beijing between 2015 and 2020

2020年北京市热环境极高风险区比例升高到12.08%,新增区域主要分布在朝阳区和海淀区外围区域、顺义区以及昌平区、平谷区、密云区的城镇地区；高风险区降至28.96%,变化区域主要位于平谷区西部、顺义区东北部和怀柔西南部；较安全区域面积比例达到29.19%,主要由密云区和房山区的临界风险区转移；临界风险区未有明显变化；理想安全区面积进一步减少,仅占全市面积的3.14%。

表6 2015—2020年北京市城市热环境风险等级空间转移矩阵

分别选取斑块数量(NP)、形状指数(LSI)和聚合度指数(AI)从数量特征、形态特征和结构特征3个方面分析2015—2020年热环境风险时空变化特征[5, 26](表7)。对于极高风险区,斑块数量有所增加(主要是零星斑块的出现),但整体形状逐渐规则,聚合程度也有所提高；高风险区和临界风险区呈现一致特征,斑块数量明显增加,斑块形状也更为复杂,聚合程度均有所降低；理想安全区和较安全区斑块数量明显减少,整体形状明显变得规则并且聚合程度有所提升。

表7 2015—2020年北京市城市热环境风险等级斑块景观格局指数

NP, 斑块数量 number of patches; LSI, 形状指数 landscape shape index; AI, 聚合度指数aggregation index

4 讨论与结论

本文构建了基于MARKOV-CA的城市热环境风险预测模型,预测2015—2020年城市热环境风险时空格局及其变化特征。所构建的城市热环境风险预测模型以土地利用类型与热力等级的空间匹配性作为CA模型因子进行迭代分析,预测未来热环境空间格局,并显示出较好的热力等级模拟效果,2015年模拟情况与真实情况的Kappa系数为0.73,该模型在一定程度可替代复杂的解析模拟模型。在此基础上,所构建的基于现状风险与预测风险的城市热环境风险划分标准,从时间过程角度实现热环境风险空间评估,具有较强的实践性。本文研究结果如下：2020年北京市热环境极高风险区域将进一步扩大,所占比例达到12.08%,有向东西方向延伸的趋势并且空间聚合程度有所提高；高风险区所占比例较2015年下降至28.96%,但斑块数量明显增加,斑块形状更为复杂,空间破碎化程度有所增加；理想安全区域面积逐步减小,多出现于北京市西部边缘地区。基于上述研究结果进一步提出城市热环境空间调控措施：在空间上,城市建设用地与极高风险区和高风险区具有极高的空间相关性,因此优化多种土地利用类型空间配置,通过科学设置生态用地分割连片的城市建设用地对于防范城市热环境风险具有重要意义；从管理上,严格控制城市建设用地增长数量,划定城市增长边界,对不同热环境风险等级区域采取不同的发展控制手段,在热环境安全较好区域可进行适度的城市用地挖潜,在热环境极高风险区严禁继续高强度开发并发展有效降温措施。

本研究所构建的城市热环境风险预测模型的最大优势是可基于土地利用类型与热力等级的空间关系快速预测未来城市热环境风险空间格局,可为热环境调控政策提供判定依据和决策支持,为城市生态评估提供可靠的评价方法,有利于辅助城市规划中各类限制界线的划定。模型具有较高的模拟和预测精度,但是仍具有较高的提升空间。一方面,所使用的MODIS地表温度产品和土地利用数据空间分辨率均为1 km,这在很大程度上影响预测精度。同时,本研究仅使用土地利用一级分类建立其与热环境的空间耦合规律。但是土地利用类型内部结构在一定程度上也会影响地表能量平衡过程,如林地盖度影响地表粗糙度、土壤湿度影响地表反照率和能量平衡过程、建筑形态及空间布局影响城市通风环境等,并最终影响城市热环境。因此,应进一步提高地表温度产品和土地利用数据的空间分辨率,深入挖掘地表温度与土地利用/覆被空间关系特征规律；另一方面,地形等自然因素会影响城市热环境。尽管本研究所使用的CA所具有的空间模拟优势在一定程度上能够减小空间异质性所带来模拟精度的降低,但是人为热排放等社会经济因素对于热环境的影响应进一步考虑,从而提高热环境风险预测的科学性。