R code to simulate data and run Stan:
library(rstan)
rstan_options(auto_write = TRUE)
options(mc.cores = parallel::detectCores())
# Create data
set.seed(12)
n <- 1000
y <- rnorm(n, mean = 3.5, sd = 1.19)
c1k <- c(-Inf, 1.5, ci, k - 0.5, Inf)
ycat <- as.numeric(cut(y, c1k))
d.stan <- list(N=n, K=max(ycat), y=ycat)
a <- list(mu = 3.5, logs = log(1.19), c = c(2.11, 3.68))
initl <- list(a,a,a)
fit <- stan(file = "mufit.stan",
data = d.stan, iter = 2000, chains = 3, init = initl,
control = list(adapt_delta = 0.95))
Stan model:
data {
int<lower=0> N; // # obs
int<lower=2> K; // # scale levels
int<lower=1,upper=K> y[N]; // scale response
}
parameters {
ordered[K-3] c; // cuts for latent normal
real mu; // latent normal parameters
real logs;
}
transformed parameters{
real s = exp(logs);
}
model {
vector[K] theta; // prob vector for each level of response
mu ~ normal(0.5*(K+1), 3);
logs ~ normal(0, 1.2);
for (i in 1:(K-3)){
c[i] ~ normal(1.5 + i, 1);
}
for (n in 1:N) {
theta[1] = Phi(1.5);
theta[2] = Phi(c[1])-Phi(1.5);
theta[3] = Phi(c[2])-Phi(c[1]);
theta[4] = Phi(4.5)-Phi(c[2]);
theta[5] = 1-Phi(4.5);
y[n] ~ categorical(theta);
}
}
I am trying to implement the model given in Analyzing Ordinal Data with Metric Models: What Could Possibly Go Wrong?.
Despite the fairly informative priors and initial values, the posteriors for mu and s are quite wide. The model doesn’t take long to run, finishes without any complaints, chains look ok, Rhat seems fine… no clue what the issue might be here.
