15  Classic least squares

References: Angela Montanari (2017), Chapter 3. Gardini A. (2000).

15.1 Assumptions

Let’s start from the classic assumptions for a linear model. The working hypothesis are:

1. The linear model approximate the conditional expectation, i.e. $\mathbb{E}\{y_i\mid {\mathbf{x}}_1,\dots ,{\mathbf{x}}_n\}=\mathbb{E}\{y_i\mid \mathbf{X}\}={\mathbf{x}}_i^{\mathrm{\top }}\mathbf{b}$ for $i=1,\dots n$.
2. The conditional variance of the response variable $y$ is constant, i.e. $\mathbb{V}\{y_i\mid {\mathbf{x}}_1,\dots ,{\mathbf{x}}_n\}=\mathbb{V}\{y_i\mid \mathbf{X}\}={\sigma }_e^2$ with $0<{\sigma }_e^2<\mathrm{\infty }$.
3. The conditional covariance of the response variable $y$ is zero, i.e. $\mathbb{C}v\{y_i,y_j\mid {\mathbf{x}}_1,\dots ,{\mathbf{x}}_n\}=\mathbb{C}v\{y_i,y_j\mid \mathbf{X}\}=0$ with $i\ne j$ and $i,j\in \{1,\dots ,n\}$.

Equivalently the formulation in terms of the stochastic component reads

1. $y_i={\mathbf{x}}_i^{\mathrm{\top }}\mathbf{b}+e_i$ for $i=1,\dots n$.
2. The residuals and the regressors are uncorrelated, i.e. $\mathbb{E}\{e_i\mid {\mathbf{x}}_1,\dots ,{\mathbf{x}}_n\}=\mathbb{E}\{e_i\mid \mathbf{X}\}=\mathbf{0}$.
3. The conditional variance of the residuals is constant, i.e. $\mathbb{V}\{e_i\mid {\mathbf{x}}_1,\dots ,{\mathbf{x}}_n\}=\mathbb{V}\{e_i\mid \mathbf{X}\}={\sigma }_e^2$ with $0<{\sigma }_e^2<\mathrm{\infty }$.
4. The conditional covariance of the residuals is zero, i.e. $\mathbb{C}v\{e_i,e_j\mid {\mathbf{x}}_1,\dots ,{\mathbf{x}}_n\}=\mathbb{C}v\{e_i,e_j\mid \mathbf{X}\}=0$ with $i\ne j$ and $i,j\in \{1,\dots ,n\}$.

Hence, in this setup the error terms $e$ are assumed to be independent and identically distributed with equal variance ${\sigma }_e$. Thus, the general expression of the covariance matrix in Equation 14.16 reduces to $\mathbf{\Sigma }={\sigma }_e^2{\mathbf{\text{I}}}_n$.

15.2 Estimator of $b$

