The Spatial-Temporal Differentiation of Aerosol Optical Properties and Types in the Beijing–Tianjin–Hebei Region Based on the Ecological Functional Zones

Abstract

Atmospheric aerosol pollution has seriously affected the ecological environment and human health in recent years. There are great differences in aerosol optical properties and types due to the influence of environmental conditions, meteorology, industrial and agricultural activities, and other factors of each ecological functional zones. Using MODIS aerosol products (including MCD19A2 and MOD04_3K), this study discussed the temporal and spatial distribution of aerosol optical depth (AOD), Ångström wavelength index (AE) and aerosol types in the Beijing–Tianjin–Hebei region (BTH region) based on the ecological functional zones from 2015 to 2020. The results showed as follows: (1) The AOD in BTH region showed an obviously spatial pattern of low in the north and high in the south, while the spatial pattern of AE was opposite to that of AOD. In addition, the dominant aerosol type of the north part was clean aerosol, the dominant aerosol type of the middle part was biomass burning or urban-industrial aerosol, while that of the other part was mixed aerosols. (2) The seasonal changes of AOD and AE indexes in each ecological functional area had obvious seasonal changes, and the AOD and AE values were highest in summer. At the same time, the proportion of biomass combustion or urban industrial aerosol was the highest in summer. (3) The ecological functional areas with fewer human activities were dominated by clean aerosols, with lower AOD and higher AE value. The ecological functional areas dominated by cities were dominated by mixed aerosols, with higher AOD. The ecological functional areas dominated by agriculture and heavy industry were dominated by biomass combustion or urban industrial aerosols, with the largest AOD. (4) Compared with 2015, the average AOD of each ecological functional area decreased significantly to 2020, and biomass combustion or urban industrial aerosols changed to mixed aerosols.

1. Introduction

Atmospheric aerosol is the suspension of solid or liquid particles in the atmosphere derived from natural sources (such as sand and dust) and artificial sources (such as industrial waste gas emissions and coal combustion) [1]. Atmospheric aerosols can affect the ecological environment and human health. For instance, they can disturb the Earth’s radiation balance, change precipitation rate, and affect water cycles [2,3,4]. More importantly, aerosol particles can seriously threaten human health by easily spreading a variety of toxic and harmful substances [5].

Accurately monitoring the temporal and spatial distribution of atmospheric aerosols is the key issue to explore the ecological, environmental, and human health effects of aerosols. The methods mainly include ground monitoring, satellite remote sensing inversion, and model simulation. Among them, satellite remote sensing has become an economical and feasible method for aerosol monitoring by obtaining comprehensive coverage and high spatial resolution of aerosol optical property data. Among them, obtaining comprehensive coverage and high spatial resolution of aerosol optical property data by satellite remote sensing inversion has become an economical and feasible method for aerosol monitoring [6,7,8]. Aerosol optical depth (AOD), Ångström wavelength index and aerosol type are three important indicators in reflecting aerosol content, aerosol particle size, and aerosol pollution source [9,10].

The ecological function zone is divided according to the characteristics of ecological environment, the sensitivity of ecological environment, and the differences and similarities in ecological service functions [11]. Therefore, there are great differences in the natural ecological environment, human activity, and industrial structure in each ecological function zones. This also leads to the obvious spatial-temporal differentiation of aerosol optical properties and types in each ecological function zone. However, the study of aerosol pollution has mainly been on the spatial-temporal distribution and driving factors of aerosol optical depth from the administrative district or urban agglomeration [12]. There are few studies on monitoring aerosol pollution based on ecological function zones. Exploring the characteristics of aerosol optical properties and types in different ecological functional zones could reveal the cause of air pollutants and help to verify the implementation effect of the air pollution control policy. In summary, it is necessary to analyze the spatial-temporal variances in atmospheric aerosol optical properties and types in different ecological functional zones.

