In [334]:
%matplotlib inline
import warnings
warnings.filterwarnings('ignore')
import datetime
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import rcdefaults
import statsmodels.api as sm
from statsmodels.tsa.stattools import acf,pacf
from statsmodels.tsa.stattools import adfuller
from statsmodels.tsa.arima_model import ARMA
rcdefaults()
In [280]:
def test_stationarity(timeseries):
#Perform Dickey-Fuller test:
print ('Results of Dickey-Fuller Test:')
dftest = adfuller(timeseries, autolag='AIC')
dfoutput = pd.Series(dftest[0:4], index=['Test Statistic','p-value','#Lags Used','Number of Observations Used'])
for key,value in dftest[4].items():
dfoutput['Critical Value (%s)'%key] = value
print (dfoutput)
In [281]:
df_milk = pd.read_csv('milk_production.csv')
In [282]:
df_milk
Out[282]:
In [283]:
rng = pd.date_range(df_milk.Month[0], periods=len(df_milk), freq='M')
In [284]:
rng
Out[284]:
In [285]:
ts_milk = pd.Series(df_milk.MilkProduction.values, index=rng)
In [286]:
ts_milk
Out[286]:
In [287]:
plt.figure(figsize=(9,6))
plt.plot(ts_milk, color = 'blue')
plt.show()
In [288]:
res = sm.tsa.seasonal_decompose(ts_milk)
resplot = res.plot()
plt.show()
In [289]:
test_stationarity(ts_milk)
In [290]:
ts_milk_fd = ts_milk - ts_milk.shift()
ts_milk_fd = ts_milk_fd[np.isnan(ts_milk_fd) == False]
test_stationarity(ts_milk_fd)
In [291]:
ts_milk_log = np.log(ts_milk)
test_stationarity(ts_milk_log)
In [292]:
plt.figure(figsize=(9,6))
plt.plot(ts_milk_fd, color = 'blue')
plt.show()
In [293]:
res = sm.tsa.seasonal_decompose(ts_milk_fd)
resplot = res.plot()
plt.show()
In [294]:
f, (ax1, ax2, ax3) = plt.subplots(3,1,figsize=(24,18))
ax1.plot(ts_milk_fd)
sm.graphics.tsa.plot_acf(ts_milk_fd.values.squeeze(), lags=40, ax = ax2)
sm.graphics.tsa.plot_pacf(ts_milk_fd.values.squeeze(), lags=40, ax = ax3)
Out[294]:
In [295]:
ts_milk_fd_train = ts_milk_fd[:datetime.datetime(1975, 9, 30)]
In [296]:
ts_milk_fd_test = ts_milk_fd[datetime.datetime(1975, 10, 31):]
In [297]:
model = ARMA(ts_milk_fd_train, order=(2,0))
results_AR = model.fit()
forecast = pd.Series(results_AR.forecast(3)[0], index = ts_milk_fd_test.index)
In [298]:
plt.figure(figsize=(10,8))
plt.plot(ts_milk_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.show()
In [299]:
plt.figure(figsize=(10,8))
plt.plot(ts_milk_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.xlim(datetime.datetime(1974,1,1), datetime.datetime(1976,3,1))
plt.show()
In [337]:
model = ARMA(ts_milk_fd_train, order=(0,2))
results_AR = model.fit()
forecast = pd.Series(results_AR.forecast(3)[0], index = ts_milk_fd_test.index)
In [338]:
plt.figure(figsize=(10,8))
plt.plot(ts_milk_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.show()
In [339]:
plt.figure(figsize=(10,8))
plt.plot(ts_milk_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.xlim(datetime.datetime(1974,1,1), datetime.datetime(1976,3,1))
plt.show()
In [340]:
model = ARMA(ts_milk_fd_train, order=(4,4))
results_AR = model.fit()
forecast = pd.Series(results_AR.forecast(3)[0], index = ts_milk_fd_test.index)
In [341]:
plt.figure(figsize=(10,8))
plt.plot(ts_milk_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.show()
In [342]:
plt.figure(figsize=(10,8))
plt.plot(ts_milk_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.xlim(datetime.datetime(1974,1,1), datetime.datetime(1976,3,1))
plt.show()
In [306]:
plt.figure(figsize=(10,8))
resid = ts_milk_fd_train - results_AR.fittedvalues
plt.plot(resid)
plt.show()
In [307]:
f, (ax1, ax2, ax3) = plt.subplots(3,1,figsize=(24,18))
ax1.plot(resid)
sm.graphics.tsa.plot_acf(resid.values.squeeze(), lags=40, ax = ax2)
sm.graphics.tsa.plot_pacf(resid.values.squeeze(), lags=40, ax = ax3)
Out[307]:
In [308]:
df_car = pd.read_csv('carsales.csv')
In [309]:
df_car
Out[309]:
In [310]:
rng = pd.date_range(df_car.Month[0], periods=len(df_car), freq='M')
In [311]:
rng
Out[311]:
In [312]:
ts_car = pd.Series(df_car.CarSales.values, index=rng)
In [313]:
ts_car
Out[313]:
In [343]:
plt.figure(figsize=(9,6))
plt.plot(ts_car, color = 'blue')
plt.show()
In [315]:
res = sm.tsa.seasonal_decompose(ts_car)
resplot = res.plot()
plt.show()
In [316]:
test_stationarity(ts_car)
In [317]:
ts_car_fd = ts_car - ts_car.shift()
ts_car_fd = ts_car_fd[np.isnan(ts_car_fd) == False]
test_stationarity(ts_car_fd)
In [318]:
plt.figure(figsize=(9,6))
plt.plot(ts_car_fd, color = 'blue')
plt.show()
In [319]:
res = sm.tsa.seasonal_decompose(ts_car_fd)
resplot = res.plot()
plt.show()
In [320]:
f, (ax1, ax2, ax3) = plt.subplots(3,1,figsize=(24,18))
ax1.plot(ts_car_fd)
sm.graphics.tsa.plot_acf(ts_car_fd.values.squeeze(), lags=40, ax = ax2)
sm.graphics.tsa.plot_pacf(ts_car_fd.values.squeeze(), lags=40, ax = ax3)
Out[320]:
In [321]:
ts_car_fd_train = ts_car_fd[:datetime.datetime(1968, 9, 30)]
In [322]:
ts_car_fd_test = ts_car_fd[datetime.datetime(1968, 10, 31):]
In [323]:
model = ARMA(ts_car_fd_train, order=(4,0))
results_AR = model.fit()
forecast = pd.Series(results_AR.forecast(3)[0], index = ts_car_fd_test.index)
In [324]:
plt.figure(figsize=(10,8))
plt.plot(ts_car_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.show()
In [325]:
plt.figure(figsize=(10,8))
plt.plot(ts_car_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.xlim(datetime.datetime(1968,1,1), datetime.datetime(1969,3,1))
plt.show()
In [326]:
model = ARMA(ts_car_fd_train, order=(0,4))
results_AR = model.fit()
forecast = pd.Series(results_AR.forecast(3)[0], index = ts_car_fd_test.index)
In [327]:
plt.figure(figsize=(10,8))
plt.plot(ts_car_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.show()
In [328]:
plt.figure(figsize=(10,8))
plt.plot(ts_car_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.xlim(datetime.datetime(1968,1,1), datetime.datetime(1969,3,1))
plt.show()
In [329]:
model = ARMA(ts_car_fd_train, order=(3,3))
results_AR = model.fit()
forecast = pd.Series(results_AR.forecast(3)[0], index = ts_car_fd_test.index)
In [330]:
plt.figure(figsize=(10,8))
plt.plot(ts_car_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.show()
In [331]:
plt.figure(figsize=(10,8))
plt.plot(ts_car_fd)
plt.plot(results_AR.fittedvalues, color='red')
plt.plot(forecast, color = 'green')
plt.xlim(datetime.datetime(1968,1,1), datetime.datetime(1969,3,1))
plt.show()
In [332]:
plt.figure(figsize=(10,8))
resid = ts_car_fd_train - results_AR.fittedvalues
plt.plot(resid)
plt.show()
In [333]:
f, (ax1, ax2, ax3) = plt.subplots(3,1,figsize=(24,18))
ax1.plot(resid)
sm.graphics.tsa.plot_acf(resid.values.squeeze(), lags=40, ax = ax2)
sm.graphics.tsa.plot_pacf(resid.values.squeeze(), lags=40, ax = ax3)
Out[333]:
In [ ]: