Year: 2020

Abstract:
We consider off-policy evaluation (OPE) in continuous action domains, such as dynamic pricing and personalized dose finding. In OPE, one aims to learn the value under a new policy using historical data generated by a different behavior policy. Most existing works on OPE focus on discrete action domains. To handle continuous action space, we develop a brand-new deep jump Q-evaluation method for OPE. The key ingredient of our method lies in adaptively discretizing the action space using deep jump Q-learning. This allows us to apply existing OPE methods in discrete domains to handle continuous actions. Our method is further justified by theoretical results, synthetic and real datasets.

heading: 1 INTRODUCTION
text: Individualization proposes to leverage omni-channel data to meet individual needs. Individualized decision making plays a vital role in a wide variety of applications. Examples include customized pricing strategy in economics (Qiang & Bayati, 2016; Turvey, 2017), individualized treatment regime in medicine (Chakraborty, 2013; Collins & Varmus, 2015), personalized recommendation system in marketing (McInerney et al., 2018; Fong et al., 2018), etc. Prior to adopting any decision rule in practice, it is crucial to know the impact of implementing such a policy. In many applications, it is risky to run a policy online to estimate its value (see, e.g., Li et al., 2011). Off-policy evaluation (OPE) thus attracts a lot of attention by learning the policy value offline using logged historical data.
Despite the popularity of developing OPE methods with a finite set of actions (see e.g., Dudı́k et al., 2011; 2014; Swaminathan et al., 2017; Wang et al., 2017), less attention has been paid to continuous action domains, such as dynamic pricing (den Boer & Keskin, 2020) and personalized dose finding (Chen et al., 2016). Recently, a few OPE methods have been proposed to handle continuous actions (Kallus & Zhou, 2018; Sondhi et al., 2020; Colangelo & Lee, 2020). All these methods rely on the use of a kernel function to extend the inverse probability weighting (IPW) or doubly robust (DR) approaches developed in discrete action domains. They suffer from three limitations. First, the validity of these methods requires the conditional mean of the reward given the feature-action pair to be a smooth function over the action space. This assumption could be violated in applications such as dynamic pricing, where the expected demand for a product has jump discontinuities as a function of the charged price (den Boer & Keskin, 2020). Second, the value estimator could be sensitive to the choice of the bandwidth parameter in the kernel function. It remains challenging to select this hyperparameter. Kallus & Zhou (2018) proposed to tune this parameter by minimizing the mean squared error of the resulting value estimator. However, their method is extremely computationally intensive in moderate or high-dimensional feature space; see Section 5 for details. Third, these kernel-based methods typically use a single bandwidth parameter. This is sub-optimal in cases where the second-order derivative of the conditional mean function has an abrupt change in the action space; see the toy example in Section 3.1 for details.
To address these limitations, we develop a deep jump Q-evaluation (DJQE) method by integrating multi-scale change point detection (see e.g., Fryzlewicz, 2014), deep learning (LeCun et al., 2015) and OPE in discrete action domains. The key ingredient of our method lies in adaptively discretizing the action space using deep jump Q-learning. This allows us to apply IPW or DR methods to handle continuous actions. It is worth mentioning that our method does not require kernel bandwidth selection. Theoretically, we show it allows the conditional mean to be either a continuous or piecewise function of the action (Theorems 1 and 2) and converges faster than kernel-based OPE (Theorem 3). Empirically, we show it outperforms state-of-the-art OPE methods in synthetic and real datasets.

heading: 2 PRELIMINARIES
text: We first formulate the OPE problem. We next discuss the kernel-based OPE methods and multi-scale change point detection, since our proposal is closely related to them.

heading: 2.1 OFF-POLICY EVALUATION
text: The observed datasets can be summarized into {(Xi, Ai, Yi)}1≤i≤n where Oi = (Xi, Ai, Yi) denotes the feature-action-reward triplet for the ith subject and n denotes the total sample size. We assume these data triplets are independent copies of some population variables (X,A, Y ). Let X and A denote the feature and action space, respectively. We focus on the setting where A is one-dimensional, as in dynamic pricing and personalized dose finding. A deterministic policy π : X → A determines the action to be assigned given the observed feature. We use b to denote the behavior policy that generates the observed data. Specifically, b(•|x) denotes the probability density or mass function of A given X = x, depending on whether A is continuous or not. Define the expected reward function conditional on the feature-action pair as
Q(x, a) = E{Y |X = x,A = a}. We refer to this function as the Q-function, to be consistent with the literature on developing individualized treatment regime (Murphy, 2003).
As standard in the OPE and the causal inference literature (see e.g., Chen et al., 2016), we assume the stable unit treatment value assumption (SUTVA), no unmeasured confounders assumption, and the positivity assumption are satisfied. These assumptions guarantee that a policy’s value is estimable from the observed data. Specifically, for a given target policy π, its value can be represented by
V (π) = E{Q(X,π(X))}. The goal of the OPE is to learn the value V (π) based on the observed data.

heading: 2.2 KERNEL-BASED OPE
text: For discrete action, Zhang et al. (2012) and Dudı́k et al. (2011) proposed a DR estimator of V (π) by
1
n n∑ i=1 ψ(Oi, π, Q̂, b̂) = 1 n n∑ i=1
[ Q̂(Xi, π(Xi)) +
I(Ai = π(Xi)) b̂(Ai|Xi) {Yi − Q̂(Xi, π(Xi))}
] , (1)
where I denotes the indicator function, Q̂ and b̂ denote some estimators for the Q-function and the behavior policy. The second term b̂−1(Ai|Xi)I(Ai = π(Xi)){Yi − Q̂(Xi, π(Xi))} inside the bracket corresponds to an augmentation term. Its expectation equals zero when Q̂ = Q. The purpose of adding this term is to offer additional protection against potential model misspecification of the Q-function. Such an estimator is doubly-robust in the sense that its consistency relies on either Q̂ or b̂ to be correctly specified. By setting Q̂ = 0, equation 1 is reduced to the IPW estimator.
In continuous action domains, the indicator function I(Ai = π(Xi)) equals zero almost surely. Consequently, naively applying equation 1 yields the plug-in estimator ∑n i=1 Q̂(Xi, π(Xi))/n. To address this concern, the kernel-based OPE proposed to replace the indicator function in equation 1 with a kernel function K{(Ai − π(Xi))/h} with some bandwidth parameter h, i.e.,
1
n n∑ i=1 ψh(Oi, π, Q̂, b̂) = 1 n n∑ i=1
[ Q̂(Xi, π(Xi)) +
K{(Ai − π(Xi))/h} b̂(Ai|Xi) {Yi − Q̂(Xi, π(Xi))}
] .
The bandwidth h represents a trade-off. The variance of the resulting value estimator decays with h. Yet, its bias increases with h. More specifically, it follows from Theorem 1 of Kallus & Zhou (2018) that the leading term of the bias is equal to
h2 ∫ u2K(u)du
2 E
( ∂2Q(X, a)
∂a2
∣∣∣∣ a=π(X) ) . (2)
To ensure the term in 2 decays to zero as h goes to 0, it requires the expected second derivative of the Q-function to exist, and thus Q(x, a) needs to be a smooth function of a. However, as commented in the introduction, this assumption could be violated in applications such as dynamic pricing.

heading: 2.3 MULTI-SCALE CHANGE POINT DETECTION
text: The change point analysis considers an ordered sequence of data, Y1:n = {Y1, · · · , Yn}, with unknown change point locations, τ = {τ1, · · · , τK} for some unknown integer K. Here, τi is an integer between 1 and n− 1 inclusive, and satisfies τi < τj for i < j. These change points split the data intoK+1 segment. Within each segments, the expected response is a constant function (see the left panel of Figure 1 for details). A number of methods have been proposed on estimating change points (see for example, Boysen et al., 2009; Killick et al., 2012; Frick et al., 2014; Fryzlewicz, 2014, and the references therein), by minimizing a penalized objective function:
arg min τ,K
( 1
n K+1∑ i=1 [ C{Y(τi−1+1):τi} ] + γK ) ,
where C is a cost function that measures the goodness-of-the-fit of the constant function within each segment and γK penalizes the number of change points with some regularization parameter γ. We remark that all the above cited works focused on models without features. Our proposal goes beyond these works in that we consider models with features and use deep neural networks (DNN) to capture the complex relationship between the response and features.

heading: 3 DEEP JUMP Q-EVALUATION
text: In section 3.1, we use a toy example to demonstrate the limitation of kernel-based methods. We present the main idea of our algorithm in Section 3.2. Details are given in Section 3.3.

heading: 3.1 TOY EXAMPLE
text: As discussed in the introduction, existing kernel-based OPE methods use a single bandwidth to construct the value estimator. Ideally, the bandwidth h in the kernel K{(Ai− π(Xi))/h} shall vary with π(Xi) to improve the accuracy of the value estimator. To elaborate this, consider the Q-function Q(x, a) = 10 max{a2 − 0.25, 0} log(x + 2) for any x, a ∈ [0, 1]. By definition, the Q-function is smooth over the entire feature-action space. However, it has different “patterns” when the action belongs to different intervals. Specifically, for a ∈ [0, 0.5], Q(x, a) is constant as a function of a. For a ∈ (0.5, 1], Q(x, a) depends quadratically in a. See the middle panel of Figure 1 for details.
Consider the target policy π(x) = x. We decompose the value V (π) into V (1)(π) +V (2)(π) where
V (1)(π) = EQ(X,π(X))I(π(X) ≤ 0.5) and V (2)(π) = EQ(X,π(X))I(π(X) > 0.5).
Similarly, denote the corresponding kernel-based value estimators by
V̂ (1) h (π) =
1
n n∑ i=1 ψh(Oi, π, Q̂, b̂)I(π(Xi) ≤ 0.5) and V̂ (2)h (π) = 1 n n∑ i=1 ψh(Oi, π, Q̂, b̂)I(π(Xi) > 0.5).
Since Q(x, a) is a constant function of a ∈ [0, 0.5], its second-order derivative ∂2Q(x, a)/∂a2 equals zero. In view of 2, when π(x) ≤ 0.5, the bias of V̂ (1)h (π) will be small even with a sufficiently large h. As such, a large h is preferred to reduce the variance of V̂ (1)h (π). When π(x) > 0.5, a small h is preferred to reduce the bias of V̂ (2)h (π). See Table 1 for details where we report the bias and standard deviation of V̂ (1)h (π) and V̂ (2) h (π) with two different bandwidths.
Due to the use of a single bandwidth, the kernel-based estimator suffers from either a large bias or a large variance. To overcome this limitation, we propose to adaptively discretize the action space into a union of disjoint intervals such that within each interval I, the Q-function {Q(x, a) : a ∈ I} can be well-approximated by some functionQI(x) that is constant in a ∈ I. Based on the discretization, one can apply IPW or DR to evaluate the value. The advantage of adaptive discretization is illustrated in the right panel of Figure 1. When a ≤ 0.5, the Q-function is constant in a. It is likely that our procedure will not further split the interval [0, 0.5]. Consequently, the corresponding DR estimator for V (1)(π) will not suffer from large variance. When a > 0.5, our procedure will split (0.5, 1] into a series of sub-intervals, approximating Q by a step function. This guarantees the resulting DR estimator for V (2)(π) will not suffer from large bias. Consequently, the proposed value estimator achieves a smaller mean squared error than kernel-based estimators. See Table 1 for details.

