Finished slide deck for w06 HPO BO.

w06_hpo_bo/images/intro_images/blackbox_HPO.tex

@@ -0,0 +1,21 @@
+\newcommand{\myblackbox}{\fcolorbox{black}{black}{
+    \minipage[t]{\dimexpr0.111\linewidth-2\fboxsep-2\fboxrule\relax}
+        ~~~\\
+        ~~~\\
+        ~~~\\
+    \endminipage}}
+
+
+    	\begin{tikzpicture}
+\tikzstyle{every node}=[draw,fill=white,minimum width=0cm,thin]
+\tikzstyle{every path}=[-latex,ultra thick]
+\node (A) [draw=white]{{\Huge{$\conf$}}};
+\node (B) [right=14mm of A,draw=white] {\myblackbox{}};
+\node (C) [right=14mm of B,draw=white] {{\Huge{$f(\conf)$}}};
+%\node (D) [below=7mm of B, align=center, fill=black!10] {\large{Bayesian}\\\large{optimization}};
+
+\draw ($(A.east)+(0.2,0.0)$) -- ($(B.west)+(-0.2,0.0)$);
+\draw ($(B.east)+(0.2,0.0)$) -- ($(C.west)+(-0.2,0.0)$);
+%\draw ($(C.south)+(0.0,-0.2)$) -| ++(0.0,0.0) |- ($(D.east)+(0.2,0.0)$);
+%\draw ($(D.west)+(-0.2,0.0)$) |- ++(0.0,0.0) -| ($(A.south)+(0.0,-0.2)$);
+\end{tikzpicture}

w06_hpo_bo/scripts/bo_configurations.py

@@ -6,6 +6,7 @@
 colors = dict({
     'observations': 'black',
     'highlighted_observations': 'green',
+    'new_observation': 'red',
     'current_incumbent': 'red',
     'highlighted_point': 'red',
     'gp_mean': '#0F028A',
@@ -16,6 +17,8 @@
     'envelope_max_opacity': 0.8,
     'minor_tick_highlight': 'red',
     'acq_func_fill': 'lightblue',
+    'acq_func_intro': 'seagreen',
+    'acq_func_intro_fill': 'mediumaquamarine'
 })
 
 # Various parameters for plotting required by our own code
@@ -59,4 +62,17 @@ def f(x):
     "x": (2, 9),
     "gp_y": (-3, 3),
     "acq_y": (0, 5),
+    "x_intro": (2, 13),
+    "y_intro": (0, 15)
+
+}
+
+zorders = {
+    'annotations_low': 20,
+    'zone_of_imp': 30,
+    'annotations_normal': 40,
+    'datapoints': 41,
+    'incumbent': 42,
+    'annotations_high': 60,
+    'legend': 100
 }

w06_hpo_bo/scripts/bo_intro.py

@@ -1,42 +1,221 @@
-from sklearn.gaussian_process import GaussianProcessRegressor
-from sklearn.gaussian_process.kernels import RBF
-from bo_intro_utils import *
+import warnings
+warnings.filterwarnings('ignore')
+import argparse
+import logging
+from functools import partial
+
+import numpy as np
+from scipy.optimize import minimize
+from sklearn.gaussian_process import GaussianProcessRegressor as GPR
+from sklearn.gaussian_process.kernels import Matern
+
 from matplotlib import pyplot as plt
-plt.style.use(['ggplot', 'seaborn-talk'])
 
