In this paper, we consider estimating parameters defined as the (unique) solution to the following d-dimensional estimating equation:

where are observed random variables with distribution unknown nuisance functions, and 0 is the zero vector in . We hope to estimate (independent and identically distributed (i.i.d) draws from the distribution P. As we will show, estimating equations of the form above are prevalent in efficient estimation in causal inference and missing data problems, with quantile treatment effect (QTE) estimation (Section 1.1) as a prominent example, among many others (Section 1.2).

One important feature of Eq. (1) is that the nuisance ) depends on the parameters to be estimated, which raises several challenges and causes existing methods to be unstable and computationally burdensome. Specifically, we could potentially use the observed data to estimate the nuisances ) and then solve a sample analogue of Eq. (1) based on the estimated nuisances in order to estimate , possibly using cross-fitting (Chernozhukov et al. 2018a). However, this requires estimating the nuisance functions of Z, and then solving for the root of an estimated function. For example, when estimating QTE (see Section 1.1), this involves estimating a whole conditional cumulative distribution function, or equivalently, infinitely many binary probability regressions (one for each threshold). This can be very unstable, especially in causal inference with observational data where typically a large number of covariates need to be conditioned on to remove confounding. Although one may discretize the space of Θand estimate ) only for finitely many , this can still be compu- tationally burdensome when the discrete grid is large, and the resulting estimator can be sensitive to the discretization scheme.

In this paper, we propose a localized debiased machine learning (LDML) approach that only requires estimating value, without estimating it for all possible values or ad-hoc discretized values of . Importantly, our estimator is asymptotically equivalent to an oracle estimator that knows the whole continuum of nuisance function other words, asymptotically, our method does not incur any loss even though it only estimates the nuisance function at a single value. Moreover, estimating this far simpler nuisance reduces to standard classification and regression tasks, i.e., fitting conditional expectations (regression) and conditional binary probabilities (classification), for which many machine learning methods exist. In particular, our approach will be shown to be largely insensitive to how these conditional expectation functions are estimated, so we may directly use off-the-shelf machine learning methods and treat them as black-box regression or classification algorithms (e.g., random forests, gradient boosting, neural networks). Therefore, our proposed method notably enables practical and efficient estimation using time-tested machine learning methods to solve Eq. (1).

In comparison, existing approaches for debiased and efficient estimation with black-box nuisance estimators either focus on settings where nuisances do not depend on target parameters or treat nuisances as abstract objects so that one must estimate a continuum of nuisances when applying to Eq. (1) thus precluding the use of standard machine-learning algorithms for regression and classification (Robins et al. 2008, Zheng and van der Laan 2011, Robins et al. 2013, Chernozhukov et al. 2018a, Bravo et al. 2020). Similarly, existing works specifically on the efficient estimation of QTEs either apply similar debiased approaches using a continuum of nuisances (Belloni et al. 2017, D´ıaz 2017) or use specific non-black-box nuisance estimators like polynomial sieves and local polynomial kernel regression and make explicit smoothness restrictions (Firpo 2007, Fr¨olich and Melly 2013). (We provide an extensive literature review in Section 7.) Compared to these works, our proposal is fully generic, flexible, and machine-learning driven in that it handles many important examples that fit into Eq. (1), as we review in the next two sections; it permits the use of flexible black-box nuisance estimators, since we only require lax rate conditions that are in fact more lax than some previous results; and these black boxes may be standard machine-learning methods for regression and classification, since whenever ) for any single is a conditional expectation, such as in all of the examples we review in the next section, our method only ever has to fit very few conditional expectations.

1.1 Motivating Example: Quantile Treatment Effects

A primary motivation of considering the setting of Eq. (1) is the estimation of QTE. In this case, we consider a population P of units, each associated with some baseline covariates , two potential outcomes for each of two possible treatments, and a treatment indicator Since both potential outcomes are included in this description, we refer to P as the complete-data distribution. We are interested in the -quantile of (assuming existence and uniqueness) for 1). And, similarly, we are interested in the quantile of Y (0) and in the difference of the quantiles, known as the quantile treatment effect (QTE), but these estimation questions are analogous so for brevity we focus just on -quantile of Y (1) (see also Remark 2). Compared to the average outcome and the average treatment effect (ATE), the quantile of outcomes and the QTE provide a more robust assessment of the effects of treatment and are very important quantities in program evaluation.

We do not observe the potential outcomes but instead only the realized factual outcome corresponding to the assigned treatment, Y = Y (T). Therefore, we only observe Z = (X, T, Y ), whose distribution P is given by coarsening P via Y = Y (T). Ignorable treatment assignment with respect to X assumes that , no unobserved confounders) and overlap assumes that P(T = 1 | X) > 0, and these together ensure that is identifiable from observations of Z. Specifically, a straightforward identification is given by the so-called inverse propensity weighting (IPW) equation:

Here estimating the propensity score function amounts to learning a conditional probabil- ity function from a binary response, for which many standard machine learning methods exist. Once we construct an estimator ˆ, we can obtain the standard IPW estimator ˆ

the particular method used to construct ˆand its convergence rate can be slowed down by that of ˆ, prohibiting the use of general nonparametric machine learning methods and potentially leading to unstable estimates.

Instead, one can alternatively obtain the following estimating equation from the efficient influence function for , Tsiatis 2006):

An important feature of the above is that it satisfies a property known as Neyman orthogonality: the moment derivatives with respect to the nuisances at This means that the estimating equation is robust to small perturbations in the nuisances so that estimation errors therein contribute only to higher-order error terms in the final estimate of In particular, Chernozhukov et al. (2018a) recently proposed to leverage Neyman orthogonality to enable the use of plug-in machine learning estimates of the nuisances. Their proposal, called debiased machine learning (DML), is as follows: split the data randomly into then for each k = 1, . . . , K, use all but the fold to construct nuisance estimates ˆfinally solve the empirical estimating equation to obtain the estimator ˆ. They prove that as long as the estimates ˆconverge to faster than , the estimate ˆwill have similar behavior to the oracle estimate that solves

, solving the empirical estimating equation using the true nuisance functions. As a result, the estimate ˆis asymptotically normal and semiparametrically efficient. Since, apart from the mild rate requirement on ˆ, no metric entropy condi- tions are assumed, this allows one to successfully use machine learning methods to learn nuisances and achieve asymptotically normal and efficient estimation.

The problem with this approach for estimating quantiles of outcomes (similarly, QTEs), however, is that it requires the estimation of a very complex nuisance function: conditional cumulative distribution function of a real-valued outcome, potentially conditioned on high-dimensional covariates. While certainly nonparametric methods for estimating conditional distributions exist such as kernel estimators, this learning problem is much harder to do in a flexi-ble, blackbox, machine-learning manner, compared to just estimating a single regression function. This indeed stands in stark contrast to the estimation of ATEs, where applying DML requires a far simpler nuisance function given by the regression of outcome on covariates and treatment, E [Y | X, T], for which a long list of practice-proven machine learning methods can be directly and successfully applied. The key difference is that the nuisance function in ATE estimation does not depend on the estimand and can therefore be estimated in an independent manner whereas the nuisance function in QTE estimation does depend on the estimand. This issue makes DML, despite its theoretical benefits, untenable in practice for the important task of QTE estimation.

The primary goal of this paper can be understood as extending DML to effectively tackle the case where nuisances depend on the estimand by alleviating this dependence via localization. In particular, this will enable efficient estimation of important quantities such as QTEs in the presence of high-dimensional nuisances by using and debiasing black-box machine learning methods for the standard regression task.

The basic idea as it applies to the estimation of the quantile of outcomes, which we will generalize and analyze thoroughly, is as follows. While perhaps inefficient, ˆrelies only on estimating a binary regression . This is amenable to machine learning approaches but may have a slow convergence rate in general. Despite its slow rate, this rough initial guess can sufficiently localize our nuisance estimation and it may suffice to only estimate , the nuisance evaluated at just a single initial estimate of . Then we use this estimated nuisance at this initial estimate of in place of ) when solving the empirical estimating equation for . For estimating the quantiles, this means we only have to regress the binary response ˆas fixed. In particular, we propose a special three-way data splitting procedure that debiases such plug-in nuisance estimates in order to obtain an estimate for with near-oracle performance.

1.2 Estimating Equations with Incomplete Data under Ignorable Treatment Assignment or Using Instrumental Variables

More generally, we can consider parameters () defined as the solution to the following esti- mating equation on the (unavailable) complete data:

for some given functions ). Quantile is one example of this. Another example is conditional value at risk (CVaR) of outcomes: the cumulative distribution function of Y (1), that is, the expectation of Y (1) conditioned on being above the -quantile (again, assuming uniqueness). CVaR is also known as expected shortfall, a popular risk measure widely used in risk management and optimization (Rockafellar and Uryasev 2002). Again, we may consider the CVaR of Y (0) and the differences of CVaRs analogously. Letting

Yet another example is the expectile, a measure for asymmetric risk (Newey and Powell 1987). The -expectile of Y (1) is defined by the following asymmetric least squares problem:

whose first-order condition gives an estimating equation for complete data:

Under ignorable treatment assignment and overlap, a general-purpose Neyman-orthogonal estimating equation for the estimand () defined by Eq. (4) is given by

Alternatively, instead of assuming ignorable treatment assignment, we may have access to an instrumental variable (IV). We considier a binary IV denoted as and assume that it satisfies identification conditions in Imbens and Angrist (1994) (namely, for potential treatments T(w) and potential outcomes Y (t, w), we have exclusion ), exogeneity (1), relevance P(T(1) = 1) > P(T(0) = 1), and monotonicity (0)). We seek to use observations of Z = (X, W, T, Y ) to estimate local parameters defined by the following estimating equation conditionally on the subpopulation of compliers (i.e., T(1) > T(0)):

For example, specializing Eq. (8) to the functions ) in Eq. (5) gives the local quantile and CVaR, which in turn gives the local QTE (LQTE). In Appendix A, we present the efficient estimating equations for these local parameters and show they also satisfy Neyman orthogonality and involve some estimand-dependent nuisance functions

In all examples above, the nuisance ) depends on the estimand. This occurs whether estimating quantiles, CVaR, or expectiles (more generally, whenever ) is not linear in and whether the identification is via ignorable treatment assignment, ignorable coarsening, or valid IV. And, in such cases, learning is practically difficult, which may involve learning a whole conditional distribution function or a whole continuum of conditional expectation functions given potentially high-dimensional covariates.

