gpmp.kernel covariance functions¶

Covariance functions accept gpmp.num arrays and covariance-parameter vectors. They can be passed directly to gpmp.core.Model. gpmp.parameter objects are optional display objects and are not part of the covariance-function call signature.

The covariance functions documented here define stationary anisotropic covariances of the form

\[k_\theta(x, y) = \sigma^2 r_\eta(h_\rho(x, y)),\]

where \(\sigma^2\) is the process variance, \(\rho_j\) are coordinate lengthscales, \(h_\rho\) is the scaled Euclidean distance, and \(r_\eta\) is a correlation kernel. This is the covariance function used in kriging / Gaussian-process interpolation; see Chilès and Delfiner [1], Stein [13] for covariance choice, regularity, and lengthscales in spatial interpolation.

Correlation kernels and covariance functions¶

Correlation kernels take a scaled distance h and return a correlation value. They do not include the variance parameter. The public correlation kernels are exponential_kernel, squared_exponential_kernel, matern12_kernel, matern32_kernel, matern52_kernel, maternp_kernel, and matern_kernel.

Covariance functions take point arrays and a covariance-parameter vector. They include the process variance and the coordinate scaling. The covariance functions documented here are squared_exponential_covariance, maternp_covariance, and matern_covariance.

Scaled-distance convention¶

For anisotropic covariance functions, covparam stores loginvrho_j = -log(rho_j). For points \(x, y \in \mathbb{R}^d\), GPmp uses

\[h_\rho(x, y) = \left( \sum_{j=0}^{d-1} \left(\exp(\mathrm{loginvrho}_j)(x_j - y_j)\right)^2 \right)^{1/2} = \left( \sum_{j=0}^{d-1} \left(\frac{x_j - y_j}{\rho_j}\right)^2 \right)^{1/2}.\]

The scaled squared distance is

\[q_\rho(x, y) = h_\rho(x, y)^2.\]

Thus a larger stored value loginvrho_j means a shorter lengthscale in coordinate j. The same convention is used by the fixed-regularity Matérn, continuous-regularity Matérn, and squared-exponential covariance functions.

Self-covariance matrix functions add a small numerical nugget \(10 \sigma^2 \epsilon_{\mathrm{mach}}\) on the diagonal. This is a linear-algebra safeguard, not an observation-noise model.

Squared-exponential scaling¶

The squared-exponential kernel follows the same scaling convention as the Matérn kernels. If h is the scaled distance, GPmp uses

\[r_{\mathrm{SE}}(h) = \exp(-h^2).\]

Equivalently,

\[k_\theta(x, y) = \sigma^2 \exp\{-q_\rho(x, y)\}.\]

GPmp does not use the alternative exp(-0.5 * h**2) convention. This choice makes the squared-exponential kernel consistent with the continuous-regularity Matérn implementation, whose large-nu limit is exp(-h**2) under the GPmp scaling.

Fixed-regularity Matérn covariance¶

maternp_covariance uses the sigma2_rho convention

[log(sigma2), -log(rho_0), ..., -log(rho_{d-1})].

The regularity is fixed by the integer p with \(\nu = p + 1/2\). GPmp evaluates the half-integer Matérn correlation as

\[r_p(h) = \exp(-2\sqrt{\nu}h) \frac{\Gamma(p+1)}{\Gamma(2p+1)} \sum_{i=0}^{p} \frac{(p+i)!}{i!(p-i)!} (4\sqrt{\nu}h)^{p-i}, \qquad \nu = p + \frac{1}{2}.\]

The closed forms matern12_kernel, matern32_kernel, and matern52_kernel are provided for the common cases p = 0, p = 1, and p = 2:

\[\begin{split}r_{1/2}(h) &= \exp(-\sqrt{2}h), \\ r_{3/2}(h) &= (1 + \sqrt{6}h)\exp(-\sqrt{6}h), \\ r_{5/2}(h) &= \left(1 + \sqrt{10}h + \frac{10h^2}{3}\right)\exp(-\sqrt{10}h).\end{split}\]

The covariance is then

\[k_\theta(x, y) = \sigma^2 r_p(h_\rho(x, y)).\]

