inla.stack {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:

- [n_l] effect blocks

- [n_k] effects

- [n_i] data values

- [n_jl] effect size for block l

- [n_j] sum_l n_jl total effect size

Input:

- [`data`] p vectors, each of length n_i

- [`A`] matrices of size n_i by n_jl

- [`effects`] 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} )

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:

- `data` =(y^1, …, y^p, \tilde{x}^1, …, \tilde{x}^{n_k})

- `A` =\tilde{A}

- `data.names` List of data names, length p

- `effect.names` List of effect names, length n_k

- `n.data` Data length, n_i

- `index` List indexed by `tag`s, each element indexing into i=1, …, n_i

### Functions

- `inla.stack.remove.unused`: Remove unused entries from an existing stack

- `inla.stack.compress`: Compress an existing stack by removing duplicates

- `inla.stack`: Shorthand for inla.stack.join and inla.stack.sum

- `inla.stack.sum`: Create data stack as a sum of predictors

- `inla.stack.join`: Join two or more data stacks

- `inla.stack.index`: Extract tagged indices

- `inla.stack.LHS`: Extract data associated with the "left hand side" of the model (e.g. the data itself, `Ntrials`, `link`, `E`)

- `inla.stack.RHS`: Extract data associated with the "right hand side" of the model (all the covariates/predictors)

- `inla.stack.data`: Extract data for an inla call, and optionally join with other variables

- `inla.stack.A`: Extract the "A matrix" for control.predictor

### 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 18.07.12 Index]