PaperlayerView exact revision on arXiv

Source-linked AI summary

A survey of Bayesian Network structure learning

Neville K. Kitson, Anthony C. Constantinou, Zhigao Guo, Yang Liu, Kiattikun Chobtham

arXiv:2109.11415v2cs.LGcs.AI

TL;DR

Learning Bayesian Network structure remains difficult, particularly for causal modelling, and inconsistent evaluation makes algorithm comparisons and choices challenging. This paper reviews 74 structure-learning algorithms and evaluation approaches, finding diverse methods and no consensus on the best algorithm for a given real-world problem.

  • Problem

    Learning Bayesian Network structure, especially under causal assumptions, remains challenging, while evaluation lacks an agreed process for determining algorithm effectiveness.

    Source: §6.1 · p. 55, §7 · p. 71

  • Method

    The paper comprehensively reviews 74 Bayesian Network structure-learning algorithms, their evaluation methods, comparative performance, noise handling, and incorporation of expert knowledge.

    Source: §7 · pp. 71–72

  • Results

    The review highlights inconsistent conclusions about which algorithm is best for a given problem and inconsistencies in the comparative literature.

    Source: §6.1 · p. 55, §7 · p. 72

  • Takeaways & Limitations

    The literature offers diverse structure-learning approaches but little consensus on the most appropriate algorithm for any given real-world context.

    Source: §7 · p. 72

    Comparisons are limited by the scarcity of real-world datasets with expert-validated causal graphs and by numerous confounding factors affecting algorithmic evaluation.

    Source: §6.2.3 · p. 61, §7 · p. 72

Abstract

from arXiv · show

Bayesian Networks (BNs) have become increasingly popular over the last few decades as a tool for reasoning under uncertainty in fields as diverse as medicine, biology, epidemiology, economics and the social sciences. This is especially true in real-world areas where we seek to answer complex questions based on hypothetical evidence to determine actions for intervention. However, determining the graphical structure of a BN remains a major challenge, especially when modelling a problem under causal assumptions. Solutions to this problem include the automated discovery of BN graphs from data, constructing them based on expert knowledge, or a combination of the two. This paper provides a comprehensive review of combinatoric algorithms proposed for learning BN structure from data, describing 74 algorithms including prototypical, well-established and state-of-the-art approaches. The basic approach of each algorithm is described in consistent terms, and the similarities and differences between them highlighted. Methods of evaluating algorithms and their comparative performance are discussed including the consistency of claims made in the literature. Approaches for dealing with data noise in real-world datasets and incorporating expert knowledge into the learning process are also covered.

1. INTRODUCTION

Bayesian networks support reasoning beyond prediction, including interventions and counterfactuals under causal assumptions, but accurately learning their graphical structure remains difficult. This survey reviews combinatoric algorithms for learning BN structure from data, alongside evaluation, practical considerations, and stated coverage boundaries.

Motivation

  • Black-box predictive models are inadequate for transparent decision making, motivating causal representations that support intervention and counterfactual reasoning.Pearl and Mackenzie (2018) distinguish associational, interventional, and counterfactual levels; causal Bayesian networks can address all three under causal assumptions.Source: §1 · p. 2

Bayesian networks

  • Bayesian networks compactly represent probabilistic dependencies among large collections of variables and have been applied across biology, healthcare, engineering, and environmental modelling.Examples include protein and gene networks, psychosis and cancer care, fault diagnosis, and air pollution modelling.Source: §1 · p. 2

Survey scope

  • The survey focuses on combinatoric algorithms that search or constrain the finite space of possible BN graphs, because recovering structure is harder than estimating parameters.It covers intuitive descriptions ranging from early algorithms that remain competitive to recent advances, while allowing graphs to be specified by experts in other settings.Source: §1 · p. 2

Coverage boundaries

  • The survey omits or briefly treats functional causal models, continuous optimization, mixed observational-experimental and cyclic-graph learning, and time-series structure learning covered by other reviews.Related surveys include Glymour et al. (2019), Vowels et al. (2021), Zanga et al. (2022), Moraffah et al. (2021), and Noguiera et al. (2022).Source: §1 · pp. 2–3

Survey organization

  • The review covers constraint-based, score-based, hybrid, and some non-combinatorial approaches, followed by evaluation, noise handling, expert knowledge, software, and algorithm-selection guidance.Figure 6 summarizes the evolution of the covered structure-learning algorithms.Source: §1 · p. 3

2. PRELIMINARIES

This section defines Bayesian Networks as DAGs with parameters that encode relationships and conditional independences, distinguishing observational prediction from causal intervention reasoning. It introduces d-separation, equivalence classes, BN assumptions, data types, and the computational challenge of structure learning.

BN representation

  • A Bayesian Network consists of a directed acyclic graph over variables and parameters specifying the strength and shape of their relationships.The graph contains nodes and directed edges, with an edge A⟶B representing a direct conditional relationship or, causally, a direct cause.Source: §2 · p. 3

Observational and causal interpretation

  • A BN compactly represents conditional independences for observational prediction, whereas a causal BN supports reasoning about interventions through its DAG.Directed edges define parent-child relationships between variables.Source: §2 · pp. 3–4

Causal classes

  • DAGs with identical conditional-independence relations belong to the same Markov Equivalence Class, while collider configurations create explaining-away dependencies called v-structures.In a v-structure, the collider is independent of its parents marginally but becomes dependent when conditioned on.Source: §2 · p. 4

Graphical separation

  • D-separation determines conditional independence by checking whether a conditioning set blocks every path, yielding Sepsets when variables are separated.A path is blocked when a collider and its descendants are excluded from the conditioning set, or when a non-collider is included; unblocked paths are active.Source: §2 · pp. 4–5

