Uncertainties in SOA simulations due to meteorological uncertainties in Mexico City during MILAGRO-2006 field campaign

The purpose of the present study is to investigate the uncertainties in simulating secondary organic aerosol (SOA) in Mexico City metropolitan area (MCMA) due to meteorological initial uncertainties using the WRF-CHEM model through ensemble simulations. The simulated periods (24 and 29 March 2006) represent two typical meteorological episodes (“Convection-South” and “Convection-North”, respectively) in the Mexico City basin during the MILAGRO2006 field campaign. The organic aerosols are simulated using a non-traditional SOA model including the volatility basis-set modeling method and the contributions from glyoxal and methylglyoxal. Model results demonstrate that uncertainties in meteorological initial conditions have significant impacts on SOA simulations, including the peak time concentrations, the horizontal distributions, and the temporal variations. The ensemble spread of the simulated peak SOA at T0 can reach up to 4.0 μg m −3 during the daytime, which is around 35 % of the ensemble mean. Both the basin wide wind speed and the convergence area affect the magnitude and the location of the simulated SOA concentrations inside the Mexico City basin. The wind speed, especially during the previous midnight and the following early morning, influences the magnitude of the peak SOA concentration through ventilation. The surface horizontal convergence zone generally determines the area with high SOA concentrations. The magnitude of the ensemble spreads may vary with different meteorological episodes but the ratio of the ensemble spread to mean does not change significantly.


