1. 引言
随着工业和经济的快速发展以及城市人口数量的增长,对能源的不断消耗使得排放出来的大量污染物严重污染了我们赖以生存的环境。人类对自然环境产生的破坏主要体现在土地、森林、水源、以及空气等方面 [1],其中对人类的生存发展影响范围最广的就是大气污染问题。由于环境质量恶化所造成的危害严重影响着人们的身体健康和生活方式,有研究表明PM2.5的污染特征对人类疾病的产生具有较高的相关性 [2]。空气污染问题的防治问题迫在眉睫。大气污染问题已越来越引起人们和各国政府的重视 [3],由国务院召开的全国环境保护大会多次在北京召开,汇集众人的智慧为中国的生态环境治理提供解决方案,不断推动生态文明建设,解决生态环境问题,构建美丽家园。
我国的空气质量监测起步较晚,但是发展较为迅速,从环境保护领导小组、环境保护局、国务院环境保护委员会,然后将国家环境保护局作为国务院环境保护委员会的办事机构。国家在保护环境的重任面前,不断地克服困难,坚决打胜这场污染防治攻坚战。想要为空气污染防治提出更为有效的解决办法,对空气质量监测提出更高的要求 [4],为制定解决策略提供理论基础。
根据《环境空气质量标准》(GB3095-2012),用于衡量空气质量的常规大气污染物共有六种,分别为二氧化硫(SO2)、二氧化氮(NO2)、粒径小于10 μm的颗粒物(PM10)、粒径小于2.5 μm的颗粒物(PM2.5)、臭氧(O3)、一氧化碳(CO)等 [5]。污染防治实践表明,建立空气质量预报模型,提前获知可能发生的大气污染过程并采取相应控制措施,是减少大气污染对人体健康和环境等造成的危害 [6],提高环境空气质量的有效方法之一。
近些年来,国内外对于空气质量预测的研究多集中于空气污染物的浓度预测。Liu [7] 提出了基于回归模型和支持向量机的空气质量预测模型,验证了组合预测模型的的预测效果较好。Zamani [8] 利用了深度学习等方法预测来了PM2.5的数值,结果表明,相比于深度学习和随机森林方法,XGBoost的预测误差最小。Bhat [9] 提出了多变量回归与双变量线性模型对印度克什米拉阳地区的PM2.5数值,经过分析对比,双变量线性模型的预测效果更好。Ma [10] 提出了基于网格重要等级的XGBoost的空气质量预测模型,预测结果的R2分数可以达到0.9以上。
目前对于大气的质量预报的研究方法主要分为两种,第一种方法是统计预报的方法,该方法的研究方法是通过大量地监测到的空气质量数据。通过统计学的分析方法,对收集到的数据进行分析,得到污染物浓度与气象条件之间的关系,通过建立预测模型得到空气质量的预报 [11]。另外一种方法是研究大气的动力学理论,通过对复杂的大气物理、化学变化的分析,建立大气污染物浓度在空气中的传输扩散数值模型,通过较为复杂的计算过程,借助于计算机完成运算,得到多种大气污染物在空气中的动态分布。
2. 监测点A单日AQI值数学模型计算
2.1. 数据预处理
提取2020年7月23日0时至2021年7月15日实测数据进行预处理,异常值采用剔除和拉格朗日插值的方法进行处理,用均值法将样本数据从逐小时处理为逐日数据,
1) 数据异常情形
关于一次预报数据:预报工作中,服务器受外接电源长时间停电等情况影响,导致部分运行日期的一次预报数据缺失。
关于实测数据:1) 因监测站点设备调试、维护等原因,实测数据在连续时间内存在部分或全部缺失的情况;2) 受监测站点及其附近某些偶然因素的影响,实测数据在某个小时(某天)的数值偏离数据正常分布;本题提供的监测气象指标共计五项(温度、湿度、气压、风向、风速),因不同监测站点使用设备存在差异,部分气象指标在某些监测站点无法获取。
2) 异常数据值处理
异常数据可能有多个来源,如数据本身、数据存储过程或者数据转换过程。由于异常数据会影响特征,也会影响最后的模型结果,因此对数据进行预处理十分必要。异常值的处理方法有多种,如删除记录、视为缺失值、平均值修正、不处理等。异常值如何处理,需要视具体应用背景分析而定。
2.2. 计算空气质量指数
根据《环境空气质量指数(AQI)技术规定(试行)》(HJ633-2012),空气质量指数(AQI)可用于判别空气质量等级。
首先需得到各项污染物的空气质量分指数(IAQI),其计算公式如下:
(1)
式中
为污染物P的空气质量分指数,
为污染物P的质量浓度值,
、
分别为与
相近的污染物浓度限值的高位值与低位值,
、
分别为与
、
对应的空气质量分指数。
空气质量指数(AQI)取各分指数中的最大值,即
(2)
式中,
为各污染物项目的分指数。在本节中,对于AQI的计算仅涉及六种污染物,因此计算公式如下:
(3)
3. 基于污染物浓度影响程度的气象条件分类
3.1. 问题分析
在某一地区,当假设污染量排放的浓度恒定不变的情况下,有些天气情况会直接影响该地区的AQI值。利用给定的监测点A的数据,分析各天气参数对该地区污染浓度的影响程度,并且对气象条件进行合理的分类,同时阐述各天气条件对该地区污染物浓度的影响特征。
分析天气条件对污染物浓度的影响,本节采用的是偏相关分析以及灰色关联分析,通过偏相关分析的控制变量的方法,不考虑天气条件之间的协同作用,单考虑单一因素对污染物浓度的影响,得到各天气条件与各种污染物的相关性,通过相关性对气象条件进行粗分类,分为对污染物浓度有正相关性或负相关两类。接下来采用灰色关联度分析模型,通过比较气象条件与各污染物浓度的关联度,以此来比较各气象条件对某一污染物浓度的影响程度强弱。为清晰化表述,设计如图1的流程图展示:
3.2. 偏相关性分析模型
首先利用附件给的数据计算样本的偏相关系数,当控制第三个变量来研究另外的变量之间的相关性时,一阶相关系数的表达式为
(4)
在实际应用中,本节的控制变量有四个,通过迭代法求取偏相关系数,迭代法认为简单的相关系数为0阶偏相关系数,任何n阶的偏相关系数都可以用3个(
)阶偏相关系数计算得到。
对变量间是否有净相关进行推断,采用t统计量作为偏相关分析的检验计量。其数学表达式如下定义:
(5)
其中r为计算的偏相关系数,n为样本数量,q为阶数。
通过对给定的数据进行分析,可以得到如图2所示的气象条件与污染物浓度的相关性,部分相关系数如图2所示。
通过对数据进行分析,可以得到各气象条件与污染物浓度之间的相关性关系,将各天气条件对污染物浓度的相关性作图,可以得到不同气象条件对于大气污染物浓度的影响程度。
由表可知对二氧化硫而言,大气压、感热通量、潜热通量、长波辐射、近地2米温度与二氧化硫浓度成正相关关系,地表温度、比湿、近地10米风向、边界层高度与二氧化硫浓度呈负相关关系,且显著性符合要求。对二氧化氮而言,近地两米温度、湿度、云量、感热通量、潜热通量与二氧化氮浓度成正比,比湿、近地10米风速与二氧化氮含量浓度呈强负相关性,且显著性复合要求。对PM10和PM2.5而言,地表温度、比湿、近地10米风速、雨量、长波辐射对粒径小的颗粒物有减弱作用,呈强负相关性,
注:黄色为在0.01水平上显著(双尾),绿色为在0.05水平上显著(双尾)
Figure 2. Correlation coefficient between meteorological conditions and pollutant concentration
图2. 气象条件和污染物浓度相关性系数
而近地2米温度、湿度、大气压 、潜热通量对PM10和PM2.5有促进作用,这些指数高的情况下会加重污染物浓度。
对臭氧而言,地表温度、湿度、近地10米风速、边界层高度、大气压、长波辐射与臭氧含量成正相关,且相关系数较大,而近地2米温度、感热通量对臭氧浓度呈负相关。对一氧化碳而言,近地2米温度、湿度、感热通量与一氧化碳浓度呈正相关性,且相关系数值较大,相关性显著,比湿、近地2米风速与一氧化碳浓度呈强负相关性。
3.3. 散点图分析模型
本节以监测点A逐小时污染物浓度与气象一次数据作为样本,将大气污染物浓度作为对比变量,将重要的气象条件变量与空气污染物浓度的关系用散点图进行分析,观察气象条件变量与空气污染物浓度之间是否存在一定关系,对指标进行了相关性分析。如图3所示是气象条件与二氧化硫浓度分布的散点图的部分结果。
从图3中可以看出,雨量、感热通量与二氧化硫浓度的散点图呈现y轴垂直的式样,这表示雨量、感热通量的变化并不能直接影响二氧化硫的浓度变化,从图中其他的点上也不能够观察出指数、对数、幂次的非线性关系。但我们可以看出散点集中于坐标轴左下方,说明二氧化硫浓度的分布集中于在雨量和感热通量较小环境下,而当雨量与感热通量较大时,二氧化硫浓度的散点分布在靠近横轴的区域。
3.4. 灰色关联度分类模型
首先参考数列的选取,以大气污染物浓度作为参考数列,记为X0,比较数列的选取,以气象条件作为比较数列,记为Xi。在附件中给出的气象参数中,由于CO的浓度与其他的污染物浓度的单位不同,需要对给出的数据进行无量纲化处理。这里采用极差值法对数据进行无量纲处理,以此来消除量纲给分析结果带来的影响。极差值法消除量纲的方式如下:
令
(6)
接下里利用关联度系数的计算公式计算关联度系数,公式为:
(7)
为
和
在第k个点的绝对误差,
为两级最小值,
为分辨率,
,在解决该问题时选取
为0.5,
与分辨率成反比,
越大分辨率越小,反之则越大。

