Investigating the response of China's surface ozone concentration to the future changes of multiple factors

Climate change and associated human response are supposed to greatly alter surface ozone (O₃), an air pollutant generated through photochemical reactions involving both anthropogenic and biogenic precursors. However, a comprehensive evaluation of China's O₃ response to these multiple changes has been lacking. We present a modeling framework under Shared Socioeconomic Pathways (SSP2-4.5), incorporating future changes in local and foreign anthropogenic emissions, meteorological conditions, and biogenic volatile organic compound (BVOC) emissions. From the 2020s to 2060s, daily maximum 8 h average (MDA8) O₃ concentration is simulated to decline by 7.7 ppb in the warm season (April–September) and 1.1 ppb in the non-warm season (October–March) over the country, with a substantial reduction in exceedances of national O₃ standards. Notably, O₃ decreases are more pronounced in developed regions such as Beijing–Tianjin–Hebei (BTH), the Yangtze River Delta (YRD), and the Pearl River Delta (PRD) during the warm season, with reductions of 9.7, 14.8, and 12.5 ppb, respectively. Conversely, in the non-warm season, the MDA8 O₃ in BTH and YRD will increase by 5.5 and 3.3 ppb, partly attributed to reduced NO_x emissions and thereby a weakened titration effect. O₃ pollution will thus expand into the non-warm season in the future. Sensitivity analyses reveal that local emission change will predominantly influence future O₃ distribution and magnitude, with contributions from other factors within ±25 %. Furthermore, the joint impact of multiple factors on O₃ reduction will be larger than the sum of individual factors, due to changes in the O₃ formation regime. This study highlights the necessity of region-specific emission control strategies to mitigate potential O₃ increases during the non-warm season and under the climate penalty.

1 Introduction

Surface ozone (O₃) is a secondary air pollutant generated by photochemical reactions in the presence of two main kinds of precursors, NO_x (NO_x= NO + NO₂) and volatile organic compounds (VOCs). It has been reported to be a non-negligible threat to both human health and crop yield and also a short-lived climate forcer with warming effect (Finlayson-Pitts and Pitts, 1997; Jerrett et al., 2009; Avnery et al., 2011; Lelieveld et al., 2015; von Schneidemesser et al., 2015; Tai and Val Martin, 2017; Yin et al., 2020; Feng et al., 2022; Niu et al., 2022). Given the abundant emissions of anthropogenic NO_x and VOCs, China suffered from extremely high and continuously increasing O₃ pollution from 2013 to 2019, with the peak season daily maximum 8 h average (MDA8) O₃ concentration over 95 µg m⁻³. The rising trend has been reversed since 2020, along with the national annual NO_x and non-methane VOC (NMVOC) emissions reduced by 28.3 % and 3.8 %, respectively, during 2013–2020 (Zheng et al., 2018; Xiao et al., 2022; Liu et al., 2023a; Wang et al., 2023). However, current O₃ concentration over China is still much higher than the global air quality guidelines (60 µg m⁻³ for the averaged peak season MDA8 O₃, WHO, 2021). This presents a great challenge for the country to meet the criteria for public heath welfare in the future (Feng et al., 2023; Jiang et al., 2023).

In addition to anthropogenic drivers, studies also addressed the roles of meteorological factors, biogenic VOC (BVOC) emissions, and transboundary transport of pollutants in O₃ enhancement in China (Monks et al., 2015; Lu et al., 2019; Cao et al., 2022; Wang et al., 2022; Weng et al., 2022; Jiang et al., 2023). Meteorological factors, including temperature, humidity, and wind, influence the chemical reactions associated with O₃ production and elimination and the transportation of O₃ precursors (Gong and Liao, 2019). The changed meteorology was estimated to enhance the summer MDA8 O₃ concentration by 1.4 ppb yr⁻¹ during 2013–2019 in the North China Plain (NCP), nearly half of the overall O₃ growth of 3.3 ppb yr⁻¹ (Li et al., 2020a). BVOCs refer to VOCs emitted from terrestrial ecosystems and possess high reactivity in atmospheric chemical processes, mainly including isoprene, monoterpenes, and sesquiterpene (Wu et al., 2020b). Cao et al. (2022) reported that BVOC emissions in summer 2018 were enhanced by 8.6 ppb MDA8 O₃ averaged over China, with the highest contribution over 30 ppb in Southern China. Moreover, O₃ and its precursors can be transported over long distances, and transboundary foreign anthropogenic emissions have been estimated to contribute 2–11 ppb to near-surface O₃ in China (Ni et al., 2018; Han et al., 2019).

In the context of future global change, substantial but uncertain changes will occur in economy, climate, and land cover. According to the Sixth Assessment Report of Intergovernmental Panel on Climate Change (IPCC AR6 report), the global surface temperature will increase 0.2–3.7 °C till 2100 under different scenarios compared to 2015 (IPCC, 2021, 2022). As a result, unfavorable meteorological extremes, such as high temperature extremes and ecological drought events, will be more frequent and intense (Hong et al., 2019; Porter and Heald, 2019; IPCC, 2021), leading to the deterioration of air quality, called the climate penalty. To conquer the climate change and resulting air quality deterioration, a series of measures will be implemented, such as accelerating the transition to clean energy, upgrading industrial production technologies, and strengthening pollution control measures. Attributable to these changes, annual mean surface O₃ in east Asia was projected to change by −13.9–6.1 ppb till 2100 compared to 2015 (IPCC, 2021). However, limited by the coarse resolution of earth system models and the lack of consideration of the regional measures for reducing air pollution and carbon emissions, global estimation is insufficient for understanding how O₃ pollution in China will respond to the complex future change.

