# Automatic Variable Selection for Partially Linear Functional Additive Model and Its Application to the Tecator Data Set

Authors: Yuping Hu; Sanying Feng; Liugen Xue
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

We introduce a new partially linear functional additive model, and we consider the problem of variable selection for this model. Based on the functional principal components method and the centered spline basis function approximation, a new variable selection procedure is proposed by using the smooth-threshold estimating equation (SEE). The proposed procedure automatically eliminates inactive predictors by setting the corresponding parameters to be zero and simultaneously estimates the nonzero regression coefficients by solving the SEE. The approach avoids the convex optimization problem, and it is flexible and easy to implement. We establish the asymptotic properties of the resulting estimators under some regularity conditions. We apply the proposed procedure to analyze a real data set: the Tecator data set.


---

## Body

## 1. Introduction

Functional data may be viewed as realization of observed stochastic processes, and it is commonly encountered in many fields of applied sciences, such as econometrics, biomedical studies, and physics experiment. The Tecator data set is collected by the Tecator company and is publicly available athttp://lib.stat.cmu.edu/datasets/tecator. This data set consists of 215 meat samples. The measurements were made through a spectrometer named the Tecator Infratec Food and Feed Analyzer, and the spectral curves were recorded at wavelengths ranging from 850 nm to 1050 nm. For each meat sample, the data consists of a 100 channel spectrum of absorbances as well as the contents of moisture (water), fat, and protein. The three contents of fat, protein, and moisture (water), measured in percentages, are determined by analytic chemistry. We aim to predict the fat content of a meat sample. In this paper, we propose a new partially linear functional additive model and apply the SEE procedure to analyze the Tecator data set.With the development of computer technology, much progress has been made on developing methodologies for analyzing functional data by many researchers, like Ramsay and Silverman [1], Cardot, Ferraty, and Sarda [2], Lian and Li [3], Fan, James, and Radchenko [4], Feng and Xue [5], Yu, Zhang, and Du [6], Zhou, Du, and Sun [7], among others. Regression models play a major role in the functional data analysis. The most widely used regression model is the following functional linear model: (1)Y=∫TβtXtdt+ε,where Y is a scalar response, functional predictor X(t) is a smooth and square-integrable random function defined on a compact domain T=[0,1] for simplicity, β(t) is the square-integrable regression parameter function, and ε is a random error, which is independent of X(t). A commonly adopted approach for fitting model is the basis expansion; that is, X(t)=μ(t)+∑k=1∞ξkϕk(t), where μ(t)=E{X(t)}. Model (1) is then transformed to a linear form with the coefficients ξk: EY∣X=a0+∑k=1∞ξkbk, where a0=∫01β(t)μ(t)dt and bk=∫01β(t)ϕk(t)dt. The basis function set {ϕk(t)} can be either predetermined (e.g., Fourier basis, wavelets, or B-splines basis) or data-driven. One convenient choice for data-driven basis is the eigenbasis of the autocovariance operator of X(t), in which case the random coefficients {ξk} are called the functional principal component (FPC) scores. The FPC scores have zero means and variances equal to the corresponding eigenvalues {λk,k=1,2,···}. We focus on the FPC representation of the functional regression throughout this paper.Müller and Yao [8] relaxed the linearity assumption and proposed a functional additive model (FAM). This leads to a more widely applicable and flexible framework for the functional regression models. In the case of scalar response, the linear structure is replaced by the sum of nonlinear functional components; that is, (2)EY∣X=a0+∑k=1∞fkξk,where fk(ξk) are unknown smooth functions. Model (2) was fitted by estimating {ξk} by the functional principal component analysis, and by estimating {fk} by the local polynomial smoothing in Müller and Yao [8]. Zhu, Yao, and Zhang [9] proposed a new regularization framework for the structure estimation of FAM in the context of reproducing kernel Hilbert spaces. The selection was achieved by the penalized least squares using a penalty which encourages the sparse structure of the additive components, and the rate of convergence was investigated.However, in many real world problems, it is common to collect information on a large number of nonfunctional predictors. How to incorporate scalar predictors into the functional regression and perform model selection are important issues. In this paper, we combine the linear model with the functional additive model and introduce a new partially linear functional additive model (PLFAM).Traditional real-value additive models were studied in Stone [10], Wang and Yang [11], Huang, Horowitz and Wei [12], and Zhao and Xue [13]. When the explanatory variables are of functional nature, Ferraty and Vieu [14] used a two-step procedure to estimate an additive model with two functional predictors. Fan, James, and Radchenko [4] suggested a new penalized least squares method to fit the nonlinear functional additive model. This method can efficiently fit the high-dimensional functional models while simultaneously performing variable selection to identify the relevant predictors. Febrero-Bande and Gonzalez-Manteiga [15] extended the ideas of the generalized additive models with multivariate data to the functional data covariates. The proposed algorithm was a modified version of the local scoring and backfitting algorithm that allows for the nonparametric estimation of the link function.In the last decades, variable selection has received substantial attention, which has been a very important topic in regression analysis. Generally speaking, most of the variable selection procedures are based on penalized estimation based on some penalty functions, like Lasso penalty [16], SCAD penalty [17], Adaptive Lasso [18], and so on. However, these penalty functions have a singularity at zero such that these penalized estimation procedures require convex optimization, while adding the burden of computation. To overcome this problem, Ueki [19] developed a new variable selection procedure called the smooth-threshold estimating equation that can automatically eliminate irrelevant parameters by setting them as zero. The method has also been successfully applied to a large class of models; for example, Lai, Wang and Lian [20] explored the generalized the estimation equation (GEE) estimation and the smooth-threshold generalized estimation equation (SGEE) variable selection for single-index model with clustered data. Li et al. [21] considered the SGEE variable selection for the generalized linear model with longitudinal data. Tian, Xue, and Xu [22] proposed a smooth-threshold estimating equation variable selection for varying coefficient models with longitudinal data.As we know, functional regression models have been widely applied to engineering problems. For example, Escabias, Aguilera, and Valderrama [23] used functional logistic regression to deal with the environmental problem, which is to estimate the risk of drought in a specific zone from time evolution of temperatures. Sonja, Branimir, and DraDen [24] dealt with tool wear in milling process and the prediction of its behaviour by utilizing functional data analysis (FDA) methodology. Pokhrel and Tsokos [25] applied functional data analysis techniques to model age-specific brain cancer mortality trend and forecast entire age-specific functions using exponential smoothing state-space models.In this article, we propose a new functional regression model and consider the variable selection problem for this model; then we apply the proposed procedure to analyze the Tecator data set. Motivated by the idea in Ueki [19], based on the functional principal components analysis and the centered spline basis function approximation, an automatic variable selection procedure is proposed using the smooth-threshold estimating equation. The proposed procedure automatically eliminates the irrelevant parameters in the model, while estimating the nonzero regression coefficients. Our approach can be easily implemented without solving any convex optimization problems, and it reduces the burden of computation. The proposed method shares some of the desired features including the oracle property. The proposed smooth-threshold estimating equation approach is flexible and easy to implement. Finally, the proposed method is applied to analyze a real data set: the Tecator data set. The validity of the partially linear functional additive model and the SEE method are confirmed.The rest of this paper is organized as follows. In Section2, we propose a variable selection procedure for PLFAM and study the asymptotic properties under some regularity conditions. In Section 3, we give the computation of the estimators as well as the choice of the tuning parameters. In Section 4, we apply the proposed method to analyze the Tecator data set. Concluding remarks are presented in Section 5. The technical proofs of all asymptotic results are provided in the Appendix.


