近年来，随着全球能量消耗，污染物排放的加剧，人们对空气质量指数(Air Quality Index，AQI)的关注程度与日俱增，一方面关注影响AQI的污染源，减少污染排放；另一方面研究对AQI的预测，以积极做好健康防护。

在AQI预测方面，多种先进智能算法被用于谋求精准预测。高帅等^[1]提出将支持向量机(SVM)和飞蛾扑火优化(MFO)相结合的算法(MFO-SVM)，张梦瑶等^[2]基于改进加权马尔科夫链，李博群等^[3]通过建立模糊时间序列模型，对一些地区的AQI进行实证分析和检验，获得了较为理想的预测结果。此外，基于小波分析^[4]、主成分分析^[5]、时间序列分析^[6]的方法也被用于AQI预测。在诸多智能算法中，神经网络作为一个高度复杂的非线性动力学系统，也经常被用于AQI预测。例如文献[7-10]所建立的AQI预报模型都是基于神经网络完成的。然而，他们对神经网络输入的处理大多基于AQI的主要污染源，忽略了次要污染源的贡献；而若将AQI时间序列本身作为神经网络输入，虽可包含所有污染源数据信息，但同时会降低神经网络的泛化能力和预测精度，因此，基于自回归时间序列输入的神经网络在AQI预测方面的研究还不多见。

现使用2013—2018年徐州地区的AQI时间序列，引入自回归输入模式，采用遗传算法(Genetic Algorithm，GA)对神经网络初始权值进行最优化筛选，利用BP算法对网络参数反向微调，并对AQI的预测结果进行了验证。

1 数据预处理及网络输入端构建 1.1 数据预处理

从江苏省气象局网站收集整理徐州地区2013年12月—2018年11月共1 814个AQI作为原始数据。原始AQI数据检验观测图见图 1。

由图 1可见，时间序列存在2类“奇异”点：一类为脱离主数据群，距离主数据群较远点，如AQI值为500的一点，距离主数据群最大值直线位移超过150；另一类是AQI值为“0”的点，该类点意味着空气质量出奇得好，没有一点污染。对于第一类奇异点，原数据组仅有1个数据，从数据预测的经验判断，该点破坏了原数据组的变化规律，会降低预测模型精度，故予以删除；对于第二类奇异点，尽管该类点距离主数据群很近，如图 1所示数据序列号为1600附近的点，但基于生活经验的常识，将“AQI为0”视为不可能事件，故将原数据组的此类14个点也予以删除。为保证原数据组的时序完整性，所删除的15个点采用cubic插值补齐。

插值补齐后数据组数据范围为[26, 343]，极差为317，仍较大，若直接作为神经网络的输入，容易造成系统震荡而出现数据欠拟合、无法收敛、预测精度降低、迭代过程中系统震荡过大导致不稳定等结果。采用神经网络常规数据处理方法，将原数据序列归一化到[-1, 1]，后发现数据由于过于密集而无法反映原数据组的变化趋势，故不考虑神经网络对数据输入处理的常规归一化方法。

建立神经网络的输入为一维时间序列，量纲统一，不受常规神经网络因多维数据输入的量纲不统一，必须进行归一的限制。因此，考虑将原数据序列映射到某一区间[x₁, x₂]上，该区间[x₁, x₂]应达到使所建立神经网络的预测值与真实值的相对误差最小。取网络预测的均方误差(Mean Square Error，MSE)为目标函数，建立对最优区间搜索的表达式：

$ \mathop {\min R}\limits_{\left[ {{x_1}, {x_2}} \right]} = {\rm{MSE}}\left( Y \right) $

(1)

式中：R——区间[x₁, x₂]内的均方根误差值，Y——所建立的神经网络。

上式的求解问题是个解空间搜索区域极大的最优极值求解问题，可采用粒子群(PSO)、蚁群等多种仿生智能算法求解。现讨论的核心是GA-BP神经网络的AQI预测性能，并非最优极值求解，对上式的求解问题做如下简化：

