Source-linked AI summary
Bayesian network learning with cutting planes
James Cussens
TL;DR
Exact Bayesian-network structure learning is difficult because the optimisation problem is NP-hard, while heuristic search does not guarantee an optimum. The paper formulates learning as an integer program, adds acyclicity constraints as cutting planes using SCIP, and reports exact learning on networks with up to 60 variables when parent sets are bounded. The approach is presented as particularly fast, although larger problems can be slow and the evaluation focuses on search rather than closeness to the true generating networks.
Problem
BN structure learning is an NP-hard optimisation problem, and heuristic search provides no guarantee that an optimal network has been found.
Method
The paper encodes Bayesian-network learning as an integer program in SCIP and adds cluster-based and Gomory cutting planes during branch-and-cut solving.
Results
The method performs exact learning on Bayesian networks with up to 60 variables when parent-set sizes are limited.
Takeaways & Limitations
Integer programming with dynamically generated cutting planes provides a fast approach to exact Bayesian-network learning under bounded parent-set sizes.
Takeaways & Limitations
The approach can solve larger problems slowly, and it evaluates search rather than closeness to the true data-generating networks.
Abstract
from arXiv · showhide
The problem of learning the structure of Bayesian networks from complete discrete data with a limit on parent set size is considered. Learning is cast explicitly as an optimisation problem where the goal is to find a BN structure which maximises log marginal likelihood (BDe score). Integer programming, specifically the SCIP framework, is used to solve this optimisation problem. Acyclicity constraints are added to the integer program (IP) during solving in the form of cutting planes. Finding good cutting planes is the key to the success of the approach -the search for such cutting planes is effected using a sub-IP. Results show that this is a particularly fast method for exact BN learning.
1 Introduction
Bayesian-network learning is an NP-hard optimisation problem, motivating exact methods that encode the search in integer programming. This paper uses SCIP and cutting planes to pursue exact learning under a parent-set-size limit.
- The formulation limits every DAG vertex to at most two parents.
- Dynamic programming has enabled exact learning up to around 30 vertices.
- The paper uses integer programming through SCIP for exact Bayesian-network learning and differs from earlier work in its cutting-plane search.
- BN learning selects a network maximizing a score for how well candidate structures fit observed data and prior knowledge.
- NP-hardness has led much prior work toward heuristic search without guarantees of finding an optimal Bayesian network.
2 Scoring Bayesian networks
The BDe score decomposes into local scores determined by each variable’s parent set, enabling preprocessing and integer-program encoding. Cluster-based constraints provide a tighter way to exclude cyclic digraphs.
- The BDe score is a log marginal likelihood using Dirichlet priors over Bayesian-network parameters.
- For fixed data, the total score is a sum of local scores, each determined entirely by a variable’s parent set.
- The method precomputes local scores for candidate parent sets and limits parent-set size to keep this preprocessing feasible.
- Cluster-based constraints exploit the fact that every vertex subset in a DAG contains a vertex without parents inside that subset.
- A local-score dominance rule prunes parent sets that are supersets of a better-scoring subset for the same variable.
3 Encoding the BN learning problem as an integer program
The integer program represents parent-set choices with family variables and enforces acyclicity through generalized k-cluster constraints. These constraints support cutting-plane separation and capture structural properties of feasible DAGs.
- The generalized constraint uses topological-order reasoning: the earliest k vertices in any cluster have fewer than k parents within that cluster.
- Family variables are binary indicators for whether a candidate parent set is selected as a variable’s parents.
- k-cluster-based constraints generalize 1-cluster constraints, with the latter recovered when k = 1.
- Convexity constraints ensure that each variable chooses exactly one parent set, possibly empty.
- Surplus variables measure how far a k-cluster constraint’s left-hand side exceeds its lower bound.
- Surplus values for clusters of sizes two and three determine a DAG’s Markov equivalence class.
- Good cutting planes should cut deeply, be mutually orthogonal, and remain suitably parallel to the objective function.
4 Cutting planes
The method iteratively solves LP relaxations and adds DAG-specific cutting planes, using a sub-IP to find cluster constraints that remove fractional or cyclic solutions. Gomory cuts provide a fallback when cluster-based cuts cannot be found.
- 4 Cutting planes: The search focuses on high-efficacy cuts, while collecting multiple sub-IP solutions because their joint benefit can exceed that of a single best-looking cut.The broader criteria also include cut orthogonality and parallelism with the objective function.
- 4 Cutting planes: The initial IP contains only convexity constraints, and its LP relaxation assigns values in [0, 1] to family variables.The relaxation is solved with the simplex algorithm to obtain an LP solution.
- 4 Cutting planes: A sub-IP searches for a cluster minimizing the numerator of the efficacy expression, with feasible solutions encoding vertices whose candidate parent sets avoid the cluster.The sub-IP uses binary J(W →u) variables and enforces J(W →u) = 1 iff u is in the cluster and W ∩C = ∅.
- 4 Cutting planes: Because LP solutions commonly contain cycles, cluster-based constraints are added during separation to remove the current solution and create a new relaxation.A cutting plane removes an LP solution from the feasible region by snipping off a polytope vertex.
- 4 Cutting planes: The algorithm repeatedly solves the current LP and adds cutting planes until none can be found or the LP solution corresponds to a DAG.Because an LP score is an upper bound on the optimal BN score, a suitable DAG can certify optimality.
- 4 Cutting planes: Gomory cuts are added only when no cluster-based constraints are found, helping the cutting-plane process avoid getting stuck at a fractional LP solution.They are generally shallow but separate the current solution and produce a different relaxation for further DAG-specific cuts.
- 4 Cutting planes: For integral cyclic solutions, the sub-IP quickly finds one cutting plane; for fractional solutions, it runs more slowly but can find more than 100 cuts.The integral case is equivalent to finding a cycle in a manner similar to depth-first search.
5 Results
The SCIP implementation applies cluster-based and Gomory cuts throughout branch-and-cut, with experiments measuring time to optimal Bayesian networks across datasets and parent-set limits. Delayed Gomory cuts were especially important on the largest problems, while the best combination of cluster constraints remained dataset-dependent.
- 5 Results: The implementation uses SCIP with a custom constraint handler and CPLEX 12 as the LP solver on a single-threaded Linux system.Experiments used a 64-bit Linux kernel, a dual 3GHz CPU, and 3.8Gb of RAM.
- 5 Results: Only 1-cluster constraints are searched directly; corresponding 2-cluster constraints are added afterward because a dedicated k = 2 sub-IP was not cost-effective.The paper reports that the time required for the k = 2 sub-IP was not compensated by the effectiveness of its cuts.
- 5 Results: DAG-specific cutting planes are described as key to the success of the IP approach, but solving sometimes proceeds to branch-and-bound when no further cuts separate the LP solution.Cuts are generated at the root and at other nodes of the branch-and-bound tree.
- 5 Results: 7.4 hours with 1-cluster constraints versus almost 12 hours with both 1- and 2-cluster constraints solved carpo10000, with neither option consistently better.The paper suggests that a more careful combination of the two constraint types may be preferable.
- 5 Results: The experiments vary dataset size, parent-set limits, and the number of family variables remaining after pruning.Table 1 defines n as the number of variables, m as the parent-set limit, and Families as the remaining family-variable count.
- 5 Results: Table 2 reports seconds to find optimal BNs using both 1- and 2-cluster constraints, with a 100-second sub-IP limit for marked problems.The table’s time columns are the main reported results, and the time limit applies to carpo100 and carpo10000.
- 4.2 Gomory cuts: 169 seconds with Gomory cuts solved hailfinder10000, while the no-Gomory run had not found an optimum after 85075 seconds and exhausted RAM.For hailfinder10000, 302 non-redundant Gomory cuts reduced the upper bound to −4.970794×10^5 before branching.
6 Ongoing work
The paper reports several ongoing-work extensions to its cutting-plane approach, but the tested elaborations did not consistently improve solving and some increased cost. Exact search also remains focused on optimization rather than closeness to the data-generating BN, while larger problems can be slow.
- Pre-adding constraints: Pre-adding k-cluster constraints led to slower solving and increased memory requirements.The constraints were added directly to the IP instead of being found as cutting planes.
- Primal heuristics: Two BN-specific primal heuristics were tested, but both found BNs less useful than those obtained by rounding.One generated random variable orders; the other derived an ordering from the current LP solution.
- Scaling: The approach has solved larger problems, albeit sometimes very slowly.The paper notes that SCIP must resort to branching on bigger problems.
- Scope: The evaluation focuses exclusively on the search component and does not assess closeness to the true data-generating BNs.That closeness depends on the score, Dirichlet priors, and amount of data.
- Markov equivalence: The approach does not exploit Markov equivalence classes, whose members have the same score under reasonable scoring functions.Markov equivalent BNs share an undirected skeleton and immoralities.
8 Conclusions and future work
The paper concludes that integer programming can perform exact Bayesian-network learning for datasets with many variables when parent-set sizes are bounded. Auxiliary line-variable branching sometimes substantially improves solving time, although it usually slows the search.
- Conclusions: 60 BN variables can be handled exactly on the carpo datasets when parent-set sizes are limited.The paper identifies this as its main finding.
- Conclusions: 6526 seconds and 84 seconds were required for alarm10000 and alarm, respectively, with line-variable branching.These times were better than the corresponding results reported in Tables 2 and 3.
- Conclusions: Branching on auxiliary line variables usually slowed solving, but it improved results for alarm10000 and alarm.The variables represent edges in the undirected skeleton and divide Markov equivalence classes.
7 Related work
The work builds on exact Bayesian-network learning methods, especially Jaakkola et al., while using SCIP to apply cluster constraints as cutting planes. It reports comparable or faster solving times on several shared datasets and discusses both the motivation and limits of exact single-model selection.
- Exact learning: Dynamic programming has supported exact learning for networks of around 30 vertices, while other applications may impose artificial parent-set limits.The paper identifies column generation as a possible way to overcome this limitation.
- Related IP methods: The paper uses the key 1-cluster constraints from Jaakkola et al. as cutting planes within the SCIP integer-programming framework.Unlike Jaakkola et al., the constraints are applied in a conventional existing IP framework, with a different search for good clusters.
- Motivation: An optimal BN may be much more probable than any sub-optimal BN, motivating exact selection despite concerns about returning only one model.The paper attributes the spiky landscape particularly to large datasets and notes that BDe corresponds to posterior probability under a uniform structure prior.
- Empirical comparison: 9 seconds, 1.3 minutes, and 2.55 minutes were reported for phen, wdbc, and alarm, respectively, versus Jaakkola et al.'s 8 seconds, 7 minutes, and 21 minutes.The comparison uses the same local scores.
- Exactness: Exact search can provide a guarantee that any yet-to-be-found BN scores no better than those already found.This guarantee can support a bound on the probability mass found so far.