A Case Study of Investigating Secondary Organic Aerosol Formation Pathways in Beijing using an Observation-based SOA Box Model

Current modeling studies have underestimated secondary organic aerosol (SOA) levels in China to a larger degree than over Europe and the United States. In this study, we investigated the SOA formation pathways in urban Beijing for the period of November 7–8, 2014, using an observation-constrained box model in which the multigenerational oxidation processes of volatile organic compounds (VOCs) and intermediate VOCs (IVOCs) and the chemical aging of semi-volatile primary organic aerosols (POAs) were taken into account. The results demonstrated that the SOA formation rate was 30.3 μg m day in Beijing during the 2-day study period. The contributions of VOCs, IVOCs, and POAs to the SOA levels were 14%, 82%, and 4%, respectively. IVOC contributions were on a scale similar to the magnitude of underestimation in a previous study. The uncertainty analysis showed that SOA levels during the study period were 55.4– 102.4 μg m (the 25 and 75 percentiles of the sensitivity simulations). The contribution of IVOCs to the SOA formation was dominant compared with that of VOCs and POAs. A more precise IVOC oxidation mechanism can thus improve the performance of the SOA model in China.


INTRODUCTION
SOA has a significant effect on visibility, atmospheric radiative balance and human health (Liu et al., 2014).Secondary organic aerosols (SOA) have been reported to account for 31%-71% of the OA levels in Beijing (Guo et al., 2014a;Huang et al., 2014;Xu et al., 2015;Sun et al., 2016).Recent studies reported that the SOA concentration in Beijing was more than 17 µg m -3 , which was 3-8 times higher than the SOA concentration in other countries (Zhang et al., 2007).The large discrepancies between SOA simulations and observations are attributable to the limited understanding of the complicated chemical and physical processes underlying SOA formation (Hallquist et al., 2009).Traditional two-product mechanisms underestimate the annual SOA concentrations in China by as much as 75% (Jiang et al., 2012).Recently, a volatility basis set (VBS) model was developed using the latest chamber data (Donahue et al., 2006), and this model was expected to narrow the gap between simulations and measurements.Lin et al. (2015) demonstrated that the VBS model only accounted for 30%-36% of the observed OA levels in Beijing; nevertheless, it was a significant improvement over the two-product approach.The performance of the VBS model was different in Switzerland, where the simulated SOA concentration was 73% of the measured concentration (Bergstrom et al., 2012).
Uncertainties in the emissions of SOA precursors are considered a key cause of the underestimations of SOA concentration in China.Miller et al. (2016) compared satellite observations and simulations and reported that aromatics were likely underestimated by a factor of 2 over the Pearl River Delta.Furthermore, bottom-up emission inventories indicated that the uncertainties of VOCs and organic carbon emissions in China were ±(68-78)% and ±(258-271)%, respectively (Kurokawa et al., 2013;Li et al., 2017).Therefore, uncertainties in the emission inventory hinder the investigation of SOA formation pathways in China using three-dimensional air quality models.
In this study, we developed a VBS approach-based box model constrained by measured data in Beijing to minimize the effects of emission inventory-induced uncertainties.The newly reported SOA formation mechanisms, including primary organic aerosol aging and IVOC oxidation, were incorporated in the VBS scheme.The contributions from different precursors and formation pathways were quantified, which have not been evaluated in previous SOA simulations in China.Finally, a series of sensitivity simulations were conducted to assess the effects of uncertainties on the estimated SOA levels in Beijing.

Measurement Techniques
The observation site is located at the Institute of Atmospheric Physics, Chinese Academy of Sciences (116.4°E, 39.9°N), which is near major roads and is representative of the urban area in Beijing.Submicron aerosol particles were measured on-line on the roof of a two-story building (at an approximate height of 8 m off the ground level) by using an Aerosol Chemical Speciation Monitor (Ng et al., 2011), and CO, NO, O 3 , and SO 2 were collected using gas analyzers (Thermo Scientific).The hourly concentrations of gaseous HNO 3 and HONO were measured using the Gas and Aerosol Collector-Ion Chromatography system designed by Dong et al. (2012).A newly developed gas chromatography-mass spectrometry-flame ionization detection system was used for on-line measurements of VOCs, including alkenes, alkanes, cycloalkanes, and aromatic hydrocarbons (Wang et al., 2014).
In addition, temperature, pressure, relative humidity, and wind were measured at a meteorological tower at the study site, and the Planetary Boundary Layer (PBL) height was observed using a dual-wavelength (1064 and 532 nm) depolarization lidar (Sugimoto et al., 2000(Sugimoto et al., , 2001) ) developed by the National Institute for Environmental Studies, Japan.