be the dimensions of , respectively, where -matrix-valued function with entry in position (i, j) and is its evaluation at -matrix-valued function with entry in position ()) to denote its largest singular value. We let ] for measurable sets A denote probabilities and expectations with respect to for measurable functions f denote expectations with respect to Z alone, while we let ) denote expectations with respect to Z and the data. Thus, if ˆdepends on the data, ) remains a function of the data while is a number. We let denote the empirical expectation: measurable function f. Moreover, for vector-valued function )). For any , we denote the open ball centered at x with radius 0 and a probability measure Q, we denote For a set of functions F, we define the covering number ) as the minimal number N of functions such that sup. For positive deterministic sequence and random variable sequence ) means for any 0, there exists M > 0 such that lim sup

We next present our methodology, first motivating the localization technique, and then explicitly stating our meta-algorithm.

2.1 Motivation

Ideally, if the nuisances were both known, then Eq. (1) suggests that could be estimated by solving the following estimating equation:

Under standard regularity conditions for Z-estimation (van der Vaart 1998), the resulting oracle estimator ˜that solves Eq. (9) is asymptotically linear (and hence-consistent and asymptotically normal):

Furthermore, if )) is the semiparametrically efficient influence function for also achieves the efficiency lower bound, that is, has minimal asymptotic variance among all regular estimators (van der Vaart 1998).

Since are actually unknown, the oracle estimator ˜is of course infeasible. Instead, we must estimate the nuisance functions. A direct application of DML would require us to learn the whole functions . That is, in order to solve Eq. (9) we would need to estimate infinitely many nuisance functions,

To avoid the daunting task of estimating infinitely many nuisances, we will instead attempt to target the following alternative oracle estimating equation

Although Eq. (11) appears very similar to Eq. (9), it only involves , as opposed to the infinitely many possible values for . In other words, among the whole family of nuisances is relevant for Eq. (11). This formulation considerably reduces the need of nuisance estimation: now we only need to estimate functions only of Z.

Figure 1: Sketch of the LDML estimation procedure and a possible initial guess estimator. Squiggly arrows “” denote estimation. Plain arrows “” denote plugging in.

The (infeasible) estimators that solve each of Eqs. (9) and (11) have the same leading asymptotic behavior as long as the respective associated Jacobian matrices coincide, as posited by the following assumption.

Assumption 1 (Invariant Jacobian)

In Appendix F, we provide a general sufficient condition for Assumption 1 in terms of a Fr´echet-differentiation variant of the Neyman orthogonality condition (Assumption 2 condition vii.). But it is easy to directly show that this invariant Jacobian assumption holds for estimating equations with incomplete data presented in Section 1.2. In particular, the estimating equation in Eq. (7) satisfies that

which does not depend on ) at all. Thus whether fixing or not, the Jacobian matrix of the estimating equation remains the same. This means that solving Eq. (9) or Eq. (11) will have the same asymptotic behavior. Both, however, are infeasible since they involve unknown nuisances. Nonetheless, Eq. (11) motivates a new algorithm that eschews estimating

2.2 The LDML Meta-Algorithm

Motivated by the new (infeasible) estimating equation in Eq. (11), we propose to estimate by the following (feasible) three-way sample splitting method, which we term localized debiased machine learning (LDML). The algorithm has two parts: three-way-cross-fold nuisance estimation and solving the estimating equation.

We start by discussing how we estimate the nuisances that we will then plug into Eq. (11).

Definition 1 (3-way-cross-fold nuisance estimation)

1. Randomly permute the data indices and let 1, . . . , K be a random even K-fold split of the data.

2. For k = 1, . . . , K:

For illustration the first iteration of step 2 above is sketched in Fig. 1(a) along with the plugging of estimated nuisances into the estimating equation (see Definitions 2 and 5 below).

Notice that since are disjoint, ) is a fixed, nonrandom function with respect to the data . That is, the nuisance estimation task in step 2b appears as estimating a rather than the estimation of

A natural question is, what might be a reasonable initial estimator. In the examples given in Sections 1.1 and 1.2, we can use an IPW estimate for ˆ(see Fig. 1(b) and Definition 4).

Given these nuisance estimates, we can obtain the LDML estimator for by approximately solving the average of the estimate of Eq. (11) in each fold.

Definition 2 (LDML). We let the estimator ˆbe given by (approximately) solving

In fact, we allow for an approximate least-squares solution, which is useful if the empirical estimating equation has no exact solution. Namely, we let ˆbe any satisfying

In Appendix D Definition 5, we give an alternative LDML estimator obtained by averaging solutions to Eq. (11) estimated in each fold separately. These two LDML estimators are asymptotically equivalent and all results in this paper apply to both, thus we focus on Definition 2 in the main text. Moreover, both estimators depend on the random splitting in Definition 1. To reduce the variance from this, we may aggregate estimates from multiple different sample splitting realizations. See Appendix E for a detailed discussion.

In this section, we provide the sufficient conditions that guarantee the proposed estimator ˆDefinitions 2 and 5 to be consistent and asymptotically normal. In particular, although the proposed estimator relies on plug-in nuisance estimators, it is asymptotically equivalent to the infeasible estimator based on Eq. (9) with true nuisances, that is, it satisfies Eq. (10). While some of our conditions are analogous to those in Chernozhukov et al. (2018a), some are not and our proof takes a different approach that enables weaker conditions for convergence rates of the nuisance estimators. Our asymptotic normality results may be stated uniformly over a sequence of models data generating distribution . Our first set of assumptions ensure that is reasonably identified by the given estimating equation for all . We also assume that our estimating equation satisfies the Neyman orthogonality condition with respect to a nuisance realization set that contains the nuisance estimates ˆ) with high probability. Note the set consists of pairs of functions of the data . Therefore, we denote members of the set as () is simply understood as a symbol representing of some fixed function of Z alone.

Assumption 2 (Regularity of Estimating Equations). Assume there exist positive constants such that the following conditions hold for all

vi. The nuisance realization set contains the true nuisance parameters (Moreover, the parameter space Θ is bounded and for each (, the function class is suitably measurable and its uniform covering entropy satisfies the following condition: for positive constants a, measurable envelope for that satisfies

(η

Assumption 2 conditions i.–v. constitute standard identification and regularity conditions for Z-estimation (with uniform guarantees; see also Remark 1 below). Assumption 2 condition vi. requires that is a well-estimable function of set of nuisances. Importantly, while it imposes a metric entropy condition on , this condition does not impose metric entropy conditions on our nuisance estimators, so flexible machine learning nuisance estimators are allowed. This assumption is very mild as Θ is finite-dimensional, so it can be ensured by some continuity and compactness condition. Finally, Assumption 2 condition vii. is the Neyman orthogonality condition (Cher- nozhukov et al. 2018a). We will show how these conditions are ensured in the incomplete data setting in Section 1.2.

Our second set of assumptions involve conditions on our nuisance estimators.

Assumption 3 (Nuisance Estimation Conditions)

i. For some sequence of constants ∆, the nuisance estimates (ˆto the realization set with probability at least 1

Here our condition on ) differs from the counterpart condition in Chernozhukov et al. (2018a), which also leads to a different proof strategy. Our condition and proof generally requires weaker conditions for convergence rates of nuisance estimators. See the discussions in Appendix H for more details. Moreover, the constants ∆are all prespecified and do not depend on any particular instance P.

Our key result in this paper is the following theorem, which shows that the asymptotic distribution of our estimator is identical to the (infeasible) oracle estimator solving the estimating equation in Eq. (9) with known nuisances.

Theorem 1 (Asymptotic Behavior of LDML). Assume Assumptions 1 to 3 hold with

Then the estimator ˆgiven in Definition 2 is asymptotically linear and converges to a Gaussian distribution uniformly over

The conditions in Eq. (16) and are fairly mild because Assumption 2 condition vi. requires 0) and Assumption 3 condition ii. requires

Remark 1 (Uniform vs non-uniform convergence). To obtain a non-uniform convergence result, we need only need set as a constant singleton in Theorem 1. In this case, much of Assumption 2 simplifies: the existence of the constants is trivial, the non-singularity of enough for to exist, and being in the interior of Θ is enough for to exist. Further, we can relax condition iv. by allowing to be zero (in which case we rephrase the asymptotic normality in Theorem 1 by putting Σ on the right-hand side of the limit rather than inverting it). Uniformity, however, is important in practice. Without uniformity, for any given sample size N there may always exist some bad instance such that the normal approximation suggested by the convergence is inaccurate (Kasy 2019).

In the previous section we established the asymptotic normality of the LDML estimator under lax conditions. This suggests that if we can estimate its asymptotic variance, then we can easily construct confidence intervals on . In this section we provide a variance estimator and prove its consistency, resulting in asymptotically calibrated confidence intervals. For DML, Chernozhukov et al. (2018a) provides variance estimates only for estimating functions that are linear in already excludes estimand-dependent nuisances. Our results are therefore notable both for handling nonlinear and non-differentiable estimating equations and for handling estimand-dependent nuisances.

Definition 3 (LDML variance estimator)from Definition 2 and ˆJ, set

We next establish the consistency of Σ, which relies on the following assumption.

Assumption 4. Assume that and that for some

Here, Eq. (17) implies continuity )) in terms of norm in the range space. Note that this condition can be satisfied even if differentiable. For example, in the estimation of QTEs, the efficient estimating equation in Eq. (3) involves the indicator function ], so the map )) is obviously not differentiable. However, the condition in Eq. (17) amounts to

In Assumption 5, we will assume that for a positive constant condition above follows if the cumulative distribution function of Y (1) is smooth enough, so that

Under Assumption 4, we now show that the variance estimator in Definition 3 is consistent and it leads to asymptotically valid confidence intervals.

Theorem 2. Assume the assumptions in Theorem 1 and Assumption 4. Then,

Given some , the confidence interval

In Definition 3, we assumed that we have a consistent estimator ˆ. How to construct such an estimator depends on the problem. When )) is differentiable, an estimator may easily be constructed as follows:

However, the estimating equation for QTE is not differentiable. Thus we rely on deriving the form of and estimate it directly, which we discuss in detail in Remark 4.

With finite sample, the variance of the LDML estimator also depends on the uncertain sample splitting in Definition 1. This uncertainty can be additionally accounted for when multiple sample splitting realizations are used, which we discuss in Appendix E.

Remark 2 (Estimating and Conducting Inference on Treatment Effects). Suppose we have two sets of parameters, , each identified by its own estimating equation, are interested in estimating the difference, . For example, quantile and/or CVaR of Y (0), Y (1), respectively, and we are interested in the QTE and/or CVaR treatment effect. To do this, we can concatenate the two estimating equations and augment them with the additional equation with LDML is equivalent to applying LDML to each of and letting ˆbe the difference of the estimates ˆ, where we may use the same data and the same folds for the two LDML procedures. For QTE and for other estimating equations with incomplete data, we can even share the nuisance estimates of the propensity score (in the below equation). The variance estimate one would derive for ˆfrom the augmented estimating equations is equivalent to

In this section, we apply our method and theory to general estimating equations with incomplete data presented in Eq. (4), which subsumes the estimation of QTEs, quantile of potential outcomes, CVaR treatment effect, CVaR of potential outcomes, expectile treatment effect, and expectile of potential outcomes. We will proceed to further specialize these results to quantile and CVaR estimation, deferring the case of expectiles to the appendix (Appendix B). We also defer the case of using IVs to estimate the solution to local estimating equations, such as those that describe the LQTE, to appendix (Appendix A).

As motivated in Section 1.1, under unconfoundedness, there is a very natural initial estimator: the IPW estimator. As we will show, the LDML estimate for this problem using the IPW initial estimator can be computed using just blackbox algorithms for (possibly binary) regression, which is the standard supervised learning task in machine learning. And, under lax conditions, the estimate is efficient, asymptotically normal, and amenable to inference.

Recall that is defined by the complete-data estimating equations in Eq. (4), namely, is identified from the incomplete-data observations Z = (X, T, Y ) where Y = Y (T). In particular, Eq. (7) provides a Neyman-orthogonal estimating equation identifying . For better interpretability, we give our nuisances names: we denote [For estimating parameters corresponding to Y (1), our estimandindependent nuisance is the propensity score ), and our estimand-dependent nuisance is ). The case for Y (0) is symmetric; and it also need the symmet- ric ignorability and overlap assumptions for identifiability: Treatment effects (e.g., QTEs) can be estimated by differences of estimates, where we can use the same data, the same fold splits, and the same estimates of for both treatments (see Remark 2).

This problem also admits a simpler but unstable (i.e., non-orthogonal) estimating equation using IPW, which suggests a possible initial estimator, using 2 in Definition 1:

Definition 4 (IPW Initial Estimator)as in Definition 1, use only the data in to construct a propensity score estimator ˆbe given by solving the following estimating equation (or, its least squares solution up to approximation error of

This procedure is illustrated in Fig. 1(b). Note that, given a fixed are conditional expectations of observable variables given X. Thus, in this setting, the whole LDML estimate using the IPW initial estimate can be computed given just blackbox algorithms for (possibly binary) regression.

5.1 Theoretical Analysis

We first study the LDML estimate for estimating equations with incomplete data by leveraging our general theory in Theorem 1. To this end, we assume a strong form of the overlap condition and specify the convergence rates of the initial estimator and nuisance estimators used. We consider a generic treatment level in these two assumptions.

Assumption 5 (Strong Overlap). Assume that there exists a positive constant for any almost surely.

Assumption 6 (Nuisance Estimation Rates). Assume that for any : condition i. of Assumption 3 holds for a sequence of constants ∆; with probability at least 1 for almost all realizations of X, and

The following theorem establishes that the asymptotic distribution of our proposed estimator is similar to the (infeasible) one that solves the semiparametric efficient estimating equation in Eq. (7) with known nuisances. This theorem is proved by verifying conditions in Theorem 1, namely Assumptions 1 to 3.

Theorem 3 (LDML for Estimating Equations with Incomplete Data). Fix t = 1 and let the estimator ˆbe given by applying Definition 2 to the estimating equation in Eq. (7). Suppose Assumptions 5 and 6 hold and that there exist positive constants such that for any the following conditions hold:

i. Conditions i. (with ), ii., v. (with ), and vi. (with ) of Assumption 2 and condition iii. of Assumption 3 for the estimating equation in Eq. (7).

ii. For is differentiable at any in a compact set Θ, and each component of its gradient is -Lipschitz continuous at . Moreover, for any with

σ

Then ˆsatisfies the conclusion of Theorem 1 for given in Eq. (7).

An analogous result for the estimating equations involving Y (0) holds when we change t = 1 to t = 0 everywhere in Theorem 3. See Remark 2 regarding estimation of the difference of the parameters (i.e., the treatment effects) and inference thereon.

In Theorem 3, conditions ii. and iii. guarantee the identification conditions iii. and iv. of Assumption 2. Condition iv. enables exchange of integration, which together with conditions v., vi., and vii. imply the rate condition ii. of Assumption 3. Note condition vii. permits nonparametric rates for nuisance estimators. Focusing on the order in the sample size and up to polylog factors, the condition allows for Note the first two restrictions are on products, permitting a trade-off between rates for different nuisances (see also Appendix H).

Remark 3 (Rate Conditions with IPW Initial Estimator). In Appendix C, we prove that if the propensity nuisance estimators used to construct the IPW initial estimators (Definition 4) also have convergence rate , then the initial estimators’ convergence rates satisfy that In this case, we are essentially imposing ): condition vii. of Theorem 3 requires ), so unless is somehow even faster than , we must need both and

5.2 Quantile and CVaR

Now we consider estimating quantile and (possibly) CVaR based on the semiparametrically efficient estimating equation in Eq. (3). Instantiating Eq. (7) for the simultaneous estimation of quantile and CVaR and rearranging, we obtain the following estimating equation:

where

We use ) to denote the conditional and unconditional cumulative distribution function of Y (t), respectively: for any The following proposition gives the asymptotic behavior of our proposed estimators for the quantile and CVaR of Y (1). This conclusion is proved by verifying all conditions in Theorem 3. Analogous conclusions also hold for Y (0) when all assumptions hold for t = 0 instead of t = 1.

Proposition 1 (LDML for Quantile and CVaR). Fix t = 1 and Let the estimator ˆbe given by applying Definition 2 to the estimating function in Eq. (18). Suppose Assumptions 5 and 6 hold and there exist positive constants , such that for any , the following conditions hold:

i. Conditions i. (with ), ii., v. (with ) of Assumption 2, condition iii. of Assumption 3, and condition vii. of Theorem 3 for the estimating function in Eq. (18) and the corresponding nuisance estimators.

ii. is twice differentiable with derivatives satisfying 0 . Moreover,

iii. At any is twice differentiable almost surely with first two order derivatives that satisfy almost surely.

Remark 4 (Estimating ) for Variance Estimation). If we want to conduct inference on the quantile or QTE using our method from Section 4, we need to estimate ). We only need to do this consistently, regardless of rate, in order to get the correct asymptotic coverage. One simple approach is to use cross-fitted IPW kernel density estimation at ˆ

where ) is a kernel function such as 0 is a bandwidth. Under Assumption 5, would be the optimal bandwidth. While this together with any consistent estimate ˆsuffices for asymptotic coverage, the estimate may be unstable. It is therefore recommended to use self-normalization by dividing the above by and to potentially clip propensities.

We first study the behavior of LDML in a simulation study. We then demonstrate its use in estimating the QTE of 401(k) eligibility on net financial assets. In Appendix A, we additionally consider estimating the LQTE of 401(k) participation using eligibility as IV.

6.1 Simulation Study

First, we consider a simulation study to compare the performance of LDML estimates to benchmarks. We consider estimating as the second tertile of Y (1) from incomplete data. The distri- bution P is as follows. First, we draw 20-dimensional covariates X from the uniform distribution on [0. Then, we draw T from Bernoulli(Φ(3(1 ))), where Φ is the standard normal cumulative distribution function, and we draw ). We only observe Y (1) when T = 1.

We consider estimating using four different methods. First, we consider LDML applied to the efficient estimating equation (Eq. (3)) with estimated using 2-fold cross- fitted IPW with random-forest-estimated propensities, and ˆ) similarly estimated by random forests. Second, we consider K = 5-fold cross-fitted IPW with random-forest-estimated propensities. Third, we consider DML with K = 5 and the estimand-dependent nuisance estimated using a discretization approach similar to the suggestion of Belloni et al. (2018): for 100 marginal quantile of ) using random forests; then apply DML with the restricted discretized estimand range We refer to this method as DML-D for discretized. Fourth, we consider DML with K = 5 and where the estimand-dependent nuisance is estimated using an approach similar to Meinshausen (2006), Bertsimas and Kallus (2014): namely, fit a random forest regression to the out-of-fold data regression trees : support(set ˆ. We refer to this method as DML-F for forest. For each method, we run it three times with new random fold splits (with the same data) and take the median of the three results to be the estimate.

For each of n = 100, 200, . . . , 25600, we consider 75 replications of drawing a dataset of size n and constructing each of the above four estimates. We plot the mean-squared error of each method and n over the 75 replications in Fig. 2(a). The shaded regions show plus/minus one standard error of this as the sample mean of 75 squared errors. We clearly see that LDML offers significant improvements over the other methods when we use flexible machine learning methods to tackle estimand-dependent nuisances.

In Fig. 2(b), we additionally report the coverage of the true parameter using the standard error estimate proposed in Remark 4 and Appendix E. Namely, for each of the three random runs of each method, we take the sample standard deviation of the estimated influence function evaluated at the final estimand and with the cross-fitted nuisances and divide it by times an estimate

Figure 2: Results for the simulation study. Shaded regions denote plus/minus one standard error for estimated mean-squared error or mean coverage computed over 75 replications.

of ) given by cross-fitted IPW kernel density estimation at the estimand. (We do the same for the IPW estimate for the sake of comparison but note IPW’s asymptotic variance may also depend on the propensity-estimation variance, unlike LDML and DML.) We take the median of these standard errors over the three runs and add to it the standard deviation of estimands over the three runs divided by3. Then we consider the 95% confidence interval given by the estimand plus/minus 1.96 of this estimated standard error. Figure 2(b) shows the sample mean coverage of over the 75 replications, and the shaded region shows plus/minus one standard error of this sample mean. LDML offers good, calibrated coverage (the 100% coverage for some n can be attributed to only observing 75 replications with 95% success probability each). The other methods have poor coverage, which may be attributed to significant bias so that confidence intervals based only on standard errors of the sample average would undercover and even get worse as samples grow and standard errors shrink relative to bias. In particular, the IPW estimate’s convergence directly depends on that of the random-forest-estimated propensities, so convergence may be slower than and/or the true standard error may be far larger than that of the cross-fitted sample average. Similarly, using DML to estimate the control-variate term using a discretization or a forest need not converge, so ultimately their convergence and standard errors may be similar to IPW’s, again leading to underconvering.

