Package {SLGP}

SLGP: A package for spatially dependent probability distributions

Description

The SLGP package provides tools for fitting, summarising, predicting, simulating, updating, and visualising Spatial Logistic Gaussian Processes models (SLGP). SLGP models define flexible conditional distributions through finite-rank Gaussian process representations and can be fitted by MAP, Laplace approximation, or MCMC

Main user interface

The recommended workflow is based on standard R generics:

• slgp: fit or initialise an SLGP model.

• summary: summarise a fitted SLGP object.

• plot: visualise fitted conditional densities.

• predict: compute densities, CDFs, quantiles, or moments.

• simulate: draw conditional samples.

• update: refit an existing SLGP object.

• coef, formula, and nobs: standard extractors.

Low-level routines

The package also retains lower-level computational routines such as predictSLGP_newNode, predictSLGP_cdf, predictSLGP_quantiles, predictSLGP_moments, sampleSLGP, and retrainSLGP. These functions are primarily intended for advanced use and backward compatibility. For ordinary use, prefer the standard methods listed above.

Typical workflow

A typical analysis consists of:

1. Fitting an SLGP model with slgp.

2. Inspecting the fitted model with summary and plot.

3. Computing predictive quantities with predict.

4. Drawing samples from the predictive distributions with simulate.

5. Updating the fitted model with update, if required.

Author(s)

Maintainer: Athénaïs Gautier athenais.gautier@onera.fr

References

Gautier, Athénaïs (2023). "Modelling and Predicting Distribution-Valued Fields with Applications to Inversion Under Uncertainty." Thesis, Universität Bern, Bern. See the thesis online at https://boristheses.unibe.ch/4377/

Examples

set.seed(1)
d <- data.frame(
  x = rep(seq(0, 1, length.out = 6), each = 5)
)
d$y <- rnorm(nrow(d), mean = sin(2 * pi * d$x), sd = 0.2)

fit <- slgp(
  y ~ x,
  data = d,
  method = "MAP",
  basisFunctionsUsed = "RFF",
  predictorsLower = 0,
  predictorsUpper = 1,
  responseRange = range(d$y),
  opts_BasisFun = list(nFreq = 20, MatParam = 5 / 2),
  seed = 1
)

summary(fit)

plot(fit, draw = c("mean", 1:5))

The SLGP S4 Class: Spatial Logistic Gaussian Process Model

Description

This S4 class represents a Spatial Logistic Gaussian Process (SLGP) model, designed for modelling conditional or spatially dependent probability distributions. It encapsulates all necessary components for training, sampling, and prediction, including the basis function setup, learned coefficients, and fitted hyperparameters.

Slots

formula

A formula specifying the model structure and covariates.

data

A data.frame containing the observations used to train the model.

responseName

A character string specifying the name of the response variable.

covariateName

A character vector specifying the names of the covariates.

responseRange

A numeric vector of length 2 indicating the lower and upper bounds of the response.

predictorsRange

A list containing:

• predictorsLower: lower bounds of the covariates;

• predictorsUpper: upper bounds of the covariates.

method

A character string indicating the training method used: one of {"MCMC", "MAP", "Laplace", "none"}.

p

An integer indicating the number of basis functions used.

basisFunctionsUsed

A character string specifying the type of basis functions used: "inducing points", "RFF", "Discrete FF", "filling FF", or "custom cosines".

opts_BasisFun

A list of additional options used to configure the basis functions.

BasisFunParam

A list containing the computed parameters of the basis functions, e.g., Fourier frequencies or interpolation weights.

coefficients

A matrix of coefficients for the finite-rank Gaussian process. Each row corresponds to a realization of the latent field: Z(x, t) = \sum_{i=1}^p \epsilon_i f_i(x, t) .

hyperparams

A list of hyperparameters, including:

• sigma: numeric signal standard deviation;

• lengthscale: a vector of lengthscales for each input dimension.

trend

A function that returns the trend of the transformed GP (not to be estimated).

