冀中南平原区近地面PM2.5遥感反演

王诗瑶, 李夫星,2, 王 晨, 赵旭雅, 玛伊热·买买提

(1.河北师范大学 地理科学学院,河北 石家庄 050024; 2.河北省环境变化遥感技术识别创新中心,河北 石家庄 050024)

大气细颗粒物PM2.5严重影响环境空气质量并危及人类身体健康,已成为衡量全球可持续发展目标(SDGs)的一项重要指标[1-3].很多流行性病学研究表明,PM2.5暴露与人体呼吸系统及心血管系统疾病的发生率及死亡率密切相关[4-6].因此,监测PM2.5质量浓度时空分布及变化特征尤为重要.地面监测方法虽然能够获取准确的浓度信息,但是难以开展大尺度空间范围内的面状监测[7].大气数值模型虽然能够给出大尺度网格化PM2.5质量浓度模拟数据,但是对近地面颗粒物理化性质及污染物排放清单要求较高,且参数众多,严重影响其可操作性[8-9].

卫星遥感反演是一种基于PM2.5质量浓度与AOD产品的相关性及物理属性的相似性,并以此构建相关性模型来模拟PM2.5质量浓度的方法.该方法克服了地面监测及数值模拟分别在空间连续性及资料齐全性等方面的不足,被广泛应用于近地面大气颗粒物浓度反演研究[10-12].根据核心方法的差异性,卫星遥感反演方法可以分为随机统计模型和物理机理性模型.随机性统计模型在AOD-PM2.5关系中增加了随机效应因子的设计,如线性混合效应模型(linear mixed effect model,LMEM)[13],能够有效捕捉AOD-PM2.5关系的时间异质性,提高了模型预测精度.近年来,机器学习模型也被广泛应用于近地面PM2.5反演,如随机森林[14]、深度学习[15]和极端梯度提升[16]等,这些模型不仅显著提升了预测精度,而且大幅度提高了AOD有效值覆盖度.然而,现有高级统计模型及机器学习模型均需要大量的地基监测数据用以拟合模型参数并验证其可靠性.由于地基监测站点在乡村和山区分布极少,主要集中在大中都市地区,因此,模型训练样本不可能全面系统地表征所有场景,其典型性和代表性对模型反演精度具有重要影响.很多研究表明,当验证数据来自相同时期的训练样本时,数据验证结果较好,反之则比较差[17].由于气溶胶与大气细颗粒物具有相近物理机制,学者们利用直接相关模型探索了基于卫星遥感反演的AOD估算近地面细颗粒物浓度的可能性,并在世界各地取得了一定效果,如美国[12]、加拿大[18]、欧洲[19]和中国[20].然而卫星遥感反演AOD为整层大气柱中气溶胶消光系数在垂直方向上的积分,并不能完全表征近地表大气细颗粒物的消光能力.此外,近地面气象要素对细颗粒物理化性质也会产生影响,如相对湿度[21].利用直接相关模型估算高精度的地表大气细颗粒物浓度普遍较难[3].因此,学者们通过对AOD进行垂直和湿度订正,发展了基于物理机制的经验反演模型,有效提高了模型预测地表细颗粒物浓度的精度[22].但云覆盖等因素的影响会导致AOD产品时空覆盖度较低,大大影响近地面PM2.5质量浓度反演精度.因此,如何提高AOD产品有效覆盖度是亟待解决的问题[23].虽然垂直校正及湿度校正等中间参数的估算会引入误差,并最终影响到地表大气细颗粒物浓度的反演精度,但基于气溶胶物理动力机制的确定性经验反演模型在优化反演思路、模型及提高反演精度上仍然有广泛的改进空间,而且在颗粒物浓度历史重建及时间预测方面仍具有十分重要的优势.

本文中,笔者以冀中南地区作为研究区(见图1),基于MAIAC AOD等遥感数据产品,在结合地面PM2.5、气象数据、AERONET AOD、土地利用等观测数据的基础上,通过修复云覆盖等因素导致的AOD数据缺失,改进物理机理性模型,发展适用于冀中南平原区地面PM2.5遥感反演模型.该研究不仅可以实现利用遥感手段定量估算近地面PM2.5质量浓度,而且可以为PM2.5对区域大气污染及流行病学影响预测研究提供科学支撑.

图1 研究区海拔及监测站点示意图Fig.1 The Elevation and Monitoring Stations in the Study Area

1 研究区概况

冀中南平原是华北平原的重要组成部分,总面积约8万km2,地跨石家庄、邢台、衡水、保定、邯郸、沧州和廊坊7市,是黄骅港及天津港的重要经济腹地.冀中南平原区属于典型的温带大陆性季风气候,夏季炎热多雨,冬季寒冷干燥.近30 a冀中南平原区城市化和工业化发展迅速,导致该地区成为中国大气污染最严重的区域之一.

