Fertilization-driven pulses of atmospheric nitrogen dioxide complicate air pollution in early spring over the North China Plain

Atmospheric nitrogen dioxide (NO₂) has shown periodic conspicuous pulses in the tropospheric column in March over the North China Plain during the past two decades. However, these repetitive pulses have never been reported, and their underlying causes remain unclear. Here, we present robust evidence to demonstrate that agricultural fertilization drives the early spring NO₂ column increase. The fertilization-driven soil NO_x (= NO+NO₂) emissions, comparable to anthropogenic sources, exert complicated influences on regional air quality. They significantly reduce nocturnal and diurnal O₃ concentrations in agricultural areas in early spring, distinct from the scenarios in summer, but increase fine particulate matter (PM_2.5) concentrations via strongly enhancing nitrate aerosol formation. The impact also extends to urban areas, approximately half that of agricultural areas. These findings have increasing implications for coordinated control of PM_2.5 and O₃ under global warming. We thus suggest that reducing NO_x emissions in croplands is essential to achieve better air quality in agricultural countries and regions.

1 Introduction

Nitrogen oxide (NO_x = nitric oxide (NO) + nitrogen dioxide (NO₂)) is a major air pollutant in the troposphere and a key precursor to ozone (O₃) and fine particulate matter (PM_2.5) due to its photochemical properties (Seinfeld and Pandis, 2006; Zhang et al., 2015). It is also a short-lived climate forcer regulated by both China and the United States (IPCC, 2023). Understanding the NO_x budget is crucial for addressing these issues. Globally, atmospheric NO_x is mainly produced by fossil fuel combustion (Crippa et al., 2023; Janssens-Maenhout et al., 2015; Yan et al., 2005), with smaller contributions from wildfires and lightning (Bauwens et al., 2020; Murray et al., 2012). Additionally, soil generates substantial NO_x through nitrification and denitrification processes (Bouwman et al., 2002; Cárdenas et al., 1993; Davidson, 1992; Yan et al., 2005), particularly after the rewetting of dry soils (Galbally and Roy, 1978; Huber et al., 2020; Yienger and Levy, 1995). On a regional scale, soil NO_x emissions may even exceed those from fossil fuel sources in summer (Almaraz et al., 2018; Sha et al., 2021). Model- and satellite-based studies estimate that global annual soil NO emissions, with the largest contributor being cultivated croplands, range from 9 to 27 Tg N (Hudman et al., 2012; Steinkamp and Lawrence, 2011; Vinken et al., 2014; Yan et al., 2005), accounting for about 15 % of total NO_x emissions (Hudman et al., 2012). This wide range is due to the complex response of soil NO_x emissions to driving factors like fertilization, temperature, and soil moisture (Huber et al., 2020; Oikawa et al., 2015), making accurate estimation challenging.

The emission rates from fertilized croplands are 1 to 2 orders of magnitude higher than those from nearby grasslands and forest soils (Almaraz et al., 2018; Anderson and Levine, 1987; Guo et al., 2020; Yienger and Levy, 1995). Recent studies show significant NO_x emissions from croplands post-fertilization, exceeding pre-fertilization rates by an order of magnitude (Almaraz et al., 2018; Hickman et al., 2017; Laville et al., 2011; Liu et al., 2005; Oikawa et al., 2015; Zhao et al., 2015). Despite this robust evidence of strong NO_x emissions from agricultural fertilization, the lack of extensive in situ measurements hinders accurate estimation of these emissions and their environmental impacts. Additionally, the effect of agricultural fertilization on air quality has not received sufficient global attention, although some pioneering studies have pointed out the implications for air quality since the 1990s (Davidson et al., 1998; Hall et al., 1996). In recent years, studies have reported that agricultural soil emissions significantly increase atmospheric NO_x levels (Almaraz et al., 2018; Hickman et al., 2017; Huang et al., 2018; Oikawa et al., 2015) and enhance O₃ formation in summer in California (Oikawa et al., 2015) or during the growing season of crops in sub-Saharan Africa (Hickman et al., 2017; Huang et al., 2018).

The North China Plain (NCP) is one of the major grain-producing regions in China. Winter wheat–maize double cropping is a typical rotation system mainly practiced in this region (Liu et al., 2003; Zhu et al., 1994). China has been the world's largest consumer of N fertilizer since 2000 (Liu et al., 2013), with annual usage peaking at approximately 31.2 Tg N in 2014 (Yu et al., 2022). About half of this fertilizer is lost to the environment (Liu et al., 2013), indicating a significant potential source for NO_x emissions from China's croplands. The agricultural management in the NCP has been known for incorporating high fertilization rates with excessive N fertilization (Sun et al., 2022; Vitousek et al., 2009; Zhao et al., 2006). Thus, this region is the largest consumer of agricultural fertilizer N in China (Yu et al., 2022) and has shown substantial soil NO_x emissions (Liu et al., 2010; Tang et al., 2020; Zhang et al., 2011). The emissions significantly increase ambient NO_x levels and enhance O₃ formation in summer (Huang et al., 2023; Lu et al., 2021; Wang et al., 2022a). These concerns typically focus on the warm season when higher temperatures favor NO_x emissions from soils. However, frequent agricultural activities and N-fertilizer use also occur during transitional seasons, and how periodic agricultural fertilization affects soil NO_x emission and regional air quality remains unclear.

In this study, we present a pulse of atmospheric NO₂ column in early spring during the past two decades over the NCP. However, this phenomenon has not been previously reported in this region. Combining agricultural fertilization records, surface NO₂ and NH₃ observations, long-term satellite observations of NO₂ and NH₃, and a flexible scheme of soil NO_x emission, we successfully explain the underlying cause for the NO₂ column peaks using a regional atmospheric transport model online coupled with chemistry and further assess the impacts of the pulsing NO_x emission on regional air quality.

