Unifying Gaussian LWF and AMP Chain Graphs to Model Interference

1 Motivation

Graphical models are among the most studied and used formalisms for causal inference. Some of the reasons of their success are that they make explicit and testable causal assumptions, and there exist algorithms for causal effect identification, counterfactual reasoning, mediation analysis, and model identification from non-experimental data [16], [20], [23]. However, in causal inference in general and graphical models in particular, one assumes more often than not that there is no interference, i. e. an intervention has no effect on units other than those to which the intervention was administered [28]. This may be unrealistic in some domains. For instance, vaccinating the mother against a disease may have a protective effect on her child, and vice versa. A notable exception is the work by Ogburn and VanderWeele [11], who distinguish three causal mechanisms that may give rise to interference and show how to model them with directed and acyclic graphs (DAGs). In this paper, we focus on what the authors call interference by contagion: One individual’s treatment does not affect another individual’s outcome directly but via the first individual’s outcome. Ogburn and VanderWeele [11] also argue that interference by contagion typically involves feedback among different individuals’ outcomes over time and, thus, it may be modeled by a DAG over the random variables of interest instantiated at different time points. Sometimes, however, the variables are observed at one time point only, e. g. at the end of a season. The observed variables may then be modeled by a Lauritzen-Wermuth-Frydenberg chain graph, or LWF CG for short [5], [7]: Directed edges represent direct causal effects, and undirected edges represents causal effects due to interference. Some notable papers on LWF CGs for modeling interference by contagion are Ogburn et al. [12], Shpitser [21],^[1] Shpitser et al. [22], and Tchetgen et al. [27]. For instance, the previous mother-child example may be modeled with the LWF CG in Figure 1 (a), where V1 and V2 represent the dose of the vaccine administered to the mother and the child, and D1 and D2 represent the severity of the disease. This paper only deals with Gaussian distributions, which implies that the relations between the random variables are linear. Note that the edge D1−D2 represents a symmetric relationship but has some features of a causal relationship, in the sense that a change in the severity of the disease for the mother causes a change in severity for the child and vice versa. This seems to suggest that we can interpret the undirected edges in LWF CGs as feedback loops. That is, that every undirected edge can be replaced by two directed edges in opposite directions. However, this is not correct in general, as explained at length by Lauritzen and Richardson [8].

Figure 1

Models of the mother-child example: (a, b) LWF CGs, and (c) UCG.

In some domains, we may need to model both interference and non-interference. This is the problem that we address in this paper. In our previous example, for instance, the mother and the child may have a gene that makes them healthy carriers, i. e. the higher the expression level of the gene the less the severity of the disease but the load of the disease agent (e. g., a virus) remains unaltered and, thus, so does the risk of infecting the other. We may model this situation with the LWF CG in Figure 1 (b), where G1 and G2 represent the expression level of the healthy carrier gene in the mother and the child. However, this model is not correct: In reality, the mother’s healthy carrier gene protects her but has no protective effect on the child, and vice versa. This non-interference relation is not represented by the model. In other words, LWF CGs are not expressive enough to model both interference relations (e. g., intervening on V1 must have an effect on D2) and non-interference relations (e. g., intervening on G1 must have no effect on D2). To remedy this problem, we propose to combine LWF CGs with Andersson-Madigan-Perlman chain graphs, or AMP CGs for short [1]. We call these new models unified chain graphs (UCGs). As we will show, it is possible to describe how UCGs exactly model interference in the Gaussian case. Works such as Ogburn et al. [12], Shpitser [21], Shpitser et al. [22], and Tchetgen et al. [27] use LWF CGs to model interference and compute some causal effect of interest. However, they do not describe how interference is exactly modeled.

The rest of the paper is organized as follows. Section 2 introduces some notation and definitions. Section 3 defines UCGs formally. Sections 4 and 5 define global, local and pairwise Markov properties for UCGs and prove their equivalence. Section 6 proposes an algorithm for maximum likelihood parameter estimation for UCGs, and reports experimental results. Section 7 considers causal inference in UCGs. Section 8 discusses identifiability of LWF and AMP CGs. Section 9 closes the paper with some discussion. The formal proofs of all the results are contained in Appendix A.

