Spatial and temporal changes in SO₂ regimes over China in the recent decade and the driving mechanism

The spatial and temporal changes in SO₂ regimes over China during 2005 to 2016 and their associated driving mechanism are investigated based on a state-of-the-art retrieval dataset. Climatological SO₂ exhibits pronounced seasonal and regional variations, with higher loadings in wintertime and two prominent maxima centered in the North China Plain and the Cheng-Yu District. In the last decade, overall SO₂ decreasing trends have been reported nationwide, with spatially varying downward rates according to a general rule – the higher the SO₂ loading, the more significant the decrease. However, such decline is in fact not monotonic, but instead four distinct temporal regimes can be identified by empirical orthogonal function analysis. After an initial rise at the beginning, SO₂ in China undergoes two sharp drops in the periods 2007–2008 and 2014–2016, amid which 5-year moderate rebounding is sustained. Despite spatially coherent behaviors, different mechanisms are tied to North China and South China. In North China, the same four regimes are detected in the time series of emission that is expected to drive the regime of atmospheric SO₂, with a percentage of explained variance amounting to 81 %. Out of total emission, those from the industrial sector dominate SO₂ variation throughout the whole period, while the role of household emission remains uncertain. In contrast to North China, SO₂ emissions in South China exhibit a continuous descending tendency, due to the coordinated cuts of industrial and household emissions. As a result, the role of emissions only makes up about 45 % of the SO₂ variation, primarily owing to the decoupled pathways of emission and atmospheric content during 2009 to 2013 when the emissions continue to decline but atmospheric content witnesses a rebound. Unfavorable meteorological conditions, including deficient precipitation, weaker wind speed and increased static stability, outweigh the effect of decreasing emissions and thus give rise to the rebound of SO₂ during 2009 to 2013.

1 Introduction

In the recent decade, air pollution has persistently plagued China, especially in leading economic and densely populated areas (Chan and Yao, 2008; Ma et al., 2012; Chai et al., 2014). In China, environmental protection agencies identify six major pollutants of concern, including sulfur dioxide (SO₂), nitrogen dioxide (NO₂), ozone (O₃), carbon monoxide (CO), fine particulate matter (PM_2.5) and coarse particulate matter (PM₁₀). Then, values of the six pollutants are transformed into a single number called the Air Quality Index (AQI) for effective communication of air quality status and the corresponding health impact (MEPC, 2012).

SO₂ is one of the six major pollutants in China (Ren et al., 2017). It is harmful to human health, affecting lung function, worsening asthma attacks and aggravating existing heart disease (WHO, 2018). It also leads to the acidification of the atmosphere, and the formed sulfate aerosol is one of the most important components of fine particles in cities (Meng et al., 2009). Overall, SO₂ is a key influencing factor for atmospheric pollution, and it poses great threats to life, property and environment (Wang et al., 2014).

Compared to airborne and ground-based remote sensing, satellite platforms permit near-global coverage on a continuing and repetitive basis, enabling quick and large-scale estimation of pollution patterns (Yu et al., 2010). Since the world's first weather satellite TIROS-I launched in 1960, satellites have become a crucial part of Earth's observations and practical applications (Yu et al., 2010). Till now, SO₂ has been measured globally by several operational satellite instruments, such as OMPS (Zhang et al., 2017), GOME-2 (Munro et al., 2006; Rix et al., 2012) and the Ozone Monitoring Instrument (OMI) (Lee et al., 2011; Li et al., 2013; Theys et al., 2015).

