Data is frequently used for various decision making which involves causal evidence seeking (or causal query) in observational data [25, 20, 38, 40]. For example, a biologist is planning for a biological experiment, and she has a collection of data from previous experiments by other researchers (those experiments do not need to have the same objective as hers) [6, 19]. She can plan the experiment on the data by “manipulating” the variable which is to be modified in the planned experiment, i.e. the treatment variable to see if the outcome is changed as expected. This type of query is useful in many real-world applications. When a public policy is in the design phase, existing data can be used to assess the policy to find out whether the policy would produce the desired outcome and what major factors would affect the outcome. In a marketing campaign, similar queries can be made to a dataset to help plan a successful campaign.

In this paper, we consider the problem that a user queries a dataset for a causal effect, i.e. how much the interested outcome Y will change due to a change of the treatment W. The answer depends on what other variables affect W and Y , called confounders or adjustment variables. Users do not know the adjustment variables in most cases. For example, when a biologist studies how a gene causes cancer, she may not know (or only partially know) the other genes which are causal factors of the cancer and regulate this gene under study. So, the adjustment variables are a part of the answer to seek too and the query answer will be a set of pairs each containing a set of possible adjustment variables (called an adjustment set in the paper) and the corresponding causal effect. To confirm which pair corresponds to the real-world mechanism, a user can design an experiment by physically controlling the variables in an adjustment set while manipulating the treatment variable, or pick up one or multiple returned adjustment sets consistent with her knowledge of the system as plausible answers.

From data, especially data with hidden variables [5, 27, 39], it may not be possible to find an exact adjustment set. The number of returned possible adjustment sets can be very large, which will make the query answer unusable. For example, a biologist has only a limited budget to try a few experiments, and she will not be able to try many possible adjustment sets. An advertisement company can try a limited number of A/B tests to make a promotion decision. Another constraint is the size of an adjustment set. A long adjustment set containing many variables will complicate the design of an experiment and make it hard to explain experimental results.

The causal query problem studied in this paper is related to adjustment variable search and causal effect estimation in data. However, existing works on causal effect estimation often have an adjustment set given [3, 14, 15] or do not explicitly indicate the adjustment set for the causal effect [16, 21, 31]. Hence they are not useful to solve our problem.

Graphical causal modelling is a principled approach to adjustment set search in data and causal effect estimation. Pearl has proposed the back-door criterion for identifying an adjustment set based on a given causal DAG (direct acyclic graph), and the do calculus for deriving identifiable causal effects from the given causal DAG and data [25]. The back-door criterion has become the main principle for identifying adjustment sets.

When there are no hidden variables or more precisely the causal sufficiency assumption is satisfied, i.e. all common causes are observed [32], we can learn from data a Markov equivalence class of causal DAGs for causal effect estimation. IDA is an algorithm for estimating causal effects directly from data [23]. As many equivalent causal DAGs can be learned from a dataset, IDA returns a multiset of causal effects with different adjustment sets blocking the back-door paths into the treatment variable W. Recently, Johansson et al. utilised deep learning to learn balanced representations for counterfactual inference [16]. Furthermore, Shalit et al. gave a meaningful and intuitive error bound to guide deep neural networks for estimating individual causal effects [31]. However, the algorithms based on

deep learning do not provide an explicit adjustment variable set.

When there are hidden variables, LV-IDA [24] extends IDA based to a MAG (Maximal Ancestral Graph) which represents a causal structure with hidden variables. LV-IDA also returns a multiset of causal effects. Furthermore, LV-IDA has low efficiency for high-dimensional or large data. Hyttinen et al. presented an ASP constraint solver for causal inference directly in data without assuming causal sufficiency [13], but it is impractical for high-dimensional data due to the fact that the query-based techniques need to perform a complete search for solutions. Entner et al. proposed a data-driven approach to estimating causal effect from data with hidden variables [10], but with impractically high time complexity. Louizos et al. developed a new neural network latent variables model to estimate population causal effects [21], but the model relies on a correct proxy variable estimation and does not provide an explicit adjustment set.

In this paper, we aim to develop a practical solution for efficient and accurate causal query in data without assuming causal suffi-ciency. The contributions of the work are:

• We have developed the theorem for finding a superset of adjustment variables for unbiased causal effect estimation based on the local causal structure learning in a MAG under a realistic pre-treatment assumption. The theorem supports efficient data-driven causal effect estimations without assuming causal sufficiency in data. To our best knowledge, they are the first theorem to support such a search.