2 Preliminaries

In this paper, set operations like union, intersection and difference have the same precedence. When several of them appear in an expression, they are evaluated left to right unless parentheses are used to indicate a different order. Unless otherwise stated, the graphs in this paper are defined over a finite set of nodes V. Each node represents a random variable. We assume that the random variables are jointly normally distributed. The graphs contain at most one edge between any pair of nodes. The edge may be undirected or directed. We consider two types of directed edges: Solid (→) and dashed (⇢).

The parents of a set of nodes X in a graph G is the set Pa(X)={A|A→BorA⇢Bis inGwithB∈X}. The neighbors of X is the set Ne(X)={A|A−Bis inGwithB∈X}. The adjacents of X is the set Ad(X)={A|A—∘BorB—∘Ais inGwithB∈X} where —∘ means →, ⇢ or −. The non-descendants of X are Nd(X)={B|neitherA—∘⋯→⋯—∘BnorA—∘⋯⇢⋯—∘Bis inGwithA∈X}. Moreover, the subgraph of G induced by X is denoted by GX. A route between a node V1 and a node Vn in G is a sequence of (not necessarily distinct) nodes V1,…,Vn such that Vi∈Ad(Vi+1) for all 1≤i<n. Moreover, Vj,…,Vj+k and Vj+k,…,Vj with 1≤j≤n−k are called subroutes of the route. The route is called undirected if Vi−Vi+1 for all 1≤i<n. If Vn=V1, then the route is called a cycle. A cycle is called semidirected if it is of the form V1→V2—∘⋯—∘Vn or V1⇢V2—∘⋯—∘Vn. A route of distinct nodes is called a path. A chain graph (CG) is a graph with (possibly) directed and undirected edges, and without semidirected cycles. A set of nodes of a CG G is connected if there exists an undirected route in G between every pair of nodes in the set. A chain component of G is a maximal connected set. Note that the chain components of G can be sorted topologically, i. e. for every edge A→B or A⇢B in G, the component containing A precedes the component containing B. The chain components of G are denoted by Cc(G). CGs without dashed directed edges are known as LWF CGs, whereas CGs without solid directed edges are known as AMP CGs.

We now recall the interpretation of LWF CGs.^[2] Given a route ρ in a LWF CG G, a section of ρ is a maximal undirected subroute of ρ. A section C1−⋯−Cn of ρ is called a collider section if A→C1−⋯−Cn←B is a subroute of ρ. We say that ρ is Z-open with Z⊆V when (i) all the collider sections in ρ have some node in Z, and (ii) all the nodes that are outside the collider sections in ρ are outside Z. We now recall the interpretation of AMP CGs.^[3] Given a route ρ in an AMP CG G, C is a collider node in ρ if ρ has a subroute A⇢C⇠B or A⇢C−B. We say that ρ is Z-open with Z⊆V when (i) all the collider nodes in ρ are in Z, and (ii) all the non-collider nodes in ρ are outside Z. Let X, Y and Z denote three disjoint subsets of V. When there is no Z-open route in a LWF or AMP CG G between a node in X and a node in Y, we say that X is separated from Y given Z in G and denote it as X⊥GY|Z. Moreover, we represent by X⊥pY|Z that X and Y are conditionally independent given Z in a probability distribution p. We say that p satisfies the global Markov property with respect to G when X⊥pY|Z for all X,Y,Z⊆V such that X⊥GY|Z. When X⊥pY|Z if and only if X⊥GY|Z, we say that p is faithful to G. Two LWF CGs or AMP CGs are Markov equivalent if the set of distributions that satisfy the global Markov property with respect to each CG is the same. Characterizations of Markov equivalence exist. See Frydenberg [5, Theorem 5.6] for LWF CGs, and Andersson et al. [1, Theorem 5] for AMP CGs. The following theorem gives an additional characterization.

