Coverage for geometric_kernels / kernels / matern_kernel.py: 80%
95 statements
« prev ^ index » next coverage.py v7.13.0, created at 2025-12-15 13:41 +0000
1"""
2Provides :class:`MaternGeometricKernel`,the geometric Matérn kernel---with
3the heat kernel as a special case---that just works.
5It wraps around different kernels and feature maps, dispatching on the space.
7Unless you know exactly what you are doing, use :class:`MaternGeometricKernel`.
8"""
10from plum import dispatch, overload
12from geometric_kernels.feature_maps import (
13 DeterministicFeatureMapCompact,
14 HodgeDeterministicFeatureMapCompact,
15 RandomPhaseFeatureMapCompact,
16 RandomPhaseFeatureMapNoncompact,
17 RejectionSamplingFeatureMapHyperbolic,
18 RejectionSamplingFeatureMapSPD,
19)
20from geometric_kernels.kernels.base import BaseGeometricKernel
21from geometric_kernels.kernels.feature_map import MaternFeatureMapKernel
22from geometric_kernels.kernels.hodge_compositional import MaternHodgeCompositionalKernel
23from geometric_kernels.kernels.karhunen_loeve import MaternKarhunenLoeveKernel
24from geometric_kernels.kernels.matern_kernel_hamming_graph import (
25 MaternKernelHammingGraph,
26)
27from geometric_kernels.spaces import (
28 CompactMatrixLieGroup,
29 DiscreteSpectrumSpace,
30 Graph,
31 GraphEdges,
32 HammingGraph,
33 HodgeDiscreteSpectrumSpace,
34 Hyperbolic,
35 HypercubeGraph,
36 Hypersphere,
37 Mesh,
38 NoncompactSymmetricSpace,
39 Space,
40 SymmetricPositiveDefiniteMatrices,
41)
44def default_feature_map(
45 *, space: Space = None, num: int = None, kernel: BaseGeometricKernel = None
46):
47 """
48 Constructs the default feature map for the specified space or kernel.
50 :param space:
51 A space to construct the feature map on. If provided, kernel must
52 either be omitted or set to None.
53 :param kernel:
54 A kernel to construct the feature map from. If provided, `space` and
55 `num` must either be omitted or set to None.
56 :param num:
57 Controls the number of features (dimensionality of the feature
58 map). If omitted or set to None, the default value for each
59 respective space is used. Must only be provided when
60 constructing a feature map on a space (not from a kernel).
62 :return:
63 Callable which is the respective feature map.
64 """
65 if kernel is not None:
66 if space is not None or num is not None:
67 raise ValueError(
68 "When kernel is provided, space and num must be omitted or set to None"
69 )
70 return feature_map_from_kernel(kernel) # type: ignore[call-overload]
71 elif space is not None:
72 if num is None:
73 num = default_num(space) # type: ignore[call-overload]
74 return feature_map_from_space(space, num) # type: ignore[call-overload]
75 else:
76 raise ValueError(
77 "Either kernel or space must be provided and be different from None"
78 )
81@overload
82def feature_map_from_kernel(kernel: MaternKarhunenLoeveKernel):
83 if isinstance(kernel.space, (CompactMatrixLieGroup, HammingGraph)):
84 # Because `CompactMatrixLieGroup` does not currently support explicit
85 # eigenfunction computation (they only support addition theorem).
86 return RandomPhaseFeatureMapCompact(
87 kernel.space,
88 kernel.num_levels,
89 MaternGeometricKernel._DEFAULT_NUM_RANDOM_PHASES,
90 )
91 else:
92 return DeterministicFeatureMapCompact(kernel.space, kernel.num_levels)
95@overload
96def feature_map_from_kernel(kernel: MaternHodgeCompositionalKernel):
97 return HodgeDeterministicFeatureMapCompact(kernel.space, kernel.num_levels)
100@overload
101def feature_map_from_kernel(kernel: MaternFeatureMapKernel):
102 return kernel.feature_map
105@dispatch
106def feature_map_from_kernel(kernel: BaseGeometricKernel):
107 """
108 Return the default feature map for the specified kernel `kernel`.
110 :param kernel:
111 A kernel to construct the feature map from.
113 :return:
114 A feature map.
115 :rtype: feature_maps.FeatureMap
117 .. note::
118 This function is organized as an abstract dispatcher plus a set of
119 @overload-decorated implementations, one for each type of kernels.
121 When followed by an "empty" @dispatch-decorated function of the same
122 name, plum-dispatch changes the default behavior of the `@overload`
123 decorator, allowing the implementations inside the preceding
124 @overload-decorated functions. This is opposed to the standard behavior
125 when @overload-decorated functions can only provide type signature,
126 while the general implementation should be contained in the function
127 of the same name without an `@overload` decorator.
129 The trick is taken from https://beartype.github.io/plum/integration.html.
131 .. note::
132 For dispatching to work, the empty @dispatch-decorated function should
133 follow (not precede) the @overload-decorated implementations in the code.
134 """
135 raise NotImplementedError(
136 "feature_map_from_kernel is not implemented for the kernel of type %s."
137 % str(type(kernel))
138 )
141@overload
142def feature_map_from_space(space: DiscreteSpectrumSpace, num: int):
143 if isinstance(space, (CompactMatrixLieGroup, HammingGraph)):
144 return RandomPhaseFeatureMapCompact(
145 space, num, MaternGeometricKernel._DEFAULT_NUM_RANDOM_PHASES
146 )
147 elif isinstance(space, Hypersphere):
148 num_computed_levels = space.num_computed_levels
149 if num_computed_levels > 0:
150 return DeterministicFeatureMapCompact(space, min(num, num_computed_levels))
151 else:
152 return RandomPhaseFeatureMapCompact(
153 space, num, MaternGeometricKernel._DEFAULT_NUM_RANDOM_PHASES
154 )
155 else:
157 return DeterministicFeatureMapCompact(space, num)
160@overload
161def feature_map_from_space(space: NoncompactSymmetricSpace, num: int):
162 if isinstance(space, Hyperbolic):
163 return RejectionSamplingFeatureMapHyperbolic(space, num)
164 elif isinstance(space, SymmetricPositiveDefiniteMatrices):
165 return RejectionSamplingFeatureMapSPD(space, num)
166 else:
167 return RandomPhaseFeatureMapNoncompact(space, num)
170@dispatch
171def feature_map_from_space(space: Space, num: int):
172 """
173 Return the default feature map for the specified space `space` and
174 approximation level `num`.
176 :param space:
177 A space to construct the feature map on.
178 :param num:
179 Approximation level.
181 :return:
182 A feature map.
183 :rtype: feature_maps.FeatureMap
185 .. note::
186 This function is organized as an abstract dispatcher plus a set of
187 @overload-decorated implementations, one for each type of spaces.
189 When followed by an "empty" @dispatch-decorated function of the same
190 name, plum-dispatch changes the default behavior of the `@overload`
191 decorator, allowing the implementations inside the preceding
192 @overload-decorated functions. This is opposed to the standard behavior
193 when @overload-decorated functions can only provide type signature,
194 while the general implementation should be contained in the function
195 of the same name without an `@overload` decorator.
197 The trick is taken from https://beartype.github.io/plum/integration.html.
199 .. note::
200 For dispatching to work, the empty @dispatch-decorated function should
201 follow (not precede) the @overload-decorated implementations in the code.
202 """
203 raise NotImplementedError(
204 "feature_map_from_space is not implemented for the space of type %s."
205 % str(type(space))
206 )
209@overload
210def default_num(space: DiscreteSpectrumSpace) -> int:
211 if isinstance(space, CompactMatrixLieGroup):
212 return MaternGeometricKernel._DEFAULT_NUM_LEVELS_LIE_GROUP
213 elif isinstance(space, (Graph, Mesh)):
214 return min(
215 MaternGeometricKernel._DEFAULT_NUM_EIGENFUNCTIONS, space.num_vertices
216 )
217 elif isinstance(space, GraphEdges):
218 return min(MaternGeometricKernel._DEFAULT_NUM_EIGENFUNCTIONS, space.num_edges)
219 elif isinstance(space, (HypercubeGraph, HammingGraph)):
220 return min(MaternGeometricKernel._DEFAULT_NUM_LEVELS, space.dim + 1)
221 else:
222 return MaternGeometricKernel._DEFAULT_NUM_LEVELS
225@overload
226def default_num(space: NoncompactSymmetricSpace) -> int:
227 return MaternGeometricKernel._DEFAULT_NUM_RANDOM_PHASES
230@dispatch
231def default_num(space: Space) -> int:
232 """
233 Return the default approximation level for the `space`.
235 :param space:
236 A space.
238 :return:
239 The default approximation level.
241 .. note::
242 This function is organized as an abstract dispatcher plus a set of
243 @overload-decorated implementations, one for each type of spaces.
245 When followed by an "empty" @dispatch-decorated function of the same
246 name, plum-dispatch changes the default behavior of the `@overload`
247 decorator, allowing the implementations inside the preceding
248 @overload-decorated functions. This is opposed to the standard behavior
249 when @overload-decorated functions can only provide type signature,
250 while the general implementation should be contained in the function
251 of the same name without an `@overload` decorator.
253 The trick is taken from https://beartype.github.io/plum/integration.html.
255 .. note::
256 For dispatching to work, the empty @dispatch-decorated function should
257 follow (not precede) the @overload-decorated implementations in the code.
258 """
259 raise NotImplementedError(
260 "default_num is not implemented for the space of type %s." % str(type(space))
261 )
264class MaternGeometricKernel:
265 """
266 This class represents a Matérn geometric kernel that "just works". Unless
267 you really know what you are doing, you should always use this kernel class.
269 Upon creation, unpacks into a specific geometric kernel based on the
270 provided space, and, optionally, the associated (approximate) feature map.
271 """
273 _DEFAULT_NUM_EIGENFUNCTIONS = 1000
274 _DEFAULT_NUM_LEVELS = 25
275 _DEFAULT_NUM_LEVELS_LIE_GROUP = 20
276 _DEFAULT_NUM_RANDOM_PHASES = 3000
278 def __new__(
279 cls,
280 space: Space,
281 num: int = None,
282 normalize: bool = True,
283 return_feature_map: bool = False,
284 **kwargs,
285 ):
286 r"""
287 Construct a kernel and (if `return_feature_map` is `True`) a
288 feature map on `space`.
290 .. note::
291 See :doc:`this page </theory/feature_maps>` for a brief
292 introduction into feature maps.
294 :param space:
295 Space to construct the kernel on.
297 :param num:
298 If provided, controls the "order of approximation" of the kernel.
299 For the discrete spectrum spaces, this means the number of "levels"
300 that go into the truncated series that defines the kernel (for
301 example, these are unique eigenvalues for the
302 :class:`~.spaces.Hypersphere` or eigenvalues with repetitions for
303 the :class:`~.spaces.Graph` or for the :class:`~.spaces.Mesh`).
304 For the non-compact symmetric spaces
305 (:class:`~.spaces.NoncompactSymmetricSpace`), this is the number
306 of random phases used to construct the kernel.
308 If num=None, we use a (hopefully) reasonable default, which is
309 space-dependent.
311 :param normalize:
312 Normalize the kernel (and the feature map). If normalize=True,
313 then either $k(x, x) = 1$ for all $x \in X$, where $X$ is the
314 `space`, or $\int_X k(x, x) d x = 1$, depending on the space.
316 Defaults to True.
318 .. note::
319 For many kernel methods, $k(\cdot, \cdot)$ and
320 $a k(\cdot, \cdot)$ are indistinguishable, whatever the
321 positive constant $a$ is. For these, it makes sense to use
322 normalize=False to save up some computational overhead. For
323 others, like for the Gaussian process regression, the
324 normalization of the kernel might be important. In
325 these cases, you will typically want to set normalize=True.
327 :param return_feature_map:
328 If `True`, return a feature map (needed e.g. for efficient sampling
329 from Gaussian processes) along with the kernel.
331 Default is False.
333 :param ``**kwargs``:
334 Any additional keyword arguments to be passed to the kernel
335 (like `key`).
337 .. note::
338 For non-compact symmetric spaces, like :class:`~.spaces.Hyperbolic`
339 or :class:`~.spaces.SymmetricPositiveDefiniteMatrices`, the `key`
340 **must** be provided in ``**kwargs``.
341 """
343 kernel: BaseGeometricKernel
344 if isinstance(space, DiscreteSpectrumSpace):
345 num = num or default_num(space)
346 if isinstance(space, HodgeDiscreteSpectrumSpace):
347 kernel = MaternHodgeCompositionalKernel(space, num, normalize=normalize)
348 elif isinstance(space, (HypercubeGraph, HammingGraph)):
349 kernel = MaternKernelHammingGraph(space, num, normalize=normalize)
350 else:
351 kernel = MaternKarhunenLoeveKernel(space, num, normalize=normalize)
352 if return_feature_map:
353 feature_map = default_feature_map(kernel=kernel)
355 elif isinstance(space, NoncompactSymmetricSpace):
356 num = num or default_num(space)
357 if "key" in kwargs:
358 key = kwargs["key"]
359 else:
360 raise ValueError(
361 (
362 "MaternGeometricKernel for %s requires mandatory"
363 " keyword argument 'key' which is a random number"
364 " generator specific to the backend used"
365 )
366 % str(type(space))
367 )
368 feature_map = default_feature_map(space=space, num=num)
369 kernel = MaternFeatureMapKernel(
370 space, feature_map, key, normalize=normalize
371 )
372 else:
373 raise NotImplementedError
375 if return_feature_map:
376 return kernel, feature_map
377 else:
378 return kernel