# State-Visitation Fairness in Average-Reward MDPs

Ganesh Ghalme111Author names appear alphabetically. Vineet Nair∗ Technion
Israel Institute of Technology Technion Israel Institute of Technology
<EMAIL_ADDRESS><EMAIL_ADDRESS>Vishakha Patil∗ Yilun
Zhou∗ Indian Institute of Science Massachusetts Institute of Technology
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Fairness has emerged as an important concern in automated decision-making in
recent years, especially when these decisions affect human welfare. In this
work, we study fairness in temporally extended decision-making settings,
specifically those formulated as Markov Decision Processes (MDPs). Our
proposed notion of fairness ensures that each state’s long-term visitation
frequency is more than a specified fraction. In an _average_ -reward MDP
(AMDP) setting, we formulate the problem as a bilinear saddle point program
and, for a generative model, solve it using a Stochastic Mirror Descent (SMD)
based algorithm. The proposed solution guarantees a _simultaneous
approximation_ on the expected average-reward and the long-term state-
visitation frequency. We validate our theoretical results with experiments on
synthetic data.

##  1 Introduction

Algorithms are being increasingly used to make important decisions that impact
human lives. While algorithmic decision-making frameworks offer increased
efficiency, speed, and scalability, the potential bias/unfairness that may be
inherently induced by them when taking high-stakes decisions has raised
several concerns. For instance, recidivism prediction [DF18], loan and credit
lending [BCGH96], online advertising [ASB+19], and recommendation systems
[YH17] are a few examples where studies have shown that the traditional
algorithms may be unfair towards certain demographics of the population. These
concerns have led to a surge in the research efforts aimed at ensuring
fairness in algorithmic decision-making frameworks [BHN18].

Machine learning-based algorithms are at the core of many of these algorithmic
decision-making frameworks. A large body of work in fair machine learning has
focused on either one-shot settings such as classification [KMR17, DHP+12], or
sequential but static settings such as multi-armed bandits where reward
distributions are stationary [CKSV19, ZL20] (also see Section 7). However, in
many real-world applications, the algorithm’s decisions may affect the future
states and influence future rewards of the system. The study of fairness in
such temporally-extended decision-making settings, which can be modeled using
the _reinforcement learning_ framework [SB18], is still in its infancy.

We introduce the problem of state-visitation fairness in Markov Decision
Processes (MDPs). Informally, our fairness notion requires that each state of
the MDP be visited some minimum (pre-specified) fraction of times. In
particular, a policy $\pi$ with stationary distribution
$\boldsymbol{\nu}^{\pi}$, is called _fair_ if $\nu^{\pi}_{s}\geq\rho_{s}$ for
every state $s$, where $\rho_{s}\in[0,1)$ specifies the fairness constraints
and is given as input to the algorithm. We note that such a _quota-based_
notion of fairness is natural in resource allocation settings where an
algorithm is required to provide a minimum quantum of the resource being
allocated to an entity (captured in our case by the state-visitation
frequency). Our notion of fairness is adapted from the recent works on
fairness in the Multi-Armed Bandits (MAB) [LLJ19, PGNN20]. Since MABs are a
restricted class of MDPs222The MAB problem is equivalent to a single-state
MDP., extending the notions of fairness from the the MAB setting to the more
general MDP setting is both interesting and challenging. Our work makes an
important contribution in this direction. Next, we motivate our fairness
notion with an example.

Consider a taxi service with multiple taxi stations spread across the city. A
passenger can only be picked up from one of the stations but can be dropped
off anywhere in the city. After the drop-off, the taxi driver decides which
station to visit next to pick up the next customer. The fare/tip received by
the driver when a passenger is picked up from a particular station corresponds
to the reward obtained in that state (station). This reward may depend on
factors such as proximity of the station to tourist spots and types of
businesses in the vicinity. In such a domain, a recommendation policy would
like to maximize the average-reward while ensuring accessibility of the
service at each station. In this context, our proposed algorithm computes a
policy that ensures – under mild feasibility assumptions – a minimum (pre-
specified) state-visitation guarantee for each state under a stationary
distribution. That is, the customers near each station are ensured certain
minimum access to the taxi service.

We study the problem of maximizing the average-reward while ensuring state-
visitation fairness. In many practical resource-allocation problems, where our
fairness notion is natural, the long-term rewards are as important as the
immediate ones. For example, a taxi driver in the above example is likely to
prefer the average-reward maximizing policy over one that maximizes the
discounted reward.

###  1.1 Main results

We formulate the problem of computing an average-reward maximizing policy
under fair state-visitation constraint as a linear program. Our main
contribution is a quantifiable guarantee on fairness and optimality of the
policy returned in a generative model333Under the assumption that a strictly
feasible solution exists.. In particular, we show that the policy returned by
our proposed algorithm (Algorithm 1) achieves a simultaneous approximation
guarantee on the average-reward and the state-visitation frequency of every
state. Algorithm 1 uses the classical Stochastic Mirror Descent (SMD)
framework to achieve the desired approximation guarantee (see Theorems 1 and
2). The running time of our algorithm depends on the desired approximation
guarantee, the mixing time of the MDP, and the fairness constraint.

We remark here that [JS20] also give a similar SMD based algorithm to compute
an approximate solution for the constrained AMDP problem _without_ any
objective function (see Section 5 of [JS20]). The problem of computing a fair
and optimal policy is an instance of the constrained AMDP problem with the
linear objective function. The constraint matrix in our case is diagonal, and
this structure enables us to achieve simultaneous guarantees on feasibility
and optimality. In our case, the structure of the constraint matrix is used to
determine the appropriate domain bounds for the primal and dual variables (see
Lemma 1). To the best of our knowledge, this work is the first attempt at
studying a constrained AMDP with a quantifiable average-reward guarantee and
constraint satisfaction guarantee in a generative model setup.

###  1.2 Related work

The study of fairness in reinforcement learning (RL) was initiated in
[JJK+17], where they extend the _meritocratic fairness_ notion in [JKMR16]
from the MAB setting to the RL setting and call a policy fair if, with high
probability, an action with a lower long-term reward is not favored over an
action with a higher long-term reward. In another work, [WBT19] modeled the
agents as states of an MDP and study demographic parity with respect to reward
to the agents. The problem of off-policy policy selection in RL with fairness
constraints, which are similar in spirit to the constraints in proposed in
[JJK+17], is studied in [DTB18]. In [ZTL+20], how algorithmic decisions impact
the evolution of feature space of the underlying population when the decision-
making process is modeled as an MDP is studied. Finally, in [CMPZ20] causal
modeling of dynamical systems to address fairness considerations is performed.

The AMDP problem has been extensively studied in the literature [Mah96, KS98,
BT02]. If the model is known and the MDP is unichain, then in [Alt99] and in
[Put14] it was shown that an optimal policy can be computed by solving a
linear program. For a generative model, in [Wan17] and [JS20] the SMD approach
[NJLS09, CJST19] was recently used to compute an approximately optimal policy.
Additionally, in [JS20] an algorithm to compute a feasible policy for
constrained AMDPs is also provided. Additional works on fairness in machine
learning, constrained MDPs, and safe-RL are discussed in Section 7.

###  1.3 Organization

In Section 2, we introduce the Fair-AMDP problem. In Section 3, we formulate
the problem as a linear program and provide the necessary background on the
SMD framework, which is used to give an algorithm for the generative model.
The proposed algorithm and its corresponding theoretical guarantees are given
in Sections 4 and 5, respectively. We experimentally validate our proposed
algorithm in Section 6. Finally, we conclude our paper with a discussion on
our work and possible future directions in Section 8.

Throughout this work, we use $[n]$ (where $n\in\mathbb{N}$) to denote the set
$\\{1,2,\ldots,n\\}$. We use $\Delta^{n}$ to denote the $n-1$ dimensional
probability simplex, i.e.,
$\Delta^{n}=\\{(x_{1},x_{2},\ldots,x_{n})\in\mathbb{R}^{n}~{}|~{}\sum_{i\in[n]}x_{i}=1\text{
and }x_{i}\geq 0~{}\forall i\in[n]\\}$.

##  2 Fair-AMDP Model

A discrete Markov Decision Process (MDP) is a sequential decision-making
framework denoted by the tuple
$\langle\mathcal{S},\mathcal{A},\mathbf{\Gamma},\boldsymbol{r},\boldsymbol{u}\rangle$.
At each discrete time $t\geq 1$, $s_{t}\in\mathcal{S}$ denotes the state of
the MDP at time $t$. A decision-maker takes an action $a_{t}$ from the action-
set444In general, the action set associated with each state can be different
and our results can be extended to this setting. $\mathcal{A}$, receives a
finite non-negative reward $r_{t}=r_{s_{t},a_{t}}$, and the MDP transitions to
a state $s_{t+1}$.555In general, the reward received at time $t$ can depend on
$s_{t},s_{t+1}$ and $a_{t}$ but we may assume without loss of generality that
rewards depend only on $s_{t}$ and $a_{t}$ . The transition to the next state
$s_{t+1}$ is governed by the transition probability function
$\mathbf{\Gamma}:\mathcal{S}\times\mathcal{A}\rightarrow\Delta^{|\mathcal{S}|}$
with the initial state $s_{1}$ sampled from the initial-state distribution
$\boldsymbol{u}\in\Delta^{|\mathcal{S}|}$. Throughout the paper, we consider a
finite MDP with $\mathcal{S}=\\{1,2,\ldots,n\\}$ and
$\mathcal{A}=\\{1,2,\ldots,m\\}$. For convenience, let $\ell=nm$ denote the
total number of state-action pairs. The entries of the vectors in
$\mathbb{R}^{\ell}$ are indexed by state-action pairs $(s,a)$, where $(s,a)$
indexes the $((s-1)m+a)$-th entry of $[\ell]$. Our main algorithm (Algorithm
1) works on a _generative model_ , where the agent has access to an oracle
that, given a state-action pair $(s,a)$, samples the next state $s^{\prime}$
according to the transition probability function $\mathbf{\Gamma}$ and returns
$s^{\prime}$ and the corresponding reward $r_{s,a}$. Following standard
convention, we use bold lower-case letters for vectors and bold upper-case
letters for matrices. For convenience, we provide a table explaining the
notations used in this work in Appendix A

Let $\pi$ be a stochastic policy with $\pi_{s,a}$ denoting the probability
with which action $a$ is taken in state $s$. Each policy $\pi$ induces a
stationary distribution over the state space denoted by
$\boldsymbol{\nu}^{\pi}\in\Delta^{n}$. One of the popular optimization
problems under an MDP framework [Alt99, Put14] is to find a policy that
maximizes the long-term expected average-reward given by

$\lim_{T\rightarrow\infty}\frac{1}{T}\cdot\sum_{t=1}^{T}\mathbb{E}_{\pi}\left[r_{t}\right].$
(1)

It is known that under ergodicity condition – where the Markov chain induced
by any policy is ergodic (Assumption 1) – the long-term state-action
visitation frequency is well-defined and given as $\nu_{s}^{\pi}\pi_{s,a}.$
The reader is referred to Chapter 8 of [Put14] for additional details.

######  Assumption 1.

A given MDP instance
$\langle\mathcal{S},\mathcal{A},\mathbf{\Gamma},\boldsymbol{r},\boldsymbol{u}\rangle$
satisfies the ergodicity condition666 Note that, even under a weaker unichain
condition (where each induced Markov chain has a single recurrent class and a
few transient states) the long-term state-action visitation frequency is well-
defined [Alt99, LP17]. .

The ergodicity condition enables us to formulate the problem of finding an
optimal policy for an AMDP problem as a linear program (Section 3.1). We
consider a constrained linear program where the constraints are in terms of
the minimum state-visitation frequency. In particular, we study the Fair-AMDP
problem with the following notion of fairness.

######  Definition 1.

Let $\boldsymbol{\rho}\in[0,1]^{n}$ such that $\sum_{s}\rho_{s}\leq 1$. Then,
a policy $\pi$ is called $\boldsymbol{\rho}$-fair if
$\nu_{s}^{\pi}\geq\rho_{s}$ for all $s\in[n]$.

A Fair-AMDP instance is denoted as
$\langle\mathcal{S},\mathcal{A},\mathbf{\Gamma},\boldsymbol{r},\boldsymbol{u},\boldsymbol{\rho}\rangle$.
We note that a Fair-AMDP instance may not even have a feasible policy. For
example, consider an AMDP instance with $n=2,m=1$ and
$\mathbf{\Gamma}=(1-\alpha,\alpha;1-\alpha,\alpha)$. Here, the stationary
distribution for the unique policy $\pi$ is
$\boldsymbol{\nu}^{\pi}=(1-\alpha,\alpha)$. Observe that, in this simple
example there does not exist a $\boldsymbol{\rho}$-fair policy for any
$\boldsymbol{\rho}$ with $\rho_{2}>\alpha$.

The objective of this paper is to find an optimal $\boldsymbol{\rho}$-fair
policy. Unfortunately, the present techniques used to compute even an
approximately feasible solution for constrained AMDPs in the generative model
require prior knowledge that the AMDP has a strictly feasible solution for the
specified constraints. For a Fair-AMDP instance with $\boldsymbol{\rho}$ such
that $\rho_{s}<1/n$ for all $s$, the strict feasibility can be ensured by
assuming the existence of a costly but fair-action $a^{*}$ in all the states
of the AMDP. The following assumption lets us state our results without
explicitly referrring to strict feasibility requirement for a large class of
$\boldsymbol{\rho}$’s.

######  Assumption 2.

There is a zero-reward fair-action $a^{*}$ available at each state that makes
the problem strictly feasible. One choice of $a^{*}$ is such that
$\mathbf{\Gamma}((s,a^{\star}),s^{\prime})=1/n$.

We restrict our attention to $\boldsymbol{\rho}$ such that $\rho_{s}<1/n$ for
all $s\in\mathcal{S}$ unless otherwise specified. We make the following
additional remarks regarding the above assumption.

i) Given a Fair-AMDP with action $a^{*}$, a policy $\pi$ that chooses action
$a^{*}$ in all the states $s$ has $\boldsymbol{\nu}^{\pi}=(1/n\ldots 1/n)$ and
hence $\boldsymbol{\nu}^{\pi}>\boldsymbol{\rho}$ implying that the Fair-AMDP
is strictly feasible. Further, if the Fair-AMDP instance is guaranteed to have
a feasible solution without using the fair-action then there is an optimal
policy that has zero probability of choosing action $a^{*}$ at any state $s$.
In the taxi-service example from Section 1, a fair-action can be an action
that uniformly distributes a taxi across all stations.

ii) For a particular $\boldsymbol{\rho}$, one may relax the assumption to the
following: the fair-action $a^{*}$ is such that
$\mathbf{\Gamma}((s,a^{*}),s^{\prime})>\max_{s}\rho_{s}$ for all
$s,s^{\prime}\in[n]$. Such an assumption is seemingly necessary to compute a
fair-algorithm in the generative model where there is only stochastic access
to the transition probability matrix via state-action queries.

##  3 Solution Approach

In this section, we summarize the approach we use to solve the Fair-AMDP
problem. In Section 3.1, we give a linear program to solve the problem. In
Section 3.2, we state the corresponding bilinear saddle point problem using
the Lagrangian of the linear program defined in Section 3.1. In Section 3.3,
we give a brief overview of the stochastic mirror descent approach used to
approximately solve the bilinear saddle point problem.

###  3.1 Linear Program for Fair-AMDP

Before delving into the linear programming approach, we define a few notations
for an easier exposition. Henceforth, the transition probability function
$\mathbf{\Gamma}$ is an $\ell\times n$ matrix whose $((s,a),s^{\prime})$-th
entry is the probability of transitioning from state $s$ to state $s^{\prime}$
when an action $a$ is taken in the current state $s$. Corresponding to a
policy $\pi$, define the policy matrix $\Pi\in\mathbb{R}^{\ell\times n}$: the
$((s,a),s^{\prime})$-th entry of $\Pi$ is equal to $\pi_{s,a}$ if
$s=s^{\prime}$ and $0$ otherwise. Define
$\mathbf{\Gamma}^{\pi}:=\Pi^{T}\mathbf{\Gamma}$ to be the transition
probability matrix of the Markov chain induced by $\pi$. Let
$\hat{\mathbf{I}}\in\mathbb{R}^{\ell\times n}$ be such that
$\hat{\mathbf{I}}((s,a),s^{\prime})=1$ if $s=s^{\prime}$ and $0$ otherwise.
The famed Bellman-equation for the (unconstrained) AMDP yields the following
primal/dual linear program [Put14].