+import bo_intro_utils as boplot
+from bo_configurations import *
+
+
+SEED = None
+INIT_X_PRESENTATION = [4.5, 10]
+NUM_ACQ_OPTS = 50 # Number of times the acquisition function is optimized while looking for the next x to sample.
+TOGGLE_PRINT = True
+
+
+def initialize_dataset(initial_design, init=None):
+    """
+    Initialize some data to start fitting the GP on.
+    :param initial_design: Method for initializing the GP, choice between 'uniform', 'random', and 'presentation'
+    :param init: Number of datapoints to initialize with, if relevant
+    :return:
+    """
+
+    # sample initial query points
+    if initial_design == 'uniform':
+        x = np.linspace(bounds["x_intro"][0], bounds["x_intro"][1], init).reshape(-1, 1).tolist()
+    elif initial_design == 'random':
+        x = np.random.uniform(bounds["x_intro"][0], bounds["x_intro"][1], init).reshape(-1, 1).tolist()
+    elif initial_design == 'presentation':
+        x = np.array(INIT_X_PRESENTATION).reshape(-1, 1).tolist()
+
+    # get corresponding response values
+    y = list(map(f, x))
+
+    return x, y
+
+
+def run_bo(acquisition, max_iter, initial_design, acq_add, init=None):
+    """
+    BO
+    :param acquisition: type of acquisition function to be used
+    :param max_iter: max number of function calls
+    :param seed: seed used to keep experiments reproducible
+    :param initial_design: Method for initializing the GP, choice between 'uniform', 'random', and 'presentation'
+    :param acq_add: additional parameteres for acquisition function (e.g. kappa for LCB)
+    :param init: Number of datapoints to initialize GP with.
+    :return: all evaluated points.
+    """
+
+    logging.debug("Running BO with Acquisition Function {0}, maximum iterations {1}, initial design {2}, "
+                  "acq_add {3} and init {4}".format(acquisition, max_iter, initial_design, acq_add, init))
+    x, y = initialize_dataset(initial_design=initial_design, init=init)
+    logging.debug("Initialized dataset with:\nsamples {0}\nObservations {1}".format(x, y))
+
+    for i in range(1, max_iter):  # BO loop
+        logging.debug('Sample #%d' % (i))
+
+        # Fit GP to the currently available dataset
+        gp = GPR(kernel=Matern())
+        logging.debug("Fitting GP to\nx: {}\ny:{}".format(x, y))
+
+
+        #Add data to fit a more stable GP
+        x_add = np.add(x, 0.05)
+        x_fit = np.append(x, x_add).reshape(-1, 1)
+        y_fit = y + list(map(f, x_add))
+        gp.fit(x_fit, y_fit)  # fit the model
+
+        # ----------Plotting calls---------------
+        fig, ax1 = plt.subplots(1, 1, squeeze=True)
+        fig.tight_layout()
+        ax1.set_xlim(bounds["x_intro"])
+        ax1.set_ylim(bounds["y_intro"])
+        ax1.set_yticks([])
+        ax1.grid()
+
+        boplot.plot_objective_function(ax=ax1, translation=10)
+
+        annotate = False
+        if i == 1:
+             annotate = True
+
+        boplot.plot_gp(model=gp, confidence_intervals=[1.0, 2.0, 3.0], ax=ax1, custom_x=x, annotate=annotate,
+                       translation=10)
+
+        mark_incumbent = False
+
+        boplot.mark_observations(X_=x, Y_=y, ax=ax1, mark_incumbent=mark_incumbent, highlight_datapoint=len(y)-1)
+
+        # # Partially initialize the acquisition function to work with the fmin interface
+        # # (only the x parameter is not specified)
+        acqui = partial(acquisition, model=gp, eta=min(y), add=acq_add)
+
+        annotate = False
+        if i == 1:
+             annotate = True
+
+        boplot.plot_acquisition_function(acquisition, min(y), gp, acq_add, invert=True, ax=ax1, annotate=annotate,
+                                         scaling=30)
+
+        # optimize acquisition function, repeat 10 times, use best result
+        x_ = None
+        y_ = 10000
+        # Feel free to adjust the hyperparameters
+        for j in range(NUM_ACQ_OPTS):
+            opt_res = minimize(acqui, np.random.uniform(bounds["x_intro"][0], bounds["x_intro"][1]), method="L-BFGS-B", bounds=[(bounds["x_intro"][0], bounds["x_intro"][1])])
+            if opt_res.fun[0] < y_:
+                x_ = opt_res.x
+                y_ = opt_res.fun[0]
+
+
+        # Update dataset with new observation
+        x.append(x_)
+        y.append(f(x_))
+        logging.info("After {0}. loop iteration".format(i))
+        logging.info("x: {0:.3E}, y: {1:.3E}".format(x_[0], y_))
+
+        if i==1:
+             annotate_x = INIT_X_PRESENTATION[0]
+             ax1.annotate("Observation", xy=(annotate_x, f([annotate_x])+ 10), xytext=(annotate_x - 1, f([annotate_x]) + 14),
+                          arrowprops={'arrowstyle': 'fancy'}, zorder=10, fontsize='x-large')
+             annotate_x = INIT_X_PRESENTATION[1]
+             ax1.annotate("Objective function", xy=(annotate_x + 1, f([annotate_x + 1])+ 10), xytext=(annotate_x + 0.25, f([annotate_x + 1])+ 7),
+                          arrowprops={'arrowstyle': 'fancy'}, zorder=10, fontsize='x-large')
+        #
+        if i==2:
+             ax1.annotate("New observation", xy=(new_observation, f(new_observation)+ 10), xytext=(new_observation - 1, f(new_observation) + 6),
+                          arrowprops={'arrowstyle': 'fancy'}, zorder=19, fontsize='x-large')
+
+        ax1.set_xlabel(labels['xlabel'])
+
+        new_observation = x_
+
+        if TOGGLE_PRINT:
+            plt.savefig("plot_{}.pdf".format(i), dpi='figure',bbox_inches = 'tight')
+        else:
+            plt.show()
+        # ---------------------------------------
+
+    return y
 
 
-# Initialize Gaussian Process
-kernel = 2.0 * RBF(length_scale=1.0, length_scale_bounds=(1e-1, 10.0))
-gp = GaussianProcessRegressor(kernel=kernel)
+def main(num_evals, init_size, repetitions, initial_design, acq_add, acquisition):
+    for i in range(repetitions):
+        bo_res_1 = run_bo(max_iter=num_evals, init=init_size, initial_design=initial_design, acquisition=acquisition, acq_add=acq_add)
 
