import requests
import pandas as pd
import math
import datetime
import urllib.request
import json
import time
import numpy as np
import matplotlib.pyplot as plt
from prophet import Prophet
import seaborn as sns

from sklearn.metrics import mean_squared_error, mean_absolute_error

pd.options.mode.chained_assignment = None  # default='warn'

import warnings
warnings.filterwarnings('ignore')

plt.style.use('ggplot')
plt.style.use('fivethirtyeight')

#hide api ket in text file so its not public
app_token = open('api_key.txt', 'r').read()
#app_token

#hide api token
#use request to pull 311 data from Jan 01, 2020 until today / can slice data in pandas later
limit = 500000
app_token = open('api_key.txt', 'r').read()
uri = f"https://data.buffalony.gov/resource/whkc-e5vr.json?$limit={limit}&$$app_token={app_token}&$where=open_date>'2020-01-10T12:00:00'"
r = requests.get(uri)
print('Status code ',r.status_code)
print('Number of rows returned ',len(r.json()))
print('Endoced URI with params ',r.url)
new_json = r.json()
#new_json

#make pandas df from  json / check head
df=pd.DataFrame(new_json)
print(df.shape)
df.head(2)

#save original Buffalo 311 Data to csv
df.to_csv('Buffalo311Data.csv', index = None, header=True)

#add count column with a value of one to rows to count calls by frequency
df['count'] = 1
#df.head()

#format open date to pd datetime to get ready to parse column
df['time'] = pd.to_datetime(df['open_date'])
#df.info() #check to make sure Dtype is datetime64

#take a look at teh columns I need
#parsing by hour will provide 24 times more training data than doing it by day
df1 = df[['time', 'count']]
df1.head(3)

#use date time to parse by year
df1['time'] = pd.to_datetime(df1['time']).dt.strftime('%Y-%m-%d %H')
df1.head()

#group by hour, count and reset index
df2 = df1.groupby(['time']).count().reset_index()

df2.info()

#make sure time is pd datetime after groupby. seems to create new object after groupby
df2['time'] = pd.to_datetime(df2['time'])

#set time as index
df2 = df2.set_index('time')

#df2.head()

#make df by day instead of hour for some charts
#ultimatly decided to use hourly because I feel it would help for staffing and letting people go home early 
df_day = df1.copy()

df_day['time'] = pd.to_datetime(df_day['time']).dt.strftime('%Y-%m-%d')
df_day = df_day.groupby(['time']).count().reset_index()
df_day['time'] = pd.to_datetime(df_day['time'])
df_day = df_day.set_index('time')

#df_day.head()

#make graph of calls by day
color_pal = sns.color_palette()
df_day.plot(style='.',
         figsize=(10,5),
         ms=1,
         color = color_pal[0],
         title='311 Calls By Day')
plt.show()

#make graph of call by hour
color_pal = sns.color_palette()
df2.plot(style='.',
         figsize=(10,5),
         ms=1,
         color = color_pal[0],
         title='311 Calls By Hour')
plt.show()

#create seasons for boxplot graph showing dayofweek and season
from pandas.api.types import CategoricalDtype

cat_type = CategoricalDtype(categories=['Monday', 'Tuesday',
                                       'Wednesday', 'Thursday',
                                       'Friday', 'Saturday', 'Sunday'],
                           ordered=True)




def create_features(df, label=None):
    """
    Creates time series features from datetime index.
    """
    df = df.copy()
    df['date'] = df.index
    df['dayofweek'] = df['date'].dt.dayofweek
    df['weekday'] = df['date'].dt.day_name()
    df['weekday'] = df['weekday'].astype(cat_type)
    df['quarter'] = df['date'].dt.quarter
    df['month'] = df['date'].dt.month
    df['year'] = df['date'].dt.year
    df['dayofyear'] = df['date'].dt.dayofyear
    df['dayofmonth'] = df['date'].dt.day
    df['weekofyear'] = df['date'].dt.isocalendar().week
    #df['weekofyear'] = df['date'].dt.weekofyear
    df['date_offset'] = (df.date.dt.month*100 + df.date.dt.day - 320)%1300
    
    df['season'] = pd.cut(df['date_offset'], [0, 300, 600, 900, 1300],
                         labels=['Spring', 'Summer', 'Fall', 'Winter'])


    X = df[['dayofweek','quarter','month','year',
           'dayofyear','dayofmonth','weekofyear', 'weekday', 'season']]
    if label:
        y = df[label]
        return X, y
    return X

