











\section{Introduction}
Reinforcement learning (RL) has demonstrated success in a wide variety of applications such as AlphaGo~\cite{DBLP:journals/nature/SilverSSAHGHBLB17},  recommendor systems~\cite{DBLP:journals/corr/abs-1810-12027}, and robotics~\cite{hoel2019combining,chang2020sound,DBLP:conf/icml/HaarnojaZAL18,DBLP:journals/corr/LillicrapHPHETS15}, albiet with some concerns regarding their sensitivity to hyperparameters~\cite{DBLP:conf/aaai/0002IBPPM18}. Using such hyperparameter sensitive algorithms might lead to an RL agent learning sub-optimal policies, which gets further exacerbated when there is a mismatch between the training environment and the intended target environment. 

A typical approach to robot learning involves training and optimization of policy in a software based simulator of the same robot~\cite{levine2016end, DBLP:journals/corr/AndrychowiczWRS17}. A robot trained in such a manner tends to perform sub-optimally in the intended target environment~\cite{DBLP:conf/icra/PengAZA18}. The primary reason behind this sub-optimal performance that these simulators fail to capture the nuances and the real world environment accurately, which is often referred to as the Sim2Real gap~\cite{DBLP:conf/wafr/TzengDHFALSD16}. These nuances either require strong assumptions (e.g., modelling contact forces) or are extremely difficult to model (e.g., human-robot interaction). The resulting sub-optimality is not only hard to characterize but it also might lead to potentially hazardous scenarios in high risk environments~\cite{li2018learning}.

To account for such model mismatches, Robust RL algorithms learn a policy which is robust to synthetic perturbations in the simulator~\cite{DBLP:conf/iros/TobinFRSZA17, DBLP:journals/corr/TzengDHFPLSD15}. These synthetic perturbations are created under the assumption that the discrepancy between the simulator and the target domain is due to an misidentification of simulator parameters. These simulators are often governed by a model with a lower complexity than the one in the target domain or struggle to capture real-world disturbances~\cite{ma2018improved}. Applying robust RL techniques to account for significant modeling error or higher complexity altogether will likely not be enough to ensure a successful domain transfer without extensive real world validation~\cite{DBLP:journals/corr/abs-2012-03806}. We believe that one of the reasons such policies fail to work in practice is due to a fundamental misunderstanding and lack of formalism of model mismatch and policy evaluation across model or environment changes. For example, using robot learning to manipulate and solve a Rubik's cube required more than two years of research on developing a state-of-the-art robust learning algorithms as compared to six months it took them to learn the same task in simulation~\cite{DBLP:journals/corr/abs-1910-07113}.

We begin our work with an assumption that it is possible to collect data of agent behaviors from the target environment. In this paper, we use a small amount of data from the target domain to learn a mismatch function that modifies simulator performance to be more reflective of true performance in deployment. We quantify this mismatch function as the as a transition ratio estimation problem between the simulator and the target domain, both of which are modelled as Markov Decision Processes (MDP). In order to estimate this transition ratio, we further pose this estimation as an optimizer of a convex risk minimization algorithm that utilizes trajectory samples from both the target and the simulator environments. We then use this estimated ratio along with the simulator to evaluate a robot agent's performance in the target domain. 

We demonstrate that our estimator tends to approximate the true transition ratio with the accuracy error decreasing with the rate of $n^{-1/4}$, where $n$ is the number of transition samples used to estimate the ratio. We also illustrate the effectiveness of our algorithm to approximate the returns for a short time horizon by proving an exponential dependence of our estimator on the time horizon.

We present the following contributions:
(1) We introduce off-environment evaluation, which learns to incorporate model mismatch into the evaluation of a fixed policy across different environments; 
(2)  We use results from convex risk minimization to efficiently learn the model mismatch as the ratio of transition probabilities across two environments; 
(3) Using tools from statistical learning theory, we bound the error rate of our estimator as a function of the number of samples from the target domain; and 
(4) We demonstrate the utility of our approach on three task settings, a large discrete example (gridworld), and continuous environments (cartpole, reacher and kinova).

This paper is organized as follows. Section \ref{sec:litreview} reviews relevant literature. Section \ref{sec:prelim} introduces necessary preliminaries to pose our problem. Section \ref{sec: OEE} discusses our proposed method of off-environment policy evaluation. Section \ref{sec: error} proves the theoretical error bounds of our algorithm. We demonstrate our results by starting with a discrete-state, discrete-action MDP and building up to a real robot manipulator in Section \ref{sec: expt} and conclude our paper in Section \ref{sec: conclude}. 
\section{Related Work}
\label{sec:litreview}
In this section, we review relevant literature that considers model mismatch, covariate shift, and off-policy evaluation  which inspired our formulation of off-environment evaluation.
\subsection{Model Mismatch and Covariate Shift}
Model mismatch between two domains is well studied in machine learning, often in context of estimating covariate shift~\cite{ DBLP:journals/jmlr/AgarwalLS11, DBLP:journals/jmlr/SugiyamaKM07}. Covariate shift occurs when the training and the test data are sampled from different distributions. The goal of the above work is to learn a classifier that generalizes well on the test distribution while being learned entirely on the training dataset. 

The impact of covariate shift is often mitigated by estimating the distributional mismatch factor between the training and the testing distributions $\beta(x)$~\cite{conf/aaai/ReddiPS15}. The factor is then used to re-weight samples on the training dataset such that higher weights are assigned to data points which have a higher likelihood of being in the test distribution than the training distribution. 
Our work extends this idea of using a distributional mismatch factor to estimate mismatch between two Markov Decision Processes (MDP), not just distributions. 
\subsection{Off-Policy Evaluation}
Approximating performance by a target policy given access to data from a different behavior policy falls under the class of Off Policy Evaluation (OPE) in reinforcement learning~\cite{DBLP:journals/corr/abs-1911-06854}. OPE literature is generally divided into three class of methods, inverse propensity scoring~\cite{Precup00eligibilitytraces, DBLP:journals/corr/abs-1806-01347}, direct methods ~\cite{pmlr-v97-le19a, DBLP:journals/corr/abs-1802-03493}, and hybrid methods ~\cite{pmlr-v48-thomasa16,DBLP:journals/corr/JiangL15}. 

We are inspired by a subset of direct methods that aims to minimize a convex risk criteria to estimate the average state-action density ratios between the behavior and the target policy based on the data from behavior policy and the functional form of target policy returns~\cite{DBLP:conf/nips/NachumCD019}. The density ratio is later used to re-weight the reward function to estimate returns using data available from the behavior policy. Similar to the direct set of methods, our work takes advantage of a convex risk minimization algorithm to estimate density ratios between the data but unlike these set of methods we use convex risk minimization to estimate the transition ratio between two different MDP environments. 

\section{Preliminaries} \label{sec:prelim}

In this section, we introduce our notation for MDPs and reinforcement learning and present the preliminaries for importance sampling and the associated estimators.

\subsection{MDPs and Reinforcement Learning}

RL problems are modelled as a Markov Decision Process (MDP), which are represented as a 5-tuple $\mathcal{M} = \langle\mathcal{S}, \mathcal{A}, \mathcal{P}, r, S_0 \rangle$, where $\mathcal{S}$ is the state space, $\mathcal{A}$ is the set of available actions to the agent, $\mathcal{P}$ gives us the underlying transition probability of the next state given the current state and action taken in the environment, and $S_0$ is the distribution from which initial states $s_0$ are sampled. 
We define a stationary policy $\pi: \mathcal{S} \rightarrow \Delta(\mathcal{A})$ as a probability distribution over actions for each states. 
At every time step $t$, the environment evolves to state $s_{t+1} \sim \mathcal{P}(\cdot| s_t, a_t)$ from state $s_{t}$ by taking an action $a_t$ according to a policy $\pi(\cdot| s_t)$. We define value of a policy by the function $J_{\mathcal{P}}(\pi)$ which calculates the average cumulative discounted sum of rewards $r: \mathcal{S} \times \mathcal{A} \rightarrow \mathbb{R}$ over an infinite time horizon:
\begin{align}
\begin{split}
    & J_{\mathcal{P}}(\pi) = \mathbb{E}\left[\sum_{t=0}^\infty \gamma^t r(s_t, a_t)\right] \\
     & s_0 \sim S_0 \quad \forall t, a_t \sim \pi(.|s_t), s_{t+1} \sim \mathcal{P}(.|s_t, a_t) 
     \end{split}
