Climate of the Earth system

Prof. Dr. Markus Meier
Leibniz Institute for Baltic Sea Research Warnemünde (IOW)
E-Mail: markus.meier@io-warnemuende.de

Spectrum#

  • = the Fourier transform of the auto-covariance function of the time series (presents the variance per frequency)

  • math. definition: Let \(\mathbf{X_t}\) be an ergodic weakly stationary stochastic process with auto-covariance function \(\gamma(\tau), ~\tau=0,\pm1,...~\). Then the spectrum (or power spectrum) \(\Gamma\) of \(\mathbf{X_t}\) is the Fourier transform \(\cal{F}\) of the auto-covariance function \(\gamma\):

(33)#\[\begin{split}\begin{align*} \Gamma(\omega) &= \cal{F}\{\gamma\}(\omega)\\ &= \sum_{\tau=-\infty}^{\infty} \gamma(\tau)e^{-2\pi i\tau\omega},~~~~~\forall \omega \in \left[-\frac{1}{2},~\frac{1}{2} \right] \end{align*}\end{split}\]
  • note that the largest frequency that a time series with time step of 1.0 can resolve is \(\omega=\frac{1}{2}\)

  • using Euler’s formular

\[e^{i\phi} = cos(\phi) + i~sin(\phi)\]
  • we can rewrite (33) to

\[\Gamma(\omega) = \gamma(0) + 2\sum_{\tau=1}^{\infty}\gamma(\tau)cos(2\pi\tau\omega)\]
  • the spectrum can be interpreted as the covariance between the auto-correlation function and the cosine function at different frequencies

  • characteristics of the spectrum:

    1. spectrum is continuous and differentiable in \([-\frac{1}{2},~\frac{1}{2}]\) (unlike the discrete time series or auto-correlation function)

    2. spectrum describes the distribution of variance across time scales: \(Var(\mathbf{X_t}) = \gamma(0) = 2\int_0^{\Sigma}\Gamma(\omega)~dw\), \(\Gamma(\omega)\) is a variance per frequency

    3. \(\gamma(\tau) = \int_{-\frac{1}{2}}^{\frac{1}{2}}\Gamma(\omega)e^{2\pi i\tau\omega}~dw\)

    4. \(\left. \frac{d}{d\omega}\Gamma(\omega)\right|_{\omega=0} =0\), spectrum must be flat at long time scales for stationary processes

    5. \(\Gamma(\omega) \propto \mathbf{X}^2~~\Rightarrow~~\sigma(\Gamma(\omega))\propto E(\Gamma(\omega))\), statistical uncertainty is proportional to expectation value (large for peaks), \(E(\mathbf{X})=k\), \(Var(\mathbf{X})=2k\) with k as the number of degrees of freedom

../_images/L16_1_SST.PNG
../_images/L16_2_spectrum_AR.PNG
  • for spectral variance per frequency on a log-linear scale as in the top row of Fig. 2, the area is the total variance. drawback: theoretical models follow power laws

../_images/L16_3_spectrum_elnino.PNG
  • peak at frequency of about \(0.25yr^{-1}\), note the different representations of that peak.

../_images/L16_4_spectrum_ts.PNG

Spectra of AR(p)-processes#

  • the spectrum of an AR(p)-process with process parameters \(\{\alpha_1,...,\alpha_p \}\) and noise variance \(Var(\mathbf{Z_t})= \sigma_Z^2\) is:

(34)#\[\Gamma(\omega) = \frac{\sigma_Z^2}{\left| 1-\sum_{k=1}^p \alpha_ke^{-2\pi ik\omega} \right|^2}\]
  • spectrum of a white noise process AR(0): variances are equally distributed on all frequencies:

\[\Gamma(\omega) = \sigma_Z^2\]

Spectrum of an AR(1)-process:#

  • for \(p=1\), (34) can be rewritten as:

(35)#\[\Gamma(\omega) = \frac{\sigma_Z^2}{\left| 1-\alpha_1e^{-2\pi i\omega} \right|^2} = \frac{\sigma_Z^2}{1 + \alpha_1^2-2\alpha_1cos(2\pi\omega)}\]
  • for small values \(\omega \in [0,~\frac{1}{2}]\) we can use the Taylor-expansion

(36)#\[cos(2\pi\omega) = 1- \frac{(2\pi\omega)^2}{2!} + ...\]
  • inserting the expansion (36) into (35) gives us the expression:

\[\Gamma(\omega) \approx \frac{\sigma_Z^2}{1 + \alpha_1^2-2\alpha_1(1-\frac{(2\pi\omega)^2}{2!})} = \frac{\sigma_Z^2}{1 + \alpha_1^2-2\alpha_1 +\alpha_1(2\pi\omega)^2} = \frac{c_1\sigma_Z^2}{c_2+\omega^2}\]
  • we can see that the spectrum follows a linear gradient of \(-2\) in log-log scale for \(\omega >>0\)

  • for \(\omega \in ]0,~\frac{1}{2}[\) the spectrum has no extrema because:

\[\frac{d}{d\omega}\Gamma(\omega) = -2\alpha_1\Gamma_1(\omega)^2sin(2\pi\omega) \neq 0\]

What’s \(\Gamma_1?\)

  • for \(\alpha_1>0\) we have a maximum (plateau) at \(\omega=0\), i.e. red noise processes!

../_images/L16_5_spectrum_AR1.PNG
../_images/L16_6_0.PNG
../_images/L16_6.PNG
  • fitting the AR(1)-process to a time series: AR(1)-processes or red noise are often chosen as a null hypothesis. AR(1)-processes are well defined by \(\sigma_{xy}\) and the lag-1 correlation because

\[\sigma_Z^2 = (1-\alpha_1^2)\sigma_{\mathbf{X_t}}^2 = (1-\rho_1^2)\sigma_{\mathbf{X_t}}^2\]
../_images/L16_7_temperatures_series.PNG
../_images/L16_8_spectrum_temp.PNG
../_images/L16_9_ts_24h.PNG
../_images/L16_10_spectrum_24h.PNG
../_images/L16_11_spectrum_AR1.PNG

Spectrum of an AR(2)-process:#

  • for \(p=2\), (34) can be rewritten as:

(37)#\[\begin{split}\Gamma(\omega) = \frac{\sigma_Z^2}{1+\alpha-1^2 + \alpha_2^2 -2g(\omega)},\\\end{split}\]
  • with:

\[g(\omega) = \alpha_1(1-\alpha_2)cos(2\pi\omega)+\alpha_2cos(4\pi\omega)\]
../_images/L16_12_0.PNG
../_images/L16_12_spectrum_AR2.PNG

Estimating the spectra (periodogram)#

  • let \(\{ \mathbf{X_1}, \mathbf{X_2}, ..., \mathbf{X_T} \}\) with t odd be a time series:

\[\mathbf{X_t} = A_0 + \sum_{k=1}^{\frac{T-1}{2}} \left( a_kcos\left(\frac{2\pi kt}{T}\right) + b_ksin\left(\frac{2\pi kt}{T}\right) \right)\]
  • remember that th spectrum is a continous function of frequency and for the periodogram \(a_k\) and \(b_k\) are discrete

\[\begin{split}\begin{align*} A_0 &= \hat{\mu} = \frac{1}{T} \sum_{t=1}^T\mathbf{X_t}\\ a_j &= \frac{2}{T} \sum_{t=1}^T\mathbf{X_t}cos(2\pi \omega_jt)\\ b_j &= \frac{2}{T} \sum_{t=1}^T\mathbf{X_t}sin(2\pi \omega_jt)\\ \omega_i &= \frac{i}{T},~~~i=1,2,...,\frac{T-1}{2} \end{align*}\end{split}\]
  • the variance is given as:

\[Var(\mathbf{X_t}) = \frac{1}{T} \sum_{t=1}^T(\mathbf{X_t}-\mathbf{\bar{X}})^2 = \frac{1}{2} \sum_{k=1}^{\frac{T-1}{2}}(a_k^2+b_k^2)\]
  • the elements \((a_k^2+b_k^2)\) are the periodogram of the finite time series \(\{ \mathbf{X_1}, \mathbf{X_2}, ..., \mathbf{X_T} \}\):

\[I_{Tk} = \frac{T}{4}(a_k^2+b_k^2)\]
  • bad characteristics:

    1. \(I_{Tk} \propto \xi^2(2)\) relatively wide distribution, strongly skewed to larger values with peak at zero

    2. Uncertainty of the spectrum coefficients is independent of the length of the time series

../_images/L16_13_spectrum_sine.PNG
  • best estimates:

    1. Divide the time series into \(m\) chunks of length \(M=\frac{T}{m}\).

    2. Compute a periodogram \(I_{Tk}^{(l)}\) for \(k=0,1,..,\frac{M}{2}\) from each chunk \(l=1,2,..,m\)

    3. Estimate the spectrum by averaging the periodograms:

\[\hat{\Gamma}(\omega_j)=\frac{1}{m}\sum_{l=1}^m I_{Tk}^{(l)},\]
  • estimator of the spectrum \(\propto \xi^2(2m)\). estimate at each frequency is representative of a special bandwidth of \(\propto \frac{1}{M}\)

../_images/L16_14_periodogram.PNG