Chemical Analysis of Surface ‐ Level Ozone Exceedances during the 2015 Pan American Games

: Surface ‐ level ozone (O 3 ) continues to be a significant health risk in the Greater Toronto Hamilton Area (GTHA) of Canada even though precursor emissions in the area have decreased significantly over the past two decades. In July 2015, Environment and Climate Change Canada (ECCC) led an intensive field study coincident with Toronto hosting the 2015 Pan American Games. During the field study, the daily 1 ‐ h maximum O 3 standard (80 ppbv) was exceeded twice at a measurement site in North Toronto, once on July 12 and again on July 28. In this study, ECCC’s 2.5 ‐ km configuration of the Global Environmental Multi ‐ scale (GEM) meteorological model was combined with the Modelling Air ‐ quality and CHemistry (MACH) on ‐ line atmospheric chemistry model and the Town Energy Balance (TEB) urban surface parameterization to create a new urban air quality modelling system. In general, the model results showed that the nested 2.5 ‐ km grid ‐ spaced urban air quality model performed better in statistical scores compared to the piloting 10 ‐ km grid ‐ spaced GEM ‐ MACH model without TEB. Model analyses were performed with GEM ‐ MACH ‐ TEB for the two exceedance periods. The local meteorology for both cases consisted of light winds with the highest O 3 predictions situated along lake ‐ breeze fronts. For the July 28 case, O 3 production sensitivity analysis along the trajectory of the lake ‐ breeze circulation showed that the region of most efficient O 3 production occurred in the updraft region of the lake ‐ breeze front, as the precursors to O 3 formation underwent vertical mixing. In this updraft region, the ozone production switches from volatile organic compound (VOC) ‐ sensitive to NOx ‐ sensitive, and the local net O 3 production rate reaches a maximum. This transition in the chemical regime is a previously unidentified factor for why O 3 surface ‐ level mixing ratios maximize along the lake ‐ breeze front. For the July 12 case, differences between the model and observed Lake Ontario water temperature and the strength of lake ‐ breeze opposing wind flow play a role in differences in the timing of the lake ‐ breeze, which impacts the predicted location of the O 3 maximum north of Toronto. and M.M.; software, S.R., A.A., R.M. ‐ A. and S.L.; supervision, C.S., P.M., S.B. and J.R. B.; validation, C.S. and S.R.; visualization, C.S.; writing—original draft, C.S.; writing—review and editing, C.S., S.R., J.Z., M.M., P.M., S.L., S.B., D.S. and J.R.B.


Introduction
A large fraction of the world population now lives in urban areas where they are exposed to harmful air pollutants. A recent study in the U.S. showed that urban areas have up to twice the burden of asthma as compared to rural areas [1]. In addition to respiratory effects, there are cardiovascular, metabolic function, pregnancy, neurologic and psychiatric outcomes associated with human exposure to air pollution [2][3][4]. All of these effects lead to more emergency room visits, hospitalizations, absenteeism and higher health care costs. Another recent study found that paediatric asthma incidence is associated with exposure to traffic-related air pollution [5]. They estimate that 1-in-5 new cases of childhood asthma are linked to on-road mobile emissions in a large North American city. All of these recent studies highlight how urban air pollution continues to negatively impact human health.
While the mix of pollutants most responsible for health effects is still unknown, nitrogen dioxide (NO2) has been used as a proxy for traffic-related air pollution in epidemiological studies because of the large NOx emission (sum of oxides of nitrogen, NO+NO2) from transportation sources and its short atmospheric lifetime [6]. The short NOx lifetime results in large concentration gradients between urban areas and regional background conditions. Ozone (O3) is another important component of surface-level air pollution. The Canadian Air Quality Health Index (AQHI) has two terms related to the acute human health impacts from the exposure to NO2 and O3 gases, as well as a term for exposure to particulate matter [7]. Odd oxygen (Ox) is defined as the sum of O3 and NO2. Ox has a longer lifetime in the atmosphere and can impact ecosystems on regional scales [8,9]. The spatial gradient in Ox between urban and rural areas is much weaker than NOx because O3 is not emitted from localized urban sources. Rather, O3 is a secondary pollutant, which is largely dependent on chemical production and depositional loss processes. Exposure of vegetation to ground-level O3 also leads to billions of lost revenue from lower crop yields worldwide [10].
Surface-level O3 mixing ratios have been a challenge to control because Ox production in the lower atmosphere is a non-linear function of its precursor gases. The two classes of important precursor gases to O3 formation are volatile organic compounds (VOCs) and oxides of nitrogen (NOx). The precursor gases also have a wide range of emission sources, both anthropogenic and biogenic. The mixing and dilution of these precursor gases, as they transport away from their emission sources, is highly dependent on wind speed and atmospheric vertical stability. Net production of O3 occurs in the troposphere when VOCs react with oxidants to produce peroxy radicals (RO2), which, in turn, react with NO to form NO2 [11]. In the presence of sunlight, the NO2 will photolyze to form an oxygen atom that can quickly react with the abundant molecular oxygen in the troposphere to form O3. Collectively, the availability of precursor emissions, local chemical production and dispersion (i.e., low wind speeds and vertical stability) are the key processes that impact whether O3 mixing ratios reach unhealthy levels within and downwind of urban areas. High concentrations of ground-level O3 continue to be observed in the Great Lakes Region despite significant reductions in the precursor emission [12,13]. The exceedances of the national standards for O3 in both the U.S. and Canada occur most often on hot, sunny days, downwind of the major urban areas and frequently near shorelines with well-developed lake-breeze circulations [14][15][16]. Foley et al. was one of the first comprehensive studies performed between 1994 and 2003 to measure O3 over Lake Michigan [17]. They found that, under southwest winds, precursor emissions could be transported out over the lake where they can photochemically react to form O3 and then can be swept on-shore in an afternoon lake-breeze. Cleary et al. measured O3 over Southern Lake Michigan using ferry transects across and also with a shoreline-based differential optical absorption spectrometer (DOAS) instrument [18]. Ferry O3 observations across the lake were on average 3.8 ppbv higher in O3 than the long-path DOAS. Comparisons of the ferry measurements to community multiscale air quality (CMAQ) model predictions showed significant model over-predictions over the lake with the bias generally larger to the east side of the lake (i.e., more chemically aged air masses). More recently, McNider et al. used the WRF (weather research forecast) and the off-line CMAQ model to study the effect of mismatches between meteorological models and air-quality models using the off-line paradigm, on O3 over-prediction over Lake Michigan. They also showed that the use of satellite-derived skin temperatures improves model predictions for temperature. Both of these factors were shown to improve their simulations of the lake-breeze circulation and atmospheric chemistry. They concluded that having a land surface scheme that can adjust moisture based on the surface energy budget is important in reducing model temperature errors [19].
Loughner et al. studied the impact of the Chesapeake Bay breeze on surface air quality in the Baltimore area [20]. On one day, they found that the convergence zone of the bay-breeze resulted in the pollution from Baltimore being transported to the upper boundary layer. Goldberg et al. found that O3 was measured at higher concentrations over Chesapeake Bay compared to a nearby upwind site [21]. During the 10-day field study, the measured O3 exceeded the U.S. 8-h standard on 4 different days. It was hypothesized that lower dry deposition rates, lower boundary layer heights and high photolysis rates due to clear skies over the water all contributed to the higher O3 levels.
Southern Ontario is bounded by three of the five North American Great Lakes: Lake Huron, Lake Erie, and Lake Ontario. Sills et al. characterized lake-breeze frequency and impact on air pollution during the Border Air Quality and Meteorology (BAQS-Met) study in the Detroit/Windsor area [22]. They noted that the cool, shallow inflow to the urban boundary layer and the clear skies associated with lake-breeze days were two important factors in the concentration of precursor emissions and photochemical O3 formation ( Figure 1). Sills et al. also compared an observation-based mesoscale analyses with a high-resolution weather forecast model. The weather forecast model accurately predicted the occurrence and type of lake-breezes, but selected cases showed differences in the timing and in-land penetration of the fronts. Brook et al. summarized air quality observations from the BAQS-Met study [23]. They concluded that the high-pressure systems associated with well-developed lake-breeze circulations could last several days, resulting in urban emissions being concentrated in the same recirculating air mass. They found that lake-driven meteorology enhances secondary pollution production and that local anthropogenic emissions have their impact closer to urban areas than would occur in the absence of the lake effects. They concluded that the presence of the large lakes causes an enhancement in the impact of local pollutant emissions on O3.
Makar et al. showed with a chemical transport model that lake-breeze convergence zones that formed in the afternoon downwind of Detroit/Windsor confined pollutants in narrow bands horizontally (few tens of kilometers) and were locations of rapid ozone formation, though the chemistry explanation for ozone enhancement along the front was not developed in that paper [15].
Hastie et al. was one of the first studies to report the impact of a lake-breeze on the Toronto air pollution plume [24]. They observed a 30 ppbv O3 increase in a 1-h increment at a downwind rural town associated with the passage of a lake-breeze front from Lake Ontario. Geddes et al. reported the observed trend in the daily maximum 8-h average O3 for eight Toronto sites and showed a near zero change in this average mixing ratio between 2000 and 2007 [25]. They argued that competing effects of: (1) a small positive change in O3 (due to NO emission deceases over the time period and a shift in O3 photo-stationary state from NO2 to O3) coupled with (2) small decreases in VOC emissions (in turn decreasing the magnitude of O3 production) resulted in near zero change in O3 mixing ratio over the 8-year period. More recently, Pugliese et al. extended the O3 trend for Toronto over a longer time period (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012) and found a small O3 decrease at urban and suburban sites, although decreases were not statistically significant (except for the northern suburban town of Newmarket) [26]. They also examined O3 observations for Toronto in 2012, which was a year with a large number of O3 exceedances. They found that year 2012 was an anomaly due to a higher-than-average number of days with light winds and winds from the west to the southeast, which cause air to be transported from other southern urban locations. There was also a higher frequency of days with stagnant wind conditions and accompanying higher-than-average frequency of clear sky days. The role of oxygenated volatile organic compounds (OVOCs) was identified as a knowledge gap due to the lack of measurements, as was the lack of data in the vertical for periods of O3 exceedances.
Wentworth et al. assessed the impact of lake-breeze circulations on pollutants in the Greater Toronto Area using a multi-year observation data set [16]. They found that the average daily 1-h O3 maximum was 42-49% higher on lake-breeze days. Case study analysis found that sites frequently located within the lake-breeze circulation experienced much higher O3 than sites just outside the circulation, even though all of the sites remained within the same synoptic-scale environment that dictates long-range transport. Morning NOx levels were greater on lake-breeze days, probably due to the stagnant, vertically stable conditions associated with lake-breeze days. They concluded that there is likely O3-enriched air over the lake, which is transported onshore in the late morning, resulting in the high O3 events for cities along the shoreline. Possible mechanisms were identified for high O3 mixing ratios over the lake including: (1) reduced mixing height over the lake and near the shoreline, (2) reduced dry deposition to the lake surface, (3) a night-time land breeze that fumigates the lake surface with pollutants from the city followed by rapid O3 production in the morning, (4) subsidence of O3-rich air from the upper boundary layer over the lake and (5) reduced cloud cover and increased solar radiation reaching the surface over the lake. Wentworth et al. (2015) concluded that models need to represent lake-breeze circulations accurately to improve air quality predictions in the region [16].
In preparation for the July 2015 Pan American Games, held in the Greater Toronto Hamilton Area (GTHA), Environment and Climate Change Canada (ECCC) developed a suite of forecast demonstration products using the GEM (global environmental multi-scale) numerical weather prediction model as the core to the forecast systems. The demonstration forecasts included high spatial resolution forecasts for heat indices, severe weather and Lake Ontario wave heights and air quality forecasts made with the GEM-modelling air quality and chemistry (MACH) chemical transport model. The real-time, high-resolution versions of GEM showed improved forecast skill compared to the ECCC operational forecast products at the time [27].
The 2015 Pan American Games science demonstration also included an intensive meteorology and air quality field study that took place in July 2015. The summer of 2015 had two cases in the GTHA when O3 exceeded the Ontario daily 1-h standard of 80 ppbv: 12 July and 28 July. Hence, these periods provide a valuable data set for testing meteorology and air quality model developments. The first goal of this paper is to assess the performance of the high-resolution 2.5-km GEM-MACH-town energy balance (TEB) model for predicting O3 and NO2 within the GTHA during the Pan American Games study. The special version of the GEM-MACH air quality model used for this study includes the town energy balance module (TEB) for urban surfaces and a high-resolution surface data assimilation scheme for updating both land and lake temperatures, as well as water vapour and surface energy fluxes at 1-h intervals. The second goal of this paper is to examine, in detail, model performance for the two O3 exceedance periods. For the 28 July case period, a chemical coordinate analysis was performed to diagnose the chemical regime for O3 production along the trajectory of the lake-breeze circulation. An important factor was uncovered, which contributes to the local enhancement in O3 mixing ratio along the lake-breeze front. It was found that the vertical mixing associated with the lake-breeze front resulted in a shift in O3 production from VOC-to-NOx sensitive conditions. O3 production maximizes during this chemical transition and, thus, contributes to the observed O3 surface maximum along the lake-breeze front.

Evaluation Datasets for the 2015 Pan American Games Period
The Greater Toronto Hamilton Area (GTHA) is a region that includes both the Greater Toronto Area, located on the northwestern shore, and the city of Hamilton, located on the western shore of Lake Ontario ( Figure 2). As of 2018, the population of the GTHA approached 6.8 million people and is projected to increase to 10.2 million by 2046 [28]. NOx and VOC emissions in the GTHA are dominated by transportation, residential and commercial sources, while Hamilton also has steel production industries producing point source emissions of SO2 and NOx. Under prevailing summertime southwesterly winds, the GHTA is also downwind of several major U.S. cities (Chicago, Detroit and Cleveland), as well as pollution sources from the industrial Ohio Valley. Figure 2 illustrates the 11 air quality (AQ) measurement stations active within the GTHA in 2015. These stations are part of ECCC's National Air Pollution Surveillance (NAPS) monitoring network. One additional AQ measurement site was deployed on Toronto Island (shown on Figure 2) to support forecasting for the Pan American Games water sports. In total, there are 37 Canadian AQ stations within the 2.5-km model domain ( Figure 3). These NAPS stations continuously measure the mixing ratio of the gases O3, NO, NO2 and SO2. In this study, the focus of the evaluation is on O3, NO2, NOx (sum of NO + NO2) and Ox (sum of O3 + NO2). In the NAPS network, O3 is measured by UV absorption detection [29]. NO is measured by chemiluminescence detection [30]. NO2 is measured by heated molybdenum conversion of NO2 to NO [30,31]. An alternating inlet measures NO and then NOx. The NO2 is measured from the difference between NOx and NO.
The ECCC mobile AQ laboratory (CRUISER) was deployed and made air quality measurements for several days of the Pan American Games period [27]. The CRUISER measured carbon monoxide (CO) with a UV absorption technique [32]. NO2 and NO were measured with a selective photolytic NO2-to-NO converter, followed by NO detection by chemiluminescence [30].
Routine surface meteorological measurements were available for sites operated permanently by the Meteorological Service of Canada (MSC). Over 40 additional compact weather stations were added in transects running north-south in Toronto [27]. The weather station data, low-level radar and GOES-13 visible satellite cloud observations were combined to form a mesoscale meteorological analysis [33]. The mesoscale analysis was used in this study to diagnose the location of lake-breeze fronts for the air quality case studies.

