# Model-Free Robust Average-Reward Reinforcement Learning

Yue Wang  Alvaro Velasquez  George Atia  Ashley Prater-Bennette  Shaofeng Zou

###### Abstract

Robust Markov decision processes (MDPs) address the challenge of model
uncertainty by optimizing the worst-case performance over an uncertainty set
of MDPs. In this paper, we focus on the robust average-reward MDPs under the
model-free setting. We first theoretically characterize the structure of
solutions to the robust average-reward Bellman equation, which is essential
for our later convergence analysis. We then design two model-free algorithms,
robust relative value iteration (RVI) TD and robust RVI Q-learning, and
theoretically prove their convergence to the optimal solution. We provide
several widely used uncertainty sets as examples, including those defined by
the contamination model, total variation, Chi-squared divergence, Kullback-
Leibler (KL) divergence and Wasserstein distance.

Robust Bellman equation,

##  1 Introduction

Reinforcement learning (RL) has demonstrated remarkable success in
applications like robotics, finance, and games. However, RL often suffers from
a severe performance degradation when deployed in real-world environments,
which is often due to the model mismatch between the training and testing
environments. Such a model mismatch could be a result of non-stationary
conditions, modeling errors between simulated and real systems, external
perturbations and adversarial attacks. To address this challenge, a framework
of robust MDPs and robust RL was developed in (Nilim & El Ghaoui, 2004;
Iyengar, 2005; Bagnell et al., 2001), where an uncertainty set of MDPs is
constructed, and a pessimistic approach is adopted to optimize the worst-case
performance over the uncertainty set. Such a minimax approach guarantees
performance for all MDPs within the uncertainty set, making it robust to model
mismatch.

Two performance criteria are commonly used for infinite-horizon MDPs: 1) the
discounted-reward setting, where the reward is discounted exponentially with
time; and 2) the average-reward setting, where the long-term average-reward
over time is of interest. For systems that operate for an extended period of
time, e.g., queue control, inventory management in supply chains, or
communication networks, it is more important to optimize the average-reward
since policies obtained from the discounted-reward setting may be myopic and
have poor long-term performance (Kazemi et al., 2022). However, much of the
work on robust MDPs has focused on the discounted-reward setting, and robust
average-reward MDPs remain largely unexplored.

Robust MDPs under the two criteria are fundamentally different. Compared to
the discounted-reward setting, the average-reward setting places equal weight
on immediate and future rewards, thus depends on the limiting state-action
frequencies of the underlying MDPs. This makes it more challenging to
investigate than the discounted-reward setting. For example, to establish a
contraction in the discounted-reward setting, it generally suffices to have a
discount factor strictly less than one, whereas in the average-reward setting,
convergence depends on the structure of the MDP. Studies on robust average-
reward MDPs are quite limited in the literature, e.g., (Tewari & Bartlett,
2007; Lim et al., 2013; Wang et al., 2023), and are model-based, where the
uncertainty set is fully known to the learner. However, the more practical
model-free setting, where only samples from the nominal MDP (the centroid of
the uncertainty set) can be observed, has yet to be explored. Algorithms and
fundamental convergence guarantees in this setting have not been established
yet, which is the focus of this paper.

###  1.1 Challenges and Contributions

In this paper, we develop a fundamental characterization of solutions to the
robust average-reward Bellman equation, design model-free algorithms with
provable optimality guarantees, and construct an unbiased estimator of the
robust average-reward Bellman operator for various types of uncertainty sets.
Our results fill the gap in model-free algorithms for robust average-reward
RL, and are fundamentally different from existing studies on robust discounted
RL and robust and non-robust average-reward MDPs. In the following, we
summarize the major challenges and our key contributions.

We characterize the fundamental structure of solutions to the robust average-
reward Bellman equation. Value-based approaches typically converge to a
solution to the (robust) Bellman equation. For example, TD and Q-learning
converge to the unique fixed-point that yields an optimal solution to the
Bellman equation in the discounted setting (Sutton & Barto, 2018). However, in
the (robust) average-reward setting, solutions to the (robust) average-reward
Bellman equation may not be unique (Puterman, 1994; Wang et al., 2023), and
the fundamental structure of the solution space usually plays an important
role in proving convergence (e.g., (Wan et al., 2021; Abounadi et al., 2001;
Tsitsiklis & Van Roy, 1999)). In the non-robust average-reward setting, the
solution to the average-reward Bellman equation is unique up to an additive
constant (Puterman, 1994). This result is largely due to the nice linear
structure of the non-robust average-reward Bellman equation. However, in the
robust setting, where the worst-case performance is of interest, the robust
average-reward Bellman equation is non-linear, which poses a significant
challenge in characterizing its solution. We demonstrate that the worst-case
transition kernel may not be unique and then show that any solution to the
robust average-reward Bellman equation is the relative value function for some
worst-case transition kernel (up to an additive constant). We note that this
structure is of key importance later in our convergence analysis.

We develop model-free approaches for robust policy evaluation and optimal
control. We then design model-free algorithms for robust average-reward RL. To
ensure stability of the iterates, we adopt a similar idea as in the non-robust
average-reward setting (Puterman, 1994; Bertsekas, 2011), which is to subtract
an offset function and design robust RVI TD and Q-learning algorithms. Unlike
the convergence proof for model-free algorithms for non-robust average-reward
MDPs, e.g., (Wan et al., 2021; Abounadi et al., 2001) and robust discounted
MDPs, e.g., (Wang & Zou, 2021; Liu et al., 2022), where the Bellman operator
is either linear or a contraction, our robust average-reward Bellman operator
is neither linear nor a contraction, which makes the convergence analysis
challenging. Nevertheless, based on the solution structure of the robust
average-reward Bellman equation discussed above, we prove that our algorithms
converge to solutions to the robust average-reward Bellman equations.
Specifically, for robust RVI TD, its output converges to the worst-case
average-reward; and for the robust RVI Q-learning, the greedy policy w.r.t.
the Q-function converges to an optimal robust policy.

We construct unbiased estimators with bounded variance for robust Bellman
operators under various uncertainty set models. One popular approach in RL is
bootstrapping, i.e., use $\mathbf{H}Q$ as an estimate of the true
$Q$-function, where $\mathbf{H}$ is the (robust) Bellman operator. In the non-
robust setting, an unbiased estimate of $\mathbf{H}$ can be easily derived
because $\mathbf{H}$ is linear in the transition kernel. This linearity,
however, does not hold in the robust setting since we minimize over the
transition kernel to account for the worst-case and, therefore, directly
plugging in the empirical transition kernel results in a biased estimate. In
this paper, we employ the multi-level Monte-Carlo method (Blanchet & Glynn,
2015) and construct unbiased estimators for five popular uncertainty models,
including those defined by contamination models, total variation, Chi-squared
divergence, Kullback-Leibler (KL) divergence, and Wasserstein distance. We
then prove that our estimator is unbiased and has a bounded variance. These
properties play an important role in establishing the convergence of our
algorithms.

###  1.2 Related Work

Robust average-reward MDPs. Studies on robust average-reward MDPs are quite
limited in the literature. Model-based robust average-reward MDPs were first
studied in (Tewari & Bartlett, 2007), where the authors showed the existence
of the Blackwell optimal policy for a specific uncertainty set of bounded
parameters. Results for general uncertainty sets were developed in recent work
(Wang et al., 2023), following a fundamentally different approach from (Tewari
& Bartlett, 2007). However, these two methods are model-based. The work of
(Lim et al., 2013) designed a model-free algorithm for the uncertainty set
defined by total variation and characterized its regret bound. However, their
method cannot be generalized to other uncertainty sets. In this paper, we
design the first model-free method for general uncertainty sets and provide
fundamental insights into the model-free robust average-reward RL problem.

Robust discounted-reward MDPs. Model-based methods for robust discounted MDPs
were studied in (Iyengar, 2005; Nilim & El Ghaoui, 2004; Bagnell et al., 2001;
Satia & Lave Jr, 1973; Wiesemann et al., 2013; Lim & Autef, 2019; Xu & Mannor,
2010; Yu & Xu, 2015; Lim et al., 2013; Tamar et al., 2014; Neufeld & Sester,
2022), where the uncertainty set is assumed to be known, and the problem can
be solved using robust dynamic programming. Later, the studies were
generalized to the model-free setting (Roy et al., 2017; Badrinath & Kalathil,
2021; Wang & Zou, 2021, 2022; Tessler et al., 2019; Liu et al., 2022; Zhou et
al., 2021; Yang et al., 2021; Panaganti & Kalathil, 2021; Goyal & Grand-
Clement, 2018; Kaufman & Schaefer, 2013; Ho et al., 2018, 2021). Our work
focuses on the average-reward setting. First, the robust average-reward
Bellman operator does not come with a simple contraction property like the one
for the discounted setting, and the average-reward depends on the limiting
behavior of the underlying MDP. Moreover, the average-reward Bellman function
is a function of two variables: the average-reward and the relative value
function, whereas its discounted-reward counterpart is a function of only the
value function. These challenges make the robust average-reward problem more
intricate. Furthermore, existing studies mostly focus on a certain type of
uncertainty set, e.g., contamination (Wang & Zou, 2021) and total variation
(Panaganti et al., 2022); in this paper, our method can be used to solve a
wide range of uncertainty sets.

Non-robust average-reward MDPs. Early contributions to non-robust average-
reward MDPs include a fundamental characterization of the problem and model-
based methods (Puterman, 1994; Bertsekas, 2011). Model-free methods in the
tabular setting, e.g., RVI Q-learning (Abounadi et al., 2001) and differential
Q-learning (Wan et al., 2021; Wan & Sutton, 2022), were developed recently and
are both shown to converge to the optimal average-reward. There is also work
on average-reward RL with function approximation, e.g., (Zhang et al., 2021b;
Tsitsiklis & Van Roy, 1999; Zhang et al., 2021a; Yu & Bertsekas, 2009). In
this paper, we focus on the robust setting, where the key challenge lies in
the non-linearity of the robust average-reward Bellman equation, whereas it is
linear in the non-robust setting.

##  2 Preliminaries and Problem Formulation

Average-reward MDPs. An MDP $(\mathcal{S},\mathcal{A},\mathsf{P},r)$ is
specified by: a state space $\mathcal{S}$, an action space $\mathcal{A}$, a
transition kernel
$\mathsf{P}=\left\\{\mathsf{P}^{a}_{s}\in\Delta(\mathcal{S}),a\in\mathcal{A},s\in\mathcal{S}\right\\}$111$\Delta(\mathcal{S})$
denotes the $(|\mathcal{S}|-1)$-dimensional probability simplex on
$\mathcal{S}$. , where $\mathsf{P}^{a}_{s}$ is the distribution of the next
state over $\mathcal{S}$ upon taking action $a$ in state $s$, and a reward
function $r:\mathcal{S}\times\mathcal{A}\to[0,1]$. At each time step $t$, the
agent at state $s_{t}$ takes an action $a_{t}$, the environment then
transitions to the next state $s_{t+1}\sim\mathsf{P}^{a_{t}}_{s_{t}}$, and
provides a reward signal $r_{t}\in[0,1]$.

A stationary policy $\pi:\mathcal{S}\to\Delta(\mathcal{A})$ maps a state to a
distribution over $\mathcal{A}$, following which the agent takes action $a$ at
state $s$ with probability $\pi(a|s)$. Under a transition kernel $\mathsf{P}$,
the average-reward of $\pi$ starting from $s\in\mathcal{S}$ is defined as

$\displaystyle
g_{\mathsf{P}}^{\pi}(s)\triangleq\lim_{T\to\infty}\mathbb{E}_{\pi,\mathsf{P}}\bigg{[}\frac{1}{T}\sum^{T-1}_{n=0}r_{t}|S_{0}=s\bigg{]}.$
(1)

The relative value function is defined to measure the cumulative difference
between the reward and $g^{\pi}_{\mathsf{P}}$:

$\displaystyle
V^{\pi}_{\mathsf{P}}(s)\triangleq\mathbb{E}_{\pi,\mathsf{P}}\bigg{[}\sum^{\infty}_{t=0}(r_{t}-g^{\pi}_{\mathsf{P}})|S_{0}=s\bigg{]}.$
(2)

Then $(g^{\pi}_{\mathsf{P}},V^{\pi}_{\mathsf{P}})$ satisfies the following
Bellman equation (Puterman, 1994):

$\displaystyle
V^{\pi}_{\mathsf{P}}(s)=\mathbb{E}_{\pi,\mathsf{P}}\bigg{[}r(s,A)-g^{\pi}_{\mathsf{P}}(s)+\sum_{s^{\prime}\in\mathcal{S}}p^{A}_{s,s^{\prime}}V^{\pi}_{\mathsf{P}}(s^{\prime})\bigg{]}.$
(3)

Robust average-reward MDPs. For robust MDPs, the transition kernel is assumed
to be in some uncertainty set $\mathcal{P}$. At each time step, the
environment transits to the next state according to an arbitrary transition
kernel $\mathsf{P}\in\mathcal{P}$. In this paper, we focus on the
$(s,a)$-rectangular uncertainty set (Nilim & El Ghaoui, 2004; Iyengar, 2005),
i.e., $\mathcal{P}=\bigotimes_{s,a}\mathcal{P}^{a}_{s}$, where
$\mathcal{P}^{a}_{s}\subseteq\Delta(\mathcal{S})$. Popular uncertainty sets
include those defined by the contamination model (Huber, 1965; Wang & Zou,
2022), total variation (Lim et al., 2013), Chi-squared divergence (Iyengar,
2005), Kullback-Leibler (KL) divergence (Hu & Hong, 2013) and Wasserstein
distance (Gao & Kleywegt, 2022). We will investigate these uncertainty sets in
detail in Section 5.

We investigate the worst-case average-reward over the uncertainty set of MDPs.
Specifically, define the robust average-reward of a policy $\pi$ as

$\displaystyle
g^{\pi}_{\mathcal{P}}(s)\triangleq\min_{\kappa\in\bigotimes_{n\geq
0}\mathcal{P}}\lim_{T\to\infty}\mathbb{E}_{\pi,\kappa}\left[\frac{1}{T}\sum^{T-1}_{t=0}r_{t}|S_{0}=s\right],$
(4)

where $\kappa=(\mathsf{P}_{0},\mathsf{P}_{1}...)\in\bigotimes_{n\geq
0}\mathcal{P}$. It was shown in (Wang et al., 2023) that the worst case under
the time-varying model is equivalent to the one under the stationary model:

$\displaystyle
g^{\pi}_{\mathcal{P}}(s)=\min_{\mathsf{P}\in\mathcal{P}}\lim_{T\to\infty}\mathbb{E}_{\pi,\mathsf{P}}\left[\frac{1}{T}\sum^{T-1}_{t=0}r_{t}|S_{0}=s\right].$
(5)

Therefore, we limit our focus to the stationary model. We refer to the
minimizers of (5) as the worst-case transition kernels for the policy $\pi$,
and denote the set of all possible worst-case transition kernels by
$\Omega^{\pi}_{g}$, i.e.,
$\Omega^{\pi}_{g}\triangleq\\{\mathsf{P}\in\mathcal{P}:g^{\pi}_{\mathsf{P}}=g^{\pi}_{\mathcal{P}}\\}$.
As shown in Example A.1 in the appendix, the worst-case transition kernel may
not be unique.

We focus on the model-free setting, where only samples from the nominal MDP
(centroid of the uncertainty set) are available. We investigate two problems:
1) given a policy $\pi$, estimate its robust average-reward
$g^{\pi}_{\mathcal{P}}$, and 2) find an optimal robust policy that optimizes
the robust average-reward:

$\displaystyle\max_{\pi}g^{\pi}_{\mathcal{P}}(s),\text{ for any
}s\in\mathcal{S}.$ (6)

We denote by
$g^{*}_{\mathcal{P}}(s)\triangleq\max_{\pi}g^{\pi}_{\mathcal{P}}(s)$ the
optimal robust average-reward.

##  3 Robust RVI TD for Policy Evaluation

In this section, we study the problem of policy evaluation, which aims to
estimate the robust average-reward $g^{\pi}_{\mathcal{P}}$ for a fixed policy
$\pi$.

For technical convenience, we make the following assumption to guarantee that
the average-reward is independent of the initial state (Abounadi et al., 2001;
Wan et al., 2021; Zhang et al., 2021a; Zhang & Ross, 2021; Chen et al., 2022;
Wang et al., 2023).

######  Assumption 1.

The Markov chain induced by $\pi$ is a unichain for all
$\mathsf{P}\in\mathcal{P}$.

In general, the average-reward depends on the initial state. For example,
imagine a policy that induces a multichain in an MDP with two closed
communicating classes. A learning algorithm would be able to learn the
average-reward for each communicating class; however, the average-rewards for
the two classes may be different. To remove this complexity, it is common and
convenient to rule out this possibility. Under Assumption 1, the average-
reward w.r.t. any $\mathsf{P}\in\mathcal{P}$ is identical for any start state,
i.e., $g^{\pi}_{\mathsf{P}}(s)=g^{\pi}_{\mathsf{P}}(s^{\prime}),\forall
s,s^{\prime}\in\mathcal{S}$.

We first revisit the robust average-reward Bellman equation in (Wang et al.,
2023), and further characterize the structure of its solutions. For
$V\in\mathbb{R}^{|\mathcal{S}|}$, denote by
$\sigma_{\mathcal{P}^{a}_{s}}(V)\triangleq\min_{p\in\mathcal{P}^{a}_{s}}pV$,
$\mathsf{P}_{V}(s,a)\triangleq\arg\min_{p\in\mathcal{P}^{a}_{s}}pV$ and let
$\mathsf{P}_{V}=\\{\mathsf{P}_{V}(s,a),s\in\mathcal{S},a\in\mathcal{A}\\}$.

######  Theorem 3.1 (Robust Bellman equation).

If $(g,V)$ is a solution to the robust Bellman equation

$\displaystyle
V(s)=\sum_{a}\pi(a|s)(r(s,a)-g+\sigma_{\mathcal{P}^{a}_{s}}(V)),\forall s,$
(7)

then 1) $g=g^{\pi}_{\mathcal{P}}$ (Wang et al., 2023); 2)
$\mathsf{P}_{V}\in\Omega^{\pi}_{g}$; 3) $V=V^{\pi}_{\mathsf{P}_{V}}+ce$ for
some $c\in\mathbb{R}$, where $e$ denotes the vector
$(1,1,...,1)\in\mathbb{R}^{|\mathcal{S}|}$.

The robust Bellman equation in (7) was initially developed in (Wang et al.,
2023). The authors also defined a robust relative value function:

$\displaystyle
V_{\mathcal{P}}^{\pi}(s)\triangleq\min_{\mathsf{P}\in\mathcal{P}}\mathbb{E}_{\pi,\mathsf{P}}\left[\sum^{\infty}_{n=0}(r_{n}-g^{\pi}_{\mathcal{P}})|S_{0}=s\right],$
(8)

which is the worst-case relative value function, and showed that
$(g^{\pi}_{\mathcal{P}},V_{\mathcal{P}}^{\pi})$ is a solution to (7). However,
conversely, it may not be the case that any solution $(g,V)$ to (7) can be
written as $V=V_{\mathcal{P}}^{\pi}+ce$ for some $c\in\mathbb{R}$. This is in
contrast to the results for the non-robust average-reward setting, where the
set of solutions to the non-robust average-reward Bellman equation can be
written as
$\\{(g^{\pi}_{\mathsf{P}},V^{\pi}_{\mathsf{P}}+ce):c\in\mathbb{R}\\}$
(Puterman, 1994).

In Theorem 3.1, we show that for any solution $(g,V)$ to (7), the transition
kernel $P_{V}\in\Omega^{\pi}_{g}$, i.e., it is a worst-case transition kernel
for $g^{\pi}_{\mathcal{P}}$. Moreover, any solution to (7) can be written as
$V=V^{\pi}_{\mathsf{P}_{V}}+ce$ for some constant $c$. This is different from
the non-robust setting, as $V^{\pi}_{\mathsf{P}_{V}}$ actually depends on $V$.
As will be seen later, this result is crucial to establish the convergence of
our robust RVI TD.

Theorem 3.1 also reveals the fundamental difference between the robust and the
non-robust average-reward settings. Under the non-robust setting, the solution
set to the Bellman equation can be written as
$\\{(g^{\pi}_{\mathsf{P}},V^{\pi}_{\mathsf{P}}+ce):c\in\mathbb{R}\\}$. The
solution is uniquely determined by the transition kernel (up to some constant
vector $ke$). In contrast, in the robust setting, the robust Bellman equation
is no longer linear. Any solution $V$ to (7) is a relative value function
w.r.t. some worst-case transition kernel $\mathsf{P}\in\Omega^{\pi}_{g}$ (up
to some additive constant vector), i.e.,
$V\in\\{V^{\pi}_{\mathsf{P}}+ce:\mathsf{P}\in\Omega^{\pi}_{g},c\in\mathbb{R}\\}$.

A natural question that arises is whether, for any
$\mathsf{P}\in\Omega^{\pi}_{g}$,
$(g^{\pi}_{\mathcal{P}},V^{\pi}_{\mathsf{P}})$ is a solution to (7)? Lemma 3.1
refutes this.

######  Lemma 3.1.

There exists a robust MDP such that for some $\mathsf{P}\in\Omega^{\pi}_{g}$,
$(g^{\pi}_{\mathcal{P}},V^{\pi}_{\mathsf{P}})$ is not a solution to (7).

