From 11d088cf903191254b53a47a6647e5512a8511ec Mon Sep 17 00:00:00 2001 From: Nicolas Date: Fri, 18 Jan 2013 15:14:23 +0000 Subject: [PATCH] added ARD flag to exponential --- GPy/kern/__init__.py | 2 +- GPy/kern/constructors.py | 22 +++------------- GPy/kern/exponential.py | 54 +++++++++++++++++++++++++--------------- 3 files changed, 39 insertions(+), 39 deletions(-) diff --git a/GPy/kern/__init__.py b/GPy/kern/__init__.py index cd893bac..4a36d6d0 100644 --- a/GPy/kern/__init__.py +++ b/GPy/kern/__init__.py @@ -2,5 +2,5 @@ # Licensed under the BSD 3-clause license (see LICENSE.txt) -from constructors import rbf, Matern32, Matern52, exponential, linear, white, bias, finite_dimensional, rbf_ARD, spline, Brownian, linear_ARD, rbf_sympy, sympykern +from constructors import rbf, Matern32, Matern52, exponential, linear, white, bias, finite_dimensional, spline, Brownian, linear_ARD, rbf_sympy, sympykern from kern import kern diff --git a/GPy/kern/constructors.py b/GPy/kern/constructors.py index e1304f11..5f676d9b 100644 --- a/GPy/kern/constructors.py +++ b/GPy/kern/constructors.py @@ -36,20 +36,6 @@ def rbf(D,variance=1., lengthscale=None,ARD=False): part = rbfpart(D,variance,lengthscale,ARD) return kern(D, [part]) -def rbf_ARD(D,variance=1., lengthscales=None): - """ - Construct an RBF kernel with Automatic Relevance Determination (ARD) - - :param D: dimensionality of the kernel, obligatory - :type D: int - :param variance: the variance of the kernel - :type variance: float - :param lengthscales: the lengthscales of the kernel - :type lengthscales: None|np.ndarray - """ - part = rbf_ARD_part(D,variance,lengthscales) - return kern(D, [part]) - def linear(D,lengthscales=None): """ Construct a linear kernel. @@ -86,7 +72,7 @@ def white(D,variance=1.): part = whitepart(D,variance) return kern(D, [part]) -def exponential(D,variance=1., lengthscales=None): +def exponential(D,variance=1., lengthscale=None, ARD=False): """ Construct a exponential kernel. @@ -96,10 +82,10 @@ def exponential(D,variance=1., lengthscales=None): variance (float) lengthscales (np.ndarray) """ - part = exponentialpart(D,variance, lengthscales) + part = exponentialpart(D,variance, lengthscale, ARD) return kern(D, [part]) -def Matern32(D,variance=1., lengthscales=None): +def Matern32(D,variance=1., lengthscale=None, ARD=False): """ Construct a Matern 3/2 kernel. @@ -109,7 +95,7 @@ def Matern32(D,variance=1., lengthscales=None): variance (float) lengthscales (np.ndarray) """ - part = Matern32part(D,variance, lengthscales) + part = Matern32part(D,variance, lengthscale, ARD) return kern(D, [part]) def Matern52(D,variance=1., lengthscales=None): diff --git a/GPy/kern/exponential.py b/GPy/kern/exponential.py index 2df6a958..21bb8398 100644 --- a/GPy/kern/exponential.py +++ b/GPy/kern/exponential.py @@ -24,37 +24,46 @@ class exponential(kernpart): :rtype: kernel object """ - def __init__(self,D,variance=1.,lengthscales=None): + def __init__(self,D,variance=1.,lengthscale=None,ARD=False): self.D = D - if lengthscales is not None: - assert lengthscales.shape==(self.D,) + self.ARD = ARD + if ARD == False: + self.Nparam = 2 + self.name = 'exp' + if lengthscale is not None: + assert lengthscale.shape == (1,) + else: + lengthscale = np.ones(1) else: - lengthscales = np.ones(self.D) - self.Nparam = self.D + 1 - self.name = 'exp' - self._set_params(np.hstack((variance,lengthscales))) + self.Nparam = self.D + 1 + self.name = 'exp_ARD' + if lengthscale is not None: + assert lengthscale.shape == (self.D,) + else: + lengthscale = np.ones(self.D) + self._set_params(np.hstack((variance,lengthscale))) def _get_params(self): """return the value of the parameters.""" - return np.hstack((self.variance,self.lengthscales)) + return np.hstack((self.variance,self.lengthscale)) def _set_params(self,x): """set the value of the parameters.""" assert x.size==(self.D+1) self.variance = x[0] - self.lengthscales = x[1:] + self.lengthscale = x[1:] def _get_param_names(self): """return parameter names.""" - if self.D==1: + if self.Nparam==2: return ['variance','lengthscale'] else: - return ['variance']+['lengthscale_%i'%i for i in range(self.lengthscales.size)] + return ['variance']+['lengthscale_%i'%i for i in range(self.lengthscale.size)] def K(self,X,X2,target): """Compute the covariance matrix between X and X2.""" if X2 is None: X2 = X - dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscales),-1)) + dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1)) np.add(self.variance*np.exp(-dist), target,target) def Kdiag(self,X,target): @@ -64,13 +73,18 @@ class exponential(kernpart): def dK_dtheta(self,partial,X,X2,target): """derivative of the covariance matrix with respect to the parameters.""" if X2 is None: X2 = X - dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscales),-1)) + dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1)) invdist = 1./np.where(dist!=0.,dist,np.inf) - dist2M = np.square(X[:,None,:]-X2[None,:,:])/self.lengthscales**3 + dist2M = np.square(X[:,None,:]-X2[None,:,:])/self.lengthscale**3 dvar = np.exp(-dist) - dl = self.variance*dvar[:,:,None]*dist2M*invdist[:,:,None] target[0] += np.sum(dvar*partial) - target[1:] += (dl*partial[:,:,None]).sum(0).sum(0) + if self.ARD == True: + dl = self.variance*dvar[:,:,None]*dist2M*invdist[:,:,None] + target[1:] += (dl*partial[:,:,None]).sum(0).sum(0) + else: + dl = self.variance*dvar*dist2M.sum(-1)*invdist + target[1] += np.sum(dl*partial) + #foo def dKdiag_dtheta(self,partial,X,target): """derivative of the diagonal of the covariance matrix with respect to the parameters.""" @@ -80,8 +94,8 @@ class exponential(kernpart): def dK_dX(self,partial,X,X2,target): """derivative of the covariance matrix with respect to X.""" if X2 is None: X2 = X - dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscales),-1))[:,:,None] - ddist_dX = (X[:,None,:]-X2[None,:,:])/self.lengthscales**2/np.where(dist!=0.,dist,np.inf) + dist = np.sqrt(np.sum(np.square((X[:,None,:]-X2[None,:,:])/self.lengthscale),-1))[:,:,None] + ddist_dX = (X[:,None,:]-X2[None,:,:])/self.lengthscale**2/np.where(dist!=0.,dist,np.inf) dK_dX = - np.transpose(self.variance*np.exp(-dist)*ddist_dX,(1,0,2)) target += np.sum(dK_dX*partial.T[:,:,None],0) @@ -101,14 +115,14 @@ class exponential(kernpart): """ assert self.D == 1 def L(x,i): - return(1./self.lengthscales*F[i](x) + F1[i](x)) + return(1./self.lengthscale*F[i](x) + F1[i](x)) n = F.shape[0] G = np.zeros((n,n)) for i in range(n): for j in range(i,n): G[i,j] = G[j,i] = integrate.quad(lambda x : L(x,i)*L(x,j),lower,upper)[0] Flower = np.array([f(lower) for f in F])[:,None] - return(self.lengthscales/2./self.variance * G + 1./self.variance * np.dot(Flower,Flower.T)) + return(self.lengthscale/2./self.variance * G + 1./self.variance * np.dot(Flower,Flower.T))