Enhanced sulfate formation through SO2+NO2 heterogeneous reactions during heavy winter haze in the Yangtze River Delta region, China

Rapid sulfate formation is recognized as key characteristics of severe winter haze in China. However, air quality 15 models tend to underestimate sulfate formation during heavy haze periods and heterogeneous formation pathways have been proposed as promising mechanisms to reduce gaps between observation and model simulation. In this study, we implemented a reactive SO2 uptake mechanism through the SO2 + NO2 heterogeneous reactions in the Comprehensive Air Quality Model with extensions (CAMx) to improve simulation of sulfate formation in the Yangtze River Delta (YRD) region. Parameterization of the SO2 + NO2 heterogeneous reactions is based on observations in Beijing and considered both impact 20 of relative humidity and aerosol pH on sulfate formation. Ammonia is reported to be critical for the formation of secondary inorganic aerosols. Estimation of ammonia emissions is usually associated with large uncertainties and model tends to underestimate ammonia concentrations substantially. Sensitivity tests were conducted to evaluate the importance of the SO2 + NO2 heterogeneous reactions as well as ammonia emissions on modelled sulfate concentrations during a period with several heavy haze episodes in the YRD region. Base case model results show large underestimation of sulfate 25 concentrations by 36 % under polluted conditions in the YRD region. Adding the SO2 + NO2 heterogeneous reactions or doubling ammonia emissions alone leads to slight model improvement (~6 %) on simulated sulfate concentrations in the YRD region. However, model performance significantly improved when both the SO2 + NO2 heterogeneous reactions and doubled ammonia emissions were included in the simulation: predicted sulfate concentrations during polluted periods increased from 23.1 μg m -3 in the base scenario to 29.1 μg m -3 (representing an increase of 26 %). Aerosol pH is crucial for 30 the SO2 + NO2 heterogeneous reactions and our calculated aerosol pH is always acidic and increased by 0.7 with doubled ammonia emissions. Modelling results also show that this reactive SO2 uptake mechanism enhanced sulfate simulations by 1 to 5 μg m -3 for the majority of eastern and central part of China, with more than 20 μg m -3 increase of sulfate concentrations

over the north-eastern plain. These findings suggest that the SO 2 + NO 2 heterogeneous reactions could be important for sulfate formation in the YRD region as well as other parts of China. More studies are needed to improve the parameterization of the SO 2 + NO 2 heterogeneous reactions based on local data further evaluate this mechanism in other regions. In addition, ammonia emissions were found to be a key driving variable of the spatial patterns of sulfate enhancement due to the new pathway. Substantial efforts are needed to improve the accuracy of ammonia emissions 5 inventory.

Introduction
Rapid sulfate (SO 4 2-) formation has been reported to be key characteristics of severe winter haze in China. However, most air quality models tend to underestimate sulfate formation during severe winter haze episodes in China because standard SO 2 oxidation pathways, including gas-phase chemistry (i.e. oxidized by hydroxyl radical OH) and aqueous-phase chemistry (i.e. 10 oxidized by ozone (O 3 ), hydrogen peroxide (H 2 O 2 )) are suppressed by weak photochemical activity and low ozone concentrations (Quan et al., 2014)). Meanwhile, analysis of severe haze events in China show enhanced secondary inorganic aerosols, especially sulfate concentrations. For example, Quan et al. (2014) found that observed sulfate accounted for 13 % of PM 2.5 (particulate matter with dynamic equivalent diameter less than 2.5 µm) on normal clean days and increased to 25 % on haze days during the infamous 2013 January Beijing haze period. For the same haze episode, Cheng et al. (2016) used 15 concentration ratios of sulfate to sulfur dioxide ([SO 4 2-]/[SO 2 ]) to diagnose sulfate production rate; this ratio increased with PM 2.5 levels and was 6 times higher under the most polluted conditions than normal conditions. Most current air quality models (e.g. CMAQ, GEOS-Chem, WRF-Chem, CAMx), which only include the traditional gaseous-or aqueous-phase mechanisms for sulfate formation, do not show very good model performances of sulfate concentrations against observations during haze periods in China Zheng et al., 2015;Gao et al., 2016aGao et al., , 2016bLi et al., 2015). The under-20 prediction of sulfate concentrations could be related to uncertainties of the emissions inventory, bias of simulated meteorological fields, and/or some missing sulfate formation mechanisms that are not included in the current models.
Heterogeneous sulfate production chemistry has been proposed by several studies to explain the high concentrations and rapid formation of sulfate during haze episode in China (e.g. He et al., 2014;Wang et al., 2014;Zheng et al., 2015;Cheng et al., 2016;Guo et al., 2017). He et al. (2014) suggested a synergistic effect between NO 2 and SO 2 on the 25 surface of mineral dust (i.e. mineral oxides) as an important source of sulfate in China and emphasized the essential role of O 2 involved in this process. More generally, heterogeneous loss of SO 2 on aerosol surfaces (not limited to mineral dust) or deliquescent aerosols is discussed by many studies, although the exact underlying mechanism is still unknown (e.g. Zheng et al., 2015). For this kind of heterogeneous reactions, the sulfate production rate has been parameterized as a pseudo first-order reaction with respect to the gaseous SO 2 concentration with the SO 2 reactive uptake coefficient (γ) on 30 aerosol surfaces being the key parameter. This uptake coefficient, representing the probability that a SO 2 gas molecule colliding with an aerosol surface results in sulfate formation, is reported to be heavily dependent on relative humidity (RH) (Zheng et al., 2015;Wang et al., 2016). Parameterized reactive uptake of SO 2 has been implemented in several current air quality models, including GEOS-Chem, WRF-Chem, CMAQ and CAMx, and generally improved model performance of sulfate concentrations during haze episodes in China (e.g. Wang et al., 2014;Zheng et al., 2015;Gao et al., 2016b). Two more recent papers, Wang et al. (2016) and Cheng et al. (2016) further suggested that reaction between NO 2 and SO 2 in aerosol water may contribute substantially to sulfate formation during haze events in China. Both studies emphasized the 5 importance of higher aerosol pH (5.4-6.2 reported by Cheng et al. (2016) and 6.0-7.6 by Wang et al. (2016)) sustained by abundant gas-phase ammonia (NH 3 ) during haze periods being an essential precondition for this mechanism. However, the near-neutralized aerosol pH that facilities SO 2 oxidation by NO 2 is questioned by Guo et al. (2017) who concluded from a thermodynamic analysis that aerosol pH was always acidic (4.  regardless of the ambient NH 3 concentrations and that the NO 2 -mediated oxidation of SO 2 was unlikely to be important in China or any other region of the world. Guo et al. (2017) 10 pointed out that within low pH ranges (up to 4.5), SO 2 oxidation catalyzed by transition metal (i.e. Fe(III) and Mn(II)) might become a dominant sulfate formation pathway in aerosol water and suggested it as an alternative to SO 2 + NO 2 reactive uptake as being a potential sulfate contributor under haze conditions. Similar conclusion is also made from a most recent work by Shao et al. (2019), who implemented four heterogeneous sulfate formation mechanisms in GEOS-Chem and assessed model performance using sulfate oxygen isotopes data in Beijing, who found that SO 2 oxidation catalyzed by 15 transition metal ion (TMI) to be the dominant sulfate formation mechanism. On the contrary, another slightly earlier study by Ye et al. (2018) concluded SO 2 oxidation by H 2 O 2 was the dominant pathway based on observations of atmospheric H 2 O 2 concentrations in Beijing. Song et al. (2018b) suggested the heterogeneous hydroxymethanesulfonate (HMS) chemistry being a potentially important contributor to heavy haze pollution in northern China. Hung et al. (2018) reported the interfacial SO 2 oxidation on the surface of aqueous micro-droplets as a potential pathway to explain fast conversion of SO 2 20 to sulfate.
To investigate whether the SO 2 + NO 2 reactions in aerosol water could help better predict the enhanced sulfate formation during haze periods in the Yangtze River Delta (YRD) region, we implemented a parameterized SO 2 + NO 2 reactive uptake mechanism in the Comprehensive Air Quality Model with Extensions (CAMx), which is a widely used air quality model in China (e.g. Wang et al., 2009;Huang et al., 2012;Li et al., 2013Li et al., , 2015Jia et al., 2017;etc.). Our parameterization 25 specifically incorporated RH and aerosol pH dependencies derived from measurement data during the 2015 Beijing haze event . Although the RH dependency of the SO 2 uptake rate has already been implemented in previous studies (e.g. Zheng et al., 2015;Wang et al., 2014), the effect of aerosol pH has not been explicitly included in most of the previous modelling studies, except for a most recent study by Shao et al (2019), who also considered aerosol pH in their model parameterization.. While most of the previous studies were trying to improve model predictions in the northern part of 30 China, especially the Beijing-Tianjin-Hebei region (e.g. Gao et al., 2016b;Zheng et al., 2015), this work is one of the few studies that focus on the Yangtze River Delta region, which has also suffered from severe haze problems in recent years due to urban expansion and industrialization (e.g. Li et al., 2011;M. Wang et al., 2015;Xu et al., 2016;Ming et al., 2017). In addition to the SO 2 + NO 2 heterogeneous reactions, we also investigated model sensitivities to ammonia emissions, which have been reported to be crucial for the formation of secondary inorganic aerosols and large uncertainties exist with current ammonia emission inventory (Huang et al., 2011;Fu et al., 2013).

Current sulfate formation pathways in CAMx
In this study, CAMx version 6.40 (Ramboll Environ, 2016) was used as the base model to simulate sulfate formation. Table  5 1 lists the sulfate formation pathways that are currently considered implemented in standard CAMx source code. In addition to the traditional SO 2 oxidation by OH in the gas phase and O 3 , H 2 O 2 , and O 2 (catalyzed by Fe(III)/Mn(II)) in cloud droplets, sulfate formation through reactions with methyl hydroperoxide and other organic hydroperoxides (MHP) as well as peracetic and other organic peracids (PAA) in the aqueous phase is also included. For heterogeneous formation pathway, the SO 2 + NO 2 reaction is currently considered and implemented as pseudo gas phase reaction with the rate parameterization based on 10 results from Zheng et al (2015), where the key parameter (i.e. gamma) is bounded between a lower and upper limit and changes linearly in response to RH. This relatively simple parameterization of SO 2 + NO 2 heterogeneous reaction has been included in many previous studies, e.g. Y. Wang et al. (2014), B. Zheng et al. (2015, etc.

SO 2 + NO 2 mechanism in CAMx
In this study, we implemented the SO 2 + NO 2 reactive uptake mechanism in CAMx version 6.40 (Ramboll Environ, 2016) as 15 a pseudo gas-phase reaction:

→
(1) Since the vapor pressure of sulfuric acid is very low, we assumed all sulfuric acid partitions to the aerosol phase. The rate constant k het is related to the reactive uptake coefficient γ for SO 2 as follows: where ̅ is the mean molecular speed (m/s), and S is the aerosol surface area concentration (m 2 /m 3 ). Based on the observations during the Chinese haze events , this uptake coefficient γ depends on aerosol pH, RH, and 20 NO 2 concentration. Therefore, we assumed a functional form of γ as the product of each of these dependencies: where k 0 (ppm -1 ) is the RH-dependent parameter; NO 2 (g) is the NO 2 gas concentration; d f is the pH-dependent distribution factor of SO 2 , i.e. the ratio of SO 2 concentration in the aqueous-phase to the gaseous-phase and is calculated as Eq. (4) in the model: where H eff is the effective Henry's low constant of SO 2 (M atm -1 ), R is the universal gas constant (L atm mol -1 K -1 ), T is air temperature (K) and w L is the aerosol water content (μg m -3 ).We used the data in Table S2 and S5 of Wang et al. (2016) to back calculate the RH dependency of k 0 under clean (observed sulfate concentration less than 10 µg m -3 ), transition (sulfate between 10 and 20 µg m -3 ), and polluted (sulfate more than 20 µg m -3 ) conditions during Beijing 2015 episodes. Aerosol pH was calculated using the ISORROPIA thermodynamic equilibrium model implemented in CAMx assuming a metastable 5 aerosol liquid phase which is an appropriate assumption for most ambient conditions including the Chinese haze events (Guo et al. 2017). Wang et al (2016) only reported NOx (not NO 2 ) concentrations in Beijing during the 2015 haze event. We simply assumed a NO 2 /NOx ratio of 0.5. Inserting NO 2 concentrations, γ values from Wang et al. (2016), and calculated aerosol pH from ISORROPIA into Eq. 3, we obtained the expression of k 0 depending upon RH as follows (parameters for k 0 calculation is shown in Table S1): 10 (5) Due to lack of observation data at high RH values, we set a constant k 0 value when RH increases from 56% and up. This would lead to underestimated sulfate formation due to the SO 2 + NO 2 heterogeneous reactions at high RH values, which is a favorable condition for the heterogeneous sulfate production. In addition, the differences of aerosol hygroscopicity in Beijing vs. Shanghai could add more uncertainties in the dependency of k 0 on RH. Reported values of hygroscopicity parameter ĸ were 0.25~0.31 for Shanghai (Ye et al., 2011;, which are higher than values reported for Beijing (0.14~0.24;Massling 15 et al., 2009). Therefore, our results represent a relatively conservative estimation of sulfate formation. The rate constant k het of SO 2 + NO 2 is formulated as: ̅ SO 2 lifetime (in hr) associated with the SO 2 + NO 2 reactive uptake mechanism is calculated as: Figure 1 shows the SO 2 lifetime as a function of aerosol pH for clean, transition, and polluted conditions, with other variables kept constant. The SO 2 lifetime shortens as aerosol pH becomes more neutralized, indicating faster conversion of 20 SO 2 to sulfate by SO 2 + NO 2 reactive uptake on aerosol. For pH within 2 to 7, one unit increase in aerosol pH shortens SO 2 lifetime by about one order of magnitude. The blue, orange, and red symbols in Figure 1 correspond to the clean, transition, and polluted conditions during Beijing 2015 based on data in Table S1. As shown in Figure 1, the aerosol pH values calculated by ISORROPIA are 5.5 (for clean conditions) and 4.1-4.2 (for transition and polluted conditions), all lower than the values (7.6) reported by Wang et al. (2016). As noted by Guo et al. (2017), it is important to make a consistent 25 assumption for aerosol state (i.e., metastable) in deriving and implementing the parameterization for reactive uptake. A most recent paper by Song et al (2018a) identified coding errors with the ISORROPIA model, which resulted unrealistic pH values of 7.7 using the standard ISORROPIA model with the stable state assumption in previous studies. Nevertheless, our results are not compromised by this coding error because the metastable assumption was chosen in our ISORROPIA calculation.

Model configuration 5
Two versions of the Comprehensive Air Quality Model with extensions (CAMx) modified based on the original version 6.40 (Ramboll Environ, 2016) were used in this study: one with the SO 2 + NO 2 heterogeneous reactions (described in Section 2.1) and one without (forcing k het equals to zero). The modeling domain consists of three nested grids ( Figure  Boundary conditions for D01 were generated from the Model for OZone And Related chemical Tracers (MOZART) global 15 chemistry model (Emmons et al., 2010). The Carbon Bond 6 (CB6) mechanism (Yarwood et al., 2010) was used for the gas phase chemistry and the static two-mode coarse/fine (CF) scheme was used to represent particle size distribution. The Zhang dry deposition (Zhang et al. 2003) and default wet deposition scheme was used to for removal processes. Anthropogenic emissions for areas outside the YRD region were from the Multi-resolution Emission Inventory for China (MEIC, http://www.meicmodel.org/). For emissions within the YRD region, an YRD-specific emission inventory (Huang et al., 2011;20 Li et al., 2011) was updated to year 2014 and utilized in this study. This YRD-specific emission inventory includes emissions from sources of combustion, industry, mobile and residential. Primary sulfate emissions over the 4km domain are estimated to be 994 tons day -1 for December 2013 (accounting for 14.8% of primary PM 2.5 ) with dense emissions from Shanghai and southern Jiangsu province (see Figure S1 for spatial distribution). At the SAES site, primary sulfate emissions were estimated to 757 kg per month (only accounting for 1.0% of primary PM 2.5 ). Biogenic emissions were simulated using 25 the Model of Emissions of Gases and Aerosols from Nature (MEGAN, version 2.1, Guenther et al. 2012) based on the WRF simulation results. The modeling episode is December 2013, during which several heavy haze events with hourly PM 2.5 concentration higher than 500 µg m -3 were observed in the YRD region.
Four simulations with identical model configuration and input data including meteorology, initial/boundary conditions, and emission inventory (except ammonia emissions) were conducted using the above two different CAMx versions: 30 -noHet (base case): simulation based on CAMx version without the SO 2 + NO 2 heterogeneous reactions (this is also our base case). Note that this CAMx version differs from the distributed CAMx v6.40 in that we removed the original heterogeneous sulfate formation reaction with NO 2 which only included a simple parameterization based on RH (ref. reaction No.7 in Table 1) in the distributed version. This is done on purpose to quantify the influence of the newly parameterized SO 2 + NO 2 heterogeneous reactions on sulfate formation.
-Het: simulation based on CAMx with the SO 2 + NO 2 heterogeneous reactions. Other model configurations were identical to scenario noHet.
-noHet_2NH 3 : CAMx version and model configurations were same as scenario noHet except ammonia emissions were 5 doubled for the 4 km domain.
-Het_2NH 3 : CAMx version and model configurations were same as scenario Het but ammonia emissions were doubled for the 4 km domain.
We first ran CAMx for 36 km/12 km domains with two-way nested; for the 4 km domain, we used boundary conditions extracted from the 12 km model outputs and conducted the above four scenarios. Fourteen vertical layers were used 10 extending from the surface to 100 mb. In addition to default CAMx outputs, we modified the source code to generate additional diagnostic variables (e.g. aerosol pH, RH, and k het ) to evaluate the SO 2 + NO 2 heterogeneous reactions.

Observations
Hourly observations of ozone, SO 2 , NO 2 , PM 2.5 and its components including sulfate, nitrate, ammonium, organic carbon (OC), and elemental carbon (EC) are available between 1 December and 29 December 2013 at a monitor site located at the 15 center of the urban area of Shanghai (referred as SAES site, 31.1695 º N, 121.4305 º E, Figure 3). Hourly PM 2.5 observations are also available at another 23 monitor sites across the YRD region ( Figure 3; see locations in Table S2). During this period, YRD region experienced relative clean days as well as several heavy haze episodes with peak PM 2.5 exceeding 600 µg m -3 during a most heavily polluted period of December 5 th to 7 th . At the SAES site, maximum hourly PM 2.5 concentration reached 540.3 µg m -3 on December 6 th with a monthly average of 118.7 µg m -3 . We followed the method in Wang et al (2016) 20 to divide the period into clean (observed sulfate <10 µg m -3 ), transition (10-20 µg m -3 ), and polluted (>20 µg m -3 ) periods based on observed hourly sulfate concentration at the SAES site. Compared with clean period, all PM species increased by more than 3 times (sulfate, nitrate and ammonium (SNA) increase by more than 5 times) during polluted period as indicated by the enhancement ratio (calculated as the ratio of average concentrations during the polluted period divided by those during the clean period). In terms of fraction of PM 2.5 , SNA increased from 44 % during clean period to 69 % during 25 polluted period while carbonaceous aerosols (OC and EC) decreased from 32 % to 24 %. This is consistent with the commonly observed characteristics of winter haze periods in China reported by many previous studies (e.g. Wang et al., 2014;Zheng et al., 2015b;Cheng et al., 2016) that SNA is playing a more important role during the heavy haze periods.
Average sulfate concentration of clean, transition and polluted periods was 6.7, 14.2, and 36.1 µg m -3 , respectively, accounting for 17-23 % of PM 2.5 ( Figure S2). 30 Observations of ambient ammonia concentrations are also available at the SAES site; however, the quality of measurements is questionable. Therefore, we used ammonia observations from a similar urban site nearby (referred as FDU site, ~15 km north from the SAES site, 31.3005 º N, 120.9778º E, Figure 3) for analysis in this study. Observations at the FDU site have been discussed by S.  and demonstrated data reliability. Diurnal NH 3 concentrations at the FDU site during our modeling period showed a weak bimodal pattern with an average of 7.3 ppb (ranging 1.6-25 ppb) during this period ( Figure S3). This two-peak diurnal variation is caused by vehicle emissions and evolution of the boundary layer (S. . In summary, observations for gases species (except NH 3 ) and PM species at the SAES site and NH 3 at the FDU site were used for model validation in this study. 5

Statistical metrics for model validation
For WRF and CAMx model performance valuation, mean bias (MB), normalized mean bias (NMB), and index of agreement (IOA) were used in this study. Calculations of these selected metrics are shown below: where P j and O j are predicted and observed hourly concentrations or values, respectively. N is the number of paired model and observation data. is the average concentration/value of observations. IOA ranges from 0 to 1 with 1 indicating perfect 10 agreement between model and observation.

WRF results evaluation
Model performance of WRF results is generally acceptable in this study. Table S3 summarizes  Temperature and relative humidity were well reproduced with NMB and NME within 37% and 41%, respectively; IOA values are above 0.8. Wind speed is overestimated with a MB of 1.5 m s -1 at Pudong and 0.5 m s -1 at Hongqiao station; NMB of predicted wind direction at the two stations is -36% and -27%, respectively. Comparisons of hourly observed and simulated relative humidity, wind speed and temperature at these two stations suggest reasonable model results in terms of 20 temporal variations ( Figure S4). Overall, the WRF simulated results are acceptable to be used in subsequent CAMx simulations.

CAMx base scenario (noHet) evaluation
Figure 4 depicts the time series of simulated and observed concentrations for sulfate and PM 2.5 during 1 to 29 December 2013 at SAES site (see Figure S5 in Supplemental Information for other species). Overall, the model is successful in capturing the temporal variations of ozone and PM species with IOA values above 0.5 (Table S4). Nevertheless, model tends to systematically underestimate all gaseous and PM species with NMB values ranging from -5% for NO 2 to -68% for NH 3 . 5 This could be partially explained by the higher simulated wind speeds compared with observed values, especially at Pudong station where the observed average wind speed during the modeling period was 4.5 m s -1 while simulated wind speed was 6.0 m s -1 , representing an overprediction by 33%. For sulfate, the model captured the day-to-day sulfate variations reasonably well with an overall MB of -2.8 µg m -3 and IOA of 0.80. For clean and transition periods, model showed slight over-prediction with MB of 1.1 and 0.5 µg m -3 (Table S5)  well captured by the model. For sites located in southern Jiangsu and southern Zhejiang province, the model showed favorable agreement with the observations. Underestimations existed for sites located in the northern part of Jiangsu and Zhejiang province. MB across all 24 monitoring sites ranged from as low as -90.4 µg m -3 (site in north Jiangsu province) to slight overestimation of 11.4 µg m -3 (site in south Zhejiang province); corresponding NMB ranged from -46 % to 16 % (Table S2). 25

Simulated sulfate concentrations at SAES site
Four scenarios -noHet, Het, noHet_2NH 3 and Het_2NH 3 were conducted to evaluate the impact of the SO 2 + NO 2 heterogeneous reactions and ammonia emissions on sulfate simulation. We first analyzed the modeled sulfate results at the SAES site; then we discussed the spatial patterns over the YRD region. Similar discussions of nitrate, ammonium and PM 2.5 are included in the supplemental information. Table 2 shows the average sulfate concentration for different scenarios by  30 clean, transition, and polluted periods; corresponding scatter plots are shown in Figure 6. A complete summary of statistical metrics for each scenario/period is presented in Table S5.

Impact of SO 2 + NO 2 heterogeneous reactions (noHet vs. Het)
As shown in Figure 6, simulated sulfate concentrations compared well with observations under clean and transition conditions in the noHet scenario with over-prediction by 16 % and 4 %, respectively. By contrast, large under-prediction of sulfate concentration existed during polluted periods (MB of -13.0 µg m -3 , NMB of -36 %). Adding the SO 2 + NO 2 heterogeneous reactions showed small enhancement on sulfate formation, reducing the overall NMB from -16 % to -12 %. If 5 only polluted periods are considered, simulated sulfate concentrations increased from 23.1 to 24.6 µg m -3 with the heterogeneous reactions, corresponding to an increase by 6.5 %. Thus even with the SO 2 + NO 2 heterogeneous reactions, model was still under-predicting sulfate concentrations on heavy haze days with a NMB of -32 %. This is because aerosol pH was always acidic (pH < 3; this will be discussed in the following section) and the SO 2 + NO 2 heterogeneous reactions were not being appreciable within this pH range (Figure 1). Model performances for clean and transition periods were 10 slightly compromised with the SO 2 + NO 2 heterogeneous reactions since the base scenario was already overestimating sulfate concentrations.

Impact of NH 3 emissions (noHet vs. noHet_2NH 3 )
Being the dominant base gas in the atmosphere, ammonia plays an essential role in the formation of secondary inorganic aerosols and estimation of ammonia emissions is usually associated with large uncertainties (e.g. Huang et al., 2011;Fu et al., 15 2013). With the base case ammonia emissions, NH 3 concentration was under-predicted by 3.0 ppb (NMB of -60 %). With doubled ammonia emissions, ammonia concentration was over-predicted by 1.7 ppb with NMB of 34 % but the MB of the total ammonia (NH 3 + ammonium) concentrations were reduced from -6.9 µg m -3 (NMB of -36%) in the base case scenario to -1.9 µg m -3 (NMB of -10%). NMB of sulfate concentrations during polluted period is -32 %, which is similar to the enhancement caused by that of the SO 2 + NO 2 heterogeneous reactions. Clearly, doubling ammonia emissions is not enough 20 to close the gap between observed and simulated sulfate concentrations during heavy haze periods. We performed additional sensitivity tests with even higher ammonia emissions and found that 10 times ammonia emissions would be needed to achieve an average sulfate concentration (33.2 µg m -3 ) that is comparable with observation (36.1 µg m -3 ) under polluted conditions (with no heterogeneous reactions). However, in that case, model performance of ammonia is significantly compromised with over-prediction by 32.3 ppb. These results indicate that the uncertainties associated with the ammonia 25 emissions are not enough to fully explain the under-prediction of sulfate formation during heavy haze periods in the YRD region.

Impact of both (noHet vs. Het_2NH 3 )
A fourth scenario (Het_2NH 3 ) with the SO 2 + NO 2 heterogeneous reactions as well as doubled ammonia emissions gave the best model performance of sulfate concentrations with an overall MB of -0.2 µg m -3 (NMB of -1 %, Figure 6). During 30 polluted periods, average sulfate concentration was predicted to be 29.1 µg m -3 (representing an increase of 26% from the base case) and NMB was reduced from -36 % in the base scenario to -19 % in the Het_2NH 3 scenario. Maximum sulfate concentration simulated under scenario Het_2NH 3 was 97.2 µg m -3 , which compared well with the observed maximum of 93.4 µg m -3 at the SAES site. With doubled ammonia emissions, the heterogeneous reactions were playing an increasing important role in sulfate formation by boosting average sulfate concentrations from 24.5 (noHet_2NH 3 ) to 29.1 µg m -3 (Het_2NH 3 ) under polluted conditions, representing an increase by 19 %. This is because with more ammonia available, aerosol pH was elevated by ~0.7 pushing it closer towards the actual pH (as discussed more in section 3.3) and the rate of the heterogeneous reactions is positively correlated with aerosol pH (Figure 1), therefore leading to the best model performances from the Het_2NH 3 scenario. These results indicate that the SO 2 + NO 2 heterogeneous reactions as well as sufficient 5 ammonia emissions are both needed to greatly improve model simulation of sulfate formation under polluted conditions. However, it is to mention that model performance under clean and transition periods got compromised most under scenario Het_2NH 3 . Figure S6 shows a Q-Q plot of modeled versus observed sulfate concentrations for the four scenarios. Underestimations of sulfate concentrations become noticeable around 35 µg m -3 in all scenarios and between 35 to 55 µg m -3 , there appears to be 10 a systematical low bias in predicted sulfate concentrations that neither doubled ammonia emissions nor the heterogeneous reactions or both could stimulate notable sulfate formation. Scenario Het_2NH 3 gives the best model performance with an overall MB of -0.2 µg m -3 but still underpredicts sulfate formation under heavy haze periods by -19 %. This could be related to still biased ammonia emissions, less direct emissions of sulfate and/or SO 2 , and/or missing of other sulfate formation pathways that needs further investigation. For example, Shao et al. (2019) included heterogeneous sulfate formation via 15 oxidations by O 3 , H 2 O 2 , and Fe(III)/Mn(II), in addition to the aqueous phase reactions and concluded that the metalcatalyzed reactions dominated the heterogeneous sulfate formation. These heterogeneous reactions were not included in the current study and could lead to some underestimated of sulfate formation. As mentioned above, the parameterization of the k 0 values is relatively conservative at high RH conditions, which are favorable for sulfate formation. In addition, reported aerosol hygroscopicity Bias in meteorology could also play some roles here as we are seeing systematically underestimation 20 of all gaseous and PM species. Another explanation is that the SO 2 + NO 2 heterogeneous reactions implemented in this study were parameterized based on observations in Beijing but the simulation is performed over the YRD region. It would be ideal to use local observations for model parameterization in future studies.

Sulfate formation budget
To gain a closer look at the sulfate formation via different pathways (e.g. gas phase vs. aqueous phase vs. heterogeneous 25 phase, Table 1), we constructed a sulfate formation budget in a similar manner as Shao et al. (2019). Figure 7 shows the relative contribution of primary sulfate emissions as well as individual sulfate formation pathway to the total sulfate concentrations at the SAES site under different conditions. Overall, primary sulfate emissions and secondary formation accounted for half of the total sulfate concentrations. Of the secondary sulfate, gas-phase reactions always dominated secondary sulfate formation, with relatively consistent contribution around 38~39% under different conditions. As pollution 30 developed, contribution from secondary formation exceeded that of primary emissions, accounting for 60% of total sulfate abundances under polluted conditions. In contrast to the relatively consistent contribution from the gas-phase formation, both aqueous and heterogeneous sulfate formation doubled from clean to polluted periods, with relative contribution increased from 4.1% to 9.4% for the former and from 5.0% to 12.6% for the latter.
If we exclude the contribution of primary sulfate emissions (i.e. smaller pie chart in Figure 7), the absolute sulfate formation via the gas-phase reactions more than doubled from clean (1.59 μg m -3 ) to polluted (3.61 μg m -3 ) periods; however, the relative contribution from gas-phase formation among all formation pathways dropped from 80.9% to 63.3% as pollution developed. Sulfate formation from all aqueous phase reactions increased from 0.17 μg m -3 under clean conditions to 0.89 μg m -3 under polluted conditions, corresponding to an increase of relative contribution from 8.6% to 15.6%. Under all 5 conditions, aqueous oxidation due to MHP and PAA is negligible, with less than 1% of sulfate contribution. The rest three aqueous pathways in turn dominated aqueous sulfate formation depending on the specific condition. For instance, under  Shao et al. (2019) and other studies cited in the paper but the overall magnitudes are well comparable. We realize that assuming constant Fe and Mn mass fraction is a simplification and latest CAMx version has the option to treat Fe and Mn as primary species. However, using this option would put even more burden on the emission inventory to have accurate source speciation profiles for different source sectors. Nevertheless, although this Fe(III)/Mn(II) catalyzed pathways stands out among all aqueous pathways under polluted conditions, the relative contribution (8.4%) is only about one third of that from the SO 2 + NO 2 heterogeneous reactions (21.1%). As for the SO 2 + NO 2 heterogeneous reactions, its contribution to sulfate formation doubled from 10.5% (0.21 μg m -3 ) under clean conditions to 21.1% (1.2 μg m -3 ) under polluted conditions. Under all conditions, the relative contribution of the SO 2 + NO 2 heterogeneous reactions exceeds the sum of all aqueous pathways, indicating the importance of heterogeneous oxidation 5 pathways compared to aqueous pathways.

Sulfate formation under selected episodes
We further selected four heavy haze episodes (EP1-EP4) with observed sulfate concentrations continuously exceeding 30 µg/m 3 (as highlighted in Figure 4) at the SAES site. These episodes lasted from 9 hours (EP2) to as long as 37 hours (EP1) with episode average sulfate concentrations are all above 50 µg m -3 ( Figure S7) except for EP3 (36.2 µg m -3 ) ( Table S6). 10 Maximum hourly sulfate concentrations ranged from 48.6 µg m -3 for EP3 to 93.4 µg m -3 for EP2. The averaged molar sulfate and SO 2 ratio ([SO 4 2-]/[SO 2 ]) for EP1 and EP2 are higher (0.52 and 0.70, respectively) than that for EP3 (0.17) and EP4 (0.19). In the base case scenario, sulfate formation was significantly underestimated for all four episodes with NMB ranging from -39% to -72%. Figure S8 shows the sulfate formation budget for the four episodes of the base case scenario. The gas phase oxidation pathway was the dominant contributor, accounting for 52% (EP2) to 79% (EP3) of total secondary sulfate 15 formation, followed by the SO 2 + NO 2 heterogeneous reactions with contributions of 20% ~ 39%. For EP1 and EP2, the Fe/Mn-catalyzed oxidation pathway contributed ~10% of sulfate formation but were negligible for the other two episodes. It is interesting to note that for all selected episodes except EP3, sulfate formation was enhanced in scenario Het_2NH 3 by 10.4 to 14.6 µg m -3 while EP3 only exhibits minimal increase of modeled sulfate concentrations by only 0.8 µg m -3 . We performed additional sensitivity tests and found that even with 10 times ammonia emissions, modeled sulfate concentration 20 during EP3 is enhanced by only 2.3 µg m -3 , which is still much lower compared to the observed values. We suspect that other factors, for example, meteorology might be biased during EP3 and lead to the underpredicted sulfate concentrations.
For instance, we looked at the model performance of WRF predictions for individual episode. All four episodes had some over-prediction of wind speeds with NMB ranging from 4% of EP2 to as much as 43% of EP3. Clearly, the large overprediction of wind speeds during EP3 contributed partially to the underestimated sulfate concentrations by the model. 25 Another potential cause for sulfate underprediction could be failure to capture episodic primary sulfate emissions during EP3.
When EP3 is excluded, modeled sulfate concentrations during heavy pollution episodes are greatly enhanced from 33.5 µg m -3 in the base scenario to 46.2 µg m -3 in scenario Het_2NH 3 (increase by 38 %), due to the combined influences of the SO 2 + NO 2 heterogeneous reactions and doubled ammonia emissions.

Observed and predicted aerosol pH at the SAES site 30
Aerosol pH, which is calculated from ISORROPIA either based on observations or within CAMx, is crucial for the heterogeneous SO 2 + NO 2 reactions to be effective. Observation-based aerosol pH was calculated using forward metastable mode by ISORROPIA to be consistent with CAMx ISORROPIA configuration. Figure 8 shows the distribution of observation-based and modeled aerosol pH at the SAES site by scenario/period. As indicated by both observation-based and modeled pH values, aerosols become more acidic as pollution develops. This is consistent with the higher SO 2 concentrations observed under polluted conditions ( Figure S9). For observation-based values, aerosol pH dropped by 35% from clean to polluted conditions while modeled aerosol pH dropped by 13~17% under different scenarios. As also shown by Figure 8, observation-based aerosol pH values are consistently higher than modeled values for all scenarios. Averaged 5 observed-based pH value during clean, transition, and polluted period is 5.5, 4.7, and 3.6, while corresponding value for base scenario (noHet) is 2.8, 2.6 and 2.3, each representing an underestimation by 48%, 45% and 34%. Maximum aerosol pH reached 5.0, 4.4, and 3.8 under clean, transition and polluted periods in the base scenario in contrast to observation-based values of 7.7, 6.5, and 5.3. Adding the SO 2 + NO 2 heterogeneous reactions causes small decrease (0.03~0.07) in predicted aerosol pH. The discrepancies between observation-based and model-based aerosol pH values might be due to significant 10 underprediction of NH 3 and ammonium concentrations. Therefore, when NH 3 emissions are doubled, modeled aerosol pH increases by ~0.7 to 3.0-3.5 and underestimation of aerosol pH for scenario noHet_2NH 3 is reduced to 36% during clean periods and 15% during polluted periods. Maximum aerosol pH during clean, transition and polluted periods is 5.7, 5.1, and 4.2 under scenario noHet_2NH 3 . Again, adding the SO 2 + NO 2 reactions on top of doubled NH 3 emissions slightly decreases the aerosol pH by 0.03-0.12, with stronger reduction associated with more enhancement of sulfate formation. Both 15 observation-based and model-based aerosol pH values at the SAES site indicate that aerosol pH is acidic, which is lower than the more neutralized values reported in previous studies for the  to 6.2 reported by Cheng et al. (2016) and 6.0 to 7.6 by Wang et al. (2016), the latter was later found to be associated with a coding bug in ISORROPIA). This difference might be due to lower ammonia levels in Shanghai compared with Beijing (S. . However, even when ammonia emissions are increased by 10 times, maximum modeled aerosol pH value is 4.8 under 20 polluted condition, which is still lower than the values reported for north China. Our results indicated that the aerosol pH at the SAES site tends to be moderately acidic regardless of the ambient ammonia concentrations. However, the acidity of aerosols in China still remains to be a vigorous debate. For example, Shi et al. (2017) reported a wide range of pH values between 0.33 and 13.6, depending on the source contributions. Xie et al. (2019) found that the predicted particulate pH values increased from moderate acid to near neutral with the increase of nitrate to sulfate molar ratio. 25 Figure 9 shows the spatial distribution of monthly mean sulfate, SO 2 , ammonia concentrations, and aerosol pH simulated in the base case and the differences between base case and other three sensitivity runs in the YRD region. Similar plots of nitrate, ammonium, and PM 2.5 are shown in Figure S10. Overall, impacts of the heterogeneous reactions and ammonia emissions over the YRD region are similar to that observed at the SAES site. With the SO 2 + NO 2 heterogeneous reactions 30 only, predicted monthly mean sulfate concentrations show ubiquitous increase of 0.1-5 µg m -3 across the domain with larger increase observed in the north and northwest directions. Regions with relative higher increase of predicted sulfate concentrations closely track regions with relatively high aerosol pH and high ammonia concentrations. Aerosol pH decreases slightly because more SO 2 is pulled into the aerosol phase. For nitrate concentrations ( Figure S10), however, the heterogeneous reactions lead to increase in the northwest region but decrease for the rest of the YRD region and magnitudes of changes in in both directions are within 1 µg m -3 . Predicted ammonium concentrations show less than 1 µg m -3 increase over the majority of the domain. Domain average PM 2.5 concentrations increased by 1.2 µg m -3 with spatial patterns similar to sulfate. 5

Spatial impact of the SO 2 + NO 2 heterogeneous reactions and ammonia emissions
With doubled ammonia emissions, predictions of all three inorganic PM species are enhanced with most profound impacts observed for nitrate (Figure 8 and Figure S10). Uniform increase across the YRD region is observed for predicted sulfate concentrations; for nitrate and ammonium, increase of predicted concentrations is more significant towards the south.
Aerosol pH is also elevated (on average by 0.3) with more ammonia available. In south Anhui and south Zhejiang provinces, 10 elevation of aerosol pH exceeds one unit. Areas with larger pH increase are also areas with relatively lower pH values in the base scenario, indicating that aerosol pH responds nonlinearly to changes in ammonia emissions.  China Plain increased from less than 20 µg m -3 during the base case scenario to more than exceed 30 µg m -3 in the Het scenario. For other regions, including the North China Plan and Sichuan Basin that show relative high sulfate concentrations in the base case scenario, sulfate concentrations were increased by 5-10 µg m -3 due to the implementation of the reactive SO 2 uptake mechanism. The spatial pattern of sulfate enhancement generally follows that of predicted ammonia concentrations, once again suggesting the important role of ammonia emissions for this mechanism. Future studies and local sulfate 30 observations are needed to evaluate this mechanism for other parts of China, especially for Northeast China Plain.
In this study, we implemented a new parameterization of the SO 2 + NO 2 heterogeneous reactions based on observations in Beijing to improve model simulation of sulfate formation under heavy haze conditions in the YRD region. Unlike previous studies that only considered the influence of relative humidity on sulfate formation, we also included the impact of aerosol pH in our parameterization. Four CAMx sensitivity runs were conducted to evaluate the importance of the SO 2 + NO 2 5 heterogeneous reactions as well as ammonia emissions on simulated sulfate concentrations in the YRD region. Base case simulation showed reasonable model performance of sulfate with an overall MB of -2.7 µg m -3 but significantly underpredicted sulfate concentrations by 36 % during polluted conditions. Implementation of the SO 2 + NO 2 heterogeneous reactions alone showed slight improvement of sulfate simulation (increase by 6.5 %) under polluted conditions due to acidic aerosol pH. Ammonia concentrations were significantly underestimated by the model. Doubling ammonia emissions alone 10 exhibited a similar impact (sulfate increase by 5.6 %) with that of the SO 2 + NO 2 heterogeneous reactions alone.
Nevertheless, aerosol pH increased by 0.7 with doubled ammonia emissions, which enabled the SO 2 + NO 2 heterogeneous reactions to become effective. Thus, in a fourth scenario where both the SO 2 + NO 2 heterogeneous reactions and doubled ammonia emissions were considered, simulated sulfate concentrations during polluted periods increased from 23.1 µg m -3 in the base case to 29.1 µg m -3 , representing an increase by 26 %. Results for sulfate simulations over entire China shows that 15 for some parts of China, especially the Northeast China Plain, implementing the SO 2 + NO 2 heterogeneous reactions could lead to as much as 20 µg m -3 increase of sulfate concentrations and the spatial pattern of sulfate enhancement follows closely to that of ammonia concentrations. These findings suggest that the SO 2 + NO 2 heterogeneous reactions could be important for sulfate formation under heavy haze periods and aerosol pH (in other words, ammonia emissions) is crucial in this process.
However, under-prediction of sulfate concentration still exists (by 20 %) in the YRD region under polluted conditions even 20 with the SO 2 + NO 2 heterogeneous reactions and doubled ammonia emissions, which urges further efforts to better constrain the parameterization of the SO 2 + NO 2 heterogeneous reactions using local data and to improve the accuracy of ammonia emissions inventory.
Date and code availability. All data and modified CAMx code is available upon request from the corresponding authors. 25 Competing interest. The authors declare that they have no conflict of interest.
Special issue statement. This article is part of the special issue -Multiphase chemistry of secondary aerosol formation under severe haze‖. It is not associated with a conference. Figure 1: SO 2 lifetime (in hr -1 ) due to SO 2 + NO 2 reactive uptake mechanism as a function of aerosol pH under clean, transition, and polluted conditions. Values of relative humidity, temperature, and NO 2 (g) concentrations are based on values in Table S1.  Table S2).