Hi,

I am trying to fit the following model in Stan. I have some idea, but I think it is incomplete. I will try to define the model below and then write the incomplete syntax I have now. I appreciate any input and suggestion. These models are described in papers before for reference.

https://doi.org/10.3102/1076998614559412

https://doi.org/10.3102/10769986031002181

The logarithm of the response times, Y_{ij} is the response time for the ith person on the jth item, are assumed to be normally distributed with a mean of \phi_j*(\lambda_j - \zeta_i) and variance of \sigma_j^2.

\phi_j : item specific time discrimination parameter

\lambda_j : item specific time intensity parameter

\zeta_i : individual specific latent working speed parameter

\sigma_j^2 : item specific measurement error variance

For model identification, suppose we fix the mean and variance of working speed parameter to 0 and 1. So, I assume \zeta_i ~ N(0,1).

In these papers, they also assume that the item specific parameters, \phi and \lambda, are correlated, and a bivariate normal distribution is assumed to describe the relationship.

Below is my effort to fit this model, but I know it is incomplete and probably wrong. I am struggling to translate this model to Stan language. Especially, the covariance matrix among parameters and how to specify priors for them (normal inverse wishart as suggested in one of these papers).

```
rt <- '
data {
int <lower=0> I; // number of individuals
int <lower=0> J; // number of items
int <lower=0,upper=1> RT[I,J]; // matrix of response times
}
parameters {
vector[J] lambda; //item specific time intensity parameter
vector <lower=0> [J] phi; //item speficic time discrimination parameter
vector <lower=0> [J] sigma; // item specific residual standard deviation,
vector[I] zeta; // person specific latent working speed
vector[2] mu; // mean vector for lambda and phi
matrix[2,2] cov; // covariance matrix between lambda and phi
}
model {
zeta ~ normal(0,1);
//mu and cov priors??????
phi ~ normal(mu[1],cov[1,1]);
lambda ~ normal(mu[2],cov[2,2]);
sigma ~ inv_gamma(5,5); // would this be appropriate??
for (i in 1:I){
for (j in 1:J){
real p = phi[j]*(zeta[i]-lambda[j]);
RT[i,j] ~ normal(p,sigma[j]);
}
}
}
'
```