WGMwgm QEqe EPep PMSpms BECbec DEde rmkRemark assumptionAssumption
zproofofProof of

[]

[]

[]

# Semiparametric inference for partially linear regressions with Box-Cox
transformation

Daniel Becker  PhD-Student in Economics, University of Bonn, Germany,
Adenauerallee 24-26, 53113 Bonn(e-mail<EMAIL_ADDRESS> Alois Kneip
Statistic Professor, University of Bonn, Germany, Adenauerallee 24-26, 53113
Bonn(e-mail<EMAIL_ADDRESS> Valentin Patilea  Statistic Professor,
CREST(Ensai), France, Campus de Ker-Lann, Rue Blaise Pascal - BP 37203
(e-mail<EMAIL_ADDRESS>
###### Abstract

In this paper, a semiparametric partially linear model in the spirit of
Robinson (1988) with Box-Cox transformed dependent variable is studied.
Transformation regression models are widely used in applied econometrics to
avoid misspecification. In addition, a partially linear semiparametric model
is an intermediate strategy that tries to balance advantages and disadvantages
of a fully parametric model and nonparametric models. A combination of
transformation and partially linear semiparametric model is, thus, a natural
strategy. The model parameters are estimated by a semiparametric extension of
the so called smooth minimum distance (SmoothMD) approach proposed by Lavergne
and Patilea (2013). SmoothMD is suitable for models defined by conditional
moment conditions and allows the variance of the error terms to depend on the
covariates. In addition, here we allow for infinite-dimension nuisance
parameters. The asymptotic behavior of the new SmoothMD estimator is studied
under general conditions and new inference methods are proposed. A simulation
experiment illustrates the performance of the methods for finite samples.

###### keywords:

Semiparametric partially linear model, Nonparametric kernel estimators, Root
N-consistent estimation, Conditional estimating equations, Hypothesis testing

##  1 Introduction

The data consists of independent copies of a response variable $Y$ and a
random covariate vector
$\left(\boldsymbol{X}^{T},\boldsymbol{Z}^{T}\right)^{T}\in\mathbb{R}^{p}\times\mathbb{R}^{q}$.111Herein,
vectors are column matrices and for any matrix $\boldsymbol{A}$,
$\boldsymbol{A}^{T}$ denotes its transpose. We adopt a semiparametric approach
to model the dependence of $Y$ on $X$ and $Z$. Starting with Robinson (1988),
the use of partial linear models of the form
$Y=\boldsymbol{X}^{T}\boldsymbol{\beta}+m(\boldsymbol{Z})+\varepsilon$ for
some unknown, nonparametric function $m$ has become popular in this context.
For an overview consider Härdle et al. (2000) and Li and Racine (2007).

In many important economic applications, the response variable $Y$ is
positive, i.e. $P(Y>0)=1$. It is then often questionable to assume (partial)
linear regression models. A commonly used remedy is to apply a partial linear
model to a suitable transformation of the response variable. In econometric
practice the dependent variable is then frequently $\log$-transformed (see for
example Acemoglu et al. (2001) and Autor and Handel (2013)). But usually no
substantial knowledge exists ensuring that this specific transformation leads
to the correct model. The transformation proposed by Box and Cox (1964)

$\displaystyle
T(Y,\lambda)=\begin{cases}\frac{Y^{\lambda}-1}{\lambda}&,\lambda\neq 0\\\
\log(Y)&,\lambda=0\end{cases}$

offers much more flexibility. Specifying the transformation up to a parameter
and estimating the parameter together with the regression coefficients, leads
to more reliable results at the cost of having to estimate only one additional
parameter. In addition, the common $\log$-transformation is nested in the Box-
Cox transformation and, thus, can be confirmed by a statistical test.

This motivates the approach adopted in this paper. To model the relationship
between a positive response variable and the covariate vectors, we consider a
_transformation partially linear_ mean regression model given by

$T(Y,\lambda)=\boldsymbol{X}^{T}\boldsymbol{\beta}+m(\boldsymbol{Z})+\varepsilon,$
(1)

where $m(\cdot)$ is an unknown function and

$\displaystyle E[\varepsilon\mid\boldsymbol{X},\boldsymbol{Z}]=0.$ (2)

The true values $\lambda_{0}$ and $\boldsymbol{\beta}_{0}$ of the parameters
$\lambda>0$ and $\boldsymbol{\beta}\in\mathbb{R}^{p}$ are unknown and have to
be estimated from an i.i.d. sample
$(Y_{i},\boldsymbol{X}_{i}^{T},\boldsymbol{Z}_{i}^{T})$, $i=1,\dots,n$, of
$(Y,\boldsymbol{X}^{T},\boldsymbol{Z}^{T})$. We impose no further assumption
on the conditional distribution of $\varepsilon$. In particular, we allow for
heteroscedasticity of unknown form. The vector $\boldsymbol{Z}$ contains only
continuous variables, but the components of $\boldsymbol{X}$ need not be
continuous.

We wish to note that the Box-Cox transformation is widely used in
applications, and is discussed in various textbooks, e.g. Amemiya (1985);
Greene (2003); Horowitz (1998); Showalter (1994); Wooldridge (1992).
Furthermore, there exist several empirical studies that employ the Box-Cox
transformation. See, for instance, Berndt et al. (1993); Heckman and Polachek
(1974) or Keane et al. (1988). For an overview of the Box-Cox transformation,
consider Horowitz (1998) and Sakia (1992). The semiparametric partially linear
specification of the conditional mean of the response seems to be quite
appealing as it allows a linear dependency on a subvector $\boldsymbol{X}$ of
covariates, which could include discrete variables, and meanwhile allows a
nonparametric additive effect of the covariates $\boldsymbol{Z}$. These
features could help practitioners faced with a large cross-sectional data set
with independent observations including many candidate explanatory variables,
who, on the basis of economic theory or past experience with similar data,
feel able to parameterize only some of them.

Similar to other semiparametric models, the major challenge is to develop
$\sqrt{n}$-consistent estimators and corresponding inference procedures for
the parameters $(\lambda_{0},\boldsymbol{\beta}_{0}^{T})$. This requires a
careful methodological development, since to our knowledge, there is no
established procedure which can readily be applied under our general setup.
Despite its popularity, even in a purely parametric framework, estimation and
inference in a Box-Cox transformation model is a difficult statistical problem
and is usually based on quite restrictive assumptions on the conditional law
of the response. See for instance chapter 5 of Horowitz (1998) for an
illuminating discussion. In particular, ordinary least squares estimation of
$(\lambda_{0},\boldsymbol{\beta}_{0}^{T})$ may lead to inconsistent results,
and more sophisticated procedures, such as the nonlinear two-stage least
squares (NL2SLS) estimator introduced by Amemiya and Powell (1981), have to be
applied. The problem becomes even more complex in the case of the
semiparametric regression (1) where one only assumes the minimal
identification condition (2).

In order to motivate our procedure, let us first consider the special case
that the true value $\lambda_{0}$ is known a priori. With
$Y_{0}=T(Y,\lambda_{0})$, we then arrive at a standard partial linear model
$Y_{0}=\boldsymbol{X}^{T}\boldsymbol{\beta}+m(\boldsymbol{Z})+\varepsilon$.
Consequently,
$m(Z)=E(Y_{0}|Z)-E[\boldsymbol{X}\mid\boldsymbol{Z}])^{T}\boldsymbol{\beta}_{0}$
and

$\displaystyle
0=E[\varepsilon|\boldsymbol{X},\boldsymbol{Z}]=E\bigg{[}Y_{0}-E[Y_{0}\mid\boldsymbol{Z}]-(\boldsymbol{X}-E[\boldsymbol{X}\mid\boldsymbol{Z}])^{T}\boldsymbol{\beta}_{0}\bigg{|}X,Z\bigg{]}.$

The basic idea of Robinson (1988) now consists in disentangling the
nonparametric estimation of unknown functions and the parametric estimation of
the coefficient vector $\boldsymbol{\beta}$. In a first step Nadaraya-Watson
kernel estimators are used for nonparametric estimation of the functions
$E[Y_{0}\mid\boldsymbol{Z}]$ and $E[\boldsymbol{X}\mid\boldsymbol{Z}]$.
Plugging in these nonparametric function estimates, an OLS regression of
$Y_{0}-E[Y_{0}\mid\boldsymbol{Z}]$ on
$\boldsymbol{X}-E[\boldsymbol{X}\mid\boldsymbol{Z}]$ then leads to an
estimator $\widehat{\boldsymbol{\beta}}$. If $q<4$ and standard nonparametric
estimators based on second order kernels are applied, then, assuming
homoscedasticity and suitable bandwidth sequences, Robinson showed that
$\widehat{\boldsymbol{\beta}}$ is a $\sqrt{n}$-consistent, asymptotically
normally distributed and efficient estimator of $\boldsymbol{\beta}_{0}$. If
the dimensionality of $Z$ is larger, i.e. $q\geq 4$, then
$\sqrt{n}$-consistent estimators can only be achieved by using higher order
kernels.

The quite straightforward way to build efficient estimators made the partially
linear model quite a popular. Versions of this model have also been studied by
Engle et al. (1986); Heckman (1986); Shiller (1984) and Wahba (1984). In order
to avoid the trimming introduced by Robinson (1988) to ensure that the
estimate of the density of $\boldsymbol{Z}$, $f_{z}(\boldsymbol{Z})$, stays
away from zero, Li (1996) considered as starting point the unfeasible OLS
regression of $(Y_{0}-E[Y_{0}\mid\boldsymbol{Z}])f_{z}(\boldsymbol{Z})$ on
$(\boldsymbol{X}-E[\boldsymbol{X}\mid\boldsymbol{Z}])f_{z}(\boldsymbol{Z})$.
Premultiplying by the density of $\boldsymbol{Z}$ does not break the
consistency of the unfeasible OLS estimator since
$E[f_{z}(\boldsymbol{Z})\varepsilon\mid\boldsymbol{X},\boldsymbol{Z}]=f_{z}(\boldsymbol{Z})E[\varepsilon\mid\boldsymbol{X},\boldsymbol{Z}]=0.$
Next, Li (1996) proposed to build OLS estimates using standard kernel
estimators instead of the unfeasible response and covariates. This new
estimator is still $\sqrt{n}$-consistent and asymptotically normally
distributed. Moreover, Li (1996) relaxed the condition on the bandwidth with
the consequence that the smoothing requires higher order kernels only if the
dimension of $\boldsymbol{Z}$ is larger than 5, instead of larger than 3 as
required in Robinson (1988).

Let us now return to the general model (1) with unknown parameter $\lambda$.
Adopting Li’s idea of premultiplying with the density $f_{z}$ of
$\boldsymbol{Z}$, the conditional moment condition
$E[\varepsilon|\boldsymbol{X},\boldsymbol{Z}]=0$ leads to

$E\left(\bigg{(}T(Y,\lambda)-E[T(Y,\lambda)\mid\boldsymbol{Z}]-(\boldsymbol{X}-E[\boldsymbol{X}\mid\boldsymbol{Z}])^{T}\boldsymbol{\beta}\bigg{)}f_{z}(\boldsymbol{Z})\bigg{|}X,Z\right)=0\quad\Longleftrightarrow\quad\lambda=\lambda_{0},\boldsymbol{\beta}=\boldsymbol{\beta}_{0}.$
(3)

It is now immediately seen that the unknown, additional parameter $\lambda$
introduces a major complication. Unlike the standard partial linear model,
there is no way to disentangle nonparametric estimation of unknown functions
and parametric estimation of the coefficient vector
$(\lambda_{0},\boldsymbol{\beta}_{0}^{T})$. The reason is that
$E[T(Y,\lambda)\mid\boldsymbol{Z}]$ depends on $\lambda$. Indeed, there does
not seem to exist a straightforward generalization of Robinson’s approach
which is able to cope with model (1).

On the other hand, by (3), our model belongs to the large class of models
identified by conditional moment restrictions. Methodologically however, we
have to deal with the obvious facts that a) the model is highly nonlinear in
$\lambda$ and b) (3) incorporates an infinite dimensional nuisance parameter
$\boldsymbol{\eta}_{\lambda}$ consisting of the functions $f_{z}(z)$,
$E[T(Y,\lambda)\mid\boldsymbol{Z}=z]$ and
$E[\boldsymbol{X}\mid\boldsymbol{Z}=z]$. Even if $\boldsymbol{\eta}_{\lambda}$
were known a priori, any use of the generalized method of moments (GMM) runs
into the problem that the conditional moment restrictions identifying our
model imply an infinite number of unconditional moment restrictions, since the
conditioning variables have a support with infinite cardinality. But GMM
relies only on a finite number of instruments and, thus, in general,
consistency of GMM requires regularization and additional assumptions. See
Dominguez and Lobato (2004). This problem has already been pointed out for the
Box-Cox transformation by Foster et al. (2001) and Shin (2008) in the linear
case. See also Horowitz (1998).

More recent work explicitly focuses on regularization techniques in order to
account for complex conditional moment conditions. Some methods rely on
increasing the number of considered unconditional estimating equations (or
instruments) with the sample size, such as the sieve minimum distance (SMD)
approach of Ai and Chen (2003), or generalizations of GMM and empirical
likelihood (EL) by Donald et al. (2003) and Hjort et al. (2009). Carrasco and
Florens (2000) use a regularization approach to generalize the GMM approach to
a continuum of estimating equations. Other EL-type estimators use
nonparametric smoothing to estimate conditional equations, such as Antoine et
al. (2007), Kitamura et al. (2004), and Smith (2007a, b). All these approaches
share one common feature. The estimators’ sensitivity to the user-chosen
parameter (number of estimating equations, regularization parameter, or
smoothing parameter) remains largely unknown.

In this paper, we rely upon the SmoothMD approach proposed by Lavergne and
Patilea (2013) to estimate $(\lambda_{0},\boldsymbol{\beta}_{0}^{T})$. Roughly
speaking, SmoothMD can be seen as a new technique to translate conditional
moment conditions into unconditional ones which can be approximated by sample
averages. Although the method involves some tuning parameters, an attractive
feature consists of the fact that a practical choice is quite uncritical,
since asymptotic results can be established for a wide range of possible
values of these tuning parameters (including values independent of the sample
size). SmoothMD thus bridges a gap between Dominguez and Lobato’s method,
which does not require a user-chosen parameter, and the competing SMD
estimator and EL and GMM-type methods that rely on smoothing with restrictive
conditions on the choice of smoothing parameters. Furthermore, although
Lavergne and Patilea (2013) rely on a more standard setup. We will show that
this technique can be well adapted to deal with complex functional nuisance
parameters.

In order to explain the methodology, we introduce some abbreviations in order
to simplify the lengthy expressions in (3). Let
$\boldsymbol{W}=\left(\boldsymbol{X}^{T},\boldsymbol{Z}^{T}\right)^{T}\in\mathbb{R}^{p}\times\mathbb{R}^{q}$
and $\boldsymbol{U}=\left(Y,\boldsymbol{W}^{T}\right)^{T}$. Moreover, set
$\boldsymbol{\theta}=\left(\lambda,\boldsymbol{\beta}^{T}\right)^{T}$, and for
a real value $\gamma$ define

$g(\boldsymbol{U};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda})=\left(T(Y,\lambda)-E[T(Y,\lambda)\mid\boldsymbol{Z}]-(\boldsymbol{X}-E[\boldsymbol{X}\mid\boldsymbol{Z}])^{T}\boldsymbol{\beta}\right)f_{z}(\boldsymbol{Z}).$
(4)

Recall that $\boldsymbol{\eta}_{\lambda}$ is an infinite-dimensional nuisance
parameter defined by

$\boldsymbol{\eta}_{\lambda}=\boldsymbol{\eta}_{\lambda}(z)=(f_{z}(z),E[T(Y,\lambda)\mid\boldsymbol{Z}=z],E[\boldsymbol{X}\mid\boldsymbol{Z}=z]^{T}).$

Condition (3) is then equivalent to requiring

$E\left(g(\boldsymbol{U};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda})\big{|}W\right)=0\quad\Longleftrightarrow\quad\boldsymbol{\theta}=\left(\lambda,\boldsymbol{\beta}^{T}\right)^{T}=\left(\lambda_{0},\boldsymbol{\beta}_{0}^{T}\right)^{T}=:\boldsymbol{\theta}_{0}.$
(5)

SmoothMD is based on the following insight: Let $\boldsymbol{U}_{1}$ and
$\boldsymbol{U}_{2}$ be two independent copies of $\boldsymbol{U}$ with
corresponding subvectors $\boldsymbol{W}_{1}$ and $\boldsymbol{W}_{2}$. For
any symmetric function $\omega(\cdot)$ of $\boldsymbol{W}$ with positive
Fourier transform, we than have that

$E\left(g(\boldsymbol{U};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda})\big{|}W\right)=0\quad\text{if
and only if}\quad
Q(\boldsymbol{\theta})=E[g(\boldsymbol{U}_{1};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda,1})g(\boldsymbol{U}_{2};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda,2})\omega(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})]=0.$

Here $\boldsymbol{\eta}_{\lambda,1}$ and $\boldsymbol{\eta}_{\lambda,2}$ are
the vector of nuisance functions corresponding to $\boldsymbol{Z}_{1}$ and
$\boldsymbol{Z}_{2}$. Indeed, in Lemma 2.1 of Section 2.1, it will be shown
that

$\displaystyle Q(\boldsymbol{\theta})=\begin{cases}\hskip 9.95863pt0&,\text{
if
}\boldsymbol{\theta}=\left(\lambda,\boldsymbol{\beta}^{T}\right)^{T}=\left(\lambda_{0},\boldsymbol{\beta_{0}}^{T}\right)^{T},\\\
>0&,\text{ else}.\end{cases}$ (6)

Lavergne and Patilea (2013) list several possible ways to define
$\omega(\cdot)$, but throughout this paper we will rely on the simple choice
$\omega(\boldsymbol{W}):=\exp\left\\{-\boldsymbol{W}^{T}\boldsymbol{D}\boldsymbol{W}\right\\}$,
where $\boldsymbol{D}$ is a diagonal matrix whose positive diagonal elements
$d_{1},\dots,d_{p+q}$ represent user selected tuning parameters. A sensible
choice consists of using the standard deviations of the components of the
vectors $(\boldsymbol{X}_{i}^{T},\boldsymbol{Z}_{i}^{T})^{T}$.

If the functions which define $\boldsymbol{\eta}_{\lambda}$ were known, then a
data-based estimator of the unconditional moment $Q(\boldsymbol{\theta})$
could be obtained by the sample averages

$Q(\boldsymbol{\theta})=\frac{1}{n^{2}}\sum_{1\leq i,j\leq
n}g(\boldsymbol{U}_{i};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda,i})g(\boldsymbol{U}_{j};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda,j})\omega(\boldsymbol{W}_{i}-\boldsymbol{W}_{j}),$
(7)

where $\boldsymbol{\eta}_{\lambda,i}$ and $\boldsymbol{\eta}_{\lambda,j}$ are
the nuisance parameter values corresponding to $\boldsymbol{Z}_{i}$ and
$\boldsymbol{Z}_{j}$, respectively. The average scheme proposed by Lavergne
and Patilea (2013) relies on leaving out diagonal elements with $i=j$. In our
setup, inclusion of these diagonal terms provides more stable and reliable
estimators. We will show that the resulting bias is asymptotically negligible.
Note that by definition of $\omega(\cdot)$

