## Frequency a1Frequency a1

Figure 5.11. Monsoon spectrum, influence of timescale errors. Focus is on two portions (a, b) of the original spectrum of the proxy record Q5 (Fig. 5.8). Errors from the age determination (Fleitmann et al. 2003: Table S1 therein) were used for timescale resampling (Section 4.1.7.2) of a piecewise linear age-depth model. The 90% bootstrap bound for the red noise (the lower of the grey lines in a and b), obtained from B = 10,000 simulations, is higher than the 90% bound (the lowest of the four black, smooth lines in a and b) obtained from ignoring timescale errors. Also shown (the upper of the grey lines in a and b) is the increase in the 99.9% red-noise bound (relative to the uppermost of the four black, smooth lines in a and b), obtained from B = 100,000 simulations. The frequency uncertainties (horizontal bars) due to timescale errors (Algorithm 5.6, B = 2000), expressed as standard errors, sef/, are compared with the half of the 6-dB bandwidth, Bs/2.

It would be premature for an analysis of the Sun-monsoon relation to stop at this point. Three lines should be explored. First, the relation can be further investigated, using the same data sets, in the time domain by means of bandpass or harmonic filtering (Section 5.2.4.3). Second, the climate physics of the Sun-monsoon link can be considered. This has been done by Kodera (2004), who explained a positive correlation between solar activity and Indian monsoon strength via a weakening of the Brewer-Dobson circulation in the stratosphere. However, this was established on measurement data from 1958-1999, and the feasibility of this or other mechanisms on longer timescales is still elusive. Third, other records of Holocene monsoon variations need to be analysed. For example, Neff et al. (2001) analysed a ¿18O record from a stalagmite from another cave than where Q5 is from, finding monsoon peaks at Tperiod = 10.7 a, 226 a and 1018 a. Combining this evidence with the information from Q5 in a new multiple test should raise the overall statistical significance. A synopsis of evidence pro and contra the Sun-monsoon hypothesis in a multiple statistical test, with timescale errors taken into account, is a major task awaiting to be done.

### 5.3 Background material

Overviews of spectral analysis have been given plentiful. Accessible ones include the classic textbooks by Priestley (1981) and Percival and Walden (1993). The latter takes another route to the definition of the power spectrum than given here in Section 5.1; it also focuses on multitaper methods. John W. Tukey's work on time series and spectrum analysis, done from the 1950s to the 1990s (and summarized by Brillinger (2002)) contains useful material for the practitioner. Kay and Marple (1981) and Bendat and Piersol (1986) offer perspectives from the engineering side. A review of quadratic spectrum estimators, including multitaper methods, was given by Mullis and Scharf (1991). Reviews of spectrum analysis written by geoscientists include MacDonald (1989), Park (1992) and Ghil et al. (2002). Stratigraphy is a geoscientific sub-field dealing with archived temporal changes in lithic or biotic units such as sediments. Such changes are often cyclical (Einsele et al. 1991) and quantifiable by means of spectral analysis (Weedon 2003). This area, cyclostratigraphy, deals also with series, {z(i),x(i)}, from the depth domain (Fig. 1.13). A prominent example is the identification of Mi-lankovitch cycles in sequences such as limestone-shale from the distant geologic past (Schwarzacher 1964, 1975, 1991, 1993, 1994). There is a discussion whether we are now living in the Anthropocene (Crutzen 2002; Crutzen and Steffen 2003; Zalasiewicz et al. 2008).

