diff --git a/basic_gp_example.py b/basic_gp_example.py index 2f53fca..f82a59e 100644 --- a/basic_gp_example.py +++ b/basic_gp_example.py @@ -1,5 +1,6 @@ """ This code will run an offline example as-is without any input arguments. Simply run the file in python. +Use the OptOnFly=True flag to optimise on the fly In order to use it on a specific machine: 1. Change the *importlib.import_module('machine_interfaces.machine_interface_example')* to your machine interface name. @@ -9,6 +10,8 @@ """ import numpy as np +import os +import datetime import importlib import time import matplotlib.pyplot as plt @@ -19,8 +22,8 @@ from modules.bayes_optimization import BayesOpt, negUCB, negExpImprove from modules.OnlineGP import OGP -saveResultsQ = False - +saveResultsQ = True +OptOnFly = False # use this to switch between using the hyperparameters from the params file, and optimising these on the fly mi_module = importlib.import_module('machine_interfaces.machine_interface_example') @@ -36,6 +39,7 @@ start_point = 0*scan_params['start_point'] #if start_point is set to None, the optimizer will start from the current device settings. mi = mi_module.machine_interface(dev_ids = dev_ids, start_point = start_point) #an isotropic n-dimensional gaussian with amplitude=1, centered at the origin, plus gaussian background noise with std dev = 0.1 + #create the gp ndim = len(dev_ids) # GP parameters @@ -45,44 +49,48 @@ hyperparams = {'precisionMatrix': gp_precisionmat, 'amplitude_covar': gp_amp, 'noise_variance': gp_noise_variance} gp = OGP(ndim, hyperparams) -#create the bayesian optimizer that will use the gp as the model to optimize the machine -opt = BayesOpt(gp, mi, acq_func="UCB", start_dev_vals = mi.x, dev_ids = dev_ids) -opt.ucb_params = scan_params['ucb_params'] #set the acquisition function parameters -print('ucb_params',opt.ucb_params) +#create the bayesian optimizer that will use the gp as the model to optimize the machine + +if OptOnFly: + # optimize_kernel_on_the_fly is the iteration number to start optimize the kernel's hyperparmaters. If None, no optimization of the hypers during BO. + optimize_kernel_on_the_fly = 2 + opt = BayesOpt(gp, mi, acq_func="UCB", start_dev_vals=mi.x, dev_ids=dev_ids, + optimize_kernel_on_the_fly=optimize_kernel_on_the_fly) + # opt.ucb_params = scan_params['ucb_params'] # set the acquisition function parameters +else: + optimize_kernel_on_the_fly = None + opt = BayesOpt(gp, mi, acq_func="UCB", start_dev_vals = mi.x, dev_ids = dev_ids) + opt.ucb_params = scan_params['ucb_params'] #set the acquisition function parameters #run the gp search for some number of steps Obj_state_s=[] -optimize_kernel_on_the_fly = None #optimize_kernel_on_the_fly is the iteration number to start optimize the kernel's hyperparmaters. If None, no optimization of the hypers during BO. - -Niter = 10 +Niter = 100 for i in range(Niter): clear_output(wait=True) print ('iteration =', i) print ('current position:', mi.x, 'current objective value:', mi.getState()[1]) Obj_state_s.append(mi.getState()[1][0]) - - f = plt.figure(figsize=(20,3)) - ax = f.add_subplot(121) - ax2 = f.add_subplot(122) - ax.set_ylabel('Input controls',fontsize=12) - ax.set_xlabel('Iteration',fontsize=12) - for x, label in zip(opt.X_obs.T, opt.dev_ids): - ax.plot(x,'.-',label = label) - ax.legend() - ax2.set_ylabel('Objective',fontsize=12) - ax2.set_xlabel('Iteration',fontsize=12) - ax2.plot(Obj_state_s,'.-') - plt.show(); - if optimize_kernel_on_the_fly is not None: if i > optimize_kernel_on_the_fly: - opt.optimize_kernel_hyperparameters() - + opt.optimize_kernel_hyperparameters() opt.OptIter() - time.sleep(acquisition_delay) - + time.sleep(acquisition_delay) + +f = plt.figure(figsize=(20, 3)) +ax = f.add_subplot(121) +ax2 = f.add_subplot(122) +ax.set_ylabel('Input controls', fontsize=12) +ax.set_xlabel('Iteration', fontsize=12) +for x, label in zip(opt.X_obs.T, opt.dev_ids): + ax.plot(x, '.-', label=label) +ax.legend() +ax2.set_ylabel('Objective', fontsize=12) +ax2.set_xlabel('Iteration', fontsize=12) +ax2.plot(Obj_state_s, '.-') +plt.show() + #save results if desired if saveResultsQ == True: timestr = datetime.now().strftime('%Y-%m-%d-%H%M%S') diff --git a/modules/bayes_optimization.py b/modules/bayes_optimization.py index 1b76664..fa54186 100644 --- a/modules/bayes_optimization.py +++ b/modules/bayes_optimization.py @@ -583,6 +583,6 @@ def negUCB(x_new, model, ndim, nsteps, nu = 1., delta = 1.): tau = 2.*np.log(nsteps**(0.5*ndim+2.)*(np.pi**2.)/3./delta) GPUCB = y_mean + np.sqrt(nu * tau * y_var) - return -GPUCB + return -GPUCB[0] diff --git a/modules/parallelstuff.py b/modules/parallelstuff.py index c0570b3..0c8bc14 100644 --- a/modules/parallelstuff.py +++ b/modules/parallelstuff.py @@ -7,6 +7,7 @@ import numpy as np import multiprocessing as mp import copy +import time # handle 'IOError: [Errno 4] Interrupted system call' errors from multiprocessing.Queue.get #https://stackoverflow.com/questions/14136195/what-is-the-proper-way-to-handle-in-python-ioerror-errno-4-interrupted-syst @@ -34,7 +35,7 @@ def mworker(f,x0,fargs,margs,out_q): #print 'worker: margs = ',margs res = minimize(f, x0, args = fargs, **margs) #return [res.x, res.fun] - out_q.put([[res.x, res.fun[0][0]]]) + out_q.put([[res.x, res.fun[0]]]) # parallelize minimizations using different starting positions using multiprocessing, scipy.optimize.minimize def parallelminimize(f,x0s,fargs,margs,v0best=None,relative_bounds=None): @@ -337,7 +338,7 @@ def parallelmap2(f,fargslist,hostlist=None): def eworker(f,x,fargs,out_q): # worker invoked in a process puts the results in the output queue out_q res = f(x, *fargs) - out_q.put(np.hstack((x, res[0][0]))) + out_q.put(np.hstack((x, res[0]))) # eval function over a range of initial points neval and return the nkeep lowest function evals def parallelgridsearch(f,x0,lengths,fargs,neval,nkeep):