X, y = create_features(df_day, label='count')

features_and_target = pd.concat([X, y], axis=1)

features_and_target.head()

fig, ax = plt.subplots(figsize=(10,5))
sns.boxplot(data=features_and_target.dropna(),
           x='weekday',
           y='count',
           hue='season',
           ax = ax,
           linewidth=1,
           palette='YlGnBu')

ax.set_title('Number of 311 Calls by Day of Week')
ax.set_xlabel('Day of Week')
ax.set_ylabel('311 Call Count')
ax.legend(bbox_to_anchor=(1,1))
plt.show()

fig, ax = plt.subplots(figsize=(10,5))
sns.boxplot(data=features_and_target.dropna(),
           x='weekofyear',
           y='count',
           hue='season',
           ax = ax,
           linewidth=1,
           palette='YlGnBu')

ax.set_title('Number Of 311 Calls By Week Of Year')
ax.set_xlabel('Week')
ax.set_ylabel('311 Call Count')
plt.xticks(rotation = 90) # rotates x-axis by 45 degrees
ax.legend(bbox_to_anchor=(1,1))
plt.show()

split_date = '1-jan-2023'
df_train = df2.loc[df2.index <= split_date].copy()
df_test = df2.loc[df2.index > split_date].copy()

# plot train and test so you can see where the split is

df_test \
    .rename(columns={'count': 'TEST SET'}) \
    .join(df_train.rename(columns={'count': 'TRAINING SET'}),
         how='outer') \
    .plot(figsize=(15,5), title='Buffalo 311 Calls By Hour', style='.')
plt.show()

#format data for prophet model using ds and y
df_train_prophet = df_train.reset_index()\
    .rename(columns={'time': 'ds',
                    'count': 'y'})

df_train_prophet.head()

%%time
model = Prophet()
model.fit(df_train_prophet)

%%time
# create test dat with same column names
df_test_prophet = df_test.reset_index()\
    .rename(columns={'time': 'ds',
                    'count': 'y'})

df_test_fcst = model.predict(df_test_prophet)

df_test_fcst.head()

fig, ax = plt.subplots(figsize=(10,5))
fig = model.plot(df_test_fcst, ax= ax)
plt.xticks(rotation = 45) # rotates x-axis by 45 degrees
ax.set_title('Prophet Forecast')
plt.show()

fig = model.plot_components(df_test_fcst)
plt.show()

#plot teh forecast with the actual
f, ax = plt.subplots(figsize=(15,5))
ax.scatter(df_test.index, df_test['count'], color ='r')
fig = model.plot(df_test_fcst, ax= ax)

#zoomed in plot of forecast vs actual
import datetime

f, ax = plt.subplots(figsize=(10,5))
ax.scatter(df_test.index, df_test['count'], color ='r')
fig = model.plot(df_test_fcst, ax= ax)
# ax.set_xbound(lower='01-01-2015',
#              upper='02-01-2015')
ax.set_xlim([datetime.date(2023, 1, 1), datetime.date(2023, 2, 1)])
ax.set_ylim(0,110)
plot = plt.suptitle('January 2023 Forecast vs Actuals')

#zoomed in plot of forecast vs actual
import datetime

f, ax = plt.subplots(figsize=(10,5))
ax.scatter(df_test.index, df_test['count'], color ='r')
fig = model.plot(df_test_fcst, ax= ax)
# ax.set_xbound(lower='01-01-2015',
#              upper='02-01-2015')
ax.set_xlim([datetime.date(2023, 1, 1), datetime.date(2023, 1, 7)])
ax.set_ylim(0,100)
plot = plt.suptitle('First Week Of January 2023 Forecast vs Actuals')

#zoomed in plot of forecast vs actual
import datetime

f, ax = plt.subplots(figsize=(10,5))
ax.scatter(df_test.index, df_test['count'], color ='r')
fig = model.plot(df_test_fcst, ax= ax)
# ax.set_xbound(lower='01-01-2015',
#              upper='02-01-2015')
ax.set_xlim([datetime.date(2023, 4, 1), datetime.date(2023, 4, 7)])
ax.set_ylim(0,70)
plot = plt.suptitle('First Week Of April 2023 Forecast vs Actuals')