Figure 1 Domain overview. (a) Geographic location of the NCP, which is predominantly characterized by plains at an elevation of less than 100 m and is known for being a major agricultural zone. (b) Extensive cultivated croplands are distributed in the NCP, marked by the orange shading, while urban areas are marked by red shading. The graphic markers denote locations of field observation sites, among which the blue circle represents the Fengqiu cropland ecological station, Chinese Academy of Sciences, with a long-term record of agricultural fertilization, and the triangle represents the rural Xianghe station with ambient NH₃ measurements. The agricultural areas in orange within the NCP are defined as croplands at an altitude of less than 100 m, and the urban areas in red are defined as built-up areas within the NCP. Publisher's remark: please note that the above figure contains disputed territories.

Table 1 Model configuration for the simulation domain, meteorological schemes, initial and boundary conditions, and emission inventories.

2 Materials and methods

2.1 Model and configurations

The Weather Research and Forecasting model fully coupled with atmospheric chemistry (WRF-Chem, version 3.6.1) that we used is a modified model by Li et al. (2011a, b, 2012, 2010) and Feng et al. (2021), in which we implement the BDSNP mechanism by Hudman et al. (2012) to calculate soil NO_x emissions related to agricultural fertilization and their influences on regional air quality in the NCP. The model is configured with grid spacing of 6 km × 6 km (240 × 280 grid cells) with the center at 38° N, 116° E (Fig. 1). Thirty-five vertical levels are employed in the stretched vertical grid with spacing ranging from 50 m near the surface to 500 m at 2.5 km and 1 km above 14 km. Meteorological initial and boundary conditions use the National Centers for Environmental Prediction (NCEP) FNL 1° × 1° analysis data, and the chemical initial and boundary conditions are interpolated from the CAM-Chem 6 h output (Buchholz et al., 2019; Emmons et al., 2020). The non-soil emission inventory is developed by Zhang et al. (2009), and the biogenic emissions are calculated online using the Model of Emissions of Gases and Aerosol from Nature (MEGAN) model (Guenther et al., 2006). Specifically, monthly ammonia (NH₃) emissions are incorporated from a high-resolution NH₃ emission inventory developed by Huang et al. (2012), which includes emissions from fertilizer application, livestock, and other sources. The model spin-up time is 2 d (Table 1).

2.2 Soil NO_X emission scheme

We implement a soil NO_x emission scheme, the Berkeley–Dalhousie Soil NO Parameterization (BDSNP) by Hudman et al. (2012), into the WRF-Chem model. The scheme comprehensively considers various factors, including available soil nitrogen content (N_avail, ng N m⁻²) from fertilizer application and nitrogen deposition, in which the soil NO_x emission (E_soil, $\mathrm{ng}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$) is a function of N_avail, climate, and edaphic conditions:

$$\begin{array}{}\text{(1)}& E_{\text{soil}}=A_{\text{biome}}^{\prime }\left(N_{\text{avail}}\right)\times f\left(T\right)\times g\left(\mathit{\theta }\right)\times P\left(l_{\text{dry}}\right)\end{array}$$

where N_avail is available soil nitrogen mass, and $A_{\text{biome}}^{\prime }$ ($\mathrm{ng}\mathrm{N}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$) represents the biome-dependent emission factor. f(T) (dimensionless) and g(θ) (dimensionless) are parameters regulated by soil temperature and moisture, respectively. P(l_dry) (dimensionless) denotes the pulsed soil emission from wetting of dry soils. The product of f(T) and g(θ) is calculated as follows:

$$\begin{array}{}\text{(2)}& f\left(T\right)\times g\left(\mathit{\theta }\right)=e^{\mathrm{0.103}T}\times a\mathit{\theta }{e}^{-b{\mathit{\theta }}^{\mathrm{2}}}\end{array}$$

where T ($\mathrm{0}\le T\le$ 30 °C) is soil temperature and θ ($\mathrm{0}\le \mathit{\theta }\le \mathrm{1}$, dimensionless) is water-filled pore space, defined as the ratio of the volumetric soil moisture content to the porosity. According to laboratory and field measurements (Hudman et al., 2012), the constants a and b are determined so that g(θ) maximizes when θ = 0.2 for arid soils and θ = 0.3 elsewhere.

The pulsing term P(l_dry), following Yan et al. (2005), describes the magnitude of the peak flux relative to the pre-wetting flux, which is parameterized as:

$$\begin{array}{}\text{(3)}& P\left(l_{\text{dry}}\right)=\left[\mathrm{13.01}\ln \right(l_{\text{dry}})-\mathrm{53.6}]\times e^{-ct}\end{array}$$