Lemma 3.1 implies that the solution set to (7) is a subset of
$\\{V^{\pi}_{\mathsf{P}}+ce,\mathsf{P}\in\Omega^{\pi}_{g},c\in\mathbb{R}\\}$.
Note that explicit characterization of the solution set to (7) is challenging
due to its non-linear structure; however, result 3 in Theorem 3.1 suffices for
the convergence proof (see appendix for proofs).

Motivated by the robust Bellman equation in (7), we propose a model-free
robust RVI TD algorithm in Algorithm 1, where $\hat{\mathbf{T}}$ and function
$f$ will be discussed later.

Algorithm 1 Robust RVI TD

Input: $V_{0},\alpha_{n},n=0,1,...,N-1$

1: for $n=0,1,...,N-1$ do

2: for all $s\in\mathcal{S}$ do

3: $V_{n+1}(s)\leftarrow
V_{n}(s)+\alpha_{n}(\hat{\mathbf{T}}V_{n}(s)-f(V_{n})-V_{n}(s))$

4: end for

5: end for

Note that (7) can be written as $V=\mathbf{T}V-g$, where $\mathbf{T}$ is the
robust average-reward Bellman operator. Since in the model-free setting
$\mathcal{P}$ is unknown, in Algorithm 1, we construct $\hat{\mathbf{T}}V$ as
an estimate of $\mathbf{T}V$ satisfying

$\displaystyle\mathbb{E}[\hat{\mathbf{T}}V]=\mathbf{T}V,\quad\text{Var}[\hat{\mathbf{T}}V(s)]\leq
C(1+\|V\|^{2}),$ (9)

for some constant $C>0$. In this paper, if not specified, $\|\cdot\|$ denotes
the infinity norm $\|\cdot\|_{\infty}$.

It is challenging to construct such $\hat{\mathbf{T}}$ as $\mathbf{T}$ is non-
linear in the nominal transition kernel from which samples are generated. In
Section 5, we will present in detail how to construct such $\hat{\mathbf{T}}$
for various uncertainty set models.

To make the iterates stable, we follow the idea of RVI in the non-robust
setting and introduce an offset function $f$ satisfying the following
assumption (Puterman, 1994).

######  Assumption 2.

$f:\mathbb{R}^{|\mathcal{S}|}\to\mathbb{R}$ is $L_{f}$-Lipschitz and satisfies
$f(e)=1,f(x+ce)=f(x)+c,f(cx)=cf(x)$, $\forall c\in\mathbb{R}$.

Assumption 2 can be easily satisfied, e.g., $f(V)=V(s_{0})$ for some reference
state $s_{0}\in\mathcal{S}$, and $f(V)=\frac{\sum_{s}V(s)}{|\mathcal{S}|}$
(Abounadi et al., 2001). Compared with the discounted setting, $f$ is critical
here. As we discussed above, in the average-reward setting, the solution to
the Bellman equation $V+ce$ can be arbitrarily large because $c$ can be any
real number. This may lead to a non-convergent sequence $V_{n}$ (see, e.g.,
example 8.5.2 of (Puterman, 1994)). Hence, a function $f$ is introduced to
"offset" $V_{n}$ and keep the iterates stable. Also, $f(V_{n})$ serves as an
estimator of the average-reward $g^{\pi}_{\mathcal{P}}$, as we shall see
later.

We then assume the Robbins-Monro condition on the stepsize, and further show
the convergence of robust RVI TD.

######  Assumption 3.

The stepsize $\\{\alpha_{n}\\}_{n=0}^{\infty}$ satisfies the Robbins-Monro
condition, i.e.,
$\sum^{\infty}_{n=0}\alpha_{n}=\infty,\quad\sum^{\infty}_{n=0}\alpha^{2}_{n}<\infty.$

######  Theorem 3.2 (Convergence of robust RVI TD).

Under Assumptions 1,2,3, and if $\hat{\mathbf{T}}$ satisfies (9), then almost
surely, $(f(V_{n}),V_{n})$ converges to a solution to (7) which may depend on
the initialization.

The result implies that $f(V_{n})\to g^{\pi}_{\mathcal{P}}$ a.s., which means
our robust RVI TD converges to the worst-case average-reward for the given
policy $\pi$.

######  Remark 3.1.

Our robust RVI TD algorithm is shown to converge to a solution to (7). This
model-free result is the same as the model-based results in (Wang et al.,
2023). Though it was shown in (Wang et al., 2023) that
$(g^{\pi}_{\mathcal{P}},V_{\mathcal{P}}^{\pi})$ (defined in (8)) is a solution
to (7), it is not guaranteed that the convergence is to
$(g^{\pi}_{\mathcal{P}},V_{\mathcal{P}}^{\pi}+ce)$ for some $c$.

######  Remark 3.2.

As we discussed following the statement of Theorem 3.1, result 3) of Theorem
3.1 is crucial to the convergence proof of Theorem 3.2. Specifically, it is
necessary in order to characterize the equilibrium of the associated ODE, and
thus the limit of the iterates $f(V_{n})\to g^{\pi}_{\mathcal{P}}$.

##  4 Robust RVI Q-Learning for Control

In this section, we study the problem of optimal control, which aims to find a
policy that optimizes the robust average-reward:
$\pi^{*}=\arg\max_{\pi}g^{\pi}_{\mathcal{P}}.$

Similar to Assumption 1, we make the following assumption to guarantee that
the average-reward is independent of the initial state (Abounadi et al., 2001;
Li et al., 2022).

######  Assumption 4.

The Markov chain induced by any $\mathsf{P}\in\mathcal{P}$ and any $\pi$ is a
unichain.

We first revisit the optimal robust Bellman equation in (Wang et al., 2023)
and further present a characterization of its solutions. Consider a
$Q$-function $Q:\mathcal{S}\times\mathcal{A}\rightarrow\mathbb{R}$, and define
$V_{Q}(s)=\max_{a}Q(s,a),\forall s\in\mathcal{S}$.

######  Lemma 4.1 (Optimal robust Bellman equation).

If $(g,Q)$ is a solution to the optimal robust Bellman equation

$\displaystyle Q(s,a)=r(s,a)-g+\sigma_{\mathcal{P}^{a}_{s}}(V_{Q}),\forall
s,a,$ (10)

then 1) $g=g^{*}_{\mathcal{P}}$ (Wang et al., 2023); 2) the greedy policy
w.r.t. $Q$: $\pi_{Q}(s)=\arg\max_{a}Q(s,a)$ is an optimal robust policy (Wang
et al., 2023); 3) $V_{Q}=V^{\pi_{Q}}_{\mathsf{P}}+ce$ for some
$\mathsf{P}\in\Omega^{\pi_{Q}}_{g},c\in\mathbb{R}$.

According to result 2 in Lemma 4.1, finding a solution to (10) is sufficient
to get the optimal robust average-reward and to derive the optimal robust
policy. We note that a complete characterization of the solution set to (10)
can be obtained similarly to the one in result 3 of Theorem 3.1. Here, we only
provide its structure to simplify the presentation and avoid cumbersome
notation. This structure as outlined in Lemma 4.1 is sufficient for our
convergence analysis.

We hence present the following model-free robust RVI Q-learning algorithm.

Algorithm 2 Robust RVI Q-learning

Input: $Q_{0},\alpha_{n}$

1: for $n=0,...,N-1$ do

2: for all $s\in\mathcal{S},a\in\mathcal{A}$ do

3: $Q_{n+1}(s,a)\leftarrow
Q_{n}(s,a)+\alpha_{n}\big{(}\hat{\mathbf{H}}Q_{n}(s,a)-f(Q_{n})-Q_{n}(s,a)\big{)}$

4: end for

5: end for

Similar to the robust RVI TD algorithm, denote the optimal robust Bellman
operator by ${\mathbf{H}}Q(s,a)\triangleq
r(s,a)+\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})$, and we construct an estimate
$\hat{\mathbf{H}}$ such that for some finite constant $C$,

$\displaystyle\mathbb{E}[\hat{\mathbf{H}}Q]$
$\displaystyle=\mathbf{H}Q,\quad\text{Var}[\hat{\mathbf{H}}Q(s,a)]\leq
C(1+\|Q\|^{2}).$ (11)

In Section 5, we will present in detail how to construct such
$\hat{\mathbf{H}}$ for various uncertainty set models.

The following theorem shows the convergence of the robust RVI Q-learning to
the optimal robust average-reward $g^{*}_{\mathcal{P}}$ and the optimal robust
policy $\pi^{*}$.

######  Theorem 4.1 (Convergence of robust RVI Q-learning).

Under Assumptions 2, 3 and 4, and if $\hat{\mathbf{H}}$ satisfies (11), then
almost surely, $(f(Q_{n}),Q_{n})$ converges to a solution to (10), i.e.,
$f(Q_{n})$ converges to $g^{*}_{\mathcal{P}}$, and the greedy policy
$\pi_{Q_{n}}(s)\triangleq\arg\max_{a}Q_{n}(s,a)$ converges to an optimal
robust average-reward $\pi^{*}$.

##  5 Case Studies

In the previous two sections, we showed that if an unbiased estimator with
bounded variance is available for the robust Bellman operator, then both
robust algorithms proposed converge to the optimum. In this section, we
present the design of these estimators for various uncertainty set models.

The major challenge in designing the estimated operators satisfying (9) and
(11) lies in estimating the support function $\sigma_{\mathcal{P}_{s}^{a}}(V)$
unbiasedly with bounded variance using samples from the nominal transition
kernel. However, the nominal transition kernel $\mathsf{P}$ in general is
different from the worst-case transition kernel, and the straightforward
estimator is shown to be biased. Specifically, if we centered at the empirical
transition kernel $\hat{\mathsf{P}}$ and construct the uncertainty set
$\hat{\mathcal{P}}$, then the estimator is biased:
$\mathbb{E}[\sigma_{{\hat{\mathcal{P}}}_{s}^{a}}(V)]\neq\sigma_{\mathcal{P}_{s}^{a}}(V).$

We consider several widely-used uncertainty models including the contamination
model, the total variation model, the Chi-square model, the KL-divergence
model, and the Wasserstein distance model. We show that our estimators are
unbiased and have bounded variance in the following theorem. We will present
the design in later sections.

######  Theorem 5.1.

For each uncertainty set, denote its corresponding estimators by
$\hat{\mathbf{T}}$ and $\hat{\mathbf{H}}$ as in Sections 5.1 and 5.2. Then,
there exists some constant $C$, such that (9) and (11) hold.

In the following sections, we construct an operator
$\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ to estimate the support function
$\sigma_{\mathcal{P}^{a}_{s}}$, $\forall s\in\mathcal{S},a\in\mathcal{A}$ for
each uncertainty set. We further define the estimated robust Bellman operators
as
$\hat{\mathbf{T}}V(s)\triangleq\sum_{a}\pi(a|s)(r(s,a)+\hat{\sigma}_{\mathcal{P}^{a}_{s}}(V))$
and $\hat{\mathbf{H}}Q(s,a)\triangleq
r(s,a)+\hat{\sigma}_{\mathcal{P}^{a}_{s}}(V_{Q}).$

###  5.1 Linear Model: Contamination Uncertainty Set

The $\delta$-contamination uncertainty set is
$\mathcal{P}^{a}_{s}=\\{(1-\delta)\mathsf{P}^{a}_{s}+\delta
q:q\in\Delta(\mathcal{S})\\},$ where $0<\delta<1$ is the radius. Under this
uncertainty set, the support function can be computed as
$\sigma_{\mathcal{P}^{a}_{s}}(V)=(1-\delta)\mathsf{P}^{a}_{s}V+\delta\min_{s}V(s),$
and this is linear in the nominal transition kernel $\mathsf{P}^{a}_{s}$. We
hence use the transition to the subsequent state to construct our estimator:

$\displaystyle\hat{\sigma}_{\mathcal{P}^{a}_{s}}(V)\triangleq(1-\delta)\gamma
V(s^{\prime})+\delta\min_{x}V(x),$ (12)

where $s^{\prime}$ is a subsequent state sample after $(s,a)$.

###  5.2 Non-Linear Models

Unlike the $\delta$-contamination model, most uncertainty sets result in a
non-linear support function of the nominal transition kernel. We will employ
the approach of multi-level Monte-Carlo which is widely used in quantile
estimation under stochastic environments (Blanchet & Glynn, 2015; Blanchet et
al., 2019; Wang & Wang, 2022) to construct an unbiased estimator with bounded
variance.

For any $s,a$, we first generate $N$ according to a geometric distribution
with parameter $\Psi\in(0,1)$. Then, we take action $a$ at state $s$ for
$2^{N+1}$ times, and observe $r(s,a)$ and the subsequent state
$\\{s_{i}^{\prime}\\},i=1,...,2^{N+1}$. We divide these $2^{N+1}$ samples into
two groups: samples with odd indices, and samples with even indices. We then
individually calculate the empirical distribution of $s^{\prime}$ using the
even-index samples, odd-index ones, all the samples, and the first sample:
$\hat{\mathsf{P}}^{a,E}_{s,N+1}=\frac{1}{2^{N}}\sum_{i=1}^{2^{N}}\mathbbm{1}_{s_{2i}^{\prime}},\quad\hat{\mathsf{P}}^{a,O}_{s,N+1}=\frac{1}{2^{N}}\sum_{i=1}^{2^{N}}\mathbbm{1}_{s_{2i-1}^{\prime}},$
$\hat{\mathsf{P}}^{a}_{s,N+1}=\frac{1}{2^{N+1}}\sum_{i=1}^{2^{N+1}}\mathbbm{1}_{s_{i}^{\prime}},\quad\hat{\mathsf{P}}^{a,1}_{s,N+1}=\mathbbm{1}_{s_{1}^{\prime}}.$
Then, we use these estimated transition kernels as nominal kernels to
construct four estimated uncertainty sets
$\hat{\mathcal{P}}^{a,E}_{s,N+1},\hat{\mathcal{P}}^{a,O}_{s,N+1},\hat{\mathcal{P}}^{a}_{s,N+1},\hat{\mathcal{P}}^{a,1}_{s,N+1}$.
The multi-level estimator is then defined as

$\displaystyle\hat{\sigma}_{\mathcal{P}^{a}_{s}}(V)\triangleq\sigma_{\hat{\mathcal{P}}^{a,1}_{s,N+1}}(V)+\frac{\Delta_{N}(V)}{p_{N}},$
(13)

where $p_{N}=\Psi(1-\Psi)^{N}$ and

$\displaystyle\Delta_{N}(V)\triangleq\sigma_{\hat{\mathcal{P}}^{a}_{s,N+1}}(V)-\frac{\sigma_{\hat{\mathcal{P}}^{a,E}_{s,N+1}}(V)+\sigma_{\hat{\mathcal{P}}^{a,O}_{s,N+1}}(V)}{2}.$

We note that in previous results of the multi-level Monte-Carlo estimator
(Blanchet & Glynn, 2015; Blanchet et al., 2019; Wang & Wang, 2022), several
assumptions are needed to show that the estimator is unbiased. These
assumptions, however, do not hold in our cases. For example, the function
$\sigma_{\mathcal{P}}(V)$ is not continuously differentiable. Hence, their
analysis cannot be directly applied here.

We then present four examples of non-linear uncertainty sets. Under each
example, a solution to the support function $\sigma_{\mathcal{P}}(V)$ is
given, and by plugging it into (13) the unbiased estimator can then be
constructed. More details can be found in Appendix D and Appendix E.

Total Variation Uncertainty Set. The total variation uncertainty set is
$\mathcal{P}^{a}_{s}=\\{q\in\Delta(|\mathcal{S}|):\frac{1}{2}\|q-\mathsf{P}^{a}_{s}\|_{1}\leq\delta\\},$
and the support function can be computed using its dual function (Iyengar,
2005):

$\displaystyle\sigma_{\mathcal{P}^{a}_{s}}(V)$ $\displaystyle=\max_{\mu\geq
0}\big{(}\mathsf{P}^{a}_{s}(V-\mu)-\delta\text{Span}(V-\mu)\big{)}.$ (14)

Chi-square Uncertainty Set. The Chi-square uncertainty set is
$\mathcal{P}^{a}_{s}=\\{q\in\Delta(|\mathcal{S}|):d_{c}(\mathsf{P}^{a}_{s},q)\leq\delta\\},$
where $d_{c}(q,p)=\sum_{s}\frac{(p(s)-q(s))^{2}}{p(s)}$. Its support function
can be computed using its dual function (Iyengar, 2005):

$\displaystyle\sigma_{\mathcal{P}^{a}_{s}}(V)=\max_{\mu\geq
0}\big{(}\mathsf{P}^{a}_{s}(V-\mu)-\sqrt{\delta\text{Var}_{\mathsf{P}^{a}_{s}}(V-\mu)}\big{)}.$
(15)

Kullback–Leibler (KL) Divergence Uncertainty Set. The KL-divergence between
two distributions $p,q$ is defined as
$D_{KL}(q||p)=\sum_{s}q(s)\log\frac{q(s)}{p(s)},$ and the uncertainty set
defined via KL divergence is

$\displaystyle\mathcal{P}_{s}^{a}=\left\\{q:D_{KL}(q||\mathsf{P}_{s}^{a})\leq\delta\right\\},\forall
s\in\mathcal{S},a\in\mathcal{A}.$ (16)

Its support function can be efficiently solved using the duality result in (Hu
& Hong, 2013):

$\displaystyle\sigma_{\mathcal{P}^{a}_{s}}(V)=-\min_{\alpha\geq
0}\left(\delta\alpha+\alpha\log\left(\mathbb{E}_{\mathsf{P}^{a}_{s}}\big{[}e^{\frac{-V}{\alpha}}\big{]}\right)\right).$
(17)

The above estimator for the KL-divergence uncertainty set has also been
developed in (Liu et al., 2022) for robust discounted MDPs. Its extension to
our average-reward setting is similar.

Wasserstein Distance Uncertainty Sets. Consider the metric space
$(\mathcal{S},d)$ by defining some distance metric $d$. For some parameter
$l\in[1,\infty)$ and two distributions $p,q\in\Delta(\mathcal{S})$, define the
$l$-Wasserstein distance between them as
$W_{l}(q,p)=\inf_{\mu\in\Gamma(p,q)}\|d\|_{\mu,l}$, where $\Gamma(p,q)$
denotes the distributions over $\mathcal{S}\times\mathcal{S}$ with marginal
distributions $p,q$, and
$\|d\|_{\mu,l}=\big{(}\mathbb{E}_{(X,Y)\sim\mu}\big{[}d(X,Y)^{l}\big{]}\big{)}^{1/l}$.
The Wasserstein distance uncertainty set is then defined as

$\displaystyle\mathcal{P}^{a}_{s}=\left\\{q\in\Delta(|\mathcal{S}|):W_{l}(\mathsf{P}^{a}_{s},q)\leq\delta\right\\}.$
(18)

To solve the support function w.r.t. the Wasserstein distance set, we first
prove the following duality lemma.

######  Lemma 5.1.

It holds that

$\displaystyle\sigma_{\mathcal{P}^{a}_{s}}(V)=\sup_{\lambda\geq
0}\left(-\lambda\delta^{l}+\mathbb{E}_{\mathsf{P}^{a}_{s}}\big{[}\inf_{y}\big{(}V(y)+\lambda
d(S,y)^{l}\big{)}\big{]}\right).$ (19)

Thus, the support function can be solved using its dual form, and the
estimator can then be constructed following (13).

##  6 Experiments

We numerically verify our previous convergence results and demonstrate the
robustness of our algorithms. Additional experiments can be found in Appendix
G.

###  6.1 Convergence of Robust RVI TD and Q-Learning

We first verify the convergence of our robust RVI TD and Q-learning algorithms
under a Garnet problem $\mathcal{G}(30,20)$ (Archibald et al., 1995). There
are $30$ states and $20$ actions. The nominal transition kernel
$\mathsf{P}=\\{\mathsf{P}^{a}_{s},s\in\mathcal{S},a\in\mathcal{A}\\}$ is
randomly generated by a normal distribution:
$\mathsf{P}^{a}_{s}\sim\mathcal{N}(1,\sigma^{a}_{s})$ and then normalized, and
the reward function $r(s,a)\sim\mathcal{N}(1,\mu^{a}_{s})$, where
$\mu^{a}_{s},\sigma^{a}_{s}\sim\textbf{Uniform}[0,100]$. We set $\delta=0.4$,
$\alpha_{n}=0.01$, $f(V)=\frac{\sum_{s}V(s)}{|\mathcal{S}|}$ and
$f(Q)=\frac{\sum_{s,a}Q(s,a)}{|\mathcal{S}||\mathcal{A}|}$. Due to the space
limit, we only show the results under the Chi-square and Wasserstein Distance
models. The results under the other three uncertainty sets are presented in
Appendix G.

For policy evaluation, we evaluate the robust average-reward of the uniform
policy $\pi(a|s)=\frac{1}{|\mathcal{A}|}$. We implement our robust RVI TD
algorithm under different uncertainty models. We run the algorithm
independently for $30$ times and plot the average value of $f(V)$ over all
$30$ trajectories. We also plot the 95th and 5th percentiles of the 30 curves
as the upper and lower envelopes of the curves. To compare, we plot the true
robust average-reward computed using the model-based robust value iteration
method in (Wang et al., 2023). It can be seen from the results in Figure 1
that our robust RVI TD algorithm converges to the true robust average-reward
value.

Figure 1: Robust RVI TD Algorithm.