Proposition 15.1 ($\mathbf{\text{Ordinary Least Squares (OLS) estimator}}$)
The ordinary least squares estimator (OLS) is the function ${\text{Q}}^{\text{OLS}}$ that minimize the sum of the squared residuals and return an estimate of the true parameter $\mathbf{b}$, i.e.  $$\begin{array}{}\text{(15.1)}& {\text{Q}}^{\text{OLS}}(\mathbf{b})=\mathbf{e}(\mathbf{b}{)}^{\mathrm{\top }}\mathbf{e}(\mathbf{b})\text{.}\end{array}$$ Formally, the OLS estimator is the solution of the following minimization problem, i.e.  $$\begin{array}{}\text{(15.2)}& {\mathbf{b}}^{\text{OLS}}=\underset{\mathbf{b}\in {\mathrm{\Theta }}_{\mathbf{b}}}{\text{argmin}}\{{\text{Q}}^{\text{OLS}}(\mathbf{b})\}\text{.}\end{array}$$ Notably, if $\mathbf{X}$ is non-singular one obtain an analytic expression, i.e.
$$\begin{array}{}\text{(15.3)}& {\mathbf{b}}^{\text{OLS}}=({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{y}\text{.}\end{array}$$ Equivalently, it is possible to express Equation 15.3 in terms of the covariance matrix of the $\mathbf{X}$ and $\mathbf{y}$, i.e.  $$\begin{array}{}\text{(15.4)}& {\mathbf{b}}^{\text{OLS}}=\mathbb{C}v\{\mathbf{X}{\}}^{-1}\mathbb{C}v\{\mathbf{X},\mathbf{Y}\}=\frac{\mathbb{C}v\{\mathbf{X},\mathbf{Y}\}}{\mathbb{V}\{\mathbf{X}\}}\text{.}\end{array}$$

Singularity of $\mathbf{X}$

Note that the solution is available if and only if $\mathbf{X}$ is non-singular. Hence, the columns should not be linearly dependent. In fact, one of the $k$-variables can be written as a linear combination of the others, then the determinant of the matrix ${\mathbf{X}}^{\mathrm{\top }}\mathbf{X}$ is zero and the inversion is not possible. Moreover, to have that $rank({\mathbf{X}}^{\mathrm{\top }}\mathbf{X})=k$ it is necessary that the number of observations have to be greater or equal than the number of regressors, i.e. $n\ge k$.

Proof: Proposition 15.1

Proof. Developing the product of the residuals in Equation 15.2: $$\begin{array}{}\text{(15.5)}& \begin{array}{rl}{\text{Q}}^{\text{OLS}}(\mathbf{b})& =\mathbf{e}(\mathbf{b}{)}^{\mathrm{\top }}\mathbf{e}(\mathbf{b})=\\ & =(\mathbf{y}-\mathbf{X}\mathbf{b}{)}^{\mathrm{\top }}(\mathbf{y}-\mathbf{X}\mathbf{b})=\\ & ={\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-(\mathbf{b}{)}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}\mathbf{y}-{\mathbf{y}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}+(\mathbf{b}{)}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}=\\ & ={\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-2(\mathbf{b}{)}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}\mathbf{y}+(\mathbf{b}{)}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}\end{array}\end{array}$$ To find the minimum, let’s compute the derivative of ${\text{Q}}^{\text{OLS}}$ with respect to $\mathbf{b}$, set it equal to zero and solve for $\mathbf{b}={\mathbf{b}}^{\text{OLS}}$, i.e. $$\begin{array}{rl}& \frac{d{\text{Q}}^{\text{OLS}}(\mathbf{b})}{d\mathbf{b}}=-2{\mathbf{X}}^{\mathrm{\top }}\mathbf{y}+2{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}=\mathbf{0}\\ & ⟹{\mathbf{X}}^{\mathrm{\top }}\mathbf{y}={\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}\\ & ⟹{\mathbf{b}}^{\text{OLS}}=({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{y}\end{array}$$ To establish if the above solution corresponds also to a global minimum, one must check the sign of the second derivative, i.e.  $$\frac{d^2{\text{Q}}^{\text{OLS}}(\mathbf{b})}{d\mathbf{b}d{\mathbf{b}}^{\mathrm{\top }}}=2{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}>0\text{,}$$ that in this case being always positive denotes a global minimum. An alternative derivation of this estimator, as in Equation 15.4, is obtained by substituting $$\frac{1}{n}\sum _{i=1}^n{\mathbf{x}}_i{\mathbf{x}}_i^{\mathrm{\top }}=\mathbb{V}\{\mathbf{X}\}\text{,}\frac{1}{n}\sum _{i=1}^n{\mathbf{x}}_i{y}_i^{\mathrm{\top }}=\mathbb{C}v\{\mathbf{X},\mathbf{Y}\}\text{,}$$ in Equation 15.3.

Intercept estimate

If in the data matrix $\mathbf{X}$ was included a column with ones, then the intercept parameter is obtained from Equation 15.3 or Equation 15.4. However, if it was not included, it is computed as: $${\alpha }^{\text{OLS}}=\mathbb{E}\{\mathbf{Y}\}-{\mathbf{b}}^{\text{OLS}}\mathbb{E}\{\mathbf{X}\}\text{.}$$

15.2.1 Projection matrices

Substituting the OLS solution (Equation 15.3) in Equation 14.12 we obtain the matrix $\mathbf{H}$, that project the vector $\mathbf{y}$ on the sub space of ${\mathbb{R}}^n$ generated by the matrix of the regressors $\mathbf{X}$, i.e.  $$\begin{array}{}\text{(15.6)}& \mathbf{H}=\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\text{,}\end{array}$$

The projection matrix $\mathbf{H}$ satisfies the following three properties, i.e. 

1. $\mathbf{H}$ is an $n\times n$ symmetric matrix.
2. $\mathbf{H}\mathbf{H}=\mathbf{H}$ is idempotent.
3. $\mathbf{H}\mathbf{X}=\mathbf{X}$.

Substituting the OLS solution (Equation 15.3) in the residuals (Equation 14.15) we obtain another projection matrix $\mathbf{M}$, that projects the vector $\mathbf{y}$ on the orthogonal sub-space with respect to the sub-space generated by the matrix of the regressors $\mathbf{X}$, i.e.  $$\begin{array}{}\text{(15.7)}& \mathbf{M}={\mathbf{\text{I}}}_n-\mathbf{H}\end{array}$$

The projection matrix $\mathbf{M}$ satisfies the following 3 properties, i.e. 

1. $\mathbf{M}$ is $n\times n$ and symmetric.
2. $\mathbf{M}\mathbf{M}=\mathbf{H}$ is idempotent.
3. $\mathbf{M}\mathbf{X}=\mathbf{0}$.

By definition $\mathbf{M}$ and $\mathbf{H}$ are orthogonal, i.e. $\mathbf{H}\mathbf{M}=\mathbf{0}$. Hence, the fitted values defined as $\hat{\mathbf{y}}=\mathbf{H}\mathbf{y}$ are the projection of the empiric values on the sub-space generated by $\mathbf{X}$. Symmetrically, the fitted residuals $\hat{\mathbf{e}}=\mathbf{M}\mathbf{y}$ are the projection of the empiric values on the sub-space orthogonal to the sub-space generated by $\mathbf{X}$.

Projection matrices

Proof. Let’s consider the property 2. of $\mathbf{H}$, i.e. $$\begin{array}{rl}\mathbf{H}\mathbf{H}& =(\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }})(\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }})=\\ & =({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}=\\ & =\mathbf{H}\end{array}$$ Let’s consider the property 3. of $\mathbf{H}$, i.e. $$\mathbf{H}\mathbf{X}=(\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }})\mathbf{X}=\mathbf{X}\text{.}$$ Let’s consider the property 2. of $\mathbf{M}$, i.e. $$\begin{array}{rl}\mathbf{M}\mathbf{M}& =({\mathbf{\text{I}}}_n-\mathbf{H})({\mathbf{\text{I}}}_n-\mathbf{H})=\\ & ={\mathbf{\text{I}}}_n-\mathbf{H}=\\ & =\mathbf{M}\end{array}$$ Let’s consider the property 3. of $\mathbf{M}$, i.e. $$\begin{array}{rl}\mathbf{M}\mathbf{X}& =({\mathbf{\text{I}}}_n-\mathbf{H})\mathbf{X}=\\ & =({\mathbf{\text{I}}}_n-\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }})\mathbf{X}=\\ & =\mathbf{X}-\mathbf{X}=\mathbf{0}\end{array}$$ Finally, let’s prove the orthogonality between $\mathbf{M}$ and $\mathbf{H}$, i.e.  $$\mathbf{H}\mathbf{M}=\mathbf{H}({\mathbf{\text{I}}}_n-\mathbf{H})=\mathbf{H}-\mathbf{H}=\mathbf{0}\text{.}$$

