Information-Theoretic Approximation to Causal Models

1 Introduction

Detecting causal relationships from data is a significant issue in many disciplines. The understanding of causal relations between variables can help to understand how a system behaves under intervention, stabilize future predictions, and has many other important implications. Identifying causal links (causal discovery) from observed data alone is only possible with further assumptions and/or additional data. Despite the various causal discovery methods available, the problem of finding the causal structure between two random variables remains notoriously hard. In this paper, we use additional data and assume a very natural principle to solve that task. Our work is based on a mathematical framework proposed by Pearl Pea09 , that formalizes causality and causal relations. It introduces causal models that represent an (unknown) underlying data generation mechanism responsible for the distribution of the sampled data PJS17 . We include sampled data from situations (environments) where interventions took place together with samples from pure observations. Recent developments in that direction revealed promising results PBM16 ; HPM18 , but often these methods are conservative, leading to situations where no direction is preferred. This paper focuses on the bivariate discrete case and is based on a natural and weak principle. The principle of independent mechanism assumes that the data generating mechanism is independent of the data that is feed into such a mechanism. From this principle, we derive an invariance relation that states that it does not matter if we observe an effect due to an observation of its cause or due to an intervention on its cause. Distributions that are generated by an underlying causal model fulfill these invariance relations. If $X$ and $Y$

are discrete random variables, then we can characterize the support set of joint distributions that fulfill these relations by embedding the distributions from observational and interventional samples into a higher dimensional space and creating a link between them. That means we first embed the empirical distributions into a higher dimensional space and then find the best approximation of this embedding to the probability distributions that are compatible with the invariance principle such that the relative entropy between them minimizes. We call this approach an information-theoretic approximation to causal models (IACM) since the relative entropy can be interpreted as an error telling us how much a finite sample deviates from a sample that comes from an assumed underlying causal model. It turns out that solving this optimization problem is equivalent to solving a linear optimization problem, which ends up in an efficient algorithm. We use IACM to formulate a causal discovery algorithm that infers the causal direction of two random variables. For this, we approximate to a causal model were

X

causes

Y

and to a model were

Y

causes

X

. We prefer the direction that has lower relative entropy. With respective preprocessing, this can also be applied to continuous data.

If we additionally assume that the underlying causal model is monotonic w.r.t. $X$ or $Y$ , then we can include this assumption into the support set characterization used by our approach. In the case of binary random variables, an approximation to a monotonic causal model enables us to calculate probabilities about how necessary, sufficient, or necessary and sufficient a cause is for an effect as defined in Pea09 . We will use this as a strength of a causal link and include it in our causal discovery algorithm.

For the rest of this paper, we assume that we have two random variables $X$ and $Y$ that attain values in finite ranges $X_{X}$ and $X_{Y}$ , respectively. The contribution of this paper is twofold. The first contribution is an approximation of (empirical) distributions to a set of distributions that is compatible with an invariance condition induced by an assumed causal model. The second contribution is a method for causal discovery based on this approximation procedure. This method can also be applied if we have observed data from $X$ and $Y$ that are heterogeneous and continuous. In experiments, we were able to verify the strength of our causal discovery approach, especially in the case that we have discrete ranges with low cardinality, we outperformed alternative state-of-the-art methods.

The paper is organized as follows. Section 2 introduces causal models and formulates the invariance statement. In Section 3 we present an information-theoretic approximation of distributions to one that is generated by causal models. We derive the theoretic foundation, illustrate the results for the binary case, and formulate the approximation algorithm. Section 4 shows applications of the approximation algorithm. In particular, the calculation of probabilities for causes and the application to causal discovery. Section 5 describes experiments to verify our approach and we conclude in Section 6.

2 Causal Models

We describe causal relations in the form of a directed graph $G = (V, E)$ with a finite vertex set $V$ and a set of directed edges $E \subset V \times V$ . A directed edge from $u \in V$ to $v \in V$

is an ordered pair

(u, v)

and often represented as an arrow between vertices, e.g.

u \to v

. For a directed edge

(u, v)

the vertex

u

is a parent of

v

and

v

is a child of

u

. The set of parents of a vertex

u

is denoted by

{P A}_{u}

. We only consider directed graphs that have no cycles and call them directed acyclic graphs (DAGs). In a DAG we interpret the vertices as random variables

V = {X_{1}, \dots, X_{n}}

and a directed edge

(X_{i}, X_{j})

as a causal link between

X_{i}

and

X_{j}

. We say that

X_{i}

is a direct cause of

X_{j}

and

X_{j}

is a direct effect of

X_{i}

. We further specify causal links by introducing functional relations between them.

Definition 1

A structural causal model (SCM) is a tuple $C := (S, P_{N})$ where $S$ is a collection of $d$ structural assignments

X_{j} := f_{j} ({P A}_{j}, N_{j}), j = 1, \dots, d,

where ${P A}_{j} \subseteq {X_{1}, \dots, X_{d}} ∖ {X_{j}}$ are the parents of $X_{j}$ and $P_{N} = P_{N_{1}, \dots, N_{d}}$ is a joint distribution over the noise variables $N_{j}$ that are assumed to be jointly independent.

