Using Gaussian Process model to make inference?

Using PyMC5, I know:

  1. How to perform Gaussian Process regression.
  2. How to construct priors and likelihood, then make sampling to get posterior.

But I’m not sure, in PyMC5, how can I use the GP regression model to construct the likelihood. To illustrate, I shall use the example in Introductory Overview of PyMC — PyMC dev documentation.

with pm.Model() as inference_model:
    alpha = pm.Normal("alpha", mu=0, sigma=10)
    beta = pm.Normal("beta", mu=0, sigma=10, shape=2)
    sigma = pm.HalfNormal("sigma", sigma=1)

    mu = alpha + beta[0] * X1 + beta[1] * X2 # REPLACE WITH GP MODEL???

    Y_obs = pm.Normal("Y_obs", mu=mu, sigma=sigma, observed=Y)
    idata = pm.sample()

In this example, mu is calculated explicitly from an analytic expression. But what if we must use the Gaussian Process model to calculate mu? \mu = \mathcal{GP}(X)

The GP may be modeled as following:

with pm.Model() as model_gp:
    ls = pm.Gamma('ls', 2, 1, shape=3)
    cov =, ls=ls)
    gp =
    f_pred = gp.marginal_likelihood('f_pred', X=X_train, y=y_train, sigma=sigma)
    gp_data = pm.sample()

I had posted another relevant question at `predictt` in PyMC5?, explaining that I was able to do it in PyMC3, but not in PyMC5. But there was no reply. So I try to ask it again in a more general way here. I apologize if this is undesirable.

For the marginal likelihood GP, the likelihood is always Gaussian, and is wrapped up in the marginal_likelihood function. I’m not exactly sure what the question is here. Are you looking to generate samples from the GP? If so, you’d want to look at using the conditional method after fitting the model. For example:

with gp_model:
  f_pred = gp.conditional('f_pred', X_New)
  pm.sample_posterior_predictive(trace, extend_inferencedata=True)
1 Like