We then consider policy optimization. We run our robust RVI Q-learning
independently for 30 times. The curves in Figure 2 show the average value of
$f(Q)$ over 30 trajectories, and the upper/lower envelopes are the 95/5
percentiles. We also plot the optimal robust average-reward
$g^{*}_{\mathcal{P}}$ computed by the model-based RVI method in (Wang et al.,
2023). Our robust RVI Q-learning converges to the optimal robust average-
reward $g^{*}_{\mathcal{P}}$ under each uncertainty set, which verifies our
theoretical results.

Figure 2: Robust RVI Q-Learning Algorithm.

###  6.2 Robustness of Robust RVI Q-Learning

We then demonstrate the robustness of our robust RVI Q-learning by showing
that our method achieves a higher average-reward when there is model deviation
between training and evaluation.

####  6.2.1 Recycling Robot

We first consider the recycling robot problem (Example 3.3 (Sutton & Barto,
2018)). A mobile robot running on a rechargeable battery aims to collect empty
soda cans. It has 2 battery levels: low and high. The robot can either 1)
search for empty cans; 2) remain stationary and wait for someone to bring it a
can; 3) go back to its home base to recharge. Under low (high) battery level,
the robot finds an empty can with probabilities $\alpha$ ($\beta$), and
remains at the same battery level. If the robot goes out to search but finds
nothing, it will run out of its battery and can only be carried back by human.
More details can be found in (Sutton & Barto, 2018).

In this experiment, the uncertainty lies in the probabilities $\alpha,\beta$
of finding an empty can if the robot chooses the action ‘search’. We set
$\delta=0.4$ and implement our algorithms and vanilla Q-learning under the
nominal environment ($\alpha=\beta=0.5$) with stepsize $0.01$. To show the
difference among the policies that the algorithms learned, we plot the
difference of $Q$ values at low battery level in Figure 3(a). In the low
battery level, the robust algorithms find conservative policies which choose
to wait instead of search, whereas the vanilla Q-learning finds a policy that
chooses to search. To test the robustness of the obtained policies, we
evaluate the average-reward of the learned policies in perturbed environments.
Specifically, let $x$ denote the amplitude of the perturbation. Then, we
estimate the worst performance of the two policies over the testing
uncertainty set $(0.5-x,0.5+x)$, and plot them in Figure 3(b). It can be seen
that when the perturbation is small, the true worst-case kernels (w.r.t.
$\delta$ during training) are far from the testing environment, and hence the
vanilla Q-learning has a higher reward; however, as the perturbation level
becomes larger, the testing environment gets closer to the worst-case kernels,
and then our robust algorithms perform better. It can be seen that the
performance of Q-learning decreases rapidly while our robust algorithm is
stable and outperforms the non-robust Q-learning. This implies that our
algorithm is robust to the model uncertainty.

(a) $Q(\text{low,wait})$-$Q(\text{low,search})$

(b) Perturbed environment

Figure 3: Recycling Robot.

####  6.2.2 Inventory Control Problem

We now consider the supply chain problem (Giannoccaro & Pontrandolfo, 2002;
Kemmer et al., 2018; Liu et al., 2022). At the beginning of each month, the
manager of a warehouse inspects the current inventory of a product. Based on
the current stock, the manager decides whether or not to order additional
stock from a supplier. During this month, if the customer demand is satisfied,
the warehouse can make a sale and obtain profits; but if not, the warehouse
will obtain a penalty associated with being unable to satisfy customer demand
for the product. The warehouse also needs to pay the holding cost for the
remaining stock and new items ordered. The goal is to maximize the average
profit.

We let $s_{t}$ denote the inventory at the beginning of the $t$-th month,
$D_{t}$ be a random demand during this month, and $a_{t}$ be the number of
units ordered by the manager. We assume that $D_{t}$ follows some distribution
and is independent over time. When the agent takes action $a_{t}$, the order
cost is $a_{t}$, and the holding cost is $3\cdot(s_{t}+a_{t})$. If the demand
$D_{t}\leq s_{t}+a_{t}$, then selling the item brings $5\cdot D_{t}$ in total;
but if the demand $D_{t}>s_{t}+a_{t}$, then there will not be any sale and a
penalty of $-15$ will be received. We set $\mathcal{S}=\\{0,1,...,16\\}$ and
$\mathcal{A}=\\{0,...,8\\}$.

We first set $\delta=0.4$ and $\alpha_{t}=0.01$, and implement our algorithms
and vanilla Q-learning under the nominal environment where
$D_{t}\sim\textbf{Uniform}(0,16)$ is generated following the uniform
distribution. To verify the robustness, we test the obtained policies under
different perturbed environments. More specifically, we perturb the
distribution of the demand to $D_{t}\sim U_{(m,b)}$, where

$U_{(m,b)}(x)=\left\\{\begin{array}[]{lr}\frac{1}{|\mathcal{S}|}+b\frac{|\mathcal{S}|-2}{2|\mathcal{S}|},\text{
if }x\in\\{m,m+1\\},&\\\ \frac{1-b}{|\mathcal{S}|},\text{
else}.&\end{array}\right.$

The results are plotted in Figure 4. We first fix $m=0$ and plot the
performance under different values of $b$ in Figure 4(a), then we fix $b=0.25$
and plot the performance under different values of $m$ in Figure 4(b).

As the results show, when $b$ is small, i.e., the perturbation of the
environment is small, the non-robust Q-learning obtains higher reward than our
robust methods; as $b$ becomes larger, the performance of the non-robust
method decreases rapidly, while our robust methods are more robust and
outperform the non-robust one. When $b$ is fixed, our robust methods
outperform the non-robust Q-learning, which also demonstrates the robustness
of our methods.

(a) $m=0$

(b) $b=0.25$

Figure 4: Inventory Control.

##  7 Conclusion

In this paper, we developed the first model-free algorithms with provable
convergence and optimality guarantee for robust average-reward RL under a
broad range of uncertainty set models. We characterized the fundamental
structure of solutions to the robust average-reward Bellman equation, which is
crucial for the convergence analysis. We designed model-free robust algorithms
based on the ideas of relative value iteration for non-robust average-reward
MDPs and the robust average-reward Bellman equation. We developed concrete
solutions to five popular uncertainty sets, where we generalized the idea of
multi-level Monte-Carlo and constructed an unbiased estimate of the non-linear
robust average-reward Bellman operator.

##  8 Acknowledgments

This work is supported by the National Science Foundation under Grants
CCF-2106560, CCF-2007783, CCF-2106339, and CCF-1552497. This material is based
upon work supported under the AI Research Institutes program by National
Science Foundation and the Institute of Education Sciences, U.S. Department of
Education through Award # 2229873 - National AI Institute for Exceptional
Education. Any opinions, findings and conclusions or recommendations expressed
in this material are those of the author(s) and do not necessarily reflect the
views of the National Science Foundation, the Institute of Education Sciences,
or the U.S. Department of Education.

## References

  * Abounadi et al. (2001) Abounadi, J., Bertsekas, D., and Borkar, V. S.  Learning algorithms for Markov decision processes with average cost.  _SIAM Journal on Control and Optimization_ , 40(3):681–698, 2001.
  * Archibald et al. (1995) Archibald, T., McKinnon, K., and Thomas, L.  On the generation of Markov decision processes.  _Journal of the Operational Research Society_ , 46(3):354–361, 1995.
  * Badrinath & Kalathil (2021) Badrinath, K. P. and Kalathil, D.  Robust reinforcement learning using least squares policy iteration with provable performance guarantees.  In _Proc. International Conference on Machine Learning (ICML)_ , pp. 511–520. PMLR, 2021.
  * Bagnell et al. (2001) Bagnell, J. A., Ng, A. Y., and Schneider, J. G.  Solving uncertain Markov decision processes.  2001\.
  * Bertsekas (2011) Bertsekas, D. P.  Dynamic Programming and Optimal Control 3rd edition, volume II.  _Belmont, MA: Athena Scientific_ , 2011.
  * Blanchet & Glynn (2015) Blanchet, J. H. and Glynn, P. W.  Unbiased Monte Carlo for optimization and functions of expectations via multi-level randomization.  In _2015 Winter Simulation Conference (WSC)_ , pp. 3656–3667. IEEE, 2015.
  * Blanchet et al. (2019) Blanchet, J. H., Glynn, P. W., and Pei, Y.  Unbiased multilevel Monte Carlo: Stochastic optimization, steady-state simulation, quantiles, and other applications.  _arXiv preprint arXiv:1904.09929_ , 2019.
  * Borkar (2009) Borkar, V. S.  _Stochastic approximation: a dynamical systems viewpoint_ , volume 48.  Springer, 2009.
  * Borkar & Soumyanatha (1997) Borkar, V. S. and Soumyanatha, K.  An analog scheme for fixed point computation. i. theory.  _IEEE Transactions on Circuits and Systems I: Fundamental Theory and Applications_ , 44(4):351–355, 1997.
  * Brockman et al. (2016) Brockman, G., Cheung, V., Pettersson, L., Schneider, J., Schulman, J., Tang, J., and Zaremba, W.  OpenAI Gym.  _arXiv preprint arXiv:1606.01540_ , 2016.
  * Chen et al. (2022) Chen, L., Jain, R., and Luo, H.  Learning infinite-horizon average-reward Markov decision processes with constraints.  _arXiv preprint arXiv:2202.00150_ , 2022.
  * Gao & Kleywegt (2022) Gao, R. and Kleywegt, A.  Distributionally robust stochastic optimization with Wasserstein distance.  _Mathematics of Operations Research_ , 2022.
  * Giannoccaro & Pontrandolfo (2002) Giannoccaro, I. and Pontrandolfo, P.  Inventory management in supply chains: a reinforcement learning approach.  _International Journal of Production Economics_ , 78(2):153–161, 2002.
  * Goyal & Grand-Clement (2018) Goyal, V. and Grand-Clement, J.  Robust Markov decision process: Beyond rectangularity.  _arXiv preprint arXiv:1811.00215_ , 2018.
  * Ho et al. (2018) Ho, C. P., Petrik, M., and Wiesemann, W.  Fast Bellman updates for robust MDPs.  In _Proc. International Conference on Machine Learning (ICML)_ , pp. 1979–1988. PMLR, 2018.
  * Ho et al. (2021) Ho, C. P., Petrik, M., and Wiesemann, W.  Partial policy iteration for L1-robust Markov decision processes.  _Journal of Machine Learning Research_ , 22(275):1–46, 2021.
  * Hu & Hong (2013) Hu, Z. and Hong, L. J.  Kullback-Leibler divergence constrained distributionally robust optimization.  _Available at Optimization Online_ , pp. 1695–1724, 2013.
  * Huber (1965) Huber, P. J.  A robust version of the probability ratio test.  _Ann. Math. Statist._ , 36:1753–1758, 1965.
  * Iyengar (2005) Iyengar, G. N.  Robust dynamic programming.  _Mathematics of Operations Research_ , 30(2):257–280, 2005.
  * Kaufman & Schaefer (2013) Kaufman, D. L. and Schaefer, A. J.  Robust modified policy iteration.  _INFORMS Journal on Computing_ , 25(3):396–410, 2013.
  * Kazemi et al. (2022) Kazemi, M., Perez, M., Somenzi, F., Soudjani, S., Trivedi, A., and Velasquez, A.  Translating omega-regular specifications to average objectives for model-free reinforcement learning.  In _Proceedings of the 21st International Conference on Autonomous Agents and Multiagent Systems_ , pp. 732–741, 2022.
  * Kemmer et al. (2018) Kemmer, L., von Kleist, H., de Rochebouët, D., Tziortziotis, N., and Read, J.  Reinforcement learning for supply chain optimization.  In _European Workshop on Reinforcement Learning_ , volume 14, 2018\.
  * Li et al. (2022) Li, T., Wu, F., and Lan, G.  Stochastic first-order methods for average-reward Markov decision processes.  _arXiv preprint arXiv:2205.05800_ , 2022.
  * Lim & Autef (2019) Lim, S. H. and Autef, A.  Kernel-based reinforcement learning in robust Markov decision processes.  In _Proc. International Conference on Machine Learning (ICML)_ , pp. 3973–3981. PMLR, 2019.
  * Lim et al. (2013) Lim, S. H., Xu, H., and Mannor, S.  Reinforcement learning in robust Markov decision processes.  In _Proc. Advances in Neural Information Processing Systems (NIPS)_ , pp. 701–709, 2013.
  * Liu et al. (2022) Liu, Z., Bai, Q., Blanchet, J., Dong, P., Xu, W., Zhou, Z., and Zhou, Z.  Distributionally robust $Q$-learning.  In _Proc. International Conference on Machine Learning (ICML)_ , pp. 13623–13643. PMLR, 2022.
  * Neufeld & Sester (2022) Neufeld, A. and Sester, J.  Robust ${Q}$-learning algorithm for Markov decision processes under Wasserstein uncertainty.  _arXiv preprint arXiv:2210.00898_ , 2022.
  * Nilim & El Ghaoui (2004) Nilim, A. and El Ghaoui, L.  Robustness in Markov decision problems with uncertain transition matrices.  In _Proc. Advances in Neural Information Processing Systems (NIPS)_ , pp. 839–846, 2004.
  * Panaganti & Kalathil (2021) Panaganti, K. and Kalathil, D.  Sample complexity of robust reinforcement learning with a generative model.  _arXiv preprint arXiv:2112.01506_ , 2021.
  * Panaganti et al. (2022) Panaganti, K., Xu, Z., Kalathil, D., and Ghavamzadeh, M.  Robust reinforcement learning using offline data.  _arXiv preprint arXiv:2208.05129_ , 2022.
  * Puterman (1994) Puterman, M. L.  Markov decision processes: Discrete stochastic dynamic programming, 1994\.
  * Roy et al. (2017) Roy, A., Xu, H., and Pokutta, S.  Reinforcement learning under model mismatch.  In _Proc. Advances in Neural Information Processing Systems (NIPS)_ , pp. 3046–3055, 2017.
  * Satia & Lave Jr (1973) Satia, J. K. and Lave Jr, R. E.  Markovian decision processes with uncertain transition probabilities.  _Operations Research_ , 21(3):728–740, 1973.
  * Sutton & Barto (2018) Sutton, R. S. and Barto, A. G.  _Reinforcement Learning: An Introduction_.  The MIT Press, Cambridge, Massachusetts, 2018.
  * Tamar et al. (2014) Tamar, A., Mannor, S., and Xu, H.  Scaling up robust MDPs using function approximation.  In _Proc. International Conference on Machine Learning (ICML)_ , pp. 181–189. PMLR, 2014.
  * Tessler et al. (2019) Tessler, C., Efroni, Y., and Mannor, S.  Action robust reinforcement learning and applications in continuous control.  In _Proc. International Conference on Machine Learning (ICML)_ , pp. 6215–6224. PMLR, 2019.
  * Tewari & Bartlett (2007) Tewari, A. and Bartlett, P. L.  Bounded parameter Markov decision processes with average reward criterion.  In _International Conference on Computational Learning Theory_ , pp. 263–277. Springer, 2007.
  * Tsitsiklis & Van Roy (1999) Tsitsiklis, J. N. and Van Roy, B.  Average cost temporal-difference learning.  _Automatica_ , 35(11):1799–1808, 1999.
  * Wan & Sutton (2022) Wan, Y. and Sutton, R. S.  On convergence of average-reward off-policy control algorithms in weakly-communicating MDPs.  _arXiv preprint arXiv:2209.15141_ , 2022.
  * Wan et al. (2021) Wan, Y., Naik, A., and Sutton, R. S.  Learning and planning in average-reward Markov decision processes.  In _Proc. International Conference on Machine Learning (ICML)_ , pp. 10653–10662. PMLR, 2021.
  * Wang & Wang (2022) Wang, G. and Wang, T.  Unbiased multilevel Monte Carlo methods for intractable distributions: Mlmc meets mcmc.  _arXiv preprint arXiv:2204.04808_ , 2022.
  * Wang & Zou (2021) Wang, Y. and Zou, S.  Online robust reinforcement learning with model uncertainty.  In _Proc. Advances in Neural Information Processing Systems (NeurIPS)_ , 2021.
  * Wang & Zou (2022) Wang, Y. and Zou, S.  Policy gradient method for robust reinforcement learning.  In _Proc. International Conference on Machine Learning (ICML)_ , volume 162, pp. 23484–23526. PMLR, 2022.
  * Wang et al. (2023) Wang, Y., Velasquez, A., Atia, G., Prater-Bennette, A., and Zou, S.  Robust average-reward Markov decision processes.  In _Proc. Conference on Artificial Intelligence (AAAI)_ , 2023.
  * Wiesemann et al. (2013) Wiesemann, W., Kuhn, D., and Rustem, B.  Robust Markov decision processes.  _Mathematics of Operations Research_ , 38(1):153–183, 2013.
  * Xu & Mannor (2010) Xu, H. and Mannor, S.  Distributionally robust Markov decision processes.  In _Proc. Advances in Neural Information Processing Systems (NIPS)_ , pp. 2505–2513, 2010.
  * Yang et al. (2021) Yang, W., Zhang, L., and Zhang, Z.  Towards theoretical understandings of robust Markov decision processes: Sample complexity and asymptotics.  _arXiv preprint arXiv:2105.03863_ , 2021.
  * Yu & Bertsekas (2009) Yu, H. and Bertsekas, D. P.  Convergence results for some temporal difference methods based on least squares.  _IEEE Transactions on Automatic Control_ , 54(7):1515–1531, 2009.
  * Yu & Xu (2015) Yu, P. and Xu, H.  Distributionally robust counterpart in Markov decision processes.  _IEEE Transactions on Automatic Control_ , 61(9):2538–2543, 2015.
  * Zhang et al. (2021a) Zhang, S., Wan, Y., Sutton, R. S., and Whiteson, S.  Average-reward off-policy policy evaluation with function approximation.  In _Proc. International Conference on Machine Learning (ICML)_ , pp. 12578–12588. PMLR, 2021a.
  * Zhang et al. (2021b) Zhang, S., Zhang, Z., and Maguluri, S. T.  Finite sample analysis of average-reward TD learning and $Q$-learning.  In _Proc. Advances in Neural Information Processing Systems (NeurIPS)_ , volume 34, pp. 1230–1242, 2021b.
  * Zhang & Ross (2021) Zhang, Y. and Ross, K. W.  On-policy deep reinforcement learning for the average-reward criterion.  In _Proc. International Conference on Machine Learning (ICML)_ , pp. 12535–12545. PMLR, 2021.
  * Zhou et al. (2021) Zhou, Z., Bai, Q., Zhou, Z., Qiu, L., Blanchet, J., and Glynn, P.  Finite-sample regret bound for distributionally robust offline tabular reinforcement learning.  In _Proc. International Conference on Artifical Intelligence and Statistics (AISTATS)_ , pp. 3331–3339. PMLR, 2021.

##  Appendix A Proof of Lemma 3.1

We construct the following example.

######  Example A.1.

Consider an MDP with 3 states (1,2,3) and only one action $a$, and set a
$(s,a)$-rectangular uncertainty set
$\mathcal{P}=\mathcal{P}^{a}_{1}\bigotimes\mathcal{P}^{a}_{2}\bigotimes\mathcal{P}^{a}_{3}$
where
$\mathcal{P}^{a}_{1}=\\{{\mathsf{P}^{a}_{1}}_{1},{\mathsf{P}^{a}_{1}}_{2}\\}$,
$\mathcal{P}^{a}_{2}=\\{(0,0,1)^{\top}\\}$ and
$\mathcal{P}^{a}_{3}=\\{(0,1,0)^{\top}\\}$, where
${\mathsf{P}^{a}_{1}}_{1}=(0,1,0)^{\top},{\mathsf{P}^{a}_{1}}_{2}=(0,0,1)^{\top}$.
Hence, the uncertainty set contains two transition kernels
$\mathcal{P}=\\{\mathsf{P}_{1},\mathsf{P}_{2}\\}$. The reward of each state is
set to be $r=(r_{1},r_{2},r_{3})$. The only stationary policy $\pi$ in this
example is $\pi(i)=a,\forall i$.

Note that this robust MDP is a unichain and hence satisfies Assumption 1 with
$g^{\pi}_{\mathsf{P}_{1}}(1)=g^{\pi}_{\mathsf{P}_{1}}(2)=g^{\pi}_{\mathsf{P}_{1}}(3),g^{\pi}_{\mathsf{P}_{2}}(1)=g^{\pi}_{\mathsf{P}_{2}}(2)=g^{\pi}_{\mathsf{P}_{2}}(3)$.

Under both transition kernels $\mathsf{P}_{1},\mathsf{P}_{2}$, the average-
reward are identical:
$g^{\pi}_{\mathsf{P}_{1}}=g^{\pi}_{\mathsf{P}_{2}}=0.5r_{2}+0.5r_{3}$. Hence,
both $\mathsf{P}_{1},\mathsf{P}_{2}$ are the worst-case transition kernels.

According to Section A.5 of (Puterman, 1994), the relative value functions
w.r.t. $\mathsf{P}_{1},\mathsf{P}_{2}$ can be computed as

$V^{\pi}_{\mathsf{P}_{1}}=\bigg{(}r_{1}-\frac{1}{4}r_{2}-\frac{3}{4}r_{3},\frac{1}{4}r_{2}-\frac{1}{4}r_{3},-\frac{1}{4}r_{2}+\frac{1}{4}r_{3}\bigg{)}^{\top},$
$V^{\pi}_{\mathsf{P}_{2}}=\bigg{(}r_{1}-\frac{3}{4}r_{2}-\frac{1}{4}r_{3},\frac{1}{4}r_{2}-\frac{1}{4}r_{3},-\frac{1}{4}r_{2}+\frac{1}{4}r_{3}\bigg{)}^{\top}.$

When $r_{3}>r_{2}$, only $V^{\pi}_{\mathsf{P}_{1}}$ is the solution to (7);
and when $r_{2}>r_{3}$, only $V^{\pi}_{\mathsf{P}_{2}}$ is the solution to
(7). Hence, this proves Lemma 3.1 and implies that not any relative value
function w.r.t. a worst-case transition kernel is a solution to (7).

##  Appendix B Robust RVI TD Method for Policy Evaluation

We define the following notation:

$\displaystyle r_{\pi}(s)$ $\displaystyle\triangleq\sum_{a}\pi(a|s)r(s,a),$
$\displaystyle\sigma_{\mathcal{P}_{s}}(V)$
$\displaystyle\triangleq\sum_{a}\pi(a|s)\sigma_{\mathcal{P}^{a}_{s}}(V),$
$\displaystyle\sigma_{\mathcal{P}}(V)$
$\displaystyle\triangleq(\sigma_{\mathcal{P}_{s_{1}}}(V),\sigma_{\mathcal{P}_{s_{2}}}(V),...,\sigma_{\mathcal{P}_{s_{|\mathcal{S}|}}}(V))\in\mathbb{R}^{|\mathcal{S}|}.$

###  B.1 Proof of Theorem 3.1

######  Theorem B.1 (Restatement of Theorem 3.1).

If $(g,V)$ is a solution to the robust Bellman equation

$\displaystyle
V(s)=\sum_{a}\pi(a|s)(r(s,a)-g+\sigma_{\mathcal{P}^{a}_{s}}(V)),\forall s,$
(20)

then 1) $g=g^{\pi}_{\mathcal{P}}$ (Wang et al., 2023); 2)
$\mathsf{P}_{V}\in\Omega^{\pi}_{g}$; 3) $V=V^{\pi}_{\mathsf{P}_{V}}+ce$ for
some $c\in\mathbb{R}$.

