findfrequency spec.ar equivalent in Python

Question

There is a very useful function in R called findfrequency on the forecast package that returns the period of the dominant frequency of a time series. More info on the function from the author can be found here: https://robjhyndman.com/hyndsight/tscharacteristics/

I want to implement something equivalent in Python and I am having trouble with the functions that should be equal to the spec.ar R function that is inside findfrequency.

The function starts from detrending the series which is easily done with x = statsmodels.tsa.tsatools.detrend(myTs, order=1, axis=0). Now that I have the residuals I would like to do in Python the equivalent of the spec.ar function in R that first fits an AR model to x (or uses the existing fit) and computes (and by default plots) the spectral density of the fitted model.

I have not found anything similar so I am doig each step at a time, first the AR and then the spec estimation. I am using the Airpassengers time series and I am not able to get the same results on R and Python for the AR order or coefficients.

My R code:

x <- AirPassengers
x <- residuals(tslm(x ~ trend))
ARmodel <- ar(x)
ARmodel

I get that 15 is the selected order for my autoregressive model.

My Python Code:

import statsmodels.api as sm
dataPeriodic = pd.read_csv('AirPassengers.csv')
tsPeriodic = dataPeriodic.iloc[:,1]
x = statsmodels.tsa.tsatools.detrend(tsPeriodic, order=1, axis=0)
n = x.shape[0]
est_order = sm.tsa.AR(x).select_order(maxlag=20, ic='aic', trend='nc')
print(est_order)

Here I get a very different result with an order selected that is equal to 10 instead of 15 and I have to specify the upper limit of the lag search with the maxlag parameter..

I have tried with the tsa.AutoReg without success, I get another different result.

So, is there a way to fit an AR model in the same way that R does ? Something similiar to spec.ar or even something similar to the findfrequency function ? I am quite confused by the big diferences the 'same' methods can output in the two languages.

Aristophanes

Hi, did you ever found a way to make it work for python? – Hotfoot

@Hotfoot so its given that Python doesn't have a simple function like R to do this, right? – Angary

Exactly, regardless i tried the package in R and while it works for some timeseries, i got false positives for some timeseries that clearly didn't have a frequency. Very hard to find a reliable way to do this programmatically. Simple Autocorrelation it's not enough either, please let me know if you stumble upon something useful – Hotfoot

@Hotfoot I agree. I ended up using Pandas as described in my answer. – Angary

Answer 1

Closest I could find in Python for findfrequency of the R forecast package was by using pandas.infer_freq like this:

>>> import pandas as pd
>>> ts_data = pd.read_csv("ts_data.csv")
>>> pd.infer_freq(ts_data.index.values)
4

Angary