Model Description
Supplemental Fig. S1 presents the conceptual diagram of the SOA scheme components in this study.VOCs (aromatics, terpenes and isoprene), IVOCs and POAs are SOA precursors.Their oxidation reactions with OH, O 3 and NO 3 produced four mixed sets of semi-volatile organic compounds (SVOCs), which were grouped on the basis of their effective saturation concentrations as follows: 1, 10, 100 and 1000 µg m -3 at 298 K. SVOCs generated from anthropogenic precursors is allowed to react with OH in the gas phase at a rate constant of 2 × 10 -11 cm 3 molecule -1 s -1 (Donahue et al., 2013), with each reaction resulting in a shift of the compound to the next lower volatility bin.
POAs traditionally treated as non-volatile in the model (Lane et al., 2008;Jiang et al., 2012;Lin et al., 2015).In this study, the POAs referred to here as primary SVOCs (P-SVOCs), which is regarded as semi-volatile and the vapor-phase portion can undergo photochemical oxidation (Robinson et al., 2007).The volatility distribution factors of POAs emission are shown in Base_vol (Table 1), which are taken from experimental results (May et al., 2013a, b, c).The vapor phase P-SVOCs reacts with OH radicals at a rate constant of 4 × 10 -11 cm 3 molecule -1 s -1 and the oxidation products are represented as a mixture of P-SVOCs and SVOCs in the next lower volatility bins (Table S1).
IVOCs refer to organic compounds which have saturation concentration between 10 3 and 10 6 µg m -3 , and are usually missing in most of the emission inventories.Recent studies indicate that POAs emissions are accompanied by the emissions of IVOCs, and the additional IVOCs emission is assumed to be 1.5 times of the POAs emission (Robinson et al., 2007).In this study, IVOCs emission was put into the bin of 10 4 µg m -3 saturation concentrations, and the rate constant of IVOCs with OH radicals is 4 × 10 -11 cm 3 molecule -1 s -1 .SOA mass yields for the VOCs and IVOCs precursors and volatility bins used in this box model listed in Table S2.
The equilibrium gas-particle partitioning of an organic i with effective gas-phase saturation mass concentration C i * (µg m -3 ) over an SOA solution is: where the partitioning coefficient ξ i is the mass ratio of particle-phase i to its total mass concentration (in both gas and particle phases), and C OA is the total mass concentration of the absorbing organic material (Donahue et al., 2006).May et al. (2003a, b, c).
For POAs and IVOCs emissions, we used the enthalpy of vaporization (ΔH vap ) values from May et al. (2013) for biomass burning and Ranjan et al. (2012) for anthropogenic emissions.This method has been used in previous studies (Murphy and Pandis, 2009).POAs includes 5-volatility bins ranging from 10 -1 to 10 3 µg m -3 saturation mass concentration (C i * ), which roughly covers the volatility range of semivolatile organic compounds (SVOCs).IVOCs are put into the sixth bin with 10000 µg m -3 C i * .
In this study, the method constraining box model by observations is taken from previous studies like Kanaya et al. (2009) and Li et al. (2011 and2015).The observed O 3 , CO, SO 2 , NO, NO 2 , HONO, HNO 3 , VOCs, temperature, Boundary layer (BL) height were averaged or interpolated with a time resolution of 5 min as constraints of the model.We assumed that the values of these observed parameters in this 5-time step were constant and input the box model.The 5-min averaged photolysis rated (J-value) was calculated by an accurate radiative transfer model (TUV4.5),which has been validated on atmospheric oxidation capacity in Beijing (Li et al., 2011).In details, the BL height was derived from the NIES two-wavelength aerosol Mie-lidar (532 and 1064 nm).The 5-min NO 2 observations were measured by an Aerodyne CAPS NO 2 monitor (Ge et al., 2013).This instrument determines NO 2 by directly measuring optical absorption of NO 2 at 450 nm in the blue region of electromagnetic spectrum and avoids overestimation of NO 2 measured by the photolysis-chemiluminescence monitors (Alpert et al., 1997).The hourly NO observation was from a commercial instrument (Model 42CTL, Thermo Electron Co.).As precursors of OH, concentrations of gaseous HNO 3 and HONO were measured using a Gas and Aerosol Collector with an Ion-exchange Chromatography (GAC-IC) system designed by Dong et al. (2012).CBMZ on gas-chemistry employed 133 reactions for 53 species and has been described by Zaveri et al. (1999).VOCs from GC-MS/FID (Wang et al., 2014) was mapped into CBMZ lumped mechanism according to suggestions by Yarwood et al. (2005) in the CAMx (http://www.camx.com/)user manual as shown in Table 2. POAs and SOA fractions were determined based on positive matrix factorization analysis of Aerosol Chemical Speciation Monitor data and are taken from Sun et al. (2016a, b).The hydrocarbon-like OA from heavy-duty vehicle emission, cooking OA, and biomass burning were considered as POAs which were used to constrain the box model.The emission rates of VOCs, NO x , and POAs were obtained from Zhang et al. (2009).
In physical processes, we only included the emissions and dry/wet deposition by assuming to be uniformly mixed in the mixing layer.This scheme and assumption have been widely used in the observation-based box-model (Zhang and Carmichael, 1999;Kanaya et al., 2009;Li et al., 2011Li et al., , 2015)).The lidar observations showed that pollutants have been mixed in the mixing layer.The observed wind velocity was only 1-2 m s -1 in the study period, which suggested that horizontal transport likely had no significant impact on concentrations in 5-min time step.
The initial time was set to 00:00 local standard time (LST), and each 1-day calculation was repeated three times to spin up the diurnal variations in the concentration of

Observations of VOCs
Fig. 1 presents the observed hourly benzene, toluene, xylene, and isoprene concentrations in Beijing.In general, these VOC species varied diurnally, with a daytime minimum and a midnight maximum.The peak values of benzene, toluene, xylene, and isoprene were 2.5, 6.3, 3.4, and 0.24 ppb, respectively, at 04:00 on November 8, 2014.The decreased concentration at 12:00 was caused by a deeper afternoon boundary layer and enhanced photochemical reactions, which converted more VOCs into semi-VOCs.Compared to other observations, the concentration of single ring aromatics in Beijing is more than twice that measured in Los Angeles and Paris (Borbon et al., 2013), and the concentration of isoprene in Beijing was usually below 0.2 ppb, which is only less than one tenth of that in southeastern United States (Xiong et al., 2015).This difference indicates that anthropogenic VOCs may be a more important role in SOA formation in Beijing than in other cities.Fig. 2 presents the observed concentrations of POAs, SOA, NO x and BL height at Beijing in the study period.Observed POAs in Beijing reached the maximum at 00:00-02:00 and the minimum around 12:00.This POAs diurnal variation presents an anticorrelation with BL height.At 23:00 on November 8, BL height was only 400 m and POAs concentrations even exceeded 50 µg m -3 at this moment.

