How The Linear Regression Got Its Bowtie.

Suppose we want to investigate the relationship between height($y$, cm) and weight($x$, kg) in a population. We have a hypothesis : any given observation of height is generated by the following data generating function (DGP):

Say we have 5 people that are the exact same weight. Do we expect them to have the exact same height ? with the above definition we say no. The error, $\epsilon$, contributes to all the randomness in $y$ — there will be a distribution around any given $x$ with a local mean defined as a lineIt is a simple claim on the face of it, but it is actually a remarkably bold assumption. We assume that the relationship between $x$ and $y$ is a simple line. Think of all the factors that influence an individuals height at a given weight — genetics, childhood nutrition, developmental disorders — the world is a complex place. However, we just say : doesn't matter, all of that complexity is basically noise compared to the simple linear relationship..

Let's fix the DGP and $X$ as follows:

n  = 100
sigma = 15
X = np.random.normal(75,10,n)
beta_0 = 1
beta_1 = 2
eps = np.random.normal(0,sigma, n)
Y = beta_0 + beta_1 * X + eps

With an ordinary linear regression (OLS) we model the following:

Notice that $\mathbb{E}$ is a single number, not a distribution We are asking what is the mean of Y at each value of $x$?. The noise hasn't disappeared, it's still in $Y$. We assumed $\epsilon \sim \mathcal{N}(0,\sigma^2)$ so a mean-zero noise term drops out under expectation.. When we plug any given $x$ into a linear regression, we are asking the question : What is the expected value of the distribution of $y$ around this $x$ ?

Suppose we have a fitted linear regression and we want to know $y$ at $x = 85$. We can see this from sampling from the DGP.

The blue line is the single number provided by the model. It has collapsed all of the local uncertainty into its prediction of the mean.

Prediction intervals - Two pieces of uncertainty.

Uncertainty around any single line

Usually, we only have one sample of data —$(x_i,y_i)$ pairs — from the population to work with. With this sample, we fit one best-fit line. Let's ignore boostrapping for now, since we assume we have enough samples and assumptions for linear regression are met.

Figure 3 shows the variation of our observed examples around the fitted line. What is the scatter of points around the fitted line ? Recall that (by definition) we say that anything not captured by our model is noise with variation $\sigma$. We can estimate this randomness as : ⊕We divide by $n - p - 1 = n-2$, where $p$ is the number of covariates. It is a generalization of Bessel's correction [1]. The formula is undefinewd for $n\leq 2$. With 2 points, the line passes perfectly through both lines. Both residuals are 0 and $\frac{0}{0}$ has no meaning. This is the math saying : you can't estimate noise from a perfect fit. You need at least one point of slack — $n \geq 3$ — for $\hat{\sigma}$ to be defined at all.

Even if we knew exactly what the data generating line was, we would still have some randomness around it (that is our assumption of the data generating process) Here we assume that the model is correctly specified. For this example, this is exactly aleatoric uncertainty because we defined the DGP as a linear function with some noise. A nuance on this is discussed later.. Note that $\sigma$ is a fixed property of the world (as seen by our model) while $\hat{\sigma}$ is our estimate — it becomes more reliable , converging to $\sigma$, as we get more data (increasing $n$). The scatter we see around the model is our model-based definition of aleatoric uncertaintyUncertainty as a function of the inherent randomness in the world. It's irreducible - No amount of additional data eliminates it..

Why can't we just slap $\pm \hat{\sigma}$ around a point prediction and call it a day? $\hat{\sigma}$ does capture randomness — our estimate of how much the points scatter around the line — but notice, there's no $x$ involved. The model assumes the noise is the same size at every $x$This is the homoscedasticity of errors assumptions. That is, the model doesn't systematically perform worse for some parts of the input space. so $\hat{\sigma}$ is a single scalar that summarizes the vertical scatter at every $x$. However, notice that our fitted line $\mathbb{E}[y|x]$ is built from a finite sample. The estimates of $\{\beta_0, \beta_1\}$ could have been different with a different sample from the same random process. That uncertainty does depend on $x$ but is absent in $\hat{\sigma}$. This uncertainty about the line itself, which depends on which $x$ you happen to sample, is what's missing.

Uncertainty about the line itself

What about the line itself ? Since we have the full data generating process, let's fit a few separate lines to 6 independent samples of data.

Figure 4 shows how the variance between samples effects the line. All of these are plausible lines generated from the same DGP.

Now, if we got the prediction for $x = 85$ from each of these models (fitted over many sample datasets), we get the following distribution of $\mathbb{E}[y|x = 85]$.

This is the uncertainty we have about the model itself. Specifically, the model parameters because we have fixed the family of the model to be a linear one. It is the standard error of the conditional meanRecall that the linear regression gives $\mathbb{E}[y|x]$.. It quantifies how much a prediction for a given input could vary over many possible model fits.

Let's look at the terms inside the square root, first:

This term measures the distance of the new example from the sample mean of $x$ relative to the total spread of $x$. This is essentially a surprise score. If $x_{\text{new}}$ is very far from the mean but the overall spread of $x$ is wide (large denominator), it's not very surprising. A wider spread of training $x$'s puts data across more of the input space, which pins down the line better. The closer a given input is to the mean, the more confident you can be about the prediction there.

