MCMC: linear model inference

Bayesian inference using MCMC

Why is MCMC so useful for Bayesian inference?

Bayes’ Theorem tells us how to use data y to change a prior into a posterior: p(\(\theta\)|y) = p(\(\theta\)) p(y|\(\theta\)) / p(y). There is usually no problem specifying a prior p(\(\theta\)) and a likelihood function L(\(\theta\)) = p(y|\(\theta\)). But we usually do not know what p(y) is. So we cannot write down the formula for the posterior. But we can still use MCMC because the troublesome p(y) terms cancel out in the ratio p(\(\theta\)a|y) / p(\(\theta\)b|y) for any pair of points \(\theta\)a and \(\theta\)b in parameter space. The ratio of posterior probabilities for \(\theta\)a and \(\theta\)b is thus the ratio of “prior times likelihood” for \(\theta\)a and \(\theta\)b:

\[p(\theta a|y) / p(\theta b|y) = p(\theta a) p(y|\theta a) / (p(\theta b) p(y|\theta b))\]

Because that posterior ratio is usually known for any pair of points, MCMC can be used to sample from posterior distributions.

Bayesian linear regression

We begin by getting our data, defining our model, prior and likelihood function. We then run the MCMC inference and analyse the outputs.

You might like to experiment changing the:

  1. Data uncertainty (standard deviation)
  2. Prior ranges
  3. MCMC chain starting point
  4. MCMC chain length
  5. Size of the proposal step

Rerun the MCMC observe any changes. This will help increase your understanding of how MCMC works.

data (y)

There are three data points. The first and second columns are “x” and “y” respectively. The third column sets the uncertainty in the y-data as a standard deviation which will be used in the log likelihood function.

linear model

The model has two parameters slope and intercept.

prior \((\theta)\)

Both parameters have a prior log uniform distribution (see log( dunif(…)) below). For a uniform distribution we set maximum and minimum ranges for each parameter.

log likelihood function

We choose a log Gaussian distribution for the likelihood.

MCMC chain length and start

chainLength is the size of the MCMC chain. pStart is the MCMC chain starting point.

run MCMC (Metropolis)

The code below runs the MCMC inference using the Metropolis algorithm. The size of the proposal step length is given by stepLength.

posterior summary statistics

Calculate and print out some summary statistics of the MCMC chain.

plot the MCMC parameter chain

model predictions

Remove the burn in and rerun the model. Plot the mean of the posterior model predictions (black) against the model prediction for the MAP value (red). The 5th and 95th quantiles of the posterior predictions are shown by the grey shaded area.