Adaptive-to-model hybrid of tests for regressions

TL;DR

This paper introduces an adaptive hybrid testing method for regressions that combines the strengths of nonparametric and empirical process-based tests, improving sensitivity and distribution tractability.

Contribution

It proposes a novel adaptive-to-model hybrid test that inherits the advantages of both test types and enhances nonparametric tests with empirical process properties.

Findings

The hybrid test effectively detects various alternative models.

Numerical studies confirm the test's superior power.

The method is adaptable to different data types.

Abstract

In model checking for regressions, nonparametric estimation-based tests usually have tractable limiting null distributions and are sensitive to oscillating alternative models, but suffer from the curse of dimensionality. In contrast, empirical process-based tests can, at the fastest possible rate, detect local alternatives distinct from the null model, but is less sensitive to oscillating alternative models and with intractable limiting null distributions. It has long been an issue on how to construct a test that can fully inherit the merits of these two types of tests and avoid the shortcomings. We in this paper propose a generic adaptive-to-model hybrid of moment and conditional moment-based test to achieve this goal. Further, a significant feature of the method is to make nonparametric estimation-based tests, under the alternatives, also share the merits of existing empirical process-based tests. This methodology can be readily applied to other kinds of data and constructing other hybrids. As a by-product in sufficient dimension reduction field, the estimation of residual-related central subspace is used to indicate the underlying models for model adaptation. A systematic study is devoted to showing when alternative models can be indicated and when cannot. This estimation is of its own interest and can be applied to the problems with other kinds of data. Numerical studies are conducted to verify the powerfulness of the proposed test.

Tables

Table 1. Table 1: Empirical sizes and powers with X ∼ N ​ ( 0 , Σ ~ ) similar-to 𝑋 𝑁 0 ~ Σ X\sim N(0,\tilde{\Sigma}) in Study 1 1 1

$p=2$	$a$	$T^{Zh}$	$T^n$	$T^{SMQ}$	$T^{GWZ}$	$T^B$
	0	0.062	0.067	0.050	0.067	0.049
	0.2	0.077	0.205	0.136	0.088	0.172
	0.4	0.197	0.631	0.426	0.211	0.598
	0.6	0.438	0.944	0.755	0.486	0.956
	0.8	0.743	0.999	0.938	0.833	1.000
	1	0.934	1.000	0.999	0.985	1.000
$p=8$	$a$	$T^{Zh}$	$T^n$	$T^{SMQ}$	$T^{GWZ}$	$T^B$
	0	0.018	0.070	0.041	0.059	0.009
	0.2	0.024	0.105	0.110	0.086	0.025
	0.4	0.025	0.268	0.291	0.101	0.198
	0.6	0.047	0.534	0.595	0.159	0.649
	0.8	0.053	0.779	0.801	0.243	0.943
	1	0.075	0.896	0.915	0.384	0.997

Table 2. Table 2: Empirical sizes and powers in Study 3 3 3

$a$	$T^{Zh}$	$T_n$	$T^{SMQ}$	$T^{GWZ}$	$T^B$
0	0.013	0.066	0.016	0.069	0.002
0.2	0.006	0.155	0.035	0.072	0.005
0.4	0.020	0.445	0.099	0.180	0.043
0.6	0.020	0.661	0.109	0.439	0.130
0.8	0.021	0.780	0.194	0.734	0.254
1	0.016	0.798	0.240	0.895	0.403

Equations

$$Y=G(X,\epsilon ),$$

$$Y=g(X,\theta )+\epsilon ,$$

$$H_0:Y=g(X,{\theta }_0)+\epsilon ,\text{for some }{\theta }_0\in \Theta ,$$

$$H_1:Y=G(X,\epsilon )\ne g(X,\theta )+\epsilon ,\text{for all }\theta \in \Theta .$$

$$V=E\{\eta w(X)[E(\eta w(X))I(q=0)+E(\eta |X)\tilde{w}(X)I(q>0)]\},$$

$$V_n=\left|{\left[\frac{1}{n}\sum _{j=1}^n{\hat{\eta }}_{j}w(x_j)\right]}^{2}I(\hat{q}=0)+\frac{1}{n(n-1)}\sum _{j=1}^n\sum _{k\ne j,k=1}^n{\hat{\eta }}_j{\hat{\eta }}_k\frac{1}{h^p}{K}_h(x_j-x_k)I(\hat{q}>0)\right|,$$

$$V_0=\frac{1}{n}\sum _{j=1}^n{\hat{\eta }}_{j}w(x_j)\text{and }{V}_1=\frac{1}{n(n-1)}\sum _{j=1}^n\sum _{k\ne j,k=1}^n{\hat{\eta }}_j{\hat{\eta }}_k\frac{1}{h^p}{K}_h(x_j-x_k).$$

$${\sigma }_0^2={\sigma }_{\epsilon }^2[E(w^2)-E(\dot{g}w{)}^{\top }E(\dot{g}{\dot{g}}^{\top }{)}^{-1}E(\dot{g}w)],$$

$${\hat{\sigma }}_0^2=\left(\frac{1}{n}\sum _{j=1}^n{\hat{\eta }}_j^2\right)\left[\frac{1}{n}\sum _{j=1}^n{w}_j^2-(\frac{1}{n}\sum _{j=1}^n{\dot{g}}_{j,n}{w}_j{)}^{\top }(\frac{1}{n}\sum _{j=1}^n{\dot{g}}_{j,n}{\dot{g}}_{j,n}^{\top }{)}^{-1}(\frac{1}{n}\sum _{j=1}^n{\dot{g}}_{j,n}{w}_j)\right].$$

$$T_n=\frac{n{V}_n}{{\hat{\sigma }}_0^2}=\frac{n{V}_0^{2}I(\hat{q}=0)}{{\hat{\sigma }}_0^2}+\frac{n{V}_{1}I(\hat{q}>0)}{{\hat{\sigma }}_0^2}.$$

$$M=\int E[X\exp (it\eta )]E^H[X\exp (it\eta )]d{F}_{\eta }(t),$$

$${\hat{m}}_n(t)=\frac{1}{n}\sum _{j=1}^n{x}_j\exp (it{\hat{\eta }}_j),$$

$${\hat{M}}_n=\frac{1}{n}\sum _{j=1}^n{\hat{m}}_n({\hat{\eta }}_j){\hat{m}}_n^H({\hat{\eta }}_j).$$

$${\hat{\theta }}_n-{\theta }^{\ast }=G^{-1}\frac{1}{n}\sum _{j=1}^n{\eta }_j\dot{g}(x_j,{\theta }^{\ast })+o_p(\frac{1}{\sqrt{n}}),$$

$${\hat{s}}_j^{\ast }=\frac{{\hat{s}}_j^2+c_{1n}}{{\hat{s}}_{j+1}^2+c_{1n}}-1\text{ and }{\hat{r}}_j=\frac{{\hat{s}}_{j+1}^{\ast }+c_{2n}}{{\hat{s}}_j^{\ast }+c_{2n}},$$

