# Assessing Spatiotemporal Characteristics of Urban PM2.5 Using Fractal Dimensions and Wavelet Analysis

Authors: Ligang Wang; Haizhong Zhang; Liang Mao; Shan Li; Hao Wu
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

Due to rapid urbanization and industrialization, atmospheric fine particulate matter (PM2.5) has become a primary urban pollutant, seriously affecting air quality and resident health. Existing airborne remote sensing and ground sensor monitoring can efficiently collect PM2.5 data. It is urgent yet challenging to fully use these two monitoring modes to analyze the spatial distribution and dynamic changes of PM2.5. This paper proposes a method to analyze the spatiotemporal characteristics of urban PM2.5 concentration using airborne and ground monitoring data. The method utilizes the boundary dimension and the radius dimension to describe the spatial distribution characteristics of PM2.5 and then adopts wavelet analysis to detect the fluctuation periods of PM2.5 concentration. Case study was performed based on the moderate-resolution imaging spectroradiometer (MODIS) data and the daily ground monitoring concentration of PM2.5 from 2014 to 2017 in Beijing. The decrease results of boundary dimension indicate that the changes in PM2.5 concentration tend to be slower over time, and the fluctuation amplitude gradually becomes smaller. The increase results closer to 2.0 in the radius dimensions suggest that the PM2.5 distribution becomes more uniform, indicating that the pollution control in Beijing had achieved initial success. The frequency results of wavelet analysis reveal that PM2.5 concentration in Beijing is mainly subject to periodic variation of 190 days. These findings can help us to gain deeper insight into the complex spatiotemporal characteristics of urban PM2.5.


---

## Body

## 1. Introduction

With the acceleration of global urbanization and industrialization, urban air pollution problems have become increasingly serious [1, 2]. Atmospheric fine particulate matter (PM2.5), as the primary pollutant in many cities, seriously affects environmental air quality and human health [3, 4]. Rapid and effective monitoring of urban PM2.5 concentration, as well as its spatiotemporal characteristics analysis, are important for preventing urban air pollution and ensuring urban environment safety and resident health [2, 5].At present, the methods for obtaining the characteristics of urban PM2.5 concentration mainly include airborne remote sensing and ground sensor monitoring [6–10]. Airborne platforms have the advantage of a wide monitoring range, fast information acquisition, and low cost, which have been widely used in many countries and regions around the world. Aerosol optical depth (AOD) data retrieved through satellite data inversion are highly correlated with ground PM2.5 concentration so as to achieve continuous spatial monitoring. Thus, it has been extensively applied in estimating the ground PM2.5 concentration and further reveals the spatiotemporal variation of PM2.5 concentration in urban landscape [11–14]. Although airborne remote sensing can provide AOD data with higher spatial resolution and continuous spatial coverage, PM2.5 monitoring at ground platforms can provide results with greater accuracy [15]. Ground platforms refer to monitoring stations established on the ground that collect pollution concentration data from nearby small areas to provide data support for air pollution monitoring, with the advantage of high precision and robustness. Nowadays, China-certified methods for the ground monitoring of PM2.5 consist of the gravimetric, β-ray, and tapered element oscillating microbalance methods. The global Aerosol Robotic Network (AERONET) jointly established by the United States and France is home to more than 500 monitoring stations in North America, South America, Europe, and Asia. The AOD data measured by AERONET ground monitoring stations are usually taken as reference for evaluating the accuracy of the AOD data retrieved through remote sensing images. Since airborne platforms and ground platforms have already been well developed for the acquisition of PM2.5 data, further exploration of the spatiotemporal distribution characteristics of PM2.5 concentration from a large amount of data needs to be urgently addressed.To date, studies on the spatiotemporal characteristics of urban PM2.5 have mainly focused on data collection and processing. For instance, many popular methods such as the deep learning [16], geographical weighted regression [17], machine learning [18], artificial neural network [19], and improved generalised regression neural network (GRNN) model [20] have been utilized for estimating PM2.5 concentration from AOD data. While the exploration of data changes is still limited to traditional statistical methods. In particular, there are few case studies providing in-depth analysis from the temporal and spatial dimensions. Yan et al. [21] made use of spatial autocorrelation analysis to explore the spatial variation of PM2.5 in the Beijing-Tianjin-Hebei region. Chen et al. [22] performed Spearman’s rank analysis and complete ensemble empirical mode decomposition with adaptive noise to study the spatiotemporal distribution characteristics of PM2.5 and the impact of meteorological factors on PM2.5 in Nanjing. Wang et al. [23] employed multitemporal and multispatial scale statistical analysis to reveal the pattern of real-time monitoring PM2.5 in 338 Chinese cities from 2014 to 2017. It is clear that the variations of urban PM2.5 have both temporal and spatial characteristics, which are similar to the spatiotemporal pattern of temperature field, humidity field, thermal field, and heat island effect in geographical studies. Therefore, the fractal dimension method widely applied to explore geographical features offers great potential in PM2.5 studies. Weng [24] studied the spatial variation of surface radiant temperature caused by the thermal behavior of different land-use types and landscape patterns by calculating the boundary dimension of transects in the distribution map of surface radiant temperature. Wu et al. [25] employed the boundary, radius, and information entropy dimensions to study the spatial pattern of urban land use. As for temporal characteristics analysis, the wavelet analysis method can reflect the detailed, local variation characteristics of time series [26] and thus can be used to analyze the time-frequency variation of PM2.5. Chakraborty and Okaya [27] noted that wavelet transform is useful to decompose signals and identify frequency features and employed it to exam the characteristics of seismic periods. Chen et al. [28] proposed a time-frequency analysis method using synchrosqueezing wavelet transform, which can accurately identify variable frequency components and further detect high-frequency attenuation anomalies and deep weak signals. Despite these significant contributions to identifying the frequency features by wavelet analysis, no substantial research has yet been attempted to investigate the temporal pattern of PM2.5 concentration. After all, the temporal records of PM2.5 concentration exhibit a high degree of complexity due to the influence of various meteorological factors and human activities.The objective of this paper is to propose a method to analyze the spatiotemporal characteristics of urban PM2.5 concentration using airborne and ground monitoring data. This method utilizes the boundary dimension and the radius dimension to describe the spatial distribution characteristics of PM2.5 and then adopts wavelet analysis to detect the fluctuation periods of PM2.5 concentration. It will be of great significance for preventing urban air pollution and ensuring urban environment safety and resident health. The paper is organized as follows. Section2 provides the detailed descriptions of the AOD retrieval and spatiotemporal analysis methods. Section 3 introduces the study area and data source. Section 4 and 5 report experimental results and analyses for spatial and temporal characteristic of PM2.5, respectively. Finally, some key conclusions drawn from the work as well as the limitations are presented in Section 6.


## 2. Methodology

### 2.1. Remote Sensing Inversion of PM2.5 Based on Airborne Platforms

Aerosol, as one of the main components of PM2.5, has become a research focus for air pollution. AOD reflects the optical properties of aerosol and is a common indicator of aerosol concentration. It is possible to revert satellite remote sensing data to AOD data and then establish a linear relationship between the AOD data and the PM2.5 concentration so as to obtain the distribution of PM2.5 concentration in the region of interest. To date, many efforts across the globe have been directed at retrieving AOD data based on remote sensing images, with great progress made in the development of aerosol inversion algorithm. This study employed the V5.2 algorithm to perform AOD inversion as follows.


#### 2.1.1. AOD Retrieval Method

Levy et al. [29] proposed V5.2 algorithm extending the dark pixel algorithm [30] by using the scattering angle θ and the vegetation index NDVISWIR as auxiliary data to improve the relationship between the apparent reflectance at 2.1 μm and the surface reflectance in the visible bands (blue and red). The V5.2 algorithm was employed for generating MODIS AOD products [31], which could accurately calculate the surface reflectance in the red and blue bands.


#### 2.1.2. Correlation Correction between AOD and PM2.5 Concentration

PM2.5 concentration values can be retrieved based on an inversion of the relationship between the AOD data and measured PM2.5 concentration values. It is noticed that the direct correlation between AOD and PM2.5 concentration is subject to the vertical distribution of AOD, humidity, chemical composition, and other factors. Therefore, it is necessary to perform vertical correction on satellite-retrieved AOD and humidity correction on ground-measured PM2.5 concentration to increase the inversion accuracy of PM2.5 concentration.(1) Vertical Correction. AOD is the integral of the extinction coefficient from the ground to the top of the atmosphere in the vertical direction, whereas the reported PM2.5 concentration only represents the concentration at the ground monitoring stations. Thus, it is necessary to introduce the concept of aerosol scale height for vertical correction on AOD to obtain near-surface extinction coefficients so that the correlation between AOD and PM2.5 concentration can be improved. The correction expression proposed by Li et al. [32] is(1)MAOD=τaHi,where MAOD is an AOD value after the vertical correction, τa is an original AOD value retrieved through remote sensing inversion, and Hi is the aerosol scale height. The Peterson model [33] quantifies the relationship between AOD, visibility, and aerosol scale height as(2)τa=Hi×3.0V−0.0146,where τa refers to AOD and V is the ground visibility in kilometers. Given that AOD is significantly different in different months, it was corrected separately for four seasons. According to equation (2), the aerosol scale height of the four seasons was calculated to be 2.094 km, 2.306 km, 1.413 km, and 2.107 km, respectively, and it was subsequently substituted in equation (1) to obtain the vertically corrected AOD.(2) Humidity Correction. Air samples should be dried before PM2.5 concentration measurement, whereas AOD inversion is performed under natural conditions. The presence of a large amount of water vapor in the air will affect inversion accuracy, as aerosol particles will become larger in volume owing to the absorption of water vapour. Thus, it is mandatory to perform humidity correction on PM2.5 mass concentration [34].(3)fRH=11.0−RH/100,where fRH is a humidity correction factor and RH refers to the relative humidity. The multiplication of PM2.5 concentration by fRH yields the corrected PM2.5 concentration.


