GPflow - GP classification with 1-dim Linear kernel fits poorly for 2 dimension data

Question

Following the issue #1435, I have an additional question to how to use GPflow.

I replicate the issue in an additional kernel: https://github.com/avalonhse/BayesNotebook/blob/master/Issue_2_GPFlow_Linear_Classification.ipynb

My purpose is fitting an additive kernel to a 2-dimensional data (squared exponential in dimension 1 and linear kernel in dimension 2). Following the instruction of #1435, I have been successfully fitting the model with kernel gpflow.kernels.Linear(variance= 0.1).

Linear kernel

However, when I use the kernel gpflow.kernels.Linear(active_dims=1,variance= 0.01) as I original planned, the model is not fitted. I used the GPy with same kernel as a reference then the result looks reasonable.

1-dim GPFlow kernel

import numpy as np
X = np.array([[ 9.96578428, 60.],[ 9.96578428, 40.],[ 9.96578428, 20.],
       [10.96578428, 30.],[11.96578428, 40.],[12.96578428, 50.],
       [12.96578428, 70.],[8.96578428, 30. ],[ 7.96578428, 40.],
       [ 6.96578428, 50.],[ 6.96578428, 30.],[ 6.96578428, 10.],
       [11.4655664 , 71.],[ 8.56605404, 63.],[12.41574177, 69.],
       [10.61562964, 48.],[ 7.61470984, 51.],[ 9.31514956, 45.]])
Y = np.array([[1., 1., 0., 0., 0., 0., 0., 0., 0., 1., 1., 0., 1., 1., 0., 0., 1., 0.]]).T

# plotting
import matplotlib.pyplot as plt
import matplotlib
%matplotlib inline
def plot(X,Y):
    mask = Y[:, 0] == 1

    plt.figure(figsize=(6, 6))
    plt.plot(X[mask, 0], X[mask, 1], "oC0", mew=0, alpha=0.5)

    plt.ylim(-10, 100)
    plt.xlim(5, 15)

    _ = plt.plot(X[np.logical_not(mask), 0], X[np.logical_not(mask), 1], "oC1", mew=0, alpha=0.5)

plot(X,Y)

# Evaluate real function and the predicted probability
res = 500
xx, yy = np.meshgrid(np.linspace(5, 15, res),
                     np.linspace(- 10, 120, res))
Xplot =  np.vstack((xx.flatten(), yy.flatten())).T

# Code followed the Notebook : https://gpflow.readthedocs.io/en/develop/notebooks/basics/classification.html

import tensorflow as tf
import tensorflow_probability as tfp
import gpflow
from gpflow.utilities import print_summary, set_trainable, to_default_float
gpflow.config.set_default_summary_fmt("notebook")

def testGPFlow(k):
    m = gpflow.models.VGP(
        (X, Y), 
        kernel= k,
        likelihood=gpflow.likelihoods.Bernoulli()
    )
    print("\n ########### Model before optimzation ########### \n")
    print_summary(m)

    print("\n ########### Model after optimzation ########### \n")
    opt = gpflow.optimizers.Scipy()
    res = opt.minimize(
        m.training_loss, variables=m.trainable_variables, options=dict(maxiter=2500), method="L-BFGS-B"
    )        

    print(' Message: ' + str(res.message) + '\n Status = ' + str(res.status) + '\n Number of iterations = ' + str(res.nit))
    print_summary(m)

    means, _ = m.predict_y(Xplot)  # here we only care about the mean
    y_prob = means.numpy().reshape(*xx.shape)

    print("Fitting model using GPFlow")
    plot(X,Y)

    _ = plt.contour(
        xx,
        yy,
        y_prob,
        [0.5],  # plot the p=0.5 contour line only
        colors="k",
        linewidths=1.8,
        zorder=100,
    )
k = gpflow.kernels.Linear(active_dims=[1],variance= 0.01)
testGPFlow(k)

k = gpflow.kernels.Linear(variance= 1)
testGPFlow(k)

The GPy code is for reference only to suggest how a fitted model should be. I am aware that GPy and GPflow use different methods. My question is why GPflow model does not fit when I specify the Linear kernel in 1 dimension.

Answer 1

Thanks for posting this question, Hoang, and for using GPflow.

When you specify input_dim in Gpy, you are telling the algorithm to act on two dimensions. Active_dims in GPflow behaves differently. It specifies which dimensions you want the kernel to act on. 'active_dims = 1' is telling GPflow to apply your linear kernel to only the y dimension.

Since you want your kernel to act on both x and y dimensions, you should specify active_dims = [0,1] rather than just 'active_dims = 1.' When I run your code with this fix, I get a result identical to GPy's result: