# Distribution of OLS predictions

Suppose: $$y = X\beta + \varepsilon$$, with $$\varepsilon \sim (0, \Omega) \Rightarrow y|X \sim (X \beta, \Omega)$$
$$\hat{\beta}_{ols} = (X'X)^{-1}X'Y = \beta +(X'X)^{-1}X'\varepsilon \sim (\beta, \Sigma)$$ with $$\Sigma \equiv (X'X)^{-1} (X'\Omega X) (X'X)^{-1}$$ Assume $$\hat{\Omega}\to^p \Omega$$, let $$\hat{\Sigma} \equiv (X'X)^{-1} (X'\hat{\Omega} X) (X'X)^{-1}$$.
Under mild conditions on $$\varepsilon$$, the asymptotic distribution $$\hat{\beta}_{ols} \sim_a N(\beta, \hat{\Sigma})$$.

Question: what is the asymptotic distribution of $$\hat{y}_{ols}\equiv X\hat{\beta}_{ols}$$?
I think: $$\hat{y}_{ols}\equiv X\hat{\beta}_{ols} \sim_a N(X\beta, X\hat{\Sigma} X')$$

Some authors online say: $$\hat{y} \sim N(X\hat{\beta}, \hat{\Omega})$$, just plugging in estimates $$\hat{\beta}, \hat{\Omega}$$.
Update: I think they mean for $$x_{new}$$ the best prediction for the distribution of
$$y|x_{new}$$ is $$N(x_{new}\hat{\beta}, \hat{\Omega})$$.
But why not $$N(x_{new}\hat{\beta}, x_{new}\hat{\Sigma}x_{new}')$$?

We agree $$E[y|x_{new}] = x_{new}\hat{\beta}$$, the difference might be
$$V[y|x_{new}]=\hat{\Omega}$$ vs $$V[E[y|x_{new}]]=x_{new}\hat{\Sigma}x_{new}'$$?

If we assume sperical errors & finite sample normality:
$$\varepsilon \sim N(0, \sigma^2 I_N) \Rightarrow \Sigma \equiv \sigma^2 (X'X)^{-1} \Rightarrow \hat{y}_{ols} \sim N(X\beta, \sigma^2 X (X'X)^{-1} X')$$.

Update 2: perhaps we are computing the distribution of different things.
Them: $$y|X \sim N(X \beta, \Omega)$$, estimate $$\hat{\beta}, \hat{\Omega}$$, plug-in: $$y|x_{new} \sim N(x_{new} \hat{\beta}, \hat{\Omega})$$.
Me: $$\hat{E}[y|x_{new}] = x_{new} \hat{\beta} \sim N(x_{new} \beta, x_{new}\hat{\Sigma}x_{new}')$$

Their estimand = the predicted distribution of $$y|x_{new}$$
My estimand = the distribution of the point prediction (conditional expectation) $$\hat{E}[y|x_{new}]$$
In general standard errors from $$\hat{\Omega}$$ are larger than $$x_{new}\hat{\Sigma}x_{new}'$$.
I was wrong when I wrote: Some authors online say: $$\hat{y} \sim N(X\hat{\beta}, \hat{\Omega})$$
I should have wrote: Some authors online say: $$\hat{y|X} \sim N(X\hat{\beta}, \hat{\Omega})$$

Consider a random sample $$\{x_i\}_{i=1}^n \sim N(\mu, \sigma^2)$$ w/ $$\sigma^2$$ known.
$$\hat{\mu}_{n} \equiv \bar{x}_{n} \sim N(\mu, \frac{\sigma^2}{n})$$
What is your best prediction for new data from same DGP $$x_{new}$$?
Predicted distribution: $$x_{new}\sim N(\hat{\mu}_{n}, \sigma^2)$$
Distribution of point prediction: $$E[x_{new}] = \hat{\mu}_{n} \sim N(\mu, \frac{\sigma^2}{n})$$
Question: how do I compute standard errors around the predicted distribution $$x_{new}$$ in the space of probability distributions?

Julia code to illustrate:

using MLJ, LinearAlgebra; using MLJ: matrix
X, y =  @load_boston; x =[ones(506) matrix(X)]; n,p =size(x);
β̂ = x\y; ŷ = x*β̂; ε̂ = y-ŷ;
s² = (ε̂'ε̂)/(n-p); Σ̂ = s²*(x'x)^-1;

√s²                 # their se(y|x) = 4.79
sqrt.(diag(x*Σ̂*x')) # my se(ŷ) =[.61 .49 .51 ...]


Summary:
Suppose: $$y = X\beta + \varepsilon$$, with $$\varepsilon \sim (0, \Omega) \Rightarrow y|X \sim (X \beta, \Omega)$$
$$\hat{\beta}_{ols} = (X'X)^{-1}X'Y \sim (\beta, \Sigma)$$ with $$\Sigma \equiv (X'X)^{-1} (X'\Omega X) (X'X)^{-1}$$
Assume $$\hat{\Omega}\to^p \Omega$$, let $$\hat{\Sigma} \equiv (X'X)^{-1} (X'\hat{\Omega} X) (X'X)^{-1}$$.

$$\hat{\beta}_{ols} \sim_a N(\beta, \hat{\Sigma})$$.
$$\hat{E}[y|x]=\hat{y}_{ols}\equiv X\hat{\beta}_{ols} \sim_a N(X\beta, X\hat{\Sigma} X')$$ (sample mean dist around truth)
$$\hat{y|x} \sim_a (X\hat{\beta}_{ols}, \hat{\Omega})$$
Note, from CLT we know the asymptotic distribution of $$\hat{E}[y|x]$$, but not $$\hat{y|x}$$.
If we make stronger assumptions, eg the finite sample distribution $$\varepsilon \sim N(0, \Omega)$$, then we can conclude $$\hat{y|x} \sim N(X\hat{\beta}_{ols}, \hat{\Omega})$$.

• Perhaps this means that $\Omega=X\Sigma X’$. I confess that I’ve forgotten this fact (if it’s even true), but can you prove it?
– Dave
Commented May 29, 2020 at 3:34
• Can you link to one of the authors online who think the variance is $\Omega$. I agree with you that they're wrong, but if you want to understand why, we'd need to see more. Commented May 29, 2020 at 3:42
• @ThomasLumley I prefer not to post direct links here. Btw, do you agree w/ my logic on the distribution of $\hat{y}$? Commented May 29, 2020 at 3:45
• Yes. I do. But I expect the people who think the variance is $\Omega$ are trying to do something else, rather than just getting this calculation wrong. Commented May 29, 2020 at 3:48
• @ThomasLumley Discussion here: github.com/alan-turing-institute/MLJ.jl/issues/… Commented May 29, 2020 at 3:57