A comparison of four techniques for separating different time scales in atmospheric variables

Atmospheric Environment 37 (2003) 313–325

A comparison of four techniques for separating different time scales in atmospheric variables Christian Hogrefea,*, Somaraju Vempatyb, S. Trivikrama Raoa, P. Steven Porterc a

Atmospheric Sciences Research Center, University at Albany, 251 Fuller Road, Albany, NY 12203, USA b Department of Mathematics, SRKR Engineering College, Bhimavaram, India c Department of Civil Engineering, University of Idaho, Idaho Falls, ID 83402, USA Received 9 April 2002; accepted 17 October 2002

Abstract In this paper, four methods for spectrally decomposing time series of atmospheric variables are compared. Two of these methods have been previously applied to the analysis of time series of atmospheric variables, while the others are being applied for the first time. This paper focuses on the practical applications of time scale separation techniques rather than on an in-depth comparison of the mathematical features of the filtering techniques. The performance of the above filtering methods is illustrated and evaluated using both simulated and observed ozone time series data. The adaptive window Fourier transform filter is shown to extract fluctuations of known frequency as cleanly as the Morlet wavelet and, therefore, is a useful new tool for time–frequency analyses of atmospheric variables. Simulation results indicate that all four of these filters provide qualitatively similar results when used to extract the energy in five frequency bands of particular interest in time series of atmospheric variables. However, differences can exist when different filters are used to study the temporal variations of the extracted components. In addition, it is shown that all filters are able to capture the year-to-year fluctuations in the magnitudes of individual components. Such analysis can be used to discern the time scales that cause long-term changes in pollutant concentrations. As no single filter performs best in separating the various time scales, great care has to be taken to match the filter characteristics with the objectives of a given analysis. r 2002 Elsevier Science Ltd. All rights reserved. Keywords: Time series analysis; Spectral decomposition; Wavelets; Ozone; Trends

1. Introduction Fourier analysis of time series of atmospheric variables reveals spectral energy over a wide range of scales, indicating the complex inter-play among different dynamical forcing mechanisms operating on different time scales. Such forcings include processes with a fixed frequency like the annual rotation of the earth around the sun or the daily rotation of the earth around its axis as well as processes that occur over a range of periods *Corresponding author. Bureau of Air Research, NYSDEC, 625 Broadway, Albany, NY 12233-3259, USA. Tel.: +1-518402-8402; fax: +1-518-402-9035. E-mail address: [email protected] (C. Hogrefe).

such as the evolution and decay of convective systems with typical periods of 30 min to several hours or the movement of synoptic-scale weather systems with characteristic time scales of 2–10 days. The strengths of these forcing mechanisms often vary in time, making the observed time series non-stationary, limiting the usefulness of traditional Fourier analysis. For a more complete study of such non-stationary time series, it is necessary to use time–frequency techniques that allow for the extraction of time series in certain frequency bands or removal of certain frequencies from time series data. In this study, we apply the Kolmogorov–Zurbenko filter (KZ) (Rao and Zurbenko, 1994), the adaptive window Fourier transform (AWFT) (Zurbenko and

1352-2310/02/$ - see front matter r 2002 Elsevier Science Ltd. All rights reserved. PII: S 1 3 5 2 - 2 3 1 0 ( 0 2 ) 0 0 8 9 7 - X

314

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

Porter, 1998), the wavelet transform (Torrence and Compo, 1998), and the elliptic filter (Rabiner and Gold, 1975; Poularika, 1998) to spectrally decompose simulated and observed ozone time series data. While the KZ and wavelet filters have been used extensively in the study of atmospheric variables (Rao et al., 1995; Milanchus et al., 1998; Sebald et al., 2000; Lau and Weng, 1995; Torrence and Compo, 1998; Yan et al., 2001), the elliptic filter and the AWFT are being applied to this problem for the first time. Thus, the object of this paper is to illustrate the application of a new tool for time–frequency analysis (i.e., AWFT), and to compare the relative performances of the four filtering techniques listed above in analyzing commonly encountered air quality problems. Specifically, we compare their ability to estimate the time series of fixed-frequency oscillations with potentially time-varying amplitude, to estimate the relative contribution of individual time scales to the overall process energy, to reproduce the temporal variations of simulated time series, and to extract longterm variations in the magnitudes of different temporal components. We also discuss important practical considerations such as computational flexibility and the ability of the different filters to work in a missing data environment. This paper focuses on the practical applications of time scale separation techniques rather than on an in-depth mathematical comparison. In Section 2, we introduce the scale separation methods considered here. A comparison of their separation characteristics and discussion of the results in the context of air quality management is presented in Section 3. The results are summarized in Section 4.

