ST308 Bayesian Inference

Class 10


Philipp Sterzinger

p.sterzinger@lse.ac.uk

04 April 2025

Recap

Multilevel logistic regression

Consider the following model

\[ y_{ij} \sim \text{Bernoulli}(\pi_{ij}), \quad \log\left(\frac{\pi_{ij}}{1-\pi_{ij}}\right) = \beta_0 + a_i + x_{ij}^\top \beta + z_{ij}^\top b_i \]

\[ \begin{aligned} \beta_0 & \sim \pi_{\beta_0}(\beta) \\ \beta & \sim \pi_\beta(\beta) \\ (a_i, b_i^\top)^\top &\overset{\textrm{ind}}{\sim} \mathcal{N}(0, \Sigma) \\ \Sigma &\sim \pi_{\Sigma}(\Sigma) \end{aligned} \]

stan_glmer

Key inputs:

  • formula: Like stan_lmer
  • family: Like stan_glm
  • prior: Prior on \(\beta\)
  • prior_intercept: Prior on \(\beta_0\)
  • prior_covariance: Prior on \(\Sigma\)
  • prior_aux: depends on family, none for logistic regression

Prior predictive distribution

Multilevel linear model

\[ y_{ij} = \beta_0 + a_i + x_{ij}^\top \beta + z_{ij}^\top b_i + \varepsilon_{ij}, \quad (j = 1,\ldots, n_i; i = 1,\ldots, K) \]

  • \(\theta \in \Theta\) are the model parameters

  • \(f(y \mid \theta, X)\) is the joint likelihood of \(y = (y_1, \ldots, y_n)^\top\) given \(\theta\), \(X\)

Joint distribution of responses \(y\) according to model and prior choice, given covariates \(X\):

\[ \textrm{pd}(y \mid X) = \int_{\Theta} f(y \mid \theta, X) \pi(\theta) d\theta \]

Prior predictive distribution

In rstanrm done using stan_glmer(..., PD = TRUE) and posterior_predict(stanreg_object)

Algorithm

For \(i\) in \(1: \textrm{N\_samples}\) do:

  1. Draw \(\beta_0, \beta, \sigma, \Sigma\)
  2. Draw \((a_i, b_i^\top)^\top \sim \mathcal{N}(0, \Sigma)\)
  3. Draw \(y_{ij}\mid_{\beta_0, \beta, a_i, b_i, X} \sim \mathcal{N}(\beta_0 + a_i + x_{ij}^\top \beta + z_{ij}^\top b_i, \sigma^2)\)

Output:

\(\textrm{N\_samples} \times n\) matrix, where each row is a sample of \(y\) given observed covariates \(X\), priors, model

Plot marginal density of observations \(y\) and simulated samples to check if values of pd are sensible (pp_check from bayesplot)

Posterior predictive distribution

Same as prior predictive distribution but with posterior distribution instead of prior

  • distribution of responses according to model, prior choice, \(X\) after observing data \((y, X)\)
  • sanity check: is the model & prior choice useful for replicating the observed data?

Examples

Multilevel logistic regression

Code
library(rstanarm) 

RadonData <- read.csv('radon.csv', header = TRUE)
RadonData$county_num <- as.numeric(as.factor(RadonData$county))
RadonData$y <- RadonData$logRadon>median(RadonData$logRadon)

bmod2 <- stan_glmer(y ~ floor + (1 | county_num), 
                    data = RadonData,
                    family = binomial(link = "logit"))
Code
prior_summary(bmod2)
Priors for model 'bmod2' 
------
Intercept (after predictors centered)
 ~ normal(location = 0, scale = 2.5)

Coefficients
  Specified prior:
    ~ normal(location = 0, scale = 2.5)
  Adjusted prior:
    ~ normal(location = 0, scale = 6.7)

Covariance
 ~ decov(reg. = 1, conc. = 1, shape = 1, scale = 1)
------
See help('prior_summary.stanreg') for more details

Multilevel logistic regression

Code
summary(bmod2)

