Journal of Machine Learning Research 21 (2020) 1-63

Submitted 10/19; Revised 7/20; Published 7/20

Double Reinforcement Learning for Efficient
Off-Policy Evaluation in Markov Decision Processes
Nathan Kallus
kallus@cornell.edu
Department of Operations Research and Information Engineering and Cornell Tech
Cornell University
New York, NY 10044, USA
Masatoshi Uehara∗

mu223@cornell.edu

Department of Computer Science and Cornell Tech
Cornell University
New York, NY 10044, USA

Editor: George Konidaris

Abstract
Off-policy evaluation (OPE) in reinforcement learning allows one to evaluate novel decision
policies without needing to conduct exploration, which is often costly or otherwise infeasible.
We consider for the first time the semiparametric efficiency limits of OPE in Markov decision
processes (MDPs), where actions, rewards, and states are memoryless. We show existing
OPE estimators may fail to be efficient in this setting. We develop a new estimator basedon cross-fold estimation of q-functions and marginalized density ratios, which we term
double reinforcement learning (DRL). We show that DRL is efficient when both components
are estimated at fourth-root rates and is also doubly robust when only one component is
consistent. We investigate these properties empirically and demonstrate the performance
benefits due to harnessing memorylessness.
Keywords: Off-policy evaluation, Markov decision processes, Semiparametric efficiency,
Double machine learning

1. Introduction
Off-policy evaluation (OPE) is the problem of estimating mean rewards of a given policy(target policy) for a sequential decision-making problem using data generated by the logof another policy (behavior policy). OPE is a key problem in reinforcement learning (RL)
(Jiang and Li, 2016; Li et al., 2015; Liu et al., 2018; Mahmood et al., 2014; Munos et al.,2016; Precup et al., 2000; Thomas and Brunskill, 2016) and it finds applications as varied
as healthcare (Murphy, 2003) and education (Mandel et al., 2014). Because data can be
scarce, it is crucial to use all available data efficiently, while at the same time using flexible,
nonparametric estimators that avoid misspecification error.
In this paper, our goal is to obtain an estimator for policy value with minimal asymptotic
mean squared error under nonparametric models for the sequential decision process andbehavior policy, that is, achieving the semiparametric efficiency bound (Bickel et al., 1998).
Toward that end, we explore the efficiency bound and efficient influence function of the∗. Corresponding authorc 2020 Nathan Kallus and Masatoshi Uehara.
License: CC-BY 4.0, see https://creativecommons.org/licenses/by/4.0/. Attribution requirements are provided
at http://jmlr.org/papers/v21/19-827.html.

Kallus and Uehara

s0

a0

r0

s1

a1

r1

s2

Figure 1: M1 : Non-Markov decision process (NMDP)

s0

a0

r0

s1

a1

r1

s2

Figure 2: M2 : Markov decision process (MDP)

EffBd(M2 )2Ω(T )=

EffBd(M1 )=EffBd(M2b )=EffBd(M1b )>EffBd(M2q )>EffBd(M1q )

>

=

>

O(T 2 )

Figure 3: Relationship between the semiparametric efficiency bounds in each model, which
lower bound achievable mean-squared error. M1 , M2 are, respectively, NMDP
and MDP with unknown behavior policy. M1b , M2b are with known behavior
policy. M1q , M2q are with parametric assumptions on q-functions. Inequalities
are generically strict (Theorem 4), and the MDP bound is generally polynomial inhorizon length T while the NMDP bound is generally exponential (see Theorem 5).

target policy value under two models: non-Markov decision processes (NMDP) and Markov
decision processes (MDP). The two models are illustrated in Figs. 1 and 2 and defined
precisely in Section 1.1. While much work has studied efficient estimation under M1 (Dudiket al., 2014; Jiang and Li, 2016; Kallus and Uehara, 2019; Thomas and Brunskill, 2016),
work on M2 has been restricted to the parametric, finite-state-finite-action case (Jiang andLi, 2016) and no globally efficient estimators have been proposed. The two models areclearly nested and indeed we obtain that the efficiency bounds are generally strictly ordered(see Fig. 3). In other words, if we correctly leverage the Markov property, we can obtain
OPE estimators that are more efficient than existing ones. This is quite important, given
the practical difficulty of evaluation in long horizons (see, e.g., Gottesman et al., 2019) andgiven that many RL problems are Markovian. In particular, our results show the NMDP
efficiency bound is generally exponential in horizon length so that estimators that target
the NMDP model necessarily suffer from the curse of horizon, a phenomenon previously
2Double Reinforcement Learning in Markov Decision Processes

identified only for specific estimators. In contrast, the MDP efficiency bound, which weachieve, is generally polynomial in horizon.
We propose the Double Reinforcement Learning (DRL) estimator, which is given byby cross-fold estimation and plug-in of the q- and density ratio functions into the efficient
influence function for each model, which we derive for the first time here. The name
DRL is inspired by the Double Machine Learning estimation procedure of Chernozhukov
et al. (2018), which we leverage, and by our simultaneous use of two learning procedures:
learning of q-functions and of density ratios. We show that DRL achieves the semiparametricefficiency bound globally even when these nuisances are estimated at slow fourth-root rates
and without restricting to Donsker or bounded entropy classes, enabling the use of flexible
machine learning method for the nuisance estimation in the spirit of Chernozhukov et al.
(2018); Zheng and van der Laan (2011). Further, we show that DRL is consistent even if
only some of the nuisances are consistently estimated, known as double robustness. To thebest of our knowledge, this is the first proposed estimator shown to be globally efficient for
OPE in MDPs.
The organization of the paper is as follows. In Section 1.1, we define the OPE problem
and our models. In Section 1.2 we summarize semiparametric inference theory and in
Section 1.3 we review the literature on OPE. In Section 2, we calculate the efficient influencefunctions and efficiency bounds in each of our models. In Section 3, we propose the DRL
estimator and prove its efficiency and robustness in each model, while also reviewing theinefficiency of other estimators. In Section 4, we discuss how to estimate q-functions in
an off-policy manner to be used in DRL as well as the efficiency bound under parametricassumptions on the q-function. In Section 5, we demonstrate the benefits of DRL empirically.
A preliminary version of this work appeared as Kallus and Uehara (2020).
1.1. Problem Setup
A (potentially) non-Markov decision process (NMDP) is given by a sequence of state andaction spaces St , At for t = 0, 1, . . . , T , an initial state distribution Ps0 (s0 ), transition
probabilities Pst (st | Hat−1 ) for t = 1, . . . , T , and emission probabilities Prt (rt | Hat )for t = 0, . . . , T , where Hat = (s0 , a0 , . . . , st , at ) is the state-action history up to at . A(non-anticipatory) policy is a sequence of action probabilities πt (at | Hst ), where Hst =(s0 , a0 , . . . , at−1 , st ) is the state-action history up to st . Together, an NMDP and a policy
define a joint distribution over trajectories H = (s0 , a0 , r0 , s1 , a1 , r1 , . . . , sT , aT , rT ), givenby the product Ps0 (s0 )π0 (a0 | Hs0 )Pr0 (r0 | Ha0 )Ps1 (s1 | Ha0 ) · · · PrT (rT | HaT ). The
dependence structure of such a distribution is illustrated in Fig. 1. We denote this distribution
by Pπ and expectations in this distribution by Eπ to highlight the dependence on π.
A (time-varying) Markov decision process (MDP) is an NMDP where transitions andemissions depend only on the recent state and action and the time index t, Pst (st | Hat−1 ) =Pst (st | st−1 , at−1 ) and Prt (rt | Hat ) = Prt (rt | st , at ), and where we restrict to policies
that depend only on the recent state, πt (at | Hst ) = πt (at | st ). MDPs have the important
property that they are memoryless: given st , the trajectory starting at st is independent of
the past trajectory, so that st fully captures the current state of the system. This imposesa stricter dependence structure, which is illustrated in Fig. 2. In particular, connections
between variables with different time indices occurs only via st .
3

Kallus and Uehara

Our ultimate goal is to estimate the average cumulative reward of a policy,
" T #
X
π
rt .
ρ = Eπt=0

The quality and value functions (q- and v-functions) are defined as the following conditional
averages of the cumulative reward to go, respectively:
" T
#
" T
#
X
X
qt (Hat ) = Eπ
rk | Hat ,vt (Hst ) = Eπ
rk | Hst = Eπ [qt (Hat ) | Hst ] .
k=t

k=t

Note that the very last expectation is taken only over at ∼ πt (at | Hst ). For MDPs, wehave qt (Hat ) = qt (st , at ) and vt (Hst ) = vt (st ) = Eπ [qt (st , at ) | st ], where again the last
expectation is taken only over at ∼ πt (at | st ). For brevity, we define the random variables
qt = qt (Hat ), vt = vt (Hst ).
The off-policy evaluation (OPE) problem is to estimate the average cumulative reward
of a given (known) target evaluation policy, π e , using n observations of trajectories D ={H(1) , . . . , H(n) } independently generated by the distribution Pπb induced by using anotherpolicy, π b , in the same decision process. This latter policy, π b , is called the behavior policy
and it may be known or unknown.
A model for the data generating process Pπ of D is given by the set of products
Ps0 (s0 )π0b (a0 | Hs0 )Pr0 (r0 | Ha0 )Ps1 (s1 | Ha0 ) · · · PrT (rT | HaT ) over some possible values foreach probability distribution in the product We let M1 denote the nonparametric model
where each distribution is unknown and free. We let M1b denote the submodel of M1
where π b is known and fixed. We let M1q denote any submodel of M1 where the functions
qt (Hat ) are restricted parametrically for t = 0, . . . , T . We let M2 , M2b , M2q denote thecorresponding models where both the decision process and the behavior policy are restrictedeto be Markovian. Since π e is given, the parameter of interest, ρπ , is a function of just thepart that specifies the decision process (initial state, transition, and emission probabilities).
To streamline notation, when no subscript is denoted, all expectations E[·] and variancesvar[·] are taken with respect to the behavior policy, π b . At the same time, all v- andq-functions are for the target policy, π e . The Lp -norm is defined as kgkp = E[|g|p ]1/p . For
any function of trajectories, we define its empirical average asP
En [f (H)] = n−1 ni=1 f (H(i) ).
We denote the density ratio at time t between the target and behavior policy byηt (Hat ) =πte (at | Hst ).
πtb (at | Hst )We denote the cumulative density ratio up to time t and the marginal density ratio at time
t by, respectively,λt (Hat ) =tYηt (Hak ),µt (st , at ) =k=0pπte (st , at ),pπb (st , at )t

4

Double Reinforcement Learning in Markov Decision Processes

where pπt (st , at ) denotes the marginal distribution of st , at under Pπ . Note that under M2 ,ηt (Hat ) = ηt (at , st ). Again, for brevity we define the variables ηt = ηt (Hat ), λt = λt (Hat ),µt = µt (st , at ).
We will often assume the following:
Assumption 1 (Sequential overlap). The density ratios ηt , µt satisfy 0 ≤ ηt ≤ C, 0 ≤ µt ≤C 0 for all t = 0, . . . , T .
Assumption 2 (Bounded rewards). The reward rt satisfies 0 ≤ rt ≤ Rmax for all t =0, . . . , T .
Assumption 1 requires that every action supported by the evaluation policy is also
supported by the behavior policy, else the evaluation policy may induce state-action combinations that we cannot possibly ever see in the data. The assumption is standard in causal
inference. Assumption 2 focuses on bounded rewards, which are common in reinforcement
learning. Both assumptions can be relaxed to Lp -norm bounds on the above variables insteadof boundedness (see Remark 9).
1.2. Summary of Semiparametric Inference
We briefly review semiparametric inference theory as it pertains to the relevance of ourresults. We provide a more complete review in Appendix B.1, while providing an accessible
casual introduction here sufficient for the reader to understand the nature of our efficiency
results. For a complete textbook presentation, we refer the reader to Bickel et al. (1998);
Kosorok (2008); Tsiatis (2006); van der Laan and Robins (2003); van der Vaart (1998).
Suppose we have a model M for the generating process of the iid data H(1) , . . . , H(n) ,that is, a (potentially nonparametric) set of possible distributions for H(i) that also contains
the true distribution F ∈ M that generated the data. Consider a (scalar) parameter ofinterest R : M → R. Given an estimator R̂ (or rather a sequence of estimators), its limiting√
law is the distribution limit of n(R̂ − R(F )), and the asymptotic mean-squared error(AMSE) is the second moment of the limiting law, which in turn lower bounds the scaled
limit of the mean-squared error (MSE), lim nE[(R̂ − R(F ))2 ], by the portmanteau lemma.
Every gradient of R(·) at F ∈ M (for paths in the model M) is an F -measurable (scalar)random variable, that is, φ(H) with H ∼ F for some function φ(·). Each such function iscalled an influence function, and the influence function φeff (·) with the smallest L2 norm is
is called the efficient influence function because
EffBd(M) = EH∼F φ2eff (H)bounds below the AMSE of any estimator that is regular with respect to the model
M.1 Regular estimators are roughly those that have risk that is invariant to vanishing
perturbations to the data generating process F (that remain inside the model M), which
is a desirable property else the estimator may be unreasonably sensitive to undetectable
1. Note that EffBd(M) depends on the estimand R(·), the model M, and the instance F in that model.
We emphasize foremost the dependence on the model to highlight the differences when we change themodel from NMDP to MDP, while our estimand is always the target policy value.

5Kallus and Ueharachanges.2 Essentially, regular estimators with respect to M are those that would work for
estimating R(F ) for any instance F ∈ M. Thus, this lower bound applies per-instance forany estimator that, in a sense, works in the model M. Note we have EffBd(M0 ) ≤ EffBd(M)whenever F ∈ M0 ⊂ M, and that these may be different even though F ∈ M0 , as the setof estimators that work in M0 is potentially larger. For example, the lower bound in the
NMDP model is still larger than (or equal to) the bound in MDP model, even if considered
at a specific instance that happens to be an MDP.
There are several further interpretations of this lower bound. By the portmanteau lemma,
the lower bound on limiting law also means that EffBd(M) lower bounds the limit of the
MSE for any regular R̂, namely
lim inf nE[(R̂ − R(F ))2 ] ≥ EffBd(M).
n→∞Moreover, standard results (e.g., van der Vaart, 1998, Thm. 25.21) establish that the lower
bound also applies to all estimators (not just regular ones) in a local minimax fashion: for√any estimator, n times the worst-case MSE in a 1/ n-sized M-neighborhood around F has
a limit infimum of at least EffBd(M). Here the ambient model M is relevant in determining
the bound as the local worst-case neighborhoods are restricted to remain inside the model.
When M is a fully parametric model the semiparametric efficiency bound is actually the sameas the Cramér-Rao bound. In fact, the semiparametric efficiency bound corresponds to thesupremum of the Cramér-Rao bounds over all regular parametric submodels F ∈ Mpara ⊂ M.
Thus, it also describes the best-achievable behavior by nonparametric estimators that workin every parametric submodel.
In these senses, EffBd(M), known as the semiparametric efficiency bound, lower bounds
the achievable MSE in estimating R on the model M. If we can find an estimator whose
limiting law has zero mean and variance EffBd(M) then it must have the smallest-possible(asymptotic) MSE, and such estimators are known as (asymptotically) efficient. Moreover,
all efficient regular estimators must satisfy
√

n

1 X
φeff (H(i) ) + op (1),n(R̂ − R(F )) = √ni=1that is, they are asymptotically linear with efficient influence function φeff . Thissuggests an
1 Pn
estimation strategy: try to approximate ψ̂(H) ≈ φeff (H)+R(F ) and use R̂ = n i=1 ψ̂(H(i) ).
Done appropriately, this can provide an efficient estimate. Therefore, deriving the efficient
influence function is important both for computing the semiparametric efficiency bound andfor coming up with good estimators.
Note that the efficiency notion of optimality is different from (non-local, non-asymptotic)
minimax optimality (e.g., as used by Wang et al., 2017 for non-sequential OPE with H = 0),
which considers the worst-case behavior against a large, fixed class of instances rather thanconsidering behavior locally at F . While efficiency provides a much finer analysis of the
precise behavior at the particular instance generating the data, it is also only asymptotic.
Nonetheless, we will also establish finite-sample guarantees for our efficient estimators, where
the leading terms will be controlled by EffBd(M).
2. See Appendix B.1 and van der Vaart, 1998, Ch. 25 for precise definitions of path derivatives and regularestimators.

6

Double Reinforcement Learning in Markov Decision Processes

