Towards Characterising Bayesian Network Models under Selection

Real-life statistical samples are often not from the population of interest but from a subpopulation with fixed values for a set of selection variables. The most prominent example is case–control studies. Suppose that $S$ is a boolean variable indicating whether an individual is included (“selected”) in the study. Then a case–control sample with variables $\mathbf{V}$ may be assumed to be a random sample from the conditional distribution over $\mathbf{V}$ given $S=\text{true}$ (Spirtes et al., 2000). We refer to conditioning on fixed values of a set of variables as selection. Selection may create “spurious” dependencies, that is, dependencies that are not present in the general population, as the following example (Pearl, 2009) shows.

This phenomenon is known as selection bias, Berkson’s paradox (Berkson, 1946), and the explaining away effect (Kim and Pearl, 1983) (in the example above, one explanation for admission renders the other one less likely).

Applying algorithms that assume a random sample in their input to selection-biased data may lead to incorrect output. Thus, several approaches have been developed to deal with selection bias in various tasks. The motivation behind this work is to improve causal structure learning from selection-biased data.

Elucidating causal relationships is of utmost importance in science. A causal directed acyclic graph (DAG) represents the direct causal relationships between a set of variables (“direct” meaning not through other variables in the set). In the absence of hidden common causes, causal feedback loops, and selection bias, the causal DAG also represents the conditional independences (CIs) between the variables based on the so-called Markov condition (Neapolitan, 2004); the causal Bayesian network (BN) model is the set of distributions that satisfy the CI constraints encoded by the causal DAG. By performing hypothesis tests of CI on a sample from the probability distribution over the variables, constraint-based structure learning algorithms such as PC (Spirtes et al., 1999) can learn the causal BN model, which amounts to learning (features of) the causal DAG. In the presence of selection bias (and/or hidden common causes), however, the probability distribution over the variables may no longer be in the causal BN model and not all constraints on the distribution are CI constraints; in that case, the PC algorithm is not appropriate.

Some previous work focussed on learning causal BN models from selection-biased samples by correcting for the selection bias. Cooper (2000) devised a Bayesian method in which the biased sample is treated as a random sample with missing values for a known or unknown number of unsampled individuals, and the likelihood of the data is computed by summing over all possibilities for the missing values and the number of unsampled individuals, if unknown. This approach is computationally intractable in all but the smallest examples. In addition, it requires knowledge of the non-random-sampling process. Borboudakis and Tsamardinos (2015) devised a CI test for case–control samples with categorical variables, characterised potentially spurious links when learning the skeleton of the causal DAG using a test for random samples, and proposed the use of their specialised test on these links as a post-processing approach to removing spurious links. The drawbacks of this approach is that it is not applicable to general selection-biased samples, the joint distribution over the selection variables needs to be known, and the specialised tests are less powerful than the ones for random samples.

A more general approach to deal with the problem of selection bias is to characterise (supermodels of) BN models under selection and design algorithms for learning those models.

An ordinary Markov model is defined by the CI constraints in a BN model under marginalisation and/or selection (Shpitser et al., 2014). A maximal ancestral graph (MAG) (Richardson and Spirtes, 2002) is a graphical representation of an ordinary Markov model that can be learned using an algorithm such as FCI (Spirtes et al., 1999; Zhang, 2008). A MAG also represents ancestral relationships in the DAG of a BN model; in the case of a causal BN model, these are (indirect) causal relationships.

There are, however, additional, non-CI constraints in a BN model under marginalisation and/or selection. A nested Markov model (Richardson et al., 2017) is defined by the CI constraints and the so-called Verma constraints (Robins, 1986) in a BN model under marginalisation. These constraints may be represented by either a marginal directed acyclic graph (mDAG) (Evans, 2016) or an acyclic directed mixed graph (ADMG) (Richardson et al., 2017) and comprise all the equality constraints in the margin of a BN model over a set of categorical variables (Evans, 2018), although there there also exist inequality constraints such as Bell’s inequalities (see, for example, Wolfe et al., 2016, for details). In contrast to MAGs, mDAGs and ADMGs represent direct causal relationships, and are therefore more expressive. In the case of selection, there exist non-CI equality constraints as well. Lauritzen (1999) showed that there exist non-CI factorisation constraints. Evans and Didelez (2015) showed the existence of non-CI, non-factorisation constraints in a certain categorical BN model under conditioning; here we show that those constraints are also constraints in another BN model under selection and explicitly identify the sole constraint in the case of binary variables.

In this work, some initial results are provided that simplify the problem of characterising BN models under selection. For the sake of simplicity, the results are presented for categorical variables, but they can be easily generalised to general state spaces where a joint density with respect to a product measure exists. All proofs can be found in the appendix of this paper.