where l_dry (hours) represents the length of the antecedent dry period, and c (c = 0.068 h⁻¹) is a constant rate denoting the rise/fall time of the pulse. Fertilizer application data are interpolated from the global gridded chemical fertilizer and manure application inventory at 0.5° × 0.5° (Potter et al., 2010; Yan et al., 2005). The chemical and manure fertilizers are obtained from the International Fertilizer Association (IFA) and the Food and Agriculture Organization of the United Nations (FAO). The Chinese chemical fertilizer application (straight N application) from IFA is about 19.6 Tg N a⁻¹ for 2000, quite close to the amount of 19.9 Tg N a⁻¹ for 2020 from the China Statistical Yearbook (https://www.stats.gov.cn/sj/ndsj/2021/indexch.htm, last access: 10 September 2025). More details of the scheme are found in related studies elsewhere (Hudman et al., 2012; Lu et al., 2021). It should be noted that we use the default BDSNP fertilizer and global emission assumptions, and these assumptions (primarily the $\mathrm{75}/\mathrm{25}$ split and not applying the topdressing in a single application) may not accurately reflect the fertilizer applications in China. The BDSNP tuning, which was done for the year 2000 agricultural emissions estimated by Stehfest and Bouwman (2006), is close to recent estimates of global emissions, possibly on the high side (Gong et al., 2025; Wang et al., 2022b).

2.3 Emission inventories

We employ two emission inventories in this study: the Hemispheric Transport of Air Pollution Version 3 (HTAP v3, 2005–2018) emission inventory, which includes soil and non-soil emissions (Li et al., 2017a), and the Multi-resolution Emission Inventory for China (MEIC v1.3, 2007–2018), which has no soil emissions (Li et al., 2017a). In the HTAP inventory, the non-soil emission inventory includes energy, industry, ground transport, residential, waste, shipping, and aviation sources, with a spatial resolution of 0.1° × 0.1° and a temporal resolution of one month. Agricultural emissions are involved in the latest HTAP v3 inventory, which includes soil NO_x emissions (Crippa et al., 2023). Nevertheless, the soil emissions in this inventory are calculated using the traditional “bottom-up” method (Kurokawa and Ohara, 2020), rather than estimated by a process-based emission module. The monthly MEIC emission inventory, with a spatial resolution of 0.25° × 0.25°, is incorporated in parallel with the HTAP emission inventory. Here, we focus on NO_x and NH₃ emissions from croplands with fertilization and adopt the NH₃ emission inventory by Huang et al. (2012) because it explicitly distinguishes NH₃ produced by agricultural fertilization from other NH₃ sources.

2.4 Air pollutant measurements

Satellite-derived tropospheric NO₂ columns are from Ozone Monitoring Instrument (OMI) hosted by the Aura satellite, which was launched by the National Aeronautics and Space Administration (NASA). The Level-3 product, where pixel-level data of good quality are binned and “averaged” into 0.25° × 0.25° grids, was retrieved and analyzed in the present study. The satellite operates in a sun-synchronous polar orbit and has a local overpass time of around 13:45 LT (local time) in North China. The dataset is for all atmospheric conditions and for sky conditions with cloud fraction less than 30 % (https://cmr.earthdata.nasa.gov/search/concepts/C1266136111-GES_DISC.html, last access: 10 September 2025). The dataset has a spatial resolution of 13 km × 24 km, with a temporal coverage of 2005–2022 (Lamsal et al., 2021). Note that the number of pixels included in NO₂ retrievals changes over time because of the increase in the number of pixels affected by the row anomaly issue, making the data unsuitable for trend analysis and possibly introducing uncertainty in seasonal averages. The Level-2 product of NH₃ columns is employed, which is provided by the Space Administration and the Infrared Atmospheric Sounding Interferometer (IASI) hosted on the MetOp series of satellites. The satellite also operates in a sun-synchronous polar orbit and has a local overpass time of around 09:30 am and 09:30 pm in North China (twice a day). We construct a 0.25° × 0.25° mesh grid and calculate the average of the NH₃ columns from IASI within each grid cell (Clarisse et al., 2023). Low-quality satellite data are filtered out due to the interference of clouds. To cover all the domain (Fig. 1), the data used in this study are averaged into 7 d mean datasets of NO₂ and NH₃ columns with a non-overlapping 7 d window during 2007–2021. The data are interpolated into the model grids using bilinear interpolation.

Figure 2 NO₂ column pulses in March, NH₃ column variation, and NO_x emissions from fossil fuel and soil sources over the NCP. (a) Long-term variation of the  7 d mean tropospheric NO₂ column observed by OMI during 2005–2022 (green) and NH₃ column retrieved from IASI during 2007–2022 (pink). Intersections of the gray bars and the green lines denote that a sub-peak of the NO₂ column occurred in each March, and the short bars represent the timing record for agricultural fertilization at the Fengqiu station in the NCP, with the red ones indicating the fertilization period in early spring. (b) Monthly mean NO_x emission rates with ± 1σ standard deviation (SD) in two sets of anthropogenic emission inventories, the HTAP v3 (2005–2018, orange) and MEIC v1.3 (2008–2017, red). (c) Same as (b), but for NO_X emission rates from soils in the HTAP v3 inventory (2005–2018).

Ambient surface NO₂, O₃, and PM_2.5 mass concentrations at 141 sites in the NCP are obtained from the China National Environmental Monitoring Centre (CNEMC, Fig. S1 in the Supplement). These in situ measurements are performed by the Thermo Scientific™ ambient particulate monitor and gas analyzers, in which NO₂ and O₃ are measured by the Model 42i Chemiluminescence NO–NO₂–NO_x Analyzer and the Model 49i UV Photometric Ozone Analyzer, respectively. PM_2.5 is measured by the Model 5030 Synchronized Hybrid Ambient, Real-time Particulate (SHARP) Monitor, which uses proprietary digital filtering to continuously calibrate mass to obtain an accurate, precise, and real-time mass concentration. The sampling time is 1 min for these monitoring devices. Agricultural NH₃ concentration is monitored by a Picarro analyzer based on the principle of cavity ring-down spectroscopy (CRDS) at the rural Xianghe station (Fig. 1), with a sampling frequency of 1 Hz. Hourly data are derived by averaging the high-frequency measurements.

3 Results and discussion

3.1 Satellite-retrieved NO₂ column pulses

During the past two decades, 7 d mean tropospheric column of NO₂ measured by the Ozone Monitoring Instrument (OMI) in the NCP exhibits a significant temporal variation, with the magnitude varying from less than 100 µmol m⁻² to more than 680 µmol m⁻² (Fig. 2a). The annual cycle is highly prominent and its seasonal variation is remarkable, with significantly higher levels in cold seasons than those in warm seasons. Throughout the year, the pattern of the NO₂ column looks like a rhinoceros horn, which is characterized by a major peak in winter and multiple noticeable sub-peaks in other seasons. These sub-peaks often occur at fixed times, such as in March, June, and October, the most noticeable of which is in March, with the highest magnitude (Figs. 2a and S1). We examine the monthly variation in anthropogenic NO_x emission rates over the NCP in global and regional emission inventories and find that the monthly variation is more evident in the regional emission inventory, with a significantly higher emission than that in the global emission inventory. Nevertheless, neither of them reveals any discernible sub-peaks of NO_x emission rates from fossil fuel combustion to coincide with the sub-peaks of the NO₂ column (Fig. 2b).

As for soil NO_x emissions, they are absent in the regional emission inventory, while in the global emission inventory, soil NO_x emissions fluctuate slightly on a monthly scale, far less than those from fossil fuel combustion, constituting less than 2 % of the total (Fig. 2c). Similar to NO_x emissions from fossil fuel combustion, there are no evident sub-peaks of soil NO_x emission to keep pace with the atmospheric NO₂ column. Soil NO_x emissions are at an even lower level in March, significantly less than emissions in adjacent months. Therefore, the known emission inventory fails to explain the occurrence of these sub-peaks. On the other hand, we compute a pollution accumulation index (PAI, Sect. S1 in the Supplement), the product of boundary layer height and wind speed, to semi-quantitatively assess the influence of atmospheric dispersion conditions on the NO₂ column. The NO₂ column seems to be somewhat dependent on the PAI, yet the noticeable discrepancies between the timing of the sub-peaks and the PAI are insufficient to account for the occurrence of each sub-peak (Fig. S3). Additionally, it is observed that the daily soil temperature was consistently higher than 0 °C during March 2020 and in the 10 days before (Fig. S4). Therefore, the sub-peak of the NO₂ column is not expected to originate from soil thaw either.

3.2 Linkage between NO₂ column pulses and agricultural fertilization

What causes these regular NO₂ sub-peaks that have occurred over the NCP during the past two decades? Measurements on ammonia (NH₃) columns also present similar pulses to those of NO₂ column during the same period, in spite of some differences in the long-term trend (Fig. 2a). Although the main peaks of NH₃ column occur in June, the sub-peaks of NH₃ column in March may provide favorable evidence that these NO₂ column sub-peaks are connected to agricultural activities because atmospheric NH₃ largely originates from fertilizer application in agriculture (Crippa et al., 2023; Huang et al., 2012; Li et al., 2017a). We note that the two time series do not have similar seasonal dynamics, which may raise some questions about whether the March NO₂ peaks represent fertilizer pulses. This discrepancy could be due to several factors: (1) contributions from fossil fuel-derived NO_x, (2) elevated background soil NO_x emissions during the “spring thaw” period, (3) differences in fertilizer type, as fertilizers vary in their potential for ammonia volatilization, and (4) variation in application methods, e.g., banding or deep soil placement vs. broadcasting, especially when the latter is done without incorporation.

Another key evidence is that the occurrence of each sub-peak of the NO₂ column is highly consistent with the record of agricultural fertilization at the Fengqiu cropland ecological station in the NCP during the past decades (Figs. 1b, 2a, and S1). The wheat–maize double-cropping system is predominant in the NCP, where the agricultural activities are strongly dependent on the lunar calendar. For winter wheat, the planting date ranges from early to mid-October (after maize harvest). Fertilization is generally divided into three stages: (1) pre-planting during late September–early October; (2) jointing stage during mid-March–early April; (3) grain filling during late April for high-yield fields. The planting date of summer maize ranges from early to mid-June (after wheat harvest), and the stages of fertilization include: (1) at planting during early June; (2) V6–V8 stage during early July; (3) tasseling stage during late July for high-yield fields. The agricultural fertilization is closely associated with three solar terms, i.e., Waking of Insects in March (the 3rd solar term), Grain in Beard in June (the 9th solar term), and Cold Dew in October (the 17th solar term). During these pulses of the NO₂ column, we found that the pulse in March is more pronounced than those in June and October, because March is the season for large-scale cultivation in the NCP, accompanied by more land preparation and fertilization. Therefore, we hypothesize that the NO₂ column pulse in March is possibly caused by fertilized croplands that accelerate NO_x emissions from agricultural soils. Field campaigns have measured a high NO emission rate of 266.3 $\mathrm{g}{\mathrm{km}}^{-\mathrm{2}}{\mathrm{h}}^{-\mathrm{1}}$ in croplands after fertilization and irrigation in autumn in eastern China (Tang et al., 2020; Tian et al., 2020) and also in other regions (Hickman et al., 2017; Huang et al., 2018; Huber et al., 2020), suggesting that agricultural fertilization is likely a significant source of atmospheric NO_x in major agricultural countries like China.

Figure 3 NO_x emissions from agricultural fertilization and resultant NO₂ column during February–April 2020 over the NCP. (a) Calculated NO_x emission rate from croplands with N-fertilizer application in the model. The black curve represents daily variation in the NO_x emission rate around the fertilization, and the blue, red, and green dashed lines correspond to the trends of NO_x emission rates in croplands during the pre-fertilization, fertilization, and post-fertilization periods, respectively. (b) Observed and simulated NO₂ column. The gray histogram represents the NO₂ column observed by satellite (OMI). The green and pink lines represent the simulated NO₂ column with and without soil NO_x emissions, and the black line shows the difference between them. The model well replicates the rapid increase in the observed NO₂ column by considering soil NO_x emissions from agricultural fertilization.

3.3 Soil NO_x emission mechanism

To examine the role of soil NO_x emissions from agricultural fertilization in the pulses of the atmospheric NO₂ column, we introduce a flexible soil NO_x emission module and NH₃ emission into the WRF-Chem model and perform two simulation experiments that include and exclude soil emissions, respectively (Table 1). Noticeably, there is an acceleration in the release of soil NO_x, and the daily mean emission rate increases from 155.6 to 438.3 $\mathrm{g}{\mathrm{km}}^{-\mathrm{2}}{\mathrm{h}}^{-\mathrm{1}}$ during the simulation period (Fig. 3a). In particular, the NO_x emission rate during the post-fertilization phase is significantly higher than that of other phases, consistent with the accelerated soil NO_x release observed in agricultural areas in California after fertilization (Oikawa et al., 2015). On average, the simulated NO_x emission rate in March is 312.9 $\mathrm{g}{\mathrm{km}}^{-\mathrm{2}}{\mathrm{h}}^{-\mathrm{1}}$, between the measured 113.6 $\mathrm{g}{\mathrm{km}}^{-\mathrm{2}}{\mathrm{h}}^{-\mathrm{1}}$ in November in eastern China (Tang et al., 2020) and 988.2 $\mathrm{g}{\mathrm{km}}^{-\mathrm{2}}{\mathrm{h}}^{-\mathrm{1}}$ in September in California (Oikawa et al., 2015), suggesting the rationality of the soil NO_x emission mechanism in the model. We should also note that the BDSNP scheme would be heavily biased toward having a fertilizer-induced emission peak only at planting/green-up.

We evaluate the model performance against the satellite-derived NO₂ column. Consequently, the modified model perfectly replicates the sudden increase in NO₂ column linked to agricultural fertilization, while the conventional WRF-Chem model fails to capture the observed NO₂ column pulse in March due to the lack of the adopted soil NO_x emission mechanism (Fig. 3b). For example, when soil NO_x emission caused by agricultural fertilization is considered, the simulated NO₂ column rapidly increases to the peak in March, matching well with the satellite observation. However, without the contribution of agricultural fertilization, the NO₂ column seems to exhibit a weak upward trend but not a significant one. Comparing these two scenarios, a substantial NO_x emission from N-fertilizer input in croplands leads to an increase in the NO₂ column by 1 to 1.5 times.

Figure 4 Simulated vs. measured surface pollutants averaged over the monitoring sites of the NCP (Fig. S1) during February–April 2020. (a–c) Temporal variations in surface NO₂ (a), O₃ (b), and PM_2.5 (c) mass concentrations. The blue curves denote the model calculations, and the black dots denote in situ measurements. Model biases are shown in the central upper position of each figure.

We also validate the modified model performance on temporal variations of routine surface pollutant measurements (NO₂, O₃, and PM_2.5) associated with NO_x emissions at the CNEMC sites throughout the simulation period (Fig. 4). Although there are some discrepancies between the simulations and observations, e.g., overestimates occur in mid-February for NO₂ and PM_2.5 levels, the model generally reproduces hourly variations in each pollutant reasonably well. The indices of agreement (IOAs, Sect. S2) between the simulated and observed near-surface concentrations of NO₂, O₃, and PM_2.5 are 0.90, 0.91, and 0.88, respectively, and the normalized mean biases (NMBs) for these pollutants are within 10 %.

We still cannot ignore the discrepancies between the model results and observations. These biases may largely originate from the soil NO_x emission mechanism. The fertilization dates in the BDSNP mechanism are determined by the beginning and end of the growing season, which are derived from the MODIS Land Cover Dynamics product (MCD12Q2) averaged over the years from 2001 to 2004 (Hudman et al., 2012). This may be quite different from practices in 2020, the year we simulated in this study. We use the default assumption in the mechanism that 75 % of fertilizer is added on the green-up day, with the remaining 25 % applied constantly throughout the rest of the season (Hudman et al., 2012). We should note that the $\mathrm{75}/\mathrm{25}$ treatment is an estimate for local farming practice in Mexico in 1994 (Matson et al., 1998), which probably introduces biases in other regions. Typical practice would place more fertilizer during topdressing and less at planting. For example, a $\mathrm{20}/\mathrm{80}$ split was used by GGCMI, the largest global crop modeling exercise (Jägermeyr et al., 2021), and a $\mathrm{40}/\mathrm{60}$ split was used for maize production in Northeast China (Zheng et al., 2023). It is also worth noting that the BDSNP mechanism treats topdressing as a series of applications spread out over several weeks or months rather than as a single event, which could further influence the modeling results. Despite the uncertainties, all of these significant improvements in the modified model we used suggest that soil NO_x emission from agricultural fertilization would exert a crucial influence on regional air quality.

Figure 5 Vertical profiles for impacts of soil emissions on gas pollutants during March 2020 over the NCP. (a)  Difference in NO₂ concentration with and without the influence of soil NO_x emission from agricultural fertilization at various heights in the near-surface layers. (b) Similar to (a), but for NH₃ concentration.

Figure 6 Contribution of soil emissions from agricultural fertilization to surface NO₂ and NH₃ during February through April 2020 over the NCP. (a, b) Change in surface NO₂ concentration with (green) and without (pink) soil NO_x emissions from agricultural fertilization, and the black line in (a) represents the observed surface NO₂ concentration. (c, d) Same as (a, b), but for NH₃. The error bar in (b) and (d) denotes ± 1σ. NO₂ observations are averaged over the 141 monitoring stations in the study area from the CNEMC network. NH₃ observations are from the rural Xianghe station (Fig. 1). According to in situ measurements of NO₂ and NH₃, the units for NO₂ and NH₃ concentrations are µg m⁻³ and ppb, respectively.

Furthermore, we examine the ability of the model to simulate the ground-level NO₂ mass concentration and NH₃ volume concentration when soil NO_x rapidly releases after fertilization. The reason is that the influences of soil emissions on atmospheric NO₂ and NH₃ concentrations are confined to the near-surface layers below 1 km, and the influences diminish rapidly as altitude increases (Fig. 5). This indicates that the impact of the soil emissions is primarily concentrated near the ground surface. With soil emissions included or not in the model, we compare the simulated NO₂ and NH₃ concentrations with near-surface observations (Fig. 6a and b). When there are no soil NO_x emissions from agricultural fertilization, the simulated NO₂ concentration is significantly lower than the observed concentration by 9.4 µg m⁻³ during February through April 2020. While considering these emissions, the mean bias (MB, Sect. S2) between the simulation and the observation decreases to 2.6 µg m⁻³, and the IOA also increases from 0.48 to 0.78. Similarly, the simulated NH₃ concentration is in good agreement with the observed concentration when the soil NH₃ emission related to agricultural fertilization is involved, e.g., the MB decreases from −12.0 to −4.4 ppb, and the IOA increases from 0.56 to 0.64 (Fig. 6c and d). It is important to acknowledge the limitation posed by the absence of direct comparisons with flux measurements of NO_x emissions from soils, due to the unavailability of such data. The simulated NO_x emission flux from the BDSNP scheme cannot be well examined, which may introduce uncertainties to the predicted emission rates and mixing ratios in the atmosphere.

Figure 7 Direct impacts of soil emissions from agricultural fertilization on (a) surface NO₂ and (b) NH₃ during March 2020 over the NCP. (a, b) Spatial distributions of changes in surface NO₂ and NH₃ concentrations due to fertilization-related soil emissions. According to in situ measurements of NO₂ and NH₃, the units for NO₂ and NH₃ concentrations are µg m⁻³ and ppb, respectively.

3.4 Significance of soil NO_x emissions from agricultural fertilization for air quality

We perform a model experiment that excludes the soil sources of NO_x and NH₃ in the study domain to examine the impacts of soil emissions on regional air quality. The model results are compared to the benchmark scenario with soil sources involved to examine these impacts. Agricultural fertilization directly leads to substantial increases in atmospheric NO_x and NH₃ concentrations. According to the spatial correlation between land use and NO₂ concentration, NO₂ concentrations increase by more than 15 µg m⁻³ over agricultural areas, with the maximal increments occurring in the densely cultivated southern region of the NCP, exceeding 40 µg m⁻³ (Figs. 1b and 7a). While in urban areas, the increase in NO₂ concentration is mostly below 10 µg m⁻³, significantly lower than that in agricultural areas. This indicates that the influence of local emissions originating from agricultural fertilization on air quality primarily concentrates in agricultural areas. Nonetheless, the influence extends to surrounding areas through atmospheric transport, leading to an inhomogeneous increase of NO₂ concentrations across the NCP. The spatial distribution of the increased NH₃ concentration is highly similar to that of the increased NO₂ concentration, but some differences exist in the southeast of the NCP. It should be noted that NH₃ emission in the model is from Huang et al. (2012), a separate monthly emission inventory. The emission rates of NH₃ in the southeast of the NCP are lower than those of NO from the BDSNP scheme. This indicates nonnegligible discrepancies in the derived emissions between these two approaches, which deserve more in-depth studies.

Figure 8 Complex impacts of agricultural fertilization on O₃, PM_2.5, and OH during March 2020 over the NCP. (a) Changes in mass concentrations of O₃, PM_2.5, aerosol constituents, i.e., nitrate, ammonium, sulfate, and secondary organics, and OH radical due to soil NO_x emission from agricultural fertilization in agricultural (green) and urban (pink) areas. The error bar denotes ± 1σ. (b) Same as (a), but for percentage changes.

A substantial amount of reactive nitrogen from agricultural fertilization suddenly enters the atmosphere and further affects air quality via photochemical reactions and aerosol chemical transformations profoundly (Seinfeld and Pandis, 2006; Wu et al., 2020). Our results reveal that the NO_x emission induced by N fertilization significantly suppresses the early spring O₃ production in the NCP, which varies remarkably with the land use, approximately twice as strong in agricultural areas as in urban areas. For instance, in agricultural areas, the emission in croplands reduces nocturnal and diurnal O₃ by 30.1 ± 6.5 µg m⁻³ (37.5 ± 8.1 %) and 15.0 ± 3.7 µg m⁻³ (18.7 ± 4.6 %), respectively, while in urban areas, the corresponding O₃ reductions are 15.6 ± 4.7 µg m⁻³ (15.6 ± 4.7 %) and 9.7 ± 3.2 µg m⁻³ (10.6 ± 3.4 %), respectively (Fig. 8). Based on the diurnal cycle of the change in O₃ concentrations (Δ[O₃]), we also find that the nighttime O₃ reduction is much higher than the daytime reduction (Fig. 9). The Δ[O₃] caused by agricultural fertilization is linearly and negatively correlated with the change in NO₂ concentration (Δ[NO₂]) (Fig. 10a–d), and the negative correlation is more pronounced at night (r < −0.99 and p < 0.001 for both the agricultural and the urban areas, Fig. 10a and b). This suggests that the O₃ concentration strongly depends on the change in NO_x levels in the NCP during early spring. Continuous agricultural NO_x (mainly NO) emissions inhibit the O₃ formation, which indicates that the NO_x abundance is excessive for O₃ formation over the NCP. Huang et al. (2023) similarly reported that variations in O₃ concentrations were inversely related to the changes in soil NO_x emissions in the NCP. Lu et al. (2021) also reported a NO_x-saturated O₃ formation regime in the NCP. On the other hand, a negative correlation between Δ[NO₂] and the change in daytime OH radical (Δ[OH]_day) suggests that the Δ[NO₂] also moderately regulates Δ[OH]_day (r = −0.50 for agricultural areas and r = −0.43 for urban areas, p < 0.001, Fig. 10e and f) through decreasing O₃ levels and reactions of NO₂ with OH radical. Both OH radical and O₃ are critical oxidants in the atmosphere, and the decrease caused by the excessive NO_x emission from agricultural fertilization weakens atmospheric oxidizing capacity (AOC) (Feng et al., 2021). The decreased AOC can further slow down the oxidation processes in homogeneous and heterogeneous reactions, which is unfavorable for the formation of secondary aerosols. We note that soil nitrous acid (HONO) emissions, which are not included in these modeling experiments, can also perturb atmospheric chemistry and the AOC (Feng et al., 2022; Tan et al., 2023) via the provision of NO and OH through photolysis. The emission rate of HONO from soil is much less than that of NO_x in the NCP (Tan et al., 2023), which increases daytime O₃ and OH concentrations slightly during summer (Feng et al., 2022; Tan et al., 2023). However, the influence in springtime still remains to be elucidated.

Figure 9 Secondary impact of soil NO_x emissions from agricultural fertilization on surface O₃ during March 2020 over the NCP. Diurnal cycles of changes in surface O₃ concentrations due to fertilization-related soil emissions over croplands (blue) and urban areas (orange) in the NCP. The blue and orange shadings show ± 1σ of the data.

Figure 10 Changes in surface NO₂ and related photochemical products. (a–f) Correlation between Δ[NO₂] and Δ[O₃] or Δ[OH] during March 2020 over the NCP. (a–d) Change in O₃ concentration is strongly dependent on the change in NO₂ concentration due to agricultural fertilization in both agricultural (a, c) and urban (b, d) areas, and the dependence is more pronounced at night, i.e., correlation coefficient r = −0.997 (r = −0.994) at night and r = −0.991 (r = −0.987) during the day in agricultural (urban) areas. (e, f) Change in daytime OH radical is also significantly influenced by the change in NO₂ concentration in both agricultural (e) and urban (f) areas.

Interestingly, these findings regarding the impacts of soil NO_x emission on O₃ formation in spring are different from previous studies that revealed that agricultural NO_x emissions enhance the O₃ formation in summer over the NCP (Tan et al., 2023; Wang et al., 2022a) and Northeast China (Shen et al., 2023) and in the Imperial Valley, California (Oikawa et al., 2015). Similar scenarios are also reported during the growing season of crops in sub-Saharan Africa (Hickman et al., 2017; Huang et al., 2018). This is largely attributed to the sensitivity of O₃ to its precursors. During early spring, a large amount of agricultural NO_x (mainly NO) emission causes a NO titration effect during daytime, decreasing O₃ concentrations, when O₃ chemistry is under the VOC-sensitive (NO_x-saturated) or the transitional regimes (Fig. S6) (Sillman, 1995). In contrast, O₃ formation chemistry in summer shifts from VOC-sensitive to NO_x-sensitive (Sha et al., 2021; Wang et al., 2022a). In this scenario, O₃ production is primarily controlled by NO_x emissions, meaning that the O₃ concentration increases with rising NO_x levels. This seasonal difference in O₃ sensitivity to its precursors highlights a seasonally dependent response of O₃ production to agricultural fertilization.

We also quantify the impact of agricultural fertilization on PM_2.5 concentrations. The NCP is characterized by an excess of NH₃, in which nitrate formation is highly sensitive to NO₂ concentration and AOC due to NO₂ oxidation to ${\mathrm{NO}}_{\mathrm{3}}^{-}$ via gas-phase and heterogeneous reactions (Feng et al., 2018; Fu et al., 2020; Liu et al., 2019; Wen et al., 2018). As atmospheric NO₂ and NH₃ concentrations rapidly increase due to emissions from fertilized croplands, nitrate aerosol in agricultural (urban) areas rises by 4.7 (4.0) µg m⁻³, corresponding to the increased percentage of 53.2 % (52.3 %), while ammonium aerosol rises by 1.3 (1.1) µg m⁻³ in agricultural (urban) areas, with an increased percentage of 27.7 % (29.4 %) (Fig. 8). However, sulfate aerosol shows a slight decrease in both agricultural and urban areas (Fig. 8a). The reason is that extra NO_x emissions from agricultural fertilization enhance nitrate formation but lower AOC, which hinders sulfate formation. Similar to sulfate aerosol, secondary organic aerosol (SOA) also has a slight reduction (Fig. 8a). The formation of SOA greatly depends on the AOC level, so decreased AOC due to NO_x emission from agricultural fertilization does not favor the conversion of organic precursors, such as VOCs and semi-volatile primary organic aerosols, into SOA.

Figure 11 Simulated changes in surface NO₂ and related aerosol-chemistry products during March 2020 over the NCP. (a) Comparison of NO₂ conversion to nitrate aerosol (${\mathrm{NO}}_{\mathrm{3}}^{-}$) formation between agricultural and urban areas, i.e., $\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\right]/\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{2}}\right]$ = 0.20 in urban areas and $\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\right]/\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{2}}\right]$ = 0.13 in agricultural areas, indicates that the change in nitrate in urban areas is more sensitive to the change in NO₂ concentration and that NO₂ conversion to ${\mathrm{NO}}_{\mathrm{3}}^{-}$ is more efficient. (b) NO₂ conversion to PM_2.5 formation is similar to (a), because nitrate aerosol is the most affected among the various aerosol constituents. The blue and orange colors correspond to the agricultural and urban areas, respectively.

