环境科学  2022, Vol. 43 Issue (5): 2262-2273   PDF    
基于时空混合效应模型的京津冀PM2.5浓度变化模拟
范丽行1,2, 杨晓辉1,3, 宋春杰1,2, 李梦诗1,2, 段继福1,2, 王卫1,2, 李夫星1,4, 李伟妙1,2     
1. 河北师范大学地理科学学院, 石家庄 050024;
2. 河北省环境演变与生态建设实验室, 石家庄 050024;
3. 河北省科学院地理科学研究所, 河北省地理信息开发应用工程技术研究中心, 石家庄 050011;
4. 河北省环境变化遥感识别技术创新中心, 石家庄 050024
摘要: 为揭示京津冀地区高精度PM2.5的时空分布特征, 以空间分辨率为1 km的MAIAC AOD数据为主要预测因子, 以气象数据、植被指数、夜间灯光数、人口密度和海拔数据作为辅助因子, 构建了一种新的时空混合效应模型(STLME), 在拟合最优次区域划分方案基础上对京津冀地区PM2.5浓度进行预测分析.结果表明, 基于STLME模型的ρ(PM2.5)预测精度高于传统的线性混合效应模型(LME), 其十折交叉验证(CV)R为0.91, 明显高于LME模型的0.87, 说明STLME模型在同时校正PM2.5-AOD关系的时空异质性方面具有优势.最优次区域划分方案识别出PM2.5-AOD关系的空间差异, 并结合缓冲区平滑方法, 提高了STLME模型预测精度.京津冀PM2.5浓度时空变化差异显著, 高值区主要分布在以石家庄、邢台和邯郸为中心的河北南部, 低值区则位于燕山-太行山区; 冬季PM2.5污染最严重, 其次是秋季和春季, 夏季污染最轻.STLME模型提供的高精度PM2.5浓度时空分布预测结果, 为京津冀地区与PM2.5污染相关的健康风险评估提供了科学依据, 也为大气污染源识别提供了科学参考.
关键词: PM2.5      多角度大气校正算法的气溶胶光学厚度(MAIAC AOD)      时空混合效应模型(STLME)      时空差异      京津冀地区(BTH)     
Modeling of PM2.5 Concentrations in the Beijing-Tianjin-Hebei Region Using a Space-time Linear Mixed Effects Model
FAN Li-hang1,2 , YANG Xiao-hui1,3 , SONG Chun-jie1,2 , LI Meng-shi1,2 , DUAN Ji-fu1,2 , WANG Wei1,2 , LI Fu-xing1,4 , LI Wei-miao1,2     
1. College of Geography Science, Hebei Normal University, Shijiazhuang 050024, China;
2. Hebei Laboratory of Environmental Evolution and Ecological Construction, Shijiazhuang 050024, China;
3. Hebei Engineering Research Center for Geographic Information Application, Institute of Geographical Sciences, Hebei Academy of Sciences, Shijiazhuang 050011, China;
4. Hebei Remote Sensing Technology Identification Innovation Center for Environmental Change, Shijiazhuang 050024, China
Abstract: To reveal the spatiotemporal distribution characteristics of high-precision PM2.5 concentrations in the Beijing-Tianjin-Hebei (BTH) region, a space-time linear mixed effects (STLME) model was developed in this study. The MAIAC AOD at a 1 km spatial resolution and the meteorological material, vegetation index, light quantity at night, population density, and altitude data were employed as the main and auxiliary predictive factors in the STLME model, respectively, to estimate the ground-level PM2.5 concentrations on the BTH region by optimizing the sub-regional division scheme. The results indicated that the STLME model with the CV R2 valued at 0.91 outperformed traditional linear-mixed effects (LME) with a CV R2 of 0.87, indicating the superiority of the STLME model in simultaneously correcting the spatiotemporal heterogeneity of the PM2.5-AOD relationship. The optimal sub-region partitioning scheme identified the spatial difference in the PM2.5-AOD relationship and, combined with the buffer smoothing method, improved the prediction accuracy of the STLME model. The PM2.5 levels in the BTH region exhibited strong spatiotemporal variations. The areas with higher PM2.5 concentrations were mainly located in the southern Hebei province centered in the Shijiazhuang, Xingtai, and Handan cities, whereas the Yanshan-Taihangshan mountainous areas were the regions with lower values. In addition, the most heavily polluted season was winter, followed by autumn and spring, and summer was the cleanest season. The spatiotemporal distribution prediction results of high-precision PM2.5 concentrations provided by the STLME model provide a scientific basis for the health risk assessment of PM2.5 pollution in the BTH region and also provide a scientific reference for the identification of air pollution sources.
Key words: PM2.5      MAIAC AOD      space-time linear mixed effects model (STLME)      spatiotemporal difference      Beijing-Tianjin-Hebei region (BTH)     

大气细颗粒物(PM2.5)严重影响区域环境空气质量并危害人类健康[1].流行病学研究表明, 长期暴露于PM2.5污染环境下会增加心血管和呼吸道疾病的发病率及死亡率[2~5].因此, 精确获取时空解析的PM2.5暴露特征对公共卫生研究及区域环境质量评价极为重要[6~9].