With the aid of satellite data, in the past decade, various attempts have been made to explore the variation of SO₂ loadings in China. Lu et al. (2010) report that total SO₂ emissions in China have increased by 53 % from 2000 to 2006, followed by a growth rate slowdown and the start of a decrease. Li et al. (2010), Yan et al. (2014), and Zhang et al. (2012) all highlight the prominent reduction of SO₂ during 2007 and 2008, as a consequence of the widespread deployment of flue-gas desulfurization and the strict control strategy implemented for preparation of the 2008 Olympic Games. Throughout the past decade, 90 % of the locations in China have shown a decline in SO₂ emissions, as highlighted by Koukouli et al. (2016). Such widespread declines are ascribed to effective air quality regulations enforced in China (van der A et al., 2017). Furthermore, Krotkov et al. (2016) and C. Li et al. (2017) both compared the sulfurous pollution in China and India, and pointed out their opposite trajectories. Since 2007, emissions in China have declined while those in India have increased substantially. Nowadays, India is overtaking China as the world's largest emitter of anthropogenic SO₂. In addition, several studies conducted analyses on SO₂ in sub-regions of China, for example Jin et al. (2016), Lin et al. (2012), Wang et al. (2015) and Su et al. (2011). All these studies contributed to a better understanding of SO₂ changes in China. However, there are still key issues to be addressed. First, with the pace of considerable progress made on SO₂ retrieval, updated data products are now available to accurately derive recent SO₂ variations in China. Second, although the general decreasing tendency has been revealed, the specific spatial and temporal regimes remain unclear. Does the SO₂ decrease monotonically, or is there a complicated oscillation? How similar/different are SO₂ variations in different parts of the country? Third, there is more to be learned about the driving mechanisms that govern SO₂ variations. Previous studies have mainly focused on the impact from amounts of emission. However, the SO₂ content is not only dependent on emissions but also on atmospheric conditions. Therefore, how large is the influence of atmospheric variability on the variation of SO₂?

The overall goal of this study is to quantify the spatial and temporal changes in SO₂ regimes over China in the last decade and to disclose the driving mechanism, based on a new generation of SO₂ retrieval datasets (Theys et al., 2015). Figure S1 in the Supplement labels the provinces of China. The paper is organized as follows. Section 2 describes the new SO₂ product, and emission inventories and atmospheric data are introduced. In Sect. 3, we evaluate the general patterns of SO₂, including mean distribution, seasonality and long-term trends. Subsequently, Sect. 4 identifies the specific regimes of SO₂ variability and the associated driving mechanisms. Finally, concluding remarks and future directions are presented in Sects. 5 and 6, respectively.

2 Data

2.1 SO₂ VCD retrievals

The OMI is one of four sensors onboard the Aura satellite launched in July 2004 (Levelt et al., 2006). In recent years, the Belgian Institute for Space Aeronomy (BIRA) and cooperators have developed an advanced differential optical absorption spectroscopy (DOAS) scheme to improve the retrieval accuracy of SO₂ in the troposphere. A SO₂ vertical column product is generated based on the algorithm applied to OMI-measured radiance spectra (Theys et al., 2015). The retrieval scheme is a based on a DOAS approach, including three steps: (1) a spectral fit in the 312–326 nm range (other fitting windows are used for volcanic scenarios but are not relevant for this study), (2) a background correction for possible bias on retrieved SO₂ slant columns, (3) a conversion into SO₂ vertical columns through radiative transfer air mass factor calculation, accounting for the SO₂ profile shape (from the IMAGES chemistry-transport model), geometry, surface reflectance and clouds.

Compared to the BRD OMI NASA SO₂ product, the BIRA retrievals proved to be better in terms of both noise level and accuracy. The BIRA product is also fully characterized (errors, averaging kernels, etc.). The improved OMI PCA SO₂ product of NASA shows similar performance and long-term trends to the BIRA product. The BIRA SO₂ product has been validated in China with long-term MAX-DOAS data (Theys et al., 2015; Wang et al., 2017).

The dataset is made available on a 0.25 and 0.25^∘ regular latitude–longitude grid over the rectangular domain 70–140^∘ E, 10–60^∘ N, and covers the period of 2005 to 2016 at a monthly interval. In addition, a cloud screening is applied to remove measurements with a cloud fraction of more than 30 %. Other details can be found in Theys et al. (2015).

Given that missing values are often presented in the satellite-retrieved product due to the limitations of retrieval algorithms in adverse environments, it is necessary to evaluate the availability of monthly SO₂ data relative to the entire period. As mapped in Fig. S2, there appears to be a substantial fraction of data gaps in West and Northeast China, especially in the winter half-year. This can be attributed to snow cover surfaces and high solar zenith angles, which invalidate the measurability. As a result, it may be problematic when sampling West and Northeast China. In contrast, the completeness across eastern parts of China is generally more than 80 % regardless of the season, sufficient for inferring the spatial and temporal structures. In what follows, the analysis is mainly confined to East China to avoid issues related to missing data.

