Ordinal response with monotonic predictor

ssalimi · February 19, 2020, 11:52am

Please also provide the following information in addition to your question:

Operating System: Windows
brms Version: GitHub version

I just recently learned about the monotonic effect of the ordinal predictor. In the manuscript, the examples are for a linear response.
I need to apply it to an ordinal outcome in repeated data and wrote the below syntax

getprior<-get_prior( bodn ~ mo(chf)+1|id,family=cumulative(),
                    data =bodn)

fit_bodn<-brm(
  formula = bodn~ mo(chf)+1|idno,family=cumulative(),
  data =bodn,prior=getprior)

But, I just got the results for intercepts not for the monotonic effect.


Group-Level Effects: 
~idno (Number of levels: 898) 
                              Estimate   Est.Error   l-95%     u-95%  Rhat 
sd(Intercept)          4.44             0.15      4.16        4.75   1.00      
sd(sevchf)              1.06            12.93     0.36        41.91 1.00   
cor(chf)                -0.01              0.58      -0.95       0.95   1.00       

Population-Level Effects: 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept[1]    -12.26      0.57   -13.43   -11.21 1.00     1399     2691
Intercept[2]     -8.76      0.32    -9.41    -8.17 1.00     1066     2270
Intercept[3]     -6.14      0.22    -6.57    -5.72 1.00      703     1945
Intercept[4]     -3.91      0.18    -4.26    -3.56 1.01      548     1452
Intercept[5]     -2.10      0.17    -2.44    -1.77 1.01      495     1282
Intercept[6]     -0.29      0.16    -0.61     0.03 1.01      454     1107
Intercept[7]      1.22      0.17     0.90     1.56 1.01      472     1185
Intercept[8]      2.84      0.18     2.50     3.20 1.01      546     1116
Intercept[9]      4.42      0.19     4.06     4.79 1.01      662     1762
Intercept[10]     5.97      0.22     5.55     6.39 1.01      814     1714
Intercept[11]     7.72      0.26     7.22     8.23 1.00     1126     2207
Intercept[12]     9.82      0.34     9.16    10.48 1.00     1610     2240
Intercept[13]    11.90      0.48    10.99    12.87 1.00     2381     2930
Intercept[14]    14.86      0.98    13.21    17.02 1.00     4245     3347

Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

What I am doing wrong here? How should I model the ordinal response with monotonic ordinal predictor n repeated data?
Thanks

torkar · February 19, 2020, 11:58am

Should you have mo(chf) + (1| id)? Not sure how much the parenthesis matter though…

What is the output from running: prior_summary(fit_bodn) on your model?

ssalimi · February 19, 2020, 12:12pm

It seems I cannot get a reasonable prior by get_prior of the model.
the

prior_summary (fit_bodn) 
student_t(3, 0, 10)

just for intercept and the same for sd for id
It did not provide prior for the predictor.
Not sure why,
I also used

prior <- prior(normal(2.1, 10), class = "b") +
  prior(normal(1.1, 10), class = "Intercept")+   prior(dirichlet(1, 1), class = "simo", coef = "mochf1")

Not sure why I don’t get the right prior by get_prior!!! I checked with adding parenthesis as well.

torkar · February 19, 2020, 12:14pm

Strange. Can you paste us the output from this:

get_prior(bodn ~ mo(chf) + (1 | id), family = cumulative, data = bodn)

ssalimi · February 19, 2020, 12:24pm

getprior
                 prior     class               coef 
1                              b                                               
2                              b    mochf                            
3  student_t(3, 0, 10) Intercept                                               
4                      Intercept                  1                            
5                      Intercept                 10                            
6                      Intercept                 11                            
7                      Intercept                 12                            
8                      Intercept                 13                            
9                      Intercept                 14                            
10                     Intercept                  2                            
11                     Intercept                  3                            
12                     Intercept                  4                            
13                     Intercept                  5                            
14                     Intercept                  6                            
15                     Intercept                  7                            
16                     Intercept                  8                            
17                     Intercept                  9                            
18 student_t(3, 0, 10) sd                                               
19                               sd                          idno                      
20                               sd          Intercept  idno                      
21                               simo      mochf1

torkar · February 19, 2020, 12:28pm

