Air quality prediction is a complex system engineering. How to fully consider the impact of meteorological, spatial and temporal factors on air quality is the core problem. To address this central conundrum, in an elaborate encoder-decoder architecture, we propose a new air quality prediction method based on multi-granularity spatiotemporal graph network. At the encoder, firstly, we use multi granularity graph and the well-known HYSPLIT model to build spatial relationship and dynamic edge relationship between nodes, respectively, while meteorological, temporal and topographic characteristics are used to build node features and LSTM (Long Short Term Memory) is used to learn the time-series relationship of pollutant concentration. At the decoder, secondly, we use the attention mechanism LSTM for decoding and forecasting of pollutant concentration. The proposed model is capable of tracking different influences on prediction resulting from the changes of air quality. On a project-based dataset, we validate the effectiveness of the proposed model and examine its abilities of capturing both fine-grained and long-term influences in pollutant process. We also compare the proposed model with the state-of-the-art air quality forecasting methods on the dataset of Yangtze River Delta city group, the experimental results show the appealing performance of our model over competitive baselines.
Air quality which is closely related to human public health, has been a common research hotspot focused by scholars all over the world. At present, many air quality monitoring stations (stations for short) have been built in major cities to monitor the concentration of air pollutants (PM2.5, PM10, O3, etc.) and meteorological parameters (temperature, pressure, wind speed, wind direction, humidity, etc.). However, these stations only can monitor real-time air quality, and fail to provide air quality prediction (AQP) and auxiliary support for urban intelligent decision-making or activity planning. How to construct an AQP model using a large amount of historical monitoring data has become a hot research topic in the field of environment engineering. Unfortunately, AQP is an extremely complex system engineering. On the one hand, air quality is related to pollutant emission, which is a type of time sequence and has periodicity; On the other hand, there exist physical and chemical changes of pollutants in the air, such as diffusion and deposition, which are greatly affected by meteorological and geographical locations; Finally, air quality also has certain probability, such as unexpected pollution leakage events will lead to a sharp decline in air quality.
The commonly used types of AQP include
mechanism model (MM) and machine learning (ML)
methods. The MM method [
With the development of deep learning, AQP
methods based on ML have attracted more attention [
In contrast, graph-based models naturally sidestep
the above issue since they shape the concentration
values into graph nodes and keep their original
distributions in graph structures. Because graph can
construct non-Euclidean entity distribution, it can
capture spatial relationships well. Thus, to compensate
for the lack of spatial relationship learning in the above
methods, the methods based on GNN (Graph Neural
Network) are applied to AQP. Qi et al. [
R=ra,1 a N as a and set
of N
we set regions ， L = la,1 a N as a location set of N regions ， Sa = sa,i ,1 i Sa as a set of stations in region ra ， a = la,i ,1 i Sa as a location set of stations in region ra . Where, la,i is composed of longitude and latitude of a station, and la is the mean value of la,i a .
features: quarter, month, week, hour and holiday. We
represent the quarter by using one-hot encoding
adopted for 4-bit binary representation, as shown in
Table 1.
of vehicle emissions by building a graph structure traffic
network, and then predict traffic pollution emissions.
Xu et al [
The above methods based on GNN use graph
structure to effectively construct the spatial
relationship among stations, but they fail to fully
construct the edges in the graph structure. For example,
literature [
In this paper, we use LSTM to learn the cycle and season of pollutant concentration, and use GNN to learn the spatial relationship of pollutants among stations, and thus building a multi granularity spatiotemporal graph neural network model, called MGST−GNN. In order to better capture the mutation of pollutant concentration and its impact on air diffusion and deposition, the mechanism model HYSLPLIT is used to dynamically construct the adjacency matrix and edge relation of the spatiotemporal graph.
The main contributions of this paper include: •We innovatively propose a dynamic spatiotemporal graph model combining mechanism model and graph neural network. The adjacency matrix and edge weight vector of dynamic graph are constructed based on the simulation results of diffusion, transport and deposition of polluted air mass by mechanism model, so that the architecture learns the spatial influence relationship among multi granularity stations.
•We propose to add time characteristic attributes of quarter, month, week, hour and holiday to each node in the encoder, use LSTM based on attention mechanism for temporal learning in the decoder to enhance MGST-GNN.
Similarly, the month is represented with 12-bit binary by one-hot encoding; week is represented with 7-bit binary by one-hot encoding; Hour is represented with 24-bit binary by one-hot coding; Holiday is represented with 2-bit binary by one-hot coding (0 means non-holidays, 1 means holiday).
Definition 3 Geomorphic Feature (GF): GF contains the topography and land usage information of a station, which consider the altitude and five land usage categories: residential area, park, mountain, water (river or pool), and industry. Where, the altitude is divided into four categories: very high (more than 1300m), high (1000-1300m), medium (500-1000m) and low (less than 500m), and the corresponding categories Spring season summertime
autumn wintertime
Spring summer fall
winter
1
0
0
0
0
1
0
0
0
0
1
0
0
0
0
1
are represented by the numbers 1, 2, 3 and 4; The land
usage type is determined by the number of major land
use types within 10km around the station. For example,
the altitude of sa,i is 800 meters, and there is one
residential area, two park, one mountain, three pool,
and two industrial facilities within the perception radius
10 kilometers of sa,i . Thus, the GF vector gfa,i can be
represented as [
ra and stations sa,i at time slot t is represented as a vector humidity, rainfall, wind speed, wind direction, and air pressure.
wdat and wdat ,i , including temperature,
pollutant concentration data includes the concentrations of six major pollutants such as PM2.5, PM10, SO2, CO2, CO and O3. Among them, PM2.5 and O3 are the most concerned at present. Therefore, the later experiments focus on the concentration forecasting of these two pollutants.
Definition 6 HYSPLIT: HYSPLIT [
In this paper, HYSPLIT is used to establish the source-receptor relationship among nodes by trajectory analysis. When HYSPLIT is used for trajectory analysis, it only needs to input the meteorological data of the simulation area and the coordinate information of the initial point of the simulation. The meteorological data can be downloaded from the official website of Air Resources Laboratory (ARL) 2 for free. Figure 1 describes an example of trajectory analysis using HYSPLIT. The example takes Beijing Center (116 ° 20 ', 39 °56') as the starting point and 10:00 on May 6, 2022 as the starting time to predict the air mass trajectory in the next 48 hours. It can be seen from the figure that the location and time of each track can be obtained by HYSPLIT. 2 https://nomads.ncep.noaa.gov/pub/data/nccf/com/hysplit
AQP task: Given region locations L , station locations , geomorphic feature gf a,i , in hours of pollutant concentration data con , in + out hours of weather data wd , and the AQP task aims to forecast the pollutant concentrations of stations for the next out hours, where in denotes the length of historical time window and out denotes forecasting horizon.
Existing AQP methods based on GNN [
We propose a multi-granularity spatiotemporal graph
network (MGST_GNN) by referring to reference [
Taking the AQP of all stations in a city as an example, all stations in the city are clustered according to their longitude and latitude information to obtain N categories, that is, N regions are clustered. Therefore, each region can be used as a node in the global graph, and each station in the region can be used as a node in the partial subgraph. A graph composed of region-level nodes is called a global graph, while a graph composed of station-level nodes is called a partial subgraph. From this, a multi-granularity graph network can be constructed.
Attribute definitions of graph node To the node in the global graph, its attribute is the pollutant concentration value of the node region. To the node in partial subgraph, its attribute is composed of the pollutant concentration value of the node station, GF and TF.
Edge and attribute definitions of graphs
Wind speed, wind direction, rainfall and other
meteorological data have a decisive influence on the
horizontal transmission of pollutants. To make use of
this domain knowledge, we adopt the air quality model
HYSPLIT to learn the relationship among nodes and
build dynamic connection edges. Taking the edge
calculation of global graph as an example, the city is
divided into grids according to the clustered regions, so
that each region falls into a unique grid. For each region
(node v a ), HYSPLIT is used to calculate all trajectories
starting from v a and stepping in hours within the
next time t. Track the time when each trajectories
crosses the grid where the node is located, and record
the number of pass trajectories and crossing time of the
grid from other nodes except v a , and dynamically
construct the global connection edge and its attribute
vector. Repeat these steps until each station is analyzed
as a starting point. For example, taking v1 as the
starting point, using HYSPLIT to calculate that there are
5 trajectories passing through v 2 in the following 48
hours, and the time interval is 1 hour, 3 hours, 4 hours,
8 hours and 48 hours in the future, then the edge
relationship between node v1 and v 2 is formed,
and its attribute vector is: e12
=[
MGST_GNN can model spatial dependencies with different granularities by message passing mechanism, which mainly includes
message aggregation and description updating. Message aggregation is shown as follows:
M a = ( xs , xa , ea,s )
sN(a) ma (M a ) (1) (2) where ， M a denotes the set of all the messages passed to node va ; xa is the attribute of node va ; xs is the attribute of a neighbor node; ea,s is the edge attribute; (a) denotes the neighbor node set of node va (Dynamically determined by
HYSPLIT trajectory)； ma is the aggregation vector of node va ; （） denotes aggregate function.
Description updating is shown below： xa 1 ( ma , xa ) , in global graph xa,i 2 ( ma,i , xa,i , ua ) , in partial graph (3) where，1and 2 denote update function; ma,i is the aggregation vector of node va,i in the a - th partial subgraph; xa,i is the attribute of node va,i in the a - th partial subgraph; ua is the global attribute of a - th
partial subgraph; the implemented using different FNNs. 1and 2 can be
The specific message aggregation method is from the partial subgraph to the global graph. The information transmitted from the partial subgraph includes the current and historical pollutant concentration values, and the average method is adopted to aggregate the information transmitted from the partial subgraph: conat = mean(conat,i ), 1 i Sa where, conat denotes the aggregate value at time t of node va
in the global graph; denotes the pollutant concentration value at time t of node va,i ; Sa
is the number of nodes of the a - th partial subgraph. To node in the global graph, message aggregation is used to calculate the global representation of pollutant concentration in each time slot, thus forming a series of sequences: cont− in+1, cont− in+2 ,, cont | 1 a N a a a These sequences are then fed into the global LSTM to learn the current and historical representations.
The specific description updating method is to update the partial subgraph using the global graph. The information transmitted from the global graph includes the historical pollutant concentration values of all (4) conat,i regions. At each time slot, the output of the global graph is fed into the FNN to obtain a downward update vector, which is used to update the global attribute u of the partial subgraph. The u include meteorological parameters and downward update vector, which are used in the partial subgraph by message passing. Therefore, the nodes in the station-level partial subgraph can make use of the historical information of neighbor nodes in the region-level global graph.
In order to capture the temporal dependence, we adopt an encoder-decoder architecture, as shown in figure 1. In each time slot, node attributes in the partial subgraph form a sequence Xa,i = xt−in +1, xt−in +2, , xat,i |1 a N,1 i Sa . That is, each a,i a,i time slot in the historical has a corresponding node attribute. The LSTM in the encoder uses X a,i as the input and the final state of the LSTM as the input of the decoder. The input of the LSTM in the decoder includes not only the output of the encoder, but also the node attribute in the partial subgraph. The output of LSTM in the decoder is used as the input of FNN, and FNN outputs the predicted pollutant concentration value in the future out time slot.
To make better use of temporal characteristics, temporal attention mechanism is introduced in the decoding stage to learn the dynamic temporal correlation between future time and historical time. Give the hidden state ht−1 and cell state ct−1 of the LSTM in the decoder at time t −1 , then at time t , the attention weight of the hidden state ht output by the encoder is calculated as follows: tt = vT tanh (W h +Uht + b ) t−1; ct−1 tt =
exp ( tt ) exp ( tt ) (5) (6) t where， t is the attention weight; v,b,W ,U are the parameters to be learned. Through the above formula (5) and (6), the attention weight of all historical hidden states in the encoder can be calculated, and then, the hidden state ht is weighted and summed to obtain the time context vector c ： ct = Tt=1ttht (7) The output result ˆ ot −1 at time t −1 of decoder, the meteorological data wdt at time t , the time feature tf t and the time context vector c t are connected as the input for the LSTM of the decoder at time t , and it is used to update the hidden state ht : ht = LSTM ( ht−1,oˆt−1;wdt ; tf t ; ct ) (8)
Jinan database (JN)：JN is a city-level data set collected by us in projects. Jinan is located in the middle of Shandong Province, China. There are 130 air monitoring stations in Jinan. Each station outputs the concentration values of pollutants (PM2.5, PM10, SO2, CO2, Co, O3) and meteorological parameters (rainfall, surface pressure, temperature, humidity, wind speed and wind direction) every hour. We collected the historical monitoring data of 130 stations in Jinan from January 1st, 2019, to January 1st, 2022, as the training and test set.
We divide the 130 stations into 13 regions, and the global graph consists of 13 regions. Each region is used as a global node in the global graph. The stations in each region form a partial subgraph, and the stations in the region are the nodes of the partial subgraph.
Yangtze River Delta city group database (YRD): The city group contains ten cities: Shanghai, Hangzhou, Suzhou, Ningbo, Shaoxing, Jiaxing, Wuxi, Zhoushan, Nantong, and Huzhou. We used air pollution prediction system3 to collect historical pollutant concentration values and meteorological parameters of corresponding stations, and the time span was from January 1th, 2019, to December 31th, 2022. Therefore, we collected 3 years of historical monitoring data as the training and test set. Each city in the Yangtze River Delta city group 3 http://airprediction.urban-computing.com is a global graph node, and each station in the city is a node of the corresponding partial subgraph.
Geographic features and weather forecast data are obtained by: Geographic features are collected from the map engine of AMAP4.The perception radius is set to 1000 m.
Weather forecast data are collected from the Air Resources Laboratory （ ARL ） 5 .This website can download the meteorological forecast data for the following 26 days at most, with an accuracy of 0.25 0.25 degrees, which can be updated four times a day.
We split the dataset into training data, validation data,
and test data by the radio of 0.7:0.1:0.2. We choose
Adam [
We implement our method by PyTorch [
We introduce two metrics: mean absolute error (MAE) and symmetric mean absolute percentage error (SMAPE) to evaluate the performances of methods. In the experiments, we utilize previous 48-hour observations and select the results of 1 hour, 6 hours, 12 hours, 18 hours, 24 hours, 36 hours and 48 hours ahead forecasting to report. All experiments are repeated 5 times to avoid contingency. Figure 4 depicts the predicted PM 2.5 and O3 concentrations and their corresponding real values of a station in Jinan city from 00:00 on January 2nd, 2021 to 23:00 on January 3rd, 2021 for consecutive 48 hours. From figure 4, we can see that it is easier to predict O3 concentration than PM2.5, because O3 has more periodic regularity and stability than PM2.5. Therefore, the subsequent experiments were carried out with the forecasting of PM2.5 concentration. 4 https://lbs.amap.com/api/webservice/guide/api/search/ 5 https://nomads.ncep.noaa.gov/pub/data/
To verify the effectiveness of multiple factors, we compare MGST_GNN with three variants, each of which removes one kind of factors. Specifically, MGST_GNN w/o wdf removes meteorological forecast data; tf MGST_GNN w/o removes time feature and MGST_GNN w/o gf removes geomorphic feature.
The performances of MGST_GNN and its variants are given in Table 2, as can be seen from the table:
MGST_GNN outperforms the other three variants, indicating that all factors can improve the performance of AQP. The rank of the effectiveness of factors is WD > TF> GF. The result shows that air quality is mostly impacted by weather conditions.
To explore the effectiveness of different components, we compare MGST_GNN with three variants:
MGST_GNN w/o multi granularity, which removes the global graph and the corresponding interactions, i.e., the region representation is removed from the global attributes of the partial subgraphs;
MGST_GNN w/o HYSPLIT, which removes dynamic edge weight vector by HYSPLIT and use vector of wind direction and distance instead;
MGST_GNN w/o taLSTM, which removes temporal attention mechanism based LSTM at the decoder, and directly use LSTM instead.
The performances of MGST_GNN and its variants are given in Table 3, as can be seen from the table:
MGST_GNN outperforms MGST_GNN w/o multi granularity in all metrics. The result indicates that the air quality of adjacent regions is beneficial, which can be used to model the diffusion processes of air pollutants from adjacent regions.
MGST_GNN outperforms MGST_GNN w/o HYSPLIT. The result indicates that compared with taking wind direction and distance as the edge weight vector, using HYSPLIT to dynamically adjust the weights of edges is a
the results of factor evaluation based on JN database method Metric 1h 6h MGST_GNN SMMAAEPE 60..1057 10.51.204 w/o wdf SMMAAEPE 07..0076 01.71.117 w/o tf MAE 6.76 17.16
SMAPE 0.07 0.11 w/o gf SMMAAEPE 06..0771 01.61.178
To verify the advanced nature of our method, we
compare MGST_GNN with the methods HighAir[
The performances of MGST_GNN and other prediction methods are given in Table 4 and 5, as can be seen from the table:
MGST_GNN outperforms GC-DCRNN, GC-LSTM and ST_UNet, especially in long-term forecasting. It indicates that a multi-granularity structure can model spatial dependencies more effectively than a flat structure. The reason is that the multi granularity graph not only considers the local impact of neighboring stations on the prediction station, but also the global impact of different regions on the prediction station, which makes the mining of spatial relationships more sufficient.
MGST_GNN outperforms HighAir. This is because in the aspect of spatial relationship learning, MGST_GNN uses the professional HYSPLIT model to build edge
Model comparison results based on JN database weight vector. It comprehensively uses meteorological and topographical conditions to calculate the influence relationship among nodes and the specific influence time, making the construction of spatial relationship more accurate and delicate. HighAir simply uses wind direction and distance to construct rough edge weight vector. In terms of temporal relation learning, MGST_GNN adds time feature, and uses the attention mechanism based LSTM in decoder. HighAir simply uses LSTM to build temporal relationships.
MGST_GNN outperforms PM2.5-GNN. This is because MGST_GNN uses a multi-granularity spatiotemporal graph network and combines meteorological data, topographic features, time features and professional models to construct node and edge attributes. However, PM2.5-GNN is only a single-granularity spatiotemporal graph network. In the construction of edge attributes, although they make use of domain knowledge, they simply list part of the influencing factors, which is still insufficient. When building node attributes, they only use meteorological data. method MGST_GNN HighAir PM2.5-GNN GC-DCRNN GC-LSTM ST-UNet
Metric MAE SMAPE MAE SMAPE MAE SMAPE MAE SMAPE MAE SMAPE MAE SMAPE 1h
MAE SMAPE
Model comparison results based on YRD database method Metric 1h 6h MGST_GNN MAE 7.16 17.29 SMAPE 0.08 0.12
MAE 7.12 17.83 HighAir SMAPE 0.08 0.12
MAE 7.15 18.09 PM2.5-GNN SMAPE 0.08 0.13
MAE 7.56 18.56 GC-DCRNN
SMAPE 0.08 0.13
MAE 8.13 18.97 GC-LSTM SMAPE 0.09 0.13
MAE 7.22 18.69 ST-UNet SMAPE 0.08 0.13
MAE 8.15 19.21 STA-LSTM SMAPE 0.09 0.14
In this paper, the influence of meteorological, spatial and temporal factors on AQP is fully considered, and an encoder-decoder architecture based on multi-granularity spatiotemporal graph network is proposed to predict pollutant concentration over a long period of time. Compared with the existing models, the striking characteristic of this paper is that the meteorological, spatial terrain and time factors are considered comprehensively through the professional air quality model, while other models take the influencing factors as the splitting parameter input.
That is, the model in this paper integrates the advantages of mechanism model and machine learning, and namely it is a comprehensive model. The experimental results show that the proposed model is of progressiveness and has good applicability.
The authors gratefully acknowledge the NOAA Air Resources Laboratory (ARL) for the provision of the HYSPLIT transport and dispersion model and/or READY website (https://www.ready.noaa.gov) used in this publication.
A.1 DATA AUGMENTATION AND EXPERIMENTATION Aiming at the difficulty of collecting training samples for air quality prediction, we use the method of data augmentation to reconstruct samples. Due to the distribution difference between real weather data and forecast weather data, Gaussian noise is introduced to the meteorological data in the sample to enhance the data and improve the generalization ability of the model. In addition, due to the influence of monitoring instruments, environment and other factors, the pollutant concentration values collected may be biased. Therefore, the pollutant concentration values of samples are disturbed up and down by 1 metric to enhance the number of samples and improve the model generalization ability. Therefore, Gaussian noise is introduced into meteorological data and random perturbation is for pollutant concentration data to enhance the samples by three times. The specific augmentation methods are shown in Table 5. Finally, 26136 samples were obtained from JN and YRD respectively.
In order to verify the data augmentation effect, we conduct a comparison experiment between data augmentation and non-data augmentation, as shown in Table 6, where w/o DA indicates that data augmentation technology is not used. As seen from Table 6: using our data augmentation method, the effect is effectively improved.
A.2 ABLATION STUDY AND EXPERIMENTATION In this paper, we are the first to use the professional model HYSPLIT to build the graph dynamically. Therefore, this appendix section will demonstrate the effectiveness of using HYSPLIT to construction graph dynamically through ablation experiments. Specifically, the ablation experiment was conducted based on the current advanced spatiotemporal graph neural network HighAir and PM2.5-GNN. In HighAir, it uses the distance among nodes to statically construct the edge of the graph, and uses the wind direction information between nodes to calculate the edge attribute vector. Therefore, we use HYSPLIT instead of the graph construction method in HighAir. In PM2.5-GNN, it uses the distance between nodes and the altitude of the position to statically construct the edge of the graph, and uses the parameter of domain knowledge between nodes to calculate the edge attribute vector. Therefore, we use HYSPLIT instead of the graph construction method in PM2.5-GNN. The experimental results are shown in Table 7, where HighAir_HYSPLIT and PM2.5-GNN_HYSPLIT respectively indicate that HYSPLIT is used to replace the original graph construction. As can be seen from the table, HYSPLIT builds dynamic graphs better than static graphs of HighAir and PM2.5-GNN.
In order to further verify the effectiveness of using HYSPLIT to dynamically construct graph, especially its advantages for air quality prediction in complex scenarios such as abrupt change in pollutant concentration. We develop a Dataset-mini, where we focus on heating season (November to February). Dataset-mini is more challenging for two reasons. Firstly, during winters, heating emissions can dramatically increase the frequency of air pollution occurrence. Secondly, the direction of prevailing wind is north or northwest, which contributes to pollutant’s long-distance transport from North China to South China. The results in Table 9 show that using HYSPLIT to dynamically construct graph can significantly improve the accuracy of the model's air quality prediction on sudden changes in pollution and regional impacts caused by strong winds. This method eliminates the construction process of specially designed auxiliary network to learn the edges of graph and provides a new way for the construction of graph neural network, which can be easily extended to other spatiotemporal forecasting tasks. For example, in water quality prediction, a professional hydrodynamic model (MIKE) can be used to dynamically construct the graph structure, so as to better learn the influence of water quality in different regions on the prediction points.