# true value vs predicted value
np.sqrt(mean_squared_error(y_true=df_test['count'],
                          y_pred=df_test_fcst['yhat']))

mean_absolute_error(y_true=df_test['count'],
                   y_pred=df_test_fcst['yhat'])

# set up holiday dataframe
from pandas.tseries.holiday import USFederalHolidayCalendar as calendar

cal = calendar()


holidays = cal.holidays(start=df2.index.min(),
                        end=df2.index.max(),
                        return_name=True)
holiday_df = pd.DataFrame(data=holidays,
                          columns=['holiday'])
holiday_df = holiday_df.reset_index().rename(columns={'index':'ds'})

holiday_df.head()

holiday_df['holiday'].value_counts()

import time

#run Prophet model with holiday df as a parameter for holiday
#%%time

model_with_holidays = Prophet(holidays=holiday_df)
model_with_holidays.fit(df_train_prophet)

# predict on training set with model
df_test_fcst_with_hols = model_with_holidays.predict(df=df_test_prophet)

df_test_fcst_with_hols.head()

model_with_holidays.plot_components(
    df_test_fcst_with_hols)

#zoomed in plot of forecast vs actual

# red dots show actual vs projected line in blue
import datetime

f, ax = plt.subplots(figsize=(10,5))
ax.scatter(df_test.index, df_test['count'], color ='r')
fig = model.plot(df_test_fcst_with_hols, ax= ax)
# ax.set_xbound(lower='07-01-2015',
#              upper='07-07-2015')
ax.set_xlim([datetime.date(2023, 1, 1), datetime.date(2023, 5, 1)])
ax.set_ylim(0,110)
plt.xticks(rotation = 45) # rotates x-axis by 45 degrees
plot = plt.suptitle('January to April Forecast vs Actuals')

#zoomed in plot of forecast vs actual
import datetime

f, ax = plt.subplots(figsize=(10,5))
ax.scatter(df_test.index, df_test['count'], color ='r')
fig = model.plot(df_test_fcst_with_hols, ax= ax)
# ax.set_xbound(lower='07-01-2015',
#              upper='07-07-2015')
ax.set_xlim([datetime.date(2023, 4, 1), datetime.date(2023, 4, 7)])
ax.set_ylim(0,80)
plot = plt.suptitle('First Week Of April Forecast vs Actuals With Holidays')

# true value vs predicted value
np.sqrt(mean_squared_error(y_true=df_test['count'],
                          y_pred=df_test_fcst_with_hols['yhat']))

mean_absolute_error(y_true=df_test['count'],
                   y_pred=df_test_fcst_with_hols['yhat'])

future = model.make_future_dataframe(periods=365*24, freq='h', include_history=False)

forecast = model_with_holidays.predict(future)

forecast.head()

fig, ax = plt.subplots(figsize=(10,5))
fig = model.plot(forecast, ax= ax)
ax.set_title('Prophet Forecast')
plt.show()

#zoomed in plot of forecast
import datetime

fig, ax = plt.subplots(figsize=(10,5))
fig = model.plot(forecast, ax= ax)
ax.set_xlim([datetime.date(2023, 1, 1), datetime.date(2023, 12, 31)])
ax.set_ylim(0,100)
plt.xticks(rotation = 45) # rotates x-axis by 45 degrees
ax.set_title('Prophet Forecast')
plt.show()

#zoomed in plot of forecast for upcoming week; first week of may
import datetime

fig, ax = plt.subplots(figsize=(10,5))
fig = model.plot(forecast, ax= ax)
ax.set_xlim([datetime.date(2023, 6, 1), datetime.date(2023, 6, 7)])
ax.set_ylim(0,80)
ax.set_title('Prophet Forecast First Week Of June')
plt.show()

#zoomed in plot of forecast for upcoming Friday; first firday of june
import datetime

fig, ax = plt.subplots(figsize=(10,5))
fig = model.plot(forecast, ax= ax)
ax.set_xlim([datetime.date(2023, 6, 2), datetime.date(2023, 6, 3)])
ax.set_ylim(0,70)
ax.set_title('Prophet Forecast First Friday Of June 2023 By Hour')
plt.show()

split_date = '1-jan-2023'
df_train2 = df2.loc[df2.index <= split_date].copy()
df_test2 = df2.loc[df2.index > split_date].copy()