In general, due to the NO_x emission from agricultural fertilization, PM_2.5 concentration increases by 5.1 (4.5) µg m⁻³ (Fig. 8a), corresponding to a percentage change of 12.9 % (13.9 %) over agricultural (urban) areas in the NCP (Fig. 8b). There is no significant difference in PM_2.5 increments between agricultural and urban areas. Nitrate aerosol is primarily responsible for the increased PM_2.5, accounting for 92.2 % in agricultural areas and 88.9 % in urban areas. Our results also indicate that changes in PM_2.5 and nitrate in urban areas are more sensitive to the change in NO₂ concentration. For instance, the ratios of nitrate change to NO₂ change ($\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\right]/\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{2}}\right]$ = 0.20) and PM_2.5 change to NO₂ change ($\mathrm{\Delta }\left[{\mathrm{PM}}_{\mathrm{2.5}}\right]/\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{2}}\right]$ = 0.24) in urban areas are both higher than those in agricultural areas ($\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{3}}^{-}\right]/\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{2}}\right]$ = 0.13 and $\mathrm{\Delta }\left[{\mathrm{PM}}_{\mathrm{2.5}}\right]/\mathrm{\Delta }\left[{\mathrm{NO}}_{\mathrm{2}}\right]$ = 0.15, Fig. 11a and b), indicating that the conversion of NO₂ to nitrate aerosol is more efficient in urban areas. Consequently, the increased percentages of PM_2.5 and ammonium aerosol in urban areas are higher than those in agricultural areas (Fig. 8b). Additionally, the ongoing stringent control measures on emission sources significantly reduce anthropogenic emissions in urban areas. Thus, the impact of agricultural fertilization on urban air quality is becoming more pronounced (Fig. S7). Since soil NO_x emission is sensitive to soil temperature, as global warming is ongoing, routine events like agricultural fertilization will continue to have amplified impact on air quality with the joint help of atmospheric dispersion/transport and chemical transformation processes (Bennetzen et al., 2016; Ma et al., 2022; Tubiello et al., 2013). These impacts are not confined to agricultural areas alone but extend to surrounding cities.