Model Info:
 function:     stan_glmer
 family:       binomial [logit]
 formula:      y ~ floor + (1 | county_num)
 algorithm:    sampling
 sample:       4000 (posterior sample size)
 priors:       see help('prior_summary')
 observations: 919
 groups:       county_num (85)

Estimates:
                                            mean   sd   10%   50%   90%
(Intercept)                                0.5    0.1  0.3   0.5   0.7 
floor                                     -1.5    0.2 -1.8  -1.5  -1.2 
b[(Intercept) county_num:1]               -1.0    0.7 -2.0  -1.0  -0.1 
b[(Intercept) county_num:2]               -1.2    0.3 -1.7  -1.2  -0.8 
b[(Intercept) county_num:3]                0.4    0.7 -0.5   0.4   1.4 
b[(Intercept) county_num:4]                0.7    0.6 -0.1   0.7   1.5 
b[(Intercept) county_num:5]               -0.5    0.7 -1.4  -0.5   0.3 
b[(Intercept) county_num:6]                0.1    0.7 -0.8   0.1   1.0 
b[(Intercept) county_num:7]                1.3    0.6  0.6   1.3   2.0 
b[(Intercept) county_num:8]                0.6    0.7 -0.3   0.6   1.5 
b[(Intercept) county_num:9]               -0.5    0.5 -1.2  -0.5   0.2 
b[(Intercept) county_num:10]               0.5    0.7 -0.4   0.5   1.4 
b[(Intercept) county_num:11]               0.4    0.7 -0.4   0.4   1.3 
b[(Intercept) county_num:12]               0.3    0.7 -0.6   0.2   1.2 
b[(Intercept) county_num:13]              -0.3    0.6 -1.1  -0.3   0.5 
b[(Intercept) county_num:14]               0.2    0.5 -0.5   0.1   0.8 
b[(Intercept) county_num:15]              -0.4    0.7 -1.3  -0.4   0.5 
b[(Intercept) county_num:16]              -0.7    0.8 -1.8  -0.7   0.3 
b[(Intercept) county_num:17]               0.2    0.7 -0.6   0.3   1.1 
b[(Intercept) county_num:18]              -0.4    0.5 -1.0  -0.4   0.3 
b[(Intercept) county_num:19]              -0.3    0.3 -0.6  -0.3   0.1 
b[(Intercept) county_num:20]               0.6    0.8 -0.3   0.6   1.6 
b[(Intercept) county_num:21]               0.9    0.6  0.1   0.9   1.7 
b[(Intercept) county_num:22]              -0.5    0.7 -1.4  -0.5   0.3 
b[(Intercept) county_num:23]              -0.6    0.8 -1.6  -0.5   0.5 
b[(Intercept) county_num:24]               0.7    0.6 -0.1   0.7   1.5 
b[(Intercept) county_num:25]               0.7    0.5  0.1   0.7   1.4 
b[(Intercept) county_num:26]              -0.2    0.2 -0.5  -0.2   0.1 
b[(Intercept) county_num:27]               0.0    0.6 -0.8   0.0   0.8 
b[(Intercept) county_num:28]              -0.1    0.7 -0.9  -0.1   0.8 
b[(Intercept) county_num:29]              -1.0    0.8 -1.9  -0.9   0.0 
b[(Intercept) county_num:30]              -1.3    0.5 -2.0  -1.2  -0.6 
b[(Intercept) county_num:31]               0.4    0.7 -0.5   0.4   1.3 
b[(Intercept) county_num:32]              -0.2    0.7 -1.1  -0.2   0.6 
b[(Intercept) county_num:33]               0.8    0.8 -0.2   0.7   1.8 
b[(Intercept) county_num:34]               0.4    0.8 -0.5   0.4   1.4 
b[(Intercept) county_num:35]              -1.2    0.7 -2.1  -1.2  -0.3 
b[(Intercept) county_num:36]               0.7    0.8 -0.3   0.6   1.7 
b[(Intercept) county_num:37]              -1.3    0.6 -2.1  -1.3  -0.6 
b[(Intercept) county_num:38]               0.1    0.7 -0.8   0.1   0.9 
b[(Intercept) county_num:39]               0.6    0.7 -0.3   0.5   1.4 
b[(Intercept) county_num:40]               0.4    0.7 -0.4   0.4   1.3 
b[(Intercept) county_num:41]               1.3    0.7  0.4   1.2   2.2 
b[(Intercept) county_num:42]               0.3    0.8 -0.8   0.3   1.3 
b[(Intercept) county_num:43]              -0.3    0.6 -1.0  -0.3   0.5 
b[(Intercept) county_num:44]               0.0    0.6 -0.8   0.0   0.8 
b[(Intercept) county_num:45]              -0.2    0.5 -0.8  -0.2   0.5 
b[(Intercept) county_num:46]              -0.5    0.6 -1.3  -0.5   0.4 
b[(Intercept) county_num:47]               0.1    0.8 -0.9   0.1   1.0 
b[(Intercept) county_num:48]              -0.7    0.6 -1.4  -0.7   0.1 
b[(Intercept) county_num:49]               0.1    0.5 -0.5   0.1   0.8 
b[(Intercept) county_num:50]               0.3    0.8 -0.8   0.3   1.4 
b[(Intercept) county_num:51]               0.8    0.7 -0.2   0.7   1.7 
b[(Intercept) county_num:52]               0.6    0.8 -0.3   0.6   1.6 
b[(Intercept) county_num:53]              -0.3    0.7 -1.2  -0.3   0.7 
b[(Intercept) county_num:54]              -0.4    0.4 -0.9  -0.3   0.2 
b[(Intercept) county_num:55]               0.0    0.6 -0.7   0.0   0.7 
b[(Intercept) county_num:56]               0.4    0.7 -0.5   0.4   1.4 
b[(Intercept) county_num:57]              -0.9    0.6 -1.7  -0.9  -0.1 
b[(Intercept) county_num:58]               0.4    0.7 -0.5   0.4   1.3 
b[(Intercept) county_num:59]               0.1    0.7 -0.8   0.1   0.9 
b[(Intercept) county_num:60]              -0.1    0.8 -1.1  -0.1   0.8 
b[(Intercept) county_num:61]              -0.8    0.4 -1.3  -0.8  -0.3 
b[(Intercept) county_num:62]               0.5    0.6 -0.3   0.5   1.3 
b[(Intercept) county_num:63]               0.3    0.7 -0.6   0.2   1.2 
b[(Intercept) county_num:64]               0.7    0.6  0.0   0.7   1.5 
b[(Intercept) county_num:65]              -0.1    0.8 -1.1  -0.2   0.8 
b[(Intercept) county_num:66]               0.0    0.5 -0.6   0.0   0.6 
b[(Intercept) county_num:67]               0.1    0.5 -0.5   0.1   0.7 
b[(Intercept) county_num:68]              -1.0    0.6 -1.7  -0.9  -0.2 
b[(Intercept) county_num:69]               0.3    0.7 -0.6   0.2   1.1 
b[(Intercept) county_num:70]              -1.7    0.3 -2.1  -1.7  -1.4 
b[(Intercept) county_num:71]               0.0    0.4 -0.6   0.0   0.5 
b[(Intercept) county_num:72]               0.2    0.6 -0.5   0.2   0.9 
b[(Intercept) county_num:73]               0.5    0.8 -0.5   0.4   1.5 
b[(Intercept) county_num:74]              -0.7    0.7 -1.6  -0.7   0.2 
b[(Intercept) county_num:75]               0.3    0.7 -0.6   0.3   1.2 
b[(Intercept) county_num:76]               0.4    0.7 -0.5   0.4   1.3 
b[(Intercept) county_num:77]               0.4    0.6 -0.4   0.3   1.1 
b[(Intercept) county_num:78]              -0.2    0.7 -1.0  -0.2   0.7 
b[(Intercept) county_num:79]              -0.5    0.7 -1.5  -0.5   0.4 
b[(Intercept) county_num:80]              -0.3    0.3 -0.7  -0.3   0.1 
b[(Intercept) county_num:81]               1.0    0.8  0.0   0.9   2.0 
b[(Intercept) county_num:82]               0.3    0.8 -0.8   0.3   1.3 
b[(Intercept) county_num:83]               0.5    0.5 -0.2   0.5   1.1 
b[(Intercept) county_num:84]               0.3    0.5 -0.3   0.3   1.0 
b[(Intercept) county_num:85]              -0.1    0.8 -1.1  -0.1   0.8 
Sigma[county_num:(Intercept),(Intercept)]  0.8    0.3  0.5   0.8   1.2 