UC-LP (D)
$\displaystyle\underset{\boldsymbol{x}\in\Delta^{\ell}}{\max}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}\boldsymbol{x}^{T}\boldsymbol{r}$
$\displaystyle\text{subject to
}~{}~{}(\hat{\mathbf{I}}-\mathbf{\Gamma})^{T}\boldsymbol{x}=\boldsymbol{0}$
UC-LP (P)
$\displaystyle\min_{\boldsymbol{\lambda}>0,\beta}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}\beta$
$\displaystyle\text{subject to
}~{}~{}(\hat{\mathbf{I}}-\mathbf{\Gamma})\boldsymbol{\lambda}+\boldsymbol{r}\leq\beta\cdot\boldsymbol{1}$

We focus our attention on UC-LP (D). Corresponding to every feasible
$\boldsymbol{x}$ of UC-LP (D), there is a feasible policy $\pi$ such that, for
each state-action pair $(s,a)$,
$\pi_{s,a}=x_{s,a}/\sum_{a^{\prime}\in[m]}x_{s,a^{\prime}}$. It is easy to
verify the following: a) For
$\boldsymbol{\nu}^{\pi}=(\sum_{a}x_{1,a},\ldots,\sum_{a}x_{n,a})$, we have
$\boldsymbol{x}=\Pi\boldsymbol{\nu}^{\pi}$, b) $\boldsymbol{\nu}^{\pi}$ is the
stationary distribution corresponding to the constructed policy $\pi$, as
$(\hat{\mathbf{I}}-\mathbf{\Gamma})^{T}\boldsymbol{x}=(\mathbf{I}-\mathbf{\Gamma}^{\pi})^{T}\boldsymbol{\nu}^{\pi}=\boldsymbol{0}$,
and c) $\hat{\mathbf{I}}\boldsymbol{x}=\boldsymbol{\nu}^{\pi}$. In particular,
from (c) we conclude that $x_{s,a}$ is the average state-action visitation
probability of the state-action pair $(s,a)$. Here, it is important to note
that Assumption 1 ensures that $\boldsymbol{\nu}^{\pi}>\boldsymbol{0}$, for
every such $\pi$.

From (a), (b), and (c) above, it follows that the desired fairness guarantee
can be achieved by ensuring that $\boldsymbol{x}$ satisfies
$\sum_{a}x_{s,a}/\rho_{s}\geq 1$ for all $s$ such that $\rho_{s}>0$. We
assume, without loss of generality, that $\rho_{s}>0$ for all $s$. We state
the Fair-LP primal/dual below, where $\mathbf{C}\in\mathbb{R}^{n\times\ell}$
is such that the $(s^{\prime},(s,a))$-entry is $1/\rho_{s}$ if $s=s^{\prime}$
and $0$ otherwise.

Fair-LP (D)
$\displaystyle\underset{\boldsymbol{x}\in\Delta^{\ell}}{\max}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}\boldsymbol{x}^{T}\boldsymbol{r}$
$\displaystyle\text{subject to
}~{}~{}(\hat{\mathbf{I}}-\mathbf{\Gamma})^{T}\boldsymbol{x}=\boldsymbol{0},~{}~{}\mathbf{C}\boldsymbol{x}\geq\boldsymbol{1}$
Fair-LP (P)
$\displaystyle\min_{\boldsymbol{\lambda}>\boldsymbol{0},\boldsymbol{\mu}>\boldsymbol{0},\beta}~{}~{}~{}~{}\beta-\sum_{s}\mu_{s}$
$\displaystyle\text{subject to
}~{}~{}(\hat{\mathbf{I}}-\mathbf{\Gamma})\boldsymbol{\lambda}+\boldsymbol{r}+\mathbf{C}^{T}\boldsymbol{\mu}\leq\beta\cdot\boldsymbol{1}$

It is worth noting here that Assumption 2 guarantees the feasibility of Fair-
LP primal and dual. Next, we define an approximate solution to the Fair-AMDP
problem.

######  Definition 2.

Let $\boldsymbol{x}^{*}$ be the optimal solution to Fair-LP (D). Then, we say
that a policy $\pi$ is $(\varepsilon_{1},\varepsilon_{2})$-approximate for the
Fair-AMDP problem if
$\boldsymbol{r}^{T}(\Pi\boldsymbol{\nu}^{\pi})\geq\boldsymbol{r}^{T}\boldsymbol{x}^{*}-\varepsilon_{1}$
and $\mathbf{C}\boldsymbol{x}\geq\boldsymbol{1}-\varepsilon_{2}$.

Proposition 1 follows easily from points (a), (b), and (c) mentioned above.

######  Proposition 1.

Let
$\langle\mathcal{S},\mathcal{A},\mathbf{\Gamma},\boldsymbol{r},\boldsymbol{u},\boldsymbol{\rho}\rangle$
be a Fair-AMDP instance and $\pi$ be a policy with $\boldsymbol{\nu}^{\pi}$ as
the induced stationary distribution. Further, let
$\boldsymbol{x}=\Pi\boldsymbol{\nu}^{\pi}$. Then, $\pi$ is a
$\boldsymbol{\rho}$-fair policy if and only if $\boldsymbol{x}$ is a feasible
solution of Fair-LP.

Proposition 1 establishes that if $\mathbf{\Gamma}$ and $\boldsymbol{r}$ are
known, we can compute an optimal $\boldsymbol{\rho}$-fair policy by solving
Fair-LP. However, when the model parameters are unknown or partially known,
even the problem of verifying the feasibility seems to be difficult. For
instance, if $\alpha$ is not known in our example from Section 2, then without
the fair-action, it is impossible to determine whether a policy is
$\boldsymbol{\rho}$-fair.

Let $(\boldsymbol{\lambda}^{*},\boldsymbol{\mu}^{*},\beta^{*})$ be an
(optimal) solution of Fair-LP (P). First, we note that
$\boldsymbol{\lambda}^{*}$ is not unique, as for every $c\in\mathbb{R}$,
$(\boldsymbol{\lambda}^{*}+c\cdot\boldsymbol{1},\boldsymbol{\mu}^{*},\beta^{*})$
is also an optimal solution. Hence, we may assume without loss of generality
that $\boldsymbol{\lambda}^{*}$ is orthogonal to the stationary distribution
$\boldsymbol{\nu}^{\pi^{*}}$ of the optimal policy $\pi^{*}$ of a given Fair-
AMDP instance. We conclude this section by stating two important lemmas
regarding the nature of $\boldsymbol{\mu}^{*}$ and $\boldsymbol{\lambda}^{*}$.
The presence of a strictly feasible solution is used in Lemma 1, which is in
turn is used to prove Lemma 2. The proofs for Lemmas 1 and 2 are in Section 9.

######  Lemma 1.

Let $\boldsymbol{\lambda}^{*},\boldsymbol{\mu}^{*},\beta^{*}$ be the solution
to the Fair-LP (P). Then $\mu_{s}^{*}\leq\frac{n\rho_{s}}{1-n\rho_{s}}$ for
all $s\in[n]$.

Remark: Suppose we had a priory knowledge that there exists a strictly
feasible policy $\pi$, that is a policy $\pi$ with stationary distribution
$\boldsymbol{\nu}$ such that $\nu_{s}>\frac{1}{\rho_{s}}$. Then it can be
shown that $\mu_{s}^{*}(\frac{\nu_{s}}{\rho_{s}}-1)\leq
1-(\Pi\boldsymbol{\nu})^{T}\boldsymbol{r}$, where
$(\Pi\boldsymbol{\nu})^{T}\boldsymbol{r}$ is the reward of the strictly
feasible policy $\pi$ and is at least $0$.

The bound on $||\boldsymbol{\lambda}^{*}||_{\infty}$ in Lemma 2 is given in
terms of the mixing time of an ergodic MDP as defined below.

######  Definition 3.

(Mixing Time) The mixing time of a given MDP
$\langle\mathcal{S},\mathcal{A},\mathbf{\Gamma},\boldsymbol{r},\boldsymbol{u}\rangle$
is given by $t_{mix}=\max_{\pi}t_{\pi}$ where, $t_{\pi}=\arg\min_{t\geq
1}\left[\max_{\boldsymbol{u}}\left[||(\mathbf{\Gamma}^{\pi^{T}})^{t}\boldsymbol{u}-\boldsymbol{\nu}^{\pi}||_{1}\leq
1/2\right]\right].$

The mixing time $t_{\text{mix}}$ of an ergodic MDP captures how fast the
Markov chain induced by any policy converges to its corresponding stationary
distribution. Note that our proposed algorithm’s convergence time depends on
the mixing time of the MDP via the following lemma.

######  Lemma 2.

Let $\boldsymbol{\lambda}^{*},\boldsymbol{\mu}^{*},\beta^{*}$ be the solution
to the Fair-LP (P). Then $||\boldsymbol{\lambda}^{*}||_{\infty}\leq
M:=2t_{\text{mix}}(1+d_{\boldsymbol{\rho}})$, where
$d_{\boldsymbol{\rho}}=\max_{s}\frac{n}{1-\rho_{s}n}$.

The bound on $||\boldsymbol{\lambda}^{*}||_{\infty}$ crucially restricts the
search space of the primal variable $\boldsymbol{\lambda}$ in the stochastic
mirror descent (SMD) approach (see Sections 3.2 and 3.3) and enables us to
give convergence guarantee for the proposed algorithm.

###  3.2 Minimax Formulation

We formulate Fair-AMDP as a bilinear saddle point problem using the Lagrangian
of Fair-LP,

$h(\boldsymbol{x},\boldsymbol{\lambda},\boldsymbol{\mu},\beta)=\boldsymbol{r}^{T}\boldsymbol{x}+\boldsymbol{\lambda}^{T}(\hat{\mathbf{I}}-\mathbf{\Gamma})^{T}\boldsymbol{x}+\boldsymbol{\mu}^{T}(\mathbf{C}\boldsymbol{x}-\boldsymbol{1})+\beta(1-\boldsymbol{1}^{T}\boldsymbol{x}).$

Let
$\boldsymbol{x}^{*},\boldsymbol{\lambda}^{*},\boldsymbol{\mu}^{*},\beta^{*}$
be a solution to

$\min_{\boldsymbol{\lambda}\in\mathbb{R}^{n},\boldsymbol{\mu}\geq\boldsymbol{0},\beta}~{}~{}\max_{\boldsymbol{x}\geq\boldsymbol{0}}~{}~{}h(\boldsymbol{x},\boldsymbol{\lambda},\boldsymbol{\mu},\beta)~{}.$
(2)

From Assumption 2, we have that Fair-LP problem has a feasible solution.
Hence, using convex optimization theory, we have that $\boldsymbol{x}^{*}$ is
the solution to Fair-LP (D), and
$\boldsymbol{\lambda}^{*},\boldsymbol{\mu}^{*}$ and $\beta^{*}$ is the
solution to Fair-LP (P). Moreover, from the KKT optimality conditions,
$(\boldsymbol{\mu}^{*})^{T}(\mathbf{C}\boldsymbol{x}^{*}-\boldsymbol{1})=0$.

Let
$\Delta_{\boldsymbol{\rho}}^{\ell}:=\\{\boldsymbol{x}\in\Delta^{\ell}\mid\sum_{a}x_{s,a}\geq\rho_{s}~{}~{}\text{for
all }s\\}$, and
$\mathbb{B}_{2M}^{n}:=\\{\boldsymbol{\lambda}~{}\mid~{}||\boldsymbol{\lambda}||_{\infty}\leq
2M\\}$, where $M$ is as defined in Lemma 2. Then, for every
$\boldsymbol{x}^{*}\in\Delta_{\boldsymbol{\rho}}^{\ell}$, note the following:
a) $\beta(1-\boldsymbol{1}^{T}\boldsymbol{x})=0$, b)
$\boldsymbol{0}\in\arg\min_{\boldsymbol{\mu}\geq\boldsymbol{0}}\boldsymbol{\mu}^{T}(\mathbf{C}\boldsymbol{x}-\boldsymbol{1})$.
Hence, using Lemma 2, and points (a) and (b) it follows that
$(\boldsymbol{x}^{*},\boldsymbol{\lambda}^{*},\boldsymbol{\mu}^{*},\beta^{*})$
is a solution to the optimization problem in Equation 2 if and only if
$(\boldsymbol{x}^{*},\boldsymbol{\lambda}^{*})$ is the solution to the
following problem

$\min_{\boldsymbol{\lambda}\in\mathbb{B}_{2M}^{n}}\ \
\max_{\boldsymbol{x}\in\Delta_{\boldsymbol{\rho}}^{\ell}}f(\boldsymbol{x},\boldsymbol{\lambda})=\boldsymbol{r}^{T}\boldsymbol{x}+\boldsymbol{\lambda}^{T}(\hat{\mathbf{I}}-\mathbf{\Gamma})^{T}\boldsymbol{x}$
(3)

Further,
$h(\boldsymbol{x}^{*},\boldsymbol{\lambda}^{*},\boldsymbol{\mu}^{*},\beta^{*})=f(\boldsymbol{x}^{*},\boldsymbol{\lambda}^{*})=\boldsymbol{r}^{T}\boldsymbol{x}^{*}$.
Though $|\boldsymbol{\lambda}^{*}|_{\infty}\leq M$, increasing the domain of
$\boldsymbol{\lambda}$ from $M$ to $2M$ helps in the proof of Theorem 2 which
provides simultaneous approximation guarantees on the optimality and fairness
of the policy that is constructed by Algorithm 1.

Next, we define the Gap function, which quantifies the closeness of a given
feasible solution $(\boldsymbol{x},\boldsymbol{\lambda})$ to the optimal
solution.

######  Definition 4.

(Gap function) The gap function
$\texttt{Gap}:\Delta_{\boldsymbol{\rho}}^{\ell}\times\mathbbm{R}^{n}\rightarrow\mathbbm{R}_{+}$
is defined as,

$\texttt{Gap}(\boldsymbol{x},\boldsymbol{\lambda})=\max_{\boldsymbol{x}^{\prime}\in\Delta_{\boldsymbol{\rho}}^{\ell}}f(\boldsymbol{x}^{\prime},\boldsymbol{\lambda})-\min_{\boldsymbol{\lambda}^{\prime}\in\mathbbm{B}_{2M}^{n}}f(\boldsymbol{x},\boldsymbol{\lambda}^{\prime})$
(4)

It is easy to see that the Gap function is non-negative and
$\texttt{Gap}(\boldsymbol{x}^{\star},\boldsymbol{\lambda}^{\star})=0$. We say
that $(\boldsymbol{x},\boldsymbol{\lambda})$ is an $\varepsilon$-approximate
solution to Equation 3 if
$\texttt{Gap}(\boldsymbol{x},\boldsymbol{\lambda})\leq\varepsilon$.

###  3.3 Stochastic Mirror Descent (SMD) Framework

In this section, we state the SMD framework and then briefly describe how we
use it to compute an (expected) $\varepsilon$-approximate solution to the
optimization problem in Equation 3. This is done using the ghost iterate
technique to compute the approximate solution using only stochastic query
access to $\mathbf{\Gamma}$ and $\boldsymbol{r}$. We begin with a few useful
definitions.

######  Definition 5.

(Strong Convexity) Let $\mathcal{X}\subset\mathbb{R}^{n}$ be a convex set. A
differentiable function $R:\mathcal{X}\rightarrow\mathbbm{R}$ is said to be
$\alpha$-strongly convex with respect to norm $||.||$ if
$R(\boldsymbol{y})\geq R(\boldsymbol{x})+\langle\nabla
R(\boldsymbol{x}),\boldsymbol{y}-\boldsymbol{x}\rangle+\frac{\alpha}{2}||\boldsymbol{y}-\boldsymbol{x}||^{2}$
for all $\boldsymbol{x},\boldsymbol{y}\in\mathcal{X}$.

