I got some questions about the Bayesian regression:

  1. Given a standard regression as $y = \beta_0 + \beta_1 x + \varepsilon$. If I want to change this into a Bayesian regression, do I need prior distributions both for $\beta_0$ and $\beta_1$ (or doesn't it work this way)?

  2. In standard regression one would try to minimize the residuals to get single values for $\beta_0$ and $\beta_1$. How is this done in Bayes regression?

I really struggle a lot here:

$$ \text{posterior} = \text{prior} \times \text{likelihood} $$

Likelihood comes from the current dataset (so it's my regression parameter but not as a single value but as a likelihood distribution, right?). Prior comes from a previous research (let's say). So I got this equation:

$$ y = \beta_1 x + \varepsilon $$

with $\beta_1$ being my likelihood or posterior (or is this just totally wrong)?

I simply can't understand how the standard regression transforms into a Bayes one.


The simple linear regression model

$$ y_i = \alpha + \beta x_i + \varepsilon $$

can be written in terms of the probabilistic model behind it

$$ \mu_i = \alpha + \beta x_i \\ y_i \sim \mathcal{N}(\mu_i, \sigma) $$

i.e. dependent variable $Y$ follows normal distribution parametrized by mean $\mu_i$, that is a linear function of $X$ parametrized by $\alpha,\beta$, and by standard deviation $\sigma$. If you estimate such a model using ordinary least squares, you do not have to bother about the probabilistic formulation, because you are searching for optimal values of $\alpha,\beta$ parameters by minimizing the squared errors of fitted values to predicted values. On another hand, you could estimate such model using maximum likelihood estimation, where you would be looking for optimal values of parameters by maximizing the likelihood function

$$ \DeclareMathOperator*{\argmax}{arg\,max} \argmax_{\alpha,\,\beta,\,\sigma} \prod_{i=1}^n \mathcal{N}(y_i; \alpha + \beta x_i, \sigma) $$

where $\mathcal{N}$ is a density function of normal distribution evaluated at $y_i$ points, parametrized by means $\alpha + \beta x_i$ and standard deviation $\sigma$.

In the Bayesian approach instead of maximizing the likelihood function alone, we would assume prior distributions for the parameters and use the Bayes theorem

$$ \text{posterior} \propto \text{likelihood} \times \text{prior} $$

The likelihood function is the same as above, but what changes is that you assume some prior distributions for the estimated parameters $\alpha,\beta,\sigma$ and include them into the equation

$$ \underbrace{f(\alpha,\beta,\sigma\mid Y,X)}_{\text{posterior}} \propto \underbrace{\prod_{i=1}^n \mathcal{N}(y_i\mid \alpha + \beta x_i, \sigma)}_{\text{likelihood}} \; \underbrace{f_{\alpha}(\alpha) \, f_{\beta}(\beta) \, f_{\sigma}(\sigma)}_{\text{priors}} $$

"What distributions?" is a different question, since there is an unlimited number of choices. For $\alpha,\beta$ parameters you could, for example, assume normal distributions parametrized by some hyperparameters, or $t$-distribution if you want to assume heavier tails, or uniform distribution if you do not want to make many assumptions, but you want to assume that the parameters can be a priori "anything in the given range", etc. For $\sigma$ you need to assume some prior distribution that is bounded to be greater than zero since standard deviation needs to be positive. This may lead to the model formulation as illustrated below by John K. Kruschke.

Bayesian linear regression model formulation. The beta parameters all have normal priors, sigma has uniform prior, and normal likelihood is used.

(source: http://www.indiana.edu/~kruschke/BMLR/)

While in the maximum likelihood you were looking for a single optimal value for each of the parameters, in the Bayesian approach by applying the Bayes theorem you obtain the posterior distribution of the parameters. The final estimate will depend on the information that comes from your data and from your priors, but the more information is contained in your data, the less influential are priors.

Notice that when using uniform priors, they take form $f(\theta) \propto 1$ after dropping the normalizing constants. This makes Bayes theorem proportional to the likelihood function alone, so the posterior distribution will reach its maximum at exactly the same point as the maximum likelihood estimate. What follows, the estimate under uniform priors will be the same as by using ordinary least squares since minimizing the squared errors corresponds to maximizing the normal likelihood.

To estimate a model in the Bayesian approach in some cases you can use conjugate priors, so the posterior distribution is directly available (see example here). However, in the vast majority of cases, posterior distribution will not be directly available and you will have to use Markov Chain Monte Carlo methods for estimating the model (check this example of using Metropolis-Hastings algorithm to estimate parameters of linear regression). Finally, if you are only interested in point estimates of parameters, you could use maximum a posteriori estimation, i.e.

$$ \argmax_{\alpha,\,\beta,\,\sigma} f(\alpha,\beta,\sigma\mid Y,X) $$

For a more detailed description of logistic regression, you can check the Bayesian logit model - intuitive explanation? thread.

For learning more you could check the following books:

Kruschke, J. (2014). Doing Bayesian Data Analysis: A Tutorial with R, JAGS, and Stan. Academic Press.

Gelman, A., Carlin, J. B., Stern, H. S., and Rubin, D. B. (2004). Bayesian data analysis. Chapman & Hall/CRC.

  • 6
    $\begingroup$ +1 Given the way the question is stated, I'd maybe emphasize a bit more this philosophical difference: In ordinary least squares and maximum likelihood estimation, we are starting with the question "What are the best values for $\beta_i$ (perhaps for later use)?", whereas in the full Bayesian approach, we start with the question "What can we say about the unknown values $\beta_i$?" and then maybe proceed to using the maximum a posteriori or posterior mean if a point estimate is needed. $\endgroup$
    – JiK
    Dec 22 '16 at 11:32
  • 3
    $\begingroup$ +1. One more thing that might be useful to point out to clarify the relationship between Bayesian and OLS approaches is that OLS can be understood as posterior mean under a flat prior (at least as far as I understand). Would be great if you could elaborate on that a bit in your answer. $\endgroup$
    – amoeba
    Dec 22 '16 at 15:07
  • $\begingroup$ @amoeba it's a good point, I'll think about it. But on another hand, I don't want to make the answer overtly long, so there is a point in going to details. $\endgroup$
    – Tim
    Dec 22 '16 at 16:08
  • 1
    $\begingroup$ @amoeba FYI, I added a brief comment on that. $\endgroup$
    – Tim
    Dec 28 '16 at 11:52

Given a data set $D = (x_1,y_1), \ldots, (x_N,y_N)$ where $x \in \mathbb{R}^d, y \in \mathbb{R}$, a Bayesian Linear Regression models the problem in the following way:

Prior: $$w \sim \mathcal{N}(0, \sigma_w^2 I_d)$$

$w$ is vector $(w_1, \ldots, w_d)^T$, so the previous distribution is a multivariate Gaussian; and $I_d$ is the $d\times d$ identity matrix.

Likelihood: $$Y_i \sim \mathcal{N}(w^T x_i, \sigma^2)$$

We assume that $Y_i \perp Y_j | w, i \neq j$

For now we'll use the precision instead of the variance, $a = 1/\sigma^2$, and $b = 1/\sigma_w^2$. We'll also assume that $a,b$ are known.

The prior can be stated as $$p(w) \propto \exp \Big\{ -\frac{b}{2} w^t w \Big\}$$

And the likelihood $$p(D|w) \propto \exp \Big\{ -\frac{a}{2} (y-Aw)^T (y-Aw) \Big\}$$

where $y = (y_1,\ldots,y_N)^T$ and $A$ is a $n\times d$ matrix where the i-th row is $x_i^T$.

Then the posterior is $$p(w|D) \propto p(D|w) p(w)$$

After many calculations we discover that

$$p(w|D) \sim \mathcal{N}(w | \mu, \Lambda^{-1})$$

where ($\Lambda$ is the precision matrix)

$$\Lambda = a A^T A + b I_d $$ $$\mu = a \Lambda^{-1} A^T y$$

Notice that $\mu$ is equal to the $w_{MAP}$ of the regular linear regression, this is because for the Gaussian, the mean is equal to the mode.

Also, we can make some algebra over $\mu$ and get the following equality ($\Lambda = aA^TA+bI_d$):

$$\mu = (A^T A + \frac{b}{a} I_d)^{-1} A^T y$$

and compare with $w_{MLE}$:

$$w_{MLE} = (A^T A)^{-1} A^T y$$

The extra expression in $\mu$ corresponds to the prior. This is similar to the expression for the Ridge regression, for the special case when $\lambda = \frac{b}{a}$. Ridge regression is more general because the technique can choose improper priors (in the Bayesian perspective).

For the predictive posterior distribution:

$$p(y|x,D) = \int p(y|x,D,w) p(w|x,D) dw = \int p(y|x,w) p(w|D) dw$$

it is possible to calculate that

$$y|x,D \sim \mathcal{N}(\mu^Tx, \frac{1}{a} + x^T \Lambda^{-1}x)$$

Reference: Lunn et al. The BUGS Book

For using a MCMC tool like JAGS/Stan check Kruschke's Doing Bayesian Data Analysis

  • $\begingroup$ Thank you jpneto. i feel that this is a great answer, but i dont understand it yet because of a lack in math-knowledge. But i will definitely read it again after gaining some math-skills $\endgroup$ Dec 21 '16 at 18:53
  • 1
    $\begingroup$ This is very nice, but the assumption that the precision is known is a bit uncommon. Isn't it much more common to assume an inverse Gamma distribution for the variance, i.e., a Gamma distribution for the precision? $\endgroup$
    – DeltaIV
    Dec 22 '16 at 10:57
  • $\begingroup$ +1. Can you comment a bit more on "Ridge regression is more general because the technique can choose improper priors"? I don't get it. I thought RR = Gaussian (proper) prior on $w$. $\endgroup$
    – amoeba
    Dec 23 '16 at 8:05
  • $\begingroup$ @amoeba: The Gaussian prior is $w \sim N(0,\lambda^{-1} I_d)$ but $\lambda$ can be zero which results on an improper prior, ie, it results in the MLE. $\endgroup$
    – jpneto
    Dec 23 '16 at 9:25
  • 1
    $\begingroup$ @DeltaIV: sure, when we have uncertainty about a parameter we can model that with a prior. The assumption of known precision is to make it easier to find an analytic solution. Usually, those analytic solutions are not possible and we must use approximations, like MCMC or some variational technique. $\endgroup$
    – jpneto
    Dec 23 '16 at 9:29

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.