-
-
Notifications
You must be signed in to change notification settings - Fork 477
/
Copy pathice.stan
55 lines (50 loc) · 1.23 KB
/
ice.stan
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
// http://www.mrc-bsu.cam.ac.uk/bugs/winbugs/Vol2.pdf
// Page 29: Ice: non-parametric smoothing in an age-cohort model
// The model is the same specified here as in JAGS example,
// but this example's JAGS version is different from that of
// WinBUGS because JAGS does not support some feature.
// TODO:
// Maybe we should check if we could specify
// the model the same as WinBUGS?
// status: the results are farely cose to those
// from JAGS
data {
int<lower=0> N;
int<lower=0> Nage;
int<lower=0> K;
array[N] int year;
array[N] int cases;
array[N] int age;
array[N] int pyr;
real alpha1;
}
parameters {
array[Nage - 1] real alpha;
array[K] real beta;
real<lower=0, upper=1> sigma;
}
model {
vector[N] r;
sigma ~ uniform(0, 1);
for (k in 1 : 2) {
beta[k] ~ normal(0, sigma * 1E3);
}
for (k in 3 : K) {
beta[k] ~ normal(2 * beta[k - 1] - beta[k - 2], sigma);
}
alpha ~ normal(0, 1000);
for (i in 1 : N) {
if (age[i] == 1) {
r[i] = alpha1 + log(pyr[i]) + beta[year[i]];
} else {
r[i] = alpha[age[i] - 1] + log(pyr[i]) + beta[year[i]];
}
}
cases ~ poisson_log(r);
}
generated quantities {
array[K] real logRR;
for (k in 1 : K) {
logRR[k] = beta[k] - beta[5];
}
}