Introduction
Meteorological, emissions, and air quality models are the key components of photochemical air quality simulation models (PAQSM).Uncertainties associated with PAQSM are varied and complex and they interact both within and across models (Fine et al., 2003).Meteorological condition simulation is critical for understanding the formation, transformation, diffusion, transport, and removal of the pollutants.Dabberdt et al. (2004) have listed the meteorological research needs for improved air quality forecasting, one of which is to provide model uncertainty information through ensemble prediction capabilities and quantify uncertainties and feedbacks between meteorological and air quality modeling components.Past studies on photochemical sensitivity to meteorological uncertainty mainly include Monte Carlo simulations (Hanna et al., 2001;Beekmann and Derognat, 2003;Irwin et al., 1987;Stuart et al., 1996;Bergin et al., 1999;Dabberdt and Miller, 2000) and adjoint sensitivity studies (Menut, 2003).The ensemble approaches have also been utilized in photochemical modeling by using different models (Galmarini et al., 2004a, b;McKeen et al., 2005), photochemical reactions (Delle Monache and Stull, 2003), emission scenarios (Delle Monache et al., 2006), and physical parameterizations (Mallet and Sportisse, 2006).In general, the ensemble means performed better than most individual models.
N. Bei et al.: Uncertainties in SOA simulations due to meteorological uncertainties Zhang et al. (2007a) have showed large uncertainties in the ozone (O 3 ) prediction in Houston and surrounding areas due to meteorological initial uncertainties through both meteorological and photochemical ensemble forecasts.Bei et al. (2010) have further investigated ozone predictabilities due to meteorological uncertainties in Mexico City Basin using ensemble forecasts.They found that the largest unpredictability in O 3 simulations was attributed to the increasing uncertainties in meteorological fields during peak O 3 period, and the impacts of wind speeds and PBL height on O 3 simulations are more straightforward.The ensemble spreads of simulated O 3 also vary with different PBL schemes and meteorological episodes.These works have demonstrated the importance of accurate representation of meteorological conditions in the air pollution studies in urban areas.
Atmospheric aerosols pose serious health risks and exert an important radiative forcing on climate.Organic aerosols (OA), accounting for 20-90 % of the total fine particulate mass in the atmosphere (Zhang et al., 2007b), comprise primary OA (POA) that is directly emitted into the atmosphere in particulate form, and secondary OA (SOA) which is formed from chemically processed gaseous organic precursors.Recent field studies have shown that the traditional semi-empirical 2-product parameterization significantly underestimates the measured SOA mass concentrations in urban and remote regions (e.g., de Gouw et al., 2009;Zhang et al., 2006).New SOA formation mechanisms have been suggested to close the gap in SOA mass concentrations between measurements and models, including the update of aromatic SOA yields (Ng et al., 2007), the SOA formation from dicarbonyl compounds (e.g., Zhao et al., 2006;Volkamer et al., 2007), and the formation of SOA from primary semivolatile and intermediate volatility species (Robinson et al., 2007;Grieshop et al., 2009).
Numerous model studies have evaluated these SOA formation mechanisms, with considerably large uncertainties compared with the field measurements.Tsimpidi et al. (2010) have evaluated the effects of the semi-volatile nature of primary organic emissions and photochemical aging of primary and secondary organics on OA levels during MCMA-2003(Molina et al., 2007) using a modified 3-D chemical transport model (CTM); on average, the model overestimates considerably the observed SOA concentrations during daytime.During the MILAGRO-2006 field campaign (Molina et al., 2010), several model studies have shown that the mechanism of Robinson et al. (2007) still underestimates the SOA observation in the urban area of Mexico City (Hodzic et al., 2010;Li et al., 2011a;Tsimpidi et al., 2011;Shrivastava et al., 2011).In addition, in the modeling study of Hodzic et al. (2010), the mechanism of Grieshop et al. (2009) significantly overestimates the SOA observations in Mexico City.Li et al. (2011a) have improved the simulation by including the contributions from dicarbonyl compounds, but the model results still fail to close the gap between the measurements and the model in Mexico City.Although the evalua-tion of the new SOA formation mechanisms using CTMs is considerably influenced by the uncertainties from measurements, emissions, aging of semi-volatile and intermediate volatile organic compounds, and contributions from background transport (Li et al., 2011a), few studies have considered the key role of meteorological conditions in the assessment of the SOA mechanism, especially when the measurements of SOA are confined in one or several supersites, which potentially constitutes one of the largest uncertainties in the evaluation of the SOA formation mechanism.
The purpose of this study is to investigate the uncertainties in simulating SOA in the Mexico City basin due to meteorological initial uncertainties based on the measurements obtained during MILAGRO-2006 field campaign (Molina et al., 2010).The impacts of meteorological uncertainties on SOA simulations are investigated through ensemble simulations using state-of-the-art meteorological and photochemical prediction models for two selected days (24 and 29 March 2006), which represent two of the typical meteorological episodes "Convection-South", and "Convection-North" in O 3 predictions in the Mexico City basin during MILAGRO-2006(de Foy et al., 2008).The methodology and experimental designs are presented in Sect. 2. The synoptic situations of the selected days are overviewed in Sect.3. The control ensemble forecasts are introduced in Sect. 4. The ensemble simulations on other day and the ensemble forecasts with different initialization method are presented in Sects.5 and 6, respectively; the summary and conclusions are given in Sect.7.

