kernels with rosen6d example

Author: bartzbeielstein

notebooks/00_spotPython_tests.ipynb

@@ -7,7 +7,7 @@ build-backend = "setuptools.build_meta"
 
 [project]
 name = "spotpython"
-version = "0.35.3"
+version = "0.36.0"
 authors = [
   { name="T. Bartz-Beielstein", email="tbb@bartzundbartz.de" }
 ]

notebooks/spot_aquisition_random_rosenbrock_2d.ipynb

@@ -49,11 +49,12 @@
     get_ith_hyperparameter_name_from_fun_control,
 )
 import plotly.graph_objects as go
-from typing import Callable
+from typing import Callable, Tuple
 from spotpython.utils.numpy2json import NumpyEncoder
 from spotpython.utils.file import load_result
 from spotpython.surrogate.plot import plotkd
 
+
 # Setting up the backend to use QtAgg
 # matplotlib.use("TkAgg")
 # matplotlib.use("Agg")
@@ -452,14 +453,14 @@ def _set_additional_attributes(self) -> None:
             matplotlib.use("TkAgg")
         self.verbosity = self.fun_control["verbosity"]
         self.acquisition_failure_strategy = self.fun_control["acquisition_failure_strategy"]
-        self.kernel = self.fun_control["kernel"]
-        self.kernel_params = self.fun_control["kernel_params"]
 
         # Surrogate control attributes:
         self.max_surrogate_points = self.surrogate_control["max_surrogate_points"]
         self.use_nystrom = self.surrogate_control["use_nystrom"]
         self.nystrom_m = self.surrogate_control["nystrom_m"]
         self.nystrom_seed = self.surrogate_control["nystrom_seed"]
+        self.kernel = self.surrogate_control["kernel"]
+        self.kernel_params = self.surrogate_control["kernel_params"]
 
         # Internal attributes:
         self.X = None

pyproject.toml

@@ -0,0 +1,98 @@
+import numpy as np
+from sklearn.gaussian_process.kernels import Kernel, NormalizedKernelMixin, Hyperparameter
+
+
+def _correlation(kernel, D, kernel_params=None):
+    """
+    Dispatches to the selected kernel function.
+    Args:
+        kernel: Kernel type (str, callable, or sklearn-compatible kernel object)
+        D (np.ndarray): Distance matrix.
+        kernel_params (dict): Parameters for the kernel (optional).
+    Returns:
+        np.ndarray: Correlation matrix.
+    """
+    kernel_params = kernel_params or {}
+    # If kernel is a sklearn-compatible kernel object
+    if isinstance(kernel, Kernel):
+        # These expect X, Y, not D, so we call on D as before for compatibility
+        return kernel(D)
+    elif callable(kernel):
+        return kernel(D, **kernel_params)
+    elif kernel == "gauss":
+        return np.exp(-D)
+    elif kernel == "matern":
+        nu = kernel_params.get("nu", 2.5)
+        if nu == 0.5:
+            return np.exp(-np.sqrt(D))
+        elif nu == 1.5:
+            sqrt3D = np.sqrt(3.0 * D)
+            return (1.0 + sqrt3D) * np.exp(-sqrt3D)
+        elif nu == 2.5:
+            sqrt5D = np.sqrt(5.0 * D)
+            return (1.0 + sqrt5D + (5.0 / 3.0) * D) * np.exp(-sqrt5D)
+        else:
+            return np.exp(-D)
+    elif kernel == "exp":
+        return np.exp(-np.sqrt(D))
+    elif kernel == "cubic":
+        return 1.0 - D**3
+    elif kernel == "linear":
+        return 1.0 - D
+    elif kernel == "rq":
+        alpha = kernel_params.get("alpha", 1.0)
+        return (1.0 + D / (2.0 * alpha)) ** (-alpha)
+    elif kernel == "poly":
+        degree = kernel_params.get("degree", 2)
+        return (1.0 + D) ** degree
+    else:
+        raise ValueError(f"Unknown kernel: {kernel}")
+
+
+# Example: Custom sklearn-compatible RBF kernel using _correlation
+class CustomRBF(NormalizedKernelMixin, Kernel):
+    def __init__(self, length_scale=1.0, length_scale_bounds=(1e-5, 1e5)):
+        self.length_scale = length_scale
+        self.length_scale_bounds = length_scale_bounds
+
+    @property
+    def hyperparameter_length_scale(self):
+        return Hyperparameter("length_scale", "numeric", self.length_scale_bounds)
+
+    def __call__(self, D, eval_gradient=False):
+        # D is assumed to be the squared distance matrix
+        K = _correlation("gauss", D)
+        if eval_gradient:
+            # Gradient not implemented
+            return K, np.empty((K.shape[0], K.shape[1], 0))
+        return K
+
+    def diag(self, X):
+        return np.ones(X.shape[0])
+
+    def is_stationary(self):
+        return True
+
+
+class CustomMatern(NormalizedKernelMixin, Kernel):
+    def __init__(self, length_scale=1.0, nu=1.5, length_scale_bounds=(1e-5, 1e5)):
+        self.length_scale = length_scale
+        self.nu = nu
+        self.length_scale_bounds = length_scale_bounds
+
+    @property
+    def hyperparameter_length_scale(self):
+        return Hyperparameter("length_scale", "numeric", self.length_scale_bounds)
+
+    def __call__(self, D, eval_gradient=False):
+        K = _correlation("matern", D, {"nu": self.nu})
+        if eval_gradient:
+            # Gradient not implemented
+            return K, np.empty((K.shape[0], K.shape[1], 0))
+        return K
+
+    def diag(self, X):
+        return np.ones(X.shape[0])
+
+    def is_stationary(self):
+        return True