There have been some studies on how the above-mentioned changes will affect future O₃ levels in China. Hong et al. (2019) reported that the 1 h maximum O₃ in April to September will be enhanced by 2–8 ppb within large areas of China under the RCP4.5 (Representative Concentration Pathway 4.5, van Vuuren et al., 2011) scenario from the 2010s to 2050s. Under high-forcing scenarios, Li et al. (2023) projected the climate-driven O₃ concentration in the 2100s and found that O₃ concentration in southeastern China would increase 5 %–20 % compared to the 2020s using a machine learning method. A warming climate should enhance the O₃ level, given the increasing frequency of atmospheric stagnation and heat waves (Hong et al., 2019; Wang et al., 2022; Gao et al., 2023; Li et al., 2023). The effect of anthropogenic emission change on China's O₃ level has been estimated by studies under different scenarios. Zhu and Liao (2016) applied global emission estimates under RCPs and found that the maximum growth of annual mean O₃ would be 6–12 ppb during 2000–2050 under different scenarios. Using the Dynamic Projection model for Emissions in China (DPEC) that more effectively includes local information on energy transition and emission controls (Cheng et al., 2021b), Xu et al. (2022) reported that the joint impact of climate change and emission reduction would reduce the annual MDA8 O₃ concentration to 63.0 µg m⁻³ under ambitious scenario of carbon neutrality. Biogenic emission change is another factor influencing future O₃ (Chen et al., 2009; Andersson and Engardt, 2010; Harper and Unger, 2018; Wang et al., 2020). Liu et al. (2019) predicted a 24 % growth of BVOC emissions driven by climate change under RCP8.5 from 2015 to the 2050s, resulting in a variation of daily 1 h maximum O₃ concentration ranging from −10.0 to 19.7 ppb across different regions in China.

Limitations exist in current studies, which prevent comprehensive assessment and understanding of the joint impacts of future changes of multiple factors on China's O₃ pollution. Firstly, the above estimations mainly focused on the influence of future changes on summertime or annual average O₃ concentration. As China's O₃ pollution has been reported to spread into spring and fall, it is of great importance to separate the impacts on warm-season (April to September, the 6 months with heaviest O₃ pollution for most of China, Liu et al., 2023a) and non-warm-season O₃ (October to March), considering the diverse air pollution sources and O₃ formation sensitivity to precursors for different seasons (Li et al., 2021; Wang et al., 2023). Recent studies have suggested diverse effects of future emission change on O₃ evolution for difference seasons in China (Hou et al., 2023; Liu et al., 2023b). In addition, the rising frequency of extreme weathers and declining anthropogenic emissions will further influence the possibility of extreme O₃ events, which has been scarcely discussed. Secondly, to restrain global warming, China has made a national commitment to achieve “carbon neutrality” by 2060 (Shi et al., 2021) and accordingly launched a series of energy and climate action plans to reduce greenhouse gas emissions. These actions will also cause substantial reductions in air pollutant emissions but have not been fully included in existing predictions of global emissions (Tong et al., 2020; Cheng et al., 2021b). Large bias will then be caused in the simulation of anthropogenic-induced future changes of air quality, with a less realistic estimate of the local emission path (Cheng et al., 2021a). Due to a decline in emissions that is probably faster in China but slower in surrounding countries in the future, the contributions of transboundary emissions to China's O₃ could be greatly changed and have not yet been fully considered (Hou et al., 2023). Thirdly, BVOC emissions will be affected not only by meteorological factors but also by land use and land cover change (Penuelas and Staudt, 2010; Szogs et al., 2017; Wang et al., 2021a). Future land management will change due to socio-economic development and necessary actions as part of the climate change response, and the changed shares of forest, cropland, and grassland will alter the magnitude and distribution of BVOC emissions and thereby affect O₃ concentration (Hurtt et al., 2020; Liao et al., 2020; Liu et al., 2022). Finally, the existing evaluations were conducted separately for individual influencing factors, with diverse methods and data. The interactions between different factors were seldom included in existing analyses, and the relative contributions of multiple factors were difficult to be evaluated or compared. Relevant studies have been conducted in developed countries (Gonzalez-Abraham et al., 2015) and are still lacking in China.

In this study, we evaluate the complex influence of future changes of multiple factors on surface O₃ concentration in China within a uniform framework. The evaluation is conducted from the perspectives of seasonal, regional, and extreme events of O₃ pollution. Four factors are included in the analyses, i.e., meteorological conditions, local anthropogenic emissions, BVOC emissions, and anthropogenic emissions from surrounding foreign countries. The analyses are conducted based on a series of sensitivity experiments in numerical modeling of future air quality, and up-to-date input data from multiple sources are utilized in the model (see details in next section). We provide a comprehensive perspective on the spatiotemporal change in China's O₃ pollution till the 2060s, under a moderate SSP (Shared Socioeconomic Pathway) scenario, SSP2 (Riahi et al., 2017), and a mid-range mitigation scenario, RCP4.5, a scenario with medium socio-economic development, with radiative forcing of 4.5 W m⁻² nominally by 2100 (Meinshausen et al., 2020). The outcomes highlight the regional and seasonal heterogeneity of O₃ pollution risks driven by complex future changes of multiple factors and support strategy design of O₃ pollution alleviation with specific principles, targets, and action pathways.

Figure 1 The modeling domain and geographical definitions (denoted by colors) of this study. Boundaries of the three regions, including BTH (Beijing–Tianjin–Hebei), the YRD (Yangtze River Delta), and the PRD (Pearl River Delta), are marked by dark-grey lines.

2 Data and methods

2.1 Main framework and research domain

The simulation framework incorporates the Weather Research and Forecasting model (WRF, version 3.7.1) to generate hourly meteorological fields, the Model of Emissions of Gases and Aerosols from Nature (MEGAN, version 2.1) to calculate gridded BVOC emissions, and the Community Multiscale Air Quality model (CMAQ, version 5.2) to simulate O₃ concentration. BVOC emission calculations and air quality simulations are driven by meteorological fields of 2018–2022 (the 2020s, representing the current situation) and 2058–2062 (the 2060s, representing the future situation). All simulation results are averaged over a period of 5 years to mitigate the influence of interannual variability of meteorology. The modeling domain, as for WRF, MEGAN, and CMAQ, covers east Asia and portions of other Asian regions, including north Asia, central Asia, south Asia, and southeast Asia (Fig. 1). It applies the Lambert conformal conic projection, centered at 34° N, 110° E, and the horizontal resolution is 27 km × 27 km, with 303×203 grids. The target area, the Chinese mainland, includes 31 provincial-level administrative regions (excluding Hong Kong SAR, Macao SAR, and Taiwan). Eight geographical regions are defined, and locations of the three regions with dense population and relatively heavy air pollution are also shown in Fig. 1, namely BTH (Beijing–Tianjin–Hebei), the YRD (Yangtze River Delta), and the PRD (Pearl River Delta).