BN assumptions

  • BNs rely on Markov and Minimality conditions, linking d-separation to the joint distribution and requiring every graph edge to be necessary.The Markov Condition permits decomposition of the joint distribution according to parent sets, while the Minimality Condition forbids removing an edge without implying an absent independence.Source: §2 · p. 5

Structure learning

  • Structure learning determines graph structure from data, expert knowledge, or both, but is NP-hard because the DAG search space grows super-exponentially with variables.Constraint-based methods use conditional-independence tests, score-based methods maximise an objective over graphs, and hybrid methods combine both approaches; both major classes remain NP-hard under favourable restrictions.Source: §2 · pp. 6–7, Table 1 · p. 6

3. CONSTRAINT-BASED LEARNING

Constraint-based learning uses conditional-independence tests to construct graphs consistent with observed relationships, typically returning an equivalence class rather than a single DAG. The section distinguishes CPDAG representations, faithfulness assumptions, and global, local, and latent-variable algorithm families.

Core approach

  • Constraint-based algorithms use conditional-independence tests to construct a graph consistent with the data, returning the DAG equivalence class rather than a single DAG.This approach commonly assumes causal faithfulness and causal sufficiency when interpreted causally.Source: §3 · p. 7

Equivalence-class representations

  • Two DAGs are equivalent when they share adjacencies and v-structures; their equivalence class is represented by a CPDAG with invariant directed and reversible undirected edges.PDAGs encode adjacencies and v-structures, while completing implicit orientations yields the CPDAG.Source: §3 · p. 7

Faithfulness assumptions

  • Faithfulness requires every distributional independence to be reflected by the DAG, whereas unfaithfulness permits independence caused by particular CPT values despite d-connected variables.The unfaithful example demonstrates causation without association; constraint-based methods often assume the DAG is a perfect map of the distribution.Source: §3 · p. 8

Algorithm families and outputs

  • The reviewed algorithms include global learners assuming causal sufficiency, local learners producing structures that can be merged, and methods accommodating latent variables through different graph representations.Table 2 summarizes algorithm scope, output type, and key assumptions; local methods may output Markov Blankets or integrated CPDAGs after software integration.Source: §3 · p. 8, Table 2 · p. 10

3.1. Determining Conditional Independence

Conditional-independence tests determine whether node pairs are independent given a conditioning set, but thresholding and finite samples can produce errors that propagate through constraint-based structure discovery. The section describes common discrete and Gaussian tests, their decision rules, and their use as either binary decisions or association measures.

Determining Conditional Independence

  • CI tests assess whether A and B are independent given S, including S = ∅, but sample-dependent errors can remove edges prematurely and induce later incorrect edges.The tests identify relationships in the dataset, which may differ from those in the true distribution; errors are more likely with smaller samples.Source: §3.1 · p. 11
  • For discrete BNs, G2, χ2, and mutual information are common tests, while Gaussian BNs commonly use Fisher’s z-test.G2 and χ2 test conditional independence as a null hypothesis: p-values below α, typically 0.05, imply dependence, whereas p-values above α imply conditional independence.Source: §3.1 · pp. 11–12
  • Degrees of freedom depend on variable cardinalities and, under the stated definition, assume every possible value combination is present; absent combinations require a heuristic adjustment.Sprites et al. (2000) suggest reducing df by 1 for every value combination where N_abs = 0.Source: §3.1 · pp. 11–12
  • Mutual information indicates conditional independence when MI(A, B | S) = 0, and in practice assumes independence when it falls below a chosen threshold ε.The threshold may be set arbitrarily, such as 0.01 (Cheng at al., 1997), or estimated by comparing predictive accuracy; G2 differs from mutual information by a scaling factor.Source: §3.1 · p. 12
  • Fisher’s z-test evaluates whether the partial correlation between Gaussian variables, conditioned on S, is zero using a transformed statistic and its normal-distribution p-value.Partial correlations can be computed recursively with conditioning sets of increasing size. CI statistics may also quantify association, such as stronger dependence corresponding to higher mutual information.Source: §3.1 · pp. 12–13

3.2. Global Discovery Algorithms

Global discovery algorithms learn Bayesian-network structure as a whole, generally assuming causal sufficiency, meaning no latent variables confound measured relationships. The section traces the progression from exhaustive but inefficient SGS to algorithms that improve efficiency, robustness to faithfulness violations, order dependence, sample limitations, or high-order conditional-independence tests.

3.2.1. SGS Algorithm

  • SGS learns a Markov equivalence class through adjacency, v-structure, and orientation-propagation phases, but its exhaustive first phase requires n(n−1)⋅2^n−3 CI tests in the worst case.SGS is relatively stable because early CI-test errors are not highly amplified within adjacency decisions, although adjacency and v-structure errors can propagate into orientation errors.Source: §3.2.1 · pp. 13–14

3.2.2. The PC algorithm

  • PC improves SGS efficiency by testing conditioning sets in increasing size and restricting them to nodes adjacent to each endpoint, but it is less stable because mistaken edge removals can affect later tests.PC reuses adjacency-phase Sepsets for v-structure identification and applies Meek Rules; it is polynomial when node degree is bounded.Source: §3.2.2 · pp. 14–15

3.2.3. The Conservative PC (CPC) algorithm

  • CPC weakens PC’s faithfulness requirements by separating adjacency- and orientation-faithfulness, using unanimous Sepset evidence to mark v-structures and thereby producing fewer erroneous orientations.On a sample of 1,000, simulations found CPC only slightly slower than PC.Source: §3.2.3 · p. 15

