Doubly Robust Bias Reduction in Infinite Horizon Off-Policy Estimation

Ziyang Tang, Yihao Feng, Lihong Li, Dengyong Zhou, Qiang Liu

Introduction

A key problem in reinforcement learning (RL) (Sutton & Barto 1998) is off-policy policy evaluation: given a fixed target policy of interest, estimating the average reward garnered by an agent that follows the policy, by only using data collected from different behavior policies. This problem is widely encountered in many real-life applications (Murphy et al. 2001; Li et al. 2011; Bottou et al. 2013; Thomas et al. 2017, e.g.,), where online experiments are expensive and high-quality simulators are difficult to build. It also serves as a key algorithmic component of off-policy policy optimization (Dudík et al. 2011; Jiang & Li 2016; Thomas & Brunskill 2016; Liu et al. 2019, e.g.,).

There are two major families of approaches for policy evaluation. The first approach is to build a simulator that mimics the reward and next-state transitions of the real environment (Fonteneau et al. 2013, e.g.,). While straightforward, this approach strongly relies on the model assumptions in building the simulator, which may invalidate evaluation results. The second approach is to use importance sampling to correct the sampling bias in off-policy data, so that an (almost) unbiased estimator can be obtained (Liu 2001; Strehl et al. 2010; Bottou et al. 2013). A major limitation, however, is that importance sampling can become inaccurate due to high variance. In particular, most existing IS-based estimators compute the weight as the product of the importance ratios of many steps in the trajectory, causing excessively high variance for problems with long or infinite horizon, yielding a curse of horizon (Liu et al. 2018a).

Recently, Liu et al. 2018a proposes a new estimator for infinite-horizon off-policy evaluation, which presents significant advantages to standard importance sampling methods. Their method directly estimates the density ratio between the state stationary distributions of the target and behavior policies, instead of the trajectories, thus avoiding exponential blowup of variance in the horizon. While Liu et al. 2018a’s method shows much promise by significantly reducing the variance, in practice, it may suffer from high bias due to the error or model misspecficiation when estimating the density ratio function.

In this paper, we develop a doubly robust estimator for infinite horizon off-policy estimation, by integrating Liu et al. 2018a’s method with information from an additional value function estimation. This significantly reduces the bias of Liu et al. 2018a’s method once either the density ratio, or the value function estimation is accurate (hence doubly robust). Since Liu et al. 2018a’s method already promises low variance, our additional bias reduction allows us to achieve significantly better accuracy for practical problems.

Technically, our new bias reduction method provides a new angle of double robustness for off-policy evaluation, orthogonal to the existing literature of doubly robust policy evaluation that solely devotes to variance reduction (Jiang & Li 2016; Thomas & Brunskill 2016; Farajtabar et al. 2018), mostly based on the idea of control variates (Asmussen & Glynn 2007, e.g.). Our double robustness for bias reduction is significantly different, and instead yields an intriguing connection with the fundamental primal-dual relations between the density (ratio) functions and value functions (Bertsekas 2000; Puterman 2014, e.g.,). Our new perspective may allow us to motivate new algorithms for more efficient policy evaluation, and develop unified frameworks for combining these two types of double robustness in future work.

Background

Let M=⟨S,A,r,T,μ0⟩M=\langle\mathcal{S},\mathcal{A},r,{\boldsymbol{T}},\mu_{0}\rangle be a Markov decision process (MDP) with state space S\mathcal{S}, action space A\mathcal{A}, reward function rr, transition probability function T{\boldsymbol{T}}, and initial-state distribution μ0\mu_{0}. A policy π\pi maps states to distributions over A\mathcal{A}, with π(a∣s)\pi(a|s) being the probability of selecting aa given ss. The average discounted reward for π\pi, with a given discount γ∈(0,1)\gamma\in(0,1) For average case when γ=1\gamma=1, the definition of RπR^{\pi} is the same. However, the definition of value function is different. We will assume γ<1\gamma<1 throughout our main paper for simplicity; for average case check appendix B for more details., is defined as

where τ={st,at,rt}0≤t≤T\tau=\{s_{t},a_{t},r_{t}\}_{0\leq t\leq T} is a trajectory with states, actions, and rewards collected by following policy π\pi in the MDP, starting from s0∼μ0s_{0}\sim\mu_{0}. Given a set of nn trajectories, D={st(i),at(i),rt(i)}1≤i≤n,0≤t≤T\mathcal{D}=\{s_{t}^{(i)},a_{t}^{(i)},r_{t}^{(i)}\}_{1\leq i\leq n,0\leq t\leq T}, collected under a behavior policy π0(a∣s)\pi_{0}(a|s), the off-policy evaluation problem aims to estimate the average discounted reward RπR^{\pi} for another target policy π(a∣s)\pi(a|s).

