MLE for Linear Regression with Gaussian Errors

Hello! Welcome to your next lesson in our journey through the foundations of machine learning.

In our last two lessons, we derived the Ordinary Least Squares (OLS) solution for linear regression. First, we used multivariate calculus to find the parameters that minimize the sum of squared residuals. Then, we used linear algebra to interpret this as an orthogonal projection, yielding the same closed-form Normal Equations. Both approaches are powerful, but they answer how to find the solution, not why the "sum of squared errors" is the right criterion to minimize in the first place.

Today, we will address that "why". Our learning outcome is to derive the probabilistic interpretation of linear regression using Maximum Likelihood Estimation (MLE) with Gaussian errors. By framing the problem statistically, you'll see that the OLS solution is not just a convenient computational trick but a direct consequence of a plausible assumption about the nature of noise in our data. This provides a deep and satisfying justification for the methods we've learned so far.

This lesson builds on your knowledge of MLE from the probability module and connects it directly to the regression problem.

1. The Probabilistic Model: Assuming Gaussian Errors

So far, we have treated the error term ${ϵ}_i$ in our model $y_i={\mathbf{x}}_i^T\mathbf{\beta }+{ϵ}_i$ as just a "residual". To build a probabilistic model, we must make an explicit assumption about where these errors come from.

The most common and foundational assumption in linear regression is that the errors are Independent and Identically Distributed (i.i.d.) according to a Gaussian (Normal) distribution with a mean of 0 and a constant variance ${\sigma }^2$. We write this as:

$${ϵ}_i\sim \mathcal{N}(0,{\sigma }^2)$$

This assumption implies three things:

1. Zero Mean: The errors are not systematic. On average, the model's predictions are correct.
2. Constant Variance (Homoscedasticity): The spread of the errors is the same for all values of the features ${\mathbf{x}}_i$.
3. Normality: The errors follow a bell-shaped curve. Small errors are common, and very large errors are rare.
4. Independence: The error for one observation provides no information about the error for another observation.

Given this assumption about the errors, what can we say about the distribution of our target variable, $y_i$? Since $y_i$ is the sum of a deterministic term (${\mathbf{x}}_i^T\mathbf{\beta }$) and a Gaussian random variable (${ϵ}_i$), $y_i$ itself is a Gaussian random variable. Its distribution is:

$$y_i\mid {\mathbf{x}}_i;\mathbf{\beta },{\sigma }^2\sim \mathcal{N}({\mathbf{x}}_i^T\mathbf{\beta },{\sigma }^2)$$

This is the heart of the probabilistic interpretation: for any given input ${\mathbf{x}}_i$, the model doesn't predict a single value for $y_i$, but a full probability distribution centered at the regression line's prediction, ${\mathbf{x}}_i^T\mathbf{\beta }$.

This video provides a very clear explanation of this idea, starting with a simple univariate case.

The image below, which you may find familiar, perfectly illustrates this. The regression line gives the mean ($\mu =Y_{est}$) of a Gaussian distribution at each point, and the actual observations (Y) are samples from these distributions.

2. Constructing the Likelihood Function

Now that we have a probabilistic model, we can use the principle of Maximum Likelihood Estimation (MLE) to find the best parameters. Recall from the statistics module that MLE seeks to find the parameter values that maximize the probability (or probability density) of observing the data we have actually collected.

Our parameters are $\mathbf{\theta }=(\mathbf{\beta },{\sigma }^2)$. Our data is the set of $n$ pairs $({\mathbf{x}}_i,y_i)$. The likelihood function, $L(\mathbf{\beta },{\sigma }^2)$, is the joint probability density of observing all the $y_i$ values, given the $\mathbf{X}$ matrix and the parameters.

The probability density for a single observation $y_i$ is given by the PDF of the normal distribution we defined above:

$$P(y_i\mid {\mathbf{x}}_i;\mathbf{\beta },{\sigma }^2)=\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp \left(-\frac{(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2}{2{\sigma }^2}\right)$$

Because we assumed our observations are independent, the joint probability density is the product of the individual densities:

$$L(\mathbf{\beta },{\sigma }^2\mid \mathbf{X},\mathbf{y})=\prod _{i=1}^{n}P(y_i\mid {\mathbf{x}}_i;\mathbf{\beta },{\sigma }^2)=\prod _{i=1}^n\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp \left(-\frac{(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2}{2{\sigma }^2}\right)$$

This is our likelihood function. Our goal is to find the $\mathbf{\beta }$ and ${\sigma }^2$ that maximize this expression.

3. Maximizing the Log-Likelihood

Directly maximizing the likelihood function, which is a product of many terms, is often difficult. It's much easier to maximize its natural logarithm, the log-likelihood function $\ell (\mathbf{\beta },{\sigma }^2)=\log L(\mathbf{\beta },{\sigma }^2)$. Since the logarithm is a monotonically increasing function, the parameters that maximize the log-likelihood are the same as those that maximize the likelihood.

Let's derive the log-likelihood:

$$\ell (\mathbf{\beta },{\sigma }^2)=\log \left(\prod _{i=1}^n\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp \left(-\frac{(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2}{2{\sigma }^2}\right)\right)$$

The log of a product is the sum of the logs:

$$\ell (\mathbf{\beta },{\sigma }^2)=\sum _{i=1}^n\log \left(\frac{1}{\sqrt{2\pi {\sigma }^2}}\exp \left(-\frac{(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2}{2{\sigma }^2}\right)\right)$$

Using the properties $\log (ab)=\log (a)+\log (b)$ and $\log (e^x)=x$:

$$\ell (\mathbf{\beta },{\sigma }^2)=\sum _{i=1}^n\left[\log \left(\frac{1}{\sqrt{2\pi {\sigma }^2}}\right)-\frac{(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2}{2{\sigma }^2}\right]$$

Simplifying the terms:

$$\ell (\mathbf{\beta },{\sigma }^2)=\sum _{i=1}^n\left[-\frac{1}{2}\log (2\pi {\sigma }^2)-\frac{1}{2{\sigma }^2}(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2\right]$$ $$\ell (\mathbf{\beta },{\sigma }^2)=-\frac{n}{2}\log (2\pi )-\frac{n}{2}\log ({\sigma }^2)-\frac{1}{2{\sigma }^2}\sum _{i=1}^n(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2$$

This is our final expression for the log-likelihood function.

The Connection to Least Squares

Now, let's find the parameter vector $\mathbf{\beta }$ that maximizes this function. We want to solve:

$${\hat{\mathbf{\beta }}}_{MLE}=\arg \underset{\mathbf{\beta }}{\max }\left(-\frac{n}{2}\log (2\pi )-\frac{n}{2}\log ({\sigma }^2)-\frac{1}{2{\sigma }^2}\sum _{i=1}^n(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2\right)$$

Notice that the first two terms do not depend on $\mathbf{\beta }$. Therefore, maximizing the entire expression is equivalent to maximizing just the third term.

$${\hat{\mathbf{\beta }}}_{MLE}=\arg \underset{\mathbf{\beta }}{\max }\left(-\frac{1}{2{\sigma }^2}\sum _{i=1}^n(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2\right)$$

Since ${\sigma }^2$ is positive, the term $-1/(2{\sigma }^2)$ is a negative constant. Maximizing a function that is a negative constant multiplied by some value is equivalent to minimizing that value. Therefore:

$${\hat{\mathbf{\beta }}}_{MLE}=\arg \underset{\mathbf{\beta }}{\min }\sum _{i=1}^n(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2$$

This is precisely the cost function for Ordinary Least Squares!

This proves that the OLS estimator is identical to the Maximum Likelihood Estimator under the assumption of i.i.d. Gaussian errors.

$${\hat{\mathbf{\beta }}}_{MLE}={\hat{\mathbf{\beta }}}_{OLS}=({\mathbf{X}}^T\mathbf{X}{)}^{-1}{\mathbf{X}}^T\mathbf{y}$$

For a rigorous, step-by-step walkthrough of this derivation in a formal style, the following resource is excellent.

Test your understanding!

The log-likelihood function also depends on the variance, ${\sigma }^2$. Using the final expression for the log-likelihood and your knowledge of calculus, can you derive the Maximum Likelihood Estimate for the variance, ${\hat{\sigma }}_{MLE}^2$?

Hint: Take the partial derivative of $\ell (\mathbf{\beta },{\sigma }^2)$ with respect to ${\sigma }^2$, set it to zero, and solve for ${\sigma }^2$. Remember to substitute ${\hat{\mathbf{\beta }}}_{MLE}$ for $\mathbf{\beta }$ in the final expression.

Show answer

Our log-likelihood function is:

$$\ell (\mathbf{\beta },{\sigma }^2)=-\frac{n}{2}\log (2\pi )-\frac{n}{2}\log ({\sigma }^2)-\frac{1}{2{\sigma }^2}\sum _{i=1}^n(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2$$

Take the partial derivative with respect to ${\sigma }^2$:

$$\frac{\partial \ell }{\partial ({\sigma }^2)}=0-\frac{n}{2}\cdot \frac{1}{{\sigma }^2}-\left(-\frac{1}{2({\sigma }^2{)}^2}\right)\sum _{i=1}^n(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2$$

Set the derivative to zero:

$$-\frac{n}{2{\sigma }^2}+\frac{1}{2({\sigma }^2{)}^2}\sum _{i=1}^n(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2=0$$

Multiply by $2({\sigma }^2{)}^2$ to clear the denominators:

$$-n{\sigma }^2+\sum _{i=1}^n(y_i-{\mathbf{x}}_i^T\mathbf{\beta }{)}^2=0$$

Solving for ${\sigma }^2$ and substituting ${\hat{\mathbf{\beta }}}_{MLE}$ gives the MLE for the variance:

$${\hat{\sigma }}_{MLE}^2=\frac{1}{n}\sum _{i=1}^n(y_i-{\mathbf{x}}_i^T{\hat{\mathbf{\beta }}}_{MLE}{)}^2$$

This is the average of the squared residuals, which is an intuitive estimate for the error variance.

Conclusion

In this lesson, we have uncovered the deep statistical foundation of Ordinary Least Squares. By moving from a purely procedural approach to a probabilistic one, we have shown that the familiar method of minimizing squared errors is not arbitrary, but is in fact the optimal approach under the widely-used assumption of normally distributed errors.

Key Takeaways:

• Linear regression can be interpreted as a probabilistic model by assuming the errors ${ϵ}_i$ are independent and identically distributed (i.i.d.) from a Gaussian distribution, $\mathcal{N}(0,{\sigma }^2)$.
• This assumption implies that each observation $y_i$ is drawn from a normal distribution centered on the regression line: $y_i\sim \mathcal{N}({\mathbf{x}}_i^T\mathbf{\beta },{\sigma }^2)$.
• The principle of Maximum Likelihood Estimation (MLE) finds the parameters that make our observed data most probable.
• Maximizing the log-likelihood function for this model with respect to the regression coefficients $\mathbf{\beta }$ is mathematically equivalent to minimizing the sum of squared residuals.
• Thus, the OLS estimator is the Maximum Likelihood Estimator for a linear model with Gaussian errors, providing a powerful statistical justification for its use. MLE also provides a principled way to estimate the error variance ${\sigma }^2$.

Preview of the Next Lesson:
We have now established the OLS solution from three different perspectives: optimization, geometry, and probability. With this solid foundation, we are ready to analyze the statistical properties of the OLS estimator. In the next lesson, we will begin our exploration of these properties by proving the unbiasedness of the OLS estimators, a crucial concept for understanding the reliability of our model's parameters. We will also introduce the set of "classical linear model assumptions" that are necessary for these properties to hold.