3.2.4. The Very Conservative SGS (VCSGS) algorithm

  • Very Conservative SGS applies minimality and triangle faithfulness to weaken faithfulness assumptions further, but its efficiency and implementation remained unresolved.The approach can identify other faithfulness violations under this weakened condition.Source: §3.2.4 · p. 15

3.2.5. The PC-Stable Algorithm

  • PC-Stable removes PC’s dependence on CI-test order by freezing adjacencies at each conditioning-set size, using majority-rule Sepset voting, and marking orientation conflicts bidirectionally.In high-dimensional simulations, PC-Stable achieved lower SHD and much smaller ordering variance than PC, while running 3–13% slower; PC and PC-Stable behaved similarly in low dimensions.Source: §3.2.5 · pp. 15–17

3.2.6. PC-MAX algorithm

  • PC-MAX selects the highest-p-value Sepset, parallelizes adjacency and v-structure phases, and combines PC-Stable’s order-independence strategies; it improved arc orientation over PC and PC-Stable on a 1,000-variable BN and learned 20,000-variable graphs in under five minutes.The reported score-based Fast Greedy Equivalence Search was faster and more accurate still.Source: §3.2.6 · p. 17

3.2.7. Three-Phase Dependency Algorithm (TPDA)

  • TPDA reduces testing to at most O(n^4) by using mutual-information-based information flow and a drafting, thickening, and thinning process, while requiring monotone-faithfulness for soundness.It then orients edges using SGS phases and was reported to have similar accuracy and efficiency to PC.Source: §3.2.7 · pp. 17–18

3.2.8. Recursive Autonomy Identification (RAI) algorithm

  • RAI recursively identifies autonomous subgraphs after lower-order CI tests, enabling independent learning that avoids many costly high-order tests; when decomposition is impossible, it behaves like PC.Across commonly evaluated BNs, the authors reported higher structural and predictive accuracy than contemporaneous algorithms including PC, fewer CI tests, and shorter runtimes.Source: §3.2.8 · p. 19

3.3. Local Discovery Algorithms

Local discovery algorithms learn each node’s neighbourhood or Markov Blanket separately, then merge local structures into a global skeleton for BN structure learning. Markov Blanket methods reduce conditional-independence testing, while parents-and-children methods target skeleton learning more directly.

3.3. Local Discovery Algorithms

  • Local algorithms learn each node’s parents-and-children set or Markov Blanket separately before merging these local structures into a global skeleton and orienting edges.Under faithfulness, a node’s Markov Blanket contains its parents, children, and spouses, shielding it from all other variables.Source: §3.3 · p. 19
  • Local structures should be symmetric, but conditional-independence test errors can create conflicts; the AND-rule includes an edge only when both nodes select each other.More sophisticated symmetry-correction approaches are also possible.Source: §3.3 · p. 20

3.3.1. Markov Blanket algorithms

  • GS (Margaritis and Thrun, 2000) grows each Markov Blanket using association-ordered tests, shrinks unnecessary nodes, then completes adjacency and performs v-structure and orientation phases.For largest Markov Blanket size b, its complexity is O(n^2 + nb^22^b), becoming exponential when b≈n but O(n^2) for sparse networks; Margaritis (2003) reported similar adjacency performance to PC and better edge orientation.Source: §3.3.1 · p. 20
  • IAMB and Inter-IAMB use conditional mutual information to improve Markov Blanket discovery, handling synthetic networks up to 1,000 nodes with better accuracy than GS.Fast-IAMB used χ2-based grouped additions, achieving similar accuracy while reducing execution time by 18-32% versus IAMB and 28-48% versus Inter-IAMB.Source: §3.3.1 · p. 20

3.3.2. Parents-and-children algorithms

  • MMPC, HITON-PC, and SI-HITON-PC learn parents-and-children sets directly, using prioritized inclusion, conditional-independence-based elimination, and different strategies for interleaving both operations.Parents-and-children sets are more directly useful for skeleton learning than Markov Blankets because spouses can be removed from the latter.Source: §3.3.2 · p. 21

3.4. Algorithms assuming the existence of latent variables

When latent confounders or selection variables may exist, observed-variable DAGs cannot represent all relevant dependencies, motivating ancestral graphs and PAG-based causal discovery. The section reviews FCI and variants that trade computational cost, assumptions, variable-type support, and graph informativeness.

3.4. Algorithms assuming the existence of latent variables

  • Latent confounders are usually handled by learning observed-variable graphs that accommodate unmeasured causes, because explicitly adding unknown latent variables can make learning intractable.DAGs cannot express some resulting dependence patterns; ancestral graphs provide a richer representation.Source: §3.4 · p. 21

3.4.1. Ancestral Graphs

  • Ancestral graphs represent observed-variable relationships under latent confounding and selection bias, are closed under conditioning and marginalisation, and encode independencies through m-separation.They extend DAGs while retaining DAGs as a special case; their edge types include directed, bidirected, and undirected edges.Source: §3.4.1 · pp. 21–22
  • Maximal ancestral graphs ensure every missing edge corresponds to a conditional independence, while PAGs represent equivalence classes of MAGs with invariant or variant endpoints.Constraint-based algorithms accounting for latent confounders and selection variables generally output PAGs, whose endpoint semantics differ from CPDAGs.Source: §3.4.1 · pp. 22–23

3.4.2. Fast Causal Inference (FCI) algorithm

  • FCI learns latent-variable causal structure by first estimating adjacencies, then refining them with Possible-D-Sep conditioning sets before repeated v-structure identification and complex orientation.Zhang (2008b) defined eleven orientation rules intended to produce a sound and complete PAG as sample size N→∞; FCI-Stable removes node-order dependence and is a common benchmark.Source: §3.4.2 · pp. 23–24

