運用機器學習方法預測空氣中臭氧濃度

蔡旺華

（ 福建省環(huán)境信息中心，福建福州 350003 ）

臭氧（O3）是一種淡藍色氣體，是天然大氣的重要微量成分,約90%的臭氧存在于平流層,僅有10%左右的臭氧分布在對流層中[1]。臭氧是《環(huán)境空氣質(zhì)量標準》（GB3095—2012）中包含光化學煙霧的標志性污染物，高濃度近地面臭氧對人體健康有較大的危害[2]。城市中的臭氧主要由氮氧化物（NOx）、揮發(fā)性有機物（VOCs）在適宜條件下作用產(chǎn)生[3],主要來源于機動車、各類工廠和自然植被。隨著經(jīng)濟活動的加劇和機動車保有量的激增，臭氧污染問題漸趨嚴重[4]，2017年福建省空氣質(zhì)量的首要污染源即為臭氧。如何準確預報臭氧濃度成為一個重要課題。

常見的大氣污染物濃度預測方法主要有兩類：基于數(shù)值模式的預測方法[5,6]和基于統(tǒng)計模型的預測方法[7-10]。數(shù)值模式的預測方法是以大氣動力學理論為基礎，在給定的氣象條件、污染源排放清單以及初始邊界條件下，通過一套復雜的偏微分方程組模擬大氣污染物在實際大氣中的各種物理化學過程，預報污染物濃度動態(tài)分布和變化趨勢?；诮y(tǒng)計模式的預測方法以空氣質(zhì)量監(jiān)測數(shù)據(jù)、氣象監(jiān)測數(shù)據(jù)等多種類型數(shù)據(jù)為基礎，利用統(tǒng)計方法，建立預測模型，實現(xiàn)大氣污染物濃度預測。傳統(tǒng)的統(tǒng)計模型預測僅考慮單一大氣污染物濃度變化，利用ARΙMA、Holt-Winters等時間序列預測模型進行濃度預測[11]，預測精度較低且預測時效較短。后有方法提出利用大氣污染物濃度和溫度、濕度等氣象因素進行多元線性回歸，預測精度相比ARΙMA等時間序列預測模型要準。近年來，隨著機器學習算法的不斷發(fā)展，支持向量機、人工神經(jīng)網(wǎng)絡等機器學習算法被不斷引入大氣污染物濃度預測任務中來[12,13]，這些方法對多源數(shù)據(jù)進行綜合考慮，預測效果有所提升。

相比基于數(shù)值模式的預測方法，基于統(tǒng)計模型的預測方法具有以下兩個優(yōu)點：一是此類方法的輸入采用真實監(jiān)測值，比采用估算值的數(shù)值模式更加準確[14]；二是此類方法計算速度快，通常在秒級即可完成計算，而數(shù)值模式運算消耗時間動輒以小時計。但是，目前國內(nèi)常見的基于統(tǒng)計模型的大氣污染物濃度預測方法一般只對空氣質(zhì)量監(jiān)測數(shù)據(jù)、氣象數(shù)據(jù)進行回歸分析。由于輸入特征單一，未考慮到引起大氣污染物濃度變化的可能原因，如點源、移動源排放變化等，從而導致了現(xiàn)有此類方法往往表現(xiàn)得不夠理想。

本文提出了一種結(jié)合臭氧濃度變化機理的機器學習方法用于臭氧濃度預測，是福建省生態(tài)環(huán)境大數(shù)據(jù)平臺建設中的大數(shù)據(jù)應用之一。福建省生態(tài)環(huán)境大數(shù)據(jù)平臺由福建省環(huán)保廳于2015年9月啟動建設，在建設過程中，正值2017年9月3日至5日金磚國家領導人廈門會晤。為此，在8月31日至9月9日，運用本文提出的方法對廈門市臭氧濃度進行了較為準確的預測，為保障空氣質(zhì)量提供了支持。

1 臭氧濃度預測模型

1.1 臭氧預測模型機制

