# How to derive the least square estimator for multiple linear regression?

In the simple linear regression case $y=\beta_0+\beta_1x$, you can derive the least square estimator $\hat\beta_1=\frac{\sum(x_i-\bar x)(y_i-\bar y)}{\sum(x_i-\bar x)^2}$ such that you don't have to know $\hat\beta_0$ to estimate $\hat\beta_1$

Suppose I have $y=\beta_1x_1+\beta_2x_2$, how do I derive $\hat\beta_1$ without estimating $\hat\beta_2$? or is this not possible?

• You can omit one of the variables and still obtain an unbiased estimate of the other if they are independent. Mar 17, 2014 at 0:10
• Oct 26, 2021 at 0:16

The derivation in matrix notation

Starting from $y= Xb +\epsilon$, which really is just the same as

$\begin{bmatrix} y_{1} \\ y_{2} \\ \vdots \\ y_{N} \end{bmatrix} = \begin{bmatrix} x_{11} & x_{12} & \cdots & x_{1K} \\ x_{21} & x_{22} & \cdots & x_{2K} \\ \vdots & \ddots & \ddots & \vdots \\ x_{N1} & x_{N2} & \cdots & x_{NK} \end{bmatrix} * \begin{bmatrix} b_{1} \\ b_{2} \\ \vdots \\ b_{K} \end{bmatrix} + \begin{bmatrix} \epsilon_{1} \\ \epsilon_{2} \\ \vdots \\ \epsilon_{N} \end{bmatrix}$

it all comes down to minimzing $e'e$:

$\epsilon'\epsilon = \begin{bmatrix} e_{1} & e_{2} & \cdots & e_{N} \\ \end{bmatrix} \begin{bmatrix} e_{1} \\ e_{2} \\ \vdots \\ e_{N} \end{bmatrix} = \sum_{i=1}^{N}e_{i}^{2}$

So minimizing $e'e'$ gives us:

$min_{b}$ $e'e = (y-Xb)'(y-Xb)$

$min_{b}$ $e'e = y'y - 2b'X'y + b'X'Xb$

