Hello Stan users and exports,
I have a problem as the tittle describes.
My model is a little complex so I paste the parts that cause error.
data{
int<lower=0> y[3,18];
real log_max; // log(1000)
//etc...
}
parameters{
vector<lower = 0>[2] sigma;
vector[2] U;
matrix[2,2] B;
matrix<upper=log_max>[2,18]log_y_2;
//etc...
}
transformed parameters{
simplex[3] theta[18];
matrix<upper=log_max>[3,18]log_y;
for(t in 1:18){
log_y[1,t]=log_y_2[1,t];
log_y[2,t]=log_y_2[2,t];
if ((exp(log_max)-exp(log_y_2[1,t])-exp(log_y_2[2,t]))<25)
log_y[3,t]=log(25);
else
log_y[3,t]=log(exp(log_max)-exp(log_y_2[1,t])-exp(log_y_2[2,t]));
}
for (t in 1:18){
theta[t]=softmax(exp(log_y[,t]));
}
}
model{
for (i in 1:18) {
y[,i] ~ multinomial(theta[i,]);
}
for(t in 2:18){
target+= multi_normal_lpdf(log_y_2[,t] |U+B*log_y_2[,t-1],diag_matrix(sigma));
//etc...
}
After many attempts, I found that log_y_3[3,t] in transformed parameters is the problem. If I write
if ((exp(log_max)-exp(log_y_2[1,t])-exp(log_y_2[2,t]))<25)
log_y[3,t]=log(25);
else
log_y[3,t]=log(exp(log_max)-exp(log_y_2[1,t])-exp(log_y_2[2,t]));
there would be a log(0) error.
However, if I write
if ((exp(log_max)-exp(log_y_2[1,t])-exp(log_y_2[2,t]))<25)
log_y[3,t]=log(25);
else
log_y[3,t]=log(50);
It will work.
But if change 50 to 1000 or 500, it will report an error once again.
Actually, I need log_y[3,t]=log(exp(log_max)-exp(log_y_2[1,t])-exp(log_y_2[2,t]));. Because log_y is a matrix that elements in each column follow multinomial distribution. How can I successfully achieve this goal?
Thank you in advance