## 2. Methodology and Main Results

LetY be a real-valued random variable, Z be a p-dimensional vector of random variables, and X(t) be a zero mean and square-integrable random function defined on interval [0,1]. Assume that {Yi,Zi,Xi(·)}i=1n are independent identically and distributed realizations of the pair {Y,Z,X(·)}. Denote by ξi,∞=(ξi1,ξi2,···)T the FPC scores sequence of Xi(t), which is associated with eigenvalues {λ1,λ2,···} with λ1≥λ2≥···≥0.


### 2.1. Smooth-Threshold Estimating Equation

For the convenience of model regularization, we would like to restrict the predictor variables to take values in[0,1] without loss of generality. This is achieved by taking a transformation of the FPC scores through a cumulative distribution function (CDF) Ψ:R→[0,1] for all {ξik}. We take the normal CDF denoted by Ψ(·,λk), with zero mean and variance λk. It is easy to see that, if ξiks follow normal distribution, the normal CDF leads to uniformly distributed transformed variables on [0,1].Denoting the transformed variable of{ξik} by {ζik}, i.e., ζik=Ψ(ξik,λk), and denoting ζi,∞=(ζi1,ζi2,···)T, we propose a partially linear functional additive model as follows: (3)Yi=ZiTβ+∑k=1∞f0kζik+εi,where εi is independent error with zero mean and variance σ2, β is a p×1 vector of unknown regression coefficients, and f0(ζi∞)=∑k=1∞f0k(ζik) is a smooth function. To ensure identifiability, we assume that Ef0k(ζik)=0. In this paper, we assume that the PLFAM is the sparse structure of the underlying true model, and this assumption is critical in the context of functional data analysis. It means that the number of important functional additive components that contribute to the response is finite, but not necessarily restricted to the leading terms. In particular, we denote by I the index set of the important FPC scores and assume that I<∞, where · denotes the cardinality of a set. In other words, there is a sufficiently large q such that I⊆{1,···,q}, which implies that fk≡0 as long as k>q. Model (3) is thus equivalent to (4)Yi=ZiTβ+∑k=1qf0kζik+εi.We replace f0k(ζik) with its basis function approximations. Let B(ζik)=(B1(ζik),···,BL(ζik))T be the centralized B-spline basis functions with the order of M, where L=K+M, and K is the number of interior knots. Thus, f0k(ζik) can be approximated by (5)f0kζik≈BζikTγk,k=1,···,q. Substituting this into model (4), we have (6)Yi≈ZiTβ+∑k=1qBζikTγk+εi.LetPi=(ZiT,WiT), WiT=(B(ζi1)T,···,B(ζiq)T), and θ=(θ1,···,θs)T=(βT,γ1T,···,γqT)T, s=p+qL. Due to the unknown ζik, we substitute the estimator of ζik by ζ^ik=Ψ(ξ^ik,λ^k), where ξ^ik and λ^k are the estimators of ξik and λk, respectively. Let P^i=(ZiT,W^iT), W^iT=(B(ζ^i1)T,···,B(ζ^iq)T), and define the following estimating function of θ(7)Uθ=∑i=1nP^iTYi-P^iθ.Motivated by the idea of Ueki [19], we propose the following smooth-threshold estimating equation: (8)Is-ΔUθ+Δθ=0,where s=p+qL, Is is the s-dimensional identity matrix, and Δ is the diagonal matrix, i.e., Δ=diag(Δ1,Δ2), where Δ1=diag{δ11,···,δ1p} is a p×p diagonal matrix, and Δ2=diag{δ21,···,δ21,···,δ2q,···,δ2q} is a qL×qL diagonal matrix. Note that δ1j=1 reduces to βj=0 for j=1,···,p, and δ2k=1 reduces to γk=0, that is, f0k(·)=0, for k=1,···,q. Therefore, (8) can yield a sparse solution. Unfortunately, we cannot directly obtain the estimator of θ by solving the smooth-threshold estimating equation (8). The reason is that (8) involves the unknown parameters δ1j and δ2k, which need to be chosen using some data-driven criteria.For the choice ofδ1j and δ2k, Ueki [19] suggested that δ1j and δ2k may be determined by the data, δ^1j=min(1,ν/β^j01+τ) and δ^2k=min(1,ν/γ^k(0)1+τ), with initial estimators β^j(0) and γ^k(0), respectively. The initial estimators β^j(0) and γ^k(0) are the solutions of the estimating equation U(θ)=0. Note that this choice involves two tuning parameters (ν,τ). In Section 3, we will propose a BIC-type criterion to select the tuning parameters. Replacing Δ in (8) by Δ^ with diagonal elements δ^1j and δ^2k, the smooth-threshold estimating equation (8) becomes (9)Is-Δ^Uθ+Δ^θ=0.The solution of (9) denoted by θ^=(θ^1,⋯,θ^p,θ^p+1,⋯,θ^s)T=(β^T,γ^1T,···,γ^qT)T is called the SEE estimator. Thus, β^ is the estimator of β, and Bk(ζ^k)Tγ^k is the estimator of f0k(ζk). For the convenience of the notations, we define A1={j:δ^1j≠1} as the set of indices of nonzero regression coefficient estimators, and define A2={k:δ^2k≠1} as the set of nonzero function estimators. A1c={j:δ^1j=1} and A2c={k:δ^2k=1} are the sets of zero regression coefficient estimators and zero function estimators, respectively. We further denote A=A1∪A2, and Ac=A1c∪A2c. The SEE estimator is defined in the following form:(10)θ^A=∑i=1nP^i,ATP^i,A-G^A-1∑i=1nP^i,ATYi,θ^Ac=0, Where P^i,A=(Zi,A1T,W^i,A2T), G^A=(IA-Δ^A)-1Δ^A, Δ^A=diag{Δ^1,A1,Δ^2,A2}, and IA is the identity matrix with the same dimension of Δ^A.


