Larger values of the posterior standard deviation at the mesh nodes?

[Vera Tibot]

Hi,

I am running a spatial error model but I getting larger values of the posterior standard deviation at the mesh nodes compared to their surroundings (see attached file). Can please someone explain to me this behaviour? 

If I increase the prior for the range, this effect is not as evident, although still present. 

The code used is the following:

response_var <- 'y'
spde <- inla.spde2.pcmatern(mesh = mesh,  alpha = 2, constr=TRUE, prior.range=c(0.01,0.1), prior.sigma=c(0.2,0.1))
s.index <- inla.spde.make.index(name = "spatial.field",n.spde = spde$n.spde)

#Create data structure
A.train <- inla.spde.make.A(mesh = mesh, loc = coordinates(data))
stack.train <- inla.stack(data  = list(response_var = data[[response_var]]),
                                A = list(A.train, 1),
                                effects = list(c(s.index, list(Intercept = 1)),
                                                data.frame(data) %>%                      
                                                as.list()),
                                tag = "train")

#Create data structure for prediction
A.pred <- inla.spde.make.A(mesh = mesh, loc = coordinates(data.grid))
stack.pred <- inla.stack(data  = list(response_var = NA),
                                A = list(A.pred, 1),
                                effects = list(c(s.index, list(Intercept = 1)),
                                                data.frame(data.grid) %>%                      
                                                as.list()),
                                tag = "pred")

#Join stack
stack.join <- inla.stack(stack.train, stack.pred)

#Fit model
ff <- as.formula(paste('response_var', paste(c('-1', 'Intercept', 'f(spatial.field, model = spde)'), collapse=" + "),sep='~'))
model <- inla(ff, data = inla.stack.data(stack.join, spde = spde),
    family = 'nbinomial',
    control.predictor = list(A = inla.stack.A(stack.join), compute = TRUE,link = 1),
    control.compute = list(cpo = TRUE, dic = TRUE), verbose = TRUE)

Thank you,

V.

[Finn Lindgren]

Hi,

this is due to the conditionally deterministic behaviour inside each triangle.
For visualisation, one can evaluate the field and uncertainty on each
mesh node, and only then interpolate to a grid; linear interpolation
of the standard deviation gives a kind of "upper bound".
This can be done with inla.mesh.project on the "sd" values at each mesh node.

The INLA::meshbuilder() "shiny" tool has mesh quality assessment tools
that one can use to visually assess whether the triangles are small
enough to resolve a field with a given correlation range. You picture
is what typically appears when the triangles are too large to capture
the small scale variation in the model.

Finn

> --
> You received this message because you are subscribed to the Google Groups "R-inla discussion group" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to r-inla-discussion...@googlegroups.com.
> To view this discussion on the web, visit https://groups.google.com/d/msgid/r-inla-discussion-group/8a938ec0-8f38-421b-bffe-03efc0b5c591n%40googlegroups.com.



--
Finn Lindgren
email: finn.l...@gmail.com