src/spotpython/spot/spot.py

@@ -8,53 +8,7 @@
 from numpy import linspace, append
 from scipy.spatial.distance import cdist, pdist, squareform
 from spotpython.surrogate.plot import plotkd
-
-# --- Kernel functions ---
-
-
-def gauss_kernel(D):
-    """Gaussian (RBF) kernel: exp(-D)"""
-    return np.exp(-D)
-
-
-def matern_kernel(D, nu=2.5):
-    """Matern kernel (default nu=2.5, smooth)."""
-    if nu == 0.5:
-        return np.exp(-np.sqrt(D))
-    elif nu == 1.5:
-        sqrt3D = np.sqrt(3.0 * D)
-        return (1.0 + sqrt3D) * np.exp(-sqrt3D)
-    elif nu == 2.5:
-        sqrt5D = np.sqrt(5.0 * D)
-        return (1.0 + sqrt5D + (5.0 / 3.0) * D) * np.exp(-sqrt5D)
-    else:
-        # Fallback to Gaussian for unsupported nu
-        return np.exp(-D)
-
-
-def exponential_kernel(D):
-    """Exponential kernel: exp(-sqrt(D))"""
-    return np.exp(-np.sqrt(D))
-
-
-def cubic_kernel(D):
-    """Cubic kernel: 1 - D^3"""
-    return 1.0 - D**3
-
-
-def linear_kernel(D):
-    """Linear kernel: 1 - D"""
-    return 1.0 - D
-
-
-def rational_quadratic_kernel(D, alpha=1.0):
-    """Rational Quadratic kernel: (1 + D/(2*alpha))^(-alpha)"""
-    return (1.0 + D / (2.0 * alpha)) ** (-alpha)
-
-
-def poly_kernel(D, degree=2):
-    """Polynomial kernel: (1 + D)^degree"""
-    return (1.0 + D) ** degree
+from spotpython.surrogate.kernels import _correlation  # Only import _correlation
 
 
 # --- The New Kriging Class with Nyström Approximation as introduced in v0.34.0 ---
