Interventions, Where and How? Experimental Design for Causal Models at Scale

Panagiotis Tigas, Yashas Annadani, Andrew Jesson, Bernhard Schölkopf, Yarin Gal, Stefan Bauer

Introduction

What is the structure of the protein-signaling network derived from a single cell? How do different habits influence the presence of disease? Such questions refer to causal effects in complex systems governed by nonlinear, noisy processes. On most occasions, passive observation of such systems is insufficient to uncover the real cause-effect relationship and costly experimentation is required to disambiguate between competing hypotheses. As such, the design of experiments is of significant interest; an efficient experimentation protocol helps reduce the costs involved in experimentation while aiding the process of producing knowledge through the (closed-loop / policy-driven) scientific method (Fig. 1).

In the language of causality [Pearl, 2009], the causal relationships are represented qualitatively by a directed acyclic graph (DAG), where the nodes correspond to different variables of the system of study and the edges represent the flow of information between the variables. The abstraction of DAGs allows us to represent the space of possible explanations (hypotheses) for the observations at hand. Representing such hypotheses as Bayesian probabilities (beliefs) allows us to formalize the problem of the scientific method as one of Bayesian inference, where the goal is to estimate the posterior distribution $p(\text{DAGs}\mid\text{Observations})$ . A posterior distribution over the DAGs allows us to employ information-theoretic acquisition functions that guide experimentation towards the most informative variables for disambiguating between competing hypotheses. Such design procedures belong to the field of Bayesian Optimal Experimental Design [Lindley, 1956] for Causal Discovery (BOECD) [Tong and Koller, 2001, Murphy, 2001].

In the Bayesian Optimal Experimental Design (BOED) [Lindley, 1956] framework, one seeks the experiment that maximizes the expected information gain about some parameter(s) of interest. In causal discovery, an experiment takes the form of a causal intervention, and the parameters of interest are the Structural Causal Model (SCM) and its associated DAG.

An intervention in a causal model refers to the variable (or target) we manipulate and the value (or strength) at which we set the variable. Hence, the design space in the case of learning causal models is the set of all subsets of the intervention targets and the possibly countably infinite set of intervention values of the chosen targets. The intervention value encapsulates important semantics in many causal inference applications. For instance, in medical applications, an intervention can correspond to the administration of different drugs and the intervention value takes the form of a dosage level for each drug. Even though the appropriate choice of this value is crucial for identifying the underlying causal model, existing work on active causal discovery focuses exclusively on selecting the intervention target [Agrawal et al., 2019, Cho et al., 2016]. There, the intervention value is generally some arbitrary fixed value (like 0) which is suboptimal (see Fig. 2a). Hence, a holistic treatment of selecting the intervention value and the target in the general case of nonlinear causal models has been missing. We present a Bayesian experimental design method (CBED - pronounced “seabed”) to acquire optimal intervention targets and values by performing Bayesian optimization.

Additionally, some settings call for the selection of a batch of interventions. The problem of batched interventions is computationally expensive as it requires evaluating all possible combinations of interventions. We extend CBED to the batch setting and propose two different batching strategies for tractable, Bayes optimal acquisition of both intervention targets and values. The first strategy — Greedy-CBED — builds up the intervention set greedily. A greedy heuristic is still near-optimal due to submodularity properties of mutual information [Krause and Guestrin, 2012, Agrawal et al., 2019, Kirsch et al., 2019]. The second strategy — Soft-CBED — constructs a set of interventions by stochastic sampling from a finite set of candidates, thereby significantly increasing computational efficiency while recovering the DAG structure and the parameters of the SCM as fast as the greedy strategy. This strategy is well suited for resource-constrained settings.

Throughout this work, we make the following standard assumptions for causal discovery [Peters et al., 2017]:

There are no hidden confounders, and all the random variables of interest are observable.

There is a finite number of observational/ interventional samples available.

The structural causal model has nonlinear conditional expectations with additive Gaussian noise.

Each intervention is atomic and applied to a single target of the SCM.