3.4.3. Really Fast Causal Inference (RFCI) algorithm

  • RFCI reduced conditional-independence testing and was around 250 times faster than FCI on some synthetic sparse high-dimensional graphs with n=500 latent-variable nodes.For a large class of graphs it produced the same PAG as FCI; when results differed, RFCI produced extra edges and therefore weaker edge meaning.Source: §3.4.3 · p. 24

3.4.4. Conservative Fast Causal Inference (CFCI) algorithm

  • CFCI weakens faithfulness by retaining only unambiguous v-structures, producing fewer arrowheads and extra PAG edges while keeping orientations closer to FCI than RFCI.The resulting PAGs had similar numbers of additional edges relative to FCI as RFCI-produced PAGs.Source: §3.4.4 · p. 24

3.4.5. Fast Causal Inference Plus (FCI+) algorithm

  • FCI+ retains ancestor-aware conditioning while efficiently identifying non-parent Sepset cases, requiring O(n^2(d+2)) CI tests, or O(PC^2), under bounded degree.The authors argued that non-parent Sepsets are relatively rare, so practical performance may approach PC despite the worst-case quadratic overhead.Source: §3.4.5 · p. 24

3.4.6. Mixed variable types – MGM-FCI-MAX algorithm

  • MGM-FCI-MAX supports mixed continuous and discrete variables through regression-based tests and parallelized phases, achieving about 30% savings with six rather than one processor core.Using an MGM graph as input yielded significantly reduced runtimes on 500-variable networks and a better precision–recall balance than CFCI and FCI.Source: §3.4.6 · p. 25

4. SCORE-BASED LEARNING

Score-based BN structure learning combines a graph-search strategy with an objective function, aiming to find high-scoring graphs despite the vast search space. Algorithms are organised mainly as approximate or exact methods according to whether they guarantee the highest-scoring graph.

Overview

  • Score-based learning uses a search strategy and objective function to evaluate candidate graphs, but exhaustive search is feasible only for problems with few variables.The search strategy determines how the graph space is traversed and the objective function scores each explored graph.Source: §4 · p. 25

Approximate and exact algorithms

  • Approximate algorithms do not guarantee the highest-scoring graph, whereas exact algorithms provide that guarantee through sound search-space pruning.Heuristic pruning may discard an optimal solution, while sound pruning preserves the optimal solution in the remaining search space.Source: §4 · pp. 25, 28

Search strategy

  • Search strategies differ by search space, traversal method, and pruning, with a greedy acyclic arc-addition procedure providing a simple approximate example.The procedure starts from an empty graph and repeatedly adds the arc that most improves the score without creating a cycle, stopping when no improving addition remains.Source: §4 · p. 25

Search spaces

  • Approximate methods explore DAG space, equivalence-class space, or node space, where node ordering constrains possible parent relationships.A node ordering is a topological ordering in which a node can only have parents higher in the ordering; DAGs and orderings can correspond to multiple alternatives.Source: §4 · pp. 25–26

Outputs and scoring

  • Approximate searches typically output a locally optimal DAG, exact searches a globally optimal DAG, and equivalence-class searches a CPDAG.Algorithm choice and score choice can be independent, although historical tables report the scores preferred when algorithms were introduced.Source: §4 · p. 28

4.1. Objective functions

BN structure-learning objectives comprise Bayesian scores, which incorporate fit and prior beliefs, and information-theoretic scores, which penalize complexity to reduce overfitting. Decomposability improves search efficiency, while score equivalence means observational data cannot distinguish DAGs within an equivalence class.

4.1. Objective functions

  • Decomposable scores recompute only nodes affected by graph changes, greatly improving computational efficiency; virtually all structure-learning algorithms use them.Scores may be Bayesian or information-theoretic, with the latter explicitly balancing goodness of fit against model complexity.Source: §4.1 · p. 28
  • Score-equivalent objectives assign identical scores to DAGs in the same equivalence class, so a returned DAG should be treated as one representative.The particular representative can depend on dataset artifacts such as variable-name ordering or dataset variable order.Source: §4.1 · p. 28

4.1.1. Bayesian scores

  • Bayesian scores assign posterior probabilities to graph structures using data and prior beliefs, enabling Bayesian Model Averaging over features such as specific arcs.For discrete data, common variants include K2, BDe, and BDeu; general BD scores require individually specified prior counts, making them impractical.Source: §4.1.1 · pp. 28–29
  • Combining observational and interventional data was most effective for identifying causally related nodes, while BDs with a marginal uniform prior outperformed BDeu with a structural uniform prior in all sixty tested combinations.BDeu is sensitive to equivalent sample size, and 11/15 real-world datasets excluded ESS = 1 from their robustness intervals; BDs is not score equivalent with missing parental combinations.Source: §4.1.1 · p. 30

4.1.2. Information-theoretic scores

  • Information-theoretic scores balance log-likelihood fit with graph complexity; BIC is widely used and can learn the true network faster than AIC, BDeu, and fNML.AIC applies a softer free-parameter penalty than BIC, while fNML performs well on small datasets and qNML is score equivalent; MIT reported better structural accuracy and data fitting than K2, BIC, and BDeu.Source: §4.1.2 · pp. 31–32

4.2. Approximate score-based algorithms

Approximate score-based algorithms seek high-scoring Bayesian-network structures without generally guaranteeing the optimum, although some are classically or high-dimensionally consistent. The section covers searches over DAGs, equivalence classes, node orderings, and ancestral graphs, including sampling and strategies for escaping local maxima.