$E\left[Q_{n}\left(\boldsymbol{\theta}\right)\right]=Q(\boldsymbol{\theta})+\frac{1}{n}E\left(g(\boldsymbol{U};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda})^{2}\right)=Q(\boldsymbol{\theta})+O(n^{-1}).$

Under model (1), the nuisance parameter $\boldsymbol{\eta}_{\lambda}$ is
unknown, and $Q_{n}(\cdot)$ cannot be directly computed. We therefore use
kernel estimation to determine nonparametric estimators
$\widehat{\boldsymbol{\eta}}_{\lambda}$. This then leads to a feasible version
$\widehat{Q}_{n}(\cdot)$. More precisely, our estimation procedure can be
described as follows. For each $\lambda$, we define the map

$\boldsymbol{\beta}\mapsto\widehat{Q}_{n}\left(\left(\lambda,\boldsymbol{\beta}^{T}\right)^{T}\right)=\frac{1}{n^{2}}\sum_{1\leq
i,j\leq
n}g(\boldsymbol{U}_{i};\boldsymbol{\theta},\widehat{\boldsymbol{\eta}}_{\lambda,i})g(\boldsymbol{U}_{j};\boldsymbol{\theta},\widehat{\boldsymbol{\eta}}_{\lambda,j})\omega(\boldsymbol{W}_{i}-\boldsymbol{W}_{j}),$

which is quadratic with an explicit unique minimum
$\widehat{\boldsymbol{\beta}}(\lambda)$. Thus, we define a profile SmoothMD
estimator of $\lambda_{0}$ as

$\widehat{\lambda}=\arg\min_{\lambda}\widehat{Q}_{n}\left(\left(\lambda,\widehat{\boldsymbol{\beta}}(\lambda)^{T}\right)^{T}\right),$

and, with at hand the estimate $\widehat{\lambda}$, we eventually calculate
$\widehat{\boldsymbol{\beta}}(\widehat{\lambda})$, the semiparametric SmoothMD
estimate of $\boldsymbol{\beta}_{0}$. In a final step, the function $m(\cdot)$
in model (1) can be estimated from the residuals
$\hat{\epsilon}_{i}=T(Y_{i},\widehat{\lambda})-\boldsymbol{X}_{i}^{T}\widehat{\boldsymbol{\beta}}(\widehat{\lambda})$
by using any established smoothing procedure.

Details of the method are described in Section 2. Under mild regularity
conditions it is then shown that our procedure leads to consistent estimators.
Using standard kernel estimators based on second order kernels and suitable
bandwidth sequences for nonparametric function estimation, we then establish
$\sqrt{n}$-consistency and asymptotic normality, provided that the dimension
of $\boldsymbol{Z}$ is $q<4$. Corresponding test procedures are described in
Section 4. All theoretical results are derived uniformly for all possible
choices in a compact set for the $d_{1},\dots,d_{p+q}$ used in (6). This
provides theoretical grounds for a sample-based choice of
$d_{1},\dots,d_{p+q}$, such as the sample standard deviations of the
components of the vectors
$(\boldsymbol{X}_{i}^{T},\boldsymbol{Z}_{i}^{T})^{T}$.

We wish to note that the generality of our approach implies that the method
may be used as a powerful tool to check parametric models. For example, in
addition to verifying a log-transformation to the response variable, one may
check linearity assumptions. The latter may be done by comparing the outputs
of the parametric model with the results of a semiparametric analysis, where
some of the regressors enter the model nonparametrically and define a
corresponding vector $Z$. This is exemplified by our real data application in
Section 4.

The remainder of the paper is organized as follows. In Section 2, we present
our new estimation method and establish identification of the model
parameters. Theoretical properties of the estimators are derived in Section 3,
while in Section 4, we investigate a distance-metric procedure for testing
restrictions on parameters. In Section 5, we study the finite sample behavior
by a simulation study and apply the estimator to a real data sample. Our
estimator performs well in our experiments and our tests yield accurate levels
and good power in moderate samples. Finally, in Section 6 we formulate few
conclusions and discuss the extension of our approach to higher-dimension
vectors $\boldsymbol{Z}$ using higher-order kernels, as well as efficiency
aspects. The proofs are left to the Appendix.

##  2 The semiparametric SmoothMD approach

In this section we formally define our semiparametric estimator. First, we
investigate two issues. On the one hand, we prove identification of the true
value $\boldsymbol{\theta}_{0}=(\lambda_{0},\boldsymbol{\beta}_{0}^{T})^{T}$
of the parameter of interest. Next, we discuss the recommendation appearing in
the literature for normalizing the response variable. This issue is specific
to the Box-Cox transformation, though similar problems occur with other
families of transformations. Finally, we define our semiparametric SmoothMD
estimator.

Before proceeding with this plan, let us slightly modify the definition of the
conditional moment equation. The function
$g(\boldsymbol{U};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda})$ defined in
(4) has zero-mean if $\boldsymbol{\theta}=\boldsymbol{\theta}_{0}$, but there
is no reason to expect its sample version to be centered, as is the case when
$\lambda_{0}$ is given and one uses least squares; see Li (1996). We therefore
propose to introduce an intercept and hereafter replace
$g(\boldsymbol{U};\boldsymbol{\theta},\boldsymbol{\eta}_{\lambda})$ with

$g(\boldsymbol{U};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda})=\left(T(Y,\lambda)-E[T(Y,\lambda)\mid\boldsymbol{Z}]-(\boldsymbol{X}-E[\boldsymbol{X}\mid\boldsymbol{Z}])^{T}\boldsymbol{\beta}\right)f_{z}(\boldsymbol{Z})-\gamma.$
(8)

The true value of $\gamma$ is known to be $\gamma_{0}=0$, but this intercept
slightly improves the results with finite samples, while it does not introduce
any additional theoretical or computational complexity.

We use the following notation throughout the remainder of the paper. For
$d_{l},d_{c}\geq 1,$ let $\mathbb{R}^{d_{l}\times d_{c}}$ denote the set of
$d_{l}\times d_{c}-$ matrices with real elements. Let $\boldsymbol{1}_{d_{l}}$
(resp. $\boldsymbol{0}_{d_{l}}$) denote the vector with all components equal
to 1 (resp. 0), $\boldsymbol{0}_{d_{l}\times d_{c}}$ the $d_{l}\times
d_{c}-$null matrix and $\boldsymbol{I}_{d_{l}\times d_{c}}$ the identity
matrix with dimension $d_{l}\times d_{c}$. For a matrix $\boldsymbol{A}$,
$\lVert\boldsymbol{A}\rVert$ is the Frobenius norm and
$\lVert\boldsymbol{A}\rVert_{\rm{Sp}}$ the spectral norm. Below,
$\boldsymbol{D}={\rm diag}(\boldsymbol{d})$ is some positive definite diagonal
matrix with $\boldsymbol{d}\in\mathcal{D}\subset\mathbb{R}^{p+q}_{+}$ being a
diagonal vector with strictly positive components. Herein, $\mathcal{D}$ is a
compact set and our asymptotic results are derived uniformly with respect to
$\boldsymbol{d}\in\mathcal{D}$.

###  2.1 Identification

Let $-\infty<\lambda_{\rm min}<\lambda_{0}<\lambda_{\rm max}<\infty$, with
$\lambda_{\min}<0$ and $\lambda_{\max}>0$. For any $\lambda\in[\lambda_{\rm
min},\lambda_{\rm max}]$, let

$(\gamma(\lambda),\boldsymbol{\beta}(\lambda)^{T})^{T}=\arg\min_{\gamma\in\mathbb{R},\boldsymbol{\beta}\in\mathbb{R}^{p}}E[g(\boldsymbol{U}_{1};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda,1})g(\boldsymbol{U}_{2};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda,2})\omega(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})],$
(9)

with
$g(\boldsymbol{U}_{1};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda,1})$
and
$g(\boldsymbol{U}_{2};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda,2})$
being independent copies of
$g(\boldsymbol{U};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda})$
defined in equation (8) with
$\boldsymbol{U}=\left(Y,\boldsymbol{W}^{T}\right)^{T}$,
$\boldsymbol{W}=\left(\boldsymbol{X}^{T},\boldsymbol{Z}^{T}\right)^{T}$,
$\boldsymbol{\theta}=(\lambda,\boldsymbol{\beta}^{T})^{T}$ and
$\boldsymbol{\eta}_{\lambda,k}=\boldsymbol{\eta}_{\lambda}(\boldsymbol{Z}_{k})$,
$k=1,2,$ where

$\boldsymbol{\eta}_{\lambda}(\boldsymbol{z})=(f_{z}(\boldsymbol{z}),E[T(Y,\lambda)\mid\boldsymbol{Z}=\boldsymbol{z}\;],E[\boldsymbol{X}\mid\boldsymbol{Z}=\boldsymbol{z}\;]^{T})^{T}.$

_Data Generating Process_

  1. 1.

The observations
$\left(Y_{i},\boldsymbol{X}_{i}^{T},\boldsymbol{Z}_{i}^{T}\right)^{T}$ ,
$1\leq i\leq n$, are i.i.d. copies of
$\left(Y,\boldsymbol{X}^{T},\boldsymbol{Z}^{T}\right)^{T}\in\mathbb{R}\times\mathbb{R}^{p}\times\mathbb{R}^{q}$.
Moreover, there exists a constant $c>0$ such that $\mathbb{P}(Y>c)=1$.

  2. 2.

The covariate vector $\boldsymbol{Z}$ admits a bounded density in
$\mathbb{R}^{q}$. The covariate vector $\boldsymbol{X}$ is split into two
subvectors $\boldsymbol{X}_{c}\in\mathbb{R}^{p_{c}}$ and
$\boldsymbol{X}_{d}\in\mathbb{R}^{p_{d}}$ with $0\leq p_{c},p_{d}\leq p$ and
$p_{c}+p_{d}=p$. The subvector $\boldsymbol{X}_{c}$ admits a bounded density
in $\mathbb{R}^{p_{c}}$. The subvector $\boldsymbol{X}_{d}$ takes values in a
finite set.

  3. 3.

The diagonal of the matrix $\boldsymbol{D}$ belongs to the $(p+q)-$dimension
cube $\mathcal{D}=[d_{L},d_{U}]^{p+q}$, with some fixed
$0<d_{L}<d_{U}<\infty$.

The assumption that the discrete components of $\boldsymbol{X}$ take values in
a finite set, is a technical condition that simplifies the proofs without
significant restriction of the generality of the applications.

_Identification_

  1. 1.

$E\left[\|\boldsymbol{X}\|^{2}\right]<\infty$,
$E\left[\|\boldsymbol{Z}\|^{2}\right]<\infty$, and
$Var\left[\boldsymbol{X}-E[\boldsymbol{X}|\boldsymbol{Z}]\right]$ has full
rank.

  2. 2.

The true value $\boldsymbol{\beta}_{0,c}\in\mathbb{R}^{p_{c}}$ of the
subvector of coefficients corresponding to $\boldsymbol{X}_{c}$ is not equal
to $\boldsymbol{0}_{p_{c}}$.

  3. 3.

The continuous random subvector $\boldsymbol{X}_{c}$ is such that, for any
$\boldsymbol{b}\in\mathbb{R}^{p_{c}}$,
$\boldsymbol{b}\neq\boldsymbol{0}_{p_{c}}$, the variable
$\boldsymbol{X}^{T}_{c}\boldsymbol{b}$ is continuous with the support equal to
the whole real line.

  4. 4.

Whenever $\lambda\neq\lambda_{0}$, for any $\boldsymbol{z}$ in the support of
$\boldsymbol{Z}$ and $\boldsymbol{x}_{d}$ in the support of the discrete
subvector $\boldsymbol{X}_{d}$, the set of values of the map
$\boldsymbol{x}_{c}\mapsto
E\left[T(Y,\lambda)-T(Y,\lambda_{0})\mid\boldsymbol{X}_{c}=\boldsymbol{x}_{c},\boldsymbol{X}_{d}=\boldsymbol{x}_{d},\boldsymbol{Z}=\boldsymbol{z}\right]$,
$\boldsymbol{x}_{c}\in\mathbb{R}^{p_{c}},$ is unbounded.

  5. 5.

$E\left[Y^{2C_{\lambda}}\right]<\infty$, where
$C_{\lambda}=\max(|\lambda_{min}|,\lambda_{max})<\infty$.

Note that $Var\left[(\boldsymbol{X}^{T},\boldsymbol{Z}^{T})^{T}\right]$
necessarily has full rank, by Assumption 2.1.1 and the fact that
$\boldsymbol{Z}$ admits a density. The complete justification of this
statement is given in the Appendix.

With all this in hand, we can now state the following identification result.

###### Lemma 2.1

Suppose that Assumptions 2.1 and 2.1 hold true. Let $\gamma(\lambda)$ and
$\boldsymbol{\beta}(\lambda)$ be defined as in (9). Then,
$\gamma(\lambda_{0})=0$ and
$\boldsymbol{\beta}(\lambda_{0})=\boldsymbol{\beta}_{0}$ and

$\displaystyle\mathbb{P}\left(E\left[(T(Y,\lambda)-E[T(Y,\lambda)\mid\boldsymbol{Z}])f_{z}(\boldsymbol{Z})-\gamma-(\boldsymbol{X}-E[\boldsymbol{X}\mid\boldsymbol{Z}])^{T}\boldsymbol{\beta}f_{z}(\boldsymbol{Z})\mid\boldsymbol{X},\boldsymbol{Z}\right]=0\right)<1,$

for all $\gamma\in\mathbb{R}$ and
$\boldsymbol{\theta}=(\lambda,\boldsymbol{\beta}^{T})^{T}\in[\lambda_{\rm
min},\lambda_{\rm max}]\times\mathbb{R}^{p}$ such that
$(\gamma,\boldsymbol{\theta}^{T})^{T}\neq(0,\boldsymbol{\theta}_{0}^{T})^{T}$.
Moreover, for any $\varepsilon>0$,

$\inf_{|\lambda-\lambda_{0}|\geq\varepsilon}\;\inf_{\boldsymbol{d}\in\mathcal{D}}E\left[g\left(\boldsymbol{U}_{1};(\lambda,\boldsymbol{\beta}(\lambda)^{T})^{T},\gamma(\lambda),\boldsymbol{\eta}_{\lambda,1}\right)g\left(\boldsymbol{U}_{2};(\lambda,\boldsymbol{\beta}(\lambda)^{T})^{T},\gamma(\lambda),\boldsymbol{\eta}_{\lambda,2}\right)\right.\\\
\times\left.\exp\left\\{-(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})^{T}\boldsymbol{D}(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})\right\\}\right]>0.$
(10)

###  2.2 Box-Cox transformation and standardized responses

Let us note that

$\lim_{\lambda\uparrow\infty}\frac{y^{\lambda}-1}{\lambda}=0\quad\text{if
$0<y<1$}\qquad
and\qquad\lim_{\lambda\downarrow-\infty}\frac{y^{\lambda}-1}{\lambda}=0\quad\text{if
$y>1$}.$

In classical estimation approaches for parametric regression models with Box-
Cox transformed response, this is likely to induce instability for the
estimation of the parameter $\lambda$. See, e.g., Khazzoom (1989), Powell
(1996) and Showalter (1994) for a discussion of this well-known issue. In
order to avoid such problems, the common recommendation is to standardize the
response by some constant, say $s$, such that

$\mathbb{P}\left(Y/s<1\right)>0\qquad\text{and}\qquad\mathbb{P}\left(Y/s>1\right)>0.$

The constant $s$ could be, for instance, the mean of $Y$ or the geometric mean
of $Y$.222The geometric mean is defined as $G(Y)=\exp\\{E(\log(Y)\\}$ and the
sample counterpart is $G_{n}=\prod\limits_{i=1}^{n}Y_{i}^{1/n}$. With finite
samples, the practitioner would first estimate such a constant using the
sample, and next would normalize the responses. The same type of problems
might occur in our semiparametric extension of the Box-Cox transformation
model. For this reason, we will replace our function
$g(\boldsymbol{U};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda})$ by
a family of functions
$s^{-\lambda}g(\boldsymbol{U};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda})$
also indexed by $s$ that we shall let depend on the sample. This change of the
family of functions is equivalent to changing $Y$ to $Y/s$ in the definition
(8), and a rescaling of the parameters $\boldsymbol{\beta}$ and $\gamma$. By
the profiling-based construction of our SmoothMD estimator, the replacement of
the response $Y$ by $Y/s$ matters only for computing $\widehat{\lambda}$.
Clearly, the identifiability property established in Lemma 2.1 is preserved.
In the remainder of the paper, we provide asymptotic results that are uniform
with respect to $s$ in order to allow for a data-driven choice of $s$, such as
for instance, the sample geometric mean of the response.

###  2.3 The estimator

Given an independent sample
$\left(Y_{1},\boldsymbol{X}^{T}_{1},\boldsymbol{Z}^{T}_{1}\right)^{T},\ldots,\left(Y_{n},\boldsymbol{X}^{T}_{n},\boldsymbol{Z}^{T}_{n}\right)^{T}$
from
$\left(Y,\boldsymbol{X}^{T},\boldsymbol{Z}^{T}\right)^{T}\in\mathbb{R}\times\mathbb{R}^{p+q}$,
let us define

$\widehat{\mathbb{Y}}_{n}(\lambda)=\left((T(Y_{1},\lambda)-\widehat{E}[T(Y_{1},\lambda)\mid\boldsymbol{Z}_{1}])\widehat{f}_{z}(\boldsymbol{Z}_{1}),\ldots,(T(Y_{n},\lambda)-\widehat{E}[T(Y_{n},\lambda)\mid\boldsymbol{Z}_{n}])\widehat{f}_{z}(\boldsymbol{Z}_{n})\right)^{T}\in\mathbb{R}^{n},$

and

$\widehat{\mathbb{X}}_{n}=\left((\boldsymbol{X}_{1}-\widehat{E}[\boldsymbol{X}_{1}\mid\boldsymbol{Z}_{1}])\widehat{f}_{z}(\boldsymbol{Z}_{1}),\ldots,(\boldsymbol{X}_{n}-\widehat{E}[\boldsymbol{X}_{n}\mid\boldsymbol{Z}_{n}])\widehat{f}_{z}(\boldsymbol{Z}_{n})\right)^{T}\in\mathbb{R}^{n\times
p}.$

For $1\leq i\leq n$,
$\widehat{\boldsymbol{\eta}}_{\lambda,i}=(\widehat{f}_{z}(\boldsymbol{Z}_{i}),\widehat{E}[T(Y_{i},\lambda)\mid\boldsymbol{Z}_{i}],\widehat{E}[\boldsymbol{X}_{i}\mid\boldsymbol{Z}_{i}]^{T})^{T}$
are nonparametric kernel estimates of
$\boldsymbol{\eta}_{\lambda,i}=(f_{z}(\boldsymbol{Z}_{i}),E[T(Y_{i},\lambda)\mid\boldsymbol{Z}_{i}],E[\boldsymbol{X}_{i}\mid\boldsymbol{Z}_{i}]^{T})^{T}$.
More precisely,