6.2 Effect of 401(k) Eligibility on Net Financial Assets

Next we consider an empirical case study to demonstrate the estimation of QTE using LDML in practice and with a variety of machine learning nuisance estimators. We use data from Cher- nozhukov and Hansen (2004) to estimate the QTEs of 401(k) retirement plan eligibility on net financial assets. Eligibility for 401(k) (here considered the treatment, T) is not randomly assigned, but is argued in Chernozhukov and Hansen (2004) to be ignorable conditioned on certain covariates: age, income, family size, years of education, marital status, two-earner household status, availability of defined benefit pension plan to household, IRA participation, and home ownership status. Net

Table 1: The QTE of 401(k) eligibility in thousand dollars (and standard error) estimated by LDML using different regression methods, and raw unadjusted differences of quantiles.

financial assets (the outcome, Y ) are defined as the sum of IRA and 401(k) balances, bank accounts, and other interest-earning accounts and assets minus non-mortgage debt. While Chernozhukov and Hansen (2004) considered controlling for these in a low-dimensional linear specification, it is not clear whether such is sufficient to account for all confounding. Consequently, Belloni et al. (2017) considered including higher-order terms and interactions, but needed to theoretically construct a continuum of LASSO estimates and may not be able to use generic black-box regression methods. Finally, Chernozhukov et al. (2018a) considered using generic machine learning methods, but only tackled ATE estimation.

In contrast, we will use LDML to estimate and conduct inference on the QTEs of 401(k) eligibility on net assets using a variety of flexible black-box regression methods. First, to understand the effect of different choices in the application of LDML to the problem, we consider estimating the 25%, 50%, and 75% QTE while varying K in {5, 15, 25} and varying the nuisance estimators. We consider estimating both propensity score and conditional cumulative distribution of: boosting (using R package gbm), LASSO (using R package hdm), and a one-hidden-layer neural network (using R package nnet). For LASSO, we use a 275-dimensional expansion of the covariates by considering higher-order terms and interactions. In each instantiation of LDML, we construct folds so to ensure a balanced distribution of treated and untreated units, we let we use the IPW initial estimator for ˆ, we normalize propensity weights to have mean 1 within each treatment group, we use estimates given by solving the grand-average estimating equation as in Definition 2, and for variance estimation we estimate using IPW kernel density estimation as in Remark 4. The solution to the LDML-estimated empirical estimating equation must occur at an observed outcome and that we can find the solution using binary search after sorting the data along outcomes. We re-randomize the fold construction and repeat each instantiation 100 times. We then remove the outlying 2.5% from each end and report ˆas in Appendix E. The resulting estimates and standard errors are shown in Table 1. The estimates appear overall roughly stable across methods and K.

Next, we consider estimating a range of QTEs. We focus on nuisance estimation using LASSO and fix K = 15. We then estimate the 10%, 11%, . . . , 89%, and 90% quantiles and QTEs. We plot the resulting LDML estimates with 90% confidence intervals in Fig. 3 and compare these to the raw unadjusted marginal quantiles within each treatment group.

Figure 3: LDML estimates of a range of quantiles and QTEs with confidence 90% intervals and comparison to raw unadjusted marginal quantiles by treatment group.

Semiparametric Estimation, Neyman orthogonality, and Debiased Machine Learning. Our work is closely related to the classical semiparametric estimation literature on constructing -consistent and asymptotically normal estimators for low dimensional target parameters in the presence of infinitely dimensional nuisances, typically estimated by conventional nonparametric estimators such as kernel or series estimators (e.g., Newey 1990, 1994, Newey et al. 1998, Ibragimov and Hasminskii 1981, Levit 1976, Bickel et al. 1998, Bickel 1982, Robinson 1988, var der Vaart 1991, Andrews 1994, Robins and Rotnitzky 1995, Linton 1996, Chen et al. 2003, Ai and Chen 2003, van der Laan and Rose 2011, Ai and Chen 2012). Our work builds on the Neyman orthogonality condition introduced by Neyman (1959)). This condition plays a critical role in many works that go beyond such nonparametric estimators, such as targeted learning (e.g., van der Laan and Rose 2011, Van der Laan and Rose 2018), inference for coefficients in high dimensional linear models (e.g., Belloni et al. 2016, 2014c, Zhang and Zhang 2014, Van de Geer et al. 2014, Javanmard and Montanari 2014, Chernozhukov et al. 2015, Ning et al. 2017), and semiparametric estimation with nuisances that involve high dimensional covariates (e.g., Belloni et al. 2017, Smucler et al. 2019, Chernozhukov et al. 2018b, Farrell 2015, Belloni et al. 2014a,b, Bradic et al. 2019, Bravo et al. 2020).

Chernozhukov et al. (2018a) further advocate the use of cross-fitting in addition to orthogonal estimating equations, so that the traditional Donsker assumption on nuisance estimators can be relaxed, and a broad array of black-box machine learning algorithms can be used instead. They refer to this generic approach as DML, which provides a principled framework to estimate low-dimensional target parameters with strong asymptotic guarantees when leveraging modern machine learning methods in nuisance estimation. Similar forms of sample splitting and cross-fitting have also appeared in Klaassen (1987), Zheng and van der Laan (2011), Fan et al. (2012), Bickel (1982), Robins et al. (2013), Schick (1986), Robins et al. (2008, 2017). Since the DML framework was introduced, numerous works have applied it in many different problems, such as heterogeneous treatment effect estimation (Kennedy 2020, Nie and Wager 2017, Curth et al. 2020, Semenova and Chernozhukov 2020, Oprescu et al. 2019, Fan et al. 2020), causal effects of continuous treatments (Colangelo and Lee 2020, Oprescu et al. 2019), instrumental variable estimation (Singh and Sun 2019, Syrgkanis et al. 2019), partial identification (Bonvini and Kennedy 2019, Kallus et al. 2019, Semenova 2017, Yadlowsky et al. 2018), difference-in-difference models (Lu et al. 2019, Chang 2020, Zimmert 2018), off-policy evaluation (Kallus and Uehara 2020, Demirer et al. 2019, Zhou et al. 2018, Athey and Wager 2017), generalized method of moments (Chernozhukov et al. 2016, Belloni et al. 2018), improved machine learning nuisance estimation (Farrell et al. 2018, Cui and Tchetgen 2019), statistical learning with nuisances (Foster and Syrgkanis 2019), causal inference with surrogate observations (Kallus and Mao 2020), linear functional estimation (Chernozhukov et al. 2018d,c, Bradic et al. 2019), etc. Our work complements this line of research by proposing a simple but effective way to handle estimand-dependent nuisances. This type of nuisances frequently appears in efficient estimation of complex causal effects such as QTEs, and applying DML directly would require estimating a continuum of nuisances, which is challenging in practice.

Efficient estimation of (L)QTE. Firpo (2007) first considered efficient estimation of QTE and proposed an IPW estimator based on propensity scores estimated by a logistic sieve estimator. Under strong smoothness conditions, this IPW estimator is-consistent and achieves the semi-parametric efficiency bound. Fr¨olich and Melly (2013) consider a weighted estimator for LQTE with weights estimated by local linear regressions using high-order kernels and show that their estimator is also semiparametrically efficient. Although these purely weighted methods bypass the estimation of nuisances that depend on the estimand, their favorable behavior is restricted to certain nonparametric weight estimators and strong smoothness requirements. D´ıaz (2017) proposed a Targeted Minimum Loss Estimator (TMLE) estimator for efficient QTE estimation. Built on the efficient influence function with nuisances that depends on the quantile itself, this estimator requires estimating a whole conditional cumulative distribution function, which as discussed may be very challenging in practice using flexible machine learning methods. Belloni et al. (2017) similarly consider efficient estimation of LQTE with high-dimensional covariates by using a Neyman-orthogonal estimating equation and discretizing a continuum of LASSO estimators for the estimand-dependent nuisance. In contrast, our proposed estimator can leverage a wide variety of flexible machine learning methods for the standard regression task to estimate nuisances, since we require estimating conditional cumulative distribution function only at a single point, which amounts to a binary regression problem.

Estimand-dependent nuisances. Besides (local) quantiles and CVaR, many efficient estimation problems involve nuisances that depends on the estimand (e.g., Tsiatis 2006, Chen et al. 2005). Previous approaches estimate the whole continuum of the estimand-dependent nuisances either by positing simple parametric model for conditional distributions (Tsiatis 2006, Chap 10), using sieve estimators (Chen et al. 2005), or discretizing a hypothetical continuum of regression estimators (Belloni et al. 2017). In contrast, our proposed method obviates the need to estimate infinitely many nuisances by fitting nuisances only at a preliminary estimate of the parameter of interest. This idea was briefly mentioned by Robins et al. (1994), focusing on parametric models for nuisance estimation. Our paper rigorously develops this approach and admits flexible machine learning methods for estimating nuisances that depend on the estimand.

In many causal inference and missing data settings, the efficient influence function involves nuisances that depend on the estimand of interest. A key example provided was that of QTE under ignorable treatment assignment and LQTE estimation using an IV, where in both cases the efficient influence function depends on the conditional cumulative distribution function evaluated at the quantile of interest. This structure, common to many other important problems, makes the application of existing debiased machine learning methods difficult in practice. In quantile estimation, it requires we learn the whole conditional cumulative distribution function. To avoid this difficulty, we proposed the LDML approach, which localized the nuisance estimation step to an initial rough guess of the estimand. This was motivated by the fact that in many applications, the oracle estimating equation is asymptotically equivalent to one where the nuisance is evaluated at the true parameter value, which our localization approach targets. Assuming only standard identification conditions, Neyman orthogonality, and lax rate conditions on our nuisance estimates, we proved the LDML enjoys the same favorable asymptotics as the oracle estimator that solves the estimating equation with the true nuisance functions. This newly enables the practical efficient estimation of important quantities such as QTEs using machine learning.

