# Semiparametric Approach to Estimation of Marginal and Quantile Effects

Seong-ho Lee1, Yanyuan Ma1, and Elvezio Ronchetti2
1Department of Statistics, Pennsylvania State University, USA
2Research Center for Statistics and GSEM, University of Geneva, Switzerland

###### Abstract

We consider a semiparametric generalized linear model and study estimation of
both marginal and quantile effects in this model. We propose an approximate
maximum likelihood estimator, and rigorously establish the consistency, the
asymptotic normality, and the semiparametric efficiency of our method in both
the marginal effect and the quantile effect estimation. Simulation studies are
conducted to illustrate the finite sample performance, and we apply the new
tool to analyze a Swiss non-labor income data and discover a new interesting
predictor.

Keywords: generalized linear model, marginal effect, quantile effect,
misspecification, robustness, semiparametric efficiency.

##  1 Introduction

Generalized linear models (GLMs) (McCullagh & Nelder, 1989) are arguably the
most frequently used models in statistical research and application. With
response variable $Y\in\mathbb{R}$ and covariates ${\bf X}\in\mathbb{R}^{p}$,
GLMs have the familiar form

$\displaystyle f_{Y\mid{\bf X}}(y,{\bf
x},{\bm{\beta}},\phi)=\exp\left\\{\frac{y{\bm{\beta}}^{\rm T}{\bf
x}-b({\bm{\beta}}^{\rm T}{\bf x})}{a(\phi)}+c(y,\phi)\right\\},$ (1)

where ${\bm{\beta}}\in\mathbb{R}^{p}$ is an unknown regression coefficient
vector, and $\phi$ is typically an unknown scalar parameter, if it is present.
Here, the functions $a(\cdot),b(\cdot)$, and $c(\cdot)$ have pre-specified
forms, and various choices of these functions lead to different GLMs
extensively studied in the statistical literature.

While a large body of the classic literature is mainly concerned with the
estimation and the inference of the model parameters ${\bm{\beta}}$ and $\phi$
in (1), an equally important target of research is to estimate marginal
effect, say ${\bm{\xi}}$, defined as

$\displaystyle{\bm{\xi}}\equiv E\left\\{\frac{\partial E(Y\mid{\bf
X})}{\partial{\bf X}}\right\\}.$ (2)

The marginal effect is a quantity of important meaning in several fields, in
particular econometrics, and is used to capture the average rate of change of
the regression mean with respect to the covariates (Greene, 2000).

A complementary, less known but important quantity is quantile effect, defined
as

$\displaystyle{\bm{\eta}}_{\tau}\equiv E\left\\{\frac{\partial
Q_{\tau}(Y\mid{\bf X})}{\partial{\bf X}}\right\\}$ (3)

for a continuous response variable $Y$, where $\tau\in(0,1)$ is the quantile
level and $Q_{\tau}(Y\mid{\bf x})$ is the $\tau$th conditional quantile of $Y$
given ${\bf x}$. Here, we restrict our attention to the continuous response
case only because the definition of conditional quantile in the discrete case
is not universally agreed upon in the literature (Parzen, 2004). The quantile
effect measures the average rate of change in the conditional quantiles, and
can be equivalently written as
$E\\{{\bm{\beta}}Q_{\tau}^{\prime}(Y\mid{\bm{\beta}}^{\rm T}{\bf X})\\}$,
where $Q_{\tau}^{\prime}(Y\mid{\bm{\beta}}^{\rm T}{\bf x})$ is the derivative
of $Q_{\tau}(Y\mid{\bm{\beta}}^{\rm T}{\bf x})$ with respect to
${\bm{\beta}}^{\rm T}{\bf x}$. Obviously, our consideration of quantile effect
is a natural consequence of considering quantile regression, which is a key
statistical tool widely used as an alternative to the traditional mean
regression models. Instead of modeling the conditional mean, the quantile
regression estimates the relationship between the covariates and the response
through the conditional quantile $Q_{\tau}(Y\mid{\bm{\beta}}^{\rm T}{\bf x})$,
and considering its corresponding rate of change directly leads to the
quantile effect defined in (3). Comparing (2) and (3), it is clear that
${\bm{\xi}}$ and ${\bm{\eta}}_{\tau}$ are in fact marginal mean effect and
marginal quantile effect respectively, while we shorten the names to marginal
effect and quantile effect.

Marginal effect has direct interpretation as measures of risk. For instance,
in logistic regression, especially in the medical literature, it is customary
to misinterpret odds ratios as measures of risk. Instead, marginal and
quantile effects represent changes in the probability of the occurrence of a
binary event with respect to given risk factors and have direct interpretation
on the probability scale. This point has been stressed in the medical
literature, e.g., by Norton et al. (2019). Marginal effect also captures how
strongly the mean treatment outcome relies on different covariates on average.
This further leads to the understanding on how average treatment effect is
affected by covariates. Treatment effect is directly linked to causal
inference which is currently one of the trending topics in Statistics, as
further confirmed by the 2021 Nobel Prize in Economics based on the seminal
paper by Imbens & Angrist (1994). In fact, estimation of marginal effect has
been studied by many earlier works in the statistics literature, where the
name “average mean derivatives estimation” was adopted. For example, Härdle &
Stoker (1989) systematically studied the estimation problem, and Newey &
Stoker (1993) studied semiparametric efficiency properties of weighted average
mean derivative estimators. Likewise, quantile effect, formerly termed as
average quantile derivatives, has also been proposed and studied as an
alternative to the marginal effect by Chaudhuri et al. (1997). Following these
earlier works, different approaches to estimating marginal effect under
various settings have been proposed, see, e.g., Hristache et al. (2001),
Cattaneo et al. (2010), and Cattaneo et al. (2013). The estimation of marginal
effect is also quite widespread in medical economics concerning health
outcomes. For instance, a two-year study performed at the University of
Chicago focuses on marginal effect on expenditures and length of stays in
hospitals, rather than on parameter estimation of GLM; see Basu & Rathouz
(2005). For the same study, Manning et al. (2005) estimate marginal effect on
inpatient expenditures using a parametric family of distributions, a three
parameter generalized Gamma distribution, which is an extension of the GLMs.
Marginal effect is also frequently the target of study in economics.
Parametric estimation of marginal effect in microeconometrics concerning the
female labor force participation for 872 women from Switzerland can be found
in Gerfin (1996) and Kleiber & Zeileis (2008). We will use this example as a
benchmark to compare our nonparametric estimation in Section 5. In a GLM, the
marginal effect happens to be identical to the regression coefficient
${\bm{\beta}}$ multiplied by average conditional variance, i.e.
${\bm{\xi}}={\bm{\beta}}E\left\\{\hbox{var}(Y\mid{\bm{\beta}}^{\rm T}{\bf
X})\right\\}$. Similarly for the quantile effect, we have
${\bm{\eta}}_{\tau}={\bm{\beta}}E\left\\{Q_{\tau}^{\prime}(Y\mid{\bm{\beta}}^{\rm
T}{\bf X})\right\\}$. This allows for their straightforward estimation.
Standard statistical packages provide direct estimation of marginal effect
under GLMs, including the margins command in Stata and the package mfx in R
(Fernihough, 2019).

However, despite its central role in the statistical literature and its wide
application, GLM has its restrictions. Indeed, GLM prespecifies the three
functions $a(\cdot),b(\cdot)$, and $c(\cdot)$ hence is purely parametric,
therefore the model is susceptible to model misspecification. In addition, to
facilitate computation, the choices of the three functions are often out of
convenience, which further increases the chance of mis-modeling. To take into
account the possibly large impact of deviation from the distributional
assumptions on the parameter estimates and on the corresponding inference
based on GLMs, a broad stream of the literature in the past decades has
focused on robust methods for the estimation of the parameters; see e.g.
Stefanski et al. (1986), Künsch et al. (1989), and Cantoni & Ronchetti (2001)
for a general class of optimally bounded influence estimators and tests for
the parameters of GLMs. The corresponding formulas for the estimation of
marginal and quantile effects remain the same as the parametric case, but with
the robust estimates of the parameters.

Here we go one step further in relaxing the parametric form of GLMs and
consider a semiparametric generalized linear model (sGLM), where we assume

$\displaystyle f_{Y\mid{\bf X}}(y,{\bf
x},{\bm{\beta}},c)=\frac{\exp\left\\{y{\bm{\beta}}^{\rm T}{\bf
x}+c(y)\right\\}}{\int\exp\left\\{y{\bm{\beta}}^{\rm T}{\bf
x}+c(y)\right\\}d\mu(y)},$ (4)

and leave $c(\cdot)$ unspecified. Here and in the following text, $\mu(\cdot)$
denotes the Lebesgue measure for a continuous variable $y$ and the counting
measure for a discrete $y$. For identifiability, we fix $c(0)=0$ and for
convenience, we assume each component of ${\bf x}$ to be centered so $E({\bf
X})={\bf 0}$. Our focus is to estimate marginal and quantile effects of (4),
i.e. ${\bm{\xi}}$ defined in (2) and ${\bm{\eta}}_{\tau}$ in (3).

Model (4) is not entirely new, and has been proposed and studied in Luo & Tsai
(2012) and Huang & Rathouz (2012) as a semiparametric proportional likelihood
ratio model. Ning et al. (2017) studied high dimensional issue in (4), while
Lin et al. (2021) studied the estimation of the parameter ${\bm{\beta}}$ in
(4). Different from these existing works, our main focus is not in the model
parameter ${\bm{\beta}}$, but in the marginal effect ${\bm{\xi}}$ and the
quantile effect ${\bm{\eta}}_{\tau}$.

In this context, we develop an approximate likelihood procedure to estimate
both the model parameters ${\bm{\beta}}$ and $c(\cdot)$ in (4), which allows
us to subsequently estimate the marginal effect ${\bm{\xi}}$ and the quantile
effect ${\bm{\eta}}_{\tau}$. We show that despite of the apparent difference
of the estimation variability results from the theoretical efficiency bounds,
our estimators are in fact semiparametricly efficient for both ${\bm{\xi}}$
and ${\bm{\eta}}_{\tau}$. We point out that our method differs from Lin et al.
(2021) primarily in how to estimate the nonparametric component $c(y)$, which
in turn has significant influence on estimation of the model parameters and
functionals. Due to the special structure of model (4), kernel estimation for
$c(y)$ is difficult to implement without giving up the model information and
resorting to model-free estimation such as the Nadaraya-Watson estimator,
hence the estimation of $c(y)$ in Lin et al. (2021) is not ideal. This leads
to efficiency loss of the ${\bm{\beta}}$ estimation in finite samples, and
greatly affects the efficiency of the ${\bm{\xi}}$ and ${\bm{\eta}}_{\tau}$
estimation. The efficiency gain from our method is demonstrated by simulation
experiments in Section 4.

Specifically, our main contributions to the literature are the following.
First, we define and consider estimating marginal effect and quantile effect
in a semiparametric model, which have been studied by many earlier works in
Statistics under various settings and are receiving more attention due to the
recent interest in mean treatment effect or mean outcome estimation. Second,
we discover that in this specific problem, kernel estimation may not be a good
approach to estimate the function $c(y)$ while B-spline based estimation fits
naturally and has a clear advantage. This is not usually the case in many
semiparametric estimation problems. Third, we provide a clear and transparent
estimation procedure, which is computationally simple and enjoys the advantage
of convex optimization. Fourth, we establish the optimality of the resulting
marginal and quantile effect estimators as functionals of the model
parameters. This optimality is not a straightforward result to obtain. It
requires keen statistical insight to link two seemingly unrelevant results and
the right mathematical tools to achieve it.

The rest of the article is organized as the following. In Section 2, we
provide the main methodological results of our work and describe the
estimation procedures for both ${\bm{\xi}}$ and ${\bm{\eta}}_{\tau}$. In
Section 3, we establish the asymptotic properties of our estimators, compare
them to the theoretical efficiency bounds, and show that they in fact reach
the efficiency bounds. Simulation studies are conducted in Section 4 to
illustrate the finite sample performance of our method in comparison with
existing methods. An interesting dataset concerning Swiss non-labor income is
analyzed in Section 5, where our nonparametric analysis reveals new insights
concerning model choice, marginal effect, and significant covariates that
would remain hidden if using a strict parametric approach. We conclude our
paper with a discussion in Section 6. All the technical details are relegated
to the supplementary material.

##  2 Methodology

###  2.1 Efficiency Bound of Marginal Effect Estimation

Before engaging ourselves in the estimation and inference of the marginal
effect ${\bm{\xi}}$ associated with the model (4), we would like first to
understand the limit of our endeavor by establishing the efficiency bound of
the ${\bm{\xi}}$ estimation. Let $v({\bm{\beta}}^{\rm T}{\bf x})\equiv
E[\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf x})\\}^{2}\mid{\bm{\beta}}^{\rm T}{\bf
x}]$. The derivation detailed in Section S.1 leads to the efficient influence
function as

$\displaystyle\bm{\phi}_{\rm eff}(y,{\bf x})={\bm{\beta}}v({\bm{\beta}}^{\rm
T}{\bf x})-{\bm{\beta}}E\\{v({\bm{\beta}}^{\rm T}{\bf
X})\\}+{\bm{\beta}}y^{2}+{\bf M}{\bf x}y+{\bf a}(y)-E\\{{\bm{\beta}}Y^{2}+{\bf
M}{\bf x}Y+{\bf a}(Y)\mid{\bf x}\\},$

where ${\bf a}(y)$ satisfies

$\displaystyle E[E\\{{\bf a}(Y)\mid{\bf X}\\}\mid y]-{\bf a}(y)$
$\displaystyle=$ $\displaystyle 2{\bm{\beta}}E[yE(Y\mid{\bf
X})-E\\{YE(Y\mid{\bf X})\mid{\bf X}\\}\mid y]+{\bf M}E[{\bf X}\\{y-E(Y\mid{\bf
X})\\}\mid y].$

and

$\displaystyle{\bf M}$ $\displaystyle=$
$\displaystyle\left(E\\{v({\bm{\beta}}^{\rm T}{\bf X})\\}{\bf
I}-E\left[2{\bm{\beta}}{\bf X}^{\rm T}Yv({\bm{\beta}}^{\rm T}{\bf X})+{\bf
a}(Y)\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf X})\\}{\bf X}^{\rm T}\right]\right)$
$\displaystyle\times\left[E\\{{\bf X}{\bf X}^{\rm T}v({\bm{\beta}}^{\rm T}{\bf
X})\\}\right]^{-1},$

Obviously, the variance of the efficient influence function, i.e.
$\hbox{var}\\{\bm{\phi}_{\rm eff}(Y,{\bf X})\\}$, is the efficiency bound in
estimating ${\bm{\xi}}$.

###  2.2 Efficiency Bound of Quantile Effect Estimation

In Section S.2, we further derive the efficient influence function of
${\bm{\eta}}_{\tau}$ under (4). Now for notational brevity, let $q(\nu)\equiv
Q_{\tau}(Y\mid\nu)$, $q^{\prime}(\nu)\equiv Q^{\prime}_{\tau}(Y\mid\nu)$,
$q^{\prime\prime}(\nu)\equiv Q^{\prime\prime}_{\tau}(Y\mid\nu)$,
$\epsilon(\nu)\equiv\tau-I\\{Y<q(\nu)\\}$, and
$\epsilon^{\prime}(\nu)\equiv-\delta\\{q(\nu)-Y\\}q^{\prime}(\nu)$. The
efficient influence function turns out to be

$\displaystyle\bm{\phi}_{\rm eff}(y,{\bf
x})={\bm{\beta}}q^{\prime}({\bm{\beta}}^{\rm T}{\bf
x})-{\bm{\beta}}E\\{q^{\prime}({\bm{\beta}}^{\rm T}{\bf X})\\}+{\bf M}_{1}{\bf
x}y+{\bf a}(y)-E\\{{\bf M}_{1}{\bf x}Y+{\bf a}(Y)|{\bf x}\\},$

where ${\bf a}(y)$ satisfies

$\displaystyle E[E\\{{\bf a}(Y)\mid{\bf X}\\}\mid y]-{\bf
a}(y)=-{\bm{\beta}}E\\{r(y,{\bm{\beta}}^{\rm T}{\bf X})|y\\}+{\bf M}_{1}E[{\bf
X}\\{y-E(Y|{\bf X})\\}|y],$

in which

$\displaystyle
r(Y,\nu)\equiv\frac{\epsilon(\nu)Y+\epsilon^{\prime}(\nu)-\epsilon(\nu)[q(\nu)+q^{\prime}(\nu)\nu+q^{\prime}(\nu)c^{\prime}\\{q(\nu)\\}]}{f\\{q(\nu),\nu\\}},$

and

$\displaystyle{\bf M}_{1}$ $\displaystyle=$
$\displaystyle(E\\{q^{\prime}({\bm{\beta}}^{\rm T}{\bf X})\\}{\bf
I}+{\bm{\beta}}E\\{{\bf X}^{\rm T}q^{\prime\prime}({\bm{\beta}}^{\rm T}{\bf
X})\\}-E[{\bf a}(Y){\bf X}^{\rm T}\\{Y-E(Y|{\bf X})\\}])$
$\displaystyle\times[E\\{{\bf X}{\bf X}^{\rm T}v({\bm{\beta}}^{\rm T}{\bf
X})\\}]^{-1}.$

Similar to the case for estimating ${\bm{\xi}}$, the variance of the efficient
influence function is the efficiency bound for estimating
${\bm{\eta}}_{\tau}$.

###  2.3 Estimation Procedure

We now propose an estimation procedure under the model (4). We first consider
the case where the response $Y$ is distributed continuously. In such case, we
approximate $c(\cdot)$ of the conditional density $f_{Y\mid{\bf X}}(y,{\bf
x},{\bm{\beta}},c)$ by a B-spline curve ${\bf B}(\cdot)^{\rm T}{\bm{\gamma}}$,
where ${\bf B}(\cdot)\equiv\\{B_{1}(\cdot),\dots,B_{m}(\cdot)\\}^{\rm T}$ is a
B-spline basis vector and
${\bm{\gamma}}\equiv(\gamma_{1},\dots,\gamma_{m})^{\rm T}$ is an unknown
coefficient of the bases. Since we assume $c(0)=0$ for identifiability, we
ignore the first B-spline basis $B_{1}(\cdot)$ which corresponds to the
intercept of the curve, by letting $\gamma_{1}=0$. Therefore, we replace
$c(y)$ in (4) by ${\bf B}(y)^{\rm T}{\bm{\gamma}}$ with $\gamma_{1}=0$, and
estimate ${\bm{\beta}}$ and ${\bm{\gamma}}$ through maximizing the approximate
loglikelihood

$\displaystyle
l({\bm{\beta}},{\bm{\gamma}})\equiv\sum_{i=1}^{n}\left[y_{i}{\bm{\beta}}^{\rm
T}{\bf x}_{i}+{\bf B}(y_{i})^{\rm
T}{\bm{\gamma}}-\hbox{log}\int\exp\left\\{y{\bm{\beta}}^{\rm T}{\bf
x}_{i}+{\bf B}(y)^{\rm T}{\bm{\gamma}}\right\\}d\mu(y)\right].$ (5)

Note that the loglikelihood $l({\bm{\beta}},{\bm{\gamma}})$ is a concave
function, hence the optimizer is unique and can be readily obtained using off-
the-shelf convex optimizers. In our implementation, we used the built-in R
function optim for optimization, yet found the computation fast with
satisfactory performance. For example, in our real data analysis where the
dataset contains $871$ observations and $7$ predictors, it takes less than 10
seconds to obtain the optimizer within the relative tolerance of $10^{-6}$.

Once we obtain $\widehat{\bm{\beta}}$ and $\widehat{c}(\cdot)\equiv{\bf
B}(\cdot)^{\rm T}\widehat{\bm{\gamma}}$, we can easily estimate the marginal
effect ${\bm{\xi}}$ using (2) through

$\displaystyle\widehat{\bm{\xi}}$ $\displaystyle\equiv$
$\displaystyle\widehat{\bm{\beta}}\widehat{E}\left\\{\widehat{\hbox{var}}(Y\mid{\bm{\beta}}^{\rm
T}{\bf X})\right\\}$ $\displaystyle=$
$\displaystyle\widehat{\bm{\beta}}\,n^{-1}\sum_{i=1}^{n}\left[\int
y^{2}f_{Y\mid{\bf X}}(y,{\bf
x}_{i},\widehat{\bm{\beta}},\widehat{c})d\mu(y)-\left\\{\int yf_{Y\mid{\bf
X}}(y,{\bf
x}_{i},\widehat{\bm{\beta}},\widehat{c})d\mu(y)\right\\}^{2}\right].$

Similarly, we estimate the quantile effect ${\bm{\eta}}_{\tau}$ using (3)
through

$\displaystyle\widehat{\bm{\eta}}_{\tau}$ $\displaystyle\equiv$
$\displaystyle\widehat{\bm{\beta}}\widehat{E}\left\\{\widehat{Q}^{\prime}_{\tau}(Y\mid{\bm{\beta}}^{\rm
T}{\bf X})\right\\}$ $\displaystyle=$
$\displaystyle\widehat{\bm{\beta}}\,n^{-1}\sum_{i=1}^{n}\frac{\tau\int_{0}^{1}yf_{Y\mid{\bf
X}}(y,{\bf
x}_{i},\widehat{\bm{\beta}},\widehat{c})d\mu(y)-\int_{0}^{q_{i}}yf_{Y\mid{\bf
X}}(y,{\bf x}_{i},\widehat{\bm{\beta}},\widehat{c})d\mu(y)}{f_{Y\mid{\bf
X}}(q_{i},{\bf
x}_{i},\widehat{\bm{\beta}},\widehat{c})}\bigg{|}_{q_{i}=\widehat{Q}_{\tau}(Y\mid{\bm{\beta}}^{\rm
T}{\bf x}_{i})},$

where $\widehat{Q}_{\tau}(Y\mid{\bm{\beta}}^{\rm T}{\bf x}_{i}),i=1,\dots,n$,
the estimated conditional $\tau$th quantiles of $Y$ given ${\bf x}_{i}$, are
obtained by solving for $q_{i}$ from

$\displaystyle\int_{0}^{q_{i}}f_{Y\mid{\bf X}}(y,{\bf
x}_{i},\widehat{\bm{\beta}},\widehat{c})d\mu(y)=\tau.$

When the response $Y$ is a categorical variable taking values $\\{0,...,m\\}$,
then the goal is to estimate ${\bm{\beta}}$ and $c(0),\dots,c(m)$. In this
case, (4) is a purely parametric model and we can proceed with maximum
likelihood estimation (MLE). An alternative way of viewing this is that we
replace $c(\cdot)$ by ${\bf B}(\cdot)^{\rm T}{\bm{\gamma}}$ where ${\bf
B}(\cdot)\equiv\\{I(\cdot=0),\dots,I(\cdot=m)\\}^{\rm T}$ with $\gamma_{1}=0$,
and maximize the approximate loglikelihood given in (5) with respect to
${\bm{\beta}}$ and ${\bm{\gamma}}$. Having obtained $\widehat{\bm{\beta}}$ and
$\widehat{c}(\cdot)$, we use

$\displaystyle\hbox{pr}(Y=y|{\bf
x};{\bm{\beta}},c)=\frac{\exp\\{y{\bm{\beta}}^{\rm T}{\bf
x}+c(y)\\}}{\sum_{y=0}^{m}\exp\\{y{\bm{\beta}}^{\rm T}{\bf x}+c(y)\\}},$

$y=0,\dots,m$ to estimate the marginal effect ${\bm{\xi}}$ by

$\displaystyle\widehat{\bm{\xi}}\equiv\widehat{\bm{\beta}}\,n^{-1}\sum_{i=1}^{n}\left[\sum_{y=0}^{m}y^{2}\,\hbox{pr}(Y=y|{\bf
x}_{i};\widehat{\bm{\beta}},\widehat{c})-\left\\{\sum_{y=0}^{m}y\,\hbox{pr}(Y=y|{\bf
x}_{i};\widehat{\bm{\beta}},\widehat{c})\right\\}^{2}\right].$

Because there is no generally accepted unique way of defining quantiles for
discrete data, we do not further study the quantile effect estimation in this
case.

##  3 Theoretical Properties

We now establish the theoretical properties of our proposed estimators for
both the marginal effect and the quantile effect in (4).

###  3.1 Continuous Response

First we analyze the properties of our estimators under the continuous
response case. We simplify $d\mu(y)$ as $dy$ below. To set notation, let

$\displaystyle f_{Y\mid{\bf X}}^{*}(y,{\bf x},{\bm{\beta}},{\bm{\gamma}})$
$\displaystyle\equiv$ $\displaystyle\frac{\exp\\{y{\bm{\beta}}^{\rm T}{\bf
x}+{\bf B}(y)^{\rm T}{\bm{\gamma}}\\}}{\int\exp\\{y{\bm{\beta}}^{\rm T}{\bf
x}+{\bf B}(y)^{\rm T}{\bm{\gamma}}\\}dy},$ $\displaystyle E^{*}\left\\{{\bf
g}(Y)|{\bf x},{\bm{\beta}},{\bm{\gamma}}\right\\}$ $\displaystyle\equiv$
$\displaystyle\int{\bf g}(y)f_{Y\mid{\bf X}}^{*}(y,{\bf
x},{\bm{\beta}},{\bm{\gamma}})dy,$ $\displaystyle\hbox{var}^{*}\left(Y|{\bf
x},{\bm{\beta}},{\bm{\gamma}}\right)$ $\displaystyle\equiv$ $\displaystyle
E^{*}(Y^{2}|{\bf x},{\bm{\beta}},{\bm{\gamma}})-\left\\{E^{*}(Y|{\bf
x},{\bm{\beta}},{\bm{\gamma}})\right\\}^{2},$
$\displaystyle\hbox{cov}^{*}\left\\{{\bf g}(Y),{\bf h}(Y)|{\bf
x},{\bm{\beta}},{\bm{\gamma}}\right\\}$ $\displaystyle\equiv$ $\displaystyle
E^{*}\left\\{{\bf g}(Y){\bf h}(Y)^{\rm T}|{\bf
x},{\bm{\beta}},{\bm{\gamma}}\right\\}$ $\displaystyle-E^{*}\left\\{{\bf
g}(Y)|{\bf x},{\bm{\beta}},{\bm{\gamma}}\right\\}E^{*}\left\\{{\bf h}(Y)^{\rm
T}|{\bf x},{\bm{\beta}},{\bm{\gamma}}\right\\},$

where we use $E^{*}$ instead of $E$ to emphasize that the expectation is
computed under the approximate model, which has the same form as (4), but with
$c(\cdot)$ replaced by ${\bf B}^{\rm T}(\cdot){\bm{\gamma}}$.

For real numbers $a_{n}$ and $b_{n}$, $a_{n}\asymp b_{n}$ denotes
$a_{n}=O(b_{n})$ and $b_{n}=O(a_{n})$ simultaneously. Similarly, for random
variables $A_{n}$ and $B_{n}$, $A_{n}\asymp_{p}B_{n}$ denotes
$A_{n}=O_{p}(B_{n})$ and $B_{n}=O_{p}(A_{n})$ simultaneously. For a vector
${\bf a}=(a_{1},\dots,a_{d})^{\rm T}\in\mathbb{R}^{d}$, we denote the
$l_{p}$-norm of ${\bf a}$ as $\|{\bf
a}\|_{p}\equiv(|a_{1}|^{p}+\dots+|a_{d}|^{p})^{1/p}$, $1\leq p\leq\infty$. For
a matrix ${\bf A}\in\mathbb{R}^{r\times c}$, we denote the induced
$l_{p}$-norm of ${\bf A}$ as $\|{\bf A}\|_{p}\equiv\sup_{{\bf
u}\in\mathbb{R}^{c},\|{\bf u}\|_{p}=1}\|{\bf A}{\bf u}\|_{p}$, $1\leq
p\leq\infty$. For a function $g(\cdot)$ in the $L^{2}$ space, we denote its
$L_{p}$-norm as $\lVert
g(\cdot)\rVert_{p}\equiv\\{\int_{0}^{1}|g(y)|^{p}dy\\}^{1/p}$. We denote the
set of the $q$th order smooth functions as $C^{q}([0,1])\equiv\\{g:g^{(q)}\in
C([0,1])\\}$.

To facilitate the theoretical derivation, we view the estimation procedure
described in Section 2.3 alternatively as a profile procedure. Specifically,
treating ${\bm{\beta}}$ as a fixed parameter, estimate $c(\cdot)$ by the
spline curve $\widehat{c}(\cdot,{\bm{\beta}})\equiv{\bf B}(\cdot)^{\rm
T}\widehat{\bm{\gamma}}({\bm{\beta}})$ via maximizing the approximate
loglikelihood

$\displaystyle
l({\bm{\beta}},{\bm{\gamma}})=\sum_{i=1}^{n}\left(y_{i}{\bm{\beta}}^{\rm
T}{\bf x}_{i}+{\bf B}(y_{i})^{\rm
T}{\bm{\gamma}}-\hbox{log}\left[\int_{0}^{1}\exp\\{y{\bm{\beta}}^{\rm T}{\bf
x}_{i}+{\bf B}(y)^{\rm T}{\bm{\gamma}}\\}dy\right]\right)$

with respect to ${\bm{\gamma}}$. Then estimate ${\bm{\beta}}$ by maximizing

$\displaystyle\widehat{l}({\bm{\beta}})\equiv\sum_{i=1}^{n}\left(y_{i}{\bm{\beta}}^{\rm
T}{\bf
x}_{i}+\widehat{c}(y_{i},{\bm{\beta}})-\hbox{log}\left[\int_{0}^{1}\exp\\{y{\bm{\beta}}^{\rm
T}{\bf x}_{i}+\widehat{c}(y,{\bm{\beta}})\\}dy\right]\right).$

We point out that profiling and performing maximization jointly with respect
to ${\bm{\beta}}$ and ${\bm{\gamma}}$ yield the same result.

In the following, we first aim at proving the convergence property of
$\widehat{c}(y,{\bm{\beta}})$ when ${\bm{\beta}}={\bm{\beta}}_{0}$, where we
let ${\bm{\beta}}_{0}$ denote the regression coefficient of the true
conditional density.

We assume the following regularity conditions.

  1. (C1)

The function $c(\cdot)\in C^{q}([0,1])$ where $q\geq 1$ and $c(0)=0$. The true
conditional density of $Y$ given ${\bf X}$, $f_{Y|{\bf X}}(y|{\bf x})$, has a
compact support $[0,1]$, is positive and bounded on its support. The marginal
density of ${\bf X}$, $f_{{\bf X}}({\bf x})$, has compact support
$\mathcal{X}$ and is bounded on its support.

  2. (C2)

The spline order $r\geq q$.

  3. (C3)

Define the knots
$t_{-r+1}=\dots=t_{0}=0<t_{1}<\dots<t_{N}<1=t_{N+1}=\dots=t_{N+r}$ where N is
the number of interior knots and [0,1] is divided into $(N+1)$ subintervals.
Let $m=N+r$. N satisfies $N\xrightarrow{}\infty$,
$N^{-q}n^{1/2}\xrightarrow{}0$, and
$N^{-1}n(\hbox{log}n)^{-1}\xrightarrow{}\infty$ as $n\xrightarrow{}\infty$.

  4. (C4)

Let $h_{p}$ be the distance between the $p$th and $(p+1)$th interior knots,
$h=\max_{r\leq p\leq N+r}h_{p}$ and $h^{\prime}=\min_{r\leq p\leq N+r}h_{p}$.
There exists a constant $C_{h}$ such that $0<C_{h}<\infty$ and
$h/h^{\prime}<C_{h}$. Therefore, $h\asymp N^{-1}$ and $h^{\prime}\asymp
N^{-1}$.

  5. (C5)

(De Boor, 1978) Under Conditions (C1)-(C4), there exists a spline coefficient
${\bm{\gamma}}_{0}$ with the first component $\gamma_{01}=0$ such that

$\displaystyle\sup_{y\in[0,1]}\left|{\bf B}(y)^{\rm
T}{\bm{\gamma}}_{0}-c(y)\right|=O(h^{q}).$

  6. (C6)

${\bf\Sigma}_{22}\equiv E[\hbox{var}\\{{\bf B}(Y)|{\bf X}\\}]$ is invertible
and $\|{\bf\Sigma}_{22}\|_{2}^{-1}{\bf\Sigma}_{22}$ has all eigenvalues
bounded above $C_{l}$, where $C_{l}>0$ is a constant.

Condition (C1) imposes the smoothness of the functional component $c(\cdot)$,
and the boundedness of the densities involved in the model, which are standard
requirements. The compact support requirement of the densities is also a
standard requirement in the B-spline literature. Condition (C2) requires that
the order of the B-spline basis is sufficiently large for the B-spline curve
to converge to the true function fast enough. We further assume that there are
an appropriate number of interior knots by Condition (C3), and that the knots
are uniformly distributed in the asymptotic sense by Condition (C4). Lastly,
we point out that Condition (C5) does not further impose any additional
requirement. It is a direct result given Conditions (C1)-(C4), and is only
stated as a condition for convenience.

######  Proposition 1.

Under Conditions (C1)-(C6),
$\|\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\|_{2}=O_{p}\\{(nh)^{-1/2}\\}$
and

$\displaystyle\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}$
$\displaystyle=$ $\displaystyle{\bf\Sigma}_{22}^{-1}n^{-1}\sum_{i=1}^{n}[{\bf
B}(y_{i})-E\\{{\bf B}(Y)|{\bf x}_{i}\\}]+{\bf r}_{1},$

where $\|{\bf r}_{1}\|_{2}=o_{p}\\{(nh)^{-1/2}\\}$. Furthermore,

$\displaystyle\sup_{y\in[0,1]}|\widehat{c}(y,{\bm{\beta}}_{0})-c(y)|=O_{p}\\{(nh)^{-1/2}+h^{q}\\}.$

Proposition 1 states that given the true regression coefficient
${\bm{\beta}}_{0}$, $\widehat{\bm{\gamma}}$ converges to the B-spline basis
coefficient ${\bm{\gamma}}_{0}$ at the nonparametric convergence rate
$(nh)^{-1/2}$, hence the B-spline curve approximates $c(\cdot)$. We further
obtain that the estimator $\widehat{c}(\cdot,{\bm{\beta}}_{0})$ converges to
the true function $c(\cdot)$ at the $(nh)^{-1/2}$ rate as well. Based on
Proposition 1, below we further establish the asymptotic property of the
estimator $\widehat{\bm{\beta}}$. Note ${\bm{\beta}}$ is not our research
interest, hence this property is stated as a by-product. We first impose one
additional regularity condition, which is a standard requirement.

  1. (C7)

The expectation of the conditional covariance of $\\{{\bf X}^{\rm T}Y,{\bf
B}^{\rm T}(Y)\\}^{\rm T}$ given ${\bf X}$, i.e.

$\displaystyle{\bf\Sigma}\equiv\begin{bmatrix}{\bf\Sigma}_{11}&{\bf\Sigma}_{12}\\\
{\bf\Sigma}_{21}&{\bf\Sigma}_{22}\end{bmatrix}\equiv E\begin{bmatrix}{\bf
X}{\bf X}^{\rm T}\hbox{var}(Y|{\bf X})&{\bf X}\hbox{cov}\\{Y,{\bf B}(Y)|{\bf
X}\\}\\\ \hbox{cov}\\{{\bf B}(Y),Y|{\bf X}\\}{\bf X}^{\rm T}&\hbox{var}\\{{\bf
B}(Y)|{\bf X}\\}\end{bmatrix},$

is invertible.

It is easy to see that the conditional covariance of $\\{{\bf X}^{\rm T}Y,{\bf
B}^{\rm T}(Y)\\}^{\rm T}$ given ${\bf X}$ is positive semidefinite. Thus,
${\bf\Sigma}$ is also positive semidefinite. Condition (C7) further requires
that ${\bf\Sigma}$ is positive definite, which is very mild. Note that this
guarantees ${\bf\Sigma}_{11}$ and the Schur complement of ${\bf\Sigma}_{22}$,
${\bf\Sigma}^{*}\equiv{\bf\Sigma}_{11}-{\bf\Sigma}_{12}{\bf\Sigma}_{22}^{-1}{\bf\Sigma}_{21}$,
are both positive definite.

######  Proposition 2.

Under Conditions (C1)-(C7),
$\|\widehat{\bm{\beta}}-{\bm{\beta}}_{0}\|_{2}=O_{p}(n^{-1/2})$ and

$\displaystyle\widehat{\bm{\beta}}-{\bm{\beta}}_{0}$ $\displaystyle=$
$\displaystyle{\bf\Sigma}^{*-1}n^{-1}\sum_{i=1}^{n}{\bf
x}_{i}\\{y_{i}-E(Y|{\bf x}_{i})\\}$
$\displaystyle-{\bf\Sigma}^{*-1}{\bf\Sigma}_{12}{\bf\Sigma}_{22}^{-1}n^{-1}\sum_{i=1}^{n}[{\bf
B}(y_{i})-E\\{{\bf B}(Y)|{\bf x}_{i}\\}]+{\bf r}_{2},$

where $\|{\bf r}_{2}\|_{2}=o_{p}(n^{-1/2})$. Furthermore,

$\displaystyle{\bf\Sigma}^{*1/2}\sqrt{n}(\widehat{\bm{\beta}}-{\bm{\beta}}_{0})\to
N({\bf 0}_{p},{\bf I}_{p})$

in distribution as $n\to\infty$.

Proposition 2 establishes how the regression coefficient estimator
$\widehat{\bm{\beta}}$ is asymptotically distributed, and allows to perform
inference since the asymptotic variance of $\widehat{\bm{\beta}}$ is estimable
based on $\widehat{\bm{\beta}}$ and $f_{Y\mid{\bf X}}(y,{\bf
x},\widehat{\bm{\beta}},\widehat{c})$. Below, we study how the estimator of
the marginal effect ${\bm{\xi}}$ is asymptotically distributed in Theorem 1.
We point out that the asymptotic variance of $\widehat{\bm{\xi}}$ is estimable
as well so that we can perform inference on the marginal effect.

######  Theorem 1.

Let ${\bf\Sigma}_{\bm{\xi}}\equiv{\bf A}{\bf\Sigma}^{-1}{\bf A}^{\rm
T}+{\bm{\beta}}_{0}{\bm{\beta}}_{0}^{\rm T}\hbox{var}\\{\hbox{var}(Y|{\bf
X})\\}$, where ${\bf A}\equiv[{\bf A}_{1},{\bf A}_{2}]$ and

$\displaystyle{\bf A}_{1}$ $\displaystyle\equiv$ $\displaystyle
E\\{\hbox{var}(Y|{\bf X})\\}{\bf I}+{\bm{\beta}}_{0}E[\\{Y-E(Y|{\bf
X})\\}^{3}{\bf X}^{\rm T}],$ $\displaystyle{\bf A}_{2}$ $\displaystyle\equiv$
$\displaystyle{\bm{\beta}}_{0}E(\\{Y-E(Y|{\bf X})\\}^{2}[{\bf B}(Y)-E\\{{\bf
B}(Y)|{\bf X}\\}])^{\rm T}.$

Under Conditions (C1)-(C7),

$\displaystyle{\bf\Sigma}_{\bm{\xi}}^{-1/2}\sqrt{n}(\widehat{\bm{\xi}}-{\bm{\xi}}_{0})\to
N({\bf 0}_{p},{\bf I}_{p})$

in distribution as $n\to\infty$.

To facilitate the analysis of the properties of the quantile effect estimator
$\widehat{\bm{\eta}}_{\tau}$, we assume two additional regularity conditions.

  1. (C8)

$c^{\prime\prime}(\cdot)$ is bounded on $[0,1]$.

  2. (C9)

The spline order $r\geq 2$.

Note that for an arbitrary twice differentiable function $g(\cdot)$,
$\|g^{\prime}(\cdot)\|_{\infty}\leq
2(\|g(\cdot)\|_{\infty}\|g^{\prime\prime}(\cdot)\|_{\infty})^{1/2}$ by the
Landau-Kolmogorov inequality. Setting $g(\cdot)={\bf B}(\cdot)^{\rm
T}{\bm{\gamma}}_{0}-c(\cdot)$, Conditions (C5), (C8), and (C9) guarantee
$\|{\bf B}^{\prime}(\cdot)^{\rm
T}{\bm{\gamma}}_{0}-c^{\prime}(\cdot)\|_{\infty}=O(h^{q/2})$, i.e. ${\bf
B}^{\prime}(\cdot)^{\rm T}{\bm{\gamma}}_{0}$ converges to $c^{\prime}(\cdot)$
uniformly at the rate $O(h^{q/2})$. Theorem 2 below provides the asymptotic
normality of the quantile effect estimator $\widehat{\bm{\eta}}_{\tau}$.

######  Theorem 2.

Let ${\bf\Sigma}_{{\bm{\eta}}_{\tau}}\equiv{\bf C}{\bf\Sigma}^{-1}{\bf C}^{\rm
T}+{\bm{\beta}}_{0}{\bm{\beta}}_{0}^{\rm T}\hbox{var}\\{q^{\prime}({\bf
X}^{\rm T}{\bm{\beta}}_{0})\\}$, where ${\bf C}\equiv[{\bf C}_{1},{\bf
C}_{2}]$ and

$\displaystyle{\bf C}_{1}$ $\displaystyle\equiv$ $\displaystyle
E\\{q^{\prime}({\bf X}^{\rm T}{\bm{\beta}}_{0})\\}{\bf
I}+{\bm{\beta}}_{0}E\Bigg{[}{\bf X}^{\rm T}\Bigg{\\{}\frac{E([\tau-I\\{Y\leq
q({\bf X}^{\rm T}{\bm{\beta}}_{0})\\}]Y^{2}|{\bf X})}{f_{Y\mid{\bf
X}}\\{q({\bf X}^{\rm T}{\bm{\beta}}_{0})|{\bf X}\\}}$
$\displaystyle-2q^{\prime}({\bf X}^{\rm T}{\bm{\beta}}_{0})q({\bf X}^{\rm
T}{\bm{\beta}}_{0})-\\{q^{\prime}({\bf X}^{\rm T}{\bm{\beta}}_{0})\\}^{2}[{\bf
X}^{\rm T}{\bm{\beta}}_{0}+c^{\prime}\\{q({\bf X}^{\rm
T}{\bm{\beta}}_{0})\\}]\\}],$ $\displaystyle{\bf C}_{2}$ $\displaystyle\equiv$
$\displaystyle{\bm{\beta}}_{0}E\Bigg{\\{}\frac{E([\tau-I\\{Y\leq q({\bf
X}^{\rm T}{\bm{\beta}}_{0})\\}]Y{\bf B}^{\rm T}(Y)|{\bf X})}{f_{Y\mid{\bf
X}}\\{q({\bf X}^{\rm T}{\bm{\beta}}_{0})|{\bf X}\\}}$
$\displaystyle-q^{\prime}({\bf X}^{\rm T}{\bm{\beta}}_{0}){\bf B}^{\rm
T}\\{q({\bf X}^{\rm T}{\bm{\beta}}_{0})\\}-\frac{E([\tau-I\\{Y\leq q({\bf
X}^{\rm T}{\bm{\beta}}_{0})\\}]{\bf B}^{\rm T}(Y)|{\bf X})}{f_{Y\mid{\bf
X}}\\{q({\bf X}^{\rm T}{\bm{\beta}}_{0})|{\bf X}\\}}$
$\displaystyle\times[q({\bf X}^{\rm T}{\bm{\beta}}_{0})+q^{\prime}({\bf
X}^{\rm T}{\bm{\beta}}_{0}){\bf X}^{\rm T}{\bm{\beta}}_{0}+q^{\prime}({\bf
X}^{\rm T}{\bm{\beta}}_{0})c^{\prime}\\{q({\bf X}^{\rm
T}{\bm{\beta}}_{0})\\}]\\}.$

Under Conditions (C1)-(C9),

$\displaystyle{\bf\Sigma}_{{\bm{\eta}}_{\tau}}^{-1/2}\sqrt{n}(\widehat{\bm{\eta}}_{\tau}-{\bm{\eta}}_{\tau
0})\to N({\bf 0}_{p},{\bf I}_{p})$

in distribution as $n\to\infty$.

The asymptotic properties we established, especially the estimation variances
of $\widehat{\bm{\xi}}$ and $\widehat{\bm{\eta}}_{\tau}$, have very different
forms from the efficiency bounds we derived in Sections 2.1 and 2.2.
Nevertheless, closer inspection, together with some basic but less frequently
adopted linear algebra tools reveal that these two sets of results have much
closer connections, and the efficiency bounds in estimating both ${\bm{\xi}}$
and ${\bm{\eta}}_{\tau}$ are actually reached by our B-spline based
approximate maximum likelihood estimators. As a by-product, we also state the
efficiency property of $\widehat{\bm{\beta}}$ as a proposition, even though
our interest is not in ${\bm{\beta}}$.

######  Proposition 3.

Under Conditions (C1)-(C7), the approximate maximum likelihood estimator
$\widehat{\bm{\beta}}$ is efficient.

######  Theorem 3.

Under Conditions (C1)-(C7), the estimator $\widehat{\bm{\xi}}$ based on the
approximate maximum likelihood estimator $\widehat{\bm{\beta}}$ and
$\widehat{c}(\cdot)$ is efficient.

######  Theorem 4.

Under Conditions (C1)-(C9), the estimator $\widehat{\bm{\eta}}_{\tau}$ based
on the approximate maximum likelihood estimator $\widehat{\bm{\beta}}$ and
$\widehat{c}(\cdot)$ is efficient.

###  3.2 Discrete Response

We now analyze the properties of our estimators under the discrete response
case. For notational simplicity, we denote
${\bm{\theta}}\equiv({\bm{\beta}}^{\rm T},{\bm{\gamma}}^{\rm T})^{\rm T}$,
${\bf p}({\bf x},{\bm{\theta}})\equiv\\{\hbox{pr}(Y=1\mid{\bf
x};{\bm{\beta}},{\bf c}),\dots,\hbox{pr}(Y=m\mid{\bf x};{\bm{\beta}},{\bf
c})\\}^{\rm T}$, and ${\bf p}_{k}({\bf
x},{\bm{\theta}})\equiv\\{1^{k}\hbox{pr}(Y=1\mid{\bf x};{\bm{\beta}},{\bf
c}),\dots,m^{k}\hbox{pr}(Y=m\mid{\bf x};{\bm{\beta}},{\bf c})\\}^{\rm T}$.
Also in this section, we use ${\bf B}(y)$ to denote a vector of indicator
functions, i.e., ${\bf B}(y)\equiv\\{I(y=1),\dots,I(y=m)\\}^{\rm T}$. Note
that we allow $m$ to grow with the sample size $n$. Here we present
theoretical results when $m$ grows to infinity. Note that in the finite $m$
case, the analysis can be done easily through incorporating the classical
maximum likelihood approach. We first list a set of regularity conditions.

  1. (D1)

The true conditional mass function of $Y$ given ${\bf X}$,
$\hbox{pr}(Y=y\mid{\bf x};{\bm{\theta}}_{0})$, has a support set
$\\{0,\dots,m\\}$. $E\\{\hbox{pr}(Y=0\mid{\bf X};{\bm{\theta}}_{0})\\}\leq
1-\delta$ for some constant $0<\delta<1$ and $E(Y^{4}\mid{\bf x})$ is bounded.
The marginal density of ${\bf X}$, $f_{{\bf X}}({\bf x})$, has compact support
$\mathcal{X}$ and is bounded on its support.

  2. (D2)

There exist constants $L_{k}$ such that $|E(Y^{k}\mid{\bf
x},{\bm{\theta}}^{*})-E(Y^{k}\mid{\bf x},{\bm{\theta}}_{0})|\leq
L_{k}\|{\bm{\theta}}^{*}-{\bm{\theta}}_{0}\|_{2}$ for $k=1,2,3$.

  3. (D3)

$m\to\infty$, $n^{-1}m^{3}\to 0$, and $E\\{\|{\bf p}({\bf
X},{\bm{\theta}}_{0})\|_{2}^{2}\\}\to 0$ as $n\to\infty$.

  4. (D4)

${\bf\Sigma}_{22}\equiv E[\hbox{var}\\{{\bf B}(Y)\mid{\bf X}\\}]$ is
invertible. $\|{\bf\Sigma}_{22}\|_{2}^{-1}{\bf\Sigma}_{22}$ has all
eigenvalues bounded above a constant $C_{l}>0$.

  5. (D5)

The expectation of the conditional covariance of $\\{{\bf X}^{\rm T}Y,{\bf
B}^{\rm T}(Y)\\}^{\rm T}$ given ${\bf X}$, i.e.

$\displaystyle{\bf\Sigma}\equiv\begin{bmatrix}{\bf\Sigma}_{11}&{\bf\Sigma}_{12}\\\
{\bf\Sigma}_{21}&{\bf\Sigma}_{22}\end{bmatrix}\equiv E\begin{bmatrix}{\bf
X}{\bf X}^{\rm T}\hbox{var}(Y\mid{\bf X})&{\bf X}\hbox{cov}\\{Y,{\bf
B}(Y)\mid{\bf X}\\}\\\ \hbox{cov}\\{{\bf B}(Y),Y\mid{\bf X}\\}{\bf X}^{\rm
T}&\hbox{var}\\{{\bf B}(Y)\mid{\bf X}\\}\end{bmatrix},$

is invertible.

Conditions (D1) and (D2) require boundedness and Lipschitz continuity on the
conditional moments of $Y$, which are standard requirements. Condition (D3)
requires $m$ to tend to infinity at the rate slower than $n^{1/3}$, and the
mass function not to concentrate on a finite subset of the support. Further,
since ${\bf\Sigma}_{22}$ and ${\bf\Sigma}$ are positive semidefinite by their
definitions, the invertibility imposed by Conditions (D4) and (D5) is very
mild. Lastly, we point out some results on ${\bf\Sigma}_{22}$ in Remark 1
below.

######  Remark 1.

Note that the sum of the eigenvalues of ${\bf\Sigma}_{22}$ is of constant
order by Conditions (D1) and (D3), because
$\hbox{trace}({\bf\Sigma}_{22})=E\left\\{1-\hbox{pr}(Y=0\mid{\bf
X};{\bm{\theta}}_{0})-\|{\bf p}({\bf X},\theta_{0})\|_{2}^{2}\right\\}$. Thus
we get $\|{\bf\Sigma}_{22}\|_{2}\asymp m^{-1}$ and
$\|{\bf\Sigma}_{22}^{-1}\|_{2}\asymp m$ because the eigenvalues are of the
same order by Condition (D4).

We now state the convergence rate and the asymptotic properties of the
estimators of the model parameters ${\bm{\beta}}$ and ${\bm{\gamma}}$. Our
analysis shows that the estimator of the regression coefficient
$\widehat{\bm{\beta}}$ achieves the parametric convergence rate under the
regularity conditions stated above. Further, we establish the asymptotic
distribution of $\widehat{\bm{\beta}}$ as a by-product, by which one can
perform inference on the regression coefficient. We formally state the results
in Proposition 4 below.

######  Proposition 4.

Under Conditions (D1)-(D5),
$\|\widehat{\bm{\beta}}-{\bm{\beta}}_{0}\|_{2}=O_{p}(n^{-1/2})$,
$\|\widehat{\bm{\gamma}}-{\bm{\gamma}}_{0}\|_{2}=O_{p}(n^{-1/2}m^{1/2})$, and

$\displaystyle\begin{bmatrix}\widehat{\bm{\beta}}-{\bm{\beta}}_{0}\\\
\widehat{\bm{\gamma}}-{\bm{\gamma}}_{0}\end{bmatrix}={\bf\Sigma}^{-1}n^{-1}\sum_{i=1}^{n}\begin{bmatrix}{\bf
x}_{i}\\{y_{i}-E(Y\mid{\bf x}_{i})\\}\\\ {\bf B}(y_{i})-E\\{{\bf B}(Y)\mid{\bf
x}_{i}\\}\end{bmatrix}+\begin{bmatrix}{\bf r}_{1}\\\ {\bf
r}_{2}\end{bmatrix},$

where $\|{\bf r}_{1}\|_{2}=o_{p}(n^{-1/2})$ and $\|{\bf
r}_{2}\|_{2}=o_{p}(n^{-1/2}m^{1/2})$. Furthermore, let
${\bf\Sigma}_{\bm{\beta}}\equiv({\bf\Sigma}_{11}-{\bf\Sigma}_{12}{\bf\Sigma}_{22}^{-1}{\bf\Sigma}_{21})^{-1}$,
then

$\displaystyle{\bf\Sigma}_{\bm{\beta}}^{-1/2}\sqrt{n}(\widehat{\bm{\beta}}-{\bm{\beta}}_{0})\to
N({\bf 0}_{p},{\bf I}_{p})$

in distribution as $n\to\infty$.

Based on Proposition 4, we further establish a theoretical result on the
marginal effect estimator, which states the convergence rate and the
asymptotic distribution of $\widehat{\bm{\xi}}$. Theorem 5 shows that
$\widehat{\bm{\xi}}$, the functional of both the parametric and the
nonparametric components, achieves the parametric convergence rate. We also
provide the closed form of the asymptotic variance of $\widehat{\bm{\xi}}$,
which can be used to infer the marginal effect of a population.

######  Theorem 5.

Let ${\bf\Sigma}_{\bm{\xi}}\equiv{\bf A}{\bf\Sigma}^{-1}{\bf A}^{\rm
T}+{\bm{\beta}}_{0}{\bm{\beta}}_{0}^{\rm T}\hbox{var}\\{\hbox{var}(Y\mid{\bf
X})\\}$, where ${\bf A}\equiv[{\bf A}_{1},{\bf A}_{2}]$ and

$\displaystyle{\bf A}_{1}$ $\displaystyle\equiv$ $\displaystyle
E\\{\hbox{var}(Y\mid{\bf X})\\}{\bf I}+{\bm{\beta}}_{0}E[\\{Y-E(Y\mid{\bf
X})\\}^{3}{\bf X}^{\rm T}],$ $\displaystyle{\bf A}_{2}$ $\displaystyle\equiv$
$\displaystyle{\bm{\beta}}_{0}E(\\{Y-E(Y\mid{\bf X})\\}^{2}[{\bf
B}(Y)-E\\{{\bf B}(Y)\mid{\bf X}\\}])^{\rm T}.$

Under Conditions (D1)-(D5),

$\displaystyle{\bf\Sigma}_{\bm{\xi}}^{-1/2}\sqrt{n}(\widehat{\bm{\xi}}-{\bm{\xi}}_{0})\to
N({\bf 0}_{p},{\bf I}_{p})$

in distribution as $n\to\infty$.

We now show that our proposed estimators $\widehat{\bm{\beta}}$ and
$\widehat{\bm{\xi}}$ for the discrete case also achieve the efficiency bounds.
Although the asymptotic variances of $\widehat{\bm{\beta}}$ and
$\widehat{\bm{\xi}}$ established in Proposition 4 and in Theorem 5 appear very
different from the efficiency bounds derived in Appendix S.9 and Section S.1,
our analysis shows that these two seemingly different variance structures are
actually identical. Below, we formally state the efficiency of
$\widehat{\bm{\xi}}$ as Theorem 6, and that of $\widehat{\bm{\beta}}$ as
Proposition 5.

######  Proposition 5.

Under Conditions (D1)-(D5), the approximate maximum likelihood estimator
$\widehat{\bm{\beta}}$ is efficient.

######  Theorem 6.

Under Conditions (D1)-(D5), the estimator $\widehat{\bm{\xi}}$ based on the
approximate maximum likelihood estimator $\widehat{\bm{\beta}}$ and
$\widehat{c}(\cdot)$ is efficient.

##  4 Simulation Experiments

We conduct simulation studies to investigate the finite sample performance of
the proposed methods. We consider the case where the response follows a GLM or
a truncated GLM. All results are based on 1000 replicates with sample size
$n=1000$. For comparison, we implemented three different estimators, our
proposed approximate maximum likelihood estimator (aMLE), the pairwise
marginal likelihood estimator (pMLE) by Lin et al. (2021), and the maximum
likelihood estimator (MLE) under the non-truncated regression model. Note that
MLE is the most efficient estimator if the response $Y$ is not truncated, but
wrongly specifies the distribution when the response $Y$ is actually
truncated.

For our proposed aMLE for the continuous response, we used the cubic B-spline
basis with the number of interior knots equal to the smallest integer larger
than $0.7n^{1/5}$, i.e., $N=\lceil 0.7n^{1/5}\rceil$, where the knots are the
quantiles of $\\{y_{i}:i=1,\dots,n\\}$ of length $(N+2)$ whose levels are
evenly spaced in $[0,1]$. On the other hand, because Lin et al. (2021) only
studied the estimation of ${\bm{\beta}}$ and $c(\cdot)$, we implemented the
estimators $\widehat{\bm{\xi}}$ and $\widehat{\bm{\eta}}_{\tau}$ of pMLE in
the same manner as for aMLE, based on the pMLE estimated ${\bm{\beta}}$ and
$c(\cdot)$.

We report the average of the absolute bias, the sample standard error
$\sigma_{\text{sim}}$, the average of the asymptotic standard error
$\widehat{\sigma}_{\text{est}}$, and the empirical coverage of the estimated
confidence interval at 95% confidence level (CI). The
$\widehat{\sigma}_{\text{est}}$ and CI of pMLE are omitted because Lin et al.
(2021) did not provide them.

###  4.1 Normal Distribution

We first examine the normal regression model. A three dimensional covariate
vector ${\bf X}_{i}$ was independently drawn from a multivariate normal
distribution with mean ${\bf 0}$ and covariance $\Sigma=(\sigma_{kl})$ where
$\sigma_{kl}=0.1^{|k-l|},k=1,\dots,3,l=1,\dots,3$. Then for the truncated
case, we further generated a response $Y_{i}$ independently from a truncated
normal distribution on $[a,b]$, which has the density

$\displaystyle
f_{Y}(y;\theta_{i},\sigma,a,b)=\frac{\frac{1}{\sigma}\phi(\frac{y-\theta_{i}}{\sigma})}{\Phi(\frac{b-\theta_{i}}{\sigma})-\Phi(\frac{a-\theta_{i}}{\sigma})},$

where $\theta_{i}={\bm{\beta}}^{\rm T}{\bf x}_{i},\sigma=1,a=-5,b=5$,
$\phi(\cdot)$ and $\Phi(\cdot)$ are the probability density function (pdf) and
the cumulative distribution function (cdf) of the standard normal
distribution. We set ${\bm{\beta}}=(1,2,3)^{\rm T}$ in the simulation
procedure. For the non-truncated case, we generated replicates with the same
parameters but without truncation. The non-truncated simulation design is
identical to that of Lin et al. (2021).

To evaluate the performance under the truncated case, we illustrate in Figure
1(a) how the estimated $c(\cdot)$ performed by our method while fixing
$c(-5)=0$. As can be seen from the plot, $\widehat{c}(\cdot)$ approximated
$c(\cdot)$ with satisfactory bias and variance at sample size 1000\. This
reflects the theoretical properties described in Lemma 6 of the supplementary
material. In addition, in Table 1, we illustrate the estimation properties for
${\bm{\beta}}$ and ${\bm{\xi}}$ via aMLE, pMLE and MLE. The aMLE and pMLE
methods were numerically almost identical in terms of both bias and standard
error, reflecting the theoretical properties established in Proposition 3 as
well as those in Theorem 1 of Lin et al. (2021). In terms of inference of
aMLE, the estimated standard error was very close to the sample standard
error, and the coverage rate of the confidence intervals was close to the
nominal level $95\%$. These indicate that our asymptotic properties are
already useful at sample size $n=1000$ in this model. In contrast, the
inference results of MLE were very bad, with the coverage rate almost 0. This
is a direct consequence of the estimation bias caused by model
misspecification. Table S.1 in the supplementary material provides the
estimation results of ${\bm{\eta}}_{\tau}$ at the quantile levels
$\tau=0.05,0.25,0.5,0.75$ and 0.95 respectively. Again, the estimation and
inference properties of aMLE were satisfactory, suggesting that the properties
described in Theorem 2 is reflected in this model for $n=1000$. Also, MLE
performed poorly in estimating ${\bm{\eta}}_{\tau}$ especially at the low or
high quantile level.

For the non-truncated case, the corresponding results are shown in Figure 1(b)
and Tables 2 and S.2 in the supplementary material. We can see from Figure
1(b) that $\widehat{c}(\cdot)$ estimated the curve $c(\cdot)$ sufficiently
well, even though the response is infinitely supported hence our compact
support assumption is violated. We also note the small biases of the
estimators $\widehat{\bm{\beta}}$, $\widehat{\bm{\xi}}$, and
$\widehat{\bm{\eta}}_{\tau}$ from Tables 2 and S.2. Moreover, interestingly,
our estimators of ${\bm{\xi}}$ and ${\bm{\eta}}_{\tau}$ performed as well as
MLE in this simulation setting, even when $\tau$ is near 0 or 1. In terms of
inference, the empirical coverage of the estimated confidence interval by our
method was close to the nominal level. These seem to suggest an empirical
robustness property of our aMLE against the compact support assumption.
Lastly, we point out that our method aMLE was more efficient than pMLE in this
example in terms of estimating ${\bm{\beta}},{\bm{\xi}}$ as well as the
quantile effect ${\bm{\eta}}_{\tau}$ at various quantile levels $\tau$.

###  4.2 Gamma Distribution

Next, we consider the situation when $Y$ given ${\bf x}$ has a gamma
distribution. We first generated a covariate vector ${\bf
X}_{i}=(X_{1i},X_{2i})^{\rm T}$, where $X_{ki},k=1,2$ are independently and
identically distributed (iid) as a uniform random variable on $[0.5,1]$. Then
for the truncated case, a response $Y_{i}$ was generated independently from a
truncated gamma distribution with the pdf

$\displaystyle f_{Y|{\bf
X}}(y;\alpha,\theta_{i},b)=\frac{y^{\alpha-1}e^{-y/\theta_{i}}}{\Gamma(\alpha)\theta_{i}^{\alpha}}\left\\{\int_{0}^{b}\frac{y^{\alpha-1}e^{-y/\theta_{i}}}{\Gamma(\alpha)\theta_{i}^{\alpha}}\right\\}^{-1},$

where $\alpha=5$, $\alpha\theta_{i}=1/{\bm{\beta}}^{\rm T}{\bf x}_{i}$, $b=2$,
and ${\bm{\beta}}=(0.5,1)^{\rm T}$. For the non-truncated case, we simply
carried out the same data generation mechanism while applying the usual gamma
distribution with the same parameters.

For the truncated case, Figure 1(c) illustrates the performance of the
estimator $\widehat{c}(\cdot)$ with $c(2)=0$. The estimation had very small
bias on most part of the support except when $y$ is close to 0 due to the
boundary effect. Indeed, the true curve $c(y)=(\alpha-1)\hbox{log}(y)$ is
unbounded near 0 and very few observations are available since the density
converges to 0\. Nevertheless, we can see from Tables 3 and 4 that aMLE
estimated ${\bm{\beta}}$, ${\bm{\xi}}$, and ${\bm{\eta}}_{\tau}$ well with
small bias, and its inference was also sufficiently precise with a good match
between the sample variance and the estimated variance, and the 95% confidence
intervals had coverage close to the nominal level. In contrast, MLE showed
bias and did not perform well in general, and deteriorated further when $\tau$
is near 1. In our application, we also find that pMLE was computationally
unstable in this setting and did not lead to reasonable results.

For the case when the response is not truncated, Figure 1(d) shows that
similar performance was observed as in the truncated gamma regression case. In
addition, Tables 5 and S.3 in the supplementary material suggest that aMLE
still provided good results in estimating the corresponding parameters and can
be used for reliable inference. In this situation, aMLE still outperformed
pMLE numerically even though both should be inconsistent in theory. Because
MLE assumes a fully parametric model which happens to be correct, it had the
best performance among all three methods as we expected.

###  4.3 Bernoulli Distribution

Our first simulation for the discrete response was carried out under a
conditional Bernoulli distribution with the intention to investigate the
performance of aMLE in the discrete response case. A three dimensional
covariate vector ${\bf X}_{i}$ was generated independently from a multivariate
normal distribution with mean ${\bf 0}$ and covariance $\Sigma=(\sigma_{kl})$
where $\sigma_{kl}=0.1^{|k-l|},k=1,\dots,3,l=1,\dots,3$. Then a binary
response $Y_{i}$ was drawn independently from a Bernoulli distribution with a
success rate $1/\\{1+\exp(-{\bm{\beta}}^{\rm T}{\bf x}_{i})\\}$. We set
${\bm{\beta}}=(-0.5,0.5,1)^{\rm T}$.

The performances of the estimators $\widehat{\bm{\beta}}$ and
$\widehat{\bm{\xi}}$ by the three methods are given in Table 6. All methods
performed similarly to each other, because the problem is fully parametric as
explained in Section 2.3.

###  4.4 Poisson and Negative Binomial Distributions

We now conduct simulation studies for the case where $Y$ is discrete with
infinite support. Specifically, $Y$ has Poisson and negative binomial
distributions respectively. We first generated a covariate vector ${\bf
X}_{i}=(X_{1i},X_{2i})^{\rm T}$, where $X_{ki},k=1,2$ are independently and
identically distributed (iid) uniform random variables on $[0.5,1]$. Then for
the Poisson regression model, a response $Y_{i}$ was generated from a Poisson
distribution with the rate $\theta_{i}=\exp({\bm{\beta}}^{\rm T}{\bf x}_{i})$
where ${\bm{\beta}}=(0,1)^{\rm T}$. For the negative binomial case, we
generated $Y_{i}$ independently from a negative binomial distribution with the
mass function

$\displaystyle\hbox{pr}(Y=y;r,\theta_{i})={y+r-1\choose
r-1}\theta_{i}^{y}(1-\theta_{i})^{r},$

where $r=2$, $\theta_{i}=\exp({\bm{\beta}}^{\rm T}{\bf x}_{i})$, and
${\bm{\beta}}=(0,-1)^{\rm T}$. In estimation, MLE assumes
$\theta_{i}=\exp({\bm{\beta}}^{\rm T}{\bf x}_{i})$ and $c(y)=-\hbox{log}(y!)$
hence results in the most efficient estimator under the Poission regression
among the three methods, but wrongly specifies the structure when the
conditional distribution is the negative binomial.

Tables 7 and 8 show the ${\bm{\beta}}$ and ${\bm{\xi}}$ estimation results
under the Poisson and the negative binomial regression models respectively. In
both cases, we found that aMLE estimated ${\bm{\xi}}$ more efficiently than
pMLE even though our method aMLE showed similar performance to pMLE in
estimating ${\bm{\beta}}$. Especially for the estimation of $\xi_{2}$ in Table
8, we can see that aMLE was about twice as efficient as pMLE in terms of both
bias and standard error. This supports that our estimation procedure is more
suitable for estimating ${\bm{\xi}}$ in finite sample situations than the
method proposed by Lin et al. (2021). Also, considering that the empirical
coverage of our approximate confidence interval reached the nominal level of
95%, the results confirm the asymptotic properties of $\widehat{\bm{\beta}}$
and $\widehat{\bm{\xi}}$ proposed in Proposition 4 and Theorem 5. On the other
hand, MLE in Table 7 indeed had the best performance in terms of both bias and
variance, benefiting from the perfectly correct parametric model setting.
However, we can observe from Table 8 that the estimators from MLE were
severely biased when $c(y)$ was misspecified.

##  5 Swiss Non-labor Income Data Analysis

We now analyze a data set concerning non-labor income situation in
Switzerland. The data set is publicly available from the R-package AER
(Kleiber & Zeileis, 2008), and consists of information collected from 871
married women randomly drawn from the representative health survey for
Switzerland (SOMIPOPS) in 1981 with income not abnormally low.

Our goal is to estimate the marginal effect and the quantile effect of
covariates on non-labor income (such as husband’s income). Thus, the non-labor
income in the log scale is used as response $Y$. All the other variables in
the dataset are considered as explanatory variables, and they include
Participation (taking the value 1 if the individual participated in the labor
market, 0 otherwise), Age (age in years), Age2 (squared age in years then
divided by 10), Education (years of formal education), Foreign (taking the
value 1 if the individual is a permanent foreign resident, 0 otherwise),
Youngkids (number of young children), Oldkids (number of older children). The
age category for children was decided by whether the age is under 7. For more
detailed description, see Gerfin (1996).

We implemented our method aMLE. For comparison, we also implemented the normal
and the gamma regression. We fitted the models without any transformation on
or interaction between the covariates, while the number of knots of the
B-spline basis is chosen in the same way as explained in Section 4.

In Figure 2, we can see that the fitted distribution from aMLE does not
resemble that of the normal regression, whose $c(y)=-y^{2}/2\sigma^{2}$ with
$\sigma>0$. It also does not resemble that of the gamma distribution whose
$c(y)=(\alpha-1)\hbox{log}(y)$ with $\alpha>0$. In Table 9, we compared the
information criteria AIC and BIC of the three models, and found that our aMLE
on (4) has the lowest AIC and BIC compared to the normal and the gamma
regression model, hence is the most suitable modeling choice. This suggests
that our model in (4) is the most suitable to use in analyzing the swiss non-
labor income data.

Table 10 shows the inference result in estimating the marginal effect. All
except one explanatory variables have similar levels of significance based on
all three methods. The exception is “Foreign”, which was selected as a
significant variable at 5% level in the analysis based on aMLE, but was
considered non-significant in both the normal and the gamma regression. In
other words, only by using (4) in combination with aMLE, we can conclude that
being a permanent foreign resident has a negative effect on the non-labor
income in Switzerland.

In Figure 3, we further illustrate the estimated $\widehat{\bm{\eta}}_{\tau}$
at the quantile level $\tau=0.05$, 0.25, 0.5, 0.75 and 0.95. We fixed the
significance level at 5%. The p-values of $\widehat{\bm{\eta}}_{\tau}$ were
largely similar to those of $\widehat{\bm{\xi}}$. However, as we can see from
the plot, the magnitude of the quantile effect increases when $\tau$ is near 0
or 1. That is, our model suggests that when the non-labor income level is
relatively low or high, the effect of the covariates on the non-labor income
is more extreme.

##  6 Discussion

We have proposed a B-spline based approximate maximum likelihood estimation
procedure to estimate both the marginal effect and the quantile effect in a
semiparametric generalized linear model. Compared to the classical GLM, our
model is more flexible hence less susceptible to model misspecification. The
estimators we proposed are shown to reach the semiparametric efficiency
bounds, hence are optimal.

Following the spirit of GLMs, we have worked with a linear summary of
covariates ${\bm{\beta}}^{\rm T}{\bf x}$ in our work. It is easy to see that
we can replace ${\bm{\beta}}^{\rm T}{\bf x}$ by a more general form $m({\bf
x},{\bm{\beta}})$, where $m$ is a known function which can be nonlinear. All
our procedures can be carried through while replacing ${\bf x}$ by $\partial
m({\bf x},{\bm{\beta}})/\partial{\bm{\beta}}$, where we only need to ensure
identifiability, sufficient smoothness, and boundedness of the corresponding
quantities to facilitate the almost identical derivation to the linear case.

An interesting but difficult further extension of our work worth mentioning is
the possibility of allowing non-compactly supported distribution of $Y$.
Although truncation is routinely done in practice to bypass the complexity, we
found it very difficult to rigorously extend our theoretical analysis to
encompass this scenario. We suspect that more fundamental work is needed,
possibly in a model simpler than the semiparametric GLM proposed in (4).

## References

  * (1)
  * Basu & Rathouz (2005) Basu, A. & Rathouz, P. J. (2005), ‘Estimating marginal and incremental effects on health outcomes using flexible link and variance function models’, Biostatistics 6, 93–109.
  * Bickel et al. (1998) Bickel, P. J., Klaassen, C. A. J., Ritov, Y. & Wellner, J. A. (1998), Efficient and Adaptive Estimation for Semiparametric Models, Springer, New York.
  * Cantoni & Ronchetti (2001) Cantoni, E. & Ronchetti, E. (2001), ‘Robust inference for generalized linear models’, Journal of the American Statistical Association 96, 1022–1030.
  * Cattaneo et al. (2010) Cattaneo, M. D., Crump, R. K. & Jansson, M. (2010), ‘Robust data-driven inference for density-weighted average derivatives’, Journal of the American Statistical Association 105(491), 1070–1083.
  * Cattaneo et al. (2013) Cattaneo, M. D., Crump, R. K. & Jansson, M. (2013), ‘Generalized jackknife estimators of weighted average derivatives’, Journal of the American Statistical Association 108(504), 1243–1256.
  * Chaudhuri et al. (1997) Chaudhuri, P., Doksum, K. & Samarov, A. (1997), ‘On average derivative quantile regression’, The Annals of Statistics 25(2), 715–744.
  * De Boor (1978) De Boor, C. (1978), A Practical Guide to Splines, Vol. 27, Springer, New York.
  * DeVore & Lorentz (1993) DeVore, R. A. & Lorentz, G. G. (1993), Constructive Approximation, Vol. 303, Springer, New York.
  * Fernihough (2019) Fernihough, A. (2019), ‘Marginal effects for generalized linear models: The mfx package for R’.
  * Gerfin (1996) Gerfin, M. (1996), ‘Parametric and semi-parametric estimation of the binary response model of labour market participation’, Journal of Applied Econometrics 11, 321–339.
  * Greene (2000) Greene, W. H. (2000), Econometric analysis, 4th edn, Prentice Hall, New Jersey.
  * Härdle & Stoker (1989) Härdle, W. & Stoker, T. M. (1989), ‘Investigating smooth multiple regression by the method of average derivatives’, Journal of the American Statistical Association 84(408), 986–995.
  * Hristache et al. (2001) Hristache, M., Juditsky, A. & Spokoiny, V. (2001), ‘Direct estimation of the index coefficient in a single-index model’, Annals of Statistics pp. 595–623.
  * Huang & Rathouz (2012) Huang, A. & Rathouz, P. J. (2012), ‘Proportional likelihood ratio models for mean regression’, Biometrika 99, 223–229.
  * Imbens & Angrist (1994) Imbens, G. W. & Angrist, J. D. (1994), ‘Identification and estimation of local average treatment effects’, Econometrica: Journal of the Econometric Society 62(2), 467–475.
  * Jiang et al. (2018) Jiang, F., Ma, Y. & Carroll, R. J. (2018), ‘A spline-assisted semiparametric approach to non-parametric measurement error models’, arXiv preprint arXiv:1804.00793 .
  * Kleiber & Zeileis (2008) Kleiber, C. & Zeileis, A. (2008), Applied econometrics with R, Springer, New York.
  * Künsch et al. (1989) Künsch, H., Stefanski, L. & Carroll, R. J. (1989), ‘Conditionally unbiased bounded-infuence estimation in general regression models with applications to generalized linear models’, Journal of the American Statistical Association 84, 460–466.
  * Lin et al. (2021) Lin, L., Liu, L., Cui, X. & Wang, K. (2021), ‘A generalized semiparametric regression and its efficient estimation’, Scandinavian Journal of Statistics 48, 1–24.
  * Luo & Tsai (2012) Luo, X. & Tsai, W. Y. (2012), ‘A proportional likelihood ratio model’, Biometrika 99, 211–222.
  * Manning et al. (2005) Manning, W. G., Basu, A. & Mullahy, J. (2005), ‘Generalized modeling approaches to risk adjustment of skewed outcomes data’, Journal of Health Economics 24, 465–488.
  * McCullagh & Nelder (1989) McCullagh, P. & Nelder, J. A. (1989), Generalized Linear Models, 2nd edn, Chapman & Hall, London.
  * Newey & Stoker (1993) Newey, W. K. & Stoker, T. M. (1993), ‘Efficiency of weighted average derivative estimators and index models’, Econometrica: Journal of the Econometric Society pp. 1199–1223.
  * Ning et al. (2017) Ning, Y., Zhao, T., Liu, H. et al. (2017), ‘A likelihood ratio framework for high-dimensional semiparametric regression’, Annals of Statistics 45, 2299–2327.
  * Norton et al. (2019) Norton, E. C., Dowd, B. E. & Maciejewski, M. L. (2019), ‘Marginal effects—quantifying the effect of changes in risk factors in logistic regression models’, Journal of the American Medical Association 321, 1304–1305.
  * Parzen (2004) Parzen, E. (2004), ‘Quantile probability and statistical data modeling’, Statistical Science 19, 652–662.
  * Stefanski et al. (1986) Stefanski, L., Carroll, R. J. & Ruppert, D. (1986), ‘Optimally bounded score functions for generalized linear models with applications to logistic regression’, Biometrika 73, 413–424.
  * Tsiatis (2006) Tsiatis, A. A. (2006), Semiparametric Theory and Missing Data, Springer, New York.

Figure 1: $c(\cdot)$ estimation results. Red: the true curve $c(\cdot)$;
Black: the median curve of $\widehat{c}(\cdot)$; Filled curves: the 2.5% and
97.5% quantiles of $\widehat{c}(\cdot)$.

(a) truncated normal

(b) normal

(c) truncated gamma

(d) gamma

Figure 2: $c(\cdot)$ estimation in the Swiss non-labor income data. Black: the curve $\widehat{c}(\cdot)$; Filled curves: the estimated pointwise confidence band of $c(\cdot)$. Figure 3: The result of ${\bm{\eta}}_{\tau}$ estimation in the Swiss non-labor income data. Table 1: ${\bm{\beta}}$ and ${\bm{\xi}}$ estimation results under the truncated normal distribution. |  $|$bias$|$ | $\sigma_{\text{sim}}$ | $\hat{\sigma}_{\text{est}}$ | C.I.
---|---|---|---|---
| aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE
$\beta_{1}$ | .049 | .049 | .263 | .061 | .061 | .035 | .060 | - | .035 | .947 | - | .000
$\beta_{2}$ | .083 | .084 | .530 | .103 | .103 | .041 | .101 | - | .036 | .951 | - | .000
$\beta_{3}$ | .118 | .118 | .793 | .145 | .145 | .044 | .145 | - | .035 | .949 | - | .000
$\xi_{1}$ | .022 | .023 | .028 | .028 | .029 | .035 | .029 | - | .035 | .948 | - | .951
$\xi_{2}$ | .028 | .027 | .033 | .035 | .034 | .041 | .035 | - | .036 | .941 | - | .912
$\xi_{3}$ | .032 | .031 | .036 | .039 | .038 | .044 | .043 | - | .035 | .964 | - | .881
Table 2: ${\bm{\beta}}$ and ${\bm{\xi}}$ estimation results under the normal distribution. |  $|$bias$|$ | $\sigma_{\text{sim}}$ | $\hat{\sigma}_{\text{est}}$ | C.I.
---|---|---|---|---
| aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE
$\beta_{1}$ | .046 | .049 | .026 | .057 | .064 | .032 | .055 | - | .032 | .938 | - | .954
$\beta_{2}$ | .081 | .088 | .025 | .099 | .114 | .032 | .096 | - | .032 | .946 | - | .957
$\beta_{3}$ | .116 | .127 | .025 | .143 | .166 | .032 | .139 | - | .032 | .940 | - | .941
$\xi_{1}$ | .026 | .029 | .026 | .032 | .038 | .032 | .032 | - | .032 | .953 | - | .954
$\xi_{2}$ | .025 | .030 | .025 | .032 | .038 | .032 | .035 | - | .032 | .961 | - | .957
$\xi_{3}$ | .026 | .030 | .025 | .032 | .038 | .032 | .038 | - | .032 | .964 | - | .941
Table 3: ${\bm{\beta}}$ and ${\bm{\xi}}$ estimation results under the truncated gamma distribution. |  $|$bias$|$ | $\sigma_{\text{sim}}$ | $\hat{\sigma}_{\text{est}}$ | C.I.
---|---|---|---|---
| aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE
$\beta_{1}$ | .097 | .114 | .074 | .124 | .141 | .074 | .122 | - | .076 | .944 | - | .889
$\beta_{2}$ | .101 | .117 | .062 | .128 | .145 | .077 | .128 | - | .077 | .942 | - | .945
$\xi_{1}$ | .064 | .284 | .108 | .081 | .154 | .059 | .080 | - | .061 | .942 | - | .589
$\xi_{2}$ | .064 | .523 | .117 | .080 | .126 | .061 | .079 | - | .061 | .947 | - | .531
Table 4: ${\bm{\eta}}_{\tau}$ estimation results under the truncated gamma distribution. $\tau$ |  | Method |  $|$bias$|$ | $\sigma_{\rm sim}$ | $\widehat{\sigma}_{\rm est}$ | C.I.
---|---|---|---|---|---|---
0.05 | $\eta_{1}$ | aMLE | 0.032 | 0.041 | 0.041 | 0.945
| pMLE | 0.120 | 0.077 | - | -
| MLE | 0.030 | 0.025 | 0.025 | 0.825
$\eta_{2}$ | aMLE | 0.038 | 0.048 | 0.048 | 0.945
| pMLE | 0.217 | 0.090 | - | -
| MLE | 0.023 | 0.026 | 0.027 | 0.943
0.25 | $\eta_{1}$ | aMLE | 0.052 | 0.066 | 0.065 | 0.942
| pMLE | 0.395 | 0.163 | - | -
| MLE | 0.048 | 0.041 | 0.041 | 0.826
$\eta_{2}$ | aMLE | 0.056 | 0.071 | 0.070 | 0.938
| pMLE | 0.750 | 0.127 | - | -
| MLE | 0.036 | 0.042 | 0.044 | 0.947
0.50 | $\eta_{1}$ | aMLE | 0.068 | 0.086 | 0.085 | 0.943
| pMLE | 0.436 | 0.192 | - | -
| MLE | 0.071 | 0.055 | 0.056 | 0.799
$\eta_{2}$ | aMLE | 0.071 | 0.089 | 0.089 | 0.939
| pMLE | 0.824 | 0.136 | - | -
| MLE | 0.055 | 0.058 | 0.059 | 0.919
0.75 | $\eta_{1}$ | aMLE | 0.083 | 0.106 | 0.105 | 0.950
| pMLE | 0.296 | 0.181 | - | -
| MLE | 0.118 | 0.073 | 0.074 | 0.669
$\eta_{2}$ | aMLE | 0.087 | 0.109 | 0.107 | 0.945
| pMLE | 0.529 | 0.144 | - | -
| MLE | 0.112 | 0.077 | 0.079 | 0.761
0.95 | $\eta_{1}$ | aMLE | 0.089 | 0.113 | 0.111 | 0.936
| pMLE | 0.120 | 0.109 | - | -
| MLE | 0.336 | 0.106 | 0.107 | 0.113
$\eta_{2}$ | aMLE | 0.089 | 0.112 | 0.113 | 0.943
| pMLE | 0.220 | 0.141 | - | -
| MLE | 0.496 | 0.111 | 0.115 | 0.008
Table 5: ${\bm{\beta}}$ and ${\bm{\xi}}$ estimation results under the gamma distribution. |  $|$bias$|$ | $\sigma_{\text{sim}}$ | $\hat{\sigma}_{\text{est}}$ | C.I.
---|---|---|---|---
| aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE
$\beta_{1}$ | .088 | .372 | .061 | .112 | .578 | .077 | .111 | - | .076 | .948 | - | .954
$\beta_{2}$ | .093 | .372 | .064 | .119 | .541 | .080 | .119 | - | .078 | .944 | - | .943
$\xi_{1}$ | .073 | 1.040 | .051 | .093 | .986 | .064 | .092 | - | .065 | .941 | - | .951
$\xi_{2}$ | .079 | 1.551 | .055 | .099 | .782 | .068 | .098 | - | .065 | .944 | - | .942
Table 6: ${\bm{\beta}}$ and ${\bm{\xi}}$ estimation results under the Bernoulli distribution. |  $|$bias$|$ | $\sigma_{\text{sim}}$ | $\hat{\sigma}_{\text{est}}$ | C.I.
---|---|---|---|---
| aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE
$\beta_{1}$ | .089 | .089 | .088 | .111 | .111 | .111 | .109 | - | .109 | .950 | - | .949
$\beta_{2}$ | .088 | .088 | .088 | .112 | .112 | .111 | .110 | - | .110 | .947 | - | .946
$\beta_{3}$ | .100 | .099 | .099 | .123 | .123 | .123 | .124 | - | .124 | .956 | - | .954
$\xi_{1}$ | .015 | .015 | .015 | .019 | .019 | .019 | .019 | - | .019 | .950 | - | .949
$\xi_{2}$ | .015 | .015 | .015 | .020 | .020 | .020 | .019 | - | .019 | .945 | - | .944
$\xi_{3}$ | .014 | .014 | .014 | .017 | .017 | .017 | .017 | - | .017 | .947 | - | .948
Table 7: ${\bm{\beta}}$ and ${\bm{\xi}}$ estimation results under the Poisson distribution. |  $|$bias$|$ | $\sigma_{\text{sim}}$ | $\hat{\sigma}_{\text{est}}$ | C.I.
---|---|---|---|---
| aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE
$\beta_{1}$ | .121 | .116 | .088 | .154 | .146 | .110 | .151 | - | .109 | .939 | - | .944
$\beta_{2}$ | .130 | .130 | .087 | .162 | .157 | .109 | .158 | - | .106 | .947 | - | .938
$\xi_{1}$ | .258 | .268 | .187 | .327 | .337 | .236 | .320 | - | .233 | .938 | - | .944
$\xi_{2}$ | .267 | .296 | .195 | .334 | .365 | .243 | .326 | - | .237 | .940 | - | .934
Table 8: ${\bm{\beta}}$ and ${\bm{\xi}}$ estimation results under the negative binomial distribution. |  $|$bias$|$ | $\sigma_{\text{sim}}$ | $\hat{\sigma}_{\text{est}}$ | C.I.
---|---|---|---|---
| aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE | aMLE | pMLE | MLE
$\beta_{1}$ | .088 | .075 | 1.388 | .112 | .095 | .171 | .113 | - | .114 | .956 | - | .000
$\beta_{2}$ | .103 | .158 | .383 | .130 | .115 | .176 | .128 | - | .120 | .951 | - | .202
$\xi_{1}$ | .336 | .395 | 2.535 | .426 | .502 | .340 | .430 | - | .223 | .958 | - | .000
$\xi_{2}$ | .424 | .886 | 2.708 | .528 | 1.008 | .322 | .518 | - | .217 | .947 | - | .000
Table 9: AIC and BIC in the Swiss non-labor income data. Method | AIC | BIC
---|---|---
aMLE | 538.800 | 600.806
Normal | 667.201 | 710.128
Gamma | 657.090 | 700.016
Table 10: ${\bm{\xi}}$ estimation results in the Swiss non-labor income data. “*” indicates the significance of the corresponding predictor at $5\%$ significance level. Variable | $\widehat{\xi}_{\rm{aMLE}}$ | p-value | $\widehat{\xi}_{\rm{Normal}}$ | p-value | $\widehat{\xi}_{\rm{Gamma}}$ | p-value
---|---|---|---|---|---|---
Participation | -0.133 |  $<$1e-3* | -0.130 |  $<$1e-3* | -0.130 |  $<$1e-3*
Age | 0.070 |  $<$1e-3* | 0.065 |  $<$1e-3* | 0.065 |  $<$1e-3*
Age2 | -0.008 |  $<$1e-3* | -0.008 |  $<$1e-3* | -0.008 |  $<$1e-3*
Education | 0.042 |  $<$1e-3* | 0.042 |  $<$1e-3* | 0.042 |  $<$1e-3*
Youngkids | 0.010 | 0.695 | 0.011 | 0.657 | 0.011 | 0.654
Oldkids | 0.019 | 0.147 | 0.023 | 0.091 | 0.022 | 0.095
Foreign | -0.079 | 0.017* | -0.060 | 0.064 | -0.062 | 0.053

## SUPPLEMENTARY MATERIAL

###  S.1 Derivation of the Efficiency Bound of Marginal Effect Estimation

Consider an arbitrary parametric submodel

$\displaystyle f_{{\bf X},Y}(y,{\bf x},{\bm{\theta}})=f_{{\bf X}}({\bf
x},{\bm{\alpha}})\frac{\exp\left\\{y{\bm{\beta}}^{\rm T}{\bf
x}+c(y,{\bm{\gamma}})\right\\}}{\int\exp\left\\{y{\bm{\beta}}^{\rm T}{\bf
x}+c(y,{\bm{\gamma}})\right\\}d\mu(y)},$ (S.1)

where ${\bm{\theta}}=({\bm{\alpha}}^{\rm T},{\bm{\beta}}^{\rm
T},{\bm{\gamma}}^{\rm T})^{\rm T}$. We can verify that the score functions
associated with an arbitrary ${\bm{\theta}}$ are ${\bf S}_{\bm{\alpha}}({\bf
x})={\bf a}({\bf x})$, where ${\bf a}({\bf x})$ can be any function that
satisfies $E\\{{\bf a}({\bf X})\\}={\bf 0}$, ${\bf S}_{\bm{\beta}}(y,{\bf
x})={\bf x}\\{y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf x})\\}$, and ${\bf
S}_{\bm{\gamma}}(y,{\bf x})={\bf a}(y)-E\\{{\bf a}(Y)\mid{\bm{\beta}}^{\rm
T}{\bf x}\\}$, where ${\bf a}(y)$ can be any function. We can verify that
${\bm{\xi}}={\bm{\beta}}E\\{v({\bm{\beta}}^{\rm T}{\bf X})\\}$ and

$\displaystyle\frac{\partial{\bm{\xi}}}{\partial{\bm{\alpha}}^{\rm T}}$
$\displaystyle=$ $\displaystyle{\bm{\beta}}E\\{v({\bm{\beta}}^{\rm T}{\bf
X}){\bf a}^{\rm T}({\bf X})\\},$
$\displaystyle\frac{\partial{\bm{\xi}}}{\partial{\bm{\beta}}^{\rm T}}$
$\displaystyle=$ $\displaystyle E\\{v({\bm{\beta}}^{\rm T}{\bf X})\\}{\bf
I}+{\bm{\beta}}E[{\bf X}^{\rm T}\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf
X})\\}^{3}],$
$\displaystyle\frac{\partial{\bm{\xi}}}{\partial{\bm{\gamma}}^{\rm T}}$
$\displaystyle=$ $\displaystyle{\bm{\beta}}E\left([{\bf a}(Y)-E\\{{\bf
a}(Y)\mid{\bm{\beta}}^{\rm T}{\bf X}\\}]^{\rm T}\\{Y-E(Y\mid{\bm{\beta}}^{\rm
T}{\bf X})\\}^{2}\right).$

Thus, following Bickel et al. (1998) and Tsiatis (2006), a possible influence
function is

$\displaystyle\bm{\phi}(y,{\bf x})$ $\displaystyle=$
$\displaystyle{\bm{\beta}}v({\bm{\beta}}^{\rm T}{\bf
x})-{\bm{\beta}}E\\{v({\bm{\beta}}^{\rm T}{\bf
X})\\}+{\bm{\beta}}\\{y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf x})\\}^{2}$
$\displaystyle-{\bm{\beta}}E[\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf
x})\\}^{2}\mid{\bm{\beta}}^{\rm T}{\bf x}]+{\bf B}E\\{v({\bm{\beta}}^{\rm
T}{\bf X})\\}[{\bf b}(y,{\bf x})-E\\{{\bf b}(Y,{\bf x})\mid{\bf x}\\}],$

where

$\displaystyle{\bf B}=\left(E[\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf X})\\}{\bf
b}(Y,{\bf X}){\bf X}^{\rm T}]\right)^{-1},$ (S.2)

and ${\bf b}(y,{\bf x})$ is such that

$\displaystyle E\\{{\bf b}(y,{\bf X})\mid y\\}=E[E\\{{\bf b}(Y,{\bf
X})\mid{\bf X}\\}\mid y].$ (S.3)

We can verify that $\partial{\bm{\xi}}/\partial{\bm{\theta}}^{\rm
T}=E(\bm{\phi}{\bf S}_{\bm{\theta}}^{\rm T})$ where ${\bf
S}_{\bm{\theta}}=({\bf S}_{\bm{\alpha}}^{\rm T},{\bf S}_{\bm{\beta}}^{\rm
T},{\bf S}_{\bm{\gamma}}^{\rm T})^{\rm T}$. Now, the tangent space is $\cal
T=\cal T_{\bm{\alpha}}\oplus(\cal T_{\bm{\beta}}+\cal T_{\bm{\gamma}})$ where

$\displaystyle\cal T_{\bm{\alpha}}$ $\displaystyle=$ $\displaystyle[{\bf
a}({\bf x}):E\\{{\bf a}({\bf X})\\}={\bf 0}],$ $\displaystyle\cal
T_{\bm{\beta}}$ $\displaystyle=$ $\displaystyle[{\bf M}{\bf
x}\\{y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf x})\\}:\forall{\bf M}],$
$\displaystyle\cal T_{\bm{\gamma}}$ $\displaystyle=$ $\displaystyle[{\bf
a}(y)-E\\{{\bf a}(Y)\mid{\bm{\beta}}^{\rm T}{\bf x}\\}:\forall{\bf a}(y)].$

Let

$\displaystyle{\bf M}$ $\displaystyle=$ $\displaystyle E\left({\bf
B}E\\{v({\bm{\beta}}^{\rm T}{\bf X})\\}{\bf b}(Y,{\bf X}){\bf X}^{\rm
T}\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf X})\\}-2{\bm{\beta}}{\bf X}^{\rm
T}Yv({\bm{\beta}}^{\rm T}{\bf X})\right.$ (S.4) $\displaystyle\left.-{\bf
a}(Y)\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf X})\\}{\bf X}^{\rm
T}\right)\left[E\\{{\bf X}{\bf X}^{\rm T}v({\bm{\beta}}^{\rm T}{\bf
X})\\}\right]^{-1},$

where ${\bf a}(y)$ is such that

$\displaystyle E[E\\{{\bf a}(Y)\mid{\bf X}\\}\mid y]-{\bf a}(y)$ (S.5)
$\displaystyle=$ $\displaystyle 2{\bm{\beta}}E[yE(Y\mid{\bf
X})-E\\{YE(Y\mid{\bf X})\mid{\bf X}\\}\mid y]+{\bf M}E[{\bf X}\\{y-E(Y\mid{\bf
X})\\}\mid y].$

Then the efficient influence function is

$\displaystyle\bm{\phi}_{\rm eff}(y,{\bf x})$ $\displaystyle=$
$\displaystyle{\bm{\beta}}v({\bm{\beta}}^{\rm T}{\bf
x})-{\bm{\beta}}E\\{v({\bm{\beta}}^{\rm T}{\bf
X})\\}+{\bm{\beta}}\\{y^{2}-E(Y^{2}\mid{\bm{\beta}}^{\rm T}{\bf x})\\}$
$\displaystyle+{\bf M}{\bf x}\\{y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf x})\\}+{\bf
a}(y)-E\\{{\bf a}(Y)\mid{\bm{\beta}}^{\rm T}{\bf x}\\},$

where ${\bf a}(y)$ satisfies (S.5), ${\bf M}$ is given in (S.4), ${\bf
b}(y,{\bf x})$ satisfies (S.3) and ${\bf B}$ is given in (S.2).

Note that once we have ${\bf a}(y)$, we can let ${\bf b}(y,{\bf
x})=2{\bm{\beta}}yE(Y\mid{\bm{\beta}}^{\rm T}{\bf x})+{\bf M}{\bf x}y+{\bf
a}(y)$ and it satisfies (S.3). Using this specific ${\bf b}(y,{\bf x})$
function, we obtain that ${\bf M}$ and ${\bf B}$ need to satisfy

$\displaystyle{\bf B}=\left(E[\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf
X})\\}\\{2{\bm{\beta}}YE(Y\mid{\bm{\beta}}^{\rm T}{\bf X})+{\bf M}{\bf
X}Y+{\bf a}(Y)\\}{\bf X}^{\rm T}]\right)^{-1}$

and

$\displaystyle{\bf M}E\\{{\bf X}{\bf X}^{\rm T}v({\bm{\beta}}^{\rm T}{\bf
X})\\}=E\\{v({\bm{\beta}}^{\rm T}{\bf X})\\}{\bf I}-E\left[2{\bm{\beta}}{\bf
X}^{\rm T}Yv({\bm{\beta}}^{\rm T}{\bf X})+{\bf
a}(Y)\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf X})\\}{\bf X}^{\rm T}\right].$

This leads to

$\displaystyle{\bf M}$ $\displaystyle=$
$\displaystyle\left(E\\{v({\bm{\beta}}^{\rm T}{\bf X})\\}{\bf
I}-E\left[2{\bm{\beta}}{\bf X}^{\rm T}Yv({\bm{\beta}}^{\rm T}{\bf X})+{\bf
a}(Y)\\{Y-E(Y\mid{\bm{\beta}}^{\rm T}{\bf X})\\}{\bf X}^{\rm T}\right]\right)$
(S.6) $\displaystyle\times\left[E\\{{\bf X}{\bf X}^{\rm T}v({\bm{\beta}}^{\rm
T}{\bf X})\\}\right]^{-1},$

where ${\bf a}(y)$ satisfies (S.5).

Gathering the above derivations and results, we obtain the summary description
of the efficient influence function as

$\displaystyle\bm{\phi}_{\rm eff}(y,{\bf x})={\bm{\beta}}v({\bm{\beta}}^{\rm
T}{\bf x})-{\bm{\beta}}E\\{v({\bm{\beta}}^{\rm T}{\bf
X})\\}+{\bm{\beta}}y^{2}+{\bf M}{\bf x}y+{\bf a}(y)-E\\{{\bm{\beta}}Y^{2}+{\bf
M}{\bf x}Y+{\bf a}(Y)\mid{\bf x}\\},$

where ${\bf a}(y)$ satisfies (S.5) and ${\bf M}$ is given in (S.6). Obviously,
the variance of the efficient influence function, i.e.
$\hbox{var}\\{\bm{\phi}_{\rm eff}(Y,{\bf X})\\}$, is the efficiency bound in
estimating ${\bm{\xi}}$.

###  S.2 Derivation of the Efficiency Bound of Quantile Effect Estimation

Recall that for an arbitrary parametric submodel in (S.1), we have already
derived the corresponding score functions ${\bf S}_{\bm{\alpha}}({\bf x})$,
${\bf S}_{\bm{\beta}}(y,{\bf x})$ and ${\bf S}_{\bm{\gamma}}(y,{\bf x})$. Now
for notational brevity, let $\nu\equiv{\bm{\beta}}^{\rm T}{\bf x}$,
$q(\nu)\equiv Q_{\tau}(Y\mid\nu)$, $\epsilon\equiv\tau-I\\{Y<q(\nu)\\}$, and
$\epsilon^{\prime}_{\nu}=-\delta\\{q(\nu)-Y\\}q^{\prime}(\nu)$. Write
$f(y,\nu)\equiv f_{Y\mid{\bf X}}(y,\nu)$. Using the quantile definition, we
further have

$\displaystyle E(\epsilon\mid\nu)$ $\displaystyle=$ $\displaystyle 0,$
$\displaystyle E(\epsilon Y\mid\nu)$ $\displaystyle=$ $\displaystyle
f\\{q(\nu),\nu\\}q^{\prime}(\nu),$ $\displaystyle\frac{E(\epsilon
Y^{2}\mid\nu)}{f\\{q(\nu),\nu\\}}$ $\displaystyle=$ $\displaystyle
2q(\nu)q^{\prime}(\nu)+q^{\prime}(\nu)^{2}[\nu+c^{\prime}\\{q(\nu)\\}]+q^{\prime\prime}(\nu),$

which can be verified based on

$\displaystyle\tau$ $\displaystyle=$ $\displaystyle E[I\\{Y\leq
q(\nu)\\}\mid\nu],$ $\displaystyle q^{\prime}(\nu)$ $\displaystyle=$
$\displaystyle\frac{E\left([\tau-I\\{Y<q(\nu)\\}]Y\mid\nu\right)}{f\\{q(\nu),\nu\\}},$
$\displaystyle q^{\prime\prime}(\nu)$ $\displaystyle=$
$\displaystyle\frac{E\left([\tau-I\\{Y<q(\nu)\\}]Y^{2}\mid\nu\right)}{f\\{q(\nu),\nu\\}}-q(\nu)q^{\prime}(\nu)-q^{\prime}(\nu)[q(\nu)+q^{\prime}(\nu)\nu+q^{\prime}(\nu)c^{\prime}\\{q(\nu)\\}].$

We also have

$\displaystyle\frac{\partial q(\nu)}{\partial{\bm{\gamma}}}$ $\displaystyle=$
$\displaystyle\frac{E\left\\{\epsilon{\bf
a}(Y)\mid\nu\right\\}}{f\\{q(\nu),\nu\\}},$ $\displaystyle\frac{\partial
q^{\prime}(v)}{\partial{\bm{\gamma}}}$ $\displaystyle=$
$\displaystyle\frac{E\left\\{\epsilon Y{\bf
a}(Y)\mid\nu\right\\}}{f\\{q(\nu),\nu\\}}-{\bf
a}\\{q(\nu)\\}q^{\prime}(\nu)-\frac{\partial
q(\nu)}{\partial{\bm{\gamma}}}\left[q(\nu)+q^{\prime}(\nu)\nu+q^{\prime}(\nu)c^{\prime}\\{q(\nu)\\}\right].$

Note that ${\bm{\eta}}_{\tau}={\bm{\beta}}E\\{q^{\prime}(\nu)\\}$. We can
verify that

$\displaystyle\frac{\partial{\bm{\eta}}_{\tau}}{\partial{\bm{\alpha}}^{\rm
T}}$ $\displaystyle=$ $\displaystyle{\bm{\beta}}E\\{q^{\prime}(\nu){\bf
a}^{\rm T}({\bf X})\\},$
$\displaystyle\frac{\partial{\bm{\eta}}_{\tau}}{\partial{\bm{\beta}}^{\rm T}}$
$\displaystyle=$ $\displaystyle E\\{q^{\prime}(\nu)\\}{\bf
I}+{\bm{\beta}}E[r(Y,\nu){\bf X}^{\rm T}\\{Y-E(Y\mid\nu)\\}],$
$\displaystyle\frac{\partial{\bm{\eta}}_{\tau}}{\partial{\bm{\gamma}}^{\rm
T}}$ $\displaystyle=$ $\displaystyle{\bm{\beta}}E\left(r(Y,\nu)[{\bf a}^{\rm
T}(Y)-E\\{{\bf a}^{\rm T}(Y)\mid{\bf X}\\}]\right),$

where

$\displaystyle r(Y,\nu)\equiv\frac{\epsilon
Y+\epsilon^{\prime}_{\nu}-\epsilon[q(\nu)+q^{\prime}(\nu)\nu+q^{\prime}(\nu)c^{\prime}\\{q(\nu)\\}]}{f\\{q(\nu),\nu\\}}.$
(S.7)

Hence, a possible influence function is

$\displaystyle\bm{\phi}(y,{\bf x})$ $\displaystyle=$
$\displaystyle{\bm{\beta}}q^{\prime}(\nu)-{\bm{\beta}}E\\{q^{\prime}(\nu)\\}+{\bm{\beta}}[r(Y,\nu)-E\\{r(Y,\nu)\mid\nu\\}]$
$\displaystyle+E\\{q^{\prime}(\nu)\\}{\bf B}[{\bf b}(y,{\bf x})-E\\{{\bf
b}(Y,{\bf x})\mid{\bf x}\\}],$

where

$\displaystyle{\bf B}\equiv\left(E[{\bf b}(Y,{\bf X}){\bf X}^{\rm
T}\\{Y-E(Y\mid\nu)\\}]\right)^{-1},$ (S.8)

and ${\bf b}(y,{\bf x})$ is such that

$\displaystyle E\\{{\bf b}(y,{\bf X})\mid y\\}=E[E\\{{\bf b}(Y,{\bf
X})\mid{\bf X}\\}\mid y].$ (S.9)

We can verify that $\partial{\bm{\eta}}_{\tau}/\partial{\bm{\theta}}^{\rm
T}=E(\bm{\phi}{\bf S}_{\bm{\theta}}^{\rm T})$. Note that we have derived the
tangent space $\cal T$ in Section 2.1. Let $\bm{\rho}(y,{\bf x})\equiv
E\\{q^{\prime}(\nu)\\}{\bf B}{\bf b}(y,{\bf x})+{\bm{\beta}}r(y,\nu)-{\bf
a}(y)$, $v(\nu)\equiv E(Y^{2}\mid\nu)-\\{E(Y\mid\nu)\\}^{2}$, and

$\displaystyle{\bf M}_{1}\equiv E\left[\bm{\rho}(Y,{\bf X}){\bf X}^{\rm
T}\\{Y-E(Y\mid\nu)\\}\right]\left[E\\{{\bf X}{\bf X}^{\rm
T}v(\nu)\\}\right]^{-1},$ (S.10)

where ${\bf a}(y)$ is such that

$\displaystyle E[E\\{{\bf a}(Y)\mid{\bf X}\\}\mid y]-{\bf
a}(y)=-{\bm{\beta}}E\\{r(y,\nu)|y\\}+{\bf M}_{1}E[{\bf X}\\{y-E(Y|{\bf
X})\\}|y].$ (S.11)

Then the efficient influence function is

$\displaystyle\bm{\phi}_{\rm eff}(y,{\bf
x})={\bm{\beta}}q^{\prime}(\nu)-{\bm{\beta}}E\\{q^{\prime}(\nu)\\}+{\bf
M}_{1}{\bf x}y+{\bf a}(y)-E\\{{\bf M}_{1}{\bf x}Y+{\bf a}(Y)|{\bf x}\\},$
(S.12)

where ${\bf a}(y)$ satisfies (S.11), ${\bf M}_{1}$ is given in (S.10), ${\bf
b}(y,{\bf x})$ satisfies (S.9) and ${\bf B}$ is given in (S.8). Note that once
we have ${\bf a}(y)$, we can let ${\bf b}(y,{\bf x})={\bm{\beta}}r(y,\nu)-{\bf
M}_{1}{\bf x}y-{\bf a}(y)$ and it satisfies (S.9). Using this specific ${\bf
b}(y,{\bf x})$ function, we obtain that ${\bf M}_{1}$ and ${\bf B}$ need to
satisfy

$\displaystyle{\bf B}=\left(E[\\{{\bm{\beta}}r(Y,\nu)-{\bf M}_{1}{\bf X}Y-{\bf
a}(Y)\\}{\bf X}^{\rm T}\\{Y-E(Y\mid\nu)\\}]\right)^{-1}$

and

$\displaystyle{\bf M}_{1}E\\{{\bf X}{\bf X}^{\rm
T}v(\nu)\\}=E\\{q^{\prime}(\nu)\\}{\bf I}+{\bm{\beta}}E\\{{\bf X}^{\rm
T}q^{\prime\prime}(\nu)\\}-E[{\bf a}(Y){\bf X}^{\rm T}\\{Y-E(Y|{\bf X})\\}].$

This leads to

$\displaystyle{\bf M}_{1}=(E\\{q^{\prime}(\nu)\\}{\bf I}+{\bm{\beta}}E\\{{\bf
X}^{\rm T}q^{\prime\prime}(\nu)\\}-E[{\bf a}(Y){\bf X}^{\rm T}\\{Y-E(Y|{\bf
X})\\}])[E\\{{\bf X}{\bf X}^{\rm T}v(\nu)\\}]^{-1},$ (S.13)

where ${\bf a}(y)$ satisfies (S.11).

Thus, in summary, the efficient influence function for estimating
${\bm{\eta}}_{\tau}$ is given in (S.12), where ${\bf a}(y)$ satisfies (S.11),
$r(Y,\nu)$ is given in (S.7), and ${\bf M}_{1}$ is given in (S.13).

###  S.3 Lemmas

We first prove several lemmas that will be useful.

######  Lemma 1.

Let $U\equiv\\{{\bf u}\in\mathbb{R}^{m}:\|{\bf u}\|_{2}=1\\}$ and ${\bf u}\in
U$. Under Conditions (C2)-(C4), $\|{\bf B}(\cdot)^{\rm T}{\bf u}\|_{2}\asymp
h^{1/2}$, $\|{\bf B}(\cdot)^{\rm T}{\bf u}\|_{1}=O(h^{1/2})$ and $\|{\bf
B}(\cdot)^{\rm T}{\bf e}_{k}\|_{1}\asymp h$.

###### Proof.

We recall Lemma 1 in the supplement of Jiang et al. (2018), which is a direct
result from Theorem 5.4.2 on page 145 of DeVore & Lorentz (1993). For each
spline $\sum_{k=1}^{m}\gamma_{k}B_{k}(y)$ and $1\leq p\leq\infty$, there
exists a constant $C_{r}>0$ such that

$\displaystyle
C_{r}\lVert{\bm{\gamma}}^{\prime}\rVert_{p}\leq\left\lVert\sum_{k=1}^{m}\gamma_{k}B_{k}(\cdot)\right\rVert_{p}\leq\lVert{\bm{\gamma}}^{\prime}\rVert_{p},$

where
${\bm{\gamma}}^{\prime}\equiv\\{\gamma_{k}\\{(t_{k}-t_{k-r})/r\\}^{1/p},k=1,\dots,m\\}^{\rm
T}$. This implies $\|{\bf B}(\cdot)^{\rm T}{\bf u}\|_{p}\asymp\|{\bf
u}^{\prime}\|_{p}$ where

$\displaystyle\|{\bf
u}^{\prime}\|_{p}=\left[\sum_{k=1}^{m}\left\\{u_{k}\left(\frac{t_{k}-t_{k-r}}{r}\right)^{1/p}\right\\}^{p}\right]^{1/p}\asymp
h^{1/p}\left(\sum_{k=1}^{m}u_{k}^{p}\right)^{1/p}=h^{1/p}\|{\bf u}\|_{p}$

by Conditions (C3) and (C4). Hence we have $\|{\bf B}(\cdot)^{\rm T}{\bf
u}\|_{2}\asymp h^{1/2}$. In addition, note that $\|{\bf
u}\|_{1}=O(m^{1/2})\asymp O(h^{-1/2})$ under Conditions (C3) and (C4), which
leads to

$\displaystyle\|{\bf B}(\cdot)^{\rm T}{\bf u}\|_{1}\asymp h\times
O(h^{-1/2})=O(h^{1/2}).$

Since $\|{\bf e}_{k}\|_{1}=1$, we get $\|{\bf B}(\cdot)^{\rm T}{\bf
e}_{k}\|_{1}\asymp h$. ∎

######  Lemma 2.

Let $U\equiv\\{{\bf u}\in\mathbb{R}^{m}:\|{\bf u}\|_{2}=1\\}$ and ${\bf u}\in
U$. Under Conditions (C1)-(C4),

$\displaystyle E[\\{{\bf B}(Y)^{\rm T}{\bf u}\\}^{2}|{\bf x}]$
$\displaystyle\asymp$ $\displaystyle h,$ $\displaystyle E\\{{\bf B}(Y)^{\rm
T}{\bf u}|{\bf x}\\}$ $\displaystyle=$ $\displaystyle O(h^{1/2}),$
$\displaystyle E\\{{\bf B}(Y)^{\rm T}{\bf e}_{k}|{\bf x}\\}$
$\displaystyle\asymp$ $\displaystyle h,$ $\displaystyle\hbox{cov}\\{Y,{\bf
B}(Y)^{\rm T}{\bf u}|{\bf x}\\}$ $\displaystyle=$ $\displaystyle O(h^{1/2})$

uniformly in ${\bf x}$. For fixed ${\bm{\beta}}$ and ${\bm{\gamma}}$, we have
the same result under the density $f_{Y\mid{\bf X}}^{*}(y|{\bf
x},{\bm{\beta}},{\bm{\gamma}})$.

###### Proof.

The result for $f_{Y\mid{\bf X}}^{*}(y|{\bf x},{\bm{\beta}},{\bm{\gamma}})$
can be obtained similarly so we omit the proof. We write $E[\\{{\bf B}(Y)^{\rm
T}{\bf u}\\}^{2}|{\bf x}]$ as

$\displaystyle E[\\{{\bf B}(Y)^{\rm T}{\bf u}\\}^{2}|{\bf x}]$
$\displaystyle=$ $\displaystyle\int_{0}^{1}\left\\{{\bf B}(y)^{\rm T}{\bf
u}\right\\}^{2}f_{Y\mid{\bf X}}(y|{\bf x})dy$ $\displaystyle=$
$\displaystyle\left\lVert{\bf B}(\cdot)^{\rm T}{\bf u}\\{f_{Y\mid{\bf
X}}(\cdot|{\bf x})\\}^{1/2}\right\rVert_{2}^{2}.$

Note that under Condition (C1), there exist constants $0<c_{f}\leq
C_{f}<\infty$ such that

$\displaystyle c_{f}\leq\inf_{y\in[0,1],{\bf x}\in\mathcal{X}}f_{Y\mid{\bf
X}}(y|{\bf x})\leq\sup_{y\in[0,1],{\bf x}\in\mathcal{X}}f_{Y\mid{\bf
X}}(y|{\bf x})\leq C_{f},$

which implies

$\displaystyle c_{f}\left\lVert{\bf B}(\cdot)^{\rm T}{\bf
u}\right\rVert_{2}^{2}$ $\displaystyle\leq$ $\displaystyle\inf_{{\bf
x}\in\mathcal{X}}\left\lVert{\bf B}(\cdot)^{\rm T}{\bf u}\\{f_{Y\mid{\bf
X}}(\cdot|{\bf x})\\}^{1/2}\right\rVert_{2}^{2}$ $\displaystyle\leq$
$\displaystyle\sup_{{\bf x}\in\mathcal{X}}\left\lVert{\bf B}(\cdot)^{\rm
T}{\bf u}\\{f_{Y\mid{\bf X}}(\cdot|{\bf x})\\}^{1/2}\right\rVert_{2}^{2}$
$\displaystyle\leq$ $\displaystyle C_{f}\left\lVert{\bf B}(\cdot)^{\rm T}{\bf
u}\right\rVert_{2}^{2}.$

Then by Lemma 1, we have $E[\\{{\bf B}(Y)^{\rm T}{\bf u}\\}^{2}|{\bf x}]\asymp
h$. Further, $E\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf x}\\}$ satisfies

$\displaystyle E\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf x}\\}=\int_{0}^{1}{\bf
B}(y)^{\rm T}{\bf u}f_{Y\mid{\bf X}}(y|{\bf x})dy\leq C_{f}\left\lVert{\bf
B}(\cdot)^{\rm T}{\bf u}\right\rVert_{1}.$

By Lemma 1, the above implies $E\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}\\}=O(h^{1/2})$. Similarly, we get $E\\{{\bf B}(Y)^{\rm T}{\bf e}_{k}|{\bf
x}\\}\asymp h$. Now, $\hbox{cov}\\{Y,{\bf B}(Y)^{\rm T}{\bf u}|{\bf x}\\}$
satisfies

$\displaystyle\hbox{cov}\\{Y,{\bf B}(Y)^{\rm T}{\bf u}|{\bf x}\\}$
$\displaystyle=$ $\displaystyle\int_{0}^{1}\\{y-E(Y|{\bf x})\\}{\bf B}(y)^{\rm
T}{\bf u}f_{Y\mid{\bf X}}(y|{\bf x})dy$ $\displaystyle\leq$
$\displaystyle\|{\bf B}(\cdot)^{\rm T}{\bf u}\|_{1}\sup_{y\in[0,1],{\bf
x}\in\mathcal{X}}|y-E(Y|{\bf x})|f_{Y\mid{\bf X}}(y|{\bf x})$ $\displaystyle=$
$\displaystyle O(h^{1/2}).$

The last equality holds because $\sup_{y\in[0,1],{\bf
x}\in\mathcal{X}}|y-E(Y|{\bf x})|f_{Y\mid{\bf X}}(y|{\bf x})$ is bounded by
Condition (C1) and $\|{\bf B}(\cdot)^{\rm T}{\bf u}\|_{1}=O(h^{1/2})$ by Lemma
1. ∎

######  Lemma 3.

Let $U\equiv\\{{\bf u}\in\mathbb{R}^{m}:\|{\bf u}\|_{2}=1\\}$ and ${\bf u}\in
U$. Under Conditions (C1)-(C4),

$\displaystyle\sup_{{\bf u}\in U}\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}\\}\asymp h$

uniformly in ${\bf x}$. For fixed ${\bm{\beta}}$ and ${\bm{\gamma}}$, we have
$\sup_{{\bf u}\in U}\hbox{var}^{*}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x},{\bm{\beta}},{\bm{\gamma}}\\}\asymp h$ uniformly in ${\bf x}$.

###### Proof.

We omit the proof for $\hbox{var}^{*}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x},{\bm{\beta}},{\bm{\gamma}}\\}$ because we can obtain the result by a
similar way. The variance $\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf x}\\}$
satisfies

$\displaystyle\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf x}\\}$
$\displaystyle=$ $\displaystyle E[\\{{\bf B}(Y)^{\rm T}{\bf u}\\}^{2}|{\bf
x}]-[E\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf x}\\}]^{2}$ $\displaystyle=$
$\displaystyle h-O(h)$ $\displaystyle=$ $\displaystyle O(h)$

uniformly in ${\bf x}$ by Lemma 2. Similarly, we get $\hbox{var}\\{{\bf
B}(Y)^{\rm T}{\bf e}_{k}|{\bf x}\\}\asymp h-h^{2}\asymp h$ by Lemma 2. Then
since ${\bf e}_{k}\in U$ and $U$ is compact, we get

$\displaystyle\sup_{{\bf u}\in U}\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}\\}\asymp h$

uniformly in ${\bf x}$ under Condition (C1). ∎

######  Lemma 4.

Let $g_{1}(\cdot)$ be a function such that $\|g_{1}(\cdot)\|_{1}<\infty$.
Under Conditions (C1)-(C5),

$\displaystyle\sup_{y\in[0,1],{\bf x}\in\mathcal{X}}|f_{Y\mid{\bf
X}}^{*}(y,{\bf x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})-f_{Y\mid{\bf X}}(y|{\bf
x})|$ $\displaystyle=$ $\displaystyle O(h^{q}),$ $\displaystyle\sup_{{\bf
x}\in\mathcal{X}}|E^{*}\\{g_{1}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}-E\\{g_{1}(Y)|{\bf x}\\}|$
$\displaystyle=$ $\displaystyle O(h^{q})\|g_{1}(\cdot)\|_{1}.$

In addition, if $g_{2}(\cdot)$ satisfies $\|g_{2}(\cdot)\|_{1}<\infty$ and
$\|g_{1}(\cdot)g_{2}(\cdot)\|_{1}<\infty$,

$\displaystyle\sup_{{\bf
x}\in\mathcal{X}}|\hbox{cov}^{*}\\{g_{1}(Y),g_{2}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}-\hbox{cov}\\{g_{1}(Y),g_{2}(Y)|{\bf
x}\\}|$ $\displaystyle=$ $\displaystyle
O(h^{q})\\{\|g_{1}(\cdot)g_{2}(\cdot)\|_{1}+\|g_{1}(\cdot)\|_{1}\|g_{2}(\cdot)\|_{1}\\}.$

###### Proof.

Conditions (C1)-(C5) implies there exists a constant $0<C_{f}^{*}<\infty$ such
that

$\displaystyle\sup_{y\in[0,1],{\bf x}\in\mathcal{X}}|f_{Y\mid{\bf
X}}^{*}(y,{\bf x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})-f_{Y\mid{\bf X}}(y|{\bf
x})|\leq C_{f}^{*}h^{q}.$

This implies

$\displaystyle\sup_{{\bf x}\in\mathcal{X}}|E^{*}\\{g_{1}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}-E\\{g_{1}(Y)|{\bf x}\\}|$
$\displaystyle=$ $\displaystyle\sup_{{\bf
x}\in\mathcal{X}}\left|\int_{0}^{1}g_{1}(y)\\{f_{Y\mid{\bf X}}^{*}(y,{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})-f_{Y\mid{\bf X}}(y|{\bf
x})\\}dy\right|$ $\displaystyle\leq$
$\displaystyle\|g_{1}(\cdot)\|_{1}\sup_{y\in[0,1],{\bf
x}\in\mathcal{X}}|f_{Y\mid{\bf X}}^{*}(y,{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})-f_{Y\mid{\bf X}}(y|{\bf x})|$
$\displaystyle\leq$ $\displaystyle C_{f}^{*}h^{q}\|g_{1}(\cdot)\|_{1}$
$\displaystyle=$ $\displaystyle O(h^{q})\|g_{1}(\cdot)\|_{1}.$

Similarly, we have

$\displaystyle\sup_{{\bf x}\in\mathcal{X}}|E^{*}\\{g_{1}(Y)g_{2}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}-E\\{g_{1}(Y)g_{2}(Y)|{\bf
x}\\}|=O(h^{q})\|g_{1}(\cdot)g_{2}(\cdot)\|_{1}.$

Now, note that under Condition (C1), there exists a constant $0<C_{f}<\infty$
such that

$\displaystyle\sup_{y\in[0,1],{\bf x}\in\mathcal{X}}f_{Y\mid{\bf X}}(y|{\bf
x})\leq C_{f}.$

Then we get

$\displaystyle|E^{*}\\{g_{1}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}E^{*}\\{g_{2}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}-E\\{g_{1}(Y)|{\bf
x}\\}E\\{g_{2}(Y)|{\bf x}\\}|$ $\displaystyle=$
$\displaystyle|O(h^{q})E\\{g_{1}(Y)|{\bf
x}\\}\|g_{2}(\cdot)\|_{1}+O(h^{q})\|g_{1}(\cdot)\|_{1}E\\{g_{2}(Y)|{\bf
x}\\}+O(h^{2q})\|g_{1}(\cdot)\|_{1}\|g_{2}(\cdot)\|_{1}|$ $\displaystyle\leq$
$\displaystyle
O(h^{q})\\{2C_{f}+O(h^{q})\\}\|g_{1}(\cdot)\|_{1}\|g_{2}(\cdot)\|_{1}$

uniformly in ${\bf x}$ by Condition (C1). Therefore,

$\displaystyle\sup_{{\bf
x}\in\mathcal{X}}|\hbox{cov}^{*}\\{g_{1}(Y),g_{2}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}-\hbox{cov}\\{g_{1}(Y),g_{2}(Y)|{\bf
x}\\}|$ $\displaystyle\leq$ $\displaystyle\sup_{{\bf
x}\in\mathcal{X}}|E^{*}\\{g_{1}(Y)g_{2}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}-E\\{g_{1}(Y)g_{2}(Y)|{\bf x}\\}|$
$\displaystyle+\sup_{{\bf x}\in\mathcal{X}}|E^{*}\\{g_{1}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}E^{*}\\{g_{2}(Y)|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}-E\\{g_{1}(Y)|{\bf
x}\\}E\\{g_{2}(Y)|{\bf x}\\}|$ $\displaystyle=$ $\displaystyle
O(h^{q})\\{\|g_{1}(\cdot)g_{2}(\cdot)\|_{1}+\|g_{1}(\cdot)\|_{1}\|g_{2}(\cdot)\|_{1}\\}.$

∎

######  Lemma 5.

Let $g_{1}(\cdot)$ be a function which satisfies
$\|g_{1}(\cdot)\|_{1}<\infty$. Let ${\bm{\beta}}^{*},{\bm{\gamma}}^{*}$
satisfy $\|{\bm{\beta}}^{*}-{\bm{\beta}}_{0}\|_{2}=o_{p}(1)$ and
$\|{\bm{\gamma}}^{*}-{\bm{\gamma}}_{0}\|_{2}=o_{p}(1)$. Under Conditions
(C1)-(C5),

$\displaystyle\sup_{y\in[0,1],{\bf x}\in\mathcal{X}}|f_{Y\mid{\bf
X}}^{*}(y,{\bf x},{\bm{\beta}}^{*},{\bm{\gamma}}^{*})-f_{Y\mid{\bf X}}(y|{\bf
x})|$ $\displaystyle=$ $\displaystyle o_{p}(1),$ $\displaystyle\sup_{{\bf
x}\in\mathcal{X}}|E^{*}\\{g_{1}(Y)|{\bf
x},{\bm{\beta}}^{*},{\bm{\gamma}}^{*}\\}-E\\{g_{1}(Y)|{\bf x}\\}|$
$\displaystyle=$ $\displaystyle o_{p}(1)\|g_{1}(\cdot)\|_{1}.$

In addition, if $g_{2}(\cdot)$ satisfies $\|g_{2}(\cdot)\|_{1}<\infty$ and
$\|g_{1}(\cdot)g_{2}(\cdot)\|_{1}<\infty$, then

$\displaystyle\sup_{{\bf
x}\in\mathcal{X}}|\hbox{cov}^{*}\\{g_{1}(Y),g_{2}(Y)|{\bf
x},{\bm{\beta}}^{*},{\bm{\gamma}}^{*}\\}-\hbox{cov}\\{g_{1}(Y),g_{2}(Y)|{\bf
x}\\}|$ $\displaystyle=$ $\displaystyle
o_{p}(1)\\{\|g_{1}(\cdot)g_{2}(\cdot)\|_{1}+\|g_{1}(\cdot)\|_{1}\|g_{2}(\cdot)\|_{1}\\}.$

###### Proof.

The result can be obtained similarly to Lemma 4 by noting that under
Conditions (C1)-(C5), we have

$\displaystyle\sup_{y\in[0,1],{\bf x}\in\mathcal{X}}|f_{Y\mid{\bf
X}}^{*}(y,{\bf x},{\bm{\beta}}^{*},{\bm{\gamma}}^{*})-f_{Y\mid{\bf X}}(y|{\bf
x})|$ $\displaystyle\leq$ $\displaystyle\sup_{y\in[0,1],{\bf
x}\in\mathcal{X}}|f_{Y\mid{\bf X}}^{*}(y,{\bf
x},{\bm{\beta}}^{*},{\bm{\gamma}}^{*})-f_{Y\mid{\bf X}}^{*}(y,{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})|$ $\displaystyle+\sup_{y\in[0,1],{\bf
x}\in\mathcal{X}}|f_{Y\mid{\bf X}}^{*}(y,{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})-f_{Y\mid{\bf X}}(y|{\bf x})|$
$\displaystyle=$ $\displaystyle o_{p}(1).$

∎

###  S.4 Proof of Proposition 1

First, we will show that
$\|\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\|_{2}=O_{p}\\{(nh)^{-1/2}\\}$.
Notice that the Hessian of the loglikelihood with respect to ${\bm{\gamma}}$

$\displaystyle\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}=-\sum_{i=1}^{n}\hbox{cov}^{*}\\{{\bf B}(Y),{\bf B}(Y)|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}\\}$

is negative definite. Hence, it suffices to show the existence of a local
maximizer inside a ball centered at ${\bm{\gamma}}_{0}$ with radius of order
$(nh)^{-1/2}$. We will show this by proving for any $\epsilon>0$, there exists
a constant $C>0$ such that

$\displaystyle\hbox{pr}\left[l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})>\sup_{\|{\bf
v}\|_{2}=C}l\\{{\bm{\beta}}_{0},{\bm{\gamma}}_{0}+(nh)^{-1/2}{\bf
v}\\}\right]\geq 1-2\epsilon$ (S.14)

for sufficiently large $n$.

Let ${\bf v}\in\mathbb{R}^{m}$ be an arbitrary vector with $\|{\bf v}\|_{2}=C$
and ${\bm{\gamma}}^{*}\equiv\alpha\\{{\bm{\gamma}}_{0}+(nh)^{-1/2}{\bf
v}\\}+(1-\alpha){\bm{\gamma}}_{0}$ for some $\alpha\in(0,1)$. By the Taylor
expansion, we have

$\displaystyle l\\{{\bm{\beta}}_{0},{\bm{\gamma}}_{0}+(nh)^{-1/2}{\bf
v}\\}-l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})=(nh)^{-1/2}\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\gamma}}^{\rm T}}{\bf
v}+\frac{1}{2}(nh)^{-1}{\bf v}^{\rm
T}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}{\bf v}.$ (S.15)

We first analyze $\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})/\partial{\bm{\gamma}}$. Let ${\bf u}\in
U$ where $U\equiv\\{{\bf u}\in\mathbb{R}^{m}:\|{\bf u}\|_{2}=1\\}$, then

$\displaystyle n^{-1}\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\gamma}}^{\rm T}}{\bf u}$
$\displaystyle=$ $\displaystyle n^{-1}\sum_{i=1}^{n}[{\bf B}(y_{i})^{\rm
T}{\bf u}-E^{*}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}]$ (S.16) $\displaystyle=$
$\displaystyle n^{-1}\sum_{i=1}^{n}[{\bf B}(y_{i})^{\rm T}{\bf u}-E\\{{\bf
B}(Y)^{\rm T}{\bf u}|{\bf x}_{i}\\}]$
$\displaystyle+n^{-1}\sum_{i=1}^{n}[E\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}_{i}\\}-E^{*}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}].$

The first term of (S.16) satisfies

$\displaystyle E\left\\{\left(n^{-1}\sum_{i=1}^{n}\left[{\bf B}(Y_{i})^{\rm
T}{\bf u}-E\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
X}_{i}\\}\right]\right)^{2}\right\\}$ $\displaystyle=$ $\displaystyle
n^{-1}E\left(\left[{\bf B}(Y)^{\rm T}{\bf u}-E\\{{\bf B}(Y)^{\rm T}{\bf
u}|{\bf X}\\}\right]^{2}\right)$ $\displaystyle=$ $\displaystyle
n^{-1}E[\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf X}\\}].$

Furthermore, by Lemma 3 we have $\sup_{{\bf u}\in U}E[\hbox{var}\\{{\bf
B}(Y)^{\rm T}{\bf u}|{\bf X}\\}]\asymp h$, which implies

$\displaystyle\left\|n^{-1}\sum_{i=1}^{n}[{\bf B}(y_{i})-E\\{{\bf B}(Y)|{\bf
x}_{i}\\}]\right\|_{2}$ $\displaystyle=$ $\displaystyle\sup_{{\bf u}\in
U}\left|n^{-1}\sum_{i=1}^{n}[{\bf B}(y_{i})^{\rm T}{\bf u}-E\\{{\bf B}(Y)^{\rm
T}{\bf u}|{\bf x}_{i}\\}]\right|$ (S.17) $\displaystyle\asymp_{p}$
$\displaystyle n^{-1/2}h^{1/2}.$

Now, the second term of (S.16) satisfies

$\displaystyle\left|n^{-1}\sum_{i=1}^{n}[E\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}_{i}\\}-E^{*}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}]\right|$ $\displaystyle\leq$
$\displaystyle n^{-1}\sum_{i=1}^{n}\sup_{{\bf x}\in\mathcal{X}}|E\\{{\bf
B}(Y)^{\rm T}{\bf u}|{\bf x}\\}-E^{*}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}|$ $\displaystyle=$ $\displaystyle
O(h^{q})\|{\bf B}(\cdot)^{\rm T}{\bf u}\|_{1}$ $\displaystyle=$ $\displaystyle
O(h^{q+1/2})$

by Lemma 1 and 4. Hence, we get

$\displaystyle\left\|n^{-1}\sum_{i=1}^{n}[E\\{{\bf B}(Y)|{\bf
x}_{i}\\}-E^{*}\\{{\bf B}(Y)|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}_{0}\\}]\right\|_{2}$ $\displaystyle=$
$\displaystyle O(h^{q+1/2})$ (S.18) $\displaystyle=$ $\displaystyle
o(n^{-1/2}h^{1/2})$

under Conditions (C3) and (C4). Therefore, by (S.16), (S.17) and (S.18), we
get

$\displaystyle\left\|\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\gamma}}}\right\|_{2}\asymp_{p}(nh)^{1/2},$

hence for any $\epsilon>0$, there exists a constant $0<C_{1}<\infty$ such that

$\displaystyle\hbox{pr}\left\\{\left\|\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\gamma}}}\right\|_{2}\leq
C_{1}(nh)^{1/2}\right\\}\geq 1-\epsilon.$ (S.19)

Next, to analyze $\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})/\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}$, we have

$\displaystyle{\bf u}^{\rm
T}\left\\{-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\\}{\bf u}$ $\displaystyle=$ $\displaystyle
n^{-1}\sum_{i=1}^{n}\hbox{var}^{*}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}^{*}\\}$ (S.20) $\displaystyle=$
$\displaystyle n^{-1}\sum_{i=1}^{n}[\hbox{var}^{*}\\{{\bf B}(Y)^{\rm T}{\bf
u}|{\bf x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}^{*}\\}-\hbox{var}\\{{\bf
B}(Y)^{\rm T}{\bf u}|{\bf x}_{i}\\}]$
$\displaystyle+\left(n^{-1}\sum_{i=1}^{n}\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf
u}|{\bf x}_{i}\\}-E[\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf X}\\}]\right)$
$\displaystyle+E[\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf X}\\}].$

Note that
$\|{\bm{\gamma}}^{*}-{\bm{\gamma}}_{0}\|_{2}=O\\{(nh)^{-1/2}\\}=o(1)$ by
construction. Then by Lemma 1 and 5, the first term satisfies

$\displaystyle\left|n^{-1}\sum_{i=1}^{n}[\hbox{var}^{*}\\{{\bf B}(Y)^{\rm
T}{\bf u}|{\bf x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}^{*}\\}-\hbox{var}\\{{\bf
B}(Y)^{\rm T}{\bf u}|{\bf x}_{i}\\}]\right|$ (S.21) $\displaystyle\leq$
$\displaystyle n^{-1}\sum_{i=1}^{n}\sup_{{\bf
x}\in\mathcal{X}}|\hbox{var}^{*}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}^{*}\\}-\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf
u}|{\bf x}\\}|$ $\displaystyle=$ $\displaystyle o_{p}(1)\left[\|\\{{\bf
B}(\cdot)^{\rm T}{\bf u}\\}^{2}\|_{1}+\|{\bf B}(\cdot)^{\rm T}{\bf
u}\|_{1}^{2}\right]$ $\displaystyle=$ $\displaystyle o_{p}(1)\left\\{\|{\bf
B}(\cdot)^{\rm T}{\bf u}\|_{2}^{2}+\|{\bf B}(\cdot)^{\rm T}{\bf
u}\|_{1}^{2}\right\\}$ $\displaystyle=$ $\displaystyle o_{p}(h).$

For the second term of (S.20), first let $w_{i}\equiv\hbox{var}\\{{\bf
B}(Y)^{\rm T}{\bf u}|{\bf x}_{i}\\}$. Note that by Lemma 3, there exists a
constant $C$ so that $0\leq w_{i}\leq Ch$ and $0\leq E(W_{i})\leq Ch$ for all
$i=1,\dots,n$. Therefore, there exists a constant $C_{w}$, where
$0<C_{w}<\infty$, such that $|w_{i}-E(W_{i})|\leq C_{w}h$ for all
$i=1,\dots,n$. This also directly implies $\hbox{var}(W_{i})\leq
C_{w}^{2}h^{2}$. Then, by the Bernstein’s inequality, we have

$\displaystyle\hbox{pr}\left\\{n^{-1}\sum_{i=1}^{n}w_{i}-E(W_{i})\geq\epsilon\right\\}$
$\displaystyle\leq$
$\displaystyle\exp\left\\{-\frac{n^{2}\epsilon^{2}/2}{\sum_{i=1}^{n}\hbox{var}(W_{i})+C_{w}hn\epsilon/3}\right\\}$
$\displaystyle\leq$
$\displaystyle\exp\left\\{-\frac{n^{2}\epsilon^{2}/2}{nC_{w}^{2}h^{2}+C_{w}hn\epsilon/3}\right\\}$
$\displaystyle\to$ $\displaystyle 0$

for $\epsilon=h(\hbox{log}n)^{-1}$ as $n\to\infty$. Note that $\epsilon=o(h)$.
Thus, we get

$\displaystyle n^{-1}\sum_{i=1}^{n}w_{i}-E(W_{i})$ $\displaystyle=$
$\displaystyle n^{-1}\sum_{i=1}^{n}\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf
x}_{i}\\}-E[\hbox{var}\\{{\bf B}(Y)^{\rm T}{\bf u}|{\bf X}\\}]$ (S.22)
$\displaystyle=$ $\displaystyle o_{p}(h).$

Now for the last term of (S.20), we have

$\displaystyle\|{\bf\Sigma}_{22}\|_{2}$ $\displaystyle=$
$\displaystyle\sup_{{\bf u}\in U}{\bf u}^{\rm T}{\bf\Sigma}_{22}{\bf u}$
(S.23) $\displaystyle=$ $\displaystyle\sup_{{\bf u}\in U}E[\hbox{var}\\{{\bf
B}(Y)^{\rm T}{\bf u}|{\bf X}\\}]$ $\displaystyle\asymp$ $\displaystyle h$

by Lemma 3. Hence Condition (C6) implies that ${\bf\Sigma}_{22}$ has all
eigenvalues of order $h$. Moreover, from (S.20), (S.21), (S.22) and (S.23), we
have

$\displaystyle{\bf u}^{\rm
T}\left\\{-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\\}{\bf u}={\bf u}^{\rm T}{\bf\Sigma}_{22}{\bf u}+o_{p}(h)$ (S.24)

uniformly for all ${\bf u}$. Thus, all eigenvalues of
$\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})/\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}$ are negative and of order $nh$ with probability approaching 1, which
implies there exists a constant $0<C_{2}<\infty$ such that

$\displaystyle\hbox{pr}\left\\{{\bf v}^{\rm
T}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}{\bf v}\leq-C^{2}C_{2}nh\right\\}\geq 1-\epsilon.$ (S.25)

Combining (S.15), (S.19), and (S.25), we get

$\displaystyle l\\{{\bm{\beta}}_{0},{\bm{\gamma}}_{0}+(nh)^{-1/2}{\bf
v}\\}-l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})$ $\displaystyle\leq$
$\displaystyle
C(nh)^{-1/2}C_{1}(nh)^{1/2}+\frac{1}{2}(nh)^{-1}(-C^{2}C_{2}nh)$
$\displaystyle=$ $\displaystyle C\left(C_{1}-\frac{C_{2}}{2}C\right)$
$\displaystyle<$ $\displaystyle 0$

with probability at least $1-2\epsilon$ when $C>2C_{1}/C_{2}$. This proves
(S.14), hence we have
$\|\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\|_{2}=O_{p}\\{(nh)^{-1/2}\\}$.

Now we analyze the asymptotic behavior of
$\widehat{\bm{\gamma}}({\bm{\beta}}_{0})$. Note that
$\widehat{\bm{\gamma}}({\bm{\beta}}_{0})$ maximizes
$l({\bm{\beta}}_{0},{\bm{\gamma}})$, then by using the Taylor expansion with
${\bm{\gamma}}^{*}$ on the line connecting
$\widehat{\bm{\gamma}}({\bm{\beta}}_{0})$ and ${\bm{\gamma}}_{0}$, we get

$\displaystyle{\bf 0}$ $\displaystyle=$ $\displaystyle\frac{\partial
l\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}}{\partial{\bm{\gamma}}}$
$\displaystyle=$ $\displaystyle\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\gamma}}}+\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\\{\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\\}.$

This leads to

$\displaystyle\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}$
$\displaystyle=$
$\displaystyle\left\\{-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\\}^{-1}n^{-1}\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\gamma}}}$
$\displaystyle=$ $\displaystyle{\bf\Sigma}_{22}^{-1}n^{-1}\sum_{i=1}^{n}[{\bf
B}(y_{i})-E\\{{\bf B}(Y)|{\bf x}_{i}\\}]+{\bf r}_{1},$

where

$\displaystyle{\bf r}_{1}$ $\displaystyle\equiv$
$\displaystyle\left[\left\\{-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\\}^{-1}-{\bf\Sigma}_{22}^{-1}\right]n^{-1}\sum_{i=1}^{n}[{\bf
B}(y_{i})-E\\{{\bf B}(Y)|{\bf x}_{i}\\}]$
$\displaystyle+\left\\{-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\\}^{-1}\left(n^{-1}\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\gamma}}}-n^{-1}\sum_{i=1}^{n}[{\bf
B}(y_{i})-E\\{{\bf B}(Y)|{\bf x}_{i}\\}]\right).$

Note that

$\displaystyle\left\|\left\\{-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\\}^{-1}-{\bf\Sigma}_{22}^{-1}\right\|_{2}$ $\displaystyle=$
$\displaystyle\left\|\left\\{-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\\}^{-1}\left\\{{\bf\Sigma}_{22}+n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\\}{\bf\Sigma}_{22}^{-1}\right\|_{2}$ $\displaystyle\leq$
$\displaystyle\left\|\left\\{-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\\}^{-1}\right\|_{2}\left\|{\bf\Sigma}_{22}+n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\|_{2}\left\|{\bf\Sigma}_{22}^{-1}\right\|_{2}$
$\displaystyle\asymp_{p}$ $\displaystyle h^{-1}\times o_{p}(h)\times h^{-1}$
$\displaystyle=$ $\displaystyle o_{p}(h^{-1}),$

by (S.23), (S.24) and Condition (C6). Furthermore, by (S.16) and (S.18), we
have

$\displaystyle\left\|n^{-1}\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\gamma}}}-n^{-1}\sum_{i=1}^{n}[{\bf
B}(y_{i})-E\\{{\bf B}(Y)|{\bf x}_{i}\\}]\right\|_{2}=o_{p}(n^{-1/2}h^{1/2}).$

Hence ${\bf r}_{1}$ satisfies

$\displaystyle\|{\bf r}_{1}\|_{2}$ $\displaystyle\leq$ $\displaystyle
o_{p}(h^{-1})n^{-1/2}h^{1/2}+h^{-1}o_{p}(n^{-1/2}h^{1/2})$ $\displaystyle=$
$\displaystyle o_{p}\\{(nh)^{-1/2}\\}$

by (S.17) and (S.24). Then by Condition (C5), this leads to

$\displaystyle\sup_{y\in[0,1]}|\widehat{c}(y,{\bm{\beta}}_{0})-c(y)|$
$\displaystyle\leq$ $\displaystyle\sup_{y\in[0,1]}|{\bf B}(y)^{\rm
T}\\{\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\\}|+\sup_{y\in[0,1]}|{\bf
B}(y)^{\rm T}{\bm{\gamma}}_{0}-c(y)|$ $\displaystyle\leq$
$\displaystyle\sup_{y\in[0,1]}\|{\bf
B}(y)\|_{2}\left(\|{\bf\Sigma}_{22}^{-1}\|_{2}\left\|n^{-1}\sum_{i=1}^{n}[{\bf
B}(y_{i})-E\\{{\bf B}(Y)|{\bf x}_{i}\\}]\right\|_{2}+\|{\bf
r}_{1}\|_{2}\right)+O(h^{q})$ $\displaystyle=$ $\displaystyle
O(1)O(h^{-1})O_{p}(n^{-1/2}h^{1/2})+o_{p}\\{(nh)^{-1/2}\\}+O(h^{q})$
$\displaystyle=$ $\displaystyle O_{p}\\{(nh)^{-1/2}+h^{q}\\}.$

∎

###  S.5 Proof of Proposition 2

First we will show for any $\epsilon>0$, there exists a constant $C>0$ such
that

$\displaystyle\hbox{pr}\left[l\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}>\sup_{\|{\bf
v}\|_{2}=C}l\\{{\bm{\beta}}_{0}+n^{-1/2}{\bf
v},\widehat{\bm{\gamma}}({\bm{\beta}}_{0}+n^{-1/2}{\bf v})\\}\right]\geq
1-3\epsilon$ (S.26)

