2_pricing_model.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "6d432863-8993-41c8-b879-68b5061f90c2",
   "metadata": {},
   "outputs": [],
   "source": [
    "from datetime import datetime\n",
    "import sys\n",
    "\n",
    "import geohashlite\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib.dates as mdates\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import plotly.express as px\n",
    "import seaborn as sns\n",
    "from sklearn.model_selection import TimeSeriesSplit\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from xgboost import XGBRegressor\n",
    "\n",
    "pd.set_option(\"display.max_columns\", None)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "1f2048f3-47b6-4147-bc96-d4b310b49dd1",
   "metadata": {},
   "outputs": [],
   "source": [
    "trip = pd.read_csv(\"data/Taxi_Trips_-_2021.csv\")\n",
    "\n",
    "# parse data\n",
    "trip['Trip Start Timestamp'] = trip['Trip Start Timestamp'].map(lambda x: datetime.strptime(x, \"%m/%d/%Y %I:%M:%S %p\") if pd.notnull(x) else None)\n",
    "trip['Trip End Timestamp'] = trip['Trip End Timestamp'].map(lambda x: datetime.strptime(x, \"%m/%d/%Y %I:%M:%S %p\") if pd.notnull(x) else None)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "086b5fed-9961-4ddf-9b13-c1a778a23293",
   "metadata": {},
   "source": [
    "# Data Cleaning "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "b758429c-9b15-47d8-8c0e-6a99785f6fb8",
   "metadata": {},
   "outputs": [],
   "source": [
    "def filter_outlier(df, colname, percentile=0.02, clip=(-sys.maxsize, sys.maxsize)):\n",
    "    \"\"\"\n",
    "    clip: (lower bound, upper bound) to clip the values\n",
    "    \n",
    "    Filter out outliers\n",
    "    Lower outer fence = Q1 - k * IQR\n",
    "    Upper outer fence = Q3 + k * IQR\n",
    "    \"\"\"\n",
    "    k = 6\n",
    "    value = df[colname]\n",
    "    q1, q3 = np.nanpercentile(value, [25, 75])\n",
    "    iqr = q3-q1\n",
    "    lower_outer_fence = max(clip[0], q1 - k * iqr)\n",
    "    upper_outer_fence = min(clip[1], q3 + k * iqr)\n",
    "    \n",
    "    print(f\"Trimming extreme values in column `{colname}`: below {round(lower_outer_fence, 2)}, or larger than {round(upper_outer_fence, 2)}\")\n",
    "    out_df = df.loc[(df[colname] < upper_outer_fence) & \n",
    "                    (df[colname] > lower_outer_fence)].copy()\n",
    "    return out_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "1b7555dc-62ae-44e6-b097-c8890adb4258",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Trimming extreme values in column `Fare`: below 0, or larger than 170.0\n",
      "Trimming extreme values in column `Trip Miles`: below 0, or larger than 65.24\n",
      "Trimming extreme values in column `Trip Seconds`: below 0, or larger than 8608.0\n"
     ]
    }
   ],
   "source": [
    "# remove rows w/o fare labels \n",
    "trip_cleaned = trip.loc[~trip.Fare.isnull()].copy() \n",
    "\n",
    "# filter extreme outlier (high tendency of data entry error)\n",
    "trip_cleaned = filter_outlier(trip_cleaned, \"Fare\", clip=(0, sys.maxsize))\n",
    "trip_cleaned = filter_outlier(trip_cleaned, \"Trip Miles\", clip=(0, sys.maxsize)) # clip left\n",
    "trip_cleaned = filter_outlier(trip_cleaned, \"Trip Seconds\", clip=(0, sys.maxsize)) # clip left"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "65fc9086-7777-4412-8c56-f4caea42137b",
   "metadata": {},
   "source": [
    "***"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3ea688f4-4812-41ee-afb8-56007851228c",
   "metadata": {},
   "source": [
    "# Feature Engineering"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "b0bf788c-310f-4ccd-be1f-8a07087ce8e9",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_geohash(lat, lon, precision=7):\n",
    "    \"\"\"\n",
    "    Convert latitude, longitude to geohash\n",
    "    \"\"\"\n",
    "    if pd.notnull(lat) and pd.notnull(lat):\n",
    "        return geohashlite.encode(lat,lon ,precision)\n",
    "    else:\n",
    "        return None"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "86157cbc-9e7d-476c-be75-377574b2ef91",
   "metadata": {},
   "source": [
    "## Feature: Demand score\n",
    "Reminder: Ensure no data leak (aggregated data score should not appear in test period)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "13ae3745-42fe-402c-a0fe-c1087ee11206",
   "metadata": {},
   "outputs": [],
   "source": [
    "def gen_demand_score(df):\n",
    "    \"\"\"\n",
    "    Generate the table for demand_score(day_of_week, hour, geohash)\n",
    "        > higher score means higher demand\n",
    "    \"\"\"\n",
    "    print(\"Generating demand score ...\")\n",
    "    \n",
    "    REQUIRED_COLS = [\"Trip ID\", \"Pickup Centroid Latitude\", \"Pickup Centroid Longitude\", \"Trip Start Timestamp\"]\n",
    "    assert all(c in df.columns for c in REQUIRED_COLS), f\"Did not have all the required columns: {REQUIRED_COLS}\"\n",
    "    \n",
    "    # transform table to get the necessary metadata\n",
    "    df = df.loc[:, REQUIRED_COLS].dropna().copy()\n",
    "    df['day_of_week'] = df['Trip Start Timestamp'].map(lambda x: int(datetime.strftime(x,\"%w\")))\n",
    "    df['hour'] = df['Trip Start Timestamp'].map(lambda x: int(datetime.strftime(x,\"%-H\")))\n",
    "    df['geohash'] = df.apply(lambda x: get_geohash(x['Pickup Centroid Latitude'], x['Pickup Centroid Longitude']),axis=1)\n",
    "    \n",
    "    # aggregate to get scores \n",
    "    pickup_demand_score = df.groupby([\"day_of_week\", \"hour\", \"geohash\"])['Trip ID'].count().reset_index(name=\"num_trips\") # \n",
    "    pickup_demand_score['demand_score'] = (pickup_demand_score.num_trips - pickup_demand_score.num_trips.mean()) / pickup_demand_score.num_trips.std() # normalise\n",
    "    \n",
    "    return pickup_demand_score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "a99892ae-843c-4cb7-b6af-62ee2c5453fb",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Generating demand score ...\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>day_of_week</th>\n",
       "      <th>hour</th>\n",
       "      <th>geohash</th>\n",
       "      <th>num_trips</th>\n",
       "      <th>demand_score</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1304</th>\n",
       "      <td>0</td>\n",
       "      <td>19</td>\n",
       "      <td>dp3qz6r</td>\n",
       "      <td>276</td>\n",
       "      <td>10.487691</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1379</th>\n",
       "      <td>0</td>\n",
       "      <td>20</td>\n",
       "      <td>dp3qz6r</td>\n",
       "      <td>271</td>\n",
       "      <td>10.288883</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1222</th>\n",
       "      <td>0</td>\n",
       "      <td>18</td>\n",
       "      <td>dp3qz6r</td>\n",
       "      <td>261</td>\n",
       "      <td>9.891267</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10079</th>\n",
       "      <td>5</td>\n",
       "      <td>16</td>\n",
       "      <td>dp3wmge</td>\n",
       "      <td>258</td>\n",
       "      <td>9.771982</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10176</th>\n",
       "      <td>5</td>\n",
       "      <td>17</td>\n",
       "      <td>dp3wmge</td>\n",
       "      <td>254</td>\n",
       "      <td>9.612936</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "       day_of_week  hour  geohash  num_trips  demand_score\n",
       "1304             0    19  dp3qz6r        276     10.487691\n",
       "1379             0    20  dp3qz6r        271     10.288883\n",
       "1222             0    18  dp3qz6r        261      9.891267\n",
       "10079            5    16  dp3wmge        258      9.771982\n",
       "10176            5    17  dp3wmge        254      9.612936"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# sanity check\n",
    "gen_demand_score(trip_cleaned.sample(frac=0.05)).sort_values(\"demand_score\", ascending=False).head()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "675a3429-a943-4469-84a1-778112c7ca9f",
   "metadata": {},
   "source": [
    "## Feature: Supply score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "d670882c-5faf-47fb-97c4-53ffa79b4baa",
   "metadata": {},
   "outputs": [],
   "source": [
    "def gen_supply_score(df):\n",
    "    \"\"\"\n",
    "    Generate the table for supply_score(day_of_week, hour, geohash)\n",
    "        > higher means supply is higher\n",
    "    \"\"\"\n",
    "    print(\"Generating supply score ...\")\n",
    "    \n",
    "    REQUIRED_COLS = [\"Taxi ID\", \"Dropoff Centroid Latitude\", \"Dropoff Centroid Longitude\", \"Trip Start Timestamp\"]\n",
    "    assert all(c in df.columns for c in REQUIRED_COLS), f\"Did not have all the required columns: {REQUIRED_COLS}\"\n",
    "    \n",
    "    # transform table to get the necessary metadata\n",
    "    df = df.loc[:, [\"Taxi ID\", \"Dropoff Centroid Latitude\", \"Dropoff Centroid Longitude\", \"Trip Start Timestamp\"]].dropna().copy()\n",
    "    df['day_of_week'] = df['Trip Start Timestamp'].map(lambda x: int(datetime.strftime(x,\"%w\")))\n",
    "    df['hour'] = df['Trip Start Timestamp'].map(lambda x: int(datetime.strftime(x,\"%-H\")))\n",
    "    df['geohash'] = df.apply(lambda x: get_geohash(x['Dropoff Centroid Latitude'], x['Dropoff Centroid Longitude']),axis=1)\n",
    "    \n",
    "    # aggregate to get scores \n",
    "    pickup_supply_score = df.groupby([\"day_of_week\", \"hour\", \"geohash\"])['Taxi ID'].agg(n_drivers=lambda x: len(set(x))).reset_index()\n",
    "    pickup_supply_score['supply_score'] = (pickup_supply_score.n_drivers - pickup_supply_score.n_drivers.mean()) / pickup_supply_score.n_drivers.std() # normalise\n",
    "    \n",
    "    return pickup_supply_score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "b0ff2fc2-eba3-416f-b0a3-f2f28eb92191",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Generating supply score ...\n",
      "   hour  day_of_week  geohash  n_drivers  supply_score\n",
      "0     0            6  dp3wvp8         94      5.026376\n",
      "1     1            6  dp3wvp8         78      4.071988\n",
      "2     2            6  dp3wvp8         62      3.117600\n",
      "3     3            6  dp3wvp8         32      1.328123\n",
      "4     4            6  dp3wvp8         27      1.029877\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "Text(0.5, 1.0, 'Supply Score')"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbIAAADvCAYAAACXFFGBAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy89olMNAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA/SUlEQVR4nO3dd3gU1frA8e+mF9JJJw0IhNBN6KggoSpNkSIIKMUrQYn8FMQGXAsKXrqicpUi1YbIVToYWiB0KaEEAgRIo6T33fn9gazGUBLYZHaT9/M888jMnJl5Zxz2Zc6cOUejKIqCEEIIYaLM1A5ACCGEeBiSyIQQQpg0SWRCCCFMmiQyIYQQJk0SmRBCCJMmiUwIIYRJk0QmhBDCpEkiE0IIYdIkkQkhhDBpksiEqESLFy9Go9Fw4cIFtUMRosqQRCaqjGPHjtGvXz8CAgKwsbHB19eXzp07M2/ePLVDq1Dr1q3j8ccfx8PDAzs7O2rXrk3//v3ZsGGD2qEJUSkkkYkqYc+ePYSHh3P06FFGjRrF/PnzGTlyJGZmZsyZM0ft8CrMp59+Sq9evdBoNEyaNIlZs2bxzDPPcPbsWVatWqV2eEJUCgu1AxDCED788EOcnJzYv38/zs7OJdalpqaqE1QFKy4u5v3336dz585s2rSp1PrKPG+dTkdhYSE2NjaVdkwhbpMnMlElnDt3joYNG5ZKYgAeHh76P1+4cAGNRsPixYtLldNoNEyZMkU/P2XKFDQaDadOnaJ///44Ojri5ubGuHHjyM/PL7Xt2LFjWb58OfXr18fGxoawsDB27Nhxz7iHDRtGzZo1KSoqKrWuS5cu1K9f/67bXrt2jczMTNq1a3fH9X8/b4D8/HymTJlCvXr1sLGxwdvbm6effppz587py+Tk5PB///d/+Pn5YW1tTf369fn000/55yAZfz/fhg0bYm1tra/KvHLlCi+++CKenp5YW1vTsGFDvvnmm3teByEehiQyUSUEBARw8OBBjh8/bvB99+/fn/z8fKZNm0aPHj2YO3cuo0ePLlUuOjqaqKgohgwZwr///W+uX79Ot27d7hnT888/z/Xr19m4cWOJ5cnJyWzbto0hQ4bcdVsPDw9sbW1Zt24dN27cuOc5aLVannrqKaZOnUpYWBj/+c9/GDduHBkZGfr4FEWhV69ezJo1i27dujFz5kzq16/PG2+8wfjx40vtc9u2bbz22msMGDCAOXPmEBgYSEpKCq1bt2bLli2MHTuWOXPmULduXUaMGMHs2bPvGaMQD0wRogrYtGmTYm5urpibmytt2rRRJkyYoGzcuFEpLCwsUS4hIUEBlEWLFpXaB6BMnjxZPz958mQFUHr16lWi3JgxYxRAOXr0aIltAeXAgQP6ZRcvXlRsbGyUvn376pctWrRIAZSEhARFURRFq9UqtWrVUgYMGFDiGDNnzlQ0Go1y/vz5e573e++9pwCKvb290r17d+XDDz9UDh48WKrcN998owDKzJkzS63T6XSKoijKzz//rADKBx98UGJ9v379FI1Go8THx5c4XzMzM+XEiRMlyo4YMULx9vZWrl27VmL5wIEDFScnJyU3N/ee5yPEg5BEJqqM2NhYpW/fvoqdnZ0+sbi7uytr167Vl3mQRLZx48YS5eLi4hRAmTZtWolt27RpU2qfAwYMUOzs7JTi4mJFUUonMkVRlIkTJyq2trZKZmamfllYWJjSrl27Mp33ihUrlPbt2ytmZmb6827evLly8uRJfZknn3xSqVmzplJUVHTX/YwePVoxNzcvEYeiKEpMTIwCKPPmzStxvh07dixRTqfTKc7Ozsro0aOVtLS0EtPt8961a1eZzkmI8pCqRVFltGjRgp9++ombN28SGxvLpEmTyMrKol+/fpw8efKB9xscHFxivk6dOpiZmZX6Fuyf5QDq1atHbm4uaWlpd93/0KFDycvLY82aNQCcPn2agwcP8vzzz5cpvkGDBrFz505u3rzJpk2beO655zh8+DA9e/bUv8s7d+4c9evXx8Li7u27Ll68iI+PDw4ODiWWN2jQQL/+74KCgkrMp6WlkZ6ezldffYW7u3uJ6YUXXgCqbsMboS5ptSiqHCsrK1q0aEGLFi2oV68eL7zwAt9//z2TJ09Go9HccRutVlvm/d9tHw8qNDSUsLAwli1bxtChQ1m2bBlWVlb079+/XPtxdHSkc+fOdO7cGUtLS5YsWcK+fft4/PHHDRrvbba2tiXmdTodAEOGDGHYsGF33KZJkyYVEouo3iSRiSotPDwcgKSkJABcXFwASE9PL1Hun08bf3f27NkSTx/x8fHodDoCAwNLlfunM2fOYGdnh7u7+z3jHDp0KOPHjycpKYkVK1bw5JNP6mN9EOHh4SxZskR/3nXq1GHfvn0UFRVhaWl5x20CAgLYsmULWVlZJZ7KTp06pV9/L+7u7jg4OKDVaomIiHjg2IUoL6laFFXC9u3bSzURB/jtt98A9M3YHR0dqVmzZqlm8Z9//vld9/3ZZ5+VmL/dU0j37t1LLI+JieHQoUP6+cTERNauXUuXLl0wNze/Z/yDBg1Co9Ewbtw4zp8/f8/Wirfl5uYSExNzx3Xr168H/jrvZ555hmvXrjF//vxSZW9ftx49eqDVakuVmTVrFhqNptT5/pO5uTnPPPMMP/744x1bat6relWIhyFPZKJKeOWVV8jNzaVv376EhIRQWFjInj17WL16NYGBgfp3NAAjR47k448/ZuTIkYSHh7Njxw7OnDlz130nJCTQq1cvunXrRkxMDMuWLeO5556jadOmJco1atSIrl278uqrr2Jtba1PjlOnTr1v/O7u7nTr1o3vv/8eZ2dnnnzyyftuk5ubS9u2bWndujXdunXDz8+P9PR0fv75Z3bu3EmfPn1o3rw5cOuJb+nSpYwfP57Y2FgeffRRcnJy2LJlC2PGjKF379707NmTjh078vbbb3PhwgWaNm3Kpk2bWLt2LVFRUdSpU+e+MX388cds376dVq1aMWrUKEJDQ7lx4waHDh1iy5Yt9/1MQIgHonZrEyEMYf369cqLL76ohISEKDVq1FCsrKyUunXrKq+88oqSkpJSomxubq4yYsQIxcnJSXFwcFD69++vpKam3rXV4smTJ5V+/fopDg4OiouLizJ27FglLy+vxD4BJTIyUlm2bJkSHBysWFtbK82bN1e2b99eotydWi3e9t133ymAMnr06DKdc1FRkbJw4UKlT58+SkBAgGJtba3Y2dkpzZs3V2bMmKEUFBSUOu+3335bCQoKUiwtLRUvLy+lX79+yrlz5/RlsrKylNdee03x8fFRLC0tleDgYGXGjBn6Jvr/PN87SUlJUSIjIxU/Pz/9cTp16qR89dVXZTovIcpLoyh3qI8RQjBlyhSmTp1KWloaNWvWvGdZjUZDZGTkHavuymrt2rX06dOHHTt28Oijjz7wfoSobuQdmRBGYuHChdSuXZv27durHYoQJkXekQmhslWrVvHHH3/w66+/MmfOHIM37xeiqpNEJoTKBg0aRI0aNRgxYgRjxoxROxwhTI68IxNCCGHS5B2ZEEIIkyaJTAghhEkzundkOp2Oq1ev4uDgIC+9hRCiGlMUhaysLHx8fDAzu/tzl9ElsqtXr+Ln56d2GEIIIYxEYmIitWrVuut6o0tktzsrTUxMxNHRUeVohBBCqCUzMxM/P79SQwv9U7kT2Y4dO5gxYwYHDx4kKSmJNWvW0KdPH/16RVGYPHkyCxcuJD09nXbt2rFgwYI7jtV0J7erEx0dHSWRCSGEuO9rpnI39sjJyaFp06alegS/bfr06cydO5cvvviCffv2YW9vT9euXfUD/AkhhBCGVO4nsu7du991OAdFUZg9ezbvvPMOvXv3BmDp0qV4enry888/M3DgwIeLVgghhPgHgza/T0hIIDk5ucSgek5OTrRq1equ4yYVFBSQmZlZYhJCGKdr2QXEnLuOVif9KAjjYdDGHsnJyQB4enqWWO7p6alf90/Tpk0r03hNf6coCsXFxeUanl5UP5aWlvcd0FLcX36Rli1xKfx06ArRZ9LQ6hQeq+fO3IHNcLazUjs8IdRvtThp0iTGjx+vn7/dSuVuCgsLSUpKIjc3tzLCEyZMo9FQq1YtatSooXYoJkenUzhw8SY/HbrMr8eSyMov1q8zN9Ow40waPefv4ssh4YT6SKMsoS6DJjIvLy8AUlJS8Pb21i9PSUmhWbNmd9zG2toaa2vrMu1fp9ORkJCAubk5Pj4+WFlZyUfT4o4URSEtLY3Lly8THBwsT2ZllHAthzWHLrPmyBUSb+Tpl/s42dCnuS9PP+JLYbHCS8sOkHgjj6cX7OaTZ5rQu5mvilGL6s6giSwoKAgvLy+2bt2qT1yZmZns27ePl19++aH3X1hYiE6nw8/PDzs7u4fen6ja3N3duXDhAkVFRZLI7iE9t5B1fyTx06HLHL6Url9ub2VOj8be9H3El9ZBbpiZ/fWPxnVj2/PqqiPsOJPGuFVHOJqYwaQeIViaS693ovKVO5FlZ2cTHx+vn09ISODIkSO4urri7+9PVFQUH3zwAcHBwQQFBfHuu+/i4+NT4luzh3WvrkqEuE2e1u+usFjH9tOp/HToMttOpVKkvdV4w0wDjwa78/QjvnQJ9cLW6s7/AHC2s2LR8BbM2nyG+dvj+WZ3AieuZjD/uUdwdyhbDYsQhlLuRHbgwAE6duyon7/9fmvYsGEsXryYCRMmkJOTw+jRo0lPT6d9+/Zs2LABGxsbw0UthHggRVodc7ee5du9F0nPLdIvD/V25OlHfOnV1AcPx7L9XTU30/B61/o08nXi9e+Psi/hBj3n7eKL58No5udcQWcgRGlGNx5ZZmYmTk5OZGRklOrZIz8/n4SEBIKCgiQxivuS+6WktKwCxq44xL6EGwB4OFjTp7kvfZv70sD74RpsxKdmM/rbA5xPy8HK3Ix/927IwJb+hghbVGP3ygd/J3V01cSFCxfQaDQcOXJE7VCECo4kptNz3i72JdyghrUFcwY2I2ZSJ97q0eChkxhAXY8arI1sR5dQTwq1Ot786RiTfjpGQbF8IiMqniQyIaq4VbGX6P9FDMmZ+dR2t+fnyHb0buaLuZlh3yE62FjyxZAw3uhaH40GVsZeYuBXe0nOkO7pRMWSRCaMXmFhodohmKSCYi2TfjrGmz8do1Cro0uoJ2sj21HXo+K+qzMz0xDZsS6LhrfAydaSw5fSeWreLmL/rM4UoiKYfCJTFIXcwmJVpvK8Xvzhhx9o3Lgxtra2uLm5ERERQU5ODh06dCAqKqpE2T59+jB8+HD9fGBgIO+//z6DBg3C3t4eX1/fUp02azQaFixYQPfu3bG1taV27dr88MMPd71mdevW5dNPPy2x/MiRI2g0mhKtUu+2/ZQpU/D398fa2hofHx9effVV/fqCggImTpyIn58f1tbW1K1bl6+//lq/Pjo6mpYtW2JtbY23tzdvvvkmxcV/fXDboUMHxo4dS1RUFDVr1qRr164AHD9+nO7du1OjRg08PT15/vnnuXbt2j1jra6SM/IZ+NVeVsZeQqOBN7rW54shYTjYWFbK8TvU92Dd2PaEeDlwLbuA5xbuZdHuhHL9nRGirFTv2eNh5RVpCX1voyrHPvnvrthZ3f8SJiUlMWjQIKZPn07fvn3Jyspi586d5fpLPWPGDN566y2mTp3Kxo0bGTduHPXq1aNz5876Mu+++y4ff/wxc+bM4dtvv2XgwIEcO3aMBg0alNiXRqPhxRdfZNGiRbz++uv65YsWLeKxxx6jbt2694zlxx9/ZNasWaxatYqGDRuSnJzM0aNH9euHDh1KTEwMc+fOpWnTpiQkJOgTzpUrV+jRowfDhw9n6dKlnDp1ilGjRmFjY8OUKVP0+1iyZAkvv/wyu3fvBiA9PZ0nnniCkSNHMmvWLPLy8pg4cSL9+/dn27ZtZb6O1cG+89eJXHGIa9mFONlaMmdgMzrU96j0OPzd7PhpTFve/PEYvxy9ytR1Jzl2OYMP+za+a7N+IR6EyScyU5CUlERxcTFPP/00AQEBADRu3Lhc+2jXrh1vvvkmAPXq1WP37t3MmjWrRCJ79tlnGTlyJADvv/8+mzdvZt68eXz++eel9jd8+HDee+89YmNjadmyJUVFRaxYsaLUU9qdXLp0CS8vLyIiIrC0tMTf35+WLVsCcObMGb777js2b96s7zy6du3a+m0///xz/Pz8mD9/PhqNhpCQEK5evcrEiRN577339N8IBgcHM336dP12H3zwAc2bN+ejjz7SL/vmm2/w8/PjzJkz1KtXr8zXsqpSFIXFey7w4a9xFOsUQrwc+Or5cPzd1Os8wM7qVsOSpn7OfPRbHD8dvtVf43Ot/BncKgAvJ2lNKh6eyScyW0tzTv67q2rHLoumTZvSqVMnGjduTNeuXenSpQv9+vXDxcWlzMdq06ZNqfnZs2fft8zdWin6+Pjw5JNP8s0339CyZUvWrVtHQUEBzz777H1jefbZZ5k9eza1a9emW7du9OjRg549e2JhYcGRI0cwNzfn8ccfv+O2cXFxtGnTpsTHyu3atSM7O5vLly/j73+ryXZYWFiJ7Y4ePcr27dvv2G/iuXPnqn0iyyvU8taaY6w5fAWAXk19+PiZxmWqMahoGo2GEe2DCPV25PXvj3IlPY952+JZ8Ps5ujbyYnjbQMIDXOQDdvHA1L/LH5JGozGKv6z3Ym5uzubNm9mzZw+bNm1i3rx5vP322+zbtw8zM7NSVYxFRUV32ZNhjRw5kueff55Zs2axaNEiBgwYUKauv/z8/Dh9+jRbtmxh8+bNjBkzhhkzZhAdHY2tra1BYrO3ty8xn52dTc+ePfnkk09Klf17v57GLq9Qy5X0XBJv5pGRW0RgTXvqezo8VFVb4o1cXvr2ICeTMjE30/BWjwa82C7Q6BJDmzpuRL/RgU0nU1i85wKxCTf49Y8kfv0jiVBvR4a3DaRXMx9syvgPRCFuM+4MUIVoNBratWtHu3bteO+99wgICGDNmjW4u7uTlJSkL6fVajl+/HiJ3lMA9u7dW2r+n+++9u7dy9ChQ0vMN2/e/K4x9ejRA3t7exYsWMCGDRvYsWNHmc/H1taWnj170rNnTyIjIwkJCeHYsWM0btwYnU5HdHR0iXHpbmvQoAE//vgjiqLof2h3796Ng4MDtWrVuuvxHnnkEX788UcCAwOxsDDe2zYjr4jLN3O5cjOPK+l5XLmZx+Xbf07P40ZO6RaYZhoIrGlPA29HGng53PqvtyPeTjb3TUY7zqTx6qrDpOcW4WZvxfznHqFNHbeKOr2HZmFuRo/G3vRo7M3Jq5ksjbnAmsNXOJmUyYQf/2Da+jgGtvRnSOsAfJ0N848iUfUZ7y9CFbJv3z62bt1Kly5d8PDwYN++faSlpdGgQQPs7e0ZP348v/76K3Xq1GHmzJmkp6eX2sfu3buZPn06ffr0YfPmzXz//ff8+uuvJcp8//33hIeH0759e5YvX05sbGyJ1oL/ZG5uzvDhw5k0aRLBwcGlqibvZvHixWi1Wlq1aoWdnR3Lli3D1taWgIAA3NzcGDZsGC+++KK+scfFixdJTU2lf//+jBkzhtmzZ/PKK68wduxYTp8+zeTJkxk/fvw9+9CMjIxk4cKFDBo0iAkTJuDq6kp8fDyrVq3iv//9ryqdAu89f53fjiWVSFpZBcX33c7B2gJfF1ucbC05l5bDtewCzqflcD4th1//+OsfNU62ljTwdiDEy5HQP5NbsGcNbCzNURSFBdHn+HTjaXQKNK3lxIIhYfiY0I9/qI8jHz/ThIndQlh9IJFvYy5yJT2PBb+f48voc3QJ9WJY20Ba13Y1uqdLYVwkkVUCR0dHduzYwezZs8nMzCQgIID//Oc/dO/enaKiIo4ePcrQoUOxsLDgtddeK/U0BvB///d/HDhwgKlTp+Lo6MjMmTP1zdJvmzp1KqtWrWLMmDF4e3uzcuVKQkND7xnbiBEj+Oijj3jhhRfKfD7Ozs58/PHHjB8/Hq1WS+PGjVm3bh1ubreeBBYsWMBbb73FmDFjuH79Ov7+/rz11lsA+Pr68ttvv/HGG2/QtGlTXF1dGTFiBO+88849j+nj48Pu3buZOHEiXbp0oaCggICAALp166ZKJ9LrjyUxduXhO46U7GZvha+LLb7Ot6ZaLrb4utjdmv8zgf1dWlYBcUmZ+ulUchbxqdlk5BWx9/wN9p7/6xssczMNQTXtqWFtwZHEdAAGhPsxtXdDk62Sc7G34l+P12HUo7XZEpfCkj0X2HPuOhtOJLPhRDIhXg4MbRNI3+a+0tpR3JH0tWgCAgMDiYqKKvW92d9pNBrWrFlT7lEGdu7cSadOnUhMTCw1srepq6j7ZWtcCi99e5BinUJEA086hrjrE5aPs61B3tkWFGuJT80mLimrRJK7+beOfi3NNUzp1ZDnWvpXuSeW08lZLIm5wJpDV8grutXNlZOtJQNa+PF86wD8XGUYp+qgrH0tyhNZNVVQUEBaWhpTpkzh2WefrXJJrKLsPJvGy8sOUaxT6NnUh9kDmhm8qycAawtzGvo40dDHSb9MURRSswo4mZTJhWs5tK7tZpB+Eo1RfS8HPurbmIldQ/j+YCJLYy5y6UYuX+04z8Kd54lo4MnwtoG0reNW5ZK4KD9JZNXUypUrGTFiBM2aNWPp0qUl1i1fvpyXXnrpjtsFBARw4sSJygjR6Ow9f51RSw9QqNXRtaEnM/s3rZAkdjcajQZPRxs8HW2gfqUdVlVOdpaMfLQ2L7QLYvupVJbEXGDn2WtsPpnC5pMpBHvUYGjbQJ5u7ou9tfycVVdStShKycrKIiUl5Y7rLC0t9R91GztD3i8HL95k6Nf7yCnU0rG+O188H4a1hbyvUUN8ajZLYy7ww8HL5BbeqnZ0sLGgf7gfQ9sEEOBmf589CFNR1qpFSWSiyjLU/XL8SgaDFu4lK7+YdnXd+HpYC5NtWFGVZOYX8cOByyyNucCF67kAaDTQsb4Hw9sG0r5uTcwq8YlZGF6VfkdmZLlXGClD3CenkjMZ8vU+svKLaRnoysKh4ZLEjISjjSUvtg9ieNtAos+msWTPBX4/nca2U6lsO5VKbXd7hrUJ5JmwWtSQascqzaT+71pa3mq2nJuba7AeJETVdXv4lwf9xiw+NZsh/91Hem4Rzfyc+Xp4uNH3IlMdmZlp6Fjfg471PTifls3SmIv8cPAy59NymPzLCWZsPE2/sFoMbOlHLRc77K3MpYFIFWNSVYtwqwPe9PR0PDw8sLOzkxtS3JFOp+Pq1av6To3Le59cvJ5D/y9jSMksINTbkZWjWuNkVzlDoIiHl11QzE+HLrNkzwXOpeWUWGemuTUIqIONBY42ljjaWuBgY4nj7WW2ljj+Y52rvRUhXg7ye1PJquQ7MrhVXZScnHzH3i+E+DszMzOCgoKwsrIq13ZX0vPo/0UMV9LzqOdZg1Wj2+BqX759COOgKAq74q+xePcFdsZfo7BY98D76tHYi/mDHpH3bpWoyiay27RabaV1ritMk5WVVbl7/UjJzKf/lzFcvJ5L7Zr2rHqpNR4O0rCoKlAUhYJiHZl5RWTmF5OZX0RWfvGf83/9OevPdX//8/m0HIp1Cm90rU9kx3uP1ycMR7XGHlqtlilTprBs2TKSk5Px8fFh+PDhvPPOOwZ9LDc3N1elfz1Rdd0eyfji9Vz8XG1ZPqqVJLEqRKPRYGNpjo2lOR7l/I589f5LTPzxGJ9uOk1jXyceq+deMUGKB2LwTuo++eQTFixYwPz584mLi+OTTz5h+vTpzJs3z9CHEsJgbuYUMuS/+ziXloOPkw0rRrbG20kaFIlbBrTwZ2ALPxQFxq06zOWbuWqHJP7G4Ilsz5499O7dmyeffJLAwED69etHly5diI2NNfShhDCIjLwihn4Ty6nkLDwcrFk+qrX05SdKmdKrIY19nbiZW8SY5YfI/7MPSKE+gyeytm3bsnXrVs6cOQPcGtl3165ddO/e/Y7lCwoKyMzMLDEJUVmyC4p5YVEsx65k4GpvxfKRrQiqKT1DiNJsLM1ZMOQRnO0s+eNyBlPXVc+u2oyRwRPZm2++ycCBAwkJCcHS0pLmzZsTFRXF4MGD71h+2rRpODk56Sc/Pz9DhyREKXmFWjafTGHYN7EcupSOk60ly0a0ItjTQe3QhBGr5WLH3IHN0WhgZWwiq/dfUjskQQW0Wly1ahVvvPEGM2bMoGHDhhw5coSoqChmzpzJsGHDSpUvKCigoKBAP5+ZmYmfn999W6kIUV6pWflsi0tlS1wKO89eo+DPptgO1hYsG9mKpn7O6gYoTMb8bWf5dNMZrCzM+PFfbWlcy+n+G4lyU635vZ+fH2+++SaRkZH6ZR988AHLli3j1KlT992+rIELcT+KonAmJZstcbd6Sr89EOVtvs62dA71ZHArf3kSE+Wi0ymM/vYAW+JS8XW25X+vtMdFvjU0ONWa3+fm5pb6dsfc3Byd7sE/RBSirIq0OvYn3GBzXApb4lJIvJFXYn3TWk5ENPAkItRTemoQD8zMTMN/+jej1/xdXLyey7jVR1g0vEWlDusj/mLwRNazZ08+/PBD/P39adiwIYcPH2bmzJm8+OKLhj6UEMCtVofRZ9LYcjKF7adTycov1q+zsjCjfd2aRDTwpFMDj1tjeQlhAE62lnwxJIy+n+9mx5k05mw5w/gu1WSgOCNj8KrFrKws3n33XdasWUNqaio+Pj4MGjSI9957r0xdBUnVoiiPDceTeGXlYYq0f93GbvZWPBHiQUSoJ48G15SOfkWF+vnwFaJWHwHg62HhdGpg2NHWFUWptjUHVb6LKiEA+n8ZQ2zCDQLd7OjWyJvOoR4083ORKh5Rqab8coLFey7gYGPBurHtCXzITzgURWFLXCrzt8dzPi2br4e1oGWQq4GiNR1VejwyIeBWK8T9F24AsHxUa3ydpScOoY63ejTg2JUMDl68yb+WHWTNmHbYWpW/Cz2dTmH98WTmbTvLqeQs/fIxyw/x66vtpWr8Lgz+HZkQlWXjiRQUBZr6OUsSE6qysjDj88GPULOGNaeSs3hrzbFyDexarNXx8+ErdJm9g8gVhziVnIW9lTn/erwOIV4OXMsuYMzyQw/Ve39VJolMmKz1x5IA6NHIS+VIhABPRxvmP9ccczMNaw5f4du9F++7TZFWx3cHEomYGU3U6iPEp2bjYGPBq52C2f3mE7zZPYQvhoThYGPBwYs3+ei3uEo4E9MjiUyYpOvZBew9fx2A7o28VY5GiFta13ZjUvcQAP697iQHL964Y7mCYi3L9l6kw4zfmfDDH1y4nouLnSVvdK3P7jefYHznejjb3WocF1jTnln9mwGweM8F1h65UinnYkrkHZkwSZtOpqBToJGvI/5u0sGvMB4j2gdx+FI6vx5LYszyQ/zvlUdxd7AGIL9Iy8rYS3wZfZ7kzHwAatawZvRjQQxuFYC99Z1/kiNCPXnlibrM2xbPmz8eo76XAyFe0hjuNklkwiT99me1ojyNCWOj0Wj4pF8TTqdkEZ+azdgVh/hqaDirYi+xcGcC17Jvdcnn5WjDS4/XZlBLf2ws798wJCqiHkcS09l59hr/+vYgv7zSHkcby4o+HZMgze+FybmZU0j4h1vQ6hS2v95BeqsXRik+NZve83eRU6jFytyMQu2thhq+zraM6ViHfmG1sLYoX8vGmzmFPDVvF1fS8+gc6smXQ8Iwq8KfmpQ1H8g7MmFyNseloNUpNPB2lCQmjFZdjxp8+mxTAAq1OgLd7Jjerwm/v9GBwa0Cyp3EAFzsrVgw5BGsLMzYfDKFBdHnDB22SZKqRWFypLWiMBXdG3vzzfBwCot1RDTwxML84Z8dmtRy5v3eDZn44zH+s+k0TWo58WiwuwGiNV3yRCZMSkZeEbvirwG3fiSEMHZPhHjSrZG3QZLYbQNa+DOwhR86BV5deZgr6Xn336gKk0QmTMrWuBSKtAr1PGtQ16OG2uEIoZopvRrS2NeJm7lFvLzsIPlFWrVDUo0kMmFSfjuWDEhrRSFsLM1ZMOQRnO0s+eNyBlPXnVA7JNVIIhMmIyu/iB1n0wDoIdWKQlDLxY65A5uj0cDK2ERW77+kdkiqkEQmTMa2U6kUFuuo7W5PPU+pVhQC4LF67vxf53oAvLv2BMcuZ6gcUeWTRCZMxvo/qxV7NPKutuMzCXEnYzrUJaKBB4XFOv617CA3cwofep+pmfkm895Nmt8Lk5BTUMz206kAdG8sze6F+DszMw3/6d+MXvN3cfF6LuNWH2HR8BZlHpdPURQuXM8lNuE6+xJusO/8Da6k5+Fka8mEbvUZ1MLfqD+8lkQmTMLvp9MoKNYR4GZHqLf0+CLEPznZWvLFkDD6fr6bHWfSmLPlDOO71L9jWUVROJua/WfSuk5swg1SswpKlcvIK+LtNcf57sBlPuzTiEa+ThV9Gg9EEpkwCb8d/6tvRalWFOLOGng7Mu3pxry2+ihzt8XT1M+ZTg080eoU4pIy2Zdwg9iE6+y/cJMb/6h+tDI3o5mfMy2DXGlV25Wmfs78dPAy/9l0hqOJ6fSav4vnWwcwvkt9nGyNq49H6WtRGL28Qi1hH2wmt1DLL2Pb0aSWs9ohCWHUJq89zpKYizjYWBAe4MKBCzfJKiguUcbG0oywABdaBrrRqrYrzfyc79h5cWpmPh/+FsfaI1eBW731v/NkA3o386nwf1SWNR9IIhNGb8PxZP617CC+zrbsmthRnsiEuI/CYh0Dv4rh0KV0/bIa1haEB7rQKsiNlkGuNPZ1wsqi7O399sRf4521xzmflgNA69qufNCnEXU9HAwdvl5Z80GFVC1euXKFiRMnsn79enJzc6lbty6LFi0iPDy8Ig4nqrj1f1Yr9mjsJUlMiDKwsjDjy+fD+fz3eHydbWkV5Eaoj2OZG3/cSdu6NVk/7lH+uzOBedvOsvf8DbrN3smox2rzyhN1sbNS702VwY988+ZN2rVrR8eOHVm/fj3u7u6cPXsWFxcXQx9KVAP5RVq2xt1urSgfQQtRVu4O1kzu2dCg+7S2MCeyY116NfVh6roTbIlLZcHv5/jlyFUm9wylc6inKv/YNHgi++STT/Dz82PRokX6ZUFBQYY+jKgmdp29RnZBMd5ONjSTd2NCGAU/Vzv+O6wFm0+mMOWXE1xJz2P0twfpFOLBlF4N8XOt3FHbDf5B9C+//EJ4eDjPPvssHh4eNG/enIULF961fEFBAZmZmSUmIW673VqxWyMvo/6ORYjqqHOoJ5vHP8aYDnWwNNew9VQqETOjmb/tLAXFlfcxtcET2fnz51mwYAHBwcFs3LiRl19+mVdffZUlS5bcsfy0adNwcnLST35+foYOSZiowmIdm0+mANK3ohDGys7KggndQlg/7lHa1HajoFjHp5vO0H32TvZfuFEpMRi81aKVlRXh4eHs2bNHv+zVV19l//79xMTElCpfUFBAQcFfH+JlZmbi5+cnrRYF20+n8sKi/Xg4WLN3Uid5IhPCyCmKwi9Hr/L+/+K4ll3AshGtaB9c84H3p1qrRW9vb0JDQ0ssa9CgAT/++OMdy1tbW2NtbW3oMEQVcHskaKlWFMI0aDQaejfzpWOIBxuPJz9UEisPg1cttmvXjtOnT5dYdubMGQICAgx9KFGFFWl1bPqzWlHGHhPCtDjaWPJseOW9JjJ4InvttdfYu3cvH330EfHx8axYsYKvvvqKyMhIQx9KVGF7z18nPbcIN3srWga5qh2OEMKIGTyRtWjRgjVr1rBy5UoaNWrE+++/z+zZsxk8eLChDyWqsNsjQXdt5PVQH3EKIaq+CvkU+6mnnuKpp56qiF2LaqBYq2PTib/GHhNCiHuRgTWF0Ym9cIPrOYW42FnSqrZUKwoh7k0SmTA6t0eC7hLqhaW53KJCiHuTXwlhVLQ6hQ1/VivKSNBCiLKQRCaMysGLN0nLKsDRxoK2dSrnGxQhhGmTRCaMym9/fgTdOdSrXGMlCSGqL/mlEEZDp1PYcPzP1opSrSiEKKMqm8gy84sqtfdl8fAOJ6aTnJlPDWuLSuvaRghh+qpkIssuKGbYN7GMXnqQ/CJJZqbidt+KEQ08sLYwVzkaIYSpqJKJ7GxKFqeSsog+k8aIJfvJLSxWOyRxH4qisP747daK8hG0EKLsqmQia+7vwpIXW2JvZc7u+OsMX7Sf7AJJZsbsj8sZXEnPw87KnMfruasdjhDChFTJRAbQMsiVpSNa4WBtQWzCDYZ9E0tmfpHaYYm7uD0S9BMhHthYSrWiEKLsqmwiAwgLcGHZyFY42lhw8OJNnv86loxcSWbGRlEUfW8eMhK0EKK8qnQiA2jq58yKUa1xsbPkaGI6z/13LzdzCtUOS/zNiauZXLqRi42lGR3qS7WiEKJ8qnwiA2jk68TK0a2pWcOKE1czGbRwL9eyC9QOS/xp/d+qFe2sKmRABiFEFVYtEhlAiJcjq0a3xt3BmlPJWQz6ai+pmflqh1Xt5Rdp+fnwVUBGghZCPJhqk8gA6no4sHp0a7wcbTibms3Ar/aSnCHJTE2Ldl/gSnoeno7WRDTwVDscIYQJqlaJDKC2ew2+e6kNvs62nL+WQ/8vY7h8M1ftsKqltKwCPtseD8CEriHYWklrRSFE+VW7RAbg72bH6pda4+9qx6UbuQz4ci+Xrksyq2wzN58hu6CYJrWc6NvcV+1whBAmqlomMoBaLreSWVBNe66k5zHgqxgSruWoHVa1cSo5k9X7LwHwzpOhmJlpVI5ICGGqqm0iA/B2smX16NbU9ahBUkY+A76MIT41S+2wqjxFUfjw1zh0yq1e7lsGuaodkhDChFV4Ivv444/RaDRERUVV9KEeiIejDatGtybEy4HUrAIGfrWX08mSzCrS9tOp7Dx7DStzM97s1kDtcIQQJq5CE9n+/fv58ssvadKkSUUe5qHVrGHNylGtaejjyLXsQgZ+FcOJqxlqh1UlFWl1fPBrHAAvtA/E381O5YiEEKauwhJZdnY2gwcPZuHChbi4uFTUYQzGxd6KFSNb07SWEzdzi3hu4T7+uJyudlhVzvK9FzmfloObvRWRHeuqHY4QogqosEQWGRnJk08+SURExD3LFRQUkJmZWWJSi5OdJd+ObMUj/s5k5BUxeOE+Dl68qVo8VU1GbhGzt54F4LXO9XC0sVQ5IiFEVVAhiWzVqlUcOnSIadOm3bfstGnTcHJy0k9+fn4VEVKZOdpYsnREK1oGuZJVUMzQr/cRm3BD1Ziqijlbz5KeW0Q9zxoMbKHu/2chRNVh8ESWmJjIuHHjWL58OTY2NvctP2nSJDIyMvRTYmKioUMqtxrWFix+oQVt67iRU6hl2Dex7Dl3Te2wTNr5tGyWxlwAbjW3tzCv1g1mhRAGZPBfk4MHD5KamsojjzyChYUFFhYWREdHM3fuXCwsLNBqtSXKW1tb4+joWGIyBnZWFnwzvAWP1XMnr0jLC4v2s+NMmtphmayPfjtFsU6hY313HpOBM4UQBmTwRNapUyeOHTvGkSNH9FN4eDiDBw/myJEjmJubTjdENpbmfPV8GJ1CPCgo1jFyyQG2nUpROyyTsyf+GlviUjA30/D2k9LcXghhWAZPZA4ODjRq1KjEZG9vj5ubG40aNTL04SqcjaU5C4aE0bWhJ4VaHS99e5CNJ5LVDstkaHUK7//Z3H5IK3/qejioHJEQoqqRFxVlYGVhxvznHuGpJt4UaRUilx/i1z+S1A7LJPxwMJG4pEwcbSyIiqindjhCiCqoUkYx/P333yvjMBXK0tyM2QOaYWluxprDV3hl5SGKtM3oI53d3lV2QTEzNp4B4NVOwbjYW6kckRCiKpInsnKwMDfj02eb8mxYLXQKvPbdEb4/oH4rS2O14Pd4rmUXEOhmx9A2gWqHI4SooiSRlZO5mYZPnmnCc638URR444c/WLHvktphGZ3LN3NZuDMBgEk9GmBlIbeaEKJiyK/LAzAz0/Bhn0YMbxsIwFtrjum/kRK3fLLhNIXFOlrXdqVLqIz8LISoOJLIHpBGo2Fyz1BGPRoEwHtrT/DfnedVjso4HLx4k3VHr6LRwLtPhaLRyFhjQoiKI4nsIWg0Gt7q0YDIjnUA+ODXOD7/PV7lqNSl0ym8/7+TADwbVouGPk4qRySEqOokkT0kjUbD613qExURDMD0DaeZs+WsylGpZ90fVzmSmI6dlTmvd6mvdjhCiGpAEpkBaDQaoiLq8UbXWz/cs7ac4ZtdCSpHVfnyCrV8sv4UAGM61MHD8f59bQohxMOSRGZAkR3rMqHbrWQ2c/MZrmcXqBxR5frvzvNczcjHx8mGkY/WVjscIUQ1IYnMwP71WB0a+TqSXVDMvG3V531ZamY+C6LPATCxewg2lqbTp6YQwrRJIjMwMzMNb3W/1THusr0XuXAtR+WIKsenm06TW6ilmZ8zvZr6qB2OEKIakURWAdrWrUmH+u4U6xRmbDytdjgVSlEUPv89nu8PXgakub0QovJJIqsgb3YPQaOBX48lcfjSTbXDqRDZBcW8vOwQ0zecRlFgZPsgwgJc1A5LCFHNSCKrICFejvR7pBYA0347haIoKkdkWOfSsunz2W42nEjG0lzDR30b885ToWqHJYSohiSRVaDxXephbWFG7IUbbIlLVTscg9l8MoU+83cTn5qNp6M1q19qw3Ot/NUOSwhRTUkiq0DeTraMaH+rC6uP18dRrNWpHNHD0eoUZm46zailB8gqKKZloCvrXmnPI/5SnSiEUI8ksgr2rw51cLGz5FxaDt8duKx2OA8sI7eIEUv2M/fPTwqGtw1k+ahWeDjIR89CCHVJIqtgjjaWvNrpVvdVs7acIaegWOWIyu9Ucia9PtvF76fTsLYwY2b/pkzp1RBLc7l9hBDqk1+iSjC4VQD+rnakZRXw352m1XXVL0ev0vezPVy8nkstF1t+fLktT//ZiEUIIYyBJLJKYGVhpu+66ssd50jLMv6uq4q1Oj789SSvrjxMXpGWR4Nrsm5sexr5Sm/2QgjjIomskjzZ2JumtZzILdQyZ+sZtcO5p+vZBTz/dax+hOeXO9Rh8QstcbG3UjkyIYQoTRJZJdFoNEzqcavrqpWxiZxLy1Y5ojv743I6PeftIub8deyszFkw+BEmdgvB3Ex66xBCGCeDJ7Jp06bRokULHBwc8PDwoE+fPpw+XbW7aSqr1rXdiGjggVanMH3DKbXDKeW7A4n0+yKGqxn5BNW0Z21kO7o39lY7LCGEuCeDJ7Lo6GgiIyPZu3cvmzdvpqioiC5dupCTUz06z72fid1CMNPAxhMpHLhwQ+1w9OZvO8uEH/6gsFhHRAMP1o5tR7Cng9phCSHEfWmUCu47KS0tDQ8PD6Kjo3nsscfuWz4zMxMnJycyMjJwdHSsyNBUM+mnP1gZm8gj/s78+HJb1TvZXb7vIm+vOQ7Aq52CieoUjJlUJQohVFbWfFDh78gyMjIAcHV1veP6goICMjMzS0xVXVREPWwtzTl0KZ2NJ5JVjWX9sSTe/flWEnvlibqM71xPkpgQwqRUaCLT6XRERUXRrl07GjVqdMcy06ZNw8nJST/5+flVZEhGwdPRhlGP3uq66pMNpylSqeuqPeeuMW7VEXQKDGrpz/jO9VSJQwghHkaFJrLIyEiOHz/OqlWr7lpm0qRJZGRk6KfExMSKDMlojH68Dm72ViRcy2Fl7KVKP/7xKxmMXnqQQq2Org09+aBPI9WrOIUQ4kFUWCIbO3Ys//vf/9i+fTu1at29Jwhra2scHR1LTNVBDWsLoiJudV01Z8tZsvKLKu3YF67lMHxRLNkFxbQKcmXOwObSvF4IYbIMnsgURWHs2LGsWbOGbdu2ERQUZOhDVBkDW/pTu6Y913MK+WrH+Uo5ZmpWPkO/ieVadiENvB1ZOCwcG0vzSjm2EEJUBIMnssjISJYtW8aKFStwcHAgOTmZ5ORk8vLyDH0ok2dpbsaEbiEALNx5npTM/Ao9XmZ+EcO+2c+lG7n4u9qx5MUWONpYVugxhRCiohk8kS1YsICMjAw6dOiAt7e3flq9erWhD1UldG3oSViAC/lFOmZtrriuq/KLtIxacoC4pExq1rDm2xEtZQgWIUSVUCFVi3eahg8fbuhDVQkajYa3etx6KvvuQCJnUrIMfgytTmHcqsPsS7hBDWsLFr/QggA3e4MfRwgh1CB9LRqBsABXujX0QqfAJ+sN23WVoii88/MxNp5IwcrcjIVDw6UHeyFElSKJzEhM6FYfczMNW0+lEnPuusH2+59NZ1gZm4iZBuYOakabOm4G27cQQhgDSWRGorZ7DZ5r6Q/Aq6sO88H/TrL3/HWKH+Jj6cW7E5i/PR6AD/o0plsj6QBYCFH1VHhfi+VVHfpavJu0rAJ6z9/F1Yy/Wi862VryRIgHEQ08eaxeTRzK2Mpw7ZErjFt1BID/61yPVzoFV0TIQghRYcqaDySRGZmcgmJ2nk1j88lUtp1K4WbuXx9KW5praF3bjc6hnnRq4Imvs+0d97HjTBojluynSKswrE0AU3o1lF47hBAmRxJZFaDVKRy6dJMtJ1PYfDKF89dKDoUT6u1IRKgnXUI9aejjiEaj4UhiOs8t3EtuoZanmngzd2Bz6QRYCGGSJJFVQefSstlyMoUtcSkcvHgT3d/+z3k52tAxxIMNx5O4mVvEo8E1+XpYC6ws5DWoEMI0SSKr4q5nF7D9dBpbTqaw42wauYVa/bomtZxYMao1NawtVIxQCCEeTlnzgfzSmSi3Gtb0C6tFv7Ba5BdpiTl/nS0nU8jML2ZKz1BJYkKIakN+7aoAG0tzOtb3oGN9D7VDEUKISicvUIQQQpg0SWRCCCFMmiQyIYQQJk0SmRBCCJNmdI09bn8NkJmZqXIkQggh1HQ7D9zvKzGjS2RZWbfG4/Lz81M5EiGEEMYgKysLJ6e7Dz9ldB9E63Q6rl69ioODw0P1D5iZmYmfnx+JiYnyYXUZyTUrH7le5SfXrHyq+/VSFIWsrCx8fHwwM7v7mzCjeyIzMzOjVq1aBtufo6NjtbwBHoZcs/KR61V+cs3Kpzpfr3s9id0mjT2EEEKYNElkQgghTFqVTWTW1tZMnjwZa2trtUMxGXLNykeuV/nJNSsfuV5lY3SNPYQQQojyqLJPZEIIIaoHSWRCCCFMmiQyIYQQJk0SmRBCCJNWJRPZZ599RmBgIDY2NrRq1YrY2Fi1QzJaU6ZMQaPRlJhCQkLUDsuo7Nixg549e+Lj44NGo+Hnn38usV5RFN577z28vb2xtbUlIiKCs2fPqhOsEbjf9Ro+fHipe65bt27qBGskpk2bRosWLXBwcMDDw4M+ffpw+vTpEmXy8/OJjIzEzc2NGjVq8Mwzz5CSkqJSxMalyiWy1atXM378eCZPnsyhQ4do2rQpXbt2JTU1Ve3QjFbDhg1JSkrST7t27VI7JKOSk5ND06ZN+eyzz+64fvr06cydO5cvvviCffv2YW9vT9euXcnPz6/kSI3D/a4XQLdu3UrccytXrqzECI1PdHQ0kZGR7N27l82bN1NUVESXLl3IycnRl3nttddYt24d33//PdHR0Vy9epWnn35axaiNiFLFtGzZUomMjNTPa7VaxcfHR5k2bZqKURmvyZMnK02bNlU7DJMBKGvWrNHP63Q6xcvLS5kxY4Z+WXp6umJtba2sXLlShQiNyz+vl6IoyrBhw5TevXurEo+pSE1NVQAlOjpaUZRb95SlpaXy/fff68vExcUpgBITE6NWmEajSj2RFRYWcvDgQSIiIvTLzMzMiIiIICYmRsXIjNvZs2fx8fGhdu3aDB48mEuXLqkdkslISEggOTm5xD3n5OREq1at5J67h99//x0PDw/q16/Pyy+/zPXr19UOyahkZGQA4OrqCsDBgwcpKioqcZ+FhITg7+8v9xlVrGrx2rVraLVaPD09Syz39PQkOTlZpaiMW6tWrVi8eDEbNmxgwYIFJCQk8Oijj+qH0xH3dvu+knuu7Lp168bSpUvZunUrn3zyCdHR0XTv3h2tVqt2aEZBp9MRFRVFu3btaNSoEXDrPrOyssLZ2blEWbnPbjG63u9F5erevbv+z02aNKFVq1YEBATw3XffMWLECBUjE1XVwIED9X9u3LgxTZo0oU6dOvz+++906tRJxciMQ2RkJMePH5d31eVQpZ7Iatasibm5eamWPCkpKXh5eakUlWlxdnamXr16xMfHqx2KSbh9X8k99+Bq165NzZo15Z4Dxo4dy//+9z+2b99eYjgrLy8vCgsLSU9PL1Fe7rNbqlQis7KyIiwsjK1bt+qX6XQ6tm7dSps2bVSMzHRkZ2dz7tw5vL291Q7FJAQFBeHl5VXinsvMzGTfvn1yz5XR5cuXuX79erW+5xRFYezYsaxZs4Zt27YRFBRUYn1YWBiWlpYl7rPTp09z6dIluc+oglWL48ePZ9iwYYSHh9OyZUtmz55NTk4OL7zwgtqhGaXXX3+dnj17EhAQwNWrV5k8eTLm5uYMGjRI7dCMRnZ2domnhYSEBI4cOYKrqyv+/v5ERUXxwQcfEBwcTFBQEO+++y4+Pj706dNHvaBVdK/r5erqytSpU3nmmWfw8vLi3LlzTJgwgbp169K1a1cVo1ZXZGQkK1asYO3atTg4OOjfezk5OWFra4uTkxMjRoxg/PjxuLq64ujoyCuvvEKbNm1o3bq1ytEbAbWbTVaEefPmKf7+/oqVlZXSsmVLZe/evWqHZLQGDBigeHt7K1ZWVoqvr68yYMAAJT4+Xu2wjMr27dsVoNQ0bNgwRVFuNcF/9913FU9PT8Xa2lrp1KmTcvr0aXWDVtG9rldubq7SpUsXxd3dXbG0tFQCAgKUUaNGKcnJyWqHrao7XS9AWbRokb5MXl6eMmbMGMXFxUWxs7NT+vbtqyQlJakXtBGRYVyEEEKYtCr1jkwIIUT1I4lMCCGESZNEJoQQwqRJIhNCCGHSJJEJIYQwaZLIhBBCmDRJZEIIIUyaJDIhhBAmTRKZEBWoQ4cOREVFqR2GEFWaJDIhhBAmTRKZEFVMYWGh2iEIUakkkQlRwXQ6HRMmTMDV1RUvLy+mTJmiX3fp0iV69+5NjRo1cHR0pH///iXGNhs+fHipXvSjoqLo0KGDfr5Dhw6MHTuWqKgoatasWa17kRfVkyQyISrYkiVLsLe3Z9++fUyfPp1///vfbN68GZ1OR+/evblx4wbR0dFs3ryZ8+fPM2DAgAc6hpWVFbt37+aLL76ogLMQwnhVufHIhDA2TZo0YfLkyQAEBwczf/58/QCJx44dIyEhAT8/PwCWLl1Kw4YN2b9/Py1atCjzMYKDg5k+fbrhgxfCBMgTmRAVrEmTJiXmvb29SU1NJS4uDj8/P30SAwgNDcXZ2Zm4uLhyHSMsLMwgsQphiiSRCVHBLC0tS8xrNBp0Ol2ZtjUzM+OfQwYWFRWVKmdvb//gAQph4iSRCaGSBg0akJiYSGJion7ZyZMnSU9PJzQ0FAB3d3eSkpJKbHfkyJHKDFMIoyeJTAiVRERE0LhxYwYPHsyhQ4eIjY1l6NChPP7444SHhwPwxBNPcODAAZYuXcrZs2eZPHkyx48fVzlyIYyLJDIhVKLRaFi7di0uLi489thjREREULt2bVavXq0v07VrV959910mTJhAixYtyMrKYujQoSpGLYTx0Sj/rIAXQgghTIg8kQkhhDBpksiEEEKYNElkQgghTJokMiGEECZNEpkQQgiTJolMCCGESZNEJoQQwqRJIhNCCGHSJJEJIYQwaZLIhBBCmDRJZEIIIUza/wOzzPO9KH7OiQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 500x200 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# sanity check\n",
    "temp = gen_supply_score(trip_cleaned.sample(frac=0.05)).groupby(\"hour\").max().reset_index()\n",
    "print(temp.head())\n",
    "\n",
    "temp['hour'] = temp.hour.astype(int)\n",
    "plt.rcParams[\"figure.figsize\"] = (5,2)\n",
    "temp.sort_values(\"hour\").plot.line(x='hour', y='supply_score')\n",
    "plt.title(\"Supply Score\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0301033a-9c24-49c1-9873-ae17c1c66aba",
   "metadata": {},
   "source": [
    "Observation:\n",
    "* Scores seems reasonable - peaked during working hour and dropped during midnight"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0c2c4008-51e2-4329-a063-30c02ce4efca",
   "metadata": {
    "tags": []
   },
   "source": [
    "## Feature: Traffic score \n",
    "Reminder: Ensure no data leak (aggregated data score should not appear in test period)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "299e3e60-d754-4a4c-affe-9ba2ebbc9806",
   "metadata": {},
   "outputs": [],
   "source": [
    "def gen_traffic_score(df):\n",
    "    \"\"\"\n",
    "    Generate the table for traffic_score(day_of_week, hour)\n",
    "        > higher means traffic is less congested\n",
    "        \n",
    "    P/S: Extension work - call Maps API to get the route from pickup to dropoff (so that traffic can be scored on a roadID level)\n",
    "    \"\"\"\n",
    "    print(\"Generating traffic score ...\")\n",
    "    \n",
    "    REQUIRED_COLS = [\"Trip ID\", \"Trip Start Timestamp\", \"Trip Seconds\", \"Trip Miles\"]\n",
    "    assert all(c in df.columns for c in REQUIRED_COLS), f\"Did not have all the required columns: {REQUIRED_COLS}\"\n",
    "    \n",
    "    # transform table to get the necessary metadata\n",
    "    df = df.loc[:, REQUIRED_COLS].dropna().copy()\n",
    "    df = df.loc[df['Trip Seconds']!=0, ]\n",
    "    df['day_of_week'] = df['Trip Start Timestamp'].map(lambda x: int(datetime.strftime(x,\"%w\")))\n",
    "    df['hour'] = df['Trip Start Timestamp'].map(lambda x: int(datetime.strftime(x,\"%-H\")))\n",
    "    df['avg_trip_speed'] = df['Trip Miles'] / df['Trip Seconds']\n",
    "\n",
    "    \n",
    "    # aggregate to get scores \n",
    "    traffic_score = df.groupby([\"day_of_week\", \"hour\"])\\\n",
    "    .agg(n_samples = (\"avg_trip_speed\", \"count\"), traffic_score = (\"avg_trip_speed\", \"mean\")).reset_index()\n",
    "    \n",
    "    # filter out states w/o sufficient samples\n",
    "    traffic_score = traffic_score.loc[traffic_score.n_samples >= 30, traffic_score.columns!=\"n_samples\"].reset_index(drop=True)\n",
    "    \n",
    "    return traffic_score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "58a618c1-d72a-4cc4-b0fc-12b941110bff",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Generating traffic score ...\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>day_of_week</th>\n",
       "      <th>hour</th>\n",
       "      <th>traffic_score</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>149</th>\n",
       "      <td>6</td>\n",
       "      <td>5</td>\n",
       "      <td>0.036008</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>46</th>\n",
       "      <td>1</td>\n",
       "      <td>22</td>\n",
       "      <td>0.020970</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>95</th>\n",
       "      <td>3</td>\n",
       "      <td>23</td>\n",
       "      <td>0.015931</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>121</th>\n",
       "      <td>5</td>\n",
       "      <td>1</td>\n",
       "      <td>0.014406</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>84</th>\n",
       "      <td>3</td>\n",
       "      <td>12</td>\n",
       "      <td>0.012402</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "     day_of_week  hour  traffic_score\n",
       "149            6     5       0.036008\n",
       "46             1    22       0.020970\n",
       "95             3    23       0.015931\n",
       "121            5     1       0.014406\n",
       "84             3    12       0.012402"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# sanity check\n",
    "gen_traffic_score(trip_cleaned.sample(frac=0.05)).sort_values(\"traffic_score\", ascending=False).head()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c64c7bee-a1ec-4869-9b34-0c82e3eeabe0",
   "metadata": {},
   "source": [
    "****"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e1e6a1ff-e2ce-453c-ae59-74c6fa4e3d0c",
   "metadata": {},
   "source": [
    "# Modelling"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "df6f71ba-b09b-448a-b469-172b11d0c4b6",
   "metadata": {},
   "outputs": [],
   "source": [
    "def add_trip_score(df, demand_score, supply_score, traffic_score):\n",
    "    \"\"\"\n",
    "    Add demand, supply, traffic scores to a booking info.\n",
    "    \"\"\"\n",
    "    \n",
    "    BOOKING_COLS = [\"Trip Start Timestamp\", \"Pickup Centroid Latitude\", \"Pickup Centroid Longitude\", \"Dropoff Centroid Latitude\", \"Dropoff Centroid Longitude\"]\n",
    "    assert all(c in df.columns for c in BOOKING_COLS), f\"Did not have all the required columns: {BOOKING_COLS}\"\n",
    "    \n",
    "    # transform to get required metadata\n",
    "    df = df.copy()\n",
    "    df['day_of_week'] = df['Trip Start Timestamp'].map(lambda x: int(datetime.strftime(x,\"%w\")))\n",
    "    df['hour'] = df['Trip Start Timestamp'].map(lambda x: int(datetime.strftime(x,\"%-H\")))\n",
    "    df['pickup_geohash'] = df.apply(lambda x: get_geohash(x['Pickup Centroid Latitude'], x['Pickup Centroid Longitude']),axis=1)\n",
    "    df['dropoff_geohash'] = df.apply(lambda x: get_geohash(x['Dropoff Centroid Latitude'], x['Dropoff Centroid Longitude']),axis=1)\n",
    "    \n",
    "    # left join to add scores\n",
    "    df = pd.merge(df, demand_score, \n",
    "                  left_on=[\"day_of_week\", \"hour\", \"pickup_geohash\"], right_on=[\"day_of_week\", \"hour\", \"geohash\"], how=\"left\")\n",
    "    df = pd.merge(df, supply_score, \n",
    "                  left_on=[\"day_of_week\", \"hour\", \"pickup_geohash\"], right_on=[\"day_of_week\", \"hour\", \"geohash\"], how=\"left\")\n",
    "    df = pd.merge(df, traffic_score, \n",
    "                  on=[\"day_of_week\", \"hour\"], how=\"left\")\n",
    "    \n",
    "    # impute missing values \n",
    "    mean_demand_score = demand_score['demand_score'].mean()\n",
    "    mean_supply_score = supply_score['supply_score'].mean()\n",
    "    mean_traffic_score = traffic_score['traffic_score'].mean()\n",
    "    \n",
    "    df['demand_score'] = df.demand_score.fillna(mean_demand_score)\n",
    "    df['supply_score'] = df.supply_score.fillna(mean_demand_score)\n",
    "    df['traffic_score'] = df.traffic_score.fillna(mean_demand_score)\n",
    "    \n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "5e4d75df-2434-43e0-953a-fc1cb25ffa42",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "---- Train-test for CV Fold 1 -----\n",
      "Generating demand score ...\n",
      "Generating supply score ...\n",
      "Generating traffic score ...\n",
      "---- Train-test for CV Fold 2 -----\n",
      "Generating demand score ...\n",
      "Generating supply score ...\n",
      "Generating traffic score ...\n",
      "---- Train-test for CV Fold 3 -----\n",
      "Generating demand score ...\n",
      "Generating supply score ...\n",
      "Generating traffic score ...\n",
      "---- Train-test for CV Fold 4 -----\n",
      "Generating demand score ...\n",
      "Generating supply score ...\n",
      "Generating traffic score ...\n",
      "---- Train-test for CV Fold 5 -----\n",
      "Generating demand score ...\n",
      "Generating supply score ...\n",
      "Generating traffic score ...\n"
     ]
    }
   ],
   "source": [
    "feat_cols_baseline = [\"Trip Miles\", \"day_of_week\", \"hour\"]\n",
    "feat_cols_v1 = [\"Trip Miles\",\"day_of_week\", \"hour\", \"demand_score\", \"supply_score\", \"traffic_score\"]\n",
    "\n",
    "cv_mse_baseline = []\n",
    "cv_mse_v1 = []\n",
    "\n",
    "tscv = TimeSeriesSplit(n_splits=5, )\n",
    "trip_cleaned = trip_cleaned.sort_values(\"Trip Start Timestamp\").reset_index(drop=True)\n",
    "\n",
    "# training-testing for 5 CV folds \n",
    "for i, (train_index, test_index) in enumerate(tscv.split(trip_cleaned)):\n",
    "    print(f\"---- Train-test for CV Fold {i+1} -----\")\n",
    "    ############### Training ###############\n",
    "    train = trip_cleaned.iloc[train_index, :]\n",
    "    \n",
    "    # create score table\n",
    "    demand_score_table = gen_demand_score(train)\n",
    "    supply_score_table = gen_supply_score(train)\n",
    "    traffic_score_table = gen_traffic_score(train)\n",
    "    \n",
    "    # add demand, supply, traffic score\n",
    "    train = add_trip_score(train, demand_score_table, supply_score_table, traffic_score_table)\n",
    "    \n",
    "    # model training\n",
    "    X_train = train.loc[:, :]\n",
    "    y_train = train.loc[:, \"Fare\"]\n",
    "    \n",
    "    # fit\n",
    "    # [v0] baseline \n",
    "    model_baseline = XGBRegressor()\n",
    "    model_baseline.fit(X_train.loc[:, feat_cols_baseline], y_train)\n",
    "    \n",
    "    # [v1] (with supply, demand, traffic scores)\n",
    "    model_v1 = XGBRegressor()\n",
    "    model_v1.fit(X_train.loc[:, feat_cols_v1], y_train)\n",
    "    \n",
    "    ############### Testing ############### \n",
    "    test = trip_cleaned.iloc[test_index, :]\n",
    "    \n",
    "    # add demand, supply, traffic score\n",
    "    test = add_trip_score(test, demand_score_table, supply_score_table, traffic_score_table)\n",
    "    \n",
    "    # model testing\n",
    "    X_test = test.loc[:, :]\n",
    "    y_test = test.loc[:, \"Fare\"]\n",
    "    \n",
    "    y_pred_baseline = model_baseline.predict(X_test.loc[:, feat_cols_baseline])\n",
    "    y_pred_v1 = model_v1.predict(X_test.loc[:, feat_cols_v1])\n",
    "    \n",
    "    cv_mse_baseline.append(mean_squared_error(y_test, y_pred_baseline))\n",
    "    cv_mse_v1.append(mean_squared_error(y_test, y_pred_v1))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6950b86d-0c4b-4f1f-85fa-920f713eb96e",
   "metadata": {},
   "source": [
    "# Results & Error Analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "b183da56-54a5-48b8-8863-6823b5977dfa",
   "metadata": {},
   "outputs": [],
   "source": [
    "res_df = pd.DataFrame()\n",
    "res_df['baseline'] = cv_mse_baseline\n",
    "res_df['v1'] = cv_mse_v1\n",
    "res_df['cv_fold'] = list(range(1, res_df.shape[0] + 1))\n",
    "\n",
    "res_df = pd.melt(res_df, id_vars=['cv_fold'], value_vars=['baseline', 'v1']).rename({\"value\": \"mse\"}, axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "id": "a40f726a-2c18-427f-af7c-4815f6d85cd1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhIAAAE8CAYAAACVc1hkAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy89olMNAAAACXBIWXMAAA9hAAAPYQGoP6dpAABF0klEQVR4nO3deXxM1/8/8NdkmewLspMQoSIaS1GCxL6EWiotgoqlKLFVVaSlEtTS2velRaVSS22tFkUtsUQJWkvliwYhicSSlSwy5/eHX+7HSMLMmGRmeD0fjzy45557znvu3Jm8c+6598qEEAJEREREGjDSdQBERERkuJhIEBERkcaYSBAREZHGmEgQERGRxphIEBERkcaYSBAREZHGmEgQERGRxphIEBERkcaYSBAREZHGmEioQSaTISIiQu3tbty4AZlMhvXr12s9plcRFRUFb29vmJqawt7eXtfhEOmlb775Bt7e3lAoFFKZpt8F+ubw4cOQyWQ4fPjwK7fVqlUrvP32268elB5o2rQpJk6cqOswDIbBJRLr16+HTCaDTCbDsWPHiq0XQsDd3R0ymQzvvfeeDiLUXNGHuujH1NQU1atXx4ABA/Dff/9pta8rV65g4MCB8PLywpo1a7B69Wqttv+miYiIgEwmw71793QdCpUiKSkJEREROH/+vMrbZGZmYs6cOQgLC4ORkcF9XRq8mTNnYufOncXKT5w4gYiICKSnp5dJv2FhYVi2bBlSUlLKpP3XjcF+MszNzREdHV2s/MiRI7h9+zbMzMx0EJV2jBkzBlFRUVi9ejW6dOmCzZs3o3HjxkhKStJaH4cPH4ZCocCiRYswcOBA9OrVS2ttE+mjpKQkREZGqpVIrF27Fk+ePEFwcHDZBaZDAQEBePz4MQICAnQdSolelEhERkaWWSLRvXt32NraYvny5WXS/uvGYBOJzp07Y+vWrXjy5IlSeXR0NBo2bAgXFxcdRfbq/P390b9/fwwaNAhLlizB3Llz8eDBA/zwww+v3HZOTg4AIDU1FQC0ekrj0aNHWmuLXo1CoUBubq6uwzB469atQ7du3WBubq7rUDRW9JkviZGREczNzTna8v8VfYcZGRnhgw8+wIYNG8DnWr6cwR49wcHBuH//Pvbv3y+V5efn4+eff0bfvn1L3CYnJwefffYZ3N3dYWZmhlq1amHu3LnFDpS8vDx8+umncHR0hI2NDbp164bbt2+X2OadO3cwePBgODs7w8zMDHXq1MHatWu190IBtGnTBgCQkJAgle3Zswf+/v6wsrKCjY0NunTpgkuXLiltN3DgQFhbW+P69evo3LkzbGxs0K9fP1SrVg1Tp04FADg6OhY737t8+XLUqVMHZmZmcHNzQ2hoaLHMv+h8aFxcHAICAmBpaYkvvvhCmg8yd+5cLFu2DNWrV4elpSU6dOiAxMRECCEwffp0VKlSBRYWFujevTsePHig1PauXbvQpUsXuLm5wczMDF5eXpg+fToKCwtLjOHy5cto3bo1LC0tUblyZXzzzTfF9mFubi4iIiLw1ltvwdzcHK6urujZsyeuX78u1VEoFFi4cCHq1KkDc3NzODs7Y/jw4Xj48KHqb1YJ8f3zzz9o2bIlLC0tUaNGDfz8888Ano6eNWnSBBYWFqhVqxYOHDigtH3R6ZIrV66gV69esLW1RaVKlTB27NhiSYJMJsOoUaOwceNG6b3bu3cvAODcuXMIDAyEra0trK2t0bZtW8TGxkrbnjlzBjKZrMREdd++fZDJZNi9e7dUpsoxX3SabsuWLYiMjETlypVhY2ODDz74ABkZGcjLy8O4cePg5OQEa2trDBo0CHl5ecX6//HHH9GwYUNYWFigYsWK6NOnDxITE0vczy86Dg4fPozGjRsDAAYNGiSdPnzRvKWEhAT8888/aNeuXal1nvWy/Zyeng5jY2MsXrxYKrt37x6MjIxQqVIlpe+hESNGFPtj6NSpU+jUqRPs7OxgaWmJli1b4vjx40p1io6Zy5cvo2/fvqhQoQJatGhRaswlzZG4evUqgoKC4OLiAnNzc1SpUgV9+vRBRkaGSvshLi4OzZo1g4WFBTw9PbFy5cpidfLy8jB16lTUqFEDZmZmcHd3x8SJE5WOAZlMhpycHPzwww/S+zVw4EBERETg888/BwB4enpK627cuCFtq85x8/x3WJH27dvj5s2bao1gvbGEgVm3bp0AIE6fPi2aNWsmPvroI2ndzp07hZGRkbhz546oWrWq6NKli7ROoVCINm3aCJlMJj7++GOxdOlS0bVrVwFAjBs3TqmP/v37CwCib9++YunSpaJnz56ibt26AoCYOnWqVC8lJUVUqVJFuLu7i2nTpokVK1aIbt26CQBiwYIFUr2EhAQBQKxbt+6Fr+3QoUMCgNi6datS+a5duwQAMWnSJCGEEBs2bBAymUx06tRJLFmyRMyZM0dUq1ZN2Nvbi4SEBGm7kJAQYWZmJry8vERISIhYuXKl2LBhg9ixY4d4//33BQCxYsUKERUVJf7++28hhBBTp04VAES7du3EkiVLxKhRo4SxsbFo3LixyM/Pl9pu2bKlcHFxEY6OjmL06NFi1apVYufOndJrrV+/vvDx8RHz588XkydPFnK5XDRt2lR88cUXolmzZmLx4sVizJgxQiaTiUGDBim93h49eohevXqJb7/9VqxYsUJ8+OGHAoCYMGGCUr2WLVsKNzc34e7uLsaOHSuWL18u2rRpIwCI33//Xar35MkT0bZtWwFA9OnTRyxdulTMmjVLtGnTRuzcuVOq9/HHHwsTExMxdOhQsXLlShEWFiasrKyKvfaSFO23tLS0EuP7/PPPxZIlS4SPj48wNjYWmzZtEi4uLiIiIkIsXLhQVK5cWdjZ2YnMzMxibfr6+oquXbuKpUuXSsfms8e9EEIAELVr1xaOjo4iMjJSLFu2TJw7d05cvHhRWFlZCVdXVzF9+nQxe/Zs4enpKczMzERsbKy0ffXq1UXnzp2Lva5BgwaJChUqSK9f1WO+6FiuX7++8PPzU3q/+/TpI/r27SsCAwPFsmXLxEcffSQAiMjISKW+Z8yYIWQymejdu7dYvny5iIyMFA4ODqJatWri4cOHah0HKSkpYtq0aQKAGDZsmIiKihJRUVHi+vXrpb6nP/74owAg/vnnn2Lrnv8uUHU/161bVwQFBUnLO3bsEEZGRgKAuHjxolRep04d8cEHH0jLBw8eFHK5XPj5+Yl58+aJBQsWiLp16wq5XC5OnTol1Ss6Znx8fET37t3F8uXLxbJly0p9jUXv06FDh4QQQuTl5QlPT0/h5uYmZsyYIb777jsRGRkpGjduLG7cuFFqO0L8731wcnISo0aNEosXLxYtWrQQAMT3338v1SssLBQdOnQQlpaWYty4cWLVqlVi1KhRwsTERHTv3l2qFxUVJczMzIS/v7/0fp04cUL8/fffIjg4WDrmitZlZ2cLIdQ7bkr6Dity+/ZtAUAsWbLkha+bhDDoRGLp0qXCxsZGPHr0SAghxIcffihat24thBDFEomdO3cKAGLGjBlK7X3wwQdCJpOJa9euCSGEOH/+vAAgRo4cqVSvb9++xb48hgwZIlxdXcW9e/eU6vbp00fY2dlJcambSKxdu1akpaWJpKQk8dtvv4lq1aoJmUwmTp8+LbKysoS9vb0YOnSo0rYpKSnCzs5OqTwkJEQpAXlWSb/4UlNThVwuFx06dBCFhYVS+dKlS6W4irRs2VIAECtXrlRqt+i1Ojo6ivT0dKk8PDxcABD16tUTBQUFUnlwcLCQy+UiNzdXKivab88aPny4sLS0VKpXFMOGDRuksry8POHi4qL0Zb127VoBQMyfP79YuwqFQgghRExMjAAgNm7cqLR+7969JZY/r7REAoCIjo6Wyq5cuSIACCMjI6VfMPv27St2jBS12a1bN6W+Ro4cKQBIyZ8QQmrz0qVLSnV79Ogh5HK50i/MpKQkYWNjIwICAqSy8PBwYWpqKh48eCCV5eXlCXt7ezF48GCpTNVjvuhYfvvtt5WSsODgYCGTyURgYKDS9n5+fqJq1arS8o0bN4SxsbH4+uuvlepduHBBmJiYKJWrehycPn1apc9hkcmTJwsAIisrq9i6578LVN3PoaGhwtnZWVoeP368CAgIEE5OTmLFihVCCCHu378vZDKZWLRokRDi6TFas2ZN0bFjR+l4FeLp58TT01O0b99eKis6ZoKDg1V6jc8nEufOnSvxjxlVFL0P8+bNk8ry8vJE/fr1hZOTk3QcREVFCSMjIxETE6O0/cqVKwUAcfz4canMyspKhISEFOvr22+/FQCU/nASQrPj5vnvsGfJ5XIxYsSIl772N53BntoAgF69euHx48fYvXs3srKysHv37lJPa/z+++8wNjbGmDFjlMo/++wzCCGwZ88eqR6AYvXGjRuntCyEwLZt29C1a1cIIXDv3j3pp2PHjsjIyMDZs2c1el2DBw+Go6Mj3Nzc0KVLF2l4r1GjRti/fz/S09MRHBys1KexsTGaNGmCQ4cOFWtvxIgRKvV74MAB5OfnY9y4cUrnTIcOHQpbW1v89ttvSvXNzMwwaNCgEtv68MMPYWdnJy03adIEANC/f3+YmJgolefn5+POnTtSmYWFhfT/rKws3Lt3D/7+/nj06BGuXLmi1I+1tTX69+8vLcvlcrz77rtKV7ls27YNDg4OGD16dLE4ZTIZAGDr1q2ws7ND+/btlfZrw4YNYW1tXeJ+VYW1tTX69OkjLdeqVQv29vaoXbu2tE+K9gOAEq/OCQ0NVVoueh1Fx2qRli1bwsfHR1ouLCzEH3/8gR49eqB69epSuaurK/r27Ytjx44hMzMTANC7d28UFBRg+/btUr0//vgD6enp6N27NwDNjvkBAwbA1NRU6XUKITB48GClek2aNEFiYqI052n79u1QKBTo1auXUj8uLi6oWbNmsfdDleNAXffv34eJiQmsra1fWE+d/ezv74+7d+8iPj4eABATE4OAgAD4+/sjJiYGAHDs2DEIIeDv7w8AOH/+PK5evYq+ffvi/v370r7IyclB27ZtcfToUaVLUwHgk08+0eg1F31m9+3bp9GcJxMTEwwfPlxalsvlGD58OFJTUxEXFwfg6Wetdu3a8Pb2Vnpvi07havpZA9Q/bl70HQYAFSpU4JVYKjB5eRX95ejoiHbt2iE6OhqPHj1CYWEhPvjggxLr3rx5E25ubrCxsVEqr127trS+6F8jIyN4eXkp1atVq5bSclpaGtLT07F69epSL50smtCorq+++gr+/v4wNjaGg4MDateuLf3yvXr1KoD/zZt4nq2trdKyiYkJqlSpolK/Rfvg+dcql8tRvXp1aX2RypUrQy6Xl9iWh4eH0nLRF5S7u3uJ5c/OQ7h06RImT56MP//8U/oCLvL8edoqVapIyUCRChUq4J9//pGWr1+/jlq1aiklMM+7evUqMjIy4OTkVOJ6Td/LkuKzs7NTaT8UqVmzptKyl5cXjIyMlM4JA0/PFz8rLS0Njx49KvZ+Ak+Pe4VCgcTERNSpUwf16tWDt7c3Nm/ejCFDhgAANm/eDAcHB+lY0+SYV+c4UCgUyMjIQKVKlXD16lUIIYq99iLPJieAasdBWVFnPxclBzExMahSpQrOnTuHGTNmwNHREXPnzpXW2draol69egD+95kPCQkpNYaMjAxUqFBBWn7+WFCVp6cnxo8fj/nz52Pjxo3w9/dHt27d0L9/f6U/DErj5uYGKysrpbK33noLwNP76TRt2hRXr17Fv//+C0dHxxLb0PSzBkDt4+ZF32HA0+T5+eOKijPoRAIA+vbti6FDhyIlJQWBgYHldmOlor8A+vfvX+oHvG7duhq17evrW+oEr6J+o6KiSrwy5flflmZmZmU2I/vZkYPnGRsbq1Uu/v9Es/T0dLRs2RK2traYNm0avLy8YG5ujrNnzyIsLKzYX14va09VCoUCTk5O2LhxY4nrS/vSexlN98OLlPbF9qL3QxW9e/fG119/jXv37sHGxga//PILgoODpWNKk2Ne09evUCggk8mwZ8+eEus+P0qgrePgWZUqVcKTJ0+QlZVV7A8QTbm5ucHT0xNHjx5FtWrVIISAn58fHB0dMXbsWNy8eRMxMTFo1qyZ9Lkt2u/ffvst6tevX2K7z++PVzkW5s2bh4EDB2LXrl34448/MGbMGMyaNQuxsbEq/1HyIgqFAr6+vpg/f36J659PMtVtW53j5mX7KT09HQ4ODhrH86Yw+ETi/fffx/DhwxEbG4vNmzeXWq9q1ao4cOBAsS+FoqHyqlWrSv8qFArpr9giRUORRYqu6CgsLFR5Vrc2FI2UODk5ab3fon0QHx+vNESbn5+PhISEcnmdhw8fxv3797F9+3ala9ufvWJFXV5eXjh16hQKCgqK/UXybJ0DBw6gefPmr/wLWduuXr2q9BfmtWvXoFAoUK1atRdu5+joCEtLy2LHLvD0uDcyMlL60u7duzciIyOxbds2ODs7IzMzU+m0THke815eXhBCwNPTU/qL9lWp+5elt7c3gKfH3ov+KFB3P/v7++Po0aPw9PRE/fr1YWNjg3r16sHOzg579+7F2bNnERkZKdUv+szb2tqW23eNr68vfH19MXnyZJw4cQLNmzfHypUrMWPGjBdul5SUhJycHKVRif/7v/8DAOl49fLywt9//422bdu+9D0pbX1p5do8bu7cuYP8/Hxp1JpKZ9BzJICnGeaKFSsQERGBrl27llqvc+fOKCwsxNKlS5XKFyxYAJlMhsDAQACQ/n32Ei0AWLhwodKysbExgoKCsG3bNly8eLFYf2lpaZq8nJfq2LEjbG1tMXPmTBQUFGi133bt2kEul2Px4sVKf8l9//33yMjIQJcuXTRuW1VFf0U8239+fv4r3RgmKCgI9+7dK/beP9tPr169UFhYiOnTpxer8+TJkzK78Y0qli1bprS8ZMkSAP87VktjbGyMDh06YNeuXUqnQe7evYvo6Gi0aNFC6VRY7dq14evri82bN2Pz5s1wdXVVSubK85jv2bMnjI2NERkZWWxUQQiB+/fvq91m0S83Vd9LPz8/AE8vj30Rdfezv78/bty4gc2bN0unOoyMjNCsWTPMnz8fBQUFUjkANGzYEF5eXpg7dy6ys7OL9a/N/Z6ZmVns3jy+vr4wMjIq8fLc5z158gSrVq2SlvPz87Fq1So4OjqiYcOGAJ5+1u7cuYM1a9YU2/7x48dK972wsrIq8f0q7b3U5nFTNKejWbNmKm/zpjL4EQngxecOi3Tt2hWtW7fGl19+iRs3bqBevXr4448/sGvXLowbN07K+uvXr4/g4GAsX74cGRkZaNasGQ4ePIhr164Va3P27Nk4dOgQmjRpgqFDh8LHxwcPHjzA2bNnceDAgWL3R9AGW1tbrFixAh999BHeeecd9OnTB46Ojrh16xZ+++03NG/evMRfmKpwdHREeHg4IiMj0alTJ3Tr1g3x8fFYvnw5GjdurDSZraw0a9YMFSpUQEhICMaMGQOZTIaoqKhXGqIeMGAANmzYgPHjx+Ovv/6Cv78/cnJycODAAYwcORLdu3dHy5YtMXz4cMyaNQvnz59Hhw4dYGpqiqtXr2Lr1q1YtGhRqfNvylpCQgK6deuGTp064eTJk/jxxx/Rt29f6Rz6i8yYMQP79+9HixYtMHLkSJiYmGDVqlXIy8sr8X4bvXv3xldffQVzc3MMGTKk2Gmx8jrmvby8MGPGDISHh+PGjRvo0aMHbGxskJCQgB07dmDYsGGYMGGC2m3a29tj5cqVsLGxgZWVFZo0aVLqfILq1avj7bffxoEDB4pNDn2eOvu5KEmIj4/HzJkzpfKAgADs2bMHZmZm0j0vgKdJxnfffYfAwEDUqVMHgwYNQuXKlXHnzh0cOnQItra2+PXXX9XaF6X5888/MWrUKHz44Yd466238OTJE0RFRUlJ5Mu4ublhzpw5uHHjBt566y1s3rwZ58+fx+rVq6XRwI8++ghbtmzBJ598gkOHDqF58+YoLCzElStXsGXLFuzbtw+NGjUC8DSJOnDgAObPny+dFmrSpImUlHz55Zfo06cPTE1N0bVrV60eN/v374eHhwcaNGig4d58g5TX5SHa8uzlny/y/OWfQgiRlZUlPv30U+Hm5iZMTU1FzZo1xbfffqt0SZUQQjx+/FiMGTNGVKpUSVhZWYmuXbuKxMTEYpd8CSHE3bt3RWhoqHB3dxempqbCxcVFtG3bVqxevVqq86r3kSitbseOHYWdnZ0wNzcXXl5eYuDAgeLMmTNSnZCQEGFlZVXi9iVdrlhk6dKlwtvbW5iamgpnZ2cxYsQIpeuvhXh66VSdOnWKbVv0Wr/99luVXltJ7+fx48dF06ZNhYWFhXBzcxMTJ06ULo8sukztRTGEhIQoXUooxNNL5b788kvh6ekpvU8ffPBBsfsIrF69WjRs2FBYWFgIGxsb4evrKyZOnCiSkpKK9fOs0i7/LCm+ko5NIZ5eUhgaGlqszcuXL4sPPvhA2NjYiAoVKohRo0aJx48fv3DbZ509e1Z07NhRWFtbC0tLS9G6dWtx4sSJEutevXpVABAAxLFjx0qso8oxr877/exrff543LZtm2jRooWwsrISVlZWwtvbW4SGhor4+HipjjrHwa5du4SPj48wMTFR6TM5f/58YW1tXeyS5JK+C9TZz05OTgKAuHv3rlR27NgxAUD4+/uXuM25c+dEz549RaVKlYSZmZmoWrWq6NWrlzh48KBU50Wf65I8f/nnf//9JwYPHiy8vLyEubm5qFixomjdurU4cODAS9sqeh/OnDkj/Pz8hLm5uahatapYunRpsbr5+flizpw5ok6dOsLMzExUqFBBNGzYUERGRoqMjAyp3pUrV0RAQICwsLAQAJQuBZ0+fbqoXLmydB+OZy8FfZXjRoin97pwdXUVkydPfunrJiFkQvD+n0T6KCIiApGRkUhLS+OELx3JyMhA9erV8c0330hXs9Drb+fOnejbty+uX78OV1dXXYej9wx+jgQRUVmxs7PDxIkT8e233xa7YoheX3PmzMGoUaOYRKjotZgjQURUVsLCwhAWFqbrMKgcnTx5UtchGBSOSBAREZHGOEeCiIiINMYRCSIiItIYEwkiIiLS2Gs/2VKhUCApKQk2NjZ8+AoR0RtKCIGsrCy4ubmV2fOH3lSvfSKRlJT0Sg+BISKi10diYqJWHj5G//PaJxJFD+hKTEws9ohtIiJ6M2RmZsLd3V1rT3Kl/3ntE4mi0xm2trZMJIiI3nA8xa19PFFEREREGmMiQURERBpjIkFEREQae+3nSBARkX4pLCxEQUFBufaZn5+PqlWrIj8/H7m5ueXatyEyNjaGiYmJSnNKXvtbZGdmZsLOzg4ZGRmcbElEpGPZ2dm4ffs2yvtXj0KhQGJiItzd3XkfCRVZWlrC1dUVcrn8hfU4IkFEROWisLAQt2/fhqWlJRwdHcv1CorCwkI8fvwY1apVg7Gxcbn1a4iEEMjPz0daWhoSEhJQs2bNFyZfTCSIiKhcFBQUQAgBR0dHWFhYlGvfhYWFAABzc3MmEiqwsLCAqakpbt68ifz8fJibm5dal+M7RERUrngvB8Og6ikgjkiQRo4dO6azvlu0aKGzvomISBlHJIiIiEhjHJEgKkccySEqrjw/FydPnlT7s9CqVSvUr18fCxcuLJugXmLgwIFIT0/Hzp079SKe5zGRICIiMiDbt2+HqamprsOQMJEgIiIyIBUrVtR1CEqYSIDDzURE9GJPnjzBqFGjEBUVBVNTU4wYMQLTpk2DTCZDVFQUFi1ahPj4eFhZWaFNmzZYuHAhnJycAAAPHz7EqFGj8McffyA7OxtVqlTBF198gUGDBgEAEhMT8dlnn+GPP/6AkZER/P39sWjRIlSrVq3EWJ4/tVGtWjUMGzYM165dw9atW1GhQgVMnjwZw4YNk7ZRtw91cLIlERHRS/zwww8wMTHBX3/9hUWLFmH+/Pn47rvvADy9P8b06dPx999/Y+fOnbhx4wYGDhwobTtlyhRcvnwZe/bswb///osVK1bAwcFB2rZjx46wsbFBTEwMjh8/Dmtra3Tq1An5+fkqxzdv3jw0atQI586dw8iRIzFixAjEx8drtY/ScESCiIjoJdzd3bFgwQLIZDLUqlULFy5cwIIFCzB06FAMHjxYqle9enUsXrwYjRs3RnZ2NqytrXHr1i00aNAAjRo1AgClUYDNmzdDoVDgu+++k+6vsW7dOtjb2+Pw4cPo0KGDSvF17twZI0eOBACEhYVhwYIFOHToEGrVqqW1PkrDEQkiIqKXaNq0qdKNtPz8/HD16lUUFhYiLi4OXbt2hYeHB2xsbNCyZUsAwK1btwAAI0aMwKZNm1C/fn1MnDgRJ06ckNr5+++/ce3aNdjY2MDa2hrW1taoWLEicnNzcf36dZXjq1u3rvR/mUwGFxcXpKamarWP0nBEgoiISEO5ubno2LEjOnbsiI0bN8LR0RG3bt1Cx44dpdMGgYGBuHnzJn7//Xfs378fbdu2RWhoKObOnYvs7Gw0bNgQGzduLNa2o6OjynE8fxWHTCaDQqEAAK31URomEkRERC9x6tQppeXY2FjUrFkTV65cwf379zF79my4u7sDAM6cOVNse0dHR4SEhCAkJAT+/v74/PPPMXfuXLzzzjvYvHkznJycyuwJ1WXdB09tEBERvcStW7cwfvx4xMfH46effsKSJUswduxYeHh4QC6XY8mSJfjvv//wyy+/YPr06UrbfvXVV9i1axeuXbuGS5cuYffu3ahduzYAoF+/fnBwcED37t0RExODhIQEHD58GGPGjMHt27e1EntZ98ERCSIi0qnyuAy+sLAQ586dQ4MGDTTafsCAAXj8+DHeffddGBsbY+zYsRg2bBhkMhnWr1+PL774AosXL8Y777yDuXPnolu3btK2crkc4eHhuHHjBiwsLODv749NmzYBACwtLXH06FGEhYWhZ8+eyMrKQuXKldG2bVutjR6UdR8yIYTQQpx6KzMzE3Z2dsjIyCh1h/E+EurjPtMM9xu9yXJzc5GQkABPT88XPpa6LDybSPAx4qpR9f3iqQ0iIiLSGBMJIiIi0hgTCSIiItIYEwkiIiLSGK/aICK9xgmqRPqNIxJERESkMbVGJNLT07Fjxw7ExMTg5s2bePToERwdHdGgQQN07NgRzZo1K6s4iYiISA+pNCKRlJSEjz/+GK6urpgxYwYeP36M+vXro23btqhSpQoOHTqE9u3bw8fHB5s3b9YokNmzZ0Mmk2HcuHFSWW5uLkJDQ1GpUiVYW1sjKCgId+/e1ah9IiIi0j6VRiQaNGiAkJAQxMXFwcfHp8Q6jx8/xs6dO7Fw4UIkJiZiwoQJKgdx+vRprFq1SunpZQDw6aef4rfffsPWrVthZ2eHUaNGoWfPnjh+/LjKbRMREVHZUSmRuHz5MipVqvTCOhYWFggODkZwcDDu37+vcgDZ2dno168f1qxZgxkzZkjlGRkZ+P777xEdHY02bdoAePr89Nq1ayM2NhZNmzZVuQ8iItJfc5K7lk9HlYGDqUCY66/l098bQqVEoiiJKCgowPDhwzFlyhR4enq+tL4qQkND0aVLF7Rr104pkYiLi0NBQQHatWsnlXl7e8PDwwMnT54sNZHIy8tDXl6etJyZmalyLERErwte7VI+kpOT8dlnn+HMmTO4du0axowZg4ULF+o6rHKl1lUbpqam2LZtm9Y637RpE86ePYtZs2YVW5eSkgK5XA57e3ulcmdnZ6SkpJTa5qxZs2BnZyf9FD3WlYiISNvy8vLg6OiIyZMno169eroORyfUvvyzR48e2Llz5yt3nJiYiLFjx2Ljxo1afXhLeHg4MjIypJ/ExESttU1ERG+W1atXw83NDQqFQqm8e/fuGDx4MKpVq4ZFixZhwIABsLOz01GUuqX2Dalq1qyJadOm4fjx42jYsCGsrKyU1o8ZM0alduLi4pCamop33nlHKissLMTRo0exdOlS7Nu3D/n5+UhPT1calbh79y5cXFxKbdfMzAxmZmbqvSgiIqISfPjhhxg9ejQOHTqEtm3bAgAePHiAvXv34vfff9dxdPpB7UTi+++/h729PeLi4hAXF6e0TiaTqZxItG3bFhcuXFAqGzRoELy9vREWFgZ3d3eYmpri4MGDCAoKAgDEx8fj1q1b8PPzUzdsIiIitVWoUAGBgYGIjo6WEomff/4ZDg4OaN26tY6j0w9qJxIJCQla6djGxgZvv/22UpmVlRUqVaoklQ8ZMgTjx49HxYoVYWtri9GjR8PPz49XbBARUbnp168fhg4diuXLl8PMzAwbN25Enz59YGTEm0MDr/isDSEEgKcjEWVhwYIFMDIyQlBQEPLy8tCxY0csX768TPoiIiIqSdeuXSGEwG+//YbGjRsjJiYGCxYs0HVYekOjdGrDhg3w9fWFhYUFLCwsULduXURFRb1yMIcPH1a6bMbc3BzLli3DgwcPkJOTg+3bt79wfgQREZG2mZubo2fPnti4cSN++ukn1KpVS2l+35tO7RGJ+fPnY8qUKRg1ahSaN28O4On1yp988gnu3buHTz/9VOtBEhER6VK/fv3w3nvv4dKlS+jfv7/SuvPnzwN4eoPFtLQ0nD9/HnK5vNQ7Qb9u1E4klixZghUrVmDAgAFSWbdu3VCnTh1EREQwkSAiIrWUx50mCwsLce7cOTRo0ECj7du0aYOKFSsiPj4effv2VVr3bJtxcXGIjo5G1apVcePGjVcJ2WConUgkJyeX+JTPZs2aITk5WStBERER6RMjIyMkJSWVuK5ovuCbSu05EjVq1MCWLVuKlW/evBk1a9bUSlBERERkGNQekYiMjETv3r1x9OhRaY7E8ePHcfDgwRITDCIiInp9qT0iERQUhL/++gsODg7YuXMndu7cCQcHB/z11194//33yyJGIiIi0lNqjUg8+/TPH3/8saxiIiIiIgOh06d/EhHRm+dNn5xoKFR9n3T29E8iInqzGBsbAwDy8/N1HAmp4tGjRwCeDiK8iM6e/klERG8WExMTWFpaIi0tDaampuX6rIrCwkIAQG5urpTQUMmEEHj06BFSU1Nhb2//0v2ls6d/EhHRm0Umk8HV1RUJCQm4efNmufatUChw79493Lhxgw/bUpG9vb1Kj6VQK5EQQuDw4cNwcnKChYWFxsHR/8xJ7qqzvsvjbnJERM+Sy+WoWbNmuZ/eyM7ORpcuXXDmzBlYW1uXa9+GyNTUVOWRG7UTiZo1a+LSpUu8+RTpDJMvIsNmZGQEc3Pzcu0zPz8fN2/ehFwuL/e+X3dqje8YGRmhZs2auH//flnFQ0RERAZE7RNFs2fPxueff46LFy+WRTxERERkQNSebDlgwAA8evQI9erVg1wuLzZX4sGDB1oLjoiIiPSb2onEwoULyyAMIiIiMkRqJxIhISFlEQcREREZII0upr1+/TomT56M4OBgpKamAgD27NmDS5cuaTU4IiIi0m9qJxJHjhyBr68vTp06he3btyM7OxsA8Pfff2Pq1KlaD5CIiIj0l9qJxKRJkzBjxgzs378fcrlcKm/Tpg1iY2O1GhwRERHpN7UTiQsXLuD9998vVu7k5IR79+5pJSgiIiIyDGonEvb29khOTi5Wfu7cOVSuXFkrQREREZFhUPuqjT59+iAsLAxbt26FTCaDQqHA8ePHMWHCBAwYMKAsYiQiLeCtxYmoLKg9IjFz5kx4e3vD3d0d2dnZ8PHxQUBAAJo1a4bJkyeXRYxERESkp9QekZDL5VizZg2++uorXLhwAdnZ2WjQoAEf4kVERPQGUjuRKOLu7g53d3dtxkJEREQGRqMbUhEREREBTCSIiIjoFTCRICIiIo0xkSAiIiKNaZRIxMTEoH///vDz88OdO3cAAFFRUTh27JhWgyMiIiL9pnYisW3bNnTs2BEWFhY4d+4c8vLyAAAZGRmYOXOm1gMkIiIi/aV2IjFjxgysXLkSa9asgampqVTevHlznD17VqvBERERkX5TO5GIj49HQEBAsXI7Ozukp6drIyYiIiIyEGonEi4uLrh27Vqx8mPHjqF69epaCYqIiIgMg9qJxNChQzF27FicOnUKMpkMSUlJ2LhxIyZMmIARI0aURYxERESkp9S+RfakSZOgUCjQtm1bPHr0CAEBATAzM8OECRMwevTosoiRiIiI9JRaIxKFhYWIiYlBaGgoHjx4gIsXLyI2NhZpaWmYPn262p2vWLECdevWha2tLWxtbeHn54c9e/ZI63NzcxEaGopKlSrB2toaQUFBuHv3rtr9EBERUdlQK5EwNjZGhw4d8PDhQ8jlcvj4+ODdd9+FtbW1Rp1XqVIFs2fPRlxcHM6cOYM2bdqge/fuuHTpEgDg008/xa+//oqtW7fiyJEjSEpKQs+ePTXqi4iIiLRP7VMbb7/9Nv777z94enq+cuddu3ZVWv7666+xYsUKxMbGokqVKvj+++8RHR2NNm3aAADWrVuH2rVrIzY2Fk2bNn3l/omIiOjVaHQfiQkTJmD37t1ITk5GZmam0o+mCgsLsWnTJuTk5MDPzw9xcXEoKChAu3btpDre3t7w8PDAyZMnS20nLy9PazERERHRi6k9ItG5c2cAQLdu3SCTyaRyIQRkMhkKCwvVau/ChQvw8/NDbm4urK2tsWPHDvj4+OD8+fOQy+Wwt7dXqu/s7IyUlJRS25s1axYiIyPVioGIqCRzkru+vFIZCXP9VWd9E6lD7UTi0KFDWg2gVq1aOH/+PDIyMvDzzz8jJCQER44c0bi98PBwjB8/XlrOzMyEu7u7NkIlIiKi56idSLRs2VKrAcjlctSoUQMA0LBhQ5w+fRqLFi1C7969kZ+fj/T0dKVRibt378LFxaXU9szMzGBmZqbVGImIiKhkaicSRR49eoRbt24hPz9fqbxu3bqvFJBCoUBeXh4aNmwIU1NTHDx4EEFBQQCe3p771q1b8PPze6U+iIiISDvUTiTS0tIwaNAgpfs9PEudORLh4eEIDAyEh4cHsrKyEB0djcOHD2Pfvn2ws7PDkCFDMH78eFSsWBG2trYYPXo0/Pz8eMUGEZEe49ySN4vaicS4ceOQnp6OU6dOoVWrVtixYwfu3r2LGTNmYN68eWq1lZqaigEDBiA5ORl2dnaoW7cu9u3bh/bt2wMAFixYACMjIwQFBSEvLw8dO3bE8uXL1Q2ZiIiIyojaicSff/6JXbt2oVGjRjAyMkLVqlXRvn172NraYtasWejSpYvKbX3//fcvXG9ubo5ly5Zh2bJl6oZJRERE5UDt+0jk5OTAyckJAFChQgWkpaUBAHx9fXH27FntRkdERER6Te1EolatWoiPjwcA1KtXD6tWrcKdO3ewcuVKuLq6aj1AIiIi0l9qn9oYO3YskpOTAQBTp05Fp06dsHHjRsjlcqxfv17b8REREZEeUzuR6N+/v/T/hg0b4ubNm7hy5Qo8PDzg4OCg1eCIiIhIv2l8H4kilpaWeOedd7QRCxERERkYtROJwYMHv3D92rVrNQ6GiIiIDIvaicTDhw+VlgsKCnDx4kWkp6dLj/smIiKiN4PaicSOHTuKlSkUCowYMQJeXl5aCYqIiIgMg9qXf5bYiJERxo8fjwULFmijOSIiIjIQWkkkAOD69et48uSJtpojIiIiA6D2qY3x48crLQshkJycjN9++w0hISFaC4yIiIj0n9qJxLlz55SWjYyM4OjoiHnz5r30ig4iIiJ6vaidSBw6dKgs4iAiIiIDpLU5EkRERPTmUXtEokGDBpDJZCrV5dNAiYiIXm9qJxKdOnXC8uXL4ePjAz8/PwBAbGwsLl26hBEjRsDCwkLrQRIREZF+UjuRSEtLw5gxYzB9+nSl8qlTpyIxMZG3yCYiInqDqD1HYuvWrRgwYECx8v79+2Pbtm1aCYqIiIgMg9qJhIWFBY4fP16s/Pjx4zA3N9dKUERERGQY1D61MW7cOIwYMQJnz57Fu+++CwA4deoU1q5diylTpmg9QCIiItJfaicSkyZNQvXq1bFo0SL8+OOPAIDatWtj3bp16NWrl9YDJCIiIv2ldiIBAL169WLSQEREROrPkUhMTMTt27el5b/++gvjxo3D6tWrtRoYERER6T+1E4m+fftKt8lOSUlBu3bt8Ndff+HLL7/EtGnTtB4gERER6S+1E4mLFy9Kkyy3bNkCX19fnDhxAhs3bsT69eu1HR8RERHpMbUTiYKCApiZmQEADhw4gG7dugEAvL29kZycrN3oiIiISK+pnUjUqVMHK1euRExMDPbv349OnToBAJKSklCpUiWtB0hERET6S+1EYs6cOVi1ahVatWqF4OBg1KtXDwDwyy+/SKc8iIiI6M2g9uWfrVq1wr1795CZmYkKFSpI5cOGDYOlpaVWgyMiIiL9ptF9JIyNjZWSCACoVq2aNuIhIiIiA6L2qQ0iIiKiIkwkiIiISGNMJIiIiEhjTCSIiIhIYxpNtjx48CAOHjyI1NRUKBQKpXVr167VSmBERESk/9ROJCIjIzFt2jQ0atQIrq6ukMlkZREXERERGQC1E4mVK1di/fr1+Oijj8oiHiIiIjIgas+RyM/PR7NmzbTS+axZs9C4cWPY2NjAyckJPXr0QHx8vFKd3NxchIaGolKlSrC2tkZQUBDu3r2rlf6JiIjo1aidSHz88ceIjo7WSudHjhxBaGgoYmNjsX//fhQUFKBDhw7IycmR6nz66af49ddfsXXrVhw5cgRJSUno2bOnVvonIiKiV6P2qY3c3FysXr0aBw4cQN26dWFqaqq0fv78+Sq3tXfvXqXl9evXw8nJCXFxcQgICEBGRga+//57REdHo02bNgCAdevWoXbt2oiNjUXTpk3VDZ+IiIi0SO1E4p9//kH9+vUBABcvXlRa96oTLzMyMgAAFStWBADExcWhoKAA7dq1k+p4e3vDw8MDJ0+eLDGRyMvLQ15enrScmZn5SjERERFR6dROJA4dOlQWcUChUGDcuHFo3rw53n77bQBASkoK5HI57O3tleo6OzsjJSWlxHZmzZqFyMjIMomRiIiIlOnNDalCQ0Nx8eJFbNq06ZXaCQ8PR0ZGhvSTmJiopQiJiIjoeRrdkOrMmTPYsmULbt26hfz8fKV127dvV7u9UaNGYffu3Th69CiqVKkilbu4uCA/Px/p6elKoxJ3796Fi4tLiW2ZmZnBzMxM7RiIiIhIfWqPSGzatAnNmjXDv//+ix07dqCgoACXLl3Cn3/+CTs7O7XaEkJg1KhR2LFjB/788094enoqrW/YsCFMTU1x8OBBqSw+Ph63bt2Cn5+fuqETERGRlqk9IjFz5kwsWLAAoaGhsLGxwaJFi+Dp6Ynhw4fD1dVVrbZCQ0MRHR2NXbt2wcbGRpr3YGdnBwsLC9jZ2WHIkCEYP348KlasCFtbW4wePRp+fn68YoOIiEgPqD0icf36dXTp0gUAIJfLkZOTA5lMhk8//RSrV69Wq60VK1YgIyMDrVq1gqurq/SzefNmqc6CBQvw3nvvISgoCAEBAXBxcdHo9AkRERFpn9ojEhUqVEBWVhYAoHLlyrh48SJ8fX2Rnp6OR48eqdWWEOKldczNzbFs2TIsW7ZM3VCJiIiojKmdSAQEBGD//v3w9fXFhx9+iLFjx+LPP//E/v370bZt27KIkYiIiPSU2onE0qVLkZubCwD48ssvYWpqihMnTiAoKAiTJ0/WeoBERESkv9ROJIruOgkARkZGmDRpklYDIiIiIsOh0Q2prl+/jsmTJyM4OBipqakAgD179uDSpUtaDY6IiIj0m9qJxJEjR+Dr64tTp05h+/btyM7OBgD8/fffmDp1qtYDJCIiIv2ldiIxadIkzJgxA/v374dcLpfK27Rpg9jYWK0GR0RERPpN7UTiwoULeP/994uVOzk54d69e1oJioiIiAyD2omEvb09kpOTi5WfO3cOlStX1kpQREREZBjUTiT69OmDsLAwpKSkQCaTQaFQ4Pjx45gwYQIGDBhQFjESERGRnlI7kZg5cya8vb3h7u6O7Oxs+Pj4ICAgAM2aNeN9JIiIiN4wat9HQi6XY82aNZgyZQouXryI7OxsNGjQADVr1iyL+IiIiEiPqZ1IFPHw8ICHh4c2YyEiIiIDo3YiIYTAzz//jEOHDiE1NRUKhUJpPZ/MSURE9OZQO5EYN24cVq1ahdatW8PZ2Rkymaws4iIiIiIDoHYiERUVhe3bt6Nz585lEQ8REREZELWv2rCzs0P16tXLIhYiIiIyMGonEhEREYiMjMTjx4/LIh4iIiIyIGqf2ujVqxd++uknODk5oVq1ajA1NVVaf/bsWa0FR0RERPpN7UQiJCQEcXFx6N+/PydbEhERveHUTiR+++037Nu3Dy1atCiLeIiIiMiAqD1Hwt3dHba2tmURCxERERkYtROJefPmYeLEibhx40YZhENERESGRO1TG/3798ejR4/g5eUFS0vLYpMtHzx4oLXgiIiISL+pnUgsXLiwDMIgIiIiQ6TRVRtEREREgIpzJHJyctRqVN36REREZJhUSiRq1KiB2bNnIzk5udQ6Qgjs378fgYGBWLx4sdYCJCIiIv2l0qmNw4cP44svvkBERATq1auHRo0awc3NDebm5nj48CEuX76MkydPwsTEBOHh4Rg+fHhZx01ERER6QKVEolatWti2bRtu3bqFrVu3IiYmBidOnMDjx4/h4OCABg0aYM2aNQgMDISxsXFZx0xERER6Qq3Jlh4eHvjss8/w2WeflVU8REREZEDUviEVERERUREmEkRERKQxJhJERESkMSYSREREpDEmEkRERKQxlROJb775Bo8fP5aWjx8/jry8PGk5KysLI0eO1G50REREpNdUTiTCw8ORlZUlLQcGBuLOnTvS8qNHj7Bq1SrtRkdERER6TeVEQgjxwmUiIiJ68+h0jsTRo0fRtWtXuLm5QSaTYefOnUrrhRD46quv4OrqCgsLC7Rr1w5Xr17VTbBERERUjE4TiZycHNSrVw/Lli0rcf0333yDxYsXY+XKlTh16hSsrKzQsWNH5ObmlnOkREREVBK1bpH93XffwdraGgDw5MkTrF+/Hg4ODgCgNH9CVYGBgQgMDCxxnRACCxcuxOTJk9G9e3cAwIYNG+Ds7IydO3eiT58+avdHRERE2qVyIuHh4YE1a9ZIyy4uLoiKiipWR1sSEhKQkpKCdu3aSWV2dnZo0qQJTp48WWoikZeXp3Q1SWZmptZiIiIiImUqJxI3btwowzCKS0lJAQA4OzsrlTs7O0vrSjJr1ixERkaWaWxERET01Gt3Q6rw8HBkZGRIP4mJiboOiYiI6LWlciJx8uRJ7N69W6lsw4YN8PT0hJOTE4YNG6Z0SuFVubi4AADu3r2rVH737l1pXUnMzMxga2ur9ENERERlQ+VEYtq0abh06ZK0fOHCBQwZMgTt2rXDpEmT8Ouvv2LWrFlaC8zT0xMuLi44ePCgVJaZmYlTp07Bz89Pa/0QERGR5lSeI3H+/HlMnz5dWt60aROaNGkiTcB0d3fH1KlTERERoXLn2dnZuHbtmrSckJCA8+fPo2LFivDw8MC4ceMwY8YM1KxZE56enpgyZQrc3NzQo0cPlfsgIiKisqNyIvHw4UOliY9HjhxRunSzcePGas9HOHPmDFq3bi0tjx8/HgAQEhKC9evXY+LEicjJycGwYcOQnp6OFi1aYO/evTA3N1erHyIiIiobKicSzs7OSEhIgLu7O/Lz83H27FmlqyOysrJgamqqVuetWrV64a22ZTIZpk2bhmnTpqnVLhEREZUPledIdO7cGZMmTUJMTAzCw8NhaWkJf39/af0///wDLy+vMgmSiIiI9JPKIxLTp09Hz5490bJlS1hbW+OHH36AXC6X1q9duxYdOnQokyCJiIhIP6mcSDg4OODo0aPIyMiAtbU1jI2NldZv3bpVun02ERERvRnUetYG8PQ21SWpWLHiKwdDREREhkXlRGLw4MEq1Vu7dq3GwRAREZFhUTmRWL9+PapWrYoGDRq88EoLIiIienOonEiMGDECP/30ExISEjBo0CD079+fpzOIiIjecCpf/rls2TIkJydj4sSJ+PXXX+Hu7o5evXph3759HKEgIiJ6Q6n19E8zMzMEBwdj//79uHz5MurUqYORI0eiWrVqyM7OLqsYiYiISE9p/BhxIyMjyGQyCCFQWFiozZiIiIjIQKiVSOTl5eGnn35C+/bt8dZbb+HChQtYunQpbt26xXtIEBERvYFUnmw5cuRIbNq0Ce7u7hg8eDB++uknODg4lGVsREREpOdUTiRWrlwJDw8PVK9eHUeOHMGRI0dKrLd9+3atBUdERET6TeVEYsCAAZDJZGUZCxERERkYtW5IRURERPQsja/aICIiImIiQURERBpjIkFEREQaYyJBREREGmMiQURERBpjIkFEREQaYyJBREREGmMiQURERBpjIkFEREQaYyJBREREGmMiQURERBpjIkFEREQaYyJBREREGmMiQURERBpjIkFEREQaYyJBREREGmMiQURERBpjIkFEREQaYyJBREREGmMiQURERBpjIkFEREQaYyJBREREGmMiQURERBoziERi2bJlqFatGszNzdGkSRP89ddfug6JiIiIYACJxObNmzF+/HhMnToVZ8+eRb169dCxY0ekpqbqOjQiIqI3nt4nEvPnz8fQoUMxaNAg+Pj4YOXKlbC0tMTatWt1HRoREdEbz0TXAbxIfn4+4uLiEB4eLpUZGRmhXbt2OHnyZInb5OXlIS8vT1rOyMgAAGRmZpbaT05OjpYiVl9uVoHO+s60Kn2fvAz3mWa439THfaYZ7rfnyv//7wAhRHmG80bQ60Ti3r17KCwshLOzs1K5s7Mzrly5UuI2s2bNQmRkZLFyd3f3MonRkEXATtchGBzuM81wv6mP+0wzL9tvWVlZsLPjvtUmvU4kNBEeHo7x48dLywqFAg8ePEClSpUgk8l0GFlxmZmZcHd3R2JiImxtbXUdjkHgPtMM95v6uM80o6/7TQiBrKwsuLm56TqU145eJxIODg4wNjbG3bt3lcrv3r0LFxeXErcxMzODmZmZUpm9vX1ZhagVtra2evWBMwTcZ5rhflMf95lm9HG/cSSibOj1ZEu5XI6GDRvi4MGDUplCocDBgwfh5+enw8iIiIgI0PMRCQAYP348QkJC0KhRI7z77rtYuHAhcnJyMGjQIF2HRkRE9MbT+0Sid+/eSEtLw1dffYWUlBTUr18fe/fuLTYB0xCZmZlh6tSpxU7FUOm4zzTD/aY+7jPNcL+9eWSC18IQERGRhvR6jgQRERHpNyYSREREpDEmEkRERKQxJhJERESkMSYSOnD06FF07doVbm5ukMlk2Llzp65D0nuzZs1C48aNYWNjAycnJ/To0QPx8fG6DkuvrVixAnXr1pVuDOTn54c9e/boOiyDMnv2bMhkMowbN07Xoei1iIgIyGQypR9vb29dh0XlhImEDuTk5KBevXpYtmyZrkMxGEeOHEFoaChiY2Oxf/9+FBQUoEOHDjp9MJG+q1KlCmbPno24uDicOXMGbdq0Qffu3XHp0iVdh2YQTp8+jVWrVqFu3bq6DsUg1KlTB8nJydLPsWPHdB0SlRO9v4/E6ygwMBCBgYG6DsOg7N27V2l5/fr1cHJyQlxcHAICAnQUlX7r2rWr0vLXX3+NFStWIDY2FnXq1NFRVIYhOzsb/fr1w5o1azBjxgxdh2MQTExMSn10Ab3eOCJBBqno8fAVK1bUcSSGobCwEJs2bUJOTg5vL6+C0NBQdOnSBe3atdN1KAbj6tWrcHNzQ/Xq1dGvXz/cunVL1yFROeGIBBkchUKBcePGoXnz5nj77bd1HY5eu3DhAvz8/JCbmwtra2vs2LEDPj4+ug5Lr23atAlnz57F6dOndR2KwWjSpAnWr1+PWrVqITk5GZGRkfD398fFixdhY2Oj6/CojDGRIIMTGhqKixcv8hysCmrVqoXz588jIyMDP//8M0JCQnDkyBEmE6VITEzE2LFjsX//fpibm+s6HIPx7KnaunXrokmTJqhatSq2bNmCIUOG6DAyKg9MJMigjBo1Crt378bRo0dRpUoVXYej9+RyOWrUqAEAaNiwIU6fPo1FixZh1apVOo5MP8XFxSE1NRXvvPOOVFZYWIijR49i6dKlyMvLg7GxsQ4jNAz29vZ46623cO3aNV2HQuWAiQQZBCEERo8ejR07duDw4cPw9PTUdUgGSaFQIC8vT9dh6K22bdviwoULSmWDBg2Ct7c3wsLCmESoKDs7G9evX8dHH32k61CoHDCR0IHs7GylTD0hIQHnz59HxYoV4eHhocPI9FdoaCiio6Oxa9cu2NjYICUlBQBgZ2cHCwsLHUenn8LDwxEYGAgPDw9kZWUhOjoahw8fxr59+3Qdmt6ysbEpNu/GysoKlSpV4nycF5gwYQK6du2KqlWrIikpCVOnToWxsTGCg4N1HRqVAyYSOnDmzBm0bt1aWh4/fjwAICQkBOvXr9dRVPptxYoVAIBWrVopla9btw4DBw4s/4AMQGpqKgYMGIDk5GTY2dmhbt262LdvH9q3b6/r0Og1c/v2bQQHB+P+/ftwdHREixYtEBsbC0dHR12HRuWAjxEnIiIijfE+EkRERKQxJhJERESkMSYSREREpDEmEkRERKQxJhJERESkMSYSREREpDEmEkRERKQxJhJERESkMSYSRG+AR48eISgoCLa2tpDJZEhPT3/pNuvXr4e9vf0L60RERKB+/fpaiZGIDBNvkU30Bvjhhx8QExODEydOwMHBAXZ2droOiYheE0wkiN4A169fR+3atfngKSLSOp7aINIjCoUC33zzDWrUqAEzMzN4eHjg66+/RrNmzRAWFqZUNy0tDaampjh69OgL22zVqhXmzZuHo0ePQiaTSQ8+e/jwIQYMGIAKFSrA0tISgYGBuHr16gvbmj17NpydnWFjY4MhQ4YgNzf3lV4vERk+JhJEeiQ8PByzZ8/GlClTcPnyZURHR8PZ2Rn9+vXDpk2b8Owz9jZv3gw3Nzf4+/u/sM3t27dj6NCh8PPzQ3JyMrZv3w4AGDhwIM6cOYNffvkFJ0+ehBACnTt3RkFBQYntbNmyBREREZg5cybOnDkDV1dXLF++XHsvnogMkyAivZCZmSnMzMzEmjVriq1LTU0VJiYm4ujRo1KZn5+fCAsLU6ntsWPHipYtW0rL//d//ycAiOPHj0tl9+7dExYWFmLLli1CCCHWrVsn7OzslPobOXKkUrtNmjQR9erVUykGIno9cUSCSE/8+++/yMvLQ9u2bYutc3R0RIcOHbBx40YAQEJCAk6ePIl+/fpp3JeJiQmaNGkilVWqVAm1atXCv//+W+o2z9YHAD8/P436J6LXBxMJIj1hYWHxwvX9+vXDzz//jIKCAkRHR8PX1xe+vr7lFB0RUcmYSBDpiZo1a8LCwgIHDx4scX337t2Rm5uLvXv3Ijo6WuPRCACoXbs2njx5glOnTkll9+/fR3x8PHx8fErd5tn6ABAbG6txDET0euDln0R6wtzcHGFhYZg4cSLkcjmaN2+OtLQ0XLp0CUOGDIGVlRV69OiBKVOm4N9//0VwcLDGfdWsWRPdu3fH0KFDsWrVKtjY2GDSpEmoXLkyunfvXuI2Y8eOxcCBA9GoUSM0b94cGzduxKVLl1C9enWN4yAiw8dEgkiPTJkyBSYmJvjqq6+QlJQEV1dXfPLJJ9L6fv36oXPnzggICICHh8cr9bVu3TqMHTsW7733HvLz8xEQEIDff/8dpqamJdbv3bs3rl+/jokTJyI3NxdBQUEYMWIE9u3b90pxEJFhkwnxzPVkRERERGrgHAkiIiLSGBMJIgMXExMDa2vrUn+IiMoST20QGbjHjx/jzp07pa6vUaNGOUZDRG8aJhJERESkMZ7aICIiIo0xkSAiIiKNMZEgIiIijTGRICIiIo0xkSAiIiKNMZEgIiIijTGRICIiIo39P2Tu7JxSrIZTAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 500x300 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "palette = {\n",
    "    'baseline': '#C3C3C3',\n",
    "    'v1':'#83FF33'\n",
    "}\n",
    "\n",
    "plt.rcParams[\"figure.figsize\"] = (5,3)\n",
    "sns.barplot(res_df, x='cv_fold', y='mse', hue=\"variable\", palette=palette)\n",
    "plt.title(\"Model Performance Improvement (lower is better)\")\n",
    "plt.ylabel(\"MSE (mean square error)\")\n",
    "plt.legend(bbox_to_anchor= (0.9, 1.0))\n",
    "plt.savefig(\"figure/model_performance.png\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1683d319-6261-4c01-af75-dfe76116daaa",
   "metadata": {},
   "source": [
    "**Observation**\n",
    "1. v1 model is consistently better than baseline\n",
    "2. Downward trend as data increase"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eec917ae-59ba-416d-8737-86e47135956c",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}