From 47e4026141b0712777eda3713b150f43d2756c11 Mon Sep 17 00:00:00 2001 From: Max Zwiessele Date: Fri, 28 Feb 2014 16:18:47 +0000 Subject: [PATCH] hierarchy edits. adding removing parameters withing hierarchy --- GPy/core/model.py | 4 +- GPy/core/parameterization/index_operations.py | 12 +++--- GPy/core/parameterization/param.py | 43 +++++++++++++------ GPy/core/parameterization/parameter_core.py | 42 ++++++++++-------- GPy/core/parameterization/parameterized.py | 35 ++++++++++----- GPy/examples/dimensionality_reduction.py | 10 ++--- .../latent_function_inference/var_dtc.py | 5 +-- GPy/kern/_src/kern.py | 10 +++-- GPy/models/sparse_gp_regression.py | 4 +- GPy/testing/parameterized_tests.py | 3 +- GPy/util/caching.py | 2 +- 11 files changed, 106 insertions(+), 64 deletions(-) diff --git a/GPy/core/model.py b/GPy/core/model.py index 0925a199..6fd80d76 100644 --- a/GPy/core/model.py +++ b/GPy/core/model.py @@ -265,7 +265,7 @@ class Model(Parameterized): and numerical gradients is within of unity. """ - x = self._get_params_transformed().copy() + x = self._get_params_transformed() if not verbose: # make sure only to test the selected parameters @@ -283,7 +283,7 @@ class Model(Parameterized): return # just check the global ratio - dx = np.zeros_like(x) + dx = np.zeros(x.shape()) dx[transformed_index] = step * np.sign(np.random.uniform(-1, 1, transformed_index.size)) # evaulate around the point x diff --git a/GPy/core/parameterization/index_operations.py b/GPy/core/parameterization/index_operations.py index b5399741..6450c41c 100644 --- a/GPy/core/parameterization/index_operations.py +++ b/GPy/core/parameterization/index_operations.py @@ -194,9 +194,13 @@ class ParameterIndexOperationsView(object): def shift_right(self, start, size): - raise NotImplementedError, 'Shifting only supported in original ParamIndexOperations' - + self._param_index_ops.shift_right(start+self._offset, size) + def shift_left(self, start, size): + self._param_index_ops.shift_left(start+self._offset, size) + self._offset -= size + self._size -= size + def clear(self): for i, ind in self.items(): self._param_index_ops.remove(i, ind+self._offset) @@ -232,9 +236,7 @@ class ParameterIndexOperationsView(object): def __getitem__(self, prop): ind = self._filter_index(self._param_index_ops[prop]) - if ind.size > 0: - return ind - raise KeyError, prop + return ind def __str__(self, *args, **kwargs): import pprint diff --git a/GPy/core/parameterization/param.py b/GPy/core/parameterization/param.py index 14cba600..d52442d1 100644 --- a/GPy/core/parameterization/param.py +++ b/GPy/core/parameterization/param.py @@ -269,7 +269,7 @@ class Param(OptimizationHandlable, ObservableArray, Gradcheckable): return [t._short() for t in self._tied_to_] or [''] def __repr__(self, *args, **kwargs): name = "\033[1m{x:s}\033[0;0m:\n".format( - x=self.hirarchy_name()) + x=self.hierarchy_name()) return name + super(Param, self).__repr__(*args, **kwargs) def _ties_for(self, rav_index): # size = sum(p.size for p in self._tied_to_) @@ -303,12 +303,12 @@ class Param(OptimizationHandlable, ObservableArray, Gradcheckable): gen = map(lambda x: " ".join(map(str, x)), gen) 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.hirarchy_name())) + return reduce(lambda a, b:max(a, len("{x:=.{0}g}".format(__precision__, x=b))), self.flat, len(self.hierarchy_name())) def _max_len_index(self, ind): return reduce(lambda a, b:max(a, len(str(b))), ind, len(__index_name__)) def _short(self): # short string to print - name = self.hirarchy_name() + name = self.hierarchy_name() if self._realsize_ < 2: return name ind = self._indices() @@ -331,8 +331,8 @@ class Param(OptimizationHandlable, ObservableArray, Gradcheckable): if lp is None: lp = self._max_len_names(prirs, __tie_name__) sep = '-' header_format = " {i:{5}^{2}s} | \033[1m{x:{5}^{1}s}\033[0;0m | {c:{5}^{0}s} | {p:{5}^{4}s} | {t:{5}^{3}s}" - if only_name: header = header_format.format(lc, lx, li, lt, lp, ' ', x=self.hirarchy_name(), c=sep*lc, i=sep*li, t=sep*lt, p=sep*lp) # nice header for printing - else: header = header_format.format(lc, lx, li, lt, lp, ' ', x=self.hirarchy_name(), c=__constraints_name__, i=__index_name__, t=__tie_name__, p=__priors_name__) # nice header for printing + if only_name: header = header_format.format(lc, lx, li, lt, lp, ' ', x=self.hierarchy_name(), c=sep*lc, i=sep*li, t=sep*lt, p=sep*lp) # nice header for printing + else: header = header_format.format(lc, lx, li, lt, lp, ' ', x=self.hierarchy_name(), c=__constraints_name__, i=__index_name__, t=__tie_name__, p=__priors_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} | {p:^{5}s} | {t:^{4}s} ".format(lc, lx, __precision__, li, lt, lp, x=x, c=" ".join(map(str, c)), p=" ".join(map(str, p)), t=(t or ''), i=i) for i, x, c, t, p in itertools.izip(indices, vals, constr_matrix, ties, prirs)]) # return all the constraints with right indices # except: return super(Param, self).__str__() @@ -356,6 +356,21 @@ class ParamConcatenation(object): 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:])] + + parents = dict() + for p in self.params: + if p.has_parent(): + parent = p._direct_parent_ + level = 0 + while parent is not None: + if parent in parents: + parents[parent] = max(level, parents[parent]) + else: + parents[parent] = level + level += 1 + parent = parent._direct_parent_ + import operator + self.parents = map(lambda x: x[0], sorted(parents.iteritems(), key=operator.itemgetter(1))) #=========================================================================== # Get/set items, enable broadcasting #=========================================================================== @@ -369,24 +384,26 @@ class ParamConcatenation(object): val = val._vals() ind = numpy.zeros(sum(self._param_sizes), dtype=bool); ind[s] = True; vals = self._vals(); vals[s] = val; del val - [numpy.place(p, ind[ps], vals[ps]) and update and p._notify_observers() + [numpy.place(p, ind[ps], vals[ps]) for p, ps in zip(self.params, self._param_slices_)] + if update: + self.update_all_params() def _vals(self): return numpy.hstack([p._get_params() for p in self.params]) #=========================================================================== # parameter operations: #=========================================================================== def update_all_params(self): - for p in self.params: - p._notify_observers() - + for par in self.parents: + par._notify_observers(-numpy.inf) + def constrain(self, constraint, warning=True): - [param.constrain(constraint, update=False) for param in self.params] + [param.constrain(constraint, trigger_parent=False) for param in self.params] self.update_all_params() constrain.__doc__ = Param.constrain.__doc__ def constrain_positive(self, warning=True): - [param.constrain_positive(warning, update=False) for param in self.params] + [param.constrain_positive(warning, trigger_parent=False) for param in self.params] self.update_all_params() constrain_positive.__doc__ = Param.constrain_positive.__doc__ @@ -396,12 +413,12 @@ class ParamConcatenation(object): fix = constrain_fixed def constrain_negative(self, warning=True): - [param.constrain_negative(warning, update=False) for param in self.params] + [param.constrain_negative(warning, trigger_parent=False) for param in self.params] self.update_all_params() constrain_negative.__doc__ = Param.constrain_negative.__doc__ def constrain_bounded(self, lower, upper, warning=True): - [param.constrain_bounded(lower, upper, warning, update=False) for param in self.params] + [param.constrain_bounded(lower, upper, warning, trigger_parent=False) for param in self.params] self.update_all_params() constrain_bounded.__doc__ = Param.constrain_bounded.__doc__ diff --git a/GPy/core/parameterization/parameter_core.py b/GPy/core/parameterization/parameter_core.py index b4483b4d..4b1b16e0 100644 --- a/GPy/core/parameterization/parameter_core.py +++ b/GPy/core/parameterization/parameter_core.py @@ -6,6 +6,11 @@ import heapq __updated__ = '2013-12-16' +class HierarchyError(Exception): + """ + Gets thrown when something is wrong with the parameter hierarchy + """ + def adjust_name_for_printing(name): if name is not None: return name.replace(" ", "_").replace(".", "_").replace("-", "").replace("+", "").replace("!", "").replace("*", "").replace("/", "") @@ -114,11 +119,11 @@ class Nameable(Parentable): self._name = name if self.has_parent(): self._direct_parent_._name_changed(self, from_name) - def hirarchy_name(self, adjust_for_printing=True): + def hierarchy_name(self, adjust_for_printing=True): if adjust_for_printing: adjust = lambda x: adjust_name_for_printing(x) else: adjust = lambda x: x if self.has_parent(): - return self._direct_parent_.hirarchy_name() + "." + adjust(self.name) + return self._direct_parent_.hierarchy_name() + "." + adjust(self.name) return adjust(self.name) @@ -175,7 +180,7 @@ class Constrainable(Nameable, Indexable): #=========================================================================== # Fixing Parameters: #=========================================================================== - def constrain_fixed(self, value=None, warning=True): + def constrain_fixed(self, value=None, warning=True, trigger_parent=True): """ Constrain this paramter to be fixed to the current value it carries. @@ -183,7 +188,7 @@ class Constrainable(Nameable, Indexable): """ if value is not None: self[:] = value - self.constrain(__fixed__, warning=warning) + self.constrain(__fixed__, warning=warning, trigger_parent=trigger_parent) rav_i = self._highest_parent_._raveled_index_for(self) self._highest_parent_._set_fixed(rav_i) fix = constrain_fixed @@ -224,7 +229,7 @@ class Constrainable(Nameable, Indexable): #=========================================================================== # Prior Operations #=========================================================================== - def set_prior(self, prior, warning=True): + def set_prior(self, prior, warning=True, trigger_parent=True): repriorized = self.unset_priors() self._add_to_index_operations(self.priors, repriorized, prior, warning) @@ -252,7 +257,7 @@ class Constrainable(Nameable, Indexable): # Constrain operations -> done #=========================================================================== - def constrain(self, transform, warning=True): + def constrain(self, transform, warning=True, trigger_parent=True): """ :param transform: the :py:class:`GPy.core.transformations.Transformation` to constrain the this parameter to. @@ -262,7 +267,7 @@ class Constrainable(Nameable, Indexable): :py:class:`GPy.core.transformations.Transformation`. """ if isinstance(transform, Transformation): - self._set_params(transform.initialize(self._get_params()), trigger_parent=True) + self._set_params(transform.initialize(self._get_params()), trigger_parent=trigger_parent) reconstrained = self.unconstrain() self._add_to_index_operations(self.constraints, reconstrained, transform, warning) @@ -275,30 +280,30 @@ class Constrainable(Nameable, Indexable): """ return self._remove_from_index_operations(self.constraints, transforms) - def constrain_positive(self, warning=True): + def constrain_positive(self, warning=True, trigger_parent=True): """ :param warning: print a warning if re-constraining parameters. Constrain this parameter to the default positive constraint. """ - self.constrain(Logexp(), warning=warning) + self.constrain(Logexp(), warning=warning, trigger_parent=trigger_parent) - def constrain_negative(self, warning=True): + def constrain_negative(self, warning=True, trigger_parent=True): """ :param warning: print a warning if re-constraining parameters. Constrain this parameter to the default negative constraint. """ - self.constrain(NegativeLogexp(), warning=warning) + self.constrain(NegativeLogexp(), warning=warning, trigger_parent=trigger_parent) - def constrain_bounded(self, lower, upper, warning=True): + def constrain_bounded(self, lower, upper, warning=True, trigger_parent=True): """ :param lower, upper: the limits to bound this parameter to :param warning: print a warning if re-constraining parameters. Constrain this parameter to lie within the given range. """ - self.constrain(Logistic(lower, upper), warning=warning) + self.constrain(Logistic(lower, upper), warning=warning, trigger_parent=trigger_parent) def unconstrain_positive(self): """ @@ -359,6 +364,9 @@ class OptimizationHandlable(Constrainable, Observable): # inverse apply transformations for parameters and set the resulting parameters self._set_params(self._untransform_params(p)) + def _size_transformed(self): + return self.size - self.constraints[__fixed__].size + def _untransform_params(self, p): p = p.copy() if self._has_fixes(): tmp = self._get_params(); tmp[self._fixes_] = p; p = tmp; del tmp @@ -373,13 +381,13 @@ class OptimizationHandlable(Constrainable, Observable): def _set_params(self, params, trigger_parent=True): # don't overwrite this anymore! - raise NotImplementedError, "This needs to be implemented seperately" + raise NotImplementedError, "This needs to be implemented in Param and Parametrizable" #=========================================================================== # Optimization handles: #=========================================================================== def _get_param_names(self): - n = np.array([p.hirarchy_name() + '[' + str(i) + ']' for p in self.flattened_parameters for i in p._indices()]) + n = np.array([p.hierarchy_name() + '[' + str(i) + ']' for p in self.flattened_parameters for i in p._indices()]) return n def _get_param_names_transformed(self): n = self._get_param_names() @@ -398,7 +406,7 @@ class OptimizationHandlable(Constrainable, Observable): import numpy as np # first take care of all parameters (from N(0,1)) # x = self._get_params_transformed() - x = np.random.randn(self.size_transformed) + x = np.random.randn(self._size_transformed()) x = self._untransform_params(x) # now draw from prior where possible [np.put(x, ind, p.rvs(ind.size)) for p, ind in self.priors.iteritems() if not p is None] @@ -435,7 +443,7 @@ class Parameterizable(OptimizationHandlable): if pname in self._added_names_: del self.__dict__[pname] self._add_parameter_name(param) - else: + elif pname not in dir(self): self.__dict__[pname] = param self._added_names_.add(pname) diff --git a/GPy/core/parameterization/parameterized.py b/GPy/core/parameterization/parameterized.py index 0093c6f3..56d785c3 100644 --- a/GPy/core/parameterization/parameterized.py +++ b/GPy/core/parameterization/parameterized.py @@ -7,7 +7,7 @@ import cPickle import itertools from re import compile, _pattern_type from param import ParamConcatenation -from parameter_core import Constrainable, Pickleable, Parentable, Observable, Parameterizable, adjust_name_for_printing, Gradcheckable +from parameter_core import Pickleable, Parameterizable, adjust_name_for_printing, Gradcheckable from transformations import __fixed__ from array_core import ParamList @@ -105,6 +105,14 @@ class Parameterized(Parameterizable, Pickleable, Gradcheckable): self.remove_parameter(param) self.add_parameter(param, index) elif param not in self._parameters_: + if param.has_parent(): + parent = param._direct_parent_ + while parent is not None: + if parent is self: + from GPy.core.parameterization.parameter_core import HierarchyError + raise HierarchyError, "You cannot add a parameter twice into the hirarchy" + parent = parent._direct_parent_ + param._direct_parent_.remove_parameter(param) # make sure the size is set if index is None: self.constraints.update(param.constraints, self.size) @@ -117,13 +125,16 @@ class Parameterized(Parameterizable, Pickleable, Gradcheckable): self.constraints.update(param.constraints, start) self.priors.update(param.priors, start) self._parameters_.insert(index, param) + param.add_observer(self, self._pass_through_notify_observers, -np.inf) + self.size += param.size + + self._connect_parameters() + self._notify_parent_change() + self._connect_fixes() else: raise RuntimeError, """Parameter exists already added and no copy made""" - self._connect_parameters() - self._notify_parent_change() - self._connect_fixes() def add_parameters(self, *parameters): @@ -146,12 +157,19 @@ class Parameterized(Parameterizable, Pickleable, Gradcheckable): del self._parameters_[param._parent_index_] param._disconnect_parent() - param.remove_observer(self, self._notify_parameters_changed) + param.remove_observer(self, self._pass_through_notify_observers) self.constraints.shift_left(start, param.size) + self._connect_fixes() self._connect_parameters() self._notify_parent_change() + parent = self._direct_parent_ + while parent is not None: + parent._connect_fixes() + parent._connect_parameters() + parent._notify_parent_change() + parent = parent._direct_parent_ def _connect_parameters(self): # connect parameterlist to this parameterized object @@ -351,7 +369,7 @@ class Parameterized(Parameterizable, Pickleable, Gradcheckable): # Printing: #=========================================================================== def _short(self): - return self.hirarchy_name() + return self.hierarchy_name() @property def flattened_parameters(self): return [xi for x in self._parameters_ for xi in x.flattened_parameters] @@ -359,11 +377,6 @@ class Parameterized(Parameterizable, Pickleable, Gradcheckable): def _parameter_sizes_(self): return [x.size for x in self._parameters_] @property - def size_transformed(self): - if self._has_fixes(): - return sum(self._fixes_) - return self.size - @property def parameter_shapes(self): return [xi for x in self._parameters_ for xi in x.parameter_shapes] @property diff --git a/GPy/examples/dimensionality_reduction.py b/GPy/examples/dimensionality_reduction.py index 2044f08d..9ebb54a2 100644 --- a/GPy/examples/dimensionality_reduction.py +++ b/GPy/examples/dimensionality_reduction.py @@ -187,10 +187,10 @@ def _simulate_sincos(D1, D2, D3, N, num_inducing, Q, plot_sim=False): _np.random.seed(1234) x = _np.linspace(0, 4 * _np.pi, N)[:, None] - s1 = _np.vectorize(lambda x: -_np.sin(_np.exp(x))) + s1 = _np.vectorize(lambda x: _np.sin(x)) s2 = _np.vectorize(lambda x: _np.cos(x)**2) s3 = _np.vectorize(lambda x:-_np.exp(-_np.cos(2 * x))) - sS = _np.vectorize(lambda x: x*_np.sin(x)) + sS = _np.vectorize(lambda x: _np.cos(x)) s1 = s1(x) s2 = s2(x) @@ -202,7 +202,7 @@ def _simulate_sincos(D1, D2, D3, N, num_inducing, Q, plot_sim=False): s3 -= s3.mean(); s3 /= s3.std(0) sS -= sS.mean(); sS /= sS.std(0) - S1 = _np.hstack([s1, s2, sS]) + S1 = _np.hstack([s1, sS]) S2 = _np.hstack([s2, s3, sS]) S3 = _np.hstack([s3, sS]) @@ -270,7 +270,7 @@ def bgplvm_simulation(optimize=True, verbose=1, from GPy import kern from GPy.models import BayesianGPLVM - D1, D2, D3, N, num_inducing, Q = 13, 5, 8, 45, 5, 9 + D1, D2, D3, N, num_inducing, Q = 13, 5, 8, 45, 3, 9 _, _, Ylist = _simulate_sincos(D1, D2, D3, N, num_inducing, Q, plot_sim) Y = Ylist[0] k = kern.Linear(Q, ARD=True)# + kern.white(Q, _np.exp(-2)) # + kern.bias(Q) @@ -294,7 +294,7 @@ def bgplvm_simulation_missing_data(optimize=True, verbose=1, from GPy.models import BayesianGPLVM from GPy.inference.latent_function_inference.var_dtc import VarDTCMissingData - D1, D2, D3, N, num_inducing, Q = 13, 5, 8, 45, 5, 9 + D1, D2, D3, N, num_inducing, Q = 13, 5, 8, 45, 7, 9 _, _, Ylist = _simulate_sincos(D1, D2, D3, N, num_inducing, Q, plot_sim) Y = Ylist[0] k = kern.Linear(Q, ARD=True)# + kern.white(Q, _np.exp(-2)) # + kern.bias(Q) diff --git a/GPy/inference/latent_function_inference/var_dtc.py b/GPy/inference/latent_function_inference/var_dtc.py index 2b7ca7ad..64707298 100644 --- a/GPy/inference/latent_function_inference/var_dtc.py +++ b/GPy/inference/latent_function_inference/var_dtc.py @@ -60,8 +60,7 @@ class VarDTC(object): _, output_dim = Y.shape #see whether we've got a different noise variance for each datum - beta = 1./max(1e-6, np.squeeze(likelihood.variance)) - + beta = 1./np.fmax(likelihood.variance, 1e-6) # VVT_factor is a matrix such that tdot(VVT_factor) = VVT...this is for efficiency! #self.YYTfactor = self.get_YYTfactor(Y) #VVT_factor = self.get_VVTfactor(self.YYTfactor, beta) @@ -214,7 +213,7 @@ class VarDTCMissingData(object): psi2_all = None Ys, traces = self._Y(Y) - beta_all = 1./max(1e-6, likelihood.variance) + beta_all = 1./np.fmax(likelihood.variance, 1e-6) het_noise = beta_all.size != 1 import itertools diff --git a/GPy/kern/_src/kern.py b/GPy/kern/_src/kern.py index eb3291e0..14e6ae49 100644 --- a/GPy/kern/_src/kern.py +++ b/GPy/kern/_src/kern.py @@ -112,10 +112,12 @@ class Kern(Parameterized): """ assert isinstance(other, Kern), "only kernels can be added to kernels..." from add import Add - return Add([self, other], tensor) - - def __call__(self, X, X2=None): - return self.K(X, X2) + kernels = [] + if not tensor and isinstance(self, Add): kernels.extend(self._parameters_) + else: kernels.append(self) + if not tensor and isinstance(other, Add): kernels.extend(other._parameters_) + else: kernels.append(other) + return Add(kernels, tensor) def __mul__(self, other): """ Here we overload the '*' operator. See self.prod for more information""" diff --git a/GPy/models/sparse_gp_regression.py b/GPy/models/sparse_gp_regression.py index 6a76df3f..4980d26a 100644 --- a/GPy/models/sparse_gp_regression.py +++ b/GPy/models/sparse_gp_regression.py @@ -8,7 +8,7 @@ from .. import likelihoods from .. import kern from ..inference.latent_function_inference import VarDTC from ..util.misc import param_to_array -from ..core.parameterization.variational import VariationalPosterior +from ..core.parameterization.variational import NormalPosterior class SparseGPRegression(SparseGP): """ @@ -47,7 +47,7 @@ class SparseGPRegression(SparseGP): likelihood = likelihoods.Gaussian() if not (X_variance is None): - X = VariationalPosterior(X,X_variance) + X = NormalPosterior(X,X_variance) SparseGP.__init__(self, X, Y, Z, kernel, likelihood, inference_method=VarDTC()) diff --git a/GPy/testing/parameterized_tests.py b/GPy/testing/parameterized_tests.py index 6f13d294..0da3f3ae 100644 --- a/GPy/testing/parameterized_tests.py +++ b/GPy/testing/parameterized_tests.py @@ -6,6 +6,7 @@ Created on Feb 13, 2014 import unittest import GPy import numpy as np +from GPy.core.parameterization.parameter_core import HierarchyError class Test(unittest.TestCase): @@ -65,7 +66,7 @@ class Test(unittest.TestCase): self.assertListEqual(self.test1.constraints[Logexp()].tolist(), [0,1]) def test_add_parameter_already_in_hirarchy(self): - self.test1.add_parameter(self.white._parameters_[0]) + self.assertRaises(HierarchyError, self.test1.add_parameter, self.white._parameters_[0]) def test_default_constraints(self): self.assertIs(self.rbf.variance.constraints._param_index_ops, self.rbf.constraints._param_index_ops) diff --git a/GPy/util/caching.py b/GPy/util/caching.py index 250efe11..a2434c80 100644 --- a/GPy/util/caching.py +++ b/GPy/util/caching.py @@ -40,7 +40,7 @@ class Cacher(object): return self.operation(*args) # TODO: WARNING !!! Cache OFFSWITCH !!! WARNING - # return self.operation(*args) + return self.operation(*args) #if the result is cached, return the cached computation state = [all(a is b for a, b in itertools.izip_longest(args, cached_i)) for cached_i in self.cached_inputs]