### 2.2. Asymptotic Properties

We first introduce some notations. Letf0(·) and β0 denote the true values of f(·) and β, respectively. γ0 is the spline coefficient vector from the spline approximation to f0(·). Denote θ0=(β0T,γ0T)T, A10={j:β0j≠0},A10c={j:β0j=0}, A20={k:f0k(·)≠0}, and A20c={k:f0k(·)=0}. We also denote A0=A10∪A20 and A0c=A10c∪A20c. The norm of f(·) is defined as f(x)=(∫01f2(x)dx)1/2.For convenience and simplicity, letC denote a positive constant that may be different at each appearance throughout this paper. We list some regularity conditions that are used in this paper.(C1) The transformation functionΨ(ξ,λk) is differentiable at ξ and λk, and satisfies the fact that ∂Ψ(ξ,λk)/∂ξ≤Cλkc0 and ∂Ψ(ξ,λk)/∂λk≤Cλkc0ξ for some positive constants C and negative constants c0.(C2) The second derivativeX(2)(t) is continuous on [0,1] with probability 1 and ∫E[X(k)(t)]4dt<∞ with probability 1, for k=0,2.(C3) The spline regression parameter is identified; that is, there is a uniqueθ0=(β0T,γ01T,⋯,γ0qT)T∈S, where the parameter space S is compact.(C4) The inner knots{ci,i=1,⋯,K} satisfy (11)max1≤i≤Khi+1-hi=oK-1,hmin1≤i≤Khi≤C0, where hi=ci-ci-1, h=max1≤i≤Khi, c0=0, and cK+1=1.(C5)f0k(·) is the rth continuously differentiable on (0,1), where r≥2.The following theorem gives the consistency of our proposed estimators.Theorem 1.
Suppose that regularity conditions (C1)-(C5) hold. For any positiveν and τ, nr/(2r+1)ν→0 and nr(1+τ)/(2r+1)ν→∞, as n→∞, and if K=O(n1/(2r+1)), one has (12)f^0kζ^k-f0kζk=Opn-r/2r+1,k=1,⋯,q.In the following theorem, we will show that such consistent estimators enjoy the sparsity property.Theorem 2.
Under the regularity conditions of Theorem1, one has(i)
P(A1c=A10c)→1,(ii)
P(A2c=A20c)→1.Next, we will show that the estimators of the nonzero coefficients for the parametric components have the same asymptotic distribution as that based on the correct model.Theorem 3.
Suppose that the regularity conditions of Theorem1 hold; as n→∞, one has, (13)Dβ^A1-1/2β^A1-β0A10⟶LN0,IA1,where D(β^A1) is defined in (A.23) in the Appendix, A1 is the number of A1, IA1 is the identity matrix, and “→L” represents the convergence in distribution.


