MGP

.travis.yml

@@ -8,7 +8,7 @@ python:
 cache: pip
 install:
   - pip install -U pip wheel
-  - pip install tensorflow==1.0.1
+  - pip install tensorflow==1.3.0rc0
   - pip install --process-dependency-links .
   - pip install .[test]
   - pip install codecov

GPflowOpt/__init__.py

@@ -20,4 +20,6 @@
 from . import transforms
 from . import scaling
 from . import objective
+from . import models
 from . import pareto
+from . import models

GPflowOpt/acquisition/acquisition.py

@@ -14,8 +14,10 @@
 
 from ..scaling import DataScaler
 from ..domain import UnitCube
+from ..models import ModelWrapper
 
 from GPflow.param import Parameterized, AutoFlow, ParamList
+from GPflow.model import Model
 from GPflow import settings
 
 import numpy as np
@@ -48,7 +50,9 @@ def __init__(self, models=[], optimize_restarts=5):
         :param optimize_restarts: number of optimization restarts to use when training the models
         """
         super(Acquisition, self).__init__()
-        self._models = ParamList([DataScaler(m) for m in np.atleast_1d(models).tolist()])
+        models = np.atleast_1d(models)
+        assert all(isinstance(model, (Model, ModelWrapper))for model in models)
+        self._models = ParamList([DataScaler(m) for m in models])
         self._default_params = list(map(lambda m: m.get_free_state(), self._models))
 
         assert (optimize_restarts >= 0)

GPflowOpt/acquisition/ei.py

@@ -57,7 +57,6 @@ def __init__(self, model):
         :param model: GPflow model (single output) representing our belief of the objective
         """
         super(ExpectedImprovement, self).__init__(model)
-        assert (isinstance(model, Model))
         self.fmin = DataHolder(np.zeros(1))
         self.setup()
 

GPflowOpt/models.py

