# Remodeling and Estimation for Sparse Partially Linear Regression Models

Authors: Yunhui Zeng; Xiuli Wang; Lu Lin
Journal: Abstract and Applied Analysis
 (2013)

---

## Abstract

When the dimension of covariates in the regression model is high, one usually uses a submodel as a working model that contains significant variables. But it may be highly biased and the resulting estimator of the parameter of interest may be very poor when the coefficients of removed variables are not exactly zero. In this paper, based on the selected submodel, we introduce a two-stage remodeling method to get the consistent estimator for the parameter of interest. More precisely, in the first stage, by a multistep adjustment, we reconstruct an unbiased model based on the correlation information between the covariates; in the second stage, we further reduce the adjusted model by a semiparametric variable selection method and get a new estimator of the parameter of interest simultaneously. Its convergence rate and asymptotic normality are also obtained. The simulation results further illustrate that the new estimator outperforms those obtained by the submodel and the full model in the sense of mean square errors of point estimation and mean square prediction errors of model prediction.


---

## Body

## 1. Introduction

Consider the following partially linear regression model:(1)Y=βTX+γTZ+f(U)+ε,
where Y is a scalar response variable, X and Z are, respectively, p-dimensional and q-dimensional continuous-valued covariates with p being finite and p≪q, β is the parameter vector of interest and γ is the nuisance parameter vector which is supposed to be sparse in the sense that ∥γ∥2 is small, f(·) is an unknown function satisfying Ef(U)=0 for identification, ε is the random error satisfying E(ε∣X,Z,U)=0. For simplicity, we assume that U is univariate. Let (Yi,Xi,Zi,Ui), i=1,…,n, be i.i.d. observations of (Y,X,Z,U) obtained from the above model.A feature of the model is that the parametric part contains both the parameter vector of interest and nuisance parameter vector. The reason for this coefficient separation is as follows. In practice we often use such a model to distinguish the main treatment variables of interest from the state variables. For instance, in a clinical trial,X consists of treatment variables and can be easily controlled, Z is a vector of many clinical variables, such as patient ages and body weights. The variables in Z may have an impact on Y but are not of primary interest and the effects may be small. In order to make up for potentially nonnegligible effects on the response Y, the nuisance covariate Z are introduced into model (1); see Shen et al. [1]. Model (1) contains all relevant covariates and in this paper we call it full model.The purpose of this paper is to estimateβ, the parameter of interest, when γTZ is removed from the model. The main idea is remodeling based on the following working model:
(2)Y=βTX+f(U)+η.
As is known, E(η∣X=x,U=u) is a nonzero function if γTE(Z∣X,U)≠0, which relies on two elements, one is E(Z∣X,U), related with the correlation between the covariates of Z and (X,U), the other is γ, determined by the nuisance parameter in the removed part. Thus the least squares estimator based on model (2) may be inconsistent. In the following, we will make use of the above two elements. Specifically, in the first stage, we shall construct a remodeled model by a multistep-adjustment to correct the submodel bias based on the correlation information between the covariates. This adjustment is motivated by Gai et al. [2]. In the paper, they proposed a nonparametric adjustment by adding a univariate nonparametric estimation to the working model (2), and it can dramatically reduce the bias of the working model. But this only holds in a subset of the covariates, although the subset may be fairly large. In order to obtain a globally unbiased working model for linear regression model, Zeng et al. [3] adjusted the working model by multiple steps. Because only those variables in Z correlated with (X,U) may have impact on estimation of β, in each step a univariate nonparametric part was added to the working model and consequently a globally unbiased working model was obtained.However, when many components ofZ are correlated with (X,U), the number of nonparametric functions added in the above working model is large. Such a model is improper in practice. Thus, in the second stage, we further simplify the above adjusted model by a semiparametric variable selection procedure proposed by Zhao and Xue [4]. Their method can select significant parametric and nonparametric components simultaneously under sparsity condition for semiparametric varying coefficient partially linear models. The relevant papers include Fan and Li [5], Wang et al. [6, 7], among others. After two-stage remodeling, the final model is conditionally unbiased. Based on this model, the estimation and model prediction are significantly improved.The rest of this paper is organized as follows. In Section2, a multistep adjustment and remodeled models are firstly proposed, then the models are further simplified via the semiparametric SCAD variable selection procedure. A new estimator of the parameter of interest based on the simplified model is derived, its convergence rate and asymptotic normality are also obtained. Simulations are given in Section 3. A short conclusion and some remarks are contained in Section 4. Some regular conditions and theoretical proofs are presented in the appendix.


## 2. New Estimator for the Parameter of Interest

In this paper, we suppose that covariateZ has zero mean, p is finite and p≪q, E(ε∣X,Z,U)=0 and Var(ε∣X,Z,U)=σ2. We also assume that covariates X and U and parameter β are prespecified, so that the submodel (2) is a fixed model.


### 2.1. Multistep-Adjustment by Correlation

In this subsection, we first adjust the submodel to be conditionally unbiased by a multistep-adjustment.WhenZ is normally distributed, the principal component analysis (PCA) method will be used. Let ΣZ be the covariance matrix of Z, then there exists an orthogonal q×q matrix Q such that QΣZQT=Λ, where Λ is the diagonal matrix diag(λ1,λ2,…,λq) with λ1≥λ2≥⋯≥λq≥0 being eigenvalues of ΣZ. Denote QT=(τ1,τ2,…,τq) and Z~(j)=τjTZ.WhenZ is centered but nonnormally distributed, we shall apply independent component analysis (ICA) method. Assume that Z is generated by a nonlinear combination of independent components Z~(i), that is Z=F(Z~), where F(·) is an unknown nonlinear mapping from Rq to Rq, Z~ is an unknown random vector with independent components. By imposing some constraints on the nonlinear mixing mapping F or the independent components Z~(i), the independent components Z~(i) can be properly estimated. See Simas Filho and Seixas [8] for an overview of the main statistical principles and some algorithms for estimating the independent components. For simplicity, in this paper we suppose that Z=(Z(1),…,Z(q))T with Z(l)=Σj=1qFlj(Z~(j)), l=1,…,q, and Flj(·) are scalar functions.In the above two cases,Z~(j)'s are independent of each other. Set K0 to be the size of set M0={j:E(Z~(j)∣X,U)≠0,1≤j≤q}. Without loss of generality, let M0={1,…,K0}.We construct the following adjusted model:(3)Y=βTX+∑j=1K0‍gj(Z~(j))+f(U)+ζK0,
where gj(Z~(j))=E(Y-βTX-f(U)∣Z~(j))=γTE(Z∣Z~(j)), j=1,…,K0 and ζK0=Y-βTX-g1(Z~(1))-⋯-gK0(Z~(K0))-f(U). The model (3) is based on Z's population and depends on the distributions of X, U and Z. It is easy to see that model (3) is conditionally unbiased, that is, E(ζK0∣X,U,Z~(j),1≤j≤K0)=0.The adjusted model (3) is an additive partially linear model, in which βTX is the parametric part, f(U) and gj(Z~(j)), j=1,…,K0, are the nonparametric parts and ζK0 is the random error. Compared with the submodel (2), the nonparametric parts gj(Z~(j)), j=1,…,K0, may be regarded as bias-corrected terms for the random error η. For centered Z, E(gj(Z~(j)))=0, j=1,…,K0, the nonparametric components g1(Z~(1)),…,gK0(Z~(K0)) can be properly identified. In fact, centered Z can be relaxed to any Z such that satisfies γTE(Z)=0.WhenZ is centered and normally distributed, the nonparametric parts gj(Z~(j))=αjτjTZ=αjZ~(j), j=1,…,K0. So the multistep adjusted model (3) is really a partially linear model
(4)Y=βTX+αTZ~K0+f(U)+ζK0
with α=(α1,…,αK0)T and Z~K0=(Z~(1),…,Z~(K0))T. Specially, when f(U)≡0, the full model is a linear model, the multistep adjusted model is also a linear model
(5)Y=βTX+αTZ~K0+ζK0.
But when the variables in Z are not jointly normal, the nonparametric parts gj can be highly nonlinear, which are similar to the results of marginal regression; see Fan et al. [9].