4.2.1. Approximate search of DAG space

  • Enhanced K2 nearly always produced the highest-scoring graphs among approximate alternatives and was generally fastest, while GOBNILP guaranteed the optimum; structural accuracy was not reported.The evaluation covered well-known networks with 5–441 variables.Source: §4.2.1 · pp. 34–35
  • DAG-space methods use greedy changes, restarts, tabu search, neighbourhood expansion, biological heuristics, or MCMC; Hybrid MCMC, REV-MCMC, and broader Gibbs proposals improve mixing or accuracy over earlier samplers.Approximate methods may return local maxima, but some are asymptotically correct as N→∞ or high-dimensionally consistent as N→∞ and n→∞.Source: §4.2 · p. 33, §4.2.1 · pp. 34–35

4.2.2. Approximate search of equivalence class space

  • Equivalence-class searches reduce the space relative to DAG searches: GES is classically consistent, GIES increasingly outperformed GES with more intervention datasets or targets, and FGES-Merge exceeded all but one DREAM5 competitor.GIES supports multiple datasets sharing variables and an underlying causal model but having different intervention targets.Source: §4.2.2 · pp. 36–37

4.2.3. Approximate search of node-ordering space

  • Node-ordering searches exploit smaller, acyclic spaces; WINASOBS scored higher than IINOBS and MINOBS except when tied with MINOBS under a one-hour limit, while Partition-MCMC removes ordering bias by scoring partitions.Order-based MCMC can mix better than DAG sampling, but node-order priors can bias DAG posteriors, especially with small samples.Source: §4.2.3 · pp. 38–40

4.2.4. Approximate search of ancestral graph space

  • GSMAG extends approximate score-based search to ancestral graphs with latent variables, achieving slightly better recall than FCI and CFCI but worse precision than CFCI and considerably slower execution.Its decomposable Gaussian MAG score recomputes only c-components affected by an edge change.Source: §4.2.4 · p. 40

4.3. Exact score-based algorithms

Exact score-based algorithms guarantee the highest-scoring DAG for the observed dataset by solving constrained combinatorial optimization problems with pruning, although score and data limitations can separate that DAG from ground truth. The survey covers exact search over node-ordering, DAG, and equivalence-class spaces, including methods for faster search and plausible-graph sets.

4.3. Exact score-based algorithms

  • Exact algorithms guarantee the highest-scoring DAG for the input dataset, using pruning within an acyclic constrained optimization problem, but the result may not match ground truth because of data limitations and score bias.They guarantee optimality for the observed data regardless of sample size, unlike asymptotic guarantees; the score may favor sparser or denser graphs.Source: §4.3 · pp. 40–41

4.3.1. Exact search of node-ordering space

  • Dynamic programming and A* search the node-ordering lattice by decomposing DAGs through sink nodes; A* uses admissible best-first bounds and was several times faster than dynamic programming and much faster than Branch & Bound.Source: §4.3.1 · pp. 41–42

4.3.2. Exact search of DAG space

  • Branch & Bound and ILP methods search DAG space while preserving optimality through parent-set caches, cycle handling, integer variables, and acyclicity constraints.Pruning reduced cache size by approximately 10^2 to 10^7 times and, in the most dramatic cases, reduced the search space by hundreds of orders of magnitude; removing low-scoring parent sets can improve runtime but makes the method approximate.Source: §4.3.2 · pp. 43–45
  • Exact DAG-space algorithms also support Bayesian model averaging, while GOBNILP-DEV returns all DAGs within a fraction ε of the optimal score as a principled set of plausible graphs.The returned graphs satisfy (1 − ε) · score(G_OPT) ≤ score(G) ≤ score(G_OPT).Source: §4.3.2 · p. 45

4.3.3. Exact search of equivalence class space

  • The EC Tree uses A* search over CPDAG equivalence classes and was always faster than an earlier dynamic-programming method, occasionally by orders of magnitude.Across eleven networks with 14–23 variables, highest-scoring equivalence classes represented from around 3.7 DAGs per CPDAG to hundreds or thousands.Source: §4.3.3 · p. 45

5. HYBRID LEARNING AND OTHER APPROACHES

This section reviews hybrid and other approaches to Bayesian-network structure learning, including methods that combine constraint-based, score-based, sampling, ordering, or ancestral-graph search with functional causal models. It also covers continuous-optimization formulations that represent network structure through real-valued adjacency matrices and differentiable acyclicity constraints.

5.1.1. Restrict/Maximise in DAG Space

  • Restrict/Maximise hybrids first constrain a graph search space with conditional-independence information, then optimize a score within that restricted space; MMHC and H2PC exemplify this pattern in DAG space.MMHC outperformed GES, PC, and TPDA and handled 1,000 variables, while H2PC was more accurate than MMHC across 10 networks up to 1,836 variables but about 10 times slower at large sample sizes.Source: §5 · p. 45, §5.1.1 · p. 48
  • DAG-space hybrids include interleaved restriction and maximization, approximate or exact score search, continuous-variable extensions, and evolutionary population search.COS used exact dynamic programming in a denser MMPC-restricted space, PCB used partial-correlation tests beyond Gaussian variables, and HEA maintained population diversity through individualized significance levels.Source: §5.1.1 · pp. 46, 48

5.1.2. Restrict/Maximise in equivalence class space

  • ARGES relaxes skeleton restrictions by temporarily adding shielding edges around v-structures, with theoretical support for both classical and high-dimensional consistency.The added edges disappear during learning and are absent from the final CPDAG.Source: §5.1.2 · pp. 48–49

5.1.3. Restrict/Sampling in ordering space

  • Restricted MCMC methods sample node-ordering or partition spaces while allowing an additional parent beyond the PC-restricted space, reducing MAP-search complexity from n^k to 2^k.PC combined with Order-MCMC achieved better SHD than GES or PC on random graphs, while Viinikka et. al. (2020) improved efficiency and ancestral-pair accuracy over Partition-MCMC.Source: §5.1.3 · p. 49

