Hello,
Thanks for the answer for all of questions.
I faced some error, which I dont understand why it happens.
I am buildling a joint likelihood model using "Poisson" and "Binomial" distirbution.
Regarding Ntrials, I have put 1 values for data stack for Poisson and associated values for data stack for Binomial.
I have checked Ntrials has no NA value.
But, INLA says
"NA's in argument 'Ntrials' are not allowed"
Could you please see my code below?
Thanks a lot!
Best
Jun-Sik
----------------------------------------
### Data A for Poisson distribution
mydata = ssubset_db_pk_env_dem
num = mydata$num
st.group = mydata$mon
m.time = st.group %>% max
spde_spatial_index <- inla.spde.make.index("spatial",
n.spde = spde$n.spde,
n.group = m.time
)
coords = mydata[c("x", "y")] %>% as.matrix
A_inci <- inla.spde.make.A(mesh = mesh_new, loc = coords, group = st.group )
stack_A <- inla.stack(
data = list(y = cbind(num, NA),
, Ntrials = cbind(rep(1, length(num)),1),
A = list(A_inci, 1),
effects = list(
list(spatial_A = spde_spatial_index$spatial,
spatial_A.group = spde_spatial_index$spatial.group),
list(intercept_A = rep(1, length(num)))
),
tag = "A"
)
### Data B for binomial distribution
st.group = mydata2$mon
n.group = mydata2 %>% pull(mon) %>% unique %>%length
coords = mydata2[c("X", "Y")] %>% as.matrix
A_sero <- inla.spde.make.A(mesh = mesh_new,
loc = coords,
group = st.group)
num = mydata2 %>% pull(num)
tot_num = mydata2 %>% pull(tot_num)
stack_B <- inla.stack(
data = list(y = cbind(NA, num), Ntrials = cbind(NA, tot_num)),
A = list(A_sero, 1),
effects = list(
list(spatial_B = spde_spatial_index$spatial,
spatial_B.group = spde_spatial_index$spatial.group),
list(intercept_B = rep(1, length(num)))
),
tag = "B"
)
formula = y ~ -1 +
intercept_A +
intercept_B +
f(spatial_A, model=spde,
group=spatial_A.group,
control.group=list(model="ar1")) +
f(spatial_B, copy = "
spatial_A",
hyper = list(beta = list(fixed = FALSE)),
group = spatial_B.group )
stack_AB= inla.stack(stack_A, stack_B)
dat <- inla.stack.data(stack_AB)
model_joint <- inla(
formula,
data = dat,
Ntrials = dat$Ntrials
family = c("poisson", "binomial"),
control.predictor = list(A = inla.stack.A(stack_AB), compute = TRUE),
control.compute = list(
dic = TRUE,
waic = TRUE,
cpo = TRUE,
return.marginals.predictor = TRUE,
config = TRUE
),
control.inla = list(
lincomb.derived.correlation.matrix = TRUE,
strategy = "gaussian",
int.strategy = "eb"
),
control.fixed = list(
correlation.matrix = TRUE
)
)