p <- get_prior(bodn ~ mo(chf) + (1 | id), family = cumulative, data = bodn)
p$prior[1] <- "normal(0,1)" # sets a N(0,1) on \beta
p$prior[3] <- "normal(0,2.5)" # for each cutpoint
p$prior[18] <- "cauchy(0,2.5)" # for sd
p$prior[21] <- "dirichlet(2,2)" # sets weakly regularizing dirichlet for 1-3 Likert scale response

Then use prior = p when you call brms. Of course, at least do pp_check() when you sample with sample_prior="only" to see that the priors are sane.

ssalimi · February 19, 2020, 12:44pm

I try the model again and get back here to update.
-Thanks

ssalimi · February 19, 2020, 12:50pm

The problem is I do not get the similar pattern as you wrote:
The response is as blow:

 p$prior[1]
[1] ""
> p$prior[3] 
[1] "student_t(3, 0, 10)"
> p$prior[1] 
[1] ""
> p$prior[3] 
[1] "student_t(3, 0, 10)"
> p$prior[18]
[1] "student_t(3, 0, 10)"
> p$prior[21] 
[1] ""

:(
I can send some of the data, if it helps to understand why!

torkar · February 19, 2020, 12:52pm

I set the first row in p to “normal(0,1)”. What happens when you cut and paste exactly what I have written above (I changed things to better fit your scenario) and then simply write p in the terminal?

ssalimi · February 19, 2020, 2:28pm

I ran the model using the above priors. The good thing is that I got simplex results. However, the Rhats are not good.
here is the results

Data: bodi (Number of observations: 2622) 
Samples: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
         total post-warmup samples = 4000

Group-Level Effects: 
~idno (Number of levels: 898) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     4.13      0.14     3.86     4.41 1.01      457     1116

Population-Level Effects: 
                  Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept[1]        -11.12      0.50   -12.11   -10.17 1.01      949     1834
Intercept[2]         -8.07      0.30    -8.67    -7.51 1.01      703     1816
Intercept[3]         -5.67      0.21    -6.10    -5.28 1.02      422     1057
Intercept[4]         -3.56      0.18    -3.92    -3.20 1.03      262      803
Intercept[5]         -1.82      0.17    -2.15    -1.50 1.04      219      514
Intercept[6]         -0.06      0.16    -0.38     0.24 1.04      125      498
Intercept[7]          1.40      0.16     1.08     1.71 1.04      130      560
Intercept[8]          2.97      0.17     2.64     3.30 1.03      172      687
Intercept[9]          4.49      0.19     4.13     4.87 1.03      259      983
Intercept[10]         5.97      0.21     5.57     6.40 1.02      373     1267
Intercept[11]         7.60      0.25     7.14     8.08 1.02      554     1719
Intercept[12]         9.49      0.31     8.92    10.10 1.01      650     2034
Intercept[13]        11.21      0.41    10.43    12.04 1.01     1412     2826
Intercept[14]        13.22      0.66    12.04    14.59 1.00     2611     2885
       mochf          0.50      0.08     0.34     0.65 1.00      915     1622

Simplex Parameters: 
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
      mochf1[1]        0.16      0.09     0.03     0.36 1.00     2083     2314
      mochf1[2]        0.84      0.09     0.64     0.97 1.00     2083     2314

Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

The good news is I got the simplex result but not good Rhat
Should I change the priors’ values?
Thanks

ssalimi · February 19, 2020, 4:03pm

I modified the prior for cut-points to (0,5) and reran the syntax and it worked well with reasonable Rhat and convergence.

Group-Level Effects: 
~idno (Number of levels: 898) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     4.32      0.15     4.04     4.62 1.00      711     1152

Population-Level Effects: 
                  Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept[1]        -11.74      0.54   -12.87   -10.72 1.00     1560     2403
Intercept[2]         -8.35      0.31    -8.97    -7.74 1.00     1307     2258
Intercept[3]         -5.80      0.23    -6.25    -5.36 1.00      859     1851
Intercept[4]         -3.61      0.19    -3.98    -3.23 1.01      639     1427
Intercept[5]         -1.81      0.18    -2.16    -1.46 1.01      560     1154
Intercept[6]         -0.01      0.17    -0.34     0.33 1.01      528     1297
Intercept[7]          1.50      0.17     1.17     1.85 1.01      529     1133
Intercept[8]          3.12      0.18     2.78     3.49 1.01      578     1077
Intercept[9]          4.70      0.20     4.32     5.08 1.01      677     1627
Intercept[10]         6.25      0.22     5.82     6.69 1.01      830     1734
Intercept[11]         7.98      0.26     7.48     8.50 1.01     1071     2221
Intercept[12]        10.02      0.33     9.38    10.67 1.00     1736     2649
Intercept[13]        12.00      0.45    11.16    12.90 1.00     2733     2877
Intercept[14]        14.64      0.87    13.10    16.56 1.00     4483     2673
        mochf         0.52      0.08     0.35     0.68 1.00     2247     2494

Simplex Parameters: 
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
          mochf1[1]     0.17      0.09     0.03     0.36 1.00     3518     3121
          mochf1[2]     0.83      0.09     0.64     0.97 1.00     3518     3121

Thanks for your help and guidance,

torkar · February 19, 2020, 9:26pm

Hi,

I’d still be careful. Have you run the model with sample_prior = "only" and used pp_check() to see how the priors look like?

Make sure to check that \widehat{R} < 1.01, that neff is more than 0.1-0.2, and that the trace plots look alright. It seems you have a Bulk_ESS which is close to 0.1 (if you ran with four chains each with 2000 iterations).

ssalimi · February 20, 2020, 6:40am

Blockquote I’d still be careful. Have you run the model with sample_prior = "only" and used pp_check() to see how the priors look like?

Blockquote
We expect the weak priors and therefore a weak posterior coefficient, right?

ssalimi · February 20, 2020, 7:02am

Blockquote brms

Ordinal response with monotonic predictor

Interfaces brms

](Profile - ssalimi - The Stan Forums)

