Temporal and spatial analysis of ozone concentrations in Europe based on timescale decomposition and a multi-clustering approach

. Air quality measures that were implemented in Europe in the 1990s resulted in reductions of ozone precursor concentrations. In this study, the effect of these reductions on ozone is investigated by analyzing surface measurements of this pollutant for the time period between 2000 and 2015. Using a nonparametric timescale decomposition methodology, the long-term, seasonal and short-term variation in ozone observations were extracted. A clustering algorithm was applied to the different timescale variations, leading to a classiﬁcation of sites across Europe based on the temporal characteristics of ozone. The clustering based on the long-term variation resulted in a site-type classiﬁcation, while a regional classiﬁcation was obtained based on the seasonal and short-term variations. Long-term trends of deseasonalized mean and meteo-adjusted peak ozone concentrations were calculated across large parts of Europe for the time period 2000–2015. A

Abstract. Air quality measures that were implemented in Europe in the 1990s resulted in reductions of ozone precursor concentrations. In this study, the effect of these reductions on ozone is investigated by analyzing surface measurements of this pollutant for the time period between 2000 and 2015. Using a nonparametric timescale decomposition methodology, the long-term, seasonal and short-term variation in ozone observations were extracted. A clustering algorithm was applied to the different timescale variations, leading to a classification of sites across Europe based on the temporal characteristics of ozone. The clustering based on the long-term variation resulted in a site-type classification, while a regional classification was obtained based on the seasonal and short-term variations. Long-term trends of deseasonalized mean and meteo-adjusted peak ozone concentrations were calculated across large parts of Europe for the time period 2000-2015. A multidimensional scheme was used for a detailed trend analysis, based on the identified clusters, which reflect precursor emissions and meteorological influence either on the inter-annual or the short-term timescale. Decreasing mean ozone concentrations at rural sites and increasing or stabilizing at urban sites were observed. At the same time, downward trends for peak ozone concentrations were detected for all site types. In addition, a reduction of the amplitude in the seasonal cycle of ozone and a shift in the occurrence of the seasonal maximum towards earlier time of the year were observed. Finally, a reduced sensitivity of ozone to temperature was identified. It was concluded that long-term trends of mean and peak ozone concentrations are mostly controlled by precursor emissions changes, while seasonal cycle trends and changes in the sensitivity of ozone to temperature are among other factors driven by regional climatic conditions.

Introduction
Tropospheric ozone (O 3 ), together with particulate matter and nitrogen dioxide (NO 2 ), is one of the most troublesome air pollutants in Europe (EEA, 2016). A total of 17 000 premature deaths every year are attributed to excess O 3 exposure, without any sign of reduction in number of fatalities (EEA, 2016). In terms of impact on ecosystems, elevated concentrations of tropospheric O 3 are responsible for damaging agricultural production and forests mainly by reducing their growth rate. In addition, tropospheric O 3 acts as a greenhouse gas with an estimated globally averaged radiative forcing of 0.4 ± 0.2 W m −2 (IPCC, 2013). In the 1990s, emission control measures on O 3 precursors, namely nitrogen oxides (NO x = NO + NO 2 ) and volatile organic compounds (VOCs), were implemented in order to regulate air pollution. As a result, concentrations of NO x and VOCs have significantly declined in Europe (EEA, 2017;Colette et al., 2011;Guerreiro et al., 2014;Henschel et al., 2015). NO x emissions in particular declined in Europe by 48 % between 199048 % between and 201548 % between (EEA, 2017. Surprisingly, O 3 concentrations have not decreased as expected (Oltmans et al., 2013;Colette et al., 2018). Mean O 3 concentrations either remained stable or even increased in rural background areas from the 1990s and even until the mid-2000s in many European countries (e.g., Boleti et al., 2018;Munir et al., 2013;Paoletti et al., 2014;Querol et al., 2016;Anttila and Tuovinen, 2009). At urban sites an increase in mean O 3 has been observed; in some cases, an increase has been found at both rural and urban sites with larger upward trends observed at urban compared to the rural sites (Paoletti et al., 2014;Querol et al., 2016;Anttila and Tuovinen, 2009). However, a change in the trend has been observed after the mid-2000s, when mean O 3 concentrations have started to decline (Boleti et al., 2018;Munir et al., 2013). On the other hand, maximum O 3 concentrations decreased continuously from the 1990s until present (Paoletti et al., 2014), except for the traffic-loaded environments (Boleti et al., 2019). Downward trends of different metrics for peak O 3 have been found at many sites across Europe (Fleming et al., 2018). However, the high year-to-year variability of O 3 tends to mask the longterm changes, leading to a large fraction of sites with nonsignificant trends. Several studies based on either observations or climate models have shown that anthropogenic emissions can affect O 3 concentrations across continents (Dentener et al., 2010;Wild and Akimoto, 2001;Lin et al., 2017). The increase in background O 3 in Europe has been associated with increasing stratospheric O 3 contribution (Ordóñez et al., 2007), as well as increased hemispheric transport of O 3 and its precursors.
A shift in the seasonal cycle of O 3 has been observed in northern midlatitudes, i.e., the peak concentrations are now observed earlier in the year compared to previous decades, with a rate of 3-6 d per decade. (Parrish et al., 2013). This shift is attributed to increasing emissions of O 3 precursors in developing countries that led to an equatorward redistribution of precursors in previous decades (Zhang et al., 2016). Negative trends of the 95th percentile of O 3 and positive trends for the 5th percentile have been detected across Europe (Yan et al., 2018). Simultaneous decrease in maximum concentrations in summer and increase in winter indicate a decrease in amplitude in the seasonal variation in O 3 , probably as a result of the regulations in the 1990s (Simon et al., 2015).
O 3 variations are largely governed by climate and weather variability (Yan et al., 2018). Temperature especially influences O 3 concentrations in the troposphere, mainly by increasing the rates of several chemical reactions and increasing emissions of biogenic VOCs with increasing temperature (Sillman and Samson, 1995). Thermal decomposition of peroxyacyl nitrates (PANs) at high air temperature conditions results in elevated O 3 concentrations (Dawson et al., 2007). Indeed, extreme O 3 concentrations in central Europe are mainly associated with high temperatures (Otero et al., 2016). However, there are indications that the relationship of O 3 to temperature has changed in the last 20 years. For instance, in the US, O 3 climate penalty -defined as the slope of the O 3 versus temperature relationship -dropped from 3.2 ppbv • C −1 before 2002 to 2.2 ppbv • C −1 after 2002 as a result of NO x emission reductions (Bloomer et al., 2009). Additionally, Colette et al. (2015), based on chemistrytransport and climate-chemistry model projections, assessed the impact of climate change on the climate penalty and found that over European land surfaces summer O 3 change is [0.44; 0.64] and [0.99; 1.50] ppbv (95 % confidence interval) for the 2041-2070 and 2071-2100 time periods, respectively.
At different locations, O 3 may show a different temporal evolution due to a variety of factors, such as local pollution, topography, influence of nearby sources, or even transboundary transport of O 3 and its precursors. In addition, meteorological conditions can vary amongst different locations within large regions such as Europe, affecting O 3 concentrations in various ways. O 3 trend studies in the past have tried to tackle this issue, mainly by using clustering techniques to categorize European measurement sites based on different O 3 metrics (e.g., Henne et al., 2010). For instance, a site-type classification representing O 3 regimes between 2007 and 2010 was obtained by Lyapina et al. (2016) using mean seasonal and diurnal variations. In addition, a geographical categorization reflecting the synoptic meteorological influence on O 3 variation between 1998 and 2012 was obtained by Carro-Calvo et al. (2017). To tackle low spatial representation of urban and rural sites across large domains, i.e., midlatitude North America, western Europe and East Asia, Chang et al. (2017) obtained a latitude-dependent site classification with lower concentrations in western and northern Europe and higher concentrations in southern Europe. These studies indicate that the selected metric used to characterize O 3 in clustering leads to site classifications that represent different aspects of O 3 variability.
In the current study, a multidimensional clustering method that captures several influencing factors for the long-term trend of O 3 is presented. The temporal and spatial evolution of O 3 concentrations between 2000 and 2015 is studied using data provided by the European Environmental Agency (EEA). Mean O 3 concentrations are decomposed into the underlying frequencies based on a nonparametric timescale decomposition method to obtain the long-term (LT), seasonal (S) and short-term (W ) variations. The multidimensional clustering approach is applied to the distinct frequency signals LT(t) and S(t) extracted from the observations.
In addition, long-term trends of deseasonalized daily mean O 3 and meteo-adjusted peak O 3 concentrations are calculated. Through deseasonalization and meteo-adjustment, a significant fraction of the meteorologically driven variability of O 3 is excluded from the observations, and uncertainty in the trend estimation is reduced by a large factor. Intersections of site groups, i.e., LT(t) and S(t) clusters, are employed to guide the study of O 3 long-term trends. Furthermore, changes in the amplitude and phase of the seasonal variability of O 3 are explored based on the S(t) signal obtained by the timescale decomposition methodology. Finally, long-term changes in the relationship between O 3 and tem-perature are estimated and discussed for the different site environments and regions in Europe.

Data
Data for O 3 surface measurements are provided by the EEA (Air Quality e-Reporting) in an hourly resolution for the period between 2000 and 2015. In this study, only time series with a maximum of 15 % of missing values and a maximum of 120 consecutive days with missing values are used for the whole period of measurements, leaving the study with 291 sites across the European domain (Fig. 1) Meteorological variables are extracted from the ERA-Interim dataset on a 1 • grid at the location (longitudelatitude-altitude) of each station and in 3-hourly intervals. The variables considered for the meteo-adjustment of the peak O 3 metrics are temperature (K), specific humidity (g kg −1 ), surface pressure (hPa), boundary layer height (m), convective available potential energy (CAPE, J kg −1 ), eastwest surface stress (N s m 2 ) and north-south surface stress (N s m 2 ).
The present trend analysis focuses on (a) the deseasonalized daily mean and MDA8 O 3 and (b) the meteoadjusted MTDM and 4-MDA8 concentrations. The analysis of changes in the seasonal cycle of O 3 across Europe is based on the daily mean O 3 concentrations.

Timescale decomposition of daily mean and MDA8 O 3
Timescale decomposition refers to decomposition of the O 3 time series into the relevant underlying frequencies: where O 3 (t) is the daily mean and MDA8 O 3 time series, LT(t) the long-term variation, S(t) the seasonal variation, W (t) the short-term variation, and E(t) the remainder of the decomposition. Timescale decomposition in this study is performed with a nonparametric method, called the ensemble empirical mode decomposition (EEMD, Huang et al., 1998;Huang and Wu, 2008;Wu and Huang, 2009), which is considered a powerful method for decomposing O 3 time series (Boleti et al., 2018). The method detects the hidden frequencies in the time series based merely on the data and yields the so-called intrinsic mode functions (IMFs); each IMF represents one distinct frequency in the signal: where y(t) is the input data, c j is the different IMFs, n is the number of the IMFs and the remainder of the time series is the LT(t) of the input data. By adding together the IMFs with frequencies between around 40 d and 3 years we obtain the seasonal variation in O 3 (S(t) = c 7 (t) + . . . + c 10 (t)), and by adding the frequencies that are smaller than 40 d the short-term variation is acquired (W (t) = c 1 (t) + . . . + c 6 (t)). A more detailed discussion on the choice of the IMFs for the seasonal and short-term variations can be found in the study by Boleti et al. (2018).

Cluster analysis of O 3 variations
Cluster analysis is referred to pattern recognition in highdimensional data. The main idea is to represent n objects by identifying k groups based on levels of similarity. Objects in the same group must have the highest level of similarity, while objects from different groups must have a low level of similarity (Jain, 2010). The partitioning around medoids (PAM) clustering algorithm is used in this study. It is based on k means (MacQueen, 1967;Hartigan and Wong, 1979), which is a widely used clustering technique (e.g., Lyapina et al., 2016). PAM is more robust than k means and less sensitive to outliers because it uses medoids (actual points in the data) instead of centroids (usually artificial points) (R Development Core Team, 2017). It works as follows: first, a set of n high-dimensional objects (measurement sites) is clustered into a set of k clusters. Initially, k clusters are generated randomly, and the empirical means m k of the Euclidean distance between their data points are calculated. Then, each data point is assigned to its nearest cluster center (centroid). Centroids are iteratively updated by taking the medoid of all data points assigned to their clusters. The squared error (ε) between the m k and the points in the cluster (x i ) is calculated as follows: Each centroid defines one of the clusters and each data point is assigned to its nearest centroid, and the iterative process is terminated when the is minimized.
For identification of the clusters, the LT(t), S(t) and W (t) of the daily mean and MDA8 O 3 were used as input time series in the PAM algorithm. A sufficient number of clusters must be defined in order to capture dominant behaviors such that redundant information is avoided but at the same time not overlooking important characteristics. To identify the optimal number of clusters, the PAM algorithm is iteratively executed for a range of k values (number of clusters) and the average sum of (SSE) is calculated for each iteration, i.e., each k.
The number of clusters with the largest reduction in SSE is considered the most representative. Eventually, the choice of the ideal number of clusters results from a combination of the SSE approach and interpretability of the obtained clusters. In addition, a silhouette width (S w ) analysis is performed to assess the goodness of the clustering (Rousseeuw, 1987). More details about the number of clusters, the goodness of the clustering and the S w are provided in the Supplement (Sects. S1 and S2).

Daily mean and MDA8 O 3 long-term trend analysis
Meteorological adjustment is essential for calculation of robust O 3 long-term trends. Thus, daily mean and MDA8 O 3 observations are deseasonalized by subtracting the S(t) obtained with the EEMD from the observations (Boleti et al., 2018): where y d (t) is the deseasonalized time series and y(t) the observations. Through deseasonalization, observations are adjusted for the effect of meteorology on the inter-annual timescale. Theil-Sen trends (Theil, 1950;Sen, 1968) are then calculated based on monthly mean deseasonalized concentrations y d (t) for the period 2000-2015. The 95 % confidence interval of the trend is obtained by bootstrapping. The Theil-Sen trends were estimated using the openair library in R (R Development Core Team, 2017).

Peak O 3 concentrations long-term trend analysis
Trend analysis of peak O 3 metrics is performed for the MTDM and the 4-MDA8 O 3 , based on a meteo-adjustment approach following Boleti et al. (2019). A different approach for meteorological adjustment was used for the peak O 3 than for the daily mean and MDA8; deseasonalization is not meaningful for peak O 3 because peak O 3 events are temporally localized. Thus, daily maximum and MDA8 O 3 observations were linked to the available meteorological variables through generalized additive models (GAMs, Hastie and Tibshirani, 1990;Wood, 2006). The models are fitted for the warm season May-September, as by definition the MTDM refers to this period of the year. GAMs are instances of generalized linear models in which the model is specified as a sum of smooth functions of the covariates. A GAM can be described as follows: where O 3 (t) stands for the O 3 time series observations (daily maximum and MDA8), α is the intercept, s i is the smooth functions (thin plates splines) of the numeric meteorological variables M i and n denotes the number of the numeric meteorological variables in the GAM. The temporal trend is represented through the smooth function s 0 (t), where t is the time variable expressed by the Julian day. Finally, stands for the residuals of the model. For the GAMs, the following meteorological variables were used: the daily maximum tempera- ture, daily mean specific humidity, daily mean surface pressure, daily maximum boundary layer height, morning mean CAPE, daily mean east-west surface stress and daily mean north-south surface stress, as well as a time variable, i.e., the Julian day. The above explanatory meteorological variables are the ones that were most often selected at the Swiss sites by the meteorological variable selection performed by Boleti et al. (2019). The GAMs were estimated with the mgcv library in R (R Development Core Team, 2017). The meteo-adjusted daily maximum and MDA8 O 3 concentrations were calculated similar to the approach applied by Barmpadimos et al. (2011): where α is the intercept of the model, s 0 (t) the time variable as Julian day and (t) the residuals. The meteo-adjusted MTDM and 4-MDA8 concentrations were estimated based on the meteo-adjusted values (O 3 adj (t)) on the same days as they were identified before the meteo-adjustment. Eventually, meteo-adjusted trends were calculated with the Theil-Sen trend estimator applied on the O 3 metadj (t).

O 3 seasonal cycle trend analysis
The S(t) signal extracted with the EEMD captures the meteorologically driven O 3 variation on yearly to multiyear timescales and is more representative compared to parametric fitting approaches (Boleti et al., 2018). Here, changes in the daily mean S(t) of O 3 throughout the studied period are identified as follows: the maximum and minimum O 3 value, as well as the day when the maximum O 3 occurred in each year, are identified in the S(t), referred to here as S max , S min and S DoM , respectively (Fig. 2). A Theil-Sen trend estimator for each of the S max (t), S min (t) and S DoM (t) is applied for each site cluster, representing the long-term temporal evolution of the amplitude and phase of S(t).

Relationship between O 3 and temperature
The relationship between O 3 and temperature is studied for the warm season between May and September. A linear regression model between daily maximum O 3 concentrations and daily maximum temperature is applied for each year throughout the studied period 2000-2015 as follows: where O 3 (t) is the time series of the daily maximum O 3 , T (t) the time series of the daily maximum temperature and n is the number of years. β 0i is the intercept and β 1i (t) the parameter describing the linear effect of temperature on O 3 . Then, a linear model is applied to β 1i (t) over all years for each site cluster to identify the long-term trend of the slope between O 3 and temperature maximum values. In addition, a linear regression model is applied on the daily maximum O 3 concentrations against binned temperature ranges and in three consequent time periods (2000-2005, 2005-2010 and 2010-2015).

Cluster analysis
Here, we present the results of the daily mean LT(t) and S(t) clustering; results for the W (t) clustering and the cluster analysis based on the MDA8 are provided in the Supplement (Sect. S3). The clusters based on daily mean and MDA8 O 3 are rather similar; thus, the choice of clusters based on these two metrics does not affect the conclusions. The daily mean O 3 clusters depict the main influencing factors for O 3 trends, i.e., proximity to emission sources and meteorological conditions. Therefore, it is appropriate to study long-term trends based on the O 3 daily mean clusters. Application of the clustering algorithm to the LT(t) leads to a site-type classification, which largely reflects the proximity to emission sources of O 3 precursors. S(t) and W (t) clustering leads to a regional site classification, which reflects the importance of the climate on the annual cycle of O 3 . It is observed that a few sites have a negative S w , which means that these sites are assigned to a certain cluster, although they do not really fit into any of the identified clusters (see Supplement, Sects. S2 an S5). Nevertheless, the sites with negative S w were not excluded from the data analysis as they do not have a noticeable influence on the results. O 3 concentrations often increase with distance from emission sources of NO x . Thus, LT(t) clustering leads to identification of site groups with similar types of environment in terms of proximity to precursor emissions and mean O 3 concentrations, which are indicative of multiannual changes in the O 3 time series (Boleti et al., 2018). This measurementbased classification can be more informative than reported station types, for example, there are rural sites with nearby pollution sources or even sites with surroundings that might have changed dramatically with time. In the following section, clusters obtained from analysis of daily mean O 3 are presented and discussed in detail; the clusters derived from MDA8 O 3 are similar and presented in the Supplement (Sect. S3).
Cluster analysis of the long-term variation LT(t) resulted in four clusters that mainly differ in the daily mean O 3 concentration levels: cluster 1 includes sites that are marked in the Air Quality e-Reporting data repository as being of rural site type and sites that are mostly located at relatively high altitudes (on average 800 to 1200 m). The sites in cluster 1 show the highest O 3 concentrations, as illustrated in Fig. 3. The high mean O 3 concentrations indicate that the sites in cluster 1 represent background situations with minor influence from nearby emissions of man-made O 3 precursors. This cluster is therefore denoted as background cluster (BAC). Cluster 2 includes mostly rural sites that are located at lower altitudes of around 300-600 m and is therefore labeled as rural cluster (RUR, Fig. 3). The sites in cluster 3 are also located at low altitudes (around 100 to 300 m) and represent rural, suburban and urban site types in similar numbers. The sites in cluster 3 seem to be influenced by nearby man-made emissions of O 3 precursors such as NO x , leading to lower mean O 3 concentrations compared to the sites in the RUR cluster. Cluster 3 consists of moderately polluted sites and is denoted as cluster MOD. Finally, cluster 4 consists mostly of urban and suburban sites showing the lowest daily mean O 3 concentrations, likely due to the proximity to sources of NO x emissions and the resulting enhanced depletion of O 3 through reaction with NO. Consequently, cluster 4 is denoted as the highly polluted cluster (HIG).
The LT(t) signal, as derived from the daily mean (Fig. 3) and MDA8 O 3 (Fig. S9 in the Supplement) observations, increases for BAC and RUR until around the beginning of the 2000s and decreases afterwards. For the MOD and HIG clusters the same pattern was observed, but the decrease starts much later than at the rural sites, i.e., around the end of the 2000s. Especially at the HIG sites a level-off is mostly observed after 2010 rather than a decrease. Similar temporal evolution with inflection points in the LT(t) has been observed in the study by Boleti et al. (2018) which was focused on trends of average O 3 concentrations in Switzerland.
Clusters derived from the daily mean S(t) show a regional representation most likely due to the influence of the regional climate conditions and the annual cycle of O 3 . The following five clusters were obtained from the daily mean S(t) (Fig. 4): 1. "CentralNorth" comprised of the northern and eastern parts of Germany; the Netherlands; and some eastern sites in Czech Republic, Poland, and Austria; 2. "CentralSouth" covers most of Austria, Switzerland, and some sites in southwestern Germany;  (Stohl, 2002;Dentener et al., 2010). Thus, the influence of background O 3 , i.e., O 3 inflow from North America and East Asia, is high at these sites (Derwent et al., 2004. In all clusters (   by spatial interpolation, leading to a larger and regular geographical coverage compared to the available observations. Compared to Carro-Calvo et al. (2017), similar geographical clusters were identified here, except for the Iberian Peninsula, eastern Europe, northern Scandinavia and the Balkan states that do not appear as separate clusters in our analysis. This is most probably due to the small number of observational sites in the above regions. In contrast to our study, gridded MDA8 O 3 concentrations during summer have exclusively been used for the cluster analysis by Carro-Calvo et al. (2017), and therefore in conditions wherein the correlation of O 3 and meteorological variables such as temperature is typically strongest. In addition, the present study results in spatial classification by utilizing the seasonal variation, while Carro-Calvo et al. (2017) have used normalized anomalies. Four site-type clusters were found based on the LT(t) in this study similar to Lyapina et al. (2016) based on absolute mixing ratios of O 3 variations, which identified five site-type clusters ranging from urban traffic (as in the HIG cluster here) to rural background environments (equivalent to RUR). Similar site classifications are obtained because the LT(t) signal of this study and the mean seasonal and diurnal profiles of Lyapina et al. (2016) both capture the O 3 concentration levels, distinguishing specific pollution regimes.

Trends of daily mean O 3 concentrations
The daily mean LT(t) and S(t) clusters identified in Sect. 4.1 are used for assessment of the temporal trends for the different site types and geographical locations. Overall, decreasing daily O 3 means are found for rural sites, while there is a tendency for increasing O 3 in more polluted urban environments (Fig. 5). MDA8 trends are similar to the daily mean trends and are shown in Sect. S4. At 64 % of all sites significant trends (p value < 0.05) were identified for the daily mean O 3 ; 61 % of the significant trends were negative and Figure 5. Box plots of deseasonalized daily mean O 3 trends for the LT(t) and S(t) clusters. LT(t) clusters represent a site-type classification, while the S(t) clusters represent a geographical classification that is influenced by regional climate conditions. Boxes include 25th to 75th percentiles with the line indicating the median value; whiskers extend to 1.5 times the interquartile range. 39 % positive. Most rural sites -BAC and RUR -experienced decreasing daily mean O 3 concentrations in all regions, as expected following the NO x and VOC reductions in Europe (Fig. 5). At the MOD and HIG sites a leveling off or increase is observed, respectively, especially in the Central-North, West and North clusters. At the HIG sites the positive trends could partly be explained by the increase in NO 2 to NO x ratio, originating from the proliferation of diesel vehicles, which have increased in the European car fleet (EEA, 2009). In addition, the strong reduction of NO x concentrations that led to less titration of O 3 by NO could also explain the positive trends at urban and suburban sites. The late inflection point at urban sites (LT(t) of HIG cluster in Fig. 3) could be an additional effect of the reduced titration of O 3 , which leads to positive trends at the HIG sites. Flat trends at central European sites might partially be explained by the increasing influence of North American and Asian emissions, which have counterbalanced the decrease in European NO x and VOC concentrations (Derwent et al., 2018;Yan et al., 2018).
In agreement with our results, significant decreases in daytime average and summertime mean of MDA8 O 3 at European rural sites, as well as small and nonsignificant downward trends of MDA8 at urban sites, have been found previously for the time period 2000-2014 (Chang et al., 2017). Similarly, in a report by EEA (2016) it was found that between 2000 and 2014 annual mean O 3 and annual mean MDA8 O 3 have been decreasing at rural background sites, while at more polluted sites influenced by nearby man-made precursor emissions upward trends have been detected.
O 3 trends at sites in the North cluster indicate changes in background O 3 , especially in the RUR clusters that are mostly free from local emissions. Here, decreasing trends of daily mean O 3 were found at RUR sites, while in the MOD and HIG the trends are slightly increasing. It is interesting to compare the trends in the North cluster with the temporal evolution of O 3 at Mace Head (a remote station in northwestern Ireland), which is representative for inflow of background O 3 into Europe. For this reason, we estimated the LT(t) variation in MDA8 O 3 and the Theil-Sen trend for the site at Mace Head (Fig. 6) to compare with the MDA8 O 3 trend identified by Derwent et al. (2018). Here, an inflection point was

Trends of peak O 3 concentrations
Peak O 3 concentrations in summertime are attributed to photochemical production during this time of the year, and the spring maximum in remote locations is linked to stratospheric influx and hemisphere-wide photochemical production during that season (Holton et al., 1995;Monks, 2000). In this study, significant negative meteo-adjusted MTDM trends were observed at 62 % of the sites (ratio of sites with negative trends and number of all sites), while without meteoadjustment significant negative trends were identified at only 19 % of the sites (4-MDA8 trends are shown in the Sect. S4). The higher number of sites with significant trends after the meteo-adjustment indicates the advantage of using meteoadjusted observations in the trend estimation. This argument is supported in the study by Fleming et al. (2018), where significant negative trends of the fourth-highest MDA8 O 3 between 2000 and 2014 have been detected at only 18 % of the studied sites across Europe, while at a large proportion of sites either weak negative to weak positive or no trends at all were found. The nonsignificant trends have been attributed by Fleming et al. (2018) to the influence of meteorology, which is not considered in their trend estimation.
Trends of meteo-adjusted MTDM are discussed here for the daily mean O 3 LT(t) and S(t) clusters. MTDM decreased for all site types and regions during the studied period 2000-2015. However, in the RUR cluster MTDM showed the strongest decrease of all LT(t) clusters (Fig. 7). Interestingly, in the BAC cluster (especially the West cluster) the decrease in MTDM was not so pronounced, likely due to the increase in hemispheric transport of O 3 in Europe (Derwent et al., 2007;Vingarzan, 2004). The same pattern was observed at  . Box plots of meteo-adjusted MTDM trends for the daily mean O 3 LT(t) and S(t) clusters. LT(t) clusters represent a sitetype classification, while the S(t) clusters a represent a geographical classification that is influenced by regional climate conditions. the high alpine site of Jungfraujoch, which is representative for European continental background O 3 concentrations (Balzani-Lööv et al., 2008;Boleti et al., 2019). Also, HIG sites in the CentralSouth cluster showed a slightly smaller decrease in MTDM compared to the other regions, possibly due to industrialization in the neighboring areas of eastern Europe (Vestreng et al., 2009). Nevertheless, in order to estimate the reasons and quantify the exact influence of the above factors on the trends, dedicated modeling studies are needed.
Our results are in line with a modeling sensitivity study, which found negative trends of the 95th percentile of O 3 concentrations at European rural background sites for the period 1995-2014 (Yan et al., 2018). For the shorter time period between 1995 and 2005, downward trends of measured MTDM have been observed in most parts of Europe as well (in the range [−0.12, −0.55] ppb yr −1 ), with the highest decrease in the Czech Republic, UK, and the Netherlands (on average −1 ppb yr −1 ) and a very small (nearly flat) decrease in Switzerland (EEA, 2009). In our study, measured MTDM trends (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015) for these regions are in the same range, i.e., the average decrease was estimated between −0.28 and −0.55 ppb yr −1 . Smaller trends in Switzerland and central and northeastern Germany have been observed by EEA (2009), which agrees with our result for the Central-South cluster sites that showed the smallest average trend of all clusters. The flat trend in Switzerland during 1995-2005 is probably linked to the disproportional decrease in NO x and VOCs until the beginning of the 2000s, when an inflection point was observed at most polluted sites (Boleti et al., 2018). In Germany, a mixed behavior was observed by the EEA (2009), with the northeastern part showing a stronger decrease and the central and northeastern region showing a smaller decrease. Similar to our differentiation between the clusters, average MTDM trends within the CentralNorth cluster (with northern and northeastern Germany included) were higher, and within the CentralSouth cluster (covering central parts of Germany) these trends were lower.

O 3 seasonal cycle trends
Analysis of S(t) allows for studying the characteristics of the annual cycle of O 3 without influence of short-term meteorological phenomena and long-term variations. Here, the trends of S max (t), S min (t) and S DoM (t) are presented for the five regions identified based on S(t), as calculated from daily mean O 3 . In Fig. 8 the average annual O 3 cycles are shown; it is clear that in the PoValley cluster the day of maximum O 3 occurs in summer (June-July), while in the North cluster it occurs around spring (late March-April). This agrees with previous studies that have reported that both the highest average O 3 concentrations and extreme O 3 episodes tend to occur over the central and southern parts of Europe during summer, while over northern Europe they occur during spring (Schnell et al., 2015;Ordóñez et al., 2017). A declining amplitude of S(t) and a simultaneous phase shift towards an earlier time in the year can be observed for the 2000 to 2015 period (Table 2).
More specifically, an overall decrease in O 3 S max (t) by around 0.05-0.18 ppb yr −1 and a simultaneous increase in S min (t) with a rate of around 0.25 ppb yr −1 was observed for all S(t) clusters (Fig. 9). However, in the North cluster the decrease in the S max (t) was very small and nonsignificant, probably due to the pronounced influence of background O 3 at these sites. The most pronounced shortening of S(t) am- plitude can be seen at the PoValley sites, where the downward trend of peak O 3 is largest (Fig. 7). The increase in the S min (t) may be partially due to the decreased titration of O 3 after reductions of NO x emissions and is probably due to the increased influx of O 3 towards northern and northwestern Europe and more cyclonic activity in the North Atlantic during winter as well (Pausata et al., 2012). A decrease in the 95th percentile and an increase in the 5th percentile of O 3 for the period 1995-2014 has been also identified in the EMEP network (rural background sites) (Yan et al., 2018). Lower summertime peaks as a result of decreased photochemical production and higher O 3 concentrations during the winter months due to decreased O 3 titration have been found in European air masses between 1987 and 2012  as well.
The trend of O 3 S DoM (t) is negative for all regions, i.e., the occurrence of the day of maximum O 3 has shifted to earlier days within the year with a rate of −0.47 to −1.35 d yr −1 (Table 2, Fig. 10). The observed shift of the day of seasonal maximum might be linked to the increase in emissions in East Asia. The associated strong photochemical reaction rates, convection, and NO x sensitivity in the tropics and subtropics (Derwent et al., 2008;West et al., 2009;Fry et al., 2012;Gupta and Cicerone, 1998) have probably contributed to increased transport of air pollution to midlatitudes and northern latitudes (Zhang et al., 2016). In addition, changes in meteorological factors may have affected the seasonal variation in O 3 . For instance, a similar behavior with an earlier onset of the summer date (the calendar day on which the daily circulation-temperature relationship switches sign) has been observed by Cassou and Cattiaux (2016) using observational data, while Peña-Ortiz et al. (2015) have found that the summer period is extending with a rate of around 2.4 d per decade based on gridded temperature data in Europe. The positive phase of the NAO leads to increased O 3 concentrations in Europe through enhanced transport of O 3 and precursors across the North Atlantic from North America to Europe Creilson et al. (2003). An increase in baseline O 3 Table 2. Linear trends of S max , S min and S DoM for the daily mean O 3 S(t) clusters during 2000-2015: * * indicates a highly significant trend (p value < 0.01), and -indicates a nonsignificant trend (p value > 0.05). related to the prevailing positive NAO index -and the associated westerly flow and intercontinental transport -during the 1990s and the beginning of the 2000s is probably a factor contributing to the increase in the winter S min O 3 values (Pausata et al., 2012). The enhanced hemispheric transport of air pollutants from North America to Europe is related to increased transport through frontal systems as well (Creilson et al., 2003;Eckhardt et al., 2003). Increased O 3 in winter and spring (but not in summer) might lead to a shifting from a pronounced maximum in late summer to a broader springsummer peak (Cooper et al., 2014). At the West cluster sites, a slightly stronger shift of the S DoMax was observed compared to other clusters, while at the North cluster sites the decrease was the smallest. The early spring maximum at the North cluster sites in April (Fig. 8) can be explained by elevated NO x that is released from PAN and alkyl nitrates that are produced during winter at northern latitudes (Brice et al., 1984;Bloomer et al., 2010).

O 3 and temperature relationship
The O 3 sensitivity to temperature is a useful metric for validation of precursor reduction scenarios and emission inventories in chemistry-transport models (Oikonomakis et al., 2018). Here, we present the long-term trends of the relationship between the daily maximum O 3 concentrations and daily maximum temperature during the warm season from May to September between 2000 and 2015. Daily maximum O 3 and temperature are chosen in order to represent peak O 3 concentrations formed during the considered days. Decreasing sensitivity of O 3 with respect to temperature was observed during the considered time period in all regions (Fig. 11). Figure 11 shows the decreasing slopes of linear regression lines of maximum O 3 against temperature for successive year groups. The decrease is consistent for all calculated regional clusters except for the North cluster. For most regions in Europe a significant downward trend of the slope of around 0.04-0.05 ppb K −1 yr −1 was found (Table 3). At PoValley sites the decrease was more pronounced (−0.083 ppb K −1 yr −1 ) because of large reductions of precursor concentrations in this region, which is characterized by high industrial emissions. Note that the average correlation between O 3 and temperature in that cluster is the highest compared to the other regions. At the North cluster sites the weakest correlation of O 3 to temperature was observed and the trend is nonsignificant. This is expected because at these high latitudes the mean temperature is lower compared to other regions in Europe (Fig. 11); thus, photochemical production of O 3 is weak during the time when O 3 typically reaches its maximum concentration. In addition, in these northern regions the variability of temperature is lower compared to the central and southern parts of Europe, while O 3 concentrations are more influenced by intercontinental transport mechanisms.
In relation to the LT(t) clusters, it was observed that the higher the pollution burden of the site is, the stronger the trend of O 3 to temperature slope will be ( Table 4). As shown here, the HIG and MOD sites have higher trends compared to the clusters BAC and RUR. Our results are in line with a box model study that tested the O 3 -temperature relationship under different NO x level scenarios (Coates et al., 2016). Coates et al. (2016) have shown that at high NO x conditions O 3 increases more strongly with temperature, while the increase is less pronounced when moving to lower NO x conditions. Consequently, regional O 3 production has mainly decreased at the most polluted locations, due to considerable reductions of precursor emissions.

Conclusions
In this study, a classification of 291 sites across Europe was performed for the time period 2000-2015. The clustering algorithm applied on the long-term changes LT(t) and the seasonal cycle S(t) of daily mean O 3 resulted in a site-type and geographical site classification, respectively. Such a twodimensional site classification scheme provides an innovative approach for O 3 trends studies in large spatial domains and can be of significant use in model evaluation studies (e.g., Otero et al., 2018). Our approach captures several features of O 3 variations, i.e., pollution level from the LT(t) clustering and influence of the regional climate conditions from the S(t) clustering, and presents a unifying perspective on past studies that report different site-type labels based on cluster analysis using different metrics of O 3 concentrations. The two-dimensional approach offers a tool beyond traditional clustering studies to study the factors that affect O 3 trends. In addition, O 3 time series analysis is complex and bene- Table 3. Linear trends of the O 3 -temperature slope (based on daily maximum values) for the daily mean O 3 S(t) clusters for the time period 2000-2015: * * indicates a highly significant trend (p value < 0.01), and -indicates a nonsignificant trend (p value > 0.05).  fits from being characterized in different ways, as well as grouped based on several features. The regional differentiation is hampered by sparse or missing measurement sites in some regions, e.g., eastern Europe or the Balkan peninsula. However, in the last few years the number and spatial distribution of sites with longer and more dense measurements has improved. A trend analysis of deseasonalized mean O 3 concentrations and meteo-adjusted peak O 3 concentrations was implemented for the considered sites. By using LT(t) and S(t) clusters, patterns of O 3 long-term trends across Europe were investigated, based on the multidimensional site classification scheme. Long-term trends of deseasonalized daily mean O 3 are decreasing at rural sites, while at suburban and urban sites they are either stable or slightly increasing. Positive or flat trends indicate that reduction of precursors has been less effective in reducing O 3 concentrations in heavily polluted environments. On the other hand, downward trends in peak O 3 concentrations were observed in all regions as a result of precursor emissions reductions. However, peak O 3 has been decreasing with the smallest rate at higher-altitude sites, especially in the western Europe, possibly due to the influence of background O 3 imported from North America and East Asia.

Daily mean O 3 S(t) cluster Trend (ppb K
The analysis of S(t) extrema revealed a decrease in summertime maxima and an increase in wintertime minima, pointing to a decreasing amplitude of the seasonal cycle of O 3 . At the same time, the occurrence of the day of maximum has shifted from summer to spring months with a rate of around −0.5 to −1.3 d yr −1 . Changes in the S(t) might be attributed to precursor reductions in Europe, changing weather patterns in the North Atlantic and increases in emissions in southern East Asia.
Finally, the sensitivity of O 3 to temperature has weakened since 2000 with a rate of around 0.04 ppb K −1 yr −1 , i.e., formation of O 3 became weaker at high temperature conditions, which can be attributed to the decrease in NO x concentrations. The trend of the sensitivity differs across sites that are influenced by different meteorological conditions. Data availability. The data used for this study can be made available upon request.
Author contributions. EB carried out the analysis. EB conceived the idea and methodology and wrote the manuscript with CH and ST. SKG contributed to data preparation and analysis. ASHP consulted regarding the methodological approach and interpretation of results.