Additionally, we assume that interventions are planned and executed in batches of size $\mathcal{B}$ , with a fixed budget of total interventions given by $\texttt{Number of Batches}\times\mathcal{B}$ . We also assume that the underlying graph is sparse, as is the case in all the real-world settings [Bengio et al., 2019, Schmidt et al., 2007]. Experimental design is preferable in sparse graph settings as the number of informative intervention targets and values would be significantly less compared to dense graphs. Many nodes corresponding to a sparse graph would have very less probability of having parent sets, and hence preforming experiments with a random policy is not maximally informative. Finally, we are interested in recovering the full graph $\mathbf{G}$ with a small number of batches. As with all causal inference tasks, the assumptions that we make above have to be carefully verified for the application of interest.

We show that our methods, Greedy-CBED and Soft-CBED, perform better than the state-of-the-art active causal discovery baselines in linear and nonlinear SCM settings. In addition, our approach achieves superior results in the real-world inspired nonlinear dataset, DREAM [Greenfield et al., 2010].

Background

A causal Bayesian network (CBN) is the pair $(\mathbf{g},P)$ such that for any $\mathbf{W}\subset\mathbf{V}$ ,

From the data generative mechanism point of view, the DAG $\mathbf{g}$ on $\mathbf{X_{V}}$ matches a set of structural equations:

where $f_{i}$ ’s are (potentially nonlinear) causal mechanisms that remain invariant when intervening on any variable $X_{j}\neq X_{i}$ . $\epsilon_{i}$ ’s are exogenous noise variables with arbitrary distribution that are mutually independent, i.e $\epsilon_{i}\mathchoice{\mathrel{\hbox to0.0pt{$ \displaystyle\perp $\hss}\mkern 2.0mu{\displaystyle\perp}}}{\mathrel{\hbox to0.0pt{$ \textstyle\perp $\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox to0.0pt{$ \scriptstyle\perp $\hss}\mkern 2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox to0.0pt{$ \scriptscriptstyle\perp $\hss}\mkern 2.0mu{\scriptscriptstyle\perp}}}\epsilon_{j}\forall i\neq j$ . (1) represents the conditional distributions in a Causal Bayesian Network and can additionally reveal the effect of interventions if the mechanisms are known [Peters et al., 2017, Pearl, 2009]. These equations together form the structural causal model (SCM), with an associated DAG $\mathbf{g}$ . Though the mechanisms $f$ can be nonparametric in the general case, we assume that there exists a parametric approximation to these mechanisms with parameters $\bm{\gamma}\in\Gamma$ . In the case of linear SCMs, $\bm{\gamma}$ corresponds to the weights of the edges in $E$ . In the nonlinear case, they could represent the parameters of a nonlinear function that parameterizes the mean of a Gaussian distribution.

A common form of (1) corresponds to Gaussian additive noise models (ANM)ANM’s can have noise variables that are non-Gaussian as well, but we restrict our exposition to the Gaussian case.:

An ANM is fully specified by a a DAG $\mathbf{g}$ , mechanisms, $f(\cdot;\bm{\gamma})=\begin{bmatrix}f_{1}(\cdot;\gamma_{1}),\dots,f_{d}(\cdot;\gamma_{d})\end{bmatrix}$ , parameterized by $\bm{\gamma}=\begin{bmatrix}\gamma_{1},\dots,\gamma_{d}\end{bmatrix}$ , and variances, $\sigma^{2}=\begin{bmatrix}\sigma^{2}_{1},\dots,\sigma^{2}_{d}\end{bmatrix}$ . For notational brevity, henceforth we denote $\bm{\theta}=(\bm{\gamma},\sigma^{2})$ and all the parameters of interest with $\bm{\phi}=(\mathbf{g},\bm{\theta})$ .

A common assumption in causal inference is that causal relations are known qualitatively and can be represented by a DAG. While this qualitative information can be obtained from domain knowledge in some scenarios, it’s infeasible in most applications. The goal of causal discovery is to recover the SCM and the associated DAG, given a dataset $\mathcal{D}$ . In general, without further assumptions about the nature of mechanisms $f$ (e.g., linear vs. nonlinear), the true SCM may not be identifiable [Peters et al., 2012] from observational data alone. This non-identifiability is because there could be multiple DAGs (and hence multiple factorizations of $P(\mathbf{X_{V}})$ ) which explain the data equally well. Such DAGs are said to be Markov Equivalent. Interventions can improve identifiability. In addition to identifiability issues, estimating the functional relationships between nodes using finite data is another source of uncertainty. Bayesian parameter estimation over the unknown SCM provides a principled way to quantify these uncertainties and obtain a posterior distribution over the SCM given observational data. An experimenter can then use the knowledge encoded by the posterior to design informative experiments that efficiently acquire interventional data to resolve unknown edge orientations and functional uncertainty.

The key challenge in performing Bayesian inference jointly over SCMs and DAGs is that the space of DAGs is discrete and superexponential in the number of variables [Peters et al., 2017]. However, recent techniques based on variational inference [Annadani et al., 2021, Lorch et al., 2021, Cundy et al., 2021] provide a tractable and scalable way of performing posterior inference of these parameters. Given a tractable distribution $q_{\psi}(\bm{\phi})$ which approximates the posterior $p(\bm{\phi}\mid\mathcal{D})$ , variational inference maximizes a lower bound on the (log-) evidence:

The key idea in these techniques is the way the variational family $\Psi$ for DAGs is parameterized. The variational family for the Variational Causal Network (VCN) method [Annadani et al., 2021] is an autoregressive Bernoulli distribution over the adjacency matrix. They further enforce the acyclicity constraint [Zheng et al., 2018] through the prior. BCD-Nets [Cundy et al., 2021] consider a distribution over node orderings through a Boltzmann distribution and perform inference with Gumbel-Sinkhorn [Mena et al., 2018] operator. DiBS [Lorch et al., 2021] consider latent variables over entries of adjacency matrix and perform inference over these latent variables using SVGD [Liu and Wang, 2016]. We demonstrate empirical results in BOECD using the DiBS model in this work because it is easily extendable to nonlinear SCMs.

Bayesian Optimal Experimental Design (BOED) [Lindley, 1956, Chaloner and Verdinelli, 1995] is an information theoretic approach to the problem of selecting the optimal experiment to estimate any parameter $\theta$ . For BOED, the utility of the experiment $\xi$ is the mutual information (MI) between the observation $\mathbf{y}$ and $\theta$ :

A common setting, called static, fixed or batch design, is to optimize $\mathcal{B}$ designs $\{\xi_{1},\dots,\xi_{\mathcal{B}}\}$ at the same time. The designs are then executed and the experimental outcomes are collected to update the model parameters in a Bayesian fashion.

Method

The above objective considers taking arg max over not just the discrete set of intervention targets $j\in\mathbf{V}$ , but also over the uncountable set of intervention values $v\subset\mathcal{X}_{j}$ . While the existing works in BOECD consider only the design of intervention targets to limit the complexity [Tong and Koller, 2001, Murphy, 2001, Agrawal et al., 2019], our approach tackles both the problems. We first outline the methodology for a single design and in Section 3.2 demonstrate how to extend this single design to a batch setting.

To maximize the objective in Equation 3, we need to (1) estimate MI for candidate interventions and (2) maximize the estimated MI by optimizing over the domain of intervention value for every candidate interventional target.

As mutual information is intractable, there are various ways to estimate it depending on whether we can sample from the posterior and whether the likelihood can be evaluated [Foster et al., 2020, Poole et al., 2019, Houlsby et al., 2011]. Since the models we consider allow both posterior sampling and likelihood evaluation, it suffices to obtain an estimator which requires only likelihood evaluation and Monte Carlo approximations of the expectations. To do so, we derive an estimator similar to Bayesian Active Learning by Disagreement (BALD) [Houlsby et al., 2011], which considers MI as a difference of conditional entropies over the outcomes $\mathbf{Y}$ :

As shown in (3), maximizing the objective is achieved not only by selecting the intervention target but also by setting the appropriate intervention value. Although optimizing the intervention target is tractable (discrete and finite number of nodes to select from), selecting the value to intervene is usually intractable since they are continuous. For any given target node $j$ , MI is a nonlinear function over $v\in\mathcal{X}_{j}$ (See Fig 2) and hence solving with gradient ascent techniques only yields a local maximum. Given that MI is expensive to evaluate, we treat MI for a given target node $j$ as a black-box function and obtain its maximum using Bayesian Optimization (BO) [Kushner, 1964, Zhilinskas, 1975, Močkus, 1975]. BO seeks to find the maximum of this function $max_{v\in\mathcal{X}_{j}}\mathcal{I}(\{(j,v)\})$ over the entire set $\mathcal{X}_{j}$ with as few evaluations as possible. See appendix E for details.

2 Batch Design

Computing the optimal solution $\mathcal{I}(\bm{\Xi}^{*})$ is computationally infeasible. However, as the conditional mutual information is submodular and non-decreasing (see Appendix B.4 for proof), we can derive a simple greedy algorithm (Algorithm 1) that can achieve at least a $(1-1/e)\approx 0.64$ approximation of the optimal solution [Krause and Guestrin, 2012, Nemhauser et al., 1978]. We denote this strategy as Greedy-CBED.

Although the greedy algorithm is tractable, it requires $O(\mathcal{B}d)$ instances of GP-UCB. Kirsch et al. show that a soft top-k selection strategy performs similarly to the greedy algorithm, reducing the computation requirements to $O(d)$ runs of GP-UCB. To achieve this, we construct a finite set of candidate intervention target-value pairs by keeping all the $T$ evaluations of GP-UCB for each node $j=\{1,\dots,d\}$ . Therefore, for $d$ nodes, our candidate set is comprised of $d\times T$ experiments. We score each experiment in this candidate set using the MI estimate. We then sample without replacement $\mathcal{B}$ times proportionally to the softmax of the MI scores (Algorithm 2). We denote this strategy as Soft-CBED.

3 Comparison with existing active causal discovery methods

We outline how our approach compares with two main existing active causal discovery methods.

The estimator of MI used in ABCD is based on weighted importance sampling. However, for the specific choice of the importance sampling weights used in ABCD, their MI estimator ends up with the same approximation as in our method (see Appendix B.5). Nevertheless, ABCD does not select intervention values but suboptimally sets them to a fixed value. In addition, our proposed Soft-CBED is a faster and more efficient batch strategy, especially when values also have to be acquired. From this perspective, our approach is an extension of ABCD with nonlinear assumptions, value acquisition, and a soft top-k batching strategy.

AIT is an F-score-based intervention target acquisition strategy. Although it is not a BOECD method, we prove here that it can be viewed as a Monte Carlo estimate of the approximation to MI when the outcomes $\mathbf{Y}$ are Gaussian. Nevertheless, AIT does not select intervention values like ABCD and does not have a batch strategy.

Let $\mathbf{Y}$ be a Gaussian random variable. Then the discrepancy score of Scherrer et al. is a Monte Carlo estimate of an approximation to mutual information (Eq. (3.1)). See Appendix B.6 for proof.

Related Work

Early efforts of using Bayesian Optimal Experimental Design for Causal Discovery (BOECD) can be found in the works of Murphy and Tong and Koller . However, these approaches deal with simple settings like limiting the graphs to topologically ordered structures, intervening sequentially, linear models, and discrete variables.

In Cho et al. and Ness et al. , BOECD was applied for learning biological networks structure. BOECD was also explored in Greenewald et al. under the assumption that undirected edges of the graph always forms a tree. More recently, ABCD framework [Agrawal et al., 2019] extended the work of Murphy and Tong and Koller in the setting where interventions can be applied in batches with continuous variables. To achieve this, they (approximately) solve the submodular problem of maximizing the batched mutual information between interventions (experiments), outcomes, and observational data, given a DAG. DAG hypotheses are sampled using DAG-bootstrap [Friedman et al., 2013]. Our work differs from ABCD in a few ways: we work with both linear and nonlinear SCMs by using state-of-the-art posterior models over DAGs [Lorch et al., 2021], we apply BO to select the value to intervene with, but we also prepare the batch using softBALD [Kirsch et al., 2021] which is significantly faster than the greedy approximation of ABCD method.

In von Kügelgen et al. the authors proposed the use of Gaussian Processes to model the posterior over DAGs and then use BO to identify the value to intervene with, however, this method was not shown to be scalable for larger than bivariate graphs since they rely on multi-dimensional Gaussian Processes for modeling the conditional distributions.

A new body of work has emerged in the field of differentiable causal discovery, where the problem of finding the structure, usually from observational data, is solved with gradient ascent and functional approximators, like neural networks [Zheng et al., 2018, Ke et al., 2019, Brouillard et al., 2020, Bengio et al., 2019]. In recent works [Cundy et al., 2021, Lorch et al., 2021, Annadani et al., 2021], the authors proposed a variational approximation of the posterior over the DAGs which allowed for modeling a distribution rather than a point estimate of the DAG that best explains the observational data $\mathcal{D}$ . Such work can be used to replace DAG-bootstrap [Friedman et al., 2013], allowing for the modeling of posterior distributions with greater support.

Besides the BOECD-based approaches, a few active causal learning works have been proposed [He and Geng, 2008, Gamella and Heinze-Deml, 2020, Scherrer et al., 2021, Shanmugam et al., 2015, Squires et al., 2020, Kocaoglu et al., 2017]. Active ICP [Gamella and Heinze-Deml, 2020] uses ICP [Peters et al., 2016] for causal learning while using an active policy to select the target, however, this work is not applicable in the setting where the full graph needs to be recovered. In Zhang et al. , the authors propose an active learning method to the problem of identifying the interventions that push a dynamical causal network towards a desired state. A few approaches tackle the problem of actively acquiring interventional data to orient edges of a skeletal graph [Shanmugam et al., 2015, Squires et al., 2020, Kocaoglu et al., 2017]. Closer to our proposal belongs AIT [Scherrer et al., 2021], which uses a neural network-based posterior model over the graphs but evaluates the F-score to select the interventions.

Experiments

We evaluate the performance of our method on synthetic and real-world causal experimental design problems and a range of baselines. We aim to investigate the following aspects empirically: (1) competitiveness of the overall proposed strategies of Greedy-CBED and Soft-CBED at scale ( $50$ nodes) on synthetic datasets; (2) performance of the value acquisition strategy based on GP-UCB; and (3) performance of the proposed approach on a real-world inspired dataset.

Random baseline acquires interventional targets at random.

active intervention targeting (AIT) [Scherrer et al., 2021] uses an f-score based acquisition strategy to select the intervention targets. See appendix B.6 for more details. Since the original proposed approach does not consider a batch setting, we introduce a variant that augments AIT with the proposed soft batching, as described in section 3.2.

These are the Monte Carlo estimates of MI, as described in section 3. CBED selects a single intervention (target and value) that maximizes the MI and this intervention is applied for the whole batch. In Greedy-CBED , the batch is built up in a greedy fashion selecting the target, value pairs one at a time (Algorithm 1). Soft-CBED is sampling (target, value) pairs proportionally to the MI scores to select a batch, as described in section 3.2 and Algorithm 2.

2 Value Selection Strategies

Fixed: This value selection strategy assumes setting the value of the intervention to a fixed value. In the experiments, we fixed the value to . Sample-Dist: This value selection strategy samples from the support of the observational data. GP-UCB: This strategy uses the proposed GP-UCB Bayesian optimization strategy to select the value that maximizes MI.

3 Tasks

In this setting, we generate Erdős–Rényi [Erdős and Rényi, 1959] (ER) and Scale-Free (SF) graphs [Barabási and Albert, 1999] of size 20 and 50. For linear SCMs, we sample the edge weights $\gamma$ uniformly at random. For the nonlinear SCM, we parameterize each variable to be a Gaussian whose mean is a nonlinear function of its parents. We model the nonlinear function with a neural network. In all settings, we set noise variance $\sigma^{2}=0.1$ . For both types of graphs, we set the expected number of edges per vertex to 1. We provide more details about the experiments in appendix D.1.

The DREAM family of benchmarks [Greenfield et al., 2010] are designed to evaluate causal discovery algorithms of the regulatory networks of a single cell. A set of ODEs and SDEs generates the dataset, simulating the reverse-engineered networks of single cells. We use GeneNetWeaver [Schaffter et al., 2011] to simulate the steady-state wind-type expression and single-gene knockouts. Refer to appendix D.2 for the exact settings.

4 Metrics

AUROC: The area under the receiver operating characteristic curve of the binary classification task of predicting the presence/ absence of all edges.

AUPRC: The area under the precision-recall curve of the binary classification task of predicting the presence/ absence of all edges.

Results

Next, we examine the importance of having a value selection strategy for active causal discovery. We use the MI estimator in Equation 3.1; moreover, we test the proposed GP-UCB with two heuristics - the fixed value strategy and sampling values from the support. As we can see in Figure 3(b), selecting the value using GP-UCB clearly benefits the causal discovery process. We expect this finding as the mutual information is not constant with respect to the intervened value. To make this point clear, we demonstrate in the appendix G the influence of the value in a simple two variables graph. In addition, we note that naively sampling from the support of the observed dataset performs worse than fixing the value to 0. We hypothesize that this is due to lower epistemic uncertainty in the high density regions of the support, hinting that these regions might be less informative.

In order to further understand how the soft batch strategy compares with other batch selection strategies, we compare the results of Soft-CBED with Greedy-CBED and CBED. We observe (Figure 3(c)) that Greedy-CBED and Soft-CBED give very similar results overall. While Greedy-CBED is optimal under certain conditions [Kirsch et al., 2019], Soft-CBED remains competitive and has the advantage that the batch can be selected in a one-shot manner. This is also evident from the runtime performance of both these batching strategies in Table 1. Both these batch selection strategies perform significantly better than selecting one intervention target/value pair, and executing them $\mathcal{B}$ times (CBED).

Summary and Conclusions

This paper studies the problem of efficiently selecting the Bayes optimal experiments to discover causal models. Our proposed framework simultaneously answers the questions of where and how to intervene in a batched setting. We present a Bayesian optimization strategy to acquire interventional targets and values. Further, we propose two different batching strategies: one based on greedy selection and the other based on soft top-k selection. The proposed methodology for selecting intervention target-value pairs in a batched setting provides superior performance over the state-of-the-art for causal models up to $50D$ variables. We validate this using synthetic datasets and using real-world inspired datasets of single-cell regulatory networks, showing the potential impact on areas like biology and other experimental sciences.

Acknowledgments and Disclosure of Funding

We would like to thank Nino Scherrer, Tom Rainforth, Desi R. Ivanova and all anonymous reviewers for sharing their valuable feedback and insights. Panagiotis Tigas is supported by the UK EPSRC CDT in Autonomous Intelligent Machines and Systems (grant reference EP/L015897/1). We are grateful for compute from the Berzelius Cluster and the Swedish National Supercomputer Centre.

References

Checklist

Do the main claims made in the abstract and introduction accurately reflect the paper’s contributions and scope? [Yes]

Did you describe the limitations of your work? [Yes] See end of Section 1. Our limitations arise from the fact that the proposed methodology and conclusion hold when the assumptions laid out are satisfied.

Did you discuss any potential negative societal impacts of your work? [Yes] Negative societal impact is discussed in appendix A and will be added to the extra page of the camera ready after the reviewing process.

Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yes]

