Deriving Linear Regression Formula

One day, during a boring lecture, I started thinking about linear regression, its usefulness and how to derive/think of some of its formulations, I decided to create a single write-up where I try to make sense of it all.

Linear regression, at is heart is solving a, seemingly, simple problem by asking a question: 'Can we model the relationship between all points, using a \(y = mx + b\) function?'. In simpler terms: 'Is there a linear function approximating our data?'. Turns out the answer is 'yes', but to understand why, let's first write a simple formula to capture our question. \[y_i = \beta_1 + x_i \beta_2\]

We know that our points (measurements taken; \(y_i\), \(x_i\)), so it's only a matter of finding \(\beta_{1, 2}\). I really want to make it easier for my brain to process all this and see the bigger picture, so I'll use three arbitrary points to see what's going on and how would I go about finding those \(\beta_{1, 2}\) values. Let's say the points are \((2;1)\), \((3;6)\) and \((2;4)\). We can make the following system of equations: \[\begin{cases}1 = \beta_1 + 2\beta_2 \\ 6 = \beta_1 + 3 \beta_2 \\ 4 = \beta_1 + 2 \beta_2 \end{cases}\]

Three points were specifically chosen to show, that this system is overfit, we have three points but only two unknowns. How can we deal with this? There's a nice 'trick' (sort of) that can help us, we introduce an 'error' (residual) value to each one of these like so: \[\begin{cases}1 = \beta_1 + 2\beta_2 + r_1 \\ 6 = \beta_1 + 3 \beta_2 + r_2 \\ 4 = \beta_1 + 2 \beta_2 + r_3 \end{cases}\]

This is very helpful, it may seem as if we introduced more complexity, but now we can treat this system of equations as functions of \(f_i(\beta_1, \beta_2) = r_i\). This will help us 'approximate' a linear function we're interested in. Okay how do we do that? By minimizing the 'error' (residuals). There's sadly no real insight into this (knowing that we need to minimize a special function we've just created), a lot of genius in math comes from very clever ingenuity and pattern recognition. You do get better at it, but you just have to keep going.

Okay, but, how do we minimize all of this (taking all these points and functions into account)? Well, we already know we can represent this system of equations as functions of \(\beta_{1, 2}\): \[\begin{cases} f_1(\beta_1, \beta_2) = r_1 = 1 - \beta_1 - 2\beta_2 \\ f_2(\beta_1, \beta_2) = r_2 = 6 - \beta_1 - 3\beta_2 \\ f_3(\beta_1, \beta_2) = r_3 = 4 - \beta_1 - 2\beta_2 \end{cases}\]

Now, instead of thinking of these as individual functions, we can get their sum and minimize that instead: \[\sum_{i=1}^N f_i(\beta_1, \beta_2)\]

This makes a lot of sense, since we want all \(r_i\)'s to be as close to zero as possible. That's almost it, there's one small thing we need take care of, 'negative error' (residuals). Currently it is very much possible for one of the \(f_i(\beta_1, \beta_2)\) functions to take on a negative value, there's a very common trick to deal with it though and it's called Least Squares Method used in all sorts of regression models. This simply squares the sum we're interested in, making 'negative error' impossible, there's also a small added bonus as well! The smaller our \(r_i\) is the smaller the square will be and vice-versa. This leaves us with the following: \[g(\beta_1, \beta_2) = \sum_{i=1}^N [f_i(\beta_1, \beta_2)]^2\]

Minimization in this and many other cases is making use of derivatives and equating them to zero, so let's take some: \[g(\beta_1, \beta_2) = \sum_{i=1}^N [f_i(\beta_1, \beta_2)]^2 = \sum_{i=1}^N [y_i - \beta_1 - x_i\beta_2]^2\] \[\frac{\partial g}{\partial \beta_1} = -2\sum_{i=1}^N [y_i - \beta_1 - x_i\beta_2]\] \[\frac{\partial g}{\partial \beta_2} = -2\sum_{i=1}^N [y_i - \beta_1 - x_i\beta_2]x_i\]

Someone more astute might think that solving this doesn't necessarily mean we have a minimum, as the derivative can only tell us where the function plateaus. This is not an issue here however, as the second derivative will always be bigger than zero in both cases, so we know for a fact that we'll always get a minimum. Solving by plugging in the arbitrary points from the beginning gives us the solution: \[\begin{cases} -2\sum_{i=1}^N [y_i - \beta_1 - x_i\beta_2] = 6\beta_1 + 14\beta_2 - 22 \\ -2\sum_{i=1}^N [y_i - \beta_1 - x_i\beta_2]x_i = 14\beta_1 + 34\beta_2 - 56 \end{cases}\] \[\begin{cases} 6\beta_1 + 14\beta_2 - 22 = 0 \\ 14\beta_1 + 34\beta_2 - 56 = 0 \end{cases}\]

Solving for and plugging the values \(\beta_1\) and \(\beta_2\) into the function from the beginning we get: \[y = -\frac{9}{2} + \frac{7}{2}x\]

Having worked through a real example, we can now begin to formulate something more general, something that would work for an arbitrary number of points compressed into a single formula. Let's first expand out sums a little bit and use common sum rules: \[\frac{\partial g}{\partial \beta_1} = -2\left[\sum_{i=1}^N y_i - \sum_{i=1}^N \beta_1 - \beta_2 \sum_{i=1}^N x_i\right] = 2\left[N \beta_1 + \beta_2 \sum_{i=1}^N x_i - \sum_{i=1}^N y_i\right]\] \[\frac{\partial g}{\partial \beta_2} = -2\left[\sum_{i=1}^N y_i x_i - \beta_1 \sum_{i=1}^N x_i - \beta_2 \sum_{i=1}^N x_i^2 \right] = 2\left[\beta_1 \sum_{i=1}^N x_i + \beta_2 \sum_{i=1}^N x_i^2 - \sum_{i=1}^N y_i x_i\right]\]

After that we just repeat the previous step, equating derivatives to zero and manipulate the equation a little: \[\sum_{i=1}^N y_i = N \beta_1 + \beta_2 \sum_{i=1}^N x_i\] \[\sum_{i=1}^N y_i x_i = \beta_1 \sum_{i=1}^N x_i + \beta_2 \sum_{i=1}^N x_i^2\]

You'd think that we'll need to do some weird manipulations in order to get \(\beta_1\) and \(\beta_2\), but we can actually just simplify these relationships with matrices: \[\begin{bmatrix} \sum_{i=1}^N y_i \\ \sum_{i=1}^N y_i x_i \end{bmatrix} = \begin{bmatrix} \beta_1 \\ \beta_2 \end{bmatrix} \begin{bmatrix} N & \sum_{i=1}^N x_i \\ \sum_{i=1}^N x_i & \sum_{i=1}^N x_i^2 \end{bmatrix}\] \[\begin{bmatrix} \beta_1 \\ \beta_2 \end{bmatrix} = \begin{bmatrix} \sum_{i=1}^N y_i \\ \sum_{i=1}^N y_i x_i \end{bmatrix} \begin{bmatrix} N & \sum_{i=1}^N x_i \\ \sum_{i=1}^N x_i & \sum_{i=1}^N x_i^2 \end{bmatrix}^{-1}\]

That's it! Something like this is very easy to plug into a program/calculator with a linear-algebra library support. It was very nice doodling this derivation, I always encourage people to re-discover things, especially in math/engineering. Walking the the same 'pligrimage' (so to speak) can be very helpful in gaining additional insight and developing problem-solving patterns.