### 2.2. Spatial Characteristics Analysis of PM2.5 Based on Fractal Dimensions

Fractal dimensions have been extensively used in analyzing the classification objects with similarity or statistical self-similarity, reflecting the effectiveness of complex shapes in occupying the space and serving as a metric of the irregularity of complex shapes. Given that the changes in urban PM2.5 mass concentration have a certain degree of self-similarity, fractal dimensions were employed to describe the complexity of PM2.5 concentration changes after the spatial distribution of PM2.5 was revealed by the inversion of airborne remote sensing images.


#### 2.2.1. Boundary Dimension

The boundary dimension was introduced to exam the city boundary complexity by Batty and Longley [35]. Letting L be a fractal curve, one obtains(4)L1=L2δ1−Dwhere D is the boundary dimension of L, δ is the scale of L, L1 is the Euclidean length of L, and L2 is the Hausdorff length of L. If δ is the step size of the divider gauge, L1 represents the length of L measured by the gauge with a spacing of δ, and L1/δ represents the number of steps in measuring L based on the step size of δ, as expressed in the following equation:(5)Nδ=L1δ.Combining equations (4) and (5), taking logarithm of both sides and taking limit, one obtains the following equation:(6)D=limδ⟶0lnNδ−lnδ,where D is the boundary dimension to be determined. By changing the step size δi, a series of lengths Niδi can be obtained. By taking the logarithm of δi and Niδi, fitting the data points lnδi,lnNiδi to a straight line, and obtaining its slope, D is thus determined.


#### 2.2.2. Radius Dimension

The radius dimension was first proposed by Frankhauser [36] and is usually used in the spatial structure analysis of urban land use. For instance, Wu et al. [37] proposed the radius fractal dimension to quantify the spatial variation of different land-use types around the heat center. If a circle of radius r is made at the center of the city, and the land-use area in the circle is Sr, and then Sr and r satisfy the following relationship:(7)lnSr=lnη+Dlnr.In equation (7), η is a constant and D is the radius dimension, which reflects the density decay of the land-use type from the city center to the surroundings. When D<2, it indicates that the spatial distribution of a given land-use type gradually decreases from the center to the surroundings, and its changes are nonlinear and nonuniform, with the decay rate becoming increasingly high; when D=2, it indicates that the spatial distribution changes uniformly from the center to the surroundings; when D>2, it indicates that the spatial distribution increases gradually from the center to the surroundings. A linear relationship between lnSri and lnri is established, and its slope is regarded as the radius dimension.


### 2.3. Wavelet-Transform-Based Temporal Characteristics Analysis of PM2.5

The PM2.5 concentration obtained by urban ground platforms is subject to a certain degree of temporal periodic variation. However, owing to the large fluctuation of PM2.5 concentration and the interference of noise, the variation trend can only be judged empirically but fail to gain deeper insight into the multiscale characteristics of PM2.5 concentration variation. Wavelet analysis is multiresolution, which can be used to analyze PM2.5 time series quickly in an intuitive manner. Its good time-frequency characteristics can be extended to the two dimensions of time and frequency, allowing the detailed analysis at multiple scales and eliminate the effects of various accidental factors, which makes it possible to have more accurate analysis of the periodicity and regularity of PM2.5 concentration.Wavelet transform is another effective time-frequency analysis method. It is multiresolution, which overcomes the shortcoming of Fourier transform with a fixed-size window. It is widely used in many fields such as signal filtering, image processing, medical diagnosis, and seismic exploration. Wavelet transform [38] refers to the superposition of wavelets at different scales and displacements, which is performed by selecting a certain scale and translation range and then performing an internal integration of the signals at different scales. For the signal ft∈L2R, its continuous wavelet transform is(8)Wfa,b=1a∫−∞+∞ftψ∗t−badt,where Wfa,b is the wavelet transform coefficient, ψ∗t−b/a is the complex conjugate of t−b/a, a is the scale factor, b is the translation factor, and t is the time. At the same temporal scale, if the wavelet coefficient is calculated to be greater than zero, it indicates that the analyzed data sequence has a relatively large value at the point; if it is less than zero, the value at the point is relatively small; if it is equal to zero, it indicates that the point is an abrupt change point, and the larger the value fluctuation scale of the data, the larger the absolute value of wavelet coefficient at the point.All the wavelet coefficients obtained through a wavelet transform of the time series are squared and then integrated with respect to the translation factorb in the corresponding time domain to obtain wavelet variance, which can be plotted against different time scales to obtain a wavelet variance diagram where the change intensity of the analyzed data is disclosed, with the corresponding time scale of maximum curve value being the major change period. Both Baubechies (db) wavelet and Morlet wavelet were used to process the time series of PM2.5 concentration in order to gain a deep insight into its periodic variation characteristics.


## 2.1. Remote Sensing Inversion of PM2.5 Based on Airborne Platforms

Aerosol, as one of the main components of PM2.5, has become a research focus for air pollution. AOD reflects the optical properties of aerosol and is a common indicator of aerosol concentration. It is possible to revert satellite remote sensing data to AOD data and then establish a linear relationship between the AOD data and the PM2.5 concentration so as to obtain the distribution of PM2.5 concentration in the region of interest. To date, many efforts across the globe have been directed at retrieving AOD data based on remote sensing images, with great progress made in the development of aerosol inversion algorithm. This study employed the V5.2 algorithm to perform AOD inversion as follows.


### 2.1.1. AOD Retrieval Method

Levy et al. [29] proposed V5.2 algorithm extending the dark pixel algorithm [30] by using the scattering angle θ and the vegetation index NDVISWIR as auxiliary data to improve the relationship between the apparent reflectance at 2.1 μm and the surface reflectance in the visible bands (blue and red). The V5.2 algorithm was employed for generating MODIS AOD products [31], which could accurately calculate the surface reflectance in the red and blue bands.


### 2.1.2. Correlation Correction between AOD and PM2.5 Concentration

PM2.5 concentration values can be retrieved based on an inversion of the relationship between the AOD data and measured PM2.5 concentration values. It is noticed that the direct correlation between AOD and PM2.5 concentration is subject to the vertical distribution of AOD, humidity, chemical composition, and other factors. Therefore, it is necessary to perform vertical correction on satellite-retrieved AOD and humidity correction on ground-measured PM2.5 concentration to increase the inversion accuracy of PM2.5 concentration.(1) Vertical Correction. AOD is the integral of the extinction coefficient from the ground to the top of the atmosphere in the vertical direction, whereas the reported PM2.5 concentration only represents the concentration at the ground monitoring stations. Thus, it is necessary to introduce the concept of aerosol scale height for vertical correction on AOD to obtain near-surface extinction coefficients so that the correlation between AOD and PM2.5 concentration can be improved. The correction expression proposed by Li et al. [32] is(1)MAOD=τaHi,where MAOD is an AOD value after the vertical correction, τa is an original AOD value retrieved through remote sensing inversion, and Hi is the aerosol scale height. The Peterson model [33] quantifies the relationship between AOD, visibility, and aerosol scale height as(2)τa=Hi×3.0V−0.0146,where τa refers to AOD and V is the ground visibility in kilometers. Given that AOD is significantly different in different months, it was corrected separately for four seasons. According to equation (2), the aerosol scale height of the four seasons was calculated to be 2.094 km, 2.306 km, 1.413 km, and 2.107 km, respectively, and it was subsequently substituted in equation (1) to obtain the vertically corrected AOD.(2) Humidity Correction. Air samples should be dried before PM2.5 concentration measurement, whereas AOD inversion is performed under natural conditions. The presence of a large amount of water vapor in the air will affect inversion accuracy, as aerosol particles will become larger in volume owing to the absorption of water vapour. Thus, it is mandatory to perform humidity correction on PM2.5 mass concentration [34].(3)fRH=11.0−RH/100,where fRH is a humidity correction factor and RH refers to the relative humidity. The multiplication of PM2.5 concentration by fRH yields the corrected PM2.5 concentration.


## 2.1.1. AOD Retrieval Method

Levy et al. [29] proposed V5.2 algorithm extending the dark pixel algorithm [30] by using the scattering angle θ and the vegetation index NDVISWIR as auxiliary data to improve the relationship between the apparent reflectance at 2.1 μm and the surface reflectance in the visible bands (blue and red). The V5.2 algorithm was employed for generating MODIS AOD products [31], which could accurately calculate the surface reflectance in the red and blue bands.