$$\hat{q}=\{\begin{array}{ll}0,& \text{ if r^j>τ for all j∈{1,…,p} },\\ \arg \underset{1\le j\le p}{\max }\left\{j:{\hat{r}}_j\le \tau \right\}& \end{array},$$

$$T_n\to {\chi }_1^2$$

$$T_n/n\to \mu >0.$$

$$Y=g(X,{\theta }_0)+{\delta }_n\ell (X)+\epsilon ,$$

$${\hat{\theta }}_n-{\theta }_0=G^{-1}\frac{1}{n}\sum _{i=1}^n{\epsilon }_i{\dot{g}}_i+{\delta }_n{G}^{-1}E(\dot{g}\ell )+o_p(\frac{1}{\sqrt{n}}).$$

$$T^{GWZ}=n{h}^{{\hat{q}}_1/2}\frac{1}{n(n-1)}\sum _{j=1}^n\sum _{k\ne j,k=1}^n{\hat{\eta }}_j{\hat{\eta }}_k\frac{1}{h^{{\hat{q}}_1}}{K}_h({\hat{B}}_1^{\top }(x_j-x_k)),$$

$$Y=0.25\exp (2{\alpha }^{\top }X)+a\sin ({\beta }^{\top }X)+\epsilon .$$

$$Y={\alpha }^{\top }X+0.8({\alpha }^{\top }X{)}^2+a\tanh ({\beta }^{\top }X)+\epsilon .$$

$$Y=0.25\exp (2X_1)+a\{0.5X_2^3+\cos (X_3)+X_4-|X_5|+\tanh (0.6\pi X_6)+X_7{X}_8\}+\epsilon .$$

$$Y=c_1{X}_1+c_2{X}_2^2+c_3{X}_3^3+c_4{X}_4{X}_5+c_5\sin (X_6)+a\{0.2X_1^2+0.3X_2^3\}+\epsilon .$$

$$E(Y|X)$$

$$E\left\{\underset{\theta \in \Theta }{\sup }{g}^2(X,\theta )\right\}<\infty ,$$

$$E\left\{\underset{\theta \in \Theta }{\sup }\Vert \frac{\partial g(X,\theta )}{\partial \theta }\frac{\partial g(X,\theta )}{\partial {\theta }^{\top }}\Vert \right\}<\infty ,$$

$$E\left\{\underset{\theta \in \Theta }{\sup }\Vert [Y-g(X,\theta ){]}^2\frac{\partial g(X,\theta )}{\partial \theta }\frac{\partial g(X,\theta )}{\partial {\theta }^{\top }}\Vert \right\}<\infty ,$$

$$E\left\{\underset{\theta \in \Theta }{\sup }\Vert [Y-g(X,\theta ){]}^2\frac{{\partial }^{2}g(X,\theta )}{\partial \theta \partial {\theta }^{\top }}\Vert \right\}<\infty .$$

$${\theta }^{\ast }=\arg \underset{\theta \in \Theta }{\inf }E\{Y-g(X,\theta ){\}}^2.$$

Full text

Keywords: High-frequency/oscillating model; Local smoothing method; Model adaptation; Model checking

1 Introduction

Suppose that $X\in R^{p}$ is the $p$-dimensional explanatory vector and $Y\in R$ the scalar response variable. Without further information from the data, the underlying regression model can be written as

[TABLE]

where $G(\cdot)$ is a general nonparametric smooth function and ${\varepsilon}$ is the error term independent of $X$. Such a model is too flexible to be interpreted, therefore in practical applications, a parametric model is often preferred:

[TABLE]

where $g(\cdot,\cdot)$ is a given function up to a parameter $\theta\in R^{d}$. A natural issue is to check whether this model is an adequate depiction of the relationship between $X$ and $Y$, since a wrongly specified model will cause unreliable statistical inferences. In this paper, we focus on this homogeneous null model. For more general heterogeneous models, we will give a discussion in Section 6.

Many efforts have been devoted to this model checking problem since 1980’s with a number of tests available. To demonstrate the necessity of revisiting this issue, we first give a brief review of the pros and cons of existing methodologies and what are the difficulties any existing test cannot well handle. In the literature, two classes of the most popularly used methods in this area are nonparametric estimation-based tests which are the so-called local smoothing tests. Early references of local smoothing tests include Härdle and Mammen (1993), Zheng (1996), Fan and Li (1996) and Koul and Ni (2004), among others. Hart (1997) is a comprehensive reference for the early research in this field. The most appealing feature of local smoothing tests is their sensitivity to high-frequency/oscillating alternative models. However, the typical drawback is the slow rates of convergence, which causes low power performance when the dimension is high. For instance, the test statistic in Zheng (1996), multiplied by $nh^{p/2}$, converges to its weak limit where $n$ is the sample size and $h\to 0$ is the bandwidth in the Nadaraya–Watson kernel estimation of the conditional moment $E(Y-g(X,\theta){|}X)$. The rate $nh^{p/2}$ can be much slower than $n$ when the dimension $p$ is even moderate and then causes that they can only detect the local alternatives distinct from the null at the rate slower than $n^{-1/2}h^{-p/4}$. Thus, the dimension of $X$ severely worsens the power performance (see Guo et al. (2016)). To alleviate this difficulty, Lavergne and Patilea (2012) used the projection method that can detect the local alternatives at the rate of order $n^{-1/2}h^{-1/4}$ in theory. Their test involves the approximation for high-dimensional integral, which is computationally expensive (see Zhu et al. (1995)). Further, when the approximation is not sufficiently good or the dimension of $X$ is high, the test cannot maintain the significance level well. Under a dimension reduction structure, Guo et al. (2016) suggested an adaptive-to-model approach that can also reach the rate of order $n^{-1/2}h^{-1/4}$, without that computational issue. However, both are still (actually this is the case for any local smoothing test in the literature) impossible to reach the fastest possible rate of order $1/\sqrt{n}$ to detect local alternatives.

Unlike local smoothing tests, global smoothing test statistics are typically the averages of weighted sums of residuals. They are called the global smoothing tests as averaging in the CUSUM nature is a global smoothing step. The most significant advantage of global smoothing tests is that they can detect local alternatives at the fastest possible rate of $n^{-1/2}$ and converge to its weak limits at this rate (or $1/n$ if the quadratic form is used). However, the main shortcomings of global smoothing methods are 1) the CUSUM structure often weakens their ability to detect oscillating alternatives; 2) the intractability of the limiting null distributions (see Stute et al. (1998); Dette et al. (2007)) requires the assistance from Monte Carlo approximation/resampling technique for critical value determination unless the dimension $p$ is $1$ (Stute, 1997) or a directional test is used (Stute and Zhu, 2002), which is time consuming.

Above expounds the well-known pros and cons of these two types of tests showing that no methodology in the literature can simultaneously inherit the respective advantages and avoid the shortcomings. In summary, having a test that has the following three features is desirable:

• $f_{1}$)

