Hi there,
I obtained a predictor 3 times per individual. Then I obtained an outcome once for each individual (N = 120).
I want to estimate the “true” value for the predictor for each individual using partial pooling. Then I want to use this estimate to predict the outcome. So the model looks like this:
Let the predictor be denoted by X and the outcome by y.
X_{obs} \sim normal(\mu_{Xtrue}, \sigma_{X})
\mu_{Xtrue} = \alpha_{subject}
\alpha_{subject} \sim normal(\bar{a}, \sigma_{subject})
\bar{a} \sim normal(0, 1)
\sigma_{X} \sim exponential(1)
\sigma_{subject} \sim exponential(1)
y \sim normal(\mu, \sigma)
\mu \sim \alpha_y + X_{true} \times \beta
\alpha_y \sim normal(0, 1)
\beta \sim normal(0, 0.5)
\sigma \sim exponential(1)
Obviously the vectors for the first model have 360 values, because for each individual there are three values. However, my X_{true} estimate should be the pooled subject mean distributions (1 per subject, thus 120). The outcome vector has 120 values. My question is, how can I model this in stan. If I write the code the way I thought it should be I get an error saying that there cannot be different vector lengths. It only works if I put for each subject the same outcome value 3 times. But then, the model assumes that the outcome has been measured three times, which was not the case.
How do I have to specify this model. Here was the attempt where I put the outcome values 3x per subject but that is not the way I need it. I must code it somehow with vector y as 120 and that for the final model 120 “true” estimated predictor values are used
data{
vector[360] y;
vector[360] predictor_obs;
int id[360];
}
parameters{
vector[120] a_predictor;
real a_predictor_bar;
real<lower=0> sigma_predictor_a;
real<lower=0> sigma_predictor;
real a_y;
real b;
real<lower=0> sigma;
}
model{
vector[360] mu;
vector[360] predictor_true;
sigma ~ exponential( 1 );
b ~ normal( 0 , 0.5 );
a_y ~ normal( 0 , 1 );
sigma_predictor ~ exponential( 1 );
sigma_predictor_a ~ exponential( 1 );
a_predictor_bar ~ normal( 0 , 1 );
a_predictor ~ normal( a_predictor_bar , sigma_predictor_a );
for ( i in 1:360 ) {
predictor_true[i] = a_predictor[id[i]];
}
predictor_obs ~ normal( predictor_true , sigma_predictor );
for ( i in 1:360 ) {
mu[i] = a_y + b * predictor_true[i];
}
y ~ normal( mu , sigma );
}