Interactive comment on “ Spectral representation of the annual cycle in the climate change signal ” by T

In the manuscript “Spectral representation of the annual cycle in the climate change signal” T. Bosshard, S. Korlarski, T. Ewen and C. Schär suggest a parametric model based on harmonic functions to describe the annual cycle in temperature and precipitation data. Compared to a moving average, this approach yields a smoother estimate of the annual cycle, is less sensitive to sampling artefacts and allows a more robust estimation of the change signal of control and scenario runs. A nice illustration of the sampling problem for moving averages is given on the basis of synthetically generated precipitation series. The order of harmonic functions is chosen by cross validation on a set of observations. For a large set of various GCM-RCM chains, the parameters


Introduction
Impacts of climate change on the hydrological cycle are both of high scientific interest as well as of high relevance for society as a whole.The former is due to the intimate coupling of the hydrological cycle and the climate system (Allen and Ingram, 2002;Wentz et al., 2007;Wild et al., 2008) while the latter is based on the manifold interactions between the anthroposphere and the hydrosphere (Kundzewicz et al., 2007).
Correspondence to: T. Bosshard (thomas.bosshard@env.ethz.ch)Hydrological impact studies focussing on runoff often use statistically post-processed global climate model (GCM) or regional climate model (RCM) data to drive a hydrological model (Hay et al., 2000;Leung et al., 2004;Wood et al., 2004;Kay et al., 2006;Buytaert et al., 2010).For this purpose, various statistical post-processing methods have been developed (see e.g.Fowler et al., 2007or Maraun et al., 2010 for comprehensive reviews).All these methods are based on statistical relationships that bridge the spatial and temporal gaps between observations and modelled data, and attempt to correct for climate model biases.Most of the available methods focus on the hydrometeorological variables temperature and precipitation (abbreviated as T and P respectively in the remainder of this article) and usually include some representation of the annual cycle.
Natural variability, both on interannual as well as intraannual time scales, impairs parameter estimates of the statistical post-processing methods.The range of the natural variability can be assessed using e.g.resampling techniques.Prudhomme and Davies (2009) and Wood et al. (2004), for example, resampled observed time series to estimate the range of natural variability of the climate change signal.Cross-validation has also been used to test the robustness of the parameter estimates to interannual variability (Terink et al., 2010;Schmidli et al., 2007;Widmann et al., 2003).However, only a few studies focussing on hydrological impacts have looked in detail at the intraannual variability of the parameters.Smoothing by averaging over seasons (see e.g.Schmidli et al., 2007) or months (see e.g.Middelkoop et al., 2001;Kleinn et al., 2005) is a common practise.An appropriate representation of the seasonal cycle, however, is not straightforward.On the one hand, the optimal choice of the averaging period is dependent on the magnitude of the Published by Copernicus Publications on behalf of the European Geosciences Union.
natural variability, the spatial averaging, and the length of the data records.The stronger the natural variability, the smaller the spatial averaging area and the shorter the data record is, the wider the averaging window has to be chosen in order to reduce the effects of natural variability on the parameter estimates.On the other hand, hydrological impact modellers are interested in an accurate representation of the annual cycle and therefore prefer as narrow averaging windows as feasible.The optimal solution is thus not trivial to find and depends on the region and application under consideration.Despite its importance, the discussion of how to optimally represent the annual cycle in climate change scenarios is often neglected in recent impact modelling.In many cases, the averaging window width is mentioned without specific justification (e.g.Cameron et al., 2000;Jasper et al., 2004;Graham et al., 2007).
This paper elaborates on the representation of the annual cycle in the climate change signal within the delta change post-processing methodology.We chose the delta change method because of its simplicity, but the results appear relevant for more sophisticated methods as well.The delta change method has been used for hydrological impact studies ever since GCM data became available, and it is still used nowadays (Gleick, 1986;Hay et al., 2000;Prudhomme et al., 2002;Lenderink et al., 2007).More sophisticated combinations of the delta change approach and weather generators have also been developed (Kilsby et al., 2007).It is noteworthy that Gleick (1986) already stressed the importance of representing the climate change throughout the annual cycle since seasonal changes tend to cancel each other out in the annual average.
Here, we test the influence of sampling variability on the annual cycle of the climate change signal by using moving averages (MA) of different window widths.As an alternative to the MA, we present a spectral approach to estimate the climate change signal.Our analysis is carried out at observational station sites in Switzerland.The spectral estimation produces smoother annual cycles of the climate change signals than MAs.
The paper is structured as follows: in Sect.2, we present the data used for the study.Section 3 introduces the delta change method and the estimation methods for the annual cycle.In Sect.4, we study the effects of sampling variability on the estimation of the annual cycle using a stochastic rainfall generator.Section 5 presents the estimation of the annual cycle of the climate change signal using a spectral smoothing method and a comparison to estimates using MAs of 31 days window width.The methodology is then systematically explored at Swiss station sites.Section 6 summarises the findings and discusses their relevance for climate impact studies.The geographical focus of our study is Switzerland.Throughout this paper, the estimation of the climate change signal is based on GCM-RCM data interpolated to station locations of MeteoSwiss (see Fig. 1).All the stations provide T and P data with at least daily resolution for the period 1980-2009.We used the four nearest gridpoints and the inverse distance weighting interpolation algorithm to spatially interpolate the GCM-RCM data to station locations.It should also be noted that any height correction is redundant since constant correction terms cancel each other in the delta change approach.Also, for simplicity, we neglect leap days in the data, unless stated otherwise.
For the stochastic rainfall generator experiments, two subsets of precipitation stations are used to estimate the rainfall generator parameters.These subsets are indicated by blue and red dots in Fig. 1 (see also Sect.4).In addition, we used long-term data series from 26 stations with records going back to 1900 from the climate monitoring network of MeteoSwiss to constrain the harmonic smoothing model.sible.The optimal solution is thus not trivial to find and depends on the region and application under consideration.Despite its importance, the discussion of how to optimally represent the annual cycle in climate change scenarios is often neglected in recent impact modelling.In many cases, the averaging window width is mentioned without specific justification (e.g.Cameron et al., 2000;Jasper et al., 2004;Graham et al., 2007).This paper elaborates on the representation of the annual cycle in the climate change signal within the delta change post-processing methodology.We chose the delta change method because of its simplicity, but the results appear relevant for more sophisticated methods as well.The delta change method has been used for hydrological impact studies ever since GCM data became available, and it is still used nowadays (Gleick, 1986;Hay et al., 2000;Prudhomme et al., 2002;Lenderink et al., 2007).More sophisticated combinations of the delta change approach and weather generators have also been developed (Kilsby et al., 2007).It is noteworthy that Gleick (1986) already stressed the importance of representing the climate change throughout the annual cycle since seasonal changes tend to cancel each other out in the annual average.

Stars mark these stations in
Here, we test the influence of sampling variability on the annual cycle of the climate change signal by using moving averages (MA) of different window widths.As an alternative to the MA, we present a spectral approach to estimate the climate change signal.Our analysis is carried out at observational station sites in Switzerland.The spectral estimation produces smoother annual cycles of the climate change signals than MAs.The paper is structured as follows: In section 2, we present the data used for the study.Section 3 introduces the delta change method and the estimation methods for the annual cycle.In section 4, we study the effects of sampling variability on the estimation of the annual cycle using a stochastic rainfall generator.Section 5 presents the estimation of the annual cycle of the climate change signal using a spectral smoothing method and a comparison to estimates using MAs of 31 days window width.The methodology is then systematically explored at Swiss station sites.Section 6 summarises the findings and discusses their relevance for climate impact studies.

Data
We

The delta change method
In the delta change method, observational records are scaled according to a climate change signal.The climate change signal is derived from climate model data as the change between a scenario period (SCE) and a control period (CTL).As a result of the scaling, the spatio-temporal patterns as well as the correlations between the variables closely follow the observed records.Thus, the delta change method is considered a robust method to generate climate impact scenarios (Graham et al., 2007).
In this study, we applied the delta change method at station sites for the SCE periods 2021-2050 and 2070-2099, both relative to the CTL period 1980-2009.At each station i, for each GCM-RCM j and for each day d in the year, we estimate the mean annual cycle of the variable of interest and denote it with X CTL i,j (d) for the control and X SCE i,j (d) for the scenario period where X stands for either T or P .The delta change method then derives an additive ( X add i,j (d)) and a multiplicative ( X mult i,j (d)) climate change signal for T and P , respectively, according to .
(2) Let X CTL i,obs (y,d) denote the continuous observational time series at station sites in the CTL period 1980-2009.Here, y represents the years in the CTL period.In the delta change method, all observational time steps in the CTL period belonging to the same day d in the year are scaled with the corresponding climate change value.Again, one commonly uses an additive or multiplicative scaling for T and P , respectively: (y,d) = P CTL i,obs (y,d) × P mult i,j (d). (4) Equations (1-4) reveal that a key issue in the delta change approach is the estimation of the climatological annual cycle in a predefined period.In fact, the delta change approach states nothing but how the climatological annual cycle changes in the transition of the climate state from CTL to SCE periods according to climate model simulations.If one fails to estimate a robust annual cycle in either the CTL or the SCE period, one will fail to project an accurate change of the annual cycle.

Estimation of the climatological annual cycle
It is not possible to derive the true climatological annual cycle of any variable but only an estimate thereof, due to the natural variability and the limited duration of observed or simulated data records.The uncertainty of the estimate might be represented by a stochastic component.Ideally, the estimated climatological annual cycle should be robust, not depend on the stochastic components in the time series, and preserve the amplitude of the annual cycle.Often, the optimisation of these criteria is a trade-off, and it is not trivial to choose an optimal method to estimate the climatological annual cycle.
In this study, we used MAs and a spectral approach as an alternative to the MA for the estimation of the climatological annual cycle of T and P .In the MA approach, the terms X(d) CTL mod and X(d) SCE mod in Eqs.
(1) and (2) become where y s and y e denote the start and end year of the chosen period and n stands for the number of days before and after the day d in each year y.days.The larger the n, the smaller the effect of the natural variability on the estimate of the climatological annual cycle.However, the amplitude of the annual cycle is more strongly damped for larger n.
In the spectral approach, we investigated a spectral reconstruction of the climatological annual cycle by a superposition of harmonics with the base period P 0 = 365 days as ) The superscript k indicates the order of the harmonic components and H is the maximum order retained.The coefficients a k i,j and b k i,j are estimated using the discrete Fourier transform (see e.g.von Storch and Zwiers, 1999) from the daily time series of the GCM-RCM j at station site i in the CTL and SCE period.
For GCM-RCMs using the Gregorian calendar, the base period P 0 is set to 365.25 days to account for leap years (Narapusetty et al., 2009).The HadCM3Q0 and HadCM3Q3 driven RCMs have a 360 days calendar.For these GCM-RCMs, we set P 0 to 360 days.Having estimated the coefficients of the harmonic smoothing model at each station site and for each GCM-RCM, we scale the different lengths of the annual cycle to fit 365 days by choosing P 0 = 365 days in the reconstruction of X i,j (d) as in Eq. ( 6).
In the spectral framework of harmonics, the choice of the maximum order H is the only free parameter.The larger H is, the more the details of the annual cycle can be resolved, but the more vulnerable the spectral model becomes to influences of natural variability and overfitting.