### 2.2. Model Simplification

When the most of the features in the full model are correlated, thenK0 is very large and even is close to q. In this case, the adjusted model (3) is improper in practice, so we shall use the group SCAD regression procedure, proposed by Wang et al. [6], and the semiparametric variable selection procedure, proposed by Zhao and Xue [4], to further simplify the model.Lets=|ℳ*| with ℳ*={1≤j≤K0:E(gj(Z~(j)))2>0}, and assume that the model (3) is sparse, that is, s is small. We define the semiparametric penalized least squares as
(6)F(β,g(·),f(·))=∑i=1n‍{Yi-βTXi-∑j=1K0‍gj(Z~i(j))-f(Ui)}2+n∑j=1K0‍pλj(∥gj(Z~(j))∥),
where ∥gj(Z~(j))∥=(E(gj(Z~(j)))2)1/2, and pλ(·) is the SCAD penalty function with λ being a tuning parameter defined as
(7)pλ′(w)=λ{I(w≤λ)+(aλ-w)+(a-1)λI(w>λ)},
with a>2, w>0 and pλ(0)=0. In (6), g(·) denotes the set {gj(Z~(j)),j=1,…,K0}. Because gj are nonparametric functions, thus they cannot be directly applied for minimization. Here we will replace f(·) and g(·) by basis function approximations. For 1≤j≤K0, let {Ψjk,k=1,…,L} be orthogonal basis functions satisfying
(8)E(ΨjkΨjl)≡∫supp‍Ψjk(Z~(j))Ψjl(Z~(j))rj(Z~(j))dZ=δkl={0,k≠l;1,k=l,
where rj(Z~(j)) is the density function of Z~(j). Similarly, let {Ψ0k,k=1,…,L} be orthogonal basis functions satisfying the above condition which is only replaced by the support and density function of U. Denote Ψj(Z~(j))=(Ψj1(Z~(j)),…,ΨjL(Z~(j)))T, Ψ0(U)=(Ψ01(U),…,Ψ0L(U))T. Then gj(Z~(j)) and f(U) can be approximated by
(9)gj(Z~(j))≈θjTΨj(Z~(j)),f(U)≈νTΨ0(U).
Denote ∥θj∥2=(θjTθj)1/2, invoking that E(Ψj(Z~(j))·ΨjT(Z~(j)))=IL the identity matrix, we get
(10)F(β,θ,ν)=∑i=1n‍{Yi-βTXi-θTΨi-νTΨ0i}2+n∑j=1K0‍pλj(∥θj∥2),
where θ=(θ1T,…,θK0T)T, Ψi≡Ψ(Z~i)=Vec(Ψ1(Z~i(1)),…,ΨK0(Z~i(K0))), Ψ0i≡Ψ0(Ui).Denote byβ^, θ^=(θ^1T,…,θ^K0T)T and ν^ the least squares estimators based on the penalized function (10), that is (β^,θ^,ν^)=argminβ∈Rp,θj∈RL,ν∈RLF(β,θ,ν). Let g^j≡g^j(Z~(j))=θ^jTΨj(Z~(j)) and f^≡f^(U)=ν^TΨ0(U), then g^j is an estimator of gj(Z~(j)), f^ is an estimator of f(U).Letℳ^n={1≤j≤K0:θ^j≠0} and Kn=|ℳ^n|. For simplicity, we assume that ℳ*={1,2,…,s} and ℳ^n={1,2,…,Kn}. So we get the following simplified working model
(11)Y=βTX+∑j=1Kn‍gj(Z~(j))+f(U)+ζKn,
where gj(Z~(j))=E(γTZ∣Z~(j)), j=1,…,Kn and ζKn=Y-βTX-g1(Z~(1))-⋯-gKn(Z~(Kn))-f(U). Under the assumption of sparsity, the model (11) contains all of significant nonparametric functions and fully utilizes both the correlation information of covariates and the model sparsity on nuisance covariate.IfZ is centered and normally distributed with covariance matrix ΣZ=Iq the identity matrix, then τj=ej, j=1,…,q, where ej denotes the unit vector with 1 at position j, and α is sparse with αj=γTτj=γj. So the model (4) is sparse. For model (5), the special case of model (4), we can apply the SCAD penalty method proposed by Fan and Li [5] to select variables in Z~K0 and estimate parameters α and β simultaneously. The selected covariate and the corresponding parameter are denoted by Z~Kn and αKn, the resulting parameter estimators are denoted by α^Kn and β^, respectively. Finally, we can use the simplified model
(12)Y=βTX+αKnTZ~Kn+ζKn
for model prediction. Under the condition of sparsity, its model size is much smaller than those of the multistep adjusted model (5) and the full model (1).


### 2.3. Asymptotic Property of Point Estimator

Letβ0, θ0, ν0, and gj0(·), f0(·) be the true values of β, θ, ν, and gj(·), f(·), respectively, in model (3). Without loss of generality, we assume that gj0(Z~(j))=0, j=s+1,…,K0, and gj0(Z~(j)), j=1,…,s, are all nonzero components.We suppose thatgj(Z~(j)), j=1,…,K0 can be expressed as ∑k=1∞‍θjkΨjk(Z~(j)) and f(U) can be expressed as ∑k=1∞‍νkΨ0k(U), θj and ν belong to the Sobolev ellipsoid S(r,M)={θ:∑k=1∞‍θk2k2r≤M,M>0,r>0}.The following theorem gives the consistency of the penalized SCAD estimators.Theorem 1.
Suppose that the regularity conditions (C1)–(C5) in the appendix hold and the number of termsL=Op(n1/(2r+1)). Then,(i)
∥β^-β0∥=Op(n-r/(2r+1)+an),(ii)
∥g^j(·)-gj0(·)∥=Op(n-r/(2r+1)+an),j=1,…,K0,(iii)
∥f^(·)-f0(·)∥=Op(n-r/(2r+1)+an),
wherean=maxj{|pλj′(∥θj0∥2)|:θj0≠0}.From the last paragraph of Section 2.2 we know that, for linear regression model and normally distributedZ, the multistep adjusted model (5) is a linear model. By orthogonal basis functions, such as power series, we have r=∞, then ∥β^-β0∥=Op(n-1/2+an), implying the estimator β^ has the same convergence rate as that of the SCAD estimator in Fan and Li [5].Theorem 2.
Suppose that the regularity conditions (C1)–(C6) in the appendix hold and the number of termsL=Op(n1/(2r+1)). Let λmax=maxj{λj} and λmin=minj{λj}. If λmax→0 and nr/(2r+1)λmin→∞ as n→∞, then, with probability tending to 1, g^j(·)≡0, j=s+1,…,K0.Remark 3.
By Remark  1 of Fan and Li [5], we have that, if λmax→0 as n→∞, then an→0. Hence from Theorems 1 and 2, by choosing proper tuning parameters, the variable selection method is consistent and the estimators of nonparametric components achieve the optimal convergence rate as if the subset of true zero coefficients was already known; see Stone [10].Letθ*=(θ1T,…,θsT)T be the nonzero components of θ, corresponding covariates are denoted by Ψi*,i=1,…,n. In addition, let
(13)Σ=1σK02{E(Ψ˘0XT)E(XXT)-E(XΨ*T)E-1(Ψ*Ψ*T)E(Ψ*XT)-E(XΨ˘0T)E-1(Ψ˘0Ψ˘0T)E(Ψ˘0XT)},
where σK02=Var(ζK0i) for homoscedastic case, Ψ˘0=Ψ0-E(Ψ0Ψ*T)E-1(Ψ*Ψ*T)Ψ*.Theorem 4.
Suppose that the regularity conditions (C1)–(C6) in the appendix hold and the number of termsL=Op(n1/(2r+1)). If Σ is invertible, then
(14)n(β^-β0)→𝒟N(0,Σ-1),
where “→𝒟” denotes the convergence in distribution.Remark 5.
From Theorems1 and 4, it can be found that the penalized estimators have the oracle property. Furthermore, the estimator of the parameter of interest has the same asymptotic distribution as that based on the correct submodel.


### 2.4. Some Issues on Implementation

In the adjusted model (4), τj,j=1,…,K0 are used. When the population distribution is not available, they need to be approximated by estimators. When Z is normally distributed and eigenvalues rj,j=1,…,q of the covariance matrix ΣZ are different from each other, then n(uj-τj) is asymptotically N(0,Vj) with Vj=∑l≠jq‍(rjrl/(rj-rl)2)τlτlT, where uj is the jth eigenvector of S=(1/(n-1))∑i=1n‍(Zi-Z-)(Zi-Z-)T with Z-=(1/n)∑i=1n‍Zi; see Anderson [11]. For the case when the population size is large and comparable with the sample size, if the covariance matrix is sparse, we can use the method in Rütimann and Bühlmann [12] or Cai and Liu [13] to estimate the covariance matrix. So we can use uj to approximate τj. When τj in model (4) are replaced by these consistent estimators, one can see that the approximation error can be neglected without changing the asymptotic property.The nonparametric partsgl(Z~(l)) in the adjusted model depend on the univariate variable Z~(l), for l=1,…,K0. So it needs to choose the steps K0 firstly. In real implementation, we compute all the q multiple correlation coefficients of Z~(l) (l=1,…,q) with X and U. Then we choose the components R={Z~(l):|mcorr(Z~(l),(X,U))|≥δ,l=1,…,q} for given small number δ>0, where mcorr(u,V) denotes the multicorrelation coefficient between u and V and can be approximated by its sample form; see Anderson [11].There are some tuning parameters needing to choose in order to implement the two-stage remodeling procedure. Fan and Li [5] showed that the SCAD penalty with a=3.7 performs well in a variety of situations. Hence, we use their suggestion throughout this paper. We still need to choose the positive integer L for basis functions and the tuning parameter λj of the penalty functions. Similar to the adaptive lasso of Zou [14], we suggest taking λj=λ/∥θ^j(0)∥2, where θ^j(0) is initial estimator of θj by using ordinary least squares method based on the first term in (10). So the two remaining parameters L and λ can be selected simultaneously using the leave-one-out CV or GCV method; see Zhao and Xue [4] for more details.


## 2.1. Multistep-Adjustment by Correlation

In this subsection, we first adjust the submodel to be conditionally unbiased by a multistep-adjustment.WhenZ is normally distributed, the principal component analysis (PCA) method will be used. Let ΣZ be the covariance matrix of Z, then there exists an orthogonal q×q matrix Q such that QΣZQT=Λ, where Λ is the diagonal matrix diag(λ1,λ2,…,λq) with λ1≥λ2≥⋯≥λq≥0 being eigenvalues of ΣZ. Denote QT=(τ1,τ2,…,τq) and Z~(j)=τjTZ.WhenZ is centered but nonnormally distributed, we shall apply independent component analysis (ICA) method. Assume that Z is generated by a nonlinear combination of independent components Z~(i), that is Z=F(Z~), where F(·) is an unknown nonlinear mapping from Rq to Rq, Z~ is an unknown random vector with independent components. By imposing some constraints on the nonlinear mixing mapping F or the independent components Z~(i), the independent components Z~(i) can be properly estimated. See Simas Filho and Seixas [8] for an overview of the main statistical principles and some algorithms for estimating the independent components. For simplicity, in this paper we suppose that Z=(Z(1),…,Z(q))T with Z(l)=Σj=1qFlj(Z~(j)), l=1,…,q, and Flj(·) are scalar functions.In the above two cases,Z~(j)'s are independent of each other. Set K0 to be the size of set M0={j:E(Z~(j)∣X,U)≠0,1≤j≤q}. Without loss of generality, let M0={1,…,K0}.We construct the following adjusted model:(3)Y=βTX+∑j=1K0‍gj(Z~(j))+f(U)+ζK0,
where gj(Z~(j))=E(Y-βTX-f(U)∣Z~(j))=γTE(Z∣Z~(j)), j=1,…,K0 and ζK0=Y-βTX-g1(Z~(1))-⋯-gK0(Z~(K0))-f(U). The model (3) is based on Z's population and depends on the distributions of X, U and Z. It is easy to see that model (3) is conditionally unbiased, that is, E(ζK0∣X,U,Z~(j),1≤j≤K0)=0.The adjusted model (3) is an additive partially linear model, in which βTX is the parametric part, f(U) and gj(Z~(j)), j=1,…,K0, are the nonparametric parts and ζK0 is the random error. Compared with the submodel (2), the nonparametric parts gj(Z~(j)), j=1,…,K0, may be regarded as bias-corrected terms for the random error η. For centered Z, E(gj(Z~(j)))=0, j=1,…,K0, the nonparametric components g1(Z~(1)),…,gK0(Z~(K0)) can be properly identified. In fact, centered Z can be relaxed to any Z such that satisfies γTE(Z)=0.WhenZ is centered and normally distributed, the nonparametric parts gj(Z~(j))=αjτjTZ=αjZ~(j), j=1,…,K0. So the multistep adjusted model (3) is really a partially linear model
(4)Y=βTX+αTZ~K0+f(U)+ζK0
with α=(α1,…,αK0)T and Z~K0=(Z~(1),…,Z~(K0))T. Specially, when f(U)≡0, the full model is a linear model, the multistep adjusted model is also a linear model
(5)Y=βTX+αTZ~K0+ζK0.
But when the variables in Z are not jointly normal, the nonparametric parts gj can be highly nonlinear, which are similar to the results of marginal regression; see Fan et al. [9].


## 2.2. Model Simplification

When the most of the features in the full model are correlated, thenK0 is very large and even is close to q. In this case, the adjusted model (3) is improper in practice, so we shall use the group SCAD regression procedure, proposed by Wang et al. [6], and the semiparametric variable selection procedure, proposed by Zhao and Xue [4], to further simplify the model.Lets=|ℳ*| with ℳ*={1≤j≤K0:E(gj(Z~(j)))2>0}, and assume that the model (3) is sparse, that is, s is small. We define the semiparametric penalized least squares as
(6)F(β,g(·),f(·))=∑i=1n‍{Yi-βTXi-∑j=1K0‍gj(Z~i(j))-f(Ui)}2+n∑j=1K0‍pλj(∥gj(Z~(j))∥),
where ∥gj(Z~(j))∥=(E(gj(Z~(j)))2)1/2, and pλ(·) is the SCAD penalty function with λ being a tuning parameter defined as
(7)pλ′(w)=λ{I(w≤λ)+(aλ-w)+(a-1)λI(w>λ)},
with a>2, w>0 and pλ(0)=0. In (6), g(·) denotes the set {gj(Z~(j)),j=1,…,K0}. Because gj are nonparametric functions, thus they cannot be directly applied for minimization. Here we will replace f(·) and g(·) by basis function approximations. For 1≤j≤K0, let {Ψjk,k=1,…,L} be orthogonal basis functions satisfying
(8)E(ΨjkΨjl)≡∫supp‍Ψjk(Z~(j))Ψjl(Z~(j))rj(Z~(j))dZ=δkl={0,k≠l;1,k=l,
where rj(Z~(j)) is the density function of Z~(j). Similarly, let {Ψ0k,k=1,…,L} be orthogonal basis functions satisfying the above condition which is only replaced by the support and density function of U. Denote Ψj(Z~(j))=(Ψj1(Z~(j)),…,ΨjL(Z~(j)))T, Ψ0(U)=(Ψ01(U),…,Ψ0L(U))T. Then gj(Z~(j)) and f(U) can be approximated by
(9)gj(Z~(j))≈θjTΨj(Z~(j)),f(U)≈νTΨ0(U).
Denote ∥θj∥2=(θjTθj)1/2, invoking that E(Ψj(Z~(j))·ΨjT(Z~(j)))=IL the identity matrix, we get
(10)F(β,θ,ν)=∑i=1n‍{Yi-βTXi-θTΨi-νTΨ0i}2+n∑j=1K0‍pλj(∥θj∥2),
where θ=(θ1T,…,θK0T)T, Ψi≡Ψ(Z~i)=Vec(Ψ1(Z~i(1)),…,ΨK0(Z~i(K0))), Ψ0i≡Ψ0(Ui).Denote byβ^, θ^=(θ^1T,…,θ^K0T)T and ν^ the least squares estimators based on the penalized function (10), that is (β^,θ^,ν^)=argminβ∈Rp,θj∈RL,ν∈RLF(β,θ,ν). Let g^j≡g^j(Z~(j))=θ^jTΨj(Z~(j)) and f^≡f^(U)=ν^TΨ0(U), then g^j is an estimator of gj(Z~(j)), f^ is an estimator of f(U).Letℳ^n={1≤j≤K0:θ^j≠0} and Kn=|ℳ^n|. For simplicity, we assume that ℳ*={1,2,…,s} and ℳ^n={1,2,…,Kn}. So we get the following simplified working model
(11)Y=βTX+∑j=1Kn‍gj(Z~(j))+f(U)+ζKn,
where gj(Z~(j))=E(γTZ∣Z~(j)), j=1,…,Kn and ζKn=Y-βTX-g1(Z~(1))-⋯-gKn(Z~(Kn))-f(U). Under the assumption of sparsity, the model (11) contains all of significant nonparametric functions and fully utilizes both the correlation information of covariates and the model sparsity on nuisance covariate.IfZ is centered and normally distributed with covariance matrix ΣZ=Iq the identity matrix, then τj=ej, j=1,…,q, where ej denotes the unit vector with 1 at position j, and α is sparse with αj=γTτj=γj. So the model (4) is sparse. For model (5), the special case of model (4), we can apply the SCAD penalty method proposed by Fan and Li [5] to select variables in Z~K0 and estimate parameters α and β simultaneously. The selected covariate and the corresponding parameter are denoted by Z~Kn and αKn, the resulting parameter estimators are denoted by α^Kn and β^, respectively. Finally, we can use the simplified model
(12)Y=βTX+αKnTZ~Kn+ζKn
for model prediction. Under the condition of sparsity, its model size is much smaller than those of the multistep adjusted model (5) and the full model (1).


## 2.3. Asymptotic Property of Point Estimator

Letβ0, θ0, ν0, and gj0(·), f0(·) be the true values of β, θ, ν, and gj(·), f(·), respectively, in model (3). Without loss of generality, we assume that gj0(Z~(j))=0, j=s+1,…,K0, and gj0(Z~(j)), j=1,…,s, are all nonzero components.We suppose thatgj(Z~(j)), j=1,…,K0 can be expressed as ∑k=1∞‍θjkΨjk(Z~(j)) and f(U) can be expressed as ∑k=1∞‍νkΨ0k(U), θj and ν belong to the Sobolev ellipsoid S(r,M)={θ:∑k=1∞‍θk2k2r≤M,M>0,r>0}.The following theorem gives the consistency of the penalized SCAD estimators.Theorem 1.
Suppose that the regularity conditions (C1)–(C5) in the appendix hold and the number of termsL=Op(n1/(2r+1)). Then,(i)
∥β^-β0∥=Op(n-r/(2r+1)+an),(ii)
∥g^j(·)-gj0(·)∥=Op(n-r/(2r+1)+an),j=1,…,K0,(iii)
∥f^(·)-f0(·)∥=Op(n-r/(2r+1)+an),
wherean=maxj{|pλj′(∥θj0∥2)|:θj0≠0}.From the last paragraph of Section 2.2 we know that, for linear regression model and normally distributedZ, the multistep adjusted model (5) is a linear model. By orthogonal basis functions, such as power series, we have r=∞, then ∥β^-β0∥=Op(n-1/2+an), implying the estimator β^ has the same convergence rate as that of the SCAD estimator in Fan and Li [5].Theorem 2.
Suppose that the regularity conditions (C1)–(C6) in the appendix hold and the number of termsL=Op(n1/(2r+1)). Let λmax=maxj{λj} and λmin=minj{λj}. If λmax→0 and nr/(2r+1)λmin→∞ as n→∞, then, with probability tending to 1, g^j(·)≡0, j=s+1,…,K0.Remark 3.
By Remark  1 of Fan and Li [5], we have that, if λmax→0 as n→∞, then an→0. Hence from Theorems 1 and 2, by choosing proper tuning parameters, the variable selection method is consistent and the estimators of nonparametric components achieve the optimal convergence rate as if the subset of true zero coefficients was already known; see Stone [10].Letθ*=(θ1T,…,θsT)T be the nonzero components of θ, corresponding covariates are denoted by Ψi*,i=1,…,n. In addition, let
(13)Σ=1σK02{E(Ψ˘0XT)E(XXT)-E(XΨ*T)E-1(Ψ*Ψ*T)E(Ψ*XT)-E(XΨ˘0T)E-1(Ψ˘0Ψ˘0T)E(Ψ˘0XT)},
where σK02=Var(ζK0i) for homoscedastic case, Ψ˘0=Ψ0-E(Ψ0Ψ*T)E-1(Ψ*Ψ*T)Ψ*.Theorem 4.
Suppose that the regularity conditions (C1)–(C6) in the appendix hold and the number of termsL=Op(n1/(2r+1)). If Σ is invertible, then
(14)n(β^-β0)→𝒟N(0,Σ-1),
where “→𝒟” denotes the convergence in distribution.Remark 5.
From Theorems1 and 4, it can be found that the penalized estimators have the oracle property. Furthermore, the estimator of the parameter of interest has the same asymptotic distribution as that based on the correct submodel.


## 2.4. Some Issues on Implementation

In the adjusted model (4), τj,j=1,…,K0 are used. When the population distribution is not available, they need to be approximated by estimators. When Z is normally distributed and eigenvalues rj,j=1,…,q of the covariance matrix ΣZ are different from each other, then n(uj-τj) is asymptotically N(0,Vj) with Vj=∑l≠jq‍(rjrl/(rj-rl)2)τlτlT, where uj is the jth eigenvector of S=(1/(n-1))∑i=1n‍(Zi-Z-)(Zi-Z-)T with Z-=(1/n)∑i=1n‍Zi; see Anderson [11]. For the case when the population size is large and comparable with the sample size, if the covariance matrix is sparse, we can use the method in Rütimann and Bühlmann [12] or Cai and Liu [13] to estimate the covariance matrix. So we can use uj to approximate τj. When τj in model (4) are replaced by these consistent estimators, one can see that the approximation error can be neglected without changing the asymptotic property.The nonparametric partsgl(Z~(l)) in the adjusted model depend on the univariate variable Z~(l), for l=1,…,K0. So it needs to choose the steps K0 firstly. In real implementation, we compute all the q multiple correlation coefficients of Z~(l) (l=1,…,q) with X and U. Then we choose the components R={Z~(l):|mcorr(Z~(l),(X,U))|≥δ,l=1,…,q} for given small number δ>0, where mcorr(u,V) denotes the multicorrelation coefficient between u and V and can be approximated by its sample form; see Anderson [11].There are some tuning parameters needing to choose in order to implement the two-stage remodeling procedure. Fan and Li [5] showed that the SCAD penalty with a=3.7 performs well in a variety of situations. Hence, we use their suggestion throughout this paper. We still need to choose the positive integer L for basis functions and the tuning parameter λj of the penalty functions. Similar to the adaptive lasso of Zou [14], we suggest taking λj=λ/∥θ^j(0)∥2, where θ^j(0) is initial estimator of θj by using ordinary least squares method based on the first term in (10). So the two remaining parameters L and λ can be selected simultaneously using the leave-one-out CV or GCV method; see Zhao and Xue [4] for more details.


## 3. Simulation Studies

In this section, we investigate the behavior of the newly proposed method by simulation studies.


### 3.1. Linear Model with Normally Distributed Covariates

The dimensions of the full model (1) and the submodel (2) are chosen to be 100 and 5, respectively. We set β=(0.5,3.5,2.5,1.5,4.0)T and γ=(γ1,γ2,055T)T, where γ2~Unif[-0.5,0.5]30, a 30-dimensional uniform distribution on [-0.5,0.5]30, and γ1 is chosen in the following ways:Case (I). γ1~Unif[0.5,1.0]10.Case (II). γ1=(1.0,1.0,1.0,1.5,1.5,1.5,2.0,2.0,2.0,2.0).We assume that(XT,ZT)T~N((15T,040T,155T)T,ΣΣT), where
(15)Σ=(σij),σij=σji={1.0,j=i,i=1,…,p+q;c,j=i+p,i=1,3,…,q;0,otherwise,
with c=0.5 or c=0.8. The error term ε is assumed to be normally distributed as N(0,0.32).Here we denote the submodel (2) as model (I), the multistep adjusted linear model (5) as model (II), the two-stage model (12) as model (III), and the full model (1) as model (IV). We compare mean square errors (MSEs) of the new two-stage estimator β^TS based on model (III) with the estimator β^S based on model (I), the multistep estimator β^M based on model (II), the SCAD estimator β^SCAD and the least squares estimator β^F based on model (IV). We also compare mean square prediction errors (MSPEs) of the above mentioned models with corresponding estimators.The data are simulated from the full model (1) with sample size n=100 and simulation times m=1000. We use the sample-based PCA approximations to substitute τj's. The parameter a in the SCAD penalty function is set to be 3.7 and λ is selected by leave-one-out CV method.Table1 reports the MSEs of point estimators on the parameter β and the MSPEs of model predictions. From the table, we have the following findings: (1)  β^F has the largest MSEs and β^S takes the second place, nearly all the new estimator β^TS has the smallest MSEs. (2) When c=0.5, the MSEs of β^SCAD are smaller than those of β^M, while when c=0.8 they are larger than those of β^M. These show that if the correlation between the covariates is strong, the MSEs of β^SCAD are larger than those of β^M, the multistep-adjustment is necessary, so the estimations and model predictions based on two-stage model are significantly improved. (3) In case (I) and (II) the simulation results have the similar performance. (4) Similar to the trend of the MSEs of the five estimators, the MSPE of the two-stage adjusted model is the smallest among the mentioned five models.Table 1
MSEs on the parameterβ and MSPEs of the two-stage adjusted linear model (12) compared with the submodel, the SCAD-penalized model, the multistep adjusted model and the full model.
No.
Item
β
^
S
β
^
SCAD
β
^
M
β
^
TS
β
^
F
0.3079
0.0457
0.0660
0.0571
1.6105
×
1
0
3
0.1763
0.0206
0.0346
0.0176
1.0940
×
1
0
3
Case (I)
MSEs
0.1396
0.0481
0.0631
0.0461
4.2049
×
1
0
3
c
=
0.5
0.1870
0.0196
0.0349
0.0186
5.0183
×
1
0
3
0.1131
0.0517
0.0609
0.0430
6.2615
×
1
0
3
MSPEs
3.4780
1.1896
1.6512
1.0679
3.0499
×
1
0
2
0.1568
0.6191
0.0934
0.0826
1.2494
×
1
0
3
0.6239
0.1060
0.0090
0.0083
1.0456
×
1
0
2
Case (I)
MSEs
0.8829
0.8173
0.0895
0.1039
2.6368
×
1
0
2
c
=
0.8
0.5882
0.0919
0.0107
0.0100
7.6452
×
1
0
1
1.0799
0.9829
0.0961
0.0929
1.1610
×
1
0
3
MSPEs
4.7930
2.6700
0.8354
0.7771
1.3223
×
1
0
2
0.4272
0.0660
0.0849
0.0557
4.3002
×
1
0
2
0.6371
0.0318
0.0499
0.0295
3.7893
×
1
0
3
Case (II)
MSEs
0.4560
0.0715
0.0927
0.0588
1.2784
×
1
0
3
c
=
0.5
0.5926
0.0306
0.0491
0.0287
6.7354
×
1
0
3
0.9052
0.0734
0.0874
0.0583
2.5047
×
1
0
2
MSPEs
6.8634
1.5096
2.0780
1.2077
5.0464
×
1
0
3
0.6764
0.4263
0.1212
0.0960
1.3904
×
1
0
3
0.9721
0.1060
0.0107
0.0102
4.0743
×
1
0
2
Case (II)
MSEs
0.6242
0.4756
0.1146
0.1003
1.0498
×
1
0
3
c
=
0.8
1.0282
0.0954
0.0112
0.0098
5.6031
×
1
0
2
1.3420
0.5474
0.1341
0.1124
9.9632
×
1
0
2
MSPEs
7.9928
2.1165
0.9514
0.8469
2.3110
×
1
0
2In summary, Table1 indicates that the two-stage adjusted linear model (12) performs much better than the full model, and better than the submodel, the SCAD-penalized model and the multistep adjusted model.


### 3.2. Partially Linear Model with Nonnormally Distributed Covariates

The dimensions of the linear part in the full model (1) and the submodel (2) are chosen to be 50 and 5, respectively. We set β=(0.5,3.5,2.5,1.5,4.0)T, γ=(γ1,γ2,025T)T, f(u)=u2*sin(3u), where γ1=(0.5,0.1,0.8,0.2,0.5,0.2,0.6,0.5,0.1,0.9),γ2~Unif[-0.3,0.3]10, a 10-dimensional uniform distribution on [-0.3,0.3].We assume that the covariates are distributed in the following two ways.Case (I). (XT,ZT,U)T~t(051T,ΣΣT), a 51-dimensional student distribution with degree of freedom df=5, where
(16)Σ=(σij),σij=σji={1.0,j=i,i=1,…,p+q+1;0.95,j=i+p,i=1,2,…,q+1;0.9,j=i+p-2,i=1,2,…,q+3;0,otherwise.Case (II). X=(1/(1+c))(W1+cV), Z=(Z1T,Z2T,Z3T,Z4T)T with Z1=(1/(1+c))(W2+cV), Z2=W3, Z3=(1/(1+c))(W4+cV), Z4=W5, U=W5(1), where W1,W2,W3,W4~Unif[-1.0,1.0]5, W5~Unif[-1.0,1.0]30, V~Unif[-1.0,1.0]5, uniform distributions on [-1.0,1.0] and constant c=0.1. All W1, W2, W3, W4, W5, and V are independent.The error termε is assumed to be normally distributed as N(0,0.32).Here we denote the submodel (2) as model (I)′, the multistep adjusted additive partially linear model (3) as model (II)′, the two-stage model (11) as model (III)′ and the full model (1) as model (IV)′. We compare mean square errors (MSEs) of the new two-stage estimator β^TS based on model (III)′ with the estimator β^S based on model (I)′, the estimator β^M based on model (II)′ and the least squares estimator β^F based on model (IV)′. We also compare the mean average square errors (MASEs) of the nonparametric estimators of f(·) and the mean square prediction errors (MSPEs) of different models with corresponding estimators.The data are simulated from the full model (1) with sample size n=100 and simulation times m=500. We use the sample-based approximations of ICA, see Hyvärinen and Oja [15]. The parameter a in the SCAD penalty function is set to be 3.7, the number L and the parameter λ is selected by GCV method. We use the standard Fourier orthogonal basis as the basis functions.Table2 reports the MSEs of point estimators on the parameter β, the MASEs of f(·) and the MSPEs of model predictions. From the table, we have the following results: (1)  β^F has the largest MSEs, its MSEs are much larger than the MSEs of the other estimators, and the new estimator β^TS always has the smallest MSEs. (2) The MASEs of f(·) have similar trend to the MSEs of the four estimators, while the differences are not very noticeable. (3) Similar to the MSEs of the estimators, the MSPEs of the two-stage adjusted model are the smallest among the four models. (4) In Case (II), the simulation results of models (I)′, (II)′ and (III)′, perform a little better than those in Case (I) because of the correlation structure among the covariates.Table 2
MSEs on the parameterβ, MASEs of f(·) and MSPEs of the two-stage adjusted model (11) compared with the submodel, multistep adjusted model and the full model.
No.
Item
β
^
S
β
^
M
β
^
TS
β
^
F
0.4352
5.0403
0.3267
2.9753
×
1
0
1
0.6859
1.2820
×
1
0
1
0.3328
1.4593
×
1
0
1
MSEs
1.1152
8.1542
0.3723
1.4391
×
1
0
1
Case (I)
1.8489
7.2055
1.3194
2.4036
×
1
0
1
3.3079
1.6144
×
1
0
1
1.9989
4.8575
×
1
0
1
MASEs
3.0887
5.9814
3.0175
3.0633
MSPEs
4.6047
7.0331
×
1
0
1
3.5536
3.9648
0.0377
0.6144
0.0191
—1
0.0449
1.0876
0.0305
—
MSEs
0.0332
3.7510
0.0246
—
Case (II)
0.0396
0.4324
0.0238
—
0.0512
1.1995
0.0335
—
MASEs
0.4722
0.5220
0.4126
0.4380
MSPEs
0.9221
9.3068
0.8053
—
1“—” denotes the algorithm collapsed and returned no value.In summary, Table2 indicates that the two-stage adjusted model (11) performs much better than the full model and the multistep adjusted model, and better than the submodel.


## 3.1. Linear Model with Normally Distributed Covariates

The dimensions of the full model (1) and the submodel (2) are chosen to be 100 and 5, respectively. We set β=(0.5,3.5,2.5,1.5,4.0)T and γ=(γ1,γ2,055T)T, where γ2~Unif[-0.5,0.5]30, a 30-dimensional uniform distribution on [-0.5,0.5]30, and γ1 is chosen in the following ways:Case (I). γ1~Unif[0.5,1.0]10.Case (II). γ1=(1.0,1.0,1.0,1.5,1.5,1.5,2.0,2.0,2.0,2.0).We assume that(XT,ZT)T~N((15T,040T,155T)T,ΣΣT), where
(15)Σ=(σij),σij=σji={1.0,j=i,i=1,…,p+q;c,j=i+p,i=1,3,…,q;0,otherwise,
with c=0.5 or c=0.8. The error term ε is assumed to be normally distributed as N(0,0.32).Here we denote the submodel (2) as model (I), the multistep adjusted linear model (5) as model (II), the two-stage model (12) as model (III), and the full model (1) as model (IV). We compare mean square errors (MSEs) of the new two-stage estimator β^TS based on model (III) with the estimator β^S based on model (I), the multistep estimator β^M based on model (II), the SCAD estimator β^SCAD and the least squares estimator β^F based on model (IV). We also compare mean square prediction errors (MSPEs) of the above mentioned models with corresponding estimators.The data are simulated from the full model (1) with sample size n=100 and simulation times m=1000. We use the sample-based PCA approximations to substitute τj's. The parameter a in the SCAD penalty function is set to be 3.7 and λ is selected by leave-one-out CV method.Table1 reports the MSEs of point estimators on the parameter β and the MSPEs of model predictions. From the table, we have the following findings: (1)  β^F has the largest MSEs and β^S takes the second place, nearly all the new estimator β^TS has the smallest MSEs. (2) When c=0.5, the MSEs of β^SCAD are smaller than those of β^M, while when c=0.8 they are larger than those of β^M. These show that if the correlation between the covariates is strong, the MSEs of β^SCAD are larger than those of β^M, the multistep-adjustment is necessary, so the estimations and model predictions based on two-stage model are significantly improved. (3) In case (I) and (II) the simulation results have the similar performance. (4) Similar to the trend of the MSEs of the five estimators, the MSPE of the two-stage adjusted model is the smallest among the mentioned five models.Table 1
MSEs on the parameterβ and MSPEs of the two-stage adjusted linear model (12) compared with the submodel, the SCAD-penalized model, the multistep adjusted model and the full model.
No.
Item
β
^
S
β
^
SCAD
β
^
M
β
^
TS
β
^
F
0.3079
0.0457
0.0660
0.0571
1.6105
×
1
0
3
0.1763
0.0206
0.0346
0.0176
1.0940
×
1
0
3
Case (I)
MSEs
0.1396
0.0481
0.0631
0.0461
4.2049
×
1
0
3
c
=
0.5
0.1870
0.0196
0.0349
0.0186
5.0183
×
1
0
3
0.1131
0.0517
0.0609
0.0430
6.2615
×
1
0
3
MSPEs
3.4780
1.1896
1.6512
1.0679
3.0499
×
1
0
2
0.1568
0.6191
0.0934
0.0826
1.2494
×
1
0
3
0.6239
0.1060
0.0090
0.0083
1.0456
×
1
0
2
Case (I)
MSEs
0.8829
0.8173
0.0895
0.1039
2.6368
×
1
0
2
c
=
0.8
0.5882
0.0919
0.0107
0.0100
7.6452
×
1
0
1
1.0799
0.9829
0.0961
0.0929
1.1610
×
1
0
3
MSPEs
4.7930
2.6700
0.8354
0.7771
1.3223
×
1
0
2
0.4272
0.0660
0.0849
0.0557
4.3002
×
1
0
2
0.6371
0.0318
0.0499
0.0295
3.7893
×
1
0
3
Case (II)
MSEs
0.4560
0.0715
0.0927
0.0588
1.2784
×
1
0
3
c
=
0.5
0.5926
0.0306
0.0491
0.0287
6.7354
×
1
0
3
0.9052
0.0734
0.0874
0.0583
2.5047
×
1
0
2
MSPEs
6.8634
1.5096
2.0780
1.2077
5.0464
×
1
0
3
0.6764
0.4263
0.1212
0.0960
1.3904
×
1
0
3
0.9721
0.1060
0.0107
0.0102
4.0743
×
1
0
2
Case (II)
MSEs
0.6242
0.4756
0.1146
0.1003
1.0498
×
1
0
3
c
=
0.8
1.0282
0.0954
0.0112
0.0098
5.6031
×
1
0
2
1.3420
0.5474
0.1341
0.1124
9.9632
×
1
0
2
MSPEs
7.9928
2.1165
0.9514
0.8469
2.3110
×
1
0
2In summary, Table1 indicates that the two-stage adjusted linear model (12) performs much better than the full model, and better than the submodel, the SCAD-penalized model and the multistep adjusted model.


## 3.2. Partially Linear Model with Nonnormally Distributed Covariates

The dimensions of the linear part in the full model (1) and the submodel (2) are chosen to be 50 and 5, respectively. We set β=(0.5,3.5,2.5,1.5,4.0)T, γ=(γ1,γ2,025T)T, f(u)=u2*sin(3u), where γ1=(0.5,0.1,0.8,0.2,0.5,0.2,0.6,0.5,0.1,0.9),γ2~Unif[-0.3,0.3]10, a 10-dimensional uniform distribution on [-0.3,0.3].We assume that the covariates are distributed in the following two ways.Case (I). (XT,ZT,U)T~t(051T,ΣΣT), a 51-dimensional student distribution with degree of freedom df=5, where
(16)Σ=(σij),σij=σji={1.0,j=i,i=1,…,p+q+1;0.95,j=i+p,i=1,2,…,q+1;0.9,j=i+p-2,i=1,2,…,q+3;0,otherwise.Case (II). X=(1/(1+c))(W1+cV), Z=(Z1T,Z2T,Z3T,Z4T)T with Z1=(1/(1+c))(W2+cV), Z2=W3, Z3=(1/(1+c))(W4+cV), Z4=W5, U=W5(1), where W1,W2,W3,W4~Unif[-1.0,1.0]5, W5~Unif[-1.0,1.0]30, V~Unif[-1.0,1.0]5, uniform distributions on [-1.0,1.0] and constant c=0.1. All W1, W2, W3, W4, W5, and V are independent.The error termε is assumed to be normally distributed as N(0,0.32).Here we denote the submodel (2) as model (I)′, the multistep adjusted additive partially linear model (3) as model (II)′, the two-stage model (11) as model (III)′ and the full model (1) as model (IV)′. We compare mean square errors (MSEs) of the new two-stage estimator β^TS based on model (III)′ with the estimator β^S based on model (I)′, the estimator β^M based on model (II)′ and the least squares estimator β^F based on model (IV)′. We also compare the mean average square errors (MASEs) of the nonparametric estimators of f(·) and the mean square prediction errors (MSPEs) of different models with corresponding estimators.The data are simulated from the full model (1) with sample size n=100 and simulation times m=500. We use the sample-based approximations of ICA, see Hyvärinen and Oja [15]. The parameter a in the SCAD penalty function is set to be 3.7, the number L and the parameter λ is selected by GCV method. We use the standard Fourier orthogonal basis as the basis functions.Table2 reports the MSEs of point estimators on the parameter β, the MASEs of f(·) and the MSPEs of model predictions. From the table, we have the following results: (1)  β^F has the largest MSEs, its MSEs are much larger than the MSEs of the other estimators, and the new estimator β^TS always has the smallest MSEs. (2) The MASEs of f(·) have similar trend to the MSEs of the four estimators, while the differences are not very noticeable. (3) Similar to the MSEs of the estimators, the MSPEs of the two-stage adjusted model are the smallest among the four models. (4) In Case (II), the simulation results of models (I)′, (II)′ and (III)′, perform a little better than those in Case (I) because of the correlation structure among the covariates.Table 2
MSEs on the parameterβ, MASEs of f(·) and MSPEs of the two-stage adjusted model (11) compared with the submodel, multistep adjusted model and the full model.
No.
Item
β
^
S
β
^
M
β
^
TS
β
^
F
0.4352
5.0403
0.3267
2.9753
×
1
0
1
0.6859
1.2820
×
1
0
1
0.3328
1.4593
×
1
0
1
MSEs
1.1152
8.1542
0.3723
1.4391
×
1
0
1
Case (I)
1.8489
7.2055
1.3194
2.4036
×
1
0
1
3.3079
1.6144
×
1
0
1
1.9989
4.8575
×
1
0
1
MASEs
3.0887
5.9814
3.0175
3.0633
MSPEs
4.6047
7.0331
×
1
0
1
3.5536
3.9648
0.0377
0.6144
0.0191
—1
0.0449
1.0876
0.0305
—
MSEs
0.0332
3.7510
0.0246
—
Case (II)
0.0396
0.4324
0.0238
—
0.0512
1.1995
0.0335
—
MASEs
0.4722
0.5220
0.4126
0.4380
MSPEs
0.9221
9.3068
0.8053
—
1“—” denotes the algorithm collapsed and returned no value.In summary, Table2 indicates that the two-stage adjusted model (11) performs much better than the full model and the multistep adjusted model, and better than the submodel.


## 4. Some Remarks

In this paper, the main objective is to consistently estimate the parameter of interestβ. When estimating the parameter of interest, its bias is mainly determined by the relevant variables, and its variance may be impacted by other variables. Because variable selection much relies on the sparsity of the parameter, when we directly consider the partially linear model, some irrelevant variables with nonzero coefficients may be selected in the final model. This may affect the estimation of the parameter β on its efficiency and stability. Thus based on the prespecified submodel, a two-stage remodeling method is proposed. In the new remodeling procedure, the correlation among the covariates (X,Z) and the sparsity of the regression structure are fully used. So the final model is sufficiently simplified and conditionally unbiased. Based on the simplified model, the estimation and model prediction are significantly improved. Generally, after the first stage the adjusted model is an additive partially linear model. Therefore, the remodeling method can be applied to partially linear regression model with linear regression model as a special case.From the remodeling procedure, we can see that it can be directly applied to additive partially linear model, in which the nonparametric functionf(U) has component-wise additive form. As for general partially linear model with multivariate nonparametric function, we should resort to multivariate nonparametric estimation method. If the dimension of covariate U is high, it may be faced with “the curse of dimensionality”.In the procedure of model simplification, orthogonal series estimation method is used. This is only for technical convenience, because the semiparametric penalized least squares (6) can be easily transformed into parametric penalized least squares (10) and then the theoretic results are obtained. Although other nonparametric methods such as kernel and spline can be used without any essential difficulty, they can not directly achieve this goal. Compared with kernel method, it is somewhat difficult for series method to establish the asymptotic normality result for the nonparametric component f(U) under primitive conditions.



---

*Source: 687151-2013-02-06.xml*