######  Definition 6.

(Distance Generating Function777A distance-generating function is also called
as a Bregman divergence or a prox-function. ) Let
$R:\mathcal{X}\rightarrow\mathbb{R}$ be a continuously differentiable,
strongly convex, real-valued function on a convex set
$\mathcal{X}\subset\mathbbm{R}^{n}$. For any
$\boldsymbol{x},\boldsymbol{y}\in\mathcal{X}$, the distance from point
$\boldsymbol{x}$ to $\boldsymbol{y}$ is given by
$V_{\boldsymbol{x}}(\boldsymbol{y})=R(\boldsymbol{y})-R(\boldsymbol{x})-\langle\nabla
R(\boldsymbol{x}),\boldsymbol{y}-\boldsymbol{x}\rangle.$

The SMD framework is a stochastic approximation approach of finding a solution
to a stochastic convex program [NJLS09]. In particular, the SMD algorithm is a
special type of Stochastic Gradient Descent (SGD) algorithm, where, in
contrast to SGD, the updates are computed in a _mirrored_ space.

The SMD framework provides an iterative procedure to select points from a
convex space $\mathcal{X}$ with stochastic query access to underlying
parameters. First, a suitable strongly convex regularizer $R$, catering to the
geometry of $\mathcal{X}$ is designed. Then, at each iterate $t$, an unbiased
estimator of gradient $\boldsymbol{g}$ given by $\tilde{\boldsymbol{g}}_{t}$
and a step size $\eta_{t}$ is computed. This estimator, along with the
previous iterate value $\boldsymbol{x}_{t}$ and step size $\eta_{t}$, is used
to compute the next iterate as given below:

$\boldsymbol{x}_{t+1}=\arg\min_{x\in\mathcal{X}}~{}~{}\langle\eta_{t}\tilde{\boldsymbol{g}}_{t},\boldsymbol{x}\rangle+V_{\boldsymbol{x}_{t}}(\boldsymbol{x})~{}.$
(5)

Since the regularizer is strictly convex and differentiable, the SMD update
can be equivalently written in the (mirrored) gradient space as follows

$\nabla R(\boldsymbol{x}_{t+1})=\nabla
R(\boldsymbol{x}_{t})-\eta_{t}\tilde{\boldsymbol{g}}_{t}~{}.$

The strong convexity of the regularizer $R$ implies the uniqueness of the
mapping of $\nabla R(\cdot)$. For more details, the reader is referred to
[NJLS09, CJST19].

In this paper, we use the above SMD framework to obtain an approximate
solution to the optimization problem defined in Equation 3. Note that, in our
problem, there are two convex spaces: $\Delta_{\boldsymbol{\rho}}^{\ell}$
corresponding to $\boldsymbol{x}$, and $\mathbbm{B}_{2M}^{\ell}$ corresponding
to $\boldsymbol{\lambda}$. First, we choose strongly convex regularizers for
the respective spaces and initialize the step sizes $\eta^{\boldsymbol{x}}$
and $\eta^{\boldsymbol{\lambda}}$, respectively888For all iterates the step
sizes remain the same.. Then, at each iterate, we compute bounded estimators
of the gradients $\boldsymbol{g}^{\boldsymbol{x}}$ and
$\boldsymbol{g}^{\boldsymbol{\lambda}}$ given by
$\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$ and
$\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}$, respectively. Here,
$\boldsymbol{g}^{\boldsymbol{x}}$ and $\boldsymbol{g}^{\boldsymbol{\lambda}}$
are the gradients of $f$ with respect to $\boldsymbol{x}$ and
$\boldsymbol{\lambda}$, respectively. Finally, we use the SMD update rule
given in Equation 5 to compute the iterate value for the respective spaces. In
this paper, we consider bounded gradient estimates of the gradients as defined
below.

######  Definition 7.

[JS20] Given the following properties on mean, scale, and variance of an
estimator $\tilde{\boldsymbol{g}}$ of the gradient $\boldsymbol{g}$:
$(i)$ Unbiasedness: $\mathbbm{E}[\tilde{\boldsymbol{g}}]=\boldsymbol{g}$
$(ii)$ Bounded maximum entry: $||\tilde{\boldsymbol{g}}||_{\infty}\leq c$ with
probability $1$
$(iii)$ Bounded second moment $\mathbbm{E}[||\tilde{\boldsymbol{g}}||^{2}]\leq
v$
we call $\tilde{\boldsymbol{g}}$ a $(v,||.||)$-bounded estimator of
$\boldsymbol{g}$ if it satisfies $(i)$ and $(iii)$, and call it a
$(c,v,||.||_{\Delta_{\boldsymbol{\rho}}^{\ell}})$-bounded estimator if it
satisfies $(i)$–$(iii)$ with local norm $||.||_{\boldsymbol{y}}$ for all
$\boldsymbol{y}\in\Delta_{\boldsymbol{\rho}}^{\ell}$.

In Section 4, we present our algorithm and the bounded gradient estimators and
regularizers corresponding to the two convex spaces (mentioned above) used in
the algorithm.

##  4 Algorithm and Estimators

In Section 3, we formulated Fair-AMDP as a bilinear saddle point problem using
the Lagrangian function of Fair-LP. We also showed how the SMD framework can
be used to compute an expected approximate solution to this bilinear problem.
Given input $\varepsilon$, Algorithm 1 in Lines 1–5 computes
$(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})$ which is
an expected $\varepsilon$-approximate solution to the optimization problem in
Equation 3, that is
$\mathbb{E}[\texttt{Gap}(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})]\leq\varepsilon$.
The details of the regularizers and estimators are provided in the next
paragraph. Then, in Line 6, the algorithm uses $\boldsymbol{x}^{\varepsilon}$
to compute $\pi^{\varepsilon}$. We remark that it is not immediately clear
that $\pi^{\varepsilon}$ satisfies the approximation guarantees (in
expectation) as mentioned in Definition 2 even though
$\boldsymbol{x}^{\varepsilon}$ is $\varepsilon$-approximate solution for the
corresponding bilinear saddle point problem. Also, note that Algorithm 1 is
oblivious to the presence of the fair-action.

Now we provide the details regarding the SMD framework used in Lines 1–5. We
use normalized entropic regularizer for $\Delta_{\boldsymbol{\rho}}^{\ell}$,
and $||.||_{2}^{2}$ regularizer for $\mathbbm{B}_{2M}^{\ell}$. See Table 1 for
the description of regularizers and their corresponding distance generating
functions999Note that in Algorithm 1
$V_{\boldsymbol{x}^{\prime}}(\boldsymbol{x})$ is identified with
$V_{\boldsymbol{x}_{t}}(\boldsymbol{x})$..

| $R$ | $V_{\boldsymbol{x}^{\prime}}(\boldsymbol{x})$
---|---|---
$\boldsymbol{x}$ \- space | $\sum_{s,a}x_{s,a}\log(x_{s,a})$ | $\sum_{s,a}x^{\prime}_{s,a}\log(x^{\prime}_{s,a}/x_{s,a})$
$\boldsymbol{\lambda}$ \- space | $\frac{1}{2}||\boldsymbol{\lambda}||_{2}^{2}$ | $\frac{1}{2}||\boldsymbol{\lambda}-\boldsymbol{\lambda}^{\prime}||_{2}^{2}$
Table 1: SMD parameter setting

The gradients of $f$ with respect $\boldsymbol{x}$, and $\boldsymbol{\lambda}$
are as follows:

$\displaystyle\boldsymbol{g}^{\boldsymbol{x}}$
$\displaystyle=(\mathbf{\Gamma}-\hat{\mathbf{I}})\boldsymbol{\lambda}-\boldsymbol{r}$
$\displaystyle\boldsymbol{g}^{\boldsymbol{\lambda}}$
$\displaystyle=(\hat{\mathbf{I}}-\mathbf{\Gamma})^{T}\boldsymbol{x}$

Note that since we maximize with respect to $\boldsymbol{x}$, we take the
negative of the gradient with respect to $\boldsymbol{x}$. We now explain how
the estimates $\tilde{\boldsymbol{g}}^{\boldsymbol{x}}_{t}$ and
$\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}_{t}$ of
$\boldsymbol{g}^{\boldsymbol{x}}$ and $\boldsymbol{g}^{\boldsymbol{\lambda}}$,
respectively, are computed (Lines 2 and 3). We note that
$(s,a)\sim\boldsymbol{x}_{t}$ denotes the state-action pair sampled from the
distribution defined by
$\boldsymbol{x}_{t}\in\Delta_{\boldsymbol{\rho}}^{\ell}$.

a) The estimator $\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$ of
$\boldsymbol{g}^{\boldsymbol{x}}$ in Line 4 is computed as follows:

$\displaystyle(s,a)$
$\displaystyle\sim\left[\frac{1}{\ell}\right],\,\,\,\,s^{\prime}\sim\mathbf{\Gamma}((s,a),s^{\prime})$
(6)
$\displaystyle\tilde{\boldsymbol{g}}^{\boldsymbol{x}}(\boldsymbol{x},\boldsymbol{\lambda})$
$\displaystyle=\ell(\lambda_{s^{\prime}}-\lambda_{s}-r_{s,a})\boldsymbol{e}_{s,a},$
(7)

where $\lambda_{s}$ denotes the $s$-th entry of $\boldsymbol{\lambda}$, and
$\boldsymbol{e}_{s,a}$ is the unit vector in $\mathbb{R}^{\ell}$. Finally,
$\tilde{\boldsymbol{g}}^{\boldsymbol{x}}_{t}=\tilde{\boldsymbol{g}}^{\boldsymbol{x}}(\boldsymbol{x}_{t-1},\boldsymbol{\lambda}_{t-1})$.

b) The estimator $\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}$ for
$\boldsymbol{g}^{\boldsymbol{\boldsymbol{\lambda}}}$ in Line 4 is computed as
follows: For an $\boldsymbol{x}\in\Delta^{\ell}_{\boldsymbol{\rho}}$,

$\displaystyle(s,a)$
$\displaystyle\sim\boldsymbol{x},\,\,\,\,s^{\prime}\sim\mathbf{\Gamma}((s,a),s^{\prime})$
(8)
$\displaystyle\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}(\boldsymbol{x},\boldsymbol{\lambda})$
$\displaystyle=\boldsymbol{e}_{s}-\boldsymbol{e}_{s^{\prime}}~{},$ (9)

where $\boldsymbol{e}_{s}$ is the unit vector in $\mathbb{R}^{n}$. Finally,
$\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}_{t}=\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}(\boldsymbol{x}_{t-1},\boldsymbol{\lambda}_{t-1})$.

In Lemmas 3 and 4, we show that $\tilde{\boldsymbol{g}}^{\boldsymbol{x}}_{t}$
and $\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}_{t}$ are bounded estimators
as desired (see Definition 7). The proofs of both the lemmas are in Section 9.

1:Input: Desired accuracy $\varepsilon$.

2:Output: An $\varepsilon$-approx policy $\pi^{\varepsilon}$.

3:Parameter: $\eta^{\boldsymbol{\lambda}}\leq\frac{\varepsilon}{16}$,
$\eta^{\boldsymbol{x}}\leq\frac{\varepsilon}{8\ell(24M^{2}+1)})$, total
iteration number
$T\geq\max(\frac{8\log\ell}{\eta^{\boldsymbol{x}}\varepsilon},\frac{32M^{2}n}{\eta^{\boldsymbol{\lambda}}\varepsilon})$.

14:for  _$t=1,\ldots,T$_ do

2

  * •

Get $\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}_{t}$ as a
$((4M+1)\ell,(24M^{2}+1)\ell,||.||_{\Delta^{\ell}})$-bounded estimator of
$\boldsymbol{g}^{\boldsymbol{\lambda}}(\boldsymbol{x}_{t},\boldsymbol{\lambda}_{t})$.

  * •

Get $\tilde{\boldsymbol{g}}^{\boldsymbol{x}}_{t}$ as a $(2,||.||_{2})$-bounded
estimator of
$\boldsymbol{g}^{\boldsymbol{x}}(\boldsymbol{x}_{t},\boldsymbol{\lambda}_{t})$.

  * •

Let
$\boldsymbol{x}_{t+1}=\arg\min_{\boldsymbol{x}\in\Delta_{\boldsymbol{\rho}}^{\ell}}\langle\eta^{\boldsymbol{x}}\tilde{\boldsymbol{g}}^{\boldsymbol{x}}_{t},\boldsymbol{x}\rangle+V_{\boldsymbol{x}_{t}}(\boldsymbol{x})$.

  * •

Let
$\boldsymbol{\lambda}_{t+1}=\arg\min_{\boldsymbol{\lambda}\in\mathbb{B}_{2M}^{n}}\langle\eta^{\boldsymbol{\lambda}}\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}_{t},\boldsymbol{\lambda}\rangle+V_{\boldsymbol{\lambda}_{t}}(\boldsymbol{\lambda})$.

3 end for

5:Let
$(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})=\frac{1}{T}\sum_{t\in[T]}(\boldsymbol{x}_{t},\boldsymbol{\lambda}_{t})$.

6:Return $\pi^{\varepsilon}$ with
$\pi^{\varepsilon}(a|s)=x_{s,a}/\sum_{a^{\prime}\in[m]}x_{s,a^{\prime}}$.

Algorithm 1 Fair State-Visitation

######  Lemma 3.

The estimator
$\tilde{\boldsymbol{g}}^{\boldsymbol{x}}(\boldsymbol{x},\boldsymbol{\lambda})$
as constructed above is a
$((4M+1)\ell,(24M^{2}+1)\ell,||.||_{\Delta^{\ell}})$-bounded estimator.

######  Lemma 4.

The estimator
$\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}(\boldsymbol{x},\boldsymbol{\lambda})$
as constructed above is a $(2,||.||_{2})$ bounded estimator.

We note that if $\boldsymbol{\rho}=\boldsymbol{0}$, then our algorithm is the
same as that in [JS20] for the unconstrained AMDP. In Section 5, we first show
in Theorem 1 that
$(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\epsilon})$ is such that
$\mathbb{E}\texttt{Gap}(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})\leq\varepsilon$.
Then, in Theorem 2 we show that the policy $\pi^{\varepsilon}$ constructed at
Line 6 is indeed $(3\varepsilon,\varepsilon)$-approximate in expectation.
Theorem 2, which shows that Algorithm 1 outputs a policy with simultaneous
approximation guarantees on fairness and optimality, is the main contribution
of our work.

##  5 Theoretical Results

We show in Theorem 1 that
$(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})$ computed
by Algorithm 1 at Line 5 satisfies
$\mathbb{E}\texttt{Gap}(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})\leq\varepsilon$.

######  Theorem 1.

Given a state-visitation fairness vector $\boldsymbol{\rho}\in[0,1/n)^{n}$,
desired accuracy $\varepsilon>0$, and bounded estimators
$\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$ and
$\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}$ as given in Lemmas 3 and
4, Algorithm 1 with step sizes
$\eta^{\boldsymbol{x}}\leq\frac{\varepsilon}{8\ell(24M^{2}+1)}$ and
$\eta^{\boldsymbol{\lambda}}\leq\varepsilon/16$, at the end of
$T\geq\max(\frac{8\log\ell}{\eta^{\boldsymbol{x}}\varepsilon},\frac{32M^{2}n}{\eta^{\boldsymbol{\lambda}}\varepsilon})$
rounds, at Line 5 computes
$(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})$ such that
$\mathbb{E}\texttt{Gap}(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})\leq\varepsilon$
.

The proof of Theorem 1 is in Section 10. We note that the expectation here is
over the random bits used by the algorithm and the stochastic transitions made
by the MDP to the next state $s^{\prime}$ given state action pair $(s,a)$. We
make additional remarks regarding Theorem 1 below.

  * •

Note that the step size $\eta^{\boldsymbol{x}}$ and time horizon $T$ depend on
the mixing time of an MDP. The problem of estimating (or providing an upper
bound) on the mixing time of an MDP is an active research topic [Pau15,
ZCKM20]. We rely on these techniques to estimate the mixing time used to
compute step size and time horizon.

  * •