logPost

A numeric value representing the (unnormalized) log-posterior of the model. Currently available only for MAP and Laplace-trained models.

Check basis function parameters

Description

Checks and completes the parameter list for a given basis function type.

Usage

check_basisfun_opts(basisFunctionsUsed, dimension, opts_BasisFun = list())

Arguments

Value

A completed list of options specific to the chosen basis function.

Extract the coefficient draws of an SLGP model

Description

Extract the coefficient draws of an SLGP model

Usage

## S4 method for signature 'SLGP'
coef(object, ...)

Arguments

Value

The matrix of finite-rank GP coefficients (draws in rows, basis functions in columns).

Computes the Euclidean distance between rows of two matrices

Description

Computes the Euclidean distance between rows of two matrices

Usage

crossdist(x, y)

Arguments

Value

Euclidean distance between rows of x and y

Evaluate basis functions at given locations.

Description

Evaluates all basis functions defined by a parameter list at new locations.

Usage

evaluate_basis_functions(parameters, X, lengthscale)

Arguments

Value

A matrix of basis function values.

Model formula of an SLGP object

Description

Model formula of an SLGP object

Usage

## S4 method for signature 'SLGP'
formula(x, ...)

Arguments

Value

The model formula.

Initialize basis function parameters

Description

Initializes the parameter list needed for a basis function.

Usage

initialize_basisfun(
  basisFunctionsUsed,
  dimension,
  lengthscale,
  opts_BasisFun = list()
)

Arguments

Value

A list of initialized basis function parameters.

Initialize parameters basis functions based on Random Fourier Features

Description

Draws parameters for standard RFF approximating a Matérn kernel.

Usage

initialize_basisfun_RFF(dimension, nFreq, MatParam = 5/2, lengthscale)

Arguments

Value

List with frequency, offset, and coefficient parameters.

Initialize discrete Fourier features

Description

Generates basis using discrete cosine/sine terms for each input dimension.

Usage

initialize_basisfun_discreteFF(dimension, maxOrdert, maxOrderx)

Arguments

Value

List with frequency, offset, and coefficient parameters.

Initialize space-filling Random Fourier Features

Description

Initializes RFF parameters with LHS-optimized frequency directions.

Usage

initialize_basisfun_fillingRFF(
  dimension,
  nFreq,
  MatParam = 5/2,
  lengthscale,
  seed = 0
)

Arguments

Value

List with frequency, offset, and coefficient parameters.

Initialize parameters for inducing-point basis functions

Description

Computes kernel matrix and its decompositions for use in inducing-point basis functions.

Usage

initialize_basisfun_inducingpt(
  dimension,
  kernel = "Mat52",
  lengthscale,
  pointscoord = NULL,
  numberPoints = NULL
)

Arguments

Value

List with kernel square root and inverse root matrices, and scaled coordinates.

Number of observations used to fit an SLGP model

Description

Number of observations used to fit an SLGP model

Usage

## S4 method for signature 'SLGP'
nobs(object, ...)

Arguments

Value

Integer number of observations.

normalize_data: Normalize data to the range [0, 1]

Description

Scales the response and covariates of a dataset to the unit interval [0,1]. This normalization is required before applying SLGP methods. If range bounds are not provided, they are computed from the data.

Usage

normalize_data(
  data,
  predictorNames,
  responseName,
  predictorsUpper = NULL,
  predictorsLower = NULL,
  responseRange = NULL
)

Arguments

Value

A normalized data frame with the same column structure as data, with values scaled to [0,1].

Plot the estimated density field of an SLGP model

Description

Visualises the conditional densities implied by a fitted (or prior) SLGP-class model. For a single covariate, the densities are drawn as response curves at a sequence of covariate slices.

Usage

## S4 method for signature 'SLGP,missing'
plot(
  x,
  y,
  newdata = NULL,
  n_slices = 8,
  n_response = 101,
  interpolateBasisFun = "WNN",
  draw = "mean",
  panels = TRUE,
  discrete = FALSE,
  ...
)

