Ordered_logistic_lpmf

Starting from the example of Bürkner & Vuorre (Ordinal Regression Modells in Psychology: A Tutorial)page 12, I use ``stancode(fit_sc1)``` to retrieve the following model:

fit_sc1 <- brm(formula = rating ~ 1+ belief, data = stemcell, family = cumulative())

to understand what’s going on in Ordinal Models.

Messing with the code I changed

for (n in 1:N) {
    target += ordered_logistic_lpmf(Y[n]| mu[n], temp_Intercept);
  }

to

 target += ordered_logistic_lpmf(Y| mu, temp_Intercept);

to try the vectorization.With the vectorization I get 10x slower code, 3951 divergences, \hat{R} \gt 4 and low ESS.
I’m totally OK with non-vectorized code… but what is happening here and what am I missing here?

This is due to a bug in the vectorised code (my fault actually!). For that combination of inputs, the derivatives aren’t being calculated properly. There’s more context on the issue here.

This bug has been fixed in the 2.20 release. In the interim, you can also fix this yourself.

First find the ordered_logistic header file with the command:

system.file("include", "stan", "math", "prim", "mat", "prob", 
            "ordered_logistic_lpmf.hpp", 
            package = "StanHeaders")

If you open up that file, on line 139 you’ll see:

        ops_partials.edge2_.partials_vec_[n](0) = d;

Change that to:

        ops_partials.edge2_.partials_vec_[n](0) += d;

i.e., change “=” to “+=”

This code was overwriting the derivatives w.r.t to the vector of cutpoints, rather than accumulating them for each observation. Let me know if that works for you or not

3 Likes

It works like a charms! And moreover… it’s 3x faster than the unvectorized code. Thanks @andrjohns!

2 Likes

Thanks! I have a graded response model that suddenly started acting strange over summer. Good to know why!