## 2.1.2. Correlation Correction between AOD and PM2.5 Concentration

PM2.5 concentration values can be retrieved based on an inversion of the relationship between the AOD data and measured PM2.5 concentration values. It is noticed that the direct correlation between AOD and PM2.5 concentration is subject to the vertical distribution of AOD, humidity, chemical composition, and other factors. Therefore, it is necessary to perform vertical correction on satellite-retrieved AOD and humidity correction on ground-measured PM2.5 concentration to increase the inversion accuracy of PM2.5 concentration.(1) Vertical Correction. AOD is the integral of the extinction coefficient from the ground to the top of the atmosphere in the vertical direction, whereas the reported PM2.5 concentration only represents the concentration at the ground monitoring stations. Thus, it is necessary to introduce the concept of aerosol scale height for vertical correction on AOD to obtain near-surface extinction coefficients so that the correlation between AOD and PM2.5 concentration can be improved. The correction expression proposed by Li et al. [32] is(1)MAOD=τaHi,where MAOD is an AOD value after the vertical correction, τa is an original AOD value retrieved through remote sensing inversion, and Hi is the aerosol scale height. The Peterson model [33] quantifies the relationship between AOD, visibility, and aerosol scale height as(2)τa=Hi×3.0V−0.0146,where τa refers to AOD and V is the ground visibility in kilometers. Given that AOD is significantly different in different months, it was corrected separately for four seasons. According to equation (2), the aerosol scale height of the four seasons was calculated to be 2.094 km, 2.306 km, 1.413 km, and 2.107 km, respectively, and it was subsequently substituted in equation (1) to obtain the vertically corrected AOD.(2) Humidity Correction. Air samples should be dried before PM2.5 concentration measurement, whereas AOD inversion is performed under natural conditions. The presence of a large amount of water vapor in the air will affect inversion accuracy, as aerosol particles will become larger in volume owing to the absorption of water vapour. Thus, it is mandatory to perform humidity correction on PM2.5 mass concentration [34].(3)fRH=11.0−RH/100,where fRH is a humidity correction factor and RH refers to the relative humidity. The multiplication of PM2.5 concentration by fRH yields the corrected PM2.5 concentration.


## 2.2. Spatial Characteristics Analysis of PM2.5 Based on Fractal Dimensions

Fractal dimensions have been extensively used in analyzing the classification objects with similarity or statistical self-similarity, reflecting the effectiveness of complex shapes in occupying the space and serving as a metric of the irregularity of complex shapes. Given that the changes in urban PM2.5 mass concentration have a certain degree of self-similarity, fractal dimensions were employed to describe the complexity of PM2.5 concentration changes after the spatial distribution of PM2.5 was revealed by the inversion of airborne remote sensing images.


### 2.2.1. Boundary Dimension

The boundary dimension was introduced to exam the city boundary complexity by Batty and Longley [35]. Letting L be a fractal curve, one obtains(4)L1=L2δ1−Dwhere D is the boundary dimension of L, δ is the scale of L, L1 is the Euclidean length of L, and L2 is the Hausdorff length of L. If δ is the step size of the divider gauge, L1 represents the length of L measured by the gauge with a spacing of δ, and L1/δ represents the number of steps in measuring L based on the step size of δ, as expressed in the following equation:(5)Nδ=L1δ.Combining equations (4) and (5), taking logarithm of both sides and taking limit, one obtains the following equation:(6)D=limδ⟶0lnNδ−lnδ,where D is the boundary dimension to be determined. By changing the step size δi, a series of lengths Niδi can be obtained. By taking the logarithm of δi and Niδi, fitting the data points lnδi,lnNiδi to a straight line, and obtaining its slope, D is thus determined.


### 2.2.2. Radius Dimension

The radius dimension was first proposed by Frankhauser [36] and is usually used in the spatial structure analysis of urban land use. For instance, Wu et al. [37] proposed the radius fractal dimension to quantify the spatial variation of different land-use types around the heat center. If a circle of radius r is made at the center of the city, and the land-use area in the circle is Sr, and then Sr and r satisfy the following relationship:(7)lnSr=lnη+Dlnr.In equation (7), η is a constant and D is the radius dimension, which reflects the density decay of the land-use type from the city center to the surroundings. When D<2, it indicates that the spatial distribution of a given land-use type gradually decreases from the center to the surroundings, and its changes are nonlinear and nonuniform, with the decay rate becoming increasingly high; when D=2, it indicates that the spatial distribution changes uniformly from the center to the surroundings; when D>2, it indicates that the spatial distribution increases gradually from the center to the surroundings. A linear relationship between lnSri and lnri is established, and its slope is regarded as the radius dimension.


## 2.2.1. Boundary Dimension

The boundary dimension was introduced to exam the city boundary complexity by Batty and Longley [35]. Letting L be a fractal curve, one obtains(4)L1=L2δ1−Dwhere D is the boundary dimension of L, δ is the scale of L, L1 is the Euclidean length of L, and L2 is the Hausdorff length of L. If δ is the step size of the divider gauge, L1 represents the length of L measured by the gauge with a spacing of δ, and L1/δ represents the number of steps in measuring L based on the step size of δ, as expressed in the following equation:(5)Nδ=L1δ.Combining equations (4) and (5), taking logarithm of both sides and taking limit, one obtains the following equation:(6)D=limδ⟶0lnNδ−lnδ,where D is the boundary dimension to be determined. By changing the step size δi, a series of lengths Niδi can be obtained. By taking the logarithm of δi and Niδi, fitting the data points lnδi,lnNiδi to a straight line, and obtaining its slope, D is thus determined.


## 2.2.2. Radius Dimension

The radius dimension was first proposed by Frankhauser [36] and is usually used in the spatial structure analysis of urban land use. For instance, Wu et al. [37] proposed the radius fractal dimension to quantify the spatial variation of different land-use types around the heat center. If a circle of radius r is made at the center of the city, and the land-use area in the circle is Sr, and then Sr and r satisfy the following relationship:(7)lnSr=lnη+Dlnr.In equation (7), η is a constant and D is the radius dimension, which reflects the density decay of the land-use type from the city center to the surroundings. When D<2, it indicates that the spatial distribution of a given land-use type gradually decreases from the center to the surroundings, and its changes are nonlinear and nonuniform, with the decay rate becoming increasingly high; when D=2, it indicates that the spatial distribution changes uniformly from the center to the surroundings; when D>2, it indicates that the spatial distribution increases gradually from the center to the surroundings. A linear relationship between lnSri and lnri is established, and its slope is regarded as the radius dimension.


## 2.3. Wavelet-Transform-Based Temporal Characteristics Analysis of PM2.5

The PM2.5 concentration obtained by urban ground platforms is subject to a certain degree of temporal periodic variation. However, owing to the large fluctuation of PM2.5 concentration and the interference of noise, the variation trend can only be judged empirically but fail to gain deeper insight into the multiscale characteristics of PM2.5 concentration variation. Wavelet analysis is multiresolution, which can be used to analyze PM2.5 time series quickly in an intuitive manner. Its good time-frequency characteristics can be extended to the two dimensions of time and frequency, allowing the detailed analysis at multiple scales and eliminate the effects of various accidental factors, which makes it possible to have more accurate analysis of the periodicity and regularity of PM2.5 concentration.Wavelet transform is another effective time-frequency analysis method. It is multiresolution, which overcomes the shortcoming of Fourier transform with a fixed-size window. It is widely used in many fields such as signal filtering, image processing, medical diagnosis, and seismic exploration. Wavelet transform [38] refers to the superposition of wavelets at different scales and displacements, which is performed by selecting a certain scale and translation range and then performing an internal integration of the signals at different scales. For the signal ft∈L2R, its continuous wavelet transform is(8)Wfa,b=1a∫−∞+∞ftψ∗t−badt,where Wfa,b is the wavelet transform coefficient, ψ∗t−b/a is the complex conjugate of t−b/a, a is the scale factor, b is the translation factor, and t is the time. At the same temporal scale, if the wavelet coefficient is calculated to be greater than zero, it indicates that the analyzed data sequence has a relatively large value at the point; if it is less than zero, the value at the point is relatively small; if it is equal to zero, it indicates that the point is an abrupt change point, and the larger the value fluctuation scale of the data, the larger the absolute value of wavelet coefficient at the point.All the wavelet coefficients obtained through a wavelet transform of the time series are squared and then integrated with respect to the translation factorb in the corresponding time domain to obtain wavelet variance, which can be plotted against different time scales to obtain a wavelet variance diagram where the change intensity of the analyzed data is disclosed, with the corresponding time scale of maximum curve value being the major change period. Both Baubechies (db) wavelet and Morlet wavelet were used to process the time series of PM2.5 concentration in order to gain a deep insight into its periodic variation characteristics.


## 3. Study Area and Data

### 3.1. Study Area

Beijing is located between 115°25′ and 117°30′ E and 39°28′–41°05′ N, with a total area of approximately 1,6807.8 km2, in the northeastern part of the North China Plain and borders the Bohai Bay. Surrounded by mountains on three sides, Beijing terrain on the whole is characterized by high elevation in the northwest and low elevation in the southeast. Beijing is located in the northern temperate zone, under a warm, semihumid, and semiarid monsoon continental climate with distinct seasons.


