EconBase
← Back to paper

Confidence Interval for Off-Policy Evaluation from Dependent Samples via Bandit Algorithm: Approach from Standardized Martingales

Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.

45,702 characters · 16 sections · 50 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Confidence Interval for Off-Policy Evaluation from Dependent Samples via Bandit Algorithm: Approach from Standardized Martingales

abstractThis study addresses the problem of off-policy evaluation (OPE) from dependent samples obtained via the bandit algorithm. The goal of OPE is to evaluate a new policy using historical data obtained from behavior policies generated by the bandit algorithm. Because the bandit algorithm updates the policy based on past observations, the samples are not independent and identically distributed (i.i.d.). However, several existing methods for OPE do not take this issue into account and are based on the assumption that samples are i.i.d. In this study, we address this problem by constructing an estimator from a standardized martingale difference sequence. To standardize the sequence, we consider using \emph{evaluation data} or \emph{sample splitting} with \emph{two-step estimation}. This technique produces an estimator with asymptotic normality without restricting a class of behavior policies. In an experiment, the proposed estimator performs better than existing methods, which assume that the behavior policy converges to a time-invariant policy.

Introduction

The multi-armed bandit (MAB) problem is one of the sequential decision-making problems, which is applied in various applications, such as ad-optimization, personalized medicine, search engines, and recommendation systems. In the MAB problem, we choose an action at a period and observe the outcome. To reduce the cost of searching for better policies, the evaluation of a new policy using historical data obtained from past trials has gathered great attention kdd2009_ads, www2010_cb, AtheySusan2017EPL. This framework is called off-policy evaluation (OPE) dudik2011doubly,wang2017optimal,narita2018,pmlr-v97-bibaut19a,Kallus2019IntrinsicallyES,Oberst2019. Although several methods for OPE have been proposed, existing studies often presume that the samples are independent and identically distributed (i.i.d.). However, in the MAB problem algorithm, the policy is usually updated based on past observations, and the samples are not i.i.d. due to the policy updating process. In this case, the consistency and asymptotic normality of the existing methods are not guaranteed. In particular, the asymptotic normality is critical because it guarantees the $\sqrt{n}$-consistency and is needed for the confidence interval in hypothesis testing to determine whether the new policy is better than the existing policy. Thus, the motivation for establishing a novel method for performing OPE from dependent samples is strong.

Several previous studies consider OPE from dependent samples Laan2008TheCA,1911.02768,KallusNathan2019EBtC,Kato2020. The strategies for deriving asymptotic normality can be classified into three approaches. In the first approach, under the assumption that the policy used in past trials converges to the time-invariant policy in probability, the asymptotic normality is derived using the theories related to martingales Laan2008TheCA,1911.02768,Kato2020. In the second approach, the policy is assumed to be batch updated, where, although the policy is updated using past observations, the sample size under a fixed policy is sufficient Hahn2011. In the third approach, both the stationarity and the conditions of mixingales KosorokMichaelR2008ItEP, which requires the independence of time-separated samples, are assumed KallusNathan2019EBtC. However, existing methods have several drawbacks. With the first approach, the assumption that the policy converges to a time-invariant policy sometimes does not hold. For example, if the algorithms of the MAB problem change during past trials, it is difficult to justify this assumption. The second and third approaches cannot be applied in cases where the policy is sequentially updated, and the conditions of mixingales do not hold, respectively.

To overcome these drawbacks, in this study, we propose a novel estimator with the asymptotic normality, which does not assume policy convergence, the existence of batches, and the independence of time-separated samples. To accomplish this goal, we refocus on the first approach. In the first approach, we use a martingale difference sequence (MDS) for constructing an estimator. In this approach, the convergence of the policy used in past trials to the time-invariant policy in probability is assumed for using central limit theorem (CLT) of an MDS, which requires that the variance of an MDS becomes asymptotically constant. In this paper, instead of assuming the convergence of the policy, we consider constructing an MDS with asymptotically constant variance by standardization using an estimator of the variance. To implement the standardization, we assume access to evaluation (test) data or conduct sample splitting. Through this technique, we can guarantee the asymptotic normality of the proposed estimator with fewer assumptions.

\paragraph{Contributions and Organization of this Paper:} This paper has three main contributions. First, this paper provides a theoretical solution for causal inference from dependent samples obtained via bandit feedback, which is more realistic than the situation that has been often considered in the existing work. Second, the proposed estimator achieves the asymptotic normality with fewer assumptions compared with existing solutions. Third, the estimator also achieves a lower mean squared error (MSE) in experiments using benchmark datasets.

\paragraph{Additional Related Work:} There are various studies for OPE under the assumption that samples are i.i.d. dudik2011doubly,wang2017optimal,narita2018,pmlr-v97-bibaut19a,Kallus2019IntrinsicallyES,Oberst2019, but there are fewer studies focus on the case in which samples are not i.i.d. When the policy converges, Laan2008TheCA, 1911.02768, and Kato2020 proposed using an estimator based on an MDS. When including the covariates for OPE, Laan2008TheCA mainly suggested using targeted maximum likelihood estimation. On the other hand, 1911.02768 and Kato2020 proposed using an estimator based on the doubly robust or augmented IPW estimator. KallusNathan2019EBtC dealt with the application of reinforcement learning and proposed an estimator constructed from dependent samples. In their study, they proposed time cross-fitting, which is a variant of cross-fitting of double/debiased machine learning ChernozhukovVictor2018Dmlf. The theoretical guarantee is based on mixingale theory. For obtaining the asymptotic normality, conducting the standardization is also proposed by Luedtke2016. Luedtke2016 tries to solve a similar problem in the batch policy updating. In addition to this difference of policy updating, there are the following two points different from us. First, in their solution, they require estimators of $f^*$ and $e^*$ to satisfy some conditions like the Donsker condition. However, in our estimator, we do not impose such restrictions on estimators of $f^*$ and $e^*$ except for the boundedness. This difference is critical because, in practice, as ChernozhukovVictor2018Dmlf pointed out, the estimators of $f^*$ and $e^*$ usually do not satisfy the Donsker condition when we use modern methods such as the Lasso and Random forest regression. Second, for obtaining an estimator with the asymptotic normality under sequential policy updating, we propose using evaluation data or sample splitting.

