<jupyter_start><jupyter_text>As expected, there is a strong daily seasonality to this time series. The strongest correlations are 1 day lag and 7 day lag (weekly component)<jupyter_code>from statsmodels.tsa.seasonal import seasonal_decompose
decomposition = seasonal_decompose(midtc_df, freq=4*24) #daily  
fig = plt.figure()  
fig = decomposition.plot()  
fig.set_size_inches(15, 8)
#there is also a weekly component
#inspect residuals
plot_acf(decomposition.resid.dropna())

#too much autocorrelation in the residuals ; we have neglected the weekly component<jupyter_output><empty_output><jupyter_text>There is still a seasonal component left in the trend. Try weekly seasonality<jupyter_code>from statsmodels.tsa.seasonal import seasonal_decompose
weekly_decomposition = seasonal_decompose(midtc_df, freq=7*4*24) #weekly  
fig = plt.figure()  
fig = weekly_decomposition.plot()  
fig.set_size_inches(15, 8)

#weekly component 
#examine the residuals
print( "Peak at {}".format(50+np.argmax(acf(weekly_decomposition.resid.dropna(), nlags=8*24*4)[50:])))
#we are now neglecting daily seasonality
mt_pacf = pacf(midtc_df, nlags=8*24*4) #8 days

#significance level of pacf with N samples

N = len(midtc_df)
sig = (np.exp(2*1.96/np.sqrt(N-3)-1)/(np.exp(2*1.96/np.sqrt(N-3)+1)))

plt.scatter(np.arange(24*4+1), mt_pacf[:24*4+1])
plt.axhline(y=sig, color='r', linestyle='-')
plt.axhline(y=-sig, color='r', linestyle='-')
#find significant components
np.argwhere(abs(mt_pacf[:24*7]) > sig)
#%pdb on
# PACF indicates a model AR(1-4) with seasonality of day. Potential for MA component given damped sin effect
# ACF indicates additional seasonality component of week

from statsmodels.tsa.statespace.sarimax  import SARIMAX
#grid search using AIC
models = []
results = []

seasonal_options = [(1,0,0,96), (0,0,1,96), (1,0,1,96)] # order 7 (weekly) doesn't converge (singular matrix)

for p in range(1,5):
    for q in range(0,5):
        for s in range(3):
            print(p,q,s)
            try:
                models.append(SARIMAX(midtc_df.values, order=(p,0,q), seasonal_order=seasonal_options[s]))
                results.append(models[-1].fit())
            except Exception as e:
                print(e)
                results.append(None)

results[0].summary()
#pick best model based on AIC & Q stat
from matplotlib import pyplot as plt
aic = pd.Series([x.aic if x else 50000 for x in results])
aic = aic.fillna(50000)
np.argmin(aic)
ax = aic.plot()
ax.set_ylim((26000,30000))
# pick simplest model with AIC below 24000
best_order = (1,0,1)
best_seasonal = seasonal_options[2] # (1,0,1,96)<jupyter_output><empty_output>
StarcoderdataJupyterStructuredCleanDedup