Fit Diagnostics:
           mean   sd   10%   50%   90%
mean_PPD 0.5    0.0  0.5   0.5   0.5  

The mean_ppd is the sample average posterior predictive distribution of the outcome variable (for details see help('summary.stanreg')).

MCMC diagnostics
                                          mcse Rhat n_eff
(Intercept)                               0.0  1.0  1745 
floor                                     0.0  1.0  4874 
b[(Intercept) county_num:1]               0.0  1.0  4073 
b[(Intercept) county_num:2]               0.0  1.0  4375 
b[(Intercept) county_num:3]               0.0  1.0  4954 
b[(Intercept) county_num:4]               0.0  1.0  5137 
b[(Intercept) county_num:5]               0.0  1.0  5325 
b[(Intercept) county_num:6]               0.0  1.0  5401 
b[(Intercept) county_num:7]               0.0  1.0  4205 
b[(Intercept) county_num:8]               0.0  1.0  5692 
b[(Intercept) county_num:9]               0.0  1.0  4925 
b[(Intercept) county_num:10]              0.0  1.0  5438 
b[(Intercept) county_num:11]              0.0  1.0  5554 
b[(Intercept) county_num:12]              0.0  1.0  5836 
b[(Intercept) county_num:13]              0.0  1.0  5204 
b[(Intercept) county_num:14]              0.0  1.0  4735 
b[(Intercept) county_num:15]              0.0  1.0  4295 
b[(Intercept) county_num:16]              0.0  1.0  4819 
b[(Intercept) county_num:17]              0.0  1.0  5209 
b[(Intercept) county_num:18]              0.0  1.0  5123 
b[(Intercept) county_num:19]              0.0  1.0  3505 
b[(Intercept) county_num:20]              0.0  1.0  4308 
b[(Intercept) county_num:21]              0.0  1.0  4916 
b[(Intercept) county_num:22]              0.0  1.0  5283 
b[(Intercept) county_num:23]              0.0  1.0  4394 
b[(Intercept) county_num:24]              0.0  1.0  5705 
b[(Intercept) county_num:25]              0.0  1.0  4161 
b[(Intercept) county_num:26]              0.0  1.0  2892 
b[(Intercept) county_num:27]              0.0  1.0  5249 
b[(Intercept) county_num:28]              0.0  1.0  5967 
b[(Intercept) county_num:29]              0.0  1.0  4505 
b[(Intercept) county_num:30]              0.0  1.0  4371 
b[(Intercept) county_num:31]              0.0  1.0  5745 
b[(Intercept) county_num:32]              0.0  1.0  5741 
b[(Intercept) county_num:33]              0.0  1.0  5895 
b[(Intercept) county_num:34]              0.0  1.0  5201 
b[(Intercept) county_num:35]              0.0  1.0  4272 
b[(Intercept) county_num:36]              0.0  1.0  4228 
b[(Intercept) county_num:37]              0.0  1.0  4876 
b[(Intercept) county_num:38]              0.0  1.0  4620 
b[(Intercept) county_num:39]              0.0  1.0  4836 
b[(Intercept) county_num:40]              0.0  1.0  5380 
b[(Intercept) county_num:41]              0.0  1.0  4250 
b[(Intercept) county_num:42]              0.0  1.0  4414 
b[(Intercept) county_num:43]              0.0  1.0  5636 
b[(Intercept) county_num:44]              0.0  1.0  5888 
b[(Intercept) county_num:45]              0.0  1.0  5165 
b[(Intercept) county_num:46]              0.0  1.0  5675 
b[(Intercept) county_num:47]              0.0  1.0  5597 
b[(Intercept) county_num:48]              0.0  1.0  5485 
b[(Intercept) county_num:49]              0.0  1.0  5439 
b[(Intercept) county_num:50]              0.0  1.0  4481 
b[(Intercept) county_num:51]              0.0  1.0  4472 
b[(Intercept) county_num:52]              0.0  1.0  4893 
b[(Intercept) county_num:53]              0.0  1.0  5658 
b[(Intercept) county_num:54]              0.0  1.0  4472 
b[(Intercept) county_num:55]              0.0  1.0  4771 
b[(Intercept) county_num:56]              0.0  1.0  5350 
b[(Intercept) county_num:57]              0.0  1.0  5390 
b[(Intercept) county_num:58]              0.0  1.0  4712 
b[(Intercept) county_num:59]              0.0  1.0  5005 
b[(Intercept) county_num:60]              0.0  1.0  5322 
b[(Intercept) county_num:61]              0.0  1.0  4597 
b[(Intercept) county_num:62]              0.0  1.0  5375 
b[(Intercept) county_num:63]              0.0  1.0  5858 
b[(Intercept) county_num:64]              0.0  1.0  5096 
b[(Intercept) county_num:65]              0.0  1.0  5559 
b[(Intercept) county_num:66]              0.0  1.0  4725 
b[(Intercept) county_num:67]              0.0  1.0  6318 
b[(Intercept) county_num:68]              0.0  1.0  4951 
b[(Intercept) county_num:69]              0.0  1.0  5740 
b[(Intercept) county_num:70]              0.0  1.0  2973 
b[(Intercept) county_num:71]              0.0  1.0  4528 
b[(Intercept) county_num:72]              0.0  1.0  5484 
b[(Intercept) county_num:73]              0.0  1.0  4704 
b[(Intercept) county_num:74]              0.0  1.0  5422 
b[(Intercept) county_num:75]              0.0  1.0  5074 
b[(Intercept) county_num:76]              0.0  1.0  5004 
b[(Intercept) county_num:77]              0.0  1.0  5234 
b[(Intercept) county_num:78]              0.0  1.0  6184 
b[(Intercept) county_num:79]              0.0  1.0  5329 
b[(Intercept) county_num:80]              0.0  1.0  3965 
b[(Intercept) county_num:81]              0.0  1.0  4378 
b[(Intercept) county_num:82]              0.0  1.0  5400 
b[(Intercept) county_num:83]              0.0  1.0  4916 
b[(Intercept) county_num:84]              0.0  1.0  5714 
b[(Intercept) county_num:85]              0.0  1.0  5222 
Sigma[county_num:(Intercept),(Intercept)] 0.0  1.0  1268 
mean_PPD                                  0.0  1.0  4342 
log-posterior                             0.3  1.0   911 