Problem Setting

In this section, we introduce our problem setting.

Date Generating Process

Let $A_t$ be the action taking variable in $\mathcal{A}=\{1,2,\dots,K\}$, $X_t$ be the covariate observed by the decision maker when choosing an action, and $\mathcal{X}$ be the domain of covariate. Let us denote a random variable of a reward at a period $t$ as a function $Y_t:\mathcal{A}\to\mathbb{R}$. Let $\mathbbm{1}[\cdot]$ be an indicator function. In this paper, we have access to a set of historical data, $\mathcal{S}_T=\{(X_t, A_t, Y_t)\}^{T}_{t=1}$, with the following data generating process (DGP):

align[align omitted — 114 chars of source]

where $Y_t = \sum^K_{a=1}\mathbbm{1}[A_t=a]Y_t(a)$, $p(x)$ denote the density of the covariate $X_t$, $p_t(a\mid x, \Omega_{t-1})$ denote the probability of assigning an action $A_t$ conditioned on a covariate $X_t$, $p(y\mid a, x)$ denote the density of an outcome $Y_t$ conditioned on $A_t$ and $X_t$, and $\Omega_{t-1}\in \mathcal{M}_{t-1}$ denotes the history defined as $\Omega_{t-1}=\{X_{t-1}, A_{t-1}, Y_{t-1}, \dots, X_{1}, A_1, Y_{1}\}$ with the domain $\mathcal{M}_{t-1}$. We assume that $p(x)$ and $p(y\mid a, x)$ are invariant across periods, but $p_t(a\mid x, \Omega_{t-1})$ can take different value across periods. Besides, we allow the decision maker to change $p_t(a\mid x, \Omega_{t-1})$ based on past observations. In this case, samples $\big\{(X_t, A_t, Y_t)\big\}^{T}_{t=1}$ are correlated across periods, i.e., samples are not i.i.d. Here, we introduce behavior policies, which determine the probability $p_t(a\mid x, \Omega_{t-1})$. Let a behavior policy $\pi_t:\mathcal{A}\times\mathcal{X}\times\mathcal{M}_{t-1}\to(0,1)$ in the $t$-th period be a function of a covariate $X_t$, an action $A_t$, and a history $\Omega_{t-1}$. In this study, we assume that $\pi_t(a\mid x, \Omega_{t-1}) = p_t(a\mid x, \Omega_{t-1})$.

remark[Evaluation Data] In Section (ref), we also assume that there is access to evaluation data, $\mathcal{E}_N = \big\{X_i\big\}^{N}_{i=1}$. In Remark (ref), we relax this assumption by introducing the sample splitting.

Off-Policy Evaluation

Under the DGP defined in the previous subsection, we consider estimating the value of an evaluation policy using samples obtained under the behavior policies. Let an evaluation policy $\pi^\mathrm{e}:\mathcal{A}\times\mathcal{X}\to(0,1)$ be a function of a covariate $X_t$ and an action $A_t$; this evaluation policy can be considered as the probability of taking the action $A_t$ conditional on the covariate $X_t$. We are interested in estimating the expected reward from a given pre-specified evaluation policy $\pi^{\mathrm{e}}(a \mid x)$. Then, we define the expected reward under $\pi^\mathrm{e}$ as $R(\pi^\mathrm{e}) := \mathbb{E}\left[\sum^{K}_{a=1}\pi^\mathrm{e}(a, X_t)Y_t(a)\right]$. We also denote $R(\pi^\mathrm{e})$ as $\theta_0$. The goal of this study is to estimate $R(\pi^{\mathrm{e}})$ using dependent samples under a sequentially-updated policy. To identify $R(\pi^{\mathrm{e}})$, we make the following assumptions on the policies and the outcomes.

assumptionThe behavior policies $\pi_t(a\mid \cdot, \Omega_{t-1})$ and the evaluation policy $\pi^\mathrm{e}(a\mid \cdot)$ are known\footnote{The proposed method requires the function of the polices, not only $\pi_t(a\mid X_t, \Omega_{t-1})$ for a specific $X_t$.}. Additionally, the evaluation policy is deterministic.
assumptionThere exist $C_1$ and $C_2$ such that $\frac{\pi^\mathrm{e}(a\mid x)}{\pi_t(a\mid x, \Omega_{t-1})}\leq C_1$ and $|Y_t| \leq C_2$.

The deterministic evaluation policy looks restrictive. However, the optimal policy that maximizes the expected reward is a deterministic policy when the policy is not restricted. Therefore, when we find the optimal policy from a set of policies, such as off-policy learning ZhaoYingqi2012EITR,KitagawaToru2018WSBT,ZhouZhengyuan2018OMPL,Chernozhukov2019, this assumption can be accommodated. Thus, restricting the policy class to the deterministic policy class would be reasonable in practice. In fact, most of the existing methods follow this way.

