This is a R package to implement certain spatial and spatio-temporal
models taking use to the cgeneric
interface in the INLA
package. This interface is a way to implement models by writing
C
code to build the precision matrix compiling it so that
INLA can use it internally.
You can install the current CRAN version of INLAspacetime:
install.packages("INLAspacetime")
You can install the latest version of INLAspacetime from GitHub with
## install.packages("remotes")
::install_github("eliaskrainski/INLAspacetime", build_vignettes=TRUE) remotes
some of the models presented in https://www.idescat.cat/sort/sort481/48.1.1.Lindgren-etal.pdf
the barrier model proposed in https://doi.org/10.1016/j.spasta.2019.01.002
This is a basic example which fit a spacetime model for some fake data. The model fitting using inlabru facilitates coding.
set.seed(1)
<- 5
n <- data.frame(
dataf s1 = runif(n, -1, 1),
s2 = runif(n, -1, 1),
time = runif(n, 1, 4),
y = rnorm(n, 0, 1))
str(dataf)
#> 'data.frame': 5 obs. of 4 variables:
#> $ s1 : num -0.469 -0.256 0.146 0.816 -0.597
#> $ s2 : num 0.797 0.889 0.322 0.258 -0.876
#> $ time: num 1.62 1.53 3.06 2.15 3.31
#> $ y : num -0.00577 2.40465 0.76359 -0.79901 -1.14766
Loading the packages:
library(INLA)
#> Loading required package: Matrix
#> Loading required package: sp
#> This is INLA_24.06.04 built 2024-06-03 20:33:29 UTC.
#> - See www.r-inla.org/contact-us for how to get help.
#> - List available models/likelihoods/etc with inla.list.models()
#> - Use inla.doc(<NAME>) to access documentation
library(INLAspacetime)
#> Loading required package: fmesher
library(inlabru)
Define spatial and temporal discretization meshes
<- inla.mesh.2d(
smesh loc = cbind(0,0),
max.edge = 5,
offset = 2)
<- inla.mesh.1d(
tmesh loc = 0:5)
Define the spacetime model object to be used
<- stModel.define(
stmodel smesh = smesh, ## spatial mesh
tmesh = tmesh, ## temporal mesh
model = '121', ## model, see the paper
control.priors = list(
prs = c(1, 0.1), ## P(spatial range < 1) = 0.1
prt = c(5, 0), ## temporal range fixed to 5
psigma = c(1, 0.1) ## P(sigma > 1) = 0.1
) )
Define the data model: the linear predictor terms
<- ~ 1 +
linpred field(list(space = cbind(s1, s2),
time = time),
model = stmodel)
Setting the likelihood
<- list(
ctrlf hyper = list(
prec = list(
initial = 10,
fixed = TRUE)
)
)<- like(
datalike formula = y ~ .,
family = "gaussian",
control.family = ctrlf,
data=dataf)
Fitting
<-
result bru(
components = linpred,
datalike,options = list(
control.inla = list(
int.strategy = "eb"
),verbose = !TRUE)
)
Summary of the model parameters
$summary.fixed
result#> mean sd 0.025quant 0.5quant 0.975quant mode kld
#> Intercept 0.5264782 3.500849 -6.33506 0.5264782 7.388016 0.5264782 0
$summary.hyperpar
result#> mean sd 0.025quant 0.5quant 0.975quant mode
#> Theta1 for field 1.190361 0.4867876 0.3624381 1.153754 2.255723 0.972674
#> Theta2 for field 1.435282 0.1709783 1.1034628 1.433661 1.776664 1.426789
Please check it out at the Tutorials