Inference for an Algorithmic Fairness-Accuracy Frontier

Inference for an Algorithmic Fairness-Accuracy Frontier

Introduction

Algorithms are increasingly used in many aspects of life, often to guide or support high-stake decisions, for example by predicting job performance, re-offense risk, loan default, college success, or patient health. These predictions feed, respectively, into the determination of who should be hired; which defendants should receive bail; who should be granted a loan; which students should be admitted to college; and which patients to treat. Yet, a growing body of literature documents that algorithms may exhibit bias against legally protected groups, both in their predictive accuracy and in the decisions they lead to ang:lar:mat:kir16, arn:dob:hul21,obe:pow:vog:mul19,ber:jei:jab:kea:rot21. The bias may arise, for example, due to the choice of labels the algorithm is trained on, the objective function that the algorithm optimizes, the training procedure, and many other factors involved in the design of the algorithm cow:tuc20.

Designing an algorithm often entails a trade-off between making it more fair, i.e., less likely to disproportionately harm a protected class, and more accurate, e.g., better at assigning treatment to those who benefit from it and withholding it from those who do not. As a result, improving fairness often comes at the cost of accuracy. Regulators, policymakers, algorithm designers, and actors affected by algorithmic predictions all have an interest in assessing various aspects of this trade-off.

We provide a set of tools for estimation of and statistical inference on a fairness-accuracy (FA) frontier recently characterized by \citet*[LLMO henceforth]{lia:lu:mu:oku24}, where fairness is measured by the gap between group-specific expected losses. The theoretical analysis in \citetalias{lia:lu:mu:oku24} assumes perfect knowledge of the population distribution of the observable variables and formalizes the trade-off between accuracy and fairness, shading light on how to use properties of the data distribution to determine whether it is optimal for the designer of the algorithm to exclude certain inputs from use. However, in practice, regulators and policymakers typically have access to only finite data. Hence, statistical inference tools are crucial for analyzing properties of algorithms and for their regulation.

We put forward a consistent estimator for \citetalias{lia:lu:mu:oku24}'s FA-frontier and derive its asymptotic distribution. For each point on the FA-frontier, we characterize an algorithm that achieves it. We then develop a method to test hypotheses such as: Is it optimal to fully exclude group identity from use in an algorithm? Does a particular algorithm lead to group-specific expected losses that are on the FA-frontier? How far from the fairest point on the FA-frontier are the group-specific expected losses associated with a given algorithm?