at the rate of $1/\sqrt{n}$, the constructed test statistic converges to a tractable weak limit such that the critical value can be easily determined without the assistance of Monte Carlo approximation/resampling technique;

• $f_{2}$)

the test can be sensitive to oscillating alternatives as local smoothing tests, but less influenced by dimensionality;

• $f_{3}$)

more than the omnibus property, at the fastest possible rate of order $1/\sqrt{n}$, the test can detect local alternatives distinct from the null model, and thus is more powerful than any existing local smoothing tests in theory.

To the best of our knowledge, any single test in the literature cannot share all these appealing features. To arrive the above goals, the test must have both local and global smoothing component. But due to their different convergence rates, any simple convex combination does not work for this purpose. We will see this clearly in Sections 2 and 4.

In this paper, we attempt to attack this longstanding problem. To be specific, we propose an adaptive-to-model hybrid of tests that is a combination of moment-based test and nonparametric conditional moment-based test component. This hybrid can automatically indicate the underlying model such that it can decide which component works and then fully inherit the merits of local and global smoothing tests described in the above features $f_{1}$)-$f_{3}$). Model adaptation is achieved through an indicative dimension of a residual-based central subspace for the underlying model in the sense of sufficient dimension reduction. Thus, under the null hypothesis, it derives a zero-dimensional projection such that the test becomes a simple moment-based test for critical value determination. Under the alternatives, the test automatically becomes a nonparametric conditional moment-based omnibus test. More interestingly and importantly, the special construction of the test makes the conditional moment-based component achieve the fastest possible rate of convergence. This is the most significant contribution of the proposed method as any existing nonparametric estimation-based test is impossible to have such a feature. The test is based on the simplest moment test and a typical nonparametric conditional moment test proposed by Zheng (1996). In effect, we can choose any nonparametric estimated-based test in the hybrid. The asymptotic properties could be similar. We will have some more discussion in Section 6.

For model adaptation, we consider a residual-based central subspace with indicative dimension. This concept is slightly different from that in sufficient dimension reduction field (Li, 1991; Cook, 1998) because this dimension can indicate the underlying model such that the proposed hybrid can do so. As a by-product, we will propose a target matrix and suggest a criterion to define an estimator of the indicative dimension. The relevant properties show that the estimated dimension can adapt to the underlying model even when the local alternative models approximate to the null model at the optimal rate that is as close to $n^{-1/2}$ as possible in a certain sense. It improves existing results for model adaptation to the underlying model through the dimension indication at the rate of order slower than $n^{-1/2}h^{-1/4}$, see Guo et al. (2016) even when a dimension reduction structure is assumed. The desirable asymptotic properties of the proposed hybrid test are derived under the null, local and global alternative hypothesis. The resulting test can be very different from the adaptive-to-model one in Guo et al. (2016) which cannot change its local smoothing nature and is impossible to inherit the merits of global smoothing tests.

It is worthwhile to mention that this generic methodology is ready to apply to any pair of tests to obtain the hybrid. Also, it can be applied to other kinds of data such as measurement error data, panel data and functional data. We will have a brief discussion in Section 6.

The materials in this paper are organized as follows. In Section 2, we present the hypotheses and the test statistic construction. A target matrix and a criterion for estimating indicative dimension are suggested in Section 3. The various rates under the null, local and global alternative hypothesis for dimension indication are systematically studied and the optimal indicative rate is derived. Section 4 contains the asymptotic properties of the test statistic. Numerical performances of the test under different models are examined by various experiments in Section 5. Section 6 includes some discussions about the main limitations and possible extension to heteroscedastic models. The regularity conditions of the theorems are listed in the Appendix A. The proofs for the theorems and details of dimension estimation are presented in the supplementary materials111See Supplementary materials..

2 The test statistic construction

2.1 A brief review of sufficient dimension reduction

Suppose there are two random column vectors $Z_{1}\in R^{p_{1}}$ and $Z_{2}\in R^{p_{2}}$. If there exists a matrix $B\in R^{p_{1}\times q}$, such that, Z_{2}\rotatebox[origin={c}]{90.0}{\models}Z_{1}{|}B^{\top}Z_{1}. Then the column space $Span\{B\}$ of $B$ is called a sufficient dimension reduction subspace of $Z_{2}$ with respect to $Z_{1}$. The intersection of all the dimension reduction subspaces is called the central subspace and denoted as $S_{Z_{2}{|}Z_{1}}$. The dimension of the central subspace is denoted as $dim(S_{Z_{2}{|}Z_{1}})$. If $B$ is the central subspace, then its column dimension $q=dim(S_{Z_{2}{|}Z_{1}})$. When the real working dimension $q$ is smaller than the original dimension $p_{1}$ of $Z_{1}$ and dimension reduction can be achieved in terms of identifying the matrix $B$ and its dimension $q$. See Li (1991) and Cook (1998) for more details. We will use these notations during the test construction next.

2.2 The hypotheses

Based on (1) and (2), the null hypothesis we concern about is

[TABLE]

where $\Theta$ is a subset of $R^{d}$ and {\varepsilon}\rotatebox[origin={c}]{90.0}{\models}X is the error term. The alternative hypothesis is

[TABLE]

Let $\theta^{\ast}=\arg\min_{\theta\in\Theta}E\{Y-g(X,\theta)\}^{2}$ and $\eta=Y-g(X,\theta^{\ast})$ provided that the involved second order moments exist. Under the conditions in Appendix, the minimizer $\theta^{\ast}$ is identifiable and can be consistently estimated. Under the null hypothesis, we have $\eta={\varepsilon}$, $\theta^{\ast}=\theta_{0}$ and then $dim(S_{\eta{|}X})=0$ by the independence between $\eta$ and $X$. While under the alternative hypothesis, the residual $\eta=G(X,{\varepsilon})-g(X,\theta^{\ast})-{\varepsilon}=\Delta(X,{\varepsilon})$, which leads to $dim(S_{\eta{|}X})>0$ since $\Delta(X,{\varepsilon})$ is a non-constant function of $X$. Thus, $dim(S_{\eta{|}X})$ is respectively equal to [math] and larger than [math] under the null and alternative hypothesis. This inspires us to construct a test that can fully use the information provided by this indicative dimension. In the following, we implement our idea.

2.3 The hybrid of tests

We in this subsection describe our idea in terms of using two simple tests as the components in the hybrid. Following description essentially provides a general framework that can be used to develop other hybrids.

Suppose $\{(x_{1},y_{1}),\dots,(x_{n},y_{n})\}$ is an available sample, where $(x_{i},y_{i}),i=1,\dots,n$ are independent and identically distributed. For the sake of illustration, the unknown parameter $\theta$ is estimated by the least squares method. Define $\hat{\theta}_{n}=\arg\min_{\theta\in\Theta}\sum_{j=1}^{n}[y_{j}-g(x_{j},\theta)]^{2}$ and $\hat{\eta}_{j}=y_{j}-g(x_{j},\hat{\theta}_{n})$.