Analysis of synthetically generated precipitation time series
Let's assume we could sample two 30 yr long precipitation time series from a stationary climate and derive the annual cycle of the precipitation change between the two time series.Stationary here means that the mean climate state is the same in both samples, but the two realisations are modulated by natural variability.Since we know that the climate is stationary by assumption, the asymptotic solution of the precipitation change (expressed as a ratio) should equal one representing no precipitation change.Any deviation from one, e.g. the occurrence of an annual cycle, is solely caused by sampling variability, and does not contain any climate signal.
Here, we investigate the sampling artefacts in the annual cycle of the precipitation changes in an idealised stochastic setup using a rainfall generator with stationary parameters.The setup consists of four experiments which are representative for precipitation series of individual station sites (STATION) and regions (REGION) both for the northern (CHN) and southern (CHS) parts of Switzerland, in order to study two distinct climates at station and regional scale.Following Wilks and Wilby (1999), we employed a first order Markov chain rainfall generator with the precipitation intensity being modelled by a two-parameter gamma distribution.Let p dw and p ww be the transition probabilities from a dry to wet and a wet to wet day, respectively.Given realisations of the uniform random number r 1 on the unit interval, the precipitation occurence Y (t) is modelled as where 1 stands for a wet and 0 for a dry day.On wet days, the precipitation intensity I (t) is sampled from a gamma probability density function according to where α and β are parameters of the gamma distribution and is the gamma function.In this case, the synthetic precipitation time series P synth is derived as For the single-station experiments CHN STATION and CHS STATION , we derived the parameters p dw , p ww , α and β from the observed daily precipitation records in the period 1980-2009 at the stations Bern (BER) and Lugano (LUG) respectively.The parameters were estimated for each season separately.At the transition from one season to the other, the parameter set is changed but the wet/dry state from the last day of the previous season is taken for the continuation of the Markov chain.
For the experiments representing regional precipitation, we first selected all stations in a radius r around BER (r = 60 km) and LUG (r = 40 km) that have less than 10 % missing values in the period 1980-2009 and calculated the mean daily precipitation time series therefrom.Remaining missing data were ignored in the averaging.The selected stations are indicated by red and blue dots in Fig. 1.We are aware that this averaging does not follow any spatial interpolation standards.The procedure suffices, though, to analyse the effect of spatial averaging on the fluctuations of the climate change signal.Table 2 lists the seasonal parameter settings for each of the four experiments.
For each experiment, we generated 100 realisations of a daily precipitation time series with a length of 30 yr.Subsequently, we randomly chose 500 pairs out of the 100 realisations and calculated the multiplicative precipitation change signal by MAs with window widths of 15, 31, 61 and 91 days.the 91d MA.Thus, the effects of sampling variability presented below are transferable to these common averaging intervals as well.

Synth ∆P [-]
The results are shown in Fig. 2. Since both time series of each pair have been generated with the same rainfall generator settings, the asymptotic solution is one, indicating no change.The dots in Fig. 2 display one randomly chosen realisation of the synthetically generated climate change signal P synth using different MA window widths.The 15d MA estimate shows large fluctuations in every experiment.The wider the MA window becomes, the smaller the fluctuations get.The 31d MA, corresponding to a monthly resolution, is a standard averaging window length in many impact studies.In the 31d MA estimates, the amplitudes of the P synth fluctuations are typically in the order of 0.2, but spikes can be as large as 1.3 or 0.7 as in the case of CHS STATION .The grey bands depict the 10th-90th % quantile range of the 500 realisations.
Comparison of the upper and lower panels in Fig. 2 shows that the spatial averaging does not reduce the band width of the 10th-90th % quantile range substantially.In the CHN experiments, spatial averaging reduces the width of the 31d MA band averaged over the year from 0.87-1.15 to 0.89-1.12whereas in the CHS cases, the width is reduced from 0.82-1.22 to 0.84-1.19.
The results indicate that on station scale as well as on regional scales relevant for hydrological impact studies, the fluctuations of P arising from sampling variability alone have to be considered in interpretations of climate change signals.The exact range of the sampling variability is dependent on the averaging window width, the spatial scale, the region of interest and the length of the climate records.For 31d MAs, our analysis shows for representative climate regions of Switzerland, that P values in the range of 0.8 to 1.2 could be solely caused by sampling variability and do not necessarily contain a climate change signal.Furthermore,  1).Results of five 30 year periods are shown in different colours.The MSE CV have been normalised by the mean MSE CV for display reasons.In the case of T , the MSE CV of HO 0 is much larger compared to higher order harmonic models and is therefore not shown.
these fluctuations.However, the maximal order of the harmonic smoothing model (see Eq. ( 6)) needs to be chosen.In this section we first define the optimal order of the harmonic smoothing model for T and P .We then present a qualitative comparison between the harmonic and the 31d MA estimates of the climate change signal at station sites, since monthly averaging periods are often employed in climate impact studies.We also discuss the limitations of the spectral smoothing methodology in detail.Finally, we show the climate change signals of T and P estimated by harmonic smoothing for 10 GCM-RCM chains at station sites in Switzerland.

Estimation of the optimal harmonic model
We use long-term observational station records of the Swiss climate monitoring network to constrain the maximum order of the harmonic smoothing model which is then applied to GCM-RCM series.This approach implicitly assumes that signal components from GCM-RCM time series having a higher frequency than the optimal harmonic order are considered as being noise.We use a cross-validation technique to specify the harmonic order for the annual cycle that optimally represents the time series.The methodology is described in detail in Narapusetty et al. (2009).Here, we give only a brief introduction and present our specific setup.We extracted 30 year time slices from 25 temperature and 26 precipitation station records with daily resolution, and split them into ten blocks of three year lengths.Five different 30 year time windows are analysed in order to test the robustness of the results with respect to decadal variability.At each station and for each order of the harmonic smoothing model, we carry out a 10-fold cross-validation by calibrating the harmonic model on 9 of the 10 blocks and validating it on the remaining block.The goal of the cross-validation is to estimate the harmonic model that has the lowest estimated prediction error (EPE).The EPE is a measure of the model error in an independent data set that was not used for calibration.It therefore penalises models that are overfitting the data (see, e. g., Chapter 7 in Hastie et al., 2009).We use the mean squared error (MSE) as a measure for the EPE and call it the cross-validated MSE (MSE CV ).The MSE CV is optimal for normally distributed and independent residuals.P time series however show strongly non-normal residual distributions.This might cause the estimation of the optimal model to be biased.We therefore carried out the cross-validation on the original data, on root-transformed and on Box-Cox transformed data (Wilks, 2006, pp. 43) in order to test how sensitive the results are with respect to the distributional characteristics of the P data.Since the Box-Cox transformation only works on positive definite variables, we replaced zeros in the P data by 0.0001.
Figure 3 shows the results of the cross-validation.For T , the MSE CV drops to a low level at the harmonic order (HO) of 2 and remains on this low level up to HO 8. Within this plateau, the differences between the models in terms of the MSE CV are small.Depending on the analysis period, the order with the lowest MSE CV varies between HO 2 and HO 8.
The MSE CV starts to consistently increase again for higher orders than HO 8.For P , the MSE CV has a minimum at HO 2, but the difference to HO 3 is very small.This result is robust for different analysis periods.We only show the results of the Box-Cox transformed data.The results based on the original and root-transformed data are very similar so are not shown here.
The above analysis yields different optimal harmonic orders for T and P .However, as the two atmospheric variables are linked through dynamical and thermodynamical processes, the optimal order for both variables should preferably be the same.We thus chose HO 3 as the optimal order for T and P .With a higher joint HO, we would accept higher-frequency precipitation fluctuations that could stem from natural variability rather than climate change.

