NGC_HT_Day2Day_ByCP_ERR_WAVG_2021.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "import plotly.express as px\n",
    "import plotly.graph_objects as go\n",
    "from plotly.subplots import make_subplots\n",
    "from plotly.offline import init_notebook_mode\n",
    "init_notebook_mode(connected=False)\n",
    "pd.set_option('display.max_colwidth', 100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "from sklearn.model_selection import KFold, train_test_split, StratifiedKFold\n",
    "from sklearn.ensemble import RandomForestRegressor, AdaBoostRegressor, ExtraTreesRegressor\n",
    "from sklearn.metrics import accuracy_score, mean_absolute_error, mean_squared_error\n",
    "from sklearn.linear_model import Ridge, Lasso, ElasticNet\n",
    "from sklearn.preprocessing import StandardScaler, MinMaxScaler, RobustScaler\n",
    "from evaluation_metrics import compute_metrics, compute_metrics_csv, mean_absolute_percentage_error, symetric_mean_absolute_percentage_error"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "from statsmodels.tsa.stattools import acf, pacf, ccf, ccovf\n",
    "from statsmodels.tsa.seasonal import STL\n",
    "from statsmodels.graphics.tsaplots import plot_acf, plot_pacf \n",
    "from statsmodels.tsa.stattools import adfuller, kpss\n",
    "from scipy.stats import boxcox, yeojohnson\n",
    "from scipy.spatial.distance import cosine, euclidean\n",
    "import rstl"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "from scipy.stats import boxcox, yeojohnson"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "from tqdm.notebook import trange, tqdm\n",
    "from datetime import datetime\n",
    "import itertools\n",
    "from numpy.polynomial import Polynomial as P"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# from lightgbm import LGBMRegressor"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "import pickle"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "import ruptures as rpt\n",
    "# from dtaidistance import dtw"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "from river import tree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df = pd.read_csv('./Data/ppnet_metar_v8_2021.csv', sep=';', index_col=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df = df[df.Year < 2021].copy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_temp_forecast = pd.read_csv('./Data/yrno_temperature.csv', sep=';', index_col=0, decimal=',')\n",
    "df_temp_forecast.index = pd.DatetimeIndex(pd.Series(df_temp_forecast.index).apply(lambda x: datetime.strptime(x, '%d.%m.%Y %H:%M')))\n",
    "df_temp_forecast = df_temp_forecast.iloc[:, 1].groupby(df_temp_forecast.index).mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_el = pd.read_csv('./Data/ppnet_metar_v8_electricity.csv', sep=';', index_col=0).Load_ffill\n",
    "df_el.index = pd.DatetimeIndex(df_el.index)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df.loc[:, 'TestSet'] = 0\n",
    "df.loc[df.Year == 2021, 'TestSet'] = 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "px.line(y=df['Consumption'], x=df.index, color=df.Year)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_gas = pd.read_csv('./Data/ppnet_metar_v8_2021.csv', sep=';', index_col=0)\n",
    "df_gas.index = pd.DatetimeIndex(df_gas.index)\n",
    "df_gas['Month'] = df_gas.index.month_name()\n",
    "df_gas['Day'] = df_gas.index.day_name()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_el = pd.read_csv('./Data/raw_electricity_2021.csv', sep=',', index_col=0)\n",
    "df_el.index = pd.DatetimeIndex(df_el.index)\n",
    "df_el['Month'] = df_el.index.month_name()\n",
    "df_el['Day'] = df_el.index.day_name()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_el"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df['Temperature_YRNO'] = pd.Series(df_temp_forecast)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "temp_err = (df.Temperature - df.Temperature_YRNO).dropna()\n",
    "temp_err.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df.isna().sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df.describe(exclude='O').columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# df['Consumption'] = np.cbrt(df['Consumption'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "FORECAST_HORIZON = 24"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df.loc[df.index <= '2014-12-31 23:00:00', 'Temperature_YRNO'] = df.loc[df.index <= '2014-12-31 23:00:00', 'Temperature']\n",
    "df.loc[:, 'Temperature_YRNO'] = df.loc[:, 'Temperature_YRNO'].interpolate()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X, y = df.loc[:, df.columns != 'Consumption'].copy(), df.Consumption.copy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X.columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X.isna().sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "y.isna().sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X = X.loc[:, ['Year', 'Month', 'Day', 'Hour', 'Day_of_week', 'Before_holiday', 'Holiday', 'Temperature', 'Cena', 'Temperature_YRNO', 'Wind_speed', 'Humidity']]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X['Cena'] = X['Cena'].bfill()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X = X.rename({'Cena' : 'Price'}, axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "def create_lagged_dataset(X, y, past_lags=24, future_temp_forecast_lags=24):\n",
    "    X_in, y_out = X.copy().drop(['Year'], axis=1), pd.DataFrame(columns=[f'H{x}' for x in range(1, 25)], index=y.index)\n",
    "    \n",
    "    for i, x in enumerate(y_out.columns):\n",
    "        y_out[f'H{i+1}'] = y.shift(-i-1)\n",
    "    \n",
    "    lagged_columns = ['Temperature', 'Price', 'Temperature_YRNO', 'Wind_speed', 'Humidity']\n",
    "    for c in lagged_columns:\n",
    "        for i in range(1, past_lags+1):\n",
    "            X_in[f'LAG_{c}_{i}'] = X[c].shift(i)\n",
    "    \n",
    "    for i in range(1, future_temp_forecast_lags + 1):\n",
    "        c = 'Temperature_YRNO'\n",
    "        X_in[f'LAG_{c}_{i}'] = X[c].shift(-i)\n",
    "        \n",
    "    for i in range(1, past_lags + 1):\n",
    "        c = 'Consumption'\n",
    "        X_in[f'LAG_{c}_{i}'] = y.shift(i)\n",
    "    \n",
    "    return (X_in, y_out)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "past_lags, future_temp_forecast_lags = 72, 24\n",
    "X_in, y_out = create_lagged_dataset(X, y, past_lags, future_temp_forecast_lags)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X_in.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "y_out.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Set the predefined setting for predefined number of change points"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "LOW, MED, HIGH = 3, 1, 0.5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "cp_params = {'l1': {111: 722870.3350949581, 62: 2132352.798169759, 46: 4274373.643209655},\n",
    " 'l2': {144: 14563484775.012444,\n",
    "  95: 59636233165.946365,\n",
    "  49: 223876329036.64392,\n",
    "  47: 244205309454.86548},\n",
    " 'rbf': {88: 15.142857142857142, 66: 25.24489795918367, 47: 45.44897959183673}}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "min_size = 24 * 7\n",
    "jump = 24\n",
    "model = 'l2'\n",
    "pen = cp_params[model][47]*HIGH\n",
    "algo = rpt.Pelt(model=model, min_size=min_size, jump=jump)\n",
    "algo.fit(y.values)\n",
    "my_bkps = algo.predict(pen=pen)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "def plot_bkps(data, my_bkps, title=None, ret=False):\n",
    "    fig = go.Figure()\n",
    "    fig.add_trace(go.Scatter(x=data.index, y=data.values,\n",
    "                            mode='lines',\n",
    "                            line_color='#333C83'))\n",
    "    # fig = px.line(x=data.index, y=data.values, title=title)\n",
    "    clr_selection = 'green'\n",
    "    color_switch = lambda x: 'blue' if x != 'blue' else 'green'\n",
    "\n",
    "    for idx, cp in enumerate(my_bkps):\n",
    "        if cp >= len(data):\n",
    "            break\n",
    "\n",
    "        fig.add_vline(x=data.index[cp], line_width=3, line_dash=\"dash\", line_color=\"red\")\n",
    "\n",
    "        if idx < len(my_bkps) - 2:\n",
    "            clr_selection = color_switch(clr_selection)\n",
    "            fig.add_vrect(x0=data.index[cp], x1=data.index[my_bkps[idx+1]], line_width=0, fillcolor=clr_selection, opacity=0.1)\n",
    "    \n",
    "    # fig.data = fig.data[::-1]\n",
    "    fig.update_layout(title=title)\n",
    "    if ret:\n",
    "        return fig\n",
    "    fig.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "plot_bkps(y, my_bkps)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X_in_nona, y_out_nona = X_in.iloc[past_lags:-future_temp_forecast_lags, :], y_out.iloc[past_lags:-future_temp_forecast_lags, :]\n",
    "X_in_nona.index = pd.DatetimeIndex(X_in_nona.index)\n",
    "y_out_nona.index = pd.DatetimeIndex(y_out_nona.index)\n",
    "X_in_nona_midnight, y_out_nona_midnight = X_in_nona[X_in_nona.index.hour == 0], y_out_nona[y_out_nona.index.hour == 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X_in_nona_midnight.index.quarter"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "len(my_bkps)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "X.index = pd.DatetimeIndex(X.index)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# bkps_index = X.index[my_bkps[:-1]]\n",
    "# bkps_index = list(bkps_index[bkps_index.year == 2013]) + list(bkps_index[bkps_index.year == 2013] + pd.offsets.DateOffset(years=1))\n",
    "# bkps_index.insert(0, bkps_index[-1] + pd.offsets.DateOffset(years=-2))\n",
    "# bkps_index = pd.DatetimeIndex(bkps_index)\n",
    "# bkps_index"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "bkps_index = X.index[my_bkps[:-1]]\n",
    "bkps_index = list(bkps_index[bkps_index.year == 2013])\n",
    "bkps_index = pd.DatetimeIndex(bkps_index)\n",
    "bkps_base = bkps_index[bkps_index.year == 2013]\n",
    "\n",
    "bkps_index = list(bkps_index)\n",
    "year_offset = 8\n",
    "\n",
    "for i in range(1, year_offset):\n",
    "    bkps_index += list(bkps_base + pd.offsets.DateOffset(years=i))\n",
    "    \n",
    "bkps_index.insert(0, bkps_index[-1] + pd.offsets.DateOffset(years=-year_offset))\n",
    "bkps_index = pd.DatetimeIndex(bkps_index)\n",
    "bkps_index"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "models_dict = {}\n",
    "\n",
    "for i in bkps_index:\n",
    "    \n",
    "    if i not in models_dict:\n",
    "        prev = i+ pd.offsets.DateOffset(years=-1)\n",
    "        if prev in models_dict:\n",
    "            models_dict[i] = models_dict[prev]\n",
    "        else:\n",
    "            models_dict[i] = {c: tree.HoeffdingTreeRegressor(grace_period=7, leaf_prediction='adaptive', model_selector_decay=0.2, tau=0.5) for c in y_out.columns}\n",
    "\n",
    "# for i in bkps_index[1:8]:\n",
    "#     print(i)\n",
    "#     models_dict[i] = {c: tree.HoeffdingTreeRegressor(grace_period=7, leaf_prediction='adaptive', model_selector_decay=0.2, tie_threshold=0.5) for c in y_out.columns}\n",
    "    \n",
    "# for i in bkps_index[8:]:\n",
    "#     print(i+ pd.offsets.DateOffset(years=-1))\n",
    "#     models_dict[i+ pd.offsets.DateOffset(years=-1)] = {c: tree.HoeffdingTreeRegressor(grace_period=7, leaf_prediction='adaptive', model_selector_decay=0.2, tie_threshold=0.5) for c in y_out.columns}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "models_dict.keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "win_width = 7\n",
    "cp_ranges = []\n",
    "\n",
    "for idx, i in enumerate(bkps_index[8:]):\n",
    "    cp_ranges.append((i - pd.offsets.DateOffset(days=win_width), i + pd.offsets.DateOffset(days=win_width), bkps_index[idx+8-1], i))\n",
    "\n",
    "cp_ranges"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "y_pred_arr1 = np.zeros(24)\n",
    "y_pred_arr2 = np.ones(24)\n",
    "np.average((y_pred_arr1, y_pred_arr2), axis=0, weights=[0.1, 0.5]).tolist()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_forecast = pd.DataFrame(columns=y_out_nona_midnight.columns, index=y_out_nona_midnight.index).fillna(0)\n",
    "limit = len(y_out_nona_midnight)\n",
    "skip = 0\n",
    "idx = 0\n",
    "last_model = None\n",
    "model_order = 0\n",
    "win_learn = True\n",
    "prev_err1, prev_err2 = 0, 1\n",
    "\n",
    "for inp, out in tqdm(zip(X_in_nona_midnight[skip:limit].iterrows(), y_out_nona_midnight[skip:limit].iterrows()), total=limit):\n",
    "    y_pred_arr1 = np.zeros(len(y_out_nona.columns))\n",
    "    y_pred_arr2 = np.zeros(len(y_out_nona.columns))\n",
    "    y_pred_arr = np.zeros(len(y_out_nona.columns))\n",
    "    q = inp[0].quarter\n",
    "    year = inp[0].year\n",
    "    model_selection = list(models_dict.values())[0]\n",
    "    \n",
    "    curr_date = inp[0]\n",
    "    for k in reversed(list(models_dict.keys())):\n",
    "        if curr_date >= k:\n",
    "            model_selection = models_dict[k]\n",
    "            model_selection_date = k\n",
    "            break\n",
    "    \n",
    "    if last_model != model_selection_date:\n",
    "        print(model_selection_date)\n",
    "        model_order += 1\n",
    "    \n",
    "    win_mode = None\n",
    "    if year > 2013:\n",
    "        for x in cp_ranges:\n",
    "            l, r = x[0], x[1]\n",
    "            if inp[0] >= l and inp[0] <= r:\n",
    "                win_mode = x\n",
    "                break\n",
    "    \n",
    "    if year == 2013 or win_mode is None:\n",
    "        for k, c in enumerate(y_out_nona.columns):\n",
    "            m = model_selection[c]\n",
    "            y_pred_arr[k] =  m.predict_one(inp[1].to_dict())\n",
    "            m.learn_one(inp[1].to_dict(), out[1][k])\n",
    "        df_forecast.iloc[idx, :] = y_pred_arr\n",
    "    else:\n",
    "        print('WMode ',win_mode[1])\n",
    "        _, rb, m1i, m2i = win_mode\n",
    "        m1s, m2s = models_dict[m1i], models_dict[m2i]\n",
    "        for k, c in enumerate(y_out_nona.columns):\n",
    "            m1, m2 = m1s[c], m2s[c]\n",
    "            # y_pred_arr1[k] = np.average((m1.predict_one(inp[1].to_dict()), m2.predict_one(inp[1].to_dict())), axis=0).tolist()\n",
    "            y_pred_arr1[k] = m1.predict_one(inp[1].to_dict())\n",
    "            y_pred_arr2[k] = m2.predict_one(inp[1].to_dict())\n",
    "            \n",
    "        if win_learn:\n",
    "            for k, c in enumerate(y_out_nona.columns):\n",
    "                m1.learn_one(inp[1].to_dict(), out[1][k])\n",
    "                m2.learn_one(inp[1].to_dict(), out[1][k])\n",
    "                # if prev_err1 < prev_err2:\n",
    "                #     m1.learn_one(inp[1].to_dict(), out[1][k])\n",
    "                # else:\n",
    "                #     m2.learn_one(inp[1].to_dict(), out[1][k])\n",
    "                    \n",
    "        # print(f'{prev_err1} :  {prev_err2} -> using model {1 if prev_err1 < prev_err2 else 2}')\n",
    "        sum_err = prev_err1 + prev_err2\n",
    "        df_forecast.iloc[idx, :] = np.average((y_pred_arr1, y_pred_arr2), axis=0, weights=[1-prev_err1/sum_err, 1-prev_err2/sum_err]).tolist()\n",
    "        prev_err1 = np.sum(np.abs(out[1] - y_pred_arr1))\n",
    "        prev_err2 = np.sum(np.abs(out[1] - y_pred_arr2))\n",
    "        if inp[0] == rb:\n",
    "            prev_err1 = 0\n",
    "            prev_err2 = 1\n",
    "            \n",
    "        \n",
    "    # print(df_forecast.iloc[idx, :].sum())\n",
    "    idx += 1\n",
    "    last_model = model_selection_date"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_forecast.dropna()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "y_out_selection = y[(y.index >= '2013-01-04 00:00:00') & (y.index <= '2020-12-30 23:00:00')]\n",
    "y_out_selection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_forecast_res = df_forecast[df_forecast.index <= '2020-12-30']\n",
    "df_forecast_res = df_forecast_res.stack().reset_index(name='y_pred')\n",
    "df_forecast_res"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "y_out_selection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_forecast_res = df_forecast[df_forecast.index <= '2020-12-30']\n",
    "df_forecast_res = df_forecast_res.stack().reset_index(name='y_pred')\n",
    "df_forecast_res.index = y_out_selection.index[:]\n",
    "\n",
    "df_forecast_res = df_forecast_res.shift()\n",
    "df_forecast_res = df_forecast_res.dropna()\n",
    "df_forecast_res = df_forecast_res.drop(['level_0', 'level_1'], axis=1)\n",
    "df_forecast_res['y_true'] = y_out_selection.iloc[1:]\n",
    "df_forecast_res"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_forecast_res.to_csv('Results/NGC_HT_Day2Day_ByCP_HIGH_ERR_WAVG_2021.csv', sep=',')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## IQR correction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# q3, q1 = np.percentile(df_forecast_res.y_pred, [75 ,25])\n",
    "# iqr = q3 - q1\n",
    "# qt, qb = q3 + 1.5*iqr, q1 - 1.5*iqr\n",
    "\n",
    "\n",
    "# df_forecast_res['y_pred'] = df_forecast_res['y_pred'].apply(lambda x: np.nan if x > qt or x < qb else x)\n",
    "# df_forecast_res['y_pred'] = df_forecast_res['y_pred'].fillna(method='ffill')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "compute_metrics(df_forecast_res[df_forecast_res.index >= '2013-01-05 00:00:00'].dropna())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "compute_metrics(df_forecast_res[df_forecast_res.index >= '2014-01-05 00:00:00'].dropna())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "compute_metrics(df_forecast_res[df_forecast_res.index >= '2014-01-05 00:00:00'].dropna()).columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_pred_all = pd.DataFrame(columns=['MAE', 'MSE', 'MAPE', 'SMAPE', 'R2', 'WAPE'], index=[x for x in range(2014, 2021)])\n",
    "\n",
    "for year in range(2014, 2021):\n",
    "    tmp = compute_metrics(df_forecast_res[(df_forecast_res.index >= f'{year}-01-01 00:00:00') & (df_forecast_res.index <= f'{year}-12-31 00:00:00')].dropna()).values\n",
    "    df_pred_all.loc[year, :] = tmp\n",
    "df_pred_all"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_forecast_res_s = df_forecast_res[df_forecast_res.index >= '2013-01-05 00:00:00'].stack().reset_index().rename({'level_1': 'Type', 0: 'Value'}, axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_forecast_res_s"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "df_err = df_forecast_res[df_forecast_res.index >= '2013-01-05 00:00:00'].copy()\n",
    "df_err.loc[:, 'AE'] = (df_err.y_true - df_err.y_pred).abs()\n",
    "df_err.loc[:, 'SE'] = (df_err.y_true - df_err.y_pred) ** 2\n",
    "df_err.loc[:, 'SAPE'] = ((df_err.y_pred - df_err.y_true).abs() / ((df_err.y_true.abs() + df_err.y_pred.abs())/2)) * 100\n",
    "df_err"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "def plot_pred_bkps(data, my_bkps, from_dt=None, title=None):\n",
    "    if from_dt:\n",
    "        data = data[data.index >= from_dt]\n",
    "    my_bkps = [x for x in my_bkps if x < len(data)]\n",
    "    \n",
    "    df_err = data.copy()\n",
    "    df_err.loc[:, 'AE'] = (df_err.y_true - df_err.y_pred).abs()\n",
    "    df_err.loc[:, 'SE'] = (df_err.y_true - df_err.y_pred) ** 2\n",
    "    df_err.loc[:, 'SAPE'] = ((df_err.y_pred - df_err.y_true).abs() / ((df_err.y_true.abs() + df_err.y_pred.abs())/2)) * 100\n",
    "    \n",
    "    fig = go.Figure()\n",
    "    fig.add_trace(go.Scatter(x=data.index, y=data.y_true,\n",
    "                            mode='lines',\n",
    "                            line_color='#333C83', name='y_true'))\n",
    "    \n",
    "    fig.add_trace(go.Scatter(x=data.index, y=data.y_pred,\n",
    "                            mode='lines',\n",
    "                            line_color='red', name='y_pred'))\n",
    "    clr_selection = 'green'\n",
    "    color_switch = lambda x: 'blue' if x != 'blue' else 'green'\n",
    "    \n",
    "    for idx, cp in enumerate(my_bkps):\n",
    "        if cp >= len(data):\n",
    "            break\n",
    "\n",
    "        fig.add_vline(x=data.index[cp], line_width=3, line_dash=\"dash\", line_color=\"red\")\n",
    "\n",
    "        if idx < len(my_bkps) - 2:\n",
    "            clr_selection = color_switch(clr_selection)\n",
    "            fig.add_vrect(x0=data.index[cp], x1=data.index[my_bkps[idx+1]], line_width=0, fillcolor=clr_selection, opacity=0.1)\n",
    "    \n",
    "    fig.update_layout(title=title)\n",
    "    fig.show()\n",
    "    \n",
    "    for roll in [24, 24*7]:\n",
    "        fig = go.Figure()\n",
    "        fig.add_trace(go.Scatter(x=data.index, y=df_err['SAPE'].rolling(roll).mean(),\n",
    "                                mode='lines',\n",
    "                                line_color='#333C83', name=f'SAPE rolling({roll})'))\n",
    "\n",
    "        for idx, cp in enumerate(my_bkps):\n",
    "                if cp >= len(data):\n",
    "                    break\n",
    "\n",
    "                fig.add_vline(x=data.index[cp], line_width=3, line_dash=\"dash\", line_color=\"red\")\n",
    "\n",
    "        fig.update_layout(title=f'SAPE rolling({roll})')\n",
    "        fig.show()\n",
    "    \n",
    "    for x in reversed(['AE', 'SE', 'SAPE']):\n",
    "        fig = go.Figure()\n",
    "        fig.add_trace(go.Scatter(x=data.index, y=df_err[x],\n",
    "                                mode='lines',\n",
    "                                line_color='#333C83', name=x))\n",
    "\n",
    "        for idx, cp in enumerate(my_bkps):\n",
    "            if cp >= len(data):\n",
    "                break\n",
    "\n",
    "            fig.add_vline(x=data.index[cp], line_width=3, line_dash=\"dash\", line_color=\"red\")\n",
    "        fig.update_layout(title=x)\n",
    "        fig.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "plot_pred_bkps(df_forecast_res, my_bkps, '2013-01-05 00:00:00')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "plot_pred_bkps(df_forecast_res, my_bkps, '2014-01-01 00:00:00')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.0rc1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}