Finally, let X, Y, W and Z be disjoint subsets of V. A probability distribution p that satisfies the following five properties is called graphoid: Symmetry X⊥pY|Z⇒Y⊥pX|Z, decomposition X⊥pY∪W|Z⇒X⊥pY|Z, weak union X⊥pY∪W|Z⇒X⊥pY|Z∪W, contraction X⊥pY|Z∪W∧X⊥pW|Z⇒X⊥pY∪W|Z, and intersection X⊥pY|Z∪W∧X⊥pW|Z∪Y⇒X⊥pY∪W|Z. If p also satisfies the following property, then it is called compositional graphoid: Composition X⊥pY|Z∧X⊥pW|Z⇒X⊥pY∪W|Z. It is known that every Gaussian distribution is a compositional graphoid [26, Lemma 2.1, Proposition 2.1 and Corollary 2.4].

3 Unified chain graphs

In this section, we introduce our causal models to represent both interference and non-interference relationships. Let p denote a Gaussian distribution that satisfies the global Markov property with respect to a LWF or AMP CG G. Then,

because K⊥GK1∪⋯∪Kn∖Pa(K)|Pa(K) where K1,…,Kn denote the chain components of G that precede K in an arbitrary topological ordering of the components. Assume without loss of generality that p has zero mean vector. Moreover, let Σ and Ω denote respectively the covariance and precision matrices of p(K,Pa(K)). Each matrix is then of dimension (|K|+|Pa(K)|)×(|K|+|Pa(K)|). Bishop [2, Section 2.3.1] shows that the conditional distribution p(K|Pa(K)) is Gaussian with covariance matrix ΛK and the mean vector is a linear function of Pa(K) with coefficients βK. Then, βK is of dimension |K|×|Pa(K)|, and ΛK is of dimension |K|×|K|. Moreover, βK and ΛK can be expressed in terms of Σ and Ω as follows:

with

and

where ΣK,Pa(K) denotes the submatrix of Σ with rows K and columns Pa(K), and ΣPa(K),Pa(K)−1 denotes the inverse of ΣPa(K),Pa(K). Moreover, (ΛK−1)i,j=0 for all i,j∈K such that i−j is not in G, because i⊥Gj|K∖{i,j}∪Pa(K) [7, Proposition 5.2].

Furthermore, if G is an AMP CG, then (βK)i,j=0 for all i∈K and j∈Pa(K)∖Pa(i), because i⊥GPa(K)∖Pa(i)|Pa(i). Therefore, the directed edges of an AMP CG are suitable for representing non-interference. On the other hand, the directed edges of a LWF CG are suitable for representing interference. To see it, note that if G is a LWF CG, then Ωi,j=0 for all i∈K and j∈Pa(K)∖Pa(i), because i⊥Gj|K∖{i}∪Pa(K)∖{j}. However, (βK)i,j is not identically zero. For instance, if G={1→3−4←2} and K={3,4}, then

because (ΩK,Pa(K))4,1=0 as 1→4 is not in G. In other words, if there is a path in a LWF CG G from j to i through nodes in K, then (βK)i,j is not identically zero. Actually, (βK)i,j can be written as a sum of path weights over all such paths. This follows directly from Equation 3 and the result by Jones and West [6, Theorem 1], who show that (ΩK,K−1)i,l can be written as a sum of path weights over all the paths in G between l and i through some (but not necessarily all) nodes in K. Specifically,

