ST308 Bayesian Inference

Class 9

Philipp Sterzinger

p.sterzinger@lse.ac.uk

28 March 2025

Recap

Normal linear model

Code

library(ggplot2)

set.seed(123)
n <- 100 

xs <- cbind(rep(1,n), rnorm(n)) 
betas <- 2 * rnorm(2) 

ys <- xs %*% betas + rnorm(n, sd = 0.25) 

dat <- data.frame(ys, xs) 
dat <- sort_by(dat, xs[,2]) 

ggplot(dat, aes(x = X2, y = ys)) + 
    geom_point(shape = 1) +                                    
    geom_abline(data = data.frame(t(betas)), 
                aes(intercept = X1, slope = X2), 
                color = "steelblue") +              
    labs(title = "Normal Linear Model", x = "x", y = "y") +    
    theme_minimal()

\[ y_i = \beta_0 + x_i^\top \beta + \varepsilon_j, \quad \varepsilon_j \mid_{x_i} \sim \mathcal{N}(0, \sigma^2) \]

The effects \(\beta_0, \beta\) are the same for all observations \(i = 1,\ldots,n\)

Frequentist

treat \(\beta_0, \beta, \sigma^2\) as fixed
estimate \(\beta\) using ML/LS with lm

Bayesian

treat \(\beta_0, \beta, \sigma^2\) as random variables
approximate \(\pi(\beta_0, \beta, \sigma^2 \mid y, X)\) using MCMC with stan_glm

Normal linear model

Interactions

model different intercepts / slopes according to different levels of a factor variable
let \(\delta_i\) be a dummy variable that is \(1\) if factor variable is at certain level, \(0\) else
No interaction: \(x_i = x_{i1}\), Interaction: \(x_i = (\delta_i, x_{1i}, \delta_i \times x_{1i})^\top\)

mod3 <- lm(logRadon ~ floor * county, data = RadonData)
mod3 <- lm(logRadon ~ floor + county + floor:county, data = RadonData)

Code

set.seed(123)
n <- 100

xs <- rnorm(n)
factor <- runif(n) < 0.5
X_matrix <- cbind(1, factor, factor * xs, xs)
betas <- -1.5 * rnorm(4)
ys <- X_matrix %*% betas + rnorm(n, sd = 0.5)


dat <- data.frame(
  ys = as.vector(ys),
  Intercept = X_matrix[, 1],
  Factor = factor,
  Factor_X = factor * xs,
  X = xs
)


ggplot(dat, aes(x = X, y = ys, color = Factor)) + 
  geom_point(shape = 1, size = 2) +
  geom_abline(aes(intercept = betas[1] + betas[2], slope = betas[3] + betas[4]),
              color = "steelblue", linewidth = 1.2, linetype = "solid") +
  geom_abline(aes(intercept = betas[1], slope = betas[4]),
              color = "tomato", linewidth = 1.2, linetype = "solid") +
  labs(title = "Normal Linear Model w/ intercations",
       x = "x",
       y = "y",
       color = "Factor") +
  theme_minimal() +
  scale_color_manual(values = c("tomato", "steelblue"))

Hierarchichal model

Avoid covariate inflation associated with interactions, make the group-level effects implicit

\[ 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) \]

\(y_{ij}\) is \(j\)th observation in group \(i\)
\(x_{ij}\), \(z_{ij}\) are observed covariates
\(\beta_0, \beta\) are effects that are the same for all groups
\(a_i\), \(b_i\) are effects that vary amongst groups

\[ (a_i, b_i)^\top \sim \mathcal{N}(0, \Sigma), \quad \varepsilon_{ij} \sim \mathcal{N}(0,\sigma^2) \]

Hierarchichal model

Frequentist approach:

Write down likelihood of observations given unobserved random effects \(f(y_{ij}\mid \beta_0, \beta, a_i, b_i)\) and “integrate out” unobserved effects
Involves high-dimensional integrals, that are often untractable, need for approximations…

Bayesian approach:

Treat \(a_i, b_i\) as model parameters with known prior, add prior on hyperparameter \(\Sigma\)
Specify priors on \(\beta_0, \beta, \sigma^2, \Sigma\)

Case study

`stan_lmer`

For the hierarchichal normal linear model we use the stan_lmer from the rstanarm R package.

If we want to use a general GLM with hierarchichal structure, we use the stan_glmerfunction.

Make sure to check out the vignette on hierarchichal models

The functions are Bayesian equivalents to the lmer, glmer functions of the lme4 R package, so check these out for further info on syntax etc. The section on model formulas in the documentation might be of interest

`stan_lmer`

Key inputs:

formula: y ~ fexp + (grp_exp | group)
- can have multiple independent groups, e.g. y ~ fexp + (grp_exp_1 | group_1) + (grp_exp_2 | group_2)
- nested groups, e.g. y ~ fexp + (grp_exp | group_1/group_2) equivalent to y ~ fexp + (grp_exp | group_1) + (grp_exp | group_1:group_2)
- || for independent group-level effects
prior: Prior on fixed effects \(\beta\)
prior_intercept: Prior on \(\beta_0\)
prior_aux: Prior on \(\sigma^2\)
prior_covariance: Prior on \(\Sigma\)

Radon case study

EDA

Code

RadonData = read.csv('radon.csv',header=TRUE)
summary(RadonData)

    county              floor           logRadon      
 Length:919         Min.   :0.0000   Min.   :-2.3026  
 Class :character   1st Qu.:0.0000   1st Qu.: 0.6931  
 Mode  :character   Median :0.0000   Median : 1.3083  
                    Mean   :0.1665   Mean   : 1.2648  
                    3rd Qu.:0.0000   3rd Qu.: 1.8083  
                    Max.   :1.0000   Max.   : 3.8774

Radon case study

Floor variable

Code

library("RColorBrewer")
colours<-brewer.pal(n = 12, name = "Paired")

boxplot(logRadon ~ floor, data=RadonData,col=colours[c(2,6)],names= c('Basement','Floor'))

Radon case study

County variable

Code

RadonData$county_num = as.numeric(as.factor(RadonData$county))
#table(RadonData$county)

boxplot(logRadon ~ county_num, data = RadonData, xlab = 'Counties')

Radon case study

Linear model `floor`

Code

mod1 = lm(logRadon ~ floor,data = RadonData)
summary(mod1)


Call:
lm(formula = logRadon ~ floor, data = RadonData)

Residuals:
     Min       1Q   Median       3Q      Max 
-3.07857 -0.52950  0.02388  0.53970  2.51502 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  1.36241    0.02855  47.725   <2e-16 ***
floor       -0.58642    0.06996  -8.382   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.7901 on 917 degrees of freedom
Multiple R-squared:  0.07116,   Adjusted R-squared:  0.07015 
F-statistic: 70.25 on 1 and 917 DF,  p-value: < 2.2e-16