An interesting future direction is to consider a uniform estimation way over the quantile we consider the setting in which is fixed to emphasize our main point. This might be possible by combining recent technique in quantile regression (Ota et al. 2019, Bradic and Kolar 2017); however, the rigorous result is left as future research.

Chunrong Ai and Xiaohong Chen. Efficient estimation of models with conditional moment restrictions containing unknown functions. Econometrica, 71(6):1795–1843, 2003.

Chunrong Ai and Xiaohong Chen. Semiparametric efficiency bound for models of sequential moment restric- tions containing unknown functions. Journal of Econometrics, 170:442–457, 2012.

Donald W. K. Andrews. Asymptotics for semiparametric econometric models via stochastic equicontinuity. Econometrica, 62(1):43–72, 1994.

Susan Athey and Stefan Wager. Efficient policy learning. arXiv preprint arXiv:1702.02896, 2017.

Alexandre Belloni, Victor Chernozhukov, and Christian Hansen. High-dimensional methods and inference on structural and treatment effects. Journal of Economic Perspectives, 28(2):29–50, 2014a.

Alexandre Belloni, Victor Chernozhukov, and Christian Hansen. Inference on treatment effects after selection among high-dimensional controls. The Review of Economic Studies, 81(2):608–650, 2014b.

Alexandre Belloni, Victor Chernozhukov, and Lie Wang. Pivotal estimation via square-root lasso in non- parametric regression. The Annals of Statistics, 42(2):757–788, 2014c.

Alexandre Belloni, Victor Chernozhukov, and Ying Wei. Post-selection inference for generalized linear models with many controls. Journal of Business & Economic Statistics, 34(4):606–619, 2016.

Alexandre Belloni, Victor Chernozhukov, Ivan Fern´andez-Val, and Christian Hansen. Program evaluation and causal inference with high-dimensional data. Econometrica, 85(1):233–298, 2017.

Alexandre Belloni, Victor Chernozhukov, Denis Chetverikov, Christian Hansen, and Kengo Kato. High- dimensional econometrics and regularized gmm. arXiv preprint arXiv: 1806.01888, 2018.

Dimitris Bertsimas and Nathan Kallus. From predictive to prescriptive analytics. arXiv preprint arXiv:1402.5481, 2014.

P. J. Bickel. On adaptive estimation. Annals of Statistics, 10(3):647–671, 1982.

P. J. Bickel, C. A. J. Klaassen, Y. Ritov, and J. A. Wellner. Efficient and Adaptive Estimation for Semiparametric Models. Springer, 1998.

Matteo Bonvini and Edward H Kennedy. Sensitivity analysis via the proportion of unmeasured confounding. arXiv preprint arXiv:1912.02793, 2019.

Jelena Bradic and Mladen Kolar. Uniform inference for high-dimensional quantile regression: linear func- tionals and regression rank scores. arXiv preprint arXiv:1702.06209, 2017.

Jelena Bradic, Victor Chernozhukov, Whitney K. Newey, and Yinchu Zhu. Minimax semiparametric learning with approximate sparsity. arXiv preprint arXiv:1912.12213, 2019.

Jelena Bradic, Stefan Wager, and Yinchu Zhu. Sparsity double robust inference of average treatment effects. arXiv preprint arXiv:1905.00744, 2019.

Francesco Bravo, Juan Carlos Escanciano, and Ingrid Van Keilegom. Two-step semiparametric empirical likelihood inference. The Annals of Statistics, 48(1):1–26, 2020.

Neng-Chieh Chang. Double/debiased machine learning for difference-in-differences models. The Econometrics Journal, 23(2):177–191, 2020.

Xiaohong Chen, Oliver Linton, and Ingrid Van Keilegom. Estimation of semiparametric models when the criterion function is not smooth. LSE Research Online Documents on Economics, 2003.

Xiaohong Chen, Han Hong, and Elie Tamer. Measurement error models with auxiliary data. The Review of Economic Studies, 72(2):343–366, 2005.

Victor Chernozhukov and Christian Hansen. The effects of 401(k) participation on the wealth distribution: an instrumental quantile regression analysis. Review of Economics and statistics, 86(3):735–751, 2004.

Victor Chernozhukov, Christian Hansen, and Martin Spindler. Valid post-selection and post-regularization inference: An elementary, general approach. 2015.

Victor Chernozhukov, Juan Carlos Escanciano, Hidehiko Ichimura, Whitney K Newey, and James M Robins. Locally robust semiparametric estimation. arXiv preprint arXiv:1608.00033, 2016.

Victor Chernozhukov, Denis Chetverikov, Mert Demirer, Esther Duflo, Christian Hansen, Whitney Newey, and James Robins. Double/debiased machine learning for treatment and structural parameters. Econometrics Journal, 21:C1–C68, 2018a.

Victor Chernozhukov, Denis Nekipelov, Vira Semenova, and Vasilis Syrgkanis. Plug-in regularized estimation of high-dimensional parameters in nonlinear semiparametric models. arXiv preprint arXiv:1806.04823, 2018b.

Victor Chernozhukov, Whitney Newey, James Robins, and Rahul Singh. Double/de-biased machine learning of global and local parameters using regularized riesz representers. arXiv preprint arXiv:1802.08667, 2018c.

Victor Chernozhukov, Whitney K Newey, and Rahul Singh. Learning l2 continuous regression functionals via regularized riesz representers. arXiv preprint arXiv:1809.05224, 8, 2018d.

Kyle Colangelo and Ying-Ying Lee. Double debiased machine learning nonparametric inference with contin- uous treatments. arXiv preprint arXiv:2004.03036, 2020.

Yifan Cui and Eric Tchetgen Tchetgen. Bias-aware model selection for machine learning of doubly robust functionals. arXiv preprint arXiv:1911.02029, 2019.

Alicia Curth, Ahmed M Alaa, and Mihaela van der Schaar. Semiparametric estimation and inference on struc- tural target functions using machine learning and influence functions. arXiv preprint arXiv:2008.06461, 2020.

Mert Demirer, Vasilis Syrgkanis, Greg Lewis, and Victor Chernozhukov. Semi-parametric efficient policy learning with continuous actions. arXiv preprint arXiv:1905.10116, 2019.

Iv´an D´ıaz. Efficient estimation of quantiles in missing data models. Journal of Statistical Planning and Inference, 190:39–51, 2017.

Jianqing Fan, Shaojun Guo, and Ning Hao. Variance estimation using refitted cross-validation in ultrahigh dimensional regression. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 74 (1):37–65, 2012.

Qingliang Fan, Yu-Chin Hsu, Robert P Lieli, and Yichong Zhang. Estimation of conditional average treat- ment effects with high-dimensional data. Journal of Business & Economic Statistics, (just-accepted): 1–39, 2020.

Max H Farrell. Robust inference on average treatment effects with possibly more covariates than observations. Journal of Econometrics, 189(1):1–23, 2015.

Max H Farrell, Tengyuan Liang, and Sanjog Misra. Deep neural networks for estimation and inference. arXiv preprint arXiv:1809.09953, 2018.

Sergio Firpo. Efficient semiparametric estimation of quantile treatment effects. Econometrica, 75:259–276, 2007.

Dylan J Foster and Vasilis Syrgkanis. Orthogonal statistical learning. arXiv preprint arXiv:1901.09036, 2019.

Markus Fr¨olich and Blaise St´ephane Melly. Unconditional quantile treatment effects under endogeneity. Journal of Business & Economic Statistics, 31(3):346–357, 2013.

I. A. Ibragimov and R. Z. Hasminskii. Statistical estimation : asymptotic theory. 1981.

Guido W Imbens and Joshua D Angrist. Identification and estimation of local average treatment effects. Econometrica, pages 467–475, 1994.

Adel Javanmard and Andrea Montanari. Confidence intervals and hypothesis testing for high-dimensional regression. The Journal of Machine Learning Research, 15(1):2869–2909, 2014.

Nathan Kallus and Xiaojie Mao. On the role of surrogates in the efficient estimation of treatment effects with limited outcome data. arXiv preprint arXiv:2003.12408, 2020.

Nathan Kallus and Masatoshi Uehara. Double reinforcement learning for efficient off-policy evaluation in markov decision processes. Journal of Machine Learning Research, 21(167):1–63, 2020.

Nathan Kallus, Xiaojie Mao, and Angela Zhou. Assessing algorithmic fairness with unobserved protected class using data combination. arXiv preprint arXiv:1906.00285, 2019.

Maximilian Kasy. Uniformity and the delta method. Journal of Econometric Methods, 8(1):1–19, 2019.

Edward H Kennedy. Optimal doubly robust estimation of heterogeneous causal effects. arXiv preprint arXiv:2004.14497, 2020.

Chris AJ Klaassen. Consistent estimation of the influence function of locally asymptotically linear estimators. The Annals of Statistics, pages 1548–1562, 1987.

B. Ya. Levit. On the efficiency of a class of non-parametric estimates. Theory of Probability and Its Applications, 20(4):723–740, 1976.

Oliver B. Linton. Edgeworth approximation for minpin estimators in semiparametric regression models. Econometric Theory, 12(1):30–60, 1996.

Chen Lu, Xinkun Nie, and Stefan Wager. Robust nonparametric difference-in-differences estimation. arXiv preprint arXiv:1905.11622, 2019.

Nicolai Meinshausen. Quantile regression forests. Journal of Machine Learning Research, 7(Jun):983–999, 2006.

Whitney K. Newey. Semiparametric efficiency bounds. Journal of Applied Econometrics, 5(2):99–135, 1990.

Whitney K. Newey. The asymptotic variance of semiparametric estimators. Econometrica, 62(6):1349–1382, 1994.

Whitney K Newey and James L Powell. Asymmetric least squares estimation and testing. Econometrica: Journal of the Econometric Society, pages 819–847, 1987.

Whitney K. Newey, Fushing Hsieh, and James Robins. Undersmoothing and bias corrected functional estimation. 1998.

Jerzy Neyman. Optimal asymptotic tests of composite statistical hypotheses. Probability and Statistics, pages 416–44, 1959.

Xinkun Nie and Stefan Wager. Quasi-oracle estimation of heterogeneous treatment effects. arXiv preprint arXiv:1712.04912, 2017.

Yang Ning, Han Liu, et al. A general theory of hypothesis tests and confidence regions for sparse high dimensional models. The Annals of Statistics, 45(1):158–195, 2017.