We consider an SCM as a model for a data generating process PJS17 . This enables us to model a system in an observational state and under perturbations at the same time. An SCM defines a unique distribution $P_{X}^{C}$ over the variables $X = (X_{1}, \dots, X_{d})$ . Perfect interventions are done by replacing an assignment in an SCM. Given an SCM $C$ we can replace the assignment for $X_{k}$ by $X_{k} := ~ f ({~ P A}_{k}, {~ N}_{k})$ . The distribution of that new SCM $~ C$ is denoted by $P_{X}^{~ C} =: P_{X}^{C; d o (X_{k} := ~ f ({~ P A}_{k}, {~ N}_{k}))}$ and called intervention distribution PJS17 ; Pea09 . When modeling causality, we assume the principle of independent mechanism. Roughly speaking, this principle states that a change in a variable does not change the underlying causal mechanism, see PJS17 . Formally for an SCM, this means that a change in a child variable $X$ will not change the mechanism $f$ that is responsible to obtain an effect from $X$ . From this principle the following invariance statement follows:

p^{C} (x_{j} | x_{{P A}_{j}}) = p^{C; d o (X_{k} := x)} (x_{j} | x_{{P A}_{j}}),

(1)

where $p^{C} (x_{j} | x_{{P A}_{j}})$ is the conditional density of $P_{X_{j} | X_{{P A}_{j}} = x_{{P A}_{j}}}^{C}$ evaluated at $x_{j}$ for some $k \neq j$ . Informally, this means that, if $X_{k}$ is a cause of $X_{j}$ , then it doesn’t matter if we observe $x_{j}$ when $x$ is present due to an observation of $X_{k}$ or when $x$ is present due to an intervention on $X_{k}$ .

3 Approximation to Causal Models

3.1 The General Case

Given two random variables $X, Y$ with finite ranges $X_{X}, X_{Y}$ , and data from observations of $X, Y$ as well as from interventions on $X$ or $Y$ .¹¹1We can also relax the assumption of having interventional data and assume that the data are heterogeneous and show a rich diversity. Alternatively, we can say that we have data of $X$ and $Y$ from different environments, where each environment belongs to a different intervention on $X$ or $Y$ . We further assume that the data from different interventions are independent of each other, and that there is no confounding variable. In practical applications, the interventional data can be obtained from experiments or more implicitly from heterogeneous data. Condition (1), is in general, not fulfilled by empirical distributions obtained from such data. We derive a method that enables us to find a joint probability distribution of $X$ and $Y$ that fulfill (1) and is closest to a given empirical distribution in an information-theoretic sense.

Without loss of generality we assume that the intervention took place on $X$ with values in $X_{X}$ and $X \to Y$ . In the following we assume that the elements of $X_{X}$ are in a fixed order. We summarize $X := (X, (X_{a})_{a \in X_{X}})$ , $Y := (Y, (Y_{a})_{a \in X_{X}})$ , where $X, Y$ are the observed data and $(X_{a})_{a \in X_{X}}, (Y_{a})_{a \in X_{X}}$ the interventional data. We define $V := {X, Y} \cup ⋃_{a \in X_{X}} Y_{a}$ that takes values in $XV:=XX×XY×\huge××a∈XXXYa$ and with $P_{V}$ we denote the joint distribution over $V$ . The space of probability distributions on $X_{V}$ is denoted by $P (X_{V})$ and for $A \subset V$ the marginalization of $P \in P (X_{V})$ is defined by $π_{A} : P (X_{V}) \to P (X_{A})$ with $π_{A} (P) (x) := \sum_{y \in X_{V ∖ A}} P (y, x)$ , where $x \in X_{A}$ and . The next Lemma gives us a characterization of distributions that fulfill (1).

Lemma 1

The set of joint probability distributions for $X, Y, Y_{a}$ , for all $a \in X_{X}$ which fulfill the consistency condition (1) is called $M_{C}$ and given as

	$M_{C}$	$=$	${P \in P (X_{V}) ∣ π_{X, Y, Y_{a}} P (a, y_{a}, {¯ ¯ ¯ y}_{a}) = π_{X, Y, Y_{a}} P (a, {¯ ¯ ¯ y}_{a}, y_{a}) = 0$
			$\forall y_{a} \in X_{Y_{a}}, a \in X_{X}},$

where ${¯ ¯ ¯ y}_{a} \in X_{Y_{a}} ∖ {y_{a}}$ .

The proof is given in the Supplementary Material. The support of $M_{C}$ is therefore given by

s u p p (M_{C}) = X_{V} ∖ ⋃ y \in X_{Y}, y_{x_{i}} \in X_{Y_{x_{i}}}, y = y_{x_{i}}, x_{i} \in X_{X}, i \in {1, \dots, | X_{X} |} x_{i} \times y \times X_{Y_{x_{1}}} \times \dots \times X_{Y_{x_{i - 1}}} \times y_{x_{i}} \times X_{Y_{x_{i + 1}}} \times \dots \times X_{Y_{x_{| X_{X} |}}} .

Given observational and interventional samples of $X$ and $Y$ and its corresponding empirical distributions $P_{X, Y}$ , $P_{Y_{a}}$ for $a \in X_{X}$ we try to find a distribution $^P \in M_{C}$ such that

π_{X, Y}^P = P_{X, Y}, a n d π_{Y_{a}}^P = P_{Y_{a}} f o r a \in X_{X} .

(2)

