From 5214c3c1ac46d60b1818e614394e10106e117bc8 Mon Sep 17 00:00:00 2001 From: Neil Lawrence Date: Wed, 19 Feb 2014 19:39:24 -0500 Subject: [PATCH] Adding update_gradients to sympy.py. --- GPy/kern/parts/rbf.py | 2 +- GPy/kern/parts/sympykern.py | 199 +++++++++++++++++++++--------------- 2 files changed, 115 insertions(+), 86 deletions(-) diff --git a/GPy/kern/parts/rbf.py b/GPy/kern/parts/rbf.py index 027aa382..8811b74a 100644 --- a/GPy/kern/parts/rbf.py +++ b/GPy/kern/parts/rbf.py @@ -109,7 +109,7 @@ class RBF(Kernpart): self.lengthscale.gradient = self._dL_dlengthscales_via_K(dL_dK, X, None) else: self.lengthscale.gradient = (self.variance / self.lengthscale) * np.sum(self._K_dvar * self._K_dist2 * dL_dK) - +b def update_gradients_sparse(self, dL_dKmm, dL_dKnm, dL_dKdiag, X, Z): #contributions from Kdiag self.variance.gradient = np.sum(dL_dKdiag) diff --git a/GPy/kern/parts/sympykern.py b/GPy/kern/parts/sympykern.py index a5bb7b1d..52813ecd 100644 --- a/GPy/kern/parts/sympykern.py +++ b/GPy/kern/parts/sympykern.py @@ -26,6 +26,8 @@ import ast from kernpart import Kernpart from ...core.parameterization import Param from ...core.parameterization.transformations import Logexp +# TODO have this set up in a set up file! +user_code_storage = tempfile.gettempdir() class spkern(Kernpart): """ @@ -61,13 +63,12 @@ class spkern(Kernpart): assert all([x.name=='x_%i'%i for i,x in enumerate(self._sp_x)]) assert all([z.name=='z_%i'%i for i,z in enumerate(self._sp_z)]) assert len(self._sp_x)==len(self._sp_z) - assert len(self._sp_x)==input_dim + x_dim=len(self._sp_x) # If it is a multi-output covariance, add an input for indexing the outputs. - self._real_input_dim = self.input_dim - if output_dim > 1: - self.input_dim += 1 - assert self.input_dim == input_dim + self._real_input_dim = x_dim + # Check input dim is number of xs + 1 if output_dim is >1 + assert self.input_dim == x_dim + int(output_dim > 1) self.output_dim = output_dim # extract parameter names from the covariance @@ -113,7 +114,6 @@ class spkern(Kernpart): #setattr(self, theta.name, val) setattr(self, theta.name, Param(theta.name, val, None)) self.add_parameters(getattr(self, theta.name)) - self.parameters_changed() # initializes cache #deal with param #self._set_params(self._get_params()) @@ -139,13 +139,15 @@ class spkern(Kernpart): extra_compile_args = [] self.weave_kwargs = { - 'support_code':self._function_code, - 'include_dirs':[tempfile.gettempdir(), os.path.join(current_dir,'parts/')], - 'headers':['"sympy_helpers.h"'], - 'sources':[os.path.join(current_dir,"parts/sympy_helpers.cpp")], + 'support_code': None, #self._function_code, + 'include_dirs':[user_code_storage, os.path.join(current_dir,'parts/')], + 'headers':['"sympy_helpers.h"', '"'+self.name+'.h"'], + 'sources':[os.path.join(current_dir,"parts/sympy_helpers.cpp"), os.path.join(user_code_storage, self.name+'.cpp')], 'extra_compile_args':extra_compile_args, 'extra_link_args':['-lgomp'], 'verbose':True} + self.parameters_changed() # initializes caches + def __add__(self,other): return spkern(self._sp_k+other._sp_k) @@ -177,31 +179,39 @@ class spkern(Kernpart): # Use weave to compute the underlying functions. if weave_available: # put the header file where we can find it - f = file(os.path.join(tempfile.gettempdir(), self.name + '.h'),'w') + f = file(os.path.join(user_code_storage, self.name + '.h'),'w') f.write(self._function_header) f.close() - - # Substitute any known derivatives which sympy doesn't compute - self._function_code = re.sub('DiracDelta\(.+?,.+?\)','0.0',self._function_code) + + if weave_available: + # Substitute any known derivatives which sympy doesn't compute + self._function_code = re.sub('DiracDelta\(.+?,.+?\)','0.0',self._function_code) + # put the cpp file in user code storage (defaults to temp file location) + f = file(os.path.join(user_code_storage, self.name + '.cpp'),'w') + else: + # put the python file in user code storage + f = file(os.path.join(user_code_storage, self.name + '.py'),'w') + f.write(self._function_code) + f.close() if weave_available: # arg_list will store the arguments required for the C code. - arg_list = (["X2(i, %s)"%x.name[2:] for x in self._sp_x] + input_arg_list = (["X2(i, %s)"%x.name[2:] for x in self._sp_x] + ["Z2(j, %s)"%z.name[2:] for z in self._sp_z]) # for multiple outputs reverse argument list is also required if self.output_dim>1: - reverse_arg_list = list(arg_list) - reverse_arg_list.reverse() + reverse_input_arg_list = list(input_arg_list) + reverse_input_arg_list.reverse() # This gives the parameters for the arg list. param_arg_list = [shared_params.name for shared_params in self._sp_theta] - arg_list += param_arg_list + arg_list = input_arg_list + param_arg_list precompute_list=[] if self.output_dim > 1: - reverse_arg_list+=list(param_arg_list) + reverse_arg_list= reverse_input_arg_list + list(param_arg_list) # For multiple outputs, also need the split parameters. split_param_arg_list = ["%s1(%s)"%(theta.name[:-2].upper(),index) for index in ['ii', 'jj'] for theta in self._sp_theta_i] split_param_reverse_arg_list = ["%s1(%s)"%(theta.name[:-2].upper(),index) for index in ['jj', 'ii'] for theta in self._sp_theta_i] @@ -218,9 +228,9 @@ class spkern(Kernpart): # Any precomputations will be done here eventually. self._precompute = \ - """ - // Precompute code would go here. It will be called when parameters are updated. - """ + """ + // Precompute code would go here. It will be called when parameters are updated. + """ # Here's the code to do the looping for K self._K_code =\ @@ -229,11 +239,11 @@ class spkern(Kernpart): // Code for computing the covariance function. int i; int j; - int N = target_array->dimensions[0]; + int n = target_array->dimensions[0]; int num_inducing = target_array->dimensions[1]; int input_dim = X_array->dimensions[1]; //#pragma omp parallel for private(j) - for (i=0;idimensions[0]; + int n = target_array->dimensions[0]; int input_dim = X_array->dimensions[1]; //#pragma omp parallel for - for (i=0;i1: grad_func_list += c_define_output_indices - grad_func_list += [' '*16 + 'TARGET1(%i+ii) += partial[i*num_inducing+j]*dk_d%s(%s);'%(self.num_shared_params+i*self.output_dim, theta.name, arg_string) for i, theta in enumerate(self._sp_theta_i)] - grad_func_list += [' '*16 + 'TARGET1(%i+jj) += partial[i*num_inducing+j]*dk_d%s(%s);'%(self.num_shared_params+i*self.output_dim, theta.name, reverse_arg_string) for i, theta in enumerate(self._sp_theta_i)] - grad_func_list += ([' '*16 + 'TARGET1(%i) += partial[i*num_inducing+j]*dk_d%s(%s);'%(i,theta.name,arg_string) for i,theta in enumerate(self._sp_theta)]) + grad_func_list += [' '*16 + 'TARGET1(%i+ii) += PARTIAL2(i, j)*dk_d%s(%s);'%(self.num_shared_params+i*self.output_dim, theta.name, arg_string) for i, theta in enumerate(self._sp_theta_i)] + grad_func_list += [' '*16 + 'TARGET1(%i+jj) += PARTIAL2(i, j)*dk_d%s(%s);'%(self.num_shared_params+i*self.output_dim, theta.name, reverse_arg_string) for i, theta in enumerate(self._sp_theta_i)] + grad_func_list += ([' '*16 + 'TARGET1(%i) += PARTIAL2(i, j)*dk_d%s(%s);'%(i,theta.name,arg_string) for i,theta in enumerate(self._sp_theta)]) grad_func_string = '\n'.join(grad_func_list) self._dK_dtheta_code =\ @@ -285,17 +296,18 @@ class spkern(Kernpart): // Code for computing gradient of covariance with respect to parameters. int i; int j; - int N = partial_array->dimensions[0]; + int n = partial_array->dimensions[0]; int num_inducing = partial_array->dimensions[1]; int input_dim = X_array->dimensions[1]; //#pragma omp parallel for private(j) - for (i=0;idimensions[0]; + int n = partial_array->dimensions[0]; int input_dim = X_array->dimensions[1]; - for (i=0;idimensions[0]; + int n = partial_array->dimensions[0]; int num_inducing = partial_array->dimensions[1]; int input_dim = X_array->dimensions[1]; //#pragma omp parallel for private(j) - for (i=0;idimensions[0]; + int n = partial_array->dimensions[0]; int input_dim = X_array->dimensions[1]; - for (int i=0;i1: - # for i, split_params in enumerate(self._split_theta_names): - # start = self.num_shared_params + i*self.output_dim - # end = self.num_shared_params + (i+1)*self.output_dim - # setattr(self, split_params, param[start:end]) + #---------------------------------------# + # Precomputations # + #---------------------------------------# - - # def _get_params(self): - # params = np.zeros(0) - # for shared_params in self._sp_theta: - # params = np.hstack((params, getattr(self, shared_params.name))) - # if self.output_dim>1: - # for split_params in self._split_theta_names: - # params = np.hstack((params, getattr(self, split_params).flatten())) - # return params - - # def _get_param_names(self): - # if self.output_dim>1: - # return [x.name for x in self._sp_theta] + [x.name[:-2] + str(i) for x in self._sp_theta_i for i in range(self.output_dim)] - # else: - # return [x.name for x in self._sp_theta] + def _K_computations(self, X, Z): + if Z is None: + self._generate_inline(self._precompute, X) + else: + self._generate_inline(self._precompute, X, Z=Z)