To achieve the optimal rate of convergence and the tractability of the limiting null distribution, we can simply use the sum of weighted residuals: $\sum_{j=1}^{n}\hat{\eta}_{j}W(x_{j})/n$ with some weight function $W(\cdot)$. This is because under the null, $E(\eta W(X))=0,$ and $\sum_{j=1}^{n}\hat{\eta}_{j}W(x_{j})/\sqrt{n}$ converges to a normal distribution $N(0,\sigma^{2})$ where $\sigma^{2}$ is easy to obtain. This simple test can fulfill feature $f_{1}$) with good performance on the critical value determination and significance level maintenance because it involves fewer unknowns in the limiting variance $\sigma^{2}$, thus its estimation error would be less. However, no researcher would simply use this moment-based test in practice as it is obviously not omnibus and not powerful at all for general alternative models. Further, notice that the conditional moment of $\eta$ given $X$ has the property: for a nonnegative weight function $\tilde{w}(X)$, $E(E(\eta{|}X)\tilde{w}(X))^{2}=0$ and $>0$ under the null and alternatives respectively. It can be estimated by, say, the Nadaraya-Watson kernel estimation to define a local smoothing test that is omnibus and sensitive to oscillating alternative models. This leads to a typical nonparametric estimation-based test, see Zheng (1996). Then feature $f_{2}$) is achieved. But as we mentioned before, it has slow rate of convergence and suffers from the curse of dimensionality, see Guo et al. (2016) for more details. Both have some very obvious shortcomings and none of these tests can achieve feature $f_{3}$). But we still stick to these two tests to see how to construct an adaptive-to-model hybrid that shares all appealing features $f_{1}$) – $f_{3}$). More importantly, the hybrid can make, under the alternatives, the above nonparametric conditional moment-based component also share the properties of global smoothing tests satisfying feature $f_{3}$). This is a somewhat surprising property that will be stated later.

For this mission, we need assistance from other information. Recall that $q=dim(S_{\eta{|}X})$ is [math] under the null and is greater than 0 under the alternatives. Consider a hybrid of weighted moment and conditional weighted moment of the residual $\eta$ in the following format:

[TABLE]

where $w(\cdot)$ and $\tilde{w}(\cdot)$ are two weight functions to be determined later. Thus, under the null, the term $V=E^{2}[\eta w(X)]=E^{2}(\eta)E^{2}(w(X))=0$ and under the alternatives $V=E(E^{2}(\eta{|}X)\tilde{w}(X){w}(X))=E\{[E(\Delta(X,{\varepsilon}){|}X)]^{2}f(X)\}>0$ when we choose $\tilde{w}=f/w$ where $f$ is the density function of $X$. Theoretically, the two weight functions can be relatively arbitrarily chosen, but we do have some preferences in practice for this choice, see, e.g. Zheng (1996). When the first component is estimated by the average of weighted residuals and the second component is estimated by a nonparametric method, the constructed test can satisfy features $f_{1}$) and $f_{2}$). Therefore, at the sample level, an estimator of $V$ in (5) is defined as

[TABLE]

where $K_{h}(\cdot)=K(\cdot/h)$, $K(\cdot)$ is the kernel function, $h$ is the corresponding bandwidth, and $\hat{q}$ is the estimated indicative dimension. Here taking the absolute value is just to keep the non-negativity of $V$ at the sample level. To get a concise expression, denote

[TABLE]

Thus $V_{n}=|V_{0}^{2}I(\hat{q}=0)+V_{1}I(\hat{q}>0)|$.

It is clear that to satisfy both features $f_{1}$) and $f_{2}$), the estimator $\hat{q}$ must also be [math] and $>0$ with a probability going to $1$ respectively under the null and alternatives. We will discuss this property in Section 3 via sufficient dimension reduction. But a very important thing is about the choice of standardizing constant for $V_{n}$ in (6) such that the test can also share feature $f_{3})$. A standardizing constant $b_{n}$ must be used such that $T_{n}=b_{n}V_{n}$ has a tractable limiting null distribution. As under the null $V_{n}=V_{0}^{2}$, $b_{n}=n/\hat{\sigma}^{2}_{0}$ should be used where $\hat{\sigma}^{2}_{0}$ is an estimator of the limiting variance such that $b_{n}V_{0}^{2}$ can be asymptotically chi-square. Under the alternatives this $b_{n}$ diverges to infinity faster than $nh^{p/2}$ (Zheng, 1996) or $nh^{1/2}$ (Lavergne and Patilea, 2012; Guo et al., 2016). Therefore, $b_{n}V_{n}=b_{n}V_{1}$ can be much powerful than Zheng (1996)’s test $nh^{p/2}V_{1}$ or the tests $nh^{1/2}V_{1}$ if $V_{1}$ is the test in Lavergne and Patilea (2012) or Guo et al. (2016). This is very unique merit of this adaptive-to-model hybrid. This discussion gives us the idea about how the constructed test could fulfill feature $f_{3}$). We will give more details about the asymptotic properties of the test in Section 4.

We now determine which estimator of $\sigma_{0}^{2}$ we should use. As $V_{n}=V_{0}^{2}$ under the null where $V_{0}$ has simple structure and its limiting variance $\sigma_{0}^{2}$ has fewer unknowns, we then use it. For the simplicity of notation, let $E(w^{2})=E(w(X)^{2})$, $E(\dot{g}w)=E(\dot{g}(X,\theta^{\ast})w(X))$ and $E(\dot{g}\dot{g}^{\top})=E(\dot{g}(X,\theta^{\ast})\dot{g}(X,\theta^{\ast})^{\top})$. Under the null, we have

[TABLE]

where $\sigma_{\varepsilon}^{2}=Var({\varepsilon})$. Write $w_{j}=w(x_{j})$ and $\dot{g}_{j,n}=\dot{g}(x_{j},\hat{\theta}_{n})$. A consistent estimator of $\sigma_{0}^{2}$ is

[TABLE]

Then, the resulting test statistic is defined as

[TABLE]

The weight function is chosen as $w(x)=c\exp(-\|x\|)$ for some constant $c>0$ to be specified in the simulations, where $\|\cdot\|$ stands for the Euclidean norm.

Remark 1**.**

*We comment on the construction in two aspects.

1). Any nonzero function is applicable as the weight function theoretically, but empirically, we find that a “small” weight helps enhance the power performance in finite sample scenarios. In theory, a natural question is whether there is an optimal choice of weight function. However, the optimality here is hard to define. For instance, if the optimality is on the magnitude of limiting variance of $V_{0}$, the standardization removes the scale. If the optimality is about the power performance of the test, it relates to the issue about how to construct a most powerful test in a certain sense, which is beyond the scope of this paper. Thus, we do not discuss this issue in more detail in this paper.

