运用机器学习方法预测空气中臭氧浓度

蔡旺华

（ 福建省环境信息中心，福建福州 350003 ）

臭氧（O3）是一种淡蓝色气体，是天然大气的重要微量成分,约90%的臭氧存在于平流层,仅有10%左右的臭氧分布在对流层中[1]。臭氧是《环境空气质量标准》（GB3095—2012）中包含光化学烟雾的标志性污染物，高浓度近地面臭氧对人体健康有较大的危害[2]。城市中的臭氧主要由氮氧化物（NOx）、挥发性有机物（VOCs）在适宜条件下作用产生[3],主要来源于机动车、各类工厂和自然植被。随着经济活动的加剧和机动车保有量的激增，臭氧污染问题渐趋严重[4]，2017年福建省空气质量的首要污染源即为臭氧。如何准确预报臭氧浓度成为一个重要课题。

常见的大气污染物浓度预测方法主要有两类：基于数值模式的预测方法[5,6]和基于统计模型的预测方法[7-10]。数值模式的预测方法是以大气动力学理论为基础，在给定的气象条件、污染源排放清单以及初始边界条件下，通过一套复杂的偏微分方程组模拟大气污染物在实际大气中的各种物理化学过程，预报污染物浓度动态分布和变化趋势。基于统计模式的预测方法以空气质量监测数据、气象监测数据等多种类型数据为基础，利用统计方法，建立预测模型，实现大气污染物浓度预测。传统的统计模型预测仅考虑单一大气污染物浓度变化，利用ARΙMA、Holt-Winters等时间序列预测模型进行浓度预测[11]，预测精度较低且预测时效较短。后有方法提出利用大气污染物浓度和温度、湿度等气象因素进行多元线性回归，预测精度相比ARΙMA等时间序列预测模型要准。近年来，随着机器学习算法的不断发展，支持向量机、人工神经网络等机器学习算法被不断引入大气污染物浓度预测任务中来[12,13]，这些方法对多源数据进行综合考虑，预测效果有所提升。

相比基于数值模式的预测方法，基于统计模型的预测方法具有以下两个优点：一是此类方法的输入采用真实监测值，比采用估算值的数值模式更加准确[14]；二是此类方法计算速度快，通常在秒级即可完成计算，而数值模式运算消耗时间动辄以小时计。但是，目前国内常见的基于统计模型的大气污染物浓度预测方法一般只对空气质量监测数据、气象数据进行回归分析。由于输入特征单一，未考虑到引起大气污染物浓度变化的可能原因，如点源、移动源排放变化等，从而导致了现有此类方法往往表现得不够理想。

本文提出了一种结合臭氧浓度变化机理的机器学习方法用于臭氧浓度预测，是福建省生态环境大数据平台建设中的大数据应用之一。福建省生态环境大数据平台由福建省环保厅于2015年9月启动建设，在建设过程中，正值2017年9月3日至5日金砖国家领导人厦门会晤。为此，在8月31日至9月9日，运用本文提出的方法对厦门市臭氧浓度进行了较为准确的预测，为保障空气质量提供了支持。

1 臭氧浓度预测模型

1.1 臭氧预测模型机制

臭氧浓度变化可以简单分为生成、分解和扩散三部分[15,16]。其中，生成部分主要考虑三类污染源（天然源、移动源和点源）的影响。因为在一段时期内，可以认为天然源的影响基本稳定，无人工干预下的移动源的影响呈现较规律的周期性变化，所以此方法主要考虑相对变化较大的点源对臭氧浓度变化的影响，从而进行臭氧浓度预测；分解和扩散的部分主要考虑湿度、风速等气象条件对臭氧浓度变化的影响[17]。因此本方法采用的数据主要有以下四类：①空气质量监测数据；②周边点源（以下本文所提“点源”均指有在线监测废气排放的污染源）排放数据；③气象监测数据；④预报时刻的气象预测数据。根据各类数据对臭氧浓度变化的影响构造特征，预测24小时后的臭氧浓度。