The lower bound on $T$ is proportional to $M$, which becomes very large as
$\rho_{s}$ approaches $1/n$. Hence, if $\rho_{s}$’s are close to $1/n$ then
the search space of $\boldsymbol{\lambda}$ in Algorithm 1 increases and it
converges slowly. Note that, if $\rho_{s}\leq 1/2n$ then
$d_{\boldsymbol{\rho}}\leq\max_{s}\frac{1}{\rho_{s}}$ then our algorithm
yields similar convergence time (in order) as that of the constrained AMDP
algorithm by [JS20] which only computes a feasible solution without reward
maximization. This bound, without Assumption 2, can be obtained from a
stationary distribution of (any) strictly feasible policy.

Next, in Theorem 2 we show that the the policy $\pi^{\varepsilon}$ computed by
Algorithm 1 at Line 6 is $(3\varepsilon,\varepsilon)$-approximate.

######  Theorem 2.

Given a state-visitation fairness vector
$\boldsymbol{\rho}\in[0,\frac{1}{n})^{n}$, desired accuracy $\varepsilon>0$
and bounded estimators $\tilde{g}^{\boldsymbol{x}}_{t}$ and
$\tilde{g}^{\boldsymbol{\lambda}}_{t}$ as given in Lemmas 3 and 4, Algorithm 1
with step sizes and $T$ as in Theorem 1, at Line 8 computes a policy
$\pi^{\varepsilon}$ which is $(3\varepsilon,\varepsilon)$-approximate in
expectation.

The proof of Theorem 2 is in Section 11. As stated earlier, Algorithm 1 is
oblivious to the presence of fair-action. Assumption 2 is not necessary if
strict feasibility (even for $\boldsymbol{\rho}\in[0,1]^{n}$ such that
$\sum_{s}\rho_{s}<1$) is known a priori. The convergence time $T$ in this case
would depend on the reward and the stationary distribution induced by the
strictly feasible policy.

##  6 Simulations

In this section, we evaluate our proposed algorithm and its several properties
in experiments. Figure 1 presents the MDP transition dynamics specification
with three states, $s_{0},s_{1},s_{2}$, and two actions, $a_{0},a_{1}$. The
reward of taking $a_{0}$ on $s_{0}$ is 1, and 0.1 otherwise.

$s_{0}$$s_{1}$$s_{2}$ Figure 1: The MDP specification, with three states and
two actions. On each state, action $a_{0}$ takes the blue transition with
probability of 90%, and takes the yellow transition with probability of 10%;
action $a_{1}$ has the converse effect.

It is easy to see the optimal unconstrained policy is a deterministic one with
$\pi(s_{0})=a_{0}$, $\pi(s_{1})=a_{1}$, and $\pi(s_{2})=a_{0}$. This policy
yields an average-reward of 0.526, with the three states $s_{0},s_{1}$, and
$s_{2}$ being visited 47.4%, 43.5%, and 9.1% of times respectively. Also,
since the MDP is feasible for $\boldsymbol{\rho}\leq 1/3$, we do not use the
fair-action. Recall, $\boldsymbol{x}^{\varepsilon}$ is the average of
$\boldsymbol{x}_{t}$, $t\in[T]$. Further, we use $\boldsymbol{x}$ to denote
the state-action frequency vector of policy $\pi^{\varepsilon}$ derived from
$\boldsymbol{x}^{\varepsilon}$ (Step 8, Algorithm 1). Note that the reward of
$\pi^{\varepsilon}$ is $\boldsymbol{r}^{T}\boldsymbol{x}$.

###  6.1 Convergence of SMD

We first empirically demonstrate the convergence of SMD. We choose a fairness
constraint $\rho=[0.1,0.1,0.25]$, which is violated by the unconstrained
optimal policy. We implemented Algorithm 1 to solve this problem with $M=100$
and $\eta^{\boldsymbol{x}}=\eta^{\boldsymbol{\lambda}}=0.01$. Figure 2 shows
the average-reward and state-visitation frequencies for policies obtained
after different number of gradient descent steps, up to a maximum of 20,000.
The results are averaged over 100 runs, with the shaded regions representing
the standard deviations. We also see that the reward
$\boldsymbol{r}^{T}\boldsymbol{x}$, approaches the fair optimal value and the
visitation frequency of $s_{2}$ approaches the required value (25%). We note
that $\boldsymbol{r}^{T}\boldsymbol{x}$ can cross the optimal value as
$\boldsymbol{x}$ is only approximately fair. In Figure 2, we also plot
$\boldsymbol{r}^{T}\boldsymbol{x}^{\varepsilon}$ along with
$\boldsymbol{r}^{T}\boldsymbol{x}$ for different values of $T$

Figure 2: Reward (left) and state-visitation frequency (right) of a policy
learned after $T$ steps, averaged over 100 runs.

###  6.2 Gap Function

The theorems prove that the difference between
$\boldsymbol{r}^{T}\boldsymbol{x}$ and the optimal reward decreases with the
value of the the gap function
Gap$(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})$. Figure
3 plots the reward difference as the gap decreases.

Figure 3: Gap function value (blue) and reward difference (orange) learned
after $T$ steps, averaged over 100 runs.

###  6.3 Varying the constraint $\boldsymbol{\rho}$

In this section, we vary $\rho_{2}\in\\{0.1,0.15,0.2,0.25,0.3\\}$. For each
$\boldsymbol{\rho}$ parameter, we run SMD 10 times, with 10,000 steps each.
Figure 4 plots the reward and state-visitation frequency on $s_{2}$ for
different $\boldsymbol{\rho}$. The colors blue, orange, green, red, and purple
indicate increasing values of $\rho_{2}$. As expected, all fairness
constraints are satisfied, but as $\rho_{2}$ increases, the optimal rewards
decreases. This decrease captures the cost of fairness as the fairness
constraints become more stringent.

Figure 4: Reward (left) and visitation frequency of $s_{2}$ (right) of
policies learned after $T$ steps, averaged over 10 runs, for each
$\rho_{2}\in\\{0.1,0.15,0.2,0.25,0.3\\}$. The colors of blue, orange, green,
red, and purple indicate increasing values of $\rho_{2}$. Dashed lines
indicate optimal rewards and visitation frequency constraints.

##  7 Additional Related Work

Fairness considerations in machine learning algorithms such as classification,
regression, recommendation systems, and online learning have been extensively
studied. Some active threads of research in this area include mathematical
modeling of various fairness notions [DHP+12, FFM+15, ADW19], proposing fair
algorithms [ZVRG17, DOBD+18], evaluating algorithms in terms of their fairness
guarantees [CDG18, GJKR18], characterizing the loss in efficiency and the
trade-off between different fairness notions[KMR16, Cho17], studying short-
term and long-term impacts of fairness constraints [ZTL+20, HNG19], and so on.

Besides fairness, uncertainty is an integral part of real-life decision-
making. Ensuring fair decision-making in the presence of uncertainty
introduces many non-trivial challenges such as variance in decision-making
[Eth20], feasibility and implementation issues [ZS15], verification of the
fairness guarantees [GJVS20, TBKN20], and so on. In many situations, it is not
clear when one should attribute the unfairness in the system to the
algorithms’ bias instead of the uncertainty in the system itself [ZB18].
Further, the ubiquity of online-decision making applications such as online
advertising and crowd-sourcing has fueled the research interest in studying
the fairness considerations in sequential rather than one-shot decision-
making. Even though algorithmic fairness in static decision-making settings
such as classification and regression has been extensively studied, the study
of fairness in sequential decision-making tasks under uncertainty is fairly
recent [ZL20]. The multi-armed bandit (MAB) problem, which is a classical
model used to capture online decision-making under uncertainty, has recently
been studied with several well-motivated notions of fairness, including
meritocratic fairness [JKMR16], quota-based fairness [LLJ19, CKSV19, PGNN20,
RBPS20], and Nash Social Welfare [HMS20].

Another relevant direction of research in RL of _safe_ RL (see [GF15] for a
comprehensive survey). The goal is safe RL is to learn policies that maximize
the expected reward while respecting certain safety constraints during
learning/deployment. The approach to safe RL that is closest to our work is
the one that uses _constrained criterion_. This setting aims to maximize the
expected reward while ensuring that some pre-specified constraints are
satisfied. Generally, these constraints are a function of the reward [Gei06,
TDCM12].

##  8 Conclusion

Discussion: In this paper, we studied the Fair-AMDP problem and proposed an
SMD-based algorithm that outputs a policy with simultaneous approximation
guarantee on average-reward and state-visitation fairness. First, we note that
our results in this work can be extended to ensure state-action visitation
fairness: for a $\boldsymbol{\rho}\in[0,1]^{\ell}$ a fair policy $\pi$ with
stationary distribution $\nu^{\pi}$ satisfies
$\Pi\cdot\nu^{\pi}\geq\boldsymbol{\rho}$, provided there is a strictly
feasible solution.

In Section 3.2, we determined the domain for the primal and dual variables for
approximately solving the bilinear saddle point problem using the SMD
framework. In particular, we restricted our domain of $\boldsymbol{x}$ to
$\Delta_{\boldsymbol{\rho}}^{\ell}$, and we removed $\boldsymbol{\mu}$
variables from the optimization problem. This in turn demands that in
Algorithm 1 at step 4, $\boldsymbol{x}_{t}$ is computed by projecting to
$\Delta_{\boldsymbol{\rho}}^{\ell}$. This approach can also be adapted for
computing a feasible solution for an AMDP with arbitrary linear constraints
and no objective function as done by [JS20].

Future Work: As noted earlier, the bound on the time complexity of our
algorithm becomes large as $\rho_{s}\rightarrow 1/n$. An immediate future
direction is to give an algorithm or a better analysis whose time complexity
has better dependence on $\boldsymbol{\rho}$. Additionally, observe that the
introduction of fair-action ensures the existence of a strictly feasible
solution when $\rho_{s}<1/n$ for all $s$. It is interesting to see whether
there is a fair action that ensures the same for a broader class of fairness
vectors $\boldsymbol{\rho}$.

We note that a convex program formulation for discounted-reward MDPs with
state-visitation fairness constraints is not immediate. Hence, the approach
discussed in this paper does not follow immediately, and we need novel
techniques to solve the state-visitation fairness problem for discounted-
reward MDPs.

The introduction of fairness constraints may lead to a loss in the maximum
average-reward that can be achieved. The difference in the average-reward of
an unconstrained optimal policy and that of an optimal policy that satisfies
fairness captures the _price of fairness_ in the system. The study of price of
fairness in the Fair-AMDP setting is an interesting future direction.

Another extension of our work is to study policies that provide both an upper
and a lower bound on the state-visitation frequency. However, it is not clear
if one can obtain such a guarantee using currently known techniques. Further,
the study of algorithms for finding policies that guarantee fairness in terms
of the relative state-visitation frequencies remains an interesting open
direction.

##  9 Proofs of Lemmas in Sections 3 and 4

Before we proceed to the proofs of the lemmas, we state the following lemma
(Lemma 5 in [JS20]) here without a proof. This lemma is used in the proof of
Lemma 2 and Theorem 2.

######  Lemma 5.

Give a mixing AMDP with mixing time $t_{\text{mix}}$, a policy $\pi$, and its
transition probability matrix $\mathbf{\Gamma}^{\pi}\in\mathbb{R}^{n\times n}$
and stationary distribution $\boldsymbol{\nu}^{\pi}$, the following holds:

$||(\mathbf{I}-\mathbf{\Gamma}^{\pi}+\boldsymbol{1}\boldsymbol{\nu}^{T})^{-1}||_{\infty}\leq
2t_{\text{mix}}~{}.$

###  9.1 Proof of Lemma 1

See 1

###### Proof.

Since $\boldsymbol{\lambda}^{*}$, $\boldsymbol{\mu}^{*}$, and $\beta^{*}$
satisfies the Fair-LP primal, we have

$(\hat{\mathbf{I}}-\mathbf{\Gamma})\boldsymbol{\lambda}^{*}+\boldsymbol{r}+\mathbf{C}^{T}\boldsymbol{\mu}^{*}\leq\beta^{*}\cdot\boldsymbol{1}~{}.$

Consider the policy $\pi^{f}$ which only pulls the fair-action at every point.
Then $\nu^{f}=(\frac{1}{n}\ldots\frac{1}{n})$ is the stationary distribution
corresponding to $\pi^{f}$. Since the reward for pulling the control action in
any state is zero, $(\Pi^{f}\boldsymbol{\nu}^{f})^{T}\boldsymbol{r}=0$.
Multiplying the above equation by $(\Pi^{f}\boldsymbol{\nu}^{f})^{T}$, we have

$\sum_{s}\frac{1}{n\rho_{s}}\mu_{i}^{*}\leq\beta^{*}$

Further, using strong duality we have
$\boldsymbol{r}^{T}\boldsymbol{x}^{*}+\sum_{s}\mu_{s}^{*}=\beta^{*}$. Hence,

$\sum_{s}\frac{1}{n\rho_{s}}\mu_{s}^{*}\leq\boldsymbol{r}^{T}\boldsymbol{x}^{*}+\sum_{s}\mu_{s}^{*}~{}.$

Since $\boldsymbol{r}\in[0,1]^{\ell}$, and $\boldsymbol{x}\in\Delta^{\ell}$,

$\sum_{s}(\frac{1-n\rho_{s}}{n\rho_{s}})\mu_{s}^{*}\leq\boldsymbol{r}^{T}\boldsymbol{x}^{*}\leq
1~{}.$

This implies for all $s$, $\mu_{s}^{*}\leq\frac{n\rho_{s}}{1-n\rho_{s}}$
assuming $\rho_{s}<\frac{1}{n}$. ∎

Remark: Suppose we had a priory knowledge that there exists a strictly
feasible policy $\pi$, that is a policy $\pi$ with stationary distribution
$\boldsymbol{\nu}$ such that $\nu_{s}>\frac{1}{\rho_{s}}$. Then using the same
arguments as above we obtain $\mu_{s}^{*}(\frac{\nu_{s}}{\rho_{s}}-1)\leq
1-(\Pi\nu)^{T}\boldsymbol{r}$, where $(\Pi\boldsymbol{\nu})^{T}\boldsymbol{r}$
is the reward of the strictly feasible policy $\pi$ and is at least $0$.

###  9.2 Proof of Lemma 2

See 2

###### Proof.

Since $\boldsymbol{\lambda}^{*},\boldsymbol{\mu}^{*}$, and $\beta^{*}$ is the
solution to the Fair-LP primal,

$(\hat{\mathbf{I}}-\mathbf{\Gamma})\boldsymbol{\lambda}^{*}+\mathbf{C}^{T}\boldsymbol{\mu}^{*}+\boldsymbol{r}\leq\beta^{*}\boldsymbol{1}~{}.$
(10)

Let $\pi^{*}$ be the optimal policy corresponding to $\boldsymbol{x}^{*}$, and
let $\Pi^{*}\in\mathbb{R}^{\ell\times n}$ and $\boldsymbol{\nu}^{*}$ be its
corresponding matrix and stationary distribution respectively. Note that
$\mathbf{\Gamma}^{*}=(\Pi^{*})^{T}\cdot\mathbf{\Gamma}$ is the probability
transition matrix corresponding to the Markov chain induced by $\pi^{*}$.
Multiplying Equation 10 by $(\Pi^{*})^{T}$ from the left, and using the KKT
condition and that $\nu^{*}_{s}>0$ for all $s$ (as the MDP is ergodic;
Assumption 1) we have

$(\mathbf{I}-\mathbf{\Gamma}^{*})\boldsymbol{\lambda}^{*}+\mathbf{D}_{\boldsymbol{\rho}}\boldsymbol{\mu}^{*}+(\Pi^{*})^{T}\boldsymbol{r}=\beta^{*}\cdot\boldsymbol{1}~{}.$
(11)