remark[Existing Methods for OPE] We review three types of standard estimators of $R(\pi^{\mathrm{e}})$ under the case where $p_1(a\mid x)=p_2(a\mid x)=\cdots=p_T(a\mid x)=p(a\mid x)$ in the DGP defined in ((ref)). The first estimator is an inverse probability weighting (IPW) estimator given by $\frac{1}{T}\sum^T_{t=1}\sum^K_{a=1}\frac{\pi^\mathrm{e}(a\mid X_t)\mathbbm{1}[A_t=a]Y_t}{p(a\mid X_t)}$ Horvitz1952,rubin87,hirano2003efficient. Although this estimator is unbiased when the behavior policy is known, it has a high variance. The second estimator is a direct method (DM) estimator $\frac{1}{T}\sum^T_{t=1}\sum^K_{a=1}\hat{f}_{T}(a, X_t)$, where $\hat{f}_{T}(a, X_t)$ is an estimator of $f^*(a, X_t)$ HahnJinyong1998OtRo. This estimator is known to be weak against model misspecification for $f^*(a, X_t)$. The third estimator is a doubly robust estimator defined as $\frac{1}{T}\sum^T_{t=1}\sum^K_{a=1}\left(\frac{\pi^\mathrm{e}(a\mid X_t)\mathbbm{1}[A_t=a]\big(Y_t - \hat{f}_{T}(a, X_t)\big)}{p(a\mid X_t)} + \pi^\mathrm{e}(a\mid X_t)\hat{f}_{T}(a, X_t)\right)$ robins94,ChernozhukovVictor2018Dmlf. Under certain conditions, it is known that this estimator achieves the lower bound (a.k.a semiparametric lower bound), which is the lower bound of the asymptotic MSE of OPE, among regular $\sqrt{T}$-consistent estimators VaartA.W.vander1998As.
remark[Unconfoundedness] Existing methods often make an assumption called unconfoundedness, that is, they assume that the outcomes $(Y_t(0), Y_t(1), \dots, Y_t(K))$ and the action $A_t$ are conditionally independent on $X_t$ for identification. In the DGP, we can obtain a similar result from the assumption that an action $A_t$ is chosen according to probability $p_t(a\mid x, \Omega_{t-1})$.
remark[Semiparametric Lower Bound] The lower bound of the variance can be defined for an estimator under some posited models of the DGP. If this posited model is a parametric model, then the lower bound is equal to the Cram\'{e}r-Rao lower bound. When this posited model is a non- or semiparametric model, the corresponding lower bound can still be defined bickel98. narita2018 shows that the semiparametric lower bound of the DGP ((ref)) under $p_1(a\mid x, \Omega_{0})=p_2(a\mid x, \Omega_{1})=\cdots=p_T(a\mid x, \Omega_{T-1})=p(a\mid x)$ is $\mathbb{E}\left[\sum^{K}_{a=1}\left\{\frac{\big(\pi^\mathrm{e}(a\mid X_t)\big)^2v^*(a, X_t)}{p(a\mid X_t)} + \Big(\pi^\mathrm{e}(a\mid X_t)f^*(a, X_t) - \theta_0\Big)^2\right\}\right]$.

\paragraph{Notations:} Let $a$ be an action in $\mathcal{A}$. We denote $\mathbb{E}[Y_t(a)\mid x]$, $\mathbb{E}[Y^2_t(a)\mid x]$, and $\mathrm{Var}(Y_t(a)\mid x)$ as $f^*(a, x)$, $e^*(a, x)$, and $v^*(a, x)$, respectively. Let $\mathcal{F}$ be the domain of $f^*(a, x)$. Let $\hat{f}_{t}(a, x)$, $\hat e_{t}(a, x)$, and $\hat{\theta}_t$ be estimators of $f^*(a, x)$, $e^*(a, x)$, and $\theta_0$ constructed from $\Omega_{t}$, respectively. Let $\mathcal{N}(\mu, \mathrm{var})$ be the normal distribution with mean $\mu$ and variance $\mathrm{var}$.

Policy Updating

In the DGP defined in Section (ref), we allow the behavior policy $\pi_t(a \mid x, \Omega_{t-1})$ to be updated based on past observations. For example, when applying some algorithms proposed in the MAB problem, we usually optimize the probability to maximize the cumulative reward using past observations. When the behavior policy $\pi_t(a \mid x, \Omega_{t-1})$ depends on past observations, the samples are not i.i.d. In this case, we cannot apply the standard methods for OPE, which assume that the samples are i.i.d. To clarify this problem, we classify the patterns of policy updating into two scenarios: sequential policy updating and batch policy updating. In sequential policy updating, the policy can be updated each period. In batch policy updating, the policy is updated after observing some samples Hahn2011,narita2018. In this study, we consider a problem using sequential policy updating.

OPE under a Converging Policy

In this section, we review OPE when a policy converges to the time-invariant policy in probability. As Laan2008TheCA, 1911.02768 and Kato2020 proposed, in such a situation, we can construct an estimator with the asymptotic normality by using an MDS hall2014martingale. To conduct OPE when the behavior policy is known, we can use the following score function $\psi:\Theta\times\mathcal{X}\times\mathcal{A}\times\mathbb{R}\times\Pi_t\times\Pi\times\mathcal{F}\to\mathbb{R}$ ChernozhukovVictor2018Dmlf:

align[align omitted — 262 chars of source]

where $\Theta$, $\Pi_t$, and $\Pi$ are spaces of a parameter $\theta$, a behavior policies $\pi_t$, and an evaluation policy $\pi^\mathrm{e}$. Then, for the true policy value $\theta_0$, the score function with $\theta = \theta_0$, $\psi_t\big(\theta_0; X_t, A_t, Y_t, \pi_t, \pi^\mathrm{e}, \hat{f}_{t - 1}\big)$, satisfies $\mathbb{E}[\psi_t\big(\theta_0; X_t, A_t, Y_t, \pi_t, \pi^\mathrm{e}, \hat{f}_{t - 1}\big)] = 0$. As we explain in Remark (ref), there are several other candidate estimators, including the IPW estimator. However, the estimator based on the score function defined in ((ref)) achieves the lowest asymptotic variance when $\pi_t$ is fixed among periods bickel98,scharfstein99. Using $\big\{\psi_t\big(\theta_0; X_t, A_t, Y_t, \pi_t, \pi^\mathrm{e}, \hat{f}_{t - 1}\big)\big\}^T_{t=1}$, we consider constructing an MDS $\big\{z_t\big\}^T_{t=1}$ that satisfies the conditions of the following CLT for an MDS.