5.1.4. Restrict/Maximise in ancestral graph space

  • Ancestral-graph hybrids extend restriction and scoring to latent-variable settings: M3HC and CCHM generally improved accuracy over FCI-family baselines, while CCHM used do-calculus to orient remaining edges.These methods assume linear Gaussian Bayesian networks and were generally slower than FCI and GFCI.Source: §5.1.4 · p. 49

5.2.1. Other hybrid approaches which search in DAG space

  • Other hybrids search alternative spaces or correct local-discovery errors: SaiyanH guarantees a connected DAG skeleton, GSP improves low-sample accuracy, GSPo targets latent variables, and GFCI combines GES with CI-based pruning.Symmetry Correction improved structural accuracy and fit in most cases across 16 Bayesian networks, especially when used during hybrid restriction.Source: §5.2.1 · p. 50, §5.1.5 · p. 50, §5.2.2 · p. 51

5.2.4. Other hybrid approaches which search in ancestral graph space

  • Bayesian scoring of independence constraints supports uncertain ancestral-graph learning: BCCD, ASP-based optimization, COmbINE, RFCI-BSC, and intervention-aware mFGS-BS convert unreliable or conflicting constraints into probabilistic decisions.RFCI-BSC generally improved structural accuracy over RFCI through repeated stochastic PAG learning and model averaging, although adjacency accuracy worsened at the largest variable count.Source: §5.2.4 · pp. 51–52
  • mFGS-BS outperforms COmbINE, GFCI, and RFCI-BSC when learning PAGs from observational and interventional datasets without assuming causal sufficiency.It combines FGES-derived CPDAGs, arc-frequency probabilities, majority-voting orientations, and do-calculus before resolving cycles and conflicts.Source: §5.2.4 · p. 52

5.3.2. Continuous Optimisation

  • Continuous optimisation treats the adjacency matrix as real-valued, converting structure learning into an equality-constrained problem closer to mainstream machine learning optimisation.This representation allows powerful off-the-shelf optimisers and contrasts with binary adjacency matrices used by combinatoric approaches.Source: §5.3.2 · pp. 53–54
  • NOTEARS minimises data-fitting error plus an L1 regulariser, which penalises complex DAGs and makes weighted edges encode both structure and relationship strength.The weighted adjacency representation naturally supports continuous variables and can extend to ordinal discrete variables through logistic regression.Source: §5.3.2 · p. 54
  • NOTEARS encodes acyclicity with a continuous, differentiable trace-based constraint, enabling augmented-Lagrangian optimisation while remaining approximate because the problem is non-convex.The constraint equals zero for DAGs and derives from the fact that adjacency-matrix traces count directed cycles; local minima remain possible.Source: §5.3.2 · p. 54
  • NOTEARS was more accurate than FGES on denser Erdős–Rényi and scale-free graphs across Gaussian, Exponential, and Gumbel noise, while matching FGES on sparser graphs.The comparison used 10, 20, 50, or 100 variables and sample sizes N=20 or N=1000.Source: §5.3.2 · pp. 54–55

6. PRACTICAL CONSIDERATIONS · 6.1. Algorithm Evaluation

The practical considerations section covers evaluation, performance, noise, expert knowledge, resources, and algorithm-use guidance for real-world BN structure learning. Its evaluation subsection emphasizes that no universal assessment process exists, so metrics must reflect the learning class, data-generation process, and intended application.

6. PRACTICAL CONSIDERATIONS

  • The practical-considerations section addresses algorithm evaluation and performance, data noise, expert knowledge, software and datasets, and guidance for choosing and using BN structure-learning algorithms.These topics are presented as considerations for applying structure-learning methods in real-world settings.Source: §6 · p. 55

6.1. Algorithm Evaluation

  • No agreed evaluation process exists, making algorithm effectiveness difficult to compare and producing conflicting conclusions about which algorithm is best.Evaluation criteria depend partly on the algorithm’s learning class and the data-generation process.Source: §6.1 · p. 55
  • Graphical evaluations compare learned and ground-truth graphs, whereas inference-based evaluations assess how well learned distributions fit observed distributions.On real data, the true graph is usually unknown, so most applications rely on inference-based evaluation, except when a knowledge-based graph is available.Source: §6.1 · p. 55
  • Evaluation standards should match the learning objective: causal constraint-based methods are mainly assessed structurally, while score-based methods can be judged by the highest-scoring graph.Using inference-based scores for constraint-based causal discovery can be inappropriate because those methods do not optimize a global graph score.Source: §6.1 · pp. 55–56
  • Precision, Recall, F1, SHD, SID, and balanced scoring provide complementary measures of graph discrepancy, causal-effect disagreement, and confusion-matrix imbalance.SHD counts edge insertions, deletions, and reversals; SID compares causal intervention statements and is generally asymmetric between two DAGs.Source: §6.1 · pp. 56–57
  • BIC is widely used for inference-based evaluation, but score equivalence can obscure DAG differences and higher scores or predictive accuracy need not indicate more accurate causal structure.This problem is amplified by limited and noisy data and reduced by larger, cleaner datasets, although inference-based evaluations remain useful.Source: §6.1 · pp. 57–58
  • Structural and inference metrics scale with graph size, and inference metrics also scale with sample size, complicating comparisons across networks and datasets.Scaled variants include SHD divided by the true graph’s arc count and Log-Likelihood divided by the product of variable count and sample size; BSF does not scale with variable count.Source: §6.1 · p. 58

6.2. Algorithm Performance · 6.3. Noise

Algorithm performance comparisons do not identify a universally optimal method: accuracy and efficiency depend on algorithm class, data, implementation, objectives, and evaluation choices. Noise further challenges structure learning, but specialized methods address missing data and measurement error under specific assumptions and conditions.

