Unconfused Optimum¶

import pandas as pd 
import expAscribe as ea
priori = pd.read_csv('./data/Wetlab/partial_links.csv',header=None)
ea.array_to_gml(priori.to_numpy(),'./data/Wetlab/influence_path.gml',["binding energy","pocket volume","bottleneck radius","length", "enzymatic activity"])

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from bayes_opt import BayesianOptimization, UtilityFunction
from matplotlib import cm
from PIL import Image
import io
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import RBF, ConstantKernel as C

'''
data_new = pd.DataFrame({
    'x': np.random.uniform(-2, 10, 100),
    'y': np.random.uniform(-2, 10, 100),
    'z': np.sin(np.random.uniform(-2, 10, 100)) + np.cos(np.random.uniform(-2, 10, 100))
})
'''

X_train = data_new[['x', 'y']].values
y_train = data_new['z'].values


kernel = C(1.5, (1e-2, 1e2)) * RBF(0.5, (1e-2, 1e2))
gp = GaussianProcessRegressor(kernel=kernel, n_restarts_optimizer=10)
gp.fit(X_train, y_train)

def target_function(x, y):
    X_test = np.array([[x, y]])
    return gp.predict(X_test)[0]

x = np.linspace(0, 2000, 100)
y = np.linspace(-10, 25, 100)
X, Y = np.meshgrid(x, y)
grid = np.c_[X.ravel(), Y.ravel()]  

Z = np.array([target_function(xi, yi) for xi, yi in grid])
Z = Z.reshape(X.shape)  

optimizer = BayesianOptimization(
    f=target_function,
    pbounds={'x': (0, 2000), 'y': (0, 20)},
    random_state=42,
    verbose=1
)
acq_function = UtilityFunction(kind="ei", kappa=5)

optimizer.register({'x': 1000, 'y': 6}, target_function(1000, 6))
optimizer.register({'x': 1200, 'y': 5}, target_function(1200, 5))

def posterior(optimizer, x_obs, y_obs, grid):
    optimizer._gp.fit(x_obs, y_obs)
    mu, sigma = optimizer._gp.predict(grid, return_std=True)
    return mu, sigma

def plot_gp(optimizer, X, Y, Z, step=0):
    fig, axs = plt.subplots(2, 2, figsize=(18, 14))
    grid = np.c_[X.ravel(), Y.ravel()]
    
    x_obs = np.array([[res["params"]["x"], res["params"]["y"]] for res in optimizer.res])
    y_obs = np.array([res["target"] for res in optimizer.res])
    
    mu, sigma = posterior(optimizer, x_obs, y_obs, grid)
    mu = mu.reshape(X.shape)
    sigma = sigma.reshape(X.shape)

    c1 = axs[0, 0].contourf(X, Y, Z, cmap='rainbow')
    axs[0, 0].scatter(x_obs[:, 0], x_obs[:, 1], color='grey', label='Observations')
    axs[0, 0].set_title("Target Function")
    fig.colorbar(c1, ax=axs[0, 0])

    c2 = axs[0, 1].contourf(X, Y, mu, cmap='coolwarm')
    axs[0, 1].scatter(x_obs[:, 0], x_obs[:, 1], color='grey', label='Observations')
    axs[0, 1].set_title(f"Gaussian Process Predicted Mean After {step} Steps")
    fig.colorbar(c2, ax=axs[0, 1])

    c3 = axs[1, 0].contourf(X, Y, sigma, cmap='coolwarm')
    axs[1, 0].scatter(x_obs[:, 0], x_obs[:, 1], color='grey', label='Observations')
    axs[1, 0].set_title(f"Gaussian Process Variance After {step} Steps")
    fig.colorbar(c3, ax=axs[1, 0])

    utility_function = UtilityFunction(kind="ei", kappa=5, xi=0)
    utility = utility_function.utility(grid, optimizer._gp, 0).reshape(X.shape)
    c4 = axs[1, 1].contourf(X, Y, utility, cmap='inferno')
    axs[1, 1].scatter(x_obs[:, 0], x_obs[:, 1], color='grey', label='Observations')
    axs[1, 1].set_title(f"Acquisition Function (EI) After {step} Steps")
    fig.colorbar(c4, ax=axs[1, 1])

    plt.tight_layout()
    buf = io.BytesIO()
    plt.savefig(buf, format='png')
    plt.close(fig)
    buf.seek(0)
    
    return buf

images = []
optimizer.maximize(init_points=5, n_iter=25)

for step in range(1, 25):
    optimizer.maximize(init_points=0, n_iter=1, acquisition_function=acq_function)
    buf = plot_gp(optimizer, X, Y, Z, step)
    image = Image.open(buf)
    images.append(image)

gif_path = 'optimization_process.gif'
images[0].save(
    gif_path, save_all=True, append_images=images[1:], duration=500, loop=0
)