EconBase
← Back to paper

Fractional order statistic approximation for nonparametric conditional quantile inference

The exact contents of citations.db main_text.text for this paper — one flattened LaTeX string, title through conclusion, appendix excluded, unmodified except for removing email addresses. This is what our citation measures are computed over.

88,937 characters

Supplemental appendix: Fractional order statistic approximation for nonparametric conditional quantile inference


\maketitle


\begin{abstract}
\iftoggle{SUPPLEMENTAL}{This supplement includes longer proofs and additional details.
}{
Using and extending fractional order statistic theory, we characterize the $O(n^{-1})$ coverage probability error of the previously proposed confidence intervals for population quantiles using $L$-statistics as endpoints in \citet{Hutson1999}.  We derive an analytic expression for the $n^{-1}$ term, which may be used to calibrate the nominal coverage level to get $O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright)$ coverage error.  Asymptotic power is shown to be optimal.  Using kernel smoothing, we propose a related method for nonparametric inference on conditional quantiles.  This new method compares favorably with asymptotic normality and bootstrap methods in theory and in simulations.  Code is provided for both unconditional and conditional inference.

{\sc JEL classification}:
                          C21

{\sc Keywords}: Dirichlet, high-order accuracy, inference-optimal bandwidth, kernel smoothing.
}
\end{abstract}

\vfill