(1) 考虑计算成本，固定x₁=0，将[x₁, x₂]双边求解简化为[0, x₂]区间的单边求解，其中0 < x₂ < 343；

(2) 基于降低所要建立神经网络的复杂性，Y使用一般BP神经网络。

基于上述简化，为确定x₂的较优解，将使用一般神经网络在MATLAB中编程求解。首先需要确定神经网络一维时间序列的输入端节点数。

1.2 自回归网络输入端构建

由于AQI数据为时间序列，不符合神经网络监督式训练模式，引入自回归模型框架，在输入端对原数据序列进行转换，以符合网络结构需求。用样本总体计算不同自回归阶数所对应的网络均方误差及测试样本相关度检验值见表 1。

由表 1可见，选择不同自回归阶数的MSE基本稳定在0.07左右，但相关度检验值存在较大差异。当自回归阶数为10时，MSE比阶数15时高0.000 8，但阶数15对于测试数据的相关度检验明显低于阶数10，选择自回归阶数为10。即以连续10 d的AQI数值为输入端，以第11 d的AQI作为标签数据，以此类推确立网络样本数据共1 804组。使用上述数据，建立输入端为10的一般BP神经网络，在MATLAB中编程，使用(1)式作为目标函数，搜索较优的x₂值为224，故确定原时间序列的输入映射区间为[0, 224]。图 2为插值并映射到[0, 224]后的数据序列。

2 GA-BP神经网络预测模型

(1) 由上文可知设定输入端节点为10；隐含层节点数目通过试错法选定，网络迭代开始前初始化各连接赋权ω_ij¹、ω_ij²及阈值θ，其中θ赋予[-0.2，0.2]的随机值，i=1, 2……10；j=1，2……1 804。

(2) 确定隐含层及输出层的传递函数为tansig函数并归一化样本到范围[-1, 1]。

(3) 将各隐含层权值进行遗传编码计算适应度并进行选择、交叉、变异操作，当种群适应度到达目标值时保留此时网络隐含层权值，设置其适应度函数与隐含层输出为：

$ fun = 1/\sum\limits_{i = 1}^n {{\rm{abs}}\left( {{d_i} - {\xi _i}} \right)} $

(2)

$ {\xi _i} = {\rm{tansig}}\left( {\sum\nolimits_{i = 1}^m {\omega _{ij}^1{x_{ij}} - {\theta _i}} } \right) $

(3)

式中：fun——种群适应度；d_i——第i个网络输出层节点的期望输出；ξ_i——第i个网络输出层节点的仿真输出；ξ_j——第j个隐含层节点输出；ω_ij¹——第i个输入层节点到第j个隐含层节点的权值；x_ij——第i个输入端节点到第j个隐含层节点的样本值；θ_i——第i个输入层节点阈值。

(4) 计算

$ {\sigma _i} = \left( {{d_i} - {\xi _i}} \right){\xi _i}\left( {1 - {\xi _i}} \right) $

(4)

式中：σ_i——第i个输出层节点校正误差；d_i——第i个网络输出层节点的期望输出；ξ_i——第i个网络输出层节点的仿真输出。

(5) 利用ω_ij²，ξ_i，ξ_j计算隐含层的校正误差：

$ {\sigma _{ij}} = {\xi _i}\left( {{d_i} - {\xi _i}} \right){\xi _j}\omega _{ij}^2。$

(5)

式中：σ_ij——第i个输入节点到第j个隐含层节点的修正误差；ξ_i——第i个网络输出层节点的仿真输出；ξ_j——第j个隐含层节点输出；ω_ij²——第i个隐含层节点到第j个输出层节点的权值。

(6) 若其误差收敛于规定精度ε时进入停止学习，反之利用ω_ij²，ξ_i，ξ_j和θ计算下一次学习中隐含层和输出层之间新的连接权值和神经元阈值，并对各个参数进行反向微调：

