Effect of Global Atmospheric Datasets in Modeling Meteorology and Air Quality in the Andean Region of Ecuador

Several components such as initial (IC) and boundary conditions (BC), generated from global atmospheric datasets for studies at the regional and mesoscale levels, play a role in numerical modeling. The analysis datasets include operationally available observations during the time of global meteorology forecast whereas reanalysis products, though available later, include more observations. Thus, it is expected that the IC and BC generated by reanalysis will be of higher quality, thereby improving the modeling results. To generate the IC and BC of a domain with a high spatial resolution (1 km) that covers the city of Cuenca, an urban area located in the Andean region of Ecuador, we employed two analysis (GFS, FNL) and three reanalysis (NCEPR2, ERA-Interim, ERA5) products. We used the Eulerian Weather Research & Forecasting model with Chemistry (WRF-ChemV3.2) to simulate meteorological and air quality variables. FNL and GFS were best fit for modeling both meteorological and air quality variables. Likewise, it also suggests their use in generating the IC and BC for modeling purposes in the Andean region of Ecuador. Currently, a few observations from the equatorial Andean region are incorporated into global atmospheric datasets. Furthermore, the atmospheric processes in this region are particularly complex and have been less studied. These limitations, which appear to be more prevalent for reanalysis products, have an impact on the generation, and thus the quality of the information stored by the atmospheric datasets for the equatorial Andean zone. Reanalysis products are not always the best choice for modeling in this area.


INTRODUCTION
Atmospheric modeling involves a variety of components.For mesoscale and regional studies, numerical models work with resolutions of a few km, requiring the generation of initial (IC) and boundary conditions (BC) from global atmospheric datasets.
The Global Forecast System (GFS), produced by the National Centers for Environmental Prediction (NCEP) (NCEI, 2021), is one of these datasets.It is generated by assimilation of atmospheric observations available at the time of generating the IC for global weather forecasting using the GFS model.
Another NCEP product is the FNL Global Operational Analysis (Final) (NCEP/NWS/NOAA/DOC, 2000), which is generated using similar weather forecast model and assimilation scheme as the GFS.Despite FNL having a delay of 60 to 90 min, it incorporates approximately 10% more observations than GFS records.
When compared to analysis products, retrospective analysis or reanalysis atmospheric datasets are generated later, thus ingesting an additional amount of records that were not operationally available at the time of generating the meteorological forecast.As a result, reanalysis datasets incorporate observations and cover previous time periods (Parker, 2016).The National Centers

Emission Inventory of Cuenca
Cuenca has a population of approximately 640 000 people.Its terrain is complex, with altitudes ranging from 1000 to 4000 masl (Fig. 1).The most recent emission inventory was conducted in 2014 (EMOV EP, 2016), and the results showed that on-road traffic was the primary source of most atmospheric pollutants (94.9% carbon monoxide (CO), 71.2% nitrogen oxides (NOx), 42.4% fine particulate matter (PM2.5),39.6 % of non-methane volatile organic compounds (NMVOC)).Other relevant sources of NMVOC were the use of solvents (29.7%), and the vegetation (19.5%).The industrial sector was the most significant source (60.1%) of sulfur dioxide (SO2).A thermoelectric plant located NE of the urban area generated 35.1% of SO2 and 18.5% of NOx emissions.Furthermore, it was reported that several artisan brick producers located NW of the urban area generated 38.5% of total PM2.5 emissions.The spatial distribution of annual PM2.5 emissions is depicted in Fig. 1.Emission inventories are typically the components with the highest level of uncertainty among the elements that participate in air quality simulation.This feature is especially important in all of Ecuador's emissions inventories prepared to date (MAE, 2014a(MAE, , 2014b)).The high levels of uncertainty are primarily due to a lack of national emission factors, which have necessitated the use of emission factors from international literature.In the case of Cuenca, the quality of the corresponding emission maps can be tested in an atmospheric model to reproduce both the observed meteorology and air quality in order to develop a new inventory with less uncertainty.To accomplish this, the best atmospheric dataset for generating the IC and BC must be identified.