@@ -0,0 +1,167 @@
+# Copyright 2017 Joachim van der Herten, Nicolas Knudde
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+from .tf_wraps import rowwise_gradients
+
+from GPflow.param import Parameterized, AutoFlow
+from GPflow.model import Model, GPModel
+from GPflow.likelihoods import Gaussian
+from GPflow import settings
+
+import tensorflow as tf
+
+float_type = settings.dtypes.float_type
+
+
+class ModelWrapper(Parameterized):
+    """
+    Class for fast implementation of a wrapper for models defined in GPflow. Once wrapped, all lookups for attributes
+    which are not found in the wrapper class are automatically forwarded to the wrapped model.
+
+    To influence the I/O of methods on the wrapped class, simply implement the method in the wrapper and call the
+    appropriate methods on the wrapped class. Specific logic is included to make sure that if AutoFlow methods are
+    influenced following this pattern, the original AF storage (if existing) is unaffected and a new storage is added
+    to the subclass.
+    """
+    def __init__(self, model):
+        """
+        :param model: model to be wrapped
+        """
+        super(ModelWrapper, self).__init__()
+
+        assert isinstance(model, (Model, ModelWrapper))
+        #: Wrapped model
+        self.wrapped = model
+
+    def __getattr__(self, item):
+        """
+        If an attribute is not found in this class, it is searched in the wrapped model
+        """
+        # Exception for AF storages, if a method with the same name exists in this class, do not find the cache
+        # in the wrapped model.
+        if item.endswith('_AF_storage'):
+            method = item[1:].rstrip('_AF_storage')
+            if method in dir(self):
+                raise AttributeError("{0} has no attribute {1}".format(self.__class__.__name__, item))
+        return getattr(self.wrapped, item)
+
+    def __setattr__(self, key, value):
+        """
+        1) If setting :attr:`wrapped` attribute, point parent to this object (the datascaler).
+        2) If setting the recompilation attribute, always do this on the wrapped class.
+        """
+        if key is 'wrapped':
+            object.__setattr__(self, key, value)
+            value.__setattr__('_parent', self)
+            return
+
+        try:
+            # If attribute is in this object, set it. Test by using getattribute instead of hasattr to avoid lookup in
+            # wrapped object.
+            self.__getattribute__(key)
+            super(ModelWrapper, self).__setattr__(key, value)
+        except AttributeError:
+            # Attribute is not in wrapper.
+            # In case no wrapped object is set yet (e.g. constructor), set in wrapper.
+            if 'wrapped' not in self.__dict__:
+                super(ModelWrapper, self).__setattr__(key, value)
+                return
+
+            if hasattr(self, key):
+                # Now use hasattr, we know getattribute already failed so if it returns true, it must be in the wrapped
+                # object. Hasattr is called on self instead of self.wrapped to account for the different handling of
+                # AF storages.
+                # Prefer setting the attribute in the wrapped object if exists.
+                setattr(self.wrapped, key, value)
+            else:
+                #  If not, set in wrapper nonetheless.
+                super(ModelWrapper, self).__setattr__(key, value)
+
+    def __eq__(self, other):
+        return self.wrapped == other
+
+    def __str__(self, prepend=''):
+        return self.wrapped.__str__(prepend)
+
+
+class MGP(ModelWrapper):
+    """
+    Marginalisation of the hyperparameters during evaluation time using a Laplace Approximation
+    Key reference:
+
+    ::
+
+       @article{Garnett:2013,
+          title={Active learning of linear embeddings for Gaussian processes},
+          author={Garnett, Roman and Osborne, Michael A and Hennig, Philipp},
+          journal={arXiv preprint arXiv:1310.6740},
+          year={2013}
+        }
+    """
+
+    def __init__(self, model):
+        assert isinstance(model, GPModel), "Object has to be a GP model"
+        assert isinstance(model.likelihood, Gaussian), "Likelihood has to be Gaussian"
+        super(MGP, self).__init__(model)
+
+    def build_predict(self, fmean, fvar, theta):
+        h = tf.hessians(self.build_likelihood() + self.build_prior(), theta)[0]
+        L = tf.cholesky(-h)
+
+        N = tf.shape(fmean)[0]
+        D = tf.shape(fmean)[1]
+
+        fmeanf = tf.reshape(fmean, [N * D, 1])      # N*D x 1
+        fvarf = tf.reshape(fvar, [N * D, 1])        # N*D x 1
+
+        Dfmean = rowwise_gradients(fmeanf, theta)   # N*D x k
+        Dfvar = rowwise_gradients(fvarf, theta)     # N*D x k
+
+        tmp1 = tf.transpose(tf.matrix_triangular_solve(L, tf.transpose(Dfmean)))    # N*D x k
+        tmp2 = tf.transpose(tf.matrix_triangular_solve(L, tf.transpose(Dfvar)))     # N*D x k
+        return fmean, 4 / 3 * fvar + tf.reshape(tf.reduce_sum(tf.square(tmp1), axis=1), [N, D]) \
+               + 1 / 3 / (fvar + 1E-3) * tf.reshape(tf.reduce_sum(tf.square(tmp2), axis=1), [N, D])
+
+    @AutoFlow((float_type, [None, None]))
+    def predict_f(self, Xnew):
+        """
+        Compute the mean and variance of the latent function(s) at the points
+        Xnew.
+        """
+        theta = self._predict_f_AF_storage['free_vars']
+        fmean, fvar = self.wrapped.build_predict(Xnew)
+        return self.build_predict(fmean, fvar, theta)
+
+    @AutoFlow((float_type, [None, None]))
+    def predict_y(self, Xnew):
+        """
+        Compute the mean and variance of held-out data at the points Xnew
+        """
+        theta = self._predict_y_AF_storage['free_vars']
+        pred_f_mean, pred_f_var = self.wrapped.build_predict(Xnew)
+        fmean, fvar = self.wrapped.likelihood.predict_mean_and_var(pred_f_mean, pred_f_var)
+        return self.build_predict(fmean, fvar, theta)
+
+    @AutoFlow((float_type, [None, None]), (float_type, [None, None]))
+    def predict_density(self, Xnew, Ynew):
+        """
+        Compute the (log) density of the data Ynew at the points Xnew
+
+        Note that this computes the log density of the data individually,
+        ignoring correlations between them. The result is a matrix the same
+        shape as Ynew containing the log densities.
+        """
+        theta = self._predict_density_AF_storage['free_vars']
+        pred_f_mean, pred_f_var = self.wrapped.build_predict(Xnew)
+        pred_f_mean, pred_f_var = self.build_predict(pred_f_mean, pred_f_var, theta)
+        return self.likelihood.predict_density(pred_f_mean, pred_f_var, Ynew)

GPflowOpt/scaling.py

@@ -12,17 +12,18 @@
 # See the License for the specific language governing permissions and
 # limitations under the License.
 
-from GPflow.param import DataHolder, AutoFlow, Parameterized
-from GPflow.model import Model, GPModel
+from GPflow.param import DataHolder, AutoFlow
+from GPflow.model import GPModel
 from GPflow import settings
 import numpy as np
 from .transforms import LinearTransform, DataTransform
 from .domain import UnitCube
+from .models import ModelWrapper
 
 float_type = settings.dtypes.float_type
 
 
-class DataScaler(GPModel):
+class DataScaler(ModelWrapper):
     """
     Model-wrapping class, primarily intended to assure the data in GPflow models is scaled. One DataScaler wraps one
     GPflow model, and can scale the input as well as the output data. By default, if any kind of object attribute
@@ -59,13 +60,8 @@ def __init__(self, model, domain=None, normalize_Y=False):
         :param normalize_Y: (default: False) enable automatic scaling of output values to zero mean and unit
          variance.
         """