Comparison between Simulations and Observations
OH radicals are key oxidants that determine the transformation rate from gaseous precursors to SOA.Unfortunately, OH was rarely observed in Beijing because  it was very hard.In this study, we employed the calculated OH value constrained by its mostly precursors' observations (NO x , VOCs, HCHO, HONO, O 3 ).Similar calculation was presented by Kanaya et al. (2009) in Mt.Tai, and Liu et al (2012) in Beijing, who found that this didn't bring significant errors.Although there were few OH observations in Beijing, we collected OH "observations" in August 2007 (Liu et al., 2012) and October 13-17, 2000 (Ren et al., 2004), in Beijing.Fig. 3 showed the comparison between observed OH in August, October and calculated OH in November in this study.We found that simulated OH was a little less than observation in October, because the driving energy of OH production (UV-B radiation) decreased from October to November rather than the model itself.The same reasons resulted in the 80% decrease of observed OH from August to October.In conclusion, the simulated OH reasonably reproduced the OH in the study period.
IVOC concentrations are rarely measured worldwide.Nevertheless, we compared the ratio of simulated IVOCs to single-ring aromatics in Beijing with a few observations in other countries.In this study, the mean simulated IVOC mass concentration was 48.3 µg m -3 , which was nearly two times that of the measured concentration of single-ring aromatics.The IVOC/single-ring aromatics ratio is comparable and slightly higher than that observed in California (1.0-1.2) during 2010 (Zhao et al., 2014).

SOA Formation Rates in Beijing
Fig. 4(a) shows the mean contribution to the net SOA production from three formation pathways in base simulation during daytime.In this 2-day simulation, the total daytime (08:00-17:00 LST) SOA net production was 60.6 µg m -3 .SOA production by the oxidation of IVOCs (SOA IV ) was the dominant pathway, which accounted for 82% (49.7 µg m -3 )  of the total SOA concentration.The magnitude of this contribution is similar to that of the underestimations reported in a previous study (Lin et al., 2015).The traditional SOA formation pathway from the oxidation of aromatics (SOA V ) only contributed 14% (8.5 µg m -3 ), and the aging of the gas-phase POAs (SOA P ) only accounted for 4% (2.4 µg m -3 ).The hourly net SOA production had a single-peak distribution pattern as shown in Fig. 4(b).After 08:00, the hourly SOA formation increased rapidly to a maximum of 5.5 µg m -3 at 11:00, which subsequently decreased to a minimum of 0.5 µg m -3 at night (not shown in Fig. 4(b)).The hourly SOA V , SOA IV , and SOA P levels ranged from 0.1 to 0.9, 0.3 to 5.3, and < 0.1 to 0.3 µg m -3 , respectively.The average SOA levels during the early afternoon (12:00-15:00 LST) were > 3.3 µg m -3 h -1 , which is four times of those reported in California (Zhao et al., 2014) because of the abundant SOA precursors present in Beijing.