proposition[CLT for an MDS, GVK126800421, Proposition 7.9, p. 194] Let $\{z_t\}^\infty_{t=1}$ be an MDS. Suppose that (a) $\mathbb{E}\big[z^2_t\big] = \nu^2_t > 0$ with $\big(1/T\big) \sum^T_{t=1}\nu^2_t\to \nu^2 > 0$; (b) $\mathbb{E}\big[|z_t|^r\big] < \infty$ for some $r>2$; (c) $\big(1/T\big)\sum^T_{t=1}z^2_t\xrightarrow{\mathrm{p}} \nu^2$. Then $\sqrt{T}\frac{1}{T}\sum^T_{t=1}z_t\xrightarrow{\mathrm{d}}\mathcal{N}(0, \nu^2)$.

When $Y_t$, $\pi^\mathrm{e}/\pi_t$, and $f_{t - 1}$ are bounded, we can easily show that the condition (b) holds (Assumption (ref)). The remaining task is to show that conditions (a) and (c) hold.

Kato2020 constructed $\ddot{z}_t = \psi_t\big(\theta_0; X_t, A_t, Y_t, \pi_t, \pi^\mathrm{e}, \hat{f}_{t - 1})$, where $\left\{\ddot{z}_t\right\}^T_{t=1}$ is an MDS. Then, by assuming $\pi_t(a\mid x, \Omega_{t-1})\xrightarrow{\mathrm{p}}\alpha(a\mid x)$ and some regularity conditions, Kato2020 proved that the sequence $\left\{\ddot{z}_t\right\}^T_{t=1}$ satisfies the conditions (a) and (c) with $\nu^2 = \mathbb{E}\left[\left(\sum^{K}_{a=1}\left\{\frac{\pi^\mathrm{e}(a\mid x)\mathbbm{1}[A_t=a]\big\{Y_t - f^*(a, X_t)\big\}}{\alpha(a\mid X_t)} + \pi^\mathrm{e}(a\mid X_t)f^*(a, X_t)\right\} - \theta_0\right)^2\right]$. Thus, by the assumption $\pi_t(a\mid x, \Omega_{t-1})\xrightarrow{\mathrm{p}}\alpha(a\mid x)$, the variance of $h_t$ converges to time-invariant value, and this property is the motivation of the assumption proposed by Kato2020. Then, under $\pi_t(a\mid x, \Omega_{t-1})\xrightarrow{\mathrm{p}}\alpha(a\mid x)$, the estimator proposed defined as $\hat{\theta}^{\mathrm{A2IPW}}_T = \frac{1}{T}\sum^T_{t=1}\sum^{K}_{a=1}\left\{\frac{\pi^\mathrm{e}(a\mid X_t)\mathbbm{1}[A_t=a]\big\{Y_t - \hat{f}_{t-1}(a, X_t)\big\}}{\pi_{t}(a\mid X_t, \Omega_{t - 1})} + \pi^\mathrm{e}(a\mid x)\hat{f}_{t-1}(a, X_t)\right\}$ can be proved to be consistent and has the asymptotic normality 1911.02768,Kato2020\footnote{They assume the almost sure convergence of the policy for showing the mean convergence of the variance, but we can relax it to the convergence in probability by the boundedness of the MDS loeve1977probability.}. We refer to this estimator as Adaptive AIPW (A2IPW). We can also define the corresponding IPW estimator and refer to it as adaptive IPW (AdaIPW) estimator.

remark[Efficiency of A2IPW] The asymptotic variance of the A2IPW estimator matches the semiparametric lower bound under a time-invariant policy $\alpha$ Kato2020.

OPE from Dependent Samples under Sequential Policy Updating

In this section, we discuss a strategy for constructing an estimator without assuming the convergence of the policy. Then, we show a method for constructing an estimator with the asymptotic normality.

Strategy for OPE from Dependent Samples under Sequential Policy Updating

However, there are various applications in which we cannot assume that $\pi_t(a\mid x, \Omega_{t-1})\xrightarrow{\mathrm{p}}\alpha(a\mid x)$. Therefore, we have a strong motivation for constructing an estimator with asymptotic normality without assuming that $\pi_t(a\mid x, \Omega_{t-1})\xrightarrow{\mathrm{p}}\alpha(a\mid x)$. Our strategy is to construct an MDS $\{\tilde{z}_t\}^T_{t=1}$ satisfying conditions (a) and (c) by standardization of $\ddot{z}_t$ as $\tilde{z}_t = \left(\sqrt{\mathrm{Var}\big(\ddot{z}_t\mid \Omega_{t-1}\big)} \right)^{-1}\ddot{z}_t$, where $\mathrm{Var}\left(\ddot{z}_t\mid \Omega_{t-1}\right)$ is the conditional variance of $\ddot{z}_t$. The standardization means that the variance of $\tilde{z}_t$ becomes $1$. Then, we can easily show that the conditions (a) and (b) hold. Thus, if we know the variance $\mathrm{Var}\big(\ddot{z}_t\mid \Omega_{t-1}\big)$, we can construct an MDS with constant variance by the standardization.

