Fair Supervised Learning with A Simple Random Sampler of Sensitive Attributes

Jinwon Sohn Qifan Song Guang Lin

Purdue University Purdue University Purdue University

###### Abstract

As the data-driven decision process becomes dominating for industrial
applications, fairness-aware machine learning arouses great attention in
various areas. This work proposes fairness penalties learned by neural
networks with a simple random sampler of sensitive attributes for non-
discriminatory supervised learning. In contrast to many existing works that
critically rely on the discreteness of sensitive attributes and response
variables, the proposed penalty is able to handle versatile formats of the
sensitive attributes, so it is more extensively applicable in practice than
many existing algorithms. This penalty enables us to build a computationally
efficient group-level in-processing fairness-aware training framework.
Empirical evidence shows that our framework enjoys better utility and fairness
measures on popular benchmark data sets than competing methods. We also
theoretically characterize estimation errors and loss of utility of the
proposed neural-penalized risk minimization problem.

##  1 Introduction

Algorithmic fairness has been a growing research area as the prediction-based
decision process becomes more and more prevalent. The legal examples include
the US Equal Credit Opportunity Act, the European Union’s General Data
Protection Regulation, and the Fair Credit Reporting Act, to name a few. On
the academic side, Mhasawade et al., (2021) discussed the importance and
challenges of algorithmic fairness in public health. Kozodoi et al., (2022)
studied various mitigation strategies for the credit scoring application. De-
Arteaga et al., (2022) gave an overview of potential areas demanding fairness-
aware business analytics with intriguing real-world examples such as dynamic
pricing, distribution of vaccines, job applications, and so forth. Similar
motivating examples can be easily found in other fields as well, such as
education (Loukina et al.,, 2019), finance (Das et al.,, 2021), mortgage
lending (Lee and Floridi,, 2021), computational medicine (Xu et al.,, 2022).

In pursuit of the rising demand for mitigating societal bias in decision-
making processes, there have been versatile approaches, which are mainly
categorized into pre-processing, in-processing, and post-processing (Barocas
et al.,, 2017; Caton and Haas,, 2020; Pagano et al.,, 2023; Xian et al.,,
2023). Pre-processing includes relabeling, reweighting, or resampling of data
instances to ease possible discrimination of a learned model (Kamiran and
Calders,, 2012). Synthesizing data to be fair also belongs to this class as
well (Xu et al.,, 2018; Sattigeri et al.,, 2019; van Breugel et al.,, 2021).
In the field of in-processing, optimization with fairness constraints has been
one of the main branches (Agarwal et al.,, 2018; Komiyama et al.,, 2018;
Agarwal et al.,, 2019; Zafar et al.,, 2019; Scutari et al.,, 2022; Jung et
al.,, 2023). Comparably, model-based fairness control has been also
extensively studied, which regularizes the degree of discrimination through an
auxiliary model (Xie et al.,, 2017; Beutel et al.,, 2017; Zhang et al.,, 2018;
Adel et al.,, 2019; Mary et al.,, 2019; Lee et al.,, 2022). Following the
model-agnostic spirit, post-processing modifies the model training outcome
such that the reporting values are non-discriminatory while minimizing the
loss of utility (Hardt et al.,, 2016; Pleiss et al.,, 2017; Zeng et al.,,
2022). Out of such mitigation categories, there have been versatile studies
imposing fairness in clustering (Wang et al.,, 2023), feature selection
(Quinzan et al.,, 2022), reinforcement learning (Deng et al.,, 2022), to name
a few.

There are three key notions for pursuing better algorithmic fairness:
independence, separation, and sufficiency (Barocas et al.,, 2017; Caton and
Haas,, 2020; Pagano et al.,, 2023). Throughout this paper, we denote $X$ as
the covariate(s), $Y$ as a response, $A$ as the sensitive random variable(s)
to be protected (e.g., race or gender), and $h$ as the scoring model that
makes a data-driven decision. Let $P(\cdot)$ denote the
(joint/marginal/conditional) distribution based on the context. The
independence condition requires that $A$ is independent of $h(X)$, i.e.,
$A\perp h(X)$; statistical parity is an example of independence. Separation
and sufficiency, on the other hand, are defined based on the conditional
independence structure, i.e., $A\perp h(X)|Y$ and $A\perp Y|h(X)$
respectively. For instance, if two groups in a binary $A$ have the same false-
positive and false-negative errors in a binary classification problem, the
circumstance achieves the balanced error rate (Chouldechova,, 2017) or
equalized odds (Hardt et al.,, 2016), implying separation. For sufficiency,
evaluating the calibration of the model across sensitive factors is a
representative example (Chouldechova,, 2017; Barocas et al.,, 2017).

In this work, we propose a model-based in-processing framework that flexibly
captures and controls underlying societal discrimination realized through the
output of $h$. Our contributions are as follows:

First, we devise a fairness penalty by a neural network denoted as $D$ that
leverages a simple random sampler of sensitive variables for independence to
measure the discrepancy between $P(h(X),A)$ and $P(h(X))P(A)$. The penalty can
accommodate various scenarios even when $A$ is a mix of continuous and
discrete variables, in which many prior works fail as exhibited in Table 1.

Second, the penalty is further extended to embrace separation while inheriting
all the remarkable properties shown in the independence case.

Third, we mathematically quantify the underlying mechanism of the proposed
fairness-controlled supervised learning. Specifically, we derive the upper
bounds of estimation error and loss of utility with the imposed penalty on the
basis of statistical learning theory.

##  2 Related Works

##### Density Matching

Achieving fairness through density matching has been steadily implemented.
Quadrianto and Sharmanska, (2017) designed a privileged learning using a
variant of a support vector machine on which a fairness constraint via maximum
mean discrepancy is imposed. Cho et al., (2020) performed matching
distributions between different sensitive groups via a differentiable kernel
density estimation technique. Li et al., (2021) adopted the loss function used
in the generative adversarial network (Goodfellow et al.,, 2014) as an extra
penalty to facilitate the scoring distributions $h$ conditional on different
$A$ values to be indistinguishable. These density-matching-based
regularizations work well but are usually limited to the case when $A$ and $Y$
are discrete or binary. That is because their backbone optimization structures
essentially depend on _sub-groups comparison_. As an example of independence
for binary $A$, the loss function relies on two subsets of data, corresponding
to $A=1$ and $A=0$ respectively, to measure the distributional difference
between sub-groups (i.e., $P(h(X)|A=1)$ vs $P(h(X)|A=0)$). Obviously, this
technique is not applicable to continuous $A$. Romano et al., (2020) suggested
a different training scheme, which uses an extra generator that generates $A$
conditionally on $Y$ to ensure separation. Although this approach technically
bypasses the sub-group comparison, critically relies on the quality of the
conditional generator. If $A$ is complex, e.g., a mix of categorical and
continuous variables, finding a good generator that captures the true
distribution could be an extremely difficult task (see, e.g., Xu et al.,
2019b, ; Kotelnikov et al.,, 2023).

##### HGR penalty and beyond

To tackle such limitations in the literature, diverse approaches have been
proposed. Mary et al., (2019) first utilized Hirschfeld-Gebelein-Rényi (HGR)
maximal correlation which can capture the nonlinear correlation between $A$
and $h(X)$ regardless of the variable type. Grari et al., (2019) and Lee et
al., (2022) approximated the HGR maximal correlation via neural networks, to
achieve independence and separation. These approaches conduct min-max
optimization where the auxiliary neural networks approximate the HGR by
maximization and the predictive model $h(\cdot)$ is trained to minimize the
approximated HGR. The auxiliary neural nets in the training process take
variable $A$ as the input instead of fitting $A$, hence the multivariate
nature of $A$ does not matter anymore. Besides, Du et al., (2021) recently
proposed a way of neutralizing the hidden layer of a neural-net model, which
can also handle multivariate $A$; Scutari et al., (2022) developed a fair
(generalized) linear model that places the ridge penalty to adjust the
violation of discrimination for any kinds of $A$ and $Y$.

Table 1: Summary of applicability of in-processing methods: “NN” denotes whether or not neural network models are supported. “Ind.” and “Sep.” imply independence and separation. Abbreviations “D.”, “C.”, and “M.” imply discrete, continuous, and mixed (continuous and discrete) type; Abbreviations “$A$” and “$Y$” mean sensitive and outcome variables. “$\triangle$” means the framework does not explicitly control specific fairness metrics. “$\square$” means the method is technically applicable but requires a non-trivial modification. “B.” implies the framework only supports a binary type.  Methods | NN | Ind. | Sep. | D. $A$ | C. $A$ | M. $A$ | D. $Y$ | C. $Y$
---|---|---|---|---|---|---|---|---
Edwards and Storkey, (2015) | ✓ | ✓ | ✗ | B | ✗ | ✗ | ✓ | ✗
Agarwal et al., (2018) | ✓ | ✓ | ✓ | ✓ | ✗ | ✗ | ✓ | ✗
Adel et al., (2019) | ✓ | ✓ | ✓ | B | ✗ | ✗ | ✓ | ✗
Cho et al., (2020) | ✓ | ✓ | ✓ | ✓ | ✗ | ✗ | ✓ | ✗
Romano et al., (2020) | ✓ | ✗ | ✓ | ✓ | ✓ | $\square$ | ✓ | ✓
Li et al., (2021) | ✓ | ✓ | ✓ | B | ✗ | ✗ | ✓ | ✗
Du et al., (2021) | ✓ | $\triangle$ | $\triangle$ | ✓ | ✓ | ✓ | ✓ | ✗
Scutari et al., (2022) | ✗ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓
Lee et al., (2022) | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓
Ours | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓ | ✓

##  3 Identifying Discrimination by $D$

###  3.1 Penalty with Resampled $A$

Focusing on independence, this section describes the main idea which
generalizes to any dimension or type of sensitive variables. The
methodological novelty originates from the use of a simple random sampler of
$A$ to quantify the degree of fairness via a discriminative neural network
$D$. Given a sufficient network capacity, the learned $D$ can capture the
violation of statistical independence, so it is used as a penalty (or
adversarial) network for the risk-minimization problem of $h$ such that the
resulting $h$ becomes fair.

To start with, let’s formally define $X\in{\cal X}\subset\mathbb{R}^{p}$ and
$A\in{\cal A}\subset\mathbb{R}^{l}$ as the multivariate non-sensitive and
sensitive random variable respectively, and $Y\in{\cal Y}\subset\mathbb{R}$ as
a univariate outcome variable. Given a scoring model $h:{\cal
X}\rightarrow{\cal S}\subset\mathbb{R}$ and a loss function $L:{\cal
Y}\times{\cal S}\rightarrow{\mathbb{R}}$, we denote the risk function as
$R(h)={\bf E}_{Y,X}[L(Y,h(X))]$ where ${\bf E}_{Y,X}$ means the expectation
w.r.t. the joint distribution of ($X$,$Y$). Note $h$ does not explicitly
depend on $A$. This structure conceptually carries on fairness via blinding
and does not input $A$ in the inference phase (Quadrianto and Sharmanska,,
2017; Zafar et al.,, 2019). In order to evaluate independence in general
settings (i.e., beyond binary $A$ and $Y$), we generalize the existing notion
of statistical parity as follows.

######  Definition 1 (Generalized Statistical Parity (GSP)).

Let $X$ and $A$ be the non-sensitive and sensitive random variables. The
scoring function $h:{\cal X}\rightarrow{\cal S}$ achieves the generalized
statistical parity if $P(h(X)|A=a)\overset{d}{=}P(h(X))$ for all $a\in{\cal
A}$.

GSP is more general than the classical statistical parity (SP); if ${\cal
A}=\\{0,1\\}$ for the binary classification problem, GSP implies SP, i.e.,
$P(h(X)>\tau|A=1){=}P(h(X)>\tau|A=0)$ for any classification threshold $\tau$.

Now, we define the penalty that promotes fairness of $h$. Inspired by the
noise-contrastive loss (Gutmann and Hyvärinen,, 2012; Goodfellow et al.,,
2014), with $D:{\cal S}\times{\cal A}\rightarrow(0,1)$, let’s define
$R_{F}(h):=\sup_{D}~{}R_{F}(h;D)$, where

$\displaystyle R_{F}(h;D)$ $\displaystyle={\bf E}_{X,A}[\log D(h(X),A)]$
$\displaystyle+{\bf E}_{X,A^{\prime}}[\log(1-D(h(X),A^{\prime}))],$

where $A^{\prime}$ shares the same distribution with $A$ but is independent to
$(X,A)$. We remark that $A^{\prime}$ can be easily obtained by applying the
simple random sampling $A^{\prime}$ from $P(A)$. Then the following
proposition theoretically validates that $R_{F}(h)$ captures the discrepancy
between $P(h(X)|A)$ and $P(h(X))$.

######  Proposition 1.

Let $p_{A}$, $p_{h(X)}$, $p_{h(X)|A}$, and $p_{h(X),A}$ be the marginal
densities of $A$ and $h(X)$, the conditional density of $h(X)$ given $A$
respectively, and the joint density of $h(X)$ and $A$. Denote
$D^{*}=\underset{D}{\arg}\max~{}R_{F}(h;D)$. Then, for all $s\in{\cal S}$ and
$a\in{\cal A}$,

$\displaystyle\dfrac{D^{*}(s,a)}{1-D^{*}(s,a)}=\dfrac{p_{h(X),A}(s,a)}{p_{h(X)}(s)p_{A}(a)}=\dfrac{p_{h(X)|A}(s|a)}{p_{h(X)}(s)},$

The proof is shown in Supplementary 7.1. This proposition provides the
theoretical justifications for the use of $R_{F}(h;D)$ as a GSP controller.
Following the argument in Theorem 1 of Goodfellow et al., (2014), we observe
that $R_{F}(h)$ can be interpreted by the Jensen-Shannon divergence
$J(\cdot,\cdot)$, i.e., $R_{F}(h;D^{*})~{}=~{}2J(P(h(X),A),P(h(X))P(A))-2\log
2$, and it implies $p_{h(X),A}(s,a)=p_{h(X)}(s)p_{A}(a)$ for all $s$ and $a$
if $J=0$, which underpins that $h$ accomplishes GSP at the minimum of
$R_{F}(h)$. Thus, we can formulate a fairness-aware optimization for $h$ by
placing the extra penalty $R_{F}(h)$ to the (discriminatory) risk-minimization
problem, i.e.,

