{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Entropy Search\n",
    "### Written by Simon Bartels, Max Planck Institute for Intelligent Systems, and Andrei Paleyes, Amazon.com\n",
    "This notebook demonstrates how to use Entropy Search (ES) in GPyOpt and compares it to Expected Improvement (EI). For details on ES have a look at the original paper:\n",
    "\n",
    "Hennig and C. J. Schuler. Entropy search for information-efficient global optimization. Journal of Machine Learning Research, 13, 2012"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import GPy\n",
    "import GPyOpt\n",
    "from GPyOpt.models.gpmodel import GPModel\n",
    "from GPyOpt.core.task.space import Design_space, bounds_to_space\n",
    "from GPyOpt.util.mcmc_sampler import AffineInvariantEnsembleSampler\n",
    "from GPyOpt.acquisitions.ES import AcquisitionEntropySearch\n",
    "from GPyOpt.acquisitions.EI import AcquisitionEI\n",
    "\n",
    "import matplotlib as mpl\n",
    "mpl.use('Agg')\n",
    "\n",
    "#configure plotting\n",
    "%matplotlib inline\n",
    "%config InlineBackend.figure_format = 'svg'\n",
    "\n",
    "import matplotlib;matplotlib.rcParams['figure.figsize'] = (8,5)\n",
    "import matplotlib;matplotlib.rcParams['text.usetex'] = True\n",
    "import matplotlib;matplotlib.rcParams['font.size'] = 16\n",
    "import matplotlib;matplotlib.rcParams['font.family'] = 'serif'\n",
    "from matplotlib import pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Toy Problem\n",
    "The following toy problem demonstrates the possible advantage Entropy Search can have over Expected Improvement. The observations are chosen in a way such that EI will evaluate at the minimum whose location is pretty clear. Entropy Search on the other hand exhibits a more explorative behavior. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "X = np.array([[-1], [1], [2]])\n",
    "y = 2 * -np.array([[.1], [.5], [.5]])\n",
    "bounds = [(-5, 5)]\n",
    "input_dim = X.shape[1]\n",
    "\n",
    "kern = GPy.kern.RBF(input_dim, variance=1., lengthscale=1.)\n",
    "model = GPModel(kern, noise_var=1e-3, max_iters=0, optimize_restarts=0)\n",
    "\n",
    "model.updateModel(X, y, None, None)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plot of Data-Set, Model and Acquisition Functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Xs = np.arange(bounds[0][0], bounds[0][1], 0.01).reshape(-1, 1)\n",
    "ys, vs = model.predict(Xs)\n",
    "\n",
    "plt.fill_between(np.ndarray.flatten(Xs), \n",
    "                 np.ndarray.flatten(ys+np.sqrt(vs)), \n",
    "                 np.ndarray.flatten(ys-np.sqrt(vs)), alpha=0.1)\n",
    "plt.plot(Xs, ys, color='b')\n",
    "plt.plot(X, y, 'x')\n",
    "\n",
    "space = Design_space(bounds_to_space(bounds))\n",
    "def normalize(vs):\n",
    "    return (vs - min(vs))/(max(vs - min(vs)))\n",
    "sampler = AffineInvariantEnsembleSampler(space)\n",
    "\n",
    "ei = AcquisitionEI(model, space)\n",
    "vei = normalize(ei.acquisition_function(Xs))\n",
    "\n",
    "es = AcquisitionEntropySearch(model, space, sampler)\n",
    "ves = normalize(es.acquisition_function(Xs))\n",
    "\n",
    "# plot Expected Improvement again\n",
    "plt.plot(Xs, ves, color='r')\n",
    "# plot Entropy Search values\n",
    "plt.plot(Xs, vei, color='g')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Expected Improvement (green line) suggests to evaluate in the location of the minimum (around 1.9). In contrast, Entropy Search (red line) is more explorative, preferring points near 4 and -4. Evaluating the minimum location would not bring much insight.\n",
    "\n",
    "## Comparison on the Branin function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# --- Function to optimize\n",
    "func  = GPyOpt.objective_examples.experiments2d.branin()\n",
    "func.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's define the necessary objects for `ModularBayesianOptimization`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "objective = GPyOpt.core.task.SingleObjective(func.f)\n",
    "space = GPyOpt.Design_space(space =[{'name': 'var_1', 'type': 'continuous', 'domain': (-5,10)},\n",
    "                                    {'name': 'var_2', 'type': 'continuous', 'domain': (1,15)}])\n",
    "acquisition_optimizer = GPyOpt.optimization.AcquisitionOptimizer(space)\n",
    "initial_design = GPyOpt.experiment_design.initial_design('random', space, 5)\n",
    "max_iter = 10"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First run Expected Improvement."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ei_model = GPyOpt.models.GPModel(optimize_restarts=5,verbose=False)\n",
    "ei = AcquisitionEI(ei_model, space, optimizer=acquisition_optimizer)\n",
    "ei_evaluator = GPyOpt.core.evaluators.Sequential(ei)\n",
    "bo_ei = GPyOpt.methods.ModularBayesianOptimization(ei_model, space, objective, ei, ei_evaluator, initial_design)\n",
    "bo_ei.run_optimization(max_iter = max_iter)\n",
    "bo_ei.plot_acquisition()\n",
    "bo_ei.plot_convergence()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And now run Entropy Search."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "es_model = GPyOpt.models.GPModel(optimize_restarts=5,verbose=False)\n",
    "ei = AcquisitionEI(es_model, space, optimizer=acquisition_optimizer)\n",
    "proposal_function = lambda x : np.clip(np.log(ei._compute_acq(x)), 0., np.PINF)\n",
    "sampler = AffineInvariantEnsembleSampler(space)\n",
    "es = AcquisitionEntropySearch(es_model, space, sampler, optimizer=acquisition_optimizer, num_representer_points=10, \n",
    "                   burn_in_steps=10, num_samples=100, proposal_function = proposal_function)\n",
    "es_evaluator = GPyOpt.core.evaluators.Sequential(es)\n",
    "bo_es = GPyOpt.methods.ModularBayesianOptimization(es_model, space, objective, es, es_evaluator, initial_design)\n",
    "bo_es.run_optimization(max_iter = max_iter)\n",
    "bo_es.plot_acquisition()\n",
    "bo_es.plot_convergence()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Let's plot the locations where Entropy Search (circles) and Expected Improvement (crosses) evaluated."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bounds = func.bounds\n",
    "x1 = np.linspace(bounds[0][0], bounds[0][1], 100)\n",
    "x2 = np.linspace(bounds[1][0], bounds[1][1], 100)\n",
    "X1, X2 = np.meshgrid(x1, x2)\n",
    "X = np.hstack((X1.reshape(100*100,1),X2.reshape(100*100,1)))\n",
    "Y = func.f(X)\n",
    "\n",
    "plt.figure()    \n",
    "plt.contourf(X1, X2, Y.reshape((100,100)),100)\n",
    "plt.plot(np.array(func.min)[:,0], np.array(func.min)[:,1], 'w.', markersize=20, label=u'Observations')\n",
    "plt.colorbar()\n",
    "plt.plot(ei_model.model.X[:, 0],ei_model.model.X[:, 1], 'o')\n",
    "plt.plot(es_model.model.X[:, 0],es_model.model.X[:, 1], 'x')\n",
    "plt.xlabel('X1')\n",
    "plt.ylabel('X2')\n",
    "plt.title(func.name)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "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.6.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}