<jupyter_start><jupyter_text>This notebook generates plots for Gaussian Process Regression (GPR) for one dimension data. The code refers [Sklearn library](https://scikit-learn.org/stable/auto_examples/gaussian_process/plot_gpr_noisy_targets.htmlsphx-glr-auto-examples-gaussian-process-plot-gpr-noisy-targets-py)<jupyter_code>import numpy as np
from matplotlib import pyplot as plt
import matplotlib as mpl
from matplotlib import cm, rcParams
%matplotlib inline

from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import RBF, ConstantKernel as C

np.random.seed(1)

plot_path = '../../paper/plots/Main_Paper/'


def f(x):
    """The function to predict."""
    return np.sin(x) + 0.05 * (x -2.5)
    # return 2 * np.sin(x) # with C = 2
    #return - np.poly1d([1,2,3,4], True)(0.3 * (x+ 2)) 
    
def gpr_plot(X, y, x, rec_size = 1, design_dist = 5, title = 'GPR', legend_flag = True):
    """
    X: x train
    y: y train
    x: x test
    rec_size: ucb rec size
    design_dist: x[::design_dist] indicates design space
    title: for plot
    """

    # Instantiate a Gaussian Process model
    # kernel = C(1.0, (1e-3, 1e3)) * RBF(0.6, (1e-2, 1e2))
    # use fixed hyperparameters
    kernel = C(2) * RBF(1)
    gp = GaussianProcessRegressor(kernel=kernel, alpha=0.1,
                                  n_restarts_optimizer=0,#)
                                 optimizer=None)
    
    # Fit to data using Maximum Likelihood Estimation of the parameters
    gp.fit(X, y)
    print(gp.kernel_)

    # Make the prediction on the meshed x-axis (ask for MSE as well)
    y_pred, sigma = gp.predict(x, return_std=True)
    
    # recommendations
    rec_idx = np.argsort(y_pred[::design_dist] + 1.9600 * sigma[::design_dist])[-rec_size:]
    rec = x_test[::design_dist][rec_idx]
    
    # Plot the function, the prediction and the 95% confidence interval based on
    # the MSE
    mpl.style.use('seaborn')

    mpl.rcParams['axes.labelsize'] = 18
    mpl.rcParams['xtick.labelsize'] = 16
    mpl.rcParams['ytick.labelsize'] = 16
    mpl.rcParams['legend.fontsize'] = 16
    mpl.rcParams['axes.titlesize'] = 20

    fig, axes = plt.subplots(1, 1, figsize=(8, 6), sharex=False)
    ymin = -5
    ymax = 5
    axes.vlines(rec, ymin, ymax, color = 'r', linestyle = 'dashed', label = 'Recommendations')
    axes.plot(x, f(x), 'g:', label='True function')
    axes.plot(X, y, 'g.', markersize=10, label='Observations')

    axes.plot(x, y_pred, '-', label='Prediction')
    axes.fill(np.concatenate([x, x[::-1]]),
             np.concatenate([y_pred - 1.9600 * sigma,
                            (y_pred + 1.9600 * sigma)[::-1]]),
             alpha=.3, ec='None', label='Confidence interval')
    axes.plot(x, y_pred + 1.96 * sigma, color = 'orange')
    axes.scatter(x[::design_dist], (y_pred + 1.96 * sigma)[::design_dist], marker = '.', color = 'orange', label = 'ucb scores')
    axes.set_xlabel('$x$')
    axes.set_ylabel('$f(x)$')
    axes.set_ylim(ymin, ymax)
    axes.set_xlim(-2,10)
#     plt.legend(loc='upper left')
    plt.title(title)
    if legend_flag:
        lgd = plt.legend(bbox_to_anchor=(1.05, 1.0), loc='upper left')
        plt.savefig(plot_path + title.replace(' ', '_').replace(':','')+'.pdf', bbox_extra_artists=(lgd,), bbox_inches='tight')
    else:
        plt.savefig(plot_path + title.replace(' ', '_').replace(':','')+'.pdf', bbox_inches='tight')
    return y_pred, sigma, rec_idx
# x_test = np.atleast_2d(np.linspace(0, 10, 101)).T
x_test = np.atleast_2d(np.linspace(-2, 10, 121)).T
design_dist = 5
x_design = x_test.reshape(-1,)[::design_dist]
noise = np.random.normal(0, 0.05)

# init
# X = np.atleast_2d([-1, 1, 5., 6., 9]).T
X = np.atleast_2d([-1, 0.5, 5., 6., 9.5]).T
# Observations
y = f(X).ravel() + noise
x_test.reshape(-1,)[::5]<jupyter_output><empty_output><jupyter_text>In the following, we show the batch recommendations based on the upper confidence bound algorithm.We consider the batch size 2, with 5 initial observations.The design (recommendation) space is 24 uniformly distributed points in the range [-2,10], i.e. {-2, -1.5, -1, ..., 9.5, 10}.   One naive approach is to train GPR on sequences with known TIR labels and recommend sequences in design space with top $n = 2$ UCB scores.<jupyter_code>top2_pred, top2_sigma, top2_rec = gpr_plot(X, y,x_test, rec_size = 2, title = 'Recommend 2 Data Points: Top 2 UCB')<jupyter_output>1.41**2 * RBF(length_scale=1)<jupyter_text>The GP-BUCB approach is designed to reduce the similarity between the recommendations in a batch.  It recommends data point sequentially by assuming have seen the previous recommended points and updating the UCB score with the updated predicted SD.   In our toy example, GP-BUCB firstly recommend the data point with maximum UCB scores based on the predictions over initial 5 observations. As shown in the following first plot.  Then we add the recommended data point (x = 2) into training data with the predicted mean as label (note it is not the true label, i.e. observation), and recommend the second data point based on the new predictions.   As we can see, since we assume we have observed x = 2, then the new predicted variance of the data points in design space around x =2  decreases, so instead of recommending a similar data point x = 2.5, we recommend the data point which is in another area and potentially having high labels (x = 7.5).     By using BUCB, we rule out the sequences with high similarities in the same batch and increase the exploration efficiency.<jupyter_code># GPBUCB
design_dist = 5

# rec1
y_pred, sigma, rec = gpr_plot(X,y, x_test, title = 'GP-BUCB: Recommend 1st Data Point', legend_flag = False)

# rec2
obs_y = y_pred[::design_dist][rec] 
batch_new_y = np.append(y, obs_y)
batch_new_X = np.append(X.reshape(-1,), x_test[::design_dist][rec]).reshape(-1,1)

batch_new_pred, batch_new_sigma, batch_rec = gpr_plot(batch_new_X,batch_new_y,x_test, title = 'GP-BUCB: Recommend 2nd Data Point', legend_flag=False)<jupyter_output>1.41**2 * RBF(length_scale=1)
1.41**2 * RBF(length_scale=1)
{'fraction_non_alphanumeric': 0.13299358868979122, 'fraction_numerical': 0.028933092224231464, 'mean_word_length': 3.7199379363847944, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 11, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '3', 'dir': 'jupyter-structured-clean-dedup', 'id': '138814', 'license': 'Apache-2.0', 'n_tokens_mistral': 2258, 'n_tokens_neox': 2087, 'n_words': 614, 'path': 'synbio_rbs/demo/method_demo.ipynb', 'repo_name': 'mholowko/Solaris', 'size': '6083'}