<jupyter_start><jupyter_text>Bayesian optimization (BO) - find the minimum value of 1 variable functionReference: https://datachemeng.com/bayesianoptimization/ http://krasserm.github.io/2018/03/21/bayesian-optimization/ https://orbi.uliege.be/bitstream/2268/226433/1/PyData%202017_%20Bayesian%20optimization%20with%20Scikit-Optimize.pdf https://github.com/Ma-sa-ue/practice/blob/master/machine%20learning(python)/bayeisan_optimization.ipynb<jupyter_code>import numpy as np
import pandas as pd
import GPy
import GPyOpt
from numpy.random import seed
import matplotlib.pyplot as plt
from IPython.display import clear_output

%matplotlib inline
# Definition of a black-box function to understand the concept
def f_y(x) :
    y = -np.abs(x * np.sin(x * np.pi/60.0) + 20)
    return y
# Draw black-box function
x = np.arange(0, 180, 0.1)
y = f_y(x)
plt.plot(x, y, c = 'r')
plt.xlabel('x')
plt.ylabel('y')

plt.show()
# Make dataframe x, y
df_xy = pd.DataFrame()
df_xy['x'] = x
df_xy['y'] = f_y(x)
print(df_xy.shape)
df_xy.head()
# Check x, y at the minimum y position
min_y = np.argmin(np.array(df_xy['y']))
print('Minimum "y" index = ', min_y)
df_xy.loc[min_y]
# Random x generation
init_exp_num = 10

np.random.seed(5)
rnd_idx_initial = [np.random.randint(0, len(df_xy)) for i in range(init_exp_num)]
print(rnd_idx_initial)
# The initial dataframe of x, y 
df_xy_initial = df_xy.loc[rnd_idx_initial]
df_xy_initial
# Determine explanatory variables X and target variable y
X = df_xy_initial['x']
y = df_xy_initial['y']
# Make the initial dataset to train the first model even though sample size is very small
X_train_init = np.array(X)
X_test  = np.array(X)
y_train_init = np.array(y)
y_test  = np.array(y)

print('The initial experiment x, y')
print('x = ', X_train_init)
print('y = ', y_train_init)<jupyter_output>The initial experiment x, y
x =  [ 86.7  20.6 172.5  99.8 114.2 142.4   7.3 103.2  74.  176.7]
y =  [ -65.40897858  -38.15485111  -86.01289208  -66.94713835  -14.15045848
 -151.27331279  -22.72281081  -59.51696665  -29.51566487  -50.37987202]<jupyter_text>Make GPR models<jupyter_code># Dictionary of kernel
kernel_dict = {1:GPy.kern.RBF(1),
               2:GPy.kern.RBF(1) + GPy.kern.Bias(1) + GPy.kern.Linear(1),
               3:GPy.kern.PeriodicMatern32(1) * GPy.kern.RBF(1),
               4:GPy.kern.Exponential(1),
               5:GPy.kern.Exponential(1) * GPy.kern.RBF(1) + GPy.kern.White(1)}
# Select a suitable GPR kernel
kernel = kernel_dict[5]
def MakeModel(X, y, kern):
    model = GPy.models.GPRegression(X.reshape(-1, 1), y.reshape(-1, 1), kernel = kern)
    model.optimize()
    
    #print(model)
    model.plot()
    plt.show()
    
    return model<jupyter_output><empty_output><jupyter_text>GPR model for the initial experiment points<jupyter_code>loop_model = MakeModel(X_train_init, y_train_init, kernel)<jupyter_output><empty_output><jupyter_text>Adaptive DOE stage from here<jupyter_code># Bounds : the range of variables
bounds = [
    {'name': 'x',
     'type': 'continuous',
     'domain': (0.0, 180.0)}
]

# Iteration to find the next point
Bayesian_opt_iter = 1 
# The function that want to mimimize (or maximize by multiplied by -1)

def my_f(X):
    y = loop_model.predict(X)[0]   # return only the mean value
    return y
# Bayesian optimization for adaptive DOE

myProblem = GPyOpt.methods.BayesianOptimization(my_f, bounds)
myProblem.run_optimization(Bayesian_opt_iter)
    
X_opt = myProblem.x_opt
y_opt = f_y(X_opt)
X_new = np.append(X_train_init, X_opt)
y_new = np.append(y_train_init, y_opt)
# Interation for adaptive DOE

iter_BO = 40

for i in range(iter_BO - 1):
    #print('Sampling number = ', i)
    loop_model = MakeModel(X_new, y_new, kernel)
    myProblem = GPyOpt.methods.BayesianOptimization(my_f, bounds)
    myProblem.run_optimization(Bayesian_opt_iter)
    
    X_opt = myProblem.x_opt
    y_opt = f_y(X_opt)
    X_new = np.append(X_new, X_opt)
    y_new = np.append(y_new, y_opt)
    clear_output(True)

print('Sampling number = ', iter_BO)
loop_model = MakeModel(X_new, y_new, kernel)
myProblem = GPyOpt.methods.BayesianOptimization(my_f, bounds)
myProblem.run_optimization(Bayesian_opt_iter)

X_opt = myProblem.x_opt
y_opt = f_y(X_opt)
X_new = np.append(X_new, X_opt)
y_new = np.append(y_new, y_opt)
# Comparison with the hypothetical curve
x_o = np.arange(0, 180, 0.1)
y_o = f_y(x_o)
y_p  = loop_model.predict(x_o.reshape(-1, 1))[0][:,0]
y_sd = loop_model.predict(x_o.reshape(-1, 1))[1][:,0]

plt.plot(x_o, y_o, c = 'r', label = 'original')
plt.plot(x_o, y_p, c = 'g', label = 'prediction', linestyle = '--')

plt.fill_between(x_o, y_p+y_sd, y_p-y_sd, color = '#888888', alpha=0.4)

plt.scatter(df_xy_initial.iloc[:, 0], df_xy_initial.iloc[:, 1], c = 'b', marker = 'x', s = 100)
plt.scatter(X_new[3:], y_new[3:], marker = '^', c = 'k', s = 100)
plt.ylim(-200, 50)

plt.legend(loc = 'best')
plt.show()<jupyter_output><empty_output><jupyter_text>Additional experiments finished, so finally find the minimum value by BO<jupyter_code>np.random.seed(777)
optimizer = GPyOpt.methods.BayesianOptimization(
    f = my_f,
    domain = bounds,
    acquisition_type ='MPI',
    acquisition_par = 0.1,
    exact_eval=True
)

optimizer.run_optimization(max_iter = 30)
print('minimum x = ', optimizer.x_opt)
print('minimum y = ', optimizer.fx_opt)
optimizer.plot_acquisition()
optimizer.plot_convergence()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13840445269016696, 'fraction_numerical': 0.048794063079777364, 'mean_word_length': 4.473096446700508, 'pattern_counts': {'":': 0, '<': 20, '<?xml version=': 0, '>': 20, 'https://': 3, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '575166', 'license': 'MIT', 'n_tokens_mistral': 2360, 'n_tokens_neox': 2165, 'n_words': 454, 'path': 'Adaptive_DOE_1_variables.ipynb', 'repo_name': 'dropbyakiba/Bayesian-Optimization', 'size': '5390'}