1.3. Summary of Literature on OPE
OPE is a central problem in both RL and in closely related problems such as dynamic
treatment regimes (DTR; Murphy et al., 2001). While the NMDP model M1 is commonly
the one assumed in the causal inference literature in the context of marginal structuralmodel estimation (Robins, 2000; Robins et al., 2000) and DTRs (Chakraborty and Moodie,
2013; Murphy et al., 2001; Zhang et al., 2013),3 in RL one often assumes that the MDP
model M2 holds. Nonetheless, with some exceptions that we review below, OPE methods
in RL have largely not leveraged the additional independence structure of M2 to improve
estimation, and in particular, the effect of this structure on efficiency has not previously
been studied and no efficient evaluation method has been proposed.
Methods for OPE can be roughly categorized into three types. The first approach isthe direct method (DM), wherein we directly estimate the q-function and use it to directly
estimate the value of the target evaluation policy. For example, one can use model-based
estimates (Mannor et al., 2007) or estimate the q-function directly using fitted LSTDQ or
more general q-iteration (Antos et al., 2008; Lagoudakis and Parr, 2004; Le et al., 2019) (wefurther review estimation of q-functions in Section 4. Once we have an estimate q̂0 of theefirst q-function, the DM estimate is simply ρ̂πDM = En [Eπe [q̂0 (s0 , a0 ) | s0 ]], where the innerexpectation is simply over a0 ∼ π e (· | s0 ) and is thus computable as a sum or integral over a
known measure and the outer expectation is simply an average over the n observations ofs0 . Recall we define all q-functions to be with respect to π e . For DM, we can leverage thestructure of M2 by simply restricting q-functions to be Markovian. However, DM can failto be efficient even under M1 unless q-functions are parametric (and correctly specified) orextremely smooth (as shown by Hahn, 1998 but only in the T = 0 case). DM is also notrobust in that, if q-functions are inconsistently estimated, the estimate will be inconsistent.
The second approach is importance sampling (IS), which averages the data weighted by
the density ratio of the evaluation and behavior policies. Given estimates λ̂t of the cumulative
density ratios
λ̂t = λt if the behavior policy is known), the IS estimate
is simplyhP (or, letting
i
h P
i
e
T
T
π
ρ̂IS = En
λ̂
r
.
(An
alternative
buthigher-varianceISestimator
is
E
λ̂
r
n
T
t=0 t tt=0 t .)
When behavior policy is known, IS is unbiased and consistent, but its variance tends tobe large due to extreme weights. In particular, under M1 , IS with λ̂t = λt is known to
be inefficient (Hirano et al., 2003), which implies it must be inefficient under M2 as well.
P
En [λ̂t rt ]A common variant of IS is the self-normalized estimate Tt=0
(Swaminathan and
En [λ̂t ]Joachims, 2015), which trades off some bias for variance but does not make IS efficient.
The third approach is the doubly robust (DR) method, which combines DM and IS
and is given by adding the estimated q-function as a control variate (Dudik et al., 2014;
πe
Jiang
hPandLi, 2016; Scharfstein et al.,i 1999). The DR estimate has the form ρ̂DR =T
En
t=0 λ̂t (rt − q̂t ) + λ̂t−1 Eπ e [q̂t |st ] .
DR is colloquially known to be efficient under M1 but no precise result is available.
When state and action spaces are finite, the model M1 is necessarily parametric, and,3. OPE is equivalent to estimating the total treatment effect of a DTR in a causal inference setting. Although
we do not explicitly use counterfactual notation (either potential outcomes or do-calculus), if we assume
the usual sequential ignorability conditions (Ertefaie and Strawderman, 2018; Luckett et al., 2018; Murphy
et al., 2001), the estimands we consider are the same and our results immediately apply.

7

Kallus and Uehara

under this parametric model, Jiang and Li (2016) study the Cramér-Rao lower bound andobserve that an infeasible DR estimator that uses oracle nuisance values instead of estimates,q̂t = qt and λ̂t = λt , would achieve the bound. For completeness, we derive precisely themore general semiparametric efficiency bound under M1 (Theorem 1) and show that two(feasible) variants of the standard DR estimate are semiparametrically efficient, either using
sample splitting with a rate condition (Theorem 6) or without sample splitting with a
Donsker condition (Theorem 9). Jiang and Li (2016) also study parametric Cramér-Rao
lower bounds under finite action and state space in the MDP model M2 , but no efficientestimator, whether parametric or nonparametric, has been proposed. See also Remarks 2
and 5. There is a significant gap to deriving the semiparametric bound, which generalizes
these results to more general action and state spaces and nonparametric models. More
importantly, our derivation yields the efficient influence function, which provides a way toconstruct an efficient estimator under M2 .
Many variations of DR have been proposed. Thomas and Brunskill (2016) propose botha self-normalized variant of DR and a variant blending DR with DM when density ratios
are extreme. Farajtabar et al. (2018) propose to optimize the choice of q̂t to minimize avariance estimate for DR rather than use a plug-in value. Kallus and Uehara (2019) propose
a DR estimator that achieves local efficiency, has certain stability properties enjoyed by
self-normalized IS, and at the same time is guaranteed to have asymptotic MSE that is neverworse than both DR, IS, and self-normalized IS.
However, all of the aforementioned IS and DR estimators do not exploit MDP structure
and, in particular, will fail to be efficient under M2 . Recently, in the finite-state-spacesetting Xie et al. (2019) studied an IS-type estimator that exploits MDP structure by
replacing density ratios with marginalized density ratios, estimated by a recursive formula.
However, this estimator is also not efficient, even in the finite tabular setting. Remark 4 ofXie et al. (2019) points out the inefficiency of the estimator proposed therein.

2. Semiparametric Inference for Off-Policy Evaluatione

In this section, we derive the efficiency bounds and efficient influence functions for ρπ underthe models M1 , M1b , M2 , and M2b . Recall that the former two models are NMDP and
the latter two are MDP and that the second and fourth assume a known behavior policy.
2.1. Semiparametric Efficiency in Non-Markov Decision Processes
First, we consider the NMDP models M1 and M1b . We do this mostly for the sake ofcompleteness since, while the influence function we derive below for the NMDP model
appears as a central object in the structure of various previously proposed doubly-robust
OPE estimators for RL (e.g., among others, Dudik et al., 2014; Farajtabar et al., 2018; Jiang
and Li, 2016; Kallus and Uehara, 2019; Thomas and Brunskill, 2016), these do not establish
it as the efficient influence function in the NMDP model or derive the semiparametricefficiency bound, with the exception of the concurrent Bibaut et al. (2019). We note that incontrast, the influence function we derive for the MDP model in the next section appears to
be novel and leads to new, more efficient estimators.

8

Double Reinforcement Learning in Markov Decision ProcesseseTheorem 1 (Efficiency bound under M1 ). The efficient influence function of ρπ under
M1 is
T
X

πe
1
φM
+
eff (H) = −ρ(λt (rt − qt ) + λt−1 vt ) .

(1)

t=0

The semiparametric efficiency bound under M1 is
EffBd(M1 ) = var(v0 ) +

T
X



E λ2t var (rt + vt+1 | Hat ) ,

(2)

t=0

where vT +1 = 0.
Under M1b , the efficient influence function and bound are the same.
Note that we do not assume Assumptions 1 and 2 in the above. The quantity EffBd(M1 )may or may not be finite. An infinite efficiency bound would imply the impossibility of√consistent n estimation. Below in Theorem 5 we show how to bound EffBd(M1 ) under
Assumptions 1 and 2.
Remark 1. The efficient influence function and bound are both the same whether we know
the behavior policy or not. Intuitively, this happens because the estimand ρπe does not infact depend on behavior policy part of the data generating distribution, Pπb , but only on thedecision process parts (initial state, transition, and emission probabilities). This phenomenon
mirrors the situation with knowledge of the propensity score in average treatment effectestimation in causal inference noted by Hahn (1998).
Remark 2. When the action and state spaces are discrete, M1 is necessarily a parametricmodel. In this discrete-space parametric model and with rt = 0 for t ≤ T − 1, Theorem 2 ofe
Jiang and Li (2016) derives the Cramér-Rao lower bound for estimating ρπ . Because thesemiparametric efficiency bound is the same as the Cramér-Rao lower bound for parametricmodels, the bound coincides with ours in this special discrete setting. Theorem 1 and therelated result in Bibaut et al. (2019) are more general, establishing the limit on estimationin non-discrete, nonparametric settings and, moreover, establishes that the efficient influencefunction coincides with the structure of many doubly-robust OPE estimators used in RL(see references above).
1
Remark 3. The efficient influence function φM
eff has the oft-noted doubly robust structure.
Specifically,
#
" T
#
" T
h
i
X
X
e
1
ρπ + E φM
λt rt + E(−λt qt + λt−1 vt )eff (H) = E

t=0

|

t=0

{z

=ρ

}

πe

= E [v0 ] + E| {z }
|
=ρπe

9

" T
X

|

{z

}

=0

#
λt (rt − qt + vt+1 ) .

t=0

{z

=0

}

Kallus and Uehara

The first term in each line corresponds to a sequential IS estimator and direct method (DM),respectively. The second term in each line is a control variate, which remain mean zero evenif we plug in different (i.e., wrong) q- and v-functions or density ratios, respectively. In thissense, it is sufficient to estimate only one part of these for consistent OPE. We will leverage
this in Theorem 10 to achieve double robustness for DRL.
Remark 4 (Dependence on Rewards and the Completely Unrestricted Model). Ourdefinition of NMDP does not allow for transitions, rewards, and policies to depend on
past rewards, as is indeed standard in RL. Nonetheless, we can easily also consider thealternative (and larger) model where we do allow these dependences, M1r . Defining
Jaj = (s0 , a0 , r0 , · · · , rj−1 , sj , aj ), and similarly Jrj and Jsj , the joint distribution over trajectories H = (s0 , . . . , rT ) in this model is given by the product Ps0 (s0 )π0 (a0 | Js0 )Pr0 (r0 |Ja0 )Ps1 (s1 | Jr0 )π1 (a0 | Js1 ) · · · PrT (rT | JaT ), where each density is free and unrestricted.
This in fact means that the model includes any joint density on the trajectory H and iscompletely unrestricted because such a sequential factorization can always be genericallydone to any density, that is, the model M1r is the model containing all joint P
densities over
the trajectory H. Redefining q- and v-functions in this model as qt (Jat ) = E[ Tj=t rj | Jat ]P
and vt (Jst ) = E[ Tj=t rj | Jst ], respectively, we can compute the efficient influence functionand efficiency bound in this model in a similar way to Theorem 1.
eTheorem 2 (Efficiency bounder under M1r ). The efficient influence function of ρπ underM1r is
πe

−ρ

+

T
X

(λt (rt − qt ) + λt−1 vt ) .

t=0

The semiparametric efficiency bound under M1r isvar(v0 ) +T
X

E[λ2t var(rt + vt+1 |Jat )].

t=02.2. Semiparametric Efficiency in Markov Decision ProcesseseNext, we derive the efficiency bound and efficient influence function for ρπ under the models
M2 and M2b , i.e., when restricting to MDP structure. To our knowledge, not only havethese never before been derived, the influence function we derive has also not appeared inany existing OPE estimators in RL. We recall that under M2 , we have qt = qt (st , at ) and
vt = vt (st ).
e

Theorem 3 (Efficiency bound under M2 ). The efficient influence function of ρπ under
M2 is
e

π
2
φM
+
eff (H) = −ρTX(µt (rt − qt ) + µt−1 vt ) .

(3)

t=0

The semiparametric efficiency bound under M2 is
EffBd(M2 ) = var(v0 ) +

T
X



E µ2t var (rt + vt+1 | st , at ) ,

t=0

10

(4)

Double Reinforcement Learning in Markov Decision Processes

where vT +1 = 0.
Under M2b , the efficient influence function and bound are the same.
Again, we do not assume Assumptions 1 and 2 in the above. Below in Theorem 5 weshow how to bound EffBd(M2 ) under Assumptions 1 and 2.
Remark 5. Again, when the action and state spaces are discrete, M2 is necessarily aparametric model. In this discrete-space parametric model and with rt = 0 for t ≤ T − 1,
Theorem 3 of Jiang and Li (2016) derives the Cramér-Rao lower bound, which must (and does)
coincide with ours in this setting. Again, our result is more general, covering nonparametricmodels and estimators, and, importantly, derives the efficient influence function, which weewill use to construct the first globally efficient estimator for ρπ under M2 .
Remark 6. The difference between the efficient influence functions in the NMDP and MDP
M2
1
models, φM
eff and φeff , is that (a) the cumulative density ratio λt is replaced with themarginalized density ratio µt and (b) that q- and v-functions only depend on recent state
and action rather than full past trajectory. Note that the latter difference is slightly hidden
M2
1
in our notation: in φM
eff , qt refers to qt (Hat ), while in φeff , qt refers to the much simpler
qt (st , at ).
Although the efficient influence function in Theorem 3 is derived de-novo in the proof,
which is the most direct route to a rigorous derivation, we can also use the geometry of
influence functions to understand the result relative to Theorem 1. The efficient influencefunction is always given by projecting the influence function of any regular asymptotic linear
1estimator onto the tangent space (Tsiatis, 2006, Thm. 4.3). Under M2 , the function φMeff (H)
from Theorem 1 can be shown to still be an influence function of some regular asymptotic
linear estimator in M2 . Projecting it onto the tangent space in M2 , where we have imposed
the independence of past and future trajectories given intermediate state, can be seen to
exactly correspond to the above marginalization over the past trajectory, explaining this2
structure of φMeff (H).
2
Remark 7. The efficient influence function φMeff (H) also has a doubly robust structure.
Specifically,

πe

ρ

" T
#
" T
#
i
X
X
M2
+ E φeff(H) = Eµt rt + E(−µt qt + µt−1 vt )
h

t=0

|

t=0

{z

}

πe

=ρ

= E [v0 ] + E| {z }
|
=ρπe

" T
X

|

{z

}

=0

#
µt (rt − qt + vt+1 ) .

t=0

{z

=0

}

The first term on the first line corresponds to the Marginalized Importance Sampling (MIS)
estimator (Xie et al., 2019). The first term on the second line corresponds to the DM
estimator. The second term on each line corresponds to control variate terms. We willleverage this in Theorem 16 to achieve double robustness for DRL.
11

Kallus and Uehara

By comparing the efficiency bounds of Theorem 1 and Theorem 3 and using Jensen’s
inequality, we can see that the Markov assumption reduces the efficiency bound, usuallystrictly so.
Theorem 4. If Pπb ∈ M2 (i.e., the underlying distribution is an MDP), then
EffBd(M2 ) ≤ EffBd(M1 ).
Moreover, the inequality is strict if there exists t ≤ T such that both λt−1 and rt−1 + vt are
not constant given st−1 , at−1 .
Beyond being sorted, we can actually see that EffBd(MDP) is generally polynomial
in T while EffBd(NMDP) is generally exponential in T . This shows that the curse ofhorizon is inevitable in NMDP. While previously it was just shown to be a limitation
specifically of IS and DR estimators (Liu et al., 2018; Xie et al., 2019), this result shows
it must plague any estimator that targets the NMDP model and that it is insurmountablewithout leveraging additional structure that further narrows the model. That EffBd(MDP)
is generally polynomial in T shows that we can potentially overcome this by efficiently
leveraging MDP structure, which is exactly what our novel DRL estimator will do.
Theorem 5. Under Assumptions 1 and 2,
2
EffBd(MDP) ≤ C 0 Rmax
(T + 1)2 ,2
EffBd(NMDP) ≤ C T +1 Rmax
(T + 1)2 .
2 ) thenIf Eπe [log(ηt )] ≥ Cmin and Eπe [log(var(rt + vt+1 | Hat ))] ≥ log(Vmin
T +1 2
EffBd(NMDP) ≥ Cmin
Vmin .


Note that Eπe [log(ηt )] = Eπe KL(πte (· | st ) || πtb (· | st )) is the KL divergence between the
distributions over actions induced by π e and π b , averaged over the states visited by π e , andthat Eπe [log(var(rt + vt+1 | Hat ))] = E [log(var(rt + vt+1 | Hat ))].
The lower bound on EffBd(NMDP) shows that the curse of horizon is inevitable. The
condition on ηt simply means that the evaluation and behavior policies are not becoming
arbitrarily similar as t grows (on-policy evaluation does not suffer from curse of horizon).
The condition on rt + vt+1 essentially ensures that rewards are not trivially constant. The
upper bound on EffBd(MDP) shows that, in contrast, variance that is polynomial in T is
possible in the MDP model.
Remark 8 (Consistency of EffBd(NMDP) and EffBd(NMDP)). NMDP models may betrivially transformed into MDP models by letting the state variable be the whole history
Hst . Then, the trajectory becomes {Hs0 , a0 , r0 , Hs1 , a1 , r1 , · · · } and the efficiency boundunder this transformed MDP matches the efficiency bound under the original NMDP sinceµt (Hst , at ) =tYpπte (Hst , at )ηk (Hak ) = λt (Hat ).
=pπb (Hst , at )
k=0

t

12

Double Reinforcement Learning in Markov Decision Processes

3. Efficient Estimation Using Double Reinforcement Learning
In this section, we construct the DRL estimator and then study its properties in the various
models. In particular, we show that DRL is globally efficient under very mild assumptions.
In the NMDP model, these assumptions are generally weaker than needed for efficiency ofprevious estimators. In the MDP model, this provides the first globally efficiency estimatorfor OPE. We further show that DRL enjoys certain double robustness properties when somenuisances are inconsistently estimated.
DRL is a meta-estimator; it takes in as input estimators for q-functions and density
ratios and combines them in a particular manner that ensures efficiency even when the
input estimators may not be well behaved. This is achieved by following the cross-foldsample-splitting strategy developed by Chernozhukov et al. (2018); Klaassen (1987); Zheng
and van der Laan (2011). We proceed by presenting DRL and its properties in each setting(NMDP and MDP). In the NMDP setting, DRL amounts to the cross-fold version of the RL
OPE doubly robust estimator, which was proposed in the experiments of Jiang and Li (2016,
Section 6.1) but not analyzed. In the MDP setting, DRL is the first semiparametricallyefficient and doubly robust estimator.
Throughout this section we assume that Assumptions 1 and 2 hold.
3.1. Double Reinforcement Learning for NMDPs
Given a learning algorithm to estimate the q-function q(Hat ) and cumulative density ratio
function λt (Hat ), DRL with K-fold sample splitting (K ≥ 2) for NMDPs proceeds as follows:
1. Randomly permute the data indices and let Dj = {d(j − 1)n/Ke + 1, . . . , djn/Ke}
for j = 1, . . . , K. Let ji be the fold containing observation i so that i ∈ Dji (namely,
ji = 1 + b(i − 1)K/nc).
(j)

(j)

2. For j = 1, . . . , K, construct estimators λ̂t (Hat ) and q̂t (Hat ) based on the training
data given by all trajectories excluding those in Dj , that is, {1, . . . , n}\Dj .
3. Leteρ̂πDRL(M1 ) =n T 
1 X X (ji ) (i)  (i)(j )λ̂t (Hat ) rt − q̂t i (Ha(i)t )n
i=1 t=0

Z
(ji )
(ji )
(i) 0
e 0
(i)
(i)
q̂t ((Hst , at ))dπt (at | Hst ) .
+ λ̂t−1 (Hat−1 )
a0t

(i)

(i)

Here (Hst , a0t ) represents the trajectory given by appending a0t to Hst ; note this differs(i)from Hat . Note further that the integral becomes a simple sum when π e has finite
support over actions (e.g., if π e is deterministic or if there are finitely many actions).
e

π from1
In other words, we approximate the efficient influence function φMeff (H) + ρ
Theorem 1 by replacing the unknown q- and density ratio functions with estimates thereof
and we take empirical averages of this approximation over the data, where for each data
point we use q- and density ratio function estimates based only on the half-sample that does
not contain the data point. While K = 2 is sufficient to achieve efficiency, larger K allows

13Kallus and Uehara

nuisances to be fit on more data and may prove practically successful. Note also that wetake only a single random K-fold split and that is enough to achieve our results below. In
practice, repeating the above process over several splits of the data and taking the averageof resulting DRL estimates can only reduce the variance without increasing bias.
This estimator has several desirable properties. To state them, we assume the following
conditions for the estimators, reflecting Assumptions 1 and 2:
(j)

(j)

Assumption 3 (Bounded estimators). lim supn→∞ kλ̂t k∞ < ∞, lim supn→∞ kq̂t k∞ < ∞for each 0 ≤ t ≤ T , 1 ≤ j ≤ K.
Assumption 3 provides that the estimators are eventually almost surely bounded. Thisto be expected when their target estimands are bounded, and, in particular, will necessarilybe the case for many non-parametric estimators such as random forests and kernel regression.
(j)(j)
And, per Assumptions 1 and 2, we have kλt k∞ ≤ C t+1 and kqt k∞ < (T + 1 − t)Rmax .
Note, however, we need not assume the same bound applies to the estimates; any finite bound(independent of n) suffices. For the rest of this subsection we will assume Assumption 3
hold.
We first prove that DRL achieves the semiparametric efficiency bound, even if each√nuisance estimator has a slow, nonparametric convergence rate (sub- n).
(j)Theorem 6 (Efficiency of ρ̂DRL(M1 ) under M1 ). Suppose kλ̂t(j)

− λt k2 kq̂t

− qt k2 =(j)(j)
op (n−1/2 ), kλ̂t − λt k2 = op (1), kq̂t − qt k2 = op (1) for 0 ≤ t ≤ T, 1 ≤ j ≤ K.

Then,

the estimator ρ̂DRL(M1 ) achieves the semiparametric efficiency bound under M1 .
Remark 9. Assumptions 1 to 3 posited L∞ bounds on density ratios, rewards, and theirestimates. These assumptions are standard in both reinforcement learning and causal
inference. It is possible to relax these to Lp bounds at the cost of requiring stronger
convergence on nuisance estimates above, requiring L2/(1−1/p) convergence instead of the L2
convergence above. Since L2 convergence in estimation is usually the standard convergencemode considered and standard results can be invoked to ensure such rates, and similarly
L∞ bounds on density ratios and rewards are also standard, we focus our analysis on thismost common case of the assumptions in order to avoid a cumbersome presentation.
Remark 10. There are two important points to make about this result. First, we have notassumed a Donsker condition (van der Vaart, 1998) on the class of estimators λ̂t and q̂t . This
is why this type of sample splitting estimator is called a double machine learning: the only
required condition is a convergence rate condition at a nonparametric rate, allowing the useof complex machine learning estimators, for which one cannot verify the Donsker condition(Chernozhukov et al., 2018). In fact, many adaptive or high-dimensional estimators failto satisfy Donsker conditions (Dı́az, 2019). Eschewing such conditions allows us to use
such estimators as the highly adaptive LASSO (Benkeser and van der Laan, 2016), càdlàgfunction estimators in very high dimensions (Bibaut and van der Laan, 2019b), and random
forests (Wager and Walther, 2016) as long as their convergence rates are ensured undercertain conditions. Benkeser and van der Laan (2016); Bibaut and van der Laan (2019b) inparticular establish op (n−1/4 ) rates, which are compatible with our assumptions. Second,
relative to the efficient influence function, which is defined in terms of the true q-function andcumulative density ratio, there is no inflation in DRL’s asymptotic variance due to plugging
14

Double Reinforcement Learning in Markov Decision Processes

in estimated nuisance functions. This is due to the doubly robust structure of efficient
influence function so that the estimation errors multiply and drop out of the first-ordervariance terms. This is in contrast to inefficient importance sampling estimators, as we will
see in Theorem 20.
In addition to efficiency, we can also establish finite-sample guarantees for DRL, where
the leading term is controlled by the efficient variance.
Theorem 7. Suppose that for some C1 , C2 , for every n, with probability at least 1 − δ, we(j)(j)
simultaneously have that kλ̂t − λt k22 ≤ κ1 = C1 (n−2α1 + log(2KT /δ)/n), kq̂t − qt k22 ≤ κ2 =(j)
C2 (n−2α2 + log(2KT /δ)/n) ∀t ≤ T, ∀j ≤ K. Suppose moreover that 0 ≤ λ̂t ≤ C t , 0 ≤(j)q̂t ≤ (T + 1 − t)Rmax for ∀t ≤ T . Then, for every n, with probability at least 1 − 7δ, wehaver
2 log(2/δ)Effbd(M1 )πeπeρ̂DRL(M1 ) − ρ ≤ns√2log(2/δ)T 2 (T Rmax C T +1 κ1 κ2 + κ1 T 2 Rmax
+ κ2 C 2(T +1) )
+ Q1
n
+ Q2

√
log(2/δ)T Rmax C T +1
+ Q3 T κ1 κ2 ,nwhere Q1 , Q2 , Q3 are constants not depending on δ, T, Rmax , C, n, C1 , C2 .
Notice that, if α1 >q0, α2 > 0, α1 + α2 > 1/2 as in Theorem 6, then the leading term

1)in the above is exactly 2 log(2/δ)Effbd(M(with no additional constant factor), while thenother terms are of strictly smaller order in n.
The rate assumptions in Theorem 7 are standard finite-sample estimation guarantees.
For example, if estimators for nuisances are obtained by empirical risk minimization methods
based on L2 -loss, then the results of Bartlett et al. (2005) apply. In this cases, the rates α1
and α2 would be determined by the local Rademacher complexity of the posited function
classes. The number 2KT in log(2KT /δ) comes from the fact that there are 2KT nuisanceestimators. The boundedness assumptions on the estimates reflect the bounds on theestimands, per Assumptions 1 and 2.
Often in RL, the behavior policy is known and need not be estimated. That is, we can(j)let λ̂t = λ. In this case, as an immediate corollary of Theorem 6, we have a much weakercondition for semiparametric efficiency: just that we estimate the q-function consistently,
without a rate.

(j)Corollary 8 (Efficiency of ρ̂DRL(M1 ) under M1b ). Suppose λ̂t = λt and kq̂t − qt k2 = op (1)for 0 ≤ t ≤ T, 1 ≤ j ≤ K. Then, the estimator ρ̂DRL(M1 ) achieves the semiparametricefficiency bound under M1b .
Without sample splitting, we have to assume a Donsker condition for the class of
estimators in order to control a stochastic equicontinutiy term (see, e.g., van der Vaart, 1998,
Lemma 19.24). Although this is more restrictive, for completeness, we also include a theorem
establishing the semiparametric efficiency of the standard plug-in doubly robust estimator for
15

Kallus and Uehara

NMDPs (Jiang and Li, 2016) when assuming the Donsker condition for in-sample-estimatednuisance functions, since this result was never precisely established before. We do notrecommend this estimator due to its restrictive requirements on nuisance estimators.
Theorem 9 (Efficiency
hP without sample splitting). Letiλ̂t , q̂t be estimators based on DeTπ
and let ρ̂DR = En
t=0 λ̂t (rt − q̂t ) − λ̂t−1 Eπ e [q̂t | Hst ] . Suppose kλ̂t − λt k2 kq̂t − qt k2 =op (n−1/2 ), kλ̂t − λt k2 = op (1), kq̂t − qt k2 = op (1) for 0 ≤ t ≤ T and that q̂t , λ̂t belong to ae
Donsker class. Then, the estimator ρ̂πDR achieves the semiparametric efficiency bound underM1 .
Thus, in M1 , in comparison to the standard doubly robust estimator, DRL enjoys
efficiency under milder conditions. To our knowledge, Theorems 6 and 9 are the first results
precisely showing semiparametric efficiency for any OPE estimator.
In addition to efficiency, DRL enjoys a double robustness guarantee (Rotnitzky and
Vansteelandt, 2014; Rotnitzky et al., 2019). Specifically, if at least just one model is correctly√
specified, then the DRL is estimator is still n-consistent.
√(j)Theorem 10 (Double robustness ( n-consistency)). Suppose kλ̂t − λ†t k2 = Op (n−α1 ) and(j)kq̂t − qt† k2 = Op (n−α2 ) for 0 ≤ t ≤ T, 1 ≤ j ≤ K. If, for each 0 ≤ t ≤ T , either λ†t = λtand α1 ≥ 1/2, α2 > 0 or qt† = qt and α2 ≥ 1/2, α1 > 0, then the estimator ρ̂DRL(M1 ) is√en-consistent around ρπ .
(j)

In particular, if the behavior policy is known so that λ̂t = λt , we can always ensure the√(j)estimator is n-consistent (an example is the IS estimator, which has q̂t = qt† = 0).
For consistency without a rate, it is sufficient for one nuisance to be consistent without arate.
Corollary 11 (Double robustness (Consistency)). Suppose kλ̂t − λ†t k2 = op (1) and kq̂t −qt† k2 = op (1) for 0 ≤ t ≤ T, 1 ≤ j ≤ K. If, for each 0 ≤ t ≤ T , either λ†t = λt or qt† = qt ,ethen the estimator ρ̂DRL(M1 ) is consistent around ρπ .
(j)

(j)

A remaining question is when can we get nonparametric estimators achieving the necessary
rates for the q- and density ratio functions. We discuss estimating q-functions in Section 4.
Regarding the density ratio, λk , if the behavior policy is known then the density ratio is
known and if the behavior policy is unknown it must be estimated. Any estimator satisfying
the slow rate conditions would suffice. Q
(j)
b,(j)
b,(j)
For example, we may let λ̂k = kt=0 πte /π̂t , where π̂t
is some nonparametric
(j)
b,(j)
regression estimator. Then λ̂k would enjoy the same rates as π̂t :b,(j)Lemma 12. Suppose π̂tand πtb are uniformly bounded by some constant below and that
b,(j)
(j)
b
−α
kπ̂t
− πt k2 = op (n ). Then, kλ̂t − λt k2 = op (n−α ).
Rates for π̂ b,(j) can be obtained from standard results for nonparametric regression, suchas for kernel and sieve estimators (Newey and Mcfadden, 1994; Stone, 1994) or nonparametricestimators suited for high dimensions and non-smooth models (Bibaut and van der Laan,
2019a; Imaizumi and Fukumizu, 2018; Khosravi et al., 2019).
16

Double Reinforcement Learning in Markov Decision Processes

Alternatively, parametric models can be used for qt and (if behavior policy is unknown)
λt . Then, under standard regularity conditions, using MLE and other parametric regression(j)estimators for behavior policy would yield kλ̂t − λ†t k2 = Op (n−1/2 ), where λ†t = λt if themodel is well-specified. Similarly, in Section 4, we discuss how using parametric q-models(j)yields kq̂t − qt† k2 = Op (n−1/2 ). If both models are correctly specified then Theorem 6
immediately implies DRL achieves the efficiency bound. When using parametric models,
this is sometimes termed local efficiency (i.e., local to the specific parametric model). If only√
one model is correctly specified then Theorem 10 ensures the estimator is still n-consistent.
3.2. Double Reinforcement Learning for MDPs
Leveraging our derivation of the efficient influence function for M2 in Theorem 3, we cansimilarly construct our DRL estimator for MDPs. Given a learning algorithm to estimate
the q-function qt (st , at ) and marginal density ratio function µt (st , at ), DRL with K-foldsample splitting (K ≥ 2) for MDPs proceeds as follows:
1. Randomly permute the data indices and let Dj = {d(j − 1)n/Ke + 1, . . . , djn/Ke}
for j = 1, . . . , K. Let ji be the fold containing observation i so that i ∈ Dji (namely,
ji = 1 + b(i − 1)K/nc).
(j)

(j)

2. For j = 1, . . . , K, construct estimators µ̂t (st , at ) and q̂t (st , at ) based on the training
data given by all trajectories excluding those in Dj , that is, {1, . . . , n}\Dj .
3. Leteρ̂πDRL(M2 ) =n T 
1 X X (ji ) (i) (i)  (i)(j )µ̂t (st , at ) rt − q̂t i (st (i) , at (i) )n
i=1 t=0

Z
(ji )
(ji )
(i)
(i)
(i) 0
e 0
(i)
+ µ̂t−1 (st−1 , at−1 )q̂t (st , at )dπt (at | st ) .
a0t

(i)

Again, note the difference between the dummy a0t and the data at , and that the
integral becomes a simple sum when π e has finite support over actions (e.g., if π e isdeterministic or if there are finitely many actions).
eπ2Again, what we have done is approximate the efficient influence function φMeff (H) + ρfrom Theorem 3 and taken its empirical average over the data, where for each data point
we use q- and marginal density ratio function estimates based only on the half-sample that
does not contain the data point. Again, taking one split suffices for our results, but one canrepeat the above over many splits and take averages without deterioration.
Since both estimators are approximating their respective influence function as we deriveeein Theorems 1 and 3, the differences between ρ̂πDRL(M1 ) and ρ̂πDRL(M2 ) , as noted in Remark 6,
is (a) λt is replaced with µt and (b) q- and v-functions only depend on recent state and
e
(j)
(j)
action rather than full past trajectory. Again notice that in ρ̂πDRL(M1 ) , q̂t refers to q̂t (Hat ),
e

(j)

(j)

while in ρ̂πDRL(M2 ) , q̂t refers to the much simpler q̂t (st , at ).
Again, to establish the properties of DRL for MDPs, we assume the following conditions
reflecting Assumptions 1 and 2:
17

Kallus and Uehara

(j)

(j)

Assumption 4 (Bounded estimators). lim supn→∞ kµ̂t k∞ < ∞, lim supn→∞ kq̂t k∞ < ∞for each 0 ≤ t ≤ T , 1 ≤ j ≤ K.
Like Assumption 3, we expect Assumption 4 to hold under Assumptions 1 and 2 as the(j)(j)
latter provides that kµt k∞ ≤ C 0 , kqt k ≤ (T + 1 − t)Rmax . And, for many non-parametricestimators, Assumption 4 will necessarily hold. Again, we need not have that the same
bounds hold for the estimators as hold for the estimands; any finite bound will do. For the
rest of this subsection we will assume Assumption 4 hold.
The following result establishes that DRL is the first efficient OPE estimator for MDPs.
In fact, it is efficient even if each nuisance estimator has a slow, nonparametric convergence√
rate (sub- n). Moreover, as before, we make no restrictive Donsker assumption; theonly required condition is the convergence rate condition.. This result leverages our novel
derivation of the efficient influence function in Theorem 3 and the structure of the influencefunction, which ensures no variance inflation due to estimating the nuisance functions.
(j)− µt k2 kq̂t(j)− qt k 2 =(j)(j)
op (n−1/2 ), kµ̂t − µt k2 = op (1), kq̂t − qt k2 = op (1) for 0 ≤ t ≤ T, 1 ≤ j ≤ K.

Then, theTheorem 13 (Efficiency of ρ̂DRL(M2 ) under M2 ). Suppose kµ̂testimator ρ̂DRL(M2 ) achieves the semiparametric efficiency bound under M2 .
Remark 11 (Example cases). We note a few specific cases of Theorem 13.
• Tabular case: Suppose the state and action spaces are finite: |St | , |At | < ∞. Thenboth µt and qt are parametric functions with parameters given by their values at each(st , at ) pair. They can therefore be easily estimated at Op (n−1/2 ) rates, ensuring theabove rate conditions are easily satisfied. For example, we can use simple frequency
estimators that simply take sample averages within each (st , at ) bin (Li and Racine,
2007, Chapter 3). Other examples and additional detail are given in Sections 3.3 and 4.
• Finite state space, known behavior policy: Suppose now only |St | < ∞ while At canbe continuous and that ηt is known. Then µt (st , at ) = ηt (st , at )wt (st ) and wt (st ) =E[λt−1 | st ] is a parametric function easily estimated at Op (n−1/2 ) rates using afrequency estimator or using a recursive estimator as in Xie et al. (2019) (more detailin Section 3.3). It therefore suffices for q-function estimators
PT to have errors that areop (1), i.e., only consistency is needed without a rate. Since k=t rk has finite variance(it is bounded) and qt is its regression on (st , at ), Theorems 4.2, 5.1, 6.1, 10.3, and 16.1
of Györfi et al. (2006) establish that this can be done with any of histogram estimates,
kernel regression, k-nearest neighbor, sieve regression, or neural networks of growing
width, respectively. (These provide L2 convergence of L2 errors, which is stronger than
the in-probability convergence of L2 errors we require.)• Nonparametric case: In the fully nonparametric case, our nuisance estimators may
converge more slowly than Op (n−1/2 ). Our result nonetheless accommodates such
lower rates and, crucially, does not impose strong metric entropy conditions that would
exclude flexible machine learning estimators. We discuss in greater detail how one
might estimate the nuisance functions in Sections 3.3 and 4.
As before, we can also obtain a finite-sample guarantee for DRL in M2 with a leading
constant controlled by the asymptotic variance, which in this case is efficient. If we can say
18

Double Reinforcement Learning in Markov Decision Processes

C1 , C2 do not depend on C T +1 , this gives a finite sample result with a sample complexity,
2 T 2.
which depends only C 0 not on C T +1 , noting Effbd(M2 ) is bounded by C 0 Rmax
Theorem 14. Suppose that for some C1 , C2 , for every n, with probability at least 1 − δ, we(j)(j)
simultaneously have that kµ̂t − µt k22 ≤ κ1 = C1 (n−2α1 + log(2KT /δ)/n), kq̂t − λt k22 ≤κ2 = C2 (n−2α2 + log(2KT /δ)/n) ∀t ≤ T, 1 ≤ j ≤ K. Suppose moreover that 0 ≤ µ̂t ≤C 0 , 0 ≤ q̂t ≤ (T + 1 − t)Rmax for t ≤ T . Then, for every n, with probability at least 1 − 7δ,
we haver
2 log(2/δ)Effbd(M2 )πeπeρ̂DRL(M2 ) − ρ ≤nr√2
log(2/δ)T 2 (T Rmax C 0 κ1 κ2 + κ1 T 2 Rmax
+ κ2 {C 0 }2 )+ Q1n√log(2/δ)T Rmax C 0+ Q2+ Q3 T κ1 κ2 ,nwhere Q1 , Q2 , Q3 are constants not depending on δ, T, Rmax , C 0 , n, C1 , C2 , C.
As before,
q if α1 > 0, α2 > 0, α1 +α2 > 1/2, then the leading term in n in the above bound2)
is exactly 2 log(2/δ)Effbd(M(with no constant factor). Note that, while C1 (embeddedninside κ1 ) does not appear in the leading n term, ensuring the conditions of Theorem 14 may(j)require C1 to depend on T , depending on which estimate µ̂t for µt is used; see Remark 17.
The DRL estimator in M2 can also achieve efficiency without sample splitting (i.e.,with adaptive in-sample estimation of nuisances) if we impose an Donsker condition on theestimated nuisances.

sample splitting). Let µ̂t , q̂t be estimators
based on D andTheorem 15 (Efficiency without
hP
i
e
T
π
let ρ̂DRL(M2 ), adaptive = En
t=0 (µ̂t (rt − q̂t ) − µ̂t−1 Eπ e [q̂t | Hst ]) . Suppose kµ̂t −µt k2 kq̂t −qt k2 = op (n−1/2 ), kµ̂t − µt k2 = op (1), kq̂t − qt k2 = op (1) for 0 ≤ t ≤ T and that q̂t , µ̂t belongeto a Donsker class. Then, the estimator ρ̂πDRL(M2 ), adaptive achieves the semiparametricefficiency bound under M2 .
In addition to efficiency, DRL enjoys a double robustness guarantee in M2 .
√(j)
Theorem 16 (Double robustness ( n-consistency)). Suppose kµ̂t − µ†t k2 = Op (n−α1 ) and(j)kq̂t − qt† k2 = Op (n−α2 ) for 0 ≤ t ≤ T, 1 ≤ j ≤ K. If, for each 0 ≤ t ≤ T , either µ†t = µtand α1 ≥ 1/2, α2 > 0 or qt† = qt and α2 ≥ 1/2, α1 > 0, then the estimator ρ̂DRL(M2 ) is√en-consistent around ρπ .
Again, we obtain consistency without a rate even if just one nuisance is consistent without
a rate.
Corollary 17 (Double robustness (consistency)). Suppose kµ̂t − µ†t k2 = op (1) and kq̂t −qt† k2 = op (1) for 0 ≤ t ≤ T, 1 ≤ j ≤ K. If, for each 0 ≤ t ≤ T , either µ†t = µt or qt† = qt ,ethen the estimator ρ̂DRL(M2 ) is consistent around ρπ .
(j)

19

(j)

Kallus and Uehara√
Remark 12. When the behavior policy is known, the estimator ρ̂DRL(M1 ) is still nconsistent under M2 even without smoothness conditions on µt because M2 is includedin M1 so that Theorem 10 applies. On the other hand, in the nonparametric setting, theestimator ρ̂DRL(M2 ) requires some smoothness conditions even if the behavior policy is
known because µt must still be estimated. In this sense, when the behavior policy is known,
ρ̂DRL(M1 ) is more robust than ρ̂DRL(M2 ) under M2 but its asymptotic variance is bigger,
and generally strictly so.
A reaming question is how to estimate the nuisances at the necessary rates. We discuss
q-function estimation in Section 4. For estimating µk , one can leverage the following
relationship to reduce it to a regression problem:
where wt (st ) = E[λt−1 | st ].

µt (st , at ) = ηt (st , at )wt (st ),(5)

Thus, for example, when the behavior policy is known, we need only estimate wt , whichamounts to regressing λt−1 on st . So, in particular, if wt (st ) belongs to a Hölder classwith smoothness α and st has dimension ds , estimating wt with a sieve-type estimator ŵt(j)(j)
based on the loss function (λt−1 − wt (st ))2 and letting µ̂t (st , at ) = ηt (st , at )ŵt (st ) will(j)give a convergence rate kµ̂t (st , at ) − µt (st , at )k2 = Op (n−α/(α+dst ) ) (Chen, 2007). Whenthe behavior policy is unknown, it can be first estimated to construct λ̂t and we can repeat
the above replacing λt with λ̂t . In particular, there will be no deterioration in rate if πtb
also belongs to a Hölder class with smoothness α and if we further split each Dj , estimate
πbt as in Theorem 12 on one half, and plug it in to estimate wt on the other half. Further
strategies for estimating µt are discussed in Section 3.3 below.
In the special case where we use parametric models for µt and qt , under some regularity
conditions, parametric estimators will generally satisfy kµ̂t − µ†t k2 = Op (n−1/2 ) and kq̂t −qt† k2 = Op (n−1/2 ), where qt† = qt and µ†t = µt if the models are well-specified. (See Section 4
regarding estimating the q-function). Thus, if both models are correctly specified, thenTheorem 13 yields local efficiency. If only one model is correctly specified, Theorem 16 yields
double robustness.
3.3. Estimating Marginalized Density Ratios and the Inefficiency of
Marginalized Importance Sampling
In this section we discuss strategies for estimating µt and also show that doing OPE
estimation using only marginalized density ratios, as recently proposed, leads to inefficientevaluation in M2 .
Given an estimate µ̂t of µt , the Marginalized Importance Sampling (MIS) is given by
" T
#
X
πe
ρ̂MIS = Enµ̂t rt ,(6)
t=0

which resembles the IS estimator but where we replace λ̂t with µ̂t . Note that µt = E[λt | st , at ],i.e., the marginalization of λt over Hat−1 ; hence the name MIS. Equation (6) can also beseen as DRL without sample splitting when we let q̂t = 0.
MIS is a generic meta-estimator that depends on a particular estimate µ̂t . As we will
see in this section, the asymptotic variance of MIS depends on how we estimate µ̂t ; this is20

Double Reinforcement Learning in Markov Decision Processes

unlike DRL, whose asymptotic variance is insensitive to this choice as long as it satisfied
lax rate conditions. To study MIS, we therefore study different instantiations of Eq. (6)for different µt estimators. Xie et al. (2019) provides one possible estimate in the special
finite-state case. We will next study general-case estimators µ̂t based on regression. See
also Remark 17 regarding different µt estimators. Note that since µt = ηt wt , when behavior
policy is known we can estimate µt as µ̂t = ηt ŵt . We focus on two cases: when ŵ is estimatedusing a histogram by averaging λt−1 by state over a finite state space and a nonparametricextension.
e

Theorem 18 (Asymptotic variance of ρ̂πMIS with finite state space). Suppose |St | < ∞ for0 ≤ t ≤ T . Let(i)i=1 I[st = st ]λt−1
.
Pn
(i)
I[s
=
s
]
t
t
i=1

Pn
ŵt (st ) =
e

(7)
e

Then ρ̂πMIS is consistent and asymptotically normal (CAN) around ρπ and its asymptotic
MSE is" T#Xvarµt rt + (λt−1 − wt )Eπe [rt | st ] .
(8)
t=0

For the proof of Theorem 18, we use an argument based on the theory of U -statistics(van der Vaart, 1998, Ch. 12) in order to rephrase
the MIS estimator with ŵ as in Eq. (7)Pein an asymptotically linear form: ρ̂πMIS = En [ Tt=0 µt rt + (λt−1 − wt )Eπe [rt |st ]] + op (n−1/2 ).
eThis influence function is different from the efficient influence function; therefore, ρ̂πMIS
with histogram nuisance estimators is not efficient (the efficient influence function is unique).
In fact, we can confirm this fact by calculating and comparing the variances.
Theorem 19. If Pπb ∈ M2 (i.e., the underlying distribution is an MDP), Eq. (8) is greaterthan or equal to EffBd(M2 ). The difference isvar[v0 ] +

T
−1
X



E (wt − λt−1 )2 var (ηt vt (st+1 )|st ) .

t=0

We now turn to the nonparametric case, where we first consider a sieve-type extensionof the ŵ estimator.
e

Theorem 20 (Asymptotic variance of ρ̂πMIS with nonparametric wt estimate). Suppose
E[(λt − µt )q ] < ∞ for some q > 1. Letŵt (st ) = arg min En [(wt (st ) − λt−1 )2 ],wt (st )∈Λα
ds

(9)

t

where Λαds is the space of Hölder functions with smoothness α and the dimension dst .
teeAssume wt ∈ Λαds , α/(2α + dst ) > 1/4. Then the estimator ρ̂πMIS is CAN around ρπ and itstasymptotic MSE is equal to Eq. (8).
21Kallus and Uehara

Remark 13. The estimator Eq. (9) is over an infinite-dimensional function space. It canbe replaced with a finite-dimensional approximation Λαn such that Λαn → Λαds . Following
Example 1(b) in Section 8 (Shen, 1997), it can be shown that this will lead to the same
asymptotic MSE as in Eq. (8) and not change the conclusion of Theorem 20.
Remark 14. When the action and sample space is continuous, the histogram estimator in
Eq. (7) can also easily be extended to a kernel estimator:(i)
i=1 Kh (st − st )λt−1,Pn(i)
i=1 Kh (st − st )Pnŵt (st ) =(10)where Kh is a kernel with a bandwidth h.
The smoothness condition in Theorem 20 ensures we can estimate wt at fourth-root
rates using Eq. (9). Following Newey and Mcfadden (1994) and utilizing a high-order kernel,
we can obtain similar fourth-root rates for Eq. (10) and a similar variance result for MIS.
Unlike Eq. (9), we cannot invoke a Donsker condition to prove a stochastic equicontinuity
condition. However, it is still possible to show this directly based on a V-statistics theory(see Chapter 8 of Newey and Mcfadden, 1994).
Finally, we also consider estimating µt directly and nonparametrically using the relationµt (st , at ) = E[λt | st , at ].

(11)

A sieve-type regression estimator for µt is then constructed asµ̂t (st , at ) =arg minµt (st ,at )∈Λα
ds +dat

En [(µt (st , at ) − λt )2 ].

(12)

t

e

Theorem 21 (Asymptotic variance of ρ̂πMIS with nonparametric µt estimate). Suppose
E[(λt − µt )q ] < ∞ for some q > 1. Let µ̂t be as in Eq. (12). Assume µt ∈ Λαds +da ,
t
t
hP
i
e
T
π
α/(2α + dst + dat ) > 1/4. Then the estimator En
and its
t=0 µ̂t rt is CAN around ρ
asymptotic MSE is equal to
" T
#
X
var
µt rt + (λt − µt )E[rt | st , at ] .
(13)
t=0

Remark 15. While both estimators for µt in Theorems 20 and 21 achieve fourth-root rates
under the respective conditions, the resulting asymptotic variances in Eqs. (8) and (13) are
different and generally incomparable. Both are inefficient, but which is larger is problemdependent. Note that, in contrast, the asymptotic variance of DRL (Theorem 13) is thesame (and is efficient) regardless of which way is used to estimate µt as long as we havethe necessary rate. When the behavior policy is known, using Eq. (9) may be better thanEq. (12) when estimating µt nonparametrically because the smoothness condition is weaker
and the convergence rate is faster (since dst < dat + dst ). However, when using parametricmodels, the rates are the same (under correct specification) and sometimes it is easier tomodel µt (st , at ) rather than wt (st ), as we do in Section 5.1.
22

Double Reinforcement Learning in Markov Decision Processes

Remark 16 (The Inefficiency of MIS). Theorems 18, 20 and 21 each study the MIS estimatoreρ̂πMIS in Eq. (6) but with different estimator for the nuisance µt = ηt wt . As noted above,
unlike DRL, the variance of the MIS estimator actually depends on the way this nuisance isestimated. And, in each case, the MIS estimator was inefficient. In the finite-state-spacesetting with known behavior policy, Xie et al. (2019) propose another, different MIS estimator
based on estimating the MDP transition kernel; but per their Remark 4 it is also inefficient.
(In contrast, Remark 11 shows that DRL is efficient in the finite-state-space setting without
requiring any smoothness conditions.) This does not immediately imply the MIS estimatoris always inefficient, as it may depend on how µt is estimated, but semiparametric theory
strongly suggests there is reason to believe that MIS would in general be inefficient.
One natural question that sheds light on this is how woulda hypothetical MIS estimatorP
perform with oracle values for µt . In fact, the variance of Tt=0 µt rt is in general incomparable
to EffBd(M2 ), that is, it may be smaller or larger depending on the particular instance.
This may surprising but is not contradictory since one can in fact prove that no regularestimator(let alone an efficient one) in either M2 or M2b could ever have P
the formPT−1/2
En [ t=0 µt rt ] + op (n), that is, asymptotically linear with influence function Tt=0 µt rt .
PTeThis is because t=0 µt rt is not a gradient of ρπ under either M2 or M2b (see Theorem 24).
PT
This is in stark contrast to IS: t=0 λt rt is always a valid influence function under either
M1b or M2b since we know its empirical average always gives an unbiasedP
linear estimator(not just asymptotically). Indeed, we similarly have that the variance of Tt=0 µt rt is also
PP
incomparable to Tt=0 λt rt . The function Tt=0 µt rt is an influence function (a gradient)under the MDP model with known transition kernel, but that is a very restrictive andunrealistic model.
One interesting specific case is the fully tabular setting (finite state and action spaces).
Since our paper was posted, the more recent Yin and Wang (2020) considered a “modified”
version of the estimator of Xie et al. (2019) in order to obtain efficiency under the tabularcase of M2b . By simple algebra, the estimator of Yin and Wang (2020), which is defined as
(i) R
(i)
(i)
1 Pn PTe
i=1
t=0 ŵt (st ) rt P̂rt (rt |st , at )π (at |st )d(rt , at ),
n


R
Q
e (a |s )P̂ (s |sπ,a)d(Hat−1 ),where ŵt (st ) = p̂ b1(st ) P̂st (st |st−1 , at−1 ) t−1
s
k
k
k
k−1
k−1
k
k=0
k
πt

and where P̂st , P̂rt , p̂πb are each an empirical frequency (histogram) estimator, can in factt
be rewritten simply as
T Z
X
t=0

rt P̂rt (rt |st , at )

t 
Y


πke (ak |sk )P̂sk (sk |sk−1 , ak−1 ) d(Hat , rt ).

k=0

This is essentially a model-based OPE estimator, where we first fit all MDP parameters
and then explicitly integrate with respect to the resulting estimated trajectory densityefunction in order to compute the expectation ρπ . This is also often called the G-formulain the causal inference literature (Hernan and Robins, 2019; Robins, 1986). In the tabular
setting, the efficiency of this estimator is immediate since it is exactly the (parametric)
MLE estimate for this setting, which is well known to achieve the Cramér-Rao lower bound,and the Cramér-Rao lower bound is the efficiency bound in the tabular setting since the23Kallus and Uehara

model is parametric. In the case of continuous state and/or action spaces, the simple
extension of replacing P̂sk+1 (sk+1 |sk , ak ), P̂rt (rt |st , at ) with some nonparametric conditional
density estimators would have poor performance since the nonparametric density estimationis unstable and would significantly inflate the variance. Alternatively, if we extend theestimator by instead estimating µt or wt nonparametrically, the above already argues why
we expect this would generally be inefficient.
Similarly to the model-based approach, in the tabular case, our results in the next section
have shown that also simple DM estimates based on q-function estimation are also efficient,since in the tabular case q-functions are parametric. Moreover, DRL was already shown tobe efficient in the the tabular MDP setting with any parametric µt estimator as they all√
have Op (1/ n) convergence in this setting, and the particular choice of estimator does notaffect this (see also Remark 11). Hence, DRL was the first efficient OPE estimator, both ingeneral and in the tabular MDP setting in particular.
Remark 17 (Other estimators for µt ). Xie et al. (2019); Yin and Wang (2020) may in fact
both offer alternative estimators for wt and hences µt (in their respective settings) that may
be used in DRL and either will ensure efficiency for DRL (see Remark 11). In particular,
since λt may have variance growing exponentially in T , this may affect the variance of
estimates of µt based on the regression of it on st or on st , at , as studied above. Although
this will not appear in the leading term of the variance of DRL and will not affect efficiency,
it may still be a concern. Developing and analyzing alternative estimators for wt and/or µt
may be fruitful future work. For example, still other possible estimation approaches for wt
include a fitted w-iteration: start with ŵ0 ≡ 1, regress ŵt−1 ηt on st using any supervised
regression method to obtain ŵt , and repeat.

4. Estimating the q-function and Efficiency Under M1q , M2q
In this section, we discuss the estimation of q-functions in an off-policy manner, parametricallyor nonparametrically, which can be plugged into our estimators, ρ̂DRL(M1 ) , ρ̂DRL(M2 ) . On
the way, we also derive the semiparametric efficiency bound when we impose parametricrestrictions on q-functions, i.e., the models M1q , M2q .
To do this, we will leverage a recursive definition of the q-functions (Bertsekas, 2012).
Under M1 , the following recursion equation holds:qt = E rt + Eπe qt+1 | Hst+1 | Hat .
(14)
Under M2 , we can further replace Hst+1 with st+1 and Hat with (st , at ) in the above.
The recursion in Eq. (14) can equivalently be written as a set of conditional moment
equations satisfied by the q-functions:mt (Hat ; {q1 , . . . , qT }) = 0 ∀t ≤ T,
(15)
 0

 0


0
0
where mt (Hat ; q1 , . . . , qT ) = E rt + Eπe qt+1 (Hat +1 ) | Hst+1 − qt (Hat ) | Hat .
This formulation of the q-function in terms of conditional moment equations, along withe
the observation that ρπ = E [Eπe [q0 (s0 , a0 ) | s0 ]] is determined by the q-function, allows us
both to estimate the q-function efficiently, either parametrically and nonparametrically, andto characterize the efficiency bounds under M1q and M2q . We start with the latter.
24

Double Reinforcement Learning in Markov Decision Processes

4.1. Efficiency Bounds Under M1q , M2q
In this section we consider the models where we restrict q-functions parametrically:
M1q = {Pπb ∈ M1 : ∃βt∗ ∈ Θβt , qt (Hat ) = qt (Hat ; βt∗ ) ∀t ≤ T } ,M2q = {Pπb ∈ M2 : ∃βt∗ ∈ Θβt , qt (st , at ) = qt (st , at ; βt∗ ) ∀t ≤ T } ,where qt (Hat ; βt ) or qt (st , at ; βt ) is some parametric model for the q-function at time t that iscontinuously differentiable with respect to the parameter βt , Θβt is some compact parameter
space, and βt∗ is the true parameter, which is assumed to lie in the interior of Θβt . For
brevity we define vt (Hst ; βt ) = Eπe [qt (Hat ; βt ) | Hst ] and similarly vt (st ; βt ).
Under M1q , M2q , Eq. (14) can be rephrased as as a set of conditional moment restrictions
on the parameter β defined by β = (β1> , · · · βT> )> . In particular, overloading notation andletting mt (Hat ; β) = mt (Hat ; {q1 (·, β1 ), . . . , qT (·, βT )}), we have that β is defined by theset of conditional moment equations mt (Hat ; β) = 0 ∀t ≤ T . This observation is key inestablishing the following result.
Theorem 22 (Efficiency bound under M1q , M2q ). Define eq,t = rt + vt+1 − qt
At = Ct−1 + Ct−1 Bt At+1 Bt> Ct−1 ,hi∗
Bt = E ∇βt qt (Hat ; βt∗ )var(eq,t | Hat )−1 ∇>
v
(H
;
β
)
st+1
t+1 ,
βt+1 t+1ih∗Ct = E ∇βt qt (Hat ; βt∗ )var(eq,t | Hat )−1 ∇>βt qt (Hat ; βt ) ,hi−1∗
AT = E ∇βT qT (HaT ; βT∗ )var(eq,T | HaT )−1 ∇>,βT qT (HaT ; βT )∗B−1 = E[η0 (s0 , a0 )∇>β0 q0 (s0 , a0 ; β0 )].

Then>
EffBd(M1q ) = var (v0 ) + B−1 A0 B−1.

Moreover, the efficiency bound for estimating βt is At .
Finally, the corresponding efficiency bounds under M2q are given by replacing Hst+1
with st+1 and Hat with (st , at ) everywhere in the above.
Remark 18. When T = 1, EffBd(M1q ) above is equal to>var[v0 (s0 )] + B−1 A0 B−1
,

where

∗ −1
A0 = E[∇β0 q0 (s0 , a0 ; β0∗ )var[r0 | s0 , a0 ]−1 ∇>β q0 (s0 , a0 ; β0 )] ,The Matrix Cauchy-Schwarz Inequality (Tripathi, 1999) immediately shows that this is
upper bounded by EffBd(M1 ), as is also implied by M1q ⊂ M1 albeit less directly.
4.2. Parametric Estimation of q-functionse

Next, we consider an estimation method for βt and ρπ . Given the above observations, anatural way to estimate β is by solving the following set of conditional moment equations
given by mt (Hat ; β) = 0 ∀t ≤ T . For example, one approach when the q-model is a linear
model specified as βt> φt (Hat ) for some dφt -dimensional feature expansion φt is to choose β̂25

Kallus and Uehara2
P PdφtEn rt + vt+1 (Hst+1 ; βt+1 ) − qt (Hat ; βt ) φti (Hat ) , which correto minimize Tt=0 i=1
sponds exactly to backward-recursive ordinary least squares. That is, first rT is regressed onφt (HaT ) to obtain β̂T , then qT (HaT ; β̂T ) is averaged over πTe (aT | HsT ) to obtain v̂T , thenrT −1 + v̂T is regressed on φt (HaT −1 ) to obtain β̂T −1 , and so on.
Although such an estimator can achieve the rate Op (n−1/2 ) under correct specificationand standard conditions for M -estimators, it might not yield an efficient estimator foreβ or for ρπ . When the q-model is linear as above, this can be easily solved by instead
applying any efficient variant of the generalized method of moments (GMM), such as twostep GMM (Hansen, 1982; Hansen et al., 1996), to the set of moment equations given by
mt (Hat ; β)φti (Hat ) = 0 ∀t ≤ T, i ≤ dφt . This is almost the same as the above backwardrecursive ordinary least squares but with an optimal weighting of the different moment
conditions in the sum above.
When the q-model may be nonlinear, we can obtain an efficient estimator by instead
applying the method of Hahn (1997) to our set of conditional moment equations. Specifically,
we can consider the set of T mn moment equations E [mt (Hat ; β)φti (Hat )] = 0 ∀t ≤ T, i ≤ mn ,where φt1 (Hat ), φt2 (Hat ), . . . is a basis expansion of the L2 -space and mn → ∞ as n → ∞.
Then, applying any efficient variant of GMM to this set of moment conditions will yield an
efficient estimator β̂ of β.
In all of the above, replacing Hst+1 with st+1 and Hat with (st , at ), the same techniques
can be applied in M2 . In either case, once we have an efficient estimate β̂ of β, an efficienteestimate for ρπ , achieving
the isemiparametric efficiency bound in the appropriate model, ishgiven by ρ̂DM = En v0 (s0 ; β̂0 ) .
Remark 19 (Tabular setting). Consider a tabular case. Then, by treating qt (st , at ) =β > φ(at , st ), where φ(a, s) = (I(a = a†1 , s = s†1 ), · · · , I(a = a†|At | , s = s†|St | )) and a†i and s†i are
the elements of the finite At and St , we can observe that Effbd(M2q ) = Effbd(M2 ) by somealgebra. This result is natural since M2q = M2 in the tabular setting.
4.3. Nonparametric Estimation of q-functions
The above observation in Eq. (14) that q-functions satisfy a set of conditional moment
equations also lends itself to nonparametric estimation of the q-functions. In this section webriefly review how one approach to this, following the application of the method of Ai andChen (2012) to this set of conditional moment equation, can obtain the necessary fourth-root
rates for use in DRL.
The estimator {q̂t }Tt=0 is constructed as the following sieve minimum distance estimator:{q̂t }Tt=0 ∈ arg minTX

qt ∈Λt,n ∀t≤T t=0hi
En m̂t (Hat ; qt )Σ̂−1
m̂
(H
;
q
)
t
at t ,twhere m̂t (Hat ; qt ) is a nonparameric estimator for mt (Hat ; qt ), Σ̂t is a nonparametric estimator for var (eq,t | Hat ), and Λk,n is a sequence of approximation space whose union ∪∞n=1 Λt,n
is dense in some infinite dimensional space Λt . Alternatively, in M2 , we replace Hat with(st , at ) in the above.
26

Double Reinforcement Learning in Markov Decision Processes

Table 1: Experiment from Section 5.1: RMSE (and standard errors).
Settingnρ̂ISρ̂DRL(M1 )

ρ̂DM

ρ̂MIS

ρ̂DRL(M2 )

(1)

1500
3000
4500

42.4 (12.4)
20.4 (3.1)
20.2 (3.1)

36.1 (16.8)
7.8 (0.8)
6.6 (0.75)

0.70 (0.002)
0.50 (0.001)
0.43 (0.001)

40.8 (12.5)
20.8 (2.8)
21.5 (3.5)

0.70 (0.002)
0.50 (0.001)
0.43 (0.001)

(2)

1500
3000
4500

42.4 (12.4)
20.4 (2.5)
20.2 (3.1)

77.6 (29.1)
36.6 (6.9)
34.4 (9.6)

10.8 (0.002)
10.8 (0.001)
10.8 (0.001)

40.8 (12.5)
20.8 (2.8)
21.5 (3.5)

10.3 (3.5)
6.0 (0.6)
5.5 (2.0)

(3)

1500
3000
4500

42.4 (12.4)
20.4 (3.1)
20.2 (3.1)

36.1 (16.8)
7.8 (0.8)
6.6 (0.75)

0.70 (0.002)
0.50 (0.001)
0.43 (0.001)

87.7 (25.5)
37.3 (3.2)
53.5 (15.1)

0.73 (0.03)
0.51 (0.002)
0.44 (0.005)

Ai and Chen (2003) prove that applying the above with appropriate nonparametricestimators, under some smoothness conditions, we can obtain kq̂t − qt kF,t = op (n−1/4 ), wherek · kF,t is the Fisher metric, which in our setting of Eq. (15) is defined askg(Hat )k2F,t = E[var(eq,t | Hat )g 2 + var(eq,t−1 | Hat−1 )Eπe [g(Hat ) | Hst ]2 ].
We omit the details and refer the interested reader to Ai and Chen (2003). We only provethat this norm is in fact equivalent to the L2 -norm under mild conditions.
Lemma 23. Suppose var[eq,t | Hat ] and var[eq,t−1 | Hat−1 ] are bounded away from zero.
Then, k · kF,k and k · k2 are equivalent norms.
This means that, under the appropriate conditions, the estimator q̂ obtains the rate
op (n−1/4 ) in terms of L2 -norm, as necessary for Theorems 6, 10, 13 and 16.

5. Experiments
We now turn to an empirical study of OPE and DRL. First, we construct a simulation
to investigate the effect of using memorylessness on estimation variance as well as theeffect of double robustness on model specification sensitivity. Then, we study comparative
performance of different OPE estimators in two standard OpenAI Gym tasks.
Replication code for all experiments is available at http://github.com/CausalML/
DoubleReinforcementLearningMDP.
5.1. The Effects of Leveraging Memorylessness and of Double Robustness
In this section we consider an MDP with a horizon of T = 30, binary actions, univariate
continuous state, initial state distribution p(s0 ) ∼ N (0.5, 0.2), transition probabilities
Pt (st+1 | st , at ) ∼ N (s + 0.3a − 0.15, 0.2). The target and behavior policies we consider are,respectively,
π e (a | s) ∼ Bernoulli(pe ), pe = 0.2/(1 + exp(−0.1s)) + 0.2U, U ∼ Uniform[0, 1]
π b (a | s) ∼ Bernoulli(pb ), pb = 0.9/(1 + exp(−0.1s)) + 0.1U, U ∼ Uniform[0, 1].
27

Kallus and Uehara

We assume the behavior policy is known. Note that this setting is an MDP and belongs toM2 .
We compare five estimators: ρ̂IS , ρ̂DRL(M1 ) , ρ̂DM , ρ̂MIS , ρ̂DRL(M2 ) when nuisancefunctions qt (s, a) and µt (s) are estimated parametrically. We consider three settings:(1) Both models correct: qt (st , at ) = β1t st + β2t st at + β3t , µt (st , at ) = β4t st + β5t st at + β6t .
(2) Only µ-model correct: qt (st , at ) = β1t s2t + β2t s2t at + β3t , µt (st , at ) = β4t st + β5t st at + β6t .
(3) Only q-model correct: qt (st , at ) = β1t st + β2t st at + β3t , µt (st , at ) = β4t s2t + β5t s2t at + β6t .
Note that in the above, the “correct” models are in fact not exactly correct because Eπe [at | st ]is actually nonlinear in st , but it is very nearly linear in the space of observed st values (forexample, best linear fit for Eπe [at | st ] has an L2 distance 3 × 10−5 on [0, 1], which spans±2.5 standard deviations for s0 ). We therefore treat them as correctly specified.
In all cases, to estimate q-models we use backward-recursive ordinary least squares as in
Section 4.2. To estimate µ-models we use ordinary least squares regression on λt (which isassumed known) as in Eq. (11).
For each n = 1500, 3000, 4500, we consider 50000 Monte Carlo replications. In each
replication, we estimate the q- and µ-models as above and compute, for each setting, each ofρ̂IS , ρ̂DRL(M1 ) , ρ̂DM , ρ̂MIS , ρ̂DRL(M2 ) . We report the RMSE of each estimator in each setting(and the standard error) in Table 1.
Our first immediate observation is that ρ̂DRL(M2 ) nearly dominates all other estimators,
achieving similar or better performance in every setting and sample size. In particular,
in settings (1) and (3), where the q-model is correct, it has performance similar to ρ̂DM .
Note that in settings (1) and (3), ρ̂DM is efficient for M2q per Section 4.2 (or almost so;
it would be efficient if we used efficient GMM instead of one-step GMM). In setting (1),ρ̂DRL(M2 ) is locally efficient, while in setting (3), it is only doubly robust and performs
almost imperceptibly worse than the efficient ρ̂DM .
In setting (2), where the q-model is incorrect, ρ̂DM is inconsistent and ρ̂DRL(M2 ) handily
outperforms it. In the same setting (2), the consistent ρ̂IS and ρ̂MIS also outperform theinconsistent ρ̂DM but not by as much as ρ̂DRL(M2 ) . While ρ̂DRL(M1 ) is doubly robust insetting (2) guaranteeing consistency, unlike the case of ρ̂DRL(M2 ) , the combination of large(unmarginalized) cumulative density ratios and a misspecified q-model leads to still worse
performance in the sample sizes tested.
Generally, ρ̂IS , ρ̂MIS , and ρ̂DRL(M1 ) all have high RMSE due to the significant mismatch
between the behavior and target policies so that cumulative density ratios are very large
and only marginalizing them without also using a q-model helps only a little. In settings (1)and (2), where the µ-model is correct, ρ̂MIS improves on ρ̂IS only slightly, while in setting(3), where µ-model is incorrect, it performs significantly worse. This highlights the potential
danger of misspecifying µ-models compared to the robustness of importance sampling withknown behavior policy (see also Remark 12).
While both ρ̂IS and ρ̂DRL(M1 ) remain consistent throughout all settings, they are outperformed by the also-consistent ρ̂DRL(M2 ) , which leverages the MDP structure of M2 andexhibits local efficiency in setting (1) and doubly robustness in settings (2) and (3).
28Double Reinforcement Learning in Markov Decision Processes

Table 2: Cliff Walking: RMSE (and standard errors)
Size

ρ̂IS

ρ̂DRL(M1 )

ρ̂DM

ρ̂MIS

ρ̂DRL(M2 )

500
1000
1500

18.8 (7.67)
7.99 (0.89)
7.64 (1.63)

3.78(1.14)
0.28 (0.026)
0.098 (0.013)

2.63 (0.01)
1.27 (0.002)
1.01 (0.001)

12.8 (4.96)
5.92 (0.78)
5.55 (1.10)

1.44 (0.29)
0.22 (0.34)
0.075 (0.008)

Table 3: Mountain Car: RMSE (and standard errors)nρ̂ISρ̂DRL(M1 )

ρ̂DM

ρ̂MIS

ρ̂DRL(M2 )

500
1000
1500

6.85 (0.13)
4.73 (0.07)
3.41 (0.04)

3.72 (0.08)
2.12 (0.04)
1.82 (0.02)

4.30 (0.05)
3.40 (0.008)
3.30 (0.008)

6.82 (0.12)
4.83 (0.06)
3.40 (0.05)

3.53 (0.12)
2.07 (0.04)
1.69 (0.03)

5.2. Investigating Performance in RL Tasks: Cliff Walking and Mountain Car
We next compare the same OPE estimators using nonparametric nuisance estimation in two
standard RL settings included in OpenAI Gym (Brockman et al., 2016): Cliff Walking andMountain Car. For further detail on each setting, see Appendix C.
First, we used q-learning to learn an optimal policy for the MDP and define it as π d .
Then we generate the dataset from the behavior policy π b = (1 − α)π d + απ u where π u is a
uniform random policy and α = 0.8. We define the target policy similarly but with α = 0.9.
Again, we assume the behavior policy is known. Note that this πd is fixed in each setting.
We estimate all µ-functions by first estimating w-functions and using Eq. (5). For Cliff
Walking, we use a histogram estimator for w as in Eq. (7). For Mountain Car, we use
a kernel estimator for w as in Eq. (10). We use the Epanechnikov kernel and choose an
optimal bandwidth based on an L2 -risk criterion for t = 1; we then use this bandwidth
for all other t values as well for simplicity. For q-functions, we use backward-recursive
regression
as in Section 4.2. For Cliff-Walking, we use a histogram model, q(s, a; β) =P>
sj ,ak ∈S,A βjk I[sj = s, ak = a]. For Mountain-Car, we use the mode q(s, a; β) = β φ(s, a)where φ(s, a) is a 400-dimensional feature vector based on a radial basis function, generatedusing the RBFSampler method of scikit-learn based on Rahimi and Recht (2008).
We again compare ρ̂IS , ρ̂DRL(M1 ) , ρ̂DM , ρ̂MIS , ρ̂DRL(M2 ) . In each setting we consider
varying evaluation dataset sizes and for each consider 1000 replications. We report the
RMSE of each estimator in each setting (and the standard error) in Tables 2 and 3.
We again find that the performance of ρ̂DRL(M2 ) is superior to all other estimators
in either setting. This is especially true in Cliff Walking. The estimator ρ̂DRL(M2 ) also
improves upon ρ̂IS and ρ̂DM but not as much as ρ̂DRL(M2 ) . The estimator ρ̂MIS offers a slight
improvement over ρ̂IS , but is still outperformed by ρ̂DRL(M2 ) , ρ̂DRL(M1 ) , and ρ̂DM . That
the improvement of ρ̂MIS over ρ̂IS and the overall improvements of ρ̂DRL(M2 ) is starker in
Cliff Walking than in Mountain Car may be attributable to the difficulty of learning wt
nonparametrically in a continuous state space.
29

Kallus and Uehara

6. Conclusions
We established the semiparametric efficiency bounds and efficient influence functions for OPE
under either NMDP or MDP model, which quantify how fast one could hope to estimate
policy value. While in the NMDP case, the influence function we derived has appeared
frequently in OPE estimators, in the MDP case, the influence function is novel and has notappeared in existing estimators. Our results also suggested how one could construct efficientestimators. We used this to develop DRL, which used our newly derived efficient influencefunction, with nuisances estimated in a cross-fold manner. This ensured efficiency under very
weak and mostly agnostic conditions on the nuisance estimation method used. Notably, DRL
is the first efficient OPE estimator for MDPs. In addition, DRL enjoyed double robustness
properties. This efficiency and robustness translated to better performance in experiments.

Acknowledgments
The authors thank Nan Jiang and Yu-Xiang Wang for helpful discussions. This material isbased upon work supported by the National Science Foundation under Grant No. 1846210.

References
Chunrong Ai and Xiaohong Chen. Efficient estimation of models with conditional moment
restrictions containing unknown functions. Econometrica, 71:1795–1843, 2003.
Chunrong Ai and Xiaohong Chen. Semiparametric efficiency bound for models of sequential
moment restrictions containing unknown functions. Journal of Econometrics, 170:442–457,
2012.
András Antos, Csaba Szepesvári, and Rémi Munos. Learning near-optimal policies withbellman-residual minimization based fitted policy iteration and a single sample path.
Machine Learning, 71:89–129, 2008.
Peter L. Bartlett, Olivier Bousquet, and Shahar Mendelson. Local rademacher complexities.
The Annals of Statistics, 33:1497–1537, 2005.
David Benkeser and Mark van der Laan. The highly adaptive lasso estimator. In 2016 IEEE
International Conference on Data Science and Advanced Analytics (DSAA), volume 2016,
pages 689–696. IEEE, 2016.
Dimitri P Bertsekas. Dynamic programming and optimal control. Athena Scientific optimization and computation series. Athena Scientific, Belmont, Mass, 4th ed. edition,
2012.
Aurélien Bibaut and Mark van der Laan. Fast rates for empirical risk minimization over
cdlg functions with bounded sectional variation norm. arXiv preprint arXiv:1907.09244,
2019a.
Aurelien Bibaut, Ivana Malenica, Nikos Vlassis, and Mark van der Laan. More efficient
off-policy evaluation through regularized targeted learning. In Proceedings of the 36th30

Double Reinforcement Learning in Markov Decision Processes

International Conference on Machine Learning, volume 97 of Proceedings of Machine
Learning Research, pages 654–663, 2019.
Aurélien F Bibaut and Mark J van der Laan. Fast rates for empirical risk minimization over c\adl\ag functions with bounded sectional variation norm. arXiv preprint
arXiv:1907.09244, 2019b.
P. J. Bickel, C. A. J. Klaassen, Y. Ritov, and J. A. Wellner. Efficient and Adaptive Estimationfor Semiparametric Models. Springer, 1998.
Clive G. Bowsher and Peter S. Swain. Identifying sources of variation and the flow of
information in biochemical networks. Proceedings of the National Academy of Sciences,
109, 2012.
G. Brockman, V. Cheung, L. Pettersson, J. Schneider, J. Schulman, J. Tang, , and W Zaremba.
Openai gym. arXiv preprint arXiv:1606.01540, 2016.
Bibhas Chakraborty and EE Moodie. Statistical methods for dynamic treatment regimes.
Springer, 2013.
Gary Chamberlain. Comment: Sequential moment restrictions in panel data. Journal ofBusiness & Economic Statistics, 10:20–26, 1992.
Xiaohong Chen. Chapter 76 large sample sieve estimation of semi-nonparametric models.
Handbook of Econometrics, 6:5549–5632, 2007.
Victor Chernozhukov, Denis Chetverikov, Mert Demirer, Esther Duflo, Christian Hansen,
Whitney Newey, and James Robins. Double/debiased machine learning for treatment andstructural parameters. Econometrics Journal, 21:C1–C68, 2018.
Iván Dı́az. Machine learning in the estimation of causal effects: targeted minimum loss-based
estimation and double/debiased machine learning. Biostatistics, 2019.
Miroslav Dudik, Dumitru Erhan, John Langford, and Lihong Li. Doubly robust policy
evaluation and optimization. Statistical Science, 29:485–511, 2014.
Ashkan Ertefaie and Robert L Strawderman. Constructing dynamic treatment regimes overindefinite time horizons. Biometrika, 105:963–977, 2018.
M. Farajtabar, Y. Chow, and M. Ghavamzadeh. More robust doubly robust off-policy
evaluation. In Proceedings of the 35th International Conference on Machine Learning,
pages 1447–1456, 2018.
Omer Gottesman, Fredrik Johansson, Matthieu Komorowski, Aldo Faisal, David Sontag,
Finale Doshi-Velez, and Leo Anthony Celi. Guidelines for reinforcement learning in
healthcare. Nat Med, 25:16–18, 2019.
László Györfi, Michael Kohler, Adam Krzyzak, and Harro Walk. A distribution-free theoryof nonparametric regression. Springer Science & Business Media, 2006.
31

Kallus and Uehara

J Hahn. Efficient estimation of panel data models with sequential moment restrictions.
Journal of Econometrics, 79:1–21, 1997.
Jinyong Hahn. On the role of the propensity score in efficient semiparametric estimation of
average treatment effects. Econometrica, 66:315–331, 1998.
Lars Peter Hansen. Large sample properties of generalized method of moments estimators.
Econometrica, 50:1029–1054, 1982.
Lars Peter Hansen, John Heaton, and Amir Yaron. Finite-sample properties of somealternative gmm estimators. Journal of Business & Economic Statistics, 14:262–280, 1996.
M.A. Hernan and J.M. Robins. Causal Inference. Boca Raton: Chapman & Hall/CRC,
2019.
Keisuke Hirano, Guido W Imbens, and Geert Ridder. Efficient estimation of averagetreatment effects using the estimated propensity score. Econometrica, 71:1161–1189, 2003.
Masaaki Imaizumi and Kenji Fukumizu. Deep neural networks learn non-smooth functions
effectively. arXiv preprint arXiv:1802.04474, 2018.
N. Jiang and L. Li. Doubly robust off-policy value evaluation for reinforcement learning. In
Proceedings of the 33rd International Conference on International Conference on Machine
Learning, pages 652–661, 2016.
Nathan Kallus and Masatoshi Uehara. Intrinsically efficient, stable, and bounded off-policy
evaluation for reinforcement learning. In Advances in Neural Information Processing
Systems (NeurIPS), 2019.
Nathan Kallus and Masatoshi Uehara. Double reinforcement learning for efficient and robust
off-policy evaluation. In Proceedings of the 37th International Conference on Machine
Learning, 2020.
Khashayar Khosravi, Greg Lewis, and Vasilis Syrgkanis. Non-parametric inference adaptive
to intrinsic dimension. arXiv preprint arXiv:1901.03719, 2019.
Chris A. J. Klaassen. Consistent estimation of the influence function of locally asymptoticallylinear estimators. The Annals of Statistics, 15:1548–1562, 1987.
Michael R Kosorok. Introduction to Empirical Processes and Semiparametric Inference.
Springer Series in Statistics. Springer New York, New York, NY, 2008.
Michail Lagoudakis and Ronald Parr. Least-squares policy iteration. Journal of Machine
Learning Research, 4:1107–1149, 2004.
Hoang Le, Cameron Voloshin, and Yisong Yue. Batch policy learning under constraints. In
International Conference on Machine Learning, pages 3703–3712, 2019.
L. Li, R. Munos, and C. Szepesvari. Toward minimax off-policy value estimation. In
Proceedings of the 18th International Conference on Artificial Intelligence and Statistics,
pages 608–616, 2015.
32

Double Reinforcement Learning in Markov Decision Processes

Qi Li and Jeffrey Scott Racine. Nonparametric econometrics : theory and practice. Princeton
University Press, Princeton, N.J., 2007.
Qiang Liu, Lihong Li, Ziyang Tang, and Dengyong Zhou. Breaking the curse of horizon:
Infinite-horizon off-policy estimation. In Advances in Neural Information Processing
Systems 31, pages 5356–5366. 2018.
Daniel J. Luckett, Eric B. Laber, Anna R. Kahkoska, David M. Maahs, Elizabeth MayerDavis, and Michael R. Kosorok. Estimating dynamic treatment regimes in mobile healthusing v-learning. Journal of the American Statistical Association, pages 1–34, 2018.
A. Rupam Mahmood, Hado P van Hasselt, and Richard S Sutton. Weighted importancesampling for off-policy learning with linear function approximation. In Advances in Neural
Information Processing Systems 27, pages 3014–3022. 2014.
T. Mandel, Y. Liu, S. Levine, E. Brunskill, and Z Popovic. Off-policy evaluation acrossrepresentations with applications to educational games. In Proceedings of the 13th
International Conference on Autonomous Agentsand Multi-agent Systems, pages 1077–
1084, 2014.
Shie Mannor, Duncan Simester, Peng Sun, and John N Tsitsiklis. Bias and varianceapproximation in value function estimates. Management Science, 53(2):308–322, 2007.
Remi Munos, Tom Stepleton, Anna Harutyunyan, and Marc Bellemare. Safe and efficient
off-policy reinforcement learning. In Advances in Neural Information Processing Systems
29, pages 1054–1062. 2016.
S. A. Murphy. Optimal dynamic treatment regimes. Journal of the Royal Statistical Society:
Series B (Statistical Methodology), 65:331–355, 2003.
S A Murphy, M J van der Laan, and J M Robins. Marginal mean models for dynamic
regimes. Journal of the American Statistical Association, 96:1410–1423, 2001.
W. K. Newey and D. L Mcfadden. Large sample estimation and hypothesis testing. Handbookof Econometrics, IV:2113–2245, 1994.
D. Precup, R. Sutton, and S Singh. Eligibility traces for off-policy policy evaluation. In
Proceedings of the 17th International Conference on Machine Learning, pages 759–766,
2000.
Ali Rahimi and Benjamin Recht. Random features for large-scale kernel machines. In
Advances in Neural Information Processing Systems 20, pages 1177–1184. 2008.
James Robins. A new approach to causal inference in mortality studies with a sustained
exposure periodapplication to control of the healthy worker survivor effect. Mathematicalmodelling, 7(9-12):1393–1512, 1986.
James M Robins. Marginal structural models versus structural nested models as tools forcausal inference. In Statistical models in epidemiology, the environment, and clinical trials,
pages 95–133. Springer, 2000.
33

Kallus and Uehara

James M Robins, Miguel Angel Hernán, and Babette Brumback. Marginal structural models
and causal inference in epidemiology. Epidemiology, 11:551, 2000.
Andrea Rotnitzky and Stijn Vansteelandt. Double-robust methods. In Handbook of missing
data methodology. In Handbooks of Modern Statistical Methods, pages 185–212. Chapman
and Hall/CRC, 2014.
Andrea Rotnitzky, Ezequiel Smucler, and James Robins. Characterization of parameterswith a mixed bias property. arXiv preprint arXiv:1509.02556, 2019.
D. Scharfstein, A. Rotnizky, and J. M. Robins. Adjusting for nonignorable dropout usingsemi-parametric models. Journal of the American Statistical Association, 94:1096–1146,
1999.
Xiaotong Shen. On methods of sieves and penalization. The Annals of Statistics, 25:2555–2591, 1997.
Charles J. Stone. Rejoinder: The use of polynomial splines and their tensor products in
multivariate function estimation. The Annals of Statistics, 22:179–184, 1994.
Richard S Sutton. Reinforcement learning : an introduction. MIT Press, Cambridge, Mass.,
2018.
Adith Swaminathan and Thorsten Joachims. The self-normalized estimator for counterfactual
learning. In Advances in Neural Information Processing Systems 28, pages 3231–3239.
2015.
P. Thomas and E. Brunskill. Data-efficient off-policy policy evaluation for reinforcement
learning. In Proceedings of the 33rd International Conference on Machine Learning, pages
2139–2148, 2016.
G Tripathi. A matrix extension of the Cauchy-Schwarz inequality. Economics Letters, 63:
1–3, 1999.
Anastasios A Tsiatis. Semiparametric Theory and Missing Data. Springer Series in Statistics.
Springer New York, New York, NY, 2006.
Mark J. van der Laan and James M Robins. Unified Methods for Censored Longitudinal
Data and Causality. Springer Series in Statistics,. Springer New York, New York, NY,
2003.
A. W van der Vaart. On differentiable functionals. Ann. Statist., 19:178–204, 03 1991.
A. W. van der Vaart. Asymptotic statistics. Cambridge University Press, Cambridge, UK,1998.
A. W. van der Vaart. Semiparametric Statistics. Lecture Notes in Mathematics ; 1781.
Springer Berlin Heidelberg, Berlin, Heidelberg, 2002.
Karel Vermeulen. Semiparametric efficiency. Gent, Faculteit WetenschappenVakgroep
Toegepaste Wiskunde en Informatica, 2010.
34

Double Reinforcement Learning in Markov Decision Processes

Stefan Wager and Guenther Walther. Adaptive concentration of regression trees, withapplication to random forests. arXiv preprint arXiv:1503.06388, 2016.
Yu-Xiang Wang, Alekh Agarwal, and Miroslav Dudik. Optimal and adaptive off-policy
evaluation in contextual bandits. In Proceedings of the 34th International Conference on
Machine Learning, pages 3589–3597, 2017.
Tengyang Xie, Yifei Ma, and Yu-Xiang Wang. Towards optimal off-policy evaluation forreinforcement learning with marginalized importance sampling. In Advances in Neural
Information Processing Systems 32, pages 9665–9675. 2019.
Ming Yin and Yu-Xiang Wang. Asymptotically efficient off-policy evaluation for tabular
reinforcement learning. In Proceedings of the 23nd International Workshop on Artificial
Intelligence and Statistics (To appear), 2020.
Baqun Zhang, Anastasios A. Tsiatis, Eric B. Laber, and Marie Davidian. Robust estimationof optimal dynamic treatment regimes for sequential treatment decisions. Biometrika, 100:681–694, 2013.
Wenjing Zheng and Mark J van der Laan. Cross-validated targeted minimum-loss-based
estimation. In Targeted Learning: Causal Inference for Observational and Experimental
Data, Springer Series in Statistics, pages 459–474. Springer New York, New York, NY,
2011.

35

Kallus and Uehara

Appendix A. Notation
We first summarize the notation we use in Table 4 and the abbreviations we use in Table 5.
Notice in particular that, following empirical process theory literature, in the proofs we also
use P to denote expectations (interchangeably with E).
Table 4: Notation∇β
rt , st , at ,Jrt , Jst , Jat
Hst , Hat
πt (at |Hst ), πt (at |st )πte , πtb
ρπvt = vt (Hst ), vt (st )qt = qt (Hat ), qt (st , at )
λt
µt
ηt
Λ
M
M1 , M1b , M1q
M2 , M2b , M2q
C, Rmax
Q
L(A|B)
k · kpEπ [·], Pπ
E[·], P
En [·], Pnnj
Enj , Pnj
Gn
Asmse[·], var[·]N (a, b)Uni[a, b]
An = op (an )An = Op (an )ΛαdDifferentiation with respect to β
Reward, state, action at t
History up to time rt , st , at , including reward variables
History up to time st , at , excluding reward variables
Policy in NMDP and MDP case, respectively
Target and behavior policies at t, respectively
P
Policy value, Eπ [ Tt=0 rt ]Value function at t, in M1 , M2 respectively
q-function at t, in M1 , M2 respectively
Q
Cumulative density ratio tk=0 πte /πtb
Marginal density ratio E[λt | st , at ]Instantaneous density ratio πte /πtb
Tangent space
A model for the data generating distribution
NMDP model with unknown behavior policy,
known behavior policy, and parametric q-function, respectively
MDP model with unknown behavior policy,
known behavior policy, and parametric q-function, respectively
Upper bound of density ratio and reward, respectively
Projection of A onto B
Direct sum
Lp -norm E[f p ]1/p
Inequality up to constantExpectation with respect to a sample from a policy πSame as above for π = π b
Empirical expectation (based on sample from a behavior policy)
The size of Dj
Empirical expectation on Dj√
Empirical process n(Pn − P)
Asymptotic variance, variance
Normal distribution with mean a and variance b
Uniform distribution on [a, b]
The term An /an converges to zero in probability
The term An /an is bounded in probability
Hölder space with smoothness α with a dimension d

36

Double Reinforcement Learning in Markov Decision Processes

Table 5: Abbreviations
Non-Markov Decision Process
Markov Decision Process
Reinforcement Learning
Contextual Bandit
Off policy Evaluation
Maximum Likelihood Estimation
Regular and Asymptotic Linear
Consistent and Asymptotically Normal
Mean Squared Error

NMDP
MDP
RL
CB
OPE
MLE
RAL
CAN
MSE

Appendix B. Proofs
Before going into details of the proof, we summarize definitions and proofs to derive asemiparametric lower bound. As we mentioned in Section 1.2, for a complete and rigorous
treatment, refer to Bickel et al. (1998); van der Laan and Robins (2003); van der Vaart(2002). Additional accessible treatments are also given in (Bibaut and van der Laan, 2019a;
Tsiatis, 2006; Vermeulen, 2010)
B.1. Semiparametric theory
We overload notation on Section 1.2. We denote the all of the history {H(i) }ni=1 as Hn , theestimand as R(F ) : M → R and the estimator as R̂ : Hn → R. First, we introduce somedefinitions.
Definition 1 (One-dimensional submodel and its score function). A one-dimensional submodel of M that passes through F at 0 is a subset of M of the form {F :  ∈ [−a, a]} forsome small a > 0 s.t. F=0 = F . The score of the submodel F at θ = 0 is defined ass(H) =log(dF /dµ)(H)|=0 .
d

Definition 2 (Tangent space). The tangent space of a model M at F denoted by TM (F ) is
the linear closure of the set of score functions of the all one-dimensional submodels regarding
M that pass through F .
Definition 3 (Influence function of estimators). An estimator R̂(Hn ) is asymptoticallylinear with influence function (IF) ψ(H) if
√

n
√
1 X
n(R̂(H ) − R(F )) = √ψ(H(i) ) + op (1/ n).
n
n

i=1

Definition 4 (Pathwise differentiability). A functional R(F ) is pathwise differentiable at
F w.r.t the model M (or w.r.t the tangent space TM (F )) if there exists a function DF (H)such that for all submodels {F : } in M satisfying F=0 = F anddR(F )
|=0 = E[DF (H)s(H)],d37

Kallus and Uehara

where s(H) is a corresponding score function for F . The function DF (H) is called a gradientof R(F ) at F w.r.t the model M. The efficient IF (EIF) of R(F ) w.r.t the model M iscalled a canonical gradient D̃F (H), which is the unique gradient of R(F ) at F w.r.t themodel M that belongs to the tangent space TM (F ).
Next, we define regular estimators. Regular estimators means estimators whose limiting
distribution is insensitive to local changes to the data generating process. It excludes awell-known Hodge estimator. Here, we denote a submodel with some score function g in agiven tangent space TM (F ) as {Ft,g ; t ∈ [−a, a]}.
Definition 5 (Regular estimators). An estimator sequence Tn is called regular at F forR(F ) w.r.t the model M (or w.r.t the tangent space TM (F )), if there exists a probability
measure L such that√n{Tn − R(F1/√n,g )}d(F1/√n,g )→L, for every g ∈ TM (F ).

The following three theorems imply that influence functions of the estimators R̂(F ) forR(F ) and gradients of R(F ) correspond to each other, and how to construct an efficientestimator. These theorems are based on Theorem 3.1 (van der Vaart, 1991).
Theorem 24 (Influence functions are gradients). Under certain regularity conditions, forP ∈ M, suppose R̂(Hn ) is a regular estimator of R(F ) w.r.t the model M, and that it isasymptotically linear with influence function DF (H). Then, R(F ) is pathwise differentiable
at F w.r.t M and DF (H) is a gradient of R(F ) at F w.r.t M.
Theorem 25 (Gradients are influence functions). Under certain regularity conditions, if a
DF (H) is a gradient of R(F ) at F w.r.t the model M, there exists an asymptotically linear
estimator of R(F ) with influence function DF (H), which is regular w.r.t the model M.
Corollary 26 (Characterization of efficient influence functions). The efficient influencefunction is the projection of any gradient onto the tangent space TM (F ).
Note that gradients w.r.t the model M are not unique if the model M is not a fully
nonparametric model. If the underlying model is fully nonparametric model, the gradient is
unique.
Strategy to calculate the EIF With the abovementioned definitions and theorems in
mind, our general strategy to compute efficient influence functions is as follows.
1. Calculate some gradient DF (H) (a candidate of EIF) of the target functional R(F )w.r.t M
2. Calculate the the tangent space TM (F ) at F
3. Show that some candidate of EIF is orthogonal to the orthogonal tangent space, i.e.,the candidate of EIF lies in the tangent space. Then, this implies that a candidate of
EIF is actually the EIF.
The other common strategy is calculating some gradient and projection it onto TM (F ).
38

Double Reinforcement Learning in Markov Decision Processes

Optimalites The efficiency bound has the following interpretations. First, the efficiency
bound is the lower bound in a local asymptotic minimax sense (van der Vaart, 1998,
Thm. 25.20).
Theorem 27 (Local Asymptotic Minimax theorem). Let R(F ) be pathwise diffentiable at F
w.r.t the model M with EIF D̃F (H). If TM (F ) is a convex cone, for any estimator sequenceR̂(Hn ), and subconvex loss function l : R → [0, ∞),Z√n√sup lim sup EF1/√n,g [l[ n{R̂(H ) − R(F1/ n,g )}]] ≥ l(u)dN (0, varF [D̃F (H)])(u),In→∞ g∈I

where the first supremum is taken over all finite subsets I of the tangent set.
Corollary 28. Under the same assumptions of Theorem 27,
Z√n
inf lim inf
sup EQ [l[ n{R̂(H ) − R(Q)}]] ≥ l(u)dN (0, varF [D̃F (H)])(u),
δ>0 n→∞ kQ−F kT ≤δwhere k · kT is a total variation distance.
Other different type of optimality is seen in the following theorem. The following theorem
state that an asymptotic variance of every regular estimator sequence R̂(Hn ) with limiting
distribution L is bounded below E[D̃F2 (H)] (van der Vaart, 1998, Thm. 25.21).
Theorem 29 (Convolution theorem). Let R(F ) be pathwise differentiable at F w.r.t themodel M with EIF D̃F (H). Let R̂(Hn ) be a regular estimator sequence at F w.r.t the tangent
space TM (F ) with limiting distribution L. Then, if the tangent space TM (F ) is a cone, then,
the termZ
u2 dL(u) − E[D̃F2 (H)]is non-negative.
B.2. Proof
Proof of Theorem 1.
Efficient influence function under M1 . The entire regular (regular model as definedin Chapter 7 van der Vaart, 1998) parametric submodel under M1 is{pθ (s0 )pθ (a0 |s0 )pθ (r0 |Ha0 )pθ (s1 |Ha0 )pθ (a1 |Hs1 )pθ (r1 |Ha1 ) · · · pθ (rT |HaT )},where it matches with a true pdf at θ = 0.
The score function of the model M1 is decomposed asg(JsT ) =TX

∇ log pθ (sk |Hak−1 ) +

k=0

=

T
X
k=0

T
X

∇ log pθ (ak |Hsk ) +k=0

gsk |Hak−1 +

T
X
k=0

gak |Hsk +

T
X
k=0

T
X
k=0

39

grk |Hak .

∇ log pθ (rk |Hak )Kallus and Uehara

We first calculate a derivative for the target functional w.r.t the model M1 . Note that
this derivative is not unique. We have
" T #X
∇θ Eπe
rt
t=0

"Z
= ∇θ

T
X

rt

t=0

=

T
X

( T
Y

)

#

pθ (sk |Hrk−1 )pπe (ak |Hsk )pθ (rk |Hak ) dµ(JsT )k=0{Eπe [{Eπe (rc |s0 ) − Eπe (rc )}gs0 ] + Eπe [{rc − Eπe (rc |Hac )}grc |Hac ]

c=0

"

T
X

"

+ Eπ eEπe#rt |Hsc+1 − Eπe

t=c+1

=

T
X

T
X

"

#!
rt |Hac#
gsc+1 |Hac }t=c+1{Eπe [{Eπe (rc |s0 ) − Eπe (rc )}g(JsT +1 )] + Eπe [{rc − Eπe (rc |Hac )}g(JsT )]

c=0

"

"

+ Eπ e

Eπ e

T
X

#
rt |Hsc +1 − Eπe

t=c+1

"
=E

πe

−ρ

+

T
X

T
X

"

#!
rt |Hac#g(JsT ) }t=c+1(
λc rc − λcTXEπe (rt |Hac ) + λc−1

t=c

c=0

T
X

)#
Eπe (rt |Hsc )!
g(JsT ) .

t=c

This concludes that the following function is a derivative:
(
)
T
T
T
X
X
X
e
π
1
ρM
+
λc rc − λcEπe (rt |Hac ) − λc−1
Eπe (rt |Hsc ) .
eff = −ρ
t=c

c=0

(16)

t=c

Next, we show that this derivative is the efficient influence function. In order to show
this, we calculate the tangent space of model M1 . The tangent space of the model M1 isthe product space:
M
M M
(At
Bt
Ct ),
0≤t≤T

At = {q(st , Hat−1 ); E[q(st , Hat−1 )|Hat−1 ] = 0, q ∈ L2 },Bt = {q(at , Hst ); E[q(at , Hst )|Hst ] = 0, q ∈ L2 },Ct = {q(rt , Hat ); E[q(rt , Hat )|Hat ] = 0, q ∈ L2 }.

The orthogonal space of the tangent space is the product of
M
M M
Ct0 )
(A0t
Bt0
0≤t≤T

40

(17)

Double Reinforcement Learning in Markov Decision Processes

such that
A0t

M

At = A00t , A00t = {q(Jst ); E[q(Jst )|Jrt−1 ] = 0, q ∈ L2 },

Bt0

M

Bt = Bt00 , Bt00 = {q(Jat ); E[q(Jat )|Jst ] = 0, q ∈ L2 },Ct0

M

Ct = Ct00 , Ct00 = {q(Jrt ); E[q(Jrt )|Jat ] = 0, q ∈ L2 }.

More specifically, we have the following lemma.
Lemma 30. The orthogonal tangent space is represented as
A0t = q(Jst ) − E[q(Jst )|Hst ]; E[q(Jst )|Jrt−1 ] = 0, q ∈ L2 ,
Bt0 = q(Jat ) − E[q(Jat )|Hat ]; E[q(Jat )|Jst ] = 0, q ∈ L2 ,

Ct0 = q(Jrt ) − E[q(Jrt )|Hat , rt ]; E[q(Jrt )|Jat ] = 0, q ∈ L2 .
Proof. We give a proof for A0t . Regarding the other cases, it is proved similarly. First, fromthe definition of the conditional expectation, A0t and At are orthogonal. Thus, what we haveto prove is E[q(Jst )|Hst ] is included in At . This is proved as follows:
E[E[q(Jst )|Hst ]|Hat−1 ] = E[q(Jst )|Hat−1 ] = E[E[q(Jst )|Jrt−1 ]|Hat−1 ] = 0.
If we can prove that the influence function Eq. (16) is orthogonal to the orthogonaltangent space Eq. (17), we can see that the above derivative is actually the efficient influencefunction under the model M1 . This fact is shown as follows.
+1
Lemma 31. The derivative Eq. (16) is orthogonal to {A0t }Tt=0, {Bt00 }Tt=0 , {Ct0 }Tt=0

Proof. The influence function is orthogonal to A0k : for t(Jsk ) ∈ A0k
"(
E

πe

−ρ

+

T
X

(
λc (Hac )rc −λc (Hac )=E

Eπe [rt |Hac ] − λc−1 (Hac−1 )

t=c

c=0

"( T
X

T
X

λc (Hac )rc − λk−1

c=k

T
X

T
X

))
Eπe [rt |Hsc ]