Estimation via Value Function

where PπV(s)\mathcal{P}^{\pi}V(s) is the average of the next value function given the current state ss and policy π\pi (check appendix A.1 for details).

The value function and the expected reward RπR^{\pi} is related in the following straightforward way

where the expectation is with respect to the distribution μ0(s)\mu_{0}(s) of the initial states s0s_{0} at time tt. Therefore, given an approximation V^\widehat{V} of VπV^{\pi}, and samples D0:={s0(i)}1≤i≤n0\mathcal{D}_{0}:=\{s_{0}^{(i)}\}_{1\leq i\leq n_{0}} drawn from μ0(s)\mu_{0}(s), we can estimate RπR^{\pi} by

Note that this estimator is off-policy in nature, since it requires no samples from the target policy π\pi.

Estimation via State Density Function

Denote dπ,t(⋅)d_{\pi,t}(\cdot) as average visitation of sts_{t} in time step tt. The state density function, or the discounted average visitation, is defined as:

where (1−γ)(1-\gamma) can be viewed as the normalization factor introduced by ∑t=0∞γt\sum_{t=0}^{\infty}\gamma^{t}.

Similar to Bellman equation for value function, the state density function can also be viewed as a fixed point to the following recursive equation (Liu et al. 2018a, Lemma 3):

The operator Tπ\mathcal{T}^{\pi} is an adjoint operator of Pπ\mathcal{P}^{\pi} used in (1). See Appendix A.1 for discussion.

If the density function dπd_{\pi} is known, it provides an alternative way for estimating the expected reward RπR^{\pi}, by noting that

We can see that both density function dπd_{\pi} and value function VπV^{\pi} can be used to estimate the expected reward RπR^{\pi}. We clarify the connection in detail in Appendix A.1.

Off-Policy State Visitation Importance Sampling

Equation (4) can not be directly used for off-policy estimation, since it involves expectation under the behavior policy π\pi. Liu et al. 2018a addressed this problem by introducing a change of measures via importance sampling:

where wπ/π0(s)w_{\pi/\pi_{0}}(s) is the density ratio function of dπd_{\pi} and dπ0d_{\pi_{0}}.

Given an approximation w^\widehat{w} of wπ/π0w_{\pi/\pi_{0}}, and samples D={st(i),at(i),rt(i)}1≤i≤n,0≤t≤T\mathcal{D}=\{s_{t}^{(i)},a_{t}^{(i)},r_{t}^{(i)}\}_{1\leq i\leq n,0\leq t\leq T} collected from policy π0\pi_{0}, we can estimate RπR^{\pi} as:

where ZZ is the normalized constant of the importance weights.

Doubly Robust Estimator

Doubly robust estimator is first proposed into reinforcement learning community to solve contextual bandit problem by Dudík et al. 2011 as an estimator combining inverse propensity score (IPS) estimator and direct method (DM) estimator.

Jiang & Li 2016 introduces the idea of doubly robust estimator into off-policy evaluation in reinforcement learning. It incorporates an approximate value function as a control variate to reduce the variance of importance sampling estimator. Inspired by previous works, we propose a new doubly robust estimator based on our infinite horizon off policy estimator R^SISπ\widehat{R}^{\pi}_{\text{SIS}}.

The value-based estimator R^VALπ[V^]\widehat{R}^{\pi}_{\text{VAL}}[\widehat{V}] and density-ratio-based estimator R^SISπ[w^]\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}] are expected to be accurate when V^\widehat{V} and w^\widehat{w} are accurate, respectively. Our goal is to combine their advantages, obtaining a doubly robust estimator that is accurate once either V^\widehat{V} or w^\widehat{w} or is accurate.

To simplify the problem, it is useful to exam the limit of infinite samples, with which R^VALπ[V^]\widehat{R}^{\pi}_{\text{VAL}}[\widehat{V}] and R^SISπ[w^]\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}] converge to the following limit of expectations:

Here and throughout this work, we assume V^\widehat{V} and w^\widehat{w} are fixed pre-defined approximations, and only consider the randomness from the data D\mathcal{D}.

A first observation is that we expect to have rπ≈V^−γPπV^r^{\pi}\approx\widehat{V}-\gamma\mathcal{P}^{\pi}\widehat{V} by Bellman equation (1), when V^\hat{V} approximates the true value VπV^{\pi}. Plugging this into RSISπ[w^]R^{\pi}_{\text{SIS}}[\widehat{w}] in Equation (7), we obtain the following “bridge estimator”, which incorporates information from both w^\widehat{w} and V^\widehat{V}.