$\widehat{f}_{z}(\boldsymbol{Z}_{i})=\frac{1}{nh^{q}}\sum_{j=1}^{n}K\left(\frac{\boldsymbol{Z}_{i}-\boldsymbol{Z}_{j}}{h}\right),\quad\widehat{E}[T(Y_{i},\lambda)\mid\boldsymbol{Z}_{i}]\widehat{f}_{z}(\boldsymbol{Z}_{i})=\frac{1}{nh^{q}}\sum_{j=1}^{n}T(Y_{j},\lambda)K\left(\frac{\boldsymbol{Z}_{i}-\boldsymbol{Z}_{j}}{h}\right),$

and

$\widehat{E}[\boldsymbol{X}_{i}\mid\boldsymbol{Z}_{i}]\widehat{f}_{z}(\boldsymbol{Z}_{i})=\frac{1}{nh^{q}}\sum_{j=1}^{n}\boldsymbol{X}_{j}K\left(\frac{\boldsymbol{Z}_{i}-\boldsymbol{Z}_{j}}{h}\right).$

Here $K(\cdot)$ is a multivariate kernel function and $h$ is the bandwidth.
Let $\boldsymbol{\Omega}_{n}$ be the $n\times n-$ symmetric matrix with
elements

$\boldsymbol{\Omega}_{n,ij}=\exp\\{-(\boldsymbol{X}_{i}^{T}-\boldsymbol{X}_{j}^{T},\boldsymbol{Z}_{i}^{T}-\boldsymbol{Z}_{j}^{T})\boldsymbol{D}(\boldsymbol{X}_{i}-\boldsymbol{X}_{j},\boldsymbol{Z}_{i}-\boldsymbol{Z}_{j})\\},\qquad
1\leq i,j\leq n.$

Typically, the components of the vector $\boldsymbol{d}$ defining the diagonal
matrix $\boldsymbol{D}$, could be taken as proportional to the standard
deviation of the components of the vectors
$(\boldsymbol{X}_{i}^{T},\boldsymbol{Z}_{i}^{T})^{T}$. The definition of
$\boldsymbol{\Omega}_{n,i,j}$ allows also to take into account, discrete
components of $\boldsymbol{X}$. For finite support discrete covariates, one
could set some large value for the corresponding diagonal element of
$\boldsymbol{D}$, which in practice would be equivalent to an indicator of the
event that the observations $i$ and $j$ have the same value for that
covariate.

We can now define, for any $\lambda$, the estimates of
$(\gamma(\lambda),\boldsymbol{\beta}(\lambda)^{T})^{T}\in\mathbb{R}^{1+p}$
introduced in (9). For any $s>0$, let

$\widehat{Q}_{n}\left(\left(\lambda,{\boldsymbol{\beta}}^{T}\right)^{T},{\gamma};s\right)=n^{-2}s^{-2\lambda}\left(\widehat{\mathbb{Y}}_{n}(\lambda)-\gamma\boldsymbol{1}_{n}-\widehat{\mathbb{X}}_{n}\boldsymbol{\beta}\right)^{T}\boldsymbol{\Omega}_{n}\left(\widehat{\mathbb{Y}}_{n}(\lambda)-\gamma\boldsymbol{1}_{n}-\widehat{\mathbb{X}}_{n}\boldsymbol{\beta}\right).$

For fixed $s$ and $\lambda$, consider the generalized least-squares problem

$\min_{\gamma,\boldsymbol{\beta}}\widehat{Q}_{n}\left(\left(\lambda,{\boldsymbol{\beta}}^{T}\right)^{T},{\gamma};s\right).$
(11)

The solution of this problem does not depend on $s^{-\lambda}$ and has the
form of standard generalized least-squares estimators:

$\widehat{\gamma}(\lambda,{\boldsymbol{\beta}}(\lambda))=\frac{1}{\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{1}_{n}}\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\left(\widehat{\mathbb{Y}}_{n}(\lambda)-\widehat{\mathbb{X}}_{n}{\boldsymbol{\beta}}(\lambda)\right),$

and

$\widehat{\boldsymbol{\beta}}(\lambda)=\left(\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\right)^{-1}\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{Y}}_{n}(\lambda),$

with

${\mathbb{D}}_{n}=\boldsymbol{\Omega}_{n}-\frac{1}{\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{1}_{n}}\boldsymbol{\Omega}_{n}\boldsymbol{1}_{n}\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\in\mathbb{R}^{n\times
n}.$ (12)

Next, plugging
$(\widehat{\gamma}(\lambda,\widehat{\boldsymbol{\beta}}(\lambda)),\widehat{\boldsymbol{\beta}}(\lambda)^{T})^{T}$
into the problem (11), for a given $s$, we define the SmoothMD estimator of
$\lambda_{0}$ as

$\widehat{\lambda}=\widehat{\lambda}(s)=\arg\min_{\lambda\in\Lambda}s^{-\lambda}\widehat{\mathbb{Y}}_{n}(\lambda)^{T}\;\widehat{\mathbb{B}}_{n}\;s^{-\lambda}\widehat{\mathbb{Y}}_{n}(\lambda),$
(13)

with

$\widehat{\mathbb{B}}_{n}={\mathbb{D}}_{n}-{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\left(\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\right)^{-1}\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\in\mathbb{R}^{n\times
n}.$

Note that, by construction,

${\mathbb{D}}_{n}\boldsymbol{1}_{n}=\widehat{\mathbb{B}}_{n}\boldsymbol{1}_{n}=\boldsymbol{0}_{n}\quad\text{
and
}\quad\widehat{\mathbb{B}}_{n}\widehat{\mathbb{X}}_{n}=\boldsymbol{0}_{n\times
p}.$

Finally, the SmoothMD estimator of $\boldsymbol{\beta}_{0}$ is
$\widehat{\boldsymbol{\beta}}(\widehat{\lambda})$. We close this section by
showing that our estimator is well-defined.

###### Lemma 2.2

If Assumptions 2.1.3 and 2.1 hold true, then, for each $n\geq 1$,

  1. 1.

the matrices $\boldsymbol{\Omega}_{n}$ and
$\widehat{\mathbb{X}}_{n}^{T}\mathbb{D}_{n}\widehat{\mathbb{X}}_{n}$ are
positive definite with probability 1. In particular,
$\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{1}_{n}>0$ and
$\widehat{\mathbb{X}}_{n}^{T}\mathbb{D}_{n}\widehat{\mathbb{X}}_{n}$ is
invertible with probability 1.

  2. 2.

the matrix $\widehat{\mathbb{B}}_{n}$ is positive semi-definite with
probability 1.

The matrix ${\mathbb{D}}_{n}$ is defined in equation (12) and has dimension
$n\times n$. Therefore, one might imagine that it becomes difficult to work
with this matrix when the sample size is large. This is not the case because
it is not necessary to estimate the matrix ${\mathbb{D}}_{n}$ itself. It
suffices to compute $\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}$ and
$\widehat{\mathbb{Y}}_{n}(\lambda)^{T}{\mathbb{D}}_{n}$ to be able to
calculate $\widehat{\boldsymbol{\beta}}(\lambda)$ and $\widehat{\lambda}$. In
chapter 5 we shall show that the estimator can be easily applied even for
$n>100,000$.

##  3 Consistency and asymptotic normality

The estimator $\widehat{\lambda}(s)$ depends on $s$, a value that in practice
could be calculated from the sample, such as for instance, the sample
geometric mean. For this reason, our asymptotic results are stated uniformly
with respect to $s$. Our asymptotic results are also stated uniformly with
respect to the diagonal of the matrix $\boldsymbol{D}$. This ensures that we
can use a data driven estimate of $\boldsymbol{D}$ proportional to the
empirical standard deviation of $\boldsymbol{X}$ or $\boldsymbol{Z}$.

Let’s introduce some more notation: for each $\lambda\in\Lambda$, let

${\mathbb{Y}}_{n}(\lambda)=\left((T(Y_{1},\lambda)-{E}[T(Y_{1},\lambda)\mid\boldsymbol{Z}_{1}]){f}_{z}(\boldsymbol{Z}_{1}),\ldots,(T(Y_{n},\lambda)-{E}[T(Y_{n},\lambda)\mid\boldsymbol{Z}_{n}]){f}_{z}(\boldsymbol{Z}_{n})\right)^{T}\in\mathbb{R}^{n},$

and

${\mathbb{X}}_{n}=\left((\boldsymbol{X}_{1}-{E}[\boldsymbol{X}_{1}\mid\boldsymbol{Z}_{1}]){f}_{z}(\boldsymbol{Z}_{1}),\ldots,(\boldsymbol{X}_{n}-{E}[\boldsymbol{X}_{n}\mid\boldsymbol{Z}_{n}]){f}_{z}(\boldsymbol{Z}_{n})\right)^{T}\in\mathbb{R}^{n\times
p}.$

Moreover,

${\mathbb{B}}_{n}={\mathbb{D}}_{n}-{\mathbb{D}}_{n}{\mathbb{X}}_{n}\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\in\mathbb{R}^{n\times
n},$

with ${\mathbb{D}}_{n}$ defined in equation (12). Again, by construction
${\mathbb{B}}_{n}\boldsymbol{1}_{n}=\boldsymbol{0}_{n}$ and
${\mathbb{B}}_{n}{\mathbb{X}}_{n}=\boldsymbol{0}_{n\times p}$.

_Consistency_

  1. 1.

The functions $f_{z}(\cdot)$, $(mf_{z})(\cdot)$,
$E[\|\boldsymbol{X}\|^{2}\mid\boldsymbol{Z}=\cdot\;]f_{z}(\cdot)$ and
$\sup_{\lambda\in\Lambda}(\partial^{2}/\partial\lambda^{2})E[T(Y,\lambda)\mid\boldsymbol{Z}=\cdot\;]f_{z}(\cdot)$
have Hölder continuous partial derivatives of order four.

  2. 2.

The kernel $K(\cdot)$ is the product of $q$ univariate kernel functions
$\widetilde{K}$ of bounded variation. Moreover, $\widetilde{K}$ is a symmetric
function with integral equal to one and $\int t^{2}\widetilde{K}(t)dt<\infty$.

  3. 3.

The bandwidth $h$ belongs to a range
$\mathcal{H}_{c,n}=[c_{min}n^{-\alpha},c_{max}n^{-\alpha}],$ with
$0<\alpha<1/q$ and $c_{min}$, $c_{max}$ positive constants.

With all this in hand, we can now state the consistency of our estimator.

######  Theorem 3.1 (Consistency)

Assume that Assumptions 2.1, 2.1 and 3 hold true. Let $s_{0}$ be some
normalizing value such that $\mathbb{P}\left(Y/s_{0}<1\right)>0$ and
$\mathbb{P}\left(Y/s_{0}>1\right)>0$ and let $S_{n}$ be an arbitrary
$o_{\mathbb{P}}(1)$ neighborhood of $s_{0}$. Then

$\sup_{h\in\mathcal{H}_{c,n}}\sup_{s\in
S_{n}}\sup_{\boldsymbol{d}\in\mathcal{D}}\left|\widehat{\lambda}-\lambda_{0}\right|=o_{\mathbb{P}}(1)\quad\text{
and }\quad\sup_{h\in\mathcal{H}_{c,n}}\sup_{s\in
S_{n}}\sup_{\boldsymbol{d}\in\mathcal{D}}\left\|\widehat{\boldsymbol{\beta}}(\widehat{\lambda})-\boldsymbol{\beta}_{0}\right\|=o_{\mathbb{P}}(1).$

In Theorem 3.1 we require that $h\in\mathcal{H}_{c,n}$. This implies that
$nh^{q}\rightarrow\infty$ and $h\rightarrow 0$ for $n\rightarrow\infty$ which
is in line with Robinson (1988) and Li (1996).

Next, we prove asymptotic normality for our estimator. For this purpose, we
first derive the asymptotic linear representation of $\widehat{\lambda}$ and
and $\widehat{\boldsymbol{\beta}}(\widehat{\lambda})$ from which the
$\sqrt{n}-$asymptotic normality follows. In the following result, we show that
$\widehat{\lambda}$ and $\widehat{\boldsymbol{\beta}}(\widehat{\lambda})$ are
asymptotically not equivalent to the infeasible estimators of $\lambda_{0}$
and $\boldsymbol{\beta}_{0}$, one would obtain when the infinite-dimensional
parameter $\boldsymbol{\eta}_{\lambda}$ is given and the intercept $\gamma$ is
equal to 0. This is in contrast to the results of Li (1996) and Robinson
(1988). The reason is that they can use the fact that
$E[{\mathbb{X}}_{n,i}\mid\boldsymbol{Z}_{i}]=0$ when controlling higher order
terms. In our case, we weight the observations by $\boldsymbol{\Omega}_{n,ij}$
and, in general,
$E[{\mathbb{X}}_{n,i}\boldsymbol{\Omega}_{n,ij}\mid\boldsymbol{Z}_{i}]\neq 0$
for $i\neq j$. This is also the reason why we need to ask for $q<4$ instead of
$q<6$ as in Li (1996). Therefore, we require that $h\in\mathcal{H}_{sc,n}$,
where $\mathcal{H}_{sc,n}=[c_{min}n^{-\alpha},c_{max}n^{-\alpha}]$, with
$\alpha\in(1/4,1/q)$. This is the small price we have to pay for getting what,
to the best of our knowledge, is the first consistent estimation procedure for
the semiparametric transformation model we investigate. However, as discussed
in section 6, one could still apply our approach with $q\geq 4$, provided one
uses higher-order kernels.

_Asymptotic Normality_

  1. 1.

$E\left[m(\boldsymbol{Z})^{4}\right]<\infty$ and
$Var\left[\frac{\partial}{\partial\lambda}{T(Y,\lambda_{0})}\right]>0$.

  2. 2.

$E\left[\varepsilon^{4}\right]<\infty$ and
$E\left[\varepsilon^{2}\mid\boldsymbol{X},\boldsymbol{Z}\right]=\sigma^{2}(\boldsymbol{X},\boldsymbol{Z})$
is in $L^{1}\cap L^{2}$.

  3. 3.

The bandwidth $h$ belongs to
$\mathcal{H}_{sc,n}=[c_{min}n^{-\alpha},c_{max}n^{-\alpha}],$ with
$\alpha\in(1/4,1/q)$ and $c_{min}$, $c_{max}>0$.

The results are again obtained uniformly with respect to the elements on the
diagonal of the matrix $\boldsymbol{D}$ that determines
$\boldsymbol{\Omega}_{n}$ and with respect to the scaling factor $s$ that
could be used for numerical stability, as mentioned in Section 2.2. In
addition, let $K_{h}(\cdot)=h^{-q}K(\cdot/h)$ and, for any $1\leq i,j\leq n$,
let

$K_{h,ij}=K_{h}(\boldsymbol{Z}_{i}-\boldsymbol{Z}_{j}).$

######  Proposition 3.1 (Asymptotic representation)

Assume that the conditions of Theorem 3.1 and Assumption 3 hold true. Then,
uniformly with respect to $h\in\mathcal{H}_{sc,n}$,
$\boldsymbol{d}\in\mathcal{D}$ and $s\in S_{n}$,

$\widehat{\lambda}-\lambda_{0}=-\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}\;{\mathbb{B}}_{n}\;\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]^{-1}\;\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}\;{\mathbb{B}}_{n}\left[(\boldsymbol{\varepsilon
f_{z}})_{n}-\left(\boldsymbol{\widehat{\varepsilon}_{|z}\widehat{f}_{z}}\right)_{n}\right]+o_{\mathbb{P}}(n^{-1/2})=O_{\mathbb{P}}(n^{-1/2}),$

and

$\widehat{\boldsymbol{\beta}}(\widehat{\lambda})-\boldsymbol{\beta}_{0}=\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\left[(\boldsymbol{\varepsilon
f_{z}})_{n}-\left(\boldsymbol{\widehat{\varepsilon}_{|z}\widehat{f}_{z}}\right)_{n}+\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\left(\widehat{\lambda}-\lambda_{0}\right)\right]+o_{\mathbb{P}}(n^{-1/2})=O_{\mathbb{P}}(n^{-1/2}),$

where $(\boldsymbol{\varepsilon
f_{z}})_{n}=(\varepsilon_{1}{f}_{z}(\boldsymbol{Z}_{1}),\ldots,\varepsilon_{n}{f}_{z}(\boldsymbol{Z}_{n}))^{T}$
and
$\left(\boldsymbol{\widehat{\varepsilon}_{|z}\widehat{f}_{z}}\right)_{n}=\left(\frac{1}{n}\sum\limits_{k=1,k\neq
1}^{n}\varepsilon_{k}K_{h,1k},\ldots,\frac{1}{n}\sum\limits_{k=1,k\neq
n}^{n}\varepsilon_{k}K_{h,nk}\right)^{T}$.

Note that the asymptotic representation of $\widehat{\lambda}$ does not depend
on $s_{0}$, i.e. the choice of $s_{0}$ does not influence the asymptotic
behavior of $\widehat{\lambda}$. This result is in line with the result of
Powell (1996).

Next, we state the asymptotic normality of our estimator. We use the notation
$\boldsymbol{\Omega}_{n,i,j}(\boldsymbol{d})=\boldsymbol{\Omega}_{n,i,j}$ and

${\mathbb{D}}_{n}(\boldsymbol{d})=\boldsymbol{\Omega}_{n}(\boldsymbol{d})-\frac{1}{\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}(\boldsymbol{d})\boldsymbol{1}_{n}}\boldsymbol{\Omega}_{n}(\boldsymbol{d})\boldsymbol{1}_{n}\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}(\boldsymbol{d}),$

to make the dependence of $\boldsymbol{\Omega}_{n}$ on $\boldsymbol{d}$
explicit. Note that with

$\displaystyle\boldsymbol{\Omega}^{X}_{n,ij}(\boldsymbol{d})=\boldsymbol{\Omega}^{X}_{n,ij}$
$\displaystyle=\exp\\{-(\boldsymbol{X}_{i}-\boldsymbol{X}_{j})^{T}{\rm
diag}(d_{1},\ldots,d_{p})(\boldsymbol{X}_{i}-\boldsymbol{X}_{j})\\}\quad\qquad\text{and}$
$\displaystyle\boldsymbol{\Omega}^{Z}_{n,ij}(\boldsymbol{d})=\boldsymbol{\Omega}^{Z}_{n,ij}$
$\displaystyle=\exp\\{-(\boldsymbol{Z}_{i}-\boldsymbol{Z}_{j})^{T}{\rm
diag}(d_{p+1},\ldots,d_{p+q})(\boldsymbol{Z}_{i}-\boldsymbol{Z}_{j})\\},\quad
1\leq i,j\leq n,$

$\boldsymbol{\Omega}_{n,ij}(\boldsymbol{d})=\boldsymbol{\Omega}^{X}_{n,ij}(\boldsymbol{d})\boldsymbol{\Omega}^{Z}_{n,ij}(\boldsymbol{d})$.
Furthermore, we define, for $1\leq i\leq n$,