Miruna Oprescu, Vasilis Syrgkanis, and Zhiwei Steven Wu. Orthogonal random forest for causal inference. In International Conference on Machine Learning, pages 4932–4941. PMLR, 2019.

Hirofumi Ota, Kengo Kato, and Satoshi Hara. Quantile regression approach to conditional mode estimation. Electronic Journal of Statistics, 13(2):3120–3160, 2019.

James Robins, Lingling Li, Eric Tchetgen, Aad van der Vaart, et al. Higher order influence functions and minimax estimation of nonlinear functionals. In Probability and statistics: essays in honor of David A. Freedman, pages 335–421. Institute of Mathematical Statistics, 2008.

James M. Robins and Andrea Rotnitzky. Semiparametric efficiency in multivariate regression models with missing data. Journal of the American Statistical Association, 90(429):122–129, 1995.

James M Robins, Andrea Rotnitzky, and Lue Ping Zhao. Estimation of regression coefficients when some regressors are not always observed. Journal of the American statistical Association, 89(427):846–866, 1994.

James M Robins, Peng Zhang, Rajeev Ayyagari, Roger Logan, Eric Tchetgen Tchetgen, Lingling Li, Thomas Lumley, Aad van der Vaart, HEI Health Review Committee, et al. New statistical approaches to semiparametric regression with application to air pollution research. Research report (Health Effects Institute), (175):3, 2013.

James M Robins, Lingling Li, Rajarshi Mukherjee, Eric Tchetgen Tchetgen, Aad van der Vaart, et al. Minimax estimation of a functional on a structured high-dimensional model. The Annals of Statistics, 45(5):1951–1987, 2017.

Peter M Robinson. Root-n-consistent semiparametric regression. Econometrica, 56(4):931–954, 1988.

R Tyrrell Rockafellar and Stanislav Uryasev. Conditional value-at-risk for general loss distributions. Journal of banking & finance, 26(7):1443–1471, 2002.

Anton Schick. On asymptotically efficient estimation in semiparametric models. The Annals of Statistics, pages 1139–1151, 1986.

Vira Semenova. Machine learning for set-identified linear models. arXiv preprint arXiv:1712.10024, 2017.

Vira Semenova and Victor Chernozhukov. Debiased machine learning of conditional average treatment effects and other causal functions. The Econometrics Journal, 2020.

Rahul Singh and Liyang Sun. De-biased machine learning for compliers. arXiv preprint arXiv:1909.05244, 2019.

Ezequiel Smucler, Andrea Rotnitzky, and James M Robins. A unifying approach for doubly-robust regularized estimation of causal contrasts. arXiv preprint arXiv:1904.03737, 2019.

Vasilis Syrgkanis, Victor Lei, Miruna Oprescu, Maggie Hei, Keith Battocchi, and Greg Lewis. Machine learn- ing estimation of heterogeneous treatment effects with instruments. In Advances in Neural Information Processing Systems, pages 15193–15202, 2019.

Anastasios. Tsiatis. Semiparametric Theory and Missing Data. Springer, New York, 2006.

Sara Van de Geer, Peter B¨uhlmann, Ya’acov Ritov, Ruben Dezeure, et al. On asymptotically optimal confidence regions and tests for high-dimensional models. The Annals of Statistics, 42(3):1166–1202, 2014.

Mark J. van der Laan and Sherri Rose. Targeted Learning: Causal Inference for Observational and Experimental Data. 2011.

Mark J Van der Laan and Sherri Rose. Targeted learning in data science. Springer, 2018.

A. W. van der Vaart. Asymptotic Statistics. Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, 1998. doi: 10.1017/CBO9780511802256.

Aad var der Vaart. On differentiable functionals. Annals of Statistics, 19(1):178–204, 1991.

Steve Yadlowsky, Hongseok Namkoong, Sanjay Basu, John Duchi, and Lu Tian. Bounds on the conditional and average treatment effect with unobserved confounding factors. arXiv preprint arXiv:1808.09521, 2018.

Cun-Hui Zhang and Stephanie S Zhang. Confidence intervals for low dimensional parameters in high dimen- sional linear models. Journal of the Royal Statistical Society: Series B: Statistical Methodology, pages 217–242, 2014.

Wenjing Zheng and Mark J van der Laan. Cross-validated targeted minimum-loss-based estimation. In Targeted Learning, pages 459–474. Springer, 2011.

Zhengyuan Zhou, Susan Athey, and Stefan Wager. Offline multi-action policy learning: Generalization and optimization. arXiv preprint arXiv:1810.04778, 2018.

Michael Zimmert. Efficient difference-in-differences estimation with high-dimensional common trend con- founding. arXiv preprint arXiv:1809.01643, 2018.

In Section 1.2, we mention that without the ignorability assumption, we can rely on an instrumental variable to identify local parameters, namely, solutions ) to the following local estimating equation:

We assume standard instrumental variable identification conditions: for potential treatments T(w) and potential outcomes Y (t, w), we have exclusion ), exogeneity (1), relevance P(T(1) = 1) > P(T(0) = 1), and monotonicity (0). Following Belloni et al. (2017), a Neyman orthogonal estimating equation for is given by

where

with nuisance functions

Here the second estimating equation denoted by the following auxiliary parameter

By redefining ˜, the estimating equation becomes

which apparently fits into our general framework in Eq. (1). Therefore, we can directly apply our LDML algorithm in Section 2.2 to estimate the local parameters ). We can also use the theory in Sections 3 and 4 to analyze the asymptotic distribution of the resulting estimators and estimate their asymptotic variances.

A.1 Estimating Local Quantiles

In particular, we take the local quantile estimation as an example, namely, the solution local estimating equation in Eq. (19) with

Its orthogonal estimating equation involves the following nuisance functions:

For better readability, we denote the event of being a complier, i.e., T(1) > T(0), as C, the nuisance functions as ˜. We fit estimators for the nuisance functions based on the sample-splitting scheme given in Definition 1, which we denote as ˆ˜ˆ˜)) respectively for k = 1, . . . , K. Finally, we obtain the estimator ˆempirical estimating equations in Definition 2 or Definition 5, specialized to Eqs. (20) and (23).

We next assume a strong form of the overlap and relevance assumptions and specify the convergence rates of the initial estimator and nuisance estimators. We again consider a generic treatment level in these two assumptions.

Assumption 7 (Strong Overlap and Relevance Assumptions). Assume that there exists a positive constant such that for any holds almost surely, and

Assumption 8 (Nuisance Estimation Rates)

In the following theorem, we derive the asymptotic distribution of the local quantile estimator, which is proved by verifying all assumptions in Theorem 1.

Proposition 2 (LDML for Local Quantile)be a compact set where given in Assumption 7. Let (ˆbe the LDML estimator given in either Definition 2 or Definition 5, specialized to Eqs. (20) and (23). Suppose that there exist constants such that the following conditions hold for any instance

i. Conditions i. (with ) and condition iii. of Assumption 3 for the estimating equation in Eqs. (20) and (23).

ii. For any , the distribution function of Y (t) for compliers, denoted as twice continuously differentible. Its first two order derivatives that

iv. For any , the conditional distribution of Y (t) given X, T (w) = 1, denoted as , is twice differentiable almost surely with

first two order derivatives that satisfy almost surely.

Then (ˆsatisfies the conclusion of Theorem 1 for given in Eq. (20) and

In particular, the local quantile estimator ˆis asymptotically linear with the following influence function:

where are given in Eq. (20). Analogous con- clusion for local quantiles of Y (0) holds when all assumptions above hold for t = 0.

A.2 Effect of 401(k) Participation on Net Financial Assets

Next, we estimate the effect of 401(k) participation on net assets. Participation in a 401(k) plan (here considered the treatment, T) is not randomly assigned: individuals with a preference for saving may save more in non-retirement accounts than others whether they were to participate in retirement savings or not. There may be many other confounding factors, such as the possibility of higher financial acumen of savers leading to higher net worth otherwise. It is unlikely that we can control for all these factors using observable covariates. Instead, we rely on instrumenting on eligibility since, as argued in Section 6.2, eligibility is ignorable given covariates. Additionally, one cannot participate if one is ineligible, ensuring monotonicity, and some eligible individuals do

Table 2: The LQTE of 401(k) participation in thousand dollars (and standard error) estimated by LDML using different regression methods, and raw unadjusted differences of marginal quantiles by eligibility.

Figure 4: LDML estimates of a range of local quantiles and LQTEs with confidence 90% intervals and comparison to raw unadjusted marginal quantiles by treatment group.

participate, ensuring relevance. Assuming that eligibility cannot affect net assets except through its effect on participation, we have that eligibility for a 401(k) (here considered as W) is valid IV. We can therefore use it to estimate local quantiles by and LQTEs of 401(k) Participation on the population of individuals that would participate if eligible.

We use LDML applied to the Neyman orthogonal estimating equation Eqs. (20) and (23). Again, we consider the impact of different choices in the application of LDML. We repeat the same specifi-cation as above, using each possible nuisance estimator to fit the conditional probabilities Eqs. (21) and (24). We display the results for the 25%, 50%, and 75% quantiles while varying K and the nuisance estimators in Table 2. The qualitative results regarding the stability of LDML across methods and K remain the same. Then, focusing as before on nuisance estimation using LASSO and on K = 15, we also estimate a range of local quantiles and QTEs, which we plot along with 90% confidence intervals in Fig. 3. Again, we compare to the raw unadjusted marginal quantiles within each treatment group.

We can also apply our method and analysis to estimating the (1), as defined in Eq. (6). Instantiating Eq. (7) for expectiles and rearranging, we get the following efficient estimating function from incomplete data:

The next result gives the asymptotic behavior of LDML applied to these equations.

Proposition 3. Fix t = 1 and let the estimator ˆbe given by applying either Definition 2 or Definition 5 to the estimating function in Eq. (25). Suppose Assumptions 5 and 6 hold and there exist positive constants , such that for any , the following conditions hold:

i. Conditions i. (with ), ii., v. (with ) of Assumption 2, condition iii. of Assumption 3, and condition vii. of Theorem 3 for the estimating function in Eq. (25) and the corresponding nuisance estimators.

ii. is continuous at . Moreover, for any such that given in Eq. (6).

iii. At any is almost surely differentiable with first-order derivative , and second-order derivative ˙that satisfies almost surely;

Then ˆsatisfies the conclusion of Theorem 1 for given in Eq. (25) and . Analogous conclusion for expectile of Y (0) holds when all assumptions above hold for t = 0.

When constructing confidence intervals, we only need to estimate ) to estimate be easily estimated by the inverse propensity reweighted estimator

Alternatively, it can be estimated by an imputation estimator based on ˆor a LDML estimator that uses both ˆ(see Remark 4).

In this part, we show that the IPW initial estimator given in Definition 4 can satisfy the conditions on ˆin Assumption 6.

Proposition 4 (IPW Initial Estimator Rate). Fix t = 1 and let the initial estimator ˆconstructed according to Definition 4 for k = 1, . . . , K. Assume the following (for t = 1):

i. For each satisfies the same conditions as for ˆAssumption 6.

iii. There exists a nuisance realization set Πthat contains the true propensity score the propensity score estimators ˆwith at least probability 1 . Moreover, any satisfies that

iv. For each , the function class 1is suitably measurable and its uniform covering entropy satisfies the following condition: for positive constants is a measurable envelope for There exists a positive constant such that for any

Then there exists a constant c that only depends on pre-specified constants in the conditions above such that with probability 1

In Remark 3, we discuss the corresponding rate conditions on other nuisance estimators when using the IPW initial estimator, based on the conclusion in Proposition 4.

In Definition 2, we construct an LDML estimator by first averaging estimates of the equation in Eq. (11) over all folds and then solving the grand-average equation approximately. Below we provide an alternative LDML estimator that first solves the estimate of Eq. (11) from each fold separately and then averages these solutions.

Definition 5 (LDML2). For k = 1, . . . , K, construct ˆby (approximately) solving

In fact, we allow for an approximate least-squares solution, which is useful if the empirical estimating equation has no exact solution. Namely, we let ˆbe any satisfying

Then, we let the final estimator be

We can easily follow previous proofs for the LDML estimator in Definition 2 to show that Theorems 1 to 3 and Proposition 1 also apply to ˆin Definition 5, provided that in Eq. (27) is , condition iii. in Assumption 3). For example, we demonstrate this at the end of the proof for Theorem 1. Thus the two LDML estimators in Definition 2 and Definition 5 are asymptotically equivalent.