臭氧濃度變化可以簡單分為生成、分解和擴散三部分[15,16]。其中，生成部分主要考慮三類污染源（天然源、移動源和點源）的影響。因為在一段時期內(nèi)，可以認為天然源的影響基本穩(wěn)定，無人工干預下的移動源的影響呈現(xiàn)較規(guī)律的周期性變化，所以此方法主要考慮相對變化較大的點源對臭氧濃度變化的影響，從而進行臭氧濃度預測；分解和擴散的部分主要考慮濕度、風速等氣象條件對臭氧濃度變化的影響[17]。因此本方法采用的數(shù)據(jù)主要有以下四類：①空氣質(zhì)量監(jiān)測數(shù)據(jù)；②周邊點源（以下本文所提“點源”均指有在線監(jiān)測廢氣排放的污染源）排放數(shù)據(jù)；③氣象監(jiān)測數(shù)據(jù)；④預報時刻的氣象預測數(shù)據(jù)。根據(jù)各類數(shù)據(jù)對臭氧濃度變化的影響構(gòu)造特征，預測24小時后的臭氧濃度。

因為各個大氣監(jiān)測站所處位置不同，各種條件一般會有較大差異，臭氧濃度變化規(guī)律也可能不同，所以對每個大氣自動監(jiān)測站單獨訓練臭氧濃度預測模型。

該方法的流程如圖1所示。首先獲取上文提及的空氣質(zhì)量監(jiān)測數(shù)據(jù)、周邊點源排放數(shù)據(jù)、氣象監(jiān)測數(shù)據(jù)和預報時刻的氣象數(shù)據(jù)等四類數(shù)據(jù)，然后進行數(shù)據(jù)預處理、特征抽取步驟，獲取訓練數(shù)據(jù)集和測試數(shù)據(jù)集。然后利用訓練數(shù)據(jù)集和機器學習的XGBoost算法進行模型訓練，得到預測模型，再利用預測模型和測試數(shù)據(jù)集得到預測臭氧濃度。

圖1 臭氧濃度預測模型流程

圖2 臭氧濃度預測模型架構(gòu)

臭氧濃度預測模型架構(gòu)如圖2所示，對圖1的特征抽取和模型訓練兩部分進行展開。由圖2可見，模型輸入特征基于起報時刻t前h天的點源數(shù)據(jù)（圖1中h=5，即t前120小時）、臭氧濃度監(jiān)測數(shù)據(jù)和氣象監(jiān)測數(shù)據(jù)以及預報時刻氣象數(shù)據(jù)構(gòu)造。模型目標為預報時刻（t+24）的臭氧濃度相較于起報時刻（t）臭氧濃度的變化值。采用變化值為目標是為了減小臭氧濃度周期性變化以及一些潛在未知因素的影響?；赬GBoost算法，可以獲得預報時刻相較起報時刻的臭氧濃度變化預測值，再結(jié)合起報時刻的臭氧濃度，即可獲得預報時刻的臭氧濃度預測值。

1.2 特征抽取

本文提出的算法包括三部分特征：點源特征、臭氧濃度和氣象特征，以下依次進行說明?？紤]到各個大氣監(jiān)測站點所處的位置不同，該算法針對單個大氣監(jiān)測站點設計。

點源特征按照點源地理位置分布抽取，如圖3所示，五星點為目標大氣監(jiān)測站點，黑色圓點為點源。針對每個大氣監(jiān)測站點，考慮與其所屬城市接壤的所有城市內(nèi)的點源，構(gòu)造如圖2的區(qū)域，并將其劃分為若干個長寬均為M的小方格。針對每個小方格，計算其中心點和當前環(huán)境監(jiān)測站點的相對方向和距離（圖2中所示的叉劃線），同時統(tǒng)計其所包含的點源，則所有點源均明確隸屬于某個方格。對前述已劃定的每個小方格，計算起報時刻t前h天內(nèi)的機制特征，以每24小時為1個時間段，每個時間段抽取二維特征，分別是污染物擴散規(guī)律特征和污染物傳播時間特征。對每個小方格每個時間段，若主要風向有利于污染物到當前關(guān)注的環(huán)境監(jiān)測站點的輸送，即可以根據(jù)平均風速和點源排放量計算該二維特征。當風向不易于污染物輸送或者方格內(nèi)暫無點源排放時，該二維特征均為0。