$\displaystyle\\!\\!\\!\\!\\!\\!\boldsymbol{\tau}_{i}(\boldsymbol{d}):=\Bigg{(}\\!\left(\\!\frac{\partial}{\partial\lambda}\mathbb{Y}_{n,i}-\frac{1}{E\left[\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}(\boldsymbol{d})\boldsymbol{1}_{n}\right]}E\left[\frac{\partial}{\partial\lambda}\mathbb{Y}_{n}^{T}\boldsymbol{\Omega}_{n}(\boldsymbol{d})\boldsymbol{1}_{n}\right]\right),-\left(\\!\mathbb{X}_{n,i}^{T}-\frac{1}{E\left[\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}(\boldsymbol{d})\boldsymbol{1}_{n}\right]}E\left[\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}(\boldsymbol{d}){\mathbb{X}}_{n}\right]\right)\\!\Bigg{)}^{\\!T}\\!,$

where
$\frac{\partial}{\partial\lambda}\mathbb{Y}_{n,i}(\lambda)=\left(\frac{\partial}{\partial\lambda}T(Y_{i},\lambda)-{E}[\frac{\partial}{\partial\lambda}T(Y_{i},\lambda)\mid\boldsymbol{Z}_{i}]\right){f}_{z}(\boldsymbol{Z}_{i})$
and
$\mathbb{X}_{n,i}=(\boldsymbol{X}_{i}-{E}[\boldsymbol{X}_{i}\mid\boldsymbol{Z}_{i}]){f}_{z}(\boldsymbol{Z}_{i})$.
In addition, let

$\displaystyle\boldsymbol{\Phi}^{X}_{n,ij}(\boldsymbol{d})=\boldsymbol{\Omega}^{X}_{n,ij}(\boldsymbol{d})-E\left[\boldsymbol{\Omega}^{X}_{n,ik}(\boldsymbol{d})\mid\boldsymbol{X}_{i}\right].$

Finally, for any vector $\boldsymbol{a}$, we denote $\boldsymbol{a}^{\otimes
2}=\boldsymbol{a}\boldsymbol{a}^{T}$. With all this in hand, we can state the
following result.

######  Theorem 3.2 (Asymptotic normality)

Assume that the conditions of Proposition 3.1 hold true. Then, uniformly with
respect to $h\in\mathcal{H}_{sc,n}$, $\boldsymbol{d}\in\mathcal{D}$ and $s\in
S_{n}$,

$\displaystyle\sqrt{n}\left((\widehat{\lambda},\widehat{\boldsymbol{\beta}}(\widehat{\lambda})^{T})^{T}-(\lambda_{0},\boldsymbol{\beta}_{0}^{T})^{T}\right)+\boldsymbol{V}(\boldsymbol{d})^{-1}\left(\frac{1}{\sqrt{n}}\sum_{j=1}^{n}\varepsilon_{j}{f}_{z}(\boldsymbol{Z}_{j})E\left[\boldsymbol{\tau}_{i}(\boldsymbol{d})\,\boldsymbol{\Omega}_{n,ij}^{Z}(\boldsymbol{d})\boldsymbol{\Phi}^{X}_{n,ij}(\boldsymbol{d})\mid\boldsymbol{X}_{j},\boldsymbol{Z}_{j}\right]\right)=o_{\mathbb{P}}\left(1\right),$

where

$\boldsymbol{V}(\boldsymbol{d})=\underset{n\rightarrow\infty}{\lim}\begin{pmatrix}E\left[n^{-2}\frac{\partial}{\partial\lambda}\mathbb{Y}_{n}(\lambda_{0})^{T}{\mathbb{D}}_{n}(\boldsymbol{d})\frac{\partial}{\partial\lambda}\mathbb{Y}_{n}(\lambda_{0})\right]&-E\left[n^{-2}\frac{\partial}{\partial\lambda}\mathbb{Y}_{n}(\lambda_{0})^{T}{\mathbb{D}}_{n}(\boldsymbol{d})\mathbb{X}_{n}\right]\\\
-E\left[n^{-2}\mathbb{X}_{n}^{T}{\mathbb{D}}_{n}(\boldsymbol{d})\frac{\partial}{\partial\lambda}\mathbb{Y}_{n}(\lambda_{0})\right]&E\left[n^{-2}\mathbb{X}_{n}^{T}{\mathbb{D}}_{n}(\boldsymbol{d})\mathbb{X}_{n}\right]\end{pmatrix}.$

As a consequence,
$\sqrt{n}\left((\widehat{\lambda},\widehat{\boldsymbol{\beta}}(\widehat{\lambda})^{T})^{T}-(\lambda_{0},\boldsymbol{\beta}_{0}^{T})^{T}\right)$
converges in distribution to a $(p+1)-$dimension centered Gaussian vector with
variance
$\boldsymbol{V}(\boldsymbol{d})^{-1}\boldsymbol{\Delta}(\boldsymbol{d})\boldsymbol{V}(\boldsymbol{d})^{-1}$
where

$\boldsymbol{\Delta}(\boldsymbol{d})=E\left\\{Var\left[\varepsilon_{j}\mid\boldsymbol{X}_{j},\boldsymbol{Z}_{j}\right]{f}^{2}_{z}(\boldsymbol{Z}_{j})\left(E\left[\boldsymbol{\tau}_{i}(\boldsymbol{d})\,\boldsymbol{\Omega}_{n,ij}^{Z}(\boldsymbol{d})\boldsymbol{\Phi}^{X}_{n,ij}(\boldsymbol{d})\mid\boldsymbol{X}_{j},\boldsymbol{Z}_{j}\right]\right)^{\otimes
2}\right\\}.$

If $\boldsymbol{\eta}_{\lambda}$ were known, which corresponds to the case
studied by Lavergne and Patilea (2013),
$\boldsymbol{\Phi}^{X}_{n,ij}(\boldsymbol{d})$ should be replaced by
$\boldsymbol{\Omega}^{X}_{n,ij}(\boldsymbol{d})$ in the expression of
$\boldsymbol{\Delta}(\boldsymbol{d})$.

We can estimate the covariance matrix by
$\widehat{\boldsymbol{V}}(\boldsymbol{d})^{-1}\widehat{\boldsymbol{\Delta}}(\boldsymbol{d})\widehat{\boldsymbol{V}}(\boldsymbol{d})^{-1}$,
where

$\displaystyle\widehat{\boldsymbol{V}}(\boldsymbol{d})$
$\displaystyle=\begin{pmatrix}n^{-2}\frac{\partial}{\partial\lambda}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda})^{T}{\mathbb{D}}_{n}(\boldsymbol{d})\frac{\partial}{\partial\lambda}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda})&-n^{-2}\frac{\partial}{\partial\lambda}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda})^{T}{\mathbb{D}}_{n}(\boldsymbol{d})\widehat{\mathbb{X}}_{n}\\\
-n^{-2}\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}(\boldsymbol{d})\frac{\partial}{\partial\lambda}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda})&n^{-2}\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}(\boldsymbol{d})\widehat{\mathbb{X}}_{n}\end{pmatrix}$
(14) and $\displaystyle\widehat{\boldsymbol{\Delta}}(\boldsymbol{d})$
$\displaystyle=n^{-3}\left(\frac{\partial}{\partial\lambda}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda}),-\widehat{\mathbb{X}}_{n}\right)^{T}{\mathbb{D}}_{n,inf}(\boldsymbol{d})\widehat{\boldsymbol{\Phi}}_{n}(\boldsymbol{d})\widehat{\boldsymbol{\Sigma}}_{n}\widehat{\boldsymbol{\Phi}}^{T}_{n}(\boldsymbol{d}){\mathbb{D}}^{T}_{n,inf}(\boldsymbol{d})\left(\frac{\partial}{\partial\lambda}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda}),-\widehat{\mathbb{X}}_{n}\right).$

Here, $\widehat{\boldsymbol{\Phi}}^{X}_{n}$ and
$\widehat{\boldsymbol{\Phi}}_{n}$ are the $n\times n-$ symmetric matrices with
elements

$\displaystyle\widehat{\boldsymbol{\Phi}}^{X}_{n,ij}(\boldsymbol{d})$
$\displaystyle=\boldsymbol{\Omega}^{X}_{n,ij}(\boldsymbol{d})-\frac{1}{n}\sum\limits_{k=1}^{n}\boldsymbol{\Omega}^{X}_{n,ik}(\boldsymbol{d}),\quad
1\leq i,j\leq n$
$\displaystyle\widehat{\boldsymbol{\Phi}}_{n,ij}(\boldsymbol{d})$
$\displaystyle=\widehat{\boldsymbol{\Phi}}^{X}_{n,ij}(\boldsymbol{d})\boldsymbol{\Omega}^{Z}_{n,ij}(\boldsymbol{d}),\hskip
49.50795pt1\leq i,j\leq n$ $\displaystyle
and\quad{\mathbb{D}}_{n,inf}(\boldsymbol{d})$
$\displaystyle=\boldsymbol{I}_{n\times
n}-\frac{1}{\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}(\boldsymbol{d})\boldsymbol{1}_{n}}\boldsymbol{\Omega}_{n}(\boldsymbol{d})\boldsymbol{1}_{n}\boldsymbol{1}_{n}^{T}.$

$\widehat{\boldsymbol{\Sigma}}_{n}=$
diag$\left(\widehat{Var}\left[\varepsilon_{1}f_{z}(\boldsymbol{Z}_{1})\mid\boldsymbol{X}_{1},\boldsymbol{Z}_{1}\right],\ldots,\widehat{Var}\left[\varepsilon_{n}{f}_{z}(\boldsymbol{Z}_{n})\mid\boldsymbol{X}_{n},\boldsymbol{Z}_{n}\right]\right)$
is an estimator of
$\rm{diag}\Big{(}Var\left[\varepsilon_{1}f_{z}(\boldsymbol{Z}_{1})\mid\boldsymbol{X}_{1},\boldsymbol{Z}_{1}\right],\\\
\ldots,Var\left[\varepsilon_{n}{f}_{z}(\boldsymbol{Z}_{n})\mid\boldsymbol{X}_{n},\boldsymbol{Z}_{n}\right]\Big{)}$.
One can use a nonparametric estimator for the conditional variance or
alternatively, use an estimate of the error terms to approximate the
conditional variance in the spirit of the Eiker-White variance estimator.
Consistency of the above estimators is straightforward to establish.

##  4 Testing based on SmoothMD for parameter restrictions

In Section 3 we established consistency and asymptotic normality of our
estimator. The asymptotic behavior of our estimator is not influenced by the
standardization with $s$ but the asymptotic variance is affected by the
estimation of $\boldsymbol{\eta}_{\lambda}$. In addition, the behavior of our
estimator is, even asymptotically, influenced by the vectors
$\boldsymbol{d}_{1}$ and $\boldsymbol{d}_{2}$. When developing a test theory,
we should take that influence into account in order to get reliable results.
That’s what we do in the following.

###  4.1 Testing the transformation parameter

When it comes to testing parameter restrictions in the semiparametric
partially linear regression model with Box-Cox transformation we might be
mainly interested in testing if $\lambda$ is zero or not and if the components
of $\boldsymbol{\beta}$ are zero. However, we shall consider here a more
general approach to allow for more complex hypotheses as well. We separate the
discussion into two parts. In the first part, we consider only restrictions
for $\lambda$ and in the second part we consider restrictions for
$\boldsymbol{\beta}$ with and without restricting $\lambda$.

Suppose we want to test the restriction for $\lambda$ given by

$\displaystyle H_{0}:\lambda_{0}=\lambda_{R}.$ (15)

In order to test this restriction, we can use the distance metric statistic
proposed by Lavergne and Patilea (2013). Adapted to our case and for testing
(15), we consider the distance

$\displaystyle
DM_{\lambda}=\frac{1}{n}\widehat{\mathbb{Y}}_{n}(\lambda_{R})^{T}\;\widehat{\mathbb{B}}_{n}\widehat{\mathbb{Y}}_{n}(\lambda_{R})-\frac{1}{n}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda})^{T}\;\widehat{\mathbb{B}}_{n}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda}).$

The distance metric is based on the object that is minimized to get the
estimate for $\lambda$, see equation (13). However, the test statistic is not
standardized by $s^{-\lambda}$ as we only need this for the estimation of
$\lambda$.

Let

$\displaystyle\boldsymbol{A}_{n}=\begin{pmatrix}\frac{\partial}{\partial\lambda}\mathbb{Y}_{n}(\lambda_{0})^{T}\\\
-\mathbb{X}_{n}^{T}\end{pmatrix}{\mathbb{D}}_{n}\left((\boldsymbol{\varepsilon
f_{z}})_{n}-\left(\boldsymbol{\widehat{\varepsilon}_{|z}\widehat{f}_{z}}\right)_{n}\right).$

We can therefore now state the following Proposition.

###### Proposition 4.1

Assume that the conditions of Proposition 3.1 hold true. Then, uniformly with
respect to $h\in\mathcal{H}_{sc,n}$, $\boldsymbol{d}\in\mathcal{D}$ and $s\in
S_{n}$,

$DM_{\lambda}-(1,\boldsymbol{0}_{p}^{T})\boldsymbol{V}(\boldsymbol{d})^{-1}n^{-3/2}\boldsymbol{A}_{n}n^{-3/2}\boldsymbol{A}_{n}^{T}\boldsymbol{V}(\boldsymbol{d})^{-1}(1,\boldsymbol{0}_{p}^{T})^{T}E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]=o_{\mathbb{P}}(1),$

under $H_{0}$ and $\mathbb{P}(n^{-1}DM_{\lambda}>c)\rightarrow 1$ for any
$c>0$ if $H_{0}$ does not hold.

The process
$(1,\boldsymbol{0}_{p}^{T})\boldsymbol{V}(\boldsymbol{d})^{-1}n^{-3/2}\boldsymbol{A}_{n}n^{-3/2}\boldsymbol{A}_{n}^{T}\boldsymbol{V}(\boldsymbol{d})^{-1}(1,\boldsymbol{0}_{p}^{T})^{T}E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]$
is asymptotically tight and for each $\boldsymbol{d}$ behaves asymptotically
as a chi-square times
$(1,\boldsymbol{0}_{p}^{T})\boldsymbol{V}(\boldsymbol{d})^{-1}\boldsymbol{\Delta}(\boldsymbol{d})\boldsymbol{V}(\boldsymbol{d})^{-1}(1,\boldsymbol{0}_{p}^{T})^{T}E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]$,
see Johnson et al. (1995). The distribution of the distance metric statistic
therefore is, in general, non-pivotal. Determining critical values requires
the estimation of
$(1,\boldsymbol{0}_{p}^{T})\boldsymbol{V}(\boldsymbol{d})^{-1}\boldsymbol{\Delta}(\boldsymbol{d})\boldsymbol{V}(\boldsymbol{d})^{-1}(1,\boldsymbol{0}_{p}^{T})^{T}E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]$,
which can rely on the estimators stated in (14).

###  4.2 Testing the slope coefficients

In the next part, we consider restrictions for $\boldsymbol{\beta}$. Suppose
we want to test $r$ linear restrictions for $\boldsymbol{\beta}$ given by

$\displaystyle H_{0}:\boldsymbol{R}\boldsymbol{\beta}_{0}=\boldsymbol{c},$
(16)

where $\boldsymbol{R}$ is a $r\times p-$ matrix of full rank and
$\boldsymbol{c}\in\mathbb{R}^{r}$. In order to test the restrictions, we need
to find the restricted estimators for $\boldsymbol{\beta}_{0}$,
$\widehat{\boldsymbol{\beta}}_{R}(\lambda)$, and $\lambda_{0}$,
$\widehat{\lambda}_{R}$. We minimize

$\displaystyle
n^{-2}s^{-2\lambda}\left(\widehat{\mathbb{Y}}_{n}(\lambda)-\widehat{\mathbb{X}}_{n}\boldsymbol{\beta}\right)^{T}\boldsymbol{\mathbb{D}}_{n}\left(\widehat{\mathbb{Y}}_{n}(\lambda)-\widehat{\mathbb{X}}_{n}\boldsymbol{\beta}\right)\quad
s.t.\quad\boldsymbol{R}\boldsymbol{\beta}=\boldsymbol{c},$

with respect to $\boldsymbol{\beta}$ and get that

$\displaystyle\widehat{\boldsymbol{\beta}}_{R}(\lambda)=\widehat{\boldsymbol{\beta}}(\lambda)-\left(\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\right)^{-1}\boldsymbol{R}^{T}\left(\boldsymbol{R}\left(\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\right)^{-1}\boldsymbol{R}^{T}\right)^{-1}\left(\boldsymbol{R}\widehat{\boldsymbol{\beta}}(\lambda)-\boldsymbol{c}\right).$

The restricted estimator for $\lambda_{0}$ is then given by

$\widehat{\lambda}_{R}=\widehat{\lambda}_{R}(s)=\arg\min_{\lambda\in\Lambda}s^{-\lambda}\left(\widehat{\mathbb{Y}}_{n}(\lambda)-\widehat{\mathbb{X}}_{n}\widehat{\boldsymbol{\beta}}_{R}(\lambda)\right)^{T}{\mathbb{D}}_{n}\,s^{-\lambda}\left(\widehat{\mathbb{Y}}_{n}(\lambda)-\widehat{\mathbb{X}}_{n}\widehat{\boldsymbol{\beta}}_{R}(\lambda)\right).$
(17)

With all the estimators in hand, we can now define our distance metric
statistic for testing (16).

$\displaystyle
DM_{\boldsymbol{\beta}}=\frac{1}{n}\left(\widehat{\mathbb{Y}}_{n}(\widehat{\lambda}_{R})-\widehat{\mathbb{X}}_{n}\widehat{\boldsymbol{\beta}}_{R}(\widehat{\lambda}_{R})\right)^{T}{\mathbb{D}}_{n}\,\left(\widehat{\mathbb{Y}}_{n}(\widehat{\lambda}_{R})-\widehat{\mathbb{X}}_{n}\widehat{\boldsymbol{\beta}}_{R}(\widehat{\lambda}_{R})\right)-\frac{1}{n}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda})^{T}\;\widehat{\mathbb{B}}_{n}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda}).$

The distance metric is based on the object that is minimized to get the
restricted estimate for $\lambda$, see equation (17). Once again the test
statistic is not standardized by $s^{-\lambda}$ as we only need this for the
estimation of $\lambda$. Let,

${\mathbb{B}}_{n,R}={\mathbb{B}}_{n}+{\mathbb{D}}_{n}{\mathbb{X}}_{n}\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}\boldsymbol{R}^{T}\left(\boldsymbol{R}\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}\boldsymbol{R}^{T}\right)^{-1}\boldsymbol{R}\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n},$

and

$\boldsymbol{V}_{R}(\boldsymbol{d})=\left(E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n,R}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]^{-1},E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n,R}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]^{-1}E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}{\mathbb{B}}_{n}^{+}\right]\right),$

where

${\mathbb{B}}_{n}^{+}=\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}-\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}\boldsymbol{R}^{T}\left(\boldsymbol{R}\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}\boldsymbol{R}^{T}\right)^{-1}\boldsymbol{R}\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}.$

We can therefore now state the following proposition.

###### Proposition 4.2

Assume that the conditions of Proposition 3.1 hold true. Then, uniformly with
respect to $h\in\mathcal{H}_{sc,n}$, $\boldsymbol{d}\in\mathcal{D}$ and $s\in
S_{n}$,