2. Database and time scale separation methods 2.1. Simulated data A simulated time series of known spectral composition is used to compare the performance of the different filters. The simulated time series consists of 20 years (175,200 h) of hourly data, containing sine waves at all 87,600 Fourier frequencies. The superposition of the sine waves at all Fourier frequencies creates a time series with the spectral properties of white noise and a standard deviation of 2.39. Two larger-amplitude sine waves corresponding to the diurnal and the annual cycle were added to this time series. The amplitude of the DU and annual cycle was chosen to be 1. 2.2. Observed data Time series of hourly ozone observations at Bayonne, NJ and Pasadena, CA were retrieved from the US EPA’s AIRS database.

2.3. KZ filter The KZ filter is a low-pass filter implemented via iteration of a simple moving average (Rao and Zurbenko, 1994; Eskridge et al., 1997). The length of the simple moving average, m; and the number of iterations, k; determine the shape of the energy transfer function (Eskridge et al., 1997). By taking the difference between two KZ filters with different parameters m and k (and thus different transmission characteristics), a bandpass filter is created. The KZ filter can be applied to time series containing missing data without the need for data restoration. Additional details on the KZ filter can be found in Rao et al. (1997), Hogrefe et al. (2000), and Porter et al. (2001).

2.4. Elliptic filter The elliptic filter is a recursive frequency domain filter. The advantage of the elliptic filters over other recursive filters such as the Butterworth or Chebyshev filter is that it provides given performance specifications with the lowest filter order. It, therefore, is often the filter of choice in electrical engineering applications (Poularika, 1998). Given an input time series X ðtÞ; the filtered time series Y ðtÞ can be written as

Y ðtÞ ¼

N X n¼0

an X ðt  nÞ þ

N X

bm Y ðt  mÞ:

ð1Þ

m¼1

The coefficients ai and bi are termed the recursion coefficients, and N is the order of the elliptic filter. The values of these variables depend on the desired spectral response properties of the elliptic filter and can be determined using a set of functions from the MATLAB signal processing toolbox (The Mathworks Inc., 2001). Specifically, once choices are made for the maximum allowable ripple in the pass band, the desired attenuation in the stop band relative to the pass band, and the desired cut-off frequency or frequencies (depending on whether the filter is employed in a low-pass, high-pass or band-pass mode), the order N of the elliptic filter can be determined by a call to the MATLAB function ellipord. Subsequently, the recursion coefficients ai and bi can be determined by a call to the MATLAB function ellip (The Mathworks Inc., 2001). For the analysis presented in this paper, a pass band ripple of 0.1 dB and an attenuation of the stop band of 40 dB below the peak of the pass band were chosen, resulting in a filter order of 4. Another feature of the elliptic filter is that the phase of the filtered signal is shifted relative to the original signal. This is illustrated in Section 3.3.

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

315

2.5. Wavelets

2.6. Adaptive window Fourier transform (AWFT)