圖3 點源網(wǎng)格劃分特征抽取示意圖

臭氧濃度特征采用起報時刻t的臭氧小時平均濃度。

氣象特征通過預報時刻t+24前若干個小時的氣象因素構(gòu)造，因為濕度和風速對臭氧的擴散和分解至關(guān)重要，所以本算法中采用其預測時刻前4小時、前6小時和前24小時的最大值、最小值和平均值作為特征。

將上述三種特征直接拼接構(gòu)成輸入特征，對應的輸出結(jié)果為預測時刻t+24與起報時刻t的臭氧濃度差值，之后利用XGBoost算法擬合針對每個大氣監(jiān)測站點的模型即可。

1.3 XGBoost算法

本文采用的XGBoost[18]（eXtreme Gradient Boosting）是一種梯度提升算法擴展。梯度提升算法是一種用于回歸和分類問題的機器學習技術(shù)，以多個弱預測模型集成的形式產(chǎn)生預測模型。通常會選擇樹這一弱預測模型作為基學習器，例如，本文所提出預測模型采用回歸樹作為基學習器。當然，XGBoost基學習器也可以是其他分類器。

假設XGBoost一共生成了k棵樹，對任一輸入x，每棵樹都有一個輸出fx(x)，則該XGBoost模型的輸出為在迭代過程中，利用梯度下降思想，以已生成樹為基礎，以最優(yōu)化目標函數(shù)為總目標，迭代生成新樹。每次迭代學習的目標是減少已生成的模型的累計結(jié)果的損失值。不同于傳統(tǒng)梯度提升決策樹（gradient boosting decision tree，GBDT），XGBoost對 代價函數(shù)進行二階泰勒展開，利用了一階導數(shù)和二階導數(shù)。另外，XGBoost的代價函數(shù)中包含正則項，這相當于給復雜模型以一定的懲罰，模型復雜度得以控制，所以模型不易出現(xiàn)過擬合問題。

XGBoost算法能夠利用CPU多線程并行，具有運行速度快，準確度較高，不易過擬合等優(yōu)點[19]，算法普適性好，在多種應用上表現(xiàn)良好，但未見在臭氧濃度預測的文獻。

2 模型預測

2.1 數(shù)據(jù)準備

該方法采用的各類數(shù)據(jù)源，除測試過程所需的氣象預報數(shù)據(jù)外，均來源于福建省生態(tài)環(huán)境大數(shù)據(jù)平臺。該方法涉及的數(shù)據(jù)源包括以下四類：大氣監(jiān)測站的空氣質(zhì)量監(jiān)測數(shù)據(jù)、大氣監(jiān)測站的氣象監(jiān)測數(shù)據(jù)、點源（在線監(jiān)測的廢氣排放污染源）數(shù)據(jù)和氣象預報數(shù)據(jù)。大氣監(jiān)測站的空氣質(zhì)量監(jiān)測數(shù)據(jù)格式參見表1，以小時平均深度代表每個大氣監(jiān)測站點在每個整點時刻的臭氧濃度；大氣監(jiān)測站的氣象監(jiān)測數(shù)據(jù)和氣象預報數(shù)據(jù)的數(shù)據(jù)格式相同，參見表2，記錄每個大氣監(jiān)測站點每個時刻的溫度、濕度、氣壓、風速和風向等5個監(jiān)測指標；點源監(jiān)測數(shù)據(jù)格式參見表3，記錄每個點源小時排放監(jiān)測值的NOx折算濃度。

表1 大氣監(jiān)測站的空氣質(zhì)量監(jiān)測數(shù)據(jù)

表2 大氣監(jiān)測站的氣象監(jiān)測數(shù)據(jù)

表3 點源（在線監(jiān)測的廢氣排放污染源）數(shù)據(jù)

