inla.stack.remove.unused {INLA}R Documentation

Data stacking for advanced INLA models

Description

Functions for combining data, effects and observation matrices into inla.stack objects, and extracting information from such objects.

Usage

inla.stack.remove.unused(stack)

inla.stack.compress(stack, remove.unused = TRUE)

inla.stack(..., compress = TRUE, remove.unused = TRUE)

inla.stack.sum(
  data,
  A,
  effects,
  tag = "",
  compress = TRUE,
  remove.unused = TRUE
)

inla.stack.join(..., compress = TRUE, remove.unused = TRUE)

inla.stack.index(stack, tag)

inla.stack.LHS(stack)

inla.stack.RHS(stack)

inla.stack.data(stack, ...)

inla.stack.A(stack)

Arguments

stack

A inla.data.stack object, created by a call to inla.stack, inla.stack.sum, or inla.stack.join.

remove.unused

If TRUE, compress the model by removing rows of effects corresponding to all-zero columns in the A matrix (and removing those columns).

...

For inla.stack.join, two or more data stacks of class inla.data.stack, created by a call to inla.stack, inla.stack.sum, or inla.stack.join. For inla.stack.data, a list of variables to be joined with the data list.

compress

If TRUE, compress the model by removing duplicated rows of effects, replacing the corresponding A-matrix columns with a single column containing the sum.

data

A list or codedata.frame of named data vectors. Scalars are expanded to match the number of rows in the A matrices, or any non-scalar data vectors. An error is given if the input is inconsistent.

A

A list of observation matrices. Scalars are expanded to diagonal matrices matching the effect vector lengths. An error is given if the input is inconsistent or ambiguous.

effects

A collection of effects/predictors. Each list element corresponds to an observation matrix, and must either be a single vector, a list of vectors, or a data.frame. Single-element effect vectors are expanded to vectors matching the number of columns in the corresponding A matrix. An error is given if the input is inconsistent or ombiguous.

tag

A string specifying a tag for later identification.

Details

For models with a single effects collection, the outer list container for A and effects may be omitted.

Component size definitions:

Input:

data

(y1,…,y2) p vectors, each of length n_i

A

(A1,…,A2) matrices of size n_i by n_jl

effects

((x_[1,1],…,x_[n_k,1]),…(x_[1,n_l],…,x_[n_k,n_l])) collections of effect vectors of length n_jl

predictor(y^1, …, y^p) ~ sum_{l=1}^{n_l} A^l sum_{k=1}^{n_k} g(k, x^{k,l}) = tilde{A} sum_{k=1}^{n_k} g(k, tilde{x}^k)

where

tilde{A} = cbind( A^1, ..., A^{n_l} )

and

tilde{x}^k = rbind( x^{k,1}, ..., x^{k,n_l} )

and for each block l, any missing x^{k,l} is replaced by an NA vector.

Value

A data stack of class inla.data.stack. Elements:

Functions

Functions

See Also

inla.spde.make.A(), inla.spde.make.index()

Examples


n <- 200
loc <- matrix(runif(n * 2), n, 2)
mesh <- inla.mesh.2d(
    loc.domain = loc,
    max.edge = c(0.05, 0.2)
)
proj.obs <- inla.mesh.projector(mesh, loc = loc)
proj.pred <- inla.mesh.projector(mesh, loc = mesh$loc)
spde <- inla.spde2.pcmatern(mesh,
    prior.range = c(0.01, 0.01),
    prior.sigma = c(10, 0.01)
)

covar <- rnorm(n)
field <- inla.qsample(n = 1, Q = inla.spde.precision(spde, theta = log(c(0.5, 1))))[, 1]
y <- 2 * covar + inla.mesh.project(proj.obs, field)

A.obs <- inla.spde.make.A(mesh, loc = loc)
A.pred <- inla.spde.make.A(mesh, loc = proj.pred$loc)
stack.obs <-
    inla.stack(
        data = list(y = y),
        A = list(A.obs, 1),
        effects = list(c(
            list(Intercept = 1),
            inla.spde.make.index("spatial", spde$n.spde)
        ),
        covar = covar
        ),
        tag = "obs"
    )
stack.pred <-
    inla.stack(
        data = list(y = NA),
        A = list(A.pred),
        effects = list(c(
            list(Intercept = 1),
            inla.spde.make.index("spatial", mesh$n)
        )),
        tag = "pred"
    )
stack <- inla.stack(stack.obs, stack.pred)

formula <- y ~ -1 + Intercept + covar + f(spatial, model = spde)
result1 <- inla(formula,
    data = inla.stack.data(stack.obs, spde = spde),
    family = "gaussian",
    control.predictor = list(
        A = inla.stack.A(stack.obs),
        compute = TRUE
    )
)

plot(y, result1$summary.fitted.values[inla.stack.index(stack.obs, "obs")$data, "mean"],
    main = "Observations vs posterior predicted values at the data locations"
)

result2 <- inla(formula,
    data = inla.stack.data(stack, spde = spde),
    family = "gaussian",
    control.predictor = list(
        A = inla.stack.A(stack),
        compute = TRUE
    )
)

field.pred <- inla.mesh.project(
    proj.pred,
    result2$summary.fitted.values[inla.stack.index(stack, "pred")$data, "mean"]
)
field.pred.sd <- inla.mesh.project(
    proj.pred,
    result2$summary.fitted.values[inla.stack.index(stack, "pred")$data, "sd"]
)

plot(field, field.pred, main = "True vs predicted field")
abline(0, 1)
image(inla.mesh.project(mesh,
    field = field,
    dims = c(200, 200)
),
main = "True field"
)
image(inla.mesh.project(mesh,
    field = field.pred,
    dims = c(200, 200)
),
main = "Posterior field mean"
)
image(inla.mesh.project(mesh,
    field = field.pred.sd,
    dims = c(200, 200)
),
main = "Prediction standard deviation"
)
plot(field, (field.pred - field) / 1,
    main = "True field vs standardised prediction residuals"
)

[Package INLA version 21.11.22 Index]