We can always find a joint distribution $^P \in P (X_{V})$ such that (2) holds, since the distributions $P_{X Y}, P_{Y_{a}}$ for all $a \in X_{X}$ are independent to each other. Although this does not guarantee $^P \in M_{C}$ , we can try to find a distribution in $M_{C}$ that has minimal relative entropy to $^P$ . This minimal relative entropy can be interpreted as an approximation error to an assumed causal model. The relative entropy or Kullback-Leibler divergence (KL-divergence) between two distributions $P, Q \in P (X_{V})$ is defined as follows

D (P | | Q) := ⎧ ⎨ ⎩ \begin{matrix} \sum_{x \in X_{V}} P (x) log (\frac{P (x)}{Q (x)}), & i f s u p p (Q) \supseteq s u p p (P), \infty, & e l s e . \end{matrix}

We use the convention that $0 log \frac{0}{q} = 0$ for $q > 0$ , see also CT91 ; Kak99 . This leads to:

min \begin{matrix} ^P \in P (X_{V}), π_{X, Y}^P = P_{X Y}, π_{Y_{a}}^P = P_{Y_{a}}, a \in X_{X} \end{matrix} min ~ P \in M_{C} D (~ P | |^P) .

(3)

That is a nonlinear min-min optimization problem with linear constraints. It turns out that in our situation, the problem simplifies to a linear optimization problem.

Proposition 1

The optimization problem (3) simplifies to the following linear optimization problem

max \begin{matrix} ^P \in P (X_{V}), π_{X, Y}^P = P_{X Y}, π_{Y_{a}}^P = P_{Y_{a}}, a \in X_{X} \end{matrix} S (^P),

with $S (^P) := \sum_{z \in s u p p (M_{C})}^P (z)$ .

The proof is given in the Supplementary Material and an application of the Lagrangian multiplier method. The statements of Proposition 1 holds also for any other support set characterization rather than $M_{C}$ . The global approximation error is given by $D (~ P | |^P) = - log (S (^P))$ . Inspired by JMZLZDSS12 ; DJMZSZS10 ; JBGS13 and by the intuition that the information loss from $π_{X}^P$ to $π_{Y}^P$ should be smaller than the other way round (due to an assumed mechanism from $X$ to $Y$ ) the quantity $D (π_{X}^P | | π_{Y}^P)$ could also be seen as a kind of approximation error to the causal model $X \to Y$ .

3.2 The Binary Case

To illustrate our approach, we consider the binary case. That means $X_{X} = X_{Y} = {0, 1}$ , and $V = {X, Y, Y_{0}, Y_{1}}$ . The set of consistent probability distributions is characterized by

M_{C} = {P \in P (X_{V}) | P_{0010} = P_{0011} = P_{0100} = P_{0101} = P_{1001} = P_{1011} = P_{1100} = P_{1110} = 0}

and therefore $s u p p (M_{C}) = {0000, 0001, 0110, 0111, 1000, 1010, 1101, 1111}$ . A probability distribution $^P \in P (X_{V})$

is a non-negative vector with

16

elements that sums up to

1

. We encode the conditions (2) into a contraint matrix

C

that takes the following form

C = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ \begin{matrix} 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 0 & 1 & 0 & 1 & 0 & 1 & 0 & 1 & 0 & 1 & 0 & 1 & 0 & 1 & 0 & 1 0 & 0 & 1 & 1 & 0 & 0 & 1 & 1 & 0 & 0 & 1 & 1 & 0 & 0 & 1 & 1 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 1 & 1 & 1 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 1 & 1 & 1 & 0 & 0 & 0 & 0 0 & 0 & 0 & 0 & 1 & 1 & 1 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \end{matrix} ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠,

and into a corresponding right-hand side

	$c$	$=$	$(1, P_{Y_{0}} (Y_{0} = 1), P_{Y_{1}} (Y_{1} = 1), P_{X, Y} (X = 1, Y = 1), P_{X, Y} (X = 1, Y = 0),$
			$P_{X, Y} (X = 0, Y = 1)) .$

The non-negativity can be encoded in an identity matrix

1_{16}

of length

16

and a zero vector

0_{16}

of length

16

as the right-hand side. A probability distribution

^P

that solves (2) is then a solution to the following linear optimization problem

min S (^P) s . t . C \cdot^P = c a n d 1_{16} \cdot^P \geq 0_{16} .

The proof of Proposition 1 tells us that a distribution $~ P$ that fulfill condition (1) and is as close as possible to $^P$ in an information-theoretic sense can be obtained by the following re-weighting of $^P$

~ P (x) := \frac{^P (x)}{S (^P)}, i f x \in s u p p (M_{C}) a n d ~ P (x) := 0, i f x \notin s u p p (M_{C}) .

3.3 Implementation

The procedure in Subsection 3.2 can be generalized for arbitrary finite ranges of $X$ and $Y$ . The pseudo-code of the algorithm is shown in Algorithm 1. The size of the finite ranges is denoted by $b_{x} := | X_{X} |$ and $b_{y} := | X_{Y} |$ . We further assume that we have for every $x \in X_{X}$ interventional data available. Therefore, the constraint matrix $C$ has dimension $b_{x} (2 b_{y} - 1) \times b_{x} b_{y}^{b_{x} + 1}$ . The first row of $C$ contains $1$ at each column, the following $b_{x} (b_{y} - 1)$ rows contain the support patterns of $P_{Y_{a}}$ and the final $b_{x} b_{y} - 1$ rows contain the support pattern of $P_{X Y}$ . The function getConstraintDistribution prepares the right hand side of $C$ accordingly. Since we assumed that the intervention took place on $X$ the underlying assumed causal model is $X \to Y$ . Note that $S (.)$ is depending on $M_{C}$ .