For each parameter, mcse is Monte Carlo standard error, n_eff is a crude measure of effective sample size, and Rhat is the potential scale reduction factor on split chains (at convergence Rhat=1).

Multilevel logistic regression

Code
ci95 <- posterior_interval(bmod2, prob = 0.95)
round(ci95, 2)
                                           2.5% 97.5%
(Intercept)                                0.24  0.83
floor                                     -1.96 -1.07
b[(Intercept) county_num:1]               -2.60  0.28
b[(Intercept) county_num:2]               -1.92 -0.59
b[(Intercept) county_num:3]               -1.05  1.85
b[(Intercept) county_num:4]               -0.43  1.96
b[(Intercept) county_num:5]               -2.02  0.86
b[(Intercept) county_num:6]               -1.37  1.49
b[(Intercept) county_num:7]                0.20  2.49
b[(Intercept) county_num:8]               -0.81  2.01
b[(Intercept) county_num:9]               -1.62  0.52
b[(Intercept) county_num:10]              -0.74  1.83
b[(Intercept) county_num:11]              -0.86  1.75
b[(Intercept) county_num:12]              -1.04  1.62
b[(Intercept) county_num:13]              -1.47  0.92
b[(Intercept) county_num:14]              -0.83  1.14
b[(Intercept) county_num:15]              -1.89  0.96
b[(Intercept) county_num:16]              -2.36  0.76
b[(Intercept) county_num:17]              -1.17  1.57
b[(Intercept) county_num:18]              -1.42  0.66
b[(Intercept) county_num:19]              -0.81  0.29
b[(Intercept) county_num:20]              -0.81  2.23
b[(Intercept) county_num:21]              -0.24  2.19
b[(Intercept) county_num:22]              -1.89  0.70
b[(Intercept) county_num:23]              -2.23  0.95
b[(Intercept) county_num:24]              -0.44  1.89
b[(Intercept) county_num:25]              -0.28  1.78
b[(Intercept) county_num:26]              -0.71  0.25
b[(Intercept) county_num:27]              -1.27  1.26
b[(Intercept) county_num:28]              -1.35  1.22
b[(Intercept) county_num:29]              -2.49  0.43
b[(Intercept) county_num:30]              -2.39 -0.25
b[(Intercept) county_num:31]              -0.86  1.75
b[(Intercept) county_num:32]              -1.61  1.07
b[(Intercept) county_num:33]              -0.65  2.35
b[(Intercept) county_num:34]              -1.05  1.96
b[(Intercept) county_num:35]              -2.72  0.08
b[(Intercept) county_num:36]              -0.87  2.21
b[(Intercept) county_num:37]              -2.68 -0.19
b[(Intercept) county_num:38]              -1.26  1.42
b[(Intercept) county_num:39]              -0.70  1.95
b[(Intercept) county_num:40]              -0.91  1.79
b[(Intercept) county_num:41]              -0.04  2.75
b[(Intercept) county_num:42]              -1.30  1.97
b[(Intercept) county_num:43]              -1.46  0.83
b[(Intercept) county_num:44]              -1.20  1.22
b[(Intercept) county_num:45]              -1.16  0.81
b[(Intercept) county_num:46]              -1.71  0.77
b[(Intercept) county_num:47]              -1.44  1.55
b[(Intercept) county_num:48]              -1.88  0.47
b[(Intercept) county_num:49]              -0.83  1.19
b[(Intercept) county_num:50]              -1.32  1.96
b[(Intercept) county_num:51]              -0.63  2.31
b[(Intercept) county_num:52]              -0.81  2.24
b[(Intercept) county_num:53]              -1.74  1.17
b[(Intercept) county_num:54]              -1.19  0.48
b[(Intercept) county_num:55]              -1.17  1.18
b[(Intercept) county_num:56]              -1.04  1.95
b[(Intercept) county_num:57]              -2.30  0.29
b[(Intercept) county_num:58]              -0.91  1.90
b[(Intercept) county_num:59]              -1.29  1.45
b[(Intercept) county_num:60]              -1.67  1.34
b[(Intercept) county_num:61]              -1.54 -0.10
b[(Intercept) county_num:62]              -0.70  1.90
b[(Intercept) county_num:63]              -1.12  1.66
b[(Intercept) county_num:64]              -0.36  1.90
b[(Intercept) county_num:65]              -1.69  1.35
b[(Intercept) county_num:66]              -0.98  0.96
b[(Intercept) county_num:67]              -0.85  1.04
b[(Intercept) county_num:68]              -2.19  0.20
b[(Intercept) county_num:69]              -1.07  1.69
b[(Intercept) county_num:70]              -2.27 -1.24
b[(Intercept) county_num:71]              -0.81  0.76
b[(Intercept) county_num:72]              -0.82  1.37
b[(Intercept) county_num:73]              -1.08  2.13
b[(Intercept) county_num:74]              -2.04  0.64
b[(Intercept) county_num:75]              -1.20  1.66
b[(Intercept) county_num:76]              -0.89  1.82
b[(Intercept) county_num:77]              -0.79  1.57
b[(Intercept) county_num:78]              -1.54  1.11
b[(Intercept) county_num:79]              -1.96  0.79
b[(Intercept) county_num:80]              -0.89  0.35
b[(Intercept) county_num:81]              -0.50  2.61
b[(Intercept) county_num:82]              -1.31  1.97
b[(Intercept) county_num:83]              -0.52  1.54
b[(Intercept) county_num:84]              -0.65  1.33
b[(Intercept) county_num:85]              -1.66  1.40
Sigma[county_num:(Intercept),(Intercept)]  0.38  1.49
Code
#launch_shinystan(bmod2, ppd = FALSE)