6.2. Algorithm Performance

  • Comparative studies support a benchmark set spanning PC-family, hill-climbing or Tabu, GES-family, and hybrid algorithms, while sampling and genetic methods are rarely included because they are nondeterministic.The diversity of algorithms and evaluation designs makes definitive context-specific recommendations difficult.Source: §6.2 · p. 58, §6.2.1 · p. 59

6.2.1. Accuracy comparisons between algorithms

  • Tabu often achieves the strongest accuracy, but performance varies by data type and comparison: it excels on discrete networks yet can produce overly dense Gaussian graphs.Tabu was most accurate in 14/20 discrete cases, including all 10 large-sample cases, but constraint-based and hybrid algorithms were more accurate for Gaussian BNs.Source: §6.2.1 · pp. 59–60

6.2.2. Efficiency comparisons between algorithms

  • Efficiency results conflict across studies: exact score-based methods are much slower, whereas approximate methods scale to tens of thousands or even one million variables with suitable computing.FGES was reported to learn 30,000-variable networks in minutes and a one-million-variable Gaussian network on a supercomputer, but its speed relative to Tabu differed across studies.Source: §6.2.2 · p. 60

6.2.3. Other factors affecting performance

  • Accuracy and runtime depend strongly on sample size, objective functions, conditional-independence tests, hyperparameters, dimensionality, and data realism; synthetic noise reduced graph accuracy by up to 37%.Increasing sample size generally improves accuracy, while CI-test choice produced over three orders of magnitude runtime variation for CPC.Source: §6.2.3 · p. 61
  • Algorithm rankings also depend on evaluation metrics, with studies rarely testing causal inferential performance or generalization to unseen data.Constantinou et al. (2021b) reported different rankings under SHD versus F1, while Malone et al. (2015) evaluated test-dataset generalizability comparatively unusually.Source: §6.2.3 · pp. 61–62

6.3. Noise

  • Noise in observed data primarily appears as missing data or measurement error, each of which can alter the relationships available for structure learning.Missingness may be MCAR, MAR, or MNAR, while measurement error makes measured values differ from their underlying true values.Source: §6.3 · p. 62, §6.3.1 · p. 62, §6.3.2 · p. 63

6.3.1. Missing Data

  • For missing data, Structural EM is general but computationally inefficient, whereas IPW- and deletion-based methods extend constraint-based and score-based learning under particular assumptions.HC-aIPW outperformed Structural EM when data were missing at random and missing not at random; test-wise deletion was shown sound for FCI and RFCI.Source: §6.3.1 · pp. 62–63

6.3.2. Measurement Error

  • Measurement error can worsen or distort learned structure, although larger samples, identifiable-model assumptions, or post-processing can mitigate its effects.FGES accuracy remained relatively high at sample size 5,000 even when 80% of variance came from measurement error; Zhang et al. (2018) established identifiability only under stated assumptions, and SED generally improved F1 and SHD.Source: §6.3.2 · pp. 63–64

6.4. Knowledge

The paper reviews soft constraints, hard constraints, and active learning as ways to incorporate expert knowledge into Bayesian-network structure learning. These approaches either influence learning before it begins, enforce structural requirements, or query experts during learning.

6.4.1. Soft Constraints

  • Soft constraints incorporate expert knowledge through priors on selected arcs, ancestral relationships, graph characteristics, or initial structures rather than specifying every possible graph.Bayesian objective functions use parameter and graph priors, but practical methods place priors on selected features or broader graph properties.Source: §6.4.1 · p. 64, §6.4 · p. 64
  • Priors favoring sparser graphs produced more accurate graphs, especially with low sample sizes, while modest ancestral knowledge reduced SHD error by around 15 to 20%.The sparsity result was observed with GOBNILP on moderately sized networks and Tabu on larger networks; the ancestral result used knowledge involving 12 of 37 nodes and a swap-equivalent operator.Source: §6.4.1 · p. 65

6.4.2. Hard Constraints

  • Hard constraints require or prohibit arcs, impose ordering or ancestral relationships, and can restrict node parent counts or parameter effects during structure learning.They have been integrated into hill-climbing, PC, Branch & Bound, ordering-based, EC Tree, and Ancestral Constraint Ordering Graph algorithms.Source: §6.4.2 · pp. 65–66
  • Around 10 ancestral constraints oriented around 30% of variant endpoints, with stronger propagation in PDAGs than PAGs, while ancestral constraints improved causal paths but could add extraneous edges.MINOBSx was more robust than CaMML because it always satisfied specified constraints; required arcs had the greatest effect in a review of several algorithms, whereas prohibiting arcs had little effect.Source: §6.4.2 · pp. 66–67

6.4.3. Active Learning

  • Active learning improved structural accuracy by querying experts about uncertain relationships, reducing SHD by around a quarter while requiring fewer judgments when queries occurred during MCMC learning.With N=50, post-processing typically used 10 to 16 queries, compared with around 5 at N=1,000; expert input also modestly improved data fitting.Source: §6.4.3 · p. 67
  • A hybrid active-learning method improved SHD by around 10% while reducing expert queries from 13-15 at N=100 to 1-4 at N=5000.The method queried experts about Markov Blankets or edges when expected information gain exceeded a threshold, with slightly better results when both learning phases used knowledge.Source: §6.4.3 · p. 67

6.5. Structure Learning Software Packages and Datasets · 6.6. Guidelines on choosing and using structure learning algorithms.

The paper catalogs software packages and evaluation datasets, then offers practical guidance for selecting, evaluating, and interpreting structure-learning algorithms. It emphasizes that no algorithm is universally best and recommends context-sensitive comparison, model averaging, and domain expertise.

6.5. Structure Learning Software Packages and Datasets

  • Open-source packages and commonly used datasets are cataloged, focusing on algorithms described in the paper rather than providing exhaustive lists.Table 7 identifies relevant software packages, while Table 8 identifies repositories of networks and datasets used to evaluate structure-learning algorithms.Source: §6.5 · pp. 67–68, Table 7 · p. 68, Table 8 · p. 68
  • Evaluation uses known reference graphs when available, but real-world datasets require model-selection scores, such as BIC, or predictive comparisons with other machine-learning approaches.Synthetic data can be generated from known or randomly generated graphs, whereas practical datasets generally lack a reference causal structure.Source: §6.5 · p. 69

6.6. Guidelines on choosing and using structure learning algorithms.

  • Algorithm choice depends partly on variable types, because implementations range from separate continuous or discrete support to restricted or unrestricted mixed networks.Support often reflects available conditional-independence tests or scores rather than a fundamental algorithmic property.Source: §6.6 · pp. 69–70
  • Only Tetrad is identified as supporting unrestricted mixtures of continuous and discrete variables through Conditional Gaussian, Degenerate Gaussian, or CI likelihood-test methods.Conditional Linear Gaussian networks allow mixed variables but restrict continuous variables from being parents of discrete variables.Source: §6.6 · p. 70
  • Exact score-based algorithms may suit low- to medium-dimensional settings up to around one hundred variables, but in-degree limits and causal complications can reduce their suitability.Latent confounding, selection bias, measurement error, and dense or scale-free graphs may limit accuracy or tractability advantages.Source: §6.6 · p. 70
  • Algorithms accounting for latent confounders and selection variables are available in Tetrad and Pcalg, yet simulations found non-accounting algorithms performed as well under latent confounding.Tetrad provides FCI, CGCI, RFCI, FCI-Max, and GIES; Pcalg provides FCI+.Source: §6.6 · p. 70
  • Practitioners should distinguish single-graph outputs from graph-equivalence representations and use knowledge or interventions to orient edges when causal inference requires a defined DAG.Constraint-based methods commonly return CPDAGs or PAGs, while hybrid and score-based methods typically return one DAG representing an equivalence class.Source: §6.6 · p. 70
  • Because rankings vary with metrics and problem characteristics, the authors recommend sensitivity analysis, bootstrapping, multiple suitable algorithms, model averaging, and domain-expert input.Model averaging can identify reliably recovered features, while sensitivity studies assess algorithm, score, CI-test, hyperparameter, column-order, and sample-size effects.Source: §6.6 · pp. 70–71

7. CONCLUDING REMARKS

The paper reviews three decades of Bayesian Network structure-learning algorithms, their assumptions, evaluation, practical limitations, and incorporation of expert knowledge. It concludes that the field remains vibrant but lacks consensus and faces major challenges in real-world validation and fair comparison.

Review scope

  • The review covers 24 constraint-based, 22 score-based, and 15 hybrid algorithms, using consistent terms to explain their principles, assumptions, similarities, differences, and evolution.It also describes pioneering and current state-of-the-art algorithms.Source: §7 · p. 71

Algorithm taxonomy

  • The paper organizes constraint-based, score-based, and hybrid methods by their tests, search spaces, guarantees, and combinations of restriction and maximization phases.Constraint-based methods include approaches with and without causal sufficiency; score-based methods include exact, approximate, and sampling algorithms.Source: §7 · pp. 71–72

Evaluation

  • Comparative evidence shows inconsistent conclusions about the best algorithms, while hill-climbing, GES, and PC remain competitive in many reviews.The paper examines graphical-discovery and inference metrics alongside accuracy and runtime in large-scale comparisons.Source: §7 · p. 72

Limitations

  • Theoretical accuracy guarantees often depend on unrealistic assumptions, and learned-graph quality varies with algorithms, tests or scores, hyperparameters, data, graphs, and evaluation metrics.These interacting factors make algorithm comparisons problematic, especially when evaluations rely on synthetic clean data.Source: §7 · p. 72

Noise robustness

  • Real-world noise such as latent variables, missing data, and measurement or discretisation error can substantially affect learning performance.The review therefore includes algorithm modifications designed to account for different forms of noise.Source: §7 · p. 72

Open questions

  • The review covers hard and soft expert constraints and active learning, while identifying open questions about algorithm choice, validated causal datasets, knowledge and intervention integration, and comparison confounders.The listed confounders include sample size, distributional assumptions, faithfulness, linearity, scores or CI tests, hyperparameters, dimensionality, noise, and evaluation metrics.Source: §7 · p. 72

APPENDIX A: GLOSSARY OF SYMBOLS

This glossary defines the notation used for Bayesian Network graphs, probabilities, datasets, statistical tests, evaluation metrics, and stochastic complexity.

Bayesian Network and data notation

  • Core Bayesian Network notation covers graphs, parameters, nodes, variables, probabilities, edges, separating sets, parents, descendants, Markov blankets, and datasets.It also defines node counts, graph counts, missing edges, data instances, sample size, and parameter-count notation.Source: Appendix A · p. 83

Statistical and evaluation notation

  • Statistical and evaluation symbols specify independence and conditioning, significance and threshold levels, degrees of freedom, test statistics, mutual information, precision, recall, and confusion-matrix outcomes.The glossary also defines partial correlation and indexing conventions for node states, parental configurations, and conditional-probability counts.Source: Appendix A · p. 83

Complexity and prior notation

  • The glossary defines r as the stochastic complexity of N values for a discrete variable with r states.Prior-parameter notation N′_ijk and θ′_ijk applies the corresponding data-count and conditional-probability symbols to prior beliefs.Source: Appendix A · p. 83

Related papers

Loading 2109.11415v2…