The exact contents of citations.db main_text.text for this paper — one flattened LaTeX string, title through conclusion, appendix excluded, unmodified except for removing email addresses. This is what our citation measures are computed over.
100,906 characters
A One-Covariate-at-a-Time Multiple Testing Approach to Variable Selection in Additive Models
\title{A One-Covariate-at-a-Time Multiple Testing Approach to Variable Selection in
Additive Models \thanks{We would like to thank the editor, Esfandiar Maasoumi,
an associate editor, and two referees for their helpful comments, which have
helped improve the paper substantially. Su gratefully acknowledges the support
from NSFC under the grant No. 72133002. Zhang gratefully acknowledges the
financial support from NSFC (Projects No.71973141 and No.71873033). The
computation code used in this paper is available upon request or can be found
on the second author's personal website. E-mail: [email removed] (L.
Su), [email removed] (T. Yang), [email removed] (Y. Zhang), and
[email removed] (Q. Zhou).}}
\author{Liangjun Su$^{a},$ Thomas Tao Yang$^{b}$, Yonghui Zhang$^{c},$ and Qiankun
Zhou$^{d}$\\$^{a}$ School of Economics and Management, Tsinghua University, China \\$^{b}$ Research School of Economics, Australian National University, Australia\\$^{c}$ School of Economics, Renmin University of China, China\\$^{d}$ Department of Economics, Louisiana State University, USA}
\maketitle
\begin{abstract}
This paper proposes a One-Covariate-at-a-time Multiple Testing (OCMT) approach
to choose significant variables in high-dimensional nonparametric additive
regression models. Similarly to \cite{Chudik_K_P}, we consider the statistical
significance of individual nonparametric additive components one at a time and
take into account the multiple testing nature of the problem. Both one-stage
and multiple-stage procedures are considered. The former works well in terms
of the true positive rate only if the net effects of all signals are strong
enough; the latter helps to pick up hidden signals that have weak net effects.
Simulations demonstrate the good finite sample performance of the proposed
procedures. As an empirical illustration, we use the OCMT procedure to a
dataset extracted from the Longitudinal Survey on Rural Urban Migration in
China. We find that our procedure works well in terms of out-of-sample
forecast root mean square errors, compared with competing methods such as
adaptive group Lasso (AGLASSO).
\vspace{2mm}
\begin{description}
\item[Keywords:] One covariate at a time, multiple testing, model selection,
high dimensionality, nonparametric, additive model.
\item[\emph{JEL classification:}] C12, C14, C21, C52.
\end{description}
\end{abstract}
\section{Introduction}
Variable selection has been playing a pivotal role in econometrics and
statistics for statistical learning and scientific discoveries. Notable early
contributions include \cite{Akaike1973}, \cite{Akaike1974}, and
\cite{Schwartz1978}. These authors suggested a unified approach to model
selection, viz., choosing a parameter vector by minimizing the conventional
criterion function plus an $L_{0}$ penalty to penalize the model size.
However, these approaches are not feasible in high-dimensional settings. The
seminal work of \cite{Tibshirani1996}\ addressed this important issue by
substituting the $L_{0}$ penalty with an $L_{1}$ penalty, and it has sparked
extensive studies in both statistics and econometrics. Important contributions
in the statistics literature include \cite{FanLi2001}, \cite{ZhouHastie2005},
\cite{Zou2006}, \cite{FanLv2008}, \cite{ZouLi2008}, \cite{Zhang2010},
\cite{FanFengSong}, \cite{FanLv2013}, and \cite{FanTang2013}. Important
contributions in the econometric literature include \cite{BelloniEtal2012},
\cite{BelloniChernozhukov2013}, \cite{BelloniChernozhukovHansen2014},
\cite{BelloniEtal2017},\ \cite{ChernozhukovEtal2018}, and \cite{Chudik_K_P}.
For a comprehensive review, see \cite{FanLiZhangZou2020}.
In this paper, we propose a multiple testing approach to variable selection
for high-dimensional nonparametric additive models. Recently \cite{Chudik_K_P}
(CKP hereafter) have proposed a One-Covariate-at-a-time Multiple Testing
(OCMT) approach for linear regression models. CKP suggest regressing the
dependent variable on each independent variable separately, retaining only
those variables that exhibit a high correlation with the dependent variable.
This strategy is often referred to as the \textquotedblleft\textit{screening
approach}\textquotedblright. CKP's main contributions are twofold. First, they
propose a criterion for variable selection by controlling the probability of
choosing all the signals (or with some pseudo-signals) in the model. Second,
unlike the usual screening approach that may miss some important
\textquotedblleft hidden\textquotedblright\ signals whose net effects on the
dependent variable are small, the CKP's OCMT procedure is able to pick up
hidden signals with very high probability. The OCMT procedure has been applied
in various applications; see, e.g., \cite{Kozbur2020}, \cite{Chudik_P_S}, and
\cite{AhmendPesaran2022}. In particular, \cite{Kozbur2020} considers a
testing-based forward model selection (TBFMS) procedure in linear regression
models that inductively selects covariates to add predictive power into a
working statistical model. But this latter paper mainly focuses on the error
bound and shows that the proposed procedure is able to achieve estimation
rates matching those of Least Absolute Shrinkage and Selection Operator
(Lasso) and post-Lasso. Furthermore, \cite{Sharifvaghefi2023} extends OCMT to
cases with many highly correlated covariates and allows the number of pseudo
signals to grow at the same rate as the sample size.
Our paper contributes to the literature by extending the CKP's\ OCMT approach
from parametric models to nonparametric additive models. Like CKP, we estimate
the net effect of each variable on the dependent variable one by one, possibly
with some preselected variables. The selected variables are those whose net
effects exceed some threshold value, specified to ensure the probability of
selecting all the signals is very high. The statistics constructed in this
paper are much more complicated than the $t$-statistics in CKP and might not
even exhibit a well-behaved limiting distribution. In addition to
investigating a different model, our paper differs from that of CKP in some
other important aspects. First, we generalize the definition of hidden signals
(in Table 2 in Section \ref{SEC:non_1st}). Technical details are updated
accordingly. Second, CKP chose tuning parameters similar to a Bonferroni
correction of the cumulative distribution function of the standard normal.
This choice implicitly requires a certain degree of approximation of the
standard normal to the $t$-statistic distribution in their study. Instead, we
select the tuning parameters using the classic Bayesian information criterion
(BIC). Third, we add an adaptive group-Lasso-based post-OCMT step to eliminate
pseudo-signals that cannot be eliminated with very high probability in CKP.
This step adds very little computation burden because the dimension is reduced
dramatically before the last-step estimation. This additional step is not
needed in theory, but it aims at eliminating the pseudo-signals and thereby
enhances the out-of-sample forecasting performance in practice.
One competing method to ours in the literature is the adaptive group Lasso
(AGLASSO) proposed by \cite{HuangEtal2010}. The AGLASSO adds some adaptive
penalty term to the usual least squares loss function in the spirit of
\cite{Zou2006}. Even though we also use AGLASSO to eliminate the
pseudo-signals after the OCMT procedure, our approach allows much faster
computation and provides more reliable estimates than their approach.
We consider various setups in the simulation studies and compare the above
post-OCMT AGLASSO procedure with that based on the OCMT\ alone or the AGLASSO
procedure of \cite{HuangEtal2010} alone. We find that the former one generally
outperforms the latter two significantly. We apply our method on a dataset
from the Longitudinal Survey on Rural Urban Migration in China (RUMiC) and the
empirical results also demonstrate the excellent performance of our procedure
in finite samples.
The remainder of the paper is structured as follows. In the next section, we
illustrate our approach through a single-stage procedure that is silent to
hidden signals. In Section \ref{SEC:multi}, we present the more powerful
multiple-stage procedure. We investigate the finite sample properties of our
procedure through Monte Carlo experiments in Section \ref{SEC:MC} and an
empirical application in Section \ref{SEC:application}. We conclude the paper
in Section \ref{SEC:conclusion}. The proofs of all propositions and theorems
in the paper are relegated to Appendix \ref{APP:real_main_proof}. The online
supplement contains some additional technical materials that include the
proofs of the technical lemmas in Appendix \ref{APP:mainproof} and some
additional results in the simulation and application. To facilitate reading,
We present our procedure in detail for practitioners in Appendix
\ref{SEC:procedure} and the idea of the proofs in Appendices
\ref{APP:tech_one_stage} and \ref{SEC:multi_tech}.
\textit{Notation.} For a generic real matrix $\boldsymbol{A=}\left\{
a_{ij}\right\} \boldsymbol{,}$ let $\left\Vert \boldsymbol{A}\right\Vert
=\left[ \lambda_{\max}\left( \boldsymbol{A}^{\prime}\boldsymbol{A}\right)
\right] ^{1/2}$ denote the spectral norm and $\left\Vert \boldsymbol{A}
\right\Vert _{\infty}=\max_{ij}\left\vert a_{ij}\right\vert .$ When
$\boldsymbol{A}$ is symmetric, $\lambda_{\max}\left( \boldsymbol{A}\right)
\ $and $\lambda_{\min}\left( \boldsymbol{A}\right) \ $denote its maximum and
minimum eigenvalues, respectively. For vector $\boldsymbol{x},$\ $\left\Vert
\boldsymbol{x}\right\Vert $ denotes its Euclidean norm. For the deterministic
series $\left\{ a_{n},b_{n}\right\} _{n=1}^{\infty}$, we denote
$a_{n}\propto b_{n}$ if $0<C_{1}\leq\lim\inf_{n\rightarrow\infty}\left\vert
a_{n}/b_{n}\right\vert \leq\lim\sup_{n\rightarrow\infty}\left\vert a_{n}
/b_{n}\right\vert \leq C_{2}<\infty$ for some constants $C_{1}$ and $C_{2},$
$a_{n}\lesssim b_{n}$ if $\lim\sup_{n\rightarrow\infty}\left\vert a_{n}
/b_{n}\right\vert \leq C<\infty$ for some $C$ that does not depend on $n,$
$a_{n}\gtrsim b_{n}$ if $b_{n}\lesssim a_{n},$ $a_{n}\ll b_{n}$ if
$a_{n}=o\left( b_{n}\right) ,$ and $a_{n}\gg b_{n}$ if $b_{n}\ll a_{n}.$
$\mathcal{A}^{c}$ denotes the complement of the set $\mathcal{A}$.
$\overset{P}{\rightarrow}$ denotes convergence in probability. $C$ and $M$
denote some positive constants that may vary from line to line.
\section{The Model and One-Stage Procedure}
In this section we give the model and definitions of various versions of
signals and noises. Then we will provide the one-stage procedure for variable
selection, present the basic assumptions and study the asymptotic properties
of our one-stage procedure.
\subsection{The Model}
Recently, CKP proposed a powerful multiple-stage procedure for linear models.
Our paper aims to generalize the results of CKP to the nonparametric additive
models. The notations and technical details in CKP are already quite tedious,
and they would be even more so in this paper. To facilitate the exposition, we
start with a simple case where there are no pre-determined variables, and we
conduct only one-stage multiple testing. We will present the more powerful
multiple-stage procedure and show its validity in Section \ref{SEC:multi}.
Suppose the model is
\begin{equation}
Y=f^{\ast}\left( X_{1},X_{2},\ldots,X_{p^{\ast}}\right) +\varepsilon,
\label{EQ:model}
\end{equation}
where $Y$ is the dependent variable, $X_{1},$ $X_{2},\ldots,$ and $X_{p^{\ast
}}$ are random independent variables, $\varepsilon$ is an unobserved error
term, and $f^{\ast}$ is an unknown smooth function. Even though we have only
$p^{\ast}$ signal variables, namely $X_{1},$ $X_{2},\ldots,$ $X_{p^{\ast}}$,
that should be included into the regression model in (\ref{EQ:model}), we do
not know this truth before the data reveal the fact. The realistic situation
is that the $p^{\ast}$ signals are contained in a set $\mathcal{S}
_{n}=\left\{ X_{j},j=1,2,\ldots,p_{n}\right\} ,$ where $p_{n}$ can be much
larger than the sample size $n$ and $p_{n}\propto n^{B_{p}}$ for some
$B_{p}>0$. $\mathcal{S}_{n}$ is the set for all candidate variables, and we
refer to it as the \textit{active }set. We assume that $E\left(
\varepsilon|X_{1},X_{2},\ldots,X_{p^{\ast}},X_{p^{\ast}+1},\ldots,X_{p_{n}
}\right) =0.$ The target of the model selection is to pick up those signals
among $\mathcal{S}_{n}.$
To avoid the curse of dimensionality in nonparametric estimation, we impose
the additive structure on $f^{\ast}$, that is,
\begin{equation}
f^{\ast}\left( X_{1},X_{2},\ldots,X_{p^{\ast}}\right) =\sum_{j=1}^{p^{\ast}
}f_{j}^{\ast}\left( X_{j}\right) . \label{EQ:model1a}
\end{equation}
We allow the additive components\textbf{ }$f_{l}^{\ast}\left( X_{l}\right) $
to change with $n$ but we suppress the dependence of $f_{l}^{\ast}\left(
X_{l}\right) $ on $n$ for notional convenience. Obviously, the individual
functions $f_{j}^{\ast},$ $j=1,\ldots,p^{\ast},$ cannot be identified without
certain suitable normalizations. We impose the normalization by assuming that
$E[f_{j}^{\ast}\left( X_{j}\right) ]=0$ for each $j$ and rewrite the model
as
\begin{equation}
Y=\mu+\sum_{j=1}^{p^{\ast}}f_{j}^{\ast}\left( X_{j}\right) +\varepsilon,
\label{EQ:model1b}
\end{equation}
where $\mu=E\left( Y\right) .$ Since $\mu$ can be estimated by the sample
mean of the $Y$ variable at the usual $\sqrt{n}$-rate and this estimation does
not affect the estimation of the nonparametric additive component, for the
simplicity of presentation, we assume $\mu=0$ below.
Since we will focus on the $B$-spline-based nonparametric theory, it is
standard to assume compact support for each regressor (see, e.g.,
\cite{Horowitz_Mammen2014}, \cite{Chen_handbook} and \cite{HuangEtal2010}).
Without loss of generality, we assume that the support for $X_{j}$ is
$[0,1\mathcal{]}$ for all $j$. Let $\alpha_{1}$ be a non-negative integer,
$\alpha_{2}\in(0,1],$ and $d=\alpha_{1}+\alpha_{2}.$\ Denote the class of all
$\alpha_{1}$ times continuously differentiable real-valued functions on
$[0,1]$ by $C^{\alpha_{1}}\left( [0,1]\right) .$\ Define $d$-th smooth
real-valued functions on $[0,1]$ as
\begin{equation}
\Lambda^{d}\left( [0,1]\right) =\left\{ h\in C^{\alpha_{1}}\left(
[0,1]\right) :\left\vert h^{\left( \alpha_{1}\right) }\left( t_{1}\right)
-h^{\left( \alpha_{1}\right) }\left( t_{2}\right) \right\vert \leq
C\left\vert t_{1}-t_{2}\right\vert ^{\alpha_{2}}\right\} . \label{EQ:model1c}
\end{equation}
For notational simplicity, we will restrict our attention to the case where
$f_{j}$'s are smooth enough and belong to $\Lambda^{d}\left( [0,1]\right) $
with $d>1$. The case of different smoothness parameters only complicates the
notation but does not bring in any new insight.
\subsection{Signals, Hidden Signals, and Noises}
The idea of one-stage procedure is that we estimate the impact of $X_{l}$ on
$Y,$ $l=1,2,...,p_{n},$ one by one. So we run $p_{n}$ estimations in total and
will keep those variables that are significant enough. Since in each
regression we only have one covariate, we are not estimating $f_{l}^{\ast}$
when the explanatory variable is $X_{l}$. Instead, we are estimating the
conditional expectation $f_{l}\left( X_{l}\right) \equiv E\left(
Y|X_{l}\right) .$ We define the \textit{net} \textit{impact} (or the net
effect) of $X_{l}$ on $Y$ as
\[
\theta_{l}\equiv\left\{ E\left[ f_{l}\left( X_{l}\right) ^{2}\right]
\right\} ^{1/2}=\left\{ E\left[ \left( \sum_{j=1}^{p^{\ast}}\sigma
_{lj}\right) ^{2}\right] \right\} ^{1/2},
\]
where $\sigma_{lj}=E[f_{j}^{\ast}\left( X_{j}\right) |X_{l}].$ Since we
allow\textbf{ }$f_{l}^{\ast}\left( X_{l}\right) $ to change with $n,$
$\theta_{l}$ might change with $n$ as well. But we suppress its dependence on
$n$ for notational convenience. Here, $\sigma_{lj}$ plays the role of the
scaled covariance between $X_{j}$ and $X_{l}$ for the linear model considered
by CKP. To better understand the connection between the net impact defined
here and in CKP, we refer the readers to the results in equations
(\ref{EQ:fnl&Ul}) and (\ref{EQ:fnl_appro}), where we approximate $f_{l}\left(
X_{l}\right) $ with certain linear functions $f_{nl}\left( X_{l}\right) .$
Ignoring the bias from the approximation, one can see the connection more clearly.
Obviously, $\theta_{l}$ can be 0 or close to 0 for signals, and $\theta_{l}$
can be nonzero or large for non-signals. As in CKP, we also have four
possibilities as tabulated in Table \ref{Table1}. Cases (I) and (IV) in Table
1 are desirable cases. Case (III) happens when some non-signals are not
independent of the signals.\ The hidden signals defined in Case (II) are rare
in the linear case, and it is also rare in the nonparametric case. We
generalize the definition of hidden signals to Table 2 in the next section
where $\theta_{l}$ is non-zero but small relative to the sample size.
\begin{table}[ptb]
\caption{The original definitions of signals and noises from CKP}
\label{Table1}
\centering{}\centering{ }
\begin{tabular}
[c]{l|cc}\hline\hline
& $\theta_{l}\neq0$ & $\theta_{l}=0$\\\hline
$\left\{ E[f_{l}^{\ast}\left( X_{l}\right) ^{2}]\right\} ^{1/2}\neq0$ &
\multicolumn{1}{|l}{(I) Signals with nonzero net effect} &
\multicolumn{1}{l}{(II) Hidden signals}\\
$\left\{ E[f_{l}^{\ast}\left( X_{l}\right) ^{2}]\right\} ^{1/2}=0$ &
\multicolumn{1}{|l}{(III) Pseudo-signals} & \multicolumn{1}{l}{(IV) Noise
variables}\\\hline\hline
\end{tabular}
\end{table}
As we shall see, our one-stage procedure is silent on picking up hidden
signals and eliminating pseudo-signals. For the hidden signals, we will
propose a multiple-stage procedure in Section \ref{SEC:multi}\ that can
effectively pick them up. To eliminate the pseudo-signals, as a post-procedure
we propose to re-estimate the model using adaptive group Lasso in Section
\ref{SEC:Pseudo}.
We assume that there are $p^{\ast\ast}$ pseudo-signals. Without loss of
generality, we denote them to be
\[
\left\{ X_{p^{\ast}+1},X_{p^{\ast}+2},\ldots,X_{p^{\ast}+p^{\ast\ast}
}\right\} .
\]
Below we focus on the one-stage procedure in this section and postpone the
multi-stage case to the next section.
\subsection{The Test Statistic and One-Stage Procedure\label{SEC:1stage}}
Suppose that we have $n$ observations $\left\{ (y_{i},x_{1i},...,x_{p_{n}
i})\right\} _{i=1}^{n}$ that are drawn from the distribution of
$(Y,X_{1},...,X_{p_{n}}).$ The data are given in a $n\times\left(
p_{n}+1\right) $ matrix
\[
\left( \boldsymbol{y},\boldsymbol{x}_{1},\boldsymbol{x}_{2},\ldots
,\boldsymbol{x}_{p_{n}}\right)
\]
where $\boldsymbol{y}=\left( y_{1},y_{2},\ldots,y_{n}\right) ^{\prime}$ and
$\boldsymbol{x}_{l}=\left( x_{l1},x_{l2},\ldots,x_{ln}\right) ^{\prime}$ for
$l=1,...,p_{n}$. We propose to choose finite order (e.g., cubic) B-spline
basis functions $\left\{ \psi_{l}\left( x\right) \right\} _{l=1}^{m_{n}}
$\ on $\left[ 0,1\right] $\ to approximate the unknown functions $f_{l}$'s.
B-splines are piecewise-defined polynomial functions that can be used to
construct curves and surfaces in numerical analysis. They offer a flexible way
to model and control the shape of these curves and surfaces. A B-spline of
order $n$ is a piecewise-defined polynomial function of degree $n-1$.
B-splines of order one are piecewise constant functions, B-splines of order
two are piecewise linear functions, B-splines of order three are piecewise
quadratic functions, and so on. The points at which different polynomial
pieces connect are called knots, and the knot vector specifies where these
knots are. For the detailed definition and properties of B-spline bases, see
\cite{Stone} and \cite{deBoor}. We list some properties of the B-spline basis
functions in Lemma \ref{LE:rank}. Other popular basis functions include
polynomials, trigonometric polynomials, splines, and orthogonal wavelets. We
refer the readers to \cite{Chen_handbook} for a nice review about sieve estimation.
Since we normalize $E[f_{j}^{\ast}\left( X_{j}\right) ]=0,$ we similarly
normalize the basis as
\[
\phi_{jl}\left( x\right) =\psi_{j}\left( x\right) -n^{-1}\sum_{i=1}
^{n}\psi_{j}\left( x_{li}\right) .
\]
This is a standard practice; see, e.g., \cite{HuangEtal2010}. For notational
simplicity, we will write $\phi_{j}\left( x\right) $ for $\phi_{jl}\left(
x\right) $. Let $P^{m_{n}}\left( x\right) =\left[ \phi_{1}\left(
x\right) ,\phi_{2}\left( x\right) ,\ldots,\phi_{m_{n}}\left( x\right)
\right] ^{\prime},$ an $m_{n}\times1$\ vector. Define
\begin{equation}
\boldsymbol{\beta}_{l}=\left\{ E\left[ P^{m_{n}}\left( X_{l}\right)
P^{m_{n}}\left( X_{l}\right) ^{\prime}\right] \right\} ^{-1}E\left[
P^{m_{n}}\left( X_{l}\right) Y\right] \text{ and }U_{l}=Y-P^{m_{n}}\left(
X_{l}\right) ^{\prime}\boldsymbol{\beta}_{l}, \label{EQ:betanl}
\end{equation}
which are the population coefficient and the error term in the regression of
$Y$ on $P^{m_{n}}\left( X_{l}\right) $. Note that we suppress the dependence
of $\boldsymbol{\beta}_{l}$ on the sample size $n$.
For the one-stage procedure, we conduct the regression of $Y$ on $P^{m_{n}
}\left( X_{l}\right) ,$ $l=1,2,\ldots,p_{n},$ one by one. Let $\mathbb{X}
_{li}=P^{m_{n}}\left( x_{li}\right) $ be the approximating function basis at
the $i$th observation for $X_{l}$.$\ $Let $\mathbb{X}_{l}=(\mathbb{X}
_{l1},\mathbb{X}_{l2},\ldots,\mathbb{X}_{ln})^{\prime}$ be the $n\times m_{n}$
\textquotedblleft design\textquotedblright\ matrix for $X_{l}.$ For $X_{l},$
we regress $\boldsymbol{y}$ on $\mathbb{X}_{l}\ $to obtain
\[
\boldsymbol{\hat{\beta}}_{l}=\left( \mathbb{X}_{l}^{\prime}\mathbb{X}
_{l}\right) ^{-1}\mathbb{X}_{l}^{\prime}\boldsymbol{y}.
\]
We construct the test statistic as\footnote{As a referee has noted, one can
define an alternative test statistic $\tilde{\mathcal{X}}_{l}=n\hat
{\boldsymbol{\beta}}^{\prime}\hat{\boldsymbol{\beta}}$, which also works under
certain rank conditions (see, e.g., Assumption \ref{A:full_rank2} below). We
opted for $\hat{\mathcal{X}}_{l}$ for two reasons. First, $\hat{\mathcal{X}
}_{l}$ resembles the usual chi-squared statistic under conditional
homoskedasticity. Because we do not want to model the conditional
heteroskedasticity of unknown form, we cannot take into account conditional
heteroskedasticity explicitly in constructing the test statistic
$\hat{\mathcal{X}}_{l}.$ Despite this, our asymptotic theory allows for
conditional heteroskedasticity in the error term. Second and more importantly,
$\hat{\mathcal{X}}_{l}$ is scale-free whereas $\tilde{\mathcal{X}}_{l}$ is
not. The latter makes it very challenging to choose the range to search the
constant $C$ in $\varsigma_{n}$ defined below.}
\begin{equation}
\mathcal{\hat{X}}_{l}=\boldsymbol{\hat{\beta}}_{l}^{\prime}\left( \hat
{\sigma}_{l}^{-2}\mathbb{X}_{l}^{\prime}\mathbb{X}_{l}\right)
\boldsymbol{\hat{\beta}}_{l}, \label{EQ:chil_definition}
\end{equation}
where $\hat{\sigma}_{l}^{2}=n^{-1}\sum_{i=1}^{n}\hat{u}_{li}^{2}$ and $\hat
{u}_{li}$ is the residual from the above regression. Then we define the
first-stage OCMT selection indictor as
\begin{equation}
\widehat{\mathcal{J}}_{l}=\mathbf{1}\left( \mathcal{\hat{X}}_{l}
>\varsigma_{n}\right) \text{ for }l=1,2,\ldots,p_{n}, \label{EQ:select1stage}
\end{equation}
where $\mathbf{1}\left( \cdot\right) $ is the usual indicator function and
$\varsigma_{n}$ is a threshold value.
For the linear model in CKP with one covariate at a time, $\mathcal{\hat{X}
}_{l}$ is asymptotically $\chi^{2}\left( 1\right) $ under conditional
homoskedasticity, and one can follow their lead to consider threshold values
for the associated $t$-statistics based on the adjusted normal critical
values. Nevertheless, such a result is not available in our framework due to
the divergent dimension of regressors in the sieve estimation. In addition,
the potential presence of conditional heteroskedasticity greatly complicates
our analysis too. What we really need is to show that $\mathcal{\hat{X}}_{l}$
behaves distinctly for signals and noises so that a suitable choice of the
threshold value $\varsigma_{n}$ can help us separate the signals from the
noises. For these reasons, we do not associate $\mathcal{\hat{X}}_{l}$ with
any asymptotic distribution. Instead, we will set $\varsigma_{n}\propto
\kappa_{n}\log\left( m_{n}\right) m_{n}$ for a positive series $\kappa_{n}$
that diverges to infinity slowly as in Assumption \ref{A:xi_n}. For more
details, see the remark on Assumption \ref{A:xi_n} in the next subsection.
\subsection{Basic Assumptions}
To study the asymptotic properties of the one-stage procedure, we impose the
following assumptions.
\begin{assumption}
\label{A:iid}$\left\{ y_{i},x_{1i},x_{2i},\ldots,x_{p_{n}i}\right\}
_{i=1}^{n}$ are independent and identically distributed (i.i.d.) across $i;$
$E(\varepsilon|X_{1},X_{2},\ldots,X_{p^{\ast}},X_{p^{\ast}+1},\ldots,$
$X_{p_{n}})$ $=0.$
\end{assumption}
\begin{assumption}
\label{A:p} $p^{\ast}\ $is a positive integer that does not vary with
$n.$\ $p^{\ast\ast}\lesssim n^{B_{p^{\ast\ast}}}\ $and $p_{n}\propto n^{B_{p}
}$ for some $B_{p}>B_{p^{\ast\ast}}\geq0.$
\end{assumption}
\begin{assumption}
\label{A:supp}The support for $X_{l}$ is $\left[ 0,1\right] ,$
$l=1,\ldots,p_{n}$. The density function for $X_{l}$ is bounded and bounded
away from $0.$
\end{assumption}
\begin{assumption}
\label{A:epsilon}$\Pr\left( \left\vert \varepsilon\right\vert >t\right) \leq
C_{1}\exp\left( -C_{2}t^{s}\right) $ holds for all $t>0$ and some $s,$
$C_{1},$ $C_{2}>0.$
\end{assumption}
\begin{assumption}
\label{A:fl}$f_{l}\left( \cdot\right) =E\left( Y|X_{l}=\cdot\right)
\in\Lambda^{d}\left( [0,1]\right) \ $with $d>1$ for $l=1,\ldots,p_{n}$.
\end{assumption}
\begin{assumption}
\label{A:tech}$\left\vert f_{j}^{\ast}\left( \cdot\right) \right\vert ,$
$j=1,...,p^{\ast},$ are uniformly bounded. $E\left( \varepsilon^{2}
|X_{1},X_{2},...,X_{p_{n}}\right) $ is uniformly bounded almost surely.
\end{assumption}
\begin{assumption}
\label{A:mn}$m_{n}\propto n^{B_{m}}$ with $1/\left( 1+2d\right) <B_{m}<1/3.$
\end{assumption}
\begin{assumption}
\label{A:xi_n}$\varsigma_{n}\propto\kappa_{n}\log\left( m_{n}\right)
m_{n}\ $for some $\kappa_{n}$ such that $\kappa_{n}>0,$ $\kappa_{n}
\rightarrow\infty$, and $\kappa_{n}=O\left[ \left( \log n\right)
^{\epsilon}\right] $ for some small positive $\epsilon$,\ as $n\rightarrow
\infty$.$\medskip$
\end{assumption}
Assumption \ref*{A:iid} imposes an i.i.d. condition on the observations and a
standard conditional moment restriction. The extension to weakly dependent
observations is possible but left for future research. The generalization to
independently non-identically distributed (i.n.i.d.) case is straightforward
because the main inequalities in Lemmas A.1 and A.2 allow for i.n.i.d.
observations. We keep using the i.i.d. assumption for notational convenience.
In Assumption \ref*{A:p}, we assume that the number of signals is fixed; the
number of pseudo-signals is allowed to increase as $n$ increases but at a
slower rate than that of the total candidate variables. It is also possible to
allow $p^{\ast}$ to diverge to infinity (see Section \ref{SEC:divergeP}).
Assumption \ref*{A:supp} restricts the support of $X$ to be $\left[
0,1\right] .$ This is a very common condition for nonparametric additive
models; see, e.g., \cite{Li2000} and \cite{Horowitz_Mammen2014}. Assumption
\ref*{A:epsilon} imposes some tail conditions $\varepsilon,$ which is also
assumed in CKP but weaker than the commonly used sub-exponential condition in
the variable selection literature and the one for the adaptive group Lasso in
\cite{HuangEtal2010}. The primary use of this condition is to derive
\textit{probability} bound for errors.
Assumption \ref*{A:fl} imposes fairly weak smooth condition on $f_{l}\left(
x\right) $, which is weaker than the commonly used condition $d\geq2$. Note
that $d>1$ is needed for Assumption \ref*{A:mn}. Assumption \ref*{A:tech} is a
technical assumption needed to simplify the proof.\ Specifically, the
boundedness of $f_{j}^{\ast}$ implies $U_{l}$ defined in equation
(\ref{EQ:fnl&Ul}) has the same tail behavior as $\varepsilon$. The bounded
conditional second moment of $\varepsilon$\ is to ensure some nice properties
of $U_{l}\phi_{j}\left( X_{l}\right) $ and this assumption is also common in
the sieve literature (see, e.g., \cite{Newey1997}).\ This assumption is mild
given that we assume all $X_{l}$'s have compact support and the tails of
$\varepsilon$ decay exponentially fast. Assumption \ref*{A:mn} imposes
conditions on $m_{n}.$ First, we need $B_{m}<1/3$ such that $nm_{n}
^{-3}\rightarrow\infty.$ The last condition is necessary for $\left\Vert
n^{-1}\mathbb{X}_{l}^{\prime}\mathbb{X}_{l}\right\Vert \propto m_{n}^{-1}$ to
hold with very high probability. To see why, note that Lemma \ref{LE:rank} in
Appendix \ref{APP:tech_lemmas} suggests that $\left\Vert E\left[ P^{m_{n}
}\left( X_{l}\right) P^{m_{n}}\left( X_{l}\right) ^{\prime}\right]
\right\Vert \propto m_{n}^{-1}.$ We need $m_{n}^{-1}\gg\left( m_{n}/n\right)
^{1/2},$ or equivalently, $nm_{n}^{-3}\rightarrow\infty,$ in order to ensure
that $n^{-1}\mathbb{X}_{l}^{\prime}\mathbb{X}_{l}$ is close to $E\left[
P^{m_{n}}\left( X_{l}\right) P^{m_{n}}\left( X_{l}\right) ^{\prime
}\right] $ with very high probability. Second, we need $B_{m}>1/\left(
1+2d\right) $ to ensure that the approximation bias is asymptotically
negligible in comparison with the asymptotic variance term: $m_{n}^{-d}
\ll\left( m_{n}/n\right) ^{1/2}$.
Assumption \ref*{A:xi_n} imposes conditions on the threshold value
$\varsigma_{n}$ that ensures the separability of the signals from noises. If
the true value of $\boldsymbol{\beta}_{l}$ is $\boldsymbol{0}$, $\mathcal{\hat
{X}}_{l}$ in equation (\ref{EQ:chil_definition}) is $O_{P}\left(
m_{n}\right) .$ In this case, to ensure $\widehat{\mathcal{J}}_{l}=0$ with
very high probability, we can take $\varsigma_{n}\propto\kappa_{n}\log\left(
m_{n}\right) m_{n}$ and lose some power up to $\kappa_{n}\log\left(
m_{n}\right) .$ Here, $\kappa_{n}$ can be any series diverging to infinity
slowly, e.g., $\left[ \log\left( m_{n}\right) \right] ^{\epsilon}$ for
some small $\epsilon>0.$ The loss of the power to some degree is inevitable
because of the nature of the multiple testing procedure when the number of
tests goes to infinity. In contrast, CKP choose their threshold by the
Bonferroni correction of the standard normal. We choose not to follow them
because of the following reasons. First, given the divergence of $m_{n},$
$\mathcal{\hat{X}}_{l}$ does not converge to a chi-square distribution
asymptotically even in the homoskedastic case so that we cannot use chi-square
distribution to approximate the finite sample distribution of $\mathcal{\hat
{X}}_{l}.$ Under conditional heteroskedasticity, $\mathcal{\hat{X}}_{l}$ does
not converge to a chi-square distribution even if $m_{n}$ is held fixed. So
our procedure does not rely on the chi-square approximation. Second, even if
we can do the approximation, the cumulative density function (CDF) of a
chi-square distribution is very complicated. We do not have a rate for the
inverse of its CDF evaluated at a certain rate (e.g., $n^{-C}$) like the case
of normal CDF. For these reasons, we do our selection based on the asymptotic
results. Specifically, we will take $\varsigma_{n}=C\kappa_{n}\log\left(
m_{n}\right) m_{n}$ for some $\kappa_{n}$ and choose the value of $C$ by the
classic BIC. The details can be found in Appendix \ref{SEC:procedure}. This
$\varsigma_{n}$ helps separate the signals from the noises with very high
probability, as demonstrated in the next section.\ Recall that in Assumption
\ref*{A:p} we assume that $p_{n}$ go to infinity at a polynomial rate of $n,$
same as $m_{n}.$ This ensures that $\log\left( m_{n}\right) \propto\log
p_{n}.$ Therefore, theoretically we only need to put $\log\left(
m_{n}\right) $ in $\varsigma_{n}$ to have a control of $p^{\ast\ast}$ and
$p_{n}$ for the TPR, FDR and FPR defined after Proposition \ref{TH:main1}
below. We postpone the discussion on the comparison of technical conditions
required for our procedure and those required for the AGLASSO to Section 3.2.
\subsection{The Asymptotic Properties of $\mathcal{\hat{X}}_{l}$ and the
One-Stage Procedure\label{SEC:non_1st}}
We present the first theoretical result in this paper. It derives the
probability bounds for the \textquotedblleft Type-I\textquotedblright\ and
\textquotedblleft Type-II\textquotedblright\ errors.
\begin{proposition}
\label{TH:main1}Suppose that $Y$ is given by equation (\ref{EQ:model}) and
Assumptions \ref*{A:iid} $-$ \ref*{A:xi_n} hold.
(i) If $\theta_{l}\lesssim\log\left( m_{n}\right) ^{1/2}\left(
m_{n}/n\right) ^{1/2},$ then for sufficiently large $n$ we have
\begin{align*}
\Pr\left( \mathcal{\hat{X}}_{l}\geq\varsigma_{n}\right) & \leq\exp\left(
-C_{1}m_{n}^{-1}\varsigma_{n}+\log m_{n}\right) +C_{2}\exp\left(
-C_{3}n^{C_{4}}\right) \\
& \leq n^{-M}+C_{2}\exp\left( -C_{3}n^{C_{4}}\right)
\end{align*}
for any fixed constant $M>0$ and some positive constants $C_{1},$ $C_{2},$
$C_{3},$ and $C_{4}.$
(ii) If $\theta_{l}\gtrsim\kappa_{n}\log\left( m_{n}\right) ^{1/2}\left(
m_{n}/n\right) ^{1/2}$ with $\kappa_{n}$ specified in Assumption
\ref*{A:xi_n}$,$ then for sufficiently large $n$ we have
\[
\Pr\left( \mathcal{\hat{X}}_{l}\geq\varsigma_{n}\right) \geq1-n^{-M}
-C_{5}\exp\left( -C_{6}n^{C_{7}}\right)
\]
for any fixed constant $M>0$ and some positive constants $C_{5},$ $C_{6},$ and
$C_{7}.$
\end{proposition}
The proof of Proposition \ref{TH:main1} is tedious. We provide some technical
details in Appendix \ref{APP:tech_one_stage} before we formally prove the
proposition in Appendix \ref{APP:real_main_proof}. An implication of
Proposition \ref{TH:main1} is that for the well-chosen threshold value
$\varsigma_{n},$ the above one-stage procedure can separate the signals with
$\theta_{l}\gtrsim\kappa_{n}\log\left( m_{n}\right) ^{1/2}\left(
m_{n}/n\right) ^{1/2}$ from the noises with $\theta_{l}=0.$ Of course, we may
have some intermediate case where $0<\theta_{l}\lesssim\log\left(
m_{n}\right) ^{1/2}\left( m_{n}/n\right) ^{1/2}$ for which the above
procedure fails to do so. Note that the convergence rate for each additive
term in \cite{Stone} is $\left( m_{n}/n\right) ^{1/2}.$ Therefore, our
procedure loses some power up to the order of $\log\left( n\right) $, a
common scenario in the variable selection literature.
Following the literature and CKP, we define the true positive rates (TPR) and
the false positive rates (FPR) respectively as
\begin{align*}
\text{TPR}_{n} & =\frac{\sum_{l=1}^{p_{n}}\mathbf{1}\left(
\widehat{\mathcal{J}}_{l}=1\text{ and }\left\{ E\left[ f_{l}^{\ast}\left(
X_{l}\right) ^{2}\right] \right\} ^{1/2}\neq0\right) }{\sum_{l=1}^{p_{n}
}\mathbf{1}\left( \left\{ E\left[ f_{l}^{\ast}\left( X_{l}\right)
^{2}\right] \right\} ^{1/2}\neq0\right) },\text{ and}\\
\text{ FPR}_{n} & =\frac{\sum_{l=1}^{p_{n}}\mathbf{1}\left(
\widehat{\mathcal{J}}_{l}=1\text{ and }\left\{ E\left[ f_{l}^{\ast}\left(
X_{l}\right) ^{2}\right] \right\} ^{1/2}=0\right) }{\sum_{l=1}^{p_{n}
}\mathbf{1}\left( \left\{ E\left[ f_{l}^{\ast}\left( X_{l}\right)
^{2}\right] \right\} ^{1/2}=0\right) }.
\end{align*}
Based on the test statistic defined in equation (\ref{EQ:chil_definition}) and
its property developed in Proposition \ref{TH:main1}, we introduce the
generalized definitions of signals and noises in Table \ref{Table2}. With this
definition of hidden signals, our result in Section \ref{SEC:multi} provides
theoretical justification for the necessity of a multi-stage procedure capable
of detecting signals not identified in the first stage, yet with non-zero net
effects.\begin{table}[ptb]
\caption{The generalized definitions of signals and noises}
\label{Table2}
\centering{}\centering{ }
\begin{tabular}
[c]{l|ll}\hline\hline
& $\theta_{l}\gg\log\left( m_{n}\right) ^{1/2}\left( m_{n}/n\right)
^{1/2}$ & $\theta_{l}\lesssim\log\left( m_{n}\right) ^{1/2}\left(
m_{n}/n\right) ^{1/2}$\\\hline
$\left\{ E[f_{l}^{\ast}\left( X_{l}\right) ^{2}]\right\} ^{1/2}\neq0$ &
(I) Signals & (II) Hidden signals\\
$\left\{ E[f_{l}^{\ast}\left( X_{l}\right) ^{2}]\right\} ^{1/2}=0$ & (III)
Pseudo-signals & (IV) Noise variables\\\hline\hline
\end{tabular}
\end{table}
With a little abuse of notation, we continue to use $p^{\ast}$ and
$p^{\ast\ast}$ to denote the number of signals and pseudo-signals at the
sample level, and assume that they satisfy Assumption \ref*{A:p}. Adopting the
definitions in Table 2, we define the false discovery rates (FDR) as
\[
\text{FDR}_{n}=\frac{\sum_{l=1}^{p_{n}}\mathbf{1}\left( \widehat{\mathcal{J}
}_{l}=1\text{, }\left\{ E\left[ f_{l}^{\ast}\left( X_{l}\right)
^{2}\right] \right\} ^{1/2}=0,\text{ and }\theta_{l}\lesssim\log\left(
m_{n}\right) ^{1/2}\left( m_{n}/n\right) ^{1/2}\right) }{\sum_{l=1}
^{p_{n}}\widehat{\mathcal{J}}_{l}+1}.
\]
Apparently, the definitions in Table 2 generalize the definitions in Table 1.
The one-stage procedure is valid in terms of TPR, only if the net effects of
all signals are strong enough. Specifically, we need the following assumption.
\begin{assumption}
\label{A;no_hidden} There are no hidden signals. That is, $\theta_{j}
\gtrsim\kappa_{n}\log\left( m_{n}\right) ^{1/2}\left( m_{n}/n\right)
^{1/2}$ for some slowly divergent series $\kappa_{n}$ as in Assumption
\ref{A:xi_n} for all $j=1,2,...,p^{\ast}.\smallskip$
\end{assumption}
We note the definition of no hidden signals in the above assumption is
equivalent to the one in Table 2, given the way $\kappa_{n}$ is defined in
Assumption 8. If Assumption \ref*{A;no_hidden} fails to hold, we need the
multiple-stage procedure to pick up the hidden signals.
Based on the results in Proposition \ref{TH:main1}, we present the results for
TPR$_{n}$, FPR$_{n}$, and FDR$_{n}$ in the following theorem.
\begin{theorem}
\label{TH:main2}Suppose that $Y$ is given by equation (\ref{EQ:model}) and
Assumptions \ref*{A:iid} $-$ \ref*{A;no_hidden}\ hold. Then after some large
$n,$
(i) $E\left( \text{\emph{TPR}}_{n}\right) \geq1-C_{1}\exp\left(
-C_{2}n^{C_{3}}\right) $ for some positive constants $C_{1},C_{2}$, and
$C_{3};$
(ii) $E\left( \text{\emph{FPR}}_{n}\right) \leq p^{\ast\ast}/\left(
p_{n}-p^{\ast}\right) +C_{4}n^{-M}+C_{5}\exp\left( -C_{6}n^{C_{7}}\right) $
for any fixed positive large constant $M$ and some positive constants $C_{4},$
$C_{5},C_{6},$ and $C_{7};$
(iii) \emph{FDR}$_{n}\overset{P}{\rightarrow}0.$
\end{theorem}
Theorem \ref{TH:main2} implies that all of TPR$_{n},$ FPR$_{n}$ and FDR$_{n}$
can be well controlled provided we assume away hidden signals at the sample
level. Note that Theorem \ref{TH:main2}(i)--(ii) focuses on the asymptotic
properties of TPR$_{n}$ and FPR$_{n}$ while the last part of Theorem
\ref{TH:main2} reveals that the false discovery rate is asymptotically
vanishing in large samples.
In the next section, we turn to the multiple-stage procedure that does not
rely on Assumption \ref*{A;no_hidden}.
\section{The Multiple-Stage Procedure\label{SEC:multi}}
In this section we propose a multiple-stage procedure to select variables for
the nonparametric additive models.
\subsection{The Test Statistic with Pre-selected Variables}
As mentioned above, we may not identify a signal $X_{l}$ whose net effect
satisfies $\theta_{l}\lesssim\kappa_{n}\log\left( m_{n}\right) ^{1/2}\left(
m_{n}/n\right) ^{1/2},$ even in the case where the \textit{marginal} effect
of $X_{l}$ on $Y,$ namely, $\left\{ E\left[ f_{l}^{\ast}\left(
X_{l}\right) ^{2}\right] \right\} ^{1/2},$ is large enough.\footnote{The
marginal effect defined here is slightly different from that in the
econometrics literature. For example, the marginal effect of $X_{l}$ on $Y$ is
defined as $\beta_{l}$ in the CKP's linear model: $Y=\beta_{0}+\sum
_{l=1}^{p^{\ast}}X_{l}\beta_{l}+\varepsilon,$ but it refers to $X_{l}\beta
_{l}$.in this paper.} Consequently, Theorem \ref{TH:main2} does not hold
without Assumption \ref*{A;no_hidden} which assumes away small sample hidden
signals. In contrast, Lemma \ref{LE:hidden} in Appendix \ref{APP:tech_lemmas}
underpins the result that as long as $\left\{ E\left[ f_{l}^{\ast}\left(
X_{l}\right) ^{2}\right] \right\} ^{1/2}\gtrsim\kappa_{n}\log\left(
m_{n}\right) ^{1/2}\left( m_{n}/n\right) ^{1/2}$ for some slowly divergent
series $\kappa_{n}$ and some full rank condition holds, the net effect of
$X_{l}$ will be strong enough to be picked up at certain stage of a
multiple-stage procedure. A by-product of this Lemma is the existence of at
least one signal for Stage 1 (see the second part of this Lemma).
To introduce the multiple-stage procedure, we need some extra notations.
Suppose at certain stage after stage 1, we have pre-selected $\iota_{n}$
variables from the active set $\mathcal{S}_{n}=\left\{ X_{j},\text{
}j=1,...,p_{n}\right\} $ based on some selection procedure to be described
below. To avoid confusion, we denote these pre-selected variables as
$Z_{1},\ldots,Z_{\iota_{n}}$. Let $\boldsymbol{Z}\equiv\left( Z_{1}
,\ldots,Z_{\iota_{n}}\right) ^{\prime}\ $and $P^{m_{n}}\left( \boldsymbol{Z}
\right) =(P^{m_{n}}\left( Z_{1}\right) ^{\prime},\ldots,$ $P^{m_{n}}\left(
Z_{\iota_{n}}\right) ^{\prime})^{\prime}.$ Note that $P^{m_{n}}\left(
\boldsymbol{Z}\right) $ is an $\iota_{n}m_{n}\times1$ vector for
$\boldsymbol{Z}.$ At the next stage, we consider the nonparametric additive
regression of $Y$ on $\boldsymbol{Z}$ and an $X_{l}$ that has not been
selected so far and we do this one by one for all $p_{n}-\iota_{n}$
non-selected variables $X_{l}.$ We define the impact of $X_{l}$ on $Y$ after
controlling $\boldsymbol{Z}$ as
\begin{align*}
\theta_{l,\boldsymbol{Z}} & \equiv\left\{ E\left\{ E\left[ \left.
Y-P^{m_{n}}\left( \boldsymbol{Z}\right) ^{\prime}\Phi_{\boldsymbol{Z}}
^{-1}E\left[ P^{m_{n}}\left( \boldsymbol{Z}\right) Y\right] \right\vert
X_{l}\right] \right\} ^{2}\right\} ^{1/2}\\
& =\left\{ E\left\{ \sum_{j=1}^{p^{\ast}}E\left[ \left. f_{j}^{\ast
}\left( X_{j}\right) -P^{m_{n}}\left( \boldsymbol{Z}\right) ^{\prime}
\Phi_{\boldsymbol{Z}}^{-1}E\left[ P^{m_{n}}\left( \boldsymbol{Z}\right)
f_{j}^{\ast}\left( X_{j}\right) \right] \right\vert X_{l}\right] \right\}
^{2}\right\} ^{1/2}\\
& =\left\{ E\left[ \left( \sum_{j=1}^{p^{\ast}}\mu_{lj,\boldsymbol{Z}
}\right) ^{2}\right] \right\} ^{1/2},
\end{align*}
where $\Phi_{\boldsymbol{Z}}\equiv E[P^{m_{n}}\left( \boldsymbol{Z}\right)
P^{m_{n}}\left( \boldsymbol{Z}\right) ^{\prime}]$ and $\mu
_{lj,\boldsymbol{Z}}\equiv E\left\{ \left. f_{j}^{\ast}\left( X_{j}\right)
-P^{m_{n}}\left( \boldsymbol{Z}\right) ^{\prime}\Phi_{\boldsymbol{Z}}
^{-1}E\left[ P^{m_{n}}\left( \boldsymbol{Z}\right) f_{j}^{\ast}\left(
X_{j}\right) \right] \right\vert X_{l}\right\} $ denotes the effect of
$X_{l}$ on $f_{j}^{\ast}\left( X_{j}\right) $ after controlling the effects
of $\boldsymbol{Z}$. Apparently, we suppress the dependence of $\theta
_{l,\boldsymbol{Z}}$ on the sample size $n$.
At the sample level, let $\mathbb{Z}$ and $\mathbb{X}_{l}$ denote the $n\times
m_{n}\iota_{n}$ and $n\times m_{n}$ \textquotedblleft design
matrices\textquotedblright\ for $\boldsymbol{Z}$ and $X_{l},$ respectively.
That is,
\begin{equation}
\mathbb{Z}\equiv\left( \mathbb{Z}_{1},\ldots,\mathbb{Z}_{\iota_{n}}\right)
\text{ and }\mathbb{X}_{l}\equiv\left( \mathbb{X}_{l1},\ldots,\mathbb{X}
_{ln}\right) ^{\prime}, \label{EQ:design}
\end{equation}
where $\mathbb{Z}_{l}=\left( \mathbb{Z}_{l1},\ldots,\mathbb{Z}_{ln}\right)
^{\prime}\ $is a $n\times m_{n}$ matrix, $\mathbb{Z}_{li}=P^{m_{n}}\left(
z_{li}\right) \ $and $\mathbb{X}_{li}=P^{m_{n}}\left( x_{li}\right) .$
Define $M_{\mathbb{Z}}\equiv I_{n}-\mathbb{Z}\left( \mathbb{Z}^{\prime
}\mathbb{Z}\right) ^{-1}\mathbb{Z}^{\prime}.$ By the result of partitioned
regressions, the coefficient of $P^{m_{n}}\left( X_{l}\right) \ $is
estimated by
\[
\boldsymbol{\hat{\beta}}_{l,\boldsymbol{Z}}=\left( \mathbb{X}_{l}^{\prime
}M_{\mathbb{Z}}\mathbb{X}_{l}\right) ^{-1}\mathbb{X}_{l}^{\prime
}M_{\mathbb{Z}}\boldsymbol{y}.
\]
To determine whether $X_{l}$ should be treated as a signal variable, we
propose the following test statistic
\begin{equation}
\mathcal{\hat{X}}_{l,\boldsymbol{Z}}=\boldsymbol{\hat{\beta}}
_{l,\boldsymbol{Z}}^{\prime}\left( \hat{\sigma}_{l,\boldsymbol{Z}}
^{-2}\mathbb{X}_{l}^{\prime}M_{\mathbb{Z}}\mathbb{X}_{l}\right)
\boldsymbol{\hat{\beta}}_{l,\boldsymbol{Z}}=\left( \boldsymbol{y}^{\prime
}M_{\mathbb{Z}}\mathbb{X}_{l}\right) \left( \hat{\sigma}_{l,\boldsymbol{Z}
}^{2}\mathbb{X}_{l}^{\prime}M_{\mathbb{Z}}\mathbb{X}_{l}\right) ^{-1}\left(
\mathbb{X}_{l}^{\prime}M_{\mathbb{Z}}\boldsymbol{y}\right) \label{EQ:XlZ}
\end{equation}
where $\hat{\sigma}_{l,\boldsymbol{Z}}^{2}=n^{-1}\sum_{i=1}^{n}\hat
{\varepsilon}_{li}^{2}$ and $\hat{\varepsilon}_{li}$ is the residual from the
regression $\boldsymbol{y}$ on $\left( \mathbb{Z},\mathbb{X}_{l}\right) $.
We will study the asymptotic properties of $\mathcal{\hat{X}}
_{l,\boldsymbol{Z}}$ in the next subsection which lay down the foundation for
our multiple-stage procedure.
\subsection{The Asymptotic Properties of $\mathcal{\hat{X}}_{l,\boldsymbol{Z}
}$}
To study the asymptotic properties of $\mathcal{\hat{X}}_{l,\boldsymbol{Z}},$
we impose the following technical conditions.\smallskip
\noindent\textbf{Assumption} \textbf{\ref*{A:p}'} $p^{\ast}\ $is a positive
integer that does not vary with $n.$\ $p^{\ast\ast}\lesssim n^{B_{p^{\ast\ast
}}}\ $and $p_{n}\propto n^{B_{p}}$ for some $B_{p}>B_{p^{\ast\ast}}\geq0.$
Further, $B_{p^{\ast\ast}}<\left( 1-3B_{m}\right) /2.$\smallskip
\noindent\textbf{Assumption} \textbf{\ref*{A:fl}'} $f_{l}\left( \cdot\right)
=E\left( Y|X_{l}=\cdot\right) \in\Lambda^{d}\left( [0,1]\right) \ $with
$d>1$ for $l=1,\ldots,p_{n}$. $f_{j}^{\ast}\in\Lambda^{d}\left( [0,1]\right)
\ $with $d>1$ for $j=1,\ldots,p^{\ast}$.\smallskip
\noindent\textbf{Assumption} \textbf{\ref*{A;no_hidden}'} $\left\{
E[f_{j}^{\ast}\left( X_{j}\right) ^{2}]\right\} ^{1/2}\gtrsim\kappa_{n}
\log\left( m_{n}\right) ^{1/2}\left( m_{n}/n\right) ^{1/2}$ for some
slowly divergent series $\kappa_{n}$ as in Assumption \ref{A:xi_n}\ and for
$j=1,\ldots,p^{\ast}$.\smallskip
\begin{assumption}
\label{A:full_rank2}Let $\boldsymbol{X}_{1}^{p^{\ast}+p^{\ast\ast}}
\equiv\left( X_{1},\ldots,X_{p^{\ast}},X_{p^{\ast}+1},\ldots,X_{p^{\ast
}+p^{\ast\ast}}\right) ^{\prime},$ the vector of all signals and pseudo
signals. Similarly, let $P^{m_{n}}(\boldsymbol{X}_{1}^{p^{\ast}+p^{\ast\ast}
})\equiv\left[ P^{m_{n}}\left( X_{1}\right) ^{\prime},\ldots,P^{m_{n}
}\left( X_{p^{\ast}+p^{\ast\ast}}\right) ^{\prime}\right] ^{\prime}\ $and
$\Phi_{\boldsymbol{X}_{1}^{p^{\ast}+p^{\ast\ast}}}\equiv E[P^{m_{n}
}(\boldsymbol{X}_{1}^{p^{\ast}+p^{\ast\ast}})$ $P^{m_{n}}(\boldsymbol{X}
_{1}^{p^{\ast}+p^{\ast\ast}})^{\prime}].$ Assume that
\[
B_{X1}m_{n}^{-1}\leq\lambda_{\min}\left( \Phi_{\boldsymbol{X}_{1}^{p^{\ast
}+p^{\ast\ast}}}\right) \leq\lambda_{\max}\left( \Phi_{\boldsymbol{X}
_{1}^{p^{\ast}+p^{\ast\ast}}}\right) \leq B_{X2}m_{n}^{-1}
\]
holds for some positive constants $B_{X1}$ and $B_{X2},$ and it also holds
when $\boldsymbol{X}_{1}^{p^{\ast}+p^{\ast\ast}}$ is augmented by an arbitrary
element $X_{l}$ with $l>p^{\ast}+p^{\ast\ast}.$
\end{assumption}
\begin{assumption}
\label{A:noisevariable}If the pre-selected variables $\boldsymbol{Z}$ are
either signals or pseudo-signals, the effect of noise variable defined in
Table 2 on $Y$ is still weak. That is, $\theta_{l,\boldsymbol{Z}}\lesssim
\log\left( m_{n}\right) ^{1/2}\left( m_{n}/n\right) ^{1/2}$ for all ${l}
${$=p^{\ast}+p^{\ast\ast}+1,\ldots,p_{n},$} when all variables in
$\boldsymbol{Z}$ are either signals or pseudo-signals.\smallskip
\end{assumption}
Assumption \ref*{A:p}' strengthens Assumption \ref*{A:p} by adding one more
condition on $B_{p^{\ast\ast}}$. It is imposed to ensure the good property of
$\boldsymbol{\tilde{u}}_{l,\boldsymbol{Z}}$ defined in equation
(\ref{EQ:utidelz}) (see Lemma \ref{LE:error_bound}). This condition can be
very restrictive on $B_{p^{\ast\ast}}.$ For example, when $B_{m}=1/4,$ this
condition implies $B_{p^{\ast\ast}}<1/8$. Assumption \ref*{A:fl}' strengthens
Assumption \ref*{A:fl} so that equation (\ref{EQ:appro}) in Appendix A.2
holds. As remarked at the beginning of last subsection, we do not impose
Assumption \ref*{A;no_hidden} for the multiple-stage procedure, as long as the
marginal effect of the signal is not too weak as imposed in Assumption
\ref*{A;no_hidden}'. This assumption seems inevitable. It is analogous to
Assumption 6 in CKP for the linear regression model and similar to the
so-called `beta-min' condition that is commonly assumed in the penalized
regression literature (see, e.g., Chapter 7.4 of \cite{Buhlmann2011}).
Assumption \ref*{A:full_rank2} is the common rank condition for nonparametric
additive regressions. We also require it to hold for the case when we add one
noise variable. This seems inevitable because we will run the regression with
regressors being all signals and pseudo-signals plus one noise variable in the
multiple-stage procedure with very high probability.\ Assumption
\ref*{A:noisevariable} is also inevitable and implicitly imposed in CKP for
the linear regression models\emph{.}
It is worth mentioning that Assumption \ref*{A:full_rank2} plays a similar
role to the \textquotedblleft restrictive eigenvalues\textquotedblright
\ condition in \cite{BickelRitovTsyvakov2009} and \cite{BelloniEtal2012}. The
restrictive-eigenvalues\ condition requires certain full rank conditions on
all possible design matrices composed of a certain number of covariates. In
contrast, our OCMT only requires full rank conditions on the design matrices
composed of signals and pseudo-signals, and permits arbitrary correlations
among the noise variables. Obviously, these two sets of conditions are
non-nested. In addition, \cite{HuangEtal2010} impose essentially the same set
of assumptions except the rank conditions discussed here. The main advantage
of the OCMT is that it does not require any numerical min-search of an
objective function, can be computed much faster, and deliver more reliable results.
The following proposition presents the probability bounds for the
\textquotedblleft Type-I\textquotedblright\ and \textquotedblleft
Type-II\textquotedblright\ errors when we have some pre-selected variables.
\begin{proposition}
\label{TH:main1_ms}Suppose that $Y$ is given by equation (\ref{EQ:model}),
Assumptions \ref*{A:iid}, \ref*{A:p}', \ref*{A:supp}, \ref*{A:epsilon},
\ref*{A:fl}', \ref*{A:tech}, \ref*{A:mn}, \ref*{A:xi_n}, and
\ref*{A:full_rank2} hold, and the pre-selected variables $\boldsymbol{Z}$ are
either signals or pseudo-signals.
(i) If $\theta_{l,\boldsymbol{Z}}\lesssim\log\left( m_{n}\right) ^{1/2}
m_{n}^{1/2}n^{-1/2},$ then
\begin{align}
\Pr\left( \mathcal{\hat{X}}_{l,\boldsymbol{Z}}\geq\varsigma_{n}\right) &
\leq\exp\left( -C_{1}m_{n}^{-1}\varsigma_{n}+\log m_{n}\right) +C_{2}
\exp\left( -C_{3}n^{C_{4}}\right) \label{EQ:xl_FDR}\\
& \leq n^{-M}+C_{2}\exp\left( -C_{3}n^{C_{4}}\right) \nonumber
\end{align}
for any fixed large constant $M>0$ and some positive constants $C_{1}
,C_{2},C_{3},$ and $C_{4},$ after some large $n.$
(ii) If $\theta_{l,\boldsymbol{Z}}\gtrsim\kappa_{n}\log\left( m_{n}\right)
^{1/2}m_{n}^{1/2}n^{-1/2}$ with $\kappa_{n}$ specified in Assumption
\ref{A:xi_n}$,$ then
\[
\Pr\left( \mathcal{\hat{X}}_{l,\boldsymbol{Z}}\geq\varsigma_{n}\right)
\geq1-n^{-M}-C_{5}\exp\left( -C_{6}n^{C_{7}}\right)
\]
for any fixed large constant $M>0$ and some positive constants $C_{5},C_{6},$
and $C_{7},$ after some large $n.$
\end{proposition}
The proof of Proposition \ref{TH:main1_ms} is rather tedious. We provide some
technical discussions on the proof in Appendix \ref{SEC:multi_tech} before we
formally prove it in Appendix \ref{APP:real_main_proof}.
Like Proposition \ref{TH:main1}, Proposition \ref{TH:main1_ms} implies that
for the well-chosen threshold value $\varsigma_{n},$ the use of the test
statistic $\mathcal{\hat{X}}_{l,\boldsymbol{Z}}$ helps to separate variables
with large value of $\theta_{l,\boldsymbol{Z}}$ from those with small value of
$\theta_{l,\boldsymbol{Z}}.$ This observation will be used in our
multiple-stage procedure to select all signal variables.
\subsection{The Multiple-Stage Procedure\label{SEC:multi_procedure}}
We present the multiple-stage procedure as follows.
We conduct the first-stage selection as in Section \ref{SEC:1stage} by
constructing the test statistic $\mathcal{\hat{X}}_{l}$ as in equation
(\ref{EQ:chil_definition}) and using the threshold value $\varsigma_{n}$ that
satisfies the condition in Assumption \ref*{A:xi_n}. We re-label the selection
indicator in equation (\ref{EQ:select1stage}) as
\[
\widehat{\mathcal{J}}_{l,\left( 1\right) }=\mathbf{1}\left( \mathcal{\hat
{X}}_{l}>\varsigma_{n}\right) \text{ for }l=1,2,\ldots,p_{n}.
\]
We collect all the variables selected in stage 1 into the vector
$\boldsymbol{Z}_{\left( 1\right) }$, and denote the index set of the
selected variables by $S_{\left( 1\right) }.$ For the second stage, we
denote the index set of the active variables in stage 2 by $\Psi_{\left(
2\right) }$ where $\Psi_{\left( 2\right) }=\left\{ 1,2,\ldots
,p_{n}\right\} \backslash S_{\left( 1\right) }$. In the second stage, we
regress $Y$ on $P^{m_{n}}\left( X_{l}\right) $ with $P^{m_{n}}\left(
\boldsymbol{Z}_{\left( 1\right) }\right) $ as pre-selected variables one by
one for $l\in\Psi_{\left( 2\right) }.$ We construct the test statistic
$\mathcal{\hat{X}}_{l,\boldsymbol{Z}_{\left( 1\right) }}$ as in equation
(\ref{EQ:XlZ}). We select the variable $X_{l}$ if $\widehat{\mathcal{J}
}_{l,\left( 2\right) }=1,$ where
\[
\widehat{\mathcal{J}}_{l,\left( 2\right) }=\mathbf{1}\left( \mathcal{\hat
{X}}_{l,\boldsymbol{Z}_{\left( 1\right) }}>\varsigma_{n}\right) \text{ for
}l\in\Psi_{\left( 2\right) }.
\]
We add all the variables selected in stage 2 into the set of variables
selected in stage 1 as a new vector, and we denote it as $\boldsymbol{Z}
_{\left( 2\right) }.$ We denote the index set of the selected variables
($\boldsymbol{Z}_{\left( 2\right) }$) by $S_{\left( 2\right) }$ and the
index set of the active variables for stage 3 as $\Psi_{\left( 3\right) },$
where $\Psi_{\left( 3\right) }=\left\{ 1,2,\ldots,p_{n}\right\} \backslash
S_{\left( 2\right) }.$ And so on and so forth. For stage $k,$ we denote the
pre-selected variables as $\boldsymbol{Z}_{\left( k-1\right) },$ and the
index set of the active variables as $\Psi_{\left( k\right) }.$ Then we
regress $Y$ on $P^{m_{n}}\left( X_{l}\right) $ with $P^{m_{n}}\left(
\boldsymbol{Z}_{\left( k-1\right) }\right) $ as pre-selected variables one
by one for $l\in\Psi_{\left( k\right) }.$ We construct the test statistic
$\mathcal{\hat{X}}_{l,\boldsymbol{Z}_{\left( k-1\right) }}$ as in equation
(\ref{EQ:XlZ}). We select the variable $X_{l}$ if $\widehat{\mathcal{J}
}_{l,\left( k\right) }=1,$ where
\[
\widehat{\mathcal{J}}_{l,\left( k\right) }=\mathbf{1}\left( \mathcal{\hat
{X}}_{l,\boldsymbol{Z}_{\left( k-1\right) }}>\varsigma_{n}\right) \text{
for }l\in\Psi_{\left( k\right) }.
\]
We add all the variables selected in stage $k$ into the set of variables
selected in stage $k-1$ as a new vector, and we denote it as $\boldsymbol{Z}
_{\left( k\right) }.$ We stop the procedure at a stage in which no new
variables are selected. We denote the stage, in which one or more variables
are selected but no new variables are selected after that, as $\hat{k}_{s}$.
So the OCMT procedures stops after stage $\hat{k}_{s}.$ The selection
indicator for variable $X_{l}$ of the OCMT procedure is defined as follows
\begin{equation}
\widehat{\mathcal{J}}_{l}=\sum_{k=1}^{\hat{k}_{s}}\widehat{\mathcal{J}
}_{l,\left( k\right) }. \label{EQ:updatedJ}
\end{equation}
By construction, $\widehat{\mathcal{J}}_{l}$ is either 1 or 0. It takes value
$1$ if $X_{l}$ is selected in the OCMT\ procedure and 0 otherwise.
The following theorem mainly studies the asymptotic properties of the
multiple-stage procedure in terms of TPR, FPR\ and FDR.
\begin{theorem}
\label{TH:stoppingFDR}Suppose\ that Assumptions \ref*{A:iid}, \ref*{A:p}',
\ref*{A:supp}, \ref*{A:epsilon}, \ref*{A:fl}', \ref*{A:tech}, \ref*{A:mn},
\ref*{A:xi_n}, \ref*{A;no_hidden}', \ref*{A:full_rank2}, and
\ref*{A:noisevariable} hold. Then after some large $n,$
(i) $\Pr\left( \hat{k}_{s}>p^{\ast}\right) \leq n^{-M_{6}}+C_{19}\exp\left(
-C_{20}n^{C_{21}}\right) $ for some fixed large positive number $M_{6}$ and
some positive constants $C_{19},C_{20},$ and $C_{21};$
(ii) $E\left( \text{\emph{TPR}}_{n}\right) \geq1-C_{1}\exp\left(
-C_{2}n^{C_{3}}\right) $ for some positive constants $C_{1},C_{2}$, and
$C_{3};$
(iii) $E\left( \text{\emph{FPR}}_{n}\right) \leq p^{\ast\ast}/\left(
p_{n}-p^{\ast}\right) +C_{4}n^{-M}$ for some positive $C_{4}$ and any fixed
positive large constant $M;$
(iv) \emph{FDR}$_{n}\overset{P}{\rightarrow}0.$
\end{theorem}
Theorem \ref{TH:stoppingFDR}(i) implies the OCMT procedure can terminate at
step $p^{\ast}$ with very high probability. Theorem \ref{TH:stoppingFDR}
(ii)-(iv) implies that all of TPR$_{n},$ FPR$_{n}$ and FDR$_{n}$ can be well
controlled. Of course, when $n\rightarrow\infty,$ we need to conduct at most
$p^{\ast}$-stage procedure to determine all signals and eliminate all noise
variables, with very high probability.
\subsection{Dealing with Pseudo-signals\label{SEC:Pseudo}}
Because of the nature of the OCMT procedure, pseudo-signals cannot be excluded
from the selection list with high probability. To eliminate pseudo-signals, we
propose to employ the adaptive group Lasso after the OCMT. Then by the
properties of the adaptive group Lasso, the event that the pseudo-signals are
excluded and the signals are kept occurs with probability approaching one
(w.p.a.1). We present this result in Theorem \ref{TH:AGLasso} below. We denote
the set of the variables selected by the OCMT procedure as $\hat
{S}_{\text{OCMT}}\equiv S_{\left( \hat{k}_{s}\right) },$ collect them into a
vector $\boldsymbol{Z}_{\text{OCMT}},$ and denote its dimension as $\hat
{p}_{\text{OCMT}}$. Similarly, let $\mathbb{Z}_{\text{OCMT}}$ ($n\times\hat
{p}_{\text{OCMT}}m_{n}$) denote the design matrix for $\boldsymbol{Z}
_{\text{OCMT}}$ as in equation (\ref{EQ:design}). The post-OCMT adaptive group
Lasso procedure goes as follows:
\begin{enumerate}
\item Obtain the group Lasso estimator by searching $\boldsymbol{\beta}_{n}$
($\hat{p}_{\text{OCMT}}m_{n}\times1$)\ to minimize
\[
L_{n1}\left( \boldsymbol{\beta}_{n},\lambda_{n1}\right) =\left\Vert
\boldsymbol{y}-\mathbb{Z}_{\text{OCMT}}\boldsymbol{\beta}_{n}\right\Vert
^{2}+\lambda_{n1}\sum_{j=1}^{\hat{p}_{\text{OCMT}}}\left\Vert
\boldsymbol{\beta}_{nj}\right\Vert ,
\]
where $\lambda_{n1}$ is a positive tuning parameter, $\boldsymbol{\beta}
_{n}=(\boldsymbol{\beta}_{n1}^{\prime},...,\boldsymbol{\beta}_{n\hat
{p}_{\text{OCMT}}}^{\prime})^{\prime},$ $\boldsymbol{\beta}_{nj}$ is a
$m_{n}\times1$ vector of coefficients of the B-spline basis for the $j$-th
element in $\boldsymbol{Z}_{\text{OCMT}}.$ Denote the above estimator as
$\boldsymbol{\tilde{\beta}}_{n}.$
\item The adaptive group Lasso estimator is obtained by searching
$\boldsymbol{\beta}_{n}$ to minimize
\[
L_{n2}\left( \boldsymbol{\beta}_{n},\lambda_{n2}\right) =\left\Vert
\boldsymbol{y}-\mathbb{Z}_{\text{OCMT}}\boldsymbol{\beta}_{n}\right\Vert
^{2}+\lambda_{n2}\sum_{j=1}^{\hat{p}_{\text{OCMT}}}\frac{1}{\left\Vert
\boldsymbol{\tilde{\beta}}_{nj}\right\Vert }\left\Vert \boldsymbol{\beta}
_{nj}\right\Vert ,
\]
where we use the convention that $0/0=0$. Denote the above estimator as
$\boldsymbol{\hat{\beta}}_{n}.$
\end{enumerate}
The post-selection estimation proceeds as follows. Denote the selected
regressor from the above procedure as $\boldsymbol{Z}_{\text{AGLASSO}},$ and
similarly denote its B-spline basis and design matrix as $P^{m_{n}}\left(
\boldsymbol{Z}_{\text{AGLASSO}}\right) $ and $\mathbb{Z}_{\text{AGLASSO}},$
respectively$.$ The post selection estimator is the OLS estimator of
regressing $\boldsymbol{y}$ on $\mathbb{Z}_{\text{AGLASSO}}$, which is
\[
\boldsymbol{\hat{\beta}}_{\text{post}}=\left( \mathbb{Z}_{\text{AGLASSO}
}^{\prime}\mathbb{Z}_{\text{AGLASSO}}\right) ^{-1}\mathbb{Z}_{\text{AGLASSO}
}^{\prime}\boldsymbol{y.}
\]
The final fitted model is
\[
P^{m_{n}}\left( \boldsymbol{Z}_{\text{AGLASSO}}\right) ^{\prime
}\boldsymbol{\hat{\beta}}_{\text{post}}.
\]
\begin{theorem}
\label{TH:AGLasso} Suppose\ that Assumptions \ref*{A:iid}, \ref*{A:p}',
\ref*{A:supp}, \ref*{A:epsilon}, \ref*{A:fl}', \ref*{A:tech}, \ref*{A:mn},
\ref*{A:xi_n}, \ref*{A;no_hidden}', \ref*{A:full_rank2}, and
\ref*{A:noisevariable} hold. Further, $\lambda_{n1}\geq C\sqrt{n\log\left(
p^{\ast\ast}m_{n}\right) }$ for a sufficient large $C$, $\lambda_{n1}\ll
\sqrt{n/m_{n}},$ and $m_{n}^{1/2}\log\left( p^{\ast\ast}m_{n}\right)
\ll\lambda_{n2}\ll nm_{n}^{-1/4}.$ Then
(i) All signal variables are kept and all pseudo signals or noise variables
are eliminated w.p.a.1;
(ii) The post OCMT estimation error satisfies
\[
P^{m_{n}}\left( \boldsymbol{Z}_{\text{AGLASSO}}\right) ^{\prime
}\boldsymbol{\hat{\beta}}_{\text{post}}-
{\displaystyle\sum_{j=1}^{p^{\ast}}}
f_{j}^{\ast}\left( X_{j}\right) =O_{P}\left( \left( m_{n}/n\right)
^{1/2}\right) .
\]
\end{theorem}
The above theorem is almost the same as that in \cite{HuangEtal2010} including
the requirements on $\lambda_{n1}$ and $\lambda_{n2}$, with the exception that
the procedure starts with the covariates post OCMT. Consequently, we only need
to show that the adaptive group Lasso procedure is still valid in the
post-OCMT situation, and the rest follows immediately from
\cite{HuangEtal2010}. In particular, under Assumption \ref*{A:mn}, the biases
of the post OCMT estimators of the nonparametric additive components are
asymptotically negligible so that their mean square errors (MSEs) are
dominated by their asymptotic variances that are of order $O\left(
m_{n}/n\right) ,$ which explains the result in Theorem \ref{TH:AGLasso}(ii).
Our procedure enjoys the same theoretical property as the adaptive group
Lasso. After the OCMT, the dimension of candidate variables is reduced
dramatically. Thus the additional computation burden applying the post-OCMT
adaptive group Lasso can be almost ignored. We note that the main advantage of
our procedure is fast and reliable computation, which delivers better small
sample performance as shown from our simulation studies. An implication of the
above theorem is that the post-OCMT adaptive group procedure improves over the
post-selection estimation results in CKP (Theorem 2 in CKP) because
pseudo-signals are now eliminated with very high probability.
\subsection{Diverging $p^{\ast}$\label{SEC:divergeP}}
Allowing a diverging $p^{\ast}$ (number of true signals) is possible. The only
additional technical condition apart from the restriction on the speed of
$p^{\ast}$ is that $\sum_{j=1}^{p^{\ast}}f_{j}^{\ast}\left( X_{j}\right) $
is uniformly bounded. Note that this condition naturally holds for a fixed
$p^{\ast}$ due to the boundedness of $f_{j}^{\ast}$. The main reason for the
requirement of this condition is technical: the uniform boundedness of
$\sum_{j=1}^{p^{\ast}}f_{j}^{\ast}\left( X_{j}\right) $ ensures that $U_{l}
$\ defined in equation (\ref{EQ:betanl}) also satisfies the exponential
decayed tail condition,\footnote{For details, see the proof of Lemma
\ref{LE:error_var}.} and the tail property is necessary to apply the main
inequalities to obtain the probability bounds. The uniform boundedness of
$\sum_{j=1}^{p^{\ast}}f_{j}^{\ast}\left( X_{j}\right) $ was also imposed in
\cite{FanFengSong}$.$ Propositions \ref{TH:main1} and \ref{TH:main1_ms} are on
individual $X_{l},$\ but we do need Assumption \ref*{A:p}\textquotedblright
\ on $p^{\ast}$ so that Proposition \ref{TH:main1_ms} holds. It is to ensure
that we can have precise estimation on a diverging design matrix.
\noindent\textbf{Assumption} \textbf{\ref*{A:p}\textquotedblright} $p^{\ast
}\lesssim n^{B_{p^{\ast}}}$ for some $B_{p^{\ast}}\geq0.$\ $p^{\ast\ast
}\lesssim n^{B_{p^{\ast\ast}}}\ $and $p_{n}\propto n^{B_{p}}$ for some
$B_{p}>B_{p^{\ast\ast}},B_{p^{\ast}}\geq0.$ Further, $B_{p^{\ast}}
+B_{p^{\ast\ast}}<\left( 1-3B_{m}\right) /2.$
We present the main results in the following theorem.
\begin{theorem}
\label{TH:stoppingFDR2}Suppose\ that Assumptions \ref*{A:iid}, \ref*{A:p}
\textquotedblright, \ref*{A:supp}, \ref*{A:epsilon}, \ref*{A:fl}',
\ref*{A:tech}, \ref*{A:mn}, \ref*{A:xi_n}, \ref*{A;no_hidden}',
\ref*{A:full_rank2}, and \ref*{A:noisevariable} hold. In addition, assume that
$\sum_{j=1}^{p^{\ast}}f_{j}^{\ast}\left( X_{j}\right) $ is uniformly
bounded. Then, the results in Theorem \ref{TH:stoppingFDR} continue to hold.
\end{theorem}
In the next section, we investigate the small sample performance of our
procedure by means of Monte Carlo experiment.
\section{Monte Carlo Simulations\label{SEC:MC}}
To investigate the finite-sample performance of our procedure, we conduct
Monte Carlo experiments in this section.
\subsection{Simulation Design}
Following \cite{HuangEtal2010}, we consider the following data generating
processes (DGPs). In what follows, we assume
\begin{align*}
f_{1}\left( x\right) & =x;\text{\ }f_{2}\left( x\right) =(2x-1)^{2}
;\text{\ }f_{3}\left( x\right) =\frac{\sin\left( 2\pi x\right) }
{2-\sin\left( 2\pi x\right) };\text{ and}\\
f_{4}\left( x\right) & =0.1\sin\left( 2\pi x\right) +0.2\cos\left( 2\pi
x\right) +0.3\sin\left( 2\pi x\right) ^{2}+0.4\cos\left( 2\pi x\right)
^{3}+0.5\sin\left( 2\pi x\right) ^{3}.
\end{align*}
For the errors, we assume $\varepsilon\sim$ i.i.d. $N\left( 0,1\right) $ for
DGPs 1--6 and 9--10. In DGPs 7--8, we check the impact of heteroskedastic
errors on our methods. Specifically, we add a heteroskedastic error to the
simplest and the most complicated designs in DGPs 1--6 to form DGPs 7 and 8,
respectively. In DGPs 9--10, we consider the case with additive components of
binary variables that mimic the application in Section \ref{SEC:application}.
\textbf{DGP 1: Four independent signals only}. $Y$ is generated as follows:
\begin{equation}
Y=2.55f_{1}\left( X_{1}\right) +2.57f_{2}\left( X_{2}\right)
+1.68f_{3}\left( X_{3}\right) +f_{4}\left( X_{4}\right) +\varepsilon.
\label{model.yi_simulation_D1}
\end{equation}
Note the coefficients before $f_{i}$'s\ are set to make each signal have the
same strength in terms of variance\ for independent uniform $X_{1}
,\ldots,X_{4}$. The covariates are generated as follows$:$
\[
X_{j}=W_{j}\text{ for }j=1,\ldots,4,\text{ and }X_{j}=\frac{W_{j}+U_{1}}
{2}\text{for }j\geq5,
\]
where $W_{j},$ $j=1,\ldots,p_{n}$, and $U_{1}$ are all independent draws from
$U(0,1)$. Thus, $p^{\ast}=4$ and $p^{\ast\ast}=0$ for DGP 1.\ Define the
Signal-to-noise ratio to be $r_{sn}=\frac{\mathtt{sd}\left( f\right)
}{\mathtt{sd}\left( \varepsilon\right) }$, and $r_{sn}=1.5$ for DGP 1.
\textbf{DGP 2: Four independent signals and two pseudo-signals}. $Y$ is
generated from equation (\ref{model.yi_simulation_D1}). The covariates are
generated as follows$:$
\begin{align*}
X_{j} & =W_{j}\text{ for }j=1,\ldots,4,\text{ }X_{5}=\frac{4X_{1}+U_{1}}
{5}\text{, }X_{6}=\frac{4X_{2}+U_{2}}{5},\text{ and }\\
X_{j} & =\frac{W_{j-2}+U_{3}}{2}\text{ for\ }j\geq7,
\end{align*}
where $W_{j},$ $j=1,\ldots,p_{n}-2$, $U_{1},$ and $U_{2}$\ are all independent
draws from $U(0,1)$. Thus, $p^{\ast}=4$ and $p^{\ast\ast}=2$ for DGP 2.
\textbf{DGP 3: Four signals, and one hidden signal}. $Y$ is generated from
\begin{equation}
Y=2.55f_{1}\left( X_{1}\right) +2.57f_{2}\left( X_{2}\right)
+1.68f_{3}\left( X_{3}\right) +f_{4}\left( X_{4}\right) +f_{5}\left(
X_{5}\right) +\varepsilon, \label{EQ:hidden_model}
\end{equation}
where $f_{5}\left( X_{5}\right) =-\mathbb{E}\left[ 2.55f_{1}\left(
X_{1}\right) +2.57f_{2}\left( X_{2}\right) +1.68f_{3}\left( X_{3}\right)
+f_{4}\left( X_{4}\right) |X_{5}\right] .$ The covariates are generated as
follows $:$
\begin{align*}
X_{j} & =W_{j}\text{ for }j=1,2,\text{ }X_{j}=\frac{W_{j}+U_{1}}{2}\text{
for }j=3,4,\text{ }X_{5}=U_{1}\text{, and}\\
X_{j} & =\frac{W_{j-1}+U_{2}}{2}\text{ for }j\geq6,
\end{align*}
where $W_{j},$ $j=1,\ldots,p_{n}-1,$ $U_{1},$ and $U_{2}$\ are independent
draws from $U(0,1)$. Then, $p^{\ast}=5$ and $p^{\ast\ast}=0$ for DGP 3, and
the fifth signal is hidden by our definition.\footnote{By the distribution of
covariates,
\begin{align*}
f_{5}\left( x\right) & \approx0.97\pi-1.2\cos(\pi x)+\sin\left( \pi
x\right) +0.6861\pi\arctan[2\left( \tan(\pi x/2)-1\right) /\sqrt{3}]\\
& -0.6861\pi\arctan[2\left( \tan(\pi x/2+\pi/2)-1\right) /\sqrt
{3}]+0.2778\cos^{3}(\pi x)-0.2222\sin^{3}(\pi x).
\end{align*}
}
\textbf{DGP 4: Four signals, two pseudo-signals, and one hidden signal}. $Y$
is generated from equation (\ref{EQ:hidden_model}). The covariates are
generated as follows $:$
\begin{align*}
X_{j} & =W_{j}\ \text{for }j=1,2,\text{ }X_{j}=\frac{W_{j}+U_{1}}
{2}\ \text{for }j=3,4,\text{ }X_{5}=U_{1},\text{\ }X_{6}=\frac{4X_{1}+U_{2}
}{5}\text{,}\\
\text{ }X_{7} & =\frac{4X_{2}+U_{3}}{5},\text{ and }X_{j}=\frac
{W_{j-3}+U_{3}}{2}\ \text{for }j\geq8,
\end{align*}
where $W_{j},$ $j=1,\ldots,p_{n}-3,$ $U_{1},$ $U_{2},$ and $U_{3}$\ are
independent draws from $U(0,1)$. Then, $p^{\ast}=5$ and $p^{\ast\ast}=2$ for
DGP 4, and the fifth signal is a hidden signal.
\textbf{DGP 5: Four correlated signals.} $Y$ is generated from equation
(\ref{model.yi_simulation_D1}). The covariates are generated as follows$:$
\begin{equation}
X_{j}=\frac{W_{j}+U_{1}}{2}\text{ for\ }j=1,\ldots,4,\text{ and }X_{j}
=\frac{W_{j}+U_{2}}{2}\ \text{for\ }j\geq5,\nonumber
\end{equation}
where $W_{j},$ $j=1,\ldots,p_{n}$, $U_{1},$ and $U_{2}$\ are independent draws
from $U(0,1)$. Thus, four signals are correlated with each other, and
$p^{\ast}=4$ and $p^{\ast\ast}=0$ for DGP 5.
\textbf{DGP 6: Four signals, many pseudo-signals, and one hidden signal}. $Y$
is generated from equation (\ref{EQ:hidden_model}). The covariates are
generated as follows $:$
\begin{align*}
X_{j} & =W_{j}\text{ for }j=1,2,\text{ }X_{j}=\frac{W_{j}+U_{1}}
{2}\ \text{for }j=3,4,\text{ }X_{5}=U_{1},\\
X_{j} & =\frac{4X_{1}+\left( j-5\right) W_{j-1}}{j-1}\text{ for
}j=6,10,14,18,\ldots,\\
X_{j} & =\frac{4X_{2}+\left( j-5\right) W_{j-1}}{j-1}\text{ for
}j=7,11,15,19,\ldots,\\
X_{j} & =\frac{4X_{3}+\left( j-5\right) W_{j-1}}{j-1}\text{ for
}j=8,12,16,20,\ldots,\text{ and}\\
X_{j} & =\frac{4X_{4}+\left( j-5\right) W_{j-1}}{j-1}\text{ for
}j=9,13,17,21,\ldots
\end{align*}
where $W_{j},$ $j=1,\ldots,p_{n}-1,$ and $U_{1}$\ are independent draws from
$U(0,1)$. Then, $p^{\ast}=5$ with one hidden signal for DGP 6.
\textbf{DGP 7: Four independent signals with heteroskedastic errors}. $Y$ is
generated from equation (\ref{model.yi_simulation_D1}) with the same
covariates as in DGP 1. We assume that conditioning on $X,$ $\varepsilon$ is
normal with mean 0 and variance $0.436[1+\left( X_{1}+X_{2}+X_{3}
+X_{4}\right) /4]^{2},$ and the unconditional variance of $\varepsilon$ is
approximately 1.
\textbf{DGP 8: Four signals, many pseudo-signals, and one hidden signal with
heteroskedastic errors}. $Y$ is generated from equation (\ref{EQ:hidden_model}
) with the same covariates as in DGP 6. We assume that conditioning on $X,$
$\varepsilon$ is normal with mean 0 and variance $0.436[1+\left( X_{1}
+X_{2}+X_{3}+X_{4}\right) /4]^{2}$.
\textbf{DGP 9: Four independent signals with some binary variables. }To mimic
the application, we consider the situation with some binary covariates. Note
that any function of a binary covariate can at most take two values. Without
loss of generality, we focus on the case in which those binary covariates
enter the model linearly. It is easy to see that our theoretical results
continue to hold in the presence of some linear additive components with the
main difference that they do not exhibit any approximation bias.\footnote{In
this case, the test statistics for the additive linear terms are the squares
of t-statistics, and the inequality continues to hold by setting, for example,
$\varsigma_{n}\propto\left[ \log p_{n}\right] ^{1.1}$ for the case when
$p_{n}$ is a polynomial of $n.$} We set the threshold as $\varsigma
_{n}=C\left[ \log p_{n}\right] ^{1.1}$. $Y$ is generated from
\[
Y=2.57f_{2}\left( X_{1}\right) +1.68f_{3}\left( X_{2}\right)
+1.47X_{3}+1.47X_{4}+\varepsilon,
\]
where $X_{j}=W_{j},$ $j=1,2,$ $X_{j}=V_{j-2},$ $j=3,4,$ $W_{1}$ and $W_{2}$
are independent $U(0,1),$ and $V_{1}$ and $V_{2}$ are independent Bernoulli
random variables with equal chances of taking value 0 or 1. The remaining
covariates are generated as follows:
\[
X_{j}=\frac{W_{j-2}+U_{1}}{2}\text{ for }j=5,6,...,\frac{p_{n}}{2},\text{ and
}X_{j}=V_{j-\frac{p_{n}}{2}+2}\text{ for }j=\frac{p_{n}}{2}+1,...,p_{n},
\]
where $W_{j},$ $j=5,6,...,\frac{p_{n}}{2}-2,$ and $U_{1},$ are independent
$U(0,1),$ $V_{j},$ $j=3,...,\frac{p_{n}}{2}+2$ are distributed the same as
$V_{1}.$ All $W$s$,V$s, and $U$ are independent of each other.
\textbf{DGP 10: Four signals with one hidden signal in the presence of some
binary variables. }For this DGP, $Y$ is generated from
\[
Y=2.57f_{2}\left( X_{1}\right) +1.5X_{2}+1.5X_{3}-X_{4}+\varepsilon,
\]
where $X_{1}=W_{1}$ is a $U(0,1),$ $X_{j}=V_{j-1},$ $j=2,3,4,$ are Bernoulli
random variables with equal chances of taking 0 or 1, and Corr$\left(
V_{1},V_{3}\right) =$Corr$\left( V_{1},V_{2}\right) =$Corr$\left(
V_{2},V_{3}\right) =1/3$. $X_{j},$ $j=5,6,...,p_{n}$ are the same as those in
DGP 9.\footnote{An example of $V_{1},V_{2},$ and $V_{3}$ is that
$V_{j}=\mathbf{1(}\tilde{U}_{j}+\tilde{U}_{4}>1)$ for $j=1,2,$ and $3,$ where
$\tilde{U}_{1},...,\tilde{U}_{4}$ are independent $U\left( 0,1\right) .$}
Some simple calculation implies that $X_{4}$ is a hidden signal.
\subsection{Tuning Parameters\label{SEC:tuning}}
The key tuning parameter in this study is the threshold $\varsigma_{n}.$\ The
CKP's Bonferroni correction strategy does not perform consistently well for
our case. We do the following instead. We set
\[
\varsigma_{n}=Cm_{n}\left[ \left( \log p_{n}\right) ^{1.1}+\left( \log
m_{n}\right) ^{1.1}\right] \text{ and }\varsigma_{n}=C\left( \log
p_{n}\right) ^{1.1}
\]
for continuous variables and binary variables, respectively. This satisfies
Assumption \ref{A:xi_n} when, in addition, Assumptions \ref{A:p} and
\ref{A:mn} hold (both $m_{n}$ and $p_{n}$ grow at the polynomial rate of $n$).
In the small samples, we set this $\varsigma_{n}$ to have control over both
$m_{n}$ and $p_{n}$. CKP suggest using a larger threshold, $\varsigma
_{n}^{\ast}$, for subsequent stages, to improve the finite sample performance.
We follow their lead and set a $\varsigma_{n}^{\ast}$ larger than
$\varsigma_{n}$\ for subsequent stages. That is, we replace $\varsigma_{n}$
with $\varsigma_{n}^{\ast}$ for $\widehat{\mathcal{J}}_{l,\left( 2\right)
},\widehat{\mathcal{J}}_{l,\left( 3\right) },...,\widehat{\mathcal{J}
}_{l,\left( k\right) }$ in Section \ref{SEC:multi_procedure}. The reason is
that the chance of including a noise variable increases quickly for subsequent
stages, and we need a larger $\varsigma_{n}^{\ast}$ for the OCMT to conclude
more easily. We set $\varsigma_{n}^{\ast}=3\varsigma_{n}$ for continuous
variables. Note that CKP take the threshold for later stages to be twice the
threshold for the first stage, and their test is based on $t$-statistics.
Since ours is based on chi-squared statistics, equivalently, we should set
$\varsigma_{n}^{\ast}$ to be $4\varsigma_{n}.$ However, some small-scale
experiments suggest that setting $\varsigma_{n}^{\ast}=3\varsigma_{n}$ can
yield better small sample performance, probably because our model is
nonparametric and our test statistic varies more than the parametric
counterpart. Note that this change does not affect our theoretical results
because $\varsigma_{n}^{\ast}$ is proportional to $\varsigma_{n}$. For binary
additive components, we continue to follow CKP and set $\varsigma_{n}^{\ast
}=4\varsigma_{n}$ because they enter the model linearly.
Another important tuning parameter is $m_{n},$ which is critical for the sieve
estimation. The optimal choice of $m_{n}$ has been studied extensively in the
literature. Popular ways to choose the value of $m_{n}$ include cross
validation, Akaike information criterion (AIC), and BIC. We refer the readers
to \cite{Chen_handbook} and \cite{Hansen2014} for a review on this important
issue. For $m_{n},$ we simply set $m_{n}=\left\lfloor n^{1/4}\right\rfloor
+1,$ where $\left\lfloor \cdot\right\rfloor $ is the floor operator$.$ This
$m_{n}$ satisfies Assumption \ref{A:mn}, if $d>3/2$. Other choices of sieve
terms such as $m_{n}=\left\lfloor n^{1/4}\right\rfloor +2$\ are also
considered in the simulations, and they yield similar results.
For the $C$ in $\varsigma_{n}=Cm_{n}\left[ \left( \log p_{n}\right)
^{1.1}+\left( \log m_{n}\right) ^{1.1}\right] $ or $\varsigma_{n}=C\left(
\log p_{n}\right) ^{1.1},$\ we test $C$ in the range of $0.5$ to $2.5$,
specifically, $0.5,0.6,...,2.5$. We determine the value of $C$ by minimizing
the following BIC:
\begin{equation}
\text{BIC}\left( C\right) =n\log\left[ \text{RSS}\left( C\right)
/n\right] +\left( \text{number of selected variables}\right) \cdot
\log\left( n\right) , \label{EQ:BIC}
\end{equation}
where RSS$\left( C\right) $ denotes the residual sum of squares from the
post-OCMT ordinary least squares regression of $Y$ on $P^{m_{n}}\left[
\boldsymbol{Z}_{\left( \hat{k}_{s}\right) }\right] $, where $\boldsymbol{Z}
_{\left( \hat{k}_{s}\right) }$ are the variables selected by the OCMT with
$C$ in use. The tuning parameters for the adaptive group Lasso\ are selected
as in Section \ref{SEC:tuningForLasso}. Another popular way to choose $C$ is
cross validation. As seen from the results, the BIC works well for our
procedure. In light of this, we will not pursue the procedure with cross validation.
For an easy reference, we present the implementation details in Appendix
\ref{SEC:procedure}.
\subsection{Estimators Compared\label{SEC:tuningForLasso}}
For comparison, we consider the adaptive group Lasso by \cite{HuangEtal2010},
which is designed for component selection in the nonparametric additive model.
It is a two-step approach---the first step is the usual group Lasso and the
second is the adaptive group Lasso with initial estimates from the first step.
We select tuning parameters $\lambda_{n1}$ in step 1 and $\lambda_{n2}$ in
step 2 ($\lambda_{n1}$ and $\lambda_{n2}$ are in the notations of
\cite{HuangEtal2010}) by BIC as well. Specifically, we set $\lambda
_{n1}=\lambda_{j}$ with
\[
\lambda_{j}=\exp\left\{ \log\left( \lambda_{\max}\right) +\left[
\log\left( \lambda_{\min}\right) -\log\left( \lambda_{\max}\right)
\right] \frac{j}{30}\right\}
\]
for $j=0,1,...,30,$ where $\lambda_{\min}=\max\left\{ 0.05,10^{-5}\left\Vert
\boldsymbol{y}\right\Vert \right\} $ and $\lambda_{\max}=0.5\left\Vert
\boldsymbol{y}\right\Vert .$ We calculate BIC$_{j}$ based on the estimation
using $\lambda_{n1}=\lambda_{j},$ and we select the $\hat{\lambda}_{n1}$ that
minimizes the BIC$_{j}$ among $j=0,1,...,30.$ We set the estimates in the
first step as the estimates using $\hat{\lambda}_{n1}.$ For the second step,
we set $\lambda_{n2}=\lambda_{j}$ for $j=0,1,...,30$ with the initial
estimates as the estimates from the first step using $\hat{\lambda}_{n1}$. We
again calculate BIC$_{j}$ based on the estimation using $\lambda_{n2}
=\lambda_{j},$ and we select the $\hat{\lambda}_{n2}$ that minimizes the
BIC$_{j}$ among $j=0,1,...,30.$ The final estimates are the adaptive group
Lasso estimates using $\hat{\lambda}_{n2}$. Note that the variables not
selected in the first step are not included in the second step, because the
penalty for those variables is infinity in the second step. We find the
solutions in both steps through the block coordinate descent algorithm (i.e.,
the \textquotedblleft shooting\textquotedblright\ algorithm).\ For the
algorithm details, see \cite{WuLange2008}.\footnote{One implementation in
MATLAB can be found at:
https://publish.illinois.edu/xiaohuichen/code/group-lasso-shooting/.}
In addition, we compare our method with the random forest regression and
bagging, both of which are commonly-used machine learning methods. Since there
is no variable selection criterion in random forest, we focus on the
comparison of out-sample forecasting performance.
\subsection{Estimation Results}
We consider the combinations of $n=200$ or $400$ and $p_{n}=100$, $200$, or
$1,000$ for each DGP. All results are based on 1,000 replications. We report
the results for five different methods. The first and second methods are our
post-OCMT procedure (Steps 1--6 in Appendix \ref{SEC:procedure}, denoted as
\textquotedblleft POST--OCMT\textquotedblright) and OCMT procedure (Steps 1--
5 in Appendix \ref{SEC:procedure}, denoted as \textquotedblleft
OCMT\textquotedblright), respectively. The third method is the adaptive group
Lasso (denoted as \textquotedblleft AGLASSO\textquotedblright), and the last
two methods are bagging (denoted as \textquotedblleft
BAGGING\textquotedblright) and random forest (denoted as \textquotedblleft
RF\textquotedblright), respectively.
Following the literature, we report the mean number of variables selected
(NV), the true positive rates (TPR), the false positive rates (FPR), the false
discovery rates (FDR), the percentage of correct selection (CS) for the first
three methods\footnote{That is, we precisely uncover the true model in
(\ref{EQ:model1b}) using all the true signals, but not any pseudo-signals or
noise variables.}, and the out-of-sample root mean squared forecast errors
(RMSFE) for all the five methods.\footnote{The forecasts are based on the
post-selection estimates from each method. In each replication, we
additionally independently generate 200 observations. Then, we calculate the
root mean square errors (RMSE) of the difference between the forecasts and $Y$
for the new observations. RMSFE is the average of those RMSE for 1,000
replications.} We report the average number of stages (STEP) for our OCMT
procedure only.
To save space, we report the results in the main body of this paper for only
DGPs 1 and 6 in Tables \ref{tabledgp1} and \ref{tabledgp6}, respectively.
These two designs correspond to the simplest and the most complicated designs
in the simulation, respectively. Results for DGPs 2--5 and 7--10 are provided
in Tables \ref{tabledgp2} to \ref{tabledgp10} in Appendix \ref{APP:tables}. To
showcase the advantage of the post-OCMT procedure compared to the OCMT
procedure only, we also consider a linear DGP (DGP 11) in the online
Appendix\textbf{ }\ref{APP:tables}\textbf{ }with results reported in Table
\ref{tabledgp11}.
\begin{table}[ptb]
\caption{DGP 1}
\label{tabledgp1}
\centering{}\centering{ } \resizebox{!}{0.67\textwidth}{
\begin{tabular}
[c]{l|ccccccc}\hline\hline
\multicolumn{8}{c}{Panel 1: $n=200$, $p=100$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 4.0080 & 0.9998 & 0.0001 & 0.0015 & 0.9910 & - & 1.0718\\
OCMT & 4.0180 & 0.9998 & 0.0002 & 0.0031 & 0.9820 & 1.0260 & 1.0886\\
AGLASSO & 4.0410 & 0.9990 & 0.0005 & 0.0075 & 0.9550 & - & 1.0753\\
BAGGING & - & - & - & - & - & - & 1.4403 \\
RF & - & - & - & - & - & - & 1.4846 \\\hline
\multicolumn{8}{c}{Panel 2: $n=200$, $p=200$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 4.0150 & 0.9998 & 0.0001 & 0.0027 & 0.9830 & - & 1.0740\\
OCMT & 4.0200 & 0.9998 & 0.0001 & 0.0035 & 0.9790 & 1.0390 & 1.0939\\
AGLASSO & 4.0510 & 0.9950 & 0.0004 & 0.0114 & 0.9330 & - & 1.0808\\
BAGGING & - & - & - & - & - & - & 1.4690 \\
RF & - & - & - & - & - & - & 1.5160 \\\hline
\multicolumn{8}{c}{Panel 3: $n=200$, $p=1000$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 4.0010 & 0.9952 & 0.0000 & 0.0034 & 0.9640 & - & 1.0773\\
OCMT & 4.0200 & 0.9952 & 0.0000 & 0.0062 & 0.9510 & 1.0880 & 1.1173\\
AGLASSO & 3.9370 & 0.9620 & 0.0001 & 0.0145 & 0.8770 & - & 1.0993\\
BAGGING & - & - & - & - & - & - & 1.5281 \\
RF & - & - & - & - & - & - & 1.5847 \\\hline
\multicolumn{8}{c}{Panel 4: $n=400$, $p=100$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 4.0010 & 1.0000 & 0.0000 & 0.0002 & 0.9990 & - & 1.0470\\
OCMT & 4.0010 & 1.0000 & 0.0000 & 0.0002 & 0.9990 & 1.0000 & 1.0470\\
AGLASSO & 4.0080 & 0.9990 & 0.0001 & 0.0020 & 0.9870 & - & 1.0531\\
BAGGING & - & - & - & - & - & - & 1.3234 \\
RF & - & - & - & - & - & - & 1.3648 \\\hline
\multicolumn{8}{c}{Panel 5: $n=400$, $p=200$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 4.0010 & 1.0000 & 0.0000 & 0.0002 & 0.9990 & - & 1.0477\\
OCMT & 4.0010 & 1.0000 & 0.0000 & 0.0002 & 0.9990 & 1.0000 & 1.0477\\
AGLASSO & 4.0210 & 1.0000 & 0.0001 & 0.0035 & 0.9790 & - & 1.0520\\
BAGGING & - & - & - & - & - & - & 1.3536 \\
RF & - & - & - & - & - & - & 1.3984 \\\hline
\multicolumn{8}{c}{Panel 6: $n=400$, $p=1000$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 4.0040 & 1.0000 & 0.0000 & 0.0007 & 0.9960 & - & 1.0519\\
OCMT & 4.0040 & 1.0000 & 0.0000 & 0.0007 & 0.9960 & 1.0000 & 1.0519\\
AGALSSO & 4.0313 & 1.0000 & 0.0000 & 0.0052 & 0.9688 & - &
1.0524\\
BAGGING & - & - & - & - & - & - & 1.4120 \\
RF & - & - & - & - & - & - & 1.4574 \\\hline\hline
\end{tabular}}\end{table}{ }
\begin{table}[ptb]
\caption{DGP 6}
\label{tabledgp6}
\centering{}\centering{ } \resizebox{!}{0.67\textwidth}{
\begin{tabular}
[c]{l|ccccccc}\hline\hline
\multicolumn{8}{c}{Panel 1: $n=200$, $p=100$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 3.6300 & 0.6628 & 0.0033 & 0.0650 & 0.1770 & - & 1.2296\\
OCMT & 4.8450 & 0.6834 & 0.0149 & 0.2163 & 0.0000 & 1.2650 & 1.4003\\
AGLASSO & 2.9840 & 0.5522 & 0.0023 & 0.0432 & 0.0040 & - & 1.3670\\
BAGGING & - & - & - & - & - & - & 1.4522 \\
RF & - & - & - & - & - & - & 1.4645 \\\hline
\multicolumn{8}{c}{Panel 2: $n=200$, $p=200$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 3.5600 & 0.6438 & 0.0017 & 0.0682 & 0.1290 & - & 1.2396\\
OCMT & 4.6700 & 0.6642 & 0.0069 & 0.2107 & 0.0000 & 1.2220 & 1.3861\\
AGLASSO & 2.6340 & 0.4866 & 0.0010 & 0.0417 & 0.0000 & - & 1.3695\\
BAGGING & - & - & - & - & - & - & 1.4855 \\
RF & - & - & - & - & - & - & 1.4969 \\\hline
\multicolumn{8}{c}{Panel 3: $n=200$, $p=1000$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 3.2550 & 0.5892 & 0.0003 & 0.0669 & 0.0440 & - & 1.2674\\
OCMT & 4.2310 & 0.6088 & 0.0012 & 0.2040 & 0.0000 & 1.0880 & 1.3311\\
AGLASSO & 1.9780 & 0.3678 & 0.0001 & 0.0306 & 0.0000 & - & 1.4304\\
BAGGING & - & - & - & - & - & - & 1.5395 \\
RF & - & - & - & - & - & - & 1.5539 \\\hline
\multicolumn{8}{c}{Panel 4: $n=400$, $p=100$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 4.8840 & 0.9448 & 0.0017 & 0.0246 & 0.6960 & - & 1.0659\\
OCMT & 7.4660 & 0.9534 & 0.0281 & 0.3056 & 0.0000 & 1.8380 & 1.4894\\
AGLASSO & 4.3480 & 0.8276 & 0.0022 & 0.0331 & 0.1060 & - & 1.1294\\
BAGGING & - & - & - & - & - & - & 1.3502 \\
RF & - & - & - & - & - & - & 1.3767 \\\hline
\multicolumn{8}{c}{Panel 5: $n=400$, $p=200$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 4.8870 & 0.9488 & 0.0007 & 0.0213 & 0.7080 & - & 1.0723\\
OCMT & 7.4120 & 0.9542 & 0.0135 & 0.3019 & 0.0000 & 1.8360 & 1.5018\\
AGLASSO & 4.2240 & 0.8122 & 0.0008 & 0.0270 & 0.0630 & - & 1.1410\\
BAGGING & - & - & - & - & - & - & 1.3848 \\
RF & - & - & - & - & - & - & 1.4118 \\\hline
\multicolumn{8}{c}{Panel 6: $n=400$, $p=1000$}\\\hline
& NV & TPR & FPR & FDR & CS & STEP & RMSFE\\\hline
POST-OCMT & 4.8380 & 0.9356 & 0.0002 & 0.0242 & 0.6740 & - & 1.0784\\
OCMT & 7.4140 & 0.9398 & 0.0027 & 0.3082 & 0.0000 & 1.8030 & 1.4863\\
AGLASSO & 4.1340 & 0.7916 & 0.0002 & 0.0299 & 0.0110 & - &
1.1509\\
BAGGING & - & - & - & - & - & - & 1.4477 \\
RF & - & - & - & - & - & - & 1.4689 \\\hline\hline
\end{tabular}}\end{table}
We can make several observations. First, POST--OCMT performs the best among
the three methods in most cases. POST--OCMT outperforms OCMT due to its
ability to eliminate pseudo-signals or noise variables.\ POST--OCMT even
outperforms OCMT slightly for DGPs without any pseudo-signals, especially for
$n=200$. This result confirms the necessity of conducting the additional
post-OCMT step. Moreover, the additional computation time for POST--OCMT
compared with OCMT can be almost ignored because the number of candidate
variables after OCMT is very small. Second, AGLASSO performs almost the same
as our procedures for DGP 2. Note that DGP 2 contains four independent
signals, two pseudo-signals, and no hidden signals. Then, this result is not
surprising, because Lasso performs very well for independent signals, and the
biggest challenge for our procedure is the possible presence of
pseudo-signals. When the signals are correlated in DGP 4, AGLASSO performs
less well and is outperformed by POST--OCMT. Third, the performance of all
methods improves as $n$ increases from $200$ to $400$. Fourth, OCMT is very
successful at picking up hidden signals, especially for $n=400$; see, for
example, the results for DGPs 3, 4, 6, and 8. Fifth, our methods perform well
in the presence of heteroskedastic errors for DGP 7 and 8. Sixth, the CS of
POST--OCMT performs well even for the complicated DGPs 6 and 8 at $n=400$,
whereas AGLASSO performs poorly in terms of CS for these two DGPs. Seventh,
the number of stages for OCMT basically confirms our theoretical findings. For
example, for DGP 1, the mean number of stages is slightly more than 1 for
$n=200$, and is 1 for $n=400$. In the presence of hidden signals, the mean
number of stages is approximately 2 for $n=400$ for DGPs\ 3, 4, 6, 8, and 10.
Note the mean number of stages is approximately 2 for $n=400$ for DGP\ 5 with
correlated signals and no hidden signals. The reason is that the correlation
makes the net effect of $X_{1}$ on $Y$ very small and $X_{1}$ almost behaves
like a \textquotedblleft hidden signal\textquotedblright. This result also
confirms the necessity of using multiple stages instead of a single-stage
procedure. Eighth, our method also works well for DGPs 9 and 10 that mimic the
application. An additional remark is that our procedure can be implemented
fast and is much faster than AGLASSO. For example, when $p=1,000,$ our
procedure took less than half a minute for one replication on average, whereas
the AGLASSO took hours for one replication. Finally, the first three
procedures have smaller RMSFE than Bagging and RF in almost all scenarios and
thus have better out-of-sample forecasting performance.\footnote{Note that
Bagging has slightly better out-of-sample forecasting performance than RF.
This is reasonable given that our DGPs have a finite number of signal
variables. In RF, many decision splits do not improve predictive accuracy
because they rely solely on noise variables to generate the trees.}\textbf{
}This is because they explicitly utilize the information (in the form of an
additive function) underlying the DGPs while BAGGING and RF do not. In
addition, POST--OCMT delivers the best out-of-sample forecasting performance
among all procedures.
To summarize, our methods perform well in small samples, and we view it as a
useful alternative to existing methods in the literature.
\section{An Application\label{SEC:application}}
In this section, we apply our method to a dataset extracted from
RUMiC.\footnote{RUMiC consists of three parts: the Urban Household Survey, the
Rural Household Survey, and the Migrant Household Survey. A group of
researchers at the Australian National University, the University of
Queensland, and the Beijing Normal University initiated this survey. The
Institute for the Study of Labor (IZA) supported it and provides the
Scientific Use Files. RUMiC had financial support from the Australian Research
Council, the Australian Agency for International Development, the Ford
Foundation, IZA, and the Chinese Foundation of Social Sciences. More
information on the survey can be found at
https://datasets.iza.org/dataset/58/longitudinal-survey-on-rural-urban-migration-in-china.}
The survey studied immigrants or workers moving from the rural areas of China
to its big cities. The survey asked interviewees (immigrants or workers) a
wide range of questions. For the detailed design of the survey and other
information, including on the construction of each variable, see the survey
website and \cite{Gongetal2008}. Currently, the 2008 wave data are publicly available.
Economic reforms since the late 1970s have brought significant changes to
China's economy. The government began relaxing its policy on population
mobility in the early 1980s. Gradually, peasants were allowed to leave
villages and work in big cities to earn higher incomes. Most migrant workers
may leave their spouses, children, or parents behind in their hometowns\, who
may need their financial support. This situation results in monetary
transfers, that is, remittances, from migrant workers to their family. In the
context of migration, family, and economic development, remittances are not
only an income source for recipients, but also reflect intrafamilial
relationships. Remittances clearly represent a dimension of family ties and
demonstrate high degrees of interaction between migrants and families at home.
In addition, remittances from the rural migrant workers also contribute
significantly to China's agricultural productivity (c.f.,
\cite{Rozelle_at_al.1999}). For these reasons, it has long been of interest to
model remittances to families or relatives in the hometown; see \cite{Li2001}
and \cite{Cai2003}, among others.
In this application, we take remittance as the dependent variable ($Y$); we
focus on the dataset from Guangdong Province (Guangzhou, Dongguan, and
Shenzhen cities) in the 2008 survey wave, and keep 78 covariates\ from the
dataset.\footnote{We keep covariates with relatively fewer missing
observations.} After dropping observations with missing information, the
number of observations is 456. We provide the definitions of the dependent
variable and covariates, and the associated summary statistics, in Tables
\ref{table:app_def} and \ref{table:summary}, respectively. We report the
original labels of all covariates in the survey in the first column of Table
\ref{table:app_def} for reference. Among the 78 covariates, there are some
continuous variables with most observations as 0 (Panel C in Table
\ref{table:app_def}), and some discrete variables with very limited support
(Panel D in Table \ref{table:app_def}). For those variables, the design matrix
of the sieves generated are either singular or close to singular. For this
reason, we add those variables linearly into the model and treat them the same
as\ dummy variables (Panel E in Table \ref{table:app_def}) for modeling.
Consequently, the way we fit the dataset resembles the approach we used for
DGPs 9 and 10 in the simulation. We explain the reason our theoretical results
continue to hold in this situation in the simulation section (DGP 9). We take
natural logarithms for the $Y$ and continuous $X$ variables to offset the
effect of outliers; otherwise, the forecast can easily take some\ extreme values.
We randomly select 400 observations as the training sample, and the remaining
56 observations as the test sample. The number of sieve terms is set as
$m_{n}=\left\lfloor 400^{1/4}\right\rfloor +1=5$ for the continuous variables
in Panel B of Table \ref{table:app_def}). We set $\varsigma_{n}=Cm_{n}\left[
\left( \log p_{n}\right) ^{1.1}+\left( \log m_{n}\right) ^{1.1}\right] $
and $\varsigma_{n}^{\ast}=3\varsigma_{n}$ for continuous variables, and
$\varsigma_{n}=C\left( \log p_{n}\right) ^{1.1}$ and $\varsigma_{n}^{\ast
}=4\varsigma_{n}$ for terms entering the model linearly (see Section
\ref{SEC:tuning} for the reason). We set $C$ in the range of $0.5$ to $2.5$,
specifically, $0.5,0.6,...,2.5$, and choose $C$ to minimize the BIC for the
model selection$.$ The competing methods are the group Lasso (labelled as
GLASSO) and the adaptive group Lasso (AGLASSO). The tuning parameters for
GLASSO and AGLASSO\ are selected as in Section \ref{SEC:tuningForLasso}. We
evaluate the performance of all methods based on the RMSFE of the test dataset
using the fitted models from different methods. We independently repeat the
above procedure 100 times.
\begin{table}[tbh]
\caption{Performance in Terms of RMSFE (100 Cases), benchmark: AGLASSO}
\label{table:ForecaseRMSE}
\centering{}\centering { }
\begin{tabular}
[c]{ccccccccccc}\hline\hline
\multicolumn{3}{c}{OCMT} & & \multicolumn{3}{c}{POST-OCMT} & &
\multicolumn{3}{c}{GLASSO}\\\cline{1-3}\cline{5-7}\cline{9-11}
Better & Same & Worse & & Better & Same & Worse & & Better & Same & Worse\\
77 & 17 & 6 & & 76 & 21 & 3 & & 5 & 88 & 7\\\hline\hline
\end{tabular}
\end{table}
\begin{table}[tbh]
\caption{Average RMSFE ratio, benchmark: AGLASSO}
\label{table:RMSEratio}
\centering{}\centering { }
\begin{tabular}
[c]{ccccc}\hline\hline
One Stage & OCMT & POST-OCMT & GLASSO & AGLASSO\\\hline
0.873 & 0.817 & 0.814 & 1.005 & 1\\\hline\hline
\end{tabular}
\end{table}
We report the results in terms of the out-of-sample RMSFE in Tables
\ref{table:ForecaseRMSE} and \ref{table:RMSEratio}, with AGLASSO as the
benchmark. Table \ref{table:ForecaseRMSE} shows that OCMT and POST-OCMT
outperform AGLASSO in the majority of cases. Of course, our methods do not
outperform AGLASSO all the time. To highlight the necessity of the multiple
stages, we also report the results of the one-stage procedure (we selected the
tuning parameter also by minimizing the BIC) in Table \ref{table:RMSEratio}.
The OCMT stops at the second stage for all 100 cases. We normalize the average
RMSFE of AGLASSO as 1.\ Notably, the average RMSFE of our methods are lower
than that of AGLASSO. The OCMT also improves the RMSFE over the one-stage
procedure, possibly owing to some hidden signals uncovered by our
multiple-stage procedure. We report the frequencies of variables (out of 100
cases) selected by all methods in Table \ref{table:varselected}. It appears
that G102 (monthly income) along with G133 (gifts to others, including
parents) and G137 (education cost for left-behind children) contribute most to
the model, as shown by all methods in general. We note that AGLASSO tends to
select more variables than our methods, on average, which was also CKP's
finding in their application. The \textquotedblleft
over-fitting\textquotedblright\ is the main cause of the relative inferior
performance of AGLASSO.
\section{Conclusion\label{SEC:conclusion}}
In this paper, we examine the one-covariate-at-a-time multiple testing
approach to model selection in additive models. The properties of the TPR,
FPR, and FDR of our approach are established based on some asymptotic
probability bounds of Type-I and II errors. The simulation experiments and one
application on the RUMiC dataset showcase excellent small-sample properties of
our methods. Just as stated by CKP for linear models, we view our approach as
a useful alternative to the model selection methods for additive models in the literature.