Figure 3. Scatter diagram of partial meteorological conditions and SO2 concentration distribution
图3. 部分气象条件与SO2浓度分布散点图
最后再根据关联度计算公式来计算其关联度,关联度计算公式如下所示:
(8)
式中ri为xi对x0的关联度。
根据灰色关联度分析,得到相应的关联度数排名结果如图4所示。
以气象条件与污染物浓度的灰色关联度数为0.5作为刻度,将气象条件分为两类,分别是对二氧化硫影响较大的气象条件和对二氧化硫影响不大的气象条件。各气象因子对污染物浓度的关联程度排序为:雨量 > 感热通量 > 潜热通量 > 短波辐射 > 地面太阳能辐射 > 边界层高度 > 近地10米风速 > 云量 > 大气压 > 比湿 > 地表温度 > 湿度 > 长波辐射 > 近地2米温度。由图可知大气压、云量、近地10米风向、近地10米风速、边界层高度、地面太阳能辐射、短波辐射、潜热通量、感热通量、雨量对二氧化硫浓度影响较大,且雨量相较于其余气象条件变量,影响效果更大。近地2米温度、长波辐射、湿度、地表温度、比湿对二氧化硫的浓度的影响效果并不明显,所以将其归为第二类。
3.5. 距离相关验证分析
在统计学中,皮尔逊相关系数是用来度量两个变量之间的线性相关,其值在−1和1之间。可以用来判断数据样本是否集中分布在同一条线上,用来衡量定距变量之间是否存在特定的线性关系。而样本