Arguments

Details

This method draws with base graphics so that the package gains a plot method without taking on a hard ggplot2 dependency. It returns, the computed prediction data.frame invisibly, so users who prefer ggplot2 can build their own figure from it.

Value

The prediction data.frame, invisibly.

Examples

set.seed(1)
d <- data.frame(x = rep(seq(0, 1, length.out = 6), each = 5))
d$y <- rnorm(nrow(d), mean = sin(2 * pi * d$x), sd = 0.2)

fit <- slgp(
  y ~ x, data = d, method = "Laplace", basisFunctionsUsed = "RFF",
  predictorsLower = 0, predictorsUpper = 1, responseRange = range(d$y),
  opts_BasisFun = list(nFreq = 20, MatParam = 5 / 2),
  seed = 1, opts=list(ndraws=5)
)

plot(fit, draw = "mean")
plot(fit, draw = c("mean", 1:5))

pre_comput_NN: Precompute quantities for SLGP basis evaluation with nearest-neighbour interpolation

Description

Computes intermediate quantities for evaluating SLGP basis functions using Nearest neighbour (NN) interpolation over a regular grid in the normalized domain.

Usage

pre_comput_NN(
  normalizedData,
  predictorNames,
  responseName,
  nIntegral = 101,
  nDiscret = 101
)

Arguments

Value

A list of intermediate quantities used in SLGP evaluation:

• nodes: grid of response × covariates,

• indNodesToIntegral: response bin indices,

• indSamplesToNodes: sample-to-node index mapping,

• weightSamplesToNodes: equal weights for NN interpolation.

pre_comput_WNN: Precompute quantities for SLGP basis evaluation with weighted nearest-neighbours

Description

Computes intermediate quantities for evaluating basis functions via weighted nearest-neighbour (WNN) interpolation on a discretised grid.

Usage

pre_comput_WNN(
  normalizedData,
  predictorNames,
  responseName,
  nIntegral = 101,
  nDiscret = 101,
  mode = "pdf"
)

Arguments

Value

A list of intermediate quantities:

• nodes: all evaluation points in response × covariates grid,

• indNodesToIntegral: indices to map nodes to response bins,

• indSamplesToNodes: index mapping from samples to grid nodes,

• weightSamplesToNodes: interpolation weights using inverse distance.

pre_comput_nothing: Precompute quantities for SLGP basis evaluation without interpolation

Description

Computes intermediate quantities for evaluating basis functions when no interpolation is used. Basis functions are evaluated at the exact covariate and response grid locations.

Usage

pre_comput_nothing(
  normalizedData,
  predictorNames,
  responseName,
  nIntegral = 101
)

Arguments

Value

A list of intermediate quantities used in SLGP basis function computation:

• nodes: all points where basis functions are evaluated,

• indNodesToIntegral: index mapping nodes to response bins,

• indSamplesToNodes: index mapping observations to nodes,

• indSamplesToPredictor: index mapping observations to unique predictors,

• weightSamplesToNodes: interpolation weights (equal to 1 here).

Predict from a fitted SLGP model

Description

Single entry point for all distributional summaries produced by an SLGP-class model. The type argument selects what is returned and dispatches to the corresponding workhorse function:

"density"

conditional density values, newdata must contain both the response and the covariate columns.

"cdf"

conditional CDF values, newdata as above.

"quantiles"

conditional quantiles at probabilities, newdata contains covariate columns only.

"moments"

conditional moments of order(s) power, newdata contains covariate columns only.

Usage

## S4 method for signature 'SLGP'
predict(
  object,
  newdata,
  type = c("density", "cdf", "quantiles", "moments"),
  probs = NULL,
  power = NULL,
  centered = FALSE,
  interpolateBasisFun = "WNN",
  nIntegral = 101,
  nDiscret = 101,
  discrete = FALSE,
  ...
)

Arguments

Value

A data.frame as returned by the corresponding predictSLGP_* function.