## 2.1. Smooth-Threshold Estimating Equation

For the convenience of model regularization, we would like to restrict the predictor variables to take values in[0,1] without loss of generality. This is achieved by taking a transformation of the FPC scores through a cumulative distribution function (CDF) Ψ:R→[0,1] for all {ξik}. We take the normal CDF denoted by Ψ(·,λk), with zero mean and variance λk. It is easy to see that, if ξiks follow normal distribution, the normal CDF leads to uniformly distributed transformed variables on [0,1].Denoting the transformed variable of{ξik} by {ζik}, i.e., ζik=Ψ(ξik,λk), and denoting ζi,∞=(ζi1,ζi2,···)T, we propose a partially linear functional additive model as follows: (3)Yi=ZiTβ+∑k=1∞f0kζik+εi,where εi is independent error with zero mean and variance σ2, β is a p×1 vector of unknown regression coefficients, and f0(ζi∞)=∑k=1∞f0k(ζik) is a smooth function. To ensure identifiability, we assume that Ef0k(ζik)=0. In this paper, we assume that the PLFAM is the sparse structure of the underlying true model, and this assumption is critical in the context of functional data analysis. It means that the number of important functional additive components that contribute to the response is finite, but not necessarily restricted to the leading terms. In particular, we denote by I the index set of the important FPC scores and assume that I<∞, where · denotes the cardinality of a set. In other words, there is a sufficiently large q such that I⊆{1,···,q}, which implies that fk≡0 as long as k>q. Model (3) is thus equivalent to (4)Yi=ZiTβ+∑k=1qf0kζik+εi.We replace f0k(ζik) with its basis function approximations. Let B(ζik)=(B1(ζik),···,BL(ζik))T be the centralized B-spline basis functions with the order of M, where L=K+M, and K is the number of interior knots. Thus, f0k(ζik) can be approximated by (5)f0kζik≈BζikTγk,k=1,···,q. Substituting this into model (4), we have (6)Yi≈ZiTβ+∑k=1qBζikTγk+εi.LetPi=(ZiT,WiT), WiT=(B(ζi1)T,···,B(ζiq)T), and θ=(θ1,···,θs)T=(βT,γ1T,···,γqT)T, s=p+qL. Due to the unknown ζik, we substitute the estimator of ζik by ζ^ik=Ψ(ξ^ik,λ^k), where ξ^ik and λ^k are the estimators of ξik and λk, respectively. Let P^i=(ZiT,W^iT), W^iT=(B(ζ^i1)T,···,B(ζ^iq)T), and define the following estimating function of θ(7)Uθ=∑i=1nP^iTYi-P^iθ.Motivated by the idea of Ueki [19], we propose the following smooth-threshold estimating equation: (8)Is-ΔUθ+Δθ=0,where s=p+qL, Is is the s-dimensional identity matrix, and Δ is the diagonal matrix, i.e., Δ=diag(Δ1,Δ2), where Δ1=diag{δ11,···,δ1p} is a p×p diagonal matrix, and Δ2=diag{δ21,···,δ21,···,δ2q,···,δ2q} is a qL×qL diagonal matrix. Note that δ1j=1 reduces to βj=0 for j=1,···,p, and δ2k=1 reduces to γk=0, that is, f0k(·)=0, for k=1,···,q. Therefore, (8) can yield a sparse solution. Unfortunately, we cannot directly obtain the estimator of θ by solving the smooth-threshold estimating equation (8). The reason is that (8) involves the unknown parameters δ1j and δ2k, which need to be chosen using some data-driven criteria.For the choice ofδ1j and δ2k, Ueki [19] suggested that δ1j and δ2k may be determined by the data, δ^1j=min(1,ν/β^j01+τ) and δ^2k=min(1,ν/γ^k(0)1+τ), with initial estimators β^j(0) and γ^k(0), respectively. The initial estimators β^j(0) and γ^k(0) are the solutions of the estimating equation U(θ)=0. Note that this choice involves two tuning parameters (ν,τ). In Section 3, we will propose a BIC-type criterion to select the tuning parameters. Replacing Δ in (8) by Δ^ with diagonal elements δ^1j and δ^2k, the smooth-threshold estimating equation (8) becomes (9)Is-Δ^Uθ+Δ^θ=0.The solution of (9) denoted by θ^=(θ^1,⋯,θ^p,θ^p+1,⋯,θ^s)T=(β^T,γ^1T,···,γ^qT)T is called the SEE estimator. Thus, β^ is the estimator of β, and Bk(ζ^k)Tγ^k is the estimator of f0k(ζk). For the convenience of the notations, we define A1={j:δ^1j≠1} as the set of indices of nonzero regression coefficient estimators, and define A2={k:δ^2k≠1} as the set of nonzero function estimators. A1c={j:δ^1j=1} and A2c={k:δ^2k=1} are the sets of zero regression coefficient estimators and zero function estimators, respectively. We further denote A=A1∪A2, and Ac=A1c∪A2c. The SEE estimator is defined in the following form:(10)θ^A=∑i=1nP^i,ATP^i,A-G^A-1∑i=1nP^i,ATYi,θ^Ac=0, Where P^i,A=(Zi,A1T,W^i,A2T), G^A=(IA-Δ^A)-1Δ^A, Δ^A=diag{Δ^1,A1,Δ^2,A2}, and IA is the identity matrix with the same dimension of Δ^A.


