8
$\begingroup$

In the Surrogate Time Series (Schreiber, Schmitz) paper, the authors claim that surrogates for a second order stationary time series can be generated by taking the Fourier Transform of the series, multiplying random phases to the coefficients, and then transforming back. This procedure should preserve the autocorrelation function. I am trying to implement this in python using numpy.

import numpy
import pandas
import matplotlib.pyplot as plt

numpy.random.seed(0)
ts = numpy.random.normal(0, 1, 1000)

pandas.tools.plotting.autocorrelation_plot(ts)
plt.ylim([-0.1,0.1])
plt.title('Autocorrelation function of random time series')

Random time series autocorrelation

The attempt at the mentioned procedure:

ts_fourier = numpy.fft.fft(ts)
ts_fourier_new = [x*numpy.exp(numpy.random.uniform(0,2*numpy.pi)*1.0j) for x in ts_fourier]
new_ts = numpy.fft.ifft(ts_fourier_new)

pandas.tools.plotting.autocorrelation_plot(new_ts)
plt.ylim([-0.1,0.1])
plt.title('Autocorrelation function of surrogate time series')

Surrogate Autocorrelation

As can be seen, the autocorrelation function of the surrogate is not identical to the original time series. What is the mistake in my implementation?

$\endgroup$
1
  • 1
    $\begingroup$ Tip: You can avoid using Python loops (which cost time) in the phase shuffling by using Numpy’s array arithmetics: Just replace the respective line with ts_fourier_new = numpy.exp(numpy.random.uniform(0,numpy.pi,1000)*1.0j)*ts_fourier . $\endgroup$
    – Wrzlprmft
    Mar 28, 2016 at 14:43

1 Answer 1

5
$\begingroup$

You need to use the Fourier transform (and inverse transform) for real time series, i.e., rfft and irfft, respectively. This way you ensure that your surrogate is real. You can do this by replacing the respective lines of your code with the following:

ts_fourier  = numpy.fft.rfft(ts)
random_phases = numpy.exp(numpy.random.uniform(0,numpy.pi,len(ts)/2+1)*1.0j)
ts_fourier_new = ts_fourier*random_phases
new_ts = numpy.fft.irfft(ts_fourier_new)

(Alternatively, you can take the regular Fourier transform and shuffle the phases in an antisymmetric manner.)

The remaining deviations originate from the finiteness of your time series: The identity of the autocorrelation functions is based on the fact that the original time series and the surrogate have per construction the same power spectrum, which in turn is linked to the autocorrelation function via the Wiener–Khinchin theorem. However the latter only holds in approximation for finite time series as it is actually about processes.

This effect is somewhat alleviated if you look at time series with “actual” prominent frequency components and not just noise.

$\endgroup$
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.