Deriving MSE($\hat{\beta}$) under Linear regression

Question

I was able to derive the MSE, but there's a part of the derivation which I don't really get. Here's what I got:

Facts:

• $\mathbb{E}(\hat{\beta})=\hat{\beta}\space$ (unbiased estimator)
• $\text{Cov}(\hat{\beta})= \sigma^2[(X^TX)^{-1}] $

By definition,

$$MSE = \mathbb{E}[||\hat{\beta}-\beta||^2] $$ $$= \space \mathbb{E}[(\hat{\beta}-\beta)^T(\hat{\beta}-\beta)]$$ Since $\hat{\beta}$ is unbiased, $$ \boldsymbol{= \text{tr}[\text{Cov}(\hat{\beta})]} *$$ $$= \text{tr}[\sigma^2(X^T X)^{-1}]$$ $\sigma^2$ is a scalar so it can be factored out, $$= \sigma^2 tr[(X^T X)^{-1}] $$

My confusion is line $*$, where I'm not sure how we got equation $ \text{tr}[\text{Cov}(\hat{\beta})] $. Here's what I understand so far:

Through Bias-Variance Decomposition, $$ \space \mathbb{E}[(\hat{\beta}-\beta)^T(\hat{\beta}-\beta)]=\mathbb{V}(\hat{\beta})+[\mathbb{E}(\hat{\beta})-\beta]^T[\mathbb{E}(\hat{\beta})-\beta]$$ Our estimator is unbiased so $\mathbb{E}(\hat{\beta})-\beta= 0$. Thus, $$\mathbb{E}[(\hat{\beta}-\beta)^T(\hat{\beta}-\beta)]=\mathbb{V}(\hat{\beta})$$

1. First, $\mathbb{V}(\hat{\beta})$ is supposed to be a scalar which doesn't really make sense to me, which brings me to my next question...
2. I assume that $V(\hat{\beta}) = \text{tr}[\text{Cov}(\hat{\beta})]$. Why is that?

Answer 1

This uses a decomposition of the expected value of the squared-norm

This MSE result is a particular application of a decomposition of the expected value of the squared-norm of a random vector. Start with the fact that the squared-norm of a vector is the sum of squares of its elements, so for any vector $\mathbf{Y} = (Y_1,...,Y_k)$ you have:

$$||\mathbf{Y}||^2 = \sum_{i=1}^k Y_i^2.$$

Taking the expectation of the squared norm then gives:$^\dagger$

$$\begin{align} \mathbb{E}(||\mathbf{Y}||^2) &= \mathbb{E} \bigg( \sum_{i=1}^k Y_i^2 \bigg) \\[6pt] &= \sum_{i=1}^k \mathbb{E}(Y_i^2) \\[6pt] &= \sum_{i=1}^k [\mathbb{E}(Y_i)^2 + \mathbb{V}(Y_i)] \\[6pt] &= \sum_{i=1}^k \mathbb{E}(Y_i)^2 + \sum_{i=1}^k \mathbb{V}(Y_i) \\[12pt] &= ||\mathbb{E}(\mathbf{Y})||^2 + \text{tr} (\mathbb{V}(\mathbf{Y})). \\[6pt] \end{align}$$

This gives you a general decomposition of the expected value of a squared-norm, which is split into the squared-norm of the expectation vector plus the trace of the variance matrix. (Note here that the variance matrix $\mathbb{V}(\mathbf{Y})$ is a $k \times k$ square matrix, which is only a scalar if the random vector contains only a single element.) The decomposition allows you to compute the expected value of the squared-norm of the random vector directly from the first two moments of that random vector.

The result you are seeing for the MSE of an unbiased estimator is just an application of this general result using the substitution $\mathbf{Y} = \hat{\boldsymbol{\beta}} - \boldsymbol{\beta}$. You have:

$$\begin{align} \text{MSE} &= \mathbb{E}(||\hat{\boldsymbol{\beta}} - \boldsymbol{\beta}||^2) \\[6pt] &= ||\mathbb{E}(\hat{\boldsymbol{\beta}} - \boldsymbol{\beta})||^2 + \text{tr} (\mathbb{V}(\hat{\boldsymbol{\beta}} - \boldsymbol{\beta})) \\[6pt] &= ||\mathbf{0}||^2 + \text{tr} (\mathbb{V}(\hat{\boldsymbol{\beta}})) \\[6pt] &= \text{tr} (\mathbb{V}(\hat{\boldsymbol{\beta}})) \\[6pt] &= \text{tr} (\sigma^2 (\mathbf{x}^\text{T} \mathbf{x})^{-1}) \\[6pt] &= \sigma^2 \text{tr} ((\mathbf{x}^\text{T} \mathbf{x})^{-1}) \\[6pt] &= \sigma^2 \sum_{i=1}^k [(\mathbf{x}^\text{T} \mathbf{x})^{-1}]_{i,i}. \\[6pt] \end{align}$$

As to your specific questions at the end: (1) the variance matrix $\mathbb{V}(\hat{\boldsymbol{\beta}})$ is a $k \times k$ square matrix with $k$ being the number of elements in the estimator $\hat{\boldsymbol{\beta}}$ (and it is only a scalar if $k=1$); (2) No, that equation is not true; the trace of the variance emerges from the above mathematics since we end up taking the sum of variances of the individual elements of the random vector.

---

$^\dagger$ I'm using slightly different notation to you: my $\mathbb{V}(\hat{\boldsymbol{\beta}})$ is your $\text{Cov}(\hat{\beta})$ (refering to the variance matrix, sometimes also called the covariance matrix).

Answer 2

When the parameter $\boldsymbol \theta\in\mathbb R^p, $ one is supposed to work with the matrix-valued squared error loss function $$\mathrm L(\boldsymbol \theta,\delta(\mathbf X))=(\delta(\mathbf X)- \boldsymbol \theta)(\delta(\mathbf X)- \boldsymbol \theta)^\top; $$ clearly the corresponding risk is $$\mathrm{MSE}(\boldsymbol \theta,\delta)=\mathsf E\left[(\delta(\mathbf X)- \boldsymbol \theta)(\delta(\mathbf X)- \boldsymbol \theta)^\top\right].$$

---

Denote by $\mathcal U(\boldsymbol\beta) $ the class of linear unbiased estimators of $\boldsymbol \beta$ in the linear model $\mathbf y=\mathbf X\boldsymbol \beta+\boldsymbol\varepsilon.$

If $g(\mathbf y)= z(\mathbf y) + \hat{\boldsymbol \beta},$ where $z(\mathbf y)\in\mathcal U(\boldsymbol 0)$ and $\hat{\boldsymbol\beta}\in\mathcal U(\boldsymbol\beta)$ is the least square estimator, then $g(\mathbf y) \in\mathcal U(\boldsymbol\beta) $ and it can be shown $$\mathsf{Cov}[g(\mathbf y) ]= \mathsf{Cov}\left(\hat{\boldsymbol \beta}\right) +\mathsf{Cov}[z(\mathbf y) ].$$

Using this identity, one can deduce for any $\tilde{\boldsymbol \beta}\in\mathcal U(\boldsymbol\beta)$ that $$\mathrm{MSE}\left(\boldsymbol\beta,\tilde{\boldsymbol\beta}\right)=\mathsf{Cov}\left(\tilde{\boldsymbol \beta}\right)=\sigma^2\mathbf B+\sigma^2\left(\mathbf X^\top\mathbf X\right)^{-1},$$ where $\bf B$ is a $p\times p$ symmetric nonnegative definite matrix.

In fact, $$\mathrm{MSE}\left(\boldsymbol\beta,\hat{\boldsymbol\beta}\right)=\mathsf{Cov}\left(\hat{\boldsymbol \beta}\right)=\sigma^2\left(\mathbf X^\top\mathbf X\right)^{-1}.$$

The real-valued unweighted squared error risk (OP's mean squared error) corresponding to the matrix-valued risk is given by $\rho\left(\boldsymbol\beta, \tilde{\boldsymbol\beta}\right) :=\operatorname{tr}\left[\mathrm{MSE}\left(\boldsymbol\beta,\tilde{\boldsymbol\beta}\right)\right];$ this would be for $\hat{\boldsymbol\beta},~\sigma^2\operatorname{tr}\left(\mathbf X^\top\mathbf X\right) ^{-1},$ which is, unsurprisingly, real-valued.

--

Reference:

Linear Regression, Jürgen Groß, Springer-Verlag, $2003, $ sec. $1.2.7, ~2.3.2.$