In Equation 11, $\mathbf{D}_{\boldsymbol{\rho}}$ is the $n\times n$ diagonal
matrix with its $s$-th entry being $\frac{1}{\rho_{s}}$. It is easy to see
that $(\Pi^{*})^{T}\mathbf{C}^{T}=\mathbf{D}_{\boldsymbol{\rho}}$, and
$(\Pi^{*})^{T}\boldsymbol{1}=\boldsymbol{1}$. Denote
$(\Pi^{*})^{T}\boldsymbol{r}$ as $\boldsymbol{r}^{*}$. Since
$\langle\boldsymbol{\lambda}^{*},\boldsymbol{\nu}^{*}\rangle=0$, Equation 11
can be rewritten as follows

$(\mathbf{I}-\mathbf{\Gamma}^{*}+\boldsymbol{1}(\boldsymbol{\nu}^{*})^{T})\boldsymbol{\lambda}^{*}+\mathbf{D}_{\boldsymbol{\rho}}\boldsymbol{\mu}^{*}+\boldsymbol{r}^{*}=\beta^{*}\cdot\boldsymbol{1}~{}.$

This implies

$\displaystyle\boldsymbol{\lambda}^{*}$
$\displaystyle=~{}(\mathbf{I}-\mathbf{\Gamma}^{*}+\boldsymbol{1}(\boldsymbol{\nu}^{*})^{T})^{-1}(\beta^{*}\cdot\boldsymbol{1}-\mathbf{D}_{\boldsymbol{\rho}}\boldsymbol{\mu}^{*}-\boldsymbol{r}^{*})$
$\displaystyle||\boldsymbol{\lambda}^{*}||_{\infty}$
$\displaystyle\leq~{}||(\mathbf{I}-\mathbf{\Gamma}^{*}+\boldsymbol{1}(\boldsymbol{\nu}^{*})^{T})^{-1}||_{\infty}||\beta^{*}\cdot\boldsymbol{1}-\mathbf{D}_{\boldsymbol{\rho}}\boldsymbol{\mu}^{*}-\boldsymbol{r}^{*}||_{\infty}$
$\displaystyle\leq~{}2t_{\text{mix}}||\beta^{*}\cdot\boldsymbol{1}-\mathbf{D}_{\boldsymbol{\rho}}\boldsymbol{\mu}^{*}-\boldsymbol{r}^{*}||_{\infty}~{}\leq~{}2t_{\text{mix}}(1+\max_{s}\frac{n}{1-\rho_{s}n})$
(12)

The first inequality in Equation 12 follows from Lemma 5, and the second
follows by showing that
$||\beta^{*}\cdot\boldsymbol{1}-\mathbf{D}_{\boldsymbol{\rho}}\boldsymbol{\mu}^{*}-\Pi^{T}\boldsymbol{r}||_{\infty}\leq(1+\max_{s}\frac{n}{1-\rho_{s}n})$,
which we argue next. By strong duality we have
$\beta^{*}-\sum_{s}\mu_{s}^{*}=\boldsymbol{r}^{T}\boldsymbol{x}^{*}$. Since
$\boldsymbol{r}\in[0,1]^{\ell}$ and $\boldsymbol{x}^{*}\in\Delta^{\ell}$,
$0\leq\beta^{*}\leq 1+\sum_{s}\mu_{s}^{*}$. Hence from Lemma 1, we have for
any $s$

$-1-\sum_{s}\frac{n\rho_{s}}{1-n\rho_{s}}\leq\beta^{*}-\frac{\mu^{*}_{s}}{\rho_{s}}-r^{*}_{s}\leq
1+\max_{s}\frac{n}{1-n\rho_{s}}\ ~{}.$

Now as $\sum_{s}\rho_{s}<1$, $1+\sum_{s}\frac{n\rho_{s}}{1-n\rho_{s}}\leq
1+\max_{s}\frac{n}{1-n\rho_{s}})$. Hence,
$||\beta^{*}\cdot\boldsymbol{1}-\mathbf{D}_{\boldsymbol{\rho}}\boldsymbol{\mu}^{*}-\boldsymbol{r}^{*}||_{\infty}\leq(1+d_{\boldsymbol{\rho}})$.
∎

###  9.3 Proof of Lemma 3

See 3

###### Proof.

From the definition of
$\tilde{\boldsymbol{g}}^{\boldsymbol{x}}(\boldsymbol{x},\boldsymbol{\lambda})$
it follows that

$\displaystyle\mathbb{E}[\tilde{\boldsymbol{g}}^{\boldsymbol{x}}]$
$\displaystyle=\sum_{s,s^{\prime}\in[n],a\in[m]}\mathbf{\Gamma}((s,a),s)(\lambda_{s^{\prime}}-\lambda_{s}-r_{s,a})e_{s,a}$
$\displaystyle=~{}~{}(\mathbf{\Gamma}-\hat{\mathbf{I}})\boldsymbol{\lambda}-\boldsymbol{r}~{}.$

Also as $||\boldsymbol{\lambda}||_{\infty}\leq 2M$ and
$\boldsymbol{r}\in[0,1]^{\ell}$,
$|\lambda_{s^{\prime}}-\lambda_{s}-r_{s,a}|\leq 4M+1$. Hence,
$||\tilde{\boldsymbol{g}}^{\boldsymbol{x}}_{t}||_{{\infty}}\leq(4M+1)\ell$.
Next, we prove the bound on the second moment. For any
$\boldsymbol{x}^{\prime}\in\Delta^{\ell}$

$\displaystyle\mathbb{E}[||\tilde{\boldsymbol{g}}^{\boldsymbol{x}}||_{\boldsymbol{x}^{\prime}}^{2}]$
$\displaystyle\leq\sum_{s,a}\frac{1}{\ell}x^{\prime}_{s,a}(24M^{2}+1)\ell^{2}~{}.$

∎

###  9.4 Proof of Lemma 4

See 4

###### Proof.

First we show that $\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}$ is a
$(2,||.||_{2})$ bounded estimator. From the definition of
$\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}(\boldsymbol{x},\boldsymbol{\lambda})$
it follows that

$\displaystyle\mathbb{E}[\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}]$
$\displaystyle=\sum_{s,s^{\prime}\in[n],a\in[m]}x_{s,a}\mathbf{\Gamma}((s,a),s^{\prime})(e_{s}-e_{s^{\prime}})$
$\displaystyle=(\hat{\mathbf{I}}-\mathbf{\Gamma})^{T}\boldsymbol{x}~{}.$

For the bound on the second-moment, observe that
$||\tilde{\boldsymbol{g}}^{\boldsymbol{\lambda}}(\boldsymbol{x},\boldsymbol{\lambda})||_{2}^{2}\leq
2$ for any $\boldsymbol{x},\boldsymbol{\lambda}$. ∎

##  10  Proof of Theorem 1

See 1

First, in Lemma 6, we show that the Gap function can be upper bounded in terms
of the sequences of gradients
($\\{\boldsymbol{g}_{t}^{\boldsymbol{x}}\\}_{t=1}^{T}$ and
$\\{\boldsymbol{g}_{t}^{\boldsymbol{\lambda}}\\}_{t=1}^{T}$) and the sequence
of choices of algorithm ($\\{\boldsymbol{x}_{t}\\}_{t=1}^{T}$ and
$\\{\boldsymbol{\lambda}_{t}\\}_{t=1}^{T}$) as follows.

######  Lemma 6.

Let $\\{\boldsymbol{x}_{t},\boldsymbol{\lambda}_{t}\\}_{t=1}^{T}$ be a
sequence of the iterate values of Algorithm 1 and
$(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon}):=\frac{1}{T}\sum_{t=1}^{T}(\boldsymbol{x}_{t},\boldsymbol{\lambda}_{t})$.
Then we have

$\small\texttt{Gap}(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon})\leq\frac{1}{T}\sup_{\boldsymbol{x},\boldsymbol{\lambda}}\sum_{t=1}^{T}[\langle\boldsymbol{g}_{t}^{\boldsymbol{x}},\boldsymbol{x}_{t}-\boldsymbol{x}_{T+1}\rangle+\langle\boldsymbol{g}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{T+1}\rangle]$

Next, we bound the individual terms in the statement of Lemma 6 using
following supporting Lemmas.

######  Lemma 7.

[Lemma 12 [CJST19]] Let $\mathcal{X}$ be a non-empty compact convex and $R$ be
a $1$-convex regularizer with respect to $||.||$. Further, let
$\mathcal{X}^{\star}$ be the dual space of $\mathcal{X}$ and the dual norm of
a vector $\boldsymbol{\gamma}\in\mathcal{X}^{\star}$ be defined as
$||\boldsymbol{\gamma}||_{\star}:=\sup_{||\boldsymbol{\lambda}||\leq
1}\langle\boldsymbol{\gamma},\boldsymbol{\lambda}\rangle$. Given a sequence
$\\{\boldsymbol{\gamma}_{t}\\}_{t=1}^{T}\in\mathcal{X}^{\star}$, let the
sequence $\\{\boldsymbol{\lambda}_{t}\\}_{t=1}^{T}\in\mathcal{X}$ be defined
as

$\boldsymbol{\lambda}_{t}=\arg\min\limits_{\boldsymbol{\lambda}\in\mathcal{X}}\langle\boldsymbol{\gamma}_{t-1},\boldsymbol{\lambda}\rangle+V_{\boldsymbol{\lambda}_{t-1}}(\boldsymbol{\lambda}).$

Then, for any time instant $T\geq 2$, for every choice of
$\boldsymbol{\lambda}_{0}\in\mathcal{X}$ we have,

$\sum_{t=1}^{T}\langle\boldsymbol{\gamma}_{t},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{T+1}\rangle\leq
V_{\boldsymbol{\lambda}_{0}}(\boldsymbol{\lambda}_{T+1})+1/2\sum_{t=1}^{T}||\boldsymbol{\gamma}_{t}||_{\star}^{2}$

######  Lemma 8.

[Lemma 13 [CJST19]] Let $R$ be an entropic regularizer defined over
$\Delta_{\boldsymbol{\rho}}^{\ell}$. Also, let a sequence
$\\{\boldsymbol{\delta}_{t}\\}_{t=1}^{T}\in\mathbbm{R}_{+}^{\ell}$ be such
that $\delta_{ti}\leq 1.79$ for all $t\in[T]$ and $i\in[\ell]$. Let the
sequence $\\{\boldsymbol{x}_{t}\\}_{t=1}^{T}\in\mathcal{X}$ is generated as

$\boldsymbol{x}_{t}=\arg\min_{\boldsymbol{x}\in\Delta_{\boldsymbol{\rho}}^{\ell}}\langle\boldsymbol{\delta}_{t-1},\boldsymbol{x}\rangle+V_{\boldsymbol{x}_{t-1}}(\boldsymbol{x})$
with $\boldsymbol{\delta}_{0}=\boldsymbol{0}$ and
$\boldsymbol{x}_{0}\in\Delta_{\boldsymbol{\rho}}^{\ell}$. Then, the KL
divergence $V_{\boldsymbol{x}}(\boldsymbol{x}^{{}^{\prime}})$ satisfies

$\sum_{t=1}^{T}\langle\boldsymbol{\delta}_{t},\boldsymbol{x}_{t}-\boldsymbol{x}_{T+1}\rangle\leq
V_{\boldsymbol{x}_{0}}(\boldsymbol{x}_{T+1})+\sum_{t=1}^{T}||-\boldsymbol{\delta}_{t}||_{\boldsymbol{x}_{t}}^{2}.$
(13)

Note that the results of Lemma 7 (and Lemma 13) holds for any arbitrary
sequence $\\{\boldsymbol{\gamma}_{t}\\}_{t=1}^{T}$ (and,
$\\{\boldsymbol{\delta}_{t}\\}_{t=1}^{T}$). We are now ready to prove the
Theorem.

###### Proof of Theorem 1.

First note that by the choice of $\eta^{\boldsymbol{x}}$ we have
$||\eta^{\boldsymbol{x}}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\infty}\leq
1/2$. Hence, we invoke Lemma 13 with
$\boldsymbol{\delta}_{t}=\eta^{\boldsymbol{x}}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$
and have

$\sum_{t=1}^{T}\langle\eta^{\boldsymbol{x}}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}},\boldsymbol{x}_{t}-\boldsymbol{x}_{T+1}\rangle\leq
V_{\boldsymbol{x}_{0}}(\boldsymbol{x}_{T+1})+\eta^{\boldsymbol{x}2}\sum_{t=1}^{T}||\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\boldsymbol{x}_{t}}^{2}.$
(14)

Similarly, using Lemma 7 with
$\boldsymbol{\gamma}_{t}=\eta^{\boldsymbol{\lambda}}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}$
we have

$\sum_{t=1}^{T}\langle\eta^{\boldsymbol{\lambda}}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{T+1}\rangle\leq
V_{\boldsymbol{\lambda}_{0}}(\boldsymbol{\lambda}_{T+1})+\frac{\eta^{\boldsymbol{\lambda}2}}{2}\sum_{t=1}^{T}||\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}||_{2}^{2}.$
(15)

Let
$\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}}:=\boldsymbol{g}_{t}^{\boldsymbol{x}}-\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$
and
$\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}:=\boldsymbol{g}_{t}^{\boldsymbol{\lambda}}-\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}$
and define the sequence

$\displaystyle\hat{\boldsymbol{x}}_{1}=\boldsymbol{x}_{1},\hskip 20.0pt$
$\displaystyle\hat{\boldsymbol{x}}_{t+1}=\arg\min_{\boldsymbol{x}\in\Delta_{\boldsymbol{\rho}}^{\ell}}\langle\eta^{\boldsymbol{x}}\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}},\boldsymbol{x}\rangle+V_{\hat{\boldsymbol{x}}_{t}}(\boldsymbol{x})$
$\displaystyle\hat{\boldsymbol{\lambda}}_{1}=\boldsymbol{\lambda}_{1},\hskip
20.0pt$
$\displaystyle\hat{\boldsymbol{\lambda}}_{t+1}=\arg\min_{\boldsymbol{\lambda}\in\mathbbm{B}_{2M}^{\ell}}\langle\eta^{\boldsymbol{\lambda}}\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}\rangle+V_{\hat{\boldsymbol{\lambda}}_{t}}(\boldsymbol{\lambda})$

Notice that,

$\displaystyle||\eta^{\boldsymbol{x}}\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\infty}$
$\displaystyle\underset{(i)}{\leq}||\eta^{\boldsymbol{x}}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\infty}+||\eta^{\boldsymbol{x}}\boldsymbol{g}_{t}^{\boldsymbol{x}}||_{\infty}\underset{(ii)}{\leq}||\eta^{\boldsymbol{x}}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\infty}+||\eta^{\boldsymbol{x}}\mathbb{E}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\infty}\underset{(iii)}{\leq}||\eta^{\boldsymbol{x}}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\infty}+\mathbb{E}||\eta^{\boldsymbol{x}}\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\infty}\leq
1.$

$(i)$ follows from the definition of
$\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$, $(ii)$ follows from fact that
$\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$ is an unbiased estimator and
finally $(iii)$ is due to Jensen’s inequality. Invoke Lemma 13 by setting
$\boldsymbol{\delta}_{t}=-\eta^{\boldsymbol{x}}\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$
we obtain,

$\displaystyle\sum_{t=1}^{T}\langle\eta^{\boldsymbol{x}}\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}},\hat{\boldsymbol{x}}_{t}-\boldsymbol{x}_{T+1}\rangle\leq
V_{\hat{\boldsymbol{x}}_{0}}(\boldsymbol{x}_{T+1})+\eta^{\boldsymbol{x}2}\sum_{t=1}^{T}||\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\hat{\boldsymbol{x}}_{t}}^{2}$
(16) Similarly, for $\boldsymbol{\lambda}$-space, using Lemma 7 with
$\boldsymbol{\gamma}_{t}=\eta^{\boldsymbol{\lambda}}\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}$,
we have
$\displaystyle\sum_{t=1}^{T}\langle\eta^{\boldsymbol{\lambda}}\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}},\hat{\boldsymbol{\lambda}}_{t}-\boldsymbol{\lambda}_{T+1}\rangle\leq
V_{\hat{\boldsymbol{\lambda}}_{0}}(\boldsymbol{\lambda}_{T+1})+\frac{\eta^{\boldsymbol{\lambda}2}}{2}\sum_{t=1}^{T}||\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}||_{2}^{2}$
(17)

Since,
$\boldsymbol{g}_{t}^{\boldsymbol{x}}=\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}+\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$
and
$\boldsymbol{g}_{t}^{\boldsymbol{\lambda}}=\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}+\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}$
we can write,

$\displaystyle\sum_{t=1}^{T}[\langle\boldsymbol{g}_{t}^{\boldsymbol{x}},\boldsymbol{x}_{t}-\boldsymbol{x}_{T+1}\rangle+\langle\boldsymbol{g}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{T+1}\rangle]$
$\displaystyle=$
$\displaystyle\sum_{t=1}^{T}[\langle\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}},\boldsymbol{x}_{t}-\boldsymbol{x}_{T+1}\rangle+\langle\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{T+1}\rangle]+\sum_{t=1}^{T}[\langle\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}},\boldsymbol{x}_{t}-\boldsymbol{x}_{T+1}\rangle+\langle\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{T+1}\rangle]$
$\displaystyle\leq$
$\displaystyle\frac{1}{\eta^{\boldsymbol{x}}}V_{\boldsymbol{x}_{0}}(\boldsymbol{x}_{T+1})+\frac{1}{\eta^{\boldsymbol{x}}}V_{\hat{\boldsymbol{x}}_{0}}(\hat{\boldsymbol{x}}_{T+1})+\sum_{t=1}^{T}\eta^{\boldsymbol{x}}\Big{[}||\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\boldsymbol{x}_{t}}^{2}+||\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}}||_{\hat{\boldsymbol{x}}_{t}}^{2}\Big{]}+\frac{1}{\eta^{\boldsymbol{x}}}\sum_{t=1}^{T}\langle\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}},\boldsymbol{x}_{t}-\hat{\boldsymbol{x}}_{t}\rangle$
$\displaystyle+\frac{1}{\eta^{\boldsymbol{\lambda}}}V_{\boldsymbol{\lambda}_{0}}(\boldsymbol{\lambda}_{T+1})+\frac{1}{\eta^{\boldsymbol{\lambda}}}V_{\hat{\boldsymbol{\lambda}}_{0}}(\hat{\boldsymbol{\lambda}}_{T+1})+\sum_{t=1}^{T}\frac{\eta^{\boldsymbol{\lambda}}}{2}\Big{[}||\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}||_{2}^{2}+||\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}||_{2}^{2}\Big{]}+\frac{1}{\eta^{\boldsymbol{\lambda}}}\sum_{t=1}^{T}\langle\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}_{t}-\hat{\boldsymbol{\lambda}}_{t}\rangle$

Note that
$\mathbbm{E}[||\hat{\boldsymbol{g}}^{\boldsymbol{x}}||_{\hat{\boldsymbol{x}}_{t}}^{2}]\leq\mathbbm{E}[||\tilde{\boldsymbol{g}}^{\boldsymbol{x}}||_{\hat{\boldsymbol{x}}_{t}}^{2}]\leq
v^{\boldsymbol{x}}$ because
$\mathbbm{E}[(X-\mathbbm{E}[X])^{2}]\leq\mathbbm{E}[X]^{2}$ and property
$(iii)$ of the bounded-estimator $\tilde{\boldsymbol{g}}^{\boldsymbol{x}}$
respectively. First, taking the supremum over
$(\boldsymbol{x},\boldsymbol{\lambda})$ and then taking expectation we obtain,

$\displaystyle\frac{1}{T}\mathbb{E}\sup_{\boldsymbol{x},\boldsymbol{\lambda}}$
$\displaystyle\Bigg{[}\sum_{t=1}^{T}[\langle\boldsymbol{g}_{t}^{\boldsymbol{x}},\boldsymbol{x}_{t}-\boldsymbol{x}\rangle+\langle\boldsymbol{g}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}\rangle]\Bigg{]}$
$\displaystyle\underset{(i)}{\leq}\sup_{\boldsymbol{x}}\frac{2}{\eta^{\boldsymbol{x}}T}V_{\boldsymbol{x}_{0}}(\boldsymbol{x})+2\eta^{\boldsymbol{x}}v^{\boldsymbol{x}}+\sup_{\boldsymbol{\lambda}}\frac{2}{\eta^{\boldsymbol{\lambda}}T}V_{\boldsymbol{\lambda}_{0}}(\boldsymbol{\lambda})+\eta^{\boldsymbol{\lambda}}v^{\boldsymbol{\lambda}}$
$\displaystyle\underset{(ii)}{\leq}\frac{2\log(\ell)}{\eta^{\boldsymbol{x}}T}+2\eta^{\boldsymbol{x}}v^{\boldsymbol{x}}+\frac{4\ell
M^{2}}{\eta^{\boldsymbol{\lambda}}T}+\eta^{\boldsymbol{\lambda}}v^{\boldsymbol{\lambda}}$
$\displaystyle\underset{(iii)}{\leq}\varepsilon$ (18)

Here, $(i)$ follows from the fact that
$\mathbb{E}[\langle\hat{\boldsymbol{g}}_{t}^{\boldsymbol{x}},\boldsymbol{x}_{t}-\hat{\boldsymbol{x}}_{t}\rangle]=\mathbb{E}[\langle\hat{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}_{t}-\hat{\boldsymbol{\lambda}}_{t}\rangle]=0$.
This is true because $\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{x}}$ and
$\tilde{\boldsymbol{g}}_{t}^{\boldsymbol{\lambda}}$ are unbiased estimators.
Next, $(ii)$ follows from the fact that the KL divergence over $\Delta^{\ell}$
(and also over $\Delta_{\boldsymbol{\rho}}^{\ell}$) is upper bounded by
$\log(\ell)$ i.e. $V_{\boldsymbol{x}_{0}}(\boldsymbol{x})\leq\log(\ell)$ for
any $\boldsymbol{x}_{0},\boldsymbol{x}\in\Delta^{\ell}$ and
$V_{\boldsymbol{\lambda}_{0}}(\boldsymbol{\lambda})=\frac{1}{2}||\boldsymbol{\lambda}_{0}-\boldsymbol{\lambda}||_{2}^{2}\leq
2\ell M^{2}$. Finally, $(iii)$ follows directly from the choice of the
parameters. Finally use Lemma 6 and Equation 18 to get the desired upper bound
on the expected duality gap i.e.
$\mathbb{E}\big{[}\texttt{Gap}((\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\varepsilon}))\big{]}\leq\varepsilon$.
∎

We now present the proof of supporting Lemmas.

###### Proof of Lemma 6.

For notational brevity denote
$\mathcal{U}:=\Delta_{\boldsymbol{\rho}}^{\ell}\times\mathbbm{B}_{2M}^{n}$ and
$\boldsymbol{u}\in\mathcal{U}$ such that $\boldsymbol{u}_{x}=\boldsymbol{x}$
and $\boldsymbol{u}_{\lambda}=\boldsymbol{\lambda}$. Similarly we denote
$\boldsymbol{u}^{\varepsilon}=(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{\epsilon})$
and $\boldsymbol{u}_{t}=(\boldsymbol{x}_{t},\boldsymbol{\lambda}_{t})$. First,
define a function

$\texttt{gap}(\boldsymbol{u};\boldsymbol{v}):=f(\boldsymbol{v}_{x},\boldsymbol{u}_{\lambda})-f(\boldsymbol{u}_{x},\boldsymbol{v}_{\lambda}).$

Note that $\texttt{gap}(\boldsymbol{u};\boldsymbol{u})=0$ for any
$\boldsymbol{u}\in\mathcal{U}$ and
$\texttt{Gap}(\boldsymbol{u})=\sup_{\boldsymbol{v}\in\mathcal{U}}\texttt{gap}(\boldsymbol{u};\boldsymbol{v})$.
Further, note that, for every $\boldsymbol{u}\in\mathcal{U}$ the function
$\texttt{gap}(\boldsymbol{u};\boldsymbol{v})$ is concave in $\boldsymbol{v}$,
hence we have,

$\texttt{gap}(\boldsymbol{u};\boldsymbol{v})\leq\texttt{gap}(\boldsymbol{u};\boldsymbol{u})+\langle\nabla_{\boldsymbol{v}}\texttt{gap}(\boldsymbol{u};\boldsymbol{u}),\boldsymbol{v}-\boldsymbol{u}\rangle=\langle\boldsymbol{g}^{\boldsymbol{u}},\boldsymbol{u}-\boldsymbol{v}\rangle$
(19)

Note here that
$\boldsymbol{g}^{\boldsymbol{u}}=-\nabla_{\boldsymbol{v}}\texttt{gap}(\boldsymbol{u};\boldsymbol{u})$.
Next, note that $\texttt{gap}(\boldsymbol{u};\boldsymbol{v})$ is convex in
$\boldsymbol{u}$ for any $\boldsymbol{v}\in\mathcal{U}$ hence we have,

$\displaystyle\texttt{gap}(\boldsymbol{u}^{\varepsilon};\boldsymbol{u})$
$\displaystyle\leq\frac{1}{T}\sum_{t=1}^{T}\texttt{gap}(\boldsymbol{u}_{t};\boldsymbol{u})$
$\displaystyle\leq\frac{1}{T}\sum_{t=1}^{T}\langle\boldsymbol{g}_{t}^{\boldsymbol{u}},\boldsymbol{u}_{t}-\boldsymbol{u}\rangle$
$\displaystyle=\frac{1}{T}\sum_{t=1}^{T}\bigg{[}\langle\boldsymbol{g}_{t}^{\boldsymbol{x}},\boldsymbol{x}_{t}-\boldsymbol{x}\rangle+\langle\boldsymbol{g}_{t}^{\boldsymbol{\lambda}},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}\rangle\bigg{]}$
(20)

Furthermore, the space $\mathcal{U}$ is closed and the function
$\texttt{gap}(.;.)$ is continuous in both the arguments we have that the
supremum is attained. Since the Equation 20 holds for every
$\boldsymbol{u}\in\mathcal{U}$ it also holds for the supremum. This completes
the proof. ∎

###### Proof of Lemma 7.

We begin with a simple observation for a distance generating function.

######  Observation 1.

For any $\boldsymbol{x},\boldsymbol{y},\boldsymbol{z}\in\mathcal{X}$ we have,
$-\langle\nabla
V_{\boldsymbol{x}}(\boldsymbol{y}),\boldsymbol{y}-\boldsymbol{z}\rangle=V_{\boldsymbol{x}}(\boldsymbol{z})-V_{\boldsymbol{y}}(\boldsymbol{z})-V_{\boldsymbol{x}}(\boldsymbol{y})$

The proof of the observation follows directly from the definition of $V$. In
what follows, we prove lemmas used in establishing the proof of Theorem 1. Let
$\boldsymbol{\lambda}_{0}\in\mathcal{X}$ and
$\boldsymbol{\gamma}_{0}=\boldsymbol{0}$.

From the first order optimality condition and convexity of $\mathcal{X}$ we
have,

$\langle\boldsymbol{\gamma}_{t-1}+\nabla
V_{\boldsymbol{\lambda}_{t-1}}(\boldsymbol{\lambda}_{t}),\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{T+1}\rangle\leq
0$ (21)

Use Equation 21 and Observation 1 with
$\boldsymbol{\lambda}:=\boldsymbol{\lambda}_{T+1}$ to get

$\displaystyle\sum_{t=1}^{T}\langle\boldsymbol{\gamma}_{t-1},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}\rangle$
$\displaystyle\leq-\sum_{t=1}^{T}\langle\nabla
V_{\boldsymbol{\lambda}_{t-1}}(\boldsymbol{\lambda}_{t}),\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda})\rangle=\sum_{t=1}^{T}V_{\boldsymbol{\lambda}_{t-1}}(\boldsymbol{\lambda})-V_{\boldsymbol{\lambda}_{t}}(\boldsymbol{\lambda})-V_{\boldsymbol{\lambda}_{t-1}}(\boldsymbol{\lambda}_{t})$
$\displaystyle=V_{\boldsymbol{\lambda}_{0}}(\boldsymbol{\lambda})+\sum_{t=1}^{T-1}V_{\boldsymbol{\lambda}_{t}}(\boldsymbol{\lambda})-\sum_{t=1}^{T}V_{\boldsymbol{\lambda}_{t}}(\boldsymbol{\lambda})-\sum_{t=0}^{T-1}V_{\boldsymbol{\lambda}_{t}}(\boldsymbol{\lambda}_{t+1})=V_{\boldsymbol{\lambda}_{0}}(\boldsymbol{\lambda})-\sum_{t=0}^{T}V_{\boldsymbol{\lambda}_{t}}(\boldsymbol{\lambda}_{t+1})$

Now simplify the LHS of the above equation as follows.

$\displaystyle\sum_{t=1}^{T}\langle\boldsymbol{\gamma}_{t-1},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}\rangle=\sum_{t=1}^{T-1}\langle\boldsymbol{\gamma}_{t},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}\rangle-\sum_{t=1}^{T-1}\langle\boldsymbol{\gamma}_{t},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{t+1}\rangle=\sum_{t=1}^{T}\langle\boldsymbol{\gamma}_{t},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}\rangle-\sum_{t=1}^{T}\langle\boldsymbol{\gamma}_{t},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{t+1}\rangle$

The first equality follows from the fact that
$\boldsymbol{\gamma}_{0}=\boldsymbol{0}$ and the second follows from
$\boldsymbol{\lambda}=\boldsymbol{\lambda}_{T+1}$. Thus, we have

$\displaystyle\sum_{t=1}^{T}\langle\boldsymbol{\gamma}_{t},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}\rangle\leq$
$\displaystyle
V_{\boldsymbol{\lambda}_{0}}(\boldsymbol{\lambda})+\sum_{t=0}^{T}\\{\langle\boldsymbol{\gamma}_{t},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{t+1}\rangle-
V_{\boldsymbol{\lambda}_{t}}(\boldsymbol{\lambda}_{t+1})\\}.$ (22)

Furthermore, the following inequalities hold for every iteration $t$,

$\langle\gamma_{t},\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{t+1}\rangle\underset{(i)}{\leq}||\boldsymbol{\gamma}_{t}||_{\star}||\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{t+1}||\underset{(ii)}{\leq}\frac{1}{2}||\boldsymbol{\gamma}_{t}||_{\star}^{2}+\frac{1}{2}||\boldsymbol{\lambda}_{t}-\boldsymbol{\lambda}_{t+1}||^{2}\underset{(iii)}{\leq}\frac{1}{2}||\boldsymbol{\gamma}_{t}||_{\star}^{2}+V_{\boldsymbol{\lambda}_{t}}(\boldsymbol{\lambda}_{t+1}).$

In the above expression, $(i)$ follows from Hölder’s inequality, $(ii)$
follows from the AM-GM inequality and $(iii)$ follows from the strong
convexity of the underlying distance-generating function. Using above
inequality with Equation 22 completes the proof of the lemma. ∎

###### Proof of Lemma 13.

Following the same steps to simplify the sum as in Lemma 7, Equation 22 we
obtain

$\displaystyle\sum_{t=1}^{T}\langle\boldsymbol{\delta}_{t},\boldsymbol{x}_{t}-\boldsymbol{x}_{T+1}\rangle$
$\displaystyle\leq
V_{\boldsymbol{x}_{0}}(\boldsymbol{x}_{T+1})+\sum_{t=0}^{T}\\{\langle\boldsymbol{\delta}_{t},\boldsymbol{x}_{t}-\boldsymbol{x}_{t+1}\rangle-
V_{\boldsymbol{x}_{t}}(\boldsymbol{x}_{t+1})\\}$
$\displaystyle=V_{\boldsymbol{x}_{0}}(\boldsymbol{x}_{T+1})+\sum_{t=0}^{T}\\{\langle-\boldsymbol{\delta}_{t},\boldsymbol{x}_{t+1}-\boldsymbol{x}_{t}\rangle-
V_{\boldsymbol{x}_{t}}(\boldsymbol{x}_{t+1})\\}$ (23)

We now provide the upper bound on the summation term in the RHS of Equation
23. Recall that the Fenchel conjugate $R^{\star}$ at point $\boldsymbol{y}$ is
given as

