Ordinal Outcome and Ordinal Predictor in brms

@paul.buerkner: I have been trying to model an ordinal outcome with multiple ordinal predictors that I was introduced to brms Tutorial and Ordinal Model in Psychology and Modelling monotonic effects of ordinal predictors in Bayesian regression models by @avehtari
Before I apply the models to my own data and to be pedagogic on solving problems I tried to reproduce the Cumulative model (CN) using opinions about funding stem cell research data as below:

set.seed(1)
stemcell=read.csv("GSS2006.csv",check.names=FALSE)

stemcell$belief<- mapvalues(stemcell$fund, from = c("moderate","fundamentalist","liberal"), to = c(1, 2,3))
stemcell$rating <- mapvalues(stemcell$scresrch, c("definitely should not fund such research", "probably should not fund such research", "probably should fund such research","definitely should fund such research"),to=c(1,2,3,4))

stemcell$rating<-as.ordered(stemcell$rating)
stemcell$belief<-as.ordered(stemcell$belief)

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

summary(fit_sc1)


Population-Level Effects: 
             Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept[1]    -0.54      0.05    -0.63    -0.45 1.00     4306     3113
Intercept[2]    -0.25      0.04    -0.33    -0.16 1.00     4841     3708
Intercept[3]     1.04      0.05     0.94     1.15 1.00     5392     3247
belief.L        -0.09      0.06    -0.22     0.04 1.00     4614     3150
belief.Q         0.34      0.07     0.21     0.47 1.00     5082     3116

1- I do not get the same results as published in the paper:

Population-Level Effects:
$  ##                          Estimate  Est.Error   l-95% CI   u-95% CI   Eff.Sample   Rhat
    ## Intercept[1]            -1.25      0.08          -1.42         -1.10           2681         1.00
    ## Intercept[2]            -0.64       0.07          -0.78         -0.50          3629         1.00
    ## Intercept[3]             0.57       0.07          0.43          0.71           3461          1.00 
## belieffundamentalist   -0.24     0.09           -0.43        -0.06            3420         1.00
 ## beliefliberal                0.31      0.09           0.13          0.50             3381        1.00

2- I used the monotonic effect of the ordinal predictor on the same data and the results are so different and I wonder how we interpret these two different results and how to understand which one is a better model.

fit_sc2<-brm(
  formula = rating ~ 1 +mo( belief),
  data = stemcell,
  family = cumulative("probit")
)

summary(fit_sc2)

Population-Level Effects: 
             Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept[1]    -0.71      0.08    -0.86    -0.56 1.00     1975     1657
Intercept[2]    -0.42      0.07    -0.56    -0.28 1.00     2040     1919
Intercept[3]     0.86      0.08     0.70     1.01 1.00     2751     2284
mobelief        -0.14      0.05    -0.23    -0.05 1.00     2081     1846

Simplex Parameters: 
             Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
mobelief1[1]     0.91      0.11     0.67     1.00 1.00     2598     2105
mobelief1[2]     0.09      0.11     0.00     0.33 1.00     2598     2105

3- Model comparison using loo(fit_sc1,fit_sc2)

```Model comparisons:
              elpd_diff   se_diff
fit_sc1     0.0          0.0   
fit_sc2     -7.4         4.1 

I am not sure which model I should go with.
In summary: I cannot reproduce the same results, using the same syntax and dataset, and I do not know why using “monotonic” causes such differences and which model is better.
Thanks in advance for your help.
Shabnam

             Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept[1]    -1.27      0.06    -1.39    -1.16 1.00     2692     2662
Intercept[2]    -0.66      0.05    -0.75    -0.57 1.00     4228     3306
Intercept[3]     0.55      0.05     0.46     0.64 1.00     4982     3133
belief.L           0.22      0.07     0.09     0.35 1.00     3531     2892
belief.Q          0.33      0.07     0.20     0.46 1.00     4361     2966

I finally could reproduce it. The problem was coding the “belief”!

I ran the model with belief as monotonic predictor using the same data and the results is:

Population-Level Effects: 
             Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept[1]    -1.11      0.07    -1.24    -0.98 1.00     2642     2771
Intercept[2]    -0.51      0.06    -0.61    -0.40 1.00     3592     2785
Intercept[3]     0.70      0.06     0.59     0.81 1.00     3676     3284
mobelief         0.22      0.04     0.13     0.30 1.00     2775     2590

Simplex Parameters: 
             Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
mobelief1[1]     0.07      0.06     0.00     0.22 1.00     2900     1900
mobelief1[2]     0.93      0.06     0.78     1.00 1.00     2900     1900

I assume the difference between two results is because the belief is not a monotonic predictor?
Any hint on this is appreciated. Thanks.

Great you could solve the reproducibility issue. You seem to have used the polynomial contrasts for the belief in the first. This is not ideal I think as the results are hard to interprete. I would try out to compare dummy coding to the monotonic effects for better interpretation of the parameters.

1 Like

@paul.buerkner
What I understood from monotonic paper, we can use ordinal predictor coded for example as 1,2,3,4 and use “mo”. Let’s say one of the Linkert responses, do I need to create dummy variables to interpret the estimate? I am always interested to compare 2, 3, or 4 to 1.
If I use the dummy variables then no need to use syntax “mo”?
I got a bit confused.
Thanks

@paul.buerkner
In the tutorial paper, in the cumulative model, we do not consider the ordinal predictor as a monotonic effect, that is why to interpret the estimate you suggests to create dummy variables and put them into the model. But in the monotonic paper, we can keep ordinal variable coded as, for example, 1 to 4 and we use “mo”. Did I get it right?

You did get it right I think. And indeed, if you want to interprete in terms of dummy variables, use dummy coding instead of monotonic effects.

Thank you.