###### Proof.

1). The robust Bellman equation in (20) can be rewritten as

$\displaystyle g+V(s)-r_{\pi}(s)=\sigma_{\mathcal{P}_{s}}(V),\forall
s\in\mathcal{S}.$ (21)

From the definition, it follows that

$\displaystyle\sigma_{\mathcal{P}_{s}}(V)=\sum_{a}\pi(a|s)\min_{\mathsf{P}^{a}_{s}\in\mathcal{P}^{a}_{s}}\mathsf{P}^{a}_{s}V.$
(22)

Hence, for any transition kernel
$\mathsf{P}=(\mathsf{P}^{a}_{s})\in\bigotimes_{s,a}\mathcal{P}^{a}_{s}$,

$\displaystyle g+V(s)-r_{\pi}(s)-\sum_{a}\pi(a|s)\mathsf{P}^{a}_{s}V\leq
0,\forall s.$ (23)

It can be further rewritten in matrix form as:

$\displaystyle ge\leq r_{\pi}+(\mathsf{P}^{\pi}-I)V,$ (24)

where $\mathsf{P}^{\pi}$ is the state transition matrix induced by $\pi$ and
$\mathsf{P}$, i.e., the $s$-th row of $\mathsf{P}^{\pi}$ is

$\displaystyle\sum_{a}\pi(a|s)\mathsf{P}^{a}_{s}.$ (25)

Note that $\mathsf{P}^{\pi}$ has non-negative components since it is a
transition matrix. Multiplying by $\mathsf{P}^{\pi}$ on both sides, we have
that

$\displaystyle\mathsf{P}^{\pi}ge=ge$
$\displaystyle\leq\mathsf{P}^{\pi}r_{\pi}+\mathsf{P}^{\pi}(\mathsf{P}^{\pi}-I)V,$
$\displaystyle ge$
$\displaystyle\leq(\mathsf{P}^{\pi})^{2}r_{\pi}+(\mathsf{P}^{\pi})^{2}(\mathsf{P}^{\pi}-I)V,$
$\displaystyle...$ $\displaystyle ge$
$\displaystyle\leq(\mathsf{P}^{\pi})^{n-1}r_{\pi}+(\mathsf{P}^{\pi})^{n-1}(\mathsf{P}^{\pi}-I)V.$
(26)

Now, by summing up all these inequalities in (24) and (B.1), we have that

$\displaystyle
nge\leq\sum^{n-1}_{i=0}(\mathsf{P}^{\pi})^{i}r_{\pi}+((\mathsf{P}^{\pi})^{n}-I)V,$
(27)

and hence,

$\displaystyle
ge\leq\frac{\sum^{n-1}_{i=0}(\mathsf{P}^{\pi})^{i}r_{\pi}}{n}+\frac{((\mathsf{P}^{\pi})^{n}-I)V}{n}.$
(28)

Let $n\to\infty$, and we have that

$\displaystyle ge$
$\displaystyle\leq\lim_{n\to\infty}\frac{\sum^{n-1}_{i=0}(\mathsf{P}^{\pi})^{i}r_{\pi}}{n}+\lim_{n\to\infty}\frac{((\mathsf{P}^{\pi})^{n}-I)V}{n}$
$\displaystyle=g^{\pi}_{\mathsf{P}}e,$ (29)

where the last inequality is from the definition of $g^{\pi}_{\mathsf{P}}$ and
the fact that $\lim_{n\to\infty}\frac{((\mathsf{P}^{\pi})^{n}-I)V}{n}=0$.
Hence, $g\leq g^{\pi}_{\mathsf{P}}$ for any
$\mathsf{P}\in\bigotimes_{s,a}\mathcal{P}^{a}_{s}$.

Consider the worst-case transition kernel $\mathsf{P}_{V}$ of $V$. The robust
Bellman equation can be equivalently rewritten as

$\displaystyle ge=r_{\pi}-V+\mathsf{P}_{V}^{\pi}V.$ (30)

This means that $(g,V)$ is a solution to the non-robust Bellman equation for
transition kernel $\mathsf{P}_{V}$ and policy $\pi$:

$\displaystyle xe=r_{\pi}-Y+\mathsf{P}_{V}^{\pi}Y.$ (31)

Thus, by Thm 8.2.6 from (Puterman, 1994),

$\displaystyle g$ $\displaystyle=g^{\pi}_{\mathsf{P}_{V}},$ (32)
$\displaystyle V$ $\displaystyle=V^{\pi}_{\mathsf{P}_{V}}+ce,\text{ for some
}c\in\mathbb{R}.$ (33)

However, note that

$\displaystyle g^{\pi}_{\mathsf{P}_{V}}=g\leq
g^{\pi}_{\mathcal{P}}=\min_{\mathsf{P}\in\mathcal{P}}g^{\pi}_{\mathsf{P}}\leq
g^{\pi}_{\mathsf{P}_{V}},$ (34)

thus,

$\displaystyle g^{\pi}_{\mathsf{P}_{V}}=g=g^{\pi}_{\mathcal{P}}.$ (35)

2). From (35),

$\displaystyle g^{\pi}_{\mathsf{P}_{V}}=g^{\pi}_{\mathcal{P}}\>.$ (36)

It then follows from the definition of $\Omega^{\pi}_{g}$ that
$\mathsf{P}_{V}\in\Omega^{\pi}_{g}$.

3). Since $(g,V)$ is a solution to the non-robust Bellman equation

$\displaystyle xe=r_{\pi}-Y+\mathsf{P}_{V}^{\pi}Y,$ (37)

the claim then follows from Theorem 8.2.6 in (Puterman, 1994). ∎

###  B.2 Proof of Theorem 3.2

######  Theorem B.2.

(Restatement of Theorem 3.2) Under Assumptions 1,2,3, $(f(V_{n}),V_{n})$
converges to a (possible sample path dependent) solution to (7) a.s..

We first show the stability of the robust RVI TD algorithm in the following
lemma.

######  Lemma B.1.

Algorithm 1 remains bounded during the update, i.e.,

$\displaystyle\sup_{n}\|V_{n}\|<\infty,a.s..$ (38)

###### Proof.

Denote by

$\displaystyle h(V)\triangleq r_{\pi}+\sigma_{\mathcal{P}}(V)-f(V)e-V.$ (39)

Then the update of robust RVI TD can be rewritten as

$\displaystyle V_{n+1}=V_{n}+\alpha_{n}(h(V_{n})+M_{n+1}),$ (40)

where $M_{n+1}\triangleq\hat{\mathbf{T}}V_{n}-r_{\pi}-\sigma_{\mathcal{P}}(V)$
is the noise term.

Further, define the limit function $h_{\infty}$:

$\displaystyle h_{\infty}(V)\triangleq\lim_{c\to\infty}\frac{h(cV)}{c}.$ (41)

Then, from $\sigma_{\mathcal{P}^{a}_{s}}(cV)=c\sigma_{\mathcal{P}^{a}_{s}}(V)$
and $f(cV)=cf(V)$, it follows that

$\displaystyle
h_{\infty}(V)=\lim_{c\to\infty}\frac{r_{\pi}}{c}+\sigma_{\mathcal{P}}(V)-f(V)e-V=\sigma_{\mathcal{P}}(V)-f(V)e-V.$
(42)

According to Section 2.1 and Section 3.2 of (Borkar, 2009), it suffices to
verify the following assumptions:

(1). $h$ is Lipschitz;

(2). Stepsize $\alpha_{n}$ satisfies Assumption 3;

(3). Denoting by $\mathcal{F}_{n}$ the $\sigma$-algebra generated by
$V_{0},M_{1},...,M_{n}$, then $\mathbb{E}[M_{n+1}|\mathcal{F}_{n}]=0$,
$\mathbb{E}[\|M_{n+1}\|^{2}|\mathcal{F}_{n}]\leq K(1+\|V_{n}\|^{2})$ for some
constant $K>0$.

(4). $h_{\infty}$ has the origin as its unique globally asymptotically stable
equilibrium.

First, note that

$\displaystyle\|h(V_{1})-h(V_{2})\|$
$\displaystyle=\max_{s}\left|\sum_{a}\pi(a|s)(\sigma_{\mathcal{P}^{a}_{s}}(V_{1})-\sigma_{\mathcal{P}^{a}_{s}}(V_{2}))-(f(V_{1})-f(V_{2}))-(V_{1}(s)-V_{2}(s))\right|$
$\displaystyle\leq\max_{s}\left\\{\left|\sum_{a}\pi(a|s)(\sigma_{\mathcal{P}^{a}_{s}}(V_{1})-\sigma_{\mathcal{P}^{a}_{s}}(V_{2}))\right|+\left|(f(V_{1})-f(V_{2}))\right|+\left|(V_{1}(s)-V_{2}(s))\right|\right\\}$
$\displaystyle\leq(2+L_{f})\|V_{1}-V_{2}\|,$ (43)

where the last inequality follows from the fact that the support function
$\sigma_{\mathcal{P}}(\cdot)$ is $1$-Lipschitz and the assumptions on $f$ in
Assumption 2. Thus, $h$ is Lipschitz, which verifies (1).

It is straightforward that (3) is satisfied if
$\mathbb{E}[\hat{\mathbf{T}}V_{n}|\mathcal{F}_{n}]=r_{\pi}+\sigma_{\mathcal{P}}(V_{n})$
and $\text{Var}[\hat{\mathbf{T}}V_{n}|\mathcal{F}_{n}]\leq
K(1+\|V_{n}\|^{2})$. As discussed in Section 3, we assume the existence of an
unbiased estimator $\hat{\mathbf{T}}$ with bounded variance here, and we will
construct the estimator in Section 5.

Then, it suffices to verify condition (4), i.e., to show that the ODE

$\displaystyle\dot{x}(t)=h_{\infty}(x(t))$ (44)

has $0$ as its unique globally asymptotically stable equilibrium.

Define an operator
${\mathbf{T}}_{0}(V)(s)\triangleq\sum_{a}\pi(a|s)\sigma_{\mathcal{P}^{a}_{s}}(V)$.
Then, any equilibrium $W$ of (44) satisfies

$\displaystyle{\mathbf{T}}_{0}(W)-f(W)e-W=0.$ (45)

This equation can be further rewritten as a set of equations:

$\left\\{\begin{aligned} W=&{\mathbf{T}}_{0}(W)-ge,\\\
g=&f(W).\end{aligned}\right.$ (46)

The equation in (46) is the robust Bellman equation for a zero-reward robust
MDP. Hence, from Theorem 3.1, any solution $(g,W)$ to (46) satisfies

$\displaystyle g=g^{\pi}_{\mathcal{P}},W=V^{\pi}_{\mathsf{P}}+ce,$ (47)

where $V^{\pi}_{\mathsf{P}}$ is the relative value function w.r.t. some worst-
case transition kernel $\mathsf{P}$ (i.e.,
$g^{\pi}_{\mathsf{P}}=\min_{\mathsf{P}\in\mathcal{P}}g^{\pi}_{\mathsf{P}}$),
and some $c\in\mathbb{R}$.

Hence, any equilibrium of (44) satisfies

$\displaystyle W=V^{\pi}_{\mathsf{P}}+ce,f(W)=g^{\pi}_{\mathcal{P}}.$ (48)

However, note that this robust Bellman equation is for a zero-reward robust
MDP, hence for any $\mathsf{P}$,

$\displaystyle
g^{\pi}_{\mathsf{P}}=\lim_{T\to\infty}\mathbb{E}_{\mathsf{P}}\left[\sum^{T-1}_{t=0}\frac{r_{t}}{T}\right]=0,$
(49) $\displaystyle
V^{\pi}_{\mathsf{P}}=\mathbb{E}_{\mathsf{P}}\left[\sum^{\infty}_{t=0}(r_{t}-g^{\pi}_{\mathsf{P}})\right]=0,$
(50)

thus $g^{\pi}_{\mathcal{P}}=0$ and $W=ce$ for some $c\in\mathbb{R}$. From
(48), it follows that $f(W)=f(ce)=0$, for any equilibrium $W$. From Assumption
2, we have that $f(ce)=cf(e)=c=0$. This further implies that

$\displaystyle W=V^{\pi}_{\mathcal{P}}+ce=0.$ (51)

Thus, the only equilibrium of (44) is $0$.

We then show that $0$ is globally asymptotically stable. Recall that the zero-
reward robust Bellman operator

$\displaystyle{\mathbf{T}}_{0}V(s)=\sum_{a}\pi(a|s)(\sigma_{\mathcal{P}^{a}_{s}}(V)).$
(52)

We further introduce two operators:

$\displaystyle{\mathbf{T}}^{\prime}_{0}V$
$\displaystyle\triangleq{\mathbf{T}}_{0}V-f(V)e,$ (53)
$\displaystyle\tilde{\mathbf{T}}_{0}V$
$\displaystyle\triangleq{\mathbf{T}}_{0}V-g^{\pi}_{\mathcal{P}}e.$ (54)

Note that in the zero-reward robust MDP, $g^{\pi}_{\mathcal{P}}=0$ and
$\tilde{\mathbf{T}}_{0}={\mathbf{T}}_{0}$, but we introduce this notation for
future use.

Consider the ODEs w.r.t. these two operators:

$\displaystyle\dot{x}$ $\displaystyle={\mathbf{T}}_{0}^{\prime}x-x,$ (55)
$\displaystyle\dot{y}$ $\displaystyle=\tilde{\mathbf{T}}_{0}y-y.$ (56)

First, it can be easily shown that both ${\mathbf{T}}^{\prime}_{0}$ and
$\tilde{\mathbf{T}}_{0}$ are Lipschitz with constants $1+L_{f}$ and $1$,
respectively. Hence, both two ODEs are well-posed. Also, it can be seen that
(55) is the same as the ODE in (44).

Since the second equation (56) is a non-expansion (Lipschitz with parameter no
larger than $1$), Theorem 3.1 of (Borkar & Soumyanatha, 1997) implies that any
solution $y(t)$ to (56) converges to the set of equilibrium points, i.e.,

$\displaystyle y(t)\to\left\\{W:W=\tilde{\mathbf{T}}_{0}W\right\\},a.s..$ (57)

Similar to the discussion for ${\mathbf{T}}_{0}$, our Theorem 3.1 implies that
the set of equilibrium points of (56) is $\\{W=ce:c\in\mathbb{R}\\}$. Hence,
for any solution $y(t)$ to (56), $y(t)\to ce$ for some constant $k$ that may
depend on the initial value of $y(t)$.

Now, consider the solution $x(t)$ to (55). According to Lemma F.1 (note that
${\mathbf{T}}_{0}$ here is a special case of ${\mathbf{T}}$ in Lemma F.1 with
$r=0$), if the solutions $x(t),y(t)$ have the same initial value $x(0)=y(0)$,
then

$\displaystyle x(t)=y(t)+r(t)e,$ (58)

where $r(t)$ is a solution to
$\dot{r}(t)=-r(t)+g^{\pi}_{\mathcal{P}}-f(y(t)),r(0)=0$.

Note that the solution $r(t)$ with $r(0)=0$ can be written as

$\displaystyle r(t)=\int^{t}_{0}e^{-(t-s)}(g^{\pi}_{\mathcal{P}}-f(y(s)))ds$
(59)

by variation of constants formula (Abounadi et al., 2001). If we denote the
limit of $y(t)$ by $y^{*}=ce$, then
$\lim_{t\to\infty}r(t)=g^{\pi}_{\mathcal{P}}-f(y^{*})$ (Lemma B.4 in (Wan et
al., 2021), Theorem 3.4 in (Abounadi et al., 2001)). Hence, $x(t)=y(t)+r(t)e$
converges to $y^{*}+(g^{\pi}_{\mathcal{P}}-f(y^{*}))e$, i.e.,

$\displaystyle x(t)\to ce-f(ce)e=0.$ (60)

Hence, any solution $x(t)$ to (55) converges to $0$, which is its unique
equilibrium. This thus implies that 0 is the unique globally asymptotically
stable equilibrium. Together with Theorem 3.7 in (Borkar, 2009), it further
implies the boundedness of $V_{n}$, which completes the proof. ∎

We can readily prove Theorem B.2.

###### Proof.

In Lemma B.1, we have shown that

$\displaystyle\sup_{n}\|V_{n}\|<\infty,a.s..$ (61)

Thus, we have verified that conditions (A1-A3) and (A5) in (Borkar, 2009) are
satisfied. Lemma 2.1 in (Borkar, 2009) thus implies that it suffices to study
the solution to the ODE $\dot{x}(t)=h(x(t))$.

For the robust Bellman operator $\mathbf{T}V=r_{\pi}+\sigma_{\mathcal{P}}(V)$,
define

$\displaystyle{\mathbf{T}}^{\prime}V$
$\displaystyle\triangleq{\mathbf{T}}V-f(V)e,$ (62)
$\displaystyle\tilde{\mathbf{T}}V$
$\displaystyle\triangleq{\mathbf{T}}V-g^{\pi}_{\mathcal{P}}e.$ (63)

From Lemma F.1, we know that if $x(t),y(t)$ are the solutions to equations

$\displaystyle\dot{x}$ $\displaystyle={\mathbf{T}}^{\prime}x-x,$ (64)
$\displaystyle\dot{y}$ $\displaystyle=\tilde{\mathbf{T}}y-y,$ (65)

with the same initial value $x(0)=y(0)$, then

$\displaystyle x(t)=y(t)+r(t)e,$ (66)

where $r(t)$ satisfies

$\displaystyle\dot{r}(t)=-r(t)+g^{\pi}_{\mathcal{P}}-f(y(t)),r(0)=0.$ (67)

Thus, by the variation of constants formula,

$\displaystyle r(t)=\int^{t}_{0}e^{-(t-s)}(g^{\pi}_{\mathcal{P}}-f(y(s)))ds.$
(68)

Note that $\tilde{\mathbf{T}}$ is also non-expansive, hence $y(t)$ converges
to some equilibrium of (65) (Theorem 3.1 of (Borkar & Soumyanatha, 1997)). The
set of equilibrium points of (65) can be characterized as

$\displaystyle\left\\{W:\tilde{\mathbf{T}}W=W\right\\}=\left\\{W:W=TW-g^{\pi}_{\mathcal{P}}e\right\\}=\left\\{W:W(s)=\sum_{a}\pi(a|s)(r(s,a)-g^{\pi}_{\mathcal{P}}+\sigma_{\mathcal{P}^{a}_{s}}(W)),\forall
s\in\mathcal{S}\right\\}.$ (69)

From Theorem 3.1, any equilibrium of (65) can be rewritten as

$\displaystyle W=V^{\pi}_{\mathsf{P}}+ce,\text{ for some
}\mathsf{P}\in\Omega^{\pi}_{g},c\in\mathbb{R}.$ (70)

Thus, $y(t)$ converges to an equilibrium denoted by $y^{*}$:

$\displaystyle y(t)\to y^{*}\triangleq V^{\pi}_{\mathsf{P}^{*}}+c^{*}e,\text{
for some }\mathsf{P}^{*}\in\Omega^{\pi}_{g},c^{*}\in\mathbb{R}.$ (71)

Similar to Lemma B.1, it can be showed that $r(t)\to
g^{\pi}_{\mathcal{P}}-f(y^{*})$ (Lemma B.4 in (Wan et al., 2021), Theorem 3.4
in (Abounadi et al., 2001)). This further implies that

$\displaystyle x(t)\to
y^{*}+(g^{\pi}_{\mathcal{P}}-f(y^{*}))e=V^{\pi}_{\mathsf{P}^{*}}+(c^{*}+g^{\pi}_{\mathcal{P}}-f(y^{*}))e,$
(72)

and we denote $m^{*}=c^{*}+g^{\pi}_{\mathcal{P}}-f(y^{*})$. Moreover, since
$f$ is continuous (because it is Lipschitz), we have that

$\displaystyle f(x(t))$ $\displaystyle\to
f(V^{\pi}_{\mathsf{P}^{*}}+(c^{*}+g^{\pi}_{\mathcal{P}}-f(y^{*}))e)$
$\displaystyle=f(V^{\pi}_{\mathsf{P}^{*}})+c^{*}+g^{\pi}_{\mathcal{P}}-f(y^{*})$
$\displaystyle=f(V^{\pi}_{\mathsf{P}^{*}})+c^{*}+g^{\pi}_{\mathcal{P}}-f(V^{\pi}_{\mathsf{P}^{*}}+c^{*}e)$
$\displaystyle=f(V^{\pi}_{\mathsf{P}^{*}})+c^{*}+g^{\pi}_{\mathcal{P}}-f(V^{\pi}_{\mathsf{P}^{*}})-c^{*}$
$\displaystyle=g^{\pi}_{\mathcal{P}}.$ (73)

Hence, we show that

$\displaystyle x(t)$ $\displaystyle\to V^{\pi}_{\mathsf{P}^{*}}+m^{*}e,$ (74)
$\displaystyle f(x(t))$ $\displaystyle\to g^{\pi}_{\mathcal{P}}.$ (75)

Following Lemma 2.1 from (Borkar, 2009), we conclude that a.s.,

$\displaystyle V_{n}$ $\displaystyle\to V^{\pi}_{\mathsf{P}^{*}}+m^{*}e,$ (76)
$\displaystyle f(V_{n})$ $\displaystyle\to g^{\pi}_{\mathcal{P}},$ (77)

which completes the proof. ∎

##  Appendix C Robust RVI Q-Learning

###  C.1 Proof of Lemma 4.1

Part of the following theorem is proved in (Wang et al., 2023), but we include
the proof for completeness.

######  Theorem C.1 (Restatement of Lemma 4.1).

If $(g,Q)$ is a solution to the optimal robust Bellman equation

$\displaystyle Q(s,a)=r(s,a)-g+\sigma_{\mathcal{P}^{a}_{s}}(V_{Q}),\forall
s,a,$ (78)

then 1) $g=g^{*}_{\mathcal{P}}$ (Wang et al., 2023); 2) the greedy policy
w.r.t. $Q$: $\pi_{Q}(s)=\arg\max_{a}Q(s,a)$ is an optimal robust policy (Wang
et al., 2023); 3) $V_{Q}=V^{\pi_{Q}}_{\mathsf{P}}+ce$ for some
$\mathsf{P}\in\Omega^{\pi_{Q}}_{g},c\in\mathbb{R}$.