Comparison between the moving average and the spectral estimation of the climate change signal
Based on the results in section 5.1, we use a third order harmonic model (HO 3) to estimate the annual cycle of T and P in the CTL and SCE periods and compare it to 31d MA estimates.We expect the HO 3 estimates to be characterised by smoother annual cycles and smaller peaks in the annual In the case of T , the MSE CV of HO 0 is much larger compared to higher order harmonic models and is therefore not shown.
the spikes within the annual cycles of P call for estimation methods that produce smoother climatological annual cycles than MAs.

Analysis of the climate change signal from regional climate models at Swiss station sites
The stochastic analysis in Sect. 4 revealed that for variables having similar characteristics as P , like e.g. a clustering of events and a heavily skewed intensity distribution, estimates of the climate change signal using MAs are prone to substantial artificial fluctuations caused by natural variability.In particular, such fluctuations lead to an impaired representation of the minima and maxima in the annual cycle of the climate change signal.Harmonic smoothing is able to filter these fluctuations.However, the maximal order of the harmonic smoothing model (see Eq. 6) needs to be chosen.In this section we first define the optimal order of the harmonic smoothing model for T and P .We then present a qualitative comparison between the harmonic and the 31d MA estimates of the climate change signal at station sites, since monthly averaging periods are often employed in climate impact studies.We also discuss the limitations of the spectral smoothing methodology in detail.Finally, we show the climate change signals of T and P estimated by harmonic smoothing for 10 GCM-RCM chains at station sites in Switzerland.

Estimation of the optimal harmonic model
We use long-term observational station records of the Swiss climate monitoring network to constrain the maximum order of the harmonic smoothing model which is then applied to GCM-RCM series.This approach implicitly assumes that signal components from GCM-RCM time series having a higher frequency than the optimal harmonic order are considered as being noise.We use a cross-validation technique to specify the harmonic order for the annual cycle that optimally represents the time series.The methodology is described in detail in Narapusetty et al. (2009).Here, we give only a brief introduction and present our specific setup.We extracted 30 yr time slices from 25 temperature and 26 precipitation station records with daily resolution, and split them into ten blocks of three year lengths.Five different 30 yr time windows are analysed in order to test the robustness of the results with respect to decadal variability.At each station and for each order of the harmonic smoothing model, we carry out a 10-fold cross-validation by calibrating the harmonic model on 9 of the 10 blocks and validating it on the remaining block.The goal of the cross-validation is to estimate the harmonic model that has the lowest estimated prediction error (EPE).The EPE is a measure of the model error in an independent data set that was not used for calibration.It therefore penalises models that are overfitting the data (see, e.g.Chapter 7 in Hastie et al., 2009).We use the mean squared error (MSE) as a measure for the EPE and call it the cross-validated MSE (MSE CV ).The MSE CV is optimal for normally distributed and independent residuals.P time series however show strongly non-normal residual distributions.This might cause the estimation of the optimal model to be biased.We therefore carried out the cross-validation on the original data, on root-transformed and on Box-Cox transformed data (Wilks, 2006, pp. 43) in order to test how sensitive the results are with respect to the distributional characteristics of the P data.Since the Box-Cox transformation only works on positive definite variables, we replaced zeros in the P data by 0.0001.
Figure 3 shows the results of the cross-validation.For T , the MSE CV drops to a low level at the harmonic order (HO) of 2 and remains on this low level up to HO 8. Within this plateau, the differences between the models in terms of the MSE CV are small.Depending on the analysis period, the order with the lowest MSE CV varies between HO 2 and HO 8.The MSE CV starts to consistently increase again for higher orders than HO 8.For P , the MSE CV has a minimum at HO 2, but the difference to HO 3 is very small.This result is robust for different analysis periods.We only show the results of the Box-Cox transformed data.The results based on the original and root-transformed data are very similar so are not shown here.cycle than 31d MA estimates.For illustration, Fig. 4 displays annual cycles of T and P at the two station sites BER and LUG as modelled by ETHZ-HadCM3Q0-CLM in the CTL period and the SCE period 2070-2099 (upper panels) as well as the climate change signal (lower panels).These two stations and the selected model chain represent typical results.
In the case of T , the annual cycle in the CTL period is well captured at both stations, biases of up to 2 K arise for individual months.The fluctuations in the 31d MA estimate of ∆T have a time scale of typically one month.The amplitudes of these fluctuations are in the order of 0.5-1 K.
The HO 3 estimate treats these fluctuations as noise and results in a smooth annual cycle of ∆T .
In the CTL period, the depicted precipitation shows a large bias in winter on the northern side of the Alps (BER), whereas in southern Switzerland (LUG), the GCM-RCM is able to reproduce the two precipitation peaks in the annual cycle but has a biased amplitude.Such biases are not uncommon in regions of complex topography.For a detailed evaluation of the ENSEMBLES GCM-RCMs, we refer to Klein Tank et al. ( 2009) and references therein.The estimates of the climatological annual cycle using a 31d MA show high  frequency fluctuations in the CTL and SCE periods, which are amplified in the annual cycle of ∆P due to the division of SCE by CTL values.A spurious amplification can be seen at the station LUG in mid October, when a decrease of P in the CTL period and a rapid increase in the SCE period occur, leading to a spike in ∆P .The HO 3 estimate is not influenced by such high frequency fluctuations and results in a smooth annual cycle.
Figure 5 shows further examples of strong fluctuations in the 31d MA around the spectrally smoothed annual cycle ∆P at different station sites and for various model chains.The fluctuations of the 31d MA estimates relative to the spectrally smoothed annual cycles are in the same order of magnitude as the climate change signal.
In climate studies, an important figure is the ensemble mean of the climate change signal.Due to the averaging effect, it is expected that the fluctuations in the climate change signal of the ensemble mean are smaller than for the individual GCM-RCM.Fig. 6 shows the ensemble mean of the ∆P signals of the individual GCM-RCMs as estimated by the spectral method and by a 31d MA.The fluctuations of the 31d MA estimate around the HO 3 estimate are still substantial and are often in the range of the climate change signal itself.Thus, smoothing might also be necessary for the ensemble mean climate change signal.For applications to impact models such as, for e. g., hydrological models, we recommend not to use the ensemble mean climate change signal, but rather the individual GCM-RCMs ∆P and to derive the ensemble mean at the end of the entire impact modelling chain.Impact models are usually non-linear and thus do not yield the same results whether the averaging over the ensemble is done at the GCM-RCM stage or at the end of the impact modelling chain.