Radon case study

Linear model `floor + county`

Code

mod2 = lm(logRadon~floor+county,data=RadonData)
summary(mod2)


Call:
lm(formula = logRadon ~ floor + county, data = RadonData)

Residuals:
    Min      1Q  Median      3Q     Max 
-3.1794 -0.4449 -0.0065  0.4406  2.6083 

Coefficients:
                        Estimate Std. Error t value Pr(>|t|)    
(Intercept)              0.88724    0.36370   2.439  0.01492 *  
floor                   -0.68923    0.07062  -9.760  < 2e-16 ***
countyANOKA              0.04337    0.37723   0.115  0.90850    
countyBECKER             0.66225    0.55568   1.192  0.23369    
countyBELTRAMI           0.69991    0.45595   1.535  0.12515    
countyBENTON             0.56730    0.51374   1.104  0.26980    
countyBIG STONE          0.64955    0.55519   1.170  0.24235    
countyBLUE EARTH         1.13764    0.41198   2.761  0.00588 ** 
countyBROWN              1.10882    0.51404   2.157  0.03129 *  
countyCARLTON            0.15844    0.42996   0.368  0.71260    
countyCARVER             0.67877    0.46931   1.446  0.14847    
countyCASS               0.54085    0.48770   1.109  0.26776    
countyCHIPPEWA           0.86861    0.51404   1.690  0.09145 .  
countyCHISAGO            0.19748    0.46931   0.421  0.67402    
countyCLAY               1.11870    0.41192   2.716  0.00675 ** 
countyCLEARWATER         0.48106    0.51404   0.936  0.34963    
countyCOOK              -0.16970    0.62945  -0.270  0.78753    
countyCOTTONWOOD         0.37977    0.51495   0.737  0.46104    
countyCROW WING          0.27276    0.41947   0.650  0.51570    
countyDAKOTA             0.48478    0.37486   1.293  0.19628    
countyDODGE              0.92958    0.55519   1.674  0.09443 .  
countyDOUGLAS            0.86332    0.43671   1.977  0.04838 *  
countyFARIBAULT         -0.10877    0.46902  -0.232  0.81666    
countyFILLMORE           0.53209    0.62945   0.845  0.39817    
countyFREEBORN           1.22567    0.43660   2.807  0.00511 ** 
countyGOODHUE            1.07757    0.41198   2.616  0.00907 ** 
countyHENNEPIN           0.50555    0.37027   1.365  0.17250    
countyHOUSTON            0.89971    0.46902   1.918  0.05541 .  
countyHUBBARD            0.38868    0.48800   0.796  0.42599    
countyISANTI             0.20364    0.55519   0.367  0.71386    
countyITASCA             0.08265    0.42458   0.195  0.84571    
countyJACKSON            1.14926    0.48770   2.356  0.01868 *  
countyKANABEC            0.38164    0.51404   0.742  0.45804    
countyKANDIYOHI          1.18859    0.51404   2.312  0.02101 *  
countyKITTSON            0.72957    0.55568   1.313  0.18957    
countyKOOCHICHING       -0.01811    0.45595  -0.040  0.96832    
countyLAC QUI PARLE      2.06359    0.62945   3.278  0.00109 ** 
countyLAKE              -0.40058    0.43671  -0.917  0.35926    
countyLAKE OF THE WOODS  0.98897    0.51404   1.924  0.05471 .  
countyLE SUEUR           0.87629    0.48739   1.798  0.07255 .  
countyLINCOLN            1.43469    0.51374   2.793  0.00535 ** 
countyLYON               1.09182    0.44500   2.454  0.01435 *  
countyMAHNOMEN           0.49905    0.81249   0.614  0.53923    
countyMARSHALL           0.75144    0.43713   1.719  0.08598 .  
countyMARTIN             0.21958    0.45545   0.482  0.62985    
countyMCLEOD             0.43240    0.41544   1.041  0.29825    
countyMEEKER             0.35916    0.48770   0.736  0.46167    
countyMILLE LACS         0.08089    0.62945   0.129  0.89778    
countyMORRISON           0.29351    0.43671   0.672  0.50171    
countyMOWER              0.83955    0.41547   2.021  0.04363 *  
countyMURRAY             1.61419    0.81249   1.987  0.04728 *  
countyNICOLLET           1.29064    0.51404   2.511  0.01224 *  
countyNOBLES             1.05532    0.55519   1.901  0.05767 .  
countyNORMAN             0.39763    0.55494   0.717  0.47386    
countyOLMSTED            0.45389    0.39368   1.153  0.24927    
countyOTTER TAIL         0.75231    0.44500   1.691  0.09129 .  
countyPENNINGTON         0.30346    0.55568   0.546  0.58514    
countyPINE              -0.07443    0.46902  -0.159  0.87396    
countyPIPESTONE          0.99059    0.51374   1.928  0.05417 .  
countyPOLK               0.85212    0.51404   1.658  0.09776 .  
countyPOPE               0.41963    0.62945   0.667  0.50518    
countyRAMSEY             0.30923    0.38546   0.802  0.42264    
countyREDWOOD            1.11084    0.48739   2.279  0.02291 *  
countyRENVILLE           0.79973    0.55494   1.441  0.14992    
countyRICE               0.97679    0.42436   2.302  0.02159 *  
countyROCK               0.44811    0.62945   0.712  0.47673    
countyROSEAU             0.79505    0.41253   1.927  0.05429 .  
countySCOTT              0.93542    0.41544   2.252  0.02460 *  
countySHERBURNE          0.23784    0.44526   0.534  0.59337    
countySIBLEY             0.38762    0.51404   0.754  0.45102    
countyST LOUIS           0.03582    0.36957   0.097  0.92281    
countySTEARNS            0.63089    0.39131   1.612  0.10728    
countySTEELE             0.71674    0.43019   1.666  0.09607 .  
countySTEVENS            0.92172    0.62945   1.464  0.14348    
countySWIFT              0.13949    0.51404   0.271  0.78618    
countyTODD               0.85192    0.55494   1.535  0.12512    
countyTRAVERSE           1.13301    0.51374   2.205  0.02770 *  
countyWABASHA            0.95532    0.45545   2.098  0.03625 *  
countyWADENA             0.43390    0.48749   0.890  0.37369    
countyWASECA            -0.18110    0.51374  -0.353  0.72455    
countyWASHINGTON         0.47614    0.37887   1.257  0.20920    
countyWATONWAN           1.81287    0.55568   3.262  0.00115 ** 
countyWILKIN             1.35347    0.81249   1.666  0.09612 .  
countyWINONA             0.76904    0.41542   1.851  0.06449 .  
countyWRIGHT             0.77923    0.41560   1.875  0.06115 .  
countyYELLOW MEDICINE    0.32956    0.62945   0.524  0.60071    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.7265 on 833 degrees of freedom
Multiple R-squared:  0.2865,    Adjusted R-squared:  0.2137 
F-statistic: 3.936 on 85 and 833 DF,  p-value: < 2.2e-16