Air Quality Network from Cuenca
To monitor meteorology and air quality in 2014, the city had an automatic station in the historic center (MUN station, Fig. 1).In addition, the monitoring network has approximately 20 passive stations that measure the monthly mean concentrations of NO2 and O3.The availability of a unique automatic station in 2014 for monitoring the air quality in the urban area of Cuenca is consistent with the minimum requirement established by the European Directive about air quality and cleaner air for Europe (EU, 2008), which indicate to assess compliance the PM2.5 exposure, it is required at least one sampling point per million inhabitants.After 2014, the air quality network improved, and currently, it has six automatic stations for measuring PM2.5 and PM1 (Parra et al., 2022).The municipality of Cuenca operates the air quality network, which has been accredited by the competent national authority, in accordance with the methods established in Ecuadorian regulations.
Between 2008 and 2019, mean annual NO2 concentrations higher than the World Health Organization (WHO) 2006 guideline (40 µg m -3 ) (WHO, 2006) were measured in BCB and VEG (Fig. 1), two passive microscale stations located in street canyons with high traffic of gasoline and diesel vehicles (mainly buses) (Parra and Espinoza, 2020).
PM2.5 concentrations (mean value of 24 h) were higher than the WHO guideline of 2006 (25 µg m -3 ) for fourteen days between 2012 and 2019.In four years of this period, the PM2.5 annual mean concentration was higher than the WHO guideline value (10 µg m -3 ) proposed in 2006.
O3 concentrations may be higher than the WHO guideline value (100 g m -3 , maximum average of 8 h) on specific days in September, due in part to high levels of solar radiation arriving perpendicular to the Equator during this month (Parra and Espinoza, 2020), which promotes the photochemical production of O3.
The WHO recently updated its guidelines for air quality, recommending stricter concentrations of particulate matter, O3, NO2, SO2, and CO (WHO, 2021).In 2020 in Cuenca, the annual mean concentrations of NO2 (19.7 µg m -3 ) and PM2.5 (8.5 µg m -3 ) were higher than the new recommended values (10.0 µg m -3 for NO2, and 5.0 µg m -3 for PM2.5).Furthermore, higher daily mean concentrations of NO2 (25.0 µg m -3 ) and PM2.5 (15.0 µg m -3 ) were recorded over 91 and 49 days, respectively, at the new guide values.The records highlight the importance of better understanding the quality problems, improving the emissions inventory, and defining a suitable approach for modeling meteorology and air quality.

Objectives
Although the results of the simulations performed thus far have been promising, it is necessary to assess the influence of other components in order to select those that contribute to improving the modeling of meteorology and air quality in Cuenca.Given that atmospheric reanalysis datasets contain more observations than analysis products, the IC and BC generated by reanalysis are expected to improve modeled results.Thus, the primary goals of the contribution are meant to address the following questions: • Are reanalysis datasets better to model meteorology and air quality in the Andean region of Ecuador? • What atmospheric dataset provides the best IC and BC to model meteorology and air quality in this region?• What is the recommended atmospheric dataset for atmospheric modeling as a quality indicator activity for a new emissions inventory?