Multilevel logistic regression

Activity 1

Fit another version of the model above using random effect on the slope and/or putting different priors on the remaining models parameters

Prior predictive distribution

\[ \textrm{height}_{ij} = \beta_0 + \beta_1 \textrm{sex}_{ij} + a_i + \varepsilon_{ij} \,, \]

where the indexing variable \(i\) is the nationality of an individual,

Need priors on \(\beta_0, \beta, \sigma_y, \sigma_r\)

Code
set.seed(123)
N <- 200 

nations <- c("A", "B", "D", "E")
nationality <- sample(nations, N, replace = TRUE) 

sigma_re <- 10
r_es <- rnorm(length(nations), mean = 0, sd = sigma_re)
r_intercepts <- r_es[as.numeric(as.factor(nationality))]

sex <- as.numeric(runif(N) < .5)

sex_eff <- 10 
base_intercept <- 160
sigma_err <- 2.5

height <- base_intercept + sex_eff * sex + r_intercepts + rnorm(N, mean = 0, sd = sigma_err) 

dat <- data.frame(height, sex, nationality) 

summary(dat) 
     height           sex       nationality       
 Min.   :146.7   Min.   :0.00   Length:200        
 1st Qu.:157.4   1st Qu.:0.00   Class :character  
 Median :163.4   Median :1.00   Mode  :character  
 Mean   :163.1   Mean   :0.54                     
 3rd Qu.:168.2   3rd Qu.:1.00                     
 Max.   :177.7   Max.   :1.00                     

