Deconfounding Reinforcement Learning in Observational Settings
Chaochao Lu Bernhard Schölkopf José Miguel Hernández-Lobato








Deconfounding Reinforcement Learning 
in Observational Settings
    Chaochao Lu cl641@cam.ac.uk 

       Department of Engineering

       University of Cambridge 

       Cambridge, United Kingdom 

 
       Department of Empirical Inference 

       Max Planck Institute for Intelligent Systems 

       Tübingen, Germany
       Bernhard Schölkopf bs@tuebingen.mpg.de

       Department of Empirical Inference 

       Max Planck Institute for Intelligent Systems 

       Tübingen, Germany
       José Miguel Hernández-Lobato jmh233@cam.ac.uk 

       Department of Engineering

       University of Cambridge 

       Cambridge, United Kingdom 

 
       Microsoft Research 

       Cambridge, United Kingdom 

 
       Alan Turing Institute 

       London, United Kingdom 
       
    December 30, 2023
==================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================




We propose a general formulation for addressing reinforcement learning (RL) problems in settings with observational data. That is, we consider the problem of learning good policies solely from historical data in which unobserved factors (confounders) affect both observed actions and rewards. Our formulation allows us to extend a representative RL algorithm, the Actor-Critic method, to its deconfounding variant, with the methodology for this extension being easily applied to other RL algorithms. 
In addition to this, we develop a new benchmark for evaluating deconfounding RL algorithms by modifying the OpenAI Gym environments and the MNIST dataset. Using this benchmark, we demonstrate that the proposed algorithms are superior to traditional RL methods in confounded environments with observational data. To the best of our knowledge, this is the first time that confounders are taken into consideration for addressing full RL problems with observational data. Code is available at https://github.com/CausalRL/DRL.




§ INTRODUCTION

In recent years, much attention has been devoted to the development of reinforcement learning (RL) algorithms with the goal of improving treatment policies in healthcare <cit.>. RL algorithms have been proposed to infer better decision-making strategies for mechanical ventilation <cit.>, sepsis management <cit.>, and treatment of schizophrenia <cit.>. In healthcare, a common practice is to focus on the observational setting in which we learn policies solely from historical data produced by real environments, instead of learning policies by actively taking actions in the traditional RL setting. The reason for this is that we do not wish to experiment with patients' lives without evidence that the proposed treatment strategy is better than the current practice <cit.>.[Another example is in finance where, due to costs in terms of time and money, it is often impractical to evaluate trading strategies by actually buying and selling stocks in the market.] As pointed out by <cit.>, RL also has advantages over other machine learning algorithms even in the observational setting, especially in two situations: when the optimal treatment strategy is unclear in medical literature <cit.>, and when training examples do not represent optimal behavior. 

Another approach that has been explored and used is causal inference <cit.>. Causal inference is a fundamental notion in science and plays an increasingly important role in healthcare and medicine <cit.>. Especially in the observational setting, causal inference offers a powerful tool to deal with confounding that occurs when a hidden variable influences both the treatment and the outcome <cit.>. It is widely acknowledged that confounding is the most crucial aspect of inferring the effect of the treatment on the outcome from observational data <cit.>, because not taking confounding into account might lead to choosing a wrong treatment. However, most approaches in causal inference to tackling confounding are restricted to fixed treatments which do not vary over time <cit.>. When treatments are allowed to vary over time, the space of combinational treatment strategies increases exponentially <cit.>. In this case, it is still unclear how to apply causal inference to dealing with confounding in such sequential data where treatments changes over time. 

On the basis of the discussion above, in this paper we attempt to combine advantages of RL and causal inference to cope with an important family of RL problems in the observational setting, that is, learning good policies solely from the historical data produced by real environments with confounding bias. This type of problem will become increasingly common in future RL research with the burgeoning development in healthcare and medicine. To the best of our knowledge, however, little work has been done in this promising area of integrating RL with causal inference <cit.>. In contrast, confounders have been extensively studied in epidemiology, sociology, and economics. Take for example the widespread kidney stones in which the size of the kidney stone is a confounding factor affecting both the treatment and the recovery <cit.>. Correcting for the confounding effect of the size of the kidney stone is crucial for determining how to choose an effective treatment. Similarly, in RL, if unobserved potential confounders exist, they affect both actions and rewards as an agent interacts with environments, and eventually influence the policy to be optimized. 

Let us for a moment stick with the example of the kidney stones and assume that physicians need to take a series of steps to treat this disease. We further assume that, during the course of treatment, physicians have natural predilections of the treatment choices as a function of the stone size. That is, it is more likely for them to choose some treatment when patients have large stones and to choose another treatment when patients have small stones. In the observational setting, given only such a set of historical records about physicians' treatments and outcomes on patients, it is extremely challenging or even impossible to learn an optimal treatment policy due to the existence of the confounder (i.e., the size of kidney stones). 

To this end, we present a general formulation for addressing RL problems with observational data, namely deconfounding reinforcement learning (DRL). More specifically, given several common confounding assumptions, we first estimate a latent-variable model from observational data. Under some conditions for identification, through the latent-variable model we can simultaneously discover the latent confounders and infer how they affect actions and reward. Then the confounders in the model can be adjusted for, using the causal language developed by <cit.>, and finally we optimize the policy based on the resulting deconfounding model. 

On the basis of the proposed formulation, we extend one popular RL algorithm, the Actor-Critic method, to its corresponding deconfounding variant. Note that our procedure for obtaining this deconfounding variant can be easily applied to other RL algorithms. Due to lack of datasets in this respect, we revise the classic control toolkit in OpenAI Gym <cit.>, making it a benchmark for comparison of DRL algorithms. In addition, we also devise a confounding version of the MNIST dataset <cit.> to verify the performance of our causal model. Finally, we conduct extensive experiments to demonstrate the superiority of the proposed formulation compared to traditional RL algorithms in confounded environments with observational data.  

To sum up, our contributions in this paper are as follows:


  * We propose a general formulation for addressing RL problems in confounded environments with observational data, namely deconfounding reinforcement learning (DRL);

  * We present the deconfounding variant of Actor-Critic methods, obtained through a methodology that can be easily applied to other RL methods;

  * We develop a benchmark for DRL by revising the toolkit for classic control in OpenAI Gym <cit.> and by devising a version of the MNIST dataset with confounders <cit.>;

  * We perform a comprehensive comparison of our DRL algorithm with its vanilla version, showing that the proposed approach has an advantage in confounded environments.

  * To the best of our knowledge, this is the first attempt to build a bridge between confounding and the full RL problem. This is one of few research papers aiming at understanding the connections between causal inference and the full RL setting.




§ BACKGROUND


In this section, we briefly review confounding in causal inference. We recommend Pearl's excellent monograph for further reading <cit.>.



 §.§ Simpson's Paradox

Let us begin with an example of one of the most famous paradoxes in statistics: Simpson's Paradox. Consider the previously mentioned kidney stones  <cit.>. We collect electronic patient records to investigate the effectiveness of two treatments against kidney stones. Although the overall probability of recovery is higher for patients who took treatment b, treatment a performs better than treatment b both on patients with small kidney stones and with large kidney stones. More precisely, we have 

    p(R=1|T=b)    > p(R=1|T=a);     but
    
      p(R=1|T=b,Z=l)    < p(R=1|T=a,Z=l), 
    
    p(R=1|T=b,Z=s)    < p(R=1|T=a,Z=s);

where Z is the size (large or small) of the stone, T the treatment, and R the recovery (all binary). How do we cope with this change in conclusions? Which treatment do you prefer if you had kidney stones? Does treatment b cause recovery? As described in the following section, the answers to these questions depend on the causal relationship between treatment, recovery, and the size of the kidney stone.





 §.§ Confounding

An intuitive explanation for this kidney stone example of Simpson's paradox is that larger stones are more severe than small stones and are much more likely to be treated with treatment a, resulting in that treatment a looks worse than treatment b overall. We assume that Figure <ref> depicts the true underlying causal diagram of the kidney stone example. Confounding occurs because the size of kidney stones influences both treatment and recovery. The size of the kidney stones is called a confounder. The term “confounding” originally meant “mixing” in English, which describes that the true causal effect T → R is “mixed” with the spurious correlation between T and R induced by the fork T ← Z → R <cit.>. In other words, we will not be able to disentangle the true effect of T on R from the spurious effect if we do not have data on Z. If we have measurements of Z or can indirectly estimate Z, it is easy to deconfound the true and spurious effects. To this end, we can adjust for Z by averaging the effect of T on R in each subgroup of Z (i.e., different size groups in the case of kidney stones).



 §.§ Do-Operator and Adjustment Criterion
 
From the viewpoint of causal inference, we can use the language of intervention, namely the do-operator, to describe when confounding happens. In fact, in the example of the kidney stones, what we are interested in is how these two treatments compare when we force all patients to take treatment a or b, rather than which treatment has a higher recovery rate given only the observational patient records. Mathematically, we focus on the true effect p(R=1|do(T=a)) (the intervention distribution obtained when patients are forced to take treatment a) instead of the conditional p(R=1|T=a) (the observational distribution obtained when patients are observed to take treatment a). Therefore, as described previously, confounding can be naturally described as the discrepancy between p(R|T) and p(R|do(T)).

The do-operator can be executed in two common ways: by Randomized Controlled Trials (RCTs) <cit.> and by adjustment formulas (i.e., Back-door criterion and Front-door criterion) <cit.>. RCTs are the gold standard but sometimes limited by practical considerations (e.g., safety, laws, ethics, physically infeasibility, etc.). The Back-door and Front-door criteria require knowledge of the causal diagram, which typically means that causal assumptions are thus provided in advance. According to the Back-door criterion, in the kidney stone example, we can immediately obtain 

    p(R=1|do(T=a)) = ∑_z∈{l, s} p(R=1|T=a, Z = z)p(Z=z).

Beyond the two adjustment formulas, the do-calculus provides a syntactic method of deriving claims about interventions <cit.>. It consists of three rules which can be repeatedly applied to simplify the expression for an interventional distribution. With the do-calculus we can calculate intervention probabilities from observation probabilities.



 §.§ Proxy Variables for Confounding


If confounders can be measured, then they can be adjusted for through the methods discussed in Section <ref>. However, in most cases where confounders are hidden or unmeasured, without further assumptions, it is impossible to estimate the effect of the intervention on the outcome. A common practice is then to leverage observed proxy variables that contain information about the unobserved confounders <cit.>. However, using proxy variables to correctly recover causal effects requires strict mathematical assumptions <cit.> and, in practice, we do not know whether or not the available data meets those assumptions. Hence, we decide to follow <cit.> and, instead of using proxy variables, we estimate a latent-variable model in which we simultaneously discover the latent confounders and infer how they affect treatment and outcome. 



§ DECONFOUNDING REINFORCEMENT LEARNING


In this section, we will formally introduce deconfounding reinforcement learning (DRL). Generally speaking, DRL consists of two steps: learning a deconfounding model shown in Figure <ref> and optimizing a policy based on the learned deconfounding model. The main idea in step 1 is to simultaneously discover hidden confounders and infer causal effects through the estimation of a latent-variable model. More specifically, we first discuss the time-independent confounding assumption in Section <ref> and then, based on this assumption, we formalize the deconfounding model in Section <ref>. Section <ref> talks about the problem of identification in our model, which is a central issue in causal inference. After that, we present details about how to learn the proposed model via variational inference in Section <ref>. Step 2 is provided in Section <ref> where we describe how to design the deconfounding actor-critic method. The proposed approach is straightforward to apply to other RL algorithms in the same manner. Finally, the training procedure 
for our deconfounding actor-critic method is presented in Section <ref>.



 §.§ Causal Assumptions
 
Without loss of generality, we assume that there exists a common confounder in the model, denoted by u in Figure <ref>, which is time-independent across a number of episodes. This assumption is sufficiently general to apply to various RL tasks across domains. For example, in personalized medicine, socio-economic status can be a confounder affecting both the medication strategy a patient has access to and the patient's general health <cit.>. In this case, socio-economic status is time-independent for each patient during the course of treatment. Another example is in agriculture, where soil fertility may serve as a time-independent confounder affecting both the application of fertilizer and the yield of each plot of land <cit.>. Finally, in the stock market example, government policy may also act as a time-independent confounder. 





 §.§ The Model
 
Given our causal assumptions, we first fit a generative model to a sequence of observational data: observations, actions, and rewards, where actions and rewards are confounded by one or several unknown factors, as shown in Figure <ref>. Formally, let x⃗=(x_1, …, x_T), a⃗=(a_1, …, a_T-1), r⃗=(r_2, …, r_T+1), z⃗=(z_1, …, z_T) be the sequences of observations, actions, rewards, and corresponding latent states, respectively. As mentioned previously, the confounder is denoted by u, and it is worth noting that here u may stand for more than one confounder in which multiple confounders are seen as a whole represented by u (i.e., the confounder u can be a vector). We assume that x_t ∈ℝ^D_x, a_t ∈ℝ^D_a, r_t ∈ℝ^D_r, z_t ∈ℝ^D_z, and u ∈ℝ^D_u, where D_z ≪ D_x. The generative model for DRL is then given by

    p(z_t) = ∏_j=1^D_z𝒩(z_tj|0, 1);      
    p(u) = ∏_j=1^D_u𝒩(u_j|0, 1); 
       p(x_t | z_t, u) = 𝒩(x_t | μ̂_̂t̂^̂x̂, σ̂_̂t̂^̂x̂^2);
      μ̂_̂t̂^̂x̂ = f_1(z_t,u),   σ̂_̂t̂^̂x̂^2 = f_2(z_t,u); 
       p(a_t | z_t, u) = 𝒩(a_t | μ̂_̂t̂^̂â, σ̂_̂t̂^̂â^2);
      μ̂_̂t̂^̂â = f_3(z_t, u),   σ̂_̂t̂^̂â^2 = f_4(z_t, u); 
       p(r_t+1 | z_t, a_t, u) = 𝒩(r_t+1 | μ̂_̂t̂^̂r̂, σ̂_̂t̂^̂r̂^2);
      μ̂_̂t̂^̂r̂ = f_5(z_t, a_t, u),   σ̂_̂t̂^̂r̂^2 = f_6(z_t, a_t, u); 
       p(z_t | z_t-1, a_t-1) = 𝒩(z_t | μ̂_̂t̂^̂ẑ, σ̂_̂t̂^̂ẑ^2);
      μ̂_̂t̂^̂ẑ = f_7(z_t-1, a_t-1),   σ̂_̂t̂^̂ẑ^2 = f_8(z_t-1, a_t-1);

where we have parametrized each probability distribution as a Gaussian with mean and variance given by nonlinear functions f_k represented by neural networks with parameters θ_k for k=1, …, 8. Note that, in our case, Equation (<ref>) is not necessary when learning the model, because the data used in our experiments are generated from a policy depending only on the confounder (see Section <ref>), that is, a_t does not depend on z_t in our data. In addition, in this case z_t is not viewed as a confounder of a_t and r_t+1, for z_t does not influence a_t. This also provides one reason why we do not need to adjust for z_t. However, in some other cases such as when applied to medical data, the strategies of treatment from physicians definitely contain valuable information about z_t and a_t (i.e., z_t does influence a_t), and therefore the arrow between them is necessary when learning the model.



 §.§ Identification of Causal Effect
 
The key component of our method that allows us to address problems with confounders is the computation of the reward according to the model from Figure <ref>. To be more precise, assuming that an agent standing at state z_t performs action a_t=𝔞, we do not use the conditional p(r_t | z_t, a_t=𝔞) as predictor for the reward as traditional RL methods do. Instead, we use the do-operator described in Section <ref> to compute the reward[To give some intuition, take for example the case of kidney stones. One prefers treatment b when considering the overall probability of recovery. By contrast, one chooses treatment a when considering the recovery rate when conditioning on each possible kidney stone size. The optimal treatment is treatment a in this example. The do-operator in Equation (<ref>) allows us to obtain the correct solution in  this type of problems with confounders.  By contrast, traditional RL methods will fail in such settings since they consider only conditional probabilities (e.g., the overall probability of recovery) instead of interventional probabilities, as the do-operator does.]:

    p(r_t | z_t, do(a_t=𝔞))
       = ∫_u p(r_t | z_t, do(a_t=𝔞), u)p(u|z_t, do(a_t=𝔞)) du 
       =∫_u p(r_t | z_t, a_t=𝔞, u) p(u) du,

where Equation (<ref>) is obtained by applying the rules of do-calculus to the causal graph in Figure <ref> <cit.>. We also can use the back-door criterion to directly obtain Equation (<ref>). Equation (<ref>) shows that p(r_t | z_t, do(a_t=𝔞)) can be identified from the joint distribution p(u, z⃗, x⃗, a⃗, r⃗), because, through Equation (<ref>), the interventional probability p(r_t | z_t, do(a_t=𝔞)) is converted to observation probabilities with regard to p(u, z⃗, x⃗, a⃗, r⃗). In other words, if we can recover p(u, z⃗, x⃗, a⃗, r⃗) then we also can recover p(r_t | z_t, do(a_t=𝔞)). 

Now the problem is reduced to whether or not we can estimate p(u, z⃗, x⃗, a⃗, r⃗) from observations of (x⃗, a⃗, r⃗). Fortunately, it is possible, because a number of works have shown that one can use the knowledge inferred from the joint distribution between proxy variables and confounders to adjust for the hidden confounders <cit.>. Here, we focus only on one possible case presented in Figure 1 of <cit.> (see Appendix <ref>), because their result can be directly used to show that the joint distribution p(u, z⃗, x⃗, a⃗, r⃗) can be approximately recovered solely from the observations (x⃗, a⃗, r⃗). The main idea is that our model can be factorized into two types of 4-tuple components, each of which can be proved using the same method in <cit.>. More precisely, in our model, at each time step the 4-tuple (z_t, x_t, a_t, r_t+1) is exactly the same as the 4-tuple (Z, X, t, y) shown in Figure 1 of <cit.>. Since it has been proved in <cit.> that the joint distribution p(Z, X, t, y) can be recovered from observations of (X, t, y), p(z_t, x_t, a_t, r_t+1) can be also approximately recovered only from the observations (x_t, a_t, r_t+1). Likewise, the joint distribution over the other 4-tuple (u, x_t, a_t, r_t+1) can be recovered from (x_t, a_t, r_t+1) in the same manner (see Appendix <ref>). Applying this rule repeatedly to the sequential data, we are finally able to approximately recover p(u, z⃗, x⃗, a⃗, r⃗) solely from the observations (x⃗, a⃗, r⃗). Hence, in the present paper we may reasonably assume that the joint distribution p(u, z⃗, x⃗, a⃗, r⃗) can be approximately recovered solely from the observations (x⃗, a⃗, r⃗). 

We estimate the joint distribution p(u, z⃗, x⃗, a⃗, r⃗) using Variational Auto-Encoders (VAEs) <cit.>, which can represent a very large class of distributions with latent variables and are easily implemented by solving an optimization problem <cit.>. However, VAEs do not guarantee that the true model can be identified through learning <cit.>. 
For example, it might occur that different model parameters could equally fit the same observed data as long as the values of the latent confounders are adjusted accordingly <cit.>. Despite this, identification is possible in both our model and the one described by <cit.> under specific conditions  <cit.>. In addition to this, in our model, each of the z_t and u have at least T proxy variables {x_t}_t=1,…,T, where z_t can be approximately viewed as hidden states in a hidden Markov model. In this case, <cit.> show that such a model can be identified if the latent confounders z_t and u are assumed to be categorical. However, in practice there is no guarantee for this to be the case. Hence, we prefer to use VAEs, which make substantially weaker assumptions about the data generating process and the structure of the hidden confounders <cit.>. Despite the lack of general identifiability results for VAEs, we empirically show that our approach is useful for learning a better policy in the presence of confounding.  

In the model from Figure <ref> there are two types of confounders: the time-independent confounder u and the time-dependent confounders {z_t}_t=1,…,T, playing different roles in the model. We refer to the former as “nuisance confounders” and the latter as “advantage confounders”. The nuisance confounder u will negatively affect the whole course of treatment and, therefore, should be adjusted for. In the example of kidney stones, the existence of the confounder (i.e., the size of stones) will lead to a wrong treatment if not adjusted for. By contrast, the advantage confounders {z_t} act as state variables and, in principle, they are not supposed to be adjusted for because in RL we aim to take optimal actions when conditioning on the current state[Another reason why we do not need to adjust for z_t has been discussed in Section <ref>.]. 



 §.§ Learning
 
The nonlinear functions parametrized by neural networks make inference intractable. Because of this, we learn the parameters of the model θ_k by using variational inference together with an inference model, a neural network which approximates the intractable posterior distributions <cit.>. We now review  how to learn a simple latent variable model p_θ(x|z) using variational inference.
In this simple model x stands for the observational variables and z for the latent variables.
By using the variational principle, we introduce an approximate posterior distribution q_ϕ(z|x) to obtain the following lower bound on the model's marginal likelihood:

    log p_θ(x) ≥q_ϕ(z|x)𝔼[log p_θ(x | z)] - KL(q_ϕ(z|x)||p_θ(z)),

