Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
158,644 characters · 18 sections · 50 citation commands
Econometric Inference with Machine-Learned Proxies: Partial Identification via Data Combination
\address{Johns Hopkins University} \email{[email removed]}
\dedicatory{Johns Hopkins University}
There is an increasing prevalence in economics and the social sciences of leveraging complex, unstructured data---such as text and imagery---to construct variables central to empirical analysis. In these contexts, machine learning (ML) algorithms are frequently employed to map these raw inputs onto proxies for latent constructs that would otherwise be inherently unobservable or prohibitively costly to quantify at scale. Such approaches have been instrumental in deriving measures of media slant groseclose_measure_2005,gentzkow_what_2010, firm-level political risk hassan_firm-level_2019, local economic activity hu_illuminating_2022, remote-work status in job postings hansen_remote_2023, and air pollution exposure sager_clean_2025. For broader overviews, see mullainathan_machine_2017 and gentzkow_text_2019, along with the references therein.
This procedure is typically implemented in two steps. First, in an upstream step, a prediction rule $g$ is trained on a labeled dataset to map unstructured or high-dimensional data $X$ to a target variable $Z$. Second, in a downstream step, a researcher studies an economic model with observed covariates $W$, latent vector $Z$, and parameter $\theta$. Although $Z$ is not observed in the downstream step, the researcher does observe the unstructured inputs $X$, which allows the researcher to apply the prediction rule from the upstream step and construct the proxy $\hat{Z}=g(X)$. While we distinguish between these two steps for clarity, the distinction is conceptual, since the upstream ML training and the downstream empirical estimation may be conducted by the same researcher.
While ML-generated proxies have substantially expanded the scope of empirical research, incorporating them into downstream economic models poses a fundamental challenge. A naive plug-in approach that treats the ML-generated proxy $\hat{Z}$ as if it were the true latent variable $Z$ ignores both measurement error and the generated-regressor problem, and can therefore lead to biased estimation and invalid inference. This issue is compounded by two additional features of modern ML applications. First, the prediction rules used in practice are often highly complex, making it analytically difficult to characterize the statistical properties of $\hat{Z}$. In particular, convergence rates for $\hat{Z}$ are often unavailable, and in some settings even consistency of $\hat{Z}$ for $Z$ is unclear. Second, because the unstructured input $X$ may contain rich information not only about $Z$ but also about the observed covariates $W$, the resulting measurement error $Z-\hat{Z}$ is generally nonclassical: it may depend on $Z$, remain correlated with $W$ even conditional on $Z$, and even be endogenous to the downstream economic model.
In this paper, we propose a new framework to address these challenges by leveraging an auxiliary validation sample. We consider a setting in which the researcher has access to two distinct datasets. The first is a downstream sample containing observed covariates $W$, unstructured inputs $X$, and the proxy $\hat{Z}$ constructed from $X$ using the upstream prediction rule $g$. The second is an auxiliary validation sample containing observations on $Z$ and $X$, and therefore also on the proxy $\hat{Z}=g(X)$. Such validation data often arise naturally in practice, for example as a held-out sample from the upstream training data used to assess the predictive performance of the rule $g$. Importantly, we do not require the validation sample to contain $W$, nor do we assume any direct correspondence or individual-level matching between the observations in the validation and downstream samples.
The key idea of this paper is to reframe the role of $\hat{Z}$. Rather than treating $\hat{Z}$ as a literal but noisy substitute for $Z$, we view it as a low-dimensional summary of the raw unstructured input $X$ that links the validation and downstream samples. Because the validation sample contains joint observations on $(Z,\hat{Z})$, it allows us to learn the conditional distribution of $Z$ given $\hat{Z}$. Because $\hat{Z}$ is also observed in the downstream sample, this information can then be carried over to the downstream setting: for each observation $(W_i,\hat{Z}_i)$, although the realization of $Z_i$ remains unobserved, the validation sample informs the conditional distribution of $Z_i$ given $\hat{Z}_i$. In this sense, $\hat{Z}$ serves not merely as a noisy proxy for $Z$, but as the key linking variable that allows information about $Z$ from the validation sample to be incorporated into the downstream economic analysis.
We then develop an identification strategy based on this idea. A natural way to implement the linking role of $\hat{Z}$ is through an optimal transport (OT) characterization conditional on $\hat{Z}$, as in fan_partial_2025. However, such a conditional OT strategy can be difficult to implement in practice, especially when $\hat{Z}$ is continuous or high-dimensional, because it requires estimating the conditional distribution of $Z$ given $\hat{Z}$ and solving a separate OT problem for each value in the support of $\hat{Z}$. To overcome this difficulty, we instead develop an unconditional OT characterization based on the decoupling idea of li_finite_2025. Our approach replaces the continuum of conditional transport problems with a single OT problem formulated in terms of the unconditional distribution of $(W,\hat{Z})$ and the unconditional distribution of $(Z,\hat{Z})$, thereby rendering the method tractable in practice. Importantly, this unconditional characterization remains sharp for the identified set of $\theta$: given the available data and maintained assumptions, the resulting bounds cannot be further tightened.
Finally, we develop an inference procedure for testing whether a candidate value of $\theta$ belongs to the identified set. Inference is challenging in this setting because OT problems typically exhibit nonstandard asymptotic behavior. A common workaround is to replace the OT problem with an entropy-regularized counterpart, but that approach is not suitable here: a fixed entropy penalty introduces regularization bias that invalidates the underlying OT characterization, whereas accommodating that bias through worst-case bounds leads to excessively conservative inference. We instead base our procedure on the Kantorovich dual formulation together with a sieve approximation to the relevant dual function space. The number of sieve terms is allowed to increase with the sample size, so that the approximation error, and hence its effect on the test statistic, vanishes asymptotically. To maintain tractability, the resulting procedure avoids resampling methods such as the bootstrap or subsampling. Instead, by applying sample splitting and cross-fitting, we obtain critical values directly from standard normal quantiles. This yields a tractable inference procedure in practice, as illustrated in our simulation exercises.
Our approach differs from the existing literature. In addressing the measurement error embedded in $\hat{Z}$, prior work has largely proceeded along two lines:
One line imposes structural restrictions on the measurement error or asymptotic requirements on the proxy itself. A standard example is a conditional independence assumption, under which the proxy $\hat{Z}$ is independent of auxiliary measurements or covariates conditional on the true latent variable $Z$ hu_instrumental_2008, hu_illuminating_2022. Related identifying assumptions also appear in recent work by rambachan_program_2024 when a joint sample of $(W,Z,\hat{Z})$ is unavailable. A different version of this approach assumes that the measurement error in $\hat{Z}$ vanishes at a rate comparable to sampling error battaglia_inference_2025. In the present setting, however, these assumptions can be restrictive. Conditional independence is difficult to justify when the unstructured input $X$ contains information that may remain correlated with other covariates in the economic model even after conditioning on $Z$, while rate conditions are hard to verify when $\hat{Z}$ is generated by complex ML procedures whose asymptotic properties are analytically intractable.
A second line assumes access to a complete validation sample containing joint observations on $(W,Z,\hat{Z})$ that is sufficient for point identification of $\theta_0$ on its own chen_unified_2000, chen_measurement_2005, allon_machine_2023, angelopoulos_prediction-powered_2023, angelopoulos_ppi_2024, zhang_debiasing_2023, zrnic_cross-prediction-powered_2023, miao_task-agnostic_2024, rambachan_program_2024, kluger_prediction-powered_2025, sanford_adversarial_2025. Within this framework, the downstream sample of $(W,\hat{Z})$ serves mainly to improve estimation efficiency. Although this approach avoids both structural assumptions on the measurement error and formal asymptotic requirements on the ML procedure, it imposes substantially stronger data requirements. In particular, this approach requires either that the downstream model not depend on the covariates $W$, or that the validation sample contain precisely the covariates $W$ needed for the downstream analysis. In practice, such complete validation samples are not always available: upstream ML researchers typically train prediction rules and release held-out validation data on $(Z,\hat{Z})$ without knowing which covariates $W$ will matter for a particular downstream application. Conversely, downstream researchers often face prohibitive costs or practical barriers to collecting the true values of $Z$ needed to obtain a joint sample of $(W,Z,\hat{Z})$.
The first contribution of this paper is to provide an alternative framework that avoids the limitations of existing approaches. Unlike the first line of the literature, our framework does not impose structural assumptions on the measurement error $Z-\hat{Z}$, nor does it require formal statistical guarantees, such as convergence rates or consistency, for the proxy $\hat{Z}$. Instead, we rely on validation data to learn the joint distribution of the true latent variable $Z$ and the proxy $\hat{Z}$. In this respect, our framework is closer in spirit to the second line of the literature. The key difference is that we do not require a complete validation sample that, by itself, is sufficient for consistent estimation of the parameter of interest. The trade-off is that our framework generally delivers partial rather than point identification for the model parameter. The informativeness of our partial-identification framework depends on the quality of the proxy. When the ML-generated proxy is highly accurate, the identified set for the parameter would be tight. In the extreme case in which $\hat{Z}=Z$, the parameter is point identified within our framework. At the same time, the validity of our approach does not depend on the accuracy of the proxy. When $\hat{Z}$ is only a crude measure of $Z$, the resulting identified set remains valid, although it may be wide and therefore only weakly informative about the parameter. This feature faithfully reflects the informational content of a poor proxy: when the proxy contains little reliable information about the latent variable, one should not expect the data to deliver sharp conclusions about the parameter of interest.
Beyond bypassing restrictive structural or asymptotic assumptions and strong data requirements, our framework offers several practical advantages. Because we treat $\hat{Z}$ not as a plug-in substitute for $Z$, but rather as a linking variable between the validation and downstream samples, the framework naturally accommodates settings in which $Z$ and $\hat{Z}$ take values in different spaces. For example, if $Z$ is a discrete group identity, $\hat{Z}$ need not itself be a categorical label. It may instead take the form of a vector of predicted probabilities, a record of the model's top-two guesses, a likelihood ranking, or some other statistical summary. When available, such richer representations of $\hat{Z}$ can carry more information about $Z$ and therefore lead to a more informative identified set. This flexibility also provides a simple way to combine proxies produced by multiple competing ML procedures: one can define a multidimensional $\hat{Z}$ whose components collect the outputs of the different procedures. More broadly, our perspective suggests that, in empirical economics, the quality of an ML procedure should be judged not solely by its out-of-sample predictive accuracy, but by the extent to which its output preserves the economically relevant information contained in the unstructured input $X$ for downstream moment conditions.
As a second contribution, the identification strategy developed in this paper applies more broadly to general data combination problems cross_regressions_2002,ridder_chapter_2007 and is therefore of independent interest. Because our framework recasts the proxy as a variable that links the downstream and validation samples, our identification strategy carries over directly to classical data combination environments. Such problems are pervasive in empirical work and remain of substantial econometric interest. Recent applications include assessing algorithmic fairness kallus_assessing_2022, measuring intergenerational mobility santavirta_name-based_2024, estimating long-term treatment effects athey_surrogate_2025, obradovic_identification_2026, and combining stated and revealed preferences meango_combining_2025.
Recent advances in this literature have used OT tools to characterize data combination problems hwang_bounding_2025, dhaultfoeuille_partially_2025, fan_partial_2025. The paper most closely related to ours is fan_partial_2025, who develop a conditional OT characterization for the same class of general moment models studied here. However, implementing their approach requires solving a separate transport problem for each realization of the overlap variable shared across the two datasets---in our setting, $\hat{Z}$. When this linking variable has large support or includes continuous components, the conditional approach can become computationally prohibitive. The unconditional OT characterization developed in this paper complements their result by providing an alternative formulation. Rather than solving a continuum of conditional OT problems and estimating the corresponding conditional distributions, we solve a single OT problem defined over unconditional distributions. The trade-off is that the resulting unconditional OT problem is higher-dimensional than each conditional problem in fan_partial_2025. We view this as a favorable trade-off when the support of $\hat{Z}$ is large or when $\hat{Z}$ contains continuous components.
The remainder of the paper is organized as follows. Section (ref) introduces the analytical framework and discusses several empirical applications that fit into this setup. Section (ref) establishes the identification results when an auxiliary validation sample containing both the proxy and the true latent target variable is available, while Section (ref) develops the corresponding statistical inference procedures. Section (ref) demonstrates the finite-sample performance of the proposed methods through Monte Carlo simulations. Section (ref) explores an extension of the benchmark model to settings in which the validation sample is only available to identify the marginal distributions of subvectors of $(Z, \hat{Z})$. Finally, Section (ref) concludes.
Throughout the paper, capital letters (e.g., $Z$, $\hat{Z}$, $W$, and $S$) denote random variables or random vectors, and lower-case letters (e.g., $z$, $\hat{z}$, $w$, and $s$) denote their corresponding realizations. All random variables are defined on a common complete probability space. Vectors are written as column vectors, and the superscript $\top$ denotes transpose. For vectors in finite-dimensional Euclidean spaces, $\left\Vert \cdot\right\Vert $ denotes the Euclidean norm. We write $\mathbb{E}_F[\cdot]$ for expectation under the distribution $F$; when no ambiguity arises, we simply write $\mathbb{E}[\cdot]$. For a distribution $F$, we write $\mathcal{L}^1(F)$ for the space of $F$-integrable functions. We use $\coloneqq$ to denote definitions. Unless stated otherwise, all equalities and inequalities involving random variables are understood to hold almost surely.
We study a two-stage empirical workflow consisting of an upstream machine-learning (ML) prediction step and a downstream econometric inference step. In the downstream step, the applied researcher seeks to learn about a finite-dimensional parameter $\theta_0\in\Theta$, characterized by the population moment conditions
where $W$ is a vector of observed covariates, $Z$ is a vector of latent target variables that are not observed in the downstream dataset, and $q(\cdot)$ is a known moment function. We assume that $W$ and $Z$ take values in finite-dimensional Euclidean spaces.
Because $Z$ is unobserved, the downstream researcher relies on a proxy generated by an upstream prediction rule,
where $X$ denotes observed predictors. The prediction rule $g(\cdot)$ is produced upstream, using separate training data on $(X,Z)$. In this paper, $X$ may consist of high-dimensional structured covariates or unstructured inputs such as text or images.
We impose minimal structure on the upstream learning step. The training procedure is left completely unrestricted, and the prediction rule $g(\cdot)$ could be misspecified for the conditional distribution of $Z$ given $X$. In particular, our analysis does not rely on any consistency, rate-of-convergence properties, or asymptotic distribution results of the estimated predictor. Dispensing with such asymptotic requirements is particularly appealing for empirical work, as formal convergence analysis are often intractable for the complex, off-the-shelf machine learning algorithms deployed in practice. For expositional simplicity in this section, we assume that the proxy $\hat{Z}$ and the latent target $Z$ share the same dimensionality. Nevertheless, our methodology readily accommodates settings in which $\hat{Z}$ and $Z$ take values in spaces of differing dimensions, as we discuss later in Section (ref).
Our framework distinguishes the informativeness of the proxy from the validity of the econometric procedure. The predictive accuracy of $g(\cdot)$ governs the informativeness of the resulting bounds. Intuitively, stronger predictive performance yields tighter restrictions linking $Z$ and $\hat Z$ and therefore a smaller identified set for $\theta_0$ in the downstream problem. By contrast, the validity of our approach does not hinge on $g(\cdot)$ being highly accurate or statistically consistent. Instead, validity is ensured whenever the downstream researcher has access to an auxiliary validation sample containing joint observations of $(Z,\hat Z)$. The joint distribution of $(Z,\hat Z)$ quantifies both the predictive performance of $g(\cdot)$ and the nature of the discrepancy between the proxy $\hat Z$ and the latent target variable $Z$. This data requirement is often modest. For example, the validation sample may be a subset of labeled observations held out from the data used to train $g(\cdot)$, or it may be obtained through a separate validation exercise after $g(\cdot)$ is trained. In the next section, we formalize the assumptions governing the validation sample and show how to combine this information with the downstream observations and moment conditions to construct sharp bounds on $\theta_0$ that account for the proxy-induced measurement error.
To summarize, in this two-stage empirical setting, the downstream researcher observes a sample drawn from the joint distribution of $(W,X)$ but lacks direct observations of the target variable $Z$. Consequently, the sample analogs of the moment conditions in (ref) cannot be directly evaluated. This motivate us to develop an approach that leverages the observed proxy $\hat{Z}$, alongside the joint distribution of $(Z, \hat{Z})$ recovered from the validation sample, to bound the moments and partially identify $\theta_0$.
We conclude this section with several examples illustrating how empirical applications map into the framework of downstream econometric model (ref) and the upstream prediction (ref).
In this section, we study partial identification of the structural parameter $\theta_0$ when an auxiliary validation sample is available. We assume that this sample contains joint observations of the latent target variable $Z$, the ML-generated proxy $\hat{Z}$, and, potentially, a low-dimensional vector of observable characteristics $h(X)$ extracted from the raw predictors $X$. This data environment arises naturally when an upstream data provider releases a held-out sample of $(Z,X)$ together with the trained prediction rule, allowing the downstream researcher to compute both $\hat{Z} = g(X)$ and $h(X)$ and hence to obtain validation data on $(Z, \hat{Z},h(X))$.
Let $S \coloneqq h(X)$ denote a low-dimensional attribute of $X$. Its role is to capture observable strata along which the relationship between $Z$ and $\hat{Z}$ may vary. For example, in Example (ref), $S$ might indicate the geographic region or time period associated with the remote-sensing data; in Example (ref), it might indicate the language or industry of the job posting. At one extreme, $S$ may simply be a constant and therefore carry no additional information. More informative choices of $S$, however, could tighten the identified set by exploiting variation in the conditional joint distribution of $(Z,\hat{Z})$ across subpopulations indexed by $S$.
Throughout the identification analysis, we treat the prediction rule $g(\cdot)$ and the feature map $h(\cdot)$ as fixed and deterministic functions. This treatment is justified under the standard assumption that the upstream training sample used to train $g(\cdot)$ is independent of both the validation sample and the downstream empirical sample. Under this assumption, the population analysis can be interpreted conditional on the realized $(g,h)$, so that the only randomness in $\hat{Z}=g(X)$ and $S=h(X)$ arises from the randomness in $X$.
A useful feature of our setup is that the validation sample need not contain the downstream covariates $W$. This is empirically important. Upstream data collection is typically organized around obtaining labeled pairs $(Z,X)$ for prediction, whereas the covariates $W$ are application-specific and often unavailable to the upstream provider. Moreover, from the downstream researcher's perspective, obtaining ground-truth measurements of $Z$ in the downstream sample is often costly or infeasible. We therefore assume the validation data to contain $(Z,\hat{Z},S)$ but not $W$.
The downstream and validation samples are informative about the same latent target only if they can be viewed as marginals of a common population distribution. We formalize this requirement as a compatibility condition. Let $F$ denote the population distribution of the downstream observables $(W,\hat Z,S)$, where $\hat Z=g(X)$ and $S=h(X)$ are computed using the fixed maps $(g,h)$. Let $G$ denote the population distribution of the validation observables $(Z,\hat Z,S)$.
Assumption (ref) requires the downstream and validation samples to agree on the overlap variables $(\hat Z,S)$. In practice, however, a validation sample constructed as a held-out subset of an upstream training dataset may exhibit covariate shift relative to the downstream population of interest, in which case the implied distributions of $(\hat Z,S)$ need not coincide across the two samples. When the support of the validation sample covers that of the downstream sample, such discrepancies can often be addressed using standard reweighting methods, such as inverse probability weighting, to align the distribution of $X$ and, hence, that of $(\hat Z,S)$ in the validation sample with that in the downstream population. We abstract from this issue here to focus on the core identification problem, and leave a formal treatment of reweighting to future work.
Under Assumption (ref), the identified set for $\theta_0$ is defined as the set of parameters consistent with the observed marginal distributions and the population moment conditions. Formally, let $\mathcal{H}(F,G)$ denote the set of all joint distributions over $(W,Z,\hat Z,S)$ whose marginals on $(W,\hat Z,S)$ and $(Z,\hat Z,S)$ are $F$ and $G$, respectively.
Characterizing the identified set $\Theta_I(F, G)$ computationally is a non-trivial task. Recently, fan_partial_2025 proposed a sharp identification method that relies on conditional optimal transport. However, when the proxy $\hat{Z}$ or the stratifying variables $S$ contain continuous components, their approach requires solving a continuum of optimal transport problems---specifically, one for every realization of $(\hat{Z}, S)$ in the joint support. Solving those conditional optimal transport problems or approximating them by fine discretization could be burdensome in practice.
We develop an alternative identification method that reduces the characterization to an unconditional optimal transport problem while preserving sharpness. The key device, adapted from li_finite_2025, is to introduce auxiliary copies of the overlap variables and move the exact-matching restrictions from the coupling class into moment conditions. Specifically, let $(\hat{Z}',S')$ be auxiliary random vectors and consider the augmented random vector
For the purpose of this representation, we view the downstream sample as providing observations on $(W,\hat{Z},S)$ with distribution $F$, and the validation sample as providing observations on $(Z,\hat{Z}',S')$ with distribution $G$. Let $\mathcal{H}'(F,G)$ denote the unconditional Fr\'{e}chet class of all joint distributions $H'$ over $(W,Z,\hat{Z},S,\hat{Z}',S')$ whose marginal distribution on $(W,\hat{Z},S)$ is $F$ and whose marginal distribution on $(Z,\hat{Z}',S')$ is $G$.
The advantage of this construction is that the almost-sure restrictions $\hat{Z}=\hat{Z}'$ and $S=S'$, which are implicit in $\mathcal{H}(F,G)$ and give rise to a computationally demanding conditional coupling problem, are no longer imposed directly in $\mathcal{H}'(F,G)$. Instead, we reintroduce these exact-matching requirements through additional moment restrictions. This yields the following equivalent representation of the identified set.
Because equations (ref) and (ref) force $\hat{Z}' = \hat{Z}$ and $S' = S$ almost surely under $H'$, any valid distribution in $\mathcal{H}'(F, G)$ that satisfies (ref) and (ref) collapses to an element in $\mathcal{H}(F, G)$, establishing that $\Theta'_I(F, G) = \Theta_I(F, G)$. While theoretically equivalent, this augmented representation successfully isolates the binding coupling constraints into the moment restrictions, allowing us to characterize the identified set using purely unconditional optimal transport.
To state this formally, let $\tilde{q}(W, Z, \hat{Z}, S, \hat{Z}', S'; \theta)$ denote the augmented moment vector stacking the structural moments and the squared proxy differences:
The system of equations (ref)--(ref) can be written compactly as the single vector moment condition $\mathbb{E}_{H'}[\tilde{q}(\cdot; \theta)] = 0$.
Let $d_q$ be the dimension of $q$, and let $\mathbb{B} \subset \mathbb{R}^{d_q + d_{\hat{Z}} + d_S}$ be any compact convex set containing an open neighborhood of the origin (for example, the closed unit ball $\mathbb{B} = \{\lambda : \left\Vert \lambda\right\Vert \le 1\}$). By standard properties of moments, $\mathbb{E}_{H'}[\tilde{q}(\cdot; \theta)] = 0$ holds if and only if
Therefore, for any valid parameter $\theta \in \Theta'_I(F, G)$, there must exist an $H' \in \mathcal{H}'(F, G)$ such that the supremum over $\lambda$ is exactly zero. Consequently, taking the minimum over all joint distributions in $\mathcal{H}'(F, G)$ yields:
Because the supremum of a minimum is always bounded above by the minimum of a supremum (i.e., the weak duality), swapping the order of the operators therefore establishes a necessary condition for identification:
Our main identification theorem establishes that this relationship is, in fact, not only necessary but also sufficient under mild regularity conditions.
Let $\mathcal{X}_d$ and $\mathcal{X}_v$ denote the supports of the downstream observables $(W, \hat{Z}, S)$ and the validation observables $(Z, \hat{Z}, S)$, respectively. Let $\mathcal{W}$ and $\mathcal{Z}$ denote the support of $W$ and $Z$ respectively.
Assumption (ref)(ref) ensures $\mathcal{X}_d$ and $\mathcal{X}_v$ are polish spaces so that probability measures $F$ and $G$ are Radon measures. Assumption (ref)(ref) imposes continuity, which automatically holds when $(w, z)$ are discrete. Finally, Assumptions (ref)(ref) and (ref) provide the necessary envelope bounds to guarantee that the expected moment conditions are finite and continuous with respect to the weak topology on the space of joint distributions.
Under these regularity conditions, we can establish the following result, the proof of which is in Appendix (ref).
Theorem (ref) provides a sharp characterization of the identified set $\Theta_I(F, G)$. For a fixed pair $(\theta,\lambda)$, the inner minimization problem in Theorem (ref) is a standard unconditional optimal transport problem between the two observed marginal distributions $F$ and $G$, with the cost function given by
As we demonstrate in Section (ref), by exploiting the Kantorovich duality of this optimal transport problem, the resulting max-min characterization can be implemented as a convex optimization problem. This greatly simplifies computation and provides a tractable foundation for estimation and inference.
The identification analysis we developed here changes the role of the proxy $\hat{Z}$. Rather than treating $\hat{Z}$ primarily as a plug-in substitute for $Z$, the theorem treats it as a bridge linking the downstream and validation samples. Accordingly, the analysis shifts attention away from the approximation error $Z-\hat{Z}$ itself and toward the compatibility of the two samples through the common variables $(\hat{Z},S)$. This perspective is what allows our results to avoid structural assumptions on the prediction error and to accommodate highly general upstream prediction rules. For the applied researcher, the implication is straightforward: one need not restrict attention to ML methods for which a full theoretical analysis is available, but may instead choose the method that is most effective in the empirical setting at hand.
This perspective also suggests a different criterion for evaluating upstream ML methods. In standard ML practice, predictors are typically judged by how accurately they approximate $Z$. From the standpoint of downstream identification, however, a useful predictor is one that preserves as much of the information in $X$ about $Z$ as possible. In the extreme case in which the conditional distribution of $Z$ given $X$ depends on $X$ only through $(\hat{Z},S)$, the pair $(\hat{Z},S)$ is as informative as $X$ itself for the downstream problem: replacing $X$ with $(\hat{Z},S)$ entails no loss of information about $Z$. A predictor $\hat{Z}$ that closely approximates $Z$ is therefore sufficient, but not necessary, for $\hat{Z}$ to serve as an effective dimension-reduction device.
An immediate implication of this viewpoint is that $Z$ and $\hat{Z}$ need not take values in the same space. For example, when $Z$ is a binary label, $\hat{Z}$ may be the predicted probability that $Z=1$ rather than a hard binary classification. Likewise, when $Z$ is a multinomial group label, $\hat{Z}$ may record the model's top two predicted classes, or even a vector of class scores, rather than only its single best guess. These richer forms of $\hat{Z}$ retain information that would be lost if the proxy were collapsed to a plug-in estimate of $Z$, and Theorem (ref) shows how such information can be incorporated into the downstream analysis.
More generally, viewing $\hat{Z}$ as a dimension-reduction device suggests a natural way to combine predictions from multiple ML methods. Suppose, for example, that two ML prediction rules produce proxies $\hat{Z}_1$ and $\hat{Z}_2$ respectively. One can then form the proxy $\hat{Z}$ as $\hat{Z} = (\hat{Z}_1,\hat{Z}_2)$, and carry out the analysis based on this combined proxy. This can be accommodated because $\hat{Z}$ need not reside in the same space as $Z$. In this way, our framework can exploit complementary information from multiple prediction methods without requiring the researcher to commit to a single proxy ex ante.
This section develops inference for testing whether a candidate parameter value $\theta$ belongs to the identified set $\Theta_I(F,G)$, building on the characterization in Theorem (ref). This inference problem is nonstandard for two reasons. First, the criterion is defined through an optimal-transport problem, and its sample analog exhibits nonstandard asymptotic behavior. Second, the identified-set characterization involves an outer max--min operation, which further complicates inference.
We address these challenges in two steps. First, we rewrite the max--min problem in (ref) using the Kantorovich dual representation, so that both $\lambda$ and the dual functions are jointly determined by a convex optimization problem. We then approximate the infinite-dimensional class of dual functions by sieve spaces, thereby reducing the problem to a finite-dimensional convex program. Details are provided in Subsection (ref). Building on this representation, we develop a testing procedure based on sample splitting and cross-fitting; see Subsection (ref). A key practical advantage of this approach is that the test can be calibrated using an asymptotically pivotal upper bound on the distribution of the test statistic, so that valid critical values can be obtained from standard reference distributions without bootstrap or other simulation-based methods.
Recall that $\mathcal{X}_d$ and $\mathcal{X}_v$ denote the supports of $(W,\hat{Z},S)$ and $(Z,\hat{Z}',S')$, respectively. Let $\mathcal{X}_{w}$ and $\mathcal{X}_{z}$ denote the suport of $W$ and $Z$ respectively. And, let $\mathcal{C}(\mathcal{X}_v)$ denote the class of real-valued continuous functions on $\mathcal{X}_v$ that are integrable with respect to $G$. The following lemma shows that, under Assumption (ref), the max--min criterion in (ref) admits a Kantorovich dual representation.
Combining Lemma (ref) with Theorem (ref) yields the following equivalent characterization of the identified set: \[ \theta\in\Theta_I(F,G) \qquad\text{if and only if}\qquad \mathscr D(F,G;\theta)\le 0. \] Compared with (ref), the dual representation in (ref) is computationally more convenient, provided that the inner minimization over $\mathcal X_v$ can be evaluated efficiently; see Remark (ref). In particular, it transforms the original max--min problem into a joint maximization over $(\lambda,\psi)$. Moreover, the objective function in (ref) is jointly concave in $(\lambda,\psi)$, since the term inside the $\inf\{\dots\}$ is affine in $(\lambda,\psi)$ for each $(z,\hat z',s')$, and the pointwise infimum of affine functions is concave.
To obtain a computationally implementable finite-dimensional approximation, we replace the infinite-dimensional class $\mathcal{C}(\mathcal{X}_v)$ with a sieve space. Let $\mathcal{S}_K(\mathcal{X}_v)$ be a linear sieve space spanned by $K$ known basis functions: \[ \mathcal{S}_K(\mathcal{X}_v) \coloneqq \left\{ \psi:\ \psi(z,\hat z',s') = \beta^{\top}\varphi(z,\hat z',s') \text{ for some } \beta\in\mathbb{R}^K \right\}, \] where $\varphi(z,\hat z',s') \coloneqq \big(\varphi_1(z,\hat z',s'),\dots,\varphi_K(z,\hat z',s')\big)^{\top}$ is a vector of known basis functions. Restricting the dual function $\psi$ in (ref) to this finite-dimensional sieve space yields the sieve approximation
Equation (ref) reduces the problem to a finite-dimensional concave maximization problem over $(\lambda,\beta)$, which is the key to the computational tractability of our inference procedure. The next lemma establishes the relation between $\mathscr{D}(F, G;\theta)$ and $\mathscr{D}_K(F, G;\theta)$.
Lemma (ref) has two important implications for inference. First, (ref) shows that the computationally tractable condition
is a valid necessary condition for $\theta\in\Theta_I(F,G)$. Our inference procedure therefore proceeds by testing (ref). Because (ref) may be weaker than (ref) for fixed $K$, the resulting test is valid but potentially conservative. Second, (ref) shows that the gap between the finite-dimensional criterion $\mathscr{D}_K(F,G;\theta)$ and the population criterion $\mathscr{D}(F,G;\theta)$ vanishes as $K\to\infty$, provided that the underlying sieve space can approximate any continuous function arbitrarily well.
One possible route to inference based on (ref) would be to study the empirical process associated with the objective function in $\mathscr D_K(F,G;\theta)$ and establish a uniform Gaussian approximation over $(\lambda,\beta)$. Such an approach would rely on Donsker-type conditions for the relevant class of objective functions. While feasible in principle, it has two drawbacks. First, the required empirical-process conditions may be restrictive and difficult to verify. Second, the resulting limiting distribution is generally non-pivotal, so critical values typically must be obtained by bootstrap or other simulation-based methods, which can be computationally burdensome in practice.
We therefore develop an alternative inference procedure based on sample splitting and cross-fitting. The primary advantage of this approach is that it places minimal regularity assumptions on the moment function $q(\cdot)$ and yields a test statistic that can be bounded using a pivotal, analytical critical value, thereby avoiding bootstrapping or other simulation methods. The theoretical trade-off is that, because the exact joint asymptotic distribution of the cross-fitted statistics remains unknown, the critical value must account for the least favorable dependence structure. While this minimax bounding renders the resulting inference procedure potentially conservative in finite samples, it guarantees correct asymptotic size control without sacrificing computational tractability.
Assume that, conditional on the trained prediction rule $g$, the downstream sample \[ \{(W_i,\hat Z_i,S_i): i=1,\ldots,n_d\} \] is i.i.d.\ from $F$, the validation sample \[ \{(Z_j,\hat Z'_j,S'_j): j=1,\ldots,n_v\} \] is i.i.d.\ from $G$, and the two samples are mutually independent. Define
Randomly partition each sample into two approximately equal folds, indexed by $m\in\{1,2\}$. Let $\mathcal I^d_m$ and $\mathcal I^v_m$ denote the index sets for fold $m$ in the downstream and validation samples, respectively.
For an arbitrary candidate parameter $\theta \in \Theta$, fix a sequence of basis dimensions $K_n\to\infty$ and a sequence of bounds $c_n\to\infty$. The proposed inference procedure consists of the following steps:
\noindentStep 1. Using only the downstream observations in fold $1$, compute an estimator $\hat{\Omega}_{1,n}$ of the covariance matrix of the vector $q(W, \hat{Z};\theta)$ and an estimator $\hat{\Lambda}_{1,n}$ of the covaraince matrix of the vector $(\hat{Z}^\top, S^\top)^\top$. Define
and \[ \mathbb B_{1,n} \coloneqq \bigl\{\hat\Sigma_{1,n}^{-1/2}\gamma:\ \|\gamma\|\le 1\bigr\}, \] where $\hat\Sigma_{1,n}^{-1/2}$ denotes an inverse square root (or generalized inverse square root, if needed) for $\hat{\Sigma}_{1,n}$. This step provides an empirical normalization of the components of the augmented moment vector $\tilde{q}$.
\noindentStep 2. Using only the downstream and validation observations in fold $1$, solve the empirical analogue of (ref):
where $\mathbb C_n \coloneqq \{\beta\in\mathbb R^{K_n}:\ \|\beta\|\le c_n\}$ and $\varphi(\cdot)$ denotes the vector of $K_n$ sieve basis functions defined on $\mathcal X_v$. Restricting the coefficients to the compact set $\mathbb C_n$ explicitly ensures a finite empirical maximizer for $\hat{\beta}_{1,n}$ always exists and controls the complexity of the sieve approximation.
\noindentStep 3. Hold $(\hat\lambda_{1,n},\hat\beta_{1,n})$ fixed and evaluate the criterion on fold $2$. Define
The resulting held-out estimate of the sieve dual criterion is \[ \widehat{\mathscr D}_2(\theta) \coloneqq \frac{1}{|\mathcal I^d_2|}\sum_{i\in\mathcal I^d_2}\hat u_i + \frac{1}{|\mathcal I^v_2|}\sum_{j\in\mathcal I^v_2}\hat v_j. \]
\noindentStep 4. Conditional on fold $1$, the quantity $\widehat{\mathscr D}_2(\theta)$ is the sum of two independent sample averages, one from the downstream sample and one from the validation sample. Its conditional variance can therefore be estimated by \[ \hat V_2(\theta) \coloneqq \frac{\hat\sigma_{u,2}^2}{|\mathcal I^d_2|} + \frac{\hat\sigma_{v,2}^2}{|\mathcal I^v_2|}, \] where $\hat\sigma_{u,2}^2$ and $\hat\sigma_{v,2}^2$ are the sample variances of $\{\hat u_i:i\in\mathcal I^d_2\}$ and $\{\hat v_j:j\in\mathcal I^v_2\}$, respectively. Since $\hat V_2(\theta)$ shrinks asymptotically at rate $1/\underline{n}$, $\widehat{\mathscr D}_2(\theta)$ should be scaled by $\sqrt{\underline{n}}$. To guard against degeneracy when $(\hat\lambda_{1,n},\hat\beta_{1,n})=(0,0)$, let $\epsilon>0$ be a small fixed constant. The fold-2 test statistic is then
\noindentStep 5. Repeat Steps 1--4 after swapping the roles of the two folds. This yields a second statistic,
where $\widehat{\mathscr D}_1(\theta)$ and $\hat V_1(\theta)$ are computed on fold $1$ using the optimizer $(\hat\lambda_{2,n},\hat\beta_{2,n})$ obtained from fold $2$. We then aggregate the two fold-specific statistics by \[ T(\theta)\coloneqq \max\{T_1(\theta),T_2(\theta)\}. \] The associated $p$-value is \[ p(\theta) \coloneqq \min\bigl\{1,\ 2\bigl[1-\Phi\bigl(T(\theta)\bigr)\bigr]\bigr\}, \] where $\Phi(\cdot)$ denotes the standard normal distribution function. We reject the null hypothesis \[ H_0:\ \theta\in\Theta_I(F,G) \] at nominal significance level $\alpha$ whenever $p(\theta)<\alpha$.
Some remarks are in order.
To establish asymptotic size control formally, we impose the following primitive regularity conditions.
Assumption (ref) consists of mild regularity conditions. Assumption (ref)(ref) requires only that both the downstream and validation sample sizes increase to infinity and that no fold becomes asymptotically negligible. In particular, it places no restriction on the relative magnitudes of $n_d$ and $n_v$. This flexibility is important in applications where the validation sample may be either much smaller or much larger than the downstream sample. Assumptions (ref)(ref)--(ref) are standard regularity conditions. In Assumption (ref)(ref), the uniform bound on the sieve basis functions is largely a normalization: by rescaling the basis, one can always impose a common sup-norm bound, such as $C=1$. Equation (ref) is a rate condition, which restricts the growth of the effective complexity of the sieve space. Intuitively, because the sieve coefficients are constrained by $\|\beta\|\le c_n$, the relevant envelope of the sieve approximation on $\mathcal X_v$ is of order $c_n\sqrt{K_n}$, and (ref) requires this envelope to grow slowly relative to sampling variation.
Theorem (ref) establishes pointwise asymptotic size control for the proposed procedure. This result can be strengthened to uniform size control over suitably restricted classes of pairs $(F,G)$, for example, classes with uniformly bounded supports and covariance matrices whose eigenvalues are uniformly bounded away from zero.
Because the critical value is calibrated against the least favorable joint distribution of $(T_1(\theta),T_2(\theta))$ consistent with their marginal standard normal limits, the resulting procedure may be conservative, so the inequality in (ref) can be strict. This conservativeness is the price of maintaining computational tractability under weak regularity conditions. Finally, both the proposed method and Theorem (ref) concern inference on the full parameter vector. A conservative procedure for subvector inference can be obtained by projecting the full-vector confidence region onto the subvector of interest. A sharper treatment of subvector inference is left for future work.
In this section, we conduct Monte Carlo simulations to assess the finite-sample performance of the proposed partial-identification and inference procedures. These simulations are designed to illustrate several key features of the methodology. First, we examine the size control of our procedure and compare it with that of a naive plug-in approach that ignores measurement error in the ML-generated proxy. Second, we illustrate the power of the proposed cross-fitted test and how it varies with sample size. Third, we show how incorporating a stratifying variable $S$ can tighten the identified bounds on the structural parameters. Fourth, we compare the results obtained using discrete and continuous proxies $\hat{Z}$. We also examine the role played by the number of sieve basis functions. This final exercise further illustrates how our framework accommodates settings in which the latent target $Z$ and its proxy $\hat{Z}$ take values in different spaces. The next subsection describes the data-generating process used throughout these simulations.
We consider a design motivated by Example (ref), in which the downstream researcher is interested in the regression model
where $Z$ is a binary regressor that is unobserved in the downstream sample. Instead, the researcher observes a machine-learned proxy $\hat Z = g(X)$ constructed from a high-dimensional predictor vector $X$. The data-generating process (DGP) is designed so that $(Z,C)$ are exogenous in the structural equation, whereas $X$ contains both exogenous and endogenous components. Consequently, although the true regressor $Z$ is exogenous, the proxy $\hat Z$ may be endogenous. In the downstream sample, the researcher observes $(Y,C,X)$. For simplicity, we take $C$ to be scalar and set the true parameter vector to $\theta_0 = (1,1)^\top$.
Specifically, let $(V,U,\xi,\nu)$ be mutually independent random variables, each distributed as $N(0,1)$. We then generate the remaining variables as follows:
By construction, $(Z,C)$ are exogenous in the downstream regression. At the same time, because $U$ enters both $\varepsilon$ and the endogenous block $X^{en}$, a flexible predictor $\hat Z = g(X)$ may inherit endogeneity through its dependence on $X^{en}$.
The prediction rule $g(\cdot)$ is estimated from an i.i.d.\ upstream training sample of $(Z,X)$ of size $n_{\text{train}}=5000$. This sample is generated independently of both the downstream sample and the auxiliary validation sample used to assess the proxy. We simulate the upstream training sample only once and keep both the realized sample and the resulting fitted rule $g(\cdot)$ fixed across all Monte Carlo replications. This design mirrors the theoretical framework in Sections (ref) and (ref), in which the downstream econometric analysis is conducted conditional on a pre-trained prediction rule.
In the simulations, we estimate $g(\cdot)$ using an $\ell_1$-penalized logistic regression (logistic LASSO). For a given penalty parameter $\gamma>0$, the estimator solves
where $\|\tau_1\|_1=\sum_j |\tau_{1,j}|$ is the $\ell_1$ norm of the slope vector, and $\Lambda(\cdot)$ denotes the standard logistic cumulative distribution function,
To choose $\gamma$, we perform five-fold cross-validation over a dense grid of candidate values. Specifically, we randomly partition the training indices into five equal-sized folds, denoted by $\mathcal I_1,\dots,\mathcal I_5$. For each candidate $\gamma$ and each fold $k\in\{1,\dots,5\}$, we fit the model on the remaining four folds to obtain fold-specific estimates $\big(\hat{\tau}_0^{(-k)}(\gamma),\hat{\tau}_1^{(-k)}(\gamma)\big)$ and evaluate predictive performance on the held-out fold $k$. We then select the tuning parameter
where \[ \hat{P}_i^{(-k)}(\gamma) = \Lambda\big(\hat{\tau}_0^{(-k)}(\gamma)+X_i^\top \hat{\tau}_1^{(-k)}(\gamma)\big) \] is the predicted probability for observation $i$ when fold $k$ is held out.
After selecting $\gamma^*$, we refit the logistic LASSO on the full upstream training sample to obtain the final estimates $\hat{\tau}_0=\hat{\tau}_0(\gamma^*)$ and $\hat{\tau}_1=\hat{\tau}_1(\gamma^*)$. The resulting predicted probability for an observation with covariates $X$ is
In the baseline designs, we convert this score into a binary proxy via
The resulting classification accuracy, measured by $\mathbb{P}(Z=\hat{Z})$, is approximately $98.42\%$ when $\sigma_\eta(X^{ex})=s_L$, $88.57\%$ when $\sigma_\eta(X^{ex})=s_M$, and $73.12\%$ when $\sigma_\eta(X^{ex})=s_H$.
We conduct the Monte Carlo study using $10{,}000$ replications. In each replication, the downstream sample consists of $n_d$ observations on $(Y,C,\hat Z)$, while the auxiliary validation sample consists of $n_v$ observations on $(Z,\hat Z)$. The baseline simulation design does not include a stratifying variable, so $S$ is omitted throughout this subsection. The role of stratification is examined separately in Subsection (ref).
For each replication, we implement the inference procedure developed in Section (ref). Throughout the simulations, we set $c_n=\underline{n}^{2/5}$. We also experimented with the alternative choice $c_n=\underline{n}^{1/3}$ and obtained very similar results. Because $(Z,\hat Z)$ is discrete and $S$ is absent, the inner $\inf$ problem in (ref) is easy to solve numerically. Moreover, in this setting the sieve approximation is exact. We therefore take the basis vector $\varphi(Z,\hat Z)$ to be the four cell indicators corresponding to the support of $(Z,\hat Z)$:
An attractive feature of our framework is its robustness to asymmetric sample sizes in the downstream and validation datasets. To illustrate this property, we examine the rejection rate of the test under four sample-size configurations: \[ (n_d,n_v)\in\{(500,500),\ (10000,10000),\ (10000,500),\ (500,10000)\}. \] The last two designs explicitly examine the performance of the procedure under substantial sample-size asymmetry. For each configuration, we report rejection frequencies at nominal significance levels of $10\%$, $5\%$, and $1\%$ when testing the true parameter value $\theta_0=(1,1)^\top$. For reference, on a Mac with an M4 Max CPU (16 CPU cores), the case $(n_d,n_v)=(10000,10000)$ takes about 90 seconds for all $10{,}000$ replications combined, or 9 milliseconds per replication.
A second implication of the theory is that the asymptotic validity of the procedure does not depend on the predictive accuracy of the upstream ML method. To illustrate this feature, we report size results for the three levels of prediction noise, $(s_L,s_M,s_H)$, introduced in the data-generating process. Although higher prediction noise naturally widens the identified set, it should not compromise the size control of the test.
As a benchmark, we also report the simulated rejection frequencies of a conventional plug-in $F$-test based on the ordinary least squares regression of $Y$ on $(\hat Z,C)$. This naive procedure treats $\hat Z$ as if it were observed without error and therefore ignores both the prediction error in $\hat Z$ and the endogeneity that $\hat Z$ may inherit from the endogenous components of $X$.
Table (ref) reports the results. Two conclusions emerge. First, our procedure controls size well across all sample-size configurations and all levels of prediction noise. Although conservative by construction, the conservativeness is mild in our simulation. For example, when $(n_d,n_v)=(10000,10000)$, the empirical rejection frequency at the $5\%$ level ranges from $4.2\%$ to $4.5\%$, depending on the prediction-noise design, while at the $1\%$ level it ranges from $0.9\%$ to $1.0\%$. Second, the naive plug-in $F$-test fails to control size. It exhibits substantial over-rejection whenever prediction noise is moderate or high, for all sample-size configurations. Even in the low-noise design, where the prediction accuracy is approximately $98.42\%$, over-rejection is still noticeable when $n_d=10{,}000$.
Figure (ref) visualizes the finite-sample power of the proposed procedure through heatmaps over the two-dimensional parameter space $(\theta_1,\theta_2)$. For each design, we evaluate the test on a grid of $10{,}000$ parameter values and repeat the experiment over $500$ Monte Carlo replications. At each grid point, we record the fraction of replications in which the null hypothesis is rejected at the $5\%$ significance level. The resulting heatmap therefore summarizes the empirical rejection probability across the parameter space. The purple end of the color scale corresponds to a rejection frequency of $0\%$, meaning that the point lies in the $95\%$ confidence set in every replication, whereas the yellow end corresponds to a rejection frequency of $100\%$, meaning that the point is excluded from the $95\%$ confidence set in every replication.
Figure (ref) reports confidence sets for nine different sample-size combinations, arranged in a $3\times 3$ panel. Across rows, the downstream sample size $n_d$ increases from $500$ to $2000$, while across columns the validation sample size $n_v$ increases from $500$ to $2000$. As the figure shows, increasing either $n_d$ or $n_v$ sharpens the rejection surface and yields a more informative confidence set. The gains in power are most pronounced, however, when both sample sizes increase simultaneously. This pattern is consistent with our asymptotic theory, under which the precision of the procedure is governed by $\underline{n}=\min(n_d,n_v)$.
To illustrate how predictive performance affects the informativeness of our procedure, Figure (ref) plots the empirical $95\%$ confidence sets under the three prediction-noise designs $s\in\{s_L,s_M,s_H\}$. The heatmaps are constructed in the same manner as in Figure (ref) with sample sizes fixed at $(n_d, n_v) = (1000, 1000)$. As expected, lower prediction noise and hence better predictive performance yield more informative inference, producing smaller confidence sets.
\thispagestyle{plain}
\thispagestyle{plain}
To illustrate the role of stratification, we modify the baseline design by allowing the scale of the latent shock in (ref) to vary with an observable component of $X$. Specifically, we let
where $X_1$ denotes the first component of $X$.\footnote{To preserve the exogeneity of $Z$ in the structural equation, we take the first component of $X$ to belong to the exogenous block in the data-generating process. This information is used only to define the stratification variable and is not otherwise exploited in the inference procedure.} We then define the stratification variable by
Under this design, the variance of the prediction noise varies across the two strata indexed by $S$. The setup is meant to mimic empirical environments in which the latent target variable $Z$ is easier to predict for some subpopulations than for others.
We compare inference with and without incorporating $S$ under two configurations of $(s_1,s_2)$. We begin with the benchmark case $(s_1,s_2)=(1,1)$. In this design, the prediction noise is homoskedastic, so stratification does not provide any additional information for identification. Accordingly, incorporating $S$ should not tighten the identified set for $\theta$. This is what we see in the first row of Figure (ref): adding $S$ delivers essentially no gain. If anything, because conditioning on $S$ increases the dimensionality of the problem, the finite-sample performance is slightly worse when stratification is included.
We next consider a highly heterogeneous design with $(s_1,s_2)=(1,20)$. In this case, conditional on $S=0$, the covariates $X$ are less informative about $Z$, so the classification accuracy of $\hat Z$ in that stratum is much lower. As shown in the lower panels of Figure (ref), incorporating stratification yields noticeably tighter inference. Although $S$ does not improve the predictive performance of $\hat Z$ itself, it identifies the subpopulations in which the proxy is informative and those in which it is not. More broadly, this exercise illustrates a central theme of our framework. The pair $(\hat Z,S)$ serves as a bridge that transfers information about $Z$ from the validation sample to the downstream moment conditions. As a result, even when $S$ does not directly enhance prediction, it can still improve inference by helping to characterize the conditional distribution of $Z$ given $(\hat Z,S)$.
To illustrate how the inference procedure applies to a continuous proxy and to assess the effect of sieve approximation, we modify the baseline design by replacing the binary proxy $\hat{Z}$ in (ref) with the continuous proxy
Because $Z$ remains binary whereas $\hat{Z}$ is continuous, this specification also shows that the procedure can accommodate settings in which $Z$ and $\hat{Z}$ take values in different spaces.
Since $\hat{Z}$ is continuous, the basis functions in (ref) are no longer appropriate. We therefore use the following order-$d$ sieve basis, which has dimension $K = 2(d + 1)$:
This basis can be viewed as the tensor product of the indicator basis $(\mathbbm{1}(Z = 1), \mathbbm{1}(Z = 0))$ and the polynomial basis $(1, \hat{Z}, \hat{Z}^2, \ldots, \hat{Z}^{d})$. We set $c_n = (\min\{n_d, n_v\})^{2/5}$.
Table (ref) reports empirical rejection probabilities for sieve orders $d = 1,2,\ldots,5$ across different sample-size configurations. In this exercise, the prediction noise is fixed at the homoskedastic medium level $s_M = 1.0$. As Lemma (ref) suggests, replacing the infinite-dimensional dual space with a finite-dimensional sieve approximation preserves the validity of the underlying necessary condition. Accordingly, the proposed procedure exhibits good size control across all polynomial orders. At the same time, the polynomial order of the sieve affects the degree of conservativeness. When the sieve order is low, the sieve space may approximate the optimal dual function $\psi$ in (ref) only coarsely, so the resulting approximation error can be non-negligible. As discussed in Section (ref), this approximation error leads to a gap between the actual dual $\mathscr{D}(F, G;\theta)$ and its finite-dimensional approximation $\mathscr{D}_K(F, G;\theta)$, therefore making the test more conservative. This pattern is visible in Table (ref), especially when the polynomial order is $1$.
As the polynomial order increases, the sieve approximation improves and the test becomes less conservative. This is borne out in the simulation results: moving from a first-order to a second-order sieve noticeably reduces conservativeness, and the third-order and the fourth-order sieve yield some further improvements. However, we observe diminishing gain as the order of polynomial increases. This suggests that, once the sieve is sufficiently flexible, the remaining conservativeness is driven primarily by the same force that underlies the results in Table (ref). Specifically, our critical values are calibrated using the least favorable joint distribution consistent with the known standard normal marginals of the fold-specific statistics. This calibration yields a tractable pivotal critical value under relatively weak regularity conditions, but it also introduces an inherent degree of conservativeness.
Finally, we compare the informativeness of the continuous proxy in (ref) with that of the binary proxy in (ref). As discussed in Remark (ref) and Section (ref), the proxy $\hat{Z}$ is used as a bridge linking the downstream and validation samples. The more information $\hat{Z}$ preserves about the relationship between the high-dimensional input $X$ and the latent target $Z$, the tighter the resulting identified set and confidence set should be. In the present design, the binary proxy is obtained from the continuous proxy by comparing it with $0.5$ and is therefore a coarsening of it. The continuous proxy thus contains more information than the binary proxy, suggesting that it should yield tighter confidence sets.
To examine this implication, we simulate confidence sets for both the binary proxy and the continuous proxy, the latter using polynomial sieve bases of orders $d=1,\ldots,5$. Figure (ref) reports the resulting confidence sets for the case $n_d=n_v=1000$, while Figure (ref) reports the corresponding results for the larger sample size $n_d=n_v=5000$. In both exercises, the prediction noise is fixed at the homoskedastic medium level $s_M=1.0$. As in Figure (ref), we evaluate the test on a grid of $10{,}000$ candidate parameter values and repeat the exercise over $500$ Monte Carlo replications. The results confirm the conjecture. Even with a first-order polynomial sieve, the continuous proxy yields tighter confidence sets than the binary proxy. This gain becomes more pronounced at the larger sample size. Moreover, the gains from using higher-order polynomial sieves are also more evident when the sample size is larger, since the reduction in sieve approximation error then more clearly outweighs the cost of estimating a higher-dimensional set of sieve coefficients.
\thispagestyle{plain}
\thispagestyle{plain}
We now consider an extension in which $\hat{Z}$ is vector-valued, but the available validation data do not jointly observe the full vector $(Z,\hat{Z},S)$. Instead, the researcher may have access to several validation samples, each of which contains information only on a subvector of $(Z,\hat{Z},S)$. In this case, the validation data identify a collection of lower-dimensional marginal distributions rather than the full joint distribution of $(Z,\hat{Z},S)$. As illustrated in Example (ref), this situation can arise when $\hat{Z}$ consists of several ML-generated proxies, each constructed separately rather than estimated jointly as a single vector-valued predictor. We show that this setting can be accommodated by a straightforward extension of the framework developed above.
Formally, partition $(Z,\hat{Z},S)$ into $L$ subvectors, \[ (Z_1,\hat{Z}_1,S_1),\ldots,(Z_L,\hat{Z}_L,S_L). \] Let $G_l$ denote the population distribution of the validation observables $(Z_l,\hat{Z}_l,S_l)$ for $l=1,\ldots,L$. The compatibility condition in Assumption (ref) is then replaced by the following.
\begin{assumption1prime}[Compatibility of Multiple Marginal Distributions] Let $H_0$ denote the true joint distribution of $(W,Z,\hat Z,S)$. Assume that (i) the marginal distribution of $(W,\hat Z,S)$ under $H_0$ is $F$, and (ii) for each $l=1,\ldots,L$, the marginal distribution of $(Z_l,\hat Z_l,S_l)$ under $H_0$ is $G_l$. \end{assumption1prime}
Define $\mathcal{H}_L\!\left(F,(G_l)_{l=1}^L\right)$ to be the set of all joint distributions $H$ over $(W,Z,\hat Z,S)$ whose marginal distribution on $(W,\hat Z,S)$ is $F$ and whose marginal distribution on $(Z_l,\hat Z_l,S_l)$ is $G_l$ for each $l=1,\ldots,L$. This leads to the following notion of identified set.
As in the baseline case, we use the same decoupling device. Let $(\hat Z',S')$ denote auxiliary random vectors, and for each $l=1,\ldots,L$, let $(\hat Z_l',S_l')$ denote the subvector corresponding to $(\hat Z_l,S_l)$. Define $\mathcal{H}_L'\!\left(F,(G_l)_{l=1}^L\right)$ to be the set of all joint distributions over $(W,Z,\hat Z,S,\hat Z',S')$ whose marginal distribution on $(W,\hat Z,S)$ is $F$ and whose marginal distribution on $(Z_l,\hat Z_l',S_l')$ is $G_l$ for each $l=1,\ldots,L$. The regularity conditions in Assumption (ref) then needs to be replaced by the following assumption.
\begin{assumption2prime}[Regularity conditions] For any parameter $\theta\in \Theta$, assume that Assumptions (ref)(ref)-(ref) hold, and that there exist non-negative continuous functions $b_1 \in \mathcal{L}^1(F)$ and $b_{2,l} \in \mathcal{L}^1(G_l)$ for each $l=1,...,L$ such that the structural moment function is bounded by an additively separable envelope: \[ \|q(w, z; \theta)\| \le b_1(w) + \sum_{l=1}^L b_{2,l}(z) \quad \text{for all } w \in \mathcal{W} \text{ and } z \in \mathcal{Z}. \] \end{assumption2prime}
We then obtain the following counterpart to Theorem (ref).
The inner minimization problem in (ref) is now a multi-marginal optimal transport problem. For each $l=1,\ldots,L$, let $\mathcal X_{v,l}$ denote the support of $(Z_l,\hat Z_l,S_l)$. This problem admits the dual representation
Consequently, the inequality in (ref) is equivalent to
To obtain a computationally tractable approximation, we again replace each infinite-dimensional dual function class with a sieve space. For simplicity, suppose that each sieve has dimension $K$, and define \[ \mathcal S_{l,K}(\mathcal X_{v,l}) \coloneqq \left\{ \psi:\ \psi(z_l,\hat z_l',s_l')=\beta_l^\top \varphi_l(z_l,\hat z_l',s_l') \text{ for some } \beta_l\in\mathbb R^K \right\}, \] where $\varphi_l(z_l,\hat z_l',s_l')$ is a vector of $K$ known basis functions on $\mathcal X_{v,l}$. Restricting each dual function $\psi_l$ in (ref) to the corresponding sieve space $\mathcal S_{l,K}(\mathcal X_{v,l})$ yields the finite-dimensional approximation
The inference procedure developed in Section (ref) extends straightforwardly to this setting. In essence, one applies the same sample-splitting and cross-fitting strategy as before. One fold is used to solve for the optimizer $(\hat{\lambda},\hat{\beta}_1,\ldots,\hat{\beta}_L)$ in the multi-marginal sieve problem, while the other fold is used to evaluate the corresponding test statistic given $(\hat{\lambda},\hat{\beta}_1,\ldots,\hat{\beta}_L)$. The roles of the two folds are then reversed to obtain a second fold-specific statistic, and the two are aggregated using their maximum. Because these modifications are largely mechanical and introduce no new conceptual issues, we omit the details for brevity.
The central conceptual innovation of our approach is to treat the ML-generated proxy not as a naive plug-in substitute for the unobserved latent variable, but as a dimension-reduction device that empirically bridges the downstream sample and an auxiliary validation sample. To operationalize this idea, we establish new identification results for data combination and propose a novel, resampling-free inference procedure based on sample splitting and cross-fitting, both of which are of independent econometric interest.
For the applied researcher, this perspective fully decouples the validity of econometric inference from the statistical consistency of the upstream ML algorithm. Researchers are free to deploy the most flexible prediction algorithms available---whether they output binary classifications, continuous predicted probabilities, or multiple competing scores---confident that the resulting econometric bounds will remain asymptotically valid. Consequently, empirical practitioners can select whichever ML procedure proves most effective in practice for their specific application, unconstrained by the need for formal statistical guarantees from the prediction step.
For researchers developing upstream ML methods, our framework motivates an alternative criterion for evaluating predictive algorithms. When the ultimate goal is downstream empirical analysis, an optimal prediction rule need not minimize the discrepancy between the prediction and the true target variable. Instead, it should serve as an effective dimension-reduction tool that preserves as much information as possible from the raw, unstructured data regarding the latent target. Although not formally explored here, the design of ML algorithms explicitly tailored for this information-preservation objective represents a promising direction for future research.
Our framework also leaves several important theoretical questions open for future research. First, our baseline identification analysis presumes that the downstream target population and the upstream validation sample exhibit compatible marginal distributions over the overlap variables. In empirical practice, however, the validation data---often a held-out subset of an upstream training corpus---may suffer from covariate shift relative to the downstream population of interest. Formally integrating reweighting techniques, such as inverse probability weighting, into the optimal transport characterization to explicitly account for sample selection and distribution shift remains a highly relevant theoretical extension. Second, the cross-fitted inference procedure developed in this paper provides valid tests for the full structural parameter vector. While applied researchers can construct valid, conservative confidence sets for individual parameters of interest (e.g., a specific treatment effect) via projection, this approach sacrifices statistical power. Developing a localized, profiling-based inference theory that achieves sharp, non-conservative bounds for subvectors within this optimal transport framework is a challenging but vital next step.
More broadly, incorporating ML-generated proxies into empirical economic analysis is an important and increasingly relevant topic. By combining tools from optimal transport, partial identification, and cross-fitting, the framework developed here offers a tractable and versatile way to do so. It allows researchers to use complex machine-learned measures while remaining agnostic about the fine details of the upstream prediction step under mild data requirement. We hope that this perspective proves useful both for applied work that relies on ML-based measurement and for future research at the intersection of econometrics and machine learning.
\printbibliography