If you are including theoretical results…

Did you state the full set of assumptions of all theoretical results? [Yes] The assumptions are laid out in Introduction as well as in Theorem 3.1.

Did you include complete proofs of all theoretical results? [Yes] The proofs are in Appendix.

Did you include the code, data, and instructions needed to reproduce the main experimental results (either in the supplemental material or as a URL)? [Yes]

Did you specify all the training details (e.g., data splits, hyperparameters, how they were chosen)? [Yes] See Appendix D.1, D.2 and H

Did you report error bars (e.g., with respect to the random seed after running experiments multiple times)? [Yes] See Figure 3 and 4.

Did you include the total amount of compute and the type of resources used (e.g., type of GPUs, internal cluster, or cloud provider)? [Yes] Please check appendix K for details.

If you are using existing assets (e.g., code, data, models) or curating/releasing new assets…

If your work uses existing assets, did you cite the creators? [Yes]

Did you mention the license of the assets? [Yes] Please check appendix L for details.

Did you include any new assets either in the supplemental material or as a URL? [No]

Did you discuss whether and how consent was obtained from people whose data you’re using/curating? [N/A] Data are simulated.

Did you discuss whether the data you are using/curating contains personally identifiable information or offensive content? [No] Data are simulated.

If you used crowdsourcing or conducted research with human subjects…