remark[Lindeberg Condition] The proposed strategy can be interpreted using the Lindeberg condition. For all $\epsilon > 0$ and an MDS $\{z_t\}^T_{t=1}$ such that $\sum^T_{t=1}\mathbb{E}\big[z^2_t \mid \Omega_{t-1}\big]/\sum^T_{t=1}\mathbb{E}\big[z^2_t\big] \xrightarrow{\mathrm{p}} 1$, the Lindeberg condition is $\left(\sum^T_{t=1}\mathbb{E}\big[z^2_t\big]\right)^{-1}\sum^T_{t=1}\mathbb{E}\left[z^2_t\mathbbm{1}\left[|z_t|\geq \epsilon \sqrt{\sum^T_{t=1}\mathbb{E}\big[z^2_t\big]}\right]\right]\xrightarrow{\mathrm{p}} 0$ as $T\to\infty$ brown1971. The first condition $\sum^T_{t=1}\mathbb{E}\big[z^2_t \mid \Omega_{t-1}\big]/\sum^T_{t=1}\mathbb{E}\big[z^2_t\big] \xrightarrow{\mathrm{p}} 1$ restricts the class of an MDS to an MDS with asymptotically time-invariant variance. Here, we note that $\sigma^2_t = \mathbb{E}\big[z^2_t \mid \Omega_{t-1}\big]$. The strategy of this paper can be intuitively interpreted as a transformation of the original sequence $\ddot{z}_t$ into such an MDS with asymptotically time-invariant variance.

OPE with a Weighted Average Estimator: Case with Known Variance

Then, we propose an estimator based on the MDS $\big\{\tilde{z}_t\big\}^T_{t=1}$ and refer to the estimator as a weighted Average Adaptive Augmented IPW (A3IPW) estimator. For $\big\{\hat{f}_t\big\}^{T-1}_{t=0}$ and $\sigma^2_t = \mathrm{Var}\big(\ddot{z}_t\mid \Omega_{t-1}\big)$, we define the A3IPW estimator as $\hat{\theta}^{\mathrm{A3IPW}}_T = \left(\sum^T_{t=1}\frac{1}{\sqrt{\sigma^2_t}}\right)^{-1}\sum^T_{t=1}\frac{1}{\sqrt{\sigma^2_t}}\sum^{K}_{a=1}\left\{\frac{\pi^\mathrm{e}(a\mid X_t)\mathbbm{1}[A_t=a]\big\{Y_t - \hat{f}_{t-1}(a, X_t)\big\}}{\pi_{t}(a\mid X_t, \Omega_{t - 1})} + \pi^\mathrm{e}(a\mid x)\hat{f}_{t-1}(a, X_t)\right\}$. Then, we can show that the estimator has the asymptotic normality as follows.

theoremSuppose that there exists $C_2$ such that $|\hat{f}_t|\leq C_2$. Then, under Assumption (ref), $\left(\sum^T_{t=1}\frac{1}{\sqrt{\sigma^2_t}}\right)\sqrt{T}\left(\hat{\theta}^{\mathrm{A3IPW}}_{T} - \theta_0\right) \xrightarrow{\mathrm{d}} \mathcal{N}\left(0, 1\right)$ as $T\to\infty$.
proofFor $\tilde{z}_t =\frac{1}{\sqrt{\sigma^2_t}}\left(\sum^{K}_{a=1}\left\{\frac{\pi^\mathrm{e}(a\mid X_t)\mathbbm{1}[A_t=a]\big\{y - \hat{f}_{t-1}(a, X_t)\big\}}{\pi_{t}(a\mid X_t, \Omega_{t - 1})} + \pi^\mathrm{e}(a\mid x)\hat{f}_{t-1}(a, X_t)\right\} - \theta_0\right)$, $\big\{\tilde{z}_t\big\}^{T}_{t=1}$ is an MDS because $\mathbb{E}\left[\tilde{z}_t\mid \Omega_{t-1}\right] = \frac{1}{\sqrt{\sigma^2_t}}\mathbb{E}\left[\ddot{z}_t\mid \Omega_{t-1}\right] = 0$. Because $\mathrm{Var}\left(\ddot{z}_t/\sqrt{\sigma^2_t}\right) = \mathbb{E}\big[\mathbb{E}[(\ddot{z}_t - \mathbb{E}[\ddot{z}])^2 \mid \Omega_{t-1}]/\sigma^2_t\big] = \mathbb{E}\big[\mathbb{E}[\ddot{z}^2_t \mid \Omega_{t-1}]/\sigma^2_t\big] = \mathbb{E}\big[\mathrm{Var}\big(\ddot{z}_t \mid \Omega_{t-1}\big)/\sigma^2_t\big] =1$, we have $\big(1/T\big) \sum^T_{t=1}\mathrm{Var}\left(\tilde{z}_t\right)= 1$. Here, we used $\mathbb{E}[\ddot{z}_t \mid \Omega_{t-1}]=\mathbb{E}[\ddot{z}_t]=0$. Besides, we can show that $\tilde{z}_t$ is bounded from assumptions, and $\big(1/T\big)\sum^T_{t=1}\tilde{z}^2_t \xrightarrow{p} 1$ from the weak law of large numbers of Proposition (ref) in Appendix (ref). Then, we can apply CLT for an MDS.

OPE with a Weighted Average Estimator: Case with Unknown Variance

In general, the variance $\sigma^2_t$ is unknown and needed to be replaced with an estimator. However, we cannot use the samples used for constructing $\pi_t$ and $\hat{f}_{t-1}$ to estimate $\sigma^2_t$ due to the dependency. For this reason, we assume access to a dataset only with covariates, $\mathcal{E}_N = \big\{X_i\big\}^{N}_{i=1}$, which is independent from the historical data $\mathcal{S}_T$. Because $\mathcal{E}_N$ is independent from $\mathcal{S}_T = \big\{(X_t, A_t, Y_t)\big\}^{T}_{t=1}$, we can approximate $\sigma^2_t$ by applying the law of large numbers to the sample average $g'_{t, N}$, defined as

