Downscaling of air pollutants in Europe using uEMEP_v6

. The air quality downscaling model uEMEP and its combination with the EMEP MSC-W chemical transport model are used here to achieve high-resolution air quality modeling at street level in Europe. By using publicly available proxy data, this uEMEP/EMEP modelling system is applied to calculate annual mean NO 2 , PM 2 . 5 , PM 10 and O 3 concentrations for all of Europe down to 100 m resolution and is validated against all available Airbase monitoring stations in Europe at 25 m resolution. Downscaling is carried out on annual mean concentrations, requiring special attention to non-linear processes, such 5 as NO 2 chemistry, where frequency distributions are applied to better represent the non-linear NO 2 chemistry. The downscaling shows signiﬁcant improvement in NO 2 concentrations where spatial correlation has been doubled for most countries and bias reduced from -46% to -18% for all stations in Europe. The downscaling of PM 2 . 5 and PM 10 does not show improvement in spatial correlation but does reduce the overall bias in the European calculations from -21% to -11% and from -39% to -30% for PM 2 . 5 and PM 10 respectively. There is improved spatial correlation in most countries after downscaling of O 3 , and a reduced 10 positive bias of O 3 concentrations from +16% to +11%. Sensitivity tests in Norway show that improvements in the emission and emission proxy data used for the downscaling can signiﬁcantly improve both the NO 2 and PM results. The downscaling development opens the way for improved exposure estimates, improved assessment of emissions as well as detailed calculations of source contributions to exceedances in a consistent way for all of Europe at high resolution.