虽然地面站点监测可以获取较为准确的PM2.5信息, 但站点稀疏且空间分布不均匀仍难以满足对区域面状空间分布的监测需求[10].相较而言, 空气质量模型可提供覆盖全域的格网数据, 但其对污染物排放清单、污染物理化性质等资料的齐全性要求较高, 致使模型不确定性增大, 影响了其模拟精度[11].基于卫星遥感的气溶胶光学厚度(aerosol optical depth, AOD)为整层大气柱中气溶胶消光系数在垂直方向上的积分, 其与大气细颗粒物具有相近物理机制[12].因此, 学者们利用不同卫星传感器AOD产品在世界各地开展了近地面PM2.5浓度估算研究[13~16].目前, 神经网络[17, 18]、支持向量机[19]和随机森林[20]等人工智能模型可以通过增加隐层节点和层的数量来解释不同变量之间复杂的非线性关系, 明显提高了PM2.5浓度时空分布的模拟精度, 但是过多的隐层节点和层会导致模型过度拟合问题.在统计模型方面, 早期研究主要通过简单线性回归模型[21~24]来探索PM2.5-AOD之间的关系, 未考虑二者关系的时空变化特性[25].近年来, 很多研究采用一些高级统计模型来提高PM2.5-AOD相关关系, 如线性混合效应(linear mixed effects, LME)[26~30]、地理加权回归(geographical weighted regression, GWR)[10]、两阶段(two stage)[31~33]、时间地理加权回归(geographically and temporally weighted regression, GTWR)[34]和实时结构自适应(timely structure adaptive, TSA)[35]模型等.与经典统计模型相比, 这些模型能够对PM2.5-AOD关系进行局部校正, 从而达到较高的预测精度.进一步改进这些统计模型来提高PM2.5浓度时空分布预测精度仍然具有挑战性.

早期的混合效应模型(LME)通过引入时间随机效应来矫正PM2.5-AOD关系的时间差异性, 证明了良好的性能[26, 27].然而上述基于时间随机效应的混合效应模型忽略了PM2.5-AOD关系的空间差异性[25], 特别是在研究区域范围较大时, 上述空间差异性不容忽视.鉴于此, 本文采用时空混合效应(space-time linear mixed effects, STLME)模型[36], 该模型在时间随机效应基础上加入空间随机效应来同时校正PM2.5-AOD关系的时空异质性, 弥补了该类模型对空间差异性表达的缺陷, 为近地面PM2.5浓度反演提供了新的研究方法.

为了提高空间分辨率和覆盖度, 以便更好地研究城市内部或较小区域PM2.5浓度的暴露评估.本文引入Lyapustin等[13]开发的多角度大气校正算法(MAIAC)生产的空间分辨率为1km的AOD产品, 并辅以PM2.5气象和土地利用数据等协变量构建一种新的STLME模型, 预测2017年京津冀地区PM2.5浓度的时空分布, 以期为该地区健康风险评估和大气污染源识别提供科学依据和参考.

1 材料与方法 1.1 研究区概况

京津冀区域陆地面积为21.8万km2, 人口约1.127亿(2018年), 位于华北平原北部, 北靠燕山山脉, 西倚太行山, 包括北京市、天津市和河北省的11个地级市(图 1), 拥有得天独厚的地理位置.该区域属于温带大陆性季风气候, 夏季受海洋水汽影响, 高温多雨; 冬季受内陆干冷空气影响, 寒冷干燥.作为全国经济最发达的地区之一[37], 京津冀地区近40年快速发展的城市化和工业化, 使其得成为了中国大气污染最严重的区域之一[38].

图 1 研究区地理位置和空气质量监测站点及其2017年PM2.5年均值分布 Fig. 1 Geographical location and air quality monitoring stations in the study area and distribution of annual mean PM2.5 in 2017

1.2 数据获取 1.2.1 PM2.5监测数据

