Package 'adlaplaceFem'

Title: Finite Element Approximations for Adlaplace
Description: Finite-element approximations for hierarchical models with the adlaplace Laplace approximation framework. Provides tensor-product B-spline FEM ingredients for Matern Gaussian random fields (mass, stiffness, and higher-order Grams), including hierarchical locally refined B-spline bases, design and Gram matrix construction on regular grids, optional terra SpatRaster / SpatVector input, AD tape builders for random_fem_ssq_* / random_fem_det_*, and a near-monotone mGP smoother (mgp) with a predictive SD (psd) interface.
Authors: Patrick Brown [aut, cre, cph], Ziang Zhang [aut]
Maintainer: Patrick Brown <[email protected]>
License: GPL (>= 2)
Version: 0.7.3
Built: 2026-10-16 21:42:57 UTC
Source: https://github.com/eborgnine/adlaplace

Help Index


Tensor-product B-spline FEM ingredients for a 2D Matern GRF

Description

Builds the design matrix of tensor-product B-splines at evaluation coordinates together with sparse mass, stiffness, and higher-order Gram matrices used in the SPDE FEM precision for operator orders alpha = 2 and 3.

Usage

fem_bspline(coords, knots, degree = 2L, ...)

## S4 method for signature 'list'
fem_bspline(coords, knots, degree = 2L, ...)

## S4 method for signature 'data.frame'
fem_bspline(coords, knots, degree = 2L, ...)

Arguments

coords

Evaluation locations. For a data.frame (e.g. sites_eval from expand.grid), rows are point coordinates with columns x and y. For a plain list(x, y) of axis sequences, the full tensor grid is formed with expand.grid. With terra, SpatRaster / SpatVector are also accepted.

knots

Knot-line positions on each axis: list(x = seq(...), y = seq(...)) of unique increasing breakpoints (domain endpoints included). These are expanded to open B-spline knot vectors. With terra, a SpatRaster whose grid defines knot lines, or an "hb_knots" object from hb_knots(), is also accepted.

degree

B-spline degree; must be ⁠>= 2⁠. Default 2 (alpha = 2). Use degree = 3 when assembling G3 for alpha = 3.

...

Passed to methods (unused for list / data.frame).

Details

Coefficient ordering is column-major (vec) over the n_x by n_y array of tensor basis indices: basis ⁠(i, j)⁠ maps to column i + (j - 1) * n_x.

Value

A list with:

A

dgCMatrix design, n by n_x * n_y

C, G, G2

sparse Grams

G3

sparse Gram or NULL if degree < 3

degree, knots, n_basis

metadata (n_basis = c(n_x, n_y))

Examples

sites_list <- list(x = seq(0, 1, by = 0.25), y = seq(0, 1, by = 0.25))
sites_eval <- do.call(expand.grid, sites_list)
fem <- fem_bspline(sites_eval, sites_list, degree = 2L)

One-dimensional linear B-spline FEM (mass and stiffness)

Description

One-dimensional linear B-spline FEM (mass and stiffness)

Usage

fem_bspline_1d(knots, degree = 1L)

Arguments

knots

Strictly increasing knot positions, including the domain ends.

degree

B-spline degree. The Matérn-1/2 term uses degree 1.

Value

List with open knots, degree, sparse Grams C and G, and n basis functions.


Symbolic LDL pattern for a sparse FEM precision structure

Description

Builds the permutation and unit-lower pattern required by adlaplace's chol_update for a precision with the same sparsity as Q_structure (typically the union of FEM Gram patterns).

Usage

fem_chol_pattern(Q_structure)

Arguments

Q_structure

Symmetric sparse matrix with the structural nonzeros of Q(theta) (numeric values only need to make the matrix positive definite for symbolic analysis; diagonals are inflated if needed).

Value

A list with perm (0-based), L1 (dtCMatrix / unit lower), and perm_inv (0-based), suitable for the chol slot of a ⁠random_fem_*⁠ precision payload.


Sparse LDL log-determinant of a FEM precision at given (range, sd)

Description

Exercises the same chol_update_csc path used by random_fem_det_*. Useful for diagnosing non-finite joint densities.

Usage

fem_logdet_debug(precision, range, sd, log_scale = NULL)

Arguments

precision

List from fem_precision_payload.

range, sd

Practical range and field SD (natural or log scale).

log_scale

Logical length-1 or length-2: whether range/sd are already on the log scale. Default FALSE (NULL).

Value

List with logdet, coef, D, first_bad, n_bad, n, alpha.


Assemble FEM Matern precision Q2 or Q3 from Grams

Description

Assemble FEM Matern precision Q2 or Q3 from Grams

Usage

fem_precision(kappa, tau, C, G, G2 = NULL, G3 = NULL, alpha = 2L)

Arguments

kappa, tau

Positive SPDE parameters.

C, G, G2, G3

Sparse Grams from fem_bspline().

alpha

1 (1D Matern nu = 1/2: ⁠tau^2 (kappa^2 C + G)⁠), 2, or 3.

Value

Sparse dgCMatrix precision.


Build ⁠random_fem_ssq_*⁠ / ⁠random_fem_det_*⁠ precision payload for adlaplace

Description

Build ⁠random_fem_ssq_*⁠ / ⁠random_fem_det_*⁠ precision payload for adlaplace

Usage

fem_precision_payload(fem, alpha = 2L)

Arguments

fem

Result of fem_bspline() (or list with C, G, G2, optional G3).

alpha

1, 2, or 3. Alpha 1 uses C and G only and stores a zero G2 so the payload layout stays the same.

Value

List for density_data@precision: Grams, chol pattern, and Q CSC coefficients aligned for on-tape assembly.


Structural nonzero pattern for Q2 or Q3 from Grams

Description

Structural nonzero pattern for Q2 or Q3 from Grams

Usage

fem_Q_structure(C, G, G2, G3 = NULL)

Arguments

C, G, G2, G3

Sparse Grams from fem_bspline(); G3 may be NULL for alpha = 2.

Value

Symmetric sparse matrix with the union of Gram nonzero patterns.


Truncated hierarchical B-spline space

Description

Expands hb_knots() breakpoint lines to open knot vectors and builds the THB embedding (active sets and sparse map S) for a B-spline degree.

Usage

hb_basis(knots, degree = 2L)

hb_summary(hb)

Arguments

knots

An "hb_knots" object from hb_knots().

degree

B-spline degree (must be ⁠>= 2⁠).

hb

An "hb_basis" object from hb_basis().

Value

An "hb_basis" object used by fem_bspline() and matern().

Functions

  • hb_summary(): Per-level active basis counts and total hierarchical dof.


Hierarchical tensor-product B-spline knot specifications

Description

Builds nested knot-line locations from a coarse raster box and optional inner refinement extents. Knot lines do not depend on B-spline degree; pass the result to hb_basis() or matern() to form the THB space.

Usage

hb_knots(outer, inner = NULL, fact = 2L)

Arguments

outer

Level-0 box: list(xmin, xmax, ymin, ymax, resolution), or a terra SpatRaster (converted to that list). resolution is dx or c(dx, dy).

inner

NULL (outer only), one extent, or a list of levels. An extent is c(xmin, xmax, ymin, ymax), a terra SpatExtent, or list(xmin, xmax, ymin, ymax).

fact

Integer ⁠>= 2⁠, recycled to length(inner). fact[i] is the cell-size ratio from the previous level, so level-i resolution is outer$resolution / prod(fact[1:i]).

Details

Each element of inner is one refinement level: a single extent, or a list of extents for sibling boxes at that level. list(ext) and list(list(ext)) are the same one-level, one-box specification. list(ext1, ext2) is two levels; siblings are list(list(ext1, ext2)). Extents are snapped to the nested cell grid at that level (spacing outer$resolution / prod(fact[1:i])) and clipped to the parent level.

Value

An "hb_knots" list, one element per level, each with knots and rasters (xmin, xmax, ymin, ymax, resolution).


Slots of a FEM quadratic form inside the joint inverse Hessian

Description

density_data method for hinv_trace_index. One 0-based index into H_inv@x per stored upper-triangle entry of Q, in that CSC order. -1L marks a pair the symbolic inverse does not store. Returns NULL when shard is not a random_fem_ssq_* density.

Usage

## S4 method for signature 'density_data'
hinv_trace_index(shard, H_inv)

Arguments

shard

A density_data shard.

H_inv

Symbolic upper triangle of the joint inverse Hessian.

Value

An integer vector, or NULL when shard is not a random_fem_ssq_* density.

See Also

hinv_trace_index


Posterior mean, SD, and simulations of a matern field on a raster

Description

Evaluates the fitted Matern GRF on eval_grid under the Laplace approximation for the random effects (outer parameters fixed at the MLE). Returns a SpatRaster with layers mean, sd, and optionally sim1 ... simn.

Usage

matern_est(fit, eval_grid, n = 0L, cores = 1L)

Arguments

fit

An "adlaplace_fit" from adlaplace::adlaplace() that includes exactly one matern() term.

eval_grid

A terra SpatRaster whose cell centers are evaluation sites.

n

Number of posterior simulations of the field. Default 0 (mean and sd only).

cores

Number of PSOCK workers used to build the design matrix. Each worker evaluates a contiguous block of cells, including A_fine %*% S for a hierarchical basis. Default 1 (no cluster). The posterior mean, SD, and simulations are formed on the calling process after the blocks are row-bound.

Value

A SpatRaster with n + 2 layers: mean, sd, and sim1...simn when n > 0.


Covariance of a FEM Matérn field at data locations

Description

Evaluates the B-spline design A at the geometry column in data and returns ⁠A Q^{-1} A'⁠. Q is the SPDE precision at the practical range sqrt(8 * nu) / kappa and marginal SD sd, with nu = shape.

Usage

matern_fem_covariance(term, data, range, sd)

Arguments

term

A matern term with cached FEM Grams.

data

Data frame containing the geometry column named in term.

range

Practical Matérn range, in the same units as the coordinates.

sd

Marginal standard deviation of the Matérn field.

Value

A dense n by n covariance matrix, with n = nrow(data).


Matern FEM spatial random field term

Description

Formula term for a 2D Matern GRF via tensor-product B-spline FEM, wired to adlaplaceFem kernels random_fem_ssq_2 / random_fem_det_2 (shape nu = 1) or random_fem_ssq_3 / random_fem_det_3 (nu = 2). The quadratic form is a random shard; the log-determinant is a parameters companion via adlaplace::extra_density(). The first argument is a column name in data (like adlaplace::iid()): that column holds observation locations as a 2-column matrix, WKT text, or HEX WKB points. Knot lines are supplied at construction; the observation design matrix is built later in design() from data.

Usage

matern(
  x,
  knots,
  shape = 1L,
  init = c(1, 1),
  lower = c(1e-09, 1e-09),
  upper = c(Inf, Inf),
  parscale = c(1, 1),
  log = TRUE
)

## S4 method for signature 'matern'
extra_density(term)

## S4 method for signature 'matern'
design(term, data)

## S4 method for signature 'matern'
precision(term, data)

## S4 method for signature 'matern'
random_info(term, data)

## S4 method for signature 'matern'
theta_info(term)

## S4 method for signature 'matern'
beta_info(term, data)

Arguments

x

Name of the geometry / coordinate column in data.

knots

Knot lines as list(x = ..., y = ...), a terra SpatRaster (extent endpoints plus interior cell centers), or an "hb_knots" object from hb_knots().

shape

Matern smoothness nu (1 or 2). Default 1L (SPDE alpha = 2, quadratic B-splines). Use 2 for nu = 2 / cubic B-splines.

init

Initial values for ⁠(range, sd)⁠ on the natural scale, where range = sqrt(8*nu)/kappa is the practical range (same units as coordinates; matches geostatsp/INLA) and sd is the marginal SD of the Matern field. Internally converted to SPDE ⁠(kappa, tau)⁠.

lower, upper

Bounds for ⁠(range, sd)⁠.

parscale

Optimization scales for ⁠(range, sd)⁠.

log

Whether each theta is optimized on the log scale (recycled to ⁠(range, sd)⁠).

term

A matern term.

data

Model data frame (passed by adlaplace generics).

Details

In 2D, B-spline degree and SPDE order alpha equal shape + 1.

Value

A matern term object.

Functions

  • extra_density(matern): Companion parameters density for the FEM log-determinant shard.

  • design(matern): Design matrix of B-spline basis evaluations at coordinates in data[[term@name]].

  • precision(matern): FEM precision payload for ⁠random_fem_*⁠.

  • random_info(matern): Random-effect metadata (one row per basis weight).

  • theta_info(matern): Theta info for practical range = sqrt(8*nu)/kappa and field sd (log scale when term@log is TRUE).

  • beta_info(matern): Beta info (none).


Boundary (null-space slope) companion of an mgp term

Description

One-column design given by the normalized Box-Cox transform of x−x0x - x_0. Default sd = Inf is a flat prior (no precision shard); a finite boundary_sd uses random_diagonal with known SD.

Usage

## S4 method for signature 'mgp_boundary'
design(term, data)

## S4 method for signature 'mgp_boundary'
precision(term, data)

## S4 method for signature 'mgp_boundary'
random_info(term, data)

Arguments

term

An mgp_boundary term.

data

Model data frame.

Functions

  • design(mgp_boundary): Normalized Box-Cox design column.

  • precision(mgp_boundary): Diagonal precision, or NULL when sd = Inf.

  • random_info(mgp_boundary): One random coefficient.

Slots

sd

Prior SD for the boundary coefficient (Inf = flat).

a

Curvature (same as the parent mgp).

shift

Internal shift cs=x0c_s = x_0.


Near-monotone mGP FEM term

Description

Cubic B-spline FEM for the near-monotone Gaussian process Laf=f′′−αf′=σξL_a f = f'' - \alpha f' = \sigma\xi with α(x)=−1/(a x)\alpha(x) = -1/(a\,x). The null space is spanned by {1,m(x)}\{1, m(x)\} where mm is the normalized Box-Cox transform of xx; the FEM weights drop the first two bases so f(x0)=f′(x0)=0f(x_0) = f'(x_0) = 0, and the slope column is returned separately as an mgp_boundary term (like adlaplace::iwp() returning rpoly).

Usage

mgp(
  x,
  a = 2,
  knots,
  ref_value = "left",
  boundary_sd = Inf,
  n_quad = 10,
  init = 0.02,
  lower = 1e-09,
  upper = Inf,
  parscale = 1,
  log = TRUE
)

## S4 method for signature 'mgp'
design(term, data)

## S4 method for signature 'mgp'
precision(term, data)

## S4 method for signature 'mgp'
theta_info(term)

## S4 method for signature 'mgp'
beta_info(term, data)

## S4 method for signature 'mgp'
random_info(term, data)

## S4 method for signature 'mgp'
extra_density(term)

## S4 method for signature 'mgp'
psd(term, h, x = NULL, sd = 1, ...)

Arguments

x

Location(s) on the original scale (default ref_value).

a

Curvature; 2 gives x\sqrt{x}, 1 gives ⁠\log⁠, Inf gives a linear base model (IWP2).

knots

Vector of knot locations on the original scale (required), same convention as adlaplace::iwp(). Must be strictly positive.

ref_value

Reference location: "left" (default), "middle", "right", or a numeric value inside the knot range (snapped to the nearest knot, like adlaplace::iwp()).

boundary_sd

Prior SD for the boundary column (Inf = flat).

n_quad

Gauss-Legendre nodes per knot span.

init, lower, upper, parscale, log

SD hyperparameter controls.

term

An mgp term.

data

Model data frame.

h

Step size(s).

sd

SD hyperparameter scale (default 1).

...

Unused; for S4 generic compatibility.

Details

Knot specification matches adlaplace::iwp(): pass original-scale knot positions (including the domain ends). Knots must be strictly positive. The reference is snapped to the nearest knot. Bases on each side of ref_value are built from the knots on that side.

The prior density uses adlaplace's random_mult kernel with a fixed precision Q (one SD theta). Use adlaplace::psd() / adlaplace::psd_prior() for a PC prior on the predictive SD.

FEM construction follows Ziang Zhang's BayesGP near-monotone FEM, with Gauss-Legendre quadrature and splines::splineDesign() (no fda dependency at runtime).

Value