$\displaystyle DM_{\boldsymbol{\beta}}$ $\displaystyle-$ $\displaystyle
n^{-3/2}\boldsymbol{A}_{n}^{T}\Bigg{(}\left(\boldsymbol{0}_{p\times
1},\boldsymbol{I}_{p\times
p}\right)^{T}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\boldsymbol{R}^{T}\left(\boldsymbol{R}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\boldsymbol{R}^{T}\right)^{-1}\boldsymbol{R}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\left(\boldsymbol{0}_{p\times
1},\boldsymbol{I}_{p\times p}\right)$
$\displaystyle\quad-\boldsymbol{V}_{R}(\boldsymbol{d})^{T}\boldsymbol{V}_{R}(\boldsymbol{d})E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n,R}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]$
$\displaystyle\quad+\boldsymbol{V}(\boldsymbol{d})^{-1}(1,\boldsymbol{0}_{p}^{T})^{T}(1,\boldsymbol{0}_{p}^{T})\boldsymbol{V}(\boldsymbol{d})^{-1}E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]\Bigg{)}\boldsymbol{A}_{n}n^{-3/2}=o_{\mathbb{P}}(1),$

under $H_{0}$ and $\mathbb{P}(n^{-1}DM_{\boldsymbol{\beta}}>c)\rightarrow 1$
for any $c>0$ if $H_{0}$ does not hold.

The process in Proposition 4.2 is asymptotically tight and for each
$\boldsymbol{d}$ behaves asymptotically as a weighted sum of $p+1-r$
independent chi-squares, where the weights are the positive eigenvalues of

$\displaystyle\left(\boldsymbol{0}_{p\times 1},\boldsymbol{I}_{p\times
p}\right)^{T}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\boldsymbol{R}^{T}\left(\boldsymbol{R}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\boldsymbol{R}^{T}\right)^{-1}\boldsymbol{R}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\left(\boldsymbol{0}_{p\times
1},\boldsymbol{I}_{p\times p}\right)\boldsymbol{\Delta}(\boldsymbol{d})$
$\displaystyle\quad-\boldsymbol{V}_{R}(\boldsymbol{d})^{T}\boldsymbol{V}_{R}(\boldsymbol{d})\boldsymbol{\Delta}(\boldsymbol{d},\boldsymbol{d})E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n,R}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]$
$\displaystyle\quad+\boldsymbol{V}(\boldsymbol{d})^{-1}(1,\boldsymbol{0}_{p}^{T})^{T}(1,\boldsymbol{0}_{p}^{T})\boldsymbol{V}(\boldsymbol{d})^{-1}\boldsymbol{\Delta}(\boldsymbol{d})E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right],$

see Johnson et al. (1995). Determining critical values requires the estimation
of the last expression. We can use the estimators stated in (14) and for all
other components, we simply replace the unknown expressions by their sample
version, e.g. estimate ${\mathbb{B}}_{n,R}$ by

$\widehat{\mathbb{B}}_{n}+{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\left(\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\right)^{-1}\boldsymbol{R}^{T}\left(\boldsymbol{R}\left(\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\right)^{-1}\boldsymbol{R}^{T}\right)^{-1}\boldsymbol{R}\left(\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\right)^{-1}\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}.$

###  4.3 Testing the transformation parameter and the slope coefficients

Finally, we consider the combined restrictions for $\boldsymbol{\beta}$ and
$\lambda$. Suppose we want to test

$\displaystyle H_{0}:\boldsymbol{R}\boldsymbol{\beta}_{0}=\boldsymbol{c}\quad
and\quad\lambda_{0}=\lambda_{R}.$

In contrast to the hypothesis stated in (16) we do not need to estimate
$\widehat{\lambda}_{R}$. Therefore, the distance metric statistic is for this
case given by

$\displaystyle
DM_{\boldsymbol{\beta},\lambda}=\frac{1}{n}\left(\widehat{\mathbb{Y}}_{n}(\lambda_{R})-\widehat{\mathbb{X}}_{n}\widehat{\boldsymbol{\beta}}_{R}(\lambda_{R})\right)^{T}{\mathbb{D}}_{n}\,\left(\widehat{\mathbb{Y}}_{n}(\lambda_{R})-\widehat{\mathbb{X}}_{n}\widehat{\boldsymbol{\beta}}_{R}(\lambda_{R})\right)-\frac{1}{n}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda})^{T}\;\widehat{\mathbb{B}}_{n}\widehat{\mathbb{Y}}_{n}(\widehat{\lambda}).$

In addition, $DM_{\boldsymbol{\beta},\lambda}$ does not converge to the same
expression as $DM_{\boldsymbol{\beta}}$ as $\lambda_{R}$ is fixed. Therefore,
we state the following proposition.

###### Proposition 4.3

Assume that the conditions of Proposition 3.1 hold true. Then, uniformly with
respect to $h\in\mathcal{H}_{sc,n}$, $\boldsymbol{d}\in\mathcal{D}$ and $s\in
S_{n}$,

$\displaystyle DM_{\boldsymbol{\beta},\lambda}$ $\displaystyle-$
$\displaystyle
n^{-3/2}\boldsymbol{A}_{n}^{T}\Bigg{(}\left(\boldsymbol{0}_{p\times
1},\boldsymbol{I}_{p\times
p}\right)^{T}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\boldsymbol{R}^{T}\left(\boldsymbol{R}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\boldsymbol{R}^{T}\right)^{-1}\boldsymbol{R}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\left(\boldsymbol{0}_{p\times
1},\boldsymbol{I}_{p\times p}\right)$
$\displaystyle\quad+\boldsymbol{V}(\boldsymbol{d})^{-1}(1,\boldsymbol{0}_{p}^{T})^{T}(1,\boldsymbol{0}_{p}^{T})\boldsymbol{V}(\boldsymbol{d})^{-1}E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right]\Bigg{)}\boldsymbol{A}_{n}n^{-3/2}=o_{\mathbb{P}}(1),$

under $H_{0}$ and
$\mathbb{P}(n^{-1}DM_{\boldsymbol{\beta},\lambda}>c)\rightarrow 1$ for any
$c>0$ if $H_{0}$ does not hold.

The process in Proposition 4.3 is asymptotically tight and for each
$\boldsymbol{d}$ behaves asymptotically as a weighted sum of $p-r$ independent
chi-squares, where the weights are the positive eigenvalues of

$\displaystyle\left(\boldsymbol{0}_{p\times 1},\boldsymbol{I}_{p\times
p}\right)^{T}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\boldsymbol{R}^{T}\left(\boldsymbol{R}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\boldsymbol{R}^{T}\right)^{-1}\boldsymbol{R}E\left[{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right]^{-1}\left(\boldsymbol{0}_{p\times
1},\boldsymbol{I}_{p\times p}\right)\boldsymbol{\Delta}(\boldsymbol{d})$
$\displaystyle\quad+\boldsymbol{V}(\boldsymbol{d})^{-1}(1,\boldsymbol{0}_{p}^{T})^{T}(1,\boldsymbol{0}_{p}^{T})\boldsymbol{V}(\boldsymbol{d})^{-1}\boldsymbol{\Delta}(\boldsymbol{d})E\left[\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})^{T}{\mathbb{B}}_{n}\frac{\partial}{\partial\lambda}{\mathbb{Y}}_{n}(\lambda_{0})\right],$

see Johnson et al. (1995). Determining critical values requires the estimation
of the last expression. We can use the estimators stated in (14) and for all
other components, we simply replace the unknown expressions by their sample
version.

The Propositions of Section 4 are also valid if the unknown parameters
$\lambda$ and $\boldsymbol{\beta}$ are estimated without the intercept
nuisance parameter $\gamma$. In that case ${\mathbb{D}}_{n}$ is replaced by
$\boldsymbol{\Omega}_{n}$ in the statements. Moreover, when estimating the
unknown variance ${\mathbb{D}}_{n,inf}$ has to be replaced by
$\boldsymbol{I}_{n\times n}$.

##  5 Small sample study and real data application

In this section we consider the small sample behavior of our estimator. We
conduct several simulation experiments to consider bias and standard deviation
for the estimated parameters. In addition, we conduct hypothesis tests as
discussed in Section 4.

In the proof of Theorem 3.2 it was established that

$\displaystyle\left((\widehat{\lambda},\widehat{\boldsymbol{\beta}}(\widehat{\lambda})^{T})^{T}-(\lambda_{0},\boldsymbol{\beta}_{0}^{T})^{T}\right)$
$\displaystyle=-\boldsymbol{V}(\boldsymbol{d})^{-1}\Bigg{(}\frac{1}{n}\sum_{j=1}^{n}\varepsilon_{j}{f}_{z}(\boldsymbol{Z}_{j})E\left[\boldsymbol{\tau}_{i}(\boldsymbol{d})\,\boldsymbol{\Omega}_{n,ij}(\boldsymbol{d})\mid\boldsymbol{X}_{j},\boldsymbol{Z}_{j}\right]$
$\displaystyle\quad-\frac{1}{n}\sum_{k=1}^{n}\varepsilon_{k}f_{z}(\boldsymbol{Z}_{k})E\left[\boldsymbol{\tau}_{i}(\boldsymbol{d})\boldsymbol{\Omega}_{n,ik}^{Z}\boldsymbol{\Omega}_{n,ij}^{X}\mid\boldsymbol{Z}_{k}\right]\Bigg{)}+o_{\mathbb{P}}\left(n^{-1/2}\right).$

The second sum in this asymptotic representation is due to the estimation of
$\boldsymbol{\eta}_{\lambda}$. In order to propose a simpler procedure, in our
simulation experiments we also investigated what happens when one neglects the
second part in the asymptotic representation, that is abusively consider
$E\left[\boldsymbol{\tau}_{i}(\boldsymbol{d})\boldsymbol{\Omega}_{n,ik}^{Z}\boldsymbol{\Omega}_{n,ij}^{X}\mid\boldsymbol{Z}_{k}\right]=\boldsymbol{0}_{p+1}$.
The estimator is labeled SmoothMD* in this section. The reason for this
investigation is that
$E\left[\boldsymbol{\tau}_{i}(\boldsymbol{d})\boldsymbol{\Omega}_{n,ij}^{Z}\boldsymbol{\Omega}_{n,ij}^{X}\mid\boldsymbol{Z}_{k}\right]$
is indeed null. It appears that considering an intercept $\gamma$ compensates
for this _ad-hoc_ simplification and allows reasonably accurate results to be
obtained.333Note that when estimating the model without constant one has to
replace $\boldsymbol{\tau}_{i}(\boldsymbol{d})$ by
$\widetilde{\boldsymbol{\tau}}_{i}(\boldsymbol{d}):=\left(\frac{\partial}{\partial\lambda}\mathbb{Y}_{n,i},-\mathbb{X}_{n,i}^{T}\right)^{T},$
but
$E\left[\widetilde{\boldsymbol{\tau}}_{i}(\boldsymbol{d})\boldsymbol{\Omega}_{n,ij}^{Z}\boldsymbol{\Omega}_{n,ij}^{X}\mid\boldsymbol{Z}_{k}\right]\neq\boldsymbol{0}_{p+1}$.

###  5.1 Simulation setup

During the simulation, we consider four different models. The models are given
by

  1. Model 1:

$T(Y,\lambda_{0})=X\beta_{0}+m(Z)+\varepsilon$,
$m(Z)=\frac{\exp\\{Z\\}}{1+\exp\\{Z\\}}+\frac{1}{3}$ with $Z\sim N(1,1)$,
$\lambda_{0}=0$ and $\beta_{0}=1$, $X=-\frac{2}{3}Z+u$ with $u\sim N(0,1)$ and
$\varepsilon=\sqrt{\frac{1+X^{2}}{2}}\;\widetilde{u}$ with $\widetilde{u}\sim
N\left(0,\frac{1}{13}\right)$.

  2. Model 2:

$T(Y,\lambda_{0})=X\beta_{0}+m(Z)+\varepsilon$,
$m(Z)=\frac{\exp\\{Z\\}}{1+\exp\\{Z\\}}+3$ with $Z\sim N(1,1)$,
$\lambda_{0}=0.5$ and $\beta_{0}=1$, $X=-\frac{2}{3}Z+u$ with $u\sim N(0,1)$
and $\varepsilon\sim N\left(0,\frac{1}{9}\right)$.

  3. Model 3:

$T(Y,\lambda_{0})=X\beta_{0}+m(Z)+\varepsilon$,
$m(Z)=\frac{\exp\\{Z\\}}{1+\exp\\{Z\\}}-1$ with $Z\sim U(-3,-1)$,
$\lambda_{0}=-1$ and $\beta_{0}=1$, $X=\frac{2}{3}Z+u$ with $u\sim U(-1,1)$
and $\varepsilon\sim U\left(-\sqrt{1/9},\sqrt{1/9}\right)$.

  4. Model 4:

$T(Y,\lambda_{0})=X_{1}\beta_{10}+\boldsymbol{X}_{2}\boldsymbol{\beta}_{20}+m(Z_{1},Z_{2})+\varepsilon$,
$m(Z_{1},Z_{2})=\frac{1}{3}+Z_{1}+Z_{2}+Z_{1}Z_{2}$ with $Z_{1},Z_{2}\sim
N(0,1)$, $\lambda_{0}=0$, $\beta_{10}=1$,
$X_{1}=-\frac{1}{3}\left(Z_{1}+Z_{2}\right)+u$ with $u\sim N(0,1)$,
$X_{2,l}\overset{i.i.d.}{\sim}Ber(0.2)$ and
$\beta_{2,l}\overset{i.i.d.}{\sim}U(-1,1)$ for $l=1,\ldots,30$ ,
$\varepsilon\sim N\left(0,\frac{1}{9}\right)$.

The main difference in the models is the transformation parameter $\lambda$.
Model 1 and Model 4 have $\lambda_{0}=0$, whereas Model 2 has
$\lambda_{0}=0.5$ and Model 3 $\lambda_{0}=-1$. To ensure that $Y>0$ in Model
3 we draw the random variables from uniform distributions. In all other models
positivity of $Y$ is ensured as well. Model 1 has heteroskedastic error terms
which are captured by the developed theory. Model 4 contains 30 dummy
variables, $\boldsymbol{X}_{2}$, which take the value 1 with probability
$20\%$. In addition, Model 4 has the same structure as the model of the
application we consider in section 5.3.

The estimators are computed by employing a normal kernel for $K(\cdot)$.
$\boldsymbol{Z}$ is standardized componentwise by the corresponding standard
deviations and $h\propto n^{-1/3.5}$. This bandwidth choice satisfies the
assumptions of Theorem 3.2. The components of $\boldsymbol{d}$ defining the
diagonal matrix $\boldsymbol{D}$ in $\boldsymbol{\Omega}_{n}$, are set equal
to the componentwise standard deviations of $X$ and $\boldsymbol{Z}$ when $X$
is continuous. In the case of the dummy variables $\boldsymbol{X}_{2}$, an
indicator of the event that the observations have the same value, is employed.
For Model 4, we ensure in the simulations that for every observation there
exists at least 4 observations with the same dummy variable combination.

In the estimation, we define a grid for values of $\lambda$ that are
considered during the optimization. This optimization grid for $\lambda$ is
given in our simulation by the grid $[\lambda_{0}-0.8,\lambda_{0}+0.8]$ with
step size $0.001$. We minimize
$G_{n}^{-\lambda}\widehat{\mathbb{Y}}_{n}(\lambda)^{T}\;\widehat{\mathbb{B}}_{n}\;G_{n}^{-\lambda}\widehat{\mathbb{Y}}_{n}(\lambda)$
over the defined grid to get $\widehat{\lambda}$ and
$\boldsymbol{\widehat{\beta}}(\widehat{\lambda})$, where
$G_{n}=\prod\limits_{i=1}^{n}Y_{i}^{1/n}$ is the sample geometric mean.

In the simulation, we compare the proposed estimator where $\gamma$ is
employed with the estimator that does not use $\gamma$. Both estimators
converge asymptotically to a normal distribution. The only difference is that
we have to replace $\boldsymbol{\tau}_{i}(\boldsymbol{d})$ by
$\widetilde{\boldsymbol{\tau}}_{i}(\boldsymbol{d}):=\left(\frac{\partial}{\partial\lambda}\mathbb{Y}_{n,i},-\mathbb{X}_{n,i}^{T}\right)^{T}$
in the case of the estimator without $\gamma$. It is therefore interesting to
compare both estimators.

We consider the bias and standard deviation of the estimators as well as the
power and size of the distance metric statistics proposed in Section 4. In
addition, we test by a simple Z-Test, if the estimated parameters are
significantly different from the true value. Therefore, we employ the variance
estimator stated in equation (14), and the necessary adjustments for the
estimator without $\gamma$ are replacing $\mathbb{D}_{n,inf}$ by
$\boldsymbol{I}_{n\times n}$ and $\mathbb{D}_{n}$ by
$\boldsymbol{\Omega}_{n}$. To estimate the error variance, we employ the
Eiker-White variance estimator. In order to see the influence of the estimated
$\boldsymbol{\eta}_{\lambda}$ on the variance, we consider all tests also
without taking the estimation error of $\boldsymbol{\eta}_{\lambda}$ into
account. Therefore, we replace $\widehat{\boldsymbol{\Phi}}_{n}$ by
$\boldsymbol{\Omega}_{n}$ in the variance estimator. As mentioned before, we
label this estimator SmoothMD*.

In addition, the Nonlinear two-stage Least Squares (NL2SLS) estimator for the
Box-Cox model introduced by Amemiya and Powell (1981) is considered as a
competitor. In order to be able to employ this estimator, it is assumed that
the function $m(\cdot)$ is known and, thus, $m(\boldsymbol{Z})$ can be added
as additional regressor. The instruments are given therefore by
$\boldsymbol{V}_{i}=(1,\boldsymbol{X}_{i},\boldsymbol{X}_{i}^{2},m(\boldsymbol{Z}_{i}),m(\boldsymbol{Z}_{i})^{2})$.
We consider the Z-Test for the NL2SLS estimator as well where we employ the
Eiker-White variance estimator again.

###  5.2 Simulation results

Table 1 states the results for the bias and standard deviation for $\lambda$
and $\beta$ in Model 1. All three estimators have comparable results for the
bias and the bias decreases with sample size for $\beta$, whereas it is the
lowest for $n=500$ for the SmoothMD estimators in the case of $\lambda$.
Surprisingly, the standard deviation is also comparable for all three
estimators even though $m(\cdot)$ is given for the NL2SLS estimator.

Table 1: Bias and Standard Deviation of the estimators for $\lambda$ and
$\beta$ in Model 1.

| $s$ | Bias | St. dev.
---|---|---|---
(lr)3-5 (lr)6-8 $n$ |  | 250 | 500 | 1000 | 250 | 500 | 1000
$\lambda$ estimator |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $0.003$ | $0.0001$ | $0.001$ | $0.042$ | $0.03$ | $0.021$
SmoothMD without $\gamma$ | $G_{n}$ | $0.002$ | $0.0001$ | $0.001$ | $0.041$ | $0.029$ | $0.02$
NL2SLS | $G_{n}$ | $-0.003$ | $-0.001$ | $0.0001$ | $0.042$ | $0.029$ | $0.02$
$\beta$ estimator |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $-0.001$ | $-0.001$ | $0.0004$ | $0.036$ | $0.025$ | $0.017$
SmoothMD without $\gamma$ | $G_{n}$ | $-0.001$ | $-0.001$ | $0.0004$ | $0.035$ | $0.024$ | $0.017$
NL2SLS | $G_{n}$ | $-0.002$ | $-0.001$ | $0.0002$ | $0.035$ | $0.024$ | $0.016$
|  |  |  |  |  |  |

