For example, the INLA code for the "seeds" example (my of my favorites) is like this:
library(INLA)
data(Seeds)
formula = r ~ x1*x2+f(plate,model="iid")
mod.seeds = inla(formula,data=Seeds,family="binomial",Ntrials=n)
## improved estimation of the hyperparameters
h.seeds = inla.hyperpar(mod.seeds)
while the Stan code is like this:
data {
int
int
int
real x1[I];
real x2[I];
}
parameters {
real alpha0;
real alpha1;
real alpha12;
real alpha2;
real
real b[I];
}
transformed parameters {
real
sigma <- 1.0 / sqrt(tau);
}
model {
alpha0 ~ normal(0.0,1.0E3);
alpha1 ~ normal(0.0,1.0E3);
alpha2 ~ normal(0.0,1.0E3);
alpha12 ~ normal(0.0,1.0E3);
tau ~ gamma(1.0E-3,1.0E-3);
b ~ normal(0.0, sigma);
for (i in 1:I) {
n[i] ~ binomial(N[i], inv_logit(alpha0
+ alpha1 * x1[i]
+ alpha2 * x2[i]
+ alpha12 * x1[i] * x2[i]
+ b[i]) );
}
}
INLA has a more R-like syntax and is much more terse, but Stan is much more flexible and can handle very complicated models that INLA may have troubles with.
No comments:
Post a Comment