-        # model sanity checks
-        assert (model is not None)
-        assert (isinstance(model, GPModel))
-        self._parent = None
-
-        # Wrap model
-        self.wrapped = model
+        # model sanity checks, slightly stronger conditions than the wrapper
+        super(DataScaler, self).__init__(model)
 
         # Initial configuration of the datascaler
         n_inputs = model.X.shape[1]
@@ -74,34 +70,8 @@ def __init__(self, model, domain=None, normalize_Y=False):
         self._normalize_Y = normalize_Y
         self._output_transform = LinearTransform(np.ones(n_outputs), np.zeros(n_outputs))
 
-        # The assignments in the constructor of GPModel take care of initial re-scaling of model data.
-        super(DataScaler, self).__init__(model.X.value, model.Y.value, None, None, 1, name=model.name+"_datascaler")
-        del self.kern
-        del self.mean_function
-        del self.likelihood
-
-    def __getattr__(self, item):
-        """
-        If an attribute is not found in this class, it is searched in the wrapped model
-        """
-        return self.wrapped.__getattribute__(item)
-
-    def __setattr__(self, key, value):
-        """
-        If setting :attr:`wrapped` attribute, point parent to this object (the datascaler)
-        """
-        if key is 'wrapped':
-            object.__setattr__(self, key, value)
-            value.__setattr__('_parent', self)
-            return
-
-        super(DataScaler, self).__setattr__(key, value)
-
-    def __eq__(self, other):
-        return self.wrapped == other
-
-    def __str__(self, prepend=''):
-        return self.wrapped.__str__(prepend)
+        self.X = model.X.value
+        self.Y = model.Y.value
 
     @property
     def input_transform(self):
@@ -216,6 +186,20 @@ def build_predict(self, Xnew, full_cov=False):
         f, var = self.wrapped.build_predict(self.input_transform.build_forward(Xnew), full_cov=full_cov)
         return self.output_transform.build_backward(f), self.output_transform.build_backward_variance(var)
 
+    @AutoFlow((float_type, [None, None]))
+    def predict_f(self, Xnew):
+        """
+        Compute the mean and variance of held-out data at the points Xnew
+        """
+        return self.build_predict(Xnew)
+
+    @AutoFlow((float_type, [None, None]))
+    def predict_f_full_cov(self, Xnew):
+        """
+        Compute the mean and variance of held-out data at the points Xnew
+        """
+        return self.build_predict(Xnew, full_cov=True)
+
     @AutoFlow((float_type, [None, None]))
     def predict_y(self, Xnew):
         """
@@ -230,6 +214,6 @@ def predict_density(self, Xnew, Ynew):
         """
         Compute the (log) density of the data Ynew at the points Xnew
         """
-        mu, var = self.build_predict(Xnew)
+        mu, var = self.wrapped.build_predict(self.input_transform.build_forward(Xnew))
         Ys = self.output_transform.build_forward(Ynew)
         return self.likelihood.predict_density(mu, var, Ys)

GPflowOpt/tf_wraps.py

@@ -0,0 +1,42 @@
+# Copyright 2017 Joachim van der Herten, Nicolas Knudde
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+import tensorflow as tf
+from GPflow import settings
+
+float_type = settings.dtypes.float_type
+
+
+def rowwise_gradients(Y, X):
+    """
+    For a 2D Tensor Y, compute the derivative of each columns w.r.t  a 2D tensor X.
+
+    This is done with while_loop, because of a known incompatibility between map_fn and gradients.
+    """
+    num_rows = tf.shape(Y)[0]
+    num_feat = tf.shape(X)[0]
+
+    def body(old_grads, row):
+        g = tf.expand_dims(tf.gradients(Y[row], X)[0], axis=0)
+        new_grads = tf.concat([old_grads, g], axis=0)
+        return new_grads, row + 1
+
+    def cond(_, row):
+        return tf.less(row, num_rows)
+
+    shape_invariants = [tf.TensorShape([None, None]), tf.TensorShape([])]
+    grads, _ = tf.while_loop(cond, body, [tf.zeros([0, num_feat], float_type), tf.constant(0)],
+                             shape_invariants=shape_invariants)
+
+    return grads

doc/source/interfaces.rst

@@ -36,3 +36,11 @@ Transform
    :special-members:
 .. autoclass:: GPflowOpt.transforms.DataTransform
    :special-members:
+
+ModelWrapper
+------------
+.. automodule:: GPflowOpt.models
+   :special-members:
+.. autoclass:: GPflowOpt.models.ModelWrapper
+   :members:
+   :special-members: