Inﬂuence of meteorology on PM 10 trends and variability in Switzerland from 1991 to 2008

. Measurements of airborne particles with aerodynamic diameter of 10 µm or less (PM 10 ) and meteorological observations are available from 13 stations distributed throughout Switzerland and representing different site types. The effect of all available meteorological variables on PM 10 concentrations was estimated using Generalized Additive Models. Data from each season were treated separately. The most important variables affecting PM 10 concentrations in winter, autumn and spring were wind gust, the precipitation rate of the previous day, the precipitation rate of the current day and the boundary layer depth. In summer, the most important variables were wind gust, Julian day and afternoon temperature. In addition, temperature was important in winter. A “weekend effect” was identiﬁed due to the selection of variable “day of the week” for some stations. Thursday contributes to an increase of 13% whereas Sunday contributes to a reduction of 12% of PM 10 concentrations compared to Monday on average over 9 stations for the yearly data. The estimated effects of meteorological variables were removed from the measured PM 10 values to obtain the PM 10 variability and trends due to other factors and processes, mainly PM 10 emissions and formation of secondary PM 10 due to trace gas emissions. After applying this process, the PM 10 variability was much lower, especially in winter where the ratio of adjusted over measured


Introduction
Airborne particles with aerodynamic diameter of 10 µm or less (PM 10 ) have important negative consequences on human health (Zemp et al., 1999;Bayer-Oglesby et al., 2005). International organizations and governments worldwide have tried to reduce emissions and human exposure to PM 10 by issuing guidelines and directives and by taking a number of legislative measures. In Switzerland yearly average PM 10 concentrations should not exceed 20 µg m −3 and a daily average value of 50 µg m −3 should not be exceeded more than once a year (Swiss Federal Council, 2010a). As of 2008 the European Union has introduced legislation for PM 2.5 as well (European Parliament and the Council, 2008). As a result of the implementation of air quality legislation, air pollution levels (including particulate matter) in most developed countries have declined after 1970(Darlington et al., 1997. Further implementation and assessment of mitigation policies require high quality long-term datasets providing quantitative air quality information at a regional level. When such datasets exist, a detailed investigation of PM 10 trends is of great value both on a scientific and on a regulatory level. An important aspect to be considered when studying trends of air pollution is that the trends do not only depend on the emissions of certain pollutants but also on the Published by Copernicus Publications on behalf of the European Geosciences Union. meteorological conditions in the area of interest (Gomiscek et al., 2004). Hooyberghs et al. (2005) reported in a study for Belgian urban centers that day-to-day PM 10 concentrations depend more on meteorology than on anthropogenic emissions. Pryor and Barthelmie (1996) emphasized that, one should first isolate and remove the effect of weather conditions before investigating the influence of anthropogenic activities on air pollution. To do so, knowledge of the dependence of the concentrations of certain pollutants on various meteorological variables is required. Only a few studies exist on this wide research field. For example, Easter and Peters (1994) demonstrated the non-linear effect of temperature and relative humidity changes on sulphate aerosol. However, temperature and relative humidity represent only a few of the meteorological variables that affect pollutant concentrations.
Although this study does not focus on PM 10 chemical composition, it is important to acknowledge that PM 10 is a complex mixture of chemical compounds whose behaviour depends strongly on the atmospheric conditions. In a paper that summarises data from a number of European sites, Putaud et al. (2004) report that for the four Swiss sites included in their study (Bern, Zurich, Basel and Chaumont), the most important compounds in terms of their contribution to total PM 10 mass are black carbon, organic matter, mineral dust, ammonium, nitrate and sulphate, while a considerable unaccounted fraction is present as well. The relative contribution of different species to the total PM 10 mass varies between different site types. The kerbside site at Bern has a relatively large fraction of black carbon and a low fraction of sulphate and ammonium. Mineral dust fraction was large at the kerbside station and at the rural station at Chaumont whereas it was low at the urban background stations at Zurich and Basel. Ammonia and nitrate have the largest fraction at urban background sites and a low fraction at the rural site. Further important consituents of PM 10 in Switzerland are organic matter and the unaccounted fraction which do not vary significantly between different site types and seasons. Information on the chemical composition of organic aerosol for different sites is given by Lanz et al. (2010), who show that the oxygenated secondary fraction is the dominating organic fraction. For PM 1 in general, larger differences are found between the seasons compared to station types. A further chemical characterisation of PM 10 at Swiss sites including water soluble ions and trace elements is found in .
Several methodologies for modeling of PM 10 trends and for prediction of PM 10 concentrations are reported in the literature. Zolghadri and Cazaurang (2006) developed an adaptive nonlinear model for the prediction of PM 10 in the Bordeaux area. However, the methodology for the development of the model has strong empirical elements, which make it strictly valid only for the area that it was developed for. Another approach for statistical modeling of PM 10 involves neural networks (Hooyberghs et al., 2005), which demonstrated a fairly good predictive skill of PM 10 for Belgium. Ordóñez  Table 1. Sidebar indicates height above mean seal level. et al. (2005) used an extensive meteorological dataset to carry out a statistical analysis of the effect of weather on ozone trends in Switzerland (see also Keller et al., 2008). In this study however, it was assumed that ozone has a linear dependence on each meteorological variable. In a more recent analysis of the effect of meteorology on PM 10 trends in the Netherlands (Velders and Matthijsen, 2009) it was assumed that the long term trend was only a result of emissions and that meteorology only affects the variability of PM 10 concentrations after removing the trends. However, it is shown later that changes in meteorology can also be responsible for long term changes in PM 10 concentrations.
To carry out an adequate investigation of PM 10 trends, one needs a good spatial and temporal coverage of the region in question (Pryor and Barthelmie, 1996). PM 10 measurements within the Swiss National Air Pollution Monitoring Network (NABEL) are available from 13 stations since 1991, providing a rare opportunity to study the PM 10 trends during a time period with marked changes in air pollution levels in developed countries. The stations represent various site types and altitudes (see Table 1). Details on the dataset used in this study are given in Sect. 2. The effect of a range of available meteorological variables on PM 10 was estimated by Generalized Additive Models (GAMs). PM 10 trends in each season were investigated separately. A description of the methodology is found in Sect. 3. The results and the conclusions are presented in Sects. 4 and 5 respectively.

Data
PM 10 concentrations used in this study were measured at 13 sites of the Swiss National Air Pollution Monitoring Network NABEL (Empa, 2010). The names and details of the considered sites are listed in Table 1 and their locations are shown in Fig. 1. The PM 10 data were obtained using high-volume samplers (DHA-80, Digitel) and gravimetric determination of the collected PM mass. From 1991 until the end of 2005 glass fiber filters (Ederol 227/1/60) were used. Since then PM 10 is collected on quartz fiber filters (Whatman QMA). Until 1997 the high volume samplers were equipped with inlets for total suspended particles (TSP). The TSP values from 1991 to 1997 were converted into PM 10 by applying a factor that was derived separately for all site types during one year parallel PM 10 measurements. For 1997, this factor is 0.63 for the urban kerbside site in Bern and between 0.81 and 0.90 for the other stations (Empa, 2010). The correction process for the TSP data ensures continuity of the TSP and PM 10 timeseries. The main assumption involved in this correction is that the PM 10 to TSP ratio at each station was approximately constant throught the 1991-1997 period and equal to the 1997 value. This is considered to be reasonable. Scatterplots of TSP measurements before and after scaling versus parallel PM 10 measurements for year 1997 and for 8 stations are provided as supplementary material.
In 2001, the frequency of PM 10 sampling for the gravimetric method was reduced from daily sampling to PM 10 collection every fourth day. Since then, the gravimetric PM 10 measurements are complemented by beta attenuation monitors (FH62-IR, Thermo ESM Andersen) and later on by TEOM FDMS (Thermo Environmental) instruments. However, the PM 10 data obtained with the PM 10 monitors were corrected to be consistent with gravimetrically determined PM 10 for the days where gravimetrically determined PM 10 were not available. The correction takes into account the previous and next gravimetric PM 10 observation using the formula put forward by Gehrig et al. (2005). The correction of the beta monitor and TEOM data introduces some random noise but no bias and therefore no breaks or inconsistencies in the PM 10 timeseries are caused by the change of the sampling period from 1 to 4 days. The observational errors of the aforementioned measurement methods at the 95% level of confidence are 1.2 µg m −3 for the gravimetric method, 1.4 µg m −3 for the beta monitor and 2.2 µg m −3 for the TEOM instrument.
At some of the sites meteorological measurements are available from the NABEL network. The remaining sites are collocated with meteorological measurements stations operated by the Swiss weather service (MeteoSwiss) and meteorological data were provided from MeteoSwiss. The meteorological data measured include wind speed, wind direction, precipitation, global radiation, net radiation, air temperature and relative humidity. NABEL and MeteoSwiss use the same instruments for their meteorological measurements.
Wind speed and wind direction are measured using a rotating cup anemometer and a weather vane located typically at 10 m above ground. Precipitation is measured by means of a tipping bucket rain gauge with a built-in heating system for measuring snow and ice. Further instrumentation installed at the sites includes a pyranometer for measuring global radiation and a pyranometer-pyrgeometer coupled system for the measurement of net radiation. Temperature is measured using a thermocouple and dewpoint is measured by a chilled mirror hygrometer. A detailed description of the instrumental setup for meteorology measurements is given in Empa (2010).
The resulting air pollution and meteorology datasets go through an automatic and manual evaluation to ensure that certain criteria are fulfilled. Namely, the data should be within a plausible range as well as showing plausible variability. Moreover, the data should reproduce to a reasonable extent the expected daily, monthly and yearly variations. Certain measured quantities which are related to each other are checked for consistency. Whenever possible, NABEL measurements are compared to nearby or similar stations with the expectation of similar values.
In addition to the daily average, the afternoon average of temperature, global irradiance, water mixing ratio and wind speed are calculated and used as meteorological variables as well. "Afternoon" values are the result of averaging hourly values between 12:00 and 18:00 (UTC+1). Such variables are hereafter denoted with the letter "a" at the beginning of their name. Averages for the morning (06:00-12:00 UTC+1) and evening (18:00-24:00 UTC+1) periods were calculated as well. These are denoted with "m" and "e" respectively.
Wind speed measurements were complemented by calculating wind gust. The maximum wind speed values of 10min intervals were recorded and the resulting six values every hour were averaged to obtain hourly wind gust. A further 24-h averaging was done to obtain the daily data.
It is known that the atmospheric boundary layer height and the strength of the temperature inversion at the top of the boundary layer play an important role in determining the levels of PM 10 by mixing and deposition mechanisms. As a proxy for the vertical gradient of the potential temperature, the difference of an approximate potential temperature between two nearby meteorological stations with large altitude difference was used (Ordóñez et al., 2005). Different pairs of stations were used for the areas north of the Alps (Swiss Plateau) and south of the Alps. A further measure of the vertical mixing in the boundary layer is the convective available potential energy (CAPE) which is estimated using meteorological soundings from the rural station at Payerne. A simple parcel method as described in Seibert et al. (2000) with a surface excess temperature of 0.1 K was used to calculate convective boundary layer depths from the same soundings.
A relationship between the transatlantic transport of ozone and the North Atlantic Oscillation (NAO) index has been demonstrated in the last few years in the literature (Creilson et al., 2003;Eckhardt et al., 2003). In order to investigate a possible influence on PM 10 concentrations in Switzerland, a monthly station-based NAO index (Jones et al., 1997) was used in the analysis.
The categorical variable "synoptic group" as calculated in the Alpine Weather Statistics AWS (Wanner et al., 1998) was used for representation of the impact of the general synoptic weather conditions on PM 10 in Switzerland. Moreover, variables "front" and "days since front" from the AWS were also included in this trend analysis to account for a possible current or past presence of a frontal passage. Another variable which records the history of the weather conditions is the amount of precipitation on the previous day, denoted here as "yesterday precipitation". Further details on the meteorological variables which were constructed and used in this study can be found in Ordóñez et al. (2005). For each station and season, only those meteorological variables with less than 30% missing values were included in the data analysis procedure described below.

Description of the statistical model
In order to properly investigate the relationship between PM 10 concentrations and meteorological variables one has to consider that these relationships can be non-linear. A way to model such non-linear relationships is to use Generalized Additive Models GAMs (Hastie and Tibshirani, 1990;Wood, 2006). GAM models were constructed for each station and season separately for the period from 1991-2008. For comparison purposes GAMs were run for the yearly data of all stations as well. All calculations were carried out using the statistical software R (R Development Core Team, 2008). The process yielded a total of 65 GAMs. The general formula of these models reads: .. : B ij denote categorical variables. Index i denotes the kind of categorical variable (e.g. synoptic group, front, day of the week etc.) and index j denotes the category. For example, categorical variable synoptic group has eight possible categories each one representing a different synoptic group. Therefore, if synoptic group is denoted by B 1j , then B 1j has eight possible values ranging from B 11 to B 18 . B ij is equal to 1 when the day in question is classified under category B ij and 0 otherwise. b ij is the corresponding coefficient. ε: error term A logarithmic link function was deemed suitable for GAMs in (1) because it improved the characteristics of the model residuals (see Sect. 3.2).
Functions s i are thin plate regression splines (see Wood, 2008 and references therein), which were estimated using the mgcv library of R. These functions are the result of a minimization procedure where smoothness is ensured by an appropriately selected penalization term. However, unlike a full thin plate spline, the basis is truncated to obtain a low rank smoother which is relatively fast to compute.
GAMs can be applied to a wide variety of data sets independently of location (including atmospheric data) under the assumption that the target variable has an additive relationship to the explanatory variables.

Model selection
A forward elimination procedure was used for the selection of the model covariates. This procedure is a modified version of the variable selection process put forward by Jackson et al. (2009). The employed algorithm includes the following 6 steps: Step 1 -The first variable was selected: All available explanatory variables were substituted one after the other into the model to construct GAMs with one variable. A version of the Akaike Information Criterion AIC (Akaike, 1974) was calculated for each model and the GAM with the lowest AIC was selected (see Eq. 3 for definition of AIC).
Step 2 -The next variable was added: The remaining explanatory variables were substituted to the existing model one after the other. The GAM with the lowest AIC was selected.
Step 3 -The choice from the first step was confirmed: The variable which was chosen in the first step was removed and each one of the remaining variables was added to the model. If this processes resulted in a model with lower AIC than the model selected in step 2, then the model from the current step was chosen.
Step 4 -Collinearity was tested: Collinearity occurs when two model covariates have a linear relationship. This makes the individual contribution of each variable difficult to discern, introduces redundancy and makes the model excessively sensitive to the data (Jackson et al., 2009). The Variance Inflation Factor VIF (Freund and Wilson, 1998) was calculated as a measure of collinearity of each covariate with each one of the other variables in the model. VIF was calculated according to the relationship where R 2 j is the coefficient of determination of the linear regression of the j -th variable in the model against the other variables. If a value of VIF = 2.5 was exceeded for a pair of covariates, the covariate whose addition to the model caused the smaller decrease of AIC was removed. The threshold of 2.5 for VIF corresponds to a coefficient of determination of 0.6 which was considered to represent a maximum acceptable degree of collinearity. Other, typically less strict cut-off values for VIF can be found in the literature. Montgomery and Peck (1992) suggest VIF = 5 and Kutner et al. (2004) suggest VIF = 10.
Step 5 -Steps 2 to 4 repeated: The aforementioned steps were repeated until the addition of an extra covariate resulted in an increase of AIC.
Step 6 -Julian day was added to the model: The PM 10 values do not only depend on weather conditions but also on a number of different processes, which were not included in the model such as emissions of PM 10 and precursors of secondary PM 10 . One additional covariate, the time variable Julian Day (here defined as the number of days since 1 January 1960) was included in the pool of available covariates to capture PM 10 variations due to any processes which were not explicitly accounted for. If this variable was not automatically selected in steps 1-5 of the variable selection algorithm, it was forcibly added in the end as an extra linear term. This was necessary as non-inclusion of this variable would return models whose error term instead of being random would potentially contain some interannual structure.
As a final step, it was checked that the residuals of the selected model fulfilled the criteria of normality, linearity and homoscedasticity (Wilks, 2006).
AIC was estimated using the relationship where L: the maximized value of the likelihood function for the estimated model edf: effective degrees of freedom of the model k: the penalty parameter Details on the calculation of AIC for GAMs can be found in Tibshirani et al. (1990) and Wood (2008). The second term on the right-hand side of Eq. (3) represents the severity of model penalization. Using the value k = 24 for the parameter k yields the classic AIC (Akaike, 1974) whereas k = log(n), n being the number of observations, yields the Bayesian Information Criterion BIC (Schwarz, 1978), which penalizes additional model terms more heavily. Nevertheless, both the use of AIC and BIC resulted in models with too many variables. A value of k = 15 was chosen here to sufficiently limit the number of selected variables in the models.
A similar model selection algorithm was used in Jackson et al. (2009). In that study however, the selection of the best model in each step of the algorithm was not based on AIC, but on deviance explained. The advantage of using the AIC in the model selection process is that, unlike simple use of deviance explained, the effective degrees of freedom of the model are taken into account. The effective degrees of freedom of a regression spline increase as wiggliness of the spline increases. Hence, using relationship (3) in the model selection process does not only control the number of variables in the model but also the wiggliness of the variables to be selected. This method has the disadvantage that it usually yields GAMs with lower deviance explained because the model with lowest AIC chosen in each step of the algorithm does not necessarily coincide with the model with the highest deviance explained. It was deemed however appropriate to choose a methodology which results in models with conceptually more robust behavior.

Adjustment of PM 10 concentrations for the effect of meteorology
Once the GAMs were constructed, PM 10 values for the period in question were adjusted for meteorology (i.e. the effect of meteorology was removed). The adjusted PM 10 values are the model estimates of the contribution of the Julian day plus the model intercept and the random model error. This reads: The adjusted PM 10 values include PM 10 variability that is not the result of meteorology but of other factors such as the anthropogenic emissions. To evaluate the model performance, all GAMs were developed using a random sample of 90% of the available observations. The remaining 10% of the observations of the explanatory variables were then substituted into the produced models to obtain modeled PM 10 concentrations. These were then compared with the corresponding PM 10 measurements (Sect. 4.4).

Important explanatory variables and their relationship to PM 10
As shown in Eq. (1), the contribution of each explanatory variable to the logarithm of PM 10 was calculated as an additive factor. Exponentiation of (1) yields a relationship between PM 10 concentrations and the selected explanatory variables which are then expressed as multiplicative factors (hereafter referred to as "PM 10 factors"), which contribute to an increase of PM 10 if they are greater than 1 and a decrease of PM 10 if they are less than 1. In this section the relationship between PM 10 concentrations and each one of these estimated factors is presented. The exponentiated output was preferred because presenting the results in terms of PM 10 concentrations was deemed more intuitive than a presentation in terms of the logarithm of PM 10 . Table 2 shows the frequency of occurrence of each explanatory variable in the GAMs corresponding to all 13 stations for each season and the whole year. A significant part of the available explanatory variables were chosen either infrequently or never. That is because they either did not have a firm relationship to PM 10 or because other variables with a similar relationship to PM 10 were preferred in the model selection process. The time variable "weekday" was only chosen once in winter whereas it was more frequently selected in summer, autumn and spring (4, 3 and 3 times respectively). All of the stations where weekday was selected are urban/street or rural/highway stations, which are more strongly affected by traffic emissions. Figure 2 shows the PM 10 factor of the "day of the week" variable along with the standard error estimates at all stations where it was selected for summer.
These PM 10 factors are expressed relative to the effect of Monday, which is taken as equal to unity. Therefore, there is no standard error for Monday. Note the increasing tendency of the influence on PM 10 during the weekdays. The maximum influence occurs on Thursday which contributes to an increase of 9% on PM 10 concentrations compared to Monday on average over the 4 stations shown in Fig. 2 for summer. Sunday contributes to a decrease of 17% compared to Monday for the same stations. "Day of the week" was selected 9 times for the yearly data and exhibited a similar weekly variation. Thursday contributes to an increase of 13% and Sunday ambient PM 10 concentrations are decrased by 12% compared to Monday on average over the 9 stations selected for the yearly data. Assuming that the primary aerosol emissions are similar during all the weekdays, the increase during weekdays might be attributed, on the one hand to the production of secondary aerosols and on the other to a possible buildup of aerosols in consecutive days with high emissions. A similar weekly PM 10 pattern was found by Zidek et al. (2002) for the city of Vancouver. The explanatory variables that were chosen for more than half (at least 7) of the models are summarized in Table 3. Figure 3 shows the PM 10 factors of the most important explanatory variables. To provide results from a representative sample of stations, the stations were separated into five site types: "urban street", "urban background", "suburban", "rural lowland" (which includes the two "rural highway" stations) and "rural mountain" (see Table 1). For each variable in Table 3, its relationship with PM 10 for one station from each site type is shown in Fig. 3. Each row of plots in this figure refers to the same variable and each column refers to the same site type. If the variable in question was selected at more than one station from the same type, the station where the variable was selected first is shown. If the variable appeared at two stations of the same type in the same order, the station where the variable caused the largest reduction of AIC is shown. In the plots in Fig. 3, the shaded areas indicate two standard errors above and below the best estimates. The values of the smooth PM 10 functions before exponentiation had the property of summing up to zero. A "rug" was added to the horizontal axis to show the distribution of the values of each explanatory variable across its range.
Wind gust is the most important variable in all seasons except summer where it is the second most important. Wind gust depends on average wind speed and on turbulence kinetic energy. Wind speed was also included as an explanatory variable; however, it was usually not preferred in the variable selection process. This is probably because wind gust carries information about turbulence in the atmosphere. Thus, including wind gust not only accounts for transport of PM 10 due to the mean wind but also for turbulent mixing. Note that it would not be possible to include both wind speed and wind gust in the GAMs because they are highly correlated and one of them would always be removed in step 4 of the algorithm (see Sect. 3). Hooyberghs et al. (2005) show a negative linear relationship between wind speed and PM 10 concentrations. The relationship found here between wind gust and PM 10 is approximately exponential (Fig. 3).
The second most important variable, chosen in spring, autumn and winter is CBL depth. CBL depth affects PM 10 concentrations by limiting the volume of turbulent mixing of PM 10 . Figure 3 shows that an increase in CBL depth results in a decrease of PM 10 as expected. Non-linear relationships are encountered at many stations. This could be attributed to the fact that CBL depth was calculated of data from Payerne, a station in the Swiss plateau north of the Alps whereas the actual boundary layer characteristics of the stations included in this study are considerably affected by the presence of cities, mountains, valleys, lakes and other topographic particularities. Moreover, the important secondary aerosol production in Switzerland (Lanz et al., 2007 does not have a linear dependence on the aerosol precursors. Unlike all other stations, at the rural mountain station Rigi, the relationship between PM 10 and CBL depth in autumn is positive.
Given that Rigi is located above 1000 m, the levels of PM 10 are affected by whether the station is in the boundary layer or above, rather than by the boundary layer depth per se.
Two other important variables that appear in spring, autumn and winter are daily precipitation and yesterday precipitation. Precipitation reduces considerably PM 10 concentrations by wet deposition. In addition, if it is frontal precipitation, some polluted boundary layer air is replaced by clean air from aloft. As shown in Fig. 3, PM 10 is particularly sensitive to low values of precipitation. PM 10 levels reduce significantly for daily precipitation up to around 5 mm whereas for larger values of precipitation its effect diminishes. This feature would not be possible to depict by a linear relationship. For precipitation above 10 mm the effect on PM 10 appears to be positive or strongly fluctuating at many stations. However, care should be taken when interpreting the results in regimes where data are sparse because the smooth functions have a large uncertainty in such regimes. Rugplots in Fig. 3 show that data are sparse above 10 mm of hourly precipitation and that the standard error becomes considerably large. The frequent selection of yesterday precipitation in the models indicates that, when investigating air pollution levels, it is not only important to take into account the current weather conditions but also the recent history of the weather. Current day precipitation and yesterday precipitation were not selected in summer where convective precipitation is encountered more frequently (Beniston et al., 1994). Convective precipitation typically lasts for a short time and exhibits a large spatial variability allowing for the presence of potentially high PM 10 concentrations before the rain and some accumulation of new PM 10 during the same day after the rain.
An important variable in spring and winter is net irradiance. Notably, this variable has opposite effects on PM 10 during winter and spring. During spring, increasing values of net irradiance are associated with a strong decrease of PM 10 concentrations. A possible explanation is that strong net irradiance contributes to large sensible heat flux which in turn drives the formation of deeper atmospheric boundary layers. However, in winter the relationship between PM 10 concentrations and net irradiance is reversed with larger positive values of net irradiance contributing to an increase of PM 10 . During winter, anticyclonic situations often create well-marked drainage flows from the Alps towards the Swiss Plateau. The cold air at the surface creates strong inversions, which in turn lead to radiation fog, high inversion fog (Hochnebel) and high levels of air pollution (Wanner and Hertig, 1984). While low negative values of net radiation occur with clear skies, the presence of Hochnebel introduces a positive term of long wave radiation in the net irradiance balance, which increases net irradiance. It is possible that through this mechanism, relatively high values of net irradiance lead to high concentrations of PM 10 .
A further important explanatory variable for spring is relative humidity. Ambient aerosol mass concentrations of hygroscopic aerosol depend on relative humidity. On the one hand, such a dependency should not be detected by filter measurements. On the other hand, it is possible to see a dependency on relative humidity with beta attenuation monitors. For example, at low-temperature and high-relative humidity ambient conditions, air sucked into the measurement station warms-up to room temperature which in turn leads to a reduction of relative humidity and evaporation of water from hygroscopic aerosols. The relationship between PM 10 and relative humidity is not the same for all stations (Fig. 3). Namely, there is a negative relationship for all sites except Lugano, Magadino (not shown) and Chaumont. These three sites exhibit a positive relationship for relative humidity below 60% and a negative relationship for relative humidity above 80% approximately. These relationships could be affected by the regional climate of these three sites. Southern Alps (where Lugano and Magadino are located) are known to be influenced by North Foehn episodes where clean air with low relative humidity descends from the Alps (Weber and Prévôt, 2002). North Foehn days amount to 9.7%, 7.6% and 4.2% of the days in March, April and May respectively. Chaumont PM 10 daily concentrations depend on whether the site is within the boundary layer or not. Polluted boundary layer air after ascending to elevated sites has increased relative humidity. In spite of the aforementioned information, the relationship between PM 10 concentrations and relative humidity in spring remains unclear and requires further investigation.
Variable "afternoon temperature" was frequently selected in summer and a positive relationship between PM 10 and afternoon temperature was identified at all stations (Fig. 3). Temperature has a particularly strong influence on PM 10 in summer since a temperature increase from 10 to 35 • C increases PM 10 by a factor of 4. This relationship can be the result of secondary aerosol production during warm summer days. The emissions of biogenic Volatile Organic Compounds (VOC), which amount to an estimated 24% of the total VOC emissions in Switzerland (Andreani-Aksoyoglu and Keller, 1995;Andreani-Aksoyoglu et al., 2008), depend on a number of meteorological variables which include temperature (Guenther, 1997). Besides, biogenic VOC include compounds with high yield, such as monoterpenes. Further explanations include good weather conditions (e.g. absence of frontal passages that could otherwise reduce PM 10 concentrations) and enhanced resuspension of soil and road dust in warm days during summer (Vardoulakis and Kassomenos, 2008).
Temperature was selected 8 times in winter. Similarly to afternoon temperature in summer, temperature in winter has a marked effect on PM 10 concentrations. Nevertheless, the relationship this time is inversed; low temperatures contribute to high PM 10 concentrations by condensation of volatile compounds and possibly because of increased wood burning for space heating. Moreover, low temperature here can act as a proxy for stable boundary layer and low vertical mixing.
Most of the important variables selected for the whole year GAMs were also selected in the other seasons and their smooth functions are similar (not shown). The smooth function for temperature is an exception as the opposite response of PM 10 to low and high temperatures is combined to give a negative relationship for low temperatures, positive relationship for high temperatures and a stationary point in between. Variable "days since front" that was only selected for whole year GAMs shows a decrease of PM 10 for "days since front" less than 1 and an increase for "days since front" greater than 2. This indicates removal of PM 10 by enhanced transport and wet deposition in the presence of a front and accumulation of PM 10 after the frontal passage.
Julian day was selected often in summer implying that less PM 10 variability is explained by meteorology during this season compared to winter. Factors that are not included in meteorology such as changes in anthropogenic emissions appear to have an important role in summer. The Julian day smooth functions provide a useful tool for the investigation of the meteorologically adjusted PM 10 trends because they are able to detect non linear evolution of PM 10 concentrations with time. The Julian day smooth functions are shown in Fig. 4. The stations shown for summer are selected in the same way as in Fig. 3 except that the stations are now split into six categories based on local emission characteristics. Namely, rural highway stations (Haerkingen and Sion) make a separate category. In the other seasons, all stations where Julian day was selected by the variable selection algorithm are shown.
All Julian day smooth functions shown in Fig. 4 have an overall decreasing trend which is interpreted to be at least to some extent the result of reduced PM 10 emissions due to the implementation of emission reduction measures and technological improvements (e.g. vehicles with more efficient engines and catalytic converters). The change in mostly traffic and industrial precursor emissions (NO x , SO 2 , VOCs) of secondary aerosols is very likely responsible for the trend as well (Vestreng et al., 2007(Vestreng et al., , 2008Stemmler et al., 2005). The reduction of the PM 10 factors is more pronounced in most of the 1990s whereas it is not as steep or even increases in the late 1990s and early 2000s. This evolution of PM 10 is consistent to the GAINS model calculations of PM 10 emissions and to the emissions of NO x , non-methane volatile organic compounds, SO 2 and NH 3 as estimated by the European Monitoring and Evaluation Programme (Monks et al., 2009). It is not clear why at some stations PM 10 concentrations increased in the early 2000s. In a source apportionment study which was carried out for London, UK similar trends were identified and the increase of PM 10 in the aforementioned period was attributed to increased traffic emissions (Fuller and Green, 2006). In Switzerland, the share of diesel vehicles has been increasing in the last two decades and the increase in the period between 2000 and 2009 has been faster than the increase in the period between 1990 and 1999 (Fig. 5). This is particularly true for diesel passenger cars (Fig. 5) whose share was increasing by 0.099% yr −1 on average in the 1990-1999 decade whereas in the 2000-2009 decade it is increasing by 1.4% yr −1 as estimated by linear regression of diesel passenger car share versus year for the two decades respectively. A further local maximum appears in 1995 or 1996 at some stations in summer and spring.
To identify an overall trend in PM 10 emissions, GAMs were created for yearly data and the corresponding Julian day smooth functions were estimated. This is shown in the last row of plots in Fig. 4. Strong seasonal features such as the 1995 maximum in the spring Lausanne PM 10 factor are largely smoothed out in the yearly trend. Also in the yearly data there is an overall decrease of PM 10 that is interrupted by two local maxima around 1996 and 2003. Two urban background stations, Zurich and Lugano (not shown) did not show a maximum in the 1990s. To better quantify and summarize the PM 10 trends for all stations and seasons, an analysis of the PM 10 trends by linear regression of the daily averages of PM 10 values before and after meteorological adjustment is presented in Sect. 4.3. Considering the generality of the selected sets of explanatory variables and their relationship to PM 10 , they are based on data from the NABEL measurement network and they are best applicable to the aforementioned sites. Nevertheless, they represent universal mechanisms of dispersion, dilution and formation of PM 10 and are therefore expected to be suitable for other locations with similar climatic characteristics. However, testing is required before using statistical models at conditions different from those they were initially developed for. Considering the applicability of the proposed GAMs to other species, it is limited because deposition, resuspension and the response of various species to some meteorological variables such as temperature, solar irradiance and relative humidity depends strongly on the examined species.

Model performance
To evaluate the GAMs, the deviance explained (Wood, 2006) was used. Deviance explained quantifies which part of the observed variance is reproduced by the GAM. As shown in Table 4, winter was the season with the highest deviance explained reaching on average 65%. Summer had an almost equally high deviance explained of 64% and at the same time a relatively small number of 6.4 selected explanatory variables. In spring an average of 8.2 selected covariates could explain on average only 56% of the observed variance, which indicates that the model has some difficulty in explaining the PM 10 concentrations by means of the variables that are available for description of the meteorological variability. The yearly GAM runs yielded models with a high 65% deviance explained on average and a high number of selected explanatory variables (11.1). It should be noted that the deviance explained varied considerably from station to station. For example, in winter where the highest average deviance explained was encountered, the rural mountain stations of Chaumont and Rigi had rather low values of deviance explained (38% and 48% respectively) indicating that for the special weather conditions in mountainous areas the available meteorological variables are not capable of describing the effect of meteorology on PM 10 concentrations. Moreover, due to the sites' elevation, PM 10 concentrations are not only affected by local meteorology and emissions but by long-range transboundary transport as well.
Another way to assess model performance is to estimate the extent to which the year-to-year variability of meteorologically adjusted values is reduced. For this purpose, the mean squared error MSE (Wilks, 2006) was calculated according to the formula where m is the number of data points used for the linear regression of the daily average PM 10 concentrations versus time (see Sect. 4.3) and R i are the model residuals defined as the difference of the predicted by the linear model minus the observed daily PM 10 . MSE was calculated for the observed and adjusted PM 10 values. The ratio of adjusted over measured MSE is tabulated in Table 4. The PM 10 variability after meteorological adjustment in winter is only 27% of the observed variability. In spring, despite the relatively low deviance explained of the model estimates, the MSE ratio is 38% indicating a strong reduction of the variability of the PM 10 values after meteorological adjustment. A smaller reduction of variability is observed in summer and autumn. The variability of the raw PM 10 concentrations in those seasons is already rather low, as indicated by the respective MSE meas values.
As already mentioned, 14-18 years of data were available for all stations. Although this is considered to be a sufficient amount of data to construct reliable GAMs, it is expected that such an analysis can be carried out with less data too. For example, data from 1994 until 2008 are available for Zurich station. A series of GAMs that were constructed using winter data since 1996, 1998, 2000, 2002 and 2004 showed that the  resulting GAMs were quite consistent for the datasets starting in 2002 or for longer datasets i. e. for a minimum of seven years of data. Even fewer years of measurements are needed for analyzing full year datasets.

PM 10 trends before and after meteorological adjustment
The trends of the PM 10 concentrations in the period from 1991 to 2008 were estimated by linear regression of the daily average PM 10 concentrations versus time. This was done for PM 10 measurements and for the PM 10 values that were adjusted for meteorology using relationship (4). The regression equations read daily average measured PM 10 = c 1 + c 2 Julian day + ε daily average adjusted PM 10 = d 1 + d 2 Julian day + ε The estimated slopes c 2 ,d 2 of the regression lines were multiplied by 365 to represent yearly PM 10 changes. The results of this process are summarized in Fig. 6 for all seasons and all stations. One should keep in mind that investigating the linear trends is a simplification because, as already discussed, the evolution of PM 10 concentrations during the period in question is often non-linear, especially in winter. All slopes are shown with their 95% confidence intervals. In addition to the results for each station, the last value called "AVG" in each plot indicates the slope averaged over all stations. The confidence interval of the average was calculated by error propagation (e.g., Bevington and Robinson, 2003). This reads where δ(PM 1 0) avg is the average upper or lower confidence interval, N = 13 is the number of stations and δ(PM 10 ) i is the upper or lower confidence interval of the i-th station. Figure 6 shows that PM 10 trends are declining at almost all stations and seasons except the winter trend at the rural stations of Chaumont and Magadino. The confidence intervals of the measured and adjusted values largely overlap making the existence of a change in the PM 10 trends due to the adjustment in many cases uncertain. However, it is discussed later that for some stations the meteorological adjustment had a considerable effect. In addition, the slopes of the adjusted trends have smaller uncertainty, which implies that the adjusted values can be used to identify significant trends from shorter time series.
Winter is the only season where the meteorological adjustment increased (in absolute terms) the slope of the PM 10 trend line on average and at most stations. The largest changes were observed at the rural highway stations Haerkingen and Sion. Rural highway stations provide ideal locations for examining air pollution trends and the effect of weather because on the one hand the meteorological observations, not being affected by local urban features, are more accurate and representative and on the other hand strong local emissions of PM 10 make it easier to detect appreciable changes in PM 10 concentrations. Considerable reduction of the slope was also observed in Magadino which is a rural station but its high annual average PM 10 concentrations indicate that it is affected by relatively strong regional PM 10 emissions and production. The stations with the largest yearly reduction of PM 10 in winter (and in all other seasons) were the urban street stations of Bern and Lausanne. The fact that urban street and rural highway stations have the fastest decreasing PM 10 concentrations indicates a reduction in traffic emissions of primary PM 10 . Therefore, regulatory action and technology advances regarding vehicle emissions have had a considerable success in reducing PM 10 concentrations.
The slopes of the PM 10 trend lines after adjustment for meteorology are around zero at Lugano, Magadino and Rigi and positive at Chaumont. Since this behavior is only encountered in winter, it is reasonable to assume that the emission sector which caused this neutral trend is space heating. A particular space heating source to be considered is woodburning since recent measurements (Szidat et al., 2007;Lanz et al., 2008Lanz et al., , 2010Sandradewi et al., 2008) demonstrate that wood-burning is a major source of particulate matter in winter both in the Swiss plateau and in the Alpine valleys. Unlike traffic emissions of PM 10 and precursors, wood burning units have not reduced in the last two decades in Switzerland, their total number remaining almost stable at around 700 000 wood burning appliances (Primas et al., 2009). It should be noted that the PM 10 trends in Lugano, Magadino, Rigi and Chaumont do not only depend on the meteorology and the emissions in Switzerland but they are considerably affected by neighboring countries too. Air pollution levels in the area of Ticino south of the Alps (including Lugano and Magadino) are considerably affected by the Milan area (Prévôt et al., 1997;Andreani-Aksoyoglu et al., 2004). Rigi and Chaumont are located at relatively high altitudes and are therefore affected by long-range transport. In addition, Chaumont is located close to the western border of Switzerland. Considering the predominant southwestern wind direction in the area (MeteoSwiss, 2010a), it is expected that PM 10 levels at this station are rather affected by emissions in France than in Switzerland.
The finding that meteorological adjustment of PM 10 trends in winter makes these trends more negative implies that the weather during winter in the 1990 and 2000 decades in Switzerland has changed in such a way that favors the accumulation of PM 10 . A possible explanation for this feature is the evolution of the North Atlantic Oscillation (NAO) index, which has reduced considerably during the 1990 decade (Hurrell and Deser, 2010). High values of NAO index are linked to mild and wet winters whereas low values of NAO are linked to cold and dry winters which tend to increase levels of PM 10 . NAO was one of the available explanatory variables which were not selected from the variable selection algorithm (see Table 2). This is attributed to the fact that the weather conditions associated with various NAO modes were described more directly by other explanatory variables. The yearly trends shown in Fig. 6 have some similarities to the winter trends. The yearly trends before and after meteorological adjustment are tabulated in Table 5. On average, a small reduction of the PM 10 trend after adjustment for meteorology is observed. The stations which actually have a considerable reduction of their PM 10 trends after adjustment for meteorology are the rural highway stations of Haerkingen and Sion where the adjusted trends are lower by 0.12 and 0.19 µg m −3 yr −1 respectively (Table 5). This probably reflects the aforementioned winter behavior of the PM 10 trends in those stations because in the other seasons the meteorological adjustment, although narrowed considerably the confidence intervals, only caused little change in the best estimates of the slope. Velders and Matthijsen (2009) Table 5). Stedman (2002) reported decreasing PM 10 trends for 8 city centre stations in the UK, these trends ranging between −0.9 µg m −3 yr −1 and −1.9 µg m −3 yr −1 for data between 1992 or 1993 until 2000. These trends are also comparable (although somewhat larger) to the trends found in Swiss cities. In a further related study, Spindler et al. (2004) report a negative trend of roughly −1.5 µg m −3 yr −1 based on measurements at some stations in the area of west Saxony, Germany, although an exponential instead of a linear relationship was actually chosen to describe the decline of PM 10 in that area.
No literature discussing PM 10 trends on the basis of seasonal data is known to the authors. The average values of the adjusted trend per season found in this study are −0.54 µg m −3 yr −1 , −0.58 µg m −3 yr −1 , −0.52 µg m −3 yr −1 , −0.35 µg m −3 yr −1 and −0.53 µg m −3 yr −1 for spring, summer, autumn, winter and yearly data respectively.
In autumn the adjusted trend at the urban street station Lausanne is less decreasing. A similar effect of reduced slope after meteorological adjustment is also evident in the stations of Bern and Magadino in summer and Lausanne in spring. Lausanne and the other urban street station, Bern have the steepest trends among all stations in all seasons implying considerable PM 10 reductions from the traffic sector.
The importance of the aforementioned long-term PM 10 changes at each site depends on the average level of PM 10 concentrations at the site in question. In order to put the PM 10 trends in this context, the adjusted trends were calculated as the percentage of the mean PM 10 observations during the period in question (see Table 1). The results are tabulated in Table 5. The urban street stations at Bern and Lausanne apart from having the largest trends, they also have the largest percentage trends amounting to 3.1 and 3.3% yr −1 respectively. The average percentage trend over all stations is −1.9% yr −1 with a range of −1.2 to −3.3% yr −1 .
As already noted, visual inspection of Fig. 4 reveals that the adjusted PM 10 trends do not have a uniform behavior throughout the period between 1991 and 2008. A further examination of the adjusted PM 10 trends was carried out by separating the PM 10 GAM estimates before 2000 and from 2000 onwards and estimating the line fitted by linear regression of PM 10 concentrations versus time in each period. At most stations the trends before and after 2000 differ by less than 0.4 µg m −3 yr −1 , in absolute terms (Table 5). Bern and Lausanne however exhibit considerably larger differences; The PM 10 trend at Bern increased from −3.4 µg m −3 yr −1 in 1991-1999 to −1.1 µg m −3 yr −1 in 2000-2008. The slope at Lausanne similarly increased from −2.8 µg m −3 yr −1 in 1994-1999 to −0.4 µg m −3 yr −1 in 2000-2008. In addition to the slopes of the linear regressions, it is useful to see the effect of the meteorological adjustment on the actual PM 10 concentrations. Figure 7 shows the yearly averaged PM 10 concentrations before and after the adjustment for meteorology and their trend lines. One station from each one of the six categories mentioned in Sect. 4.1 is included. The stations with the largest changes in their PM 10 trend after meteorological adjustment were preferred.
Winter has the highest PM 10 variability and the highest reduction of PM 10 variability after the meteorological adjustment (see MSE values in Table 4). For example, as shown in Fig. 7, extreme high PM 10 concentrations in winter seasons 1992, 1997, 2003 and 2006 in Bern were heavily reduced after meteorological adjustment. Similar effects are seen in other seasons too (e.g. spring 1996 in Lausanne station). Examination of all plots in Fig. 7 (and the remaining plots which are not included in this figure) shows that GAMs have a good skill in removing PM 10 variability even in cases where PM 10 concentrations were unusually high because of weather conditions.

Model validation
The models described above were developed using a random 90% of the available observations. The remaining 10% of the observations were used to validate the models obtained from the 90% of the observations. In this process, the values of the explanatory variables which were chosen in the GAMs were substituted into the GAM formulas to give predictions of PM 10 values. These predictions were then compared to the corresponding PM 10 measurements.
A graphic comparison between model predictions and observations is done by examination of scatter plots (Fig. 8).
A one-to-one line was added to facilitate comparison to the ideal model. The factor-of-two scatter is indicated by dashed lines. Some points are missing from these plots either because of missing PM 10 measurements or missing meteorological observations in the selected 10% of the days. The stations included in Fig. 8 are the same as in Fig. 7.
To validate quantitavely the GAMs, the fraction of daily GAM predictions which lie between 50% and 200% of the corresponding measurements was calculated. This is the so-called "factor of 2" (FAC2) statistic (Chang and Hanna, 2004) and reads: FAC2 = fraction of data that satisfy 0.5 ≤ PM 10 predicted PM 10 observed ≤ 2.0 (9) It follows that for an ideal model, FAC2 = 1. Data points which lie in the area within the dashed lines in the scatterplots in Fig. 8 satisfy Eq. (9). Visual inspection of the scatterplots in Fig. 8 reveals a very good predictive skill of the GAMs since the vast majority of the data fall within the factor of two area marked by the dash lines. However, in most stations in summer, autumn and spring, some underprediction of PM 10 values is observed for high (above 30-40 µg m −3 ) PM 10 measurements. Such deviation can occur in cases where irregular emissions such as, for example, agricultural burning or fireworks impact the measurement site.
The FAC2 values are very high, often approaching unity. Their averages over all stations are 0.90, 0.98, 0.92, 0.93 and 0.93 for spring, summer, autumn, winter and yearly data respectively. The FAC2 values for some individual stations are included in Fig. 8. FAC2 values at the mountainous sites at Rigi and Chaumont, although still satisfactory, are somewhat lower (usually between 80%-90%).
The linear dependence of the predicted versus the measured PM 10 concentration was estimated using the Pearson correlation coefficient. A relatively high linear dependence was found with the average correlation coefficient values over all stations being 0.62, 0.78, 0.72, 0.72 and 0.72 for spring, summer, autumn, winter and yearly data respectively.

Conclusions
The trends of PM 10 and the effect of meteorology were investigated for 13 air quality stations in Switzerland for the period between 1991 and 2008. This was done by constructing Generalised Additive Models (GAMs) which incorporate non-parametric relationships between PM 10 and meteorological and time variables. These relationships were based on quality-checked meteorological and PM 10 long-term observations including a wide spectrum of meteorological variables with good spatial and temporal coverage.
Most of the estimated functions provided fairly intuitive answers with regards to which direction PM 10 concentrations change when a change of a certain meteorological variable occurs. For example, it was found that PM 10 concentrations have negative relationships with wind gust, yesterday precipitation and convective boundary layer depth whereas they have positive relationships with afternoon sunshine duration and afternoon temperature. The contribution of this study was to provide further insight in those relationships by quantifying them and accounting for non-linearities. For instance, a small precipitation rate strongly reduces PM 10 concentrations whereas a further increase in the precipitation rate has no further considerable effect. Another interesting finding is the importance of each meteorological variable in the different seasons. Apart from wind gust which was one of the most important variables in all seasons, large interseasonal differences were found with regards to the preferred variables. The convective boundary layer depth was important in winter, spring and autumn while PM 10 concentrations appear to be not so sensitive to changes in the boundary layer depth during summer where deep convective boundary layers form relatively often. The rate of precipitation on the current and the previous day was frequently selected for winter, autumn and spring and not so frequently selected for summer where convective precipitation is more frequent.
The variability of PM 10 unrelated to meteorology (represented by the Julian day variable) was found to be important in summer and not so important in spring, winter and autumn. Examination of the Julian day smooth functions showed a non-linear trend in PM 10 concentrations. A relatively fast yearly reduction of PM 10 was observed for the years until 2000 and a slow reduction from 2000 on, especially at the urban street stations Bern and Lausanne.
The GAMs showed a good ability to reduce the observed PM 10 variability in winter. The reduction of the PM 10 variability after the adjustment for meteorology was lower in the other seasons, which anyway had a rather low variability of the raw PM 10 values.
The strong effect of meteorology in PM 10 concentrations was also reflected on the PM 10 trends in winter. Removal of the effect of meteorology resulted in a more negative trend at most stations. In almost all seasons and stations, there was a downward trend of the adjusted PM 10 concentrations during the period in question. The seasonal trends averaged over all stations amount to −0.54 µg m −3 yr −1 , −0.58 µg m −3 yr −1 , −0.52 µg m −3 yr −1 , −0.35 µg m −3 yr −1 and −0.53 µg m −3 yr −1 for spring, summer, autumn, winter and yearly data, respectively. These trends correspond to respective percentage changes of −2.0% yr −1 , −2.2% yr −1 , −1.8% yr −1 , −1.1% yr −1 and −1.9% yr −1 .
The reduction of PM 10 after adjustment for meteorology demonstrates that air pollution abatement policies in Switzerland have been effective to a certain extent. This also applies to EU countries surrounding Switzerland which affect air pollution levels in Swizterland through transboundary pollution. Air pollution abatement measures in Switzerland include the adoption of cleaner road vehicle technologies (Euro 1, 2, 3, and 4 emission standards), promotion of public transport and improvements in the rail and road infrastructure (Berger et al., 2009). The continuation of regulatory action towards the reduction of air pollution is expected to contribute to a further decrease of PM 10 in the future. For instance, as of 2011 only Euro 5 vehicles will be imported in Switzerland (Swiss Federal Council, 2010b). An important feature of the Euro 5 technology is that passenger cars and light duty vehicles using diesel have to be equipped with particle traps. As a large fraction of particulate matter is due to secondary aerosols from precursor gases like NO x , SO 2 and NH 3 , future trends will depend a lot on the future emissions on these precursors. Especially for NO x there is still a large reduction potential. Further recent measures include stricter regulation for emissions from construction equipment and stricter regulation for wood burning emissions (Swiss Federal Council, 2010a). The road infrastructure is planned to develop further so that traffic jams are reduced (Swiss Federal Council, 2010c). Some relief of the road network is also expected by the planned increase of the capacity and speed of rail transport (Swiss Federal Council, 2010c). Nevertheless, the actual PM 10 trends in the future will also be affected by the increasing demand for transport and other particle and precursor emitting activities (Swiss Federal Council, 2010c).
The prognostic skill of the model was examined by model validation. The FAC2 statistic for different seasons ranged between 0.90 and 0.93 whereas the Pearson correlation coefficient ranged between 0.62 and 0.78 on average over all stations. Therefore, GAMs have considerable potential as a cheap tool for air quality forecasting, provided that a good weather forecast is available.