## 2.2. Asymptotic Properties

We first introduce some notations. Letf0(·) and β0 denote the true values of f(·) and β, respectively. γ0 is the spline coefficient vector from the spline approximation to f0(·). Denote θ0=(β0T,γ0T)T, A10={j:β0j≠0},A10c={j:β0j=0}, A20={k:f0k(·)≠0}, and A20c={k:f0k(·)=0}. We also denote A0=A10∪A20 and A0c=A10c∪A20c. The norm of f(·) is defined as f(x)=(∫01f2(x)dx)1/2.For convenience and simplicity, letC denote a positive constant that may be different at each appearance throughout this paper. We list some regularity conditions that are used in this paper.(C1) The transformation functionΨ(ξ,λk) is differentiable at ξ and λk, and satisfies the fact that ∂Ψ(ξ,λk)/∂ξ≤Cλkc0 and ∂Ψ(ξ,λk)/∂λk≤Cλkc0ξ for some positive constants C and negative constants c0.(C2) The second derivativeX(2)(t) is continuous on [0,1] with probability 1 and ∫E[X(k)(t)]4dt<∞ with probability 1, for k=0,2.(C3) The spline regression parameter is identified; that is, there is a uniqueθ0=(β0T,γ01T,⋯,γ0qT)T∈S, where the parameter space S is compact.(C4) The inner knots{ci,i=1,⋯,K} satisfy (11)max1≤i≤Khi+1-hi=oK-1,hmin1≤i≤Khi≤C0, where hi=ci-ci-1, h=max1≤i≤Khi, c0=0, and cK+1=1.(C5)f0k(·) is the rth continuously differentiable on (0,1), where r≥2.The following theorem gives the consistency of our proposed estimators.Theorem 1.
Suppose that regularity conditions (C1)-(C5) hold. For any positiveν and τ, nr/(2r+1)ν→0 and nr(1+τ)/(2r+1)ν→∞, as n→∞, and if K=O(n1/(2r+1)), one has (12)f^0kζ^k-f0kζk=Opn-r/2r+1,k=1,⋯,q.In the following theorem, we will show that such consistent estimators enjoy the sparsity property.Theorem 2.
Under the regularity conditions of Theorem1, one has(i)
P(A1c=A10c)→1,(ii)
P(A2c=A20c)→1.Next, we will show that the estimators of the nonzero coefficients for the parametric components have the same asymptotic distribution as that based on the correct model.Theorem 3.
Suppose that the regularity conditions of Theorem1 hold; as n→∞, one has, (13)Dβ^A1-1/2β^A1-β0A10⟶LN0,IA1,where D(β^A1) is defined in (A.23) in the Appendix, A1 is the number of A1, IA1 is the identity matrix, and “→L” represents the convergence in distribution.