該方法采用2016年9月1日至2017年8月30日一年的上述三種數(shù)據(jù)源，針對廈門市溪東、洪文、鼓浪嶼和湖里中學四個大氣監(jiān)測站點各自的訓練模型。利用金磚領導人廈門會晤前后共10天（2017年8月31日至9月9日）的數(shù)據(jù)和對應的氣象預報數(shù)據(jù)進行模型測試。訓練數(shù)據(jù)和測試數(shù)據(jù)量如表4所示，因數(shù)據(jù)缺失或異常，各站點的訓練、測試數(shù)據(jù)條數(shù)可能不同。

表4 各站點數(shù)據(jù)量統(tǒng)計

2.2 模型訓練和參數(shù)設置

對每個大氣監(jiān)測站點，在訓練集上使用10折交叉驗證方法進行模型訓練和參數(shù)調(diào)優(yōu)，最后在測試集上測試結(jié)果。所謂10折交叉驗證，即將數(shù)據(jù)集劃分為10個子集，依次取其中9個子集為訓練集，剩余1個子集為測試集，最終以在10個測試集上的平均精度評價模型。

本次實驗在訓練XGBoost模型時設置了三個參數(shù)：樹的深度（max_depth）為4、學習率（learning_rate）為0.1和樹的個數(shù)（n_estimators）為500，其余參數(shù)使用默認值，具體設置方法可參見XGBoost文檔[20]。根據(jù)以上參數(shù)在四個站點的訓練集上訓練各個站點的XGBoost模型，并在測試集上進行結(jié)果驗證。

3 預測結(jié)果分析

福建省廈門市東臨臺灣海峽，與泉州、漳州和龍巖市毗鄰。廈門市有溪東、洪文、鼓浪嶼和湖里中學四個大氣自動監(jiān)測站，其中溪東是背景點，受城市污染干擾較少。

以上所提出的結(jié)合臭氧生成機理和機器學習方法的臭氧預測方法對廈門市8月31日至9月9日四個大氣自動監(jiān)測站的臭氧1小時平均濃度進行預測，實驗結(jié)果如圖4～圖7所示，其中氣象數(shù)據(jù)采用氣象預報數(shù)據(jù)。

圖4 溪東站臭氧濃度實測值、預測值對比

圖5 洪文站臭氧濃度實測值、預測值對比

圖6 鼓浪嶼站臭氧濃度實測值、預測值對比

圖7 湖里中學站臭氧濃度實測值、預測值對比

在圖4～圖7中，實線表示各站點實測值，虛線表示模型預報值。每個圖的上半部分為臭氧1小時平均濃度值對比（一級限值是160μg/m3,二級限值是 200 μg/m3），下半部分為臭氧8小時滑動平均值對比（一級限值是100μg/m3,二級 限 值 是 160 μg/m3）。 可 以 看出，本文提出的方法對臭氧濃度的變化趨勢捕捉較為準確，除了比較準確地表現(xiàn)臭氧濃度的日周期性變化，同時對峰值和低谷處能進行較為有效的捕捉和刻畫。此外，該方法的預報濃度和等級準確率都較高。

臭氧濃度預測值與實測值對比如表5所示，其中對比了臭氧小時濃度的絕對誤差最大值、絕對誤差最小值、平均絕對誤差和平均相對誤差，臭氧八小時濃度滑動平均值的絕對誤差最大值、絕對誤差最小值、平均絕對誤差和平均相對誤差。

從表5中可以看出，除背景監(jiān)測站點溪東外，其余三個大氣監(jiān)測站點的臭氧小時濃度的平均相對誤差在21.8%～26.3%，平均絕對誤差在14.21～15.97μg/m3；臭氧八小時濃度滑動平均值的平均相對誤差為19.0%～23.4%，平均絕對誤差為13.21～15.01μg/m3。四個站點臭氧小時濃度平均值和八小時濃度滑動平均值的絕對誤差最小值均小于0.2μg/m3。以上結(jié)果顯示，該方法有較高的預報準確率。此外，表5的結(jié)果顯示溪東站點的平均相對誤差較大，是因為溪東站點存在小時臭氧濃度極小的時刻，如9月9日1:00-7:00這7個時刻的臭氧濃度均在10μg/m3以下，因而導致該時刻的臭氧小時濃度平均相對誤差極大。

