EDIT This issue appears to have resulted from rstan not having up-to-date stan syntax because of the restrictions that CRAN places on package updates. A more up-to-date rstan package appears to be available at the mc-stan repo, but does not yet have binaries for the most recent R version (4.3).
I am not sure exactly where this question goes, but haven't found an answer anywhere else, so here goes.
I need to learn stan to run some CJS mark-recapture models, and have been working from the stan user guide here to do so. The collective Cormack-Jolly-Seber model works fine, and is fairly straightforward in terms of the data going into it, but I really need to get a version of an individual CJS working. The problem is that the guide does a very poor job of describing the data that goes into the model and there is no data provided. As a result, when I run the example model with real data, it keeps throwing back errors that suggest something is wrong with the input data.
Here is the data block, which requires a capture history matrix y of I rows and T columns:
data {
int<lower=2> T;
int<lower=0> I;
int<lower=0,upper=1> y[I, T];
}
From the parameters, I know that it also requires phi and p vectors of lengths T-1 and T, respectively.
parameters {
vector<lower=0,upper=1>[T-1] phi;
vector<lower=0,upper=1>[T] p;
}
So far this is consistent with the collective CJS model and shouldn't be causing issues. I have also included values for T and I because of their inclusion in the data block. So so far my input data looks like this (I am using rstan, so this is a toy version of what it looks like in r):
history.array<-list(y=matrix(c(1,0,0,1,0,0,0,0,1,0,0,1,
0,0,0,0,1,1,0,1,1,0,0,1,
1,1,0,0,1,1,0,0,0,0,0,0,
1,1,0,0,1,0,0,1), nrow=11),
phi=rep(.5,4-1),
p=rep(.5,4),
T = as.integer(4),
I = as.integer(11))
And here is the full of the stan code from the guide (which I have saved as "individual_CJS.stan").
data {
int<lower=2> T;
int<lower=0> I;
int<lower=0, upper=1> y[I,T];
}
functions {
int first_capture(int[] y_i) {
for (k in 1:size(y_i))
if (y_i[k])
return k;
return 0;
}
int last_capture(int[] y_i) {
for (k_rev in 0:(size(y_i) - 1)) {
int k;
k = size(y_i) - k_rev;
if (y_i[k])
return k;
}
return 0;
}
vector prob_uncaptured(int T, vector p, vector phi) {
vector[T] chi;
chi[T] = 1.0;
for (t in 1:(T - 1)) {
int t_curr;
int t_next;
t_curr = T - t;
t_next = t_curr + 1;
chi[t_curr] = (1 - phi[t_curr])
+ phi[t_curr]
* (1 - p[t_next])
* chi[t_next];
}
return chi;
}
}
transformed data {
int<lower=0,upper=T> first[I];
int<lower=0,upper=T> last[I];
vector<lower=0,upper=I>[T] n_captured;
for (i in 1:I)
first[i] = first_capture(y[i]);
for (i in 1:I)
last[i] = last_capture(y[i]);
n_captured = rep_vector(0, T);
for (t in 1:T)
for (i in 1:I)
if (y[i, t])
n_captured[t] = n_captured[t] + 1;
}
parameters {
vector<lower=0,upper=1>[T-1] phi;
vector<lower=0,upper=1>[T] p;
}
transformed parameters {
vector<lower=0,upper=1>[T] chi;
chi = prob_uncaptured(T,p,phi);
}
model {
for (i in 1:I) {
if (first[i] > 0) {
for (t in (first[i]+1):last[i]) {
1 ~ bernoulli(phi[t - 1]);
y[i, t] ~ bernoulli(p[t]);
}
1 ~ bernoulli(chi[last[i]]);
}
}
}
generated quantities {
real beta;
vector<lower=0>[T] pop;
beta = phi[T - 1] * p[T];
pop = n_captured ./ p;
pop[1] = -1;
}
But when I try to run the program with these inputs:
warmups <- 1000
total_iterations <- 2000
max_treedepth <- 10
n_chains <- 4
n_cores <- 4
adapt_delta <- 0.95
test_ind_CJS<-stan(
file = "individual_CJS.stan",
chains = n_chains,
warmup = warmups,
data = history.array,
iter = total_iterations,
cores = n_cores,
refresh = 250,
control = list(max_treedepth = max_treedepth,
adapt_delta = adapt_delta)
),
I get this error:
PARSER FAILED TO PARSE INPUT COMPLETELY
STOPPED AT LINE 10:
functions {
int first_capture(int[] y_i) {
for (k in 1:size(y_i))
if (y_i[k])
return k;
return 0;
}.......*the rest of the stan code is printed here*.......
Error in stanc(file = file, model_code = model_code, model_name = model_name, :
failed to parse Stan model 'individual_CJS' due to the above error.
So that first_capture function seems to need something, but I can't figure out what. y_i, shouldn't need to be input if it is just the ith row of matrix y. What am I missing here? thanks in advance