0
$\begingroup$

My objective is to detect all kinds of seasonalities and their time periods that are present in a timeseries waveform.

I'm currently using the following dataset: https://www.kaggle.com/rakannimer/air-passengers

At the moment, I've tried the following approaches:

1) Use of FFT:

import pandas as pd
import numpy as np
from statsmodels.tsa.seasonal import seasonal_decompose
 
#https://www.kaggle.com/rakannimer/air-passengers
df=pd.read_csv('AirPassengers.csv')
 
df.head()

frequency_eval_max = 100
A_signal_rfft = scipy.fft.rfft(df['#Passengers'], n=frequency_eval_max)
n = np.shape(A_signal_rfft)[0] # np.size(t)
frequencies_rel = len(A_signal_fft)/frequency_eval_max * np.linspace(0,1,int(n))

fig=plt.figure(3, figsize=(15,6))
plt.clf()
plt.plot(frequencies_rel, np.abs(A_signal_rfft), lw=1.0, c='paleturquoise')
plt.stem(frequencies_rel, np.abs(A_signal_rfft))
plt.xlabel("frequency")
plt.ylabel("amplitude")

This results in the following plot: enter image description here

But it doesn't result in anything conclusive or comprehensible.

Ideally I wish to see the peaks representing daily, weekly, monthly and yearly seasonality.

Could anyone point out what am I doing wrong?

2) Autocorrelation:

from pandas.plotting import autocorrelation_plot
plt.rcParams.update({'figure.figsize':(10,6), 'figure.dpi':120})
autocorrelation_plot(df['#Passengers'].tolist())

After doing which I get a plot like the following: enter image description here

But how do I read this plot and how can I derive the presence of the various seasonalities and their periods from this?

3) SLT Decomposition Algorithm

df.set_index('Month',inplace=True)
df.index=pd.to_datetime(df.index)
#drop null values
df.dropna(inplace=True)
df.plot()

result=seasonal_decompose(df['#Passengers'], model='multiplicable', period=12)

result.seasonal.plot()

This gives the following plot: enter image description here

But here I can only see one kind of seasonality.

So how do we detect all the types of seasonalities and their time periods that are present using this method?


Hence, I've tried 3 different approaches but they seem either erroneous or incomplete.

Could anyone please help me out with the most effective approach (even apart from the ones I've tried) to detect all kinds of seasonalities and their time periods for any given timeseries data?



EDIT: As suggested by Max, I tried out FFT (with 500 samples transformed in FFT) & Autocorrelation on a different dataset having more no. of data points resulting in the following. (Note: this is a square wave having daily and weekly periodicity)

The following is my code for the same:

  1. Generation of the waveform:
import json
import sys, os
import numpy as np
import pandas as pd
import glob
import pickle

from statsmodels.tsa.stattools import adfuller, acf, pacf
from scipy.signal import find_peaks, square
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
import matplotlib.pyplot as plt

#GENERATION OF A FUNCTION WITH DUAL SEASONALITY & NOISE

def white_noise(mu, sigma, num_pts):
    """ Function to generate Gaussian Normal Noise
    Args:
        sigma: std value
        num_pts: no of points
        mu: mean value

    Returns:
        generated Gaussian Normal Noise
    """
    
    noise = np.random.normal(mu, sigma, num_pts)
    return noise

def signal_line_plot(input_signal: pd.Series, title: str = "", y_label: str = "Signal"):
    """ Function to plot a time series signal
    Args:
        input_signal: time series signal that you want to plot
        title: title on plot
        y_label: label of the signal being plotted
        
    Returns:
        signal plot
    """
    
    plt.plot(input_signal)
    plt.title(title)
    plt.ylabel(y_label)
    plt.show()

t_week = np.linspace(1,480, 480)
t_weekend=np.linspace(1,192,192)
T=96 #Time Period
x_weekday = 10*square(2*np.pi*t_week/T, duty=0.7)+10 + white_noise(0, 1,480)
x_weekend = 2*square(2*np.pi*t_weekend/T, duty=0.7)+2 + white_noise(0,1,192)
x_daily_weekly = np.concatenate((x_weekday, x_weekend)) 
x_daily_weekly_long = np.concatenate((x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly,x_daily_weekly))
signal_line_plot(x_daily_weekly_long)
signal_line_plot(x_daily_weekly_long[0:1000])

#x_daily_weekly_long is the final waveform
  1. FFT:
frequency_eval_max = 500
A_signal_rfft = scipy.fft.rfft(x_daily_weekly_long, n=frequency_eval_max)
n = np.shape(A_signal_rfft)[0] # np.size(t)
frequencies_rel = len(A_signal_fft)/frequency_eval_max * np.linspace(0,1,int(n))

fig=plt.figure(3, figsize=(15,6))
plt.clf()
plt.plot(frequencies_rel, np.abs(A_signal_rfft), lw=1.0, c='paleturquoise')
plt.stem(frequencies_rel, np.abs(A_signal_rfft))
plt.xlabel("frequency")
plt.ylabel("amplitude")

enter image description here 3) Autocorrelation:

from pandas.plotting import autocorrelation_plot
plt.rcParams.update({'figure.figsize':(10,6), 'figure.dpi':120})
autocorrelation_plot(x_daily_weekly_long)

enter image description here

$\endgroup$
1

1 Answer 1

1
$\begingroup$

Your dataset is sampled on a monthly basis. The shortest seasonality you will possibly be able to see is $2$ months. Daily, weekly and monthly are impossible.

While your FFT plot is somewhat clumsy to interpret due to the fact, that you choose to transform 100 samples and to scale the plot considering the DC part (the first FFT bin, $0\text{Hz}$), some things can be garnered from it. The peak around $0.25$ represents the period of $2$ years length, the one at about $0.5$ is the year long period. The peak around $0.95$ stands for $8$ months. The autocorrelation plot clearly depicts the year long period. That's where the peaks at multiples of 12 come from, so does the STL plot (as you concluded yourself. $8$ months is visible here, too.

EDIT

Your synthesized signal has a sampling rate of 96/day (one sample every quarter hour). The length of the signal is 6780 samples, which translates to 10 weeks. The maximum frequency represented by the FFT plot is 48/day. The frequency axis of the FFT plot has then a bin-to-bin distance of $$\frac{48}{(6780/2)}\cdot\frac{1}{\text{day}} \approx 0,014\cdot\frac{1}{\text{day}}$$

This means, that at sample $1/0,014\approx71$ there should be a peak representing the daily periodicity.

Usually, you would go the other way around: for each peak, you calculate the period its bin stands for. For example the high peak in the 5th bin stands for the two week period. $\frac{1}{0,014\cdot5}\approx 14$

Like this, you can interpret each peak in the FFT.

For the autocorrelation it is analog: you look for peaks and dips, and since it is in the time domain, you can read the period directly from the position. E.g. a peak at sample 96 represents a daily period, since the sampling frequency is 96/day.

If you want a somewhat smoother FFT, try sine signals instead of square waves, as these tend to be "noisy" in the frequency domain (a rect function transforms to an si function).

$\endgroup$
14
  • $\begingroup$ I'd like to know, how were you able to determine the time duration (of 2 years, 8 months, etc.) based on the frequency values from the FFT plot? Sure, I can increase the samples to say a higher value like 5000 or so - but upon doing so I'm getting the peaks at pretty much the same frequency values. So even after having more samples in the FFT, how would it be possible to determine the daily, weekly or monthly periodicities? $\endgroup$
    – EnigmAI
    Commented Mar 29, 2022 at 12:18
  • $\begingroup$ The data sample period is one month, so your sampling frequency is 1/month. You took 100 samples and calculated the FFT. That yields 50 frequency bins for the right side. (rfft) By Nyquist, the highest frequency component in your signal is half the sampling frequency, so in your case 0.5/month. Combining those numbers, your frequency bins are spaced by 0.005/month. To find the bin representing 12 months, you solve $n=1/(0,005 \cdot 12$ which yields (rounded) 17. So the 17th bin stands for a period of one year. More FFT samples will do no good. In the raw data is no information on days or weeks $\endgroup$
    – Max
    Commented Mar 29, 2022 at 12:32
  • $\begingroup$ Alright, I got the calculation. But coming to my original issue, is making use of an FFT the right way to calculate the periodicity of a timeseries data in the 1st place? Based on your calculation, it comes out to be the 200th (=1/0.005*1) bin representing 1 month, which isn't a part of this plot itself. So how do we identify the time-period from the FFT plot, which is 1 month for this data? $\endgroup$
    – EnigmAI
    Commented Mar 29, 2022 at 14:05
  • $\begingroup$ A period of one month cannot be evaluated. The shortest possible period is two months. Read about Nyquists sampling theorem. And yes, the FFT is an apt tool to identify periodic structures in data. Your data is just coarsely resolved. If you had data for each day, you could find periodicities as short as 2 days. $\endgroup$
    – Max
    Commented Mar 29, 2022 at 14:10
  • $\begingroup$ As you suggested, I tried out FFT and Autocorrelation on a different dataset having more no. of points. I've added its details in the question itself. Can you please now help me interpret the two plots? $\endgroup$
    – EnigmAI
    Commented Mar 29, 2022 at 18:53

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.