Downscaling of Regional Air Quality Model using Gaussian Plume Model and Random Forest Regression

1. Introduction

$P{M}_{10}$ is a major air pollutant that can significantly impact human health, including respiratory and cardiovascular diseases [1]. Therefore, accurate and detailed information on $P{M}_{10}$ concentrations is essential for assessing and managing the impacts of air pollution.

In recent years, advances in remote sensing and satellite data have enabled high-resolution $P{M}_{10}$ concentration maps at regional and global scales [2,3]. These maps can provide valuable information to policymakers, researchers, and the general public. They help to plan and evaluate efforts to improve air quality and protect human health. However, due to inherent remote sensing limitations, they remain a temporal snapshot (or a mosaic of snapshots), sometimes misinterpreted as monthly or annual average concentrations. Still, the demand for high-resolution concentration maps needs additional means [4].

Using mathematical equations, air quality models can simulate the transport and dispersion of $P{M}_{10}$ particles in the atmosphere. These models use meteorological data, emission inventories, and other inputs to estimate $P{M}_{10}$ concentrations at high spatial and temporal resolution. This method can generate $P{M}_{10}$ concentration maps from regional to global scales. While air quality models can be accurate, they require significant computational resources and depend on the quality of input data. There are several classes of air quality models[5], including CTM (Chemical Transport Models), which are based on mass and momentum conservation equations and are the most universal and general. Since pollutant transport is highly linked to meteorological conditions, CTMs are often run together with meteorological models (online approach) or demand the results of a meteorological model (offline approach).

Gaussian plume models are analytical and widely used due to their low computational demands [6]. These models provide an analytical solution to the pollutant transport equation, assuming one of the atmosphere stability classes. They are usually applied for a single stack to assess its environmental impact [7].

Another rapidly emerging approach for obtaining air quality concentration maps is the application of machine learning techniques (also known as data-driven modelling) [8]. In this case, the user is only responsible for defining the input data set (like meteorological parameters, land use, and emission inventory data) and target variables (usually concentrations). These models act as black-box models, trained using supervised learning processes. Depending on the output type - they can be classified as classification or regression models. The most popular ones include neural networks (focusing on LSTM networks [9]), random forests [10] or spatial kriging algorithms [11].

Data fusion methods can combine data from different sources [12], such as ground-based monitoring, satellite remote sensing, and air quality models, to generate high-resolution $P{M}_{10}$ concentration maps. These methods use statistical and machine learning techniques to merge the data and estimate $P{M}_{10}$ concentrations at locations where data are missing or incomplete. While data fusion methods can improve the accuracy and spatial resolution of $P{M}_{10}$ concentration maps, benefiting from multiple approaches, they must be applied cautiously as they often produce non-physical results [13,14].

This paper uses Random Forest to describe and evaluate a data fusion technique from two sources (the regional GEM-AQ model[15] and the local Gaussian plume model). We discuss several approaches to using random forests in a diverse area of southern Poland. The proposed approach was tested on a regional domain, routinely modelled with a coarse-resolution regional model. The proposed approach can also act as a form of air quality model results downscaling. This region was chosen due to data availability and reported air-quality issues in the past. Observation and emission data from the while year 2021 were used.

2. Data and Methods

In this study, we used the GEM-AQ model 24 h forecast from the operational run and the national air quality monitoring network observations. The study period covered the year 2021.

2.1. Study Area

The study area is located in southern Poland, covering an area of around 5300 $k{m}^2$ in Silesian and Lesser Poland voivodships (Figure 1A). The area is populated with almost 1.7 mln people. The largest cities include Bielsko Biała (170 000 inhabitants), Rybnik (140 000 inhabitants) and Jastrzębie-Zdrój (91 000 inhabitants). The northern part of the study area covers the upper Vistula valley with a high urbanization level (Figure 2B). In contrast, the southern part reaches the Carpathian mountains, which limit air mass exchange (Figure 1B). Temperature inversion is frequently observed, especially in foothill valleys in the winter period. This fact limits the boundary layer mixing and contributes to poor air quality [16]. As a consequence, cities within the study area suffer from poor air quality due to high $P{M}_{10}$ concentration[17]

2.2. The GEM-AQ Model

The GEM-AQ is a semi-Lagrangian chemical weather model in which air quality processes (chemistry and aerosols) and tropospheric chemistry are implemented online in the operational weather prediction model, the Global Environmental Multiscale (GEM) [18] model, which was developed at Environment Canada. The gas-phase chemistry mechanism used in the GEM-AQ model is based on a modified version of the Acid Deposition and Oxidants Model (ADOM) [19], where additional reaction in the free troposphere was included [15].

The GEM-AQ model is set up to perform calculations using 28 vertical layers, out of which the lower 21 layers are classified as the troposphere.

Emission data from the polish national emission inventory drive emission sources within the model. These data are based on annual reporting obligations the facilities’ owners fulfil. Annual emissions are transformed into monthly emission rates using the weighting factor from annual emission profiles. Emission profiles are assigned to so-called SNAP categories [20].