where operator Pπ\mathcal{P}^{\pi} is defined in Bellman equation (1). The corresponding empirical estimator is defined by

where Z1=∑i=1n∑t=0T−1γtw^(st(i))Z_{1}=\sum_{i=1}^{n}\sum_{t=0}^{T-1}\gamma^{t}\widehat{w}(s_{t}^{(i)}) and Z2=∑i=1n∑t=0T−1γt+1w^(st(i))βπ/π0(at(i)∣st(i))Z_{2}=\sum_{i=1}^{n}\sum_{t=0}^{T-1}\gamma^{t+1}\widehat{w}(s_{t}^{(i)})\beta_{\pi/\pi_{0}}(a_{t}^{(i)}|s_{t}^{(i)}) are self-normalized constant of important weights each empirical estimation.

However, directly estimating RπR^{\pi} using the bridge estimator R^bridgeπ[V^,w^]\widehat{R}^{\pi}_{\text{bridge}}[\widehat{V},\widehat{w}] yields a poor estimation, because it includes the errors from both w^\widehat{w} and V^\widehat{V} and is in some sense “doubly worse”. However, we can construct our “doubly robust” estimator by “canceling R^bridgeπ[V^,w^]\widehat{R}^{\pi}_{\text{bridge}}[\widehat{V},\widehat{w}] out from RSISπ[w^]R^{\pi}_{\text{SIS}}[\widehat{w}] and RVALπ[V^]R^{\pi}_{\text{VAL}}[\widehat{V}]”.

The double robustness of R^DRπ[V^,w^]\widehat{R}^{\pi}_{\text{DR}}[\widehat{V},\widehat{w}] is reflected in the following key theorem, which shows that it is accurate once either V^\widehat{V} or w^\widehat{w} is accurate.

Let RDRπ[V^,w^]:=lim⁡n0,n,T→∞R^DRπ[V^,w^]R^{\pi}_{\text{DR}}[\widehat{V},\widehat{w}]:=\lim_{n_{0},n,T\to\infty}\widehat{R}^{\pi}_{\text{DR}}[\widehat{V},\widehat{w}] be the limit of R^DRπ\widehat{R}^{\pi}_{\text{DR}} when it has infinite samples. Following the definition above, we have

where εV^\varepsilon_{\widehat{V}} and εw^\varepsilon_{\widehat{w}} are errors of V^\widehat{V} and w^\widehat{w}, respective, defined as follows

The error εw^\varepsilon_{\widehat{w}} of w^\widehat{w} is measured by the difference with the true density ratio dπ(s)/dπ0(s)d_{\pi}(s)/d_{\pi_{0}}(s), and the error εV^\varepsilon_{\widehat{V}} of V^\widehat{V} is measured using the Bellman residual.

If V^\widehat{V} is exact (V^≡Vπ\widehat{V}\equiv V^{\pi}), we have εV^≡0\varepsilon_{\widehat{V}}\equiv 0; if w^\widehat{w} is exact (w^≡dπ/dπ0\widehat{w}\equiv d_{\pi}/d_{\pi_{0}}), we have εw^≡0\varepsilon_{\widehat{w}}\equiv 0. Therefore, our estimator is consistent (i.e., lim⁡n,n0→∞R^DRπ[V^,w^]=Rπ\lim_{n,n_{0}\to\infty}\widehat{R}^{\pi}_{\text{DR}}[\widehat{V},\widehat{w}]=R^{\pi}) if either V^\widehat{V} or w^\widehat{w} are exact. In comparison, R^SISπ[w^]\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}] and R^VALπ[V^]\widehat{R}^{\pi}_{\text{VAL}}[\widehat{V}] are sensitive to the error of w^\widehat{w} and V^\widehat{V}, respectively. We have

Variance Analysis

Different from the bias reduction, the doubly robust estimator does not guarantee to the reduce the variance over R^SISπ[w^]\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}] or R^VALπ[V^]\widehat{R}^{\pi}_{\text{VAL}}[\widehat{V}]. However, as we show in the following result, the variance of R^DRπ[V^,w^]\widehat{R}^{\pi}_{\text{DR}}[\hat{V},\hat{w}] is controlled by R^SISπ[w^]\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}] and R^VALπ[V^]\widehat{R}^{\pi}_{\text{VAL}}[\widehat{V}], both of which are already relatively small by the design of both methods. In addition, our method can have significant reduction of variance over R^SISπ[w^]\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}] when V^≈V\hat{V}\approx V, which can have much larger variance than R^VALπ[V^]\widehat{R}^{\pi}_{\text{VAL}}[\widehat{V}].