## 3. Issues in Practical Implementation

### 3.1. Computational Algorithm

Since the transformed FPC scores{ζi} cannot be observed, we first need to estimate the FPC scores before the estimation and selection of f0k(·). In what follows, we propose the algorithm to implement the estimation procedure.Step 1.
Apply the functional principal component analysis (Zhu, Yao, and Zhang, 2014) to estimate the FPC scores{ξi1,···,ξiq} of Xi(t) denoted by {ξ^i1,···,ξ^iq}. Then we can obtain the transformed variables ζ^ik=Ψ(ξ^ik,λ^k), where λ^k is the estimated eigenvalue, and q is chosen to explain nearly 100% of the total variation.Step 2.
Calculate the initial estimateθ^0 of θ by solving the estimating equation U(θ)=0.Step 3.
Choose the tuning parametersν and τ by the BIC-type criterion in the next subsection.Step 4.
Solve the smooth-threshold estimating equations (9) and update the estimator of θ as follows:(14)θ^A=∑i=1nP^i,ATP^i,A-G^A-1∑i=1nP^i,ATYi,θ^Ac=0.


### 3.2. Selection of the Tuning Parameters

Following Ueki [19], we minimize the following BIC-type criterion to choose the tuning parameters (ν,τ). (15)BICν,τ=∑i=1nYi-P^iθ^ν,τTYi-P^iθ^ν,τ+log⁡n·dfν,τ,where θ^ν,τ is the SEE estimator for given (ν,τ), and dfν,τ is the number of nonzero coefficients θ^=(β^T,γ^1,…,γ^q)T.


