Summarize numerical conditioning of a fitted GP — gp_numerical

Reports numerical quantities separately from statistical fit diagnostics. In particular, observation-noise variance and numerical jitter remain distinct.

Usage

gp_numerical_diagnostics(
  model,
  condition_tolerance = sqrt(.Machine$double.eps)
)

Arguments

model: A fitted gaussianprocesses_model.
condition_tolerance: Positive reciprocal-condition threshold used only to flag a potentially ill-conditioned stabilized covariance matrix.

Value

A named list of numerical diagnostics, including:

observation_noise_variance: The statistical noise variance.
numerical_jitter: The diagonal jitter \(\delta\) that the factorization of the training covariance \(C\) needed.
covariance_scale: The mean of the diagonal of \(C\).
relative_jitter: numerical_jitter / covariance_scale.
cholesky_attempts: Factorization attempts, 1 if no jitter was needed.
reciprocal_condition_number, condition_number: The 1-norm condition estimate of \(C + \delta I\) and its reciprocal.
potentially_ill_conditioned: Whether the reciprocal condition number is below condition_tolerance.
observation_types: A data frame with one row per observation type ("value", or "derivative x1", ... for derivative observations): the number of observations, the mean diagonal of their block of \(C\), and the reciprocal condition number of that block with the jitter. Derivative blocks scale with \(1 / \ell^2\), so their conditioning can differ from that of the values.

Stability

Experimental: this interface may change in a minor release, and every change is listed in NEWS. See gaussianprocesses-package for the policy.

Examples

# Duplicated inputs without observation noise need numerical jitter, which is
# reported separately from the (zero) noise variance.
model <- fit_gp(c(0, 0, 1), c(1, 1, 0), kernel = rbf_kernel(), noise_variance = 0)

diagnostics <- gp_numerical_diagnostics(model)
diagnostics[c(
  "observation_noise_variance",
  "numerical_jitter",
  "relative_jitter",
  "cholesky_attempts"
)]
#> $observation_noise_variance
#> [1] 0
#> 
#> $numerical_jitter
#> [1] 1e-10
#> 
#> $relative_jitter
#> [1] 1e-10
#> 
#> $cholesky_attempts
#> [1] 2
#>