Dear Stan community,
d.Rdata (10.2 KB)
Intercept_model.stan (1.9 KB)
model_fitting.rmd (1.4 KB)
I am fitting an intercept model to my data, but I encountered confusing fitting problems, and I would like to ask for help.
In my data, the outcome/dependent variables are participants’ skin conductance responses (SCR). Please see attached files, where I attached the Stan model, an example dataset with 10 participants, and a RMarkdown file in which using the rstan to do the model fitting.
There are three things that should be noted:
(1) Due to the experimental protocol, half of the SCR are missing (i.e., NaN), and I don’t want to delete them, because in my experiment the trial sequences are meaningful.
(2) Since SCR (i.e., the DV) can only take positive values, so in my model, I employed a truncated Gaussian likelihood.
(3) My model includes varying/random intercepts at the participant level.
Following is my Stan model, which is syntactically correct:
data {
int<lower=0> N_obs; // number of observed data points
int<lower=0> N_mis; // number of missing data points
int<lower=1> N_all; // number of total data points (i.e., observed + missing)
int<lower=1, upper=N_all> ii_obs[N_obs]; // indexes of the observed data
int<lower=1, upper=N_all> ii_mis[N_mis]; // indexes of the missing data
int<lower=1> S; // number of participants
int<lower=1, upper=S> ID[N_all]; // participant index
array[N_obs] real<lower=0> SCR_obs; // observed SCR data
}
parameters {
// missing data
array[N_mis] real<lower=0> SCR_mis; // missing SCR data
// Hierarchical mu parameters (at the population-level)
real mu_b0; // intercept
// sigma parameters
real<lower=0> sigma; // for the truncated Gaussian likelihood
real<lower=0> sigma_b0; // from the hierarchical distribution of intercept
// Normal distributions for non-center parameterization, i.e., Matt trick
real z_b0[S];
}
transformed parameters {
// data
array[N_all] real<lower=0> SCR;
SCR[ii_obs] = SCR_obs;
SCR[ii_mis] = SCR_mis;
// participant-level parameters
real b0[S];
for (s in 1:S) {
b0[s] = mu_b0 + z_b0[s]*sigma_b0;
}
}
model {
mu_b0 ~ normal(0, 10);
sigma ~ gamma(2, 2);
sigma_b0 ~ gamma(2, 2);
z_b0 ~ normal(0, 1);
for (n in 1:N_all) {
target += normal_lpdf(SCR[n] | b0[ID[n]], sigma) -
normal_lccdf(0 | b0[ID[n]], sigma);
}
}
generated quantities {
vector[N_all] log_lik;
{for (n in 1:N_all) {
log_lik[n] = normal_lpdf(SCR[n] | b0[ID[n]], sigma) -
normal_lccdf(0 | b0[ID[n]], sigma);
}
}
}
When I fitted this model, I first tried one chain with 5,00 warmup and 1,000 iteration, and adapt_delta = 0.95. But the model fitting failed and I got the following warnings:
I have no idea what’s going wrong. Might anyone know what’s going wrong? I am immeasurably grateful for any suggestions/recommendations/comments.
I hope I included all relevant information and this is the right place to ask a question like this.
Thank you very much!!!
Operating System: Windows 10
rstan Version: 2.26.13
Stan Version: 2.26.1
All the best,