Modeling Approach
Using the Eulerian Weather Research & Forecasting with Chemistry model (WRF-ChemV3.2),we generated the IC and BC for modeling meteorological and air quality variables in Cuenca using GFS, FNL, NCEPR2, ERA-Interim, and ERA5.WRF is a non-hydrostatic 3-D model widely used for meteorological research and forecasting (NCAR, 2021), which allows modeling of chemical transport of air pollutants to be included at the same time.A master domain and three nested subdomains (Fig. 1) were used to run the meteorological simulations.The nested grid ratio between domains was 3 to 1.The internal subdomain is a grid of 100 × 82 cells, 1 km on each side, that covers the territory of Cuenca (Fig. 1(c)).We used 35 vertical levels, up to the top pressure of 50 hPa (about 20 km altitude), being comparable with the vertical resolution of ERA-Interim and ERA5 datasets (38 vertical levels) and more detailed than GFS (27), FNL ( 27), and NCEP2 (18) vertical resolution.The option of chemical transport of pollutants was activated for the internal subdomain, allowing the reading of hourly emissions previously elaborated from the 2014 emission inventory.Based on previous modeling studies in Cuenca (EMOV EP, 2016;Parra, 2018b), BC for chemical species were assumed, using 85% of the default concentrations from WRF-Chem V3.2.Due to the O3 levels, the study period run from 4 to 27 September 2014.Furthermore, the activity of on-road traffic and other emission sources is typical of other months.On-road traffic emissions for September were estimated from the total yearly emissions through the fuel consumption ratio during this month to the total annual amount.The daily traffic emissions were estimated based on the traffic intensities during weekdays and weekends.Hourly on-road traffic emissions were defined from typical daily profiles provided by the Municipality of Cuenca.NMVOC emissions from vegetation were estimated from modeled maps of temperature and solar radiation for each day of the study period.The Carbon Bond Mechanism -Z (CBMZ) (Zaveri and Peters, 1999) and the Model for Simulation Aerosol Chemistry and Interaction (MOSAIC) were used to speciate hourly emissions (Zaveri et al., 2008).The main characteristics of the atmospheric bases used in this study are shown in Table S1.Table S2 shows the main parameters and options chosen for the WRF-Chem model.The option to use the online approach, which considers the direct effects of aerosols on meteorology, was activated.
We made the simulations for three-day periods, using six h as spin-up time for each one, and activating the option for using the default WRF-ChemV3.2's idealized profile to initialize chemistry.
The study covers 24 days, with a proportion of 80% of the monthly period, which is consistent with the required proportion (at least 75%) when aggregating data and calculating statistical parameters, according to the European Directive about ambient air quality and cleaner air for Europe (EU, 2008).

Modeling Performance
The meteorology model's performance was defined using the MUN station's hourly records of temperature, wind speed, and wind direction.For temperature and wind speed, we use the following metrics: where GE is the gross error, MB is the mean bias, IOA is index of agreement, RMSE is the root mean square error, N is the number of values, Pm is the mean modeled value, Om is the mean observed value, Pi is the modeled value, and Oi is the record.Wind direction was assessed through GE and MB.Wind direction also was evaluated by comparing the modeled wind roses to the one obtained from the records.
To evaluate the performance of modeling short-term air quality, the CO (maximum mean in 1 h and 8 h per day), PM2.5 (mean in 24 h), and O3 (maximum mean in 8 h per day) records were used.These periods are consistent with national air quality legislation and the WHO guidelines Furthermore, we calculated the percentages of meteorological and air quality records captured by modeling.We consider that a record was "positive modeled or captured" if the maximum deviation between the observed and modeled value agreed with the accuracy proposed by the European modeling quality Directive (EEA, 2011) (Table S3).A percentage of 100% means that all records were captured.This approach was applied in previous studies to assess the modeling performance in Cuenca (Parra, 2021(Parra, , 2018b)).
For long-term air quality records (monthly means of NO2 and O3), modeling performance was defined as the percentage of passive stations with a maximum difference of 30%.The percentage of 100% indicates that the monthly mean concentrations in all passive stations that measure each pollutant were properly modeled.