Methodology and experimental descriptions
The Advanced Research WRF (ARW) v3.2 (Skamarock et al., 2008) is used in meteorological deterministic and ensemble forecasts.The model simulations adopt horizontal resolution of 12 km and 35 sigma levels in the vertical direction with the grid size of 259 × 160 (Fig. 1a).The WRF model is initialized at 00:00 UTC and integrated for 30 h for all the selected days.The National Centers for Environmental Prediction final operational global gridded analysis (NCEP-FNL) is used to produce the initial and boundary conditions for the reference deterministic forecast.The physical process parameterization schemes used in the reference deterministic forecasts include the Grell-Devenyi ensemble scheme for cumulus scheme (Grell and Devenyi, 2002), the WRF Single Moment (WSM) three-class microphysics (Hong et al., 2004), and Mellor-Yamada-Janjic (MYJ) TKE scheme (Janjic, 2002) for PBL processes.The physical process parameterizations, particularly the PBL parameterization, play an important role in the air quality simulation.We have performed sensitivity studies to investigate the impact of different PBL schemes on ozone and aerosol simulations and found that the MYJ TKE PBL scheme yields more reasonable results than the other PBL schemes in the WRF model  compared to the observations.However, it is worth mentioning that the MYJ TKE PBL scheme is appropriate in the simulations in Mexico City, but might not work well in other megacities due to different meteorological situations, topography, land use, etc.
The ensemble initialization method is similar to the one employed in our previous study on ozone predictability due to meteorological uncertainty (Bei et al., 2010).The initial ensemble is generated with the WRF-3DVAR (Barker et al., 2004) using Background Error Statistics (BES) option cv5.A set of random control vectors with a normal distribution was generated.A control increment vector is then transformed back to model space via an empirical orthogonal functions (EOF) transform, a recursive filter, and physical transformation via balance equation.The perturbed variables include the horizontal wind components, potential temperature, perturbation pressure, and mixing ratio of water vapor, whose error statistics are defined by the domain specific climatological background error covariance that are derived from one-month simulations in the same domain using the NMC method (Parrish and Derber, 1992).Other prognostic variables such as vertical velocity (w) and mixing ratios of cloud water (q c ), rainwater (q r ), snow (q s ) and graupel (q g ) are not perturbed.The perturbations generated through this method are random and balanced noises, and their magnitudes are also small compared to the typical sounding observational and analysis errors (Nielsen-Gammon et al., 2007).The boundary conditions are perturbed in the same manner as the initial ensemble.Figure 2a-d show the vertical distribution of the initial ensemble spread, the average spread is 0.7-1.7 m s −1 for horizontal winds (u, v), 0.4-1.7 K for temperature (T ), 0-0.75 hPa for pressure (p), and 0-1.9 kg kg −1 for the water vapor mixing ratio (q).The 12-km meteorological ensemble simulations are then used to drive a 30-member 3-km photochemical ensemble simulation (97 × 97 grids) using the WRF-CHEM model.
We have also investigated the ensemble simulations using a "climatological ensemble initialization method" (hereafter referred to as "climatological method") in which dynamically consistent initial and boundary conditions are statistically sampled from a seasonal meteorological data set (Aksoy et al., 2005(Aksoy et al., , 2006;;Zhang et al., 2007a).To represent the springtime climatological statistics, a data set for the period of 1 February to 15 May 2006 is generated from NCEP-FNL 1 • × 1 • reanalysis data.30 ensemble perturbations are randomly selected from this climatological data set.Similarly, boundary conditions for each ensemble member are generated from the same data set beginning at the randomly selected initial time of the given member, and extending for the same length of time as the control run.Deviations of the initial and boundary condition data for each member from the climatological mean for the entire period are then scaled down to 20 % to reduce the ensemble spread to below typical observation error magnitudes (Kalnay, 2003) and added to the unperturbed initial and boundary conditions derived directly from the NCEP-FNL analyses valid at 00:00 UTC, 29 March, which are used for the 12-km domain ensemble simulation.Figure 2e-f show the vertical distribution of the initial ensemble spread, the average spread is 0.4-2.2m s −1 for horizontal winds (u, v), 0.5-0.65 K for temperature (T ), 0-0.17 hPa for pressure (p), and 0-0.45 g kg −1 for the water vapor mixing ratio (q).Except that the initial ensemble spread of u component is slightly bigger than the WRF-3DVAR method (see Fig. 2a-d), the initial ensemble spreads of other variables are smaller than their typical magnitudes of observation error.
In the present study, a specific version of the WRF-CHEM model (Grell et al., 2005) is used for photochemical ensemble simulations.The version of the WRF-CHEM model is developed by Li et al. (2010Li et al. ( , 2011aLi et al. ( , b, 2012) ) at the Molina Center for Energy and the Environment, with a new flexible gas phase chemical module which can be utilized in different chemical mechanisms, including CBIV, RADM2, and  SAPRC.The gas-phase chemistry differential equations are solved by an Eulerian backward Gauss-Seidel iterative technique.The short-lived species, such as OH and O( 1 D), are assumed to be in steady state.The solution is iterated until all species are within 0.1 % of their previous iterative values.For the aerosol simulations, the CMAQ (version 4.6) aerosol module developed by EPA, which is designed to be an efficient and economical depiction of aerosol dynamics in the atmosphere, is incorporated in the WRF-CHEM model (Binkowski and Roselle, 2003).The photolysis rates are calculated using the FTUV (Tie et al., 2003;Li et al., 2005).
The secondary organic aerosol (SOA) formation is simulated using a non-traditional SOA model including the volatility basis-set modeling method in which primary organic components are assumed to be semi-volatile and photochemically reactive and are distributed in logarithmically spaced volatility bins (Li et al., 2011a).The partitioning of semi-volatile organic species is calculated using the algorithm suggested by Koo et al. (2003), in which the bulk gas and particle phases are in equilibrium and all condensable organics form a pseudo-ideal solution (Odum et al., 1996).Nine surrogate species with saturation concentrations from 10 −2 to 10 6 µg m −3 at room temperature are used for the primary organic aerosol (POA) components following the approach of Shrivastava et al. (2008).The SOA formation from each anthropogenic or biogenic precursor is predicted using four semi-volatile organic compounds whose effective saturation concentrations at 298 K are 1, 10, 100, and 1000 µg m −3 , respectively.The NO x -dependent SOA yields from anthropogenic and biogenic precursors are included (Lane et al., 2008), and the oxidation hypothesis of semivolatile and intermediate volatile organic compounds by Grieshop et al. (2009) is used.The contributions of glyoxal and methylglyoxal are also considered in the study.Detailed description about the volatility basis-set approach can be found in Li et al. (2011a).
The emission inventory used in this study is developed at the Molina Center by Lei et al. (2012), including fossil fuel combustion (mobile, area and point sources) and open burning of biomass and trash.The biogenic emissions are calculated using the on-line MEGAN model (Model of Emissions of Gases and Aerosols from Nature) developed by Guenther et al. (2006), in order to consider the variations of biogenic emissions due to the temperature change in ensemble simulations.The method proposed by Tsimpidi et al. (2010) is applied to redistribute the POA emissions.The chemical initial and boundary conditions for the WRF-CHEM model simulations are interpolated from MOZART 3-h output (Horowitz et al., 2003).Considering that we mainly concentrate on the effects caused by changes in the meteorological fields, the initial and boundary conditions for chemical fields and the emission inventory are the same for all ensemble experiments.
Both meteorological and photochemical ensemble simulations are conducted on two selected days (24 and 29 March 2006).We choose 29 March as a control ensemble run (CTRL), and a detailed analysis is presented on this day.The physical process parameterization schemes used in the CTRL run are the same as those used in the reference deterministic forecast.
The ensemble simulation results are compared to the Aerosol Mass Spectrometry (AMS) observations analyzed using the Positive Matrix Factorization (PMF) technique at an urban background site (T0) and a suburban background site (T1) in Mexico City.The T0 monitoring station is located in the northwestern part of the basin of Mexico City, influenced by road traffic emissions (300 m from four major roads surrounding it), domestic and residential emissions, and also potentially influenced by local industrial emissions and from the Tula industrial area (60 km to the north-northwest, in the Hidalgo State).T1 supersite is located around 50 km to the north of Mexico City, in an area isolated from major urban agglomerations but close to small populated agglomerations, and around 500 m from the closest road.

Synoptic overview
The two days selected (29 and 24 March) in the study represent two different meteorological episode types in Mexico City, which are defined as "Convection-North" and "Convection-South" in de Foy et al. (2008).29 March is classified as "Convection-North", which represents northerly wind aloft and rain in the northern part of the basin.24 March is classified as "Convection-South", which represents southerly wind aloft and rain in the southern part of the basin.At 500 hPa, the dominant winds over the Mexico City basin are northerly (Fig. 3a) and southerly (Fig. 3b) for these two days, respectively.At 700 hPa (Fig. 3c-d), the anti-cyclones on both days lead to subsidence over Mexico City basin.At the surface (Fig. 3e-f), there are convergences in the Mexico City basin on both days but at different location, which leads to different location of the precipitation on these two days.In addition, on 29 March, the southerly wind is slightly stronger, causing the precipitation occurred in the northern part of the basin.On 24 March, the northerly wind is slightly stronger, causing the precipitation occurred in the southern part of the basin.
According to the flow type, "Convection-North" and "Convection-South" days can be classified into "O 3 -South" and "O 3 -North" episode types, respectively, except that the convective activity prevents the formation of a clean convergence zone sweeping through the basin in late afternoon (de Foy et al., 2005(de Foy et al., , 2008)).Considering that "O 3 -South" and "O 3 -North" episodes dominate the MILAGRO-2006 field campaign period (de Foy et al, 2008), the two days we have simulated represent most of the meteorological situations during the campaign.

Control ensemble simulations
A detailed analysis of the ensemble simulations on 29 March is presented in this section; the other ensemble simulation on 24 March will be presented in the next section.[POA] against observations at T0 on 29 March for the ensemble mean, best, minimal, and maximal member, and the reference deterministic forecast.In Fig. 5, overall, the ensemble mean and best member exhibit better performance than the reference deterministic forecast.Further analysis on this best member will be presented in the next subsection.

Overview of the control ensemble performance
We have also compared the ensemble simulations and the observation at T1 site, located in the northwest of the Mexico City basin and used as a suburban background site during MILAGRO (Fig. 6).The ensemble means fail to capture the peak time concentrations (including both timing and magnitude) in the morning at T1, but still performs better than the reference simulation.One of the possible reasons for the deviation between measurements and the model at T1 is that T1 is in the downwind of the MCMA, and is frequently influenced by the transport from MCMA determined generally by local meteorological conditions.Therefore, the meteorological uncertainties cause more sensitivities of the simulated [SOA] and [POA] at T1 than at T0. Figure 7 shows the evolution of the ensemble mean of the surface [SOA] distributions along with the ensemble mean wind vectors in the MCMA and the surrounding area simulated by the CTRL ensemble.From 00:00 to 06:00 CDT (Fig. 6a), the ensemble mean of the [SOA] is low within 37 Figure 6 Fig. 6.Same as Fig. 4, but for T1 site.
the urban area of the Mexico City basin due to the lack of photochemical activities.After 06:00 CDT, the [SOA] inside the basin increases due to the nearly calm wind inside the basin and the increasing photochemical activities.The high ensemble mean [SOA] first occurs around 09:00 CDT (Fig. 6b) as observed, with the maximum concentration area located nearly in the southwest of the basin.Along with the development of the mixing layer height and the increase of the surface divergence due to the topography, the [SOA] decreases from 09:00 to 12:00 CDT (Fig. 6c).From 12:00 to 15:00 CDT (Fig. 6d), in association with the increase of northerly wind, gap wind in the southeast and the downhill wind in the west, south and east edge of the basin, a convergence zone is formed in the southwest of the Mexico City basin, leading to the increase of the simulated [SOA] and the second occurrence of the maximum ensemble mean [SOA] inside the basin around 17:00 CDT (not shown), which is 2-h later than the observation from T0. From 18:00 to 21:00 CDT (Fig. 6e-f), the high ensemble mean [SOA] area moves southward along with the increased northwesterly winds inside the basin and decreases again due to the weakened convergence and photochemical activities.

Uncertainties in SOA and POA simulations
Although the initial meteorological uncertainties are smaller than typical observational and analysis errors, our control ensemble simulations demonstrate that large uncertainties still exist in [SOA] and [POA] simulations, especially during the peak time periods (see Fig. 4).
To illustrate the discrepancy between different ensemble members, we have chosen two ensemble members: EN-34 and EN-30, which represent the highest and lowest [SOA] at T0, respectively.Figure 8   distributions in Mexico City basin are principally attributed to the striking discrepancies in the surface winds between these two extreme members.Before the peak [SOA] time (Fig. 7a-b), the weaker southerly surface winds in the northern basin transports less precursors outside of the Mexico City basin in EN-34, while the stronger southerly surface winds in the northern basin transports more precursors outside of the Mexico City basin in EN-30.At the peak time (Fig. 7c-d), the winds from south, west, and northwest around the basin are all stronger in EN-34 but weaker in EN-30, thus more [SOA] and its precursors accumulate inside the basin in EN-34 but less [SOA] and its precursors accumulate inside the basin in EN-30.After the peak time (Fig. 7e-f), the [SOA] is higher in EN-34 due to the stronger convergence inside the basin in EN-34, while the [SOA] is lower in EN-30 due to the weaker convergence inside the basin in EN-30.The high [SOA] area is basically consistent with the convergence zone.Figure 4 shows that the ensemble mean [SOA] is consistently better than the reference deterministic simulation.We have also identified the best member (EN-14) that agrees well with the observed [SOA] at T0, including the amplitude and timing.Figure 9 shows the evolution of the simulated surface [SOA] distributions from the reference forecast and the best member (EN-14), demonstrating the large differences of the location and movement of high [SOA] area between the above two simulations.At the first peak time (Fig. 8a-b: 10:00 CDT), due to the stronger uphill wind along the south edge of the basin, the simulated high [SOA] area is located along the west and south edge of the basin in the reference simulation.While the simulated high [SOA] of the best member is located in the southwest of the basin because of the stronger uphill wind along the west and the east edge, and the weaker uphill wind along the south edge of the basin.At the second peak time (Fig. 8c-d: 15:00 CDT), the northwesterly wind in the northwest of the basin and the southerly wind along the south edge are both much stronger in the reference simulation but much weaker in the best member simulation, which cause the [SOA] plume to move to the south edge (southwest) of the basin in the reference simulation but to the southwest of the basin in the best member simulation.At the 17:00 CDT (Fig. 8e-f), the high [SOA] areas in reference simulation and best member are generally consistent with the location of the convergence line (area) inside the basin.The differences in SOA simulations are mainly attributed to the differences in the simulated horizontal convergence distributions.The discrepancies in the basin scale wind and the related SOA simulations between our reference deterministic run and the best member are still manifest, indicating the importance of considering the uncertainties in current photochemical simulations due to the meteorological fields.

Potential impact of meteorological fields on chemical simulations
We have discussed the role of meteorological conditions in photochemical simulation in our previous studies (Bei et al., 2008(Bei et al., , 2010) ) and in the previous subsections in the present study.Here, we attempt to summarize the role of the wind speed and the horizontal convergence in the basin wide area in simulating SOA and POA. Figure 10 shows the time evolution of the basin-wide domain averaged (hereafter domainaveraged, averaged domain shown in Fig. 1b) surface wind speed and convergence.During 06:00 to 09:00 CDT, the domain-averaged convergence/divergence is relatively small but the domain-averaged wind speed is the lowest on that day, which are both favorable for the accumulation of the pollutants and the formation of the first peak [SOA] and [POA].During 12:00 to 15:00 CDT, the domain-averaged wind speed slightly increases but the convergence in the basin wide reaches its maximum value during the day.Therefore, during the two peak times of [SOA], the meteorological condition in the basin wide is favorable for the formation of high [SOA].In addition, from midnight to the morning time of the following day (00:00 to 12:00 CDT, Fig. 9b), two extreme members EN-30 and EN-34 have maximum and minimum wind speed in the basin wide, respectively, showing that the ventilation inside the basin from the midnight to the following morning is also important to the formation of the first [SOA] peak.During the two peak times (09:00 and 17:00 CDT) of the simulated [SOA], two extreme members EN-30 and EN-34 also have maximum divergence and convergence in the basin wide (Fig. 9a), respectively, indicating that the convergence is also crucial to the formations of both [SOA] peak.

Ensemble simulations on 24 March 2006
In order to explore the impacts of meteorological uncertainties on SOA predictability under different meteorological conditions, we have further conducted ensemble forecasts on 24 March that represent the "Convection-South" episode in Mexico City basin (de Foy et al., 2008), using the same models and ensemble initialization method as the control ensemble simulation.port from the source region, which are more sensitive to the meteorological uncertainties.
Overall, the uncertainties in the SOA and POA simulations due to the initial meteorological errors are comparable in magnitude and also significant to the ensemble mean on these two days.

Ensemble simulations with other initialization method
Initial perturbation is a key problem in ensemble forecasting.Saito et al. (2011) compared the five initial perturbation methods for the mesoscale ensemble prediction for the Beijing 2008 Olympics Research and Development Project and showed the considerable impact of different initial perturbation methods on the mesoscale ensemble prediction.We have also investigated the ensemble simulations on 29 March using another ensemble initialization method named as "climatological method" (see Sect. 2).
Figure 13 shows the temporal evolution of the ensemble mean and spread of the surface [POA] and [SOA] along with the reference deterministic forecast and the observations at T0 site.The ensemble mean [SOA] is slightly higher than that of WRF-3DVAR method (Fig. 4a) during the peak time, but the time evolution pattern are the same as WRF-3DVAR method.Except some differences in individual members, the ensemble mean [POA] are very similar to the results of WRF-3DVAR method (Fig. 4b).The above two initialization methods have very similar results, indicating that there is not much difference between these two methods.

Conclusions and discussions
We have investigated the uncertainties in simulating SOA due to meteorological initial uncertainties using the WRF-CHEM model through ensemble simulations in Mexico City on two selected days (24 and 29 March 2006), which represent two different meteorological episodes ("Convection-South" and "Convection-North") in the Mexico City basin.We choose 29 March 2006 as a control run to provide a detailed analysis of the model results.The control initial ensemble is generated with the WRF-3DVAR system.
In the urban area of Mexico City, the ensemble means of the [SOA] and [POA] in the control run compare reasonably well with the observations, including the sharp buildup of [SOA] and [POA] in the morning and the second peak of the [SOA] in the early afternoon.The ensemble mean is generally lower than the observations but considerably better than the reference deterministic forecast.The ensemble spread of the simulated peak [SOA] and [POA] at T0 can reach up to 4.0 (1.6) µg m −3 during the daytime, respectively, which is around 35 % and 40 % of the ensemble mean.During the two peak times of [SOA], the meteorological condition is favorable for the accumulation of pollutants and the formation of high [SOA] in the basin wide, such as the low wind speed in the morning and the strong convergence in the early afternoon.In addition, the analysis of two extreme members shows that the ventilation inside the basin from the midnight to the following early morning is also important to the formation of the [SOA] peak in the morning.However, in the suburban area of Mexico City, the ensemble means still deviate appreciably from observations, and the meteorological uncertainties result in larger sensitivities in simulating the [SOA] and [POA] because of the dominant impact of meteorological fields on the downwind transport of the urban pollution.The uncertainties in the SOA and POA simulations due to the initial meteorological errors on 24 March 2006 are comparable in magnitude to those in the control run and also significant compared to the ensemble mean.
We have also demonstrated the uncertainties in SOA simulations using another ensemble initialization method, namely "climatological method".These two initialization methods yield very similar results, indicating that there is not much difference between these two methods.However, the ensemble mean is not yet consistent with the measurements in both urban and suburban area of Mexico City, showing that meteorological initial uncertainties can only partially explain the uncertainties in the SOA simulations.Other uncertainties, such as those from meteorological model, e.g., PBL schemes (Bei et al., 2010), emission (Li et al., 2011a), SOA formation mechanism (Li et al., 2011a), and photochemical models, should be considered also in the ensemble simulation system.Future studies need to be performed to further investigate the impacts of these uncertainties on SOA simulations.
The uncertainties from the SOA formation mechanisms, emissions, and aging of semi-volatile and intermediate volatile organic compounds also contribute to the uncertainties in the SOA simulations.However, in the study of Li et al. (2011a), the difference of SOA simulations between the mechanism of Robinson et al. (2007) and Grieshop et al. (2009) is about 1 µg m −3 , the SOA contributions from glyoxal and methylglyoxal do not exceed 0.8 µg m −3 , and the assumption of the continued chemical aging of the SVOCs produced from the oxidation of anthropogenic VOCs also increases the SOA formation by about 0.7 µg m −3 in the urban area.In addition, Li et al. (2011a) have suggested that the uncertainties of SOA formation from the emission inventory are not more than 0.5 µg m −3 in Mexico City based on the comparison between modeled precursors and measurements, and if the aging process of semi-volatile and intermediate volatile organic compounds does not have feedback on the OH in the gas-phase chemistry, the SOA production is enhanced by up to 1.0 µg m −3 .All the uncertainties of SOA formation in Li et al. (2011a) is much less than those caused by the meteorological initial uncertainties, which produce up to 4.0 µg m −3 ensemble spread in the urban area of Mexico City.Hence, meteorological uncertainties constitute one of the largest uncertainties in the evaluation of the SOA formation mechanism using CTMs based on the measurements confined at one or several supersites, and also provide a potentially reasonable explanation for the large difference in recent SOA model studies (e.g., Hodzic et al., 2010;Tsimpidi et al., 2011;Shrivastava et al., 2011;Li et al., 2011a).Furthermore, the meteorological ensemble is possibly an efficient method to reduce the meteorological uncertainties in simulations of CTMs.

Fig. 1 .
Fig. 1.(a) WRF domain (red box is the WRF-CHEM domain) and (b) WRF-CHEM domain (red box indicated in a) and the observation sites for aerosol measurements in MCMA (red dot: T0, blue dot: T1).Inner box indicates the domain shown in Figs.6-8.Contours in both panels represent the terrain height with the intervals of 500 m (top) and 200 m (bottom), respectively.

Figure 3 Fig. 3 .
Figure 4a-b show the temporal evolution of the ensemble mean and the spread of the surface SOA concentrations ([SOA]) and POA concentrations ([POA]) along with the reference deterministic forecast and the observations at T0 (location shown in Fig. 1b), a supersite located near the urban center of Mexico City.The ensemble mean captures reasonably well the sharp buildup of the [SOA] and [POA] in

35 Figure 4 Fig. 4 .Fig. 5 .
Figure 4 Fig. 4. The temporal evolution of (a) the surface POA and (b) SOA concentration from each ensemble member (thin green lines), ensemble mean (bold black line), reference deterministic forecast (bold blue line), and the best member (compare to observations, bold red line) of the CTRL ensemble simulations (29 March 2006) and observations (red dots) at T0. (c) and (d) denote the ratio of the ensemble spread and the ensemble mean (brown line) for SOA and POA, respectively.
presents the horizontal distributions of the surface [SOA] along with surface winds from EN-34 and EN-30.The large variations in[SOA]

39 Figure 8 Fig. 8 .
Figure 8 Fig. 8. Horizontal distributions of [SOA] and surface winds simulated by two extreme members (which have maximum and minimum [SOA] values) along with the measurements inside the Mexico City basin for 10:00 CDT, 15:00 CDT, and 17:00 CDT on 29 March 2006.

40 Figure 9 Fig. 9 .
Figure 9 Fig. 9. Same as Fig. 8, but for the reference deterministic run and the best member.

41Figure 10
Figure 10 Fig. 10.Temporal evolution of the Mexico basin domain-wide averaged divergence (a) and wind speed (b) for 29 March 2006.
Figure 11a-b show the temporal evolutions of the simulated ensemble mean and spread and its ratio of the [SOA] and [POA] at T0 along with the observations on 24 March.Apparently, the diurnal pattern of the simulated [SOA] on 24 March is different from that on 29 March, which is caused by the different synoptic conditions of the two days.The temporal evolutions of the ensemble means of the [SOA] and [POA] generally agree with the observation at T0, but the ensemble mean [SOA] are slightly underestimated during the daytime and the ensemble mean [POA] are overestimated during the morning hours.The ensemble means are much better than the reference deterministic forecasts, especially during the peak times (around 14:00 CDT for [SOA], around 10:00 CDT for [POA]).The temporal evolution of the ratio of the ensemble spread and the ensemble mean is different from that of 29 March.The maximum ratio value can reach up to 70 % during the high [SOA] build up period.