Enhanced atmospheric oxidation and particle reductions driving changes to nitrate formation mechanisms across coastal and inland regions of north China

Nitrate (${\mathrm{NO}}_{\mathrm{3}}^{-}$) has surpassed sulfate as the dominant secondary inorganic ion, posing a significant challenge to air quality improvement measures in China. We utilized the WRF–CMAQ (Weather Research and Forecasting and Community Multiscale Air Quality) model and isotope analysis to investigate the nitrate formation mechanisms driving regional changes in inland and coastal cities in north China during the winters of 2013 and 2018. Among the nitrate formation pathways, the oxidation reaction of OH radicals with NO₂ (OH+NO₂) and the heterogeneous reaction of N₂O₅ (hetN₂O₅) were determined to be the dominant pathways (88 %–95.5 %), whereas the other pathways contributed less than 12.4 % to the total amount of nitrate formation. In inland cities, 63.7 %–85.6 % of nitrate was formed via OH+NO₂ and 8.3 %–27.7 % was formed from hetN₂O₅. In coastal cities, approximately half of the nitrate (48.2 %–56.5 %) was produced from OH+NO₂, whereas hetN₂O₅ contributed 37.0 %–45.7 % due to higher N₂O₅ concentrations and longer NO₃ radical lifetimes. Compared with that in 2013, the OH+NO₂ contribution in 2018 increased by 7.6 % in inland cities and 3.6 % in coastal cities due to the increased atmospheric oxidizing capacity. Scenario simulations indicated that a 60 % reduction in NO_x emissions led to a 4.5 % decrease in nitrate concentrations in Beijing. The reduction reached 32.4 % reduction in Qingdao. A 60 % combined reduction in NH₃, NO_x, and volatile organic compounds (VOCs) yielded 44.2 % and 60.0 % reductions in nitrate in Beijing and Qingdao, respectively, underscoring the necessity of multipollutant control strategies.

1 Introduction