2.2 Emission inventory

The SO₂ emissions at national and provincial level are collected from the China Statistical Yearbook on Environment (NBS, 2003–2015), which is compiled jointly by the National Bureau of Statistics and Ministry of Environmental Protection. It is an annual statistics publication, with industrial and household emissions listed separately. Currently, this publicly available dataset spans the period from 2003 to 2015, covering 31 provinces in China other than Taiwan, Hong Kong and Macau. Industrial emissions refer to volume of SO₂ emission from fuel burning and production processes in premises of enterprises for a given period, while household emissions are calculated on the basis of consumption of coal by households and the sulfur content of coal. Note that power generation is incorporated into industrial sources and emissions from transportation sources are not reported. This emission inventory released in the official yearbook (OYB for short) has been cited or used in several works, i.e., M. Li et al. (2017), Yan et al. (2017), and Hou et al. (2018).

Since a credible emission inventory is the key foundation of this study, the Multi-resolution Emission Inventory for China (MEIC) developed by Tsinghua University (M. Li et al., 2017; Zheng et al., 2018) has been adopted to verify the OYB inventory as well as to corroborate our findings. The MEIC is a bottom–up emission inventory model including more than 700 anthropogenic sources and then aggregated into five sectors: power, industry, residential, transportation and agriculture. Unlike the OYB estimate, emissions from power plants in MEIC are considered to be a single sector and are presented separately. Here, we use province-level emissions from 2003 to 2015, together with the monthly gridded emissions at 0.25^∘ × 0.25^∘ horizon resolution for the years 2008, 2010, 2012, 2014 and 2016. To be in line with the OYB inventory, the transportation and agriculture sectors are excluded when calculating the summed emission, and the power sector is folded into the industrial sector.

Figure 1 (a) National total SO₂ emissions estimated by OYB (solid black), MEIC (solid red), REAS (dashed blue) and Zhao (dashed orange) between 2003 and 2015. (b) Scatter diagrams and regression lines for the OYB estimate (x axis) against the other three products (y axis). (c) The province-by-province correlations between the OYB and MEIC products, with the significance levels of 0.1 and 0.05 being marked by open and filled triangles, respectively.

Figure 1 compares the OYB and MEIC emission inventories in terms of both national and regional scales. In addition, the other two candidates on national annual totals including REASv2 (Kurokawa et al., 2013) and Zhao (Xia et al., 2016) are overlaid. Figure 1a shows that considerable differences exist with regards to the magnitude among the four datasets and, in particular, OYB emissions are generally lower than those deduced from other inventories. However, their temporal variations are characterized in a very similar manner. As further illustrated in the scatter plot of OYB against the other three (Fig. 1b), highly linear clustered markers with correlations above 0.92 confirm such temporal consistency. On an even smaller regional scale, as shown in Fig. 1c, high degrees of correspondence between OYB and MEIC overwhelm the whole of East China, with most correlations exceeding the 0.05 significance level. In comparison, West China features relatively less agreement, but it is not a major concern in this study.

In short, all the datasets capture coherent temporal behaviors, despite the spread in their magnitudes. We emphasize that this study is centered on the fluctuation patterns rather than the magnitude itself. Therefore, the above evidence justifies the use of the OYB dataset in the following text. Meanwhile, in order to test whether results were robust to using a different dataset, all analyses have been repeated using the MEIC inventory.

2.3 Meteorological fields

The large-scale meteorological conditions are taken from Japanese 55-year Reanalysis (JRA-55) data, prepared by the Japan Meteorological Agency (Kobayashi et al., 2015; Harada et al., 2016). The variables analyzed include total column precipitable water, horizontal wind and temperature at pressure levels.

3 General patterns of SO₂ over China

3.1 Mean distribution