PM2.5日均浓度数据来源于全国空气质量实时发布平台(http://106.37.208.233:20035/).共获取研究区2017年79个国控监测点的PM2.5日均值数据, 剔除低于监测限值(日均值<2 μg ·m-3)的ρ(PM2.5)[32].

1.2.2 MAIAC AOD

MAIAC AOD产品来自NASA气候模拟中心(https://portal.nccs.nasa.gov/datashare/maiac/DataRelease/).采用h03v01、h03v02和h04v01这3个条带内的2017年京津冀地区AOD数据.为提高空间覆盖度, 采用最小二乘法对Terra和Aqua卫星的AOD数据分别逐月建立线性回归模型, 以此互相预测缺失的AOD值; 最后将每日Terra和Aqua卫星AOD通过计算平均值进行数据融合, 该值表征了当地10:30~13:30的MAIAC AOD平均状况[31], 增强了与PM2.5浓度日均值数据的相关性.

1.2.3 PM2.5时间-空间预测因子

PM2.5时间-空间预测因子包括气象因子、归一化植被指数和夜间灯光指数.其中气象因子包括大气边界层高度(boundary layer height, BLH)、地面气温、风速、气压和相对湿度等.其中, 大气边界层高度数据来源于第5版的戈达德地球观测系统数据同化系统的网格化气象数据产品(ftp://rain.ucis.dal.ca), 时空分辨率为0.25°×0.312 5°的小时数据; 地面气温、相对湿度、风速和气压等气象数据来源于国家青藏高原数据中心[39, 40](http://data.tpdc.ac.cn/zh-hans/data/), 时空分辨率为0.1°×0.1°的日均值.所有数据都使用双线性内插法重采样成与AOD相同分辨率的栅格数据.

归一化植被指数来源于NASA(https://modis.gsfc.nasa.gov/)提供的2017年MODIS产品, 产品序列号为MOD13A2, 时空分辨率为1 km的16 d合成图像.通过ENVI 5.3进行影像拼接、投影变换和裁剪等预处理, 得到京津冀地区的NDVI数据.

夜间灯光指数数据来源于NPP/VIIRS(National Polar-orbiting Partnership's Visible Infrared Imaging Radiometer Suite)月度合成夜间灯光数据(https://earthdata.nasa.gov/earth-observation-data/), 空间分辨率为0.5 km.原数据通过投影变换、重采样以及裁剪和消除负值等预处理, 最终得到京津冀地区的夜间灯光指数数据.

1.2.4 PM2.5空间预测因子

空间预测因子包括人口密度和海拔高度.人口数据为美国能源部橡树岭国家实验室(ORNL)(https://www.satpalda.com/product/landscan/)2010年1 km格网数据集.海拔数据为NASA(https://ladsweb.modaps.eosdis.nasa.gov/)航天飞机雷达地形测绘任务(SRTM)提供的分辨率30 m的DEM数据.

1.2.5 PM2.5和AOD数据预处理

AOD是解释PM2.5浓度最重要的因子, 为优化每日PM2.5观测值和AOD之间的关系, 在建模之前, 需要对数据进行预处理.AOD数据在反演过程中由于云污染(无法观测到地表)、水或者积雪(反射率过高无法区分地表和大气的贡献)等原因, 会导致数据异常.为了确保高质量的AOD数据, 上述被污染的数据则需要删除[41, 42].首先, 考虑到PM2.5-AOD两者正相关关系, 删除二者关系明显相冲突(AOD<第50个百分位数和PM2.5>第85个百分位数, 以及PM2.5< 第50个百分位数和AOD>第80个百分位数)的数据; 其次, 剔除每个监测站点少于30d有效记录的数据; 最后, 剔除因水面和冰雪覆盖(NDVI≤0)造成异常的AOD数据.通过上述处理, 本文删除了约4%的AOD异常数据.

1.3 研究方法 1.3.1 STLME模型

STLME模型包括固定效应和随机效应, 其中随机效应又包含时间随机效应和空间随机效应.STLME模型与流行的LME模型不同之处在于, 在时间随机效应中嵌套空间随机效应来表征PM2.5-AOD关系的空间异质性.因此STLME模型能够表征PM2.5-AOD关系在京津冀随时空变化的特性.STLME模型公式如下:

(1)

式中, PMij为研究区第i个地面监测站点第j天的PM2.5浓度监测值; aβ1为截距项和AOD项的固定效应, 分别表示固定截距和固定斜率; μjvj为第j天的随机效应, 分别表示随机截距和随机斜率; gjreghjreg分别为次级区域第j天的随机截距和随机斜率, 其中gjreghjreg嵌套在时间随机效应μjvj中; AODij为第i个地面监测站点所在网格(1 km×1 km)第j天遥感反演的气溶胶光学厚度值; x1mii网格点第m个空间预测因子的值(即人口密度、海拔); x2mij为第i个网格点第j天第m个时间-空间预测因子值(即边界层高度、地面气温、相对湿度、风速、气压、归一化植被指数和夜间灯光指数), γ为各个时空预测因子对PM2.5浓度影响的固定效应. εij是站点i在第j天的随机误差项; σ1为时间随机效应协方差矩阵, 其对角元素为σμσv; σ2为空间随机效应协方差矩阵, 其对角元素为σgσh.模型中的固定效应代表各个因子对PM2.5浓度的固定影响, 它们并不随时空变化而变化; 随机效应μjvj随时间变化而改变但并不随空间变化, 体现PM2.5-AOD关系随时间的变化; 随机效应gjreghjreg随着站点所在次区域的改变而改变, 反映了PM2.5-AOD关系在空间上的差异.

1.3.2 次区域划分方案优选与PM2.5预测结果的边界平滑

由于大气污染物排放强度、颗粒物组成成分和颗粒物浓度垂直廓线等因素存在次区域尺度上的空间差异, 从而导致PM2.5-AOD关系存在次区域尺度的空间异质性[36].因此, 次区域优选就成为STLME模型构建的重要一环, 本研究在初步选择多个次区域划分方案的基础上, 采用模型预测精度验证法进行次区域划分方案优选, 选出模型预测精度最高的次区域划分方案.

由于次区域划分方案一经确定其范围和边界就是固定的, 有研究表明[36], 如果相邻次区域之间PM2.5浓度差异较大时, 容易在相邻边界两侧产生“断崖式”差异, 这种现象与实际情况不符.为了消除这一问题, 本研究以LME模型预测结果为背景参照, 通过在次区域边界两侧设置缓冲区并计算相邻两区域的预测均值替代原单一预测值, 依据背景参照确定缓冲区的最佳宽度, 从而解决上述“断崖式”差异问题.

1.3.3 模型验证

为消除模型的过拟合问题, 采用十折交叉验证法来验证STLME模型预测精度.本文用模型预测值与实测值间的决定系数R2、均方根预测误差(RMSE)和相对预测误差(RPE)等统计指标来评估.其中R2是预测值与观测值的拟合程度, R2预测越接近于1, 表明预测精度越高; 斜率解释模型的预测偏差, 斜率越接近于1, 表明预测偏差越小; RMSE和RPE分别反映模型预测的绝对误差和相对误差, 二者数值越小, 模型精度越高.其中, RMSE和RPE的计算公式如下:

(2)
(3)

式中, yobs, i为地面各监测站点PM2.5浓度监测值; ymodel, i为各站点PM2.5浓度模型预测值; n为建模数据集的样本总数, y为PM2.5浓度的监测均值.

为探查STLME模型相比于LME模型是否提高模型精度, 本文用时间维决定系数(时间R2)和空间维决定系数(空间R2)来评估STLME模型对校正PM2.5-AOD关系的时空差异效果.其中, 时间R2是所有站点每天ΔPM2.5与预测ΔPM2.5回归拟合所得, ΔPM2.5是所有站点第j天的PM2.5浓度的平均值, 预测ΔPM2.5是所有站点第j天的预测PM2.5浓度的平均值; 空间R2为各站点全年ΔPM2.5与预测ΔPM2.5回归拟合所得.

2 结果与分析 2.1 描述性统计分析

京津冀地区2017年PM2.5地面监测站点的有效观测量为28 275个, AOD的有效观测量17 288条.顺应建模要求, 在栅格匹配过程中剔除了不满足4个分区均有有效数据条件的天数, 共删除50 d, 剩余315 d(86%的时间覆盖率)16 163条(61.86%的空间覆盖率)有效建模数据记录.

表 1为建模数据中PM2.5及AOD的描述性统计分析.京津冀地区各城市的PM2.5浓度之间的差异很大, 其中, 邯郸ρ(PM2.5)的均值最高(86.50μg ·m-3), 保定次之, 张家口最低, 空气质量最好.AOD中邯郸的均值同样最大, 邢台次之, 张家口最小(0.23).从整体来看, PM2.5分布规律与AOD空间分布大体一致, 都是南高北低.

表 1 2017年建模数据中PM2.5和AOD的描述性统计量 Table 1 Descriptive statistics of PM2.5 and AOD in modeling data in 2017

图 2为STLME模型的各变量的直方图和描述性统计.结果表明, 除地面气温外, 其余变量均近似正态分布.其中, AOD平均值为0.50, 标准差为0.49, 时空差异较大; ρ(PM2.5)的范围为2~499μg ·m-3, 标准差为47.56μg ·m-3, 说明研究区内PM2.5污染程度时空差异同样较大. ρ(PM2.5)年均值为62.20μg ·m-3, 约是国家环境空气质量(GB 3095-2012)二级标准限值(35μg ·m-3)的1.8倍, 表明京津冀地区PM2.5污染较为严重.

图 2 京津冀地区2017年建模数据各参数的直方图和描述性统计 Fig. 2 Histograms and descriptive statistics of modeling data parameters for the Beijing-Tianjin-Hebei region in 2017

2.2 模型诊断

图 3为模型拟合结果的随机误差项诊断, LME和STLME模型随机误差频率分布直方图均符合均值为0的正态分布的要求, 随机误差分布呈现随机性特征, 说明模型模拟结果可靠.

图 3 LME和STLME模型拟合结果的随机误差项诊断 Fig. 3 Random error term diagnosis of the fitting results of the LME and STLME models

表 2为STLME模型固定效应的解, 为使不同的变量具有相同的尺度, 对原始数据进行标准化处理.固定效应(固定截距和固定斜率)表示全年内所有自变量对PM2.5浓度的固定影响, 其中MAIAC AOD固定斜率标准值为0.982, 说明PM2.5-AOD两者呈正相关, 且两者相关关系显著; BLH、气压、风速和归一化植被指数等与PM2.5浓度呈负相关; 地面气温、相对湿度和海拔等与PM2.5浓度呈正相关.其中, 相对湿度、归一化植被指数、风速和地面气温等因子相关性显著.

表 2 STLME模型固定效应的解 Table 2 Fixed effects in the STLME model

2.3 模型次区域划分方案拟合及交叉验证

为了识别不同分区方案模型拟合结果差异, 确定最优分区方案, 在保持京津和河北省11个地级市行政区划完整的基础上, 共选取了4种方案.方案一:将京津冀地区按地理条件和发展水平差异划分成京津(北京、天津)、河北山地(承德、张家口)、河北内陆(石家庄、保定、衡水、廊坊、邢台、邯郸)和河北沿海(秦皇岛、唐山、沧州)这4个次区域; 方案二:按照地貌类型和发展水平差异把京津冀分成京津平原(北京、天津)、河北平原(石家庄、保定、衡水、廊坊、邢台、邯郸、皇岛、唐山、沧州)和山地(承德、张家口)这3个次区域; 方案三:按扩散条件差异将京津冀分成沿海(秦皇岛、唐山、沧州和天津)、内陆平原(石家庄、保定、衡水、廊坊、邢台、邯郸、北京)和内陆山地(承德、张家口)这3个次区域; 方案四:按照地貌类型差异将京津冀划分成山地(承德、张家口)和平原(石家庄、保定、衡水、廊坊、邢台、邯郸、北京、秦皇岛、唐山、沧州和天津)这2个次区域.从图 4可知, 4种分区方案的STLME模型拟合R2分别为0.94、0.93、0.93和0.92, RMSE分别为11.90、14.09、13.91和15.08 μg ·m-3; 经十折交叉验证(CV)后R2分别为0.91、0.90、0.90和0.90, RMSE分别为14.43、16.49、16.12和16.19 μg ·m-3.上述结果表明:不同分区方案影响STLME模型的预测精度, 其中, 方案一模拟效果最好, 精度最高, 因此本文采用方案一的次区域划分方案.

图 4 不同分区方案拟合和交叉验证结果对比 Fig. 4 Comparison of the fitting and cross validation results of different regional division schemes

2.4 LME模型和STLME模型比较

为了更深入了解STLME模型对区域因子的敏感性, 图 5为LME模型拟合与交叉验证结果.结果表明:当随机效应中只有时间效应时, LME模型拟合及交叉验证后的斜率分别是0.90和0.87, 决定系数R2分别为0.91和0.87, RMSE分别为15.59μg ·m-3和19.23μg ·m-3, RPE分别为25.92%和31.52%, 精度较高; 当加入空间随机效应后, STLME模型(图 4中方案一)的拟合及交叉验证的R2分别提高到0.94和0.91.反映出STLME模型优于LME模型, LME模型只能校正PM2.5-AOD关系的时间差异, 而STLME模型既能校正时间差异, 又能校正空间差异, 提高了模型精度.

图 5 2017年LME模型拟合及交叉验证结果 Fig. 5 Results of LME model fitting and cross validation in 2017

LME及STLME模型时空维统计结果表明(表 3), LME模型时间CV R2为0.99, 空间CV R2为0.87, 空间CV RMSE为6.20μg ·m-3; STLME模型时间CV R2为0.99, 空间CV R2为0.91, 空间CV RMSE为4.96μg ·m-3.可以看出在模型中加入区域因子后, 空间R2有所提高, 空间RMSE和时间RMSE有所降低, 说明STLME模型能够揭示PM2.5-AOD关系的时空异质性.

表 3 2017年LME及STLME模型时空维统计结果 Table 3 Spatiotemporal dimensional statistical results of LME and STLME models in 2017

为比较两种模型在站点尺度上的预测精度, 图 6为LME和STLME模型站点预测误差空间分布.其中, STLME和LME模型预测误差范围分别为-8.56~10.24和-15.19~13.08, 从中可以看出, STLME模型预测误差相对较小, 预测精度更高.并且, 在站点预测误差的空间分布上, LME模型表现为非随机性, 张家口、承德和秦皇岛等低值区为正误差, 保定和石家庄等高值区为负误差; 而STLME模型预测误差空间分布则呈现随机性特征, 说明STLME模型预测精度和稳健性更高.

图 6 LME和STLME模型预测误差对比 Fig. 6 Comparison of prediction errors of the LME and STLME model

图 7为STLME模型交叉验证后, ρ(PM2.5)的预测值和实测值在月、季和年尺度上的拟合状况.月、季和年尺度的R2都是0.95; RMSE分别是6.095、4.848和3.791μg ·m-3; RPE分别是10.19%、8.194%和6.243%.从结果上看, 随着时间尺度的增大, 模型的RMSE和RPE逐渐减小, 预测精度升高.

图 7 基于STLME模型交叉验证结果在月、季和年尺度上站点的观测值和预测值拟合 Fig. 7 Fitting results between predicted and observed PM2.5 concentration on monthly, seasonal, and yearly scales based on the results of cross-validation in the STLME model

2.5 次区域边界平滑结果

图 8(a)(对照)和图 8(b)可以看出, STLME模型出现“断崖式”差异问题主要分布在京津与河北平原、唐秦与承德的交界处.因此本文重点对上述界线两侧创建缓冲区, 缓冲区内采用相邻两个区域PM2.5浓度预测结果的算数平均值消除界线两侧“断崖式”差异.以LME模型预测结果在次区域界线两侧的变化趋势为参照, 尝试不同的缓冲区范围, 最终选择与背景参照变化趋势相匹配的宽度90km的缓冲区方案[图 8(c)], 较好展现了次区域界线两侧预测结果的渐变特征.

图 8 2017年LME模型冬季PM2.5浓度预测结果和STLME模型平滑处理前后PM2.5浓度预测结果 Fig. 8 PM2.5 estimations by the LME model, STLME model, and STLME model with smoothing treatment during winter in 2017

2.6 京津冀PM2.5浓度时空分布

图 9为基于STLME模型预测京津冀地区PM2.5浓度年均值分布, 2017年京津冀全域ρ(PM2.5)年均值为62.20μg ·m-3.低值区主要分布在北部燕山山区与西部太行山区, ρ(PM2.5)年均值在25~35μg ·m-3之间; 次低值区为北京和天津等地, ρ(PM2.5)年均值在35~50μg ·m-3之间; 而高值区主要分布在河北南部, 以石家庄、邢台和邯郸为主, ρ(PM2.5)年均值高于70μg ·m-3, 又以邯郸污染最为严重.

图 9 2017年PM2.5浓度年均值分布 Fig. 9 Spatial distribution of annual mean PM2.5 concentrations in 2017

京津冀地区PM2.5污染具有季节性特征, 春(3~5月)、夏(6~8月)、秋(9~11月)和冬季(12月~次年2月)的PM2.5浓度平均值分布如图 10.从中可以看出2017年京津冀地区PM2.5浓度季均值差异明显, 呈现“冬高夏低, 春秋过渡”的特征.冬季ρ(PM2.5)最高(74.52μg ·m-3), 原因是由于冬季取暖增加了污染物排放和不利的天气条件(逆温和沙尘天气等); 夏季污染最轻(42.72μg ·m-3), 主要源于大气扩散条件好和降水的清洁作用.

图 10 2017年PM2.5浓度季均值空间分布 Fig. 10 Spatial distribution of seasonal mean PM2.5 concentrations in 2017

3 讨论 3.1 模型中分区方案与缓冲区设置

基于多个分区方案的STLME模型优化, 不仅提高了模型的预测精度, 也识别了影响PM2.5-AOD关系区域差异的主要因子, 从而有利于京津冀地区大气污染治理的精准施策.从最优分区方案可以获得以下认识:气溶胶组成成分的区域差异是分区的重要因子, 主要区分了京津次区域与河北省各个次区域的差异, 后者仍然以煤烟型气溶胶成分为主, 导致河北平原特别是南部平原PM2.5浓度明显高于京津地区; 由大气扩散条件差异导致的边界层大气气溶胶垂直廓线结构的区域差异也是分区的重要因子, 主要区分了河北内陆和河北沿海两大次区域的差异, 其中, 河北沿海由于海陆热力性质差异导致大气扩散条件较好, 区域大气环境容量较大, PM2.5浓度明显低于河北内陆次区域; 而河北内陆恰好相反, 由于西侧太行山地的阻挡和动力作用, 该区域大气扩散条件差, 冬季常常导致区域大气污染物的堆积, 大气环境容量小, 成为该区大气污染治理难度大的主要原因之一.

STLME模型由于同时校正了PM2.5-AOD关系的时空异质性, 进一步提高了混合效应模型的预测精度.但该方法也存在一个技术性缺陷, 即容易出现次区域边界两侧PM2.5浓度预测值的“断崖式”差异, 这与实际情况相悖.本文提出的缓冲区平滑方法消除了这一缺陷, 完善了STLME模型.

3.2 模型预测精度国内外比较

与本研究使用相同或相近研究方法的研究案例中(表 4), STLME模型优于LME模型, 例如Zheng等[30]使用LME模型(CV R2为0.77)模拟京津冀地区PM2.5浓度和Ma等[43]在长江三角洲采用LME模型的研究(CV R2为0.67).本模型的性能也超过了Ma等[44]人使用的两阶段模型结果, 其R2(0.82)和CV R2(0.79)略低于STLME模型.与Kloog等[36]使用相同模型研究2003~2013年以色列PM2.5时空变化相比(CV R2=0.72), R2显著提高, 时间R2和空间R2均显著提升, 但是斜率稍有降低(CV Slope=0.99), 可能由于京津冀地区比以色列大气污染严重致使AOD非随机缺失严重, 从而使预测偏差加大.

表 4 STLME模型与相近或相同方法对PM2.5浓度时空变化研究结果对比1) Table 4 Comparison between STLME model and similar or identical methods for the spatiotemporal variation in PM2.5 concentrations

随着人工智能的发展, 许多机器学习模型被用于PM2.5浓度预测.表 5为本研究与京津冀区域内使用各种机器学习模型的研究结果对比, 其中, 最流行的机器学习模型是反向传播神经网络、人工神经网络、深度学习和随机森林等.相对于STLME模型, 人工神经网络[45, 46](CV R2=0.55)在估算PM2.5浓度方面较差, 决策树[47](CV R2=0.88)和随机森林[48](CV R2=0.86)模型表现出较好的性能, 但本研究仍然表现出更高的模拟精度.

表 5 STLME模型与机器学习模型对京津冀区域PM2.5浓度时空变化研究结果对比 Table 5 Comparison between STLME model and machine learning model on spatiotemporal variation of PM2.5 concentrations in BTH

4 结论

(1) 与传统LME模型相比, 加入空间随机效应后构建的STLME模型在近地面PM2.5浓度估算精度方面有了很大提升, CV R2由0.87提升到0.91, 该模型在校正PM2.5-AOD关系的时空异质性方面具有很大优势, 可用于PM2.5浓度预测研究.

(2) STLME模型识别出京津冀PM2.5-AOD关系的主要空间差异, 表现为京津地区、河北山地区、河北内陆平原区和河北沿海平原区4个次区域差异.其中, 京津地区与河北平原区差异的主要影响因子为产业结构与能源结构, 河北沿海与内陆平原区差异的主要影响因子为大气扩散条件.

(3) 2017年京津冀地区PM2.5浓度时空分布差异性较大, 高值区主要分布在以石家庄、邢台和邯郸为中心的河北南部地区, 低值区主要分布在燕山-太行山山区; 冬季是PM2.5污染最为严重的季节, 其次是秋季和春季, 夏季污染最轻.

参考文献
[1] Zhang Y Q, Ding Z, Xiang Q Q, et al. Short-term effects of ambient PM1 and PM2.5 air pollution on hospital admission for respiratory diseases: case-crossover evidence from Shenzhen, China[J]. International Journal of Hygiene and Environmental Health, 2020, 224. DOI:10.1016/j.ijheh.2019.11.001
[2] Lim J M, Jeong J H, Lee J H, et al. The analysis of PM2.5 and associated elements and their indoor/outdoor pollution status in an urban area[J]. Indoor Air, 2011, 21(2): 145-155. DOI:10.1111/j.1600-0668.2010.00691.x
[3] Lelieveld J, Evans J S, Fnais M, et al. The contribution of outdoor air pollution sources to premature mortality on a global scale[J]. Nature, 2015, 525(7569): 367-371. DOI:10.1038/nature15371
[4] Feng H H, Zou B. Satellite-based estimation of the aerosol forcing contribution to the global land surface temperature in the recent decade[J]. Remote Sensing of Environment, 2019, 232. DOI:10.1016/j.rse.2019.111299
[5] Peng J, Chen S, Lü H L, et al. Spatiotemporal patterns of remotely sensed PM2.5 concentration in China from 1999 to 2011[J]. Remote Sensing of Environment, 2016, 174: 109-121. DOI:10.1016/j.rse.2015.12.008
[6] He Q Q, Gu Y F, Zhang M. Spatiotemporal trends of PM2.5 concentrations in central China from 2003 to 2018 based on MAIAC-derived high-resolution data[J]. Environment International, 2020, 137. DOI:10.1016/j.envint.2020.105536
[7] Cui Y Q, Sun Q H, Liu Z G. Ambient particulate matter exposure and cardiovascular diseases: a focus on progenitor and stem cells[J]. Journal of Cellular and Molecular Medicine, 2016, 20(5): 782-793. DOI:10.1111/jcmm.12822
[8] Ma Z W, Hu X F, Huang L, et al. Estimating ground-level PM2.5 in China using satellite remote sensing[J]. Environmental Science & Technology, 2014, 48(13): 7436-7444.
[9] Zhang Y, Li Z Q. Remote sensing of atmospheric fine particulate matter (PM2.5) mass concentration near the ground from satellite observation[J]. Remote Sensing of Environment, 2015, 160: 252-262. DOI:10.1016/j.rse.2015.02.005
[10] Hu X F, Waller L A, Al-Hamdan M Z, et al. Estimating ground-level PMPM2.5concentrations in the Southeastern U.S.using geographically weighted regression[J]. Environmental Research, 2013, 121: 1-10. DOI:10.1016/j.envres.2012.11.003
[11] Song W Z, Jia H F, Huang J F, et al. A satellite-based geographically weighted regression model for regional PM2.5 estimation over the Pearl River Delta region in China[J]. Remote Sensing of Environment, 2014, 154: 1-7. DOI:10.1016/j.rse.2014.08.008
[12] Wu J S, Yao F, Li W F, et al. VIIRS-based remote sensing estimation of ground-level PM2.5 concentrations in Beijing-Tianjin-Hebei: a spatiotemporal statistical model[J]. Remote Sensing of Environment, 2016, 184: 316-328. DOI:10.1016/j.rse.2016.07.015
[13] Lyapustin A, Wang Y J, Korkin S, et al. MODIS collection 6 MAIAC algorithm[J]. Atmospheric Measurement Techniques, 2018, 11(10): 5741-5765. DOI:10.5194/amt-11-5741-2018
[14] Kloog I, Chudnovsky A A, Just A C, et al. A new hybrid spatio-temporal model for estimating daily multi-year PM2.5concentrations across Northeastern USA using high resolution aerosol optical depth data[J]. Atmospheric Environment, 2014, 95: 581-590. DOI:10.1016/j.atmosenv.2014.07.014
[15] Xie Y Y, Wang Y X, Zhang K, et al. Daily estimation of ground-level PM2.5 concentrations over Beijing using 3 km resolution MODIS AOD[J]. Environmental Science & Technology, 2015, 49(20): 12280-12288.
[16] Liang F C, Xiao Q Y, Wang Y J, et al. MAIAC-based long-term spatiotemporal trends of PM2.5 in Beijing, China[J]. Science of the Total Environment, 2018, 616-617: 1589-1598. DOI:10.1016/j.scitotenv.2017.10.155
[17] Ni X L, Cao C X, Zhou Y K, et al. Spatio-temporal pattern estimation of PM2.5 in Beijing-Tianjin-Hebei region based on MODIS AOD and meteorological data using the back propagation neural network[J]. Atmosphere, 2018, 9(3). DOI:10.3390/atmos9030105
[18] Zang L, Mao F Y, Guo J P, et al. Estimating hourly PM1 concentrations from Himawari-8 aerosol optical depth in China[J]. Environmental Pollution, 2018, 241: 654-663. DOI:10.1016/j.envpol.2018.05.100
[19] Wang P, Zhang H, Qin Z D, et al. A novel hybrid-garch model based on ARIMA and SVM for PM2.5 concentrations forecasting[J]. Atmospheric Pollution Research, 2017, 8(5): 850-860. DOI:10.1016/j.apr.2017.01.003
[20] Brokamp C, Jandarov R, Hossain M, et al. Predicting daily urban fine particulate matter concentrations using a random forest model[J]. Environmental Science & Technology, 2018, 52(7): 4173-4179.
[21] Engel-Cox J A, Holloman C H, Coutant B W, et al. Qualitative and quantitative evaluation of MODIS satellite sensor data for regional and urban scale air quality[J]. Atmospheric Environment, 2004, 38(16): 2495-2509. DOI:10.1016/j.atmosenv.2004.01.039
[22] Gupta P, Christopher S A. Particulate matter air quality assessment using integrated surface, satellite, and meteorological products: multiple regression approach[J]. Journal of Geophysical Research: Atmospheres, 2009, 114(D14). DOI:10.1029/2008JD011496
[23] 金囝囡, 杨兴川, 晏星, 等. 京津冀及周边MAIAC AOD和PM2.5质量浓度特征及相关性分析[J]. 环境科学, 2021, 42(6): 2604-2615.
[24] Xin J Y, Zhang Q, Wang L L, et al. The empirical relationship between the PM2.5 concentration and aerosol optical depth over the background of north China from 2009 to 2011[J]. Atmospheric Research, 2014, 138: 179-188. DOI:10.1016/j.atmosres.2013.11.001
[25] Gryparis A, Paciorek C J, Zeka A, et al. Measurement error caused by spatial misalignment in environmental epidemiology[J]. Biostatistics, 2009, 10(2): 258-274. DOI:10.1093/biostatistics/kxn033
[26] Lee H J, Chatfield R B, Strawa A W. Enhancing the applicability of satellite remote sensing for PM2.5 estimation using MODIS Deep Blue AOD and land use regression in California, United States[J]. Environmental Science & Technology, 2016, 50(12): 6546-6555.
[27] 郝静, 孙成, 郭兴宇, 等. 京津冀内陆平原区PM2.5浓度时空变化定量模拟[J]. 环境科学, 2018, 39(4): 1455-1465.
[28] 孙成, 王卫, 刘方田, 等. 基于线性混合效应模型的河北省PM2.5浓度时空变化模型研究[J]. 环境科学研究, 2019, 32(9): 1500-1509.
[29] Li R, Gong J H, Chen L F, et al. Estimating ground-level PM2.5 using fine-resolution satellite data in the megacity of Beijing, China[J]. Aerosol and Air Quality Research, 2015, 15(4): 1347-1356. DOI:10.4209/aaqr.2015.01.0009
[30] Zheng Y X, Zhang Q, Liu Y, et al. Estimating ground-level PM2.5 concentrations over three megalopolises in China using satellite-derived aerosol optical depth measurements[J]. Atmospheric Environment, 2016, 124: 232-242. DOI:10.1016/j.atmosenv.2015.06.046
[31] 杨晓辉, 宋春杰, 范丽行, 等. 京津冀地区高分辨率PM2.5浓度时空变化模拟与分析[J]. 环境科学, 2021, 42(9): 4083-4094.
[32] Hu X F, Waller L A, Lyapustin A, et al. Estimating ground-level PM2.5 concentrations in the Southeastern United States using MAIAC AOD retrievals and a two-stage model[J]. Remote Sensing of Environment, 2014, 140: 220-232. DOI:10.1016/j.rse.2013.08.032
[33] Yao F, Wu J S, Li W F, et al. A spatially structured adaptive two-stage model for retrieving ground-level PM2.5 concentrations from VIIRS AOD in China[J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2019, 151: 263-276. DOI:10.1016/j.isprsjprs.2019.03.011
[34] He Q Q, Huang B. Satellite-based high-resolution PM2.5 estimation over the Beijing-Tianjin-Hebei region of China using an improved geographically and temporally weighted regression model[J]. Environmental Pollution, 2018, 236: 1027-1037. DOI:10.1016/j.envpol.2018.01.053
[35] Fang X, Zou B, Liu X P, et al. Satellite-based ground PM2.5 estimation using timely structure adaptive modeling[J]. Remote Sensing of Environment, 2016, 186: 152-163. DOI:10.1016/j.rse.2016.08.027
[36] Kloog I, Sorek-Hamer M, Lyapustin A, et al. Estimating daily PM2.5 and PM10 across the complex geo-climate region of Israel using MAIAC satellite-based AOD data[J]. Atmospheric Environment, 2015, 122: 409-416. DOI:10.1016/j.atmosenv.2015.10.004
[37] 陈辉, 厉青, 李营, 等. 京津冀及周边地区PM2.5时空变化特征遥感监测分析[J]. 环境科学, 2019, 40(1): 33-43.
[38] 郝静. 京津冀地区PM2.5浓度时空变化定量模拟[D]. 石家庄: 河北师范大学, 2018.
[39] He J, Yang K, Tang W J, et al. The first high-resolution meteorological forcing dataset for land process studies over China[J]. Scientific Data, 2020, 7(1). DOI:10.1038/s41597-020-0369-y
[40] Yang K, He J, Tang W J, et al. On downward shortwave and longwave radiations over high altitude regions: observation and modeling in the Tibetan Plateau[J]. Agricultural and Forest Meteorology, 2010, 150(1): 38-46. DOI:10.1016/j.agrformet.2009.08.004
[41] De Hoogh K, Héritier H, Stafoggia M, et al. Modelling daily PM2.5 concentrations at high spatio-temporal resolution across Switzerland[J]. Environmental Pollution, 2018, 233: 1147-1154. DOI:10.1016/j.envpol.2017.10.025
[42] Stafoggia M, Schwartz J, Badaloni C, et al. Estimation of daily PM10 concentrations in Italy (2006-2012) using finely resolved satellite data, land use variables and meteorology[J]. Environment International, 2017, 99: 234-244. DOI:10.1016/j.envint.2016.11.024
[43] Ma Z W, Liu Y, Zhao Q Y, et al. Satellite-derived high resolution PM2.5concentrations in Yangtze River Delta region of China using improved linear mixed effects model[J]. Atmospheric Environment, 2016, 133: 156-164. DOI:10.1016/j.atmosenv.2016.03.040
[44] Ma Z W, Hu X F, Sayer A M, et al. Satellite-based spatiotemporal trends in PM2.5 concentrations: China, 2004-2013[J]. Environmental Health Perspectives, 2016, 124(2): 184-192. DOI:10.1289/ehp.1409481
[45] Gupta P, Christopher S A. Particulate matter air quality assessment using integrated surface, satellite, and meteorological products: 2. A neural network approach[J]. Journal of Geophysical Research-atmospheres, 2009, 114(D20). DOI:10.1029/2008JD011497
[46] Wang X P, Sun W B. Meteorological parameters and gaseous pollutant concentrations as predictors of daily continuous PM2.5 concentrations using deep neural network in Beijing-Tianjin-Hebei, China[J]. Atmospheric Environment, 2019, 211: 128-137. DOI:10.1016/j.atmosenv.2019.05.004
[47] Ding Y, Chen Z Q, Lu W F, et al. A catBoost approach with wavelet decomposition to improve satellite-derived high-resolution PM2.5estimates in Beijing-Tianjin-Hebei[J]. Atmospheric Environment, 2021, 249. DOI:10.1016/J.ATMOSENV.2021.118212
[48] Zhao C, Wang Q, Ban J, et al. Estimating the daily PM2.5 concentration in the Beijing-Tianjin-Hebei region using a random forest model with a 0.01°×0.01° spatial resolution[J]. Environment International, 2020, 134. DOI:10.1016/j.envint.2019.105297