2 数据及方法

2.1 数据来源

2.1.1 AOD产品数据

作为世界上先进的多光谱遥感传感器,MODIS共有36个观测通道,覆盖了当前主要遥感卫星的主要观测数据[24-27].MAIAC AOD产品数据空间分辨率为1 km×1 km,获取源为MODIS传感器搭载的Terra卫星(过境时间为地方时10:30)与Aqua卫星(过境时间为地方时13:30).数据通过NASA官网(https:∥portal.nccs.nasa.gov/datashare/maiac/Data Release/)下载,条带数据集h03v01,h03v02,h04v01下载覆盖京津冀地区的MAIAC AOD数据.AERONET AOD产品数据集由太阳光度计(CE318)测量获得,被广泛用于卫星遥感反演AOD验证.研究区内设有北京、邢台等4个站点,可以从AERONET官网(http:∥aeronet.gsfc.nasa.gov/)获取Level 2的数据.图2为2016年AERONET AOD与MAIAC AOD线性拟合图,决定性系数R2=0.889,样本数量N=432,均方根误差RMSPE=0.192,说明MAIAC AOD可以被用于近地面PM2.5质量浓度反演.

图2 AERONET AOD与MAIAC AOD拟合结果Fig.2 Fitting Result of AERONET AOD and MAIAC AOD

2.1.2 气象数据

2016年气象要素小时数据由河北省气象中心提供,共获取研究区内63个气象站点监测数据,包括气温、气压、大气能见度、相对湿度及水汽压等数据.气象台站每日2:00,8:00,14:00,20:00定时观测4次,为排除弱光线、辐射雾、逆温等天气因素的影响,研究使用14:00观测气象数据进行AOD缺值修复及湿度校正处理.

2.1.3 行星边界层高度(atmospheric boundary layer height,PBLH)数据

该数据来源于欧洲中期气象预报中心(European center for medium-range weather forecast,ECMWF)的气象再分析资料ERA5数据(http:∥apps.ecmwf.int/datasets/),其空间分辨率为13.875 km,时间分辨率为3 h,该数据用于AOD垂直校正处理.

将MAIAC AOD,PBLH等栅格数据与地面站点气象监测数据及PM2.5质量浓度数据进行匹配.

2.2 研究方法

2.2.1 AOD补值方法

云覆盖、亮地表等因素容易导致AOD产品数据缺失,因此基于地面气象数据,利用改进型Elterman模型(R-ERM)对MAIAC AOD产品进行补值,获取时空全覆盖AOD产品.

(1)

H1=e0.023 SV+0.031 WVP-0.402.

(2)

式中:TAOD为气溶胶光学厚度;SV为大气能见度,单位为km;λ为波长,单位为μm;H1为大气边界层内气溶胶标高,单位为km;Z为海拔,单位为km;H2为大气18 km高度以上的气溶胶标高,单位为km;WVP为水汽压,单位为hPa;v*表示每个对数直径区间内气溶胶颗粒(0.1～10 μm)的总体积,接近于常数3.

将R-ERM模型计算得到的AOD数据填补MAIAC AOD产品缺失数据,获得研究区时空全覆盖AOD产品.R-ERM模型模拟精度及有效性评价详细内容可参考文献[20].

2.2.2 垂直校正方法

AOD为气溶胶消光系数在整层大气垂直方向上的积分,表征整层大气中的气溶胶对光的削弱作用.而近地面PM2.5质量浓度表征近地面气溶胶质量浓度,因此需对AOD数据进行垂直校正,获取近地表气溶胶湿消光系数.使用PBLH对AOD进行垂直校正,校正方程为

bext=TAOD/HPBLH.

(3)

式中,bext为近地面气溶胶湿消光系数,TAOD为MAIAC AOD产品数据,HPBLH为行星边界层高度数据.

2.2.3 湿度校正方法

由于大气中的水汽会造成细颗粒物的吸湿增长,导致细颗粒物形状、粒径、消光性等理化性质发生改变.为降低大气湿度对气溶胶颗粒消光的影响,对垂直校正后的湿消光系数继续进行湿度校正,获取干消光系数,订正公式为

(4)

式中αext为干消光系数,bext为湿消光系数,WRH为相对湿度.

3 结果与分析

3.1 AOD时空变化

3.1.1 月平均值空间分布