Based on 12 years of SO₂ column data over China, Fig. 2a shows the spatial pattern of the long-term mean. Overall, SO₂ distribution is of great inhomogeneity in China, with two maximum centers: one is the North China Plain (NCP for short), and the other is Cheng-Yu (CY) district in Southwest China. In particular, the SO₂ amount in the NCP exceeds 1.2 DU. There are two essential causes responsible for high SO₂ loading in the two areas. On the one hand, the combined effect of rapid economic and industrial development as well as population growth leads to a high degree of anthropogenic SO₂ emission. Figure 2c and d show the emission strengths, defined as emitted SO₂ per unit area, in each province based on OYB and MEIC, respectively. Note that in the rest of this paper, the terms “emission” or “emission amount” always refer to “per unit area emission”. It is obvious that the two regions release above 8.0 t km⁻² SO₂ per year, which is 3 times greater than the average level of China. Although OYB exhibits a smaller magnitude of emissions than MEIC, the spatial patterns in terms of relative difference across space are generally consistent. On the other hand, as shown in Fig. 2b, either of the two regions is surrounded or partly surrounded by mountains, which makes it difficult for the pollutants to dissipate.

In contrast, over the sparsely populated western part of China, low SO₂ concentrations of less than 0.2 DU are observed, except over some provincial capitals. Since the western part of China is less affected by human activities, anthropogenic sources of SO₂ are much smaller than natural emissions including emissions from terrestrial ecosystems and oxidation of H₂S to SO₂ (Wang, 1999). Between latitude 30 and 40^∘ N, for example, the SO₂ amounts over the eastern regions (110–120^∘ E) are 6–12 times greater than western regions (80–110^∘ E). In addition, note that the low-level SO₂ columns in West China are subject to large uncertainties and the background correction is an important source of error. However, West China with weak SO₂ signals/background SO₂ is not the subject of the present work, since we mainly focus on highly polluted East China.

Besides the NCP and CY regions with the highest SO₂ loadings, this study is also interested in the Yangtze River Delta (YRD) and Pearl River Delta (PRD), the other two economic mega-urban zones in China. These four identified hotspots NCP, CY, YRD and PRD are outlined in Fig. 2a and will be specially examined in the following discussion.

3.2 Seasonal cycle

The annual total is decomposed into the seasonal cycle, as shown in Fig. 3. In East China, about 35 % of the annual total is from winter, while SO₂ in summer only accounts for 15 %; the remaining 50 % is almost equally divided into spring and autumn. Seasonal variations measured in the fractional contribution are similar within East China.

Figure 2 (a) Spatial distribution of 12-year (2005–2016) averaged SO₂ columns over China. (b) Topography of China in meters. (c, d) SO₂ emission (t km⁻²) among Chinese provinces based on OYB and MEIC.

To unveil the underlying mechanism, Fig. 4 illustrates the annual cycle of SO₂ VCDs in relation to sulfur emission, precipitable water and temperature at the four hotspots. Intensive heating during winter in North China raises sulfur release. However, emissions alone are not sufficient to explain the pronounced seasonality of SO₂. The remaining variation is associated with the seasonal change in the meteorological conditions. Temperature and humidity are cold and dry in winter due to the influence of winter monsoon, and jointly weaken the rate of oxidation and wet deposition. Thus, one expects that SO₂ molecules will have a longer lifetime and therefore will accumulate more easily. The opposite is true for summer, when chemical reaction is active and wet removal is effective. In summary, both emission and meteorological change explain the seasonality of the atmospheric SO₂ loadings.

Due to the climate transition from South China to North China, the annual range of SO₂ rises progressively from south to north. NCP has the greatest amplitude of up to 1 DU, while there is virtually no annual cycle in PRD. A larger amplitude for SO₂ cycles in NCP arises from the significantly reversed source–sink imbalance between summer and winter. In contrast, the climate in PRD is characterized by a smoother transition over the whole year and there is no heating season, which explains the insignificant seasonal variation of SO₂ in PRD. The other two regions, CY and YRD, have approximately the same amplitude of 0.6 DU, because they are on the same line of latitude.

Figure 3 Seasonal SO₂ columns over China: (a) winter, (b) spring, (c) summer and (d) autumn.

3.3 Long-term trends