Notes: For the SmoothMD estimators, $h\propto n^{-1/3.5}$. The components of
$\boldsymbol{d}$ are set equal to the componentwise standard deviations for
all variables. The grid for $\lambda$ is $[\lambda_{0}-0.8,\lambda_{0}+0.8]$.
2000 Monte Carlo samples were used for all simulations.

Table 2 states the results for the bias and standard deviation for $\lambda$,
$\beta_{1}$ and $\beta_{2}$, one representative parameter out of the 30
parameters in $\boldsymbol{\beta}_{2}$, in Model 4. All three estimators have
comparable results for the bias and the bias decreases with sample size for
all three parameters. In contrast to the results for Model 1, the standard
deviation is smaller in the case of the NL2SLS estimator for $\beta_{1}$ and
$\beta_{2}$. This result should be expected as $m(\cdot)$ is given for the
NL2SLS estimator. The standard deviations for the SmoothMD estimators with and
without $\gamma$ are as in Model 1, nearly the same.

Table 2: Bias and Standard Deviation of the estimators for $\lambda$,
$\beta_{1}$ and $\beta_{2}$ in Model 4.

| $s$ | Bias | St. dev.
---|---|---|---
(lr)3-5 (lr)6-8 $n$ |  | 250 | 500 | 1000 | 250 | 500 | 1000
$\lambda$ estimator |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $0.0001$ | $-0.0002$ | $-0.0001$ | $0.015$ | $0.01$ | $0.008$
SmoothMD without $\gamma$ | $G_{n}$ | $0.0001$ | $-0.0002$ | $-0.0001$ | $0.015$ | $0.01$ | $0.008$
NL2SLS | $G_{n}$ | $-0.0004$ | $-0.0002$ | $-0.0002$ | $0.014$ | $0.01$ | $0.005$
$\beta_{1}$ estimator |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $-0.002$ | $-0.002$ | $-0.001$ | $0.036$ | $0.023$ | $0.017$
SmoothMD without $\gamma$ | $G_{n}$ | $-0.002$ | $-0.002$ | $-0.001$ | $0.036$ | $0.023$ | $0.017$
NL2SLS | $G_{n}$ | $0.0004$ | $-0.0001$ | $-0.0001$ | $0.025$ | $0.015$ | $0.011$
$\beta_{2}$ estimator |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $0.004$ | $-0.001$ | $-0.0002$ | $0.133$ | $0.065$ | $0.042$
SmoothMD without $\gamma$ | $G_{n}$ | $0.004$ | $-0.001$ | $-0.0002$ | $0.133$ | $0.065$ | $0.042$
NL2SLS | $G_{n}$ | $-0.001$ | $0.002$ | $-0.0006$ | $0.091$ | $0.046$ | $0.028$
|  |  |  |  |  |  |

Notes: For the SmoothMD estimators, $h\propto n^{-1/3.5}$. The components of
$\boldsymbol{d}$ are set equal to the componentwise standard deviations for
all continuous variables and for the dummy variables an indicator of the event
that the observations have the same value is employed. The grid for $\lambda$
is $[\lambda_{0}-0.8,\lambda_{0}+0.8]$. 2000 Monte Carlo samples were used for
all simulations.

Table 3: Empirical Level for distance metric statistics of the estimators for
$\lambda$ and $\beta$ in Model 2.

| $s$ | 5% level | 10% level
---|---|---|---
(lr)3-5 (lr)6-8 $n$ |  | 250 | 500 | 1000 | 250 | 500 | 1000
Test for $\lambda$ |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $10.3$ | $8.15$ | $7.8$ | $12.75$ | $10.65$ | $11.45$
SmoothMD* with $\gamma$ | $G_{n}$ | $10.15$ | $7.95$ | $7.7$ | $12.75$ | $10.6$ | $11.45$
SmoothMD without $\gamma$ | $G_{n}$ | $9.55$ | $7.0$ | $6.0$ | $12.45$ | $10.85$ | $10.35$
Test for $\beta$ |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $10.4$ | $8.25$ | $7.95$ | $12.55$ | $11.1$ | $11.9$
SmoothMD* with $\gamma$ | $G_{n}$ | $10.2$ | $8.3$ | $7.7$ | $12.55$ | $11.3$ | $11.1$
SmoothMD without $\gamma$ | $G_{n}$ | $9.7$ | $6.85$ | $6.7$ | $12.55$ | $10.2$ | $11.2$
|  |  |  |  |  |  |

Notes: For the SmoothMD estimators, $h\propto n^{-1/3.5}$. The components of
$\boldsymbol{d}$ are set equal to the componentwise standard deviations for
all variables. The variances are estimated by the Eiker-White variance
estimator. For SmoothMD* the additional variance part due to the estimation of
$\boldsymbol{\eta_{\lambda}}$ is not taken into account. For SmoothMD the
additional variance part is taken into account. 2000 Monte Carlo samples were
used for all simulations.

Table 3 states the empirical level for distance metric statistics of the
estimators for $\lambda$ and $\beta$ in Model 2. Here we state the results for
the SmoothMD estimator with correctly estimated variance as well as with the
variance estimate that does not account for the estimation of
$\boldsymbol{\eta}_{\lambda}$. For all three estimators, the empirical levels
converge to the nominal levels if the sample size increases and $\beta$ seems
to need a larger sample size than $\lambda$ to get close to the nominal level.
However, in this setup, the results of the SmoothMD estimators with and
without $\gamma$ differ. In addition, it can be seen that the estimator
SmoothMD* leads to almost the same results as SmoothMD.

Table 4 states the empirical level for the Z-Tests for $\lambda$, $\beta_{1}$
and $\beta_{2}$ in Model 4. The fact that we do not consider the estimation
error has almost no influence on the results. In addition, both SmoothMD
versions lead to similar results. However, in order to get close to the
nominal level, the sample size needs to be large as only for $n=1000$ do the
SmoothMD estimators get close to the nominal level. The NL2SLS estimator gives
more convincing results for smaller sample sizes. Note that the dummy variable
coefficient $\beta_{2}$ seems to require a larger sample size than the other
two parameters to get close to the nominal level when employing the SmoothMD
estimators.

Table 4: Empirical Level for Z-Tests of the estimators for $\lambda$,
$\beta_{1}$ and $\beta_{2}$ in Model 4.

| $s$ | 5% level | 10% level
---|---|---|---
(lr)3-5 (lr)6-8 $n$ |  | 250 | 500 | 1000 | 250 | 500 | 1000
Test for $\lambda$ |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $8.4$ | $7.1$ | $5.5$ | $15.45$ | $13.5$ | $10.55$
SmoothMD* with $\gamma$ | $G_{n}$ | $8.5$ | $7.1$ | $5.45$ | $15.3$ | $13.45$ | $10.55$
SmoothMD without $\gamma$ | $G_{n}$ | $9.6$ | $8.15$ | $5.2$ | $16.45$ | $12.6$ | $11.65$
NL2SLS | $G_{n}$ | $9.6$ | $7.6$ | $6.0$ | $16.8$ | $12.5$ | $11.1$
Test for $\beta_{1}$ |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $11.8$ | $8.55$ | $6.4$ | $18.25$ | $14.05$ | $11.75$
SmoothMD* with $\gamma$ | $G_{n}$ | $11.8$ | $8.55$ | $6.45$ | $18.25$ | $14.05$ | $11.8$
SmoothMD without $\gamma$ | $G_{n}$ | $11.75$ | $8.9$ | $6.4$ | $18.35$ | $14.8$ | $12.1$
NL2SLS | $G_{n}$ | $8.25$ | $4.95$ | $5.1$ | $13.95$ | $10.4$ | $10.35$
Test for $\beta_{2}$ |  |  |  |  |  |  |
SmoothMD with $\gamma$ | $G_{n}$ | $13.6$ | $8.25$ | $7.05$ | $20.7$ | $15.65$ | $12.7$
SmoothMD* with $\gamma$ | $G_{n}$ | $13.7$ | $8.25$ | $7.05$ | $20.75$ | $15.55$ | $12.75$
SmoothMD without $\gamma$ | $G_{n}$ | $13.75$ | $8.35$ | $6.65$ | $20.6$ | $14.55$ | $12.25$
NL2SLS | $G_{n}$ | $7.65$ | $6.55$ | $4.7$ | $13.45$ | $13.05$ | $9.55$
|  |  |  |  |  |  |

Notes: For the SmoothMD estimators, $h\propto n^{-1/3.5}$. The components of
$\boldsymbol{d}$ are set equal to the componentwise standard deviations for
all continuous variables and for the dummy variables, an indicator of the
event that the observations have the same value is employed. The variances are
estimated by the Eiker-White variance estimator. For SmoothMD* the additional
variance part due to the estimation of $\boldsymbol{\eta_{\lambda}}$, is not
taken into account. For SmoothMD, the additional variance part is taken into
account. $\beta_{2}$ is one representative parameter out of the 30 parameters
in $\boldsymbol{\beta}_{2}$. 2000 Monte Carlo samples were used for all
simulations.

Figure 1: Power function of the distance metric
statistic for $\lambda$ of Model 3 with $n=250$.

Figure 2: Power function of the distance metric
statistic for $\lambda$ of Model 3 with $n=1000$.

Figures 2 and 2 state the power functions of the distance metric statistic for
$\lambda$ in Model 3 with $n=250$ and $n=1000$. In the case of $n=250$, the
power function is skewed and the power for values larger than $-1$ is small.
In addition, the power function is smaller than the nominal value at $-0.85$
and $-0.7$. For the SmoothMD estimator without $\gamma$, the power function is
larger than for the SmoothMD estimator with $\gamma$ at values larger than
$-1$. These issues disappear for the larger sample size $n=1000$.

Figures 4 and 4 state the power functions of the distance metric statistic for
$\beta$ in Model 2 with $n=250$ and $n=500$. As in Figure 2 the power function
for $n=250$ is skewed but the effect is less distinct. However, the power
function is smaller than the nominal value at $0.8$. For the SmoothMD
estimator without $\gamma$, the power function is larger than for the SmoothMD
estimator with $\gamma$ at values smaller than $1$ for both sample sizes. For
$n=500$, the skewness is less pronounced and the power function has no values
lower than the nominal value. The main conclusion from both power functions is
that the samples size should not be too small so that the tests have a
reasonable power.

Figure 3: Power function of the distance metric
statistic for $\beta$ of Model 2 with $n=250$.

Figure 4: Power function of the distance metric
statistic for $\beta$ of Model 2 with $n=500$.

Before we consider a real data application, we close the discussion with
Figures 6 and 6. The figures state the estimated $m(Z)$ for Model 1 with
$n=250$ and for Model 3 with $n=500$. For the estimation, the NW estimator
with same kernel and bandwidth as for the SmoothMD estimator was used. No
matter if the SmoothMD estimator with or without $\gamma$ is employed, the
results are very accurate. In practice one can of course use cross validation
to choose the bandwidth or employ the local linear estimator instead of the NW
estimator.

Figure 5: Estimated $m(Z)$ for Model 1 with $n=250$.

Figure 6: Estimated $m(Z)$ for Model 3
with $n=500$.

###  5.3 Real data application

We consider in this section an application of our estimator to investigate the
returns of social and cognitive skills in the labor market. For this purpose,
we apply the proposed transformation partially linear estimator to a dataset
studied in Deming (2017). In particular, we consider the regressions (4) and
(5) in TABLE I of Deming (2017) that are based on the National Longitudinal
Survey of Youth 1979 (NLSY79). NLSY79 is a nationally representative sample
taken in the US, of young people aged from 14 to 22. The survey was conducted
yearly from 1979 to 1993 and biannually from 1994 through 2012. Deming (2017)
estimates the model

$\displaystyle\log(wage_{ijt})$ $\displaystyle=\alpha+\beta_{1}\cdot
COG_{i}+\beta_{2}\cdot SS_{i}+\beta_{3}\cdot COG_{i}\times
SS_{i}+\beta_{4}\cdot NCOG_{i}$ (18)
$\displaystyle\quad+\boldsymbol{C}_{ijt}^{T}\boldsymbol{\rho}+\delta_{j}+\zeta_{t}+\varepsilon_{ijt},$

where $COG$, $SS$ and $NCOG$ denote measures of cognitive, social and
noncognitive skills. The model includes controls $\boldsymbol{C}$ for race-by-
gender indicators, indicators for region and urbanicity as well as age
(indexed by $j$) and year (indexed by $t$) fixed effects.

In his paper, Deming (2017) develops a theoretical model that is written in
levels instead of logs as in equation (18). Nevertheless, he estimates the
log-linearized model in his paper to follow standard practice in the
literature, as he argues. Results for the model in levels are stated in an
online appendix. Therefore, it makes sense to use the Box-Cox transformation
for $wage$ and estimate the transformation parameter $\lambda$ together with
the remaining model parameters to decide whether the model in logs or in
levels is more appropriate.

Furthermore, we consider an unknown functional form for cognitive and social
skills to see if the linear form, $\beta_{1}\cdot COG+\beta_{2}\cdot
SS+\beta_{3}\cdot COG\times SS$ used by Deming (2017), is reasonable. The
transformation partially linear model is, thus, given by

$\displaystyle T(wage_{ijt},\lambda)=m(COG_{i},SS_{i})+\beta\cdot
NCOG_{i}+\boldsymbol{C}_{ijt}^{T}\boldsymbol{\rho}+\delta_{j}+\zeta_{t}+\varepsilon_{ijt},$
(19)

where $m(\cdot)$ is an unknown function. The model stated in (19) is closely
related to Model 4 of the simulations in section 5.1, where $Y=wage$,
$\boldsymbol{Z}=(COG,SS)^{T}$, $X_{1}=NCOG$ and
$X_{2}=(\boldsymbol{C}^{T},1,1)^{T}$.

As proxy for cognitive skills, the Armed Forces Qualifying Test (AFQT) was
taken. Deming (2017) uses raw scores from Altonji et al. (2012) and normalizes
them to have mean 0 and standard deviation 1. The social skill measure is
constructed from the following four variables of the NLSY79:

  1. 1.

Self-reported sociability in 1981 (extremely shy, somewhat shy, somewhat
outgoing, extremely outgoing)

  2. 2.

Self-reported sociability in 1981 at age 6 (retrospective)

  3. 3.

The number of clubs in which the respondent participated in high school

  4. 4.

Participation in high school sports (yes/no).

Each variable is normalized to have mean 0 and standard deviation 1. The
social skill measure is the average of these four normalized variables (also
normalized to standard deviation 1). In addition to social and cognitive skill
measures Deming (2017) includes a noncognitive skill measure in his
regression. He uses the Rotter Locus of Control and the Rosenberg Self-Esteem
Scale as also used by Heckman et al. (2006). In the following discussion we
use these variables to estimate the models stated in (18) and (19).

Deming (2017) used a weighted log-linearized OLS estimator to estimate the
returns of cognitive and social skills on wage and excluded respondents under
the age of $23$ or who were enrolled in school. The weighting was necessary as
in each survey year of the NLSY79 a set of sampling weights was constructed.
These weights provided the researcher with an estimate of how many individuals
in the United States each respondent’s answers represented. We also employ
these weights in our analysis.

Table 5 shows the regression results. The first column, (4), provides the
results of Deming (2017) estimating equation (18). The second and third column
state the transformation partially linear estimator of equation (19) with and
without employing $\gamma$. The fourth and fifth column state the
transformation partially linear estimator of equation (19) with and without
employing $\gamma$ imposing $\lambda=0$. This is a standard partially linear
model as studied by Robinson (1988) and Li (1996).

For the inner smoothing of the estimations in column 2-5, we use a normal
kernel with $h\propto n^{-1/3.5}$, which is in line with the developed theory
and usual bandwidth choices for bivariate smoothing with the Nadaraya-Watson
estimator. The components of $\boldsymbol{d}$ defining the diagonal matrix
$\boldsymbol{D}$ in $\boldsymbol{\Omega}_{n}$ are set equal to the
componentwise standard deviations for all continuous variables. In the case of
the controls and fixed effects, an indicator of the event that the
observations have the same value is employed.

The optimization grid for $\lambda$ is given by $[-0.1,0.1]$ with step size
$0.001$.444We evaluated subsamples of the dataset before we made the final
estimation. The estimated $\lambda$’s in the subsamples are contained in the
used grid. We minimize
$G_{n}^{-\lambda}\widehat{\mathbb{Y}}_{n}(\lambda)^{T}\;\widehat{\mathbb{B}}_{n}\;G_{n}^{-\lambda}\widehat{\mathbb{Y}}_{n}(\lambda)$
over the defined grid to get $\widehat{\lambda}$ and the estimates of the
remaining coefficients, where $G_{n}$ is the sample geometric mean.

The results in the first column of Table 5 show that all Deming (2017)
estimated coefficients are significantly different from $0$. In the remaining
four columns, we cannot state parameter estimates for cognitive and social
skills and the interaction of both as these variables are contained in
$m(\cdot)$. However, the parameter estimates for noncognitive skills are
comparable to the estimate from the first column. In addition, the estimates
for $\lambda$ with and without $\gamma$ are equal and close to zero which
would imply that a log-transformation of the dependent variable is
appropriate. The estimated coefficient for noncognitive skills is
significantly different from $0$ in all SmoothMD estimations whereas both
estimates for $\lambda$ are not significantly different from $0$.

In order to check if the linear specification for cognitive and social skills
employed by Deming (2017) is reasonable, we proceed as follows. We estimate
the parameters of the transformation partially linear model as stated in (19)
to obtain the residuals

$\widehat{\varepsilon}_{ijt}=T(wage_{ijt},\widehat{\lambda})-\widehat{\beta}\cdot
NCOG_{i}-\boldsymbol{C}_{ijt}^{T}\widehat{\boldsymbol{\rho}}-\widehat{\delta}_{j}-\widehat{\zeta}_{t}.$

We now estimate the unknown function $m(\cdot)$ by smoothing
$\widehat{\varepsilon}$ with the NW estimator. In addition, we also regress
$\widehat{\varepsilon}$ on $COG$, $SS$ and $COG\times SS$. To see if the
linear specification is appropriate, we compare the MSE of the linear and
nonlinear estimates. We employ a normal density kernel for the NW estimator
and let $h\propto n^{-1/6}$.

Table 6 states the results where OLS indicates that we used the linear model
to fit the residuals. The MSE of the linear and nonlinear estimates are
identical no matter if we use the SmoothMD estimator with or without $\gamma$
to estimate the unknown model parameters. The same holds true for the SmoothMD
estimator with or without $\gamma$ where $\lambda=0$ is imposed. All results
show that the linear representation of Deming (2017) seems to be reasonable.

Table 5: Labor Market Returns to Cognitive and Social Skills in the NLSY79

