knmodels.RdIt implements the models in Knorr-Held, L. (2000) with three different constraint approaches: sum-to-zero, contrast or diagonal add.
The formula specifying the other model components, without
the spacetime interaction term. The spacetime interaction term will be added
accordly to the specification in the control.st argument. See
inla
If it is to be shown the model fitting progress. Useful if more than one interaction type is being fitted.
Named list of arguments to control the spacetime interaction. It should contain:
to be used as the index set for the main temporal effect which will be considered for the constraints when it is the case.
to be used as the index set for the main spatial effect which will be considered for the constraints when it is the case.
to be the index set for the spacetime interaction effect.
to be the graph for the spatial neighbor
structure to be used in a f() term for the main spatial random
effect term or for building the spacetime interaction model.
to
specify the spacetime interaction type. 1 to 4 corresponds to
the four interaction types in Knorr-Held, L. (2000) with all the needed
sum-to-zero constraints. 2c, 3c and 4c are the
contrast version considering the first time or space constrained to be equal
to zero. 2d, 3d and 4d are the corresponding versions
when considering the diagonal add approach.
to be the value to be added to the diagonal when using the diagonal add approach.
to specify the time point to be the reference time in the contrast parametrization.
to specify the area to be the reference for the contrast parametrization.
where additional
arguments can be passed to f() function. Specification of the
hyperparameter, fixed or random, initial value, prior and its parameters for
the spacetime interaction. See ?inla.models and look for
generic0. By default we scale it and use the PC-prior to set the
prior using the pc.prec prior with param = c(0.5, 0.5). See
documentation with ?inla.doc("pc.prec").
Arguments to be passed to the inla() function.
Environment in which to evaluate the ... arguments.
inla.knmodels returns an object of class "inla". or a
list of objects of this class if it is asked to compute more than one
interaction type at once. Note: when the model type is 2c, 3c, 4c, 2d, 3d
or 4d, it also includes linear combinations summary.
inla.knmodels.sample() to sample from
### define space domain as a grid
grid <- sp::SpatialGrid(sp::GridTopology(c(0,0), c(1, 1), c(4, 5)))
(n <- nrow(xy <- sp::coordinates(grid)))
#> [1] 20
### build a spatial neighborhood list
jj <- lapply(1:n, function(i)
which(sqrt((xy[i,1]-xy[,1])^2 + (xy[i,2]-xy[,2])^2)==1))
### build the spatial adjacency matrix
graph <- sparseMatrix(rep(1:n, sapply(jj, length)),
unlist(jj), x=1, dims=c(n, n))
### some random data at 10 time point
dat <- inla.knmodels.sample(graph, m=10, tau.t=2, tau.s=2, tau.st=3)
str(dat)
#> List of 4
#> $ time : int [1:200] 1 1 1 1 1 1 1 1 1 1 ...
#> $ space : int [1:200] 1 2 3 4 5 6 7 8 9 10 ...
#> $ spacetime: int [1:200] 1 2 3 4 5 6 7 8 9 10 ...
#> $ x :List of 8
#> ..$ t.str : num [1:10] 1.425 0.95 0.481 0.897 -1.222 ...
#> ..$ t.iid : num [1:10] -0.7769 -0.2598 1.2863 2.1116 -0.0256 ...
#> ..$ time : num [1:20] 0.542 0.462 0.783 1.349 -0.733 ...
#> ..$ s.iid : num [1:20] -1.81 -2.14 -1.56 -2.24 -1.12 ...
#> ..$ s.str : num [1:20] 0.64 0.769 1.246 -0.322 0.366 ...
#> ..$ space : num [1:40] -0.324 -0.374 0.132 -1.059 -0.219 ...
#> ..$ spacetime: num [1:200] 0.7681 0.0152 0.2515 -0.8529 0.9289 ...
#> ..$ eta : num [1:200] 0.986 0.184 0.926 -1.37 1.252 ...
sapply(dat$x, summary)
#> t.str t.iid time s.iid s.str
#> Min. -1.221710e+00 -1.97580178 -1.476350342 -2.2431731 -1.911770e+00
#> 1st Qu. -7.039139e-01 -0.64758701 -0.755964829 -1.4034535 -2.510841e-01
#> Median -8.521006e-02 -0.07563295 0.007120603 -0.1948509 8.851493e-03
#> Mean 6.425849e-16 -0.02920444 -0.005655415 -0.2861701 -1.497077e-16
#> 3rd Qu. 7.930931e-01 0.19383505 0.602150032 0.6288812 4.110711e-01
#> Max. 1.424618e+00 2.11162990 1.424617924 2.4299344 1.245502e+00
#> space spacetime eta
#> Min. -1.91176980 -1.910377e+00 -3.31059420
#> 1st Qu. -0.39325099 -2.909498e-01 -1.06841564
#> Median -0.08801616 1.775803e-02 -0.04157377
#> Mean -0.05541659 1.591609e-18 -0.12214402
#> 3rd Qu. 0.41107108 3.769689e-01 0.75948787
#> Max. 1.24550186 1.950664e+00 2.50165651
nd <- length(dat$x$eta)
dat$e <- runif(nd, 0.9, 1.1)*rgamma(n, 40, 2)
dat$y <- rpois(nd, dat$e*exp(dat$x$eta-3))
summary(dat$y)
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 0.00 0.00 1.00 1.43 2.00 12.00
### fit the type 4 considering three different approaches
tgraph <- sparseMatrix(i=c(2:10, 1:9), j=c(1:9, 2:10), x=1)
res <- inla.knmodels(y ~ f(time, model='bym2', graph=tgraph) +
f(space, model='bym2', graph=graph),
data=dat, family='poisson', E=dat$E, progress=TRUE,
control.st=list(time=time, space=space,
spacetime=spacetime, graph=graph, type=c(4, '4c')),
control.compute=list(dic=TRUE, waic=TRUE, cpo=TRUE))
#> r.size.r = 10 20
#> Added model term:
#> f(spacetime, model="generic0", constr=FALSE, Cmatrix=kronecker(R.t, R.s), extraconstr=list(A=rbind(M2,M3)[-1,], e=rep(0,20+10-1)))
#> type = 4, cpu = 11.24455
#> type = 4c, cpu = 11.17978
sapply(res, function(x)
c(dic=x$dic$dic, waic=x$waic$waic, cpo=-sum(log(x$cpo$cpo))))
#> 4 4c
#> dic 569.1414 571.8434
#> waic 575.8232 581.3414
#> cpo 796.1873 734.2524