Periodogram tests belong to the historically earliest tools in statistical time series analysis. A review, on which this paragraph is oriented, was given by Priestley (1981: Section 6.1.4 therein), see also Priestley (1997). The periodogram was not only invented by Schuster (1898), this man also devised a test for the significance of an I (fj) value based on the assumption of Gaussian white noise in Eq. (5.15) and the chi-squared distribution. See Brillinger (1975) for a rigorous description of CIs and other properties of the periodogram. Schuster (1898) applied his test to a "supposed 26 day period of meteorological phenomena" (the declination of the Earth's magnetic field at Prague, measured during 1870 with ( = 1 day; the supposition had been made by Hornstein (1871)), but found little evidence in favour of a true periodicity. Later, Schuster (1906) analysed monthly sunspot time series for the interval 1749 to 1901 and detected periodicities, the major one at Tperiod = 11.125 a. In that paper, Schuster also considered nonstationarity. Gilbert Walker, a physicist with contributions to meteorology, looked on max(I(/j))/S2, with S2 replaced by the sample variance estimator (Walker 1914), and found an asymptotic distribution. Fisher (1929) derived the exact distribution of a related test quantity for n odd, and also Hartley (1949) took Walker's test statistic, changed the denominator and derived the distribution of this re-studentized quantity. It is natural to test not only for one (max(I(/j))) but also for more harmonic components in a time series, and relevant work on this topic includes that published by Whittle (1952), Grenander and Rosenblatt (1956), Siegel (1980) and Walden (1992). A test for the number K of frequencies to include in the harmonic model (Eq. 5.15) was developed by Quinn (1989). A test for peaks in the spectrum estimated with maximum likelihood (instead of periodogram estimation) was presented by Foias et al. (1988). A serious caveat against all tests described so far in this paragraph is their assumption of a white Gaussian background noise against which to test. We assume climate processes to have rather a non-white background, that is, to exhibit a mixed spectrum (Fig. 5.3c). Statistical tests can still be constructed for mixed spectra based on analysis of I(/j)/h(/), that means, the periodogram divided by the power spectral density function of the background process. The serious practical problem here is that h(/) is usually unknown and has to be replaced by an appropriate estimate, and obtaining a background estimate requires in principle the harmonic peaks to be detected. If the background spectrum has a narrow local maximum (e.g., an AR(2) spectrum), then it may be impossible to distinguish between background maximum and periodogram peak (noise and signal). Periodogram test methods to deal with such a situation require adapted background spectrum estimation (Whittle 1952; Hannan 1960, 1961; Priestley 1962a,b). An interesting alternative to periodogram tests is Thomson's F test using the background spectrum estimated with the multitaper technique (Section 5.2.3.1). The tests described so far in this paragraph were developed under the assumption of even time spacing. There exists a test using the Lomb-Scargle periodogram for unevenly spaced time series (Scargle 1982; Horne and Baliunas 1986), which is similar to Schuster's (1898)—including the restrictive assumption of white background noise. In their review of Lomb-Scargle periodogram analysis, VanDongen et al. (1997) and Van Dongen et al. (1999) mention the permutation test by Linnell Nemec and Nemec (1985). However, the permutation resampling does not preserve redness, and, hence, also this test assumes a white background. Summarizing, periodogram tests may be useful for analysing processes with line components and little/no background noise (e.g., astronomical cycles, tides), but they have little relevance for time series from climate processes with mixed, non-white spectra. We do not share the view of Muller and MacDonald (2000: p. 56 therein) that Priestley's (1981: p. 420 therein) remark of the unusefulness of the periodogram for the estimation of continuous spectra is misleading. On page 431 of his book, Priestley defends the periodogram's usefulness for estimating line spectra (Fig. 5.3a). The question raised by Muller and MacDonald (2000) is more whether their study object, Milankovitch cycles embedded in climate noise in the form of Pleistocene ice-volume changes, should indeed be analysed by means of periodogram estimation.

Superresolution refers to (almost) purely harmonic processes with a line spectrum, where a higher frequency resolution (i.e., < Afj) than for spectrum estimation can be achieved (Thomson 1990a). Fields for application in climatology include frequency estimation and separation of tide components (Munk and Hasselmann 1964). Also Hannan and Quinn (1989) and Quinn and Hannan (2001) studied frequency separability in dependence of n, S2 and the amplitude of the sinusoidal components. The latter book contains further statistical tests and considers also nonstationarity in form of slowly changing frequencies.

Nonstationarity in the context of this chapter has something to do with a "time-dependent spectrum." The problem is that this is not well defined; some assumptions have to be met, and some variables to be introduced, to be able to speak of a "time-dependent frequency" or a "time-dependent power," as was reiterated by Priestley (1981, 1988). One assumption is that the time-dependences are slow and smooth. It is possible to erect a "nonstationary spectral analysis" on wavepackets or wavelets (Fig. 5.12) that have an oscillatory and a smoothly damped part (Priestley 1996). The estimation means to effectively compose a time series, {t(i),x(i)}™=1, using shifted (in time) and scaled (in time) versions of the "mother wavelet," ^(T). Most estimation algorithms seem to require (1) even spacing and (2) n to be a power of two. Evidently, interpolation methods can free the climate time series analyst from those two strong restrictions, but this seems to be at the expense of introducing heteroscedasticity and introducing or enhancing autocorrelation (Silverman 1999). The second effect could make tests of red-noise alternatives more difficult, the first require the analyst to reunite S(T) with Xnoise(T). It is fair to say that a systematic and wide knowledge about interpolation effects on time-dependent spectrum estimates obtained with wavelets is not yet reached (Daubechies

Time, T

Figure 5.12. Wavelet, "Mexican hat" function, ^(T) = (1 - T2)exp(-T2/2).

Time, T

Figure 5.12. Wavelet, "Mexican hat" function, ^(T) = (1 - T2)exp(-T2/2).

et al. 1999; Vidakovic 1999; Sweldens and Schroder 2000). (For the case of nonlinear wavelet regression (time-dependent mean), Hall and Turlach (1997) studied interpolation effects theoretically and by means of a Monte Carlo simulation.) Papers on the application of wavelet estimation for unevenly spaced astronomical/climatological time series include Foster (1996a,b), Scargle (1997), Witt and Schumann (2005) and Milne and Lark (2009). A recent contribution from theory is Mon-dal and Percival (in press), who propose new, unbiased estimators of the wavelet power spectrum for even spacing with missing observations, analyse their large sample properties and methods for CI construction. Mondal and Percival (in press) show also an application to annual runoff minima from the Nile for the interval from a.d. 622 to 1921. To summarize, wavelet models may offer many new insights into time-dependent climate processes, but more theoretical and simulation work needs to be done, and software tools to be developed, to understand the robustness and accuracy of results with respect to uneven spacing, aliasing and timescale errors. Another technique applicable to slowly changing time-dependent spectra is to form time intervals ("windows") and estimate the spectrum separately for the windows. For example, Berger et al. (1998) studied the stability of the Milankovitch periods of variations in the Earth's orbital geometry over the interval from 1.5 Ma ago to 0.5 Ma into the future by means of a windowed multitaper estimation. Urban et al. (2000) applied the same method to look on the ENSO history within 1840-1995 as provided by the ¿18O proxy record from a coral (d = 2 months) from the central western Pacific. The ENSO spectrum exhibits power in the range of 2.2-15 years period, more broadly and not in the form of sharp peaks, and the analysis (Urban et al. 2000) shed light on the time-frequency composition of the ENSO. Schulz et al.

(1999) used windowed Lomb-Scargle estimation with WOSA for quantifying amplitude variations of the "1500-year cycle," recorded by the unevenly spaced ¿18O record from the GISP2 ice core. The caveats against wavelet estimation regarding robustness and accuracy of results apply also to windowed spectrum estimation techniques. One should bear in mind that the various difficulties in spectral analysis are ultimately rooted in the ambition to estimate the spectrum at frequency points fj, which are O (n), on the basis of a data sample of size n. Allowing time-dependence introduces a second dimension, and estimating a quantity at O (n2) time-frequency points cannot be expected to reduce the difficulties. Genton and Hall (2007) present an interesting alternative, namely estimation of parametric models for the time-dependences of frequency and amplitude, fitted in the time domain by means of kernel functions, supported by bootstrap CIs, and applicable also to uneven spacing. A notable tool for unevenly spaced series is also period analysis using robust regression in the time domain (Oh et al. 2004), which can be combined with bootstrap resampling.

The 100-ka cycle is the dominant type of changes of global ice volume and, related, temperature and atmospheric CO2 concentrations during the late Pleistocene (Fig. 1.3). Explaining this cycle is a challenge to the Milankovitch theorists for two reasons (Raymo and Huybers 2008). This theory of how variations in Earth orbital parameters influence climate has been successful regarding changes in the obliquity (Tperiod w 40 ka; Fig. 5.4) and precession

(Tperiod ww 19 — 23 ka) bands (Imbrie et al. 1992). First, the 100-ka cycle has a distinct sawtooth shape, which is absent in the more sinusoidal orbital time series. Second, the eccentricity component (ellipse) has a peak in that frequency range, but associated with clearly less power than the obliquity or precession components have. Ideas on how to reconcile Milankovitch theory with the 100-ka cycle include some nonlinear amplification of the eccentricity component in the climate system (Imbrie et al. 1993) and combinations of obliquity and precession components into a ~ 100-ka component (Raymo 1997; Huybers and Wunsch 2005). An astronomical cause, suggested not by Milankovitch but Muller, is variations in orbital inclination (Tperiod w 95 ka), see Section 5.2.1. Non-astronomical explanations view the 1/100 ka-1 as kind of an eigenfrequency of the ice-bedrock-carbon cycle system (DeBlonde and Peltier 1991; Saltzman and Verbitsky 1993). It is difficult to distinguish among the various explanations on basis of statistical analyses by spectrum estimation because the 100-ka cycle came into existence as late as approximately 650 ka ago, as found by Mudelsee and Schulz (1997) using a windowed version of the harmonic filter (Section 5.4). This short time span means large bandwidth and fre quency resolution, A/j. Spectrum estimation for the 100-ka cycle could be possibly improved by considering non-sinusoidal (i.e., sawtooth) basis functions (Thomson 1982).

The multitaper method with jackknife resampling for CI or standard error determination has been employed in studies of various aspects of the climate system. Among them are the following. Diaz and Pulwarty (1994) analysed centuries-long proxy records of ENSO and records of potentially related variables by means of a cross-spectral analysis (spectral analysis in a bivariate setting), and Rodo et al. (1997) tested the significance of peaks in the spectra of Iberian rainfall records from 1910 to 1994. Thomson (1997) examined variations of global temperature and the logarithm of solar irradiance during nineteenth and twentieth century and calculated jackknife CIs using = 6. This low number and the resulting large sampling fluctuations could explain the discrepancies in CI length he noted between the jackknife approach and one based on the Gaussian assumption. Hinnov et al. (2002) investigated the interhemispheric relations among various time series of D-O variations over the past 100 ka with cross-spectral analysis, and Prokopenko et al. (2006) applied the same method to solar insolation time series and a sedimentary record from Lake Baikal covering the past ~ 1.8 Ma. The criticism regarding the low values and jackknife replicates (Section 5.2.3.4) is not restricted to the paper by Thomson (1997). Also the other studies mentioned so far used similar low values. The jackknife should yield more accurate results when applied to long instrumental time series, such as oceanographic (Chave et al. 1997) or seismologic (Prieto et al. 2007), for which multisegmenting is possible.

The Lomb-Scargle method with bootstrap resampling for bias correction has been utilized in various climate studies. Among the analysed archives "containing" the unevenly spaced records are the following: stalagmites that provide ¿18O and ¿13C proxy evidence about changes in precipitation and temperature on Holocene and late Pleistocene timescales (Niggemann et al. 2003; Holzkamper et al. 2004; Fleit-mann et al. 2007a); Antarctic ice cores that give methanesulfonic acid proxy evidence about changes in winter sea ice extent over the past 100 years (Abram et al. 2007); a loess section from Nebraska, absolutely dated with radiocarbon and dosimeter technologies, that informs via colour parameters and organic carbon content about drought variations on Holocene timescales (Miao et al. 2007); a Pacific sediment core that supplies nitrogen isotopic proxy evidence about nutrient concentrations for phytoplankton growth (i.e., carbon sequestration) over the past 70 ka (De Pol-Holz et al. 2007); and, finally, a sediment core from Bear Lake

(Utah-Idaho) that documents via pollen content the regional vegetation and climate history over the past 225 ka (Jimenez-Moreno et al. 2007).

Aliasing and uneven spacing. The effects of sampling a continuous-time climate process X(T) with spectrum h(f) depend on the temporal spacing of the discrete sampling points, see Priestley (1981: Section 7.1.1 therein) and Masry (1984). Even spacing bears the risk of aliasing (Section 5.2.7). Different types of uneven spacing can be distinguished.

1. An independently jittered spacing (Eq. 5.36) amounts to a "disturbed" even spacing. (This is equivalent to the timescale model given by Eqs. (4.31) and (4.33) for evenly spaced Ttrue(i).) Independent jitter with ¿2 ^ d2 still leads to aliasing effects (Akaike 1960; Shapiro and Silverman 1960; Moore and Thomson 1991). Instead of a Gaussian, also another shape may be employed for the innovation term in the equation, the jitter. A model of a jitter uniformly distributed over the interval between —d/2 and +d/2 (excluding the endpoints) respects the condition of monotonic growth for a climate archive. Beutler (1970) shows that for frequencies below 1/(2d), this jitter model leads to an alias-free spectral estimation. This paper demonstrates also that spectral estimation can lead to meaningful results even when the {t(i)} are unknown and only their rank is known, a possible situation in paleoclimate time series analysis. The independent jitter model may be applicable to climate time series when an originally even spacing is superimposed by small time uncertainties (e.g., radar measurements, which are influenced by the travel time, other instrumental observations).

2. Dependent jitter means that the innovations in the spacing equation (Eq. 5.36) are autocorrelated. Here it is more difficult than for independent jitter to obtain analytical results on second-order properties of a process such as the spectrum (Thomson and Robinson 1996). The dependent jitter model may be the norm for many climate archives such as speleothems or sediment cores (Pisias and Mix 1988). This model was also used by De Ridder et al. (2006) for modelling the shell growth of a mollusk (a climate archive).

3. Poisson sampling refers to a more irregular spacing, where the times are realizations of a homogeneous Poisson process, that is, they are uniformly distributed (Chapter 6). Then the deviation from the case of even spacing is large and spectral estimation is alias-free (Shapiro and Silverman 1960).

However, above mentioned papers on aliasing do not assume application of Lomb-Scargle spectrum estimation. Instead, they study what results when (1) the {x(i)}™=1 are assumed to be realizations of a process sampled on a discrete, evenly spaced time grid and (2) a spectrum estimation method is used that assumes even spacing. Analytical results on aliasing seem hardly to exist for unevenly spaced time series analysed with the Lomb-Scargle method. Scargle (1989) states that aliasing is then diminished. This is supported by a Monte Carlo simulation study (Press et al. 1992: Fig. 13.8.1 therein), where the sampling was Poisson, n = 100 and d = 1. The sinusoidal component with prescribed frequency 0.81, larger than 1/(2(1"), was detected with high confidence. An interesting discussion was initiated by the suggestion (Wunsch 2000) that the so-called "1500-year cycle," found in late Pleistocene and Holocene climate proxy records, is an alias of the annual cycle. Meeker et al. (2001) made it clear that, at least for the Ca record from the GISP2 ice core (Mayewski et al. 1997), interval 30 to 36 ka, the annual cycle is not preserved because of a finite sample duration (D(i); Fig. 1.13) and diffusion, D'(i) > D(i). Nevertheless, this time interval displays Dansgaard-Oeschger variations in Ca (Meeker et al. 2001: Fig. 1C therein), for which visual inspection infers a period of roughly 1500 years. This argument against aliasing was accepted by Wunsch (2001). A climatological objection against the existence of the "1500-year cycle," however, is that the weak stationar-ity assumption (time constant second-order properties) is violated: this "cycle" is restricted to this time interval (D-O events 5, 6 and 7), as was shown for the GISP2 ¿1SO record (Schulz 2002). Amazingly, there may exist not a cycle but rather a "1500-year pacing" of the onset of D-O events, that is, the onsets (during the late Pleistocene) are not always separated by ~ 1500 years but sometimes by multiples of this period (Schulz 2002; Rahmstorf 2003).

Timescale error influences. The piecewise linear age-depth model with constraint "monotonic growth," which was used for analysing the effects of timescale errors on spectrum estimates for stalagmite Q5 (Section 5.2.9), had been suggested previously (McMillan et al. 2002) for sedimentary sequences in general, where the timescale is constructed using dating points and interpolation. These authors considered also the inclusion of interpolation error models.

Sun—climate connections on timescales between those of the by everyone experienced daily and annual cycles and, on the other hand, Milankovitch cycles (Tperiod ^ 19 ka) are not based on changes in the geometry but on solar activity variations. On shorter, decadal timescales (sunspot cycle), activity variations (Woods and Lean 2007) lead to a direct climate forcing smaller than that of greenhouse gas emissions (Hansen and Lacis 1990), although inclusion of solar activity is in dispensable for climate models to reproduce the observations (Hegerl et al. 2007b). Solar irradiance changes on timescales between 11 years (sunspots) and 1500 years (Bond et al. 2001) and their role for paleo-climate where assessed by numerous authors, including Cini Castagnoli and Provenzale (1997), Hoyt and Schatten (1997) and Bard and Frank (2006). The paper by Rind (2002) considers nearly all timescales (up to the age of the Earth). Pittock (1978) took a "critical look at long-term Sun-weather relationships" and detected statistical deficits in many papers that claimed a strong Sun-climate connection. However, one may in response extend the regret for an absent comprehensive synopsis of the Sun-monsoon hypothesis (Section 5.2.9) to the Sun-climate system. Required are multiple statistical tests, which take also timescale errors into account. Regarding the spectral peak at Tperiod = 963a (Fig. 5.8), it may correpond to reported cycles at approximately 900 year period (Schulz and Paul 2002; Rimbu et al. 2004; Wanner et al. 2008), but its nature (solar or not) deserves further analysis.

Resampling in the frequency domain. The periodogram has the property that the covariance between two points, COV [I(f1),I(f2)], vanishes for f1, f2 € {1/(nd), 2/(nd),...} under some conditions (normal shape, even spacing). This led to the idea to resample (ordinary bootstrap) periodogram values (frequency domain) and not residuals (time domain). One technique based on that is to resample I (fj) locally, that means close to a frequency of interest f', in order to determine a confidence interval for the spectrum estimate, h(f'), see Paparodi-tis (2002). This technique can be applied also to spectrum estimation with tapers (Politis et al. 1992; Politis and Romano 1992b). We have written the harmonic process (Eq. 5.15) with sinus and cosinus components; an alternative notation uses terms Aj cos(2nfjT(i) + $j), where $j is the phase. Estimation of the phase over the frequency range, the phase spectrum, becomes important for climatology in bivariate settings. There one is interested in leads and lags at a certain period between two climate variables. The second resampling technique mentioned in this paragraph resamples the phase spectrum estimates, while leaving the amplitude spectrum estimates (Aj) intact. The resampled data are transformed back into the time domain and serve as surrogate bootstrap data (Nordgaard 1992; Theiler et al. 1992; Kantz and Schreiber 1997; Hidalgo 2003). This surrogate data technique has been extended into the "wavelet domain" (Angelini et al. 2005). These techniques apply to evenly spaced time series, and adapting them to the case of uneven spacing should be worth the effort. One problem with spectrum estimation then, however, is that the Lomb-Scargle periodogram does not exhibit vanishing covariances (Section 5.2.4.2).

Other spectrum estimation methods than multitaper (even spacing) or Lomb-Scargle (uneven spacing) are not recommended. The Blackman-Tukey approach (Jenkins and Watts 1968) goes via Eqs. (5.7), (5.8), (5.9) and (5.10) on the sample level and takes the autocovari-ance or autocorrelation function, truncates it at a certain lag to smooth and transforms into the frequency domain. Because of the estimation bias and variance involved therein (Chapter 2), we subscribe to Thomson's (1990a: p. 543 therein) remark that "if your spectrum estimate explicitly requires sample autocorrelations you are almost certainly doing something wrong." Multitaper and Lomb-Scargle are nonparametric methods. Parametric methods include fitting AR(p) models in the time domain and taking the fitted model in the frequency domain. Various procedures of fitting led to various names associated with such methods: Yule-Walker, Levinson-Durbin, maximum entropy or Burg's algorithm (Percival and Walden 1993: Chapter 9 therein). We do not dispute their usefulness for fitting time-series models to data, but we remain cautious regarding parametric estimation of climate spectra. In the case of even spacing, the bias and variance properties of the estimates these parametric methods produce, should be less good than those of the optimal (least-squares sense) multitaper method. In the case of uneven spacing, these methods are not available without interpolation.

Interpolation of an unevenly spaced time series to equidistance for spectrum estimation is not recommended. Besides the ambiguity which interpolation type (linear, cubic spline, Akima spline, etc.) to take, interpolation leads generally to smoothing and distortion of the data (Belcher et al. 1994). It may introduce spurious peaks, especially at higher frequencies, as was demonstrated by Horowitz (1974) or Schulz and Stattegger (1997). Interpolation may also corrupt the test of the red-noise alternative (Section 5.2.5) because of the amount of serial dependence artificially introduced.

Bispectrum. A zero-mean process X(i) can be written in the Volterra expansion, to second order, as

+ 22 j • ef(O, a2)(i - j) •eF(O,ct2) (i -j,k where the g and are parameters and EF(0, a2)(i) is a random variable with mean zero, variance ct2 and distribution function F (Stine 1997). The linear term (the first on the right-hand side) of the expansion gives a complete description if F is Gaussian. The parameters g of this term are related to the spectral density function h(/). The nonlinear term is required for describing non-Gaussian or nonlinear processes. It is related (Stine 1997) to a function h(f1, f2), the bispectrum. Estimation of the bispectrum (Subba Rao and Gabr 1984) can add to a spectral characterization of an observed process. Muller and MacDonald (1997b,c) applied bispectral analysis to evenly spaced climatological and astronomical records in order to support their hypothesis that the 100-ka cycle corresponds to variations of orbital inclination. A limitation of bispectral analysis is that currently available implementations seem to be restricted to even time spacing. Notwithstanding this, researchers studying late Pleistocene global climate changes applied bispectral analysis to time series originally unevenly spaced (Hagelberg et al. 1991; King 1996; Rutherford and D'Hondt 2000). Little knowledge seems to exist on testing noise alternatives and quantifying the robustness of bispectral estimates against timescale errors.

The trade-off between variance and bandwidth of spectrum estimators has been referred to by many authors as "Heisenberg's uncertainty principle," although the latter is a concept from quantum physics. It would be more apt to speak of "Grenander's uncertainty principle," after, for example, the paper by Grenander (1958).

## Post a comment