C \leftarrow createConstraintMatrix (b_{x}, b_{y})

c \leftarrow getConstraintDistribution (P, b_{x}, b_{y})

3:Solve LP problem:

min S (^P)

s.t.

C^P = c

and

1_{b_{x} b_{y} b_{y}^{b_{x}}}^P \geq 0_{b_{x} b_{y} b_{y}^{b_{x}}}

~ P (x) \leftarrow ⎧ ⎨ ⎩ \begin{matrix} \frac{^P (x)}{S (^P)}, & f o r x \in s u p p (M_{C}), 0, & f o r x \notin s u p p (M_{C}) . \end{matrix}

D_{X \to Y} \leftarrow - log (S (^P))

D_{X \to Y} \leftarrow D (π_{X}^P | | π_{Y}^P)

depending on the setting

6:return

~ P, D_{X \to Y}

Algorithm 1 IACM(

P

b_{x}

b_{y}

M_{C}

)

We implemented this procedure in Python and used the cvxpy

package to solve the linear program. The dimension of

C

will grow exponentially with the size of ranges for

X

and

Y

. However, it turns out that it is enough to consider low range sizes

b_{x} \leq 4, b_{y} \leq 4

. For possibly preprocessed sample data of

X

and

Y

with much higher or continuous range sizes, we apply an equal-frequency discretization based on quantiles, that is done before calculating and feeding

P

into Algorithm 1.

4 Applications

The approximation approach described in Section 3 has several applications. We describe two of them in the following subsections.

4.1 Probabilities for Causes

To measure the effect of a cause-effect relation Pearl proposed in Pea09 counterfactual statements that give information about the necessity, the sufficiency, and the necessity and sufficiency of cause-effect relations. A counterfactual statement is a do-statement in a hypothetical situation that can, in general, not observed or simulated. Formally this means we condition an SCM to an observed situation and apply a do-operator. The corresponding intervention distribution reads for example $P^{C | (X, Y) = (1, 1); d o (X = 0)} (Y = 0)$ which means the probability that $Y$ equals $0$ if $X$ would have been $0$ where indeed we observed that $X$ is $1$ and $Y$ is $1$ .

Definition 2

Let $X, Y$ be random variables in an SCM $C$ such that $X$ is a (hypothetical) cause of $Y$ and $x \in X_{X}, y \in X_{Y}$ .

The probability that $X = x$ is necessary as a cause for an effect $Y = y$ is defined as ${P N}_{x \to y} := P^{C | (X, Y) = (x, y); d o (X \in ¯ ¯ ¯ x)} (Y \in ¯ ¯ ¯ y),$ where $¯ ¯ ¯ x = X_{X} ∖ {x}$ .
The probability that $X = x$ is sufficient as a cause for an effect $Y = y$ is defined as ${P S}_{x \to y} := P^{C | (X, Y) \in (¯ ¯ ¯ x, ¯ y); d o (X = x)} (Y = y) .$
The probability that $X = x$ is necessary and sufficient as a cause for an effect $Y = y$ is defined as ${P N S}_{x \to y} := P (X = x, Y = y) {P N}_{x \to y} + P (X \in ¯ ¯ ¯ x, Y \in ¯ ¯ ¯ y) {P S}_{x \to y} .$

In general, counterfactual statements cannot be calculated from observational data and without knowing the true underlying SCM. However, Pearl identified situations in which we can exploit the presence of observational and interventional data to calculate the probabilities defined above. One such situation is when the underlying SCM is monotonic.

Definition 3

An SCM $C$ with $Y := f (X, N_{Y})$ for two random variables $X$ and $Y$ is called monotonic relative to $X$ , if and only if $f$ is monotonic in $X$ independent of $N_{Y}$ .

If $X$ and $Y$ are binary and if $Y$ is increasing monotonic relative to $X$ , then Theorem 9.2.15 in Pea09 give us

{P N}_{1 \to 1} = \frac{P (Y = 1) - P^{C; d o (x = 0)} (Y = 1)}{P (Y = 1, X = 1)},

(4)

{P S}_{1 \to 1} = \frac{P^{C; d o (x = 1)} (Y = 1) - P (Y = 1)}{P (Y = 0, X = 0)},

(5)

{P N S}_{1 \to 1} = P^{C; d o (x = 1)} (Y = 1) - P^{C; d o (x = 0)} (Y = 1) .

(6)

Similar if $Y$ is decreasing monotonic relative to $X$ , then we could also derive in the same fashion as Pearl did it the following formulas

{P N}_{0 \to 1} = \frac{P^{C; d o (x = 1)} (Y = 0) - P (Y = 0)}{P (Y = 1, X = 0)},

(7)

{P S}_{0 \to 1} = \frac{P (Y = 0) - P^{C; d o (x = 0)} (Y = 0)}{P (Y = 0, X = 1)},

(8)

{P N S}_{0 \to 1} = P^{C; d o (x = 0)} (Y = 1) - P^{C; d o (x = 1)} (Y = 1) .

(9)