图3为冀中南平原区2016年AOD月平均值空间分布图.可以看出,1月AOD高值区主要分布在冀中南平原中西部区域,约为0.9,低值区主要位于东北部区域,AOD值大都低于0.5.进入2月后,AOD值下降明显,2～5月数值大都处于0.29～0.50,且呈现相似的空间分布规律,即冀中南平原区周边区域AOD值略高,而中间区域数值相对较低.进入夏季后,研究区AOD值升高迅速,其中7月最为明显,高值区域主要分布在衡水—石家庄一线和邯郸南部地区,最高值达到1.98,低值区主要位于研究区北部及邢台地区;相较而言,6月和8月AOD值略低,大部分区域数值低于1.2,但空间分布模式均与7月相似.进入秋季后,冀中南平原区AOD值进一步降低,9～11月AOD数值为0.5～1.0,空间分布与夏季大体相似.12月AOD数值有所升高,最高值达到了1.3,高值区主要位于石家庄—衡水一线,南部邢台、邯郸等地部分地区AOD数值居中,保定、沧州、廊坊等地北部地区数值较低.

图3 AOD月平均值空间分布图Fig.3 Spatial Distribution of Monthly Mean AOD

3.1.2 季平均值与年平均值空间分布

图4为2016年冀中南平原区AOD季平均值及年平均值空间分布图.可以看出,AOD值在时间分布上整体呈夏季最高、秋冬季次之、春季最低的特点.空间分布上,春季AOD值由北向南呈由低到高的空间分布模式,邯郸南部为高值区,最高可达0.60,保定及廊坊北部地区为低值区,最低值约为0.44.与春季相比,夏季AOD值具有非常显著的升高趋势,冀中南平原区大部分区域AOD值超过了0.90,其中石家庄、保定白洋淀和邯郸南部地区数值最高,达到了1.43;邢台和邯郸北部地区数值相对较低,最低值约为0.90.秋季AOD值空间分布与夏季大体相似,但数值却大幅降低,约为0.60～1.10.冬季AOD值持续降低,最高值约为0.90,主要集中在石家庄地区,而保定和廊坊北部区域AOD值较低,最低值约为0.50.由图4可以看出,研究区内AOD整体呈现低—高—低的空间分布特征.其中,AOD高值区主要位于石家庄东部和邯郸南部地区,最高值约为0.96;石家庄东部—衡水一线存在一个稍低于石家庄的次高值区,最高值约为0.85;而AOD低值区主要集中在保定—廊坊—沧州一线北部地区和邢台东部地区,数值为0.66～0.75.

图4 AOD季平均值及年平均值空间分布图Fig.4 Spatial Distributions of Seasonal and Annual Mean AODs

3.2 垂直校正与吸湿校正

3.2.1 垂直校正

图5a为研究区地面站点PM2.5质量浓度监测数据月平均值与MAIAC AOD相关性图.可以看出,PM2.5与MAIAC AOD相关性较差,线性拟合决定性系数R2=0.000 5,主要原因是由于AOD表征整层大气中的气溶胶消光积分,而近地面PM2.5质量浓度只能表征近地表大气颗粒物浓度情况,因此很难基于AOD数据反演近地表PM2.5质量浓度.图5b为AOD经过垂直校正后得到的近地表消光系数与近地面PM2.5质量浓度相关性图.可以看出,经过垂直校正后的近地面消光系数与PM2.5质量浓度相关性明显高于AOD,决定性系数R2=0.603 5,有了显著提升.

图5 PM2.5质量浓度月平均值与AOD值(a)及湿消光系数(b)相关性Fig.5 Correlations Between Monthly Mean PM2.5 Concentration and AOD(a) and Wet Extinction Coefficient(b)

由图6可以看出,冬、春季节2者相关性较强,决定性系数R2=0.618 7和R2=0.533 2;其次为秋季,R2=0.378 2;夏季2者相关性最差,R2=0.004 3.主要原因可能是夏季大气相对湿度较大,气溶胶颗粒吸湿增长后,造成消光系数大幅度升高,从而导致近地表湿消光系数与PM2.5相关性不大;与之相反,冬春季节天气寒冷干燥,大气相对湿度较低,气溶胶颗粒吸湿增长作用不强,进而导致近地表消光系数与PM2.5相关性较强.

图6 湿消光系数季平均值与PM2.5质量浓度相关性Fig.6 Correlations Between Seasonal Mean Wet Extinction Coefficients and PM2.5 Concentrations

综上所述,虽然AOD经过垂直校正得到的近地表湿消光系数与PM2.5相关性程度有了很大的提升,但是由于气溶胶粒子的吸湿增长作用,限制了消光系数与PM2.5的相关性程度.因此,需要对近地表湿消光系数进行进一步湿度校正,获取干消光系数,进而提升近地表PM2.5质量浓度的估算精度.

3.2.2 湿度校正