$\displaystyle\underset{h}{\min}~{}R(h)+\lambda R_{F}(h),$ (1)

where $\lambda$ trades off between goodness-of-fit (utility) and the degree of
GSP; as $\lambda$ becomes larger, the solution model $h$ becomes fairer (in
terms of GSP) where $D^{*}$ gets closer to 0.5 but meanwhile undergoes the
loss of utility.

Finally, we employ an empirical min-max optimization structure since (1) is
practically intractable. The population joint density of $Y$, $X$, and $A$ is
not available, and the optimal $D^{*}$ is unknown in general and requires
maximization over $R_{F}(h;D)$. Let’s denote by
$\\{(X_{i},A_{i},Y_{i})\\}_{i=1}^{n}$ the observed data set with sample size
$n$, $\\{A_{i}^{\prime}\\}_{i=1}^{n}$ the set of resampled
$\\{A_{i}\\}_{i=1}^{n}$ via the simple random sampler, and further denote by
$\hat{R}(h)=\frac{1}{n}\sum_{i=1}^{n}L(Y_{i},h(X_{i}))$ and
$\hat{R}_{F}(h;D)=\frac{1}{n}\sum_{i=1}^{n}\left(\log
D(h(X_{i}),A_{i})+\log(1-D(h(X_{i}),A^{\prime}_{i}))\right)$ the unbiased
estimators of $R$ and $R_{F}$ respectively. Then, the empirical version of (1)
is

$\displaystyle\min_{h}\max_{D}~{}\hat{R}(h)+\lambda\hat{R}_{F}(h;D),$ (2)

where we use neural networks to model both $h$ and $D$. The optimization
problem (2) can be solved via an alternative min-max strategy (Supplementary
7.2). For every iteration, $D$ is first trained to capture the degree of
discrimination against GSP. The model $h$ is then trained to minimize the risk
plus the fairness penalty evaluated by $D$, which in turn hinders $D$ from
identifying the presence of discrimination in the next iteration. By repeating
this adversarial game between the two networks $D$ and $h$, the model $h$
settles down to an equilibrium, determined by $\lambda$, between utility and
fairness. We name $\hat{R}_{F}(h;D)$ as _simple random sampling (SRS)-based
penalty (SBP)_.

This penalty is not only adequately flexible to capture nonlinear
discriminatory dependency between $h(X)$ and $A$ but also differentiable while
most fairness criteria are usually non-differentiable in general (Zafar et
al.,, 2019; Cotter et al.,, 2019). Noteworthily, _this penalty requires
neither partitioning data sets into sub-groups nor predicting sensitive
attributes_ in contrast to a lot of previous works (Edwards and Storkey,,
2015; Xie et al.,, 2017; Beutel et al.,, 2017; Zhang et al.,, 2018; Adel et
al.,, 2019; Zhao et al.,, 2019; Li et al.,, 2021; Cho et al.,, 2020; Du et
al.,, 2021). These properties greatly expand the applicability of the proposed
method to numerous real-world problems having a mix of continuous and discrete
sensitive attributes. We also remark that the previous sampler-based work of
Romano et al., (2020) is not applicable to GSP.

###  3.2 Extension for Separation

We also devise a penalty for separation that inherits all the advantages
discussed in Section 3.1. To begin with, we define the generalized equalized
odds (GEO) which is an intuitive generalization of equalized odds (Hardt et
al.,, 2016) beyond binary outcomes.

######  Definition 2 (Generalized Equalized Odds (GEO)).

The model $h$ is called to satisfy the generalized equalized odds if

$\displaystyle p(h(X)|A=a,Y=y)=p(h(X)|Y=y),$

for all $a\in{\cal A}$ and $y\in{\cal Y}$.

To control GEO, we suggest specifying $R_{F}(h;D)=$

$\displaystyle{\bf E}_{X,A,Y}[\log D(h(X),A,Y)]$ $\displaystyle+{\bf
E}_{A^{\prime}}{\bf
E}_{X,Y}[\beta(A^{\prime},Y)\log(1-D(h(X),A^{\prime},Y))],$

with some function
$\beta:\mathcal{A}\times\mathcal{Y}\rightarrow\mathbb{R}^{+}$ and $D:{\cal
S}\times{\cal A}\times{\cal Y}\rightarrow(0,1)$, where $A^{\prime}$ is
statistically independent to $(X,A,Y)$ but $A=A^{\prime}$ in distribution, and
the penalty accompanies the next proposition.

######  Proposition 2.

Let $p_{h(X)|Y}$ be the conditional density function of $h(X)$ given $Y$ and
$p_{h(X)|A,Y}$ be of given $Y$ and $A$. For $R_{F}(h;D)$, if
$D^{*}=\underset{D}{\arg}\max~{}R_{F}(h;D)$, then $D^{*}(s,a,y;\beta)=$

$\displaystyle\dfrac{p_{h(X)|A,Y}(s|a,y)}{p_{h(X)|A,Y}(s|a,y)+\beta(a,y)p_{h(X)|Y}(s|y)\frac{p_{A^{\prime},Y}(a,y)}{p_{A,Y}(a,y)}}$

for all $s\in{\cal S}$, $a\in{\cal A}$, and $y\in{\cal Y}$, where
$p_{A^{\prime},Y}$ and $p_{A,Y}$ be the joint density functions of
$A^{\prime}$ and $Y$ and of $A$ and $Y$ respectively.

The proof appears in Supplementary 8.1. If $\beta$ cancels out the density
ratio between $p_{A^{\prime},Y}$ and $p_{A,Y}$ (i.e.,
$\beta(a,y)p_{A^{\prime},Y}(a,y)/p_{A,Y}(a,y)=1$), then by the same argument
in Proposition 1, $R_{F}(h;D^{*})$ is minimized when $h$ satisfies
$P(h(X)|A,Y)\overset{d}{=}P(h(X)|Y)$, i.e., the penalty $R_{F}(h;D^{*})$
promotes GEO. Thus, we fit a density-ratio estimator $\hat{\beta}$ by
maximizing $R_{\beta}(D_{\beta})=$

$\displaystyle{\bf E}_{A,Y}[\log D_{\beta}(A,Y)]+{\bf
E}_{A^{\prime},Y}[\log(1-D_{\beta}(A^{\prime},Y))],$

where $D_{\beta}:{\cal A}\times{\cal Y}\rightarrow(0,1)$ is modeled by a
neural network. Then
$D_{\beta}^{*}=\arg_{D_{\beta}}\max~{}R_{\beta}(D_{\beta})$ satisfies
$D_{\beta}^{*}(a,y)/(1-D_{\beta}^{*}(a,y))=p_{A,Y}(a,y)/p_{A^{\prime},Y}(a,y)$.
Given an empirical minimizer $\hat{D}_{\beta}$ of $R_{\beta}(D_{\beta})$, let
$\hat{\beta}(a,y)=\hat{D}_{\beta}(a,y)/(1-\hat{D}_{\beta}(a,y))$, then by
Proposition 2, $D^{*}(s,a,y;\hat{\beta})$ properly regularizes GEO. The
numerical algorithm to expedite GEO appears in Supplementary 8.2 on the basis
of the expression (2). The algorithm needs a pre-training phase to estimate
$\hat{\beta}$ that is leveraged as the adaptive weights for the evaluation of
the penalty $R_{F}(h;D)$. We observe that $\hat{\beta}$ performs a powerful
density-ratio estimation on toy examples in Supplementary 8.3. It is worth
mentioning that when $A$ and $Y$ are both discrete, $\hat{\beta}$ can be
easily found without neural network training, e.g., using the empirical
probability mass functions.

Proposition 2 clearly shows the fundamental differences of our simple random
sampler approach compared to the conditional sampler approach (Romano et al.,,
2020), i.e., $A^{\prime}\sim P(A|Y)$. In our case, we additionally employ
$\beta$ in exchange for sampling $A^{\prime}$ directly from $P(A)$. In
contrast, the prior work does not need $\beta$ but has to obtain the
conditional sampler $A^{\prime}=G(\epsilon,Y)$ for some generative model $G$
and random noise $\epsilon$. There may be room for debate but we believe that
estimating $\beta$ is generally more straightforward than training $G$. The
former is equivalent to a DNN binary classification problem, which in modern
practice, is hardly impacted by the dimension or characteristics of $A$ for
tabular data sets. In contrast, the latter may necessitate the use of advanced
generative models, such as generative adversarial networks (GAN) or diffusion
models which are well-known for their notorious training difficulty (Xu et
al., 2019b, ; Kotelnikov et al.,, 2023), when dealing with tabular-type
attributes of $A$. Moreover, even a successfully trained generative model may
only learn the support rather than the shape of a distribution.

######  Remark 1.

Although our approach necessitates $\beta$ for the theoretical justification,
we find that the performance of our simple random sampler-based method is
fairly robust to a poor $\beta$ estimation. Refer to Section 5 and
Supplementary 10.4 for more details.

Finally, it is worth mentioning that this penalty can also be used to obtain
fair representation (Zhao et al.,, 2019; Du et al.,, 2021). Let’s denote by
$E:{\cal X}\rightarrow{\cal E}$ the encoder which precedes $h_{E}:{\cal
E}\rightarrow{\cal S}$ with the risk $\hat{R}(h_{E}\circ E)$. Then solving
$\min_{h_{E},E}\max_{D}\hat{R}(h_{E}\circ E)+\lambda\hat{R}_{F}(E;D)$ brings a
fair encoder $E$ for either independence or separation with $D:{\cal
E}\times{\cal A}\rightarrow(0,1)$ or $D:{\cal E}\times{\cal A}\times{\cal
Y}\rightarrow(0,1)$ respectively.

##  4 Theory

In this section, we characterize the estimation error and the loss of utility
of the proposed fairness-aware optimization scheme. For simplicity of the
presentation, the analysis focuses on the solution of (2) with the GSP
penalty. Our analysis borrows some proof techniques in the literature of
generative adversarial modeling (Ji et al.,, 2021). For readers who are
interested in the details of the proof, please refer to Supplementary 9.

Let’s suppose ${\cal X}=\\{x:||x||\leq B,\,x\in\mathbb{R}^{p}\\}$, ${\cal
Y}=[0,1]$, and ${\cal A}=[0,1]^{l}$ where $||\cdot||$ denotes the Euclidean
norm. For mathematical convenience, we set $D(\cdot)=\sigma(f(\cdot))$ with
$\sigma(x)=(1+\exp(-x))^{-1}$ and consider the following neural networks:
$f(x,a)=f_{{\bf w}}(x,a)=w^{\top}_{d}\kappa_{d-1}(W_{d-1}\kappa_{d-2}(\cdots
W_{1}[x^{\top},a^{\top}]^{\top}))$ and $h(x)=h_{{\bf
v}}(x)=v^{\top}_{g}\psi_{g-1}(V_{g-1}\psi_{g-2}(\cdots V_{1}x))$ where ${\bf
w}=(W_{1},\dots,W_{d-1},w_{d})\in{\bf W}$ and ${\bf
v}=(V_{1},\dots,V_{g-1},v_{g})\in{\bf V}$. We denote by ${\cal F}$ and ${\cal
H}$ the function classes of $f$ and $h$. Also, it is assumed that ${\bf
W}=\bigotimes_{i=1}^{d-1}\\{W_{i}\in\mathbb{R}^{p_{i+1}\times
p_{i}}:||W_{i}||_{F}\leq
M_{w}(i)\\}\bigotimes\\{w_{d}\in\mathbb{R}^{p_{d}\times 1}:||w_{d}||\leq
M_{w}(d)\\}$ and ${\bf
V}=\bigotimes_{i=1}^{g-1}\\{V_{i}\in\mathbb{R}^{q_{i+1}\times
q_{i}}:||V_{i}||_{F}\leq
M_{v}(i)\\}\bigotimes\\{v_{g}\in\mathbb{R}^{q_{g}\times 1}:||v_{g}||\leq
M_{v}(g)\\}$ with constants $M_{w}(\cdot)$ and $M_{v}(\cdot)$, $p_{1}=p+l$,
$q_{1}=p$, and the Frobenius norm $||\cdot||_{F}$. This also induces the class
of $D$ by functional composition, defined as ${\cal D}:=\\{f_{{\bf w}}(h_{{\bf
v}}(x),a):{\bf w}\in{\bf W},{\bf v}\in{\bf V}\\}$

We further assume that the activation functions $\psi_{u}$ and $\kappa_{t}$
are $K_{\psi}(u)$ and $K_{\kappa}(t)$-Lipschitz for all $t=1,\dots,d-1$ and
$u=1,\dots,g-1$. ReLU and Sigmoid are examples of the 1-Lipschitz functions.
These assumptions hint $1>\gamma_{1}\geq\sigma(h_{{\bf
v}}(x))\geq\gamma_{0}>0$ and $1>\nu_{1}\geq\sigma(f_{{\bf
w}}(x,a))\geq\nu_{0}>0$ for all $x$ and $a$ where the upper and the lower
limits $(\gamma_{1},\gamma_{0},\nu_{1},\nu_{0})$ depend on the assumed bounds
on the parameter spaces (i.e., $M_{w}(\cdot)$ and $M_{v}(\cdot)$) and the
Lipschitz constants of the activation functions.

Now, let’s denote $d(h;\lambda):=R(h)+\lambda R_{F}(h)$ and
$\hat{d}(h;\lambda):=\hat{R}(h)+\lambda\hat{R}_{F}(h)$ as the population-level
target function and its empirical version respectively. The estimation error
then can be characterized by evaluating the empirical solution on the
population objective based on the Rademacher complexity.

######  Definition 3 (Rademacher Complexity).

Let ${\cal H}$ be the function class of $h$. Denote by $X_{1},\dots,X_{n}$
random samples that are independent and identically distributed (i.i.d.) to
$P_{X}$. Then the Rademacher complexity ${\cal R}({\cal H})$ is defined as

