Hi,

I have a question about model interpetation. I originally posed a question about the type of model I should run, and Solomon was very kind to help Setting priors and interpreting ordered logistic model).

I am trying to determine whether three different continents (Africa, Asia, Europe) have different numbers of life expectancies at different levels (low, med, high) (this data and question are entirely fabricated)
Essentially: does continent predict how many low/med/high life expectancies in that group?

My data:

``````#load libraries
library(brms)
library("gapminder")
library("dplyr")

#Edit data
data(gapminder)
gapminder <-gapminder
#select columns of interest
gapminder <-gapminder %>%
select(continent,lifeExp, pop)
#select continents of interest
gapminder <-gapminder %>%
filter(continent == c("Africa","Asia","Europe"))
#Make ordered factors
gapminder\$lifeExp<- cut(gapminder\$lifeExp, 3, labels=c("low","med","high"), ordered_result=TRUE)
``````

After advice, I have settled on an ordinal cumulative model with probit link:

``````mod2selection_fit <- brm(
formula = lifeExp ~ 1 + continent,
data = gapminder,
family = cumulative(probit),
prior = c(
prior(normal(0, 2), class = Intercept),
prior(normal(0, 1), class = b)))
``````

which gives:

``````Family: cumulative
Links: mu = probit; disc = identity
Formula: lifeExp ~ 1 + continent
Data: gapminder (Number of observations: 460)
Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
total post-warmup draws = 4000

Population-Level Effects:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept        0.10      0.08    -0.06     0.26 1.00     3951     2911
Intercept        1.47      0.11     1.26     1.69 1.00     2764     3112
continentAsia       1.17      0.13     0.91     1.42 1.00     2899     2919
continentEurope     3.03      0.21     2.64     3.44 1.00     2482     2857

Family Specific Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
disc     1.00      0.00     1.00     1.00   NA       NA       NA

Draws were sampled 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).
``````

My question is about interpreting this, what it means to the average person (me!)

Lets taker the continent Asia as an Example, the estimate is 1.17, and the 95%CI do not cross 0. Could I interpret this as Asia has 1.17 years more life expectancy than Africa? Or should I transform this as it is a probit link?Or is it in terms of risk? I’m trying to understand exactly what the numbers in the output are telling me.

Also, what if the 95%CI did NOT cross 0, then I couldn’t be confident in the estimate?

Finally, does anyone have any suggestions for how one would report this in a paper, is there a Bayes equivalent of the frequentist (p value, confidence intervals) reporting?

Thank you so much and Happy new year to all!

When you use the cumulative probit, the ordinal `LifeExp` values are expressed as an underlying latent variable on the standardized normal scale (mean = 0, sd = 1). In the case of your model, this means that Africa is the reference category, with latent `LifeExp` scores with a mean of zero. The `continentAsia` coefficient then suggests that Life expectancy in Asia is 1.17 standard scores higher than Africa, like a Cohen’s d of 1.17.

As to whether the 95% CI crosses zero, that just tells you whether zero is within the range of values the model considers 95% probable for a given parameter. You could think about this from an NHST-type perspective, but I’d encourage you not to. Rather, think in terms of effect sizes and the precision with which you estimate them.

But if you really want to take a p-value oriented approach, look into the `brms::hypothesis()` function. There are some examples of that approach floating around here on the forums.

2 Likes

Also, I just published a blog post on the Bayesian cumulative probit via brms at Notes on the Bayesian cumulative probit | A. Solomon Kurz. It’s not exactly a beginner’s post, but it is pretty comprehensive. The closest model in that post to your data would be `fit2`.

2 Likes

Thank you so much, this is so helpful! Especially thinking of the coefficient as something like cohens d. Can I just follow up on reporting, do you/anyone reading this know of a standard way to report bayesian model outcomes? Do people generally give the estimate, and error? If anyone knows of papers (I have googled!) that have these types of models in them that I could use as a reference, that would be super helpful!

You might also find the parameters and bayestestR summaries useful for summarizing posterior results