I'm trying to create a multilevel ordinal logistic regression model in Stan and the following code would seem to work, in the sense that Stan seems to convergence to sensible answers:

stanmodel <- '
data {
  int<lower=2> K;  // ordinal response with 4 values, 3 cutpoints
  int<lower=0> N;  // number of measurements
  int<lower=1,upper=K> y[N]; // response

  int<lower=0> Ntests;         // number of groups
  int<lower=1,upper=Ntests> tests[N];  // groups

parameters {
  // population cutpoints and associated variance.
  ordered[K-1]  CutpointsMean;
  real<lower=0> CutpointsSigma[K-1];   

  ordered[K-1]  Cutpoints[Ntests];   // ordinal response cutpoints for groups

model {

  CutpointsSigma ~ exponential(1);
  CutpointsMean  ~ normal(2, 3);

  for (i in 1:Ntests) {
    Cutpoints[i][1] ~ normal(CutpointsMean[1] , CutpointsSigma[1]);
    Cutpoints[i][2] ~ normal(CutpointsMean[2] , CutpointsSigma[2]);
    Cutpoints[i][3] ~ normal(CutpointsMean[3] , CutpointsSigma[3]);

  for (i in 1:N)
    y[i] ~ ordered_logistic(0, Cutpoints[tests[i]]);


I have removed the part relating to the covariates for clarity.

'CutpointsMean' and 'CutpointsSigma' define the population global ordinal response while Cutpoints[i][1:3] is the ordinal response for group i.

The idea is that for example the first 'cutpoint' of each group is generated from a normal distribution centered on the first 'cutpoint' of the overall population.

The second 'cutpoint' of each group is generated from a normal distribution centered on the second 'cutpoint' of the overall population and so on.

As Cutpoints[i] is an ordered vector of 3 elements, what happens when I write directly into Cutpoints[i][2] ?

Is the write operation rejected if the constraints are not satisfied or simply the entry is written in the vector and the results is sorted?

Is this the correct way of modelling a multilevel ordinal response in Stan?

  • 1
    $\begingroup$ You can't write to the elements of something declared in the parameters block; you will get a parser error. What you have looks like legal syntax, but I don't think a normal prior for the cutpoints makes much sense. At a minimum, you need a truncated normal. $\endgroup$ – Ben Goodrich Jul 31 '16 at 18:29
  • $\begingroup$ HI Ben, thanks for your reply. I have also created another model where the first cutpoint is normal and the other two cutpoints are half normal. The cutpoints are then recovered by adding the offsets back together: ordererd_cuts[1] <- cuts[1]; ordererd_cuts[2] <- cuts[1] + cuts[2]; ordererd_cuts[3] <- cuts[1] + cuts[2] + cuts[3]; where cuts[2] and cuts[3] are positive half normal $\endgroup$ – Anton Aug 1 '16 at 11:49
  • $\begingroup$ Try stan_polr from rstanarm! $\endgroup$ – Zach May 2 '18 at 1:06

For those interested, the answer to this problem requires the Ordered Inverse Transform, the complete answer can be found here:


| cite | improve this answer | |
  • $\begingroup$ Link not working $\endgroup$ – Sandu Ursu Aug 7 '19 at 17:51
  • $\begingroup$ It still works for me, I joined the stan users mailing list years ago where the post belongs and I am logged in my google account. I presume one of these factors could affect your ability to see it. $\endgroup$ – Anton Sep 1 '19 at 17:20

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.