Assume R^DRπ[V^,w^]\widehat{R}^{\pi}_{\text{DR}}[\widehat{V},\widehat{w}] is estimated based sample D0∼μ0\mathcal{D}_{0}\sim\mu_{0} and Dπ0∼dπ0\mathcal{D}_{\pi_{0}}\sim d_{\pi_{0}}, which we assume to be independent with each other. For simplicity, assume constant normalization is used in importance sampling (hence an unbiased estimator). We have

The theorem shows the variance of our doubly robust comes from two parts: the variance for value function estimation and a variance-reduced variant of R^SISπ\widehat{R}^{\pi}_{\text{SIS}}, when V^≈Vπ\widehat{V}\approx V^{\pi}. (13) shows that our variance is always larger than that R^VALπ[V^]\widehat{R}^{\pi}_{\text{VAL}}[\widehat{V}], however, it can have lower variance than R^SISπ[w^]\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}], relevant to practice. This is because the variance of R^VALπ[V^]\widehat{R}^{\pi}_{\text{VAL}}[\widehat{V}] can be very small if we can draw a lot of samples from μ0\mu_{0}, and R^SISπ[w^]\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}] may have larger variance if the variance of the density ratio w^(s)\widehat{w}(s) and wπ/π0(s)w_{\pi/\pi_{0}}(s) are large. Meanwhile, the variance of both R^VALπ[V^]\widehat{R}^{\pi}_{\text{VAL}}[\widehat{V}] and R^SISπ[w^]\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}], by their design, are already much smaller than typical trajectory-based importance sampling methods.

The fact that the variance in (13) is a sum of two terms is because of the assumption that samples from μ0\mu_{0} and dπ0d_{\pi_{0}} are independent. In practice they have dependency but it is possible to couple the samples from μ0\mu_{0} and dπ0d_{\pi_{0}} in a certain way to even decrease the variance. We leave this to future work.

Proposed Algorithm for Off-Policy Evaluation

Suppose we have already get V^\widehat{V}, an estimation of VπV^{\pi} and w^\widehat{w}, an estimation of wπ/π0w_{\pi/\pi_{0}}, we can directly use equation (11) to estimate RπR^{\pi}. A detail procedure is described in Algorithm 1.

Double Robustness and Lagrangian Duality

We reveal a surprising connection between our double robustness and Lagrangian duality. We show that our doubly robust estimator is equivalent to the Lagrangian function of primal dual formulation of policy evaluation. This connection is of its own interest, and may provide a foundation for deriving more new algorithms in future works.

We start with the following classical optimization formulation of policy evaluation (Puterman 2014):

where we find VV to maximize its average value, subject to an inequality constraint on the Bellman equation. It can be shown that the solution of (14) is achieved by the true value function VπV^{\pi}, hence yielding an true expected reward RπR^{\pi}.

Introducing a Lagrangian multiplier ρ≥0\rho\geq 0, we can derive the Lagrangian function L(V,ρ)L(V,\rho) of (14),

Comparing L(V,ρ)L(V,\rho) with our estimator R^DRπ[V^,w^]\widehat{R}^{\pi}_{\text{DR}}[\widehat{V},\widehat{w}] in (11), we can see that they are in fact equivalent in expectation.

I) Define wρ/π0(s)=ρ(s)dπ0(s)w_{\rho/\pi_{0}}(s)=\frac{\rho(s)}{d_{\pi_{0}}(s)}. We have

which suggests that L(V,ρ)L(V,\rho) is “doubly robust” in that it equals RπR^{\pi} if either V=VπV=V^{\pi} or ρ=dπ\rho=d_{\pi}.

II) The primal problem (14) forms a strong duality with the following dual problem,

where Tπ\mathcal{T}^{\pi} is defined in (3).

This shows that the dual problem is equivalent to constraint ρ\rho using the fixed point equation (3) and maximize the average reward given distribution ρ\rho. Since the unique fixed point of (3) is dπ(s)d_{\pi}(s), the solution of (16) naturally yields the true reward RπR^{\pi}, hence forming a zero duality gap with (14).

It is natural to intuitize the double robustness of the Lagrangian function. From (15), L(V,ρ)L(V,\rho) can be viewed as estimating the reward RπR^{\pi} using value function with a correction of Bellman residual (V−rπ−γPπV)(V-r^{\pi}-\gamma\mathcal{P}^{\pi}V). If V=VπV=V^{\pi}, the estimation equals the true reward and the correction equals zero. From the dual problem (16), L(V,ρ)L(V,\rho) can be viewed as estimating RπR^{\pi} using density function ρ\rho, corrected by the residual (ρ−(1−γ)μ0−γTπρ)(\rho-(1-\gamma)\mu_{0}-\gamma\mathcal{T}^{\pi}\rho). We again get the true reward if ρ=dπ\rho=d_{\pi}.