where πi,l denotes the set of paths in G between i and l through nodes in K, |ρ| denotes the number of nodes in a path ρ, ρn denotes the n-th node in ρ, and (ΩK,K)∖ρ is the matrix with the rows and columns corresponding to the nodes in ρ omitted. Moreover, the determinant of a zero-dimensional matrix is taken to be 1. As a consequence, a LWF CG G does not impose zero restrictions on βK, because K is connected by definition of chain component. Previous works on the use of LWF CGs to model interference (e. g., [12], [21], [22], [27]) focus on developing methods for computing some causal effects of interest, and do not give many details on how interference is really being modeled. In the case of Gaussian LWF CGs, Equations 3 and 5 shed some light on this question. For instance, consider extending the mother-child example in Section 1 to a group of friends. The vaccine for individual j has a protective effect on individual j developing the disease, which in its turn has a protective effect on her friends, which in its turn has a protective effect on her friends’ friends, and so on. Moreover, the protective effect also works in the reverse direction, i. e. the vaccine for the latter has a protective effect on individual j. In other words, the protective effect of the vaccine for individual j on individual i is the result of all the paths for the vaccine’s effect to reach individual i through the network of friends. This is exactly what Equations 3 and 5 tell us (assuming that the neighbors of an individual’s disease node in G are her friends’ disease nodes). Appendix B gives an alternative decomposition of the covariance in terms of path coefficients, which may give further insight into Equation 5.

The discussion above suggests a way to model both interference and non-interference by unifying LWF and AMP CGs, i. e. by allowing →, ⇢ and − edges as long as no semidirected cycle exists. We call these new models unified chain graphs (UCGs). The lack of an edge − imposes a zero restriction on the elements of ΩK,K, as in LWF and AMP CGs. The lack of an edge ⇢ in an UCG imposes a zero restriction on the elements of βK, as in AMP CGs. Finally, the lack of an edge → imposes a zero restriction on the elements of ΩK,Pa(K) but not on the elements of βK, as in LWF CGs. Therefore, edges ⇢ may be used to model non-interference, whereas edges → may be used to model interference. For instance, the mother-child example in Section 1 may be modeled with the UCG in Figure 1 (c).

Every UCG is a CG, but the opposite is not true due to the following requirement. We divide the parents of a set of nodes X in an UCG G into mothers and fathers as follows. The fathers of X are Fa(X)={B|B→Ais inGwithA∈X}. The mothers of X are Mo(X)={B|B⇢Ais inGwithA∈X}. We require that Fa(K)∩Mo(K)=∅ for all K∈Cc(G). Therefore, the lack of an edge ⇢ imposes a zero restriction on the elements of βK corresponding to the mothers, and the lack of an edge → imposes a zero restriction on the elements of ΩK,Pa(K) corresponding to the fathers. In other words, G imposes the following zero restrictions:

The reason why we require that Fa(K)∩Mo(K)=∅ is to avoid imposing contradictory zero restrictions on βK, e. g. the edge j→i excludes the edge j⇢i by definition of CG, which implies that (βK)i,j is identically zero, but j→i implies the opposite. In other words, without this constraint, UCGs would reduce to AMP CGs. The following lemma formalizes this statement.

The lemma above implies that the UCG and the AMP CG can model the same Gaussian distributions. Lastly, note that the zero restrictions associated with an UCG (i. e., Equations 6-8) induce a Gaussian distribution by Equation 1 and Bishop [2, Section 2.3.3].

4 Global Markov property

In this section, we present a separation criterion for UCGs. Given a route ρ in an UCG, C is a collider node in ρ if ρ has a subroute A⇢C⇠B or A⇢C−B or A⇢C←B. A section of ρ is a maximal undirected subroute of ρ. A section C1−⋯−Cn of ρ is called a collider section if A→C1−⋯−Cn←B is a subroute of ρ. We say that ρ is Z-open with Z⊆V when (i) all the collider nodes in ρ are in Z, (ii) all the collider sections in ρ have some node in Z, and (iii) all the nodes that are outside the collider sections in ρ and are not collider nodes in ρ are outside Z. Let X, Y and Z denote three disjoint subsets of V. When there is no Z-open route in G between a node in X and a node in Y, we say that X is separated from Y given Z in G and denote it as X⊥GY|Z. Note that this separation criterion unifies the criteria for LWF and AMP CGs reviewed in Section 3. Finally, we say that a probability distribution p satisfies the global Markov property with respect to an UCG G if X⊥pY|Z for all X,Y,Z⊆V such that X⊥GY|Z.

The next theorem states the equivalence between the global Markov property and the zero restrictions associated with an UCG.

