Exact distribution of the F-statistic under heteroskedasticity of unknown form for improved inference

Abstract

In this paper, we derive the exact finite sample distribution of the F ($=t^2$) statistic for a single linear restriction on the regression parameters. We show that the F statistic can be expressed as a ratio of quadratic forms, and therefore its exact cumulative distribution under the null hypothesis can be derived from the result of Imhof [Computing the distribution of quadratic forms in normal variables. Biometrika. 1961;48(3/4):419–426]. A numerical calculation is carried out for the exact distribution of the F statistic using various HC covariance matrix estimators, and the rejection probability under the null hypothesis (size) based on the exact distribution is examined. The results show the exact finite sample distribution is remarkably reliable, while, in comparison, the use of the F-table leads to a serious over-rejection when the sample is not large or leveraged/unbalanced.

Keywords:

• Heteroskedasticity
• finite-sample theory
• Imhof distribution

JEL Classifications:

• C1
• C2

Acknowledgements

The authors are grateful for the useful comments from the editor and three anonymous referees.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Notes

Notes: As noted in Remark  3.1, the Infeasible Imhof distribution of the F statistic is $\mathrm{Imhof}(\hat{\mathbf{\Omega }},\mathbf{N}\mathbf{\Omega },c)$ obtained from the cumulative distribution function in Equation (15(15) $Pr(F\le c)=\frac{1}{2}-\frac{1}{\pi }{\int }_0^{\mathrm{\infty }}\frac{\sin (\frac{1}{2}\sum _{j=1}^J[{\tan }^{-1}({\lambda }_{j}u)])}{u\prod _{j=1}^J(1+{\lambda }_j^2{u}^2{)}^{\frac{1}{4}}}\mathrm{d}u\equiv \mathrm{Imhof}(\hat{\mathbf{\Omega }},\mathbf{N}\mathbf{\Omega },c).$(15) ). The Benchmark distribution is computed from the cumulative distribution function in Equation (21(21) $Pr(F^{\ast }\le c)=\frac{1}{2}-\frac{1}{\pi }{\int }_0^{\mathrm{\infty }}\frac{\sin (\frac{1}{2}\sum _{j=1}^J[{\tan }^{-1}({\lambda }_j^{\ast }u)]-\frac{1}{2}cu)}{u\prod _{j=1}^J(1+{\lambda }_j^{\ast 2}{u}^2{)}^{\frac{1}{4}}}\mathrm{d}u\equiv \mathrm{Benchmark}(\mathbf{\Omega },\mathbf{L},c).$(21) ). In Tables 1–2, the Feasible Imhof distribution, $\mathrm{Imhof}(\hat{\mathbf{\Omega }},\mathbf{N}\hat{\mathbf{\Omega }},c)$ obtained with $\mathbf{N}\hat{\mathbf{\Omega }}$ estimated using the HC estimators in Section 2.

Notes: See Table 1.

Notes: See Table 1.

Notes: See Table 1.

Notes: See Table 1

Notes: Reported are the F statistics testing for each of the six different single restrictions. Four different p-values of each F statistic are computed from using the Imhof distribution, the F-table, and two different ways of the wild bootstrap methods.

1 A related literature to the heteroskedasticity-consistent (HC) covariance matrix estimation of the regression coefficients is on the estimation of the heteroskedastic error variances themselves. Rao [23] proposes the MInimum Norm Quadratic Unbiased Estimation (MINQUE) of the heteroskedastic error variances. Horn et al. [16] proposes an ‘almost unbiased’ estimator as an estimator of heteroskedastic error variance, which MacKinnon and White [5] uses to construct HC2 as appears below in Equation (7(7) $\hat{\mathbf{\Omega }}=\mathrm{diag}\left\{\frac{e_i^2}{1-h_{ii}}\right\}.$(7) ). The variance estimator of Horn et al. [16] is almost unbiased as its asymptotic bias is of order $n^{-1}$ when errors are heteroskedastic, and it is unbiased when errors are homoskedastic.

2 For applications of Imhof [15] in a separate context of studying distributional properties of estimators, see Bao et al. [24], Ullah [25], Nakamura and Nakamura [26], and Farebrother [27], among others.

3 The accurate test size of using the exact Imhof distribution would also produce more accurate confidence intervals which can be computed from inverting the test statistics using the Imhof distribution.

4 The original HC3 in [5, Equation 12] is slightly different. See also [8, p.233].

5 The notation $h_r,r=1,\dots ,m$, used in Equations (1.1) and (3.2) of [15] is ${\nu }_j,j=1,\dots ,J$ in this paper. The values of ${\delta }_r$ and x in these Imhof's equations are both zero in this paper.