By approximating empirical observational and interventional distributions to a monotonic causal model we can calculate ${P N}_{x \to y}$ , ${P S}_{x \to y}$ , and ${P N S}_{x \to y}$ . To do this, we need to further restrict the set $M_{C}$ and note that the monotonicity of $f$ implies that either $Y_{0} = 1$ and $Y_{1} = 0$ has zero probability or that $Y_{0} = 0$ and $Y_{1} = 1$ has zero probability. This means that either $P_{0110} = P_{1010} = 0$ or $P_{0001} = P_{1101} = 0$ has to hold in addition to the conditions given in $M_{C}$ . We define $M_{M_{i}} := {P \in M_{C} | P_{0110} = P_{1010} = 0}$ as the set of probability conditions with an underlying monotonic increasing data generation process and $M_{M_{d}} := {P \in M_{C} | P_{0001} = P_{1101} = 0}$ as the set of probability conditions with an underlying monotonic decreasing data generation process. An approximation in the sense of Subsection 3.1 to $M_{M_{d}}$ or $M_{M_{i}}$ instead of $M_{C}$ will only change the definition of $S (P)$ , the rest will remain the same. In order to calculate ${P N}_{x \to y}$ , ${P S}_{x \to y}$ , and ${P N S}_{x \to y}$ we approximate to $M_{M_{d}}$ and $M_{M_{d}}$ , choose the one with the least approximation error and use the formulas given above. We state the pseudo-code in Algorithm 2.

{~ P}_{i}, D_{i} \leftarrow I A C M (P, 2, 2, M_{M_{i}})

{~ P}_{d}, D_{d} \leftarrow I A C M (P, 2, 2, M_{M_{d}})

3:if

D_{i} < D_{d}

then

4: Calculate

P N, P S, P N S

using

{~ P}_{i}

and formulas (4) - (6)

5:else

6: Calculate

P N, P S, P N S

using

{~ P}_{d}

and formulas (7) - (9)

7:end if

8:return

P N

P S

P N S

Algorithm 2 CalcCausalProbabilities(

P

)

4.2 Causal Discovery

When we assume that $X \to Y$ we can test how well the given data fit that assumption and obtain an approximation error $D_{X \to Y}$ . Switching the roles of $X$ and $Y$ we get $D_{Y \to X}$ . The direction with the smallest error is the one we infer as the causal direction. If the error difference is below a small tolerance $ϵ > 0$ , we consider both directions as equal and return "no decision". If $X$ and $Y$ are binary and the error to the monotone models is smaller than to the non-monotone models, then we apply Algorithm 2 to determine $P N S$ for both directions and use this as a decision criterion for the preferred direction (the direction with the higher $P N S$ determines the direction). In general, some kind of data preprocessing and discretization before applying the causal discovery method is of advantage. In our implementation, we included several different preprocessing steps that treat $X$ and $Y$ different depending on the assumed causal direction, see Supplementary Material for more details.

{d a t a}_{X}, b_{x} \leftarrow

preprocessing of

X, Y

w.r.t.

X

{d a t a}_{Y}, b_{y} \leftarrow

preprocessing of

X, Y

w.r.t.

Y

3:if

b_{x} = b_{y} = 2

AND monotone model is preferred then

4: use CalcCausalProbabilities to get

P N S, D_{X \to Y}, D_{Y \to X}

for

X \to Y

and

Y \to X

5: If

(D_{X \to Y} - D_{Y \to X}) < ϵ

then return direction with highest

P N S

6:else

9: If

(D_{X \to Y} - D_{Y \to X}) < ϵ

then return no decision

10:end if

11:If

D_{X \to Y} < D_{Y \to X}

then return

X \to Y

else return

Y \to X

Algorithm 3 IACMDiscovery(X, Y)

5 Experiments

We test Algorithm 3 with synthetic and real-world benchmark data against alternative causal discovery methods. It runs with $b_{x} = b_{y} = 2$

, and with a heuristic preprocessing procedure described in the Supplementary Material. Depending on the input data structure, we used

D (π_{X}^P | | π_{Y}^P)

- log (S (^P))

as an approximation error, see also Supplementary Material for more details.

5.1 Pairwise Causal Discovery Methods

Among the various causal discovery approaches for continuous, discrete, nonlinear bivariate data, we select those that do not include any training of labeled cause-effect pairs to have a fair comparison. One well-known method uses additive noise models (ANM) that assume SCMs with additive noise and applies for continuous and discrete data HJMPS09 ; PJS11 . Furthermore, we select an information-geometric approach (IGCI) JMZLZDSS12 designed for continuous data and some recent methods for discrete data that are using minimal description length (CISC) BV17 , Shannon entropy (ACID) BV18 , and a compact representation of the causal mechanism (HCR) CQZZH18 . We further select conditional distribution similarity (CDS) F16 , regression error based causal inference (RECI) BJWSS18 , and (nonlinear) invariant causal prediction (nonICP, ICP) HPM18 ; PBM16 as baseline methods. For all methods we used the default parameter settings.²²2For HCR, nonICP, and ICP we use the R-packages from the references, for CISC, ACID the corresponding Python code and for ANM, IGCI, CDS, RECI the Python package causal discovery toolbox KG19 .

5.2 Synthetic Data