Prior predictive distribution

Code
library(bayesplot) 
bhmod1_prior <- stan_glmer(height ~ sex + (1 | nationality), 
                     family = gaussian(link = "identity"), 
                     data = dat,
                     prior_intercept = normal(160, .1), 
                     prior = normal(10, .1), 
                     prior_aux = normal(2.5 / 0.8),
                     prior_covariance = decov(regularization = 0.1),
                     prior_PD = TRUE)

bhmod2_prior <- stan_glmer(height ~ sex + (1 | nationality), 
                     family = gaussian(link = "identity"), 
                     data = dat,
                     prior_intercept = normal(0, 100), 
                     prior = normal(100, 0.5), 
                     prior_PD = TRUE) 

Good priors

Code
pp_check(bhmod1_prior) 

Poor priors

Code
pp_check(bhmod2_prior) 

Posterior predictive distribution

Code
library(rstanarm) 
library(bayesplot) 
bhmod1 <- stan_glmer(height ~ sex + (1 | nationality), 
                     family = gaussian(link = "identity"), 
                     data = dat,
                     prior_intercept = normal(160, .1), 
                     prior = normal(10, .1), 
                     prior_aux = normal(2.5 / 0.8),
                     prior_covariance = decov(regularization = 0.1),
                     prior_PD = FALSE)

bhmod2 <- stan_glmer(height ~ sex + (1 | nationality), 
                     family = gaussian(link = "identity"), 
                     data = dat,
                     prior_intercept = normal(0, 100), 
                     prior = normal(100, 0.5), 
                     prior_PD = FALSE) 

Good priors

Code
pp_check(bhmod1) 

Poor priors

Code
pp_check(bhmod2)