Wavelet transform methods have found wide application in various fields of science including meteorology and oceanography (e.g. Lau and Weng, 1995; FoufoulaGeorgiou and Kumar, 1995; Torrence and Compo, 1998; Yan et al., 2001). In this paper, we follow closely Torrence and Compo (1998) for the application of wavelets. We used the Fortran wavelet package developed and provided by Torrence and Compo (available online at http://paos.colorado.edu/research/wavelets/) in our analysis. In addition, we restrict our analysis to the application of the Morlet wavelet since it has been used extensively (Torrence and Compo, 1998; Yan et al., 2001). The Morlet wavelet is a plane complex wave modulated by a Gaussian. The exact shape of the Morlet wavelet is determined by the non-dimensional frequency o0 that can be modified in order to change the wavelet resolution properties as discussed below. The wavelet time series at different scales s are estimated by translating and dilating the base wavelet. The wavelet scale s can then be related to the frequency o0 : For nonorthogonal wavelets, such as the Morlet, there is no preset basis of scales at which the wavelet analysis is performed. Thus, the choice of scales is arbitrary. However, a set of scales based on fractional powers of 2 is frequently used to perform the wavelet analysis (Torrence and Compo, 1998)

Like the wavelet method, the AWFT (see Zurbenko and Porter (1998) for details) is designed to extract time– frequency information from time series. Following Zurbenko and Porter (1998), the AWFT of a time series X ðtÞ at a given frequency o is computed as follows: m X 1 Y ðtÞ ¼ ei2pot X ðt þ sÞ: ð5Þ 2m þ 1 s¼m

Sj ¼ S0 2j dj ; J ¼ dj 1 log ðN dt=S0 Þ;

j ¼ 0; 1; 2; y; J;

ð2Þ

where S0 ð¼ 2 dtÞ is the smallest resolvable scale and J determines the largest scale. The letter N represents the length of data, dt stands for the sampling interval, and dj is the spacing between the discrete scales. For the Morlet wavelet, dj should be less than about 0.5 to ensure adequate scale (frequency) resolution (Torrence and Compo, 1998). For the Morlet wavelet, the scale s can be related to the frequency o as follows: qffiffiffiffiffiffiffiffiffiffiffiffiffiffi o0 þ 2 þ o20 o¼ ; ð3Þ 4ps where o0 is the Morlet wavelet dimensionless frequency tuning parameter. The time series for a frequency band of interested can be obtained by summing the real part of the wavelet transform over the corresponding scales (frequencies). The wavelet filtered time series over a set of scales may be written as j2 RfWn ðsj Þg dj dt1=2 X Xn ¼ ; 1=2 Cd C0 ð0Þ j¼j1 s

ð4Þ

j

where Wn ðsÞ represent the wavelet coefficients, and c0 ð0Þ and Cd are constants for each wavelet function. For a Morlet wavelet, Cd is 0.776 and c0 ð0Þ has a value of p1=4 :

Similar to the KZ moving average filter, the AWFT consists of k iterations in which the input time series X ðtÞ for iteration i þ 1 is replaced by the Y ðtÞ output from iteration i: The parameter m determines the half length of the averaging window. The window length is taken to be a multiple of the period of interest: m ¼ no1 ;

ð6Þ

where n is an integer. Thus, the AWFT filter is described by three parameters m; k; and o: For a given frequency o; the other two parameters determine the time and frequency resolution properties. In addition, since the length of the averaging window is a function of the frequency o; the AWFT avoids problems of the windowed Fourier transform (WFT) pointed out by Torrence and Compo (1998). The WFT uses a sliding window of fixed length on which a Fast Fourier Transform (FFT) is performed to estimate the local spectral density at all frequencies. This leads to oversampling of high frequencies and undersampling of low frequencies. Further, the AWFT filter is also able to work in a missing data environment as discussed in Section 3.1.

3. Results and discussion In the following, we compare the different filters as applied to four commonly encountered tasks in the analysis of time series of meteorological and air quality variables, namely, extraction of time series of a signal with known and fixed frequency, determination of process energies within different spectral bands, extraction of time series of components within certain spectral bands, and analysis of temporal variations (trends) in the strengths of different frequency bands. 3.1. Extraction of fixed frequencies We compare two methods for extracting a signal of known and fixed frequency from a time series. This problem is typically encountered when one wishes to extract and study the time-varying amplitude of the 24 h diurnal or the 365 day annual oscillation embedded in time series of many atmospheric variables. Additionally, if this procedure is repeated over a range of frequencies for a given time series, it is possible to generate

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

316

time–frequency plots which depict the temporal evolution of the different fluctuations embedded in the time series. Both the Morlet wavelet and AWFT can be used for this type of analysis. For the Morlet wavelet, the desired sampling frequency must first be converted to the dimensionless scale parameter s (Eq. (3)) that can then be used to construct the appropriate wavelet (Torrence and Compo, 1998). In addition, as discussed in Section 2.5, the parameter o0 determines the sharpness of the scale (frequency) resolution of the wavelet. For the AWFT, the desired frequency is o in Eq. (5), and the time–frequency resolution properties of the AWFT are determined by the parameters m and k: Figs. 1a and b display the energy transfer functions for the Morlet wavelet and AWFT filters centered at a period of 24 h with different resolution parameters o0 and m; respectively. For the three values of o0 (4, 8, and 20) shown for the Morlet wavelet, the dimensionless scale parameter s is 15.7, 30.8, and 76.5, respectively.

This ensures that the Morlet wavelet transfer function is centered at a period of 24 h. It is evident that these transfer functions have very similar shapes, and that the sharpness of frequency separation can be controlled effectively by the parameters o0 and m: It is important to bear in mind that a sharper resolution in the frequency domain comes at the expense of decreased localization in the time domain. Both methods can be tuned to achieve the desired trade-off between resolution in the time and the frequency domain. There is one significant difference between these two methods: namely, their ability to work in a missing data environment. Most frequently, wavelet transforms are computed with an FFT, utilizing the convolution theorem for computational efficiency (Torrence and Compo, 1998). However, since the FFT method cannot handle missing data, it is necessary to generate data for periods of missed observations, by using either interpolation or sophisticated statistical models (Yan et al.,

(a)

(b)

Fig. 1. Energy transfer functions for the (a) AWFT and (b) Morlet wavelet filters centered at a period of 24 h for different tuning parameters m and o0 ; respectively.

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

2001). Obviously, such data restoration procedures conceivably alter the underlying spectral properties of the process (i.e. unwanted frequencies might be introduced depending on the procedure used for filling in the missing data) and, therefore, the results of the wavelet transform. The computational implementation of the AWFT does not require a substitution for missing data and considers only non-missing data for the summation shown in Eq. (5). Only when the percentage of missing data in a given summation window is less than a userdefinable constant pcmiss (e.g. 25%), the term in Eq. (5) is actually computed; otherwise, a missing value is returned by the subroutine. Thus, the AWFT for a given point in the time series is only computed if enough data are present in the window surrounding the point to compute the value given by Eq. (5). Since the userdefinable constant pcmiss is specified as a percentage of the averaging window length m (Eq. (6)), the estimation of AWFT for higher frequencies (shorter periods and, thus, shorter averaging window) will be more affected by the presence of missing data than the estimation of AWFT for lower frequencies. To illustrate the AWFT’s ability to handle missing data, Fig. 2a displays the simulated time series described in Section 2.1 with 60 days flagged as missing. An AWFT filtration was performed to extract the 24 h and 365 day waves from this time series; the parameters m and pcmiss were set to 10* period and 25%, respectively. It can be seen in Figs. 2b and c that the AWFT for the 24 h wave was not estimated for a certain period because of missing data, while the AWFT for the 365 day wave can be estimated for the entire time period, with no apparent effect of the 60 days of missing data. Also, note that the amplitude of the estimated 24 h wave (the actual cycle is not visible because the displayed time series contains more than 800 of these cycles) fluctuates around the ‘‘true’’ amplitude of 1 due to the presence of other surrounding frequencies in the simulated time series. As stated above, the AWFT routine automatically handles missing data, and the AWFT filtered time series might still contain periods of missing data if not enough information was present to calculate the transform. When plotting or analyzing the AWFT filtered time series, the missing data can easily be omitted since they are clearly flagged as such. Since missing data is a very common problem in air quality and meteorological time series, this might be an important consideration when deciding which filter to use. 3.2. Extraction of different frequency bands In this section, we discuss techniques for separating a time series into different spectral bands. In contrast to the previous section, the object here is to separate different frequency bands rather than separate a single

317

oscillation with a fixed frequency from the original time series. This problem typically arises when one wants to determine the relative contribution of different dynamical processes to the overall spectral energy (variance), and when some of these processes such as synoptic-scale weather variations or inter-annual variability in climate are characterized by a broad band of frequencies rather than by a fixed frequency; the elliptic and KZ filters are explicitly designed for this purpose. Wavelet transform and AWFT can be adapted to extract a band of frequencies by performing multiple wavelet/AWFT filtrations at various sampling frequencies within the desired band and then adding these filtered time series. The four filtering techniques discussed above have been applied to the simulated time series described in Section 2.1. The spectral representation of this time series is flat (equivalent to that of a white noise process) with two peaks corresponding to the DU and annual oscillations. We are interested in separating five different bands of frequencies representative of different physical processes typically present in time series of atmospheric variables. The five time scales of interest are the intraday (ID) component with periods less than 11 h, the diurnal (DU) component with periods ranging from 11 h to 2.5 days, the synoptic (SY) component with periods from 2.5 to 21 days, the seasonal (SE) component with periods covering 21 days to 2.5 years, and the long-term (LT) component with periods greater than 2.5 years variations. The ID component contains fast-acting, local-scale processes; the DU component is dominated by the effect of day and night differences on atmospheric variables; the SY component contains fluctuations related to changing weather patterns; the SE component is dominated by the effect of the annual solar cycle on atmospheric variables; and the LT component contains the low-frequency part of the signal such as inter-annual climate variability or trends. For the KZ and elliptic filter, once the four separation frequencies between the five components are specified (i.e., 11 h, 2.5, 21 days, and 2.5 years), the desired components can be estimated directly by specifying the appropriate filter parameters as described in Section 2. For the KZ filter, m values of 3, 13, 103 and 8760 and k values of 3, 5, 5, and 3, respectively, were selected. For the wavelet and AWFT estimated ID component, all wavelet and AWFT components with periods less than 11 h were added; for the DU time scale, all wavelet and AWFT components with periods between 11 h and 2.5 days were added, and so on. The wavelet analysis presented here used 34 dimensionless scales s (related to periods via Eq. (3)) ranging from 2 (21) to 92,682 (216.5). The AWFT transform was performed at the 34 frequencies corresponding to these scales for the sake of consistency. Figs. 3–6 depict the power spectra of the five components extracted with the wavelet, AWFT, KZ,

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

318

(a)

(b)

(c)

Fig. 2. (a) Simulated time series consisting of a 24 h sine wave and a 365 day sine wave with amplitudes 1 superimposed on white noise with standard deviation 2.39 with 60 days of missing data, (b) 24 h oscillation extracted using the AWFT filter, and (c) 365 day oscillation extracted using the AWFT filter.

and elliptic filters, along with the power spectrum of the original simulated time series. These power spectra were obtained by applying the respective filters and can be viewed as their empirical transfer functions when compared to the power spectrum of the original time series. The dashed vertical lines in these figures represent the four periods separating the five components of interest. The power spectra for these components estimated by ideal filters would have a box-car shape, i.e., an ideal filter would transfer all spectral energy within the frequency band of interest and none of the spectral energy outside. For example, the power spectrum of the ID component estimated by an ideal

filter would have a spectral energy of 3  105 for all periods up to 11 h and drop to 0 for higher periods. Departures of the shape of the power spectra of real filters from such an ideal shape imply that the separation between different components is not perfect (causing non-zero covariances among the separated components) and that not all spectral energies within the frequency band of interest are transferred equally. The power spectra of the ID, DU and SY bands extracted by the elliptic filter show the desirable steep drop off from pass bands to stop bands, but the SE and LT components do not (Fig. 3). It can also be seen from the flat top portion of the power spectra that the

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

319

Fig. 3. Power spectra of the original simulated time series as described in the text and the ID, DU, SY, SE, and LT components estimated by the elliptic filter. The dashed vertical lines represent the desired separation points between the five components of interest.

application of the elliptic filter did not alter the spectral energies within the frequency bands of interest. The power spectra for the KZ filter in Fig. 4 reveal a significant departure from the box-car shape, indicating the presence of covariances among the separated components and an uneven transfer of the spectral energies within the frequency bands of interest. The power spectra for the wavelet and AWFT transforms in Figs. 5 and 6 display a steep roll-off from the pass bands to the stop bands for the five spectral bands of interest. However, because the components are estimated by adding multiple time series sampled at different frequencies, not all spectral energies within the desired frequency bands are transferred equally, as evidenced by the ripple in the pass bands. The amount of ripple for

both methods can be controlled but not eliminated by changing their tuning parameters as discussed in Section 3.1. To provide a quantitative assessment of the different methods’ abilities in separating the simulated time series into the five temporal components of interest, we computed the variance of each of the estimated components as well as the total covariances between the estimated components, and compared them to their known counterparts in the simulated time series (Table 1). It should be noted that—except for the KZ filter— the estimated variances and covariances in this table were normalized to facilitate the comparison of the performance of different methods. It can be seen that both AWFT and wavelet filters closely preserve the true

320

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

Fig. 4. Power spectra of the original simulated time series as described in the text and the ID, DU, SY, SE, and LT components estimated by the KZ moving average filter. The dashed vertical lines represent the desired separation points between the five components of interest.

distribution of variance among different temporal components and have covariances less than 5%. Both elliptic and KZ filters show covariances greater than 10% among the separated components. The correlation coefficients between the estimated and true simulated time series for each component are shown in Table 2. It can be seen that the correlations generally are highest for the wavelet filter, followed by the AWFT, KZ and elliptic filters. The poor correlations of the components estimated by the elliptic filter with the ‘‘true’’ simulated time series are likely attributable to the phase shift introduced by the elliptic filter as mentioned in Section 2. Fig. 7 shows phase shift diagrams for each of the five components estimated by the elliptic filter. It is evident that the elliptic filter

introduces a frequency-dependent phase shift to each of the components. Thus, while the elliptic filter can be a useful tool to study the variance of different components, the phase shift makes a meaningful interpretation of the temporal variation of extracted components extremely difficult. It can be seen that for all methods except the KZ filter, the correlation decreases for the LT component. It is important to note that while the LT components estimated by the KZ and elliptic filters are open-ended (i.e. the transfer functions show the characteristic of a low-pass filter)—the wavelet and AWFT estimates of the LT components consist of the superposition of a finite number of periods with a fixed longest period determined by the length of the data set. In other words,

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

321

Fig. 5. Power spectra of the original simulated time series as described in the text and the ID, DU, SY, SE, and LT components estimated by the Morlet wavelet filter. The dashed vertical lines represent the desired separation points between the five components of interest.

while the KZ and elliptic filter pass very low frequencies (even though they are not entirely contained within the data set), the wavelet and AWFT LT estimates do not contain periods that are only partially contained in the data set. However, a low-pass LT component could be estimated from the AWFT and wavelet methods by taking the difference between the original time series and the sum of all wavelets/AWFT waves up to a certain frequency. These residuals would then contain all fluctuations with frequencies lower than the selected frequency, including partial oscillations and the lowfrequency components of a linear trend. Thus, these

residuals would then have properties similar to the LT component estimated by the KZ and elliptic filters. 3.3. Application to real-world data Spectral decomposition techniques have been applied in the past for the analysis of meteorological and air quality variables (Eskridge et al., 1997; Rao et al., 1997). As the results presented above illustrate, the choice of methods for a particular application may affect the outcome of the analysis. However, for applications that are rather qualitative in nature (such as variance

322

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

Fig. 6. Power spectra of the original simulated time series as described in the text and the ID, DU, SY, SE, and LT components estimated by the AWFT filter. The dashed vertical lines represent the desired separation points between the five components of interest.

Table 1 Variance and covariance (expressed in percent of total variance) for the five spectral components estimated by the four different filtering techniques from the ‘‘true’’ simulated time series described in Section 2.1

ID variance DU variance SY variance SE variance LT variance Sum of variances Total covariances Variances+covariances

True (%)

KZ (%)

Elliptic (%)

AWFT (%)

Wavelet (%)

69.7 20.2 2.4 7.7 0.007 100 0 100

57.5 17.1 2.2 7.7 0.005 84.5 15.5 100

62.9 19.6 2.6 3.3 0.06 88.5 11.5 100

72.3 17.7 1.8 8.1 0.002 99.9 0.1 100

65.5 20.6 2.4 8.2 0.003 96.7 3.3 100

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

distribution analysis or air quality model evaluation), the impact of the filters’ different operating characteristics probably is of relatively minor concern as long as the results are interpreted prudently. As shown in Table 3 for the distribution of variance analysis for the observed ozone time series at Bayonne, NJ, the use of different filters does not affect the ranking of the different estimated components in terms of their relative importance. One exception to this is the elliptic filter, which underestimates the importance of the SE component relative to the SY component. The detection of changes in atmospheric pollutant concentrations and their possible attribution to anthropogenic activities or policy measures is another area of great interest. In this context, the analysis of changes in pollutant concentration by linear trend analysis is often Table 2 Correlation coefficient between the ‘‘true’’ and estimated time series of five spectral components

R(ID Estimated, ID True) R (DU estimated, DU true) R (SY estimated, SY true) R (SE estimated, SE true) R (LT estimated, LT true)

KZ

Elliptic

AWFT

Wavelet

0.32

0.67

0.34

0.21

0.82

0.52

0.90

0.93

0.83

0.30

0.93

0.97

0.99

0.87

0.99

0.99

0.94

0.29

0.60

0.67

These components were estimated by four different filtering techniques from the ‘‘true’’ simulated time series described in Section 2.1.

323

not applicable (Reid et al., 2001; Civerolo et al., 2001; Porter et al., 2001). Instead, it may be more meaningful to examine the changes in the magnitudes of different spectral components (reflecting different dynamical processes) over time. This approach has been taken by Chan et al. (1999) and more recently by Yan et al. (2001). In their study, Yan et al. (2001) investigated long term changes in the strength of weather and SE cycles by applying the Morlet wavelet to extract time series in the SY and SE wavebands from temperature time series and then examined how the strength of these fluctuations changed over time. In a similar effort, we applied the four filtering methods to study the LT behavior of the strength of the DU component estimated from the observed ozone time series at Pasadena from 1984 to 1999. To this end, after the component time series were estimated by each of the four methods, their standard Table 3 Variance distribution (expressed in percentage of total variance) for five spectral components estimated from the four different filters from hourly ozone observations at Bayonne, NJ

ID variance DU variance SY variance SE variance LT variance Sum of variances Total covariances Variances +covariances

KZ (%) Elliptic (%)

AWFT (%)

Wavelet (%)

4.5 30.2 15.6 26 0.3 76.6 23.4 100

7.1 45.7 16.5 31.3 0.06 100.7 0.7 100

6.9 43.5 18.5 25.8 0.05 94.8 5.2 100

9.2 38.1 18.4 11.7 5.1 82.5 17.5 100

Fig. 7. Phase shift of the ID, DU, SY, SE, and LT components estimated by the elliptic filter relative to the phase of the original time series as a function of period. The dashed vertical lines represent the separation points between the five components.

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

324

deviations were determined for each individual year. Since the standard deviation is proportional to the amplitude of the signal, a plot of the standard deviations for each year for each of the components can provide information about possible changes in the strengths of these components over time and about the ability of the different filters to discern these changes. The results, shown in Fig. 8, indicate two features. First, all methods estimated the same year-to-year variation in the strength of the DU components. Second, the strength of the DU component clearly has decreased over the time period analyzed here. This means that the decrease in the observed ozone concentrations at this site (Chan et al., 1999) is in part due to the decrease in the strength of the DU component. In summary, all of the examined methods are suitable for examining long-term changes

in the magnitudes of individual temporal components, and such analysis—when performed for different components—may be more informative about the nature of long-term changes than a simple linear trend analysis.

4. Summary Four methods for spectrally decomposing time series of atmospheric variables have been compared. Two of these (the Morlet wavelet and the AWFT filters) can be applied both to extract time series at a given frequency as well as time series within a certain frequency band (like a bandpass filter). The KZ and elliptic filters are often used to extract certain frequency bands. In this

(a)

(b)

(c)

(d)

Fig. 8. Strength of the DU component in different years. The DU component was extracted by (a) the AWFT filter, (b) the Morlet wavelet filter, (c) the KZ moving average filter, and (d) the elliptic filter.

C. Hogrefe et al. / Atmospheric Environment 37 (2003) 313–325

study, the AWFT filter was introduced to the analysis of atmospheric variables; its performance characteristics have been found to be very similar to the Morlet wavelet. In addition, the application of the AWFT is conceptually simple and, more importantly, it can work in the presence of missing data. When the four filters were used in a simulation study to extract five frequency bands of interest when analyzing atmospheric variables, they all gave qualitatively similar results in terms of the variance distribution among components. However, the phase shift introduced by the elliptic filter makes a meaningful interpretation of the temporal variation of extracted components extremely difficult. In addition, all filters have been able to capture year-to-year fluctuations in the strengths of individual components. It was shown that such analysis can be used to discern the time scales contributing to long-term changes (trends) in ambient pollutant concentrations. Since no single filter performs best for all applications, great care has to be taken in choosing the filter for the intended application.

Acknowledgements This research was supported by the New York State Energy Research and Development Authority under contract no. 4914-ERTER-S99. The authors would like to thank Dr. Torrence and Dr. Compo for providing the wavelet software used in this analysis (available online at http://paos.colorado.edu/research/wavelets). Sincere thanks are extended to Dr. R.F. Henry for many helpful discussions on scale separation.

References Chan, D., Rao, S.T., Zurbenko, I.G., Porter, P.S., 1999. Linking changes in ozone to changes in emissions and meteorology. In: Brebbia, C.A., Jacobson, M., Power, H. (Eds.), Air Pollution VII. WIT Press, Southampton, UK, pp. 663–675. Civerolo, K.L., Brankov, E., Rao, S.T., Zurbenko, I.G., 2001. Assessing the impact of the acid deposition control program. Atmospheric Environment 35, 4135–4148. Eskridge, R.E., Ku, J.-Y., Rao, S.T., Porter, P.S., Zurbenko, I.G., 1997. Separating different scales of motion in time series of meteorological variables. Bulletin of the American Meteorological Society 78, 1473–1483.

325

Foufoula-Georgiou, E., Kumar, P. (Eds.), 1995. Wavelets in Geophysics. Academic Press, San Diego, CA, 373pp. Hogrefe, C., Rao, S.T., Zurbenko, I.G., Porter, P.S., 2000. Interpreting information in time series of ozone observations and model predictions relevant to regulatory policies in the eastern United States. Bulletin of the American Meteorological Society 81, 2083–2106. Lau, K.-M., Weng, H.-Y., 1995. Climate signal detection using wavelet transform: how to make a time series sing. Bulletin of the American Meteorological Society 76, 2391–2402. Milanchus, M.L., Rao, S.T., Zurbenko, I.G., 1998. Evaluating the effectiveness of ozone management efforts in the presence of meteorological variability. Journal of Air & Waste Management Association 48, 201–215. Porter, P.S., Rao, S.T., Zurbenko, I.G., Dunker, A.M., Wolff, G.T., 2001. Ozone air quality over North America: Part II—an analysis of trend detection and attribution techniques. Journal of Air & Waste Management Association 51, 283–306. Poularika, A.D., 1998. The Handbook of Formulas and Tables for Signal Processing. CRC Press, Boca Raton, FL. Rabiner, L.R., Gold, B., 1975. Theory and Application of Digital Signal Processing. Prientice Hall, Englewood Cliffs, NJ. Rao, S.T., Zurbenko, I.G., 1994. Detecting and tracking changes in ozone air quality. Journal of Air & Waste Management Association 44, 1089–1092. Rao, S.T., Zalewsky, E.E., Zurbenko, I.G., 1995. Determining temporal and spatial variations in ozone air quality. Journal of Air & Waste Management Association 45, 57–61. Rao, S.T., Zurbenko, I.G., Neagu, R., Porter, P.S., Ku, J.Y., Henry, R.F., 1997. Space and time scales in ambient ozone data. Bulletin of the American Meteorological Society 78, 2153–2166. Reid, N., Misra, P.K., Bloxam, R., Yap, D., Rao, S.T., Civerolo, K., Brankov, E., 2001. Do we understand trends in atmospheric sulfur species? Journal of Air & Waste Management Association 51, 1561–1567. Sebald, L., Treffeisen, R., Reimer, E., Hies, T., 2000. Spectral analysis of air pollutants part 2 ozone time series. Atmospheric Environment 34, 3503–3509. The Mathworks Inc., 2001. Signal Processing Toolbox, available online at http://www.mathworks.com/access/helpdesk/ help/toolbox/signal/signal.shtml. Torrence, C., Compo, G.P., 1998. A practical guide to wavelet analysis. Bulletin of the American Meteorological Society 79, 61–78. . H., Davies, T.D., Yan, Z., Jones, P.D., Moberg, A., Bergstrom, Yang, C., 2001. Recent trends in weather and seasonal cycles: an analysis of daily data from Europe and China. Journal of Geophysical Research 106, 5123–5138. Zurbenko, I.G., Porter, P.S., 1998. Construction of highresolution wavelet. Signal Processing 65, 315–327.