Calculate Credible Intervals for group means with ± 0.5 sum contrasts?

General
1

I used ± 0.5 sum contrasts for a mixed effect logistic regression with categorical binary predictors. I know how to get the coefficient estimates for each factor level from the estimated mean in the model summary, but how do I calculate the credible intervals? is there any package in R for this?

2

The credible intervals for the coefficients should also be in the model summary. If you are interested in credible interval of the difference between groups, you should be able to get it by extracting individual samples for both coefficients, then for each sample multiply the coefficients by ± 0.5 (the same way you did to create your design matrix) and add them. This will give you posterior samples of the difference and then you can use quantile to compute the credible interval of this quantity.

Does that make sense?

3

What I want to know is the coefficient and the credible interval for the reference level of a binary predictor. So if the predictor X has two levels A and B, with sum contrasts I get the average of A and B in the summary together with the CI. So if A is the reference level, to get the estimate for A alone I add the Intercept to the coefficient estimate (model estimate+ intercept= estimate for A). My question is how to get the credible interval for the estimate for A. Hope that is clearer.

4

So, with half-sum contrasts (i.e. -.5/+.5) the intercept parameter reflects the mean of the two conditions and the effect parameter reflects the difference between the two conditions. To get the posterior for a given condition, you can do this either in the Stan model in the generated quantities section as:

generated quantities{
    real condition_A = intercept - effect;
    real condition_B = intercept + effect;
}

Or, if the model was already sampled and it took a long time so you don’t want to re-sample, you can do it in R via:

intercept = extract(from_stan_sampling,par='intercept',permute=F)[[1]]
effect = extract(from_stan_sampling,par='effect',permute=F)[[1]]
condition_A = intercept-effect
condition_B = intercept+effect
quantile(condition_A,probs=c(.025,.5,.975))
quantile(condition_B,probs=c(.025,.5,.975))
1 Like
5

Thank you, this is really helpful. Does it matter if there are more predictors than just one? They are all binary and what about interactions? Is the calculation the same ? or does it change anything?

6

What I showed is the simple case of wanting to look at the condition means while accounting for just a simple single-2-level predictor. For anything more complicated you can still do things by hand but you’ll probably save yourself some headache by learning a package that does things automatically for you like tidyBayes. My own package, ezStan, has a get_condition_post() function that obtains all condition posteriors if you have a certain type of hierarchical model (see: https://www.youtube.com/watch?v=G75WrqPyQuQ&list=PLu77iLvsj_GNmWqDdX-26kZ56dCZqkzBO&index=23)