heading: 3.2 THE MAIN IDEA
text: For simplicity, we set the action space A = [0, 1]. From now on, we focus on a subset of intervals in [0, 1]. By interval we always refer to those of the form [c, d) for some 0 ≤ c < d < 1, or [c, 1] for some 0 ≤ c < 1, denoted as I. A discretization D for A is defined as a collection of mutually disjoint intervals that coversA. Let |D| denote the number of intervals inD and |I| denote the length of the interval I. We aim to identify an “optimal” discretization D̂ such that for each interval I ∈ D̂, Q(x, a) is approximately a constant function of a ∈ I. The number of intervals in D̂ represents a trade-off. If |D̂| is too large, then D̂ will contain many short intervals, the resulting IPW or DR estimator might suffer from large variance. Yet, a smaller value of |D̂| might result in a large bias. Our proposed method adaptively determines D̂ and its size |D̂| as illustrated below. To begin with, we cut the entire action space A into m initial intervals: [0, 1/m),[1/m, 2/m), . . . , [(m − 1)/m, 1]. The number m shall be sufficiently large such that the Q-function can be wellapproximated by a piecewise function on these intervals. In practice, we recommend to set the initial number of intervals m to be proportional to the sample size n. Note the set of these initial intervals is not the final partition D̂ that we recommend, but only serve as the initial candidate intervals. We next adaptively combine some of these initial intervals to form the final partition D̂. As shown in our numerical studies (see Table 5 in Appendix B for more details), the size of the final partition |D̂| is usually much less than m. More specifically, denote by B(m) as the set of discretizations D such that the end-points of each interval I ∈ D lie on the grid {j/m : j = 0, 1, · · · ,m}. We associate to each partition D ∈ B(m) a collection of functions {QI}I∈D. These functions depend only on features, not the action. They are used to produce a piecewise approximation of the Q-function such that Q(a, •) ≈ ∑ I I(a ∈ I)QI(•). We model these QI using deep neural networks, to capture the complex dependence between the response and features. When the Q-function is well-approximated, we expect the least square loss ∑ I∈D [∑n i=1 I(Ai ∈ I) { Yi −QI(Xi) }2] , will be small.
Consequently, D̂ can be estimated by solving
(D̂, {q̂I : I ∈ D̂}) = arg min (D∈B(m),{QI∈Q:I∈D}) (∑ I∈D [ 1 n n∑ i=1 I(Ai ∈ I) { Yi −QI(Xi) }2] + γ|D| ) , (3)
for some regularization parameter γ and DNN class Q. Here, the penalty term γ|D| in equation 3 controls the total number of intervals in D̂, as in multi-scale change point detection. A large γ results in few intervals in D̂ and a potential large bias of the value estimator, whereas a small γ
could procedure a large number of intervals in D̂, leading to a noisy value estimator. In practice, we use cross-validation to select the regularization parameter γ that minimizes the mean square error of the fitted Q-function. We refer to this step as deep jump Q-learning. Details of this step are given in the next section. Given D̂, one can apply IPW or DR (see equation 1) to derive the value estimates. To further reduce the bias of the value estimator, we employ a data splitting and cross-fitting strategy, which is commonly used in statistics (Romano & DiCiccio, 2019). That is, we use different subsets of data samples to learn the discretization D̂ and to construct the value estimator. A pseudocode summarizing our algorithm is given in Algorithm 1 in Appendix A. We present the details below.