Figure 5 depicts the spatial pattern of linear trends in annual and seasonal SO₂ from 2005 to 2016. Overall, apparent downward trends overwhelm most parts of East China, while West China has experienced little change. In particular, the most significant reduction occurred in the highly SO₂-polluted regions, with the decreasing rates amounting to 0.1 DU a⁻¹. This result suggests that the governments and communities in these economically developed regions have done their best to effectively control environmental pollution, including energy saving, emission cut, adjustment of energy consumption structure, shutdown of the most polluting factories, and upgrade of coal quality. Besides, enforcement of environmental protection laws is becoming more and more rigorous (van der A et al., 2017). Therefore, under collaborative efforts, the SO₂ levels in these highly developed regions with high background concentrations have been decreasing markedly in the last decade. Moreover, the pattern correlation between mean (Fig. 2a) and trends (Fig. 5a) of SO₂ reaches −0.77, implying that the downward rate over China can be summarized into a general rule – the higher the SO₂ loading, the more significant the decrease.

Figure 5b–e portray the long-term trends of SO₂ on a seasonal basis. On the one hand, every season has witnessed SO₂ reduction, with the strongest decrease occurring in winter and autumn. Consequently, it can be concluded that the SO₂ decreases in winter and autumn contribute most to the reduction of annual SO₂. On the other hand, the highly SO₂-polluted regions have experienced the most pronounced decrease across all seasons, which is consistent with the annual outcomes. It is worth pointing out that a belt of large positive values extends along 40^∘ N in winter (Fig. 5b). This feature is a known artefact related to the large solar zenith angles at high northern latitudes.

Last, we discuss the trends of the four hotspots involved. Figure 6 depicts the SO₂ columns from 2005 to 2016 as a function of year (y axis) and calendar month (x axis). The horizontal axis is the month of the year, the vertical axis is the year, and the color is the SO₂ VCD for that month and year. SO₂ VCDs exhibit a decreasing tendency during the last decade, regardless of the time of the year. Quantitatively, SO₂ in NCP, CY, YRD and PRD had undergone an overall downward trend with a rate of 0.062, 0.059, 0.046 and 0.055 DU per year, respectively.

4 Specific regimes of SO₂ variability and causes

4.1 Specific regimes of SO₂ variability

The above investigation presents SO₂ patterns and trends across China, but some elusive non-monotonic behaviors are not fully understood. In this section, we aim to detect the specific regimes of SO₂ variability and associated responsible mechanisms.

Figure 4 Annual cycle of SO₂ VCD (unit: DU, pink line), MEIC SO₂ emission (unit: t km⁻², yellow line), precipitable water (unit: kg m⁻², blue bar) and temperature at 925 hPa (unit: ^∘C, green line) for NCP (a), CY (b), YRD (c) and PRD (d).

Spatiotemporal regimes of SO₂ over China are mapped by using empirical orthogonal function (EOF) decomposition (Hannachi, 2004), which is a useful tool to reduce the data dimensionality to two dimensions. One dimension represents the spatial structure and the other the temporal dimension. Figure 7 illustrates the leading mode (top) and the corresponding principal component (PC, bottom) obtained from EOF, since only the first mode is statistically well separated. Compared to the first EOF mode explaining 36.8 % of the total variance, each of the other modes is characterized by less than 6 % contribution and is thus discarded. On the one hand, the variation of SO₂ is dominated by a spatially uniform feature with large loadings in NCP and CY, suggesting that SO₂ changes would be in the same phase but with varying amplitude across the entire region. On the other hand, the corresponding PC exhibits overall declines during 2005 to 2016. However, the result does not imply a simple continuous decrease. In fact, there appears to be a transient increase until a peak and thereafter two sharp drops occur in the periods 2007–2008 and 2014–2016, amid which SO₂ concentrations are in the process of slightly rebounding. In short, the SO₂ variability is characterized by four distinct temporal regimes.

Moreover, Fig. 8 demonstrates the time series for each province in East China, with the segment over 2009–2013 highlighted by red color. It reflects extensive common variation that goes through four stages – that is, a short-lived increasing period at the beginning, a steep drop period during 2007 to 2008, a rebound period of 2009 to 2013 and another drastic drop period during 2014 to 2016. Most importantly, it confirms that the SO₂ does not evolve in a monotonic way but shows a striking rebound during 2009 to 2013. This pattern is true throughout most of the region, with the only two exceptions of Guizhou and Guangdong provinces, which have experienced a continuous decrease since 2005.

