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