## 3.1. Computational Algorithm

Since the transformed FPC scores{ζi} cannot be observed, we first need to estimate the FPC scores before the estimation and selection of f0k(·). In what follows, we propose the algorithm to implement the estimation procedure.Step 1.
Apply the functional principal component analysis (Zhu, Yao, and Zhang, 2014) to estimate the FPC scores{ξi1,···,ξiq} of Xi(t) denoted by {ξ^i1,···,ξ^iq}. Then we can obtain the transformed variables ζ^ik=Ψ(ξ^ik,λ^k), where λ^k is the estimated eigenvalue, and q is chosen to explain nearly 100% of the total variation.Step 2.
Calculate the initial estimateθ^0 of θ by solving the estimating equation U(θ)=0.Step 3.
Choose the tuning parametersν and τ by the BIC-type criterion in the next subsection.Step 4.
Solve the smooth-threshold estimating equations (9) and update the estimator of θ as follows:(14)θ^A=∑i=1nP^i,ATP^i,A-G^A-1∑i=1nP^i,ATYi,θ^Ac=0.


## 3.2. Selection of the Tuning Parameters

Following Ueki [19], we minimize the following BIC-type criterion to choose the tuning parameters (ν,τ). (15)BICν,τ=∑i=1nYi-P^iθ^ν,τTYi-P^iθ^ν,τ+log⁡n·dfν,τ,where θ^ν,τ is the SEE estimator for given (ν,τ), and dfν,τ is the number of nonzero coefficients θ^=(β^T,γ^1,…,γ^q)T.


## 4. Application to Real Data