2.3. The Gaussian plume model

A Gaussian plume model is a widely used mathematical model for predicting the dispersion of pollutants in the atmosphere. The model assumes that the dispersion of pollutants can be approximated as a two-dimensional Gaussian distribution, which spreads out in a pattern similar to the shape of a bell curve.

The Gaussian plume model is based on the idea that a combination of atmospheric turbulence and wind patterns determines the dispersion of pollutants. The model considers factors such as the source strength and height, wind speed and direction, and atmospheric stability. The topography of the modelled area is not considered.

The Gaussian plume model is used in various applications, including air quality assessments for industrial facilities, roadway emissions, and wildfires. It is often used with other models or measurement techniques to provide a comprehensive picture of air quality in a given area.

The Gaussian plume model was implemented in a parallel Python code. The model is based on the Gaussian plume formula [21], which describes the one-hour average concentration distribution at the surface level:

$$C(x,y)=\frac{E}{2\pi {\sigma }_y{\sigma }_z\overline{u}}exp\left(-\frac{y^2}{2{\sigma }_y^2}\right)exp\left(-\frac{H^2}{2{\sigma }_z^2}\right)·1000$$

(1)

where $\overline{u}$ is the one-hour average wind velocity (assumed to be uniform over the whole computational domain), ${\sigma }_y\left(x\right)$ and ${\sigma }_z\left(x\right)$ are the standard deviations (horizontal and vertical) of plume concentration spatial distribution, which were estimated using formulas proposed by Briggs [22,23]:

$${\sigma }_y\left(x\right)=k_{11}x{(1+k_{22}x)}^{k_{33}}$$

(2)

$${\sigma }_z\left(x\right)=k_{44}x{(1+k_{55}x)}^{k_{66}}$$

(3)

Coefficients $k_{11},...,k_{66}$ depend on the atmosphere stability class (A-F).

Atmosphere stability was classified based on gradient Richardson number criteria [24]. The gradient Richardson number was estimated based on the meteorological output from the GEM-AQ model. The vertical gradients were calculated between the two lowest layers.

H is the plume rise above the surface, which is a sum of stack height $H_e$ and plume rise $dh$ calculated using the combination of Holland and CONCAWE formulas ($d{h}_H$,$d{h}_C$ respectively)[25], depending on the heat flux:

$$dh=\left\{\begin{array}{cc}d{h}_C\hfill & forQ>24000\\ d{h}_H·(24000-Q)/8000+d{h}_C·(Q-16000)/8000\hfill & forQ\in \langle 16000;24000)\\ d{h}_H\hfill & forQ<16000\end{array}\right.$$

(4)

2.4. Emission data

The Polish national emission inventory fully covers the study area. Three major emission SNAP categories from the inventory were used in gaussian modelling: domestic, industrial and transport emissions[26]. Annual emissions are transformed into monthly emission rates using the weighting factor from annual emission profiles

Traffic emissions were represented as point sources distributed across 30 meters along the road network (Figure 3A). For the uplift formula (eq. 4), we assumed fumes temperature 500K and velocity $1m/s$. Domestic emissions are based on the National Database of Topographic Objects (BDOT), a nationwide system of collecting and sharing topographic data, including vector data describing each building as the basis for national emission inventory[26]. For the gaussian model, we assumed stack height as 3.5 times the number of floors + 0.5 meters. Fumes temperature was assumed to be 400K and velocity $0.5m/s$.

Industrial emissions are based on annual reporting obligations, which the facilities’ owners fulfil to the unit responsible for managing the national emission inventory.

2.5. Surface Observations

There are nine air quality stations within the study area. Each of them measures $P{M}_{10}$ concentration with an hourly time step. One meteorological station is located in the centre of the area (Figure 1). Table 1 summarizes the annual observed air quality series. Despite the annual mean at a moderate level (30-40$ug/m^3$), the number of days with the legal threshold (50$ug/m^3$) exceeded is quite significant and covers the major part of the winter season.

2.6. Random forest

As a data fusion algorithm, Random Forest was used. Random Forest is a robust machine learning algorithm used for classification [27], and regression tasks [28]. It is a type of ensemble learning algorithm that combines multiple decision trees to improve the predictive performance of the model[29]. Each decision tree in a random forest is constructed using a different subset of the training data and a random subset of the input features. This is done to introduce diversity and reduce overfitting, as each tree is trained on a different subset of the data and features.

In this work, we trained a random forest algorithm to predict observed $P{M}_{10}$ concentrations at the observation station locations. Input features included concentrations from GEM-AQ and the Gaussian plume model. The models were trained using a 5-fold cross-validation process. We also attempted to use calendar-related variables such as the day of the week number and the month and observed meteorological parameters as additional features. The training dataset was based on time series observations from all nine observation stations.

As a second trial, we tried using datasets based on observations from a single observation station’s time series. This is because air quality observations are primarily influenced by each station’s location, and merging time series from multiple locations may not always be the best approach [30].

