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.
90,645 characters · 21 sections · 93 citation commands
Robust Design and Evaluation of Predictive Algorithms under Unobserved Confounding
{\singlespacing \thispagestyle{empty} \setcounter{page}{0}
}
Predictive algorithms inform high-stakes decisions in pretrial release, consumer lending, health care, and many other domains. Designing and evaluating such tools is often complicated by the “selective labels” problem: outcomes are observed only for units selected by past decision makers. For example, pretrial risk assessments predict the likelihood defendants would commit pretrial misconduct if they were released before trial; but we only observe whether a defendant committed misconduct if a judge decided to release them KLLLM(18), Rambachan(21), AngelovaEtAl(22). Consumer credit scores predict the likelihood of default if an applicant were granted a loan; but we only observe whether an applicant defaulted if the financial institution approved them and they accepted the loan terms BlattnerNelson(21), CostonEtAl(21)-Rashomon.\footnote{The selective labels problem arises when analyzing predictive algorithms in many other empirical settings as well, such as medical diagnosis MullainathanObermeyer(21), child protective services ChouldechovaEtAl(18), hiring LiRaymondBergman(20), education BergmanEtAl(21), and tax audits BlackEtAl(22)-TaxAudits, BattagliniEtAl(22), ElzaynEtAl(23).}
The selective labels problem creates an identification challenge under unobserved confounding--when selected and unselected units differ in ways not captured by observed covariates. Judges may release defendants based on courtroom interactions not in administrative data; loan applicants may accept offers based on competing credit options unavailable to researchers. Ignoring unobserved confounding leads to inaccurate predictions and misleading evaluations of predictive algorithms.\footnote{While there exists methods for designing and evaluating predictive algorithms in computer science from selectively observed data Coston(20)-Evaluation, CostonEtAl(21)-Rashomon, MishlerEtAl(21), mishler2021fade, guerdan2023counterfactual, they tackle this challenge assuming that the historical selection decision was unconfounded conditional on observed covariates.} Existing solutions are limited. Applied researchers rely on informal heuristics---imputing missing outcomes using proxy variables ChouldechovaEtAl(18), BlattnerNelson(21), MullainathanObermeyer(21) or applying ad-hoc adjustments to observed outcomes among selected units. Recent work uses instrumental variables but typically focuses on specific estimands lakkaraju2017selective, ArnoldDobbieHull(20)-discrimination_alg or predictive algorithms that are linear in covariates ArnoldEtAl(24).
We develop a framework for robust design and evaluation of predictive algorithms under unobserved confounding. We bound how much outcomes may differ between selected and unselected units with the same observed characteristics. These bounds formalize common empirical strategies --- including proxy outcomes and instrumental variables --- and partially identify key predictive estimands such as conditional likelihoods, mean square error, calibration, and true/false positive rates. Our estimation approach works across multiple bounding strategies and performance measures, flexibly incorporating machine learning methods.
We consider data $(X_i, D_i, Y_i)$ for $i = 1, \hdots, n$ drawn i.i.d. from a joint distribution $\mathbb{P}(\cdot)$, where $X_i$ is a vector of covariates and $D_i \in \{0, 1\}$ is a selection decision by past decision makers (e.g., judges or loan applicants). The observed outcome satisfies $Y_i = D_i \times Y_i^*$ for some true outcome $Y_i^* \in \{0, 1\}$; for example, $Y_i^*$ is whether a defendant would commit pretrial misconduct if released before trial or whether an applicant would default if granted a loan.
Researchers have two objectives. First, design predictive algorithms $s(X_i)$ that predict the probability $Y_i^* = 1$ given the covariates $X_i$ and construct data-driven decision rules. Second, evaluate existing predictive algorithms through performance measures, such as overall performance given by $\mbox{perf}(s; \beta) = \mathbb{E}[\beta_0(X_i; s) + \beta_1(X_i; s) Y_i^*]$, and class-specific performance given by $\mbox{perf}_{+}(s; \beta) = \mathbb{E}[\beta_0(X_i; s) \mid Y_i^* = 1]$ and $\mbox{perf}_{-}(s; \beta) = \mathbb{E}[\beta_0(X_i; s) \mid Y_i^* = 0]$, where $\beta_0(X_i; s), \beta_1(X_i; s)$ are researcher-specified functions. Alternative choices recover popular diagnostics such as mean square error, true positive rate, and false positive rate. Comparing these measures across sensitive attributes evaluates fairness properties MitchellEtAl(19), barocas-hardt-narayanan.
Rather than assuming unconfounded selection, we bound the confounding function $\delta(x) := \mathbb{P}(Y_i^* = 1 \mid D_i = 0, X_i = x) - \mathbb{P}(Y_i^* = 1 \mid D_i = 1, X_i = x)$ using researcher-specified bounding functions that may depend on covariates and identified nuisance parameters. The confounding function measures how much unselected units differ from selected units with the same characteristics; bounding it formalizes assumptions about unobserved confounding severity without requiring a full selection model. Under these bounds, the conditional likelihood and predictive performance measures are partially identified.
Our framework encompasses popular strategies for addressing selective labels as special cases. When researchers observe proxy outcomes correlated with $Y_i^*$ --- such as default on other credit products BlattnerNelson(21) or long-term health outcomes MullainathanObermeyer(21), ChanGentzkowYu(21) --- this implies specific bounds. When researchers have instrumental variables like randomly assigned decision makers lakkaraju2017selective, KLLLM(18), this implies alternative bounds. Our framework accommodates other plausible identifying assumptions.
We apply this framework to two objectives. First, we develop estimators for the bounds on the conditional likelihood $P(Y_i^* = 1 \mid X_i = x)$ (e.g., the pretrial misconduct rate or default rate given covariates), enabling design of predictive algorithms robust to unobserved confounding. Building on advances in causal machine learning NieWager(20), kennedy2020optimal, we propose pseudo-outcome regression estimators with sample splitting: estimate nuisance parameters on one fold, construct influence-function-based pseudo-outcomes, then regress them on covariates using machine learning on a second fold. A key innovation is an oracle inequality extending kennedy2020optimal's pointwise analysis to integrated mean square error: our feasible estimators approximate an infeasible oracle up to a smoothed, doubly-robust product of first-stage errors. This applies to any regression method satisfying a mild stability condition, including random forests, kernel methods, and series estimators. We leverage this to characterize integrated mean square error convergence rates and establish regret bounds for plug-in decision rules, connecting estimation guarantees to decision-theoretic performance.
Second, we develop estimators for bounds on the performance of a given predictive algorithm $s(\cdot)$ under unobserved confounding. For overall performance measures $\mbox{perf}(s; \beta)$ --- which nest mean square error, calibration, and precision --- we construct debiased estimators using crossfitting and debiased machine learning RobinsEtAl(08), ZhengvanderLaan(11), ChernozhukovEtAl(18)-DML. Under mild regularity conditions, these estimators are $\sqrt{n}$-consistent and jointly asymptotically normal. For class-specific performance $\mbox{perf}_{+}(s; \beta)$ and $\mbox{perf}_{-}(s; \beta)$ --- such as true and false positive rates --- we show sharp bounds correspond to optimal values of linear fractional programs with nuisances in both the objective and constraints. We estimate using two folds: estimate nuisances and bounding functions on the first fold, then solve the sample optimization program on the second fold. We establish partial double robustness: estimation error decomposes into sampling error, a doubly-robust product for objective nuisances, and constraint estimation errors entering through root mean square error. Proving this exploits the optimization structure of linear fractional programs.\footnote{Since linear fractional programs arise in several related sensitivity analysis frameworks in causal inference AronowLee(2013), MiratrixWagerZubizaretta(18), KallusMaoZhou(18)-IntervalEstimation, ZhaoSmallBhattacharya(19), KallusZhou(21), our results may be broadly applicable beyond the selective labels setting.} Our estimation strategies accommodate multiple bounding strategies, flexibly incorporating domain-specific knowledge while maintaining statistical guarantees.
As an empirical illustration, we apply our framework to administrative data on personal loan applications from a large Australian bank. First, we construct confounding-robust credit risk scores under varying confounding assumptions. The resulting algorithms produce substantially different risk rankings that could meaningfully affect credit access. Second, inspired by recent work in consumer finance BlattnerNelson(21), FusterEtAl(22)predictably, DiMaggioEtAl(22)-Fintech, we evaluate a benchmark credit score's performance across income groups. An analysis ignoring unobserved confounding suggests the score is more accurate (in mean square error) for higher-income applicants. However, we show this disparity disappears at plausible levels of unobserved confounding, illustrating how researchers can assess the robustness of fairness conclusions.
Our framework prioritizes breadth and applicability. Rather than bespoke methods for each combination of bounding strategy and performance measure, we provide a unified approach for a variety of assumptions and performance measures. We deliver practical estimation procedures that enable researchers to robustly design and evaluate predictive algorithms under unobserved confounding.
\paragraph{Related Work:} Our analysis builds on foundational work on partial identification with selective labels Manski(89), Manski(90), Manski(94)-SampleSelectionBounds. We generalize Manski(03)'s “approximate mean independence” assumption by incorporating covariates and nuisance parameters into the bounding functions, enabling flexible formalization of empirical heuristics using proxy outcomes and other strategies. See also diaz2013sensitivity, luedtke2015statistics, diaz2018sensitivity for similar bounding assumptions without covariates. We characterize identified sets for novel estimands arising in predictive algorithm evaluation.
Instrumental variable strategies are widely used to address unobserved confounding in this setting. lakkaraju2017selective and KLLLM(18) evaluate specific performance measures by returning a single point in the identified set. ArnoldDobbieHull(20)-discrimination_alg, ArnoldEtAl(24) extrapolate across instrument values using identification-at-infinity arguments. We instead derive covariate-dependent bounds following Manski(94)-SampleSelectionBounds, accommodating IVs as a special case while also encompassing other identifying strategies.\footnote{In related work, LevisEtAl(23)-CovariateAssistedIV analyze covariate-assisted Balke-Pearl bounds on the average treatment effect. semenova2023generalized generalizes the Lee(09) bounds on the average treatment effect on the always takers, which requires a binary instrument that monotonically affects selection, to incorporate covariate information. Monotonicity is particularly unlikely to hold when the instrument arises from the random assignment of decision makers. See, for example, FrandsenLefgrenLeslie(19), ArnoldDobbieHull(20)-discrimination_bail, ChanGentzkowYu(21).}
Our work relates to sensitivity analysis in causal inference, where researchers bound treatment effects under relaxations of unconfoundedness. Influential frameworks include the marginal sensitivity model Tan(06), KallusMaoZhou(18)-IntervalEstimation, ZhaoSmallBhattacharya(19), DornGuoKallus(21), KallusZhou(21), DornGuo(21), Rosenbaum's sensitivity model Rosenbaum(02), yadlowsky2018bounds, and partial $c$-dependence MastenPoirier(18), MastenPoirier(20). These models bound how unobservables affect selection propensities—assumptions that may be difficult to formulate when decisions involve human judgment. We instead bound outcome differences, which may be more intuitive in applied settings. Nonetheless, these models imply specific bounding functions in our framework, as we discuss in Remark (ref). chernozhukov2022long develop omitted variable bounds for linear functionals; this covers our overall performance measures but not class-specific performance.
We consider data $(X_i, D_i, Y_i)$ for $i = 1, \hdots, n$ drawn i.i.d. from a joint distribution $\mathbb{P}(\cdot)$, where $X_i \in \mathcal{X} \subseteq \mathbb{R}^{d}$ are covariates and $D_i \in \{0, 1\}$ is the decision. There is a “selective labels” problem: the observed outcome satisfies $Y_i = D_i \, Y_i^*$ for some true outcome $Y_i^* \in \{0, 1\}$ KLLLM(18). We assume $0 < \mathbb{P}(Y_i^* = 1) < 1$ so that we face a non-trivial prediction problem and that the decision satisfies strict overlap (i.e., there exists $\epsilon > 0$ such that $\mathbb{P}(D_i = 1 \mid X_i = x) \geq \epsilon$ with probability one).
\paragraph{Example: lending decisions} We observe data on loan applications, where $X_i$ contains information reported on the application such as income and prior credit history, $D_i$ is whether the application was granted a loan, and $Y_i = D_i \, Y_i^*$ is whether the applicant defaulted if granted the loan. We study credit scores that predict whether a new applicant would default $Y_i^* = 1$ BlattnerNelson(21), DiMaggioEtAl(22)-Fintech, FusterEtAl(22)predictably. $\blacktriangle$
\paragraph{Example: pretrial release} We observe data on defendants, where $X_i$ contains information about the defendant such as their current charge and prior conviction history, $D_i$ is whether the defendant was released before trial, and $Y_i = D_i \, Y_i^*$ is whether the defendant committed pretrial misconduct if released. We study pretrial risk tools that predict whether a new defendant would commit pretrial misconduct $Y_i^* = 1$ KLLLM(18), ArnoldDobbieHull(20)-discrimination_alg, AngelovaEtAl(22). $\blacktriangle$
\paragraph{Example: medical testing} We observe data on patients, where $X_i$ contains patient information such as vital signs at admission and prior medical history, $D_i$ is whether a doctor decided to conduct a diagnostic test such as a cardiac stress test or lab test, and $Y_i = D_i \, Y_i^*$ is whether the patient was diagnosed with a heart attack or bacterial infection. We study medical risk tools that predict whether a new patient suffered a heart attack or bacterial infection $Y_i^* = 1$ MullainathanObermeyer(21), HuangEtAl(22). $\blacktriangle$
Let $\pi_d(x) := \mathbb{P}(D_i = d \mid X_i = x)$ for $d \in \{0, 1\}$ be the propensity score and $\mu_{1}(x) := \mathbb{P}(Y_i^* = 1 \mid D_i = 1, X_i = x)$. Averages of a random variable $V_i$ are $\mathbb{E}_n[V_i] := n^{-1} \sum_{i=1}^{n} V_i$. We let $\| \cdot \|$ denote the appropriate $L_2$-norm by context, where $\|f \| = \left( \int f(v)^2 d\mathbb{P}(v) \right)^{1/2}$ for a measurable function $f(\cdot)$ and $\|v\| = \left( \sum_{j=1}^{k} v_j^2 \right)^{1/2}$ for a vector $v \in \mathbb{R}^k$.
A predictive algorithm $s(\cdot) \colon \mathcal{X} \rightarrow [0, 1]$ predicts the probability $Y_i^* = 1$ given covariates $X_i$, therefore estimating the conditional likelihood $\mu^*(x) := \mathbb{P}(Y_i^* = 1 \mid X_i = x)$. Given a predictive algorithm, researchers evaluate its performance in two ways.
First, given a predictive algorithm, researchers construct decision rules $d(\cdot) \colon \mathcal{X} \rightarrow \{0, 1\}$ by applying a threshold rule; that is, defining $d(X_i) = 1\{s(X_i) \geq \tau\}$ for some threshold $\tau \in [0, 1]$. Alternative decision rules are evaluated based on their expected payoff
where each $u_{d,y^*}(\cdot) \geq 0$ for $d,y^* \in \{0,1\}^2$ specifies the known payoff associated with a combination of decision $D_i$ and outcome $Y_i^*$ at features $X_i$. Payoffs are normalized so that $\sum_{d,y^* \in \{0,1\}^2} u_{d,y^*}(x) = 1$ for all $x \in \mathcal{X}$. In consumer lending, the profitability of approving customers that would not default $u_{1,0}(\cdot)$ may vary based on features like the requested loan size. In pretrial release, the societal benefits of releasing defendants that would not commit misconduct $u_{1,0}(\cdot)$ may vary based on the defendant's age, charge severity, or prior history of pretrial misconduct Baughman(17).
Second, given a predictive algorithm $s(\cdot)$, researchers evaluate its performance over the full population or specific groups. The target parameters are its overall performance
and its positive class and negative class performance
where $\beta_0(X_i; s), \beta_1(X_i; s)$ are researcher-specified functions that may depend on covariates and the predictive algorithm. Alternative choices of $\beta_0(X_i; s), \beta_1(X_i ; s)$ recover popular evaluations used in empirical research. We discuss several leading examples.
\paragraph{Example: mean square error and calibration} Mean square error is $\mathbb{E}[\left( s(X_i) - Y_i^* \right)^2] = \mbox{perf}(s; \beta)$ for $\beta_0(X_i; s) = s^2(X_i)$, $\beta_1(X_i; s) = 1 - 2s(X_i)$. Calibration at prediction bin $[r_1, r_2] \subseteq [0,1]$ is $\mathbb{E}[Y_i^* \mid r_1 \leq s(X_i) \leq r_2] = \mbox{perf}(s; \beta)$ for $\beta_0(X_i; s) := 0$, $\beta_1(X_i; s) := \frac{1\{r_1 \leq s(X_i) \leq r_2\}}{\mathbb{P}\left(r_1 \leq s(X_i) \leq r_2 \right)}$ assuming $\mathbb{P}(r_1 \leq s(X_i) \leq r_2) > 0$. $\blacktriangle$
\paragraph{Example: generalized true positive and false positive rates} Generalized true positive rate is $\mathbb{E}[s(X_i) \mid Y_i^* = 1] = \mbox{perf}_{+}(s; \beta)$ for $\beta_0(X_i; s) = s(X_i)$, and generalized false positive rate is $\mathbb{E}[s(X_i) \mid Y_i^* = 0] = \mbox{perf}_{-}(s; \beta)$ for $\beta_0(X_i; s) = s(X_i)$. $\blacktriangle$
\paragraph{Example: ROC curve} The true positive rate or recall at threshold $\tau$ is $\mathbb{P}\left( s(X_i) \geq \tau \mid Y_i^* = 1\right) = \mbox{perf}_{+}(s; \beta)$ for $\beta_{0}(X_i; s) = 1\{s(X_i) \geq \tau\}$. The false positive rate at threshold $\tau$ is $\mathbb{P}\left(s(X_i) \geq \tau \mid Y_i^* = 0\right) = \mbox{perf}_{-}(s; \beta)$. The Receiver Operating Characteristic (ROC) curve $\{( \mathbb{P}(s(X_i) \geq \tau \mid Y_i^* = 0 ), \mathbb{P}( s(X_i) \geq \tau \mid Y_i^* = 1 ) ) \colon \tau \in [0, 1] \}$ summarizes the predictive algorithm's ability to separate the positive class $Y_i^* = 1$ from the negative class $Y_i^* = 0$ as $\tau$ varies. $\blacktriangle$
\paragraph{Example: precision-recall curve} The precision at threshold $\tau$ is $\mathbb{P}\left( Y_i^* = 1 \mid s(X_i) \geq \tau \right) := \mbox{perf}(s; \beta)$ for $\beta_0(X_i; s) = 0$ and $\beta_1(X_i; s) = \frac{1\{ s(X_i) \geq \tau \}}{\mathbb{P}\left( s(X_i) \geq \tau \right)}$ assuming $\mathbb{P}(s(X_i) \geq \tau) > 0$. The precision-recall curve $\{ \left( \mathbb{P}\left( s(X_i) \geq \tau \mid Y_i^* = 1\right), \mathbb{P}\left( Y_i^* = 1 \mid s(X_i) \geq \tau \right) \right) \colon \tau \in [0, 1] \}$ summarizes the predictive algorithm's ability to accurately classify $Y_i^* = 1$ as $\tau$ varies. $\blacktriangle$
Returning to the earlier examples, BlattnerNelson(21), FusterEtAl(22)predictably, and DiMaggioEtAl(22)-Fintech compare mean square error, ROC curves, and precision-recall curves of alternative predictive algorithms for default likelihood in consumer finance. KLLLM(18) and MullainathanObermeyer(21) analyze predictive algorithms in pretrial release and medical diagnosis by reporting calibration at various prediction bins and ROC curves.\footnote{In analyzing pretrial risk assessments, lakkaraju2017selective and AngelovaEtAl(22) consider the “failure rate” or “counterfactual misconduct rate” $\mathbb{P}(Y_i^* = 1, s(X_i) \leq \tau)$ at threshold $\tau \in [0,1]$. This is recovered by the overall performance $\mbox{perf}(s; \beta)$ for $\beta_0(X_i; s) = 0, \beta_1(X_i;s) = 1\{s(X_i) \leq \tau\}$.}
Selective labels create an identification problem: since $Y_i^*$ is only observed when $D_i = 1$, the conditional likelihood $\mathbb{P}(Y_i^* = 1 \mid D_i = 0, X_i)$ is not point identified from the data alone. Consequently, neither are the conditional likelihood, expected payoffs for any decision rule, nor predictive performance measures defined in Section (ref).
A natural approach would be to assume unconfounded selection --- that is, $Y_i^* \perp \!\!\!\perp D_i \mid X_i$ Coston(20)-Evaluation, mishler2021fade. Under unconfoundedness, the conditional likelihood among unselected units equals that among selected units with the same covariates: $P(Y_i^* = 1 \mid D_i = 0, X_i) = P(Y_i^* = 1 \mid D_i = 1, X_i)$. However, this assumption is often implausible in practice. Selected and unselected units may differ systematically in ways not captured by observed covariates $X_i$. In consumer lending, applicants' decisions to accept loan terms may depend on competing offers from other lenders or their private assessment of their ability to repay --- information unavailable to the researcher. In pretrial release, judges may rely on courtroom interactions or details from case files that are not recorded in administrative data. Consequently, unobserved confounding --- unobservables correlated with both the outcome $Y_i^*$ and the decision $D_i$ --- is a central challenge.
Rather than imposing unconfoundedness, we partially identify the conditional likelihood, expected payoffs, and predictive performance measures by placing restrictions on how much outcomes can differ between unselected and selected units with the same covariates. These restrictions formalize researcher beliefs about the severity of unobserved confounding in a transparent way and nest popular empirical strategies as special cases. The key to our approach is bounding the confounding function $\delta(x) := \mathbb{P}(Y_i^* = 1 \mid D_i = 0, X_i = x) - \mathbb{P}(Y_i^* = 1 \mid D_i = 1, X_i = x)$, which summarizes the difference in outcomes between unselected and selected units conditional on covariates. When $\delta(x) = 0$ for all $x \in \mathcal{X}$, the selection decision is unconfounded. When $\delta(x) \neq 0$, unobserved confounding is present: units with the same observed covariates have systematically different outcomes depending on whether they were selected.
Assumption (ref) restricts how much outcomes can differ between unselected and selected units with the same observed covariates. In consumer lending, this bounds the difference in default rates between unfunded and funded applicants conditional on application characteristics. In pretrial release, this bounds the difference in misconduct rates between detained and released defendants conditional on charge and criminal history.
Assumption (ref) allows researchers to encode beliefs about unobserved confounding without committing to a specific selection model. It nests two important special cases: unconfoundedness with $\underline{\delta}(x; \eta) = \overline{\delta}(x; \eta) = 0$ for all $x \in \mathcal{X}$, and the assumption-free, worst-case bounds with $\underline{\delta}(x; \eta) = -\mathbb{P}(Y_i^* = 1 \mid D_i = 1, X_i = x)$, $\overline{\delta}(x; \eta) = 1 - \mathbb{P}(Y_i^* = 1 \mid D_i= 1, X_i = x)$ Manski(89), Manski(94)-SampleSelectionBounds. Thus Assumption (ref) interpolates between strong identifying assumptions and possibly uninformative bounds.
More broadly, Assumption (ref) follows the structure of widely-used partial identification frameworks in econometrics, most notably Manski(03)'s "approximate mean independence" assumption. It extends this by conditioning on covariates and allowing bounds to depend on identified nuisance parameters. In this sense, Assumption (ref) formalizes restrictions in terms that align with domain knowledge: researchers in these settings naturally reason about outcome differences rather than, say, how unobservables affect selection propensities. As shown in Section (ref), the bounding functions can be chosen through direct calibration based on domain expertise or derived indirectly from additional data features such as proxy outcomes or instrumental variables. Each strategy is instantiated by specifying both the functional form of the bounding functions and the collection of nuisance parameters $\eta$ (possibly infinite-dimensional) to be estimated.
Under Assumption (ref), the conditional likelihood and evaluation estimands are partially identified. We assume throughout that the chosen bounding functions yield valid conditional probabilities. Specifically, for all $\delta \in \Delta$ and $x \in \mathcal{X}$, $\mu_1(x) + \pi_0 \delta(x; \eta) \in [0, 1]$ and $\delta(x; \eta) \in [-1, 1]$. This is a mild restriction in practice; as we discuss below, for common bounding functions, these constraints can be enforced through the researcher's specification of sensitivity parameters or are satisfied automatically, as we will discuss next. If these constraints are violated and the researcher proceeds with our identification framework, the resulting bounds will be conservative (non-sharp) but may still be informative.
Since $\mu^*(x) = \mu_1(x) + \pi_0(x) \delta(x)$, the identified set is the set of all values consistent with the bounds on $\delta(x)$,
The identified sets for the performance measures, $\mathcal{H}(\mbox{perf}(s; \beta))$ and $\mathcal{H}(\mbox{perf}_{+}(s; \beta))$, and for the expected payoff of a decision rule $\mathcal{H}(U(d))$ are defined analogously. As shorthand, let $\beta_{0,i} := \beta_0(X_i; s), \beta_{1,i} = \beta_{1}(X_i; s)$ and $\underline{\delta}_i := \underline{\delta}(X_i; \eta)$, $\overline{\delta}_i := \overline{\delta}_i(X_i; \eta)$.
The remainder of the paper applies this partial identification framework to the two objectives for predictive algorithms commonly seen in applied work. Section (ref) addresses algorithm design: we estimate bounds on the conditional likelihood, thereby constructing predictive algorithms robust to unobserved confounding. Section (ref) addresses algorithm evaluation: we estimate bounds on performance measures for a given predictive algorithm. Before proceeding to these estimation results, we next discuss how researchers might choose the bounding functions in Assumption (ref).
How should the bounding functions $\underline{\delta}(\cdot; \eta)$ and $\overline{\delta}(\cdot; \eta)$ be chosen? The choice encodes assumptions about how much unselected units might differ from selected units with the same observed characteristics. In consumer lending: How much more likely are unfunded applicants to default than funded applicants with the same credit profile? In pretrial release: How much more likely are detained defendants to commit misconduct than released defendants with the same criminal history? We next show how alternative bounding functions formalize three common strategies for addressing selective labels in applied research.
A natural starting point is to assume that the unobserved outcome rate among unselected units is similar to the observed outcome rate among selected units with the same characteristics. For instance, in pretrial release, detained defendants with a given criminal history might be 1.5 times as likely to commit misconduct as released defendants with the same history. In consumer lending, unfunded applicants might be twice as likely to default as funded applicants with the same credit profile.
Observed outcome bounds formalize this intuition by specifying, for researcher-chosen constants $\underline{\Gamma}, \overline{\Gamma} > 0$ and all $x \in \mathcal{X}$,
Under this choice, the unobserved conditional probability among unselected units satisfies $\underline{\Gamma} \mu_1(x) \leq P(Y_i^* = 1 \mid D_i = 0, X_i = x) \leq \overline{\Gamma} \mu_1(x)$ --- that is, it is bounded by constant multiples of the observed conditional probability $\mu_1(x)$. Setting $\underline{\Gamma} = \overline{\Gamma} = 1$ recovers unconfoundedness. For observed outcome bounds, the probability constraints are satisfied for all choices of $\underline{\Gamma}, \overline{\Gamma}$ in a neighborhood of one under strict overlap on the decision and observed outcome probabilities being bounded away from zero and one.
This approach formalizes strategies commonly used in applied work. KLLLM(18) evaluate a constructed pretrial risk tool under varying assumptions on how much misconduct rates among detained defendants differ from released defendants (their Table 5). In consumer finance, researchers impute default rates for rejected applicants by applying multiplicative adjustments to observed default rates among accepted applicants HandHenley(93), ZengZhao(14). Observed outcome bounds flexibly encode such adjustments while allowing the imputed values to vary richly with covariates.
Sometimes researchers observe a proxy outcome $\widetilde{Y}_i \in \{0, 1\}$ for both selected and unselected units that is statistically related to the true outcome $Y_i^*$. This proxy can help discipline assumptions about unobserved confounding. For example, in consumer lending, whether an applicant defaulted on other credit products (e.g., credit cards) serves as a proxy for whether they would default on the loan of interest BlattnerNelson(21). In medical testing, long-term health outcomes or total medical spending serves as a proxy for whether a patient suffered an acute condition ObermeyerEtAl(19), ChanGentzkowYu(21), MullainathanObermeyer(21).
The key assumption is that the proxy outcome is equally informative about the true outcome for both selected and unselected units. Formally, suppose the proxy outcome satisfies, for all $x \in \mathcal{X}$,
Conditional on the covariates, the probability that the true and proxy outcome agree is the same regardless of selection status. This assumption, combined with the observed distribution of $\widetilde{Y}_i$ among unselected units, pins down specific bounding functions on $\delta(x)$.
Under mild restrictions satisfied by the proxy outcomes used in BlattnerNelson(21) and MullainathanObermeyer(21), the bounding functions simplify to
where $\gamma_1(x) := \mathbb{P}(Y_i^* = \tilde{Y}_i \mid D_i = 1, X_i = x)$ is the probability the true and proxy outcome agree among selected units and $\widetilde{\mu}_0(x) := \mathbb{P}(\tilde{Y}_i = 1 \mid D_i = 0, X_i = x)$ is the observed proxy outcome rate among unselected units. These bounding functions satisfy the probability constraints under the maintained proxy assumptions under strict overlap on the decision and observed outcome probabilities being bounded away from zero and one. Online Appendix (ref) provides the general characterization without these simplifying restrictions. To our knowledge, proxy outcome assumptions of this form have not been explored in the causal inference sensitivity analysis literature.
Sometimes researchers observe an instrumental variable $Z_i \in \mathcal{Z}$ with finite support that generates quasi-random variation in selection decisions but does not directly affect outcomes (i.e., $Y_i^*$ independent of $Z_i$ given $X_i$). For example, in pretrial release, judges randomly assigned to cases vary in their propensity to release defendants KLLLM(18), ArnoldDobbieHull(20)-discrimination_bail, Rambachan(21), AngelovaEtAl(22). The day-of-week affects whether certain tests are ordered MullainathanObermeyer(21), and state-level variation in banking regulations affects loan approval rates BlattnerNelson(21).
An instrument $Z_i$ implies specific bounding functions through a classic result of Manski(94)-SampleSelectionBounds. For any fixed instrument value $z \in \mathcal{Z}$, the confounding function satisfies, for all $x \in \mathcal{X}$, $$ \underline{\delta}_{z}(x; \eta) \leq \delta(x) \leq \overline{\delta}_{z}(x; \eta), $$ where $\underline{\delta}_{z}(x) = \left( \mathbb{E}[Y_i D_i \mid X_i = x, Z_i = z] - \mu_1(x) \right)/\pi_0(x)$ and $\overline{\delta}_z(x) := ( \pi_0(x, z) + \mathbb{E}[Y_i D_i \mid X_i = x, Z_i = z] - \mu_1(x) )/\pi_0(x)$. Our identification and estimation results apply to these bounds for any fixed $z \in \mathcal{Z}$.
Researchers can obtain sharper bounds by aggregating across all instrument values: $\max_{z \in \mathcal{Z}} \ \underline{\delta}_{z}(x; \eta) \leq \delta(x) \leq \min_{z \in \mathcal{Z}} \ \overline{\delta}_{z}(x; \eta)$. These intersection bounds Chernozhukov-IntersectionBounds are tighter but involve maximum and minimum operators, creating non-smoothness that complicates estimation. Following LevisEtAl(23)-CovariateAssistedIV, we address this challenge using smooth approximations to the intersection bounds. Online Appendix (ref) discusses both approaches: bounds based on a fixed instrument value and smooth approximations to intersection bounds.\footnote{LevisEtAl(23)-CovariateAssistedIV and semenova2023adaptive consider an alternative strategy for dealing with pointwise intersection bounds that invokes a “margin condition” on the separation of the maximal value of the intersection bounds. Under such a margin condition, existing results would apply to bounds on overall performance. Extending this approach to estimators for the conditional likelihood bounds or positive class performance bounds is not immediate. We leave the full analysis of margin conditions for future work.}
In this section, we propose estimators for the bounds on the conditional likelihood $\mu^*(\cdot)$ under Assumption (ref). We build on recent work studying the estimation of heterogeneous treatment effects in causal inference KunzelEtAl(19), NieWager(20), Kennedy(22)-towardsoptimal, and our estimators use cross-fitting and doubly-robust bias corrections. Since the bounds on the conditional likelihood are an infinite-dimensional parameter, we establish finite-sample bounds on our estimator's integrated mean square error (MSE) convergence using a novel oracle inequality for pseudo-outcome regression procedures. These convergence guarantees enable us to construct plug-in decision rules with provable regret bounds.
To simplify exposition, we develop our estimator assuming the researcher specifies observed outcome bounds (Section (ref)) and illustrate it with two folds. The analysis for multiple folds and cross-fitting is straightforward. We extend to proxy outcome bounds in Online Appendix (ref) and instrumental variable bounds in Online Appendix (ref).
\paragraph{Estimation Procedure for Observed Outcome Bounds:} Under Assumption (ref), the unobserved conditional probability among unselected units can be decomposed as $P(Y^*_i = 1 | D_i = 0, X_i = x) = \mu_1(x) + \delta(x)$, where $\mu_1(x)$ is the observed rate among selected units and $\delta(x)$ measures confounding. Our estimation strategy constructs pseudo-outcomes based on influence functions that combine observed outcomes among selected units with extrapolation to unselected units. The pseudo-outcomes correct for both observable selection (via inverse propensity weighting) and unobservable confounding (via the sensitivity parameters $\underline{\Gamma}, \overline{\Gamma}$).
We estimate the bounds on the conditional likelihood under observed outcome bounds with the following steps. Split the data into two folds and estimate the nuisance functions $\hat \eta = \left( \hat{\mu}_1(\cdot), \hat{\pi}_1(\cdot) \right)$ on the first fold. On the second fold, construct the pseudo-outcomes $\phi_{\mu,i}(\hat \eta) + \left( \overline{\Gamma} - 1 \right) \phi_{\pi \mu, i}(\hat \eta)$ for the upper bound and $\phi_{\mu,i}(\hat \eta) + \left( \underline{\Gamma} - 1 \right) \phi_{\pi \mu, i}(\hat \eta)$ for the lower bound, where $\phi_{\mu,i}(\eta) := \mu_1(X_i) + \frac{D_i}{\pi_1(X_i)} (Y_i - \mu_1(X_i))$ and $\phi_{\pi \mu, i}(\eta) := \left((1-D_i) - \pi_0(X_i)\right) \mu_1(X_i) + \frac{D_i}{\pi_1(X_i)} (Y_i - \mu_1(X_i)) \pi_0(X_i) + \pi_0(X_i) \mu_1(X_i)$ are the influence functions for $\mathbb{E}\left[ \mu_1(X_i)\right]$ and $\mathbb{E}\left[ \pi_0(X_i) \mu_1(X_i) \right]$ respectively. These can be derived using standard arguments Hines_2022. We regress the estimated pseudo-outcomes on the covariates $X_i$ using a researcher-specified nonparametric regression procedure in the second fold. This yields our estimators $\widehat{\underline{\mu}}(\cdot)$, $\widehat{\overline{\mu}}(\cdot)$.
We analyze the integrated MSE convergence rate of our proposed estimators by comparing them against an infeasible oracle nonparametric regression. This oracle has access to the true nuisance functions $\eta$ and estimates the conditional likelihood bounds by regressing the true pseudo-outcomes $\phi_{\mu,i}(\eta) + \left( \overline{\Gamma} - 1 \right) \phi_{\pi \mu, i}(\eta)$ for the upper bound and $\phi_{\mu,i}(\eta) + \left( \underline{\Gamma} - 1 \right) \phi_{\pi \mu, i}(\eta)$ for the lower bound on $X_i$ in the second fold using the same nonparametric regression procedure specified by the researcher. We denote these oracle estimators by $\widehat{\underline{\mu}}_{oracle}(\cdot), \widehat{\overline{\mu}}_{oracle}(\cdot)$.
This oracle comparison is useful because it decomposes the integrated MSE of our proposed estimators into two sources. The first is approximation error from using a finite-sample nonparametric regression procedure to estimate the infinite-dimensional conditional likelihood bounds --- this is captured by the oracle's integrated MSE. The second is estimation error from first-stage nuisance function estimation. Our next result shows that the latter only affects our proposed estimator through a smoothed, doubly-robust remainder term.
Proposition (ref) delivers two insights. First, our analysis is modular: researchers can pair choices of nuisance function estimators in the first step with any choice of nonparametric regression method (satisfying Assumption (ref)) in the second step, and the integrated MSE is characterized by the oracle rate plus a doubly-robust remainder term. Second, Proposition (ref) shows that the cost of estimating the nuisance functions in the first step only enters through the product of nuisance function errors.
To prove this result, we establish an oracle inequality on the $L_2(\mathbb{P})$-error of pseudo-outcome regression procedures (Lemma (ref)), extending Kennedy(22)-towardsoptimal's pointwise analysis. This result may be of independent interest. The $L_2(\mathbb{P})$-stability condition (Assumption (ref)) on $\widehat{\mathbb{E}}_n[\cdot \mid X_i = x]$ is quite mild in practice. It is satisfied by a variety of generic linear smoothers such as linear regression, series regression, nearest neighbor matching, random forest models, and several others (Proposition (ref)). Once again, in this sense, Proposition (ref) is agnostic to both the researcher's choice of ${\widehat{\mathbb{E}}_n[\cdot \mid X_i = x]}$ and nuisance function estimators.
Furthermore, Proposition (ref) can be applied in settings where nuisance functions satisfy additional smoothness or sparsity conditions, and for particular choices of the second-stage regression estimator and nuisance function estimators. Known results on mean-squared error convergence rates of nonparametric regression procedures can then be used to characterize $\|\widehat{\overline{\mu}}_{oracle}(\cdot) - \overline{\mu}^*(\cdot) \|$ and $\| \tilde{R}(\cdot) \|$, yielding specific rates of convergence in terms of primitives such as underlying smoothness assumptions and sample size. This follows in the spirit of recent work on optimal estimation of heterogeneous treatment effects Kennedy(22)-towardsoptimal, KennedyEtAl(23)-Minimax.
Researchers often apply a threshold rule to an estimated predictive algorithm to make decisions. For example, KLLLM(18), Rambachan(21), AngelovaEtAl(22) compare outcomes under alternative decisions rules that threshold a predictive algorithm in pretrial release. We show that plug-in decision rules based on our estimators of the conditional likelihood bounds achieve low regret, ensuring that the statistical uncertainty from estimation does not lead to substantially worse decisions than if the true conditional probabilities were known.
Since the expected payoff $U(d)$ of any decision rule $d(\cdot) \colon \mathcal{X} \rightarrow \{0, 1\}$ is partially identified (Lemma (ref)), we evaluate decision rules by comparing their worst-case expected welfare Manski(07). An immediate consequence of Lemma (ref) is that the optimal max-min decision rule $d^*(\cdot) \in \arg \max_{d(\cdot) \colon \mathcal{X} \rightarrow [0,1]} \underline{U}(d)$ is a threshold rule: $$ d^*(X_i) = 1\{ \widetilde{\mu}^*(X_i) \leq u_{1,0,i} + u_{0,0,i}\}, $$ where $\widetilde{\mu}^*(x) = (u_{1,1,i} + u_{1,0,i}) \overline{\mu}^*(x) + (u_{0,0,i} + u_{0,1,i}) \underline{\mu}^*(x)$ is a welfare-weighted average of the conditional probability bounds. We construct a feasible plug-in decision rule $\widehat{d}(X_i)$ by substituting in our nonparametric estimators $\widehat{\overline{\mu}}(\cdot)$, $\widehat{\underline{\mu}}(\cdot)$.
The regret of the plug-in decision rule is $R(\hat{d}) = \underline{U}(d^*) - \underline{U}(\widehat{d})$, which measures the welfare loss from using estimated rather than true conditional probabilities. Our next result bounds squared regret by the integrated MSE of the oracle estimators plus the smoothed, doubly-robust remainder term.
The proof proceeds in two steps. First, we show that the regret of the plug-in decision rule is bounded by the integrated MSE of our estimators for the conditional probability bounds --- this follows because the optimal decision rule is a threshold rule, and threshold rules are sensitive to estimation errors measured in $L_2(\mathbb{P})$ norm. Second, we apply Proposition (ref) to decompose this integrated MSE into the oracle rate plus the doubly-robust remainder term. Consequently, worst-case regret converges to zero whenever the oracle's integrated MSE converges to zero. Plug-in decision rules based on our nonparametric estimators are therefore robust: they achieve nearly optimal welfare despite estimation uncertainty.
In this section, we construct estimators for the bounds on the performance of a given predictive algorithm $s(\cdot)$ under Assumption (ref). This setting arises empirically in two scenarios: (i) when researchers evaluate an externally provided algorithm (e.g., a proprietary credit score or a commercial risk assessment tool), and (ii) when researchers conduct hold-out evaluation of an algorithm trained on separate data. In both cases, the algorithm $s(\cdot)$ is taken as fixed, and the goal is to assess its performance under varying assumptions about unobserved confounding. Our estimators are consistent for the bounds on the identified set as sample size grows large, and we characterize their convergence rates in terms of errors in the first-step estimation of nuisance parameters.
We first construct estimators for the bounds on overall performance $\mbox{perf}(s; \beta)$. To simplify exposition, we develop estimators for observed outcome bounds (Section (ref)). We extend to proxy outcome bounds and instrumental variable bounds in Online Appendices (ref)-(ref) respectively. Under observed outcome bounds, Lemma (ref) implies that the upper bound can be written as $$ \overline{\mbox{perf}}(s; \beta) := \mathbb{E}\left[ \beta_{0,i} + \beta_{1,i} \mu_1(X_i) + \beta_{1,i} \left( 1\{ \beta_{1,i} > 0 \} \left( \overline{\Gamma} - 1 \right) + 1\{ \beta_{1,i} \leq 0 \} (\underline{\Gamma} - 1) \right) \pi_0(X_i) \mu_1(X_i) \right]. $$ The lower bound $\underline{\mbox{perf}}(s; \beta)$ has an analogous expression. Both are linear functionals of known functions of the data and identified nuisance parameters $\eta = (\mu_1, \pi_0)$. This linearity enables us to construct debiased estimators using standard arguments based on influence functions and cross-fitting ChernozhukovEtAl(18)-DML, ChernozhukovEtAl(22)-LocallRobust, Kennedy(22)-IFreview.
\paragraph{Estimation Procedure for Observed Outcome Bounds:} For completeness, we sketch the construction of our estimators based on $K$-fold cross-fitting. We randomly split the data into $K$ disjoint folds. For each fold $k$, we estimate the nuisance functions $\widehat{\eta}_{-k}$ using all observations not in the $k$-th fold and construct $$ \overline{\mbox{perf}}_{i}(\hat \eta_{-k}) := \beta_{0,i} + \beta_{1,i} \phi_{\mu,i}(\hat \eta_{-k}) + \beta_{1,i} \left( 1\{ \beta_{1,i} > 0 \} \left( \overline{\Gamma} - 1 \right) + 1\{ \beta_{1,i} \leq 0 \} (\underline{\Gamma} - 1) \right) \phi_{\pi \mu, i}(\hat \eta_{-k}) $$ for each observation $i$ in the $k$-th fold, where $\phi_{\mu,i}(\eta)$ and $\phi_{\pi \mu, i}(\eta)$ are defined as before. We then average over all observations $\widehat{\overline{\mbox{perf}}}(s; \beta) := \mathbb{E}_n\left[ \overline{\mbox{perf}}_i(\hat \eta_{-K_i}) \right],$ where $K_i$ is the fold containing observation $i$. Our estimator for the lower bound $\widehat{\underline{\mbox{perf}}}(s;\beta)$ is defined analogously.
Our next result characterizes the convergence rate and asymptotic distribution of our estimators for the bounds on overall performance under observed outcome bounds.
The estimators converge at rate $O_{\mathbb{P}}(1/\sqrt{n})$ plus a term that captures the nuisance estimation error. Nuisance errors enter multiplicatively --- through the product of propensity score and outcome regression errors --- rather than additively. The rate condition required for the estimators to be jointly asymptotically normal is satisfied under the standard condition that all nuisance function estimators converge at a rate faster than $O_\mathbb{P}(n^{-1/4})$. This is the familiar condition from debiased machine learning, satisfied by a wide range of nonparametric regression and modern machine learning methods.
Online Appendix (ref) constructs a consistent estimator of the asymptotic covariance matrix in Proposition (ref). This implies that confidence intervals for the identified set $\mathcal{H}(\mbox{perf}(s; \beta))$ can be constructed using standard methods ImbensManski(04), Stoye(09)--PICI.
We next estimate bounds on positive class performance $\mbox{perf}_{+}(s; \beta)$. Unlike overall performance, this is a ratio of expectations, and the bounds are given by the optimal values of linear fractional programs as shown in Lemma (ref). This creates additional challenges in estimation. In this section, we will develop estimators for these bounds and establish their convergence properties.
We describe the construction of our estimators using sample-splitting across two folds; extending to cross-fitting with multiple folds is straightforward. We randomly split the data into two folds (we assume $n$ is divisible by $2$ and each fold contains $n/2$ observations for simplicity). In the first fold, we estimate the nuisance functions $\hat \eta := (\widehat{\mu}_1(\cdot), \widehat{\pi}_1(\cdot))$ and the bounding functions $(\widehat{\underline{\delta}}(\cdot), \widehat{\overline{\delta}}(\cdot))$ using one of the strategies in Section (ref). On the second fold, we construct $\phi_{\mu,i}(\hat \eta)$ and solve the sample analogue
where $\widehat{\Delta}_n := \left\{ \tilde{\delta} \colon (1 - D_i) \widehat{\underline{\delta}}(X_i) \leq \tilde{\delta}_i \leq (1 - D_i) \widehat{\overline{\delta}}(X_i) \right\}$ and $\mathbb{E}_{n}[\cdot]$ denotes the sample average over the second fold. The optimization problem in Equation (ref) can be solved efficiently using standard techniques. It can be written as a linear program by applying the Charnes-Cooper transformation CharnesCooper(62). See Online Appendix (ref) for details. The construction of this estimator involves two key ideas. First, we replace $\mu_1(\cdot)$ in the objective with its estimated influence function $\phi_{\mu,i}(\hat \eta)$. This orthogonalizes against first-stage estimation errors in nuisance functions that only enter the objective. Second, we plug in estimated bounding functions into the constraint set.
We provide an error decomposition for the positive class performance bounds. This result characterizes how errors in estimating the nuisance parameters $(\mu_1(\cdot), \pi_1(\cdot))$ and the bounding functions $(\underline{\delta}(\cdot), \overline{\delta}(\cdot))$ propagate to the final estimator.
Proposition (ref) reveals the structure of how estimation errors affect our estimator. The first term $O_{\mathbb{P}}(1/\sqrt{n})$ reflects sampling error that would arise even if the nuisance functions $(\mu_1(\cdot), \pi_1(\cdot))$ and the bounding functions $(\underline{\delta}(\cdot), \overline{\delta}(\cdot))$ were known. The second term $O_{\mathbb{P}}(R_{1,n})$ reflects errors from estimating the nuisances $(\mu_1(\cdot), \pi_1(\cdot))$ that appear in the numerator and denominator of the linear fractional program. These enter through a doubly robust product, and this behavior is inherited from the influence function structure underlying the numerator and denominator. The final term captures errors from estimating the bounding functions $(\underline{\delta}(\cdot), \overline{\delta}(\cdot))$ that define the constraint set $\Delta$. This error enters through the root mean square error. This error can be controlled using standard nonparametric regression techniques, and researchers can trade off bias and variance in familiar ways.
For this reason, Proposition (ref) establishes a “partial” double robustness: given estimated constraints, our estimator is doubly robust to errors in the objective nuisances $(\mu_1(\cdot), \pi_1(\cdot))$. Constraint estimation errors enter additively through RMSE. Achieving full double robustness --- where constraint estimation errors would also enter with product structure --- remains an open problem, which we discuss in Remark (ref). Remark (ref) describes a simple alternative approach that obtains fully doubly-robust estimators for non-sharp bounds on positive class performance.
Proposition (ref) and our estimator are practically valuable. Researchers can use state-of-the-art machine learning methods for estimating the objective nuisances $(\mu_1(\cdot), \pi_1(\cdot))$, while separately controlling the constraint estimation error through their choice of estimator for the bounding functions. The error decomposition applies to any estimator of the bounding functions with favorable $L_2(\mathbb{P})$ properties, and we discuss alternative strategies and their tradeoffs in Section (ref). This modularity is a key practical advantage of our framework.
The error decomposition in Proposition (ref) applies to any estimator of the bounding functions with favorable mean square error properties. This modularity gives researchers flexibility. We discuss two natural strategies and their relative merits, focusing on observed outcome bounds for concreteness.
The simplest approach directly substitutes nuisance estimates into the bounding function formulas. For observed outcome bounds, this means $\widehat{\underline{\delta}}(x) = (\underline{\Gamma} - 1) \widehat{\mu}_1(x)$ and $\overline{\delta}(x) = (\overline{\Gamma} - 1) \widehat{\mu}_1(x)$. The mean square error properties of the estimated bounding functions are then inherited from those of $\widehat{\mu}_1(x)$. This approach requires no additional estimation beyond the first-stage nuisances. An alternative approach is to construct influence-function-based pseudo-outcomes and regress them on covariates using the full sample (i.e., our pseudo-outcome regression procedure in Section (ref)). We develop this approach below as it may offer bias-variance advantages in practice. Ultimately, both strategies are valid under our framework.
To understand why pseudo-outcome regression may be advantageous, consider that $\widehat{\mu}_1(x)$ is typically constructed using only observations with $D_i = 1$. Consequently, standard regression procedures target the mean square error conditional on $D_i = 1$ However, the error bound in Proposition (ref) depends on the unconditional mean square error. These two objectives could lead to different bias-variance tradeoffs when the selected population ($D_i=1$) is small relative to the full population or there is covariate shift between the selected and unselected populations. In such settings, an estimator that directly targets the unconditional MSE may achieve better performance. Pseudo-outcome regression accomplishes this by constructing unbiased pseudo-outcomes using observations from both the selected and unselected populations, then regressing these pseudo-outcomes on the full sample.
\paragraph{Estimation Procedure for Observed Outcome Bounds:} To illustrate how Proposition (ref) may be applied, we estimate bounding functions using our pseudo-outcome regression procedure. This illustrates how our partial double robustness result (Proposition (ref)) can be combined with analyses of nonparametric regression estimators --- in particular, leveraging our oracle inequality for pseudo-outcome regression (Appendix (ref)).
We now randomly split the data into three folds (we assume $n$ is divisible by $3$ and each fold contains $n/3$ observations for simplicity). We construct nuisance function estimates $\hat \eta$ using the first fold. On the second fold, we construct the pseudo-outcomes $\phi_{\mu,i}(\hat \eta)$ and regress the estimated pseudo-outcomes on the covariates $X_i$ using a researcher-specified nonparametric regression procedure, yielding $\widehat{\delta}(\cdot)$. On the third fold, we solve the sample analogue in Equation (ref). Our analysis in Section (ref) immediately allows us to analyze the errors from the estimated bounding functions.
We validate the finite-sample performance of our estimators through Monte Carlo simulations in Appendix (ref). We now apply our framework to credit risk prediction, where financial institutions use scores $s(\cdot)$ to predict default risk. We observe data on past loan applications but only observe defaults $Y_i = D_i Y_i^*$ for funded applications ($D_i = 1$). It is implausible to assume unconfounded funding decisions. Applicants who accept loan offers may differ systematically from those who decline in ways that affect their default risk.\footnote{Lending institutions face two distinct inference problems. First, "policy rejects" with $\mathbb{P}(D_i = 1 \mid X_i) = 0$ (i.e., overlap violations) due to underwriting rules require extrapolation-based methods that are beyond our scope. Second, even among fundable applicants with $\mathbb{P}(D_i = 1 \mid X_i) > 0$, those accepting versus declining offers may differ in unobserved ways --- the selective labels problem we address. Extending our framework to handle policy rejects is an important direction for future research.}
We analyze 372,346 personal loan applications submitted to Commonwealth Bank of Australia from July 2017 to July 2019. These loans are repaid in monthly installments and used to purchase vehicles, refinance debt, or cover major expenses CostonEtAl(21)-Rashomon. Loan amounts ranged up to AU\$50,000 (median AU\$10,000) with a median offered interest rate of 14.9% per annum. We observe rich application-level covariates including reported income, occupation, and credit history at CommBank. We only observe whether an applicant defaulted within 5 months $Y_i^* \in \{0, 1\}$ if the application was funded $D_i = 1$. Approximately one-third of applications were funded, with a 2.0% default rate among funded loans.
We consider two exercises. First, we construct credit risk scores robust to varying confounding assumptions. The resulting algorithms produce substantially different risk rankings that could affect credit access. Second, following recent work in consumer finance BlattnerNelson(21), FusterEtAl(22)predictably, DiMaggioEtAl(22)-Fintech, we evaluate how a benchmark score's accuracy varies across income groups. While an analysis ignoring confounding suggests better performance for higher-income applicants, this disparity disappears under plausible confounding levels.
We construct robust credit risk scores using our estimator from Section (ref), exploring how varying confounding assumptions affect individual default risk predictions. We estimate the upper bound on the conditional default probability with observed outcome bounds, $\overline{\mu}(\cdot)$. We split applications into two folds: the first estimates $\pi_1(\cdot)$, $\mu_1(\cdot)$ using random forests; the second regresses estimated pseudo-outcomes on application-level covariates $X_i$ using cross-validated Lasso regression. We set $\underline{\Gamma} = 1$ and vary $\overline{\Gamma} \in \{1, 2, 3\}$.
We compare our robust estimator (with $\underline{\Gamma} = 1, \overline{\Gamma} = 2$) against a benchmark score trained only on funded applications. Figure (ref)'s left panel shows the joint distribution of their predictions. For each benchmark decile, the heatmap displays the percentage of applications falling into each decile of our robust predictions. Applications exhibit substantial reranking: among those in the benchmark's 5th decile, 18.7% shift to extreme deciles (1-3 or 8-10) under the robust score. These differences meaningfully affect credit access. For instance, if only the least risky third of applications receive funding, 10.1% would have their decision reversed when comparing benchmark versus robust scores.
Figure (ref)'s right panel examines how estimated coefficients vary with $\overline{\Gamma}$ (variable definitions in Online Appendix Table (ref)). Confounding assumptions affect different covariates heterogeneously. Some coefficients—total net income and credit bureau score—remain zero across all $\overline{\Gamma}$ values. Others, like the number of recent credit card applications, receive zero weight in the benchmark but non-zero weight for all $\overline{\Gamma}$>1. Still others, such as occupation type and maximum recent delinquency, exhibit large magnitude changes as $\overline{\Gamma}$ varies, highlighting how confounding assumptions shape which applicant characteristics the algorithm prioritizes.
We now evaluate an existing credit risk score's performance. We split applications into training and evaluation data. On the training data, we construct a benchmark score predicting 5-month default risk among funded applications using random forests. On the evaluation data, we analyze mean square error and ROC curves, following recent mortgage lending studies BlattnerNelson(21), FusterEtAl(22)predictably. Using observed outcome bounds with $\underline{\Gamma} = 1$, we assess how the benchmark's estimated performance varies as we increase $\overline{\Gamma}$, allowing unfunded applicants to have progressively higher default rates than funded applicants with similar characteristics.
Figure (ref)(a) shows how estimated MSE bounds vary with $\overline{\Gamma}$. When $\overline{\Gamma} = 1$ (the unconfounded case), both bounds equal the benchmark's MSE on funded applications (dashed red line; $0.143$). As $\overline{\Gamma}$ increases, allowing unfunded applicants higher default rates, bounds widen but remain informative. At $\overline{\Gamma} = 2$ --- permitting unfunded applicants to default at twice the rate of funded applicants with similar characteristics --- the estimated bounds are $[0.138, 0.182]$. Figure (ref)(b) examines how estimated ROC curves vary with $\overline{\Gamma}$. We summarize ROC curves using Area Under the Curve (AUC), computed via the trapezoidal rule. Under unconfoundedness ($\overline{\Gamma} = 1$), the benchmark's AUC is $0.637$. At $\overline{\Gamma} = 1.5$, AUC bounds are $[0.575, 0.705]$. As with MSE, bounds widen with $\overline{\Gamma}$ but remain informative at empirically plausible confounding levels.
We finally investigate the benchmark credit risk score's predictive disparities across income groups. We define $G_i \in \{0, 1\}$ to be whether an application is below or above the median personal income of all submitted applications. We examine how the benchmark credit risk score's mean square error and ROC curve vary across these income groups, thereby studying whether it is a “noisier” signal of 5-month default risk for lower income applications, as in BlattnerNelson(21), FusterEtAl(22)predictably.
Figure (ref)(a) compares the mean square error of the benchmark credit risk across income groups. Interestingly, for $\overline{\Gamma} = 1$, the mean square error of the benchmark credit risk score is significantly larger on applications below the median income (0.169) than those above it (0.136). This difference persists as we vary assumptions on unobserved confounding. The upper bound on the mean square error for applications above the median income is only larger than the lower bound on the mean square error for applications below the median income for values $\overline{\Gamma} \geq 1.69$. Ruling out this mean square error disparity across income groups would require that unfunded applications are at least 1.69 times as likely to default within 5 months as funded applications conditional on covariates.\footnote{Rather than examining disparity bounds across a range of $\overline{\Gamma}$ values, researchers could compute the “identification breakdown point” $\overline{\Gamma}^*$: the minimal value of $\overline{\Gamma}$ (holding $\underline{\Gamma}=1$) such that $0 \in H(\text{disp}(s;\beta))$. A natural estimator is $\widehat{\overline{\Gamma}}^*$, the minimal $\overline{\Gamma}$ such that zero is contained in the estimated $(1-\alpha)$ confidence set $\widehat{H}(\text{disp}(s;\beta))$. Since our confidence sets have valid coverage by the results in Section (ref) and Appendix (ref), $P(\widehat{\overline{\Gamma}}^* \geq \overline{\Gamma}^*) \geq 1-\alpha$.} Figure (ref)(b) also compares the ROC curves of the benchmark credit risk score across income groups. Intriguingly, there exists effectively no disparity across income groups based on the ROC curves. When $\overline{\Gamma} = 1$, the implied AUC for applications above the median income is $0.603$, whereas the implied AUC is $0.616$ for applications below the median income. Any disparity in terms of ROC curves is already ruled out at $\overline{\Gamma} = 1.25$ as the bounds on the AUC for applications above the median income are $[0.567, 0.653]$ and $[0.581, 0.654]$ for applications below the median income.
We developed a framework for robust design and evaluation of predictive algorithms under unobserved confounding. Our approach nests popular empirical strategies --- observed outcome bounds, proxy outcomes, and instrumental variables --- by bounding outcome differences between selected and unselected units (Assumption (ref)). We provide estimators that work well across bounding strategies and performance measures: conditional likelihoods for algorithm design, overall metrics like MSE and calibration, and class-specific measures like TPR and FPR. This breadth makes our framework valuable across domains where predictive algorithms face selective labels—pretrial release, consumer lending, medical diagnosis, hiring, and child welfare.
Several open questions present promising research directions. First, achieving full double robustness for positive class performance bounds remains open, though our technical results suggest possible paths forward. Second, uniform inference for breakdown frontiers would let researchers formally test how strong confounding must be to overturn conclusions MastenPoirier(20). Third, we evaluate fixed algorithms; extending our results to settings where the same data is used to train and evaluate an algorithm is an important direction. Finally, extending to multi-valued or continuous outcomes would broaden applicability, though defining class-specific performance in these settings requires care.
As algorithms increasingly drive high-stakes decisions in credit, criminal justice, and healthcare, robust evaluation under realistic confounding becomes essential. Our framework and these open questions offer pathways to strengthen algorithmic accountability where the stakes are highest.
{\singlespacing }
\subfile{Supplement.tex}