The term is zero when $x_{\text{new}} = \bar{x}$ — predicting at the center of your data is the most confident you can be. As $x_{\text{new}}$ moves away from $\bar{x}$, the term grows quadratically, but a wide spread of training $x$'s damps that growth — you've earned the right to be confident further out. This quadratic growth is the mathematical signature of the bowtie shape we see in the bands : narrow at the center of the data, flaring outward as you move toward the edges.

The $\frac{1}{n}$ term keeps the value in the square root strictly positive — even when $x_{\text{new}} = \bar{x}$ zeros out the surprise term. Suppose the noise is high (large $\hat{\sigma}$) but we have lots of data. The noise on individual points averages out across many points, and the line gets pinned down even through a noisy cloud. We can't make the world less noisy, but we can buy a more confident estimate of the line by collecting more data.

This can be seen in the following animation. For the first few data points, the fitted lineAlways important to keep in mind that the fitted line gives $\mathbb{E}[y|x]$ — the mean of the distribution of $y$ at each $x$, not a single deterministic prediction. jumps around a lot — a few points do very little to pin down the line. As new points arrive from around in the input space, the line settles towards true line.This isn't $\hat{se}_{line}$ shrinking per se. That would be best illustrated with many lines fitted with $m$ samples of size $n$, with increasing $n$. It would show the lines cluster together as $n$ increases.

So the line uncertainty has two separable knobs: $\hat{\sigma}$ — our estimate of the noise in the world — and $\frac{1}{n} + \text{surprise}$ — the information factor, which we control through sample size and where we choose to sample. Drive the information factor toward zero with enough well-spread data, and the line uncertainty vanishes. What remains is just $\sigma$ — the irreducible scatter of new points around a now-well-known line.

With $\hat{\text{se}}_{\text{line}}$, we can calculate the confidence band for a specific $x$ within which we can expect a line to fall:

where $t$ is how many standard deviations wide do I want the interval to be? (e.g., 1.98We have a finite sample estimate of $\sigma$ which in itself is noisy. A different sample would have given us a different estimate. Here, what we want are the bounds that capture some proportion (e.g., 95%) of the values $\mathbb{E}[y|x_{new}]$ can take. For that, we need a distribution. We use a $t$-distribution here because it has fatter tails, reflecting that we should be more uncertain as compared to using a normal distribution. to capture 95% of the probability mass). And now, for a range of test inputs, we can create a band around the entire line which shrinks and grows depending on the value of $x$.

This is our model-based definition of epistemic uncertaintyUncertainty due to the lack of knowledge. We could conceivable improve our model with more data (keywords being our model, more below).

Bringing the parts together

What we are after is the prediction interval — a band that tells us with some confidence where a new, unseen point will land. We need a way to combine the effects of both $\hat{\sigma}$ and $\hat{se}_{\text{line}}$ into a standard error of prediction, $\hat{se}_{\text{pred}}$.

Recall that we have two sources of uncertainty about the true $y$ of a new example: a) where the true line might actually be at $x_{\text{new}}$, and b) how a new observation will scatter around that line. Both contribute to the error of a new prediction, and both are measured in standard deviation units.

We can't just add them directly. Adding standard deviations would assume the two errors always reinforce each other — but they don't. They're independent: sometimes they push in the same direction (the line happens to be too high and the new point happens to land high, doubling up), and sometimes they push in opposite directions (the line is too high but the new point lands low, partially cancelling). Over many predictions, the average magnitude of their combined effect is less than the sum of their individual spreads. Picture the prediction error $y_{\text{new}} - \hat{y}(x_{\text{new}})$ as a quantity being pushed around by two independent random forces — one from the line being off, one from the new point's noise. The Pythagorean combination accounts for the mix of reinforcement and cancellation honestly.

To combine them properly, we work in variance units (where the cancellation/reinforcement bookkeeping works out cleanly), then translate back to a standard deviation:

and finally, to establish the bounds for a new observation, we multiply by a desired confidence level to get the final prediction interval:

A tale of two uncertainties ?

Why have I constantly reminded you that all of this is model-based ? Let's look at all the calculations in one place.

Notice how $\hat{\sigma}$ shows up everywhere and it is calculated as the error between the model and the true $y$. We labelled this as aleatoric uncertainty — irreducible uncertainty about the world — but consider the following hypothetical DGP and fitted linear model.

n  = 100
sigma = 15
X = np.random.normal(75, 10, n)
beta_0 = 5
beta_1 = 1 
eps = np.random.normal(0, sigma, n)

Y = beta_0 + beta_1 * (X - 60)**2 + eps

We are using the model to estimate $\hat{\sigma}$ and if the model is badly misspecified, the world looks much noisier than it really is. This is because, according to the model, the world is linear and any deviation from this hypothetical line is irreducible noise. However, clearly, if we had a better model (i.e., more knowledge), we could reduce this "irreducible" uncertainty.

So, the split between epistemic and aleatoric uncertainty is not as clean as it sounds. One might say it's all epistemic uncertainty and all we lack is a better model and more data. We could be missing a crucial covariate (e.g., gender) in our modelFixes structural ignorance which reduces that apparent aleatoric uncertainty. or perhaps we could do a wider samplePins down the parameters better which reduces epistemic uncertainty once we've made a commitment to a model.. I find this to be a more optimistic position than resigning parts of the world to inscrutable randomness. It says that we can do better, we just have to keep going. That, though is philosophical debate better tackled elsewhere.

Thanks for reading.