Uncertainties in the Estimated SOA Formation Rates from Different Pathways
As shown in Table 3, 39 sensitivity simulations were conducted to evaluate the uncertainties of estimated SOA formation rates.As mentioned previously, the uncertainty in the POAs emission rate was ± 271% (Kurokawa et al., 2013).Two sensitivity simulations (Table 3, S01 and S02) High_vol scenario, and IVOCs/POAs = 3.0 for transport emission S20 High_vol scenario, and IVOCs/POAs = 3.5 for transport emission S21 High_vol scenario, and IVOCs/POAs = 4.0 for transport emission S22 Low_vol scenario, and IVOCs/POAs = 2.0 for transport emission S23 Low_vol scenario, and IVOCs/POAs = 2.5 for transport emission S24 Low_vol scenario, and IVOCs/POAs = 3.0 for transport emission S25 Low_vol scenario, and IVOCs/POAs = 3.5 for transport emission S26 Low_vol scenario, and IVOCs/POAs = 4.0 for transport emission S27 POAs emission is 2.5 times that of Base S28 POAs emission is 1.5 times that of Base S29 POAs emission is 2.5 times that of Base, High_vol scenario S30 POAs emission is 1.5 times that of Base, High_vol scenario S31 POAs emission is 2.5 times that of Base, Low_vol scenario S32 POAs emission is 1.5 times that of Base, Low_vol scenario S33 SOAs mass yields under high-NO x conditions multiply by 1.2 S34 anthropogenic SOA aging reaction rate use 1 × 10 -11 cm 3 molecule -1 s -1 S35 anthropogenic SOA aging reaction rate use 4 × 10 -11 cm 3 molecule -1 s -1 S36 anthropogenic SOA aging reaction rate use 4 × 10 -12 cm 3 molecule -1 s -1 S37 the concentration of OH radical multiplied by 2 S38 the concentration of OH radical multiplied by 3 S39 the concentration of OH radical multiplied by 5 were conducted to evaluate the effects of POAs emission rate uncertainty.The POAs emission rates in S01 and S02 were 200% and 50% those in the base case, respectively.The total SOA production in S01 and S02 was 69.8 and 55.0 µg m -3 , respectively (Figs. 5(a) and 5(b)).The SOA IV formation in S01 and S02 were 6.3 µg m -3 more than that in the base case (49.7 µg m -3 ) and 2.8 µg m -3 less than that in the base case, respectively.The change in the POAs emission rates not only affected the contributions from SOA P but also affected the total SOA formation through gas-particle partitioning.
POAs volatility distribution factors, as listed in the Base_vol scenario (Table 1), were used in a three-dimensional chemical transport model (Koo et al., 2014).May et al. (2013a, b, c) reported a wide range of values of POAs volatility distribution factors.In this study, two sensitivity simulations, in which the POAs have the highest (S03, High_vol scenario in Table 1) and the lowest (S04, Low_vol scenario in Table 1) volatilities, were conducted to evaluate the effects of POAs volatility on SOA formation (Table 3).The difference in the portion of gaseous POAs only slightly changed the SOA P levels (~1 µg m -3 ) in the base case, S03,and S04 (Figs. 5(c) and 5(d)).Uncertainties in the POAs volatility distribution (May et al., 2013a, b, c) exerted a higher effect on particle-phase POAs than SOA, and the total losses of particle-phase POAs mass in the base case, S03, and S04 were 38.4, 62.3, and 19.2 µg m -3 , respectively.The volatility distribution factors can strongly affect the total OA (including particle-phase SOA and POAs) mass as well as the SOA/OA ratio.In particular, the simulated SOA/OA ratio was largely underestimated in China relative to the measured SOA/OA ratio.Jiang et al. (2012) found that the simulated SOA/OA ratio in Beijing in summertime was underestimated to a certain degree of 40% compared with the observations.During the 2-day study period, the SOA/OA ratio in base case, S03, S04 and observation was 38.2%, 39.6%, 36.9% and 42.8%, respectively.
The IVOC emission rate has not been reported in previous bottom-up emission inventories.Robinson et al. (2007) suggested a IVOC/POAs emission ratio of almost 1.5 and applied it to three-dimensional modeling studies in Europe and America (Tsimpidi et al., 2010).However, a recent study showed that IVOC emissions could be as much as 4.5 times higher than POAs emissions (Woody et al., 2015).In this study, the higher IVOC emission simulation (Table 3, S05) and no IVOC emission simulation (Table 3, S06) were performed to estimate the extreme influence of IVOCs (Table 3).The total SOA formed in S05 was 139.9 µg m -3 (Fig. 5), which was more than two times that in the base case.Without contributions from SOA IV (S06), the total SOA production was only 8.2 µg m -3 .The comparison among base case (IVOCs/POA ratio of 1.5), S05 (change the IVOCs/POAs ratio to 4.5 only for transportation sources) and S06 (no IVOCs emission) shown large discrepancy of total SOA production in Beijing.Considering the impacts of IVOCs emission uncertainties on SOA formation, additional studies on IVOCs emissions are warranted to provide additional insights into SOA formation in China.
OH concentrations and the OH aging rates of SOA precursors (4 × 10 -11 ) were the uncertainty sources in this study.Six sensitivity simulations (S34-S39, Table 3) were conducted to examine the effects of these uncertainties.The different anthropogenic SOA aging rate by OH radicals lead to the differences in the total net production of SOA.However, the relative contribution from each pathway nearly the same (Fig. 6).Fig. 7 presents the boxplots of SOA formation from three pathways based on all sensitivity simulations (including S01-S39).The median value was regarded as the SOA production rate.In this study, the 25 th and 75 th percentiles represent the varying ranges of SOA formation rates, which account for the uncertainties related to emissions and model parameters.In general, the median SOA formation of S01-S39 during this 2-day study period was 73.8 µg m -3 , and the 25 th and 75 th percentiles was 55.4 µg m -3 and 102.4 µg m -3 , respectively.The SOA formation from three pathways (SOA V , SOA IV , and SOA P ), was 7.7-9.5 µg m -3 , 46.7-90.9µg m -3 and 1.0-2.7 µg m -3 (25 th and 75 th percentiles), respectively.