We generate a set of synthetic data that are different in its structure (linear, nonlinear, discrete, non-discrete) and its range size. These synthetic data consists of observed data and data from perfect interventions. We use SCMs with additive or multiplicative noise to generate these data $X := N_{X}, Y := f (X) + N_{Y}$ or $Y := f (X) * N_{Y}$ where $N_{X}, N_{Y}$ are sampled independently from a $t$

-distribution for which the degrees of freedom are chosen randomly from

{2, \dots, 10}

and we randomly decide if we use an additive or multiplicative model. The nonlinear function

f

is randomly selected between the following functions

f_{1} (x) := max (0, x)

f_{2} (x) := sin (2 π x)

, and

f_{3} (x) := s i g n (x) \sqrt{| x |}

. The linear function is given by

f (x) := α x

, where

α

is randomly selected from the interval

[- 10, 10]

. The discrete data are generated using a

k

-bins discretization. We simulate perfect interventions on

X

by setting them to every value in the range if the range is discrete and to some randomly selected value if the range is continuous. The sample size is chosen randomly from

{100, 500, 1000}

and we run

100

simulations for each configuration. Figure 1 shows the averaged accuracy of correct inferred causal direction for linear synthetic data relative to the difference in range size

| X_{X} | - | X_{Y} |

for small range sizes. Our method performs substantially better for positive differences than all alternative approaches. A similar picture can be seen in Figure 1 for nonlinear synthetic data. Therefore, it seems that our method is well suited for situations where the range size of the cause is greater than the range size of the effect. In Figure 2 we see that our method performs also well for synthetic linear and nonlinear continuous data.

5.3 Real-World Data

As a benchmark set consisting of real-world data, we use $102$ manually labeled continuous cause-effect pairs (CEP) from different contexts MPJZS16 ; DG19 . As real-world discrete data sets, we use $4$ anonymous discrete cause-effect pairs where food intolerances cause health issues (Food)³³3This dataset, given as discrete timeseries data, has been provided by the author and the causal direction has been independently confirmed by medical tests., the Pittsburgh bridges dataset (Bridge, $4$ pairs) from DG19 as it has been used in CQZZH18 , and the Abalone dataset (Abalone, $3$

pairs) from the UCI Machine Learning Repository

DG19 . For the CEP data set IACM outperforms all other methods. Also for discrete real-world data we can see in Figure 2 that IACM successfully recovers all causal directions and can keep pace with state-of-the-art methods.

(a) Accuracies of synthetic linear and nonlinear continuous data and of continuous real-world data (CEP).

6 Discussion

In this paper, we proposed a way how empirical distributions coming from observations and experiments can be approximated to ones that follow the restrictions enforced by an assumed causal model. This can be used to calculate probabilities of causation and leads to a new causal discovery method. In our experiments, we could confirm that our approach can compete with the current state-of-the-art methods, also on real-world datasets (continuous and discrete) and without the explicit knowledge of experimental data. Especially for the discrete setting in which the range size of the cause is greater than the range size of the effect, our method has advantages compared to other approaches. Since IACM ran with small range sizes $b_{x}$ and $b_{y}$ , it seems that in many cases the essential cause-effect information can be encoded with much less information than we might have in the data. This is interesting by itself and could serve as a base for future research.

Broader Impact

As all methods for causal discovery that use observational and/or data from implicit interventions the work in this paper could help to avoid unethical experiments. Furthermore, it contributes to a more relieable detection of causal relations, since it fills a gap in the existing causal discovery landscape. Therefore, our research can help during the evaluation and design of studies with few discrete features and contribute to more solid conclusions of those studies. On the other hand there is a potential risk that our method is used for data that are not following the assumptions of this paper. This may lead to wrong causal directions and to wrong conclusions, but can be avoided by checking the assumptions on the data before applying our method. Finally, it should be noted that this article may inspire future research projects in the field.

References