It turns out that we can use the primal-dual formula when γ=1\gamma=1 to obtain the double robust estimator for the average reward case. We clarify it in appendix B.

The fact that the density function dπd_{\pi} forms a dual variable of the value function VπV^{\pi} is widely known in the optimal control and reinforcement learning literature (Bertsekas 2000; Puterman 2014; de Farias & Van Roy 2003, e.g.,), and has been leveraged in various works for policy optimization. However, it does not seem to be well exploited in the literature of off-policy policy evaluation.

Related Work

The problem of off-policy value evaluation has been studied in contextual bandits (Dudík et al. 2011; Wang et al. 2017) and more general finite horizon RL settings (Fonteneau et al. 2013; Li et al. 2015; Jiang & Li 2016; Thomas & Brunskill 2016; Liu et al. 2018b; Farajtabar et al. 2018; Xie et al. 2019). However, most of the existing works are based on importance sampling (IS) to correct the mismatch between the distribution of the whole trajectories induced by the behavior and target policies, which faces the “curse of horizon” (Liu et al. 2018a) when extended to long-horizon (or infinite-horizon) problems.

Several other works (Guo et al. 2017; Hallak & Mannor 2017; Liu et al. 2018a; Gelada & Bellemare 2019; Nachum et al. 2019) have been proposed to address the high variance issue in the long-horizon problems. Liu et al. 2018a apply importance sampling on the average visitation distribution of state-action pairs, instead of the distribution of the whole trajectories, which provides a unified approach to break “the curse of horizon”. However, they require to learn a density ratio function over the whole state-action pairs, which may induce large bias. Our work incorporates the density ratio and value function estimation, which significantly reduces the induced bias of two estimators, resulting a doubly robust estimator.

Our work is also closely related to DR techniques used in finite horizon problems (Murphy et al. 2001; Dudík et al. 2011; Jiang & Li 2016; Thomas & Brunskill 2016; Farajtabar et al. 2018), which incorporate an approximate value function as control variates to IS estimators. Different from existing DR approaches, our work is related to the well known duality between the density and the value function, which reveals the relationship between density (ratio) learning (Liu et al. 2018a) and value function learning. Based on this interesting observation, we further obtain the doubly robust estimator for estimating average reward in infinite-horizon problems.

Primal-Dual Value Learning

Primal-dual optimization techniques have been widely used for off-policy value function learning and policy optimization (Liu et al. 2015; Chen & Wang 2016; Dai et al. 2017a; Dai et al. 2017b; Feng et al. 2019). Nevertheless, the duality between density and value function has not been well explored in the literature of off policy value estimation. Our work proposes a new doubly robustness technique for off-policy value estimation, which can be naturally viewed as the Lagrangian function of the primal-dual formulation of policy evaluation, providing an alternative unified view for off policy value evaluation.

Experiment

In this section, we conduct simulation experiments on different environmental settings to compare our new doubly robust estimator with existing methods. We mainly compare with infinite horizon based estimator including state importance sampling estimator (Liu et al. 2018a) and value function estimator. We do not report results on the vanilla trajectory-based importance sampling estimators because of their significant higher variance, but we do compare with the doubly robust version induced by Thomas & Brunskill 2016 (self-normalized variant of Jiang & Li 2016). In all experiments we compare with Monte Carlo and naive average as Liu et al. 2018a suggested. The ground truth for each environment is calculated by averaging Monte Carlo estimation with a very large sample size.

We follow Liu et al. 2018a’s tabular environment Taxi, which has 20002000 states and 66 actions in total. For more experimental details, please check appendix C.1.

Figure 1(a)-(c) show results of comparison for different methods as we changing the number of trajectories. We can see that the MSE performance of value function(R^VALπ\widehat{R}^{\pi}_{\text{VAL}}) and state visitation importance sampling(R^SISπ\widehat{R}^{\pi}_{\text{SIS}}) estimators are mainly impeded by their large biases, while our method has much less bias thus it can keep decreasing as sample size increase and achieves same performance as on policy estimator. Figure 1(d) shows results if we change the horizon length. Notice that here we keep the number of samples to be the same, so if we increase our horizon length we will decrease the number of trajectories in the same time. We can see that our method alongside with all infinite horizon methods will get better result as horizon length increase. Figure 1(e)-(f) indicate the “double robustness” of our method, where our method benefits from either a better VV or a better ρ\rho.

Puck-Mountain

