# Copyright (c) 2012, GPy authors (see AUTHORS.txt). # Licensed under the BSD 3-clause license (see LICENSE.txt) from kernpart import Kernpart import numpy as np class PeriodicSpikes(Kernpart): """ periodic spikes kernel .. math:: k(r) = cos(sin(2*pi*r/p) * exp( s*cos(2*pi*r/p) - s) where r^2 = (x-y)^2 :param input_dim: the number of input dimensions :type input_dim: int (input_dim=1 is the only value currently supported) :param shape: the shape parameter s, should be >1 :type shape: float :param period: the period p :type period: float :rtype: Kernpart object """ def __init__(self,input_dim,shape=1.,period=1.): assert input_dim == 1, "For this kernel we assume input_dim=1" self.input_dim = input_dim self.num_params = 2 self.name = 'ps' self.shape = shape self.period = period def _get_params(self): return np.hstack((self.shape,self.period)) def _set_params(self,x): self.shape = x[0] self.period = x[1] def _get_param_names(self): return ['shape','period'] def K(self,X,X2,target): if X2 is None: X2 = X pcomp = 2.*np.pi*(X-X2.T)/self.period target += np.cos(np.sin(pcomp))*np.exp(self.shape*(np.cos(pcomp)-1)) def Kdiag(self,X,target): target += 1 def dK_dtheta(self,dL_dK,X,X2,target): if X2 is None: X2 = X pcomp=2.*np.pi*(X-X2.T)/self.period dshape = np.cos(np.sin(pcomp))*np.exp(self.shape*(np.cos(pcomp)-1))*(np.cos(pcomp)-1) dper = sin(sin(pcomp))*cos(pcomp)*pcomp*np.exp(self.shape*(cos(pcomp)-1))+cos(sin(pcomp))*np.exp(self.shape*(np.cos(pcomp)-1))*self.shape*pcomp*np.sin(pcomp); target[0] += np.sum(dshape*dL_dK) target[1] += np.sum(dper*dL_dK) def dKdiag_dtheta(self,dL_dKdiag,X,target): # !!! these are not implemented target[0] += np.sum(dL_dKdiag) # here self.lengthscale and self.power have no influence on Kdiag so target[1:] are unchanged def dK_dX(self,dL_dK,X,X2,target): # !!! these are not implemented """derivative of the covariance matrix with respect to X.""" # if X2 is None: # dist2 = np.square((X-X.T)/self.lengthscale) # dX = -2.*self.variance*self.power * (X-X.T)/self.lengthscale**2 * (1 + dist2/2./self.lengthscale)**(-self.power-1) # else: # dist2 = np.square((X-X2.T)/self.lengthscale) # dX = -self.variance*self.power * (X-X2.T)/self.lengthscale**2 * (1 + dist2/2./self.lengthscale)**(-self.power-1) target += np.sum(dL_dK*dX,1)[:,np.newaxis] def dKdiag_dX(self,dL_dKdiag,X,target): pass