where we have used Jensen's inequality, and KL(·||·) denotes the Kullback-Leibler divergence between two distributions
and ϕ are the parameters of the inference model q(z|x). 



  §.§.§ Variational Lower Bound

By directly applying the lower bound in Equation (<ref>) to our model, we obtain

    log p_θ(x⃗, a⃗, r⃗) 
    ≥   
    q_ϕ(z⃗, u|x⃗, a⃗, r⃗)𝔼   [log p_θ(x⃗, a⃗, r⃗ , z⃗, u)] - KL(q_ϕ(z⃗, u|x⃗, a⃗, r⃗)||p_θ(z⃗, u)) 
    
    ≡ℒ(x⃗, a⃗, r⃗;θ, ϕ).

By using the Markov property of our model, we can factorize the full joint distribution in the following way:

    p_θ(x⃗, a⃗, r⃗, z⃗, u) =     p(u)p(z_1)[∏_t=1^T p(x_t|z_t,u)p(a_t|z_t, u)p(r_t+1|z_t, a_t, u)]
       [∏_t=2^T p(z_t|z_t-1, a_t-1)].

In addition, for simplicity, we also have a similar factorization assumption in the posterior approximation for z⃗ and u:

    q_ϕ(z⃗, u|x⃗, a⃗, r⃗) = q(u|x⃗, a⃗, r⃗)q(z_1|x⃗, a⃗, r⃗)∏_t=2^T q(z_t|z_t-1, x⃗, a⃗, r⃗).

Combining equations (<ref>), (<ref>) and (<ref>) yields

    ℒ(x⃗, a⃗, r⃗;θ, ϕ)
    =   ∑_t=1^T z_t ∼ q(z_t|z_t-1, x⃗,a⃗,r⃗) 
     u ∼ q(u|x⃗,a⃗,r⃗)𝔼[log p(x_t|z_t,u)+log p(a_t|z_t, u)+log p(r_t+1|z_t, a_t, u)] 
       -KL(q(u|x⃗,a⃗,r⃗) || p(u)) 
    -KL(q(z_1|x⃗,a⃗,r⃗) || p(z_1))
       -∑_t=2^T z_t-1∼ q(z_t-1|z_t-2,x⃗,a⃗,r⃗)𝔼[KL(q(z_t|z_t-1, x⃗,a⃗,r⃗) || p(z_t|z_t-1,a_t-1))],

where we omit the subscripts θ and ϕ. A more detailed derivation can be found in Appendix <ref>. Equation (<ref>) is differentiable with respect to the model parameters θ and the inference parameters ϕ and, by using the reparametrization trick <cit.>, we can directly apply backpropagation to optimize this objective. 



  §.§.§ Inference Model
 
From the factorization in Equation (<ref>), we can see that there are two types of inference models: q(u|x⃗, a⃗, r⃗) and q(z_t|z_t-1, x⃗, a⃗, r⃗). Similar to the generative model in Section <ref>, we also parametrize both of them as Gaussian:

    q(u|x⃗, a⃗, r⃗) = 𝒩(u | μ̂_̂t̂^̂û, σ̂_̂t̂^̂û^2);
      μ̂_̂t̂^̂û = f_9(x⃗, a⃗, r⃗),   σ̂_̂t̂^̂û^2 = f_10(x⃗, a⃗, r⃗); 
       q(z_t|z_t-1, x⃗, a⃗, r⃗) = 𝒩(z_t | μ̂_̂t̂^̂𝔷̂, σ̂_̂t̂^̂𝔷̂^2);
      μ̂_̂t̂^̂𝔷̂ = f_11(z_t-1, x⃗, a⃗, r⃗),   σ̂_̂t̂^̂𝔷̂^2 = f_12(z_t-1, x⃗, a⃗, r⃗).

In fact, as shown in Equation (<ref>), q(z⃗|x⃗, a⃗, r⃗) can be factorized as the product of q(z_t|z_t-1, x⃗, a⃗, r⃗) for t=1, …, T. By using the Markov property of our model, we have

    z_t  x_1,…, x_t-1, a_1, …, a_t-2, r_2, …, r_t | z_t-1.

Therefore, we can further simplify each of these terms as follows:

    q(z_t|z_t-1, x⃗,a⃗,r⃗)=q(z_t|z_t-1,a_t-1_past, x_t,a_t,r_t+1_present,x_t+1,a_t+1,r_t+2,…,x_T,a_T,r_T+1_future).

Equation (<ref>) tells us that z_t depends on z_t-1 and all the current and future observed data in (x⃗,a⃗,r⃗). Meanwhile, the conditional independence above means that z_t-1 contains all the historical data. Therefore, it is natural to calculate z_t based on the whole sequence of data. This can be done using recurrent neural networks (RNNs). Inspired by <cit.>, we similarly choose a bi-directional LSTM <cit.> to parameterize f_11 and f_12 in Equation (<ref>). Since Equation (<ref>) has the same structure as Equation (<ref>), f_9 and f_10 are parameterized by a bi-directional LSTM as well. More details about the architecture can be found in Appendix <ref>.

Note that to generate data from the model, at any time step t, we require to know a_t and r_t+1 before inferring the distribution over z_t when conditioning on the observed x_t. Hence, we need to introduce two auxiliary distributions to perform counterfactual reasoning. This corresponds to the problem of predicting x_t+1 given an x_t unseen in the training set. To be more precise, we have

    q(a_t|x_t) = 𝒩(a_t|μ̂_̂t̂^̂𝔞̂, σ̂_̂t̂^̂𝔞̂^2);
      μ̂_̂t̂^̂𝔞̂ = f_13(x_t),   σ̂_̂t̂^̂𝔞̂^2 = f_14(x_t); 
       q(r_t+1|x_t, a_t) = 𝒩(r_t+1|μ̂_̂t̂^̂𝔯̂, σ̂_̂t̂^̂𝔯̂^2);
      μ̂_̂t̂^̂𝔯̂ = f_15(x_t, a_t),   σ̂_̂t̂^̂𝔯̂^2 = f_16(x_t, a_t),

where f_13, f_14, f_15, and f_16 are also parameterized by neural networks. To estimate the parameters of these inference models, we will add extra terms in the variational lower bound given by Equation (<ref>):

    ℒ_DRL(x⃗, a⃗, r⃗;θ, ϕ) = ℒ(x⃗, a⃗, r⃗;θ, ϕ) + ∑_t=1^T(log q(a_t|x_t) + log q(r_t+1|x_t, a_t)).






 §.§ Deconfounding RL Algorithms
 
We have now all the building blocks for our DRL algorithm. Once our model is learned from the observational data, it can be directly used as a dynamic environment like those in OpenAI Gym <cit.>. We can exploit the learned model to generate rollouts for policy evaluation and optimization. In practice, Equation (<ref>) is approximated using the Monte Carlo method as follows:

    p(r_t | z_t=𝔷, do(a_t=𝔞))=1/N∑_i=1^N p(r_t | z_t=𝔷, a_t=𝔞, u_i)   
    u_i ∼ p(u),

where N is the sample size from the prior p(u). Given observational data, we sample u instead from the approximate posterior q(u|x⃗, a⃗, r⃗) which we compute using the inference network presented in Section <ref>.

By using this deconfounding reward function, it is straightforward to extend traditional RL algorithms to their corresponding deconfounding version. In this paper, we select and extend one representative RL algorithm: the Actor-Critic method <cit.>. Nevertheless, our methodology can be used to extend other algorithms as well in a straightforward manner.




  
Deconfounding Actor-Critic Methods The actor-critic method is a policy-based RL method directly parameterizing a policy function π(a|z; θ). The goal is to reduce variance in the estimate of the policy gradient by subtracting a learned function of the state, known as a baseline, from the return. The learned value function V(z; ϕ) is commonly used as the baseline. Taking into consideration that the return is estimated by Q(z, a; ϕ_Q), we can write the gradient of the actor-critic loss function at step time t as

    ∇ J(θ) = 𝔼_π[(Q(z_t, a_t; ϕ_Q)-V(z_t; ϕ_V))∇_θlnπ(a_t|z_t; θ)],

where Q(z_t, a_t; ϕ_Q)-V(z_t; ϕ_V) is an estimate of the advantage of action a_t in state z_t. In practice, Q(z_t, a_t; ϕ_Q) is usually replaced with the one-step return, that is, r_t+1+V(z_t+1; ϕ_V). The crucial difference in deconfounding actor-critic methods is to use r_t+1∼ p(r_t+1|z_t, do(a_t)) given by Equation (<ref>), as opposed to the r_t+1∼ p(r_t+1|z_t, a_t) used by vanilla actor-critic methods.



 §.§ Training
 

As mentioned previously, DRL consists of two steps: learning a deconfounding model and optimizing a policy based on the learned deconfounding model. In step 1 we learn the model by optimizing the objective
given by Equation (<ref>). Once the deconfounding model is learned, we have an estimate of the state transition function p(z_t | z_t-1, a_t-1), as given by the model, and can also calculate the deconfounding reward function p(r_t | z_t, do(a_t)) according to Equation (<ref>). In step 2, we treat the learned deconfounding model as an RL environment like CartPole in OpenAI Gym and generate trajectories/rollouts using the estimated state transition function and deconfounding reward function. These trajectories/rollouts are then used to train the policy by following the gradient given by Equation (<ref>).



 §.§ Implementation Details
 

We used Tensorflow <cit.> for the implementation of our model and the proposed DRL algorithm. Optimization was done using Adam <cit.>. Unless stated otherwise, the setting of all hyperparameters and network architectures can be found in Appendix <ref>. 

To assess the equality of the learned model, we performed two types of tasks: reconstruction and counterfactual reasoning. The reconstructions were performed by first feeding an input sequence x⃗ into the learned inference network, then sampling from the resulting posterior distribution over z⃗ according to Equation (<ref>), and finally feeding those samples into the generative network described in Equation (<ref>) to reconstruct the original observed sequence x⃗. The counterfactual reasoning, that is, the prediction of x_t+1 given an x_t that we have not seen during training, was executed in four steps:

    
  1) Given the new x_t, we estimate a_t and r_t+1 using equations (<ref>) and (<ref>). 
  2) Once we have x_t, a_t, and r_t+1, we estimate z_t using Equation (<ref>). 
    
  3) Given the estimated z_t and another uniformly randomly selected a_t', we can directly compute z_t+1 from Equation (<ref>). 
    
  4) The final step is to reconstruct x_t+1 from z_t+1 and u according to Equation (<ref>). 

By repeating the last two steps, we can counterfactually reason out a sequence of data. 

We may also be interested in estimating the approximate posterior over the confounder u from observed data. We consider two possible scenarios for this. In the easy one we are given the observed data (x⃗, a⃗, r⃗) and we estimate the posterior by using Equation (<ref>). The more challenging scenario involves estimating the posterior from only x⃗ and no action or reward data. In this case, we follow the same steps used in the task of counterfactual reasoning. We first compute a_t and r_t+1 for t=1 to T using Equations (<ref>) and (<ref>), and then estimate the posterior through Equation (<ref>). 



§ EXPERIMENTAL RESULTS


The evaluation of any method that deals with confounders is always challenging due to a lack of real-world benchmark problems with known ground truth. Furthermore, little work has been done before on the task of deconfounding RL. All this creates difficulties in the evaluation of our algorithms and motivates us to develop several new benchmark problems by revising the MNIST dataset <cit.> and by revising two environments in OpenAI Gym <cit.>, CartPole and Pendulum.



 §.§ New Confounding Benchmarks
 

We now describe three new confounding benchmark problems, all of them including a single binary confounder.
The procedure to create these benchmarks is inspired by <cit.>, who synthesized a dataset mimicking healthcare data under realistic conditions (e.g., noisy laboratory measurements, surgeries and drugs affected by patient age and gender, etc.). The data
used is from either popular RL environments in OpenAI Gym such as Pendulum and CartPole or popular machine learning datasets such as MNIST.



  
Confounding Pendulum and CartPole
We revised two environments in OpenAI Gym <cit.>: Pendulum and CartPole. To obtain a confounding version of Pendulum, we selected 100 different screen images of Pendulum and created a synthetic problem in which actions are joint effort[More details can be found at <https://github.com/openai/gym/wiki/Pendulum-v0>.] with values between -2 and 2. We added 20% bit-flip noise to each image and then, a random policy confounded with a binary factor was used to select the action applied at each time step. This is repeated multiple times to produce a large number of 5-step sequences of noisy images. In each generated sequence, one block with three consecutive squares (2 × 2 pixels) is superimposed with the top-left corner of the images in a random starting location. These squares are intended to be analogous to seasonal flu or other ailments that a patient could exhibit, which are independent of the actions taken and which last several time steps <cit.>. The goal is to show that our model can learn long-range patterns since these play an important role in medical applications. We treat the generated sequences of images as the observations x⃗. The training, validation and test sets respectively comprise 140,000, 28,000 and 28,000 sequences of length five.

The key characteristic of confounding Pendulum is the relationship between confounder, actions, and rewards. For simplicity of notation, we denote the confounder by u, the action by a, and the reward by r. In this case, u is a binary variable mimicking socio-economic status (i.e., the rich and the poor). The range of actions a ∈ [-2,2] is grouped into two categories T_1 (i.e., 1 ≤ |a| ≤ 2) and T_2 (i.e., 0 ≤ |a| ≤ 1) representing different treatments[Note that T_1 is a better treatment than T_2 in our setting described in Appendix <ref>.]. The treatment selection depends on u, that is, u determines the probabilities of choosing T_1 and T_2. The reward is defined as

    r = r_o + r_c,

where r_o is the original reward[In fact, here r_o is a function of both a and a state, and we mention only a to emphasize that the confounder affects the action.] in Pendulum as a function of a, and r_c is the extra reward as a function of both u and a, defined by a two-component Gaussian mixture: 

    r_c ∼∑_r_c∈{R_1, R_2}π(r_c|a,u)𝒩(μ_r_c, σ^2),

with σ and μ_r_c fixed and mixing coefficient π(r_c|a,u) determined by both u and a. More details are available in Appendix <ref>. 

Obviously, in the definition above, r depends on a and u. Furthermore, u has an influence on x_t through a_t-1 and z_t, meaning that x_t contains information about u and, therefore, can be viewed as a proxy variable for the confounder. We assume that the influence between the confounder, action and reward is stochastic while generating the data. The reason for this is that, in practice, we do not have an oracle that tells us which treatment is better so we might make wrong decisions sometimes. For example, take the case of kidney stones. Even though treatment a is better than treatment b, there are still some patients choosing treatment b with positive probability in each category.

All the details about the data generation process can be found in Appendix <ref>, where a straightforward analogy is provided as well. 

The confounding CartPole problem is implemented in the exact same manner, except for the action which is now binary and can be naturally divided into two categories[More details can be found at https://github.com/openai/gym/wiki/CartPole-v0.]. Accordingly, the binary confounder determines the probabilities of choosing which of two actions, and the probabilities are set to the same value as those in the confounding Pendulum given in Appendix <ref>.  



  
Confounding MNIST We follow the same procedure to obtain a confounding MNIST problem. However, the definitions of action and the original reward term are now different. Similar to the Healing MNIST dataset <cit.>, the actions encode rotations of digit images with each entry in a⃗ satisfying (-π/4≤ a ≤ +π/4) and the actions being divided into two categories (π/8≤ |a| ≤π/4) and (0≤ |a| < π/8) according to the confounder u. The original reward term r_o is defined as the minus degree between the upright position and the position that the digit rotates to. For example, if the digit rotates to the position of 3 o'clock or 9 o'clock, then both rewards are -π/2. 









 §.§ Performance Analysis of the Deconfounding Model


We assess the performance of the deconfounding model from Figure <ref>, denoted by M_decon, and compare it with a similar alternative model that does not include the confounder u and that is denoted by M_alt. We train M_decon by optimizing Equation (<ref>) but train M_alt using a different loss function which excludes the confounder u and whose full derivation can be found in Appendix <ref>. Both models are separately trained in a minibatch manner on the training set of the confounding dataset (i.e., 140K image sequences of length five). Afterwards, following the steps depicted in Section <ref>, we use each trained model to perform the reconstruction task on the training set, and both reconstruction and counterfactual reasoning tasks on the testing set (i.e., 28K image sequences of length five).  

Figure <ref> presents a comparison of M_decon and M_alt in terms of reconstruction and counterfactual reasoning on confounding Pendulum. The second row is based on M_decon whilst the top row comes from M_alt. The results generated by the deconfounding model are superior to those produced by the model not taking into account the confounder. To be more specific, as shown in the zoom of samples on the bottom row, M_alt generates more blurry images than M_decon because, without modelling the confounder u, M_alt is forced to average over its multiple latent states, resulting in more blurry samples. 

We obtain similar results regarding the samples produced by M_alt and M_decon on the confounding MNIST dataset, as shown in Figure <ref>. Looking closely at the squares on the generated digit samples (inside a yellow box), we observe that M_alt generates non-consecutive white squares in the task of counterfactual reasoning, which does not really make sense because only consecutive patterns appear in the training set. The generated images are also more blurry. By contrast, this does not occur on the samples from M_decon, showing that our deconfounding model is able to cope with long-range patterns and describe the data better. 

Figure <ref> visualizes approximate posterior samples of the 2-dimensional confounder u. We can see that, although the prior distribution of the confounder is assumed to be a factorized standard Gaussian distribution, the model still identifies two clear clusters from the data because u is originally a binary variable. This demonstrates that our model can learn confounders even if the assumed prior is not that accurate. 



 §.§ Comparison of RL Algorithms
 

In this section, we evaluate the proposed deconfounding actor-critic (AC) method by comparing with its vanilla version on confounding Pendulum. In the vanilla AC method, given a learned M_alt, we optimize the policy by calculating the gradient presented in Equation (<ref>) on the basis of the trajectories/rollouts generated through M_alt. Equation (<ref>) involves two functions: V(z_t;ϕ_V) and π(a_t | z_t; θ), both of which are represented using neural networks and their corresponding hyperparameters can be found in Appendix <ref>. It is worth noting that, in this vanilla case, each reward r_t+1 is sampled from the conditional distribution p(r_t+1 | z_t, a_t). By contrast, the proposed deconfounding AC method uses M_decon instead. The optimizaiton of the policy is again performed using the gradient from Equation (<ref>). However, for the deconfounding AC approach, we use trajectories/rollouts generated by M_decon in which each reward r_t+1 is obtained using the interventional distribution p(r_t+1 | z_t, do(a_t)) computed according to Equation (<ref>). For completeness, we also compare with the direct AC method in which the AC method is directly trained on the training data instead of the trajectories/rollouts. 



In the training phase, we respectively run the vanilla AC and the deconfounding AC algorithms over 2000 episodes with 200 time steps each. For the fair comparison, we also run the direct AC algorithm over 2000 episodes, each of which has 200 pairs of state-action-reward randomly selected from the training data. In order to reduce non-stationarity and to decorrelate updates, the generated data is stored in an experience replay memory and then randomly sampled in a batch manner <cit.>. We summarize all the rewards produced during the rollouts in each episode and further average the summarized rewards over a window of 100 episodes to obtain a smoother curve. Figure <ref>(a) shows that
our deconfounding AC algorithm performs significantly better than the vanilla AC algorithm in the confounded environment. Here, the direct AC algorithm is not included because it does not generate rollouts.

In the testing phase, we first randomly select 100 samples from the testing set, each starting a new episode, and then use the learned policies to generate trajectories over 200 time steps as we did during training. From the resulting 100 episodes, we plot the total reward obtained by each method, see in Figure <ref>(b), and also plot the percentage of times that the optimal action T_1 is chosen in each episode, as shown in Figure <ref>(c). Figure <ref>(b) shows that our deconfounding AC obtains on average much higher reward at test time than the baselines. Figure <ref>(c) further tells us that our deconfounding AC is much more likely to choose the optimal action at each time step, whilst the vanilla AC and the direct AC make a wrong decision more than half of the times. 




§ RELATED WORK

<cit.> used deep neural networks to construct nonlinear state space models and leveraged a structured variational approximation parameterized by recurrent neural networks to mimic the posterior distribution. <cit.> reformulated RL and control problems using probabilistic inference, which allows us to bring to bear a large pool of approximate inference methods, and flexibly extend the models. <cit.> exploited continuous state-space models and deep RL to obtain improved treatment policies for septic patients from observational data. <cit.> discussed problems when evaluating RL algorithms in observational health setting. However, all the works mentioned above do not take into account confounders in the proposed models.

<cit.> attempted to learn individual-level causal effects from observational data in the non-temporal setting. They also used a variational auto-encoder to estimate the unknown confounder given a causal graph. <cit.> developed predictive models based on electronic medical records without using causal inference. <cit.> proposed a nonparametric Bayesian method to analyze clinical temporal data. <cit.> represented the treatment response curves using linear time-invariant dynamical systems. This provides a flexible approach to modeling response over time. Although the latter two works model sequential data, they both do not consider RL or causal inference.

<cit.> considered the problem of bandits with unobserved confounders. <cit.> and <cit.> further studied contextual bandits with latent confounders. <cit.> circumvented some problems caused by unobserved confounders in multi-armed bandits by using counterfactual-based decision making. <cit.> leveraged causal inference to tackle the problem of transferring knowledge across bandit agents. Unlike our method, all these works are restricted to bandits, which corresponds to a simplified RL setting without state transitions.

Last but not least, we want to mention some connections between our method and partially observable Markov decision processes (POMDPs). To the best of our knowledge, existing work on POMDPs has not considered problems with confounders yet. Apart from the confounding part, in POMDPs the observation only provides partial information about the actual state, so the agent in POMDPs does not necessarily know which actual state it is in. By contrast, for convenience in practice, we simplify this setting and assume that the observation provides all the information about the actual state the agent is in but with some noise. Hence we only need to denoise the observation to obtain the actual state. In this sense, a more related work to our model is probably the world model <cit.>, because both models used variational inference to estimate the actual state from the observation. If we treat our model as POMDPs with the confounding bias, it will make our model more complicated but it is worth exploring in the future. 

As far as we are concerned, this is the first attempt to build a bridge between confounding and the full RL problem with observational data. 



§ CONCLUSION AND FUTURE WORK


We have introduced deconfounding reinforcement learning (DRL), a general method for 
addressing reinforcement learning problems with observational data in which hidden confounders affect
observed rewards and actions.
We have used DRL to obtain deconfounding variants of actor-critic methods and showed how these new variants perform better than
the original vanilla algorithms on new confounding benchmark problems.
In the future, we will collaborate with hospitals and apply our approach to real-world medical datasets. We also hope that our work will stimulate further investigation of connections between causal inference and RL.








§ VARIATIONAL LOWER BOUND FOR M_DECON
 

    log p_θ(x⃗, a⃗, r⃗) 
    
    =   log∫_u∫_z⃗ p_θ(x⃗, a⃗, r⃗, z⃗, u) dz⃗du 
    ≥   ∫_u∫_z⃗ q_ϕ(z⃗, u|x⃗, a⃗, r⃗)logp_θ(x⃗, a⃗, r⃗, z⃗, u)/q_ϕ(z⃗, u|x⃗, a⃗, r⃗)dz⃗du 
    
    =    ∫_u∫_z⃗ q_ψ(u|x⃗, a⃗, r⃗)q_ϕ(z⃗|x⃗, a⃗, r⃗)logp_θ(x⃗, a⃗, r⃗, z⃗, u)/q_ψ(u|x⃗, a⃗, r⃗)q_ϕ(z⃗|x⃗, a⃗, r⃗)dz⃗du     (factorization assumption)
    
    =   ∫_u∫_z⃗ q_ψ(u|x⃗, a⃗, r⃗)q_ϕ(z⃗|x⃗, a⃗, r⃗) 
       logp(u)p(z_1)[∏_t=1^T p(x_t|z_t,u)p(a_t|z_t, u)p(r_t+1|z_t, a_t, u)][∏_t=2^T p(z_t|z_t-1, a_t-1)]/q_ψ(u|x⃗, a⃗, r⃗)q_ϕ(z⃗|x⃗, a⃗, r⃗)dz⃗du 
    
    =   ∫_u∫_z⃗q(u|x⃗, a⃗, r⃗)q(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗) 
       logp(u)p(z_1)∏_t=1^T p(x_t|z_t,u)p(a_t|z_t, u)p(r_t+1|z_t, a_t, u)∏_t=2^T p(z_t|z_t-1, a_t-1)/q(u|x⃗, a⃗, r⃗)q(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗)dz⃗du 
    
    =   ∑_t=1^T∫_u∫_z_1⋯∫_z_Tq(u|x⃗, a⃗, r⃗)q(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗) 
       log(p(x_t|z_t,u)p(a_t|z_t, u)p(r_t+1|z_t, a_t, u)) dz⃗du 
       +∫_u∫_z_1⋯∫_z_Tq(u|x⃗, a⃗, r⃗)q(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗) logp(u)/q(u|x⃗, a⃗, r⃗)dz⃗du 
       +∫_u∫_z_1⋯∫_z_Tq(u|x⃗, a⃗, r⃗)q(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗) logp(z_1)/q(z_1|x⃗, a⃗, r⃗)dz⃗du 
       + ∑_t=2^T∫_u∫_z_1⋯∫_z_Tq(u|x⃗, a⃗, r⃗)q(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗)logp(z_t|z_t-1, a_t-1)/q(z_t|z_t-1, x⃗, a⃗, r⃗)dz⃗du 
    
    =   ∑_t=1^T∫_u∫_z_t q(u|x⃗, a⃗, r⃗)q(z_t | z_t-1, x⃗, a⃗, r⃗)log(p(x_t|z_t,u)p(a_t|z_t, u)p(r_t+1|z_t, a_t, u)) dz_tdu 
       + ∫_u q(u | x⃗, a⃗, r⃗)logp(u)/q(u|x⃗, a⃗, r⃗)du 
       + ∫_z_1 q(z_1 | x⃗, a⃗, r⃗)logp(z_1)/q(z_1|x⃗, a⃗, r⃗)dz_1 
       + ∑_t=2^T∫_z_t-1∫_z_tq(z_t|z_t-1, x⃗, a⃗, r⃗)logp(z_t|z_t-1, a_t-1)/q(z_t|z_t-1, x⃗, a⃗, r⃗)dz_tdz_t-1


    =   ∑_t=1^T z_t ∼ q(z_t|z_t-1, x⃗,a⃗,r⃗) 
     u ∼ q(u|x⃗,a⃗,r⃗)𝔼[log p(x_t|z_t,u)+log p(a_t|z_t, u)+log p(r_t+1|z_t, a_t, u)] 
       -KL(q(u|x⃗,a⃗,r⃗) || p(u)) 
       -KL(q(z_1|x⃗,a⃗,r⃗) || p(z_1))
       -∑_t=2^T z_t-1∼ q(z_t-1|z_t-2,x⃗,a⃗,r⃗)𝔼[KL(q(z_t|z_t-1, x⃗,a⃗,r⃗) || p(z_t|z_t-1,a_t-1))].




§ VARIATIONAL LOWER BOUND FOR M_ALT
 


    log p_θ(x⃗, a⃗, r⃗) 
    
    =   log∫_z⃗ p_θ(x⃗, a⃗, r⃗, z⃗) dz⃗
    ≥   ∫_z⃗ q_ϕ(z⃗|x⃗, a⃗, r⃗)logp_θ(x⃗, a⃗, r⃗, z⃗)/q_ϕ(z⃗|x⃗, a⃗, r⃗)dz⃗
    
    =   ∫_z⃗ q_ϕ(z⃗|x⃗, a⃗, r⃗) logp(z_1)[∏_t=1^T p(x_t|z_t)p(a_t|z_t)p(r_t+1|z_t, a_t)][∏_t=2^T p(z_t|z_t-1, a_t-1)]/q_ϕ(z⃗|x⃗, a⃗, r⃗)dz⃗
    
    =   ∫_z⃗q(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗) 
       logp(z_1)∏_t=1^T p(x_t|z_t)p(a_t|z_t)p(r_t+1|z_t, a_t)∏_t=2^T p(z_t|z_t-1, a_t-1)/q(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗)dz⃗
    
    =   ∑_t=1^T∫_z_1⋯∫_z_Tq(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗) log(p(x_t|z_t)p(a_t|z_t)p(r_t+1|z_t, a_t)) dz⃗
       +∫_z_1⋯∫_z_Tq(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗) logp(z_1)/q(z_1|x⃗, a⃗, r⃗)dz⃗
       + ∑_t=2^T∫_z_1⋯∫_z_Tq(z_1|x⃗, a⃗, r⃗)⋯ q(z_T|z_T-1, x⃗, a⃗, r⃗)logp(z_t|z_t-1, a_t-1)/q(z_t|z_t-1, x⃗, a⃗, r⃗)dz⃗
    
    =   ∑_t=1^T∫_z_t q(z_t | z_t-1, x⃗, a⃗, r⃗)log(p(x_t|z_t)p(a_t|z_t)p(r_t+1|z_t, a_t)) dz_t 
       + ∫_z_1 q(z_1 | x⃗, a⃗, r⃗)logp(z_1)/q(z_1|x⃗, a⃗, r⃗)dz_1 
       + ∑_t=2^T∫_z_t-1∫_z_tq(z_t|z_t-1, x⃗, a⃗, r⃗)logp(z_t|z_t-1, a_t-1)/q(z_t|z_t-1, x⃗, a⃗, r⃗)dz_tdz_t-1
    
    =   ∑_t=1^T z_t ∼ q(z_t|z_t-1, x⃗,a⃗,r⃗)𝔼[log p(x_t|z_t)+log p(a_t|z_t)+log p(r_t+1|z_t, a_t)] 
       -KL(q(z_1|x⃗,a⃗,r⃗) || p(z_1))
       -∑_t=2^T z_t-1∼ q(z_t-1|z_t-2,x⃗,a⃗,r⃗)𝔼[KL(q(z_t|z_t-1, x⃗,a⃗,r⃗) || p(z_t|z_t-1,a_t-1))].






§ BI-DIRECTIONAL LSTM
 
In our inference model, we use a similar architecture of bi-directional LSTM to that in <cit.>. Apart from different inputs, the main difference is to introduce a_t-1 to computing the combined hidden feature as shown in the following formula.


    h_combined   = 1/4(tanh(W_z · z_t-1+b_z) + tanh(W_a · a_t-1+b_a) + h_t^left + h_t^right) 
    μ_t    = W_μ· h_combined + b_μ
    σ_t^2    = softplus(W_σ^2· h_combined + b_σ^2)




§ EXAMPLE OF A PROXY VARIABLE
 

For convenience, we directly move Figure 1 of <cit.> here, as shown in Figure <ref>.





§ ANALYSIS OF THE DECONFOUNDING MODEL
 

At each time step t, from the deconfounding model presented in Figure <ref>, we can extract two 4-tuple components, centering at z_t (Figure <ref>) and u (Figure <ref>), respectively. It is apparent to observe that both Figure <ref> and Figure <ref> share the exact same structure with Figure <ref>.







§ CONFOUNDING DATASETS
 

The probabilities between the confounder, action and reward are presented as follows, 

    p(u=0) = 0.8,    p(u=1)=0.2 ; 
       p(a=T_1 | u=0) = 0.24,    p(a=T_2 |u =0) = 0.76; 
       p(a=T_1 | u=1) = 0.77,    p(a=T_2 |u =1) = 0.23; 
       p(r_c=R_1 | a=T_1, u=0) = 0.93,    p(r_c=R_2 | a=T_1, u=0) = 0.07;  
       p(r_c=R_1 | a=T_2, u=0) = 0.87,    p(r_c=R_2 | a=T_2, u=0) = 0.13;  
       p(r_c=R_1 | a=T_1, u=1) = 0.73,    p(r_c=R_2 | a=T_1, u=1) = 0.27;  
       p(r_c=R_1 | a=T_2, u=1) = 0.69,    p(r_c=R_2 | a=T_2, u=1) = 0.31;  
       μ_r_c=R_1 = -1,   μ_r_c=R_2 = -200,   σ=2.

Note that, in our model r_c is a function of u, a_t, and z_t. We omit the notation of the state in the conditional probability of r_c because z_t is fixed at each time step, which does not affect the probability. To help readers understand this setting, we can make a straightforward analogy. More specifically, the confounder u stands for the socio-economic status where u=0 represents the poor and u=1 the rich. Treatment T_1 is more expensive but much better than Treatment T_2. R_1 means a healthier feedback/recovery than R_2. Equation (<ref>) says that the poor are much more than the rich (it is indeed in reality). Equation (<ref>) and Equation (<ref>) tell us the fact that the poor tend to choose the cheaper but worse treatment T_2 whilst the rich prefer the more expensive but better treatment T_1, which also makes sense in the real world. Finally, Equations (<ref>)-(<ref>) reveal that T_1 is a better treatment than T_2 within both poor and rich population. 



§ EXPERIMENTAL SETTINGS
 



[h]

Function     Architecture  
 
2cDeconfounding Model 

f_1, f_2     FC_512→Conv_7, 32→Conv_14, 16→Conv_28, 1→FC_784→{sigmoid, softplus} 
 
 f_3, f_4     {FC_100, FC_100}→FC_200→FC_200→FC_1→{tanh, softplus}  
  
 f_5, f_6     {FC_100, FC_100, FC_100}→FC_100→FC_100→FC_1→{sigmoid, softplus}  
 
 f_7, f_8     {FC_100, FC_100}→FC_100→FC_100→FC_50→{None, softplus}  
  
 f_9, f_10     { [Conv_5,16→Conv_5,32→Conv_5,32→FC_100], FC_100, FC_100}→FC_100→FC_100→LSTM_100, 5→FC_1→{None, softplus} 
 
 f_11, f_12      { [Conv_5,16→Conv_5,32→Conv_5,32→FC_100], FC_100, FC_100}→FC_100→FC_100→LSTM_100, 5→FC_50→{None, softplus} 
 
 f_13, f_14     Conv_5,16→Conv_5,32→Conv_5,32→FC_1→{tanh, softplus}  
 
 f_15, f_16      {[Conv_5,16→Conv_5,32→Conv_5,32→FC_100], FC_100}→FC_100→FC_1→{sigmoid, softplus} 
  
 2cDeconfounding AC 

 V(z_t;ϕ_V)     FC_300→FC_300→{None, softplus} 

 π(a_t | z_t; θ)     FC_300→FC_300→{tanh, softplus} 

 
Architectures for deconfounding reinforcement learning. Here FC_k stands for a fully-connected layer with k units, Conv_k,n for a convolution layer with n filters of size k × k, LSTM_n, t for a LSTM layer rolling out for t steps with latent size of n, {·} for the parallel operators, and [·] for the sequential operators. In default, FC_k and Conv_k,n are followed by a softplus activation layer and a batch-norm layer, which are omitted here for simplicity. Note that, in our setting, the two functions in each pair {f_2i-1, f_2i}_i=1,…,8 share the same parameters except for the last layer.




