notebooks/5_Model_Selection_and_Evaluation.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-info\">\n",
    "Section of the book chapter: <b>5.3 Model Selection, Optimization and Evaluation</b>\n",
    "</div>\n",
    "\n",
    "# 5. Model Selection and Evaluation\n",
    "\n",
    "**Table of Contents**\n",
    "\n",
    "* [5.1 Hyperparameter Optimization](#5.1-Hyperparameter-Optimization)\n",
    "* [5.2 Model Evaluation](#5.2-Model-Evaluation)\n",
    "\n",
    "**Learnings:**\n",
    "\n",
    "- how to optimize machine learning (ML) models with grid search, random search and Bayesian optimization,\n",
    "- how to evaluate ML models.\n",
    "\n",
    "\n",
    "\n",
    "### Packages"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "%config InlineBackend.figure_format = 'retina'\n",
    "\n",
    "# ignore warnings\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib as mpl\n",
    "from sklearn.ensemble import RandomForestRegressor\n",
    "\n",
    "import utils"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Read in Data\n",
    "\n",
    "**Dataset:** Felix M. Riese and Sina Keller, \"Hyperspectral benchmark dataset on soil moisture\", Dataset, Zenodo, 2018. [DOI:10.5281/zenodo.1227836](http://doi.org/10.5281/zenodo.1227836) and [GitHub](https://github.com/felixriese/hyperspectral-soilmoisture-dataset)\n",
    "\n",
    "**Introducing paper:** Felix M. Riese and Sina Keller, “Introducing a Framework of Self-Organizing Maps for Regression of Soil Moisture with Hyperspectral Data,” in IGARSS 2018 - 2018 IEEE International Geoscience and Remote Sensing Symposium, Valencia, Spain, 2018, pp. 6151-6154. [DOI:10.1109/IGARSS.2018.8517812](https://doi.org/10.1109/IGARSS.2018.8517812)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_train, X_test, y_train, y_test = utils.get_xy_split()\n",
    "\n",
    "print(X_train.shape, X_test.shape, y_train.shape, y_test.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Fix Random State"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.seed(42)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "***\n",
    "\n",
    "## 5.1 Hyperparameter Optimization\n",
    "\n",
    "Content:\n",
    "\n",
    "- [5.1.1 Grid Search](#5.1.1-Grid-Search)\n",
    "- [5.1.2 Randomized Search](#5.1.2-Randomized-Search)\n",
    "- [5.1.3 Bayesian Optimization](#5.1.3-Bayesian-Optimization)\n",
    "\n",
    "### 5.1.1 Grid Search"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# NBVAL_IGNORE_OUTPUT\n",
    "\n",
    "from sklearn.svm import SVR\n",
    "from sklearn.model_selection import GridSearchCV\n",
    "\n",
    "# example mode: support vector regressor\n",
    "model = SVR(kernel=\"rbf\")\n",
    "\n",
    "# define parameter grid to be tested\n",
    "params = {\n",
    "    \"C\": np.logspace(-4, 4, 9),\n",
    "    \"gamma\": np.logspace(-4, 4, 9)}\n",
    "\n",
    "\n",
    "# set up grid search and run it on the data\n",
    "gs = GridSearchCV(model, params)\n",
    "%timeit gs.fit(X_train, y_train)\n",
    "print(\"R2 score = {0:.2f} %\".format(gs.score(X_test, y_test)*100))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5.1.2 Randomized Search"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# NBVAL_IGNORE_OUTPUT\n",
    "\n",
    "from sklearn.svm import SVR\n",
    "from sklearn.model_selection import RandomizedSearchCV\n",
    "\n",
    "# example mode: support vector regressor\n",
    "model = SVR(kernel=\"rbf\")\n",
    "\n",
    "# define parameter grid to be tested\n",
    "params = {\n",
    "    \"C\": np.logspace(-4, 4, 9),\n",
    "    \"gamma\": np.logspace(-4, 4, 9)}\n",
    "\n",
    "# set up grid search and run it on the data\n",
    "gsr = RandomizedSearchCV(model, params, n_iter=15, refit=True)\n",
    "%timeit gsr.fit(X_train, y_train)\n",
    "print(\"R2 score = {0:.2f} %\".format(gsr.score(X_test, y_test)*100))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5.1.3 Bayesian Optimization\n",
    "\n",
    "Implementation: [github.com/fmfn/BayesianOptimization](https://github.com/fmfn/BayesianOptimization)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# NBVAL_IGNORE_OUTPUT\n",
    "\n",
    "from sklearn.svm import SVR\n",
    "from bayes_opt import BayesianOptimization\n",
    "\n",
    "# define function to be optimized\n",
    "def opt_func(C, gamma):\n",
    "    model = SVR(C=C, gamma=gamma)\n",
    "    return model.fit(X_train, y_train).score(X_test, y_test)\n",
    "\n",
    "# set bounded region of parameter space\n",
    "pbounds = {'C': (1e-5, 1e4), 'gamma': (1e-5, 1e4)}\n",
    "\n",
    "# define optimizer\n",
    "optimizer = BayesianOptimization(\n",
    "    f=opt_func,\n",
    "    pbounds=pbounds,\n",
    "    random_state=1)\n",
    "\n",
    "# optimize\n",
    "%time optimizer.maximize(init_points=2, n_iter=15)\n",
    "print(\"R2 score = {0:.2f} %\".format(optimizer.max[\"target\"]*100))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "***\n",
    "\n",
    "## 5.2 Model Evaluation\n",
    "\n",
    "Content:\n",
    "\n",
    "- [5.2.1 Generate Exemplary Data](#5.2.1-Generate-Exemplary-Data)\n",
    "- [5.2.2 Plot the Data](#5.2.2-Plot-the-Data)\n",
    "- [5.2.3 Evaluation Metrics](#5.2.3-Evaluation-Metrics)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import sklearn.metrics as me"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5.2.1 Generate Exemplary Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### generate example data\n",
    "np.random.seed(1)\n",
    "\n",
    "# define x grid\n",
    "x_grid = np.linspace(0, 10, 11)\n",
    "y_model = x_grid*0.5\n",
    "\n",
    "# define first dataset without outlier\n",
    "y1 = np.array([y + np.random.normal(scale=0.2) for y in y_model])\n",
    "\n",
    "# define second dataset with outlier\n",
    "y2 = np.copy(y1)\n",
    "y2[9] = 0.5\n",
    "\n",
    "# define third dataset with higher variance\n",
    "y3 = np.array([y + np.random.normal(scale=1.0) for y in y_model])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5.2.2 Plot the Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot example data\n",
    "fig, (ax1, ax2, ax3) = plt.subplots(1,3, figsize=(12,4))\n",
    "fontsize = 18\n",
    "titleweight = \"bold\"\n",
    "titlepad = 10\n",
    "\n",
    "scatter_label = \"Data\"\n",
    "scatter_alpha = 0.7\n",
    "scatter_s = 100\n",
    "ax1.scatter(x_grid, y1, label=scatter_label, alpha=scatter_alpha, s=scatter_s)\n",
    "ax1.set_title(\"(a) Low var.\", fontsize=fontsize, fontweight=titleweight, pad=titlepad)\n",
    "\n",
    "ax2.scatter(x_grid, y2, label=scatter_label, alpha=scatter_alpha, s=scatter_s)\n",
    "ax2.set_title(\"(b) Low var. + outlier\", fontsize=fontsize, fontweight=titleweight, pad=titlepad)\n",
    "\n",
    "ax3.scatter(x_grid, y3, label=scatter_label, alpha=scatter_alpha, s=scatter_s)\n",
    "ax3.set_title(\"(c) Higher var.\", fontsize=fontsize, fontweight=titleweight, pad=titlepad)\n",
    "\n",
    "for i, ax in enumerate([ax1, ax2, ax3]):\n",
    "    i += 1\n",
    "    \n",
    "    # red line\n",
    "    ax.plot(x_grid, y_model, label=\"Model\", c=\"tab:red\", linestyle=\"dashed\", linewidth=4, alpha=scatter_alpha)\n",
    "    \n",
    "    # x-axis cosmetics\n",
    "    ax.set_xlabel(\"x in a.u.\", fontsize=fontsize)\n",
    "    for tick in ax.xaxis.get_major_ticks():\n",
    "        tick.label.set_fontsize(fontsize) \n",
    "    \n",
    "    # y-axis cosmetics\n",
    "    if i != 1:\n",
    "        ax.set_yticklabels([])\n",
    "    else:\n",
    "        ax.set_ylabel(\"y in a.u.\", fontsize=fontsize, rotation=90)\n",
    "        for tick in ax.yaxis.get_major_ticks():\n",
    "            tick.label.set_fontsize(fontsize) \n",
    "    ax.set_xlim(-0.5, 10.5)\n",
    "    ax.set_ylim(-0.5, 6.5)\n",
    "    # ax.set_title(\"Example \"+str(i), fontsize=fontsize)\n",
    "    if i == 2:\n",
    "        ax.legend(loc=2, fontsize=fontsize*1.0, frameon=True)\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.savefig(\"plots/metrics_plot.pdf\", bbox_inches=\"tight\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5.2.3 Evaluation Metrics"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# calculating the metrics\n",
    "for i, y in enumerate([y1, y2, y3]):\n",
    "    print(\"Example\", i+1)\n",
    "    print(\"- MAE = {:.2f}\".format(me.mean_absolute_error(y_model, y)))\n",
    "    print(\"- MSE = {:.2f}\".format(me.mean_squared_error(y_model, y)))\n",
    "    print(\"- RMSE = {:.2f}\".format(np.sqrt(me.mean_squared_error(y_model, y))))\n",
    "    print(\"- R2 = {:.2f}%\".format(me.r2_score(y_model, y)*100))\n",
    "    print(\"-\"*20)"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "# print out for LaTeX table\n",
    "\n",
    "descriptions = {\n",
    "    1: \"Low variance\",\n",
    "    2: \"Low variance and one outlier\",\n",
    "    3: \"Higher variance\",}\n",
    "bold = [[False, False, False, False], [False, True, True, True], [True, False, False, False]]\n",
    "def make_bold(is_bold):\n",
    "    if is_bold:\n",
    "        return \"\\\\bfseries\"\n",
    "    return \"\"\n",
    "\n",
    "for i, y in enumerate([y1, y2, y3]):\n",
    "    print(\"{description} & {bold1} {mae:.2f} & {bold2} {mse:.2f} & {bold3} {rmse:.2f} & {bold4} {r2:.2f} \\\\\\\\\".format(\n",
    "        description=descriptions[i+1],\n",
    "        mae=me.mean_absolute_error(y_model, y),\n",
    "        mse=me.mean_squared_error(y_model, y),\n",
    "        rmse=np.sqrt(me.mean_squared_error(y_model, y)),\n",
    "        r2=me.r2_score(y_model, y)*100,\n",
    "        bold1=make_bold(bold[i][0]),\n",
    "        bold2=make_bold(bold[i][1]),\n",
    "        bold3=make_bold(bold[i][2]),\n",
    "        bold4=make_bold(bold[i][3]),))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.9.1"
  },
  "toc": {
   "nav_menu": {},
   "number_sections": false,
   "sideBar": true,
   "skip_h1_title": false,
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": "block",
   "toc_window_display": false
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}