###### Proof.

Taking the maximum on both sides of (78) w.r.t. $a$, we have that

$\displaystyle\max_{a}Q(s,a)=\max_{a}\\{r(s,a)-g+\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})\\},\forall
s\in\mathcal{S}.$ (79)

This is equivalent to

$\displaystyle
V_{Q}(s)=\max_{a}\\{r(s,a)-g+\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})\\},\forall
s\in\mathcal{S}.$ (80)

By Theorem 7 in (Wang et al., 2023), we can show that $g=g^{*}_{\mathcal{P}}$,
which proves claim (1).

Recall that $V_{Q}(s)=\max_{a}Q(s,a)$. It can be also written as

$\displaystyle V_{Q}(s)=\sum_{a}\pi_{Q}(a|s)Q(s,a).$ (81)

Here, we slightly abuse the notation of $\pi_{Q}$, and use $\pi_{Q}(s)$ and
$\pi_{Q}(a|s)$ interchangeably.

Then, the optimal robust Bellman equation in (79) can be rewritten as

$\displaystyle
Q(s,\pi_{Q}(s))=r(s,\pi_{Q}(s))-g+\sigma_{\mathcal{P}^{\pi_{Q}(s)}_{s}}\bigg{(}\sum_{a}\pi_{Q}(a|\cdot)Q(\cdot,a)\bigg{)}.$
(82)

Moreover, if we denote by $W(s)=Q(s,a)=Q(s,\pi_{Q}(s))=\max_{a}Q(s,a)$, then
the equation above is equivalent to

$\displaystyle
W(s)=\sum_{a}\pi_{Q}(a|s)(r(s,a)-g+\sigma_{\mathcal{P}^{a}_{s}}(W)).$ (83)

Therefore, $(W,g)$ is a solution to the robust Bellman equation for the policy
$\pi_{Q}$ in Theorem 3.1. By Theorem 3.1, we have that

$\displaystyle g$ $\displaystyle=g^{\pi_{Q}}_{\mathcal{P}},$ (84)
$\displaystyle W$ $\displaystyle=V^{\pi_{Q}}_{\mathsf{P}}+ce,$ (85)

for some $\mathsf{P}\in\Omega_{g}^{\pi_{Q}}$ and $c\in\mathbb{R}$.

Combining this with the claim (1) implies that $\pi_{Q}$ is an optimal robust
policy. Claims (2) and (3) are thus proved. ∎

###  C.2 Proof of Theorem 4.1

######  Lemma C.1.

If $\hat{\mathbf{H}}$ satisfies that for any
$Q,s\in\mathcal{S},a\in\mathcal{A}$,
$\mathbb{E}[\hat{\mathbf{H}}Q(s,a)]=\mathbf{H}Q(s,a)$ and
$\text{Var}(\hat{\mathbf{H}}Q(s,a))\leq C(1+\|Q\|^{2})$ for some constant $C$,
then under Assumptions 2, 4 and 3, Algorithm 2 remains bounded during the
update almost surely, i.e.,

$\displaystyle\sup_{n}\|Q_{n}\|<\infty,a.s..$ (86)

###### Proof.

Denote by

$\displaystyle h(Q)\triangleq r_{\pi}+\sigma_{\mathcal{P}}(V_{Q})-f(Q)e-Q.$
(87)

Then, the update of robust RVI Q-learning can be rewritten as

$\displaystyle Q_{n+1}=Q_{n}+\alpha_{n}(h(Q_{n})+M_{n+1}),$ (88)

where
$M_{n+1}\triangleq\hat{\mathbf{H}}Q_{n}-r_{\pi}-\sigma_{\mathcal{P}}(V_{Q})$
is the noise term.

Further, define the limit function $h_{\infty}$:

$\displaystyle h_{\infty}(Q)\triangleq\lim_{c\to\infty}\frac{h(cQ)}{c}.$ (89)

Then, note that
$\sigma_{\mathcal{P}^{a}_{s}}(V_{cQ})=\sigma_{\mathcal{P}^{a}_{s}}(cV_{Q})=c\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})$
for $c>0$ and $f(cQ)=cf(Q)$. It then follows that

$\displaystyle
h_{\infty}(Q)=\lim_{c\to\infty}\frac{r_{\pi}}{c}+\sigma_{\mathcal{P}}(V_{Q})-f(Q)e-Q=\sigma_{\mathcal{P}}(V_{Q})-f(Q)e-Q.$
(90)

Similar to the proof of Lemma B.1, it suffices to verify the following
conditions:

(1). $h$ is Lipschitz;

(2). Stepsize $\alpha_{n}$ satisfies Assumption 3;

(3). $\mathbb{E}[M_{n+1}|\mathcal{F}_{n}]=0$, and
$\mathbb{E}[\|M_{n+1}\|^{2}|\mathcal{F}_{n}]\leq K(1+\|Q_{n}\|^{2})$ for some
constant $K$.

(4). $h_{\infty}$ has the origin as its unique globally asymptotically stable
equilibrium.

Clearly, (2) and (3) can be verified similarly to Lemma B.1. We then verify
(1) and (4).

Firstly, it can be shown that

$\displaystyle|h(Q_{1})(s,a)-h(Q_{2})(s,a)|$
$\displaystyle=\big{|}\sigma_{\mathcal{P}^{a}_{s}}(V_{Q_{1}})-f(Q_{1})-Q_{1}(s,a)-\sigma_{\mathcal{P}^{a}_{s}}(V_{Q_{2}})-f(Q_{2})-Q_{2}(s,a)\big{|}$
$\displaystyle\leq\big{|}\sigma_{\mathcal{P}^{a}_{s}}(V_{Q_{1}})-\sigma_{\mathcal{P}^{a}_{s}}(V_{Q_{2}})\big{|}+|f(Q_{1})-f(Q_{2})|+|Q_{1}(s,a)-Q_{2}(s,a)|$
$\displaystyle\leq\|V_{Q_{1}}-V_{Q_{2}}\|+L_{f}\|Q_{1}-Q_{2}\|+\|Q_{1}-Q_{2}\|$
$\displaystyle\leq(2+L_{f})\|Q_{1}-Q_{2}\|,$ (91)

where the last inequality is from the fact that
$\|V_{Q_{1}}-V_{Q_{2}}\|\leq\|Q_{1}-Q_{2}\|$. This implies that $h$ is
Lipschitz.

To verify (4), note that the stability equation is

$\displaystyle\dot{X}(t)=h_{\infty}(X(t))=\sigma_{\mathcal{P}}(V_{X}(t))-f(X(t))e-X(t),$
(92)

where $V_{X}(t)$ is a $|\mathcal{S}|$-dimensional vector with
$V_{X}(t)(s)=\max_{a}X(t)(s,a).$

Any equilibrium $Q$ of the stability equation (92) satisfies that

$\displaystyle Q(s,a)=\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})-f(Q)e,$ (93)

which can be viewed as an optimal robust Bellman equation (10) with zero
reward. Hence, by Lemma 4.1, it implies that

$\displaystyle f(Q)$ $\displaystyle=g^{*}_{\mathcal{P}}=0,$ (94)
$\displaystyle V_{Q}$ $\displaystyle=V^{\pi_{Q}}_{\mathsf{P}}+ce\text{ for
some }\mathsf{P}\in\Omega^{\pi_{Q}}_{g},c\in\mathbb{R}.$ (95)

In the zero-reward MDP, we have that $V^{\pi}_{\mathsf{P}}=0$ for any
$\pi,\mathsf{P}$, thus $V_{Q}(s)=\max_{a}Q(s,a)=c$ for any $s\in\mathcal{S}$.

Note that from (93), $Q$ satisfies that

$\displaystyle
Q(s,a)=\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})=\sigma_{\mathcal{P}^{a}_{s}}(ce)=c.$
(96)

Since $f(Q)=0$, it implies that

$\displaystyle f(Q)=f(ce)=c=0.$ (97)

Therefore,

$\displaystyle c=0,$ (98) $\displaystyle Q=0.$ (99)

Thus, $0$ is the unique equilibrium of the stability equation.

We then show that $0$ is globally asymptotically stable. Define the zero-
reward optimal robust Bellman operator

$\displaystyle{\mathbf{H}}_{0}Q(s,a)=\sigma_{\mathcal{P}^{a}_{s}}(V_{Q}),$
(100)

and further introduce two operators

$\displaystyle{\mathbf{H}}^{\prime}_{0}Q(s,a)$
$\displaystyle=\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})-f(Q),$ (101)
$\displaystyle\tilde{\mathbf{H}}_{0}Q(s,a)$
$\displaystyle=\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})-g^{*}_{\mathcal{P}}.$ (102)

It is straightforward to verify that $\tilde{\mathbf{H}}_{0}$ is non-
expansive. Hence by (Borkar & Soumyanatha, 1997), the solution $y(t)$ to
equation

$\displaystyle\dot{y}=\tilde{\mathbf{H}}_{0}y-y$ (103)

converges to the set of equilibrium points

$\displaystyle\\{W:W(s,a)=\sigma_{\mathcal{P}^{a}_{s}}(V_{W})-g^{*}_{\mathcal{P}}\\},a.s..$
(104)

This again can be viewed as an optimal robust Bellman equation with zero-
reward. Hence, any equilibrium $W$ of (103) satisfies

$\displaystyle\max_{a}W(s,a)=c,\forall s.$ (105)

This together with (104) further implies that the equilibrium $W$ of (103)
satisfies

$\displaystyle
W(s,a)=\sigma_{\mathcal{P}^{a}_{s}}(V_{W})=\sigma_{\mathcal{P}^{a}_{s}}(ce)=c,$
(106)

and hence $y(t)$ converges to $\\{ce:c\in\mathbb{R}\\}$. We denote its limit
by $y^{*}=c^{*}e$.

Lemma F.6 implies the solution $x(t)$ to the ODE
$\dot{x}={\mathbf{H}}_{0}^{\prime}(x)-x$ can be decomposed as
$x(t)=y(t)+r(t)e$, where $r(t)$ satisfies
$\dot{r}(t)=-r(t)+g^{*}_{\mathcal{P}}-f(y(t)),r(0)=0$.

Then, similar to Lemma B.1, Lemma B.4 in (Wan et al., 2021) and Theorem 3.4 in
(Abounadi et al., 2001), it can be shown that $r(t)\to
g^{*}_{\mathcal{P}}-f(y(t))=-c^{*}$. Hence,

$\displaystyle x(t)\to 0,$ (107)

which proves the asymptotic stability.

Thus, we conclude that $0$ is the unique globally asymptotically stable
equilibrium of the stability equation, which implies the boundedness of
$\\{Q_{n}\\}$ together with results from Section 2.1 and 3.2 from (Borkar,
2009). ∎

######  Theorem C.2 (Restatement of Theorem 4.1).

The sequence $\\{Q_{n}\\}$ generated by Algorithm 2 converges to a solution
$Q^{*}$ to the optimal robust Bellman equation a.s., and $f(Q_{n})$ converges
to the optimal robust average-reward $g^{*}_{\mathcal{P}}$ a.s..

###### Proof.

According to Lemma 1 from (Borkar, 2009) and Theorem 3.5 from (Abounadi et
al., 2001), the sequence $\\{Q_{n}\\}$ converge to the same limit as the
solution $x(t)$ to the ODE $\dot{x}={\mathbf{H}}^{\prime}x-x$. Hence the proof
can be completed by showing convergence of $x(t)$ and $f(x(t))$.

For the optimal robust Bellman operator,

$\displaystyle{\mathbf{H}}Q(s,a)=r(s,a)+\sigma_{\mathcal{P}^{a}_{s}}(V_{Q}),$
(108)

define two operators

$\displaystyle{\mathbf{H}}^{\prime}Q$
$\displaystyle\triangleq{\mathbf{H}}Q-f(Q)e,$ (109)
$\displaystyle\tilde{\mathbf{H}}Q$
$\displaystyle\triangleq{\mathbf{H}}Q-g^{*}_{\mathcal{P}}e.$ (110)

From Lemma F.6, we know that if $x(t),y(t)$ are the solutions to equations

$\displaystyle\dot{x}$ $\displaystyle={\mathbf{H}}^{\prime}x-x,$ (111)
$\displaystyle\dot{y}$ $\displaystyle=\tilde{\mathbf{H}}y-y,$ (112)

with the same initial value $x(0)=y(0)$, then

$\displaystyle x(t)=y(t)+r(t)e,$ (113)

where $r(t)$ satisfies

$\displaystyle\dot{r}(t)=-r(t)+g^{*}_{\mathcal{P}}-f(y(t)),r(0)=0.$ (114)

It can be easily verified that $\tilde{\mathbf{H}}$ is non-expansive. Hence
$y(t)$ converges to the set of equilibrium points of of (112) (Theorem 3.1 of
(Borkar & Soumyanatha, 1997)), which can be characterized as

$\displaystyle\left\\{W:\tilde{\mathbf{H}}W=W\right\\}=\left\\{W:W={\mathbf{H}}W-g^{*}_{\mathcal{P}}e\right\\}=\left\\{W:W(s,a)=r(s,a)-g^{*}_{\mathcal{P}}+\sigma_{\mathcal{P}^{a}_{s}}(V_{W}),\forall
s,a\right\\}.$ (115)

From Lemma 4.1, any equilibrium $W$ satisfies

$\displaystyle V_{W}=V^{\pi_{W}}_{\mathsf{P}}+ce,\text{ for some
}\mathsf{P}\in\Omega^{\pi_{W}}_{g},c\in\mathbb{R},$ (116)

and $\pi_{W}$ is robust optimal. We denote the limit of $y(t)$ by $W^{*}$.

Similar to (72) to (76), it can be shown that $r(t)\to
g^{*}_{\mathcal{P}}-f(W^{*})$. This further implies that

$\displaystyle x(t)\to W^{*}+(g^{*}_{\mathcal{P}}-f(W^{*}))e\triangleq
W^{*}+m^{*}e,$ (117)

where $m^{*}=g^{*}_{\mathcal{P}}-f(W^{*})$. Note that $W^{*}+m^{*}e$ is a
solution to the optimal robust Bellman equation, hence $x(t)$ converges to a
solution to (10). Moreover, since $f$ is continuous (because it is Lipschitz),
we have that

$\displaystyle f(x(t))$ $\displaystyle\to f(W^{*}+m^{*}e)$
$\displaystyle=f(W^{*})+g^{*}_{\mathcal{P}}-f(W^{*})$
$\displaystyle=g^{*}_{\mathcal{P}}.$ (118)

This completes the proof. ∎

##  Appendix D Case Studies for Robust RVI TD

In this section, we provide the proof of the first part of Theorem 5.1, i.e.,
that $\hat{\mathbf{T}}$ is unbiased and has bounded variance under each
uncertainty model.

We first show a lemma, by which the problem can be reduced to investigating
whether $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ is unbiased and has bounded
variance.

######  Lemma D.1.

If

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]=\sigma_{\mathcal{P}^{a}_{s}}(V),\forall
s,a,$ (119)

and moreover, there exists a constant $C$, such that

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\leq
C(1+\|V\|^{2}),\forall s,a,$ (120)

then

$\displaystyle\mathbb{E}[\hat{\mathbf{T}}V(s)]=\mathbf{T}V(s),\forall s,$
(121)

and

$\displaystyle\text{Var}(\hat{\mathbf{T}}V(s))\leq|\mathcal{A}|C(1+\|V\|^{2}),\forall
s.$ (122)

###### Proof.

From the definition,
$\hat{\mathbf{T}}V(s)=\sum_{a}\pi(a|s)(r(s,a)+\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)$.
Thus,

$\displaystyle\mathbb{E}[\hat{\mathbf{T}}V(s)]$
$\displaystyle=\mathbb{E}\bigg{[}\sum_{a}\pi(a|s)(r(s,a)+\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\bigg{]}$
$\displaystyle=\sum_{a}\pi(a|s)(r(s,a)+\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V])$
$\displaystyle=\sum_{a}\pi(a|s)(r(s,a)+\sigma_{\mathcal{P}^{a}_{s}}(V))=\mathbf{T}V(s),$
(123)

which shows that $\hat{\mathbf{T}}$ is unbiased. On the other hand, we have
that

$\displaystyle\text{Var}(\hat{\mathbf{T}}V(s))$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}\sum_{a}\pi(a|s)(r(s,a)+\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)-\mathbb{E}\bigg{[}\sum_{a}\pi(a|s)(r(s,a)+\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\bigg{]}\bigg{)}^{2}\bigg{]}$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}\sum_{a}\pi(a|s)(r(s,a)+\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)-\sum_{a}\pi(a|s)(r(s,a)+\mathbb{E}\big{[}\hat{\sigma}_{\mathcal{P}^{a}_{s}}V\big{]}\bigg{)}^{2}\bigg{]}$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}\sum_{a}\pi(a|s)(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)-\mathbb{E}\big{[}\hat{\sigma}_{\mathcal{P}^{a}_{s}}V\big{]}\bigg{)}^{2}\bigg{]}$
$\displaystyle\overset{(a)}{\leq}\mathbb{E}\bigg{[}\sum_{a}\pi(a|s)(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V-\mathbb{E}\big{[}\hat{\sigma}_{\mathcal{P}^{a}_{s}}V\big{]})^{2}\bigg{]}$
$\displaystyle=\sum_{a}\pi(a|s)\mathbb{E}\bigg{[}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V-\mathbb{E}\big{[}\hat{\sigma}_{\mathcal{P}^{a}_{s}}V^{2}\big{]})^{2}\bigg{]}$
$\displaystyle{\leq}\sum_{a}\pi(a|s)\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)$
$\displaystyle\leq|\mathcal{A}|C(1+\|V\|^{2}),$ (124)

where $(a)$ is because $(\mathbb{E}[X])^{2}\leq\mathbb{E}[X^{2}]$, which
completes the proof. ∎

This lemma implies that to prove Theorem 5.1, it suffices to show that
$\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ is unbiased and has bounded variance.

###  D.1 Contamination Uncertainty Set

######  Theorem D.1.

$\hat{\mathbf{T}}$ defined in (12) is unbiased and has bounded variance.

###### Proof.

First, note that

$\displaystyle V_{n+1}(s)$
$\displaystyle=V_{n}(s)+\alpha_{n}(r(s,a)+((1-\delta)V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)-f(V_{n})-V_{n}(s))$
$\displaystyle=V_{n}(s)+\alpha_{n}({\mathbf{T}}V_{n}(s)-f(V_{n})-V_{n}(s)+M_{n}(s)),$
(125)

where

$\displaystyle
M_{n}(s)=r(s,a)+(1-\delta)V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)-{\mathbf{T}}V_{n}(s),$
(126)

and

$\displaystyle{\mathbf{T}}V_{n}(s)=\sum_{a}\pi(a|s)\bigg{(}r(s,a)+(1-\delta)\sum_{s^{\prime}}\mathsf{P}^{a}_{s,s^{\prime}}V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)\bigg{)}.$
(127)

Thus,