因为各个大气监测站所处位置不同，各种条件一般会有较大差异，臭氧浓度变化规律也可能不同，所以对每个大气自动监测站单独训练臭氧浓度预测模型。

该方法的流程如图1所示。首先获取上文提及的空气质量监测数据、周边点源排放数据、气象监测数据和预报时刻的气象数据等四类数据，然后进行数据预处理、特征抽取步骤，获取训练数据集和测试数据集。然后利用训练数据集和机器学习的XGBoost算法进行模型训练，得到预测模型，再利用预测模型和测试数据集得到预测臭氧浓度。

图1 臭氧浓度预测模型流程

图2 臭氧浓度预测模型架构

臭氧浓度预测模型架构如图2所示，对图1的特征抽取和模型训练两部分进行展开。由图2可见，模型输入特征基于起报时刻t前h天的点源数据（图1中h=5，即t前120小时）、臭氧浓度监测数据和气象监测数据以及预报时刻气象数据构造。模型目标为预报时刻（t+24）的臭氧浓度相较于起报时刻（t）臭氧浓度的变化值。采用变化值为目标是为了减小臭氧浓度周期性变化以及一些潜在未知因素的影响。基于XGBoost算法，可以获得预报时刻相较起报时刻的臭氧浓度变化预测值，再结合起报时刻的臭氧浓度，即可获得预报时刻的臭氧浓度预测值。

1.2 特征抽取

本文提出的算法包括三部分特征：点源特征、臭氧浓度和气象特征，以下依次进行说明。考虑到各个大气监测站点所处的位置不同，该算法针对单个大气监测站点设计。

点源特征按照点源地理位置分布抽取，如图3所示，五星点为目标大气监测站点，黑色圆点为点源。针对每个大气监测站点，考虑与其所属城市接壤的所有城市内的点源，构造如图2的区域，并将其划分为若干个长宽均为M的小方格。针对每个小方格，计算其中心点和当前环境监测站点的相对方向和距离（图2中所示的叉划线），同时统计其所包含的点源，则所有点源均明确隶属于某个方格。对前述已划定的每个小方格，计算起报时刻t前h天内的机制特征，以每24小时为1个时间段，每个时间段抽取二维特征，分别是污染物扩散规律特征和污染物传播时间特征。对每个小方格每个时间段，若主要风向有利于污染物到当前关注的环境监测站点的输送，即可以根据平均风速和点源排放量计算该二维特征。当风向不易于污染物输送或者方格内暂无点源排放时，该二维特征均为0。

图3 点源网格划分特征抽取示意图

臭氧浓度特征采用起报时刻t的臭氧小时平均浓度。

气象特征通过预报时刻t+24前若干个小时的气象因素构造，因为湿度和风速对臭氧的扩散和分解至关重要，所以本算法中采用其预测时刻前4小时、前6小时和前24小时的最大值、最小值和平均值作为特征。

将上述三种特征直接拼接构成输入特征，对应的输出结果为预测时刻t+24与起报时刻t的臭氧浓度差值，之后利用XGBoost算法拟合针对每个大气监测站点的模型即可。

1.3 XGBoost算法

本文采用的XGBoost[18]（eXtreme Gradient Boosting）是一种梯度提升算法扩展。梯度提升算法是一种用于回归和分类问题的机器学习技术，以多个弱预测模型集成的形式产生预测模型。通常会选择树这一弱预测模型作为基学习器，例如，本文所提出预测模型采用回归树作为基学习器。当然，XGBoost基学习器也可以是其他分类器。

假设XGBoost一共生成了k棵树，对任一输入x，每棵树都有一个输出fx(x)，则该XGBoost模型的输出为在迭代过程中，利用梯度下降思想，以已生成树为基础，以最优化目标函数为总目标，迭代生成新树。每次迭代学习的目标是减少已生成的模型的累计结果的损失值。不同于传统梯度提升决策树（gradient boosting decision tree，GBDT），XGBoost对 代价函数进行二阶泰勒展开，利用了一阶导数和二阶导数。另外，XGBoost的代价函数中包含正则项，这相当于给复杂模型以一定的惩罚，模型复杂度得以控制，所以模型不易出现过拟合问题。

XGBoost算法能够利用CPU多线程并行，具有运行速度快，准确度较高，不易过拟合等优点[19]，算法普适性好，在多种应用上表现良好，但未见在臭氧浓度预测的文献。

2 模型预测

2.1 数据准备

该方法采用的各类数据源，除测试过程所需的气象预报数据外，均来源于福建省生态环境大数据平台。该方法涉及的数据源包括以下四类：大气监测站的空气质量监测数据、大气监测站的气象监测数据、点源（在线监测的废气排放污染源）数据和气象预报数据。大气监测站的空气质量监测数据格式参见表1，以小时平均深度代表每个大气监测站点在每个整点时刻的臭氧浓度；大气监测站的气象监测数据和气象预报数据的数据格式相同，参见表2，记录每个大气监测站点每个时刻的温度、湿度、气压、风速和风向等5个监测指标；点源监测数据格式参见表3，记录每个点源小时排放监测值的NOx折算浓度。

表1 大气监测站的空气质量监测数据

表2 大气监测站的气象监测数据

表3 点源（在线监测的废气排放污染源）数据

该方法采用2016年9月1日至2017年8月30日一年的上述三种数据源，针对厦门市溪东、洪文、鼓浪屿和湖里中学四个大气监测站点各自的训练模型。利用金砖领导人厦门会晤前后共10天（2017年8月31日至9月9日）的数据和对应的气象预报数据进行模型测试。训练数据和测试数据量如表4所示，因数据缺失或异常，各站点的训练、测试数据条数可能不同。

表4 各站点数据量统计

2.2 模型训练和参数设置

对每个大气监测站点，在训练集上使用10折交叉验证方法进行模型训练和参数调优，最后在测试集上测试结果。所谓10折交叉验证，即将数据集划分为10个子集，依次取其中9个子集为训练集，剩余1个子集为测试集，最终以在10个测试集上的平均精度评价模型。

本次实验在训练XGBoost模型时设置了三个参数：树的深度（max_depth）为4、学习率（learning_rate）为0.1和树的个数（n_estimators）为500，其余参数使用默认值，具体设置方法可参见XGBoost文档[20]。根据以上参数在四个站点的训练集上训练各个站点的XGBoost模型，并在测试集上进行结果验证。

3 预测结果分析

福建省厦门市东临台湾海峡，与泉州、漳州和龙岩市毗邻。厦门市有溪东、洪文、鼓浪屿和湖里中学四个大气自动监测站，其中溪东是背景点，受城市污染干扰较少。

以上所提出的结合臭氧生成机理和机器学习方法的臭氧预测方法对厦门市8月31日至9月9日四个大气自动监测站的臭氧1小时平均浓度进行预测，实验结果如图4～图7所示，其中气象数据采用气象预报数据。

图4 溪东站臭氧浓度实测值、预测值对比

图5 洪文站臭氧浓度实测值、预测值对比

图6 鼓浪屿站臭氧浓度实测值、预测值对比

图7 湖里中学站臭氧浓度实测值、预测值对比

在图4～图7中，实线表示各站点实测值，虚线表示模型预报值。每个图的上半部分为臭氧1小时平均浓度值对比（一级限值是160µg/m3,二级限值是 200 µg/m3），下半部分为臭氧8小时滑动平均值对比（一级限值是100µg/m3,二级 限 值 是 160 µg/m3）。 可 以 看出，本文提出的方法对臭氧浓度的变化趋势捕捉较为准确，除了比较准确地表现臭氧浓度的日周期性变化，同时对峰值和低谷处能进行较为有效的捕捉和刻画。此外，该方法的预报浓度和等级准确率都较高。

