Sampling Distribution of OLS

The slope is a random variable

When you run reg y x, Stata hands back one number for the slope. It looks final. But that number is an accident of the particular sample you happened to draw. Draw a fresh sample from the same population and the line tilts to a slightly different angle; \(\hat{\beta}_1\) comes back a different number. The estimator \(\hat{\beta}_1\) is not a number — it is a recipe that turns whatever data you feed it into a slope. Because the data are random, its output is random too.

The sampling distribution of \(\hat{\beta}_1\) is the histogram you would get by running that recipe on every sample you could have drawn. Everything in this part of the course — unbiasedness, the variance formula, standard errors, \(t\)-stats — is a statement about that histogram. This page lets you build it.

Where the proof’s weights come from

Start from the slope formula you derived: \[\hat{\beta}_1 = \frac{\sum_i (x_i - \bar{x})(y_i - \bar{y})}{\sum_i (x_i - \bar{x})^2}.\]

Because \(\sum_i (x_i - \bar{x})\,\bar{y} = 0\), the \(\bar{y}\) drops out and you can write the slope as a weighted sum of the \(y\)’s: \[\hat{\beta}_1 = \sum_i w_i\, y_i, \qquad w_i = \frac{x_i - \bar{x}}{SST_x}, \qquad SST_x = \sum_i (x_i - \bar{x})^2.\]

The weights \(w_i\) are fixed the moment you know the \(X\)’s — larger for points far from \(\bar{x}\), which is why outlying \(X\) values swing the line hardest. Now substitute the true model \(y_i = \beta_0 + \beta_1 x_i + u_i\). Two of the three pieces collapse (\(\sum_i w_i = 0\) and \(\sum_i w_i x_i = 1\)), leaving the one line that drives the whole of Topic 2: \[\hat{\beta}_1 = \beta_1 + \sum_i w_i\, u_i.\]

Read it in words: your estimate equals the truth plus a weighted average of the errors. Everything follows from here.

  • Unbiasedness. Take expectations. With \(E(u\mid X) = 0\), each \(E(w_i u_i) = 0\), so \(E(\hat{\beta}_1) = \beta_1\). The estimate is centered on the truth — the histogram sits over \(\beta_1\).
  • Variance. The only random part is \(\sum_i w_i u_i\). With independent, constant-variance errors, \(\mathrm{Var}(\hat{\beta}_1) = \sigma^2 \sum_i w_i^2 = \sigma^2 / SST_x\) — the width of the histogram.

This is the same \(w_i\) that appears out of nowhere on the lecture slide. It is not a trick — it is just the slope written as a weighted sum. Once you see \(\hat{\beta}_1 = \beta_1 + \sum_i w_i u_i\), the unbiasedness and variance proofs are each one short step away.

Build it yourself

The simulator fixes a true line \(Y = 2 + 3X\) (so \(\beta_1 = 3\)), holds the \(X\)’s fixed, and redraws the errors \(u\) to generate one sample after another — exactly the “shuffling” from lecture, but you turn the knobs. The left panel shows a single sample: the true line (dashed) and the line OLS fits to that one draw (solid). They rarely coincide. The right panel collects the fitted slope from every sample into the sampling distribution.

Watch two things: the histogram centers on 3 (unbiasedness), and its width equals \(\sigma_u/\sqrt{SST_x}\) (the variance formula) — the red curve is that exact normal.

#| standalone: true
#| viewerHeight: 620

library(shiny)