$\displaystyle\mathbb{E}[M_{n}(s)]$
$\displaystyle=\mathbb{E}\big{[}r(s,a)+(1-\delta)V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)\big{]}-\sum_{a}\pi(a|s)\bigg{(}r(s,a)+(1-\delta)\sum_{s^{\prime}}\mathsf{P}^{a}_{s,s^{\prime}}V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)\bigg{)}$
$\displaystyle=\sum_{a}\pi(a|s)\bigg{(}r(s,a)+(1-\delta)\sum_{s^{\prime}}\mathsf{P}^{a}_{s,s^{\prime}}V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)\bigg{)}$
$\displaystyle\quad\quad\quad\quad\quad\quad\quad\quad\quad-\sum_{a}\pi(a|s)\bigg{(}r(s,a)+(1-\delta)\sum_{s^{\prime}}\mathsf{P}^{a}_{s,s^{\prime}}V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x))\bigg{)}$
$\displaystyle=0.$ (128)

Hence, the operator is unbiased.

We also have that

$\displaystyle\mathbb{E}[|M_{n}(s)|^{2}]$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}r(s,a)+(1-\delta)V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)-{\mathbf{T}}V_{n}(s)\bigg{)}^{2}\bigg{]}$
$\displaystyle\leq
2\mathbb{E}\bigg{[}\bigg{(}r(s,a)+(1-\delta)V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)\bigg{)}^{2}\bigg{]}+2\mathbb{E}[({\mathbf{T}}V_{n}(s))^{2}]$
$\displaystyle\overset{(a)}{\leq}8+8\|V_{n}\|^{2}$ $\displaystyle\leq
8(1+\|V_{n}\|^{2}),$ (129)

where $(a)$ is from the fact that
$\mathbb{E}\big{[}\big{(}(1-\delta)V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)\big{)}^{2}\big{]}=\mathbb{E}\big{[}\big{|}(1-\delta)V_{n}(s^{\prime})+\delta\min_{x}V_{n}(x)\big{|}^{2}\big{]}\leq\mathbb{E}\big{[}\big{(}\big{|}(1-\delta)V_{n}(s^{\prime})\big{|}+\big{|}\delta\min_{x}V_{n}(x)\big{|}\big{)}^{2}\big{]}\leq\mathbb{E}\big{[}\big{(}(1-\delta)\|V_{n}\|+\big{(}\delta\|V_{n}\|\big{)}^{2}\big{]}\leq\|V_{n}\|^{2}$.

The proof is completed. ∎

###  D.2 Total Variation Uncertainty Set

The estimator under the total variation uncertainty set can be written as

$\displaystyle\hat{\sigma}_{\mathcal{P}^{a}_{s}}(V)=\max_{\mu\geq
0}\big{(}\hat{\mathsf{P}}^{a,1}_{s,N+1}(V-\mu)-\delta\text{Span}(V-\mu)\big{)}+\frac{\Delta_{N}(V)}{p_{N}},$
(130)

where

$\displaystyle\Delta_{N}(V)=\max_{\mu\geq
0}\big{(}\hat{\mathsf{P}}^{a}_{s,N+1}(V-\mu)-\delta\text{Span}(V-\mu)\big{)}$
$\displaystyle\quad-\frac{1}{2}\max_{\mu\geq
0}\big{(}\hat{\mathsf{P}}^{a,O}_{s,N+1}(V-\mu)-\delta\text{Span}(V-\mu)\big{)}$
$\displaystyle\quad-\frac{1}{2}\max_{\mu\geq
0}\big{(}\hat{\mathsf{P}}^{a,E}_{s,N+1}(V-\mu)-\delta\text{Span}(V-\mu)\big{)}.$
(131)

######  Theorem D.2.

The estimated operator $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ defined in (130)
is unbiased, i.e.,

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]=\sigma_{\mathcal{P}^{a}_{s}}(V).$
(132)

###### Proof.

First, denote the dual function (14) by $g$:

$\displaystyle
g_{s,a}^{V}(\mu)=\mathsf{P}^{a}_{s}(V-\mu)-\delta\text{Span}(V-\mu),$ (133)

and denote its optimal solution by $\mu^{V}_{s,a}$:

$\displaystyle\mu^{V}_{s,a}=\arg\max_{\mu\geq
0}\big{(}\mathsf{P}^{a}_{s}(V-\mu)-\delta\text{Span}(V-\mu)\big{)}.$ (134)

Then, the support function
$\sigma_{\mathcal{P}^{a}_{s}}(V)=g_{s,a}^{V}(\mu^{V}_{s,a})$. Similarly,
define the empirical function

$\displaystyle{g}_{s,a,N+1}^{V}(\mu)$
$\displaystyle=\hat{\mathsf{P}}^{a}_{s,N+1}(V-\mu)-\delta\text{Span}(V-\mu),$
(135) $\displaystyle{g}_{s,a,N+1,O}^{V}(\mu)$
$\displaystyle=\hat{\mathsf{P}}^{a,O}_{s,N+1}(V-\mu)-\delta\text{Span}(V-\mu),$
(136) $\displaystyle{g}_{s,a,N+1,E}^{V}(\mu)$
$\displaystyle=\hat{\mathsf{P}}^{a,E}_{s,N+1}(V-\mu)-\delta\text{Span}(V-\mu),$
(137)

and their optimal solutions are denoted by
$\mu^{V}_{s,a,N+1},\mu^{V}_{s,a,N+1,O},\mu^{V}_{s,a,N+1,E}$. We have that

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]$
$\displaystyle=\mathbb{E}\bigg{[}\max_{\mu\geq
0}\big{(}\hat{\mathsf{P}}^{a,1}_{s,N+1}(V-\mu)-\delta\text{Span}(V-\mu)\big{)}+\frac{\Delta_{N}(V)}{p_{N}}\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\mathbb{E}\bigg{[}\frac{\Delta_{N}(V)}{p_{N}}\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}p(N=n)\mathbb{E}\bigg{[}\frac{\Delta_{N}(V)}{p_{N}}|N=n\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}[\Delta_{n}(V)]$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}\bigg{[}g^{V}_{s,a,n+1}(\mu^{V}_{s,a,n+1})-\frac{g^{V}_{s,a,n+1,O}(\mu^{V}_{s,a,n+1,O})+g^{V}_{s,a,n+1,E}(\mu^{V}_{s,a,n+1,E})}{2}\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}\bigg{[}g^{V}_{s,a,n+1}(\mu^{V}_{s,a,n+1})-g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]},$
(138)

where the last inequality is from Lemma F.2. The last equation can be further
rewritten as

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}\bigg{[}g^{V}_{s,a,n+1}(\mu^{V}_{s,a,n+1})-g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]}$
$\displaystyle=\lim_{n\to\infty}\mathbb{E}\bigg{[}g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]}.$
(139)

To show that $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ is unbiased, it suffices to
prove that

$\displaystyle\lim_{n\to\infty}\mathbb{E}\bigg{[}g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]}=g^{V}_{s,a}(\mu^{V}_{s,a}).$
(140)

For any arbitrary i.i.d. samples $\\{X_{i}\\}$ and its corresponding function
$g^{V}_{s,a,n}$, together with Lemma F.3, we have that

$\displaystyle|g^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a})|$
$\displaystyle=|\max_{0\leq\mu\leq
V+\|V\|e}g^{V}_{s,a}(\mu)-\max_{0\leq\mu\leq V+\|V\|e}g^{V}_{s,a,n}(\mu)|$
$\displaystyle\leq\max_{0\leq\mu\leq
V+\|V\|e}|g^{V}_{s,a}(\mu)-g^{V}_{s,a,n}(\mu)|$
$\displaystyle=\max_{0\leq\mu\leq
V+\|V\|e}|\mathsf{P}^{a}_{s}(V-\mu)-\delta\text{Span}(V-\mu)-\hat{\mathsf{P}}^{a}_{s,n}(V-\mu)+\delta\text{Span}(V-\mu)|$
$\displaystyle=\max_{0\leq\mu\leq
V+\|V\|e}|\mathsf{P}^{a}_{s}(V-\mu)-\hat{\mathsf{P}}^{a}_{s,n}(V-\mu)|$
$\displaystyle\leq\max_{0\leq\mu\leq
V+\|V\|e}\|V-\mu\|\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1}$
$\displaystyle\leq
3\|V\|\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1}.$ (141)

Thus, by Hoeffding’s inequality and Theorem 3.7 from (Liu et al., 2022),

$\displaystyle\mathbb{E}[|g^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a})|]\leq
3\|V\|\frac{|\mathcal{S}|^{2}\sqrt{\pi}}{2^{\frac{n+1}{2}}},$ (142)

which implies that

$\displaystyle\lim_{n\to\infty}\mathbb{E}\bigg{[}g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]}=g^{V}_{s,a}(\mu^{V}_{s,a}),$
(143)

completing the proof. ∎

######  Theorem D.3.

The estimated operator $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ defined in (130)
has bounded variance, i.e., there exists a constant $C_{0}$, such that

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\leq(1+18(1+2\delta)^{2}+2C_{0})\|V\|^{2}.$
(144)

###### Proof.

Similar to Theorem D.2, we have that

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)$
$\displaystyle=\mathbb{E}[(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)^{2}]-\sigma_{\mathcal{P}^{a}_{s}}(V)^{2}$
$\displaystyle\leq\mathbb{E}\bigg{[}\bigg{(}g^{V}_{s,a,0}(\mu^{V}_{s,a,0})+\frac{\Delta_{N}(V)}{p_{N}}\bigg{)}^{2}\bigg{]}+(\sigma_{\mathcal{P}^{a}_{s}}(V))^{2}$
$\displaystyle\leq
2\mathbb{E}\bigg{[}\bigg{(}g^{V}_{s,a,0}(\mu^{V}_{s,a,0})\bigg{)}^{2}\bigg{]}+2\mathbb{E}\bigg{[}\bigg{(}\frac{\Delta_{N}(V)}{p_{N}}\bigg{)}^{2}\bigg{]}+(\sigma_{\mathcal{P}^{a}_{s}}(V))^{2}$
$\displaystyle\leq(1+18(1+2\delta)^{2})\|V\|^{2}+2\sum^{\infty}_{i=0}\frac{\mathbb{E}[(\Delta_{i}(V))^{2}]}{p_{i}},$
(145)

where the last inequality is from Lemma F.3. For any $n\geq 1$, we have that

$\displaystyle\mathbb{E}[(\Delta_{n}(V))^{2}]$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}{g}^{V}_{s,a,n}(\mu^{V}_{s,a,n})-\frac{{g}^{V}_{s,a,n,E}(\mu^{V}_{s,a,n,E})+{g}^{V}_{s,a,n,O}(\mu^{V}_{s,a,n,O})}{2}\bigg{)}^{2}\bigg{]}$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}{g}^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a})+g^{V}_{s,a}(\mu^{V}_{s,a})-\frac{{g}^{V}_{s,a,n,E}(\mu^{V}_{s,a,n,E})+{g}^{V}_{s,a,n,O}(\mu^{V}_{s,a,n,O})}{2}\bigg{)}^{2}\bigg{]}$
$\displaystyle\leq
2\mathbb{E}[({g}^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a}))^{2}]+2\mathbb{E}\bigg{[}\bigg{(}g^{V}_{s,a}(\mu^{V}_{s,a})-\frac{{g}^{V}_{s,a,n,E}(\mu^{V}_{s,a,n,E})+{g}^{V}_{s,a,n,O}(\mu^{V}_{s,a,n,O})}{2}\bigg{)}^{2}\bigg{]}$
$\displaystyle\overset{(a)}{=}2\mathbb{E}[({g}^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a}))^{2}]+2\mathbb{E}[({g}^{V}_{s,a,n-1}(\mu^{V}_{s,a,n-1})-g^{V}_{s,a}(\mu^{V}_{s,a}))^{2}]$
$\displaystyle\leq
18\|V\|^{2}\mathbb{E}[\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1}^{2}]+18\|V\|^{2}\mathbb{E}[\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n-1}\|_{1}^{2}],$
(146)

where $(a)$ is due to Lemma F.2 and the last inequality follows a similar
argument to (D.2). Note that $p_{n}=\Psi(1-\Psi)^{n}$ for $\Psi\in(0,0.5)$,
thus similar to Theorem 3.7 of (Liu et al., 2022), we can show that there
exists a constant $C_{0}$, such that

$\displaystyle\sum^{\infty}_{i=0}\frac{\mathbb{E}[(\Delta_{i}(V))^{2}]}{p_{i}}\leq
C_{0}\|V\|^{2}.$ (147)

Thus,

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\leq(1+18(1+2\delta)^{2})\|V\|^{2}+2C_{0}\|V\|^{2}=(1+18(1+2\delta)^{2}+2C_{0})\|V\|^{2}\>.$
(148)

∎

###  D.3 Chi-Square Uncertainty Set

The estimator under the Chi-square uncertainty set can be written as

$\displaystyle\hat{\sigma}_{\mathcal{P}^{a}_{s}}V$
$\displaystyle=\max_{\mu\geq
0}\big{(}\hat{\mathsf{P}}^{a,1}_{s,N+1}(V-\mu)-\sqrt{\delta\text{Var}_{\hat{\mathsf{P}}^{a,1}_{s,N+1}}(V-\mu)}\big{)}$
$\displaystyle\quad+\frac{\Delta_{N}(V)}{p_{N}},$ (149)

where

$\displaystyle\Delta_{N}(V)=\max_{\mu\geq
0}\big{(}\mathbb{E}_{\hat{\mathsf{P}}^{a}_{s,N+1}}[V-\mu]-\sqrt{\delta\text{Var}_{\hat{\mathsf{P}}^{a}_{s,N+1}}(V-\mu)}\big{)}$
$\displaystyle\quad-\frac{1}{2}\max_{\mu\geq
0}\big{(}\mathbb{E}_{\hat{\mathsf{P}}^{a,O}_{s,N+1}}[V-\mu]-\sqrt{\delta\text{Var}_{\hat{\mathsf{P}}^{a,O}_{s,N+1}}(V-\mu)}\big{)}$
$\displaystyle\quad-\frac{1}{2}\max_{\mu\geq
0}\big{(}\mathbb{E}_{\hat{\mathsf{P}}^{a,E}_{s,N+1}}[V-\mu]-\sqrt{\delta\text{Var}_{\hat{\mathsf{P}}^{a,E}_{s,N+1}}(V-\mu)}\big{)}.$

######  Theorem D.4.

The estimated operator defined in (D.3) is unbiased, i.e.,

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]=\sigma_{\mathcal{P}^{a}_{s}}(V).$
(150)

###### Proof.

Denote the dual function (15) by $g$:

$\displaystyle
g_{s,a}^{V}(\mu)=\mathsf{P}^{a}_{s}(V-\mu)-\sqrt{\delta\text{Var}_{\mathsf{P}^{a}_{s}}(V-\mu)},$
(151)

and denote its optimal solution by $\mu^{V}_{s,a}$:

$\displaystyle\mu^{V}_{s,a}=\arg\max_{\mu\geq
0}\big{(}\mathsf{P}^{a}_{s}(V-\mu)-\sqrt{\delta\text{Var}_{\mathsf{P}^{a}_{s}}(V-\mu)}\big{)}.$
(152)

Then, the support function
$\sigma_{\mathcal{P}^{a}_{s}}(V)=g_{s,a}^{V}(\mu^{V}_{s,a})$. Similarly,
define the empirical function

$\displaystyle{g}_{s,a,N+1}^{V}(\mu)$
$\displaystyle=\hat{\mathsf{P}}^{a}_{s,N+1}(V-\mu)-\sqrt{\delta\text{Var}_{\hat{\mathsf{P}}^{a}_{s,N+1}}(V-\mu)},$
(153) $\displaystyle{g}_{s,a,N+1,O}^{V}(\mu)$
$\displaystyle=\hat{\mathsf{P}}^{a,O}_{s,N+1}(V-\mu)-\sqrt{\delta\text{Var}_{\hat{\mathsf{P}}^{a,O}_{s,N+1}}(V-\mu)},$
(154) $\displaystyle{g}_{s,a,N+1,E}^{V}(\mu)$
$\displaystyle=\hat{\mathsf{P}}^{a,E}_{s,N+1}(V-\mu)-\sqrt{\delta\text{Var}_{\hat{\mathsf{P}}^{a,E}_{s,N+1}}(V-\mu)},$
(155)

and their optimal solutions are denoted by
$\mu^{V}_{s,a,N+1},\mu^{V}_{s,a,N+1,O},\mu^{V}_{s,a,N+1,E}$. We have that

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\mathbb{E}\bigg{[}\frac{\Delta_{N}(V)}{p_{N}}\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}p(N=n)\mathbb{E}\bigg{[}\frac{\Delta_{N}(V)}{p_{N}}|N=n\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}[\Delta_{n}]$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}\bigg{[}g^{V}_{s,a,n+1}(\mu^{V}_{s,a,n+1})-\frac{g^{V}_{s,a,n+1,O}(\mu^{V}_{s,a,n+1,O})+g^{V}_{s,a,n+1,E}(\mu^{V}_{s,a,n+1,E})}{2}\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}\bigg{[}g^{V}_{s,a,n+1}(\mu^{V}_{s,a,n+1})-g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]},$
(156)

where the last inequality is from Lemma F.2. The last equation can be further
rewritten as

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\mu^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}\bigg{[}g^{V}_{s,a,n+1}(\mu^{V}_{s,a,n+1})-g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]}$
$\displaystyle=\lim_{n\to\infty}\mathbb{E}\bigg{[}g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]}.$
(157)

To show that $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ is unbiased, it suffices to
prove that

$\displaystyle\lim_{n\to\infty}\mathbb{E}\bigg{[}g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]}=g^{V}_{s,a}(\mu^{V}_{s,a}).$
(158)

For any arbitrary i.i.d. samples $\\{X_{i}\\}$ and its corresponding function
$g^{V}_{s,a,n}$, together with Lemma F.4, we have that

$\displaystyle|g^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a})|$
$\displaystyle=|\max_{0\leq\mu\leq
V+\|V\|e}g^{V}_{s,a}(\mu)-\max_{0\leq\mu\leq V+\|V\|e}g^{V}_{s,a,n}(\mu)|$
$\displaystyle\leq\max_{0\leq\mu\leq
V+\|V\|e}|g^{V}_{s,a}(\mu)-g^{V}_{s,a,n}(\mu)|$
$\displaystyle=\max_{0\leq\mu\leq
V+\|V\|e}\bigg{|}\mathsf{P}^{a}_{s}(V-\mu)-\hat{\mathsf{P}}^{a}_{s,n}(V-\mu)-\bigg{(}\sqrt{\delta\text{Var}_{{\mathsf{P}^{a}_{s}}}(V-\mu)}-\sqrt{\delta\text{Var}_{\hat{\mathsf{P}}^{a}_{s,n}}(V-\mu)}\bigg{)}\bigg{|}$
$\displaystyle\leq\max_{0\leq\mu\leq
V+\|V\|e}|\mathsf{P}^{a}_{s}(V-\mu)-\hat{\mathsf{P}}^{a}_{s,n}(V-\mu)|+\max_{0\leq\mu\leq
V+\|V\|e}\bigg{|}\bigg{(}\sqrt{\delta\text{Var}_{{\mathsf{P}^{a}_{s}}}(V-\mu)}-\sqrt{\delta\text{Var}_{\hat{\mathsf{P}}^{a}_{s,n}}(V-\mu)}\bigg{)}\bigg{|}$
$\displaystyle\overset{(a)}{\leq}\max_{0\leq\mu\leq
V+\|V\|e}\|V-\mu\|\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1}+\max_{0\leq\mu\leq
V+\|V\|e}\sqrt{|\delta\text{Var}_{{\mathsf{P}^{a}_{s}}}(V-\mu)-\delta\text{Var}_{\hat{\mathsf{P}}^{a}_{s,n}}(V-\mu)|},$
(159)

where $(a)$ is due to $|\sqrt{x}-\sqrt{y}|\leq\sqrt{|x-y|}$. Note that for any
distribution $p,q\in\Delta(|\mathcal{S}|)$ and any random variable $X$,

$\displaystyle|\text{Var}_{p}[X]-\text{Var}_{q}[X]|$
$\displaystyle=|\mathbb{E}_{p}[X^{2}]-\mathbb{E}_{p}[X]^{2}-\mathbb{E}_{q}[X^{2}]+\mathbb{E}_{q}[X]^{2}|$
$\displaystyle\leq|\mathbb{E}_{p}[X^{2}]-\mathbb{E}_{q}[X^{2}]|+|(\mathbb{E}_{p}[X]+\mathbb{E}_{q}[X])(\mathbb{E}_{p}[X]-\mathbb{E}_{q}[X])|$
$\displaystyle\leq\sup|X^{2}|\|p-q\|_{1}+2(\sup|X|)^{2}\|p-q\|_{1}.$ (160)

Hence,

$\displaystyle\sqrt{|\delta\text{Var}_{{\mathsf{P}^{a}_{s}}}(V-\mu)-\delta\text{Var}_{\hat{\mathsf{P}}^{a}_{s,n}}(V-\mu)|}$
$\displaystyle\leq\sqrt{3\delta\|V-\mu\|^{2}\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1}}.$
(161)

Thus, by Hoeffding’s inequality and Theorem 3.7 from (Liu et al., 2022),

$\displaystyle\mathbb{E}[|g^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a})|]\leq
3\|V\|\left(\frac{|\mathcal{S}|^{2}\sqrt{\pi}}{2^{\frac{n+1}{2}}}+\sqrt{\frac{3\delta|\mathcal{S}|^{2}\sqrt{\pi}}{2^{\frac{n+1}{2}}}}\right),$
(162)

which implies that

