[SciPy-User] Autocorrelation function: Convolution vs FFT

Skipper Seabold jsseabold@gmail....
Tue Jun 22 12:49:30 CDT 2010


I am trying to compute the autocorrelation via convolution and via fft
and am far from an expert in DSP.  I'm wondering if someone can spot
anything that might introduce numerical inaccuracies or if I'm stuck
with the following two being slightly different.

Generate some autocorrelated data:

import numpy as np
nobs = 150000
x = np.zeros((nobs))
for i in range(1,nobs):
    x[i] = .85 * x[i-1] + np.random.randn()

# compute ACF using convolution

x0 = x - x.mean()

# this takes a while for the big data
acf1 = np.correlate(x0,x0,'full')[nobs-1:]/nobs
acf1 /= acf1[0]


# compute ACF using FFT

Frf = np.fft.fft(x0, n=2*nobs) # zero-pad for separability
Sf = Frf * Frf.conjugate()
acf2 = np.fft.ifft(Sf)
acf2 = acf2[1:nobs+1]/nobs
acf2 /= acf2[0]
acf2 = acf2.real

np.linalg.norm(acf1-acf2, ord=2)

They are pretty close, but I would expect them to be closer than this.

np.max(acf1-acf2)
0.006581962491189159

np.min(acf1-ac2)
-0.0062705596399049799

Skipper


More information about the SciPy-User mailing list