Limitations of the methodology
The employed harmonic smoothing model uses a sharp spectral low-pass filter as it removes all harmonic components above the order of 3 and retains the orders 0 to 3 without any The above analysis yields different optimal harmonic orders for T and P .However, as the two atmospheric variables are linked through dynamical and thermodynamical processes, the optimal order for both variables should preferably be the same.We thus chose HO 3 as the optimal order for T and P .With a higher joint HO, we would accept higher-frequency precipitation fluctuations that could stem from natural variability rather than climate change.

Comparison between the moving average and the spectral estimation of the climate change signal
Based on the results in Sect.5.1, we use a third order harmonic model (HO 3) to estimate the annual cycle of T and P in the CTL and SCE periods and compare it to 31d MA estimates.We expect the HO 3 estimates to be characterised by smoother annual cycles and smaller peaks in the annual cycle than 31d MA estimates.For illustration, Fig. 4 displays annual cycles of T and P at the two station sites BER and LUG as modelled by ETHZ-HadCM3Q0-CLM in the CTL period and the SCE period 2070-2099 (upper panels) as well as the climate change signal (lower panels).These two stations and the selected model chain represent typical results.cycle than 31d MA estimates.For illustration, Fig. 4 displays annual cycles of T and P at the two station sites BER and LUG as modelled by ETHZ-HadCM3Q0-CLM in the CTL period and the SCE period 2070-2099 (upper panels) as well as the climate change signal (lower panels).These two stations and the selected model chain represent typical results.
In the case of T , the annual cycle in the CTL period is well captured at both stations, although biases of up to 2 K arise for individual months.The fluctuations in the 31d MA estimate of ∆T have a time scale of typically one month.The amplitudes of these fluctuations are in the order of 0.5-1 K.
The HO 3 estimate treats these fluctuations as noise and results in a smooth annual cycle of ∆T .
In the CTL period, the depicted precipitation shows a large bias in winter on the northern side of the Alps (BER), whereas in southern Switzerland (LUG), the GCM-RCM is able to reproduce the two precipitation peaks in the annual cycle but has a biased amplitude.Such biases are not uncommon in regions of complex topography.For a detailed evaluation of the ENSEMBLES GCM-RCMs, we refer to Klein Tank et al. ( 2009) and references therein.The estimates of the climatological annual cycle using a 31d MA show high  Thus, smoothing might also be necessary for the ensemble mean climate change signal.For applications to impact models such as, for e. g., hydrological models, we recommend not to use the ensemble mean climate change signal, but rather the individual GCM-RCMs ∆P and to derive the ensemble mean at the end of the entire impact modelling chain.Impact models are usually non-linear and thus do not yield the same results whether the averaging over the ensemble is done at the GCM-RCM stage or at the end of the impact modelling chain.

Limitations of the methodology
The employed harmonic smoothing model uses a sharp spectral low-pass filter as it removes all harmonic components above the order of 3 and retains the orders 0 to 3 without any In the case of T , the annual cycle in the CTL period is well captured at both stations, although biases of up to 2 K arise for individual months.The fluctuations in the 31d MA estimate of T have a time scale of typically one month.The amplitudes of these fluctuations are in the order of 0.5-1 K.The HO 3 estimate treats these fluctuations as noise and results in a smooth annual cycle of T .
In the CTL period, the depicted precipitation shows a large bias in winter on the northern side of the Alps (BER), whereas in southern Switzerland (LUG), the GCM-RCM is able to reproduce the two precipitation peaks in the annual cycle but has a biased amplitude.Such biases are not uncommon in regions of complex topography.For a detailed evaluation of the ENSEMBLES GCM-RCMs, we refer to Klein Tank et al. ( 2009) and references therein.The estimates of the climatological annual cycle using a 31d MA show high frequency fluctuations in the CTL and SCE periods, which are amplified in the annual cycle of P due to the division of SCE by CTL values.A spurious amplification can be seen at the station LUG in mid October, when a decrease of P in the CTL period and a rapid increase in the SCE period occur, leading to a spike in P .The HO 3 estimate is not influenced by such high frequency fluctuations and results in a smooth annual cycle.
Figure 5 shows further examples of strong fluctuations in the 31d MA around the spectrally smoothed annual cycle P at different station sites and for various model chains.The fluctuations of the 31d MA estimates relative to the spectrally smoothed annual cycles are in the same order of magnitude as the climate change signal.
In climate studies, an important figure is the ensemble mean of the climate change signal.Due to the averaging effect, it is expected that the fluctuations in the climate change signal of the ensemble mean are smaller than for an individual GCM-RCM.Figure 6  damping.Such sharp spectral filters are prone to overshootings which can occur in situations when sudden changes in a time series take place within a time scale that cannot be resolved by the spectral model.This is also known as the Gibbs phenomenon.Overshootings are critical as, e. g., negative precipitation values can occur so results need to be scanned for such overshootings.For the application to the delta change methodology, there are three types of overshooting cases: overshootings occur in the smoothed mean annual cycle of (a) observed time series, (b) climate model time series in the CTL period (and possibly also in the SCE period) and (c) climate model time series in the SCE period only.Case a) indicates that the spectral smoothing model is not appropriate for the climate in the region of interest.Case (b) indicates that the climate model is spectrally biased since its time series contain more spectral frequencies than estimated from the observed time series.In case (c), the climate model is not spectrally biased but changes its spectral characteristics substantially in the transition from the CTL to the SCE period.In case (a), one should not use the spectral model but resort to other smoothing methods such as, e. g., generalized linear models or smoothing filters that do not produce overshootings.In cases (b) and (c), one could use a less sharp spectral filter to estimate the mean annual cycle.
Figure 7 shows the application of the sharp (i. e., having a sharp cutoff) and a less sharp (i. e., having a gradual cutoff) HO 3 filter as well as the 31d MA to the time series of the DMI-ARPEGE-HIRHAM GCM-RCM at the station LUG.The inset in Fig. 7 shows the response function of the 31d MA and the two spectral filters constructed following Duchon (1979).The response function of the 31d MA is a damped sinus oscillation (black line).The sharp filter has a cutoff at HO 3 (red line).The gradual cutoff filter (green line) has a linear transition from the response value of 1 to 0. We subjectively chose a transition range of ± 2 harmonic orders around HO 3.
In the example shown, the GCM-RCM has a pronounced summer drying in the CTL period.In the SCE period, the summer drying still exists, but precipitation starts to increase a bit earlier in the year.As a result, the 31d MA shows a very high and unrealistic multiplicative climate change signal of up to 4.5 in August.The sharp cutoff HO 3 filter has problems resolving the pronounced low precipitation period in summer and produces overshootings that cause negative precipitation.The climate change signal is therefore not meaningful.The gradual cutoff HO 3 filter does not produce overshootings and results in a reasonable climate change signal.It remains up to the users of the presented methodology to choose whether GCM-RCMs that show problems such as those in cases (b) and (c) should be included in the ensemble or not.We chose not to use GCM-RCMs that have severe overshootings.In principle though, the approach of a spectral filter with a gradual cutoff would be a suitable workaround if such GCM-RCMs should be included.For brevity, we show results of the climate signal's annual cycle only at the two exemplary stations BER and LUG (see Fig. 1). Figure 8 shows each GCM-RCMs annual cycle of ∆T and ∆P respectively for both SCE periods relative to the CTL period 1980-2009.To compare the changes to the natural variability range, we resampled each station's observed precipitation record of the CTL period.
We constructed 100 realisations with a length of 30 years by resampling with replacement the years of the observed record.From the 100 realisations, we randomly chose 500 pairs and estimated the climate change signal between the pairs.The range of ± 1 standard deviation (σ) of the 500 resampled realisations is shown as a grey band in Fig. 8.  signal itself.Thus, smoothing might also be necessary for the ensemble mean climate change signal.For applications to impact models such as, for e.g.hydrological models, we recommend not to use the ensemble mean but rather the individual GCM-RCMs climate change signal, and to derive the ensemble mean at the end of the entire impact modelling chain.Impact models are usually non-linear and thus do not yield the same results whether the averaging over the ensemble is done at the GCM-RCM stage or at the end of the impact modelling chain.