$\displaystyle\lim_{n\to\infty}\mathbb{E}\bigg{[}g^{V}_{s,a,n}(\mu^{V}_{s,a,n})\bigg{]}=g^{V}_{s,a}(\mu^{V}_{s,a}),$
(163)

which completes the proof. ∎

######  Theorem D.5.

The estimated operator $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ defined in (D.3)
has bounded variance, i.e., there exists a constant $C_{0}$, such that

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\leq(1+18(1+\sqrt{2\delta})^{2}+2C_{0})\|V\|^{2}.$
(164)

###### Proof.

We have that

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)$
$\displaystyle=\mathbb{E}[(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)^{2}]-\sigma_{\mathcal{P}^{a}_{s}}(V)^{2}$
$\displaystyle\leq\mathbb{E}\bigg{[}\bigg{(}g^{V}_{s,a,0}(\mu^{V}_{s,a,0})+\frac{\Delta_{N}(V)}{p_{N}}\bigg{)}^{2}\bigg{]}+(\sigma_{\mathcal{P}^{a}_{s}}(V))^{2}$
$\displaystyle\leq
2\mathbb{E}\bigg{[}\bigg{(}g^{V}_{s,a,0}(\mu^{V}_{s,a,0})\bigg{)}^{2}\bigg{]}+2\mathbb{E}\bigg{[}\bigg{(}\frac{\Delta_{N}(V)}{p_{N}}\bigg{)}^{2}\bigg{]}+(\sigma_{\mathcal{P}^{a}_{s}}(V))^{2}$
$\displaystyle\leq(1+18(1+\sqrt{2\delta})^{2})\|V\|^{2}+2\sum^{\infty}_{i=0}\frac{\mathbb{E}[(\Delta_{i}(V))^{2}]}{p_{i}},$
(165)

where the last inequality is from Lemma F.4. For any $n\geq 1$, we have that

$\displaystyle\mathbb{E}[(\Delta_{n}(V))^{2}]$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}{g}^{V}_{s,a,n}(\mu^{V}_{s,a,n})-\frac{{g}^{V}_{s,a,n,E}(\mu^{V}_{s,a,n,E})+{g}^{V}_{s,a,n,O}(\mu^{V}_{s,a,n,O})}{2}\bigg{)}^{2}\bigg{]}$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}{g}^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a})+g^{V}_{s,a}(\mu^{V}_{s,a})-\frac{{g}^{V}_{s,a,n,E}(\mu^{V}_{s,a,n,E})+{g}^{V}_{s,a,n,O}(\mu^{V}_{s,a,n,O})}{2}\bigg{)}^{2}\bigg{]}$
$\displaystyle\leq
2\mathbb{E}[({g}^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a}))^{2}]+2\mathbb{E}\bigg{[}\bigg{(}g^{V}_{s,a}(\mu^{V}_{s,a})-\frac{{g}^{V}_{s,a,n,E}(\mu^{V}_{s,a,n,E})+{g}^{V}_{s,a,n,O}(\mu^{V}_{s,a,n,O})}{2}\bigg{)}^{2}\bigg{]}$
$\displaystyle\overset{(a)}{=}2\mathbb{E}[({g}^{V}_{s,a,n}(\mu^{V}_{s,a,n})-g^{V}_{s,a}(\mu^{V}_{s,a}))^{2}]+2\mathbb{E}[({g}^{V}_{s,a,n-1}(\mu^{V}_{s,a,n-1})-g^{V}_{s,a}(\mu^{V}_{s,a}))^{2}]$
$\displaystyle\leq
18(1+\sqrt{3\delta})^{2}\|V\|^{2}\mathbb{E}[\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1}^{2}+\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1}]$
$\displaystyle\quad+18(1+\sqrt{3\delta})^{2}\|V\|^{2}\mathbb{E}[\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n-1}\|_{1}^{2}+\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n-1}\|_{1}],$
(166)

where $(a)$ is due to Lemma F.2 and the last inequality follows a similar
argument to (D.3). Note that $p_{n}=\Psi(1-\Psi)^{n}$ for
$\Psi\in\big{(}0,1-\frac{\sqrt{2}}{2}\big{)}$. Thus, similar to Theorem 3.7 of
(Liu et al., 2022), we can show that there exists a constant $C_{0}$, such
that

$\displaystyle\sum^{\infty}_{i=0}\frac{\mathbb{E}[(\Delta_{i}(V))^{2}]}{p_{i}}\leq
C_{0}\|V\|^{2}.$ (167)

Thus,

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\leq(1+18(1+\sqrt{2\delta})^{2})\|V\|^{2}+2C_{0}\|V\|^{2}=(1+18(1+\sqrt{2\delta})^{2}+2C_{0})\|V\|^{2}.$
(168)

∎

###  D.4 KL-Divergence Uncertainty Sets

The estimator under the KL-Divergence uncertainty set can be written as

$\displaystyle\hat{\sigma}_{\mathcal{P}^{a}_{s}}V$
$\displaystyle\triangleq-\min_{\alpha\geq
0}\bigg{(}\delta\alpha+\alpha\log\left(e^{\frac{-V(s^{\prime}_{1})}{\alpha}}\right)\bigg{)}+\frac{\Delta_{N}(V)}{p_{N}},$

where

$\displaystyle\Delta_{N}(V)=-\min_{\alpha\geq
0}\left(\delta\alpha+\alpha\log\left(\mathbb{E}_{\hat{\mathsf{P}}^{a}_{s,N+1}}\big{[}e^{\frac{-V}{\alpha}}\big{]}\right)\right)$
$\displaystyle\quad+\frac{1}{2}\min_{\alpha\geq
0}\left(\delta\alpha+\alpha\log\left(\mathbb{E}_{\hat{\mathsf{P}}^{a,O}_{s,N+1}}\big{[}e^{\frac{-V}{\alpha}}\big{]}\right)\right)$
$\displaystyle\quad+\frac{1}{2}\min_{\alpha\geq
0}\left(\delta\alpha+\alpha\log\left(\mathbb{E}_{\hat{\mathsf{P}}^{a,E}_{s,N+1}}\big{[}e^{\frac{-V}{\alpha}}\big{]}\right)\right).$
(169)

######  Theorem D.6.

(Liu et al., 2022) The estimated operator $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$
is unbiased and has bounded variance, i.e., there exists a constant $C_{0}$,
such that $\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\leq
C_{0}(1+\|V\|^{2}).$

###  D.5 Wasserstein Distance Uncertainty Sets

To study the support function w.r.t. this uncertainty model, we first
introduce some notation.

######  Definition D.1.

For any function $f:\mathcal{Z}\to\mathbb{R}$, $\lambda\geq 0$ and
$x\in\mathcal{Z}$, define the regularization operator

$\displaystyle\Phi(\lambda,x)\triangleq\inf_{y\in\mathcal{Z}}(\lambda
d(x,y)^{l}+f(y)).$ (170)

The growth rate $\kappa$ of function $f$ and any distribution $q$ over
$\mathcal{Z}$ is defined as

$\displaystyle\kappa_{q}\triangleq\inf\bigg{(}\lambda\geq
0:\sum_{x\in\mathcal{Z}}q(x)\Phi(\lambda,x)>-\infty\bigg{)}.$ (171)

######  Lemma D.2.

(Gao & Kleywegt, 2022) Consider the distributional robust optimization of a
function $f$:

$\displaystyle\inf_{W_{l}(q,p)\leq\delta}\mathbb{E}_{x\sim q}[f(x)],$ (172)

and define its dual problem as

$\displaystyle\sup_{\lambda\geq
0}(-\lambda\delta^{l}+\sum_{x\in\mathcal{Z}}p(x)\inf_{y\in\mathcal{Z}}(f(y)+\lambda
d(x,y)^{l})).$ (173)

If $\kappa_{p}<\infty$, then the strong duality holds, i.e.,

$\displaystyle\inf_{W_{l}(q,p)\leq\delta}\mathbb{E}_{x\sim
q}[f(x)]=\sup_{\lambda\geq
0}(-\lambda\delta^{l}+\sum_{x\in\mathcal{Z}}p(x)\inf_{y\in\mathcal{Z}}(f(y)+\lambda
d(x,y)^{l})).$ (174)

We first verify that this strong duality holds for our support function.

######  Lemma D.3.

(Restatement of Equation 19) It holds that

$\displaystyle\sigma_{\mathcal{P}^{a}_{s}}(V)=\sup_{\lambda\geq
0}\bigg{(}-\lambda\delta^{l}+\sum_{x}\mathsf{P}^{a}_{s,x}\inf_{y}(V(y)+\lambda
d(x,y)^{l})\bigg{)}.$ (175)

###### Proof.

In our case, the regularization operator is

$\displaystyle\Phi(\lambda,x)=\inf_{s\in\mathcal{S}}(\lambda
d(s,x)^{l}+V(s)).$ (176)

Note that for any $\lambda\geq 0$,

$\displaystyle\sum_{x\in\mathcal{S}}\mathsf{P}^{a}_{s}(x)\Phi(\lambda,x)=\sum_{x\in\mathcal{S}}\mathsf{P}^{a}_{s}(x)\inf_{s\in\mathcal{S}}(\lambda
d(s,x)^{l}+V(s))\geq-\|V\|>-\infty.$ (177)

Hence, the growth rate $\kappa_{\mathsf{P}^{a}_{s}}=0<\infty$. Thus, the
strong duality holds. ∎

Then, the estimator under the Wasserstein distance uncertainty set can be
constructed as

$\displaystyle\hat{\sigma}_{\mathcal{P}^{a}_{s}}V$
$\displaystyle\triangleq\sup_{\lambda\geq
0}\bigg{(}-\lambda\delta^{l}+\inf_{y}(V(y)+\lambda
d(s_{1}^{\prime},y)^{l})\bigg{)}+\frac{\Delta_{N}(V)}{p_{N}}+r(s,a),$ (178)

where

$\displaystyle\Delta_{N}(V)$ $\displaystyle=\sup_{\lambda\geq
0}\bigg{(}-\lambda\delta^{l}+\mathbb{E}_{\hat{\mathsf{P}}^{a}_{s,N+1}}\bigg{[}\inf_{y}(V(y)+\lambda
d(S,y)^{l})\bigg{]}\bigg{)}$ $\displaystyle\quad-\sup_{\lambda\geq
0}\bigg{(}-\lambda\delta^{l}+\mathbb{E}_{\hat{\mathsf{P}}^{a,O}_{s,N+1}}\bigg{[}\inf_{y}(V(y)+\lambda
d(S,y)^{l})\bigg{]}\bigg{)}$ $\displaystyle\quad-\sup_{\lambda\geq
0}\bigg{(}-\lambda\delta^{l}+\mathbb{E}_{\hat{\mathsf{P}}^{a,E}_{s,N+1}}\bigg{[}\inf_{y}(V(y)+\lambda
d(S,y)^{l})\bigg{]}\bigg{)}.$

######  Theorem D.7.

The estimated operator defined in (178) is unbiased, i.e.,

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]=\sigma_{\mathcal{P}^{a}_{s}}(V).$
(179)

###### Proof.

Denote the dual function (19) by $g$:

$\displaystyle
g_{s,a}^{V}(\lambda)=-\lambda\delta^{l}+\mathbb{E}_{S\sim\mathsf{P}^{a}_{s}}[\inf_{x\in\mathcal{S}}(V(x)+\lambda
d(S,x)^{l})],$ (180)

and denote its optimal solution by $\lambda^{V}_{s,a}$:

$\displaystyle\lambda^{V}_{s,a}=\arg\max_{\lambda\geq
0}\bigg{(}-\lambda\delta^{l}+\mathbb{E}_{S\sim\mathsf{P}^{a}_{s}}[\inf_{x\in\mathcal{S}}(V(x)+\lambda
d(S,x)^{l})]\bigg{)}.$ (181)

Then, the support function
$\sigma_{\mathcal{P}^{a}_{s}}(V)=g_{s,a}^{V}(\lambda^{V}_{s,a})$. Similarly,
define the empirical function
${g}_{s,a,N+1}^{V},{g}_{s,a,N+1,O}^{V},{g}_{s,a,N+1,E}^{V}$, and denote their
optimal solutions by
$\lambda^{V}_{s,a,N+1},\lambda^{V}_{s,a,N+1,O},\lambda^{V}_{s,a,N+1,E}$. We
have that

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\lambda^{V}_{s,a,0})]+\mathbb{E}\bigg{[}\frac{\Delta_{N}(V)}{p_{N}}\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\lambda^{V}_{s,a,0})]+\sum^{\infty}_{n=0}p(N=n)\mathbb{E}\bigg{[}\frac{\Delta_{N}(V)}{p_{N}}|N=n\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\lambda^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}[\Delta_{n}]$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\lambda^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}\bigg{[}g^{V}_{s,a,n+1}(\lambda^{V}_{s,a,n+1})-\frac{g^{V}_{s,a,n+1,O}(\lambda^{V}_{s,a,n+1,O})+g^{V}_{s,a,n+1,E}(\lambda^{V}_{s,a,n+1,E})}{2}\bigg{]}$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\lambda^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}\bigg{[}g^{V}_{s,a,n+1}(\lambda^{V}_{s,a,n+1})-g^{V}_{s,a,n}(\lambda^{V}_{s,a,n})\bigg{]},$
(182)

where the last inequality is from Lemma F.2. The last equation can be further
rewritten as

$\displaystyle\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]$
$\displaystyle=\mathbb{E}[g^{V}_{s,a,0}(\lambda^{V}_{s,a,0})]+\sum^{\infty}_{n=0}\mathbb{E}\bigg{[}g^{V}_{s,a,n+1}(\lambda^{V}_{s,a,n+1})-g^{V}_{s,a,n}(\lambda^{V}_{s,a,n})\bigg{]}$
$\displaystyle=\lim_{n\to\infty}\mathbb{E}\bigg{[}g^{V}_{s,a,n}(\lambda^{V}_{s,a,n})\bigg{]}.$
(183)

To show that $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ is unbiased, it suffices to
prove that

$\displaystyle\lim_{n\to\infty}\mathbb{E}\bigg{[}g^{V}_{s,a,n}(\lambda^{V}_{s,a,n})\bigg{]}=g^{V}_{s,a}(\lambda^{V}_{s,a}).$
(184)

For any arbitrary i.i.d. samples $\\{X_{i}\\}$ and its corresponding function
$g^{V}_{s,a,n}$, together with Lemma F.5, we have that

$\displaystyle|g^{V}_{s,a,n}(\lambda^{V}_{s,a,n})-g^{V}_{s,a}(\lambda^{V}_{s,a})|$
$\displaystyle=|\max_{0\leq\lambda\leq\frac{2\|V\|}{\delta^{l}}}g^{V}_{s,a}(\lambda)-\max_{0\leq\lambda\leq\frac{2\|V\|}{\delta^{l}}}g^{V}_{s,a,n}(\lambda)|$
$\displaystyle\leq\max_{0\leq\lambda\leq\frac{2\|V\|}{\delta^{l}}}|g^{V}_{s,a}(\lambda)-g^{V}_{s,a,n}(\lambda)|$
$\displaystyle=\max_{0\leq\lambda\leq\frac{2\|V\|}{\delta^{l}}}\bigg{|}\mathbb{E}_{S\sim\mathsf{P}^{a}_{s}}[\inf_{x\in\mathcal{S}}(V(x)+\lambda
d(S,x)^{l})]-\mathbb{E}_{S\sim\hat{\mathsf{P}}^{a}_{s,n}}[\inf_{x\in\mathcal{S}}(V(x)+\lambda
d(S,x)^{l})]\bigg{|}$
$\displaystyle\leq\max_{0\leq\lambda\leq\frac{2\|V\|}{\delta^{l}}}\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1}\sup_{x,S\in\mathcal{S}}(|V(x)+\lambda
d(S,x)^{l}|)$
$\displaystyle\leq\bigg{(}1+\frac{2D^{l}}{\delta^{l}}\bigg{)}\|V\|\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1},$
(185)

where the last inequality is from the bound on $\lambda$ and $D$ is the
diameter of the metric space $(\mathcal{S},d)$.

By Hoeffding’s inequality and similar to the previous proofs, we have that

$\displaystyle\mathbb{E}[|g^{V}_{s,a,n}(\lambda^{V}_{s,a,n})-g^{V}_{s,a}(\lambda^{V}_{s,a})|]\leq\bigg{(}1+\frac{2D^{l}}{\delta^{l}}\bigg{)}\bigg{(}\frac{|\mathcal{S}|^{2}\sqrt{\pi}}{2^{\frac{n+1}{2}}}\bigg{)}\|V\|,$
(186)

which implies that

$\displaystyle\lim_{n\to\infty}\mathbb{E}\bigg{[}g^{V}_{s,a,n}(\lambda^{V}_{s,a,n})\bigg{]}=g^{V}_{s,a}(\lambda^{V}_{s,a}).$
(187)

This completes the proof. ∎

######  Theorem D.8.

The estimated operator $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ defined in (178)
has bounded variance, i.e., there exists a constant $C_{0}$, such that

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\leq(3+2C_{0})\|V\|^{2}.$
(188)

###### Proof.

We first have that

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)$
$\displaystyle=\mathbb{E}[(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)^{2}]-\sigma_{\mathcal{P}^{a}_{s}}(V)^{2}$
$\displaystyle\leq\mathbb{E}\bigg{[}\bigg{(}g^{V}_{s,a,0}(\lambda^{V}_{s,a,0})+\frac{\Delta_{N}(V)}{p_{N}}\bigg{)}^{2}\bigg{]}+(\sigma_{\mathcal{P}^{a}_{s}}(V))^{2}$
$\displaystyle\leq
2\mathbb{E}\bigg{[}\bigg{(}g^{V}_{s,a,0}(\lambda^{V}_{s,a,0})\bigg{)}^{2}\bigg{]}+2\mathbb{E}\bigg{[}\bigg{(}\frac{\Delta_{N}(V)}{p_{N}}\bigg{)}^{2}\bigg{]}+(\sigma_{\mathcal{P}^{a}_{s}}(V))^{2}$
$\displaystyle\leq
3\|V\|^{2}+2\sum^{\infty}_{i=0}\frac{\mathbb{E}[(\Delta_{i}(V))^{2}]}{p_{i}},$
(189)

where the last inequality is from Lemma F.4. For any $n\geq 1$, we have that

$\displaystyle\mathbb{E}[(\Delta_{n}(V))^{2}]$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}{g}^{V}_{s,a,n}(\lambda^{V}_{s,a,n})-\frac{{g}^{V}_{s,a,n,E}(\lambda^{V}_{s,a,n,E})+{g}^{V}_{s,a,n,O}(\lambda^{V}_{s,a,n,O})}{2}\bigg{)}^{2}\bigg{]}$
$\displaystyle=\mathbb{E}\bigg{[}\bigg{(}{g}^{V}_{s,a,n}(\lambda^{V}_{s,a,n})-g^{V}_{s,a}(\lambda^{V}_{s,a})+g^{V}_{s,a}(\lambda^{V}_{s,a})-\frac{{g}^{V}_{s,a,n,E}(\lambda^{V}_{s,a,n,E})+{g}^{V}_{s,a,n,O}(\lambda^{V}_{s,a,n,O})}{2}\bigg{)}^{2}\bigg{]}$
$\displaystyle\leq
2\mathbb{E}[({g}^{V}_{s,a,n}(\lambda^{V}_{s,a,n})-g^{V}_{s,a}(\lambda^{V}_{s,a}))^{2}]+2\mathbb{E}\bigg{[}\bigg{(}g^{V}_{s,a}(\lambda^{V}_{s,a})-\frac{{g}^{V}_{s,a,n,E}(\lambda^{V}_{s,a,n,E})+{g}^{V}_{s,a,n,O}(\lambda^{V}_{s,a,n,O})}{2}\bigg{)}^{2}\bigg{]}$
$\displaystyle\overset{(a)}{=}2\mathbb{E}[({g}^{V}_{s,a,n}(\lambda^{V}_{s,a,n})-g^{V}_{s,a}(\lambda^{V}_{s,a}))^{2}]+2\mathbb{E}[({g}^{V}_{s,a,n-1}(\lambda^{V}_{s,a,n-1})-g^{V}_{s,a}(\lambda^{V}_{s,a}))^{2}]$
$\displaystyle\leq
2\bigg{(}1+\frac{2D^{l}}{\delta^{l}}\bigg{)}^{2}\|V\|^{2}\mathbb{E}[\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n}\|_{1}^{2}]+2\bigg{(}1+\frac{2D^{l}}{\delta^{l}}\bigg{)}^{2}\|V\|^{2}\mathbb{E}[\|\mathsf{P}^{a}_{s}-\hat{\mathsf{P}}^{a}_{s,n-1}\|_{1}^{2}],$
(190)

where $(a)$ is due to Lemma F.2 and the last inequality follows a similar
argument to (186). Note that $p_{n}=\Psi(1-\Psi)^{n}$ for
$\Psi\in\big{(}0,0.5\big{)}$, thus similar to Theorem 3.7 of (Liu et al.,
2022), we can show that there exists a constant $C_{0}$, such that

$\displaystyle\sum^{\infty}_{i=0}\frac{\mathbb{E}[(\Delta_{i}(V))^{2}]}{p_{i}}\leq
C_{0}\|V\|^{2}.$ (191)

Thus, we have that

$\displaystyle\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\leq
3\|V\|^{2}+2C_{0}\|V\|^{2}=(3+2C_{0})\|V\|^{2}.$ (192)

∎

##  Appendix E Case Studies for Robust RVI Q-Learning