$\displaystyle{\cal R}({\cal H})={\bf E}_{X,\epsilon}\left[\sup_{h\in{\cal
H}}\left|\dfrac{1}{n}\sum_{i=1}^{n}\epsilon_{i}h(X_{i})\right|\right]$

where $\epsilon_{1},\dots,\epsilon_{n}\sim{\rm Unif}\\{-1,1\\}$ i.i.d.

The Rademacher complexity for other function classes is defined in the same
fashion. Theorem 1 holds for either the cross-entropy loss with
$y\in\\{0,1\\}$ or the mean absolute error for $y\in[0,1]$.

######  Theorem 1 (Estimation Error).

Let $\hat{h}^{*}=\arg_{h\in{\cal H}}\min\hat{d}(h;\lambda)$ and define ${\cal
L}:=\\{L(y,h_{{\bf v}}(x)):{\bf v}\in V\\}$ for $x\in{\cal X}$ and $y\in{\cal
Y}$. Then for any given $\lambda\geq 0$, an upper bound of the estimation
error is $|d(\hat{h}^{*};\lambda)-\inf_{h\in{\cal H}}d(h;\lambda)|\leq$

$\displaystyle\underbrace{4{\cal R}({\cal L})+2F_{{\bf
V},\psi,B,\gamma_{0},\gamma_{1}}\sqrt{\dfrac{\log(1/\delta)}{2n}}}_{\text{from
}\hat{R}(h)}$
$\displaystyle+\underbrace{2\lambda\left(F_{\nu_{0},\nu_{1}}{\cal R}({\cal
D})+F_{{\bf W},{\bf
V},B,\kappa,\psi,l,\nu_{0},\nu_{1}}\sqrt{\dfrac{\log(1/\delta)}{2n}}\right)}_{\text{from
}\hat{R}_{F}(h)},$

with the probability $1-3\delta$ where $F_{\nu_{0},\nu_{1}}$, $F_{{\bf
V},\psi,B,\gamma_{0},\gamma_{1}}$, and $F_{{\bf W},{\bf
V},B,\kappa,\psi,l,\nu_{0},\nu_{1}}$ are constants depending on the
architectures of neural networks $h\in{\cal H}$ and $f\in{\cal F}$ whose exact
values can be found in Supplementary 9.

Theorem 1 shows how the complexity of the neural network function classes
affects the estimation accuracy through the Rademacher complexity and
$F_{\cdot}$ constants. As sample size $n$ increases (while other settings are
fixed), the Rademacher complexity generally decreases to zero, leading to
consistency according to Theorem 1. The hyperparameter $\lambda$ controls the
balance between utility and fairness. Larger $\lambda$ enforces better
fairness at the expense of possible utility loss; the next corollary studies
this utility loss.

######  Corollary 1 (Loss of Utility).

Define $\hat{h}^{*}=\arg_{h\in{\cal H}}\min\hat{d}(h;\lambda)$,
${h}^{*}_{0}=\arg_{h\in{\cal H}}\min{d}(h;\lambda=0)$, and
${h}^{*}=\arg_{h\in{\cal H}}\min{d}(h;\lambda)$. Let’s denote
$\Delta(h_{0}^{*},h^{*}):=R_{F}(h_{0}^{*})-R_{F}(h^{*})$. Then the loss of
utility has an upper bound
$|d(\hat{h}^{*};\lambda=0)-d(h^{*}_{0};\lambda=0)|\leq$

$\displaystyle 4{\cal R}({\cal L})+2F_{{\bf
V},\psi,B,\gamma_{0},\gamma_{1}}\sqrt{\dfrac{\log(1/\delta)}{2n}}+\lambda\Delta(h_{0}^{*},h^{*}),$

with 1-$\delta$ probability.

The loss of utility can take place from two main sources. Basically, the loss
of utility could stem from sampling errors, which unfold the first two terms
in the upper bound of Corollary 1. Secondly, if the true $h_{0}^{*}$, which is
of best utility, is exposed to huge discrimination against $A$ in the
population, i.e., $h$ statistically strongly depends on $A$, then the penalty
could lead to tremendous loss of utility in an effort to remove the dependency
w.r.t. $A$ (i.e., the term $\lambda\Delta(h_{0}^{*},h^{*})$).

##  5 Simulation

The performance of supervised learning with the proposed penalties is verified
in the following three scenarios: (I) discrete outcome and sensitive
attribute, (II) discrete outcome and mixed sensitive attributes, and (III)
continuous outcome and sensitive attribute. The main text only delivers
results about GSP in Scenario I and GEO in Scenario II due to the page limit.
Results about GEO in Scenario I, GSP in Scenario II, Scenario III, and more
additional studies (e.g., fair representation) are presented in Supplementary
Section 10.

We compare the performance of our model (SRS-based penalty, SBP) against
several newly proposed competing methods: HGR (Lee et al.,, 2022) that uses
the neural-net approximation to calculate the soft HGR, the feature
neutralization (Du et al., (2021), NEU) that interpolates feature points in a
hidden layer such that the mapped space is independent of sensitive
attributes, the kernel density estimation (Cho et al., (2020), KDE) that
employs $Q$-function to approximate a distribution function, and CON (Romano
et al.,, 2020) that uses the learned conditional sampler of $A$ given $Y$. All
the above in-processing models have the trade-off parameter $\lambda$ to the
fairness penalty $L_{F}$ for fairness on top of the main loss term $L_{M}$. As
a means of matching the scale of $\lambda$ for the different models, all of
them are tested under the formula $(1-\lambda)L_{M}+\lambda L_{F}$ with 5
replicated experiments for various $0<\lambda<1$. All methods learn a neural
network $h$ with 3 hidden layers with 64 nodes. For more details about
simulation setups, please refer to Supplementary Section 10.1.

Three benchmark data sets are considered for Scenario I: Adult
Data111https://archive.ics.uci.edu/ml/datasets/adult where $Y$ is whether or
not the annual income is greater than $50K and $A$ is whether an individual is
white or non-white; Law School Admission Data
222http://www.seaphe.org/databases.php where $Y$ is whether or not an
applicant receives admission and $A$ whether an individual is white or non-
white; and Credit Card Default Data1 where $Y$ is whether or not a customer
declares default and $A$ is the gender. We refer to the work (Cho et al.,,
2020) to specify the sensitive attributes for analysis. Each data set is split
by 80% and 20% for training and validation during the training course.

For the evaluation of fairness in Scenario I, statistical parity (SP) is
assessed by calculating

SP $\displaystyle=\left|{\bf E}(\hat{Y}=1|A=1)/{\bf
E}(\hat{Y}=1|A=0)-1\right|,$

where ${\cal Y}=\\{0,1\\}$, and $\hat{Y}=I(\hat{h}^{*}(X)>\tau)$ with the
threshold $\tau$ maximizing the area under ROC curve (AUC) on the validation
set. The AUC value is also used as a measure of utility. We use the
Kolmogorov-Smirnov statistics (KS) to assess GSP, which is calculated as

KS-GSP
$\displaystyle=\sum_{a\in\\{0,1\\}}\max_{h_{x}}|\hat{P}(\hat{h}^{*}(X)\leq
h_{x}|A=a)$
$\displaystyle~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}-\hat{P}(\hat{h}^{*}(X)\leq
h_{x})|,$

where $\hat{P}$ denotes the empirical distribution.

Figure 1: (Scenario I) Pareto frontiers: the first row includes pairs of SP
and AUC, and the second row shows pairs of KS-GSP and AUC from 5 experiments
for each $\lambda$. SBP (ours) tends to be more tightly in the upper-left
corner than the competitors.

Since the model intrinsically has multi-objectives (utility and fairness) in
nature, the trained models at all iterations are evaluated over the validation
data. Then, we plot the Pareto frontier curve of (AUC versus SP/KS) for each
independent run. The Pareto frontier (Emmerich and Deutz,, 2018) is a set of
solutions that are not dominated by other pairs. For instance, a pair of AUC
and SP, e.g., $(0.5,0.5)$, is not dominated by $(0.4,0.4)$ but by $(0.6,0.4)$.
This analysis helps design fair comparison studies without designing different
or non-comparable early stopping strategies for the competitors (for example,
NEU tends to degenerate to the trivial solution, which is perfectly fair but
has no utility, for long-run iterations, so comparing the last iterates among
the competing methods is not fair).

Table 2: Averages of the 5 smallest SP/KS-GSPs whose AUCs are greater than the thresholds. Those scores are selected in the Pareto solutions appearing in Figure 1. Standard deviations are in the parentheses next to the averages. All values are rounded to the third decimal place. |  Adult (AUC $\geq$ 0.85) |  Cred. Card. (AUC $\geq$ 0.75) |  Law School. (AUC $\geq$ 0.80)
---|---|---|---
| SP ($\downarrow$) | KS-GSP ($\downarrow$) | SP ($\downarrow$) | KS-GSP ($\downarrow$) | SP ($\downarrow$) | KS-GSP ($\downarrow$)
SBP | 0.001 ($\approx 0$) | 0.014 ($\approx 0$) |  $\approx 0$ ($\approx 0)$ | 0.011 (0.001) | 0.004 (0.003) | 0.018 (0.001)
HGR | 0.070 (0.005) | 0.040 ($\approx 0$) | 0.001 (0.001) | 0.023 (0.001) | 0.096 (0.003) | 0.064 (0.002)
KDE | 0.355 (0.010) | 0.124 (0.002) | 0.013 (0.009) | 0.024 ($\approx 0)$ | 0.239 (0.013) | 0.132 (0.003)
NEU | 0.446 (0.029) | 0.154 (0.003) | 0.040 (0.010) | 0.020 (0.003) | 0.214 (0.017) | 0.098 (0.005)

Figure 1 clearly illustrates the remarkable advantages of using SBP over the
competitors. The closer the set of points is to the upper left corner of the
figure, the higher performance is verified. We observe that SBP provides more
consistent and better solutions in the sense that the Pareto frontiers are
much more tightly gathered along the trade-off path than other competing
methods. Note CON is not available for SP. Table 2 reports the fairness scores
(SP/KS) of the Pareto solutions whose AUC is above the thresholds. The table
succinctly demonstrates our superiority which is consistent with Figure 1.
Those thresholds are chosen such that all methods yield a sufficient number of
candidate SP/KS scores by referring to Figure 1. More tables with different
thresholds appear in Supplementary 10.2.2.

Table 3: Averages of the 5 smallest EO/KS-GEOs whose AUCs are greater than the thresholds. Those scores are selected by referring to Figure 4 in Supplementary 10.2. Standard deviations are in the parentheses. |  |  Adult (AUC $\geq$ 0.80) |  Cred. Card. (AUC $\geq$ 0.75) |  ACSEmpl. (AUC $\geq$ 0.75)
---|---|---|---|---
Metric | Model | Race ($\downarrow$) | Age ($\downarrow$) | Gender ($\downarrow$) | Age ($\downarrow$) | Gender ($\downarrow$) | Age ($\downarrow$)
EO | SBP | 0.047 (0.004) | 0.078 (0.002) | 0.001 ($\approx$ 0) | 0.016 (0.001) | 0.066 (0.004) | 0.209 (0.002)
CON | 0.358 (0.029) | 0.162 (0.027) | 0.006 (0.003) | 0.016 (0.004) | 0.155 (0.006) | 0.316 (0.011)
HGR | 0.273 (0.073) | 0.307 (0.002) | 0.002 (0.001) | 0.033 (0.001) | 0.193 (0.005) | 0.447 (0.002)
NEU | 0.292 (0.019) | 0.201 (0.007) | 0.027 (0.007) | 0.055 (0.005) | 0.099 (0.003) | 0.276 (0.006)
KS- GEO | SBP | 0.098 (0.003) | 0.120 (0.001) | 0.040 (0.001) | 0.036 (0.001) | 0.071 (0.001) | 0.167 (0.002)
CON | 0.177 (0.001) | 0.094 (0.007) | 0.049 (0.001) | 0.038 ($\approx 0$) | 0.100 (0.002) | 0.207 (0.003)
HGR | 0.160 (0.003) | 0.208 (0.008) | 0.041 (0.003) | 0.047 (0.001) | 0.122 (0.001) | 0.320 (0.004)
NEU | 0.156 (0.005) | 0.130 (0.001) | 0.063 (0.004) | 0.054 (0.001) | 0.074 (0.001) | 0.205 (0.002)

For Scenario II, we use an additional data set
ACSEmployment333https://github.com/socialfoundations/folktables from
California in 2018 where $Y$ is whether or not an individual is employed and
$A$ is a vector of age (continuous) and gender (discrete). We choose age
(continuous) and race (discrete) for Adult and also age (continuous) and
gender (discrete) for Credit Card Default as sensitive attributes respectively
while the output variables are the same as in Scenario I. For CON, we devise a
conditional GAN model (Mirza and Osindero,, 2014) because the trivial
estimation method for $P(A|Y)$ in the original work (Romano et al.,, 2020) is
not directly applicable when $A$ is a mix of discrete and categorical
variable. For a fair comparison, we use the same network structure for $\beta$
in SBP and the discriminator for CON. To measure fairness for discrete
attributes, we define

$\displaystyle\text{EO}=\sum_{y\in{\cal Y}}\left|\frac{{\bf
E}(\hat{Y}=1|A=1,Y=y)}{{\bf E}(\hat{Y}=1|A=0,Y=y)}-1\right|,$

and $\text{KS-GEO}=$

$\displaystyle\sum_{y\in{\cal Y},a\in{\\{1,0\\}}}\max_{h_{x}}|$
$\displaystyle\hat{P}(\hat{h}^{*}(X)\leq h_{x}|A=a,Y=y)$
$\displaystyle-\hat{P}(\hat{h}^{*}(X)\leq h_{x}|Y=y)|.$

We similarly define the fairness measures for continuous attributes; please
refer to Supplementary 10.2 to see the formal definition.

Table 4: Training times (mins) for the first 1000 iterations on A30 GPU for separation in Adult, including/excluding the pre-training time. Note that for Scenario II, CON, NEU, and SBP need a pre-training (additional NN models) course but KDE and HGR do not. Refer to Supplementary 10.1 for more details. | Without | With
---|---|---
| pre-training | pre-training
Method | Scenario I | Scenario II | Scenario II
SBP | 0.22 (0.02) | 0.59 (0.05) | 0.80 (0.07)
CON | 0.21 (0.03) | 0.40 (0.03) | 1.45 (0.09)
HGR | 0.47 (0.08) | 0.77 (0.10) | -
KDE | 1.31 (0.09) | - | -
NEU | 1.01 (0.09) | 1.46 (0.03) | 3.61 (0.09)