Implications for Three-dimensional SOA Simulations in China
In most current SOA three-dimensional simulations in China, SOA were formed only from VOCs, which were underestimated by approximately 70% compared with the observations in Beijing (Guo et al., 2014b;Lin et al., 2015).In particular, the multigeneration oxidation of VOCs was considered a possible method to improve model performance.Zhang et al. (2014) suggested that SOA mass yields under high NO x conditions may be underestimated by a factor of approximately 1.2 due to vapor losses.In this study, according to the measured VOCs, contributions from the multigeneration oxidation of VOCs were only 14%.Furthermore, corrections in mass yields only increased the predicted SOA V to 9.6 µg m -3 and did not alter the relative contribution of VOCs (15%) substantially, indicating that   The POAs mechanism was expected to largely explain discrepancies between observations and simulations and could contribute approximately 50% to the total SOA levels in Europe (Bergstrom et al., 2012); however, modeling studies in China have rarely considered this mechanism.In this study, SOA P accounted for only less than 5% of the total SOA levels in Beijing, which is much less than those reported in Europe and the United States (Hodzic et al., 2010).
In this study, IVOCs were the most important precursors of SOA.SOA production was remarkably increased through the IVOC mechanism in Beijing, which is expected to reduce the gap between simulations and observations.As mentioned previously, the SOA formed by IVOCs remain a large source of uncertainties (Fig. 7(b)), and the lack of observations in China hinder the selection of correct modeling parameters.The emission sources and reaction rates of IVOCs are required for further identification and quantification.

CONCLUSIONS AND RECOMMENDATIONS
We used an observation-constrained box model and related measurements to examine the current SOA formation in urban Beijing in autumn.In the base case, the total SOA yield during the 2-day study period was 60.6 µg m -3 , which is consistent with the yield reported for a previously recommended model configuration.The contributions from the oxidation of VOCs, IVOCs, and semi-volatile POAs (SOA V , SOA IV , and SOA P ) were 14%, 82%, and 4%, respectively.Considering the uncertainties in the emission rates and POAs volatility distribution factors, the median SOA net production of the sensitivity simulations in this study was 73.8 µg m -3 (55.4-102.4µg m -3 ).The amount of SOA formed via the three pathways (SOA V , SOA IV , and SOA P ), the SOA formation was 7.7-9.5 µg m -3 , 46.7-90.9µg m -3 , and 1.0-2.7 µg m -3 , respectively.IVOC oxidation was the major contributor to SOA production in Beijing during the study period.The present finding, which has not been reported for previous simulations in China, is likely to improve the predictive accuracy of three-dimensional SOA models in the nation.

Fig. 3 .
Fig. 3. Simulated diurnal OH concentrations in Nov. in Beijing.Also shown are observed OH in Aug. and Oct. in Beijing.The grey area shows the 20-day average diurnal variation of the OH radical with the hourly standard deviations.

Fig. 4 .
Fig. 4. (a) The averaged contribution from each pathway during daytime (08:00-17:00) and (b) the hourly net SOA production rate from different pathways.Note that SOA V , SOA IV and SOA P represents the SOA from oxidation of VOCs, IVOCs and semi-volatile POA.

Fig. 7 .
Fig. 7. Boxplots of the hourly net SOA formation production rates (µg m -3 hr -1 ) from different pathways based on sensitivity simulations including S01-S39 (red line, blue box, whiskers, and markers beyond the range of whiskers represent the median value, 25 th and 75 th percentiles, the 1.5 interquartile range, and outliers, respectively).The green solid points represent the base case.(a) SOA from oxidation of VOCs, (b) SOA from oxidation of IVOCs and (c) SOA from oxidation of semi-volatile POAs.

Table 1 .
Volatility distribution factors of POAs emissions.
b C * = 0 means nonvolatile.c Others mean other anthropogenic emission sources.d the distribution factors of POAs emissions are taken from

Table 2 .
Assignment from VOC Species to CBMZ Model Species.