$\frac{\partial(e'e)}{\partial b} = -2X'y + 2X'Xb \stackrel{!}{=} 0$

$X'Xb=X'y$

$b=(X'X)^{-1}X'y$

One last mathematical thing, the second order condition for a minimum requires that the matrix $X'X$ is positive definite. This requirement is fulfilled in case $X$ has full rank.

The more accurate derivation which goes trough all the steps in greater dept can be found under http://economictheoryblog.com/2015/02/19/ols_estimator/

• This derivation is precisely what I was searching for. NO SKIPPED STEPS. Surprising how difficult to find same. Mar 20, 2015 at 1:33
• In the matrix equation, shouldn't the second * be a +? Also, shouldn't it be $b_K$ instead of $b_N$ to get the dimensions to match? Sep 7, 2017 at 16:45
• Alexis Olson, you are right! I edited my answer. Sep 8, 2017 at 5:43

It is possible to estimate just one coefficient in a multiple regression without estimating the others.

The estimate of $\beta_1$ is obtained by removing the effects of $x_2$ from the other variables and then regressing the residuals of $y$ against the residuals of $x_1$. This is explained and illustrated How exactly does one control for other variables? and How to normalize (a) regression coefficient?. The beauty of this approach is that it requires no calculus, no linear algebra, can be visualized using just two-dimensional geometry, is numerically stable, and exploits just one fundamental idea of multiple regression: that of taking out (or "controlling for") the effects of a single variable.

In the present case the multiple regression can be done using three ordinary regression steps:

1. Regress $y$ on $x_2$ (without a constant term!). Let the fit be $y = \alpha_{y,2}x_2 + \delta$. The estimate is $$\alpha_{y,2} = \frac{\sum_i y_i x_{2i}}{\sum_i x_{2i}^2}.$$ Therefore the residuals are $$\delta = y - \alpha_{y,2}x_2.$$ Geometrically, $\delta$ is what is left of $y$ after its projection onto $x_2$ is subtracted.

2. Regress $x_1$ on $x_2$ (without a constant term). Let the fit be $x_1 = \alpha_{1,2}x_2 + \gamma$. The estimate is $$\alpha_{1,2} = \frac{\sum_i x_{1i} x_{2i}}{\sum_i x_{2i}^2}.$$ The residuals are $$\gamma = x_1 - \alpha_{1,2}x_2.$$ Geometrically, $\gamma$ is what is left of $x_1$ after its projection onto $x_2$ is subtracted.

3. Regress $\delta$ on $\gamma$ (without a constant term). The estimate is $$\hat\beta_1 = \frac{\sum_i \delta_i \gamma_i}{\sum_i \gamma_i^2}.$$ The fit will be $\delta = \hat\beta_1 \gamma + \varepsilon$. Geometrically, $\hat\beta_1$ is the component of $\delta$ (which represents $y$ with $x_2$ taken out) in the $\gamma$ direction (which represents $x_1$ with $x_2$ taken out).

Notice that $\beta_2$ has not been estimated. It easily can be recovered from what has been obtained so far (just as $\hat\beta_0$ in the ordinary regression case is easily obtained from the slope estimate $\hat\beta_1$). The $\varepsilon$ are the residuals for the bivariate regression of $y$ on $x_1$ and $x_2$.

The parallel with ordinary regression is strong: steps (1) and (2) are analogs of subtracting the means in the usual formula. If you let $x_2$ be a vector of ones, you will in fact recover the usual formula.

This generalizes in the obvious way to regression with more than two variables: to estimate $\hat\beta_1$, regress $y$ and $x_1$ separately against all the other variables, then regress their residuals against each other. At that point none of the other coefficients in the multiple regression of $y$ have yet been estimated.

• Great answer, here is a general theorem en.wikipedia.org/wiki/… Dec 6, 2015 at 14:45
• Thank you for this insightful answer - a stupid question I have as I'm learning is, how do you know that your $\hat{\beta}_1$ you write in step (3) is in fact $\hat{\beta}_1$? Like what is the "proof" that the $\hat{\beta}_1$ obtained from this 3 step process is the same $\hat{\beta}_1$ you get from the original linear model formulation? Apr 5 at 17:56
• @HavelTheGreat There are myriad ways to prove this. If you don't like the geometrical demonstration I gave, you can translate it into linear algebra in various ways. But the geometry would have been familiar to Euclid and is perfectly rigorous: see the first link in my answer. It comes down to proving "we could adjust the match to get it even closer." That follows from the Pythagorean Theorem.
– whuber
Apr 5 at 18:04
• Yeah chasing the links in your answer I'm getting there, thank you. Sorry if this is off topic, but do you have any text recommendations that go into this process in detail? I'm trying to work through some variation of the same problem but extended to multiple regression, your answers have been super helpful. Apr 5 at 18:08
• @HavelTheGreat The Tukey & Mosteller reference ("Data Analysis and Regression") bases its account of multiple regression on this sequential approach. There occasionally appear articles in places like The American Statistician about this: "Gram-Schmidt" is a good keyword to search.
– whuber
Apr 5 at 18:13

A simple derivation can be done just by using the geometric interpretation of LR.

Linear regression can be interpreted as the projection of $$Y$$ onto the column space $$X$$. Thus, the error, $$\hat{\epsilon}$$ is orthogonal to the column space of $$X$$.

Therefore, the inner product between $$X'$$ and the error must be 0, i.e.,

$$ = 0$$

$$X'y - X'X\hat{\beta} = 0$$

$$X'y = X'X\hat{\beta}$$

Which implies that,

$$(X'X)^{-1}X'y = \hat{\beta}$$.

Now the same can be done by:

(1) Projecting $$Y$$ onto $$X_2$$ (error $$\delta = Y-X_2 \hat{D}$$), $$\hat{D} = (X_2'X_2)^{-1}X_2'y$$,

(2) Projecting $$X_1$$ onto $$X_2$$ (error $$\gamma = X_1 - X_2 \hat{G}$$), $$\hat{G} = (X_2'X_2)^{-1}X_2X_1$$,

and finally,

(3) Projecting $$\delta$$ onto $$\gamma$$, $$\hat{\beta}_1$$

The ordinary least squares estimate of $\beta$ is a linear function of the response variable. Simply put, the OLS estimate of the coefficients, the $\beta$'s, can be written using only the dependent variable ($Y_i$'s) and the independent variables ($X_{ki}$'s).

To explain this fact for a general regression model, you need to understand a little linear algebra. Suppose you would like to estimate the coefficients $(\beta_0, \beta_1, ...,\beta_k)$ in a multiple regression model,

$$Y_i = \beta_0+\beta_1X_{1i}+...+\beta_kX_{ki}+\epsilon_i$$

where $\epsilon_i \overset{iid}{\sim} N(0,\sigma^2)$ for $i=1,...,n$. The design matrix $\mathbf{X}$ is a $n\times k$ matrix where each column contains the $n$ observations of the $k^{th}$ dependent variable $X_k$. You can find many explanations and derivations here of the formula used to calculate the estimated coefficients $\boldsymbol{\hat{\beta}}=(\hat{\beta}_0, \hat{\beta}_1, ..., \hat{\beta}_k)$, which is

$$\boldsymbol{\hat{\beta}}=(\mathbf{X}^\prime \mathbf{X})^{-1}\mathbf{X}^\prime \mathbf{Y}$$

assuming that the inverse $(\mathbf{X}^\prime \mathbf{X})^{-1}$ exists. The estimated coefficients are functions of the data, not of the other estimated coefficients.

• I have a follow up question, on the simple regression case, you make $y_i=\beta_0+\beta_1\bar x+\beta_1(x_i-\bar x)+e_i$ then $X$ becomes a matrix of $(1,...,1)$ and $(x_1-\bar x,...,x_n-\bar x)$, then follow through the $\hat\beta=(X'X)^(-1)X'Y$. How should I rewrite the equation in my case? Dec 18, 2012 at 8:20
• And 1 more question, does this apply to cases where $x_1$ and $x_2$ are not linear, but the model is still linear? For example the decay curve $y=\beta_1 e^{x_1t}+\beta_2 e^{x_2t}$, can I substitute the exponential with $x_1'$ and $x_2'$so it becomes my original question? Dec 18, 2012 at 8:52
• In your first comment, you can center the variable (subtract its mean from it) and use that is your independent variable. Search for "standardized regression". The formula you wrote in terms of matrices is not correct. For your second question, yes you may do that, a linear model is one that is linear in $\beta$, so as long as $y$ equal to a linear combination of $\beta$'s you are fine. Dec 18, 2012 at 9:01
• (+1). But shouldn't it be "$n \times k$ matrix" instead of $k \times n$? Dec 18, 2012 at 11:18

One small minor note on theory vs. practice. Mathematically $\beta_0, \beta_1, \beta_2 ... \beta_n$ can be estimated with the following formula:

$$\hat{\beta} = (X'X)^{-1} X'Y$$

where $X$ is the original input data and $Y$ is the variable that we want to estimate. This follows from minimizing the error. I will proove this before making a small practical point.

Let $e_i$ be the error the linear regression makes at point $i$. Then:

$$e_i = y_i - \hat{y_i}$$

The total squared error we make is now:

$$\sum_{i=1}^n e_i^2 = \sum_{i=1}^n (y_i - \hat{y_i})^2$$

Because we have a linear model we know that:

$$\hat{y_i} = \beta_0 + \beta_1 x_{1,i} + \beta_2 x_{2,i} + ... + \beta_n x_{n,i}$$

Which can be rewritten in matrix notation as:

$$\hat{Y} = X\beta$$

We know that

$$\sum_{i=1}^n e_i^2 = E'E$$

We want to minimize the total square error, such that the following expression should be as small as possible

$$E'E = (Y-\hat{Y})' (Y-\hat{Y})$$

This is equal to:

$$E'E = (Y-X\beta)' (Y-X\beta)$$

The rewriting might seem confusing but it follows from linear algebra. Notice that the matrices behave similar to variables when we are multiplying them in some regards.

We want to find the values of $\beta$ such that this expression is as small as possible. We will need to differentiate and set the derivative equal to zero. We use the chain rule here.

$$\frac{dE'E}{d\beta} = - 2 X'Y + 2 X'X\beta = 0$$

This gives:

$$X'X\beta = X'Y$$

Such that finally: $$\beta = (X'X)^{-1} X'Y$$

So mathematically we seem to have found a solution. There is one problem though, and that is that $(X'X)^{-1}$ is very hard to calculate if the matrix $X$ is very very large. This might give numerical accuracy issues. Another way to find the optimal values for $\beta$ in this situation is to use a gradient descent type of method. The function that we want to optimize is unbounded and convex so we would also use a gradient method in practice if need be.

• except that you don't actually need to compute $(X'X)^{-1}$... Dec 18, 2012 at 13:43
• valid point. one could also use the gram schmidt process, but I just wanted to remark that finding the optimal values for the $\beta$ vector can also be done numerically because of the convexity. Dec 19, 2012 at 17:44