diff --git a/GPy/core/parameterization/parameter_core.py b/GPy/core/parameterization/parameter_core.py index b99ae25a..e0dadf18 100644 --- a/GPy/core/parameterization/parameter_core.py +++ b/GPy/core/parameterization/parameter_core.py @@ -500,8 +500,9 @@ class Indexable(Nameable, Observable): #=========================================================================== def tie(self, name): + from ties_and_remappings import Tie #remove any constraints - old_const = self.constraints.properties()[:] + old_const = [c for c in self.constraints.properties() if not isinstance(c,Tie)] self.unconstrain() #see if a tie exists with that name @@ -510,14 +511,14 @@ class Indexable(Nameable, Observable): else: #create a tie object value = np.atleast_1d(self.param_array)[0]*1 - from ties_and_remappings import Tie t = Tie(value=value, name=name) #add the new tie object to the global index self._highest_parent_.ties[name] = t self._highest_parent_.add_parameter(t) + #constrain the tie as we were constrained - if len(old_const)==1: + if len(old_const)>0: t.constrain(old_const[0]) self.constraints.add(t, self._raveled_index()) diff --git a/GPy/core/parameterization/ties_and_remappings.py b/GPy/core/parameterization/ties_and_remappings.py index a4cf4e6c..61c40a8e 100644 --- a/GPy/core/parameterization/ties_and_remappings.py +++ b/GPy/core/parameterization/ties_and_remappings.py @@ -46,7 +46,7 @@ class Tie(Remapping): def callback(self, param=None, which=None): """ This gets called whenever any of the tied parameters changes. we spend - considerable effort working out whhat has changed ant to what value. + considerable effort working out what has changed and to what value. Then we store that value in self.value, and broadcast it everywhere with parameters_changed. """ @@ -54,11 +54,13 @@ class Tie(Remapping): index = self._highest_parent_.constraints[self] if len(index)==0: return # nothing to tie together, this tie exists without any tied parameters - self.value.gradient[:] = self._highest_parent_.gradient[index].sum() + self.collate_gradient() vals = self._highest_parent_.param_array[index] uvals = np.unique(vals) if len(uvals)==1: #all of the tied things are at the same value + if (self.value==uvals[0]).all(): + return # DO NOT DO ANY CHANGES IF THE TIED PART IS NOT CHANGED! self.value[...] = uvals[0] elif len(uvals)==2: #only *one* of the tied things has changed. it must be different to self.value @@ -69,7 +71,7 @@ class Tie(Remapping): raise ValueError, "something is wrong with the tieing" def parameters_changed(self): super(Tie,self).parameters_changed() - self.value.gradient[:] = self._highest_parent_.gradient[self._highest_parent_.constraints[self]].sum() + self.collate_gradient() def mapping(self): return self.value diff --git a/GPy/core/parameterization/variational.py b/GPy/core/parameterization/variational.py index b9811af8..e70b9b71 100644 --- a/GPy/core/parameterization/variational.py +++ b/GPy/core/parameterization/variational.py @@ -40,7 +40,6 @@ class SpikeAndSlabPrior(VariationalPrior): self.pi = Param('pi', pi, Logistic(1e-10,1.-1e-10)) self.variance = Param('variance',variance) self.add_parameters(self.pi) - self.group_spike_prob = False def KL_divergence(self, variational_posterior): mu = variational_posterior.mean @@ -56,11 +55,7 @@ class SpikeAndSlabPrior(VariationalPrior): S = variational_posterior.variance gamma = variational_posterior.binary_prob - if self.group_spike_prob: - gamma_grad = np.log((1-self.pi)/self.pi*gamma/(1.-gamma))+(np.square(mu)+S-np.log(S)-1.)/2. - gamma.gradient -= gamma_grad.mean(axis=0) - else: - gamma.gradient -= np.log((1-self.pi)/self.pi*gamma/(1.-gamma))+(np.square(mu)+S-np.log(S)-1.)/2. + gamma.gradient -= np.log((1-self.pi)/self.pi*gamma/(1.-gamma))+(np.square(mu)+S-np.log(S)-1.)/2. mu.gradient -= gamma*mu S.gradient -= (1. - (1. / (S))) * gamma /2. self.pi.gradient = (gamma/self.pi - (1.-gamma)/(1.-self.pi)).sum(axis=0) diff --git a/GPy/kern/_src/psi_comp/rbf_psi_gpucomp.py b/GPy/kern/_src/psi_comp/rbf_psi_gpucomp.py index 43fffcd1..92387945 100644 --- a/GPy/kern/_src/psi_comp/rbf_psi_gpucomp.py +++ b/GPy/kern/_src/psi_comp/rbf_psi_gpucomp.py @@ -242,14 +242,14 @@ gpu_code = """ class PSICOMP_RBF_GPU(PSICOMP_RBF): - def __init__(self, GPU_direct=False): + def __init__(self, threadnum=128, blocknum=15, GPU_direct=False): assert gpu_init.initSuccess, "GPU initialization failed!" self.GPU_direct = GPU_direct self.cublas_handle = gpu_init.cublas_handle self.gpuCache = None - self.threadnum = 128 - self.blocknum = 15 + self.threadnum = threadnum + self.blocknum = blocknum module = SourceModule("#define THREADNUM "+str(self.threadnum)+"\n"+gpu_code) self.g_psi1computations = module.get_function('psi1computations') self.g_psi1computations.prepare('PPdPPPPiii') diff --git a/GPy/models/ss_gplvm.py b/GPy/models/ss_gplvm.py index ab32e462..fa469b3b 100644 --- a/GPy/models/ss_gplvm.py +++ b/GPy/models/ss_gplvm.py @@ -44,15 +44,12 @@ class SSGPLVM(SparseGP): if X_variance is None: # The variance of the variational approximation (S) X_variance = np.random.uniform(0,.1,X.shape) - gamma = np.empty_like(X, order='F') # The posterior probabilities of the binary variable in the variational approximation + gamma = np.empty_like(X) # The posterior probabilities of the binary variable in the variational approximation gamma[:] = 0.5 + 0.1 * np.random.randn(X.shape[0], input_dim) gamma[gamma>1.-1e-9] = 1.-1e-9 gamma[gamma<1e-9] = 1e-9 gamma[:] = 0.5 - - if group_spike: - gamma[:] = gamma[:,0] - + if Z is None: Z = np.random.permutation(X.copy())[:num_inducing] assert Z.shape[1] == X.shape[1] @@ -73,14 +70,13 @@ class SSGPLVM(SparseGP): X = SpikeAndSlabPosterior(X, X_variance, gamma) - if group_spike: - kernel.group_spike_prob = True - self.variational_prior.group_spike_prob = True - SparseGP.__init__(self, X, Y, Z, kernel, likelihood, inference_method, name, **kwargs) self.add_parameter(self.X, index=0) self.add_parameter(self.variational_prior) + if self.group_spike: + [self.X.gamma[:,i].tie('tieGamma'+str(i)) for i in xrange(self.X.gamma.shape[1])] # Tie columns together + if mpi_comm != None: from ..util.mpi import divide_data Y_start, Y_end, Y_list = divide_data(Y.shape[0], mpi_comm)