We have mentioned before that Jones and West [6] prove that the covariance between two nodes i and j can be written as a sum of path weights over the paths between i and j in a certain undirected graph (recall Equation 5 for the details). Wright [29] proves a similar result for DAGs. The following theorem generalizes both results.

5 Block-recursive, pairwise and local Markov properties

In this section, we present block-recursive, pairwise and local Markov properties for UCGs, and prove their equivalence to the global Markov property for Gaussian distributions. Equivalence means that every Gaussian distribution that satisfies any of these properties with respect to an UCG also satisfies the global Markov property with respect to the UCG, and vice versa. The relevance of these results stems from the fact that checking whether a distribution satisfies the global Markov property can now be performed more efficiently: We do not need to check that every global separation corresponds to a statistical independence in the distribution, it suffices to do so for the local (or pairwise or block-recursive) separations, which are typically considerably fewer. This is the approach taken by most learning algorithms based on hypothesis tests, such as the PC algorithm for DAGs [23], and its posterior extensions to LWF CGs [24] and AMP CGs [14]. The results in this section pave the way for developing a similar learning algorithm for UCGs, something that we postpone to a future article.

We say that a probability distribution p satisfies the block-recursive Markov property with respect to G if for any chain component K∈Cc(G) and vertex i∈K

1. i⊥pNd(K)∖K∖Fa(K)∖Mo(i)|Fa(K)∪Mo(i),

2. i⊥pNd(K)∖K∖Fa(i)∖Mo(K)|K∖{i}∪Fa(i)∪Mo(K) and

3. i⊥pj|K∖{i,j}∪Pa(K) for all j∈K∖{i}∖Ne(i).

We say that a probability distribution p satisfies the pairwise Markov property with respect to an UCG G if for any non-adjacent vertices i and j of G with i∈K and K∈Cc(G)

1. i⊥pj|Nd(i)∖K∖{j} if j∈Nd(i)∖K∖Fa(K)∖Mo(i), and

2. i⊥pj|Nd(i)∖{i,j} if j∈Nd(i)∖{i}∖Fa(i)∖Mo(K).

We say that a probability distribution p satisfies the local Markov property with respect to an UCG G if for any vertex i of G with i∈K and K∈Cc(G)

1. i⊥pNd(i)∖K∖Fa(K)∖Mo(i)|Fa(K)∪Mo(i) and

2. i⊥pNd(i)∖{i}∖Pa(i)∖Ne(i)∖Mo(Ne(i))|Pa(i)∪Ne(i)∪Mo(Ne(i)).

Theorems 3 and 5-7 imply the following corollary, which summarizes our results.

6 Maximum likelihood parameter estimation

In this section, we introduce a procedure for computing maximum likelihood estimates (MLEs) of the parameters of an UCG G, i. e. the MLEs of the non-zero entries of (βK)K,Mo(K), ΩK,K and ΩK,Fa(K) for every chain component K of G. Specifically, we adapt the procedure proposed by Drton and Eichler [4] for computing MLEs of the parameters of an AMP CG. The procedure combines generalized least squares and iterative proportional fitting. Suppose that we have access to a data matrix D whose column vectors are the data instances and the rows are indexed by V. Moreover, let DX denote the data over the variables X⊆V, i. e. the rows of D corresponding to X. Note that thanks to Equation 1, the MLEs of the parameters corresponding to each chain component K of G can be obtained separately. Moreover, recall that K|Pa(K)∼N(βKPa(K),ΩK,K−1). Therefore, we may compute the MLEs of the parameters corresponding to the component K by iterating between computing the MLEs ΩˆK,K and ΩˆK,Fa(K) and thus (βˆK)K,Fa(K) while fixing (βˆK)K,Mo(K), and computing the MLEs (βˆK)K,Mo(K) while fixing (βˆK)K,Fa(K). Specifically, we initialize (βˆK)K,Mo(K) to zero and then iterate through the following steps until convergence:

1. Compute ΩˆK,K and ΩˆK,Fa(K) from the data DFa(K) and the residuals DK−(βˆK)K,Mo(K)DMo(K).

2. Compute (βˆK)K,Fa(K) from ΩˆK,K and ΩˆK,Fa(K).

3. Compute (βˆK)K,Mo(K) from the data DMo(K) and the residuals DK−(βˆK)K,Fa(K)DFa(K).

7 Causal inference

In this section, we show how to compute the effects of interventions in UCGs. Intervening on a set of variables X⊆V modifies the natural causal mechanism of X, as opposed to (passively) observing X. For simplicity, we only consider interventions that set X to constant values. We represent an intervention that sets X=x as do(X=x). Given a chain component K of an UCG G, we have that

as discussed at length in Section 3, or equivalently

where ϵK∼N(0,ΛK). We interpret the last equation as a structural equation, i. e. it defines the natural causal mechanism of the variables in K. Specifically, the natural causal mechanism of Vi∈K is given by the following structural equation:

where ϵi denotes an error variable representing the unmodelled causes of Vi. Then, the intervention do(Vi=a) amounts to replacing the last equation with the structural equation Vi=a, which we dub the interventional causal mechanism of Vi. The resulting system of equations represents the behavior of the phenomenon being modelled under the intervention do(Vi=a), i. e. p(V∖{Vi}|do(Vi=a)). Moreover, p(V∖{Vi}|do(Vi=a)) satisfies the global Markov property with respect to GV∖{Vi}. To see it, note that every {Vi}-open route in G that does not include Vi is a {Vi}-open route in GV∖{Vi} too. However, every {Vi}-open route in G that includes Vi must contain a subroute of the form A⇢Vi⇠B or A⇢Vi−B or A⇢Vi←B or A→V1−⋯−Vi−⋯−Vn←B. Note that all of these subroutes are part of the natural causal mechanism of Vi which has been replaced and, thus, they are inactive, i. e. they are not really {Vi}-open.^[4]

The single variable interventions described in the paragraph above can be generalized to sets by simply replacing the corresponding equation for each variable in the set.

Graphically, we can represent the natural and interventional causal mechanisms of the variables in an UCG G by adding a new parent to each variable. We denote the resulting UCG by G′. The new parent of a variable Vi is a variable FVi that has the same domain as Vi plus a state labeled idle: FVi=a represents the intervention do(Vi=a), whereas FVi=idle represents no intervention on Vi. In other words, FVi=a represents that the natural causal mechanism is inactive and the interventional causal mechanism is active, whereas FVi=idle represents the opposite. See Pearl [16, Section 3.2.2] for further details. In order to decide whether to augment G with FVi→Vi or FVi⇢Vi, we have to study the interventional causal mechanism. This is unlike previous works where the value set by the intervention is all that matters. Specifically, if FVi=idle then the interventional causal mechanism is inactive and, thus, it does not matter whether we add FVi→Vi or FVi⇢Vi. We may even decide not to add FVi at all. On the other hand, if FVi=a then the interventional causal mechanism is active and, thus, the type of arrow added matters: If the interventional causal mechanism is such that the intervention delivered may affect (i. e., interfere with) the rest of the variables in K, then we augment G with FVi→Vi, otherwise we augment it with FVi⇢Vi. We illustrate this with the mother-child example from Section 1. An intervention that immunizes the mother against the disease (i. e., do(D1=0)) may or may not protect the child (i. e., interfere with D2), depending on how the intervention is delivered: It protects the child when the interventional causal mechanism is the intake of some medication (different from the vaccination V1 but comparable), but it does not protect the child when the interventional causal mechanism is a gene therapy (different from the healthy carrier genotype G1 but comparable). Then, the former case should be modeled as FD1→D1, whereas the latter should be modelled as FD1⇢D1.