臭氧浓度预测值与实测值对比如表5所示，其中对比了臭氧小时浓度的绝对误差最大值、绝对误差最小值、平均绝对误差和平均相对误差，臭氧八小时浓度滑动平均值的绝对误差最大值、绝对误差最小值、平均绝对误差和平均相对误差。

从表5中可以看出，除背景监测站点溪东外，其余三个大气监测站点的臭氧小时浓度的平均相对误差在21.8%～26.3%，平均绝对误差在14.21～15.97µg/m3；臭氧八小时浓度滑动平均值的平均相对误差为19.0%～23.4%，平均绝对误差为13.21～15.01µg/m3。四个站点臭氧小时浓度平均值和八小时浓度滑动平均值的绝对误差最小值均小于0.2µg/m3。以上结果显示，该方法有较高的预报准确率。此外，表5的结果显示溪东站点的平均相对误差较大，是因为溪东站点存在小时臭氧浓度极小的时刻，如9月9日1:00-7:00这7个时刻的臭氧浓度均在10µg/m3以下，因而导致该时刻的臭氧小时浓度平均相对误差极大。

除了对比臭氧小时浓度平均值和臭氧八小时浓度滑动平均值之外，根据臭氧日报等级计算方法，计算并对比了每天臭氧浓度八小时浓度滑动平均最大值计算对应的空气质量分指数ΙAQΙ及其等级，统计结果参见表6。四个站点均有1天日报等级错误，日报预报等级正确率较高，均为90%。

臭氧日报ΙAQΙ值及等级参见表7。溪东站点在9月9日等级预报错误，9月6日洪文、鼓浪屿和湖里中学三个站点日报等级均错误。除预报等级错误的天数外，其余日报臭氧ΙAQΙ预测值和真实值间的差值多在10以内。

表5 厦门四站点臭氧预测值与实测值对比

表6 厦门市四站点臭氧预测结果日报等级正确率

表7 厦门市四站点日报等级及预测值

通过厦门市四个大气监测点，可计算厦门市日报等级及预测值，参见表8，仅有9月6日一天日报等级错误，其余9天内7天日报ΙAQΙ差值都在5以内，说明等级预报较准，且日报ΙAQΙ差异较小。

表8 厦门市日报等级及预测值

通过以上结果可以看出，该方法也存在一些不足。由于目前对造成臭氧浓度变化的原因仅考虑点源的排放，针对其他原因导致的臭氧浓度突变可能捕捉不到，因此会出现如限行情况变化（8月31日至9月5日厦门机动车单双号限行，外地车辆不许进入岛内）和中元节期间（闽南地区祭祀祖宗习俗）等特殊原因导致9月6日臭氧浓度实测值较高但预测值较低的情况。等级预报方面，在日报等级接近级别边界的情况下，由于臭氧小时浓度的预报不够精准，在计算日报结果（臭氧八小时均值的最大值）时，可能出现等级预报错误的现象。这些都是该方法未来的改进方向。

4 结论

本文应用机器学习方法，较为准确地预测了大气中的臭氧浓度。该方法有别于数值模式的预报方法，可扩展性强，能够直接增加其他影响臭氧浓度变化的特征，并实现快速滚动预报。本文提出的预测方法充分考虑不同站点所处的不同地理位置及其可能受到的不同周边的影响，因此，可以方便推广到其他地区的臭氧浓度预测，但需采用该地区一年或一年以上的数据集进行学习训练。同时，该方法针对臭氧浓度突变值和等级边界的预报准确性还有待进一步提高。

[1] 唐孝炎. 大气环境化学[M]. 北京: 高等教育出版社, 1990: 60-70.