2). Another issue is about choosing an estimator for the limiting variance of $T_{n}$ as $T_{n}$ involves two tests in effect. We now explain the use of $\hat{\sigma}_{0}^{2}$. Let $\Sigma_{Zh}$ be the limiting variance of $nh^{p/2}V_{1}$. From Lemma 3.3 of Zheng (1996), we can see that $\sigma_{0}^{2}$ can be smaller than $\Sigma_{Zh}$ even under the null. In other words, if we use an estimator of $\Sigma_{Zh}$ in lieu of $\sigma_{0}^{2}$, the limiting distribution under the null cannot be a standard chi-square and the values of $T_{n}$ under the alternatives can be smaller to lower the power.*

3 Estimation of the indicative dimension

In this section, we propose a method to estimate the dimension $q$. To this end, we define a target matrix with $q$ nonzero eigenvalues below.

3.1 Target matrix

As $\eta$ needs to be replaced by its estimator, the estimation of target matrix and the dimension $q$ of the corresponding central subspace becomes more complicated. We then consider the following target matrix to make the estimation relatively easier although existing methods such as sliced inverse regression (Li, 1991) and sliced average variance estimation (Cook and Weisberg, 1991) could also be used. Denote $A^{H}$ as the conjugate transpose of a matrix $A$ and define the target matrix as

[TABLE]

where $i=\sqrt{-1}$ is the complex number. For ease of illustration, assume $E(X)=0$ and write $Var(X)=\Sigma$. Let $B$ be a basis of $S_{\eta{|}X}$ and $P_{B}(\Sigma)=B(B^{\top}\Sigma B)^{-1}B^{\top}\Sigma$. The linearity condition is $E(X{|}B^{\top}X)=P_{B}^{\top}(\Sigma)X$, see condition 3.1 in Li (1991), which is widely used in sufficient dimension reduction. It holds for $X$ that follows symmetric elliptical distribution.

Lemma 1**.**

Under the linearity condition, $Span\{M\}\subseteq\Sigma S_{\eta{|}X}$. If the rank of $M$ is $q$, then $Span(M)=\Sigma S_{\eta{|}X}$.

Further, we have $q=rank(M)=dim(S_{\eta{|}X})$ since $dim(\Sigma S_{\eta{|}X})=dim(S_{\eta{|}X})$. Let $\lambda_{1}\geq\dots\geq\lambda_{p}\geq 0$ be the ordered eigenvalues of $M$ defined in (10). Then under the null hypothesis, the independence between $\eta$ and $X$ leads to $\lambda_{1}=\dots=\lambda_{p}=0$ and $q=0$. In contrast, under the alternatives, $\lambda_{1}\geq\dots\geq\lambda_{q}>0=\lambda_{q+1}=\dots=\lambda_{p}$ for some $1\leq q\leq p$.

3.2 Estimation

Consider estimating the target matrix first. Define $m(t)=E(X\exp(it\eta))$. An estimator of $m(t)$ is

[TABLE]

where $\hat{\eta}_{j}=y_{j}-g(x_{j},\hat{\theta}_{n})$. The corresponding estimator of $M$ is

[TABLE]

To study this, we first give the consistency of $\hat{\theta}_{n}$ in the following lemma.

Lemma 2**.**

Under the regularity conditions in Appendix, we have $\|\hat{\theta}_{n}-\theta^{\ast}\|=O_{p}(n^{-1/2})$ and

[TABLE]

where $G=E[\dot{g}(X,\theta^{\ast})\dot{g}^{\top}(X,\theta^{\ast})]$.

This asymptotically linear representation is a standard result for the least squares estimation, see Zheng (1996) for details.

Let $X_{i}$ be the $i$-th component of $X$ and $\|\cdot\|$ be the Frobenius norm. The following theorem states the consistency of $\hat{M}_{n}$.

Theorem 1**.**

Suppose $\max_{1\leq j\leq p}E(|X_{j}|^{4})<\infty$. Then

under the null hypothesis, $\|\hat{M}_{n}-M\|=O_{p}(n^{-1})$; 2. 2.

under the alternative hypothesis, $\|\hat{M}_{n}-M\|=O_{p}(n^{-1/2})$.

These results provide a base for estimating the indicative dimension $q$. We then suggest a criterion that is a slight modification of the thresholding double ridge ratio (TDRR, hereafter) method developed by Zhu et al. (2016).

Define the eigenvalues of $\hat{M}_{n}$ to be $\hat{\lambda}_{1}\geq\dots\geq\hat{\lambda}_{p}\geq 0$ and $\hat{s}_{j}=\hat{\lambda}_{j}/(\hat{\lambda}_{j}+1)$. We then use the following two-step ratio to determine the estimation of dimension $q$. Let

[TABLE]

where $c_{1n}$ and $c_{2n}$ are the two ridges that converge to [math] in proper rates to be selected later, such that the dimension $q$ can be identified. The criterion for the determination of $q$ is

[TABLE]

with a threshold $0<\tau<1$. Based on the rule of thumb as in Zhu et al. (2016), we also set $\tau=0.5$. Details of the modified TDRR method and corresponding proofs are shown in the supplementary materials.

Following theorem states the consistency of $\hat{q}$ to indicate the underlying model.

Theorem 2**.**

Under the conditions in Appendix, we have

if $c_{1n}\to 0$, $c_{2n}\to 0$ and $n^{2}c_{1n}c_{2n}\to\infty$, then under the null hypothesis of (3), ${\mathbb{P}}(\hat{q}=0)\to 1$; 2. 2.

if $c_{1n}\to 0$, $c_{2n}\to 0$ and $nc_{1n}c_{2n}\to\infty$, then under the alternative hypothesis of (4), ${\mathbb{P}}(\hat{q}=q>0)\to 1$;

Due to the different properties of the target matrix from those in Zhu et al. (2016), the ridges are chosen differently from the ones used in their method. Based on our experience in the numerical studies, the recommended ridges are $c_{1n}=3e^{-4}\sqrt{8}\log n/\sqrt{n}$ and $c_{2n}=4\sqrt{8}\log n/5\sqrt{n}$.

4 Asymptotic properties

4.1 Asymptotics under different hypotheses

Theorem 2 offers a very important result such that the proposed test is of the following asymptotic model adaptation property.

Lemma 3**.**

Assume the conditions in Appendix. Then as $n\to\infty$, with a probability going to $1$, under the null hypothesis of (3), $T_{n}=nV_{0}^{2}/\hat{\sigma}^{2}_{0}$ and under the alternative hypothesis of (4), $T_{n}=n|V_{1}|/\hat{\sigma}^{2}_{0}$.

From this lemma, we first give the result about the limiting null distribution.

Theorem 3**.**

Under the null hypothesis of (3) with the regularity conditions in Appendix, the test statistic satisfies

[TABLE]

in distribution where $\chi_{1}^{2}$ stands for chi-square distribution with one degree of freedom.

Now we investigate the power performance of the test under the alternative hypothesis of (4). Recall that $\eta=\Delta(X,{\varepsilon})$. For notational simplicity, rewrite $\Delta(X,{\varepsilon})$ as $\Delta$, $f(X)$ as $f$, $\dot{g}(X,\theta^{\ast})$ as $\dot{g}$ and $w(X)$ as $w$ hereafter. Denote $\Sigma_{01}=E(\Delta^{2})[E(w^{2})-E(\dot{g}w)^{\top}E(\dot{g}\dot{g}^{\top})^{-1}E(\dot{g}w)]$ and $\mu=E(\Delta^{2}f)/\Sigma_{01}$. We have the following result.