Figure 5 Spatial pattern of SO₂ linear trends (2005–2016) in annual (a) and seasonal (b–e) values.

Figure 6 SO₂ amounts from 2005 to 2016 as a function of year (y axis) and calendar month (x axis) for NCP (a), CY (b), YRD (c) and PRD (d).

4.2 Causes

In this section, we diagnose the likely mechanisms behind the observed SO₂ variability. Generally, emissions and meteorological conditions are two main factors that essentially exert influence on atmospheric pollutant load. The impact of changes in emitted SO₂, as the main driving force, is first examined. To this end, the temporal classifications of SO₂ emission for each province based on OYB and MEIC are, respectively, depicted in Fig. 9a and b, in which the red upward pointing triangle implies non-monotonic decrease with a rebound in the middle, whereas persistent decrease is denoted by the green downward pointing triangle. In North China except Henan province, both the OYB and MEIC datasets show that the emission passed its secondary peak during 2009 to 2013. In South China, however, discrepancies between OYB and MEIC emerge in some provinces, namely Jiangxi, Hunan, Guangxi and Guizhou. Even so, we are still confident enough that the majority of South China has witnessed a successive drop in emitted SO₂. In addition, an auxiliary map is presented in Fig. 9c showing the slope of the linear regression of MEIC gridded emission over years 2008, 2010 and 2012. We can see that most of North China is subject to a positive rate of change, while the opposite holds true over most of South China, which confirms the above findings. Eventually, we come to the conclusion that despite spatial uniformity in temporal-pattern classification of SO₂ VCD (Fig. 8), the temporal structure of the emission demonstrates a strong south–north contrast (Fig. 9). Therefore, it is advantageous to treat North China and South China separately, as delineated by the dotted line in Fig. 9. Regional averaged quantities are estimated as a weighted average by assigning the district area as a weight. In addition, to evade possible contaminations, we have ruled out Henan and Jiangxi provinces in OYB and Henan, Hunan, Guizhou and Guangxi provinces in MEIC.

Figure 7 The first leading EOF mode (a) and the corresponding principal component (b).

Although we divide East China into north and south blocks, the inter-regional transport cannot be neglected. Therefore, we construct an Effective Emission Index (EEI) to account for impacts from both local and remote sources. Here, we directly adopt the results obtained by Zhang et al. (2015), who divided East China into three parts, North China, Southeast China and Southwest China, and quantified the percent contributions of within-region versus inter-regional transport on sulfate concentrations. The geographical partition in their work broadly coincides with ours, with the only difference that South China is further split into two parts. Given that the ratio of Southeast China to Southwest China is 1.4, we merge the percent contributions over these two portions via simple conversions. This produces, for North China, a within-region SO₂ emission contribution of 68 %, followed by 19 % from South China and 13 % from other regions; for South China, within-region emissions provide 66 %, while transport from North China and other regions amounts to 17 % and 17 %, respectively. With these statistics, the EEI is formulated as follows.

North China	EEI${}_{\mathrm{1}}=\mathrm{0.68}+\mathrm{0.19}+\mathrm{0.13}=\mathrm{1}$
	EEI${}_m=\mathrm{0.68}\cdot \frac{N_m}{N_{\mathrm{1}}}+\mathrm{0.19}\cdot \frac{S_m}{S_{\mathrm{1}}}+\mathrm{0.13}$

South China	EEI${}_{\mathrm{1}}=\mathrm{0.17}+\mathrm{0.66}+\mathrm{0.17}=\mathrm{1}$
	EEI${}_m=\mathrm{0.17}\cdot \frac{N_m}{N_{\mathrm{1}}}+\mathrm{0.66}\cdot \frac{S_m}{S_{\mathrm{1}}}+\mathrm{0.17}$