Limitations of the methodology
The employed harmonic smoothing model uses a sharp spectral low-pass filter as it removes all harmonic components above the order of 3 and retains the orders 0 to 3 without any damping.Such sharp spectral filters are prone to overshootings which can occur in situations when sudden changes in a time series take place within a time scale that cannot be resolved by the spectral model.This is also known as the Gibbs phenomenon.Overshootings are problematic as, e.g.negative precipitation values can occur, so results need to be scanned for such overshootings.For the application to the delta change methodology, there are three types of overshooting cases: overshootings occur in the smoothed mean annual cycle of (a) observed time series, (b) climate model time series in the CTL period (and possibly also in the SCE period) and (c) climate model time series in the SCE period only.Case (a) indicates that the spectral smoothing model is not appropriate for the climate in the region of interest.Case (b) indicates that the climate model is spectrally biased since its time series contain more spectral frequencies than estimated from the observed time series.In case (c), the climate model is not spectrally biased but changes its spectral characteristics substantially in the transition from the CTL to the SCE period.In case (a), one should not use the spectral model but resort to other smoothing methods such as, e.g.generalized linear models or smoothing filters that do not produce overshootings.In cases (b) and (c), one could use a less sharp spectral filter to estimate the mean annual cycle.
Figure 7 shows the application of the sharp (i.e.having a sharp cutoff) and a less sharp (i.e.having a gradual cutoff) HO 3 filter as well as the 31d MA to the time series of the DMI-ARPEGE-HIRHAM GCM-RCM at the station indicates that the climate model is spectrally biased since its time series contain more spectral frequencies than estimated from the observed time series.In case (c), the climate model is not spectrally biased but changes its spectral characteristics substantially in the transition from the CTL to the SCE period.In case (a), one should not use the spectral model but resort to other smoothing methods such as, e. g., generalized linear models or smoothing filters that do not produce overshootings.In cases (b) and (c), one could use a less sharp spectral filter to estimate the mean annual cycle.Figure 7 shows the application of the sharp (i. e., having a sharp cutoff) and a less sharp (i. e., having a gradual cutoff) HO 3 filter as well as the 31d MA to the time series of the DMI-ARPEGE-HIRHAM GCM-RCM at the station LUG.The inset in Fig. 7 shows the response function of the 31d MA and the two spectral filters constructed following Duchon (1979).The response function of the 31d MA is a damped sinus oscillation (black line).The sharp filter has a cutoff at HO 3 (red line).The gradual cutoff filter (green line) has a linear transition from the response value of 1 to 0. We subjectively chose a transition range of ± 2 harmonic orders around HO 3. In the example shown, the GCM-RCM has a pronounced summer drying in the CTL period.In the SCE period, the summer drying still exists, but precipitation starts to increase a bit earlier in the year.As a result, the 31d MA shows a very high and unrealistic multiplicative climate change signal of up to 4.5 in August.The sharp cutoff HO 3 filter has problems resolving the pronounced low precipitation period in pairs and estimated the climate change signal between the pairs.The range of ± 1 standard deviation (σ) of the 500 resampled realisations is shown as a grey band in Fig. 8.  LUG.The inset in Fig. 7 shows the response function of the 31d MA and the two spectral filters constructed following Duchon (1979).The response function of the 31d MA is a damped sinus oscillation (black line).The sharp filter has a cutoff at HO 3 (red line).The gradual cutoff filter (green line) has a linear transition from the response value of 1 to 0. We subjectively chose a transition range of ±2 harmonic orders around HO 3.
In the example shown, the GCM-RCM has a pronounced summer drying in the CTL period.In the SCE period, the summer drying still exists, but precipitation starts to increase a bit earlier in the year.As a result, the 31d MA shows a very high and unrealistic multiplicative climate change signal of up to 4.5 in August.The sharp cutoff HO 3 filter has problems resolving the pronounced low precipitation period in summer and produces overshootings that cause negative precipitation.The climate change signal is therefore not meaningful.The gradual cutoff HO 3 filter does not produce overshootings and results in a reasonable climate change signal.
It remains up to the users of the presented methodology to choose whether GCM-RCMs that show problems such as those in cases (b) and (c) should be included in the ensemble or not.We chose not to use GCM-RCMs that have severe overshootings.In principle though, the approach of a spectral filter with a gradual cutoff would be a suitable workaround if such GCM-RCMs should be included.