Radon case study

Linear model `floor * county`

Code

mod3 = lm(logRadon~floor * county, data = RadonData)
#dim(model.matrix(mod3))
summary(mod3)


Call:
lm(formula = logRadon ~ floor * county, data = RadonData)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.9561 -0.4108  0.0000  0.4204  2.8264 

Coefficients: (25 not defined because of singularities)
                              Estimate Std. Error t value Pr(>|t|)    
(Intercept)                    0.67561    0.41411   1.631 0.103201    
floor                          0.15730    0.82823   0.190 0.849420    
countyANOKA                    0.27753    0.42660   0.651 0.515524    
countyBECKER                   0.80599    0.82823   0.973 0.330781    
countyBELTRAMI                 1.04496    0.58565   1.784 0.074769 .  
countyBENTON                   0.66784    0.58565   1.140 0.254498    
countyBIG STONE                0.86118    0.58565   1.470 0.141840    
countyBLUE EARTH               1.20169    0.46299   2.595 0.009625 ** 
countyBROWN                    1.44404    0.65477   2.205 0.027718 *  
countyCARLTON                  0.21629    0.47818   0.452 0.651168    
countyCARVER                   1.11362    0.58565   1.902 0.057604 .  
countyCASS                     0.75248    0.52382   1.437 0.151258    
countyCHIPPEWA                 1.08024    0.54782   1.972 0.048978 *  
countyCHISAGO                  0.40911    0.50719   0.807 0.420126    
countyCLAY                     1.43878    0.47216   3.047 0.002388 ** 
countyCLEARWATER               0.57715    0.65477   0.881 0.378346    
countyCOOK                     0.04193    0.65477   0.064 0.948955    
countyCOTTONWOOD               1.25591    0.82823   1.516 0.129832    
countyCROW WING                0.31157    0.47818   0.652 0.514867    
countyDAKOTA                   0.70961    0.42451   1.672 0.095010 .  
countyDODGE                    1.14122    0.58565   1.949 0.051699 .  
countyDOUGLAS                  1.04074    0.48559   2.143 0.032406 *  
countyFARIBAULT                0.05341    0.52382   0.102 0.918809    
countyFILLMORE                 0.51831    0.82823   0.626 0.531626    
countyFREEBORN                 1.47644    0.49496   2.983 0.002945 ** 
countyGOODHUE                  1.14859    0.46299   2.481 0.013321 *  
countyHENNEPIN                 0.72844    0.42067   1.732 0.083741 .  
countyHOUSTON                  1.05467    0.54782   1.925 0.054569 .  
countyHUBBARD                  0.33317    0.65477   0.509 0.611012    
countyISANTI                   0.41527    0.58565   0.709 0.478488    
countyITASCA                   0.29428    0.46718   0.630 0.528948    
countyJACKSON                  1.36089    0.52382   2.598 0.009555 ** 
countyKANABEC                  0.59327    0.54782   1.083 0.279163    
countyKANDIYOHI                1.40022    0.54782   2.556 0.010779 *  
countyKITTSON                  1.27030    0.82823   1.534 0.125499    
countyKOOCHICHING             -0.31276    0.58565  -0.534 0.593470    
countyLAC QUI PARLE            2.10321    0.82823   2.539 0.011299 *  
countyLAKE                    -0.25633    0.48559  -0.528 0.597741    
countyLAKE OF THE WOODS        1.09747    0.65477   1.676 0.094121 .  
countyLE SUEUR                 0.99185    0.54782   1.811 0.070602 .  
countyLINCOLN                  1.76358    0.58565   3.011 0.002686 ** 
countyLYON                     1.20972    0.49496   2.444 0.014745 *  
countyMAHNOMEN                 0.71068    0.82823   0.858 0.391116    
countyMARSHALL                 1.88545    0.54782   3.442 0.000609 ***
countyMARTIN                   0.58595    0.50719   1.155 0.248327    
countyMCLEOD                   0.74099    0.47818   1.550 0.121644    
countyMEEKER                   0.57079    0.52382   1.090 0.276195    
countyMILLE LACS               1.08225    0.82823   1.307 0.191703    
countyMORRISON                 0.42516    0.48559   0.876 0.381551    
countyMOWER                    1.19267    0.46718   2.553 0.010874 *  
countyMURRAY                   1.82583    0.82823   2.204 0.027783 *  
countyNICOLLET                 1.50227    0.54782   2.742 0.006243 ** 
countyNOBLES                   1.26695    0.58565   2.163 0.030821 *  
countyNORMAN                   0.64177    0.65477   0.980 0.327322    
countyOLMSTED                  0.74264    0.44409   1.672 0.094875 .  
countyOTTER TAIL               0.83166    0.52382   1.588 0.112767    
countyPENNINGTON               1.09934    0.82823   1.327 0.184788    
countyPINE                     0.10945    0.52382   0.209 0.834543    
countyPIPESTONE                1.10456    0.58565   1.886 0.059662 .  
countyPOLK                     1.43341    0.65477   2.189 0.028883 *  
countyPOPE                     0.63126    0.65477   0.964 0.335303    
countyRAMSEY                   0.43824    0.43501   1.007 0.314037    
countyREDWOOD                  0.89095    0.54782   1.626 0.104282    
countyRENVILLE                 0.62276    0.65477   0.951 0.341844    
countyRICE                     1.21427    0.47216   2.572 0.010305 *  
countyROCK                     0.65974    0.65477   1.008 0.313971    
countyROSEAU                   0.94043    0.50719   1.854 0.064089 .  
countySCOTT                    0.90470    0.47818   1.892 0.058868 .  
countySHERBURNE                0.44947    0.48559   0.926 0.354934    
countySIBLEY                   0.59926    0.54782   1.094 0.274345    
countyST LOUIS                 0.22503    0.42028   0.535 0.592503    
countySTEARNS                  0.88653    0.44271   2.003 0.045578 *  
countySTEELE                   0.92837    0.47216   1.966 0.049632 *  
countySTEVENS                  1.13335    0.65477   1.731 0.083867 .  
countySWIFT                    0.35112    0.54782   0.641 0.521750    
countyTODD                     0.68438    0.65477   1.045 0.296250    
countyTRAVERSE                 1.42212    0.58565   2.428 0.015397 *  
countyWABASHA                  1.32879    0.50719   2.620 0.008967 ** 
countyWADENA                   0.74321    0.58565   1.269 0.204808    
countyWASECA                   0.03617    0.58565   0.062 0.950766    
countyWASHINGTON               0.72212    0.42900   1.683 0.092727 .  
countyWATONWAN                 1.95628    0.82823   2.362 0.018423 *  
countyWILKIN                   1.56510    0.82823   1.890 0.059173 .  
countyWINONA                   1.28000    0.47216   2.711 0.006858 ** 
countyWRIGHT                   0.99934    0.46299   2.158 0.031201 *  
countyYELLOW MEDICINE          0.54120    0.65477   0.827 0.408753    
floor:countyANOKA             -1.23704    0.93164  -1.328 0.184634    
floor:countyBECKER            -0.74469    1.20734  -0.617 0.537546    
floor:countyBELTRAMI          -1.08002    0.99301  -1.088 0.277101    
floor:countyBENTON            -0.40213    1.17129  -0.343 0.731448    
floor:countyBIG STONE               NA         NA      NA       NA    
floor:countyBLUE EARTH         0.18649    0.99301   0.188 0.851083    
floor:countyBROWN             -1.09370    1.09565  -0.998 0.318480    
floor:countyCARLTON            0.69127    1.12143   0.616 0.537800    
floor:countyCARVER            -1.29298    1.01437  -1.275 0.202811    
floor:countyCASS                    NA         NA      NA       NA    
floor:countyCHIPPEWA                NA         NA      NA       NA    
floor:countyCHISAGO                 NA         NA      NA       NA    
floor:countyCLAY              -1.22610    0.93061  -1.318 0.188053    
floor:countyCLEARWATER        -0.61544    1.09565  -0.562 0.574471    
floor:countyCOOK                    NA         NA      NA       NA    
floor:countyCOTTONWOOD        -1.73254    1.17129  -1.479 0.139502    
floor:countyCROW WING         -0.15523    0.95636  -0.162 0.871100    
floor:countyDAKOTA            -1.05437    0.90736  -1.162 0.245587    
floor:countyDODGE                   NA         NA      NA       NA    
floor:countyDOUGLAS           -0.53865    1.12461  -0.479 0.632101    
floor:countyFARIBAULT         -0.54985    1.14164  -0.482 0.630203    
floor:countyFILLMORE          -0.39571    1.30955  -0.302 0.762601    
floor:countyFREEBORN          -1.02266    1.00831  -1.014 0.310789    
floor:countyGOODHUE            0.13772    0.99301   0.139 0.889730    
floor:countyHENNEPIN          -0.95396    0.85919  -1.110 0.267214    
floor:countyHOUSTON           -0.67651    1.03529  -0.653 0.513660    
floor:countyHUBBARD           -0.40130    1.05579  -0.380 0.703981    
floor:countyISANTI                  NA         NA      NA       NA    
floor:countyITASCA                  NA         NA      NA       NA    
floor:countyJACKSON                 NA         NA      NA       NA    
floor:countyKANABEC                 NA         NA      NA       NA    
floor:countyKANDIYOHI               NA         NA      NA       NA    
floor:countyKITTSON           -1.34018    1.20734  -1.110 0.267333    
floor:countyKOOCHICHING        0.03946    0.99301   0.040 0.968316    
floor:countyLAC QUI PARLE     -0.50250    1.30955  -0.384 0.701288    
floor:countyLAKE              -0.24011    1.12461  -0.214 0.830991    
floor:countyLAKE OF THE WOODS -0.64027    1.09565  -0.584 0.559136    
floor:countyLE SUEUR          -0.36614    1.15285  -0.318 0.750877    
floor:countyLINCOLN           -1.31556    1.17129  -1.123 0.261714    
floor:countyLYON              -0.09672    1.12869  -0.086 0.931732    
floor:countyMAHNOMEN                NA         NA      NA       NA    
floor:countyMARSHALL          -2.50680    0.95785  -2.617 0.009041 ** 
floor:countyMARTIN            -1.92968    1.13410  -1.702 0.089249 .  
floor:countyMCLEOD            -1.16164    0.93367  -1.244 0.213818    
floor:countyMEEKER                  NA         NA      NA       NA    
floor:countyMILLE LACS        -2.42598    1.30955  -1.853 0.064329 .  
floor:countyMORRISON          -0.12666    1.12461  -0.113 0.910354    
floor:countyMOWER             -1.76618    0.99497  -1.775 0.076273 .  
floor:countyMURRAY                  NA         NA      NA       NA    
floor:countyNICOLLET                NA         NA      NA       NA    
floor:countyNOBLES                  NA         NA      NA       NA    
floor:countyNORMAN            -0.94405    1.20734  -0.782 0.434497    
floor:countyOLMSTED           -1.43773    0.93978  -1.530 0.126458    
floor:countyOTTER TAIL        -0.49377    0.97997  -0.504 0.614501    
floor:countyPENNINGTON        -1.72290    1.20734  -1.427 0.153979    
floor:countyPINE              -0.68000    1.14164  -0.596 0.551594    
floor:countyPIPESTONE         -0.45587    1.17129  -0.389 0.697237    
floor:countyPOLK              -1.58583    1.09565  -1.447 0.148192    
floor:countyPOPE                    NA         NA      NA       NA    
floor:countyRAMSEY             0.03475    0.93552   0.037 0.970375    
floor:countyREDWOOD            1.31109    1.15285   1.137 0.255780    
floor:countyRENVILLE           0.31928    1.20734   0.264 0.791504    
floor:countyRICE              -1.13089    1.11888  -1.011 0.312458    
floor:countyROCK                    NA         NA      NA       NA    
floor:countyROSEAU            -0.73058    0.91434  -0.799 0.424522    
floor:countySCOTT             -0.05889    0.93367  -0.063 0.949724    
floor:countySHERBURNE               NA         NA      NA       NA    
floor:countySIBLEY                  NA         NA      NA       NA    
floor:countyST LOUIS          -0.68398    0.85045  -0.804 0.421495    
floor:countySTEARNS           -1.12157    0.91601  -1.224 0.221173    
floor:countySTEELE                  NA         NA      NA       NA    
floor:countySTEVENS                 NA         NA      NA       NA    
floor:countySWIFT                   NA         NA      NA       NA    
floor:countyTODD               0.29100    1.20734   0.241 0.809599    
floor:countyTRAVERSE          -1.15641    1.17129  -0.987 0.323806    
floor:countyWABASHA           -1.97937    1.13410  -1.745 0.081325 .  
floor:countyWADENA            -1.09073    1.05579  -1.033 0.301883    
floor:countyWASECA            -0.86908    1.17129  -0.742 0.458322    
floor:countyWASHINGTON        -1.16249    0.89521  -1.299 0.194480    
floor:countyWATONWAN          -0.74419    1.20734  -0.616 0.537822    
floor:countyWILKIN                  NA         NA      NA       NA    
floor:countyWINONA            -2.14361    0.95336  -2.248 0.024827 *  
floor:countyWRIGHT            -0.95678    1.11504  -0.858 0.391119    
floor:countyYELLOW MEDICINE         NA         NA      NA       NA    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.7173 on 774 degrees of freedom
Multiple R-squared:  0.3539,    Adjusted R-squared:  0.2337 
F-statistic: 2.944 on 144 and 774 DF,  p-value: < 2.2e-16

