I am trying to fit a set of parameters (\theta_1 and \theta_2) for subjects (s
) over different time points (t
) with a Gaussian dynamic process. My generative model looks like this:
where initial state follows:
y
is choice data for subject s in time t that is being fitted.
U
is also observed data for subject s in time t
A
, B
, C
, X
, u
, R
and V
are dynamic parameters of the Gaussian process.
I implement this process in Stan in the following way:
model {
.
.
.
for (s in 1:S) {
for (t in 1:T) {
vector[2] params_pr;
params_pr[1] = theta1[s,t];
params_pr[2] = theta2[s,t];
if (w == 1) { // initial state
to_vector(X[s,t,]) ~ normal(0,sigma_v);
}
else {
to_vector(X[s,t,]) ~ normal(diag_matrix(A) * to_vector(X[s,t-1,]) + B * to_vector(U[s, t-1,]), Q);
}
params_pr ~ normal(mu_d + C * to_vector(X[s,t,]),sigma_r);
y[s,t,] ~ bernoulli_logit(theta1[s,t] * dat[s,t,]) ./ (1 + theta2[s,t]));
}
}
}
My general chain diagnostics (bin/diagnose from cmdstan) is not bad - no hitting max tree depth, almost no divergent transitions (0 to 3 transitions), satisfactory E-BFMI and ESS.
This is my typical bin/diagnose
output:
Checking sampler transitions treedepth.
Treedepth satisfactory for all transitions.
Checking sampler transitions for divergences.
No divergent transitions found.
Checking E-BFMI - sampler transitions HMC potential energy.
E-BFMI satisfactory for all transitions.
Effective sample size satisfactory.
Parameters \theta_1 and \theta_2 are fitted just fine - stable chain traces (see examples below; upper two), R-hat is 1.00 and also my predictive checks look amazing.
Now my problem is that A
, B
, C
and X
on the other hand are a complete disaster, as you can see in the example trace plots (bottom two), which is also reflected in huge R-hats.
I tried changing the number of draws, adapt_delta, step_size, put stricter priors, implementing a non-centered reparametrization. Nothing helped.
Can this be an identifiability issue, especially given the fact that I am able to recover the original behavior data? Or maybe I am missing something simple in how I specify my model in Satn? Any suggestions are welcome!
My system:
python 3.7
cmdstanpy 0.9.63
cmdstan 2.25.0
CentOS7
gcc 9.2.0