2020, Vol. 41
2. 云南师范大学西部资源环境地理信息技术教育部工程研究中心, 昆明 650500;
3. 云南电网有限责任公司信息中心, 昆明 650217
2. GIS Technology Engineering Research Centre for West-China Resources and Environment, Ministry of Education, Yunnan Normal University, Kunming 650500, China;
3. Information Center of Yunnan Power Grid Co., Ltd., Kunming 650217, China
在过去近30年里, 中国以“高能耗、高污染”的方式实现了工业化、城市化和社会经济的快速发展, 发展过程中伴随着城市空间盲目扩张、机动车拥有量快速增加、能源、资源消耗迅速增长[1]和城市植被覆盖率快速下降[2]等现象, 间接引发了一系列严峻的包括空气污染在内的环境问题.作为最大的发展中国家, 中国正面临着严重的空气污染问题[3].近年来中国多地空气重度污染事件频发, 从2012年开始, 全国338个地级及以上城市陆续建成了PM2.5地面监测网络, 并开展了空气质量新标准监测. 2017年在全国338个地级及以上城市中, 以PM2.5为首要污染物的天数占重度污染及以上污染天数的74.2%, PM2.5已经成为我国大部分城市环境空气质量的首要污染物[4].PM2.5污染给交通运输、自然环境、人类健康和社会可持续发展带来严重影响, 因此全面准确地了解我国PM2.5污染的时空变化情况对控制PM2.5污染和确保我国的可持续发展至关重要.
连续、长时间和大范围的PM2.5浓度数据是研究PM2.5污染时空变化规律的基础, 然而我国PM2.5浓度地面监测网络建成时间较晚, 且地面监测站点覆盖范围较小, 空间分布不均匀, 可获取的传统地面监测数据通常难以满足PM2.5相关研究的需要.气溶胶光学厚度(aerosol optical depths, AOD)与PM2.5有很强的相关性[5, 6], 通过卫星遥感数据产品AOD来估算近地面PM2.5浓度为从长时间序列、大空间尺度研究PM2.5浓度提供了一种低成本且有效的方法[5~9].近年来, 国内外相关学者以全球不同的地区[9~20]为研究对象, 使用多元回归[10]、非线性模型[11]、线性混合效应[12]、两阶段空间统计模型[13]、两阶段广义相加模型[14]、地理加权回归模型[7]、化学传输模型[15]、ANN[10]、BP神经网络[12]和支持向量回归(SVR)[16]等模型, 建立了PM2.5和卫星遥感反演的AOD之间的关系模型, 研究结果表明AOD和PM2.5之间的关系存在明显的时间差异和空间差异[17, 18].有研究表明, 在AOD-PM2.5关系模型中加入气象参数可以有效提高模型估算精度[7, 13~15, 19]；Gupta等[10]使用MODIS AOD和气象数据估算了美国东南部的PM2.5浓度, 结果表明人工神经网络估算效果比简单线性回归和多元回归更好；郭建平等[12]使用中国东部地区5个大气成分站气象数据和11个PM2.5地面监测站数据, 结合2007~2008年MODIS/Terra AOD, 针对各个大气成分站分别构建了BP神经网络模型估算中国PM2.5浓度, 5个模型的拟合系数在0.4~0.6之间；Nguyen等[16]使用多元线性回归(MLR)和SVR模型估算了越南河内的细微颗粒(PM2.5和PM10)的浓度, 研究结果表明SVR的回归效果普遍优于MLR.
尽管已经有大批学者在我国PM2.5污染方面开展了大量的研究, 但这些研究大多关注于PM2.5污染严重的区域[14~20], 或直接以全国作为一个整体为研究对象[7, 20].我国的东部和西部、沿海和内陆地区在自然条件、经济水平和人口规模等方面存在着明显差异, 造成PM2.5浓度变化特征差异明显.根据自然和社会经济同质性将中国完全划分为多个子区域, 然后对PM2.5的长期时空变化特征进行估算和分析是有必要的.同时, 从大气污染防治行动计划[21]发布以来, 各地空气污染状况已经在一定程度上得到缓解, 很明显2000~2017年我国PM2.5浓度是一个阶段性动态变化的过程, 现有研究主要通过整个时间序列数据进行趋势分析[19], 而线性回归只能评估数据序列在某时间段内的整体变化, 无法客观地体现指定时间段内的逐步变化过程.
在前期研究基础上[22~24], 本研究以2013年1月至2017年12月内所有可获取的国控空气质量地面监测站的实测PM2.5浓度数据、2000年3月至2017年12月内卫星遥感反演的MODIS AOD数据以及气象参数(温度、风速、风向、边界层高度和相对湿度)为基础数据集, 综合考虑AOD与PM2.5关系的季节性和区域性差异性, 基于支持向量回归机与思维进化算法优化后的BP人工神经网络, 构建了一个全国尺度下近地面PM2.5浓度的二阶段组合估算模型.在此基础上, 从长时间和大空间尺度分析了中国PM2.5污染演变的时空变化过程, 以期为开展全球PM2.5时空演变过程和特征分析提供方法借鉴, 并为揭示PM2.5形成机制、空间异质性以及其它与PM2.5相关的研究提供数据基础.
1 材料与方法 1.1 研究区域本研究以中国为研究区域, 经纬度范围为73°~135°E、17°~54°N.国务院发展研究中心根据各省自然条件及经济的同质性于2005年提出了对大陆31个省的八大经济区域划分方式, 八大经济区域包括:东北地区(northeast region, NE)、北部沿海(northern coastal area, NC)、东部沿海(eastern coastal area, EC)、南部沿海(southern coastal area, SC)、黄河中游(middle yellow river, MYR)、长江中游(middle changjiang river, MCR)、西南地区(southwest China, SW)和大西北地区(northwest China, NW).图 1为八大经济区域、国控空气质量地面监测站点的空间分布, 2000~2017年内每年各经济区域的GRP(地区生产总值)占我国GDP(国内生产总值, 不包括港澳台地区)的占比情况(各省地区生产总值数据来自于国家数据官网[25]).可以看到, 区域GRP占比最高的区域主要为EC和NC地区；而在SC和MCR地区的GRP占比较高；在SW和MYR地区的GRP占比较低；在NE和NW地区的GRP占比最低.
|
图 1 研究区域概况 Fig. 1 Overview of study area |
本研究从大气存档和分配系统[26]获取了2000年3月至2017年12月期间内月均尺度的MOD08 M3和MYD08 M3文件(C6, 空间分辨率为1°×1°), 并选用了使用深蓝算法反演的550nm处陆地气溶胶数据产品(参数名:Deep_Blue_Aerosol_Optical_ Depth_550_Land_Mean_ Mean).
1.2.2 气象数据气象数据(如温度、风速、风向、边界层高度和相对湿度)能够有效提高AOD-PM2.5关系模型的精度, 从欧洲中期天气预报中心[27]获取了2000年3月至2017年12月的5种空间分辨率为1°×1°的月均气象数据作为本研究中PM2.5浓度估算模型的协变量, 包括2 m气温(T2M, K)、10 m风的东向分量(10U, m·s-1)、10 m风的北向分量(10 V, m·s-1)、边界层高度(BLH, m)和相对湿度(RH, %).其中, 10 U与10 V两个参数组合在一起可给出水平10 m风的速度和方向.
1.2.3 PM2.5地面监测数据我国从2013年1月陆续开始发布PM2.5地面监测数据, 本研究从全国城市空气质量实时发布平台[28]获取了2013年1月至2017年12月内的62 166 756行每小时实测数据, 共涉及1 621个国控空气质量地面监测站点, 其中, PM2.5地面监测站点分布如图 1所示, 监测站点的数量呈现明显的由东向西、由南向北逐渐递减的分布特点.
2 数据处理 2.1 二阶段PM2.5浓度组合估算模型大部分研究结果表明, 回归型支持向量机[29]和BP神经网络[29, 30]的非线性回归效果普遍优于传统的统计回归模型.BP神经网络对初始值(初始权重和阈值)设置较敏感, 通过思维进化算法(mind evolutionary computation, MEC)对其进行优化, 可有效避免其陷入局部最小[30].MEC-BP算法不受初始种群规模限制, 学习速度快, 既能满足精度要求, 也能降低网络复杂性, 具有很好地扩展性.通常在解决小样本的非线性回归问题时, ε-SVR具有许多独特优势, 但当其解决大数据集的非线性回归问题时处理速度较慢, 而BP神经网络的收敛速度较快.本研究数据样本较大, 气象参数、AOD与PM2.5之间的关系是非线性的, AOD与PM2.5关系存在明显的季节性和区域性差异性.因此, 本研究提出了一种基于ε-SVR-MEC-BP算法的近地面PM2.5浓度二阶段组合估算模型, 组合模型的训练和估算过程如下.
2.1.1 数据预处理(1) MODIS AOD提取、融合及插值 因从MOD08和MYD08文件提取的两种AOD数据的时间覆盖范围和数据缺失情况有差异, 故使用一元线性回归模型以实现MODIS AOD数据的融合[7, 9, 20], 当两种AOD数据中仅有一个数据缺失时, 可以结合线性回归模型估算另一个数据, 最终两种AOD数据的均值参与模型构建.由于云层覆盖、高表面反射率(如沙漠、积雪天气、冰川覆盖)或检索错误等原因[14], 几乎每个月的AOD数据都存在不同程度的非随机性缺失, 故使用克里金插值对AOD数据缺失率较低的春、夏和秋季进行插值, 但冬季缺失率较高, 所以不进行插值处理.经过插值处理后的AOD数据不参与模型训练, 仅在获得理想模型后估算PM2.5时使用.
(2) 气象数据和PM2.5实测数据的提取 每年的气象数据以一个三维数组的形式存储([经度][纬度][月份]), 批量提取每一种气象数据.本研究获取到的地面PM2. 5监测数据为每小时的点值, 为了确保由点到面转化后的PM2.5浓度的月均值能够更准确地表征每一个网格内的近地面PM2.5浓度, 在提取数据过程中引入了《环境空气质量标准》(GB 3095-2012):首先, 剔除PM2.5每小时实测数据的缺失值、异常值(PM2.5≤0), 以网格为单位计算每小时PM2.5浓度均值, 然后计算日均PM2.5浓度(每日多于20 h), 最后计算PM2.5浓度月均值(2月多于25 d, 其余月份多于27 d).
(3) 划分数据集 删除没有覆盖在我国陆地区域的网格(剩余1 112个)及其所对应的AOD数据和气象数据, 并完成气象数据、MODIS AOD数据和PM2.5数据的时空匹配, 得到集合S.将数据集S分别按季节和区域划分为子数据集(Sz1, Sz2, …, Sz8、St1, St2, …, St4, 其中, Sz表示按区域划分, St表示按季节划分), 其中八大经济区域划分主要是针对内陆地区31个省(不包括港澳台), 因港澳台地理位置、自然条件与南部沿海地区较接近, 故将港澳台与南部沿海地区划为一类.
2.1.2 模型构建及训练(1) 第一阶段ε-SVR分区域、分季节建模
将子数据集(Sz1, Sz2, …, Sz8、St1, St2, …, St4)内的PM2.5、MODIS AOD和气象参数(10U、10V、BLH、RH)均非空的数据用于建模.基于每一个子数据集分别训练12个ε-SVR模型, 选取了高斯径向基核函数, 并使用交叉检验用于确定ε-SVR参数.在多次实验获得理想模型后, 将12个PM2.5浓度估算序列合并, 最后得到由分区域估算的PM2.5、分季节估算的PM2.5和实际PM2.5浓度序列组成的数据集S1.
(2) 第二阶段MEC-BP神经网络建模
数据集S1内分区域估算的PM2.5和分季节估算的PM2.5作为输入变量, 实际的PM2.5浓度作为输出变量, 使用MEC-BP进行模型训练, 其中输入层节点数量为2, 输出节点数量为1, 通过多次实验得到当隐含层节点数为5时模型精度最高, 故BP神经网络的拓扑结构为2-5-1.MEC算法初始化参数如下:种群大小为200, 子种群大小为20, 优胜子种群个数为5, 临时子种群个数为5, 迭代次数为100.
2.1.3 模型精度评价及数据估算本研究使用PM2.5浓度的估算值与观测值的均方误差(MSE, μg·m-3)、绝对误差(MAE, μg·m-3)、平均绝对百分比误差(MAPE, %)、均方百分比误差(MSPE, %)、均方根误差(RMSE, μg·m-3)和决定系数R2评价模型的性能.得到理想的模型后, 将AOD和气象参数均为非空并且实测的PM2.5浓度数据为空的数据集用作输入变量以估算在我国2000-03~2017-01期间内缺失的PM2.5浓度数据.
2.2 PM2.5浓度时间序列的趋势和突变检验基于实测的和估算的PM2.5浓度数据, 将突变检验和趋势分析相结合, 以揭示近20年来我国PM2.5浓度的阶段性动态变化特征.首先, 使用Pettitt突变检验[31](Pettitt test)对PM2.5浓度时间序列做突变点分析, 如果检测到突变点, 则以突变点为分割点, 将该PM2.5时间序列分割为两个序列, 并对两个子序列分别进行Mann-Kendall(M-K)趋势检验[31].同时, 本研究也使用M-K趋势检验对完整的PM2.5浓度时间序列进行趋势检验, 以揭示研究时间范围内PM2.5浓度的整体变化特征.此外, 还计算了多个PM2.5浓度时间序列的Theil-Sen's[32]倾斜度以量化变化趋势的大小.
2.2.1 Pettitt突变检验Pettitt突变检验是一种基于秩的检验, 能检测到时间序列的平均值的突然变化.通过突变检验可找出时间序列的突变时间t0, 并计算出统计变量P, 若P小于刚给定显著性水平时(取α=0.1), 则认为检测出的突变点在统计学意义上是显著的, 已知长度为n的按年份顺序排列的PM2.5年均浓度序列X, 具体的计算公式如下:
|
(1) |
|
(2) |
若在t0年满足k(t0)=max|Sk|, k=2, 3, …, n, 则t0年即为突变点.
|
(3) |
Mann-Kendall (M-K)趋势检验是一种非参数检验, 通常用于监测环境数据、气候数据或水文数据的单调趋势, 它不要求样本遵从一定的分布, 也不受少数异常值的干扰, 它已被普遍用于环境方面时间序列的趋势检测.通过Mann-Kendall趋势检验可以检测到某个时间序列中是否存在单调的增加或减小的趋势, 以及检测到的趋势是否具有统计学意义, 即通过Mann-Kendall趋势检验可以计算出指定时间序列的测试统计量Z, 正的Z表示上升趋势, 负的Z表示下降趋势, 其中统计量Z遵循正正态分布, 当|Z|>U1-α/2时表示变化趋势具有统计显著性, 其中U1-α/2是标准正态分布的(1-α/2)分位数.
已知长度为n的按年份顺序排列的PM2.5年均浓度序列X, Mann-Kendall趋势检验的过程即先将每年的年均PM2.5 Xi(i=1, 2, …, n)作为参考点将其与其余数据点Xj (j=i+1, 2, …, n)进行比较并计算统计量Z, 具体计算公式如下:
|
(4) |
|
(5) |
|
(6) |
|
(7) |
通过Theil-Sen's趋势斜率估计可以计算出某个时间序列的趋势斜率, 即对趋势幅度进行的稳健估计, 该方法已经被广泛用于识别时间序列中趋势线的斜率.通过计算变量β(倾斜度)即可知道序列随时间增加或减小的具体值[32].
|
(8) |
式中, xi和xj分别表示在ti和tj年的PM2.5浓度(i>j).
3 结果与分析 3.1 子数据集的描述性统计基于分区域和季节划分的12个子数据集Szi(i=1, 2, …, 8)和Stj(j=1, 2, …, 4)分别计算了各子集内PM2.5、MODIS AOD、T2M、10U、10V、BLH和RH这7个变量之间的Pearson相关系数和各变量的概率分布.结果表明不同的区域和季节, PM2.5与MODIS AOD的概率分布是相似的, 且MODIS AOD与PM2.5的相关系数均大于0且通过了95%的显著性水平检验即两者之间正相关；因为季节、气象条件和地理位置存在差异, PM2.5浓度的概率分布, MODIS AOD与PM2.5之间的关系均呈现明显差异, 该结果也进一步论证了本研究分区域、分季节建模的合理性.
分区域角度, MODIS AOD与PM2.5的相关系数在0.24~0.56范围内, 其中南部沿海(SC) MODIS AOD与PM2.5的相关系数最低(r=0.24)这可能因为南部沿海地区临海面积大, 海盐颗粒物为该区域气溶胶的主要组成部分, PM2.5的占比较低；同时, 沿海地区的下浮地表特征和复杂的气溶胶类型会对中等分辨率光谱仪的算法精度产生影响, 进而加大了MODIS AOD的反演难度[33].在地形复杂、有冰川分布且沙漠分布较密集的西北地区和冬季积雪较多的东北地区AOD与PM2.5的相关系数也较低(r < 0.4)；在西南地区(SW)MODIS AOD与PM2.5浓度的相关性最高(r=0.56)；此外, 气象参数(T2M、10U、10V、BLH、RH)与PM2.5的相关系数大部分为负数, 这说明温度、风速、边界层高度和相对湿度对PM2.5的增长有一定的抑制作用.而在分季节角度, 夏季里MODIS AOD与PM2.5的相关性最高(r=0.67), 冬季相关性最低(r=0.53), 在春季和秋季相关性较高.
3.2 模型的训练和验证经过第一阶段ε-SVR分区域建模后, 各经济区域模型效果呈由沿海向内陆逐渐递减的趋势；在监测站点分布较密集的东部沿海、北部沿海、南部沿海、黄河中游和长江中游地区模型拟合效果较好；在大西北地区模型的决定系数最低且误差类指标较差, 主要原因是大西北地区PM2.5监测点分布较稀疏、地形复杂、沙漠、冰川分布较密集.经过第一阶段ε-SVR分季节建模后, 模型拟合效果在春季最差, 在夏、秋和冬季模型拟合效果较好.其中, 在MODIS AOD数据缺失率较高的冬季和东北区域模型拟合效果较好的原因主要是:缺失的MODIS AOD数据剔除了雨雪等极端天气的AOD数据, 而这部分数据恰好是与地面PM2.5浓度相关性较弱的AOD数据, 进而提升了模型精度.表 1为将样本集S1内分区域估算和分季节估算的PM2.5序列作为模型的输入变量, 实测的PM2.5序列作为模型的输出变量, 经过第二阶段建模后的完整模型的效果, 以及单独使用ELM(极限学习机)、ε-SVR、MEC-BP对整个未划分的数据集S进行模型拟合的效果.分析可知, 相比只使用ELM、ε-SVR或MEC-BP, 通过本研究提出的组合模型进行建模获得了最好的PM2.5浓度估算效果.
|
表 1 模型的具体性能评价参数及对比 Table 1 Comparison of specific performance and specific performance evaluation parameters |
为了比较经过二阶段MEC-BP建模后在各区域、各季节角度模型估算效果的改善情况, 将组合模型最终估算结果序列和实测PM2.5序列组成的数据集, 按照区域和季节进行划分, 并计算了划分后子集的各项模型性能评价指标.分别以不同的区域和季节计算各项评价指标差值(完整模型的各项评价指标分别减去第一阶段ε-SVR分区域、分季节建模的各项评价指标)情况(图 2), 可以清楚地看到, 经过第二阶段MEC-BP模型的进一步校正, 从各经济区域、各个季节和全国角度看, 结合PM2.5浓度的估算序列和实测序列计算的决定系数(ΔR2＞0)和大部分误差指标都有明显改善(图 2中绿色点数明显多于红色点).
|
①完整模型效果与一阶段ε-SVR分区域建模合并后对比；②完整模型效果与一阶段ε-SVR分季节建模合并后对比 图 2 完整模型与一阶段ε-SVR分区域、分季节建模效果对比 Fig. 2 Comparison of the performance between the complete model and the first stage of the regional and seasonal ε-SVR model |
在完成模型精度评价及验证后, 将估算的历史月均PM2.5浓度数据与实测的PM2.5浓度数据合并, 进一步从多尺度角度分析了全国PM2.5浓度的时空分布及时间变化特征.其中, 当计算2000年的年均PM2.5浓度时, 选用了2001年1月和2月的数据填补了2000年1~2月每一个格网内每一个格网内的PM2.5数据.
3.3 2000~2017年中国PM2.5浓度的时空分布特征 3.3.1 全国尺度图 3为全国2000~2017内每年的年均PM2.5浓度的冷热点分析结果叠加(置信水平为90%、95%和99%).近20年内中国PM2.5年均浓度总体呈现东高西低的空间分布特点.在北部沿海地区、东部沿海地区和长江中游地区观察到了大范围高浓度的PM2.5；最高的PM2.5浓度普遍分布在北部沿海地区的河北省和天津市, 同时在北京以北区域PM2.5污染稍轻；相比北部沿海和东部沿海, 南部沿海区域PM2.5浓度最低；较低的PM2.5浓度主要在西南地区以西, 大兴安岭地区、东北的小兴安岭地区, 西藏以东地区被观测到；在西南地区, PM2.5浓度在重庆市和四川以东的区域最高, 在广西区和贵州省较高, 在云南省和四川以西区域最低；东北地区的PM2.5浓度自西向东逐渐递减.
|
图 3 2000~2017年PM2.5年均值的冷热点分析 Fig. 3 Analysis of cold and hot points of PM2.5 annual means from 2000 to 2017 |
年均PM2.5浓度的高值(热点)主要集中于在人口密集、经济发达、工业化水平较高(如北部沿海、东部沿海、长江中游地区)以及戈壁沙漠分布密集的区域(如大西北地的新疆塔克拉玛干沙漠一带以及黄河中游的陕西、山西、甘肃和内蒙古以西区域).年均PM2.5浓度的低值(冷点)主要集中于四川、云南和西藏接壤处、以及大兴安岭地区与东北地区的小兴安岭地区.
3.3.2 八大经济区域尺度图 4为2000-03~2017-12范围内每月和每个季节各经济区域PM2.5变化曲线.分析可知大部分经济区域的月均PM2.5浓度在1月达到最大值, 在1~6月快速下降, 到6月达到最低值, 然后在6~12月快递增大；部分地区偶尔会呈现高峰值在12月到次年1月之间的双峰型；只有在大西北地区每年的月均PM2.5浓度的波动较稳定, 而在其余区域每年的月均PM2.5波动都不太一致.从季节均值上看, 各经济区域内的冬季PM2.5浓度最高, 冬季到夏季逐渐递减, 到夏季达到最低值, 然后夏季到冬季逐渐递增；在不同的地区和不同的年份, 春季和秋季的浓度大小存在明显差别.
|
图 4 2000-03~2017-12各区域的月均及季均PM2.5浓度变化曲线 Fig. 4 Variation curve of monthly and seasonal average of PM2.5 concentrations in eight economic regions from March 2000 to December 2017 |
根据《环境空气质量标准》(GB 3095-2012)规定, PM2.5年均值一级浓度限值为15 μg·m-3, 二级浓度限值为35 μg·m-3, 分析全国PM2.5年均浓度变化曲线以及每年PM2.5均值在≤15、15~35和>35 μg·m-3这3个范围内的1°×1°格网数量占全国格网数量(共1 112个)的占比变化情况可知, 2000~2017年内我国近65%的区域PM2.5浓度超过二级浓度限值, PM2.5污染形势严重(图 5).
|
图 5 2000~2017年间PM2.5浓度变化特征 Fig. 5 Changing characteristics of PM2.5 concentrations from 2000 to 2017 |
对2000~2017年内全国年均PM2.5序列进行M-K趋势检验的结果未通过95%置信水平的显著性检验, 呈不显著的增长趋势(Theil-Sen's倾斜度为0.011).通过Pettitt突变检验检测到2014年为研究时间范围内全国年均PM2.5序列的显著突变时间, 以2014年为分界线, 将全国年均PM2.5序列分为两个子序列并分别对其进行M-K趋势检验, 两个子序列的变化趋势均通过了95%置信水平的显著性检验, 2000~2014年中国年均PM2.5浓度缓慢持续增长(Theil-Sen's倾斜度为0.124), 2014~2017年中国年均PM2.5浓度快速减小(Theil-Sen's倾斜度为-1.757).检验结果如表 2所示.
|
|
表 2 PM2.5年均序列的突变检验、M-K趋势检验及Theil-Sen's倾斜度结果1) Table 2 Change-point test, M-K trend test, and Theil-Sen's slope result of PM2.5 concentrations |
3.4.2 八大经济区域尺度
图 6为2000~2017年各经济区域的年均PM2.5浓度变化, 分析结果可知, 在北部沿海、东部沿海、长江中游和黄河中游的PM2.5污染较严重, 区域年均PM2.5浓度普遍高于45 μg·m-3(超过了PM2.5二级浓度限值35 μg·m-3)；八大经济区域中, PM2.5浓度较高的区域主要为北部沿海、其次是东部沿海和长江中游, 在南部沿海和西南地区PM2.5浓度较低.在黄河中游地区和大西北地区, 沙尘暴(风扬起的灰尘和沙粒)为PM2.5的主要来源之一, 该区域沙漠、沙地分布较密集, 位于我国西北方的干旱和半干旱地区, 是亚洲尘埃的主要来源区之一, 尤其是在春季和冬季受当地沙尘暴影响较大, 导致空气中包括PM2.5在内的颗粒物浓度急剧增加, 伴随着强烈的西伯利亚/蒙古冷锋, 粒径较小的颗粒物会沿着东南向轨道移出戈壁和沙漠地区, 途经北部沿海地区的北京、天津和河北, 导致途经区域的PM2.5浓度的增加.此外, 发达的经济和密集的人口及密集的人类活动和快速的城市化也间接导致了北部沿海地区存在较高浓度的PM2.5.
|
图 6 各区域PM2.5年均值突变检验结果 Fig. 6 Results of annual mean change-point test of PM2.5 concentrations in each region |
使用Pettitt test分别对各经济区域的PM2.5年均序列进行突变检验(表 2), 结果表明, 在北部沿海、东部沿海、长江中游、东北地区和黄河中游存在显著的突变点；其中东部沿海和长江中游的突变时间为2014年, 在北部沿海和东北地区突变时间为2015年, 在黄河中游地区突变时间为2017年；大西北地区、南部沿海和西南地区检测到的突变点未通过显著性检验, 可能因为大西北地区分布有较多的沙漠和戈壁的沙尘, 植被覆盖较少, 该区域的空气污染治理是一个长期持续的过程, 而在西南地区和南部沿海地区空气条件较好所以年均PM2.5时间序列变化幅度较小, 所以这3个区域内很难存在明显的突变点.对2000~2017年内各区域的PM2.5年均序列(表 2中①)进行M-K趋势检验分析的结果表明, 只有黄河中游的年均PM2.5浓度序列呈显著下降趋势(倾斜度为-0.166), 其余7个区域的年均PM2.5浓度均的趋势变化均未通过95%置信水平的显著性检验, 变化趋势不明显.为了更好地分析各区域PM2.5浓度阶段变化特点, 以Pettitt test检测到的突变点为分割点, 将其划分为两段子序列(表 2中突变点前②和突变点后③)并分别对每个子序列其进行M-K趋势检验和Theil-Sen's倾斜度估计.
分析M-K趋势检验结果可知, 我国大部分经济区域的PM2.5年均浓度呈现明显的阶段变化, 突变点之前的子序列呈现显著的增长趋势, 突变点之后的序列呈现显著的减小趋势, 虽然西南地区和南部沿海地区的PM2.5年均浓度序列的突变点未通过显著性检验, 但是依然表现出明显的阶段性变化特征.分析Theil-Sen's倾斜度估计结果可知, 除黄河中游外, 其余7个区域在突变点之前的年均PM2.5浓度子序列的均呈显著增长趋势, 其倾斜度按东部沿海、北部沿海、长江中游、南部沿海、西南地区、东北地区和大西北的顺序递减；突变点后的年均PM2.5浓度子序列均呈显著下降趋势, 倾斜度的大小按北部沿海、东部沿海、长江中游、东北地区、南部沿海、西南地区和西北地区NW的顺序递减.各区域在突变点之后的年均PM2.5浓度序列的倾斜度明显普遍高于突变之前的子序列, 从2014年开始全国PM2.5污染情况明显改善.
3.4.3 1°×1°的格网尺度图 7(a)为每一个格网内的年均PM2.5浓度.全国范围PM2.5年均浓度呈显著下降趋势和显著增长趋势的格网分布较离散, 呈显著增长趋势的格网主要集中在北部沿海地区；而PM2.5年均浓度呈显著下降趋势的格网主要集中在黄河中游地区, 较少格网的Theil-Sen's倾斜度超过1；同时, 从格网数来看, 呈显著下降趋势的格网少于呈显著增长趋势的格网数.结合全国尺度和八大经济区域尺度的PM2.5浓度年均序列的Pettitt突变检验、M-K趋势检验和Theil-Sen's倾斜度估计分析结果可知, 2014年是全国以及PM2.5污染较严重的几个地区的PM2.5年均浓度序列从持续缓慢上升到快速下降的关键转折点.
|
图 7 PM2.5年均值的M-K趋势检验结果及其趋势斜率的估计 Fig. 7 M-K trend test result and Theil-Sen slope estimator of the annual PM2.5 concentrations |
图 7(b)为2000~2014年内每一个格网的M-K趋势检验结果.在该阶段内有近24.37%格网的年均PM2.5浓度呈显著增长趋势(Theil-Sen's倾斜度为0.1~1.9), 近7.28%格网的PM2.5年均PM2.5浓度呈显著下降趋势(Theil-Sen's倾斜度为-0.1~-1.4), 其余格网的PM2.5浓度变化趋势不显著；PM2.5显著增长的格网主要集中在北部沿海的山东省、长江中游、东部沿海、黄河中游地区的陕西省、东北地区的辽宁省等；PM2.5年均PM2.5浓度显著下降的格网主要分布在大西北地区和黄河中游的内蒙古地区.
图 7(c)为2014~2017年内每一个格网的M-K趋势检验结果.在该阶段内近17.8%的格网的PM2.5年均浓度呈显著下降趋势(Theil-Sen's倾斜度为-0.1~-16.5), 仅有近1.7%的格网PM2.5年均浓度呈显著增长趋势(Theil-Sen's倾斜度为0.1~6.7)；表明从2014年开始全国范围内的PM2.5污染状况得到快速有效的改善.
4 讨论本研究构建的二阶段PM2.5浓度组合估算模型综合考虑了AOD-PM2.5之间的关系存在着非线性、季节差异、区域差异的特点, 并结合气象参数(气温、风速、风向、边界层高度、相对湿度), 将ε-SVR与MEC-BP神经网络相组合, 提供了中国2000~2017年内分辨率为1°×1°的可靠的PM2.5月度地面浓度数据估算, 相比现有研究具有更广的时空覆盖度.在实测PM2.5浓度数据由点到面的提取过程中, 引入了《环境空气质量标准》(GB 3095-2012), 较大程度地增加了估算结果的准确度.本研究成果可为其它的地区估算大范围内的历史PM2.5浓度提供方法借鉴.
本研究构建的组合模型的决定系数R2为0.838, 均方根误差RMSE为11.512 μg·m-3, 平均绝对百分比误差MAPE为14.905%, 均方百分比误差MSPE为0.243%, 绝对误差MAE为6.476 μg·m-3, 均方误差MSE为132.519 μg·m-3, 具有误差小, 泛化高的综合估算能力.地面监测站空间分布的不均匀可能会导致PM2.5估算结果存在一定程度的误差, 大部分PM2.5监测站点都位于PM2.5污染较严重的城市地区, 如在大西北(NW)的农村地区基本无PM2.5地面监测站点分布；同时, 冬季在东北的黑龙江、吉林和黄河中游地区的大兴安岭有大面积积雪覆盖、在大西北地区的西北边缘和新疆与西南地区接壤处的横断山脉一带有大规模的冰川分布等, 这些原因都可能会导致卫星采样偏差, 进而导致MODIS AOD产生误差或缺失；在新疆和内蒙古境内分布有大面积的沙漠和沙地也可能导致卫星遥感反演的MODIS AOD值偏高, 进而导致估算的PM2.5浓度偏高.基于估算数据和实测数据, 本研究分析了2000~2017年中国PM2.5浓度的时空分布特征, 在全国和区域尺度中国PM2.5浓度的月均值从1~12月呈明显的”U”型变化特征, 从6月到次年6月呈明显的倒”U”型变化特征.中国的PM2.5污染在冬季最严重, 在夏季最轻, 在秋季和春季较轻.高浓度的PM2.5主要分布在人口密集、经济发达、工业化水平较高(如北部沿海、东部沿海、长江中游)以及戈壁沙漠密集的区域(如大西北地区、黄河中游地区), 较低浓度的PM2.5主要分布在南部沿海和西南地区, 该研究结果与大部分已有研究结果类似.综上, 本研究构建的二阶段组合估算模型的PM2.5浓度估算结果能够较好地表征2000-03~2017-12内中国月均PM2.5的时空变化过程.
此外, 本研究还利用数理统计与地理空间分析手段, 从多空间尺度揭示了以PM2.5为主要代表的中国空气质量在2000~2017年内的一个退化及改善的阶段性动态变化过程.
5 结论(1) 本研究提出的二阶段组合估算模型提供了中国2000~2017年内空间分辨率为1°×1°的月度近地面PM2.5浓度的可靠估算, 有效地弥补了中国地面监测网络在时间和空间上的空白.
(2) 本研究时间范围内, 全国尺度PM2.5年均浓度序列的整体变化趋势不显著, 但是呈明显的阶段变化特征:2000~2014年全国年均PM2.5浓度缓慢持续增长, 2014~2017年全国年均PM2.5浓度快速减小, 其中从2014年开始, 空气污染较严重的北部沿海、东部沿海和长江中游地区的PM2.5污染情况改善更显著.总体而言, 2014年是2000~2017年内中国PM2.5浓度从持续缓慢上升到快速下降的关键转折点.
(3) 本研究时间范围内全国有超过65%的区域PM2.5年均浓度超过了二级限值(35 μg·m-3), 虽然全国PM2.5污染情况有一定程度的改善, 但是空气污染形势依然严峻.
| [1] | Wang S J, Fang C L, Guan X L, et al. Urbanisation, energy consumption, and carbon dioxide emissions in China:a panel data analysis of China's provinces[J]. Applied Energy, 2014, 136: 738-749. DOI:10.1016/j.apenergy.2014.09.059 |
| [2] | Guan X L, Wei H K, Lu S S, et al. Assessment on the urbanization strategy in China:achievements, challenges and reflections[J]. Habitat International, 2018, 71: 97-109. DOI:10.1016/j.habitatint.2017.11.009 |
| [3] |
姜磊, 周海峰, 柏玲, 等. 空气质量指数(AQI)的社会经济影响因素分析——基于指数衰减效应视角[J]. 环境科学学报, 2018, 38(1): 390-398. Jiang L, Zhou H F, Bai L, et al. The analysis of socio-economic factors of air quality index (AQI) based on the perspective of the exponential decay effects[J]. Acta Scientiae Circumstantiae, 2018, 38(1): 390-398. |
| [4] |
郦嘉诚, 高庆先, 李亮, 等. 对首要污染物所揭示的京津冀环境空气质量状况的认识启迪与对策建议[J]. 环境科学研究, 2018, 31(10): 1651-1661. Li J C, Gao Q X, Li L, et al. Enlightenment and suggestions on the air quality of Beijing, Tianjin and Hebei revealed by primary pollutants[J]. Research of Environmental Sciences, 2018, 31(10): 1651-1661. |
| [5] | Kahn R, Banerjee P, McDonald D, et al. Sensitivity of multiangle imaging to aerosol optical depth and to pure-particle size distribution and composition over ocean[J]. Journal of Geophysical Research:Atmospheres, 1998, 103(D24): 32195-32213. DOI:10.1029/98JD01752 |
| [6] | Kahn R, Banerjee P, McDonald D. Sensitivity of multiangle imaging to natural mixtures of aerosols over ocean[J]. Journal of Geophysical Research:Atmospheres, 2001, 106(D16): 18219-18238. DOI:10.1029/2000JD900497 |
| [7] | 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. |
| [8] | Wang J, Christopher S A. Intercomparison between satellite-derived aerosol optical thickness and PM2.5 mass:Implications for air quality studies[J]. Geophysical Research Letters, 2003, 30(21): 2095. DOI:10.1029/2003GL018174 |
| [9] |
赵仕伟, 高晓清. 利用MODIS C6数据分析中国西北地区气溶胶光学厚度时空变化特征[J]. 环境科学, 2017, 38(7): 2637-2646. Zhao S W, Gao X Q. Analysis of spatio-temporal distribution and variation characteristics of aerosol optical depth over the northwest of China by MODIS C6 product[J]. Environmental Science, 2017, 38(7): 2637-2646. |
| [10] | 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): D20205. DOI:10.1029/2008jd011497 |
| [11] | You W, Zang Z L, Pan X B, et al. Estimating PM2.5 in Xi'an, China using aerosol optical depth:a comparison between the MODIS and MISR retrieval models[J]. Science of the Total Environment, 2015, 505: 1156-1165. DOI:10.1016/j.scitotenv.2014.11.024 |
| [12] |
郭建平, 吴业荣, 张小曳, 等. BP网络框架下MODIS气溶胶光学厚度产品估算中国东部PM2.5[J]. 环境科学, 2013, 34(3): 817-825. Guo J P, Wu Y R, Zhang X Y, et al. Estimation of PM2.5 over eastern China from MODIS aerosol optical depth using the back propagation neural network[J]. Environmental Science, 2013, 34(3): 817-825. |
| [13] | 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 |
| [14] | Liu Y, Paciorek C J, Koutrakis P. Estimating regional spatial and temporal variability of PM2.5 concentrations using satellite data, meteorology, and land use information[J]. Environmental Health Perspectives, 2009, 117(6): 886-892. DOI:10.1289/ehp.0800123 |
| [15] | Xiao Q Y, Wang Y J, Chang H H, et al. Full-coverage high-resolution daily PM2.5 estimation using MAIAC AOD in the Yangtze River Delta of China[J]. Remote Sensing of Environment, 2017, 199: 437-446. DOI:10.1016/j.rse.2017.07.023 |
| [16] | Nguyen T N T, Ta V C, Le T H, et al. Particulate matter concentration estimation from satellite aerosol and meteorological parameters: data-driven approaches[A]. In: Huynh V N, Denoeux T, Tran D H, et al (Eds.). Knowledge and Systems Engineering: Proceedings of the Fifth International Conference KSE 2013, Volume 1[M]. Cham: Springer, 2014. 351-362. |
| [17] | 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 |
| [18] | Guo J P, Xia F, Zhang Y, et al. Impact of diurnal variability and meteorological factors on the PM2.5-AOD relationship:Implications for PM2.5 remote sensing[J]. Environmental Pollution, 2017, 221: 94-104. DOI:10.1016/j.envpol.2016.11.043 |
| [19] | 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 |
| [20] | Liu M, Bi J, Ma Z. Visibility-based PM2.5 concentrations in China:1957-1964 and 1973-2014[J]. Environmental Science & Technology, 2017, 51(22): 13161-13169. |
| [21] | 国务院.大气污染防治行动计划[EB/OL]. http://www.jingbian.gov.cn/gk/zfwj/gwywj/41211.htm?from=timeline&tdsourcetag=s_pctim_aiomsg, 2013-09-10. |
| [22] |
罗毅, 邓琼飞, 杨昆, 等. 近20年来中国典型区域PM2.5时空演变过程[J]. 环境科学, 2018, 39(7): 3003-3013. Luo Y, Deng Q F, Yang K, et al. Spatial-temporal change evolution of PM2.5 in typical regions of China in recent 20 years[J]. Environmental Science, 2018, 39(7): 3003-3013. |
| [23] |
杨昆, 喻臻钰, 罗毅, 等. 湖泊表面水温预测与可视化方法研究[J]. 仪器仪表学报, 2017, 38(12): 3090-3099. Yang K, Yu Z Y, Luo Y, et al. Lake surface water temperature prediction and visualization[J]. Chinese Journal of Scientific Instrument, 2017, 38(2): 3090-3099. |
| [24] |
杨昆, 罗毅, 徐玉妃, 等. 基于无线传感器网络与GIS的蓝藻水华爆发动态监测与模拟[J]. 农业工程学报, 2016, 32(24): 197-205. Yang K, Luo Y, Xu Y F, et al. Dynamic monitoring and simulation of Cyanobacteria Bloom based on wireless sensor network and GIS[J]. Transactions of the Chinese Society of Agricultural Engineering, 2016, 32(24): 197-205. DOI:10.11975/j.issn.1002-6819.2016.24.026 |
| [25] | 国家统计局. 2000-2017年地区分省年度数据地区生产总值[EB/OL]. http://data.stats.gov.cn/easyquery.htm?cn=E0102. |
| [26] | 大气存档和分配系统. MOD08 M3和MYD08 M3文件[EB/OL]. https://ladsweb.modaps.eosdis.nasa.gov/search/. |
| [27] | 欧洲中期天气预报中心.气象数据[EB/OL]. https://apps.ecmwf.int/datasets/. |
| [28] | 全国城市空气质量实时发布平台. PM2.5地面监测数据[EB/OL]. http://106.37.208.233:20035. |
| [29] | Kavousi-Fard A, Samet H, Marzbani F. A new hybrid modified firefly algorithm and support vector regression model for accurate short term load forecasting[J]. Expert Systems with Applications, 2014, 41(13): 6047-6056. DOI:10.1016/j.eswa.2014.03.053 |
| [30] | Zhang J, Yang Y P, Zhang J M. A MEC-BP-Adaboost neural network-based color correction algorithm for color image acquisition equipments[J]. Optik, 2016, 127(2): 776-780. DOI:10.1016/j.ijleo.2015.10.120 |
| [31] |
强皓凡, 靳晓言, 刘超, 等. 基于水热耦合平衡的长江源实际蒸散变化研究[J]. 干旱区资源与环境, 2018, 32(3): 106-111. Qiang H F, Jin X Y, Liu C, et al. Variations of actual evapotranspiration in the Yangtze River headwater region based on coupled water-energy balance[J]. Journal of Arid Land Resources and Environment, 2018, 32(3): 106-111. |
| [32] | Some'e B S, Ezani A, Tabari H. Spatiotemporal trends and change point of precipitation in Iran[J]. Atmospheric Research, 2012, 113: 1-12. DOI:10.1016/j.atmosres.2012.04.016 |
| [33] | Lue Y L, Liu L Y, Hu X, et al. Characteristics and provenance of dustfall during an unusual floating dust event[J]. Atmospheric Environment, 2010, 44(29): 3477-3484. DOI:10.1016/j.atmosenv.2010.06.027 |