The proposed LDML estimator ˆin Definition 2 or Definition 5 relies on nuisance estimates based on random sample splitting (Definition 1). Although the uncertainty due to sample splitting does not affect the asymptotic theory, it may influence the finite-sample performance of the LDML estimator.

To make the results more robust to sample splitting, we may consider aggregating the estimates over different random splitting realizations. In particular, it is possible to use many other different ways of splitting data. For example, in both Definitions 2 and 5 we may average more than just K solutions or equations. For each k, we can permute over alland folds used for fitting ˆ. Or, we could even permute over all ways to split {1, . . . , K}\{k} into two. Or, we can even repeat the initial random splitting into K folds many times over and average the resulting estimates from either Definition 5 or 2, or take their median to avoid outliers, or solve the grand-mean of estimating equations. All of these procedures can provide improved finite-sample performance in practice as they can only reduce variance without affecting bias, and we do recommend these, but they have no effect on the leading asymptotic behavior, which remains the same whether you use one or more splits of the data into folds and/or one or more splits of {1, . . . , K} \ {k} into two.

With estimates from multiple random splitting realizations, we may also improve variance estimation and to account for the variance due to random splitting. In particular, letting ˆbe the parameter and variance estimates for each run of LDML for s = 1, . . . , S, we can let ˆ) be the final parameter and variance estimates. Like ˆ, the first term in reduces the variance in the estimate itself. The second term in accounts for the variance of ˆdue to random splitting. Notice that the second term vanishes as ; indeed then ˆhas no variance due to random splitting as it is fully averaged over. Because ˆare each consistent, the second term also vanishes as . Removing the factor in the second term we can instead get an estimate of the variance of each single ˆ, rather than of ˆ, accounting for random splitting. This procedure extends a similar proposal by Chernozhukov et al. (2018a) for inference in linear estimating equations.

The key condition that motivates our LDML approach is the invariant Jacobian condition in Assumption 1. Below Assumption 1, we directly show that this condition is satisfied for efficient estimating equations in incomlete data settings. In the following proposation, we provide a more general sufficient condition for the invariant Jacobian condition.

Proposition 5 (Sufficient Conditions for Invariant Jacobian). Assume that the map (is Fr´echet differentiable at (. Namely, assume that there exists a bounded linear operator , such that for any (within a small open neighborhood

Assume further that there exists C > 0 such that for any (

Then Assumption 1 is satisfied.

Here the condition in Eq. (29) is an orthogonality condition using the Fr´echet derivative, which is stronger than the Gˆateaux differentiability required in Neyman orthogonality (see Assumption 2 condition vii.). In Eq. (12), we already show that ))] for the efficient estimating function in the incomplete data setting does not depend on at all. Thus, its Fr´echet derivative with respect to trivially exists and is always 0, and therefore our Assumption 1 will be satisfied per Proposition 5.

Our work is closely related to the classical semiparametric estimation literature on constructing -consistent and asymptotically normal estimators for low dimensional target parameters in the presence of infinitely dimensional nuisances, typically estimated by conventional nonparametric estimators such as kernel and series estimators (e.g., Newey 1990, 1994, Newey et al. 1998, Ibragimov and Hasminskii 1981, Levit 1976, Bickel et al. 1998, Bickel 1982, Robinson 1988, var der Vaart 1991, Andrews 1994, Robins and Rotnitzky 1995, Linton 1996, Chen et al. 2003, van der Laan and Rose 2011, Ai and Chen 2012). Our work builds on the Neyman orthogonality condition introduced by Neyman (1959) (see Eq. (29) in Proposition 5 below). This condition plays a critical role in many works that go beyond such nonparametric estimators, such as targeted learning (e.g., van der Laan and Rose 2011, Van der Laan and Rose 2018), inference for coefficients in high dimensional linear models (e.g., Belloni et al. 2016, 2014c, Zhang and Zhang 2014, Van de Geer et al. 2014, Javanmard and Montanari 2014, Chernozhukov et al. 2015, Ning et al. 2017), and semiparametric estimation with nuisances that involve high dimensional covariates (e.g., Belloni et al. 2017, Smucler et al. 2019, Chernozhukov et al. 2018b, Farrell 2015, Belloni et al. 2014a,b, Bradic et al. 2019).

Chernozhukov et al. (2018a) further advocate the use of cross-fitting in addition to orthogonal estimating equations, so that the traditional Donsker assumption on nuisance estimators can be relaxed, and a broad array of black-box machine learning algorithms can be used instead. They refer to this generic approach as DML, which provides a principled framework to estimate low-dimensional target parameters with strong asymptotic guarantees when leveraging modern machine learning methods in nuisance estimation. Other forms of sample splitting and cross-fitting have also appeared in Klaassen (1987), Zheng and van der Laan (2011), Fan et al. (2012), Bickel (1982), Robins et al. (2013), Schick (1986), van der Vaart (1998), Robins et al. (2008, 2017). Since the DML framework was introduced, numerous works have applied it in many different problems, such as heterogeneous treatment effect estimation (Kennedy 2020, Nie and Wager 2017, Curth et al. 2020, Semenova and Chernozhukov 2020, Oprescu et al. 2019, Fan et al. 2020), causal effects of continuous treatments (Colangelo and Lee 2020, Oprescu et al. 2019), instrumental variable estimation (Singh and Sun 2019, Syrgkanis et al. 2019), partial identification (Bonvini and Kennedy 2019, Kallus et al. 2019, Semenova 2017, Yadlowsky et al. 2018), difference-in-difference models (Lu et al. 2019, Chang 2020, Zimmert 2018), off-policy evaluation (Kallus and Uehara 2020, Demirer et al. 2019, Zhou et al. 2018, Athey and Wager 2017), generalized method of moments (Chernozhukov et al. 2016, Belloni et al. 2018), improved machine learning nuisance estimation (Farrell et al. 2018, Cui and Tchetgen 2019), statistical learning with nuisances (Foster and Syrgkanis 2019), causal inference with surrogate observations (Kallus and Mao 2020), linear functional estimation (Chernozhukov et al. 2018d,c, Bradic et al. 2019), etc. Our work complements this line of research by proposing a simple but effective way to handle estimand-dependent nuisances. This type of nuisances frequently appears in efficient estimation of complex causal effects such as QTEs, and applying DML directly would require estimating a continuum of nuisances, which is challenging in practice.

Our proof of Theorem 1 and the proof of Theorem 3.3 in Chernozhukov et al. (2018a) are overall similar, but critically differ in Step II. In Step II, both proofs are based on the following decompo-

sition:

where

and ) is proved analogously in both proofs.

However, our proof and the proof in Chernozhukov et al. (2018a) assume different rate on thus

Under our condition, , then jointly considering the left hand side and

right hand side in Eq. (30) gives ), which in turn implies that thus the asserted conclusion in Theorem 1. In contrast, the counterpart condition in Chernozhukov et al. (2018a) guarantees that ) directly without needing to consider both sides of Eq. (30) jointly.

Now we use the example of estimating equation for incomplete data to show that the condition Eq. (32) in Chernozhukov et al. (2018a) generally requires stronger conditions for the convergence rates of nuisance estimators than our condition Eq. (31).

According to Eq. (44), under suitable regularity conditions,

Since Step I in the proof of Theorem 1 already proves that to guarantee our condtion. Thus our condition in Eq. (31) only requires that the product error rates to vanish faster than ), which is common in debiased machine learning for linear estimating equation (Chernozhukov et al. 2018a).

In contrast, to guarantee the condition in Chernozhukov et al. (2018a) given in Eq. (32), we need to assume that , besides the conditions on product error rates. Therefore, following the proof in Chernozhukov et al. (2018a) directly will require the propensity score to converge faster than ), no matter how fast the initial estimator ˆand the regression estimator ˆ) converge.