Examples

set.seed(1)
d <- data.frame(
  x = rep(seq(0, 1, length.out = 6), each = 5)
)
d$y <- rnorm(nrow(d), mean = sin(2 * pi * d$x), sd = 0.2)

fit <- slgp(
  y ~ x,
  data = d,
  method = "MAP",
  basisFunctionsUsed = "RFF",
  predictorsLower = 0,
  predictorsUpper = 1,
  responseRange = range(d$y),
  opts_BasisFun = list(nFreq = 20, MatParam = 5 / 2),
  seed = 1
)

## Prediction grid for density and CDF evaluations:
grid <- expand.grid(
  y = seq(min(d$y), max(d$y), length.out = 100),
  x = c(0.25, 0.75)
)

## Predict conditional densities
pred_density <- predict(fit, newdata = grid, type = "density")

## Predict conditional cumulative distribution functions
pred_cdf <- predict(fit, newdata = grid, type = "cdf")

## Prediction locations for summaries depending only on the covariates
newx <- data.frame(x = c(0.25, 0.75))

## Predict conditional quantiles
pred_q <- predict(fit, newdata = newx, type = "quantiles",
                  probs = c(0.25, 0.5, 0.75))

## Predict the first two raw moments
pred_m <- predict(fit, newdata = newx, type = "moments", power = c(1, 2))

Predict cumulative distribution values at new locations using a SLGP model

Description

predictSLGP_cdf() is deprecated; use predict(object, type = "cdf") instead.

Computes posterior predictive cumulative distribution function values at specified response and covariate locations for a fitted Spatial Logistic Gaussian Process (SLGP) model.

Usage

predictSLGP_cdf(
  SLGPmodel,
  newNodes,
  interpolateBasisFun = "WNN",
  nIntegral = 101,
  nDiscret = 101,
  discrete = FALSE
)

Arguments

Value

A data frame combining newNodes with columns named cdf_1, cdf_2, ..., containing predictive CDF evaluations for each coefficient draw stored in SLGPmodel.

See Also

predict for the recommended user interface.

Predict conditional moments

Description

predictSLGP_moments() is deprecated; use predict(object, type = "moments") instead.

Computes raw or centered moments of the posterior predictive distributions at specified covariate locations.

Usage

predictSLGP_moments(
  SLGPmodel,
  newNodes,
  power,
  centered = FALSE,
  interpolateBasisFun = "WNN",
  nIntegral = 101,
  nDiscret = 101,
  discrete = FALSE
)

Arguments

Value

A data frame containing the covariates in newNodes, a column power, and columns named mSLGP_1, mSLGP_2, ..., containing predictive moments for each coefficient draw stored in SLGPmodel.

See Also

predict for the recommended user interface.

Predict densities at new covariate locations using a given SLGP model

Description

predictSLGP_newNode() is deprecated; use predict(object, type = "density") instead.

Computes the posterior predictive probability densities at new covariate points using a fitted Spatial Logistic Gaussian Process (SLGP) model.

Usage

predictSLGP_newNode(
  SLGPmodel,
  newNodes,
  interpolateBasisFun = "WNN",
  nIntegral = 101,
  nDiscret = 101,
  normalize = TRUE,
  discrete = FALSE
)

Arguments

Value

A data frame combining newNodes with columns named pdf_1, pdf_2, ..., representing the posterior predictive density for each coefficient draw stored in SLGPmodel.

See Also

predict for the recommended user interface.

Predict conditional quantiles

Description

predictSLGP_quantiles() is deprecated; use predict(object, type = "quantiles") instead.

Computes predictive quantiles at specified covariate locations by numerical inversion of the posterior predictive CDF.

Usage

predictSLGP_quantiles(
  SLGPmodel,
  newNodes,
  probs,
  interpolateBasisFun = "WNN",
  nIntegral = 101,
  nDiscret = 101,
  discrete = FALSE
)

Arguments

Value