As a key component of fine particulate matter (PM_2.5), nitrate (${\mathrm{NO}}_{\mathrm{3}}^{-}$) exacerbates health risks (Tie et al., 2009; Sun et al., 2014) and affects the physical and chemical properties of particulates, such as their hygroscopicity, light absorption, and acidity (Cao et al., 2013; Wang et al., 2023). These properties directly influence the atmospheric radiation balance (Ramanathan and Feng, 2009; Tegen et al., 2000), atmospheric visibility, and air quality. Furthermore, nitrate serves as a key source of cloud condensation nuclei, affecting cloud formation and precipitation patterns which subsequently influence the global water cycle and climate regulation (Yu et al., 2020; Kalkavouras et al., 2019). Moreover, the photolysis process of nitrate in atmospheric boundary layers is highly active, serving as an important source of NO_x and regulating the atmospheric oxidation capacity, thus influencing the formation of secondary pollutants such as sulfate (${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$) and brown carbon (BrC) (Ye et al., 2016; Xue et al., 2019; Zheng et al., 2020; Yang et al., 2021). Therefore, nitrate is closely related to regional haze pollution occurrence (Zhai et al., 2021; Zhang et al., 2021; Xu et al., 2019a, b; Fu et al., 2020). Understanding the formation mechanism of nitrate is essential not only for advancing the atmospheric chemical processes but also for mitigating regional haze pollution.

In tropospheric atmospheres, nitrate formation primarily follows two pathways (Fig. 1). During the daytime, NO₂ is oxidized by hydroxyl radicals (OH) to produce gaseous HNO₃ (Reaction R1). Conversely, at night, O₃ oxidizes NO₂, leading to the formation of the NO₃ radical, which combines with NO₂ to form N₂O₅. This compound can subsequently be adsorbed onto particles via heterogeneous reactions, resulting in the formation of HNO₃ (Reaction R2) (Hallquist et al., 1999; Pathak et al., 2011). These dynamics highlight the intrinsic dependence of atmospheric nitrate formation on NO_x and oxidants, especially O₃ and OH radicals. The generation of OH and O₃ is intricately linked to the photochemical reactions of NO_x and volatile organic compounds (VOCs) (Atkinson, 2000). OH radicals and hydroperoxyl radicals (HO₂) are produced via the photolysis of ozone (O₃), nitrous acid (HONO), oxygenated volatile organic compounds (OVOCs), and hydrogen peroxide (H₂O₂) and via reactions between O₃ and VOCs. OH subsequently reacts with VOCs to generate organic peroxy radicals (RO₂) and HO₂, which are then recycled back to OH via their interaction with nitrogen monoxide (NO). During this cycle, NO is transformed into NO₂, which yields O₃ upon photolysis (Fu et al., 2020).

Figure 1 Formation mechanisms of nitrate in the atmosphere (Wang et al., 2023).

Nitrate formation mechanisms exhibit significant regional differences, especially in winter. Some studies based on isotope analysis have shown that the OH pathway is important for nitrate formation in inland cities. For example, in Beijing, the OH pathway contributes 66 %–92 % to winter nitrate formation, whereas the heterogeneous reaction of the N₂O₅ (hetN₂O₅) pathway contributes 8 %–34 % to nitrate (Chen et al., 2020). In Xi'an, the contribution of the hetN₂O₅ pathway to winter nitrate formation ranges from 13 % to 35 % (Wu et al., 2021). In coastal areas, which are influenced by high humidities, high sea salt levels, and the combined effects of marine emissions and air masses (Zhong et al., 2023; Athanasopoulou et al., 2008; Zhao et al., 2024), the contribution of the OH pathway is smaller than that in inland cities, whereas the hetN₂O₅ pathway plays a more significant role. For example, isotope analysis revealed that the contribution of the OH+NO₂ pathway to winter nitrate formation is 48 %–74 % in Shanghai (He et al., 2020), whereas the contributions of the OH, hetN₂O₅, and NO₃+HC pathways to winter nitrate formation are 20.2 %, 38.2 %, and 21.6 %, respectively, in Xiamen (Li et al., 2022). Broad coastal studies have revealed similar significant regional differences. Michalski et al. (2003) noted that in coastal California, approximately 90 % of winter nitrate originates from the hetN₂O₅ pathway. In the polar environment of Greenland (Kunasek et al., 2008), nearly all nitrate is formed via the ${\mathrm{NO}}_{\mathrm{3}}+\mathrm{HC}/\mathrm{DMS}$ and hetN₂O₅ pathways in winter, with contributions of 60 % and 40 %, respectively. Furthermore, the impacts of continental and marine air masses on nitrate formation mechanisms differ. In the south Yellow Sea, under the dominance of marine air masses, the contributions of the OH+NO₂, ${\mathrm{NO}}_{\mathrm{3}}+\mathrm{HC}/\mathrm{DMS}$, and hetN₂O₅ pathways to nitrate formation are 43.9 %, 22.4 %, and 33.6 %, respectively, whereas the contributions (12.6 %, 59.0 %, and 28.4 %, respectively) in the Bohai and north Yellow seas are affected by continental pollution (Zhao et al., 2024). These studies have confirmed substantial regional differences in nitrate formation mechanisms.

In addition to isotope analysis, model simulation studies have provided significant support for elucidating nitrate formation mechanisms. Global models such as GEOS-Chem have revealed significant land–ocean differences in nitrate formation pathways. Alexander et al. (2020), using the GEOS-Chem model at a global scale, reported that the OH+NO₂ reaction (41 %–42 %) and N₂O₅ hydrolysis (28 %–41 %) were the main pathways involved below 1 km altitude, with N₂O₅ hydrolysis dominating in winter in the mid- to high-latitude northern continents; however, in remote marine regions, XNO₃ hydrolysis emerged as the primary pathway. However, owing to their relatively low resolution (typically 0.5° × 0.625° or coarser) and simplified treatment of aerosol–surface heterogeneous reactions, global models clearly have limitations in capturing regional differences (Heald et al., 2012). Consequently, the average outputs from global models may not accurately reflect the nitrate formation mechanisms of specific regions (Alexander et al., 2020), leading to potential overestimations or underestimations of local nitrate formation, particularly in terms of ammonia emissions and nitrate production rates. Compared with global models, regional models, such as WRF-Chem (Weather Research and Forecasting model coupled with Chemistry) and WRF–CMAQ (Weather Research and Forecasting and Community Multiscale Air Quality model), can offer greater accuracy in capturing local chemical processes. For example, Sun et al. (2022) used the WRF–CMAQ model to study the Yangtze River Delta (YRD) and reported that during winter, the OH+NO₂ reaction contributed 69.3 % of nitrate formation in Shanghai (a coastal area) and 66.9 % in Hefei (an inland area), whereas the contribution of hetN₂O₅ in Hefei was 27.1 %. Similarly, Li et al. (2021) used the WRF-Chem model to investigate the North China Plain (NCP) and YRD and reported that near the surface, the OH+NO₂ reaction dominated nitrate formation (60 %–92 %), with hetN₂O₅ contributing 8 %–40 %. Moreover, Kim et al. (2014b) applied the WRF–CMAQ model to the Great Lakes region of North America and reported that hetN₂O₅ accounted for 57 % of nitrate formation, highlighting the regional variability. For targeted regions, especially heavily polluted urban areas or those with unique ecosystems, employing high-resolution regional simulations and more detailed chemical mechanisms allows for a more precise characterization of nitrate formation processes.

The North China Plain (NCP) region is a focal point for studying air pollution due to its severe winter haze issues, which consistently rank among the areas with the highest PM_2.5 concentrations globally (Chen et al., 2016). Before 2013, sulfate was the dominant pollutant during haze events in the NCP, accounting for a significant proportion of the total PM_2.5 (Li et al., 2019). However, following the implementation of the Clean Air Action Plan (CAAP) in 2013, which aimed to reduce PM_2.5 concentrations and improve air quality (Zhang et al., 2023a), the composition of winter haze pollution shifted from being sulfate-dominated to nitrate-dominated (Li et al., 2019; Xu et al., 2019a). Data from the Multiresolution Emission Inventory for China (MEIC) (Zheng et al., 2018) show that from 2013 to 2017, anthropogenic emissions of SO₂, NO_x, CO, and NH₃ decreased by 59 %, 21 %, 23 %, and −3 %, respectively. Correspondingly, the national annual average PM_2.5 level decreased by 30 %–50 % from 2013–2018 (Zhai et al., 2019). While these emission reduction measures led to significant declines in CO and SO₂ levels and in particulate sulfate concentrations, the reduction in particulate nitrate levels was smaller than expected. In fact, nitrate mass concentrations have increased in some major cities (Xu et al., 2019b; Shao et al., 2018; Zhou et al., 2019; Fu et al., 2020). This divergence in the responses of sulfate and nitrate to emission controls underscores the complex chemistry of nitrate formation and its precursors. Most studies have focused on nitrate formation mechanisms in inland cities such as Beijing–Tianjin–Hebei (Fan et al., 2020b; Fu et al., 2020; Yang et al., 2024; Chen et al., 2020); thus, the formation mechanisms in coastal cities remain poorly understood. The coastal regions of China experience monsoon transitions, where marine air masses are transported inland under the influence of sea and land breezes, affecting the formation of nitrate over terrestrial areas.

To understand the factors controlling nitrate concentration variations in inland and coastal cities, especially during periods dominated by different primary pollutants, we selected the winters of 2013, when the CAAP was initiated, and 2018, after the project had concluded, to explore the nitrate formation mechanisms in five inland cities and two coastal cities on the NCP via a regional multiscale model combined with isotope analysis. The primary objectives of this study were to quantify the contributions of different nitrate formation pathways, such as OH+NO₂ and hetN₂O₅, in inland and coastal cities, investigate the impacts of emission control measures on nitrate formation, and provide insights into regional differences in nitrate control.

2 Methods

2.1 Model configuration

The Community Multiscale Air Quality (CMAQ) model (version 5.3.3) was used to simulate the chemical reactions and physical processes that contributed to ${\mathrm{NO}}_{\mathrm{3}}^{-}$ and TNO₃ (${\mathrm{HNO}}_{\mathrm{3}}+{\mathrm{NO}}_{\mathrm{3}}^{-}$) formation in the NCP region. Two distinct periods were considered, i.e., one period from 1 December 2013 to 28 February 2014, representing the period when the CAAP was initiated, and one period from 1 December 2018 to 28 February 2019, representing the period after the CAAP was completed. As shown in Fig. 2, a dual-layer grid nesting method was implemented in our simulations. Domain 1 (D01) encompassed the majority of China with a 36 km horizontal resolution, whereas domain 2 (D02) covered the NCP with a 12 km horizontal resolution. The major cities within the NCP included Beijing (BJ), Tianjin (TJ), Shijiazhuang (SJZ), Jinan (JN), Zhengzhou (ZZ), Qingdao (QD), and Yantai (YT). The terrain elevation data used in Fig. 2 came from the GEBCO 2024 Grid provided by the GEBCO Bathymetric Compilation Group (2024, https://doi.org/10.5285/1c44ce99-0a0d-5f4f-e063-7086abc0ea0f). The model was vertically segmented into 14 layers, stretching from the Earth's surface to the troposphere, with the first layer height at approximately 31 m. The CB6 chemical mechanism was chosen to simulate gas-phase chemistry (Luecken et al., 2019; Yan et al., 2021) and the aerosol mechanism. Aerosol module version 7 (AERO7; Appel et al., 2021) was used to simulate particulate-phase chemistry. The initial and boundary conditions used as the inputs of D01 were provided via the GEOS-Chem global simulation results, and D02 was provided via D01. To minimize the impact of the initial conditions, we initiated the model run 7 d before the analysis period. Anthropogenic emission data were sourced from the Multiresolution Emission Inventory for China (MEIC) source emission inventory maintained by Tsinghua University (http://meicmodel.org/, last access: 17 June 2024), while biogenic emission estimates were derived from the Model of Emissions of Gases and Aerosols from Nature version 2.0.4 (MEGANv2.0.4, https://bai.ess.uci.edu/megan, last access: 3 June 2024).

Figure 2 Terrain heights of the NCP region and locations of the major cities (Beijing, Tianjin, Shijiazhuang, Jinan, Zhengzhou, Qingdao, and Yantai) within the study domain. The white dots indicate meteorological observation stations, and the blue crosses indicate air quality monitoring stations. Publisher's remark: please note that the above figure contains disputed territories.

The required meteorological data were generated with the Weather Research and Forecasting (WRF) model version 3.7, a system for predicting mesoscale weather patterns. The selection of physics options for this model conformed with the methodologies applied in our earlier studies (Liu et al., 2021; Chen et al., 2021). The WRF simulations in this study depended on the final operational global analysis (FNL) datasets, which are global reanalysis data with temporal and spatial resolutions of 6 h and 1° × 1°, respectively, provided by the U.S. Environmental Prediction Center and the National Center for Atmospheric Research (https://rda.ucar.edu/datasets/ds083-2/, last access: 12 March 2024). The WRF model outputs were subsequently processed via Model-3/CMAQ Modeling System Interface Processor (MCIP) version 4.3, ensuring compatibility with the CMAQ model format.

2.2 Process analysis

Within the CMAQ modeling framework, the integrated reaction rate (IRR) module of the process analysis (PA) tool was used to simulate the formation reactions of TNO₃ (${\mathrm{HNO}}_{\mathrm{3}}+{\mathrm{NO}}_{\mathrm{3}}^{-}$). The reaction rates of the chemical reactions at each moment were provided by the IRR module, enabling the quantitative evaluation of the target reactions. The complex chemical formation of TNO₃ involved seven reaction pathways, which were categorized into three main pathways on the basis of their significance: OH+NO₂ (Reaction R1), hetN₂O₅ (Reaction R2), and other (Reactions R3 to R7). The chemical equations for these seven production reaction pathways could be expressed as follows.

$$\begin{array}{}\text{(R1)}& {\displaystyle }{\mathrm{NO}}_{\mathrm{2}}\left(\mathrm{g}\right)+\mathrm{OH}\left(\mathrm{g}\right)\to {\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{g}\right)\text{(R2)}& {\displaystyle }{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}+{\mathrm{H}}_{\mathrm{2}}\mathrm{O}\left(\mathrm{l}\right)\to \mathrm{2}{\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{aq}\right)\text{(R3)}& {\displaystyle }{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}\left(\mathrm{g}\right)+{\mathrm{H}}_{\mathrm{2}}\mathrm{O}\left(\mathrm{g}\right)\to \mathrm{2}{\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{g}\right)\text{(R4)}& {\displaystyle }{\mathrm{NO}}_{\mathrm{3}}\left(\mathrm{g}\right)+\text{VOCs(g)}\to {\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{g}\right)\text{(R5)}& {\displaystyle }\mathrm{2}{\mathrm{NO}}_{\mathrm{2}}\left(\mathrm{g}\right)+\mathrm{2}{\mathrm{H}}_{\mathrm{2}}\mathrm{O}\left(\mathrm{l}\right)\to {\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{aq}\right)+\mathrm{HONO}\left(\mathrm{aq}\right)\text{(R6)}& {\displaystyle }{\mathrm{ClNO}}_{\mathrm{3}}\left(\mathrm{g}\right)+{\mathrm{H}}_{\mathrm{2}}\mathrm{O}\left(\mathrm{l}\right)\to {\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{aq}\right)+\mathrm{HOCl}\left(\mathrm{aq}\right)\text{(R7)}& {\displaystyle }{\mathrm{RONO}}_{\mathrm{2}}\left(\mathrm{g}\right)+{\mathrm{H}}_{\mathrm{2}}\mathrm{O}\left(\mathrm{l}\right)\to {\mathrm{HNO}}_{\mathrm{3}}\left(\mathrm{aq}\right)+\mathrm{ROH}\left(\mathrm{l}\right)\end{array}$$

2.3 Observation data

In this study, the simulation results, including meteorological parameters and atmospheric pollutants, for the major cities in the NCP were validated in depth. The meteorological parameters included the 2 m temperature (T2), relative humidity (RH), 10 m wind speed (WS10), and 10 m wind direction (WD10). The atmospheric pollutants included the maximum 8 h moving average ozone (MDA8 O₃), NO₂, PM_2.5, and PM_2.5 components.

2.3.1 Meteorological data and sources

Meteorological data, comprising T2, WS10, WD10, and RH, were sourced from the National Climatic Data Center (NCDC) of the National Oceanic and Atmospheric Administration (NOAA) (https://www.ncdc.noaa.gov/cdo-web/, last access: 17 October 2023). Data were collected from a total of 68 stations, with a temporal resolution of either 3 or 1 h.

2.3.2 Atmospheric pollutant concentrations

Data on O₃, NO₂, and PM_2.5 concentrations were retrieved from the real-time national urban air quality dataset disseminated by the China Environmental Monitoring Center (https://air.cnemc.cn:18007/, last access: 13 February 2024). The ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, ${\mathrm{NH}}_{\mathrm{4}}^{+}$, and ${\mathrm{NO}}_{\mathrm{3}}^{-}$ concentration data were sourced from our own observations in Qingdao and published data (detailed information can be found in Table S3).

In the city of Qingdao, we collected total suspended particle (TSP) samples during two winter seasons. The TSP samples were collected on preheated quartz fiber filters with a high-volume (1.05 m³ min⁻¹) aerosol sampler (Qingdao Laoshan Electronics Co., Ltd., China). The sampling site was situated on the roof of Darwin Hall (36°06^′ N, 120°33^′ E; 16 m) at the Ocean University of China, approximately 1 km from the beach. Detailed information can be found in Ding et al. (2019). Water-soluble inorganic ions, such as ${\mathrm{NO}}_{\mathrm{3}}^{-}$, ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, and ${\mathrm{NH}}_{\mathrm{4}}^{+}$, were extracted from the samples using ultrapure water (> 18.2 MΩ cm), and their concentrations were measured with a Thermo Scientific Dionex ICS-1100 IC system, as described in a previous study (Qi et al., 2020). The ionic concentrations in the TSP samples were calibrated by subtracting the concentrations in the blanks. The δ¹⁸O and δ¹⁵N values of ${\mathrm{NO}}_{\mathrm{3}}^{-}$ in the TSP samples were determined via the bacterial denitrifier method (Casciotti et al., 2002; Sigman et al., 2001) with the Gasbench-IRMS system (Delta V model, Thermo Scientific). For analytical calibration, four international reference standards were employed, namely USGS32, USGS34, USGS35, and IAEA-${\mathrm{NO}}_{\mathrm{3}}^{-}$. The standard deviations of the replicates were ±0.2 ‰ for δ¹⁵N–${\mathrm{NO}}_{\mathrm{3}}^{-}$ and ±0.5 ‰ for δ¹⁸O–${\mathrm{NO}}_{\mathrm{3}}^{-}$. The analytical procedures for determining δ¹⁵N–${\mathrm{NO}}_{\mathrm{3}}^{-}$ and δ¹⁸O–${\mathrm{NO}}_{\mathrm{3}}^{-}$ were detailed in a previous study (Luo et al., 2021). The blank contributed less than 5 % to the sample results. Additional details are provided in Sect. S1 in the Supplement.

2.4 Dual-isotopic Bayesian mixing model

The Bayesian isotope mixing model (Stable Isotope Analysis in R, SIAR) and dual-isotopic compositions (δ¹⁵N–${\mathrm{NO}}_{\mathrm{3}}^{-}$ and δ¹⁸O–${\mathrm{NO}}_{\mathrm{3}}^{-}$) determined in Qingdao and adapted from a reference study in Beijing (Zong et al., 2020; Fan et al., 2020b) were employed to estimate the contributions of OH radicals and the hetN₂O₅ pathway to particulate ${\mathrm{NO}}_{\mathrm{3}}^{-}$. Equations (S1)–(S7) for calculating the isotopic compositions of nitrogen oxides in the atmosphere are provided in the Supplement. The atmospheric δ¹⁵N–${\mathrm{NO}}_{\mathrm{3}}^{-}$ and δ¹⁸O–${\mathrm{NO}}_{\mathrm{3}}^{-}$ can be expressed by Eqs. (S1) and (S2), respectively. The endmembers of [δ¹⁵N−HNO₃]_OH, $[{\mathit{\delta }}^{\mathrm{15}}\mathrm{N}-{\mathrm{HNO}}_{\mathrm{3}}{]}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$, [δ¹⁸O−HNO₃]_OH, and $[{\mathit{\delta }}^{\mathrm{18}}\mathrm{O}-{\mathrm{HNO}}_{\mathrm{3}}{]}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ are expressed by Eqs. (S3), (S4), (S5), and (S6), respectively (Walters and Michalski, 2016). This study was focused on the seasonal scale, where observations of δ¹⁸O–${\mathrm{NO}}_{\mathrm{3}}^{-}$ did not reveal a significant diurnal difference (Feng et al., 2020; Zhang et al., 2022), as shown in Fig. S1 in the Supplement. Therefore, we did not account for the diurnal variation in δ¹⁸O–NO₂ in the calculation of endmembers (Walters et al., 2018; Zhang et al., 2025; Albertin et al., 2021, 2024). The δ¹⁵N values of tropospheric NO_x and the δ¹⁸O values of tropospheric H₂O_(g), NO_x, O₃, and •OH occurred within certain ranges, as described in Sect. S2 and Table S1 in the Supplement. Therefore, the endmembers of δ¹⁵N and δ¹⁸O for the two pathways could be estimated according to fNO₂ (the molar ratio of NO₂ and NO_x) and the isotope fractionation values of nitrogen and oxygen, i.e., ${\mathit{\alpha }}_{{\mathrm{NO}}_{\mathrm{2}}/\mathrm{NO}}$, ${\mathit{\alpha }}_{\mathrm{OH}/{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}$, and ${\mathit{\alpha }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}/{\mathrm{NO}}_{\mathrm{2}}}$. The nitrogen and oxygen isotope fractionations were temperature-dependent and could be estimated via Eq. (S7) and Table S2. On the basis of the endmembers of [δ¹⁵N−HNO₃]_OH, $[{\mathit{\delta }}^{\mathrm{15}}\mathrm{N}-{\mathrm{HNO}}_{\mathrm{3}}{]}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$, [δ¹⁸O−HNO₃]_OH, and $[{\mathit{\delta }}^{\mathrm{18}}\mathrm{O}-{\mathrm{HNO}}_{\mathrm{3}}{]}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$, the contribution of the OH radical formation pathway (γ) was estimated via the dual-isotopic Bayesian mixing model (Luo et al., 2021). The uncertainty of the Bayesian model mainly stemmed from the standard deviations in the endmember calculations for each formation pathway and the posterior probability distributions output by the model (Moore and Semmens, 2008).

2.5 Model evaluation

The following statistical indicators were used to evaluate the simulation effect: the mean deviation (MB), normalized mean deviation (NMB), normalized mean error (NME), correlation coefficient (R), root mean square error (RMSE), and index of agreement (IOA). The definitions and standards of all the statistical indicators are provided in Table S3. To evaluate model performance, Emery et al. (2001) provided benchmarks of the above statistical indicators for major meteorological parameters such as T2, WS10, and WD10, which are presented in Table 1. Similarly, Emery et al. (2017) and Huang et al. (2021) proposed benchmarks for the concentrations of major air pollutants, including PM_2.5, NO₂, MDA8 O₃, and ${\mathrm{NO}}_{\mathrm{3}}^{-}$, as comprehensively detailed in Table 2.

Table 1 Statistical performance of the modeled meteorological parameters in the NCP during the winters of 2013 and 2018 (68 sites).

^* T2 denotes the 2 m temperature. RH denotes the relative humidity. WS10 denotes the 10 m wind speed. WD10 denotes the 10 m wind direction.

Table 2 Model performance for the major air pollutants in typical cities of the NCP during the winters of 2013 and 2018.

^* The benchmarks for these pollutants, including NO₂ according to US-EPA (2007), MDA8 O₃, and PM_2.5, were proposed by Emery et al. (2017) and Huang et al. (2021).

For the comparison between the model simulations and observational data, the model grid cell corresponding to the geographic coordinates of each observation station was first identified. At each station, the observed values were compared with the simulated values from the corresponding grid cell on a point-by-point, hour-by-hour basis. These paired comparisons were then aggregated across all stations, and statistical evaluation metrics were generated from the entire dataset to assess model performance comprehensively. All data processing and comparisons were conducted on an hourly basis.

2.6 Emission reduction scenario simulation design

Notably, emission reduction scenarios were specifically designed to investigate how reductions in NH₃, NO_x, and VOC emissions affected PM_2.5 nitrate concentrations in Beijing and Qingdao, which are typical inland and coastal cities, respectively, in northern China. The simulations included single-pollutant reduction strategies for NH₃, NO_x, and VOCs and combined reduction scenarios. For each pollutant, emissions were reduced by 20 %, 40 %, and 60 %, with a focus on assessing the resulting influence on nitrate formation. The combined reduction scenarios entailed simultaneous reductions in NH₃, NO_x, and VOCs to evaluate synergistic effects. The nitrate concentration responses to these reductions were analyzed to determine the most effective strategies for controlling winter nitrate levels.

3 Results and discussion

3.1 Model evaluation

3.1.1 Evaluation of meteorological parameters

The statistical model performance for meteorological parameters such as T2, WS10, WD10, and RH in the NCP is summarized in Table 1, covering the winter months from December 2013 to February 2014 (winter 2013) and from December 2018 to February 2019 (winter 2018). The simulated T2 and RH values exhibited satisfactory reproducibility, with temperature simulations exhibiting MB values slightly higher than the recommended threshold (MB ≤ ±0.5). The simulated wind speeds during both winters were slightly overestimated, with the RMSE slightly exceeding the criterion (2 m s⁻¹; Emery et al., 2001), whereas the wind direction results fully satisfied the criterion. The observed wind speed overestimation by the WRF model could be attributed to its inability to accurately capture the impact of high aerosol loadings on shortwave radiation in winter, which likely reduced the near-surface wind speed (Tan et al., 2017; Jacobson and Kaufman, 2006). In general, the meteorological field simulations of the WRF model were reliable and could effectively reveal the changes in various meteorological elements in the NCP.

3.1.2 Evaluation of atmospheric pollutants

As indicated in Table 2, the model exhibited a favorable simulation ability for the concentrations of NO₂, MDA8 O₃, and PM_2.5 pollutant in the NCP region. Specifically, the simulations accurately captured NO₂ concentrations, with NMB values ranging from −30 % to 32 %. Notably, during the winter of 2018, the simulated NO₂ concentration in Beijing was slightly higher than the recommended range (by 2 %), whereas that in other cities remained within this range. The MDA8 O₃ simulations generally exhibited slight underestimation, whereas in Tianjin, Shijiazhuang, and Zhengzhou, the simulations slightly exceeded the benchmark. This difference was primarily due to the uncertainty in anthropogenic VOC emissions (Wang et al., 2014), making accurate ozone simulations particularly challenging (Sun et al., 2022; Yang et al., 2024). The PM_2.5 simulations were generally accurate, although the simulated concentrations in Tianjin during the winter of 2013 were slightly elevated, with an NMB of 32 %, exceeding the standard by 2 %; however, those in the other cities met the standard requirements.

Figure 3 shows a comparison between the simulated and observed ${\mathrm{NO}}_{\mathrm{3}}^{-}$ concentrations in the six NCP cities during the winters of 2013 and 2018. Yantai was excluded from the comparison because of the inadequacy of the observation data. The data sources are listed in Table S4. The comparison results indicated that the simulated ${\mathrm{NO}}_{\mathrm{3}}^{-}$ trends in the NCP were accurate during the winters of 2013 and 2018, with R values of 0.61 and 0.72, respectively, meeting the benchmark of 0.60. However, the simulation for 2013 yielded slightly overestimated values, with an MB of 1.27, corresponding to NMB = 7 % and NME = 51 %. In contrast, the simulation results for 2018 were closer to the observed values, with NMB = 5 % and NME = 51 %. Overall, the model demonstrated suitable stability and accuracy in simulating atmospheric pollutants in the NCP region, providing a solid foundation for future analysis.

Figure 3 Scatter plots of the simulated and observed daily mean ${\mathrm{NO}}_{\mathrm{3}}^{-}$ concentrations in the NCP during the winter seasons of 2013 and 2018; the ^* symbol indicates R values that are statistically significant (P < 0.05).

3.2 Variations in the concentrations of PM_2.5 and its components between 2013 and 2018

During the winter of 2013, the average PM_2.5 concentration in the NCP reached 80.52 ± 77.13 µg m⁻³, with high-concentration areas mainly concentrated in the provinces of Beijing, Tianjin, southern Hebei, Henan, Shandong, and Anhui (Fig. 4). By 2018, the average PM_2.5 concentration had decreased to 50.74 ± 49.26 µg m⁻³, representing a significant reduction of 37 ± 10 %. Previously, scholars reported similar trends in PM_2.5 reduction across China (Li et al., 2023; Zhai et al., 2019). We found that the ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ concentration decreased by 53 ± 5 %, whereas the ${\mathrm{NO}}_{\mathrm{3}}^{-}$ and ${\mathrm{NH}}_{\mathrm{4}}^{+}$ concentrations decreased by 13 ± 20 % and 33 ± 13 %, respectively (refer to Fig. S2 for city-specific changes). Our analysis revealed a pronounced difference in nitrate concentration. Specifically, inland cities experienced a more pronounced reduction in nitrate concentrations (4.6 % to 41.7 %), whereas coastal cities showed relatively small changes, with Yantai experiencing a 4.8 % decrease and Qingdao exhibiting a 3.9 % increase. The substantial reduction in the three inorganic components (${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, ${\mathrm{NH}}_{\mathrm{4}}^{+}$, and ${\mathrm{NO}}_{\mathrm{3}}^{-}$) led to significant changes in the PM_2.5 composition, with the proportion of nitrate increasing relative to that of sulfate (Fig. 5). Nitrate concentrations demonstrated significantly different patterns of reduction between coastal and inland regions, necessitating the identification of their dominant controlling factors to facilitate effective regional pollution control strategies.

Figure 4 Spatial distribution and changes in PM_2.5 and its components in the NCP during the winters of 2013 and 2018. The up arrows indicate increases, and the down arrows indicate decreases.

Figure 5 Spatial distributions of the nitrate and sulfate proportions in PM_2.5 in the NCP region during the winters of 2013 and 2018.

3.3 Rates, contributions, and diurnal variations in nitrate formation reactions in the major cities

The ${\mathrm{NO}}_{\mathrm{3}}^{-}$ concentration slightly decreased from 2013 to 2018 due to a decrease in the TNO₃ production rate (0–0.05 ppb h⁻¹), which could be attributed to the production rate of the hetN₂O₅ pathway (0–0.03 ppb h⁻¹). According to Table 3, the contribution increased in Zhengzhou from 2013 to 2018, whereas it decreased in six cities, including Beijing (−4.4 %), Tianjin (−4.6 %), Shijiazhuang (−5.9 %), Jinan (−7.8 %), Qingdao (−4.5 %), Yantai (−2.2 %), and Zhengzhou (+2.1 %). For the OH+NO₂ reaction pathway, the rate remained unchanged in all but four cities (Beijing, Jinan, Zhengzhou, and Qingdao); it varied between −0.02 and 0.01 ppb h⁻¹ in the other cities. However, the contribution of the OH+NO₂ reaction pathway to TNO₃ formation increased overall due to the decreased rate of the hetN₂O₅ reaction pathway, although the production rate remained almost unchanged. Moreover, the contribution increase was greater in inland cities than in coastal cities. Among the five inland cities (Beijing, Tianjin, Shijiazhuang, Jinan, and Zhengzhou), the OH+NO₂ reaction pathway was the main nitrate formation mechanism, with its contribution ratio ranging from 63.7 %–77.8 % in 2013 and increasing to 72.9 %–85.6 % in 2018. In contrast, in the two coastal cities (Qingdao and Yantai), the OH+NO₂ reaction pathway contributed to 48.2 %–56.5 % of nitrate formation. However, the contribution of the hetN₂O₅ pathway ranged from 37.0 %–45.7 %. Therefore, both pathways were equally important in nitrate formation in coastal cities. These findings highlighted notable differences in nitrate formation mechanisms between inland cities and coastal cities, demonstrating that nitrate formation mechanisms were strongly affected by the ocean.

Table 3 TNO₃ production rates (ppb h⁻¹) in the surface layer and contributions (%) of the major production pathways in the surface atmosphere in seven representative cities of the NCP in 2013 and 2018.

Since the OH+NO₂ and hetN₂O₅ reactions were the major pathways affecting nitrate formation, we comprehensively analyzed the diurnal variations in these two TNO₃ formation reactions in the NCP, considering the inland city of Beijing and the coastal city of Qingdao as examples. Figure 6 shows the average diurnal variations in the TNO₃ production rates via the different pathways in Qingdao and Beijing. We found that the OH+NO₂ reaction pathway exhibited similar diurnal characteristics in these two cities. Specifically, from night to the early morning (00:00 to 08:00 CST (UTC+8)), the reaction rate under this pathway remained almost zero because of the lack of photochemical processes that generated OH radicals; these processes were limited by solar radiation (Liu et al., 2020; Sun et al., 2022; Tan et al., 2021). Observational data from the Pearl River Delta region supported this result, showing that OH radicals only began to accumulate significantly after 06:00 CST (UTC+8) (Hofzumahaus et al., 2009; Lu et al., 2012). After approximately 08:00 CST (UTC+8), the reaction rate rapidly increased, reaching a peak between 11:00 and 13:00 CST (UTC+8) at approximately 0.4 ppb h⁻¹ (∼ 2.8 $\mathrm{\mu }\mathrm{g}{\mathrm{m}}^{-\mathrm{3}}{\mathrm{h}}^{-\mathrm{1}}$), resembling the findings for the NCP and Shanghai (Tan et al., 2021; Sun et al., 2022; Liu et al., 2020). The reaction rate decreased with decreasing sunlight during the afternoon and returned to low levels at night due to a decrease in the number of photochemical processes that occurred. Compared with those during the winter of 2013, the peak reaction rates of the OH pathway decreased from 0.43 to 0.42 ppb h⁻¹ and from 0.45 to 0.44 ppb h⁻¹ in 2018 in Qingdao and Beijing, respectively. The contributions of the OH+NO₂ reaction pathway to TNO₃ formation were in the ranges of 90.2 %–95.1 % and 75.8 %–82.2 % in Beijing and Qingdao, respectively, during the daytime. However, the other cities in the NCP exhibited similar diurnal variations. Therefore, the OH+NO₂ reaction pathway dominated during the daytime in both the inland and coastal cities.

Figure 6 Average diurnal variation in the TNO₃ production rate under different pathways (daytime: 06:00–18:00 CST (UTC+8); nighttime: 18:00–06:00 CST (UTC+8)): (a) Beijing 2013, (b) Qingdao 2013, (c) Beijing 2018, and (d) Qingdao 2018.

At night, the hetN₂O₅ reaction pathway notably affected TNO₃ formation. The rate of the hetN₂O₅ reaction pathway showed similar bimodal diurnal variations in Qingdao and Beijing. However, the peak values occurred at 02:00 and 20:00 CST (UTC+8) in Qingdao, with much higher peak values than those in Beijing, where the peak values occurred from 04:00–05:00 and 18:00 CST (UTC+8). The hetN₂O₅ reaction rate was relatively high at night but decreased rapidly by early morning, reaching its lowest point during the daytime (09:00 to 16:00 CST (UTC+8)). This phenomenon occurred because the N₂O₅ hydrolysis reaction depended heavily on NO₃ radicals and N₂O₅, which were only stable in weakly illuminated environments (Tie et al., 2003a, b; Atkinson et al., 2004).

High concentrations of NO₂ and O₃, precursors for N₂O₅, were considered key in this process (Dentener and Crutzen, 1993).

$$\begin{array}{}\text{(R8)}& {\displaystyle }{\mathrm{NO}}_{\mathrm{2}}+{\mathrm{O}}_{\mathrm{3}}\to {\mathrm{NO}}_{\mathrm{3}}+{\mathrm{O}}_{\mathrm{2}}\text{(R9)}& {\displaystyle }{\mathrm{NO}}_{\mathrm{3}}+{\mathrm{NO}}_{\mathrm{2}}\leftrightharpoons {\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}\end{array}$$

However, owing to high rates of photolysis during the day, NO₃ radicals and N₂O₅ could accumulate only at night (Zhao et al., 2023; Atkinson et al., 2004). In Beijing, the hetN₂O₅ reaction rate was relatively low, with a peak value of only 0.03 ppb h⁻¹ and a nighttime contribution rate of 46.9 %. In contrast, Qingdao presented a peak hetN₂O₅ reaction rate of 0.2 ppb h⁻¹ (approximately 1.4 $\mathrm{\mu }\mathrm{g}{\mathrm{m}}^{-\mathrm{3}}{\mathrm{h}}^{-\mathrm{1}}$) at night, with an average contribution of 75.2 %. This difference arose primarily due to the higher N₂O₅ concentrations and RH levels in Qingdao than in Beijing, as detailed in Sect. 3.5.3. The other inland cities exhibited diurnal variations similar to those in Beijing, whereas the coastal cities presented trends similar to those in Qingdao. Compared with those during the winter of 2013, the reaction rates of the hetN₂O₅ pathway in both cities during the winter of 2018 decreased from 0.19 to 0.18 ppb h⁻¹ in Qingdao and from 0.07 to 0.03 ppb h⁻¹ in Beijing. A comparative analysis revealed significant variations in the nitrate formation pathways between the coastal and inland areas, particularly in the hetN₂O₅ reaction rate (Fig. 7). The factors influencing these variations are discussed in Sect. 3.5.3, which specifically addresses the underlying mechanisms responsible for the differences in hetN₂O₅ reaction rates between coastal and inland regions.

Figure 7 Spatial distributions of the mean reaction rates for the nitrate formation pathways (OH+NO₂ and hetN₂O₅) in the NCP during the winters of 2013 and 2018.

3.4 Assessment of the simulated nitrate formation reactions via the isotopic method

To verify the accuracy of nitrate formation reaction simulations, we conducted a comprehensive isotope analysis in Qingdao and Beijing. The isotope analysis conducted in Qingdao and Beijing for the 2013 and 2018 data (Fig. 8) revealed that the model-simulated and isotope-based contributions of ${\mathrm{NO}}_{\mathrm{3}}^{-}$ formation pathways generally exhibited consistent trends, although some differences were observed. In 2013 in Qingdao, the CMAQ model yielded contribution values of 72.9 % (OH+NO₂) and 18.7 % (hetN₂O₅), whereas the isotope analysis suggested 62.4 ± 22.5 % and 37.6 ± 22.5 %, respectively. By 2018 in Beijing, the model yielded 89.5 % vs. an isotope-derived 76.5 ± 15.8 % for OH+NO₂ and 5.4 % vs. 23.5 ± 15.8 % for hetN₂O₅. Overall, although the differences persisted, the model captured the general trends, and the deviations fell within an acceptable range, indicating that the model simulations generally agreed with the isotope analysis. Several factors could explain these discrepancies. First, the isotope analysis carried intrinsic uncertainties, commonly approximately 10 %–30 %, due to the standard deviations (SDs) in each pathway's endmember calculations and the posterior probability distributions generated by the Bayesian model (Moore and Semmens, 2008). Second, the isotope analysis considered only two principal pathways: OH+NO₂ and hetN₂O₅. Although the other pathways only contributed to a small fraction of ${\mathrm{NO}}_{\mathrm{3}}^{-}$ formation, their omission could cause biases, especially under specific meteorological or pollution conditions. Third, the IRR module did not consider long-range transport effects, which could further contribute to the discrepancies in the isotope results. In brief, the model simulation results align reasonably well with the isotope analysis findings despite some discrepancies between the two methods.

Figure 8 Time series of the contributions of the atmospheric ${\mathrm{NO}}_{\mathrm{3}}^{-}$ formation pathways: (a, b) Beijing 2018, (c, d) Beijing 2013, and (e, f) Qingdao 2018, based on (a, c, e) dual-isotope analysis and (b, d, f) model simulations.

Notably, both the model and isotope results indicated an increasing trend in the OH+NO₂ pathway contribution in Beijing from 2013 to 2018 (Fig. 9), with the model and isotope analysis showing increases of 5.1 % and 9.7 %, respectively. This trend reflected the increased atmospheric oxidation capacity in the region. Regional comparative analysis revealed significant differences in nitrate formation pathways between Qingdao and Beijing (Fig. 9). In 2018, the Qingdao model results showed contributions of 72.9 % and 18.7 % to the total nitrate formation for the OH+NO₂ and hetN₂O₅ pathways, respectively, compared with the isotope results of 62.4 ± 22.5 % and 37.6 ± 22.5 %, respectively. In contrast, the Beijing model results indicated contributions of 89.5 % and 5.4 % to the total nitrate formation for the OH+NO₂ and hetN₂O₅ pathways, respectively, whereas the isotope results were 76.5 ± 15.8 % and 23.5 ± 15.8 %, respectively. Both methods demonstrated that the OH+NO₂ pathway dominated nitrate formation in Beijing, while the hetN₂O₅ pathway played a substantial role in Qingdao, highlighting regional differences in the mechanism of nitrate formation and resulting in different changes in the nitrate concentration.

Figure 9 Comparison of the contributions of the atmospheric ${\mathrm{NO}}_{\mathrm{3}}^{-}$ formation pathways based on the dual-isotope results and model simulations for Beijing in 2013 and 2018 and for Qingdao in 2018.

As discussed above, field observation data could include contributions from long-range transport (Walters et al., 2024); however, the CMAQ model's IRR module confined nitrate formation to individual grid cells. To quantify this discrepancy, we employed a backward trajectory analysis method incorporating temporal and spatial constraints to calculate the source regions and their contributions to the nitrate formation pathway.

The backward trajectory correction significantly altered the estimated nitrate formation pathways (Fig. 10). In Qingdao, the corrected model results revealed an average OH+NO₂ pathway contribution of 33.3 %, which was substantially lower than the isotope analysis result of 62.4 ± 22.5 %, whereas the hetN₂O₅ pathway contribution was 58.4 %, which was significantly greater than the isotope result of 37.6 ± 22.5 %. In Beijing, the corrected model in 2013 indicated an OH+NO₂ pathway contribution of 48.3 % (isotope: 66.8 ± 26.1 %) and a hetN₂O₅ pathway contribution of 38.5 % (isotope: 33.2 ± 26.1 %). In 2018, the corrected model showed an OH+NO₂ pathway contribution of 57.2 % (isotope: 76.5 ± 15.8 %) and a hetN₂O₅ pathway contribution of 35.1 % (isotope: 23.5 ± 15.8 %). Compared with the uncorrected model, the backward trajectory correction significantly improved the estimation of the hetN₂O₅ pathway contribution (Qingdao: from −18.9 % to +20.8 %; Beijing: from −18.1 % to +11.6 %) but exacerbated the discrepancy in the OH+NO₂ pathway contribution (Qingdao: from +10.5 % to −29.1 %; Beijing: from +13.0 % to −19.3 %). The effects of this correction exhibited regional differences: the improvement in the hetN₂O₅ pathway was more pronounced in Qingdao, where regional transport had a greater impact, whereas the overadjustment of the OH+NO₂ pathway was more evident in Beijing, where local processes were dominant.

Figure 10 Backward trajectories and nitrate formation pathways: (a–c) backward trajectories for (a) Beijing 2018, (b) Beijing 2013, and (c) Qingdao 2018; time series of the contributions of the atmospheric ${\mathrm{NO}}_{\mathrm{3}}^{-}$ formation pathways based on backward-trajectory-constrained model simulations for (c) Beijing 2018, (d) Beijing 2013, and (e) Qingdao 2018.

Remarkably, Walters et al. (2024) achieved relatively accurate results using the unmodified IRR module without considering transport effects, reporting an overall RMSE of 2.6 ‰ for Δ¹⁷O(HNO₃) between modeled and observed values. This finding indicated that if the trends in nitrate formation pathways were similar across regions and vertical layers, the impact of transport on the overall simulation accuracy could be limited (Walters et al., 2024). Our results highlighted that in regions with significant spatial heterogeneity in formation pathways, such as Qingdao and Beijing, incorporating transport effects through backward trajectory analysis remained crucial for improving model performance.

These results demonstrated that backward trajectory correction played a positive role in improving the estimation of the hetN₂O₅ pathway, particularly in regions significantly influenced by regional transport. However, the overadjustment of the OH+NO₂ pathway highlighted the need for further optimization of the correction method. Key areas for improvement would include refining OH radical concentration estimates and optimizing the quantification of regional transport. Currently, the overly broad consideration of transport regions could obscure local formation signals, suggesting a need for increasingly precise spatial delineation in transport modeling. These enhancements would be essential for improving the model's applicability across different regions and under varied meteorological conditions. Currently, owing to the lack of successful transport correction methodologies, we discuss the influencing factors without applying correctional procedures.

3.5 Factors influencing nitrate formation pathways

To explore the reasons for the observed coastal–inland differences in nitrate formation, we further investigated the influences of several key factors, including NH₃, NO₂, OH, O₃, RH, NO, and N₂O₅, on these pathways.

3.5.1 Influence of NH₃ on nitrate formation

The reaction of NH₃ with nitric acid (HNO₃) to form ammonium nitrate (NH₄NO₃) was considered an important process in the formation of atmospheric aerosols. In the absence of sufficient NH₃, the material would first react with H₂SO₄ to form (NH₄)₂SO₄ and then with HNO₃ to form NH₄NO₃ (Zhai et al., 2021). The availability of NH₃ was considered a critical factor governing the partitioning of HNO₃ to particulate nitrate (pNO₃). The gas ratio (GR) (Fu et al., 2020; Ansari and Pandis, 1998) was used to study whether NH₃ limits nitrate formation. The GR could be calculated as follows.

$$\text{GR}={\displaystyle \frac{\left(\left[{\mathrm{NH}}_{\mathrm{3}}\right]+\left[{\mathrm{NH}}_{\mathrm{4}}^{+}\right]\right)-\mathrm{2}\times \left[{\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}\right]}{\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\right]+\left[{\mathrm{HNO}}_{\mathrm{3}}\right]}}$$

When the GR exceeded 1, the NH₃ concentration in the atmosphere was considered sufficient. A GR value between 0 and 1 indicated NH₃-neutral conditions, whereas a value less than 0 suggested NH₃-poor conditions, thus limiting ${\mathrm{NO}}_{\mathrm{3}}^{-}$ formation due to insufficient NH₃. In the NCP region, the GR values were generally greater than 2 in 2013 and 2018 (Fig. 11), indicating that the region mainly exhibited a state of sufficient NH₃. Similar phenomena were observed by other researchers in the NCP region (Zhai et al., 2021; Xu et al., 2019c; Li et al., 2018). Therefore, in the NCP region, the formation of NH₄NO₃ was not limited by NH₃ in terms of supply.

Figure 11 Spatial distributions and interannual variations in the GR, HONO, N₂O₅, and NO_x concentrations over the North China Plain during the winters of 2013 and 2018. The percentage changes (diff) represent the relative differences between 2018 and 2013.

The GR in 2018 was greater than that in 2013, indicating an increasing surplus of NH₃. We believed that the variation in GR was caused by the significant reduction in sulfate and slight reduction in nitrate in PM_2.5 (Fig. 12a) because of the minimal changes in NH₃ emissions. This hypothesis could be verified from an emission perspective. The total ammonia emissions in China increased from 9.64 to 9.75 Tg from 2013–2015 and then gradually decreased to 9.12 Tg by 2018 (Liao et al., 2022). Additionally, there were coastal–inland differences in NH₃ surplus levels. Notably, the NH₃ surplus ratio was consistently greater than 3.4 in inland regions, whereas it ranged from 1.7 to 2.5 in coastal cities. This difference was attributed primarily to intensive industrial activities related to fertilizer production in inland areas, which resulted in increased NH₃ emissions. Data showed that NH₃ emissions from fertilizer production in the inland regions of the NCP ranged from 1.0 to 40.0 Gg yr⁻¹, whereas emissions in coastal areas were significantly lower, ranging only from 0 to 0.1 Gg yr⁻¹ (Chen et al., 2022).

Figure 12 (a) GR and (b) N₂O₅ concentrations across the North China Plain and seven major cities during the winters of 2013 and 2018.

In summary, NH₃ played a critical role in influencing particulate nitrate concentrations by affecting the gas-to-particle conversion of HNO₃ in the NCP region, although its availability was sufficient. Furthermore, variations in NH₃ emissions and distinct land–ocean transport patterns significantly affected the formation and spatial distribution of particulate nitrate.

3.5.2 Factors influencing the OH + NO₂ reaction rate

As the OH+NO₂ reaction significantly contributed to TNO₃ formation, the OH and NO₂ concentrations served as key factors. The NO₂ concentration in the NCP was 30.3 µg m⁻³ during the winter of 2013 and decreased to 24.7 µg m⁻³ during the winter of 2018, and its concentration remained several orders of magnitude greater than that of OH radicals. Figure S3 shows that the ${\mathrm{NO}}_{\mathrm{2}}/\mathrm{OH}$ molar ratio in the NCP region generally exceeded 10⁸, indicating that the NO₂ concentration was in excess under emission reduction measures. Therefore, the reaction primarily depended on the OH radical concentration. Next, we examined the impact of the OH radical concentration on the reaction rate. OH radicals were produced primarily via photochemical reactions, such as by the decomposition of HONO into OH and NO under sunlight (Song et al., 2023), which could account for 20 % to 90 % of the total primary production of OH radicals (Song et al., 2023; Xue et al., 2020; Kim et al., 2014a). Additionally, OH radicals could be produced indirectly via O₃ photolysis in the presence of water vapor (Fu et al., 2020; Kim et al., 2014a).

We compare the simulated and observed OH radical concentrations in Table 4. The simulated daily average OH radical concentrations were relatively close to the observed values. However, the peak concentrations of OH radicals were lower than those observed. Underprediction of the OH peak concentration was determined to be a common issue in current model simulations (Czader et al., 2013; Stone et al., 2012; Xue et al., 2020). Although the CMAQ model often underestimated high daytime concentrations (Czader et al., 2013), it satisfactorily captured the diurnal variation in OH radicals, likely because it underestimated HONO concentrations (Zhang et al., 2023b). The data (Fig. 13a–c) revealed that the peak value of OH radicals in Beijing significantly increased by 0.98 × 10⁶ molec. cm⁻³ in 2018 compared with that in 2013. A similar trend was observed in Qingdao, with an increase of 0.58 × 10⁶ molec. cm⁻³. The GEOS–Chem simulations (Zhang et al., 2023a) indicated a stable upward trend in OH radical concentrations across the NCP from 2014 to 2017, ranging from 0.05 to 0.17 × 10⁶ molec. cm⁻³ yr⁻¹. Our study data suggested that despite China's emission reduction measures, these efforts did not effectively control OH radical concentrations. The continuous increase in OH radicals accelerated the nitrate formation rate via the OH+NO₂ pathway.

Table 4 Observations and CMAQ model simulation results for the winter OH radical concentration.

Figure 13 Diurnal variations in OH radicals and O₃ in Beijing and Qingdao and average concentrations across seven major cities in the NCP (2013, 2018): (a) diurnal variation in OH radicals in Beijing, (b) diurnal variation in OH radicals in Qingdao, (c) average OH radical concentration in the seven cities, (d) diurnal variation in O₃ in Beijing, (e) diurnal variation in O₃ in Qingdao, and (f) average O₃ concentration in the seven cities.

One of the primary pathways for the production of OH radicals was the photolysis of nitrous acid (HONO) in the tropospheric atmosphere. This photolytic reaction can be expressed as follows.

$$\begin{array}{}\text{(R10)}& \mathrm{HONO}+hv\left(\mathit{\lambda }<\mathrm{400}\mathrm{nm}\right)\to \mathrm{OH}+\mathrm{NO}\end{array}$$

However, the data indicated that the HONO concentration decreased between 2013 and 2018 (Fig. 11), suggesting that this production pathway did not significantly contribute to the increase in the OH radical concentration.

In the atmosphere, another source of OH was the reaction between water vapor (H₂O) and excited oxygen atoms (O(¹D)). Notably, excited oxygen atoms would typically be produced via O₃ photolysis (Kim et al., 2014a; Tan et al., 2019).

$$\begin{array}{}\text{(R11)}& {\mathrm{O}}_{\mathrm{3}}+h\mathit{\nu }\left(\mathit{\lambda }<\mathrm{340}\mathrm{mm}\right)\to {\mathrm{O}}_{\mathrm{2}}+\mathrm{O}{(}^{\mathrm{1}}\mathrm{D})\end{array}$$

The excited oxygen atoms (O(¹D)) subsequently reacted with water vapor to produce OH radicals.

$$\begin{array}{}\text{(R12)}& {\mathrm{H}}_{\mathrm{2}}\mathrm{O}+\mathrm{O}{(}^{\mathrm{1}}\mathrm{D})\to \mathrm{2}\mathrm{OH}\end{array}$$

Therefore, O₃ photolysis played an indirect but crucial role in the formation process. O₃ could oxidize NO to NO₂, which then could react with OH radicals to form HNO₃, significantly influencing the OH+NO₂ pathway. As shown in Fig. 13d–f, the increase in the O₃ concentration indirectly increased the OH+NO₂ reaction rate for nitrate formation. This effect was greater in inland cities, where the increase in the O₃ concentration was greater than that in coastal cities, leading to a greater contribution of the OH+NO₂ pathway. Thus, the OH and O₃ concentrations affected the formation rate and contribution of the OH+NO₂ pathway to nitrate formation due to the excess NO₂ concentration, although the NO₂ emissions decreased. Therefore, inland cities should prioritize controlling the emissions of VOCs and NO_x. These pollutants could generate OH radicals through photochemical reactions, which would subsequently dominate nitrate formation via the OH+NO₂ pathway. By reducing precursor emissions (VOCs and NO_x), the production of OH radicals could be suppressed, thereby slowing the reaction kinetics and mitigating atmospheric nitrate accumulation.

3.5.3 Factors influencing hetN₂O₅ reaction rates

The heterogeneous N₂O₅ reaction was considered a major pathway for nitrate formation influenced by many factors, including the N₂O₅ concentration, aerosol surface area, and uptake coefficient of N₂O₅ (γN₂O₅). The value of γN₂O₅ varied widely, ranging from an extremely low value of 10⁻⁴ to a relatively high value of 0.1 (Wang et al., 2023). Experimental results indicated that environmental factors such as temperature and humidity; the composition of aerosols such as ${\mathrm{NO}}_{\mathrm{3}}^{-}$, Cl⁻, and ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$; the liquid water content in aerosols; the presence of organic compounds; and the mixing state of aerosols were critical factors influencing γN₂O₅ (Folkers et al., 2003; Mitroo et al., 2019; Thornton and Abbatt, 2005; Wahner et al., 1998; Wang et al., 2020). These factors collectively influenced the reaction rate and yield of heterogeneous N₂O₅ reactions.

With increasing concentrations of atmospheric O₃ (Fig. 13f), the generation rate of N₂O₅ significantly increased. Consequently, the average N₂O₅ concentration in the NCP was 35.7 ppt in 2013 and increased to 53.2 ppt in 2018 (Fig. 12b). However, the N₂O₅ concentrations were generally higher in ocean areas and western mountainous regions than in central inland areas. We found that this spatial distribution of N₂O₅ was caused by NO. In urban areas, especially during severe haze nighttime periods, the rapid titration of NO with O₃ led to near-zero concentrations of surface O₃ (Zang et al., 2022), thereby inhibiting the in situ generation of NO₃ radicals and N₂O₅ (Zhao et al., 2023). As shown in Fig. 11, the NO concentration exhibited a spatial distribution inversely correlated with that of N₂O₅, where higher NO concentrations corresponded to lower N₂O₅ levels and vice versa. This inverse relationship explains the low N₂O₅ concentrations in inland cities. Additionally, observations have showed that the lifetime of NO₃ radicals in marine atmospheres could reach 30 min because of the low NO_x mixing ratio in marine air masses, which was longer than 1 min in inland atmospheres (Crowley et al., 2011) and contributed to the high N₂O₅ concentrations in coastal cities. Observations in Hebei Province revealed low N₂O₅ concentrations in summer, with lifetimes ranging from 0.1 to 10 min (Tham et al., 2018). Therefore, the spatial differences in hetN₂O₅ reaction rates depended on NO and the lifetime of NO₃ radicals.

Moreover, we speculated that the decreased rate of heterogeneous N₂O₅ was influenced not only by N₂O₅ but also by other factors, such as RH and the surface area of aerosol particles, due to the interannual variations in RH and PM. Compared with 2013, the RH in the NCP decreased by 3.1 % in 2018 (Fig. S4). The sensitivity experiments revealed that both the N₂O₅ reaction rate and its contribution to total nitrates remained essentially unchanged when the humidity decreased by 10 %. This finding indicated that humidity changes were not the main factor influencing the decrease in the reaction rate of hetN₂O₅ in the NCP. The total dry surface area of the particles in the CMAQ model was calculated as the sum of the modal dry surface areas, corresponding to the particles in the Aitken, the accumulation levels, and coarse modes (Bergin et al., 2022). During the winter of 2013 in the NCP, the average surface reaction factor (SRF) value was 689 µm² cm⁻³, which decreased to 425 µm² cm⁻³ during the winter of 2018, indicating a reduction of 38.3 %. This significant reduction in the aerosol surface area likely contributed to the observed decrease in the reaction rate of heterogeneous N₂O₅. Air quality monitoring data revealed that PM₁₀ decreased by 20 % and PM_2.5 decreased by 28 % in China from 2014 to 2018 (Fan et al., 2020a). Therefore, we considered that the reduction in the PM_2.5 concentration was a major factor contributing to the decrease in the hetN₂O₅ generation rate.

In brief, hetN₂O₅ reaction rates were affected by three key factors: N₂O₅ concentrations, aerosol surface area, and the spatial distribution of NO. The decline in hetN₂O₅ reaction rates across the NCP between 2013 and 2018 was driven primarily by a substantial reduction in aerosol surface area, particularly due to decreased PM_2.5 concentrations. Furthermore, regional variations in NO concentrations and NO₃ radical lifetimes significantly influenced N₂O₅ concentrations, resulting in distinct spatial patterns. To further mitigate nitrate pollution, coastal cities should prioritize reducing particulate matter concentrations, which would reduce hetN₂O₅ reaction rates and enhance air quality.

3.6 Impact of reducing NO_x and VOC emissions on the nitrate concentration

To reduce the nitrate content in PM_2.5, effective control strategies must target key factors, such as OH radicals, O₃, and the photochemical reactions of NO_x and VOCs. However, a previous study revealed that even with a 30 % reduction in VOC and NO_x emissions, the winter nitrate concentration in PM_2.5 in the NCP decreased by only 8.6 % (Fu et al., 2020). NH₃ played a crucial role in the gas-to-particle conversion process of HNO₃ in the NCP region. Therefore, we designed single and combined pollutant reduction strategies for NH₃, NO_x, and VOCs to examine the effects of the emission of these species on nitrate in PM_2.5.

The effectiveness of emission reduction scenarios notably differed between Beijing and Qingdao (Fig. 14). In Beijing, single-pollutant reduction strategies yielded outcomes dependent on target factors. The city showed limited sensitivity to NO_x reduction, with nitrate concentrations initially increasing before showing a 4.5 % decrease in nitrate concentration until reaching a 60 % NO_x reduction threshold. A 60 % reduction in VOC emissions led to a 10.5 % decrease in nitrate concentration, whereas NH₃ reduction led to a 38.4 % decrease at the same level. In contrast, Qingdao demonstrated relatively strong responsiveness to NO_x controls, with a 60 % reduction resulting in a 32.4 % decrease in nitrate concentration. However, Qingdao showed relatively weak sensitivity to VOC controls, with a 60 % reduction in VOC emissions resulting in an 11.1 % decrease in nitrate concentration. Qingdao exhibited the strongest response to NH₃ reduction measures, where a 60 % decrease in NH₃ emissions led to a substantial 51.7 % reduction in nitrate concentration.

Figure 14 Relative changes in the PM_2.5 nitrate contribution in response to the emission reduction scenarios in (a) Beijing and (b) Qingdao, compared with the winter 2018 baseline conditions

The comprehensive reduction scenarios revealed more pronounced intercity differences. Beijing showed a relatively weak response to combined NO_x and VOC reductions, with only a 13.0 % decrease in nitrate concentration at a 60 % reduction in both pollutants. Qingdao, however, achieved a significantly greater reduction of 34.7 % under the same conditions. The most effective strategy involved simultaneous reductions in NH₃, VOCs, and NO_x, with Qingdao achieving a 60.0 % reduction and Beijing achieving a 44.2 % reduction at 60 % emission cuts. These findings demonstrated that effective nitrate control measures would require the coordinated reduction of multiple precursors. While multipollutant strategies yielded optimal results, NH₃ reductions were determined to be the most effective single-pollutant control measure for both cities.

4 Conclusion

In this study, nitrate formation mechanisms and influencing factors in seven cities across the NCP during the winters of 2013 and 2018 were investigated. Despite inherent uncertainties, the model simulations and isotopic measurements showed remarkable agreement, demonstrating the reliability of the modeling approach in elucidating nitrate formation mechanisms. Our findings revealed that nitrate formation was significantly influenced by NH₃, NO₂, OH radicals, O₃, NO, and N₂O₅, with distinct differences between regions. In inland cities, the OH+NO₂ reaction pathway dominated, contributing 63.7 %–85.6 % of nitrate formation. This phenomenon was largely driven by the increased concentrations of OH radicals and O₃, resulting in a 7.6 % greater contribution of this pathway in inland cities in 2018 than in 2013. In contrast, coastal cities presented a greater contribution of the hetN₂O₅ pathway (37.0 % to 45.7 %) because of the higher N₂O₅ concentrations and longer NO₃ radical lifetimes. High NO concentrations in inland areas facilitated O₃ titration, inhibiting N₂O₅ formation and further differentiating nitrate formation processes between the two regions.

Our emission reduction experiments revealed distinct regional responses and highlighted the need for tailored, multipollutant strategies to effectively control nitrate pollution. Among the single-pollutant control strategies, NH₃ reduction demonstrated the highest effectiveness, with a 60 % emission reduction leading to a 38.4 % decrease in nitrate concentrations in Beijing and a substantial 51.7 % reduction in Qingdao. The combined emission control strategy, which involved simultaneous reductions in multiple pollutants, was more effective than NH₃ reduction alone, achieving a 44.2 % nitrate reduction in Beijing and a 60.0 % decrease in Qingdao. This study advanced our understanding of regional variations in atmospheric nitrogen chemistry and provided essential insights for developing targeted measures to mitigate nitrate pollution and improve air quality in northern China. To further advance research in this field, it is recommended that researchers in the future should focus on developing robust methodologies to comprehensively account for the influence of atmospheric transport on nitrate formation.