除了對比臭氧小時濃度平均值和臭氧八小時濃度滑動平均值之外，根據(jù)臭氧日報等級計算方法，計算并對比了每天臭氧濃度八小時濃度滑動平均最大值計算對應的空氣質(zhì)量分指數(shù)ΙAQΙ及其等級，統(tǒng)計結(jié)果參見表6。四個站點均有1天日報等級錯誤，日報預報等級正確率較高，均為90%。

臭氧日報ΙAQΙ值及等級參見表7。溪東站點在9月9日等級預報錯誤，9月6日洪文、鼓浪嶼和湖里中學三個站點日報等級均錯誤。除預報等級錯誤的天數(shù)外，其余日報臭氧ΙAQΙ預測值和真實值間的差值多在10以內(nèi)。

表5 廈門四站點臭氧預測值與實測值對比

表6 廈門市四站點臭氧預測結(jié)果日報等級正確率

表7 廈門市四站點日報等級及預測值

通過廈門市四個大氣監(jiān)測點，可計算廈門市日報等級及預測值，參見表8，僅有9月6日一天日報等級錯誤，其余9天內(nèi)7天日報ΙAQΙ差值都在5以內(nèi)，說明等級預報較準，且日報ΙAQΙ差異較小。

表8 廈門市日報等級及預測值

通過以上結(jié)果可以看出，該方法也存在一些不足。由于目前對造成臭氧濃度變化的原因僅考慮點源的排放，針對其他原因?qū)е碌某粞鯘舛韧蛔兛赡懿蹲讲坏?，因此會出現(xiàn)如限行情況變化（8月31日至9月5日廈門機動車單雙號限行，外地車輛不許進入島內(nèi)）和中元節(jié)期間（閩南地區(qū)祭祀祖宗習俗）等特殊原因?qū)е?月6日臭氧濃度實測值較高但預測值較低的情況。等級預報方面，在日報等級接近級別邊界的情況下，由于臭氧小時濃度的預報不夠精準，在計算日報結(jié)果（臭氧八小時均值的最大值）時，可能出現(xiàn)等級預報錯誤的現(xiàn)象。這些都是該方法未來的改進方向。

4 結(jié)論

本文應用機器學習方法，較為準確地預測了大氣中的臭氧濃度。該方法有別于數(shù)值模式的預報方法，可擴展性強，能夠直接增加其他影響臭氧濃度變化的特征，并實現(xiàn)快速滾動預報。本文提出的預測方法充分考慮不同站點所處的不同地理位置及其可能受到的不同周邊的影響，因此，可以方便推廣到其他地區(qū)的臭氧濃度預測，但需采用該地區(qū)一年或一年以上的數(shù)據(jù)集進行學習訓練。同時，該方法針對臭氧濃度突變值和等級邊界的預報準確性還有待進一步提高。

[1] 唐孝炎. 大氣環(huán)境化學[M]. 北京: 高等教育出版社, 1990: 60-70.

[2] 姜峰, 荀鈺嫻. 城市臭氧濃度變化規(guī)律分析[J]. 環(huán)境保護與循環(huán)經(jīng)濟, 2015, 35（2）: 55-59.

[3] 嚴茹莎, 李莉, 安靜宇, 等. 上海市夏季臭氧生成與其前體物控制模擬研究[J]. 環(huán)境污染與防治, 2016, 38（1）: 30-35, 40-40.

[4] 王雪松, 李金龍. 北京地區(qū)臭氧源識別個例研究[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] 劉烽, 徐怡珊. 臭氧數(shù)值預報模型綜述[J]. 中國環(huán)境監(jiān)測,2017, 33（4）: 1-16.

[15] 王占山, 李云婷, 陳添, 等. 北京城區(qū)臭氧日變化特征及與前體物的相關(guān)性分析[J]. 中國環(huán)境科學, 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] 王磊, 劉端陽, 韓桂榮, 等. 南京地區(qū)近地面臭氧濃度與氣象條件關(guān)系研究[J/OL]. 環(huán)境科學學報: （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的商業(yè)銷售預測[J]. 南昌大學學報（理科版）, 2017, 41（3）: 275-281.

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