A data frame containing the covariates in newNodes, a column probs, and columns named qSLGP_1, qSLGP_2, ..., containing predictive quantiles for each coefficient draw stored in SLGPmodel.

See Also

predict for the recommended user interface.

Print a brief description of an SLGP model

Description

Compact, one-screen summary of a fitted (or prior) SLGP-class object: the model formula, the estimation method, the basis family and its rank, the response range, the data dimensions, and the key hyperparameters. No internal slots are dumped.

Usage

## S4 method for signature 'SLGP'
print(x, ...)

## S4 method for signature 'SLGP'
show(object)

Arguments

Value

x, invisibly.

Examples

set.seed(1)
d <- data.frame(
  x = rep(seq(0, 1, length.out = 6), each = 5)
)
d$y <- rnorm(nrow(d), mean = sin(2 * pi * d$x), sd = 0.2)

fit <- slgp(
  y ~ x,
  data = d,
  method = "none",
  basisFunctionsUsed = "RFF",
  predictorsLower = 0,
  predictorsUpper = 1,
  responseRange = range(d$y),
  opts_BasisFun = list(nFreq = 20, MatParam = 5 / 2),
  seed = 1
)

fit

Print method for SLGP summaries

Description

Print method for SLGP summaries

Usage

## S3 method for class 'summary.SLGP'
print(x, ...)

Arguments

Value

x, invisibly.

Refit an SLGP model

Description

retrainSLGP() is deprecated; use update(object, ...) instead.

Re-estimates an existing SLGP model, optionally with new data, under a chosen estimation method while reusing the existing basis representation and model ranges.

Usage

retrainSLGP(
  SLGPmodel,
  newdata = NULL,
  epsilonStart = NULL,
  method,
  interpolateBasisFun = "WNN",
  nIntegral = 101,
  nDiscret = 101,
  hyperparams = NULL,
  sigmaEstimationMethod = "none",
  seed = NULL,
  opts = list(),
  trend = NULL,
  discrete = FALSE,
  verbose = FALSE
)

Arguments

Value

An updated object of class SLGP-class with retrained coefficients and updated posterior information.

See Also

update for the recommended user interface.

Rosenblatt transform to multivariate Student distribution

Description

Auxiliary function that maps uniform samples in [0, 1]^d to samples from the spectral density of a Matérn kernel (i.e., a multivariate Student distribution).

Usage

rosenblatt_transform_multivarStudent(x, dimension, MatParam = 5/2)

Arguments

Value

A matrix with transformed coordinates following a multivariate Student distribution.

Draw conditional samples

Description

sampleSLGP() is deprecated; use simulate(object) instead.

Draws samples from the posterior predictive distributions at specified covariate locations using inverse transform sampling based on the estimated predictive CDF.

Usage

sampleSLGP(
  SLGPmodel,
  newX,
  n,
  interpolateBasisFun = "WNN",
  nIntegral = 101,
  nDiscret = 101,
  seed = NULL,
  discrete = FALSE
)

Arguments

Value

A data frame containing simulated responses, with the covariate columns from newX and one response column named after the response variable of SLGPmodel.

See Also

simulate for the recommended user interface.

Simulate responses from a fitted SLGP model

Description

This method provides the standard user interface for simulating from a fitted SLGP-class object.

Usage

## S4 method for signature 'SLGP'
simulate(
  object,
  nsim = 1,
  seed = NULL,
  newdata,
  interpolateBasisFun = "WNN",
  nIntegral = 101,
  nDiscret = 101,
  discrete = FALSE,
  ...
)

Arguments

Value

A data.frame of sampled responses with the covariate columns from newdata.

Examples

set.seed(1)
d <- data.frame(
  x = rep(seq(0, 1, length.out = 6), each = 5)
)
d$y <- rnorm(nrow(d), mean = sin(2 * pi * d$x), sd = 0.2)