Having characterised equality constraints in BN models under selection, it may be possible to devise a graphical representation of them. Structure-learning algorithms for selection-biased data can then be developed that make use of the non-CI information in the data and potentially enable the learning of more causal relationships than is currently possible.

In this paper, sets are in boldface (e.g., $\mathbf{S}$), $X$ is used a shortcut for the singleton $\{X\}$ in places were a set is expected, and $\mathbf{A}\dot{\cup}\mathbf{B}$ denotes the union of disjoint sets $\mathbf{A}$ and $\mathbf{B}$. Random variables are denoted by capital letters (e.g., $X$) and their values by the respective lowercase letters (e.g., $x$); fixed values are denoted by a hat (e.g., $\hat{x}$). If $\mathbf{X}$ and $\mathbf{Y}$ are random variables, $\mathbf{x}\cup\mathbf{y}$ denotes the values of $\mathbf{X}\cup\mathbf{Y}$; the same is true for set operations other than the union. Probability distributions are denoted by capital letters (e.g., $P$) and their respective probability density functions by the respective lowercase letters (e.g., $p$); $p(\mathbf{x})$ is used as a shortcut for $p(\mathbf{X}=\mathbf{x})$. If $P$ is a distribution over $\mathbf{V}$, $\mathbf{X}$, $\mathbf{Y}$, and $\mathbf{Z}$ are distinct subsets of $\mathbf{V}$, and $\mathbf{X}$ and $\mathbf{Y}$ are nonempty, the conditional independence of $\mathbf{X}$ and $\mathbf{Y}$ given $\mathbf{Z}$ in $P$ is denoted by $\mathbf{X}\Perp_{P}\mathbf{Y}\mid\mathbf{Z}$. If $P$ is a distribution over $\mathbf{O}\dot{\cup}\mathbf{H}\dot{\cup}\mathbf{C}\dot{\cup}\mathbf{S}$, then the marginal/conditional distribution of $P$ over $\mathbf{O}$ given $\mathbf{C}$ and $\mathbf{S}=\hat{\mathbf{s}}$ is denoted by $P[_{\mathbf{H}}^{\mathbf{C},\mathbf{S}=\hat{\mathbf{s}}}$. Finally, graphs are in calligraphic (e.g., $\mathcal{G}$).

A graph $\mathcal{G}$ is an ordered pair $(\mathbf{V},\mathbf{E})$ of a set of nodes $\mathbf{V}$ and a set of edges $\mathbf{E}$ that connect pairs of distinct nodes in $\mathbf{V}$. If there is an edge between nodes $X$ and $Y$ in $\mathcal{G}$, then $X$ and $Y$ are adjacent in $\mathcal{G}$. The union of the sets of nodes adjacent to nodes $\mathbf{X}$ in $\mathcal{G}$ is denoted by $\mathbf{ADJ}_{\mathcal{G}}(\mathbf{X})$. A sequence of $n\geq 2$ nodes $(X_{1},\ldots,X_{n})$ such that, for $2\leq i\leq n$, $X_{i-1}$ and $X_{i}$ are adjacent, is called a path from $X_{1}$ to $X_{n}$. $\mathcal{G}$ is connected if there is a path from every node to every other node in the graph. Let $p=(X_{1},\ldots,X_{n})$ be a path. The nodes $X_{2},\ldots,X_{n-1}$ are called interior nodes on $p$. Path $(X_{i},\ldots,X_{j})$, where $1\leq i<j\leq n$, is the subpath of $p$ from $X_{i}$ to $X_{j}$ and is denoted by $p(X_{i},Y_{j})$. If $X_{1}=X_{n}$, $p$ is a cycle; if also $X_{1},\ldots,X_{n-1}$ are distinct, $p$ is a simple cycle. A path is simple if no subpath is a cycle. A triple is a simple path with three nodes. A triple $(X,Z,Y)$ is shielded if $X$ and $Y$ are adjacent. A simple path is shielded if there is a shielded triple on the path. Let $\mathcal{G}_{1}=(\mathbf{V}_{1},\mathbf{E}_{1})$ and $\mathcal{G}_{2}=(\mathbf{V}_{2},\mathbf{E}_{2})$ be graphs. $\mathcal{G}_{1}$ is a subgraph of $\mathcal{G}_{2}$ (denoted by $\mathcal{G}_{1}\subseteq\mathcal{G}_{2}$) if $\mathbf{V}_{1}\subseteq\mathbf{V}_{2}$ and $\mathbf{E}_{1}\subseteq\mathbf{E}_{2}$. The induced subgraph of $\mathcal{G}$ over $\mathbf{A}\subseteq\mathbf{V}$ (denoted by $\mathcal{G}_{\mathbf{A}}$) is the graph with set of nodes $\mathbf{A}$ and the edges in $\mathcal{G}$ between nodes in $\mathbf{A}$.

A graph is called directed (resp. undirected) when its edges are directed (resp. undirected). A tree is a connected undirected graph without (simple) cycles. In a tree, a leaf is a node which is adjacent to a single node. If there is an edge $X\rightarrow Y$ in directed graph $\mathcal{G}$, then $X$ is a parent of $Y$ and $Y$ a child of $X$ in $\mathcal{G}$; the edge is said to be out of $X$ and into $Y$. The union of the sets of parents (resp. children) of nodes $\mathbf{X}$ in $\mathcal{G}$ is denoted by $\mathbf{PA}_{\mathcal{G}}(\mathbf{X})$ (resp. $\mathbf{CH}_{\mathcal{G}}(\mathbf{X})$). The set $X\cup\mathbf{PA}_{\mathcal{G}}(X)$ is called the family of $X$ in $\mathcal{G}$. A node without children is called a source (resp. sink). A link is an edge without regard of direction, and the skeleton of a directed graph is the undirected graph whose edges corresponds to links in the directed graph. A triple $(X,Y,Z)$ such that $X\rightarrow Z\leftarrow Y$ is a collider. A path from $X$ to $Y$ is out of (resp. into) $X$ and out of (resp. into) $Y$ if the first edge of the path is out of (resp. into) $X$ and the last edge is out of (resp. into) $Y$. A simple path from $X$ to $Y$ where all edges are directed towards $Y$ is called directed. If there is a directed path from $X$ to $Y$ or $X=Y$, then $X$ is an ancestor of $Y$ and $Y$ a descendant of $X$. The union of the sets of ancestors (resp. descendants) of nodes $\mathbf{X}$ in $\mathcal{G}$ is denoted by $\mathbf{AN}_{\mathcal{G}}(\mathbf{X})$ (resp. $\mathbf{DE}_{\mathcal{G}}(\mathbf{X})$). A simple cycle $(X_{1},\ldots,X_{n})$ is directed if for $2\leq i\leq n$, the edge between $X_{i-1}$ and $X_{i}$ is directed towards $X_{i}$. A DAG is a directed graph without directed cycles. A conditional DAG is a DAG ($\mathbf{X}\dot{\cup}\mathbf{Y}$, $\mathbf{E}$) such that the nodes in $\mathbf{Y}$ are sources (Evans, 2018). The nodes in $\mathbf{X}$ and $\mathbf{Y}$ are called random nodes and fixed nodes, respectively; fixed nodes are drawn in a rectangle. The result of fixing $\mathbf{A}\subseteq\mathbf{V}$ in DAG $\mathcal{G}=(\mathbf{V},\mathbf{E})$ (denoted by $\phi_{\mathbf{A}}(\mathcal{G})$) is the conditional DAG with random nodes $\mathbf{V}\setminus\mathbf{A}$, fixed nodes $\mathbf{A}$, the edges in $\mathcal{G}$ between nodes in $\mathbf{V}\setminus\mathbf{A}$, and the edges in $\mathcal{G}$ from nodes in $\mathbf{A}$ to nodes in $\mathbf{V}\setminus\mathbf{A}$. A partially directed acyclic graph (PDAG) is a partially directed graph without directed cycles.

Let $\mathbf{X}$ and $\mathbf{Y}$ be distinct sets of variables. A (conditional) model over $\mathbf{X}$ given $\mathbf{Y}$ is a set of (conditional) probability distributions over $\mathbf{X}$ given $\mathbf{Y}$.

Let $\mathbf{X}$ be a set of categorical variables, $\mathcal{G}$ be a DAG over $\mathbf{X}$, and $P\in\operatorname{\mathbf{BNM}}(\mathcal{G})$. Then $P$ equals the product of the conditional distributions of the nodes in $\mathcal{G}$ given their parents in $\mathcal{G}$: (Neapolitan, 2004, Theorem 1.4)

Suppose that conditional distributions of the nodes in $\mathcal{G}$ given their parents in $\mathcal{G}$ are specified. Then the product of the distributions is in $\operatorname{\mathbf{BNM}}(\mathcal{G})$ (Neapolitan, 2004, Theorem 1.5).

BN models are widely used for causal inference. A causal BN model is defined by a causal DAG, whose edges denote direct causal relationships between the variables (See Neapolitan, 2004, for the exact definition of (direct) causation used). In the absence of hidden common causes, causal feedback loops, and selection bias, the distribution over the variables in a causal DAG $\mathcal{G}$ is in $\operatorname{\mathbf{BNM}}(\mathcal{G})$ (Neapolitan, 2004). Therefore, structure learning from a sample amounts to learning direct causal relationships between the variables in the sample.

Different DAGs may impose the same factorisation on the distribution over the variables and, therefore, impose the same CI constraints and define the same BN model; these DAGs are called Markov equivalent and said to belong to the same Markov equivalence class. Two DAGs are in the same class if and only if they have the same skeleton and unshielded colliders. The skeleton and the unshielded colliders are the same within the class. A Markov equivalence class of DAGs can be represented by a completed partially directed acyclic graph (CPDAG), a PDAG with the same skeleton as the DAGs in the class and the directed edges that are present in every DAG in the class. Any DAG in the class can be obtained by orienting the undirected edges in the CPDAG, as long as no directed cycles or new unshielded colliders are created. Without further assumptions, structure learning algorithms cannot distinguish between Markov equivalent DAGs; therefore, DAGs can be learned up to Markov equivalence.

We start by formally defining conditioning and selection of models.

According the following lemma, model conditioning is commutative.

Let $\mathbf{O}$ and $\mathbf{S}=\{S_{1},\ldots,S_{n}\}$ be distinct sets of categorical variables, $\mathbf{F}=\{\mathbf{F}_{1},\ldots,\mathbf{F}_{n}\}$ be a set of inclusion-maximal subsets of $\mathbf{O}$ such that every variable in $\mathbf{O}$ is included in at least one set in $\mathbf{F}$, $\mathbf{E}_{i}$ ($1\leq i\leq n$) be the set of edges from each node in $\mathbf{F}_{i}$ to $S_{i}$, and $\mathcal{G}=(\mathbf{O}\cup\mathbf{S},\bigcup_{1\leq i\leq n}\mathbf{E}_{i})$. Lauritzen (1999) proved that $\operatorname{\mathbf{BNM}}(\mathcal{G})[^{\mathbf{S}=\hat{\mathbf{s}}}=\mathbf{HM}(\mathbf{F})$, which means that every categorical hierarchical model is a BN model under selection.

In fact, every BN model under selection is a submodel of a (not necessarily saturated) hierarchical model; consider the following example.

We refer to the constraints in the selection hierarchical model as factorisation constraints. In the end of this section, an example of a non-CI, non-factorisation constraint is given.

The results that follow simplify the characterisation problem.

Owing to the lemma below, it is sufficient to characterise BN models under selection when the selection nodes are sinks.

Based on the following lemma and due to commutativity of model conditioning (Lemma 1), only cases where the sets of parents of the selection nodes are non-nested need to be considered.

We note that Lemma 3 and 4 are dual to Lemma 1 and 2, respectively, about marginalisation, in Evans (2016).

The notion of conditional BN model is needed for the theorem that follows.

Clearly, a BN model is a conditional BN model.

Let $\mathbf{X}$ and $\mathbf{Y}$ be distinct sets of categorical variables, $\mathcal{G}$ be a conditional DAG with random nodes $\mathbf{X}$ and fixed nodes $\mathbf{Y}$, and $P\in\operatorname{\mathbf{BNM}}(\mathcal{G})$. The proof of Theorem 1.4 in Neapolitan (2004) can be easily adapted to show that $P$ equals the product of the conditional distributions of the random nodes in $\mathcal{G}$ given their parents in $\mathcal{G}$:

Suppose that conditional distributions of the nodes in $\mathcal{G}$ given their parents in $\mathcal{G}$ are specified. It is straightforward to adapt the proof of Theorem 1.5 in Neapolitan (2004) to show that the product of the distributions is in $\operatorname{\mathbf{BNM}}(\mathcal{G})$.

According to the theorem below, selection only affects the ancestors of the selection nodes; it is therefore sufficient to characterise the case where all nodes are ancestors of the selection nodes.

In general, different Markov equivalent DAGs have different ancestors of the selection nodes. When computing constraints in a model under selection, Theorem 1 may be applied to a Markov equivalent DAG with the fewest ancestors of the selection nodes. As it turns out, such a DAG contains only the compelled ancestors of the selection nodes. Let $\mathbf{G}$ be a Markov equivalence class of DAGs. $X$ is a compelled ancestor of $Y$ in $\mathbf{G}$ if $X\in\mathbf{AN}_{\mathcal{G}}(Y)$ for every $\mathcal{G}\in\mathbf{G}$. Clearly, compelled ancestry is a transitive relation. The union of the sets of compelled ancestors of $\mathbf{Y}$ in $\mathbf{G}$ is denoted by $\mathbf{compAN}_{\mathbf{G}}(\mathbf{Y})$. Note that, if $\mathbf{G}$ is the Markov equivalence class of an (unknown) causal DAG, then $\mathbf{compAN}_{\mathbf{G}}(\mathbf{Y})$ is the set of definite causes of $Y$.

The theorem that follows states that, for a Markov equivalence class of DAGs, there exists a DAG in the class such that the ancestors of a set of variables in the DAG are the compelled ancestors of the variables in the class.

The lemma below is useful to obtaining such a DAG.

The following theorem gives necessary and sufficient criteria for identifying compelled ancestors from a CPDAG, thereby eliminating the need to enumerate the DAGs in the Markov equivalence class and take the intersection of the ancestors in each DAG.

Algorithm 1 can be used for the identification of compelled ancestors based on Theorem 3. The algorithm first identifies $\mathbf{A}=\mathbf{AN}_{\mathcal{P}}(\mathbf{Y})$ and then recursively finds the interior nodes $\mathbf{U}$ on the unshielded undirected paths from each member $X$ of $\mathbf{A}$ to $\mathbf{A}\setminus\{X\}$ with no interior node in $\mathbf{A}$ in $\mathcal{P}$. The second part of the algorithm is based on the lemma below.

For each member $X$ of $\mathbf{A}$, function call $\textsc{find\_uup\_interior\_nodes}(\mathcal{U},X,\mathbf{A})$ at line 7 of Algorithm 1 implicitly performs a depth-first search (DFS) in the tree of unshielded undirected paths from $X$ to $\mathbf{V}\setminus\{X\}$ in $\mathcal{P}$, starting at $X$. When either a member of $\mathbf{A}$ or a leaf in the tree is encountered, the algorithm backtracks. When at node $Y$ during backtracking, if a member of $\mathbf{A}$ was encountered beyond $Y$ (as indicated by the boolean flag $found$), $Y$ is added to $\mathbf{U}$ (unless $Y\in\mathbf{A}$). Since identifying $\mathbf{A}$ is $O(|\mathbf{V}|)$, each DFS is $O(|\mathbf{V}|)$, and a DFS is performed for each member of $\mathbf{A}$, Algorithm 1 is $O(|\mathbf{V}|^{2})$.

Let $\mathcal{G}$ be a DAG over $\mathbf{V}$, $\mathbf{G}$ be the Markov equivalence class of $\mathcal{G}$, $\mathbf{Y}\subseteq\mathbf{V}$, and $\mathbf{Z}=\mathbf{compAN}_{\mathbf{G}}(\mathbf{Y})$. A DAG $\mathcal{G}^{\prime}$ in $\mathbf{G}$ such that $\mathbf{AN}_{\mathcal{G}^{\prime}}(\mathbf{Y})=\mathbf{Z}$ can be obtained as follows. First, obtain the CPDAG $\mathcal{P}$ of $\mathbf{G}$ (using, e.g., the algorithm of Chickering, 1995). Then identify $\mathbf{Z}$ from $\mathcal{P}$ using Algorithm 1. Finally, use Algorithm 10.3 in Neapolitan (2004), orienting edges between $Z\in\mathbf{Z}$ and $X\in\mathbf{V}\setminus\mathbf{Z}$ out of $Z$ in $\mathcal{P}$ in Step 1 of the algorithm, to orient edges in $\mathcal{P}$.

Finally, the following theorem states that, in the case of a single selection node, selection reduces to conditioning on the parents of the node in the induced subgraph over the observed nodes.

We conclude this section with an example of a non-CI, non-factorisation constraint in a categorical BN model under selection.

In this work, preliminary results were provided towards characterising the constraints imposed on a distribution in a BN model under selection. Specifically, it was shown that the only cases that need to be characterised are the ones where the selection nodes have no children, the sets of parents of the selection nodes are non-nested, and all nodes are compelled ancestors of the selection nodes. In addition, in the case of a single selection node, selection reduces to conditioning on the parents of the selection node in the induced subgraph over the observed nodes. Furthermore, an algorithm was designed for identifying compelled ancestors from the CPDAG, thereby eliminating the need to enumerate the DAGs in the Markov equivalence class. This is a useful result on its own in causal structure learning, as compelled ancestors correspond to definite causes. Finally, a non-CI, non-factorisation constraint in a BN model under selection was computed for the first time.

Future work includes further reducing the characterisation problem, interpreting the constraints, devising a graph representation of them, unifying that representation with mDAGs, characterising the equivalence classes of the unified graphs, and, ultimately, devising structure-learning algorithms for the graphs.