\end{align}
where $\gamma \in (0, 1]$ is the discount factor. 
An alternative formulation to calculate $J_{\mathcal{P}}(\pi) = \frac{1}{1 - \gamma}\mathbb{E}_{(s, a) \sim d^{\pi}}[r(s, a)]$ involves calculating expectation under the distribution $d^{\pi}_{\mathcal{P}}(s, a)$, which is the normalized discounted stationary distribution: 
\begin{align}\label{eq: stationarydist}
\begin{split}
    d_{\mathcal{P}}^\pi(s, a) &= (1-\gamma) \sum_{t=0}^\infty \gamma^t \mathcal{P}r(s_t = s, a_t = a) \\
    &  s_0 \sim S_0, \forall t, a_t \sim \pi(\cdot|s_t), s_{t+1} \sim \mathcal{P}(\cdot|s_t, a_t). 
    \end{split}
\end{align}
We similarly define the normalized discounted state-action-state distribution. We observe that the relationship between $d_{\mathcal{P}}(s, a)$  and $d_{\mathcal{P}}(s, a, s')$ can be expressed by:
\begin{align}\label{eq: statedistribution}
    d_{\mathcal{P}}^\pi(s, a, s') = d_{\mathcal{P}}^\pi(s, a) \mathcal{P}(s'|s, a).
\end{align} 
It is worth noting that this relationship becomes necessary to estimate the mismatch between two MDPs. To estimate this mismatch, we first define two MDP environments: training $MDP_{tr}$ and testing $MDP_{te}$. \newline 
\begin{align}\label{eq: MDPS}
\begin{split}
    &MDP_{tr} = \left \langle \mathcal{S}, \mathcal{A}, \mathcal{P}_{tr}, r, S_0 \right\rangle, \\
    &MDP_{te} = \left \langle \mathcal{S}, \mathcal{A}, \mathcal{P}_{te}, r, S_0 \right \rangle
    \end{split}
\end{align}
Note that both these MDPs are identical in all elements except the environment model. 

Mathematically, we wish to estimate $J_{\mathcal{P}_{te}}(\pi)$ while having full access to $MDP_{tr}$ and limited data from $MDP_{te}$. First, we wish to estimate the ratio of transition probability $\zeta = \frac{\mathcal{P}_{te}(s' | s, a)}{\mathcal{P}_{tr}(s'|s, a)}$ and use this ratio to estimate $J_{\mathcal{P}_{te}}(\pi)$ using $MDP_{tr}$ and concepts from importance sampling.

\subsection{Importance Sampling}
Importance sampling is a technique to estimate the statistical properties of a given probability distribution, given access to samples from a similar distribution. Consider two probability distributions on the same sample space $(\Omega, \mathcal{F}, P)$ and $(\Omega, \mathcal{F}, Q)$.  We want to infer a statistical property $\mathbb{E}_{x \sim P}[G(x)]$, for any $G: \Omega \rightarrow \mathbb{R}$, with only access to samples from $x \sim Q$. Under the assumption that $\forall x \in \Omega$, if $P(x) > 0 \Rightarrow Q(x) > 0$, we can re-write the expectation in the following form:
\begin{align}
\begin{split}
    E_{x \sim P}[G(x)] &= \sum_{x \in \Omega} P(x) G(x) \\
    &= E_{x \sim Q}\left[\frac{P(x)}{Q(x)} G(x)\right] \label{eq: impsamp}
    \end{split}
\end{align}
The term $\frac{P(x)}{Q(x)}$ is often referred to as the importance sampling factor $\beta(x)$. When using inverse propensity scoring methods for off policy evaluation methods, this ratio $\beta(x)$ is generally known in advance. However, in some applications, this ratio needs to be estimated based on samples from the distribution $x \sim P$ and $x' \sim Q$. In the next section, we discuss about one such method to estimate the ratio $\beta(x)$. 

\subsection{Estimating the Importance Sampling Factor} \label{sec: is}
In our approach, we wish to learn $\beta(x)$ as a function of the state space. For that purpose, we look at the entropy regularisation class of methods which re-formulates the f-divergence between these two distributions in terms of their convex dual~\cite{DBLP:journals/tit/NguyenWJ10}. Specifically, the Kullback-Leiber divergence between two distributions $P$ and $Q$ can be rewritten as a maximisation problem over a space of functions $\mathcal{G}$:
\begin{align}\label{eq: expoptimizer}
    \hat{\mathcal{D}}_{KL}(P || Q) =  \max_{g \in \mathcal{G}}\sum_{x' \in \Omega}& P(x') \log(g(x'))  \\
    &- \sum_{x \in \Omega} Q(x) g(x)  +  1 \nonumber.
\end{align}
Note that the optima of the above problem is $g^* \in \mathcal{G} = \frac{P(x)}{Q(x)}$, as desired. In practice, we only have access to $n$ samples each from $P$ and $Q$. Thus, the optimal $g^*$ is approximated by searching for $\hat{g}_n$ that minimizes the empirical loss along with a regularizer $I(g)$ to prevent over-fitting:
\begin{align}\label{eq: optimizer}
    \hat{g}_n = \arg \min_{g \in \mathcal{G}}
    \frac{1}{n}\sum_{i=1}^n g(x_i') - \frac{1}{n}\sum_{i=1}^n \log(g(x_i)) + \frac{\lambda_n}{2} I(g)^2
\end{align}
In this work, we use neural networks as the class of functions $\mathcal{G}$ to show the effectiveness of our approach. In practice, researchers have also used other different function spaces like Reproducible Kernel Hilbert Space (RKHS)~\cite{DBLP:journals/corr/abs-1810-12429}. For small scale (finite state, finite action) experiments, the function space is often represented by matrices~\cite{DBLP:conf/nips/NachumCD019}. 

\section{Off Environment Policy Evaluation}\label{sec: OEE} 
We present an algorithm to evaluate performance of an agent in the testing MDP using a training MDP along with minimal data from the testing MDP. We first show that to approximate returns using the training MDP, one requires the ratio of transition probabilities of the two distributions $\zeta$. We further use the convex optimization method to propose an algorithm to estimate this ratio. 

\subsection{Methodology}
Consider a training and a testing MDP (equation \ref{eq: MDPS}). Given only a few empirical samples from the testing MDP $\{s_t, a_t, s_{t+1}\}^{tr}_{t=1, ... ,N}$, we wish to estimate $J_{\mathcal{P}_{te}}(\pi)$ using access to environment $MDP_{tr}$. To do this, consider the following expression for $J_{\mathcal{P}_{te}}(\pi)$:
\begin{flalign}\label{eq: optima}
\begin{split}
    J_{\mathcal{P}_{te}}(\pi) &=  \mathbb{E}\left [\sum_{t=0}^{\infty} \gamma^t r(s_t, a_t)\right] \\
    = & \mathbb{E}_{s_0 \sim S_0} [ r(s_0, a_0) ...   \\
    & \quad + \mathbb{E}_{s_1 \sim \mathcal{P}_{te}(\cdot|s_0, a_0)}[\gamma r(s_1, a_1) ...
    \\ & \quad \quad +  \mathbb{E}_{s_2 \sim \mathcal{P}_{te}(\cdot|s_1, a_1)} [\gamma^2 r(s_2, a_2)...]] ]  
    \end{split}
\end{flalign}
where $s_{t+1} \sim \mathcal{P}_{te}(\cdot|s_t, a_t)$, $\forall t =0, 1, .., \infty$. The last equation comes from the tower property of the expectation. Next, we use importance sampling equation \ref{eq: impsamp} to re-write the expectation in terms of $MDP_{tr}$:
\begin{flalign}\label{eq: value_approx}
\begin{split}
  & J_{\mathcal{P}_{te}}(\pi) =  \sum_{t=0}^\infty \mathbb{E}_{s_0 \sim S_0}\\ 
  &\prod_{k=0}^{t} \mathbb{E}_{s_{k} \sim \mathcal{P}_{tr}(\cdot | s_{k-1}, a_{k-1})} \frac{\mathcal{P}_{te}(s_k | s_{k-1}, a_{k-1})}{\mathcal{P}_{tr}(s_k | s_{k-1}, a_{k-1})}\gamma^t r(s_t, a_t)   \\
  = &\prod_{k=0}^{t} \mathbb{E}_{s_{k} \sim \mathcal{P}_{tr}(\cdot | s_{k-1}, a_{k-1})}  \zeta(s_k|s_{k-1}, a_{k-1})\gamma^t r(s_t, a_t)
  \end{split}
\end{flalign}
We thus define $\zeta$ as the ratio of conditional probabilities:
\begin{align}\label{eq: zeta}
\begin{split}
    \zeta(s_k|s_{k-1}, a_{k-1}) &= \frac{\mathcal{P}_{te}(s_k | s_{k-1}, a_{k-1})}{\mathcal{P}_{tr}(s_k | s_{k-1}, a_{k-1})} 
    \end{split}
\end{align}
Using importance sampling we thus show our methodology that can assess the policy performance in the testing environment using the training environment and the ratio of transition probabilities between two MDPs.

\subsection{Estimating $\zeta$}
We showed how a convex risk minimization algorithm can estimate the ratio of two probability distributions. Now, we use these principles to estimate the ratio of two conditional probabilities (i.e., the transition probabilities for environment models). 
Using the trajectory data collected on the testing MDP, we apply convex risk minimization algorithm previously described. We wish to calculate two types of density ratios: $\frac{d^{\mathcal{D}_{te}}(s, a, s')}{d^{\mathcal{D}_{tr}}(s, a, s')}$ and $\frac{d^{\mathcal{D}_{te}}(s, a)}{d^{\mathcal{D}_{tr}}(s, a)}$. 
We calculate $\zeta$ by calculating the ratio of these two density ratios which follows simply from the relation as mentioned in equation $\ref{eq: statedistribution}$:
\begin{align}\label{eq: ratiocond}
    \zeta(s' | s, a) = \Large \left({\frac{d^{\mathcal{D}_{te}}(s, a, s')}{d^{\mathcal{D}_{tr}}(s, a, s')}}\Large\right)\bigg/ \large \left({\frac{d^{\mathcal{D}_{te}}(s, a)}{d^{\mathcal{D}_{tr}}(s, a)}} \large \right)
\end{align}
Using the $\zeta$ estimator along with with equation \ref{eq: value_approx}, we can now approximate $J_{\mathcal{P}_{te}}(\pi)$ using just the training MDP. Note that whenever we use an estimate $\hat{\zeta}_n$ to approximate $J_{\mathcal{P}_{te}}$ we use the notation $J^{\hat{\zeta}_n}_{\mathcal{P}_{te}}$.
\section{Error Analysis}\label{sec: error}
To characterize the error of our estimator, we use existing results developed for empirical risk minimization algorithms to upper-bound the error on the KL-divergence. Using first order Taylor expansion, we later use that upper-bound the error in estimating $\zeta$. We also use this estimator error to approximate the upper-bounds in evaluating $J_{\mathcal{P}_{te}}(\pi)$ using the training environment. 
\subsection{Estimation Error}
Before analysing estimator error, we make the following assumptions.
\begin{assumption}\label{ass: bounds}
For any $g \in \mathcal{G}$, we assume that $\nu \leq g(x) \leq \mu \quad \forall x \in \Omega$. Note that $\nu < 1$ and $\mu > 1$
\end{assumption}
Note that estimating $\zeta$ (equation \ref{eq: ratiocond}) requires calculating a ratio of two entities. The lower bound assumption is to make sure we do not run into numerical issues. 
\begin{assumption}\label{ass: existence}
$\exists$ $g^* \in \mathcal{G}$ such that $g^*(x) = \frac{P(x)}{Q(x)}$
\end{assumption}
For the sake of simplicity, we will also assume that the optimizer contains the true minimizer $g^*$. Under these assumptions, we analyse the error between the empirically estimated optimizer $\hat{g}_n$ and the true estimator $g^*$.\newline 
\begin{theorem} \label{thm: estimation}
Under assumption \ref{ass: bounds}, we can upper bound the estimator error between $\hat{g}_n$ and $g^*$, with probability at-least $1- \delta$ for a sufficiently large $n$ by:
\begin{flalign}\label{eq: upperbound}
\begin{split}
&  \left \| \hat{g}_n - g^* \right\|_\infty^2 \leq Ke^{D_{\infty}(P \| Q)} \\
& \left(\sqrt{\frac{1}{n}} (\mu  + \max(\log \mu, -\log \nu)) + \sqrt{\frac{2 \log\frac{1}{\delta}}{n}}  \right)
\end{split}
\end{flalign}
For the sake of simplicity, we denote the upper-bound as a function of the function class $\mathcal{G}$ and $\delta$ as $M(\mathcal{G}, \delta)$. The detailed proofs along with additional experiments can be found in the supplementary materials\footnote{ Supplementary Materials: \href{https://pulkitkatdare.web.illinois.edu/icrasupp.pdf}{https://pulkitkatdare.web.illinois.edu/icrasupp.pdf}}\footnote{Code implementation: \href{https://github.com/pulkitkatdare/offenveval}{https://github.com/pulkitkatdare/offenveval}}. 
\end{theorem}
We observe that the error decreases at the rate of $O(n^{-\frac{1}{4}})$ and scales exponentially with the $\infty$-renyi divergence. Note that in theorem \ref{thm: estimation}, we bounded the error in estimating \ref{eq: optimizer}. Using this estimator error, we now upper-bound the error in estimating $\hat{\zeta}_n$.
\begin{lemma}\label{le: zetaerror}
Applying theorem \ref{thm: estimation}, we can bound estimation error for $\hat{\zeta}_n$ as follows:
\begin{align}
    \| \hat{\zeta}_n - \zeta^* \| &\leq \frac{\mu(1 + \nu \mu)}{\nu^2} \frac{M(\mathcal{G}, \delta)}{n^{1/4}}\\ \nonumber
\end{align}
\end{lemma}
\begin{theorem}\label{th: evalerror}
Under assumptions \ref{ass: existence} and \ref{ass: bounds} and the assumption that the reward function is bounded by a constant $\|r(s,a)\|_\infty \leq R$, we can bound the error in performance predicted by our algorithm $J_{\mathcal{P}_{te}}^{\hat{\zeta}_n}(\pi)$ for a time horizon $T$ by the following expression:\\
\begin{align}
\begin{split}
    &\|J_{\mathcal{P}_{te}}(\pi)  - J_{\mathcal{P}_{te}}^{\hat{\zeta}_n}(\pi)\|_\infty ^2 \\
    &\leq \frac{TM (\mathcal{G}, \delta)^2 R^2 \gamma }{\nu \sqrt{n}} \left(\frac{1 - (T+1)\frac{\gamma}{\nu}^T + T\frac{\gamma}{\nu}^{T+1}}{(1 - \frac{\gamma}{\nu})^2} \right) 
    \end{split}
\end{align}
\end{theorem}
with probability at-least $1-\delta$.

 We observe that although the bound above converges at the rate of $n^{-1/2}$, the bound also has an exponential dependence on the time horizon. Interestingly, the variance for inverse propensity scoring also depends exponentially on time horizon~\cite{DBLP:journals/corr/abs-1911-06854}. 
To reduce the error in practice, we employ to shorter time horizons $T \simeq 100$.
\section{Experiments and Results}\label{sec: expt}
Our experiments are guided by the following questions:
\begin{enumerate}
    \item How well does the ratio of two estimated importance sampling (equation \ref{eq: ratiocond}) factors work to estimate conditional probabilities? 
    \item Can we use our estimator to validate performance? What is the quality of this validation? 
    \item Can our estimator be applied to real world robots?
\end{enumerate}
We report the average performance by training 10 of these networks, all of which are initialized with different random seeds. Our algorithm is summarized in algorithm \ref{algo1}.

\subsection{Gridworld}\label{sec: gridworld}
\begin{figure*}
\centering
\begin{minipage}{0.23\textwidth}
\subfigure[Error vs. Sample Size]{\label{fig:gridworld_1}\scalebox{0.4}{\hspace{-1.0em}\hbox{\hspace{-1.0em}\input{GridWorld/gridworld_samplesize}}}}%
\end{minipage}%
\begin{minipage}{0.23\textwidth}
\subfigure[$\delta = 0.1$]{\label{fig:gridworld_eval_1}\scalebox{0.4}{\hspace{-1.0em}\hbox{\hspace{-1.0em}\input{GridWorld/gridworld_evaluation1_3162}}}}%
\end{minipage}
\begin{minipage}{0.23\textwidth}
\subfigure[$\delta=0.5$]{\label{fig:gridworld_eval_5}\scalebox{0.4}{\hbox{\hspace{-1.0em}\input{GridWorld/gridworld_evaluation5_3162}}}}
\end{minipage}
\begin{minipage}{0.23\textwidth}
\subfigure[$\delta=0.9$]{\label{fig:gridworld_eval_9}\scalebox{0.4}{\hbox{\hspace{-1.0em}\input{GridWorld/gridworld_evaluation9_3162}}}}
\end{minipage}
\caption{\textbf{Gridworld}. (a) We plot the absolute error between the transition ratios and the estimated transition ratios averaged across all state-action-state pairs. (b, c, d) Evaluation using our approach is depicted in green, importance sampling baseline in blue, the oracle in red and the true value.}
\label{fig: gridworld_2}
\vspace{-0.7cm}
\end{figure*}
In this experiment, we show that our algorithm is able to approximate the ratio of transition probabilites. We also show that it is also possible to validate the policy performance in the test environment using the simulator and our estimator. To that end, we define $\text{Gridworld}_\varepsilon$ as a $n \times n$ sized grid space. The agent starts with a fixed initial state at the bottom left corner of the grid and can take one of the four actions \{N, S, E, W\}. The MDP is designed such that whenever the agent takes an action, (e.g., N) they move in that direction with 1-$\varepsilon$ probability and with $\varepsilon$ probability end up in any of the remaining adjacent states. An episode ends when the agent reaches the goal state in the top-right corner of the grid. 

We use $\text{Gridworld}$ with ${\varepsilon = 0.3}$ as the training MDP and $\text{Gridworld}$ with ${\varepsilon = 0.1}$ as the testing MDP. In figure \ref{fig:gridworld_1}, we plot the mean absolute error of the transition probabilities over all the state-action-state pairs along with one standard deviation. The absolute error is plotted against the sample size which was used to train the estimator. We also compare our estimation error against three different grid sizes $\{10 \times 10, 20 \times 20, 40 \times 40\}$. We observe that the estimation error decreases with sample size, while larger grids require more samples to achieve reasonable error margins. These results highlight the dependence of the estimator on the size of the state space and the need for more samples in more complex environments.

Next, we evaluate the performance of a policy using the training environment and the estimated transition probability. Using equation \ref{eq: value_approx}, we evaluate performance for a range of policies (parameterised by $\delta$) over 200 timesteps. The functional form of these policies can be written as follows, $\pi_\delta(\cdot | s) = (1-\delta) \cdot U(\{N, S, E, W\}) + \delta \cdot \pi_E(.|s)$. $U(.)$ is the policy that picks each of the four actions randomly at every state and $\pi_E(.)$ is the expert policy. We evaluate the performance for $\delta=0.1, 0.5, 0.9$ over 200 timesteps respectively. These plots are presented in figures \ref{fig:gridworld_eval_1}, \ref{fig:gridworld_eval_5} and \ref{fig:gridworld_eval_9} respectively.

In each of these figures, we compare our performance against the true value which is the actual performance in the test environment, the oracle, and the importance sampling estimator. The oracle is the estimator that knows the true value of these transition probabilities, so the only error that arises in the evaluation is the error due to averaging over finite number of samples. We also compare our method against the importance sampling (IS) estimator. IS estimators trace their origins from Off Policy Evaluation methods where we have data from a behavior policy and we wish to evaluate performance for a given target policy. In all of these figures, we notice that the importance sampling method performs better whenever the target policy is close to behavior policy, the performance deteriorates as target policy moves away from the behavior policy. We observe that our method tends to generalize well over a the range of different policies. 

It is important to note that the best performance that our estimator can ideally get is the one approximated by the oracle. One might notice here that our estimator tends to perform better than the oracle. A possible reason could be that, while optimizing, we never constrain the neural network to learn over the structure of the ratio of two probability distributions. For example, suppose the function $\hat{\beta} \in \mathcal{G}$ is the optimum function that maximizes equation \ref{eq: optimizer}. On top of optimizing the equation, $\hat{\beta}$ also needs to satisfy the constraint 
\begin{align*}
\sum_{x} \hat{\beta}(x)q(x) = 1.    
\end{align*}
We assume that this constraint will be implicitly satisfied while optimizing, but since the algorithm empirical samples from distributions to calculate $\hat{\beta}$, it may not satisfy the above constraint.
\subsection{Cartpole Environment}
We now examine the Cartpole environment, which has a continuous state space and discrete action space. Since the Carpole environment is a highly unstable environments, good sampling strategies are very important, as the estimator tends to fails to generalize if only poor data is collected.
For the Cartpole environment, the training MDP has gravity $g = 10.0m/s^2$. We evaluate our performance on a range of test MDPs with gravity $g = \{7.5, 10.0, 12.5, 15.0\} m/s^2$.
As in the gridworld environment, we demonstrate our performance over a range of policies parameterized by $\delta$, $\pi_\delta(. | s) = \delta \cdot U(|\mathcal{A}|) + (1-\delta) \cdot \pi_E(.|s)$. The expert policy is a deterministic policy trained using the cross entropy method~\cite{DBLP:conf/icml/MannorRG03}. We plot of our results for each of the target environments in figure \ref{fig: cartpole_oee} over 100 timesteps.  We compare our results against a model based baseline which learns the model using the same data that was used to learn the transition ratio. We also use an IS based OPE as a baseline to compare performance. 

We observe that our algorithm closely approximates the true performance over a range of policies in the target environment. We also observe that as divergence between the two environments increases, the performance predicted by our algorithm tends to deteriorate. In comparison, the importance sampling baseline works only well when the target policies near the data collection policy but the performance soon deteriorates after that. On the other hand, the MLE which learns the model of the environment using the target environment data, tends to perform worse.
\begin{figure*}
\centering
\begin{minipage}{0.23\textwidth}
\subfigure[gravity = $7.5 m/s^2$]{\label{fig:cartpole_75}\scalebox{0.40}{\hbox{\hspace{-1.0em}\input{Cartpole/cartpole_eval_trained7.5}}}}
\end{minipage}%
\begin{minipage}{0.23\textwidth}
\subfigure[gravity = $10.0 m/s^2$]{\label{fig:cartpole_100}\scalebox{0.40}{\hbox{\hspace{-1.0em}\input{Cartpole/cartpole_eval_trained10.0}}}}
\end{minipage}%
\begin{minipage}{0.23\textwidth}
\subfigure[gravity = $12.5 m/s^2$]{\label{fig:cartpole_125}\scalebox{0.40}{\hbox{\hspace{-1.0em}\input{Cartpole/cartpole_eval_trained12.5}}}}
\end{minipage}%
\begin{minipage}{0.23\textwidth}
\subfigure[gravity = $15.0 m/s^2$]{\label{fig:cartpole_125}\scalebox{0.40}{\hbox{\hspace{-1.0em}\input{Cartpole/cartpole_eval_trained15.0}}}}
\end{minipage}
\caption{\textbf{Cartpole Environment}. Evaluation using our algorithm is depicted in green, MLE baseline is in red, the importance sampling baseline is in blue, and the true performance is in black.}
\label{fig: cartpole_oee}
\end{figure*}
\subsection{Reacher Environment}
\begin{figure}
\begin{flushleft}
\begin{minipage}{0.23\textwidth}
\subfigure[$l_{te} = 2 l_{tr} $]{\label{fig:reacher_2}\scalebox{0.4}{\hbox{\hspace{-1.0em}\input{Reacher/reacher_eval2}}}}%
\end{minipage}
\begin{minipage}{0.23\textwidth}
\subfigure[$l_{te} = 0.5 l_{tr} $]{\label{fig:reacher_5}\scalebox{0.4}{\hbox{\hspace{-1.0em}\input{Reacher/reacher_eval3}}}}
\end{minipage}
\caption{\textbf{Reacher Environment}. Evaluation using our algorithm is depicted in green, the simulator performance in red, the MLE baseline is in blue, and the true value is in black.}
\vspace{-0.7cm}
\label{fig: oee_reacher}
\end{flushleft}
\end{figure}
We also evaluate our performance on a continuous state, continous action Reacher Environment. For these experiments pertaining to the Reacher environment, the training environment has both link lengths set to $0.1 m$. We validate the performance of our algorithm on a set of two test environments. The first target environment has twice the length of the link as the training environment. The second target environment, has half the length of the link in the training environment. 

We evaluate the performance using the training environment and the OEE estimator over a range of policies parameterized by $\epsilon$, $\pi_\delta(. | s) = \delta \cdot U(|\mathcal{A}|) + (1-\delta) \cdot \pi_E(.|s)$. These performance is reported along with their standard error in figure \ref{fig: oee_reacher}. For comparison we also reported performances predicted by the simulator, and the one predicted by a learned model that trained over the same data collected over the target environment. We observe that even when the difference in the performance between the training and the testing environment is large, our method is able to make up for that difference and predict performance accurately over a range of policies.  
\subsection{Real World Kinova Environment}
\begin{table}[!htb]
\captionsetup{singlelinecheck = false, justification=justified}
\captionof{table}{Performance as predicted on a Kinova Robot}
\normalsize
\label{tab:kinova}
\begin{center} 
\begin{tabular}{@{}lcc@{}}
\toprule 
Algorithm & Mean & Std\\
\midrule
True Value & -0.911 & 0.016 \\
OEE & -1.095 & 0.101 \\
Simulator & -1.181 & 0.170 \\
\bottomrule
\end{tabular}
\end{center}
\vspace{-0.5cm}
\end{table}
To demonstrate the practicality of our algorithm in a real world setting Sim2Real tasks, we evaluate the performance of our algorithm on 
Kinova Gen3 modular robotic arm. We use a  Kinova simulator that was designed for training RL algorithms~\cite{chang2020robot}. Using this simulator along with demonstration data collected on the Robot, we optimize equation \ref{eq: optimizer} calculate the ratio of transition probabilities. 

Using this estimated conditional probability we evaluate the performance of a proportional controller that maneuvers the robot arm from an initial point to a final goal state. To that end, we choose 5 different goal states and the cost function incurred at each of the goal states is averaged over 10 runs each. We compare the performance against the one predicted by our algorithm using the simulator and the estimated conditional probability over 50 trajectories. In addition, we also report the returns approximated by the simulator. All of these results are tabulated in Table \ref{tab:kinova}. 

We observe that the performance of the simulator as compared to the one observed on the robot is off by around 30\%. Using our approach we are able to decrease that discrepancy to around 20\%. It is important to read these results in context. We note performance improvement is important because it decreases the Sim2Real gap, which helps us learn policies that are optimal in the target environment 
\section{Conclusion and Remarks} \label{sec: conclude}
In this paper, we propose an algorithm that can evaluate performance of a deployed agent using a minimal amount of data in the target environment along with access to the training environment. We first estimate the transition probability ratio between two MDPs using convex risk minimization algorithm. We demonstrate that this ratio can be used to measure performance of an RL policy in the target environment using the training environment and the estimator. 

Our experiments demonstrate that, using data from the real world which combined with the simulator, one can indeed approximate performance of any policy. We believe that the ability to validate performance of agents in another environment can further help us learn policies in the training environment that are well performing in the target domain. However, we note that the quality of data collected over target environment does impact performance. It was observed that data collected using a sub-optimal policy tends to generalize well. While we demonstrated our method on standard simulators plus real world robot, exploring the utility of such methods on more complex systems will be explored in future work. 


\section{Supplementary Material}
\subsection{Theorem 5.1}
The Kullback-Leiber divergence between two distributions $P$ and $Q$ can be rewritten as a maximisation problem over a space of functions $\mathcal{G}$ ~\cite{DBLP:journals/tit/NguyenWJ10}:
\begin{align}\label{eq: expoptimizer}
\begin{split}
    &\hat{\mathcal{D}}_{KL}(P || Q) =\\
    -&\left (\min_{g \in \mathcal{G}} \sum_{x \in \Omega} Q(x) g(x) 
    -\sum_{x' \in \Omega}P(x') \log(g(x')) -  1\right)  \\
    &= \min_{g \in \mathcal{G}} L(g)
    \end{split}
\end{align}
Note that the optima of the above problem is $g^*= \frac{P(x)}{Q(x)}$, as desired. In practice we only have access to $n$ samples each from $P$ and $Q$, the optimal $g^*$ is approximated by searching for $\hat{g}_n$ that minimizes the empirical loss along with a regularizer $I(g)$ to prevent over-fitting:
\begin{align}\label{eq: optimizer}
    \hat{g}_n = \arg \min_{g \in \mathcal{G}}
    \frac{1}{n}\sum_{i=1}^n g(x_i') - \frac{1}{n}\sum_{i=1}^n \log(g(x_i)) + \frac{\lambda_n}{2} I(g)^2
\end{align}
For implementation purposes we make sure that the space of functions $\mathcal{G}$ satisfies the following assumptions. 
\begin{assumption}\label{ass: bounds}
For any $g \in \mathcal{G}$, we assume that $\nu \leq g(x) \leq \mu \quad \forall x \in \Omega$. Note that $\nu < 1$ and $\mu > 1$
\end{assumption}
\begin{assumption}\label{ass: existence}
$\exists$ $g^* \in \mathcal{G}$ such that $g^*(x) = \frac{P(x)}{Q(x)}$
\end{assumption}
For the sake of simplicity we are also going to assume that $\mathcal{G}$ contains the true minimizer $g^*$. 
Under these assumptions we can analyze the error between the empirically estimated optimizer $\hat{g}_n$ and the true optimizer $g^*$ as follows, 

\textbf{Theorem 5.1} Under assumption \ref{ass: bounds} and \ref{ass: existence}, we can upper bound the estimator error $\hat{g}_n$ and $g^*$, with probability at-least $1- \delta$ and for a sufficiently large $n$ by,
\begin{align}\label{eq: upperbound}
\begin{split}
&\| \hat{g}_n - g^* \|_\infty^2 \leq \\&
Ke^{D_{\infty}(P \| Q)}(\sqrt{\frac{1}{n}}  (\mu  + \max(\log(\mu), -\log(\nu)) \\ &+ \sqrt{2 \log(\frac{1}{\delta})}) \\
\end{split}
\end{align}
For the sake of simplicity we abbreviate the upper-bound as as a function of the function class $\mathcal{G}$ and $\delta$ as $M(\mathcal{G}, \delta)$. 
\subsection{Proof of Theorem 5.1}\label{sec: error}
Since, equation \ref{eq: optimizer} is an empirical form of equation \ref{eq: expoptimizer} it is an Empirical Risk Minimization algorithm. We denote the original loss with respect to a function $g \in \mathcal{G}$ as $L(g)$.  Using familiar result from learning theory (Corollary 6.1 ~\cite{DBLP:journals/corr/Hajek}) with probability at-least $1-\delta$
\begin{align}
    L(\hat{g}_n) - L(g^*) \leq 4 \mathbb{E} \mathcal{R}_n(\mathcal{G}(X^n, Y^n)) + \sqrt{\frac{2 \log(\frac{1}{\delta})}{n}} \label{eq: slt}
\end{align}
Where, $\mathcal{R}_n(\mathcal{G}(X^n, Y^n))$ is the Radamacher average such that for $\forall \varepsilon_i$ is binomial distribution distributed uniformly in $\{-1, 1\}$. 
\begin{align}
    \mathcal{R}_n(\mathcal{G}(X^n, Y^n)) &= \frac{1}{n} \mathbb{E}_{\varepsilon^n} [\sup_{g \in \mathcal{G}} |\sum_{i} \varepsilon_i (g(x_i) - \log(g(y_i)|] 
\end{align}
Using the fact that the points $x_i \sim Q$ and $y_i \sim P$ are independent from each other. Using Jensen's inequality for a convex function $\phi(x)$, $\phi(\mathbb{E}[x]) \leq \mathbb{E}[\phi(x)])$. Using the fact that $\varepsilon_i$'s are independent of each other. 
\begin{align}
\begin{split}
    &\mathcal{R}_n(\mathcal{G}(X^n, Y^n)) = \frac{1}{n} \mathbb{E}_{\varepsilon^n} [\sup_{g \in \mathcal{G}} |\sum_{i} \varepsilon_i (g(x_i) - \log(g(y_i)|] \\
    &\leq \frac{1}{n} \mathbb{E}_{\varepsilon^n} \sup_{g \in \mathcal{G}} \sum_i |\varepsilon_i g(x_i) | +  \sup_{g \in \mathcal{G}} \sum_i |\varepsilon_i \log(g(y_i)) | \\
    &\leq \frac{1}{n} (\sqrt{\sum_i g(x_i)^2 } + \sqrt{\sum_i log(g(y_i))^2 })\\
    &\leq \sqrt{\frac{1}{n}} (\mu + max(\log(\mu), -\log(\nu)) \label{eq:error_bound}
    \end{split}
\end{align}
The last set of inequality comes from assuming that $ g(x) \in [\nu, \mu]$ $\forall g \in \mathcal{G}, \forall x \in \Omega$. And the fact that $\nu \in (0, 1)$ and  $\mu \geq 1$. Using equation \ref{eq:error_bound} in equation \ref{eq: slt} we get, 
\begin{align}\label{eq: radaapprox}
\begin{split}
    &L(\hat{g}_n) - L(g^*) \leq \\
    &4 \sqrt{\frac{1}{n}} \left(\mu + \max(\log(\mu), -\log(\nu)) + \sqrt{ 2\log(\frac{1}{\delta})}\right)
    \end{split}
\end{align}
With probability at least $1-\delta$. Now, let's turn our attention to the left hand side. Before we end up doing that let's define the estimation error $\Bar{e}_n(x) = \hat{g}_n(x) - g^*(x)$. Thus, we can re-write the left hand side in terms of $\Bar{e}$
\begin{align}
\begin{split}
    L(\hat{g}_n) - L(g^*)  &= L(g^* + \Bar{e}_n) - L(g^*)\\
    &= \sum_{x \in \Omega} Q(x) \bar{e}_n -\sum_{x \in \Omega} P(x) log(\frac{\Bar{e}_n + g^*}{g^*}) \\
    &= \sum_{x \in \Omega} Q(x) (\bar{e}_n - g^* \log(1 + \frac{\bar{e}_n}{g^*})) 
    \end{split}
\end{align}
Assuming that $n$ is sufficiently large such that $|\frac{\Bar{e}_n}{g^*}| \leq 1$. We can now use second order Taylor approximation for $\log (1 + x)$ for $|x| < 1$
\begin{align}\label{eq:taylorapprox}
\begin{split}
     L(\hat{g}_n) - L(g^*)  &= \sum_{x \in \Omega} Q(x) (\Bar{e}_n - g^*( \frac{\Bar{e}_n}{g^*} - \frac{1}{2} (\frac{\Bar{e}_n}{g^*})^2)) \\
     &= \sum_{x \in \Omega} Q(x)\frac{1}{2}(\frac{\Bar{e}_n^2}{g^*})
     \end{split}
\end{align}
Combining equations \ref{eq: radaapprox} and the above equation, we can bound the error with probability at least $1-\delta$ that, 

\begin{align}
\begin{split}
\sum_{x \in \Omega} &Q(x)\frac{1}{2}(\frac{\Bar{e}_n^2}{g^*}) \leq\\
&4 \sqrt{\frac{1}{n}} \left 
(\mu + \max(\log(\mu), -\log(\nu)) + \sqrt{2 \log(\frac{1}{\delta})} \right)
\end{split}
\end{align}
Under assumption \ref{ass: bounds} and \ref{ass: existence} $\exists x^* \in \Omega$ such that $|\bar{e}_n(x^*)| =  \| \hat{g}_n - g^*\|_\infty$. Thus, the equation above can be re-written as 
\begin{align}
\begin{split}
    &\frac{1}{K }\| \Bar{e}_n \|_\infty^2 \leq  \\
    &8 g^*\sqrt{\frac{1}{n}} \left (\mu + \max(\log(\mu), -\log(\nu)) + \mu \sqrt{2 \log(\frac{1}{\delta})}\right) \\
    &8 e^{D_{\infty(P \| Q)}}\sqrt{\frac{1}{n}} \left(\mu + \max(\log(\mu), -\log(\nu)) + \mu \sqrt{2 \log(\frac{1}{\delta})}\right)
    \end{split}
\end{align}
where $Q(x^*) = \frac{1}{K}$. The last inequality comes from the fact that $g^*(x^*) = \frac{P(x)}{Q(x)} \leq \sup_{x \in \Omega} \frac{P(x)}{Q(x)} = e^{D_{\infty} (P \| Q)}$. 
This completes the proof. 
\section{Estimating $\zeta$}
Using this theorem, we now bound the error in estimation of $\zeta$
\begin{align}\label{eq: zeta}
    \zeta(s' | s, a) = \Large \left({\frac{d^{\mathcal{D}_{te}}(s, a, s')}{d^{\mathcal{D}_{tr}}(s, a, s')}}\Large\right)\bigg/ \large \left({\frac{d^{\mathcal{D}_{te}}(s, a)}{d^{\mathcal{D}_{tr}}(s, a)}} \large \right)
\end{align}
Note that $\zeta$ estimation requires taking a ratio of two such aforementioned estimations. Since, each of these estimations are error prone, it also adds to error of estimating $\zeta$. We first take note of the following lemma from numerical analysis~\cite{numericalanalysis}
\begin{lemma}\label{le: product}
    Consider the following the function $f(x_1, x_2, ..., x_n) = \Pi_{i=1}^n x_i$. Now, instead of actual values of $x_i$'s suppose we have noisy estimates $x_i' = x_i + \varepsilon_i$. Using familiar results from numerical analysis, we can bound the first order error between the noisy estimate and the actual value as follows, 
\begin{align}
\begin{split}
    &| f(x_1', x_2', x_3', ..., x_n') - f(x_1, x_2, ..., x_n)| \\ &\leq {f(x_1, x_2, ..., x_n)} \sum_{i=1}^n \frac{\varepsilon_i}{x_i} \\
    & x_i \neq 0 \quad\forall i \in \{1,2,3,..., n\} 
    \end{split}
\end{align}
\end{lemma}
\begin{proof}
Let us first calculate what $| f(x_1', x_2', x_3', ..., x_n') - f(x_1, x_2, ..., x_n)|$ looks like
\begin{align}
\begin{split}
    &| f(x_1', x_2', x_3', ..., x_n') - f(x_1, x_2, ..., x_n)| \\&= 
    \sum_{i=1}^n \left ( \prod_{j \neq i } x_j \right)  \varepsilon_i  + O(\varepsilon^2)\\
    &\leq\sum_{i=1}^n \left ( \prod_{j } x_j \right)  \frac{\varepsilon_i}{x_i} \\
        &\leq \left ( \prod_{j } x_j \right)  \sum_{i=1}^n  \frac{\varepsilon_i}{x_i}  = f(x_1, x_2, ..., x_n) \sum_{i=1}^n  \frac{\varepsilon_i}{x_i}
    \end{split}
\end{align}
This completes the proof. 
\end{proof}
Using the above lemma along with theorem \textit{5.1} we now propose the error bound for estimating $\zeta$. 
\begin{lemma}\label{le: zetaerror}
Applying theorem 5.1 and lemma \ref{le: product} we can bound estimation error for $\hat{\zeta}_n$ as follows, 
\begin{align}
    \| \hat{\zeta}_n - \zeta^* \| &\leq \frac{\mu(1 + \nu \mu)}{\nu^2} \frac{M(\mathcal{G}, \delta)}{n^{1/4}}
\end{align}
\end{lemma}
\begin{proof}
Assuming that the $\zeta$ estimation can be written in functional form as follows, 
\begin{align}
\begin{split}
    &f({\frac{d_{\mathcal{P}_{te}}^\pi(s, a, s')}{d_{\mathcal{P}_{tr}}^\pi(s, a, s')}},{\frac{d_{\mathcal{P}_{tr}}^\pi(s, a)}{d_{\mathcal{P}_{te}}^\pi(s, a)}}) \\
    &= \left({\frac{d_{\mathcal{P}_{te}}^\pi(s, a, s')}{d_{\mathcal{P}_{tr}}^\pi(s, a, s')}}\Large\right) \large \left({\frac{d_{\mathcal{P}_{tr}}^\pi(s, a)}{d_{\mathcal{P}_{te}}^\pi(s, a)}} \large \right)
    \end{split}
\end{align}
We can now apply lemma \ref{le: product} along with theorem 5.1 which completes the proof. 
\end{proof}
Using the bounds on the estimation of $\hat{\zeta}_n$. We now prove the error bound in on the evaluation of returns $ {J_{\mathcal{P}_{te}}(\pi)}$
\section{Theorem 5.2}
In this section, we wish to bound the error in evaluation using the estimator that we calculated using equation \ref{eq: zeta}.  Using this calculated $\hat{\zeta}_n$, we evaluate returns (equation \ref{eq:jzetan}). To analyse the error we again re-define these two terms. 
\begin{flalign}
   J_{\mathcal{P}_{te}}(\pi)  &=  \mathbb{E}[\sum_{t=0}^{T} \gamma^t \prod_{k=1}^t \zeta(s_{k-1}, a_{k-1}, s_k) r(s_t, a_t) ]\\
   J_{\mathcal{P}_{te}}^{\hat{\zeta}_n}(\pi)  &=  \mathbb{E}[\sum_{t=0}^{T} \gamma^t \prod_{k=1}^t \hat{\zeta}_n(s_{k-1}, a_{k-1}, s_k) r(s_t, a_t) ] \label{eq:jzetan}\\
    & a_t = \pi(s_t) \quad s_{t+1} \sim \mathcal{P}_{tr}(.|s_t, a_t) \quad s_0 \sim S_0 \nonumber
\end{flalign}
\textbf{Theorem 5.2}
Under assumptions \ref{ass: existence} and \ref{ass: bounds} and the assumption that the reward function is bounded by a constant $\|r(s,a)\|_\infty \leq R$, we can bound $|J_{\mathcal{P}_{te}}(\pi)  - J_{\mathcal{P}_{te}}^{\hat{\zeta}_n}(\pi)|^2$ with atleast $1-\delta$ probability: 
\begin{align}
\begin{split}
    &|J_{\mathcal{P}_{te}}(\pi)  - J_{\mathcal{P}_{te}}^{\hat{\zeta}_n}(\pi)|^2 \\
    &\leq \frac{TM (\mathcal{G}, \delta)^2 R^2 \gamma }{\nu \sqrt{n}} \left(\frac{1 - (T+1)\frac{\gamma}{\nu}^T + T\frac{\gamma}{\nu}^{T+1}}{(1 - \frac{\gamma}{\nu})^2} \right)
    \end{split}
\end{align}

\begin{proof}
In order to estimate error, we first use the triangular inequality to divide the left hand side into terms dependent on a particular time instant.
\begin{align}
\begin{split}
     & |J_{\mathcal{P}_{te}}(\pi)  -  J_{\mathcal{P}_{te}}^{\hat{\zeta}_n}(\pi)| ^2 \\
    &=  |\mathbb{E}[\sum_{t=1}^T \gamma^t \Pi_{k=1}^{t-1}  (\zeta^*(s_{k-1}, a_{k-1}, s_{k}) \\
    & - \hat{\zeta}_n(s_{k-1}, a_{k-1}, s_{k}))  r(s_t, a_t) ]\|_\infty^2 \\
    &\leq T (\sum_{t=1}^T \mathbb{E}[\gamma^t \Pi_{k=1}^{t-1} (\zeta^*(s_{k-1}, a_{k-1}, s_{k}) \\
    & - \hat{\zeta}_n(s_{k-1}, a_{k-1}, s_{k}))  r(s_t, a_t)])^2
    \end{split}
\end{align}
where, the last inequality follows from applying the triangular inequality. We complete the proof by using lemma \ref{le: product} along with lemma \ref{le: zetaerror} to upper-bound each of these terms as follows, 
\begin{align}
\begin{split}
   &|J_{\mathcal{P}_{te}}(\pi)  - J_{\mathcal{P}_{te}}^{\hat{\zeta}_n}(\pi)|_\infty ^2 \\
   &\leq T \sum_{t=1}^T ((\gamma \frac{\mu}{\nu})^t \frac{\mu(1 + \nu \mu)}{\nu^2} \frac{M(\mathcal{G}, \delta)}{n^{1/4}}tR)^2 \\
    &\leq (\frac{\mu^2(1 + \nu \mu)^2}{\nu^4})\frac{TM(\mathcal{G}, \delta)^2 R^2}{\sqrt{n}} \sum_{t=1}^T (\frac{\gamma \mu}{\nu})^t t  \\
    &\leq  (\frac{\mu^2(1 + \nu \mu)^2}{\nu^4})\frac{TM(\mathcal{G}, \delta)^2 R^2}{\sqrt{n}} \times \\
&    (\frac{1 - T(\frac{\gamma \mu}{\nu})^T}{1 - \frac{\gamma \mu }{\nu}} + \frac{\frac{\gamma \mu }{\nu}(1 - (\frac{\gamma \mu}{\nu})^{T-1}) }{(1 - \frac{\gamma \mu}{\nu})^2}) \nonumber \\
    \end{split}
\end{align}
where the last inequality comes from applying the sum of Arithmetic Geometric Progression (AGP) series which completes the proof.  
\end{proof}
\section{Implementation Details}
For all the experiments, we keep our class of classifiers as a 3-layered neural network with relu based activation function. The final output has a tanh layer, which is scaled to $[\nu, \mu]$ (Assumption \ref{ass: bounds}). The network was trained using data from both the training and testing MDP using a stochastic gradient descent algorithm with $10^{-5}$ learning rate. 

Using data from the testing MDP $\mathcal{D}^{te}$ along with the training MDP $\mathcal{D}^{tr}$, we train two network using the equation \ref{eq: optimizer}. One, where we sample the state-action pairs $\{(s, a)\}_{te} \sim \mathcal{D}^{te},\{(s, a)\}_{tr }\sim \mathcal{D}^{tr}$ to calculate the state-action density ratio between testing and the training MDP. Two, where we sample state-action-state tuples $\{(s, a, s')\}_{te} \sim \mathcal{D}^{te}, \{(s, a, s')\}_{tr} \sim \mathcal{D}^{tr}$ to calculate state-action-state density ratio between two MDPs. We recall the definition of state-action density and state-action-state density with respect to a policy $\pi$ below
\begin{align}\label{eq: stationarydist}
\begin{split}
    &d_{\mathcal{P}}^\pi(s, a) = (1-\gamma) \sum_{t=0}^\infty \gamma^t \mathcal{P}r(s_t = s, a_t = a) \\
    &d_{\mathcal{P}}^\pi(s, a, s') = (1-\gamma) \sum_{t=0}^\infty \gamma^t \mathcal{P}r(s_t = s, a_t = a, s_{t+1} = s') \\
   & s_0 \sim S_0 \quad a_t \sim \pi(.|s_t) \quad s_{t+1} \sim \mathcal{P}(.|s_t, a_t)
    \end{split}
\end{align}

We further exploit the relationship $
    d_{\mathcal{P}}^\pi(s, a, s') = d_{\mathcal{P}}^\pi(s, a) \mathcal{P}(s'|s, a)$  between these two trained networks to calculate $\zeta$ (equation \ref{eq: zeta})

For the purpose of experiments, we assumed that one of the simulation parameters as the training MDP and use different testing MDPs over which to test our approach. For example, for the cartpole set of experiments we keep gravity value $g = 10.0 m/s^2$ as the training MDP, the results are evaluated on a range of different testing MDPs which are parameterised by gravity $g$. 

Specifically for the cartpole set of experiments, we also add a small amount of gaussian noise with standard deviation $10^{-3}$ to keep things in line with the assumption \ref{ass: bounds} and \ref{ass: existence}. 
\subsection{Estimating $\beta$}\label{sec: beta}
\begin{figure*}[!ht]
\centering
\begin{minipage}{0.23 \textwidth}
\subfigure[$P \sim \mathcal{N}(2, 1)$  ]{\label{fig:betahat1}\scalebox{0.4}{\hbox{\input{Beta/beta_2142}}}}
\end{minipage}%
\begin{minipage}{0.23\textwidth}
\subfigure[$P \sim \mathcal{N}(3, 1)$ ]{\label{fig:betahat2}\scalebox{0.4}{\hbox{\input{Beta/beta_3142}}}}
\end{minipage}%
\begin{minipage}{0.23\textwidth}
\subfigure[$P \sim \mathcal{N}(4, 1)$]{\label{fig:betahat3}\scalebox{0.4}{\hbox{\input{Beta/beta_4142}}}}
\end{minipage}
\caption{Estimating $\beta = \frac{P(x)}{Q(x)}$ for two normal distributions with $Q \sim \mathcal{N}(4, 2)$}
\end{figure*}
We evaluate the effectiveness of the estimator to calculate the ratio of probability distributions using sampled points from two different distributions $p$ and $q$. To that end, we minimize the loss function mentioned in \ref{eq: optimizer} where the optimizer would be the estimated ratio of probability distributions. We report performance for three pairs of normal distributions $p, q$ over different sample sizes used to train for this ratio. Our results are presented in figures \ref{fig:betahat1}, and \ref{fig:betahat2} and \ref{fig:betahat3}. We observe that when the two distributions are similar (figure \ref{fig:betahat3}) the estimator learns the ratio much more quickly as compared to the case when the divergence is more (figure \ref{fig:betahat1}). For example, 2000 samples is good enough to give a good enough estimate of the ratio when the distributions are closer. However, as the divergence between two distributions increases, we require larger sample sizes to get a reasonable estimate of the density ratio.
\subsection{Archery Environment}
\begin{figure*}
\centering
\begin{minipage}{.23\textwidth}
\centering
\subfigure[$W_{te} \sim \mathcal{N}(2,1)$]{\label{fig:dartenv_2}\scalebox{0.4}{\hbox{\input{DartEnv/dart_evaluation_2142}}}}
\end{minipage}
\begin{minipage}{.23\textwidth}
\centering
  \subfigure[$W_{te} \sim \mathcal{N}(3,1)$]{\label{fig:dartenv_3}\scalebox{0.40}{\hbox{\input{DartEnv/dart_evaluation_3142}}}}
\end{minipage}%
\begin{minipage}{.23\textwidth}
\centering
\subfigure[$W_{te} \sim \mathcal{N}(4,1)$]{\label{fig:dartenv_4}\scalebox{0.4}{\hbox{\input{DartEnv/dart_evaluation_4142}}}}
\end{minipage}%
\begin{minipage}{.23\textwidth}
\centering
\subfigure[Archery \newline  Schematic]{\label{fig:dartenv_schematic}\scalebox{0.4}{\input{DartEnv/dartschematic}}}
\vspace{-7pt}
\end{minipage}
\caption{Off Environment Evaluation for Archery Environments with $MDP_{tr} \sim \mathcal{N}(4,2)$}
\label{fig: dartenv}
\end{figure*}
In this section we demonstrate the ability of our algorithm to estimate the ratio of transition probabilities. We use the archery environment to demonstrate the effectiveness of our approach over a three different pair of environments. \\\\
As shown in the figure \ref{fig:dartenv_schematic} is a schematic of the archery task. The agent has to choose an angle of attack such that it hits the bulls eye which is set at the origin. The arrow faces resistance in form of horizontal wind flow with varying speeds. In these set of experiments, we parameterize the wind with different probability distributions to create different environments. Note that this system acts as a 1-step MDP, because the distribution of final position of the arrow is dependent on the angle of attack. We show that it is possible to compute the ratio of transition probabilities using samples from the two different environment, by evaluating the performance of the target environment by using the training environment and the estimated ratio (equation 4.2).\\\\
For these experiments we used a training environment where the wind speed is sampled from a normal distribution $\mathcal{N}(4, 2)$. Using this training environment we approximate the ratio of transition probabilities for three different target environments. Each of these target environments have horizontal wind speed distributed as $\mathcal{N}(2, 1), \mathcal{N}(3, 1), \mathcal{N}(4, 1)$, respectively. To demonstrate the effectiveness of this approximation, we evaluate the average performance of the agent in the target domain using the simulator and the estimated transition probability ratio. These performance results are plotted against the true performance of the same agent in these target environments. The cost function that we use is the negative of absolute distance of the arrow from the bull's eye (which is the origin). \\\\
The evaluation results for each of these environments are shown in figures \ref{fig:dartenv_2}, \ref{fig:dartenv_3} and \ref{fig:dartenv_4} respectively. We observe that our method in all the three instances was able to approximate performance with very high accuracy. We also observe that the deviations between estimated performance and the true value is larger when the data sampled is sparse. For example, in all the three figures we observe that the error is larger around $\theta = 0.5$. 