fig, ax = plt.subplots(figsize=(15,5))
df_train2.plot(ax=ax, label='Training Data', title='Train/Test Data Split')
df_test2.plot(ax=ax, label='Testing Data')
ax.axvline(split_date, color='black', ls='--')
ax.legend(['Training Data', 'Test Data'])
plt.show()

#check out two weeks in Jan 2022 to get a better picture of the data
df_train2.loc[(df_train.index > '01-01-2022') & (df_train2.index < '01-15-2022')].plot(figsize=(15,5), title='Two Weeks in Jan 2022')
plt.show()

def create_features_2(df):
    '''
    create time series features based on index date
    '''
    df = df.copy()
    df['hour'] = df.index.hour
    df['dayofweek'] = df.index.dayofweek
    df['quarter'] = df.index.quarter
    df['month'] = df.index.month
    df['year'] = df.index.year
    df['dayofyear'] = df.index.dayofyear
    return df

df3 = create_features_2(df2)

df3.head(2)
df3.shape

fig, ax = plt.subplots(figsize=(10,8))
sns.boxplot(data = df3, x='hour', y='count')
ax.set_title('Buffalo 311 Call Count By Hour')
plt.show()

fig, ax = plt.subplots(figsize=(10,8))
sns.boxplot(data = df3, x='month', y='count', palette='Blues')
ax.set_title('Buffalo 311 Call Count By Month')
plt.show()

fig, ax = plt.subplots(figsize=(10,8))
sns.boxplot(data = df3, x='quarter', y='count', palette='Blues')
ax.set_title('Buffalo 311 Call Count By Quarter')
plt.show()

#note that  2023 data is only until 05/04/2023
fig, ax = plt.subplots(figsize=(10,8))
sns.boxplot(data = df3, x='year', y='count', palette = 'Reds')
ax.set_title('Buffalo 311 Call Count By Year')
plt.show()

fig, ax = plt.subplots(figsize=(10,8))
sns.boxplot(data = df3, x='dayofweek', y='count', palette='Greens')
ax.set_title('Buffalo 311 Call Count By Day Of Week')
plt.show()

fig, ax = plt.subplots(figsize=(10,8))
sns.boxplot(data = df3, x='dayofyear', y='count', palette='Purples')
ax.set_title('Buffalo 311 Call Count By Day Of Year')
plt.show()

import xgboost as xgb
from sklearn.metrics import mean_squared_error

df2.head()

df_train = create_features_2(df_train)
df_test = create_features_2(df_test)

print(df_train.columns)
print(df_test.columns)

features = ['hour', 'dayofweek', 'quarter', 'month', 'year', 'dayofyear']
target = 'count'

X_train = df_train[features]
y_train = df_train[target]

X_test = df_test[features]
y_test = df_test[target]

#had to lower learning rate to yield better MSE results and have the model not overfit
reg = xgb.XGBRegressor(n_estimators = 1000, early_stopping_rounds = 50,
                      learning_rate = 0.001)
reg.fit(X_train, y_train,
       eval_set = [(X_train, y_train), (X_test, y_test)],
       verbose = 75)

from xgboost import plot_importance
plot_importance(reg)
plt.show()

df_test['prediction'] = reg.predict(X_test)

df3 = df3.merge(df_test[['prediction']], how='left', left_index=True, right_index=True)

ax = df3[['count']].plot(figsize = (15,5))
df3['prediction'].plot(ax=ax, style='.')
plt.legend(['True Data', 'Predictions'])
ax.set_title('Raw Data & Predictions')
plt.show()

ax = df_test.loc[(df_test.index > '02-01-2023') & (df_test.index < '03-01-2023')]['count'] \
    .plot(figsize=(15,5), title = ' Feb 2023')
df_test.loc[(df_test.index > '02-01-2023') & (df_test.index <'03-01-2023')]['prediction'] \
    .plot(style='.')
plt.legend(['Prediction','Real Data' ])
plt.show()

score = np.sqrt(mean_squared_error(df_test['count'], df_test['prediction']))
print(f'RMSE Score on Test set: {score:0.2f}')

df_test['error'] = np.abs(df_test[target] - df_test['prediction'])

df_test['date'] = df_test.index.date

#worst projected days happended close to the blizzard
df_test.groupby('date')['error'].mean().sort_values(ascending=False).head(5)

#best projected days of forecast
df_test.groupby('date')['error'].mean().sort_values(ascending=True).head(5)

