inla.stack turning effects into NAs

[Mark Myer]

Hello, 

I'm having a basic issue in conducting a discrete space-time separable INLA SPDE as outlined in the SPDE Tutorial, Section 5.1. My example data is structured as such, with 5474 observations of a binomial response over the course of 23 weeks, along with the coordinates of each observation and the precipitation that week:

> str(example.df)

'data.frame': 5474 obs. of  5 variables:

 $ dat$WNVLabResult: int  0 0 0 0 1 0 0 0 1 0 ...

 $ dat$Prcp        : num  8.14 9.29 0 8.14 1 ...

 $ dat$MLatX      : num  4519163 4519163 4519163 4519163 4519163 ...

 $ dat$MLongY     : num  606832 606832 606832 606832 606832 ...

 $ dat$Week        : num  3 4 5 6 8 9 10 13 14 14 ...

I'm creating the mesh, spde, projector matrix, and stack using the following code and everything seems to work fine, until I try to check the contents of the stack under 'effects' and find that my precipitation variable is being read as 'NA' and is not the same length as my number of observations. 

#Create the mesh for INLA observations

loc <- cbind(dat$MLongY/1000,dat$MLatX/1000) #Divide by 1000 to get km

bnd<-inla.nonconvex.hull(loc) #Makes a nonconvex hull around the points rather than the spatial extent of the island

mesh<-inla.mesh.2d(boundary=bnd, cutoff=1, max.edge=c(2,3))

#Create the SPDE model using the pc-prior 

spde <- inla.spde2.pcmatern(

  mesh=mesh, alpha=2, ### mesh and smoothness parameter

  prior.range=c(10, 0.01), ### P(practic.range<10 km)=0.01

  prior.sigma=c(1, 0.01) ### P(sigma>1 km)=0.01

  ) 

#Create the index set for specifying number of time steps in spacetime model

iset <- inla.spde.make.index('i', n.spde=spde$n.spde, n.group=max(dat$Week))

#Create the projector matrix based on the locations and time of the observations

A <- inla.spde.make.A(mesh=mesh,

                      loc=loc,

                      group=dat$Week)

#Create the data stack. For this example we're just going to use precipitation as the only predictor

sdat <- inla.stack(tag='stdata', data=list(y=dat$WNVLabResult),

                   A=list(A,1), effects=list(iset, list(Prcp=dat$Prcp)))

my effects list is as follows, and I can't figure out why it has 2529 observations rather than 5474, and why the precipitation variable, Prcp, is being turned into 'NA':

> str(sdat$effects)

List of 5

 $ data :'data.frame': 2529 obs. of  5 variables:

  ..$ i        : int [1:2529] 236 367 613 108 113 136 173 186 197 289 ...

  ..$ i.group  : int [1:2529] 1 1 1 2 2 2 2 2 2 2 ...

  ..$ i.repl   : int [1:2529] 1 1 1 1 1 1 1 1 1 1 ...

  ..$ Intercept: num [1:2529] 1 1 1 1 1 1 1 1 1 1 ...

  ..$ Prcp     : num [1:2529] NA NA NA NA NA NA NA NA NA NA ...

 $ nrow : int 2529

 $ ncol : Named int [1:5] 1 1 1 1 1

  ..- attr(*, "names")= chr [1:5] "i" "i.group" "i.repl" "Intercept" ...

 $ names:List of 5

  ..$ i        : chr "i"

  ..$ i.group  : chr "i.group"

  ..$ i.repl   : chr "i.repl"

  ..$ Intercept: chr "Intercept"

  ..$ Prcp     : chr "Prcp"

 $ index:List of 1

  ..$ stdata: int [1:31763] NA NA NA NA NA NA NA NA NA NA ...

 - attr(*, "class")= chr "inla.data.stack.info"

> 

If any further information is needed to diagnose this issue, please let me know. 

Appreciate the help,

Mark H. Myer

myer...@gmail.com

[Finn Lindgren]

Hi,

The purpose of inla.stack is to automate the construction of input vectors with NA in the required places, so this is most likely not an error, but rather the intended behaviour.

Please see the "read this first" spde tutorial for an explanation of precisely what inla.stack is meant to do.

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 post to this group, send email to r-inla-disc...@googlegroups.com.
Visit this group at https://groups.google.com/group/r-inla-discussion-group.
For more options, visit https://groups.google.com/d/optout.

[Mark Myer]

Finn, 

I see, I thought it would be an issue when looking at the contents of the stack, but it seems that the model is finding the predictor values OK when I call them later on in the inla(). Thank you!

---Mark

To unsubscribe from this group and stop receiving emails from it, send an email to r-inla-discussion-group+unsub...@googlegroups.com.