Radon case study

Bayesian Linear Regression - single level `floor`

Code

library(rstanarm)
bmod1 = stan_glm(logRadon ~ floor, data = RadonData)


SAMPLING FOR MODEL 'continuous' NOW (CHAIN 1).
Chain 1: 
Chain 1: Gradient evaluation took 3.2e-05 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.32 seconds.
Chain 1: Adjust your expectations accordingly!
Chain 1: 
Chain 1: 
Chain 1: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 1: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 1: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 1: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 1: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 1: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 1: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 1: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 1: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 1: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 1: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 1: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 1: 
Chain 1:  Elapsed Time: 0.014 seconds (Warm-up)
Chain 1:                0.061 seconds (Sampling)
Chain 1:                0.075 seconds (Total)
Chain 1: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 2).
Chain 2: 
Chain 2: Gradient evaluation took 5e-06 seconds
Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.05 seconds.
Chain 2: Adjust your expectations accordingly!
Chain 2: 
Chain 2: 
Chain 2: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 2: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 2: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 2: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 2: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 2: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 2: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 2: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 2: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 2: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 2: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 2: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 2: 
Chain 2:  Elapsed Time: 0.014 seconds (Warm-up)
Chain 2:                0.059 seconds (Sampling)
Chain 2:                0.073 seconds (Total)
Chain 2: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 3).
Chain 3: 
Chain 3: Gradient evaluation took 4e-06 seconds
Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.04 seconds.
Chain 3: Adjust your expectations accordingly!
Chain 3: 
Chain 3: 
Chain 3: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 3: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 3: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 3: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 3: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 3: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 3: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 3: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 3: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 3: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 3: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 3: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 3: 
Chain 3:  Elapsed Time: 0.013 seconds (Warm-up)
Chain 3:                0.058 seconds (Sampling)
Chain 3:                0.071 seconds (Total)
Chain 3: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 4).
Chain 4: 
Chain 4: Gradient evaluation took 5e-06 seconds
Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.05 seconds.
Chain 4: Adjust your expectations accordingly!
Chain 4: 
Chain 4: 
Chain 4: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 4: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 4: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 4: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 4: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 4: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 4: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 4: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 4: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 4: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 4: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 4: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 4: 
Chain 4:  Elapsed Time: 0.014 seconds (Warm-up)
Chain 4:                0.06 seconds (Sampling)
Chain 4:                0.074 seconds (Total)
Chain 4:

Code

#bmod1 = stan_glm(logRadon~floor,data=RadonData, 
#                 prior = student_t(df = 7, 0, 5),
#                 prior_intercept = student_t(df = 7, 0, 5))

#prior_summary(bmod1)

#summary(bmod1) 


#ci95 <- posterior_interval(bmod1, prob = 0.95)
#round(ci95, 2)

#launch_shinystan(bmod1, ppd = FALSE)

Radon case study

Bayesian Linear Regression - single level `floor + county`

Code

bmod2 = stan_glm(logRadon ~ floor + county, data = RadonData)


SAMPLING FOR MODEL 'continuous' NOW (CHAIN 1).
Chain 1: 
Chain 1: Gradient evaluation took 1.8e-05 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.18 seconds.
Chain 1: Adjust your expectations accordingly!
Chain 1: 
Chain 1: 
Chain 1: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 1: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 1: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 1: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 1: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 1: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 1: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 1: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 1: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 1: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 1: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 1: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 1: 
Chain 1:  Elapsed Time: 0.218 seconds (Warm-up)
Chain 1:                0.384 seconds (Sampling)
Chain 1:                0.602 seconds (Total)
Chain 1: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 2).
Chain 2: 
Chain 2: Gradient evaluation took 1e-05 seconds
Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.1 seconds.
Chain 2: Adjust your expectations accordingly!
Chain 2: 
Chain 2: 
Chain 2: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 2: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 2: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 2: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 2: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 2: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 2: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 2: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 2: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 2: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 2: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 2: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 2: 
Chain 2:  Elapsed Time: 0.176 seconds (Warm-up)
Chain 2:                0.3 seconds (Sampling)
Chain 2:                0.476 seconds (Total)
Chain 2: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 3).
Chain 3: 
Chain 3: Gradient evaluation took 3.4e-05 seconds
Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.34 seconds.
Chain 3: Adjust your expectations accordingly!
Chain 3: 
Chain 3: 
Chain 3: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 3: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 3: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 3: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 3: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 3: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 3: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 3: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 3: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 3: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 3: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 3: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 3: 
Chain 3:  Elapsed Time: 0.205 seconds (Warm-up)
Chain 3:                0.296 seconds (Sampling)
Chain 3:                0.501 seconds (Total)
Chain 3: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 4).
Chain 4: 
Chain 4: Gradient evaluation took 1.5e-05 seconds
Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.15 seconds.
Chain 4: Adjust your expectations accordingly!
Chain 4: 
Chain 4: 
Chain 4: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 4: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 4: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 4: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 4: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 4: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 4: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 4: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 4: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 4: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 4: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 4: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 4: 
Chain 4:  Elapsed Time: 0.21 seconds (Warm-up)
Chain 4:                0.363 seconds (Sampling)
Chain 4:                0.573 seconds (Total)
Chain 4:

Code

#summary(bmod2) 

#ci95 <- posterior_interval(bmod2, prob = 0.95)
#round(ci95, 2)

#launch_shinystan(bmod2, ppd = FALSE)

Radon case study

Bayesian Linear Regression - single level `floor * county`

Code

bmod3 = stan_glm(logRadon ~ floor * county, data = RadonData)