Named list with the mgp term and its mgp_boundary companion.

Functions

  • design(mgp): Sparse B-spline design at data[[x]].

  • precision(mgp): Fixed random_mult precision payload.

  • theta_info(mgp): One SD theta.

  • beta_info(mgp): No fixed effects.

  • random_info(mgp): One row per FEM weight.

  • extra_density(mgp): No companion log-det density (random_mult includes it).

  • psd(mgp): Predictive SD via Green's-function quadrature.

Slots

a

Curvature parameter (Inf gives IWP2: α=0\alpha = 0, m(u)=um(u)=u).

shift

Internal shift cs=x0c_s = x_0 used in centred coordinates.

n_quad

Gauss-Legendre nodes per knot span.

fem

Cached list from the internal FEM assembler.

two_sided

Logical; both sides of ref_value are active.

See Also

adlaplace::psd(), adlaplace::psd_prior(), mgp_boundary

Examples

knots <- seq(8, 100, by = 4)
mgp_term <- mgp(x = "pm10", a = 2, knots = knots, ref_value = 30)
mgp_term

Replace data-column names in an rsmatern term

Description

Method for adlaplace::replace_vars(). Rewrites the time column, the exposure column (mult), and the grouping variable (by).

Usage

## S4 method for signature 'rsmatern'
replace_vars(x, map)

## S4 method for signature 'mgp'
replace_vars(x, map)

## S4 method for signature 'mgp_boundary'
replace_vars(x, map)

Arguments

x

An rsmatern term.

map

Named character vector: names are current column names, values are replacement column names. Names absent from map are left unchanged.

Value

The term with column references and the label updated.

See Also

adlaplace::replace_vars()


Random-slope Matérn-1/2 (Ornstein-Uhlenbeck) FEM term

Description

Stationary one-dimensional Matérn field with smoothness 1/2, the continuous analogue of AR(1), multiplied by (mult - ref_mult). Groups in by share one ⁠(range, sd)⁠ and have independent fields. The precision is the linear B-spline FEM ⁠Q = tau^2 (kappa^2 C + G)⁠ with practical range rho = 2/kappa and tau = 1 / (sd * sqrt(2 * kappa)).

Usage

rsmatern(
  x,
  mult,
  by,
  knots,
  ref_mult = 0,
  levels = NULL,
  init = c(1, 0.05),
  lower = NULL,
  upper = NULL,
  log = TRUE,
  parscale = NULL
)

## S4 method for signature 'rsmatern'
design(term, data)

## S4 method for signature 'rsmatern'
precision(term, data)

## S4 method for signature 'rsmatern'
random_info(term, data)

## S4 method for signature 'rsmatern'
theta_info(term)

## S4 method for signature 'rsmatern'
beta_info(term, data)

## S4 method for signature 'rsmatern'
extra_density(term)

Arguments

x

Time column.

mult

Exposure column. The basis row is multiplied by data[[mult]] - ref_mult.

by

Grouping column. Levels share ⁠(range, sd)⁠.

knots

Increasing knot positions on the time axis. Linear B-splines.

ref_mult

Reference exposure (default 0).

levels

Optional group levels, in the order used for the blocks.

init

Length-2 starting values ⁠(range, sd)⁠.

lower, upper

Length-2 bounds. The default lower range is the knot spacing.

log

Log-transform the two hyperparameters (default TRUE).

parscale

Optional length-2 step scale.

term

An rsmatern term.

data

Model data frame (passed by adlaplace generics).

Value

An rsmatern term.

Functions

  • design(rsmatern): Design: basis at x, times (mult - ref_mult), written into the block for by.

  • precision(rsmatern): Block-diagonal FEM precision, one block per group.

  • random_info(rsmatern): One metadata row per basis weight and group.

  • theta_info(rsmatern): Hyperparameters ⁠(range, sd)⁠.

  • beta_info(rsmatern): No fixed effects.

  • extra_density(rsmatern): Log-determinant density random_fem_det_1.

Slots

mult

Exposure column multiplied into the basis.

ref_mult

Reference subtracted from mult.

by

Grouping variable; one independent field per level.

fem

Cached 1D Grams from fem_bspline_1d().