Draw functions from a Gaussian-process prior

Draws latent function values \(f(x) \sim \mathcal{N}(m(x), K(x, x))\) at the supplied inputs. Use simulate_gp_data() for a single prior draw with added observation noise.

Usage

sample_gp_prior(
  x,
  kernel,
  mean = zero_mean(),
  n_draws = 1L,
  seed = NULL,
  initial_jitter = 0,
  fallback_jitter = 1e-10,
  jitter_multiplier = 10,
  max_attempts = 8L,
  symmetry_tolerance = sqrt(.Machine$double.eps)
)

Arguments

x

Numeric vector or matrix of inputs. Matrix rows are input points.

kernel

A Gaussian-process kernel specification, including ARD and composite kernels.

mean

Gaussian-process mean specification. Its coefficients must be fixed or have a Gaussian coefficient_prior(), whose uncertainty the draws include.

n_draws

Number of independent function draws.

seed

Optional non-negative integer seed. The caller's random state is restored afterwards.

initial_jitter

Non-negative jitter tried first, relative to the covariance scale (the mean of its diagonal).

fallback_jitter

Positive relative jitter tried after the first failed factorization.

jitter_multiplier

Multiplicative jitter escalation factor.

max_attempts

Maximum Cholesky attempts.

symmetry_tolerance

Relative tolerance for checking that the covariance matrix is symmetric.

Value

An object of class gaussianprocesses_samples. Its latent element is a matrix with one row per input point and one column per draw. It also contains the prior mean, the latent_covariance that was sampled, the numerical jitter, the relative_jitter and covariance_scale it was computed from (jitter is relative jitter times the mean diagonal of the covariance), and the cholesky_attempts.

Details

Each draw is \(m(x) + R^\top z\) with \(K(x, x) = R^\top R\) and \(z \sim \mathcal{N}(0, I)\). Prior covariances on closely spaced inputs are often numerically singular; the Cholesky factor then uses the package's jitter policy, and the jitter is reported in the result. A single draw with the same seed equals the latent values of simulate_gp_data().

Stability

Stable: from version 1.0.0 this interface changes incompatibly only in a major release, after a deprecation period. Results and options that concern an experimental model class, kernel, or argument follow that interface's tier. See gaussianprocesses-package for the policy.

Examples

x <- seq(0, 5, length.out = 100)
draws <- sample_gp_prior(x, rbf_kernel(length_scale = 0.8), n_draws = 3, seed = 1)
draws
#> Gaussian-process samples
#>   distribution: prior
#>   draws: 3 at 100 input point(s)
#>   contents: latent function draws
#>   numerical jitter: 1e-10 (1e-10 times the covariance scale 1)

matplot(x, draws$latent, type = "l", lty = 1, ylab = "f(x)")