Hi @bbbales2 - sorry to leave the declarations out of the thread… The structure of the model at each level is similar, where I am trying to define a vector of residuals at each level, and in my current model, i have also defined the offsets as vectors…
parameters {
vector[J] uj;
vector <lower=0> [P] sigma_item;
vector <lower=0> [M] fam_resid;
vector <lower=0> [Q] mod_resid;
vector[K] betai_offset;
vector[P] fammu_offset;
vector[M] modmu_offset;
vector[Q] template_mu;
real disp_horiz;
real disp_verb;
real char_m10;
real char_d3;
}
At each level of the model, there is a regression, and you’ll see here I have used the .* operator to multiply the offsets by the vector of residuals - and this is where I think my model is having difficulty.
transformed parameters{
vector[N] eta;
vector[K] betai;
vector[P] family_mu;
vector[M] model_mu;
// varying item model difficulty across item models within templates
// Item model difficulty is decomposed into a linear combination of the template difficulty
// and the effect(s) of display characteristics
// variances differ across item models
model_mu = template_mu[tt] + modmu_offset .* mod_resid[tt];
model_mu[vv] = model_mu[vv] + disp_verb;
model_mu[hh] = model_mu[hh] + disp_horiz;
// varying item family diffculty across item families within models
// decomposition of family means into a model mean and fixed effects
family_mu = model_mu[mm] + fammu_offset .* fam_resid[mm];
family_mu[multten] = family_mu[multten] + char_m10;
family_mu[threedig] = family_mu[threedig] + char_d3;
// item difficulties parameterized as random, with parent-specific means
betai = family_mu[parent] + betai_offset .* sigma_item[parent];
//log odds of a correct probability
eta = uj[peeps]-betai[item];
}
The model statement is here - and you will see that the priors for the random effects are now normal rather than cauchy, per [another thread/discussion](Cross-classified Hierarchical Model, Heterogenous Variances - #18 by bbbales2 - Thank you)!
model {
//hyperprior
betai_offset ~ normal(0,1);
fammu_offset ~ normal(0,1);
modmu_offset ~ normal(0,1);
//prior on random variable theta to scale the item difficulty parameters
uj ~ normal(0,1);
//weakly informative prior on item variance
sigma_item ~ normal(0,1);
fam_resid ~ normal(0,1);
mod_resid ~ normal(0,1);
//likelihood function
score ~ bernoulli_logit(eta);
}
I have done due diligence, I think to attempt to explore different priors, chain lengths, starting values, and estimation parameters (i.e. adapt_delta, etc.) in an effort to recover the variance parameters - so I think there is an issue with how the residuals and/or offsetts (often called “tau”) are being coded.
As I am writing this - I’m wondering if the issue may be with how I have defined/indexed the residuals… I feel as though as I have been attempting to translate my models into code, I have gotten a bit tangled, even doing the step-through process of building from simple to more complex that folks have advocated on this forum!
Thanks again, and in advance!