df2['count'].plot(kind='hist', bins=1000)

df2.query('count < 10').plot(figsize=(15,5), style='.')

df2.query('count > 50').plot(figsize=(15,5), style='.')

split_date = '1-jan-2023'
df_train3 = df2.loc[df2.index <= split_date].copy()
df_test3 = df2.loc[df2.index > split_date].copy()

fig, ax = plt.subplots(figsize=(15,5))
df_train3.plot(ax=ax, label='Training Data', title='Train/Test Data Split')
df_test3.plot(ax=ax, label='Testing Data')
ax.axvline(split_date, color='black', ls='--')
ax.legend(['Training Data', 'Test Data'])
plt.show()

from sklearn.model_selection import TimeSeriesSplit

df2.shape

print(df2.head())
print(df2.tail())

#tss = TimeSeriesSplit(n_splits=2, test_size=24*365, gap = 24)
tss = TimeSeriesSplit(n_splits=8, test_size=24*60*1, gap = 24)
df4 = df2.sort_index()

df4.head()

fig, axs = plt.subplots(8, 1, figsize=(15,15),
                            sharex = True)

fold = 0

for train_idx, val_idx in tss.split(df4):
    train = df4.iloc[train_idx]
    test = df4.iloc[val_idx]
    train['count'].plot(ax=axs[fold],
                        label='Training Set', 
                        title=f'Data Train/Test Split Fold {fold}')
    test['count'].plot(ax=axs[fold],
                       label = 'Test Set')
    axs[fold].axvline(test.index.min(), color='black', ls='--')
    fold +=1

# make function for new time features 
def create_features4(df):
    """
    Create time series features based on time series index.
    """
    df = df.copy()
    df['hour'] = df.index.hour
    df['dayofweek'] = df.index.dayofweek
    df['quarter'] = df.index.quarter
    df['month'] = df.index.month
    df['year'] = df.index.year
    df['dayofyear'] = df.index.dayofyear
    df['dayofmonth'] = df.index.day
    df['weekofyear'] = df.index.isocalendar().week
    return df

df4 = create_features4(df4)

df4.head()

# make lag feature function
def add_lags(df):
    target_map = df['count'].to_dict()
    df['lag1'] = (df.index - pd.Timedelta('182 days')).map(target_map)
    df['lag2'] = (df.index - pd.Timedelta('364 days')).map(target_map)
    df['lag3'] = (df.index - pd.Timedelta('728 days')).map(target_map)
    return df

#test lags function
#df4 = add_lags(df)

#check df after lags function test
#df4.tail()

tss = TimeSeriesSplit(n_splits=8, test_size=24*60*1, gap = 24)
df4 = df2.sort_index()

fold = 0
preds = []
scores = []
for train_idx, val_idx in tss.split(df4):
    train = df4.iloc[train_idx]
    test = df4.iloc[val_idx]
    
    train = create_features4(train)
    test = create_features4(test)
    
    train = add_lags(train)
    test = add_lags(test)
    
    FEATURES4 = ['dayofyear', 'hour', 'dayofweek', 'quarter', 'month','year',
                'lag1','lag2', 'lag3']
    TARGET4 = 'count'
    
    X_train = train[FEATURES4]
    y_train = train[TARGET4]

    X_test = test[FEATURES4]
    y_test = test[TARGET4]
    
    
    reg = xgb.XGBRegressor(base_score=0.5, booster='gbtree',
                           n_estimators=100000,
                           early_stopping_rounds=50,
                           objective='reg:linear',
                           max_depth=3,
                           learning_rate=0.1)
    reg.fit(X_train, y_train,
            eval_set=[(X_train, y_train), (X_test, y_test)],
            verbose=10000)

    y_pred = reg.predict(X_test)
    preds.append(y_pred)
    score = np.sqrt(mean_squared_error(y_test, y_pred))
    scores.append(score)
    
    

print(f'Scores across all folds {np.mean(scores):0.4f}')
print(f'Fold scores:{scores}')

#retrain all data to prepare for forecasting
#still leverage all data for forecast prediction
#switchedn_estimator to 1500 because the model still seemed like it was getting better
#also moved the learning rate down

df4 = create_features4(df4)
df4 = add_lags(df4)

FEATURES4 = ['dayofyear', 'hour', 'dayofweek', 'quarter', 'month','year',
                'lag1','lag2', 'lag3']
