fix convergence experiments to use sigmoidal invlink, and tune M

shgp/classification/experiments/convergence_experiment.py

@@ -4,7 +4,7 @@
 import tensorflow as tf
 
 from shgp.data.metadata_reinit import ReinitMetaDataset
-from shgp.data.metadata_convergence import PimaConvergenceMetaDataset
+from shgp.data.metadata_convergence import BananaConvergenceMetaDataset
 from shgp.inducing.initialisation_methods import h_reinitialise_PGPR, k_means, uniform_subsample
 from shgp.utilities.train_pgpr import train_pgpr
 from shgp.utilities.train_svgp import train_svgp
@@ -17,6 +17,11 @@
 the convergence of the ELBO using Adam vs BFGS on both k_means and hgv. This allows us to analyse 
 whether the inducing points of hgv are initialised closer than k_means to the gradient-optimised result.
 To evaluate other inducing point selection methods, change the model definitions below.
+
+Observations:
+ HGV is always better than K-means for both Adam & BFGS.
+ BFGS learns quickest initially, but Adam often overtakes at the end.
+ SVGP requires much larger M to converge to its optimal elbo than PGPR requires to converge to its optimal ELBO.
 """
 
 
@@ -26,8 +31,7 @@ def run_convergence_experiment(X, Y, M, opt_iters):
  ##############################################
 
  print("Training non-sparse model...")
- svgp, elbo_svgp = train_svgp(X, Y, M=len(X), train_iters=250, init_method=uniform_subsample)
- optimal_kernel = svgp.kernel
+ _, elbo_svgp = train_svgp(X, Y, M=len(X), train_iters=250, init_method=uniform_subsample)
  optimal = np.full(opt_iters + 1, elbo_svgp)
 
  print("Non-sparse result:")
@@ -53,7 +57,7 @@ def run_convergence_experiment(X, Y, M, opt_iters):
  ####################
 
  results_hgv_bfgs, results_hgv_adam, results_kmeans_bfgs, results_kmeans_adam = \
- test_convergence((X, Y), opt_iters, hgv_Z, kmeans_Z, optimal_kernel)
+ test_convergence((X, Y), opt_iters, hgv_Z, kmeans_Z)
 
  results = list(zip(results_hgv_bfgs, results_hgv_adam, results_kmeans_bfgs, results_kmeans_adam))
 
@@ -67,44 +71,43 @@ def run_convergence_experiment(X, Y, M, opt_iters):
 
 
 # Test the performance of a model with a given number of inducing points, M.
-def test_convergence(data, opt_iters, hgv_Z, kmeans_Z, optimal_kernel):
+def test_convergence(data, opt_iters, hgv_Z, kmeans_Z):
  print("Training HGV BFGS...")
  bfgs_opt = gpflow.optimizers.Scipy()
- results_hgv_bfgs = train_convergence_svgp(data, opt_iters, optimal_kernel, hgv_Z, bfgs_opt)
+ results_hgv_bfgs = train_convergence_svgp(data, opt_iters, hgv_Z, bfgs_opt)
  print("HGV BFGS trained: ELBO = {}".format(results_hgv_bfgs[-1]))
  print("Training HGV Adam...")
  adam_opt = tf.optimizers.Adam(beta_1=0.5, beta_2=0.5)
- results_hgv_adam = train_convergence_svgp(data, opt_iters, optimal_kernel, hgv_Z, adam_opt)
+ results_hgv_adam = train_convergence_svgp(data, opt_iters, hgv_Z, adam_opt)
  print("HGV Adam trained: ELBO = {}".format(results_hgv_adam[-1]))
 
  print("Training K-means BFGS...")
  bfgs_opt = gpflow.optimizers.Scipy()
- results_kmeans_bfgs = train_convergence_svgp(data, opt_iters, optimal_kernel, kmeans_Z, bfgs_opt)
+ results_kmeans_bfgs = train_convergence_svgp(data, opt_iters, kmeans_Z, bfgs_opt)
  print("K-means BFGS trained: ELBO = {}".format(results_kmeans_bfgs[-1]))
  print("Training K-means BFGS...")
  adam_opt = tf.optimizers.Adam(beta_1=0.5, beta_2=0.5)
- results_kmeans_adam = train_convergence_svgp(data, opt_iters, optimal_kernel, kmeans_Z, adam_opt)
+ results_kmeans_adam = train_convergence_svgp(data, opt_iters, kmeans_Z, adam_opt)
  print("K-means Adam trained: ELBO = {}".format(results_kmeans_adam[-1]))
 
  return results_hgv_bfgs, results_hgv_adam, results_kmeans_bfgs, results_kmeans_adam
 
 
-def train_convergence_svgp(data, opt_iters, kernel, Z, opt):
+def train_convergence_svgp(data, opt_iters, Z, opt):
  """
  Train an SVGP model until completion, recording the ELBO every 'elbo_period' iterations.
  If we error, keep retrying until success - this is due to a spurious Cholesky error.
  """
  model = gpflow.models.SVGP(
- kernel=kernel,
- likelihood=gpflow.likelihoods.Bernoulli(),
- inducing_variable=Z
+ kernel=gpflow.kernels.SquaredExponential(),
+ likelihood=gpflow.likelihoods.Bernoulli(tf.sigmoid),
+ inducing_variable=Z.copy()
  )
- gpflow.set_trainable(model.kernel, False) # fixed hyperparameters
- gpflow.set_trainable(model.inducing_variable, True) # optimised Z
- results = [model.elbo(data).numpy()]
+ gpflow.set_trainable(model.inducing_variable, True)
 
  # Try to run the full optimisation cycle.
  try:
+ results = [model.elbo(data).numpy()]
  if isinstance(opt, gpflow.optimizers.Scipy):
  # Optimise kernel hyperparameters, recording the ELBO after each step.
  for _ in range(opt_iters):
@@ -124,7 +127,7 @@ def train_convergence_svgp(data, opt_iters, kernel, Z, opt):
  raise error
  else:
  print("Cholesky error caught, retrying...")
- return train_convergence_svgp(data, opt_iters, kernel, Z, opt)
+ return train_convergence_svgp(data, opt_iters, Z, opt)
 
  return results
 
@@ -146,7 +149,7 @@ def plot_results(name, opt_iters, results, optimal):
  plt.xlabel('Iterations')
  # Axis limits
  plt.xlim(0, opt_iters)
- plt.xticks([0, 50, 100, 150])
+ plt.xticks(np.arange(0, opt_iters + 1, 50))
  plt.yscale('symlog')
 
  # Setup each subplot
@@ -163,7 +166,7 @@ def plot_results(name, opt_iters, results, optimal):
 
 if __name__ == '__main__':
  # Load data
- dataset = PimaConvergenceMetaDataset() # to test another dataset, just change this definition
+ dataset = BananaConvergenceMetaDataset() # to test another dataset, just change this definition
  X, Y = dataset.load_data()
 
  # Get convergence results

shgp/data/metadata_convergence.py

@@ -19,34 +19,34 @@ class ConvergenceMetaDataset:
 class BananaConvergenceMetaDataset(BananaDataset, ConvergenceMetaDataset):
  def __init__(self):
  BananaDataset.__init__(self)
- ConvergenceMetaDataset.__init__(self, 40, 150)
+ ConvergenceMetaDataset.__init__(self, 80, 250)
 
 
 class CrabsConvergenceMetaDataset(CrabsDataset, ConvergenceMetaDataset):
  def __init__(self):
  CrabsDataset.__init__(self)
- ConvergenceMetaDataset.__init__(self, 10, 200)
+ ConvergenceMetaDataset.__init__(self, 20, 1000)
 
 
 class HeartConvergenceMetaDataset(HeartDataset, ConvergenceMetaDataset):
  def __init__(self):
  HeartDataset.__init__(self)
- ConvergenceMetaDataset.__init__(self, 35, 150)
+ ConvergenceMetaDataset.__init__(self, 60, 100)
 
 
 class IonosphereConvergenceMetaDataset(IonosphereDataset, ConvergenceMetaDataset):
  def __init__(self):
  IonosphereDataset.__init__(self)
- ConvergenceMetaDataset.__init__(self, 80, 250)
+ ConvergenceMetaDataset.__init__(self, 150, 350)
 
 
 class BreastCancerConvergenceMetaDataset(BreastCancerDataset, ConvergenceMetaDataset):
  def __init__(self):
  BreastCancerDataset.__init__(self)
- ConvergenceMetaDataset.__init__(self, 50, 150)
+ ConvergenceMetaDataset.__init__(self, 150, 250)
 
 
 class PimaConvergenceMetaDataset(PimaDataset, ConvergenceMetaDataset):
  def __init__(self):
  PimaDataset.__init__(self)
- ConvergenceMetaDataset.__init__(self, 60, 150)
+ ConvergenceMetaDataset.__init__(self, 100, 150)