-Resolution of the sub-grids varies according to the application but maps are made at 100 m and calculations at monitoring sites are made at 25 m. based on the official data submissions to EMEP Centre on Emission Inventories and Projections (Pinterits et al., 2020) in 2020, in which the PM emissions from the residential combustion sector (GNFR C) are replaced by a bottom-up estimate of TNO (Denier van der Gon et al., 2020;Fagerli et al., 2020) for 2017. This TNO dataset should represent an improved estimate of residential combustion emissions of PM, accounting for condensable organics in a consistent way.
The EMEP model calculates and outputs the 'local fraction', used by the uEMEP downscaling to remove double counting 70 of emissions (Denby et al., 2020;. The local fraction is the contribution of emissions in one EMEP grid to itself and to its neighbouring grids. For this application only a 3 x 3 grid contribution region is calculated, though for other applications this can be much larger. By tagging the grid emissions in this way the local contribution from EMEP can be removed and replaced by the high-resolution uEMEP sub-grid calculation. 75 The uEMEP model is described in a recent publication (Denby et al., 2020). In that article the Norwegian forecast application of uEMEP is described where hourly downscaling using bottom-up emission inventories is carried out. For the European application calculations are made on annual mean data, creating air quality maps for Europe down to 100 m resolution and calculating concentrations at Airbase stations positions down to 25 m.

uEMEP model implementation
Downscaling is carried out in the following way. EMEP grid emissions per sector and per source are redistributed to uEMEP 80 sub-grids using the proxy emission data described in Sect. 2.3. Each sub-grid dispersion calculation uses only the emission sub-grids within the moving window emission region, see Fig. 1. For annual mean calculations dispersion is carried out using a rotationally symmetric Gaussian dispersion kernel (Denby et al., 2020), given an initial plume size and height. These parameters are provided in Table 1. The initial horizontal plume size is determined by the size of the sub-grid. The Gaussian dispersion parameters used are based on the K z dispersion methodology described in Denby et al. (2020) but adapted to the rotationally 85 symmetric dispersion kernel. The EMEP model local fraction contribution originating from within the moving window region is removed and replaced with the sub-grid dispersion calculation from uEMEP, thus avoiding double counting of emissions ( Fig. 1). For these simulations this region corresponds to 2 x 2 EMEP grids, i.e. within an area that is ± 0.1 o in both latitude and longitude. This ensures that no matter where the uEMEP calculation sub-grid is placed that the moving window region will always be covered by the 3 x 3 local fraction region. 90 Downscaling with uEMEP occurs only for primary emissions. NO 2 is calculated from NO X and O X using the same travel time parameterisation described in Denby et al. (2020) but applied to annual mean wind speeds, photo-dissociation rates and concentrations. To account for non-linearity in the NO 2 chemistry, when calculating with annual means, an additional frequency distribution correction factor is implemented, see Sect. 2.4. Annual mean downscaled O 3 is also determined using this same parameterisation. 95 To improve efficiency of the calculations, Europe is split into a number of tiles that cover the European land domain. For the 100 m resolution mapping calculations there are 1097 tiles, each of which is 100 km x 100 km. These tiling regions are shown in Fig. 2. Table 1. Initial dispersion (σz0) and emission height (hemis) for the three downscaled sources for all primary pollutants.

uEMEP proxy data
We use road data from OSM to redistribute the traffic emissions. Though the spatial coverage of OSM is very good, it does not 100 contain actual traffic data. Redistribution of the emission data is achieved by weighting the different road categories provided in OSM. The following road categories are considered: motorway, trunk, primary, secondary, tertiary, unclassified and residential.
Each is weighted relative to the other so that emissions can be redistributed and attributed to the road links. Estimates of the weights are based on the representative average daily traffic (ADT) for different road categories for Norwegian average road situations. The weighting currently employed for the redistribution is shown in Fig. 3. It is also worth noting that for major 105 roads, such as motorways, OSM often represents these as dual carriageways, i.e. as two separate road links. In these cases the weighting of a motorway will be twice that indicated here. Sensitivity tests with alternative weighting, see Sect. 5.2, show the choice of weighting does impact on the results but that the current choice provides close to optimal spatial correlation when compared to measurements.
A global population dataset from the GHS is used as the proxy for redistributing residential heating emissions. We choose the 110 highest available resolution of 9 arcsec (0.0025 o ) from the year 2015. The coordinate system is WGS84. This dataset indicates emissions we still use them as proxy data to redistribute EMEP gridded emissions, to be consistent in the methodology.
2.4 uEMEP chemistry parameterisation for annual mean NO 2 uEMEP downscales only primary pollutants. It is thus necessary to apply chemistry parameterisations to the NO X and EMEP O 3 concentrations to derive downscaled NO 2 and O 3 concentrations. Two methods for doing this are described in Denby et al. (2020). One for hourly concentrations using a weighted travel time parcel method and the other for annual means, which is 120 based on a simple empirical relationship between observed NO 2 and NO X . It is desirable to apply the model based chemistry scheme rather than an empirical scheme, however due to the non-linearity of the NO 2 chemistry the chemical scheme cannot be directly applied to annual mean concentrations.
To solve this, chemistry is not calculated on just a single annual mean value for NO X and O 3 but on a frequency distribution for these parameters that represent the variability over a year. This can be illustrated for the photostationary case where the 125 NO 2 concentrations can be derived from NO X and O X using where the concentrations are annual mean values in molecules/cm 3 , k 1 is the production reaction rate and J the photodissociation rate for NO 2 . If the frequency distribution for the three annual mean variables [NO X ], [O X ] and J/k 1 is known 130 then we can integrate over Eq. (1) using these distributions as weighting functions. An appropriate probability distribution function for the concentrations is the log-normal distribution, which can be written as where the log-normal parameters µ and σ are determined from the mean (m) and the standard deviation (s) by The frequency distribution of J is not log-normally distributed, since it is dependent on the solar zenith angle (ZA) and various other meteorological parameters, such as cloud cover and water vapour content. The EMEP model uses lookup tables 140 based on precalculated J values from the Phodis model (Jonson et al., 2000). In order to implement the frequency distribution of J in uEMEP a power law fit is made to the tabulated values. This can be written as: where C j is a constant that is normalised out when producing the normalised frequency distribution and p j = 0.28. The standard 145 deviation of k 1 , dependent on air temperature, is significantly smaller than for J so it is treated as a constant.
A new value [NO 2 ] pdf can then be determined using these frequency distributions and a correction term showing its difference from the mean is defined as When calculating in three orthogonal dimensions then it is assumed there is no correlation between the variables.
To implement this procedure the standard deviation s must be known for NO X and O X . Values for s nox and s ox have been derived from earlier model calculations for Norwegian stations. Linear regression provides robust values for s ox /m ox and 155 s nox /m nox of 0.21 and 1.14 respectively (see Fig. S5). The variability of NO X reflects the variability of the traffic emissions, for stations within the influence of traffic, and this should be generally applicable throughout Europe. 72 sites are used for the calculation. The calculation of the distribution correction is carried out numerically after calculation of the NO 2 concentrations.
Implementation of the frequency distribution for concentrations has a significant impact, with a general reduction in NO 2 concentrations compared to the annual mean calculation using Eq. 1. This reduction leads to correction terms (f no2,pdf ) of 160 between 0 to -25%. The highest corrections occur around NO X =O X , where Eq. (1) shows the most non-linear behaviour. On average for all station sites in Europe, around a -16% reduction on the initially calculated [NO 2 ] mean has been determined. In contrast to the distribution correction for concentrations, the distribution correction for J leads to an increase in NO 2 of around 6%. This is because roughly half of the frequency distribution for J is 0, i.e night time, when there is no photo-dissociation.
In Sect. 5.5 the impact of this and other chemistry schemes is further discussed. More information concerning this scheme is 165 contained in the supplementary material S1.4. In this section we present example maps that are generated from the EMEP model and uEMEP simulations. 100 km example tiles are shown in Fig. 4 -6, demonstrating the original EMEP model calculations and the downscaled maps using uEMEP for shown as circles. This is the same tile as is shown in Fig. 4.   Lefebvre et al. (2013) in Antwerp, where both Gaussian and street canyon models were applied at 15 street canyon modelling sites, showed an average street canyon modelling increment of just 11% for NO 2 . We include all available sites in this study because in Europe the 185 majority of traffic sites appear not to be street canyons, though information on this is unclear (Tarrasón et al., 2021). Also, even without including obstacles, the increased model resolution (up to 25 m) allows the concentration gradients at roadside to be better described.

NO 2
In Fig. S6 and Fig. S7 scatter plots for NO 2 are shown for each country and Europe as a whole. These results are summarised 190 in Fig. 7 where the annual mean concentration and spatial correlation are shown.
In a majority of countries the spatial correlation for NO 2 is more than doubled when implementing uEMEP. The two exceptions are Ireland (IE), where the spatial correlation hardly changes with the downscaling, and Bosnia and Herzegovina (BA), where the spatial correlation is significantly reduced. Both these countries have very few stations. The highest spatial correlation is for Poland (PL) with r 2 = 0.85.

195
It is worth noting that the average spatial correlation per country is r 2 = 0.62 which is higher than the spatial correlation when assessed for all stations in Europe (r 2 = 0.57). This indicates that a part of the variability occurs between countries and can be interpreted to reflect differences related to emission reporting from each country. If the NO X emissions from individual countries have uncorrelated bias then this will reduce the overall spatial correlation.
The relative bias is improved for all countries with the exception of Greece (EL), which is the only country with a significant 200 positive bias. Overall for Europe bias is improved from -46% for the EMEP model to -18% when using uEMEP. Of the 28 countries with 10 or more monitoring sites, 18 of these have an absolute bias less than 25% after downscaling. Turkey (TR) has the largest negative bias, after downscaling, of -45%.

PM 2.5
In Fig. S8 and Fig. S9 scatter plots for PM 2.5 are shown for each country and Europe as a whole. These results are summarised 205 in Fig. 8 where the annual mean concentration and spatial correlation are presented.
Unlike the NO 2 downscaling, there is generally no improvement in the spatial correlation when applying uEMEP for PM 2.5 .
Only 6 out of 17 countries show improved spatial correlation and overall for Europe there is a slight decrease, from r 2 = 0.49 for the EMEP model to 0.46 for uEMEP. This result is further discussed in Sect. 6.    The relative bias is however reduced for almost all countries. For Europe as a whole the relative bias went from -21% for 210 the EMEP model to -11% for uEMEP. Only the three countries Austria (AT), Sweden (SE) and Finland (FI), that had almost no bias with the the EMEP model calculation, achieve a positive bias with uEMEP.

PM 10
In Fig. S10 and Fig. S11 scatter plots for PM 10 are shown for each country and Europe as a whole. These are summarised in Fig. 9.

215
The results for PM 10 are similar to those for PM 2.5 . In this case though the majority of countries, 21 out of 27, have improved spatial correlation with the application of uEMEP. The spatial correlation for all of Europe using uEMEP is unaltered compared to the EMEP model calculation, with r 2 =0.34. This is lower than the spatial correlation found for PM 2.5 by around 0.12. As with PM 2.5 the relative bias is reduced with the uEMEP downscaling. For Europe we see the relative bias went from -39% for the EMEP model to -30% for uEMEP. Ozone is generally reduced with the downscaling due to an increase in NO X concentrations. In general for Europe we see a reduced positive bias from +16% for the EMEP model to +11% for uEMEP. Spatial correlation is also improved in 21 of the 225 24 countries. The only countries to show significant degradation in the downscaling results are Switzerland (CH) and Greece (EL). This is likely due to the overestimated NO X concentrations there (Fig. 7).

Sensitivity studies
In this Section we present the results of several sensitivity calculations using uEMEP. These include sensitivity to sub-grid resolution, traffic emission proxies, residential combustion emission proxies, sensitivity to alternative bottom up emissions in 230 Norway and sensitivity to the NO 2 chemistry scheme.

Sensitivity to resolution
When calculating concentrations at station positions a grid resolution of 25 m is used. However, when mapping all of Europe a lower resolution of 100 m is employed. In Fig. 11 we show the results of a change in resolution on the annual mean NO 2 and PM 2.5 concentrations for resolutions from 25 to 500 m. PM 2.5 r 2 due to changes in uEMEP grid resolution

Sensitivity to OSM weighting 240
In Fig. 3 the weighting imposed on the OSM road categories is shown. This weighting specifies the relative contribution of the different road categories to the redistribution of the gridded traffic emissions in uEMEP. This weighting is based on an analysis of Norwegian traffic data but it is worthwhile assessing the sensitivity of the calculated NO 2 concentrations using different weights. To assess this sensitivity a power law is applied to the weighting. For power indices greater than 1 then more weighting is applied to the major roads, for power law indices less than 1 then more weight is applied to the the minor roads.

245
The different weights for the different power law indices are shown in Fig. 12. The results of this sensitivity test, presented in terms of relative bias and spatial correlation (r 2 ), are shown in Fig. 13.
Bias is quite strongly affected by the change in weighting. Higher concentrations are calculated when more weight is given to the minor roads. This is likely because most measurement sites are not on major roads. Increasing the weighting to minor roads will generally increase the urban background levels. The spatial correlation is among the highest for the current weighting with a power index of 1. This confirms that the initial estimate, based on Norwegian traffic, reflects a good general distribution of traffic in Europe. If real traffic volume were available then the weighting would be more precise. Tests on Norwegian data, Sect. 5.4, confirm that spatial correlation is significantly improved when using real traffic data for the redistribution weighting.

Sensitivity to the residential combustion emission proxy
For the PM 2.5 calculations presented in Sect. 4.2 population density data at 0.0025 o has been used to redistribute the residential 255 combustion emissions in uEMEP. The results indicate a slightly reduced spatial correlation but also with an improved negative bias. In this section we assess the sensitivity of the redistribution proxy to a number of alternative proxies. Firstly a power law is applied to the population density data. A lower power law index will reduce the weighting towards highly populated regions. PM 2.5 r 2 for residential combustion proxy parameters Figure 14. Change in bias and spatial correlation (coefficient of determination r 2 ) for PM2.5 as a result of changes in the residential combustion proxy. A lower power law index gives less weight to the population redistribution, a higher power law indices give more weight.
'Building/population' is the building density masked by population data. See text for details.
A power law index of 0 will work as a mask, redistributing the EMEP emissions evenly to any 250 m sub-grid that contains population. As an alternative to the population data, building density data has also been extracted from the OpenStreetMap 260 dataset. This has also been placed on a 0.0025 o grid for all of Europe. Two alternatives with this proxy are tested. The first using building density as the weighting proxy and the second using building density masked with population, so that only areas with both buildings and population are used for redistribution. In addition to the alternative proxy data the sensitivity of the calculations to emission height, currently set to 15 m, is also assessed.
The results are shown in Fig. 14. Here we see that a power law of 0.25 gives slightly improved spatial correlation and that 265 the use of building density also slightly improves spatial correlation compared to population. However, none of the alternative proxies significantly improves the spatial distribution of PM 2.5 and none attain the spatial correlation of the EMEP model calculations at 0.1 o . There is a general trend for reduced negative bias to lead to reduced spatial correlation in all calculations, so when the contribution from the downscaled residential combustion increases then spatial correlation decreases. This infers that the redistribution is not improving the results.

270
In addition to the proxy sensitivity the result of the EMEP model calculation where all local EMEP model grid contributions (± 1 o ) have been removed is shown in Fig. 14. This shows firstly that around 10% of the PM 2.5 in the EMEP model comes from within this local region and that the inclusion of these emissions does add to improved spatial correlation at the EMEP model 0.1 o scale, from r 2 = 0.467 to 0.488. Here we see more clearly that while the bias is improved by downscaling the spatial correlation is not and is similar to the spatial correlation obtained from the non-local contributions. However, it is possible to 275 achieve improved spatial correlation when more appropriate downscaling proxies are used. This is presented in Sect. 5.4.

Results of improved emission data in Norway
Throughout the uEMEP downscaling simulations we used the 0.1 o country reported emission data and redistributed it using population, OpenStreetMap data and AIS shipping data as redistribution proxies. However, many countries have more detailed emission data sets, including Norway, that could be used to improve the downscaling calculations. To test the impact of more We make four separate downscaling calculations for Norway using the two emissions, 'European emissions' and 'Norwegian emissions', and the two high-resolution proxy datasets, 'European proxy downscaling' and 'Norwegian proxy downscaling', respectively. Shipping is not changed in these simulations and in this case the calculation year is 2017. Though the resolution 295 of the EMEP model in the Norwegian forecasting system is nominally 2.5 km, for these simulations we use the same 0.1 o EMEP model grid resolution. The results are shown in Fig. 15 for NO 2 , PM 2.5 and PM 10 where the relative bias (%) and spatial correlation (r 2 ) are presented.
For NO 2 in Norway the large negative bias seen in the EMEP model is almost completely removed by the use of the traffic downscaling, using either the Norwegian or European emission data. On a national level the local Norwegian (bottom up) traffic 300 NO X emissions are roughly 25% higher than the EMEP (top down) emissions. NO 2 concentrations are slightly overestimated when using the Norwegian proxy data for traffic. Spatial correlation is improved with the use of the Norwegian proxy data for traffic, compared to European emissions that use OSM data, from r 2 = 0.6 to 0.72. It is worth noting that in the complete Norwegian calculation reported in Denby et al. (2020) using hourly calculations that the spatial correlation is even higher at r 2 = 0.78, but the bias is less at -5%.

305
For PM 2.5 biases are very similar for both the European and Norwegian proxy data sets when using either the European or Norwegian emissions. The spatial correlation however is significantly higher when using the Norwegian emissions, both at grid level and after downscaling. There is significant improvement, r 2 increases from 0.37 to 0.55, when both changing European emission to Norwegian emission and changing the residential heating proxy from population (European proxy) to the MetVed model (Norwegian proxy). This indicates that improved spatial representation can be attained when both the gridded and the 310 proxy data are consistent and more representative. However, little can be improved with downscaling when the initial gridded emissions are not well distributed, even with improved proxy data. Interestingly we see the same result as reported in Sect.
4.2, that the spatial correlation is reduced when applying the European proxy data to the European emissions. These results  indicate that significant improvements can still be obtained in the downscaling if improved emissions and emission proxies are implemented.

315
For PM 10 both bias and spatial correlation are significantly improved with the implementation of the local emissions and proxies. This is to a large extent due to the improvement in the road dust emission contribution but also due to an improvement in the residential heating distribution. Spatial correlation is also significantly increased, from r 2 = 0.27 to 0.49.

Sensitivity to the NO 2 chemistry scheme
Included in uEMEP are a number of simplified NO 2 chemistry schemes, used to derive downscaled NO 2 concentrations from 320 NO X and O 3 concentrations. In the results presented so far we have used the weighted travel time parcel method, as applied and described in Denby et al. (2020), with the additional use of the frequency distribution correction described in Sect. 2.4 ("Travel time" in Fig. 16). Two additional chemistry-based schemes and two empirically-based schemes are also available.
The two alternative chemistry schemes are the photo-stationary formulation ("Photostationary" in Fig. 16) and an alternative stationary formulation ("Stationary" in Fig. 16) that also allows for deviation from the photo-stationary state (Maiheu et al.,325 2017). The first empirical scheme is the Romberg scheme  ("Romberg" in Fig. 16), also described in Denby et al. (2020), that directly converts NO X to NO 2 concentrations. The parameters for this equation have been updated by fitting to all available Airbase data for the year 2017. The other empirical formulation is the SRM scheme (Wesseling and van Velze, 2014) ("SRM" in Fig. 16) that is also based on a fit to measurement data but includes background O 3 as one of the input parameters. The advantage of the two empirical fits is that they should convert NO X to NO 2 in a manner that is consistent with 330 the observations, and as such can be applied to annual mean concentrations directly, without any correction for non-linearity.
All methods are described in Sect. S1.
In Fig. 16 we provide the results of the sensitivity tests, showing bias and spatial correlation for both NO 2 and O 3 . The three chemistry based schemes give similar results indicating that in all three cases the calculations are close to photo-stationary.
The two empirical fits also give similar results, with the largest negative bias in NO 2 given by the Romberg scheme with -25%.

335
Since the Romberg scheme is specifically designed to reflect measurements, providing the correct NO 2 /NO X ratio, it can be regarded as the closest to the measurements. The bias differences between chemistry schemes and the Romberg scheme indicate that chemistry schemes have higher concentrations of NO 2 than the Romberg scheme, thus overestimate the NO 2 contribution when applied to annual mean concentrations. This is partially due to the positive bias in the EMEP model O 3 concentrations of 16%, but this only accounts for around 4% of the additional NO 2 . Included in Fig. 16 is the annual mean calculation without 340 the frequency distribution correction ("Travel time (annual)"), showing a 10% difference in bias when compared to calculations that use this correction. Spatial correlation is also improved by using the frequency distribution methodology.

Discussion
Downscaling only applies to emissions within a limited region of ± 0.1 o surrounding each receptor sub-grid. Based on the uEMEP calculation, the local contributions to NO X are significantly larger than for PM. The different source contributions at 345 measurement sites are given in Table 2 and this shows that, on average in Europe, 58% of the NO X contributions come from traffic within this limited region. In contrast only 19% of the PM 2.5 is attributable to residential heating, the largest downscaled contribution, from inside this region.
NO 2 is well modelled with high spatial correlation for many countries, but still with a significant negative bias of -18%.
There is significant variation in bias between countries even though the methodology is consistently applied to all countries.  it is generally the case that the spatial representativeness of the uEMEP calculations is suitable for comparison with these measurements (Lefebvre et al., 2013). Variation in bias between countries is then no longer a case of a mismatch in resolution but most likely reflects bias in the national emissions. uEMEP may be used to investigate this variability between countries further and to help harmonise future emission inventories across Europe.
There is a significant difference between the results achieved for the downscaling of PM compared to NO 2 . NO 2 is dominated by traffic emissions and this is spatially well defined using OSM as a proxy. The largest contributor to PM in the downscaled sources is residential heating, with contributions of 19% and 16% for PM 2.5 and PM 10 , respectively (Table 2). This is inline with other estimates of residential combustion in Europe. Thunis et al. (2017) calculated a contribution of 13% from residential 360 combustion from primary PM 2.5 averaged over 150 European cities, without downscaling. Population is used as a downscaling proxy for the residential source, but it appears that this is not a good proxy for high-resolution emission redistribution. Though clearly residential heating emissions occur where people live there can be large variation from city to city and from urban to suburban and to rural areas as heating practices vary significantly depending on housing type and on availability of alternative heating sources. To some extent this has been taken into account in the emission inventory at 0.1 o , but the emission proxy used 365 in uEMEP is likely not consistent with the EMEP emission inventory.
The Norwegian sensitivity tests show that when consistent emissions and emission proxies are used then spatial correlation can be significantly improved. For the application of uEMEP in Europe this was not the case since each country has their own methodology for calculating gridded EMEP emissions that may or may not make use of the downscaling proxies applied in uEMEP. A more consistent approach, as applied in Norway, would be to use the same spatial redistribution proxies in both 370 the gridded EMEP emissions and the downscaling proxies. This would require additional interaction and cooperation between emission inventory developers and air quality modellers.
It is worth noting that no selection of the Airbase monitoring data was carried out. All available stations with more than 75% coverage were used. This includes mountain stations, all traffic stations as well as industrial sited stations. In comparisons with the EMEP model these types of sites are often removed. All stations were also assumed to be sited at 3 m above the surface.

375
It is quite possible that different results would be obtained if a selection of stations was carried out. This will be assessed at a later time.

Conclusions
Downscaling of annual mean concentrations from the EMEP model have been carried out for NO 2 , PM 2.5 , PM 10 , and O 3 using the uEMEP model. Downscaling redistributes EMEP gridded emission data, using suitable proxy data, to high-resolution sub- The results for NO 2 show significant improvement with a doubling of spatial correlation for most countries and a significant reduction in negative bias. For NO 2 the downscaling works very well, which is due to the fact that NO X emissions are mainly 385 attributable to traffic and these emissions are well defined spatially with the proxy data used. O 3 concentrations are decreased due to higher NO X concentrations. Both concentrations and spatial correlations of O 3 are better simulated with uEMEP.
Neither PM 2.5 nor PM 10 shows any improvement in spatial correlation with the downscaling, though the negative bias in PM concentrations is improved. The spatial distribution of PM emissions can be improved, as demonstrated for Norway, with more accurate proxy data, but emissions of PM remain difficult to quantify properly at high resolutions and will require further 390 effort. One way forward is to harmonise the proxies used for both the EMEP gridded emissions and the uEMEP downscaling.
This has been shown to improve results in Norway.
Downscaling can provide additional information concerning the contributions of local sources. This may be combined with the EMEP model source-receptor calculations to provide a more complete picture of local and long-transported contributions.
The method can lead to a better assessment of local versus regional mitigation strategies to improve air quality in Europe at 395 high resolution. It also shows good potential to be used to improve exposure estimates.
Code and data availability.  Denby et al. (2020), together with the frequency distribution correction scheme (Section 2.4). In addition to this parameterisation two other chemistry based formulas are available. The first being the photo-stationary formulation, also 5 described in Denby et al. (2020), and the second a stationary formulation that allows for deviation from the photo-stationary state that may be caused by emissions, advection gradients or additional chemistry . Two empirical formulations are also included that are based on a fit to measurement data. The first, the Romberg scheme , is already described in Denby et al. (2020) and directly converts NO X to NO 2 concentrations. The parameters for this equation have been updated by fitting to all available Airbase data for the year 2017. The second empirical scheme, the SRM 10 scheme (Wesseling and van Velze, 2014), includes background O 3 and NO 2 combined with local NO X as input parameters.
The advantage of the two empirical fits is that they should convert NO X to NO 2 in a manner that is consistent with the observations, and as such can be applied to annual mean concentrations without correcting for non-linear chemistry. The first method is already described in (Denby et al., 2020), the remaining methods are presented in the following sections.

15
The solution to the photo-stationary equilibrium of NO 2 that only considers the NO X titration of O 3 (reaction rate k 1 ) and the photo-disassociation of NO 2 (photo-disassociation rate J) is already presented in Section 2.4. This steady state solution For the case where there is equilibrium but other terms such as emissions, deposition or additional chemistry are involved then we can define a constant λ such that For λ < 1 this implies a sink of NO 2 and for λ > 1 a source. The steady state solution can then be written, as in Eq.
(1), as By simply extracting the parameters listed in Eq. (S1) from the EMEP model we can determine λ and by using Eq. (S3) we can calculate the downscaled NO 2 concentrations. At the station sites used in this study we find on average λ = 0.95 with a standard deviation of 0.35. Around 60% of the sites have λ < 1. This means there is significant divergence from the photo-stationary 30 assumption in the EMEP annual mean calculations. λ is also strongly correlated with NO 2 , increasing with increasing NO 2 , indicating that emissions are perhaps the dominating factor leading to non photo-stationary equilibrium. It should be noted that is usually significantly larger than J/k 1 so even with a significant deviation of λ from unity then the NO 2 concentrations will not be largely affected.
This parameterisation is particularly useful if NO 2 concentrations from EMEP are to be recalculated, for example when 35 using the local fraction to assess the impact on concentrations with changes in emissions in a post-processing step. It is not necessarily the case that the conditions leading to a deviation of λ from unity are equally applicable for the downscaling chemistry.

S1.2 Updated Romberg scheme
The Romberg scheme  is the simplest method for converting NO X to NO 2 . It is based on an empirical 40 fit to the formula The parameters a, b (µg/m 3 ) and c are derived by fitting to observed annual mean concentrations. The term c should in 45 some way reflect the ratio of NO 2 /NO X emissions, i.e. the asymptotic limit for large NO X . Values for these parameters were given in Denby et al. (2020), based on a fit to Norwegian measurement data. For the European application a new fit is made to annual mean concentrations from all available European observations for 2017. This is shown in Fig. S1. The fitted parameters are: a = 41.1, b = 56.4 µg/m 3 and c = 0.162. The root mean square error of this fit is 2.9 µg/m 3 , or a normalised error of around 14%.

50
When implemented in uEMEP a small adjustment is made to Eq. (S4). It is desirable that background NO 2 levels provided by EMEP are not altered by this calculation. Background NO 2 and O 3 concentrations are calculated in uEMEP after the removal of the local fraction contribution from NO X . In Denby et al. (2020) this is referred to as the non-local contribution. To achieve this we rewrite Eq. (S4) as follows where ∆[N O 2 ] bg is the difference between the EMEP calculated NO 2 non-local contribution and the Romberg calculation of the non-local contribution given by This ensures that the calculated NO 2 concentrations are unchanged when the local contribution to NO X is negligible.

S1.3 SRM scheme
The SRM scheme for NO 2 calculation is implemented in The Netherlands as part of the Standard Calculation Method for air quality (Wesseling and van Velze, 2014). This scheme requires information on background levels of O 3 and NO 2 . Whilst this information is directly available from the model calculations this is more difficult to extract from measurement data. For this reason no fitting or changes to the parameters given in Wesseling and van Velze (2014) have been made. The equation is slightly re-written to be comparable with the Romberg scheme and the uEMEP formulation as: where 70 b srm = 100 1 − c srm and c srm = 0.15 The value for c srm applied here is the traffic exhaust NO 2 /NO X emission ratio used in uEMEP, and the value for b srm (µg/m 3 ) comes directly from Wesseling and van Velze (2014). The parameters c srm and b srm are refered to as F and K respectively in Wesseling and van Velze (2014).

75
This equation clearly resembles the Romberg scheme, Eq. (S4), with the exception that a background NO 2 and a local NO X are used. Here the background ozone would then be equivalent to the parameter a srm = [O 3 ] nonlocal .

S1.4 Additional information on the frequency distribution correction scheme
A frequency distribution correction scheme for annual mean calculations of NO 2 is described in Section 2.4. This is implemented to deal with the non-linear nature of the NO 2 chemistry when calculating annual means. Here we present some of the 80 background information used to derive this scheme.
An example of the NO X , O X and J frequency distributions is shown in Fig. S2. Here we see the log-normal distribution of the concentrations and the quite different distribution of the photo-dissociation rate J. J has a frequency of close to 0.5 at J = 0, since this, over the course of a year, is the amount of time the sun is below the horizon. The example shown here is taken from a low latitude European urban background station, with moderate NO X and high O 3 .

85
An assessment of measurement and model data for 2018 was carried out for Norwegian stations in order to derive the standard deviations of the NO X and O X concentrations. This assessment was limited to available data and a more substantial assessment could be carried out on a larger set of data. Focus is on the modelled distributions since this is what needs to be reproduced in order to correct the annual mean calculation. Firstly the concentration distributions were tested for a log-normal distribution. An example for both NO X and O X , both modelled and observed is shown for the urban background station  Figure S2. Example of the frequency distribution for NOX , OX and J concentrations used in the model. This is taken from the calculation of an actual site in Southern Europe. X axis is logarithmic.
Norway. Further assessment of the log-normal distribution was carried out for all modelled data at the 72 sites modelled in Norway for that year. This assessment showed that the concentrations were very close to log-normally distributed at all sites.  The standard deviations were then derived and compared to mean values. The results for both NO X and O X at modelled stations are shown in Fig. S5. For NO X the observed values are also included since these are more readily available. The 95 relationship for NO X is very robust and will reflect the temporal variation of the NO X sources, mostly traffic, and of the meteorology. It is worth noting that the normalised standard deviation of the traffic time profile applied in the Norwegian calculations is 0.71. Additional variability will be added due to meteorology. For O X the standard deviation is less dependent on the mean concentration. Even so, the normalised standard deviation of O X is significantly smaller than that of NO X .
Finally the correlation between NO X and O X was addressed since the frequency distribution correction assumes these two 100 concentrations are not correlated. This assessment showed both negative and positive correlations at the modelling sites. Sites with high NO X concentrations showed positive correlation and background sites showed mostly negative correlation. Since NO 2 is part of both NO X and O X then this positive correlation is not surprising. Values for the correlation (r) ranged from -0.4 to +0.6. This result shows there is correlation between NO X and O X but this would be difficult to take account of in the frequency distribution correction currently implemented in the model.  Figure S5. Scatter plot of hourly standard deviation versus annual mean NOX and OX for 2018 for all the 72 modelled sites. Also included are measurement data from 41 of the sites that measured NOX with data coverage >75%. The regression slope, intercept set to 0, and the correlation are also shown in the plots.

S7
S2 Scatter plots per country Figure S6. Scatter plots of annual mean NO2 concentrations per country for 2018 calculated with uEMEP. Only countries with 10 or more stations are shown individually but all stations are included in the final EU plot. Figure S7. Scatter plots of annual mean NO2 concentrations per country for 2018 calculated with EMEP. Only countries with 10 or more