Theorem 4**.**

Given the regularity conditions in Appendix, under the alternative hypothesis of (4), in probability,

[TABLE]

Now we consider the sequence of local alternative models:

[TABLE]

where $\delta_{n}\to 0$ as $n\to\infty$ and $\ell(\cdot)$ is a non-constant function.

Following lemma states the consistency of $\hat{\theta}_{n}$ under such local alternatives.

Lemma 4**.**

Under the local alternatives in (16) and the conditions in Appendix, we have

[TABLE]

Before presenting the asymptotic result of the test, we give some results about the estimator $\hat{q}$ under the local alternatives to show why we call the dimension the indicative dimension.

Theorem 5**.**

Assume that the conditions in Appendix hold, then under the local alternative model of (16),

if $\delta_{n}=n^{-1/2}$, $c_{1n}\to 0$, $c_{2n}\to 0$ and $n^{2}c_{1n}c_{2n}\to\infty$, then ${\mathbb{P}}(\hat{q}=0)\to 1$; 2. 2.

if $\delta_{n}=n^{-\alpha}$ for some $0<\alpha<1/2$, $c_{1n}=o_{p}(\delta_{n}^{4})$, $c_{2n}\to 0$, $nc_{1n}c_{2n}/\delta_{n}^{2}\to\infty$ and $c_{1n}c_{2n}/\delta_{n}^{6}\to\infty$, then ${\mathbb{P}}(\hat{q}=1)\to 1$.

Remark 2**.**

Theorem 5 provides two pieces of important information on the model indication through $\hat{q}$. First, for any rate of $\delta_{n}$ as close to $1/\sqrt{n}$ as possible, $\hat{q}=1$ cannot be equal to the true dimension $q$, but can still well indicate the local alternatives. This rate is the fastest possible rate for such a separation in this research area as when $\delta_{n}$ goes to zero at or faster than $1/\sqrt{n}$, $\hat{q}$ loses the indication ability to the local alternatives. In other words, TDRR can reach the optimal rate $\delta_{n}$ of convergence while the only existing method proposed by Guo et al. (2016) can only identify the dimension at the rate such that $n^{-1/2}h^{-1/4}/\delta_{n}=o(1)$ even under the dimension reduction framework. We also note that these results, similarly as those in Guo et al. (2016), basically have theoretical meaning, unless we have prior information on the closeness of local alternatives to the null, we can then set those values. It deserves further study to make it useful in practice.

Based on these results, we are now in the position to discuss the asymptotic properties of the proposed test in detail. Define $\mu_{0}=E(\ell w)-E(\dot{g}w)^{\top}E(\dot{g}\dot{g}^{\top})^{-1}E(\dot{g}\ell)$, $\mu_{1}=E(\ell^{2}f)/\sigma_{0}^{2}$, $\Sigma_{1}=2\sigma_{\varepsilon}^{2}\int K^{2}(u)du\int f^{2}(x)dx$ and $Var=\Sigma_{1}/\sigma_{0}^{2}$. We state the power performance of the test in various scenarios.

Theorem 6**.**

Given conditions in Appendix with the ridges $c_{1n}$ and $c_{2n}$ satisfying the conditions in Theorem 5, if $h\to 0$ and $nh^{p}\to\infty$, then under the local alternatives of (16),

if $\delta_{n}=n^{-1/2}$, $T_{n}\to\chi_{1}^{2}(\mu_{0}^{2}/\Sigma_{0})$ in distribution; 2. 2.

if $\delta_{n}=n^{-\alpha}$, $0<\alpha<1/2$, and

1. (a)

$n^{1/2}h^{p/4}\delta_{n}\to 0$, $T_{n}/h^{-p/2}\to|N(0,Var)|$ in distribution, 2. (b)

$\delta_{n}=n^{-1/2}h^{-p/4}$, $T_{n}/h^{-p/2}\to|N(\mu_{1},Var)|$ in distribution, 3. (c)

$n^{1/2}h^{p/4}\delta_{n}\to\infty$, $T_{n}/(n\delta_{n}^{2})\to\mu_{1}$ in probability.

Remark 3**.**

Together with Theorem 4, the results are informative and important which shows that feature $f_{3}$) can be achieved. As commented in Remark 1, Zheng (1996)’s test and the tests in Lavergne and Patilea (2012) and Guo et al. (2016) diverge to infinity at the nonparametric rate of order $nh^{p/2}$ and of $nh^{1/2}$. Theorem 4 shows that although under the global alternative $T_{n}$ uses the lobal smoothing test proposed by Zheng (1996) up to the standardizing constant, it can diverge to infinity at the fastest possible rate of order $n$. Further, result 1) shows that the nonparametric estimation-based component $V_{1}$ is used, $T_{n}$ can still detect the local alternatives distinct from the null at the fastest possible rate $1/\sqrt{n}$. Result 2) gives a full picture to show that when $\delta_{n}$ converges to zero slower, the test statistic can have different weak limits at different rates. We can see that in case (b), $T_{n}$ diverges to infinity at the rate of $h^{-p/2}$ whereas Zheng (1996)’s test goes to a finite limit. All these results demonstrate that the special construction of the hybrid makes the nonparametric estimation-based component behave like a global smoothing test. This is a significant contribution of the method. The numerical studies reported later will also justify the power performance improvement over existing local smoothing tests.

5 Numerical studies

5.1 Simulations

In this section, we conduct several simulation studies to examine the performance of the proposed test through comparisons with several typical local and global smoothing tests in different scenarios.

The proposed test is written as $T_{n}$. The typical local smoothing tests are compared: the kernel estimation-based test $T^{Zh}$ proposed by Zheng (1996) and the kernel estimation-based adaptive-to-model test $T^{GWZ}$ in Guo et al. (2016). Here we make a slight modification of the Guo et al. (2016)’s test such that it can be applied to the general parametric model considered in this paper other than the parametric single-index model in Guo et al. (2016). To be precise, assume $q_{1}=dim(S_{Y{|}X})$ and $Span\{B_{1}\}=S_{Y{|}X}$ with the rank of $B_{1}$ equal to $q_{1}$. Without confusion, we still write the modified Guo et al. (2016)’s test as $T^{GWZ}$ in the following:

[TABLE]

where $\hat{q}_{1}$ and $\hat{B}_{1}$ are the corresponding estimate of $q_{1}$ and $B_{1}$. In the experiments, we use the DEE-SIR method and the BIC-type criteria used in Guo et al. (2016) to estimate the target matrix $B_{1}$ and its rank $q_{1}$. Furthermore, we also compare it with two typical global smoothing tests: the ICM test $T^{B}$ proposed by Bierens (1982) and the empirical process-based test $T^{SMQ}$ proposed by Stute et al. (1998).