Meteorology
The metrics indicate that GFS, FNL, and NCEPR2 provided the best fit to temperature measurements.As shown in Table 1, the GE indicator values were in the reference range for all datasets (reference: < 2°C, minimum: 1.3°C for FNL, maximum: 1.8°C for ERA5), whereas the MB indicator was within the ± 0.5°C boundaries for all databases, except for ERA5 and ERA-Interim.Although, on average, ERA-Interim overestimated (MB = 1.0°C) and ERA5 underestimated (MB = -0.7°C) the temperature, a strong relationship between measured and modeled values (IOA = 0.9) was established for all datasets.
When FNL, GFS, ERA-Interim, NCEPR2, and ERA5 were applied, 77.3%, 75.7%, 74.1%, 72.4%, and 64.3% of the hourly temperature records were captured, respectively (Table 3).In addition, the FNL produced the highest proportion of modeled values within the range of variation established for temperature (Fig. 2(b)), as well as matched the recorded daily temperature profile most closely (Fig. 2(f)).Although a higher MB was obtained for GFS (0.4°C) relative to FNL (0.2°C), the former model produced output that was most consistent with records.On the other hand, NCEPR2 tended to underestimate and overestimate temperatures during the early morning hours and at noon, respectively.As can be seen from Fig. 3, for 25 September 2014 (13:00 LT), the GFS and FNL modeled temperatures between 18°C and 21°C (shown in violet color) throughout the urban area, whereas values between 21°C and 24°C (red color) were computed using the reanalysis datasets.
According to the obtained metrics, all atmospheric datasets yielded acceptable modeled wind speed values.Even though all RMSE, MB, and IOA values are within their respective reference ranges, as can be seen from Table 1, the FNL produced the best fit to the data (RMSE = 0.9°C, MB = 0.2°C, IOA = 0.9), as it captured 77.6% of the recorded values (Table 3).In contrast, ERA5 captured 63.8% of the records.For all databases, RMSE (reference: < 2 m s -1 , minimum: 0.9 m s -1 for FNL, maximum: 1.3 m s -1 for ERA5) and MB (reference: < ± 0.5 m s -1 , minimum: 0.2 m s -1 for FNL, maximum: 0.3 m s -1 for the remaining datasets) were within the reference range.A strong linear relationship between the records and modeled wind speed values was obtained in all cases (IOA of 0.9 and 0.8 calculated for FNL and other models, respectively).According to the records, winds were primarily coming from the east-northeast (ENE) direction (Fig. 4(a)).However, GFS, FNL, NCEPR2, and ERA5 indicated east (E), while only the output produced by ERA-Interim unmatched the recorded pattern (Fig. 4(e)).Although FNL (64.0°) and GFS (65.0°) were better in terms of GE, none of the datasets reached values for wind direction into the benchmark range (< 30°).NCEPR2 (-0.8°) and ERA-Interim (-4.4°) got MB values into the benchmark ranges.However, the ERA-Interim MB value resulted from significant negative and positive differences between wind speed records and computed values, as shown by GE (90.1°).
The metrics indicated that GFS, FNL, and NCEPR2 provided the best fit to model temperature, wind speed, and wind direction.

Air Quality
The GFS, FNL and ERA-Interim produced the best fit to model the maximum CO daily 1-h mean concentration.It is important to note that GFS and ERA-Interim did not overestimate or underestimate (MB = 0.0 mg m -3 , Table 2) whereas NCEPR2 (MB = 1.0 mg m -3 ) and ERA5 (MB = 1.1 mg m -3 ) overestimated these concentrations.The magnitude of difference between the modeled and observed values were smaller for GFS, FNL, and ERA-Interim (RMSE = 0.6 mg m -3 ) compared to NCEPR2 (1.2 mg m -3 ) and ERA5 (1.5 = mg m -3 ).However, there was no strong linear relationship between the recorded and modeled values (r = 0.4 for GFS, FNL, and ERA5).The FNL (91.7%) and GFS (83.3%) presented the highest percentages of captured records (Table 3).
Similarly, GFS, FNL, and ERA-Interim produced the best fit for modeling the maximum CO 8-h mean concentration (Fig. 5).On the other hand, as seen from Table 2, the datasets showed a slight underestimation (MB = -0.1 mg m -3 ) in the concentrations.The magnitude of difference between the modeled and observed values were lower for GFS and FNL (RMSE = 0.2 mg m -3 ).It is also worthy of note that there was no strong linear relationship between the recorded and modeled values (r = 0.4) for FNL.For these concentrations, FNL (100.0%) and GFS (95.8%) had the highest capture percentages (Table 3).
The best options for modeling daily mean PM2.5 concentrations were GFS and ERA-Interim (Table 3, Fig. 6).On the other hand, ERA5 had the lowest performance.However, there was a weak linear correlation between the recorded and modeled values for the GFS (r = 0.2) and ERA-Interim (r = 0.3).The GFS and ERA-Interim produced the most captured records (75.0%) to model PM2.5 concentrations.ERA5 produced the least capture percentage (45.8%).
The majority of the databases had a capture percentage above 80.0% of the maximum O3 daily 8-h mean concentrations (Table 3, Fig. 7).ERA5 had the best metrics to model these concentrations.All of the datasets provided mean O3 profiles that matched the behavior of the records during the early hours of the morning, though the modeled values were higher after noon and in the afternoon (Fig. 7(f)).datasets.They were consistent, despite minor differences.The area adjacent to the NW of the urban area presented concentrations between 75.0 µg m -3 and 80.0 µg m -3 (purple color) for GFS and FNL.However, ERA-Interim and ERA5 generated values between 80.0 µg m -3 and 85.0 µg m -3 (red color).Despite the high percentages of records captured, the linear relationship between recorded and modeled values was low (r between -0.1 and 0.2, Table 3, Fig. 7).
In terms of NO2 monthly mean concentrations, GFS, FNL, and ERA-Interim captured 100.0%records from all passive stations that measure the pollutant while 73.3% was captured by ERA5 (Table 3, Fig. S2(a)).However, ERA5 captured the O3 monthly mean concentrations in 87.5% of the stations (Table 3, Fig. S2(b)).
FNL and GFS were the datasets best fit for modeling meteorology and air quality, thus suggesting their ability to generate IC and BC for atmospheric modeling purposes in the study region.

