From a7d2b7da828d7b32b551524c78a3da1f2c92cb20 Mon Sep 17 00:00:00 2001 From: Duncs <49517207+Duncan10@users.noreply.github.com> Date: Sat, 27 Aug 2022 11:45:40 +0200 Subject: [PATCH 1/4] add Wiener Velocity Kernel --- GPy/kern/__init__.py | 1 + GPy/kern/src/wiener_velocity.py | 61 +++++++++++++++++++++++++++++++++ 2 files changed, 62 insertions(+) create mode 100644 GPy/kern/src/wiener_velocity.py diff --git a/GPy/kern/__init__.py b/GPy/kern/__init__.py index 804672901..f6989a474 100644 --- a/GPy/kern/__init__.py +++ b/GPy/kern/__init__.py @@ -38,6 +38,7 @@ from .src.linear import Linear, LinearFull from .src.static import Bias, White, Fixed, WhiteHeteroscedastic, Precomputed from .src.brownian import Brownian +from .src.wiener_velocity import WienerVelocity from .src.stationary import Exponential, OU, Matern32, Matern52, ExpQuad, RatQuad, Cosine, Sinc, ExpQuadCosine from .src.mlp import MLP from .src.periodic import PeriodicExponential, PeriodicMatern32, PeriodicMatern52 diff --git a/GPy/kern/src/wiener_velocity.py b/GPy/kern/src/wiener_velocity.py new file mode 100644 index 000000000..4dcd0aadd --- /dev/null +++ b/GPy/kern/src/wiener_velocity.py @@ -0,0 +1,61 @@ +# Copyright (c) 2012, GPy authors (see AUTHORS.txt). +# Licensed under the BSD 3-clause license (see LICENSE.txt) + +from .kern import Kern +from ...core.parameterization import Param +from paramz.transformations import Logexp +import numpy as np + +class WienerVelocity(Kern): + """ + Wiener Velocity in 1D only. + + Negative times are treated as a separate (backwards!) motion. + + :param input_dim: the number of input dimensions + :type input_dim: int + :param variance: + :type variance: float + """ + def __init__(self, input_dim=1, variance=1., active_dims=None, name='WienerVelocity'): + assert input_dim==1, "Wiener velocity in 1D only" + super(WienerVelocity, self).__init__(input_dim, active_dims, name) + + self.variance = Param('variance', variance, Logexp()) + self.link_parameters(self.variance) + + def to_dict(self): + """ + Convert the object into a json serializable dictionary. + Note: It uses the private method _save_to_input_dict of the parent. + :return dict: json serializable dictionary containing the needed information to instantiate the object + """ + + input_dict = super(RBF, self)._save_to_input_dict() + input_dict["class"] = "GPy.kern.WienerVelocity" + return input_dict + + def K(self,X,X2=None): + if X2 is None: + X2 = X + return self.variance*np.where(np.sign(X)==np.sign(X2.T),(np.fmin(np.abs(X),np.abs(X2.T))**3) / 3 + np.abs(X - X2.T) * (np.fmin(np.abs(X),np.abs(X2.T))**2) / 2, 0.) + + def Kdiag(self,X): + return self.variance*np.divide(np.abs(X.flatten())**3,3) + + def update_gradients_full(self, dL_dK, X, X2=None): + if X2 is None: + X2 = X + self.variance.gradient = np.sum(dL_dK * np.where(np.sign(X)==np.sign(X2.T),np.fmin(np.abs(X),np.abs(X2.T)), 0.)) + + #def update_gradients_diag(self, dL_dKdiag, X): + #self.variance.gradient = np.dot(np.abs(X.flatten()), dL_dKdiag) + + #def gradients_X(self, dL_dK, X, X2=None): + #if X2 is None: + #return np.sum(self.variance*dL_dK*np.abs(X),1)[:,None] + #else: + #return np.sum(np.where(np.logical_and(np.abs(X) Date: Sun, 29 Jan 2023 13:13:10 +0100 Subject: [PATCH 2/4] add check_kernel_gradient_functions test case for new kernel --- GPy/kern/src/wiener_velocity.py | 59 +++++++++++++++++++-------------- GPy/testing/kernel_tests.py | 5 +++ 2 files changed, 39 insertions(+), 25 deletions(-) diff --git a/GPy/kern/src/wiener_velocity.py b/GPy/kern/src/wiener_velocity.py index 4dcd0aadd..b772a7aab 100644 --- a/GPy/kern/src/wiener_velocity.py +++ b/GPy/kern/src/wiener_velocity.py @@ -6,56 +6,65 @@ from paramz.transformations import Logexp import numpy as np + class WienerVelocity(Kern): """ - Wiener Velocity in 1D only. - + Wiener Velocity is a 1D kernel only. Negative times are treated as a separate (backwards!) motion. + The Wiener velocity kernel corresponds to a once integrated Brownian motion kernel, + as described in Solin: "Stochastic Differential Equation Methods for Spatio-Temporal Gaussian Process Regression", 2016. + URL: http://urn.fi/URN:ISBN:978-952-60-6711-7 :param input_dim: the number of input dimensions :type input_dim: int :param variance: :type variance: float """ + def __init__(self, input_dim=1, variance=1., active_dims=None, name='WienerVelocity'): - assert input_dim==1, "Wiener velocity in 1D only" + assert input_dim == 1, "Wiener velocity in 1D only" super(WienerVelocity, self).__init__(input_dim, active_dims, name) self.variance = Param('variance', variance, Logexp()) self.link_parameters(self.variance) - - def to_dict(self): - """ - Convert the object into a json serializable dictionary. - Note: It uses the private method _save_to_input_dict of the parent. - :return dict: json serializable dictionary containing the needed information to instantiate the object - """ - input_dict = super(RBF, self)._save_to_input_dict() + def to_dict(self): + input_dict = super(WienerVelocity, self)._save_to_input_dict() input_dict["class"] = "GPy.kern.WienerVelocity" + input_dict["variance"] = self.variance.values.tolist() return input_dict - def K(self,X,X2=None): + @staticmethod + def _build_from_input_dict(kernel_class, input_dict): + useGPU = input_dict.pop('useGPU', None) + return WienerVelocity(**input_dict) + + def K(self, X, X2=None): if X2 is None: X2 = X - return self.variance*np.where(np.sign(X)==np.sign(X2.T),(np.fmin(np.abs(X),np.abs(X2.T))**3) / 3 + np.abs(X - X2.T) * (np.fmin(np.abs(X),np.abs(X2.T))**2) / 2, 0.) + return (self.variance*np.where(np.sign(X) == np.sign(X2.T), (np.fmin(np.abs(X), np.abs(X2.T))**3) / + 3 + np.abs(X - X2.T) * (np.fmin(np.abs(X), np.abs(X2.T))**2) / 2, 0.)) - def Kdiag(self,X): - return self.variance*np.divide(np.abs(X.flatten())**3,3) + def Kdiag(self, X): + return self.variance*np.divide(np.abs(X.flatten())**3, 3) def update_gradients_full(self, dL_dK, X, X2=None): if X2 is None: X2 = X - self.variance.gradient = np.sum(dL_dK * np.where(np.sign(X)==np.sign(X2.T),np.fmin(np.abs(X),np.abs(X2.T)), 0.)) - - #def update_gradients_diag(self, dL_dKdiag, X): - #self.variance.gradient = np.dot(np.abs(X.flatten()), dL_dKdiag) - - #def gradients_X(self, dL_dK, X, X2=None): - #if X2 is None: - #return np.sum(self.variance*dL_dK*np.abs(X),1)[:,None] - #else: - #return np.sum(np.where(np.logical_and(np.abs(X) Date: Sun, 5 Feb 2023 21:59:51 +0100 Subject: [PATCH 3/4] removed update_gradients_diag --- GPy/kern/src/wiener_velocity.py | 15 --------------- 1 file changed, 15 deletions(-) diff --git a/GPy/kern/src/wiener_velocity.py b/GPy/kern/src/wiener_velocity.py index b772a7aab..63a6fad87 100644 --- a/GPy/kern/src/wiener_velocity.py +++ b/GPy/kern/src/wiener_velocity.py @@ -53,18 +53,3 @@ def update_gradients_full(self, dL_dK, X, X2=None): X2 = X self.variance.gradient = (np.sum(dL_dK * np.where(np.sign(X) == np.sign(X2.T), (np.fmin(np.abs(X), np.abs(X2.T))**3) / 3 + np.abs(X - X2.T) * (np.fmin(np.abs(X), np.abs(X2.T))**2) / 2, 0.))) - - def update_gradients_diag(self, dL_dKdiag, X, X2=None): - if X2 is None: - X2 = X - self.variance.gradient = (np.sum(dL_dKdiag * np.where(np.sign(X) == np.sign(X2.T), (np.fmin(np.abs(X), np.abs(X2.T))**3) / - 3 + np.abs(X - X2.T) * (np.fmin(np.abs(X), np.abs(X2.T))**2) / 2, 0.))) - - # def update_gradients_diag(self, dL_dKdiag, X): - # self.variance.gradient = np.dot(np.abs(X.flatten()), dL_dKdiag) - - # def gradients_X(self, dL_dK, X, X2=None): - # if X2 is None: - # return np.sum(self.variance*dL_dK*np.abs(X),1)[:,None] - # else: - # return np.sum(np.where(np.logical_and(np.abs(X) Date: Sun, 26 Feb 2023 00:21:47 +0100 Subject: [PATCH 4/4] insert useGPU keyword in init for model loading compatibility --- GPy/kern/src/wiener_velocity.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/GPy/kern/src/wiener_velocity.py b/GPy/kern/src/wiener_velocity.py index 63a6fad87..dcf623579 100644 --- a/GPy/kern/src/wiener_velocity.py +++ b/GPy/kern/src/wiener_velocity.py @@ -21,7 +21,7 @@ class WienerVelocity(Kern): :type variance: float """ - def __init__(self, input_dim=1, variance=1., active_dims=None, name='WienerVelocity'): + def __init__(self, input_dim=1, variance=1., active_dims=None, name='WienerVelocity', useGPU=False): assert input_dim == 1, "Wiener velocity in 1D only" super(WienerVelocity, self).__init__(input_dim, active_dims, name) @@ -36,7 +36,6 @@ def to_dict(self): @staticmethod def _build_from_input_dict(kernel_class, input_dict): - useGPU = input_dict.pop('useGPU', None) return WienerVelocity(**input_dict) def K(self, X, X2=None):