Syntax for Lognormal Response Time Model - how to specify priors for parameters with multivariate normal distribution



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.

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.

\left ( \begin{matrix} \phi_j\\\lambda_j \end{matrix} \right ) \sim N\left (\left [ \begin{matrix} \mu_\phi\\ \mu_\lambda \end{matrix} \right ],\left [\begin{matrix} \sigma_\phi^2 & \sigma_{\phi,\lambda} \\ \sigma_{\lambda,\phi} & \sigma_\lambda^2 \end{matrix} \right ] \right )

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]);


Check out the manual chapter on regression. It explains our preferred approach to multivariate priors, which doesn’t involve Wisharts, but rather independent priors on the correlations and on the scales.

You need constrained types for things like covariance. This is also explained in the manual. But why define a covariance matrix if you’re just going to ignore it and sample phi and lambda from teh diagonal?

We don’t recommend inverse gammas, but (5, 5) is probably OK.

You can vectorize that RT sampling statement. You probably need to log it first.