t=c

)

#

Eπe [rt |Hsk ] t(Jsk )

t=k

= 0.
The influence function is orthogonal to Bk00 : for t(Jak ) ∈ Bk00 ;
"(
E

πe

−ρ

+

T
X

(
λc rc −λc=E(
λc rc −

λc

=E

c=k

T
X

)
λc rc

(
−

λk

T
X

Eπe [rt |Hac ] − λc−1TX

Eπe [rt |Hsc ]#t(Jak )))
Eπe [rt |Hsc ]

t=c
T
X

))

t=c

t=c

c=k

"(( T
X

Eπe [rt |Hac ] − λc−1

t=c

c=0

"( T
X

T
X

))
Eπe [rt |Hak ]

t=k

41

#
t(Jak )#t(Jak )#t(Jsk )

Kallus and Uehara

"((
=E

λk

T
X

)
Eπe [rt |Hak ]

(
−

λk

t=k

T
X

#

))
Eπe [rt |Hak ]t(Jak ) = 0.

t=kThe influence function is orthogonal to Ck0 : for t(Jrk ) ∈ Ck0 ;
"(
πe

−ρ

E

+

T
X

(
λc (Hac )rc −λc (Hac )=E

Eπe [rt |Hac ] − λc−1 (Hac−1 )#))Eπe [rt |Hsc ]t(Jrk )#)λc (Hac )rc

T
X
t=c

t=c

c=0

"( T
X

T
X

t(Jrk )

c=k

"((
=E

λk−1
"((

=E

λk−1

T
X
t=k
T
X

#

))
Eπe [rt |Jrk ]"(t(Jrk ) = E