As mentioned, our goal is to compute or identify a causal effect p(Y|do(X=x)) with the help of a given UCG. That is, we would like to leverage the UCG to transform the causal effect of interest into an expression that contains neither the do operator nor latent variables, so that it can be estimated from observational data. Pearl’s do-calculus consists of three rules whose repeated application together with standard probability manipulations do the job for acyclic directed mixed graphs [16, Section 3.4]. We show next that the calculus carries over into UCGs. Specifically, the calculus consists of the following three rules:

1. Rule 1 (insertion/deletion of observations):

   p(Y|do(X),Z∪W)=p(Y|do(X),W)ifY⊥(GV∖X)′Z|W.

2. Rule 2 (intervention/observation exchange):

   p(Y|do(X∪Z),W)=p(Y|do(X),Z∪W)ifY⊥(GV∖X)′FZ|W∪Z.

3. Rule 3 (insertion/deletion of interventions):

   p(Y|do(X∪Z),W)=p(Y|do(X),W)ifY⊥(GV∖X)′FZ|W.

Note that the corollary above implies that the causal effect of any intervention that involves interference is identifiable in a parametrized UCG. The parameter values may be provided by an expert or estimated from data as shown in Section 6. In the latter case, all the variables in the model are assumed to be measured in the data. When the model contains latent variables, we may still perform causal effect identification via rules 1–3. It is also worth mentioning that the corollary above is generalization of causal effect identification in LWF CGs as proposed by Lauritzen and Richardson [8], which in turn is a generalization of causal effect identification in DAGs [16]. This may come as a surprise because undirected edges in UCGs represent interference relationships whereas, in the work by Lauritzen and Richardson [8], they represent dependencies in the equilibrium distribution of a dynamic system with feedback loops. However, it has been suggested that interference is nothing but dependencies in an equilibrium distribution [11], [12], [21].

8 Identifiability of LWF and AMP CGs

This section proves that identifiability of LWF and AMP CGs is possible when the error variables have equal variance. The error variable ϵA associated with a variable A∈V represents the unmodelled causes of A. In other words, given a probability distribution p that is faithful to a LWF or AMP CG G, we can identify the Markov equivalence class of G (recall Theorem 1) from p by running, for instance, the learning algorithm developed by Studený [24] for LWF CGs and by Peña [14] for AMP CGs. We prove below that we can actually identify G from p if the error variables have equal variance. We discuss this assumption at the end of the section. Our result generalizes a similar result reported by Peters and Bühlmann [18] for DAGs.

Specifically, let G denote a LWF or AMP CG. Assume that the non-zero entries of βK and ΛK−1 in Equation 2 have been selected at random. This implies that p is faithful to G with probability almost 1 by Peña [13, Theorems 1 and 2] and Levitz et al. [10, Theorem 6.1].^[5] We therefore assume faithfulness hereinafter. Moreover, we rewrite Equation 2 as

in distribution, where ϵK∼N(0,ΛK). For any Vi∈V, we have then that

in distribution, where K denotes the chain component of G that contains Vi, and ϵi denotes an error variable representing the unmodelled causes of Vi. All such error variables are jointly denoted by ϵ which is distributed according to N(0,Λ), where Λ is a block diagonal matrix whose blocks correspond to the covariance matrices ΛK for all K∈Cc(G). Moreover, we assume that the errors ϵi have equal but unknown variance λ2. Note that if the error variances are unequal and unknown but have the form Λi,i=λi2λ2 for some known ratios λi2, then we can satisfy the equal error variance assumption by rescaling each variable Vi by dividing it with λi, i. e. Vi↦Vi/λi. This implies that the linear coefficients and the errors get rescaled as βi,j↦βi,jλj/λi and ϵi↦ϵi/λi. The following lemma proves that, after the rescaling, the error covariance matrix is still positive definite and keeps all the previous (in)dependencies which implies that, after the rescaling, p is still faithful to G with probability almost 1.

We are now ready to state formally the main result of this section.