ssalimi

torkar

18h

The problem is I do not get the similar pattern as you wrote:
The response is as blow:

 p$prior[1]
[1] ""
> p$prior[3] 
[1] "student_t(3, 0, 10)"
> p$prior[1] 
[1] ""
> p$prior[3] 
[1] "student_t(3, 0, 10)"
> p$prior[18]
[1] "student_t(3, 0, 10)"
> p$prior[21] 
[1] ""

:(
I can send some of the data, if it helps to understand why!

SolutionReply

torkar

18h

I set the first row in p to “normal(0,1)”. What happens when you cut and paste exactly what I have written above (I changed things to better fit your scenario) and then simply write p in the terminal?

SolutionReply

ssalimi

17h

I ran the model using the above priors. The good thing is that I got simplex results. However, the Rhats are not good.
here is the results

Data: bodi (Number of observations: 2622) 
Samples: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
         total post-warmup samples = 4000

Group-Level Effects: 
~idno (Number of levels: 898) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     4.13      0.14     3.86     4.41 1.01      457     1116

Population-Level Effects: 
                  Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept[1]        -11.12      0.50   -12.11   -10.17 1.01      949     1834
Intercept[2]         -8.07      0.30    -8.67    -7.51 1.01      703     1816
Intercept[3]         -5.67      0.21    -6.10    -5.28 1.02      422     1057
Intercept[4]         -3.56      0.18    -3.92    -3.20 1.03      262      803
Intercept[5]         -1.82      0.17    -2.15    -1.50 1.04      219      514
Intercept[6]         -0.06      0.16    -0.38     0.24 1.04      125      498
Intercept[7]          1.40      0.16     1.08     1.71 1.04      130      560
Intercept[8]          2.97      0.17     2.64     3.30 1.03      172      687
Intercept[9]          4.49      0.19     4.13     4.87 1.03      259      983
Intercept[10]         5.97      0.21     5.57     6.40 1.02      373     1267
Intercept[11]         7.60      0.25     7.14     8.08 1.02      554     1719
Intercept[12]         9.49      0.31     8.92    10.10 1.01      650     2034
Intercept[13]        11.21      0.41    10.43    12.04 1.01     1412     2826
Intercept[14]        13.22      0.66    12.04    14.59 1.00     2611     2885
       mochf          0.50      0.08     0.34     0.65 1.00      915     1622

Simplex Parameters: 
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
      mochf1[1]        0.16      0.09     0.03     0.36 1.00     2083     2314
      mochf1[2]        0.84      0.09     0.64     0.97 1.00     2083     2314

Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

The good news is I got the simplex result but not good Rhat
Should I change the priors’ values?
Thanks

SolutionReply

ssalimi

15h

I modified the prior for cut-points to (0,5) and reran the syntax and it worked well with reasonable Rhat and convergence.

Group-Level Effects: 
~idno (Number of levels: 898) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     4.32      0.15     4.04     4.62 1.00      711     1152

Population-Level Effects: 
                  Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept[1]        -11.74      0.54   -12.87   -10.72 1.00     1560     2403
Intercept[2]         -8.35      0.31    -8.97    -7.74 1.00     1307     2258
Intercept[3]         -5.80      0.23    -6.25    -5.36 1.00      859     1851
Intercept[4]         -3.61      0.19    -3.98    -3.23 1.01      639     1427
Intercept[5]         -1.81      0.18    -2.16    -1.46 1.01      560     1154
Intercept[6]         -0.01      0.17    -0.34     0.33 1.01      528     1297
Intercept[7]          1.50      0.17     1.17     1.85 1.01      529     1133
Intercept[8]          3.12      0.18     2.78     3.49 1.01      578     1077
Intercept[9]          4.70      0.20     4.32     5.08 1.01      677     1627
Intercept[10]         6.25      0.22     5.82     6.69 1.01      830     1734
Intercept[11]         7.98      0.26     7.48     8.50 1.01     1071     2221
Intercept[12]        10.02      0.33     9.38    10.67 1.00     1736     2649
Intercept[13]        12.00      0.45    11.16    12.90 1.00     2733     2877
Intercept[14]        14.64      0.87    13.10    16.56 1.00     4483     2673
        mochf         0.52      0.08     0.35     0.68 1.00     2247     2494

Simplex Parameters: 
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
          mochf1[1]     0.17      0.09     0.03     0.36 1.00     3518     3121
          mochf1[2]     0.83      0.09     0.64     0.97 1.00     3518     3121

Thanks for your help and guidance,

Solution

2

Reply

torkar

10h

Hi,

I’d still be careful. Have you run the model with sample_prior = "only" and used pp_check() to see how the priors look like?

Make sure to check that ˆR < 1.01, that neff is more than 0.1-0.2, and that the trace plots look alright. It seems you have a Bulk_ESS which is close to 0.1 (if you ran with four chains each with 2000 iterations).

Blockquote

torkar · February 20, 2020, 7:15am

Hi,

please post the output of this:

pp_check(fit_bodn, type = "bars")

ssalimi · February 20, 2020, 7:20am

please post the output of this:

pp_check(fit_bodn, type = "bars")

PPcheck_bar

ssalimi · February 20, 2020, 7:29am

@torkar
just one person has bodn=15, is this a problem?

torkar · February 20, 2020, 7:35am

Your priors are a bit heavy on low and high numbers, while they don’t cover, e.g., category 2 very much. I’d like to see a more balance view.

I presume your pp_check() on the full model (not with priors only) is ok?

ssalimi · February 20, 2020, 8:39am

[quote=“torkar, post:18, topic:13165, full:true”]
Your priors are a bit heavy on low and high numbers, while they don’t cover, e.g., category 2 very much. I’d like to see a more balance view.

I presume your pp_check() on the full model (not with priors only) is ok?

It seems so, this is the pp_check (fit_bodn) for the full model!
What do you think? Thanks a lot for all your time. I learned many points!
bodn_chf_ppCheck

torkar · February 20, 2020, 8:48am

I think it looks ok :) I would now compare that model with a model with no monotonic effects,

bodn ~ chf + (1 | id)

and a general intercept model,

bodn ~ 1

using LOO: https://rdrr.io/cran/brms/man/loo_compare.brmsfit.html

That is however completely up to you :)

Topic		Replies	Views
Monotonic effects with ordered predictors and ordinal response variables Modeling brms	4	782	June 10, 2022
Conditional monotonicity in cumulative probit model brms ordinal-response	2	376	December 21, 2023
Interpretion of Monotonic Effects brms interpret-results	2	824	December 10, 2018
Ordinal Outcome and Ordinal Predictor in brms brms	6	1483	February 22, 2020
How to simulate monotonic effects with an ordinal outcome brms brms	3	703	October 4, 2022

Ordinal response with monotonic predictor

Ordinal response with monotonic predictor

Related topics