07c - BlackBoxOptimization.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "4pJysFMoQ2Du"
   },
   "source": [
    "*Credit*: This example is taken from the `skopt` [Bayesian Optimization tutorial](https://scikit-optimize.github.io/notebooks/bayesian-optimization.html)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 0,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "C135LrgaNh2P"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "!pip install scikit-optimize\n",
    "from skopt import gp_minimize"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "12l1n44ZNh2S"
   },
   "source": [
    "# Black-box function optimization with skopt\n",
    "\n",
    "Scikit-Optimize, or [skopt](https://scikit-optimize.github.io/), is a simple and efficient library to minimize (very) expensive and noisy black-box functions. It implements several methods for sequential model-based optimization.\n",
    "\n",
    "Alternative libraries include [Spearmint](https://github.com/HIPS/Spearmint),  [Hyperopt](http://hyperopt.github.io/hyperopt/), [Optuna](https://optuna.org) and [BoTorch](https://github.com/pytorch/botorch).\n",
    "\n",
    "Black-box algorithms do not need any knowledge of the gradient. These libraries provide algorithms that are more powerful and scale better than the *Nelder-Mead* simplex algorithm we saw in the *Constrained Optimization* notebook. Modern black-box (or sequential model-based) optimization algorithms are increasingly popular for optimizing the *hyperparameters* (user-tuned \"knobs\") of machine learning models."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "q4IljA6khhI2"
   },
   "source": [
    "Let's assume the following noisy function $f$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 0,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "uVPjKKszhhI3"
   },
   "outputs": [],
   "source": [
    "noise_level = 0.1\n",
    "\n",
    "def f(x, noise_level=noise_level):\n",
    "    return np.sin(5 * x[0]) * (1 - np.tanh(x[0] ** 2)) + np.random.randn() * noise_level"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "L44jv7DLhhI4"
   },
   "source": [
    "In `skopt`, functions $f$ are assumed to take as input a 1D vector $x$ represented as an array-like and to return a scalar $f(x)$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 0,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "gkzcS1wwhhI5"
   },
   "outputs": [],
   "source": [
    "# Plot f(x) + contours\n",
    "x = np.linspace(-2, 2, 400).reshape(-1, 1)\n",
    "fx = [f(x_i, noise_level=0.0) for x_i in x]\n",
    "plt.plot(x, fx, \"r--\", label=\"True (unknown)\")\n",
    "plt.fill(np.concatenate([x, x[::-1]]),\n",
    "         np.concatenate(([fx_i - 1.9600 * noise_level for fx_i in fx], \n",
    "                         [fx_i + 1.9600 * noise_level for fx_i in fx[::-1]])),\n",
    "         alpha=.2, fc=\"r\", ec=\"None\")\n",
    "plt.legend()\n",
    "plt.grid()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "zOdmwQGVhhI7"
   },
   "source": [
    "Note that this function is differentiable and wouldn't actually be that difficult to differentiate analytically. But let's pretend for now that it's impossible to differentiate it and this precludes gradient-based optimization.\n",
    "\n",
    "Bayesian Optimization based on Gaussian Process regression is implemented in `skopt.gp_minimize` and can be carried out as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 0,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "nR8lj4qbhhI7"
   },
   "outputs": [],
   "source": [
    "res = gp_minimize(f,                  # the function to minimize\n",
    "                  [(-2.0, 2.0)],      # the bounds on each dimension of x\n",
    "                  acq_func=\"EI\",      # the acquisition function\n",
    "                  n_calls=15,         # the number of evaluations of f \n",
    "                  n_random_starts=5,  # the number of random initialization points\n",
    "                  noise=0.1**2,       # the noise level (optional)\n",
    "                  random_state=123)   # the random seed"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "RugtPENehhI9"
   },
   "source": [
    "Accordingly, the approximated minimum is found to be:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 0,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "XRouBMV9hhI9"
   },
   "outputs": [],
   "source": [
    "print(\"x^*=%.4f, f(x^*)=%.4f\" % (res.x[0], res.fun))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "JMi5vSfzhhI_"
   },
   "source": [
    "For further inspection of the results, attributes of the res named tuple provide the following information:\n",
    "\n",
    "- `x [float]`: location of the minimum.\n",
    "- `fun [float]`: function value at the minimum.\n",
    "- `models`: surrogate models used for each iteration.\n",
    "- `x_iters [array]`: location of function evaluation for each iteration.\n",
    "- `func_vals [array]`: function value for each iteration.\n",
    "- `space [Space]`: the optimization space.\n",
    "- `specs [dict]`: parameters passed to the function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 0,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "PfBIZvdohhJA"
   },
   "outputs": [],
   "source": [
    "print(res)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "Io5PrMuRhhJB"
   },
   "source": [
    "Together these attributes can be used to visually inspect the results of the minimization, such as the convergence trace or the acquisition function at the last iteration:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 0,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "Sq1QUZrbhhJC"
   },
   "outputs": [],
   "source": [
    "from skopt.plots import plot_convergence\n",
    "plot_convergence(res);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "u3JhgSEXFC5q"
   },
   "source": [
    "## Branin-Hoo function\n",
    "\n",
    "The [Branin-Hoo function](http://www.sfu.ca/~ssurjano/branin.html) is a popular optimization benchmark. \n",
    "\n",
    "The 2-dimensional function has the form\n",
    "$$\n",
    "f(\\mathbf{x}) = a(x_2 - bx_1^2 + cx_1 - r)^2 + s(1-t)cos(x_1) + s .\\\n",
    "$$\n",
    "\n",
    "The recommended values of $a$, $b$, $c$, $r$, $s$ and $t$ are $a=1$, $b = 5.1 / (4 \\pi^2)$, $c= 5/\\pi$, $r=6$, $s=10$ and $t=1/(8\\pi)$. \n",
    "\n",
    "The function is usually evaluated on the square $x_1 \\in [-5, 10], x_2 \\in [0, 15]$.\n",
    "\n",
    "It has three global minima. Using the recommended parameters above, $f(\\mathbf{x}^*)=0.397887$ at $\\mathbf{x}^* = (-\\pi, 12.275), (\\pi, 2.275)$ and $(9.42478, 2.475)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 0,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "rXk3sRHLGhpR"
   },
   "outputs": [],
   "source": [
    "def branin(x, a=1, b=5.1 / (4 * np.pi**2), c=5. / np.pi,\n",
    "           r=6, s=10, t=1. / (8 * np.pi)):\n",
    "    \"\"\"Branin-Hoo function is defined on the square x1 ∈ [-5, 10], x2 ∈ [0, 15].\n",
    "\n",
    "    It has three minima with f(x*) = 0.397887 at x* = (-pi, 12.275),\n",
    "    (+pi, 2.275), and (9.42478, 2.475).\n",
    "\n",
    "    More details: <http://www.sfu.ca/~ssurjano/branin.html>\n",
    "    \"\"\"\n",
    "    return (a * (x[1] - b * x[0] ** 2 + c * x[0] - r) ** 2 +\n",
    "            s * (1 - t) * np.cos(x[0]) + s)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "XEela3dxGq1x"
   },
   "source": [
    "### Exercise\n",
    "\n",
    "Adapt the Bayesian Optimization with GP regression example above to solve the Branin-Hoo function. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "OVuoMxRuHjfH"
   },
   "source": [
    "### Exercise\n",
    "\n",
    "Experiment with different settings of the acquisition function (parameter `acq_func`). What is an acquisition function? You can find a great write-up on Bayesian Optimization on [Martin Krasser's blog](http://krasserm.github.io/2018/03/21/bayesian-optimization/)."
   ]
  }
 ],
 "metadata": {
  "colab": {
   "collapsed_sections": [],
   "name": "BlackBoxOptimization.ipynb",
   "provenance": [
    {
     "file_id": "1txpNobA9vwwCtR-EPoGLkBz4rwLIMHWY",
     "timestamp": 1569987035009
    },
    {
     "file_id": "1FMM-aAALAeLY-ctwPOmjeadN2eFWjsHC",
     "timestamp": 1569985457190
    },
    {
     "file_id": "1u9v8Xc7KilebO7Wti1ysIJWQmWNi8VS2",
     "timestamp": 1569981925535
    }
   ]
  },
  "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.6.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}