### 3.2. Data and Sources

This study involves three data sources, namely, MODIS remote sensing image data, ground-monitored PM2.5 mass concentration data, and other ancillary data. (1) MODIS L2 aerosol product at 3 km resolution were selected in Beijing from 2014 to 2017 (https://modis.gsfc.nasa.gov/). (2) The real-time, hourly ground-monitored PM2.5 concentration at 35 ground monitoring stations in Beijing from 2014 to 2017 were collected from the Beijing Municipal Environmental Protection Bureau (http://zx.bjmemc.com.cn), with the locations of the 35 stations shown in Figure 1(c). For each station, PM2.5 concentration was averaged every 24 hours to obtain the daily mean concentration of the station, and each station’s daily mean concentration was summed and divided by 35 to obtain the daily mean concentration of Beijing. The Beijing daily mean concentration of PM2.5 was plotted against the date over the period from 2014 to 2017 (Figure 2). (3) Ancillary data consist of AERONET ground-station aerosol data and relative humidity and ground visibility data in Beijing were collected from the AERONET and the US NOAA (https://gis.ncdc.noaa.gov).Figure 1
Schematic diagram of the study region: (a, b) Beijing location; (c) Beijing administrative division and 35 ground monitoring stations; (d) major sources of PM2.5 such as mobile sources, dust sources, coal combustion, and industrial sources in Beijing, as well as the locations of the three section lines of analysis.
(a)
(b)
(c)
(d)Figure 2
Daily mean PM2.5 concentration at ground monitoring stations in Beijing from 2014 to 2017.Existing studies show that urban PM2.5 is mainly derived from the influence of mobile sources, dust sources, coal combustion, and industrial sources. Combining both the airborne and the ground monitoring data, this study analyzed the effects of the above four types of sources on the Beijing PM2.5 concentration. Among them, mobile sources were taken into account from the four aspects, including airports, railway stations, famous scenic spots, and major highways. Dust source analysis was focused on the effect of the dust from construction processes. Coal combustion analysis was focused on four major coal-fired power plants, namely, Huaneng Beijing Thermal Power Plant, Guohua Beijing Thermal Power Plant, Gaojing Thermal Power Plant, and Shijingshan Thermal Power Plant. Industrial source analysis was performed from the aspects of petrochemical companies, automobile companies, and printing companies. Centered in Tian’anmen, three section lines were selected to analyze the spatiotemporal distribution characteristics of PM2.5 as shown in Figure1(d). Particularly, the first section line I passed sequentially through Beijing Dongfang Petrochemical Co., Ltd., Beijing Hualian Printing Co., Ltd., Jiugong Town Green-Belt Industry Development Project of Beijing Hongtu Jianye Investment Management Center, Heaven Temple, Tian’anmen, Prince Kung’s Mansion, Beijing North Railway Station, Summer Palace, Beiqi Foton Motor Co., Ltd., China Mobile Communication Group International Information Port Project, Ming Tombs, and Badaling Section of the Great Wall. The second section line II passed sequentially through Beijing Wenlinyuan Printing Co., Ltd., Beijing Lianxing Shengye Printing Co., Ltd., Nanyuan Airport, Beijing South Railway Station, Tian’anmen, Capital International Airport, Beijing Automobile Group Co., Ltd., Beijing Hyundai Motor Co., Ltd., Beijing Tairui Printing Co., Ltd., Longshan International Hotel and Taiwan Building Project of Beijing Huairou District, Beijing Bosa Auto Parts Co., Ltd., and Beijing Zhongxing Printing Co., Ltd. As for the third section line III, it passed sequentially through Beijing Tianlun Auto Parts Co., Ltd., Beijing Xishan Yijing Jiayuan Project, Gaojing Thermal Power Plant, Shijingshan Thermal Power Plant, Beijing West Railway Station, Tian’anmen, Beijing Railway Station, Guohua Beijing Thermal Power Plant, Beijing East Railway Station, Huaneng Beijing Thermal Power Plant, and Land Block Comprehensive Development Project of Tongzhou District Canal Core Area.


## 3.1. Study Area

Beijing is located between 115°25′ and 117°30′ E and 39°28′–41°05′ N, with a total area of approximately 1,6807.8 km2, in the northeastern part of the North China Plain and borders the Bohai Bay. Surrounded by mountains on three sides, Beijing terrain on the whole is characterized by high elevation in the northwest and low elevation in the southeast. Beijing is located in the northern temperate zone, under a warm, semihumid, and semiarid monsoon continental climate with distinct seasons.


## 3.2. Data and Sources

This study involves three data sources, namely, MODIS remote sensing image data, ground-monitored PM2.5 mass concentration data, and other ancillary data. (1) MODIS L2 aerosol product at 3 km resolution were selected in Beijing from 2014 to 2017 (https://modis.gsfc.nasa.gov/). (2) The real-time, hourly ground-monitored PM2.5 concentration at 35 ground monitoring stations in Beijing from 2014 to 2017 were collected from the Beijing Municipal Environmental Protection Bureau (http://zx.bjmemc.com.cn), with the locations of the 35 stations shown in Figure 1(c). For each station, PM2.5 concentration was averaged every 24 hours to obtain the daily mean concentration of the station, and each station’s daily mean concentration was summed and divided by 35 to obtain the daily mean concentration of Beijing. The Beijing daily mean concentration of PM2.5 was plotted against the date over the period from 2014 to 2017 (Figure 2). (3) Ancillary data consist of AERONET ground-station aerosol data and relative humidity and ground visibility data in Beijing were collected from the AERONET and the US NOAA (https://gis.ncdc.noaa.gov).Figure 1
Schematic diagram of the study region: (a, b) Beijing location; (c) Beijing administrative division and 35 ground monitoring stations; (d) major sources of PM2.5 such as mobile sources, dust sources, coal combustion, and industrial sources in Beijing, as well as the locations of the three section lines of analysis.
(a)
(b)
(c)
(d)Figure 2
Daily mean PM2.5 concentration at ground monitoring stations in Beijing from 2014 to 2017.Existing studies show that urban PM2.5 is mainly derived from the influence of mobile sources, dust sources, coal combustion, and industrial sources. Combining both the airborne and the ground monitoring data, this study analyzed the effects of the above four types of sources on the Beijing PM2.5 concentration. Among them, mobile sources were taken into account from the four aspects, including airports, railway stations, famous scenic spots, and major highways. Dust source analysis was focused on the effect of the dust from construction processes. Coal combustion analysis was focused on four major coal-fired power plants, namely, Huaneng Beijing Thermal Power Plant, Guohua Beijing Thermal Power Plant, Gaojing Thermal Power Plant, and Shijingshan Thermal Power Plant. Industrial source analysis was performed from the aspects of petrochemical companies, automobile companies, and printing companies. Centered in Tian’anmen, three section lines were selected to analyze the spatiotemporal distribution characteristics of PM2.5 as shown in Figure1(d). Particularly, the first section line I passed sequentially through Beijing Dongfang Petrochemical Co., Ltd., Beijing Hualian Printing Co., Ltd., Jiugong Town Green-Belt Industry Development Project of Beijing Hongtu Jianye Investment Management Center, Heaven Temple, Tian’anmen, Prince Kung’s Mansion, Beijing North Railway Station, Summer Palace, Beiqi Foton Motor Co., Ltd., China Mobile Communication Group International Information Port Project, Ming Tombs, and Badaling Section of the Great Wall. The second section line II passed sequentially through Beijing Wenlinyuan Printing Co., Ltd., Beijing Lianxing Shengye Printing Co., Ltd., Nanyuan Airport, Beijing South Railway Station, Tian’anmen, Capital International Airport, Beijing Automobile Group Co., Ltd., Beijing Hyundai Motor Co., Ltd., Beijing Tairui Printing Co., Ltd., Longshan International Hotel and Taiwan Building Project of Beijing Huairou District, Beijing Bosa Auto Parts Co., Ltd., and Beijing Zhongxing Printing Co., Ltd. As for the third section line III, it passed sequentially through Beijing Tianlun Auto Parts Co., Ltd., Beijing Xishan Yijing Jiayuan Project, Gaojing Thermal Power Plant, Shijingshan Thermal Power Plant, Beijing West Railway Station, Tian’anmen, Beijing Railway Station, Guohua Beijing Thermal Power Plant, Beijing East Railway Station, Huaneng Beijing Thermal Power Plant, and Land Block Comprehensive Development Project of Tongzhou District Canal Core Area.


## 4. Results and Discussion of the Spatial Distribution Characteristics of PM2.5

### 4.1. Results and Validation of Remote Sensing Inversion of PM2.5

AOD distribution maps were retrieved from the MODIS aerosol product using the V5.2 algorithm in Beijing from 2014 to 2017 shown in Figure3. The monitoring values of AOD and PM2.5 at Beijing’s only NOAA meteorological station (40.08° N, 116.585° E) were selected for the training and validation of the inversion model. Given that the NOAA station did not have PM2.5 concentration data, seven PM2.5 monitoring stations close to the NOAA station were selected and their PM2.5 concentration data at satellite transit times were used for kriging interpolation using ArcGIS software, with the interpolated values considered as the PM2.5 concentration data of the NOAA station. In short, a total of 419 effective monitoring data of PM2.5 were collected, of which 312 were used for training the model and 107 for validation. After the vertical and humidity correction, the coefficient of determination (R2) between AOD and PM2.5 concentration (0.8019) shown in Figure 4 indicates a strong linear relationship so that PM2.5 concentration can be accurately estimated.Figure 3
Spatial distribution of annual mean AOD in Beijing from 2014 to 2017: (a) 2014; (b) 2015; (c) 2016; (d) 2017.
(a)
(b)
(c)
(d)Figure 4
Correlation between AOD and PM2.5 concentration after vertical and humidity correction.The spatial distribution of the annual mean of PM2.5 mass concentration from 2014 to 2017 is shown in Figure5. From the temporal aspect, the spatial range of high PM2.5 concentration (PM2.5 > 200) showed a significant decreasing trend over time during the period from 2014 to 2017, and the spatial ranges with the PM2.5 concentration of 100–150 μg·m−3 and 150–200 μg·m−3 also became smaller over time to some extent. From the spatial aspect, the PM2.5 concentration in the central, southern, and eastern parts of Beijing was relatively high but low in the north and west, and the boundary between the high-concentration and low-concentration areas was generally the boundary between the mountains and the plains. Most of the central, southern, and eastern parts are plains with a high population density and many buildings in the presence of intensive industrial activities, where there is serious air pollution from automobile exhaust emissions and industrial emissions, whereas the vegetation coverage is small and the air mobility is poor owing to the blocking by the northern Yanshan Mountains and buildings. All aforementioned factors jointly contribute to the accumulation of pollutants and consequently high PM2.5 concentration. In contrast, the north and west are predominantly mountainous areas with high vegetation coverage, where the population density is low with few factories, so that there are few pollution sources and good air mobility, allowing PM2.5 concentration to be maintained in a low range.Figure 5
Spatial distribution of the annual mean PM2.5 concentration from 2014 to 2017: (a) 2014; (b) 2015; (c) 2016; (d) 2017.
(a)
(b)
(c)
(d)To validate the inversion accuracy of PM2.5 mass concentration, the PM2.5 concentration derived from satellite-based AOD with vertical and humidity correction was compared with the aforementioned 107 measured PM2.5 data. Owing to space limitations, only valid data of the first 20 days were listed in this study, as shown in Table1, which illustrates that the calculated PM2.5 concentration was relatively close to the true concentration, with a mean relative error of 20.34% and a maximum relative error of 40.94%. Despite that there are few data points with large errors, the inversion model on the whole was low in relative error and has proven reliable for PM2.5 concentration inversion in the study area.Table 1
Comparison of estimated PM2.5 concentration with measured PM2.5 concentration (unit:μg·m−3).
Date
Measurement
Estimation
Absolute error
Relative error (%)
2014.01.15
34.93
27.13
7.8
22.33
2014.01.26
95.53
73.74
21.79
22.81
2014.02.06
172.68
202.29
29.61
17.15
2014.03.18
144.16
148.45
4.29
2.98
2014.04.07
165.5
183.47
17.97
10.86
2014.04.09
66.06
79.93
13.87
20.0
2014.05.21
15.47
11.65
3.82
24.70
2014.06.04
72.14
79.27
7.13
9.88
2014.06.11
25.8
23.4
2.4
9.30
2014.06.23
15.67
12.89
2.78
17.74
2014.07.02
8.11
11.43
3.32
40.94
2014.07.09
18.44
22.35
3.91
21.20
2014.08.15
13.46
14.37
0.91
6.76
2014.08.26
9.81
5.93
3.88
39.55
2014.08.28
13.61
9.83
3.78
27.77
2014.08.29
18.94
15.95
2.99
15.79
2014.09.12
7.66
10.07
2.41
31.46
2014.09.16
44.98
59.49
14.51
32.26
2014.09.26
38.71
25.82
12.89
33.30
2014.10.08
54.33
46.71
7.62
14.03


### 4.2. Linear Spatial Characteristics of PM2.5 Concentration

The boundary dimension and coefficient of determination (R2) of the PM2.5 concentration distribution curve at the section lines I, II, and III from 2014 to 2017 are shown in Table 2, where R2 is close to 1, indicating that the regression equation had high goodness-of-fit, and the regression results were highly reliable. Due to the reason that section lines I, II, and III are 1-dimensional fractal initiators, we only discussed the linear spatial characteristics of PM2.5 concentration by the boundary dimension rather than form dimension, although there is a function between the boundary dimension and form dimension reported by Chen [39]. The boundary dimension at the section line I from 2014 to 2017 was 1.1625, 1.0724, 1.0892, and 1.0678, respectively, with a significant decrease in boundary dimension from 2014 to 2015. Although the year 2016 witnessed a slight increase in boundary dimension compared with that in 2015, the boundary dimension in 2017 was still lower than that in 2015, indicating that it showed an overall decreasing trend over the four years. The boundary dimension at the section line II during the same period was 1.249, 1.1945, 1.1424, and 1.1343, respectively, showing an annual decreasing trend. The boundary dimension at the section line III from 2014 to 2017 was 1.2487, 1.1552, 1.1976, and 1.1924, respectively, with a dramatic decrease in 2017 compared with that in 2014. The decrease in boundary dimension indicated that the changes in PM2.5 concentration on these three section lines tended to be slower over time, and the fluctuation amplitude gradually became smaller. In addition, all three section lines passed through the central area of Beijing and a large number of rural areas, whereas the PM2.5 concentration of the rural areas was relatively stable. The decrease in boundary dimension was mainly caused by the gradual decrease in PM2.5 concentration in the central area of Beijing, implying that air quality control has achieved initial success in recent years.Table 2
Boundary dimensions of PM2.5 concentration at each section line from 2014 to 2017 in Beijing.
Year
Section line
Regression equation
Boundary dimension
Coefficient of determination
2014
I
lnN = −1.1625lnδ + 5.9243
1.1625
0.99858
II
lnN = −1.249lnδ + 6.2777
1.249
0.99744
III
lnN = −1.2487lnδ + 6.386
1.2487
0.99949
2015
I
lnN = −1.0724lnδ + 5.7465
1.0724
0.99904
II
lnN = −1.1945lnδ + 5.7843
1.1945
0.99597
III
lnN = −1.1552lnδ + 5.8351
1.1552
0.99972
2016
I
lnN = −1.0892lnδ + 5.6857
1.0892
0.99862
II
lnN = −1.1424lnδ + 5.8284
1.1424
0.99958
III
lnN = −1.1976lnδ + 5.9991
1.1976
0.99722
2017
I
lnN = −1.0678lnδ + 5.6544
1.0678
0.99925
II
lnN = −1.1343lnδ + 5.8918
1.1343
0.99492
III
lnN = −1.1924lnδ + 6.1498
1.1924
0.99889The PM2.5 concentration on each section line from 2014 to 2017 was extracted for spatial variation trend analysis, with the results shown in Figure6 and elaborated below. (1) The airports and the railway stations were vehicle gathering areas, and the automobile exhaust had a significant contribution to PM2.5. The scenic spots located in the center of Beijing with a large passenger flow, such as Tian’anmen, Heaven Temple, Prince Kung’s Mansion, and Summer Palace, had high PM2.5 concentration, whereas the Ming Tombs and the Badaling section of the Great Wall are far away from the city center with a mountainous terrain, and hence, their PM2.5 concentration was low. (2) A significant amount of dust was generated during the construction process, aggravating the environmental pollution. For example, the research site of the Jiugong Town Green-Belt Industry Development Project of Beijing Hongtu Jianye Investment Management Center had the PM2.5 concentration of 151 μg·m−3, 162 μg·m−3, 155 μg·m−3, and 153 μg·m−3 during the period from 2014 to 2017, respectively, which was attributed to dust generation during construction. (3) From 2014 to 2017, the Beijing Municipal Government gradually shut down large coal-fired power plants, and consequently the concentration of PM2.5 also dropped significantly. The Gaojing Thermal Power Plant was shut down in July 2014, where the annual mean PM2.5 concentration in successive years during the period from 2014 to 2017 was 60 μg·m−3, 58 μg·m−3, 50 μg·m−3, and 46 μg·m−3, respectively, showing a gradual decreasing trend. The Shijingshan Thermal Power Plant was officially shut down in March 2015, where the annual mean PM2.5 concentration from 2014 to 2017 was 60 μg·m−3, 61 μg·m−3, 52 μg·m−3, and 47 μg·m−3, respectively, with a dramatic decrease observed in 2016. The Guohua Beijing Thermal Power Plant was shut down in March 2015, and its annual mean PM2.5 concentration decreased significantly from 2014 to 2015, as shown in Figure 6. The Huaneng Beijing Thermal Power Plant was shut down in March 2017, and the PM2.5 concentration has not yet been improved significantly. (4) Petrochemical industries, automobile manufacturing, and printing industries also contributed to PM2.5. The PM2.5 concentration at Beijing Hualian Printing Co., Ltd. in 2014, 2016, and 2017 was 147 μg·m−3, 158 μg·m−3, and 160 μg·m−3, respectively, all in the high-concentration range of PM2.5. The PM2.5 concentration at the research sites of Beijing Automobile Group Co., Ltd. and Beijing Hyundai Motor Co., Ltd. was relatively high, exhibiting an annual increasing trend.Figure 6
PM2.5 concentration variation at each section line.


### 4.3. Planar Spatial Characteristic of PM2.5 Concentration

With Tian’anmen as the center of Beijing, a series of circular buffers was generated on the distribution map of annual mean PM2.5 concentration from 2014 to 2017, with the radius increasing from 4 km to 128 km by increments of 4 km, leading to a total of 32 rings shown in Figure7. According to China’s latest ambient air quality standard (GB 3095-2012), the annual mean concentration limit of national first-class air quality in Beijing is 35 μg·m−3. The area with PM2.5 concentration greater than 35 μg·m−3 in each buffer zone was calculated for each year, followed by taking a natural logarithm of the area and the radius and plotting a scatter plot of lnSr against ln r. Linear regression was performed on the data points, and the slope of the linear fitted line was considered as the radius dimension D.Figure 7
Calculation and analysis of the planar spatial characteristics of Beijing PM2.5.Results show that the radius dimensions of Beijing from 2014 to 2017 were 1.7608, 1.8234, 1.8548, and 1.8622, respectively, as shown in Table3. The variation trend can be observed from Figure 8 that, from 2014 to 2017, the PM2.5 mass concentration in Beijing decreased from the center to the surroundings, and the PM2.5 distribution becomes more uniform, indicating that the pollution control in Beijing had achieved initial success. In addition, the radius dimension in each year from 2014 to 2017 was less than 2, indicating that, for a PM2.5 concentration above 35 μg·m−3, the gradually declining spatial distribution from the center to the surroundings was nonlinear and nonuniform.Table 3
Radius dimensions of PM2.5 concentration from 2014 to 2017 in Beijing.
Year
Regression equation
Radius dimension
Coefficient of determination
2014
lnS(r) = 1.7608lnr + 1.4643
1.7608
0.9979
2015
lnS(r) = 1.8234lnr + 1.3892
1.8234
0.9986
2016
lnS(r) = 1.8548lnr + 1.359
1.8548
0.9987
2017
lnS(r) = 1.8522lnr + 1.3739
1.8522
0.9979Figure 8
Variation trend in the radius dimension of Beijing PM2.5 concentration from 2014 to 2017.


## 4.1. Results and Validation of Remote Sensing Inversion of PM2.5

AOD distribution maps were retrieved from the MODIS aerosol product using the V5.2 algorithm in Beijing from 2014 to 2017 shown in Figure3. The monitoring values of AOD and PM2.5 at Beijing’s only NOAA meteorological station (40.08° N, 116.585° E) were selected for the training and validation of the inversion model. Given that the NOAA station did not have PM2.5 concentration data, seven PM2.5 monitoring stations close to the NOAA station were selected and their PM2.5 concentration data at satellite transit times were used for kriging interpolation using ArcGIS software, with the interpolated values considered as the PM2.5 concentration data of the NOAA station. In short, a total of 419 effective monitoring data of PM2.5 were collected, of which 312 were used for training the model and 107 for validation. After the vertical and humidity correction, the coefficient of determination (R2) between AOD and PM2.5 concentration (0.8019) shown in Figure 4 indicates a strong linear relationship so that PM2.5 concentration can be accurately estimated.Figure 3
Spatial distribution of annual mean AOD in Beijing from 2014 to 2017: (a) 2014; (b) 2015; (c) 2016; (d) 2017.
(a)
(b)
(c)
(d)Figure 4
Correlation between AOD and PM2.5 concentration after vertical and humidity correction.The spatial distribution of the annual mean of PM2.5 mass concentration from 2014 to 2017 is shown in Figure5. From the temporal aspect, the spatial range of high PM2.5 concentration (PM2.5 > 200) showed a significant decreasing trend over time during the period from 2014 to 2017, and the spatial ranges with the PM2.5 concentration of 100–150 μg·m−3 and 150–200 μg·m−3 also became smaller over time to some extent. From the spatial aspect, the PM2.5 concentration in the central, southern, and eastern parts of Beijing was relatively high but low in the north and west, and the boundary between the high-concentration and low-concentration areas was generally the boundary between the mountains and the plains. Most of the central, southern, and eastern parts are plains with a high population density and many buildings in the presence of intensive industrial activities, where there is serious air pollution from automobile exhaust emissions and industrial emissions, whereas the vegetation coverage is small and the air mobility is poor owing to the blocking by the northern Yanshan Mountains and buildings. All aforementioned factors jointly contribute to the accumulation of pollutants and consequently high PM2.5 concentration. In contrast, the north and west are predominantly mountainous areas with high vegetation coverage, where the population density is low with few factories, so that there are few pollution sources and good air mobility, allowing PM2.5 concentration to be maintained in a low range.Figure 5
Spatial distribution of the annual mean PM2.5 concentration from 2014 to 2017: (a) 2014; (b) 2015; (c) 2016; (d) 2017.
(a)
(b)
(c)
(d)To validate the inversion accuracy of PM2.5 mass concentration, the PM2.5 concentration derived from satellite-based AOD with vertical and humidity correction was compared with the aforementioned 107 measured PM2.5 data. Owing to space limitations, only valid data of the first 20 days were listed in this study, as shown in Table1, which illustrates that the calculated PM2.5 concentration was relatively close to the true concentration, with a mean relative error of 20.34% and a maximum relative error of 40.94%. Despite that there are few data points with large errors, the inversion model on the whole was low in relative error and has proven reliable for PM2.5 concentration inversion in the study area.Table 1
Comparison of estimated PM2.5 concentration with measured PM2.5 concentration (unit:μg·m−3).
Date
Measurement
Estimation
Absolute error
Relative error (%)
2014.01.15
34.93
27.13
7.8
22.33
2014.01.26
95.53
73.74
21.79
22.81
2014.02.06
172.68
202.29
29.61
17.15
2014.03.18
144.16
148.45
4.29
2.98
2014.04.07
165.5
183.47
17.97
10.86
2014.04.09
66.06
79.93
13.87
20.0
2014.05.21
15.47
11.65
3.82
24.70
2014.06.04
72.14
79.27
7.13
9.88
2014.06.11
25.8
23.4
2.4
9.30
2014.06.23
15.67
12.89
2.78
17.74
2014.07.02
8.11
11.43
3.32
40.94
2014.07.09
18.44
22.35
3.91
21.20
2014.08.15
13.46
14.37
0.91
6.76
2014.08.26
9.81
5.93
3.88
39.55
2014.08.28
13.61
9.83
3.78
27.77
2014.08.29
18.94
15.95
2.99
15.79
2014.09.12
7.66
10.07
2.41
31.46
2014.09.16
44.98
59.49
14.51
32.26
2014.09.26
38.71
25.82
12.89
33.30
2014.10.08
54.33
46.71
7.62
14.03


## 4.2. Linear Spatial Characteristics of PM2.5 Concentration

The boundary dimension and coefficient of determination (R2) of the PM2.5 concentration distribution curve at the section lines I, II, and III from 2014 to 2017 are shown in Table 2, where R2 is close to 1, indicating that the regression equation had high goodness-of-fit, and the regression results were highly reliable. Due to the reason that section lines I, II, and III are 1-dimensional fractal initiators, we only discussed the linear spatial characteristics of PM2.5 concentration by the boundary dimension rather than form dimension, although there is a function between the boundary dimension and form dimension reported by Chen [39]. The boundary dimension at the section line I from 2014 to 2017 was 1.1625, 1.0724, 1.0892, and 1.0678, respectively, with a significant decrease in boundary dimension from 2014 to 2015. Although the year 2016 witnessed a slight increase in boundary dimension compared with that in 2015, the boundary dimension in 2017 was still lower than that in 2015, indicating that it showed an overall decreasing trend over the four years. The boundary dimension at the section line II during the same period was 1.249, 1.1945, 1.1424, and 1.1343, respectively, showing an annual decreasing trend. The boundary dimension at the section line III from 2014 to 2017 was 1.2487, 1.1552, 1.1976, and 1.1924, respectively, with a dramatic decrease in 2017 compared with that in 2014. The decrease in boundary dimension indicated that the changes in PM2.5 concentration on these three section lines tended to be slower over time, and the fluctuation amplitude gradually became smaller. In addition, all three section lines passed through the central area of Beijing and a large number of rural areas, whereas the PM2.5 concentration of the rural areas was relatively stable. The decrease in boundary dimension was mainly caused by the gradual decrease in PM2.5 concentration in the central area of Beijing, implying that air quality control has achieved initial success in recent years.Table 2
Boundary dimensions of PM2.5 concentration at each section line from 2014 to 2017 in Beijing.
Year
Section line
Regression equation
Boundary dimension
Coefficient of determination
2014
I
lnN = −1.1625lnδ + 5.9243
1.1625
0.99858
II
lnN = −1.249lnδ + 6.2777
1.249
0.99744
III
lnN = −1.2487lnδ + 6.386
1.2487
0.99949
2015
I
lnN = −1.0724lnδ + 5.7465
1.0724
0.99904
II
lnN = −1.1945lnδ + 5.7843
1.1945
0.99597
III
lnN = −1.1552lnδ + 5.8351
1.1552
0.99972
2016
I
lnN = −1.0892lnδ + 5.6857
1.0892
0.99862
II
lnN = −1.1424lnδ + 5.8284
1.1424
0.99958
III
lnN = −1.1976lnδ + 5.9991
1.1976
0.99722
2017
I
lnN = −1.0678lnδ + 5.6544
1.0678
0.99925
II
lnN = −1.1343lnδ + 5.8918
1.1343
0.99492
III
lnN = −1.1924lnδ + 6.1498
1.1924
0.99889The PM2.5 concentration on each section line from 2014 to 2017 was extracted for spatial variation trend analysis, with the results shown in Figure6 and elaborated below. (1) The airports and the railway stations were vehicle gathering areas, and the automobile exhaust had a significant contribution to PM2.5. The scenic spots located in the center of Beijing with a large passenger flow, such as Tian’anmen, Heaven Temple, Prince Kung’s Mansion, and Summer Palace, had high PM2.5 concentration, whereas the Ming Tombs and the Badaling section of the Great Wall are far away from the city center with a mountainous terrain, and hence, their PM2.5 concentration was low. (2) A significant amount of dust was generated during the construction process, aggravating the environmental pollution. For example, the research site of the Jiugong Town Green-Belt Industry Development Project of Beijing Hongtu Jianye Investment Management Center had the PM2.5 concentration of 151 μg·m−3, 162 μg·m−3, 155 μg·m−3, and 153 μg·m−3 during the period from 2014 to 2017, respectively, which was attributed to dust generation during construction. (3) From 2014 to 2017, the Beijing Municipal Government gradually shut down large coal-fired power plants, and consequently the concentration of PM2.5 also dropped significantly. The Gaojing Thermal Power Plant was shut down in July 2014, where the annual mean PM2.5 concentration in successive years during the period from 2014 to 2017 was 60 μg·m−3, 58 μg·m−3, 50 μg·m−3, and 46 μg·m−3, respectively, showing a gradual decreasing trend. The Shijingshan Thermal Power Plant was officially shut down in March 2015, where the annual mean PM2.5 concentration from 2014 to 2017 was 60 μg·m−3, 61 μg·m−3, 52 μg·m−3, and 47 μg·m−3, respectively, with a dramatic decrease observed in 2016. The Guohua Beijing Thermal Power Plant was shut down in March 2015, and its annual mean PM2.5 concentration decreased significantly from 2014 to 2015, as shown in Figure 6. The Huaneng Beijing Thermal Power Plant was shut down in March 2017, and the PM2.5 concentration has not yet been improved significantly. (4) Petrochemical industries, automobile manufacturing, and printing industries also contributed to PM2.5. The PM2.5 concentration at Beijing Hualian Printing Co., Ltd. in 2014, 2016, and 2017 was 147 μg·m−3, 158 μg·m−3, and 160 μg·m−3, respectively, all in the high-concentration range of PM2.5. The PM2.5 concentration at the research sites of Beijing Automobile Group Co., Ltd. and Beijing Hyundai Motor Co., Ltd. was relatively high, exhibiting an annual increasing trend.Figure 6
PM2.5 concentration variation at each section line.


## 4.3. Planar Spatial Characteristic of PM2.5 Concentration

With Tian’anmen as the center of Beijing, a series of circular buffers was generated on the distribution map of annual mean PM2.5 concentration from 2014 to 2017, with the radius increasing from 4 km to 128 km by increments of 4 km, leading to a total of 32 rings shown in Figure7. According to China’s latest ambient air quality standard (GB 3095-2012), the annual mean concentration limit of national first-class air quality in Beijing is 35 μg·m−3. The area with PM2.5 concentration greater than 35 μg·m−3 in each buffer zone was calculated for each year, followed by taking a natural logarithm of the area and the radius and plotting a scatter plot of lnSr against ln r. Linear regression was performed on the data points, and the slope of the linear fitted line was considered as the radius dimension D.Figure 7
Calculation and analysis of the planar spatial characteristics of Beijing PM2.5.Results show that the radius dimensions of Beijing from 2014 to 2017 were 1.7608, 1.8234, 1.8548, and 1.8622, respectively, as shown in Table3. The variation trend can be observed from Figure 8 that, from 2014 to 2017, the PM2.5 mass concentration in Beijing decreased from the center to the surroundings, and the PM2.5 distribution becomes more uniform, indicating that the pollution control in Beijing had achieved initial success. In addition, the radius dimension in each year from 2014 to 2017 was less than 2, indicating that, for a PM2.5 concentration above 35 μg·m−3, the gradually declining spatial distribution from the center to the surroundings was nonlinear and nonuniform.Table 3
Radius dimensions of PM2.5 concentration from 2014 to 2017 in Beijing.
Year
Regression equation
Radius dimension
Coefficient of determination
2014
lnS(r) = 1.7608lnr + 1.4643
1.7608
0.9979
2015
lnS(r) = 1.8234lnr + 1.3892
1.8234
0.9986
2016
lnS(r) = 1.8548lnr + 1.359
1.8548
0.9987
2017
lnS(r) = 1.8522lnr + 1.3739
1.8522
0.9979Figure 8
Variation trend in the radius dimension of Beijing PM2.5 concentration from 2014 to 2017.


## 5. Results and Discussion

### 5.1. Wavelet Decomposition of the Temporal Distribution Characteristics of PM2.5 Concentration

The db wavelet transform was employed to PM2.5 concentration from 2014 to 2017. The original time series was decomposed into four layers, high-frequency coefficients, and the low-frequency coefficients of the fourth layer as shown in Figure9. The high-frequency coefficient is mainly composed of various disturbance noises and random fluctuations of abnormal mutations, reflecting the sudden changes and disturbances of PM2.5 concentration. In contrast, the low-frequency coefficient mainly consists of deterministic components, which reflects the variation characteristics of PM2.5 concentration, In general, the PM2.5 concentration was different in each year during the period from 2014 to 2017, but it could be roughly divided into three stages each year, namely, fluctuation, stationary, and re-fluctuation stages.Figure 9
Results of db6 wavelet analysis on PM2.5 concentration from 2014 to 2017. The original time series (s) was decomposed into four layers (d4 to d1) and low-frequency coefficients of the fourth layer (a4).The reconstructed low-frequency coefficients in the fourth layer of PM2.5 concentration (Figure10) show that the first stage was from January to April of each year, when the PM2.5 mass concentration showed significant fluctuation, experiencing a complex oscillation process of increasing, decreasing, re-increasing, and re-decreasing, and the maximum PM2.5 concentration decreased over time. The second stage, a stationary stage, was from May to October of each year, when the most significant feature was that the PM2.5 mass concentration was relatively stable, always at a low level. The third stage was from November to December of each year, when the PM2.5 concentration showed more significant intermonth fluctuation than the counterpart in the second stage, and the concentration increased dramatically in this stage in 2015 and 2016 compared with the second stage of the same years; moreover, this stage was similar to the first stage in showing complex cycles of increase, decrease, and re-increase in the overall variation trend of PM2.5 concentration.Figure 10
Distribution of reconstructed low-frequency coefficients in the fourth layer of PM2.5 concentration from 2014 to 2017 in Beijing.The abrupt nature of the PM2.5 time series usually indicates a status point where atmospheric pollution is quite severe. This study employed the db1 wavelet with good regularity to decompose the time series of PM2.5 concentration from 2014 to 2017 in Beijing into three layers, with the reconstruction result of the high-frequency coefficients in the first layer shown in Figure11. The points with large coefficient amplitude imply the abrupt change of PM2.5 concentration. As shown in Figure 11, there were only a few points with a large amplitude of reconstructed coefficients in the spring (March to May) and summer (June to August), which indicates of a low frequency of abrupt incidents of PM2.5 mass concentration. In contrast, there were a great number of points with a large amplitude of reconstructed coefficients in the autumn (September to November) and winter (December to next February), suggesting that these two seasons had a high frequency of abrupt incidents of PM2.5 mass concentration.Figure 11
Reconstruction of the high-frequency signals in the first layer of PM2.5 concentration from 2014 to 2017.


### 5.2. Periodicity Characteristics of the Multiscale Variation of PM2.5

To further identify the multiscale period of PM2.5 concentration change, the data were processed using the Morlet wavelet, which generated a contour map of the real part of the wavelet coefficients. The contour map represents a periodic oscillation of the PM2.5 concentration at different temporal scales. In particular, the larger (or smaller) the wavelet coefficient, the higher (or lower) the pollution concentration of atmospheric particulate matter, with 0 representing the abrupt change point. In contrast, the larger the absolute value of the wavelet coefficient, the more significant the change at the time scale.Figure12 shows that (1) the PM2.5 concentration from 2014 to 2017 in Beijing underwent changes at multiple temporal scales, i.e., 170–210 days, 110–140 days, 50–75 days, 30–50 days, and 10–20 days, namely, five change periods of different durations. However, periodic changes at different temporal scales showed a large difference in the temporal distribution characteristics. (2) The PM2.5 concentration during this period exhibited a very significant main change period of approximately 170–210 days, which ran through the entire timeline and exhibited a pattern of alternating high and low PM2.5 concentration. (3) The PM2.5 concentration had four distinct minor periods in the time domain, namely, 110–140 days, 50–75 days, 30–50 days, and 10–20 days. The four minor periods were subject to a large fluctuation in the time domain, indicating that the PM2.5 concentration series contained local features that varied significantly in minor periods.Figure 12
Contour map of the real part of the wavelet coefficients using Morlet wavelet analysis on PM2.5 concentration from 2014 to 2017 in Beijing.To explore the variation of the main periods of PM2.5 concentration in Beijing, the wavelet variance of the daily mean concentration of PM2.5 was calculated as shown in Figure13. It is possible to determine the relative intensity and major period of disturbances at different scales in the signals. The corresponding time scale of the highest peak intensity in the wavelet variance curve represents period 1, namely, the major period, whereas the others represent minor periods. The results show that there is an obvious major period of 190 days, indicating that the PM2.5 concentration in Beijing changed nearly every six months as a period. This is due to the strong winds in spring and autumn, but few winds in summer and winter in Beijing. The secondary period of 125 days also suggests that the variation in PM2.5 mass concentration in Beijing is closely season transitional. In addition, the PM2.5 concentration exhibits other relatively significant minor periods of 18 to 68 days. Given the minor period, it would be desirable to consider using the date of serious air pollution as the starting point to conduct short-term motor vehicle traffic restrictions nearly every 18 days so as to reduce the emission of atmospheric pollution sources.Figure 13
Wavelet variance curve of PM2.5 concentration from 2014 to 2017 in Beijing.


## 5.1. Wavelet Decomposition of the Temporal Distribution Characteristics of PM2.5 Concentration

The db wavelet transform was employed to PM2.5 concentration from 2014 to 2017. The original time series was decomposed into four layers, high-frequency coefficients, and the low-frequency coefficients of the fourth layer as shown in Figure9. The high-frequency coefficient is mainly composed of various disturbance noises and random fluctuations of abnormal mutations, reflecting the sudden changes and disturbances of PM2.5 concentration. In contrast, the low-frequency coefficient mainly consists of deterministic components, which reflects the variation characteristics of PM2.5 concentration, In general, the PM2.5 concentration was different in each year during the period from 2014 to 2017, but it could be roughly divided into three stages each year, namely, fluctuation, stationary, and re-fluctuation stages.Figure 9
Results of db6 wavelet analysis on PM2.5 concentration from 2014 to 2017. The original time series (s) was decomposed into four layers (d4 to d1) and low-frequency coefficients of the fourth layer (a4).The reconstructed low-frequency coefficients in the fourth layer of PM2.5 concentration (Figure10) show that the first stage was from January to April of each year, when the PM2.5 mass concentration showed significant fluctuation, experiencing a complex oscillation process of increasing, decreasing, re-increasing, and re-decreasing, and the maximum PM2.5 concentration decreased over time. The second stage, a stationary stage, was from May to October of each year, when the most significant feature was that the PM2.5 mass concentration was relatively stable, always at a low level. The third stage was from November to December of each year, when the PM2.5 concentration showed more significant intermonth fluctuation than the counterpart in the second stage, and the concentration increased dramatically in this stage in 2015 and 2016 compared with the second stage of the same years; moreover, this stage was similar to the first stage in showing complex cycles of increase, decrease, and re-increase in the overall variation trend of PM2.5 concentration.Figure 10
Distribution of reconstructed low-frequency coefficients in the fourth layer of PM2.5 concentration from 2014 to 2017 in Beijing.The abrupt nature of the PM2.5 time series usually indicates a status point where atmospheric pollution is quite severe. This study employed the db1 wavelet with good regularity to decompose the time series of PM2.5 concentration from 2014 to 2017 in Beijing into three layers, with the reconstruction result of the high-frequency coefficients in the first layer shown in Figure11. The points with large coefficient amplitude imply the abrupt change of PM2.5 concentration. As shown in Figure 11, there were only a few points with a large amplitude of reconstructed coefficients in the spring (March to May) and summer (June to August), which indicates of a low frequency of abrupt incidents of PM2.5 mass concentration. In contrast, there were a great number of points with a large amplitude of reconstructed coefficients in the autumn (September to November) and winter (December to next February), suggesting that these two seasons had a high frequency of abrupt incidents of PM2.5 mass concentration.Figure 11
Reconstruction of the high-frequency signals in the first layer of PM2.5 concentration from 2014 to 2017.


## 5.2. Periodicity Characteristics of the Multiscale Variation of PM2.5

To further identify the multiscale period of PM2.5 concentration change, the data were processed using the Morlet wavelet, which generated a contour map of the real part of the wavelet coefficients. The contour map represents a periodic oscillation of the PM2.5 concentration at different temporal scales. In particular, the larger (or smaller) the wavelet coefficient, the higher (or lower) the pollution concentration of atmospheric particulate matter, with 0 representing the abrupt change point. In contrast, the larger the absolute value of the wavelet coefficient, the more significant the change at the time scale.Figure12 shows that (1) the PM2.5 concentration from 2014 to 2017 in Beijing underwent changes at multiple temporal scales, i.e., 170–210 days, 110–140 days, 50–75 days, 30–50 days, and 10–20 days, namely, five change periods of different durations. However, periodic changes at different temporal scales showed a large difference in the temporal distribution characteristics. (2) The PM2.5 concentration during this period exhibited a very significant main change period of approximately 170–210 days, which ran through the entire timeline and exhibited a pattern of alternating high and low PM2.5 concentration. (3) The PM2.5 concentration had four distinct minor periods in the time domain, namely, 110–140 days, 50–75 days, 30–50 days, and 10–20 days. The four minor periods were subject to a large fluctuation in the time domain, indicating that the PM2.5 concentration series contained local features that varied significantly in minor periods.Figure 12
Contour map of the real part of the wavelet coefficients using Morlet wavelet analysis on PM2.5 concentration from 2014 to 2017 in Beijing.To explore the variation of the main periods of PM2.5 concentration in Beijing, the wavelet variance of the daily mean concentration of PM2.5 was calculated as shown in Figure13. It is possible to determine the relative intensity and major period of disturbances at different scales in the signals. The corresponding time scale of the highest peak intensity in the wavelet variance curve represents period 1, namely, the major period, whereas the others represent minor periods. The results show that there is an obvious major period of 190 days, indicating that the PM2.5 concentration in Beijing changed nearly every six months as a period. This is due to the strong winds in spring and autumn, but few winds in summer and winter in Beijing. The secondary period of 125 days also suggests that the variation in PM2.5 mass concentration in Beijing is closely season transitional. In addition, the PM2.5 concentration exhibits other relatively significant minor periods of 18 to 68 days. Given the minor period, it would be desirable to consider using the date of serious air pollution as the starting point to conduct short-term motor vehicle traffic restrictions nearly every 18 days so as to reduce the emission of atmospheric pollution sources.Figure 13
Wavelet variance curve of PM2.5 concentration from 2014 to 2017 in Beijing.


## 6. Conclusion

At present, many studies of urban PM2.5 mainly focus on data collection and processing. However, little work has been conducted to reveal how to take advantage of multimodal observation data to assess the spatiotemporal characteristics. This study proposes a spatiotemporal characteristics research method for urban PM2.5 concentration based on both airborne and ground monitoring mode. The proposed method was validated using the MODIS data and the daily mean PM2.5 concentration at ground monitoring stations from 2014 to 2017 in Beijing, and the following conclusions were drawn:(1)
Given the PM2.5 monitoring data of different frequencies and different granularities provided by the airborne and ground monitoring mode, this study utilized fractal dimensions and wavelet analysis to assess the spatiotemporal characteristics of urban PM2.5 concentration.(2)
It is feasible for us to use fractal dimensions to quantify the spatial distribution characteristics of PM2.5. The decrease results of boundary dimension indicate that the changes in PM2.5 concentration tend to be slower over time, and the fluctuation amplitude gradually becomes smaller. The increase results closer to 2.0 in the radius dimensions suggest that the PM2.5 distribution becomes more uniform, indicating that the pollution control in Beijing had achieved initial success.(3)
Wavelet analysis can detect different fluctuation periods of PM2.5 concentration at a long time series. The results reveal that PM2.5 concentration in Beijing is mainly subject to main periodic variation of 190 days, which was not reported before.Despite the findings of this study, further studies are required to gain deeper insight into some issues, especially the driving effect of meteorological factors on the formation of PM2.5. Therefore, in future studies, it is necessary to take more advantage of the data acquisition capability of the combined airborne and ground monitoring mode to not only obtain PM2.5 monitoring data but also design more meteorological indexes so that the meteorological factors affecting the dynamics of the spatiotemporal characteristics of PM2.5 can be simultaneously analyzed. In addition, the relationships between fractal dimension and wavelet analysis should be further discussed. Both wavelet transform and multifractal scaling can be used to characterize the multiscaling feature of spatial distributions of geographical phenomena so as to enhance the understanding of spatial variation of PM2.5.



---

*Source: 8091515-2020-02-10.xml*