Bayesian Linear Regression
The debate over statistical legitimacy between frequentists and bayesianist has raged for ages. But when it comes to simple machine learning approaches, like regression, there is a large amount of dominance by frequentists’ method to sample and analyze data based on those assumptions. As a fledgling bayesianist, I explored if linear regression is possible using the Bayesian approach and I was not disappointed. Before delving into the mathematics of the Bayesian method to address regression problems, it is important to do a quick recap on the frequentist approach to linear regression.
Frequentists’ Approach to Linear Regression
Linear regressions are often presented as well drawn problems from a uniformly distributed population. This assumptions conforms well with the central limit theorem. In brief, the classical representation of linear regression looks something like this:
Where ε and y terms are sampled from a uniform distribution:
The terms X and y represent the feature matrix and the target that we want to estimate with our model. Classical regression uses Ordinary Least Squares (OLS) where a single fit is calculated that minimizes the L2 loss function (also called residual sum) represented as follows.
Finding the β estimates (β hat) values that can minimize the loss could be computed using the following equation.
Inverting a matrix is computationally costly when the sizes are large. Some machine learning tools take different route in doing this. Scikit learn uses iterative gradient descent approach while statsmodel employs a direct method, techniques probably like Moore–Penrose inverse that avoid matrix singularity.
The flaw of this method is that it does not account for priors and updating the posterior, unlike the Bayesian statistics. In addition, in the classical method, only one fit is given as a solution to predict future cases which is unrealistic in the real world. The other shortcoming is that the β terms predicted using OLS have high variance, meaning that small changes in x can lead to high change in the coefficients, another undesirable effect.
The Bayesian Approach
Unironically, the Bayesian approach uses Bayes’ theorem while also considering the central limit theorem. The following picture shades light on what terminologies are used in describing Bayesian statistics.
When used for linear regression, the equation looks like this.
MAP stands for Maximum a posteriori estimation, which is to say maximizing the posterior probability. The direct interpretation is that what are the β estimates that can maximize the prediction, given features and target values. The conditional probability can be replaced its Bayesian equivalent.
The denominator has no β term thus can be ignored which reduces our equation to this.
Taking the log,
Now, let’s treat each term separately. Just as the frequentist approach, the conditional probability of y given β is assumed to be Gaussian. Therefore, the probability becomes:
Taking the log of the function,
Now for the second part. The prior in this case is hard to generalize. Making assumptions on βs without targets really does not make sense. But we are going to take a leap of faith here and claim that β is also a normally distributed with 0 mean and τ standard deviation as such.
Then the probability becomes:
Taking the log, the equation is changed to:
Putting it all together, our MAP becomes:
Since all the terms are negative, the ARGMAX is changed to ARGMIN with positive values. The constants can be taken out and the function is written in a vectorized form:
Again substituting σ²/τ² by some constant λ, the function now becomes:
For those familiar with classical linear regression, this is what is called Ridge regression (also called L2 regularization). Better predictions can be made by penalizing the β estimators by some constant. The Bayesian explanation of this is,when λ is increased, τ is reduced, squeezing the β normal distribution towards the center. This leads to smaller β estimations.
Get Abe Fa’s stories in your inbox
Join Medium for free to get updates from this writer.
Let’s take a few steps back and considered a different distribution instead, like a Laplace distribution instead. The β now becomes:
The probability and log of probability are now changed to:
Now the MAP function becomes a Lasso or L1 regularization:
It is important to establish the theoretical background of some regression techniques the reasoning of which is strongly linked to Bayesian approach.
Computing Bayesian Linear Regression
Now that the mathematical formulation is covered, how do we implement Bayesian linear regression. Remember the vectorized forms of β MAP estimates are summations of vectors. This means that methods used for integrating continuous functions for unknown functions can be applied here. The most common one is usually used in tandem with Bayesian statistics is the Markov Chain Monte Carlo (MCMC) method. Samples are drawn from likelihood and prior to compute the MAP distribution. Here it is important to note that unlike the classical linear regression, β is computed as a distribution, and not as a vector of coefficients.
Libraries and code samples
PyMC3 provides excellent tools for MCMC estimation. Some notable GitHub repositories that help to understand the implementation of PyMC3 are written by Michael Pyrcz (video lectures also found here) and a blog by Will Koehrsen.
Conclusion
Bayesian statistics has wide variety of applications. But in this blog, I tried to show how some of the methods used in the frequentist approach have Bayesian background. Although not described in detail, interested readers are encouraged to explore the resources provided here and other materials in relation to Bayesian linear regression. The blog might be updated with complex applications of Bayesian regression, hopefully soon. Finally, I would like to thank Ritvikmath for his excellent video on how to derive necessary equations for Bayesian linear regression.