2.2 Data sources and processing methods

We use the bias-corrected RCP4.5 output of the National Center for Atmospheric Research's Community Earth System Model (NCAR CESM) as initial and boundary conditions for WRF (Monaghan et al., 2014). A 10-year dynamic downscaling simulation for 2018–2022 and 2058–2062 is conducted. Note we do not utilize the real-time reanalysis data to drive the simulation of the 2020s, in order to minimize the systematic error between the simulation driven by real meteorological conditions (for current simulations) and climate projection (for future simulations).

The BVOC emissions are basically determined by meteorology and vegetation. The meteorological conditions are supplied by WRF. The vegetation data, including leaf area index (LAI), plant functional types (PFTs), and emission factors (EFs) of each PFT, are determined for 2020 and 2060. Gridded LAI data for 2020 are obtained from Global Land Surface Satellite product (Liang et al., 2021), and those for 2060 under the SSP2-4.5 scenario are downscaled from the daily CESM2 output of Coupled Model Intercomparison Project Phase 6 (CMIP6). PFT data for 2020 are derived from the MCD12C1 product of the Moderate-Resolution Imaging Spectroradiometer (MODIS) dataset and mapped to the 16 types required for MEGAN following Liao et al. (2020). The PFT data for 2060 in China are obtained from Liao et al. (2020) under the SSP2-4.5 scenario, with other regions maintaining 2020 levels. EFs for each PFT are taken from Guenther et al. (2012).