We demonstrate the effectiveness of the proposed method by an application for the Tecator data set. The Tecator data set contains 215 samples; each sample contains the finely chopped pure meat with different moisture, fat, and protein contents, which are measured in percentages and are determined by the analytical chemistry. The functional covariate byX(t) for each food sample consists of a 100-channel spectrum of absorbances recorded on a Tecator Infratec Food and Feed Analyzer working in the wavelength range 850-1050 nm by the near-infrared transmission (NIT) principle. In this analysis, Y is the fat content, {ξk} are functional principal component (FPC) scores of X(t), and ζk=Ψ(ξk,λk), we take the protein and the moisture content by Z1 and Z2, respectively. In order to predict the fat content of a meat sample, many models and algorithms are proposed to fit the data; see, for example, Aneiros-Pérez and Vieu [26]. In this paper, to fit the data, we consider the following partially linear functional additive model: (16)Yi=β0+β1Zi1+β2Zi2+∑k=1qf0kζik+εi.To compare the performance of different models, the sample is divided into two data sets: the training sample I1={1,⋯,165} is used to obtain the estimators of the parameter and the nonparametric function, and the testing sample I2={166,⋯,215} is used to verify the quality of prediction by the following mean squared errors of prediction (PMSE), (17)PMSE=150∑j∈I2Yj-Y^j2varI2Y. In this section, we consider the variable selection problem for model (16). We apply the proposed smooth-threshold estimating equation approach to eliminate the irrelevant parameters in the model, while estimating the nonzero regression coefficients. Steps are as follows.Step 1.
Estimate the FPC scores{ξi1,···,ξiq} of Xi(t) denoted by {ξ^i1,···,ξ^iq}. Then, we can obtain the transformed variables ζ^ik=Ψ(ξ^ik,λ^k), where λ^k is the estimated eigenvalue, and q=4 is chosen to explain nearly 100% of the total variation.Step 2.
Calculate the initial estimateθ^0 of θ by solving the estimating equation U(θ)=0.Step 3.
Choose the tuning parametersν and τ by the BIC-type criterion.Step 4.
Solve the smooth-threshold estimating equations (9).The estimated values of theβ0, β1, and β2 are β^0=95.2503, β^1=-1.0223, and β^2=-0.7093, and f02(·) is estimated to be 0. The SEE method selects two parametric components Z1, Z2 and three nonparametric components: f01(·), f03(·), and f04(·). Hence, in conclusion, we would say that the strong linear relationship between the fat content and the protein and moisture contents, and f01(·), f03(·), and f04(·), are important functional nonparametric part. Our approach can be easily implemented and reduces the burden of computation.Aneiros-Pérez and Vieu [26] proposed several additive semifunctional models, Y=μ+m1(Z1)+m2(X)+ε, Y=μ+m3(Z2)+m4(X)+ε, and Y=μ+m5(Z1)+m6(Z2)+m7(X)+ε, where mi(·) are nonparametric function for i=1,⋯,7, and μ=E(Y). To assess the performance of the proposed PLFAM and SEE method, we compare with these additive semifunctional models and the SCAD-penalized method and the Lasso-penalized method. The PMSE results are reported in Table 1. From Table 1 we can see the following.Table 1
The prediction mean squared error (PMSE) for the different models.
Models PMSE Y=μ+m1(Z1)+m2(X)+ε, 0.0242 Y=μ+m3(Z2)+m4(X)+ε, 0.0368 Y=μ+m5(Z1)+m6(Z2)+m7(X)+ε, 0.0395 Y = β 0 + β 1 Z 1 + β 2 Z 2 + ∑ k = 1 4 ‍ f 0 k ( ζ k ) + ε, (SEE Method) 0.0153 Y = β 0 + β 1 Z 1 + β 2 Z 2 + ∑ k = 1 4 ‍ f 0 k ( ζ k ) + ε, (SCAD Method) 0.0160 Y = β 0 + β 1 Z 1 + β 2 Z 2 + ∑ k = 1 4 ‍ f 0 k ( ζ k ) + ε, (LASSO Method) 0.0234( 1 ) Real explanatory variables Z1 and Z2 can be used to improve the accuracy of the prediction; it is consistent with the conclusion of Aneiros-Pérez and Vieu [26].( 2 ) Compared with these additive semifunctional models, the PLFAM has a much smaller PMSE. The proposed PLFAM performs better than other models.( 3 ) Compared to the SCAD-penalized method and the Lasso-penalized method, the PMSE of PLFAM (SEE method) is smaller than all of the results of them.These conclusions confirm the validity of the proposed PLFAM and SEE method.


## 5. Concluding Remarks

The article develops a SEE procedure for automatic variable selection in the partially linear functional additive model. The proposed procedure can identify nonzero regression coefficients significant variables from the parametric components and the nonparametric components simultaneously, and it automatically eliminates the irrelevant parameters by setting them as zero and simultaneously estimates the nonzero regression coefficients. It is noteworthy that the proposed procedure avoids the convex optimization problem, and the resulting estimator enjoys the oracle property. The application to the Tecator data set confirms the validity of the proposed PLFAM model and the SEE method.



---

*Source: 5683539-2018-08-16.xml*