Continuous-regularity Matérn covariance¶

matern_covariance uses the sigma2_nu_rho convention

[log(sigma2), log(nu), -log(rho_0), ..., -log(rho_{d-1})].

Here nu is optimized or supplied as part of the covariance-parameter vector. This convention requires an initializer and selection method that know about the extra log(nu) component.

For \(\nu > 0\), the continuous-regularity Matérn correlation is

\[r_\nu(h) = \frac{2^{1-\nu}}{\Gamma(\nu)} (2\sqrt{\nu}h)^\nu K_\nu(2\sqrt{\nu}h), \qquad r_\nu(0)=1,\]

where \(K_\nu\) is the modified Bessel function of the second kind. The covariance is

\[k_\theta(x, y) = \sigma^2 r_\nu(h_\rho(x, y)).\]

Under this scaling, \(r_\nu(h)\) tends to \(\exp(-h^2)\) as \(\nu \to \infty\), which explains the squared-exponential convention used in GPmp. The regularity parameter controls mean-square smoothness of the GP and is a central parameter in Matérn modeling [13].

The backend-dependent evaluation policy, the Torch derivatives, and the nodes, t_max, and large_nu_method arguments are documented in Numerical evaluation of the Matérn covariance. Normal model construction and parameter selection use the default large_nu_method="auto" policy.

Squared-exponential covariance summary¶

squared_exponential_covariance uses the sigma2_rho convention

[log(sigma2), -log(rho_0), ..., -log(rho_{d-1})].

Under the GPmp distance scaling, the correlation is exp(-q), where q is the scaled squared distance.

Function reference¶

exponential_kernel¶

gpmp.kernel.exponential_kernel(h)[source]

Exponential kernel.

\[k(h) = \exp(-h)\]

Parameters:: h (array_like) – Scaled distances.
Returns:: Correlation values with the same backend as h.
Return type:: array_like

squared_exponential_kernel¶

gpmp.kernel.squared_exponential_kernel(h)[source]¶

Squared-exponential kernel.

\[k(h) = \exp(-h^2)\]

The argument h is the scaled distance used by the GPmp Matérn kernels. With this convention the squared-exponential kernel is the limiting case of the continuous-regularity Matérn kernel. GPmp therefore uses exp(-h**2), not the alternative exp(-0.5 * h**2) convention.

Parameters:: h (array_like) – Scaled distances between points.
Returns:: Correlation values with the same backend as h.
Return type:: array_like

squared_exponential_covariance¶

gpmp.kernel.squared_exponential_covariance(x, y, param, pairwise=False)[source]¶

Squared-exponential covariance.

param follows the convention [log(sigma2), -log(rho_0), ..., -log(rho_{d-1})]. The correlation term is exp(-q), where q is the scaled squared distance associated with this covariance-parameter convention.

Parameters:

x (array_like, shape (nx, d)) – First set of points.
y (array_like or None) – Second set of points. If y is None or y is x, the self-covariance of x is returned.
param (array_like, shape (d + 1,)) – Covariance parameters [log(sigma2), -log(rho_0), ..., -log(rho_{d-1})].
pairwise (bool, optional) – If True, return elementwise covariances. If False, return a full covariance matrix.

Returns:

Covariance matrix or vector, depending on pairwise.

Return type:

array_like

matern12_kernel¶

gpmp.kernel.matern12_kernel(h)[source]¶

Matérn 1/2 kernel.

\[K(h) = \exp(-2\sqrt{1/2}\,h)\]

Parameters:: h (array_like) – Scaled distances.
Returns:: Correlation values with the same backend as h.
Return type:: array_like

matern32_kernel¶

gpmp.kernel.matern32_kernel(h)[source]¶

Matérn 3/2 kernel.

\[K(h) = (1 + 2\sqrt{3/2}\,h) \exp(-2\sqrt{3/2}\,h)\]

Parameters:: h (array_like) – Scaled distances.
Returns:: Correlation values with the same backend as h.
Return type:: array_like

matern52_kernel¶

gpmp.kernel.matern52_kernel(h)[source]¶

Matérn 5/2 kernel.