Table 3 contrasts all models except KDE on the three data sets. Note KDE can
not be implemented for continuous sensitive attributes. The thresholds are
chosen based on Figure 4 in Supplementary 10.2. As is consistent with the
table, the figure illustrates that SBP tends to outperform the competitors for
both sensitive variables in general. For the computing time, ours is generally
more efficient as shown in Table 4.

To intuitively explain that SBP outperforms others, we comment that HGR
involves more burdensome approximation because it has to estimate four
auxiliary neural networks in the adversarial training procedure, leading to
larger estimation errors and computing time; NEU, on the other hand, does not
explicitly control fairness metrics in its training process, leading to uneven
performance. CON is significantly defeated by SBP mainly because estimating a
good conditional generator $A|Y$ via the GAN approach is more difficult than
estimating a good $\hat{\beta}$. On the other hand, we observe that SBP for
separation is robust against a poor estimation of $\hat{\beta}$ (details in
Supplementary 10.4).

In spite of such remarkable performance, an unavoidable issue of the proposed
penalty remains in regard to the adversarial game structure between $h$ and
$D$ and the training instability of such an adversarial game. As a remedy to
this obstacle, users can consider saving snapshots of $h$ being learned over
the Pareto frontier and deploying the final model having the desired Pareto
solution for the purpose.

##  6 Conclusion

The proposed fair-ML framework which employs simple random sampling is
universally applicable to classification/regression problems across various
fairness criteria, without worrying about whether sensitive attributes are
continuous or discrete, or even a mix of them. This methodological versatility
is of great importance as multifarious communities in society pay more
attention to the unbiased data-driven decision process while having more and
more diverse variables with distinctive characteristics to be protected
together (Loukina et al.,, 2019; Lee and Floridi,, 2021; Kozodoi et al.,,
2022; De-Arteaga et al.,, 2022). In future research, the fundamental concept
of unconditionally sampling sensitive attributes will be a promising tool for
promoting fairness for other related works that heavily depend on the
discreteness, mixed or even high-dimensional sensitive variables, such as
synthesizing a fair data set (Xu et al., 2019a, ) or fair uncertainty
quantification (Lu et al.,, 2022).

## References

  * Adel et al., (2019) Adel, T., Valera, I., Ghahramani, Z., and Weller, A. (2019).  One-network adversarial fairness.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 33, pages 2412–2420.
  * Agarwal et al., (2018) Agarwal, A., Beygelzimer, A., Dudík, M., Langford, J., and Wallach, H. (2018).  A reductions approach to fair classification.  In International Conference on Machine Learning, pages 60–69. PMLR.
  * Agarwal et al., (2019) Agarwal, A., Dudík, M., and Wu, Z. S. (2019).  Fair regression: Quantitative definitions and reduction-based algorithms.  In International Conference on Machine Learning, pages 120–129. PMLR.
  * Barocas et al., (2017) Barocas, S., Hardt, M., and Narayanan, A. (2017).  Fairness in machine learning.  Nips tutorial, 1:2017.
  * Beutel et al., (2017) Beutel, A., Chen, J., Zhao, Z., and Chi, E. H. (2017).  Data decisions and theoretical implications when adversarially learning fair representations.  arXiv preprint arXiv:1707.00075.
  * Caton and Haas, (2020) Caton, S. and Haas, C. (2020).  Fairness in machine learning: A survey.  arXiv preprint arXiv:2010.04053.
  * Cho et al., (2020) Cho, J., Hwang, G., and Suh, C. (2020).  A fair classifier using kernel density estimation.  Advances in neural information processing systems, 33:15088–15099.
  * Chouldechova, (2017) Chouldechova, A. (2017).  Fair prediction with disparate impact: A study of bias in recidivism prediction instruments.  Big data, 5(2):153–163.
  * Cotter et al., (2019) Cotter, A., Jiang, H., Gupta, M. R., Wang, S., Narayan, T., You, S., and Sridharan, K. (2019).  Optimization with non-differentiable constraints with applications to fairness, recall, churn, and other goals.  J. Mach. Learn. Res., 20(172):1–59.
  * Das et al., (2021) Das, S., Donini, M., Gelman, J., Haas, K., Hardt, M., Katzman, J., Kenthapadi, K., Larroy, P., Yilmaz, P., and Zafar, B. (2021).  Fairness measures for machine learning in finance.
  * De-Arteaga et al., (2022) De-Arteaga, M., Feuerriegel, S., and Saar-Tsechansky, M. (2022).  Algorithmic fairness in business analytics: Directions for research and practice.  Production and Operations Management, 31(10):3749–3770.
  * Deng et al., (2022) Deng, Z., Sun, H., Wu, Z. S., Zhang, L., and Parkes, D. C. (2022).  Reinforcement learning with stepwise fairness constraints.  arXiv preprint arXiv:2211.03994.
  * Du et al., (2021) Du, M., Mukherjee, S., Wang, G., Tang, R., Awadallah, A., and Hu, X. (2021).  Fairness via representation neutralization.  Advances in Neural Information Processing Systems, 34:12091–12103.
  * Edwards and Storkey, (2015) Edwards, H. and Storkey, A. (2015).  Censoring representations with an adversary.  arXiv preprint arXiv:1511.05897.
  * Emmerich and Deutz, (2018) Emmerich, M. T. and Deutz, A. H. (2018).  A tutorial on multiobjective optimization: fundamentals and evolutionary methods.  Natural computing, 17:585–609.
  * Goodfellow et al., (2014) Goodfellow, I. J., Pouget-Abadie, J., Mirza, M., Xu, B., Warde-Farley, D., Ozair, S., Courville, A., and Bengio, Y. (2014).  Generative adversarial networks.
  * Grari et al., (2019) Grari, V., Ruf, B., Lamprier, S., and Detyniecki, M. (2019).  Fairness-aware neural r$\backslash$’eyni minimization for continuous features.  arXiv preprint arXiv:1911.04929.
  * Gutmann and Hyvärinen, (2012) Gutmann, M. U. and Hyvärinen, A. (2012).  Noise-contrastive estimation of unnormalized statistical models, with applications to natural image statistics.  Journal of machine learning research, 13(2).
  * Hardt et al., (2016) Hardt, M., Price, E., and Srebro, N. (2016).  Equality of opportunity in supervised learning.  Advances in neural information processing systems, 29.
  * Ji et al., (2021) Ji, K., Zhou, Y., and Liang, Y. (2021).  Understanding estimation and generalization error of generative adversarial networks.  IEEE Transactions on Information Theory, 67(5):3114–3129.
  * Jung et al., (2023) Jung, S., Park, T., Chun, S., and Moon, T. (2023).  Re-weighting based group fairness regularization via classwise robust optimization.  arXiv preprint arXiv:2303.00442.
  * Kamiran and Calders, (2012) Kamiran, F. and Calders, T. (2012).  Data preprocessing techniques for classification without discrimination.  Knowledge and information systems, 33(1):1–33.
  * Komiyama et al., (2018) Komiyama, J., Takeda, A., Honda, J., and Shimao, H. (2018).  Nonconvex optimization for regression with fairness constraints.  In International conference on machine learning, pages 2737–2746. PMLR.
  * Kotelnikov et al., (2023) Kotelnikov, A., Baranchuk, D., Rubachev, I., and Babenko, A. (2023).  Tabddpm: Modelling tabular data with diffusion models.  In International Conference on Machine Learning, pages 17564–17579. PMLR.
  * Kozodoi et al., (2022) Kozodoi, N., Jacob, J., and Lessmann, S. (2022).  Fairness in credit scoring: Assessment, implementation and profit implications.  European Journal of Operational Research, 297(3):1083–1094.
  * Lee et al., (2022) Lee, J., Bu, Y., Sattigeri, P., Panda, R., Wornell, G., Karlinsky, L., and Feris, R. (2022).  A maximal correlation approach to imposing fairness in machine learning.  In ICASSP 2022-2022 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), pages 3523–3527. IEEE.
  * Lee and Floridi, (2021) Lee, M. S. A. and Floridi, L. (2021).  Algorithmic fairness in mortgage lending: from absolute conditions to relational trade-offs.  Minds and Machines, 31(1):165–191.
  * Li et al., (2021) Li, X., Cui, Z., Wu, Y., Gu, L., and Harada, T. (2021).  Estimating and improving fairness with adversarial learning.  arXiv preprint arXiv:2103.04243.
  * Loukina et al., (2019) Loukina, A., Madnani, N., and Zechner, K. (2019).  The many dimensions of algorithmic fairness in educational applications.  In Proceedings of the Fourteenth Workshop on Innovative Use of NLP for Building Educational Applications, pages 1–10.
  * Lu et al., (2022) Lu, C., Lemay, A., Chang, K., Höbel, K., and Kalpathy-Cramer, J. (2022).  Fair conformal predictors for applications in medical imaging.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 36, pages 12008–12016.
  * Mary et al., (2019) Mary, J., Calauzenes, C., and El Karoui, N. (2019).  Fairness-aware learning for continuous attributes and treatments.  In International Conference on Machine Learning, pages 4382–4391. PMLR.
  * Mhasawade et al., (2021) Mhasawade, V., Zhao, Y., and Chunara, R. (2021).  Machine learning and algorithmic fairness in public and population health.  Nature Machine Intelligence, 3(8):659–666.
  * Mirza and Osindero, (2014) Mirza, M. and Osindero, S. (2014).  Conditional generative adversarial nets.  arXiv preprint arXiv:1411.1784.
  * Pagano et al., (2023) Pagano, T. P., Loureiro, R. B., Lisboa, F. V., Peixoto, R. M., Guimarães, G. A., Cruz, G. O., Araujo, M. M., Santos, L. L., Cruz, M. A., Oliveira, E. L., et al. (2023).  Bias and unfairness in machine learning models: A systematic review on datasets, tools, fairness metrics, and identification and mitigation methods.  Big data and cognitive computing, 7(1):15.
  * Pleiss et al., (2017) Pleiss, G., Raghavan, M., Wu, F., Kleinberg, J., and Weinberger, K. Q. (2017).  On fairness and calibration.  Advances in neural information processing systems, 30.
  * Quadrianto and Sharmanska, (2017) Quadrianto, N. and Sharmanska, V. (2017).  Recycling privileged learning and distribution matching for fairness.  Advances in neural information processing systems, 30.
  * Quinzan et al., (2022) Quinzan, F., Khanna, R., Hershcovitch, M., Cohen, S., Waddington, D. G., Friedrich, T., and Mahoney, M. W. (2022).  Fast feature selection with fairness constraints.  arXiv preprint arXiv:2202.13718.
  * Romano et al., (2020) Romano, Y., Bates, S., and Candes, E. (2020).  Achieving equalized odds by resampling sensitive attributes.  Advances in neural information processing systems, 33:361–371.
  * Sattigeri et al., (2019) Sattigeri, P., Hoffman, S. C., Chenthamarakshan, V., and Varshney, K. R. (2019).  Fairness gan: Generating datasets with fairness properties using a generative adversarial network.  IBM Journal of Research and Development, 63(4/5):3–1.
  * Scutari et al., (2022) Scutari, M., Panero, F., and Proissl, M. (2022).  Achieving fairness with a simple ridge penalty.  Statistics and Computing, 32(5):1–17.
  * van Breugel et al., (2021) van Breugel, B., Kyono, T., Berrevoets, J., and van der Schaar, M. (2021).  Decaf: Generating fair synthetic data using causally-aware generative networks.  Advances in Neural Information Processing Systems, 34:22221–22233.
  * Wang et al., (2023) Wang, J., Lu, D., Davidson, I., and Bai, Z. (2023).  Scalable spectral clustering with group fairness constraints.  In International Conference on Artificial Intelligence and Statistics, pages 6613–6629. PMLR.
  * Xian et al., (2023) Xian, R., Yin, L., and Zhao, H. (2023).  Fair and optimal classification via post-processing.
  * Xie et al., (2017) Xie, Q., Dai, Z., Du, Y., Hovy, E., and Neubig, G. (2017).  Controllable invariance through adversarial feature learning.  Advances in neural information processing systems, 30.
  * (45) Xu, D., Wu, Y., Yuan, S., Zhang, L., and Wu, X. (2019a).  Achieving causal fairness through generative adversarial networks.  In Proceedings of the Twenty-Eighth International Joint Conference on Artificial Intelligence.
  * Xu et al., (2018) Xu, D., Yuan, S., Zhang, L., and Wu, X. (2018).  Fairgan: Fairness-aware generative adversarial networks.  In 2018 IEEE International Conference on Big Data (Big Data), pages 570–575. IEEE.
  * Xu et al., (2022) Xu, J., Xiao, Y., Wang, W. H., Ning, Y., Shenkman, E. A., Bian, J., and Wang, F. (2022).  Algorithmic fairness in computational medicine.  EBioMedicine, 84:104250.
  * (48) Xu, L., Skoularidou, M., Cuesta-Infante, A., and Veeramachaneni, K. (2019b).  Modeling tabular data using conditional gan.  Advances in neural information processing systems, 32.
  * Zafar et al., (2019) Zafar, M. B., Valera, I., Gomez-Rodriguez, M., and Gummadi, K. P. (2019).  Fairness constraints: A flexible approach for fair classification.  The Journal of Machine Learning Research, 20(1):2737–2778.
  * Zeng et al., (2022) Zeng, X., Dobriban, E., and Cheng, G. (2022).  Bayes-optimal classifiers under group fairness.  arXiv preprint arXiv:2202.09724.
  * Zhang et al., (2018) Zhang, B. H., Lemoine, B., and Mitchell, M. (2018).  Mitigating unwanted biases with adversarial learning.  In Proceedings of the 2018 AAAI/ACM Conference on AI, Ethics, and Society, pages 335–340.
  * Zhao et al., (2019) Zhao, H., Coston, A., Adel, T., and Gordon, G. J. (2019).  Conditional learning of fair representations.  arXiv preprint arXiv:1910.07162.

