@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