Finally, we attempted station vs station cross-validation, as we anticipate the presence of clusters of similar stations in terms of air quality dynamics within the analyzed area. Additionally, information about outliers (i.e., stations that are different from the others) would be helpful in excluding them from the training dataset.

3. Results

3.1. Overall Performance

In order to assess the reproducibility of $P{M}_{10}$ concentration dynamics, we examined the coefficient of determination $R^2$ of the GEM-AQ and Gaussian plume models within a monthly time window (Figure 4 and Figure 5). Both models appear to perform better in winter months (October-March) than in summer (April-September). This pattern is observed at all air quality stations and can be explained by the meteorological factors driving $P{M}_{10}$ concentration, which are better reproduced in winter than in summer. Another reason might be the temporal emission profile, which is more uniform in winter (due to constant demand for heating) than in summer when the daily profile is not uniform.

Both GEM-AQ and Gaussian plume models perform better regarding daily averaged concentration. This fact is due to rapid changes in observed concentrations, which cannot be simulated by any of these models [31]. Some authors [32] explain this concentration variation by emission, which is driven by the air temperature.

Random forest performance was evaluated based on $R^2$, suitable for assessing the dynamics reproducibility. For assessing accuracy, the accuracy coefficient ($ACC$), which describes how accurate the results predicted by the random forest ($\widehat{y}$) are in comparison to observations (y):

$$ACC=\left(1-\sum \frac{|\widehat{y}-y|}{y}\right)·100\%$$

(5)

As the Table 2 reveal - reproducing hourly dynamics of $P{M}_{10}$ concentration was challenging, regardless of the extra features. Using daily averaged concentrations instead of hourly concentration has increased the $R^2$ from around 0.2 to 0.4 and accuracy from 48% to 60%. Using additional features improves model dynamics in all cases while the accuracy remains almost the same.

3.2. Temporal Comparison

We analyzed the performance of random forest models trained on data from one month only. We used 5-fold cross-validation and data from all the observation stations. The process was repeated for hourly data and daily averages. As the data from Table 3 reveal - the best performance in terms of dynamics reproduction ($R^2$) was obtained for winter months. At the same time, the accuracy (thus reproduction of the magnitude of observed concentrations) was better in the summer months. This fact can be explained by a general tendency of ensemble methods which are not very good at reproducing peak values. Also, some authors claim that a significant amount of emissions is not included in the national emission inventory [33,34].

3.3. Spatial Comparison

Finally, we analyzed if the choice of observation station location influenced the performance of the random forest. We used a 5-fold cross-validation process and data from one station at a time. The results make it possible to distinguish stations with sufficiently better performance (SlBielKossak, SlWodzGalczy) and stations with significantly worse one (MpSuchaNiesz, MpOswiecBema) - Table 4. This difference can be explained when we look at the station location. The former is located in dense urban areas with local district heating systems. In contrast, the latter is located in a single-family housing area with typical low-stack residential heating emissions.

3.4. Annual statistics

The annual statistics of the GEM-AQ model (Figure 6A, Figure 9A and Figure 8A) resemble the emission pattern (Figure 3 and Figure 2B). High concentrations are observed in the western part of the study area (Rybnik and Wodzisław cities) and the centre (Bielsko-Biała and Żywiec). Results from the Gaussian plume model highly underestimate the average $P{M}_{10}$ concentration (Figure 6B). Both the Gaussian plume model and random forest reveal a complex concentration pattern in the southern part of the study area, resulting from complex topography. The order of magnitude of the random forest results is similar to that of the GEM-AQ model. The random forest is not good at reproducing peak concentrations; thus, the percentile 90.2% for the random forest is generally lower than for GEM-AQ (Figure 8). Also, the number of days when the legal threshold of $50\mu g/m^3$ is exceeded is lower for a random forest than for GEM-AQ.

4. Conclusion

Random forest regression is a powerful and robust technique for developing non-linear regression models. As we have shown, it can be applied to obtaining high-resolution concentration maps based on regional model results. As a random forest cannot extrapolate data, its results are slightly underestimated.

The accuracy of random forest improves when applied to daily averaged values. This is likely due to the smoothing effect of averaging, thanks to which no sharp gradients must be simulated. Additional improvements to the random forest regression model can be made by using additional features. Including the day of the week and month improved both accuracy and dynamics of all random forest variants. These features act as a non-explicit temporal profile, which helps to adjust regression to emission temporal changes. Including meteorological observation (temperature and wind) as additional features, is also helpful in improving random forest regression results. However, the improvement is less significant in this case. The effect of meteorological observations is likely a way of fixing the inaccuracy of meteorological results of the GEM-AQ model, which later on affected the air quality results.

The choice of observation stations for random forest training should be made with care. Some stations tend to deliver observations which could be more challenging to replicate. On the other hand, using single observation station time series as a training target over large areas may produce results which lack universality.