We want to provide the following information on the performance. First, note that in our test, the second component $V_{1}$ is exactly the same as that in Zheng (1996) except for different standardizing constant $n/\hat{\sigma}_{0}^{2}$. Guo et al. (2016) has well demonstrated that Zheng (1996) suffers from the curse of dimensionality severely and $T^{GWZ}$ fully utilizes the dimension reduction structure with much better power performance than $T^{Zh}$ when the dimension is high. Thus, through the comparison with $T^{GWZ}$, we want to check how the hybrid can well overcome the shortcomings $T^{Zh}$ has. Second, we choose two global smoothing tests to compare as we want to see how the hybrid can share the nice features of global smoothing tests under low-frequency alternatives and can acquire higher power under high-frequency alternatives. Based on these considerations, we design the following studies.

Study 1 and Study 2 consider some nonlinear null models against high-frequency and low-frequency alternative model respectively. In Study 1, we also examine the influence of dependency among the components of $X$ to the performances of the competitors. These studies focus on the scenarios with $q=1$ under the alternatives. In Study 3, we consider a model with $q=p-1$ under the alternatives. The null models in these $3$ studies are all parametric single-index. To further check the performance of the hybrid test, in Study 4, we design a model without dimension reduction structure under the null hypothesis. Every experiment in the simulations is repeated $1,000$ times to compute the empirical sizes and powers. The significance level is set to be $0.05$.

Study 1:  Consider a model as:

[TABLE]

The parameters are $\alpha=(1,\dots,1,0,\dots,0)/\sqrt{p/2}$ and $\beta=(0,\dots,0,1,\dots,1)/\sqrt{p/2}$, where $\|\alpha\|=\|\beta\|=1$. The value $a=0$ corresponds to the null hypothesis and $a\not=0$ to the alternatives. The sample size $n=200$ and $p=2,8$ to respond low and high-dimensional scenarios respectively. We do not consider higher dimension only because when $p=8$, local smoothing tests have already clearly shown their poor power performance as Guo et al. (2016) demonstrated. We consider two different distributions of $X$ with independent and dependent components: $N(0,I_{p})$ and $N(0,\tilde{\Sigma})$ where the $(i,j)$-th element of $\tilde{\Sigma}$ is $0.5^{|i-j|}$. Empirical sizes and powers with independent $X$ are shown in Figure 1. Table 1 displays the results with $X\sim N(0,\tilde{\Sigma})$.

For $X$ with independent components, when $p=2$, the hybrid test and Bierens (1982)’s test perform best and when $p=8$, our test shows a great advantage over the other competitors. For the tests in Zheng (1996), Bierens (1982) and Stute et al. (1998), their powers are severely deteriorated by the dimension $p$. Guo et al. (2016)’s test keeps a stable performance when $p$ goes from $2$ up to $8$. This makes sense as the model structure is single-index, which is in favor of their test.

For $X\sim N(0,\tilde{\Sigma})$, again the proposed test and Bierens (1982)’s test are the winners when $p=2$. Our test is compatible with the global smoothing test in Stute et al. (1998), which works best with $p=8$. Although Bierens (1982)’s test has a higher power than Stute et al. (1998)’s, its empirical size is far away from the significance level of $0.05$. Thus, the significance level maintenance is an issue for this test.

The results with $X\sim N(0,I_{p})$ and $X\sim N(0,\tilde{\Sigma})$ suggest that our test is robust to the dependency between components of $X$. Guo et al. (2016)’s test has a better performance with $X\sim N(0,\tilde{\Sigma})$ than $X\sim N(0,I_{p})$ when $p=2$. However, its power drops quickly with increasing dimension when the components of $X$ are correlated. Stute et al. (1998)’test performs much better under $X\sim N(0,\tilde{\Sigma})$ than under $X\sim N(0,I_{p})$.

The comparison with Zheng (1996)’s test clearly shows the advantage of the hybrid even when Zheng (1996)’s test is, up to the standardizing constant, its nonparametric conditional comment component.

Study 2: Consider a nonlinear null model against low-frequency alternative models:

[TABLE]

The parameters $\alpha$ and $\beta$ follows the previous settings. $X$ and ${\varepsilon}$ are generated from $N(0,I_{p})$ and $N(0,1)$ respectively. Simulation results are presented in Figure 2.

The hybrid test performs best among the competitors in both cases with $p=2$ and $p=8$. Except for the proposed test and Guo et al. (2016)’s test, the others suffer from severe power lose due to the dimension increasing. Bierens (1982)’s test is compatible with our test when $p=2$, but deviates far from the significance level and low power when $p=8$. Zheng (1996)’s test, as discussed before, slowly converges to its weak limit when $p=8$, which may be the cause for low power. Guo et al. (2016)’s test performs much better than Zheng (1996)’s when $p=8$. Power of Stute et al. (1998)’s test declines sharply as $p$ increases up to $8$.

Study 3:  In the above two studies, the indicative dimension $q=1$ under the alternatives. In this study, we conduct experiments for models with $q=p-1>1$ under the alternatives. Let $X_{i}$ be the $i$-th component of $X$. Consider

[TABLE]

The sample size $n=200$ and dimension $p=8$. The explanatory vector $X$ and the error ${\varepsilon}$ are drawn independently from $N(0,I_{p})$ and $N(0,1)$. Table 2 presents the empirical sizes and powers.

The results show that the proposed test and Guo et al. (2016)’s test work better than the others. When the deviation from the null hypothesis is small, which corresponds to $a=0.2,0.4,0.6,0.8$, the proposed test constantly surpasses the other competitors. This is consistent with the property that our test can better detect the local alternatives. Compared with the others, Guo et al. (2016)’s test also works well, and when $a=1$, its power even suddenly jumps up to be higher than ours.

Under the null hypothesis, the models in Study 1 - Study 3 are all single-index where $dim(S_{Y{|}X})=1$ under the null hypothesis. To make the comparison more comprehensive, we consider a null model with $dim(S_{Y{|}X})>1$ in the following study.

Study 4:  Recall that $X_{i}$ represents the $i$-th component of $X$. Consider a null model without a dimension reduction structure:

[TABLE]

The unknown parameters in the null model are $(c_{1},c_{2},c_{3},c_{4},c_{5})=C/\|C\|$ where $C=(1,1/2,1/3,1,1)$. Figure 3 displays the empirical sizes and powers of the all competitors.

From the results, we can see that the proposed test has a great advantage over the competitors, especially under the scenarios with small $a=0.2,0.4,0.6$. Zheng (1996)’s test and Stute et al. (1998)’s test does not work well with this model. The power of Guo et al. (2016)’s test is also close to $1$ when $a=1$, but its ability to detect the alternatives with small $a$ is not comparable with the proposed test. Again, the empirical size of Bierens (1982)’s test is too conservative to maintain the significance level.

5.2 A real data example