(BJW $^{+}$ 18) P. Bloebaum, D. Janzing, T. Washio, S. Shimizu, and B. Schölkopf, Cause-effect inference by comparing regression errors
, International Conference on Artificial Intelligence and Statistics (2018), 900–909.
(BV17) K. Budhathoki and J. Vreeken, MDL for causal inference on discrete data, 2017 IEEE International Conference on Data Mining (ICDM) (2017), 751–756.
(BV18) , Accurate causal inference on discrete data, 2018 IEEE International Conference on Data Mining (ICDM) (2018), 881–886.
(CLRS01) T. H. Cormen, C. E. Leiserson, R. L. Rivest, and C. Stein, Introduction to algorithms, MIT Press, Cambridge, 2001.
(CQZ $^{+}$ 18) R. Cai, J. Qiao, K. Zhang, Z. Zhang, and Z. Hao, Causal discovery from discrete data using hidden compact representation, Adv Neural Inf Process Syst (2018), 2666‐2674.
(CT91) T. Cover and J. Thomas, Elements of information theory, John Wiley and Sons, 1991.
(DG19) D. Dua and C. Graff, UCI machine learning repository, 2019.
(DJM $^{+}$ 10) P. Daniušis, D. Janzing, J. Mooij, J. Zscheischler, B. Steudel, K. Zhang, and B. Schölkopf, Inferring deterministic causal relations, Proceedings of the Twenty-Sixth Conference on Uncertainty in Artificial Intelligence (2010), 143–150.
(Fon19) J. Fonollosa, Conditional distribution variability measures for causality detection, Cause Effect Pairs in Machine Learning (I. Guyon, A. Statnikov, and B. Batu, eds.), Springer, 2019.
(HDPM18) C. Heinze-Deml, J. Peters, and N. Meinshausen, Invariant causal prediction for nonlinear models, Journal of Causal Inference 6 (2018), no. 2.
(HJM $^{+}$ 09) P. Hoyer, D. Janzing, J. Mooij, J. Peters, and B. Schölkopf, Nonlinear causal discovery with additive noise models, In Neural Information Processing Systems (NIPS) (2009), 689–696.
(JBGWS13) D. Janzing, D. Balduzzi, M. Grosse-Wentrup, and B. Schölkopf, Quantifying causal influences, The Annals of Statistics 41 (2013), no. 5, 2324–2358.
(JMZ $^{+}$ 12) D. Janzing, J. Mooij, K. Zhang, J. Lemeire, J. Zscheischler, P. Daniušis, B. Steudel, and B. Schölkopf, Information-geometric approach to inferring causal directions, Artificial Intelligence 182-183 (2012), 1–31.
(Kak99) Y. Kakihara, Abstract methods in information theory, World Scientific Publishing Co. Pte. Ltd., 1999.
(KG19) D. Kalainathan and O. Goudet, Causal discovery toolbox: Uncover causal relationships in python.
(MPJ $^{+}$ 16) J. Mooij, J. Peters, D. Janzing, J. Zscheischler, and S. B., Distinguishing cause from effect using observational data: methods and benchmarks, Journal of Machine Learning Research 17 (2016), 1–102.
(PBM16) J. Peters, P. Bühlmann, and N. Meinshausen, Causal inference by using invariant prediction: identification and confidence intervals
, Journal of the Royal Statistical Society: Series B (Statistical Methodology) 78 (2016), no. 5, 947–1012.
(Pea09) J. Pearl, Causality, models, reasoning, and inference, Cambridge University Press, 2009.
(PJS11) J. Peters, D. Janzing, and B. Schölkopf, Causal inference on discrete data using additive noise models, IEEE Transactions on Pattern Analysis and Machine Intelligence 33 (2011), 2436–2450.
(PJS17) J. Peters, D. Janzing, and B. Schölkopf, Elements of causal inference, MIT Press, 2017.

Supplementary Material for Information-Theoretic Approximation to Causal Models

Appendix A Proofs

a.1 Proof of Lemma 1

Lemma 1

The set of joint probability distributions for $X, Y, Y_{a}$ , for all $a \in X_{X}$ which fulfill the consistency condition (1) is called $M_{C}$ and given as

	$M_{C}$	$=$	${P \in P (X_{V}) ∣ π_{X, Y, Y_{a}} P (a, y_{a}, {¯ ¯ ¯ y}_{a}) = π_{X, Y, Y_{a}} P (a, {¯ ¯ ¯ y}_{a}, y_{a}) = 0$
			$\forall y_{a} \in X_{Y_{a}}, a \in X_{X}},$

where ${¯ ¯ ¯ y}_{a} \in X_{Y_{a}} ∖ {y_{a}}$ .

Proof.

The consistency condition (1) implies the following relation for some $P \in P (X_{V})$ and $a \in X_{X}$

π_{X, Y} P (a, y_{a})

=

π_{X, Y, Y_{a}} P (a, y_{a}, y_{a}) .

These relation implies

π_{X, Y, Y_{a}} P (a, y_{a}, {¯ ¯ ¯ y}_{a}) = π_{X, Y, Y_{a}} P (a, {¯ ¯ ¯ y}_{a}, y_{a}) = 0, f o r a \in X_{X},

which characterizes the joint distributions that satisfy (1). $□$

a.2 Proof of Proposition 1

Proposition 1

The optimization problem (3) simplifies to the following linear optimization problem

max \begin{matrix} ^P \in P (X_{V}), π_{X, Y}^P = P_{X Y}, π_{Y_{a}}^P = P_{Y_{a}}, a \in X_{X} \end{matrix} S (^P),

with $S (^P) := \sum_{z \in s u p p (M_{C})}^P (z)$ .

Proof.

We first consider the inner minimization problem of (3) for a given joint distribution $^P \in P (X_{V})$ . This is a constrained optimization problem where the constraints in $M_{C}$ are equivalent to the equation

S (~ P) = 1,

since $~ P$ is a probability distribution. Therefore, the Lagrange functional of this minimization problem reads

Λ (~ P) := D (~ P | |^P) + λ (S (~ P) - 1),

with $λ$ as Lagrange multiplier. Using the Lagrange multiplier method we obtain explicit expressions for the approximating distribution $~ P \in M_{C}$

~ P (z) = \frac{^P (z)}{S (^P)}

for $z \in s u p p (M_{C})$ and $~ P (z) = 0$ for all $z \notin s u p p (M_{C})$ . Thus we have solved the inner minimization problem explicitly and the relative entropy simplifies to

D (~ P | |^P) = - log (S (^P)) .

Therefore, we can now optimize on the space of possible joint distributions and (3) simplifies to

max \begin{matrix} ^P \in P (X_{V}), π_{X, Y}^P = P_{X Y}, π_{Y_{a}}^P = P_{Y_{a}}, a \in X_{X} \end{matrix} log (S (^P)) .