GEM-MACH-TEB Model Description and Setup
The GEM-MACH (global environmental multiscale-modelling air quality and chemistry) numerical model is a chemical transport model composed of dynamics, physics and atmospheric chemistry modules [34][35][36][37][38], where the chemistry modules run on-line within the GEM meteorological model [39,40]. The on-line paradigm has significant advantages over off-line models, including removal of potential interpolation errors between weather forecast and air-quality model grids, avoiding the use of conflicting boundary layer parameterizations between weather and air-quality models [19] reducing the large amount of memory space required for the storage of high resolution meteorological model output files, and reducing the overall processing time required to carry out a large domain, high resolution simulation. The on-line paradigm also enables chemistry feedbacks on weather. However, in this work, the feedback of MACH predicted aerosol on the meteorological forecast was not included. Here, we wanted to keep the MACH version similar to the ECCC GEM-MACH operational version to see the differences related to the urban canopy model and the higher spatial resolution.
Numerical modelling for the Pan American Games study was undertaken using a nested GEM-MACH domain, comprised of an outer model domain at 10-km grid spacing (750 × 620 grid cells) and an inner nested model domain at 2.5-km grid spacing (502 × 402 grid cells) centered on the GTHA. Figure 3 shows the outer GEM-MACH model domain in green and the inner 2.5-km GEM-MACH-TEB domain contoured by predicted surface O3 mixing ratio. Both domains used the same GEM-MACH version (v2.3.1). Both GEM-MACH domains also had 80 vertical levels, up to a pressure of 0.1 hPa, and used a staggered vertical grid (based on Charney-Phillips discretization) with thermodynamic and chemical tracers on alternating levels from the dynamical tracers. In Figure 3, the 10-km GEM outer domain in blue provides the meteorological lateral boundary conditions for the 10-km GEM-MACH model (in green). The 10-km GEM model is used to create a meteorological analysis to initialize the 10-km GEM-MACH model and 2.5-km GEM-MACH-TEB model. The surface scheme in GEM-MACH-TEB over land is based on an advanced soil moisture and land surface temperature assimilation system (termed CALDAS, Canadian Land Data Assimilation System) [41]. Hourly lake-water temperature output, for driving the GEM-MACH-TEB lake surface, comes from a 2-km coupled lake-atmosphere model system (NEMO, Nucleus for European Modelling of the Ocean). NEMO air temperature is initialized from the 10-km GEM analysis. The 10-km GEM analysis includes a surface analysis for lake temperature [42,43].
As noted above, another recent improvement to the 2.5-km GEM-MACH model was the inclusion of the town energy balance (TEB) parameterization [44][45][46]. TEB simulates the urban heat island effect more accurately by including terms in the surface heat-flux equation for building roofs, walls and road asphalt. TEB also includes the change in optical and thermal properties for urban surfaces and the change in aerodynamic roughness associated with urban buildings. Ren et al., (2020) describes the implementation of TEB into GEM-MACH [47]. A parameterization of traffic heat flux scaled from mobile pollutant emissions (as published in [48]) was used to calculate a traffic heat input into GEM-MACH using the most recent CO and NOx mobile emission maps. For the remainder of this paper, the urbanized GEM-MACH model is referred to as GEM-MACH-TEB. Figure 4 shows a schematic of the dynamic, physics and chemistry processes represented in GEM-MACH-TEB. GEM-MACH-TEB includes a comprehensive chemistry process package that represents gas-phase chemistry, aqueous-phase chemistry and particle microphysics (nucleation, condensation, coagulation, settling and deposition). Table 1 lists the key model settings for chemistry and physics. The 2.5-km inner domain was initialized, based on an interpolation from the 10-km outer domain. The inner model was run in 24-h cycles with a 24-h pollutant output used to initialize the next model cycle. Biogenic emissions are calculated off-line with the BEISv3.09C (Biogenic Emission Inventory System) module [38] with updates for the biogenic emission rate for the boreal forest tree species based on a consideration of 3-year SPOT satellite-derived leaf area index (LAI) maps (see Figure A1 and Table A1). GEM-MACH v2.3.1 includes an improvement to the gas-phase dry deposition routine by assimilating LAI based on the 3-yr SPOT satellite-derived climatology. The lateral boundary conditions for the outer 10-km GEM-MACH domain were also updated by averaging the global MOZART4 model output for a summer season [49]. GEM-MACH-TEB was run on a Cray XC40 supercomputer with parallelization using a message passing interface (MPI) and OpenMP.

Pollutant Emission Inventories and Emission Processing for GEM-MACH-TEB
For the model simulations for July 2015, the 2015 Canadian National Pollutant Release Inventory (NPRI) was considered for point-source emissions, the 2015 Canadian Air Pollutant Emission Inventory (APEI) was used for area-source emissions and the 2013 APEI was considered for Canadian transportation emissions, both on-road and off-road [52-54]. On-road NOx emission for Ontario's Drive Clean Program Phase 1 area, which is approximately the area of the GTHA, from the 2013 APEI were 3869 tonnes and were created using the MOVES (mobile vehicle emission simulator) traffic emission model [55]. For the U.S., 2017 projected U.S. emissions were available, with mobile emissions also based on the MOVES traffic model. Note that total NOx on-road emissions reported for Ohio for 2015 were 124 kt [56,57], which is comparable to our projected 2017 on-road NOx emissions for Ohio (118 kt).
To generate 2.5-km gridded emissions for GEM-MACH-TEB, some improvements were made to the spatial surrogate fields and temporal profiles used to transform the emissions inventories into model-ready emission maps, particularly for on-road emissions. To process Canadian on-road emissions, a new set of spatial surrogates was generated using the Canadian National Road Network and population shape-files. The Canadian on-road inventory was also modified to split emissions between different road types to make better use of the new spatial surrogate fields. New diurnal profiles were used for light-duty and heavy-duty vehicles for major Canadian highways based on diesel traffic count data from a special field study [32], as described in [58]. The impact of using the new Canadian spatial surrogate fields was to redistribute the on-road mobile emissions more to the major highways and less to city streets compared to previous emissions.
In order to switch to MOVES to generate the on-road mobile emission inventories in Canada and the U.S., a new set of on-road source classification codes was introduced (2 new rural types and 2 new urban types), which required the creation of 4 new spatial surrogate fields. In addition, new off-network vehicle emission sources (from cold-starts, idling and fuel venting) were introduced with MOVES and required the creation of 11 new spatial surrogate fields (e.g., gas stations, truck stops, landfills, bus terminals and various building types). The Canadian population spatial surrogate field was also capped at 2000 people/km 2 , based on results from [59], which show no further increase in traffic emissions with population density above this threshold. The Canadian population spatial surrogate field is used to allocate local road emissions, and having this cap reduced the on-road emissions allocated to downtown Toronto.
The new source classification codes in MOVES also enabled very detailed weekly and hourly temporal emission profiles to be created and used for the U.S. emissions (based on the U.S. Vehicle Travel Information System). A recent version of the U.S. EPA SPECIATE database (v4.5) was used for VOC and NOx speciation profiles and applied to both the Canadian and U.S. mobile emissions [60]. The new off-network vehicle sources added with MOVES increase urban background emissions relative to previous inventories. The changes to the U.S. temporal profiles tend to increase the onroad emissions in the morning hours but reduce emissions in the afternoon rush hours and overnight [58].
The emissions of O3 precursor species (NOx, VOCs and CO) in Ontario and neighbouring states have decreased considerably over the past 15 years, largely due to the use of tailpipe catalytic converters on vehicles and the closure of several fossil-fuel power plants. The Supplementary includes a model emission map ( Figure S1) and summaries of the emission trends over the land decade. It is interesting to look at how O3 mixing ratios in Ontario have responded to these emission reductions [61]. Even with large emission reductions, summertime O3 continues to exceed the 1-h 80 ppbv ambient air quality criterion, (AAQC) for O3 (e.g., 6 times at North Toronto site in 2015). The O3 trend in the mean concentration for the summertime months for 39

Evaluating Pollutant Emissions for a Downtown Toronto Model Grid Cell
The CRUISER mobile laboratory was deployed on 27 July 2015 with a driving route created to spatially sample throughout one 2.5-km × 2.5-km GEM-MACH-TEB grid cell (centered at 43.66 N, −79.39 W) in downtown Toronto under light southerly winds and sunny skies. The model output was extracted along the path of the CRUISER and both model and measurement time series were averaged for the time the CRUISER was within the grid cell (15:45-17:30 UTC). The median value of a running 120 s period was used to filter the CRUISER data, thus, removing points with very high concentrations that represent fresh emissions from nearby vehicles and that would not be representative of the grid-volume concentrations [63]. The time-average of the measured values for NOx and CO was 17.6 ± 3.7 ppbv and 421 ± 12 ppbv, respectively. The modeled mixing ratios of the two species for the same period (15:45-17:30 UTC) were 19.3 ppbv and 519 ppbv, respectively, corresponding to model over-predictions of +10% and +22%. These numbers may be compared to regional AQ model evaluations in the literature [36], where North American NO2 normalized mean biases from the GEM-MACH, WRF-CHEM and WRF-CMAQ modelling platforms ranged from 8.09% to 30.20%, NO normalized mean biases ranged from +42.15 to −61.96%, and CO normalized mean biases ranged from −7.02 to −20.09%. Urban-station-only evaluations for a European domain in the same study showed larger normalized mean biases for three different WRF-CHEM configurations (NO2: −60.58% to −61.11%, NO: −81.46% to −86.67% and CO −34.57%). A recent global model CO intercomparison was published using satellite data, but this is not at the same scale as our study [64]. The models in this latter case had an average mean bias of −20%. The CO evaluation results shown here are considerably better than prior results with ECCC's predecessor off-line chemical transport model (AURAMS). For a simulation of the 2007 BAQS-Met study, AURAMS CO predictions had a normalized mean bias of +134% for a measurement comparison at an urban site (Windsor) in Southern Ontario [65]. Figure 5 shows the time series for modeled and measured hourly O3 and NO2 mixing ratio at the North Toronto measurement site for July 2015. The characteristic diurnal cycling for O3 and NO2 was evident. There were two periods in Figure 5a, indicated with arrows, when the observed hourly O3 surpassed the 80 ppbv Ontario AAQC for ozone. On 28 July at 05:00 p.m. EDT, the observed O3 reached 90 ppbv. On 12 July at 04:00 p.m. EDT (local time), O3 reached 84 ppbv. The model underpredicted the 12 July exceedance by about 20 ppbv. It did capture the observed daytime maxima (70 ppbv) on the days before and after, and also the night-time O3 minimum and NO2 maximum that occurred before the exceedance. Figure 5b also includes the time series for modelled surface-level wind speed at this station. Wind speed plays a critical role in the dilution of emissions and mixing of pollutants vertically. For the two case study periods, the modelled daytime wind speeds were quite low (<4 m/s). Conversely, afternoons with modelled wind speeds >8 m/s had the lowest ozone mixing ratios. Both case study time periods will be discussed more in Section 3.4.