Supplementary Materials: Fair Supervised Learning with A Simple Random Sampler
of Sensitive Attributes

##  7 Details of Section 3.1

###  7.1 Proof of Proposition 1

Let’s recall

$\displaystyle R_{F}(h;D)={\bf E}_{X,A}[\log D(h(X),A)]+{\bf
E}_{X,A^{\prime}}[\log(1-D(h(X),A^{\prime}))],$

where $A$ and $A^{\prime}$ are independent but identically distributed. Let
$s=h(x)$ where $x$ is a realization of $X$ and also $a$ is of $A$. The loss
function can be written

$\displaystyle R_{F}(h;D)$ $\displaystyle=\int\log
D(s,a)p(s,a)dsda+\int\log(1-D(s,a^{\prime}))p(s,a^{\prime})dsda^{\prime},$
$\displaystyle=\int\log D(s,a)p(s|a)p(a)+\log(1-D(s,a))p(s)p(a)dsda.$

By the proof of Proposition 1 in Goodfellow et al., (2014), $R_{F}(h;D)$ is
maximized at

$\displaystyle
D^{*}(s,a)=\dfrac{p(s|a)p(a)}{p(s|a)p(a)+p(s)p(a)}=\dfrac{p(s|a)}{p(s|a)+p(s)}.$

for any $s\in{\cal S}$ and $a\in{\cal A}$.

###  7.2 Algorithm of GSP penalty

Let $\\{(x_{i},a_{i},y_{i})\\}_{i=1}^{n}$ be a set of the realization of
$\\{X_{i},A_{i},Y_{i}\\}_{i=1}^{n}$. The simple random sampler of $A$ can be
easily obtained by the minibatch construction (Algorithm 1). It is allowed to
update ${\bf w}$ up to $T^{\prime}\geq 1$ times for every single update of
${\bf v}$ for better approximation for $\hat{R}_{F}(h)$, e.g.,
$T^{\prime}=1,\dots,10$. In every iteration, $D$ serves as a fairness critic
by quantifying the degree of discrimination against GSP. The model $h$ is then
trained to minimize the risk but at the same time to be debiased such that $D$
would not capture the discrimination. Note ${\bf w}$ and ${\bf v}$ are
parameters of $D$ and $h$. The algorithm to control GSP appears in Algorithm
2.

Data: Let’s ${\bf D}_{n}=\\{(x_{i},a_{i},y_{i})\\}_{i=1}^{n}$ be a set of
training data set. Set the minibatch size $n_{b}$. The subscript $(i)$ denotes
the $i$th drawn sample.

Result: ${\bf D}_{n_{b}}$

${\bf D}=\\{(x_{(i)},a_{(i)},y_{(i)})\\}_{i=1}^{n_{b}}$ is randomly drawn from
${\bf D}_{n}$.

${\bf
D}^{\prime}=\\{(x_{(i)}^{\prime},a_{(i)}^{\prime},y_{(i)}^{\prime})\\}_{i=1}^{n_{b}}$
is randomly drawn from ${\bf D}_{n}\setminus{\bf D}$.

Construct ${\bf
D}_{n_{b}}=\\{(x_{(i)},a_{(i)},y_{(i)},a_{(i)}^{\prime})\\}_{i=1}^{n_{b}}$ by
selecting $\\{a_{(i)}^{\prime}\\}_{i=1}^{n_{b}}$ from ${\bf D}^{\prime}$ and
combining it into ${\bf D}$.

Algorithm 1 Minibatch Construction (MC) at the $t$th iteration

Data: Let’s ${\bf D}_{n}=\\{(x_{i},a_{i},y_{i})\\}_{i=1}^{n}$ be a set of
training data set. Denote the $t$th iterate of model parameters by ${\bf
w}^{(t)}$ and ${\bf v}^{(t)}$ and $\alpha$ is a learning rate. Fix the number
of training iterations $T$ and $T^{\prime}$; set $t=0$; and initialize ${\bf
w}^{(0)}$ and ${\bf v}^{(0)}$.

Result: $h_{{\bf v}^{(T)}}$

while  _$t\leq T$_ do

$\\{(x_{(i)},a_{(i)},y_{(i)},a_{(i)}^{\prime})\\}_{i=1}^{n_{b}}=\text{MC}({\bf
D}_{n})$

$t^{\prime}=0$

while  _$t^{\prime}\leq T^{\prime}$_ do

$t^{\prime}=t^{\prime}+1$

$\hat{R}_{F}(h_{{\bf v}};D_{{\bf
w}^{(t)}})=\frac{1}{n_{b}}\sum_{i=1}^{n_{b}}\left(\log D_{{\bf
w}^{(t)}}(h_{{\bf v}^{(t)}}(x_{(i)}),a_{(i)})+\log(1-D_{{\bf w}^{(t)}}(h_{{\bf
v}^{(t)}}(x_{(i)}),a^{\prime}_{(i)}))\right)$

${\bf w}^{(t+1)}={\bf w}^{(t)}+\alpha\dfrac{\partial}{\partial{\bf
w}^{(t)}}\hat{R}_{F}(h_{{\bf v}^{(t)}};D_{{\bf w}^{(t)}})$

end while

$\hat{R}(h_{{\bf v}^{(t)}})=\frac{1}{n_{b}}\sum_{i=1}^{n_{b}}L(y_{(i)},h_{{\bf
v}^{(t)}}(x_{(i)}))$

${\bf v}^{(t+1)}={\bf v}^{(t)}-\alpha\dfrac{\partial}{\partial{\bf
v}^{(t)}}(\hat{R}(h_{{\bf v}^{(t)}})+\lambda\hat{R}_{F}(h_{{\bf
v}^{(t)}};D_{{\bf w}^{(t)}}))$

end while

Algorithm 2 Generalized Statistical Parity

##  8 Details of Section 3.2

###  8.1 Proof of Proposition 2

The penalty for separation is

$\displaystyle R_{F}(h;D)$ $\displaystyle={\bf E}_{X,A,Y}[\log
D(h(X),A,Y)]+{\bf E}_{A^{\prime}}{\bf
E}_{X,Y}[\beta(A^{\prime},Y)\log(1-D(h(X),A^{\prime},Y))],$
$\displaystyle=\int\log
D(s,a,y)p(s|a,y)p(a,y)dadyds+\beta(a^{\prime},y)\log(1-D(s,a^{\prime},y))p(s|y)p(a^{\prime})p(y)da^{\prime}dyds,$
$\displaystyle=\int\log
D(s,a,y)p(s|a,y)p(a,y)+\beta(a,y)\log(1-D(s,a,y))p(s|y)p(a)p(y)dadyds.$

By the same argument in Supplementary 7.1, $D^{*}=\arg_{D}\max R_{F}(h:D)$ has
the form of

$\displaystyle
D^{*}(s,a,y;\beta)=\dfrac{p(s|a,y)p(a,y)}{p(s|a,y)p(a,y)+\beta(a,y)p(s|y)p(a)p(y)}=\dfrac{p(s|a,y)}{p(s|a,y)+\beta(a,y)p(s|y)\dfrac{p(a)p(y)}{p(a,y)}}.$

###  8.2 Algorithm of GEO penalty

GEO needs a pre-training step for $D_{\beta}$ that is also modeled by a neural
network with parameter ${\bf u}$. We additionally denote by $D_{{\bf
u}^{(l)}}$ the $l$th iterate of $D_{\beta}$. With the same notation in
Supplementary 7.2, the numerical algorithm for GEO appears in Algorithm 3.

Data: Let’s $\\{(x_{i},a_{i},y_{i})\\}_{i=1}^{n}$ be a set of training data
set and obtain $\\{a_{i}^{\prime}\\}_{i=1}^{n}$ by permuting
$\\{a_{i}\\}_{i=1}^{n}$. Denote the $t$th iterate of model parameters by ${\bf
u}^{(t)}$, ${\bf w}^{(t)}$ and ${\bf v}^{(t)}$ and $\alpha$ is a learning
rate. Fix $L$, $T$, and $T^{\prime}$; set $t=0$ and $l=0$; and initialize
${\bf u}^{(0)}$, ${\bf w}^{(0)}$ and ${\bf v}^{(0)}$.

Result: $h_{{\bf v}^{(T)}}$

while  _$l\leq L$_ do

$\\{(x_{(i)},a_{(i)},y_{(i)},a_{(i)}^{\prime})\\}_{i=1}^{n_{b}}=\text{MC}({\bf
D}_{n})$

$\hat{R}_{\beta}(D_{{\bf
u}^{(l)}})=\frac{1}{n_{b}}\sum_{i=1}^{n_{b}}\left(\log D_{{\bf
u}^{(l)}}(a_{(i)},y_{(i)})+\log(1-D_{{\bf
u}^{(l)}}(a^{\prime}_{(i)},y_{(i)}))\right)$

${\bf u}^{(l+1)}={\bf u}^{(l)}+\alpha\dfrac{\partial}{\partial{\bf
u}^{(l)}}\hat{R}_{\beta}(D_{{\bf u}^{(l)}})$

end while

while  _$t\leq T$_ do

$\\{(x_{(i)},a_{(i)},y_{(i)},a_{(i)}^{\prime})\\}_{i=1}^{n_{b}}=\text{MC}({\bf
D}_{n})$

$t^{\prime}=0$

while  _$t^{\prime}\leq T^{\prime}$_ do

$t^{\prime}=t^{\prime}+1$

$\hat{R}_{F}(h_{{\bf v}};D_{{\bf
w}^{(t)}})=\frac{1}{n_{b}}\sum_{i=1}^{n_{b}}\left(\log D_{{\bf
w}^{(t)}}(h_{{\bf v}^{(t)}}(x_{(i)}),a_{(i)},y_{(i)})+\dfrac{D_{{\bf
u}^{(L)}}(a_{(i)},y_{(i)})}{1-D_{{\bf
u}^{(L)}}(a_{(i)},y_{(i)})}\log(1-D_{{\bf w}^{(t)}}(h_{{\bf
v}^{(t)}}(x_{(i)}),a^{\prime}_{(i)},y_{(i)}))\right)$

${\bf w}^{(t+1)}={\bf w}^{(t)}+\alpha\dfrac{\partial}{\partial{\bf
w}^{(t)}}\hat{R}_{F}(h_{{\bf v}^{(t)}};D_{{\bf w}^{(t)}})$

end while

$\hat{R}(h_{{\bf v}^{(t)}})=\frac{1}{n_{b}}\sum_{i=1}^{n_{b}}L(y_{(i)},h_{{\bf
v}^{(t)}}(x_{(i)}))$

${\bf v}^{(t+1)}={\bf v}^{(t)}-\alpha\dfrac{\partial}{\partial{\bf
v}^{(t)}}(\hat{R}(h_{{\bf v}^{(t)}})+\lambda\hat{R}_{F}(h_{{\bf
v}^{(t)}};D_{{\bf w}^{(t)}}))$

end while

Algorithm 3 Generalized Equalized Odds

###  8.3 Estimation performance from $D_{\beta}$ on a toy experiment

We test the estimation performance of $D_{\beta}$ for the density ratio
estimator $\beta$. Suppose $p(y|a)=P(Y=y|A=a)=(1+\exp(-a))^{-1}$ with
$P(A=a)=0.5$ for $y,a\in\\{0,1\\}$. $D_{\beta}$ is assumed to have 2 dense
hidden layers with 64 nodes and Relu activation functions, and its output
layer has a single dimension with Sigmoid activation. In the total number of
10000 training iterations with 100 minibatch size, the last iterate is
selected as the point estimate of $p(y|a)/p(y)$. Based on 5 replicated
experiments, we calculate the averages and the standard deviation of the point
estimates. As it is shown in Table 5, $\hat{\beta}$ successfully estimates the
true density ratios.

Table 5: Estimation performance for $\hat{\beta}$ in the toy example. The point estimate is found by averaging 5 outputs. Std. implies the standard deviation calculated from the 5 outputs.  |  $p(1|1)/p(1)$ |  $p(1|0)/p(1)$ |  $p(0|1)/p(0)$ |  $p(0|0)/p(0)$
---|---|---|---|---
True | 1.1877 | 0.8123 | 0.6995 | 1.3005
Point Estimate (Std.) | 1.1810 (0.0077) | 0.8126 (0.0029) | 0.7082 (0.0018) | 1.2955 (0.0048)

##  9 Details of Section 4

###  9.1 Proof of Theorem 1

Our theory specifically deals with binary cross-entropy and mean absolute
error for binary and continuous outcomes respectively. The overall proof
strategy is to check the bounded difference condition to use McDiarmid’s
inequality.

The estimation error can be defined as

$\displaystyle|d(\hat{h}^{*};\lambda)-\inf_{h\in{\cal H}}d(h;\lambda)|$ (3)

where $\hat{h}^{*}=\arg_{h}\min\hat{d}(h;\lambda)$. The estimation error is
further decomposed as follows.

$\displaystyle
d(\hat{h}^{*};\lambda)-\inf_{h}d(h;\lambda)=d(\hat{h}^{*};\lambda)-d(h^{*};\lambda)$
$\displaystyle=d(\hat{h}^{*};\lambda)-\hat{d}(\hat{h}^{*};\lambda)$ (4)
$\displaystyle\quad+\hat{d}(h^{*};\lambda)-d(h^{*};\lambda)$ (5)
$\displaystyle\quad+\hat{d}(\hat{h}^{*};\lambda)-\hat{d}(h^{*};\lambda),$ (6)

where $h^{*}=\arg_{h}\min d(h;\lambda)$, and it is trivial to see
$\eqref{thm:bdd3}\leq 0$. The first line (4) is equivalent to

$\displaystyle
d(\hat{h}^{*};\lambda)-\hat{d}(\hat{h}^{*};\lambda)=R(\hat{h}^{*})-\hat{R}(\hat{h}^{*})+\lambda(R_{F}(\hat{h}^{*})-\hat{R}_{F}(\hat{h}^{*})).$

For (5), we have

$\displaystyle\hat{d}(h^{*};\lambda)-d(h^{*};\lambda)=\hat{R}(h^{*})-R(h^{*})+\lambda(\hat{R}_{F}(h^{*})-R_{F}(h^{*})).$

Therefore, the estimation error is bounded by

$\displaystyle|d(\hat{h}^{*};\lambda)-\inf_{h}d(h;\lambda)|\leq
2\underbrace{\sup_{h}|R(h)-\hat{R}(h)|}_{\text{I}}+2\lambda\underbrace{\sup_{h}|R_{F}(h)-\hat{R}_{F}(h)|}_{\text{II}}.$
(7)

Since $\hat{R}(h)$ is an empirical risk function, we denote
$T((X_{1},Y_{1}),\dots,(X_{n},Y_{n}))=\sup_{h}|R(h)-\hat{R}(h)|$. For the
binary cross-entropy loss function, i.e.,
$L(Y_{i},h(X_{i}))=Y_{i}\log\sigma(h(X_{i}))+(1-Y_{i})\log(1-\sigma(h(X_{i})))$,
the bounded difference of the $i$th differing variable can be bounded by

$\displaystyle|T((X_{1},Y_{1}),\dots,(X_{i},Y_{i}),\dots,(X_{n},Y_{n}))$
$\displaystyle-T((X_{1},Y_{1}),\dots,(X_{i}^{\dagger},Y_{i}^{\dagger}),\dots,(X_{n},Y_{n}))|$
$\displaystyle\leq\sup_{h}|\dfrac{1}{n}(L(Y_{i},h(X_{i}))-L(Y_{i}^{\dagger},h(X_{i}^{\dagger})))|,$
$\displaystyle\leq\dfrac{1}{n}\sup_{h}|Y_{i}\log\sigma(h(X_{i}))+(1-Y_{i})\log(1-\sigma(h(X_{i})))-$
$\displaystyle
Y_{i}^{\dagger}\log\sigma(h(X_{i}^{\dagger}))-(1-Y_{i}^{\dagger})\log(1-\sigma(h(X_{i}^{\dagger})))|.$
(8)

Without loss of generality, let’s consider (i) $Y_{i}=Y_{i}^{\dagger}=1$, (ii)
$Y_{i}=Y_{i}^{\dagger}=0$, and (iii) $Y_{i}=1$ and $Y_{i}^{\dagger}=0$. Since
the sigmoid function $\sigma$ is a 1-Lipschitz function, (i) upper bounds

$\displaystyle|\log\sigma(h(X_{i}))-\log\sigma(h(X_{i}^{\dagger}))|$
$\displaystyle\leq\dfrac{1}{\gamma_{0}}|\sigma(h(X_{i}))-\sigma(h(X_{i}^{\dagger}))|,$
$\displaystyle\leq\dfrac{1}{\gamma_{0}}|h(X_{i})-h(X_{i}^{\dagger})|,$
$\displaystyle\leq\dfrac{1}{\gamma_{0}}\prod_{j=1}^{g}M_{v}(j)\prod_{k=1}^{g-1}K_{\psi}(k)\times
2B.$

The first inequality comes from the Lipschitz property of the logarithm whose
domain is bounded below by a positive constant. The Cauchy-Schwarz inequality
and the Lipschitz conditions of the activation functions lead the last
inequality. Similarly, for (ii), we have

$\displaystyle|\log(1-\sigma(h(X_{i})))-\log(1-\sigma(h(X_{i}^{\dagger})))|$
$\displaystyle\leq\dfrac{1}{1-\gamma_{1}}|\sigma(h(X_{i}))-\sigma(h(X_{i}^{\dagger}))|,$
$\displaystyle\leq\dfrac{1}{1-\gamma_{1}}\prod_{j=1}^{g}M_{v}(j)\prod_{k=1}^{g-1}K_{\psi}(k)\times
2B,$

For (iii), we obtain

$\displaystyle|\log\sigma(h(X_{i}))-\log(1-\sigma(h(X_{i}^{\dagger})))|$
$\displaystyle\leq\dfrac{1}{\min\\{\gamma_{0},1-\gamma_{1}\\}}|\sigma(h(X_{i}))+\sigma(h(X_{i}^{\dagger}))-1|\leq\dfrac{|2\gamma_{1}-1|}{\min\\{\gamma_{0},1-\gamma_{1}\\}}.$

Therefore, the upper bound (9.1) is represented as

$\displaystyle\max\left\\{\dfrac{|2\gamma_{1}-1|}{n\gamma_{0,1}},\dfrac{1}{n\gamma_{0,1}}\prod_{j=1}^{g}M_{v}(j)\prod_{k=1}^{g-1}K_{\psi}(k)\times
2B\right\\},$

where $\gamma_{0,1}=\min\\{\gamma_{0},1-\gamma_{1}\\}$.

On the one hand, if the underlying loss function has the mean absolute error,
i.e., $L(Y_{i},h(X_{i}))=|Y_{i}-h(X_{i})|$, the upper bound of the bounded
difference is

$\displaystyle|T((X_{1},Y_{1}),\dots,(X_{i},Y_{i}),\dots,(X_{n},Y_{n}))$
$\displaystyle-T((X_{1},Y_{1}),\dots,(X_{i}^{\dagger},Y_{i}^{\dagger}),\dots,(X_{n},Y_{n}))|$
$\displaystyle\leq\sup_{h}|\dfrac{1}{n}(L(Y_{i},h(X_{i}))-L(Y_{i}^{\dagger},h(X_{i}^{\dagger})))|,$
$\displaystyle\leq\dfrac{1}{n}\sup_{h}(|Y_{i}-Y_{i}^{\dagger}|+|h(X_{i})-h(X_{i}^{\dagger})|),$
$\displaystyle\leq\dfrac{1}{n}\left(1+\prod_{j=1}^{g}M_{v}(j)\prod_{k=1}^{g-1}K_{\psi}(k)\times
2B\right).$

Next, we take expectation to (I) in (7) with respect to the random samples. We
obtain

$\displaystyle{\bf E}_{X,Y}\sup_{h}|R(h)-\hat{R}(h)|$ $\displaystyle={\bf
E}_{X,Y}\sup_{h}|{\bf
E}_{\tilde{X},\tilde{Y}}\dfrac{1}{n}\sum_{i=1}^{n}L(\tilde{Y}_{i},h(\tilde{X}_{i}))-\dfrac{1}{n}\sum_{i=1}^{n}L(Y_{i},h(X_{i}))|,$
$\displaystyle\leq{\bf
E}_{X,Y,\tilde{X},\tilde{Y}}\sup_{h}|\dfrac{1}{n}\sum_{i=1}^{n}L(\tilde{Y}_{i},h(\tilde{X}_{i}))-\dfrac{1}{n}\sum_{i=1}^{n}L(Y_{i},h(X_{i}))|,$
$\displaystyle={\bf
E}_{X,Y,\tilde{X},\tilde{Y},\epsilon}\sup_{h}|\dfrac{1}{n}\sum_{i=1}^{n}\epsilon_{i}(L(\tilde{Y}_{i},h(\tilde{X}_{i}))-L(Y_{i},h(X_{i})))|,$
$\displaystyle\leq 2{\bf
E}_{X,Y,\epsilon}\sup_{h}|\dfrac{1}{n}\sum_{i=1}^{n}\epsilon_{i}L(Y_{i},h(X_{i}))|:=2{\cal
R}({\cal L})$

where $\epsilon_{i}\sim{\rm Unif}\\{-1,1\\}$ are i.i.d,
$(\tilde{X}_{i},\tilde{Y}_{i})$ are i.i.d. copies (ghost samples) of $(X,Y)$,
and ${\cal L}:=\\{L(y,h_{{\bf v}}(x)):{\bf v}\in V\\}$ is a function class of
the given loss $L$.

Therefore, McDiarmid’s inequality implies

$\displaystyle\text{I}\leq 2{\cal R}({\cal L})+F_{{\bf
V},\psi,B,\gamma_{0},\gamma_{1}}\sqrt{\dfrac{\log(1/\delta)}{2n}},$ (9)

with the 1-$\delta$ probability, where ${\cal R}$ is the Rademacher complexity
of ${\cal L}$, and the involved constant is

$\displaystyle F_{{\bf
V},\psi,B,\gamma_{0},\gamma_{1}}=\begin{cases}\max\left\\{\dfrac{|2\gamma_{1}-1|}{\gamma_{0,1}},\dfrac{1}{\gamma_{0,1}}\prod_{j=1}^{g}M_{v}(j)\prod_{k=1}^{g-1}K_{\psi}(k)\times
2B\right\\}\quad\text{if }L\text{ is cross-enstropy}\\\
\left(1+\prod_{j=1}^{g}M_{v}(j)\prod_{k=1}^{g-1}K_{\psi}(k)\times
2B\right)\quad\text{if }L\text{ is mean absolute difference}\end{cases}$

Similarly, we can upper bound (II). First, we check the bounded difference
condition that is induced by the neural penalty, i.e.,

II $\displaystyle\leq\underbrace{\sup_{h,D}[|{\bf E}_{X,A}\log
D(h(X),A)-\hat{{\bf E}}\log
D(h(X),A)|]}_{U_{1}((X_{1},A_{1}),\dots,(X_{i},A_{i}),\dots,(X_{n},A_{n}))}$
$\displaystyle+\underbrace{\sup_{h,D}[|{\bf
E}_{X,A^{\prime}}\log(1-D(h(X),A^{\prime}))-\hat{{\bf
E}}\log(1-D(h(X),A^{\prime}))|]}_{U_{2}((X_{1},A^{\prime}_{1}),\dots,(X_{i},A^{\prime}_{i}),\dots,(X_{n},A^{\prime}_{n}))}.$

Then, for any differing the $i$th coordinate,

$\displaystyle|U_{1}((X_{1},A_{1}),\dots,(X_{i},A_{i}),\dots,(X_{n},A_{n}))$
$\displaystyle-
U_{1}((X_{1},A_{1}),\dots,(X_{i}^{\dagger},A_{i}^{\dagger}),\dots,(X_{n},A_{n}))|$
$\displaystyle\leq\dfrac{1}{n}\sup_{h,D}|\log D(h(X_{i}),A_{i})-\log
D(h(X_{i}^{\dagger}),A_{i}^{\dagger})|,$
$\displaystyle\leq\dfrac{1}{n\nu_{0}}|\sigma(f(h(X_{i}),A_{i}))-\sigma(f(h(X_{i}^{\dagger}),A_{i}^{\dagger}))|$
$\displaystyle\leq\dfrac{1}{n\nu_{0}}\prod_{k=1}^{d}M_{w}(k)\prod_{j=1}^{d-1}K_{\kappa}(j)\times||[h(X_{i}),A_{i}]-[h(X_{i}^{\dagger}),A_{i}^{\dagger}]||,$
$\displaystyle\leq\dfrac{1}{n\nu_{0}}\prod_{k=1}^{d}M_{w}(k)\prod_{j=1}^{d-1}K_{\kappa}(j)\times\sqrt{l+4B^{2}\left(\prod_{k=1}^{g}M_{v}(k)\prod_{l=1}^{g-1}K_{\psi}(l)\right)^{2}}$
$\displaystyle:=\dfrac{F_{{\bf W},{\bf V},B,\kappa,\psi,l}}{n\nu_{0}}.$

Likewise, we get the similar result for $U_{2}$ as follows,

$\displaystyle|U_{2}((X_{1},A^{\prime}_{1}),\dots,(X_{i},A^{\prime}_{i}),\dots,(X_{n},A^{\prime}_{n}))$
$\displaystyle-
U_{2}((X_{1},A^{\prime}_{1}),\dots,(X_{i}^{\dagger},(A_{i}^{\prime})^{\dagger}),\dots,(X_{n},A^{\prime}_{n}))|$
$\displaystyle\leq\dfrac{F_{{\bf W},{\bf V},B,\kappa,\psi,l}}{n(1-\nu_{1})}$

In addition,

$\displaystyle{\bf
E}_{X,A}U_{1}((X_{1},A_{1}),\dots,(X_{i},A_{i}),\dots,(X_{n},A_{n}))$
$\displaystyle\leq{\bf E}_{X,A}{\bf
E}_{\tilde{X},\tilde{A}}\sup_{h,D}|\dfrac{1}{n}\sum_{i=1}^{n}\log
D(h(\tilde{X}_{i}),\tilde{A}_{i})-\log D(h(X_{i}),A_{i})|,$
$\displaystyle\leq\dfrac{1}{\nu_{0}}{\bf E}_{X,A}{\bf
E}_{\tilde{X},\tilde{A}}\sup_{h,f}|\dfrac{1}{n}\sum_{i=1}^{n}f(h(\tilde{X}_{i}),\tilde{A}_{i})-f(h(X_{i}),A_{i})|,$
$\displaystyle=\dfrac{1}{\nu_{0}}{\bf E}_{X,A}{\bf
E}_{\tilde{X},\tilde{A},\epsilon}\sup_{h,f}|\dfrac{1}{n}\sum_{i=1}^{n}\epsilon_{i}(f(h(\tilde{X}_{i}),\tilde{A}_{i})-f(h(X_{i}),A_{i}))|,$
$\displaystyle\leq\dfrac{2}{\nu_{0}}{\bf
E}_{X,A,\epsilon}\sup_{h,f}|\dfrac{1}{n}\sum_{i=1}^{n}\epsilon_{i}f(h(X_{i}),A_{i})|:=\dfrac{2}{\nu_{0}}{\cal
R}({\cal D})$

where $\epsilon_{i}\sim{\rm Unif}\\{-1,1\\}$, are i.i.d.,
$(\tilde{X}_{i},\tilde{A}_{i})$ are i.i.d. copies (ghost sample) of $(X,A)$,
and ${\cal D}=\\{f_{{\bf w}}(h_{{\bf v}}(x),a):{\bf w}\in{\bf W},{\bf
v}\in{\bf V}\\}$ is a compositional function class. Therefore, with the
probability $1-\delta$,

$\displaystyle
U_{1}((X_{1},A_{1}),\dots,(X_{i},A_{i}),\dots,(X_{n},A_{n}))\leq\dfrac{2}{\nu_{0}}{\cal
R}({\cal D})+\dfrac{1}{\nu_{0}}F_{{\bf W},{\bf
V},B,\kappa,\psi,l}\sqrt{\dfrac{\log(1/\delta)}{2n}}.$

Similarly for $U_{2}$, with $1-\delta$ probability,

$\displaystyle
U_{2}((X_{1},A^{\prime}_{1}),\dots,(X_{i},A^{\prime}_{i}),\dots,(X_{n},A^{\prime}_{n}))\leq\dfrac{2}{1-\nu_{1}}{\cal
R}({\cal D})+\dfrac{1}{(1-\nu_{1})}F_{{\bf W},{\bf
V},B,\kappa,\psi,l}\sqrt{\dfrac{\log(1/\delta)}{2n}}.$

Therefore, with $1-2\delta$ probability,

$\displaystyle\text{II}\leq
2\left(\dfrac{1}{\nu_{0}}+\dfrac{1}{1-\nu_{1}}\right){\cal R}({\cal
D})+\left(\dfrac{1+\nu_{0}-\nu_{1}}{\nu_{0}(1-\nu_{1})}\right)F_{{\bf W},{\bf
V},B,\kappa,\psi,l}\sqrt{\dfrac{\log(1/\delta)}{2n}}.$

Consequently, by combining I and II, the estimation error is bounded above by

$\displaystyle|d(\hat{h}^{*};\lambda)-$
$\displaystyle\inf_{h}d(h;\lambda)|\leq$ $\displaystyle 4{\cal R}({\cal
L})+2F_{{\bf
V},\psi,B}\sqrt{\dfrac{\log(1/\delta)}{2n}}+2\lambda\left(F_{\nu_{0},\nu_{1}}{\cal
R}({\cal D})+F_{{\bf W},{\bf
V},B,\kappa,\psi,l,\nu_{0},\nu_{1}}\sqrt{\dfrac{\log(1/\delta)}{2n}}\right)$

with $1-3\delta$ probability, where
$F_{\nu_{0},\nu_{1}}=2\left(\frac{1}{\nu_{0}}+\frac{1}{1-\nu_{1}}\right)$ and
$F_{{\bf W},{\bf
V},B,\kappa,\psi,l,\nu_{0},\nu_{1}}=\left(\frac{1+\nu_{0}-\nu_{1}}{\nu_{0}(1-\nu_{1})}\right)\times
F_{{\bf W},{\bf V},B,\kappa,\psi,l}$.

###  9.2 Proof of Corollary 1

Let $\hat{h}^{*}=\arg\min\hat{d}(h;\lambda)$, $h^{*}_{0}=\arg\min
d(h;\lambda=0)$, and $h^{*}=\arg\min d(h;\lambda)$. Let’s denote by
$\Delta(h_{0}^{*},h^{*}):=R_{F}(h_{0}^{*})-R_{F}(h^{*})$ the amount of
decision discrimination in the population with respect to the sensitive
information. Since $d(h^{*};\lambda)\leq d(h_{0}^{*};\lambda)$, we obtain

$\displaystyle d(\hat{h}^{*};\lambda=0)-\inf_{h}d(h;\lambda=0)$
$\displaystyle=R(\hat{h}^{*})-R(h_{0}^{*}),$
$\displaystyle=R(\hat{h}^{*})-R(h^{*})+R(h^{*})-R(h_{0}^{*}),$
$\displaystyle\leq R(\hat{h}^{*})-R(h^{*})+\lambda\Delta(h_{0}^{*},h^{*}).$

Similar to $\eqref{thm:bdd1}\sim\eqref{thm:bdd3}$, the loss of utility has the
following upper bound

$\displaystyle|d(\hat{h}^{*};\lambda=0)-\inf_{h}d(h;\lambda=0)|$
$\displaystyle\leq
2\underbrace{\sup_{h}|R(h)-\hat{R}(h)|}_{\text{I}}+\lambda\Delta(h_{0}^{*},h^{*}),$

where I has the same upper bound (9).

##  10 Details of Simulation

This section explains simulation setups and delivers additional empirical
studies to justify the performance of our method. The outline of this section
is as follows.

  * •

10.1 explains the simulation settings concretely;

  * •

10.2 has all other simulation results omitted in the main text such as
simulation results for (generalized) statistical parity in Scenario II and
(generalized) equalized odds in Section I, trade-off figures, and additional
tables summarizing fairness scores;

  * •

10.3 compares our model to HGR and CON in Scenario III;

  * •

10.4 discusses the estimation stability of our method when the density ratio
estimator (for separation) is poor;

  * •

10.5 illustrates our method outperforms in having fair representation.

###  10.1 Overall simulation setting

##### Target objective function

We first clarify that the notation $(1-\lambda)L_{M}+\lambda L_{F}$ in the
manuscript. For the competing methods, we refer to the notation of their
original papers. All methods below share $\hat{R}(h)$ to express an empirical
risk function.

  * •

For Ours,

$\displaystyle L_{M}=\hat{R}(h),\quad L_{F}=\hat{R}_{F}(h;D)$

  * •

CON has the same $\hat{R}(h)$ and $\hat{R}_{F}(h;D)$ but it uses
$A^{\prime}\sim P(A|Y)$.

  * •

For HGR (Lee et al.,, 2022),

$\displaystyle L_{M}$ $\displaystyle=\hat{R}(h),$ $\displaystyle L_{F}$
$\displaystyle={\text{HGR}_{\text{soft}}}(h(X),A\otimes
Y)-{\text{HGR}_{\text{soft}}}(h(X),Y)\quad\text{for EO}$ $\displaystyle L_{F}$
$\displaystyle={\text{HGR}_{\text{soft}}}(h(X),A)\quad\text{for SP}.$

Refer to the original paper to see the exact form of
${\text{HGR}_{\text{soft}}}$.

  * •

For KDE (Cho et al.,, 2020),

$\displaystyle L_{M}$ $\displaystyle=\hat{R}(h),$ $\displaystyle L_{F}$
$\displaystyle={\text{DDP}}\quad\text{for SP},$ $\displaystyle L_{F}$
$\displaystyle={\text{DEO}}\quad\text{for EO},$

Refer to the original paper to see the exact expression of DDP and DEO.

  * •

For NEU (Du et al.,, 2021),

$\displaystyle L_{M}={\cal L}_{\text{MSE}},\quad L_{F}={\cal
L}_{\text{Smooth}}.$

Refer to the original paper to see the exact expression of ${\cal
L}_{\text{MSE}}$ and ${\cal L}_{\text{Smooth}}$. It is assumed that the
sensitive information is available during the training phase. For NEU, the
model $h$ is pre-trained by minimizing $\hat{R}(h)$ with an early-stopping
procedure. Then the discriminatory $h$ is trained to minimize
$(1-\lambda)L_{M}+\lambda L_{F}$ to make it debiased.

##### Network architectures for the simulation studies

All methods have the same neural network model for $h$. More specifically,

  1. A.

Scenario I & II (binary classification with cross-entropy)

     * •

$h$ has [Dense(64)-BN-ReLU]*3-[Dense(1)-Sigmoid].

     * •

$D$ has [Dense(64)-BN-ReLU]*2-[Dense(1)-Sigmoid].

     * •

$D_{\beta}$ has [Dense(64)-ReLU]*2-[Dense(1)-Sigmoid] (for separation).

For $\text{HGR}_{\text{soft}}$, two neural networks are needed, and each of
them is set to [Dense(64)-BN-ReLU]*2-[Dense(1)]. As a result, HGR employs four
auxiliary networks in total for separation and two networks for independence.
For NEU, $h$’s first output layer is neutralized. For CON in Scenario II, the
discriminator has the same structure with $D_{\beta}$ while the generator
$G=[G_{D},G_{C}]$ with $G_{D}=G_{B}$-Dense(1)-Sigmoid and
$G_{C}=G_{B}$-Dense(1) where
$G_{B}=$[Concat(Unif(2,[0,1]),$Y$)-Dense(64)-ReLU-Dense(64)-ReLU].
Unif(2,[0.1]) denotes a 2-dimensional uniform random variable each of which is
in $[0,1]$. The binary value is determined by binning the outputs of $G_{D}$
with the threshold 0.5.

  2. B.

Scenario III (regression with mean absolute error)

     * •

$h$ has [Dense(16)-BN-ReLU]*3-[Dense(1)].

     * •

$D$ has [Dense(16)-BN-ReLU]*2-[Dense(1)-Sigmoid]

     * •

$D_{\beta}$ has [Dense(16)-ReLU]*2-[Dense(1)-Sigmoid] (for separation).

$\text{HGR}_{\text{soft}}$ has two [Dense(16)-BN-ReLU]*2-[Dense(1)]. For CON,
the generator only has a continuous output
$G=$[Concat(Unif(2,[0,1]),$Y$)-Dense(16)-ReLU-Dense(16)-Dense(1)].

##### Hyperparameter

All models experimented 5 times for each $\lambda\in\\{0.1,0.3,0.5,0.7,0.9\\}$
based on 80% training and 20% validation set. Seed numbers are specified such
that all models in the comparison are trained on the same data sets for each
$\lambda$. For optimization, we adopt stochastic gradient descent with a
learning rate of 0.005. For SBP and HGR, the number of training iterations for
the maximization part of both algorithms, e.g., the notation $T^{\prime}$ in
our work (Algorithm 1 and 2), is set to 1. The evaluation metrics are
calculated every 100 iterations. Each metric is basically expressed as a sum,
but it is equivalent to finding an average value. To see other details
including a mini-batch size or the number of epochs for all data sets, please
refer to the shared code scripts.

##### Training time comparison shown in Table 4

NEU, SBP, and CON require pre-training courses, i.e., pre-training of $h$ for
NEU, of $\hat{\beta}$ for SBP, the GAN model for CON. As a possible case of
NEU, pre-training of $h$ may need substantial training time or longer
evaluation epochs for using early stopping. In the case of CON, learning the
generator gets longer as the complexity of $A|Y$ increases. Tables 4 and 6
present the training times including/excluding the pre-training courses during
early iterations for EO and SP to the Adult data. The same pattern is observed
in other data sets, so they are omitted. All experiments, written in
Tensorflow 2.4.0, run on CentOS 7 featuring Nvidia A30 GPU and 192GB of RAM.

Table 6: Training times (mins) for the first 1000 iterations on A30 GPU for independence in Adult. Note NEU needs the pre-training course while others do not, and KDE cannot be applied to Scenario II. | Without pre-training | With pre-training
---|---|---
Method | Scenario I | Scenario II | Scenario II
SBP | 0.23 (0.02) | 0.56 (0.05) | -
HGR | 0.24 (0.02) | 0.59 (0.06) | -
KDE | 0.88 (0.07) | - | -
NEU | 1.01 (0.09) | 1.46 (0.14) | 3.61 (0.73)

###  10.2 More simulation results in Scenario I & II

For the continuous variable, we calculate

SP $\displaystyle=~{}|{\cal A}^{*}|^{-1}\sum_{a\in{\cal A}^{*}}\left|{\bf
E}(\hat{Y}=1|A\leq a)/{\bf E}(\hat{Y}=1)-1\right|,$ KS-GSP
$\displaystyle=~{}|{\cal A}^{*}|^{-1}\sum_{a\in{\cal
A}^{*}}\max_{h_{x}}|\hat{P}(\hat{h}^{*}(X)\leq h_{x}|A\leq
a)-\hat{P}(\hat{h}^{*}(X)\leq h_{x})|.$

where ${\cal A}^{*}=\\{\tilde{q}_{10},\dots,\tilde{q}_{90}\\}$ with the
$\tilde{q}_{r}:=\frac{r}{100}$th quantile of $A$. For EO and KS-GEO,

EO $\displaystyle=|{\cal A}^{*}|^{-1}\sum_{y\in{\cal Y},a\in{\cal
A}^{*}}\left|\frac{{\bf E}(\hat{Y}=1|A\leq a,Y=y)}{{\bf
E}(\hat{Y}=1|Y=y)}-1\right|,$ KS-GEO $\displaystyle=|{\cal
A}^{*}|^{-1}\sum_{y\in{\cal Y},a\in{\cal
A}^{*}}\max_{h_{x}}|\hat{P}(\hat{h}^{*}(X)\leq h_{x}|A\leq
a,Y=y)-\hat{P}(\hat{h}^{*}(X)\leq h_{x}|Y=y)|.$

Following the same manner in the manuscript, we collect the Pareto frontiers
for all 5 experiments for each $\lambda$. It is possible that each experiment
produces a set of Pareto solutions. Figures 2, 3, and 4 overlay those Pareto
frontiers for all competitors and fairness measures. In Figure 2, we observe
that CON works as well as SBP in Scenario I. This is because $P(A|Y)$ can be
easily estimated by the maximum likelihood method.

Figure 2: (Scenario I) Pareto frontiers: the first row includes pairs of EO
and AUC, and the second row shows pairs of KS-GEO and AUC from 5 experiments
for each $\lambda$. Since SBP and CON are similar, they are directly compared
in the below figure. SBP and CON illustrate better results than others in
Adult and Law School Admission as they are tightly in the upper-left corner.
Figure 3: (Scenario II) Pareto frontiers: the first and the second row
correspond to Adult, Credit Card Default, and ACSEmployment respectively from
5 experiments for each $\lambda$. SBP is superior to both NEU and HGR overall
but comparable to HGR in Adult. Note CON cannot handle statistical parity.
Figure 4: (Scenario II) Pareto frontiers: the first and the second row
correspond to Adult, Credit Card Default, and ACSEmployment respectively from
5 experiments for each $\lambda$. Remarkably, SBP (ours) outperforms the
competitors for the most part. Note KDE cannot handle continuous attributes.

####  10.2.1 Trade-off curves by differing $\lambda$

Figure 5 shows that our method can capture the trade-off between utility and
fairness in independence as $\lambda$ differs. For clear visualization,
$\lambda$ is specifically chosen as 0.1, 0.5, and 0.5, respectively. The
trade-off is also shown in separation, so it is omitted.

Figure 5: (Scenario I of SBP) Pareto frontiers of SP (and KS-GSP) and AUC by
differing $\lambda$ from 0.1 to 0.9 for all 5 experiments.

####  10.2.2 Evaluation tables with other AUC thresholds

  * •

(Scenario I-SP) Table 7 uses other thresholds compared to Table 2.

  * •

(Scenario I-EO) Tables 8 and 9 show the tables with different thresholds.

  * •

(Scenario II-SP) Table 10 show the results in SP.

  * •

(Scenario II-EO) Table 11 uses other thresholds compared to Table 3.

These extra tables show consistent results with the tables in the main text
and with their corresponding figures. In general, our method outperforms the
competing methods. The thresholds are chosen such that there are a sufficient
number of fairness scores to calculate their averages and standard deviations
from the 5 smallest scores. For instance, if the larger threshold is chosen
such as AUC $\geq 0.9$ in Adult, there are no available fairness scores in all
comparison methods.

Table 7: (Scenario I) Averages of the 5 smallest of SP/KS-GSPs whose AUCs are greater than the thresholds. Those scores are selected in the Pareto solutions appearing in Figure 1. Standard deviations are in the parentheses next to the averages. Note CON is not available for SP. |  Adult (AUC $\geq$ 0.80) |  Cred. Card. (AUC $\geq$ 0.70) |  Law School. (AUC $\geq$ 0.85)
---|---|---|---
| SP ($\downarrow$) | KS-GSP ($\downarrow$) | SP ($\downarrow$) | KS-GSP ($\downarrow$) | SP ($\downarrow$) | KS-GSP ($\downarrow$)
SBP | 0.001 ($\approx 0$) | 0.014 ($\approx 0$) |  $\approx 0$ ($\approx 0$) | 0.011 (0.001) | 0.199 (0.032) | 0.112 (0.011)
HGR | 0.070 (0.005) | 0.040 ($\approx 0$) | 0.001 (0.001) | 0.020 (0.001) | 0.123 (0.005) | 0.089 (0.001)
KDE | 0.280 (0.018) | 0.110 (0.005) | 0.005 (0.003) | 0.023 (0.001) | 0.252 (0.006) | 0.151 (0.001)
NEU | 0.372 (0.018) | 0.137 (0.002) | 0.002 (0.001) | 0.013 (0.001) | 0.385 (0.008) | 0.159 (0.006)
Table 8: (Scenario I) Averages of the 5 smallest of EO/KS-GEOs whose AUCs are greater than the thresholds. Those scores are selected in the Pareto solutions appearing in Figure 2. Standard deviations are in the parentheses next to the averages. |  Adult (AUC $\geq$ 0.85) |  Cred. Card. (AUC $\geq$ 0.75) |  Law School. (AUC $\geq$ 0.85)
---|---|---|---
| EO ($\downarrow$) | KS-GEO ($\downarrow$) | EO ($\downarrow$) | KS-GEO ($\downarrow$) | EO ($\downarrow$) | KS-GEO ($\downarrow$)
SBP | 0.009 (0.005) | 0.070 (0.002) | 0.003 (0.001) | 0.040 (0.001) | 0.166 (0.027) | 0.140 (0.005)
CON | 0.019 (0.008) | 0.066 (0.002) | 0.001 (0.001) | 0.038 (0.002) | 0.217 (0.033) | 0.190 (0.009)
HGR | 0.164 (0.031) | 0.170 (0.006) | 0.002 (0.001) | 0.047 (0.002) | 0.891 (0.080) | 0.409 (0.014)
KDE | 0.383 (0.013) | 0.172 (0.003) | 0.002 (0.001) | 0.050 (0.001) | 0.311 (0.128) | 0.262 (0.006)
NEU | 0.450 (0.028) | 0.205 (0.005) | 0.030 (0.019) | 0.046 (0.001) | 0.421 (0.022) | 0.299 (0.003)
Table 9: (Scenario I) Averages of the 5 smallest of EO/KS-GEOs whose AUCs are greater than the thresholds. Those scores are selected in the Pareto solutions appearing in Figure 2. Standard deviations are in the parentheses next to the averages. |  Adult (AUC $\geq$ 0.80) |  Cred. Card. (AUC $\geq$ 0.70) |  Law School. (AUC $\geq$ 0.80)
---|---|---|---
| EO ($\downarrow$) | KS-GEO ($\downarrow$) | EO ($\downarrow$) | KS-GEO ($\downarrow$) | EO ($\downarrow$) | KS-GEO ($\downarrow$)
SBP | 0.009 (0.005) | 0.070 (0.003) | 0.003 (0.001) | 0.040 (0.001) | 0.020 (0.010) | 0.078 (0.001)
CON | 0.016 (0.007) | 0.066 (0.002) | 0.001 (0.001) | 0.038 (0.002) | 0.014 (0.003) | 0.074 (0.002)
HGR | 0.127 (0.025) | 0.167 (0.007) | 0.002 ($\approx 0$) | 0.047 (0.002) | 0.587 (0.152) | 0.365 (0.027)
KDE | 0.282 (0.064) | 0.155 (0.009) | 0.001 (0.001) | 0.049 (0.001) | 0.194 (0.016) | 0.216 (0.007)
NEU | 0.364 (0.013) | 0.182 (0.003) | 0.016 (0.003) | 0.046 (0.001) | 0.380 (0.033) | 0.245 (0.012)
Table 10: (Scenario II) Averages of the 5 smallest of SP/KS-GSPs whose AUCs are greater than the thresholds. Those scores are selected by referring to Figure 3. Standard deviations are in the parentheses next to the averages. |  |  Adult (AUC $\geq$ 0.80) |  Cred. Card. (AUC $\geq$ 0.75) |  ACSEmpl. (AUC $\geq$ 0.75)
---|---|---|---|---
Metric | Method | Race ($\downarrow$) | Age ($\downarrow$) | Gender ($\downarrow$) | Age ($\downarrow$) | Gender ($\downarrow$) | Age ($\downarrow$)
SP | SBP | 0.002 (0.001) | 0.018 (0.002) |  $\approx 0$ ($\approx 0$) | 0.004 (0.001) | 0.001 (0.001) | 0.116 (0.002)
HGR | 0.148 (0.006) | 0.019 (0.001) |  $\approx 0$ ($\approx 0$) | 0.010 (0.001) |  $\approx 0$ ($\approx 0$) | 0.196 (0.004)
NEU | 0.338 (0.018) | 0.216 (0.010) | 0.050 (0.011) | 0.040 (0.004) |  $\approx 0$ ($\approx 0$) | 0.216 (0.008)
KS- GSP | SBP | 0.019 (0.002) | 0.015 ($\approx 0$) | 0.014 ($\approx 0$) | 0.010 ($\approx 0$) | 0.019 (0.003) | 0.078 (0.002)
HGR | 0.062 (0.005) | 0.021 (0.001) | 0.015 ($\approx 0$) | 0.012 ($\approx 0$) | 0.017 (0.002) | 0.118 (0.001)
NEU | 0.127 (0.001) | 0.088 (0.003) | 0.025 (0.003) | 0.026 ($\approx 0$) | 0.012 ($\approx 0$) | 0.140 (0.001)
Table 11: (Scenario II) Averages of the 5 smallest of EO/KS-GEOs whose AUCs are greater than the thresholds. Those scores are selected by referring to Figure 4. Standard deviations are in the parentheses next to the averages. |  |  Adult (AUC $\geq$ 0.85) |  Cred. Card. (AUC $\geq$ 0.70) |  ACSEmpl. (AUC $\geq$ 0.80)
---|---|---|---|---
Metric | Method | Race ($\downarrow$) | Age ($\downarrow$) | Gender ($\downarrow$) | Age ($\downarrow$) | Gender ($\downarrow$) | Age ($\downarrow$)
EO | SBP | 0.143 (0.013) | 0.122 (0.002) | 0.001 ($\approx 0$) | 0.016 (0.001) | 0.168 (0.004) | 0.373 (0.006)
CON | 0.425 (0.017) | 0.241 (0.002) | 0.006 (0.003) | 0.016 (0.003) | 0.224 (0.004) | 0.458 (0.031)
HGR | 0.355 (0.010) | 0.307 (0.002) | 0.002 (0.001) | 0.033 (0.001) | 0.195 (0.004) | 0.447 (0.002)
NEU | 0.385 (0.022) | 0.268 (0.012) | 0.015 (0.010) | 0.038 (0.004) | 0.143 (0.008) | 0.443 (0.017)
KS- GEO | SBP | 0.147 (0.002) | 0.124 (0.002) | 0.040 (0.001) | 0.036 (0.001) | 0.094 (0.001) | 0.247 (0.002)
CON | 0.177 (0.001) | 0.138 (0.002) | 0.049 (0.001) | 0.038 ($\approx 0$) | 0.115 (0.003) | 0.312 (0.002)
HGR | 0.162 (0.004) | 0.221 (0.001) | 0.041 (0.003) | 0.047 (0.001) | 0.122 (0.001) | 0.356 (0.005)
NEU | 0.196 (0.008) | 0.166 (0.004) | 0.053 (0.003) | 0.051 (0.001) | 0.101 (0.006) | 0.334 (0.009)

###  10.3 Simulation results (for both independence and separation) in
Scenario III

We further compare SBP with HGR and CON in Scenario III using Community and
Crime1 data. Following Lee et al., (2022), the number of violent crimes per
population and the ratio of black people in the population are used as a
continuous outcome and a sensitive variable respectively. For the evaluation
of fairness, we consider the following metrics,

SP $\displaystyle=|{\cal A}^{*}|^{-1}\sum_{a\in{\cal A}^{*}}\left|{\bf
E}(\hat{h}^{*}(X)|A\leq a)/{\bf E}(\hat{h}^{*}(X))-1\right|,$ EO
$\displaystyle={|{\cal Y}^{*}|}^{-1}\sum_{y\in{\cal Y}^{*}}{|{\cal
A}^{*}|}^{-1}\sum_{a\in{\cal A}^{*}}\left|{\bf E}(\hat{h}^{*}(X)|A\leq a,Y\leq
y)/{\bf E}(\hat{h}^{*}(X)|Y\leq y)-1\right|,$ KS-GSP $\displaystyle=|{\cal
A}^{*}|^{-1}\sum_{a\in{\cal A}^{*}}\max_{h_{x}}|\hat{P}(\hat{h}^{*}(X)\leq
h_{x}|A\leq a)-\hat{P}(\hat{h}^{*}(X)\leq h_{x})|,$ KS-GEO
$\displaystyle={|{\cal Y}^{*}|}^{-1}\sum_{y\in{\cal Y}^{*}}{|{\cal
A}^{*}|}^{-1}\sum_{a\in{\cal A}^{*}}\max_{h_{x}}|\hat{P}(\hat{h}^{*}(X)\leq
h_{x}|A\leq a,Y\leq y)-\hat{P}(\hat{h}^{*}(X)\leq h_{x}|Y\leq y)|.$

where ${\cal A}^{*}$ and ${\cal Y}^{*}$ are the sets of quantile values of $A$
and $Y$; ${\cal Y}^{*}=\\{\tilde{b}_{10},\dots,\tilde{b}_{90}\\}$ with the
$\tilde{b}_{r}:=\frac{r}{100}$th quantile of $Y$.

In Scenario III, NEU cannot be implemented because NEU requires partitioning a
data set with respect to the outcome variable, which is practically impossible
for the continuous outcome. Figure 6 and Table 12 show that SBP and HGR are
comparable in Scenario III.

Figure 6: (Scenario III) Pareto frontiers: the points are pairs of the fairness metrics and MAE. Note the bottom-left tendency implies a better trade-off.  Table 12: (Scenario III) Averages of the 5 smallest of SP/KS-GSP/EO/KS-GEOs whose MAEs are less than the thresholds. Those scores are selected in the Pareto solutions appearing in Figure 6. Standard deviations are in the parentheses next to the averages. |  Community and Crime (MAE $\leq$ 0.12)
---|---
| SP ($\downarrow$) | KS-GSP ($\downarrow$) | EO ($\downarrow$) | KS-GEO ($\downarrow$)
SBP | 0.225 (0.007) | 0.141 (0.004) | 0.083 (0.003) | 0.060 (0.002)
CON | - | - | 0.094 (0.004) | 0.062 (0.002)
HGR | 0.247 (0.002) | 0.155 (0.001) | 0.082 (0.004) | 0.069 (0.002)

###  10.4 Robust estimation against the poor estimate of $\beta(a,y)$

We also investigate the impact of $\hat{\beta}$ discussed in Section 3.2 based
on II. Interestingly, we observe that $\hat{\beta}$ hardly affects the
performance of SBP in spite of the fact that the density-ratio estimator is
required to guarantee GEO theoretically. Figure 7 compares the Pareto
frontiers of SBP in Figure and SBP with $\hat{\beta}(a,y)=1$ for all $a,y$
(SBP_NW). SBP_NW is set to having the same simulation configuration as SBP.
The figure describes that SBP and SBP_NW are almost the same in both
scenarios. The same tendency is also found in Scenario I, so it is omitted.

Figure 7: (Scenario II, the impact of $\hat{\beta}$) Pareto frontiers: the
points are pairs of EO and AUC or KS-GEO and AUC.

###  10.5 Fair Representation

We carry out additional simulation studies to verify that the neural penalty
can be used to produce fair representation in Scenario II for EO. Following
the same notation in Section 3.2, we set $E$ as [Dense(64)-BN-ReLU] and
$h_{E}$ as [Dense(64)-BN-ReLU]*2-[Dense(1)-Sigmoid], so that $h=h_{E}\circ E$.
This is the same architecture as NEU. For HGR, we calculate
$\text{HGR}_{\text{soft}}(E(X),A)$ for SP and
${\text{HGR}_{\text{soft}}}(E(X),A\otimes
Y)-{\text{HGR}_{\text{soft}}}(E(X),Y)$ for EO. When training SBP and HGR, the
number of training iterations for the maximization part is set to 5, i.e.,
$T^{\prime}=5$.

Figure 8 highlights that SBP succeeds in generating fair representation as
achieving comparable or even better performance than the competing methods.
The collected Pareto frontiers of the models on the benchmark data sets
illustrate that SBP tends to defeat others slightly. We see that the SBP can
represent the fair representation for SP as well, but it is omitted to save
the page.

Figure 8: (Scenario II, fair representation) Pareto frontiers: the points are
pairs of EO and AUC or KS-GEO and AUC.