Bayesian Analysis of Posterior Predictions from GP model

Questions
13

I didn’t realize it but yes I am doing Bayesian Optimization. I ditched the idea of working with the GPy model and built this code instead, does it seem right? I was getting a mass matrix contains zeros on the diagonal error when I was running with a small test data set.

with pm.Model() as model:
    theta = pm.Uniform('theta', lower=np.array([#lowest input values of data set]), upper=np.array([#highest input values of data set]), shape=13)
    m = pm.gp.mean.Constant(labels.mean())
    cov_func = pm.gp.cov.ExpQuad(13, ls=theta)
    gp = pm.gp.Marginal(m, cov_func=cov_func)
    sig_noise = pm.HalfNormal('sig_noise')
    y_ = gp.marginal_likelihood("y_", X=features, y=labels, noise=sig_noise, is_observed=False, observed=yobs)

Also, when I increase the data set to the actual data set, the time is extremely large. Is this normal?