Annual cycle of the climate change signal
For brevity, we show results of the climate change signal's annual cycle only at the two exemplary stations BER and LUG (see Fig. 1). Figure 8 shows each GCM-RCMs annual cycle of T and P respectively for both SCE periods relative to the CTL period 1980-2009.To compare the changes to the natural variability range, we resampled each station's observed precipitation record of the CTL period.We constructed 100 realisations with a length of 30 yr by resampling with replacement the years of the observed record.From the 100 realisations, we randomly chose 500 pairs and estimated the climate change signal between the pairs.The range of ±1 standard deviation (σ ) of the 500 resampled realisations is shown as a grey band in Fig. 8.In the case of T , the ensemble mean shows peaks in winter and summer in both SCE periods.The model spread is largest in summer, which is mainly due to a strong summer warming of HadCM3Q0-driven experiments.Generally, the T signal for both SCE periods is distinctively above the estimated natural variability range.The natural variability range of P relative to the projected P values is much larger than in the case of T .Also, the range of natural variability strongly differs from station to station.Only the decrease of P in the summer of the later scenario period as projected by a majority of the GCM-RCMs substantially exceeds the range of the natural variability.
The results also indicate that for T and to a lesser extent for P , GCM-RCMs belonging to the same GCM family show similar patterns in the annual cycle of the climate change signal.

Spatial patterns of seasonal mean changes
Figure 9 shows the mean seasonal pattern of T and P for both scenario periods.Only the results for the seasons DJF and JJA are shown since the analysis of the annual cycles showed these seasons to have stronger climate change signals than the transition seasons.
For T , the spatial pattern is homogenous across Switzerland with the exception of the Alpine ridge region in JJA that generally shows higher T than other regions in Switzerland.The strongest warming is projected for JJA.The median of all station's ensemble mean T for JJA is 1.35 K for 2021-2050 and 3.72 K for 2070-2099.At most stations, the ensemble mean T is larger than 2 times the standard deviation of the natural variability for both scenario periods.
The strongest seasonal P signal is projected for JJA in the SCE period 2070-2099 with ensemble mean P values around 0.8 in large parts of Switzerland.Southern Switzerland is projected to have the strongest decrease of summer precipitation, with P values around 0.7.For DJF, an increase of P can be expected but the strength of the signal is smaller than for JJA.In the period 2021-2050, the P values generally do not exceed the range of estimated natural variability.
www.hydrol-earth-syst-sci.net/15/2777/2011/ Hydrol.Earth Syst.Sci., 15, 2777-2788, 2011 The ensemble mean's projected seasonal changes for the SCE period 2070-2099 are consistent with the results from the PRUDENCE project (Christensen et al., 2007).In the PRUDENCE project, the ensemble mean T in the Alpine region for the SCE period 2070-2100 relative to the CTL period 1961-1990 was +2 K for winter and +4 K for summer.The estimated ensemble mean P was +10 % and −30 % for winter and summer respectively (Christensen and Christensen, 2007).

Summary and conclusions
The delta change method commonly used in climate impact modelling studies requires a representation of the climate change signal's annual cycle.This implies the estimation of the annual cycle of T and P both in the CTL and the SCE period.Using a stochastic rainfall generator, we showed that climate change signals of mean precipitation derived by MAs are strongly affected by sampling artefacts.Spatial aggregation to a region corresponding to the area of a few RCM grid cells does not reduce the effect of sampling variability on the climate change signal substantially.
Climate change signals estimated using narrow averaging intervals, such as monthly means, should thus be regarded with caution, since associated artificial peaks in the annual cycle can lead to undesirable effects when used in combination with non-linear impact models.
We used a spectral smoothing to ameliorate the effects of natural variability on artificial fluctuations in the annual cycle.Compared to 31d MA estimates, the spectral smoothing successfully filters intraannual fluctuations.In a few cases when a strong amplitude of the annual precipitation cycle is paired with a large relative precipitation change, the spectral smoothing produces overshootings.
The derived climate change signal for the ENSEMBLES GCM-RCM chains is particularly clear for the later SCE period 2070-2099.The peak in the ensembles mean's T is around 4 K.In the case of P , a pronounced decrease of summer precipitation is projected for the whole of Switzerland.In the other seasons, precipitation is projected to increase.
In this study we focussed on changes in the annual cycle of mean T and P as used in the delta change method.This is a statistical model of low complexity in the whole variety of statistical post-processing methods.More complex models, such as, e.g.bias-correction methods involving a mean and variance scaling or quantile-mapping, are possibly even more sensitive to natural variability than the delta change method.The quantification of the sensitivity, however, requires further study.The study presented here indicates that the representation of the annual cycle in any statistical post-processing or downscaling method should be addressed with care.
used daily near-surface T and P data from 10 GCM-RCM model chains provided by the ENSEMBLES project (van der Linden and Mitchell, 2009) to estimate the annual cycle of the climate change signal (see Tab. 1).All model chains use the A1B emission scenario, cover the period 1951-2099 and have a horizontal resolution of about 25 km.From the whole set of model chains available through the ENSEMBLES project, we excluded the HadCM3Q16 driven model chains as well as DMI-ARPEGE-ALADIN in our analysis, due to pronounced summer dryings that caused severe overshooting in the spectral smoothing (see section 5.3 for details).The geographical focus of our study is Switzerland.Throughout this paper, the estimation of the climate change signal is based on GCM-RCM data interpolated to station

Fig. 1 .
Fig. 1.Map with station locations for T (top) and P (bottom).Stars indicate stations belonging to the long-term Swiss climate monitoring network.Blue and red dots show the selected stations for the CHN REGION (93 stations) and CHS REGION (23 stations) experiments, respectively (see Sect. 4).
Also, note that the commonly used mean monthly climate change signals are a subsample of the 31d MA.The same similarity exists between seasonal climate change signals and Hydrol.Earth Syst.Sci., 15, 2777-2788, 2011 www.hydrol-earth-syst-sci.net/15/2777/2011/ Jan Feb Mar Apr Mai Jun Jul Aug Sep Oct Nov Dec 0