Figure 2(a)-(c) show results of comparison for different methods as we changing the number of trajectories. Similar to taxi, we find our method has much lower bias than density ratio and value function estimation, which yields a better MSE. In Figure 2(d) the performance for all infinite horizon estimator will not degenerate as horizon increases, while finite horizon method such as finite weighted horizon doubly robust will suffer from larger variance as horizon increases.

InvertedPendulum

InvertedPendulum is a pendulum that has its center of mass above its pivot point. We use the implementation of InvertedPendulum from OpenAI gym (Brockman et al. 2016), which is a continuous control task with state space in R4\mathcal{R}^{4} and we discrete the action space to be {−1,−0.3,−0.2,0,0.2,0.3,1}\{-1,-0.3,-0.2,0,0.2,0.3,1\}. More experiment details can be found in appendix C.2.

In Figure 3(a)-(c) our method again significantly reduces the bias, which yields a better MSE comparing with value and density estimation. Figure 3(d) also shows that our method consistently outperforms all other methods as the horizon increases with a fixed total timesteps.

Conclusion

In this paper, we develop a new doubly robust estimator based on the infinite horizon density ratio and off policy value estimation. Our new proposed doubly robust estimator can be accurate as long as one of the estimators are accurate, which yields a significant advantage comparing to previous estimators. Future directions include deriving more novel optimization algorithms to learn value function and density(ratio) function by using the primal dual framework.

References

Appendix A Proof

For simplicity, we define the following two operators thorough our proofs to simplify our notations.

Given a policy π\pi and the unknown environment transition T{\boldsymbol{T}}, we define Tπ\mathcal{T}^{\pi} and Pπ\mathcal{P}^{\pi} over any function f:S→Rf:\mathcal{S}\to\mathcal{R} as

Using these operator notations, we can rewrite the above two recursive equations as:

These transition operators have the following nice adjoint property.

For two function ff and gg, if the following summation is finite, we will have

Using this property we can actually using Bellman Equations to re-derive the two different way to get RπR^{\pi}.

A.2 Proof of Theorem 3.1

Using the property of the operator, we can rewrite (1−γ)μ0(s)(1-\gamma)\mu_{0}(s) using Bellman equation as dπ−γTπdπd_{\pi}-\gamma\mathcal{T}^{\pi}d_{\pi}, thus we have

and similarly if we break rπr^{\pi} as (I−γPπ)Vπ(I-\gamma\mathcal{P}^{\pi})V^{\pi}, for RSISπ[w^]R^{\pi}_{\text{SIS}}[\widehat{w}] we have:

where dw^=dπ0w^d_{\widehat{w}}=d_{\pi_{0}}\widehat{w} for short.

Compare with Rπ=∑s(I−γPπ)Vπ(s)dπ(s)R^{\pi}=\sum_{s}\left(I-\gamma\mathcal{P}^{\pi}\right)V^{\pi}(s)d_{\pi}(s), we can see the main difference between RSISπR^{\pi}_{\text{SIS}} and RVALπR^{\pi}_{\text{VAL}} with RπR^{\pi} are they replace dπd_{\pi} and dw^d_{\widehat{w}} and VπV^{\pi} with V^\widehat{V} respectively. If we add them together and minus the connection estimator, we have we will have:

where (I−γPπ)(Vπ−V^)=(I−γPπ)Vπ−(I−γPπ)V^=rπ−(I−γPπ)V^\left(I-\gamma\mathcal{P}^{\pi}\right)(V^{\pi}-\widehat{V})=\left(I-\gamma\mathcal{P}^{\pi}\right)V^{\pi}-\left(I-\gamma\mathcal{P}^{\pi}\right)\widehat{V}=r^{\pi}-\left(I-\gamma\mathcal{P}^{\pi}\right)\widehat{V}. ∎

A.3 More discussions on the Variance in Theorem 3.2

where εV^(s)=V^(s)−rπ(s)−γPπV^(s′)\varepsilon_{\widehat{V}}(s)=\widehat{V}(s)-r^{\pi}(s)-\gamma\mathcal{P}^{\pi}\widehat{V}(s^{\prime}) is the Bellman residual, δ1(s,a)=π(a∣s)π0(a∣s)r(s,a)−rπ(s)\delta_{1}(s,a)=\frac{\pi(a|s)}{\pi_{0}(a|s)}r(s,a)-r^{\pi}(s) is the randomness for action and δ2(s,a,s′)=π(a∣s)π0(a∣s)V^(s′)−PπV^(s)\delta_{2}(s,a,s^{\prime})=\frac{\pi(a|s)}{\pi_{0}(a|s)}\widehat{V}(s^{\prime})-\mathcal{P}^{\pi}\widehat{V}(s) is the randomness for transition operator over function V^\widehat{V}. Both δ1\delta_{1} and δ2\delta_{2} is zero mean if we condition over ss.