In this section, we provide the proof of the second part of Theorem 5.1, i.e.,
$\hat{\mathbf{H}}$ is bounded and unbiased under each uncertainty model. We
note that the proofs in this part can be easily derived by following the ones
in Appendix D.

We first prove a lemma necessary to the proofs in this section.

######  Lemma E.1.

It holds that

$\displaystyle\|V_{Q}\|\leq\|Q\|.$ (193)

###### Proof.

From the definition of $V_{Q}$, we have that

$\displaystyle\|V_{Q}\|=\max_{s}|V_{Q}(s)|=\max_{s}|\max_{a}Q(s,a)|\triangleq|Q(s^{*},a^{*})|.$
(194)

Clearly, $|Q(s^{*},a^{*})|\leq\max_{s,a}|Q(s,a)|$, hence

$\displaystyle\|V_{Q}\|\leq\|Q\|.$ (195)

∎

Similar to Appendix D, the propositions of $\hat{\mathbf{H}}$ can be reduced
to the ones of $\hat{\sigma}_{\mathcal{P}^{a}_{s}}$.

######  Lemma E.2.

If
$\mathbb{E}[\hat{\sigma}_{\mathcal{P}^{a}_{s}}V]=\sigma_{\mathcal{P}^{a}_{s}}(V)$,
and moreover there exists a constant $C$, such that for any $s,a$,
$\text{Var}(\hat{\sigma}_{\mathcal{P}^{a}_{s}}V)\leq C(1+\|V\|^{2})$, then
$\mathbb{E}[\hat{\mathbf{H}}Q(s,a)]=\mathbf{H}Q(s,a),$ and
$\text{Var}(\hat{\mathbf{H}}Q(s,a))\leq C(1+\|Q\|^{2}).$

###### Proof.

First, we have that

$\displaystyle\mathbb{E}[\hat{\mathbf{H}}Q(s,a)]=\mathbb{E}[r(s,a)+\hat{\sigma}_{\mathcal{P}^{a}_{s}}V_{Q}(s)]=r(s,a)+\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})=\mathbf{H}Q(s,a).$
(196)

For boundedness, note that

$\displaystyle\text{Var}(\hat{\mathbf{H}}Q(s,a))$
$\displaystyle=\mathbb{E}\bigg{[}(\hat{\mathbf{H}}Q(s,a)-\mathbf{H}Q(s,a))^{2}\bigg{]}$
$\displaystyle=\mathbb{E}\big{[}\big{(}\hat{\sigma}_{\mathcal{P}^{a}_{s}}V_{Q}(s)-\sigma_{\mathcal{P}^{a}_{s}}(V_{Q})\big{)}^{2}\big{]}$
$\displaystyle\leq C(1+\|V_{Q}\|^{2})$ $\displaystyle\leq C(1+\|Q\|^{2}),$
(197)

where the last inequality is from Lemma E.1. ∎

This implies that the problem is reduced to verifying whether
$\hat{\sigma}_{\mathcal{P}^{a}_{s}}$ is unbiased and has bounded variance,
which is identical to the results in Appendix D. We thus omit the proofs for
this part.

##  Appendix F Technical Lemmas

######  Lemma F.1.

For a robust Bellman operator $\mathbf{T}$, define

$\displaystyle{\mathbf{T}}^{\prime}V$
$\displaystyle\triangleq{\mathbf{T}}V-f(V)e,$ (198)
$\displaystyle\tilde{\mathbf{T}}V$
$\displaystyle\triangleq{\mathbf{T}}V-g^{\pi}_{\mathcal{P}}e.$ (199)

Assume that $x(t),y(t)$ are the solutions to equations

$\displaystyle\dot{x}$ $\displaystyle={\mathbf{T}}^{\prime}x-x,$ (200)
$\displaystyle\dot{y}$ $\displaystyle=\tilde{\mathbf{T}}y-y,$ (201)

with the same initial value $x(0)=y(0)=x_{0}$. Then,

$\displaystyle x(t)=y(t)+r(t)e,$ (202)

where $r(t)$ satisfies

$\displaystyle\dot{r}(t)=-r(t)+g^{\pi}_{\mathcal{P}}-f(y(t)).$ (203)

###### Proof.

Note that
${\mathbf{T}}^{\prime}V=\tilde{\mathbf{T}}V+(g^{\pi}_{\mathcal{P}}-f(V))e$,
then from the variation of constants formula, we have that

$\displaystyle x(t)$
$\displaystyle=x_{0}e^{-t}+\int^{t}_{0}e^{-(t-s)}\tilde{\mathbf{T}}(x(s))ds+\left(\int^{t}_{0}e^{-(t-s)}(g^{\pi}_{\mathcal{P}}-f(x(s)))ds\right)e,$
(204) $\displaystyle y(t)$
$\displaystyle=x_{0}e^{-t}+\int^{t}_{0}e^{-(t-s)}\tilde{\mathbf{T}}(y(s))ds.$
(205)

Hence, the maximal and minimal components of $x(t)-y(t)$ can be bounded as:

$\displaystyle\max_{i}(x_{i}(t)-y_{i}(t))$
$\displaystyle\leq\int^{t}_{0}e^{-(t-s)}\max_{i}(\tilde{\mathbf{T}}_{i}(x(s))-\tilde{\mathbf{T}}_{i}(y(s)))ds+\int^{t}_{0}e^{-(t-s)}(g^{\pi}_{\mathcal{P}}-f(x(s)))ds,$
$\displaystyle\min_{i}(x_{i}(t)-y_{i}(t))$
$\displaystyle\geq\int^{t}_{0}e^{-(t-s)}\min_{i}(\tilde{\mathbf{T}}_{i}(x(s))-\tilde{\mathbf{T}}_{i}(y(s)))ds+\int^{t}_{0}e^{-(t-s)}(g^{\pi}_{\mathcal{P}}-f(x(s)))ds.$
(206)

This hence implies that

$\displaystyle\text{Span}(x(t)-y(t))$
$\displaystyle\leq\int^{t}_{0}e^{-(t-s)}\text{Span}(\tilde{\mathbf{T}}(x(s))-\tilde{\mathbf{T}}(y(s)))ds$
$\displaystyle\leq\int^{t}_{0}e^{-(t-s)}\text{Span}(x(s)-y(s))ds,$ (207)

where the last inequality is because $\tilde{\mathbf{T}}$ is non-expansive
w.r.t. the span semi-norm (Wang et al., 2023).

Gronwall’s inequality implies that $\text{Span}(x(t)-y(t))\leq
0\cdot\int^{t}_{0}e^{-(t-s)}ds=0$ for any $t\geq 0$. However, since Span is
non-negative, then $\text{Span}(x(t)-y(t))=0$. Hence, we have that
$x(t)=y(t)+r(t)e$ for some $r(t)$ satisfying $r(0)=0$.

Also note that the differential of $r(t)$ can be written as

$\displaystyle\dot{r}(t)e$ $\displaystyle=\dot{x}(t)-\dot{y}(t)$
$\displaystyle=\tilde{\mathbf{T}}x(t)+(g^{\pi}_{\mathcal{P}}-f(x(t)))e-x(t)-\tilde{\mathbf{T}}y(t)+y(t)$
$\displaystyle=(-r(t)+g^{\pi}_{\mathcal{P}}-f(y(t)))e,$ (208)

where the last equation is because

$\displaystyle\tilde{\mathbf{T}}x(t)$
$\displaystyle=\tilde{\mathbf{T}}(y(t)+r(t)e)=\tilde{\mathbf{T}}(y(t))+r(t)e,$
(209) $\displaystyle f(x(t))$ $\displaystyle=f(y(t)+r(t)e)=f(y(t))+r(t).$
(210)

This completes the proof. ∎

######  Lemma F.2.

For any function $g:\Delta(|\mathcal{S}|)\to\mathbb{R}$, assume there are
$2^{n+1}$ i.i.d. samples $X_{i}\sim q$. Denote the empirical distributions
from samples
$\\{X_{i}:i=1,...,2^{n+1}\\},\\{X_{2i-1}:i=1,...,2^{n}\\},\\{X_{2i}:i=1,...,2^{n}\\}$
by $\hat{q}_{n+1},\hat{q}_{n+1,O},\hat{q}_{n+1,E}$. Then,

$\displaystyle\mathbb{E}[g(\hat{q}_{n+1,O})]=\mathbb{E}[g(\hat{q}_{n+1,E})]=\mathbb{E}[g(\hat{q}_{n})].$
(211)

###### Proof.

Note that

$\displaystyle\hat{q}_{n+1,O}(s)=\frac{\sum_{i=1}^{2^{n}}\mathbbm{1}_{X_{2i-1}=s}}{2^{n}},$
(212)

hence,

$\displaystyle\mathbb{E}[g(\hat{q}_{n+1,O})]$
$\displaystyle=\sum_{p=(p_{1},...,p_{|\mathcal{S}|})\in\Delta(|\mathcal{S}|)}g(p)\mathbb{P}(\hat{q}_{n+1,O}=p)$
$\displaystyle=\sum_{p=(p_{1},...,p_{|\mathcal{S}|})\in\Delta(|\mathcal{S}|)}\mathbb{P}\bigg{(}\frac{\sum^{2^{n}}_{i=1}\mathbbm{1}_{X_{2i-1}=s_{1}}}{2^{n}}=p_{1},...,\frac{\sum^{2^{n}}_{i=1}\mathbbm{1}_{X_{2i-1}=s_{|\mathcal{S}|}}}{2^{n}}=p_{|\mathcal{S}|}\bigg{|}q\bigg{)}g(p),$
(213)

where $2^{n}p_{i}\in\mathbb{N}$ and $\sum_{i=1}^{|\mathcal{S}|}p_{i}=1$. On
the other hand,

$\displaystyle\mathbb{E}[g(\hat{q}_{n})]$
$\displaystyle=\sum_{p=(p_{1},...,p_{|\mathcal{S}|})\in\Delta(|\mathcal{S}|)}g(p)\mathbb{P}(\hat{q}_{n}=p)$
$\displaystyle=\sum_{p=(p_{1},...,p_{|\mathcal{S}|})\in\Delta(|\mathcal{S}|)}\mathbb{P}\bigg{(}\frac{\sum^{2^{n}}_{i=1}\mathbbm{1}_{X_{i}=s_{1}}}{2^{n}}=p_{1},...,\frac{\sum^{2^{n}}_{i=1}\mathbbm{1}_{X_{i}=s_{|\mathcal{S}|}}}{2^{n}}=p_{|\mathcal{S}|}\bigg{|}q\bigg{)}g(p).$
(214)

Note that $X_{i}$ are i.i.d., hence,

$\displaystyle\mathbb{P}\bigg{(}\frac{\sum^{2^{n}}_{i=1}\mathbbm{1}_{X_{i}=s_{1}}}{2^{n}}=p_{1},...,\frac{\sum^{2^{n}}_{i=1}\mathbbm{1}_{X_{i}=s_{|\mathcal{S}|}}}{2^{n}}=p_{|\mathcal{S}|}\bigg{|}q\bigg{)}=\mathbb{P}\bigg{(}\frac{\sum^{2^{n}}_{i=1}\mathbbm{1}_{X_{2i-1}=s_{1}}}{2^{n}}=p_{1},...,\frac{\sum^{2^{n}}_{i=1}\mathbbm{1}_{X_{2i-1}=s_{|\mathcal{S}|}}}{2^{n}}=p_{|\mathcal{S}|}\bigg{|}q\bigg{)}.$

Thus,

$\displaystyle\mathbb{E}[g(\hat{q}_{n+1,O})]=\mathbb{E}[g(\hat{q}_{n})].$
(215)

Similarly, $\mathbb{E}[g(\hat{q}_{n+1,E})]=\mathbb{E}[g(\hat{q}_{n})]$ and
hence it completes the proof. ∎

######  Lemma F.3.

Under the total variation uncertainty model, the optimal solution and optimal
value for $g_{s,a}^{V},g_{s,a,N+1}^{V},g_{s,a,N+1,E}^{V},g_{s,a,N+1,O}^{V}$
are bounded. Specifically,

$\displaystyle\mu^{V}_{s,a},\mu^{V}_{s,a,N+1},\mu^{V}_{s,a,N+1,E},\mu^{V}_{s,a,N+1,O}$
$\displaystyle\leq V+\|V\|e,$ (216)
$\displaystyle\|\mu^{V}_{s,a}\|,\|\mu^{V}_{s,a,N+1}\|,\|\mu^{V}_{s,a,N+1,E}\|,\|\mu^{V}_{s,a,N+1,O}\|$
$\displaystyle\leq 2\|V\|,$ (217)
$\displaystyle|g_{s,a}^{V}(\mu^{V}_{s,a})|,|g_{s,a,N+1}^{V}(\mu^{V}_{s,a,N+1})|,|g_{s,a,N+1,O}^{V}(\mu^{V}_{s,a,N+1,O})|,|g_{s,a,N+1,E}^{V}(\mu^{V}_{s,a,N+1,E})|$
$\displaystyle\leq 3(1+2\delta)\|V\|.$ (218)

###### Proof.

First we show the bounds on the optimal solutions. If we denote the minimal
entry of $V$ by $w$: $w=\min_{s}V(s)$, then $W\triangleq V-we\geq 0$. Note
that,

$\displaystyle\mu^{W}_{s,a}$ $\displaystyle=\arg\max_{\mu\geq
0}\big{(}\mathsf{P}^{a}_{s}(W-\mu)-\delta\text{Span}(W-\mu)\big{)}$
$\displaystyle=\arg\max_{\mu\geq
0}\big{(}-w+\mathsf{P}^{a}_{s}(V-\mu)-\delta\text{Span}(V-\mu)\big{)},$ (219)

which is because $\text{Span}(V+ke)=\text{Span}(V)$ and
$\mathsf{P}^{a}_{s}(V+ke)=k+\mathsf{P}^{a}_{s}V$. Hence,
$\mu^{W}_{s,a}=\mu^{V}_{s,a}$. Moreover note that $W\geq 0$, hence
$\mu^{W}_{s,a}$ is bounded: $\mu^{W}_{s,a}\leq W$, this further implies that

$\displaystyle\|\mu^{V}_{s,a}\|=\|\mu^{W}_{s,a}\|\leq\|W\|\leq 2\|V\|.$ (220)

The bounds on $\mu^{V}_{s,a,N+1},\mu^{V}_{s,a,N+1,O},\mu^{V}_{s,a,N+1,E}$ can
be similarly derived.

We then consider the optimal value. Note that,

$\displaystyle g_{s,a}^{V}(\mu^{V}_{s,a})$
$\displaystyle=\mathsf{P}^{a}_{s}(V-\mu^{V}_{s,a})-\delta\text{Span}(V-\mu^{V}_{s,a})$
$\displaystyle\leq\|V\|+\|\mu^{V}_{s,a}\|+\delta|\max_{i}(V(i)-\mu^{V}_{s,a}(i))|+\delta|\min_{i}(V(i)-\mu^{V}_{s,a}(i))|$
$\displaystyle\leq 3\|V\|+2\delta(\|V\|+\|\mu^{V}_{s,a}\|)$ $\displaystyle\leq
3(1+2\delta)\|V\|.$ (221)

On the other hand,

$\displaystyle g_{s,a}^{V}(\mu^{V}_{s,a})\geq
g_{s,a}^{V}(0)=\mathsf{P}^{a}_{s}V-\delta\text{Span}(V)=\mathsf{P}^{a}_{s}V-\delta\max_{i}V(i)+\delta\min_{i}V(i).$
(222)

Denote the maximal and minimal entries of $V$ by $V(M)$ and $V(m)$, then we
have that

$\displaystyle\mathsf{P}^{a}_{s}V-\delta\max_{i}V(i)+\delta\min_{i}V(i)$
$\displaystyle=\sum_{x}\mathsf{P}^{a}_{s,x}V(x)-\delta V(M)+\delta V(m)$
$\displaystyle\geq-\|V\|-2\delta\|V\|,$ (223)

where the last inequality is from $\|V\|\geq V(i)\geq-\|V\|$ for any entry
$i$. Thus, combining (222) and (F) implies that

$\displaystyle-(1+2\delta)\|V\|\leq g_{s,a}^{V}(\mu^{V}_{s,a})\leq
3(1+2\delta)\|V\|.$ (224)

Similarly, the bounds on
$g_{s,a,N+1}^{V}(\mu^{V}_{s,a,N+1}),g_{s,a,N+1,O}^{V}(\mu^{V}_{s,a,N+1,O}),g_{s,a,N+1,E}^{V}(\mu^{V}_{s,a,N+1,E})$
can be derived. ∎

######  Lemma F.4.

Under the chi-square uncertainty model, the optimal solution and optimal value
for $g_{s,a}^{V},g_{s,a,N+1}^{V},g_{s,a,N+1,E}^{V},g_{s,a,N+1,O}^{V}$ are
bounded. Specifically,

$\displaystyle\mu^{V}_{s,a},\mu^{V}_{s,a,N+1},\mu^{V}_{s,a,N+1,E},\mu^{V}_{s,a,N+1,O}$
$\displaystyle\leq V+\|V\|e,$ (225)
$\displaystyle\|\mu^{V}_{s,a}\|,\|\mu^{V}_{s,a,N+1}\|,\|\mu^{V}_{s,a,N+1,E}\|,\|\mu^{V}_{s,a,N+1,O}\|$
$\displaystyle\leq 2\|V\|,$ (226)
$\displaystyle|g_{s,a}^{V}(\mu^{V}_{s,a})|,|g_{s,a,N+1}^{V}(\mu^{V}_{s,a,N+1})|,|g_{s,a,N+1,O}^{V}(\mu^{V}_{s,a,N+1,O})|,|g_{s,a,N+1,E}^{V}(\mu^{V}_{s,a,N+1,E})|$
$\displaystyle\leq 3(1+\sqrt{2\delta})\|V\|.$ (227)

###### Proof.

First, we show the bounds on the optimal solutions. If we denote the minimal
entry of $V$ by $w$: $w=\min_{s}V(s)$, then $W\triangleq V-we\geq 0$. Note
that,

$\displaystyle\mu^{W}_{s,a}$ $\displaystyle=\arg\max_{W\geq\mu\geq
0}\big{(}\mathsf{P}^{a}_{s}(W-\mu)-\sqrt{\delta\text{Var}_{\mathsf{P}^{a}_{s}}(W-\mu)}\big{)}$
$\displaystyle=\arg\max_{W\geq\mu\geq
0}\big{(}-w+\mathsf{P}^{a}_{s}(V-\mu)-\sqrt{\delta\text{Var}_{\mathsf{P}^{a}_{s}}(V-\mu)}\big{)},$
(228)

which is because $\text{Var}_{\mathsf{P}^{a}_{s}}(V-\mu-
we)=\text{Var}_{\mathsf{P}^{a}_{s}}(V-\mu)+\text{Var}_{\mathsf{P}^{a}_{s}}(we)-2\text{Cov}_{\mathsf{P}^{a}_{s}}(V-\mu,we)=\text{Var}_{\mathsf{P}^{a}_{s}}(V-\mu)$.
Hence $\mu^{W}_{s,a}=\mu^{V}_{s,a}$. Moreover note that $W\geq 0$, hence
$\mu^{W}_{s,a}$ is bounded: $\mu^{W}_{s,a}\leq W$, this further implies that

$\displaystyle\|\mu^{V}_{s,a}\|=\|\mu^{W}_{s,a}\|\leq\|W\|\leq 2\|V\|.$ (229)

The bounds on $\mu^{V}_{s,a,N+1},\mu^{V}_{s,a,N+1,O},\mu^{V}_{s,a,N+1,E}$ can
be similarly derived.

We then consider the optimal value. Note that,

$\displaystyle|g_{s,a}^{V}(\mu^{V}_{s,a})|$
$\displaystyle=|\mathsf{P}^{a}_{s}(V-\mu^{V}_{s,a})-\sqrt{\delta\text{Var}_{\mathsf{P}^{a}_{s}}(V-\mu^{V}_{s,a})}|$
$\displaystyle\leq\|V\|+\|\mu^{V}_{s,a}\|+\sqrt{2\delta\|V-\mu^{V}_{s,a}\|^{2}}$
$\displaystyle\leq 3\|V\|+\sqrt{2\delta}(\|V\|+\|\mu^{V}_{s,a}\|)$
$\displaystyle\leq 3(1+\sqrt{2\delta})\|V\|.$ (230)

Similarly, the bounds on
$g_{s,a,N+1}^{V}(\mu^{V}_{s,a,N+1}),g_{s,a,N+1,O}^{V}(\mu^{V}_{s,a,N+1,O}),g_{s,a,N+1,E}^{V}(\mu^{V}_{s,a,N+1,E})$
can be derived. ∎

######  Lemma F.5.

Under the Wasserstein distance uncertainty model, the optimal solution and
optimal value for
$g_{s,a}^{V},g_{s,a,N+1}^{V},g_{s,a,N+1,E}^{V},g_{s,a,N+1,O}^{V}$ are bounded.
Specifically,

$\displaystyle\lambda^{V}_{s,a},\lambda^{V}_{s,a,n},\lambda^{V}_{s,a,n,O},\lambda^{V}_{s,a,n,E}\leq\frac{2\|V\|}{\delta^{l}},$
(231)
$\displaystyle|g^{V}_{s,a}(\lambda^{V}_{s,a})|,|g^{V}_{s,a,n}(\lambda^{V}_{s,a,n})|,|g^{V}_{s,a,n,O}(\lambda^{V}_{s,a,n,O})|,|g^{V}_{s,a,n,E}(\lambda^{V}_{s,a,n,E})|\leq\|V\|.$
(232)

###### Proof.