\[K(h) = \left(1 + t + \frac{t^2}{3}\right)\exp(-t), \quad t = 2\sqrt{5/2}\,h\]

Parameters:: h (array_like) – Scaled distances.
Returns:: Correlation values with the same backend as h.
Return type:: array_like

maternp_kernel¶

gpmp.kernel.maternp_kernel(p: int, h)[source]¶

Matérn kernel with half-integer regularity \(\nu = p + 1/2\).

Using the half-integer simplification from Watson (1922) and Abramowitz and Stegun:

\[K(h) = \exp(-2\sqrt{\nu}\,h)\, \frac{\Gamma(p+1)}{\Gamma(2p+1)} \sum_{i=0}^{p} \frac{(p+i)!}{i!(p-i)!}\,(4\sqrt{\nu}h)^{\,p-i}\]

Parameters:

p (int) – Nonnegative integer with \(\nu = p+1/2\).
h (array_like) – Scaled distances.

Returns:

Correlation values with the same backend as h.

Return type:

array_like

maternp_covariance¶

gpmp.kernel.maternp_covariance(x, y, p, param, pairwise=False)[source]¶

Matérn covariance (\(\nu = p+1/2\)). Wrapper.

Parameters:

x (array_like, shape (nx, d))
y (array_like or None)
p (int)
param (array_like, shape (1 + d,))
pairwise (bool)

Return type:

array_like

matern_kernel¶

gpmp.kernel.matern_kernel(nu, h, *, nodes: int = 160, t_max: float = 20.0, large_nu_method: str = 'auto')[source]¶

Matérn kernel with continuous regularity nu.

h follows the same scaled-distance convention as maternp_kernel. Internally, the continuous Matérn primitive uses the standard sqrt(2 * nu) * r parameterization, so the squared distance is multiplied by two for compatibility with the existing GPmp half-integer kernels.

Parameters:

nu (scalar or array_like) – Positive Matérn regularity parameter.
h (scalar or array_like) – Scaled distance.
nodes (int, optional) – Base number of quadrature nodes used by the Torch backend. NumPy accepts this argument for API consistency, but it has no effect.
t_max (float, optional) – Base upper integration limit for adaptive Torch quadrature.
large_nu_method (str, optional) – Numerical evaluation policy. In both backends, the default "auto" uses a uniform expansion for scalar nu >= 70. Below this threshold, NumPy uses SciPy, while Torch selects among exact formulas, calibrated recurrences, and quadrature. "quadrature" forces the SciPy path in NumPy and quadrature in Torch. Explicit "uniform_asymptotic" requires scalar nu >= 70.

Returns:

Matérn correlation values.

Return type:

scalar or array_like

matern_covariance¶

gpmp.kernel.matern_covariance(x, y, param, pairwise=False, *, nodes: int = 160, t_max: float = 20.0, large_nu_method: str = 'auto')[source]¶

Matérn covariance with continuous regularity.

param follows the convention [log(sigma2), log(nu), -log(rho_0), ..., -log(rho_{d-1})].

Parameters:

x (array_like, shape (nx, d)) – First set of points.
y (array_like or None) – Second set of points. If y is None or y is x, the self-covariance of x is returned.
param (array_like, shape (d + 2,)) – Covariance parameters [log(sigma2), log(nu), -log(rho_0), ..., -log(rho_{d-1})].
pairwise (bool, optional) – If True, return elementwise covariances. If False, return a full covariance matrix.
nodes (int, optional) – Base number of quadrature nodes used by the Torch backend. NumPy accepts this argument for API consistency, but it has no effect.
t_max (float, optional) – Base upper integration limit for adaptive Torch quadrature.
large_nu_method (str, optional) – Numerical evaluation policy. In both backends, the default "auto" uses a uniform expansion for scalar nu >= 70. Below this threshold, NumPy uses SciPy, while Torch selects among exact formulas, calibrated recurrences, and quadrature. "quadrature" forces the SciPy path in NumPy and quadrature in Torch. Explicit "uniform_asymptotic" requires scalar nu >= 70.

Returns:

Covariance matrix or vector, depending on pairwise.

Return type:

array_like

Notes

"quadrature" selects the backend reference path used to validate the automatic large-order branch near its switching threshold.