TARGET4 = 'count'

X_all = df4[FEATURES4]
y_all = df4[TARGET4]

reg = xgb.XGBRegressor(base_score=0.5,
                       booster='gbtree',
                       #n_estimators=25000,
                       n_estimators=100000,
                       early_stopping_rounds=50,
                       objective='reg:linear',
                       max_depth=3,
                       #learning_rate=0.0001) # LEARNING FASTER SEEMS TO BOOST PERFORMANCE OF THE MODEL
                       learning_rate=0.1)

reg.fit(X_all, y_all,
       eval_set=[(X_all, y_all)],
       verbose=10000)

#find max value and set rage from last date 
FurtueStartDate = df4.index.max()
'''
df_day['time'] = pd.to_datetime(df_day['time']).dt.strftime('%Y-%m-%d')
'''
#FurtueStartDate = FurtueStartDate.dt.strftime('%Y-%m-%d')
FurtueStartDate

#can't project longer than the smallest lag
future = pd.date_range('2023-05-04', '11/02/2023' , freq='1h')
future_df = pd.DataFrame(index=future)

future_df['is_future'] = True
df4['is_future'] = False
df_and_future = pd.concat([df4, future_df])

df_and_future = create_features4(df_and_future)
df_and_future = add_lags(df_and_future)

df_and_future

#may need to go in and fill in lag Nan's with ffill
future_with_features = df_and_future.query('is_future').copy()
future_with_features

FEATURES4 = ['dayofyear', 'hour', 'dayofweek', 'quarter', 'month','year',
                'lag1','lag2', 'lag3']
TARGET4 = 'count'

#run trained model on features in future df
#save in new column called pred
future_with_features['pred'] = reg.predict(future_with_features[FEATURES4])

import matplotlib.pyplot as plt

future_with_features['pred'].plot(figsize = (10,5),
                                  #color= 'Blue',
                                  color= color_pal[0],
                                  ms=1, 
                                  lw=1,
                                  title='XGBOOST Future Predictions')
plt.show()

reg.save_model('311XGB_model.json')

reg_new = xgb.XGBRegressor()
reg_new.load_model('311XGB_model.json')

#run saved model under reg_new and plot results
future_with_features['pred'] = reg_new.predict(future_with_features[FEATURES4])

future_with_features['pred'].plot(figsize = (10,5),
                                  #color= 'Blue',
                                  color= color_pal[0],
                                  ms=1, 
                                  lw=1,
                                  title='XGBOOST Future Predictions')
plt.show()

#zoomed in plot of forecast for upcoming week; first week of may
import datetime
import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(10,5))
ax.plot(future_with_features.index, future_with_features['pred'])
#fig = future_with_features['pred'].plot(future_with_features['pred'], ax= ax)
ax.set_xlim([datetime.date(2023, 7, 1), datetime.date(2023, 7, 7)])
ax.set_ylim(0,80)
ax.set_title('XGBOOST Forecast First Week Of July')
plt.show()

#zoomed in plot of forecast for upcoming week; first week of may
import datetime
import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(10,5))
ax.plot(future_with_features.index, future_with_features['pred'])
#fig = future_with_features['pred'].plot(future_with_features['pred'], ax= ax)
ax.set_xlim([datetime.date(2023, 11, 1), datetime.date(2023, 11, 2)])
ax.set_ylim(0,80)
ax.set_title('XGBOOST Forecast November 1st By Hour')
plt.show()

from xgboost import plot_importance
plot_importance(reg_new)
plt.show()

!pip3 install pmdarima

df_day['count'].plot(figsize=(12,5))

!pip install statsmodels

#import statsmodels.api as sm
import statsmodels.api
import statsmodels.formula.api as smf

from statsmodels.tsa.stattools import adfuller

def ad_test(dataset):
    dftest = adfuller(dataset, autolag = 'AIC')
    print('1. ADF : ', dftest[0])
    print('2. P-Value : ', dftest[1])
    print('3. Num Of Lags : ', dftest[2])
    print('4. Num Of Observations Used For ADF Regression and Critical Values Calculation : ', dftest[3])
    print('5. Critical Values : ')
    for key, val in dftest[4].items():
        print('\t', key, ':',val)

# P-Value is probability and should be as low as possible
#the smaller the better P-Value
#small value means teh data set is stationary
ad_test(df_day['count'])