Figure4. Line chart of grey correlation degree
图4. 灰色关联度折线图
之间线性关系不明显,因此皮尔逊相关性分析并不适合本体的数据样本。距离分析可以分为个案间和变量间举例分析两种,分析的方法有相似性和不相似性分析两种。距离分析客服了传统的相关分析的缺点,能够分析非线性数据变量之间的关系,对于两个随机变量x,y而言,距离相关系数可以定义为
(9)
其中:
(10)
同理:
(11)
本节所测量的变量如表1所示:
最后我们采用距离相关系数对分类效果进行验证,相关系数热力图如图5所示,由图可知,近地2米温度、长波辐射、湿度、地表温度、比湿这几个变量的相关性均在0.7以上,说明对大气污染物浓度影响不明显的一类中变量相关性显著,而其余变量之间的相关性较小,符合上述分类结果,故分类效果较好。

Figure 5. Heat diagram of correlation coefficient of meteorological condition variables
图5. 气象条件变量相关系数热力图
4. 大气污染物浓度预测
4.1. 大气污染物浓度趋势时间序列分析
ARIMA模型是差分模型、自回归模型和移动平均模型的结合,差分模型是为了实数剧更加平缓。
1) 自回归模型AR
自回归模型是利用因素自身的历史数据对自身进行预测,预测只能在满足平稳性要求的情况下才能进行,平稳性检测采取ADF指标进行检验。
自回归模型表达式如下:
(12)
式中,
表示当前值,μ为常数,p为阶数,
为自相关系数,
为误差。
2) 移动平均模型MA移动平均模型是为了消除自回归模型的随机波动即误差
,其表达式如下:
(13)
如图6所示,时间序列模型的流程可总结如下:

Figure 6. Flow chart of time series model
图6. 时间序列模型流程图
本文利用SPSS软件对数据进行时间序列预测,发现污染物浓度的规律有明显的季节性,根据时间序列预测模型,本文通过两年中每个月的同一天的污染物浓度变化规律设计时间序列模型,得出三天的污染物浓度变化数值,部分数据趋势图如图7所示。

Figure7. 3-day concentration change trend diagram
图7. 三天浓度变化趋势图
对与大气污染物浓度的时间序列预测,从总体而言建立的时间序列预测模型效果较好,由于预测事件较短,预测出来的结果符合本节关于时间序列图中的大气污染物浓度趋势,在7月份夏季时部分污染物浓度处于低值,符合预期,通过观察时间序列预测的结果可以发现预测结果与历史周期规律高度一致,但是严格来说由于大气污染物浓度的变化和气象条件的变化具有时间滞后性,所以在现实应用中参考性不大,无法保证其准确性,因此本文使用BP神经网络预测模型加以参考,提高预测结果的可靠性。
4.2. 大气污染物浓度趋势BP神经网络分析
二次空气质量的预测模型,采用的是多层感知器的误差反向传播算法,也就是所说的BP神经网络,具体的学习步骤如图8所示。

Figure 8. BP neural network learning flow chart
图8. BP神经网络学习流程图
建立BP神经网络模型时,将样本通常分为训练集和测试集,训练集用来构建BP神经网络的模型结构,测试集用来分析和检验模型的实际效果,一般通常以百分之八十和百分之二十分配数据样本个数,通常在设计神经网络结构时,隐含层层数不超过两层,隐含层层数设置过多不仅会减慢其收敛速度还会对模型最后的预测精度造成一定影响,本文采取三层经典神经网路结构:输入层,隐含层,输出层。对于节点个数的选取,本文选取四个节点用来设计神经网络模型。
1) 训练数据和测试数据确定
对于附件一中的数据,取监测点A逐小时污染物浓度与气象一次预报数据作为神经网络训练数据,将监测点A逐小时污染物浓度与气象实测数据的后50行的气象数据作为参数输入到网络中,得出后50行的六种污染物的浓度值,并将此数据与A点逐小时的实测数据进行误差分析,验证二次模型的可行性。基于此设定,先将数据进行归一化处理:
(14)
2) 参数设定
本研究设计的BP神经网络图设计4个参数,分别为迭代次数、隐含层元素个数、训练目标最小误差和学习速率。由于设计数据较大,本研究设计的迭代次数为10,000,训练目标最小误差为0.000001。然后采用控制变量法分别讨论隐含层元素个数和学习速率的设定,以要求六种污染物浓度的相对误差为最小。首先设定学习速率为0.1,隐含层个数分别为20、45、60、90,然后求得在不同隐含层下各污染物的总偏差。通过对比污染物的总偏差大小和偏差均值来确定隐含层元素的个数。
从图9可以看出当隐含层元素个数为20时总偏差较为平稳且偏差均值最低。因此设置隐含层的个数为20。在此基础上,分别设定学习速率为0.1、0.4、0.7、0.9。