λk−1

T
X

)

#

Eπe [rt |Hak , rk ] t(Jrk )t=k))Eπe [rt |Hak , rk ]#E[t(Jrk )|Hak , rk ] = 0.

t=kThis concludes the proof for M1 .
Efficient influence function under M1b . Next, we show that the efficiency bound is
still the same even if we know the target policy. To show that, we derive an orthogonal
space of the tangent space of the regular parametric submodel:{pθ (s0 )p(a0 |s0 )pθ (r0 |Ha0 )pθ (s1 |Hr0 )p(a1 |Hs1 )pθ (r1 |Ha1 ) · · · pθ (rT |HaT )},where p(at |Hst ) is fixed at πtb . This is equal to
M

(A0t

M

Bt00

M

Ct0 )

(18)

0≤t≤T

This space Eq. (18) is orthogonal to the obtained efficient influence function under M1 .
Therefore, the efficient influence function under M1b is the same as the one under Mb .
Efficiency bound. We use a law of total variance (Bowsher and Swain, 2012) to computethe variance of the efficient influence function.
" T
#
X
var
(λt rt − (λt qt − λt−1 vt ))

=

t=0
"
T
+1
X

"

E var E λt−1 rt−1 +

t=0

=

T
+1
X

"

"

E var E λt−1 rt−1 +