-# Generate data and fit GP
-noise = np.random.rand(4)
-data = np.linspace(0.15, 0.9, 4)[:, np.newaxis]
-y = np.sin(data[:, 0]*12)+1/2*np.sin(data[:, 0]*11)+1/2*np.sin(data[:, 0]*23)
-gp.fit(data, y)
 
-# Plot datapoints
-plt.scatter(data[:, 0], y, c='k', marker='X', s=100, zorder=9)
-plt.xlim(0, 1)
-plt.ylim(-6, 6)
-plt.xlabel("x")
-plt.ylabel("f(x)")
-plt.savefig("plot_datapoints.pdf", format='pdf')
-plt.show()
+if __name__ == '__main__':
+    cmdline_parser = argparse.ArgumentParser('AutoMLLecture')
 
-X_ = np.linspace(0, 1, 600)
-y_mean, y_cov = gp.predict(X_[:, np.newaxis], return_cov=True)
+    cmdline_parser.add_argument('-n', '--num_func_evals',
+                                default=10,
+                                help='Number of function evaluations',
+                                type=int)
+    cmdline_parser.add_argument('-f', '--init_db_size',
+                                default=2,
+                                help='Size of the initial database',
+                                type=int)
+    cmdline_parser.add_argument('-i', '--initial_design',
+                                default="presentation",
+                                choices=['random', 'uniform', 'presentation'],
+                                help='How to choose first observations.')
+    cmdline_parser.add_argument('-v', '--verbose',
+                                default=False,
+                                help='verbosity',
+                                action='store_true')
+    cmdline_parser.add_argument('-a', '--acquisition',
+                                default='EI',
+                                choices=['LCB', 'EI', 'PI'],
+                                help='acquisition function')
+    cmdline_parser.add_argument('-s', '--seed',
+                                default=15,
+                                help='Which seed to use',
+                                required=False,
+                                type=int)
+    cmdline_parser.add_argument('-r', '--repetitions',
+                                default=1,
+                                help='Number of repeations for the experiment',
+                                required=False,
+                                type=int)
+    cmdline_parser.add_argument('-p', '--print',
+                                default=True,
+                                help='Print graphs to file instead of displaying on screen.',
+                                action='store_true')
 
+    args, unknowns = cmdline_parser.parse_known_args()
+    log_lvl = logging.DEBUG if args.verbose else logging.INFO
+    logging.basicConfig(level=log_lvl)
 
-samples = [1, 3, 10, 100, 1000]
-seed = 13
+    if unknowns:
+        logging.warning('Found unknown arguments!')
+        logging.warning(str(unknowns))
+        logging.warning('These will be ignored')
 
-# Plot samples from GP
-for i, num_sample in enumerate(samples):
-    plot_sample_gp(num_samples=num_sample, data=data, y=y, X_domain=X_, gp=gp, rnd_state=seed)
+    # init_size = max(1, int(args.num_func_evals * args.fraction_init))
+    # Seed the RNG to obtain reproducible results
+    SEED = args.seed
+    np.random.seed(SEED)
 
-# Plot GP posterior
-plot_posterior_and_density(data=data, y=y, X_=X_, y_mean=y_mean, y_cov=y_cov, gp=gp, rnd_state=seed)
+    TOGGLE_PRINT = args.print
+    if TOGGLE_PRINT:
+        boplot.enable_printing(figsize=(30, 10))
+    else:
+        boplot.enable_onscreen_display()
 
 
+    main(   num_evals=args.num_func_evals,
+            init_size=args.init_db_size,
+            repetitions=args.repetitions,
+            initial_design=args.initial_design,
+            acquisition=acquisition_functions[args.acquisition],
+            acq_add=1
+            )