Computationally efficient goodness-of-fit tests for the error distribution in nonparametric regression

REVSTAT – Statistical Journal Volume 16, Number 1, January 2018, 137–166 COMPUTATIONALLY EFFICIENT GOODNESS-OF-FIT TESTS FOR THE ERROR DISTRIBUTION IN NONPARAMETRIC REGRESSION Authors: G.I. Rivas-Mart´ ınez – Laboratorio de Sistemas de Potencia y Control, Universidad Nacional de Asunci´on, Paraguay [email protected] M.D. Jim´ enez-Gamero – Dpto. de Estad´ıstica e Investigaci´on Operativa, Universidad de Sevilla, Espa˜na [email protected] Received: April 2016 Revised: August 2016 Accepted: August 2016 Abstract: •Several procedures have been proposed for testing goodness-of-fit to the error distribution in nonparametric regression models. The null distribution of the associated test statistics is usually approximated by means of a parametric bootstrap which, under certain conditions, provides a consistent estimator. This paper considers a goodness-of-fit test whose test statistic is an L2norm of the difference between the empirical characteristic function of the residuals and a parametric estimate of the characteristic function in the null hypothesis. It is proposed to approximate the null distribution through a weighted bootstrap which also produces a consistent estimator of the null distribution but, from a computational point of view, is more efficient than the parametric bootstrap. Key-Words: •goodness-of-fit; empirical characteristic function; regression residuals; weighted bootstrap; consistency. AMS Subject Classification: •62G08, 62G09, 62G10. 138 G.I. Rivas-Mart´ınez and M.D. Jim´enez-Gamero Computationally Efficient Goodness-of-Fit Tests... 139 1. INTRODUCTION Let (X, Y ) be a bivariate random vector satisfying the general nonparametric regression model (1.1) Y=m(X) + σ(X)ε, where m(x) = E(Y|X=x) is the regression function, σ2(x) = Var(Y|X=x) is the conditional variance function and εis the regression error, which is assumed to be independent of X. Note that, by construction, E(ε) = 0 and Var(ε) = 1. The covariate Xis continuous with density function fX. The regression function, the variance function, the error distribution and that of the covariate are unknown and no parametric models are assumed for them. Because the knowledge of the error distribution will improve the statistical analysis of model (1.1), several authors have proposed tests for such distribution, that is, tests of the null hypothesis H0:F∈ F, versus the alternative H1:F /∈ F, where Fstands for the cumulative distribution function (CDF) of εand Fis a parametric family, F=F(·;θ), θ ∈Θ,Θ⊆Rp. Examples are the tests in Neumeyer et al. [17] and Heuchenne and Van Keilegom [6], which are based on comparing the empirical CDF of the residuals to a parametric estimator of the CDF under the null hypothesis. Since the equality of the CFDs can be also interpreted in terms of the associated characteristic functions (CFs), Huˇskov´a and Meintanis [11] have proposed a test for H0that is based on comparing the empirical CF of the residuals to a parametric estimator of the CF under the null hypothesis. As commented in Jim´enez-Gamero [13], it is interesting to observe that the last paper requires weaker conditions for the validity of the procedures than the ones based on the CDF. Nevertheless, in all cases the limit distribution of the proposed test statistics is unknown, even under the null distribution, because it depends on the unknown value of the parameter θ. To overcome this difficulty, these papers propose to use a parametric bootstrap (PB) for approximating the null distribution of the test statistic. Although very easy to implement, the PB can become very computationally expensive as the sample size and/or the number of unknown parameters increase. This paper studies another method for estimating the null distribution of the test statistic Tn,w(ˆ θ) in [11]. Specifically, a weighted bootstrap (WB) approximation in the sense of Burke [2] is considered (see also Zhu [23]). This method 140 G.I. Rivas-Mart´ınez and M.D. Jim´enez-Gamero has been previously suggested in Kojadinovic and Yan [15], to approximate the null distribution of goodness-of-fit (GOF) tests based on the empirical CDF, and in Jim´enez-Gamero and Kim [14], to approximate the null distribution of GOF tests based on the empirical CF (ECF), among others. Both papers assume observable independent and identically distributed (IID) data. They show that the properties of the WB are quite similar to those of the PB (it provides a consistent estimator of the null distribution and the resulting test is able to detect any alternative) but, from a computational point of view, it is more efficient. In view of the good properties of the WB in these and other papers, it is also expected to work satisfactorily for estimating the null distribution of the test statistic considered in this paper. The purpose of the current study is to investigate, both theoretically and empirically, the use of the WB for approximating the null distribution of Tn,w(ˆ θ). A main difference between the setting in this paper and the one in [14, 15] is that in our case the errors are not observable. So we replace the errors by the residuals, but the residuals are not independent. The paper is organized as follows. Section 2 describes the test statistic and explains some problems with the WB approximation. Section 3 gives a solution to the problems described in the previous section and proves the consistency of the proposed WB approximation. It also shows that the resulting test is consistent, in the sense of being able to detect any alternative. The application of the proposed WB approximation requires the estimation of certain functions appearing in the linear expansion of the parameter estimators. The estimation of such functions is dealt with in Section 4. Section 5 reports the results of some simulation experiments designed to study the finite sample performance of the proposed approximation and to compare it to the PB. From this numerical study it is concluded that both approximations behave quite closely but, from a computational point of view, the WB outperforms the PB. Section 6 concludes and outlines possible extensions of the results presented in this paper. All proofs and technical details are deferred to the last section. The following notation will be used along the paper: all vectors are column vectors; for any vector a,akdenotes its k-th coordinate and kakits Euclidean norm; the superscript Tdenotes transpose; Eθand Pθdenote expectation and probability, respectively, assuming that the data has CDF F(·;θ); P∗denotes the conditional probability law, given the data; all limits in this paper are taken when n→ ∞;L →denotes convergence in distribution; P →denotes convergence in probability; a.s. →denotes the almost sure convergence; for any complex number z=a+ ib,|z|is its modulus; an unspecified integral denotes integration over the whole real line R; for a given non-negative real-valued function wwe denote k·kwto the norm and h·,·iwto the scalar product in the Hilbert space L2(w) = {g:R→C,R|g(t)|2w(t)dt < ∞}; if Fis a CDF, then L2(F) = {g:R→C,R|g(t)|2dF(t)<∞}; for any real function f(t;θ) differentiable at t∈Rand at θ= (θ1, θ2, ..., θp)T∈Rpthe following notations will be Computationally Efficient Goodness-of-Fit Tests... 141 used: f′(t;θ) = ∂ ∂tf(t;θ), f(r)(t;θ) = ∂ ∂θr f(t;θ),1≤r≤p, ∇f(t;θ) = f(1)(t;θ), f(2)(t;θ), ..., f(p)(t;θ)T. 2. THE TEST STATISTIC Let (X1, Y1), ..., (Xn, Yn) be IID from model (1.1), that is, Yj=m(Xj) + σ(Xj)εj, 1 ≤j≤n. Since the hypothesis H0is on the common error distribution, ε1, ..., εn, and the errors are not observable, the inference must be based on the residuals, ˆεj=Yj−ˆm(Xj) ˆσ(Xj),1≤j≤n, where ˆm(·) and ˆσ(·) are estimators of m(·) and σ(·), respectively. Several choices are possible for ˆm(·) and ˆσ(·). Here, as in [11], we use the following kernel estimators for the density function fXof X, the regression function m(·) and the variance function σ2(·), ˆ fX(x) = 1 n n X j=1 Khn(Xj−x), ˆm(x) = 1 nˆ fX(x) n X j=1 Khn(Xj−x)Yj, ˆσ2(x) = 1 nˆ fX(x) n X j=1 Khn(Xj−x){Yj−ˆm(x)}2, where Khn(·) = 1 hnK(· hn), K(·) is a kernel and hnis the bandwidth, satisfying certain conditions that will be specified later. Huˇskov´a and Meintanis [11] proposed the following test for testing H0, Ψ = (1,if Tn,w(ˆ θ)≥tn,ω,α, 0,otherwise, where tn,ω,α is the 1 −αpercentile of the null distribution of Tn,ω(ˆ θ), (2.1) Tn,ω(ˆ θ) = nZ|cn(t)−c(t, ˆ θ)|2ω(t)dt =nkcn(t)−c(t, ˆ θ)k2 w, cn(t) is the ECF of the residuals, cn(t) = 1 n n X j=1 exp(itˆεj) = 1 n n X j=1 cos(tˆεj) + i 1 n n X j=1 sin(tˆεj), 142 G.I. Rivas-Mart´ınez and M.D. Jim´enez-Gamero c(t;θ) is the CF associated to F(ε;θ), that is, c(t;θ) = Eθ{exp(itε)}=R(t;θ) + iI(t;θ), ω(t) is a nonnegative function such that Rω(t)dt < ∞, which may depend on θ, and ˆ θis a consistent estimator of θsatisfying the following assumption. (A.1) Under H0,√n(ˆ θ−θ0) = 1 √n n X j=1 ψ(εj;θ0) + op(1), where θ0is the true parameter value, Eθ0{ψ(εj;θ0)}= 0 and Eθ0kψ(εj;θ0)k2<∞. Assumption (A.1) implies that, when the null hypothesis is true and θ0 denotes the true parameter value, √n(ˆ θ−θ0) is asymptotically normally distributed. This assumption is satisfied by commonly used estimators such as maximum likelihood estimators and method of moment estimators when ε1, ..., εn are observable and, in such a case, the expression of the function ψis well-known (see, for example, [1, Ch. 5]). In our setting, the errors are not observable and the expression of the function ψdiffers from the observable case. This topic will be discussed in detail in Section 4. Theorem 1 in [11] states that if ˆ θsatisfies (A.1), H0is true and θ0is the true parameter value, under certain additional conditions (assumptions (A.2)–(A.7) in Section 7), (2.2) Tn,ω(ˆ θ)L −→ kZ(t;θ0)k2 ω, where {Z(t;θ0), t ∈R}is a centered Gaussian process on L2(ω) with covariance structure of the form Covθ0{Z1(ε;t, θ0, ψ), Z1(ε;s, θ0, ψ)}, (2.3) Z1(ε;t, θ, ψ) = cos(tε) + sin(tε)−R(t;θ)−I(t;θ)−tε{R(t;θ)−I(t;θ)} −tε2−1 2{R′(t;θ) + I′(t;θ)}−ψT(ε;θ){∇R(t;θ) + ∇I(t;θ)}. Clearly, the asymptotic null distribution of Tn,ω(ˆ θ) is unknown. It depends on the hypothetical the error distribution, on the chosen estimator and the true unknown value of the parameter. In order to try to approximate the null distribution of Tn,ω(ˆ θ) we first observe that it resembles a degree-2 V-statistic, because Tn,ω(ˆ θ) = 1 n n X j=1 n X k=1 ρ(ˆεj,ˆεk;ˆ θ), with ρ(ε, z;θ) = u(ε−z)−u0(ε;θ)−u0(z;θ)+u00(θ), u0(ε;θ) =Ru(ε−z)dF(z;θ), u00(θ) = Ru(ε−z)dF(ε;θ)dF(z;θ), and u(t) = Rcos(tε)ω(ε)dε. Dehling and Mikosch [4] (see also Huˇskov´a and Janssen [10]) showed that if ε1, ..., εnare IID, ξ1, ..., ξnare IID with E(ξ1) = 0 and Var(ξ1) = 1, independent of ε1, ..., εnand Vn=1 n2P1≤j,k≤ng(εj, εk) is a degenerate degree-2 V-statistic, Computationally Efficient Goodness-of-Fit Tests... 143 then the conditional distribution, given ε1, ..., εn, of 1 nX 1≤j,k≤n g(εj, εk)ξjξk consistently estimates that of nVn. In the light of this result, since ˆεjand ˆ θ are approximations to εjand θ, respectively, one may be tempted to estimate the null distribution of Tn,ω(ˆ θ) by means of the conditional distribution, given (X1, Y1), ..., (Xn, Yn), of (2.4) W∗=1 nX 1≤j,k≤n ρ(ˆεj,ˆεk;ˆ θ)ξjξk. We will see that this approach is wrong. The next result gives the limit distribution of W∗. The required assumptions are listed in Section 7. Theorem 2.1. Suppose that kˆ θ−θ1k=op(1), for some θ1∈Θ, that assumptions (A.2)–(A.6) hold, that the first partial derivatives R(r)(t;θ),I(r)(t;θ), 1≤r≤p, exist and are continuous functions ∀θ∈U(θ1)⊆Θ, an open neighborhood of θ1, and they are bounded by functions in L2(ω),∀θ∈U(θ1), then sup x|P∗{W∗≤x}−P{W0≤x}| P −→ 0, where W0=kZ0(t;θ1)k2 ω,{Z0(t;θ1), t ∈R}is a centered Gaussian process on L2(ω)with covariance structure of the form Cov{Z0(ε;t, θ1), Z0(ε;s, θ1)}, Z0(ε;t, θ) = cos(tε) + sin(tε)−R(t;θ)−I(t;θ). From the result in Theorem 2.1 and (2.2), it is clear that the conditional distribution of W∗does not provide a consistent estimator of the null distribution of Tn,ω(ˆ θ) because replacing m(·), σ(·) and θby ˆm(·), ˆσ(·) and ˆ θ, respectively, has an impact on the asymptotic null distribution of the test statistic that is not captured by the conditional distribution of W∗. The next Section shows how to deal with this problem. Before ending this section we do some comments on the behaviour of ˆ θ under the alternative. Theorem 2.1 assumes that ˆ θhas a limit (in probability), θ1. In practice, to estimate θone proceeds as if H0were true. For example, θis usually estimated by its quasi maximum likelihood estimator, which maximizes the likelihood under the null hypothesis (with the errors replaced by the residuals). If H0is true, under certain assumptions, the resulting estimator converges to the true parameter value (see Section 4); if H0is not true, then proceeding as in White [22] for observable data, it can shown that, under certain conditions, the estimator also converges to a well-defined limit. Similar comments could be done for other estimators. 144 G.I. Rivas-Mart´ınez and M.D. Jim´enez-Gamero 3. CONSISTENCY OF THE WB APPROXIMATION If assumptions (A.1)–(A.7) hold and H0is true, from the proof of Theorem 1 in [11], it follows that (3.1) Tn,ω(ˆ θ) = T1,n,ω(θ0) + op(1), where T1,n,ω(θ) = k1 √n n X j=1 Z1(εj;t, θ, ψ)k2 ω, with Z1(ε;t, θ, ψ) as defined in (2.3). Now, from (3.1) and applying the results in [4], we get that the conditional distribution, given (X1, Y1), ..., (Xn, Yn), of T∗ 1,n,ω(θ0) = k1 √n n X j=1 Z1(εj;t, θ0, ψ)ξjk2 ω, provides a consistent estimator of the distribution of Tn,ω(ˆ θ), when H0is true. From a practical point of view, this result is useless because Z1(εj;t, θ0, ψ) depends on the non-observable error εj, on the unknown value of θ0and on the function ψ(εj;θ0), whose explicit expression is usually unknown. Suppose that kˆ θ−θ1k=op(1), for some θ1∈Θ, θ1being the true parameter value if H0is true. To overcome these difficulties we replace εjby ˆεj,θ0by ˆ θand ψ(εj;θ0) by ψn(ˆεj;ˆ θ), where ψn(·;ˆ θ) is a function of the data which approximates ψin such a way that (3.2) 1 n n X j=1 kψn(ˆεj;ˆ θ)−ψ1(εj;θ1)k2P −→ 0, with E{kψ1(ε;θ1)k2}<∞and ψ1(ε;θ1) = ψ(ε;θ1) if H0is true. The choice of ψnwill depend on ψ, that is, on the estimator of θconsidered. Section 4 studies some proposals for ψnsatisfying (3.2) for two common choices for ˆ θ: the maximum likelihood estimator and the method of moments estimator, both based on the residuals. So, the null distribution of Tn,ω(ˆ θ) is now estimated by means of the conditional distribution, given (X1, Y1), ..., (Xn, Yn), of T∗ 2,n,ω(ˆ θ) = k1 √n n X j=1 Z1(ˆεj;t, ˆ θ, ψn)ξjk2 ω. The next theorem gives the limit of the conditional distribution of T∗ 2,n,ω(ˆ θ), given (X1, Y1), ...(Xn, Yn). Theorem 3.1. Suppose that kˆ θ−θ1k=op(1), for some θ1∈Θ,θ1being the true parameter value if H0is true, and that assumptions (A.1)–(A.7) and (3.2) hold, then sup xP∗nT∗ 2,n,ω(ˆ θ)≤xo−P{T2≤x} P −→ 0, Computationally Efficient Goodness-of-Fit Tests... 145 where T2=kZ2(t;θ1)k2 ω,{Z2(t;θ1), t ∈R}is a centered Gaussian process on L2(ω)with covariance structure of the form Cov{Z1(ε;t, θ1, ψ1), Z1(ε;s, θ1, ψ1)}. The result in Theorem 3.1 is valid whether the null hypothesis H0is true or not. An immediate consequence of this fact and (2.2) is the following. Corollary 3.1. If H0is true and the assumptions in Theorem 3.1 hold, then sup xP∗nT∗ 2,n,ω(ˆ θ)≤xo−Pθ1nTn,ω(ˆ θ)≤xo P −→ 0. Let α∈(0,1) and Ψ∗=(1,if Tn,ω(ˆ θ)≥t∗ 2,n,ω,α, 0,otherwise, where t∗ 2,n,ω,α is the 1 −αpercentile of the conditional distribution of T∗ 2,n,ω(ˆ θ), or equivalently, Ψ∗= 1 if p∗≤α, where p∗=P∗nT∗ 2,n,ω(ˆ θ)≥Tn,ω(ˆ θ)obsoand Tn,ω(ˆ θ)obs is the observed value of the test statistic. The result in Corollary 3.1 states that Ψ∗is asymptotically correct, in the sense that its type I error is asymptotically equal to the nominal value α. Corollary 3.2. Suppose that H0is not true and let c(t)denote the true CF of the errors. If the assumptions in Theorem 3.1 hold and ωis such that (3.3) κ=kc(t)−c(t;θ1)k2 ω>0, then P(Ψ∗= 1) →1. Corollary 3.2 shows that, if ωis such that (3.3) holds, then the test Ψ∗is consistent in the sense of being able to asymptotically detect any (fixed) alternative. Since two distinct characteristic functions can be equal in a finite interval (Feller [5, p.506]), a general way to ensure (3.3) is to take ωpositive for almost all (with respect to the Lebesgue measure) points in R. Remark 3.1. If model (1.1) is homoscedastic, that is, if σ(x) = σ,∀x, for some unknown σ > 0, we can use the residuals ˜εj=Yj−ˆm(Xj), 1 ≤j≤n, and consider σas a parameter of the family F. In this framework, the result in Theorem 3.1 (with weaker assumptions) keeps on being true with the following simpler expression for Z1(ε;t, θ, ψ), Z1(ε;t, θ, ψ) = cos(tε)−R(t;θ) + sin(tε)−I(t;θ)−tεR(t;θ) + tεI(t;θ) −ψT(ε;θ){∇R(t;θ) + ∇I(t;θ)}. 152 G.I. Rivas-Mart´ınez and M.D. Jim´enez-Gamero distribution with 3 degrees of freedom (denoted as χ2 3) and the Student tdistribution with 5 degrees of freedom (denoted as t5). All aforementioned distributions were conveniently centered and scaled to have mean 0 and variance 1. To approximate the p-value, 1000 replications were generated for both the PB and the WB. For the WB, the raw multipliers and the centered multipliers were considered, denoted by WB1 and WB2 in the tables, respectively. The multipliers were generated from a univariate standard normal distribution. As for the asymptotic distribution (when the null hypothesis is simple, denoted as A in the tables), it is rather difficult to calculate because it coincides with that of Pj≥1λjχ2 1,j, where χ2 1,1, χ2 1,2, ... are independent chi-squared variables with one degree of freedom, the set {λj, j ≥1}are the non-null eigenvalues of the integral equation RC(t, s)Gj(t)dt =λjGj(s), with corresponding eigenfunctions {Gj(·), j ≥1},C(t, s) is the covariance kernel of Z1(ε;t) (see Remarks 3.2 and 3.3 for the expression of Z1(ε;t)), and determining the eigenvalues of an integral equation is tricky. Because of this reason, we approximated it by generating 10,000 samples of size 1000 obeying H0and calculated the test statistic at each sample, obtaining 10,000 values. The empirical CDF of these 10,000 values was taken as an approximation to the asymptotic null distribution. 1000 samples with size n= 25 were generated from each distribution and the fractions of p-values less than or equal to 0.05 and 0.1 were calculated. The experiment was repeated for n= 50,100. 5.1. Homoscedastic model The reported results correspond to the model Yj=Xj+X2 j+εj,1≤j≤n, where Xjfollows the uniform (0,1) distribution. We first considered that θis known. Since the model is homoscedastic and the null hypothesis is simple, the simplifications in Remark 3.3 can be applied. Table 1 displays the results obtained for the type I error and the power for testing normality and Table 2 for testing GOF to the Laplace distribution. Looking at these tables it can be concluded that, in terms of type I error, both the PB and the WB behave very close to the nominal levels, while the asymptotic approximation is a bit conservative, specially for testing GOF for the Laplace distribution. As for the power, the test based on the WB approximation seems to be a bit more powerful than one based on the PB. In most cases (all but alternatives βand χ2 3in Table 2) the WB approximation is also more powerful than one based on the asymptotic approximation. Tables 3 and 4 show the results when θis assumed to be unknown. In this case, the simplifications in Remark 3.1 can be applied. Looking at these tables it Computationally Efficient Goodness-of-Fit Tests... 153 can be concluded that, in terms of the type I error, as before, both the PB and the WB behave very close to the nominal levels. As for the power, for n= 25,50 in some cases the WB is more powerful than the PB, but in others cases the opposite is observed; for n= 100 the test based on the WB approximation seems to be a bit more powerful than one based on the PB. Table 1: (Homoscedastic model, simple null hypothesis) Percentage of rejections for the normality null hypothesis at the significance levels 5% (upper entry) and 10% (lower entry). n= 25 n= 50 n= 100 A PB WB1 WB2 A PB WB1 WB2 A PB WB1 WB2 N3.60 6.10 4.10 6.40 4.00 5.00 4.12 4.84 5.20 4.74 4.12 4.74 8.20 11.50 10.20 12.30 9.00 10.04 9.24 10.48 10.20 9.64 9.34 10.40 LP 25.50 36.10 57.80 64.80 45.40 56.30 86.60 88.30 76.60 77.70 98.90 99.00 35.90 48.70 70.60 74.70 57.40 68.50 90.40 91.00 83.60 83.20 99.60 99.70 LG 10.30 57.60 56.40 63.10 12.70 88.10 87.30 89.30 17.80 99.90 100.00 100.00 18.10 70.40 72.00 76.00 20.60 93.20 94.30 95.10 27.80 99.90 100.00 100.00 G18.40 33.50 45.80 52.00 36.70 61.80 87.30 89.30 71.70 90.70 100.00 100.00 30.60 46.30 62.80 67.70 49.80 74.40 94.30 95.10 81.80 96.70 100.00 100.00 β54.10 37.50 76.20 83.10 87.50 61.20 98.40 99.00 99.70 85.30 100.00 100.00 65.20 49.40 87.60 89.60 92.70 69.90 99.60 98.80 99.90 90.70 100.00 100.00 χ2 3 48.60 44.20 76.50 82.50 84.60 73.40 98.40 98.60 99.90 94.50 100.00 100.00 61.30 57.30 87.80 89.60 92.70 83.10 99.10 99.30 99.90 97.00 100.00 100.00 t5 15.50 44.50 49.10 55.00 24.50 74.00 87.30 89.20 39.30 97.50 99.90 99.90 25.00 59.50 63.00 67.90 35.40 84.70 93.70 94.90 51.10 99.50 100.00 100.00 Table 2: (Homoscedastic model, simple null hypothesis) Percentage of rejections for the Laplace null hypothesis at the significance levels 5% (upper entry) and 10% (lower entry). n= 25 n= 50 n= 100 A PB WB1 WB2 A PB WB1 WB2 A PB WB1 WB2 N3.70 22.60 17.20 19.10 4.20 42.60 38.10 39.30 8.30 69.70 68.60 69.20 8.40 30.90 25.00 27.20 9.10 51.80 48.10 50.20 14.60 77.60 78.10 78.20 LP 2.70 4.70 3.60 4.20 3.80 4.80 3.80 3.80 3.90 5.50 4.40 4.50 7.30 9.40 7.70 8.90 8.20 10.60 8.00 9.20 8.90 9.20 9.00 9.10 LG 4.20 25.60 18.90 20.60 4.70 40.60 36.90 37.50 5.90 69.90 70.00 70.70 7.50 35.00 28.50 31.20 9.30 48.80 46.60 47.50 11.90 78.30 78.30 79.00 G6.00 23.30 17.70 18.80 11.60 41.60 36.60 38.20 27.10 67.10 68.20 68.90 10.90 31.70 25.90 28.20 20.70 50.60 47.50 48.60 40.40 77.10 77.10 77.80 β35.50 12.80 13.40 15.30 78.60 19.30 30.90 32.60 99.30 36.20 66.00 66.60 48.80 19.20 21.80 24.30 86.20 27.70 43.20 44.50 99.60 46.60 74.00 75.60 χ2 3 17.50 20.00 16.00 17.80 44.40 34.00 32.60 33.80 92.30 61.50 65.60 66.50 27.20 27.60 24.00 25.70 59.10 44.00 43.70 44.60 96.90 72.00 76.20 76.70 t5 3.20 21.50 16.20 18.10 5.40 39.30 35.00 36.70 8.70 71.60 70.80 71.40 8.00 31.10 24.60 27.10 10.00 49.60 45.90 48.00 14.10 79.80 80.10 80.70 154 G.I. Rivas-Mart´ınez and M.D. Jim´enez-Gamero Table 3: (Homoscedastic model, composite null hypothesis) Percentage of rejections for the normality null hypothesis at significance levels 5% (upper entry) and 10% (lower entry). n= 25 n= 50 n= 100 PB WB1 WB2 PB WB1 WB2 PB WB1 WB2 N6.50 5.60 7.30 5.20 4.80 5.60 5.40 5.10 5.20 10.70 10.90 14.50 10.00 9.90 11.10 9.20 9.20 9.60 LP 29.90 15.30 21.50 33.60 40.40 43.50 38.30 80.50 81.60 40.50 26.20 30.60 44.10 56.30 58.70 54.70 90.40 91.00 LG 30.30 44.10 50.50 47.80 86.60 89.00 94.90 99.90 99.90 40.30 60.60 65.80 63.90 93.80 94.60 98.50 99.99 99.99 G29.10 18.50 21.50 35.70 42.50 43.50 51.80 80.50 83.60 43.50 29.20 30.60 51.30 58.30 59.70 66.10 90.40 95.90 β18.00 16.40 20.40 23.30 39.40 42.70 67.30 80.80 82.10 25.40 27.10 30.50 32.80 55.30 56.60 72.10 89.70 91.60 χ2 3 37.30 51.40 53.80 58.90 77.30 80.70 83.10 89.90 91.30 48.50 63.20 64.20 67.80 85.40 87.20 91.50 97.70 98.80 t5 40.40 14.50 21.50 52.90 38.90 42.40 76.80 82.30 83.10 58.70 28.70 31.40 69.20 53.70 56.00 88.20 89.50 90.30 Table 4: (Homoscedastic model, composite null hypothesis) Percentage of rejections for the Laplace null hypothesis at significance levels 5% (upper entry) and 10% (lower entry). n= 25 n= 50 n= 100 PB WB1 WB2 PB WB1 WB2 PB WB1 WB2 N53.20 56.90 58.80 62.80 64.40 66.20 69.30 71.40 77.20 66.30 68.20 71.10 74.50 75.40 76.60 80.60 80.90 81.20 LP 4.30 3.80 4.50 4.60 4.60 4.40 5.00 4.70 4.90 9.20 8.30 9.20 10.30 9.30 10.40 9.50 9.80 9.50 LG 52.40 48.20 50.50 60.40 58.50 60.20 74.60 77.50 78.50 65.30 62.00 65.70 72.10 71.70 73.90 90.80 93.20 93.70 G52.20 47.20 50.30 50.40 51.10 58.70 63.80 65.50 66.90 64.30 60.70 64.60 62.20 61.50 73.20 80.40 82.30 83.10 β50.50 57.00 62.90 55.80 60.60 65.60 76.40 83.50 87.70 63.60 71.50 76.50 72.30 74.20 77.30 87.70 95.60 98.80 χ2 3 37.50 67.30 70.60 41.40 78.50 80.10 43.50 88.00 88.40 51.50 79.60 82.30 54.10 91.30 93.20 59.60 97.30 98.30 t5 33.30 42.20 44.60 38.10 44.80 44.80 44.10 51.20 52.00 46.40 52.80 56.70 52.30 56.40 58.90 60.90 65.00 65.80 5.2. Heteroscedastic model The reported results correspond to the model Yj=Xj+X2 j+ (Xj+ 0.5)εj,1≤j≤n, Computationally Efficient Goodness-of-Fit Tests... 155 where Xjfollows the uniform (0,1) distribution. Since the model is heteroscedastic and the null hypothesis is simple, the simplifications in Remark 3.2 can be applied. Table 5 displays the results obtained for the type I error and the power for testing normality and Table 6 for testing GOF to the Laplace distribution. Similar conclusions to those given for Tables 1 and 2 can be also expressed in this case. Table 5: (Heteroscedastic model) Percentage of rejections for the normality null hypothesis at the significance levels 5% (upper entry) and 10% (lower entry). n= 25 n= 50 n= 100 A PB WB1 WB2 A PB WB1 WB2 A PB WB1 WB2 N4.50 6.00 5.40 6.50 4.90 5.36 4.82 5.92 4.90 5.32 5.08 5.74 10.30 10.80 10.20 12.50 10.50 10.70 10.12 11.74 9.40 10.30 10.64 11.24 LP 16.40 43.00 60.00 64.30 34.20 60.00 87.10 88.50 65.00 75.60 99.50 99.50 23.40 54.00 70.40 73.70 44.50 71.40 92.10 92.80 73.80 82.30 99.80 99.80 LG 7.40 57.60 56.40 63.10 8.50 91.90 94.70 95.20 12.60 99.80 100.00 100.00 12.10 70.40 72.00 76.00 15.00 95.90 97.40 98.10 20.20 99.90 100.00 100.00 G19.40 39.10 56.40 63.10 36.90 68.10 94.70 95.20 67.20 93.60 100.00 100.00 29.90 55.00 72.00 76.00 49.90 80.30 97.40 98.10 76.10 97.60 100.00 100.00 β43.00 16.10 57.60 63.30 86.20 77.00 99.80 99.80 99.90 95.20 100.00 100.00 56.20 26.10 70.00 74.70 92.30 86.20 100.00 100.00 100.00 97.60 100.00 100.00 χ2 3 50.90 41.60 85.50 89.10 83.00 71.30 99.70 99.70 99.20 95.70 100.00 100.00 61.80 54.80 92.50 93.90 91.00 83.10 99.90 99.90 99.70 98.80 100.00 100.00 t5 9.20 51.00 59.10 65.40 15.90 80.20 92.90 94.30 27.90 99.00 100.00 100.00 16.20 65.70 71.60 76.50 23.40 89.30 97.70 98.00 36.80 99.90 100.00 100.00 Table 6: (Heteroscedastic model) Percentage of rejections for the Laplace null hypothesis at the significance levels 5% (upper entry) and 10% (lower entry). n= 25 n= 50 n= 100 A PB WB1 WB2 A PB WB1 WB2 A PB WB1 WB2 N2.00 31.80 25.00 27.10 2.60 55.30 51.20 52.50 2.80 86.10 85.70 86.20 4.90 40.30 34.70 37.10 7.40 64.80 61.80 62.90 7.80 90.50 91.20 91.40 LP 2.10 4.60 3.70 4.60 3.00 5.70 4.00 4.40 3.60 4.40 4.00 4.40 6.80 10.00 8.00 9.60 7.30 11.50 9.20 10.20 7.80 9.10 8.40 9.00 LG 2.10 33.80 27.10 29.30 2.30 54.80 50.80 52.30 3.10 85.00 84.40 84.70 6.30 43.80 37.60 40.20 6.80 64.40 61.60 62.50 7.00 89.30 89.60 89.90 G2.10 31.30 23.50 25.50 2.80 53.90 50.20 51.50 3.00 85.30 85.10 85.60 6.70 41.10 34.40 37.10 6.80 65.10 62.70 63.70 7.50 91.10 91.10 91.50 β3.00 19.20 18.40 21.00 6.00 33.50 43.20 45.90 27.60 56.70 81.20 81.50 8.00 27.40 29.10 31.70 14.60 43.70 55.30 56.80 39.60 68.50 87.30 87.80 χ2 3 2.70 22.30 18.60 20.80 3.40 43.10 42.90 44.50 5.60 78.40 81.30 81.90 7.10 30.80 27.30 30.10 7.60 54.50 54.10 56.60 12.70 84.10 87.30 87.70 t5 2.90 30.60 22.80 24.50 3.90 56.80 53.20 53.90 4.60 84.30 83.90 84.30 6.30 41.50 33.70 38.00 6.50 66.70 64.30 65.20 9.40 90.20 90.40 90.70 156 G.I. Rivas-Mart´ınez and M.D. Jim´enez-Gamero 5.3. Time consumed Table 7 compares the PB and the WB (with raw and centered multipliers) in terms of the required CPU time. This table shows the CPU time consumed in seconds to get a p-value for testing GOF for the normal and the Laplace distributions in the homoscedastic (for both single and composite null hypothesis) and the heteroscedastic models with sample sizes n= 25,50,100,200. Looking at this table it becomes evident that the WB is more efficient than the PB, in terms of the required computing time, specially for larger sample sizes. The difference in time when using the raw and the centered multipliers is rather small. Table 7: CPU time consumed for the calculation of one p-value in seconds for testing normality and Laplace distribution for the homoscedastic model and composite null hypothesis (upper entry), the heteroscedastic model (middle entry) and the homoscedastic model and single null hypothesis (lower entry). nNormal distribution Laplace distribution PB/WB1 WB1 WB2 PB/WB1 WB1 WB2 2.72 0.71 0.74 3.49 1.00 1.01 25 7.45 0.33 0.35 7.17 0.54 0.60 4.42 0.31 0.34 5.34 0.50 0.55 5.61 0.71 0.70 7.51 1.08 1.09 50 30.88 0.17 0.22 38.15 0.26 0.25 15.63 0.19 0.19 23.68 0.28 0.25 12.15 0.84 0.86 23.40 1.11 1.12 100 52.80 0.25 0.27 74.33 0.42 0.45 30.64 0.25 0.26 64.56 0.37 0.39 27.56 1.25 1.27 76.37 1.54 1.58 200 66.19 0.59 0.62 127.80 0.83 0.83 41.14 0.56 0.58 117.51 0.78 0.76 The gain in computational efficiency of the WB over the PB stems from the fact that one does not have to re-estimate the parameters at each iteration, which slows down the process considerably. Note that in the WB the parameter θ, the regression function m(.) and the conditional variance function σ(·) are estimated only one time. For the WB approximation, once the set {mjk,1≤j≤k≤n}is computed, the WB replicates T∗1 2,n,ω(ˆ θ), ..., T∗B 2,n,ω(ˆ θ) can be calculated very rapidly. Computationally Efficient Goodness-of-Fit Tests... 157 6. CONCLUSIONS This paper proposes a WB approximation for the null distribution of a test statistic for testing GOF to the error distribution in nonparametric models. It provides a consistent estimator. The WB and the PB share this property. Nevertheless, from a computational point of view, the WB approximation is more efficient, in the sense of requiring less computation time. The numerical examples support these attributes. In addition, in cases were the asymptotic null distribution does not depend on unknown quantities, the simulations carried out declare that, for small to moderate sample sizes, the WB provides a better fit than the asymptotic distribution. To derive the results in this paper we considered certain estimators for the regression function and the conditional variance function. In addition, we assumed that the covariate was univariate. The results could be extended by considering other estimators (such as other local polynomial estimators) as well as covariates with higher dimension. The null distribution of other test statistics (for example, those based on the empirical CDF) could be similarly approximated. 7. APPENDIX 7.1. Assumptions (A.2) The weight function ωsatisfies (7.1) ω(t) = ω(−t),∀t, ω(t)≥0, ∀t, and Rt4ω(t)dt < ∞. There is no restriction in assuming that the weight function ω(t) satisfies (7.1) because otherwise by defining ω1(t) = 0.5{ω(t) + ω(−t)}, which satisfies (7.1), we have that Tn,ω(ˆ θ) = Tn,ω1(ˆ θ). (A.3)ε1, ..., εnare IID with E(ε4 j)<∞and ε1, ..., εnand X1, ..., Xnare independent. Recall that by construction we have that E(εj) = 0 and Var(εj) = 1. (A.4) (i)Xhas a compact support S. (ii)fX,mand σare twice continuously differentiable on S. (iii) infx∈SfX(x)>0 and infx∈Sσ(x)>0. 158 G.I. Rivas-Mart´ınez and M.D. Jim´enez-Gamero (A.5)nh4 n→0, nh2 n/ln n→ ∞. (A.6)Kis a twice continuously differentiable symmetric pdf with compact support. Assumptions (A.4)–(A.6) are mainly needed to guarantee the uniform consistency of the kernel estimators ˆ fX(·), ˆm(·) and ˆσ(·) for fX(·), m(·) and σ(·), respectively. (A.7) The first partial derivatives R′(t;θ), I′(t;θ), R(r)(t;θ), I(r)(t;θ), 1≤r≤p, exist and are continuous functions ∀t∈R,∀θin an open neighborhood of θ1. In addition, R′(t;θ), I′(t;θ), R(r)(t;θ), I(r)(t;θ), tR′(t;θ), tI′(t;θ), tR(r)(t;θ), tI(r)(t;θ), 1 ≤r≤p, are bounded by functions in L2(ω), ∀θin an open neighborhood of θ1. The following assumption will be used for the maximum likelihood estimator of the parameter. (A.9) The following functions exist ∀θin an open neighborhood of θ1: ur(x;θ) = ∂ ∂θrlog f(x;θ) , u1,r(x;θ) = ∂2 ∂x∂θrlog f(x;θ), u0,r,s(x;θ) = ∂2 ∂θr∂θslog f(x;θ) , u2,r(x;θ) = ∂3 ∂x2∂θrlog f(x;θ), u1,r,s(x;θ) = ∂3 ∂x∂θr∂θslog f(x;θ) , and satisfy |u1,r(a1+a2x;θ)| ≤ b1,r(x),with xb1,r(x), b1,r(x)∈L2(F) , |u0,r,s(a1+a2x;θ)| ≤ b0,r,s(x)∈L2(F) , |u2,r(a1+a2x;θ)| ≤ b2,r(x)∈L2(F) , |u1,r,s(a1+a2x;θ)| ≤ b1,r,s(x)∈L2(F) , ∀a1,a2,θ such that |a1|,|a2−1|,|θ−θ1|≤δ, for some small δ, 1≤r,s≤p. In addition, the following expectations exist: Eur(ε;θ1)us(ε;θ1), Eεu1,r(ε;θ1), 1≤r, s ≤p. The following assumption will be used for the method of moment estimator of the parameter, which assumes that under the null hypothesis, θ0=g(µ0), for some known function g= (g1, ..., gp)T,gr:Rk−1→R, 1 ≤r≤p: (A.10)gris twice continuously differentiable at a neighborhood of µF, 1≤r≤p. Computationally Efficient Goodness-of-Fit Tests... 159 7.2. Proofs We now sketch the proofs of the results stated in the previous sections, as well as some preliminary results. Along this section Mdenotes a generic positive constant taking many different values. Lemma 7.1. Suppose that assumptions (A.3)–(A.6) hold, then (a)1 nPn j=1(εj−ˆεj)2=op(1). (b)1 nPn j=1(ˆε2 j−ε2 j)2=op(1). (c)1 nPn j=1(ˆε2 j−1)2=Op(1). (d)1 nPn j=1 ˆε2 j=Op(1). Proof: First, observe that under the considered assumptions (see, for example, Masry [16]) sup x∈S|ˆm(x)−m(x)|=op(n−1/4),(7.2) sup x∈S|ˆσ(x)−σ(x)|=op(n−1/4).(7.3) The difference between the residuals and the errors can be written as follows (7.4) ˆεj−εj=εjσ(Xj)−ˆσ(Xj) ˆσ(Xj)+m(Xj)−ˆm(Xj) ˆσ(Xj). The results in (a)–(d) follow from (7.2)–(7.4). Lemma 7.2. If kˆ θ−θ1k=op(1) and (A.7) holds, then (a)kt{R′(t;ˆ θ)−R′(t;θ1)}k2 ω=op(1), kt{I′(t;ˆ θ)−I′(t;θ1)}k2 ω=op(1). (b)Rk∇R(t;ˆ θ)−∇R(t;θ1)k2ω(t)dt =op(1), Rk∇I(t;ˆ θ)−∇I(t;θ1)k2ω(t)dt =op(1). (c)kR(t;ˆ θ)−R(t;θ1)k2 ω=op(1), kI(t;ˆ θ)−I(t;θ1)k2 ω=op(1). (d)kt{R(t;ˆ θ)−R(t;θ1)}k2 ω=op(1), kt{I(t;ˆ θ)−I(t;θ1)}k2 ω=op(1). Proof: (a) From (A.7) tR′(t;θ)∈L2(ω), ∀θin a neighborhood of θ1. Since ˆ θP →θ1, the integral R{R′(t;ˆ θ)−R′(t;θ1)}2t2ω(t)dt is finite with probability tending to 1. Thus, ∀ǫ > 0,∃M=M(ǫ)>0 such that (7.5) ZR\[−M,M]{R′(t;ˆ θ)−R′(t;θ1)}2t2ω(t)dt < ǫ, 160 G.I. Rivas-Mart´ınez and M.D. Jim´enez-Gamero with probability tending to 1. tR′(t;θ) is a uniformly continuous function in [−M, M]×Bδ(θ1) = C, where Bδ(θ1) = {θ:kθ−θ1k ≤ δ}. Thus, ∀ǫ > 0, ∃ρ= ρ(ǫ)>0 such that ∀(ta, θa),(tb, θb)∈Csatisfying k(ta, θa)−(tb, θb)k< ρ, we have |t1R′(ta;θa)−t2R′(tb;θb)|< ǫ/ι, with ι=Rω(t)dt. As a consequence (7.6) ZM −M{R′(t;ˆ θ)−R′(t;θ1)}2t2ω(t)dt < ǫ, with probability tending to 1. As ǫis arbitrary, the result in (a) for the real part follows from (7.5) and (7.6). The proof for the imaginary part is parallel. (b) The proof of this part is quite similar to that of part (a). Parts (c) and (d) can be proven by applying the mean value theorem. Proof of Theorem 2.1: W∗can be expressed as W∗=W1+W2+ 2W3, where W2 3≤W1W2,W1=k1 √nPn j=1 Z0(εj;t, θ1)ξjk2 ω,W2=k1 √nPn j=1{Z0(ˆεj;t, ˆ θ) −Z0(εj;t, θ1)}ξjk2 ω. From the results in [4], sup x|P∗{W1≤x}−P{W0≤x}| a.s. −→ 0. Thus, to show the result it suffices to see that W2=op∗(1) in probability. With this aim, observe that W2can be expressed as W2=P4 j=1 Sj+Pj6=kSjk, with S2 jk ≤SjSk, 1 ≤j, k ≤4. In the proof of Theorem 3.1 it is given the expression of Sjand it is also proven that Sj=op∗(1) in probability, 1 ≤j≤4. This proves the result. Proof of Theorem 3.1: T∗ 2,n,ω(ˆ θ) can be expressed as T∗ 2,n,ω(ˆ θ) = D1+ D2+ 2D3, where D1=k1 √nPn j=1 Z2(εj;t, θ1)ξjk2 ω,D2=k1 √nPn j=1{Z2(ˆεj;t, ˆ θ)− Z2(εj;t, θ1)}ξjk2 ω,D2 3≤D1D2. From the results in [4], sup x|P∗{D1≤x}−P{T2≤x}| a.s. −→ 0. Thus, to show the result it suffices to see that D2=op∗(1) in probability. With this aim, observe that D2can be expressed as D2= 10 X j=1 Sj+X k<j Sjk, with S2 jk ≤SjSk, 1 ≤j, k ≤10, S1=k1 √nPn j=1{cos(tεj)−cos(tˆεj)}ξjk2 ω, S2=k1 √nPn j=1{sin(tεj)−sin(tˆεj)}ξjk2 ω, S3=k1 √n{R(t;ˆ θ)−R(t;θ1)}Pn j=1 ξjk2 ω, Computationally Efficient Goodness-of-Fit Tests... 161 S4=k1 √n{I(t;ˆ θ)−I(t;θ1)}Pn j=1 ξjk2 ω, S5=kt √nPn j=1{ˆεjR(t;ˆ θ)−εjR(t;θ1)}ξjk2 ω, S6=kt √nPn j=1{ˆεjI(t;ˆ θ)−εjI(t;θ1)}ξjk2 ω, S7=kt 2√nPn j=1{(ˆε2 j−1)R′(t;ˆ θ)−(ε2 j−1)R′(t;θ1)}ξjk2 ω, S8=kt 2√nPn j=1{(ˆε2 j−1)I′(t;ˆ θ)−(ε2 j−1)I′(t;θ1)}ξjk2 ω, S9=k1 √nPn j=1{ψT n(ˆεj;ˆ θ)∇R(t;ˆ θ)−ψT 1(εj;θ)∇R(t;θ1)}ξjk2 ω, S10 =k1 √nPn j=1{ψT n(ˆεj;ˆ θ)∇I(t;ˆ θ)−ψT 1(εj;θ)∇I(t;θ1)}ξjk2 ω. We will show that Sj=op∗(1) in probability, 1 ≤j≤10. By the mean value theorem, S1=1 n n X j,k=1 ξjξk(εj−ˆεj)(εk−ˆεk)Zt2sin(t∼ εj) sin(t∼ εk)ω(t)dt, where ∼ εj=αjεj+ (1 −αj)ˆεj, for some αj∈(0,1). Then, from Lemma 7.1 (a), E∗(S1)≤1 n n X j=1 (εj−ˆεj)2Zt2ω(t)dt =op(1), which implies S1=op∗(1) in probability. Analogously, S2=op∗(1) in probability. Since S3=1 √nPn j=1 ξj2kR(t;ˆ θ)−R(t;θ1)k2 ω, the central limit theorem and Lemma 7.2 (c) imply that S3=op∗(1) in probability. Analogously, S4= op∗(1) in probability. Observe that S5=S51 +S52 + 2S53, with S2 53 ≤S51S52, S51 =1 nPn j,k=1(ˆεj−εj)(ˆεk−εk)ξjξkktR(t;ˆ θ)k2 ω, S52 =1 nPn j,k=1 εjεkξjξkkt{R(t;ˆ θ)−R(t;θ1)}k2 ω. From Lemma 7.1 (a) and Assumption (A.2), it follows that E∗(S51) = op(1) and thus S51 =op∗(1), in probability. From Lemma 7.2 (d), it follows that E∗(S52) = op(1) and thus S52 =op∗(1), in probability. Therefore, S5=op∗(1), in probability. Analogously, S6=op∗(1), in probability. Observe that S7=S71 +S72 + 2S73, with S2 73 ≤S71S72, S71 =1 4 1 nPn j,k=1(ˆε2 j−1)(ˆε2 k−1)ξjξkkt{R′(t;ˆ θ)−R′(t;θ1)}k2 ω, S72 =1 4 1 nPn j,k=1(ˆε2 j−ε2 j)(ˆε2 k−ε2 k)ξjξkktR′(t;θ1)k2 ω. From Lemma 7.1 (c) and Lemma 7.2 (a), it follows that E∗(S71) = op(1) and thus S71 =op∗(1), in probability.