Model Evaluation
The evaluation of the 2.5-km GEM-MACH-TEB model predictions is described in detail in the Supplementary (Figures S2-S5, Tables S1-S3). Here, we briefly summarized the most important results. Table 2 compares the statistics for the 10-km GEM-MACH model with 2.5-km GEM-MACH-TEB for all the Canadian sites in the high-resolution 2.5-km domain (including hourly O3 data from the U.S. Air Quality System). The 10-km model is analogous to the operational model used by the Meteorological Service of Canada for AQHI forecasting. The metrics are typical of those found in the literature to assess AQ model performance [66][67][68]; namely, normalized mean bias (NMB), root mean square error (RMSE) and the correlation coefficient (R). In looking at the NOx statistics, the model NMB was +5.4% and correlation coefficient, R = 0.64. The Ox NMB was also very good with almost no bias, −0.3%, and a strong correlation coefficient, R = 0.80. In looking at the scores for NO2 and O3, there was a small O3 over-prediction (NMB=5.4%) and a larger NO2 over-prediction (NMB =1 9%). The RMSE for NO2, O3, NOx and Ox were all lower for the high-resolution 2.5-km nested model compared to the 10-km resolution GEM-MACH piloting model. Correlation coefficient scores were slightly worse with the higher spatial resolution model for NO2, NOx and Ox and slightly better for O3. Simon et al. (2012) reported RMSE values in the range, 15-20 ppbv, for hourly O3 comparisons in other model validation studies [68], which are higher than the RMSE values here. The SI section also includes a section evaluating the 2.5-km GEM-MACH-TEB model against the 10-km GEM-MACH model for all the hourly data at the 12 GTHA sites. The same conclusion was found that the higher resolution model, with the urban canopy model, outperformed the 10-km version for most statistics.  Figure 6 shows a scatterplot of the correlation between modelled and measured daily 1-h maximum O3 at the North Toronto site for July 2015. The correlation coefficient, slope and model intercept (y-axis) were R = 0.85, 0.97 and −2.2 ppbv, respectively. The SI section also included a model evaluation using the Canadian 8-h O3 metric. Overall, results here and in the SI section show that the model performance was comparable to other peer AQ models [63,[66][67][68][69][70] and met the performance criteria defined in the U.S. EPA guidance documents [71,72]. The SI section also included an evaluation of Ox predictions grouped according to three different wind conditions (southwesterly, northerly and light wind speeds from any direction) and two temporal periods (mid-morning and mid-afternoon). The selection of data for different wind conditions helps to separate the Ox data by source region and assess how the model performs for these different wind conditions. Figure S5 in the SI section summarizes these results.

Case Study Analysis for Periods of O3 Exceedance
In this section, the two periods with 1-h O3 exceedances are analyzed to understand the meteorological and chemical factors that contribute to the build-up of poor air quality in Toronto. We begin in Section 3.4.1 by characterizing the synoptic-scale setup for each case by presenting U.S. NCEP (National Center for Environmental Prediction) surface reanalysis and 500-hPa height maps for North America. This is followed by plots of 48-h air mass back trajectories, derived from our 2.5km GEM-MACH-TEB output, for the middle afternoon of each case, which allows us to characterize the origin of air masses on a regional scale. Next, Section 3.4.2 presents an analysis of pollutant time series for 4 locations, along a south-to-north transect, in the GTHA, and Section 3.4.3 presents the high-resolution pollutant surface maps, vertical cross sections, high-resolution model wind fields and the observation-derived mesoscale meteorological analysis.

Synoptic-Scale Meteorology and Back Trajectories
The U.S. NCEP archives synoptic-scale reanalysis maps that show surface pressure systems and fronts each day. Geopotential heights on the 500 hPa pressure level were plotted to show the midtroposphere jet stream. The 27 July 2015 surface map [73] has a surface high-pressure system centered over Southern Ontario. The upper-level map has a large but weak ridge building into Southern Ontario. This ridge provided for some amplification, but only small transport of surface high pressure. Figure 7 shows the surface and upper air NCEP reanalysis for 28 July 2015. A double highpressure system has developed over Southern Ontario and Northern Pennsylvania. The upper level ridge has moved somewhat eastward into Ontario. This setup should bring air from the northern Midwest U.S. over the Great Lakes region and then SE into the GTHA. The upper-level forcing on 28 July was weak, so the surface high-pressure system did not move much throughout the daytime on 28 July. Air mass back trajectories were calculated using the 2.5-km GEM-MACH-TEB output wind fields. Figure 8 shows the 48-h back trajectory ending in Toronto at the North Toronto station at a height of 175-m at 20:00 UTC 28 July 2015. Every point along the trajectory is a 1-h increment. For the trajectory, the time series was initialized with a height of 175-m. The trajectory starts over Northern Michigan and travels southeast over unpopulated regions of Southern Ontario before arriving over Lake Ontario. The lake breeze circulation then brings the air mass northward over downtown Toronto before arriving at North Toronto. The SI section compares the 2.5-km GEM-MACH-TEBbased back trajectory with the U.S. NOAA HYSPLIT (Hybrid Single-Particle Lagrangian Integrated Trajectory) model back trajectories calculated using the 12-km NAM (North American Mesoscale) meteorological data and the 3-km HRRR (high resolution rapid refresh) meteorological data. The 2.5km GEM-MACH-TEB back trajectory is consistent with the HYSPLIT trajectories, both of which also originate over northern Michigan. The NCEP reanalysis map for 11 July has a surface high-pressure system over Northwestern Ohio. A weak ridge over Michigan provides for little amplification. Figure 9 shows that the surface high-pressure system has moved to Northern Pennsylvania on 12 July. The upper-level pattern on 12 July was a weak divergent zonal flow. This setup should bring air masses from the Midwest U.S. to the GTHA. The weak synoptic conditions are ideal for mesoscale features such as lake-breezes and urban heat islands to impact the local meteorology over the GTHA.
(a) Surface Meteorology (b) 500-mb Heights  Figure 10 shows a 48-h back trajectory ending in Toronto at 20:00 UTC on 12 July 2015 based on our 2.5-km GEM-MACH-TEB wind fields. The trajectory originates over central Michigan, travels east towards Lake Ontario, passes over Lake Ontario, turns north in a lake breeze circulation, passes over downtown Toronto and finally arrives at the North Toronto station. This back trajectory is quite similar to the back trajectories calculated with the HYSPLIT model using either the 12-km NAM wind data or the 4-km HRRR data (see the SI section). Note that selecting final altitudes between 50 and 300 m did not result in significant differences in the resulting back trajectories. A 96-h HYSPLIT backward trajectory was also run for each case using the GDAS (global data assimilation system) data. Both cases were similar in that the 96-h trajectories originated from the northern Midwestern U.S. region (i.e., state of Wisconsin) and travelled over sparsely populated terrain before arriving in the GTHA. The impact from pollution sources in the Southern and Eastern U.S. are thus likely minimal. All the trajectories suggest that for both cases the enhanced pollution observed at the North Toronto site in the late afternoon is likely photochemically produced from local sources within the GTHA. Figure 11a presents the O3 time series comparison for the 27-29 July 2015 period at the North Toronto and Newmarket measurement sites. On the afternoons of 27-28 July there was a large horizontal gradient in O3 mixing ratio as illustrated by the 50 ppbv difference in measured O3 between the North Toronto and Newmarket sites (a horizontal distance apart of only 30 km). This is consistent with a local influence from Toronto emissions, as also suggested by the backward trajectory for North Toronto. The model predicted the daytime O3 maxima quite well at North Toronto on 28 July and the NO2 maximum the night before, but under-predicted the maximum on 27 July and over-predicted the maximum on 29 July.

Case Study Time Series Analysis
At the time of the 1-h O3 exceedance in North Toronto (late afternoon on 28 July), model underpredicted at the downtown site and this under-prediction continued throughout the evening. A very stable, light-wind surface layer developed overnight and O3 remained under-predicted throughout the night for the downtown location. The O3 time series between 11 and 13 July for the four sites perpendicular to the lakeshore are shown in Figure 12. The daytime O3 maxima on July 11 were predicted well at all four sites. On 12 July, the O3 observed at North Toronto reached over 80 ppbv in the late afternoon. At this time (20:00 UTC), the model under-predicted O3 by 20 ppbv. The model predicted over 80 ppbv O3 at 20:00 UTC but for the Newmarket site. The measured O3 at Newmarket reached its maximum (73 ppbv) a couple of hours later (23:00 UTC) and the modeled O3 at Newmarket was predicted well at 23:00 UTC. The model also under-predicted at downtown Toronto on the afternoon of 18:00 UTC 12 July. The predictions at Toronto Island did not show a strong bias in the afternoons on 11-13 July; however, the observed O3 maxima did not exceed 65 ppbv.

Modeled Pollutant Spatial Distributions, Vertical Cross Sections and Meteorological Analysis 28 July 2015 Case Study
Few observations are available in the literature on the vertical extent of lake-breeze circulations associated with O3 exceedance periods. In this section, the modeled spatial distribution of pollutants was described for 21:00 UTC 28 July with a focus on the vertical distribution of pollutants. The SI section described the meteorology and chemistry leading up to 21:00 UTC (Figures S10 and S11 at 12:00 UTC, Figures S12 and S13 at 15:00 UTC and Figures S14 and S15 at 18:00 UTC). At the end of this case study description, an O3 sensitivity analysis was presented for 28 July at different points along the trajectory of the lake-breeze circulation. An important chemistry explanation for the O3 maximum along a lake-breeze front would be described. Figure 13 illustrates the temperature, relative humidity, vertical wind velocity near surface and the mixing length at 21:00 UTC (05:00 p.m. local time). The temperature over land reached a maximum up to 33 °C and the land/lake temperature difference was the maximum for the day. This temperature difference drove a strong lake-breeze front that reached as far inland as Vaughan to the north (30-km from shoreline), as illustrated by the vertical wind velocity in Figure 13c. The reader is referred to Figure 2 for the locations of cities. The vertical cross section of mixing height from Toronto Island (point A) to Newmarket (point E) shows uniformly well mixed convection all the way from point 'C' to point 'E' with mixing up to 2.3-km altitude. Horizontal winds were from the south over Toronto in the boundary layer. The lake-breeze circulation resulted in subsidence of air over Lake Ontario. A shallow layer of stable air remained over the lake surface, as shown by the dark blue color in mixing length in Figure 13d. Horizontal winds above the shallow surface layer were from the north and complete the lake breeze circulation. is terrain-following and can be linked to height above surface as follows: 1.0 is the surface, 0.95 is 500m above ground level (agl), 0.90 is 1-km agl and 0.70 is 3-km agl. The wind arrows in the vertical cross section do not show all the levels near surface for clarity. The arrows shown in the four panels and vertical cross-section represent wind in horizontal plane. Note that panel (d) shows a close-up over the west end of Lake Ontario.
A Doppler LIDAR was located at the North Toronto site on the afternoon of 28 July 2015. A description of the instrument, its deployment and the data image processing can be found in [33]. The Doppler LIDAR vertical velocity profile on 28 July at 20Z shows periods of strong up-draft followed by periods of down-draft (see Figure 6 in [27]). The up-drafts reached a maximum of 2.2km above ground at 4 pm local time, which is the time of the lake-breeze passage. This was consistent with the height of the boundary layer predicted by GEM-MACH-TEB model and provided an independent validation of the model mixing height. The magnitude of the updrafts was similar (1 m/s) between Figure 6 in [27] and Figure 13c.
The primary emitted pollutants, NOx and VOCs, show a uniform, well mixed region in the updraft, just south of the lake-breeze front ( Figure 14). The NOx mixing ratio is in the 3-6 ppbv range, which is optimum for efficiently producing O3. The lumped TOLU species is in the range 1-2 ppbv throughout the updraft region. O3 has its highest concentrations in North Toronto and Richmond Hill in the region of updraft with predicted O3 mixing ratios of 90-95 ppbv. The region of maximum O3 extends all the way up to 2.3-km altitude. The mixing of the O3 with the background air was evident as the plume moves back over the lake under the northerly flow above 1-km altitude, i.e., return flow of the lake-breeze circulation. The subsidence over the lake brings cleaner air from the upper boundary layer down to the residual layer above Lake Ontario and appeared to dilute the air in the layer from 250-m to 1-km altitude above points 'B' to 'C' over Toronto Island and part of downtown Toronto (see Figure 14b). Figure 15 shows the back trajectory for 21:00 UTC 28 July 2015 for the North Toronto site. It is clear that the air traveled over the city the night before and returned over downtown Toronto in the afternoon lake-breeze circulation.   Figure 16 shows the mesoscale analysis created using data from the expanded ECCC weather station network for the study, as well as satellite and radar data. The availability of the expanded station data enables an accurate diagnosis of the location of the lake-breeze front from the observations. The observation-derived lake-breeze front at 20:00 UTC is just south of Hwy 7 and just north of the position of the CRUISER mobile vehicle location at this time. Here, the CRUISER measured up to 90 ppbv O3, consistent with the model. The red star on Figure 16 denotes the CRUISER position at 20:00 UTC. By 21:00 UTC, the lake-breeze front is located just north of Hwy 7, passing through the city of Vaughan. Figure 14c shows the model-derived lake-breeze front on the lumped TOLU species panel (c). The gradient across the front is particularly striking for aromatic species because of their urban sources and short chemical lifetimes. The modeled lake-breeze front (red line) location is only slightly north of the observation-derived front (pink line). Both diagnosed fronts are north of the North Toronto site. The good agreement provides further validation of the modeled meteorology on 28 July 2015.
Ozone Production Sensitivity along the Lake-Breeze Air Mass Trajectory Rapid O3 production occurred on the afternoon of 28 July over Toronto, particularly in Uptown and North Toronto. It is informative to study the O3 production rate sensitivity to changes in precursor concentrations, as this will shed light on the optimal emission reduction strategies to reduce O3 mixing ratios. Kleinman et al. (2001) developed an analytical expression for the sensitivity of the O3 production rate to changes in precursor concentrations as a function of the fraction of radicals (OH + HO2 + RO + RO2) that are lost by reaction with NOx [74]. Figure 17 is a plot of the two sensitivity curves. The red curve is the O3 production sensitivity to NOx changes and the green curve is the sensitivity to VOC changes. The x-axis is the fraction of radicals lost by reaction with NOx relative to the total radical loss by reactions with NOx, HO2 and RO2 (termed Ln/Q). Crossing the ratio from less than 0.5 to more than 0.5 indicates a change from NOx-sensitive to VOC-sensitive O3 production. In looking at the curve for VOC changes, the value of O3 production sensitivity is always positive. This means that the O3 production rate always increases with VOC increases. However, in looking at the curve for NOx changes, the O3 production sensitivity is positive at low Ln/Q, but turns negative at a cross-over point of Ln/Q = 0.67 and becomes increasing negative at higher Ln/Q. Figure  17 and Table 3 present the results of the chemical sensitivity analysis using the model output at different locations, times and altitudes along the lake-breeze circulation. Figure 16. Mesoscale meteorological analysis at 20:00 UTC and 21:00 UTC 28 July 2015. The purple lines are the measurement-derived lake-breeze fronts. The wind data at stations is for a 10-m height. The colour scale is radar reflectivity. The red star in the first panel is the location of the ECCC mobile AQ laboratory (CRUISER) at 20:00 UTC when it measured 90 ppbv O3. Figure 17. Sensitivity of ozone production to changes in precursor concentrations (a value of +1.0 means that a 10% increase in precursor concentration (e.g., NOx, red curve) results in a 10% increase in ozone production, a value of +0.5 means that a 10% increase in precursor concentrations results in a 5% increase in the ozone production rate, a value of −0.5 means that a 10% increase in precursor concentrations results in a −5% decrease in ozone production rate). The right axis is the sensitivity to changes in volatile organic compound (VOC) concentrations (green curve). The points represent positions along the lake-breeze circulation (A0 is surface over North Toronto, A1.3 is at 1.3-km altitude over North Toronto, A2 represents 2-km altitude over North Toronto and B2 represents 2-km altitude over Uptown Toronto). Figure 17 is the location of the CRUISER (43.80° N, 79.42° W), near the North Toronto site, when recording the 90 ppbv O3 on at 20:00 UTC 28 July. Point A0 refers to the ground level at point A. Point A1.3 is at an altitude of 1.3-km above Point A. Point B is in Uptown Toronto (see Figure  2), just south of Hwy 401 (43.73° N, 79.40° W). Points A2 and B2 are at locations A and B, respectively, each at 2-km above ground. The analysis of the model reaction rates at point A0 resulted in a Ln/Q value of 0.84. This corresponds to the VOC-sensitive regime with a dlnPO3/dln(VOC) = 0.72. At this point, a 10% decrease in VOC concentration will result in a 7.2% decrease in O3 production rate. However, in this regime, a 10% decrease in NO concentration will result in a 4.5% increase in the O3 production rate. Above point A, at an altitude of 1.3 km, the vertical mixing is most intense resulting in a decrease in NOx concentration because of both dilution and chemical loss. Point A1.3 corresponds to an Ln/Q value of 0.58, just above the transition point (at 0.5) from VOC-sensitive to NOx-sensitive conditions. This point is close to the maximum O3 production rate for a given [VOC] reactivity (Ln/Q = 0.67 is the transition point from negative to positive value). The rapid O3 production rate in this lake-breeze-front updraft can help to further explain why O3 concentration maximizes along the lakebreeze convergence zone and why the high O3 concentrations extend all the way up to 2.5-km. At point A2, the updraft has now crossed over into the NOx-limited regime. At this point, a 10% decrease in NOx will result in a 6% decrease in O3 production rate. At point B2, there was significant dilution of plume with background air from the upper boundary layer. Point B2 had an Ln/Q value of only 0.07. The air mass was very sensitive to NO changes and relatively insensitive to VOC changes. The modeled NOx mixing ratios at each point in the lake-breeze circulations is also shown in Table 3. Table 3 also includes a column for downtown Toronto at surface (43.65° N, 79.38° W). The NOx mixing ratios spanned the range from 12 ppbv downtown at surface to 0.56 ppbv uptown in the upper boundary layer. Table 3. Chemical coordinate analysis results for ozone production sensitivity to precursor concentration changes for different locations in lake-breeze circulation on 28 July 2015.  Table 3 includes two other modeled species ratio indicators. The H2O2/HNO3 ratio is a weighted measure for the O3 production sensitivity (NOx vs. VOC sensitive) accumulated along the history of the backward trajectory, i.e., how much of the O3 production in an air mass has been from VOC vs. NOx sensitive chemistry. A high H2O2/HNO3 ratio (>0.6, by mass) indicates a NOx-sensitive O3 production regime and a low H2O2/HNO3 ratio (<0.3, by mass) indicates a VOC-sensitive O3 production regime [75,76]. In Table 3, at the surface in downtown and North Toronto, the ratio indicated a transition regime. At the North Toronto location, in the updraft, the O3 production was VOC-sensitive and in the return flow at 2-km altitude over Uptown Toronto, the O3 production turned to NOx-sensitive.

Chemical Sensitivity
There has been one recent chemical transport model study that looked at O3 production sensitivity to precursor species. Whaley et al. (2015) used the GEOS-Chem global model to study the transport patterns for pollution events in Toronto from 2004-2007 [77]. They found that for the midto long-range transport cases, the O3 maximum in Toronto was sensitive to regional-scale NOx and isoprene emissions. For cases with stronger local O3 production, there were high NOx concentrations in Toronto and the local O3 production was sensitive to anthropogenic VOC emissions in Toronto. Our high-resolution model result for light wind conditions show similar precursor sensitivity as in [77].
Finally, there are few data on OVOC contributions to total VOC reactivity for urban areas due to the challenge of measuring polar OVOCs. The OVOC/VOC ratio is the fraction of the total VOC reactivity that results from OVOC model species reacting with the OH radical. Here, we present the modeled contributions of OVOCs to the total VOC reactivity. The fraction ranged from a low of 0.28 at the downtown site to a value of 0.37 in the updraft at 2-km over North Toronto. These large fractions suggest that OVOCs are important precursors to O3 production and the monitoring community should continue to improve network measurement methods for OVOCs.

July 2015 Case Study
The 2.5-km GEM-MACH-TEB model had a 'bust' forecast for the afternoon of 12 July 2015 when the measured O3 exceeded 80 ppbv at the Toronto North site. In this section, we present the modelpredicted spatial maps for both air quality and meteorology fields, compare them to the observationderived mesoscale meteorological analysis, and try to answer the question as to why the model under-performed. Figure 18 presents the modeled O3 and NO2 surface spatial distributions at 15:00 UTC, 18:00 UTC and 21:00 UTC. The black box encompassing the downtown and North Toronto is a key region of focus. This same box is shown on the meteorological mesoscale analysis (Figure 19), and the model 10-m wind field panels ( Figure S16). At 15:00 UTC, the model predicted a moderate southwesterly wind over Lake Ontario. Southwesterly winds also characterized most of the key region (black box) with the exception of the northwest corner of the box, where the model winds were more southerly, and the most southern part of the box, where a lake breeze front was just pushing on-shore. NO2 predictions were highest over the southern portion of the box due to the urban emissions. By 18:00 UTC, the model predicted the southwesterly flow for the left half of the box and moderate southerly lake breeze flow on the right half of box. The modeled lake-breeze front arrived at Newmarket, as shown by the jump in O3 mixing ratio in the time series (Figure 12), just after 18:00 UTC. The region of maximum ozone production was already north of Toronto. By 21:00 UTC, the model predicted over 90 ppbv O3 north of Newmarket as the modeled lake-breeze front reached all the way to Lake Simcoe. We next compared the model meteorology with observations to understand what caused the 'bust' O3 forecast. Figure 19 shows the observation-derived mesoscale meteorological analysis from the over 40 meteorological stations. Key stations within the box were A4T (Brampton) in the west of the box, YYZ (Pearson International Airport) in the southwest of the box, King City (L2D) and Vaughan (A2T) in the north of the box and Concord (L1F) and Downsview Park (L1E) in the center of the box. At 15:00 UTC, Brampton measured the southwesterly flow, King City measured the westerly flow and the stations in the center of the box measured the light northerly flows. Pearson International Airport also measured the light northerly flow. By 17:00 UTC, the observed lake-breeze front had moved past Downsview Park and was south of Concord station. The lake-breeze front had also passed Pearson International Airport, while Brampton and King City were still reporting a southwesterly flow. By 18:00 UTC, the observed lake-breeze front had slowed its progression and was still located south of Concord. By 19:00 UTC, the front had stalled. However, by 21:00 UTC, the front had continued its progression and was now located just south of King City. By 22:00 UTC, the front had passed King City and was over Newmarket. The passage of the front at Newmarket on the mesoscale meteorological analysis was coincident with the measured O3 maximum at Newmarket on Figure 12. Figure 19. Mesoscale meteorological analysis for 12 July 2015. The box highlights the same region as Figure 18. Wind barbs report 10-m station data. The purple front is the observation-derived lake breeze front. Figure S16 shows the detailed predicted surface wind fields at different times. These panels were used with the mesoscale meteorological analysis to create Table 4, which compares the times of lakebreeze passage from the model and from observations along the south-north line from the lakeshore. There was a 2-h difference in the timing of the lake-breeze front passing Concord (center of box on plots). By the time the front reached Newmarket, there was a 3-h difference. The consequence of this difference in timing is that the model predicted its daytime O3 maximum too far north (Newmarket) compared to the observed location of the O3 maximum (North Toronto). Figure 20 shows a 24-h airmass back trajectory arriving at Newmarket at 21:00 UTC. It is clear that the origin of the air mass earlier in the day was over Lake Ontario. A key question to answer is why does the model predict a stronger lake-breeze flow on this afternoon? The model could differ in two potential ways. The modeled regional scale wind pattern opposing the lake-breeze front could differ from observations and/or the modeled lake/land temperature contrast could differ from observations. Figure S16 shows the predicted wind fields at 16:00-17:00 UTC. The model predicted a southwesterly flow coming into the box and southwesterly flow north of Toronto, whereas the mesoscale meteorological analysis had a weak northerly flow north of Toronto at Concord (L1F) and Vaughan (A2T) at these times. Thus, in reality on the regional scale, there was a more opposing flow to the passage of lake-breeze front. Table 5 presents the modelled values for lake water (Tw) and lake air (Tair) temperatures compared to observations. There was a Watch Keeper buoy just south of Toronto Island, which measured Tw. The modeled Tw was almost 2 ˚C colder than that measured on the afternoon of 12 July. There was also a buoy in Southwestern Lake Ontario that recorded Tair over 2 °C warmer than the model predictions, and there was a boat deployed at a location in Western Lake Ontario that measured Tair 1.5 °C warmer than the model at its location. Finally, Toronto Island airport measured Tair slightly higher (0.3-1.0 °C) than model values on this afternoon. The wind speeds modeled at Toronto Island Airport (YTZ) were 7-9 knots (3.6-4.6 m/s) between 15:00 and 17:00 UTC. The wind speed observations at YTZ were 4-6 knots (2.0-3.1 m/s) between 15:00 and 17:00 UTC. Collectively, these results suggest the modeled Tw and Tair were cooler than the measurements in the western part of Lake Ontario. The surface scheme in 2.5-km GEM-MACH-TEB used the 2-km coupled NEMO atmosphere-lake model output over the Great Lakes. NEMO Tair was initialized from the 10-km GEM regional model meteorology [42,43]. The 10-km regional GEM model was initiated each cycle with a 10-km surface analysis, including lake Tw. A 6-h spin-up of the 2-km NEMO was used before extracting the lake Tw field for the driving 2.5-km GEM-MACH-TEB surface.
The satellite imagery on the mesoscale meteorological analysis maps suggests 12 July was a partly cloudy day over Lake Ontario. The heterogeneity in surface solar heating from broken cloud cover may have been a challenge for the NEMO model to represent on this afternoon. On 28 July, there were clear skies all morning and afternoon over Lake Ontario. In conclusion, model comparison with observations suggest the stronger modeled lake breeze on 12 July and the different timing of the front passage is due to a weak but opposing wind flow north of Toronto and a larger modeled difference in the land-lake temperature. Further improvements to the model system should consider a high-resolution, more frequent lake-water-temperature-analysis based on satellite observations. The recent McNider et al. study showed a strong sensitivity of air quality predictions on their satellite data-assimilated lake water temperatures [19].

Recommendations for Emission Reduction Strategies
Geddes et al. (2009) found that simultaneous reductions in observed mixing ratios of both VOCs and NOx in Toronto from 2000 to 2007 did not result in a significant change in observed O3 mixing ratios [25]. They showed using a chemical box model that this weak response was not unexpected, due to the non-linear nature of O3 production. They also showed that, for average conditions in downtown Toronto, the O3 production is VOC-limited. They recommended that VOC emission reductions would be most effective at decreasing local O3 production rates. One caveat, in this prior study, was the use of 24-h VOC canister data as no high-temporal-resolution VOC data sets were available for the photochemically active times of the day. The OVOC data were also limited. They also noted the challenge in the future in reducing VOC levels as a significant fraction are biogenic in nature and biogenic emissions will likely increase in a warming climate. Makar et al. (2010) examined the sensitivities of different VOCs towards ozone formation, noting that the products of aromatic oxidation had the largest impact on urban ozone concentrations in the Great Lakes region and suggested that a reduction in on-road and shipping VOC emissions (the largest local sources of aromatics) could be one route for reducing ozone concentrations [15].
The two 1-h exceedance periods presented here are remarkably similar in terms of the meteorological setup on the synoptic scale. Based on the back trajectories, the model suggests the production of O3 in North Toronto on these days is largely from local production. Differences in the extent of northern encroachment of the modelled lake-breeze front on each day determines which northern suburbs are affected by the maximum O3, but both cases predict rapid production of Ox up to 90 ppbv at the front. The chemical sensitivity analysis suggests that the O3 production is VOCsensitive. Thus, reductions in VOC emissions over Toronto along the path from downtown to North Toronto will be most effective in reducing the rapid O3 production rates. Indeed, conditions are NOx saturated, in that reductions in NOx emissions will likely increase O3 production rates along the case study trajectories over Toronto. However, as the Toronto air mass mixes with background air from aloft, the plume will turn from VOC-sensitive to NOx-sensitive. It seems the most prudent strategy is to reduce both urban VOC and NOx emissions, but to reduce the urban VOC emissions more aggressively so that urban O3 production is not enhanced further.

Conclusions
The first goal of this work was to evaluate 2.5-km GEM-MACH-TEB model predictions against observations for the daily maximum O3 metric in the GTHA. The 2015 Pan American Games had 4 days in the month of July with a daily maximum 8-h average (DM8A) O3 mixing ratio above the Canadian standard of 63 ppbv and 2 days with O3 above the Ontario daily 1-h standard of 80 ppbv. Overall, the model performed well, with its O3 predictions meeting the performance criteria set by the U.S. EPA. During the study, an additional measurement site was located on Toronto Island about 2-km off-shore from downtown Toronto. The model did not over-predict O3 at this location, unlike previous model studies, including previous versions of GEM-MACH. During the study, the two highest 1-h O3 measurements were observed farther inland along the lake-breeze front in the late afternoon. The two periods when O3 levels were above the 1-h standard were analyzed in detail. For 12 July, at midday, GEM-MACH-TEB predicted a strong lake-breeze front that pushed northward and that by 18:00 UTC was located north of Newmarket. However, the observationderived mesoscale analysis diagnosed the lake-breeze front location at 18:00 UTC to be about 30 km south of Newmarket, near the town of Vaughan. The northwesterly flow, observed at the A2T station in Vaughan, also disagrees with the strong modeled southerly flow at this location. By 21:00 UTC, the modelled lake-breeze front was predicted to have reached Lake Simcoe, whereas the mesoscale analysis diagnosed the lake-breeze-front location as still south of Newmarket. This case illustrates the challenge in predicting the in-land penetration of a lake-breeze front when the synoptic forcing is weak and there is a tenuous balance of mesoscale forces resulting from land/lake temperature differences, the urban heat island, and topography changes. The model did predict the fine-scale structure and resulting O3 maximum, but the structure was displaced northward. We identified modelled lake air and water temperatures that were too cool as a factor contributing to the difference in strength and timing of the lake-breeze front.
For 28 July, the coupled high-resolution meteorology/chemistry/urban canopy model was successful at predicting the location of the daytime O3 maximum because the model was able to predict the location of the lake-breeze front in the afternoon. An important factor for the O3 maximum along lake-breeze fronts was discovered when performing the O3 production sensitivity analysis at different points along the trajectory of the lake-breeze circulation. In the shallow inflow region of the lake-breeze circulation over downtown Toronto, O3 production is VOC-limited and suppressed by the high NOx conditions. However, in the updraft region of the lake-breeze frontal zone, the primary pollutant concentrations begin to dilute and the O3 production begins to transition from VOCsensitive to NOx-sensitive. In this transition zone, O3 production efficiency is a maximum and the still relatively high urban VOC and NOx absolute concentrations result in large O3 production rates and the maximum in O3 mixing ratios.
In general, the model results from the Pan American Games period showed that the highresolution urbanized GEM-MACH-TEB model performed better in statistical scores (RMSE, MB) compared to the 10-km GEM-MACH model for O3 and NO2. The 10-km model performed slightly better for the correlation coefficient (R), likely due to the challenge in modelling the correct location of more spatially resolved pollution plumes. Future work will evaluate the high-resolution urbanized model, run in near-real time, for an entire year period against ECCC's operational air quality model.

Supplementary Materials:
The following are available online at www.mdpi.com/xxx/s1, Table S1: Model evaluation for hourly data at 12 GTHA sites, Figure S1: Monthly Mean Maps, Table S2: Model evaluation for 1h O3 maximum at 4 GTHA sites, Figure S2: Correlation North Toronto, Figure S3: Correlation downtown Toronto, Table S3: Model evaluation for 8-h O3 at 4 GTHA sites, Figure S4: Ox Bias for Different Wind Directions, Figure S5: HYSPLIT back trajectory for 28 July with NAMS, Figure S6: HYSPLIT back trajectory for 28 July with HRRR, Figure S7: HYSPLIT back trajectory for 12 July with NAMS, Figure S8: HYSPLIT back trajectory for 12 July with HRRR, Figure S9 Table A1 lists the standard biogenic emission rates updated in the BEIS program [78,79]. The MEGAN v2.10 (Model of Emissions of Gases and Aerosol from Nature) emission rates were used for Birch, Larch, Balsam Poplar, Aspen and other Populus tree types. This resulted in a significant increase for 'Monoterpene' species. The dominant tree species in the Canadian boreal forest are pine and spruce trees. We reduced the standard emission rate for these tree species by a factor of 2 to be consistent with the satellite-derived LAI for this ecosystem.