1 Bayes’ Theorem, for Dummies

Imagine you have a belief about something (a coin’s fairness, a regression coefficient, whatever). Then you observe data. Bayes’ theorem tells you exactly how to update your belief given that new evidence.

\[ \text{Posterior} \;\propto\; \text{Likelihood} \times \text{Prior} \]

Formally:

\[ P(\theta \mid \text{data}) = \frac{P(\text{data} \mid \theta) \, P(\theta)}{P(\text{data})} \]

Analogy: A doctor thinks you have a 1% chance of a rare disease (prior). You test positive on a test that’s 90% accurate (likelihood). Bayes’ theorem combines those two things to tell you your actual updated probability of having the disease (posterior) — which is often much lower than 90%, because the prior matters a lot.

2 Bayesian vs. OLS — the Core Difference

OLS (Frequentist) Bayesian
Parameters Fixed, unknown constants Random variables with a probability distribution
Goal Find the single point estimate minimizing squared residuals Compute a full posterior distribution over parameters
Prior info Not used (unless you add penalties like ridge/lasso) Explicitly incorporated via a prior
Output Point estimate + confidence interval Full posterior + credible interval
Interval interpretation “95% of intervals built this way would contain the true value” “There’s a 95% probability the true value is in this interval”
With a flat/vague prior — Posterior mean \(\approx\) OLS estimate (they converge!)
With an informative prior — Estimates shrink toward the prior (regularization effect, like ridge regression)

Key insight: OLS is basically a special case of Bayesian inference — the case where you use a completely flat (uninformative) prior. Bayesian analysis generalizes it by letting you add real prior knowledge, and gives you a full distribution of plausible parameter values instead of a single number.

3 Example 1 — Beta-Binomial (Coin Flip): the Simplest Bayes Update

## [1] 0.7083333
## [1] 0.75

Notice the Bayesian estimate is pulled slightly toward the prior belief of 0.5, compared to the raw frequentist estimate. That “pull toward the prior” is the entire philosophical difference from OLS.

4 Example 2 — Bayesian Linear Regression vs. OLS

This example shows mathematically how a Bayesian posterior mean converges to the OLS estimate as the prior becomes vague, and diverges (shrinks) when the prior is informative. It uses only base R matrix algebra — no extra packages required.

4.1 OLS Estimate

##       [,1]
##   1.687471
## x 3.027159
## (Intercept)           x 
##    1.687471    3.027159

4.2 Bayesian Estimate

We assume a prior \(\beta \sim \text{Normal}(0, \tau^2 I)\) and treat \(\sigma^2\) as known for simplicity. The posterior is then also Normal, with a closed-form mean:

##       [,1]
##   1.687470
## x 3.027159
##       [,1]
##   1.270270
## x 3.089465

What this shows:

  • With a vague/weak prior, Bayesian and OLS estimates are essentially the same.
  • With a strong/informative prior, the Bayesian estimate is “pulled” toward the prior — this is mathematically identical to ridge regression (L2 regularization). Ridge regression is literally Bayesian linear regression with a \(\text{Normal}(0, \tau^2)\) prior on the coefficients.

4.3 Side-by-Side Comparison Table

OLS vs Bayesian estimates under different priors
Method Intercept Slope
OLS 1.687 3.027
Bayes (vague prior) 1.687 3.027
Bayes (strong prior) 1.270 3.089

5 A More Realistic Bayesian Workflow (MCMC via rstanarm)

In practice you rarely do the matrix algebra by hand — you use a package that runs MCMC sampling to get the full posterior distribution, which also works for non-conjugate models. This chunk is set to eval = FALSE by default since it requires installing rstanarm (which depends on Stan and can take a while to install).

## 
## SAMPLING FOR MODEL 'continuous' NOW (CHAIN 1).
## Chain 1: 
## Chain 1: Gradient evaluation took 7.4e-05 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.74 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.062 seconds (Warm-up)
## Chain 1:                0.048 seconds (Sampling)
## Chain 1:                0.11 seconds (Total)
## Chain 1: 
## 
## SAMPLING FOR MODEL 'continuous' NOW (CHAIN 2).
## Chain 2: 
## Chain 2: Gradient evaluation took 1.8e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.18 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.04 seconds (Warm-up)
## Chain 2:                0.054 seconds (Sampling)
## Chain 2:                0.094 seconds (Total)
## Chain 2: 
## 
## SAMPLING FOR MODEL 'continuous' NOW (CHAIN 3).
## Chain 3: 
## Chain 3: Gradient evaluation took 1.6e-05 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.16 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.042 seconds (Warm-up)
## Chain 3:                0.045 seconds (Sampling)
## Chain 3:                0.087 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.037 seconds (Warm-up)
## Chain 4:                0.042 seconds (Sampling)
## Chain 4:                0.079 seconds (Total)
## Chain 4:
## stan_glm
##  family:       gaussian [identity]
##  formula:      y ~ x
##  observations: 100
##  predictors:   2
## ------
##             Median MAD_SD
## (Intercept) 1.7    0.5   
## x           3.0    0.1   
## 
## Auxiliary parameter(s):
##       Median MAD_SD
## sigma 1.8    0.1   
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg
##                  2.5%    97.5%
## (Intercept) 0.7448792 2.637528
## x           2.8529457 3.204966
## sigma       1.6022757 2.115897
## 
## Call:
## lm(formula = y ~ x, data = df)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.7768 -1.0133  0.0245  1.0821  5.7248 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.68747    0.47975   3.517 0.000662 ***
## x            3.02716    0.08767  34.531  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.817 on 98 degrees of freedom
## Multiple R-squared:  0.9241, Adjusted R-squared:  0.9233 
## F-statistic:  1192 on 1 and 98 DF,  p-value: < 2.2e-16

stan_glm() gives you the full posterior distribution for each coefficient (so you can plot it, compute any probability statement you want, e.g. “P(slope > 2.5)”), whereas lm() only gives you a point estimate and a confidence interval based on repeated-sampling theory.

6 TL;DR