1
$\begingroup$

I am using the following function to estimate the Gaussian window rolling average of my timeseries. Though it works great from small size averaging windows, it crushes (or gets extremely slow) for larger averaging windows.

def norm_factor_Gauss_window(s, dt):
    
    numer         = np.arange(-3*s, 3*s+dt, dt)
    multiplic_fac = np.exp(-(numer)**2/(2*s**2))
    norm_factor   = np.sum(multiplic_fac)
    window        = len(multiplic_fac)
    
    return window,  multiplic_fac, norm_factor

# Create dataframe for MRE
aa = np.sin(np.linspace(0,2*np.pi,1000000))+0.15*np.random.rand(1000000)
df = pd.DataFrame({'x':aa})

hmany  = 10
dt     = 1      # ['seconds']
s      = hmany*dt  # Define averaging window size ['s']

# Estimate multip factor, normalizatoon factor etc
window, multiplic_fac, norm_factor= norm_factor_Gauss_window(s, dt)

# averaged timeseries
res2 =(1/norm_factor)*df.x.rolling(window, center=True).apply(lambda x: (x * multiplic_fac).sum(), raw=True, engine='numba', engine_kwargs= {'nopython': True,  'parallel': True} , args=None, kwargs=None)

#Plot
plt.plot(df.x[0:2000])
plt.plot(res2[0:2000])

I am aware that people usually speed up moving average operations using convolve(e.g., https://stackoverflow.com/questions/14313510/how-to-calculate-rolling-moving-average-using-python-numpy-scipy)

Would it be possible to use convolve here somehow to fix this issue? Also, are there any other suggestion that would help me speed up the operation for large averaging windows?

$\endgroup$
7
  • 1
    $\begingroup$ yeah, you should just use convolve instead of implementing a convolution yourself as you seem to be doing in res2=…. $\endgroup$ Commented Sep 8, 2022 at 13:11
  • 1
    $\begingroup$ This looks mostly like a simple linear filter operation with a long-ish impulse response. By far the most effective way to implement this is Overlap Add or Overlap Save algorithms (leveraging the FFT). I would also take a look at the overhead of using DataFrame objects compared to simple arrays. $\endgroup$
    – Hilmar
    Commented Sep 8, 2022 at 13:29
  • $\begingroup$ @MarcusMüller would you mind adding an answer implementing this to the code. I don't think I really understnad how convolve works to implement it myself. Thank you for your time! $\endgroup$
    – Jokerp
    Commented Sep 8, 2022 at 15:35
  • $\begingroup$ @Hilmar would scipy.signal work for such a job? $\endgroup$
    – Jokerp
    Commented Sep 8, 2022 at 15:36
  • 1
    $\begingroup$ docs.scipy.org/doc/scipy/reference/generated/… should work, also according to signal;.convolve() will automatically select this if it's better $\endgroup$
    – Hilmar
    Commented Sep 8, 2022 at 16:25

1 Answer 1

1
$\begingroup$

With the help of @Hilmar an @Marcus I was able to drastically improve the speed of this code using the following:

from scipy import signal
def norm_factor_Gauss_window(s, dt):
    
    numer         = np.arange(-3*s, 3*s+dt, dt)
    multiplic_fac = np.exp(-(numer)**2/(2*s**2))
    norm_factor   = np.sum(multiplic_fac)
    window        = len(multiplic_fac)
    
    return window,  multiplic_fac, norm_factor

# Create dataframe for MRE
aa = np.sin(np.linspace(0,2*np.pi,1000000))+0.15*np.random.rand(1000000)
df = pd.DataFrame({'x':aa})

hmany  = 10
dt     = 1      # ['seconds']
s      = hmany*dt  # Define averaging window size ['s']

# Estimate multip factor, normalizatoon factor etc
window, multiplic_fac, norm_factor= norm_factor_Gauss_window(s, dt)

# averaged timeseries

res2 = (1/norm_factor)*signal.fftconvolve(df.x.values, multiplic_fac[::-1], 'same')

#Plot
plt.plot(df.x[0:2000])
plt.plot(res2[0:2000])
$\endgroup$
1
  • $\begingroup$ You should be able to give yourself the check mark now, so please do so! $\endgroup$
    – Peter K.
    Commented Oct 8, 2022 at 19:32

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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