# MAP of Gaussian Process Classification in Tensorflow Probability

I'm attempting to implement Gaussian Process Classification learning in tensorflow-probability, but my estimator turns out to be very biased toward zero. As opposed to sklearn, I attempted to optimize the posterior log likelihood directly to find the maximum a posteriori estimator.

Is there anything I'm doing wrong in the process? Or is it in fact impossible to learn Gaussian Process Classification directly from the true posterior log likelihood?

In order to test my setup, I generate a single function $$f$$ from a GP distribution at given independent points $$X$$, and generate independent observation at each point according to $$P(Y_i = 1) = \sigma(f(X_i))$$ where $$\sigma$$ is the logistic function $$\frac{1}{1 + \exp(-z)}$$.

The following graph shows how much far from the truth my estimator tends to be:

Generating the data as follows:

import numpy as np
import tensorflow as tf
import tensorflow_probability as tfp

# GP kernel
# Independent variable points
X = np.random.uniform(-1,1, (30,1)).astype("float32")

# Generating a latent function at observable points
gp = tfp.distributions.GaussianProcess(kernel=kernel, index_points=X)
f = gp.sample(1)

with tf.Session() as sess:
f_ = sess.run(f)

# Generating observed classes (+-1)
Y = 2*tfp.distributions.Bernoulli(logits=tf.reshape(f_, (-1,))).sample(1) - 1

with tf.Session() as sess:
Y_ = sess.run(Y)


Learning the latent space using scikit-learn whice provides an estimator using the Laplace approximation

from sklearn.gaussian_process import GaussianProcessClassifier

clf = GaussianProcessClassifier()
clf.fit(X, Y_.reshape(-1))
sigma_f_sklearn = clf.predict_proba(X)[:,1]


Optimizing over the posterior distribution directly using tensorflow

def generate_posterior_log_likelihood(locations, observations):
"""
Build posterior log likelihood function with GP prior at specific locations
"""
gp_infer = tfp.distributions.GaussianProcess(
kernel=kernel, index_points=locations
)
def posterior_log_likelihood(f):
return (
gp_infer.log_prob(f) +
tf.reduce_sum(
tf.log_sigmoid(
tf.multiply(tf.cast(observations, tf.float32), f)
)
)
)
return posterior_log_likelihood

f = tf.get_variable("latent_f", shape=Y.shape,
initializer=tf.initializers.random_normal())
neg_log_likelihood = -generate_posterior_log_likelihood(X, Y_)(f)
.minimize(neg_log_likelihood)

with tf.Session() as sess:
sess.run(tf.global_variables_initializer())
prev_neg_log_likelihood = np.inf
# Optimizing until convergence
_, neg_log_likelihood_ = sess.run([optimize, neg_log_likelihood])
while prev_neg_log_likelihood - neg_log_likelihood_ > 1e-12:
prev_neg_log_likelihood = neg_log_likelihood_
_, neg_log_likelihood_ = sess.run([optimize, neg_log_likelihood])

# Calculating the points
f_map = sess.run(f)



I think your Bernoulli log prob expression is wrong (though I haven't thought through the correct expression for a {-1, 1}-valued Bernoulli). To be safe, you could try
# go back to {0, 1}