align*[align* omitted — 291 chars of source]

The consistency of $g'_{t, N}$ to the variance of $\psi_t\big(\theta_0; X_t, A_t, Y_t, \pi_t, \pi^\mathrm{e}, f^*\big)$ is showed as follows.

theoremSuppose that, for all $x\in\mathcal{X}$ and $a\in\mathcal{A}$, $\hat{f}_{t-1}(a, x)-f^*(a, x)\xrightarrow{\mathrm{p}}0$, $\hat{e}_{t-1}(a, x)-e^*(a, x)\xrightarrow{\mathrm{p}}0$ and $\hat{\theta}_t\xrightarrow{\mathrm{p}}\theta_0$ as $t\to\infty$. Then, for $\pi_t\in\Pi$, we have $g'_{t, (1-r)T} - \sigma^{*2}_t\xrightarrow{\mathrm{p}}0$ as $t\to\infty$ and $N\to\infty$, where $\sigma^{*2}_t = \mathbb{E}\left[\sum^{K}_{a=1}\left\{\frac{\big(\pi^\mathrm{e}(a\mid X_t)\big)^2\nu^*(a, X_t)}{\pi_{t}(a\mid X_t, \Omega_{t - 1})} + \Big(\pi^\mathrm{e}(a\mid X_t)f^*(a, X_t)-\theta_0\Big)^2\right\}\right]$.

In practice, we need to avoid having the estimator of $\sigma^{*2}_t$ be less than or equal to $0$. Thus, instead of $g'_{t, (1-r)T}$, we use the estimator $g_{t, N}$, defined as $g_{t, (1-r)T}=\max\big\{g'_{t, N}, \epsilon\big\}$, where $\epsilon > 0$ is the lower bound of the variance $\sigma^{*2}_t$. The estimator $g_{t, N}$ is consistent, as is $g'_{t, N}$. The value $\epsilon$ is assumed to be known, and we can use any sufficiently small value for $\epsilon$ because it is simply a technical term. The proof is shown in Appendix (ref). Then, using $g_{t, N}$, we define the Feasible A3IPW (FA3IPW) estimator $\hat{\theta}^{\mathrm{FA3IPW}}_{T}$ as $\left(\sum^{T}_{t=1}\frac{1}{\sqrt{g_{t, N}}}\right)^{-1}\sum^{T}_{t=1}\sum^{K}_{a=1}\frac{1}{\sqrt{g_{t, N}}}\left\{\frac{\pi^\mathrm{e}(a\mid X_t)\mathbbm{1}[A_t=a]\big\{Y_t - \hat{f}_{t-1}(a, X_t)\big\}}{\pi_{t}(a\mid X_t, \Omega_{t - 1})} + \pi^\mathrm{e}(a\mid x)\hat{f}_{t-1}(a, X_t)\right\}$. Similarly, we can also define the corresponding IPW estimator and refer to it as Feasible weighted Average Adaptive IPW estimator (FA2daIPW).

In some applications, we can obtain the dataset $\mathcal{E}_N = \big\{X_i\big\}^{N}_{i=1}$ as evaluation (test) data, which is the target of a new policy kato_uehara_2020. On the other hand, by introducing a sample splitting, we can relax the assumption as the following remark.