Compared with Var[R^SISπ[w^]]\text{Var}\left[\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}]\right] we have:

R^resπ[V^,w^]\widehat{R}^{\pi}_{\text{res}}[\widehat{V},\widehat{w}] can be written as

where βπ/π0(a∣s)\beta_{\pi/\pi_{0}}(a|s) is short for π(a∣s)π0(a∣s)\frac{\pi(a|s)}{\pi_{0}(a|s)}. We can break βπ/π0(a∣s)(r+γV^(s′))−V^(s)\beta_{\pi/\pi_{0}}(a|s)(r+\gamma\widehat{V}(s^{\prime}))-\widehat{V}(s) into

where εV^=V^−rπ−PπV^\varepsilon_{\widehat{V}}=\widehat{V}-r^{\pi}-\mathcal{P}^{\pi}\widehat{V} is the Bellman residual and the if we condition over ss we have the expectations for δ1\delta_{1} and δ2\delta_{2} are 0. Also notice that if we condition over ss then εV^(s)\varepsilon_{\widehat{V}}(s) become a constant thus it is independent to δ1\delta_{1} and δ2\delta_{2}. Thus we have:

For Var[R^SISπ[w^]]\text{Var}\left[\widehat{R}^{\pi}_{\text{SIS}}[\widehat{w}]\right] we have:

From the theorem we can see that the variance of residual comes from two parts, the majority part relies on the variance of ∣εV^(s)∣|\varepsilon_{\widehat{V}}(s)| is usually much smaller than rπr^{\pi} as the majority variance of state visitation importance sampling.

A.4 Proof of Theorem 4.1

We can see that the Lagrangian L(V,ρ)L(V,\rho) is actually our doubly robust estimator RDRπ[V,wρ/π0]R^{\pi}_{\text{DR}}[V,w_{\rho/\pi_{0}}].

From the last equation we can derive our dual as:

Appendix B Doubly Robust Estimator for Average Case

We start from primal dual framework to get our doubly robust estimator similar to section 4. To estimate the average reward of a given policy π\pi, we can consider solve the following linear programming:

where ρ(s)\rho(s) is the stationary distribution of states under Pπ\mathcal{P}^{\pi}, and the objective is the average reward given π\pi.

Consider the Lagrangian of above linear programming:

From Equation (21) we can get the dual formula as:

where V(s)V(s) is the value function and vˉ\bar{v} is the average reward we want to optimize.

Notice that in average case, Vπ(s)V^{\pi}(s) can be viewed as the fixed-point solution to the following Bellman equation:

Note that this explains the constraint and only if we pick vˉ=Rπ\bar{v}=R^{\pi}, we can find a VV to guarantee the constraint vˉ+V(s)−PπV(s)−rπ(s)≥0\bar{v}+V(s)-\mathcal{P}^{\pi}V(s)-r^{\pi}(s)\geq 0 holds true.

B.2 Doubly Robust Estimator

We want to build the doubly robust estimator via the lagrangian. However, the Lagrangian consist of three term ρ,V\rho,V and vˉ\bar{v}. It is counter-intuitive if we already given an estimator of vˉ≈Rπ\bar{v}\approx R^{\pi} into our estimator.

A better way to solve this problem is to remove the constraint ∑ρ(s)=1\sum\rho(s)=1, but we divide it as an self-normalization. Then our Lagrangian becomes

which can be utilized to define the doubly robust estimator for average reward.

Given a learned value function V^(s){\hat{V}}(s) for policy π\pi and an estimated density ratio w^(s)\hat{w}(s) for wπ/π0(s)w_{\pi/\pi_{0}}(s), we define

where βπ/π0(a∣s)=π(a∣s)π0(a∣s)\beta_{\pi/\pi_{0}}(a|s)=\frac{\pi(a|s)}{\pi_{0}(a|s)}.

Similarly to Theorem 3.1 we will have our double robustness for our average doubly robust estimator:

Suppose we have infinite samples and we can get

where εV^\varepsilon_{\widehat{V}} and εw^\varepsilon_{\widehat{w}} are errors of V^\widehat{V} and w^\widehat{w}, respective, defined as follows

Similar to discounted case we have RDRπ[V^,w^]=RπR^{\pi}_{\text{DR}}[\widehat{V},\widehat{w}]=R^{\pi} if either w^\widehat{w} or V^\widehat{V} is accurate.

Appendix C Experimental Details

