INLA Keeping on crashing

[Kalondepatrick]

I am learning the process of fitting a spatiotemporal model using INLA. In my pursuit, I simulated some data for Malawi, and every time I tried to fit my model with random effects, my model keep on crashing. I do not understand why this is the case, but the error that pops up indicates non-numeric argument to binary operator. I have attached a reproducible example on my Git profile (https://github.com/Kalondepatrick/spatiotemporal). To those interested to help, you will find the script that I was using in the folder scripts.

[finnlindgren]

You had some basic errors in your inla() call syntax (, instead of + after pop, incorrect idarea variable name, and extraneous Ntrials) that are fix by this:

res <- inla(cases ~pop +
            f(idate, model = 'iid') +
              f(idarea, model ='besag', graph = malawi_graph),
            family = "poisson", data = map_sf2, verbose = T)

That allows it to get further, but gets a different error:

idx[51] =  51  val = 1.09892e+11
        idx[52] =  52  val = 1.45924e+11
        idx[53] =  53  val = 3.61802e+1inla.mkl: idxval.c:841: GMRFLib_idxval_nsort_x_core: Assertion `0 == 1' failed.
Aborted (core dumped)

I don't know what that's about.

[hrue]

I think this is what you want.

res <- inla(cases ~ f(idate, model = 'iid') + f(idarea, model ='besag', graph = malawi_graph, scale.model = TRUE), E = pop, family = "poisson", data = map_sf2)

[Kalondepatrick]

Thanks all. Specifically @hrue the approach that you suggested worked for me. Just for the sake of my own learning, why do we have to specify our explanatory variable as E = x. How can the same be done if we have multiple explanatory variables i.e. population, average surface temperature, average rainfall, etc? In case there is some documentation that you can recommend to me for my own readings regarding this, please share them. Thanks

[hrue]

I guess population is an offset and not a covariate.

since mean=exp(linear.predictor), you need to give this in logscale, so

.. ~ offset(log(pop)) + ...

and

.. ~ ...

E=pop,

is the same

see the documentation

inla.doc("poisson")

On Mon, 2023-02-06 at 22:15 -0800, Patrick Ken Kalonde wrote:

Thanks all. Specifically @hrue the approach that you suggested worked for me. Just for the sake of my own learning, why do we have to specify our explanatory variable as E = x. How can the same be done if we have multiple explanatory variables i.e. population, average surface temperature, average rainfall, etc? In case there is some documentation that you can recommend to me for my own readings regarding this, please share them. Thanks — Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were mentioned.Message ID: @.***>

-- Håvard Rue Professor of Statistics Chair of the Statistics Program CEMSE Division King Abdullah University of Science and Technology Thuwal 23955-6900 Kingdom of Saudi Arabia

@.*** Office: +966 (0)12 808 0640   Mobile: +966 (0)54 470 0421 Research group: bayescomp.kaust.edu.sa   R-INLA project: www.r-inla.org Zoom: kaust.zoom.us/my/haavard.rue

--

This message and its contents, including attachments are intended solely for the original recipient. If you are not the intended recipient or have received this message in error, please notify me immediately and delete this message from your computer system. Any unauthorized use or distribution is prohibited. Please consider the environment before printing this email.

[Kalondepatrick]

@hrue thanks so much for the guidance. I have two questions related to the inquiry above.

1. Given that population has been embedded in the model as an offset, how can one incorporate a covariate, i.e., adding elevation in the model. I tried to add it as cases ~pop, but doing so INLA keeps on crashing.

2. Another one is on the possibility of using the spatiotemporal model to make predictions. I am interested in using it for predicting cases, and I am unable to find resources that has presented such instructions

[hrue]

y ~ offset(log(pop)) + elevation + ...

for prediction, see the SPDE book https://www.r-inla.org/learnmore/books

[Kalondepatrick]

Thanks so much, @hrue . I went through some of the books and noticed that in most cases where prediction is made, it focused much on Relative Risk. This makes me to wonder whether this can be expanded to the prediction of actual counts

[hrue]

if you want to predict the acutal counts you have to do that manually, either by sampling using inla.posterior.sample() and inla.posterior.sample.eval(), or by numerical integration. All depending on what your aim is.

On Fri, 2023-02-24 at 01:59 -0800, Patrick Ken Kalonde wrote:

Thanks so much, @hrue . I went through some of the books and noticed that in most cases where prediction is made, it focused much on Relative Risk. This makes me to wonder whether this can be expanded to the prediction of actual counts — Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were mentioned.Message ID: @.***>

-- Håvard Rue Professor of Statistics Chair of the Statistics Program CEMSE Division King Abdullah University of Science and Technology Thuwal 23955-6900 Kingdom of Saudi Arabia

@.*** Office: +966 (0)12 808 0640   Mobile: +966 (0)54 470 0421 Research group: bayescomp.kaust.edu.sa   R-INLA project: www.r-inla.org Zoom: kaust.zoom.us/my/haavard.rue

--

This message and its contents, including attachments are intended solely for the original recipient. If you are not the intended recipient or have received this message in error, please notify me immediately and delete this message from your computer system. Any unauthorized use or distribution is prohibited. Please consider the environment before printing this email.

[Kalondepatrick]

@hrue are there any tutorials or examples on the manual prediction of the counts

[hrue]

library(INLA)

data(Germany) g = system.file("demodata/germany.graph", package="INLA") source(system.file("demodata/Bym-map.R", package="INLA")) summary(Germany)

just make a duplicated column

Germany = cbind(Germany,region.struct=Germany$region)

Germany <- Germany[1:10, ]

g <- inla.read.graph(matrix(1, 10, 10))

form = Y ~ f(region.struct,model="besag",graph=g) + f(region,model="iid") + f(x, model="rw2")

r = inla(form, family="poisson", data=Germany, E=E, verbose = TRUE, control.predictor = list(compute = TRUE), control.compute = list(config = TRUE))

xx <- inla.posterior.sample(10000, r) sim.y <- function(E = E) { rpois(length(E), E * exp(Predictor)) } yy <- inla.posterior.sample.eval(sim.y, xx, E = Germany$E)

r$summary.fitted.values$mean[1] * Germany$E[1] mean(yy[1,])

hist(yy[1, ], prob = TRUE)

On Fri, 2023-02-24 at 04:10 -0800, Patrick Ken Kalonde wrote:

@hrue are there any tutorials or examples on the manual prediction of the counts — Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you were mentioned.Message ID: @.***>

-- Håvard Rue Professor of Statistics Chair of the Statistics Program CEMSE Division King Abdullah University of Science and Technology Thuwal 23955-6900 Kingdom of Saudi Arabia

@.*** Office: +966 (0)12 808 0640   Mobile: +966 (0)54 470 0421 Research group: bayescomp.kaust.edu.sa   R-INLA project: www.r-inla.org Zoom: kaust.zoom.us/my/haavard.rue

--

This message and its contents, including attachments are intended solely for the original recipient. If you are not the intended recipient or have received this message in error, please notify me immediately and delete this message from your computer system. Any unauthorized use or distribution is prohibited. Please consider the environment before printing this email.