SAMPLING FOR MODEL 'continuous' NOW (CHAIN 1).
Chain 1: 
Chain 1: Gradient evaluation took 2.9e-05 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.29 seconds.
Chain 1: Adjust your expectations accordingly!
Chain 1: 
Chain 1: 
Chain 1: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 1: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 1: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 1: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 1: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 1: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 1: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 1: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 1: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 1: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 1: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 1: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 1: 
Chain 1:  Elapsed Time: 1.193 seconds (Warm-up)
Chain 1:                1.253 seconds (Sampling)
Chain 1:                2.446 seconds (Total)
Chain 1: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 2).
Chain 2: 
Chain 2: Gradient evaluation took 2e-05 seconds
Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.2 seconds.
Chain 2: Adjust your expectations accordingly!
Chain 2: 
Chain 2: 
Chain 2: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 2: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 2: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 2: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 2: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 2: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 2: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 2: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 2: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 2: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 2: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 2: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 2: 
Chain 2:  Elapsed Time: 0.972 seconds (Warm-up)
Chain 2:                1.216 seconds (Sampling)
Chain 2:                2.188 seconds (Total)
Chain 2: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 3).
Chain 3: 
Chain 3: Gradient evaluation took 3.5e-05 seconds
Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.35 seconds.
Chain 3: Adjust your expectations accordingly!
Chain 3: 
Chain 3: 
Chain 3: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 3: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 3: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 3: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 3: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 3: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 3: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 3: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 3: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 3: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 3: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 3: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 3: 
Chain 3:  Elapsed Time: 1.119 seconds (Warm-up)
Chain 3:                1.32 seconds (Sampling)
Chain 3:                2.439 seconds (Total)
Chain 3: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 4).
Chain 4: 
Chain 4: Gradient evaluation took 3e-05 seconds
Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.3 seconds.
Chain 4: Adjust your expectations accordingly!
Chain 4: 
Chain 4: 
Chain 4: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 4: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 4: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 4: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 4: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 4: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 4: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 4: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 4: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 4: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 4: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 4: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 4: 
Chain 4:  Elapsed Time: 1.072 seconds (Warm-up)
Chain 4:                1.215 seconds (Sampling)
Chain 4:                2.287 seconds (Total)
Chain 4:

Code

#summary(bmod3) 

#launch_shinystan(bmod2, ppd = FALSE)

Radon case study

Bayesian Linear Regression - multi level

Code

bhmod1 = stan_lmer(logRadon ~ floor + (1 | county_num), data = RadonData)


SAMPLING FOR MODEL 'continuous' NOW (CHAIN 1).
Chain 1: 
Chain 1: Gradient evaluation took 4.9e-05 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.49 seconds.
Chain 1: Adjust your expectations accordingly!
Chain 1: 
Chain 1: 
Chain 1: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 1: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 1: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 1: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 1: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 1: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 1: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 1: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 1: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 1: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 1: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 1: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 1: 
Chain 1:  Elapsed Time: 1.368 seconds (Warm-up)
Chain 1:                0.46 seconds (Sampling)
Chain 1:                1.828 seconds (Total)
Chain 1: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 2).
Chain 2: 
Chain 2: Gradient evaluation took 3.5e-05 seconds
Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.35 seconds.
Chain 2: Adjust your expectations accordingly!
Chain 2: 
Chain 2: 
Chain 2: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 2: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 2: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 2: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 2: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 2: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 2: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 2: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 2: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 2: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 2: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 2: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 2: 
Chain 2:  Elapsed Time: 1.495 seconds (Warm-up)
Chain 2:                0.459 seconds (Sampling)
Chain 2:                1.954 seconds (Total)
Chain 2: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 3).
Chain 3: 
Chain 3: Gradient evaluation took 3.1e-05 seconds
Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.31 seconds.
Chain 3: Adjust your expectations accordingly!
Chain 3: 
Chain 3: 
Chain 3: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 3: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 3: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 3: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 3: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 3: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 3: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 3: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 3: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 3: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 3: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 3: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 3: 
Chain 3:  Elapsed Time: 1.067 seconds (Warm-up)
Chain 3:                0.523 seconds (Sampling)
Chain 3:                1.59 seconds (Total)
Chain 3: 

SAMPLING FOR MODEL 'continuous' NOW (CHAIN 4).
Chain 4: 
Chain 4: Gradient evaluation took 3.8e-05 seconds
Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.38 seconds.
Chain 4: Adjust your expectations accordingly!
Chain 4: 
Chain 4: 
Chain 4: Iteration:    1 / 2000 [  0%]  (Warmup)
Chain 4: Iteration:  200 / 2000 [ 10%]  (Warmup)
Chain 4: Iteration:  400 / 2000 [ 20%]  (Warmup)
Chain 4: Iteration:  600 / 2000 [ 30%]  (Warmup)
Chain 4: Iteration:  800 / 2000 [ 40%]  (Warmup)
Chain 4: Iteration: 1000 / 2000 [ 50%]  (Warmup)
Chain 4: Iteration: 1001 / 2000 [ 50%]  (Sampling)
Chain 4: Iteration: 1200 / 2000 [ 60%]  (Sampling)
Chain 4: Iteration: 1400 / 2000 [ 70%]  (Sampling)
Chain 4: Iteration: 1600 / 2000 [ 80%]  (Sampling)
Chain 4: Iteration: 1800 / 2000 [ 90%]  (Sampling)
Chain 4: Iteration: 2000 / 2000 [100%]  (Sampling)
Chain 4: 
Chain 4:  Elapsed Time: 1.216 seconds (Warm-up)
Chain 4:                0.472 seconds (Sampling)
Chain 4:                1.688 seconds (Total)
Chain 4:

Radon case study

Bayesian Linear Regression - multi level

Code

summary(bhmod1)


Model Info:
 function:     stan_lmer
 family:       gaussian [identity]
 formula:      logRadon ~ 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)                                1.5    0.1  1.4   1.5   1.6 