t=0

=

T
+1
X
t=0

"

"

E var E λt−1 rt−1 +

T
X
k=0
T
X

#

!#

(λk rk − {λk qk − λk−1 vk }) |Jat |Jat−1#!#(λk rk − {λk qk − λk−1 vk }) |Jat |Jat−1

k=t
T
X

!
λk rk

k=t

42

#
− {λt qt − λt−1 vt }|Jat |Jat−1#!

Double Reinforcement Learning in Markov Decision Processes

=

T
+1
X



E λ2t−1 var rt−1 + vt (Hst )|Hat−1 .

t=0Here, we used E[PT

k=t λk rk |Jak ] = λk qk .

Proof of Theorem 2. The proof is analogous to the proof of Theorem 1.
Proof of Theorem 3.
Efficient influence function under M2 .

The entire regular parametric submodel is{pθ (s0 )pθ (a0 |s0 )pθ (r0 |s0 , a0 )pθ (s1 |s0 , a0 )pθ (a1 |s1 )pθ (r1 |s1 , a1 ) · · · pθ (rT |sT , aT )}.
The score function of the parametric submodel isg(JsT ) =
=

T
X
k=0
T
X

∇θ log pθ (sk | sk−1 , ak−1 ) + ∇θ log pθ (ak+1 | sk ) + ∇θ log pθ (rk | sk , ak )gsk |sk−1 ,ak−1 +