We use an on-policy Q-learning to get a sequence of policy π0,π1,...,π19\pi_{0},\pi_{1},...,\pi_{19} as data size increases. We pick the last policy π19\pi_{19} (almost optimum) as our target policy and π18\pi_{18} as our behavior policy to guarantee that those policies are not far away. We set our discounted factor γ=0.99\gamma=0.99.

Train V^\widehat{V} and ρ^\widehat{\rho}

Separate from testing, we use a set of independent sample to first train a value function V^\widehat{V} and a density function ρ^\widehat{\rho}. Both V^\widehat{V} and ρ^\widehat{\rho} have bias due to finite sample approximation.

For how to train V^\widehat{V} and ρ^\widehat{\rho}, we choose to use Monte Carlo method to estimate V^\widehat{V} and ρ^\widehat{\rho}. We first use the finite samples to get an estimated model T^(s′∣s,a)\widehat{T}(s^{\prime}|s,a) and rewards function r^(s,a)\widehat{r}(s,a) and d^0\widehat{d}_{0}. Then we solve the following linear equation (by iteration like power method, which is actually Monte Carlo):

Estimate RπR^{\pi} Using V^\widehat{V} and ρ^\widehat{\rho}

We put V^\widehat{V} and ρ^\widehat{\rho} into the Lagrangian as equation (15) as our doubly robust estimator. For those states we haven’t visited, we set V^(s)\widehat{V}(s) and ρ^(s)\widehat{\rho}(s) as 00 and we self-normalized the ρ^\widehat{\rho} to get a fair estimation.

C.2 Continuous States Off-Policy Evaluation

We evaluate our method on two infinite horizon environments: Puck-Mountain and InvertedPendulum.

Puck-Mountain is an environment similar to Mountain-car, except that the goal of the task is to push the puck as high as possible in a local valley whose initial position is at the bottom of the valley. If the ball reaches the top sides of the valley, it will hit a roof and changes the speed to its opposite direction with half of its original speeds. The rewards was determined by the current velocity and height of the puck.

InvertedPendulum is a pendulum that has its center of mass above its pivot point. It is unstable and without additional help will fall over. We train a near optimal policy that can make the pendulum balance for a long horizon. For both behavior and target policies, we assume they are good enough to keep the pendulum balance and will never fall down until they reach the maximum timesteps. We use the implementation from OpenAI Gym (Brockman et al. 2016) and changing the dynamic by adding some additional zero mean Gaussian noise to the transition dynamic.

Behavior and Target Policies Learning

We use the open source implementation https://github.com/openai/baselines of deep Q-learning to train a 32×3232\times 32 MLP parameterized Q-function to converge. We then use the softmax policy of learned the Q-function as the target policy π\pi, which has a default temperature τ=1\tau=1. For the behavior policy π0\pi_{0}, we set a relative large temperature which encourages exploration. We set the temperature of the behavior policy τ0=1.88\tau_{0}=1.88 for Puck-Mountain and τ0=1.50\tau_{0}=1.50 for InvertedPendulum respectively.

Training of density ratio w^​(s)\hat{w}(s) and value function V^​(s)\hat{V}(s)

We use a seperate training dataset with 200 trajectories whose horizon length is 1000 to learn the density ratio w^(s)\hat{w}(s) and the value function V^(s)\hat{V}(s). For the training of density ratio, we adapt the algorithm 2 in Liu et al. 2018a to train a neural network parameterized wθ(s)w_{\theta}(s). Instead of taking the test function f(s)f(s) into an RKHS HK\mathcal{H}_{\mathcal{K}}, we parameterize the test function f(s)=fβ(s)f(s)=f_{\beta}(s) to be a neural network with parameter β\beta, and perform minimax optimization over parametr θ\theta and β\beta. A detail description can be found in Algorithm 2.

Similarly, for the training of value function, we use primal-dual based optimization methods (Dai et al. 2017b; Feng et al. 2019) to minimize the bellman residual:

where Vϕ(s)V_{\phi}(s) is the parameterized value function and fβ(s)f_{\beta}(s) is the test function. We also perform minimax over parameter ϕ\phi and β\beta. A detail description can be found in Algorithm 3.

For the network structures, we use 32×3232\times 32 feed forward neural networks to parameterize value function VϕV_{\phi} and density ratio wθ(s)w_{\theta}(s), and we use one hidden neural network with 10 units to parameterize the test function fβ(s)f_{\beta}(s). We use Adam Optimizer for all our experiments.

Estimate RπR^{\pi} using V^\widehat{V} and w^\widehat{w}

Given data samples from the policy π0\pi_{0}, We can directly use R^DRπ\widehat{R}^{\pi}_{\text{DR}} in equation (11) to estimate RπR^{\pi}.