Mean adjusted $R^2$ for linear regression with gaussian noise covariates

Question

Consider the simple regression model $y = a x + b + \varepsilon$, where $x$ is a covariate, $y$ is the observed response and $\varepsilon$ is the unobserved noise (with no distribution assumption). We can add to the model a covariate $\omega$ which is a realisation of a gaussian noise, and this obviously increases the $R^2$ coefficient of the fit.

However, is the expected adjusted $R^2$ cofficient of the augmented model the same as the adjusted $R^2$ coefficient of the original model?

This would mean that the adjusted $R^2$ is apt at eliminating the inflation of $R^2$ which is due to random noise in the covariates. The answer may be completely trivial. In any case, I'm interested to learn more about this and in references (for a mathematically mature audience) on the subject.

---

Remark. An example seems to indicate it may be the case. With R's Ozone data, the adjusted $R^2$ coefficient of the model $\texttt{O3}$ ~ $\texttt{T12}$ is found as follows:

out = with(ozone, lm(O3 ~ T12))
summary(out)$adj.r.squared
## 0.2640621

For the model $\texttt{O3}$ ~ $\texttt{T12}$ + $\texttt{noise}$, where $\texttt{noise}$ is an instance of white noise, the average adjusted $R^2$ is almost the same:

r2.aug <- function() {
    df = with(ozone, data.frame(O3, T12, noise=rnorm(50)))
    summary(lm(df$O3 ~ df$T12 + df$noise))$adj.r.squared
}
set.seed(1)
mean(replicate(10000, r2.aug()))
## 0.264054

Answer 1

Suppose we have a general multiple linear regression model (without assuming normality):

$$\begin{matrix} \boldsymbol{Y} = \boldsymbol{X}\boldsymbol{\beta} + \boldsymbol{\epsilon} & & \mathbb{E}(\boldsymbol{\epsilon}) = 0 & & \mathbb{V}(\boldsymbol{\epsilon}) = \sigma^2 \boldsymbol{I}. \end{matrix}$$

We assume that the design matrix $\boldsymbol{X}$ has an intercept term (column of ones) and $m$ other columns and is of full rank (i.e., all the regressors are linearly independent). It is a matrix with dimensions $n \times (m+1)$. Using ordinary least-squares (OLS) estimation of the coefficient parameters, the sums-of-squares can be written in terms of the annihilator matrix $\boldsymbol{A} (\boldsymbol{X}) \equiv \boldsymbol{I} - \boldsymbol{X} (\boldsymbol{X}^\text{T} \boldsymbol{X})^{-1} \boldsymbol{X}^\text{T}$ and the matrix $\boldsymbol{M} \equiv I - \frac{1}{n} \boldsymbol{1}_{n \times n}$. Specifically, it can be shown that the residual and total mean-squares are given respectively by:

$$\begin{matrix} \hat{\sigma}^2 \equiv \text{MSE}(\boldsymbol{Y, \boldsymbol{X}}) = \frac{1}{n-m-1} \boldsymbol{Y}^\text{T} \boldsymbol{A}(\boldsymbol{X}) \boldsymbol{Y} & &\hat{\sigma}_Y^2 \equiv \text{MST}(\boldsymbol{Y}) = \frac{1}{n-1} \boldsymbol{Y}^\text{T} \boldsymbol{M} \boldsymbol{Y}. \end{matrix}$$

Using the general formula for the expected value of a quadratic form we then obtain:

$$\begin{matrix} \mathbb{E}(\hat{\sigma}^2|\boldsymbol{X}) = \sigma^2 & & & & & \mathbb{E}(\hat{\sigma}_Y^2)|\boldsymbol{X}) = \boldsymbol{\beta}^\text{T} \hat{\Sigma}_X \boldsymbol{\beta} + \sigma^2, \end{matrix}$$

where $\hat{\Sigma}_X \equiv \frac{1}{n-1} \boldsymbol{X}^\text{T} \boldsymbol{M} \boldsymbol{X}$ is the sample covariance of the regressors. Now, for a model with an intercept term, the adjusted coefficient of determination is:

$$R^2_{\text{adj}} = 1 - \frac{\hat{\sigma}^2}{\hat{\sigma}_Y^2}.$$

Rearranging and taking expectations of both sides gives the general result:

$$\mathbb{E}(\hat{\sigma}_Y^2\cdot R^2_{\text{adj}} | \boldsymbol{X}) = \boldsymbol{\beta}^\text{T} \hat{\Sigma}_X \boldsymbol{\beta} .$$

In this sense, the adjusted coefficient of determination is calibrated to adjust for the number of regressors in the model.

---

Partial application to your problem: Suppose you compare the general model above to a model where you add an additional regressor to the design matrix, but (by assumption) this additional regressor has no relationship with the response vector. In both cases you have the same model equation $\boldsymbol{Y} = \boldsymbol{X}\boldsymbol{\beta} + \boldsymbol{\epsilon}$, but in the latter case your OLS estimate is different, owing to the extra regressor in the design matrix (which flows through to affect the annihilator matrix).

Since both models have the same underlying model equation we clearly have:

$$\mathbb{E}(\hat{\sigma}_Y^2 \cdot R^2_{\text{adj}, 1}) = \mathbb{E}(\hat{\sigma}_Y^2 \cdot R^2_{\text{adj}, 2}) ,$$

where the subscripts indicate the two models. This gets you a rough equivalency result, but to go further you would need to see what happens when you partition the annihilator matrix in the second regression into parts attributable to the original regressors, and a part attributable to the newly introduced noise term. If you can decompose the adjusted coefficient in this way, you might be able to go further than this.