With regarding to OLS estimators why $\hat\beta_1$ and $\hat\beta_0$ have a variance?
3 Answers
The result that we obtain from linear regression is a function of random variables, so the parameters are random variables. You can calculate variance for any random variable. The variance tells us how uncertain the estimates are (in absolutely noiseless data, or with an overfitting model, the variances would be zeros).
The data you are feeding into your OLS model are random draws from some underlying population. You could either be drawing from the joint distribution of the predictors and the outcome, or from the distribution of the outcome conditional on the predictors (so you consider the predictors fixed - this is the assumption in OLS).
In either case, your regression coefficient estimates will be random variables themselves, because they are functions of the random variable $y$. The statistical properties of your parameter estimates, like the mean and the variance, will depend on the properties of the outcomes, which is why the variance of parameter estimates depends on both the model (via the hat matrix) and the variance of the outcome.
The other answers are correct, but I think it might be helpful to simulate what is happening.
library(ggplot2)
set.seed(2023)
# Define sample size
#
N <- 100
# Define number of times to repeat the simulation
#
R <- 1000
# Fix values of x
#
x <- seq(0, 1, 1/(N - 1))
# Define conditional expected values of y as E[y|x] = 1 + 2x
#
Ey <- 1 + 2*x
B_01 <- B_02 <- B_03 <- B_11 <- B_12 <- B_13 <- rep(NA, R)
for (i in 1:R){
# Simulate iid Gaussian error terms
#
e1 <- rnorm(N, 0, 1)
e2 <- rnorm(N, 0, 2)
e3 <- rnorm(N, 0, 3)
# Define observed values of y as the sum of the expected value and the error
#
y1 <- Ey + e1
y2 <- Ey + e2
y3 <- Ey + e3
# Fit regressions and extract the estimated regression coefficients
#
L1 <- lm(y1 ~ x)
L2 <- lm(y2 ~ x)
L3 <- lm(y3 ~ x)
#
B_01[i] <- summary(L1)$coefficients[1, 1]
B_02[i] <- summary(L2)$coefficients[1, 1]
B_03[i] <- summary(L3)$coefficients[1, 1]
#
B_11[i] <- summary(L1)$coefficients[2, 1]
B_12[i] <- summary(L2)$coefficients[2, 1]
B_13[i] <- summary(L3)$coefficients[2, 1]
}
# Make a data frame of the coefficients
#
d_01 <- data.frame(
Estimate = B_01,
Coefficient = "Intercept",
Variance = "1"
)
d_02 <- data.frame(
Estimate = B_02,
Coefficient = "Intercept",
Variance = "2"
)
d_03 <- data.frame(
Estimate = B_03,
Coefficient = "Intercept",
Variance = "3"
)
d_11 <- data.frame(
Estimate = B_11,
Coefficient = "Slope",
Variance = "1"
)
d_12 <- data.frame(
Estimate = B_12,
Coefficient = "Slope",
Variance = "2"
)
d_13 <- data.frame(
Estimate = B_13,
Coefficient = "Slope",
Variance = "3"
)
d <- rbind(d_01, d_02, d_03, d_11, d_12, d_13)
# Plot
#
ggplot(d, aes(x = Estimate, fill = Variance)) +
geom_density(alpha = 0.25) +
facet_grid(~Coefficient) +
theme(legend.position="bottom")
As the variance of the error term gets larger, the estimated slope $\hat\beta_1$ and estimated intercept $\hat\beta_0$ bounce around more.
In terms of the math, below is a common way to write the OLS estimates, which may show why the coefficient estimates vary.
$$ \hat\beta_1 = \dfrac{ \overset{N}{\underset{i = 1}{\sum}}\left[ (x_i - \bar x)(y_i - \bar y) \right] }{ \overset{N}{\underset{i = 1}{\sum}}\left[ (x_i - \bar x)^2 \right] }\\ \hat\beta_0 = \bar y - \hat\beta_1\bar x $$
Since $y$ is a random variable (due to the randomness of the error term), each of these estimates is random and depends on the exact observed values of $y$. As these observed values of $y$ change, so do the estimated coefficients. If those observed values of $y$ change a lot, such as with a large error variance, then the estimated coefficients change a lot which is the exact behavior in the above simulation.
(It could also be argued that $x$ has randomness, such as in an observational study. However, this is not necessary for the discussion of why the coefficient estimates have variances, and it adds confusion without any clear benefit, so I am setting aside $x$ randomness.)