$R^{\star}(\boldsymbol{y})=\sup\\{\langle\boldsymbol{y},\boldsymbol{x}\rangle-R(\boldsymbol{x})|\boldsymbol{x}\in\Delta_{\boldsymbol{\rho}}^{\ell}\\}$

This implies,

$R^{\star}(\boldsymbol{y})\geq\langle\boldsymbol{y},\boldsymbol{x}\rangle-R(\boldsymbol{x})\hskip
10.0pt\text{for all }\boldsymbol{x}\in\Delta_{\boldsymbol{\rho}}^{\ell}\hskip
10.0pt\text{ and }\hskip 10.0ptR^{\star}(\nabla
R(\boldsymbol{x}))=\langle\nabla
R(\boldsymbol{x}),\boldsymbol{x}\rangle-R(\boldsymbol{x}).$ (24)

Let $\boldsymbol{y}:=\nabla R(\boldsymbol{x})$. From maximizing argument we
have, $\nabla R^{\star}(\boldsymbol{y})=\boldsymbol{x}$ hence we have

$\boldsymbol{x}=\nabla R^{\star}(\nabla R(\boldsymbol{x})).$ (25)

Using definition of $V_{\boldsymbol{x}}(\boldsymbol{x}^{\prime})$, Equation 24
and Equation 25 we obtain,

$\displaystyle\langle\boldsymbol{\delta},\boldsymbol{x}^{\prime}-\boldsymbol{x}\rangle-
V_{\boldsymbol{x}}(\boldsymbol{x}^{\prime})$ $\displaystyle=\langle\nabla
R(\boldsymbol{x})+\boldsymbol{\delta},\boldsymbol{x}^{\prime}\rangle-R(\boldsymbol{x}^{\prime})-[\langle\nabla
R(\boldsymbol{x}),\boldsymbol{x}\rangle-R(\boldsymbol{x})]-\langle\boldsymbol{x},\boldsymbol{\delta}\rangle$
$\displaystyle\leq R^{\star}(\nabla
R(\boldsymbol{x})+\boldsymbol{\delta})-R^{\star}(\nabla
R(\boldsymbol{x}))-\langle\nabla R^{\star}(\nabla
R(\boldsymbol{x})),\boldsymbol{\delta}\rangle$ (26)

Note that $R(.)$ is an entropic regularizer whose Fenchel dual is given by
$R^{\star}(\boldsymbol{x})=\log(\sum_{i=1}^{\ell}e^{\boldsymbol{x}_{i}})$
(See, Table 2.1, [SS11] , for reference). Hence,

$\displaystyle R^{\star}$ $\displaystyle(\nabla
R(\boldsymbol{x})+\boldsymbol{\delta})-R^{\star}(\nabla
R(\boldsymbol{x}))=\log\Big{(}\frac{\sum_{i\in[\ell]}e^{(\nabla
R(\boldsymbol{x})+\boldsymbol{\delta})_{i}}}{\sum_{i\in[\ell]}e^{(\nabla
R(\boldsymbol{x}))_{i}}}\Big{)}$
$\displaystyle\leq\log\Big{(}1+\frac{\sum_{i\in[\ell]}e^{(\nabla
R(\boldsymbol{x}))_{i}}(\delta_{i}+\delta_{i}^{2})}{\sum_{i\in[\ell]}e^{(\nabla
R(\boldsymbol{x}))_{i}}}\Big{)}$ ( $e^{a}\leq 1+a+a^{2}$ for $a\leq 1.79$)
$\displaystyle=\log\big{(}1+\langle\nabla R^{\star}(\nabla
R(\boldsymbol{x})),\boldsymbol{\delta}+\boldsymbol{\delta}^{2}\rangle\big{)}$
(Here, $\boldsymbol{\delta}^{2}$ is a vector with $i^{th}$ coordinate
$\delta_{i}^{2}$. ) $\displaystyle\leq\langle\nabla R^{\star}(\nabla
R(\boldsymbol{x})),\boldsymbol{\delta}\rangle+\langle\nabla R^{\star}(\nabla
R(\boldsymbol{x})),\boldsymbol{\delta}^{2}\rangle$ ( $\log(1+a)\leq a$)

Thus we have, $R^{\star}(\nabla
R(\boldsymbol{x})+\boldsymbol{\delta})-R^{\star}(\nabla
R(\boldsymbol{x}))-\langle\nabla R^{\star}(\nabla
R(\boldsymbol{x})),\boldsymbol{\delta}\rangle\leq\langle\nabla
R^{\star}(\nabla
R(\boldsymbol{x})),\boldsymbol{\delta}^{2}\rangle=||\boldsymbol{\delta}||_{\nabla
R^{\star}(\nabla
R(\boldsymbol{x}))}^{2}=||\boldsymbol{\delta}||_{\boldsymbol{x}}^{2}$. The
last equality follows from Equation 25. Using this in Equation 26 we have

$\langle\boldsymbol{\delta},\boldsymbol{x}^{\prime}-\boldsymbol{x}\rangle-
V_{\boldsymbol{x}}(\boldsymbol{x}^{\prime})\leq||\boldsymbol{\delta}||_{\boldsymbol{x}}^{2}.$
(27)

Finally, we complete the proof of the lemma by using the above result
(Equation 27) in Equation 23 i.e.

$\sum_{t=1}^{T}\langle\boldsymbol{\delta}_{t},\boldsymbol{x}_{t}-\boldsymbol{x}_{T+1}\rangle\leq
V_{\boldsymbol{x}_{0}}(\boldsymbol{x}_{T+1})+\sum_{t=1}^{T}||-\boldsymbol{\delta}_{t}||_{\boldsymbol{x}_{t}}^{2}$

∎

##  11 Proof of Theorem 2

See 2 In the proof of this theorem, for convenience we denote
$\pi^{\varepsilon}$ as $\pi$, and $\boldsymbol{\nu}^{\pi^{\varepsilon}}$ equal
to the stationary distribution corresponding to policy $\pi^{\varepsilon}$ as
$\boldsymbol{\nu}$. Further let $\boldsymbol{u}\in\mathbb{R}^{n}$ be such that
$u_{s}=\sum_{a\in[m]}x^{\varepsilon}_{s,a}$, and note that
$x^{\varepsilon}_{s,a}=\pi_{s,a}u_{s}$ implying that
$\boldsymbol{x}^{\varepsilon}=\Pi\cdot\boldsymbol{u}$, where
$\Pi\in\mathbb{R}^{\ell\times n}$ is the policy matrix corresponding to $\pi$.
Also, let $\boldsymbol{x}\in\mathbb{R}^{\ell}$ be such that
$\boldsymbol{x}=\Pi\mathbf{\nu}$. Finally, note that
$\mathbf{\Gamma}^{\pi}=\Pi^{T}\mathbf{\Gamma}$ is the probability transition
matrix of the Markov chain induced by the policy $\pi$. Hence,

$(\mathbf{I}-\mathbf{\Gamma}^{\pi})^{T}\boldsymbol{\nu}^{\pi}=\boldsymbol{0}~{}.$
(28)

Finally, let $\boldsymbol{x}^{*},\boldsymbol{\lambda}^{*}$ be the optimal
solution to the problem in Equation 3. We first prove a few lemmas which would
be used to prove the fairness and optimality.

######  Lemma 9.

$\mathbb{E}\left[(\boldsymbol{x}^{\varepsilon})^{T}\left((\mathbf{\Gamma}-\hat{\mathbf{I}})\boldsymbol{\lambda}^{*}-\boldsymbol{r}\right)+v^{*}\right]\leq\varepsilon$

###### Proof.

From Theorem 1, we have

$\mathbb{E}[f(\boldsymbol{x}^{*},\boldsymbol{\lambda}^{\varepsilon})-f(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{*})]\leq\varepsilon~{}.$

Substituting appropriate value of $f$ at the respective points and observing
that $(\boldsymbol{x}^{*})^{T}\boldsymbol{r}=v^{*}$ the lemma follows. ∎

######  Lemma 10.

$\mathbb{E}[\max_{\lambda\in\mathbb{B}_{2M}^{n}}\boldsymbol{u}^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})(\boldsymbol{\lambda}^{*}-\boldsymbol{\lambda})]\leq\varepsilon$

###### Proof.

From Theorem 1

$f(\boldsymbol{x}^{*},\boldsymbol{\lambda}^{\varepsilon})-\min_{\boldsymbol{\lambda}\in\mathbb{B}_{2M}^{n}}f(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda})\leq\varepsilon$
(29)

Further note that since
$(\boldsymbol{x}^{*})^{T}(\hat{\mathbf{I}}-\mathbf{\Gamma})=\boldsymbol{0}$,
we have

$f(\boldsymbol{x}^{\varepsilon},\boldsymbol{\lambda}^{*})\leq
f(\boldsymbol{x}^{*},\boldsymbol{\lambda}^{*})=f(\boldsymbol{x}^{*},\boldsymbol{\lambda}^{\varepsilon})~{}.$
(30)

Using Equations 29 and 30,

$\mathbb{E}[\max_{\lambda\in\mathbb{B}_{2M}^{n}}\boldsymbol{u}^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})(\boldsymbol{\lambda}^{*}-\boldsymbol{\lambda})]\leq\varepsilon~{}.$

∎

######  Lemma 11.

$\mathbb{E}||(\boldsymbol{u}-\boldsymbol{\nu}^{\pi})^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi}+\boldsymbol{1}(\boldsymbol{\nu}^{\pi})^{T})||_{1}\leq\frac{\varepsilon}{M}$
.

###### Proof.

First we show that
$\mathbb{E}||\boldsymbol{u}^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})||_{1}\leq\frac{\varepsilon}{M}$:

$\displaystyle
2M\mathbb{E}||\boldsymbol{u}^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})||_{1}$
$\displaystyle=\mathbb{E}\left[\max_{\boldsymbol{\lambda}\in\mathbb{B}_{2M}^{n}}\boldsymbol{u}^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})(-\boldsymbol{\lambda})\right]$
(31)
$\displaystyle=\mathbb{E}\left[\max_{\boldsymbol{\lambda}\in\mathbb{B}_{2M}^{n}}\boldsymbol{u}^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})(\boldsymbol{\lambda}^{*}-\boldsymbol{\lambda})-\boldsymbol{u}^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})\boldsymbol{\lambda}^{*}\right]$
(32)
$\displaystyle\leq\varepsilon+||\boldsymbol{\lambda}^{*}||_{\infty}\mathbb{E}||\boldsymbol{u}^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})||_{1}\leq\varepsilon+M\mathbb{E}||\boldsymbol{u}^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})||_{1}~{}.$
(33)

The last but one inequality follows from Lemma 10 and
$||\boldsymbol{\lambda}^{*}||_{\infty}\leq M$. Here, we have made use of the
fact that the box is over the range $2M$ whereas
$||\boldsymbol{\lambda}^{*}||_{\infty}\leq M$. In particular, this is the
place where enlarging the box helps (as in [JS20]). Finally, to prove the
lemma observe that
$(\boldsymbol{\nu}^{\pi})^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})=\boldsymbol{0}$
and
$||(\boldsymbol{u}-\boldsymbol{\nu}^{\pi})^{T}(\boldsymbol{1}(\boldsymbol{\nu}^{\pi})^{T})||_{1}=0$,
and hence

$\mathbb{E}||(\boldsymbol{u}-\boldsymbol{\nu}^{\pi})^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi}+\boldsymbol{1}(\boldsymbol{\nu}^{\pi})^{T})||_{1}=\mathbb{E}||(\boldsymbol{u}-\boldsymbol{\nu}^{\pi})^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})||_{1}\leq\frac{\varepsilon}{M}~{}.$

∎

######  Lemma 12.

Let $\boldsymbol{r}^{\pi}=\Pi^{T}\boldsymbol{r}$. Then
$\mathbb{E}\left[(\boldsymbol{u}-\boldsymbol{\nu}^{\pi})^{T}\boldsymbol{r}^{\pi}\right]\leq\varepsilon$
.

###### Proof.

The proof of the lemma is completed using the following sequence of equations.

$\displaystyle\mathbb{E}\left[(\boldsymbol{u}-\boldsymbol{\nu}^{\pi})^{T}\boldsymbol{r}^{\pi}\right]$
$\displaystyle=\mathbb{E}\left[(\boldsymbol{u}-\boldsymbol{\nu}^{\pi})^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi}+\boldsymbol{1}(\boldsymbol{\nu}^{\pi})^{T})(\mathbf{I}-\mathbf{\Gamma}^{\pi}+\boldsymbol{1}(\boldsymbol{\nu}^{\pi})^{T})^{-1}\boldsymbol{r}^{\pi}\right]$
(34)
$\displaystyle\leq\mathbb{E}||(\boldsymbol{u}-\boldsymbol{\nu}^{\pi})^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi}+\boldsymbol{1}(\boldsymbol{\nu}^{\pi})^{T})||_{1}||(\mathbf{I}-\mathbf{\Gamma}^{\pi}+\boldsymbol{1}(\boldsymbol{\nu}^{\pi})^{T})^{-1}\boldsymbol{r}^{\pi}||_{\infty}$
(35)
$\displaystyle\leq\mathbb{E}||(\boldsymbol{u}-\boldsymbol{\nu}^{\pi})^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi}+\boldsymbol{1}(\boldsymbol{\nu}^{\pi})^{T})||_{1}||(\mathbf{I}-\mathbf{\Gamma}^{\pi}+\boldsymbol{1}(\boldsymbol{\nu}^{\pi})^{T})^{-1}||_{\infty}||\boldsymbol{r}^{\pi}||_{\infty}$
(36) $\displaystyle\leq\frac{\varepsilon}{M}\cdot
2t_{\text{mix}}\leq\varepsilon~{}.$ (37)

The first inequality in the last line follows from Lemmas 5 and 11. The last
inequality follows by observing that $M\geq 2t_{\text{mix}}$. ∎

Proof of Fairness: Recall that $\mathbf{D}_{\boldsymbol{\rho}}$ is the
$n\times n$ diagonal matrix with its $s$-th entry being $\frac{1}{\rho_{s}}$.
It is easy to see that $\mathsf{C}\Pi=\mathbf{D}_{\boldsymbol{\rho}}$, Now we
show that
$\mathbf{D}_{\boldsymbol{\rho}}\boldsymbol{\nu}^{\pi}=\mathbf{C}\boldsymbol{x}\geq\boldsymbol{1}-\varepsilon$.
First we note that as $\boldsymbol{x}_{t}\in\Delta_{\boldsymbol{\rho}}^{\ell}$
for all $t\in[T]$,
$\boldsymbol{x}^{\varepsilon}\in\Delta^{\ell}_{\boldsymbol{\rho}}$. Hence, we
have
$\mathbf{D}_{\boldsymbol{\rho}}\boldsymbol{x}^{\varepsilon}\geq\boldsymbol{1}$.
The policy $\pi$ is $\varepsilon$-fair follows from sequence of equations
below:

$\displaystyle\mathbf{C}\boldsymbol{x}$
$\displaystyle=\mathbf{C}\boldsymbol{x}^{\varepsilon}+\mathbf{C}(\boldsymbol{x}-\boldsymbol{x}^{\varepsilon})=\mathbf{C}\boldsymbol{x}^{\varepsilon}+\mathbf{C}\Pi(\boldsymbol{\nu}^{\pi}-\boldsymbol{u})$
$\displaystyle=\mathbf{C}\boldsymbol{x}^{\varepsilon}+\mathbf{D}_{\boldsymbol{\rho}}(\mathbf{I}-(\mathbf{\Gamma}^{\pi})^{T}+\boldsymbol{\nu}^{\pi}\boldsymbol{1}^{T})^{-1}(\mathbf{I}-(\mathbf{\Gamma}^{\pi})^{T}+\boldsymbol{\nu}^{\pi}\boldsymbol{1}^{T})(\boldsymbol{\nu}^{\pi}-\boldsymbol{u})$
$\displaystyle\geq\boldsymbol{1}-||\mathbf{D}_{\boldsymbol{\rho}}(\mathbf{I}-(\mathbf{\Gamma}^{\pi})^{T}+\boldsymbol{\nu}^{\pi}\boldsymbol{1}^{T})^{-1}(\mathbf{I}-(\mathbf{\Gamma}^{\pi})^{T}+\boldsymbol{\nu}^{\pi}\boldsymbol{1}^{T})(\boldsymbol{\nu}^{\pi}-\boldsymbol{u})||_{\infty}\cdot\boldsymbol{1}$
$\displaystyle\geq\boldsymbol{1}-\left(||\mathbf{D}_{\boldsymbol{\rho}}(\mathbf{I}-(\mathbf{\Gamma}^{\pi})^{T}+\boldsymbol{\nu}^{\pi}\boldsymbol{1}^{T})^{-1}||_{\infty}||(\mathbf{I}-(\mathbf{\Gamma}^{\pi})^{T}+\boldsymbol{\nu}^{\pi}\boldsymbol{1}^{T})(\boldsymbol{\nu}^{\pi}-\boldsymbol{u})||_{1}\right)\cdot\boldsymbol{1}$
$\displaystyle\geq\boldsymbol{1}-((\max_{s}\frac{1}{\rho_{s}})\cdot
2t_{\text{mix}}\cdot\frac{\varepsilon}{M})\cdot\boldsymbol{1}\geq\boldsymbol{1}-\varepsilon\boldsymbol{1}~{}.$

