MachineLearningNotebooks/how-to-use-azureml/automated-machine-learning/forecasting-energy-demand/auto-ml-forecasting-energy-demand.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Copyright (c) Microsoft Corporation. All rights reserved.\n",
        "\n",
        "Licensed under the MIT License."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "![Impressions](https://PixelServer20190423114238.azurewebsites.net/api/impressions/MachineLearningNotebooks/how-to-use-azureml/automated-machine-learning/forecasting-energy-demand/auto-ml-forecasting-energy-demand.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# Automated Machine Learning\n",
        "_**Energy Demand Forecasting**_\n",
        "\n",
        "## Contents\n",
        "1. [Introduction](#Introduction)\n",
        "1. [Setup](#Setup)\n",
        "1. [Data](#Data)\n",
        "1. [Train](#Train)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Introduction\n",
        "In this example, we show how AutoML can be used for energy demand forecasting.\n",
        "\n",
        "Make sure you have executed the [configuration](../../../configuration.ipynb) before running this notebook.\n",
        "\n",
        "In this notebook you would see\n",
        "1. Creating an Experiment in an existing Workspace\n",
        "2. Instantiating AutoMLConfig with new task type \"forecasting\" for timeseries data training, and other timeseries related settings: for this dataset we use the basic one: \"time_column_name\" \n",
        "3. Training the Model using local compute\n",
        "4. Exploring the results\n",
        "5. Viewing the engineered names for featurized data and featurization summary for all raw features\n",
        "6. Testing the fitted model"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Setup\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "import azureml.core\n",
        "import pandas as pd\n",
        "import numpy as np\n",
        "import logging\n",
        "\n",
        "from azureml.core.workspace import Workspace\n",
        "from azureml.core.experiment import Experiment\n",
        "from azureml.train.automl import AutoMLConfig\n",
        "from matplotlib import pyplot as plt\n",
        "from sklearn.metrics import mean_absolute_error, mean_squared_error, r2_score"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "As part of the setup you have already created a <b>Workspace</b>. For AutoML you would need to create an <b>Experiment</b>. An <b>Experiment</b> is a named object in a <b>Workspace</b>, which is used to run experiments."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "ws = Workspace.from_config()\n",
        "\n",
        "# choose a name for the run history container in the workspace\n",
        "experiment_name = 'automl-energydemandforecasting'\n",
        "# project folder\n",
        "project_folder = './sample_projects/automl-local-energydemandforecasting'\n",
        "\n",
        "experiment = Experiment(ws, experiment_name)\n",
        "\n",
        "output = {}\n",
        "output['SDK version'] = azureml.core.VERSION\n",
        "output['Subscription ID'] = ws.subscription_id\n",
        "output['Workspace'] = ws.name\n",
        "output['Resource Group'] = ws.resource_group\n",
        "output['Location'] = ws.location\n",
        "output['Project Directory'] = project_folder\n",
        "output['Run History Name'] = experiment_name\n",
        "pd.set_option('display.max_colwidth', -1)\n",
        "outputDf = pd.DataFrame(data = output, index = [''])\n",
        "outputDf.T"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Data\n",
        "Read energy demanding data from file, and preview data."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = pd.read_csv(\"nyc_energy.csv\", parse_dates=['timeStamp'])\n",
        "data.head()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# let's take note of what columns means what in the data\n",
        "time_column_name = 'timeStamp'\n",
        "target_column_name = 'demand'"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Split the data into train and test sets\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "X_train = data[data[time_column_name] < '2017-02-01']\n",
        "X_test = data[data[time_column_name] >= '2017-02-01']\n",
        "y_train = X_train.pop(target_column_name).values\n",
        "y_test = X_test.pop(target_column_name).values"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Train\n",
        "\n",
        "Instantiate a AutoMLConfig object. This defines the settings and data used to run the experiment.\n",
        "\n",
        "|Property|Description|\n",
        "|-|-|\n",
        "|**task**|forecasting|\n",
        "|**primary_metric**|This is the metric that you want to optimize.<br> Forecasting supports the following primary metrics <br><i>spearman_correlation</i><br><i>normalized_root_mean_squared_error</i><br><i>r2_score</i><br><i>normalized_mean_absolute_error</i>\n",
        "|**iterations**|Number of iterations. In each iteration, Auto ML trains a specific pipeline on the given data|\n",
        "|**iteration_timeout_minutes**|Time limit in minutes for each iteration.|\n",
        "|**X**|(sparse) array-like, shape = [n_samples, n_features]|\n",
        "|**y**|(sparse) array-like, shape = [n_samples, ], targets values.|\n",
        "|**n_cross_validations**|Number of cross validation splits.|\n",
        "|**path**|Relative path to the project folder.  AutoML stores configuration files for the experiment under this folder. You can specify a new empty folder. "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "automl_settings = {\n",
        "    \"time_column_name\": time_column_name    \n",
        "}\n",
        "\n",
        "\n",
        "automl_config = AutoMLConfig(task = 'forecasting',\n",
        "                             debug_log = 'automl_nyc_energy_errors.log',\n",
        "                             primary_metric='normalized_root_mean_squared_error',\n",
        "                             iterations = 10,\n",
        "                             iteration_timeout_minutes = 5,\n",
        "                             X = X_train,\n",
        "                             y = y_train,\n",
        "                             n_cross_validations = 3,\n",
        "                             path=project_folder,\n",
        "                             verbosity = logging.INFO,\n",
        "                            **automl_settings)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Submitting the configuration will start a new run in this experiment. For local runs, the execution is synchronous. Depending on the data and number of iterations, this can run for a while. Parameters controlling concurrency may speed up the process, depending on your hardware.\n",
        "\n",
        "You will see the currently running iterations printing to the console."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "local_run = experiment.submit(automl_config, show_output=True)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "local_run"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Retrieve the Best Model\n",
        "Below we select the best pipeline from our iterations. The get_output method on automl_classifier returns the best run and the fitted model for the last fit invocation. There are overloads on get_output that allow you to retrieve the best run and fitted model for any logged metric or a particular iteration."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "best_run, fitted_model = local_run.get_output()\n",
        "fitted_model.steps"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### View the engineered names for featurized data\n",
        "Below we display the engineered feature names generated for the featurized data using the time-series featurization."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "fitted_model.named_steps['timeseriestransformer'].get_engineered_feature_names()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Test the Best Fitted Model\n",
        "\n",
        "For forecasting, we will use the `forecast` function instead of the `predict` function. There are two reasons for this.\n",
        "\n",
        "We need to pass the recent values of the target variable `y`, whereas the scikit-compatible `predict` function only takes the non-target variables `X`. In our case, the test data immediately follows the training data, and we fill the `y` variable with `NaN`. The `NaN` serves as a question mark for the forecaster to fill with the actuals. Using the forecast function will produce forecasts using the shortest possible forecast horizon. The last time at which a definite (non-NaN) value is seen is the _forecast origin_ - the last time when the value of the target is known. \n",
        "\n",
        "Using the `predict` method would result in getting predictions for EVERY horizon the forecaster can predict at. This is useful when training and evaluating the performance of the forecaster at various horizons, but the level of detail is excessive for normal use."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Replace ALL values in y_pred by NaN. \n",
        "# The forecast origin will be at the beginning of the first forecast period\n",
        "# (which is the same time as the end of the last training period).\n",
        "y_query = y_test.copy().astype(np.float)\n",
        "y_query.fill(np.nan)\n",
        "# The featurized data, aligned to y, will also be returned.\n",
        "# This contains the assumptions that were made in the forecast\n",
        "# and helps align the forecast to the original data\n",
        "y_fcst, X_trans = fitted_model.forecast(X_test, y_query)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# limit the evaluation to data where y_test has actuals\n",
        "def align_outputs(y_predicted, X_trans, X_test, y_test, predicted_column_name = 'predicted'):\n",
        "    \"\"\"\n",
        "    Demonstrates how to get the output aligned to the inputs\n",
        "    using pandas indexes. Helps understand what happened if\n",
        "    the output's shape differs from the input shape, or if\n",
        "    the data got re-sorted by time and grain during forecasting.\n",
        "    \n",
        "    Typical causes of misalignment are:\n",
        "    * we predicted some periods that were missing in actuals -> drop from eval\n",
        "    * model was asked to predict past max_horizon -> increase max horizon\n",
        "    * data at start of X_test was needed for lags -> provide previous periods\n",
        "    \"\"\"\n",
        "    df_fcst = pd.DataFrame({predicted_column_name : y_predicted})\n",
        "    # y and X outputs are aligned by forecast() function contract\n",
        "    df_fcst.index = X_trans.index\n",
        "    \n",
        "    # align original X_test to y_test    \n",
        "    X_test_full = X_test.copy()\n",
        "    X_test_full[target_column_name] = y_test\n",
        "\n",
        "    # X_test_full's does not include origin, so reset for merge\n",
        "    df_fcst.reset_index(inplace=True)\n",
        "    X_test_full = X_test_full.reset_index().drop(columns='index')\n",
        "    together = df_fcst.merge(X_test_full, how='right')\n",
        "    \n",
        "    # drop rows where prediction or actuals are nan \n",
        "    # happens because of missing actuals \n",
        "    # or at edges of time due to lags/rolling windows\n",
        "    clean = together[together[[target_column_name, predicted_column_name]].notnull().all(axis=1)]\n",
        "    return(clean)\n",
        "\n",
        "df_all = align_outputs(y_fcst, X_trans, X_test, y_test)\n",
        "df_all.head()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Looking at `X_trans` is also useful to see what featurization happened to the data."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "X_trans"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Calculate accuracy metrics\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "def MAPE(actual, pred):\n",
        "    \"\"\"\n",
        "    Calculate mean absolute percentage error.\n",
        "    Remove NA and values where actual is close to zero\n",
        "    \"\"\"\n",
        "    not_na = ~(np.isnan(actual) | np.isnan(pred))\n",
        "    not_zero = ~np.isclose(actual, 0.0)\n",
        "    actual_safe = actual[not_na & not_zero]\n",
        "    pred_safe = pred[not_na & not_zero]\n",
        "    APE = 100*np.abs((actual_safe - pred_safe)/actual_safe)\n",
        "    return np.mean(APE)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "print(\"Simple forecasting model\")\n",
        "rmse = np.sqrt(mean_squared_error(df_all[target_column_name], df_all['predicted']))\n",
        "print(\"[Test Data] \\nRoot Mean squared error: %.2f\" % rmse)\n",
        "mae = mean_absolute_error(df_all[target_column_name], df_all['predicted'])\n",
        "print('mean_absolute_error score: %.2f' % mae)\n",
        "print('MAPE: %.2f' % MAPE(df_all[target_column_name], df_all['predicted']))\n",
        "\n",
        "# Plot outputs\n",
        "%matplotlib notebook\n",
        "test_pred = plt.scatter(df_all[target_column_name], df_all['predicted'], color='b')\n",
        "test_test = plt.scatter(y_test, y_test, color='g')\n",
        "plt.legend((test_pred, test_test), ('prediction', 'truth'), loc='upper left', fontsize=8)\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The distribution looks a little heavy tailed: we underestimate the excursions of the extremes. A normal-quantile transform of the target might help, but let's first try using some past data with the lags and rolling window transforms.\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Using lags and rolling window features to improve the forecast"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "We did not use lags in the previous model specification. In effect, the prediction was the result of a simple regression on date, grain and any additional features. This is often a very good prediction as common time series patterns like seasonality and trends can be captured in this manner. Such simple regression is horizon-less: it doesn't matter how far into the future we are predicting, because we are not using past data.\n",
        "\n",
        "Now that we configured target lags, that is the previous values of the target variables, and the prediction is no longer horizon-less. We therefore must specify the `max_horizon` that the model will learn to forecast. The `target_lags` keyword specifies how far back we will construct the lags of the target variable, and the `target_rolling_window_size` specifies the size of the rolling window over which we will generate the `max`, `min` and `sum` features."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "automl_settings_lags = {\n",
        "    'time_column_name': time_column_name,\n",
        "    'target_lags': 1,\n",
        "    'target_rolling_window_size': 5,\n",
        "    # you MUST set the max_horizon when using lags and rolling windows\n",
        "    # it is optional when looking-back features are not used \n",
        "    'max_horizon': len(y_test), # only one grain\n",
        "}\n",
        "\n",
        "\n",
        "automl_config_lags = AutoMLConfig(task = 'forecasting',\n",
        "                             debug_log = 'automl_nyc_energy_errors.log',\n",
        "                             primary_metric='normalized_root_mean_squared_error',\n",
        "                             iterations = 10,\n",
        "                             iteration_timeout_minutes = 5,\n",
        "                             X = X_train,\n",
        "                             y = y_train,\n",
        "                             n_cross_validations = 3,\n",
        "                             path=project_folder,\n",
        "                             verbosity = logging.INFO,\n",
        "                            **automl_settings_lags)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "local_run_lags = experiment.submit(automl_config_lags, show_output=True)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "best_run_lags, fitted_model_lags = local_run_lags.get_output()\n",
        "y_fcst_lags, X_trans_lags = fitted_model_lags.forecast(X_test, y_query)\n",
        "df_lags = align_outputs(y_fcst_lags, X_trans_lags, X_test, y_test)\n",
        "df_lags.head()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "X_trans_lags"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "print(\"Forecasting model with lags\")\n",
        "rmse = np.sqrt(mean_squared_error(df_lags[target_column_name], df_lags['predicted']))\n",
        "print(\"[Test Data] \\nRoot Mean squared error: %.2f\" % rmse)\n",
        "mae = mean_absolute_error(df_lags[target_column_name], df_lags['predicted'])\n",
        "print('mean_absolute_error score: %.2f' % mae)\n",
        "print('MAPE: %.2f' % MAPE(df_lags[target_column_name], df_lags['predicted']))\n",
        "\n",
        "# Plot outputs\n",
        "%matplotlib notebook\n",
        "test_pred = plt.scatter(df_lags[target_column_name], df_lags['predicted'], color='b')\n",
        "test_test = plt.scatter(y_test, y_test, color='g')\n",
        "plt.legend((test_pred, test_test), ('prediction', 'truth'), loc='upper left', fontsize=8)\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### What features matter for the forecast?"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "from azureml.train.automl.automlexplainer import explain_model\n",
        "\n",
        "# feature names are everything in the transformed data except the target\n",
        "features = X_trans.columns[:-1]\n",
        "expl = explain_model(fitted_model, X_train, X_test, features = features, best_run=best_run_lags, y_train = y_train)\n",
        "# unpack the tuple\n",
        "shap_values, expected_values, feat_overall_imp, feat_names, per_class_summary, per_class_imp = expl\n",
        "best_run_lags"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Please go to the Azure Portal's best run to see the top features chart.\n",
        "\n",
        "The informative features make all sorts of intuitive sense. Temperature is a strong driver of heating and cooling demand in NYC. Apart from that, the daily life cycle, expressed by `hour`, and the weekly cycle, expressed by `wday` drives people's energy use habits."
      ]
    }
  ],
  "metadata": {
    "authors": [
      {
        "name": "xiaga, tosingli"
      }
    ],
    "kernelspec": {
      "display_name": "Python 3.6",
      "language": "python",
      "name": "python36"
    },
    "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.6.7"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 2
}