N and S denote the emission amount in North China and South China, respectively, and subscripts 1 and m the 1st and mth time node, respectively. The fundamental assumptions to derive the formula are that EEI is linearly dependent on N and S and the external contributions remain fixed (without interannual variation). For comparison purposes, we also define an Emission Index (EI) that involves a single effect from within-region emission, as written below.

North China	EI1=1
	EI${}_m=\frac{N_m}{N_{\mathrm{1}}}$

South China	EI1=1
	EI${}_m=\frac{S_m}{S_{\mathrm{1}}}$

The notations of symbols are identical to those in the EEI definition. In the case of a large scale, integrating the role of inter-regional transport does not alter the overall pattern, as proved in the following analyses.

Figure 10 presents time series and scatter plots of SO₂ VCD and emission with its variants EI and EEI, and Fig. 11 is designed to show the total emission generated by industries and households. These two figures are created based on the OYB inventory, while their counterparts obtained from the MEIC inventory are shown in Figs. S3 and S4. As shown in Fig. 10a, North China features a good correspondence between amount and either EI or EEI, with a linear correlation of 0.9. Time series of emission also indicate the existence of four distinct regimes that are likely to drive the regime of SO₂ VCD directly. This is confirmed by the scatter plot (Fig. 10b), in which the points are tightly clustered around the regression line. Based on variance analysis, emission accounts for an 81 % fraction of SO₂ VCD variation over North China. In parallel, the same procedure relying on the MEIC inventory yields nearly identical results, as shown in Fig. S3a and b. Furthermore, how large do industrial and household sectors, respectively, contribute to the total trends? Figures 11a and S4a indicate that the industrial emissions play a crucial role in SO₂ VCD variation throughout the whole period, while the influence induced by residential activity is secondary. A more in-depth comparison between OYB and MEIC shows some dissimilarity in household emission: OYB-based household emission acts to offset the industrial effect, while the opposite function is identified for the MEIC-based one. However, this does not seriously affect the major conclusion, due to the marginal impacts caused by households.

Figure 8 Temporal evolution of annual SO₂ (unit: DU) from 2005 to 2016 in each province of eastern China, with the segment over 2009–2013 highlighted by red color.

The close linear relationship observed in North China is not found in South China, since the two curves appear to become non-adherent in Fig. 10c and the points in the scatter plot Fig. 10d are widely spread around the regression line. Variance analysis suggests that only 45 % of SO₂ VCD variability is forced by emissions, suggesting that the SO₂ variations in South China cannot be explained by emission changes alone. This is mainly ascribed to the decoupled pathways of emission and SO₂ VCD during 2009 to 2013, as the emission continues to decline but SO₂ VCD witnesses a rebound. MEIC emissions also exhibit a general decreasing tendency in spite of a transient pause embedded, as shown in Fig. S3c. Moreover, Figs. 11b and S4b suggest that the cuts of industrial and household emissions collectively promote the continuous decrease in total emission in South China, which are different from that in North China. However, the emission decrease in the household sector is differently reported in the OYB and MEIC inventories, the former one showing a sudden shift, while the latter displays a gradual decrease. Anyway, one is assured that household emissions in South China have undergone a reduction, irrespective of the exact manner.

Why do decreasing emissions not cause a reduction of SO₂ VCD in South China during 2009 to 2013? To answer this question, the atmospheric conditions during 2009 to 2013 are compared with those during the rest of the years, as depicted in Fig. 12. The period 2009 to 2013 is characterized by prolonged dry conditions in South China with the precipitable water and precipitation being lower than usual (Fig. 12a), which weakens wet adsorption and scavenging. At the same time, this period is also associated with relatively weaker wind speed (Fig. 12b) and increased static stability (Fig. 12c and d), reducing the ability of the atmosphere to diffuse, leading to the accumulation of SO₂ loads. In brief, unfavorable meteorological conditions produce the observed rebound of SO₂ during 2009 to 2013, despite the continued decrease in emission.

Figure 9 (a, b) Temporal structure classification of SO₂ emission based on OYB and MEIC. Red upward pointing triangles imply non-monotonic decrease with a rebound in the middle, while monotonic decrease is denoted by green downward pointing triangles. (c) Slope of the linear regression of MEIC gridded emission over years 2008, 2010 and 2012. The black or red dotted line delimits North China and South China.