6 Note that $\mathbf{L}$ and $\mathbf{H}\equiv {\mathbf{\Omega }}^{1/2}{\mathbf{L}\mathbf{\Omega }}^{-1/2}={\mathbf{\Omega }}^{1/2}({\mathbf{N}}^{\ast }){\mathbf{\Omega }}^{1/2}$ have the same characteristic equations and the same eigenvalues, because $|\mathbf{L}-\lambda \mathbf{I}|=0$ and $|{\mathbf{\Omega }}^{1/2}||\mathbf{L}-\lambda \mathbf{I}||{\mathbf{\Omega }}^{-1/2}|=0.$ As $b^{\ast }$ is a scalar, $\mathbf{H}={\mathbf{\Omega }}^{1/2}\mathbf{X}({\mathbf{X}}^{\mathrm{\prime }}\mathbf{X}{)}^{-1}\mathbf{r}{\mathbf{r}}^{\mathrm{\prime }}({\mathbf{X}}^{\mathrm{\prime }}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{\prime }}{\mathbf{\Omega }}^{1/2}/b^{\ast }$ is symmetric and idempotent. Hence, the eigenvalues of $\mathbf{H}$ and $\mathbf{L}$ are either 1's or 0's. Further, rank$(\mathbf{H})=1$ because $\mathbf{H}=\xi {\xi }^{\mathrm{\prime }}/b^{\ast }$ with $\xi \equiv {\mathbf{\Omega }}^{1/2}\mathbf{X}({\mathbf{X}}^{\mathrm{\prime }}\mathbf{X}{)}^{-1}\mathbf{r}$ being an $n\times 1$ vector. Therefore, only one of the eigenvalues of $\mathbf{H}$ or $\mathbf{L}$, is 1 and all the other eigenvalues are 0. Then, (19(19) $\begin{array}{rl}\theta (u)& =\frac{1}{2}\sum _{j=1}^J[{\nu }_j{\tan }^{-1}({\lambda }_j^{\ast }u)]-\frac{1}{2}cu,\end{array}$(19) ) and (20(20) $\begin{array}{rl}\rho (u)& =\prod _{j=1}^J(1+{\lambda }_j^{\ast 2}{u}^2{)}^{\frac{1}{4}{\nu }_j},\end{array}$(20) ) are simplified as $\theta (u)=\frac{1}{2}{\tan }^{-1}(u)-\frac{1}{2}cu$ and $\rho (u)=(1+u^2{)}^{\frac{1}{4}}$.

7 For $\mathbf{y}\sim N(\mathbf{X}\beta \mathbf{,}\mathbf{\Omega })$, ${\mathbf{r}}^{\mathrm{\prime }}\mathbf{(}{\mathbf{X}}^{\mathrm{\prime }}\mathbf{X}{\mathbf{)}}^{\mathbf{-}\mathbf{1}}{\mathbf{X}}^{\mathrm{\prime }}\mathbf{y}\mathbf{-}{\mathbf{r}}^{\mathrm{\prime }}\beta \sim N(0,b^{\ast })$, where $b^{\ast }={\mathbf{r}}^{\mathrm{\prime }}\mathbf{(}{\mathbf{X}}^{\mathrm{\prime }}\mathbf{X}{\mathbf{)}}^{\mathbf{-}\mathbf{1}}{\mathbf{X}}^{\mathrm{\prime }}\mathbf{\Omega }\mathbf{X}\mathbf{(}{\mathbf{X}}^{\mathrm{\prime }}\mathbf{X}{\mathbf{)}}^{\mathbf{-}\mathbf{1}}\mathbf{r}$. Therefore, ${\left(\frac{{\mathbf{r}}^{\mathrm{\prime }}\mathbf{(}{\mathbf{X}}^{\mathrm{\prime }}\mathbf{X}{\mathbf{)}}^{\mathbf{-}\mathbf{1}}{\mathbf{X}}^{\mathrm{\prime }}\mathbf{y}\mathbf{-}{\mathbf{r}}^{\mathrm{\prime }}\beta }{\sqrt{b^{\ast }}}\right)}^2\sim {\chi }_1^2.$

8 Unlike in simulation where we have the Benchmark in (21(21) $Pr(F^{\ast }\le c)=\frac{1}{2}-\frac{1}{\pi }{\int }_0^{\mathrm{\infty }}\frac{\sin (\frac{1}{2}\sum _{j=1}^J[{\tan }^{-1}({\lambda }_j^{\ast }u)]-\frac{1}{2}cu)}{u\prod _{j=1}^J(1+{\lambda }_j^{\ast 2}{u}^2{)}^{\frac{1}{4}}}\mathrm{d}u\equiv \mathrm{Benchmark}(\mathbf{\Omega },\mathbf{L},c).$(21) ) which can be used to evaluate the performance of the Imhof distribution, we do not have such a benchmark criterion in empirical applications. Thus we compare the feasible Imhof p-values with the wild bootstrap p-values. Many simulation studies suggested that the wild bootstrap gives good performance (e.g. [28]), while some recent studies, e.g. [29], suggested that alternative variants of the wild bootstrap may perform quite differently. Here, we use the wild bootstrap just to verify that the Imhof procedure works.

Additional information

Funding

This work was supported by the National Natural Science Foundation of China [grant number 71801184], Natural Science Foundation of Fujian Province [grant number 2018J01114] and the UCR Academic Senate for research funds.