From d3721b76a8ec4f98932474834ca9add20e7f04e8 Mon Sep 17 00:00:00 2001 From: Max Zwiessele Date: Fri, 25 Oct 2013 15:29:04 +0100 Subject: [PATCH] changing all parameterized objects to be compatible with the new parameterization --- GPy/core/gp.py | 34 ++-- GPy/core/gp_base.py | 9 +- GPy/core/index_operations.py | 1 - GPy/core/model.py | 46 ++--- GPy/core/parameter.py | 212 +++++++++++++---------- GPy/core/parameterized.py | 282 +++++++++++++++++++++---------- GPy/examples/regression.py | 9 +- GPy/kern/kern.py | 281 ++++++++++++++---------------- GPy/kern/parts/bias.py | 10 +- GPy/kern/parts/kernpart.py | 15 +- GPy/kern/parts/linear.py | 38 ++--- GPy/kern/parts/prod.py | 50 +++--- GPy/kern/parts/rbf.py | 62 ++++--- GPy/kern/parts/rbfcos.py | 53 +++--- GPy/kern/parts/spline.py | 23 +-- GPy/kern/parts/symmetric.py | 15 +- GPy/kern/parts/white.py | 8 +- GPy/likelihoods/gaussian.py | 10 +- GPy/likelihoods/likelihood.py | 4 +- GPy/testing/kernel_tests.py | 2 +- GPy/testing/parameter_testing.py | 10 +- 21 files changed, 645 insertions(+), 529 deletions(-) diff --git a/GPy/core/gp.py b/GPy/core/gp.py index 166cac29..31210aa1 100644 --- a/GPy/core/gp.py +++ b/GPy/core/gp.py @@ -23,6 +23,7 @@ class GP(GPBase): """ def __init__(self, X, likelihood, kernel, normalize_X=False): GPBase.__init__(self, X, likelihood, kernel, normalize_X=normalize_X) + #self._set_params(self._get_params()) def getstate(self): @@ -89,21 +90,26 @@ class GP(GPBase): return (-0.5 * self.num_data * self.output_dim * np.log(2.*np.pi) - 0.5 * self.output_dim * self.K_logdet + self._model_fit_term() + self.likelihood.Z) +# def _log_likelihood_gradients(self): +# """ +# The gradient of all parameters. +# +# Note, we use the chain rule: dL_dtheta = dL_dK * d_K_dtheta +# """ +# #return np.hstack((self.kern.dK_dtheta(dL_dK=self.dL_dK, X=self.X), self.likelihood._gradients(partial=np.diag(self.dL_dK)))) +# +# if not isinstance(self.likelihood,EP): +# tmp = np.hstack((self.kern.dK_dtheta(dL_dK=self.dL_dK, X=self.X), self.likelihood._gradients(partial=np.diag(self.dL_dK)))) +# else: +# tmp = np.hstack((self.kern.dK_dtheta(dL_dK=self.dL_dK, X=self.X), self.likelihood._gradients(partial=np.diag(self.dL_dK)))) +# return tmp + + def dL_dtheta(self): + return self.kern.dK_dtheta(self.dL_dK, self.X) - def _log_likelihood_gradients(self): - """ - The gradient of all parameters. - - Note, we use the chain rule: dL_dtheta = dL_dK * d_K_dtheta - """ - #return np.hstack((self.kern.dK_dtheta(dL_dK=self.dL_dK, X=self.X), self.likelihood._gradients(partial=np.diag(self.dL_dK)))) - - if not isinstance(self.likelihood,EP): - tmp = np.hstack((self.kern.dK_dtheta(dL_dK=self.dL_dK, X=self.X), self.likelihood._gradients(partial=np.diag(self.dL_dK)))) - else: - tmp = np.hstack((self.kern.dK_dtheta(dL_dK=self.dL_dK, X=self.X), self.likelihood._gradients(partial=np.diag(self.dL_dK)))) - return tmp - + def dL_dlikelihood(self): + return self.likelihood._gradients(partial=np.diag(self.dL_dK)) + def _raw_predict(self, _Xnew, which_parts='all', full_cov=False, stop=False): """ Internal helper function for making predictions, does not account diff --git a/GPy/core/gp_base.py b/GPy/core/gp_base.py index bb892874..41c4a001 100644 --- a/GPy/core/gp_base.py +++ b/GPy/core/gp_base.py @@ -31,9 +31,9 @@ class GPBase(Model): self._Xoffset = np.zeros((1, self.input_dim)) self._Xscale = np.ones((1, self.input_dim)) - self.set_as_parameters(*self.kern._parameters_) - self.set_as_parameters(*self.likelihood._parameters_) - + self.add_parameter(self.kern, gradient=self.dL_dtheta) + self.add_parameter(self.likelihood, gradient=self.dL_dlikelihood) + # Model.__init__(self) # All leaf nodes should call self._set_params(self._get_params()) at # the end @@ -53,7 +53,8 @@ class GPBase(Model): self.likelihood, self.output_dim, self._Xoffset, - self._Xscale] + self._Xscale, + ] def setstate(self, state): self._Xscale = state.pop() diff --git a/GPy/core/index_operations.py b/GPy/core/index_operations.py index 6c342c6e..2e1fc774 100644 --- a/GPy/core/index_operations.py +++ b/GPy/core/index_operations.py @@ -34,7 +34,6 @@ class ParamDict(defaultdict): for a in self.iterkeys(): if numpy.all(a==key) and a._parent_index_==key._parent_index_: return super(ParamDict, self).__setitem__(a, value) - raise KeyError, key defaultdict.__setitem__(self, key, value) diff --git a/GPy/core/model.py b/GPy/core/model.py index 3a26fba5..0e136920 100644 --- a/GPy/core/model.py +++ b/GPy/core/model.py @@ -12,6 +12,7 @@ import numpy as np from domains import _POSITIVE, _REAL from numpy.linalg.linalg import LinAlgError from index_operations import ParameterIndexOperations +import itertools # import numdifftools as ndt class Model(Parameterized): @@ -28,8 +29,15 @@ class Model(Parameterized): def log_likelihood(self): raise NotImplementedError, "this needs to be implemented to use the model class" def _log_likelihood_gradients(self): + #def dK_d(self, param, dL_dK, X, X2) + g = np.zeros(self.size) + try: + [g.__setitem__(s, self.gradient_mapping[p]().flat) for p, s in itertools.izip(self._parameters_, self._param_slices_) if not p.is_fixed] + except KeyError: + raise KeyError, 'Gradient for {} not defined, please specify gradients for parameters to optimize'.format(p.name) + return g raise NotImplementedError, "this needs to be implemented to use the model class" - + def getstate(self): """ Get the current state of the class. @@ -153,19 +161,19 @@ class Model(Parameterized): return ret - def _transform_gradients(self, g): - x = self._get_params() - for constraint, index in self._constraints_.iteritems(): - if constraint != __fixed__: - g[index] = g[index] * constraint.gradfactor(x[index]) - #[np.put(g, i, v) for i, v in [(t[0], np.sum(g[t])) for t in self.tied_indices]] - [np.put(g, i, v) for i, v in [[i, t.sum()] for p in self._parameters_ for t,i in p._tied_to_me_.iteritems()]] -# if len(self.tied_indices) or len(self.fixed_indices): -# to_remove = np.hstack((self.fixed_indices + [t[1:] for t in self.tied_indices])) -# return np.delete(g, to_remove) -# else: - if self._fixes_ is not None: return g[self._fixes_] - return g +# def _transform_gradients(self, g): +# x = self._get_params() +# for constraint, index in self.constraints.iteritems(): +# if constraint != __fixed__: +# g[index] = g[index] * constraint.gradfactor(x[index]) +# #[np.put(g, i, v) for i, v in [(t[0], np.sum(g[t])) for t in self.tied_indices]] +# [np.put(g, i, v) for i, v in [[i, t.sum()] for p in self.flattened_parameters for t,i in p._tied_to_me_.iteritems()]] +# # if len(self.tied_indices) or len(self.fixed_indices): +# # to_remove = np.hstack((self.fixed_indices + [t[1:] for t in self.tied_indices])) +# # return np.delete(g, to_remove) +# # else: +# if self._fixes_ is not None: return g[self._fixes_] +# return g def randomize(self): """ @@ -173,15 +181,15 @@ class Model(Parameterized): Make this draw from the prior if one exists, else draw from N(0,1) """ # first take care of all parameters (from N(0,1)) - x = self._get_params_transformed() - x = np.random.randn(x.size) + #x = self._get_params_transformed() + x = np.random.randn(self.size_transformed) self._set_params_transformed(x) # now draw from prior where possible x = self._get_params() if self.priors is not None: [np.put(x, i, p.rvs(1)) for i, p in enumerate(self.priors) if not p is None] self._set_params(x) - self._set_params_transformed(self._get_params_transformed()) # makes sure all of the tied parameters get the same init (since there's only one prior object...) + #self._set_params_transformed(self._get_params_transformed()) # makes sure all of the tied parameters get the same init (since there's only one prior object...) def optimize_restarts(self, num_restarts=10, robust=False, verbose=True, parallel=False, num_processes=None, **kwargs): @@ -480,7 +488,7 @@ class Model(Parameterized): names = self._get_param_names_transformed() except NotImplementedError: names = ['Variable %i' % i for i in range(len(x))] - + import ipdb;ipdb.set_trace() # Prepare for pretty-printing header = ['Name', 'Ratio', 'Difference', 'Analytical', 'Numerical'] max_names = max([len(names[i]) for i in range(len(names))] + [len(header[0])]) @@ -522,7 +530,7 @@ class Model(Parameterized): d = '%.6f' % float(difference) g = '%.6f' % gradient ng = '%.6f' % float(numerical_gradient) - grad_string = "{0:^{c0}}|{1:^{c1}}|{2:^{c2}}|{3:^{c3}}|{4:^{c4}}".format(formatted_name, r, d, g, ng, c0=cols[0] + 9, c1=cols[1], c2=cols[2], c3=cols[3], c4=cols[4]) + grad_string = "{0:<{c0}}|{1:^{c1}}|{2:^{c2}}|{3:^{c3}}|{4:^{c4}}".format(formatted_name, r, d, g, ng, c0=cols[0] + 9, c1=cols[1], c2=cols[2], c3=cols[3], c4=cols[4]) print grad_string def input_sensitivity(self): diff --git a/GPy/core/parameter.py b/GPy/core/parameter.py index 0144733a..79660162 100644 --- a/GPy/core/parameter.py +++ b/GPy/core/parameter.py @@ -6,7 +6,7 @@ Created on 4 Sep 2013 import itertools import numpy from transformations import Logexp, NegativeLogexp, Logistic -from GPy.util.misc import fast_array_equal +from parameterized import Parentable ###### printing __constraints_name__ = "Constraint" @@ -14,9 +14,9 @@ __index_name__ = "Index" __tie_name__ = "Tied to" __precision__ = numpy.get_printoptions()['precision'] # numpy printing precision used, sublassing numpy ndarray after all __print_threshold__ = 5 -###### +###### -class Param(numpy.ndarray): +class Param(numpy.ndarray, Parentable): """ Parameter object for GPy models. @@ -40,13 +40,12 @@ class Param(numpy.ndarray): See :py:class:`GPy.core.parameterized.Parameterized` for more details. """ __array_priority__ = -numpy.inf # Never give back Param - def __new__(cls, name, input_array, gradient, *args, **kwargs): + def __new__(cls, name, input_array, *args, **kwargs): obj = numpy.atleast_1d(numpy.array(input_array)).view(cls) obj._name_ = name - obj._parent_ = None + obj._direct_parent_ = None obj._parent_index_ = None - obj._gradient_parent_ = None - obj._gradient_ = gradient + obj._highest_parent_ = None obj._current_slice_ = (slice(obj.shape[0]),) obj._realshape_ = obj.shape obj._realsize_ = obj.size @@ -62,9 +61,9 @@ class Param(numpy.ndarray): if obj is None: return self._name_ = getattr(obj, '_name_', None) self._current_slice_ = getattr(obj, '_current_slice_', None) - self._parent_ = getattr(obj, '_parent_', None) + self._direct_parent_ = getattr(obj, '_direct_parent_', None) self._parent_index_ = getattr(obj, '_parent_index_', None) - self._gradient_ = getattr(obj, '_gradient_', None) + self._highest_parent_ = getattr(obj, '_highest_parent_', None) self._tied_to_me_ = getattr(obj, '_tied_to_me_', None) self._tied_to_ = getattr(obj, '_tied_to_', None) self._realshape_ = getattr(obj, '_realshape_', None) @@ -72,7 +71,11 @@ class Param(numpy.ndarray): self._realndim_ = getattr(obj, '_realndim_', None) self._updated_ = getattr(obj, '_updated_', None) self._original_ = getattr(obj, '_original_', None) - self._gradient_parent_ = getattr(obj, '_gradient_parent_', None) + def __eq__(self, other): + return other is self + if other is self: + return True + return super(Param, self).__eq__(other) def __array_wrap__(self, out_arr, context=None): return out_arr.view(numpy.ndarray) @@ -83,10 +86,9 @@ class Param(numpy.ndarray): func, args, state = super(Param, self).__reduce__() return func, args, (state, (self._name_, - self._parent_, + self._direct_parent_, self._parent_index_, - self._gradient_parent_, - self._gradient_, + self._highest_parent_, self._current_slice_, self._realshape_, self._realsize_, @@ -106,16 +108,16 @@ class Param(numpy.ndarray): self._realsize_ = state.pop() self._realshape_ = state.pop() self._current_slice_ = state.pop() - self._gradient_ = state.pop() - self._gradient_parent_ = state.pop() + self._highest_parent_ = state.pop() self._parent_index_ = state.pop() - self._parent_ = state.pop() + self._direct_parent_ = state.pop() self._name_ = state.pop() #=========================================================================== # get/set parameters #=========================================================================== def _set_params(self, param): self.flat = param + self._fire_changed() def _get_params(self): return self.flat @property @@ -132,7 +134,10 @@ class Param(numpy.ndarray): def name(self, new_name): from_name = self.name self._name_ = new_name - self._parent_._name_changed(self, from_name) + self._direct_parent_._name_changed(self, from_name) + @property + def _parameters_(self): + return [] #=========================================================================== # Fixing Parameters: #=========================================================================== @@ -142,21 +147,15 @@ class Param(numpy.ndarray): :param warning: print a warning for overwriting constraints. """ - self._parent_._fix(self,warning) + self._highest_parent_._fix(self,warning) fix = constrain_fixed def unconstrain_fixed(self): """ This parameter will no longer be fixed. """ - self._parent_._unfix(self) + self._highest_parent_._unfix(self) unfix = unconstrain_fixed #=========================================================================== - # Convenience methods: - #=========================================================================== - @property - def is_fixed(self): - return self._parent_._is_fixed(self) - #=========================================================================== # Constrain operations -> done #=========================================================================== def constrain(self, transform, warning=True, update=True): @@ -171,10 +170,10 @@ class Param(numpy.ndarray): if self._original_: # this happens when indexing created a copy of the array self.__setitem__(slice(None), transform.initialize(self), update=False) else: - self._parent_._get_original(self).__setitem__(self._current_slice_, transform.initialize(self), update=False) - self._parent_._add_constrain(self, transform, warning) + self._direct_parent_._get_original(self).__setitem__(self._current_slice_, transform.initialize(self), update=False) + self._highest_parent_._add_constrain(self, transform, warning) if update: - self._parent_.parameters_changed() + self._highest_parent_.parameters_changed() def constrain_positive(self, warning=True): """ :param warning: print a warning if re-constraining parameters. @@ -204,7 +203,7 @@ class Param(numpy.ndarray): remove all :py:class:`GPy.core.transformations.Transformation` transformats of this parameter object. """ - self._parent_._remove_constrain(self, *transforms) + self._highest_parent_._remove_constrain(self, *transforms) def unconstrain_positive(self): """ Remove positive constraint of this parameter. @@ -239,23 +238,22 @@ class Param(numpy.ndarray): if self._original_: # this happens when indexing created a copy of the array self[:] = param else: - self._parent_._get_original(self)[self._current_slice_] = param + self._direct_parent_._get_original(self)[self._current_slice_] = param except ValueError: raise ValueError("Trying to tie {} with shape {} to {} with shape {}".format(self.name, self.shape, param.name, param.shape)) - - self._parent_._get_original(self)._tied_to_ += [param] + self._direct_parent_._get_original(self)._tied_to_ += [param] param._add_tie_listener(self) - self._parent_._set_fixed(self) -# self._parent_._add_tie(self, param) + self._highest_parent_._set_fixed(self) +# self._direct_parent_._add_tie(self, param) def untie(self, *ties): """ remove tie of this parameter to ties it was tied to. """ - [self._parent_._get_original(t)._remove_tie_listener(self) for t in self._tied_to_] - self._tied_to_ = [tied_to for tied_to in self._tied_to_ for t in tied_to._tied_to_me_ if self._parent_index_==self._parent_._get_original(t)._parent_index_] - self._parent_._set_unfixed(self) -# self._parent_._remove_tie(self, *params) + [t._direct_parent_._get_original(t)._remove_tie_listener(self) for t in self._tied_to_] + self._tied_to_ = [tied_to for tied_to in self._tied_to_ for t in tied_to._tied_to_me_ if self._parent_index_==t._direct_parent_._get_original(t)._parent_index_] + self._highest_parent_._set_unfixed(self) +# self._direct_parent_._remove_tie(self, *params) def _fire_changed(self): for tied, ind in self._tied_to_me_.iteritems(): tied._on_change(self.base, list(ind)) @@ -266,7 +264,7 @@ class Param(numpy.ndarray): if t._parent_index_ == to_remove._parent_index_: new_index = list(set(t._raveled_index()) - set(to_remove._raveled_index())) if new_index: - tmp = self._parent_._get_original(t)[numpy.unravel_index(new_index,t._realshape_)] + tmp = t._direct_parent_._get_original(t)[numpy.unravel_index(new_index,t._realshape_)] self._tied_to_me_[tmp] = self._tied_to_me_[t] del self._tied_to_me_[t] if len(self._tied_to_me_[tmp]) == 0: @@ -279,7 +277,7 @@ class Param(numpy.ndarray): if self._original_: self.__setitem__(slice(None), val[ind], update=False) else: # this happens when indexing created a copy of the array - self._parent_._get_original(self).__setitem__(self._current_slice_, val[ind], update=False) + self._direct_parent_._get_original(self).__setitem__(self._current_slice_, val[ind], update=False) self._fire_changed() self._updated_ = False #=========================================================================== @@ -291,16 +289,16 @@ class Param(numpy.ndarray): Set prior for this parameter. """ - if not hasattr(self._parent_, '_set_prior'): - raise AttributeError("Parent of type {} does not support priors".format(self._parent_.__class__)) - self._parent_._set_prior(self, prior) + if not hasattr(self._highest_parent_, '_set_prior'): + raise AttributeError("Parent of type {} does not support priors".format(self._highest_parent_.__class__)) + self._highest_parent_._set_prior(self, prior) def unset_prior(self, *priors): """ :param priors: priors to remove from this parameter Remove all priors from this parameter """ - self._parent_._remove_prior(self, *priors) + self._highest_parent_._remove_prior(self, *priors) #=========================================================================== # Array operations -> done #=========================================================================== @@ -321,7 +319,7 @@ class Param(numpy.ndarray): numpy.ndarray.__setitem__(self, s, val) self._fire_changed() if update: - self._parent_.parameters_changed() + self._highest_parent_.parameters_changed() #=========================================================================== # Index Operations: #=========================================================================== @@ -368,22 +366,48 @@ class Param(numpy.ndarray): return numpy.r_[:b] return itertools.imap(f, itertools.izip_longest(slice_index[:self._realndim_], self._realshape_, fillvalue=slice(self.size))) #=========================================================================== - # Printing -> done + # Convienience #=========================================================================== @property - def _desc(self): - if self.size <= 1: return "%f"%self - else: return str(self.shape) - @property - def _constr(self): - return ' '.join(map(lambda c: str(c[0]) if c[1].size==self._realsize_ else "{"+str(c[0])+"}", self._parent_._constraints_iter_items(self))) + def is_fixed(self): + return self._highest_parent_._is_fixed(self) def round(self, decimals=0, out=None): view = super(Param, self).round(decimals, out).view(Param) view.__array_finalize__(self) return view round.__doc__ = numpy.round.__doc__ + def _get_original(self, param): + return self + #=========================================================================== + # Printing -> done + #=========================================================================== + @property + def _description_str(self): + if self.size <= 1: return ["%f"%self] + else: return [str(self.shape)] + def _parameter_names(self, add_name): + return [self.name] + @property + def flattened_parameters(self): + return [self] + @property + def parameter_shapes(self): + return [self.shape] + @property + def _constraints_str(self): + return [' '.join(map(lambda c: str(c[0]) if c[1].size==self._realsize_ else "{"+str(c[0])+"}", self._highest_parent_._constraints_iter_items(self)))] + @property + def _ties_str(self): + return [t._short() for t in self._tied_to_] or [''] + @property + def name_hirarchical(self): + if self.has_parent(): + return self._direct_parent_.hirarchy_name()+self.name + return self.name def __repr__(self, *args, **kwargs): - return "\033[1m{x:s}\033[0;0m:\n".format(x=self.name)+super(Param, self).__repr__(*args,**kwargs) + name = "\033[1m{x:s}\033[0;0m:\n".format( + x=self.name_hirarchical) + return name + super(Param, self).__repr__(*args,**kwargs) def _ties_for(self, rav_index): ties = numpy.empty(shape=(len(self._tied_to_), numpy.size(rav_index)), dtype=Param) for i, tied_to in enumerate(self._tied_to_): @@ -395,7 +419,7 @@ class Param(numpy.ndarray): #[ties.__setitem__(i, ties[i] + [tied_to]) for i in t._raveled_index()] return map(lambda a: sum(a,[]), zip(*[[[tie.flatten()] if tx!=None else [] for tx in t] for t,tie in zip(ties,self._tied_to_)])) def _constraints_for(self, rav_index): - return self._parent_._constraints_for(self, rav_index) + return self._highest_parent_._constraints_for(self, rav_index) def _indices(self, slice_index=None): # get a int-array containing all indices in the first axis. if slice_index is None: @@ -412,17 +436,18 @@ class Param(numpy.ndarray): def _max_len_names(self, gen, header): return reduce(lambda a, b:max(a, len(b)), gen, len(header)) def _max_len_values(self): - return reduce(lambda a, b:max(a, len("{x:=.{0}G}".format(__precision__, x=b))), self.flat, len(self.name)) + return reduce(lambda a, b:max(a, len("{x:=.{0}g}".format(__precision__, x=b))), self.flat, len(self.name)) def _max_len_index(self, ind): return reduce(lambda a, b:max(a, len(str(b))), ind, len(__index_name__)) - def _short(self, slice_index=None): + def _short(self): # short string to print + name = self._direct_parent_.hirarchy_name() + self.name if self._realsize_ < 2: - return self.name - ind = self._indices(slice_index) + return name + ind = self._indices() if ind.size > 4: indstr = ','.join(map(str,ind[:2])) + "..." + ','.join(map(str,ind[-2:])) else: indstr = ','.join(map(str,ind)) - return self.name+'['+indstr+']' + return name+'['+indstr+']' def __str__(self, constr_matrix=None, indices=None, ties=None, lc=None, lx=None, li=None, lt=None): filter_ = self._current_slice_ vals = self.flat @@ -437,7 +462,7 @@ class Param(numpy.ndarray): if lt is None: lt = self._max_len_names(ties, __tie_name__) header = " {i:^{2}s} | \033[1m{x:^{1}s}\033[0;0m | {c:^{0}s} | {t:^{3}s}".format(lc,lx,li,lt, x=self.name, c=__constraints_name__, i=__index_name__, t=__tie_name__) # nice header for printing if not ties: ties = itertools.cycle(['']) - return "\n".join([header]+[" {i!s:^{3}s} | {x: >{1}.{2}G} | {c:^{0}s} | {t:^{4}s} ".format(lc,lx,__precision__,li,lt, x=x, c=" ".join(map(str,c)), t=(t or ''), i=i) for i,x,c,t in itertools.izip(indices,vals,constr_matrix,ties)]) # return all the constraints with right indices + return "\n".join([header]+[" {i!s:^{3}s} | {x: >{1}.{2}g} | {c:^{0}s} | {t:^{4}s} ".format(lc,lx,__precision__,li,lt, x=x, c=" ".join(map(str,c)), t=(t or ''), i=i) for i,x,c,t in itertools.izip(indices,vals,constr_matrix,ties)]) # return all the constraints with right indices #except: return super(Param, self).__str__() class ParamConcatenation(object): @@ -449,7 +474,12 @@ class ParamConcatenation(object): See :py:class:`GPy.core.parameter.Param` for more details on constraining. """ - self.params = params + #self.params = params + self.params = [] + for p in params: + for p in p.flattened_parameters: + if p not in self.params: + self.params.append(p) self._param_sizes = [p.size for p in self.params] startstops = numpy.cumsum([0] + self._param_sizes) self._param_slices_ = [slice(start, stop) for start,stop in zip(startstops, startstops[1:])] @@ -458,7 +488,7 @@ class ParamConcatenation(object): #=========================================================================== def __getitem__(self, s): ind = numpy.zeros(sum(self._param_sizes), dtype=bool); ind[s] = True; - params = [p.flatten()[ind[ps]] for p,ps in zip(self.params, self._param_slices_) if numpy.any(p.flat[ind[ps]])] + params = [p._get_params()[ind[ps]] for p,ps in zip(self.params, self._param_slices_) if numpy.any(p._get_params()[ind[ps]])] if len(params)==1: return params[0] return ParamConcatenation(params) def __setitem__(self, s, val, update=True): @@ -467,7 +497,7 @@ class ParamConcatenation(object): [numpy.place(p, ind[ps], vals[ps]) and p._fire_changed() for p, ps in zip(self.params, self._param_slices_)] if update: - self.params[0]._parent_.parameters_changed() + self.params[0]._highest_parent_.parameters_changed() def _vals(self): return numpy.hstack([p._get_params() for p in self.params]) #=========================================================================== @@ -505,6 +535,8 @@ class ParamConcatenation(object): def unconstrain_bounded(self, lower, upper): [param.unconstrain_bounded(lower, upper) for param in self.params] unconstrain_bounded.__doc__ = Param.unconstrain_bounded.__doc__ + def untie(self, *ties): + [param.untie(*ties) for param in self.params] __lt__ = lambda self, val: self._vals() 1: - self.tie_params(index) - for i, t in zip(prev_constr_ind, prev_constr): - self.constrain(np.where(index_param == i)[0], t) - - def _get_params(self): - return np.hstack(self._parameters_) - return np.hstack([p._get_params() for p in self.parts]) +# def _follow_constrains(self, K1, K2): +# +# # Build the array that allows to go from the initial indices of the param to the new ones +# K1_param = [] +# n = 0 +# for k1 in K1.parts: +# K1_param += [range(n, n + k1.num_params)] +# n += k1.num_params +# n = 0 +# K2_param = [] +# for k2 in K2.parts: +# K2_param += [range(K1.num_params + n, K1.num_params + n + k2.num_params)] +# n += k2.num_params +# index_param = [] +# for p1 in K1_param: +# for p2 in K2_param: +# index_param += p1 + p2 +# index_param = np.array(index_param) +# +# # Get the ties and constrains of the kernels before the multiplication +# prev_ties = K1.tied_indices + [arr + K1.num_params for arr in K2.tied_indices] +# +# prev_constr_ind = [K1.constrained_indices] + [K1.num_params + i for i in K2.constrained_indices] +# prev_constr = K1.constraints + K2.constraints +# +# # prev_constr_fix = K1.fixed_indices + [arr + K1.num_params for arr in K2.fixed_indices] +# # prev_constr_fix_values = K1.fixed_values + K2.fixed_values +# +# # follow the previous ties +# for arr in prev_ties: +# for j in arr: +# index_param[np.where(index_param == j)[0]] = arr[0] +# +# # ties and constrains +# for i in range(K1.num_params + K2.num_params): +# index = np.where(index_param == i)[0] +# if index.size > 1: +# self.tie_params(index) +# for i, t in zip(prev_constr_ind, prev_constr): +# self.constrain(np.where(index_param == i)[0], t) +# +# def _get_params(self): +# return np.hstack(self._parameters_) +# return np.hstack([p._get_params() for p in self._parameters_]) - def _set_params(self, x): - import ipdb;ipdb.set_trace() - [p._set_params(x[s]) for p, s in zip(self.parts, self.param_slices)] +# def _set_params(self, x): +# import ipdb;ipdb.set_trace() +# [p._set_params(x[s]) for p, s in zip(self._parameters_, self._param_slices_)] - def _get_param_names(self): - # this is a bit nasty: we want to distinguish between parts with the same name by appending a count - part_names = np.array([k.name for k in self.parts], dtype=np.str) - counts = [np.sum(part_names == ni) for i, ni in enumerate(part_names)] - cum_counts = [np.sum(part_names[i:] == ni) for i, ni in enumerate(part_names)] - names = [name + '_' + str(cum_count) if count > 1 else name for name, count, cum_count in zip(part_names, counts, cum_counts)] - - return sum([[name + '_' + n for n in k._get_param_names()] for name, k in zip(names, self.parts)], []) +# def _get_param_names(self): +# # this is a bit nasty: we want to distinguish between parts with the same name by appending a count +# part_names = np.array([k.name for k in self._parameters_], dtype=np.str) +# counts = [np.sum(part_names == ni) for i, ni in enumerate(part_names)] +# cum_counts = [np.sum(part_names[i:] == ni) for i, ni in enumerate(part_names)] +# names = [name + '_' + str(cum_count) if count > 1 else name for name, count, cum_count in zip(part_names, counts, cum_counts)] +# +# return sum([[name + '_' + n for n in k._get_param_names()] for name, k in zip(names, self._parameters_)], []) def K(self, X, X2=None, which_parts='all'): """ @@ -357,17 +334,17 @@ class kern(Parameterized): handles this as X2 == X. :param which_parts: a list of booleans detailing whether to include each of the part functions. By default, 'all' - indicates [True]*self.num_parts + indicates all parts """ if which_parts == 'all': - which_parts = [True] * self.num_parts + which_parts = [True] * self.size assert X.shape[1] == self.input_dim if X2 is None: target = np.zeros((X.shape[0], X.shape[0])) - [p.K(X[:, i_s], None, target=target) for p, i_s, part_i_used in zip(self.parts, self.input_slices, which_parts) if part_i_used] + [p.K(X[:, i_s], None, target=target) for p, i_s, part_i_used in zip(self._parameters_, self.input_slices, which_parts) if part_i_used] else: target = np.zeros((X.shape[0], X2.shape[0])) - [p.K(X[:, i_s], X2[:, i_s], target=target) for p, i_s, part_i_used in zip(self.parts, self.input_slices, which_parts) if part_i_used] + [p.K(X[:, i_s], X2[:, i_s], target=target) for p, i_s, part_i_used in zip(self._parameters_, self.input_slices, which_parts) if part_i_used] return target def dK_dtheta(self, dL_dK, X, X2=None): @@ -384,11 +361,11 @@ class kern(Parameterized): returns: dL_dtheta """ assert X.shape[1] == self.input_dim - target = np.zeros(self.num_params) + target = np.zeros(self.size) if X2 is None: - [p.dK_dtheta(dL_dK, X[:, i_s], None, target[ps]) for p, i_s, ps, in zip(self.parts, self.input_slices, self.param_slices)] + [p.dK_dtheta(dL_dK, X[:, i_s], None, target[ps]) for p, i_s, ps, in zip(self._parameters_, self.input_slices, self._param_slices_)] else: - [p.dK_dtheta(dL_dK, X[:, i_s], X2[:, i_s], target[ps]) for p, i_s, ps, in zip(self.parts, self.input_slices, self.param_slices)] + [p.dK_dtheta(dL_dK, X[:, i_s], X2[:, i_s], target[ps]) for p, i_s, ps, in zip(self._parameters_, self.input_slices, self._param_slices_)] return self._transform_gradients(target) @@ -404,18 +381,18 @@ class kern(Parameterized): target = np.zeros_like(X) if X2 is None: - [p.dK_dX(dL_dK, X[:, i_s], None, target[:, i_s]) for p, i_s in zip(self.parts, self.input_slices)] + [p.dK_dX(dL_dK, X[:, i_s], None, target[:, i_s]) for p, i_s in zip(self._parameters_, self.input_slices)] else: - [p.dK_dX(dL_dK, X[:, i_s], X2[:, i_s], target[:, i_s]) for p, i_s in zip(self.parts, self.input_slices)] + [p.dK_dX(dL_dK, X[:, i_s], X2[:, i_s], target[:, i_s]) for p, i_s in zip(self._parameters_, self.input_slices)] return target def Kdiag(self, X, which_parts='all'): """Compute the diagonal of the covariance function for inputs X.""" if which_parts == 'all': - which_parts = [True] * self.num_parts + which_parts = [True] * self.size assert X.shape[1] == self.input_dim target = np.zeros(X.shape[0]) - [p.Kdiag(X[:, i_s], target=target) for p, i_s, part_on in zip(self.parts, self.input_slices, which_parts) if part_on] + [p.Kdiag(X[:, i_s], target=target) for p, i_s, part_on in zip(self._parameters_, self.input_slices, which_parts) if part_on] return target def dKdiag_dtheta(self, dL_dKdiag, X): @@ -423,49 +400,49 @@ class kern(Parameterized): assert X.shape[1] == self.input_dim assert dL_dKdiag.size == X.shape[0] target = np.zeros(self.num_params) - [p.dKdiag_dtheta(dL_dKdiag, X[:, i_s], target[ps]) for p, i_s, ps in zip(self.parts, self.input_slices, self.param_slices)] + [p.dKdiag_dtheta(dL_dKdiag, X[:, i_s], target[ps]) for p, i_s, ps in zip(self._parameters_, self.input_slices, self._param_slices_)] return self._transform_gradients(target) def dKdiag_dX(self, dL_dKdiag, X): assert X.shape[1] == self.input_dim target = np.zeros_like(X) - [p.dKdiag_dX(dL_dKdiag, X[:, i_s], target[:, i_s]) for p, i_s in zip(self.parts, self.input_slices)] + [p.dKdiag_dX(dL_dKdiag, X[:, i_s], target[:, i_s]) for p, i_s in zip(self._parameters_, self.input_slices)] return target def psi0(self, Z, mu, S): target = np.zeros(mu.shape[0]) - [p.psi0(Z[:, i_s], mu[:, i_s], S[:, i_s], target) for p, i_s in zip(self.parts, self.input_slices)] + [p.psi0(Z[:, i_s], mu[:, i_s], S[:, i_s], target) for p, i_s in zip(self._parameters_, self.input_slices)] return target def dpsi0_dtheta(self, dL_dpsi0, Z, mu, S): target = np.zeros(self.num_params) - [p.dpsi0_dtheta(dL_dpsi0, Z[:, i_s], mu[:, i_s], S[:, i_s], target[ps]) for p, ps, i_s in zip(self.parts, self.param_slices, self.input_slices)] + [p.dpsi0_dtheta(dL_dpsi0, Z[:, i_s], mu[:, i_s], S[:, i_s], target[ps]) for p, ps, i_s in zip(self._parameters_, self._param_slices_, self.input_slices)] return self._transform_gradients(target) def dpsi0_dmuS(self, dL_dpsi0, Z, mu, S): target_mu, target_S = np.zeros_like(mu), np.zeros_like(S) - [p.dpsi0_dmuS(dL_dpsi0, Z[:, i_s], mu[:, i_s], S[:, i_s], target_mu[:, i_s], target_S[:, i_s]) for p, i_s in zip(self.parts, self.input_slices)] + [p.dpsi0_dmuS(dL_dpsi0, Z[:, i_s], mu[:, i_s], S[:, i_s], target_mu[:, i_s], target_S[:, i_s]) for p, i_s in zip(self._parameters_, self.input_slices)] return target_mu, target_S def psi1(self, Z, mu, S): target = np.zeros((mu.shape[0], Z.shape[0])) - [p.psi1(Z[:, i_s], mu[:, i_s], S[:, i_s], target) for p, i_s in zip(self.parts, self.input_slices)] + [p.psi1(Z[:, i_s], mu[:, i_s], S[:, i_s], target) for p, i_s in zip(self._parameters_, self.input_slices)] return target def dpsi1_dtheta(self, dL_dpsi1, Z, mu, S): target = np.zeros((self.num_params)) - [p.dpsi1_dtheta(dL_dpsi1, Z[:, i_s], mu[:, i_s], S[:, i_s], target[ps]) for p, ps, i_s in zip(self.parts, self.param_slices, self.input_slices)] + [p.dpsi1_dtheta(dL_dpsi1, Z[:, i_s], mu[:, i_s], S[:, i_s], target[ps]) for p, ps, i_s in zip(self._parameters_, self._param_slices_, self.input_slices)] return self._transform_gradients(target) def dpsi1_dZ(self, dL_dpsi1, Z, mu, S): target = np.zeros_like(Z) - [p.dpsi1_dZ(dL_dpsi1, Z[:, i_s], mu[:, i_s], S[:, i_s], target[:, i_s]) for p, i_s in zip(self.parts, self.input_slices)] + [p.dpsi1_dZ(dL_dpsi1, Z[:, i_s], mu[:, i_s], S[:, i_s], target[:, i_s]) for p, i_s in zip(self._parameters_, self.input_slices)] return target def dpsi1_dmuS(self, dL_dpsi1, Z, mu, S): """return shapes are num_samples,num_inducing,input_dim""" target_mu, target_S = np.zeros((2, mu.shape[0], mu.shape[1])) - [p.dpsi1_dmuS(dL_dpsi1, Z[:, i_s], mu[:, i_s], S[:, i_s], target_mu[:, i_s], target_S[:, i_s]) for p, i_s in zip(self.parts, self.input_slices)] + [p.dpsi1_dmuS(dL_dpsi1, Z[:, i_s], mu[:, i_s], S[:, i_s], target_mu[:, i_s], target_S[:, i_s]) for p, i_s in zip(self._parameters_, self.input_slices)] return target_mu, target_S def psi2(self, Z, mu, S): @@ -478,13 +455,13 @@ class kern(Parameterized): """ target = np.zeros((mu.shape[0], Z.shape[0], Z.shape[0])) - [p.psi2(Z[:, i_s], mu[:, i_s], S[:, i_s], target) for p, i_s in zip(self.parts, self.input_slices)] + [p.psi2(Z[:, i_s], mu[:, i_s], S[:, i_s], target) for p, i_s in zip(self._parameters_, self.input_slices)] # compute the "cross" terms # TODO: input_slices needed crossterms = 0 - for [p1, i_s1], [p2, i_s2] in itertools.combinations(zip(self.parts, self.input_slices), 2): + for [p1, i_s1], [p2, i_s2] in itertools.combinations(zip(self._parameters_, self.input_slices), 2): if i_s1 == i_s2: # TODO psi1 this must be faster/better/precached/more nice tmp1 = np.zeros((mu.shape[0], Z.shape[0])) @@ -501,14 +478,14 @@ class kern(Parameterized): def dpsi2_dtheta(self, dL_dpsi2, Z, mu, S): """Gradient of the psi2 statistics with respect to the parameters.""" target = np.zeros(self.num_params) - [p.dpsi2_dtheta(dL_dpsi2, Z[:, i_s], mu[:, i_s], S[:, i_s], target[ps]) for p, i_s, ps in zip(self.parts, self.input_slices, self.param_slices)] + [p.dpsi2_dtheta(dL_dpsi2, Z[:, i_s], mu[:, i_s], S[:, i_s], target[ps]) for p, i_s, ps in zip(self._parameters_, self.input_slices, self._param_slices_)] # compute the "cross" terms # TODO: better looping, input_slices - for i1, i2 in itertools.permutations(range(len(self.parts)), 2): - p1, p2 = self.parts[i1], self.parts[i2] + for i1, i2 in itertools.permutations(range(len(self._parameters_)), 2): + p1, p2 = self._parameters_[i1], self._parameters_[i2] # ipsl1, ipsl2 = self.input_slices[i1], self.input_slices[i2] - ps1, ps2 = self.param_slices[i1], self.param_slices[i2] + ps1, ps2 = self._param_slices_[i1], self._param_slices_[i2] tmp = np.zeros((mu.shape[0], Z.shape[0])) p1.psi1(Z, mu, S, tmp) @@ -518,12 +495,12 @@ class kern(Parameterized): def dpsi2_dZ(self, dL_dpsi2, Z, mu, S): target = np.zeros_like(Z) - [p.dpsi2_dZ(dL_dpsi2, Z[:, i_s], mu[:, i_s], S[:, i_s], target[:, i_s]) for p, i_s in zip(self.parts, self.input_slices)] + [p.dpsi2_dZ(dL_dpsi2, Z[:, i_s], mu[:, i_s], S[:, i_s], target[:, i_s]) for p, i_s in zip(self._parameters_, self.input_slices)] # target *= 2 # compute the "cross" terms # TODO: we need input_slices here. - for p1, p2 in itertools.permutations(self.parts, 2): + for p1, p2 in itertools.permutations(self._parameters_, 2): # if p1.name == 'linear' and p2.name == 'linear': # raise NotImplementedError("We don't handle linear/linear cross-terms") tmp = np.zeros((mu.shape[0], Z.shape[0])) @@ -534,11 +511,11 @@ class kern(Parameterized): def dpsi2_dmuS(self, dL_dpsi2, Z, mu, S): target_mu, target_S = np.zeros((2, mu.shape[0], mu.shape[1])) - [p.dpsi2_dmuS(dL_dpsi2, Z[:, i_s], mu[:, i_s], S[:, i_s], target_mu[:, i_s], target_S[:, i_s]) for p, i_s in zip(self.parts, self.input_slices)] + [p.dpsi2_dmuS(dL_dpsi2, Z[:, i_s], mu[:, i_s], S[:, i_s], target_mu[:, i_s], target_S[:, i_s]) for p, i_s in zip(self._parameters_, self.input_slices)] # compute the "cross" terms # TODO: we need input_slices here. - for p1, p2 in itertools.permutations(self.parts, 2): + for p1, p2 in itertools.permutations(self._parameters_, 2): # if p1.name == 'linear' and p2.name == 'linear': # raise NotImplementedError("We don't handle linear/linear cross-terms") tmp = np.zeros((mu.shape[0], Z.shape[0])) @@ -549,7 +526,7 @@ class kern(Parameterized): def plot(self, x=None, plot_limits=None, which_parts='all', resolution=None, *args, **kwargs): if which_parts == 'all': - which_parts = [True] * self.num_parts + which_parts = [True] * self.size if self.input_dim == 1: if x is None: x = np.zeros((1, 1)) diff --git a/GPy/kern/parts/bias.py b/GPy/kern/parts/bias.py index 7010db12..22e0882b 100644 --- a/GPy/kern/parts/bias.py +++ b/GPy/kern/parts/bias.py @@ -15,13 +15,9 @@ class Bias(Kernpart): :param variance: the variance of the kernel :type variance: float """ - super(Bias, self).__init__(input_dim) - self.input_dim = input_dim - self.num_params = 1 - self.name = 'bias' - - self.variance = Param(lambda: self.name+"_variance", variance, None) - self.set_as_parameters(self.variance) + super(Bias, self).__init__(input_dim, 'bias') + self.variance = Param("variance", variance, None) + self.add_parameter(self.variance) #self._set_params(np.array([variance]).flatten()) # def _get_params(self): diff --git a/GPy/kern/parts/kernpart.py b/GPy/kern/parts/kernpart.py index 461a152c..9843bf7d 100644 --- a/GPy/kern/parts/kernpart.py +++ b/GPy/kern/parts/kernpart.py @@ -1,10 +1,10 @@ # Copyright (c) 2012, GPy authors (see AUTHORS.txt). # Licensed under the BSD 3-clause license (see LICENSE.txt) -from ...core.parameter import Param #from ...core.parameterized.Parameterized import set_as_parameter +from GPy.core.parameterized import Parameterized -class Kernpart(object): - def __init__(self,input_dim): +class Kernpart(Parameterized): + def __init__(self,input_dim,name): """ The base class for a kernpart: a positive definite function which forms part of a covariance function (kernel). @@ -14,21 +14,14 @@ class Kernpart(object): Do not instantiate. """ + super(Kernpart, self).__init__(name) # the input dimensionality for the covariance self.input_dim = input_dim # the number of optimisable parameters - self.num_params = 1 # the name of the covariance function. - self.name = 'unnamed' # link to parameterized objects self._parameters_ = [] - def set_as_parameters(self, *params): - self._parameters_.extend(params) - - def parameters_changed(self): - pass - # def set_as_parameter_named(self, name, gradient, index=None, *args, **kwargs): # """ # :param names: name of parameter to set as parameter diff --git a/GPy/kern/parts/linear.py b/GPy/kern/parts/linear.py index ffcbcf5e..f36e5968 100644 --- a/GPy/kern/parts/linear.py +++ b/GPy/kern/parts/linear.py @@ -7,6 +7,7 @@ import numpy as np from ...util.linalg import tdot from ...util.misc import fast_array_equal from scipy import weave +from GPy.core.parameter import Param class Linear(Kernpart): """ @@ -26,11 +27,9 @@ class Linear(Kernpart): """ def __init__(self, input_dim, variances=None, ARD=False): - self.input_dim = input_dim + super(Linear, self).__init__(input_dim, 'linear') self.ARD = ARD if ARD == False: - self.num_params = 1 - self.name = 'linear' if variances is not None: variances = np.asarray(variances) assert variances.size == 1, "Only one variance needed for non-ARD kernel" @@ -38,32 +37,33 @@ class Linear(Kernpart): variances = np.ones(1) self._Xcache, self._X2cache = np.empty(shape=(2,)) else: - self.num_params = self.input_dim - self.name = 'linear' if variances is not None: variances = np.asarray(variances) assert variances.size == self.input_dim, "bad number of lengthscales" else: variances = np.ones(self.input_dim) - self._set_params(variances.flatten()) + + self.variances = Param('variances', variances) + self.add_parameters(self.variances) # initialize cache self._Z, self._mu, self._S = np.empty(shape=(3, 1)) - self._X, self._X2, self._params = np.empty(shape=(3, 1)) + self._X, self._X2 = np.empty(shape=(2, 1)) - def _get_params(self): - return self.variances - - def _set_params(self, x): - assert x.size == (self.num_params) - self.variances = x +# def _get_params(self): +# return self.variances +# +# def _set_params(self, x): +# assert x.size == (self.num_params) +# self.variances = x + def parameters_changed(self): self.variances2 = np.square(self.variances) - - def _get_param_names(self): - if self.num_params == 1: - return ['variance'] - else: - return ['variance_%i' % i for i in range(self.variances.size)] +# +# def _get_param_names(self): +# if self.num_params == 1: +# return ['variance'] +# else: +# return ['variance_%i' % i for i in range(self.variances.size)] def K(self, X, X2, target): if self.ARD: diff --git a/GPy/kern/parts/prod.py b/GPy/kern/parts/prod.py index 0549ea22..2569c51c 100644 --- a/GPy/kern/parts/prod.py +++ b/GPy/kern/parts/prod.py @@ -18,35 +18,37 @@ class Prod(Kernpart): """ def __init__(self,k1,k2,tensor=False): - self.num_params = k1.num_params + k2.num_params - self.name = '['+k1.name + '**' + k2.name +']' - self.k1 = k1 - self.k2 = k2 if tensor: - self.input_dim = k1.input_dim + k2.input_dim - self.slice1 = slice(0,self.k1.input_dim) - self.slice2 = slice(self.k1.input_dim,self.k1.input_dim+self.k2.input_dim) + super(Prod, self).__init__(k1.input_dim + k2.input_dim, '['+k1.name + '**' + k2.name +']') else: assert k1.input_dim == k2.input_dim, "Error: The input spaces of the kernels to sum don't have the same dimension." - self.input_dim = k1.input_dim - self.slice1 = slice(0,self.input_dim) - self.slice2 = slice(0,self.input_dim) + super(Prod, self).__init__(k1.input_dim, '['+k1.name + '*' + k2.name +']') + #self.num_params = k1.num_params + k2.num_params + self.k1 = k1 + self.k2 = k2 +# if tensor: +# self.slice1 = slice(0,self.k1.input_dim) +# self.slice2 = slice(self.k1.input_dim,self.k1.input_dim+self.k2.input_dim) +# else: +# self.slice1 = slice(0,self.input_dim) +# self.slice2 = slice(0,self.input_dim) - self._X, self._X2, self._params = np.empty(shape=(3,1)) - self._set_params(np.hstack((k1._get_params(),k2._get_params()))) + self._X, self._X2 = np.empty(shape=(2,1)) + self.add_parameters(self.k1, self.k2) +# self._set_params(np.hstack((k1._get_params(),k2._get_params()))) - def _get_params(self): - """return the value of the parameters.""" - return np.hstack((self.k1._get_params(), self.k2._get_params())) - - def _set_params(self,x): - """set the value of the parameters.""" - self.k1._set_params(x[:self.k1.num_params]) - self.k2._set_params(x[self.k1.num_params:]) - - def _get_param_names(self): - """return parameter names.""" - return [self.k1.name + '_' + param_name for param_name in self.k1._get_param_names()] + [self.k2.name + '_' + param_name for param_name in self.k2._get_param_names()] +# def _get_params(self): +# """return the value of the parameters.""" +# return np.hstack((self.k1._get_params(), self.k2._get_params())) +# +# def _set_params(self,x): +# """set the value of the parameters.""" +# self.k1._set_params(x[:self.k1.num_params]) +# self.k2._set_params(x[self.k1.num_params:]) +# +# def _get_param_names(self): +# """return parameter names.""" +# return [self.k1.name + '_' + param_name for param_name in self.k1._get_param_names()] + [self.k2.name + '_' + param_name for param_name in self.k2._get_param_names()] def K(self,X,X2,target): self._K_computations(X,X2) diff --git a/GPy/kern/parts/rbf.py b/GPy/kern/parts/rbf.py index 4912d241..df615ac9 100644 --- a/GPy/kern/parts/rbf.py +++ b/GPy/kern/parts/rbf.py @@ -33,20 +33,17 @@ class RBF(Kernpart): """ def __init__(self, input_dim, variance=1., lengthscale=None, ARD=False, name='rbf'): - super(RBF, self).__init__(input_dim) + super(RBF, self).__init__(input_dim, name) self.input_dim = input_dim - self.name = name self.ARD = ARD if not ARD: - self.num_params = 2 if lengthscale is not None: lengthscale = np.asarray(lengthscale) assert lengthscale.size == 1, "Only one lengthscale needed for non-ARD kernel" else: lengthscale = np.ones(1) else: - self.num_params = self.input_dim + 1 if lengthscale is not None: lengthscale = np.asarray(lengthscale) assert lengthscale.size == self.input_dim, "bad number of lengthscales" @@ -54,9 +51,10 @@ class RBF(Kernpart): lengthscale = np.ones(self.input_dim) #self._set_params(np.hstack((variance, lengthscale.flatten()))) - self.variance = Param(lambda: self.name+'_variance', variance, None) - self.lengthscale = Param(lambda: self.name+'_lengthscale', lengthscale, None) - self.set_as_parameters(self.variance, self.lengthscale) + self.variance = Param('variance', variance, None) + self.lengthscale = Param('lengthscale', lengthscale, None) + + self.add_parameters(self.variance, self.lengthscale) # self.set_as_parameter('variance', self.variance, None) # self.set_as_parameter('lengthscale', self.lengthscale, None) @@ -68,37 +66,35 @@ class RBF(Kernpart): self.weave_options = {'headers' : [''], 'extra_compile_args': ['-fopenmp -O3'], # -march=native'], 'extra_link_args' : ['-lgomp']} - def parameters_changed(self): self.lengthscale2 = np.square(self.lengthscale) - # reset cached results - self._X, self._X2, self._params_save = np.empty(shape=(3, 1)) - self._Z, self._mu, self._S = np.empty(shape=(3, 1)) # cached versions of Z,mu,S - - def _get_params(self): - return np.hstack((self.variance, self.lengthscale)) -# - def _set_params(self, x): - assert x.size == (self.num_params) - #self.variance = x[0] - #self.lengthscale = x[1:] - - #self.lengthscale2 = np.square(self.lengthscale) - # reset cached results #self._X, self._X2, self._params_save = np.empty(shape=(3, 1)) #self._Z, self._mu, self._S = np.empty(shape=(3, 1)) # cached versions of Z,mu,S - - def _get_param_names(self): - if self.num_params == 2: - return ['variance', 'lengthscale'] - else: - return ['variance'] + ['lengthscale_%i' % i for i in range(self.lengthscale.size)] + self._X, self._X2 = np.empty(shape=(2, 1)) + self._Z, self._mu, self._S = np.empty(shape=(3, 1)) # cached versions of Z,mu,S - def _dK_dvariance(self, model, target): - pass +# def _get_params(self): +# return np.hstack((self.variance, self.lengthscale)) +# # +# def _set_params(self, x): +# assert x.size == (self.num_params) +# #self.variance = x[0] +# #self.lengthscale = x[1:] +# +# #self.lengthscale2 = np.square(self.lengthscale) +# +# # reset cached results +# #self._X, self._X2, self._params_save = np.empty(shape=(3, 1)) +# #self._Z, self._mu, self._S = np.empty(shape=(3, 1)) # cached versions of Z,mu,S +# +# def _get_param_names(self): +# if self.num_params == 2: +# return ['variance', 'lengthscale'] +# else: +# return ['variance'] + ['lengthscale_%i' % i for i in range(self.lengthscale.size)] def K(self, X, X2, target): self._K_computations(X, X2) @@ -243,10 +239,10 @@ class RBF(Kernpart): #---------------------------------------# def _K_computations(self, X, X2): - params = self._get_params() - if not (fast_array_equal(X, self._X) and fast_array_equal(X2, self._X2) and fast_array_equal(self._params_save , params)): + #params = self._get_params() + if not (fast_array_equal(X, self._X) and fast_array_equal(X2, self._X2)):# and fast_array_equal(self._params_save , params)): self._X = X.copy() - self._params_save = params.copy() + #self._params_save = params.copy() if X2 is None: self._X2 = None X = X / self.lengthscale diff --git a/GPy/kern/parts/rbfcos.py b/GPy/kern/parts/rbfcos.py index 4d09dfdb..28630f21 100644 --- a/GPy/kern/parts/rbfcos.py +++ b/GPy/kern/parts/rbfcos.py @@ -5,6 +5,7 @@ from kernpart import Kernpart import numpy as np +from ...core.parameter import Param class RBFCos(Kernpart): def __init__(self,input_dim,variance=1.,frequencies=None,bandwidths=None,ARD=False): @@ -41,29 +42,30 @@ class RBFCos(Kernpart): else: bandwidths = np.ones(1) + self.variance = Param('variance', variance) + self.frequencies = Param('frequencies', frequencies) + self.bandwidths = Param('bandwidths', bandwidths) + #initialise cache - self._X, self._X2, self._params = np.empty(shape=(3,1)) + self._X, self._X2 = np.empty(shape=(3,1)) - self._set_params(np.hstack((variance,frequencies.flatten(),bandwidths.flatten()))) +# def _get_params(self): +# return np.hstack((self.variance,self.frequencies, self.bandwidths)) +# def _set_params(self,x): +# assert x.size==(self.num_params) +# if self.ARD: +# self.variance = x[0] +# self.frequencies = x[1:1+self.input_dim] +# self.bandwidths = x[1+self.input_dim:] +# else: +# self.variance, self.frequencies, self.bandwidths = x - def _get_params(self): - return np.hstack((self.variance,self.frequencies, self.bandwidths)) - - def _set_params(self,x): - assert x.size==(self.num_params) - if self.ARD: - self.variance = x[0] - self.frequencies = x[1:1+self.input_dim] - self.bandwidths = x[1+self.input_dim:] - else: - self.variance, self.frequencies, self.bandwidths = x - - def _get_param_names(self): - if self.num_params == 3: - return ['variance','frequency','bandwidth'] - else: - return ['variance']+['frequency_%i'%i for i in range(self.input_dim)]+['bandwidth_%i'%i for i in range(self.input_dim)] +# def _get_param_names(self): +# if self.num_params == 3: +# return ['variance','frequency','bandwidth'] +# else: +# return ['variance']+['frequency_%i'%i for i in range(self.input_dim)]+['bandwidth_%i'%i for i in range(self.input_dim)] def K(self,X,X2,target): self._K_computations(X,X2) @@ -94,6 +96,11 @@ class RBFCos(Kernpart): def dKdiag_dX(self,dL_dKdiag,X,target): pass + def parameters_changed(self): + self._rbf_part = np.exp(-2.*np.pi**2*np.sum(self._dist2*self.bandwidths,-1)) + self._cos_part = np.prod(np.cos(2.*np.pi*self._dist*self.frequencies),-1) + self._dvar = self._rbf_part*self._cos_part + def _K_computations(self,X,X2): if not (np.all(X==self._X) and np.all(X2==self._X2)): if X2 is None: X2 = X @@ -107,11 +114,3 @@ class RBFCos(Kernpart): #ensure the next section is computed: self._params = np.empty(self.num_params) - - if not np.all(self._params == self._get_params()): - self._params == self._get_params().copy() - - self._rbf_part = np.exp(-2.*np.pi**2*np.sum(self._dist2*self.bandwidths,-1)) - self._cos_part = np.prod(np.cos(2.*np.pi*self._dist*self.frequencies),-1) - self._dvar = self._rbf_part*self._cos_part - diff --git a/GPy/kern/parts/spline.py b/GPy/kern/parts/spline.py index e675420c..e59642fd 100644 --- a/GPy/kern/parts/spline.py +++ b/GPy/kern/parts/spline.py @@ -5,6 +5,7 @@ from kernpart import Kernpart import numpy as np import hashlib +from ...core.parameter import Param def theta(x): """Heaviside step function""" return np.where(x>=0.,1.,0.) @@ -25,16 +26,18 @@ class Spline(Kernpart): assert self.input_dim==1 self.num_params = 1 self.name = 'spline' - self._set_params(np.squeeze(variance)) - - def _get_params(self): - return self.variance - - def _set_params(self,x): - self.variance = x - - def _get_param_names(self): - return ['variance'] + self.variance = Param('variance', variance) + self.lengthscale = Param('lengthscale', lengthscale) + self.add_parameters(self.variance, self.lengthscale) + +# def _get_params(self): +# return self.variance +# +# def _set_params(self,x): +# self.variance = x +# +# def _get_param_names(self): +# return ['variance'] def K(self,X,X2,target): assert np.all(X>0), "Spline covariance is for +ve domain only. TODO: symmetrise" diff --git a/GPy/kern/parts/symmetric.py b/GPy/kern/parts/symmetric.py index bbdd5ac0..d47fbd9d 100644 --- a/GPy/kern/parts/symmetric.py +++ b/GPy/kern/parts/symmetric.py @@ -24,19 +24,8 @@ class Symmetric(Kernpart): self.num_params = k.num_params self.name = k.name + '_symm' self.k = k - self._set_params(k._get_params()) - - def _get_params(self): - """return the value of the parameters.""" - return self.k._get_params() - - def _set_params(self,x): - """set the value of the parameters.""" - self.k._set_params(x) - - def _get_param_names(self): - """return parameter names.""" - return self.k._get_param_names() + self.add_parameter(k) + #self._set_params(k._get_params()) def K(self,X,X2,target): """Compute the covariance matrix between X and X2.""" diff --git a/GPy/kern/parts/white.py b/GPy/kern/parts/white.py index 00422c9b..bca0918a 100644 --- a/GPy/kern/parts/white.py +++ b/GPy/kern/parts/white.py @@ -15,12 +15,10 @@ class White(Kernpart): :type variance: float """ def __init__(self,input_dim,variance=1.): - super(White, self).__init__(input_dim) + super(White, self).__init__(input_dim, 'white') self.input_dim = input_dim - self.num_params = 1 - self.name = 'white' - self.variance = Param(lambda: self.name+'_variance', variance, None) - self.set_as_parameters(self.variance) + self.variance = Param('variance', variance, None) + self.add_parameters(self.variance) # self._set_params(np.array([variance]).flatten()) self._psi1 = 0 # TODO: more elegance here diff --git a/GPy/likelihoods/gaussian.py b/GPy/likelihoods/gaussian.py index 5774537b..21b4f4d3 100644 --- a/GPy/likelihoods/gaussian.py +++ b/GPy/likelihoods/gaussian.py @@ -15,7 +15,7 @@ class Gaussian(likelihood): :type normalize: False|True """ def __init__(self, data, variance=1., normalize=False): - super(Gaussian, self).__init__() + super(Gaussian, self).__init__('gaussian') self.is_heteroscedastic = False self.num_params = 1 self.Z = 0. # a correction factor which accounts for the approximation made @@ -33,11 +33,12 @@ class Gaussian(likelihood): self.set_data(data) - self.variance = Param('noise_variance', variance, None) - self.set_as_parameters(self.variance) - + self.variance = Param('variance', variance) self._variance = variance + 1 + self.add_parameter(self.variance) + + # self._set_params(np.asarray(variance)) @@ -62,6 +63,7 @@ class Gaussian(likelihood): # return ["noise_variance"] # # def _set_params(self, x): +# self.variance = x[0] def parameters_changed(self): if np.any(self._variance != self.variance): if np.all(self.variance == 0.):#special case of zero noise diff --git a/GPy/likelihoods/likelihood.py b/GPy/likelihoods/likelihood.py index 4abf89af..d94af758 100644 --- a/GPy/likelihoods/likelihood.py +++ b/GPy/likelihoods/likelihood.py @@ -21,8 +21,8 @@ class likelihood(Parameterized): - self.V : self.precision * self.Y """ - def __init__(self): - Parameterized.__init__(self) + def __init__(self, name=None): + Parameterized.__init__(self, name) # def _get_params(self): # raise NotImplementedError diff --git a/GPy/testing/kernel_tests.py b/GPy/testing/kernel_tests.py index f64dac2b..64d057af 100644 --- a/GPy/testing/kernel_tests.py +++ b/GPy/testing/kernel_tests.py @@ -78,7 +78,7 @@ class KernelTests(unittest.TestCase): self.assertTrue(GPy.kern.kern_test(kern, verbose=verbose)) def test_heterokernel(self): - kern = GPy.kern.hetero(5, mapping=GPy.mappings.Linear(5, 1), transform=GPy.core.transformations.logexp()) + kern = GPy.kern.hetero(5, mapping=GPy.mappings.Linear(5, 1), transform=GPy.core.transformations.Logexp()) self.assertTrue(GPy.kern.kern_test(kern, verbose=verbose)) def test_mlpkernel(self): diff --git a/GPy/testing/parameter_testing.py b/GPy/testing/parameter_testing.py index 0a2a4ff2..4ab9b810 100644 --- a/GPy/testing/parameter_testing.py +++ b/GPy/testing/parameter_testing.py @@ -31,10 +31,10 @@ class Test(unittest.TestCase): self.Y = numpy.random.multivariate_normal(numpy.zeros(self.N), K + numpy.eye(self.N) * .2, self.D).T - self.bgplvm = BayesianGPLVM(Gaussian(self.Y, variance=self.noise_variance), self.Q, self.X, self.X_variance, kernel=self.kern) - self.bgplvm.ensure_default_constraints(warning=False) - self.bgplvm.tie_params("noise_variance|white_variance") - self.bgplvm.constrain_fixed("rbf_var", warning=False) +# self.bgplvm = BayesianGPLVM(Gaussian(self.Y, variance=self.noise_variance), self.Q, self.X, self.X_variance, kernel=self.kern) +# self.bgplvm.ensure_default_constraints(warning=False) +# self.bgplvm.tie_params("noise_variance|white_variance") +# self.bgplvm.constrain_fixed("rbf_var", warning=False) self.parameter = Parameterized([ Parameterized([ Param('X', self.X), @@ -110,7 +110,7 @@ class Test(unittest.TestCase): self.parameter[''].unconstrain() self.parameter.X.constrain_positive() self.parameter.X[:,numpy.s_[0::2]].unconstrain_positive() - assert(numpy.alltrue(self.parameter._constraints_.indices()[0] == numpy.r_[1:self.N*self.Q:2])) + assert(numpy.alltrue(self.parameter.constraints.indices()[0] == numpy.r_[1:self.N*self.Q:2])) def testNdarrayFunc(self): assert(numpy.alltrue(self.parameter.X * self.parameter.X == self.X * self.X))