The first inequality in the last line follows from Lemmas 5 and 11, and the
last inequality follows by observing that
$M\geq(\max_{s}\frac{1}{\rho_{s}})\cdot 2t_{\text{mix}}$.

Proof of Optimality: We have the following sequence of equations.

$\displaystyle
v^{\pi}=\mathbb{E}\left[(\boldsymbol{\nu}^{\pi})^{T}\boldsymbol{r}^{\pi}\right]$
$\displaystyle=\mathbb{E}\left[(\boldsymbol{\nu}^{\pi})^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})\boldsymbol{\lambda^{*}}+(\boldsymbol{\nu}^{\pi})^{T}\boldsymbol{r}^{\pi}\right]$
$\displaystyle=\mathbb{E}\left[(\boldsymbol{\nu}^{\pi}-\boldsymbol{u})^{T}((\mathbf{I}-\mathbf{\Gamma}^{\pi})\lambda^{*}+\boldsymbol{r}^{\pi})+\boldsymbol{u}^{T}((\mathbf{I}-\mathbf{\Gamma}^{\pi})\boldsymbol{\lambda^{*}}+\boldsymbol{r}^{\pi})\right]$
$\displaystyle\geq\mathbb{E}\left[(\boldsymbol{\nu}^{\pi}-\boldsymbol{u})^{T}(\mathbf{I}-\mathbf{\Gamma}^{\pi})\boldsymbol{\lambda^{*}}\right]+\mathbb{E}\left[(\boldsymbol{\nu}^{\pi}-\boldsymbol{u})^{T}\boldsymbol{r}^{\pi}\right]+v^{*}-\varepsilon$
$\displaystyle\geq v^{*}-3\varepsilon~{}.$

The last but one inequality follows from Lemma 9, and the last inequality
follows from Lemmas 11 and 12.

## Acknowledgements

Ganesh Ghalme is thankful for the financial support from the Israeli Ministry
of Science and Technology grant 19400214. Vineet Nair is thankful to be
supported by the European Union’s Horizon 2020 research and innovation program
under grant agreement No 682203 -ERC-[Inf-Speed-Tradeoff]. Vishakha Patil is
grateful for the support of a Google PhD Fellowship.

## References

  * [ADW19] Alekh Agarwal, Miroslav Dudík, and Zhiwei Steven Wu.  Fair regression: Quantitative definitions and reduction-based algorithms.  In Proceedings of the 36th International Conference on Machine Learning, ICML, volume 97, pages 120–129. PMLR, 2019.
  * [Alt99] Eitan Altman.  Constrained Markov decision processes, volume 7.  CRC Press, 1999.
  * [ASB+19] Muhammad Ali, Piotr Sapiezynski, Miranda Bogen, Aleksandra Korolova, Alan Mislove, and Aaron Rieke.  Discrimination through optimization: How facebook’s ad delivery can lead to biased outcomes.  Proc. ACM Hum. Comput. Interact., 3(CSCW):199:1–199:30, 2019\.
  * [BCGH96] James A Berkovec, Glenn B Canner, Stuart A Gabriel, and Timothy H Hannan.  Mortgage discrimination and FHA loan performance.  Cityscape, pages 9–24, 1996.
  * [BHN18] Solon Barocas, M. Hardt, and A. Narayanan.  Fairness and Machine Learning Limitations and Opportunities.  https://fairmlbook.org, 2018.
  * [BT02] Ronen I Brafman and Moshe Tennenholtz.  R-max-a general polynomial time algorithm for near-optimal reinforcement learning.  Journal of Machine Learning Research, 3(Oct):213–231, 2002.
  * [CDG18] Sam Corbett-Davies and Sharad Goel.  The measure and mismeasure of fairness: A critical review of fair machine learning.  arXiv preprint arXiv:1808.00023, 2018.
  * [Cho17] Alexandra Chouldechova.  Fair prediction with disparate impact: A study of bias in recidivism prediction instruments.  Big data, 5(2):153–163, 2017.
  * [CJST19] Yair Carmon, Yujia Jin, Aaron Sidford, and Kevin Tian.  Variance reduction for matrix games.  In Advances in Neural Information Processing Systems, pages 11381–11392, 2019.
  * [CKSV19] L Elisa Celis, Sayash Kapoor, Farnood Salehi, and Nisheeth Vishnoi.  Controlling polarization in personalization: An algorithmic framework.  In Proceedings of the conference on fairness, accountability, and transparency, pages 160–169, 2019.
  * [CMPZ20] Elliot Creager, David Madras, Toniann Pitassi, and Richard Zemel.  Causal modeling for fairness in dynamical systems.  In International Conference on Machine Learning, pages 2185–2195. PMLR, 2020.
  * [DF18] Julia Dressel and Hany Farid.  The accuracy, fairness, and limits of predicting recidivism.  Science advances, 4(1), 2018.
  * [DHP+12] Cynthia Dwork, Moritz Hardt, Toniann Pitassi, Omer Reingold, and Richard Zemel.  Fairness through awareness.  In Proceedings of the 3rd Innovations in Theoretical Computer Science conference, pages 214–226, 2012.
  * [DOBD+18] Michele Donini, Luca Oneto, Shai Ben-David, John S Shawe-Taylor, and Massimiliano Pontil.  Empirical risk minimization under fairness constraints.  In Advances in Neural Information Processing Systems, pages 2791–2801, 2018.
  * [DTB18] Shayan Doroudi, Philip S Thomas, and Emma Brunskill.  Importance sampling for fair policy selection.  In Proceedings of the 27th International Joint Conference on Artificial Intelligence, pages 5239–5243, 2018.
  * [Eth20] Kawin Ethayarajh.  Is your classifier actually biased? measuring fairness under uncertainty with bernstein bounds.  In Proceedings of the 58th Annual Meeting of the Association for Computational Linguistics, ACL, pages 2914–2919, 2020.
  * [FFM+15] Michael Feldman, Sorelle A Friedler, John Moeller, Carlos Scheidegger, and Suresh Venkatasubramanian.  Certifying and removing disparate impact.  In Proceedings of the 21st ACM SIGKDD international conference on knowledge discovery and data mining, pages 259–268, 2015.
  * [Gei06] Peter Geibel.  Reinforcement learning for mdps with constraints.  In European Conference on Machine Learning, pages 646–653. Springer, 2006.
  * [GF15] Javier Garcıa and Fernando Fernández.  A comprehensive survey on safe reinforcement learning.  Journal of Machine Learning Research, 16(1):1437–1480, 2015.
  * [GJKR18] Stephen Gillen, Christopher Jung, Michael Kearns, and Aaron Roth.  Online learning with an unknown fairness metric.  In Advances in neural information processing systems, pages 2600–2609, 2018.
  * [GJVS20] Philips George John, Deepak Vijaykeerthy, and Diptikalyan Saha.  Verifying individual fairness in machine learning models.  In Jonas Peters and David Sontag, editors, Proceedings of the 36th Conference on Uncertainty in Artificial Intelligence (UAI), volume 124, pages 749–758. PMLR, 2020.
  * [HMS20] Safwan Hossain, Evi Micha, and Nisarg Shah.  Fair algorithms for multi-agent multi-armed bandits.  arXiv preprint arXiv:2007.06699, 2020.
  * [HNG19] Hoda Heidari, Vedant Nanda, and Krishna P. Gummadi.  On the long-term impact of algorithmic decision policies: Effort unfairness and feature segregation through social learning.  In Proceedings of the 36th International Conference on Machine Learning, ICML, volume 97, pages 2692–2701. PMLR, 2019.
  * [JJK+17] Shahin Jabbari, Matthew Joseph, Michael Kearns, Jamie Morgenstern, and Aaron Roth.  Fairness in reinforcement learning.  In International Conference on Machine Learning, pages 1617–1626. PMLR, 2017.
  * [JKMR16] Matthew Joseph, Michael Kearns, Jamie H Morgenstern, and Aaron Roth.  Fairness in learning: Classic and contextual bandits.  Advances in neural information processing systems, 29:325–333, 2016\.
  * [JS20] Yujia Jin and Aaron Sidford.  Efficiently solving MDPs with stochastic mirror descent.  In Proceedings of the 37th International Conference on Machine Learning, pages 4890–4900. PMLR, 2020.
  * [KMR16] Jon Kleinberg, Sendhil Mullainathan, and Manish Raghavan.  Inherent trade-offs in the fair determination of risk scores.  arXiv preprint arXiv:1609.05807, 2016.
  * [KMR17] Jon M. Kleinberg, Sendhil Mullainathan, and Manish Raghavan.  Inherent trade-offs in the fair determination of risk scores.  In 8th Innovations in Theoretical Computer Science Conference, ITCS, volume 67, pages 43:1–43:23, 2017.
  * [KS98] Michael J. Kearns and Satinder P. Singh.  Near-optimal reinforcement learning in polynominal time.  In Proceedings of the Fifteenth International Conference on Machine Learning (ICML), pages 260–268. PMLR, 1998.
  * [LLJ19] Fengjiao Li, Jia Liu, and Bo Ji.  Combinatorial sleeping bandits with fairness constraints.  IEEE Transactions on Network Science and Engineering, 2019.
  * [LP17] David A Levin and Yuval Peres.  Markov chains and mixing times, volume 107.  American Mathematical Soc., 2017.
  * [Mah96] Sridhar Mahadevan.  Average reward reinforcement learning: Foundations, algorithms, and empirical results.  Machine Learning, 22(1-3):159–195, 1996.
  * [NJLS09] Arkadi Nemirovski, Anatoli Juditsky, Guanghui Lan, and Alexander Shapiro.  Robust stochastic approximation approach to stochastic programming.  SIAM Journal on optimization, 19(4):1574–1609, 2009.
  * [Pau15] Daniel Paulin.  Concentration inequalities for markov chains by marton couplings and spectral methods.  Electronic Journal of Probability, 20:32 pp., 2015.
  * [PGNN20] Vishakha Patil, Ganesh Ghalme, Vineet Nair, and Y Narahari.  Achieving fairness in the stochastic multi-armed bandit problem.  In Proceedings of the AAAI Conference on Artificial Intelligence, pages 5379–5386, 2020.
  * [Put14] Martin L Puterman.  Markov decision processes: discrete stochastic dynamic programming.  John Wiley & Sons, 2014.
  * [RBPS20] Tom Ron, Omer Ben-Porat, and Uri Shalit.  Corporate social responsibility via multi-armed bandits.  In AI for Social Good Workshop, 2020.
  * [SB18] Richard S Sutton and Andrew G Barto.  Reinforcement learning: An introduction.  MIT press, 2018.
  * [SS11] Shai Shalev-Shwartz.  Online learning and online convex optimization.  Machine Learning, 4(2):107–194, 2011.
  * [TBKN20] Bahar Taskesen, Jose Blanchet, Daniel Kuhn, and Viet Anh Nguyen.  A statistical test for probabilistic fairness, 2020.
  * [TDCM12] Aviv Tamar, Dotan Di Castro, and Shie Mannor.  Policy gradients with variance related risk criteria.  In Proceedings of the 29th International Conference on Machine Learning, pages 1651–1658, 2012.
  * [Wan17] Mengdi Wang.  Primal-dual $\pi$ learning: Sample complexity and sublinear run time for ergodic markov decision problems.  arXiv preprint arXiv:1710.06100, 2017.
  * [WBT19] Min Wen, Osbert Bastani, and Ufuk Topcu.  Fairness with dynamics.  arXiv preprint arXiv:1901.08568, 2019.
  * [YH17] Sirui Yao and Bert Huang.  Beyond parity: Fairness objectives for collaborative filtering.  In Advances in Neural Information Processing Systems, pages 2921–2930, 2017.
  * [ZB18] Junzhe Zhang and Elias Bareinboim.  Fairness in decision-making—the causal explanation formula.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 32, 2018.
  * [ZCKM20] Tom Zahavy, Alon Cohen, Haim Kaplan, and Yishay Mansour.  Unknown mixing times in apprenticeship and reinforcement learning.  In Proceedings of the 36th Conference on Uncertainty in Artificial Intelligence (UAI), volume 124, pages 430–439. PMLR, 2020.
  * [ZL20] Xueru Zhang and Mingyan Liu.  Fairness in learning-based sequential decision algorithms: A survey.  arXiv preprint arXiv:2001.04861, 2020.
  * [ZS15] Chongjie Zhang and Julie Shah.  On fairness in decision-making under uncertainty: Definitions, computation, and comparison.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 29, 2015.
  * [ZTL+20] Xueru Zhang, Ruibo Tu, Yang Liu, Mingyan Liu, Hedvig Kjellström, Kun Zhang, and Cheng Zhang.  How do fair decisions fare in long-term qualification?  In Advances in Neural Information Processing Systems 33: Annual Conference on Neural Information Processing Systems 2020, NeurIPS, 2020.
  * [ZVRG17] Muhammad Bilal Zafar, Isabel Valera, Manuel Gomez Rogriguez, and Krishna P Gummadi.  Fairness constraints: Mechanisms for fair classification.  In Artificial Intelligence and Statistics, pages 962–970. PMLR, 2017.

##  Appendix A Notation Summary

Notation | Meaning | Notation | Meaning
---|---|---|---
$\mathcal{S}$ | State space | $\boldsymbol{\mu}$ | Fair-LP (P) variable
$\mathcal{A}$ | Action space | $\boldsymbol{\lambda}$ | Fair-LP (P) variable
$\mathbf{\Gamma}$ | Transition probability matrix | $\mathbf{C}$ | Constraint matrix in LP
$\boldsymbol{r}$ | Reward vector | $\mathbf{D}_{\boldsymbol{\rho}}$ |  $n\times n$ diagonal matrix with entries as $1/\rho_{i}$
$\boldsymbol{\rho}$ | Fairness vector | $d_{\boldsymbol{\rho}}$ | $\max_{s}\frac{n}{1-\rho_{s}n}$
$n$ | Number of states $|\mathcal{S}|$ | $M$ | $2t_{\text{mix}}(1+d_{\boldsymbol{\rho}})$
$m$ | Number of actions $|\mathcal{A}|$ | $\Delta_{\boldsymbol{\rho}}^{\ell}$ | Fair simplex
$\ell$ | Number of state-action pairs $nm$ | $f(.)$ | Optimization function
$\pi$ | Policy | $\texttt{gap}(\texttt{Gap})$ | Gap function
$\Pi$ | Matrix representation of policy $\pi$ | $V_{\boldsymbol{x}}(\boldsymbol{y})$ | Distance-generating function
$\boldsymbol{\nu}^{\pi}$ | Stationary distribution of policy $\pi$ | $\boldsymbol{g}$ | Gradient vector
$t_{\text{mix}}$ | Mixing-time of MDP | $v$ | Gradient estimate moment bound
$\boldsymbol{x}$ | Fair-LP (D) variable/state-action visitation frequency | $\eta$ | Step size
|  | $T$ | Number of gradient steps

Table 2: The notation used in the paper.