• We have proposed an efficient and effective Data-drIven Causal Effect estimation (DICE) algorithm. DICE is faster than LV-IDA, the only practical data-driven causal effect estimation method with hidden variables in graphical causal modelling. Extensive experiments show that DICE works in datasets where LV-IDA fails, and produces better causal effect estimation than LV-IDA. The causal effects estimated by DICE are as accurate as those by the state-of-the-art methods using domain knowledge.

2.1 Notation and basic definitions

We use upper case letters to represent variables and bold-faced upper case letters to denote sets of variables. Let G = (V, E) be a graph, where is the set of nodes and E a set of edges. Two nodes are adjacent if there is an edge between them. The set of all the nodes adjacent to node V is denoted as Adj(V ). A path is a sequence of distinct nodes such that every pair of successive nodes are adjacent in G. A subpath of is denoted by , if there exists is a parent of and we use to denote the set of all parents of . If there exists is a spouse of we use to denote the set of all spouses of . A path from is a directed path if for any pair of adjacent nodes in the path where is a parent of is an ancestor of if there is a directed path from . The ancestor set of denoted as

A DAG is a directed graph without directed cycles (See an example of directed cycle in Figure 1 (a)). When a DAG satisfies the following Markov condition and faithfulness assumption, we can read dependency/independency in a distribution from the DAG.

Definition 2.1 (Markov condition [25]). Given a DAG G = (V, E) and P(V), the joint probability distribution of V, G satisfies the Markov condition if for probabilistically independent of all non-descendants of , given the parents of

Figure 1. (a) An example of a directed cycle, (b) An example of an almost directed cycle.

Based on the Markov condition, according to G P(V) can be factorized into: P

Definition 2.2 (Faithfulness [32]). A DAG G = (V, E) is faithful to P(V) iff every independence presenting in P(V) is entailed by G and fulfills the Markov condition. A distribution P(V) is faithful to a DAG G iff there exists a DAG G which is faithful to P(V).

A causal DAG is a DAG in which a node’s parents are interpreted as its direct causes. To learn a causal DAG from observational data, we also need to assume causal sufficiency.

Definition 2.3 (Causal sufficiency [32]). A dataset satisfies causal sufficiency if for every pair of variables , all their common causes are also in V.

When causal sufficiency is violated, ancestral graphs are used to represent data generating processes. An ancestral graph M is a mixed graph that does not contain directed cycles or almost directed cycles. A graph is called a mixed graph if it includes directed and bi-directed edges. When there exists a bi-directed edge M, there is not directed path directed cycle occurs when Figure 1 (b) is such an example. Let “” denote any allowed edge marks. For a path , a non-ending node contains collider path if each node excluding the ending nodes on is a collider.

