IMRF: Interpolated Markov Random Field models
In spaMM: Mixed-Effect Models, with or without Spatial Random Effects
multIMRF R Documentation
Interpolated Markov Random Field models
Description
spaMM can fit random-effect terms of the forms considered by Lindgren et al. (2011) or Nychka et al. (2015, 2018).
The random effects considered here all involve a multivariate Gaussian random effect over a lattice, from which the random-effect value in any spatial position is determined by interpolation of values on the lattice. IMRF stands for Interpolated Markov Random Field because the specific process considered on the lattice is currently known as a Gaussian Markov Random Field (see the Details for further information). Lindgren et al. considered irregular lattices designed to approximate of the Matern correlation model with fixed smoothness <= 2, while Nychka et al. considered regular grids.
The correlation model of Lindgren et al. (2011) can be fitted by spaMM by declaring an IMRF random effect term in the model formula, with a model argument in the right-hand side whose value is the result of INLA::inla.spde2.matern (or INLA::inla.spde2.pcmatern) for given smoothness. The spaMM functions for such a fit do not call INLA functions. Alternatively, the same model with variable smoothness can be fitted by declaring a corrFamily term whose structure is described through the MaternIMRFa function, whose respective documentations should be considered for more details. In the latter case INLA::inla.spde2.matern is called internally by spaMM. The correlation models thus defined are fitted by the same methods as other models in spaMM.
Regular lattices can also be declared as an IMRF term (with arguments distinct from model). The multIMRF syntax implements the multiresolution model of Nychka et al. Any multIMRF term in a formula is immediately converted to IMRF terms for regular grids wih different step sizes. This has distinct implications for controlling the parameters of these or other random effects in the model by init or fixed values: see Details if you need such control.
Usage
# IMRF( 1 | <coordinates>, model, nd, m, no, ce, ...)
# multIMRF( 1 | <coordinates>, levels, margin, coarse=10L,
# norm=TRUE, centered=TRUE )
Arguments
model
An inla.spde2 object as produced by INLA::inla.spde2.matern or
INLA::inla.spde2.pcmatern (see Examples below, and https://www.r-inla.org for further information).
levels
integer;
Number of levels in the hierarchy, i.e. number of component IMRFs.
margin, m
integer;
width of the margin, as a number of additional grid points on each side (applies to all levels of the hierarchy).
coarse
integer;
number of grid points (excluding the margins) per dimension for the coarsest IMRF. The number of grids steps nearly doubles with each level of the hierarchy (see Details).
nd
integer;
number of grid steps (excluding the margins) per dimension for the given IMRF.
norm, no
Boolean; whether to apply normalization (see Details), or not.
centered, ce
Boolean; whether to center the grid in all dimensions, or not.
...
Not documented, for programming purposes
Details
Formulation of the covariance models:
Gaussian Markov Random Field (MRF) and conditional autoregressive models are essentially the same thing, apart from details of specification.
adjacency and AR1 random effects can be seen as specific MRFs.
The common idea is the Markov-like property that the distribution of each element b_i of the random-effect b, given values of a few specific elements (the “neighbours” of i), is independent of other elements (i.e., of non-neighbours). The non-zero non-diagonal elements of a precision matrix characterize the neighbours.
Given the inferred vector b of values of the MRF on the lattice, the interpolation of the MRF in any focal point is of the form Ab where each row of A weights the values of b according to the position of the focal point relative to the vertices of the lattice. Following the original publications,
* for grids given by model=<inla.spde2 object>, the non-zero weights are the barycentric coordinates of the focal point in the enclosing triangle from the mesh triangulation (points from outside the mesh would have zero weights, so the predicted effect Ab=0);
* for regular grids (NULL model), the weights are computed as <Wendland function>(<scaled Euclidean distances between focal point and vertices>).
The IMRF model defines both a lattice in space, the precision matrix for a Gaussian MRF over this lattice, and the A matrix of weights. The full specification of the MRF on irregular lattices is complex. The \kappa (kappa) parameter considered by spaMM is the \kappa scale parameter considered by Lindgren et al and comparable to the \rho scale factor of the Matérn model. The \alpha argument of the INLA::inla.spde2.matern controls the smoothness of the approximated Matern model, as \alpha=\nu + d/2) where d is the dimension of the space.
Correlation models created by INLA::inla.spde2.pcmatern are handled so as to give the same correlation values as when INLA::inla.spde2.matern is used with the same mesh and alpha argument (thus, the extra functionalities of “pc”matern are ignored).
Not all options of the INLA functions may be compatible or meaningful when used with spaMM (only the effects of alpha and cutoff have been checked).
Normalization:
For the MRFs on default regular grids (missing model argument), the precision matrix is defined (up to a variance parameter) as M'M where the diagonal elements m_{ii} of M are 4+\kappa^2 and the m_{ij} for the four nearest neighbours are -1 (note that M'M involves more than these four neighbours). The precision matrix defined in this way is the inverse of an heteroscedastic covariance matrix C, but (following Nychka et al.) by default a normalization is applied so that the random effect in each data position is homoscedastic (the precision matrix for the latent effect in grid positions is not modified, but it is the A matrix of weights which is is modified). As for other random effects, the variance is further controlled by a multiplicative factor \lambda.
Without normalization, the covariance matrix of the random effect in data locations is \lambdaALL'A' (A being the above-described weight matrix, and L is a “square root” of C), and AL is the original “design matrix” of the random effect. \lambda may then be quite different from the marginal variance of the random effect, and is difficult to describe in a simple way.
For normalization, A is modified as WA where W is a diagonal matrix such that WAL is a correlation matrix (WALL'A'W' has unit diagonal); then, \lambda is the marginal variance of the random effect.
For irregular grids specified using the model argument, the precision matrix described by this object is also the inverse of an heteroscedastic covariance matrix, but here (again following original publicatiosn such as Lindgren at al. 2011) the normalization is not applied by default (and was not even an option before version 4.3.23). But for ease of presentation and interpretation, if for no other reason, the normalized model may be preferable.
Details for rectangular grids:
By default (meaning in particular that model is not used to specify a lattice defined by the INLA procedures), the IMRF lattice is rectangular (currently the only option) and is made of a core lattice, to which margins of margin steps are added on each side. The core lattice is defined as follows: in each of the two spatial dimensions, the range of axial coordinates is determined. The largest range is divided in nd-1 steps, determining nd values and step length L. The other range is divided in steps of the same length L. If it extends over (say) 2.5 L, a grid of 2 steps and 3 values is defined, and by default centered on the range (the extreme points therefore typically extend slightly beyond the grid, within the first of the additional steps defined by the margin; if not centered, the grid start from the lower coordinate of the range).
multIMRF implements multilevel IMRFs. It defines a sequence of IMRFs, with progressively finer lattices, a common \kappa value hy_kap for these IMRFs, and a single variance parameter hy_lam that determines \lambda values decreasing by a factor of 4 for successive IMRF terms. By default, each component IMRF is normalized independently as described above (as in Nychka et al. 2019), and hy_lam is the sum of the variances of these terms (e.g., if there are three levels and hy_lam=1, the successive variances are (1,1/4,1/16)/(21/16) ). The nd of the first IMRF is set to the coarse value, and its lattice is defined accordingly. If coarse=4 and margin=5, a grid of 14 coordinates is therefore defined over the largest range. In the second IMRF, the grid spacing is halved, so that new steps are defined halfway between the previous ones (yielding a grid of 27 step in the widest range). The third IMRF proceeds from the second in the same way, and so on.
To control initial or fixed values of multIMRF \kappa and variance parameters, which are hyper-parameter controlling several IMRF terms, the hyper syntax shown in the Examples should be used. hyper is a nested list whose possible elements are named "1", "2", ... referring to successive multIMRF terms in the input formula, not to successive random effect in the expanded formula with distinct IMRF terms (see Examples). But the different IMRF terms should be counted as distinct random effects when controlling other parameters (e.g., for fixing the variances of other random effects).
References
D. Nychka, S. Bandyopadhyay, D. Hammerling, F. Lindgren, S. Sain (2015)
A multiresolution gaussian process model for the analysis of large spatial datasets.
Journal of Computational and Graphical Statistics 24 (2), 579-599. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1080/10618600.2014.914946")}
D. Nychka, D. Hammerling, Mitchel. Krock, A. Wiens (2018) Modeling and emulation of nonstationary Gaussian fields. Spat. Stat. 28: 21-38. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1016/j.spasta.2018.08.006")}
Lindgren F., Rue H., Lindström J. (2011) An explicit link between Gaussian fields and Gaussian Markov random fields: the stochastic partial differential equation approach Journal of the Royal Statistical Society: Series B (Statistical Methodology), 73: 423-498. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1111/j.1467-9868.2011.00777.x")}
Examples
if (spaMM.getOption("example_maxtime")>6) {
data("blackcap") ## toy examples; but IMRF may be useful only for much larger datasets
## and when using the 'cutoff' parameter of INLA::inla.mesh.2d()
########################## Irregular lattice specified by 'model':
#
data("small_spde") ## load object of class 'inla.spde2', created and saved by :
# spd <- sp::SpatialPointsDataFrame(coords = blackcap[, c("longitude", "latitude")],
# data = blackcap)
# small_mesh <- INLA::inla.mesh.2d(loc = INLA::inla.mesh.map(sp::coordinates(spd)),
# max.n=100, # only for demonstration purposes
# max.edge = c(3, 20))
# small_spde <- INLA::inla.spde2.matern(small_mesh)
# save(small_spde, file="small_spde.RData", version=2)
#
fit_SPDE <- fitme(migStatus ~ means + IMRF(1|longitude+latitude, model=small_spde),
data=blackcap)
########################## Regular lattices:
#
#Using 'hyper' to control fixed hyper-parameters
#
(mrf <- fitme(migStatus ~ 1 + (1|pos) +
multIMRF(1|longitude+latitude,margin=5,levels=2),
data=blackcap, fixed=list(phi=1,lambda=c("1"=0.5),
hyper=list("1"=list(hy_kap=0.1,hy_lam=1)))) )
# Using 'hyper' to control initial hyper-parameters
#
(mrf <- fitme(migStatus ~ 1 + multIMRF(1|longitude+latitude,margin=5,levels=2),
data=blackcap, method="ML", fixed =list(phi=1),
init=list(hyper=list("1"=list(hy_kap=0.1,hy_lam=1)))) )
# *Independent* IMRF terms with default rectangular lattice (often giving dubious results)
#
(mrf <- fitme(migStatus ~ 1 + IMRF(1|longitude+latitude,margin=5, nd=4L)
+ IMRF(1|longitude+latitude,margin=5, nd=7L),
data=blackcap,
fixed=list(phi=1,lambda=c(1/4,1/16),
corrPars=list("1"=list(kappa=0.1),"2"=list(kappa=0.1)))))
}
spaMM documentation built on Aug. 30, 2023, 1:07 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
| multIMRF | R Documentation |
Interpolated Markov Random Field models
Description
spaMM can fit random-effect terms of the forms considered by Lindgren et al. (2011) or Nychka et al. (2015, 2018). The random effects considered here all involve a multivariate Gaussian random effect over a lattice, from which the random-effect value in any spatial position is determined by interpolation of values on the lattice. IMRF stands for Interpolated Markov Random Field because the specific process considered on the lattice is currently known as a Gaussian Markov Random Field (see the Details for further information). Lindgren et al. considered irregular lattices designed to approximate of the Matern correlation model with fixed smoothness <= 2, while Nychka et al. considered regular grids.
The correlation model of Lindgren et al. (2011) can be fitted by spaMM by declaring an IMRF random effect term in the model formula, with a model argument in the right-hand side whose value is the result of INLA::inla.spde2.matern (or INLA::inla.spde2.pcmatern) for given smoothness. The spaMM functions for such a fit do not call INLA functions. Alternatively, the same model with variable smoothness can be fitted by declaring a corrFamily term whose structure is described through the MaternIMRFa function, whose respective documentations should be considered for more details. In the latter case INLA::inla.spde2.matern is called internally by spaMM. The correlation models thus defined are fitted by the same methods as other models in spaMM.
Regular lattices can also be declared as an IMRF term (with arguments distinct from model). The multIMRF syntax implements the multiresolution model of Nychka et al. Any multIMRF term in a formula is immediately converted to IMRF terms for regular grids wih different step sizes. This has distinct implications for controlling the parameters of these or other random effects in the model by init or fixed values: see Details if you need such control.
Usage
# IMRF( 1 | <coordinates>, model, nd, m, no, ce, ...)
# multIMRF( 1 | <coordinates>, levels, margin, coarse=10L,
# norm=TRUE, centered=TRUE )
Arguments
model |
An |
levels |
integer; Number of levels in the hierarchy, i.e. number of component IMRFs. |
margin, m |
integer; width of the margin, as a number of additional grid points on each side (applies to all levels of the hierarchy). |
coarse |
integer; number of grid points (excluding the margins) per dimension for the coarsest IMRF. The number of grids steps nearly doubles with each level of the hierarchy (see Details). |
nd |
integer; number of grid steps (excluding the margins) per dimension for the given IMRF. |
norm, no |
Boolean; whether to apply normalization (see Details), or not. |
centered, ce |
Boolean; whether to center the grid in all dimensions, or not. |
... |
Not documented, for programming purposes |
Details
Formulation of the covariance models:
Gaussian Markov Random Field (MRF) and conditional autoregressive models are essentially the same thing, apart from details of specification.
adjacency and AR1 random effects can be seen as specific MRFs.
The common idea is the Markov-like property that the distribution of each element b_i of the random-effect b, given values of a few specific elements (the “neighbours” of i), is independent of other elements (i.e., of non-neighbours). The non-zero non-diagonal elements of a precision matrix characterize the neighbours.
Given the inferred vector b of values of the MRF on the lattice, the interpolation of the MRF in any focal point is of the form Ab where each row of A weights the values of b according to the position of the focal point relative to the vertices of the lattice. Following the original publications,
* for grids given by model=<inla.spde2 object>, the non-zero weights are the barycentric coordinates of the focal point in the enclosing triangle from the mesh triangulation (points from outside the mesh would have zero weights, so the predicted effect Ab=0);
* for regular grids (NULL model), the weights are computed as <Wendland function>(<scaled Euclidean distances between focal point and vertices>).
The IMRF model defines both a lattice in space, the precision matrix for a Gaussian MRF over this lattice, and the A matrix of weights. The full specification of the MRF on irregular lattices is complex. The \kappa (kappa) parameter considered by spaMM is the \kappa scale parameter considered by Lindgren et al and comparable to the \rho scale factor of the Matérn model. The \alpha argument of the INLA::inla.spde2.matern controls the smoothness of the approximated Matern model, as \alpha=\nu + d/2) where d is the dimension of the space.
Correlation models created by INLA::inla.spde2.pcmatern are handled so as to give the same correlation values as when INLA::inla.spde2.matern is used with the same mesh and alpha argument (thus, the extra functionalities of “pc”matern are ignored).
Not all options of the INLA functions may be compatible or meaningful when used with spaMM (only the effects of alpha and cutoff have been checked).
Normalization:
For the MRFs on default regular grids (missing model argument), the precision matrix is defined (up to a variance parameter) as M'M where the diagonal elements m_{ii} of M are 4+\kappa^2 and the m_{ij} for the four nearest neighbours are -1 (note that M'M involves more than these four neighbours). The precision matrix defined in this way is the inverse of an heteroscedastic covariance matrix C, but (following Nychka et al.) by default a normalization is applied so that the random effect in each data position is homoscedastic (the precision matrix for the latent effect in grid positions is not modified, but it is the A matrix of weights which is is modified). As for other random effects, the variance is further controlled by a multiplicative factor \lambda.
Without normalization, the covariance matrix of the random effect in data locations is \lambdaALL'A' (A being the above-described weight matrix, and L is a “square root” of C), and AL is the original “design matrix” of the random effect. \lambda may then be quite different from the marginal variance of the random effect, and is difficult to describe in a simple way.
For normalization, A is modified as WA where W is a diagonal matrix such that WAL is a correlation matrix (WALL'A'W' has unit diagonal); then, \lambda is the marginal variance of the random effect.
For irregular grids specified using the model argument, the precision matrix described by this object is also the inverse of an heteroscedastic covariance matrix, but here (again following original publicatiosn such as Lindgren at al. 2011) the normalization is not applied by default (and was not even an option before version 4.3.23). But for ease of presentation and interpretation, if for no other reason, the normalized model may be preferable.
Details for rectangular grids:
By default (meaning in particular that model is not used to specify a lattice defined by the INLA procedures), the IMRF lattice is rectangular (currently the only option) and is made of a core lattice, to which margins of margin steps are added on each side. The core lattice is defined as follows: in each of the two spatial dimensions, the range of axial coordinates is determined. The largest range is divided in nd-1 steps, determining nd values and step length L. The other range is divided in steps of the same length L. If it extends over (say) 2.5 L, a grid of 2 steps and 3 values is defined, and by default centered on the range (the extreme points therefore typically extend slightly beyond the grid, within the first of the additional steps defined by the margin; if not centered, the grid start from the lower coordinate of the range).
multIMRF implements multilevel IMRFs. It defines a sequence of IMRFs, with progressively finer lattices, a common \kappa value hy_kap for these IMRFs, and a single variance parameter hy_lam that determines \lambda values decreasing by a factor of 4 for successive IMRF terms. By default, each component IMRF is normalized independently as described above (as in Nychka et al. 2019), and hy_lam is the sum of the variances of these terms (e.g., if there are three levels and hy_lam=1, the successive variances are (1,1/4,1/16)/(21/16) ). The nd of the first IMRF is set to the coarse value, and its lattice is defined accordingly. If coarse=4 and margin=5, a grid of 14 coordinates is therefore defined over the largest range. In the second IMRF, the grid spacing is halved, so that new steps are defined halfway between the previous ones (yielding a grid of 27 step in the widest range). The third IMRF proceeds from the second in the same way, and so on.
To control initial or fixed values of multIMRF \kappa and variance parameters, which are hyper-parameter controlling several IMRF terms, the hyper syntax shown in the Examples should be used. hyper is a nested list whose possible elements are named "1", "2", ... referring to successive multIMRF terms in the input formula, not to successive random effect in the expanded formula with distinct IMRF terms (see Examples). But the different IMRF terms should be counted as distinct random effects when controlling other parameters (e.g., for fixing the variances of other random effects).
References
D. Nychka, S. Bandyopadhyay, D. Hammerling, F. Lindgren, S. Sain (2015) A multiresolution gaussian process model for the analysis of large spatial datasets. Journal of Computational and Graphical Statistics 24 (2), 579-599. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1080/10618600.2014.914946")}
D. Nychka, D. Hammerling, Mitchel. Krock, A. Wiens (2018) Modeling and emulation of nonstationary Gaussian fields. Spat. Stat. 28: 21-38. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1016/j.spasta.2018.08.006")}
Lindgren F., Rue H., Lindström J. (2011) An explicit link between Gaussian fields and Gaussian Markov random fields: the stochastic partial differential equation approach Journal of the Royal Statistical Society: Series B (Statistical Methodology), 73: 423-498. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1111/j.1467-9868.2011.00777.x")}
Examples
if (spaMM.getOption("example_maxtime")>6) {
data("blackcap") ## toy examples; but IMRF may be useful only for much larger datasets
## and when using the 'cutoff' parameter of INLA::inla.mesh.2d()
########################## Irregular lattice specified by 'model':
#
data("small_spde") ## load object of class 'inla.spde2', created and saved by :
# spd <- sp::SpatialPointsDataFrame(coords = blackcap[, c("longitude", "latitude")],
# data = blackcap)
# small_mesh <- INLA::inla.mesh.2d(loc = INLA::inla.mesh.map(sp::coordinates(spd)),
# max.n=100, # only for demonstration purposes
# max.edge = c(3, 20))
# small_spde <- INLA::inla.spde2.matern(small_mesh)
# save(small_spde, file="small_spde.RData", version=2)
#
fit_SPDE <- fitme(migStatus ~ means + IMRF(1|longitude+latitude, model=small_spde),
data=blackcap)
########################## Regular lattices:
#
#Using 'hyper' to control fixed hyper-parameters
#
(mrf <- fitme(migStatus ~ 1 + (1|pos) +
multIMRF(1|longitude+latitude,margin=5,levels=2),
data=blackcap, fixed=list(phi=1,lambda=c("1"=0.5),
hyper=list("1"=list(hy_kap=0.1,hy_lam=1)))) )
# Using 'hyper' to control initial hyper-parameters
#
(mrf <- fitme(migStatus ~ 1 + multIMRF(1|longitude+latitude,margin=5,levels=2),
data=blackcap, method="ML", fixed =list(phi=1),
init=list(hyper=list("1"=list(hy_kap=0.1,hy_lam=1)))) )
# *Independent* IMRF terms with default rectangular lattice (often giving dubious results)
#
(mrf <- fitme(migStatus ~ 1 + IMRF(1|longitude+latitude,margin=5, nd=4L)
+ IMRF(1|longitude+latitude,margin=5, nd=7L),
data=blackcap,
fixed=list(phi=1,lambda=c(1/4,1/16),
corrPars=list("1"=list(kappa=0.1),"2"=list(kappa=0.1)))))
}
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.