DISCUSSION AND CONCLUSIONS
We used a state-of-the-art 3-D numerical model to simulate meteorology and air quality, using atmospheric analysis and reanalysis datasets to define the IC and BC of a nested domain with a high spatial resolution (1 km) that covers the territory of an Andean equatorial city.The modeling approach examined the direct effects of aerosols on meteorology, based on the benefits of using an online approach reported in the literature.In this context, we begin with the premise that the best database should accurately model both meteorology and air quality.
The metrics indicated that the model's performance was better when the analysis datasets (FNL and GFS) were used.The lack of improvement, and more specifically, the decrease in the model's performance when using reanalysis products (NCEPR2, ERA-Interim, and ERA5), was unexpected due to the high quality of these datasets since they were made with more observations.FNL and GFS performed best in terms of modeling surface temperature, which is an important meteorological parameter that is strongly related to the height of the planetary boundary layer (PBL) and atmospheric stability.These databases provided minimum mean temperatures during the early morning hours close to the observations (Fig. 2(f)), thus suggesting that the corresponding heights of the PBL and atmospheric stability were also adequately modeled.The peak of vehicular traffic emissions in Cuenca typically occurs between 07:00 and 08:00 LT, due to stable atmospheric conditions that prevent dispersion and produce higher CO concentrations around these hours.This relationship is consistent with the best results and metric values for modeling CO concentrations in the analysis datasets.With lower PBL modeled heights (not indicated) and higher CO concentrations, NCEPR2 and ERA5 produced lower minimum temperature values than the records.This interaction explains why these datasets overestimated CO concentrations.
In the case of wind speed, despite that all datasets obtained metrics within the reference ranges, FNL performed the best and, on average, its results were more in line with the daily behavior of this parameter (not indicated).
All datasets overestimated PM2.5 concentrations, suggesting that other components should be evaluated, such as the quality of the emissions inventory and the influence on the modeling of direct and indirect effects between aerosols and meteorology.
Although ERA5 captured the majority of the short-term and long-term O3 records, the model performance with respect to the meteorological variables and other air quality variables was the lowest.On average, all datasets modeled higher hourly solar radiation levels than the records (not indicated), with ERA5 producing the most accurate results.Because solar radiation influences the photochemical production of O3, the more approximate levels of solar radiation of ERA5 explains, at least in part, why this database presented a higher performance in O3 modeling.The overestimation of solar radiation highlights the importance of evaluating other schemes, such as cumulus parameterization, O3 concentrations in IC and BC, and the emission inventory of O3 precursors.
The results suggest that instead of reanalysis datasets, FNL or GFS should be used to generate the IC and BC for modeling meteorology and air quality in the equatorial Andean region.This finding is consistent with the preliminary decrease reported when using NCEPR2 instead of GFS or FNL to model four volcanic ash emission and dispersion events in Ecuador, which occurred in December 2012, July 2013, February 2014, and August 2015 in the Tungurahua volcano and the Cotopaxi volcano, respectively (Parra, 2018c).We emphasize that these events occurred over a period of months and years, and that the modeling done on a national scope.The current contribution and the volcanic ash dispersion study results were consistent, and it was suggested that the IC and BC be generated using GFS or FNL for other regions of Ecuador and at the national level.
Some studies found similar results in other parts of the world.Table S4 summarizes the main conclusions of selected studies on the impact of global atmospheric datasets.These conclusions, as well as the results of this contribution, suggest that reanalysis products are not always the best default option.Global atmospheric datasets require dedicated evaluations to determine their use.
Although NCEPR2 includes a higher number of observations compared to GFS and FNL, this reanalysis product is made with lower horizontal (2.5°) and vertical (18 vertical layers) resolution (Table S1).The lower spatial resolution of NCEPR2 could partly explain why this reanalysis obtained a lower performance than GFS (0.5°, 27 vertical layers) or FNL (1.0°, 27 vertical layers (Wagner et al., 2018).However, when using ERA5, a fifth-generation atmospheric dataset and one of the most advanced reanalysis products with high horizontal (0.25°) and vertical (38 vertical layers) resolution, overall modeling performance decreased.
Although the study of the causes of the decline in the performance of the reanalysis datasets is beyond the scope of this manuscript, we argue that one reason corresponds to the complex conditions imposed by the Andes Mountains.The limitations of the reanalysis datasets are evident in areas with complex orography, particularly in mountainous regions and in areas where assimilation and processing schemes cannot reproduce actual atmospheric processes (WMO, 2018).
Although the national meteorological and hydrological services of the Andean countries provide a large amount of surface data, some of these observations do not meet the World Meteorological Organization's requirements (WMO, 2020).In addition, national databases are partially shared with international entities in charge of generating global atmospheric databases.On the other hand, high altitude areas are poorly instrumented (Condom et al., 2020;Kull et al., 2021), and the equatorial Andean region has poor coverage of radiosonde observations (Durre et al., 2018;Cazorla et al., 2021).These limitations imply an incomplete description of the atmospheric conditions in this region, which has an impact on the assimilation processes used to prepare global atmospheric datasets.
Aside from the few observations from the equatorial Andean region, the reduction in modeling quality when using reanalysis datasets, particularly ERA5, was unexpected.Increased records from this region are required to prepare global atmospheric datasets and improve understanding of atmospheric processes in the equatorial Andean region.The atmosphere in this region exhibits complex behavior, with unique characteristics that combine a complex topography, the influence of the Intertropical Convergence Zone, the presence of breezes from the Pacific Ocean and the Amazon region, which promote vertical convective movements, and consequently, the persistence of clouds.The atmospheric models include parameterizations and schemes that were developed and tested primarily for the northern hemisphere.These constraints are currently affecting the development of atmospheric databases.These limitations, which appear to be more prevalent for reanalysis products, have an impact on the generation, and therefore the quality of the information stored by the atmospheric datasets for the equatorial Andean zone.Reanalysis products are not always the best choice for modeling in this area.
Initial and boundary conditions are essential components for mesoscale and regional atmospheric modeling.This contribution provided information for selecting the best atmospheric datasets to model meteorology and air quality in the Andean region of Ecuador.Other parameters and components, such as the parameterization of the earth's surface and cumulus clouds, the effect of assimilation of meteorological and air quality data, the influence of spin-up time, and the indirect effects between meteorology and aerosols, must be evaluated.Complementary studies based on updated and improved emission inventories, with coverage of other months and the use of other chemical mechanisms and aerosol models, will be required in the future.
Within the surface layer, the bottom 20 to 200 m of the planetary boundary layer, the frictional drag, heat conduction, and evaporation from the surface produce significant changes in wind speed and temperature with height, defining the land-atmosphere interaction.The surface layer and components as the planetary boundary layer and land surface are typically parameterized in atmospheric modeling (Table S2).Although this study focused on the influence of IC and BC, we highlight the need to assess the effect of parameterized variables in the future.

(
WHO, 2006).Based on the recommendations ofSimon et al. (2012), we present as metrics the MB, RMSE, the fractional bias (FB), the mean normalized bias (MNB), and the correlation coefficient (r).The following equations define the last three indicators:

Table 1 .
Metrics for meteorological modeling.Green background color indicates metrics with values into the benchmark ranges.

Table 2 .
Metrics for air quality modeling.Green background color indicates metrics with the best performance values.

Table 3 .
Percentages of captured records by modeling.