A crash course in Bayesian mixed models with brms (Lesson 5)

Introduction and review

Part 5 of a practical introduction to fitting Bayesian multilevel models in R and Stan
Uses the brms package (Bayesian Regression Models using Stan)

Slides and text version of lessons are online
Fill in code in the worksheet (replace ... with code)
You can always copy and paste code from text version of lesson if you fall behind

Relationship between prior, likelihood, and posterior
How to sample from the posterior with Hamiltonian Monte Carlo
How to use posterior distributions from our models to make predictions and test hypotheses
GLMMs with non-normal response distributions and random intercepts and slopes
Models with residuals correlated in space and time
Spatial generalized additive mixed models
Beta regression for proportion data
Zero-inflated and hurdle models for data with zeros

Explore in more detail how link functions work in mixed models
Discover how every parameter in a nonlinear mixed model can have fixed and random effects

Image (c) Allison Horst

We’ve seen binary response data in lesson 2
Nominal data (for example: blood group) has statistical models but we won’t cover them today

Lots of data in ag science are on an ordinal rating scale
Common in phenotyping: capture a complex phenotype, or a quantitative one that is too time-consuming to measure directly, with a simple score
- Rate disease severity on a scale from 1 to 5
- Instead of counting hairs on a stem, rate low/medium/high

It is faster to score phenotypes on a rating scale than to physically measure the underlying variable directly
But this speed has a price: more data points are required for similarly precise estimates, than if you had a continuous variable
“Quick and dirty” way: treat the numerical ratings as continuous variables with normally distributed error
Often OK, but we have statistical models for ordinal response data!

For binary data we use a logistic model: binomial response distribution with logit link function
We’re trying to predict a probability, which has to be between 0 and 1
But our linear predictor can be any value, positive or negative
Logit, or log-odds, link function maps a probability to a scale that can have any value
- \(\text{logit } p = \log \frac{p}{1-p}\)

Logit model for binary outcome, we only need to predict single probability
Probability of success + probability of failure = 1
If we know one, we know the other
Model this process as following a binomial distribution which only has one parameter, \(\theta\)
Example: coin flips

If we have \(n\) ordered categories, we have to estimate \(n - 1\) probabilities
They all have to be between 0 and 1
Instead of a binomial distribution, use a multinomial distribution.
\(n - 1\) parameters to estimate, the probability of each outcome except for the last
Example: dice roll modeled with a multinomial distribution, parameter \(p = (\frac{1}{6}, \frac{1}{6}, \frac{1}{6}, \frac{1}{6}, \frac{1}{6}, \frac{1}{6})\)

We still use a logit link function but it’s a bit more complicated
Example: we have four categories, 1, 2, 3, 4
- Estimate probability of category 1, probability of 2 or less, and probability of 3 or less
- (Probability of 4 or less = 1)
These are cumulative probabilities so we’ll fit a cumulative logit model

\(n - 1\) link functions, one for each category except the last
Map our cumulative probabilities for each category to the linear predictor scale
If probabilities are \(p_{y \le 1}\), \(p_{y \le 2}\), and \(p_{y \le 3}\), the link functions are:
- \(\log \frac {p_{y = 1}}{1-p_{y = 1}}\)
- \(\log \frac {p_{y = 1} + p_{y = 2}}{1-(p_{y = 1} + p_{y = 2})}\)
- \(\log \frac {p_{y = 1} + p_{y = 2} + p_{y = 3}}{1-(p_{y = 1} + p_{y = 2} + p_{y = 3})}\)

karcher.turfgrass from agridat package
Effect of management style and nitrogen level on turfgrass color
Completely randomized design
- Four management levels are different ways of applying N (at the surface or injected at different depths)
- Two nitrogen levels: low and high N addition