图7为经垂直校正与湿度校正后的月均干消光系数与PM2.5质量浓度相关性图.可以看出,经过2次校正后,近地表干消光系数与PM2.5相关性程度有了明显提升,决定性系数R2=0.661 2,高于近地表湿消光系数与PM2.5相关性系数(R2=0.603 5).研究结果表明:湿度校正可有效降低大气湿度对气溶胶消光的影响,提高PM2.5质量浓度与干消光系数之间的相关性.

图7 干消光系数月平均值与PM2.5质量浓度相关性图Fig.7 Correlation Between Monthly Mean Dry Extinction Coefficients and PM2.5 Concentrations

3.3 近地面PM2.5质量浓度反演

3.3.1 PM2.5质量浓度估算值验证

图8为冀中南平原区PM2.5质量浓度估算值与监测值拟合图.可以看出,MAIAC AOD经过缺失修复、垂直校正与湿度校正后计算得到的近地面PM2.5质量浓度与站点监测值线性拟合度较好,决定性系数R2=0.700 8,说明基于MAIAC AOD遥感数据通过垂直和湿度校正能够较为准确地估算近地面PM2.5质量浓度.

图8 PM2.5质量浓度监测值与估算值拟合图Fig.8 Fitting Result Between Estimated and Observed PM2.5 Concentrations

3.3.2 冀中南平原区PM2.5质量浓度时空分布

图9给出了基于MAIAC AOD数据经垂直和湿度校正预测的冀中南平原区近地面PM2.5质量浓度月平均值.可以看出,1月和12月近地面PM2.5质量浓度最高,可达382.619 μg/m3,但高值区空间分布存在明显差异,1月高值区主要位于石家庄、衡水、邯郸南部及与山东分界区域,而12月则主要位于太行山山前平原区.2～3,9～10,11月PM2.5质量浓度相差不大,大都低于250 μg/m3,但各月PM2.5质量浓度空间分布差异性较大,2～3月高值区主要位于冀中南平原南部及以城市为中心的部分地区,9月和11月高值区主要位于石家庄地区,10月高值区则主要位于石家庄北部及廊坊南部地区.

图9 PM2.5质量浓度月平均值空间分布图Fig.9 Spatial Distribution of Monthly Mean PM2.5 Concentrations

图10为冀中南平原区PM2.5质量浓度季平均值及年平均值空间分布图.可以看出,2016年冀中南平原区PM2.5年平均值存在明显空间差异性,高值区主要分布在石家庄、定州东部及雄安新区等区域,并且PM2.5质量浓度大都超过了135 μg/m3,最高可达181.383 μg/m3;低值区主要位于沧州、衡水及邢台一带,年平均值大都低于120 μg/m3.

图10 PM2.5质量浓度季平均值及年平均值空间分布图Fig.10 Spatial Distribution of Seasonal and Annual Mean PM2.5 Concentrations

冀中南平原区PM2.5污染具有季节性特征,从图10可以看出,2016年研究区PM2.5质量浓度季平均值差异明显,呈现冬高夏低、春秋过渡的特征.冬季PM2.5质量浓度最高可达320.809 μg/m3,高值区主要位于冀中南平原区南部地区,低值区则位于东北部区域,PM2.5质量浓度也大都超过了150 μg/m3,主要原因是由于冬季取暖增加了污染物排放,且天气寒冷气流下沉,容易出现静稳天气,不利于污染物扩散.夏季PM2.5质量浓度最低,绝大部分区域质量浓度低于80 μg/m3,主要原因是夏季空气对流较强,利于污染物扩散,且降水集中于夏季,能起到清洁作用.

4 结 论

基于2016年MODIS AOD遥感数据和气象监测数据,通过对AOD遥感数据进行补值、垂直校正和湿度校正,有效提升了近地面PM2.5质量浓度估算精度.

1) 冀中南平原区AOD数值呈现夏高冬低、春秋过渡的季节分布规律;空间分布上,高值区主要位于石家庄—衡水沿线及邯郸南部地区,低值区主要分布在保定北部及邢台区域.

2) 经垂直校正和湿度校正后,AOD与PM2.5质量浓度相关性得到极大提升,决定性系数R2由0.000 5提升到0.661 2;经校正后的PM2.5质量浓度估算值与监测值拟合系数R2=0.700 8.

3) 2016年冀中南平原区PM2.5质量浓度时空变化差异显著.空间分布上,高值区主要分布在石家庄和定州东部及雄安新区等区域,低值区主要位于沧州、衡水及邢台一带.时间分布上,PM2.5质量浓度呈现冬高夏低、春秋过渡的特征,主要是由于冬季取暖增加了污染物排放且大气扩散条件差,而夏季降水集中且大气扩散能力强,有效降低了大气污染物质量浓度.