fit <- slgp(
  y ~ x,
  data = d,
  method = "MAP",
  basisFunctionsUsed = "RFF",
  predictorsLower = 0,
  predictorsUpper = 1,
  responseRange = range(d$y),
  opts_BasisFun = list(nFreq = 20, MatParam = 5 / 2),
  seed = 1
)

## Draw 10 samples from the conditional distributions at two locations
sim <- simulate(
  fit,
  newdata = data.frame(x = c(0.25, 0.75)),
  nsim = 10
)

head(sim)

Fit a Spatial Logistic Gaussian Process model

Description

Builds a finite-rank Spatial Logistic Gaussian Process (SLGP) model for conditional density estimation. The model can be fitted by MAP, MCMC, or Laplace approximation, or initialised without fitting by setting method = "none".

Usage

slgp(
  formula,
  data,
  epsilonStart = NULL,
  method,
  basisFunctionsUsed,
  interpolateBasisFun = "NN",
  nIntegral = 101,
  nDiscret = 101,
  hyperparams = NULL,
  predictorsUpper = NULL,
  predictorsLower = NULL,
  responseRange = NULL,
  sigmaEstimationMethod = "none",
  seed = NULL,
  opts_BasisFun = list(),
  BasisFunParam = NULL,
  opts = list(),
  trend = NULL,
  discrete = FALSE,
  verbose = FALSE
)

Arguments

Value

An object of S4 class SLGP-class. Standard methods are available for fitted objects, including summary, plot, predict, simulate, update, coef, formula, and nobs.

References

Gautier, Athénaïs (2023). "Modelling and Predicting Distribution-Valued Fields with Applications to Inversion Under Uncertainty." Thesis, Universität Bern, Bern. https://boristheses.unibe.ch/4377/

Examples

set.seed(1)
d <- data.frame(
  x = rep(seq(0, 1, length.out = 6), each = 5)
)
d$y <- rnorm(nrow(d), mean = sin(2 * pi * d$x), sd = 0.2)

fit <- slgp(
  y ~ x,
  data = d,
  method = "none",
  basisFunctionsUsed = "RFF",
  predictorsLower = 0,
  predictorsUpper = 1,
  responseRange = range(d$y),
  opts_BasisFun = list(nFreq = 20, MatParam = 5 / 2),
  seed = 1
)

fit
summary(fit)

Summarise a fitted SLGP model

Description

Extends print with fit diagnostics where available: the log-posterior at the mode (MAP / Laplace), the number of coefficient draws, and basic structural information. The returned object has class "summary.SLGP" and its own print method.

Usage

## S4 method for signature 'SLGP'
summary(object, ...)

Arguments

Value

An object of class "summary.SLGP" (a list), returned invisibly.

Re-fit an SLGP model under a new method or with new data

Description

Standard-generic front end to retrainSLGP: re-estimates an existing SLGP-class model, optionally with new data, under a chosen estimation method, reusing the existing basis and ranges.

Usage

## S4 method for signature 'SLGP'
update(
  object,
  method,
  newdata = NULL,
  epsilonStart = NULL,
  interpolateBasisFun = "WNN",
  nIntegral = 101,
  nDiscret = 101,
  discrete = FALSE,
  hyperparams = NULL,
  sigmaEstimationMethod = "none",
  seed = NULL,
  opts = list(),
  trend = NULL,
  verbose = FALSE,
  ...
)

Arguments

Value

An updated object of class SLGP-class.

See Also

retrainSLGP for the low-level routine.

Examples

set.seed(1)
d <- data.frame(
  x = rep(seq(0, 1, length.out = 6), each = 5)
)
d$y <- rnorm(nrow(d), mean = sin(2 * pi * d$x), sd = 0.2)

prior <- slgp(
  y ~ x,
  data = d,
  method = "none",
  basisFunctionsUsed = "RFF",
  predictorsLower = 0,
  predictorsUpper = 1,
  responseRange = range(d$y),
  opts_BasisFun = list(nFreq = 20, MatParam = 5 / 2),
  seed = 1
)

## Refit the same model structure by MAP
fit_map <- update(prior, method = "MAP")