Since $log$ is a monotone function it suffices to maximize $S (^P)$ given the constraints. But this is nothing than a linear optimization problem which can be solved by linear programming using the simplex algorithm, see, for example, CLRS01 . $□$

Appendix B Application to Timeseries Data

Algorithm 1 can also be applied when we assume that the underlying causal model has a time lag of $T$ , which is $X_{t - T} \to Y_{t}$ , and the observational and interventional data have a time order. We only have to shift the incoming data for $X_{t}$ and $Y_{t}$ so that Algorithm 1 applies to $X_{t}, Y_{t + T}$ , and has to take care that we preserve the data order during preprocessing steps. If we do not know the exact time lag we can run the approximation several times with different time lags to find the approximation with the lowest error.

Appendix C Experiments

c.1 Data Preprocessing

Algorithm 3 includes several preprocessing steps for sample data from $X$ and $Y$ that can be parametrized. The preprocessing steps split the data into obersvational data and interventional data for $X$ and $Y$ accordingly. These steps are:

discrete-cluster: discretize data from $X$ and $Y$ using KBinsDiscretizer to a number of bins. Apply a $k$
-means clustering to the discretized data with fixed number of clusters. We split the data according to the variance in the identified cluster. If
$X$ is the assumed cause, then the cluster with the lowest variance in $X$ determines the observed data and the union of the other clusters the interventional data. If $Y$ is the assumed cause, then we apply the same procedure by exchanging $X$ with $Y$ .
cluster-discrete: the steps of discrete-cluster are interchanged. The data are first clustered, then discretized, and split according the variance in the clusters.

c.2 Parameter Setting in Experiments

In the experimental runs we used the following parameter configuration to run IACM. If $| X_{X} | = 2, | X_{Y} | < 5$ or $| X_{X} | < 5, | X_{Y} | = 2$ we used no preprocessing, since the binary discretization build into IACM itself seems sufficient enough. If $| X_{X} | - | X_{Y} | < 0$ we used discrete-cluster as preprocessing and used $D (π_{X}^P | | π_{Y}^P)$ as a decision criterion for the causal direction. In the other case when $| X_{X} | - | X_{Y} | \geq 0$ we used discrete-cluster but $- log (S (^P))$ as decision criterion. For data with continuous or large range sizes we used cluster-discrete as preprocessing and the approximation error as decision criterion. The number of bins used in the discretization part is chosen such that it is a little above the range size of the data. The number of clusters is determined by a simple search procedure. We apply KMeans clustering for every number of clusters between $2$

and the number of bins to the data and determine the relative entropies between the observational data and between the interventional data resulting from the cluster split as described above. The number of clusters with the lowest sum of those relative entropies is chosen for preprocessing. Furthermore, before feeding data into preprocessing we applied a robust scaling to them that is robust w.r.t. outliers. Our experiments also show that these heuristics defining the meta parameters for IACM is by far not optimal. We left that for future research.

Information-Theoretic Approximation to Causal Models

Telling Cause from Effect using MDL-based Local and Global Regression

A Primer on Causal Analysis

Causal Inference on Discrete Data using Additive Noise Models

Causal Inference on Multivariate and Mixed-Type Data

Information in probability: Another information-theoretic proof of a finite de Finetti theorem

Causal Discovery by Kernel Intrinsic Invariance Measure

Genome-Wide Association Studies: Information Theoretic Limits of Reliable Learning

1 Introduction

2 Causal Models

Definition 1

3 Approximation to Causal Models

3.1 The General Case

Lemma 1

Proposition 1

3.2 The Binary Case

3.3 Implementation

4 Applications

4.1 Probabilities for Causes

Definition 2

Definition 3

4.2 Causal Discovery

5 Experiments

5.1 Pairwise Causal Discovery Methods

5.2 Synthetic Data

5.3 Real-World Data

6 Discussion

Broader Impact

References

Supplementary Material for Information-Theoretic Approximation to Causal Models

Appendix A Proofs

a.1 Proof of Lemma 1

a.2 Proof of Proposition 1

Appendix B Application to Timeseries Data

Appendix C Experiments

c.1 Data Preprocessing

c.2 Parameter Setting in Experiments

Information-Theoretic Approximation to Causal Models

Related Research

Telling Cause from Effect using MDL-based Local and Global Regression

A Primer on Causal Analysis

Causal Inference on Discrete Data using Additive Noise Models

Causal Inference on Multivariate and Mixed-Type Data

Information in probability: Another information-theoretic proof of a finite de Finetti theorem

Causal Discovery by Kernel Intrinsic Invariance Measure

Genome-Wide Association Studies: Information Theoretic Limits of Reliable Learning

1 Introduction

2 Causal Models

Definition 1

3 Approximation to Causal Models

3.1 The General Case

Lemma 1

Proposition 1

3.2 The Binary Case

3.3 Implementation

4 Applications

4.1 Probabilities for Causes

Definition 2

Definition 3

4.2 Causal Discovery

5 Experiments

5.1 Pairwise Causal Discovery Methods

5.2 Synthetic Data

5.3 Real-World Data

6 Discussion

Broader Impact

References

Supplementary Material for Information-Theoretic Approximation to Causal Models

Appendix A Proofs

a.1 Proof of Lemma 1

a.2 Proof of Proposition 1

Appendix B Application to Timeseries Data

Appendix C Experiments

c.1 Data Preprocessing

c.2 Parameter Setting in Experiments