from pmdarima import auto_arima
# Ignore  warnings
import warnings
warnings.filterwarnings('ignore')

df_day.head()

split_date = '1-jan-2023'
df_train_day = df_day.loc[df_day.index <= split_date].copy()
df_test_day = df_day.loc[df_day.index > split_date].copy()

fig, ax = plt.subplots(figsize=(15,5))
df_train_day.plot(ax=ax, label='Training Data', title='Train/Test Data Split BY Day')
df_test_day.plot(ax=ax, label='Testing Data')
ax.axvline(split_date, color='black', ls='--')
ax.legend(['Training Data', 'Test Data'])
plt.show()

import statsmodels.api as sm
from statsmodels.tsa.arima_model import ARIMA

import statsmodels.api
import statsmodels.formula.api as smf

stepwise_fit = auto_arima(df_day['count'],
                          trace=True,
                          supress_warnings=True)

stepwise_fit.summary()  

#use order from best order in the step above

model_day=sm.tsa.ARIMA(df_train_day['count'], order=(5,1,5))
model_day=model_day.fit()
model_day.summary()

start=len(df_train_day)
end=len(df_train_day)+len(df_test_day)-1
pred_day=model_day.predict(start=start, end=end, typ='levels')
#print(pred) # dates did not appear do use next line of code
pred_day.index=df_day.index[start:end+1]
print(pred_day)

pred_day.plot(legend=True)
df_test_day['count'].plot(legend=True)

df_day['count'].mean()

from sklearn.metrics import mean_squared_error
from math import sqrt
rmse_day=sqrt(mean_squared_error(pred_day, df_test_day['count']))
rmse_day

#Fit model on all df_day
model_day = sm.tsa.ARIMA(df_day['count'], order=(5,0,5))
model_day=model_day.fit()
df_day.tail() # predict into the future after training

index_future_dates=pd.date_range(start='2023-05-09', end='2023-12-31', freq='D')
pred_day=model_day.predict(start=len(df_day), end=len(df_day)+(236), type='levels').rename('ARIMA Predictions BY Day')
pred_day.index=index_future_dates
print(pred_day)

pred_day.info()

pred_day.plot(figsize=(12,5),legend=True)

df5 = df2.copy()
#df5 = df5.drop(['is_future', 'lag1', 'lag2', 'lag3'], axis=1)
df5.info()

df5['count'].plot(figsize=(12,5))

# P-Value is probability and should be as low as possible
#the smaller the better P-Value
#small value means teh data set is stationary
ad_test(df5['count'])

stepwise_fit5 = auto_arima(df5['count'],
                          trace=True,
                          supress_warnings=True)

stepwise_fit5.summary()    

split_date = '1-jan-2023'
df_train_day = df_day.loc[df_day.index <= split_date].copy()
df_test_day = df_day.loc[df_day.index > split_date].copy()

fig, ax = plt.subplots(figsize=(15,5))
df_train_day.plot(ax=ax, label='Training Data', title='Train/Test Data Split BY Day')
df_test_day.plot(ax=ax, label='Testing Data')
ax.axvline(split_date, color='black', ls='--')
ax.legend(['Training Data', 'Test Data'])
plt.show()

#from statsmodels.tsa.arima_model import ARIMA
model_day=sm.tsa.ARIMA(df_train_day['count'], order=(5,0,5))
model_day=model_day.fit()
model_day.summary()

start=len(df_train_day)
end=len(df_train_day)+len(df_test_day)-1
pred_day=model_day.predict(start=start, end=end, typ='levels')
#print(pred) # dates did not appear do use next line of code
pred_day.index=df_day.index[start:end+1]
print(pred_day)

pred_day.plot(legend=True)
df_test_day['count'].plot(legend=True)

df_day['count'].mean()

from sklearn.metrics import mean_squared_error
from math import sqrt
rmse_day=sqrt(mean_squared_error(pred_day, df_test_day['count']))
rmse_day

#Fit model on all df5
model_day = sm.tsa.ARIMA(df_day['count'], order=(5,0,5))
model_day=model_day.fit()
df_day.tail() # predict into the future after training

index_future_dates=pd.date_range(start='2023-05-09', end='2023-12-31', freq='D')
pred_day=model_day.predict(start=len(df_day), end=len(df_day)+(236), type='levels').rename('ARIMA Predictions BY Day')
pred_day.index=index_future_dates
print(pred_day)