heading: 3.3 THE COMPLETE ALGORITHM
text: We present the details for DJQE in this section. It consists of three steps: data splitting, deep jump Q-learning, and cross-fitting.
Step 1: Data Splitting: We divide all n samples into L subsets of equal size, where L` denotes the indices of samples in the `th subset for ` = 1, · · · ,L. Let Lc` = {1, 2, · · · , n} − L` as the complement of L`.
Step 2: Deep Jump Q-Learning: For each ` = 1, · · · ,L, we propose to apply deep jump Qlearning to compute a discretization D̂(`) and {q̂(`)I : I ∈ D̂(`)} by solving a version of equation 3 using the data subset in Lc` only. We next present the computational details for solving this optimization. Our approach is motivated by the PELT method (Killick et al., 2012) in multi-scale change point detection. Specifically, for any interval I, define q̂(`)I as the minimizer of
arg min QI∈Q
1 |Lc`| ∑ i∈Lc` I(Ai ∈ I) { QI(Xi)− Yi }2 , (4)
where |Lc`| denotes the number of samples in Lc`. Define the cost function C(`)(I) as the minimum value of the objective function 4, i.e, C(`)(I) = 1|Lc` | ∑ i∈Lc` I(Ai ∈ I) { q̂ (`) I (Xi) − Yi }2 . In our implementation, we set Q as the class of multilayer perceptrons (MLPs, see Figure 3 in Appendix A for an illustration) with L hidden layers and H nodes each hidden layer. The above optimization can be solved via the MLP regressor implementation of Pedregosa et al. (2011).
Computation of D̂(`) relies on dynamic programming (Friedrich et al., 2008). For any integer 1 ≤ v∗ < m, denote by B(m, v∗) the set consisting of all possible discretizations Dv∗ of [0, v∗/m). Set B(m,m) = B(m), we define the Bellman function as
Bell(v∗) = inf Dv∗∈B(m,v∗) ( ∑ I∈Dv∗ C(`)(I) + γ(|Dv∗ | − 1) ) , and Bell(0) = −γ.
Our algorithm recursively updates the Bellman function based on the following formula,
Bell(v∗) = min v∈Rv∗
{ Bell(v) + C(`)([v/m, v∗/m)) + γ } , ∀v∗ ≥ 1, (5)
whereRv∗ is the candidate change points list updated by
{v ∈ Rv∗−1 ∪ {v∗ − 1} : Bell(v) + C(`)([v/m, (v∗ − 1)/m)) ≤ Bell(v∗ − 1)}, (6) during each iteration with R0 = {0}. The constraint listed in 6 is important as it facilitates the computation by discarding change points not relevant to obtain the final discretization, leading to a linear computational cost (Killick et al., 2012).
To solve equation 5, we search the optimal change point location v that minimizes Bell(v∗). This requires to apply the MLP regressor to learn q̂(`)[v/m,v∗/m) and C
(`)([v/m, v∗/m)) for each v ∈ Rv∗ . Let v1 be the corresponding minimizer. We then define the change points list τ(v∗) = {v1, τ(v1)}. This procedure is iterated to compute Bell(v∗) and τ(v∗) for v∗ = 1, . . . ,m. The optimal partition D̂(`) is determined by the values stored in τ (see Algorithm 1 in Appendix A for details).
Step 3. Cross-Fitting: For each interval in the estimated optimal partition D̂(`), we estimate the generalized propensity score function Pr(A ∈ I|X = x) via the MLP regressor using the training dataset Lc`. Let b̂(`)(I|x) denote the resulting estimate. The final estimated value for V (π) is
constructed via cross-fitting, given by,
V̂ (π) = 1
n L∑ `=1 ∑ I∈D̂(`) ∑ i∈L`
[ I(Ai ∈ I)
I{π(Xi) ∈ I} b̂ (`) I (I|Xi)
{ Yi − q̂(`)I (Xi) } + I(Ai ∈ I)q̂(`)I (Xi) ] . (7)
Note the samples used to construct V̂ inside bracket are independent from those to estimate q̂(`)I , b̂ (`) I and D̂(`). This helps remove the bias induced by overfitting in the estimation of q̂(`)I , b̂ (`) I and D̂(`).

heading: 4 THEORY
text: We investigate the theoretical properties of the proposed estimator. For simplicity, let the support X = [0, 1]p. We will show our estimator is consistent when the Q-function is either a piecewise function or a continuous function of a. Specifically, consider the following two model assumptions.
Model 1 (Piece-wise constant function). Suppose Q(x, a) = ∑ I∈D0 qI(x)I(a ∈ I), ∀x ∈ X , ∀a ∈ A, (8)
for some partition D0 of [0, 1] and a collection of continuous functions (qI)I∈D0 . Model 2 (Continuous function). Suppose Q is a continuous function of a and x.
We first consider the case where the value function takes the form of equation 8. We remark that kernel-based estimators will fail under this model assumption, as it requires the second-order derivative of the Q-function to exist. Without loss of generality, assume qI1 6= qI2 for any two adjacent intervals I1, I2 ∈ D0. This guarantees that the representation in equation 8 is unique. For any partition D = {[0, τ1), [τ1, τ2), · · · , [τK , 1]}, we use J(D) to denote the set of change points {τ1, · · · , τK}. We impose the following conditions to establish our theories.
Assumption 1. The number of layers L and the number of nodes in each hidden layer H diverge with n, in that HL = O(nρ), for some constant 0 < ρ < 1/2.
Assumption 2. Functions {q̂(`)I (•)}I∈D̂(`) are uniformly bounded.
Assumption 1 is mild, as both H and L are parameters we specify. The part that HL = O(nρ) ensures that the stochastic error resulting from the parameter estimation in the MLP is negligible. Assumption 2 ensures that the optimizer would not diverge in the `∞ sense. Similar assumptions are commonly imposed in the literature to derive the convergence rates of DNN estimators (see e.g., Farrell et al., 2018). These two assumptions guarantee the uniform consistency of the DNN estimator {q̂I}I∈D̂(`) . See Lemma 1 in Appendix D for details.
Theorem 1 Suppose Model 1, Assumptions 1 and 2 hold. Suppose m diverges to infinity with n. Then, there exists some constant γ0 such that as long as 0 < γ ≤ γ0, the following events occur with probability approaching 1 (w.p.a.1),
(i) |D̂(`)| = |D0|; (ii) maxτ∈J(D0) minτ̂∈J(D̂(`)) |τ̂ − τ | = op(1). (iii) V̂ (π) = V (π) + op(1) for any policy π such that for any τ0 ∈ J(D0), Pr(π(X) ∈ [τ0 − , τ0 + ])→ 0 as → 0.
Theorem 1 establishes the properties of our method under settings where the Q(a, x) is piecewise function in a. Results in (i) imply that deep jump Q-learning correctly identifies the number of change points. Results in (ii) imply that any change point in D0 can be consistently identified. In particular, D̂(`) corresponds to a subset of {1/m, 2/m, · · · , (m − 1)/m}. With a sufficiently large m, for any true change point τ in D0, there will be a change point in D̂(`) that approaches τ . Consequently, the change point locations can be consistently estimated. To ensure the consistency of the proposed value estimator, we require that the distribution of the random variable π(X) induced by the target policy does not have point-masses at the change point locations. This condition is also mild. For instance, it automatically holds when π(X) has a density function on [0, 1].
Theorem 2 Suppose Model 2, Assumptions 1 and 2 hold. Suppose m diverges to infinity and γ decays to zero. Then we have
(i) maxI∈D̂(`) supa∈I E|q̂ (`) I (X)−Q(X, a)|2 = o(1); (ii) V̂ (π)− V (π) = op(1) for any π.
Theorem 2 establishes the properties of our method under settings where Q is continuous in a. Results in (i) imply that q̂(`)I (•) can be used to uniformly approximate Q(a, •) for any a ∈ I. The consistency of the value in (ii) thus follows.
Finally, we conduct in-depth theoretical analysis to demonstrate the advantage of our estimator. Due to space constraint, we present an informal statement below. Details are given in Appendix C.
Theorem 3 (Informal Statement) (i) When the Q-function belongs to a class of piecewise constant functions of the action a, the minimax rate of convergence of the proposed value estimator isOp(n−1/2). In contrast, the minimax convergence rate of kernel-based estimator is Op(n−1/3).
(ii) When the Q-function belongs to a class of smooth functions of a, the minimax convergence rate of our estimator is faster than kernel-based estimators if the bandwdith undersmoothes or oversmoothes the data.

heading: 5 EXPERIMENTS
text: In this section, we investigate the finite sample performance of the proposed DJQE on the synthetic and real datasets, in comparison to two kernel-based methods. The computing infrastructure used is a virtual machine in the AWS Platform with 72 processor cores and 144GB memory.

heading: 5.1 SYNTHETIC DATA
text: Synthetic data are generated from the following model:
Y |X,A ∼ N{Q(X,A), 1}, b(A|X) ∼ Unif[0, 1] and X(1), X(2), . . . , X(p) iid∼ Unif[−1, 1],
where X = [X(1), X(2), . . . , X(p)]. Consider the following different scenarios: S1: Q(x, a) = (1 + x(1))I(a < 0.35) + (x(1) − x(2))I(0.35 ≤ a < 0.65) + (1− x(2))I(a ≥ 0.65); S2: Q(x, a) = I(a < 0.25) + sin(2πx(1))I(0.25 ≤ a < 0.5){0.5 − 8(x(1) − 0.75)2}I(0.5 ≤ a < 0.75) + 0.5I(a ≥ 0.75); S3 (toy): Q(x, a) = 10 max{a2 − 0.25, 0} log(x(1) + 2); S4: Q(x, a) = 0.2(8 + 4x(1) − 2x(2) − 2x(3))− 2(1 + 0.5x(1) + 0.5x(2) − 2a)2. The Q-function is a piecewise function of a under Scenarios 1 and 2, and is continuous under Scenarios 3 (toy example considered in Section 3.1) and 4. We set the target policy to be the optimal policy that achieves the highest possible mean reward. We list the oracle mean value under the optimal policy for each scenario in the first column of Table 4 in Appendix B.
We apply the proposed DJQE and two kernel-based methods (Kallus & Zhou, 2018; Colangelo & Lee, 2020) to Scenario 1-4 with 20-dimensional covariates (p = 20) and n ∈ {50, 100, 200, 300}. For the DJQE, we select γ ∈ {0.1, 0.2, 0.3, 0.4, 0.5}n0.4 based on five-fold cross-validation. Here, we setm = n/10 to achieve a good balance between the bias and the computational cost (see Figure 4 in Appendix B for the detailed computational cost of the DJQE and the resulting bias as a function of m in Scenario 1 with n = 100). We find it extremely computationally intensive to compute the optimal bandwidth h∗ using Kallus & Zhou (2018)’s method (see the detailed comparison of computational cost under different methods based on Scenario 1 in Table 3). Thus, as suggested in Kallus & Zhou (2018), we first compute h∗ using data with sample size n0 = 50. To accommodate data with different sample sizes n, we adjust h∗ by setting h∗{n0/n}0.2. To implement Colangelo
& Lee (2020)’s estimator, we consider a list of bandwidths, given by h = cσAn−0.2 with c ∈ {0.5, 0.75, 1.0, 1.5} and σA is the sample standard deviation of the action. We then manually select the best bandwidth such that the resulting value estimator achieves the smallest mean squared error. The conditional mean value and generalized propensity score are fitted via the MLP regressor with 10 hidden layer and 10 neurons in each layer. The average estimated value and its standard deviation over 100 replicates are illustrated in Figure 2 for different methods, with detailed values reported in Table 4 in Appendix B. In addition, we provide the size of the final estimated partition under the DJQE in Table 5 in Appendix B, which is much smaller than m in most cases.
It can be seen from Figure 2 that the proposed DJQE is very efficient and outperforms all competing methods in almost all cases. We note that the proposed method performs reasonably well even when the sample size is small (n = 50). In contrast, two kernel-based methods fail to accurately estimate the value even in some cases when n = 300. Among the two kernel-based OPE approaches, we observe that the method developed by Kallus & Zhou (2018) performs better in general.

heading: 5.2 REAL DATA: PERSONALIZED DOSE FINDING
text: Warfarin is commonly used for preventing thrombosis and thromboembolism. We use the dataset provided by the International Warfarin Pharmacogenetics (Consortium, 2009) for analysis. We choose p = 81 baseline covariates considered in Kallus & Zhou (2018). This yields a total of 3964 with complete records of baseline information. The response is defined as the absolute distance between the international normalized ratio (INR, a measurement of the time it takes for the blood to clot) after the treatment and the ideal value 2.5, i.e, Y = −|INR−2.5|. We use the min-max normalization to convert the range of the dose level A into [0, 1].
To compare among different methods, we calibrate the dataset to generate simulated outcomes. Specifically, we first estimate the Q-function via a MLP regressor with 10 hidden layers and 50 neurons in each layer using the whole dataset. The goodness-of-the-fit of the fitted model under the MLP regressor is reported in Table 6 in Appendix B. We next use the fitted Q-function Q̂(X,A) to simulate the data. Given a randomly sampled feature-action pair (aj , xj) from {(A1, X1), · · · , (An, Xn)}, we set the reward rj toN{Q̂(xj , aj), σ̂2}, where σ̂ is the standard deviation of the fitted residual {Yi − Q̂(Xi, Ai)}i. Given the simulated data {(xj , aj , rj) : 1 ≤ j ≤ n}, we are interested in evaluating the optimal policy: π?(X) ≡ arg maxa∈[0,1] Q̂(X, a). The oracle value under the optimal policy is V = −0.278. We apply the DJQE on the calibrated Warfarin data, against two kernel-based methods. Due to the extremely intensive computation in Kallus & Zhou (2018), we directly apply the estimated optimal bandwidth h∗ in their real data analysis, since they used the same dataset. Biases, standard deviations, and mean squared errors of the estimated values are reported in Table 2 over 20 replicates for sample size n = 500 under different methods.
It can be observed from Table 2 that our proposed DJQE achieves much smaller mean squared error than the two kernel-based methods, when evaluating the optimal policy. Specifically, the DJQE yields bias as 0.259 with the standard deviation as 0.416, in contrast to the large bias as 0.662 with the standard deviation as 0.742 under Kallus & Zhou (2018)’s method and the bias as 0.442 with the large standard deviation as 1.164 under Colangelo & Lee (2020)’s method. Therefore, our proposed DJQE for off-policy evaluation with continuous actions works better than the kernel-based methods.

heading: 6 DISCUSSION
text: Currently, we focus on settings with a single decision point. It would be practically interesting to extend our proposal to sequential decision making. A potential drawback of our method is that it would be very computationally intensive for a large m as the runtime increases linearly in m.

heading: A MORE ON THE IMPLEMENTATION
text: We summarize our algorithm in Algorithm 1.
Global: data {(Xi, Ai, Yi)}1≤i≤n; number of initial intervals m; penalty term γ; target policy π. Local: an upper triangular matrix of cost C ∈ Rm(m+1)/2; Bellman function Bell ∈ Rm; partitions D̂;
DNN functions {q̂I , b̂I : I ∈ D̂}; a vector τ ∈ Nm; a set of candidate point listsR. Output: the value estimator for target policy V̂ (π). I. Split all n samples into L subsets as {L1, · · · ,LL}; V̂ (π)← 0; II. Initialize an even segment on the action space with m pieces: {I} = {[0, 1/m), [1/m, 2/m), . . . , [(m− 1)/m, 1]}; III. For ` = 1, · · · , L:
1. Set the training dataset as Lc` = {1, 2, · · · , n} − L`; 2. Bell(0)← −γ; D̂ = [0, 1]; τ ← Null;R(0)← {0}; 3. Collect cost function: For r = 1, . . . ,m: For l = 0, . . . , (r − 1):
(i). Let I = [l/m, r/m) if r < m else I = [l/m, 1]; (ii). Fit a MLP regressor: q̂I(·)← I(i ∈ Lc`)I(Ai ∈ I)Yi ∼ I(Ai ∈ I)MLP (Xi); (iii). Calculate the cost: C(I)← ∑ i∈Lc` I(Ai ∈ I) { q̂I(Xi)− Yi }2 ;
4. Apply the PELT method to get partitions: For v∗ = 1, . . . ,m: (i).Bell(v∗) = minv∈R(v∗){Bell(v) + C([v/m, v∗/m)) + γ}; (ii). v1 ← arg minv∈R(v∗){Bell(v) + C([v/m, v∗/m)) + γ}; (iii). τ(v∗)← {v1, τ(v1)}; (iv). R(v∗)← {v ∈ R(v∗ − 1) ∪ {v∗ − 1} : Bell(v) + C([v/m, (v∗ − 1)/m)) ≤ Bell(v∗ − 1)}; 5. Construct the DR value estimator: r ← m; l← τ [r]; While r > 0: (i) Let I = [l/m, r/m) if r < m else I = [l/m, 1]; D̂ ← D̂ ∪ I; (ii) Recall fitted MLP: q̂I(·)← I(i ∈ Lc`)I(Ai ∈ I)Yi ∼ I(Ai ∈ I)MLP (Xi); (iii) Fit propensity score: b̂I(·)← I(i ∈ Lc`)I(Ai ∈ I) ∼ I(Ai ∈ I)MLP (Xi); (iv) r ← l; l← τ(r); 6. Evaluation using testing dataset L`: V̂ (π)+ = ∑ I∈D̂ (∑ i∈L` I(Ai ∈ I) [ I{π(Xi)∈I} b̂I(I|Xi) { Yi − q̂I(Xi) } + q̂I(Xi) ]) ;
return V̂ (π)/n . Algorithm 1: Deep Jump Q-Evaluation

heading: B ADDITIONAL EXPERIMENTAL RESULTS
text: We include additional experimental results in this section.

heading: C ADDITIONAL THEORETICAL RESULTS
text: In this section, we conduct in-depth theoretical analysis to compare the minimax convergence rate of the proposed value estimator with kernel-based value estimator.
We first briefly summarize our theoretical findings. When the Q-function belongs to a class of piecewise constant functions of the action a, the proposed estimator converges at a faster rate than kernel-based estimators. Specifically, the kernel-based OPE converges at a rate of n−1/3 whereas our estimator converges at a rate of n−1/2.
When the Q-function belongs to a class of Lipschitz continuous functions of a, our estimator converges at a rate of n−1/5 whereas kernel-based estimators converge at a slower rate when the bandwdith undersmoothes or oversmoothes the data.
As we have commented, it remains extremely challenging to tune the bandwidth of kernel-based OPE in practice, since OPE is an unsupervised problem. The use of a suboptimal bandwidth would lead to a slower rate of convergence for kernel-based estimators. In contrast, we develop a supervised learning algorithm to adaptively discretize the action space. The tuning parameters in our procedure can be selected via cross-validation.
All the above findings are supported by our observations in the toy example (Section 3.1) and numerical experiments (Section 5) as well.
We next present our theoretical results. We first introduce some notations. Define the following classes of Q-functions:
Q1 = Q : Q(a, x) = ∑ I∈D0 ∑ a∈I I(a ∈ I)qI(x), max I1 6=I2 I1,I2∈D0 E|qI1(X)− qI2(X)|2 ≥ ,
|D0| ≤ C1, Q(a, ·) ∈ Λ(β,C2),∀a} ,
(9)
Q2 = { Q : sup
a1,a2,x |Q(a1, x)−Q(a2, x)| ≤ C3, Q(a, ·) ∈ Λ(β,C2),∀a
} ,
for some sufficiently small constant > 0 and some sufficiently large constants C1, C2 > 0 where Λ(β,C2) denotes the class of β-smooth functions (see e.g., Section 3.1.1, Shi et al., 2020). By definition, the first set Q1 consists of all piecewise constant functions of a with at most C1 − 1 change points. The second set Q2 contains the class of Lipschitz continuous functions of a. Finally, define the policy class Π = {π : π(X) has a density function bounded by C4} for some sufficiently large constant C4 > 0. To simplify the analysis, we assume the behavior policy is
5MSE = 1 n ∑n i=1(Yi − Ŷi)
2. See https://en.wikipedia.org/wiki/Mean_squared_ error.
6NRMSD = √ MSE
max(Y )−min(Y ) . See https://en.wikipedia.org/wiki/ Root-mean-square_deviation.
7MAE = 1 n ∑n i=1 |Yi − Ŷi|. See https://en.wikipedia.org/wiki/Mean_absolute_
error. 8NMAE = MAE
max(Y )−min(Y ) . See https://en.wikipedia.org/wiki/Root-mean-square_ deviation.
known such that b̂ = b̂(`)I = b. In addition, assume the conditional variance of Y given A and X is uniformly bounded away from zero.
Theorem 3 Suppose 2β > p. Then for any π ∈ Π and any Q ∈ Q1, with proper choice of γ and the MLP class, the minimax rate of convergence of the proposed value estimator is Op(n−1/2). In contrast, the minimax convergence rate of kernel-based estimator is Op(n−1/3).
Suppose 4β > 3p. For any Q ∈ Q2, with proper choice of γ and the MLP class, the minimax rate of convergence of the proposed value estimator is Op(n−1/5), up to some logarithmic factors. In contrast, the minimax convergence rate of kernel-based estimator is slower when h n−κ for some κ < 1/5 or κ > 3/5.
We briefly comment on the condition 2β > p. This condition allows that the MLP regressor to converge at a rate faster than n−1/4 (see e.g., Imaizumi & Fukumizu, 2019). Such a condition is commonly assumed in the literature on evaluating average treatment effects (see e.g., Chernozhukov et al., 2017; Farrell et al., 2018). In the second part of Theorem 3, we further require the MLP regressor to converge at a rate faster than n−3/10.

heading: D TECHNICAL PROOF
text: Throughout the proof, we use c, C, c0, c̄, etc., to denote some universal constants whose values are allowed to change from place to place. Let Oi = {Xi, Yi} denote the data summarized from the ith observation. For any interval I, let b(I|x) denote the integral ∫ a∈I b(a|x)da.
Proofs of Theorems 1 and 2 rely on Lemmas 1 and 2. Specifically, Lemma 1 establishes the uniform convergence rate of q̂I for any I whose length is no shorter than o(γ) and belongs to the set of intervals:
I(m) = {[i1/m, i2/m) : for some integers i1 and i2 that satisfy 0 ≤ i1 < i2 < m} ∪ {[i3/m, 1] : for some integers i3 that satisfy 0 ≤ i3 < m}.
To state this lemma, we first introduce some notations. For any such interval I, define the function qI(x) = E(Y |A ∈ I, X = x). It is immediate to see that the definition of qI here is consistent with the one defined in equation 8 for any I ⊆ D0.
Lemma 1 The following holds when either the conditions in Theorem 1 or 2 hold:
max I∈I(m),|I|≥cγ
E[|qI(X)− q̂(`)I (X)| 2|{Oi}i∈L` ] = op(1),
for any positive constant c > 0. x
Lemma 2 When either the conditions in Theorem 1 or 2 hold, minI∈D̂(`) |I| ≥ Cγ w.p.a.1 for some constant C > 0.
We first present the proofs for these two lemmas. Next we present the proofs for Theorems 1 and 2.
D.1 PROOF OF LEMMA 1
For any sufficiently small constant > 0, it suffices to show
max I∈I(m),|I|≥cγ
E[|qI(X)− q̂(`)I (X)| 2|{Oi}i∈L` ] ≤ C 2,
w.p.a.1, for some constant C > 1 whose value will be specified later. The proof is divided into two parts. In the first part, we show there exist some constants {θ∗I}I such that
max I∈I(m) max x∈X |qI(x)−MLP(x; θ∗I)| ≤ . (10)
In the second part, we show
max I∈I(m),|I|≥cγ
E[|MLP(X; θ∗I)− q̂ (`) I (X)| 2|{Oi}i∈L` ] ≤ (C − 1) 2. (11)
The proof is hence completed.
Part 1. Under the step model assumption, for any I0 ∈ D0, qI0(•) is continuous. Since the support X is compact, qI0(•) is uniformly continuous. Similarly, we can show b is uniformly continuous as well. For any I ∈ Im, we have
qI(x) =
∫ a∈I ∑ I0∈D0 b(a|x)I(a ∈ I0)qI0(x)da
b(I|x) .
It follows from the uniform continuity of qI0 , b and the positivity assumption that qI(•) is continuous for any I ∈ I(m). Similarly, under the model assumption in Theorem 2, we can show qI(•) is continuous for any I ∈ I(m) as well. Consequently, when either the conditions in Theorem 1 or 2 hold, we obtain that qI(•) is continuous for any I ∈ I(m). By Stone-Weierstrass theorem, there exists a multivariate polynomial function q∗I such that the absolute value of the residual qI−q∗I is uniformly bounded by /2. By Theorem 1 of Yarotsky (2017), there exists a feedforward neural network with a bounded number of hidden units that uniformly approximates q∗I , with the approximation error uniformly bounded by /2 in absolute value. By Lemma 1 of Farrell et al. (2018), such a feedforward network can be embedded into an MLP with a bounded number of hidden units. Since we allow H and L to diverge, such an MLP can be further embedded into an MLP with L layers and the widths of all layers being proportional to H . This yields equation 10. The proof for Part 1 is thus completed.
Part 2. We aim to show equation 11 holds. Under the boundedness assumption on Y , {qI(•)}I∈I(m) are uniformly bounded. We first observe that MLP(•, θ∗I) is a bounded function, by equation 10 and the fact that qI(•) is bounded. Next, it follows from equation 10 that
EI(A ∈ I){Y −MLP(X; θ∗I)}2
Pr(A ∈ I) =
EI(A ∈ I){Y − E(Y |A ∈ I, X)}2
Pr(A ∈ I)
+ EI(A ∈ I){E(Y |A ∈ I, X)−MLP(X; θ∗I)}2 Pr(A ∈ I) ≤ EI(A ∈ I){Y − E(Y |A ∈ I, X)} 2 Pr(A ∈ I) + 2.
(12)
By definition,∑ i∈L` I(Ai ∈ I){Yi −MLP(Xi; θ ∗ I)}2 Pr(A ∈ I) ≥ ∑ i∈L` I(Ai ∈ I){Yi −MLP(Xi; θ̂ (`) I )}2 Pr(A ∈ I) . (13)
Suppose we can show
sup I∈I(m),|I|≥cγ
θI∈Θ∗
∣∣∣∣∣ ∑ i∈L` I(Ai ∈ I){Yi −MLP(Xi; θI)} 2 |L`|Pr(A ∈ I) − EI(A ∈ I){Y −MLP(X; θI)} 2 Pr(A ∈ I) ∣∣∣∣∣ = op(1), (14)
where the set Θ∗ consists of all θI such that supx |MLP(x; θI)| ≤ M for some sufficiently large constant M such that both θ∗I and θ̂ (`) I satisfy this constraint. It thus follows from equation 12 and equation 13 that
EI(A ∈ I){Y −MLP(X; θ̂(`)I )}2 Pr(A ∈ I) ≤ EI(A ∈ I){Y − E(Y |A ∈ I, X)} 2 Pr(A ∈ I) + 2 2,
w.p.a.1. Similar to equation 12, we have
EI(A ∈ I){Y −MLP(X; θ̂(`)I )}2
Pr(A ∈ I) =
EI(A ∈ I){Y − E(Y |A ∈ I, X)}2
Pr(A ∈ I)
+ EI(A ∈ I){E(Y |A ∈ I, X)−MLP(X; θ̂(`)I )}2
Pr(A ∈ I) .
It follows that
EI(A ∈ I){E(Y |A ∈ I, X)−MLP(X; θ̂(`)I )}2
Pr(A ∈ I) =
Eb(I|X){E(Y |A ∈ I, X)−MLP(X; θ̂(`)I )}2
Pr(A ∈ I) ≤ 2 2.
Under the positivity assumption, supx b(I|x)/Pr(A ∈ I) ≥ 2(c − 1)−1 for some constant c > 1. This yields equation 11.
To complete the proof, it remains to show equation 14 holds.
Using similar arguments in Section A.2.2 of Farrell et al. (2018), for any I, we have with probability at least 1− 2 exp(−γ̄) that∣∣∣∣∣ ∑ i∈L` I(Ai ∈ I){Yi −MLP(Xi; θI)} 2 |L`|Pr(A ∈ I) − EI(A ∈ I){Y −MLP(X; θI)} 2 Pr(A ∈ I) ∣∣∣∣∣ ≤ c̄ { Mγ̄
n +
√ γ̄
n σ(I, θI) + σ(I, θI)
√ H2L2
n
( log M
σ(I, θI) + log n
)} ,
for some constant c̄ > 0, where σ2(I, θ) corresponds to the variance of 1
Pr(A ∈ I) I(A ∈ I){Y −MLP(X; θ)}2.
Under the positivity assumption, we have σ2(I, θ) ≤ O(1)|I|−1 where O(1) denotes some positive constant. Consequently, for any I whose length is greater than cγ−1, we have with probability at least 1− 2 exp(−γ̄) that∣∣∣∣∣ ∑ i∈L` I(Ai ∈ I){Yi −MLP(Xi; θI)} 2 |L`|Pr(A ∈ I) − EI(A ∈ I){Y −MLP(X; θI)} 2 Pr(A ∈ I)
∣∣∣∣∣ ≤ c̄0 { Mγ̄
n +
√ γ̄
nγ +
√ H2L2
nγ (log γ + log n)
} ,
for some constant c̄0 > 0.
Note that the number of elements in I(m) is upper bounded by (m + 1)2. Since m is proportional to n, by setting γ̄ = c̄∗ log n for some constant c̄∗ > 0, we have 1 − 2(m + 1)2 exp(−γ̄) → 1. Consequently, it follows from Bonferroni’s inequality that the following event occurs w.p.a.1 for any I ∈ I(m) such that |I| ≥ cγ,∣∣∣∣∣ ∑ i∈L` I(Ai ∈ I){Yi −MLP(Xi; θI)} 2 |L`|Pr(A ∈ I) − EI(A ∈ I){Y −MLP(X; θI)} 2 Pr(A ∈ I)
∣∣∣∣∣ ≤ O(1) { M log n
n +
√ log n
nγ +
√ H2L2
nγ log n
} ,
(15)
whereO(1) denotes some positive constant. Under the given conditions, the RHS is o(1). The proof is hence completed.
D.2 PROOF OF LEMMA 2
Consider a given interval I ∈ D̂(`). We can find some interval I ′ ∈ I(m) ∩ D̂(`) that is adjacent to I. Consequently, the interval I ∪ I ′ belongs to I(m) as well. It follows that
1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂(`)I (Xi)} 2 + 1 |L`| ∑ i∈L` I(Ai ∈ I ′){Yi − q̂(`)I′ (Xi)} 2
≤ 1 |L`| ∑ i∈L` I(Ai ∈ I ∪ I ′){Yi − q̂(`)I∪I′(Xi)} 2 − γ.
(16)
By definition, q̂(`)I∪I′ minimizes the loss ∑ i∈Li I(Ai ∈ I ∪ I
′){Yi − q̂(`)I∪I′(Xi)}2. It follows that∑ i∈L` I(Ai ∈ I ∪ I ′){Yi − q̂(`)I∪I′(Xi)} 2 ≤ ∑ i∈L` I(Ai ∈ I ∪ I ′){Yi − q̂(`)I′ (Xi)} 2.
Combining this together with equation 16 yields that
γ ≤ 1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂(`)I′ (Xi)} 2. (17)
As both Y and q̂(`)I are uniformly bounded, we have γ ≤ c|L`|−1 ∑ i∈L` I(Ai ∈ I) for some constant c > 0. By Bernstein’s inequality, we have for any t > 0 such that
1 |L`| ∑ i∈L` I(Ai ∈ I) ≤ Pr(A ∈ I) + t, (18)
with probability 1− exp(−nt2/c̄(|I|+ t)) for some constant c̄ > 0. Similar to equation 15, by setting tI = c̄0 max( √ |I|n−1 log n, n−1 log n) for some constant c̄0 > 0, we obtain w.p.a.1 that
1 |L`| ∑ i∈L` I(Ai ∈ I) ≤ Pr(A ∈ I) + tI ,
for any I ∈ I(m). Under the condition that b is continuous and A × X is compact, b is bounded. Consequently, the probability Pr(A ∈ I) is proportional to the length of I. In view of equation 18, for any I ∈ D̂(`), it shall satisfy
γ c̄1 ≤ |I|+ log n n + 2
√ |I| log n
n ≤ 2|I|+ 2log n n ,
for some constant c̄1 > 0, w.p.a.1, where the last inequality follows from Cauchy-Schwarz inequality. As γ n−1 log n, we obtain |I| ≥ γ/c̄2 for any I ∈ D̂(`) and some constant c̄2 > 0, w.p.a.1. The proof is hence completed.
D.3 PROOF OF THEOREM 1
We begin with an outline of the proof. The proof is divide into four steps. In the first step, we show
Pr(|D̂(`)| ≤ |D`|)→ 1. (19) In the second step, we show
max τ∈J(D0) min τ̂∈J(D̂(`))
|τ̂ − τ | < δmin, (20)
where δmin = minI∈D0 |I|/3. By the definition of δmin, this implies that
Pr(|D̂(`)| ≥ |D`|)→ 1. (21) Combining equation 22 together with equation 19 proves (i) in Theorem 1. This proves (ii) in Theorem 1. In the third step, we show
max τ∈J(D0) min τ̂∈J(D̂(`))
|τ̂ − τ | = op(1). (22)
In the last step, we show (iii) holds. The proof is thus completed.
We next detail the proof for each of the step.
Step 1. Assume |D0| > 1. Otherwise, equation 19 automatically holds. Consider the partition D = {[0, 1]} which consists of a single interval and a zero q-function q[0,1](x) = 0 for any x. By definition, we have∑
I∈D̂(`) ∑ i∈L` I(Ai ∈ I){Yi − q̂I(Xi)}2 + |L`|γ|D̂(`)| ≤ ∑ i∈L` Y 2i + |L`|γ.
Since Y is uniformly bounded and γ = o(1), the right-hand-side (RHS) is O(n). Consequently, we obtain
|D̂(`)| ≤ c0γ−1, (23)
for some constant c0 > 0.
Notice that∑ I∈D̂(`) ∑ i∈L` I(Ai ∈ I){Yi − q̂(`)I (Xi)} 2 ≥ ∑ I∈D̂(`) ∑ i∈L`
I(Ai ∈ I){Yi − qI(Xi)}2︸ ︷︷ ︸ η1
+ ∑ I∈D̂(`) ∑ i∈L` I(Ai ∈ I){q̂(`)I (Xi)− qI(Xi)} 2
︸ ︷︷ ︸ η2
− 2 ∑ I∈D̂(`) ∣∣∣∣∣∑ i∈L` I(Ai ∈ I){Yi − qI(Xi)}{q̂(`)I (Xi)− qI(Xi)} ∣∣∣∣∣︸ ︷︷ ︸ η3 .
(24)
We next show η2, η3 = op(1).
Consider η2 first. Under Lemma 2, we have w.p.a.1 that
η2 = ∑
I∈D̂(`),|I|≥Cγ ∑ i∈L` I(Ai ∈ I){q̂(`)I (Xi)− qI(Xi)} 2.
We decompose the RHS by∑ I∈D̂(`),|I|≥Cγ |L`|E[I(A ∈ I){q̂(`)I (X)− qI(X)} 2|{Oi}i∈Li ]
+ ∑
I∈D̂(`),|I|≥Cγ ∑ i∈L` [I(Ai ∈ I){q̂(`)I (Xi)− qI(Xi)} 2 − E{I(A ∈ I){q̂(`)I (X)− qI(X)} 2|{Oi}i∈Li}].
The first line is op(n), by Lemma 1. Using similar arguments in equation 15, we can show the second term is upper bounded by
c̄|L`| ∑
I∈D̂(`),|I|≥Cγ
{ M log n
n +
√ |I|H2L2 log n
n
} ,
for some constant c̄ > 0, w.p.a.1. In view of equation 23, the above expression can be further bounded by
c̄n c0M log nγn + ∑ I∈D̂(`) √ |I|H2L2 log n n  ≤ O(1)n { log n γn + √ H2L2 log n γn } ,
where O(1) denotes some positive constant and the second inequality follows from the CauchySchwarz inequality. This yields that η2 = op(n).
Using similar arguments, we can show that η3 = op(n). It follows that∑ I∈D̂(`) ∑ i∈L` I(Ai ∈ I){Yi − q̂(`)I (Xi)} 2 ≥ η1 + op(n). (25)
Notice that η1 = ∑ I∈D̂(`) ∑ i∈L` I(Ai ∈ I){Yi −Q(Ai, Xi) +Q(Ai, Xi)− qI(Xi)}2
= ∑ I∈D̂(`) ∑ i∈L`
I(Ai ∈ I){Yi −Q(Ai, Xi)}2︸ ︷︷ ︸ η4
+ ∑ I∈D̂(`) ∑ i∈L`
I(Ai ∈ I){Q(Ai, Xi)− qI(Xi)}2︸ ︷︷ ︸ η5
+ 2 ∑ I∈D̂(`) ∑ i∈L`
I(Ai ∈ I){Yi −Q(Ai, Xi)}{Q(Ai, Xi)− qI(Xi)}︸ ︷︷ ︸ η6 .
Using similar arguments in bounding η2 and η3, we can show η6 = op(n) and that η5 = |L`| ∑ I∈D̂(`) ∫ I Eb(a|X)|Q(X, a)− qI(X)|2da+ op(n).
The first term on the RHS is greater than c0n ∑ I∈D̂(`) ∫ I E|Q(X, a)−qI(X)|
2da for some constant c0 > 0. It follows that
η1 ≥ η4 + c0|L`| ∑ I∈D̂(`) ∫ I E|Q(X, a)− qI(X)|2da+ op(n).
Note that under the step model assumption, η4 can be rewritten as ∑ I∈D0 ∑ i∈L` I(Ai ∈ I){Yi − qI(Xi)}2. This together with equation 25, yields that∑ I∈D̂(`) ∑ i∈L` I(Ai ∈ I)(Yi − q̂(`)I (Xi)) 2
≥ ∑ I∈D0 ∑ i∈L` I(Ai ∈ I)(Yi − qI(Xi))2 + c0|L`| ∑ I∈D̂(`) ∫ I E|Q(X, a)− qI(X)|2da+ op(n),
For any integer k such that 1 ≤ k ≤ |D0|−1, let τ∗0,k be the change point location that satisfies τ∗0,k = i/m for some integer i and that |τ0,k − τ∗0,k| < m−1. Denoted by D∗ the oracle partition formed by the change point locations {τ∗0,k} |D0|−1 k=1 . Set τ ∗ 0,0 = 0, τ ∗ 0,|D0| = 1 and q ∗ [τ∗0,k−1,τ ∗ 0,k) = q[τ0,k−1,τ0,k) for 1 ≤ k ≤ |D0| − 1 and q∗[τ∗0,K−1,1] = q[τ0,K−1,1]. Let ∆k = [τ ∗ 0,k−1, τ ∗ 0,k) ∩ [τ0,k−1, τ0,k)c for 1 ≤ k ≤ |D0| − 1 and ∆|D0| = [τ∗0,|D0|−1, 1] ∩ [τ0,|D0|−1, 1] c. The length of each interval ∆k is at most m−1. It follows that∑ I∈D∗ ∑ i∈L` I(Ai ∈ I){Yi − q∗I(Xi)}2 − ∑ I∈P0 ∑ i∈L` I(Ai ∈ I){Yi − qI(Xi)}2
≤ |D0|∑ k=1 ∑ i∈Li I(Ai ∈ ∆k) { Y 2i + sup I⊆[0,1] sup x q2I(x) } .
(26)
Using similar arguments in bounding the RHS of equation 17 in the proof of Lemma 2, we can show the RHS of equation 26 is op(n), asm diverges to infinity. Combining this together with equation 25 yields ∑
I∈D̂(`) ∑ i∈L` I(Ai ∈ I){Yi − q̂(`)I (Xi)} 2 ≥ ∑ I∈D∗ ∑ i∈L` I(Ai ∈ I){Yi − q∗I(Xi)}2
+c0|L`| ∑ I∈D̂(`) ∫ I E|Q(X, a)− qI(X)|2da+ op(n).
(27)
By definition, we have ∑ I∈D̂(`) ∑ i∈L` I(Ai ∈ I)(Yi − q̂(`)I (Xi)) 2 + |L`|γ|D̂(`)|
≤ ∑ I∈P∗ ∑ i∈L` I(Ai ∈ I){Yi − q∗I(Xi)}2 + |L`|γ|D∗|.
Since |D∗| = |D0|, it follows from equation 27 that∑ I∈D̂(`) ∫ I E|Q(X, a)− qI(X)|2da ≤ γ{|D0| − |D̂(`)|}+ op(1). (28)
Note that the left-hand-side (LHR) is non-negative. As γ > 0, we obtain that |D̂(`)| ≤ |D0|, w.p.a.1. This completes the proof for the first step.
Step 2. It follows from equation 28 that∑ I∈D̂(`) ∫ I E|Q(X, a)− qI(X)|2da ≤ γ{|D0| − 1}+ op(1),
and hence ∑ I∈D̂(`) ∫ I E|Q(X, a)− qI(X)|2da ≤ 2γ{|D0| − 1}, (29)
w.p.a.1. We aim to show equation 20 holds under the event defined in equation 29. Otherwise, there exists some τ0 ∈ J(P0) such that |τ̂ − τ0| ≥ δmin, for all τ̂ ∈ J(P̂ ). Under the event defined in equation 29, we obtain that∫ τ0+δmin
τ0−δmin E|Q(X, a)− qI(X)|2da ≤ 2γ(|D0| − 1), (30)
w.p.a.1. On the other hand, sinceQ(X, a) is a constant function on [τ0−δmin, τ0) or [τ0, τ0 +δmin), we have ∫ τ0+δmin
τ0−δmin E|Q(X, a)− qI(X)|2da
≥ min q
( δminE|q[τ0−δmin,τ0)(X)− q(X)| 2 + δminE|q[τ0,τ0+δmin)(X)− q(X)| 2 2 ) ≥ δmin
2 E|q[τ0−δmin,τ0)(X)− q[τ0,τ0+δmin)(X)| 2 ≥ δminκ0 2 ,
where
κ0 ≡ min I1,I2∈P0
I1 and I2 are adjacent
E|qI1(X)− qI2(X)|22 > 0.
This apparently violates equation 30 where γ ≤ δminκ0/(4|D0−4|). equation 20 thus holds w.p.a.1.
Step 3: Combining the results obtained in the first two steps yields that |D(`)| = |D0|, w.p.a.1. It follows from equation 28 that∑
I∈D̂(`)
∫ I E|Q(X, a)− qI(X)|2da = op(1).
Using similar arguments in Step 2, we can show maxτ∈J(D0) minτ̂∈J(D̂(`)) |τ̂ − τ | < w.p.a.1, for any sufficiently small > 0. The proof is hence completed.
Step 4: We begin with some notations. For any π, we define a random policy πD̂(`) according to the partition D̂(`) as follows:
πD̂(`)(a|x) = ∑ I⊆D̂(`) I{π(x) ∈ I, a ∈ I} b(a|x) b(I|x) .
Note that ∫ 1
0 πD̂(`)(a|x)da = ∑ I⊆D̂(`) I{π(x) ∈ I} = 1 for any x. Consequently, πD̂(`) is a “valid”
random policy.
The proposed value estimator V̂ (`) is doubly-robust to V (πD̂(`)). By Lemma 1, the estimated Qfunction is consistent. Consequently, V̂ (`) is consistent to V (πD̂(`)). Since the propose estimator is a weighted average of V̂ (`), it suffices to show V (πD̂(`)) is consistent to V (π). Note that
V (πD̂(`)) = E ∫
[0,1] Q(X, a) ∑ I⊆D̂(`) I{π(X) ∈ I, a ∈ I} b(a|X) b(I|X) da
= ∑ I0∈D0 ∫ I0 EqI0(X) ∑ I⊆D̂(`) I{π(X) ∈ I}b(I ∩ I0|X) b(I|X) .
Similarly, we can show V (π) = ∑ I0∈D0 ∫ I0 EqI0(X)I{π(X) ∈ I0}. (31)
For each I0 ∈ D0, there exists an interval I ∈ D̂(`) such that |I ∪ I0|/|I| → 1. We use Î(`)0 to denote this interval. Under the given conditions, we have
V (πD̂(`)) = ∑ I0∈D0 ∫ I0 EqI0(X)I{π(X) ∈ Î (`) 0 }+ op(1).
In view of equation 31, the value difference |V (πD̂(`))− V (π)| can be upper bounded by∑ I0∈D0 Pr(π(X) ∈ I0 − Î(`)0 ) + ∑ I0∈D0 Pr(π(X) ∈ Î(`)0 − I0) + op(1).
The first two terms decay to zero under the given conditions on π(X). The proof is hence completed.
D.4 PROOF OF THEOREM 2
We need the following lemma to prove Theorem 2.
Lemma 3 Assume conditions in Theorem 2 hold. Then for any interval I ∈ I(m) with |I| γ and any interval I ′ ∈ D̂(`) with I ⊆ I ′, we have w.p.a.1 that
E|qI(X)− qI′(X)|2 = o(1), where the little-o term is uniform in I and I ′.
The rest of the proof is divided into three steps. In the first step, we show Assertion (i) in Theorem 2 holds. In the second step, we show Assertion (ii) in Theorem 2 holds. Finally, we present the proof for Lemma 3.
Step 1. Consider a sequence {dn}n such that dn → 0 and dn γ. We aim to show inf a∈I′ I′∈D̂(`) E[|Q(X, a) − q̂I′(X)|2|{Oi}i∈L` ] = op(1). By Lemma 1, it suffices to show inf a∈I′ I′∈D̂(`) E|Q(X, a)− qI′(X)|2 = o(1).
Suppose |I ′| ≥ dn. Then according to Lemma 3, we can find some I such that a ∈ I ⊆ I ′, |I| → 0, E|qI(X) − qI′(X)|2 = o(1). Since |I| → 0 and a ∈ I, it follows from the uniform continuity of the Q-function that |qI(X)−Q(X, a)| → 0. The assertion thus follows. Next, suppose |I ′| < dn. Then |I ′| → 0 as well. It follows from the uniform continuity of the Q-function that infa∈I′ |qI′(X) − Q(X, a)| → 0. The assertion thus follows. This completes the proof for the first step.
Step 2. Using similar arguments in Step 4 of the proof of Theorem 1, it suffices to show V (πD̂(`)) = V (π) + op(1). By definition,
V (πD̂(`))− V (π) = ∑ I∈D̂(`) ∫ I EQ(X, a)I(π(X) ∈ I) b(a|X) b(I|X) da− EQ(π(X), X)
= ∑ I∈D̂(`) ∫ I E{Q(X, a)−Q(π(X), X)}I(π(X) ∈ I) b(a|X) b(I|X) da.
It follows that |V (πD̂(`))− V (π)| ≤ ∑ I∈D̂(`) inf a′∈I ∫ I E|Q(a′, X)−Q(π(X), X)|I(π(X) ∈ I) b(a|X) b(I|X) da
= ∑ I∈D̂(`) inf a′∈I E|Q(a′, X)−Q(π(X), X)|I(π(X) ∈ I)
≤ inf a′,a′′∈I,I∈D̂(`)
E|Q(a′, X)−Q(a′′, X)|.
In Step 1 of the proof, we have shown that infa∈I,I∈D̂(`) E|Q(X, a) − qI(X)| 2 = o(1). It follows that infa,a′∈I,I∈D̂(`) E|Q(X, a) − Q(a ′, X)|2 = o(1) and hence infa′,a′′∈I,I∈D̂(`) E|Q(a
′, X) − Q(a′′, X)| = o(1), by Cauchy-Schwarz inequality. This completes the proof for the second step.
Step 3. For a given interval I ′ ∈ D̂(`), the set of intervals I considered in Lemma 3 can be classified into the following three categories.
Category 1: I = I ′. It is immediate to see that qI = qI′ and the assertion automatically holds. Category 2: There exists another interval I∗ ∈ I(m) that satisfies I ′ = I∗ ∪ I. Notice that the partition D̂(`)∗ = D̂(`) ∪ {I∗} ∪ I − {I ′} corresponds to another partition. By definition, we have
1 |L`| ∑ i∈L` ∑ I0∈D̂(`)∗ I(Ai ∈ I0){Yi − q̂I0(Xi)}2 + γ|D̂(`)∗|
≥ 1 |L`| ∑ i∈L` ∑ I0∈D̂(`) I(Ai ∈ I0){Yi − q̂I0(Xi)}2 + γ|D̂(`)|,
and hence 1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂I(Xi)}2 + 1 |L`| ∑ i∈L` I(Ai ∈ I∗){Yi − q̂I∗(Xi)}2
≥ 1 |L`| ∑ i∈L` I(Ai ∈ I ′){Yi − q̂I′(Xi)}2 − γ.
It follows from the definition of q̂I∗ that 1 |L`| ∑ i∈L` I(Ai ∈ I∗){Yi − q̂I∗(Xi)}2 ≤ 1 |L`| ∑ i∈L` I(Ai ∈ I∗){Yi − q̂I′(Xi)}2.
Therefore, we obtain 1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂I(Xi)}2 ≥ 1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂I′(Xi)}2 − γ. (32)
Category 3: There exist two intervals I∗, I∗∗ ∈ I(m) that satisfy I ′ = I∗ ∪I ∪I∗∗. Using similar arguments in proving equation 32, we can show that
1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂I(Xi)}2 ≥ 1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂I′(Xi)}2 − 2γ.
Hence, regardless of whether I belongs to Category 2, or it belongs to Category 3, we have 1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂I(Xi)}2 ≥ 1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂I′(Xi)}2 − 2γ. (33)
Using similar arguments in equation 15, we can show w.p.a.1 that 1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂I(Xi)}2 = E[I(A ∈ I){Y − q̂I(X)}2|{Oi}i∈L` ] + o( √ γ|I|),
1 |L`| ∑ i∈L` I(Ai ∈ I){Yi − q̂I′(Xi)}2 = E[I(A ∈ I){Y − q̂I′(X)}2|{Oi}i∈L` ] + o( √ γ|I|),
where the little-o terms are uniform in I and I ′. Combining these together with equation 33 yields E[I(A ∈ I){Y − q̂I(X)}2|{Oi}i∈L` ] ≥ E[I(A ∈ I){Y − q̂I′(X)}2|{Oi}i∈L` ]− 2γ + o( √ γ|I|),
for any I and I ′, w.p.a.1. Note that qI satisfies E[I(A ∈ I){Y − qI(X)}|X] = 0. We have E[I(A ∈ I){qI(X)− q̂I(X)}2|{Oi}i∈L` ] ≥ E[I(A ∈ I){qI(X)− q̂I′(X)}2|{Oi}i∈L` ]− 2γ + o( √ γ|I|).
Consider the first term on the RHS. Note that
E[I(A ∈ I){qI(X)− q̂I′(X)}2|{Oi}i∈L` ] = E[I(A ∈ I){qI(X)− qI′(X)}2|{Oi}i∈L` ] +E[I(A ∈ I){q̂I′(X)− qI′(X)}2|{Oi}i∈L` ]− 2E[I(A ∈ I){qI(X)− qI′(X)}{q̂I′(X)− qI′(X)}|{Oi}i∈L` ]. By Cauchy-Schwarz inequality, the last term on the RHS can be lower bounded by
−1 2 E[I(A ∈ I){qI(X)− qI′(X)}2|{Oi}i∈L` ]− 2E[I(A ∈ I){q̂I′(X)− qI′(X)}2|{Oi}i∈L` ].
It follows that
E[I(A ∈ I){qI(X)− q̂I′(X)}2|{Oi}i∈L` ] ≥ 1
2 E[I(A ∈ I){qI(X)− qI′(X)}2|{Oi}i∈L` ] −3E[I(A ∈ I){q̂I′(X)− qI′(X)}2|{Oi}i∈L` ], and hence 1 2 E[I(A ∈ I){qI(X)− qI′(X)}2|{Oi}i∈L` ]− 2γ + o( √ γ|I|) ≤ E[I(A ∈ I){qI(X)− q̂I(X)}2|{Oi}i∈L` ] +3E[I(A ∈ I){qI′(X)− q̂I′(X)}2|{Oi}i∈L` ]. By Lemma 1 and the positivity assumption, the RHS is op(|I|). Note that the little-o terms are uniform in I and I ′. As |I| γ, we obtain that
E[I(A ∈ I){qI(X)− qI′(X)}2|{Oi}i∈L` ] = op(|I|), uniformly for any I and I ′, or equivalently,
E [ b(I|X) |I| {qI(X)− qI′(X)}2|{Oi}i∈L` ] = op(1).
By the positivity assumption, we have
E[{qI(X)− qI′(X)}2|{Oi}i∈L` ] = op(1), uniformly for any I and I ′, or equivalently,
E{qI(X)− qI′(X)}2 = op(1). This yields that
E{qI(X)− qI′(X)}2 = o(1), uniformly for any I and I ′. The proof is hence completed.
D.5 PROOF OF THEOREM 3
In the first two steps, we compare the minimax convergence rate when Q ∈ Q1. In the next two steps, we compare the minimax convergence rate when Q ∈ Q2. Step 1: We provide a lower bound for the minimax convergence of kernel-based OPE whenQ ∈ Q1 in this step. Consider a piecewise constant Q-function
Q(a, x) = { 0, if a ≤ 1/2, 1, otherwise.
Apparently, we have Q ∈ Q1 when C2 ≥ 1 ≥ . Define a policy π such that the density function of π(X) equals {
4/3, if 1/4 ≤ π(x) ≤ 1/2, 2/3, else if 1/2 ≤ π(x) < 4/3, 0, otherwise.
We aim to show for such Q and π, the best possible convergence rate of kernel-based estimator is n−1/3.
We first consider its variance. Since the conditional variance of Y |A,X is uniformly bounded away from 0 and 1, similar to Theorem 1 of Colangelo & Lee (2020), we can show the variance is lower bounded by O(1)(nh)−1 where O(1) denotes some positive constant.
We next consider its bias. Since the behavior policy is know. The bias is equal to
E K{(A− π(X))/h} hb(A|X) {Y −Q(π(X), X)} = EK{(A− π(X))/h} hb(A|X) {Q(A,X)−Q(π(X), X)}
= E ∫ π(X)+h/2 π(X)−h/2 K { a− π(X) h } {I(π(X) ≤ 1/2 < a)− I(a ≤ 1/2 < π(X))}da.
Using the change of variable a = ht+ π(X), the bias equals E ∫ 1/2 −1/2 K(t){I(π(X) ≤ 1/2 < π(X) + ht)− I(π(X) + ht ≤ 1/2 < π(X))}dt.
Consider any 0 < h ≤ for some sufficiently small > 0. The bias is then equal to
4
3 ∫ 1/2 1/2− /2 ∫ 1/2 −1/2 K(t){I(a ≤ 1/2 < a+ ht)− I(a+ ht ≤ 1/2 < a)}dtda
+ 2
3 ∫ 1/2+ /2 1/2 ∫ 1/2 −1/2 K(t){I(a ≤ 1/2 < a+ ht)− I(a+ ht ≤ 1/2 < a)}dtda.
Under the symmetric condition on the kernel function, the above quantity is equal to
2
3 ∫ 1/2 1/2−h/2 ∫ 1/2 (1−2a)/2h K(t)dtda ≥ 2 3 ∫ 1/2−h/4 1/2−h/2 ∫ 1/2 (1−2a)/2h K(t)dtda
≥ 2 3 ∫ 1/2−h/4 1/2−h/2 ∫ 1/2 1/4 K(t)dtda = h 6 ∫ 1/2 1/4 K(t)dt.
Consequently, the bias is lower bounded by O(1)h where O(1) denotes some positive constant.
To summarize, the root mean squared error of kernel based estimator is lower bounded by O(1){(nh)−1/2 + h} where O(1) denotes some positive constant. The optimal choice of h that minimizes such lower bound would be of the order O(n−1/3). Consequently, the convergence rate is lower bounded by O(1)n−1/3.
Step 2: We derive an upper bound for the minimax convergence rate of our estimator when Q ∈ Q1 in this step. Using similar arguments in the proof of Lemma 1 (see Section D.1) and the proof of Theorem 1 in Imaizumi & Fukumizu (2019), we can show that with proper choices of the MLP networks, the following holds with probability at least 1 − O(n−C) for some sufficiently large constant C > 0,
max I∈I(m),|I|≥cγ
E[|qI(X)− q̂(`)I (X)| 2|{Oi}i∈L` ] ≤ C(n|I|)−2β/(2β+d) log 2 n, (34)
where β is defined in Q1 and Q2. By equation 34 and the condition that γ n−2β/(2β+p) log2 n, using similar arguments in Steps 1 and 2 of the proof of Theorem 1, we can show that the estimated change point locations are consistent and D(`) = D0 with probability 1 − O(n−C). Similar to equation 27, it follows from equation 34 that ∑
I∈D̂(`)
∫ I E|Q(X, a)− qI(X)|2da = O(n−2β/(2β+p)),
up to some logarithmic factors. This further implies that the change point locations in D̂(`) converge at a rate of O(n−2β/(2β+p)) up to some logarithmic factors.
We next establish the statistical properties of our value estimates. We first observe that when b is known, our estimator is unbiased to the inverse-propensity score weighted estimator
1
n L∑ `=1 ∑ I∈D̂(`) ∑ i∈L` I(Ai ∈ I) I(π(Xi) ∈ I) b(I|Xi) Yi. (35)
Consequently, its bias is equal to
1
L L∑ `=1 ∑ I∈D̂(`) EI(A ∈ I) I(π(X) ∈ I) b(I|X) {Q(A,X)−Q(π(X), X)}
= 1
L L∑ `=1 ∑ I∈D̂(`) E ∫ a∈I I(π(X) ∈ I) b(I|X) {Q(a,X)− qI(X)}b(a|X)da.
Since b is bounded away from zero and infinity, and that supI Pr(π(X) ∈ I)/|I| ≤ C4, the absolute value of the bias is upper bounded by
O(1) 1
L L∑ `=1 ∑ I∈D̂(`) E ∫ a∈I I(π(X) ∈ I) |I| |Q(a,X)− qI(X)|da
≤ O(1) 1 L L∑ `=1 ∑ I∈D̂(`) ∫ a∈I sup x |Q(a, x)− qI(x)|da,
where O(1) denotes some positive constant. For each I ∈ D̂(`), we could find some I0 ∈ D0 such that the lengths of I −I0 and I0−I are upper bounded by O(n−2β/(2β+p)) up to some logarithmic factors. By definition, we have
qI(x) =
∫ a∈I ∑ I0∈D0 b(a|x)I(a ∈ I0)qI0(x)da
b(I|x) .
Consequently, supx |qI(x)− qI0(x)| = O(n−2β/(2β+p)) up to some logarithmic factors. It follows that the absolute value of the bias of bias is upper bounded by O(|D(`)|n−2β/(2β+p)) up to some logarithmic factors. As |D(`)| = |D0| ≤ C1 with probability tending to 1, the absolute value of the bias is upper bounded by
Op(n −2β/(2β+p)), (36)
up to some logarithmic factors.
We next consider the variance. The variance is of the same order of magnitude of that of equation 35. Since L is finite, it suffices to consider the variance of
L n ∑ I∈D̂(`) ∑ i∈L` I(Ai ∈ I) I(π(Xi) ∈ I) b(I|Xi) Yi. (37)
The variance of equation 37 is upper bounded by
O(n−1)E  ∑ I∈D̂(`) I(A ∈ I) I(π(X) ∈ I) b(I|X)  2 = O(n−1)  ∑ I∈D̂(`) EI(A ∈ I) I(π(X) ∈ I) b2(I|X)  = O(n−1)
 ∑ I∈D̂(`) E I(π(X) ∈ I) b(I|X)  = O(n−1|D̂(`)|), where the last equality is due to that b is uniformly bounded way from zero and that Pr(π(X) ∈ I) ≤ O(1)|I| for some positive constant O(1). As |D(`)| = |D0| ≤ C1 with probability tending to 1, the variance of equation 37 is upper bounded by Op(n−1). Consequently, the variance of our value estimator is bounded by Op(n−1) as well. In addition, we note this bound and the bias bound
in equation 36 is uniform in Q ∈ Q1 ∪ Q2 and π ∈ Π. The proposed value converges at a rate of Op(n
−2β/(2β+p)) up to some logarithmic terms. When 2β/(2β + p) > 1/2, both the bias and the variance decay at a rate of n−1/2 up to some logarithmic factors. This completes the proof for the first part.
Step 3. We provide a lower bound for the minimax convergence of kernel-based OPE whenQ ∈ Q2 in this step. Similar to Step 1, we can show its variance is lower bounded by O(n−1h−1). When h n−κ for some κ > 3/5, the root mean squared error is larger than or equal to O(n−(1−κ)/2). Consider the Q-function
Q(x, a) = Ch−1K
( a− π(x)
h
) ,
for some constant C > 0. With proper choice of C, we can show that such a choice of Q-function belongs to Q2. Using similar arguments in Step 1, we can show the bias equals
EC−1 K2{(A− π(X))/h} h2b(A|X) ≥ C−1EK 2{(A− π(X))/h} h2 .
Using similar arguments in Step 1, we can show the right-hand-side is lower bounded O(1)h. When h n−κ for some κ < 1/5, the root mean squared error is larger than or equal to O(n−κ). To summarize, when h n−κ for some κ < 1/5 or κ > 3/5, kernel-based OPE converges at a rate of n−κ ∗ for some κ∗ < 1/5.
Step 4. We provide an upper bound for the minimax convergence of our estimator when Q ∈ Q2 in this step.
Consider its bias first. Using similar arguments in Step 2, the bias is upper bounded by O(L−1) ∑L `=1 ∑ I′∈D̂(`) ∫ a∈I′ supx |Q(a, x) − qI′(x)|da. Similar to Lemma 3, we can show for any interval I ∈ I(m) with |I| γ and any interval I ′ ∈ D̂(`) with I ⊆ I ′, we have w.p.a.1 that
E|qI(X)− qI′(X)|2 = O(1) √ γ
|I| , (38)
provided that γ is proportional to n−2β/(2β+p) up to some logarithmic factors. For any a, there exists an interval I whose length is proportional to γ log n that covers a. Since Q ∈ Q2, we have supx |Q(a, x)− qI(x)| = O(|I|). This together with equation 38 yields that∑
I′∈D̂(`)
∫ a∈I′ sup x |Q(a, x)− qI′(x)|da ≤ O(1) ( |I|+ √ γ |I| ) .
Set I = γ1/3, the bias is upper bounded by O(1)γ1/3. Using similar arguments in Step 2 of the proof, the standard deviation is upper bounded by ∑L `=1 L−1 √ n−1|D̂(`)|. By equation 23, It is upper bounded by O(n−1/2γ−1/2). By setting γ to be proportional to n−3/5 (this rate is achievable under the condition that 4β > 3p), we obtain the rate of n−1/5.