Note that the theorem above implies that two LWF CGs or AMP CGs that represent the same separations are not Markov equivalent under the constraint of equal error variances, i. e. Theorem 1 does not hold under this constraint. The suitability of the equal error variances assumption should be assessed on a per domain basis. However, it may not be unrealistic to assume that it holds when the variables correspond to a relatively homogeneous set of individuals. For instance, in the case of a contagious disease, the error variable represents the unmodelled causes of an individual developing the disease, e. g. environmental factors. We may assume that these factors are the same for all the individuals, given their homogeneity. Therefore, we may assume equal error variances. We conjecture that the theorem above also holds for UCGs. However, a formal proof of this result requires first a characterization of Markov equivalence for UCGs, something that we postpone to a future article.

9 Discussion

LWF and AMP CGs are typically used to represent independence models. However, they can also be used to represent causal models. For instance, LWF CGs have been shown to be suitable for representing the equilibrium distributions of dynamic systems with feedback loops [8]. AMP CGs have been shown to be suitable for representing causal linear models with additive Gaussian noise [15]. LWF CGs have been extended into segregated graphs, which have been shown to be suitable for representing causal models with interference [21]. In this paper, we have shown how to combine LWF and AMP CGs to represent causal models of domains with both interference and non-interference relationships. Moreover, we have defined global, local and pairwise Markov properties for the new models, which we have coined unified chain graphs (UCGs), and shown that these properties are equivalent for Gaussian distributions. We have also proposed and evaluated an algorithm for computing MLEs of the parameters of an UCG. Finally, we have shown how to perform causal inference in UCGs.

It is worth mentioning that we are not the first to unify LWF and AMP CGs. Lauritzen and Sadeghi [9] recently proposed a new class of graphical models that unify many existing classes, including LWF and AMP CGs. Specifically, they consider acyclic graphs with four types of edges: Directed edges, bidirected edges, and solid and dotted undirected edges. Several edges between any pair of nodes are allowed. If their graphs only contain directed and solid undirected edges, then they coincide with LWF CGs. If they only contain directed and dotted undirected edges, then they coincide with AMP CGs. The authors develop global and pairwise Markov properties for the new models, and prove their equivalence. However, the pairwise Markov property only applies to graphs that have no dotted undirected edges. So, it does not apply to AMP CGs or superclasses of it. Other differences with our work are that no local Markov property is proposed, no parameterization or parameter learning algorithm is proposed, and no causal interpretation is given. However, the main difference with our work is that their models cannot accommodate both interference and non-interference relationships, because they rely on a single type of directed edge. For instance, our mother-child example may be modeled with a graph that contains the edges V1→D1, G1→D1 and D1−D2 or D1⋯D2 or both. However, this graph cannot represent that intervening on V1 must have an effect on D2 while intervening on G1 must not, because the paths from V1 and G1 to D2 contain the same types of edges in the same order, namely a directed edge followed by a solid or dotted undirected edge. This leads us to conclude that the models proposed by Lauritzen and Sadeghi [9] do not subsume UCGs.

Finally, we would like to mention some questions that we have not studied in this paper but which we will. We plan to extend UCGs to categorical random variables. When dealing with continuous random variables, assuming that these are jointly Gaussian simplifies the problem by restricting the relations to be linear. However, in our opinion, the main simplification that the Gaussian assumption brings is that checking whether an independence holds reduces to checking whether a linear coefficient or an entry in a precision matrix is identically zero. Discrete UCGs will not enjoy this advantage. In any case, finding a suitable/amenable parameterization of discrete UCGs is of utmost importance. Moreover, Drton [3] has shown that discrete LWF CGs are smooth models but discrete AMP CGs are not. Non-smoothness implies that some standard asymptotic distribution results (e. g., normal distribution limits for MLEs and χ2-limits for likelihood ratios) may not hold for the model at hand. Therefore, we need to investigate whether non-smoothness hinders discrete UCGs from representing interference and non-interference relationships. Other questions that we plan to study are (i) characterize Markov equivalent UCGs, (ii) develop structure learning algorithms based on the local and pairwise Markov properties, (iii) extend UCGs to model confounding via bidirected edges, and (iv) make use of the linearity of the relations for causal effect identification in UCGs along the lines in Pearl [16, Chapter 5].