Outcome: (log) hourly wage | (4) | SmoothMD | SmoothMD | SmoothMD | SmoothMD
---|---|---|---|---|---
(in 2012 dollars) |  | with $\gamma$ | without $\gamma$ | with $\gamma$, $\lambda=0$ | without $\gamma$, $\lambda=0$
$\lambda$ | - | -0.007 | -0.007 | - | -
|  | [0.005] | [0.005] |  |
Cognitive skills | 0.189*** | - | - | - | -
| [0.007] |  |  |  |
Social skills | 0.043*** | - | - | - | -
| [0.006] |  |  |  |
Cognitive $\times$ Social | 0.019*** | - | - | - | -
| [0.006] |  |  |  |
Noncognitive skills | 0.048*** | 0.047*** | 0.047*** | 0.048*** | 0.048***
| [0.006] | [0.004] | [0.004] | [0.004] | [0.004]
Demographics and age/ | X | X | X | X | X
year fixed effects |  |  |  |  |
Number of Observations | 126191 | 126191 | 126191 | 126191 | 126191
|  |  |  |  |

Notes: The data source is the National Longitudinal Survey of Youth 1979
cohort (NLSY79). (4) denotes the OLS regression proposed by Deming (2017). In
all SmoothMD estimations, $h\propto n^{-1/3.5}$. The components of
$\boldsymbol{d}$ are set equal to the componentwise standard deviations for
all continuous variables and for controls and fixed effects, an indicator of
the event that the observations have the same value is employed. The grid for
$\lambda$ is $[-0.1,0.1]$ and $s=G_{n}$. Cognitive skills are measured by each
NLSY79 respondent’s score on the Armed Forces Qualifiying Test (AFQT) and are
normalized to have mean 0 and standard deviation 1. The AFQT score crosswalk
of Altonji et al. (2012) is used. Social skill is a standardized composite of
four variables, (i) sociability in childhood, (ii) sociability in adulthood,
(iii) participation in high school clubs and (iv) participation in team
sports; see the text and Deming (2017) for details on the construction of the
social skills measure. The noncognitive skills measure is the normalized
average of the Rotter and Rosenberg scores in the NLSY. The regressions also
control for race-by-gender indicator variables, age, year, census region and
urbanicity. Standard errors are in brackets and are clustered at the
individual level for (4). The remaining standard errors are estimated by the
Eiker-White variance estimator. ***$p<.01$, **$p<.05$, *$p<.1$

Table 6: MSE of estimated nonlinear part in the transformation partially
linear model

| SmoothMD | SmoothMD | SmoothMD | SmoothMD
---|---|---|---|---
| with $\gamma$ | without $\gamma$ | with $\gamma$, $\lambda=0$ | without $\gamma$, $\lambda=0$
(lr)2-3 (lr)4-5 (lr)6-7 (lr)8-9 | OLS | NW | OLS | LL | OLS | NW | OLS | LL
MSE | 0.282 | 0.277 | 0.282 | 0.277 | 0.293 | 0.288 | 0.293 | 0.288
Number of Observations | 126191 | 126191 | 126191 | 126191 | 126191 | 126191 | 126191 | 126191
|  |  |  |  |  |  |  |

Notes: For the NW estimator a normal kernel with $h\propto n^{-1/6}$ is
employed. OLS indicates that the linear model is used to fit the residuals.

In a second step, we include years of completed education as an additional
explanatory variable in the regression models. In one of his estimations
Deming (2017) controls for years of education as well. Table 7 states the
regression results for all considered models. The first column, (5), provides
the results of Deming (2017) estimating equation (18) with years of completed
education as control. The results show that all Deming (2017) estimated
coefficients are significantly different from $0$. However, the coefficients
become smaller compared to the first specification. In addition, the
coefficient of the interactive effect is only significant at the $10\%$ level.
In the remaining four columns, the parameter estimates for noncognitive skills
are comparable to the estimate from the first column. In addition, the
estimates for $\lambda$ with and without $\gamma$ are equal and close to zero
which would imply that a log-transformation of the dependent variable is
appropriate. The estimated coefficient for noncognitive skills is
significantly different from $0$ in all SmoothMD estimations whereas both
estimates for $\lambda$ are not significantly different from $0$.

Table 8 states the MSE of the estimated nonlinear part in the transformation
partially linear models. The MSE of the linear and nonlinear estimates are
identical no matter if we use the SmoothMD estimator with or without $\gamma$
to estimate the unknown model parameters. The same holds true for the SmoothMD
estimator with or without $\gamma$ where $\lambda=0$ is imposed. All results
show that the linear representation of Deming (2017) seems to be reasonable.

Table 7: Labor Market Returns to Cognitive and Social Skills in the NLSY79
controlling for education

Outcome: (log) hourly wage | (5) | SmoothMD | SmoothMD | SmoothMD | SmoothMD
---|---|---|---|---|---
(in 2012 dollars) |  | with $\gamma$ | without $\gamma$ | with $\gamma$, $\lambda=0$ | without $\gamma$, $\lambda=0$
$\lambda$ | - | 0.002 | 0.002 | - | -
|  | [0.005] | [0.005] |  |
Cognitive skills | 0.126*** | - | - | - | -
| [0.008] |  |  |  |
Social skills | 0.029*** | - | - | - | -
| [0.006] |  |  |  |
Cognitive $\times$ Social | 0.011* | - | - | - | -
| [0.006] |  |  |  |
Noncognitive skills | 0.040*** | 0.037*** | 0.037*** | 0.037*** | 0.037***
| [0.006] | [0.004] | [0.004] | [0.004] | [0.004]
Demographics and age/ | X | X | X | X | X
year fixed effects |  |  |  |  |
Years of completed education | X | X | X | X | X
Number of Observations | 126191 | 126191 | 126191 | 126191 | 126191
|  |  |  |  |

Notes: The data source is the National Longitudinal Survey of Youth 1979
cohort (NLSY79). (5) denotes the OLS regression proposed by Deming (2017). In
all SmoothMD estimations, $h\propto n^{-1/3.5}$. The components of
$\boldsymbol{d}$ are set equal to the componentwise standard deviations for
all continuous variables and for controls and fixed effects, an indicator of
the event that the observations have the same value is employed. The grid for
$\lambda$ is $[-0.1,0.1]$ and $s=G_{n}$. Cognitive skills are measured by each
NLSY79 respondent’s score on the Armed Forces Qualifiying Test (AFQT) and are
normalized to have mean 0 and standard deviation 1. The AFQT score crosswalk
of Altonji et al. (2012) is used. Social skill is a standardized composite of
four variables (i) sociability in childhood, (ii) sociability in adulthood,
(iii) participation in high school clubs and (iv) participation in team
sports; see the text and Deming (2017) for details on the construction of the
social skills measure. The noncognitive skills measure is the normalized
average of the Rotter and Rosenberg scores in the NLSY. The regressions also
control for race-by-gender indicator variables, age, year, census region,
urbanicity and years of completed education. Standard errors are in brackets
and are clustered at the individual level for (5). The remaining standard
errors are estimated by the Eiker-White variance estimator. ***$p<.01$,
**$p<.05$, *$p<.1$

Before we close the section, we plot the estimated labor market returns to
cognitive and social skills of model (19) with and without controlling for
years of completed education. The returns are estimated with the NW estimator
employing a normal kernel with $h\propto n^{-1/6}$. Figures 8 and 8 present
the results. Note that the mean was subtracted. The return increases no matter
if the social or cognitive indicator is increased. However, the cognitive
effect seems to be stronger. In addition, if we control for education, it
seems that being too social might sometimes lower the wage somewhat.
Nevertheless, both plots confirm that the linear model used by Deming (2017)
is reasonable.

Table 8: MSE of estimated nonlinear part in the transformation partially
linear model controlling for education

| SmoothMD | SmoothMD | SmoothMD | SmoothMD
---|---|---|---|---
| with $\gamma$ | without $\gamma$ | with $\gamma$, $\lambda=0$ | without $\gamma$, $\lambda=0$
(lr)2-3 (lr)4-5 (lr)6-7 (lr)8-9 | OLS | NW | OLS | LL | OLS | NW | OLS | LL
MSE | 0.288 | 0.283 | 0.288 | 0.283 | 0.284 | 0.280 | 0.284 | 0.280
Number of Observations | 126191 | 126191 | 126191 | 126191 | 126191 | 126191 | 126191 | 126191
|  |  |  |  |  |  |  |

Notes: For the NW estimator a normal kernel with $h\propto n^{-1/6}$ is
employed. OLS indicates that the linear model is used to fit the residuals.

Figure 7: Estimated Labor Market Returns to
Cognitive and Social Skills in the NLSY79.

Figure 8: Estimated Labor Market Returns to
Cognitive and Social Skills in the NLSY79
controlling for education.

##  6 Discussion

In this paper, we study the semiparametric partially linear model with Box-Cox
transformed dependent variable. We follow the SmoothMD approach introduced by
Lavergne and Patilea (2013), which is based on conditional moment
restrictions, and we extend it to the case with an infinite-dimensional
nuisance parameter. Our results are new both for the transformation regression
models and the semiparametric partially linear model. We establish model
identification, consistency as well as $\sqrt{n}$-asymptotic normality. In
addition, we proposed a distance metric statistic to test the model
parameters. A Monte Carlo experiment showed the usefulness of the proposed
estimator in finite samples. An application to a large real data sample
already studied in Labor Economics is also reported.

The SmoothMD approach is a convenient approach for nonlinear regression
models. It could be interpreted as a generalized least-squares method where
the weights are given by a suitable positive-definite matrix with each entry
representing a measure of discrepancy between a pair of covariate vectors. It
is unnecessary to localize the measure of discrepancy between the observed
covariate vectors, and this represents a significant advantage for the
practitioner who thus avoids the choice of an additional tuning parameter. One
price to pay is on the semiparametric efficiency for the estimators of the
finite dimensional parameters. However, Lavergne and Patilea (2013) showed
that a two-step procedure, where the first step allows nonparametrically
suitable weights to be estimated for building an asymptotically optimal
measure of discrepancy, allows to achieve semiparametric efficiency. An
efficient estimator would also induce a distance metric test statistic with a
usual asymptotic chi-square distribution under the null hypothesis. We expect
that the same results extend to the present framework. However, with our
semiparametric model, the two-step procedure would likely result in a
numerically unstable, complex to calibrate, inference procedure. We believe
this theoretical, and quite technical, refinement to be of little use for the
applications, and we therefore do not consider it.

Another slight drawback for using SmoothMD with a fixed weighting matrix, is
the condition that $\boldsymbol{Z}$ should be of dimension $q$ less than or
equal to 3. This restriction would not be binding in most applications.
However, if necessary, one could use higher-order kernels to diminish the bias
induced by the nonparametric estimation of the nuisance parameter, and thus
allow for larger $q$. It is noticeable that SmoothMD does not involve any
denominator and therefore the higher-order kernels would not induce the
numerical problems, due to division by zero, usually encountered in
semiparametric methods.

## Acknowledgments

This research was supported by the DFG through KN 567/5-1 and the Hausdorff
Center for Mathematics. We furthermore thank the Regional Computing Center of
the University of Cologne (RRZK) for providing computing time on the DFG-
funded High Performance Computing (HPC) system CHEOPS as well as support. This
work was started while the first author was visiting CREST-ENSAI. V. Patilea
acknowledges support from ‘Models and mathematical processing of very large
data’, a Joint Research Initiative under the aegis of Risk Foundation, with
partnership of MEDIAMETRIE and GENES, and from the Romanian Minister of
Education and Research, CNCS – UEFISCDI, project number PN-III-P4-ID-
PCE-2020-1112, within PNCDI III.

## Appendix

### Appendix A: Proofs

Indeed, for any $\boldsymbol{u}\in\mathbb{R}^{p}$ and
$\boldsymbol{v}\in\mathbb{R}^{q}$ such that
$(\boldsymbol{u}^{T},\boldsymbol{v}^{T})^{T}\neq\boldsymbol{0}_{p+q},$ we can
write

$Var\left[\boldsymbol{u}^{T}\boldsymbol{X}+\boldsymbol{v}^{T}\boldsymbol{Z}\right]=E\left[Var\left[\boldsymbol{u}^{T}(\boldsymbol{X}-E\left[\boldsymbol{X}\mid\boldsymbol{Z}\right])\mid\boldsymbol{Z}\right]\right]+Var\left[\boldsymbol{u}^{T}E\left[\boldsymbol{X}\mid\boldsymbol{Z}\right]+\boldsymbol{v}^{T}\boldsymbol{Z}\right].$

If $\boldsymbol{u}\neq\boldsymbol{0}_{p}$,

$Var\left[\boldsymbol{u}^{T}\boldsymbol{X}+\boldsymbol{v}^{T}\boldsymbol{Z}\right]\geq
E\left[Var\left[\boldsymbol{u}^{T}(\boldsymbol{X}-E\left[\boldsymbol{X}\mid\boldsymbol{Z}\right])\mid\boldsymbol{Z}\right]\right]\\\
=\boldsymbol{u}^{T}E\left[Var\left[\boldsymbol{X}-E\left[\boldsymbol{X}\mid\boldsymbol{Z}\right]\mid\boldsymbol{Z}\right]\right]\boldsymbol{u}=\boldsymbol{u}^{T}Var\left[\boldsymbol{X}-E\left[\boldsymbol{X}\mid\boldsymbol{Z}\right]\right]\boldsymbol{u}>0,$

where the last inequality is guaranteed by Assumption 2.1.1. When
$\boldsymbol{u}=\boldsymbol{0}_{p}$, we obtain

$Var\left[\boldsymbol{u}^{T}\boldsymbol{X}+\boldsymbol{v}^{T}\boldsymbol{Z}\right]=Var\left[\boldsymbol{v}^{T}\boldsymbol{Z}\right]=\boldsymbol{v}^{T}Var\left[\boldsymbol{Z}\right]\boldsymbol{v}>0,$

where the last inequality holds because $\boldsymbol{v}\neq\boldsymbol{0}_{q}$
and $Var\left[\boldsymbol{Z}\right]$ has necessarily full rank provided
$\boldsymbol{Z}$ admits a density.

By construction, $\gamma(\lambda_{0})=0$ and
$\boldsymbol{\beta}(\lambda_{0})=\boldsymbol{\beta}_{0}$. Following the lines
of Shin (2008), for any $(\gamma,\lambda,\boldsymbol{\beta}^{T})^{T}$ we note
that

$\mathbb{P}\bigg{(}E\left[\left(T(Y,\lambda)-E[T(Y,\lambda)\mid\boldsymbol{Z}]\right)f_{z}(\boldsymbol{Z})-\gamma-(\boldsymbol{X}-E[\boldsymbol{X}\mid\boldsymbol{Z}])^{T}\boldsymbol{\beta}f_{z}(\boldsymbol{Z})\mid\boldsymbol{X},\boldsymbol{Z}\right]=0\bigg{)}\\\
=\mathbb{P}\bigg{(}E\left[T(Y,\lambda)-T(Y,\lambda_{0})|\boldsymbol{X},\boldsymbol{Z}\right]-\boldsymbol{X}^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})=E[T(Y,\lambda)\mid\boldsymbol{Z}]-E[T(Y,\lambda_{0})|\boldsymbol{Z}]\\\
+\gamma
f^{-1}_{z}(\boldsymbol{Z})-E[\boldsymbol{X}\mid\boldsymbol{Z}]^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})\bigg{)}.$

Hence, it suffices to prove that the last probability could not be equal to 1
when
$(\gamma,\boldsymbol{\theta}^{T})^{T}\neq(0,\boldsymbol{\theta}_{0}^{T})^{T}$.
Note that

$\displaystyle
E[T(Y,\lambda)\mid\boldsymbol{Z}]-E[T(Y,\lambda_{0})\mid\boldsymbol{Z}]+\gamma
f^{-1}_{z}(\boldsymbol{Z})-E[\boldsymbol{X}\mid\boldsymbol{Z}]^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0}),$

does not depend on $\boldsymbol{X}$ anymore but only on $\boldsymbol{Z}$.

If $\lambda=\lambda_{0}$ the result follows immediately from the full rank
condition in Assumption 2.1.1. Indeed, by the variance decomposition formula
and Assumption 2.1.1, for any $\boldsymbol{a}\in\mathbb{R}^{p}$,
$\boldsymbol{a}\neq\boldsymbol{0}_{p}$,

$\boldsymbol{a}^{T}Var(\boldsymbol{X}-E[\boldsymbol{X}|\boldsymbol{Z}])\boldsymbol{a}=E\left[Var(\boldsymbol{a}^{T}(\boldsymbol{X}-E[\boldsymbol{X}|\boldsymbol{Z}])\mid\boldsymbol{Z})\right]>0.$

This implies

$\boldsymbol{a}^{T}Var(f_{z}(\boldsymbol{Z})\left(\boldsymbol{X}-E[\boldsymbol{X}|\boldsymbol{Z}]\right))\boldsymbol{a}=E\left[f^{2}_{z}(\boldsymbol{Z})\boldsymbol{a}^{T}Var(\left(\boldsymbol{X}-E[\boldsymbol{X}|\boldsymbol{Z}]\right)\mid\boldsymbol{Z})\boldsymbol{a}^{T}\right]\\\
=E\left[f^{2}_{z}(\boldsymbol{Z})Var(\boldsymbol{a}^{T}\left(\boldsymbol{X}-E[\boldsymbol{X}|\boldsymbol{Z}]\right)\mid\boldsymbol{Z})\right]>0.$

Thus,
$f_{z}(\boldsymbol{Z})(\boldsymbol{X}-E[\boldsymbol{X}\mid\boldsymbol{Z}]^{T})^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})$
cannot be equal to a constant almost surely, as is necessarily the case when
$\lambda=\lambda_{0}$.

Next, consider the case $\lambda\neq\lambda_{0}$. Without loss of generality,
we could assume that $\lambda>\lambda_{0}$.

  1. 1.

Consider the case $\boldsymbol{\beta}\neq\boldsymbol{\beta}_{0}$ and let
introduce the event

$\mathcal{E}=\left\\{E\left[T(Y,\lambda)-T(Y,\lambda_{0})\mid\boldsymbol{X},\boldsymbol{Z}\right]-\boldsymbol{X}^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})\right.\\\
\left.=E[T(Y,\lambda)\mid\boldsymbol{Z}]-E[T(Y,\lambda_{0})|\boldsymbol{Z}]+\gamma
f^{-1}_{z}(\boldsymbol{Z})-E[\boldsymbol{X}\mid\boldsymbol{Z}]^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})\right\\}.$

Taking conditional expectation given $\boldsymbol{Z}$ on both sides, we deduce
that on $\mathcal{E}$ we necessarily have $\gamma=0$. Thus it suffices to
investigate the probability of the larger event

$\mathcal{E}^{\prime}=\left\\{E\left[T(Y,\lambda)-T(Y,\lambda_{0})\mid\boldsymbol{X},\boldsymbol{Z}\right]-\boldsymbol{X}^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})\right.\\\
\left.=E[T(Y,\lambda)\mid\boldsymbol{Z}]-E[T(Y,\lambda_{0})|\boldsymbol{Z}]-E[\boldsymbol{X}\mid\boldsymbol{Z}]^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})\right\\}.$

