# SKTIME
# Period: quarterly, Metrics: Population
# datasets: https://www.sktime.org/en/stable/api_reference/datasets.html?highlight=datasets#

pip install sktime

pip install pmdarima

import matplotlib
import matplotlib.pyplot as plt
from matplotlib.pyplot import title
import pandas as pd
import pmdarima
import sktime
from sktime.datasets import load_airline
from sktime.utils.plotting import plot_series
import numpy as np
from sktime.utils import plotting
from statsmodels.tsa.api import *
from sktime.registry import all_estimators
from sktime.forecasting.arima import ARIMA
from sktime.forecasting.var import VAR

from sktime.forecasting.model_selection import temporal_train_test_split
from sktime.forecasting.base import ForecastingHorizon
from sktime.forecasting.theta import ThetaForecaster
from sktime.forecasting.naive import NaiveForecaster
from sktime.performance_metrics.forecasting import mean_absolute_percentage_error
from sktime.forecasting.arima import AutoARIMA
from sktime.utils.plotting import plot_series
import seaborn as sns

# 1. plot the time series, use the dataset provided by sktime library with 1 Date variable and 1 Float variable

# define the time series
y = load_airline()

# plotting the graph of the basic series
fig, ax = plotting.plot_series(y, labels=["y"],markers=[''])
ax.legend(loc='upper left', frameon=False)

# define colors for line and legend
ax.get_lines()[0].set_c("silver")
leg = ax.get_legend()
leg.legendHandles[0].set_color('silver')

y

# revise the type of series
y.dtypes

# the primary index of time series is Period
y.index

# 2. plot the forecasted y
# forecasting horizon defines how many steps/periods ahead are forecasted
# there are many ways how to define forecasting horizon

fh = np.arange(1, 37) # the first 37 steps ahead
fh

fh = np.array([2, 5])  # 2nd and 5th step ahead

fh = ForecastingHorizon(
    pd.PeriodIndex(pd.date_range("1961-01", periods=36, freq="M")), is_relative=False)
fh

cutoff = pd.Period("1960-12", freq="M") # define the beginning

fh.to_relative(cutoff)

fh.to_absolute(cutoff)

# data specification
y = load_airline()


# specify forecasting horizon
fh = np.arange(1, 37)

# specify the forecasting algorithm
forecaster = NaiveForecaster(strategy="last", sp=12)

# fit the forecaster
forecaster.fit(y)

# create predictions
y_pred = forecaster.predict(fh)

# plot predictions and past data
fig, ax = plotting.plot_series(y, y_pred, labels=["y", "y_pred"],markers=['',''])
ax.legend(loc='upper left', frameon=False)

ax.get_lines()[0].set_c("silver")
ax.get_lines()[-1].set_c("black")

leg = ax.get_legend()
leg.legendHandles[0].set_color('silver')
leg.legendHandles[-1].set_color('black')

# There are many sktime forecasters: NaiveForecaster, Thetaforecaster, 
# NaiveVariance, Arima, Exponential smoothing, Multiplex forecaster,
# TransformendTargetForecaster


# 3. Univariate Theta Forecaster

# until fit, define how many periods we want to forecast ahead, e.g. 19

y = load_airline()

fh = np.arange(1, 20)

forecaster = ThetaForecaster(sp=19)
forecaster.fit(y, fh=fh)

# coverage has values between 0 and 1, it defines the width of filling
# coverage is a special parameter of predict_interval
coverage = 0.9
y_pred_ints = forecaster.predict_interval(coverage=coverage)
y_pred_ints

from sktime.utils import plotting

# create predictions
y_pred = forecaster.predict()

fig, ax = plotting.plot_series(y, y_pred, labels=["y", "y_pred"])
ax.get_lines()[0].set_c("silver") # set color for y line
ax.get_lines()[-1].set_c("black") # set color for y_pred line
ax.legend(loc='upper left', frameon=False)

ax.fill_between(
    ax.get_lines()[-1].get_xdata(),
    y_pred_ints["Coverage"][coverage]["lower"], # lower limit for filling
    y_pred_ints["Coverage"][coverage]["upper"], # upper limit for filling
    alpha=0.1,# transparency of filling
    color='k',
    label=f"{coverage}% prediction intervals",

)

ax.get_lines()[-1].get_c()

# 4. Univariate ARIMA with deseasonaliser and detrender
# TransformedTargetForecasted has hardcoded statistical methodology to consider the season and trend as well

from sktime.datasets import load_macroeconomic
from sktime.forecasting.arima import ARIMA
from sktime.forecasting.compose import TransformedTargetForecaster
from sktime.transformations.series.detrend import Deseasonalizer, Detrender

# use the available training time series dataset and only 1 Date variable and 1 Float variable
data = load_macroeconomic()
y = pd.Series(data["unemp"], name='Unemployment')

y

# forecast 3 periods ahead
forecaster = TransformedTargetForecaster(
    [
        ("deseasonalize", Deseasonalizer(sp=12)), # get the trend for 12 periods, 12 quarters, 3 years
        ("detrend", Detrender()),
        ("forecast", ARIMA()),
    ]
)

forecaster.fit(y, fh=[1, 2, 3, 4, 5, 6, 7, 8, 9 ]) # forecast 9 periods, 9 quarters ahead
forecaster.predict_interval()

forecaster.predict_quantiles()

forecaster = Deseasonalizer(sp=12) * Detrender() * ARIMA()

forecaster.fit(y, fh=[1, 2, 3, 4, 5, 6, 7, 8, 9])
forecaster.predict_interval()

coverage = 0.9
y_pred_ints = forecaster.predict_interval(coverage=coverage)
y_pred_ints

# create unemployment rate forecasts for given periods with lower and upper limits
# see how the each further period has less clear forecasted values / variance / higher error
y_pred = forecaster.predict()

fig, ax = plotting.plot_series(y, y_pred, labels=["y", "y_pred"],markers=['',''])
ax.get_lines()[0].set_c("black") # set color for y line
ax.get_lines()[-1].set_c("silver") # set color for y_pred line
ax.legend(loc='upper left', frameon=False)

leg = ax.get_legend()
leg.legendHandles[0].set_color('black')
leg.legendHandles[-1].set_color('silver')

ax.fill_between(
    ax.get_lines()[-1].get_xdata(),
    y_pred_ints["Coverage"][coverage]["lower"], # lower limit for filling
    y_pred_ints["Coverage"][coverage]["upper"], # upper limit for filling
    alpha=0.1,# transparency of filling
    color='k',
    label=f"{coverage}% prediction intervals",

)