fbachl
commented
6 years ago There are actually three ways of doing this: fit3.spde <- lgcp(mdl3, pts3, ips = ipoints(c(0,55), 200, name = "x")) fit3a.spde <- lgcp(mdl3, pts3, domain = list(x = c(0,55))) fit3b.spde <- lgcp(mdl3, pts3, domain = list(x = mesh1D))
the difference is how the integration points are constructed: 1) Manually set the ips to wherever you want them to be. Useful for debugging but you need to know what you are doing. 2) This only states the boundaries (0 and 55) of the area to integrate over. By default 100 (i think) equidistant integration points are set. Also not really ideal. 3) Use the mesh to set the integration points exactly at the vertex locations. Numerically this is the most stable approach and best practice... But harder to explain to newbies.
The problem with 1) and 2) is what we observed as a big issue with the very first integration schemes we used during the project. Sightings pull the intensity up, integration points pull it down. If there are too few integration points and a sighting you end up with infinite intensity. To be precise: The vertices that surround a sighting should always have an integration point attached to them.
The error message you got is a ggplot problem. One of the numbers was super large and apparently ggplot doesn't like that
Summary of problem:
Sometimes “domain=“ does not work properly. For example, first chunk of code below (using “ips=“) works fine, but second chunk (using “domain=“) does not. Symptom: very different posteriors, apparently no variance on matern covariance function, and error message (shown at bottom) when trying to plot (because of Infs in prediction, I think).
Get data and set up model:
library(inlabru) data(Poisson3_1D) cd3 <- countdata3a x <- seq(0, 55, length = 50) mesh1D <- inla.mesh.1d(x, boundary = "free") xs <- seq(0, 55, length = 150) pr.df <- data.frame(x=xs)
mdl3 <- x ~ spde1D(map = x, model = inla.spde2.pcmatern(mesh1D, prior.range = c(150,0.75), prior.sigma = c(0.1,0.75))) + Intercept
Fit Using ips=
fit3.spde <- lgcp(mdl3, pts3, ips = ipoints(c(0,55), 200, name = "x")) post.range = spde.posterior(fit3.spde,name="spde1D",what="range");plot(post.range) post.variance = spde.posterior(fit3.spde,name="spde1D",what="variance");plot(post.variance) post.matern = spde.posterior(fit3.spde,name="spde1D",what="matern.covariance");plot(post.matern) pred3_spde <- predict(fit3.spde, pr.df, ~ exp(spde1D + Intercept)) plot(pred3_spde)
Fit Using domain=
fit3a.spde <- lgcp(mdl3, pts3, domain = list(x = c(0,55))) post.rangea = spde.posterior(fit3a.spde,name="spde1D",what="range");plot(post.range) post.variancea = spde.posterior(fit3a.spde,name="spde1D",what="variance");plot(post.variance) post.materna = spde.posterior(fit3a.spde,name="spde1D",what="matern.covariance");plot(post.matern) pred3a_spde <- predict(fit3a.spde, pr.df, ~ exp(spde1D + Intercept)) plot(pred3a_spde)
Error message from last command above:
Error in if (score > best$score && (!only.loose || (lmin <= dmin && lmax >= : missing value where TRUE/FALSE needed