With the rapid economic development, industrialization and urbanization, the ecological environment has been severely damaged in the Beijing–Tianjin–Hebei (BTH) region, especially atmospheric aerosol pollution. In 2013, the BTH region has been severely polluted by PM_2.5 (the daily concentration of PM_2.5 is higher than 150 μg/m³) for more than 270 days [13]. In addition, BTH region has obviously spatial heterogeneity containing a variety of land cover types, so there are naturally different ecological function zones. However, it is unclear the temporal and spatial differentiation characteristics of aerosol optical properties and types in different ecological function zones. Therefore, this study analyzes the spatial-temporal diversity analysis of aerosol attributes and types in different ecological function areas during the 13th Five Year Plan period, in order to reveal the causes and driving mechanisms of atmospheric pollution in different ecological functional zones, and then formulate reasonable and effective measures for the prevention and controlling of spatial pollution.

2. Materials and Methods

2.1. Overview of the Study Area

The Beijing–Tianjin–Hebei region (36–42° N, 112–121° E) is located in the North China Plain, the eastern coast of China. The terrain of the BTH region is high in the northwest and low in the southeast. The climate of the BTH region is typical warm temperate continental monsoon climate zone with hot, rainy summers and cool, dry winters. The BTH region could provide abundant anthropogenic aerosols due to more industrial and agricultural production and higher population density. In addition, the BTH region also receives regional transmission of aerosol particles due to its closeness to the sand source [14].

The ecological functional zones of BTH region refers to the general office of the Hebei Provincial People’s Government on the issuance of Hebei Province’s support for the construction of the BTH ecological environment. According to the notice of the “14th Five-Year Plan” of the district, the whole area of BTH is divided into the six regions (Figure 1c), following as: Ecological transition zone around Beijing and Tianjin (ET), Bashang Plateau Ecological Protection Zone (BP), Yanshan-Taihang Mountain Ecological Conservation Zone (YT), Low Plain Ecological Restoration Zone (LP), Coastal Ecological Protection Zone (CE), and Urban Core Functional Zone (UC).

2.2. Data Sources and Methodology

2.2.1. Data Sources and Preprocessing