15.3 Properties OLS

Theorem 15.1 ($\mathbf{\text{Gauss-Markov theorem}}$)
Under the Gauss-Markov hypothesis the Ordinary Least Square (OLS) estimate is $\mathbf{\text{BLUE}}$ (Best Linear Unbiased Estimator), where “best” stands for the estimator with minimum variance in the class of linear unbiased estimators of the unknown true population parameter $\mathbf{b}$. More precisely, the Gauss-Markov hypothesis are:

1. $\mathbf{y}=\mathbf{X}\mathbf{b}+\mathbf{e}$.
2. $\mathbb{E}\{\mathbf{e}\}=0$.
3. $\mathbb{E}\{\mathbf{e}{\mathbf{e}}^{\mathrm{\top }}\}={\sigma }_e^2{\mathbf{\text{I}}}_n$, i.e. omoskedasticity.
4. $\mathbf{X}$ is non-stochastic and independent from the errors for all $n$’s.

Proposition 15.2 ($\mathbf{\text{Properties OLS estimator}}$)
1. Unbiased: ${\mathbf{b}}^{\text{OLS}}$ is correct and it’s conditional expectation is equal to true parameter in population, i.e.  $$\begin{array}{}\text{(15.8)}& \mathbb{E}\{{\mathbf{b}}^{\text{OLS}}\mid \mathbf{X}\}=\mathbf{b}\text{.}\end{array}$$ 2. Linear in the sense that it can be written as a linear combination of $\mathbf{y}$ and $\mathbf{X}$, i.e. ${\mathbf{b}}^{\text{OLS}}={\mathbf{A}}_{\text{x}}\mathbf{y}$, where ${\mathbf{A}}_{\text{x}}$ do not depend on $\mathbf{y}$, i.e.
$$\begin{array}{}\text{(15.9)}& {\mathbf{b}}^{\text{OLS}}={\mathbf{A}}_{\text{x}}\mathbf{y},{\mathbf{A}}_{\text{x}}=({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\text{.}\end{array}$$ 3. Under the Gauss-Markov hypothesis (Theorem 15.1) ${\mathbf{b}}^{\text{OLS}}$ is the estimator that has the minimum variance in the class of the unbiased linear estimators of $\mathbf{b}$ and it’s variance reads: $$\begin{array}{}\text{(15.10)}& \mathbb{V}\{{\mathbf{b}}^{\text{OLS}}\mid \mathbf{X}\}={\sigma }_e^2({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}\text{.}\end{array}$$

Proof: Proposition 15.2

Proof.

1. The OLS estimator is correct: it’s expected value is computed from Equation 15.3 and substituting Equation 14.12, is equal to the true parameter in population, i.e.
   $$\begin{array}{rl}\mathbb{E}\{{\mathbf{b}}^{\text{OLS}}\mid \mathbf{X}\}& =\mathbb{E}\{({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{y}\mid \mathbf{X}\}=\\ & =\mathbb{E}\{({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}(\mathbf{X}\mathbf{b}+\mathbf{e})\mid \mathbf{X}\}=\\ & =({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}+({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbb{E}\{\mathbf{e}\mid \mathbf{X}\}=\\ & =\mathbf{b}\end{array}$$

2. In general, applying the properties of the variance operator, the variance of ${\mathbf{b}}^{\text{OLS}}$ is computed as: $$\begin{array}{rl}\mathbb{V}\{{\mathbf{b}}^{\text{OLS}}\mid \mathbf{X}\}& =\mathbb{V}\{({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{y}\mid \mathbf{X}\}=\\ & =\mathbb{V}\{({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}(\mathbf{X}\mathbf{b}+\mathbf{e})\mid \mathbf{X}\}=\\ & =\mathbb{V}\{({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}+({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{e}\mid \mathbf{X}\}=\\ & =\mathbb{V}\{\mathbf{b}+({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{e}\mid \mathbf{X}\}=\\ & =\mathbb{V}\{({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{e}\mid \mathbf{X}\}\end{array}$$ Then, since $\mathbf{X}$ is non-stochastic one can bring it outside the variance thus obtaining: $$\begin{array}{}\text{(15.11)}& \begin{array}{rl}\mathbb{V}\{{\mathbf{b}}^{\text{OLS}}\mid \mathbf{X}\}& =({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbb{V}\{\mathbf{e}\mid \mathbf{X}\}\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}=\\ & =({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbb{E}\{\mathbf{e}{\mathbf{e}}^{\mathrm{\top }}\mid \mathbf{X}\}\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}\end{array}\end{array}$$ Under the Gauss Markov hypothesis (Theorem 15.1) the conditional variance $\mathbb{V}\{\mathbf{e}\mid \mathbf{X}\}=\sigma \cdot {\mathbf{\text{I}}}_n$ and therefore the Equation 15.11 reduces to: $$\begin{array}{rl}\mathbb{V}\{{\mathbf{b}}^{\text{OLS}}\mid \mathbf{X}\}& ={\sigma }_e^2({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}=\\ & ={\sigma }_e^2({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}\end{array}$$

15.4 Variance decomposition

In a linear model, the deviance (or total variance) of the dependent variable $\mathbf{y}$ can be decomposed into the sum of the regression variance and the dispersion variance. This decomposition helps us understand how much of the total variability in the data is explained by the model and how much is due to unexplained variability (residuals).

• Total Deviance ($\mathbb{D}ev\{\mathbf{y}\}$): represents the total variability of the dependent variable $\mathbf{y}$. It is calculated as the sum of the squared difference of $y_i$ from its mean $\overline{y}$.

• Regression Deviance ($\mathbb{D}ev\mathbb{R}eg\{\mathbf{y}\}$): represents the portion of variability that is explained by the regression model. It is computed as the sum of the squared differences between the fitted values ${\hat{y}}_i$ and $\overline{y}$.

• Dispersion Deviance ($\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}$): represents the portion of variability that is not explained by the model. It is computed as the sum of the squared differences between the observed values $y_i$ and the fitted values ${\hat{y}}_i$ (Equation 14.14).

Hence, the total deviance of $\mathbf{y}$ can be decomposed as follows: $$\begin{array}{}\text{(15.12)}& \begin{array}{rl}\mathbb{D}ev\{\mathbf{y}\}& =\mathbb{D}ev\mathbb{R}eg\{\mathbf{y}\}+\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}\\ \sum _{i=1}^n(y_i-\overline{y}{)}^2& =\sum _{i=1}^n({\hat{y}}_i-\overline{y}{)}^2+\sum _{i=1}^n({\hat{y}}_i-y_i{)}^2\\ {\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-n{\overline{y}}^2& ={\mathbf{b}}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}-n{\overline{y}}^2+{\mathbf{e}}^{\mathrm{\top }}\mathbf{e}\end{array}\end{array}$$

Regression deviance

Proof. Let’s prove the expression for the regression deviance $\mathbb{D}ev\mathbb{R}eg\{\mathbf{y}\}$, i.e.  $$\begin{array}{rl}\mathbb{D}ev\mathbb{R}eg\{\mathbf{y}\}& =\mathbb{D}ev\{\mathbf{y}\}-\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}=\\ & ={\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-n{\overline{y}}^2-{\mathbf{e}}^{\mathrm{\top }}\mathbf{e}=\\ & ={\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-n{\overline{y}}^2-{(\mathbf{y}-\mathbf{X}\mathbf{b})}^{\mathrm{\top }}(\mathbf{y}-\mathbf{X}\mathbf{b})=\\ & ={\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-n{\overline{y}}^2+{\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-{\mathbf{y}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}-\mathbf{y}{\mathbf{b}}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}+{\mathbf{b}}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}=\\ & =2{\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-n{\overline{y}}^2-2{\mathbf{y}}^{\mathrm{\top }}(\mathbf{X}\mathbf{b})+{\mathbf{b}}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}=\\ & ={\mathbf{b}}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}-n{\overline{y}}^2\end{array}$$

15.4.1 Estimator of ${\sigma }_e^2$

The OLS estimator do not depend on variance of the residuals ${\sigma }_e^2$ and it is not possible to obtain in one step both the estimators. As far as we know ${\sigma }_e^2$ is the variance of the residuals of which we know the realized values on the sample $\hat{\mathbf{e}}=\{{\hat{e}}_1,{\hat{e}}_2,\dots ,{\hat{e}}_n\}$. Hence, let’s define an unbiased estimator of the population variance ${\sigma }_e^2$ as: $${\hat{s}}_e^2=\frac{{\hat{\mathbf{e}}}^{\mathrm{\top }}\hat{\mathbf{e}}}{n-k-1}\text{.}$$ In general the regression variance overestimate the true variance ${\sigma }_e^2$, i.e. $${\hat{s}}_r^2=k{\sigma }_e^2+g(\mathbf{b},\mathbf{X}),g(\mathbf{b},\mathbf{X})\ge 0\text{.}$$ Only in the special case where $b_1=b_2=\cdots =b_k=0$ in population, then $g(\mathbf{b},\mathbf{X})=0$ and also the regression variance produces a correct estimate of ${\sigma }_e^2$.

Correct estimator of ${\sigma }_e^2$

Proof. By definition, the residuals can be computed pre multiplying the matrix $\mathbf{M}$ to $\mathbf{y}$, i.e.  $$\begin{array}{rl}\hat{\mathbf{e}}& =\mathbf{y}-\hat{\mathbf{y}}=\\ & =\mathbf{y}-\mathbf{X}{\mathbf{b}}^{\text{OLS}}=\\ & =\mathbf{y}-\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{y}=\\ & =({\mathbf{\text{I}}}_n-\mathbf{H})\mathbf{y}=\\ & =\mathbf{M}\mathbf{y}\end{array}$$ Substituting $\mathbf{y}=\mathbf{X}{\mathbf{b}}^{\text{OLS}}+\mathbf{e}$, since $\mathbf{M}\mathbf{X}=0$, one obtain $$\begin{array}{rl}\hat{\mathbf{e}}& =\mathbf{M}(\mathbf{X}{\mathbf{b}}^{\text{OLS}}+\mathbf{e})=\\ & =\mathbf{M}\mathbf{X}{\mathbf{b}}^{\text{OLS}}+\mathbf{M}\mathbf{e}=\\ & =\mathbf{M}\mathbf{e}\end{array}$$ Then, being $\mathbf{M}$ symmetric and idempotent: $$\begin{array}{rl}{\hat{\mathbf{e}}}^{\mathrm{\top }}\hat{\mathbf{e}}& =(\mathbf{M}\mathbf{e}{)}^{\mathrm{\top }}(\mathbf{M}\mathbf{e})=\\ & ={\mathbf{e}}^{\mathrm{\top }}{\mathbf{M}}^{\mathrm{\top }}\mathbf{M}\mathbf{e}=\\ & ={\mathbf{e}}^{\mathrm{\top }}\mathbf{M}\mathbf{e}\end{array}$$ Thus, since ${\hat{\mathbf{e}}}^{\mathrm{\top }}\mathbf{M}\hat{\mathbf{e}}$ is a scalar, the expected value of the deviance of dispersion read $$\begin{array}{rl}\mathbb{E}\{{\hat{\mathbf{e}}}^{\mathrm{\top }}\hat{\mathbf{e}}\}& =\mathbb{E}\{{\mathbf{e}}^{\mathrm{\top }}\mathbf{M}\mathbf{e}\}=\\ & =\mathbb{E}\{\text{trace}({\mathbf{e}}^{\mathrm{\top }}\mathbf{M}\mathbf{e})\}=\\ & =\mathbb{E}\{\text{trace}(\mathbf{M}\mathbf{e}{\mathbf{e}}^{\mathrm{\top }})\}=\\ & =\text{trace}(\mathbf{M}\mathbb{E}\{\mathbf{e}{\mathbf{e}}^{\mathrm{\top }}\})=\\ & =\mathbb{E}\{\mathbf{e}{\mathbf{e}}^{\mathrm{\top }}\}\cdot \text{trace}(\mathbf{M}{\mathbf{\text{I}}}_n)=\\ & ={\sigma }_e^2\cdot \text{trace}(\mathbf{M}{\mathbf{I}}_n)=\\ & ={\sigma }_e^2\cdot \text{trace}(\mathbf{M})\end{array}$$ where the trace of the matrix $\mathbf{M}$ is: $$\begin{array}{rl}\text{trace}(\mathbf{M})& =\text{trace}({\mathbf{\text{I}}}_n-\mathbf{H})=\\ & =\text{trace}({\mathbf{\text{I}}}_n)-\text{trace}(\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\top }})=\\ & =\text{trace}({\mathbf{\text{I}}}_n)-\text{trace}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1})=\\ & =\text{trace}({\mathbf{\text{I}}}_n)-\text{trace}({\mathbf{\text{J}}}_{k+1})=\\ & =n-k-1\end{array}$$ Hence the expectation of the deviance of dispersion is equal to $$\mathbb{E}\{\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}\}={\sigma }_e^2\cdot (n-k-1)⟹{\hat{s}}_e^2=\frac{{\hat{\mathbf{e}}}^{\mathrm{\top }}\hat{\mathbf{e}}}{n-k-1}$$

The decomposition of the deviance of $\mathbf{y}$ holds true also with respect to the correspondents degrees of freedom,

Table 15.1: Deviance and variance decomposition in a multivariate linear model

Deviance	Degrees of freedom	Variance
$\mathbb{D}ev\{\mathbf{y}\}=\sum _{i=1}^n(y_i-\overline{y}{)}^2$	$n-1$	${\hat{s}}_y^2=\frac{\mathbb{D}ev\{\mathbf{y}\}}{n-1}$
$\mathbb{D}ev\mathbb{R}eg\{\mathbf{y}\}=\sum _{i=1}^n({\hat{y}}_i-\overline{y}{)}^2$	$k$	${\hat{s}}_r^2=\frac{\mathbb{D}ev\mathbb{R}eg\{\mathbf{y}\}}{k}$
$\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}=\sum _{i=1}^n({\hat{y}}_i-y_i{)}^2$	$n-k-1$	${\hat{s}}_e^2=\frac{\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}}{n-k-1}$

15.5 ${\text{R}}^2$

The ${\text{R}}^2$ statistic, also known as the coefficient of determination, is a measure used to assess the goodness of fit of a regression model. In a multivariate context, it evaluates how well the independent variables explain the variability of the dependent variable.

Definition 15.1 ($\mathbf{\text{Multivariate }}{\text{R}}^2$)
The $R^2$ represents the proportion of the variation in the dependent variable that is explained or predicted by the independent variables. Formally, it is defined as the ratio of the deviance explained by the model ($\mathbb{D}ev\mathbb{R}eg\{\mathbf{y}\}$) to the total deviance ($\mathbb{D}ev\{\mathbf{y}\}$). It can also be expressed as one minus the ratio of the residual deviance ($\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}$) to the total deviance, i.e.  $$\begin{array}{}\text{(15.13)}& {\text{R}}^2=\frac{\mathbb{D}ev\mathbb{R}eg\{\mathbf{y}\}}{\mathbb{D}ev\{\mathbf{y}\}}=1-\frac{\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}}{\mathbb{D}ev\{\mathbf{y}\}}\text{.}\end{array}$$ Using the variance decomposition (Equation 15.12), it is possible to write a multivariate version of the ${\text{R}}^2$ as: $${\text{R}}^2=\frac{{\mathbf{b}}^{\mathrm{\top }}{\mathbf{X}}^{\mathrm{\top }}\mathbf{X}\mathbf{b}-n{\overline{y}}^2}{{\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-n{\overline{y}}^2}=1-\frac{{\mathbf{e}}^{\mathrm{\top }}\mathbf{e}}{{\mathbf{y}}^{\mathrm{\top }}\mathbf{y}-n{\overline{y}}^2}\text{.}$$

The numerator represents the variance explained by the regression model, while the denominator the total variance in the dependent variable. The term ${\mathbf{e}}^{\mathrm{\top }}\mathbf{e}$ in the second expression represents the variance of the residuals, or the variance not explained by the model. A value of the ${\text{R}}^2$ close to 1 denotes that a large proportion of the variability of the dependent variable has been explained by the regression model, while a value close to 0 indicates that the model explains very little of the variability.

Variance Inflation Factor (VIF)

The elements on the diagonal of the matrix $({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}$ determine the variances while the other elements the covariances. In general the variance of the $j$-th regressor is denoted as $\mathbb{V}\{b_j\}={\sigma }_e^2{c}_{jj}$ where $c_{jj}$ is the $j$-th element on the diagonal of $({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}$. An alternative expression for the variance is: $$\mathbb{V}\{b_j\}=\frac{{\sigma }_e^2}{\mathbb{D}ev\{X_j\}}\frac{1}{1-R_{j0}^2}\text{,}$$ where $R_{j0}^2$ is the multivariate coefficient of determination on the regression of $X_j$ on the other regressors. The term $\frac{1}{1-R_{j0}^2}$ is also called Variance Inflation Factor and denoted as ${\text{VIF}}_j$.

Limitations of ${\text{R}}^2$

The ${\text{R}}^2$ statistic has some limitations. Firstly, it can be close to 1 even if the relationship between the variables is not linear. Additionally, ${\text{R}}^2$ increases whenever a new regressor is added to the model, making it unsuitable for comparing models with different numbers of regressors.

Definition 15.2 A more robust indicator that does not always increase with the addition of a new regressor is the adjusted ${\text{R}}^2$, which is computed as: $${\overline{\text{R}}}^2=1-\frac{n-1}{n-k-1}\frac{\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}}{\mathbb{D}ev\{\mathbf{y}\}}=1-\frac{{\hat{s}}_e^2}{{\hat{s}}_y^2}\text{.}$$ The ${\overline{\text{R}}}^2$ can be negative, and its value will always be less than or equal to that of ${\text{R}}^2$. Unlike ${\text{R}}^2$, the adjusted version increases only when the new explanatory variable improves the model more than would be expected simply by adding another variable.

${\text{R}}^2$

Proof. To arrive at the formulation of the adjusted ${\text{R}}^2$ let’s consider that under the null hypothesis ${\mathcal{H}}_0:b_1=b_2=\cdots =b_k$ the variance of regression ${\hat{s}}_r^2$ (Table 15.1) is a correct estimate of the variance of the residuals ${\sigma }_e^2$. Hence, under ${\mathcal{H}}_0$: $$\frac{n-1}{k}\mathbb{E}\left\{\frac{\mathbb{D}ev\mathbb{R}eg\{\mathbf{y}\}}{\mathbb{D}ev\{\mathbf{y}\}}\right\}\stackrel{\sim }{=}1\text{.}$$ This implies that the expectation of the $R^2$ is not zero (as it should be under ${\mathcal{H}}_0$) but: $$\mathbb{E}\{{\text{R}}^2\}\stackrel{\sim }{=}\frac{k}{n-1}\text{.}$$ Let’s rescale the ${\text{R}}^2$ such that when ${\mathcal{H}}_0$ holds true it is equal to zero, i.e.  $${\text{R}}_c^2={\text{R}}^2-\frac{k}{n-1}\text{.}$$ However, the specification of ${\text{R}}_c^2$ implies that when ${\text{R}}^2=1$ (perfect linear relation between $\mathbf{X}$ and $\mathbf{y}$) the value of ${\text{R}}_c^2<1$, i.e. ${\text{R}}_c^2=\frac{n-k-1}{n-1}<1$. Hence, let’s correct again the indicator such that it takes values in $[0,1]$, i.e.  $$\begin{array}{rl}{\overline{\text{R}}}^2& =({\text{R}}^2-\frac{k}{n-1})\frac{n-1}{n-k-1}=\\ & =\left(\frac{{\text{R}}^2(n-1)-k}{n-1}\right)\frac{n-1}{n-k-1}=\\ & =\frac{n-1}{n-k-1}{\text{R}}^2-\frac{k}{n-k-1}\end{array}$$ Remembering that ${\text{R}}^2$ can be rewritten as in Equation 15.13 one obtain: $$\begin{array}{rl}{\overline{\text{R}}}^2& =\frac{n-1}{n-k-1}(1-\frac{\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}}{\mathbb{D}ev\{\mathbf{y}\}})-\frac{k}{n-k-1}=\\ & =\frac{(n-1)\mathbb{D}ev\{\mathbf{y}\}-(n-1)\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}}{\mathbb{D}ev\{\mathbf{y}\}(n-k-1)}-\frac{k}{n-k-1}=\\ & =\frac{n-1}{n-k-1}-\frac{n-1}{n-k-1}\frac{\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}}{\mathbb{D}ev\{\mathbf{y}\}}-\frac{k}{n-k-1}=\\ & =1-\frac{n-1}{n-k-1}\frac{\mathbb{D}ev\mathbb{D}isp\{\mathbf{y}\}}{\mathbb{D}ev\{\mathbf{y}\}}=\\ & =1-\frac{{\hat{s}}_e^2}{{\hat{s}}_y^2}\end{array}$$

15.6 Diagnostic

Let’s consider a linear model where the residuals $\mathbf{e}$ are IID normally distributed random variables. Hence, the working hypothesis of the Gauss Markov theorem holds true.

15.6.1 $\text{t}$-test for $b_j$

A $\text{t}$-test valuates the significance of the parameter of a regressor, given the effect of the others $k-1$ regressors, by testing the null hypothesis of linear independence between $\mathbf{y}$ and ${\mathbf{X}}_j$, i.e.  $${\mathcal{H}}_0:b_j=0\text{.}$$ Under the normality assumption on the distribution of the residuals, the vector of parameter $\hat{\mathbf{b}}$ is distributed as a multivariate normal random vector, thus also the marginal distribution of ${\hat{b}}_j$ is normal. Therefore, given the expectation and variance of ${\hat{b}}_j$ one can standardize it to obtain $$\frac{{\hat{b}}_j-\mathbb{E}\{{\hat{b}}_j\}}{\sqrt{\mathbb{V}\{{\hat{b}}_j\}}}=\frac{{\hat{b}}_j-b_j}{\sqrt{{\sigma }_e^2\cdot c_{jj}}}\sim \mathcal{N}(0,1)\text{.}$$ Substituting the unknown ${\sigma }_e^2$ with it’s correct estimate ${\hat{\sigma }}_e^2$ one obtain the statistic, i.e. $$\begin{array}{}\text{(15.14)}& {\text{t}}_j=\frac{{\hat{b}}_j-b_j}{\sqrt{{\hat{\sigma }}_e^2\cdot c_{jj}}}\sim t_{n-k-1}\text{.}\end{array}$$ that is Student-$\text{t}$ distributed (Equation 35.2) with $\nu =n-k-1$ degrees of freedom. Under the null hypothesis ${\mathcal{H}}_0:{\hat{b}}_j=0$ one obtain the $\text{t}$-test statistic, i.e. $$\begin{array}{}\text{(15.15)}& {\text{t}}_j\stackrel{{\mathcal{H}}_0}{=}\frac{{\hat{b}}_j}{\sqrt{{\hat{\sigma }}_e^2\cdot c_{jj}}}\sim t_{n-k-1}\text{.}\end{array}$$

15.6.2 Confidence intervals for $b$

Under the assumption of normality, from Equation 15.14, one can build a confidence interval for $b_j$, i.e.  $$b_j={\hat{b}}_j\pm {\text{t}}_{n-k-1}^{\alpha /2}\sqrt{{\sigma }_e^2\cdot c_{jj}}\text{.}$$ where $\alpha$ is the confidence level, ${\text{t}}_{n-k-1}^{\alpha /2}$ is the quantile at level $\alpha /2$ of a Student-t distribution with $n-k-1$ and $c_{jj}$ is the $j$-th element on the diagonal of $({\mathbf{X}}^{\mathrm{\top }}\mathbf{X}{)}^{-1}$.

15.6.3 F-test for the regression

The ${\text{F}}_{\text{test}}$ evaluates the significance of the entire regression model by testing the null hypothesis of linear independence between $\mathbf{y}$ and $\mathbf{X}$, i.e.  $${\mathcal{H}}_0:b_1=b_2=\cdots =b_k=0\text{.}$$ where the only coefficient different from zero is the intercept. In this case, the test statistic reads $${\text{F}}_{\text{test}}=\frac{{\hat{s}}_r^2}{{\hat{s}}_e^2}=\frac{\mathbb{D}ev\mathbb{R}eg(\mathbf{y})\cdot (n-k-1)}{k\cdot \mathbb{D}ev\mathbb{D}isp(\mathbf{y})}\sim {\text{F}}_{k,n-k-1}\text{,}$$ that is distributed with an $\text{F}$-Fischer (Equation 35.3) with ${\nu }_1=k$ and ${\nu }_2=n-k-1$ degrees of freedom. ${\hat{s}}_r^2$ is the regression variance and ${\hat{s}}_e^2$ is the dispersion variance. By fixing a significance level $\alpha$, the null hypothesis ${\mathcal{H}}_0$ is rejected if ${\text{F}}_{\text{test}}>{\text{F}}_{k,n-k-1}^{\alpha }$. Remembering the relation between the deviance and the ${\text{R}}^2$, i.e. $\mathbb{D}ev\mathbb{R}eg(\mathbf{y})={\text{R}}^2\mathbb{D}ev(\mathbf{y})$ and $\mathbb{D}ev\mathbb{D}isp(\mathbf{y})=(1-{\text{R}}^2)\mathbb{D}ev(\mathbf{y})$, it is possible to express the $\text{F}$-test in terms of the multivariate ${\text{R}}^2$ as: $${\text{F}}_{\text{test}}=\frac{{\text{R}}^2}{1-{\text{R}}^2}\frac{n-k-1}{k}\sim {\text{F}}_{k,n-k-1}\text{.}$$

Interpretation $\text{F}$-test

If the null hypothesis $H_0$ is rejected then:

• The variability of $Y$ explained by the model is significantly greater than the residual variability.
• At least one of the $k$ regressors has a coefficient $b_k$ that is significantly different from zero in the population.

On contrary if ${\mathcal{H}}_0$ is not rejected, then the model is not adequate and there is no evidence of a linear relation between $\mathbf{y}$ and $\mathbf{X}$.

Angela Montanari. 2017. “Appunti Sui Modelli Lineari.”

Gardini A., Costa M., Cavaliere G. 2000. Econometria, Volume Primo. FrancoAngeli.