4 Conclusions and implications

The impact of soil NO_x emissions from agricultural fertilization on the atmospheric environment remains unclear worldwide (Guo et al., 2020; Huang et al., 2018; Sha et al., 2021; Shen et al., 2023). In particular, this issue has not yet received enough attention in China, where substantial N fertilizers are consumed year by year due to extensive agricultural cultivation areas (Sun et al., 2022; Vitousek et al., 2009; Zhao et al., 2006). Our results indicate that agricultural fertilization is highly responsible for the periodic pulse of the atmospheric NO₂ column in the NCP over the past two decades. A two-decade record of fertilization events at a research station and model results both provide evidence consistent with a cause-and-effect relationship. For example, the fertilization timing is found to match well with the occurrence of satellite-derived NO₂ column pulse in the region. Moreover, the model reasonably captures the regular sub-peak of the NO₂ column in March by introducing an independent module that specifically describes soil NO_x emissions from agricultural fertilization.

These additional NO_x emissions released by croplands directly lead to an elevated level of surface NO_x concentration. Consequently, the increased atmospheric NO_x concentration significantly inhibits O₃ production in early spring, distinct from the impacts in summer (Sha et al., 2021; Wang et al., 2022a), but enhances nitrate formation. For example, soil emissions linked to agricultural fertilization dramatically reduce nighttime O₃ concentrations by 30.1 and 15.0 µg m⁻³ in croplands and urban areas, respectively. During daytime, the decreased O₃ concentrations are 15.6 and 9.7 µg m⁻³, respectively. In contrast, soil emissions elevate ambient PM_2.5 concentrations by more than 4.5 µg m⁻³, accounting for 12 % of the PM_2.5 mass over the NCP in March 2020. The opposite effects are challenging for China to improve air quality, because China is the world's largest consumer of food, and food production strongly depends on N-fertilizer input. As the emission from fossil fuel combustion has been gradually decreasing, emissions from agricultural fertilization have increasing implications for air quality. We thus highlight that reducing NO_x emissions from agricultural fertilization is of great importance to air quality improvement. In China, the excessive use of N fertilizer still remains severe (Sun et al., 2022; Vitousek et al., 2009; Zhao et al., 2006), though a lot of efforts have been taken to increase the N-fertilizer efficiency and to reduce N losses from fertilizer (Li et al., 2018; Qiao et al., 2022; Yan et al., 2008). Fortunately, the consumption of N fertilizer reached its peak in 2014 in China and has been decreasing since then (Yu et al., 2022). Policymakers should manage to further reduce emissions from N-fertilizer application; for example, improving N-fertilizer efficiency and developing alternative fertilizers friendly to the environment are highly necessary. These measures will greatly minimize the adverse effects of agricultural fertilization on air quality, human health, and the ecological environment. Nevertheless, one should be aware of the limitation in the present case study that there are only three months of simulation as the basis for all of the insights into the soil NO_x emission and its influences on atmospheric chemistry and composition. More studies in terms of soil NO_x emissions, particularly during springtime, are in need to validate and generalize our model results.

