Smoothness with higher order neighborhood matrix

69 views
Skip to first unread message

Hanah

unread,
Jun 3, 2022, 11:56:10 AM6/3/22
to R-inla discussion group
Hi everybody,

I'm using a besag model to estimate the relative risks. I tried fitting with different neighborhood matrices, with the adjacency order runs up to the 1st, 3rd and 9th order (I'm expanding the neighborhood structure).  An example of the matrix with 3rd order is given below. 
  Capture.PNG

The matrix is added into the model using the "graph" argument as:
fit <- inla(dth ~  1+f(ID, model = "besag",graph=Wmat),
            data = as.data.frame(data),
            family = "poisson",
            E = data$E2020, control.predictor = list(compute = TRUE),
            control.compute = list(dic = TRUE,waic=TRUE,config=T),
            control.inla=list(strategy="gaussian"))

I expected that when expanding the neighborhood matrix, I would get smoother trend over the study region, but I somehow saw the opposite. I also tried adjusting the neighborhood matrix by either setting all the values >1 equal to 1, or  equal to (1/order number) but they give similar estimates. 

Could you kindly advise if I made mistakes and how to correctly put high-order neighborhood matrix into the model?

Thanks a lot in advance!

Best regards,
Hanah.

Helpdesk

unread,
Jun 3, 2022, 12:02:12 PM6/3/22
to Hanah, R-inla discussion group
the answer is in the SPDE paper 2011 (JRSS-B), and also in the
discussion contribution there from John Kent.

On Fri, 2022-06-03 at 08:56 -0700, Hanah wrote:
> Hi everybody,
>
> I'm using a besag model to estimate the relative risks. I tried
> fitting with different neighborhood matrices, with the adjacency order
> runs up to the 1st, 3rd and 9th order (I'm expanding the neighborhood
> structure).  An example of the matrix with 3rd order is given below. 
>   Capture.PNG
>
> The matrix is added into the model using the "graph" argument as:
> fit <- inla(dth ~  1+f(ID, model = "besag",graph=Wmat),
>             data = as.data.frame(data),
>             family = "poisson",
>             E = data$E2020, control.predictor = list(compute = TRUE),
>             control.compute = list(dic = TRUE,waic=TRUE,config=T),
>             control.inla=list(strategy="gaussian"))
>
> I expected that when expanding the neighborhood matrix, I would get
> smoother trend over the study region, but I somehow saw the opposite.
> I also tried adjusting the neighborhood matrix by either setting all
> the values >1 equal to 1, or  equal to (1/order number) but they give
> similar estimates. 
>
> Could you kindly advise if I made mistakes and how to correctly put
> high-order neighborhood matrix into the model?
>
> Thanks a lot in advance!
>
> Best regards,
> Hanah.
> --
> 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/75e49585-b0e5-41a4-90ee-894edf1a1600n%40googlegroups.com
> .

--
Håvard Rue
he...@r-inla.org

Hanah M

unread,
Jun 3, 2022, 12:46:18 PM6/3/22
to Helpdesk, R-inla discussion group
Thank you so much for the prompt response!

Best regards,
Hanah.

On Fri, 3 Jun 2022 at 18:07, Helpdesk <he...@r-inla.org> wrote:
like

>  g = inla.read.graph(system.file("demodata/germany.graph",
package="INLA"))
> Q = INLA:::inla.pc.bym.Q(g)
> Q2=Q%*%Q
> g2 = inla.read.graph(Q2)

you can input Q2 like into the generic model, but I guess you want to
scale it first to control the precision. for this you can use

> constr=list(A=matrix(1,1,n),e=0)
> Q2.scaled = INLA:::inla.scale.model(Q2, constr)
Reply all
Reply to author
Forward
0 new messages