for a sufficiently large $n$. This implies there exists a local maximizer of
$l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}$, say
$\widetilde{\bm{\beta}}$, such that
$\|\widetilde{\bm{\beta}}-{\bm{\beta}}_{0}\|_{2}=O_{p}(n^{-1/2})$. We will
further show the Hessian of
$l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}$ is negative definite
for any ${\bm{\beta}}$, hence $\widetilde{\bm{\beta}}$ is the global maximizer
of $l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}$, i.e.
$\widehat{\bm{\beta}}=\widetilde{\bm{\beta}}$, and hence
$\|\widehat{\bm{\beta}}-{\bm{\beta}}_{0}\|_{2}=O_{p}(n^{-1/2})$.

By the Taylor expansion,

$\displaystyle l\\{{\bm{\beta}}_{0}+n^{-1/2}{\bf
v},\widehat{\bm{\gamma}}({\bm{\beta}}_{0}+n^{-1/2}{\bf
v})\\}-l\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}$ (S.27)
$\displaystyle=$ $\displaystyle
n^{-1/2}\frac{dl\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}}{d{\bm{\beta}}^{\rm
T}}{\bf v}+\frac{1}{2}n^{-1}{\bf v}^{\rm
T}\frac{d^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{d{\bm{\beta}}d{\bm{\beta}}^{\rm
T}}{\bf v},$

where ${\bm{\beta}}^{*}\equiv\alpha_{1}({\bm{\beta}}_{0}+n^{-1/2}{\bf
v})+(1-\alpha_{1}){\bm{\beta}}_{0}$ for some $\alpha_{1}\in(0,1)$. Noting that
$\partial
l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}/\partial{\bm{\gamma}}={\bf
0}$ for any ${\bm{\beta}}$, we get

$\displaystyle\frac{dl\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}}{d{\bm{\beta}}}$
$\displaystyle=$ $\displaystyle\frac{\partial
l\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}}{\partial{\bm{\beta}}}+\frac{\partial\widehat{\bm{\gamma}}^{\rm
T}({\bm{\beta}}_{0})}{\partial{\bm{\beta}}}\frac{\partial
l\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}}{\partial{\bm{\gamma}}}$
$\displaystyle=$ $\displaystyle\frac{\partial
l\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}}{\partial{\bm{\beta}}}$
$\displaystyle=$ $\displaystyle\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\beta}}}+\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\\{\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\\},$

where
${\bm{\gamma}}^{*}=\alpha_{2}\widehat{\bm{\gamma}}({\bm{\beta}}_{0})+(1-\alpha_{2}){\bm{\gamma}}_{0}$
for some $\alpha_{2}\in(0,1)$. In addition,

$\displaystyle\frac{d^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{d{\bm{\beta}}d{\bm{\beta}}^{\rm
T}}$ $\displaystyle=$ $\displaystyle\frac{d}{d{\bm{\beta}}^{\rm
T}}\left[\frac{\partial
l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\beta}}}\right]\biggr{|}_{{\bm{\beta}}={\bm{\beta}}^{*}}$
$\displaystyle=$
$\displaystyle\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}}+\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\frac{\partial\widehat{\bm{\gamma}}({\bm{\beta}}^{*})}{\partial{\bm{\beta}}^{\rm
T}}.$

Now since $\partial
l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}/\partial{\bm{\gamma}}={\bf
0}$ for any ${\bm{\beta}}$,

$\displaystyle{\bf 0}$ $\displaystyle=$
$\displaystyle\frac{d}{d{\bm{\beta}}^{\rm T}}\left[\frac{\partial
l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}}\right]$
$\displaystyle=$
$\displaystyle\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\beta}}^{\rm
T}}+\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\frac{\partial\widehat{\bm{\gamma}}({\bm{\beta}})}{\partial{\bm{\beta}}^{\rm
T}}.$

Note that

$\displaystyle\frac{\partial^{2}l({\bm{\beta}},{\bm{\gamma}})}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}=-\sum_{i=1}^{n}\hbox{cov}^{*}\\{{\bf B}(Y),{\bf B}(Y)|{\bf
x}_{i},{\bm{\beta}},{\bm{\gamma}}\\}$

is negative definite. Then,

$\displaystyle\frac{\partial\widehat{\bm{\gamma}}({\bm{\beta}})}{\partial{\bm{\beta}}^{\rm
T}}=-\left[\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right]^{-1}\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\beta}}^{\rm
T}}.$ (S.28)

Hence, we can write (S.27) as

$\displaystyle l\\{{\bm{\beta}}_{0}+n^{-1/2}{\bf
v},\widehat{\bm{\gamma}}({\bm{\beta}}_{0}+n^{-1/2}{\bf
v})\\}-l\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}$
$\displaystyle=$ $\displaystyle n^{-1/2}\left[\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\beta}}}+\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\\{\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\\}\right]^{\rm
T}{\bf v}$ $\displaystyle+\frac{1}{2}{\bf v}^{\rm
T}n^{-1}\Bigg{(}\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}}-\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\left[\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right]^{-1}\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\gamma}}\partial{\bm{\beta}}^{\rm
T}}\Bigg{)}{\bf v}.$

In (S.5), we first analyze $\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})/\partial{\bm{\beta}}$. We have

$\displaystyle n^{-1}\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\beta}}}$ $\displaystyle=$
$\displaystyle n^{-1}\sum_{i=1}^{n}{\bf x}_{i}\\{y_{i}-E^{*}(Y|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})\\}$ $\displaystyle=$ $\displaystyle
n^{-1}\sum_{i=1}^{n}{\bf x}_{i}\\{y_{i}-E(Y|{\bf
x}_{i})\\}+n^{-1}\sum_{i=1}^{n}{\bf x}_{i}\\{E(Y|{\bf x}_{i})-E^{*}(Y|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})\\}.$

By Conditions (C1) and (C7), the first term satisfies

$\displaystyle\left\|n^{-1}\sum_{i=1}^{n}{\bf x}_{i}\\{y_{i}-E(Y|{\bf
x}_{i})\\}\right\|_{2}\asymp_{p}n^{-1/2}.$

For the second term, by Lemma 4 we get

$\displaystyle\left\|n^{-1}\sum_{i=1}^{n}{\bf x}_{i}\\{E(Y|{\bf
x}_{i})-E^{*}(Y|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})\\}\right\|_{2}$ $\displaystyle\leq$
$\displaystyle n^{-1}\sum_{i=1}^{n}\left\|{\bf x}_{i}\right\|_{2}\sup_{{\bf
x}\in\mathcal{X}}|E(Y|{\bf x})-E^{*}(Y|{\bf
x},{\bm{\beta}}_{0},{\bm{\gamma}}_{0})|$ $\displaystyle=$ $\displaystyle
O_{p}(h^{q}),$

where the last equality holds because $E(\|{\bf X}\|_{2})$ is bounded by
Condition (C1). Note that $h^{q}=o_{p}(n^{-1/2})$ under Conditions (C3) and
(C4). Therefore, we have

$\displaystyle\left\|\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\beta}}}\right\|_{2}\asymp_{p}n^{1/2},$

and

$\displaystyle\left\|\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\beta}}}-\sum_{i=1}^{n}{\bf
x}_{i}\\{y_{i}-E(Y|{\bf x}_{i})\\}\right\|_{2}=o_{p}(n^{1/2}).$ (S.30)

Moreover, for any $\epsilon>0$, there exists a constant $0<C_{1}<\infty$ such
that

$\displaystyle\hbox{pr}\left\\{\left\|\frac{\partial
l({\bm{\beta}}_{0},{\bm{\gamma}}_{0})}{\partial{\bm{\beta}}}\right\|_{2}\leq
C_{1}n^{1/2}\right\\}\geq 1-\epsilon.$ (S.31)

Now we analyze
$\partial^{2}l({\bm{\beta}},{\bm{\gamma}})/\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}$. Let ${\bf u}_{m}\in U_{m}\equiv\\{{\bf u}\in\mathbb{R}^{m}:\|{\bf
u}\|_{2}=1\\}$, then for fixed ${\bm{\beta}}$ and ${\bm{\gamma}}$ we have

$\displaystyle\left\|-n^{-1}\frac{\partial^{2}l({\bm{\beta}},{\bm{\gamma}})}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}{\bf u}_{m}\right\|_{2}$ $\displaystyle=$
$\displaystyle\left\|n^{-1}\sum_{i=1}^{n}{\bf x}_{i}\hbox{cov}^{*}\\{Y,{\bf
B}(Y)^{\rm T}{\bf u}_{m}|{\bf
x}_{i},{\bm{\beta}},{\bm{\gamma}}\\}\right\|_{2}$ $\displaystyle\leq$
$\displaystyle n^{-1}\sum_{i=1}^{n}\left\|{\bf x}_{i}\right\|_{2}\sup_{{\bf
x}\in\mathcal{X}}|\hbox{cov}^{*}\\{Y,{\bf B}(Y)^{\rm T}{\bf u}_{m}|{\bf
x},{\bm{\beta}},{\bm{\gamma}}\\}|$ $\displaystyle=$ $\displaystyle
O_{p}(h^{1/2})$

by Condition (C1) and Lemma 2. Then using the fact that $U_{m}$ is compact,
taking supremum with respect to ${\bf u}_{m}$ gives

$\displaystyle\left\|\frac{\partial^{2}l({\bm{\beta}},{\bm{\gamma}})}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\right\|_{2}=O_{p}(nh^{1/2}).$ (S.32)

Further recall that
$\|\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\|_{2}=O_{p}\\{(nh)^{-1/2}\\}$
by Proposition 1. This implies

$\displaystyle\left\|\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\\{\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\\}\right\|_{2}=O_{p}(n^{1/2}),$

and equivalently, there exists a constant $0<C_{2}<\infty$ such that

$\displaystyle\hbox{pr}\left\\{\left\|\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\\{\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\\}\right\|_{2}\leq
C_{2}n^{1/2}\right\\}\geq 1-\epsilon.$ (S.33)

We next show
$-n^{-1}\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})/\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}$ converges to ${\bf\Sigma}_{12}$. Note that
$\|{\bm{\gamma}}^{*}-{\bm{\gamma}}_{0}\|_{2}=O_{p}\\{(nh)^{-1/2}\\}=o_{p}(1)$.
Then by Condition (C1), Lemma 1 and 5, we get

$\displaystyle\left\|\left(-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}-{\bf\Sigma}_{12}\right){\bf u}_{m}\right\|_{2}$ (S.34) $\displaystyle\leq$
$\displaystyle\left\|n^{-1}\sum_{i=1}^{n}{\bf x}_{i}[\hbox{cov}^{*}\\{Y,{\bf
B}(Y)^{\rm T}{\bf u}_{m}|{\bf
x}_{i},{\bm{\beta}}_{0},{\bm{\gamma}}^{*}\\}-\hbox{cov}\\{Y,{\bf B}(Y)^{\rm
T}{\bf u}_{m}|{\bf x}_{i}\\}]\right\|_{2}$
$\displaystyle+\left\|n^{-1}\sum_{i=1}^{n}{\bf x}_{i}\hbox{cov}\\{Y,{\bf
B}(Y)^{\rm T}{\bf u}_{m}|{\bf x}_{i}\\}-E[{\bf X}\hbox{cov}\\{Y,{\bf
B}(Y)^{\rm T}{\bf u}_{m}|{\bf X}\\}]\right\|_{2}$ $\displaystyle\leq$
$\displaystyle\left(n^{-1}\sum_{i=1}^{n}\|{\bf
x}_{i}\|_{2}\right)o_{p}(h^{1/2})+o_{p}(n^{-1/2})$ $\displaystyle=$
$\displaystyle o_{p}(h^{1/2}).$

The last equality holds because $n^{-1/2}=o(h^{1/2})$ by Conditions (C3) and
(C4). This implies

$\displaystyle\left\|-n^{-1}\frac{\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}-{\bf\Sigma}_{12}\right\|_{2}=o_{p}(h^{1/2}).$ (S.35)

Now we show
$-n^{-1}\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}/\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}$ converges to ${\bf\Sigma}_{11}$. First we will show that
$\|\widehat{\bm{\gamma}}({\bm{\beta}}^{*})-{\bm{\gamma}}_{0}\|_{2}=O_{p}\\{(nh)^{-1/2}\\}$.
Note that since

$\displaystyle\begin{bmatrix}\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}/\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}&\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}/\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}\\\
\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}/\partial{\bm{\gamma}}\partial{\bm{\beta}}^{\rm
T}&\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}/\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}\end{bmatrix}$ $\displaystyle=$
$\displaystyle-\sum_{i=1}^{n}\begin{bmatrix}{\bf x}_{i}{\bf x}_{i}^{\rm
T}\hbox{var}^{*}\\{Y|{\bf
x}_{i},{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}&{\bf
x}_{i}\hbox{cov}^{*}\\{Y,{\bf B}(Y)|{\bf
x}_{i},{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}\\\ [{\bf
x}_{i}\hbox{cov}^{*}\\{Y,{\bf B}(Y)|{\bf
x}_{i},{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}]^{\rm
T}&\hbox{cov}^{*}\\{{\bf B}(Y),{\bf B}(Y)|{\bf
x}_{i},{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}\end{bmatrix}$
$\displaystyle=$
$\displaystyle-\sum_{i=1}^{n}\hbox{cov}^{*}\left\\{\begin{bmatrix}{\bf
x}_{i}Y\\\ {\bf B}(Y)\end{bmatrix},\begin{bmatrix}{\bf x}_{i}Y\\\ {\bf
B}(Y)\end{bmatrix}\Bigg{|}{\bf
x}_{i},{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\right\\}$

is negative definite,
$\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}/\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}$ and
$\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}/\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}$ are negative definite, and the Schur complement of
$\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}/\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}$,

$\displaystyle\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}}-\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\left[\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right]^{-1}\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\beta}}^{\rm
T}},$

is also negative definite. Let ${\bf u}_{p}\in U_{p}\equiv\\{{\bf
u}\in\mathbb{R}^{p}:\|{\bf u}\|_{2}=1\\}$, then we have

$\displaystyle 0\leq{\bf u}_{p}^{\rm
T}\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\left[-\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right]^{-1}\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\beta}}^{\rm
T}}{\bf u}_{p}<{\bf u}_{p}^{\rm
T}\left[-\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}}\right]{\bf u}_{p}.$

It is easy to see that

$\displaystyle{\bf u}_{p}^{\rm
T}\left[-\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}}\right]{\bf u}_{p}=\sum_{i=1}^{n}({\bf x}_{i}^{\rm T}{\bf
u}_{p})^{2}\hbox{var}^{*}\\{Y|{\bf
x}_{i},{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}\asymp_{p}n,$

and

$\displaystyle\left\|\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right\|_{2}$ $\displaystyle=$ $\displaystyle\sup_{{\bf u}_{m}\in
U_{m}}{\bf u}_{m}^{\rm
T}\left[-\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right]{\bf u}_{m}$ $\displaystyle=$ $\displaystyle\sup_{{\bf u}_{m}\in
U_{m}}\sum_{i=1}^{n}\hbox{var}^{*}\\{{\bf B}(Y)^{\rm T}{\bf u}_{m}|{\bf
x}_{i},{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}$
$\displaystyle\asymp_{p}$ $\displaystyle nh$

by Lemma 3. Hence we have

$\displaystyle\left\|\left[\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right]^{-1}\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\beta}}^{\rm
T}}\right\|_{2}=O_{p}(h^{-1/2}).$ (S.36)

Now, (S.28) and the Taylor expansion of
$\widehat{\bm{\gamma}}({\bm{\beta}}^{*})$ with ${\bm{\beta}}^{**}$ on the line
connecting ${\bm{\beta}}_{0}$ and ${\bm{\beta}}^{*}$ gives

$\displaystyle\widehat{\bm{\gamma}}({\bm{\beta}}^{*})-{\bm{\gamma}}_{0}$
$\displaystyle=$
$\displaystyle\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}-\left[\frac{\partial^{2}l\\{{\bm{\beta}}^{**},\widehat{\bm{\gamma}}({\bm{\beta}}^{**})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right]^{-1}\frac{\partial^{2}l\\{{\bm{\beta}}^{**},\widehat{\bm{\gamma}}({\bm{\beta}}^{**})\\}}{\partial{\bm{\gamma}}\partial{\bm{\beta}}^{\rm
T}}({\bm{\beta}}^{*}-{\bm{\beta}}_{0}).$

Recall that
$\|\widehat{\bm{\gamma}}({\bm{\beta}}_{0})-{\bm{\gamma}}_{0}\|_{2}=O_{p}\\{(nh)^{-1/2}\\}$
by Proposition 1 and $\|{\bm{\beta}}^{*}-{\bm{\beta}}_{0}\|_{2}=O(n^{-1/2})$
by construction. These imply
$\|\widehat{\bm{\gamma}}({\bm{\beta}}^{*})-{\bm{\gamma}}_{0}\|_{2}=O_{p}\\{(nh)^{-1/2}\\}$.
Therefore, we get

$\displaystyle\left|{\bf u}_{p}^{\rm
T}\left[-n^{-1}\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}}-{\bf\Sigma}_{11}\right]{\bf u}_{p}\right|$ $\displaystyle\leq$
$\displaystyle\left|n^{-1}\sum_{i=1}^{n}({\bf x}_{i}^{\rm T}{\bf
u}_{p})^{2}[\hbox{var}^{*}\\{Y|{\bf
x}_{i},{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}-\hbox{var}(Y|{\bf
x}_{i})]\right|$ $\displaystyle+\left|n^{-1}\sum_{i=1}^{n}({\bf x}_{i}^{\rm
T}{\bf u}_{p})^{2}\hbox{var}(Y|{\bf x}_{i})-E\\{({\bf X}^{\rm T}{\bf
u}_{p})^{2}\hbox{var}(Y|{\bf X})\\}\right|$ $\displaystyle=$ $\displaystyle
o_{p}(1)$

by Condition (C1) and Lemma 5, which implies

$\displaystyle\left\|-n^{-1}\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}}-{\bf\Sigma}_{11}\right\|_{2}=o_{p}(1).$ (S.37)

Below, we show
$-n^{-1}d^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}/d{\bm{\beta}}d{\bm{\beta}}^{\rm
T}$ converges to ${\bf\Sigma}^{*}$. Using similar argument to (S.34), we get

$\displaystyle\left\|-n^{-1}\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}-{\bf\Sigma}_{12}\right\|_{2}=o_{p}(h^{1/2}).$ (S.38)

In addition, by the similar proof to showing
$-n^{-1}\partial^{2}l({\bm{\beta}}_{0},{\bm{\gamma}}^{*})/\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}$ converges to ${\bf\Sigma}_{22}^{-1}$ in Proposition 1, we get

$\displaystyle\left\|\left[-n^{-1}\frac{\partial^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right]^{-1}-{\bf\Sigma}_{22}^{-1}\right\|_{2}=o_{p}(h^{-1}).$ (S.39)

Now, note that

$\displaystyle\frac{d^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{d{\bm{\beta}}d{\bm{\beta}}^{\rm
T}}=\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\beta}}\partial{\bm{\beta}}^{\rm
T}}-\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\beta}}\partial{\bm{\gamma}}^{\rm
T}}\left[\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\gamma}}^{\rm
T}}\right]^{-1}\frac{\partial^{2}l\\{{\bm{\beta}},\widehat{\bm{\gamma}}({\bm{\beta}})\\}}{\partial{\bm{\gamma}}\partial{\bm{\beta}}^{\rm
T}}.$

Combining (S.37), (S.38), and (S.39), we have

$\displaystyle\left\|-n^{-1}\frac{d^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{d{\bm{\beta}}d{\bm{\beta}}^{\rm
T}}-{\bf\Sigma}^{*}\right\|_{2}=o_{p}(1),$ (S.40)

where ${\bf\Sigma}^{*}$ is given in Condition (C7). This shows the Hessian at
${\bm{\beta}}^{*}$ is negative definite and its 2-norm is of order $n$ with
probability approaching 1. Hence, there exists a constant $0<C_{3}<\infty$
such that

$\displaystyle\hbox{pr}\left\\{{\bf v}^{\rm
T}\frac{d^{2}l\\{{\bm{\beta}}^{*},\widehat{\bm{\gamma}}({\bm{\beta}}^{*})\\}}{d{\bm{\beta}}d{\bm{\beta}}^{\rm
T}}{\bf v}\leq-C^{2}C_{3}n\right\\}\geq 1-\epsilon.$ (S.41)

Combining (S.5), (S.31), (S.33), and (S.41), we have

$\displaystyle l\\{{\bm{\beta}}_{0}+n^{-1/2}{\bf
v},\widehat{\bm{\gamma}}({\bm{\beta}}_{0}+n^{-1/2}{\bf
v})\\}-l\\{{\bm{\beta}}_{0},\widehat{\bm{\gamma}}({\bm{\beta}}_{0})\\}$