# How can I compute correct standard errors after implementing the FWL Theorem?

I am trying to implement the FWL theorem for some sample data in Stata. This theorem tells us that given a multivariate regression of the form $$y = \beta_{1}x_{1} + \beta_{2}x_{2} + \varepsilon$$, the OLS estimator obtained by regressing $$y$$ on $$\gamma_1$$ where $$\gamma_{1}$$ are the residuals from the regression of $$x_{1}$$ on $$x_{2}$$ will be the same as the OLS estimator on $$x_{1}$$ from the multivariate regression.

I have been able to produce this result for the coefficient, but how can I find the correct standard errors for the OLS estimator from the partialled out regression? I understand that the issue has to do with the degrees of freedom when computing the new standard errors and I found a useful answer here, but this example requires that the independent and the dependent variables are both partialled out. Is there a way to find the correct standard errors when you only partial out the independent variables? In case it is useful, my Stata code to compute the standard errors in question is below:

sysuse auto2, clear

*compute the multivariate regression

*partialling out
predict double resid_x2, resid

*the coefficient is the same in this case as in the multivariate regression, but the standard error is different
regress price resid_x2

*correcting for the degrees of freedom, this line should be equal to the standard
*error in the multivariate regression
display sqrt(73/71)*_se[resid_x2]

$$$$


The "problem" here is that when you only partial out the explanatory variables, it is - unlike in the case of the "full" FWL approach where, as you showed, only the degrees of freedom need to be adjusted - no longer true that the residuals of the partial FWL approach are identical to those of the full regression.

That the residuals for full partialing out are identical to those for full OLS can be seen as follows. Let $$e$$ denote the residuals of a regression of $$y$$ on all regressors, denoted $$X$$, i.e. $$y=X_1b_1+X_2b_2+e$$ Then, multiplying through with $$M_{X_2}$$ yields the regression of the residuals on each other, $$M_{X_2}y=M_{X_2}X_1b_1+M_{X_2}e,$$ where we have used $$M_{X_2}X_2=0$$ By FWL, $$b_1$$ is the estimated coefficient vector of this regression, so that $$M_{X_2}e$$ are the residuals. Now, by linear regression properties $$X_2'e=0$$, thus $$M_{X_2}e=e.$$

On the contrary, the partial approach "does not multiply through", so that equality of residuals does not result: $$M_{X_2}y-M_{X_2}X_1b_1\neq y-M_{X_2}X_1b_1$$

Hence, the "numerators" $$\hat\sigma^2$$ of the standard error formula $$\sqrt{diag(\hat\sigma^2(X_1'M_{X_2}X_1)^{-1})}$$ also no longer agree, as $$\hat\sigma^2$$ is, of course, computed as the variance of the residuals.

So to make standard errros agree, you would need to compute your estimate of $$\sigma^2$$ from the residuals of the full regression (which of course then raises the practical question of whether the partial approach still is useful).

Numerical illustration (in R - as my answer hopefully illustrates, the issue is not related to the software used):

n <- 20
x1 <- rnorm(n)
x2 <- rnorm(n)
y <- rnorm(n)

(ls <- summary(lm(y~x1+x2)))                                               # full LS
(fwl <- summary(lm(resid(lm(y~x2))~resid(lm(x1~x2)))))                     # full FWL
(partial.fwl <- summary(lm(y~resid(lm(x1~x2)))))                           # partial FWL

sqrt(n-2)/sqrt(n-3)*partial.fwl$$coefficients[2,2] # adjusted s.e.s from full FWL... all.equal(sqrt(n-2)/sqrt(n-3)*fwl$$coefficients[2,2], ls$coefficients[2,2]) # ...agree with those from full LS cbind(resid(ls), resid(fwl), resid(partial.fwl)) # residuals of partial do not agree with the other two partial.fwl$$coefficients[2,2]*ls$$sigma/partial.fwl$$sigma # adjusted s.e.s from partial FWL... all.equal(partial.fwl$$coefficients[2,2]*ls$$sigma/partial.fwl$$sigma, ls$coefficients[2,2])                                            # agree with full OLS
`