predict-energy-behavior-of-prosumers-predictions.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "_cell_guid": "b1076dfc-b9ad-4769-8c92-a6c4dae69d19",
    "_uuid": "8f2839f25d086af736a60e9eeb907d3b93b6e0e5"
   },
   "source": [
    "# Predict Energy Behavior of Prosumers - Predictions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Import libraries"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# data processing\n",
    "import numpy as np \n",
    "import pandas as pd\n",
    "from geopy import distance\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "\n",
    "import os\n",
    "\n",
    "# modeling\n",
    "from sklearn.model_selection import KFold\n",
    "from sklearn.ensemble import RandomForestRegressor\n",
    "from sklearn import svm\n",
    "import tensorflow as tf\n",
    "from tensorflow import keras\n",
    "from tensorflow.keras import layers\n",
    "\n",
    "# evaluation\n",
    "from sklearn.metrics import mean_squared_error, mean_absolute_error\n",
    "\n",
    "# visualization\n",
    "import plotly.express as px "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Import data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_dir = os.path.join(os.getcwd(), \"data\\\\\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "prosumer = pd.read_csv(os.path.join(data_dir, 'prosumer.csv'))\n",
    "weather_history = pd.read_csv(os.path.join(data_dir, 'weather_history.csv'))\n",
    "weather_forecast = pd.read_csv(os.path.join(data_dir, 'weather_forecast.csv'))\n",
    "gas_prices = pd.read_csv(os.path.join(data_dir, 'gas_prices.csv'))\n",
    "electricity_prices = pd.read_csv(os.path.join(data_dir, 'electricity_prices.csv'))\n",
    "client = pd.read_csv(os.path.join(data_dir, 'client.csv'))\n",
    "counties = pd.read_json(os.path.join(data_dir, 'county_id_to_name_map.json'), typ=\"series\")\n",
    "weather_stations = pd.read_csv(os.path.join(data_dir, 'weather_station_to_county_mapping.csv'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# how big subset of the data to include in dataset\n",
    "data_sample_size = 10000"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-12-08T11:48:43.630407Z",
     "iopub.status.busy": "2023-12-08T11:48:43.629884Z",
     "iopub.status.idle": "2023-12-08T11:48:43.637205Z",
     "shell.execute_reply": "2023-12-08T11:48:43.635463Z",
     "shell.execute_reply.started": "2023-12-08T11:48:43.630365Z"
    }
   },
   "source": [
    "## 3. Preprocess"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-12-11T08:10:35.278618Z",
     "iopub.status.busy": "2023-12-11T08:10:35.278264Z",
     "iopub.status.idle": "2023-12-11T08:10:35.284941Z",
     "shell.execute_reply": "2023-12-11T08:10:35.283414Z",
     "shell.execute_reply.started": "2023-12-11T08:10:35.278588Z"
    }
   },
   "source": [
    "### Add weather station data into weather data\n",
    "\n",
    "Add information about the closest weather station into the weather data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "weather_stations.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Number of weather stations: \", len(weather_stations))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Number stations with county info: \", len(weather_stations[~weather_stations[\"county\"].isna()]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-12-13T09:11:47.360601Z",
     "iopub.status.busy": "2023-12-13T09:11:47.360264Z",
     "iopub.status.idle": "2023-12-13T09:11:47.365502Z",
     "shell.execute_reply": "2023-12-13T09:11:47.364397Z",
     "shell.execute_reply.started": "2023-12-13T09:11:47.360574Z"
    }
   },
   "source": [
    "**We need county info to connect weather data into prosumer data. Since, let's use only stations with county info available. \n",
    "For each weather data row, let's find the closest weather station with county info available.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def find_closest_weather_station(point: pd.Series, wss: pd.DataFrame):\n",
    "    \"\"\"\n",
    "    Find closest weather station with county information to a point.\n",
    "    \"\"\"\n",
    "    \n",
    "    point_coordinates = [point.latitude, point.longitude]\n",
    "    weather_station_coordinates = [[ws.latitude, ws.longitude] for ws in wss.itertuples()]\n",
    "    \n",
    "    # calculate distances to every weather station\n",
    "    dists = [distance.distance(point_coordinates, station).km\n",
    "             for station in weather_station_coordinates]\n",
    "    # get closest station\n",
    "    closest_dist = np.min(dists)\n",
    "    closest_station = wss.iloc[np.argmin(dists), :]\n",
    "    return closest_station, closest_dist\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def process_weather_data(wd: pd.DataFrame, wss: pd.DataFrame, feature_names: list):\n",
    "    \"\"\"\n",
    "    Process weather data by\n",
    "    - Adding weather station info to weather data by finding the closest weather station to a point\n",
    "    \"\"\"\n",
    "    \n",
    "    county_names = []\n",
    "    county_ids = []\n",
    "    for row in wd.itertuples():\n",
    "        closest_station, dist = find_closest_weather_station(row, wss)\n",
    "        if dist < 30:  # if station is less than 30 km away from point\n",
    "            county_names.append(closest_station.county_name)\n",
    "            county_ids.append(closest_station.county)\n",
    "        else:\n",
    "            county_names.append(None)\n",
    "            county_ids.append(None)\n",
    "    \n",
    "    \n",
    "    # add info to weather data\n",
    "    wd_processed = wd.copy()\n",
    "    wd_processed[\"county_name\"] = county_names\n",
    "    wd_processed[\"county\"] = county_ids\n",
    "\n",
    "\n",
    "    # take only data points from weather stations where county info is available\n",
    "    wd_processed = wd_processed.dropna(subset=[\"county\"], axis=0)\n",
    "\n",
    "    # format timestamps\n",
    "    wd_processed[\"forecast_datetime\"] = wd_processed.forecast_datetime.apply(lambda x: pd.to_datetime(x).tz_localize(None).strftime(\"%Y-%m-%d %H:%M:%S\"))\n",
    "\n",
    "    # aggregate data per county, per timestamp\n",
    "    wd_processed = wd_processed.groupby([\"county\", \"forecast_datetime\"])[feature_names].mean()\n",
    "\n",
    "    return wd_processed"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2023-12-13T09:37:01.342686Z",
     "iopub.status.busy": "2023-12-13T09:37:01.342337Z",
     "iopub.status.idle": "2023-12-13T09:37:01.346109Z",
     "shell.execute_reply": "2023-12-13T09:37:01.345423Z",
     "shell.execute_reply.started": "2023-12-13T09:37:01.342644Z"
    }
   },
   "source": [
    "### Add capacity and consumption info for prosumers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# capacity and consumption info is found in client.csv-file\n",
    "client.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def process_prosumer_data(prosumers: pd.DataFrame, clients: pd.DataFrame):\n",
    "    \"\"\"\n",
    "    Process prosumer data by adding capacity and consumption info into it.\n",
    "    \"\"\"\n",
    "    prosumers_proc = prosumers.copy()\n",
    "    # add feature for date\n",
    "    prosumers_proc[\"date\"] = prosumers_proc.datetime.apply(lambda x: pd.to_datetime(x).strftime(\"%Y-%m-%d\"))\n",
    "\n",
    "    # calculate consumption and capacity per segment per date\n",
    "    cons = client.groupby([\"product_type\", \"county\", \"is_business\", \"date\"])[\"eic_count\"].sum()\n",
    "    cap = client.groupby([\"product_type\", \"county\", \"is_business\", \"date\"])[\"installed_capacity\"].sum()\n",
    "\n",
    "    # add to prosumer data\n",
    "    prosumers_proc = pd.merge(pd.merge(prosumers_proc, cons, on=[\"product_type\", \"county\", \"is_business\", \"date\"]), cap, on=[\"product_type\", \"county\", \"is_business\", \"date\"])\n",
    "\n",
    "    return prosumers_proc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def process_elec_price_data(data: pd.DataFrame):\n",
    "    \"\"\"\n",
    "    Process electricity price data.\n",
    "    \"\"\"\n",
    "\n",
    "    data_processed = data.rename(columns={\"euros_per_mwh\": \"electricity_price\"})\n",
    "    return data_processed"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def process_gas_price_data(data: pd.DataFrame):\n",
    "    \"\"\"\n",
    "    Process gas price data.\n",
    "    \"\"\"\n",
    "    data_processed = data.copy()\n",
    "    data_processed[\"avg_gas_price\"] = data.apply(lambda row: np.mean([row.lowest_price_per_mwh, row.highest_price_per_mwh]), axis=1)\n",
    "    return data_processed\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Feature selection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Make dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def make_dataset(prosumer: pd.DataFrame, \n",
    "                 weather_forecast_data: pd.DataFrame, \n",
    "                 weather_stations: pd.DataFrame,\n",
    "                 client: pd.DataFrame,\n",
    "                 electricity_prices: pd.DataFrame,\n",
    "                 gas_prices: pd.DataFrame,\n",
    "                 weather_feature_names: list,\n",
    "                 scaler\n",
    "                 ):\n",
    "    \"\"\"\n",
    "    Make dataset for predicting prosumer consumption and production.\n",
    "    \"\"\"\n",
    "\n",
    "    # process data\n",
    "    weather_data = process_weather_data(weather_forecast_data, weather_stations, weather_feature_names)\n",
    "    prosumer_data = process_prosumer_data(prosumer, client)\n",
    "    electricity_data = process_elec_price_data(electricity_prices)\n",
    "    gas_data = process_gas_price_data(gas_prices)\n",
    "\n",
    "    # combine together\n",
    "    data = pd.merge(prosumer_data, \n",
    "                    weather_data, \n",
    "                    left_on=[\"county\", \"datetime\"], \n",
    "                    right_index=True,\n",
    "                    how=\"inner\")\n",
    "\n",
    "    # add price data\n",
    "    data = pd.merge(\n",
    "                pd.merge(data, \n",
    "                         electricity_data[[\"forecast_date\", \"electricity_price\"]], \n",
    "                         left_on=\"datetime\", \n",
    "                         right_on=\"forecast_date\", \n",
    "                         how=\"left\"), \n",
    "                gas_data[[\"forecast_date\", \"avg_gas_price\"]], \n",
    "                left_on=\"date\", \n",
    "                right_on=\"forecast_date\", \n",
    "                how=\"left\")\n",
    "\n",
    "    # select features for dataset\n",
    "    feats = [\"is_consumption\", \"eic_count\", \"installed_capacity\", \"electricity_price\", \"avg_gas_price\"] + weather_feature_names + [\"target\"]\n",
    "    data = data[feats].reset_index(drop=True)\n",
    "\n",
    "    if scaler is not None:\n",
    "        target = data.target\n",
    "        data = pd.DataFrame(scaler.fit_transform(data.drop(\"target\", axis=1)))\n",
    "        data[\"target\"] = target\n",
    "\n",
    "    return data    \n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# select which weather parameters to use\n",
    "weather_feat_names = ['direct_solar_radiation','surface_solar_radiation_downwards']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# select scaler \n",
    "scaler = StandardScaler()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "weather_forecast"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create dataset\n",
    "dataset = make_dataset(prosumer.head(data_sample_size), \n",
    "                       weather_forecast.iloc[:data_sample_size, :], \n",
    "                       weather_stations, \n",
    "                       client, \n",
    "                       electricity_prices, \n",
    "                       gas_prices,\n",
    "                       [\"direct_solar_radiation\", \"surface_solar_radiation_downwards\"],\n",
    "                       scaler\n",
    "                       )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Size of the dataset: \", dataset.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Modeling\n",
    "\n",
    "Let's define three different models to be fit into training data. The fitted models are later used to predict prosumers' energy consumption and production.\n",
    "\n",
    "Let's try to find the best model out of the three by applying the K-Fold method. In this method, each model is fit into subsets of the training data at the time, and their performance is evaluated\n",
    "by predicting values for other subset of the data.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def create_rf_model():\n",
    "    \"\"\"\n",
    "    Create a Random Forest model\n",
    "    \"\"\"\n",
    "    model = RandomForestRegressor()\n",
    "    return model\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def create_svm_model():\n",
    "    \"\"\"\n",
    "    Create a Support Vector Machine model\n",
    "    \"\"\"\n",
    "\n",
    "    model = svm.SVR()\n",
    "    return model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def create_nn_model():\n",
    "    \"\"\"\n",
    "    Create a neural network model\n",
    "    \"\"\"\n",
    "    model = keras.Sequential()\n",
    "    model.add(keras.layers.Dense(64))\n",
    "    model.add(keras.layers.Dense(32))\n",
    "    model.add(keras.layers.Dense(16))\n",
    "    model.add(keras.layers.Dense(8))\n",
    "    model.add(keras.layers.Dense(1))\n",
    "\n",
    "    model.compile(optimizer=\"adam\", loss=\"mae\")\n",
    "\n",
    "    return model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###  K-fold cross validation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# split data into training and test sets\n",
    "X_train, X_test, y_train, y_test = train_test_split(dataset.drop(\"target\", axis=1), dataset.target, test_size=0.2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cross_validate(data: pd.DataFrame, models, n_splits):\n",
    "    \"\"\"\n",
    "    Cross validate different models to find out the best fitting model for the data.\n",
    "    \"\"\"\n",
    "\n",
    "    kf = KFold(n_splits=10)\n",
    "    # for storing results\n",
    "    kfold_results = {}\n",
    "    \n",
    "    # define independent and dependent varibles\n",
    "    X = dataset.iloc[:, :-1]\n",
    "    y = dataset.target\n",
    "\n",
    "    # cross validate\n",
    "    run = 1\n",
    "    for train_idx, test_idx in kf.split(X):\n",
    "        scores = {}\n",
    "        \n",
    "        # define training and test sets\n",
    "        X_train, X_test, y_train, y_test = X.iloc[train_idx, :], X.iloc[test_idx, :], y.iloc[train_idx], y.iloc[test_idx]\n",
    " \n",
    "        # fit models and evaluate performaces\n",
    "        for model_name in list(models.keys()):\n",
    "            model = models.get(model_name)\n",
    "            if model_name == \"nn\":    \n",
    "                model.fit(X_train, y_train, verbose=0, epochs=20)\n",
    "            else:\n",
    "                model.fit(X_train, y_train)\n",
    "            y_hat = model.predict(X_test)\n",
    "            scores[model_name] = mean_absolute_error(y_test, y_hat)\n",
    "            \n",
    "        # store evaluation results\n",
    "        kfold_results[\"Validation round {}\".format(run)] = scores\n",
    "        run += 1\n",
    "\n",
    "    # find out the best model by taking mean score from all splits\n",
    "    mean_scores = {model: sum(kfold_results[round_name][model] \n",
    "                              for round_name in kfold_results) / len(kfold_results) \n",
    "                   for model in kfold_results[list(kfold_results.keys())[0]]\n",
    "                  }\n",
    "    best_model_idx = np.argmin(list(mean_scores.values()))\n",
    "    best_model_name = list(models.keys())[best_model_idx]\n",
    "    best_model = models[best_model_name]\n",
    "    best_model_score = mean_scores[best_model_name]\n",
    "\n",
    "    return best_model, best_model_score\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. Fit and evaluate models"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "best_model, score = cross_validate(pd.concat([X_train, y_train], axis=1), models={\"rf\": create_rf_model(), \"smv\": create_svm_model(), \"nn\": create_nn_model()}, n_splits=10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Best model: {}\\nMean absolute error: {}\".format(best_model, score))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Predict with test data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_hat = best_model.predict(X_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 7. Visualize predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sample_size = 100  # how many data points to visualize"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "results = pd.DataFrame(data=pd.concat([pd.Series(y_test.reset_index(drop=True)).iloc[:sample_size], pd.Series(y_hat).iloc[:sample_size]]), columns=[\"target\"])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "results[\"type\"] = pd.concat([pd.Series(np.repeat(\"Ground truth\", sample_size)), pd.Series(np.repeat(\"Prediction\", sample_size))]).values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = px.line(results, \n",
    "              y=\"target\", \n",
    "              color=\"type\", \n",
    "              title=\"Predictions and true values for energy amount\",\n",
    "              labels={\"target\": \"Energy amount\", \"index\": \"Timestep\"},\n",
    "              symbol=\"type\",\n",
    "              height=600\n",
    "              )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig.update_traces(selector=({\"name\": \"Prediction\"}), opacity=.6)\n",
    "fig.update_traces(selector=({\"name\": \"Ground truth\"}), opacity=.8)\n",
    "fig.show()"
   ]
  }
 ],
 "metadata": {
  "kaggle": {
   "accelerator": "none",
   "dataSources": [
    {
     "databundleVersionId": 7230081,
     "sourceId": 57236,
     "sourceType": "competition"
    }
   ],
   "dockerImageVersionId": 30615,
   "isGpuEnabled": false,
   "isInternetEnabled": true,
   "language": "python",
   "sourceType": "notebook"
  },
  "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.10.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}