Figure 9. Comparison of total deviation of different hidden layers and total deviation of different learning rates
图9. 不同隐含层总偏差与不同学习速率总偏差对比图
结合图9可以看出,当学习速率为0.7时污染物的总偏差平稳且偏差均值为最小。综上所述,该神经网络的迭代次数为10,000,训练目标最小误差为0.000001,隐含层元素个数为20且学习速率为0.7。
在设定好的神经网络下,将附件1中的监测点A逐小时污染物浓度与气象实测2021/7/13日之前的数据输入到神经网络中进行训练,将监测点A逐小时污染物浓度与气象一次预报2021/7/13~2021/7/15的五个气象数据(温度、湿度、气压、风速、风向)作为参数输入到神经网络中,来输出7/13~7/15逐小时的六种污染物浓度。最后求取3天逐小时的污染物平均浓度作为3天逐日的污染物浓度并且算出AQI。
将两种模型得出的结果进行对比可以发现,时间序列模型的短期预测结果一般,最大相对误差较大,BP神经网络的预测结果较好,最大相对误差较小,但浓度分布均处在正常范围之内,因此可以接受,在条件允许的情况下,可综合考虑时间预测模型和BP神经网络预测模型的预测结果。
5. 粒子群算法优化
粒子群优化算法对神经网络的调整过程如下所示:1) 对参数运用网络自身进行二次优化,直到得到最佳的权值和阈值之后停止搜索。2) 粒子群替代梯度下降法对神经网络权值、阈值进行优化,直到计算出来的适应度值无法继续下降迭代才会停止。神经网络中的权值和阈值数据会被粒子群算法充当粒子的位置向量,通过寻求粒子群的最佳位置同时寻找最佳权值和阈值,再利用BP神经网络的正向传播计算粒子群算法的适应度。
算法流程图如图10所示。

Figure 10. Logical diagram of BP neural network prediction model based on particle swarm optimization
图10. 基于粒子群优化的BP神经网络预测模型逻辑图
将修正好的一次预报数据,作为训练集输入到基于粒子群优化的BP神经网络中,将监测点逐小时的实测数据的后50行的气象条件作为参数输入到网络中,得到后50行的六种污染物浓度的预测值,并与逐小时的实测数据进行对比,如下图11所示。
由图11可知,由于对一次预报数据进行了修正及对BP神经网络进行了粒子群算法的优化,优化后的模型预测精度更高,准确度更好。
6. 模型评价
如图12所示为不同模型下六种污染物浓度的总偏差,其中经过粒子群优化的BP神经网络总偏差较小且相对未优化的BP神经网络更加平稳。
为了便于比较,取均方根(RMS)为性能指标:
(15)
由表2可知,均方根值随着数据的增大而增大,且已优化的BP神经网络预测出的六种浓度的均方根值远远小于未优化的BP神经网络预测值。

Figure 11. Comparison of pollutant concentration prediction of different models
图11. 不同模型的污染物浓度预测对比图
7. 结语
1) 根据大气污染物浓度变量的取值范围,剔除不在范围的变量数据样本;对数据缺失值运用拉格朗日插值法进行填充。根据预处理后的数据,通过给定的环境空气质量分指数公式计算监测点A处各项污染物的IAQI值,再取污染物浓度计算得到的最大值作为最终单日的AQI值,最大值的污染物作为首要污染物。

Figure 12. Comparison of total deviation of pollutants in different models
图12. 不同模型的污染物总偏差对比图

Table 2. Root-mean-square table of pollutant concentration prediction of different models
表2. 不同模型的污染物浓度预测均方根表
2) 根据偏相关系数分析气象条件之间的相关关系,再利用散点图分析气象条件对污染物浓度的影响规律。利用灰色关联度分析气象条件与污染物浓度之间的关联度,以关联度的强弱将气象条件分成两类。最后通过距离相关得出不同类别气象条件之间的相关系数热力图验证分类的准确度。
3) 本文首先采用时间序列预测模型根据历史数据对7月13日至7月15日的数据进行初步预测,再利用BP神经网络模型进行污染物浓度的预测。然后将一次预报数据作为训练集输入到经过粒子群算法优化的BP神经网络中,得到三点的二次预报数据,利用该数据修正A点的二次预报数据,采用BP神经网络的预测方法得到预测值和实测值的误差,得出二次模型预测效果良好。
本文在一次预测模型的结果的基础上,研究表明BP神经网络的预测效果要大于时间序列预测效果。结合更多的空气数据对预测模型进行二次建模,二次建模的结果准确性大大优于原始的预测结果,利用粒子群算法对BP神经网络的预测模型进行优化,能够不断对参数进行更新与修正,从而实现对空气质量的准确预测,预测结果有一定的鲁棒性,对后续空气质量预测有一定的研究意义。