Anthropogenic emissions for the Chinese mainland are obtained from the Multi-resolution Emission Inventory for China (MEIC; http://meicmodel.org.cn/?page_id=560, last access: 24 March 2024) for 2020 and DPEC version 1.1 under SSP2-4.5, incorporating the best available end-of-pipe pollution control technologies, for 2060. The annual total emissions by country/region outside the Chinese mainland are obtained from the CMIP6 dataset under the SSP2-4.5 scenario (O'Neill et al., 2016; Gidden et al., 2019), and these emissions are downscaled into gridded monthly data for CMAQ simulation, based on the spatial and temporal distributions of emissions in MIX Asian emission inventory (Li et al., 2017). The speciation profiles of NMVOCs are taken from MIX as well. Figure S1 in the Supplement shows the emissions of two main precursors of O₃ by year and region. The NO_x and NMVOC emissions for the Chinese mainland were estimated to decline 58 % and 51 % from 2020 to 2060, respectively, much faster than those of surrounding areas within the modeling domain (8 % and 14 %, respectively). In particular, the NO_x emissions would decline 57 %–62 % for the three developed regions, BTH, the YRD, and the PRD, while the reductions of anthropogenic NMVOCs would vary a lot among regions (36 %, 49 %, and 60 % for BTH, the YRD, and the PRD, respectively).

Carbon Bond 2005 (CB05; Yarwood et al., 2005) is adopted as the gas-phase chemical mechanism and the sixth-generation CMAQ aerosol module AERO6 (Appel et al., 2013) as the aerosol chemistry mechanism. The initial and boundary conditions are set by default clean air conditions in CMAQ, and the first 10 d for each year is determined as the spin-up period to minimize the effects of initial and boundary conditions.

Table 1 List of simulation cases to investigate the impact of future change upon surface O₃ in China, with sensitivity experiments from the perspectives of four main influencing factors.

2.3 Simulation cases

Six cases of CMAQ simulations are conducted to investigate the impacts of future change of the four factors on O₃ concentration in China (Table 1). Cases 1 and 2 represent the current (2020s) and future (2060s) baseline, respectively, and the difference between them indicates the joint effect of the future changes of multiple factors. Cases 3–6 each apply the prediction for 2060s for one specific factor but keeps the remaining factors at current condition (2020s). Thus, the difference between each of those four cases and Case 1 indicates the impact of individual factor, including meteorological conditions (Case 3), domestic anthropogenic emissions (Case 4), BVOC emissions (Case 5), and anthropogenic emissions of surrounding countries (countries other than the Chinese mainland within the modeling domain, Case 6). Each case contains a 5-year (2018–2022 or 2058–2062) WRF–MEGAN–CMAQ simulation driven by the varying meteorological conditions for individual years, and the 5-year average of simulated O₃ concentrations is adopted for further analyses.

2.4 Model performance

To evaluate the model performance, we conduct a comparative analysis between simulations and observations for meteorological factors and O₃ concentrations, as well as an intercomparison for BVOC estimates between different studies.

We first examine the capability of downscaled CESM climate projections in capturing the meteorological conditions of the 2020s. We applied the meteorological data from the National Climate Data Center (NCDC, archived at https://quotsoft.net/air, last access: 13 December 2021) in 2020, and the statistical metrics are presented in Table S1 in the Supplement. The modeled temperature at 2 m (T2) is in good spatiotemporal agreement with the observations, with the correlation coefficient (R) of 0.96 and index of agreement (IOA) of 0.98. The relative humidity (RH) is also well predicted with R and IOA at 0.78 and 0.88, respectively. The model shows an overestimation on the wind speed by 1.41 m s⁻¹, which is also reported by Hu et al. (2022). The correlation coefficients of wind speed and direction are higher than 0.5. Overall, the modeled meteorological fields have basically captured the conditions in China and are appropriate for subsequent MEGAN and CMAQ simulations.

For BVOC emissions, we compare our estimates for the 2020s with previous studies, as summarized in Table S2. The total BVOC, isoprene, and terpene emissions in this study are estimated at 33.55, 21.08, and 3.30 Tg yr⁻¹, respectively, and are comparable to other studies. In particular, our estimate is larger than others, except for Li et al. (2020b) for isoprene, while it is smaller than others, except for Wu et al. (2020b) for terpenes. The differences between studies might result from the diverse strategies of mapping PFTs from the original satellite products and the difference between downscaled climate conditions and the real meteorological fields.

We apply the observed MDA8 O₃ concentration data from the national network of the China Ministry of Ecology and Environment (archived at https://quotsoft.net/air, last access: 2 April 2021) to evaluate CMAQ performance. As shown in Fig. S2, the simulation could capture the spatiotemporal distribution of surface MDA8 O₃ concentration for the whole country and specific O₃ pollution hot spots, e.g., BTH and eastern Sichuan Province with their surrounding areas. The statistical metrics of the comparisons between the simulated and observed monthly average MDA8 O₃ concentration of 2020s are summarized in Table S1. The normalized mean biases (NMBs) are calculated at 14.12 % and 10.90 % for the warm and non-warm season, and R values at 0.71 and 0.32, respectively. Even with a slight overestimation, the reliability of our simulation is comparable to most previous studies in China, with a better performance in the warm season (Hu et al., 2016; Lu et al., 2019; Gao et al., 2020; Yang and Zhao, 2023).

We evaluate the interannual variability within each of the 5-year simulations, based on the coefficient of variation (CV), the ratio of standard deviation to the mean of the simulated O₃ concentration. As shown in Table S3, the CVs are generally below 5 % in most cases, indicating relatively small interannual variability in O₃ concentration simulation. The results thus justify the representativeness of the 5-year averages for present and future scenarios.

Figure 2 Projected changes of the ozone-related meteorological factors, including daily maximum temperature at 2 m (T_max; a, d), relative humidity (RH; b, e), and wind speed (WS; c, f), from the 2020s to 2060s. Panels (a)–(c) represent the warm season, and panels (d)–(f) represent the non-warm season.

3 Results and discussions

3.1 Future change of meteorology and BVOC emissions

The downscaled changes in the meteorological factors from the 2020s to 2060s (SSP2-4.5 scenario) are shown in Fig. 2, including temperature, RH, and wind speed (WS). The changes are analyzed separately for April–September (warm season) and October–March (non-warm season). For the warm season, daily maximum temperature at 2 m (T_max) will increase across China with an average change of 1.0 °C, and the minimum and maximum changes are found in Tibetan Plateau at 0.1 °C and in Heilongjiang Province at 2.1 °C, respectively. The RH will decrease slightly by −0.6 % for the whole country, with the changes for most areas within the range between −3 % and 0 % except for some areas of Northwestern China, Southwestern China, and the Tibetan Plateau (see the region definitions in Fig. 1). The growing T_max and declining RH will enhance the photochemical production of O₃ and BVOC emissions. For the non-warm season, the national average growth of T_max will be smaller at 0.2 °C, and some areas in Northeastern, Northern, and Eastern China will even experience a decline ranging from −1.8 to 0 °C. The RH will change diversely across the country, ranging from −6.0 % to 6.3 %. Very limited change in WS will occur, ranging from −0.1 to 0.2 m s⁻¹ in most areas of the country. Generally, the decreasing wind speed in future east Asia could be attributed to weakened atmospheric circulation (Coumou et al., 2018; Deng et al., 2021). The increasing wind speed in the non-warm season might result from the temperature and pressure gradients between the land and adjacent oceans (Yao et al., 2019; Wu et al., 2020a). The spatial distribution of downscaled future meteorological field changes is generally in agreement with those predicted by Hong et al. (2019) and Hu et al. (2022). Some discrepancies in temperature and wind speed change of the non-warm season between studies result from the different choices of base year and parameterization schemes of WRF.

Table 2 The BVOC estimation over China and emission intensity in BTH, the YRD, and the PRD of the 2020s and 2060s, as well as the corresponding growth rates over this period.

Table 2 shows China's BVOC emissions of the 2020s and 2060s (SSP2-4.5 scenario) estimated with MEGAN, as well as the BVOC emission intensity (emissions per unit area) for the three developed regions. The emissions will increase from 33.6 Tg yr⁻¹ for the 2020s to 43.4 Tg yr⁻¹ for the 2060s. The growth rates in BTH, YRD, and PRD are predicted to be 21.4 %, 23.9 %, and 23.0 %, respectively, smaller than for the whole country (29.2 %). The spatial distributions of BVOC emissions for the 2020s and the changes from the 2020s to 2060s are shown in Fig. S3. Areas all over China will experience growth in BVOC emissions, and this will be more prominent in areas with high vegetation coverage (e.g., Southern and Southwestern China) rather than urban areas. The growth of BVOC emissions will enhance the contribution of natural sources to O₃ formation, especially along with declining anthropogenic emissions in the future (Penuelas and Llusia, 2003; Riahi et al., 2017; Gao et al., 2022).

Figure 3 Simulation and projection of seasonal average MDA8 O₃ in the 2020s (Case 1, a, b) and 2060s (Case 2, d, e) and the changes over this period (Case 2 − Case 1, c, f). Panels (a)–(c) represent the warm season, and panels (d)–(f) represent the non-warm season. Regional mean concentrations across China (CHN), BTH, the YRD, and the PRD are shown in the inset on the left.

3.2 Response of surface O₃ concentration to combined future changes

Figure 3 illustrates the spatial distributions of MDA8 O₃ concentrations for the warm and non-warm seasons of the 2020s and 2060s (SSP2-4.5 scenario), as well as the differences between the two periods. Briefly, future changes of the four factors under SSP2-4.5 are estimated to jointly reduce MDA8 O₃ by 7.7 and 1.1 ppb in the warm and non-warm season, respectively, while the O₃ responses to future changes will differ by region.

In the warm season of the 2020s (Fig. 3a), the nationwide average MDA8 O₃ concentration is simulated at 57.3 ppb, and the concentrations in BTH, the YRD, and the PRD are 73.7, 68.7, and 52.3 ppb, respectively. Hot spots of O₃ pollution, with average MDA8 O₃ over 75 ppb, are mainly located in Northern China and Sichuan Province. The pattern is predicted to persist into the 2060s (Fig. 3b), with a decline in both the severity and size of highly polluted regions. The nationwide MDA8 O₃ concentration will decline 13.4 % to 49.6 ppb, and in most areas of China it will be within the range of 37.5–67.5 ppb. The highest concentration will be lower than 75 ppb for the two hot spots of Northern China and Sichuan. BTH will remain the most O₃-polluted area in the warm season, with the O₃ concentration at 63.9 ppb (13.3 % smaller than the 2020s), while that of YRD and PRD will decrease to 53.9 (21.5 %) and 39.8 ppb (23.9 %), respectively. O₃ concentration in the developed regions will decline faster than or roughly the same as that for the whole country. The reductions in MDA8 O₃ from 2020s to 2060s will be 10–20 ppb for Northern, Eastern, Central, and Southern China and 0–10 ppb for Northeastern and Northwestern China as well as the Tibetan Plateau (Fig. 3c). Notably, some areas in Sichuan are expected to experience a substantial decline in MDA8 O₃ over 20 ppb.

O₃ concentration of the non-warm season is simulated to be much lower than that of the warm season. The 2020s average MDA8 O₃ is 48.4 ppb, ranging from 30.0 to 67.5 ppb in most areas of China (Fig. 3d). Different from the warm season in which the highest concentration is found for Northern China and Sichuan, Southern China, and Southwestern China suffer the highest O₃ level for the non-warm season. A general west–east and south–north gradient is found for MDA8 O₃, with the lowest concentration found in Northern and Northeastern China. The concentrations in BTH and YRD are simulated at 33.8 and 45.1 ppb, respectively, much lower than the concentration in the PRD (58.9 ppb). Relatively high temperature during even the non-warm season is expected to expand the O₃ pollution period in Southern China. Resulting from complex change of multiple factors, the national average MDA8 O₃ concentration in the non-warm season of 2060s will decrease slightly to 47.3 ppb under SSP2-4.5, and in most regions it will be within the range of 37.5–52.5 ppb, except for some areas in Northeastern China and the Tibetan Plateau (Fig. 3e). The MDA8 O₃ concentrations of the three developed regions will become closer at 39.3, 48.4, and 51.6 ppb for BTH, the YRD, and the PRD, respectively. As illustrated in Fig. 3f, MDA8 O₃ is predicted to increase in BTH and YRD and the surrounding areas, with the growth mostly ranging 0–15 ppb. In other areas (especially in Southern China), the concentration will decrease in the non-warm season by −15 to −5 ppb. As a result of the increased O₃ in the less polluted regions of Eastern and Northern China and decreased O₃ in the more polluted regions of Southwestern and Southern China, the 2060s regional disparity in the non-warm season O₃ pollution will get smaller compared to the 2020s (Fig. 3d and e).

Figure 4 Simulation and projection of monthly average MDA8 O₃ in the 2020s and 2060s across CHN (a), BTH (b), YRD (c), and PRD (d).

To further explore the temporal pattern of O₃ level in the future, we compare the monthly average MDA8 O₃ in the 2020s and 2060s under SSP2-4.5 for the whole country and three developed regions (Figs. 4 and S4). For the whole country (Fig. 4a), the changes of monthly average MDA8 O₃ from 2020s to 2060s are estimated to range from −3.2 to −10.7 ppb in the warm season but to be less prominent in the non-warm season (from −2.7 to 0.9 ppb). Along with the greater reduction in summertime (June, July, and August), in particular, the periods with the highest O₃ concentration will expand into spring (March) and fall (October), as presented in Fig. S4. For the three regions, a greater decline in O₃ concentration is found in the warm season, while a smaller or even a growth is found in the non-warm season. For BTH (Fig. 4b), the monthly MDA8 O₃ concentrations range between 24.7 and 88.4 ppb in the 2020s, with a clear difference between the warm and non-warm season. This pattern will remain in the 2060s, with a smaller difference between months (30.6–70.2 ppb). The temporal change pattern of the YRD is similar to that in BTH, with a decline in the warm season and growth in the non-warm season (Fig. 4c). The shift in O₃ pollution from the warm towards the non-warm season is more prominent in the PRD, the only region where O₃ concentration of all the months in the 2060s is predicted to decline (Fig. 4d). Different from BTH and the YRD, as mentioned above, higher O₃ concentrations during spring and fall and lower in summer (due to the abundant summertime precipitation and high humidity) are found for PRD in the 2020s (Gao et al., 2020; Han et al., 2020). With great O₃ decline in the warm season, the periods experiencing peak O₃ pollution are predicted in the non-warm season of the 2060s, predominantly between October and March (Fig. S4).

Figure 5 Projected changes of MDA8 O₃ from the 2020s to 2060s attributed to anthropogenic emissions from local sources (Case 3 − Case 1, a, e), meteorological conditions (Case 4 − Case 1, b, f), BVOC emissions (Case 5 − Case 1, c, g), and anthropogenic emissions from surrounding countries (Case 6 − Case 1, d, h). Panels (a)–(d) represent the warm season, and panels (e)–(h) represent the non-warm season. Regional mean changes across CHN, BTH, the YRD, and the PRD are shown in the inset on the left.

3.3 Identifying surface O₃ response to individual factors

3.3.1 Local anthropogenic emission change

Figure 5 shows the influences of changes of individual factors (local anthropogenic emissions, meteorological conditions, BVOC emissions, and anthropogenic emissions from surrounding countries) on the warm and non-warm season O₃ concentrations. Out of the four, the change of local anthropogenic emissions is predicted to be the most influential factor, resulting in a national average decline of 7.2 and 0.8 ppb for the warm and non-warm season, respectively (Fig. 5a and e). In the warm season, the emission reduction will play a positive role in reducing O₃ pollution in most areas of China, and the decrease will exceed 10 ppb across Northern, Eastern, Central, Southern, and part of Southwestern China. In the non-warm season, emission reduction will have contrasting effects on MDA8 O₃ levels in the north and south of China, enhancing MDA8 O₃ by 0–15 ppb for the former while restraining it by 0–10 ppb for the latter. In particular, the emission reduction is predicted to elevate the O₃ concentration by 5.9 and 4.0 ppb for BTH and the YRD, respectively.

Figure S5 shows the relative emission reductions from the 2020s to 2060s by region. Under the SSP2-4.5 scenario, the reductions in NO_x and VOC emissions will range from 35.6 % to 63.6 % for different regions, and VOC emission reduction will be less than that of NO_x except for the PRD. As the NO_x-limited regime for O₃ formation (i.e., O₃ is more sensitive to NO_x emission change) occurs more frequently in the warm season, while the VOC-limited regime occurs more in the non-warm season, the larger decline in NO_x emissions than VOCs should be more effective in restraining the warm season O₃ pollution, but it has less benefit or even a negative effect in the non-warm season (Sillman and He, 2002). Wintertime in the NCP and YRD has been reported under the VOC-limited regime, and the excessive NO_x emissions play an important role in removing O₃ via the titration effect (Jin and Holloway, 2015; Li et al., 2021; Wang et al., 2021b). This may explain the MDA8 O₃ increase during the non-warm season with an insufficient reduction in VOCs (35.6 % and 49.5 %) but a sharp reduction in NO_x of 53.4 % and 60.3 % for the NCP and YRD, respectively. Similarly, Hou et al. (2023) and Liu et al. (2023b) also predicted a growth of O₃ concentration in the non-warm season over BTH and the YRD under a net-zero carbon emission scenario, resulting from a weakened titration effect. Figure S6 shows the monthly variation of O₃ and odd oxygen (O_x, O_x= O₃+ NO₂, representing the real photochemical production potential of O₃ considering the titration effect) in the 2020s and 2060s. It should be noted that the growth of O₃ in the non-warm season in the 2060s for BTH and the YRD will be accompanied by minimal change of O_x, while the declines of O₃ and O_x will appear simultaneously in the warm season for the three regions and in the non-warm season for the PRD. This indicates that the growth of non-warm season O₃ in BTH and the YRD should result partly from NO_x reduction and thereby weakened NO titration, as titration is a key pathway of O₃ loss when the chemical reactivity is relatively low in winter (Gao et al., 2013; Akimoto and Tanimoto, 2022). The differentiated O₃ responses to precursor reduction between YRD and PRD were also detected during the COVID-19 breakout period. With the O₃ isopleth plots, Wang et al. (2021b) illustrated that 40 %–60 % reduction in NO_x and VOCs enhanced the O₃ formation in YRD under the VOC-limited regime but suppressed O₃ in PRD under the transitional regime (a regime that is between NO_x- and VOC-limited). Therefore, VOC emission controls should be better addressed for O₃ pollution alleviation when it expands to the non-warm season in the future.

3.3.2 Meteorological condition change

As shown in Fig. 5b and f, the influence of meteorological change exhibits different patterns for the warm and non-warm season. In the warm season, meteorological change due to global warming will play a positive role on O₃ formation in most of China, with the enhancement within 0–4 ppb, but it will reduce the O₃ level in remote areas like the Tibetan Plateau. The national average growth will be 0.3 ppb and that for the YRD, the PRD, and BTH will be 1.9, 0.7, and 0.3 ppb, respectively. The response of O₃ to meteorological change is associated with some specific variables (Hong et al., 2019). For example, the great enhancement of O₃ in YRD might be attributable to a hotter, drier, and more stable atmosphere, with a growth in T_max (over 0.6 °) and decline in RH and WS (Fig. 2). The result is similar to Hong et al. (2019), which reported a change of 2–8 ppb of daily 1 h maximum O₃ concentration for the peak season from the 2010s to 2050s under RCP4.5. In addition, the declining O₃ in Tibetan Plateau and the surrounding areas might result partly from the weakened long-range transport of peroxyacetyl nitrate (PAN; the principal NO_x reservoir) from the polluted areas (Fischer et al., 2014). Driven by the elevated temperature, PAN from relatively polluted regions will undergo stronger thermal decomposition locally and thus fail to be transported far away to the remote regions to promote O₃ formation (Liu et al., 2013; Lu et al., 2019).

The influence of meteorological change on O₃ production is predicted to be much smaller for the non-warm season, with the magnitude within ±1 ppb in most areas and nationwide average at −0.2 ppb. In the three developed regions, the changes are predicted to range from −0.4 to 0.3 ppb, with little regional difference. The limited influence might be attributable to the modest change in temperature and RH in the non-warm season.

3.3.3 BVOCs and surrounding anthropogenic emission change

Compared to domestic emissions, a change in BVOC emissions and anthropogenic emissions from surrounding countries will have less of an influence (within ±3 ppb) on surface O₃ concentration in China. BVOC change tends to enhance O₃, while foreign emission change tends to restrain it in most areas (Fig. 5).

The growing BVOC emissions due to vegetation and climate change are estimated to enhance O₃ concentration by 0–3 ppb in the most areas of China, with a larger influence of 0.6 ppb in the warm season than that of 0.3 ppb in the non-warm season across the country (Fig. 5c and g). In the warm season, relatively large growth of O₃ concentration will occur in BTH at 2.1 ppb, and in the YRD and PRD, it will be 1.5 and 1.0 ppb, respectively. The abundant NO_x emissions in BTH are expected to result in a larger O₃ concentration response to BVOC emission change than in the YRD and PRD; even the BVOC emission change in BTH will be smaller than the other two regions (Table 2). The result is in agreement with other numerical simulation experiments. Liu et al. (2019) reported a prominent O₃ enhancement even with a low BVOC emission rate under RCP8.5, in a NO_x-abundant environment. In remote areas like the Tibetan Plateau and Northeastern China, the increased BVOCs will remove O₃ due to the isoprene ozonolysis in a low-NO_x environment (Hollaway et al., 2017; Zhu et al., 2022). In general, regions with higher O₃ pollution levels and NO_x emissions will suffer more risk of O₃ growing from rising BVOC emissions in the future.

Most areas of China will benefit from the foreign emission change in terms of O₃ pollution alleviation (Fig. 5d and h). An exception is the Tibetan Plateau and its surrounding areas, which will be affected by the elevated emissions of NO_x and VOCs from south Asia under SSP2-4.5. Limited by the range of pollutant transport, greater impacts will be found for coastal and border areas and less for inland areas (Ni et al., 2018). Larger O₃ changes in the three developed regions are predicted than in the whole country, benefitting from the precursor emission reduction in east Asia and southeast Asia.

Figure 6 The relationships between the separate MDA8 O₃ changes attributed to the four factors (denoted by the name of Case 3–6) and the total changes from the 2020s to 2060s over China and the three regions. Panels (a)–(d) represent the warm season, and panels (e)–(h) represent the non-warm season. The relative contributions of the four factors to the total influence of future change are shown in the light-grey box.

3.4 The relationship between the joint and separate effects of multiple factors

Figure 6 summarizes the contributions of individual factors to the total O₃ change by region and season. Due to the nonlinear response of O₃ to multi-factor changes, the aggregated contribution of the four factors does not equal the joint contribution (i.e., there are gaps between the difference of the 2020s and 2060s and the aggregated contribution of four factors).

The varying domestic anthropogenic emissions are predicted to dominate the change of the future O₃, with a relative contribution ranging from 75 % to 117 % for different regions and seasons. The relative contributions of the other three factors are estimated to be limited within ±25 % at national and regional level. Among different regions, YRD will be more affected by climate change with the contribution of −13 % and −12 % for the warm and non-warm season, respectively, far greater than that of BTH and PRD (−6 % to 0 %). BTH will be more affected by BVOC emission change than other regions in the warm season (−21 %), while the YRD and PRD will be more affected in the non-warm season with the relative contributions of 17 % and −20 %, respectively. Little regional difference is found for the relative contributions of foreign emission change.

Figure 7 The relationships between simulated hourly NO₂ and O₃ concentrations with the lognormal fits for different regions and seasons. The colored circles, representing different cases, come from the seasonal average concentrations for each grid in the target region. The specific circles with black border represent the regional average situation, and the turning points of every fitted curve are marked by the “+” sign. The density plots of the 2020s and 2060s are shown in the inset on the left.

To better understand the regional and seasonal differences of the relative contributions of future changes to O₃ concentration, we examine the nonlinear response of O₃ to precursor change in the three developed regions. We follow Chen et al. (2021) and Schroeder et al. (2017) and conduct a fit of lognormal distribution for the relationship of modeled hourly O₃ and NO₂ concentrations, as shown in Fig. 7. The data points on the left of the turning point of the fitted curve suggest a NO_x-limited regime and on the right a VOC-limited regime, and data points around the turning point are under a transitional regime.

The O₃–NO₂ relationship from the 2020s to 2060s will be mostly influenced by the changing domestic anthropogenic emissions, indicated by the close distributions of data points and fitted curves between “EMIS” and “2060s” in Fig. 7. In the warm season, the future O₃–NO₂ relations in BTH and YRD are predicted to change greatly from a highly O₃ polluted situation with moderate NO₂ concentration to a situation with a relatively low level of NO₂ (mostly under 10 ppb) and a moderate level of O₃ (under 60 ppb). A weak VOC-limited regime appeared for the whole BTH in the 2020s, consistent with recent observation-based analysis (Chen et al., 2023; Kong et al., 2024). There is big diversity within the region, including a dense area with a strong VOC-limited regime and other areas with a transitional or NO_x-limited regime (Fig. 7a). Represented by the moving of most points from the right of the turning point to near or left of the turning point, the NO_x-limited and transitional regimes will dominate BTH in the 2060s. Compared to 2020s, the data points of 2060s are more closely distributed, indicating a reduced diversity of O₃ formation regime in the region. For YRD, most areas were under a transitional or weak VOC-limited regime in the 2020s with limited diversity within the region, and the situation in the 2060s will be similar to that of BTH (Fig. 7a and b). The shift from a weak VOC-limited regime in the 2020s to a transitional or NO_x-limited regime in the 2060s for BTH and the YRD implies the influence of emission reduction on altering the sensitivity of O₃ formation to precursors. Most areas of PRD in the 2020s are under transitional or NO_x-limited regimes, and the regime will transfer to a strong NO_x-limited one in the 2060s, with an almost positive correlation between NO₂ and O₃ in a low-NO₂ environment (Fig. 7c). In the non-warm season, O₃ and NO₂ will remain negatively correlated for BTH and YRD till the 2060s, which suggests a persistent VOC-limited regime and explains the O₃ concentration growth along with substantial precursor emission reductions. The turning points are simulated at extremely low NO₂ concentrations of 2.0 and 1.2 ppb for BTH and YRD, respectively (Fig. 7d and e). A big challenge still exists on effective emission controls to reduce the O₃ concentration in the non-warm season for the two regions. Differently, the O₃ formation sensitivity in most of the PRD will shift from a transitional regime towards a more NO_x-limited situation (Fig. 7f). A simple comparison with the O₃ evolution and its precursor emission changes in developed country provided more implications for policy. According to Chen et al. (2021) and the US Environmental Protection Agency, the northeastern United States experienced a rapid surpassing of the turning point of O₃ formation sensitivity during the 1990s–2010s, with approximately 60 % reductions in both anthropogenic NO_x and VOC emissions. In BTH, the emissions are predicted to decline 57 % and 36 % for NO_x and NMVOCs during 2020s–2060s under SSP2-4.5. Therefore, more ambitious reductions in NMVOCs will be necessary (ideally double the current projected abatement under SSP2-4.5) to accelerate the shift in the O₃ chemical regime for BTH.

The fitted curves of other three factors are similar to those of the 2020s, and the change of these factors will make little difference on NO₂ concentration but will result in moderate changes on O₃ concentration within ±2 ppb. The limited changes of climate, BVOC emissions, and foreign anthropogenic emissions will not essentially alter the O₃ formation regime but may change the O₃ production under nearly the same NO₂ concentration. Changes in individual meteorological factors are expected to easily influence the O₃ and NO₂ concentrations (Pope et al., 2015; Liu and Wang, 2020; Dewan and Lakhani, 2022). The little modeled response of NO₂ to meteorological change, except that in the non-warm season for BTH, might be attributed to the compensating effect of different variables. The limited influence of BVOCs on the O₃ formation sensitivity to precursors is consistent with Gao et al. (2022), which reported comparable empirical kinetic modeling approach (EKMA) curves with and without BVOC emissions. The transboundary O₃ pollution results from the transport of both O₃ and its precursors (mainly associated with PAN), while NO₂ is less influenced by long-range transport due to its shorter lifetime (Ni et al., 2018; Yin et al., 2022).

The change in O₃ formation regime might partly explain the finding that the joint effect of multiple factors on restraining O₃ pollution will be larger than the aggregated effects of individual factors. Under a NO_x-limited regime, O₃ is less sensitive to changing VOC emissions (e.g., BVOC emissions) than under a VOC-limited one. Therefore, the enhancement of O₃ due to BVOC emission growth in the future will be restrained with a much lower NO₂ concentration. This indicates a co-benefit of reducing the anthropogenic emissions to restrain the potential O₃ pollution elevation due to growing BVOC emissions (as a part of the climate penalty) in the future.

Figure 8 Projected annual O₃ exceedance over the east of China in the 2020s and 2060s and the exceedance changes when the four factors at 2060s level. Regional mean changes across CHN, BTH, the YRD, and the PRD are shown in the inset on the left.

3.5 Change of O₃ exceedance events over the east of China

Figure 8 shows the “O₃ exceedance events” over the east of China (mainly including Northern, Eastern, Central, and Southern China) in the 2020s and 2060s, and the changes influenced by different factors. The exceedance is defined as number of days with the MDA8 O₃ exceeding the Chinese National Air Quality Standard Grade II (160 µg m⁻³ or 81.6 ppb). The exceedance events mainly appear in the warm season (Fig. S7). Areas with frequent exceedance (over 50 d) in the 2020s were mainly located in Northern China. Much fewer exceedances are found for YRD and PRD (19.3 and 8.2 d in the 2020s, respectively). In the 2060s, the O₃ exceedance events will drop significantly. The exceedance days will be fewer than 10 d for most of the country, except for some areas in BTH which will still have more than 20 exceedance days over the year.

Domestic emission abatement will be the most important factor reducing the O₃ exceedance, particularly in Northern China. The exceedance days will be cut by 45.3, 19.1, and 8.1 d for BTH, the YRD, and the PRD, respectively, with the maximum reduction reaching 80 d within BTH and the YRD. Notably, the spatial pattern of changing O₃ exceedance due to emission reduction is different from that of changing MDA8 O₃ due to emission reduction as shown in Fig. 5a. Even the warm-season MDA8 O₃ concentration of BTH will decline only 9.7 ppb, the O₃ exceedance events will be greatly reduced, indicating that national emission controls will be especially effective in reducing serious O₃ pollution. Climate change will mainly affect Jiangsu, Anhui, Henan, and Hebei provinces, elevating the exceedance by more than 15 d in most of these areas. For the YRD and PRD, climate change will elevate the exceedance by 9.5 and 3.3 d, respectively. Some areas of BTH will benefit from climate change, with the exceedance declining 0–10 d. The influences of BVOCs and foreign emission change on exceedance days are of limited regional differences, with a growth of 5 to 15 d for the former and a decline of −5 to 0 d for the latter. The exceedances elevated by BVOC emission growth will be 6.6, 6.1, and 2.8 d for BTH, the YRD, and the PRD with the maximum reaching 19, 18, and 12 d within the region, respectively, reflecting a role of biogenic source change in future O₃ episodes that cannot be neglected.

4 Conclusions

We explore the response of China's surface O₃ concentration to the future changes of multiple factors under SSP2-4.5, based on a series of sensitivity experiments with WRF–MEGAN–CMAQ simulations. From the 2020s to 2060s, the MDA8 O₃ concentration is predicted to decline by 7.7 and 1.1 ppb in the warm and non-warm season, respectively, and the O₃ exceedances of the Chinese National Air Quality Standard (Grade II) will be largely eliminated. In the warm season, MDA8 O₃ in BTH, the YRD, and the PRD will decline by 9.7, 14.8 and 12.5 ppb, respectively, larger than the national average level. However, MDA8 O₃ will increase in BTH and the YRD in the non-warm season, attributed to the reduced NO_x emissions and thereby a weakened titration effect. The O₃ pollution will expand towards the non-warm season in the future, bringing a new challenge for policy makers to optimize the strategy of precursor emission controls based on local conditions.

Reduction of local anthropogenic emissions is estimated to dominate the spatial distribution and magnitude of future O₃ change. Meteorological variation will lead to a change of MDA8 O₃ ranging between −1 and 4 ppb for most areas in the warm season. The influences of changing BVOCs and foreign anthropogenic emissions will be within ±3 ppb, with the former elevating O₃ and the latter reducing O₃. Especially in areas with high O₃ pollution and intense NO_x emissions, the growing BVOC emissions will enhance the risk of O₃ pollution even more. The joint effect of multiple factors on restraining O₃ pollution will be larger than the aggregated effects of individual factors, which can be partly explained by the changing O₃ formation regime. A large amount of emission reduction under SSP2-4.5 will reshape the O₃ formation sensitivity to precursors. In BTH and the YRD, O₃ formation in the warm season is projected to shift from a weak VOC-limited to transitional or NO_x-limited regime, while a VOC-limited regime will still dominate in the non-warm season. In the future, O₃ will be less sensitive to BVOC change in a low-NO_x environment, along with persistent emission controls, highlighting the benefit of anthropogenic emissions abatement for mitigating the climate penalty and limiting O₃ pollution.

Limitations exist in current study. Firstly, the future climate data are taken from one single model CESM, subject to bias in the assessment of meteorological influence on O₃. Secondly, some factors that will influence future O₃ level are not included in our analyses, such as the changing CH₄ concentration, increasing soil NO_x emissions, and the stratosphere–troposphere exchange of O₃. For example, the hot spot of soil NO_x emissions in Northern China is also the region with large reduction of anthropogenic NO_x emission but relatively small decline in O₃ concentrations. Under a warmer climate, a growing trend of soil NO_x emissions is expected for the future and may thus present an additional challenge for O₃ pollution alleviation. Thirdly, there are gaps between the downscaled and realistic conditions of meteorology for the 2020s, leading to uncertainty in the O₃ simulation. Finally, the changing O₃ formation regime is presented through the relation between O₃ and NO₂ concentrations, and the mechanism how the climate penalty will influence O₃ formation under substantial reduction of anthropogenic emissions needs to be better analyzed in future studies.