k=0

T
X

gak+1 |sk +

k=0

T
X

grk |sk ,ak .

k=0

We first calculate a derivative of the target functional w.r.t the model M2 . Note that
this derivative is not only derivative. We haveTX∇θ E [
rt ]
πe

t=0

= ∇θ

Z X
T

rt

t=0

=

T
X

(T
Y

)
pθ (sk |ak−1 , sk−1 )pπke (ak |sk )pθ (rk |ak , sk ) dµ(JsT )t=0{Eπe [(Eπe [rc |s0 ] − Eπe [rc ])gs0 ] + Eπe [(rc − Eπe [rc |sc , ac ])grc |sc ,ac ]c=0"+ Eπ e

Eπe [

T
X

c=t+1
T
X

rt |sc+1 ] − Eπe [

T
X

!

#

rt |sc , ac ] gsc+1 |sc ,ac }

c=t+1




pπe (sc , ac )={E[(E[rc |s0 ] − Eπe [rc ])g] + E(rc − E[rc |sc , ac ])gpπb (sc , ac )
c=0
"
#
T
T
X
X
pπe (sc , ac )rt |sc+1 ] − E[
rt |sc , ac ])g }+E(E[pπb (sc , ac )
t=c+1
t=c+1
""
(
)#
#
T
T
T
X
pπe (sc , ac )pπe (sc , ac ) X
pπe (sc−1 , ac−1 ) X
πe= E −ρ +rc −Eπ [rt |sc , ac ] +Eπ [rt |sc ]g(JsT )
pπb (sc , ac )pπb (sc , ac ) t=c epπb (sc−1 , ac−1 ) t=c e
c=043

Kallus and Uehara

Therefore, the following function is a derivative:
)
(
T
T
T
X
pπce (sc−1 , ac−1 ) X
pπce (sc , ac ) X
pπce (sc , ac )
πce
−ρ +
rc −
Eπ [rt |sc , ac ] −Eπ [rt |sc ] .
pπcb (sc , ac )pπcb (sc , ac ) t=c epπcb (sc−1 , ac−1 ) t=c e
c=0

(19)
We will show this derivative is the efficient influence function.
In order to show this, we calculate the tangent space of model M2 . The tangent spaceof the model M2 is the product space;
M
M M
(At
Bt
Ct ),
0≤t≤T

At = {q(st , st−1 , at−1 ); E[q(st , st−1 , at−1 )|st−1 , at−1 ] = 0, q ∈ L2 },Bt = {q(at , st ); E[q(at , st )|st ] = 0, q ∈ L2 },Ct = {q(rt , st , at ); E[q(rt , st , at )|st , at ] = 0, q ∈ L2 }.

The orthogonal space of the tangent space is the product of
M M
M
(A0t
Bt0
Ct0 ),

(20)

0≤t≤T

such that
A0t

M

At = A00t , A00t = {q(Jst ); E[q(Jst )|Jrt−1 ] = 0, q ∈ L2 },

Bt0

M

Bt = Bt00 , Bt00 = {q(Jat ); E[q(Jat )|Jst ] = 0, q ∈ L2 },Ct0

M

Ct = Ct00 , Ct00 = {q(Jrt ); E[q(Jrt )|Jat ] = 0, q ∈ L2 }.

More specifically, the orthogonal tangent space is represented as
A0t = q(Jst ) − E[q(Jst )|st , at−1 , st−1 ]; E[q(Jst )|Jrt−1 ] = 0, q ∈ L2 ,
Bt0 = q(Jat ) − E[q(Jrt )|st , at ]; E[q(Jat )|Jst ] = 0, q ∈ L2 ,

Ct0 = q(Jrt ) − E[q(Jrt )|rt , st , at ]; E[q(Jrt )|Jat ] = 0, q ∈ L2 .

If we can prove that the derivative Eq. (19) is orthogonal to the orthogonal tangent space
Eq. (20), we can see that the above derivative is actually the efficient influence function
under the model M2 . This fact is shown as follows.
+1
Lemma 32. The derivative Eq. (19) is orthogonal to {A0t }Tt=0, {Bt00 }Tt=0 , {Ct0 }Tt=0 .

Proof. First, the influence function Eq. (19) is orthogonal to A0k ; for t(Jsk ) ∈ A0k
"(
)
#
T
X
E
v0 +
µt (st , at )(rt + vt+1 − qt ) t(Jsk )t=044

Double Reinforcement Learning in Markov Decision Processes

"(
=E

T
X

#

)

µt (st , at )(rt + vt+1 − qt ) t(Jsk )t=k−1= E [µk−1 (sk−1 , ak−1 )(rk−1 + vk − qk−1 )t(Jsk )]= E [µk−1 (sk−1 , ak−1 )vk t(Jsk )]= E [µk−1 (sk−1 , ak−1 )vk E[t(Jsk )|sk , ak−1 , sk−1 ]] = 0
Second, the influence function Eq. (19) is orthogonal to Bk00 ; for t(Jak ) ∈ Bk00
"(
E

v0 +

=E

T
X

#

)
µt (st , at )(rt + vt+1 − qt ) t(Jak )

t=0
"( T
X

#

)

µt (st , at )(rt + vt+1 − qt ) t(Jak ) = 0.

t=kThird, the influence function Eq. (19) is orthogonal to Ck0 ; for t(Jrk ) ∈ Ck0
"(
E

=E

v0 +

T
X

)

#

µt (st , at )(rt + vt+1 − qt ) t(Jrk )

t=0
"( T
X

)

#

µt (st , at )(rt + vt+1 − qt ) t(Jrk )t=k= E [{µk (sk , ak )(rk + vk+1 − qk )} t(Jrk )]= E [{µk (sk , ak )(E[rk + E[vk+1 |Jr,k ] − qk )} t(Jrk )]= E [{µk (sk , ak )(rk + E[vk+1 |sk ] − qk )} E[t(Jrk )|sk , ak , rk ]] = 0.
Efficient influence function under M2b . In Theorem 32, we check that the φM2
eff isorthogonal to B”t . This concludes the proof noting that the orthogonal tangent space ofM2b is
M
M
M
(A0t
B”t
Ct0 ).
0≤t≤T

Efficiency bound. To show an efficiency bound, we use a law of total variance (Bowsher
and Swain, 2012). Recall that we can also easily derive this variance form using anotherequivalent form of efficient influence function.
"#
TX
var v0 +µt (st , at )(rt + vt+1 − qt )
t=0

=

=

T
+1
X
t=0
T
+1
X
t=0

"

" "

E var E v0 +

T
X

#

##

µk (sk , ak )(rk + vk+1 − qk )|Jat |Jat−1

t=0

"

" "

E var ETX#
µk (sk , ak )(rk + vk+1 − qk )|Jat |Jat−1

k=t−1

45

##

Kallus and Uehara

=

=

=

T
+1
X
t=0
T
+1
X
t=0
T
+1
X




E var E [µt−1 (st−1 , at−1 )(rt−1 + vt − qt−1 )|Jat ] |Jat−1
E var µt−1 (st−1 , at−1 )(rt−1 + vt − qt−1 )|Jat−1E µ2t−1 (st−1 , at−1 )var (rt−1 + vt (st ))|Jat−1 .

t=0Proof of Theorem 4. From Jensen’s inequality,
T
+1
X

+1
 2
 TX


E E(λ2t−1 |st−1 , at−1 )var {rt + vt (st )|st−1 , at−1 }E λt−1 var {rt + vt (st )|st−1 , at−1 } =

t=0
T
+1
X

≥

t=0
+1

 TX


E E(λt−1 |st−1 , at−1 )2 var {rt + vt (st )|st−1 , at−1 } =E µ2t−1 var {rt + vt (st )|st−1 , at−1 } .

t=0

t=0

When λ2t−1 is not constant given st−1 , at−1 and var {rt + vt (st )|st−1 , at−1 } =6 0, the inequality
is strict.
Proof of Theorem 5. By changing the limits of summation and letting r−1 = 0, λ0 = 1, we
can write the efficiency bound under NMDP as
T
+1
X

T
+1
X






E λ2t−1 var rt−1 + vt (Hst ) | Hat−1 ≤ C T +1
E λt−1 var rt−1 + vt (Hst )|Hat−1

t=0= C T +1
= C T +1

t=0
T
+1
X
t=0
T
+1
X




Eπe varπe rt−1 + vt (Hst )|Hat−1
Eπe var rt−1 + vt (Hst )|Hat−1

t=0= C T +1 var[T+1
X

rt−1 ]

t=02≤ C T +1 (T + 1)2 Rmax.

The last equality follows by the law of total variance.
Similarly, the efficiency bound under MDP is
T
+1
X

T
+1
X
 2

0
E µt−1 var {rt−1 + vt (st ) | st−1 , at−1 } ≤ C
E [µt−1 var {rt−1 + vt (st ) | st−1 , at−1 }]

t=0

= C0

t=0
T
+1
X
t=0

46

Eπe [var {rt−1 + vt (st ) | st−1 , at−1 }]Double Reinforcement Learning in Markov Decision Processes

= C0

T
+1
X

Eπe [varπe {rt−1 + vt (st ) | st−1 , at−1 }]t=00= C var[T+1X

rt−1 ]

t=0
2≤ C 0 (T + 1)2 Rmax.

The last equality again follows by the law of total variance.
Finally, for the NMDP lower bound we have by Jensen’s inequality
T
+1
X

+1

 TX


 2

Eπe λt−1 var rt−1 + vt (Hst ) | Hat−1
E λt−1 var rt−1 + vt (Hst ) | Hat−1 =

t=0

t=0
T
+1
X

≥

≥

≥

t=0
T
+1
X
t=0
T
+1
X




exp Eπe log(λt−1 var rt−1 + vt (Hst ) | Hat−1 )exp(Eπe [log(λt−1 )] + Eπe var rt−1 + vt (Hst ) | Hat−1 ) )2
exp(t log Cmin + log Vmin
)