Definition 2.4 (m-separation [27]). In an ancestral graph G, a path is said to be m-separated by a set of nodes Z (possibly does not contain any collider which is in Z, or (2). for any collider on the path and no descendant of . Two nodes are said to be m-connected by Z in are not m-separated by Z.

If are m-separated by , the information flow from is “blocked” by Z. In this paper, we call Z a block set of the ordered pair (

Definition 2.5 (Maximal ancestral graph (MAG)). An ancestral graph M = (V, E) is called a maximal ancestral graph if every pair of non-adjacent nodes can be m-separated by a set

We present an example of causal MAG in Figure 2 (b). Figure 2 (a) is its corresponding causal DAG with the hidden variables shown. In a MAG M, if a path between is a directed path, it is called a causal path from . A non-directed path between non-causal path. In Figure 2(b), path causal path, and other paths between W and Y are non-causal paths.

Definition 2.6 (Visibility [41]). Given a MAG M, a directed edge is visible if there is a node not adjacent to , such that either there is an edge between that is into , or there is a collider path between that is into and every node in this path is a parent of . Otherwise, is said to be invisible.

Figure 2. (a) A DAG with two hidden variables . (b) The corresponding MAG of the DAG.

Figure 3. Two possible configurations of the visible edge for

Figure 3 shows two different graphical configurations where the edge is visible. In a given MAG M, a visible edge implies that there are no hidden variables between Otherwise, the edge may include hidden variables.

2.2 Adjustment set for causal effect estimation

Let W be a binary treatment variable, Y the outcome variable. The causal effect of W on Y , CE(W, Y ) is the change of Y due to a change of W. When we use the do operator [25] to represent a manipulation of a variable (e.g. set W to 1) and denote W = 1 and respectively, CE(W, Y ) is defined as:

where E is the Expectation. CE(W, Y ) can be found in an experimental setting where W is physically set as 0 and 1, respectively. Let us assume that an underlying mechanism dictates the causal relationships among variables, and CE(W, Y ) is determined by the mechanism.

When estimating the causal effect of W on Y, the effect of other variables needs to be eliminated or adjusted to obtain an unbiased estimation of the causal effect. In graphical causal modelling, the back-door criterion is the main principle for identifying a set of adjustment variables, denoted as Z in this paper.

Definition 2.7 (Back-door criterion in DAG). A set of variables Z satisfies the back-door criterion relative to (W, Y ) in a DAG G if (1) Z does not contain descendants of W; (2) Z blocks every back-door path between W and Y (i.e. the paths with an arrow pointing to W).

Given an adjustment set Z, the causal effect CE(W, Y ) can be estimated unbiasedly as follows:

When the assumption of causal sufficiency does not hold, the generalised back-door criterion in a MAG can be used to search for an appropriate adjustment set Z for estimating CE(W, Y ) [22].

Definition 2.8 (Back-door path in MAG). For the ordered pair (W, Y ) in a MAG, a path from W to Y is a back-door path if it does not have a visible edge out of W.

Definition 2.9 (Generalised back-door criterion for MAG). Given a MAG satisfies the generalised back-door criterion w.r.t. (W, Y ) in M if (1) Z does not contain descendants of W; (2) Z blocks every back-door path between W and Y .

If a set Z in M satisfies the generalised back-door criterion, relative to (W, Y ), then the unbiased estimation of CE(W, Y ) can be achieved by Eq.(2).

In some MAGs, the adjustment sets could not be found due to the causal ambiguities, and this can be judged based on the amenability of a MAG [35] as defined below.

Definition 2.10 (Amenable MAG w.r.t. (W, Y )). Let W and Y be two nodes in a MAG M = (V, E). The MAG is adjustment amenable w.r.t. is visible.

If the MAG M is not adjustment amenable w.r.t. (W, Y ), then no adjustment set Z in V \ {W, Y } can be found [35].

3.1 The Theorem

Let D be a dataset containing a binary treatment variable W, an outcome variable Y which is binary or numerical, and X, a set of all other variables of any type. We assume that the dataset is generated from an underlying causal MAG M which is adjustment amenable to the ordered pair (W, Y ). We also assume that all the variables in X are measured before applying the treatment and observing Y , indicating that variables in X are all non-descendants of W or Y in M. The pretreatment variable assumption is realistic as it reflects what normally happens in practice, that is, all the other variables are often measured before the treatment is applied and the outcome is observed. This assumption is commonly used in [9, 12, 14, 36].

In the proposed Theorem 1 below, denotes the manipulated M from which has been removed, and is still a MAG by the closure property of a MAG [27].

Theorem 1 (Adjustment set discovery via local search). Given a causal MAG M which is adjustment amenable to the ordered pair (W, Y ), there exists at least an adjustment set such that Z is a subset of shorthand of

Proof. The requirement that M is adjustment amenable relative to the ordered pair (W, Y ) is to ensure that the edge is visible and the causal effect of W on Y is identifiability. There are two types of paths between W and Y , causal paths and back-door paths. In our problem setting, pretreatment variables X are non-descendants of W and Y , and hence does not contain any causal path between W and Y .

To prove that an adjustment set Z in the MAG M is a subset of , we will show that the set Pa(Y ) \ {W} blocks all back-door paths from be any back-door path between W and Y , and let X denote the node adjacent to . Note that X is not W, for the assumption, is visible and so is not a back-door path. There are two cases as follows.

collider on and is in the set is blocked by Pa(Y ) \ {W}.

2. The edge between . Then it must be bi-directed edge: , for otherwise X would be a descendant of Y and hence a descendant of W, contradicting the assumption. This means that X is not an ancestor of Y and hence also not an ancestor of W in the MAG , there must be a collider. Let C be the collider closest to must be a descendant of X (possibly X itself), which means that C does not have a descendant in is blocked by Pa(Y ) \ {W}.

Therefore, given a MAG M, Pa(Y ) \ {W}, which is a subset of , is an adjustment set for the ordered pair (W, Y ).

3.2 Causal effect estimation

Based on the above Theorem developed, we can firstly employ a local causal structure discovery algorithm to obtain data, then our answer to a causal query is an Adjustment Set-Causal Effect Table, or ASCET, where each row contains a candidate adjustment set (a subset of ) and the causal effect estimated corresponding to the candidate adjustment set.

Example. Suppose , Table 1 shows an example ASCET where 1 denotes a variable is in a candidate adjustment set and 0 otherwise.

3.3 Removing insignificant adjustment variables

We reduce the size of ASCET by removing variables that do not sig-nificantly affect causal effect estimation. To test if a variable significantly affects causal effect estimation, we compare the cause effect estimated by including the variable with the cause effect estimated by excluding the variable in all cases. If the average difference is smaller than a threshold, the variable is insignificant. More precisely, we have the following definition.

Table 1. An example ASCET.

1 1 0 0.4

1 0 1 0.1

Table 2. The final ASCET.

Definition 3.1. For each be a candidate adjustment set containing sensitivity of X is defined as:

where is the size of the power set of and are average causal effects of W on Y with the adjustment sets respectively.

The sensitivity of variable X reflects its impact on the estimated causal effect of W on Y . When X has a low sensitivity, no significant error will be introduced if we exclude it from is removed from , and the size of ASCET is reduced.

Example. We use the example in Table 1 to illustrate how we reduce the size of ASCET. We use as an example, as the size of the power set of is 4, then

In the same way, Hence, has a low sensitivity such that we can exclude it from

and the corresponding rows in ASCET can be removed. The final result is presented as Table 2.

3.4 The proposed DICE algorithm

Algorithm 1 presents our proposed algorithm for Data-drIven Causal Effect estimation without causal sufficiency (DICE). DICE contains two parts. Part 1 (lines 1 to 13) is for finding a local causal structure learning algorithm (such as PC-Select [4] or HITON-PC [1]) is used to find the parents of W and Y . In this work, we use PC-Select, given its high accuracy, PC-Select is a local version of PC [32], an algorithm for learning a global causal structure. Then we check whether Adj(W) = Adj(W) \ {Y } is empty to determine if the causal effect of W on Y is identifiable or not from data. If Adj(W) is empty, no adjustment set can be found, and DICE returns an empty ASCET and terminates.

Part 2 of Algorithm 1 (lines 7 to 16) firstly estimates each possible causal effect when adjusting each subset of use Propensity Score Matching (PSM) [28] to estimate causal effects where propensity score is calculated by glm and matching is performed by Match in the R packages stats and Matching [30], respectively. The adjustment sets and the corresponding causal effects estimated are added to the ASCET. Next, DICE calculates the sensitivity of each variable in , and if a variable’s sensitivity is below the given threshold, DICE removes all the adjustment sets containing the variable and the corresponding causal effects from the ASCET.

Time-complexity analysis. PC-Select and PSM contribute the most to the time complexity of DICE. In the worst case, the complexity of PC-Select is [4], where p is the number of variables and n the number of samples. In our problem, the worst case is when all variables of X are causes of W and/or Y , which, however rarely occurs. In most cases, PC-Select can handle thousands of variables [4]. The time complexity of PSM is experiment in Section 4 shows that DICE scales with both p and n well.

Figure 4. The Precision, Recall and F-score of FCI & RFCI & PC-Select algorithms on two groups of synthetic datasets w.r.t. the number of variables and samples. The qualities of parent discovery by the three algorithms are very consistent. but the time efficiencies of the algorithms are quite different as shown in Figure 6

4.1 The quality of local structure learning

Before evaluating DICE, we assess the quality of local structure learning algorithm, i.e. PC-Select used in our paper versus global structure learning algorithms. The synthetic datasets are generated according to the underlying causal DAG in Fig 2 (a) with 8 pretreatment variables, 2 hidden variables, W and Y by following the procedure in [11]. The underlying causal MAG is presented in Fig 2 (b), i.e. are removed from these generated datasets.

We generate two groups of datasets to evaluate the impact of sample size and the number of variables on the quality of the local structure learning algorithm, respectively. The first group of datasets include 10 datasets , each containing 5K samples for 10 variables generated from Fig 2(b), plus 0, 10, 20, . . . , 90 random variables, respectively. The second group of datasets include 20 datasets with 10 fixed variables of the MAG, but with varying sample sizes of 5K, 10K, . . . , 95K and 100K. These datasets with hidden variables are used to evaluate the quality of the local structure learning algorithm.

We use two global structure learning algorithms, Fast causal inference (FCI) [32] and Really FCI (RFCI) [8] to learn from each of the generated datasets a PAG (Partial Ancestral Graph) which represents a Markov equivalence class of MAGs encoding the same dependence/independence relations in the data. Then we extract Adj(W) and Adj(Y ) from the learned PAG. For PC-Select, we apply it twice to a dataset, one having W as the target and the other having Y as the target to learn Adj(W) and Adj(Y ). The implementations of FCI, RFCI and PC-Select are from the R package pcalg [17], with the default parameter settings and a significance level (

Results. We perform the experiments 10 times on each of the 30 datasets. We draw the mean results (precision, recall and F-score) of the three algorithms in Figure 4. The local structure learning algorithm PC-Selects achieves similar performance as the global structure learning algorithm in Adj(W) and Adj(Y ), indicating that ad-

justment set discovery in data by local search is reliable.

4.2 Evaluating DICE with real-world datasets

We evaluate the effectiveness of DICE on three real-world datasets: Jobs training (Jobs for short) [18], IHDP [12] and Twins [2]. A brief description of the datasets is shown in Table 3.

Table 3. Summary of the real-world datasets.

Table 4. Results on Jobs. Bias (%) is the relative error comparing to the ground-truth causal effect. The estimation with the lowest bias is highlighted in each category.

DICE 1745.5 97.0% The most probable estimate PSM -947.6 201.0% CBPS 423.3 52.0% CFR 742.0 16.0% LASSO -475.3 153.6% BART -245.2 127.7% TMLE -1901.4 314.6% CF -4438.4 600.9% The best of 30 forests

Table 5. Results on IHDP w.r.t. CE and Bias.

DICE 4.48 2.5% The most probable estimate PSM 4.95 13.2% CBPS 4.56 4.3% CFR 4.86 11.0% LASSO 5.06 15.7% BART 4.69 7.3% TMLE 4.86 11.2% CF 4.25 2.9% The best of 30 forests

Jobs consists of the original LaLonde dataset (297 treated samples and 425 control samples) [18] and the Panel Study of Income Dynamics (PSID) observational group (2490 control samples) [14]. Each sample has 9 variables, including 7 pretreatment variables (which are age in years, schooling in years, indicators for black and Hispanic, marriage status, school degree, previous earnings in 1974 and 1975, and whether the 1974 earnings variable is missing), employment status with/without job training as treatment variable, and 1978 earning as the outcome variable. Because the dataset contains records of people taking part in the training only, as in [18], we estimate Average Treatment Effect on Treated (ATT) for DICE and all comparisons, against the ground truth ATT which is $886 [14].

IHDP is related to the Infant Health and Development Program (IHDP) on low-birth-weight premature infants [12]. IHDP has 24 pretreatment variables and we follow the method in [12] to generate the simulated outcome with the ground-truth of 4.38 for the causal effect (CE) by the

Table 6. Results on Twins w.r.t. CE and Bias..

The Twins dataset consists of samples from twin births in the USA between 1989 and 1991 with the birth weight less than 2000g [2]. We eliminate records with missing values and have 4821 twin pairs left. The treatment variable is birth weight: W=1 for a baby who is heavier in a twin; W=0 otherwise. The mortality after one year is the true outcome for each twin, and the ground-truth causal effect is -0.025. To simulate an observational study, according to [21], we use the following Bernoulli distribution to randomly select one of the two twins as the observation and hide the other: denotes the pre-treatment variables of the sample i, x denotes the samples set, and

There are two major objectives for the experiments.

Firstly, we will show that the true adjustment set can be found as shown in Theorem 1, and therefore the smallest bias of DICE will be very small. We compare DICE with a data-driven causal effect estimation method, LV-IDA [24]4 which uses FCI to learn a PAG and searches for an adjustment set in each MAG enumerated from the PAG. LV-IDA returns a set of causal effects, one for the adjustment set discovered from each MAG. In the experiment, we report the best estimates for both methods to show the quality of adjustment sets found.

Secondly, we will show that a causal effect estimated by DICE is comparative with the best estimates by methods using domain knowledge. Some statistical and machine learning methods are available for causal effect estimation, with the assumptions of known covariate set and unfoundedness. The typical methods include PSM, propensity score matching with logistic regression [28]; CBPS, covariate balancing propensity score [14]; CFRNET, a deep learning framework for counterfactual regression with a theoretical error bound (named as CFR) [31]; LASSO, Linear regression with the regularization -norm to predict the factual outcome [33]; BART, Bayesian additive regression tree [7], a non-linear model which has been applied for counterfactual inference [12]; TMLE, targeted maximum likelihood estimation, a doubly robust method [34]; CausalForest, Random forest regression to estimate causal effect [37] (CF for short). In the experiments, we compare DICE with all the above mentioned methods.

The three datasets will favour the above methods since the unconfoundedness assumption is satisfied. Since the datasets were obtained from well designed observational studies and the covariates were chosen by domain experts.

DICE returns ASCET, the set of adjustment set and causal effect pairs. When reporting causal effects estimated by DICE, we use the most probable value in the ASCET after removing insignificant variables with the sensitivity threshold of . We group the causal effect in the ASCET by using the bin size of (max(Y )/100). The average causal effect in the most frequent bin is used as an estimated causal effect.

The parameter setting of DICE is that the Match function of estimate is set to “ATT” for Jobs and “ATE” for IHDP and Twins as in the prior work in [12, 31]. Other parameters are set as the default. For CF, the parameter num.trees is set from 10 to 300 with the increment of 10, and the best result is reported. For PSM and CBPS, all variables are included in the adjustment variable set.

The experimental results are shown in Tables 4, 5, 6, and Figure 5. From the results, we have the following conclusions.

Firstly, DICE finds the correct adjustment set through local search. Comparing with the other methods, the lowest biases of DICE are consistently the smallest in the three datasets, and are significantly smaller than the lowest biases of LV-IDA, which searchers adjustment sets by global search.

Secondly, causal effects estimated by DICE are comparative with those by other methods which make use of domain knowledge. For the IHDP dataset, DICE achieves the smallest bias. For the Jobs dataset, DICE is ranked the 3rd based on smallest biases. For the Twin dataset, the largest absolute bias is only 0.015, and hence we consider that the estimates by all methods are similar and good.

To show that those comparison methods do not work in a dataset when the unconfoundedness assumption is not satisfied, we have generated datasets with M-bias [26]. The performance of some methods deteriorates greatly. Implementations of other methods do not work with datasets with more than 200K records. We do not show the results in this paper due to the space limitation.

4.3 Efficiency evaluation of DICE

As PC-Select is a major contributor to the complexity of DICE, in this section, we firstly evaluate the efficiency of PC-Select, and then the overall efficiency of DICE. The computations were performed on a PC with 2.6 GHz Intel Core i7 and 16GB of memory. Efficiency evaluation of the structure learning algorithms.. We record the running time of the structure learning algorithms on the synthetic datasets introduced in Section 4.1 and draw the mean running time in Figure 6. PC-Select is faster than FCI and RFCI, and comparing to FCI and RFCI, PC-Select scales with the number of variables and samples very well.

Experiments on ten standard benchmark BNs. To evaluate the efficiency of DICE, we use the datasets generated from the ten benchmark BNs from the BN repository5: SACHS, CHILD, INSURANCE, WATER, ALARM, BARLEY, HAILFINDER, HEPAR II, WIN95 and ANDES. With each BN, we choose a variable without child nodes as the outcome variable and then select one of its parents without child nodes (except the outcome) as the treatment variable. We generate ten synthetic datasets from the ten BNs with 5000 samples each using the R package bnlearn [29]. Then we hide 5% variables (only one variable is hidden for the small BNs SACHS and CHILD) which lie on a non-causal path between the treatment and outcome variables. The same parameter settings are used for DICE and LV-IDA as in Subsection 4.2.

Results. As shown in Figure 7, DICE is faster than LV-IDA across all the 10 datasets. Especially, LV-IDA did not return results for the 5 datasets generated using the larger BNs within two hours, while DICE completed in seconds (at most 50+ seconds) for all datasets.

Figure 5. Histograms of causal effects in ASCETs for Jobs, IHDP and Twins. The red line denotes the ground truth and the bin size is (max(Y )/100).

Figure 6. The runtime on two groups of synthetic datasets w.r.t. the number of variables and samples.

Figure 7. The runtime on ten BNs. Note: LV-IDA did not return results in two hours on BARLEY, HAILFINDER, HEPAR II, WIN95 and ANDES.

Causal query in data is an important means for evidence-based decision making, without relying on or being restricted by domain knowledge. However, its widespread applications are hindered by the low efficiency and unsatisfactory accuracy of existing methods. In this paper, we have developed a theorem which supports the utilisation of efficient local causal structure discovery for causal query and assures the correctness of the result of causal query obtained with local search. Based on the proposed theorem, we have developed DICE to query data for causal effects and the confounders potentially impacting the causal effects. The results returned by DICE provide decision makers not only valuable information on the effect of changing one variable in their systems, but also the awareness of the other factors which could influence the effect. Experiments with synthetic and real-world data have demonstrated the effectiveness and efficiency of DICE. It has been shown that DICE produces better causal effect estimation than LV-IDA, and works on datasets where LV-IDA fails. The causal effects estimated by DICE are as good as the state-of-the-art methods which use domain knowledge. In future, we will evaluate DICE and other methods on purely observational datasets which may include M-bias, and apply it to real-world biological applications.

We would like to thank the anonymous reviewers for their constructive suggestions. We would also like to thank Assoc. Prof. Jiji Zhang from Lingnan University for many helpful comments. This research is partially funded by ARC DP170101306 and the National Science Foundation of China (under grant 61876206). The first author is supported by China Scholarship Council.

[1] Constantin F Aliferis, Alexander Statnikov, Ioannis Tsamardinos, Subramani Mani, and Xenofon D Koutsoukos, ‘Local causal and markov blanket induction for causal discovery and feature selection for classifi-cation part i: Algorithms and empirical evaluation’, Journal of Machine Learning Research, 11(Jan), 171–234, (2010).

[2] Douglas Almond, Kenneth Y Chay, and David S Lee, ‘The costs of low birth weight’, The Quarterly Journal of Economics, 120(3), 1031–1083, (2005).

[3] Susan Athey, Guido W Imbens, and Stefan Wager, ‘Approximate residual balancing: debiased inference of average treatment effects in high dimensions’, Journal of the Royal Statistical Society: Series B (Statistical Methodology), 80(4), 597–623, (2018).

[4] Peter B¨uhlmann, Markus Kalisch, and Marloes H Maathuis, ‘Variable selection in high-dimensional linear models: partially faithful distributions and the pc-simple algorithm’, Biometrika, 97(2), 261–278, (2010).

[5] Ruichu Cai, Feng Xie, Clark Glymour, Zhifeng Hao, and Kun Zhang, ‘Triad constraints for learning causal structure of latent variables’, in NeurIPS, pp. 12863–12872, (2019).

[6] Ruichu Cai, Zhenjie Zhang, and Zhifeng Hao, ‘Causal gene identifi-cation using combinatorial v-structure search’, Neural Networks, 43, 63–71, (2013).

[7] Hugh A Chipman, Edward I George, Robert E McCulloch, et al., ‘Bart: Bayesian additive regression trees’, The Annals of Applied Statistics, 4(1), 266–298, (2010).

[8] Diego Colombo, Marloes H Maathuis, Markus Kalisch, and Thomas S Richardson, ‘Learning high-dimensional directed acyclic graphs with latent and selection variables’, The Annals of Statistics, 294–321, (2012).

[9] Xavier De Luna, Ingeborg Waernbaum, and Thomas S Richardson, ‘Covariate selection for the nonparametric estimation of an average treatment effect’, Biometrika, 98(4), 861–875, (2011).

[10] Doris Entner, Patrik Hoyer, and Peter Spirtes, ‘Data-driven covariate selection for nonparametric estimation of causal effects’, in Artificial Intelligence and Statistics, pp. 256–264, (2013).

[11] Jenny H¨aggstr¨om, ‘Data-driven confounder selection via markov and bayesian networks’, Biometrics, 74(2), 389–398, (2018).

[12] Jennifer L Hill, ‘Bayesian nonparametric modeling for causal inference’, J. Comput. Graph. Statist., 20(1), 217–240, (2011).

[13] Antti Hyttinen, Frederick Eberhardt, and Matti J¨arvisalo, ‘Do-calculus when the true graph is unknown.’, in UAI, pp. 395–404, (2015).

[14] Kosuke Imai and Marc Ratkovic, ‘Covariate balancing propensity score’, Journal of the Royal Statistical Society: Series B (Statistical Methodology), 76(1), 243–263, (2014).

[15] Guido W Imbens and Donald B Rubin, Causal inference in statistics, social, and biomedical sciences, Cambridge University Press, 2015.

[16] Fredrik Johansson, Uri Shalit, and David Sontag, ‘Learning representations for counterfactual inference’, in ICML, pp. 3020–3029, (2016).

[17] Markus Kalisch, Martin M¨achler, Diego Colombo, Marloes H Maathuis, Peter B¨uhlmann, et al., ‘Causal inference using graphical models with the r package pcalg’, Journal of Statistical Software, 47(11), 1–26, (2012).

[18] Robert J LaLonde, ‘Evaluating the econometric evaluations of training programs with experimental data’, The American economic review, 604–620, (1986).

[19] Thuc Duy Le, Lin Liu, Anna Tsykin, Gregory J Goodall, Bing Liu, Bing-Yu Sun, and Jiuyong Li, ‘Inferring microrna–mrna causal regulatory relationships from expression data’, Bioinformatics, 29(6), 765– 771, (2013).

[20] Jiuyong Li, Thuc Duy Le, Lin Liu, Jixue Liu, Zhou Jin, Bingyu Sun, and Saisai Ma, ‘From observational studies to causal rule mining’, ACM TIST, 7(2), 14, (2016).

[21] Christos Louizos, Uri Shalit, Joris M Mooij, David Sontag, Richard Zemel, and Max Welling, ‘Causal effect inference with deep latentvariable models’, in NeurIPS, pp. 6446–6456, (2017).

[22] Marloes H Maathuis, Diego Colombo, et al., ‘A generalized back-door criterion’, The Annals of Statistics, 43(3), 1060–1088, (2015).

[23] Marloes H Maathuis, Markus Kalisch, Peter B¨uhlmann, et al., ‘Estimating high-dimensional intervention effects from observational data’, The Annals of Statistics, 37(6A), 3133–3164, (2009).

[24] Daniel Malinsky and Peter Spirtes, ‘Estimating causal effects with ancestral graph markov models’, in Conference on Probabilistic Graphical Models, pp. 299–309, (2016).

[25] Judea Pearl, Causality, Cambridge university press, 2009.

[26] Judea Pearl, ‘Myth, confusion, and science in causal analysis’, Tech. Rep. R-348, (2009). Los Angeles, CA: University of California.

[27] Thomas Richardson, Peter Spirtes, et al., ‘Ancestral graph markov models’, The Annals of Statistics, 30(4), 962–1030, (2002).

[28] Paul R Rosenbaum and Donald B Rubin, ‘The central role of the propensity score in observational studies for causal effects’, Biometrika, 70(1), 41–55, (1983).

[29] Marco Scutari, ‘Learning bayesian networks with the bnlearn r package’, arXiv preprint arXiv:0908.3817, (2009).

[30] Jasjeet S Sekhon, ‘Multivariate and propensity score matching software with automated balance optimization: the matching package for r’, Journal of Statistical Software, Forthcoming, (2008).

[31] Uri Shalit, Fredrik D Johansson, and David Sontag, ‘Estimating individual treatment effect: generalization bounds and algorithms’, in Proceedings of the 34th ICML-Volume 70, pp. 3076–3085, (2017).

[32] Peter Spirtes, Clark N Glymour, Richard Scheines, David Heckerman, Christopher Meek, Gregory Cooper, and Thomas Richardson, Causation, prediction, and search, MIT press, 2000.

[33] Robert Tibshirani, ‘Regression shrinkage and selection via the lasso’, Journal of the Royal Statistical Society: Series B (Methodological), 58(1), 267–288, (1996).

[34] Mark J Van Der Laan and Daniel Rubin, ‘Targeted maximum likelihood learning’, The International Journal of Biostatistics, 2(1), (2006).

[35] Benito van der Zander, Maciej Liskiewicz, and Johannes Textor, ‘Constructing separators and adjustment sets in ancestral graphs.’, in CI@ UAI, pp. 11–24, (2014).

[36] Tyler J VanderWeele and Ilya Shpitser, ‘A new criterion for confounder selection’, Biometrics, 67(4), 1406–1413, (2011).

[37] Stefan Wager and Susan Athey, ‘Estimation and inference of heterogeneous treatment effects using random forests’, Journal of the American Statistical Association, 113(523), 1228–1242, (2018).

[38] Kui Yu, Lin Liu, and Jiuyong Li, ‘Learning markov blankets from multiple interventional data sets’, IEEE Trans. Neural Netw. Learn. Syst, (2019).

[39] Kui Yu, Lin Liu, Jiuyong Li, and Huanhuan Chen, ‘Mining markov

blankets without causal sufficiency’, IEEE Trans. Neural Netw. Learn. Syst, 29(12), 6333–6347, (2018).

[40] Kui Yu, Lin Liu, Jiuyong Li, Wei Ding, and Thuc Duy Le, ‘Multi-source causal feature selection’, IEEE Transactions on Pattern Analysis and Machine Intelligence, (2019).

[41] Jiji Zhang, ‘Causal reasoning with ancestral graphs’, Journal of Machine Learning Research, 9(Jul), 1437–1474, (2008).