Answers to the first two questions inform the regulation of algorithms and the determination of whether discrimination occurred. The law recognizes two main categories of discrimination: disparate treatment, where individuals are deliberately treated differently based on their membership in a protected class; and disparate impact, where protected classes are adversely affected disproportionately, no matter the intent kle:lud:mul:sun18,bla:spi22. Often, as part of an effort to avoid disparate treatment, algorithms are designed so that they do not take race, gender, or other sensitive attributes as input. Even class-blind algorithms, however, may lead to disparate impact. Our first test informs a fairness-minded policymaker interested in assessing whether banning group identity has the potential to mitigate disparate impact.\footnote{This question is of interest, e.g., when assessing the recent U.S. Supreme Court decision to rule out the use of race in college admissions; for an overview, see \href{https://highered.collegeboard.org/recruitment-admissions/policies-research/access-diversity/2023-scotus-decision/overview}{College Board}.}

Our second test evaluates whether a given algorithm lies on the frontier—and thus whether a less discriminatory alternative (LDA) exists. This test is relevant to both plaintiffs (e.g., job applicants) and defendants (e.g., hiring companies) in disparate impact disputes. For example, if a selection process yields disparate impact, the hiring company may invoke business necessity to justify it. The challenger must then show the existence of an LDA, i.e., a fairer algorithm that is just as accurate. If our test rejects the null that the current algorithm is on the frontier, it supports the plaintiff's claim. Conversely, if the test fails to reject, there is no statistical evidence that the hiring company can build a fairer algorithm without sacrificing accuracy, supporting the business necessity defense. When a given algorithm is not on the frontier, we characterize alternative algorithms that improve upon it in terms of accuracy or fairness (or both).

The third inferential method yields an estimator of the distance from a given algorithm to the fairest point on the frontier and constructs a confidence interval around it. This tool may interest any fairness-minded agent (e.g., a college) willing to trade some accuracy for reduced disparate impact (e.g., via affirmative action), as it provides a measure of the trade-off between promoting equity and achieving accuracy.

The key insight underlying our proposed inference methods is that since the feasible set of group-specific expected losses associated with all possible algorithms is convex, it can be fully represented by its support function. As the FA-frontier is a portion of the boundary of the feasible set, we characterize and estimate it through this support function. We express the hypotheses listed above as restrictions on the support function, yielding easy to understand test statistics that essentially rely on judicious use of the separating hyperplane theorem. Throughout our analysis, the support function serves as a unifying tool for inference on properties of algorithms and of the FA-frontier.

We provide a consistent debiased machine learning (DML) estimator of the support function and establish that it converges to a tight Gaussian process as sample size increases, building on and extending results in ber:mol08, cha:che:mol:sch18 and sem23. We show how to allow for infimum-type test statistics that are directionally-differentiable mappings of the support function, building on results of fan:san19. Earlier uses of the support function for inference in partially identified models ber:mol08, bon:mag:mau12, kai:san14, kai16, cha:che:mol:sch18,mol20,sem23 did not include tests for hypotheses such as the ones we consider. Expressing these hypotheses in terms of restrictions on the support function is one of our main contributions.

We evaluate the finite-sample performance of our inference toolkit using extensive Monte Carlo simulations. We then demonstrate its empirical value by reanalyzing algorithms for high-risk care assignment at a research hospital studied by \citet*{obe:pow:vog:mul19}. We fail to reject the hypothesis that the hospital’s status-quo algorithm admits an LDA, and document fairness and accuracy gains from several alternative algorithms on the frontier that we characterize.

Related Literature. A growing literature in computer science and statistics studies algorithmic fairness; see cho:rot18, bar:har:nar23, and Corbett24 for comprehensive overviews and open questions. Models have been developed to explain algorithmic bias by decomposing disparity sources ram:kle:lud:mul20 or incorporating taste-based discrimination and unobservables in label generation ram:rot20. Fairness has been modeled as a constraint or regularizer in the objective function that maximizes predictive accuracy dwork12,ber:hei:jab:jos:kea:mor:jam:nee:roth17 and incorporated in the preferences of a social planner that uses algorithms in their decision-making process kle:lud:mul:ram18,ram:kle:mul:lud20. In optimal policy targeting, fairness has been set as the criterion to be maximized when choosing a policy from the set of welfare-maximizing rules viv:bra23. When protected class membership is not observed in the data but proxy variables are available, data combination methods have been proposed to partially identify disparity measures kal:mao:zhou22. Yet, tests of hypotheses for properties of the trade-off between fairness and accuracy of algorithms are scant in the literature. aue:lia:tab:oku24 propose a test, which can incorporate exogenous constraints on the algorithm space, using sample splitting for the union null hypothesis that a status quo algorithm can be weakly improved in terms of both fairness and accuracy. Their test is based on finding another algorithm within a user-specified subclass, subject to the constraint that the alternative algorithm is at least as accurate as the status quo algorithm.

In contrast, our test for existence of an LDA is a one-shot test, valid across all algorithms rather than a specific subclass, that reveals if the status quo algorithm is on the FA-frontier and does not require first estimating an alternative algorithm. Characterizing the entire FA-frontier allows us to provide a comprehensive toolkit for statistical inference that can be useful for regulators to determine what algorithm design restrictions and reporting requirements to impose on entities making decisions using algorithms.

Outline. Section (ref) lays out notation and summarizes the derivation of the FA-frontier in \citetalias{lia:lu:mu:oku24}. Section (ref) characterizes the support function of interest and uses it to describe the FA-frontier. Section (ref) derives the asymptotic properties of our DML estimator for the support function. Section (ref) uses these results to obtain a consistent estimator and an asymptotically valid confidence set for the FA-frontier, and characterizes algorithms attaining points on it. Section (ref) formulates hypotheses of interest in the fairness literature as restrictions on the support function and proposes asymptotically valid tests. Section (ref) provides an estimator and inference method for the distance between the expected group-specific losses associated with a given algorithm and the fairest point on the frontier. Section (ref) presents our Monte Carlo simulations and re-evaluation of obe:pow:vog:mul19's study on a research hospital’s use of algorithms for assigning patients to a high-risk care management program. Section (ref) concludes. Our main proofs are in Appendix (ref); Appendix (ref) reports auxiliary results and extensions, and Appendix (ref) includes supplemental empirical results.

Setup

Let a population of individuals be described by an outcome $Y\in\mathcal{Y}\subset\mathbb{R}$, a binary group identity $G\in\{r,b\}$ (red or blue), and a vector of covariates $X\in\mathcal{X}\subset\mathbb{R}^{d_X}$, with the population distribution of $(Y,G,X)$ denoted $\mathbb{P}$. For example, $Y$ may denote an individual's number of active chronic illnesses in the subsequent year, $G$ may denote their race, and $X$ may include age, gender, biomarkers, comorbidity, costs and medication variables. The relation between $G$ and $X$ is left unspecified{, but $G$ is not part of $X$}. Throughout, we assume $G$ is binary, though the results extend to multiple groups. Each individual receives a binary decision $D\in\{0,1\}$, e.g., whether they are automatically enrolled in a high-risk care management program. An algorithm $a:\mathcal{X} \mapsto [0,1]$ assigns a probability distribution to $D$; e.g., the algorithm assigns each patient a health risk score in $[0,1]$, which for simplicity we take to be the only input to the {enrollment decision, and hence to coincide with the enrollment probability}. Let $\mathcal A(\mathcal{X})$ denote the set of all algorithms that map from the input space $\mathcal{X}$ to a probability distribution over $D$, and $\ell:\{0,1\}\times\mathcal{Y}\mapsto\mathbb R$ be a function that measures the loss associated with decision $d\in\{0,1\}$ for an individual with outcome $y\in\mathcal{Y}$. In the example discussed so far, the algorithm designer observes training data consisting of covariates $X$ and a binary outcome $Y$ indicating whether someone has a number of chronic illnesses exceeding a given threshold; the loss function $\ell$ may be the classification loss, $\ell(D,Y)=\mathds{1}\{D\neq Y\}$, which returns the value $1$ if the algorithm mistakenly enrolls a healthy person in the high-risk care program or fails to enroll someone who is very sick. We assume throughout that the training data is drawn from the same distribution as the population that we eventually apply the algorithm to (i.e., the subpopulation for which labels are observed is representative of the entire population).

Given an algorithm $a\in\mathcal A(\mathcal{X})$, let the population expected loss for group $g\in\{r,b\}$ be

where the expectation is taken with respect to $\mathbb{P}$. We refer to the group-specific expected loss in Eq. (ref) as group risk. Following \citetalias{lia:lu:mu:oku24}, we define a preference ordering over group risk pairs so that $e=(e_r,e_b)$ is preferred to $e'=(e_r',e_b')$, denoted $e>_{FA}e'$, if

with at least one strict inequality. As shown in \citetalias{lia:lu:mu:oku24}, all of utilitarian, Rawlsian, egalitarian, and various other preferences are consistent with this ordering. One can then define the feasible set of group risk pairs across algorithms from the class $\mathcal{A}(\mathcal{X})$ as

and the fairness-accuracy (FA) frontier as

For finite $(\mathcal{X},\mathcal{Y})$, \citetalias{lia:lu:mu:oku24} show that $\mathcal{E}\big(\mathbb{P}, \mathcal{A}(\mathcal{X})\big)$ is a closed convex set {(we extend this convexity result to general $(\mathcal{X},\mathcal{Y})$ in the proof of Proposition (ref))} and $\mathcal{F}\big(\mathbb{P}, \mathcal{A}(\mathcal{X})\big)$ is a specific portion of its boundary connecting three points: the feasible point that minimizes the risk for group $r$, denoted $R$; the feasible point that minimizes the risk for group $b$, denoted $B$, and the feasible point that minimizes the absolute difference in group risks, denoted $F$.\footnote{Ties are broken in favor of the other group's risk for $R$ or $B$. If there are multiple feasible points that minimize the absolute difference in group risks, $F$ is chosen to be the one that has the lowest risk for both groups.} Adapting \citetalias{lia:lu:mu:oku24} nomenclature, call $\big(\mathbb{P}, \mathcal{A}(\mathcal{X})\big)$ group-balanced if $\mathcal{E}\big(\mathbb{P}, \mathcal{A}(\mathcal{X})\big)$ has $R$ and $B$ such that either $R=B=F$, or $e_r < e_b$ at $R$ and $e_r > e_b$ at $B$; call $\big(\mathbb{P}, \mathcal{A}(\mathcal{X})\big)$ $r$-skewed if $e_r<e_b$ at $R$ and $e_r\leq e_b$ at $B$, and $b$-skewed if $e_r\geq e_b$ at $R$ and $e_r > e_b$ at $B$.

To ease notation, we drop the dependence of $\mathcal{E}$ and $\mathcal{F}$ on $\big(\mathbb{P}, \mathcal{A}(\mathcal{X})\big)$ unless explicitly needed. Figure (ref) illustrates these sets and key points on them under a smoothness condition stated in Assumption (ref). \citetalias{lia:lu:mu:oku24} (Theorem 1) show that the shape of $\mathcal{F}$ depends entirely on whether $\big(\mathbb{P}, \mathcal{A}(\mathcal{X})\big)$ is group-balanced or $g$-skewed. If group-balanced, $\mathcal{F}$ is the curve connecting $R$ and $B$, coinciding with the Pareto frontier (panel (a)); if $g$-skewed, $\mathcal{F}$ connects $F$ and the feasible point minimizing risk for group $g$ (panels (b) and (c) for the $r$-skewed case; omitted panels for the $b$-skewed case). \citetalias{lia:lu:mu:oku24} (Proposition 6) further show that excluding group identity as an algorithmic input is uniformly welfare-reducing under strict group balance, where $R$ and $B$ are strictly separated by the 45-degree line.

{

} Notation. We denote by $\|\cdot\|_E$, $\|\cdot\|_{L^2(\mathbb{P})}$, $\|\cdot\|_{\infty}$, respectively, the Euclidean norm, the $L^2$-norm under the probability measure $\mathbb{P}$, and the $L^\infty$-norm (or sup-norm). For a vector $\mathbf{a}$, let $\|\mathbf{a}\|_{L^2(\mathbb{P})}\equiv\bigl\|\|\mathbf{a}\|_E\bigr\|_{L^2(\mathbb{P})}$ and $\|\mathbf{a}\|_{\infty}$ be the supremum over the largest component of $\mathbf{a}$. For a matrix $\mathbf{A}$, let $\|\mathbf{A}\|_{\max}$ denote its max norm (the maximum absolute value among its entries). For two sequences $a_n$ and $b_n$, $a_n\lesssim b_n$ means $a_n\leq c\cdot b_n$ for some constant $c>0$.

Support Function Based Characterizations

We leverage the convexity of the feasible set $\mathcal{E}$ to characterize it by its support function and express the points $R$, $B$, $F$, and the FA-frontier $\mathcal{F}$ in Eq. (ref) through this support function. We begin by observing that Eq. (ref) and the law of iterated expectations yield

where we denote $\theta_d^g(X)\equiv\mathbb{E}[\ell(d,Y)\mathds{1}\{G=g\}|X]$ the (measurable) conditional expectation of $L_d^g\equiv\ell(d,Y)\mathds{1}\{G=g\}$ given $X$; $\mu_g\equiv \mathbb{P}(G=g)$ the population proportion of group $g\in\{r,b\}$; and the expectation in Eq. (ref) is taken with respect to the population marginal distribution of the covariates, $\mathbb{P}(X)$. To make sure that Eq. (ref) is well defined, we assume:

Throughout, we let $\boldsymbol{\theta}(X)\equiv[\theta_1^r(X) ~~ \theta_0^r(X) ~~ \theta_1^b(X) ~~ \theta_0^b(X)]^\intercal$ and

Support Function of the Feasible Set

Given Eqs. (ref)-(ref)-(ref), $\mathcal{E}$ can be written as

with $\operatorname{conv}(\cdot)$ the convex hull of the set in parentheses, $\Lambda(X)\equiv \operatorname{conv}\left(\{\boldsymbol{\theta}_0(X), \boldsymbol{\theta}_1(X)\}\right)$ a random interval, $\mathcal{M}$ in Eq. (ref), and $\mathbf{E}\left[\mathcal{M}\Lambda(X)\right]$ the Aumann expectation of the scaled random interval $\mathcal{M}\Lambda(X)$ mol:mol18.

As the set $\mathcal{E}$ is non-empty, compact, and convex, its support function in each direction $q=[q_1~q_2]^\intercal\in\mathbb{S}^1\equiv\{v\in\mathbb{R}^2:\, \Vert v\Vert_E=1\}$, defined as

uniquely characterizes $\mathcal{E}$ through the identity roc97

We next provide a closed-form expression for $h_{\mathcal{E}}(q)$.

The support function $h_\mathcal{E}(q)$ is our key inferential tool. One can garner the intuition behind its closed-form expression by rewriting $e_g(a)=\mathbb{E}\left[\tfrac{\theta_0^g(X)}{\mu_g}\right]+\mathbb{E}\left[a(X)\tfrac{\theta_1^g(X)-\theta_0^g(X)}{\mu_g}\right]$ and

The maximum in the above expression is achieved by the algorithm $a^{\texttt{opt}}(X;q)=\mathds{1}\{k(\boldsymbol{\theta}(X),\cMq)>0\}$, yielding Eq. (ref) upon applying the law of iterated expectations.

Best Group-Specific Points on the FA-Frontier

We next define the support set of $\mathcal{E}$ in direction $q\in\mathbb{S}^1$:

i.e., $\mathcal{S}_\mathcal{E}(q)$ is the intersection between $\mathcal{E}$ and the hyperplane with normal vector $q$ and constant $h_{\mathcal{E}}(q)$, hence collecting the extreme point(s) of $\mathcal{E}$ in direction $q$. To derive a closed-form expression for $\mathcal{S}_\mathcal{E}(q)$, we impose the following assumption:

Assumption (ref) is a margin condition that guarantees sufficient smoothness in the distribution of $\boldsymbol{\theta}_1(X)-\boldsymbol{\theta}_0(X)$ for us to show that $h_\mathcal{E}(q)$ is differentiable in $q\in\mathbb{S}^1$ and consequently $\mathcal{S}_\mathcal{E}(q)$ includes a single element in each direction $q$ sch93. We further show that $\mathcal{S}_\mathcal{E}(q)$ equals the gradient of the support function $h_{\mathcal{E}}(\cdot)$ with respect to $q$. We denote by $\mathcal{S}_\mathcal{E}(q)$ both the singleton set and its only element.

Let $\mathfrak{u}_1\equiv[-1\,\,\,0]^\intercal$ and $\mathfrak{u}_2\equiv[0\,\,\,-1]^\intercal$. The best group-specific points satisfy:

Fairest Point on the FA-Frontier

{

} Determining the coordinates of the fairest point $F$ is more laborious, as they depend on whether $\mathcal{E}$ lies entirely above, entirely below, or on top of the 45-degree line. Figure (ref) illustrates all possible locations of the feasible set $\mathcal{E}$ relative to the 45-degree line. When $\mathcal{E}$ lies entirely on one side of the 45-degree line, as depicted in panels (b) and (d), $F$ is the support set of $\mathcal{E}$, respectively, in directions $\mathfrak{u}_2-\mathfrak{u}_1=[1\,\,-1]^\intercal$ and $\mathfrak{u}_1-\mathfrak{u}_2=[-1\,~~\,1]^\intercal$:

Complications arise when $\mathcal{E}$ intersects with the 45-degree line, as depicted in panels (a), (c), and (e) of Figure (ref). In this case, the direction at which we can obtain $F$ as the support set of $\mathcal{E}$ is difficult to determine. To circumvent this challenge, we propose a different approach. We focus on the convex set that results when $\mathcal{E}$ intersects the 45-degree line:

The new set $\Tilde{\mathcal{E}}$ is depicted in panels (a), (c), and (e) of Figure (ref) as an orange line segment. In these cases, $F$ is the support set of $\Tilde{\mathcal{E}}$ in direction $\mathfrak{u}_1$ with identical values for its two coordinates. Hence,

We are left with providing an expression for $h_{\Tilde{\mathcal{E}}}(q)$. When $\mathcal{E}$ intersects $\mathcal{H}_{45}$,

where the first equality follows from roc97, and the second follows from the fact that $h_{\mathcal{H}_{45}}(p_2)$ is bounded from above only along the direction $p_2=c[1\,\,-1]^\intercal$ for any scalar $c\in\mathbb{R}$, in which case $h_{\mathcal{H}_{45}}(p_2)=0$. Importantly, the last expression in Eq. (ref) is always well defined, regardless of whether $\mathcal{E}$ intersects with $\mathcal{H}_{45}$ or not; the infimum equals a bounded scalar in the case of intersection and $-\infty$ otherwise.\footnote{To see this, let $q(c)\equivq-c[1~~-1]^\intercal$ and note that $h_\mathcal{E}(q(c))=\|q(c)\|_E\cdoth_\mathcal{E}\big(\frac{q(c)}{\|q(c)\|_E}\big)$. When $\mathcal{E}$ intersects with $\mathcal{H}_{45}$, $h_\mathcal{E}\big(\frac{q(c)}{\|q(c)\|_E}\big)$ is bounded and nonnegative along the sequences $c\to\infty$ and $c\to-\infty$, but when $\mathcal{E}$ and $\mathcal{H}_{45}$ are disjoint, it takes negative value along one of these sequences, yielding $\inf_c h_\mathcal{E}(q(c))=-\infty$.}

Support Function-Based Characterization of the FA-Frontier

We next show that the FA-frontier put forward by \citetalias{lia:lu:mu:oku24} and reproduced in our Eq. (ref) can be characterized using the support function of the feasible set $\mathcal{E}$ and that of an auxiliary set that we introduce in this subsection.

Given an algorithm $a^*\in\mathcal{A}(\mathcal{X})$ that induces the risk pair $e^*=[e_r^*, e_b^*]^\intercal\in\mathcal{E}$, let

{

} denote the set of risk allocations $e\in\mathbb{R}^2$---whether or not they are feasible---that are both weakly more accurate and weakly fairer than $e^*$. The set $\mathcal{C}(e^*)$, depicted in the two panels of Figure (ref) as the shaded green regions corresponding to two different values of $e^*$, is a closed and convex subset of $\mathbb{R}^2$. Panel (a) depicts a case where $e^*\in\mathcal{F}$, whereas panel (b) depicts a case where $e^*\notin\mathcal{F}$. The key insight from the figure, which we prove can be used to characterize the FA-frontier using the support function of $\mathcal{E}$ and that of $\mathcal{C}(e^*)$, is that for any $e^*\in\mathcal{F}$, the sets $\mathcal{E}$ and $\mathcal{C}(e^*)$ can be properly separated sch93, while in the case where $e^*\notin\mathcal{F}$ they cannot.

Remarkably, this characterization of the FA-frontier does not require knowledge of whether $\mathcal{E}$ intersects the $45^o$ line, or whether there is group balance or skew. As we show in Sections (ref)-(ref), the characterization of the FA-frontier in Eq. (ref) is very helpful for estimation and inference, and for testing whether there exists an LDA to a given algorithm.

Support Function Estimator and Its Asymptotic Distribution

In practice, the policymaker does not have perfect knowledge of $\mathbb{P}$. Hence, $h_{\mathcal{E}}(q)$ can only be estimated from a finite sample. Let a sample of size $n$, $\{(Y_i,G_i,X_i)\}_{i=1}^n$, drawn independently and identically from $\mathbb{P}$, be available. Recall Eq. (ref): $h_\mathcal{E}(q)=\mathbb{E}\bigl[(\cMq)^\intercal\mathbf{L}_0+(\cMq)^\intercal(\mathbf{L}_1-\mathbf{L}_0)\mathds{1}\{k(\boldsymbol{\theta}(X),\cMq)>0\}\bigr]$ with $k(\boldsymbol{\theta},\cMq)\equiv(\cMq)^\intercal\bigl(\boldsymbol{\theta}_1(X)-\boldsymbol{\theta}_0(X)\bigr)$, so that $\boldsymbol{\theta}(X)\equiv[\theta_1^r(X) ~ \theta_0^r(X) ~ \theta_1^b(X) ~ \theta_0^b(X)]^\intercal$ enters the expression for $h_\mathcal{E}(q)$ only through $(\theta_1^r(X)-\theta_0^r(X))$ and $(\theta_1^b(X)- \theta_0^b(X))$. We propose estimating $h_{\mathcal{E}}(q)$ by first estimating the finite dimensional parameters $\mathcal{M}$ by sample averages and the nuisance functions $\Delta\boldsymbol{\theta}(X)\equiv[(\theta_1^r(X)-\theta_0^r(X)) ~~ (\theta_1^b(X)- \theta_0^b(X))]^\intercal$ by flexible machine learning methods (allowing the complexity of the parameter space containing the estimator to grow with sample size), and then plugging their estimators, denoted $\widehat{\mathcal{M}}$ and $\widehat{\Delta\boldsymbol{\theta}}(X)$, into the sample analogue of Eq. (ref). Following the literature on debiased machine learning new94, che:che:dem:duf:han:whi:rob18, sem:che21, che:esc:ich:new:rob22, ich:new22, we show that, at the population $\mathcal{M}$, the moment in Eq. (ref) is Neyman-orthogonal and hence “insensitive” to the errors in the first-stage estimation of $\Delta\boldsymbol{\theta}$ ney79, ney95. We use sample splitting to relax the otherwise needed Donsker condition that limits the complexity of the relevant parameter space bic82, rob:li:tch:van08, rob:li:muk:tch:van17, and we account for the estimation error of $\mathcal{M}$.

Recall that $\mathbf{L}_d\equiv[L_d^r~~L_d^b]^\intercal$, with $L_d^g\equiv \ell(d,Y)\mathds{1}\{G=g\}$, and $\theta_d^g(X)=\mathbb{E}[L_d^g|X]$. Let $\DeltaL^g\equivL_1^g-L_0^g$, so that $\Delta\theta^g(X)=\mathbb{E}[\DeltaL^g|X]$. To learn the nuisance parameter $\Delta\theta^g(X)$, the effective label that, given $X$, we train machine learners to predict is $\DeltaL^g$. Let $\Theta$ denote the convex nuisance parameter space (a subset of a vector space with $L^2(\mathbb{P})$ norm) to which $\Delta\boldsymbol{\theta}=[\Delta\theta^g(X)~~\Delta\theta^b(X)]^\intercal$ belongs and $\Delta\boldsymbol{\vartheta}\equiv[\Delta\vartheta^r(X)~~\Delta\vartheta^b(X)]^\intercal$ be a generic element from $\Theta$ (to simplify notation, we drop the dependence of $\boldsymbol{\theta}$ on $X$ unless explicitly needed). For the $i$-th observation and a given $2\times2$ diagonal matrix $\mathring{\mathcal{M}}$, we note that

and using Eq. (ref) to recognize that we estimate $\Delta\boldsymbol{\theta}$ while keeping the notation as close as possible to that in Eq. (ref), we define the mapping $\zeta_i(\mathring{\mathcal{M}}q;\,\cdot\,) : \Theta \to \mathbb{R}$ as

By Proposition (ref), $h_{\mathcal{E}}(q)=\mathbb{E}[\zeta_i(\cMq;\boldsymbol{\theta})]$. We next show that, when $\mathring{\mathcal{M}}$ is fixed at the population $\mathcal{M}$, the score function $\zeta_i(\cMq;\boldsymbol{\vartheta})$ is Neyman-orthogonal at $\Delta\boldsymbol{\vartheta}=\Delta\boldsymbol{\theta}$.

Intuitively, Proposition (ref) shows that, under its maintained assumptions, the first-order mistake in the sign of $k(\boldsymbol{\theta}(X),\cMq)$ due to the estimation error in $\Delta\boldsymbol{\theta}(X)$ is negligible.

We next show that, provided $n^{1/4}$-consistent first-stage estimators are available, the estimated support function using sample splitting and cross-fitting, as described in Definition (ref) below, converges to a Gaussian process uniformly in $q\in\mathbb{S}^1$. The proof requires showing the residual in estimating the indicator functions is bounded by the quadratic rate of convergence of the nuisance parameter, which we establish under this assumption:

Assumption (ref) requires $\Delta\boldsymbol{\theta}(X)$ to be linearizable in a known set of covariates $(X_1,X_2)$ and that $\alpha^b\cdot\beta^r\neq\alpha^r\cdot\beta^b$. This assures that for each $q\in\mathbb{S}^1$, $k(\boldsymbol{\theta},\cMq)$ depends on at least one of $X_1$ or $X_2$, so that we can employ a proof technique similar to that in sem23 to show that the bias induced by errors in estimating the sign of $k(\boldsymbol{\theta},\cMq)$ is bounded by the desired quadratic $L^2$-rate of the nuisance estimator, accommodating a menu of flexible machine learners in the first stage. For example, one can adapt the method in rob88 to estimate $\Delta\theta^g$. Example machine learning methods that, under suitable conditions, are $n^{1/4}$-consistent in the $L^2$-norm include $\ell_1$-penalized methods, boosting, regression trees and forests, and neural nets che:che:dem:duf:han:whi:rob18. Assumption (ref) can be eliminated at the price of using a more restrictive class of machine learning methods, by deriving rate bounds that depend on the squared $L^\infty$-rate of the first-stage estimation bel:che:fer:han17, sem23. Alternatively, Assumption (ref) can be eliminated by smoothing the indicators che:aus:syr23, par24 at the cost of introducing an additional tuning parameter that controls the degree of smoothing.\footnote{An earlier version of this paper liu:mol24v1 obtains $\sqrt{n}$-Gaussianty of the second-stage estimator without imposing Assumption (ref), but under the classical Donsker condition that limits $\Theta$'s complexity.} Importantly, Assumption (ref) implies the margin condition in Assumption (ref) with $m=1$ whenever $(X_1,X_2)$ is continuously distributed with a bounded density and $X_{[3:d_X]}$ can be either continuous or discrete; see Assumption (ref) and Proposition (ref).

We next define cross-fitting, adapting Definition 3.2 in che:che:dem:duf:han:whi:rob18:

In Eq. (ref), to simplify notation and keep it as close as possible to that in Eq. (ref) and $h_{\mathcal{E}}(q)=\mathbb{E}[\zeta_i(\cMq;\boldsymbol{\theta})]$, we use the shorthand notation $\widehat{\Delta\boldsymbol{\theta}}(X_i)\equiv\widehat{\Delta\boldsymbol{\theta}}_k(X_i)$ for $i\in I_k$, suppress the dependence of $\widehat{h}_{\mathcal{E}}(q; \widehat{\boldsymbol{\theta}})$ on $\widehat{\mathcal{M}}$, and adapt Eqs. (ref)-(ref) to let

Our main asymptotic result shows that $\widehat{h}_{\mathcal{E}}(q; \widehat{\boldsymbol{\theta}})$ converges to a Gaussian process uniformly in the direction $q\in\mathbb{S}^1$, where the score $\zeta_i(\cMq;\boldsymbol{\vartheta})$ in Eq. (ref) evaluated at $\Delta\boldsymbol{\vartheta}=\Delta\boldsymbol{\theta}$ is the influence function that governs the part of the asymptotic distribution of $\widehat{h}_{\mathcal{E}}(q; \widehat{\boldsymbol{\theta}})$ due to $\widehat{\Delta\boldsymbol{\theta}}$, and the the remaining part is attributed to estimating $\mathcal{M}$:

Estimation and Inference for the Frontier

Estimation and Inference for the Feasible Set

We propose an estimator of the set $\mathcal{E}$ based on $\widehat{h}_\mathcal{E}(q)$ and Eq. (ref), given by

which is convex, almost surely compact, and non-empty with probability approaching one if $\mathcal{E}$ has a non-empty interior. As the Hausdorff distance between two non-empty convex and compact sets $A,B\in\mathbb R^d$, denoted $\mathbf{d}_H(A,B)$, equals the uniform distance between their support functions mol:mol18, when $\widehat{\mathcal{E}}$ is non-empty we have

By Theorem (ref) and the continuous mapping theorem, $\mathbf{d}_H(\widehat{\mathcal{E}},\mathcal{E})\xrightarrow[]{p} 0$ and $\sqrt{n}\mathbf{d}_H(\widehat{\mathcal{E}},\mathcal{E})\xrightarrow[]{d}\sup_{q\in\mathbb{S}^1}|\mathbb{G}[ \zeta_i^*(\cMq;\boldsymbol{\theta})]|$. Hence, asymptotically valid tests of hypotheses about $\mathcal{E}$, and confidence sets covering it, can be obtained as in ber:mol08.

Estimation and Inference for the FA-Frontier

As shown in \citetalias{lia:lu:mu:oku24} (Theorem 1), when $\big(\mathbb{P}, \mathcal{A}(\mathcal{X})\big)$ is group-balanced, $\mathcal{F}$ is the curve connecting $R$ and $B$ and coincides with the Pareto frontier (panel (a) in Figure (ref)). However, when $\big(\mathbb{P}, \mathcal{A}(\mathcal{X})\big)$ is $g$-skewed, $\mathcal{F}$ is the curve connecting $F$ with the feasible point that minimizes the risk for group $g$ (panels (b)-(e) in Figure (ref)). A further challenge is that, while it is simple to express $F$ through $\mathcal{S}_\mathcal{E}$ when $\mathcal{E}$ is fully contained in one of the two half-spaces defined by the 45-degree line $\mathcal{H}_{45}$, as shown in Eqs. (ref)-(ref) (panels (b) and (d) in Figure (ref)), characterizing $F$ is not as straightforward when $g$-skew occurs but $\mathcal{E}\cap\mathcal{H}_{45}\neq\emptyset$.

We therefore leverage the characterization of the FA-frontier $\mathcal{F}$ in Proposition (ref), whereby $\mathcal{F} = \left\{e\in\mathcal{E}:~\left[\max_{q\in\tilde\mathbb{S}^1}(-h_{\mathcal{C}(e)}(q)-h_\mathcal{E}(-q))\right]_-=0\right\}$, and the definition of the set $\mathcal{C}(\cdot)$ in Eq. (ref) to sidestep these difficulties. We propose to estimate $\mathcal{F}$ using

where $[\,\,\cdot\,\,]_+\equiv\max\{\,\cdot\,, 0\}$, $\kappa_n=o(\sqrt{n})$ is a sequence that diverges to infinity, and $\mathbb{B}_C\equiv\{e\in\mathbb{R}^2:\Vert e\Vert_E \le C\}$, with $\mathcal{E}\subset \mathbb{B}_C$ by Assumption (ref) and $C<\infty$ a constant pinned down by $c_1,c_2$ defined in Assumption (ref). The first maximization problem in Eq. (ref) is the sample analog to the requirement that $e\in\mathcal{E}$; the second is the sample analog to the requirement that $\mathcal{E}$ and $\mathcal{C}(e)$ can be properly separated. To implement this estimator, we derive a closed-form expression for the support function of $\mathcal{C}(e)$ in direction $q=[q_1,~q_2]^\intercal$. Only two points are “active” for evaluating $h_{\mathcal{C}(e)}(q)$: $[\min\{e_r,2e_b-e_r\},\,\,e_b]^\intercal$ and $[e_r,\,\,\min\{e_b, 2e_r-e_b\}]^\intercal$, which correspond to the points $e$ and $[e_r,\,\,2e_r-e_b]^\intercal$ (respectively, $[2e_b-e_r,\,\,e_b]^\intercal$ and $e$) when $e$ is above (respectively, below) the 45-degree line as shown in Figure (ref)-panel (a) (respectively, panel (b)). By Proposition (ref), it is without loss of generality to focus on $q\in\tilde{\mathbb{S}}^1$. Hence, the support function of $\mathcal{C}(e)$ at any $q\in\tilde{\mathbb{S}}^1$ equals

We next propose a test statistic for the null hypothesis $\mathsf{H}_0:~e\in\mathcal{F}$, against the alternative $\mathsf{H}_A:~e\notin\mathcal{F}$. Our test statistic is given by

As shown in the proof of Proposition (ref), if Eq. (ref) is satisfied and consequently $\mathcal{E}$ has no kinks, the large sample distribution of $T_n^\mathcal{F}(e)$, denoted $\psi^{\mathcal{F}}(e)$, simplifies to

where $q^*_{\mathbb{S}^1}(e)\equiv\arg\max_{q\in\mathbb{S}^1}q^\intercal e -h_\mathcal{E}(q)$. The quantiles of this distribution can be estimated by standard methods. We build a confidence set for the elements of $\mathcal{F}$ by test inversion:

where for any $\beta\in(0,1)$, $c_\beta^{\mathcal{F}}$ is the $\beta$-quantile of $\psi^{\mathcal{F}}$. When Eq. (ref) is not assumed and $\mathcal{E}$ has kinks, the expression for $\psi^{\mathcal{F}}(e)$ is more complex and given in Eq. (ref). In this case, restrictions that would guarantee uniform continuity and strict increasing properties for $\psi^{\mathcal{F}}(e)$ are harder to verify, and we follow and:shi13 to replace the confidence set in Eq. (ref) with $\mathcal{CS}_n(\mathcal{F}) = \left\{e\in\mathbb{B}_C:T_n^\mathcal{F}(e)\le c_{1-\alpha+\varsigma}^{\mathcal{F}}(e)+\varsigma\right\}$, for $\varsigma>0$ an arbitrarily small constant. In this case, the critical value can be approximated through bootstrap methods, as in Procedure (ref) in Section (ref).

We next establish consistency of our estimator $\widehat{\mathcal{F}}$ and validity of the confidence set, building respectively on che:hon:tam07 and fan:san19.

Estimation and Inference for the Pareto Frontier

The Pareto Frontier ($\mathcal{PF}$) is the lower boundary of the feasible set $\mathcal{E}$ connecting points $R$ and $B$ (see Figure (ref)). Denoting $\mathbb{Q}\equiv\left\{q\in\mathbb{S}^1:q=[\cos\gamma~~\sin\gamma]^\intercal,~\gamma\in[\pi,3/2\pi]\right\}$, we can express $\mathcal{PF}$ through the support set $\mathcal{S}_\mathcal{E}(\cdot)$, or equivalently through the support function:

where the first condition in Eq. (ref) enforces that $e\in\mathcal{E}$ and the second that it belongs to the supporting hyperplane of $\mathcal{E}$ in a direction $q\in\mathbb{Q}$. Knowing the directions that determine the support points comprising $\mathcal{PF}$ simplifies the expression for it in Eq. (ref) relative to Eq. (ref). It also allows for an estimator, put forward in Proposition (ref) below, based on the support set characterization in Eq. (ref) that is free from the tuning parameter $\kappa_n$, which instead is needed for the estimator of $\mathcal{F}$ in Eq. (ref). Let $\boldsymbol{\zeta}_{\mathcal{S},i}(\hcMq;{\widehat{\boldsymbol{\theta}}})\equiv\widehat{\mathcal{M}}\mathbf{L}_{0_i}+\widehat{\mathcal{M}}(\mathbf{L}_{1_i}-\mathbf{L}_{0_i})\cdot\mathds{1}\bigl\{k\bigl({\widehat{\boldsymbol{\theta}}(X_i)},\hcMq\bigr)>0\bigr\}$, with $k\bigl(\widehat{\boldsymbol{\theta}}(X_i),\hcMq\bigr)=q^\intercal\widehat{\mathcal{M}}\widehat{\Delta\boldsymbol{\theta}}(X_i)$ as in Eq. (ref) and\footnote{As in Eq. (ref), this expression is a shorthand for $\widehat{\mathcal{S}}_\mathcal{E}(q;\widehat{\boldsymbol{\theta}})\equiv\frac{1}{K}\sum_{k\in[K]}\left(\frac{1}{n/K}\sum_{i\in I_k}\boldsymbol{\zeta}_{\mathcal{S},i}({\widehat{\mathcal{M}}}q;\widehat{\boldsymbol{\theta}}_k)\right)$.}

Proposition (ref) delivers an estimator for $\mathcal{PF}$ based on $\widehat{\mathcal{S}}_\mathcal{E}(q;\widehat{\boldsymbol{\theta}})$ in Eq. (ref) and establishes its Hausdorff-consistency.

We next provide a method to test, for a given $e\in\mathbb{R}^2$,

using the characterization of $\mathcal{PF}$ in Eq. (ref). We first propose a test statistic and derive its asymptotic distribution in Proposition (ref) below, building on results in kai16.

If Eq. (ref) is satisfied and the set $\mathcal{E}$ has no kinks (see Remark (ref)), $Q_{\mathbb{S}^1}^*(e)$ and $Q_{\mathbb{Q}}^*(e)$ are singletons that can be consistently estimated through standard methods, and the limit distribution in Eq. (ref) coincides with that in Eq. (ref). If instead kinks are not ruled out for $q\in\mathbb{Q}$, $Q_{\mathbb{S}^1}^*(e)$ and $Q_{\mathbb{Q}}^*(e)$ may not be singletons. In this case we can consistently estimate these sets kai16 as

where $\mathbb{T}\in\{\mathbb{S}^1,\mathbb{Q}\}$ and $\kappa_n=o(\sqrt{n})$ is a sequence that diverges to infinity. We use Procedure (ref)-Step 1 to obtain a valid bootstrap-based approximation to the Gaussian process in Theorem (ref). Denote $\mathbb{P}^*$ this bootstrap distribution conditional on $\{(Y_i,G_i,X_i)\}_{i=1}^n$. Let

When $Var(\mathbf{L}_d|X)$ is positive definite, it follows that if $\alpha\in(0,0.5)$,

kai16. A confidence set that covers each $e\in\mathcal{PF}$ with asymptotic probability at least equal to $(1-\alpha)$ can be obtained by test inversion.

Algorithms Yielding a Risk Allocation on the Frontier

An algorithm designer or a regulator may wonder if one can characterize the algorithm yielding a specific point on the FA-frontier $\mathcal{F}$ or the Pareto frontier $\mathcal{PF}$. It turns out that, using our support function approach, the answer to this question is affirmative, and particularly simple for points in $\mathcal{PF}$.

Suppose we have two data sets: one for training the algorithm and the other for evaluating it. We denote the training sample as $\{(\tilde{Y}_i,\tilde{G}_i,\tilde{X}_i)\}_{i=1}^{n_1}$ and the evaluation sample as $\{(Y_j,G_j,X_j)\}_{j=1}^{n_2}$, with $\frac{n_1}{n_2}\to c$ for some positive constant $c$ and both samples drawn from the same distribution $\mathbb{P}$. Use the training data to estimate $\Delta\boldsymbol{\theta}$ via machine learning and $\mathcal{M}$ by sample averages as described in Section (ref) and denote the resulting estimators $\widehat{\Delta\boldsymbol{\theta}}_{n_1}$ and $\widehat{\mathcal{M}}_{n_1}$, where the subscript $n_1$ indicates that it is based on the training sample. Let

with $k\bigl(\widehat{\boldsymbol{\theta}}_{n_1}(X_j),\widehat{\mathcal{M}}_{n_1}q\bigr)$ as in Eq. (ref). Note that $\widehat{a}_{n_1}(X_j;q)$ only takes as input the covariates, $X_j$, and does not depend on group identity $G_j$. Consequently, individuals with the same covariates are assigned the same treatment irrespective of their group. Nonetheless, information on group identity contained in the training data is used in the prediction model (to separately estimate $\Delta\theta^g$, $g\in\{r,b\})$, thereby offering a compromise between a utilitarian perspective as in man:mul:ven23 which advocates for group-aware decisions, and proponents of group-blind decision making.

{Then, for any $q\in\mathbb{S}^1$, algorithm $\widehat{a}_{n_1}(\,\cdot\,;q)$ leads to a loss for each observation $i$ in the evaluation sample equal to $\ell(0,Y_j)+(\ell(1,Y_j)-\ell(0,Y_j))\cdot\mathds{1}\bigl\{k\bigl(\widehat{\boldsymbol{\theta}}_{n_1}(X_j),q\bigr)>0\bigr\}$. Therefore, the average losses for the $r$ group and for the $b$ group in the evaluation sample equal

} By the same argument as in the proof of Proposition (ref), it follows that

Hence, the algorithm in Eq. (ref) with $q\in\mathbb{Q}$ gives consistent estimators of points on $\mathcal{PF}$.

Algorithms that return consistent estimators of points on $\mathcal{F}$ (other than those on $\mathcal{PF}$) are harder to obtain, as one needs to determine which direction $q$ corresponds to a point $e\in\mathcal{F}$ and plug that direction in the algorithm in Eq. (ref). Suppose Eq. (ref) is satisfied and $\mathcal{E}$ has no kinks. Use the same DML construction in Section (ref) applied to the training sample to obtain an estimator of $h_\mathcal{E}(q)$, denoted $\widehat{h}_{\mathcal{E},n_1}(q;\widehat{\boldsymbol{\theta}}_{n_1})$, as in Eq. (ref), and an estimator of $\mathcal{F}$, denoted $\widehat{\mathcal{F}}_{n_1}$, as in Eq. (ref). {Select $\widehat{e}_{n_1}\in\widehat{\mathcal{F}}_{n_1}$ through a projection method detailed in the proof of Proposition (ref) so that $\Vert\widehat{e}_{n_1}-e\Vert=o_p(1)$ for some $e\in\mathcal{F}$ and let}

Let $\widehat{\mathcal{S}}_{\mathcal{E},n_2}(\widehat{q}^*_{n_1}(\widehat{e}_{n_1});\widehat{\boldsymbol{\theta}}_{n_1})$ be as in Eq. (ref) but with $\widehat{q}^*_{n_1}(\widehat{e}_{n_1})$ replacing $q$. The next proposition establishes that $\widehat{\mathcal{S}}_{\mathcal{E},n_2}(\widehat{q}^*_{n_1}(\widehat{e}_{n_1});\widehat{\boldsymbol{\theta}}_{n_1})$ is a consistent estimator of $\mathcal{S}_\mathcal{E}(q^*_{\mathbb{S}^1}(e))\in\mathcal{F}$.

Regardless of whether one aims at obtaining points in $\mathcal{PF}$ or the entire $\mathcal{F}$, the direction $q$ can be interpreted as the vector of weights that the agent choosing the algorithm puts on each group's risk. In other words, one may think of the agent as evaluating group risks according to the welfare loss function $U(e;q)\equivq_1 e_r(a)+q_2 e_b(a)$. For example, the more the agent cares about group $r$, the closer $q$ is to $\mathfrak{u}_1$. We note that $\widehat{a}_{n_1}(X;q)$ is an empirical success rule, and leave its statistical decision theory analysis to future research.

Hypothesis Testing

In this Section we propose hypothesis tests to answer the following policy questions:

(1) Should the policymaker consider banning group identity as an input to the algorithm?

(2) Is there a less discriminatory alternative (LDA) to an existing algorithm?

We show how to express the first policy question in terms of restrictions on $h_\mathcal{E}(q)$, and we leverage Proposition (ref) to do the same for the second policy question. We then establish asymptotic validity of the corresponding testing procedures.

How to Test Whether Group Identity Should be Banned

\citetalias{lia:lu:mu:oku24} (Proposition 6) show that using $X$ only instead of $(X,G)$ as algorithmic input uniformly worsens the frontier if $R$ and $B$, obtained when only $X$ is used as input, are strictly separated by the 45-degree line.\footnote{“Uniformly worsening the frontier” means that, under the preference relations defined in Eq. (ref), every point on the frontier $\mathcal{F}\big(\mathbb{P},\mathcal{A}(\mathcal{X})\big)$ is dominated by a point on the frontier $\mathcal{F}\big(\mathbb{P},\mathcal{A}(\mathcal{X}\times\{r,b\})\big)$.} We therefore aim at testing the null hypothesis that $R$ and $B$ lie weakly on the same side of the 45 degree line, i.e., that the difference in the two coordinates of $R$ has the same sign as that of $B$, against the alternative that they are strictly separated by the 45-degree line:

If the null in Eq. ((ref)) is rejected, the policymaker should not ban group identity $G$ as algorithm's input. A Type-I error amounts to the case where one concludes that there is strict group-balance, while instead weak group-skew holds. As a consequence, one does not ban $G$ as input to the algorithm, thinking that banning $G$ is uniformly welfare-reducing, while instead depending on the preferences of the designer it might not be the case. Another interpretation of this test amounts to determining, based on whether $\mathsf{H}_0$ is rejected or not, if one can justify implementing algorithms that lead to Pareto-dominated risks based on the designer's preference over fairness and accuracy. When $\mathsf{H}_0$ holds true, this justification is possible, as the frontier includes an upward-sloping segment (e.g., Panels (b)-(c) of Figure (ref)). On the other hand, when $\mathsf{H}_A$ holds true, such justification is untenable, as the frontier coincides with the Pareto frontier (e.g., Panel (a) of Figure (ref)). Hence, a Type-I error can also be interpreted as a case where one concludes that only Pareto-optimal risks should be implemented, while fairness considerations may justify Pareto-dominated risks.

{As discussed in Remark (ref), carrying out inference if we use DML to directly estimate both coordinates of the points $R$ and $B$ is difficult.} Yet, using Eqs. (ref) and (ref), we can represent these points through moment equalities and inequalities that involve $h_\mathcal{E}(q)$ only:

where in Eq. (ref), the equality constraint for $R$ restricts it to have horizontal coordinate equal to that of $\mathcal{S}_\mathcal{E}(\mathfrak{u}_1)$, as per Eq. (ref), and the continuum of inequality constraints indexed by $q\in\mathbb{S}^1$ restricts $R$ to be an element of $\mathcal{E}$. The moment constraints that define $B$ are interpreted similarly. Because the support set $\mathcal{S}_\mathcal{E}(\cdot)$ in any direction is a singleton by Proposition (ref), the moments in Eq. (ref) yield points $R$ and $B$ that coincide with Eq. (ref).

We test the null in Eq. (ref) at a given significance level $\alpha\in(0,1)$ based on our procedure to test, for given $e\in\mathbb{R}^2$, whether $e\in\mathcal{PF}$:

We note that the test in Eq. (ref) may be conservative as it is based on projection.

How to Test for the Existence of an LDA

Given an algorithm $a^*\in\mathcal{A}(\mathcal{X})$ that induces the risk pair $e^*=(e_r^*, e_b^*)\in\mathcal{E}$, call another algorithm that yields a feasible risk pair $e=(e_r,e_b)\in\mathcal{E}$ an LDA if it is at least as accurate as $a^*$ for both groups and at least as fair, with one of these inequalities strict. It follows from the characterization in Proposition (ref) and Eq. (ref) that no LDA to $a^*$ exists if and only if $\mathcal{E}$ can be properly separated from

Recall that the closed form expression for the support function of $\mathcal{C}^*\equiv\mathcal{C}(e^*)$ in direction $q=[q_1,~q_2]^\intercal\in\tilde{\mathbb{S}}^1$ is given in Eq. (ref). We then test the null hypothesis

against the alternative that $\mathsf{H}_0$ is false. Rejecting the null in Eq. (ref) means that $e^*\notin\mathcal{F}$ and there exists an LDA. We propose estimating $h_{\mathcal{C}^*}(q)$ for $q\in\tilde{\mathbb{S}}^1$ by

where for $g\in\{r,b\}$, $e_g^*$ is estimated by sample means,

and $\widehat{\mu}_g\equiv\frac{1}{n}\sum_{i=1}\mathds{1}\{G_i=g\}$. We propose the following test statistic:

$T_n^{\texttt{LDA}}$ differs from $T_n^\mathcal{F}$ in Eq. (ref) only in that $\widehat{e}^*$ is estimated in the former.

When Eq. (ref) holds and $Var(\mathbf{L}_d|X)$ is positive definite for $d\in\{0,1\},~X$-a.s., one can take $\varsigma=0$ in Proposition (ref) and the expression for $\psi^{\texttt{LDA}}$ simplifies to that in Eq. (ref). When kinks might be present and the limit distribution is not guaranteed to be continuous and strictly increasing, we take $\varsigma>0$ as in and:shi13. The derivation of $\psi^{\texttt{LDA}}$ uses the fact that $T_n^{\texttt{LDA}}$ is formed by compositions of the $\max$ and $\min$ functions in Eq. (ref) and the $\max$ function in Eq. (ref)---all of which are Hadamard directionally differentiable, as shown in fan:san19 and car:cue:alb20; since compositions preserve directional differentiability Shapiro90, an extension to the functional Delta method can be applied to Theorem (ref). However, standard bootstraps are inconsistent due to the lack of full differentiability fan:san19. As such, we leverage results in fan:san19 to approximate the distribution of $\psi^{\texttt{LDA}}$ and its quantiles via a modified multiplier bootstrap procedure detailed below, similar to that in sem23, where we denote $\phi$ any generic Hadamard directionally differentiable function, $\widehat{he^*}=\widehat{he^*}(\widehat{\boldsymbol{\theta}})\equiv[\widehat{h}_\mathcal{E}(q;\widehat{\boldsymbol{\theta}}),~\widehat{e}_r^*,~\widehat{e}_b^*]^\intercal$ the vector of estimators, and $he^*\equiv[h_\mathcal{E}(q),~e_r^*,~e_b^*]^\intercal$ the vector of truths.

The consistency of the bootstrap outlined in Procedure (ref) is stated in the following result.

The proof of Proposition (ref) shows that $\psi^{\texttt{LDA}}=\phi_{he^*}'(\mathbb{G}_{he^*})$ for a particular $\phi_{he^*}'(\cdot)$ that is the composition of the directional derivatives of the $\min$, $\max$, and $\inf$ functions that constitute $T_n^{\texttt{LDA}}$. The expression of $\phi_{he^*}'(\cdot)$, given in Eq. (ref), is complex and hence we recommend the numerical approximation approach in Step 2 of Procedure (ref). Under Proposition (ref), we can estimate $c_\beta^{\texttt{LDA}}$, the $\beta$-quantile of $\psi^{\texttt{LDA}}$, by going through the steps in Procedure (ref), where we replace $\phi(\cdot)$ by the composition of the above mentioned directionally differentiable functions accordingly. The same bootstrap procedure can be used to consistently estimate the critical values put forward to carry out inference in Section (ref). The use of the infinitesimal constant $\varsigma$ in Proposition (ref) and Procedure (ref) accounts for the possibility that the cdf of $\phi_{he^*}'(\mathbb{G}_{he^*})$ may not be continuous and increasing at $c_{1-\alpha}$.

Algorithms Yielding an LDA

{

} An algorithm designer or a regulator may wonder if one can characterize the algorithms yielding LDAs to a given algorithm $e^*$. It turns out that it is possible to do so, by combining the algorithm that we put forward in Eq. (ref) with a careful use of our characterization of $\mathcal{F}$ in Proposition (ref). We illustrate the idea for the case that Eq. (ref) holds and $\mathcal{E}$ has no kinks. For a given risk pair $e^*$ induced by an algorithm $a^*\in\mathcal{A}$ such that $e^*\notin\mathcal{F}$, by definition the set of risk pairs $e\in\mathcal{F}$ such that $e>_{FA}e^*$ are:

The set $\mathcal{F}^*$ is depicted in Figure (ref) as the orange portion of $\mathcal{F}$ that intersects with $\mathcal{C}(e^*)$. Using the same notation as in Section (ref), we can use a training sample $\{(\tilde{Y}_i,\tilde{G}_i,\tilde{X}_i)\}_{i=1}^{n_1}$ and a construction that mimics the procedure to build the estimator of $\mathcal{F}$ in Eq. (ref) to obtain a consistent estimator $\widehat{\mathcal{F}}^*_{n_1}$ of $\mathcal{F}^*$ (the consistency of this estimator can be established through the same steps as in the proof of Proposition (ref)). {Select $\widehat{e}_{n_1}\in\widehat{\mathcal{F}}_{n_1}$ through a projection method detailed in the proof of Proposition (ref) and denote by $e\in\mathcal{F}^*$ the point to which it converges.} Let $\widehat{q}^*_{n_1}(\widehat{e}_{n_1})$ be the consistent estimator of $q^*_{\mathbb{S}^1}(e)$ defined in Eq. (ref). Let $\widehat{\mathcal{S}}_{\mathcal{E},n_2}(\widehat{q}^*_{n_1}(\widehat{e}_{n_1});\widehat{\boldsymbol{\theta}}_{n_1})$ be defined as in Eq. (ref) but with $\widehat{q}^*_{n_1}(\widehat{e}_{n_1})$ replacing $q$. Then $\widehat{\mathcal{S}}_{\mathcal{E},n_2}(\widehat{q}^*_{n_1}(\widehat{e}_{n_1});\widehat{\boldsymbol{\theta}}_{n_1})$ is a consistent estimator of $\mathcal{S}_\mathcal{E}(q^*_{\mathbb{S}^1}(e))\in\mathcal{F}^*$, by the same argument used to establish Proposition (ref).

Distance to the Fairest Point

In this section, we propose a method to build a confidence interval for the distance between the risk $e^*$ induced by a given algorithm and the fairest point $F$ on the frontier, denoted $\rho(e^*,F)$, where $\rho$ is a Hadamard directionally differentiable distance function (e.g., the Euclidean distance, Manhattan distance, Chebyshev distance, etc.). This can inform a decision maker of the relative merits in promoting equity and achieving business efficiency of different algorithms by comparing the confidence intervals on their distance to $F$.

Recall $\mathcal{H}_{45}\equiv\{e\in\mathbb{R}^2:e_r=e_b\}$ denotes the 45-degree line; let $\mathcal{H}_{45}^+\equiv\{e\in\mathbb{R}^2:e_r<e_b\}$ and $\mathcal{H}_{45}^-\equiv\{e\in\mathbb{R}^2:e_r>e_b\}$ denote, respectively, the open halfspace above and below the 45-degree line. As shown in Section (ref), the coordinates of $F$ depend on whether $\mathcal{E}$ intersects with $\mathcal{H}_{45}$. When $\Tilde{\mathcal{E}}\equiv\mathcal{E}\cap\mathcal{H}_{45}\neq\emptyset$, as shown in Eqs. (ref)-(ref) we have $F=\mathfrak{u}\cdoth_{\Tilde{\mathcal{E}}}\left(\mathfrak{u}_1\right)$, with $h_{\Tilde{\mathcal{E}}}\left(\mathfrak{u}_1\right)=\inf_{c\in\mathbb{R}}h_\mathcal{E}\left(\mathfrak{u}_1(c)\right)$, $\mathfrak{u}_1(c)\equiv\mathfrak{u}_1-c[1~~-1]^\intercal$, and $\mathfrak{u}\equiv(\mathfrak{u}_1+\mathfrak{u}_2)$. Note that $h_{\Tilde{\mathcal{E}}}\left(\cdot\right)$ is a Hadamard directionally differentiable function of $h_\mathcal{E}(\cdot)$, and its composition with the distance function $\rho$ is again directionally differentiable. Hence, we use our DML estimator, Theorem (ref), and the results in fan:san19 to directly obtain the limit distribution of an estimator for $\rho(e^*, F)$. However, when $\mathcal{E}\subset \mathcal{H}_{45}^+$ (respectively, $\mathcal{E}\subset \mathcal{H}_{45}^-$), $F$ is given by the support set of $\mathcal{E}$ in direction $(\mathfrak{u}_2-\mathfrak{u}_1)$ (respectively, $(\mathfrak{u}_1-\mathfrak{u}_2)$); see Eqs. (ref)-(ref). {As discussed in Remark (ref), carrying out inference if we use DML to directly estimate both coordinates of a support point is difficult,} and therefore we use Eq. (ref) to represent $F$ through moments that involve $h_\mathcal{E}(\cdot)$ only; see Eqs. (ref)-(ref) below.

Observe that $\mathcal{E}\subset \mathcal{H}_{45}^+$ if and only if $F \in \mathcal{H}_{45}^+$, and $\mathcal{E}\subset \mathcal{H}_{45}^-$ if and only if $F \in \mathcal{H}_{45}^-$, as illustrated in Panels (b) and (d) of Figure (ref). Hence, we partition the parameter space $\mathbb{B}_C$, to which $F$ belongs Assumption (ref), into three sets: $\mathbb{B}_C^+\equiv\mathbb{B}_C\cap\mathcal{H}_{45}^+$, $\mathbb{B}_C^-\equiv\mathbb{B}_C\cap\mathcal{H}_{45}^-$, and $\mathbb{B}_C^{45}\equiv\mathbb{B}_C\cap\mathcal{H}_{45}$. We then have the following expressions for $\rho(e^*, F)$:

In each of Eqs. (ref)-(ref), the first condition pins down $\Tilde{e}$ to equal $e^*$; the second and third conditions restrict, respectively, $\tilde{F}$ to be a point on the supporting hyperplane of $\mathcal{E}$ in the appropriate direction and $\tilde{F}\in\mathcal{E}$, which together restricts $\tilde{F}$ to be the support set of $\mathcal{E}$ in this direction. We propose the following testing procedure:

Intuitively, this construction inverts a test that jointly assesses the location of $\mathcal{E}$ relative to $\mathcal{H}_{45}$ and the value of $\rho(e^*,F)$. For example, if $F\in\mathcal{H}_{45}^-$, then both $\mathcal{CS}_n^+(e^*,F)$ and $\mathcal{CS}_n^{45}(\rho(e^*,F))$ are empty with probability approaching one (recall from Section (ref) that $\inf_{c\in\mathbb{R}}h_\mathcal{E}\left(\mathfrak{u}_1(c)\right)$ is unbounded when $F\notin\mathcal{H}^{45}$). Our next result shows that Procedure (ref) delivers an asymptotically valid confidence interval.

As we show in the proof, $\psi^{45}$ is again a composition of Hadamard directionally differentiable functions that define $T_n^{\rho(\Tilde{e}, \Tilde{F})}$ in Eq. (ref). We can therefore employ the same bootstrap method detailed in Procedure (ref) to consistently estimate the quantiles $c_\beta^{45}$ for $\beta\in(0,1)$, where we replace $\phi(\cdot)$ by the composition of the directionally differentiable functions, including $\inf$ and $\rho$, that define $T_n^{45}$, whose exact expression we relegate to Appendix (ref).

Monte Carlo Experiments and Empirical Illustration

Monte Carlo Simulations

We evaluate the finite sample properties of the tests introduced in Sections (ref)-(ref) using two distinct data generating processes (DGPs). For both DGPs, the covariates $X\equiv[X_1,\dots,X_{20}]^\intercal\in\mathbb{R}^{20}$ and group identity $G$ are drawn from the following distributions:

We consider two different ways of generating the outcome $Y$:

where the group-balanced DGP is such that $X$ is informative about $Y$ in opposite directions for group $r$ ($G=1$) and group $b$ ($G=0$), but its predictive power is similar across groups. In contrast, the $r$-skewed DGP is such that $X$ is systematically more informative about the group-$r$ outcome. We take the loss function to be the classification error, $\ell(d,y)=\mathds{1}\{d\neq y\}$. We also construct a status quo algorithm $a^*$, that we fix in the simulations, by training once a logistic regression on a sample of size $10,000$, with $5,000$ observations from the balanced DGP and $5,000$ observations from the $r$-skewed DGP.\footnote{We train the logistic regression on a mixture of group-balanced and $r$-skewed data so that we test for existence of an LDA to the same algorithm $a^*$ in both DGPs. DGP-specific logistic regressions trained on $10,000$ observations drawn form that DGP for each case yield group risks $e^*$ very close to the ones plotted in Figure (ref).} {

}

Figure (ref) depicts the population feasible set $\mathcal{E}$ corresponding to each DGP (pink region), along with a $95\%$ confidence set for the frontier $\mathcal{F}$ (shaded grey region), constructed using a random sample of $10,000$ observations drawn from the respective DGPs, and the risk $e^*$ induced by the status-quo algorithm $a^*$ (green asterisk). Throughout Section (ref), we fit the nuisance parameter $\Delta\boldsymbol{\theta}$ using logit lasso and $5$-fold sample splitting. To assess the finite sample properties of the LDA test in Section (ref) and the distance-to-$F$ test in Section (ref), we test whether $e^*$ is on the frontier and whether it is at a specific distance from $F$.

Table (ref) reports the simulation results. Overall, the Monte Carlo exercise suggests good finite sample properties for our proposed tests, especially as the sample size increases. The top panel corresponds to the test for weak group skew (Section (ref)), where the third column reports the frequency with which the $95\%$ confidence set for the vector $(R, B)$ (Eq. (ref)) based on the available sample fails to cover the true value of $(R,B)$. While the $r$-skewed DGP exhibits some over-rejection at $n=1{,}000$, which is a small sample size relative to the complexity of DML estimation with 20 covariates and estimating the support function across directions, this over-rejection quickly disappears as sample size increases. On the other hand, the weak group skew test (fourth column) rejects the null of weak group skew in the balanced DGP and fails to reject it in the $r$-skewed DGP essentially with probabilities 1 and 0, respectively. This is not surprising because the simulation DGPs are far from the boundary of the null, and because the test is conservative due to the projection step. {

}

The middle panel reports the LDA test results (Section (ref)). The third (fourth) column shows the frequency with which the population point $R$ ($B$), which by definition belongs to $\mathcal{F}$, is rejected by the test of the null that it belongs to $\mathcal{F}$. The fifth (sixth) column reports the frequency with which $(R+B)/2$ (the risk $e^*$ associated with the logit algorithm $a^*$), which by construction does not belong to $\mathcal{F}$, is rejected by the test as an element of $\mathcal{F}$. While at small sample size ($n=1,000$), the test exhibits some over-rejection for the point $B$, this quickly disappears as sample size grows. At small and medium sample sizes ($n=1,000$ and $n=5,000$), the rejection probability for the false null that $(R+B)/2\in\mathcal{F}$ is low for the $r$-skewed DGP, while it is high at all sample sizes for the balanced DGP. This is justifiable in light of Figure (ref), which shows that in the $r$-skewed DGP the chord between $R$ and $B$ is close to the frontier and lies largely inside its 95% confidence set. For $n=10,000$, the test detects the false nulls in columns five and six with substantially higher probability.

The bottom panel reports the results for the distance-to-$F$ test (Section (ref)). This test is more delicate because its implementation requires solving an optimization problem to estimate $F$. Here we observe more substantial over-rejection for $n=1,000$, but the distortion gets markedly reduced as sample size increases.

Empirical Illustration

We revisit the analysis in \citet*[OPVM henceforth]{obe:pow:vog:mul19}, who analyze properties of the algorithm used by a research hospital to determine if a patient should be automatically enrolled in a high-risk care management program. The algorithm used by the research hospital aims at predicting patients' health needs based on total medical expenditures (the label on which the algorithm is trained). It produces a health risk score and automatically enrolls a patient in the high-risk care management program if that patient's risk score exceeds the $97^{\text{th}}$ percentile of all predicted scores. We reassess this algorithm through our testing procedures. To do so, we use the synthetic data made available by \citet*{li:lin:obe19} at \href{https://gitlab.com/labsysmed/dissecting-bias}{GitLab} to replicate all analyses in \citetalias{obe:pow:vog:mul19}. The data include $48,784$ patient observations, of which $5,582$ self report as Black and the others self report as White ($g\in\{\textsf{bl},\textsf{wh}\}$). The data include $149$ covariates, such as age, gender, comorbidity and medication variables, costs, and biomarkers, and provide information about each patient's number of active chronic conditions in the subsequent year, viewed as the true measure of health needs of a patient. Following \citetalias{lia:lu:mu:oku24}, we let $\ell(d,Y)=\mathds{1}\{Y\neq d\}$ be the classification loss and, unless explicitly stated otherwise, $Y_i=\mathds{1}\{$patient $i$ has 6 or more chronic conditions$\}$, with the choice of 6 driven by the fact that it is the $97^{\text{th}}$ percentile of active chronic condition numbers across patients in the sample. We use random forests with $5,000$ trees to estimate $\Delta\boldsymbol{\theta}$ as described in Section (ref).

Feasible Set Estimation and Inference for the Frontier

We report in Figure (ref) an estimate of the feasible set $\mathcal{E}$ based on Eq. (ref) using $1,000$ directions (top-left panel), along with 100 estimated supporting hyperplanes (top-right panel), where across all panels the horizontal (vertical) axis is the risk for Blacks (Whites), denoted as $e_{\textsf{bl}}$ ($e_{\textsf{wh}}$). We zoom in to show the estimated FA-frontier $\widehat\mathcal{F}$ and its $95\%$ confidence set (bottom-left panel), and zoom in further to show the best risk achievable for the Black patients (bottom-right panel, red point labeled $\textsf{BL}$), which coincides with and is overlaid by the best risk for the White patients (bottom-right panel, blue point labeled $\textsf{WH}$).

Hypothesis Testing

We next report the weak group skew test results in Figure (ref)-Panel (a), where both $\textsf{BL}$ and $\textsf{WH}$ are below the $45^\circ$ degree line and the feasible set is $\textsf{wh}$-skewed; the test fails to reject the null of weak group skew, suggesting that implementing Pareto-dominated algorithms could be justified based on the designer's preferences over fairness and accuracy. In Figure (ref)-Panel (b) we plot $\widehat{\mathcal{F}}$, along with its $95\%$ confidence set. Figure (ref)-Panel (b) also plots the estimated group risks associated with the original algorithm used by the hospital (a green asterisk) and three alternative algorithms that \citetalias{obe:pow:vog:mul19} experiment with to assess whether different label choices yield decision rules that are more accurate and fairer than the algorithm currently used by the hospital. One of these algorithms is trained to predict total cost (hollow diamond with a cross in Figure (ref)-Panel (b)), one to predict avoidable costs (filled diamond), and the other one to predict the number of active chronic conditions (hollow diamond). All our tests take into account the finite sample estimation error of the group risks induced by these four algorithms. The top panel of Table (ref) reports the values of the LDA test statistics and the associated critical values that we compute for the four algorithms considered by \citetalias{obe:pow:vog:mul19}. The hypothesis that the original algorithm yields group risks on the frontier is rejected, and so is the same hypothesis for group risks associated with the algorithm trained to predict total costs. On the other hand, we fail to reject that the algorithms trained to predict avoidable costs and the number of active chronic conditions yield group risks on the FA-frontier. Regarding the three algorithms proposed by \citetalias{obe:pow:vog:mul19}, if one were to plot the set $\mathcal{C}^*$ corresponding to the original algorithm, it would be immediate to see that one cannot reject the hypothesis that all other algorithms considered by \citetalias{obe:pow:vog:mul19} improve upon it, both in terms of fairness and accuracy.

The bottom panel of Table (ref) shows that the closest risk to $F$ is the one associated with the algorithm predicting total costs. However, all confidence intervals overlap.

{

}

Performance of Algorithms Constructed to Be on a Constrained Frontier

Finally, we compare various outcomes associated with the algorithms considered by \citetalias{obe:pow:vog:mul19} with those of decision rules resulting from the algorithms that we propose in Section (ref). To do so, we randomly split the sample into two halves, and use one half (training data) to estimate the nuisance parameter $\Delta\boldsymbol{\theta}$, and implement Eq. (ref) on the other half (evaluation data) for the following choices of $q$:

Throughout, we recognize that to compare algorithms we should enforce a global capacity constraint on the total percentage of patients assigned to the high-risk case management program, so that variation across algorithms is not confounded with a possibly more generous care program. To enforce the capacity constraint, we need to require the algorithms to satisfy $\int a(x)d\mathbb{P}_X\le \bar{a}$ for some known constant $\bar{a}$. For example, $\bar{a}=0.03$ when comparing with the algorithm currently used by the hospital which assigns only patients with risk score above the $97^\text{th}$ percentile to the high-risk care program. Recall from the discussion following Proposition (ref) that without capacity constraints, we have

When $\mathcal{A}(\mathcal{X})$ is constrained to only include algorithms such that $\int a(x)d\mathbb{P}_X\le \bar{a}$, maximization in Eq. (ref) is achieved by setting

for $\mathrm{quant}_{k(\boldsymbol{\theta}(X),\cMq)}(\alpha)$ the $\alpha$-quantile of $k(\boldsymbol{\theta}(X),\cMq)$ and $k(\boldsymbol{\theta}(X),\cMq)$ as in Eq. (ref).\footnote{In this case, the constrained feasible set is $\mathcal{E}^\texttt{co}_{\bar{a}}\equiv \left\{\bigl(e_r(a),e_b(a)\bigr)\in\mathbb{R}^2:a\in\mathcal{A}(\mathcal{X})~\text{and}~\int a(x)d\mathbb{P}_X\le \bar{a}\right\}$, a convex subset of $\mathcal{E}$, and the algorithms in Eq. (ref) is on its frontier.} In words, for our example with $\bar{a}=0.03$, high-risk care management is assigned to those patients with positive values of $k(\boldsymbol{\theta}(X),\cMq)$ that exceed its $97^\text{th}$ percentile.

The results of our first exercise are reported in Figure (ref). The top panel replicates Figure 1-(a) in \citetalias{obe:pow:vog:mul19}, and shows that at each percentile of the original algorithm's risk score (the horizontal axis), Black patients have a substantially larger number of active chronic conditions (the vertical axis) than White patients.\footnote{The top panel of Figure (ref) is not identical to Figure 1-(a) of \citetalias{obe:pow:vog:mul19} because it is plotted using the synthetic data, instead of the real data, and because the code used to generate Figure 1-(a) of \citetalias{obe:pow:vog:mul19} had an error that was later corrected after the publication of \citetalias{obe:pow:vog:mul19}; see the documentation of li:lin:obe19 for more details.} Even among patients automatically enrolled in the high-risk care management program (those with a risk score above the $97^{\text{th}}$ percentile), Black patients appear to be in worse health than White patients.

We ask whether the algorithms that we propose are able to select for treatment patients that, in fact, exhibit substantially worse health outcomes one year later. The bottom four panels in Figure (ref) plot, for each of the algorithms described at the beginning of this section, the number of active chronic conditions for (a) Black patients that our algorithm does not assign to the high-risk care program (dash-dotted, light-purple line); (b) White patients that our algorithm does not assign to the high-risk care program (solid salmon line); (c) Black patients that our algorithm assigns to the high-risk care program (dash-dotted, dark-purple line); (b) White patients that our algorithm assigns to the high-risk care program (solid orange line). For comparability with the top panel, on the horizontal axis we continue to report the risk score produced by the algorithm used by the hospital.

The main takeaway from Figure (ref) is that our four proposed algorithms, which use group identity to estimate $\Delta\boldsymbol{\theta}$ but not for treatment assignment, are successful at selecting for treatment patients who, one year later, experience substantially worse health outcomes. Moreover, with the Ralwsian, Majority, and Utilitarian algorithms, Black and White patients assigned to treatment have similar numbers of chronic conditions. Only the Egalitarian algorithm shows notable disparities for patients with hospital risk scores below the $70^{\text{th}}$ percentile. We think this result might be due to two factors: estimation of $F$ is challenging and hence the direction $\widehat{q}(\widehat{F})$ might be imprecisely estimated; and the feasible set is $\textsf{wh}$-skewed and hence the Egalitarian algorithm leads to a Pareto-dominated outcome.

Our last exercise aims at further assessing the extent to which our proposed algorithms may reduce the substantial disparities between Black and White patients in the current program screening practices documented by \citetalias{obe:pow:vog:mul19}. To illustrate the potential for improvement over the hospital's algorithm, \citetalias{obe:pow:vog:mul19} simulate “a counterfactual world with no gap in health conditional on risk” (p.3). They construct an infeasible, couterfactual algorithm that uses group identity and patients' ex-post active chronic health conditions to find the sickest Black patient with health risk score just below a threshold (the “inframarginal Black patient”) and the healthiest White patient with health risk score just above the same threshold (the “supramarginal White patient”). If the number of chronic health conditions of the inframarginal Black patient is larger than that of the supramarginal White patient, they iteratively swap them until the number of chronic health conditions of the inframarginal Black patient equals that of the supramarginal White patient. \citetalias{obe:pow:vog:mul19} find that at all risk thresholds $\alpha$ above the median, the counterfactual algorithm increases the fraction of Black patients treated. Table (ref), columns 2-3, show that across the various thresholds for automatic enrollment in the program that we consider, the fraction of Black patients would rise between 6 and 41 percentage points.\footnote{In columns 2-3, the fractions reported are based on a denominator that equals the total number of patients with risk scores above a certain percentile of risk scores, e.g., the $55^{\text{th}}$ percentile, with this threshold viewed as the threshold above which the patient is automatically enrolled in the high-risk care program, and a numerator that equals the number of Black patients with risk score above that percentile.}

We ask how much of this gap could be filled if one were to use our algorithms, which are feasible and do not rely on ex-post knowledge of the number of active chronic conditions nor on using group identity for assignment to the high-risk care program. The results, reported in Table (ref), columns 4-7,\footnote{For each algorithm, the fraction of Black patients treated at each capacity threshold is computed as the following ratio. The denominator equals the total number of patients treated under a given algorithm, e.g., the Rawlsian algorithm, subject to the constraint that at each threshold $\bar{a}\in[0.55,0.97]$ at most $100(1-\bar{a})\%$ of the evaluation sample is treated, and for each threshold $\bar{a}$ the outcome $Y$ equals the indicator of whether the number of active chronic conditions exceeds the $\bar{a}$-quantile of the distribution of the number of active chronic conditions.} show that each of our four algorithm yields an increase in the fraction of Black patients treated at all capacity thresholds. The Rawlsian and the Utilitarian algorithms yield the largest increases, ranging between 4 and 18 percentage points across the various capacity thresholds. This shows that these two feasible, easy to implement algorithms can close between $26\%-69\%$ of the gap between the algorithm that the hospital uses and the counterfactual, infeasible algorithm simulated by \citetalias{obe:pow:vog:mul19}.

{

}

We conclude by noting that results based on random forests trained using the grf package are not guaranteed to reproduce exactly across platforms, even with the same seed (all simulations and estimation are in R). This is a known feature of grf (see the \href{https://grf-labs.github.io/grf/REFERENCE.html\#forests-predict-different-values-depending-on-the-platform-even-though-the-seed-is-the-same}{reference manual}). Similarly, tests based on optimization via stochastic gradient descent implemented by the torch package are not exactly reproducible, even after seed setting, which affects the $F$ estimator and thus the distance-to-$F$ test (see the \href{https://github.com/mlverse/torch/issues/1311}{repository discussion}). To assess sensitivity of our results to these features, we repeat the entire empirical exercise 20 times with different seeds. Appendix (ref) reports results and includes a robustness check using logit lasso for estimating the nuisance function.\footnote{For the simulations in Section (ref), we expect cross-platform variability to be negligible, as results are averaged over 1000 Monte Carlo replications.} While the results exhibit some nontrivial variation across seeds (although the qualitative results are unchanged), we view this as expected, considering that the nonparametric estimation step involves 149 covariates against a total sample size of $48,784$ and a minority group of size $5,582$.

Conclusion

We provide a consistent nonparametric estimator for a theoretical fairness-accuracy frontier proposed by \citetalias{lia:lu:mu:oku24} and algorithms that attain points on this frontier. We obtain the estimator through judicious use of the separating hyperplane theorem and the support function of the (convex) feasible set of expected losses associated with all possible algorithms, a portion of whose boundary coincides with the FA-frontier. We provide a DML estimator of the support function and show it converges to a tight Gaussian process as sample size increases. We formulate important policy-relevant hypotheses that have received much attention in the fairness literature as restrictions on the support function and construct valid test statistics. We provide an estimator for the distance between a given algorithm and the fairest point on the frontier. We carry out a Monte Carlo exercise that illustrates the good finite sample properties of our method, and demonstrate its practical relevance by revisiting the empirical analysis in \citetalias{obe:pow:vog:mul19}. Our results show that the algorithm that a research hospital employs to screen patients for high-risk care is not on the frontier. Our proposed algorithms substantially improve over this status quo in terms of both fairness and accuracy.