[2] 姜峰, 荀钰娴. 城市臭氧浓度变化规律分析[J]. 环境保护与循环经济, 2015, 35（2）: 55-59.

[3] 严茹莎, 李莉, 安静宇, 等. 上海市夏季臭氧生成与其前体物控制模拟研究[J]. 环境污染与防治, 2016, 38（1）: 30-35, 40-40.

[4] 王雪松, 李金龙. 北京地区臭氧源识别个例研究[J]. 北京大学学报（自然科学版）, 2003, 39（2）: 244-253.

[5] BΙNKOWSKΙ F S, ROSELLE S J. Models-3 community Multiscale Air Quality （CMAQ） model aerosol component 1. Model description[J]. Journal of geophysical research: atmospheres, 2003,108（D6）: 4183.

[6] TΙE X X, GENG F H, PENG L, et al. Measurement and modeling of O3variability in Shanghai, China: Application of the WRF-Chem model[J]. Atmospheric environment, 2009, 43（28）: 4289-4302.

[7] THOMPSON M L, REYNOLDS J, COX L H, et al. A review of statistical methods for the meteorological adjustment of tropospheric ozone[J]. Atmospheric environment, 2001, 35（3）: 617-630.

[8] LU W Z, WANG D. Ground-level ozone prediction by support vector machine approach with a cost-sensitive classification scheme[J].Science of the total environment, 2008, 395（2-3）: 109-116.

[9] COMAN A, ΙONESCU A, CANDAU Y. Hourly ozone prediction for a 24-h horizon using neural networks[J]. Environmental modelling& software, 2008, 23（12）: 1407-1421.

[10] YΙ J, PRYBUTOK V R. A neural network model forecasting for prediction of daily maximum ozone concentration in an industrialized urban area[J]. Environmental pollution, 1996, 92（3）: 349-357.

[11] PRYBUTOK V R, YΙ J, MΙTCHELL D. Comparison of neural network models with ARΙMA and regression models for prediction of Houston’s daily maximum ozone concentrations[J]. European journal of operational research, 2000, 122（1）: 31-40.

[12] 陈俏. 支持向量机应用于大气污染物浓度预测[D]. 西安：西安: 西安科技大学, 2010.

[13] ZHENG Y, YΙ X W, LΙ M, et al. Forecasting fine-grained air quality based on big data[C]//Proceedings of the 21th ACM SΙGKDD Ιnternational Conference on Knowledge Discovery and Data Mining. New York: ACM, 2015: 2267-2276.

[14] 刘烽, 徐怡珊. 臭氧数值预报模型综述[J]. 中国环境监测,2017, 33（4）: 1-16.

[15] 王占山, 李云婷, 陈添, 等. 北京城区臭氧日变化特征及与前体物的相关性分析[J]. 中国环境科学, 2014, 34（12）: 3001-3008.

[16] GRAŠΙč B, MLAKAR P, BOŽNAR M Z. Ozone prediction based on neural networks and Gaussian processes[J]. Nuovo Cimento della Societa Ιtaliana di Fisica Sect. C., 2006, 29（6）: 651-661.

[17] 王磊, 刘端阳, 韩桂荣, 等. 南京地区近地面臭氧浓度与气象条件关系研究[J/OL]. 环境科学学报: （2017-10-12）. http://kns.cnki.net/kcms/detail/11.1843.X.20171012.1702.002.html.

[18] CHEN T Q, GUESTRΙN C. XGBoost: A scalable tree boosting system[C]//Proceedings of the 22nd ACM SΙGKDD Ιnternational Conference on Knowledge Discovery and Data Mining. New York:ACM, 2016: 785-794.

[19] 叶倩怡, 饶泓, 姬名书. 基于Xgboost的商业销售预测[J]. 南昌大学学报（理科版）, 2017, 41（3）: 275-281.

[20] DMLC. Python APΙ Reference[EB/OL]. http://xgboost.readthedocs.io/en/latest/python/python_api.html.