$ \begin{array}{l} \;\;\;\;\;\;\;\theta \left( {t + 1} \right) = \theta \left( t \right) + \gamma \left\{ {\eta \left( t \right){\sigma _i} + \alpha \left[ {\theta \left( t \right) - } \right.} \right.\\ \left. {\left. {\theta \left( {t - 1} \right)} \right]} \right\}, \eta \left( t \right) = {\eta _0}\left[ {1 - t/\left( {T + M} \right)} \right] \end{array} $

(6)

$ \begin{array}{l} \;\;\;\;\omega _{ij}^2\left( {t + 1} \right) = \omega _{ij}^2\left( t \right) + \gamma \left\{ {\eta \left( t \right){\sigma _i}{\xi _i} + \alpha \left[ {\omega _{ij}^2\left( t \right) - } \right.} \right.\\ \left. {\left. {j\left( {t - 1} \right)} \right]} \right\} \end{array} $

(7)

$ \begin{array}{l} \;\;\;\;\omega _{ij}^1\left( {t + 1} \right) = \omega _{ij}^1\left( t \right) + \gamma \left\{ {\eta \left( t \right){\sigma _{ij}}{\xi _i} + \alpha \left[ {\omega _{ij}^1\left( t \right)} \right.} \right.\\ - \left. {\left. {j\left( {t - 1} \right)} \right]} \right\} \end{array} $

(8)

式中：θ(t+1)——第t+1次更新后的阈值；θ(t)——第t次更新后的阈值；γ——学习率；η(t)——第t次更新后的梯度；σ_i——第i个输出层节点校正误差；α——动量系数，且α∈(0, 1)；θ(t－1)——第t-1次更新后的阈值；η₀——初始梯度；t——学习次数；j=1，2……1 804；T——总迭代次数；M——为防止分母为0的整数；ω_ij²(t+1)——第t+1次更新后第i个隐含层节点到第j个输出层节点的权值；ω_ij²(t)——第t次更新后第i个隐含层节点到第j个输出端的权值；ω_ij¹(t+1)——第t+1次更新后第i个输入层节点到第j个隐含层节点的权值；ω_ij¹(t)——第t次更新后第i个输入端到第j个隐含层节点的权值；σ_ij——第i个输入节点到第j个隐含层节点的修正误差；ξ_i——第i个网络输出层节点的仿真输出；ξ_j——第j个隐含层节点输出。

(7) 停止学习并保存训练好的各网络参数。算法整体流程见图 3。图中k为进化代数，K为预设进化代数，ε为网络误差，σ为网络期望误差。

3 仿真实验

为保证对网络泛化能力的测试，对于夏季，截取2018年7月前后共50 d样本作为测试集，训练样本总数1 637组；对于冬季，截取2018年11月底之前共50 d样本作为测试集，冬季测试时样本总数为1 754组。考虑到夏、冬季AQI所表现的范围并不相同，分别对夏、冬季进行连续的长期预测以展示所建立GA-BP神经网络模型的性能。

经过调试，所确定的GA-BP预测模型参数如下：输入端节点10，隐含层节点30，输出端节点1，种群为30，种群进化代数10，交叉概率0.4，变异概率0.1，学习率0.1，动量参数0.9，最大迭代次数100，在MATLAB环境下仿真输出的夏季和冬季的预测结果见图 4(a)(b)(RMSE为均方根误差)。

由图 4可见，计算夏、冬季2个预测系统所对应的MSE分别为0.052和0.046，均接近于0；其中夏季预测相对误差18.23%，仿真RMSE为14.59；冬季相对误差9.14%，仿真RMSE为11.47，故认为训练效果良好，2个系统的RMSE相差不大，但可看出冬季相对于夏季其仿真结果要更加稳定一些；2个系统虽在趋势追踪上效果良好，但对局部极大或极小值的追踪上效果相对要差一些。

4 结论

建立的网络模型能够准确预测徐州市空气质量指数的变化趋势，其中夏季预测相对误差18.23%，仿真RMSE为14.59；冬季预测相对误差9.14%，仿真RMSE为11.47。