floor                                     -0.7    0.1 -0.8  -0.7  -0.6 
b[(Intercept) county_num:1]               -0.3    0.2 -0.6  -0.3   0.0 
b[(Intercept) county_num:2]               -0.5    0.1 -0.6  -0.5  -0.4 
b[(Intercept) county_num:3]                0.0    0.3 -0.3   0.0   0.3 
b[(Intercept) county_num:4]                0.0    0.2 -0.2   0.0   0.3 
b[(Intercept) county_num:5]                0.0    0.2 -0.3   0.0   0.3 
b[(Intercept) county_num:6]                0.0    0.3 -0.3   0.0   0.3 
b[(Intercept) county_num:7]                0.4    0.2  0.2   0.4   0.6 
b[(Intercept) county_num:8]                0.2    0.2 -0.1   0.2   0.5 
b[(Intercept) county_num:9]               -0.3    0.2 -0.5  -0.3   0.0 
b[(Intercept) county_num:10]               0.0    0.2 -0.2   0.0   0.3 
b[(Intercept) county_num:11]               0.0    0.2 -0.3   0.0   0.3 
b[(Intercept) county_num:12]               0.1    0.2 -0.2   0.1   0.4 
b[(Intercept) county_num:13]              -0.2    0.2 -0.5  -0.2   0.1 
b[(Intercept) county_num:14]               0.4    0.2  0.1   0.4   0.6 
b[(Intercept) county_num:15]              -0.1    0.2 -0.4  -0.1   0.2 
b[(Intercept) county_num:16]              -0.2    0.3 -0.6  -0.2   0.1 
b[(Intercept) county_num:17]              -0.1    0.2 -0.4  -0.1   0.2 
b[(Intercept) county_num:18]              -0.2    0.2 -0.5  -0.2   0.0 
b[(Intercept) county_num:19]              -0.1    0.1 -0.2  -0.1   0.0 
b[(Intercept) county_num:20]               0.1    0.3 -0.2   0.1   0.5 
b[(Intercept) county_num:21]               0.2    0.2 -0.1   0.2   0.4 
b[(Intercept) county_num:22]              -0.4    0.2 -0.7  -0.4  -0.1 
b[(Intercept) county_num:23]               0.0    0.3 -0.4   0.0   0.3 
b[(Intercept) county_num:24]               0.4    0.2  0.1   0.4   0.6 
b[(Intercept) county_num:25]               0.3    0.2  0.1   0.3   0.5 
b[(Intercept) county_num:26]              -0.1    0.1 -0.2  -0.1   0.0 
b[(Intercept) county_num:27]               0.2    0.2 -0.1   0.1   0.4 
b[(Intercept) county_num:28]              -0.1    0.2 -0.4  -0.1   0.2 
b[(Intercept) county_num:29]              -0.1    0.3 -0.5  -0.1   0.2 
b[(Intercept) county_num:30]              -0.4    0.2 -0.6  -0.4  -0.1 
b[(Intercept) county_num:31]               0.3    0.2  0.0   0.3   0.6 
b[(Intercept) county_num:32]              -0.1    0.2 -0.4  -0.1   0.2 
b[(Intercept) county_num:33]               0.3    0.2 -0.1   0.3   0.6 
b[(Intercept) county_num:34]               0.0    0.3 -0.3   0.0   0.4 
b[(Intercept) county_num:35]              -0.4    0.2 -0.6  -0.4  -0.1 
b[(Intercept) county_num:36]               0.4    0.3  0.0   0.4   0.8 
b[(Intercept) county_num:37]              -0.6    0.2 -0.9  -0.6  -0.4 
b[(Intercept) county_num:38]               0.2    0.2 -0.1   0.2   0.5 
b[(Intercept) county_num:39]               0.1    0.2 -0.2   0.1   0.4 
b[(Intercept) county_num:40]               0.4    0.2  0.0   0.4   0.7 
b[(Intercept) county_num:41]               0.3    0.2  0.0   0.3   0.5 
b[(Intercept) county_num:42]               0.0    0.3 -0.4   0.0   0.3 
b[(Intercept) county_num:43]               0.1    0.2 -0.2   0.1   0.3 
b[(Intercept) county_num:44]              -0.2    0.2 -0.5  -0.2   0.0 
b[(Intercept) county_num:45]              -0.1    0.2 -0.4  -0.1   0.1 
b[(Intercept) county_num:46]              -0.1    0.2 -0.4  -0.1   0.2 
b[(Intercept) county_num:47]              -0.1    0.3 -0.5  -0.2   0.2 
b[(Intercept) county_num:48]              -0.2    0.2 -0.4  -0.2   0.1 
b[(Intercept) county_num:49]               0.2    0.2 -0.1   0.2   0.4 
b[(Intercept) county_num:50]               0.2    0.3 -0.2   0.2   0.6 
b[(Intercept) county_num:51]               0.3    0.3  0.0   0.3   0.6 
b[(Intercept) county_num:52]               0.2    0.3 -0.2   0.2   0.5 
b[(Intercept) county_num:53]              -0.1    0.3 -0.4  -0.1   0.2 
b[(Intercept) county_num:54]              -0.1    0.1 -0.3  -0.1   0.1 
b[(Intercept) county_num:55]               0.1    0.2 -0.2   0.1   0.3 
b[(Intercept) county_num:56]              -0.1    0.2 -0.4  -0.1   0.2 
b[(Intercept) county_num:57]              -0.4    0.2 -0.7  -0.4  -0.1 
b[(Intercept) county_num:58]               0.2    0.2 -0.1   0.2   0.5 
b[(Intercept) county_num:59]               0.1    0.2 -0.2   0.1   0.4 
b[(Intercept) county_num:60]              -0.1    0.3 -0.4  -0.1   0.3 
b[(Intercept) county_num:61]              -0.3    0.1 -0.4  -0.3  -0.1 
b[(Intercept) county_num:62]               0.2    0.2 -0.1   0.2   0.5 
b[(Intercept) county_num:63]               0.1    0.3 -0.2   0.1   0.4 
b[(Intercept) county_num:64]               0.2    0.2  0.0   0.2   0.5 
b[(Intercept) county_num:65]               0.0    0.3 -0.4   0.0   0.3 
b[(Intercept) county_num:66]               0.1    0.2 -0.1   0.1   0.3 
b[(Intercept) county_num:67]               0.2    0.2  0.0   0.2   0.4 
b[(Intercept) county_num:68]              -0.2    0.2 -0.5  -0.2   0.0 
b[(Intercept) county_num:69]              -0.1    0.2 -0.4  -0.1   0.2 
b[(Intercept) county_num:70]              -0.6    0.1 -0.7  -0.5  -0.4 
b[(Intercept) county_num:71]               0.0    0.1 -0.2   0.0   0.2 
b[(Intercept) county_num:72]               0.1    0.2 -0.2   0.1   0.3 
b[(Intercept) county_num:73]               0.1    0.3 -0.3   0.1   0.4 
b[(Intercept) county_num:74]              -0.2    0.2 -0.5  -0.2   0.1 
b[(Intercept) county_num:75]               0.1    0.3 -0.2   0.1   0.4 
b[(Intercept) county_num:76]               0.2    0.2 -0.1   0.2   0.5 
b[(Intercept) county_num:77]               0.2    0.2 -0.1   0.2   0.5 
b[(Intercept) county_num:78]              -0.1    0.2 -0.4  -0.1   0.2 
b[(Intercept) county_num:79]              -0.3    0.2 -0.7  -0.3   0.0 
b[(Intercept) county_num:80]              -0.1    0.1 -0.3  -0.1   0.0 
b[(Intercept) county_num:81]               0.4    0.3  0.1   0.4   0.8 
b[(Intercept) county_num:82]               0.1    0.3 -0.2   0.1   0.5 
b[(Intercept) county_num:83]               0.1    0.2 -0.1   0.1   0.3 
b[(Intercept) county_num:84]               0.1    0.2 -0.1   0.1   0.3 
b[(Intercept) county_num:85]              -0.1    0.3 -0.4  -0.1   0.3 
sigma                                      0.7    0.0  0.7   0.7   0.7 
Sigma[county_num:(Intercept),(Intercept)]  0.1    0.0  0.1   0.1   0.1 

Fit Diagnostics:
           mean   sd   10%   50%   90%
