Hi all
I am trying to fit a fairly large spatial-temporal model, but am having problems using the full dataset. The data have about 3000 locations and 365 time steps, and fits well, with no problems to about 60% of the data (all spatial locations but 60% of the time steps).
I found another answer here (https://groups.google.com/g/r-inla-discussion-group/c/uf2ZGh4jmWc) that suggested adding a small non-negative diagonal value in control.inla, and then using the output of that model as starting values in a second model with diagonal = 0. If I do this, then the first model will run on the full dataset, but I cannot get the second model to run. This crashes with the following error:
inla.mkl: smtp-pardiso.c:757: GMRFLib_pardiso_solve_core: Assertion `store->pstore->done_with_chol == GMRFLib_TRUE' failed.
Segmentation fault
Error in inla.inlaprogram.has.crashed() :
The inla-program exited with an error. Unless you interupted it yourself, please rerun with verbose=TRUE and check the output carefully.
If this does not help, please contact the developers at <he...@r-inla.org>.
Do you have any suggestions about how to fix this? Or where I may be going wrong? I’ve copied the code below.
Thank you!
Simon Brewer
## Priors
### BYM2
prec_u <- .5/.31
prec_alpha <- .01
phi_u <- .5
phi_alpha <- 2/3
prior_bym2 <- list(
prec = list(
prior = "pc.prec",
param = c(prec_u, prec_alpha)),
phi = list(
prior = "pc",
param = c(phi_u, phi_alpha))
)
### RW1 From INLA Germany example
rw1_u = .2/.31
rw1_alpha = .01
prior_rw1 <- list(theta = list(prior="pc.prec", param=c(rw1_u, rw1_alpha)))
## Model formula
f8 <- TESTSNEW7POS_CDC ~ scale(Pop_m) + scale(gini) + scale(Uninsured) + scale(Production) + RUCC + scale(PCP) + scale(RPL_THEMES) +
f(date1, model = "rw1", hyper = prior_rw1,
scale.model = TRUE) +
f(struct, model = "bym2", graph = g, hyper = prior_bym2,
scale.model = TRUE)
## Model w/ diagonal value
res1 <- inla(f8, data = dat2,
control.predictor = list(compute = TRUE),
control.compute = list(dic = TRUE, waic = TRUE),
control.fixed = list(prec.intercept = 1),
control.inla = list(strategy = "gaussian",
int.strategy = "eb",
diagonal = 5e-2),
verbose = TRUE, num.threads = 32)
## Model w/ diagonal value
res1 <- inla(f8, data = dat2,
control.predictor = list(compute = TRUE),
control.compute = list(dic = TRUE, waic = TRUE),
control.fixed = list(prec.intercept = 1),
control.inla = list(strategy = "gaussian",
int.strategy = "eb",
diagonal = 5e-2),
verbose = TRUE, num.threads = 32)
## Model w/ restart
res1 <- inla(f8, data = dat2,
control.predictor = list(compute = TRUE),
control.compute = list(dic = TRUE, waic = TRUE),
control.fixed = list(prec.intercept = 1),
control.mode = list(result = res1, restart = TRUE),
control.inla = list(diagonal = 0),
verbose = TRUE, num.threads = 32,
)