Dynamic Causal Bayesian Optimization
Year: 2021

Abstract:
This paper studies the problem of performing a sequence of optimal interventions in a causal dynamical system where both the target variable of interest and the inputs evolve over time. This problem arises in a variety of domains e.g. system biology and operational research. Dynamic Causal Bayesian Optimization (DCBO) brings together ideas from sequential decision making, causal inference and Gaussian process (GP) emulation. DCBO is useful in scenarios where all causal effects in a graph are changing over time. At every time step DCBO identifies a local optimal intervention by integrating both observational and past interventional data collected from the system. We give theoretical results detailing how one can transfer interventional information across time steps and define a dynamic causal GP model which can be used to quantify uncertainty and find optimal interventions in practice. We demonstrate how DCBO identifies optimal interventions faster than competing approaches in multiple settings and applications.

text: 1 Introduction Causal Dimension
Solving decision making problems in a variety of domains requires understanding of cause-effect relationships in a system. This can be obtained by experimentation. However, deciding how to intervene at every point in time is particularly complex in dynamical systems, due to the evolving nature of causal effects. For instance, companies need to decide how to allocate scarce resources across different quarters. In system biology, scientists need to identify genes to knockout at specific points in time. This paper describes a probabilistic framework that finds such optimal interventions over time.
Focusing on a specific example, consider a setting in which Yt denotes the unemployment-rate of an economy at time t, Zt is the economic growth and Xt the inflation rate. Fig. 1a depicts the causal graph [26] representing an agent’s understanding of the causal links between these variables. The agent aims to determine, at each time step t ∈ {0, 1, 2}, the optimal action to perform in order to minimize the current unemployment rate
∗Denotes equal contribution
35th Conference on Neural Information Processing Systems (NeurIPS 2021).
Yt while accounting for the intervention cost. The investigator could frame this setting as a sequence of global optimization problems and find the solutions by resorting to Causal Bayesian Optimization [CBO, 2]. CBO extends Bayesian Optimization [BO, 30] to cases in which the variable to optimize is part of a causal model where a sequence of interventions can be performed. However, CBO does not account for the system’s temporal evolution thus breaking the time dependency structure existing among variables (Fig. 1b). This will lead to sub-optimal solutions, especially in non-stationary scenarios. The same would happen when using Adaptive Bayesian Optimization [ABO, 25] (Fig. 1c) or BO (Fig. 1d). ABO captures the time dependency of the objective function but neither considers the causal structure among inputs nor their temporal evolution. BO disregards both the temporal and the causal structure. Our setting differs from both reinforcement learning (RL) and the multi-armed bandits setting (MAB). Differently from MAB we consider interventions on continuous variables where the dynamic target variable has a non-stationary interventional distribution. In addition, compared to RL, we do not model the state dynamics explicitly and allow the agent to perform a number of explorative interventions which do not change the underlying state of the system, before selecting the optimal action. We discuss these points further in §1.1.
Dynamic Causal Bayesian Optimization2, henceforth DCBO, accounts for both the causal relationships among input variables and the causality between inputs and outputs which might evolve over time. DCBO integrates CBO with dynamic Bayesian networks (DBN), offering a novel approach for decision making under uncertainty within dynamical systems. DBN [19] are commonly used in time-series modelling and carry dependence assumptions that do not imply causation. Instead, in probabilistic causal models [27], which form the basis for the CBO framework, graphs are buildt around causal information and allow us to reason about the effects of different interventions. By combining CBO with DBNs, the proposed methodology finds an optimal sequence of interventions which accounts for the causal temporal dynamics of the system. In addition, DCBO takes into account past optimal interventions and transfers this information across time, thus identifying the optimal intervention faster than competing approaches and at a lower cost. We make the following contributions:
• We formulate a new class of optimization problems called Dynamic Causal Global Optimization (DCGO) where the objective functions account for the temporal causal dynamics among variables.
• We give theoretical results demonstrating how interventional information can be transferred across time-steps depending on the topology of the causal graph.
• Exploiting our theoretical results, we solve the optimization problem with DCBO. At every time step, DCBO constructs surrogate models for different intervention sets by integrating various sources of data while accounting for past interventions.
• We analyze DCBO performance in a variety of settings comparing against CBO, ABO and BO.

heading: 1.1 Related Work
text: Dynamic Optimization Optimization in dynamic environments has been studied in the context of evolutionary algorithms [14, 16]. More recently, other optimization techniques [28, 32, 10] have been adapted to dynamic settings, see e.g. [9] for a review. Focusing on BO, the literature on dynamic settings [3, 7, 25] is limited. The dynamic BO framework closest to this work is given by Nyikosa et al. [25] and focuses on functions defined on continuous spaces that follow a more complex behaviour than a simple Markov model. ABO treats the inputs as fixed and not as random variables, thereby disregarding their temporal evolution and, more importantly, breaking their causal dependencies.
Causal Optimization Causal BO [CBO, 2] focuses instead on the causal aspect of optimization and solves the problem of finding an optimal intervention in a DAG by modelling the intervention functions with single GPs or a multi-task GP model [1]. CBO disregards the existence of a temporal evolution in both the inputs and the output variable, treating them as i.i.d. overtime. While disregarding time significantly simplifies the problem, it prevents the identification of an optimal intervention at every t.
Bandits and RL In the broader decision-making literature, causal relationships have been considered in the context of bandits [4, 20–22] and reinforcement learning [23, 8, 13, 36, 24]. In these cases, actions or arms, correspond to interventions on a causal graph where there exists complex relationships between the agent’s decisions and the received rewards. While dynamic settings have been considered in acausal bandit algorithms [5, 33, 35], causal MAB have focused on static settings. Dynamic settings are instead considered by RL algorithms and formalized through Markov decision processes (MDP).
2A Python implementation is available at: https://github.com/neildhir/DCBO.
In the current formulation, DCBO does not consider an MDP as we do not have a notion of state and therefore do not require an explicit model of its dynamics. The system is fully specified by the causal model. As in BO, we focus on identifying a set of time-indexed optimal actions rather than an optimal policy. We allow the agent to perform explorative interventions that do not lead to state transitions. More importantly, differently from both MAB and RL, we allow for the integration of both observational and interventional data. An expanded discussion on the reason why DCBO should be used and the links between DCBO, CBO, ABO and RL is included in the supplement (§8).

heading: 2 Background and Problem Statement
text: Let random variables and values be denoted by upper-case and lower-case letters respectively. Vectors are represented shown in bold. do(X = x) represents an intervention on X whose value is set to x. p(Y | X = x) represents an observational distribution and p(Y | do(X = x)) represents an interventional distribution. This is the distribution of Y obtained by intervening on X and fixing its value to x in the data generating mechanism (see Fig. 2), irrespective of the values of its parents. Evaluating p(Y | do(X = x)) requires “real” interventions while p(Y | X = x) only requires “observing” the system. DO and DI denote observational and interventional datasets respectively. Consider a structural causal model (SCM) defined in Definition 1. Definition 1. (Structural Causal Model) [27, p. 203]. A structural causal model M is a triple 〈U,V, F )〉 where U is a set of background variables (also called exogenous), that are determined by factors outside of the model. V is a set {V1, V2, . . . , V|V|} of observable variables (also called endogenous), that are determined by variables in the model (i.e., determined by variables in U ∪V). F is a set of functions {f1, f2, . . . , fn} such that each fi is a mapping from the respective domains of Ui ∪ Pa(Vi) to Vi, where Ui ⊆ U and Pa(Vi) ⊆ V \ Vi and the entire set F forms a mapping from U to V. In other words, each {fi ∈ vi ← fi(Pa(vi) , ui) | i = 1, . . . , n}, assigns a value to Vi that depends on the values of the select set of variables (Ui ∪ Pa(Vi)). M is associated to a directed acyclic graph (DAG) G, in which each node corresponds to a variable and the directed edges point from members of Pa(Vi) and Ui to Vi. We assume G to be known and leave the integration with causal discovery [15] methods for future work. Within V, we distinguish between three different types of variables: non-manipulative variables C, which cannot be modified, treatment variables X that can be set to specific values and output variable Y that represents the agent’s outcome of interest. Exploiting the rules of do-calculus [27] one can compute p(Y | do(X = x)) using observational data. This often involves evaluating intractable integrals which can be approximated by using observational data to get a Monte Carlo estimate p̂(Y | do(X = x)) ≈ p(Y | do(X = x)). These approximations will be consistent when the number of samples drawn from p(V) is large.
Causality in time One can encode the existence of causal mechanisms across time steps by explicitly representing these relationships with edges in an extended graph denoted by G0:T . For instance, the DAG in Fig. 1(a) can be seen as one of the DAGs in Fig. 1(b) propagated in time. The DAG in Fig. 1(a) captures both the causal structure existing across time steps and the causal mechanism within every “time-slice” t [19]. In order to reason about interventions that are implemented in a sequential manner, that is at time t we decide which intervention to perform in the system and so define: Definition 2. Mt is the SCM at time step t defined asMt = 〈U0:t,V0:t,F0:t〉 where 0 : t denotes the union of the corresponding variables or functions up to time t (see Fig. 2). V0:t includes X0:t = Xt, Y0:t = Yt and C0:t = Ct ∪C0:t−1. The functions in F0:t corresponding to intervened variables are replaced by constant values while the exogenous variables related to them are excluded from U0:t. Definition 3. Gt is the causal graph associated to Mt. In Gt, the incoming edges in variables intervened at 0 : t− 1 are mutilated while intervened variables are represented by deterministic nodes (squares) – see Fig. 2.
Dynamic Causal Global Optimization (DCGO) The goal of this work is to find a sequence of interventions, optimizing a target variable, at each time step, in a causal DAG. Given Gt and Mt, at every time step t, we wish to optimize Yt by intervening on a subset of the manipulative variables Xt. The optimal intervention variables X?s,t and intervention levels x ? s,t are given by:
X?s,t,x ? s,t = argmin Xs,t∈P(Xt),xs∈D(Xs,t) E[Yt | do(Xs,t = xs,t) ,1t>0 · I0:t−1] (1)
where I0:t−1 = ⋃t−1 i=0 do ( X?s,i = x ? s,i ) denotes previous interventions, 1t>0 is the indicator function and P(Xt) is the power set of Xt. D(Xs,t) represents the interventional domain of Xs,t. In the sequel we denote the previously intervened variables by IV0:t−1 = ⋃t−1 i=0 X ? s,i and implemented intervention levels by IL0:t−1 = ⋃t−1 i=0 x ? s,i. The cost of each intervention is given by cost(Xs,t,xs,t). In order to solve the problem in Eq. (1) we make the following assumptions : Assumptions 1. Denote by G(t) the causal graph including variables at time t in G0:T and let Y PTt = Pa(Yt) ∩ Y0:t−1 be the set of variables in G0:T that are both parents of Yt and targets at previous time step. Let the set Y PNTt = Pa(Yt) \Y PTt be the complement and denote by fYt(·) the functional mapping for Yt in Mt. We make the following assumptions:
1. Invariance of causal structure: G(t) = G(0),∀t > 0. 2. Additivity of fYt(·) that is Yt = fYt(Pa(Yt)) + with fYt(Pa(Yt)) = fYY (Y PTt ) + fNYY (Y PNTt )
where fYY and f NY Y are two generic unknown functions and ∼ N (0, σ2).
3. Absence of unobserved confounders in G0:T .
Assumption (3) implies the absence of unobserved confounders at every time step. For instance, this is the case in Fig. 1a. Still in this DAG, Assumption (2) implies fYt(Pa(Yt)) = f Y Y (Yt−1) + fNYY (Zt) + Yt ,∀t > 0. Finally, Assumption (1) implies the existence of the same variables at every time step and a constant orientation of the edges among them for t > 0.
Notice that Assumptions 1 imply invariance of the causal structure within each time-slice, i.e. the structure, edges and vertices, concerning the nodes with the same time index. This means that, across time steps, both the graph and the functional relationships can change. Therefore, not only can the causal effects change significantly across time steps, but also the input dimensionality of the causal functions we model, might change. For instance, in the DAG of Fig. 3(c), the target function for Y2 has dimensionality 3 and a function fYt(·) that is completely different from the one assumed for Y1 that has only two parents. We can thus model a wide variety of settings and causal effects despite this
assumption. Furthermore, even though we assume an additive structure for the functional relationship on Y , the use of GPs allow us to have flexible models with highly non-linear causal effects across different graph structures. In the causality literature, GP models are well established and have shown good performances compared to parametric linear and non-linear models (see e.g. [31, 34, 37]). The sum of GPs gives a flexible and computationally tractable model that can be used to capture highly non-linear causal effects while helping with interpretability [12, 11].

heading: 3 Methodology
text: In this section we introduce Dynamic Causal Bayesian Optimization (DCBO), a novel methodology addressing the problem in Eq. (1). We first study the correlation among objective functions for two consecutive time steps and use it to derive a recursion formula that, based on the topology of the graph, expresses the causal effects at time t as a function of previously implemented interventions (see square nodes in Fig. 2). Exploiting these results, we develop a new surrogate model for the objective functions that can be used within a CBO framework to find the optimal sequence of interventions. This model enables the integration of observational data and interventional data, collected at previous time-steps and interventional data collected at time t, thereby accelerating the identification of the current optimal intervention.

heading: 3.1 Characterization of the time structure in a DAG with time dependent variables
text: The following result provides a theoretical foundation for the dynamic causal GP model introduced later. In particular, it derives a recursion formula allowing us to express the objective function at time t as a function of the objective functions corresponding to the optimal interventions at previous time steps. The proof is given in the appendix (§2). Definition 4. Consider a DAG G0:T and the objective function E[Yt | do(Xs,t = xs,t) , I0:t−1] for a generic time step t ∈ {0, . . . , T}. Denote by Y PTt = (Pa(Yt) ∩ Y0:t−1) the parents of Yt that are targets at previous time steps and by Y PNTt = Pa(Yt) \Y PTt the remaining parents. For any Xs,t ∈ P(Xt) and IV0:t−1 ⊆ X0:t−1 we define the following sets:
• XPYs,t = Xs,t ∩ Pa(Yt) includes the variables in Xs,t that are parents of Yt.
• IPY0:t−1 = I V 0:t−1 ∩ Pa(Yt) includes the variables in IV0:t−1 that are parents of Yt.
• W ⊂ Pa(Yt) such that Pa(Yt) = (Pa(Yt) ∩ Y0:t−1) ∪XPYs,t ∪ IPY0:t−1 ∪W . W includes variables that are parents of Yt but are not targets nor intervened variables.
The values of I0:t−1, XPYs,t , I PY 0:t−1 and W will be denoted by i, x PY, iPY and w respectively. Theorem 1. Time operator. Consider a DAG G0:T and the related SCM satisfying Assumptions 1. It is possible to prove that, ∀Xs,t ∈ P(Xt), the intervention function fs,t(x) = E[Yt | do(Xs,t = x) ,1t>0 · I0:t−1] with fs,t(x) : D(Xs,t)→ R can be written as:
fs,t(x) = f Y Y (f ?) + Ep(w|do(Xs,t=x),i) [ fNYY (x PY, iPY,w) ]
(2)
where f? = {E [ Yi|do ( X?s,i = x ? s,i ) , I0:i−1 ] }Yi∈Y PTt that is the set of previously observed optimal targets that are parents of Yt. fYY denotes the function mapping Y PT t to Yt and f NY Y represents the function mapping Y PNTt to Yt. Eq. (2) reduces to Ep(w| do(Xs,t=x),i) [ fNYY (x PY, iPY,w) ]
when Yt does not depend on previous targets. This is the setting considered in CBO that can be thus seen as a particular instance of DCBO. Exploiting Assumptions (1), it is possible to further expand the second term in Eq. (2) to get the following expression. A proof is given in the supplement (§2). Corollary 1. Given Assumptions 1, we can write:
Ep(w| do(Xs,t=x),i) [ fNYY (x PY, iPY,w) ] = Ep(U0:t) [ fNYY (x PY, iPY, {C(W )}W∈W) ]
(3) where p(U0:t) is the distribution for the exogenous variables up to time t and C(W ) is given by:
C(W ) =    fW (uW ,x PW, iPW) if R = ∅ fW (uW ,x PW, iPW, r) if R ⊆ Xs,t ∪ IV0:t−1
fW (uW ,x PW, iPW, C(R)) if R 6⊆ Xs,t ∪ IV0:t−1
where fW represents the functional mapping for W in the SCM and uW is the set of exogenous variables with edges into W . xPW and iPW are the values corresponding to XPWs,t and I PW 0:t−1 which in turn represent the subset of variables in Xs,t and IV0:t−1 that are parents of W . Finally r is the value of R = Pa(W ) \(XPYs,t ∪ IPW0:t−1).
Examples for Eq. (2): For the DAG in Fig. 1(a), at time t = 1 and with IV0:t−1 = {Z0}, we have E[Y | do(Z1 = z) , I0] = fYY (y?0) + fNYY (z). Indeed in this case W = ∅, xPY = z and f? = {y?0 = E[Y0|do(Z0 = z0)]}. Still at t = 1 and with IV0:t−1 = {Z0}, the objective function for Xs,t = {X1} can be written as fYY (y?0) + Ep(z1| do(X1=x),I0) [ fNYY (z1) ] as W = {Z1}. All derivations for these expressions and alternative graphs are given in the supplement (§2).

heading: 3.2 Restricting the search space
text: The search space for the problem in Eq. (1) grows exponentially with |Xt| thus slowing down the identification of the optimal intervention when G includes more than a few nodes. Indeed, a naive approach of finding X?s,t at t = 0, . . . , T would be to explore the 2
|Xt| sets in P(Xt) at every t and keep 2|Xt| models for the objective functions. In the static setting, CBO reduces the search space by exploiting the results in [21]. In particular, it identifies a subset of variables M ⊆ P(X) worth intervening on thus reducing the size of the exploration set to 2|M|.
In our dynamic setting, the objective functions change at every time step depending on the previously implemented interventions and one would need to recompute M at every t. However, it is possible to show that, given Assumptions 1, the search space remains constant over time. Denote by Mt the set M at time t and let M0 represent the set at t = 0 which corresponds to M computed in CBO. For t > 0 it is possible to prove that: Proposition 3.1. MIS in time. If Assumptions 1 are satisfied, Mt = M0 for t > 0.

heading: 3.3 Dynamic Causal GP model
text: Here we introduce the Dynamic Causal GP model that is used as a surrogate model for the objective functions in Eq. (1). The prior parameters are constructed by exploiting the recursion in Eq. (2). At each time step t, the agent explores the sets in Mt ⊆ P(Xt) by selecting the next intervention to be the one maximizing a given acquisition function. The DCBO algorithm is shown in Algorithm 1.
Prior Surrogate Model At each time step t and for each Xs,t ∈ Mt, we place a GP prior on the objective function fs,t(x) = E[Yt|do(Xs,t = x) ,1t>0 · I0:t−1]. We construct the prior parameters exploiting the recursive expression in Eq. (2):
fs,t(x) ∼ GP(ms,t(x), ks,t(x,x′)) where ms,t(x) = E [ fYY (f ?) + Ê[fNYY (xPY, iPY,w)] ]
ks,t(x,x ′) = kRBF(x,x ′) + σs,t(x)σs,t(x ′) with
σs,t(x) = √ V[fYY (f?) + Ê [ fNYY (x PY, iPY,w) ]
and kRBF(x,x′) := exp(− ||x−x ′||2
2l2 ) represents the radial basis function kernel [29]. We have it that
Ê [ fNYY (x PY, iPY,w) ] = Êp(w| do(Xs,t=x),i) [ fNYY (x PY, iPY,w) ]
represents the expected value of fNYY (x PY, iPY,w) with respect to p(w | do(Xs,t = x) , i) which is estimated via the do-calculus using observational data. The outer expectation in ms,t(x) and the variance in σs,t(x) are computed with respect to p(fYY , f NY Y ) which is also estimated using observational data. In this work we model fYY , f NY Y and all functions in the SCM by independent GPs.
Both ms,t(x) and σs,t(x) can be equivalently written by exploiting the equivalence in Eq. (3). In both cases, this prior construction allows the integration of three different types of data: observational data, interventional data collected at time t and the optimal interventional data points collected in the past. The former is used to estimate the SCM model and p(w | do(Xs,t = x) , i) via the rules of
do-calculus. The optimal interventional data points at 0 : t− 1 determine the shift fYY (f?) while the interventional data collected at time t are used to update the prior distribution on fs,t(x). Similar prior constructions were previously considered in static settings [2, 1] where only observational and interventional data at the current time step were used. The additional shift term appears here as there exists causal dynamics in the target variables and the objective function is affected by previous decisions. Fig. 2 in the appendix shows a synthetic example in which accounting for the dynamic aspect in the prior formulation leads to a more accurate GP posterior compared to the baselines, especially when the the optimum location changes across time steps.
Algorithm 1: DCBO Data: DO, {DIs,t=0}s∈{0,...,|M0|}, G0:T , H . Result: Optimal intervention path {X?s,t,x?s,t, y?t }Tt=1 Initialise: M, DI0 and initial optimal DI? = ∅. for t = 0, . . . , T do
1. Initialise dynamic causal GP models for all Xs,t ∈Mt using DI?,t−1 if t > 0. 2. Initialise interventional dataset {DIs,t}s∈{0,...,|Mt|} for h = 1, . . . ,H do 1. Compute EIs,t(x) for each Xs,t ∈Mt. 2. Obtain (s?, α?) 3. Intervene and augment DIs=s?,t 4. Update posterior for fs=s?,t end 3. Return the optimal intervention (X?s,t,x ? s,t) 4. Append optimal interventional data DI?,t = DI?,t−1 ∪ ((X?s,t,x?s,t), y?t )
end
Likelihood Let DIs,t = (XI ,YIs,t) be the set of interventional data points collected for Xs,t with XI being a vector of intervention values and YIs,t representing the corresponding vector of observed target values. As in standard BO we assume each ys,t in YIs,t to be a noisy observation of the function fs,t(x) that is ys,t(x) = fs,t(x) + s,t with s,t ∼ N (0, σ2) for s ∈ {1, . . . , |Mt|} and t ∈ {0, . . . , T}. In compact form, the joint likelihood function for DIs,t is p(YIs,t | fs,t, σ2) = N (fs,t(XI), σ2I).
Acquisition Function Given our surrogate models at time t, the agent selects the interventions to implement by solving a Causal Bayesian Optimization problem [2]. The agent explores the sets in Mt and decides where to intervene by maximizing the Causal Expected Improvement (EI). Denote by y?t the optimal observed target value in {YIs,t}|Mt|s=1 that is the optimal observed target across all intervention sets at time t. The Causal EI is given by
EIs,t(x) = Ep(ys,t)[max(ys,t − y?t , 0)]/cost(Xs,t,xs,t). Let α1, . . . , α|Mt| be solutions of the optimization of EIs,t(x) for each set in Mt and α? := max{α1, . . . , α|Mt|}. The next best intervention to explore at time t is given by s? = argmaxs∈{1,··· ,|Mt|} αs. Therefore, the set-value pair to intervene on is (s
?, α?). At every t, the agent implement H explorative interventions in the system which are selected by maximizing the Causal EI. Once the budget H is exhausted, the agent implements what we call the decision intervention It, that is the optimal intervention found at the current time step, and move forward to a new optimization at t+ 1 carrying the information in y?0:t−1. The parameter H determines the level of exploration of the system and acts as a budget for the CBO algorithm. Its value is determined by the agent and is generally problem specific.
Posterior Surrogate Model For any set Xs,t ∈Mt, the posterior distribution p(fs,t | DIs,t) can be derived analytically via standard GP updates. p(fs,t | DIs,t) will also be a GP with parameters
ms,t(x | DIs,t) = ms,t(x) + ks,t(x,XI)[ks,t(XI ,XI) + σ2I](YIs,t −ms,t(XI) and ks,t(x,x ′ | DIs,t) = ks,t(x,x′)− ks,t(x,XI)[ks,t(XI ,XI) + σ2I]ks,t(XI ,x′).

heading: 4 Experiments
text: We evaluate the performance of DCBO in a variety of synthetic and real world settings with DAGs given in Fig. 3. We first run the algorithm for a stationary setting where both the graph structure and the SCM do not change over time (STAT.). We then consider a scenario characterised by increased observation noise (NOISY) for the manipulative variables and a settings where observational data are missing at some time steps (MISS.). Still assuming stationarity, we then test the algorithm in a DAG where there are multivariate interventions in Mt (MULTIV.). Finally, we run DCBO for a non-stationary graph where both the SCM and the DAG change over time (NONSTAT.). To conclude,
we use DCBO to optimize the unemployment rate of a closed economy (DAG in Fig. 3d, ECON.) and to find the optimal intervention in a system of ordinary differential equation modelling a real predator-prey system (DAG in Fig. 3e, ODE). We provide a discussion on the applicability of DCBO to real-world problems in §7 of the supplement together with all implementation details.
Baselines We compare against the algorithms in Fig. 1. Note that, by constructions, ABO and BO intervene on all manipulative variables while DCBO and CBO explore only Mt at every t. In addition, both DCBO and ABO reduce to CBO and BO at the first time step. We assume the availability of an observational dataset DO and set a unit intervention cost for all variables. Performance metric We run all experiments for 10 replicates and show the average convergence path at every time step. We then compute the values of a modified “gap” metric3 across time steps and with standard errors across replicates. The metric is defined as
Gt = [ y(x?s,t)− y(xinit) y? − y(xinit) + H −H(x?s,t) H ]/( 1 + H − 1 H ) (4)
where y(·) represents the evaluation of the objective function, y? is the global minimum, and xinit and x?s,t are the first and best evaluated point, respectively. The term H−H(x?s,t) H with H(x ? s,t) denoting the number of explorative trials needed to reach x?s,t captures the speed of the optimization. This term is equal to zero when the algorithm is not converged and equal to (H − 1)/H when the algorithm converges at the first trial. We have 0 ≤ Gt ≤ 1 with higher values denoting better performances. For each method we also show the average percentage of replicates where the optimal intervention set X?s,t is identified.

heading: 4.1 Synthetic Experiments
text: Stationary DAG and SCM (STAT.) We run the algorithms for the DAG in Fig. 1(a) with T = 3 and N = 10. For t > 0, DCBO converges to the optimal value faster than competing approaches (see Fig. 2 in the supplement, right panel, 3rd row). DCBO identifies the optimal intervention set in 93% of the replicates (Table 2) and reaches the highest average gap metric (Table 1). In this experiment the location of the optimum changes significantly both in terms of optimal set and intervention value when going from t = 0 to t = 1. This information is incorporated by DCBO through the prior dependency on y?0:t−1. In addition, ABO performance improves over time as it accumulates interventional data and uses them to fit the temporal dimension of the surrogate model. This benefits ABO in a stationary setting but might penalise it in non-stationary setting where the objective functions change significantly.
Noisy manipulative variables (NOISY): The benefit of using DCBO becomes more apparent when the manipulative variables observations are noisy while the evolution of the target variable is more accurately detected. In this case both the convergence of DCBO and CBO are slowed down by noisy observations which are diluting the information provided by the do-calculus making the priors less informative. However, the DCBO prior dependency on y?0:t−1 allows it to correctly identify the shift in the target variable thus improving the prior accuracy and the speed-up of the algorithm (Fig. 4).
Missing observational data (MISS.) Incorporating dynamic information in the surrogate model allows us to efficiently optimise a target variable even in setting where observational data are missing.
3This metric is a modified version of the one used in [18].
We consider the DAG in Fig. 1(a) with T = 6, N = 10 for the first three time steps and N = 0 afterwards. DCBO uses the observational distributions learned with data from the first three time steps to construct the prior for t > 3. On the contrary, CBO uses the standard prior for t > 3. In this setting DCBO consistently outperforms CBO at every time step. However, ABO performance improves over time and outperforms DCBO starting from t = 4 due to its ability to exploit all interventional data collected over time (see Fig. 3 in the supplement).
Multivariate intervention sets (MULTIV.) When the optimal intervention set is multivariate, both DCBO and CBO convergence speed worsen. For instance, for the DAG in Fig. 3a, |M| = 5 thus both CBO and DCBO will have to perform more explorative interventions before finding the optimum. At the same time, ABO and BO consider interventions only on {Wt, Xt, Zt},∀t and need to explore an even higher intervention space. The performance of all methods decrease in this case (Table 1) but DCBO still identifies the optimal intervention set in 93% of the replicates (Table 2).
Independent manipulative variables (IND.): Having to explore multiple intervention sets significantly penalises DCBO and CBO when there is no causal relationship among manipulative variables which are also the only parents of the target. This is the case for the DAG in Fig. 3b where the optimal intervention is {Xt, Zt} at every time step. In this case, exploring M and propagating uncertainty in the causal prior slows down DCBO convergence and decreases both its performance (Table 1) and capability to identify the optimal intervention set (Table 2).
Non-stationary DAG and SCM (NONSTAT.): DCBO outperforms all approaches in non-stationary settings where both the DAG and the SCM change overtime – see Fig. 3c. Indeed, DCBO can timely incorporate changes in the system via the dynamic causal prior construction while CBO, BO and ABO need to perform several interventions before accurately learning the new objective functions.

heading: 4.2 Real experiments
text: Real-World Economic data (ECON.) We use DCBO to minimize the unemployment rate Ut of a closed economy. We consider its causal relationships with economic growth (Gt), inflation rate (Rt) and fiscal policy (Tt)4. Inspired by the economic example in [17] we consider the DAG in
4The causality between economic variables is oversimplified in this example thus the results cannot be used to guide public policy and are only meant to showcase how DCBO can be used within a real application.

heading: 5 Conclusions
text: We consider the problem of finding a sequence of optimal interventions in a causal graph where causal temporal dependencies exist between variables. We propose the Dynamic Causal Bayesian Optimization (DCBO) algorithm which finds the optimal intervention at every time step by intervening in the system according to a causal acquisition function. Importantly, for each possible intervention we propose to use a surrogate model that incorporates information from previous interventions implemented in the system. This is constructed by exploiting theoretical results establishing the correlation structure among objective functions for two consecutive time steps as a function of the topology of the causal graph. We discuss the DCBO performance in a variety of setting characterized by different DAG properties and stationarity assumptions. Future work will focus on extending our theoretical results to more general DAG structures thus allowing for unobserved confounders and a changing DAG topology within each time step. In addition, we will work on combining the proposed framework with a causal discovery algorithm so as to account for uncertainty in the graph structure.

heading: Acknowledgements
text: This work was supported by the EPSRC grant EP/L016710/1, The Alan Turing Institute under EPSRC grant EP/N510129/1, the Defence and Security Programme at The Alan Turing Institute, funded by the UK Government and the Lloyds Register Foundation programme on Data Centric Engineering through the London Air Quality project. TD acknowledges support from UKRI Turing AI Fellowship (EP/V02678X/1).
5Data were downloaded from https://www.data.oecd.org/ [Accessed: 01/04/2021]. All details in the supplement.

Review summaries:
The paper proposes a model, Dynamic Causal Bayesian Optimization (DCBO), which builds upon existing models to find a sequence of interventions , optimizing the target variable at each time step in a causal system. It combines Causal Bayesian Optimization, which builds around casual information and finds the optimal intervention in a DAG without temporal evolution, with Dynamic Bayesian Network, which is used in time-series modelling and carries dependence assumptions that don’t imply causation. DCBO represents associations between features and the causality between input and outputs, both of which can change over time. The modeling therefore offers a novel approach for decision making in dynamic systems ie. find optimal intervention at time t. The paper compares the model’s results against algorithms which either account for temporal evolution or optimize the sequence of interventions of a system (as no model exists which does both). The algorithm is run on both interventional and observational data. Finally, synthetic experiments are performed, such as creating noisy variables, incorporating missing observational data, and creating non-stationary settings. 

This paper addresses the Dynamic Causal Global Optimization (DCGO) problem that aims to find both the optimal set of manipulative variables and the optimal intervention levels to the selected variables. This is a time-series optimization that determines the best way to intervene in the present based on the results of interventions made in the past, and is considered to be a formulation of a particularly important problem e.g. in the medical field. The authors make three assumptions about the DCGO problem: (1) the causal structure is invariant with respect to time,  (2) the function value of a variable Y_t is determined from the function value of its parent node, and the function value of the parent node can be expressed as the sum of the function value of the target node and that of the non-target node,  (3) the DAG does not contain any unobserved confounding factors.   Under these assumptions, it follows that the objective function of DCGO at time t can be expressed recursively using the results of interventions up to time t-1, and that the search space of DCGO is invariant with respect to time. Under the above problem setup, the proposed method DCBO is an EI-based Bayesian optimization algorithm based on the Dynamic Causal Gaussian Process model whose mean and covariance function are defined by the recursive representation of the objective function.    Experiments to validate the effectiveness of the proposed method have been conducted on both synthetic and real data. Experiments with synthetic data have been conducted in six situations: when DAGs and SCMs are stationary/non-stationary with respect to time, when the manipulative variables contain noise, when there are missing observations, when there is no causal relationship among manipulative variables,  and when there are multiple intervening variables. Two real data experiments were conducted: one on the problem of minimizing the unemployment rate in a closed economy (economics), and the other on the problem of finding an intervention that would reduce the phenomenon of a concentration of dead animals in a chemostat (biology).

This paper works with the Dynamic Causal Global Optimization problem. For a temporal structural causal model, which is a sequence of structural causal models with causality in time, meaning there are dependencies on the values of variables from previous step. The goal is to optimize the outcome variable by manipulating some observable variables (those are chosen by the model).  An algorithm called dynamic causal bayesian optimization is proposed to solve such a problem. The objective function $f_{s,t}$ is decomposed into the sum of a function depending on the optimized observed variables and another function depending on others (variables in the current time). Thus the optimization can be done step by step in time.  A dynamic causal GP model is introduced for practically implement the DCBO algorithm.  For experiments, both Synthetic and Real data are used, in several different structural causal models, to compare the behavior of DCBO, to that of CBO, ABO and BO.

In this paper, the authors extend the notion of casual Bayesian optimisation from a static to a dynamic scenario which, indeed, has many relevant applications. The paper theoretically analyses the proposed approach and run the algorithm in a set of real-world and synthetic data. This work is really interesting but I have some questions related to scalability and applicability that I hope the authors could help me understand.