Fitted model objects — gaussianprocesses

The structure, versioning, and persistence of the models returned by fit_gp(), optimize_gp(), fit_time_series_gp(), fit_sparse_gp(), fit_heteroscedastic_gp(), fit_latent_gp(), optimize_latent_gp(), and optimize_sparse_gp().

Exact models

Objects of class gaussianprocesses_model, from fit_gp(), optimize_gp(), and fit_time_series_gp(), document these fields:

x, y: Training inputs as a matrix and responses.
derivative: For derivative observations, the type of each observation: 0 for a function value and $d$ for $\partial f / \partial x_d$. NULL when every observation is a function value.
kernel, mean: The kernel and mean specifications. A polynomial mean records the centring and scaling of its inputs.
mean_fit: For estimated or marginalized mean coefficients, a list with treatment, method ("ml", "reml", "gaussian", or "vague"), the coefficients $\bar\beta$, the basis matrix basis ($H$), and coefficient_factor, an upper-triangular $G$ with $G^\top G = H^\top C^{-1} H$ plus $B^{-1}$ under a Gaussian prior. NULL for other means. Use coef() and vcov() to read the coefficients.
noise_variance, training_noise_variance, noise_structure: The noise variance as supplied, one value per training observation, and "homoscedastic" or "known_heteroscedastic".
cholesky: Upper-triangular $R$ with $R^\top R = C + \delta I$, where $C = K + \Sigma_\varepsilon$ and $\delta$ is jitter.
alpha: $\alpha = (C + \delta I)^{-1}(y - m(X))$, with $m(X) = H \bar\beta$ for estimated or marginalized coefficients.
jitter: Diagonal jitter $\delta$ that the factorization needed; usually 0. It is a multiple of the mean diagonal of $C$ (see fit_gp()).
n_observations, n_features: Training size and input dimension.
optimization: Optimizer diagnostics, for optimize_gp() models.
time: Time metadata, for fit_time_series_gp() models.
schema_version, package_version: See below.

The kernel matrix $K$ and the observation covariance $C$ are not stored, because they can be recomputed exactly: evaluate_kernel(model$kernel, model$x) is $K$, and adding diag(model$training_noise_variance) gives $C$. Storing only $R$ keeps a model with $n$ observations at about $8 n^2$ bytes. Other fields, such as centered_y, prior_mean, cholesky_attempts, and numerical_policy, are internal and may change without notice.

Sparse and heteroscedastic models

Sparse models (gaussianprocesses_sparse_model) document x, y, kernel, mean, noise_variance, inducing_points, n_inducing, method ("fitc" or "vfe"), approximation ("FITC" or "VFE"), log_marginal_likelihood (the FITC objective or the VFE bound), trace_term, and, for optimize_sparse_gp() models, optimization; their other fields are internal. Heteroscedastic models (gaussianprocesses_heteroscedastic_model) document x, y, kernel, noise_kernel, training_noise_variance, converged, iterations, and history, and contain two exact models, mean_model and noise_model. Both record schema_version and package_version.

Latent models

Latent models (gaussianprocesses_latent_model) document the fields listed in fit_latent_gp(). Like exact models, they store the Cholesky factor of $B = I + W^{1/2} K W^{1/2}$ but not $K$, and record schema_version and package_version; they were added at schema 3.

State-space models

State-space models (gaussianprocesses_state_space_model), from fit_time_series_gp() with method = "state_space", store the training data, kernel, mean, noise, log marginal likelihood, and the smoothed fitted means and variances, but no covariance matrix or factor: every prediction filters and smooths again. They were added at schema 4.

Schema versions

schema_version identifies the structure of a model, and package_version is the version of gaussianprocesses that created it. The current schema is 4: schema 2 added mean_fit, schema 3 added derivative, and schema 4 added the method and trace_term of sparse models. Schema 4 is the format of version 1.0. Models saved by versions that predate schema versions, which stored $K$ and $C$, are schema 0, and models saved before parametric means were added are schema 1. Every function that accepts a model first converts older schemas to schema 4; the converted model gives the same results, and a converted sparse model is a FITC model. A converted schema-0 model has package_version NA, because the creating version is unknown.

A model with a newer or unknown schema, or without the fields its class requires, raises an error of class gaussianprocesses_model_format_error. The condition object has fields found_schema, supported_schemas, and package_version.

Persistence

Models can be saved with saveRDS() and restored with readRDS(); a restored model gives the same predictions, likelihoods, gradients, diagnostics, and seeded draws as the original.

Every 1.x release reads models saved by any earlier version, back to 0.1.0. A 1.x release may introduce a new schema only together with an upgrade from schema 4, so models saved by version 1.0 stay readable and give the same results. A model saved by a newer version than the one reading it raises gaussianprocesses_model_format_error rather than giving wrong results. The tests read models saved by the 0.1.0 release and schema-4 models written before 1.0, and check that they give the results they were saved with. Keeping the training data with the kernel, mean, and noise also reproduces a model exactly: refitting them with the same numerical settings gives the same model.

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(-2, 2, length.out = 10)
model <- fit_gp(x, sin(x), rbf_kernel(), noise_variance = 0.01)

model$schema_version
#> [1] 4

# C is recomputed from the kernel, and R'R reproduces it.
C <- evaluate_kernel(model$kernel, model$x) +
  diag(model$training_noise_variance)
all.equal(crossprod(model$cholesky), C + diag(model$jitter, nrow(C)))
#> [1] TRUE

# A saved model gives identical predictions.
path <- tempfile(fileext = ".rds")
saveRDS(model, path)
identical(predict_gp(readRDS(path), 0.5), predict_gp(model, 0.5))
#> [1] TRUE