Review summaries:
The main contribution of this paper is a new algorithm to learn the expected reward function for a given target policy using the historical data generated by a different behavior policy in continuous action domains.  All current Offline-Policy Evaluation (OPE) methods for handling continuous action domains use a kernel function to extend Inverse Probability Weighting (IPW) or Doubly Robust (DR) approaches for discrete action domains.  The algorithm proposed in this work adaptively discretizes the action space by combining methods in multi-scale changepoint detection, multi-layer perceptron regression and OPE in discrete action domains.  The finite sample performance of the proposed method, known as Deep-Jump Q-Evaluation (DJQE), is compared to that of two kernel-based methods, one due to Kallus and Zhou (2018) and another due to Colangelo and Lee (2020), on synthetic as well as real-world data.  To generate synthetic data, four scenarios are considered, where in each case the Q-function is continuous in the action domain or is a piecewise function of the action.  In almost all of these cases, DJQE outperforms the two kernel-based methods.  Similarly, when applied to real-world Warfarin data (after calibration), DJQE outperforms the two kernel-based methods with respect to the bias, standard deviation and mean squared error, even when the sample size is small (n=50).  The average runtime of DJQE in each scenario (for synthetic or real-world data) is about 5 minutes. 

This paper proposes a new method for offline evaluation when the action space is continuous, one dimensional. This overcomes the drawbacks of the kernel based method, which cannot be applied to non-smooth Q functions and requires heavy computation to optimize the bandwidth. The proposed method can be applied to discontinuous Q functions like step functions, and achieves smaller bias. This is made possible by the adaptive jump q learning method.

This paper considers the problem of off-policy evaluation with continuous actions. The main idea is to first using multi-scale change point detection to discretize the action space and then apply traditional IPW or DR methods to estimate the value. The DJQE method is theoretically analyzed under both the cases that the Q function is either a piecewise function or a continuous function. For continuous function, it is not surprising that as the number of splits m goes to infinity as n, the estimation is consistent, while additional results in Theorem 2 also shows that for limited m, the estimator can also be shown as a uniform approximation of the Q value. Experiments consider both a toy dataset and a real problem in personalized does finding, and the results show that the DJQE method is superior than existing methods for continuous Q evaluation.