\begin{flushright}
\copyright\ 2016 by the authors.
This manuscript version is made available under the CC-BY-NC-ND 4.0 license: \url{http://creativecommons.org/licenses/by-nc-nd/4.0/}
\end{flushright}






\doublespacing

\iftoggle{SUPPLEMENTAL}{\vfill\pagebreak}{\newpage}

\section{Introduction\label{sec:intro}}

Quantiles contain information about a distribution's shape.  Complementing the mean, they capture heterogeneity, inequality, and other measures of economic interest.  Nonparametric conditional quantile models further allow arbitrary heterogeneity across regressor values.
This paper concerns nonparametric inference on quantiles and conditional quantiles.
In particular, we characterize the high-order accuracy of both \citeposs{Hutson1999} $L$-statistic-based confidence intervals (CIs) and our new conditional quantile CIs.

Conditional quantiles appear across diverse topics because they are fundamental statistical objects.
Such topics include wages \citep{Hogg1975,Chamberlain1994,Buchinsky1994},
infant birthweight \citep{Abrevaya2001}, demand for alcohol \citep{ManningEtAl1995},
and Engel curves \citep[][pp.\ 81--82]{AlanEtAl2005,Deaton1997}, which we examine in our empirical application.


We formally derive the coverage probability error (CPE) of the CIs from \citet{Hutson1999}, as well as asymptotic power of the corresponding hypothesis tests.
\citet{Hutson1999} had proposed CIs for quantiles using $L$-statistics (interpolating between order statistics) as endpoints and found they performed well, but formal proofs were lacking.
Using the analytic $n^{-1}$ term we derive in the CPE, we provide a new calibration to achieve $O\bigl(n^{-3/2}[\log(n)]^3\bigr)$ CPE, analogous to the \citet{HoLee2005a} analytic calibration of the CIs in \citet{BeranHall1993}.

The theoretical results we develop contribute to the fractional order statistic literature and provide the basis for inference on other objects of interest explored in \citet{GoldmanKaplan2014b} and \citet{Kaplan2014}.  In particular, Theorem \ref{thm:cdferror} tightly links the distributions of $L$-statistics from the observed and `ideal' (unobserved) fractional order statistic processes.
Additionally, Lemma \ref{lem:den-i} provides Dirichlet PDF and PDF derivative approximations.

High-order accuracy is important for small samples (e.g.,\ for experiments) as well as nonparametric analysis with small \emph{local} sample sizes.  For example, if $n=1024$ and there are five binary regressors, then the smallest local sample size cannot exceed $1024/2^5=32$.


For nonparametric conditional quantile inference, we apply the unconditional method to a local sample (similar to local constant kernel regression), smoothing over continuous covariates and also allowing discrete covariates.
CPE is minimized by balancing the CPE of our unconditional method and the CPE from bias due to smoothing.
We derive the optimal CPE and bandwidth rates, as well as a plug-in bandwidth when there is a single continuous covariate.

Our $L$-statistic method has theoretical and computational advantages over methods based on normality or an unsmoothed bootstrap.
The theoretical bottleneck for our approach is the need to use a uniform kernel.
Nonetheless, even if normality or bootstrap methods assume an infinitely differentiable conditional quantile function (and hypothetically fit an infinite-degree local polynomial), our CPE is still of smaller order with one or two continuous covariates.
Our method also computes more quickly than existing methods (of reasonable accuracy), handling even more challenging tasks in 10--15 seconds instead of minutes.

Recent complementary work of \citet{FanLiu2016} also concerns a ``direct method'' of nonparametric inference on conditional quantiles.  They use a limiting Gaussian process to derive first-order accuracy in a general setting, whereas we use the finite-sample Dirichlet process to achieve high-order accuracy in an iid setting.
\citet{FanLiu2016} also provide uniform (over $X$) confidence bands.
We suggest a confidence band from interpolating a growing number of joint CIs (as in \citet{HorowitzLee2012}), although it will take additional work to rigorously justify.  A different, ad hoc confidence band described in Section \ref{sec:sim} generally outperformed others in our simulations.

If applied to a local constant estimator with a uniform kernel and the same bandwidth, the \citet{FanLiu2016} approach is less accurate than ours due to the normal (instead of beta) reference distribution and integer (instead of interpolated) order statistics in their CI in equation (6).
However, with other estimators like local polynomials or that in \citet{DonaldEtAl2012}, the \citet{FanLiu2016} method is not necessarily less accurate.  One limitation of our approach is that it cannot incorporate these other estimators, whereas Assumption GI(iii) in \citet{FanLiu2016} includes any estimator that weakly converges (over a range of quantiles) to a Gaussian process with a particular structure.  We compare further in our simulations.  One open question is whether using our beta reference and interpolation can improve accuracy for the general \citet{FanLiu2016} method beyond the local constant estimator with a uniform kernel; our Lemma \ref{lem:u1u2-approx} shows this at least retains first-order accuracy.




The order statistic approach to quantile inference uses the idea of the probability integral transform, which dates back to R.\ A.\ \citet{Fisher1932}, Karl \citet{Pearson1933}, and \citet{Neyman1937}.
For continuous $X_i\stackrel{iid}{\sim}F(\cdot)$, $F(X_i)\stackrel{iid}{\sim}\textrm{Unif}(0,1)$.
Each order statistic from such an iid uniform sample has a known beta distribution for any sample size $n$.
We show that the $L$-statistic linearly interpolating consecutive order statistics also follows an approximate beta distribution, with only $O(n^{-1})$ error in CDF.
Although $O(n^{-1})$ is an asymptotic claim, the CPE of the CI using the $L$-statistic endpoint is bounded between the CPEs of the CIs using the two order statistics comprising the $L$-statistic, where one such CPE is too small and one is too big, for any sample size.
This is an advantage over methods more sensitive to asymptotic approximation error.





Many other approaches to one-sample quantile inference have been explored.
With Edgeworth expansions, \citet{HallSheather1988} and \citet{Kaplan2015} obtain two-sided $O(n^{-2/3})$ CPE.
With bootstrap, smoothing is necessary for high-order accuracy.  This increases the computational burden and requires good bandwidth selection in practice.\footnote{For example, while achieving the impressive two-sided CPE of $O(n^{-3/2})$, \citet[][p.\ 833]{PolanskySchucany1997} admit, ``If this method is to be of any practical value, a better bandwidth estimation technique will certainly be required.''} See \citet[\S1]{HoLee2005b} for a review of bootstrap methods.
Smoothed empirical likelihood \citep{ChenHall1993} also achieves nice theoretical properties, but with the same caveats.

Other order statistic-based CIs dating back to \citet{Thompson1936} are surveyed in \citet[\S7.1]{DavidNagaraja2003}.  Most closely related to \citet{Hutson1999} is \citet{BeranHall1993}.  Like \citet{Hutson1999}, \citet{BeranHall1993} linearly interpolate order statistics for CI endpoints, but with an interpolation weight based on the binomial distribution. Although their proofs use expansions of the \citet{Renyi1953} representation instead of fractional order statistic theory, their $n^{-1}$ CPE term is identical to that for \citet{Hutson1999} other than the different weight.
Prior work \citep[e.g.,][]{Bickel1967,Shorack1972} has established asymptotic normality of $L$-statistics and convergence of the sample quantile process to a Gaussian limit process, but without such high-order accuracy.

The most apparent difference between the two-sided CIs of \citet{BeranHall1993} and \citet{Hutson1999} is that the former are symmetric in the order statistic index, whereas the latter are equal-tailed.  This allows \citet{Hutson1999} to be computed further into the tails.
Additionally, our framework can be extended to CIs for interquantile ranges and two-sample quantile differences \citep{GoldmanKaplan2014b},
which has not been done in the R\'enyi representation framework.

For nonparametric conditional quantile inference, in addition to the aforementioned \citet{FanLiu2016} approach, \citet{Chaudhuri1991} derives the pointwise asymptotic normal distribution of a local polynomial estimator.
\citet{QuYoon2015} propose modified local linear estimators of the conditional quantile process that converge weakly to a Gaussian process, and they suggest using a type of bias correction that strictly enlarges a CI to deal with the first-order effect of asymptotic bias when using the MSE-optimal bandwidth rate.




Section \ref{sec:cdf-err} contains our theoretical results on fractional order statistic approximation, which are applied to unconditional quantile inference in Section \ref{sec:inf-unconditional}.
Section \ref{sec:inf-conditional} concerns our new conditional quantile inference method.
An empirical application and simulation results are in Sections \ref{sec:empirical} and \ref{sec:sim}, respectively.
Proof sketches are collected in Appendix \ref{sec:app-pfs}, while the supplemental appendix contains full proofs.
The supplemental appendix also contains details of the plug-in bandwidth calculations, as well as additional empirical and simulation results.

Notationally,
$\phi(\cdot)$ and $\Phi(\cdot)$ are respectively the standard normal PDF and CDF,
$\doteq$ should be read as ``is equal to, up to smaller-order terms'',
$\asymp$ as ``has exact (asymptotic) rate/order of'',
and $A_n=O(B_n)$ as usual.
Acronyms used are those for cumulative distribution function (CDF), confidence interval (CI), coverage probability (CP), coverage probability error (CPE), and probability density function (PDF).

















\section{Fractional order statistic theory}\label{sec:cdf-err}

In this section, we introduce notation and present our core theoretical results linking unobserved `ideal' fractional $L$-statistics with their observed counterparts.

Given an iid sample $\{X_i\}_{i=1}^n$ of draws from a continuous CDF denoted\footnote{$F$ will often be used with a random variable subscript to denote the CDF of that particular random variable.  If no subscript is present, then $F(\cdot)$ refers to the CDF of $X$.  Similarly for the PDF $f(\cdot)$.} $F(\cdot)$, interest is in $Q(p) \equiv F^{-1}(p)$ for some $p\in(0,1)$, where $Q(\cdot)$ is the quantile function.
For $u\in(0,1)$, the sample $L$-statistic commonly associated with $Q(u)$ is
\begin{equation}
\label{eqn:QXLdef}
\hat Q^L_X(u)
  \equiv (1-\epsilon) X_{n:k}
         +\epsilon X_{n:k+1} , \quad
k \equiv \lfloor u(n+1)\rfloor, \quad
\epsilon \equiv u(n+1) - k ,
\end{equation}
where $\lfloor\cdot\rfloor$ is the floor function, $\epsilon$ is the interpolation weight, and $X_{n:k}$ denotes the $k$th order statistic (i.e.,\ $k$th smallest sample value).  While $Q(u)$ is latent and nonrandom, $\hat Q^L_X(u)$ is a random variable, and $\hat Q^L_X(\cdot)$ is a stochastic process, observed for arguments in $[1/(n+1),n/(n+1)]$.

Let $\Xi_n\equiv\mathopen{}\mathclose\bgroup\originalleft\{k/(n+1)\aftergroup\egroup\originalright\}_{k=1}^n$ denote the set of quantiles corresponding to the observed order statistics.  If $u\in\Xi_n$, then no interpolation is necessary and $\hat Q^L_X(u)=X_{n:k}$.  As detailed in Section \ref{sec:inf-unconditional}, application of the probability integral transform yields exact coverage probability of a CI endpoint $X_{n:k}$ for $Q(p)$:
$ P\mathopen{}\mathclose\bgroup\originalleft(X_{n:k}<F^{-1}(p)\aftergroup\egroup\originalright) = P\mathopen{}\mathclose\bgroup\originalleft( U_{n:k}<p\aftergroup\egroup\originalright) $,
where $U_{n:k}\equiv F(X_{n:k})\sim\beta(k,n+1-k)$ is equal in distribution to the $k$th order statistic from $U_i\stackrel{iid}{\sim}\textrm{Unif}(0,1)$, $i=1,\ldots,n$ \citep[][\textbf{8.7.4}]{Wilks1962}.
However, we also care about $u\notin\Xi_n$, in which case $k$ is fractional.
To better handle such fractional order statistics, we will present a tight link between the marginal distributions of the stochastic process $\hat Q^L_X(\cdot)$ and those of the analogous `ideal' (I) process
\begin{align}
\label{eqn:Qxidef}
\tilde Q^I_X(\cdot) &\equiv F^{-1}\mathopen{}\mathclose\bgroup\originalleft(\tilde Q^I_U(\cdot)\aftergroup\egroup\originalright) ,
\end{align}
where $\tilde Q^I_U(\cdot)$ is the ideal (I) uniform (U) fractional order ``statistic'' process.  We use a tilde in $\tilde Q^I_X(\cdot)$ and $\tilde Q^I_U(\cdot)$ instead of the hat like in $\hat Q^L_X(\cdot)$ to emphasize that the former are unobserved (hence not true statistics), whereas the latter is computable from the sample data.

This $\tilde Q^I_U(\cdot)$ in \eqref{eqn:Qxidef} is a Dirichlet process \citep{Ferguson1973,Stigler1977} on the unit interval with index measure $\nu\mathopen{}\mathclose\bgroup\originalleft([0,t]\aftergroup\egroup\originalright)=(n+1)t$.  Its univariate marginals are
\begin{align}\label{eqn:QUIdist}
\tilde Q^I_U(u) = U_{n:(n+1)u} \sim \beta\bigl((n+1)u,(n+1)(1-u)\bigr) .
\end{align}
The marginal distribution of $\mathopen{}\mathclose\bgroup\originalleft(\tilde Q^I_U(u_1),\tilde Q^I_U(u_2)-\tilde Q^I_U(u_1),\ldots,\tilde Q^I_U(u_k)-\tilde Q^I_U(u_{k-1})\aftergroup\egroup\originalright)$ for $u_1<\cdots<u_k$ is Dirichlet with parameters $\mathopen{}\mathclose\bgroup\originalleft(u_1(n+1),(u_2-u_1)(n+1),\ldots,(u_k-u_{k-1})(n+1)\aftergroup\egroup\originalright)$.

For all $u\in\Xi_n$, $\tilde Q^I_X(u)$ coincides with $\hat Q^L_X(u)$; they differ only in their interpolation between these points.  Proposition \ref{prop:error-prob} shows $\tilde Q^I_X(\cdot)$ and $\hat Q^L_X(\cdot)$ to be closely linked in probability.
\begin{proposition}\label{prop:error-prob}
For any fixed $\delta>0$ and $m>0$, define
 $\mathcal{U}^\delta \equiv \{u\in(0,1) \mid \forall\,t\in(u-m,u+m), f\mathopen{}\mathclose\bgroup\originalleft(F^{-1}(t)\aftergroup\egroup\originalright) \ge \delta \}$ and
 $\mathcal{U}^\delta_n \equiv \mathcal{U}^{\delta} \cap [\frac{1}{n+1},\frac{n}{n+1}]$; then, $\underset{u \in \mathcal{U}^\delta_n}{\sup}\mathopen{}\mathclose\bgroup\originalleft| \tilde Q^I_X(u) - \hat Q^L_X(u) \aftergroup\egroup\originalright| = O_p\mathopen{}\mathclose\bgroup\originalleft(n^{-1}\log(n)\aftergroup\egroup\originalright)$.
\end{proposition}
Although Proposition \ref{prop:error-prob} motivates approximating the distribution of $\hat Q^L_X(u)$ by that of $\tilde Q^I_X(u)$, it is not relevant to high-order accuracy.  In fact, its result is achieved by any interpolation between $X_{n:k}$ and $X_{n:k+1}$, not just $\hat Q^L_X(u)$; in contrast, the high-order accuracy we establish in Theorem \ref{thm:IDEAL-single} is only possible with precise interpolations like $\hat Q^L_X(u)$.

Next, we consider marginal distributions of fixed dimension $J$.
We also consider the Gaussian approximation to the sampling distribution of fractional order statistics.
It is well known that the centered and scaled empirical process for standard uniform random variables converges to a Brownian bridge.  For standard Brownian bridge process $B(\cdot)$, we index by $u\in(0,1)$ the additional stochastic processes
\begin{align*}
\tilde Q^B_U(u) \equiv u + n^{-1/2} B(u) \qquad\textrm{and} \qquad \tilde Q^B_X(u) \equiv F^{-1}\bigl(\tilde Q^B_U(u)\bigr).
\end{align*}
The vector $\tilde Q^I_U(\mathbf u)$ has an ordered Dirichlet distribution (i.e.,\ the spacings between consecutive $\tilde Q^I_U(u_j)$ follow a joint Dirichlet distribution), while $\tilde Q^B_U(\mathbf u)$ is multivariate Gaussian.
Lemma \ref{lem:den-i} in the appendix shows the close relationship between multivariate Dirichlet and Gaussian PDFs and PDF derivatives.

Theorem \ref{thm:cdferror} shows the close distributional link among linear combinations of ideal, interpolated, and Gaussian-approximated fractional order statistics.  Specifically, for arbitrary weight vector $\boldsymbol{\psi} \in \mathbb{R}^J$, we (distributionally) approximate
\begin{align}
\label{eqn:def-L}
L^L \equiv \sum_{j=1}^J \psi_j \hat Q^L_X(u_j)
\quad&\textrm{by}\quad
L^I \equiv \sum_{j=1}^J \psi_j \tilde Q^I_X(u_j)
,\\
\notag
\quad\textrm{or alternatively } &\textrm{by}\quad
L^B \equiv \sum_{j=1}^J \psi_j \tilde Q^B_X(u_j).
\end{align}

Our assumptions for this section are now presented, followed by the main theoretical result.  Assumption \ref{a:hut-pf} ensures that the first three derivatives of the quantile function are uniformly bounded in neighborhoods of the quantiles, $u_j$, which helps bound remainder terms in the proofs.  We use \textbf{bold} for vectors and \underline{underline} for matrices.
\begin{assumption}\label{a:iid}
Sampling is iid: $X_i\stackrel{iid}{\sim} F$, $i=1,\ldots,n$.
\end{assumption}
\begin{assumption}\label{a:hut-pf}
For each quantile $u_j$, the PDF $f(\cdot)$ (corresponding to CDF $F(\cdot)$ in \ref{a:iid}) satisfies
(i)  $f(F^{-1}(u_j))>0$;
(ii) $f''(\cdot)$ is continuous in some neighborhood of $F^{-1}(u_j)$, i.e.,\ $f\in C^2\mathopen{}\mathclose\bgroup\originalleft(U_\delta\mathopen{}\mathclose\bgroup\originalleft(F^{-1}(u_j)\aftergroup\egroup\originalright)\aftergroup\egroup\originalright)$ with $U_\delta(x)$ denoting some $\delta$-neighborhood of point $x\in\mathbb R$.
\end{assumption}
\begin{theorem}\label{thm:cdferror}
Define $\underline{\mathcal{V}}$ as the $J\times J$ matrix with row $i$, column $j$ entries
$\underline{\mathcal{V}}_{i,j} = \min\{u_i,u_j\} - u_i u_j$.
Let $\underline{\mathcal{A}}$ be the $J\times J$ matrix with main diagonal entries $\underline{\mathcal{A}}_{j,j}=f\mathopen{}\mathclose\bgroup\originalleft(F^{-1}(u_j)\aftergroup\egroup\originalright)$ and zeros elsewhere, and let
\begin{gather*}
\mathcal{V}_{\boldsymbol{\psi}} \equiv \boldsymbol{\psi}'\mathopen{}\mathclose\bgroup\originalleft(\underline{\mathcal{A}}^{-1}\underline{\mathcal{V}}\,\underline{\mathcal{A}}^{-1}\aftergroup\egroup\originalright)\boldsymbol{\psi}, \quad
\mathbb{X}_0  \equiv \sum_{j=1}^J \psi_j F^{-1}(u_j).
\end{gather*}
Let Assumption \ref{a:iid} hold, and let \ref{a:hut-pf} hold at $\mathbf{\bar u}$.
Given the definitions in \eqref{eqn:QXLdef}, \eqref{eqn:Qxidef}, and \eqref{eqn:def-L}, the following results hold uniformly over $\mathbf u=\mathbf{\bar u}+o(1)$.
\begin{enumerate}\item \label{thm:cdferror-ptwise}For a given constant $K$,
\begin{align*}
P & \bigg(L^L<\mathbb{X}_{0} + n^{-1/2}K\bigg)
   -P\bigg(L^I <\mathbb{X}_{0} + n^{-1/2}K\bigg)  \\
  &= \frac{K\exp\mathopen{}\mathclose\bgroup\originalleft\{-K^2/(2 \mathcal{V}_{\boldsymbol{\psi}})\aftergroup\egroup\originalright\}}{\sqrt{2\pi \mathcal{V}_{\boldsymbol{\psi}}^3}}
     \mathopen{}\mathclose\bgroup\originalleft[\sum_{j=1}^J \mathopen{}\mathclose\bgroup\originalleft(\frac{\psi_j^2 \epsilon_j(1-\epsilon_j) }{f\mathopen{}\mathclose\bgroup\originalleft[F^{-1}(u_j)\aftergroup\egroup\originalright]^2}\aftergroup\egroup\originalright)\aftergroup\egroup\originalright] n^{-1}
    +O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright) ,
\end{align*}
where the remainder is uniform over all $K$.
\item \label{thm:cdferror-unif} Uniformly over $K$,
\begin{align*}
\sup_{K\in\mathbb R}
  &\mathopen{}\mathclose\bgroup\originalleft[P\bigg(L^L<\mathbb{X}_{0} + n^{-1/2}K\bigg)
        -P\bigg(L^I <\mathbb{X}_{0} + n^{-1/2}K\bigg)\aftergroup\egroup\originalright] \\
  &= \frac{e^{-1/2}}{\sqrt{2\pi \mathcal{V}_{\boldsymbol{\psi}}^2}}
     \mathopen{}\mathclose\bgroup\originalleft[\sum_{j=1}^J \mathopen{}\mathclose\bgroup\originalleft(\frac{\psi_j^2 \epsilon_j(1-\epsilon_j) }{f\mathopen{}\mathclose\bgroup\originalleft[F^{-1}(u_j)\aftergroup\egroup\originalright]^2}\aftergroup\egroup\originalright)\aftergroup\egroup\originalright] n^{-1}
    +O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright) , \\
    \sup_{K\in\mathbb R}
      &\mathopen{}\mathclose\bgroup\originalleft|P\bigg(L^L<\mathbb{X}_{0} + n^{-1/2}K\bigg)
            -P\bigg(L^B <\mathbb{X}_{0} + n^{-1/2}K\bigg)\aftergroup\egroup\originalright|
      = O\mathopen{}\mathclose\bgroup\originalleft(n^{-1/2}[\log(n)]^3\aftergroup\egroup\originalright) .
\end{align*}
\end{enumerate}
\end{theorem}






















\section{Quantile inference: unconditional}\label{sec:inf-unconditional}

For inference on $Q(p)$, we continue to maintain \ref{a:iid} and \ref{a:hut-pf}.
For $p\in(0,1)$ and confidence level $1-\alpha$, define $u^h(\alpha)$ and $u^l(\alpha)$ to solve
\begin{align}
\label{eqn:uhdef}
\alpha &= P\Bigl(\tilde Q^I_U\mathopen{}\mathclose\bgroup\originalleft(u^h(\alpha)\aftergroup\egroup\originalright)<p\Bigr)  ,
\quad
\alpha = P\Bigl(\tilde Q^I_U\mathopen{}\mathclose\bgroup\originalleft(u^l(\alpha)\aftergroup\egroup\originalright)>p\Bigr)  ,
\end{align}
with $\tilde Q^I_U(u)\sim\beta\bigl((n+1)u,(n+1)(1-u)\bigr)$ from \eqref{eqn:QUIdist}, parallel to (7) and (8) in \citet{Hutson1999}.

One-sided CI endpoints for $Q(p)$ are $\hat Q^L_X(u^h)$ or $\hat Q^L_X(u^l)$.
Two-sided CIs replace $\alpha$ with $\alpha/2$ in \eqref{eqn:uhdef} and use both endpoints.  This use of $\alpha/2$ yields the equal-tailed property; more generally, $t\alpha$ and $(1-t)\alpha$ can be used for $t\in(0,1)$.

Figure \ref{fig:hutson-u-example} visualizes an example.  The beta distribution's mean is $u^h$ (or $u^l$).  Decreasing $u^h$ increases the probability mass in the shaded region below $u$, while increasing $u^h$ decreases the shaded region, and vice-versa for $u^l$.  Solving \eqref{eqn:uhdef} is a simple numerical search problem.
\begin{figure}[htbp]
  \centering
  \hspace*{\fill}
  \includegraphics[clip=true,trim=5 45 25 55,width=0.45\textwidth]{Hutson_u_example_lower.pdf}
  \hfill
  \includegraphics[clip=true,trim=5 45 25 55,width=0.45\textwidth]{Hutson_u_example_upper.pdf}
  \hspace*{\fill}
  \caption{\label{fig:hutson-u-example}Example of one-sided CI endpoint determination, $n=11$, $p=0.65$, $\alpha=0.1$.  Left: $u^l$ makes the shaded region's area $P\bigl(\tilde Q^I_U(u^l)>p\bigr)=\alpha$.  Right: similarly, $u^h$ solves $P\bigl(\tilde Q^I_U(u^h)<p\bigr)=\alpha$.}
\end{figure}

Lemma \ref{lem:u1u2-approx} shows the CI endpoint indices converge to $p$ at a $n^{-1/2}$ rate and may be approximated using quantiles of a normal distribution.
\begin{lemma}\label{lem:u1u2-approx}
Let $z_{1-\alpha}$ denote the $(1-\alpha)$-quantile of a standard normal distribution, $z_{1-\alpha}\equiv \Phi^{-1}(1-\alpha)$.  From the definitions in \eqref{eqn:uhdef}, the values $u^l(\alpha)$ and $u^h(\alpha)$ can be approximated as
\begin{align*}
u^l(\alpha)
  &= p - n^{-1/2}z_{1-\alpha}\sqrt{p(1-p)} - \frac{2p-1}{6n}(z_{1-\alpha}^2+2) +O(n^{-3/2}) , \\
u^h(\alpha)
  &= p + n^{-1/2}z_{1-\alpha}\sqrt{p(1-p)} - \frac{2p-1}{6n}(z_{1-\alpha}^2+2) +O(n^{-3/2}) .
\end{align*}
\end{lemma}


For the lower one-sided CI, using \eqref{eqn:uhdef}, the $1-\alpha$ CI from \citet{Hutson1999} is
\begin{align}
\label{eqn:hutson-CI-lower}
\Bigl(-\infty,\hat Q^L_X\bigl(u^h(\alpha)\bigr)\Bigr) .
\end{align}
 Coverage probability is
\begin{align*}
P & \mathopen{}\mathclose\bgroup\originalleft\{Q(p) \in \mathopen{}\mathclose\bgroup\originalleft(-\infty, \hat Q^L_X\bigl(u^h(\alpha)\bigr)\aftergroup\egroup\originalright)\aftergroup\egroup\originalright\}
   = P\mathopen{}\mathclose\bgroup\originalleft(\hat Q^L_X\bigl(u^h(\alpha)\bigr)>Q(p)\aftergroup\egroup\originalright) \\
  &\!\!\!\!\stackrel{\textrm{Thm \ref{thm:cdferror}}}{=}
     P\mathopen{}\mathclose\bgroup\originalleft(\tilde Q^I_X\bigl(u^h(\alpha)\bigr)>Q(p)\aftergroup\egroup\originalright)
    +\frac{\epsilon_h(1-\epsilon_h)z_{1-\alpha}\exp\{-z_{1-\alpha}^2/2\}} {\sqrt{2\pi}u^h(\alpha)(1-u^h(\alpha))} n^{-1}
    +O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright)  \\
  &= 1-\alpha
    +\frac{\epsilon_h(1-\epsilon_h)z_{1-\alpha}\phi(z_{1-\alpha})}{p(1-p)}n^{-1}
    +O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright),
\end{align*}
where $\phi(\cdot)$ is the standard normal PDF and the $n^{-1}$ term is non-negative.
Similar to the \citet{HoLee2005a} calibration, we can remove the analytic $n^{-1}$ term with the calibrated CI
\begin{align}
\label{eqn:hutson-CI-lower-calibrated}
  \mathopen{}\mathclose\bgroup\originalleft(-\infty, \hat Q^L_X\mathopen{}\mathclose\bgroup\originalleft(u^h\mathopen{}\mathclose\bgroup\originalleft(\alpha+\frac{\epsilon_h(1-\epsilon_h)z_{1-\alpha}\phi(z_{1-\alpha})}{p(1-p)}n^{-1}\aftergroup\egroup\originalright)\aftergroup\egroup\originalright)\aftergroup\egroup\originalright),
\end{align}
which has CPE of order $O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright)$.  We follow convention and define $\textrm{CPE}\equiv\textrm{CP}-(1-\alpha)$, where CP is the actual coverage probability and $1-\alpha$ the desired confidence level.

By parallel argument, \citeposs{Hutson1999} uncalibrated upper one-sided and two-sided CIs also have $O(n^{-1})$ CPE, or $O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright)$ with calibration.
For the upper one-sided case, again using \eqref{eqn:uhdef}, the $1-\alpha$ Hutson CI and our calibrated CI are respectively given by
\begin{equation}
\label{eqn:hutson-CI-upper}
\mathopen{}\mathclose\bgroup\originalleft(\hat Q^L_X\bigl(u^l(\alpha)\bigr) , \infty\aftergroup\egroup\originalright) , \quad
\mathopen{}\mathclose\bgroup\originalleft(\hat Q^L_X\mathopen{}\mathclose\bgroup\originalleft(u^l\mathopen{}\mathclose\bgroup\originalleft(\alpha + \frac{\epsilon_\ell(1-\epsilon_\ell)z_{1-\alpha}\phi(z_{1-\alpha})}{p(1-p)}n^{-1}\aftergroup\egroup\originalright)\aftergroup\egroup\originalright) , \infty\aftergroup\egroup\originalright) ,
\end{equation}
and for equal-tailed two-sided CIs,
\begin{align}
\label{eqn:hutson-CI-2s}
&\mathopen{}\mathclose\bgroup\originalleft(\hat Q^L_X\mathopen{}\mathclose\bgroup\originalleft[u^l\mathopen{}\mathclose\bgroup\originalleft(\alpha/2\aftergroup\egroup\originalright)\aftergroup\egroup\originalright],\hat Q^L_X\mathopen{}\mathclose\bgroup\originalleft(u^h\mathopen{}\mathclose\bgroup\originalleft(\alpha/2\aftergroup\egroup\originalright)\aftergroup\egroup\originalright)\aftergroup\egroup\originalright)   \quad\textrm{and}   \\
\label{eqn:hutson-CI-2s-calibrated}
\begin{split}
&\Bigg(\enspace \hat Q^L_X\mathopen{}\mathclose\bgroup\originalleft(u^l\mathopen{}\mathclose\bgroup\originalleft(\frac{\alpha}{2} + \frac{\epsilon_\ell(1-\epsilon_\ell)z_{1-\alpha/2}\phi(z_{1-\alpha/2})}{p(1-p)}n^{-1}\aftergroup\egroup\originalright)\aftergroup\egroup\originalright),  \\*
  &\qquad\quad  \hat Q^L_X\mathopen{}\mathclose\bgroup\originalleft(u^h\mathopen{}\mathclose\bgroup\originalleft(\frac{\alpha}{2} + \frac{\epsilon_h(1-\epsilon_h)z_{1-\alpha/2}\phi(z_{1-\alpha/2})}{p(1-p)}n^{-1}\aftergroup\egroup\originalright)\aftergroup\egroup\originalright)\enspace\Bigg).
\end{split}
\end{align}
Without calibration, in all cases the $n^{-1}$ CPE term is non-negative (indicating over-coverage).

For relatively extreme quantiles $p$ (given $n$), the $L$-statistic method cannot be computed because the $(n+1)$th (or zeroth) order statistic is needed.  In such cases, our code uses the Edgeworth expansion-based CI in \citet{Kaplan2015}.  Alternatively, if bounds on $X$ are known a priori, they may be used in place of these ``missing'' order statistics to generate conservative CIs.  Regardless, as $n\to\infty$, the range of computable quantiles approaches $(0,1)$.

The hypothesis tests corresponding to all the foregoing CIs achieve optimal asymptotic power against local alternatives.  The sample quantile is a semiparametric efficient estimator, so it suffices to show that power is asymptotically first-order equivalent to that of the test based on asymptotic normality.
Theorem \ref{thm:IDEAL-single} collects all of our results on coverage and power.

\begin{theorem}\label{thm:IDEAL-single}
Let $z_\alpha$ denote the $\alpha$-quantile of the standard normal distribution, and let $\epsilon_h=(n+1)u^h(\alpha)-\lfloor(n+1)u^h(\alpha)\rfloor$ and $\epsilon_\ell=(n+1)u^l(\alpha)-\lfloor(n+1)u^l(\alpha)\rfloor$.
 Let Assumption \ref{a:iid} hold, and let \ref{a:hut-pf} hold at $p$.  Then, we have the following.
\begin{enumerate} \item\label{cor:IDEAL-single-1s-CP} The one-sided lower and upper CIs in \eqref{eqn:hutson-CI-lower} and \eqref{eqn:hutson-CI-upper} have coverage probability
 \[ 1-\alpha + \frac{\epsilon(1-\epsilon)z_{1-\alpha}\phi(z_{1-\alpha})}{p(1-p)}n^{-1} + O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright),\]
 with $\epsilon=\epsilon_h$ for the former and $\epsilon=\epsilon_\ell$ for the latter.
 \item\label{cor:IDEAL-single-2s-CP} The equal-tailed, two-sided CI in \eqref{eqn:hutson-CI-2s} has coverage probability
 \[ 1-\alpha + \frac{[\epsilon_h(1-\epsilon_h)+\epsilon_\ell(1-\epsilon_\ell)]z_{1-\alpha/2}\phi(z_{1-\alpha/2})}{p(1-p)}n^{-1} + O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright). \]
 \item\label{cor:IDEAL-single-calibrated} The calibrated one-sided lower, one-sided upper, and two-sided equal-tailed CIs given in \eqref{eqn:hutson-CI-lower-calibrated}, \eqref{eqn:hutson-CI-upper}, and \eqref{eqn:hutson-CI-2s-calibrated}, respectively, have $O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright)$ CPE.
 \item\label{cor:IDEAL-single-power}
 The asymptotic probabilities of excluding $D_n=Q(p)+\kappa n^{-1/2}$ from lower one-sided (l), upper one-sided (u), and equal-tailed two-sided (t) CIs (i.e.,\ asymptotic power of the corresponding hypothesis tests) are
 \begin{align*}
  \mathcal P_n^l(D_n) &\to \Phi\mathopen{}\mathclose\bgroup\originalleft(z_\alpha+S\aftergroup\egroup\originalright), \;
  \mathcal P_n^u(D_n)  \to \Phi\mathopen{}\mathclose\bgroup\originalleft(z_\alpha-S\aftergroup\egroup\originalright), \;
  \mathcal P_n^t(D_n)  \to \Phi\mathopen{}\mathclose\bgroup\originalleft(z_{\alpha/2}+S\aftergroup\egroup\originalright)
                          +\Phi\mathopen{}\mathclose\bgroup\originalleft(z_{\alpha/2}-S\aftergroup\egroup\originalright),
 \end{align*}
 where $S\equiv \kappa f(F^{-1}(p))/\sqrt{p(1-p)}$.
\end{enumerate}
\end{theorem}

The equal-tailed property of our two-sided CIs is a type of median-unbiasedness.  If $(\hat{L},\hat{H})$ is a CI for scalar $\theta$, then an equal-tailed CI is ``unbiased'' under loss function $L(\theta,\hat L,\hat H)=\max\{0,\theta-\hat H,\hat L-\theta\}$, as defined in (5) of \citet{Lehmann1951}.
This median-unbiased property may be desirable \citep[e.g.,][footnote 11]{AndrewsGuggenberger2014}, although it is different than the usual ``unbiasedness'' where a CI is the inversion of an unbiased test.
More generally, in \eqref{eqn:hutson-CI-2s}, we could replace $u^l(\alpha/2)$ and $u^h(\alpha/2)$ by $u^l(t\alpha)$ and $u^h((1-t)\alpha)$ for $t\in[0,1]$.  Different $t$ may achieve different optimal properties, which we leave to future work.






\section{Quantile inference: conditional}\label{sec:inf-conditional}

\subsection{Setup and bias}\label{sec:setup}

Let $Q_{Y|X}(u;x)$ be the conditional $u$-quantile function of scalar outcome $Y$ given conditioning vector $X\in\mathcal X\subset\mathbb R^d$, evaluated at $X=x$.
The object of interest is $Q_{Y|X}(p;x_0)$, for $p\in(0,1)$ and interior point $x_0$.
The sample $\{Y_i,X_i\}_{i=1}^n$ is drawn iid.
Without loss of generality, let $x_0=0$.

If $X$ is discrete so that $P(X=0)>0$, we can take the subsample with $X_i=0$ and compute a CI from the corresponding $Y_i$ values, using the method in Section \ref{sec:inf-unconditional}.
Even with dependence like strong mixing among the $X_i$, CPE is the same $O(n^{-1})$ from Theorem \ref{thm:IDEAL-single} as long as the subsample's $Y_i$ are independent draws from the same $Q_{Y|X}(\cdot;0)$ and $N_n\stackrel{a.s.}{\asymp}n$.

If $X$ is continuous, then $P(X_i=0)=0$, so observations with $X_i\ne0$ must be included.
If $X$ contains mixed continuous and discrete components, then we can apply our method for continuous $X$ to each subsample corresponding to each unique value of the discrete subvector of $X$.
The asymptotic rates are unaffected by the presence of discrete variables (although the finite-sample consequences may deserve more attention), so we focus on the case where all components of $X$ are continuous.

We now present definitions and assumptions, continuing the normalization $x_0=0$.

\begin{defn}[local smoothness]\label{def:smoothness}
Following \citet[pp.\ 762--3]{Chaudhuri1991}: if, in a neighborhood of the origin, function $g(\cdot)$ is continuously differentiable through order $k$, and its $k$th derivatives are uniformly H\"older continuous with exponent $\gamma\in(0,1]$, then $g(\cdot)$ has ``local smoothness'' of degree $s=k+\gamma$.
\end{defn}

\begin{assumption}\label{a:sampling}
Sampling of $(Y_i,X_i')'$ is iid, for continuous scalar $Y_i$ and continuous vector $X_i\in\mathcal X\subseteq\mathbb R^d$.  The point of interest $X=0$ is in the interior of $\mathcal X$, and the quantile of interest is $p\in(0,1)$.
\end{assumption}
\begin{assumption}\label{a:f}
The marginal density of $X$, denoted $f_X(\cdot)$, satisfies $0<f_X(0)<\infty$ and has local smoothness $s_X=k_X+\gamma_X>0$.
\end{assumption}
\begin{assumption}\label{a:Q}
For all $u$ in a neighborhood of $p$, $Q_{Y|X}(u;\cdot)$ (as a function of the second argument) has local smoothness\footnote{Our $s_Q$ corresponds to variable $p$ in \citet{Chaudhuri1991}; \citet{BhattacharyaGangopadhyay1990} use $s_Q=2$ and $d=1$.} $s_Q=k_Q+\gamma_Q>0$.
\end{assumption}
\begin{assumption}\label{a:h}
As $n\to\infty$, the bandwidth satisfies
(i) $h\to0$,
(i') $h^{b+d/2}\sqrt{n}\to0$ with $b\equiv\min\{s_Q,s_X+1,2\}$,
(ii) $nh^d/[\log(n)]^2\to\infty$.
\end{assumption}
\begin{assumption}\label{a:GK1}
For all $u$ in a neighborhood of $p$ and all $x$ in a neighborhood of the origin, $f_{Y|X}\mathopen{}\mathclose\bgroup\originalleft(Q_{Y|X}(u;x);x\aftergroup\egroup\originalright)$ is uniformly bounded away from zero.
\end{assumption}
\begin{assumption}\label{a:GK2}
For all $y$ in a neighborhood of $Q_{Y|X}(p;0)$ and all $x$ in a neighborhood of the origin, $f_{Y|X}\mathopen{}\mathclose\bgroup\originalleft(y;x\aftergroup\egroup\originalright)$ has a second derivative in its first argument ($y$) that is uniformly bounded and continuous in $y$, having local smoothness $s_Y=k_Y+\gamma_Y>2$.
\end{assumption}

Definition \ref{def:cube} refers to a window whose size depends on $h$: $C_h=[-h,h]$ if $d=1$, or more generally a hypercube as in \citet[pp.\ 763]{Chaudhuri1991}:
letting $\|\cdot\|_\infty$ denote the $L_\infty$-norm,
\begin{align}
\label{eqn:Ch}
C_h &\equiv \{x:x\in\mathbb R^d,\|x\|_\infty\le h\}, \quad
N_n \equiv \#\bigl(\{Y_i: X_i\in C_h, 1\le i\le n\}\bigr)  .
\end{align}
\begin{defn}[local sample]\label{def:cube}
Using $C_h$ and $N_n$ defined in \eqref{eqn:Ch}, the ``local sample'' consists of $Y_i$ values from observations with $X_i\in C_h\subset\mathbb R^d$,
and the ``local sample size'' is $N_n$.
Additionally, let the local quantile function $Q_{Y|X}(p;C_h)$ be the $p$-quantile of $Y$ given $X\in C_h$, satisfying
$p = P\mathopen{}\mathclose\bgroup\originalleft(Y<Q_{Y|X}(p;C_h)\mid X\in C_h\aftergroup\egroup\originalright)$; similarly define the local CDF $F_{Y|X}(\cdot;C_h)$, local PDF $f_{Y|X}(\cdot;C_h)$, and derivatives thereof.
\end{defn}

Given fixed values of $n$ and $h$, Assumption \ref{a:sampling} implies that the $Y_i$ in the local sample are independent and identically distributed,\footnote{This may be the case asymptotically even with substantial dependence, although we do not explore this point.  For example, \citet[p.\ 237]{PolonikYao2002} write, ``Only the observations with $X_t$ in a small neighbourhood of $x$ are effectively used\ldots [which] are not necessarily close with each other in the time space. Indeed, they could be regarded as asymptotically independent under appropriate conditions such as strong mixing\ldots.''}
which is needed to apply Theorem \ref{thm:IDEAL-single}.  However, they do not have the quantile function of interest, $Q_{Y|X}(\cdot;0)$, but rather the biased $Q_{Y|X}(\cdot;C_h)$.  This is like drawing a global (any $X_i$) iid sample of wages, $Y_i$, and restricting it to observations in Japan ($X\in C_h$) when our interest is only in Tokyo ($X=0$): our restricted $Y_i$ constitute an iid sample from Japan, but the $p$-quantile wage in Japan may differ from that in Tokyo.
Assumptions \ref{a:f}--\ref{a:h}(i) and \ref{a:GK2} are necessary for the calculation of this bias, $Q_{Y|X}(p;C_h)-Q_{Y|X}(p;0)$, in Lemma \ref{lem:bias}.
Assumptions \ref{a:h}(ii) and \ref{a:GK1} (and \ref{a:sampling}) ensure $N_n\stackrel{a.s.}{\to}\infty$.  Assumptions \ref{a:GK1} and \ref{a:GK2} are conditional versions of Assumptions \ref{a:hut-pf}(i) and \ref{a:hut-pf}(ii), respectively. Their uniformity ensures uniformity of the remainder term in Theorem \ref{thm:IDEAL-single}, accounting for the fact that the local sample's distribution, $F_{Y|X}(\cdot;C_h)$, changes with $n$ (through $h$ and $C_h$).

From \ref{a:h}(i), asymptotically $C_h$ is entirely contained within the neighborhoods implicit in \ref{a:f}, \ref{a:Q}, and \ref{a:GK2}.  This in turn allows us to examine only a local neighborhood around $p$ (e.g.,\ as in \ref{a:Q}) since the CI endpoints converge to the true value at a $N_n^{-1/2}$ rate.

The $X_i$ being iid helps guarantee that $N_n$ is almost surely of order $nh^d$.
The $h^d$ comes from the volume of $C_h$.
Larger $h$ lowers CPE via $N_n$ but raises CPE via bias.  This tradeoff determines the optimal rate at which $h\to0$ as $n\to\infty$.  Using Theorem \ref{thm:IDEAL-single} and additional results on CPE from bias below, we determine the optimal value of $h$.

\begin{defn}[steps to compute CI for $Q_{Y|X}(p;0)$]\label{def:method}
First, $C_h$ and $N_n$ are calculated as in Definition \ref{def:cube}.
Second, using the $Y_i$ from observations with $X_i\in C_h$, a $p$-quantile CI is constructed as in \citet{Hutson1999}.
If additional discrete conditioning variables exist, then repeat separately for each combination of discrete conditioning values.
This procedure may be repeated for any number of $x_0$.
For the bandwidth, we recommend the formulas in Section \ref{sec:opt-h}.
\end{defn}


The bias characterized in Lemma \ref{lem:bias} is the difference between these two population conditional quantiles.
\begin{lemma}\label{lem:bias}
Define $b$ as in \ref{a:h} and let $B_h\equiv Q_{Y|X}(p;C_h)-Q_{Y|X}(p;0)$.
If Assumptions \ref{a:f}, \ref{a:Q}, \ref{a:h}(i), \ref{a:GK1}, and \ref{a:GK2} hold, then the bias is of order
$|B_h| = O(h^b)$.
Defining
\begin{equation*}
\xi_p \equiv Q_{Y|X}(p;0) , \quad
F_{Y|X}^{(0,1)}(\xi_p;0) \equiv \mathopen{}\mathclose\bgroup\originalleft.\frac{\partial}{\partial x}F_{Y|X}(\xi_p;x)\aftergroup\egroup\originalright|_{x=0} , \quad
F_{Y|X}^{(0,2)}(\xi_p;0) \equiv \mathopen{}\mathclose\bgroup\originalleft.\frac{\partial^2}{\partial x^2}F_{Y|X}(\xi_p;x)\aftergroup\egroup\originalright|_{x=0} ,
\end{equation*}
with $d=1$, $k_X\ge1$, and $k_Q\ge2$, the bias is
\begin{equation}\label{eqn:Bh}
B_h
   = -h^2
      \frac{f_X(0) F_{Y|X}^{(0,2)}(\xi_p;0)
            +2 f_X'(0) F_{Y|X}^{(0,1)}(\xi_p;0)}
           {6 f_X(0) f_{Y|X}(\xi_p;0)}
           +o(h^2) .
\end{equation}
\end{lemma}
Equation \eqref{eqn:Bh} is the same as in \citet{BhattacharyaGangopadhyay1990}, who derive it using different arguments.

















\subsection{Optimal CPE order}\label{sec:opt-CPE}

The CPE-optimal bandwidth minimizes the sum of the two dominant high-order CPE terms.  It must be small enough to control the $O(h^b+N_nh^{2b})$ (two-sided) CPE from bias, but large enough to control the $O(N_n^{-1})$ CPE from applying the unconditional $L$-statistic method.
The following theorem summarizes optimal bandwidth and CPE results.

\begin{theorem}\label{thm:h-rate}
Let Assumptions \ref{a:sampling}--\ref{a:GK2} hold.
The following results are for the method in Definition \ref{def:method}.
For a one-sided CI, the bandwidth $h^*$ minimizing CPE has rate
$ h^* \asymp n^{-3/(2b+3d)} $,
corresponding to CPE of order $O(n^{-2b/(2b+3d)})$.
For a two-sided CI, the optimal bandwidth rate is $h^*\asymp n^{-1/(b+d)}$, and the optimal CPE is $O(n^{-b/(b+d)})$.
Using the calibration in Section \ref{sec:inf-unconditional}, if $p=1/2$, then the nearly (up to $\log(n)$) CPE-optimal two-sided bandwidth rate is $h^*\asymp n^{-5/(4b+5d)}$, yielding CPE of order $O\mathopen{}\mathclose\bgroup\originalleft(n^{-6b/(4b+5d)}[\log(n)]^3\aftergroup\egroup\originalright)$; if $p\ne1/2$, then $h^*\asymp n^{-3/(b+3d)}$ and CPE is $O\mathopen{}\mathclose\bgroup\originalleft( n^{-3b/(2b+6d)} [\log(n)]^3 \aftergroup\egroup\originalright)$.
The nearly CPE-optimal calibrated one-sided bandwidth rate is $h^*\asymp n^{-2/(b+2d)}$, yielding CPE of order $O\mathopen{}\mathclose\bgroup\originalleft(n^{-3b/(2b+4d)}[\log(n)]^3\aftergroup\egroup\originalright)$.
\end{theorem}







As detailed in the supplemental appendix, Theorem \ref{thm:h-rate} implies that for the most common values of dimension $d$ and most plausible values of smoothness $s_Q$, even our uncalibrated method is more accurate than inference based on asymptotic normality with a local polynomial estimator.
The same comparisons apply to basic bootstraps, which claim no refinement over asymptotic normality; in this (quantile) case, even Studentization does not improve theoretical CPE without the added complications of smoothed or $m$-out-of-$n$ bootstraps.

The only opportunity for normality to yield smaller CPE is to greatly reduce bias by using a very large local polynomial if $s_Q$ is large; our approach implicitly uses a uniform kernel, so bias reduction beyond $O(h^2)$ is impossible.  Nonetheless, our method has smaller CPE when $d=1$ or $d=2$ even if $s_Q=\infty$, and in other cases the necessary local polynomial degree may be prohibitively large given common sample sizes.



\begin{figure}[htbp]
\centering
\hfill
 \includegraphics[clip=true,trim=10 10 30 70,width=0.45\textwidth]{CPE_comp_s2.pdf}
 \hfill
 \includegraphics[clip=true,trim=10 10 30 70,width=0.45\textwidth]{CPE_comp_match.pdf}
 \hfill
 \null
 \caption{\label{fig:CPE-comp}Two-sided CPE comparison between new (``L-stat'') method and the local polynomial asymptotic normality method based on \citet{Chaudhuri1991}.  Left: with $s_Q=2$ and $s_X=1$, writing CPE as $n^\kappa$, comparison of $\kappa$ for different methods and different values of $d$.  Right: required smoothness $s_Q$ for the local polynomial normality-based CPE to match that of L-stat, as well as the corresponding number of terms in the local polynomial, for different $d$.}
\end{figure}

Figure \ref{fig:CPE-comp} (left panel) shows that if $s_Q=2$ and $s_X=1$, then the optimal CPE from asymptotic normality is always larger (worse) than our method's CPE.
As shown in the supplement, CPE with normality is nearly $O\bigl(n^{-2/(4+2d)}\bigr)$.
With $d=1$, this is $O\bigl(n^{-1/3}\bigr)$, much larger than our two-sided $O\bigl(n^{-2/3}\bigr)$.
With $d=2$, $O\bigl(n^{-1/4}\bigr)$ is larger than our $O\bigl(n^{-1/2}\bigr)$.
It remains larger for all $d$ since the bias is the same for both methods while the unconditional $L$-statistic inference is more accurate than normality.

Figure \ref{fig:CPE-comp} (right panel) shows the required amount of smoothness and local polynomial degree for asymptotic normality to match our method's CPE.
For the most common cases of $d=1$ and $d=2$, two-sided CPE with normality is larger even with infinite smoothness and a hypothetical infinite-degree polynomial.
With $d=3$, to match our CPE, normality needs $s_Q\ge12$ and a local polynomial of degree $k_Q\ge11$.
Since interaction terms are required, an $11$th-degree polynomial has $\sum_{T=d-1}^{k_Q+d-1} \binom{T}{d-1} = 364$ terms, which requires a large $N_n$ (and yet larger $n$).
As $d\to\infty$, the required number of terms in the local polynomial only grows larger and may be prohibitive in realistic finite samples.





\subsection{Plug-in bandwidth}\label{sec:opt-h}

We propose a feasible bandwidth value with the CPE-optimal rate.
To avoid recursive dependence on $\epsilon$ (the interpolation weight), we fix its value.  This does not achieve the theoretical optimum, but it remains close even in small samples and seems to work well in practice.  The CPE-optimal bandwidth value derivation is shown for $d=1$ in the supplemental appendix; a plug-in version is implemented in our code.
For reference, the plug-in bandwidth expressions are collected here.
The $\alpha$-quantile of $N(0,1)$ is again denoted $z_\alpha$.
We let $\hat B_h$ denote the estimator of bias term $B_h$; $\hat f_X$ the estimator of $f_X(x_0)$; $\hat f_X'$ the estimator of $f_X'(x_0)$; $\hat F_{Y|X}^{(0,1)}$ the estimator of $F_{Y|X}^{(0,1)}(\xi_p;x_0)$; and $\hat F_{Y|X}^{(0,2)}$ the estimator of $F_{Y|X}^{(0,2)}(\xi_p;x_0)$, with notation from Lemma \ref{lem:bias}.

When $d=1$, the following are our CPE-optimal plug-in bandwidths.
\begin{itemizecomp}
\item For one-sided inference, let
\begin{align}
\label{eqn:h-plugin-1s-pm}
\hat h_{+-}
  &= n^{-3/7}
     \mathopen{}\mathclose\bgroup\originalleft(
      \frac{z_{1-\alpha}}
           {3 \mathopen{}\mathclose\bgroup\originalleft[p(1-p)\hat f_X\aftergroup\egroup\originalright]^{1/2}
            \mathopen{}\mathclose\bgroup\originalleft[\hat f_X \hat F_{Y|X}^{(0,2)}
                   +2 \hat f_X' \hat F_{Y|X}^{(0,1)} \aftergroup\egroup\originalright]}
     \aftergroup\egroup\originalright)^{2/7}  , \\
\label{eqn:h-plugin-1s-pp}
\hat h_{++}
  &= -0.770\hat h_{+-}  .
\end{align}
  For lower one-sided inference, $\hat h_{+-}$ should be used if $\hat B_h<0$, and $\hat h_{++}$ otherwise.
  For upper one-sided inference, $\hat h_{++}$ should be used if $\hat B_h<0$, and $\hat h_{+-}$ otherwise.
\item For two-sided inference with general $p\in(0,1)$,
\begin{align}
\label{eqn:h-plugin-2s-general}
\hat h
  &= n^{-1/3} \mathopen{}\mathclose\bgroup\originalleft( \frac
       { (\hat B_h/|\hat B_h|) (1-2p)
        +\sqrt{(1-2p)^2  +4} }
       {2 \mathopen{}\mathclose\bgroup\originalleft|   \hat f_X  \hat F_{Y|X}^{(0,2)}
                +2 \hat f_X' \hat F_{Y|X}^{(0,1)} \aftergroup\egroup\originalright| }
     \aftergroup\egroup\originalright)^{1/3} ,
\end{align}
which simplifies to $\hat h = n^{-1/3} \bigl| \hat f_X \hat F_{Y|X}^{(0,2)} +2 \hat f_X' \hat F_{Y|X}^{(0,1)} \bigr|^{-1/3}$ with $p=0.5$.
\end{itemizecomp}

While we suggest the CPE-optimal bandwidths for moderate $n$, we suggest shifting toward a larger bandwidth as $n\to\infty$.  Once CPE is small over a range of bandwidths, a larger bandwidth in that range is preferable since it yields shorter CIs.
As an initial suggestion, we use a coefficient of $\max\{1,n/1000\}^{5/60}$ that keeps the CPE-optimal bandwidth for $n\le1000$ and then moves toward a $n^{-1/20}$ under-smoothing of the MSE-optimal bandwidth rate, as in \citet[p.\ 205]{FanLiu2016}.
































\section{Empirical application}\label{sec:empirical}

We present an application of our $L$-statistic inference to \citet{Engel1857} curves.  Code is available from the latter author's website, and the data are publicly available.

\citet{BanksEtAl1997} argue that a linear Engel curve is sufficient for certain categories of expenditure, while adding a quadratic term suffices for others.  Their Figure 1 shows nonparametrically estimated mean Engel curves (budget share $W$ against log total expenditure $\ln(X)$) with $95\%$ pointwise CIs at the deciles of the total expenditure distribution, using a subsample of 1980--1982 U.K.\ Family Expenditure Survey (FES) data.

We present a similar examination, but for quantile Engel curves in the 2001--2012 U.K.\ Living Costs and Food Surveys \citep{UKLCFS2012}, which is a successor to the FES.  We examine the same four categories as in the original analysis: food; fuel, light, and power (``fuel''); clothing and footwear (``clothing''); and alcohol.  We use the subsample of households with one adult male and one adult female (and possibly children) living in London or the South East, leaving $8{,}528$ observations.  Expenditure amounts are adjusted to 2012 nominal values using annual CPI data.{\interfootnotelinepenalty10000\footnote{\url{http://www.ons.gov.uk/ons/datasets-and-tables/data-selector.html?cdid=D7BT&dataset=mm23&table-id=1.1}}}

\begin{table}[htbp]
\centering
\caption{\label{tab:emp}$L$-statistic $99\%$ CIs for various unconditional quantiles ($p$) of the budget share distribution, for different categories of expenditure described in the text.}
\begin{tabular}{lccc}
\multicolumn{1}{c}{Category} & $p=0.5$ & $p=0.75$ & $p=0.9$ \\
\hline
food     & (0.1532,0.1580) & (0.2095,0.2170) & (0.2724,0.2818) \\
fuel     & (0.0275,0.0289) & (0.0447,0.0470) & (0.0692,0.0741) \\
clothing & (0.0135,0.0152) & (0.0362,0.0397) & (0.0697,0.0761) \\
alcohol  & (0.0194,0.0226) & (0.0548,0.0603) & (0.1012,0.1111) \\
\hline
\hline
\end{tabular}
\end{table}

Table \ref{tab:emp} shows unconditional $L$-statistic CIs for various quantiles of the budget share distributions for the four expenditure categories.  (Due to the large sample size, calibrated CIs are identical at the precision shown.)  These capture some population features, but the conditional quantiles are of more interest.

\begin{figure}[htbp]
 \centering
\hfill
 \includegraphics[clip=true,trim=10 40 30 70,width=0.49\textwidth]{quantile_inf_np_ex_jt_food.pdf}
 \hfill
 \includegraphics[clip=true,trim=40 40 0 70,width=0.49\textwidth]{quantile_inf_np_ex_jt_fuel.pdf}
\hfill\null
 \\
\hfill
 \includegraphics[clip=true,trim=10 5 30 80,width=0.49\textwidth]{quantile_inf_np_ex_jt_clothing.pdf}
 \hfill
 \includegraphics[clip=true,trim=40 5 0 80,width=0.49\textwidth]{quantile_inf_np_ex_jt_alc.pdf}
\hfill\null
 \caption{\label{fig:emp}Joint (over the nine expenditure levels) $90\%$ confidence intervals for quantile Engel curves: food (top left), fuel (top right), clothing (bottom left), and alcohol (bottom right).}
\end{figure}

Figure \ref{fig:emp} is comparable to Figure 1 of \citet{BanksEtAl1997} but with $90\%$ joint (over the nine expenditure levels) CIs instead of $95\%$ pointwise CIs, alongside quadratic quantile regression estimates.  (To get joint CIs, we simply use the Bonferroni adjustment and compute $1-\alpha/9$ pointwise CIs.)  Joint CIs are more intuitive for assessing the shape of a function since they jointly cover all corresponding points on the true curve with $90\%$ probability, rather than any given single point.  The CIs are interpolated only for visual convenience.  Although some of the joint CI shapes do not look quadratic at first glance, the only cases where the quadratic fit lies outside one of the intervals are for alcohol at the conditional median and clothing at the conditional upper quartile, and neither is a radical departure.  With a $90\%$ confidence level and 12 confidence sets, we would not be surprised if one or two did not cover the true quantile Engel curve completely.  Importantly, the CIs are relatively precise, too; the linear fit is rejected in 8 of 12 cases.  Altogether, this evidence suggests that the benefits of a quadratic (but not linear) approximation may outweigh the cost of approximation error.

The supplemental appendix includes a similar figure but with a nonparametric (instead of quadratic) conditional quantile estimate along with joint CIs from \citet{FanLiu2016}.

























\section{Simulation study}\label{sec:sim}

Code for our methods and simulations is available on the latter author's website.


\subsection{Unconditional simulations}\label{sec:sim-unconditional}

We compare two-sided unconditional CIs from the following methods: ``L-stat'' from Section \ref{sec:inf-unconditional}, originally in \citet{Hutson1999}; ``BH'' from \citet{BeranHall1993}; ``Norm'' using the sample quantile's asymptotic normality and kernel-estimated variance; ``K15'' from \citet{Kaplan2015}; and ``BStsym,'' a symmetric Studentized bootstrap (99 draws) with bootstrapped variance (100 draws).\footnote{Other bootstraps were consistently worse in terms of coverage: (asymmetric) Studentized bootstrap, and percentile bootstrap with and without symmetry.}


Overall, L-stat and BH have the most accurate coverage probability (CP), avoiding under-coverage while maintaining shorter length than other methods achieving at least $95\%$ CP.  Near the median, L-stat and BH are nearly identical.  Away from the median, L-stat is closer to equal-tailed and often shorter than BH.  Farther into the tails, L-stat can be computed where BH cannot.

\begin{table}[htbp]
\caption{\label{tab:sim-un1}CP and median CI length, $1-\alpha=0.95$; $n$, $p$, and distributions of $X_i$ ($F$) shown in table; $10{,}000$ replications.  ``Too high'' is the proportion of simulation draws in which the lower endpoint was above the true $F^{-1}(p)$, and ``too low'' is the proportion when the upper endpoint was below $F^{-1}(p)$.}
\centering
\begin{tabular}[c]{cccccccc}
 $n$  &   $p$   &       $F$       & Method & CP & Too low & Too high & Length \\
\hline
$ 25$ & $0.5  $ &          Normal & L-stat & 0.953 & 0.022 & 0.025 & 0.99 \\
$ 25$ & $0.5  $ &          Normal & BH     & 0.955 & 0.021 & 0.024 & 1.00 \\
$ 25$ & $0.5  $ &          Normal & Norm   & 0.942 & 0.028 & 0.030 & 1.02 \\
$ 25$ & $0.5  $ &          Normal & K15    & 0.971 & 0.014 & 0.015 & 1.19 \\
$ 25$ & $0.5  $ &          Normal & BStsym & 0.942 & 0.028 & 0.030 & 1.13 \\[2pt]
$ 25$ & $0.5  $ &         Uniform & L-stat & 0.953 & 0.022 & 0.025 & 0.37 \\
$ 25$ & $0.5  $ &         Uniform & BH     & 0.954 & 0.021 & 0.025 & 0.37 \\
$ 25$ & $0.5  $ &         Uniform & Norm   & 0.908 & 0.046 & 0.046 & 0.35 \\
$ 25$ & $0.5  $ &         Uniform & K15    & 0.963 & 0.018 & 0.020 & 0.44 \\
$ 25$ & $0.5  $ &         Uniform & BStsym & 0.937 & 0.031 & 0.032 & 0.45 \\[2pt]
$ 25$ & $0.5  $ &     Exponential & L-stat & 0.953 & 0.024 & 0.023 & 0.79 \\
$ 25$ & $0.5  $ &     Exponential & BH     & 0.954 & 0.024 & 0.022 & 0.80 \\
$ 25$ & $0.5  $ &     Exponential & Norm   & 0.924 & 0.056 & 0.020 & 0.75 \\
$ 25$ & $0.5  $ &     Exponential & K15    & 0.968 & 0.022 & 0.010 & 0.96 \\
$ 25$ & $0.5  $ &     Exponential & BStsym & 0.941 & 0.039 & 0.020 & 0.93 \\
\hline
\hline
\end{tabular}
\end{table}

Table \ref{tab:sim-un1} shows nearly exact CP for both L-stat and BH when $n=25$ and $p=0.5$.
``Norm'' can be slightly shorter, but it under-covers.
The bootstrap has only slight under-coverage, and K15 none, but their CIs are longer than L-stat's.
Additional results are in the supplemental appendix, but the qualitative points are the same.

\begin{table}[htbp]
\caption{\label{tab:sim-un2}CP and median CI length, as in Table \ref{tab:sim-un1}.}
\centering
\begin{tabular}[c]{cccccccc}
 $n$  &   $p$   &       $F$       & Method & CP & Too low & Too high & Length \\
\hline
$ 99$ & $0.037$ &          Normal & L-stat & 0.951 & 0.023 & 0.026 & 1.02 \\
$ 99$ & $0.037$ &          Normal & BH     &   NA &   NA &   NA &   NA \\
$ 99$ & $0.037$ &          Normal & Norm   & 0.925 & 0.016 & 0.059 & 0.83 \\
$ 99$ & $0.037$ &          Normal & K15    & 0.970 & 0.009 & 0.021 & 1.55 \\
$ 99$ & $0.037$ &          Normal & BStsym & 0.950 & 0.020 & 0.030 & 1.20 \\[2pt]
$ 99$ & $0.037$ &          Cauchy & L-stat & 0.950 & 0.022 & 0.028 & 39.37 \\
$ 99$ & $0.037$ &          Cauchy & BH     &   NA &   NA &   NA &   NA \\
$ 99$ & $0.037$ &          Cauchy & Norm   & 0.784 & 0.082 & 0.134 & 18.90 \\
$ 99$ & $0.037$ &          Cauchy & K15    & 0.957 & 0.002 & 0.041 & 36.55 \\
$ 99$ & $0.037$ &          Cauchy & BStsym & 0.961 & 0.002 & 0.037 & 48.77 \\[2pt]
$ 99$ & $0.037$ &         Uniform & L-stat & 0.951 & 0.024 & 0.026 & 0.07 \\
$ 99$ & $0.037$ &         Uniform & BH     &   NA &   NA &   NA &   NA \\
$ 99$ & $0.037$ &         Uniform & Norm   & 0.990 & 0.000 & 0.010 & 0.12 \\
$ 99$ & $0.037$ &         Uniform & K15    & 0.963 & 0.028 & 0.009 & 0.11 \\
$ 99$ & $0.037$ &         Uniform & BStsym & 0.924 & 0.053 & 0.022 & 0.08 \\
\hline
\hline
\end{tabular}
\end{table}

Table \ref{tab:sim-un2} shows a case in the lower tail with $n=99$ where BH cannot be computed (because it needs the zeroth order statistic).
Even then, L-stat's CP remains almost exact, and it is closest to equal-tailed.
``Norm'' under-covers for two $F$ (severely for Cauchy) and is almost twice as long as L-stat for the third.
BStsym has less under-coverage, and K15 none, but both are generally longer than L-stat.
Again, additional results are in the supplemental appendix, with similar patterns.


The supplemental appendix contains additional simulation results for $p\ne0.5$ but where BH is still computable.
L-stat and BH both attain $95\%$ CP, but L-stat is much closer to equal-tailed and is shorter.
The supplemental appendix also has results illustrating the effect of calibration.




Table \ref{tab:sim-beta-norm1} isolates the effects of using the beta distribution rather than the normal approximation, as well as the effects of interpolation.  Method ``Normal'' uses the normal approximation to determine $u^h$ and $u^l$ but still interpolates, while ``Norm/floor'' uses the normal approximation with no interpolation as in equations (5) and (6) of \citet[Ex.\ 2.1]{FanLiu2016}.

\begin{table}[htbp]
\caption{\label{tab:sim-beta-norm1}CP and median CI length, $n=19$, $Y_i\stackrel{iid}{\sim}N(0,1)$, $1-\alpha=0.90$, $1{,}000$ replications, various $p$.  In parentheses below CP are probabilities of being too low or too high, as in Table \ref{tab:sim-un1}.
Methods are described in the text.}
\centering
\begin{tabular}[c]{lcccccccc}
 && \multicolumn{3}{c}{Two-sided CP} && \multicolumn{3}{c}{} \\
 && \multicolumn{3}{c}{(Too low, Too high)} && \multicolumn{3}{c}{Median length} \\
\cline{3-5}\cline{7-9}
Method && $p=0.15$ & $p=0.25$ & $p=0.5$ && $p=0.15$ & $p=0.25$ & $p=0.5$ \\
\hline
L-stat      && 0.905 & 0.901 & 0.898 && 1.20 & 1.03 & 0.93 \\
            && (0.048,0.047) & (0.050,0.049) & (0.052,0.050) &&  &  &  \\
Normal      &&    NA & 0.926 & 0.912 &&   NA & 1.22 & 1.00 \\
            && (NA,NA) & (0.062,0.012) & (0.045,0.043) &&  &  &  \\
Norm/floor  &&    NA & 0.913 & 0.876 &&   NA & 1.47 & 0.91 \\
            && (NA,NA) & (0.083,0.004) & (0.087,0.037) &&  &  &  \\
\hline
\end{tabular}
\end{table}

Table \ref{tab:sim-beta-norm1} shows several advantages of L-stat.
First, for $p=0.15$, Normal and Norm/floor cannot even be computed (hence ``NA'') because they require the zeroth order statistic, which does not exist, whereas L-stat is computable and has nearly exact CP ($0.905$).
Second, with $p=0.25$ and $p=0.5$, the normal approximation (Normal) makes the CI needlessly longer than L-stat's CI.
Third, additionally not interpolating (Norm/floor) makes the CI even longer for $p=0.25$ but leads to under-coverage for $p=0.5$.
Fourth, whereas the L-stat CIs are almost exactly equal-tailed, the normal-based CIs are far from equal-tailed at $p=0.25$, where Norm/floor is essentially a one-sided CI.










\subsection{Conditional simulations}\label{sec:sim-conditional}


For conditional quantile inference, we compare our $L$-statistic method
(``L-stat'') with a variety of others.
Implementation details may be seen in the supplemental appendix and available code.
The first other method (``rqss'') is from the popular \texttt{quantreg} package in R \citep{R.quantreg}.
The second (``boot'') is a local cubic method following \citet{Chaudhuri1991} but with bootstrapped standard errors;
the bandwidth is L-stat's multiplied by $n^{1/12}$ to get the local cubic CPE-optimal rate.
The third (``QYg'') uses the asymptotic normality of a local linear estimator with a Gaussian kernel, using results and ideas from \citet{QuYoon2015}, although they are more concerned with uniform (over quantiles) inference; they suggest using the MSE-optimal bandwidth (Corollary 1) and a particular type of bias correction (Remark 7).
The fourth (``FLb'') is from Section 3.1 in \citet{FanLiu2016}, based on a symmetrized $k$-NN estimator using a bisquare kernel; we use the code from their simulations.\footnote{Graciously provided to us.  The code differs somewhat from the description in their text, most notably by an additional factor of $0.4$ in the bandwidth.}  Interestingly, although in principle they are just slightly undersmoothing the MSE-optimal bandwidth, their bandwidth is very close to the CPE-optimal bandwidth for the sample sizes considered.

We now write $x_0$ as the point of interest, instead of $x_0=0$; we also take $d=1$, $b=2$, and focus on two-sided inference, both pointwise (single $x_0$) and joint (over multiple $x_0$).
Joint CIs for all methods are computed using the Bonferroni approach.
Uniform bands are also examined, with L-stat, QYg, and boot relying on the adjusted critical value from the \citet{Hotelling1939} tube computations in \texttt{plot.rqss}.
Each simulation has $1{,}000$ replications unless otherwise noted.


Figure \ref{fig:FLsim} uses Model 1 from \citet[p.\ 205]{FanLiu2016}:
$Y_i=2.5+\sin(2X_i)+2\exp\bigl(-16X_i^2\bigr)+0.5\epsilon_i$,
$X_i\stackrel{iid}{\sim}N(0,1)$,
$\epsilon_i\stackrel{iid}{\sim}N(0,1)$,
$X_i\protect\mathpalette{\protect\independenT}{\perp}\epsilon_i$,
$n=500$, $p=0.5$.
The ``Direct'' method in their Table 1 is our FLb.
All methods have good pointwise CP (top left).
L-stat has the best pointwise power (top right).

\begin{figure}[thbp]
  \centering
\hfill
  \includegraphics[clip=true,trim=15 15 10 58,width=0.445\textwidth]
    {qinfnp_2015_08_26_ptCP_F001_a05_n500_p50_numx3_trimx04_nrep1000.pdf}
    \hfill
  \includegraphics[clip=true,trim=15 15 10 58,width=0.445\textwidth]
    {qinfnp_2015_08_26_ptPWR2_F001_a05_n500_p50_numx3_trimx04_nrep1000.pdf}
\hfill
  \null
\\
\hfill
  \includegraphics[clip=true,trim=15 15 10 58,width=0.445\textwidth]
    {qinfnp_2015_08_26_jtPWR_F001_a05_n500_p50_numx3_trimx04_nrep1000.pdf}
    \hfill
  \includegraphics[clip=true,trim=15 15 10 58,width=0.445\textwidth]
    {qinfnp_2015_08_27_unifPWR_F001_a05_n500_p50_numx231_trimx04_nrep1000.pdf}
\hfill
\null
  \caption{\label{fig:FLsim}Results from DGP in Model 1 of \citet{FanLiu2016}, $n=500$, $p=0.5$.  Top left: pointwise CP at $x_0\in\{0,0.75,1.5\}$, interpolated for visual ease.  Top right: pointwise power at the same $x_0$ against deviations of $\pm0.1$.  Bottom left: joint power curves.  Bottom right: uniform power curves.}
\end{figure}

Figure \ref{fig:FLsim} (bottom left) shows power curves of the hypothesis tests corresponding to the joint (over $x_0\in\{0,0.75,1.5\}$) CIs, varying $H_0$ while maintaining the same DGP.  The deviations of $Q_{Y|X}(p;x_0)$ shown on the horizontal axis are the same at each $x_0$; zero deviation implies $H_0$ is true, in which case the rejection probability is the type I error rate.
All methods have good type I error rates: L-stat's is 6.2\%, and other methods' are below the nominal 5\%.  L-stat has significantly better power, an advantage of 20--40\% at the larger deviations.
The bottom right graph in Figure \ref{fig:FLsim} is similar, but based on uniform confidence bands evaluated at $231$ different $x_0$.  Only L-stat has nearly exact type I error rate and good power.



Next, we use the simulation setup of the \texttt{rqss} vignette in \citet{R.quantreg}, which in turn came in part from \citet[\S17.5.1]{RuppertEtAl2003}.  Here, $n=400$, $p=0.5$, $d=1$, $\alpha=0.05$, and
\begin{equation}\label{eqn:DGP-rqss}
X_i \stackrel{iid}{\sim}\textrm{Unif}(0,1), \quad
Y_i = \sqrt{X_i(1-X_i)}\sin\bigl(2\pi(1+2^{-7/5})/(X_i+2^{-7/5})\bigr) +\sigma(X_i)U_i,
\end{equation}
where the $U_i$ are iid $N(0,1)$, $t_3$, Cauchy, or centered $\chi^2_3$, and $\sigma(X)=0.2$ or $\sigma(X)=0.2(1+X)$.
The conditional median function is graphed in the supplemental appendix.
Although the function as a whole is not a common shape in economics (with multiple local maxima and minima), it provides insight into different types of functions at different points.
For pointwise and joint CIs, we consider 47 equispaced points, $x_0=0.04,0.06,\ldots,0.96$; uniform confidence bands are evaluated at 231 equispaced values of $x_0$.



\begin{figure}[htbp]
  \centering
  \hfill
  \includegraphics[clip=true,trim=20 45 30 25,width=0.32\textwidth]
    {qinfnp_2015_08_16_ptCP_F010_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=51 45 30 25,width=0.30\textwidth]
    {qinfnp_2015_08_16_ptCP_F011_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=10 45 30 25,width=0.329\textwidth]
    {qinfnp_2015_08_16_jtPWR_F010_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill\null
\\
  \hfill
  \includegraphics[clip=true,trim=20 45 30 70,width=0.32\textwidth]
    {qinfnp_2015_08_16_ptCP_F012_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=51 45 30 70,width=0.30\textwidth]
    {qinfnp_2015_08_16_ptCP_F013_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=10 45 30 70,width=0.329\textwidth]
    {qinfnp_2015_08_16_jtPWR_F012_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill\null
\\
  \hfill
  \includegraphics[clip=true,trim=20 45 30 70,width=0.32\textwidth]
    {qinfnp_2015_08_16_ptCP_F014_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=51 45 30 70,width=0.30\textwidth]
    {qinfnp_2015_08_16_ptCP_F015_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=10 45 30 70,width=0.329\textwidth]
    {qinfnp_2015_08_16_jtPWR_F014_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill\null
\\
  \hfill
  \includegraphics[clip=true,trim=20 15 30 70,width=0.32\textwidth]
    {qinfnp_2015_08_16_ptCP_F016_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=51 15 30 70,width=0.30\textwidth]
    {qinfnp_2015_08_16_ptCP_F017_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=10 15 30 70,width=0.329\textwidth]
    {qinfnp_2015_08_16_jtPWR_F016_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill\null
  \caption{\label{fig:ptCP-1}Pointwise CP (first two columns) and joint power curves (third column), $1-\alpha=0.95$, $n=400$, $p=0.5$, DGP in \eqref{eqn:DGP-rqss}.  Distributions of $U_i$ are, top row to bottom row: $N(0,1)$, $t_3$, Cauchy, and centered $\chi^2_3$.  Columns 1 \& 3: $\sigma(x)=0.2$; Column 2: $\sigma(x)=(0.2)(1+x)$.}
\end{figure}

Figure \ref{fig:ptCP-1}'s first two columns show that across all eight DGPs (four error distributions, homoskedastic or heteroskedastic), L-stat has consistently accurate pointwise CP.
At the most challenging points (smallest $x_0$), L-stat can under-cover by around five percentage points.
Otherwise, CP is near $1-\alpha$ for all $x_0$ in all DGPs.

In contrast, with the exception of boot, the other methods can have significant under-coverage.
As seen in the first two columns of Figure \ref{fig:ptCP-1}, rqss has under-coverage (as low as 50--60\% CP) for $x_0$ closer to zero.
QYg has under-coverage with the $\chi^2_3$ and (especially) Cauchy.
FLb has good CP except with the Cauchy, where CP can dip below 70\%.

Figure \ref{fig:ptCP-1}'s third column shows the joint power curves.
The horizontal axis of the graphs indicates the deviation of $H_0$ from the true values.  For example, letting $\xi_{p,j}$ be the true conditional quantiles at the $j=1,\ldots,47$ values of $x_0$ (say, $x_j$), $-0.1$ deviation refers to $H_0:\{Q_{Y|X}(p;x_j)=\xi_{p,j}-0.1\textrm{ for }j=1,\ldots,47\}$ (which is false), and zero deviation means $H_0$ is true.
Our method's type I error rate is close to $\alpha$ under all four $U_i$ distributions (5.7\%, 5.8\%, 7.3\%, 6.3\%).
In contrast, other methods show size distortion under Cauchy and/or $\chi^2_3$ $U_i$; among them, boot is closest but still has $10.3\%$ type I error rate with the $\chi^2_3$.
Next-best is rqss; size distortion for FLb and QYg is more serious.
L-stat also has the steepest joint power curves among all methods.
Beyond steepness, they are also the most robust to the underlying distribution.
L-stat's type I error rate is near 5\% for all four distributions.
In contrast, boot ranges from only 1.2\% for the Cauchy, leading to worse power, up to 10.3\% for the $\chi^2_3$.

The supplemental appendix shows a comparison of hypothesis tests based on uniform confidence bands.
The results are similar to the joint power curves, but with slightly higher rejection rates all around.

\begin{figure}[thb]
  \centering
  \hfill
  \includegraphics[clip=true,trim=15 15 10 58,width=0.445\textwidth]
    {qinfnp_2015_08_16_ptPWR2_F010_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=51 15 10 58,width=0.406\textwidth]
    {qinfnp_2015_08_16_ptPWR2_F012_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill\null
\\
  \hfill
  \includegraphics[clip=true,trim=15 15 10 58,width=0.445\textwidth]
    {qinfnp_2015_08_16_ptPWR2_F014_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill
  \includegraphics[clip=true,trim=51 15 10 58,width=0.406\textwidth]
    {qinfnp_2015_08_16_ptPWR2_F016_a05_n400_p50_numx47_trimx04_nrep1000.pdf}
  \hfill\null
  \caption{\label{fig:ptPWR2-1}Pointwise power (described in text), $1-\alpha=0.95$, $n=400$, $p=0.5$, DGP from \eqref{eqn:DGP-rqss}, $\sigma(x)=0.2$. The $U_i$ are $N(0,1)$ (top left), $t_3$ (top right), Cauchy (bottom left), and centered $\chi^2_3$ (bottom right).}
\end{figure}

Figure \ref{fig:ptPWR2-1} shows pointwise power.  Specifically, for a given $x_0$, this is the proportion of simulation draws in which $Q_{Y|X}(p;x_0)-0.1$ is excluded from the CI, averaged with the corresponding proportion for $Q_{Y|X}(p;x_0)+0.1$.
L-stat generally has the best power among methods with correct CP (per first column of Figure \ref{fig:ptCP-1}).

The supplemental appendix contains results for $p=0.25$, where L-stat continues to perform well.
One additional advantage is that L-stat's joint test is nearly unbiased, whereas the other joint tests are all biased.

The supplemental appendix also shows the computational advantage of our method.  For example, with $n=10^5$ and $100$ different $x_0$, L-stat takes only $10$ seconds, whereas the local cubic bootstrap takes $141$ seconds; rqss is even slower.

Overall, the simulation results show the new L-stat method to be fast and accurate.
Besides L-stat, the only method to avoid serious under-coverage is the local cubic with bootstrapped standard errors, perhaps due to its reliance on our newly proposed CPE-optimal bandwidth.  However, L-stat consistently has better power, greater robustness across different conditional distributions, and less bias of its joint hypothesis tests.















\section{Conclusion}

We derive a uniform $O(n^{-1})$ difference between the linearly interpolated and ideal fractional order statistic distributions.
We generalize this to $L$-statistics to help justify quantile inference procedures.
In particular, this translates to $O(n^{-1})$ CPE for the quantile CIs proposed by \citet{Hutson1999}, which we improve to $O\mathopen{}\mathclose\bgroup\originalleft(n^{-3/2}[\log(n)]^3\aftergroup\egroup\originalright)$ via calibration.
We extend these results to a nonparametric conditional quantile model, with both theoretical and Monte Carlo success.  The derivation of an optimal bandwidth value (not just rate) and a fast approximation thereof are important practical advantages.

Our results can be extended to other objects of interest, such as interquantile ranges and two-sample quantile differences \citep{GoldmanKaplan2014b}, quantile marginal effects \citep{Kaplan2014}, and entire distributions \citep{GoldmanKaplan2015c}.

In ongoing work, we consider the connection with Bayesian bootstrap quantile inference, which may be a way to ``relax'' the iid assumption.
Other future work may improve finite-sample performance, e.g.,\ by smoothing over discrete covariates \citep{LiRacine2007}.



\singlespacing
\iftoggle{SUPPLEMENTAL}{}{
\bibliographystyle{chicago}
\begin{thebibliography}{}

\bibitem[\protect\citeauthoryear{Abrevaya}{Abrevaya}{2001}]{Abrevaya2001}
Abrevaya, J. (2001).
\newblock The effects of demographics and maternal behavior on the distribution
  of birth outcomes.
\newblock {\em Empirical Economics\/}~{\em 26\/}(1), 247--257.

\bibitem[\protect\citeauthoryear{Alan, Crossley, Grootendorst, and Veall}{Alan
  et~al.}{2005}]{AlanEtAl2005}
Alan, S., T.~F. Crossley, P.~Grootendorst, and M.~R. Veall (2005).
\newblock Distributional effects of `general population' prescription drug
  programs in {Canada}.
\newblock {\em Canadian Journal of Economics\/}~{\em 38\/}(1), 128--148.

\bibitem[\protect\citeauthoryear{Andrews and Guggenberger}{Andrews and
  Guggenberger}{2014}]{AndrewsGuggenberger2014}
Andrews, D. W.~K. and P.~Guggenberger (2014).
\newblock A conditional-heteroskedasticity-robust confidence interval for the
  autoregressive parameter.
\newblock {\em Review of Economics and Statistics\/}~{\em 96\/}(2), 376--381.

\bibitem[\protect\citeauthoryear{Banks, Blundell, and Lewbel}{Banks
  et~al.}{1997}]{BanksEtAl1997}
Banks, J., R.~Blundell, and A.~Lewbel (1997).
\newblock Quadratic {Engel} curves and consumer demand.
\newblock {\em Review of Economics and Statistics\/}~{\em 79\/}(4), 527--539.

\bibitem[\protect\citeauthoryear{Beran and Hall}{Beran and
  Hall}{1993}]{BeranHall1993}
Beran, R. and P.~Hall (1993).
\newblock Interpolated nonparametric prediction intervals and confidence
  intervals.
\newblock {\em Journal of the Royal Statistical Society: Series B (Statistical
  Methodology)\/}~{\em 55\/}(3), 643--652.

\bibitem[\protect\citeauthoryear{Bhattacharya and Gangopadhyay}{Bhattacharya
  and Gangopadhyay}{1990}]{BhattacharyaGangopadhyay1990}
Bhattacharya, P.~K. and A.~K. Gangopadhyay (1990).
\newblock Kernel and nearest-neighbor estimation of a conditional quantile.
\newblock {\em Annals of Statistics\/}~{\em 18\/}(3), 1400--1415.

\bibitem[\protect\citeauthoryear{Bickel}{Bickel}{1967}]{Bickel1967}
Bickel, P.~J. (1967).
\newblock Some contributions to the theory of order statistics.
\newblock In {\em Proceedings of the Fifth Berkeley Symposium on Mathematical
  Statistics and Probability, Volume 1: Statistics}. The Regents of the
  University of California.

\bibitem[\protect\citeauthoryear{Buchinsky}{Buchinsky}{1994}]{Buchinsky1994}
Buchinsky, M. (1994).
\newblock Changes in the {U.S.} wage structure 1963--1987: Application of
  quantile regression.
\newblock {\em Econometrica\/}~{\em 62\/}(2), 405--458.

\bibitem[\protect\citeauthoryear{Chamberlain}{Chamberlain}{1994}]{Chamberlain1994}
Chamberlain, G. (1994).
\newblock Quantile regression, censoring, and the structure of wages.
\newblock In {\em Advances in Econometrics: Sixth World Congress}, Volume~2,
  pp.\  171--209.

\bibitem[\protect\citeauthoryear{Chaudhuri}{Chaudhuri}{1991}]{Chaudhuri1991}
Chaudhuri, P. (1991).
\newblock Nonparametric estimates of regression quantiles and their local
  {Bahadur} representation.
\newblock {\em Annals of Statistics\/}~{\em 19\/}(2), 760--777.

\bibitem[\protect\citeauthoryear{Chen and Hall}{Chen and
  Hall}{1993}]{ChenHall1993}
Chen, S.~X. and P.~Hall (1993).
\newblock Smoothed empirical likelihood confidence intervals for quantiles.
\newblock {\em Annals of Statistics\/}~{\em 21\/}(3), 1166--1181.

\bibitem[\protect\citeauthoryear{DasGupta}{DasGupta}{2000}]{DasGupta2000}
DasGupta, A. (2000).
\newblock Best constants in {Chebyshev} inequalities with various applications.
\newblock {\em Metrika\/}~{\em 51\/}(3), 185--200.

\bibitem[\protect\citeauthoryear{David and Nagaraja}{David and
  Nagaraja}{2003}]{DavidNagaraja2003}
David, H.~A. and H.~N. Nagaraja (2003).
\newblock {\em Order Statistics\/} (3rd ed.).
\newblock New York: Wiley.

\bibitem[\protect\citeauthoryear{Deaton}{Deaton}{1997}]{Deaton1997}
Deaton, A. (1997).
\newblock {\em The analysis of household surveys: a microeconometric approach
  to development policy}.
\newblock Baltimore: The Johns Hopkins University Press.

\bibitem[\protect\citeauthoryear{Donald, Hsu, and Barrett}{Donald
  et~al.}{2012}]{DonaldEtAl2012}
Donald, S.~G., Y.-C. Hsu, and G.~F. Barrett (2012).
\newblock Incorporating covariates in the measurement of welfare and
  inequality: methods and applications.
\newblock {\em The Econometrics Journal\/}~{\em 15\/}(1), C1--C30.

\bibitem[\protect\citeauthoryear{Engel}{Engel}{1857}]{Engel1857}
Engel, E. (1857).
\newblock Die productions- und consumtionsverh\"altnisse des k\"onigreichs
  sachsen.
\newblock {\em Zeitschrift des Statistischen Bureaus des K\"oniglich
  S\"achsischen, Ministerium des Inneren\/}~{\em 8--9}, 1--54.

\bibitem[\protect\citeauthoryear{Fan, Grama, and Liu}{Fan
  et~al.}{2012}]{FanEtAl2012}
Fan, X., I.~Grama, and Q.~Liu (2012).
\newblock Hoeffding's inequality for supermartingales.
\newblock {\em Stochastic Processes and their Applications\/}~{\em 122\/}(10),
  3545--3559.

\bibitem[\protect\citeauthoryear{Fan and Liu}{Fan and Liu}{2016}]{FanLiu2016}
Fan, Y. and R.~Liu (2016).
\newblock A direct approach to inference in nonparametric and semiparametric
  quantile models.
\newblock {\em Journal of Econometrics\/}~{\em 191\/}(1), 196--216.

\bibitem[\protect\citeauthoryear{Ferguson}{Ferguson}{1973}]{Ferguson1973}
Ferguson, T.~S. (1973).
\newblock A {Bayesian} analysis of some nonparametric problems.
\newblock {\em Annals of Statistics\/}~{\em 1\/}(2), 209--230.

\bibitem[\protect\citeauthoryear{Fisher}{Fisher}{1932}]{Fisher1932}
Fisher, R.~A. (1932).
\newblock {\em Statistical Methods for Research Workers\/} (4th ed.).
\newblock Edinburg: Oliver and Boyd.

\bibitem[\protect\citeauthoryear{Goldman and Kaplan}{Goldman and
  Kaplan}{2016a}]{GoldmanKaplan2015c}
Goldman, M. and D.~M. Kaplan (2016a).
\newblock Evenly sensitive {KS-type} inference on distributions.
\newblock Working paper, available at
  \url{http://faculty.missouri.edu/~kaplandm}.

\bibitem[\protect\citeauthoryear{Goldman and Kaplan}{Goldman and
  Kaplan}{2016b}]{GoldmanKaplan2014b}
Goldman, M. and D.~M. Kaplan (2016b).
\newblock Nonparametric inference on conditional quantile differences, linear
  combinations, and vectors, using {$L$-statistics}.
\newblock Working paper, available at
  \url{http://faculty.missouri.edu/~kaplandm}.

\bibitem[\protect\citeauthoryear{Hall and Sheather}{Hall and
  Sheather}{1988}]{HallSheather1988}
Hall, P. and S.~J. Sheather (1988).
\newblock On the distribution of a {Studentized} quantile.
\newblock {\em Journal of the Royal Statistical Society: Series B (Statistical
  Methodology)\/}~{\em 50\/}(3), 381--391.

\bibitem[\protect\citeauthoryear{Ho and Lee}{Ho and Lee}{2005a}]{HoLee2005a}
Ho, Y. H.~S. and S.~M.~S. Lee (2005a).
\newblock Calibrated interpolated confidence intervals for population
  quantiles.
\newblock {\em Biometrika\/}~{\em 92\/}(1), 234--241.

\bibitem[\protect\citeauthoryear{Ho and Lee}{Ho and Lee}{2005b}]{HoLee2005b}
Ho, Y. H.~S. and S.~M.~S. Lee (2005b).
\newblock Iterated smoothed bootstrap confidence intervals for population
  quantiles.
\newblock {\em Annals of Statistics\/}~{\em 33\/}(1), 437--462.

\bibitem[\protect\citeauthoryear{Hogg}{Hogg}{1975}]{Hogg1975}
Hogg, R. (1975).
\newblock Estimates of percentile regression lines using salary data.
\newblock {\em Journal of the American Statistical Association\/}~{\em
  70\/}(349), 56--59.

\bibitem[\protect\citeauthoryear{Horowitz and Lee}{Horowitz and
  Lee}{2012}]{HorowitzLee2012}
Horowitz, J.~L. and S.~Lee (2012).
\newblock Uniform confidence bands for functions estimated nonparametrically
  with instrumental variables.
\newblock {\em Journal of Econometrics\/}~{\em 168\/}(2), 175--188.

\bibitem[\protect\citeauthoryear{Hotelling}{Hotelling}{1939}]{Hotelling1939}
Hotelling, H. (1939).
\newblock Tubes and spheres in $n$-space and a class of statistical problems.
\newblock {\em American Journal of Mathematics\/}~{\em 61}, 440--460.

\bibitem[\protect\citeauthoryear{Hutson}{Hutson}{1999}]{Hutson1999}
Hutson, A.~D. (1999).
\newblock Calculating nonparametric confidence intervals for quantiles using
  fractional order statistics.
\newblock {\em Journal of Applied Statistics\/}~{\em 26\/}(3), 343--353.

\bibitem[\protect\citeauthoryear{Jones}{Jones}{2002}]{Jones2002}
Jones, M.~C. (2002).
\newblock On fractional uniform order statistics.
\newblock {\em Statistics \& Probability Letters\/}~{\em 58\/}(1), 93--96.

\bibitem[\protect\citeauthoryear{Kaplan}{Kaplan}{2014}]{Kaplan2014}
Kaplan, D.~M. (2014).
\newblock Nonparametric inference on quantile marginal effects.
\newblock Working paper, available at
  \url{http://faculty.missouri.edu/~kaplandm}.

\bibitem[\protect\citeauthoryear{Kaplan}{Kaplan}{2015}]{Kaplan2015}
Kaplan, D.~M. (2015).
\newblock Improved quantile inference via fixed-smoothing asymptotics and
  {Edgeworth} expansion.
\newblock {\em Journal of Econometrics\/}~{\em 185\/}(1), 20--32.

\bibitem[\protect\citeauthoryear{Kaplan and Sun}{Kaplan and
  Sun}{2016}]{KaplanSun2016}
Kaplan, D.~M. and Y.~Sun (2016).
\newblock Smoothed estimating equations for instrumental variables quantile
  regression.
\newblock {\em Econometric Theory\/}~{\em XX\/}(XX), XX--XX.
\newblock Forthcoming.

\bibitem[\protect\citeauthoryear{Koenker}{Koenker}{2012}]{R.quantreg}
Koenker, R. (2012).
\newblock {\em quantreg: Quantile Regression}.
\newblock R package version 4.81.

\bibitem[\protect\citeauthoryear{Kumaraswamy}{Kumaraswamy}{1980}]{Kumaraswamy1980}
Kumaraswamy, P. (1980).
\newblock A generalized probability density function for double-bounded random
  processes.
\newblock {\em Journal of Hydrology\/}~{\em 46\/}(1--2), 79--88.

\bibitem[\protect\citeauthoryear{Lehmann}{Lehmann}{1951}]{Lehmann1951}
Lehmann, E.~L. (1951).
\newblock A general concept of unbiasedness.
\newblock {\em Annals of Mathematical Statistics\/}~{\em 22\/}(4), 587--592.

\bibitem[\protect\citeauthoryear{Li and Racine}{Li and
  Racine}{2007}]{LiRacine2007}
Li, Q. and J.~S. Racine (2007).
\newblock {\em Nonparametric econometrics: Theory and practice}.
\newblock Princeton University Press.

\bibitem[\protect\citeauthoryear{Manning, Blumberg, and Moulton}{Manning
  et~al.}{1995}]{ManningEtAl1995}
Manning, W., L.~Blumberg, and L.~Moulton (1995).
\newblock The demand for alcohol: the differential response to price.
\newblock {\em Journal of Health Economics\/}~{\em 14\/}(2), 123--148.

\bibitem[\protect\citeauthoryear{Muir}{Muir}{1960}]{Muir1960}
Muir, T. (1960).
\newblock {\em A Treatise on the Theory of Determinants}.
\newblock Dover Publications.

\bibitem[\protect\citeauthoryear{Neyman}{Neyman}{1937}]{Neyman1937}
Neyman, J. (1937).
\newblock {\guillemotright}{Smooth} test{\guillemotright} for goodness of fit.
\newblock {\em Skandinavisk Aktuarietidskrift\/}~{\em 20\/}(3--4), 149--199.

\bibitem[\protect\citeauthoryear{{Office for National Statistics and Department
  for Environment, Food and Rural Affairs}}{{Office for National Statistics and
  Department for Environment, Food and Rural Affairs}}{2012}]{UKLCFS2012}
{Office for National Statistics and Department for Environment, Food and Rural
  Affairs} (2012).
\newblock {Living Costs and Food Survey}.
\newblock {2nd Edition. Colchester, Essex: UK Data Archive.
  \url{http://dx.doi.org/10.5255/UKDA-SN-7472-2}}.

\bibitem[\protect\citeauthoryear{Pearson}{Pearson}{1933}]{Pearson1933}
Pearson, K. (1933).
\newblock On a method of determining whether a sample of size n supposed to
  have been drawn from a parent population having a known probability integral
  has probably been drawn at random.
\newblock {\em Biometrika\/}~{\em 25}, 379--410.

\bibitem[\protect\citeauthoryear{Peizer and Pratt}{Peizer and
  Pratt}{1968}]{PeizerPratt1968}
Peizer, D.~B. and J.~W. Pratt (1968).
\newblock A normal approximation for binomial, {$F$}, beta, and other common,
  related tail probabilities, {I}.
\newblock {\em Journal of the American Statistical Association\/}~{\em
  63\/}(324), 1416--1456.

\bibitem[\protect\citeauthoryear{Polansky and Schucany}{Polansky and
  Schucany}{1997}]{PolanskySchucany1997}
Polansky, A.~M. and W.~R. Schucany (1997).
\newblock Kernel smoothing to improve bootstrap confidence intervals.
\newblock {\em Journal of the Royal Statistical Society: Series B (Statistical
  Methodology)\/}~{\em 59\/}(4), 821--838.

\bibitem[\protect\citeauthoryear{Polonik and Yao}{Polonik and
  Yao}{2002}]{PolonikYao2002}
Polonik, W. and Q.~Yao (2002).
\newblock Set-indexed conditional empirical and quantile processes based on
  dependent data.
\newblock {\em Journal of Multivariate Analysis\/}~{\em 80\/}(2), 234--255.

\bibitem[\protect\citeauthoryear{Pratt}{Pratt}{1968}]{Pratt1968}
Pratt, J.~W. (1968).
\newblock A normal approximation for binomial, {$F$}, beta, and other common,
  related tail probabilities, {II}.
\newblock {\em Journal of the American Statistical Association\/}~{\em
  63\/}(324), 1457--1483.

\bibitem[\protect\citeauthoryear{Qu and Yoon}{Qu and Yoon}{2015}]{QuYoon2015}
Qu, Z. and J.~Yoon (2015).
\newblock Nonparametric estimation and inference on conditional quantile
  processes.
\newblock {\em Journal of Econometrics\/}~{\em 185\/}(1), 1--19.

\bibitem[\protect\citeauthoryear{R{\'e}nyi}{R{\'e}nyi}{1953}]{Renyi1953}
R{\'e}nyi, A. (1953).
\newblock On the theory of order statistics.
\newblock {\em Acta Mathematica Hungarica\/}~{\em 4\/}(3), 191--231.

\bibitem[\protect\citeauthoryear{Robbins}{Robbins}{1955}]{Robbins1955}
Robbins, H. (1955).
\newblock A remark on {Stirling's} formula.
\newblock {\em The American Mathematical Monthly\/}~{\em 62\/}(1), 26--29.

\bibitem[\protect\citeauthoryear{Ruppert, Wand, and Carroll}{Ruppert
  et~al.}{2003}]{RuppertEtAl2003}
Ruppert, D., M.~P. Wand, and R.~J. Carroll (2003).
\newblock {\em Semiparametric Regression}.
\newblock Cambridge Series in Statistical and Probabilistic Mathematics.
  Cambridge University Press.

\bibitem[\protect\citeauthoryear{Shorack}{Shorack}{1972}]{Shorack1972}
Shorack, G.~R. (1972).
\newblock Convergence of quantile and spacings processes with applications.
\newblock {\em Annals of Mathematical Statistics\/}~{\em 43\/}(5), 1400--1411.

\bibitem[\protect\citeauthoryear{Shorack and Wellner}{Shorack and
  Wellner}{1986}]{ShorackWellner1986}
Shorack, G.~R. and J.~A. Wellner (1986).
\newblock {\em Empirical Processes with Applications to Statistics}.
\newblock New York: John Wiley \& Sons.

\bibitem[\protect\citeauthoryear{Stigler}{Stigler}{1977}]{Stigler1977}
Stigler, S.~M. (1977).
\newblock Fractional order statistics, with applications.
\newblock {\em Journal of the American Statistical Association\/}~{\em
  72\/}(359), 544--550.

\bibitem[\protect\citeauthoryear{Thompson}{Thompson}{1936}]{Thompson1936}
Thompson, W.~R. (1936).
\newblock On confidence ranges for the median and other expectation
  distributions for populations of unknown distribution form.
\newblock {\em Annals of Mathematical Statistics\/}~{\em 7\/}(3), 122--128.

\bibitem[\protect\citeauthoryear{Wilks}{Wilks}{1962}]{Wilks1962}
Wilks, S.~S. (1962).
\newblock {\em Mathematical Statistics}.
\newblock New York: Wiley.

\end{thebibliography}

}


\doublespacing
\onehalfspacing
\singlespacing




















\iftoggle{SUPPLEMENTAL}{\pagebreak\setcounter{page}{2}}{}