@@ -258,38 +212,17 @@ def __init__(
         self.Rinv_r_ = None  # R^{-1} r (n,)
 
     # --- Kernel dispatch ---
-
     def _correlation(self, D):
         """
-        Dispatches to the selected kernel function.
+        Dispatches to the selected kernel function using the shared _correlation() from kernels.py.
 
         Args:
             D (np.ndarray): Distance matrix.
 
         Returns:
             np.ndarray: Correlation matrix.
         """
-        if callable(self.kernel):
-            return self.kernel(D, **self.kernel_params)
-        elif self.kernel == "gauss":
-            return gauss_kernel(D)
-        elif self.kernel == "matern":
-            nu = self.kernel_params.get("nu", 2.5)
-            return matern_kernel(D, nu=nu)
-        elif self.kernel == "exp":
-            return exponential_kernel(D)
-        elif self.kernel == "cubic":
-            return cubic_kernel(D)
-        elif self.kernel == "linear":
-            return linear_kernel(D)
-        elif self.kernel == "rq":
-            alpha = self.kernel_params.get("alpha", 1.0)
-            return rational_quadratic_kernel(D, alpha=alpha)
-        elif self.kernel == "poly":
-            degree = self.kernel_params.get("degree", 2)
-            return poly_kernel(D, degree=degree)
-        else:
-            raise ValueError(f"Unknown kernel: {self.kernel}")
+        return _correlation(self.kernel, D, self.kernel_params)
 
     # -------- Basis correlation construction (identical to kriging.py) --------
 

src/spotpython/surrogate/kernels.py

@@ -57,8 +57,6 @@ def fun_control_init(
     horizon=None,
     hyperdict=None,
     infill_criterion="y",
-    kernel="gauss",
-    kernel_params={},
     log_every_n_steps=50,
     log_level=50,
     lower=None,
@@ -211,11 +209,6 @@ def fun_control_init(
             For example: `spotriver.hyperdict.river_hyper_dict import RiverHyperDict`
         infill_criterion (str):
             Can be `"y"`, `"s"`, `"ei"` (negative expected improvement), or `"all"`. Default is "y".
-        kernel (str):
-            The kernel to be used by the Kriging surrogate model.
-            Can be either "gauss", "matern32", "matern52", or "exponential". Default is "gauss".
-        kernel_params (dict):
-            The parameters for the kernel function. Default is an empty dictionary.
         log_every_n_steps (int):
             Lightning: How often to log within steps. Default: 50.
         log_level (int):
@@ -404,8 +397,6 @@ def fun_control_init(
                 'horizon': 7,
                 'infill_criterion': 'y',
                 'k_folds': None,
-                'kernel': 'gauss',
-                'kernel_params': {},
                 'loss_function': None,
                 'lower': None,
                 'max_surrogate_points': 100,
@@ -498,8 +489,6 @@ def fun_control_init(
         "hyperdict": hyperdict,
         "infill_criterion": infill_criterion,
         "k_folds": 3,
-        "kernel": kernel,
-        "kernel_params": kernel_params,
         "log_every_n_steps": log_every_n_steps,
         "log_graph": False,
         "log_level": log_level,
@@ -765,6 +754,8 @@ def design_control_init(init_size=10, repeats=1) -> dict:
 
 
 def surrogate_control_init(
+    kernel="gauss",
+    kernel_params={},
     log_level: int = 50,
     method="regression",
     model_optimizer=differential_evolution,
@@ -789,6 +780,11 @@ def surrogate_control_init(
     """Initialize surrogate_control dictionary.
 
     Args:
+        kernel (str):
+            The kernel to be used by the Kriging surrogate model.
+            Can be either "gauss", "matern32", "matern52", or "exponential". Default is "gauss".
+        kernel_params (dict):
+            The parameters for the kernel function. Default is an empty dictionary.
         model_optimizer (object):
             The optimizer object used for the search on surrogate.
             Default is differential_evolution.
@@ -871,6 +867,8 @@ def surrogate_control_init(
         "min_theta": min_theta,
         "max_theta": max_theta,
         "isotropic": isotropic,
+        "kernel": kernel,
+        "kernel_params": kernel_params,
         "p_val": p_val,
         "n_p": n_p,
         "optim_p": optim_p,

src/spotpython/surrogate/kriging.py

@@ -36,13 +36,12 @@
     tolerance_x=tolerance_x,
     TENSORBOARD_CLEAN=True,
     tensorboard_log=True,
-    kernel=kernel,
-    kernel_params=kernel_params,
     selection_method=selection_method,
     min_success_rate=min_success_rate
 )
 design_control = design_control_init(init_size=init_size)
-surrogate_control_exact = surrogate_control_init(use_nystrom=use_nystrom, method=method, max_surrogate_points=max_surrogate_points, min_Lambda=min_Lambda, max_Lambda=max_Lambda, min_theta=min_theta, max_theta=max_theta, isotropic=isotropic)
+surrogate_control_exact = surrogate_control_init(use_nystrom=use_nystrom, method=method, max_surrogate_points=max_surrogate_points, min_Lambda=min_Lambda, max_Lambda=max_Lambda, min_theta=min_theta, max_theta=max_theta, isotropic=isotropic,    kernel=kernel,
+    kernel_params=kernel_params)
 
 def test_run_1():
     spot_exact_y = Spot(