Did you include the full text of instructions given to participants and screenshots, if applicable? [N/A]

Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals, if applicable? [N/A]

Did you include the estimated hourly wage paid to participants and the total amount spent on participant compensation? [N/A]

Appendix A Potential negative societal impacts

Causal Experimental Design has the potential to impact several sectors; healthcare, biology, mechanical and material engineering, computational advertisement, to name a few. As any AI powered field, it can have negative societal impact when being used by malicious actors.

Appendix B Theoretical Results

In the following lemma, we derive the mutual information over outcomes given in (3.1).

B.2 Estimating the Mutual Information over Outcomes

where $\widehat{\mathbf{y}}_{i,k}\sim p(\mathbf{y}\mid\widehat{\bm{\phi}}_{l},\{(j,v)\})$ is one of $m$ samples from the density parameterised by the $i$ th of $c_{o}$ SCMs $\widehat{\bm{\phi}}_{i}\sim p(\bm{\phi}\mid\mathcal{D})$ augmented by intervention $\{(j,v)\}$ . The likelihood of the sample $\widehat{\mathbf{y}}_{i,k}$ is then evaluated under the parameterisation of the $l$ th of $c_{in}$ additional SCMs $\widehat{\bm{\phi}}_{l}\sim p(\bm{\phi}\mid\mathcal{D})$ augmented by intervention $\{(j,v)\}$ .

by Jensen’s inequality. We can then define an unbiased estimator of this lower bound.

where $\widehat{\mathbf{y}}_{i,k}\sim p(\mathbf{y}\mid\widehat{\bm{\phi}}_{i},\{(j,v)\})$ is one of $m$ samples from the density parameterised by the $i$ th of $c_{o}$ graphs $\widehat{\bm{\phi}}_{i}\sim p(\bm{\phi}\mid\mathcal{D})$ augmented by intervention $\{(j,v)\}$ .

B.3 Monte Carlo Estimator of the Batch Mutual Information

While Equation 3.1 pertains to MI for a single design, we present here the MI estimator for the batch design.

B.4 Mutual Information Submodularity and Monotonicity Proofs

$\mathbf{I}(Y;\omega\mid X)$ is submodular.

The proof follows the structure of [Kirsch et al., 2019, Appendix A].

$\mathbf{I}(Y;\omega\mid X)$ is non-decreasing.

B.5 Relation to MI Approximation in ABCD

Here we demonstrate that though ABCD [Agrawal et al., 2019] uses an importance weighted estimate of mutual information, for the specific choice of importance weights used in ABCD, the MI estimate turns out to be the same as the one used in this work.

We note that ABCD decomposes the MI as entropy over the SCM as opposed to the entropy over outcomes used in this work.

The mutual information in (3) can be written as:

The above equation cannot be estimated from samples of $q(\bm{\phi}\mid\mathcal{D})\approx p(\bm{\phi}\mid\mathcal{D})$ since the posterior of the SCM would change when the interventional outcome $\mathbf{y}$ is conditioned on. To address this problem, ABCD [Agrawal et al., 2019] proposes to use weighted importance sampling with weights $w=p(\mathbf{y}\mid\bm{\phi},\{(j,v)\},\mathcal{D})$ and use samples from $q(\bm{\phi}\mid\mathcal{D})$ .

The weighted importance sampling estimate of entropy over SCM (10) with weights $w(\bm{\phi})$ is given by

B.5.2 Entropy Over Outcomes.

We can instead consider an alternative factorisation of (3) which would not require importance sampling and also compute entropies in the lower dimensional space of experimental outcomes, as given in Equation 3.1.

The Monte Carlo estimate of entropy over outcomes (3.1) is given by

B.5.3 Relation between Approximations with Entropy over SCM and Entropy over Outcomes

We prove below that for specific choice of importance weights $w(\bm{\phi})\coloneqq=p(\mathbf{y}\mid\bm{\phi},\{(j,v)\},\mathcal{D})$ used in ABCD, the MI approximations due to the above two factorizations are the same.

Consider the importance weighted estimate of the above equation with weights $w(\bm{\phi})$ . We can rewrite $p(\bm{\phi}\mid\mathbf{y},\{(j,v)\},\mathcal{D})$ as:

Let $\{\widehat{\bm{\phi}}_{i}\sim p(\bm{\phi}\mid\mathcal{D})\}_{i=1}^{c_{o}}$ , using (15) in (14),

Plugging the above result back in (16b) and noticing that second term in the above equation cancels with first term in (16b), we get:

B.6 Information Theoretic Interpretation of Neural Causal Models with Active Interventions

Here we provide the proof for Theorem 3.1.

The discrepancy score for a target $j$ in AIT [Scherrer et al., 2021] is given by:

where $\widehat{\mathbf{y}}_{i,j,k}$ is the interventional sample from a hypothetical intervention on node $j$ on a graph $i$ sampled from the model. $\widehat{\mu}^{j}_{i}$ is the sample mean over samples $k$ in $\widehat{\mathbf{y}}_{i,j,k}$ and $\widehat{\mu}^{j}$ is the mean over all graphs and samples.

We restate Theorem 3.1 for the sake of completeness. See 3.1

Mutual Information over outcomes is given by

Since the bounds are in the opposite direction in (20e) and (20i), we cannot obtain a single common bound to MI but instead only a rough approximation given by (20i). We can now define a Monte Carlo estimate of the above approximation:

Appendix C Models

For optimizing DiBS [Lorch et al., 2021] we used RMSProp with learning rate 0.005. Additionally, per dataset we set the following hyperparameters:

C.2 DAG Bootstrap

The DAG bootstrap bootstraps observations and interventions to infer a different causal structure per bootstrap. We used GIES as the causal inference algorithm because of the adaptation of GES on interventional data as well. In our experiments, we used the pcalg R implementation https://github.com/cran/pcalg/blob/master/R/gies.R to discover 100 graphs. Each graph can be seen as a posterior sample from $p(\mathbf{G}\mid\mathcal{D})$ . For each of the sampled graphs $G_{i}$ we compute the appropriate $\theta_{\text{MLE}}$ under linear Gaussian assumption for the conditional distributions.

Appendix D Datasets and Experiment details

In the synthetic data experiments, we focus on two types of graphs. The Erdős-Rényi and Scale Free.

We used networkxhttps://networkx.org/documentation/networkx-1.10/reference/generated/networkx.generators.random_graphs.fast_gnp_random_graph.html and method fast_gnp_random_graph [Batagelj and Brandes, 2005] to generate graphs based on the Erdős-Rényi model. We set expected number of edges per vertex to 1.

We used igraphhttps://igraph.org/python/api/latest/igraph._igraph.GraphBase.html#Barabasi package to generate the graphs. We set the expected number of edges per vertex to 1.

For all the synthetic graph experiments, we used batch size of 10 and number of iterations of 10.

D.2 DREAM Experiments

For the DREAM experiments, we used GeneNetWeaver [Schaffter et al., 2011], a simulator of gene regulatory networks, based on stochastic differential equations. This simulator was used to generate data for Dialogue for Reverse Engineering Assessments and Methods (DREAM) [Sachs et al., 2005] competition with three network inference challenges (DREAM3, DREAM4 and DREAM5). We used the GeneNetWeaver v3.1https://github.com/tschaffter/genenetweaver.

Each experiment is parametrized as an xml file describing the network topology but also the crucial parameters of the stochastic differential equation that GeneNetWeaver simulates. In our experiments, we used Ecoli1, Ecoli2, Yeast1 and Yeast2 networks for 10 and 50 nodes.

Each experiment was initialized with 100 observational data. For the observational data, we used the steady state Steady state is considered the result of the simulation of the SDE for maximum 2000 steps. of wild-type experiments. For the interventional data, we used the steady-state of knock-out experiments. Each observational or interventional sample was conducted by running the simulator with a different seed per draw.

Appendix E Bayesian Optimisation

For the Gaussian Process, we used the following hyperparameters. Matern kernel, length scale $1.0$ , length scale bounds (lower= $1e-5$ , upper= $1e5$ ), Nu (smoothness of learned function) $2.5$ . Also added $1e-6$ to the diagonal of the kernel matrix.

Appendix F Related Work

Appendix G Mutual Information per value for two Variables graph

Appendix H metrics