Fig. 2 .
Fig. 2. Annual cycles of the multiplicative precipitation change signals (Synth P ) for the four experiments CHN STATION , CHS STATION , CHN REGION , CHS REGION estimated from time series pairs generated by a stationary first order Markov chain rainfall generator.The correct asymptotic solution is 1 representing no change.The dots indicate a specific realisation whereas the grey bands show the 10th-90th quantile range of 500 realisations.The precipitation change signals were estimated by moving averages (MA) based on averaging window widths of 15, 31, 61 and 91 days over 30 yr of daily data.

Fig. 3 .
Fig.3.Mean over all station's MSE CV of harmonic models with increasing harmonic order (HO 1 to HO 12) for observed daily T (top) and Box-Cox transformed P (bottom) series at Swiss climate monitoring stations (see Fig.1).Results of five 30 yr periods are shown in different colours.The MSE CV have been normalised by the mean MSE CV for display reasons.In the case of T , the MSE CV of HO 0 is much larger compared to higher order harmonic models and is therefore not shown.

Fig. 4 .
Fig. 4. Comparison between annual cycles estimated by a 31d MA (dashed lines) and a third order harmonic smoothing (HO 3; solid lines).Annual cycles of T and P are displayed in the left and right panels respectively.Results of the model chain ETHZ-HadCM3Q0-CLM at the two exemplary stations BER (top) and LUG (bottom) are shown.In each of the four blocks, the top panels contain the annual cycles in the CTL period 1980-2009 and SCE period 2070-2099.The annual cycles of the observed records in the CTL period are added in grey.In the bottom panels, the delta change signals are shown.

Fig. 5 .
Fig. 5. Examples of ∆P as estimated by a 31d MA (dashed lines) and the spectral smoothing (solid lines) of various model chains at different station sites as indicated in each panel's title.Examples are shown for the SCE period 2070-2099.

Fig. 4 .
Fig. 4. Comparison between annual cycles estimated by a 31d MA (dashed lines) and a third order harmonic smoothing (HO 3; solid lines).Annual cycles of T and P are displayed in the left and right panels respectively.Results of the model chain ETHZ-HadCM3Q0-CLM at the two exemplary stations BER (top) and LUG (bottom) are shown.In each of the four blocks, the top panels contain the annual cycles in the CTL period 1980-2009 and SCE period 2070-2099.The annual cycles of the observed records in the CTL period are added in grey.In the bottom panels, the delta change signals are shown.

TFig. 4 .
Fig. 4. Comparison between annual cycles estimated by a 31d MA (dashed lines) and a third order harmonic smoothing (HO 3; solid lines).Annual cycles of T and P are displayed in the left and right panels respectively.Results of the model chain ETHZ-HadCM3Q0-CLM at the two exemplary stations BER (top) and LUG (bottom) are shown.In each of the four blocks, the top panels contain the annual cycles in the CTL period 1980-2009 and SCE period 2070-2099.The annual cycles of the observed records in the CTL period are added in grey.In the bottom panels, the delta change signals are shown.

Fig. 5 .
Fig. 5. Examples of ∆P as estimated by a 31d MA (dashed lines) and the spectral smoothing (solid lines) of various model chains at different station sites as indicated in each panel's title.Examples are shown for the SCE period 2070-2099.

Fig. 5 .
Fig. 5. Examples of P as estimated by a 31d MA (dashed lines) and the spectral smoothing (solid lines) of various model chains at different station sites as indicated in each panel's title.Examples are shown for the SCE period 2070-2099.

Fig. 6 .
Fig. 6.Examples of the ensemble mean ∆P as estimated by a 31d MA (dashed lines) and the HO 3 spectral smoothing (solid lines) at the station sites BER and LUG.Examples are shown for the SCE period 2070-2099.

5. 4
Climate change signal at Swiss station sites 5.4.1 Annual cycle of the climate change signal

Fig. 7 .
Fig. 7. Illustration of the overshooting problem for the example of the GCM-RCM DMI-ARPEGE-HIRHAM at the station LUG.Annual cycle of P in the CTL (top) and SCE 2070-2099 (middle) period as well as its climate change signal (bottom) estimated by use of different filters.The inset shows the response function of the spectral filters.

Fig. 6 .
Fig. 6.Examples of the ensemble mean P as estimated by a 31d MA (dashed lines) and the HO 3 spectral smoothing (solid lines) at the station sites BER and LUG.Examples are shown for the SCE period 2070-2099.

Fig. 7 .
Fig. 7. Illustration of the overshooting problem for the example of the GCM-RCM DMI-ARPEGE-HIRHAM at the station LUG.Annual cycle of P in the CTL (top) and SCE 2070-2099 (middle) period as well as its climate change signal (bottom) estimated by use of different filters.The inset shows the response function of the spectral filters.

Fig. 7 .
Fig. 7. Illustration of the overshooting problem for the example of the GCM-RCM DMI-ARPEGE-HIRHAM at the station LUG.Annual cycle of P in the CTL (top) and SCE 2070-2099 (middle) period as well as its climate change signal (bottom) as estimated by use of different filters.The inset shows the response function of the filters.

Fig. 8 .
Fig. 8. Annual cycle of T (left panel group) and P (right panel group) for the scenario period 2021-2050 (top) and 2070-2099 (bottom) at the two exemplary stations BER and LUG.Each of the 10 GCM-RCMs is shown with an individual colour.The ensemble mean is indicated by a black line.The grey band shows the range of the estimated natural variability as ±σ of the resampled T s or P s at the station site.Note the different scales in upper and lower panels.

Table 1 .
List of the employed climate model chains from the ENSEMBLES project.
We used daily near-surface T and P data from 10 GCM-RCM model chains provided by the ENSEMBLES project (van derLinden and Mitchell, 2009)to estimate the annual cycle of the climate change signal (see Table1).All model chains use the A1B emission scenario, cover the period 1951-2099 and have a horizontal resolution of about 25 km.From the whole set of model chains available through the ENSEMBLES project, we excluded the HadCM3Q16 driven model chains as well as DMI-ARPEGE-ALADIN in our analysis, due to pronounced summer dryings that caused severe overshooting in the spectral smoothing (see Sect. 5.3 for details).

Table 1 .
List of the employed climate model chains from the EN-SEMBLES project.

Table 2 .
Seasonal parameter settings of the rainfall generator for the four experiments CHN STATION , CHN REGION , CHS STATION and CHS REGION .