I.1 Proofs for Section 2

) such that (, the asserted Fr´echet differentiability and orthogonality condition imply that

This means that

I.2 Proofs for Section 3

Proof for Theorem 1. Fix any sequence that generates the observed data () and satisfies that 1. Because this sequence is chosen arbitrarily, to prove that the asserted conclusion holds uniformly over , we only need to prove

For to represent the empirical average operator based on example, Analogously, is the empirical average operator for the whole dataset, i.e., is the empirical process operator Moreover, for a given population average operator P are all derived from the underlying true distribution supress such dependence for ease of notation. Throughout the proof, we condition on the event (ˆ, which happens with at least -probability 1according to Assumption 3 condition i.. All statements involving notations in this proof depend on only constants pre-specified in Assumptions 2 and 3, and do not depend on constants specific to the instance . This should be clear from the proof, and the fact that the maximal inequality in Lemma 6.2 of Chernozhukov et al. (2018a) only depend on pre-specified parameters. Here we prove the asymptotic distribution of ˆgiven in Definition 2 first.

-probability 1 (1). Here we prove this by showing that with -probability 1

so that Assumption 2 implies that

Since the singular values of are lower bounded by 0, we can conclude that with probability 1 exceeding an instance-independent threshold.

In order to prove Eq. (33), we use the following decomposition:

Denote

Then obviously,

To bound (c), note that Eq. (14) implies

Thus

(

Therefore,

Note that Assumption 3 condition ii. implies that and the Assumption 3 condition iii. implies that

for a positive constant that only depends on specified in Assumption 2.

Then conditionally on ˆ), we can use Lemma 6.2 eq. (A.1) in Cher- nozhukov et al. (2018a) to prove that with -probability 1

which also holds unconditionally according to Lemma 6.1 of in Chernozhukov et al. (2018a). This further implies that Thus

In Step I, we proved that -probability 1 (1). Conditioned on this event, we will show that

where

Here Assumption 3 condition ii. guarantees that and the assumption that guarantees that . In step III and IV, we will further bound respectively.

Consequently, with -probability 1

and

By Assumption 2 condition v. and Markov inequality, ). Thus, with -probability 1

Plugging this back into Eq. (36) gives

Thus,

because

Now we prove the decomposition Eq. (35). Note that for any

If we apply Eq. (37) with th fold, and apply Eq. (14), then

Here

and the second order tayler expansion at

where the third equality uses the Neyman orthogonality in Assumption 2 condition vii..

Combining Eq. (38), Eq. (39) and Eq. (40) gives decomposition Eq. (35).

, we still condition on ˆ, and then apply Lemma 6.2 in Chernozhukov et al. (2018a) with function class

, . . . , d, θ

We can verify that satisfies similar entropy condition with envelope Moreover, Assumption 3 implies that

Thus conditionally on ˆ), we can use Lemma 6.2 eq. (A.1) in Cher- nozhukov et al. (2018a) to show that with -probability 1

which also holds unconditionally according to Lemma 6.1 in Chernozhukov et al. (2018a) .

Since is nonsingular with singular values bounded away from 0 by According to Assumption 2 condition i., probability 1 (1). Therefore,

Then apply the linearization Eq. (37) and taylor expansion similar to Eq. (40) with

where the last equality here holds because consequence of the special construction of

By applying step I to IV to sample estimating equation Eq. (27), we can get that for k = 1, . . . , K,

Since K is a fixed integer that does not grow with N, the equation above implies that the asserted conclusion in Theorem 1 also holds for ˆ

I.3 Proofs for Section 4

Proof of Theorem 2. We still consider data generating processes defined in the proof for Theorem 1, and define . Now we prove that

for any Then, the statement in Theorem 2 is immediately concluded. For all ), Eq. (41) follows once we have

Here, to simplify the notation, we use ˆ). Obviously we have

and we show that each term here is

We first bound . This is upper bounded as

Here, we use the fourth moment assumption in Assumption 4. From conditional Markov inequality, we have

Next, we bound . Following the proof of Theorem 3.2 (Chernozhukov et al. 2018a), we have

In addition, from the fourth moment assumption in Assumption 4

It follows from Markov inequality that

It remains to bound

Then, the first term of Eq. (42) is upper bounded with -probability 1

In the last inequality, we use Lemma 6.2 (Chernozhukov et al. 2018a). Here, the envelops exists since

for some constant C depending on . According to condition vi. in Assumption 2, it satisfies the moment condition . In addition, the metricy entropy assumption is satisfied since

Similarly, with -probability probability 1 (1), the second term of Eq. (42) can be upper bounded as follows:

In the last inequality, we use Lemma 6.2 (Chernozhukov et al. 2018a) and Assumption 3. In the end, we have

This concludes the proof.

I.4 Proofs for Section 5

Proof for Theorem 3. In this part, we prove the asymptotic distribution of our estimators corresponding to the general estimating equation Eq. (7). We prove this by verifying all conditions in the assumptions for Theorem 1.

Verifying Assumption 1.

where the second and third equality follow because

Verifying Assumption 2. We first verify conditions iii. and iv. in Assumption 2. We denote that component of condition ii., ) is Lipschitz continuous at with Lipschitz constant . So for any belongs to the open ball

By first-order Taylor expansion, for any ), there exists ) such that

where the second last inequality holds if we choose is the smallest singular value of

Moreover, the singular values are bounded between according to condition iii..

We then verify condition vii. in Assumption 2: for any (

to be the set that contains all ()) that satisfies the following conditions:

Then Assumption 6 and condition vii. in Theorem 3 guarantee that the nuisance estimates (ˆwith probability, namely, condition i. in Assumption 3 is satisfied.

Before verifying other conditions, first note that the condition vi. states that

Now we verify the condition on : for any (

Thus, the condition on is satisfied with

Next, we verify the condition on , and any ((

Finally, to verify the condition on , we note that for any any (

Thus the first-order derivative is

The second order derivative is

+P

+P

Above, we use condition iv. in Theorem 3 to ensure exchange of integration and differentiation so we can get terms

=

=

and

Thus for any

Given . Moreover, when Consequently,

Proof for Proposition 1. First note that

When ) is differentiable, )] is also differentiable by Leibnitz integral rule, with derivative

Thus

Now we prove Proposition 1 by verifying the assumptions in Theorem 3.

Verifying condition i in Theorem 3. We only need to verify that condition vi. of Assumption 2 hold. Since Θ is compact, obviously Donsker classes, so condition vi. of Assumption 2 is satisfied.

Verifying conditions ii. and iii. in Theorem 3. It is straightforward to show that invertible with the following matrix as its inverse:

Note that

Thus condition iii. in Theorem 3 is satisfied with When we estimate quantile only, then only ) matters. Then condition iii. in Theorem 3 is satisfied with

Since , it follows that each element in ) is Lipschtiz continuous at with Lipschitz constant . Moreover, for Θ such that have 2. This means that condition ii. in Theorem 3 is satisfied with . When we estimate quantile only, we only require for . Then condition ii. in Theorem 3 is satisfied with

Verifying condition iv. in Theorem 3. This condition can be verified by the following facts: for any

and

Verifying conditions v. and vi. in Theorem 3. For any (

By first-order Taylor expansion, for any , there exists ˜

Moreover, for any

Proof for Proposition 4. We follow the proof of Theorem 1 to consider any sequence of data generating process but we suppress it for ease of notation. We prove the conclusion for a generic , we denote as the empirical average operator and empirical process operator for data in the . Throughout the proof, we condition on the event that the convergence rate of propensity score estimator ˆin mean squared error is is lower bounded by , which holds with at least probability 1 according to Assumption 6. In this proof, all notations only involve pre-specified constants and not any instance-dependent constants.

We use the following decomposition analogous to that in Step I of proof for Theorem 1.

where

Given that condition vi. in Assumption 2 is satisfied for the estimating equation follow the end of step I in the proof for Theorem 1 to prove that with probability 1

for a constant c that depends on only constants in the assumptions,

Therefore, with -probability 1

In the proof of Theorem 3, we have showed that conditions ii. and iii. in Theorem 3 imply that

Therefore, with probability 1

I.5 Proofs for Appendix

Proof of Proposition 2. In this part, we prove the asymptotic distribution of our estimator ˆcorresponding to Eqs. (20) and (23). We denote prove this by verifying all conditions in the assumptions in Theorem 1.

Verifying Assumption 1. Similar to the proof of Theorem 3, we can easily show that

Verifying Assumption 2. We first verify conditions iii. and iv. in Assumption 2. We can easily derive that

and its Jacobian matrix is given by

This means that

Therefore,

and

The latter implies that

Therefore, condition iv. in Assumption 2 is satisfied with

Moreover, for any (, so we have

that entries in ) are all Lipschtiz with

a valid Lipschitz constant. Moreover, we have 2Θ such that . By following the proof of Theorem 3, we can easily

verify condition iii. in Assumption 2.

Next, we verify condition vi. in Assumption 2. For any fixed , the class only through and . Since the latter two classes are Donsker class, vi. in Assumption 2 for the function class has to be satisfied as well.

to be the set that contains all (())) that satisfies the following conditions: for w = 0, 1,

given in Eq. (45).

Then Assumption 8 and Proposition 2 condition v. ensure that the nuisance estimates (ˆwith probability, namely, condition i. in Assumption 3 is satisfied.

Before verifying other conditions, we first note that

It follows from Item iv. that for any

Next, we verify Assumption 3 condition ii..

We first verify the condition on . By following the proof of Theorem 3, we can show that for any (

The last inequality holds because

and so is ˜). This means that the condition on is satisfied with

Next, we verify the condition on . Again, by following the proof of Theorem 3, we have that for

and

where

and

Also define

Morevoer, when ˜

Moreover, we can similarly show that

provided that ˜

which verifies Assumption 3 condition ii..

Therefore, we have

This means that

Proof for Proposition 3. We only need to verify the conditions in Theorem 3.

Verifying condition i. in Theorem 3. We only need to verify that condition vi. of Assumption 2 hold. Since Θ is compact, obviously Donsker classes, condition vi. of Assumption 2 is satisfied.

Verifying condition ii. and iii. in Theorem 3. According to Eq. (6), the estimating function for complete data is given by

It follows that

Here the differentiability of )] is guaranteed by Leibniz integral rule, the con- tinuity of its derivative at is guaranteed by the continuity of ), whose singular value bounded between implies that )] is Lipschtiz continuous with Lipschitz constant maxTherefore, the constants in condition ii. and constant in iii. of Theorem 3 can be set as ensures the condition ii. of Theorem 3 with constant

Verifying condition iv. in Theorem 3. Note that for any

Thus

which trivially imply condition iv. in Theorem 3.

Verifying condition iv in Theorem 3. Again

The the asserted assumpton iv means thatany is upper bounded by

Plus, for any

and there exists ˜