The data of real estate price were collected from June 2012 to May 2013, see Yeh and Hsu (2018). This data set contains $414$ samples during this period. Let $Y$ be the price per unit area of the interested real estate. Analysis in Yeh and Hsu (2018) shows that the important attributes involve the date of transaction ($X_{1}$), the house age ($X_{2}$), the distance to the nearest metro station ($X_{3}$), the number of convenience stores within walking distance ($X_{4}$) and the geographic location ( denoted by northern latitude $X_{5}$ and eastern longitude $X_{6}$). The transaction dates are transformed into real numbers. For instance, if a house was traded in September 17, 2012, then its transaction date is presented as $2012.917$. Compared with traditional comparative approaches, regression analysis costs lower and is more reasonable since it reduces the bias caused by subjectivity. The following regression model was applied to fit this real estate price dataset:

[TABLE]

They also compared the advantages and disadvantages of various methods to appraise the real estate price and concluded that the above regression model is then preferred by industry for its interpretability and low cost in large scale evaluation. Now we conduct a test to verify the adequacy of this model. The value of the test statistic is $1.0249$ and the $p$-value is $0.3114$. Therefore the model is plausible in describing the relationship between the target real estate price and the important attributes. The residual plot in Figure 4 shows that there is no obvious nonlinear pattern in the residuals, which also implies this conclusion.

6 Discussions

In this research, we propose a hybrid test that can have model adaptation property. This is the first test in this research field, to the best of our knowledge, to inherit all main advantages of nonparametric estimation-based and empirical process-based tests in the literature. Further, this methodology is ready to be applied to other testing problems with more complicated data.

It is interestingly observed that we use two tests that are definitely not powerful to construct a very powerful test, even more powerful than any existing tests in the sense that it can handle both low-frequency and high-frequency alternative models. Thus, the most significant contribution of this method is that under the alternatives, the hybrid can make a local smoothing test share all appealing features global smoothing tests have. Note that the test under the null is just for critical value determination. Thus, we should use a test that is simple such that the significance level can be easily maintained. For the test under the alternatives, it could be any nonparametric estimation-based test. A natural question, though beyond the scope of this paper, is whether we can choose a proved powerful test to make the test more powerful. It could be possible because the hybrid can consist of any two tests and thus deserves further study. On the other hand, we also see from the results in Sections 2 and 4 that for any nonparametric estimation-based test, the asymptotic properties of the resulting hybrid test should be similar, whereas the pair of the tests are moment-based and empirical process-based should make the resulting test only satisfy features $f_{1}$ and $f_{3}$, but not feature $f_{2}$ as both are global smoothing tests.

From the proof for Theorem 5, we can find the method of dimensionality determination has the following property:   When the null model is linear, even under the local alternatives (16), the estimation $\hat{q}$ equals [math] with a probability going to $1$.  In other words, when the null model is linear, the dimension $\hat{q}$ cannot indicate the local alternatives designed in our paper. We find that this is because of the target matrix construction. For a linear null model, we have $g(X,\theta)=\theta^{\top}X$ and thus $\dot{g}=X$. Together with $m(t)=E(X\exp(it\eta))$ in Section 3.2, it will lead to $\Delta_{\ell}=E(X\ell)-E(X\dot{g}^{\top})E(\dot{g}\dot{g}^{\top})^{-1}E(\dot{g}\ell)=0$, as shown by the proof in the supplementary materials. To avoid this problem, we may use a more general function $L(X)$ in lieu of $X$ such that $\Delta_{\ell}=E(L(X)\ell)-E(L(X)\dot{g}^{\top})E(\dot{g}\dot{g}^{\top})^{-1}E(\dot{g}\ell)\not=0$. We point out this problem here, but in this paper we still use the target matrix provided in Section 3.2 for simplicity. How to choose a target matrix such that the dimensionality determination works better deserves future study.

Another issue is about the extension of the method to heteroscedastic models. The current approach has a main limitation that only homoscedastic models can be handled. From the results, we can see that the key to handling heteroscedastic null models is to define an indicative dimension and to have a method to identify it. The method currently relies on the independence between the error and covariates. For heteroscedastic models, the independence no longer holds. Thus, in a special case where $\varepsilon=b(X)\epsilon$ and $X$ is independent of $\epsilon$, we may estimate $b$ and then use $\eta/\hat{b}$ instead of $\eta$ to construct the test. For more general paradigms, it deserves further study.

SUPPLEMENTARY MATERIAL

Supplemenraty of adaptive-to-model hybrid test for regressions

Technical details of TDRR and proofs of the theorems. (.pdf file)

Appendix A Appendix

Regularity conditions

Condition 1**.**

$\{(x_{i},y_{i})\}_{i=1,\dots,n}$* are i.i.d. random samples from $(X,Y)$ in $R^{p}\times R$ and $EY^{2}<\infty$.*

Condition 2**.**

The parameter space $\Theta$ is compact and convex.

Condition 3**.**

The regression function $g(x,\theta)$ is a Borel measurable real function on $R^{p}$ for each $\theta$ and is twice continuously differentiable with respect to $\theta$ for each $x$.

Condition 4**.**

Let $\|\cdot\|$ represent the Euclidean norm.

[TABLE]

Condition 5**.**

There exists a unique minimizer $\theta^{\ast}$ such that

[TABLE]

Under the null hypothesis, $\theta^{\ast}$ is an interior point of $\Theta$.

Condition 6**.**

The matrix $E\{\dot{g}(X,\theta^{\ast})\dot{g}(X,\theta^{\ast})^{\top}\}$ is nonsingular.

Condition 7**.**

The independent random vector $X$ satisfies the linearity condition: $E(X{|}B^{\top}X)=P_{B}^{\top}(\Sigma)X$, where $B\in R^{p\times q}$ is any basis of $S_{\eta{|}X}$.

Condition 8**.**

Let $\tilde{m}_{1}(t)=E(it\exp(it{\varepsilon})X\ell)$, $\tilde{m}_{2}(t)=E(itX\dot{g}^{\top}\exp(it{\varepsilon}))G^{-1}E(\dot{g}\ell)$ and $\tilde{m}(t)=\tilde{m}_{1}(t)-\tilde{m}_{2}(t)$. The matrix $\tilde{M}=\int\tilde{m}(t)\tilde{m}^{H}(t)dF_{\varepsilon}(t)$ is nonsingular.

Condition 9**.**

The kernel function $K(\cdot)$ is a symmetric nonnegative and continuous function which is bounded and $\int K(u)du=1$.

Condition 10**.**

The bandwidth in kernel estimation satisfies $h\to 0$ and $nh^{p}\to\infty$ as $n\to 0$.

Conditions $1-6$ are commonly required for the consistency and asymptotic normality of the least squares estimation for the parameter $\theta$, see Zheng (1996). Condition $7$ is the linearity condition for the independent vector $X$ such that the target matrix $M$ is contained in the central subspace $S_{\eta{|}X}$. Condition $8$ is assumed for the target matrix such that the indicative dimension can be well identified and then can be used for Theorems 5 and 6. Conditions $9$ and $10$ are commonly used conditions on the kernel estimation and the asymptotic normality of the corresponding test statistic.