We should not ignore the uncertainties regarding the BDSNP scheme. The default $\mathrm{75}/\mathrm{25}$ split in fertilizer application may not be widely suitable for the globe. Using a $\mathrm{20}/\mathrm{80}$ split, as is commonly used in crop modeling (Jägermeyr et al., 2021), or a $\mathrm{40}/\mathrm{60}$ split, as has been reported to be common in Northeast China (Zheng et al., 2023), would be expected to result in differences in the magnitude and timing of emissions compared to the default scheme. Among other impacts, there is no canopy interception of emitted soil NO_x at planting, which would result in substantially larger emissions to the atmosphere under the default $\mathrm{75}/\mathrm{25}$ split than in a $\mathrm{20}/\mathrm{80}$ split. In addition, because BDSNP applies the 25 % topdressing application evenly over the growing season following the 75 % basal application, it is less likely to produce sizable pulses of emissions. Because fertilizer applications are kept constant and global fertilizer emissions are constrained to 1.8 Tg N, the BDSNP mechanism is unable to reproduce historical trends or capture significant interannual variability in emissions.

The discrepancies between NO₂ and NH₃ column densities suggest substantial differences in the soil emission mechanisms of NO_x and NH₃, especially after fertilization. Future studies could incorporate a dynamic bidirectional NH₃ scheme alongside the BDSNP scheme to further investigate the nature of fertilizer-induced emission pulses.