t=0
T +1
2
.
≥ VminCmin

Proof of Theorem 6. Without loss of generality, we consider the case K = 2.
φ({λ̂k }, {q̂k }) as:
T
X

Defineλ̂k rk − {λ̂k q̂k − λ̂k−1 Eπe [q̂k (Hak )|Hsk ]}.

k=0

The estimator ρ̂M1
DRL(M1 ) is given by
n1
n2
(1)
(1)
(2)
(2)
Pn1 φ({λ̂k }, {q̂k }) + Pn2 φ({λ̂k }, {q̂k }),nnwhere Pn1 is an empirical approximation based on a set of samples such that i ∈ D1 , Pn2 is
an empirical approximation based on a set of samples such that i ∈ D2 . Then, we have
p
√
e
(1)
(1)
(1)
(1)
n(Pn1 φ({λ̂k }, {q̂k }) − ρπ ) = n/n1 Gn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]
(21)
p
(1)
(1)
+ n/n1 Gn1 [φ({λk }, {qk })]
(22)
√
e
(1)
(1)
(1)
(1)
+ n(E[φ({λ̂k }, {q̂k })|{λ̂k }, {q̂k }] − ρπ ).
(23)
We analyze each term. To do that, we use the following relation:φ({λ̂k }, {q̂k }) − φ({λk }, {qk }) = D1 + D2 + D3 ,
47

where

Kallus and Uehara

D1 =

T
X
(λ̂k − λk )(−q̂k + qk ) + (λ̂k−1 − λk−1 )(v̂k − vk ),
k=0

D2 =

D3 =

T
X

λk (−q̂k + qk ) + λk−1 (v̂k − vk ),

k=0
T
X

(λ̂k − λk )(rk − qk + vk+1 ).

k=0First, we show the term Eq. (21) is op (1).
Lemma 33. The term Eq. (21) is op (1).
Proof. If we can show that for any  > 0,
lim

√

n→∞

(1)

(1)

(1)

(1)

n1 P [Pn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]

(24)
(1)

(1)

− E[φ({λ̂k }, {q̂k }) − φ({λk }, {qk })|{λ̂k }, {q̂k }] > |D2 ] = 0,
Then, by bounded convergence theorem, we would have
lim

√

n→∞

(1)

(1)

(1)

(1)

n1 P [Pn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]
(1)

(1)

− E[φ({λ̂k }, {q̂k }) − φ({λk }, {qk })|{λ̂k }, {q̂k }] > ] = 0,
yielding the statement.
To show Eq. (24), we show that the conditional mean is 0 and conditional variance isop (1). The conditional mean is
(1)

(1)

(1)

(1)

E[Pn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })|{λ̂k }, {q̂k }]
(1)

(1)

− P[φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]|D2 ] = 0.
(1)

(1)

Here, we leveraged the sample splitting construction, that is, λ̂k and q̂k
D2 . The conditional variance is
√
(1)
(1)
var[ n1 Pn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]|D2 ]only depend on

(1)

(1)

= E[E[D12 + D22 + D32 + 2D1 D2 + 2D2 D3 + 2D2 D3 |{q̂k }, {λk }]|D2 ]= op (1).
(1)(1)

Here, we used the convergence rate assumption and the relation kv̂k − vk k2 ≤ kq̂k − qk k2
arising from the fact that the former is the marginalization of the latter over πke . Then, from
Chebyshev’s inequality:
√

(1)

(1)

(1)

(1)

n1 P [Pn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })] − E[φ({λ̂k }, {q̂k })
(1)

(1)

− φ({λk }, {qk })|{λ̂k }, {q̂k }] > |D2 ]
≤

√
1
(1)
(1)
var[ n1 Pn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]|D2 ] = op (1).
2
48

Double Reinforcement Learning in Markov Decision Processes

Lemma 34. The term Eq. (23) is op (1).
Proof.
√

(1)

(1)

(1)

(1)

nE[φ({λ̂k }, {q̂k }) − E[φ({λk }, {qk })]|{λ̂k }, {q̂k }]

=

+

√

√

T
X
(1)
(1)
(1)
(1)
(1)
(1)
nE[ (λ̂k − λk )(−q̂k + qk ) + (λ̂k−1 − λk−1 )(−v̂k + vk )|{λ̂k }, {q̂k }]
k=0
T
X

nE[

(1)

(1)

(1)

(1)

λk (−q̂k + qk ) + λk−1 (v̂k − vk )|{λ̂k }, {q̂k }]

k=0
T
X
√
(1)
(1)
(1)
+ nE[ (λ̂k − λk )(rk − qk + vk+1 )|{λ̂k }, {q̂k }]

=

=

k=0
T
X

√
nE[
√

n

k=0
T
X

(1)

(1)

(1)

(1)

(1)

(λ̂k − λk )(−q̂k + qk ) + (λ̂k−1 − λk−1 )(−v̂k + vk )|{λ̂k }, {q̂k }]
(1)

(1)

O(kλ̂k − λk k2 kq̂k − qk k2 ) =

k=0

T
√ X
n
op (n−1/2 ) = op (1).
k=0

Finally, we get
p
√
e
(1)
(1)
n(Pn1 φ({λ̂k }, {q̂k }) − ρπ ) = n/n1 Gn1 [φ({λk }, {qk })] + op (1).
Therefore,√een(ρ̂πDRL(M1 ) − ρπ )
√
√
e
e
(1)
(1)
(2)
(2)
= n1 /n nφ({λ̂k }, {q̂k }) − ρπ ) + n2 /n n(Pn2 φ({λ̂k }, {q̂k }) − ρπ )pp
= n1 /nGn1 [φ({λk }, {qk })] + n2 /nGn2 [φ({λk }, {qk })] + op (1)= Gn [φ({λk }, {qk })] + op (1),econcluding the proof by showing the influence function of ρ̂πDRL(M1 ) is the efficient one.
Proof of Theorem 7. Here, a / b means there exists constant Q not depending on T, Rmax , C, n, C1 , C2such that a < Qb.
Define φ({λ̂k }, {q̂k }) as:TX

λ̂k rk − {λ̂k q̂k − λ̂k−1 Eπe [q̂k (Hak )|Hsk ]}.

k=0

The estimator ρ̂M1
DRL(M1 ) is given by
n1
n2
(1)
(1)
(2)
(2)
Pn φ({λ̂k }, {q̂k }) + Pn2 φ({λ̂k }, {q̂k }),n 1nwhere Pn1 is an empirical approximation based on a set of samples such that i ∈ D1 , Pn2 is
an empirical approximation based on a set of samples such that i ∈ D2 .
49

Kallus and Uehara

Then, we have
(1)

e

(1)

p
(1)
(1)
1/n1 Gn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]p+ 1/n1 Gn1 [φ({λk }, {qk })](Pn1 φ({λ̂k }, {q̂k }) − ρπ ) =

(25)
(26)

e
(1)
(1)
(1)
(1)
+ (E[φ({λ̂k }, {q̂k })|{λ̂k }, {q̂k }] − ρπ ).

(27)

We analyze each term. To do that, we use the following relation:φ({λ̂k }, {q̂k }) − φ({λk }, {qk }) = D1 + D2 + D3 ,
D1 =

D2 =

where

T
X

(λ̂k − λk )(−q̂k + qk ) + (λ̂k−1 − λk−1 )(v̂k − vk ),

k=0
T
X

λk (−q̂k + qk ) + λk−1 (v̂k − vk ),k=0TX(λ̂k − λk )(rk − qk + vk+1 ).
D3 =
k=0

Lemma 35. With probability 1 − 2δ, the absolute value of the term Eq. (25) is bounded bys√2log(2/δ)T 2 (T Rmax C T +1 κ1 κ2 + κ1 T 2 Rmax
+ κ2 C 2(T +1) ) log(2/δ)T Rmax C T +1+,nn
up to some constant independent of n, C, Rmax , T .
Proof. From Bernstein inequality, with probability 1 − δ,
(1)

(1)

(1)

(1)

P [|Pn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]− E[φ({λ̂k }, {q̂k }) − φ({λk }, {qk })|D1 ]| > ]|D2 ]≤ 2 exp −0.5n2(1)(1)

E[{φ({λ̂k }, {q̂k }) − φ({λk }, {qk })}2 |D2 ] + Q1 T Rmax C T +1 noting the conditional mean is 0;
(1)

(1)

E[Pn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })|D1 ]
(1)

(1)

− P[φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]|D2 ] = 0.
Here, Q1 is some constant independent of n, C, Rmax , T, δ.
With probability 1 − δ, the conditional variance is
(1)

(1)

E[{φ({λ̂k }, {q̂k }) − φ({λk }, {qk })}2 |D2 ]= E[D12 + D22 + D32 + 2D1 D2 + 2D2 D3 + 2D2 D3 |D2 ]√/ T 2 (T Rmax C T +1 κ1 κ2 + κ1 (T Rmax )2 + κ2 C 2(T +1) ).
50!
,

Double Reinforcement Learning in Markov Decision Processes

Therefore, with probability 1 − 2δ,
(1)

(1)

(1)

(1)

P [|Pn1 [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]− E[φ({λ̂k }, {q̂k }) − φ({λk }, {qk })|D1 ]| > ]|D2 ]
0.5n1 2
≤ 2 exp −.
√2
Q2 T 2 (T Rmax C T +1 κ1 κ2 + κ1 T 2 Rmax
+ κ2 C 2(T +1) ) + Q1 T Rmax C T +1 Here, Q2 is some constant independent of n, C, Rmax , T, δ. Then, by the law of total
probability and some algebra, the statement is concluded.
Lemma 36. With probability 1 − δ, the absolute value of the term Eq. (27) is bounded by√T κ1 κ2 .
up to some constant independent of n, C, Rmax , T, δ.
Proof. Here, we have
(1)
(1)
(1)
(1)
E[φ({λ̂k }, {q̂k }) − E[φ({λk }, {qk })]|{λ̂k }, {q̂k }] ≤

T
X

(1)

(1)

2kλ̂k − λk k2 kq̂k − qk k2 .

k=0√Then, with probability 1 − δ, this is bounded by T κ1 κ2 .
Combining all results so far, with probability 1 − 6δ, we have
n1
n2
(1)
(1)
(2)
(2)
Pn1 φ({λ̂k }, {q̂k }) + Pn2 φ({λ̂k }, {q̂k }) − ρπe
n
n
s
√
2
log(2/δ)T 2 (T Rmax C T +1 κ1 κ2 + κ1 T 2 Rmax
+ κ2 C 2(T +1) ) Q1 log(2/δ)T Rmax C T +1√+
< Q3 T κ1 κ2 + Q2
nn
+ Pn1 [φ({λk }, {qk })] + Pn2 [φ({λk }, {qk })] − ρπe .
Here, Q3 is some constant independent of n, C, Rmax , T . Noting |Pn [φ({λk }, {qk })] − ρπe | isbounded byr
2 log(2/δ)Effbd(M1 ) Q1 log(2/δ)T Rmax C T +1+nnwith probability 1 − δ, the statement is concluded.
Proof of Theorem 9. We define φ({λ̂k }, {q̂k }) as:TX

λ̂k rk − λ̂k−1 {η̂k q̂k − Eπe [q̂k (Hak )|Hsk ]}.

k=0e

The estimator ρ̂πDR is given by Pn φ({λ̂k }, {q̂k }). Then, we have√en(Pn φ({λ̂k }, {q̂k }) − ρπ ) = Gn [φ({λ̂k }, {q̂k }) − φ({λk }, {qk })]
51

(28)

Kallus and Uehara+ Gn [φ({λk }, {qk })]√e+ n(E[φ({λ̂k }, {q̂k })|{λ̂k }, {q̂k }] − ρπ ).

(29)
(30)

If we can prove that the term Eq. (28) is op (1), the statement is concluded as in the proof of
Theorem 6. We proceed to prove this.
First, we show φ({λ̂k }, {q̂t }) − φ({λt }, {qt }) belongs to a Donsker class. The transformation({λk }, {qk }) 7→TX

λk rk − {λk qk − λk−1 Eπe [qk (Hak )|Hsk ]}k=0

is a Lipschitz function. Therefore, by Example 19.20 in van der Vaart (1998), φ({λ̂k }, {q̂k })−φ({λk }, {qk }) is an also Donsker class. In addition, we can also show thatkφ({λ̂k }, {q̂k }) − φ({λk }, {qk })k2 = op (1),as in Lemma 33. Therefore, from Lemma 19.24 in van der Vaart (1998), the term Eq. (28)is op (1), concluding the proof.
Proof of Theorem 10. Without loss of generality, we consider the case K = 2.
We use the following doubly robust structure" T#X
λk rk − {λk qk − λk−1 Eπe (qk |Hsk )}E
k=0= E{Eπe (q0 |s0 )} + E

" T
X

#
e

λk {rk − qk + Eπe (qk |Hsk+1 )} = ρπ .

k=0Then, as in the proof of Theorem 6,
√

(1)

(1)

e

n(Pn1 φ({λ̂k }, {q̂k }) − ρπ )
p
p
(1)
(1)
= n/n1 Gn1 [φ({λ̂k }, {q̂k }) − φ({λ†k }, {qk† })] + n/n1 Gn1 [φ({λ†k }, {qk† })]
p
√
e
(1)
(1)
(1)
(1)
+ n/n1 (E[φ({λ̂k }, {q̂k })|{λ̂k }, {q̂k }] − E[φ({λ†k }, {qk† })]) + n(E[φ({λ†k }, {qk† })] − ρπ )p√e
= n/n1 Gn1 [φ({λ†k }, {qk† })] + n(E[φ({λ†k }, {qk† })] − ρπ ) + op (1).

Here, we used
p
(1)
(1)
n/n1 Gn1 [φ({λ̂k }, {q̂k }) − φ({λ†k }, {qk† })] = op (1)from Theorem 33 and
p
(1)
(1)
(1)
(1)
n/n1 (E[φ({λ̂k }, {q̂k })|{λ̂k }, {q̂k }] − E[φ({λ†k }, {qk† })]) = op (1),which we prove below as in Theorem 34.
Lemma 37.
p
(1)
(1)
(1)
(1)
n/n1 (E[φ({λ̂k }, {q̂k })|{λ̂k }, {q̂k }] − E[φ({λ†k }, {qk† })]) = op (1).
52

Double Reinforcement Learning in Markov Decision Processes

Proof. First, consider the case where λk = λ†k .
√nE[φ({λ̂k }, {q̂k }) − E[φ({λk }, {qk† })]|{λ̂k }, {q̂k }]
(1)

=

√

(1)

(1)

(1)

T
X
(1)
(1)
(1)
(1)
(1)
nE[ (λ̂k − λk )(−q̂k + qk† ) + (λ̂k−1 − λk−1 )(−v̂k + vk† )|{λ̂k }, {q̂k }]
k=0

T
X
√
(1)
(1)
(1)
+ nE[
λk (q̂k − qk† ) + λk−1 (v̂k − vk† )|{λ̂k }, {q̂k }]

+

=

√

k=0
T
X

nE[

√

†
(λ̂k − λk )(rk − qk† + vk+1)|{λ̂k }, {q̂k }]

k=0
T
X

nE[

(1)

(1)

(1)

(λ̂k − λk )(−q̂k + qk† ) + (λ̂k−1 − λk−1 )(−v̂k + vk† )|{λ̂k }, {q̂k }]
(1)

(1)

(1)

(1)

(1)

