Dear all,

I’ve recently transitioned from JAGS to STAN for the purpose of fitting large-scale categorical logit models. I’m running into a problem with the model below, in that the model fails to initialize:

Rejecting initial value:

Chain 1: Log probability evaluates to log(0), i.e. negative infinity.

Chain 1: Stan can’t start sampling from this initial value.

I’m not sure how to diagnose this - could you point me in the right direction?

Thanks,

Roberto

Code:

```
data {
int N; // N. Obs
int Nv; // N. Levels
int Nz; // N. Covariates
int Nj; // N. Choices
int V[N]; // Choices
matrix[N,Nz] Z; // Covariates
int var_id[Nz]; // Level-id per Covariate
}
parameters {
matrix[Nz,Nj] Beta; // Regression Coefficients of interest
matrix[Nz,Nj - 1] Theta_raw; // Prior Mean
vector[Nj-1] alpha_raw; // Choice-Intercept
cholesky_factor_corr[Nz] L_Omega; // off-diagonal Cholesky Prior
vector<lower=0>[Nv] L_sigma_raw; // Prior Scale
}
transformed parameters {
vector[Nj] alpha;
vector<lower=0>[Nz] L_sigma;
matrix[Nz, Nz] L_Sigma;
matrix[N,Nj] ZB;
matrix[Nz,Nj] Theta;
vector[Nj] nu[N];
for(z in 1:Nz) L_sigma[z] = L_sigma_raw[var_id[z]];
L_Sigma = diag_pre_multiply(L_sigma, L_Omega);
for (z in 1:Nz) {
Theta[z,1] = 0;
for (j in 1:(Nj-1)) {
Theta[z,j + 1] = Theta_raw[z ,j];
} }
alpha[1] = 0;
for (j in 1:(Nj-1)) {
alpha[j + 1] = alpha_raw[j];
}
ZB = Z * Beta;
for (i in 1:N) for (j in 1:Nj) nu[i][j] = alpha[j] + ZB[i,j];
}
model {
//priors
for(j in 1:(Nj-1)) Theta_raw[,j] ~ normal(0, 10);
for (v in 1:Nv) L_sigma_raw[v] ~ cauchy(0,2.5);
L_Omega ~ lkj_corr_cholesky(4);
//likelihood
for(j in 1:Nj) Beta[,j] ~ multi_normal_cholesky(Theta[,j], L_Sigma);
for (i in 1:N) V[i] ~ categorical_logit(nu[i]);
}'
```

Notes:

All Zs are standardised; var_id specifies the groups of Z that share a variance parameter for regularisation purposes; Nz is approximately 1,000; N is approximately 30,000.