The AOD data used in this paper are all from the Moderate Resolution Imaging Spectrometer (MODIS) mounted on the EOS series Aqua and Terra satellite detectors, including MCD19A2 and MOD04_3K products (https://ladsweb.modaps.eosdis.nasa.gov/, accessed on 15 January 2022). Since 2000, the aerosol products of MODIS have been widely used in aerosol spatial-temporal distribution, air pollution monitoring, and environmental pollution dynamic change monitoring [15]. Among them, the daily MCD19A2 product adopts the new MAIAC aerosol algorithm by the combination of Terra satellite and Aqua satellite with 1 km spatial resolution [16]. This study selects the dataset of Optical_Depth_550 dataset in the MCD19A2 product to obtain the 550 nm_AOD between 2015 and 2020. The MOD04_3K product generated by the Terra satellite uses the dark target (DT) aerosol algorithm with a spatial resolution of 3 km. This study obtained 470nm_AOD and 660nm_AOD using the “Corrected_Optical_Depth_Land” dataset of the MOD04_3K product in 2015 and 2020.

In this study, the AOD dataset was preprocessed to extract the dataset, filter QA value and calculate the mean band value by Envi5.3 + idl. The MCD19A2 dataset was extracted by first obtaining the multi-value data of the 550 nm band AOD that passed the quality control in the raw data, converting it to the actual value, and then getting the daily mean value 550nm_AOD. Through image stitching, projection conversion, and image cropping, the daily AOD data of the BTH region was obtained, which was used to obtain the annual and monthly average AOD values of the BTH area. The MOD04_3K data mainly used the dataset of “Corrected_Optical_Depth_Land_470nm” and “Corrected_Optical_Depth_Land_660nm”. Finally, all AOD data projection information uses the WGS1984 coordinate system.

2.2.2. Methodology

Standard Error Ellipse Analysis

In this paper, the standard error ellipse is used to analyze the distribution characteristics of AOD (AOD > 0.48) and the AE (AE < 1.2) in the BTH region. Standard error ellipse (directional distribution) is a common method proposed by Lefever [17] to study the distribution characteristics of partial or overall regional data sets. Calculate the average center of the data set distribution, i.e., the ellipse center, which is represented by the following equation.

Center coordinates:

$$SD{E}_x=\sqrt{\frac{{{\displaystyle \sum }}_{i=1}^n{\left(x_i-\overline{X}\right)}^2}{n}},SD{E}_y=\sqrt{\frac{{{\displaystyle \sum }}_{i=1}^n{\left(y_i-\overline{Y}\right)}^2}{n}}$$

(1)

where $x_i$ and $y_i$ represent the coordinates of the element point i, $\overline{X}$ and $\overline{Y}$ represent the average center of the element, and n is the element number.

The long semi-axis of the ellipse represents the direction of the data distribution, and the short semi-axis represents the data distribution range. Conversely, the directionality is less noticeable if the long and short semi-axes are closer [18]. The direction of the ellipse is represented by Equation (2).

Ellipse direction:

$$\begin{array}{c}tan\theta =\frac{A+B}{C},A=\left({\displaystyle \sum }_{i=1}^n{\tilde{x}}_i^2-{\displaystyle \sum }_{i=1}^n{\tilde{y}}_i^2\right)\\ B=\sqrt{{\left({\displaystyle {\displaystyle \sum }_{i=1}^n}{\tilde{x}}_i^2-{\displaystyle {\displaystyle \sum }_{i=1}^n}{\tilde{y}}_i^2\right)}^2+4{\left({\displaystyle {\displaystyle \sum }_{i=1}^n}{\tilde{x}}_i{\tilde{y}}_i\right)}^2},\\ C=2{\displaystyle {\displaystyle \sum }_{i=1}^n}{\tilde{x}}_i{\tilde{y}}_i\end{array}$$

(2)

where ${\tilde{x}}_i$ and ${\tilde{y}}_i$ are the difference between the mean center and the xy coordinates, respectively, and the direction represents the rotation of the long axis measured clockwise from the vertex.

When the features have a spatially normal distribution (that is, they are densest at the center of the ellipse and gradually become sparser as they expand from the center to the edge). The range of the first level standard deviation can include the centroid of about 68% of the total input elements. The second standard deviation range will cover about 95% of the total elements, while the third standard deviation range will cover about 99% of the total elements. In comparison, the third level of the normal deviation range will cover the centroids of about 99% of the total number of features [19]. The standard deviation of the x-axis and the y-axis is represented by Equation (3). The first-level standard deviation is used in this paper, and the meaning of each parameter is consistent with Equations (1) and (2).

x-axis standard deviation:

$${\sigma }_x=\sqrt{2}\sqrt{\frac{{{\displaystyle \sum }}_{i=1}^n{\left({\tilde{x}}_{i}cos\theta -{\tilde{y}}_{i}sin\theta \right)}^2}{n}}$$

y-axis standard deviation:

$${\sigma }_y=\sqrt{2}\sqrt{\frac{{{\displaystyle \sum }}_{i=1}^n{\left({\tilde{x}}_{i}sin\theta -{\tilde{y}}_{i}cos\theta \right)}^2}{n}}$$

(3)

Ångström Index

AE, a key optical parameter of aerosols, was used to quantitatively characterize aerosol size. The higher AE indicates finer aerosols particles [20]. The calculation method of the wavelength index is:

$$AE=-\frac{\ln \left[\tau \left({\lambda }_1\right)\right]/\ln \left[\tau \left({\lambda }_2\right)\right]}{\ln \left({\lambda }_{1/}{\lambda }_2\right)}$$

(4)

In the formula: $\tau \left(\lambda \right)$ is the aerosol optical thickness of the $\lambda$ band, ${\lambda }_1$ and ${\lambda }_2$ are the wavelengths of two different bands; the more significant the wavelength interval, the more realistic the $\alpha$ value is. In this paper, the AOD data of the 470 nm and 660 nm bands of the MOD04_3K data are selected, and the value of $\alpha$ is calculated using ArcGIS.

Dynamic Analysis of AOD and AE

The difference method has the advantage of being concise and intuitive in analyzing the dynamic change characteristics [21], and its formula is as follows:

$$\mathsf{\Delta }x=x_2-x_1$$

(5)

$\mathsf{\Delta }x$ is the change of AOD or AE, $x_2$ was the value of AOD or AE in 2020, and $x_1$ was the value of AOD or AE in 2015. Among them, the change of aerosol optical thickness is set as four categories, namely rising ($\mathsf{\Delta }x$ > 0), unchanged ($\mathsf{\Delta }x$ = 0), significant falling (−0.1 < $\mathsf{\Delta }x$ < 0), and falling ($\mathsf{\Delta }x$ < −0.1), AE. The index changes are set to three categories, namely rising ($\mathsf{\Delta }x$ > 0), unchanged ($\mathsf{\Delta }x$ = 0), and falling ($\mathsf{\Delta }x$ < 0), to study the dynamic characteristics of AOD and AE indices in 2015 and 2020.

Aerosol Type Analysis

Referring to the method of K. Raghavendra Kumar [22], AOD_550nm and AE_470-660nm are used to classify the types of aerosols into the following five categories: (1) Values with AOD < 0.2 and AE > 1.0 indicate the continental clean (CC) aerosols representing the clean background conditions. (2) The marine (MA) type of aerosols are recognized between AOD < 0.2 and AE < 0.9; (3) AOD > 0.3 and AE > 1.0 represent the excessively turbid atmosphere dominated by biomass burning/urban-industrial (BU) aerosol type. (4) AOD > 0.5 and AE < 0.7 are characterized by desert dust (DD) type aerosols. (5) The remaining cases not belonging to any of the above categories are identified as mixed (MX) type or undetermined aerosols.

3. Results

3.1. Spatial-Temporal Variation of Aerosol Optical Depth (AOD)

3.1.1. Spatial Pattern of Annual AOD

In general, the spatial distribution of AOD was similar in 2015 and 2020, with a low in the north and a high in the south (Figure 2a,b). In addition, the average annual AOD (0.308) in 2020 is obviously lower than that in 2015 (0.408). To further explore the variations in spatial distribution of the annual average AOD from 2015 to 2020, standard error ellipse analysis of the highly polluted areas (AOD > 0.48) was calculated (Figure 2c). It was found that the long axis and short axis of the error ellipse in 2020 (114 km, 79 km) were shorter than those in 2015 (263 km, 188 km), and the centers of the two ellipses were 100 km apart. The results indicated that the range of the high-value area in 2020 was lower than that in 2015, and the center of the high-value area (the center of the ellipse) was moving southward.

Through the analysis of AOD in different ecological function zones, the AOD value of the ET and LP zones was obviously higher than that in other zones, with an annual average value of about 0.5. Moreover, the lowest AOD was located in the BP and YT zones, with an annual average AOD value of 0.2–0.3 (

Table 1. The change (%) of annual average AOD in each ecological function zone from 2015 to 2020.

Ecological Function Zone	Annual Average AOD in 2015	Annual Average AOD in 2020	Reduction of Annual Average AOD from 2015 to 2020 (%)
ET	0.655	0.442	32.52%
BP	0.174	0.155	10.92%
YT	0.299	0.234	21.74%
LP	0.644	0.499	22.52%
CE	0.485	0.374	22.89%
UC	0.445	0.316	28.99%
BTH	0.412	0.308	25.24%

). In addition, the annual average AOD value of each ecological function zone in 2020 was lower than that in 2015, especially in ET and UC zones (Table 1). However, there was little difference in the AOD values of the BP area from 2015 to 2020.

3.1.2. Temporal Pattern of Monthly AOD

On the whole, the monthly variations in AOD in different ecological functional zones were similar, with the highest AOD value occurring in July and the lowest average AOD value occurring in February and December (Figure 3). The lowest values of monthly AOD were located in the BP and YT zones, while it was higher in ET and LP zones. Compared with 2015, the values of monthly AOD in 2020 were almost lower than that in 2015, especially in July and August.

3.2. Spatial and Temporal Distribution of Aerosol Particle Size

3.2.1. Spatial Distribution of Annual AE

There was similar spatial pattern of AE indices in the BTH region between 2015 and 2020 (i.e., a high in the north while a low in the south) (Figure 4a,b). To further explore the spatial variation of the AE, this study conducts a standard deviation ellipse analysis on the regions where the AE is less than 1.2 in 2015 and 2020 (Figure 4c). It was found that there were few changes in the values of the long axes, short axes, and ellipse centers, indicating that the spatial pattern of AE values was similar between 2020 and 2015. From the difference analysis, it could be seen that the areas with declining AE values are mainly concentrated in the central and southern parts (Figure 4d).

From the perspective of ecological function zones, the AE values in the BP and YT zones were higher than that in other ecological function zones (

Table 2. The proportion of AOD changes of ecological function zone between 2015 and 2020.

Ecological Function Zone	Annual Average AE in 2015	Annual Average AE in 2020	Reduction of Annual Average AE from 2015 to 2020 (%)
ET	1.13	1.046	−7.43%
BP	1.344	1.303	−3.05%
YT	1.329	1.356	2.03%
LP	0.971	0.965	−0.62%
CE	1.044	0.994	−4.79%
UC	1.188	1.138	−4.21%
BTH	1.215	1.209	−0.49%

). In addition, the AE values in most zones were above 1.3, indicating that the aerosol particle size of these zones was fine-grained. In addition, the AE values in the LP and CE zones were relatively low, with a value generally less than 1. However, the AE value in the YT zone was higher in 2020, which indicated that the proportion of fine particles had increased since 2015. The results showed that the proportion of coarse particles had increased in all the ecological function zones apart from the YT zone.

3.2.2. Temporal Pattern of Monthly AE

On the whole, the monthly AE values in the CE, LP, and ET zones were obviously higher in the second half of the year than that in the first half of the year for both 2015 and 2020 (Figure 5). In contrast, the monthly AE value in each ecological functional zone in the second half of 2020 was lower than that in 2015, while the first half of the year was just the opposite.

3.3. Temporal and Spatial Variation of Aerosol Types

3.3.1. Spatial Distribution of Aerosol Types

The spatial distribution of aerosol types was similar between 2015 and 2020 (Figure 6a,b). In the BP and north of YT zones, CC aerosols were the leading aerosol type; BU aerosols were the dominant aerosol type of UC and ET regions, and MX aerosols were the dominant aerosol type of LP and CE zones. In addition, the area where MX aerosols were converted to CC aerosols was mainly located in the southern YT zone, and the area where BU aerosols were converted to MX aerosols was mainly located in the UC zone (Figure 6c).

To analyze the distribution of aerosol types in each ecological functional area, the proportions of different aerosol types in 2015 and 2020 were counted and plotted (Figure 7). The results showed that the BU aerosols in the BTH region declined obviously from 37% in 2015 to 23% in 2020. In addition, the BU aerosols also declined from 88%, 53%, and 58% in 2015 to 65%, 16%, and 40% in 2020 in the ET, UC, and CE zones, respectively, whereas the MX aerosols simultaneously increased from 12% to 35%, from 45% to 70%, and from 41% to 58% in these regions. In addition, the CC aerosols prevail over the BP region, slightly increasing from 89% to 93% since 2015.

3.3.2. Monthly Pattern of Aerosol Types

Based on the ecological functional areas, we have counted the changes of aerosol types in each ecological functional area in each month (Figure 8). In the ET, YT, LP, and CE zones, the proportion of MX aerosols was higher in the first half of the year, and the proportion of BU aerosols was higher in the second half of the year. In addition, the proportion of BU aerosols in the second half of 2020 declined significantly compared with 2015, while the proportion of MX aerosols increased. In the BP area, the proportion of CC aerosols is relatively high in each month, and the proportion of CC aerosols in winter months (December, January, and February.) declined from 2015 to 2020, but the proportion of MX and MA aerosols obviously rose. Furthermore, the BU, MX, and CC aerosols are the dominant aerosol types for most months in UC region. Compared with 2015, the proportion of BU aerosol is lower, MX and CC aerosol is higher in the second half of 2020, the proportion of BU and CC aerosol declined and only the proportion of MX aerosol rose in the first half of the year.

4. Discussion

4.1. Spatial-Temporal Differentiation of Aerosol Optical Properties and Types

In general, the atmospheric aerosol pollution in the BTH region was increasing from north to south. The areas that is high AE value and low value of AOD were located on the BP and the YT zones, where the dominated aerosol type was CC type. This phenomenon indicated that the degree of aerosol pollution in the BP and YT zones was relatively slight, due to small population density, high vegetation coverage and fewer aerosol emission sources. The higher AOD value and AE being close to 1 in the UC zone implied that its dominant aerosol type was BU aerosols, indicating that atmospheric aerosol pollution was more serious in the UC zone because of a large number of aerosol sources. In addition, the LP zone seemed the major aerosol source because there had the highest AOD, AE value was about 1, and MX was the dominant aerosol type. At the same time, the value of AE where is around the CE area is low, because the dominant aerosol type is sea-salt aerosol whose particle size is larger than that of other aerosol types [23].

Compared with 2015, the decrease in AOD and increase in AE, and the conversion of the dominant aerosol type which was from biomass combustion or urban industrial aerosol to mixed aerosol and clean aerosol were occurred in 2020, especially in UC, CE, and ET zones. This might be due to a series of measures during the 13th Five-Year Plan period, such as the control of bulk coal in winter, the supervision of motor vehicles and the increasing proportion of receiving external power transmission [24].

4.2. Monthly Variation of Aerosol Optical Properties and Types

The value of AOD and AE were largest, and the dominant aerosol type of each ecological function zone was BU type in July. On the one hand, the incidence of photochemical reactions and the formation of secondary inorganic aerosols would be promoted in hot and rainy summer, and therefore the aerosol concentration were also obviously improved [25]. On the other hand, the crop straw burning in June contributes obviously to the aerosol concentration of July which reached the maximum value of the year [26]. In addition, the value of AE was small by supplying a large number of coarse particles thought frequent sandstorms in spring. Simultaneously, aerosol pollution is more severe in winter and spring which are the heating seasons in the north of China. Except for a slight rebound in the value of AOD in the winter of 2020(the specific seasonal distribution is shown in Appendix A), the atmospheric aerosol pollution in other seasons was more improving than that in the same period of 2015. This was due to the prohibition of straw burning and the comprehensive management of critical industries in summer. In addition, the degree of aerosol pollution was highly alleviated in 2020, due to the measurement of reducing coal burning.

5. Conclusions

By analyzing the spatial-temporal pattern of aerosol optical properties and types in the BTH region between 2015 and 2020, it was found that the AOD was low in the north and high in the south, but that of the AE was the opposite; the dominant aerosol type was CC type in the north, BU type in the middle, and MX type in the south. In general, the atmospheric aerosol pollution of the northern BTH region was relatively lighter than that in the central and southern region. Simultaneously, the value of both AOD and AE was the highest in July. This phenomenon might be due to the high temperature and rain in July, which has promoted the growth of hygroscopic aerosol particles and the generation of secondary aerosols from photochemical reactions. Meanwhile, straw burning was frequently occurred during this period, which would generate a large amount of black carbon aerosols and eventually increase the aerosol concentration. Compared with 2015, the AOD was obviously declined, and the dominant aerosol type was changed from BU type to MX or CC type in 2020. The significant improvement in atmospheric environment in 2020 might be due to a series of strict environmental protection measures such as the “Ten Atmospheric” measures including the treatment of loose coal in winter, increasing non-fossil energy supply, and removing production capacity of traditional heavy industries.

There was also obviously a spatial-temporal pattern of ecological functional zones. In BP and YT zones with high vegetation coverage and less human activities, the AOD was low, the AE value was high, and the dominant aerosol type was CC aerosol. In other ecological functional zones, the value of AOD was high while the dominant aerosol types were MX and BU types. From 2015 to 2020, the declining AOD, increasing AE, and the conversion from BU type to MX type in ET, UC, and CE zones. This might be caused by the implementation of severe air pollution control and prevention measures to avoid serious air pollution in winter.