remark[Sample Splitting] To relax the access to $\mathcal{E}_N$, we consider splitting $\big\{(X_t, A_t, Y_t)\big\}^{T}_{t=1}$ into $\big\{(X_t, A_t, Y_t)\big\}^{\lfloor rT\rfloor}_{t=1}$ and $\big\{(X_t, A_t, Y_t)\big\}^{T}_{t=\lfloor rT\rfloor + 1}$, where $0 < r < 1$. In the following discussion, for simplicity, we assume that $rT$ is an integer (i.e., $\lfloor rT\rfloor=rT$). The samples $\big\{X_t\big\}^{T}_{t=rT}$ and $\big\{(X_t, A_t, Y_t)\big\}^{rT}_{t=1}$ are independent. Thus, to approximate $\sigma^2_t$, we can apply the law of large numbers to the sample average $g_{t, (1-r)T} = \max\left\{g'_{t, (1-r)T}, \epsilon\right\}$, where $g'_{t, (1-r)T}$ is defined as $\frac{1}{(1-r)T} \sum^T_{s = rT} \sum^{K}_{a=1}\left\{\frac{\big(\pi^\mathrm{e}(a\mid X_s)\big)^2\big(\hat{e}_{t-1}(a, X_s) - \hat{f}^2_{t-1}(a, X_s)\big)}{\pi_{t}(a\mid X_s, \Omega_{t - 1})} + \Big(\pi^\mathrm{e}(a\mid X_s)\hat{f}_{t-1}(a, X_s) - \hat{\theta}_{t-1} \Big)^2\right\}$. More details on sample splitting are provided in Appendix (ref).

Asymptotic Properties of a FA3IPW Estimator

Let $\hat{z}_t = \sum^{T}_{t=1}\frac{1}{\sqrt{g_{t, N}}}\left(\sum^{K}_{a=1}\left\{\frac{\pi^\mathrm{e}(a\mid X_t)\mathbbm{1}[A_t=a]\big\{Y_t - \hat{f}_{t-1}(a, X_t)\big\}}{\pi_{t}(a\mid X_t, \Omega_{t - 1})} + \pi^\mathrm{e}(a\mid X_t)\hat{f}_{t-1}(a, X_t)\right\} - \theta_0\right)$. We can show that $\hat{z}_t$ is an MDS as $\mathbb{E}\left[\hat{z}_t\mid \Omega_{t-1}\right] = \frac{1}{\sqrt{g_{t, N}}}\mathbb{E}\left[\ddot{z}_t\mid \Omega_{t-1}\right] = 0$. Based on this property, the following theorems gives us the consistency and asymptotic normality of $\hat{\theta}^{\mathrm{FA3IPW}}_{T}$. First, we show the consistency of $\hat{\theta}^{\mathrm{FA3IPW}}_{T}$, i.e., $\hat{\theta}^{\mathrm{FA3IPW}}_{T} \xrightarrow{\mathrm{p}} \theta_0$.

theorem[Consistency of the Proposed Estimator] Suppose that there exist $C_2$ and $C_3$ such that $|\hat{f}_t|\leq C_2$ and $0 < 1/\sqrt{g_{t, N}} \leq C_3$. Then, under Assumptions (ref), $\hat{\theta}^{\mathrm{FA3IPW}}_{T} \xrightarrow{\mathrm{p}} \theta_0$.
proofFrom the boundedness, we can apply the law of large numbers for an MDS (Proposition (ref) in Appendix (ref)). Therefore, $\frac{1}{T}\sum^{T}_{t=1}\hat{z}_t \xrightarrow{\mathrm{p}} 0$.

The following theorem shows the asymptotic normality of $\hat{\theta}^{\mathrm{FA3IPW}}_T$. The proof is in Appendix (ref).

theorem[Asymptotic Distribution of the Proposed Estimator] Suppose that \begin{description} • there exist $C_2$ and $C_3$ such that $|\hat{f}_t|\leq C_2$ and $0 < 1/\sqrt{g_{t, N}} \leq C_3$; • Pointwise convergence in probability of $\hat{f}_t$, i.e., $\forall x \in \mathcal{X}$, $\hat{f}_t(a, x)- f^*(a, x)\xrightarrow{\mathrm{p}}0$ as $t\to \infty$; • $g_{t, N}- \sigma^{*2}_t\xrightarrow{\mathrm{p}}0$ as $t\to \infty$ and $N\to \infty$. \end{description} Then, under Assumptions (ref), $\left(\frac{1}{\sqrt{T}}\sum^{T}_{t=1}\frac{1}{\sqrt{g_{t, N}}}\right)\left(\hat{\theta}^{\mathrm{FA3IPW}}_{T} - \theta_0\right) \xrightarrow{\mathrm{d}} \mathcal{N}\left(0, 1\right)$ as $T, N\to \infty$.
remark[Difference of Assumption between the Consistency and Asymptotic Normality] The asymptotic normality requires $g_{t, N}- \sigma^{2*}_t\xrightarrow{\mathrm{p}}0$. However, we can obtain the consistency $\hat{\theta}^{\mathrm{A3IPW}}_T \xrightarrow{\mathrm{p}} \theta_0$ without assuming $g_{t, N}- \sigma^{*2}_t\xrightarrow{\mathrm{p}}0$. From this property, we can obtain the consistent estimator of the variance, $\sigma^{*2}_t$, which requires the consistent estimator of the policy value, $\theta_0$.
remark[Heuristic Stabilization] The proposed estimator requires a sequence of $\big\{f_t\big\}^{T-1}_{t=0}$. Because the early estimator of $f_t$ may be unstable due to the insufficient amount of training data, some heuristic stabilizations may improve its performance. 1911.02768 proposed a similar estimator and suggested using a weight sequence for stabilization. In the setting of Kato2020, they can change $p_t(a\mid x, \Omega_{t-1})$ and proposed adjusting the probability to stabilize their proposed estimator. In this study, we take a simple approach by suggesting discarding some samples at first and constructing the proposed estimator without them. We call this estimator the stabilized FA3IPW (SFA3IPW) estimator. In some experiments, this approach improves the estimator's performance.

Efficiency of A3IPW Estimator

In the problem setting, it is not easy to discuss an efficiency of an estimator. For simplicity, we consider the lower bound of the variance of an estimator within a class of estimators such that $\hat{\theta}_T = \sum^{T}_{t=1}w_t\tilde{\psi}_t(X_t, A_t, Y_t, \pi_t)$, where $\tilde{\psi}$ is a score function and $w_t$ is a weight satisfying $\sum^T_{t=1}w_t = 1$. Let us note that we usually use $w_t=1/T$ in the case where samples are i.i.d. First, we consider the variance of $\tilde{\psi}_t$. In OPE, the lowest variance is given as $\mathbb{E}\left[\sum^{K}_{a=1}\left\{\frac{\big(\pi^\mathrm{e}(a\mid X_t)\big)^2v^*(a, X_t)}{\pi_t(a\mid X, \Omega_{t-1})} + \Big(\pi^\mathrm{e}(a\mid X_t)f^*(a, X_t) - \theta_0\Big)^2\right\}\right]$. In the proposed estimator, each element of $\{\tilde{\psi}_t\}^{T}_{t=1}$ achieves this variance. In addition, it is well known that for a sequence of positive values $\phi_t = \mathbb{E}[\tilde{\psi}^2_t \mid \Omega_{t-1}]$, the weight $w_t =\frac{1}{\sqrt{\phi_t}}/\sum^{T}_{i=1}\frac{1}{\sqrt{\phi_i}}$ minimizes the variance of $\hat{\theta}_T$ when the covariances among the score functions are $0$. When $\{\phi_t\}^T_{t=1}$ is an MDS, the covariance is $0$. Thus, the proposed estimator uses the weights and scores that minimize the variance for an estimator with a form that satisfies $\sum^{T}_{t=1}w_t\tilde{\psi}_t(X_t, A_t, Y_t, \pi_t)$.

Main Algorithm: Two-step FA3IPW Estimation

Under the assumption of $g_{t, N}- \sigma^{*2}_t\xrightarrow{\mathrm{p}}0$, the proposed estimator $\hat{\theta}^{\mathrm{A3IPW}}_T$ has the asymptotic normality. However, to estimate $\sigma^{*2}_t = \mathrm{Var}\big(\ddot{z}_t\mid \Omega_{t-1}\big)$, a consistent estimator of $\theta_0$ is required, which is why we use $g_{t, N}$. However, we can obtain a consistent estimator of $\theta_0$ without imposing the assumption that $g_{t, N} - \sigma^{*2}_t\xrightarrow{\mathrm{p}} 0$. Based on these properties, we propose two-step estimation. First, using arbitrary appropriate values for $g_{t, N} = \tilde{g}_{t}$ for $t=1,2,\dots,T-1$, such as $\tilde{g}_{t}=1$, we obtain the initial estimates, $\left\{\tilde{\theta}_t\right\}^{T-1}_{t=1}$. Then, using $\tilde{\theta}_t \xrightarrow{\mathrm{p}} \theta_0$, we construct $g_{t, N} - \sigma^{*2}_t \xrightarrow{\mathrm{p}} 0$. Finally, we can estimate the next estimator of $\theta_0$, $\hat{\theta}^{\mathrm{TSFA3IPW}}_{T}$, with asymptotic normality from Theorem (ref). We summarize the process of two-step FA3IPW (TSFA3IPW) estimation in Algorithm (ref).

algorithm[algorithm omitted — 1,164 chars of source]
remark[Estimation of Functions \texorpdfstring{$f^*_t$}{f} and \texorpdfstring{$e^*_t$}{g}] Several nonparametric regression estimators can be proved to be consistent for dependent samples via certain bandit algorithms yang2002.
table*[table* omitted — 3,358 chars of source]

Experiments

In this section, using benchmark datasets, we demonstrate the effectiveness of the proposed FA3IPW estimator and the FA3IPW estimator with discarding some early samples for stabilization (SFA3IPW estimator in Remark (ref)). Following dudik2011doubly and Chow2018, we evaluate the proposed estimators using classification datasets by transforming them into contextual bandit data. From the LIBSVM repository, we use the satImage, pendigits, mnist, dna, letter, sensorless, connect-4, and covtype datasets \footnote{\url{https://www.csie.ntu.edu.tw/ cjlin/libsvmtools/datasets/}}. First, we create an adaptive policy $\pi^a_t$ by following a random walk and perform upper confidence bound (UCB) algorithms in the MAB problem. Then, we construct a behavior policy as a mixture of $\pi^a_t$ and the uniform random policy $\pi^u_t$, defined as $\pi_t = 0.7\pi^a_t+0.3\pi^u_t$ Kallus2019IntrinsicallyES. When we use the RW policy, the policy $\pi^a_t$ does not converge. In contrast, when we use the UCB policy to generate $\pi^a_t$, the policy converges to a time-invariant policy. To construct an evaluation policy, we create a deterministic policy $\pi_d$ by training a logistic regression classifier on historical data. Then, we construct the evaluation policy $\pi^\mathrm{e}$ as a mixture of $\pi^d$ and the uniform random policy $\pi^u$, defined as $\pi^\mathrm{e} = 0.7\pi^d+0.3\pi^u$. Through experiments, the behavior policy $\pi_t$ is assumed to be known. More details are in Appendix (ref).

We compare the MSEs of six estimators, the FA3IPW, SFA3IPW, DM, AdaIPW, A2IPW, and FA2daIPW estimators. For the FA3IPW, SFA3IPW, and FA2daIPW estimators, we applied the two-step estimation. In each experiment, we have access to historical data with $T=1000$ and evaluation data with $N=1000$. With the SFA3IPW estimator, we discard the first $1000$ samples of the historical data. When estimating $f^*$ and $e^*$, we use the Nadaraya-Watson regression (NW) estimator, which is used in yang2002. The resulting MSEs and their standard deviations (SDs) over $20$ replications of each experiment are shown in Tables (ref) and (ref). In the experimental results, when the RW policy is used, the FA3IPW and SFA3IPW estimators show preferable performances compared with the A2IPW estimator. This is because the performance of the A2IPW estimator is guaranteed only when the policy converges to a time-invariant policy. On the other hand, when the UCB policy is used, the A2IPW estimator also works well as we can expect. From the results, we recommend using the FA3IPW estimator when the policy does not converge. More importantly, in the case, we can construct confidence intervals from the FA3IPW estimator, but cannot construct it from the A2IPW estimator. In Appendix (ref), we show the experimental results with other settings.

Conclusion

This study presented solutions for causal inference from dependent samples obtained via bandit feedback. By utilizing evaluation data or sample splitting, we can standardize the variance of an MDS. Then, we can apply the CLT for an MDS to obtain the asymptotic normality of the proposed estimator. In experiments, the proposed estimators showed theoretically expected performances.

Broader Impact

In various applications, such as A/B testing and clinical trials, the problem of causal inference from dependent samples arises. Because this study promotes the use of OPE for those applications by presenting an effective solution to this problem, the ethical issues related to such applications will become more critical. For example, in clinical trials, the proposed method also should take the criteria proposed by the US Food and Drug Administration fda into account.

Besides, although this study solved an essential problem with OPE, the limitations of methods for OPE should still be recognized. The existing methods for OPE try to find causality only from data without requiring domain-specific knowledge, but they sometimes lead to inappropriate decision making by producing incorrect results. When finding causality in the real world, it might be desirable to utilize domain-specific knowledge combined with the results obtained from machine learning methods.