ui <- fluidPage(
  tags$head(tags$style(HTML("
    .stats-box {
      background: #eaf2f8; border-radius: 6px; padding: 14px;
      margin-top: 12px; font-size: 14px; line-height: 1.8;
    }
    .stats-box b { color: #2c3e50; }
  "))),

  sidebarLayout(
    sidebarPanel(
      width = 3,

      sliderInput("n", "Sample size (n):",
                  min = 5, max = 500, value = 250, step = 5),

      sliderInput("sd_u", "Error SD (sigma_u):",
                  min = 1, max = 12, value = 6, step = 0.5),

      sliderInput("sd_x", "Spread of X (SD of X):",
                  min = 1, max = 8, value = 3, step = 0.5),

      sliderInput("reps", "Number of samples:",
                  min = 100, max = 3000, value = 1000, step = 100),

      uiOutput("stats_box")
    ),

    mainPanel(
      width = 9,
      fluidRow(
        column(6, plotOutput("scatter_plot",  height = "380px")),
        column(6, plotOutput("sampling_plot", height = "380px"))
      )
    )
  )
)

server <- function(input, output, session) {

  b0 <- 2
  b1 <- 3

  # X is held fixed across repeated samples: we condition on the sample
  # values X_1, ..., X_n. Changing n or the spread of X draws a fresh design.
  xvals <- reactive({
    rnorm(input$n, mean = 0, sd = input$sd_x)
  })

  sim <- reactive({
    x    <- xvals()
    n    <- input$n
    sd_u <- input$sd_u
    reps <- input$reps

    xbar  <- mean(x)
    sst_x <- sum((x - xbar)^2)

    # Redraw only the errors each replication; recompute the OLS slope.
    slopes <- replicate(reps, {
      y <- b0 + b1 * x + rnorm(n, 0, sd_u)
      sum((x - xbar) * (y - mean(y))) / sst_x
    })

    list(x = x, slopes = slopes, n = n, sd_u = sd_u,
         sst_x = sst_x, theo_se = sd_u / sqrt(sst_x), b0 = b0, b1 = b1)
  })

  # Left: one sample, true line vs fitted line
  output$scatter_plot <- renderPlot({
    s <- sim()
    x <- s$x
    y <- s$b0 + s$b1 * x + rnorm(s$n, 0, s$sd_u)
    b1hat <- sum((x - mean(x)) * (y - mean(y))) / sum((x - mean(x))^2)
    b0hat <- mean(y) - b1hat * mean(x)

    par(mar = c(4.5, 4.5, 3, 1))
    plot(x, y, pch = 19, col = "#6c8ebf80",
         main = "One sample", xlab = "X", ylab = "Y")
    abline(a = s$b0,  b = s$b1,  col = "#2c3e50", lwd = 2.5, lty = 2)
    abline(a = b0hat, b = b1hat, col = "#e74c3c", lwd = 2.5)
    legend("topleft", bty = "n", cex = 0.9,
           legend = c(paste0("True line (slope ", s$b1, ")"),
                      paste0("Fitted line (slope ", round(b1hat, 3), ")")),
           col = c("#2c3e50", "#e74c3c"), lwd = 2.5, lty = c(2, 1))
  })

  # Right: sampling distribution of the slope
  output$sampling_plot <- renderPlot({
    s <- sim()
    par(mar = c(4.5, 4.5, 3, 1))
    hist(s$slopes, breaks = 40, probability = TRUE,
         col = "#dae8fc", border = "#6c8ebf",
         main = paste0("Sampling Distribution of the Slope (n = ", s$n, ")"),
         xlab = expression(hat(beta)[1]), ylab = "Density")

    x_seq <- seq(min(s$slopes), max(s$slopes), length.out = 300)
    lines(x_seq, dnorm(x_seq, mean = s$b1, sd = s$theo_se),
          col = "#e74c3c", lwd = 2.5)
    abline(v = s$b1, lty = 2, lwd = 2, col = "#2c3e50")

    legend("topright", bty = "n", cex = 0.9,
           legend = c("Normal approximation", "True slope"),
           col    = c("#e74c3c", "#2c3e50"),
           lwd    = c(2.5, 2),
           lty    = c(1, 2))
  })

  output$stats_box <- renderUI({
    s <- sim()
    tags$div(class = "stats-box",
      HTML(paste0(
        "<b>True slope &beta;<sub>1</sub>:</b> ", s$b1, "<br>",
        "<b>Mean of &beta;&#770;<sub>1</sub>:</b> ", round(mean(s$slopes), 4),
        " <small>(unbiasedness)</small><br>",
        "<hr style='margin:8px 0'>",
        "<b>Theoretical SE (&sigma;<sub>u</sub>/&radic;SST<sub>x</sub>):</b> ",
        round(s$theo_se, 4), "<br>",
        "<b>Observed SD of &beta;&#770;<sub>1</sub>:</b> ",
        round(sd(s$slopes), 4)
      ))
    )
  })
}

shinyApp(ui, server)

Why \(X\) is held fixed. We redraw only the errors, which is what Wooldridge means by conditioning on the sample values \(X_1, \dots, X_n\). It keeps \(SST_x\) constant, so the theoretical width \(\sigma_u/\sqrt{SST_x}\) is exact. Changing \(n\) or the spread of \(X\) draws a fresh \(X\) design.

Three levers on the variance

The formula \(\mathrm{Var}(\hat{\beta}_1) = \sigma^2 / SST_x\) has exactly three moving parts, and each is a slider:

Slider Effect on the histogram Why
Error SD \(\sigma_u\) ↑ wider noisier errors (\(\sigma^2\) up) → harder to see the slope
Sample size \(n\) ↑ narrower more points → \(SST_x\) grows → more information
Spread of \(X\) ↑ narrower \(X\)’s more spread out → \(SST_x\) grows → line better pinned

The last one is the least obvious and the most useful: you learn a slope best when your \(X\)’s are far apart. Two points at the ends of the range pin a line down better than a hundred bunched in the middle.

Estimating \(\sigma^2\): where \(n-2\) comes from

The variance formula \(\mathrm{Var}(\hat{\beta}_1) = \sigma^2 / SST_x\) contains one piece you never observe: the error variance \(\sigma^2\). To report a standard error you estimate it from the residuals \(\hat{u}_i = y_i - \hat{\beta}_0 - \hat{\beta}_1 x_i\), using

\[s^2 = \frac{SSR}{n-2}, \qquad SSR = \sum_{i=1}^{n} \hat{u}_i^{\,2}.\]

Here \(SSR\) is the sum of squared residuals, the Residual line Stata prints and the same \(SSR\) from the \(R^2\) page: the total leftover misfit after the line has done its best. Its square root \(s\) is the residual standard error, roughly the typical vertical distance from a point to the fitted line, in the units of \(y\). And \(s^2\) is the estimate of \(\sigma^2\).

The one thing that needs explaining is the \(n-2\). It is the degrees of freedom, and it counts how many residuals are actually free to vary.

Compare the ordinary sample variance, which divides by \(n-1\). There you estimate one quantity first, the mean, and the deviations are forced to satisfy one constraint, \(\sum_i (y_i - \bar{y}) = 0\), leaving \(n-1\) free. In regression you estimate two quantities, \(\hat{\beta}_0\) and \(\hat{\beta}_1\), so the residuals satisfy two constraints, which come straight from the least-squares first-order conditions:

\[\sum_{i} \hat{u}_i = 0, \qquad \sum_{i} x_i \hat{u}_i = 0.\]

Once any \(n-2\) residuals are known, these two equations pin down the remaining two. Only \(n-2\) are free, so dividing \(SSR\) by \(n-2\) (rather than \(n\)) is what makes \(s^2\) an unbiased estimator of \(\sigma^2\). Dividing by \(n\) would understate it, because the fitted line is chosen to hug the data, so the residuals come out systematically a little smaller than the true errors \(u_i\).

A quick sanity check: with \(n=2\) the line passes through both points exactly, \(SSR = 0\), and \(s^2 = 0/0\) is undefined. Two points carry no information about the noise; you need at least three before you can estimate \(\sigma^2\) at all.

Putting the pieces together, the number printed next to the slope in any regression output is

\[\widehat{\mathrm{SE}}(\hat{\beta}_1) = \sqrt{\frac{s^2}{SST_x}} = \sqrt{\frac{SSR/(n-2)}{\sum_i (x_i - \bar{x})^2}},\]

the estimated counterpart of the theoretical width \(\sigma / \sqrt{SST_x}\) from the simulation above. This is the standard error you use for a \(t\)-test or a confidence interval on \(\hat{\beta}_1\).

Things to try

  • Defaults (\(n=250\), \(\sigma_u=6\), spread \(=3\)): the mean of \(\hat{\beta}_1\) lands on ~3.00, and the observed SD matches the theoretical SE almost exactly.
  • Crank \(\sigma_u\) to 12: the histogram roughly doubles in width. That is noise, not bias — it is still centered on 3.
  • Drop \(n\) to 20: wide and lumpy. The slope is barely identified from so few points.
  • Push spread of \(X\) from 1 to 8 at fixed \(n\): the histogram tightens sharply, even though nothing about the errors changed. That is \(SST_x\) at work.
  • Any setting: the mean of \(\hat{\beta}_1\) never leaves 3. Unbiasedness does not care about \(\sigma_u\), \(n\), or the spread of \(X\) — only about \(E(u\mid X) = 0\).

From here

  • The assumptions behind all of this — SLR.1–SLR.4 for unbiasedness, SLR.5 for the variance formula — are on the Gauss-Markov page.
  • The algebra that produces \(\hat{\beta}_1 = S_{xy}/S_{xx}\) in the first place is on The Algebra Behind OLS.
  • In practice \(\sigma_u\) is unknown. Replacing it with an estimate — and the \(n-2\) degrees of freedom that makes that estimate unbiased — is what turns the theoretical SE above into the standard error Stata prints next to your slope.