Figure 10 Time series plots of SO₂ VCD and EI/EEI (a, c), and scatter plots with the regression line of SO₂ VCD and emission (b, d) for North China (a, b) and South China (c, d). Each marker in (b) and (d) corresponds to one year and one province.

5 Conclusions

In this study, the spatiotemporal variability of SO₂ columns over China and the associated driving mechanisms are examined over the past decade. Based on a state-of-the-art SO₂ retrieval dataset recently derived from the OMI, we elaborate on the characteristics of specific SO₂ regimes over China and underlying causes.

Climatological SO₂ in China has an uneven spatial distribution in space and time. East China is far more exposed to SO₂ pollution than West China, with two maxima centered in NCP and CY. From analysis of the annual cycles we conclude that 35 % of the annual totals are from winter, while SO₂ in summer only accounts for 15 %. In addition, the annual amplitude of SO₂ rises progressively from south to north.

Figure 11 Annual SO₂ emission (t km⁻²) generated by industries (upward blue bars) and households (downward red bars) in North China (a) and South China (b). Note that the y axis in a positive direction does not start at zero.

From 2005 through 2016, most of East China presents a clear decreasing tendency for SO₂, while West China experiences little change. Spatially, the decreasing rate is generally enhanced for high SO₂ loads. When computed seasonally, SO₂ reductions in winter and autumn contribute most to the reduction of annual SO₂.

Figure 12 Comparison of atmospheric conditions between the period of 2009–2013 and the other years: (a) composite difference in precipitable water (unit: kg m⁻²), (b) composite difference in wind velocity at 925 hPa (unit: m s⁻¹), (c) composite difference in static stability at 925 hPa (unit: K hPa⁻¹), and (d) averaged vertical profile of static stability over the 23–31^∘ N, 105–122^∘ N rectangle for the two episodes.

Four stages of variation are identified by EOF analysis. The first regime (2005–2006) features a transient increasing trend, the second (2007–2008) and the last (2014–2016) regimes show sharp drops, and the third regime (2009–2013) manifests itself in 5-year moderate rebounding. Although temporal regimes of SO₂ are coherent throughout the country, different driving forces are tied to North China and South China. In North China, the atmospheric SO₂ and emission vary essentially in the same way. Therefore, the atmospheric SO₂ variability is primarily associated with the emission variability, which accounts for 81 % of the total variance. Further, the emission generated by the industrial sector is largely responsible for the atmospheric SO₂ variability. The household emissions appear to remain uncertain, due to the dissimilarity between the OYB and MEIC inventories.

SO₂ emissions in South China exhibit a continuous decreasing tendency, due to the coordinated cuts of industrial and household emissions. As a result, the role of emissions only contributes 45 % of the SO₂ variation, primarily owing to the decoupled pathways of emission and atmospheric content during 2009 to 2013 when the emission continues to decline but atmospheric content witnesses a rebound. It is found that such rebound occurs in response to the joint effect of deficient precipitation, weaker wind speed and increased static stability during 2009 to 2013.

6 Future directions

As shown by this study, the spatial and temporal changes in SO₂ regimes over China in the recent decade become clear. However, there is much left to be learned about the responsible driving mechanisms. First, a major obstacle to cause-and-effect relation surveys stems from uncertainties in the current emission inventories. In this study, many facets inferred by OYB and MEIC are convergent, because we look at a large spatial scale and long-term general tendency that help filter out or attenuate some uncertainties. However, if the aim is to focus on smaller spatial or temporal scales or on specific sectors, there is still great uncertainty. To overcome these barriers, emission inventories should be further improved and more observational products should be used for comparison. Second, this work investigates the impact of emission, inter-regional transport and meteorology using purely statistical techniques, but finer-scale investigations require numerical simulations using coupled chemical-transport models. Third, the analysis presented in Sect. 4 is constrained to provincial or multi-provincial levels, due to the limitation that only continuous emission data on provinces are gathered at hand. In reality, however, either emission or atmospheric loadings can be quite inhomogeneous within the same region. Therefore, future studies should use both gridded SO₂ VCDs and gridded SO₂ emission inventories.