k=0
T
X
√
(1)
(1)
(1)
†
)|{λ̂k }, {q̂k }]+ nE[ (λ̂k − λk )(rk − qk† + vk+1

=

=

√

√

n

n

k=0
T
X

O(kλ̂k − λk k2 kq̂k − qk† k2 + kλ̂k − λk k2 )

k=0
T
X

(1)

(1)

(1)


Op (n−α1 )Op (n−α2 ) + Op (n−α1 ) = Op (1).

k=0

Next, consider the case where qk = qk† :√nE[φ({λ̂k }, {q̂k }) − E[φ({λ†k }, {qk })]|{λ̂k }, {q̂k }]

=

+

=

(1)

√

√

(1)

T
X
(1)
(1)
(1)
(1)
(1)
nE[ (λ̂k − λ†k )(−q̂k + qk ) + (λ̂k−1 − λ†k−1 )(−v̂k + vk )|{λ̂k }, {q̂k }]
k=0
T
X

nE[

√

(1)

(1)

n

k=0
T
X

λ†k (q̂k − qk ) + λ†k−1 (v̂k − vk )|{λ̂k }, {q̂k }]
(1)

(1)

(1)

O(kλ̂k − λ†k k2 kq̂k − qk k2 + kq̂k − qk k2 )
(1)

(1)

(1)

k=0
T

√ X
= n
Op (n−α1 )Op (n−α2 ) + Op (n−α2 ) = Op (1).
k=0Using the above result, we prove the statement for each case below.
λ-model is well-specified. First, consider the case when λ†k = λk :E[φ({λk }, {qk† })] = E[TX

[λk rk − {λk qk† (Hak ) − λk−1 Eπe [qk† (Hak )|sk ]]

k=0

53

Kallus and Uehara

T
X
e
= E[
λk rk ] = ρπ .
k=0

Then,
√

(1)

(1)

e

n(Pn1 φ({λ̂k }, {q̂k }) − ρπ ) =pn/n1 Gn1 [φ({λk }, {qk† })] + Op (1).

Therefore,
p
p
√
e
e
n(ρ̂πDRL(M1 ) − ρπ ) = n1 /nGn1 [φ({λk }, {qk† })] + n2 /nGn2 [φ({λk }, {qk† })] + Op (1)= Gn [φ({λk }, {qk† })] + Op (1) = Op (1),√eewhich shows ρ̂πDRL(M1 ) is n-consistent around ρπ when the model for the behavior policy
is well-specified.
q-model is well-specified. Next, consider the case where qk† = qk .
"E[φ({λ†k }, {qk })] = E

Eπe [q0 (Ha0 )|s0 ] +TX#λ†k {rk − qk (Hak ) + Eπe [qk (Hak )|sk+1 ]}k=0e= E [Eπe [q0 (Ha0 )|s0 ]] = ρπ .
Then,
√

(1)

(1)

e

n(Pn1 φ({λ̂k }, {q̂k }) − ρπ ) =pn/n1 Gn1 [φ({λ†k }, {qk })] + Op (1).

Therefore,
p
p
√
e
e
n(ρ̂πDRL(M1 ) − ρπ ) = n1 /nGn1 [φ({λ†k }, {qk })] + n2 /nGn2 [φ({λ†k }, {qk })] + Op (1)= Gn [φ({λ†k }, {qk })] + Op (1) = Op (1).
√eewhich shows ρ̂πDRL(M1 ) is n-consistent around ρπ when the model for the q-function iswell-specified.
Proof of Theorem 11. Without loss of generality, we consider the case K = 2.
Then, as in the proof of Theorem 6,
(1)

(1)

e

Pn1 φ({λ̂k }, {q̂k }) − ρπ
p
p
(1)
(1)
= 1/n1 Gn1 [φ({λ̂k }, {q̂k }) − φ({λ†k }, {qk† })] + 1/n1 Gn1 [φ({λ†k }, {qk† })]
p
e
(1)
(1)
(1)
(1)
+ 1/n1 (E[φ({λ̂k }, {q̂k })|{λ̂k }, {q̂k }] − E[φ({λ†k }, {qk† })]) + (E[φ({λ†k }, {qk† })] − ρπ )pe
= 1/n1 Gn1 [φ({λ†k }, {qk† })] + (E[φ({λ†k }, {qk† })] − ρπ ) + op (1)= (Pn1 − P)[φ({λ†k }, {qk† })] + op (1).
Here, under the assumption, we use the following equations:
p
(1)
(1)
1/n1 Gn1 [φ({λ̂k }, {q̂k }) − φ({λ†k }, {qk† })] = op (1)54Double Reinforcement Learning in Markov Decision Processes(E[φ({λ̂k }, {q̂k })|{λ̂k }, {q̂k }] − E[φ({λ†k }, {qk† })]) = op (1)
(1)

(1)

(1)

(1)

e

(E[φ({λ†k }, {qk† })] − ρπ ) = 0.
These equations are proved as in the proof of Theorem 10. Then,ee(ρ̂πDRL(M1 ) − ρπ ) = (Pn1 − P)[φ({λ†k }, {qk† })] + (Pn2 − P)[φ({λ†k }, {qk† })] + op (1)= (Pn − P)[φ({λ†k }, {qk† })] + op (1).
From the law of large numbers, the statement is concluded.
Proof of Theorem 12. Let any π̂tb be given. Due to boundedness away from zero, we have
k
Y
πe

T
i
k
i−1
k
X
Y
πte Y πte Y πte Y πte
t
≤
−
λ
−
k
π̂ b
π̂ b
π b t=0 π̂tb t=i πtb
0=t t
i=0 t=0 t t=i t
2kX≤
op (k1/π̂tb − 1/πtb k2 )
t=0
k
X
≤
op (kπ̂tb − πtb k2 )t=0
≤ op (n−α ).

!
2Proof of Theorem 13. Without loss of generality, we consider the case K = 2.
Define φ({µ̂k }, {q̂k }) as:TXµ̂k {rk − q̂k } − µ̂k−1 Eπe [q̂k (Hak )|Hsk ]}.

k=0

The estimator ρ̂DRL(M2 ) is given by
n1
n2
(1)
(1)
(2)
(2)
Pn1 φ({µ̂k }, {q̂k }) + Pn2 φ({µ̂k }, {q̂k }).
nnThen, we have
√

(1)

e

(1)

p
(1)
(1)
n/n1 Gn1 [φ({µ̂k }, {q̂k }) − φ({µk }, {qk })]p+ n/n1 Gn1 [φ({µk }, {qk })]
√
e
(1)
(1)
(1)
(1)
+ n(E[φ({µ̂k }, {q̂k })|{µ̂k }, {q̂k }] − ρπ ).

n(Pn1 φ({µ̂k }, {q̂k }) − ρπ ) =We analyze each term. To do that, we use the following relation;
(1)

(1)

φ({µ̂k }, {q̂k }) − φ({µk }, {qk }) = D1 + D2 + D3 ,
55

where

(31)
(32)
(33)

Kallus and Uehara

D1 =

T
X

(1)

(1)

(1)

(µ̂k − µk )(−q̂k + qk ) + (µ̂k−1 − µk−1 )(v̂k − vk ),

k=0

D2 =

D3 =

T
X
k=0
T
X

(1)

(1)

µk (−q̂k + qk ) + µk−1 (v̂k − vk ),(1)(µ̂k − µk )(rk − qk + vk+1 ).

k=0First, we show the term Eq. (31) is op (1).
Lemma 38. The term Eq. (31) is op (1).
Proof. If we can show that for any  > 0,
lim

√

n→∞

(1)

(1)

(1)

(1)

(1)

(1)

n1 P [Pn1 [φ({µ̂k }, {q̂k }) − φ({µk }, {qk })]

(34)
(1)

(1)

− E[φ({µ̂k }, {q̂k }) − φ({µk }, {qk })|{µ̂k }, {q̂k }] > |D2 ] = 0.
Then, by bounded convergence theorem, we have
lim

n→∞

√

(1)

(1)

(1)

(1)

n1 P [Pn1 [φ({µ̂k }, {q̂k }) − φ({µk }, {qk })]
(1)

(1)

− E[φ({µ̂k }, {q̂k }) − φ({µk }, {qk })|{µ̂k }, {q̂k }] > ] = 0,
yielding the statement.
To show Eq. (34), we show that the conditional mean is 0 and conditional variance isop (1). The conditional mean is
(1)

(1)

(1)

(1)

E[Pn1 [φ({µ̂k }, {q̂k }) − φ({µk }, {qk })|{µ̂k }, {q̂k }]−
(1)

(1)

P[φ({µ̂k }, {q̂k }) − φ({µk }, {qk })]|D2 ] = 0.
(1)

(1)

Here, we used a sample splitting construction, that is, µ̂k and q̂k only depend on D2 . The
conditional variance is
√
(1)
(1)
var[ n1 Pn1 [φ({µ̂k }, {q̂k }) − φ({µk }, {qk })]|D2 ]
(1)

(1)

=

E[E[D12 + D22 + D32 + 2D1 D2 + 2D2 D3 + 2D2 D3 |{q̂k }, {µk }]|D2 ]=op (1).
(1)(1)

Here, we used the convergence rate assumption and the relation kv̂k − vk k2 ≤ kq̂k − qk k2 .
Then, from Chebyshev’s inequality;
√

(1)

(1)

(1)

(1)

(1)

(1)

n1 P [Pn1 [φ({µ̂k }, {q̂k }) − φ({µk }, {qk })] − E[φ({µ̂k }, {q̂k })−φ({µk }, {qk })|{µ̂k }, {q̂k }] > |D2 ]
√
1
(1)
(1)
≤ 2 var[ n1 Pn1 [φ({µ̂k }, {q̂k }) − φ({µk }, {qk })]|D2 ] = op (1).
56

Double Reinforcement Learning in Markov Decision Processes

Lemma 39. The term Eq. (33) is op (1).
Proof.
√

(1)

(1)

(1)

(1)

nE[φ({µ̂k }, {q̂k }) − E[φ({µk }, {qk })]|{µ̂k }, {q̂k }]

=

+

+

=

=

T
X
√
(1)
(1)
(1)
(1)
(1)
(1)
nE[ (µ̂k − µk )(−q̂k + qk ) + (µ̂k−1 − µk−1 )(v̂k − vk )|{µ̂k }, {q̂k }]

√

√

k=0
T
X

nE[

k=0
T
X

√

k=0
T
X

nE[

n

(1)

(1)

(1)

(1)

(1)

(1)

(µ̂k − µk )(rk − qk + vk+1 )|{µ̂k }, {q̂k }]

nE[

√

(1)

µk (−q̂k + qk ) + µk−1 (v̂k − vk )|{µ̂k }, {q̂k }]

k=0
T
X

(1)

(1)

(1)

(1)

(1)

(1)

(µ̂k − µk )(−q̂k + qk ) + (µ̂k−1 − µk−1 )(−v̂k + vk )|{µ̂k }, {q̂k }]
(1)

(1)

O(kµ̂k − µk k2 kq̂k − qk k2 ) =T√ Xn
op (n−1/2 ) = op (1).
t=0

k=0

Finally, we get
√

(1)

e

(1)

n(Pn1 φ({µ̂k }, {q̂k }) − ρπ ) =pn/n1 Gn1 [φ({µk }, {qk })] + op (1).

Therefore,√een(ρ̂πDRL(M2 ) − ρπ )
√
√
e
e
(1)
(1)
(2)
(2)
= n1 /n nφ({µ̂k }, {q̂k }) − ρπ ) + n2 /n n(Pn2 φ({µ̂k }, {q̂k }) − ρπ )pp
= n1 /nGn1 [φ({µk }, {qk })] + n2 /nGn2 [φ({µk }, {qk })] + op (1)= Gn [φ({µk }, {qk })] + op (1),econcluding the proof by showing the influence function of ρ̂πDRL(M2 ) is the efficient one.
Proof of Theorem 14. Almost same as the proof of Theorem 7. The differences are replacing
λt with µt , and C T +1 with C 0
Proof of Theorem 15. We define φ({µ̂k }, {q̂k }) as:TXµ̂k {rk − q̂k } − µ̂k−1 Eπe [q̂k (Hak )|Hsk ]}.

k=0e

The estimator ρ̂πDRL(M2 ), adaptive is given by Pn φ({µ̂k }, {q̂k }). Then, we have√en(Pn φ({µ̂k }, {q̂k }) − ρπ ) = Gn [φ({µ̂k }, {q̂k }) − φ({µk }, {qk })]+ Gn [φ({µk }, {qk })]
57

(35)
(36)

Kallus and Uehara

+

√

e

n(E[φ({µ̂k }, {q̂k })|{µ̂k }, {q̂k }] − ρπ ).

(37)

If we can prove that the term Eq. (35) is op (1), the statement is concluded as in the proof of
Theorem 13. We proceed to prove this.
First, we show φ({µ̂k }, {q̂t }) − φ({µt }, {qt }) belongs to a Donsker class. The transformation({µk }, {qk }) 7→

TXµk rk − {µk qk − µk−1 Eπe [qk (Hak )|Hsk ]}k=0

is a Lipschitz function. Therefore, by Example 19.20 in van der Vaart (1998), φ({µ̂k }, {q̂k })−φ({µk }, {qk }) is an also Donsker class. In addition, we can also show thatkφ({µ̂k }, {q̂k }) − φ({µk }, {qk })k2 = op (1),as in Lemma 33. Therefore, from Lemma 19.24 in van der Vaart (1998), the term Eq. (35)is op (1), concluding the proof.
Proof of Theorem 16. Without loss of generality, we consider the case K = 2.
We use the following doubly robust structure
E

" T
X

#
µk rk − {µk qk − µk−1 Eπe (qk |sk )}k=0= E[Eπe (q0 |s0 )]] + E

" T
X

#
e

µk {rk − qk + Eπe (qk |sk+1 )} = ρπ .

k=0Then, as in the proof of Theorem 6,
√

(1)

e

(1)

n(Pn1 φ({µ̂k }, {q̂k }) − ρπ )
p
p
(1)
(1)
= n/n1 Gn1 [φ({µ̂k }, {q̂k }) − φ({µ†k }, {qk† })] + n/n1 Gn1 [φ({µ†k }, {qk† })]
p
(1)
(1)
(1)
(1)
+ n/n1 (E[φ({µ̂k }, {q̂k }); {µ̂k }, {q̂k }]√e− E[φ({µ†k }, {qk† })]) + n(E[φ({µ†k }, {qk† })] − ρπ )p√e
= n/n1 Gn1 [φ({µ†k }, {qk† })] + n(E[φ({µ†k }, {qk† })] − ρπ ) + Op (1).

We proceed by considering each case.
µ-model is well-specified. First, consider the case when µ†k = µk :E[φ({µk }, {qk† })] = E[TX[µk rk − {µk qk† (sk , ak ) − µk−1 Eπe [qk† (sk , ak )|sk ]]
k=0

T
X

= E[

e

µk rk ] = ρπ .

k=0

58

Double Reinforcement Learning in Markov Decision Processes

Then,
√

(1)

e

(1)

n(Pn1 φ({µ̂k }, {q̂k }) − ρπ ) =pn/n1 Gn1 [φ({µk }, {qk† })] + Op (1).

Therefore,
p
p
√
†
πe
πe
[φ({µ
},
{q
})]
+
n(ρ̂DRL(M
−
ρ
)
=
n
/nG
n2 /nGn2 [φ({µk }, {qk† })] + Op (1)
1
n
k
1
k
2)
= Gn [φ({µk }, {qk† })] + Op (1) = Op (1),√eewhich shows ρ̂πDRL(M2 ) is n-consistent around ρπ when the model for the µ-function iswell-specified.
q-model is well-specified. Next, consider the case where qk† = qk :#"TX
E[φ({µ†k }, {qk })] = E Eπe [q(sk , ak )|s0 ] +µ†k {rk − qk (sk , ak ) + Eπe [qk (sk , ak )|sk+1 ]}k=0πe= E [Eπe [q0 (s0 , a0 )|s0 ]] = ρ .
We have
√

(1)

e

(1)

n(Pn1 φ({µ̂k }, {q̂k }) − ρπ ) =pn/n1 Gn1 [φ({µ†k }, {qk })] + Op (1).

Therefore,
p
p
√
†
πe
πe
n(ρ̂DRL(M
−
ρ
)
=
n/n
G
n/n2 Gn2 [φ({µ†k }, {qk })] + Op (1)
[φ({µ
},
{q
})]
+
1
n
k
1
)
k
2
= Gn [φ({µ†k }, {qk })] + Op (1) = Op (1),√eewhich shows ρ̂πDRL(M2 ) is n-consistent around ρπ when the model for the q-function iswell-specified.
Proof of Theorem 17. Almost the same as the proof of Theorem 11
Proof of Theorem 18. We first prove
Pn [ŵt (st )ηt rt ] = Pn [wt (st )ηt rt ] + Pn [(λt−1 − wt (st ))Eπe [rt |st ]] + op (n−1/2 ).
Notingni
1 X h (i)I st = st λt−1 − pπte (st )
n
i=1

1
n
−1

â/b̂ = b−1

{1 − b̂nXh= op (n−1/4 ),
∞

i

(i)

= op (n−1/4 ),I st = st − pπb (st )
t

i=1

∞

(b̂ − b)}{(â − a) − a/b(b̂ − b)},we haveŵt (st ) − wt (st ) + op (n−1/2 )
59

(38)

Kallus and Uehara1=pπb (st )"(ni
1 X h (i)I st = st λt−1 − pπte (st )n)(− wt (st )

i=1

t

n
i
1 X h (i)I st = st − pπb (st )
t
n

)#
.

i=1

Then,
n

1X
(i) (i) (i)wt (st )ηt rt +
n
i=1


n
n
n
(i) (i) X

h
i
h
i
X
X
1
ηt rt
(j)
(i)
(i)
(i)
(j)
(i)
+ 2
I st = st λt−1 − wt (st )I st = st(i)np b (s ) Pn [ŵt (st )ηt rt ] =

i=1

πt

t

j=1

j=1

X
2
0.5(aij + aji ) + op (n−1/2 ),=n(n − 1)
i<j

where
aij =

h
i
(i) (i)(j)(i)
ηt rt I st = st (i)pπb (st )


(j)
(i)
λt−1 − wt (st ) .

tFrom U -statistics theory, by defining bij (H(i) , H(j) ) = 0.5(aij + aji ), we have
n

X
2
2X
bij (H(i) , H(j) ) =E[bij (H(i) , H(j) )|H(i) ] + op (n−1/2 ).
n(n − 1)
n
i<j

i=1

In addition,(i) (i)(i)(i)

E[ai,j |H(i) ] = ηt rt {wt (st ) − wt (st )} = 0,i
 (j) (j) h (j)(i)
ηt rt I st = st (i)(j)λt−1 − wt (st ) |H(i) E[aj,i |H(i) ] = E (j)pπb (st )
t

(i)
(i)
(i) (i) (i)= (λt−1 − wt )E[ηt rt |st ]Therefore, we have shown Eq. (38). Summing over t yields
" T
#
" T
#
X
X
Pn
ŵt (st )ηt rt = Pn
{wt (st )ηt rt + λt−1 Eπe [rt |st ] − wt Eπe [rt |st ]} + op (n−1/2 ),
t=0

t=0

which concludes the proof by establishing the influence function for ρ̂MIS .
Proof of Theorem 19. The difference of the influence functions belongs to the orthogonaltangent space. Therefore, the difference of variances is the variance of the difference of theinfluence functions. This is equal to"#TXvar v0 +−µt qt + µt vt+1 − {λt−1 − wt (st )}Eπe [rt |st ]t=060

Double Reinforcement Learning in Markov Decision Processes= var[v0 ] +

T
+1
X

"
E var E

" T
X

t=1

= var[v0 ] +−µk qk + µk vk+1 − {λk−1 − wk (sk )}Eπe [rk |sk ]|Jst |Jst−1

k=0

T
+1
X

"

"

E var E µt−1 vt +t=1

= var[v0 ] +

T
+1
X

"

"

E var E µt−1 vt −t=1

= var[v0 ] += var[v0 ] += var[v0 ] +

T
X
k=t
T
X

#

!#

−µk rk − {λk−1 − wk (sk )}Eπe [rk |sk ]|Jst |Jst−1!##
λk−1 Eπe [rk |sk ]|Jst |Jst−1

k=t

T
+1
X



E var E [µt−1 vt − λt−1 vt (st )|Jst ] |Jst−1

t=1
T
+1
X
t=1
T
X

!#

#



E var {µt−1 − λt−1 }vt (st )|Jst−1E {wt−1 − λt−1 }2 var (ηt−1 vt (st )|st−1 ) .

t=1

Proof of Theorem 20. We have
√

T
X
e
n{P[
ŵc ηc (sc , ac )rc ] − ρπ }c=0TX

= Gn [ŵc ηc rc −

T
X

c=0

+

√

T
X
wc ηc rc ] + Gn [wc ηc rc ]

c=0

(
n E{

T
X

c=0

)
ŵc ηc rc |{µ̂c }} − ρπec=0

= Gn [TX

wc ηc rc ] +

c=0

√

T
X
e
ŵc ηc rc |{µ̂c }} − ρπ
n E{

(

)
+ op (1)

c=0

T
T
X
X
= op (1) + Gn {wc ηc rc } + Gn [ {λc−1 − wc (sc )}Eπe [rc |sc ]].
c=0

c=0

From the second line to the third line, we used a fact that the stochastic equicontinuity
term is op (1) as in the proof of Theorem 9 because {ŵc ηc rc } belongs to a Donsker class andthe convergence rate condition holds. This fact is confirmed by the fact that a Hölder classwith α > dHsc /2 is a Donsker class (Example 19.9 in van der Vaart, 1998).
From the third line to the fourth line, we have used a result of Example 1(a) in Section 8of Shen (1997). More specifically, the functional derivative of the loss function with respect
to ŵc isg(sc ) → {wc (sc ) − λc−1 }g(sc ),and the induced metric from the loss function is L2 -metric with respect to pπb (sc ). The
functional derivative of the target function with respect to µ̂c isg(sc ) → E[g(sc )ηc rc ] = E[g(sc )Eπe [rc |sc ]],61

Kallus and Uehara

and Riesz representation of the Hilbert space with the induced L2 -metric with respect topπb (sc ) is Eπe [rc |sc ]. Therefore, from Theorem 1 in Shen (1997),
T
T
X
X
e
E[
ŵc ηc rc |{µ̂c }] − ρπ = (Pn − P)
(λc−1 − wc (sc ))Eπe [rc |sc ] + op (n−1/2 ).
c=0

c=0

Proof of Theorem 21. The proof is done as in the proof of Theorem 20. We study thefollowing drift term:
" T
#
X
E
µ̂c rc | µ̂c .
c=0Here, the functional derivative of the loss function with respect to µ̂c (sc , ac ) isg(sc , ac ) → E[g(sc , ac )(µc − λc )].
and Riesz representation of the Hilbert space with the induced L2 -metric with respect topπb (sc , ac ) is E[rc |sc , ac ]. On the other hand, the functional derivative of the target function
with respect to µ̂c isg(sc , ac ) → E[g(sc , ac )rc ] = E[g(sc , ac )E[rc |sc , ac ]],From Theorem 1 in Shen (1997),
" T
#
" T
#
X
X
E
µ̂c rc | µ̂c = Pn
(λc − µc )E[rc |sc , ac ] + op (n−1/2 ).
c=0

c=0

Proof of Theorem 22. We use the general framework developed in Chamberlain (1992) for
establishing the efficiency bounds. For the current problem, noting that the orthogonal
moment condition
E[eq,k+i eq,k ] = E[E[eq,k+i |Hak+i ]eq,k ]] = 0, (0 ≤ k < k + i ≤ T )holds, the efficiency bound for β is represented as
( T
)−1
X
>
∗
∇β mk (Hak ; β ∗ )Σ−1,
k (Hak )∇β mk (Hak ; β )k=0where Σk (Hak ) = var[eq,k |Hak ]. The statement of the theorem for the efficiency boundefor β is arrived at by algebraic simplification of the above. The efficiency bound of ρπ iscalculated similarly.
Proof of Theorem 23. Note var[eq,k |Hak ] and var[eq,k−1 |Hak−1 ] are upper and lower boundedby some constants by assumption. From Jensen’s inequality, we also know
E[g 2 ] = E[E[g 2 |Hsk ]] ≥ E[E[g|Hsk ]2 ].
This conludes that there exists some constant C1 , C2 such thatC1 kgk2 ≤ kgkF,k ≤ C2 kgk2 .
62

Double Reinforcement Learning in Markov Decision Processes

Appendix C. Additional Details from Section 5.2
Cliff Walking. This RL task is detailed in Example 6.6 in Sutton (2018). We consider aboard of size 4 × 12. The horizon was set to T = 400. Each time step incurs −1 reward
until the goal is reached, at which point it is 0, and stepping off the cliff incurs −100 reward
and a reset to the start.
Mountain Car. The RL task is as follows: a car is between two hills in the interval[−0.7, 0.5] and the agent must move back and forth to gain enough power to reach the topof the right hill. The state space comprises position and velocity. There are three discrete
actions: (1) forward, (2) backward, and (3) stay-still. The horizon was set to T = 200. The
reward for each step is −1 until the position 0.5 is reached, at which point it is 0. The state
space was continuous; thus, we obtained a 400-dimensional feature expansion using a radialbasis function kernel as mentioned.
The Policy πd . We construct the policy πd using standard q-learning (Sutton, 2018). For
Cliff Walking, we use a q-learning in a tabular manner. Regarding a Mountain Car, we use
q-learning based on the same feature expansion as above. We use 4000 sample to learn an
optimal policy.

63