mean_PPD 1.3    0.0  1.2   1.3   1.3  

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  1643 
floor                                     0.0  1.0  4979 
b[(Intercept) county_num:1]               0.0  1.0  5468 
b[(Intercept) county_num:2]               0.0  1.0  4369 
b[(Intercept) county_num:3]               0.0  1.0  5655 
b[(Intercept) county_num:4]               0.0  1.0  5608 
b[(Intercept) county_num:5]               0.0  1.0  6392 
b[(Intercept) county_num:6]               0.0  1.0  4751 
b[(Intercept) county_num:7]               0.0  1.0  4742 
b[(Intercept) county_num:8]               0.0  1.0  5596 
b[(Intercept) county_num:9]               0.0  1.0  5479 
b[(Intercept) county_num:10]              0.0  1.0  4872 
b[(Intercept) county_num:11]              0.0  1.0  6262 
b[(Intercept) county_num:12]              0.0  1.0  6035 
b[(Intercept) county_num:13]              0.0  1.0  5341 
b[(Intercept) county_num:14]              0.0  1.0  4920 
b[(Intercept) county_num:15]              0.0  1.0  5345 
b[(Intercept) county_num:16]              0.0  1.0  5090 
b[(Intercept) county_num:17]              0.0  1.0  6303 
b[(Intercept) county_num:18]              0.0  1.0  5688 
b[(Intercept) county_num:19]              0.0  1.0  3467 
b[(Intercept) county_num:20]              0.0  1.0  6435 
b[(Intercept) county_num:21]              0.0  1.0  4806 
b[(Intercept) county_num:22]              0.0  1.0  4556 
b[(Intercept) county_num:23]              0.0  1.0  5484 
b[(Intercept) county_num:24]              0.0  1.0  5352 
b[(Intercept) county_num:25]              0.0  1.0  4955 
b[(Intercept) county_num:26]              0.0  1.0  2618 
b[(Intercept) county_num:27]              0.0  1.0  4513 
b[(Intercept) county_num:28]              0.0  1.0  5315 
b[(Intercept) county_num:29]              0.0  1.0  5743 
b[(Intercept) county_num:30]              0.0  1.0  4921 
b[(Intercept) county_num:31]              0.0  1.0  4713 
b[(Intercept) county_num:32]              0.0  1.0  5686 
b[(Intercept) county_num:33]              0.0  1.0  4778 
b[(Intercept) county_num:34]              0.0  1.0  5356 
b[(Intercept) county_num:35]              0.0  1.0  4862 
b[(Intercept) county_num:36]              0.0  1.0  4608 
b[(Intercept) county_num:37]              0.0  1.0  4129 
b[(Intercept) county_num:38]              0.0  1.0  5392 
b[(Intercept) county_num:39]              0.0  1.0  5563 
b[(Intercept) county_num:40]              0.0  1.0  5458 
b[(Intercept) county_num:41]              0.0  1.0  4918 
b[(Intercept) county_num:42]              0.0  1.0  5958 
b[(Intercept) county_num:43]              0.0  1.0  5536 
b[(Intercept) county_num:44]              0.0  1.0  5127 
b[(Intercept) county_num:45]              0.0  1.0  5219 
b[(Intercept) county_num:46]              0.0  1.0  5832 
b[(Intercept) county_num:47]              0.0  1.0  5365 
b[(Intercept) county_num:48]              0.0  1.0  6184 
b[(Intercept) county_num:49]              0.0  1.0  4931 
b[(Intercept) county_num:50]              0.0  1.0  4025 
b[(Intercept) county_num:51]              0.0  1.0  4997 
b[(Intercept) county_num:52]              0.0  1.0  5085 
b[(Intercept) county_num:53]              0.0  1.0  6263 
b[(Intercept) county_num:54]              0.0  1.0  4831 
b[(Intercept) county_num:55]              0.0  1.0  5527 
b[(Intercept) county_num:56]              0.0  1.0  5498 
b[(Intercept) county_num:57]              0.0  1.0  4751 
b[(Intercept) county_num:58]              0.0  1.0  5691 
b[(Intercept) county_num:59]              0.0  1.0  5307 
b[(Intercept) county_num:60]              0.0  1.0  5311 
b[(Intercept) county_num:61]              0.0  1.0  4500 
b[(Intercept) county_num:62]              0.0  1.0  5656 
b[(Intercept) county_num:63]              0.0  1.0  6090 
b[(Intercept) county_num:64]              0.0  1.0  4792 
b[(Intercept) county_num:65]              0.0  1.0  5061 
b[(Intercept) county_num:66]              0.0  1.0  4385 
b[(Intercept) county_num:67]              0.0  1.0  5397 
b[(Intercept) county_num:68]              0.0  1.0  5155 
b[(Intercept) county_num:69]              0.0  1.0  6470 
b[(Intercept) county_num:70]              0.0  1.0  2854 
b[(Intercept) county_num:71]              0.0  1.0  4103 
b[(Intercept) county_num:72]              0.0  1.0  5844 
b[(Intercept) county_num:73]              0.0  1.0  6013 
b[(Intercept) county_num:74]              0.0  1.0  5465 
b[(Intercept) county_num:75]              0.0  1.0  6543 
b[(Intercept) county_num:76]              0.0  1.0  4655 
b[(Intercept) county_num:77]              0.0  1.0  5333 
b[(Intercept) county_num:78]              0.0  1.0  5039 
b[(Intercept) county_num:79]              0.0  1.0  4979 
b[(Intercept) county_num:80]              0.0  1.0  4215 
b[(Intercept) county_num:81]              0.0  1.0  4593 
b[(Intercept) county_num:82]              0.0  1.0  5350 
b[(Intercept) county_num:83]              0.0  1.0  5473 
b[(Intercept) county_num:84]              0.0  1.0  4971 
b[(Intercept) county_num:85]              0.0  1.0  5095 
sigma                                     0.0  1.0  4471 
Sigma[county_num:(Intercept),(Intercept)] 0.0  1.0  1242 
mean_PPD                                  0.0  1.0  4171 
log-posterior                             0.3  1.0   862 

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).

Code

#prior_summary(bhmod1)

#launch_shinystan(bhmod1, ppd = FALSE)

ST308 Bayesian Inference Class 9

Recap

Normal linear model

Frequentist

Bayesian

Normal linear model

Interactions

Hierarchichal model

Hierarchichal model

Frequentist approach:

Bayesian approach:

Case study

stan_lmer

stan_lmer

Key inputs:

Radon case study

EDA

Radon case study

Floor variable

Radon case study

County variable

Radon case study

Linear model floor

Radon case study

Linear model floor + county

Radon case study

Linear model floor * county

Radon case study

Bayesian Linear Regression - single level floor

Radon case study

Bayesian Linear Regression - single level floor + county

Radon case study

Bayesian Linear Regression - single level floor * county

Radon case study

Bayesian Linear Regression - multi level

Radon case study

Bayesian Linear Regression - multi level

ST308 Bayesian Inference

Class 9

`stan_lmer`

`stan_lmer`

Linear model `floor`

Linear model `floor + county`

Linear model `floor * county`

Bayesian Linear Regression - single level `floor`

Bayesian Linear Regression - single level `floor + county`

Bayesian Linear Regression - single level `floor * county`