Note that the right hand side equality does not depend on $\boldsymbol{X}$. We
distinguish two sub-cases. First, the case where the components of
$(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})$ corresponding to
$\boldsymbol{X}_{c}$ are equal to zero. Thus the linear combination
$\boldsymbol{X}^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})$ does not
include any of the continuous components of $\boldsymbol{X}$. In this case,
for any value of $\boldsymbol{Z}$, the support of
$\boldsymbol{X}^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})$ is finite and
independent of the value of $\boldsymbol{X}_{c}$. Then Assumption 2.1.4
guarantees that the probability of the event $\mathcal{E}^{\prime}$ could not
be equal to 1. Next, consider the case where
$\boldsymbol{X}^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})$ includes
continuous components of $\boldsymbol{X}$. In this case, by Assumption 2.1.3,
the support of the variable
$\boldsymbol{X}^{T}(\boldsymbol{\beta}-\boldsymbol{\beta}_{0})$ is the whole
real line and, by the monotonicity of $\lambda\mapsto T(y;\lambda)$ for each
value $y>0$,
$E\left[T(Y,\lambda)-T(Y,\lambda_{0})|\boldsymbol{X},\boldsymbol{Z}\right]\geq
0$ almost surely, the statement follows again.

  2. 2.

Consider the case $\gamma\neq 0$ and
$\boldsymbol{\beta}=\boldsymbol{\beta}_{0}$. In this case

$f_{z}(\boldsymbol{Z})\left(E\left[T(Y,\lambda)-T(Y,\lambda_{0})\mid\boldsymbol{X},\boldsymbol{Z}\right]-E\left[T(Y,\lambda)-T(Y,\lambda_{0})\mid\boldsymbol{Z}\right]\right)=\gamma.$

Taking expectation on both sides, we deduce that necessarily $\gamma=0$.

  3. 3.

Consider the case
$(\gamma,\boldsymbol{\beta}^{T})^{T}=(0,\boldsymbol{\beta}_{0}^{T})^{T}$. Then
necessarily

$\displaystyle
E\bigg{[}T(Y,\lambda)-T(Y,\lambda_{0})|\boldsymbol{X},\boldsymbol{Z}\bigg{]}=E[T(Y,\lambda)|\boldsymbol{Z}]-E[T(Y,\lambda_{0})|\boldsymbol{Z}]\quad\text{almost
surely}.$

Once again the right hand side does not depend on $\boldsymbol{X}$, and thus
the probability of the event $\mathcal{E}$ could not be equal to 1 because of
Assumption 2.1.3.

Therefore, the first statement follows. Consider now the second statement of
the Lemma 2.1. First, note that the maps $\lambda\mapsto\gamma(\lambda)$,
$\lambda\mapsto\beta(\lambda)$ and

$\lambda\mapsto
E\left[g\left(\boldsymbol{U}_{1};\boldsymbol{\theta}(\lambda),\gamma(\lambda),\boldsymbol{\eta}_{\lambda,1}\right)g\left(\boldsymbol{U}_{2};(\lambda,\boldsymbol{\beta}(\lambda)^{T})^{T}\\!,\gamma(\lambda),\boldsymbol{\eta}_{\lambda,2}\right)\exp\left\\{-(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})^{T}\boldsymbol{D}(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})\right\\}\right],$
(A.1)

with $\lambda\in\Lambda:=[\lambda_{\rm min},\lambda_{\rm max}]$ and
$\boldsymbol{\theta}(\lambda)=(\lambda,\boldsymbol{\beta}(\lambda)^{T})^{T}$,
are continuous. Indeed, it is quite clear that
$\inf_{\boldsymbol{d}\in\mathcal{D}}E[\omega(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})]>0$.
Next, by Lebesgue’s Dominated Convergence Theorem we have

$\gamma(\lambda)=\frac{1}{E[\omega(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})]}E\left[\left(T(Y_{1},\lambda)-E[T(Y_{1},\lambda)\mid\boldsymbol{Z}_{1}]-(\boldsymbol{X}_{1}-E[\boldsymbol{X}_{1}\mid\boldsymbol{Z}_{1}])^{T}\boldsymbol{\beta}\right)f_{z}(\boldsymbol{Z}_{1})\omega(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})\right],$

and

$\boldsymbol{\beta}(\lambda)=E\left[(\boldsymbol{X}_{1}-E[\boldsymbol{X}_{1}\mid\boldsymbol{Z}_{1}])(\boldsymbol{X}_{2}-E[\boldsymbol{X}_{2}\mid\boldsymbol{Z}_{2}])^{T}f_{z}(\boldsymbol{Z}_{1})f_{z}(\boldsymbol{Z}_{2})\omega(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})\right]^{-1}\\\
E\left[(\boldsymbol{X}_{1}-E[\boldsymbol{X}_{1}\mid\boldsymbol{Z}_{1}])f_{z}(\boldsymbol{Z}_{1})\left(\left(T(Y_{2},\lambda)-E[T(Y_{2},\lambda)\mid\boldsymbol{Z}_{2}]\right)f_{z}(\boldsymbol{Z}_{2})-\gamma(\lambda)\right)\omega(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})\right],$

which are clearly continuous. Again, by dominated convergence argument, the
continuity of the map (A.1) follows. Finally, by the same inverse Fourier
Transform argument used by Lavergne and Patilea (2013) we get that

$E\left[g\left(\boldsymbol{U}_{1};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda,1}\right)g\left(\boldsymbol{U}_{2};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda,2}\right)\omega(\boldsymbol{W}_{1}-\boldsymbol{W}_{2})\right]=\frac{\pi^{-(p+q)/2}}{\sqrt{d_{1}\cdots
d_{p+q}}}\\\
\times\int_{\mathbb{R}^{p+q}}\left|E\left[E[g\left(\boldsymbol{U};\boldsymbol{\theta},\gamma,\boldsymbol{\eta}_{\lambda}\right)\mid\boldsymbol{X},\boldsymbol{Z}]\exp\left\\{2i\pi\boldsymbol{w}^{T}\left(\boldsymbol{X}^{T},\boldsymbol{Z}^{T}\right)^{T}\right\\}\right]\right|^{2}\exp\left\\{-\boldsymbol{w}^{T}\boldsymbol{D}^{-1}\boldsymbol{w}\right\\}d\boldsymbol{w}.$

For any $\boldsymbol{x},\boldsymbol{z}$, the map $\lambda\mapsto
E[g\left(\boldsymbol{U};\boldsymbol{\theta}(\lambda),\gamma(\lambda),\boldsymbol{\eta}_{\lambda}\right)\mid\boldsymbol{X}=\boldsymbol{x},\boldsymbol{Z}=\boldsymbol{z}]$
is continuous. By Lebesgue Dominated Convergence Theorem, the map

$\lambda\mapsto\int_{\mathbb{R}^{p+q}}\left|E\left[E[g\left(\boldsymbol{U};\boldsymbol{\theta}(\lambda),\gamma(\lambda),\boldsymbol{\eta}_{\lambda}\right)\mid\boldsymbol{X},\boldsymbol{Z}]\exp\left\\{2i\pi\boldsymbol{w}^{T}\left(\boldsymbol{X}^{T},\boldsymbol{Z}^{T}\right)^{T}\right\\}\right]\right|^{2}\exp\left\\{-\boldsymbol{w}^{T}{\rm
diag}(d_{L},\ldots,d_{L})^{-1}\boldsymbol{w}\right\\}d\boldsymbol{w},$

is continuous, and thus attains its minimum on the compact set
$\Lambda\setminus[\lambda_{0}-\varepsilon,\lambda_{0}+\varepsilon]$. The
minimum value is necessarily positive. Since $(d_{1}\cdots
d_{p+q})^{-1/2}\exp\left\\{-\boldsymbol{w}^{T}\boldsymbol{D}\boldsymbol{w}\right\\}\geq
d_{U}^{-(p+q)/2}\exp\left\\{-\boldsymbol{w}^{T}{\rm
diag}(d_{L},\ldots,d_{L})^{-1}\boldsymbol{w}\right\\}$, the last statement in
the Lemma follows.

  1. 1.

First, we note that

$\mathbb{P}\left(\left\\{\boldsymbol{u}^{T}\boldsymbol{X}+\boldsymbol{v}^{T}\boldsymbol{Z}=0,\;\;\forall(\boldsymbol{u}^{T},\boldsymbol{v}^{T})^{T}\text{
with norm equal to }1\right\\}\right)=0.$ (A.2)

This is a consequence of the fact that $\boldsymbol{W}$ is not degenerated.
Given a sample
$\left(\boldsymbol{X}^{T}_{1},\boldsymbol{Z}^{T}_{1}\right)^{T},\ldots,\left(\boldsymbol{X}^{T}_{n},\boldsymbol{Z}^{T}_{n}\right)^{T}$,
and a vector $\boldsymbol{a}=(a_{1},\ldots,a_{n})\in\mathbb{R}^{n}$, using the
inverse Fourier Transform, we could write

$\boldsymbol{a}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{a}=\frac{\pi^{-(p+q)/2}}{\sqrt{d_{1}\cdots
d_{p+q}}}\int_{\mathbb{R}^{p+q}}\left|\sum_{j=1}^{n}a_{j}\exp\left\\{2i\boldsymbol{w}^{T}\left(\boldsymbol{X}^{T}_{j},\boldsymbol{Z}^{T}_{j}\right)^{T}\right\\}\right|^{2}\exp\left\\{-\boldsymbol{w}^{T}\boldsymbol{D}^{-1}\boldsymbol{w}\right\\}d\boldsymbol{w},$

where $\boldsymbol{D}={\rm diag}(d_{1},\ldots,d_{p+q})$ with
$d_{1},\ldots,d_{p+q}\in[d_{L},d_{U}]$; see Assumption 2.1.3. Then,
necessarily

$\boldsymbol{a}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{a}=0\Longrightarrow\sum_{j=1}^{n}a_{j}\exp\left\\{2i\pi\boldsymbol{w}^{T}\left(\boldsymbol{X}^{T}_{j},\boldsymbol{Z}^{T}_{j}\right)^{T}\right\\}=0,\qquad\forall\boldsymbol{w}\in\mathbb{R}^{p+q}.$
(A.3)

Equation (A.2) indicates that, with probability 1, the unique vector
$\boldsymbol{a}$ which satisfies the right-hand side of (A.3) is
$\boldsymbol{a}=\boldsymbol{0}_{n}$ . This means that, with probability 1, the
matrix $\boldsymbol{\Omega}_{n}$ is positive definite.

Next, we use the following Cauchy-Schwarz555A similar so-called Cauchy-Schwarz
inequality was proposed by Lavergne (2008). To justify the statement, it
suffices to notice that
$\boldsymbol{B}^{T}\boldsymbol{B}-\boldsymbol{B}^{T}\boldsymbol{A}(\boldsymbol{A}^{T}\boldsymbol{A})^{-1}\boldsymbol{A}^{T}\boldsymbol{B}=\boldsymbol{\Gamma}^{T}\boldsymbol{\Gamma}$
with
$\boldsymbol{\Gamma}=\boldsymbol{B}-\boldsymbol{A}(\boldsymbol{A}^{T}\boldsymbol{A})^{-1}\boldsymbol{A}^{T}\boldsymbol{B}.$
inequality for matrices: let $\boldsymbol{A}\in\mathbb{R}^{n\times p_{1}}$
such that $\boldsymbol{A}^{T}\boldsymbol{A}$ is invertible and let
$\boldsymbol{B}\in\mathbb{R}^{n\times p_{2}}$. Then

$\boldsymbol{B}^{T}\boldsymbol{B}-\boldsymbol{B}^{T}\boldsymbol{A}(\boldsymbol{A}^{T}\boldsymbol{A})^{-1}\boldsymbol{A}^{T}\boldsymbol{B}\;\text{
is positive semi-definite. }$

Moreover, the equality
$\boldsymbol{B}^{T}\boldsymbol{B}=\boldsymbol{B}^{T}\boldsymbol{A}(\boldsymbol{A}^{T}\boldsymbol{A})^{-1}\boldsymbol{A}^{T}\boldsymbol{B}$
is equivalent to the relationship
$\boldsymbol{B}=\boldsymbol{A}(\boldsymbol{A}^{T}\boldsymbol{A})^{-1}\boldsymbol{A}^{T}\boldsymbol{B}.$
For any non null vector $\boldsymbol{u}\in\mathbb{R}^{p}$, taking

$\boldsymbol{B}=\boldsymbol{\Omega}_{n}^{1/2}\quad\text{ and
}\boldsymbol{A}=\boldsymbol{\Omega}_{n}^{1/2}\boldsymbol{1}_{n},$

we deduce that $\mathbb{D}_{n}$ is positive semi-definite and thus
$\boldsymbol{u}^{T}\widehat{\mathbb{X}}_{n}^{T}\mathbb{D}_{n}\widehat{\mathbb{X}}_{n}\boldsymbol{u}\geq
0$. (Herein, $\boldsymbol{\Omega}_{n}^{1/2}$ is the positive definite square
root of $\boldsymbol{\Omega}_{n}$.) Meanwhile, by elementary matrix algebra,
we deduce that, for any $\boldsymbol{a}\in\mathbb{R}^{n}$,

$\boldsymbol{a}^{T}\mathbb{D}_{n}\boldsymbol{a}=0\quad\Leftrightarrow\quad\left(\boldsymbol{a}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{a}^{T}\right)\left(\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{1}_{n}\right)=\left(\boldsymbol{a}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{1}_{n}\right)^{2}.$

Then, the Cauchy-Schwarz inequality indicates that
$\boldsymbol{a}^{T}\mathbb{D}_{n}\boldsymbol{a}=0$ if and only if
$\boldsymbol{a}=a\boldsymbol{1}_{n}$ for some scalar $a\neq 0$. Thus,
$\boldsymbol{u}^{T}\widehat{\mathbb{X}}_{n}^{T}\mathbb{D}_{n}\widehat{\mathbb{X}}_{n}\boldsymbol{u}=0$
if and only if $\widehat{\mathbb{X}}_{n}\boldsymbol{u}=a\boldsymbol{1}_{n}$
for some $a\neq 0$. By construction,
$\boldsymbol{1}_{n}^{T}\widehat{\mathbb{X}}_{n}=\boldsymbol{0}_{n}$, and thus
necessarily $a=0$, which leads to a contradiction. Thus,
$\widehat{\mathbb{X}}_{n}^{T}\mathbb{D}_{n}\widehat{\mathbb{X}}_{n}$ is almost
surely invertible. Note that we could also write

$\mathbb{D}_{n}=\left[\boldsymbol{I}_{n\times
n}-\frac{1}{\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{1}_{n}}\boldsymbol{1}_{n}\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\right]^{T}\boldsymbol{\Omega}_{n}\left[\boldsymbol{I}_{n\times
n}-\frac{1}{\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{1}_{n}}\boldsymbol{1}_{n}\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\right]$

and deduce the positive semi-definiteness of $\mathbb{D}_{n}$ from the
positive definiteness of $\boldsymbol{\Omega}_{n}$.

  2. 2.

We could rewrite $\widehat{\mathbb{B}}_{n}$ under the form

$\widehat{\mathbb{B}}_{n}=\left[\boldsymbol{I}_{n\times
n}-\widehat{\mathbb{X}}_{n}\left(\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\right)^{-1}\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\right]^{T}\mathbb{D}_{n}\left[\boldsymbol{I}_{n\times
n}-\widehat{\mathbb{X}}_{n}\left(\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\widehat{\mathbb{X}}_{n}\right)^{-1}\widehat{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}\right]$

and deduce its positive semi-definiteness from the positive definiteness of
$\mathbb{D}_{n}$.

Recall that $\Lambda:=[\lambda_{\rm min},\lambda_{\rm max}]$ and let

$\widehat{M}_{n}(\lambda)=n^{-2}\widehat{\mathbb{Y}}_{n}(\lambda)^{T}\;\widehat{\mathbb{B}}_{n}\;\widehat{\mathbb{Y}}_{n}(\lambda)$

such that

$s^{-2\lambda}\widehat{M}_{n}(\lambda)=n^{-2}s^{-\lambda}\widehat{\mathbb{Y}}_{n}(\lambda)^{T}\;\widehat{\mathbb{B}}_{n}\;s^{-\lambda}\widehat{\mathbb{Y}}_{n}(\lambda),$

and, thus,
$\widehat{\lambda}=\arg\min\limits_{\lambda\in\Lambda}s^{-2\lambda}\widehat{M}_{n}(\lambda)$.
Next, let

$M_{n}(\lambda)=n^{-2}{\mathbb{Y}}_{n}(\lambda)^{T}\;{\mathbb{B}}_{n}\;{\mathbb{Y}}_{n}(\lambda).$

By construction,

$M_{n}(\lambda)=Q_{n}((\lambda,{\boldsymbol{\beta}}_{n}(\lambda)^{T})^{T},{\gamma}_{n}(\lambda)),$

where

$Q_{n}((\lambda,{\boldsymbol{\beta}}^{T})^{T},{\gamma})=n^{-2}\left({\mathbb{Y}}_{n}(\lambda)-\gamma\boldsymbol{1}_{n}-{\mathbb{X}}_{n}\boldsymbol{\beta}\right)^{T}\boldsymbol{\Omega}_{n}\left({\mathbb{Y}}_{n}(\lambda)-\gamma\boldsymbol{1}_{n}-{\mathbb{X}}_{n}\boldsymbol{\beta}\right)$

and

$\gamma_{n}(\lambda)=\frac{1}{\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\boldsymbol{1}_{n}}\boldsymbol{1}_{n}^{T}\boldsymbol{\Omega}_{n}\left({\mathbb{Y}}_{n}(\lambda)-{\mathbb{X}}_{n}{\boldsymbol{\beta}}(\lambda)\right)\qquad\text{
and
}\qquad{\boldsymbol{\beta}}_{n}(\lambda)=\left({\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{X}}_{n}\right)^{-1}{\mathbb{X}}_{n}^{T}{\mathbb{D}}_{n}{\mathbb{Y}}_{n}(\lambda).$

Let

$Q((\lambda,{\boldsymbol{\beta}}^{T})^{T},{\gamma})=E\left[Q_{n}((\lambda,{\boldsymbol{\beta}}^{T})^{T},{\gamma})\right],\qquad\lambda\in\Lambda,\boldsymbol{\beta}\in\mathbb{R}^{p},\gamma\in\mathbb{R}.$

Next, let $c>0$ be a lower bound of the support of $Y$. Then, necessarily
$c<s_{0}$ and we could work on the event $c\leq\inf S_{n}$, that is $s$ stays
away from zero. In order to prove the uniform consistency it will suffice to
prove

$\sup_{h\in\mathcal{H}^{c}_{n}}\sup_{\boldsymbol{d}\in\mathcal{D}}\sup_{\lambda\in\Lambda}\left|\widehat{M}_{n}(\lambda)-M_{n}(\lambda)\right|=o_{\mathbb{P}}(1),$
(A.4)
$\sup_{h\in\mathcal{H}^{c}_{n}}\sup_{\boldsymbol{d}\in\mathcal{D}}\sup_{\lambda\in\Lambda}\left|Q_{n}((\lambda,{\boldsymbol{\beta}}_{n}(\lambda)^{T})^{T},{\gamma}_{n}(\lambda))-Q_{n}((\lambda,\widehat{\boldsymbol{\beta}}(\lambda)^{T})^{T},\widehat{\gamma}(\lambda))\right|=o_{\mathbb{P}}(1),$
(A.5)