Paper deep dive
Neural Bayesian Networks for Modular and Scalable Decision-Making Under Uncertainty
Giovanni Briglia, Marco Lippi, Stefano Mariani, Franco Zambonelli
Intelligence
Status: not_run | Model: - | Prompt: - | Confidence: 0%
Entities (0)
Relation Signals (0)
No relation signals yet.
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:Bayesian Networks (BNs) are established tools for letting agents reason under uncertainty in sequential decision-making. They split a joint distribution into local conditional mechanisms arranged on a directed acyclic graph, which makes them modular, interpretable, and ready to perform probabilistic inference. Their classical formulation, however, has limitations in both expressiveness and scale: each conditional is usually stored as a table or a linear-Gaussian form, and inference evaluates the graph one node and one query at a time. In this paper we introduce Neural Bayesian Networks (NBNs), in which each conditional is a differentiable, tensor-valued mechanism, while the graph and its factorization are unchanged. NBNs retain the modularity and local mechanisms of BNs, can represent any learnable conditional distribution, and run natively on GPUs. We also release NBN, a PyTorch-based library for learning and inference in NBNs, together with a public benchmark. Our evaluation shows gains in both expressiveness and scale: NBNs fit distributions that classical BN libraries struggle to represent, and on large networks and large batches of queries they run two orders of magnitude faster thanks to GPU parallelization. NBN is available at this https URL.
Tags
Links
- Source: https://arxiv.org/abs/2610.16902v1
- Canonical: https://arxiv.org/abs/2610.16902v1
Trouble viewing inline? Open PDF directly →
Full Text
118,958 characters extracted from source content.
Expand or collapse full text
Neural Bayesian Networks for Modular and Scalable Decision-Making Under Uncertainty Giovanni Briglia University of Modena and Reggio Emilia Reggio Emilia, Italy giovanni.briglia@unimore.it Marco Lippi University of Florence Florence, Italy marco.lippi@unifi.it Stefano Mariani University of Modena and Reggio Emilia Reggio Emilia, Italy stefano.mariani@unimore.it Franco Zambonelli University of Modena and Reggio Emilia Reggio Emilia, Italy franco.zambonelli@unimore.it Abstract Bayesian Networks (BNs) are established tools for letting agents reason under uncertainty in sequential decision-making. They split a joint distribution into local conditional mechanisms arranged on a directed acyclic graph, which makes them modular, interpretable, and ready to perform probabilistic inference. Their classical formu- lation, however, has limitations in both expressiveness and scale: each conditional is usually stored as a table or a linear-Gaussian form, and inference evaluates the graph one node and one query at a time. In this paper we introduce Neural Bayesian Networks (NBNs), in which each conditional is a differentiable, tensor-valued mechanism, while the graph and its factorization are unchanged. NBNs retain the modularity and local mechanisms of BNs, can rep- resent any learnable conditional distribution, and run natively on GPUs. We also releaseNBN, a PyTorch-based library for learning and inference in NBNs, together with a public benchmark. Our evaluation shows gains in both expressiveness and scale: NBNs fit distributions that classical BN libraries struggle to represent, and on large networks and large batches of queries they run two orders of magnitude faster thanks to GPU parallelization.NBNis available at https://github.com/Giovannibriglia/NeuralBayesianNetworks. Keywords Bayesian networks; Decision Making; Neural Networks; Model- Based Agents 1 Introduction An agent that acts under uncertainty would benefit from a model of how its environment responds to its actions to reason about prob- abilistic outcomes reliably. Model-based reinforcement learning, planning, and diagnosis all can leverage a transition model that can be learned from data and queried under evidence [20,29,37]. If the dynamics of the environment are known, at least in the sense of which variables influence which, the decision-making problem can be written as a factored Markov decision process (F-MDP) [10,19,30], whose transition model is a Bayesian network (BN) over the state variables. The BN is then the tool that lets the agent reason under uncertainty: it can predict the next state, and it can explain an observation in terms of the state that produced it. This matters most when the state is high-dimensional, because the fac- tored model stays compact where a flat one does not. Classical BNs, however, were mostly designed for discrete vari- ables and small graphs [33,44]. A conditional probability table grows combinatorially with the parents and is undefined when a parent is continuous. Additionally, both exact and approximate inference evaluate the conditionals one node and one sample at a time. These limitations stem from representing each conditional as a scalar function and from evaluating the graph sequentially. They are thus inherited by mainstream BN software packages (Table 4). We propose Neural Bayesian Network (NBN) as a new tool for autonomous agents when dealing with decision-making under un- certainty. NBNs keep the structural properties of BNs (i.e., the graph) while changing the conditionals to leverage more expressive and scalable representations based on neural networks. In an NBN each local conditional푝(푥 푖 | 푥 Pa(푖) )is a differentiable, tensor-valued mechanism that takes a batch of parent configurations and returns the full conditional distribution for each row (Figure 1). The graph and the factorization are exactly those of a BN. What changes is their computational representation: the mechanisms can be tabular, parametric, neural or non-parametric; every one of them exposes density, sampling, and gradients through one unified interface; and a whole layer of the graph is evaluated in one tensor operation on a GPU, for thousands of samples and evidence rows at once. This advancement enables BNs to be used efficiently where they were not practical before: on continuous and non-Gaussian vari- ables without reasonable discretization available, on networks of thousands of nodes, and for thousands of queries answered in one call. As such, NBNs make BNs seamless to integrate with learning- based agent architectures (e.g., reinforcement learning and cogni- tive architectures) and, more generally, with any machine learning toolchain: the NBN model is a differentiable module that is trained by gradient descent and runs on the same computational accelera- tors as a neural network (e.g., a GPU), while keeping the explicit graph, modular representation, and local mechanisms factorisation of the BNs’ symbolic model. Our contributions are: (1)Representation. We define NBNs (Section 2) and show that they are semantically equivalent to classical BNs, end-to-end differentiable, and whose inference scales with the number of active nodes instead of samples or queries. (2)Library. We releaseNBN, a PyTorch implementation with tabu- lar, linear-Gaussian, neural and non-parametric mechanisms 1 arXiv:2610.16902v1 [cs.AI] 12 Oct 2026 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli 푠 1 푠 2 푎 푠 ′ 1 푠 ′ 2 푠 1 ,푠 2 ,푎푠 ′ 1 =0 푠 ′ 1 =1 0, 0, 0.92 .08 0, 0, 1.35 .65 . . . . . . . . . 퐾 1 퐾 2 퐾 푎 rows, 퐾 ′ 1 columns Classical BN scalar CPD per node: a table, orN(푤 ⊤ 푥+푏,휎 2 ) fit by counting (MLE); queried by exact inference one node, one sample, one query at a time same DAG same Î 푖 푝 푖 푠 1 푠 2 푎 푠 ′ 1 푠 ′ 2 퐵×|Pa| M 1 (· ;휙 1 ) table· linear-Gaussian MDN· kNN· series 퐵 distr. log_prob, sample rsample,∇ 휙 one call per graph layer Neural BN tensor-valued differentiable mechanism per node, any family fit by gradient descent (or closed form); queried by VE / LW / AVI one tensor op per graph layer: all samples and all queries at once, on the GPU Figure 1: From a classical BN to an NBN. Left: each node stores a scalar conditional, evaluated one node and one configuration at a time. Right: each node stores a tensor-valued differentiable mechanismM 푖 that maps a batch of퐵parent configurations to a batch of conditional distributions. The DAG and the factorization are unchanged. and with exact, sampling-based and amortized inference en- gines, all batched on the GPU. (3)Benchmark. We release a benchmark for BN learning and in- ference and compareNBNwith state-of-the-art BN frameworks on discrete, Gaussian, and non-Gaussian networks of different dimensions (Section 3). Although NBNs do not yet win all comparisons, the benchmark con- firms their potential with regards to scalability and expressiveness, while also pointing to further assessment and development. The paper is organised as follows. Section 2 introduces NBNs, contrasts them with BNs, and define their properties. Section 3 presents the design of the benchmark and Section 4 its results. Section 5 reviews related works, and Section 6 concludes. 2 Neural Bayesian Networks A Neural Bayesian Network (NBN) is a BN whose conditional dis- tributions are learnable modules. The graph acts as the architecture of the model: each node hosts one module, its mechanism, which receives the values of the node’s parents and returns the distribu- tion of the node (Figure 1). Learning fits the mechanisms to data and inference runs them along the graph, both on batches, as in a neural network. Similarly to a BN, an NBN over variables푋= (푋 1 , . . .,푋 푑 )is a directed acyclic graph (DAG)G= (푉,퐸)with푉= 1, . . .,푑 and one conditional distribution푝 푖 (푥 푖 | 푥 Pa(푖) )per node, where Pa(푖)are the parents of푖inG[33,44]. The joint distribution is the product of the conditionals, 푝(푥 1 , . . .,푥 푑 )= 푑 Ö 푖=1 푝 푖 푥 푖 | 푥 Pa(푖) ,(1) so each variable is independent of its non-descendants given its parents. Differently from a BN, where a conditional is a conditional probability table (CPT) or a linear-Gaussian form, an NBN admits ar- bitrary representations of a conditional distribution: tabular entries, parametric and non-parametric expressions, and neural network based approximations, possibly mixed in the same graph. Additionally, in a BN, parameter learning estimates the condition- als one at a time using the log-likelihood Í 푛 Í 푖 log푝 푖 (푥 (푛) 푖 | 푥 (푛) Pa(푖) ) based on counting for a CPT and on the least squares method for a linear-Gaussian form [22]. In an NBN, instead, every conditional is differentiable, so the same log-likelihood can also be maximized by gradient descent on batches of data—for all the conditionals at once, and for those that admit no closed-form estimate, too. Also inference proceeds one variable at a time in a BN, by computing the posterior푝(푥 푄 | 푥 퐸 )of query variables푄given observed evidence 푥 퐸 using variable elimination at a cost exponential in the treewidth of the graph [13]. In an NBN, instead, each conditional is evaluated once for a whole batch of samples and evidence assignments, so that many queries are answered in a single pass over the graph, by exact elimination where the treewidth allows it and by sampling elsewhere. As we will show in the following, BNs are special case of NBNs when the closed-form methods are feasible. NBNis written in PyTorch [42]. The network is annn.Module. A user builds it from a DAG and a per-node mechanism assignment, callsfitonce, and thenquery. The benchmark, its configuration, and the complete results are released with the library. 1 2.1 Motivation BNs are a common “world model” for agents that decide under uncertainty: they are the transition models of factored Markov decision processes in planning and reinforcement learning [10,19, 30], and the diagnostic models used in troubleshooting, risk analysis, and healthcare [11, 34, 55]. Consider a decision-making agent that has to deal with a Markov decision process (MDP) formalised as(S,A,푃,푅,훾)with state푠 ∈ S, action푎 ∈ Aand transition kernel푃(푠 ′ | 푠,푎). The MDP treats the state as a single object, so푃is one distribution over the joint next state for every pair(푠,푎), and its size grows exponentially with the number of state variables. If we know which current variables each next-state variable depends on, we can move to a factored MDP [10, 30], which is much lighter (Figure 2): 푃(푠 ′ | 푠,푎)= 푛 Ö 푗=1 푝 푠 ′ 푗 | Pa(푠 ′ 푗 ) ,(2) where the parentsPa(푠 ′ 푗 )of each next-state variable are a few of the current state and action variables. This is a two-slice BN over 푋= (푠,푎,푠 ′ ). The same formulation applies to agents that face non-sequential scenarios, such as while doing diagnosis, where the BN has a single slice. The agent has to deal with two problems: P1 Learning. GivenGand a datasetDof푁samples, estimate every conditional of Eq. 1, for discrete and for continuous vari- ables, without fixing the shape of the conditional in advance. 1 https://github.com/Giovannibriglia/NeuralBayesianNetworks 2 Neural Bayesian Networks 푠 푎 푠 ′ MDP one kernel 푃(푠 ′ | 푠,푎) over the joint state known structure 푠 1 푠 2 푠 3 푎 푠 ′ 1 푠 ′ 2 푠 ′ 3 time푡 time푡+1 Factored MDP Î 푗 푝 푠 ′ 푗 | Pa(푠 ′ 푗 ) one small conditional per variable Figure 2: From an MDP to a factored MDP. When it is known which variables influence which, the transition kernel over the joint state becomes a two-slice BN with one conditional per next-state variable. P2Inference. Given the learned model, compute푝(푥 푄 | 푥 퐸 )for query variables푄and evidence푥 퐸 , for many evidence assign- ments at once. In this work we assume that the structure is given and focus on parameter learning and inference, while leaving structure learning and the sampling of trajectories for planning to future work. More precisely, we make the standard assumptions of parameter learning in BNs. (A1) The graphGis given and the data distribution factor- izes according to it. (A2) The data are complete, so every variable is observed in every sample, and the samples are independent and identically distributed. (A3) Each variable is a discrete variable with finitely many values or a real-valued scalar. All statements below are for a BN over 푑 variables and cover both cases. Classical BNs solve P1 and P2 only in part, because of three limitations of scalar conditionals evaluated one at a time. (L1) Rep- resentation, which limits P1. A CPT for푋 푖 has퐾 푖 Î 푗∈Pa(푖) 퐾 푗 entries, with퐾 푗 the number of values of푋 푗 , and is undefined for a continu- ous parent. Discretizing a parent introduces a binning parameter and loses information, and the linear-Gaussian alternative cannot represent multi-modal or heteroscedastic conditionals. Both fix the shape of the outcome uncertainty before seeing the data: a table needs the continuous outcome discretized, a Gaussian has one mode and a spread that does not depend on the parents. (L2) Sequential evaluation, which limits P1 and P2. A likelihood evaluation on one sample visits the푑nodes in order, and each visit is a scalar lookup or a scalar density. Nothing is shared across samples or across queries, so푆samples for퐵evidence assignments costΘ(푑·푆·퐵)scalar oper- ations. (L3) Exact inference, which limits P2. Variable elimination on a network of treewidth푤costs푂(푑 푐 푤 )with푐being the maximum cardinality [13], and it exists only for discrete and linear-Gaussian networks. Sampling avoids the exponent but inherits (L2). A solution to P1 therefore needs conditionals that are not re- stricted to tables and linear-Gaussian forms (L1) and that are fit on batches of data (L2). A solution to P2 needs an evaluation whose cost does not grow with푆 · 퐵(L2), with exact inference where the treewidth allows it and batched sampling elsewhere (L3). As described in Section 2.2, NBNs provide all of these properties. 2.2 Formal definition and Properties LetGbe a DAG over푋 1 , . . .,푋 푑 with푋 푖 ∈ X 푖 . An NBN replaces each conditional by a parametrised mapping 푝 휙 푖 푥 푖 | 푥 Pa(푖) =M 푖 푥 Pa(푖) ; 휙 푖 ,(3) whereM 푖 takes a batch푋 Pa(푖) ∈ R 퐵×|Pa(푖)| of parent configurations and returns, for every row, a distribution overX 푖 with a density log푝 휙 푖 , a sampler, and where it exists a reparametrised sampler. The joint is still Eq. 1: 푝 휙 (푋)= 푑 Ö 푖=1 푝 휙 푖 푋 푖 | 푋 Pa(푖) .(4) Two design choices are worth emphasising here. First,M 푖 returns the full conditional, not a sample or a point estimate, so learning and inference reach every mechanism through one interface. Second,휙 푖 are tensors on which autograd works, so a table, a linear-Gaussian form and a neural network are instances of the same class and can be mixed in one graph. All nodes at the same depth ofGare evaluated in one batched call. With this representation, NBNs support the following properties. Proposition 2.1 (Semantic eqivalence). Let푝be a BN onG and푝 휙 an NBN onGwith푝 휙 푖 (푥 푖 | 푥 Pa(푖) )= 푝 푖 (푥 푖 | 푥 Pa(푖) )for all푖 and all arguments. Then푝 휙 = 푝, and푝 휙 satisfies the global Markov property ofG. Proof. Both joints are the product in Eq. 1 of the same factors. The Markov property depends only on the factorization [33]. □ Every classical BN is an NBN with tabular or linear-Gaussian mechanisms, and everything that follows from the factorization (e.g., conditional independences) holds for NBNs. Proposition 2.2 (Differentiability). If everylog푝 휙 푖 is differ- entiable in휙 푖 , the log-likelihoodℓ(휙;푥)= Í 푖 log푝 휙 푖 (푥 푖 | 푥 Pa(푖) ) is differentiable in휙,∇ 휙 푖 ℓdepends only on the term of node푖, and the whole gradient is one backward pass. Proof. ℓ is a sum of terms depending on one휙 푖 ; reverse-mode differentiation of the sum evaluates all terms in a single pass. □ Parameter learning (P1) is therefore gradient descent on−ℓ, per node since the terms are independent (closed form where the family admits it), or jointly when a downstream differentiable objective couples the nodes. Proposition 2.3 (Cost of inference). Consider a query푝 휙 (푥 푄 | 푥 퐸 )answered by ancestral simulation with푆samples and퐵evidence assignments, and let퐴= Anc(푄 ∪ 퐸)be the active nodes. An NBN issues푂(depth(G[퐴]))sequential tensor operations, each of cost 푐(M,퐵·푆), so the wall-clock time is푂(depth(G[퐴])·푐(M,퐵푆)); a classical implementation issuesΘ(|퐴|· 퐵· 푆)scalar operations. For a fixed family, 푐(M,퐵푆) is constant in 퐵푆 until the device saturates. Proof.Nodes outside퐴are barren and pruned. Within퐴, all nodes at the same depth have their parents available and are evalu- ated in one call over the퐵·푆rows; there aredepth(G[퐴])such levels. The scalar count is the number of (node, sample, row) triplets.□ 3 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli The claim is not that NBNs change the dependence on the graph, but that for inference (P2) they remove dependence on푆and퐵up to device saturation, so one or a thousand queries cost the same—on a GPU the saturation point is far, on a CPU it is near. These gains do not come at the cost of sacrificing the usual BN properties that are behind their adoption: modularity is preserved since each mechanism has its own parameters and can be inspected, refit, or replaced without touching the others; the explicit structure stands still (in the DAG), which a designer can read and edit; and the “one model for every query” workflow is also there, since a fitted NBN answers 푝(푥 푄 | 푥 퐸 ) for any 푄 and 퐸 just like BNs. Other routes to expressive and scalable probabilistic models give up at least one of them: see Section 5. 2.3 Learning and Inference Methods Learning mechanisms differ in how much of the conditional is fixed in advance. Parametric mechanisms fix its form and have a closed- form fit: a probability table for a discrete child, by counting under a Dirichlet prior [22], and the linear-Gaussian formN(푊푥+푏,휎 2 )for a continuous one [35], by ridge regression [23]. The neural mecha- nism fixes the output form but not its dependence on the parents: a Mixture Density Network (MDN) [8] maps the parents with a multilayer perceptron to the parameters of a Gaussian mixture, fit by gradient descent on−ℓ. It is still parametric, but linear in the number of parents and multi-modal. Non-parametric mechanisms fix no form and grow with the data: KNN is a local mixture over the푘=⌈ √ 푁⌉training rows with the closest parents [53], and Flex- Code [25] expands the conditional density in a cosine basis, one regression per coefficient. Three inference engines answer푝 휙 (푥 푄 | 푥 퐸 ). Variable Elimina- tion (VE) [40,57] is exact and discrete-only: it drops the nodes that cannot affect the query and sums out the others by multiplying their tables in one batched tensor contraction, so퐵evidence rows cost one pass. Being exponential in the induced width, it is our exact reference, not a scalable method. Likelihood Weighting (LW) [15,51] simulates푆samples with the evidence clamped and aver- ages the query node, each sample weighted by the probability of the evidence; it needs only a sampler and a density per mechanism, so it works with every family. Amortized Variational Inference (AVI) [18,36] avoids the samples LW wastes when the evidence is un- likely: a recognition network푞 휓 (푥 −퐸 | 푥 퐸 ), trained once per model on the evidence lower bound [32], proposes samples close to the posterior for every query at constant cost; it falls back to LW when the proposal is worse than LW’s own Monte-Carlo noise. 3 Benchmark Design To the best of our knowledge, a fully reproducible benchmark for BNs has not been provided so far. Thebnlearnrepository [50] offers domains and their structures, and it is the usual source of test networks, but the benchmark configuration is unspecified: which queries to ask, how to obtain the ground truth, which metrics to report, or how to count a method that fails or does not finish is not established. In this way, results are hardly comparable thus advance- ments difficult to identify. Hence, we built one and release it with NBN. A benchmark is a configuration file that fixes the problems, the queries, the methods and the budgets, and the harness writes one row per problem, seed, method, query and metric, failures in- cluded. Its problems come from two sources: (i) the real domains ofbnlearn, and (i) synthetic domains that we generate. Synthetic networks are needed to define whatever network structure (varying nodes and edges configuration) with whatever underlying condi- tional probabilities (discrete, continuous gaussian and continuous non-gaussian), and evaluate sample efficiency, scalability and accu- racy across whatever BN. We use the benchmark onNBNto answer three research questions: RQ1.Where the classical BN assumptions hold, doesNBNmatch the classical BN libraries? RQ2. Where they do not hold, how much accuracy is gained? RQ3.Does the cost model of Proposition 2.3 show up as wall-clock time when many queries are answered at once? The numerical settings that follow are fields of the configuration file. We report here those needed to read the results, and all the others in Appendix A. Problems. A problem is a triplet of family, number of nodes푛, and seed. The family is the type of the conditionals of the true model: discrete, where every conditional is a CPT; linear-Gaussian; or non-Gaussian, where every conditional is a set of components forming a Gaussian mixture (that a linear-Gaussian model cannot represent correctly). The seed draws a sparse random DAG over the푛nodes and then the parameters of every conditional. From the true model we draw a training set (20 480 rows), a test set that supplies the evidence values of the queries, and a reference set that supplies the ground truth. Every problem is repeated over 5 seeds. The real networks are the 13 networks of Section 4.4. We define a cell as the combination of one method on one problem. Queries. A query has two components: a target variable, and an evidence set, which pairs each observed variable with its value. The evidence values are always copied from a row of the test set, so the evidence is a configuration that the true model produces. Three selectors choose targets and evidence variables, one per purpose: (1)Random (accuracy): the target and the evidence variables are drawn uniformly for the problem, and the query is asked for the evidence values of 128 test rows; it measures the typical error over evidence values. (2) Topological (real networks): targets are chosen by their role in the graph, among the most connected nodes (hubs), the nodes whose removal disconnects the graph (cut nodes), and the deep nodes with few descendants (terminals). Each target is asked in the prediction direction, with evidence among its ancestors, or in the diagnosis direction, with evidence among its descendants. (3) Heaviest (scalability and speed): only the top-ranked hubs, cut nodes, and terminals are queried, each in both directions, once with the evidence observed (full) and once with it marginalized (empty). Reference and metrics. For a discrete target the reference poste- rior ˆ 푝 is its empirical distribution among the reference samples that agree with the evidence, and the metric is the total variation (TV) between ˆ 푝and the answer푞,TV( ˆ 푝,푞)= 1 2 Í 푥 푄 | ˆ 푝(푥 푄 )−푞(푥 푄 )|[54]. For a continuous target the reference is a sample from푝(푥 푄 | 푥 퐸 ) obtained by likelihood weighting on the true model, and the metric 4 Neural Bayesian Networks is the 1-Wasserstein distance (W1) between the empirical distribu- tion ˆ 휇of this sample and the empirical distribution휈of the samples returned by the method,푊 1 ( ˆ 휇,휈)= inf 훾∈Π( ˆ 휇,휈) E (푥,푥 ′ )∼훾 |푥 − 푥 ′ |, whereΠ( ˆ 휇,휈)is the set of couplings of ˆ 휇and휈[16]. In both cases the reference is a finite sample and not the exact posterior, so the measured accuracy has a floor that all methods share. On the dis- crete networks of the accuracy scenario it is a TV of about 0.06 at 10 nodes and 0.08 from 50 nodes on, which are the values exact engines reach in Figure 3. We use a sampled reference because an exact one is not available for every family and size: variable elimination on the true model is exponential in the treewidth, and non-Gaussian conditionals have no closed-form posterior. Com- puting exact references when tractable is left for future work. The main text reports this accuracy metric (TV or W1) for the accuracy, scalability, and real-network scenarios, and the time per query for the speed batching scenario, measured as the wall-clock time of one call divided by the number of answered queries. Complete fit and query times are in Appendix B. Budgets and failures. Every cell has a budget of 1200 s for the fit, a separate budget of 1200 s for its queries, and a memory cap: 80% of RAM (∼100 GB on our server, that is a more favourable scenarios for cpu-based libraries e.g.pgmpy.), and the 12 GB of the GPU. It ends with a status:ok,timeout,oomornot_supported. In the figures a bar is the interquartile mean (IQM) over the seeds a method solved, as [2] suggests;푓/푛is printed when푓of푛seeds failed, and a missing bar indicates the failure status. Scenarios. Each scenario is one configuration. Accuracy fixes the query load (128 random queries) and varies the family and the size, 푛 ∈ 10,50,100,500,1000: on the discrete and the linear-Gaussian family the classical model is correct, which isolates RQ1, and on the non-Gaussian family it is misspecified, which isolates RQ2. Scalability pushes the network size to푛 ∈ 10,10 2 ,10 3 ,10 4 on the heaviest queries and records where each method stops. Batching fixes the network size (푛=50) and varies only the number퐵 ∈ 1,128,512,1024of queries answered in one call, which is the variable of Proposition 2.3 (RQ3). Real networks repeat the accuracy scenario on published structures with topological queries. Learning scores the fitted conditionals without any query and addresses P1. Methods and hardware. The comparison baselines arepgmpywith maximum-likelihood and with Bayesian CPTs, both queried by its variable elimination;pgmpylinear-Gaussian, which conditions the joint Gaussian in closed form;pomegranate, which is discrete only and uses its built-in loopy belief propagation; andPyro, where the BN is written as a program with smoothed tables or linear- Gaussian conditionals and queried by importance sampling from the prior with 50 particles (the largest number its sequential sam- pler affords within the budgets). On continuous networkspgmpy andPyroalways fit linear-Gaussian conditionals, so on the non- Gaussian family they are the misspecified classical model. NBNruns every applicable pair of a mechanism and an engine, with 2048 particles for the sampling engines (16 in the batching scenario, which measures time only). To keep the figures readable, each panel shows all applicable external baselines plus the bestNBN methods of the family, selected by accuracy. Appendix B reports every method on every scenario. All runs use one server with a 20-core CPU, 126 GB of RAM and one GPU with 12 GB of memory; NBNandpomegranaterun on the GPU,pgmpyon the CPU, and Pyro on the CPU except in the batching scenario. 4 Experimental Results Section 2 makes two claims: an NBN keeps what a BN offers, while enhancing expressiveness of conditionals and scaling up inference performance. The results support both, within limits that the bench- mark makes visible. In terms of our research questions: RQ1.No loss where the classical model is right (Proposition 2.1). On discrete networks, exact inference inNBNreturns the same TV aspgmpyat every size, and on real networks it matches or improves it on seven out of eight; learning is identical as well. The exception is inference on linear-Gaussian networks, where the closed form ofpgmpyis more accurate than our sampled posterior (W1 0.039 against 0.074 at 1000 nodes). RQ2.Expressiveness pays where it is wrong. On non-Gaussian net- works the MDN has between one third and 60% of the W1 of a linear-Gaussian model, and achieves better held-out likelihood: where the outcome has several modes or a state- dependent spread, the linear-Gaussian family misstates the uncertainty while the MDN learns its shape. RQ3.The cost model holds on the GPU (Proposition 2.3). With the sampling engines the time per query falls by two orders of magnitude between one query and a batch of 1024, down to tens of microseconds. Whilepgmpyachieves top accuracy wherever one of its two fami- lies is the right model, it cannot represent non-Gaussian condition- als and answers one query at a time.pomegranateleverages the GPU, but its loopy belief propagation returns undefined posteriors or exceeds budget on all networks but the smallest.Pyroruns a generic sequential importance sampler, which the budget limits to 50 particles: it has the highest error on all synthetic families and times out from 500 nodes onwards. The performance ofNBNcan still be improved, and the benchmark shows where. The sampling engines pay a Monte-Carlo error where a closed form exists, the MDN exceeds the fit budget between 100 and 500 nodes, and the nearest-neighbour mechanism exhausts the GPU memory at 500. Below,NBNmethods are named by their mechanisms and engines: cat+VE is categorical table queried by variable elimination; likewise for linear-Gaussian (LG) and MDN mechanisms, and for likelihood weighting (LW) and amortized variational inference (AVI). 4.1 Inference accuracy (RQ1, RQ2) Figure 3 shows the accuracy of the 128 queries for each family. Discrete. cat+VE andpgmpyreturn the same TV at every size because they run exact inference on the same maximum-likelihood tables: this is Proposition 2.1, measured. Their common value is the floor of the reference (Section 3). The sampling engines ofNBNstay within 0.01 of it, so on discrete networks approximation costs little. Linear-Gaussian. This is the one family where a baseline is ahead: pgmpyconditions the joint Gaussian in closed form, whileNBNsam- ples. The gap grows with the size because, in deeper networks, fewer prior samples agree with the evidence: the median effec- tive sample size falls from all 2048 particles at 10 nodes to 13% of 5 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli discrete pgmpy Bayes+VE pgmpy MLE+VE pomegranate pyro IS NBN cat+VE NBN cat+LW 10501005001000 n (nodes) 0.05 0.1 0.2 TV ↓ 1/5 2/5 nannantimeouttimeouttimeouttimeout continuous (Gaussian) pgmpy LG pyro IS NBN LG+LW NBN kNN+LW 10501005001000 n (nodes) 0.05 0.1 0.2 0.5 W1 ↓ timeouttimeoutoomoom continuous (non-Gaussian) pgmpy LG pyro IS NBN MDN+LW NBN kNN+LW NBN LG+LW 10501005001000 n (nodes) 0.05 0.1 0.2 timeouttimeouttimeouttimeoutoomoom Figure 3: Inference accuracy against network size on synthetic networks (TV for discrete, W1 for continuous, lower is better). All applicable external baselines are shown, plus the selectedNBNmethods per family. A missing bar is a method that failed on every seed, with its status printed; 푓/푛 marks 푓 failed seeds out of 푛. them at 1000. The error is due to the engine, not the model, and a closed-form engine for all-Gaussian networks would remove it. Non-Gaussian. Here the model is decisive: every linear-Gaussian method has the same error, inNBNor inpgmpyand whatever the engine, because the family cannot represent the data. The MDN has one third to 60% of that error and kNN one half to three quarters. The price is paid at fit time (Table 13 in Appendix B.1): the MDN exceeds the fit budget at 500 nodes, where kNN runs out of memory, so beyond 100 nodes only the misspecified models remain. Table 1: Accuracy on the heaviest queries with full evidence against network size, up to 10 000 nodes (IQM over seeds, lower is better).(푓/푛)marks푓failed seeds out of푛. Bold: fewest failed seeds, then best value. Times in Appendix B.2. Method푛=10 푛=100 푛=1000 푛=10000 Discrete, TV pgmpy Bayes+VE0.04910.05780.0645 0.0323 (1/5) pgmpy MLE+VE0.04380.04060.05820.0484 pomegranate0.0276 (2/5)nannantimeout Pyro IS0.1200.1880.1510.141 (1/5) NBN cat+VE0.04350.0476 0.0516 0.0376 (2/5) NBN cat+LW0.04950.06000.06190.0586 Linear-Gaussian, W1 pgmpy LG0.03030.0553 0.03780.0466 Pyro IS0.1380.3300.3380.226 (1/5) NBN LG+LW0.02910.05260.06520.0560 NBN MDN+LW0.03260.0534 timeouttimeout kNN+LW0.05570.102oomoom Non-Gaussian, W1 pgmpy LG0.1150.142 0.08300.101 Pyro IS0.2320.1930.1750.209 (1/5) NBN LG+LW0.1250.1370.09170.0992 NBN MDN+LW0.0504 0.0414 (1/5) timeouttimeout NBN kNN+LW0.07550.0779oomoom 4.2 Scalability Table 1 pushes the heaviest queries to 10 000 nodes. Accuracy stays in the same range at every size; what this scenario reveals is why each method stops and how long a query takes. Exact inference is bounded by the elimination itself.pgmpycompletes every seed but needs tens of seconds per query from 1000 nodes; cat+VE answers in 0.2 to 0.4 s, two orders of magnitude faster, and exceeds the 12 GB of the GPU on two of five seeds at 10 000 nodes. The sampling engines have no such wall: cat+LW and LG+LW finish every seed at 10 000 nodes, where LG+LW is also faster than the closed form ofpgmpy (2.6 s against 17.5 s per query), being the fastest up to 1000 nodes. The flexible mechanisms stop earlier due to the cost of learning: at 1000 nodes the MDN exceeds the fit budget and kNN the memory. 4.3 Batched queries (RQ3) Proposition 2.3 postulates that the cost of a batch does not grow with its size until the device saturates, so the time per query should fall as 1/퐵. We fix 50-node networks and ask each heaviest query for the evidence values of퐵different rows (Figure 4). The sampling engines follow the prediction closely up to퐵=128, where a query costs about one hundredth of a single one, and then flattens as the GPU saturates: at퐵=1024 the gain is 456×for cat+AVI and 227× for LG+AVI, down to 0.09 and 0.05 ms per query. Exact cat+VE is the fastest method on one query (2.8 ms) but gains only 4×: every factor of the elimination carries the batch axis, so its work grows with퐵and batching removes only the per-call overhead. Among the baselines onlypomegranatebatches, and at퐵=1024 it is 60 times slower than cat+AVI;pgmpyanswers one query at a time. The gain is in the batch, not in the latency of one query: on a single continuous query the closed form ofpgmpytakes 2.6 ms andNBN 12 ms. An agent that asks one question at a time is better served by the classical code; one that scores many candidate actions or particles at once is the case NBNs are built for. 6 Neural Bayesian Networks discrete pgmpy Bayes+VE pgmpy MLE+VE pomegranate pyro IS NBN cat+VE NBN cat+AVI 11285121024 B (queries per batch) 10 −4 10 −3 10 −2 10 −1 Time per query (s) 4/5 1/5 timeouttimeout continuous (Gaussian) pgmpy LG pyro IS NBN LG+AVI NBN FlexCode+LW 11285121024 B (queries per batch) 10 −4 10 −3 10 −2 timeout continuous (non-Gaussian) pgmpy LG pyro IS NBN LG+AVI NBN FlexCode+LW 11285121024 B (queries per batch) 10 −4 10 −3 10 −2 timeout Figure 4: Time per query vs. batch size퐵(log scale) on 50-node networks, full evidence. Non-batching methods are at퐵=1 only. 4.4 Real networks Table 2 repeats the accuracy comparison on thirteenbnlearnnet- works: eight discrete ones, from 8 to 1041 nodes, two Gaussian ones, and three conditional linear-Gaussian (CLG) ones, which mix discrete and continuous variables. The conclusions of the synthetic scenarios carry over. On the discrete networks cat+VE has the lowest TV on five, including the two largest, and is within 0.001 of the best on two more;pgmpyis better onmildewonly. Amortized inference suffers more here than on synthetic networks (up to four times the TV of exact inference), so the exact engine is the one to use whenever the treewidth allows it. On the Gaussian networksNBNis again close to the closed form of pgmpy. The CLG networks are where coverage matters:pgmpyand Table 2: Accuracy onbnlearnnetworks, grouped by kind and ordered by node count푛(IQM over seeds; TV for discrete, W1 for Gaussian and CLG, lower is better). Each block names the method of each library; – denotes no applicable engine. pgmpywith Bayesian estimates and the otherNBNengines are in Appendix B.4. Network (푛) pgmpy pomegr. Pyro NBN NBN Discrete, TVMLE+VEBPIScat+VE cat+AVI asia (8)0.003430.003800.03360.003540.00478 sachs (11)0.01280.01220.08180.01200.0162 mildew (35)0.0325nan0.2870.04650.0703 barley (48)0.0275nan0.2150.02400.0401 hailfinder (56)0.02030.02080.1360.02060.0266 diabetes (413)0.0393 (1/5)nan timeout0.03430.128 link (724)0.00197 timeout timeout 0.001740.00612 munin (1041)0.00464 timeout timeout 0.004010.0105 Gaussian, W1LG–ISLG+AVI kNN+LW ecoli70 (46)0.0418–0.2450.05410.111 arth150 (107)0.0210–0.1070.02120.0418 CLG, W1–IS MDN+LW kNN+LW healthcare (7)–10662523.4 sangiovese (15)–0.05830.02420.0280 mehra (24)– 3.74·10 7 1.89·10 6 2.41·10 6 pomegranatehave no applicable engine, and the bestNBNmecha- nism is belowPyroon all three. Which mechanism is best, however, depends on the network: the MDN wins onsangioveseandmehra and has six times the error ofPyroonhealthcare, where kNN wins. The choice of the mechanism is a modelling decision, and the benchmark is the tool to make it. 4.5 Parameter learning and sample efficiency The learning scenarios score the fitted conditionals with no query: by the parameter-recovery TV on discrete networks, the distance be- tween learned and true rows of the CPTs, and by held-out negative log-likelihood (NLL) per node on continuous ones (Appendix A). Table 3: Sample efficiency on 50-node networks against the number of training rows푛 tr (IQM over seeds, lower is better): parameter-recovery TV on discrete networks, held-out NLL per node on continuous ones. timeout: the 1200 s fit budget was exceeded on every seed. Fit times and the network-size sweep are in Appendices B.5 and B.6. Method푛 tr =4096 푛 tr =16384 푛 tr =32768 푛 tr =65536 Discrete, parameter-recovery TV pgmpy Bayes0.09770.04740.03350.0236 pgmpy MLE0.09790.04740.03360.0236 pomegranate0.09790.04740.03360.0236 Pyro0.08470.04530.03280.0233 NBN cat0.09790.04740.03360.0236 NBN cat (Bayes) 0.08470.04530.03280.0233 Gaussian, NLL per node pgmpy LG0.9430.9420.9440.945 NBN LG0.935 (1/5)0.9420.9440.945 NBN MDN0.9530.9460.958 (1/5)timeout NBN kNN1.010.9880.9810.966 (1/5) Non-Gaussian, NLL per node pgmpy LG1.261.271.271.27 NBN LG1.261.271.271.27 NBN MDN1.201.191.19timeout NBN kNN1.251.231.221.21 7 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli When the model is correct,NBNlearns what competitors do, as the estimator is the same: the tabular mechanism equals the max-likelihood tables ofpgmpyandpomegranateat all sizes, its smoothed variant equals the tables ofPyro, and the linear-Gaussian mechanism equals the least squares ofpgmpy(Table 3; Appen- dix B.5). It is also faster: a 1000-node network is fitted in 2.1 s (tables) and 2.6 s (linear-Gaussian), against 9.8 s and 4.0 s forpgmpy. The accuracy drop when data is Gaussian is negligible: the NLL of the MDN is within 0.004 of the exact estimator. When the model is wrong, more data does not repair it (Table 3). On 50-node non-Gaussian networks the NLL of the linear-Gaussian models stays at 1.26–1.27 from 4096 to 65 536 rows, while the MDN is at 1.20 with 4096 rows already and kNN improves steadily to 1.21. Flexibility is thus sample-efficient but not compute-efficient: the fit time of the MDN grows linearly with the rows and nodes, exceeding the budget at 65 536 rows and from 100 nodes (Appendix B.6). 5 Related Works NBNs extend classical BNs in two directions, the expressiveness of the conditionals and the efficiency and scale of learning and inference; we review the works that pursue each of them. Expressiveness. Kernel density estimation was used for continu- ous conditionals in Bayesian classifiers [28], and non-parametric BNs based on kernels or copulas followed [21,24,38]. Semiparamet- ric BNs [4,5] mix linear-Gaussian and kernel conditionals in one network; surveys of inference with continuous variables are [39,48]. Neural structural causal models replace each structural equation with a neural density, typically a normalizing flow [17,26,31,43,56]. Closest to ours in name, the Neural Bayesian Network Understudy [46] trains a single Neural network to approximate the conditional probabilities of a given BN, with its independences encoded through training. These works fix one family of conditionals and one in- ference procedure. At the other extreme, deep generative models such as normalizing flows [41,47] and neural world models [20] are expressive but hide the dependencies in the weights of one network: there is no local conditional to inspect or replace, and the posterior given arbitrary evidence is not directly available. An NBN changes the representation rather than adding a family. With tensor- valued mechanisms the BN is kept exactly (Proposition 2.1); tabular, linear-Gaussian, neural and non-parametric [25] mechanisms can be mixed in one graph and trained by gradient, also jointly with a downstream objective (Proposition 2.2); and any mechanism that can sample and evaluate a density is queried by the same batched engine, with no family-specific algorithm (Proposition 2.3). Efficiency and scalability. Exact inference has been parallelized on multi-core CPUs [27] and on GPUs [58], and written as tensor contractions [40], but only for discrete networks and at a cost that stays exponential in the treewidth. Probabilistic circuits [12,45] keep inference tractable and batched, but trade the per-variable conditionals for latent sum and product units. When exact infer- ence is too expensive, a BN is queried by sampling. Likelihood weighting (Section 2.3) is simple and general, but its samples are drawn without looking at the evidence, so most of them are wasted when the evidence is unlikely. Amortized inference [18, 36] trains a neural network to propose samples that already agree with the Table 4: Bayesian-network and probabilistic-modeling li- braries. Disc.: discrete; LG: linear-Gaussian; NG: continuous non-Gaussian. Batch: the API accepts many queries in one call. GPU: BN queries run in parallel on the accelerator. Diff.: differentiable conditionals.do: native interventions.✓sup- ported,• partial,✗ absent. Library FamiliesComputeCapabilities Disc. LG NGBatch GPUDiff. do pgmpy [3]✓✗✓✗ • pyAgrum [14]✓✗✓✗ bnlearn [50] ✓✗✓✗ CausalNex [6]✓✗✓✗✓ pomegranate [49]✓ •✓ •✗ PyMC [1]✓•✓•✗ Pyro [7]✓ • DoWhy [52]✓ •✗✓ NBN (ours)✓ evidence and reuses it for every later query. We do not propose a new estimator: NBNs change how these estimators are executed. All the nodes at the same depth of the graph are evaluated in one tensor operation over every particle and every evidence assignment, so the cost of a query stops growing with the number of samples and of queries (Proposition 2.3). Libraries. Table 4 summarizes NBNs competitors. Classical BN libraries (pgmpy[3],pyAgrum[14],bnlearn[50],CausalNex[6]) implement tabular or linear-Gaussian conditionals with exact in- ference; they run on the CPU and accept several queries in one call but execute them one by one, which are the limitations (L1) and (L2) of Section 2.1.pomegranate[49] is PyTorch-based and runs on the GPU, but its BN inference is limited to discrete networks and to loopy belief propagation, which is not guaranteed to converge on graphs with loops. Probabilistic programming systems (Pyro[7], PyMC[1]) support any distribution and GPU execution, but they do not expose the BN as an object: the graph, the per-node condi- tionals, and the inference procedure over them have to be written by hand for every model, and a generic sampler cannot eliminate variables exactly, prune the nodes that a query does not need, or evaluate a level of the graph at once.DoWhy[9,52] estimates the causal effect of a treatment on an outcome and does not answer posterior queries over arbitrary variables. To the best of our knowledge, no existing library combines an explicit DAG with differentiable conditionals of arbitrary family and GPU-batched learning and inference, which is the gapNBN fills. Our benchmark compares against one representative per class: pgmpyfor the classical libraries,pomegranateas the PyTorch-based one, and Pyro as the probabilistic programming system. 6 Conclusions and Future Works Neural Bayesian Networks (NBNs) are an extension to Bayesian networks that combine explicit graphical structure with learn- able, tensor-valued conditional mechanisms. NBNs preserve the probabilistic factorization and modularity of BNs while bringing heterogeneous conditional models into a common framework for learning and GPU-batched inference. The accompanying PyTorch 8 Neural Bayesian Networks library and reproducible benchmark provide a basis for evaluat- ing the trade-offs between expressiveness, accuracy, and computa- tional cost. NBNs provide to autonomous agents a tool for decision- making under uncertainty that combines structured world models with the flexibility and computational tools of modern machine learning. Explicit dependencies and separately replaceable mecha- nisms support adaptation of individual model components, while batched inference serve agents that evaluate many candidate ac- tions or hypotheses. Evaluation of NBNs within agent decision making tasks (e.g. causal reinforcement learning), removing the assumption of having the graph available, evaluating sampling, and considering CPU execution along GPU are amongst the future works we are considering. References [1]Oriol Abril-Pla, Virgile Andreani, Colin Carroll, Larry Dong, Christopher J Fon- nesbeck, Maxim Kochurov, Ravin Kumar, Junpeng Lao, Christian C Luhmann, Osvaldo A Martin, et al.2023. PyMC: a modern, and comprehensive probabilistic programming framework in Python. PeerJ Computer Science 9 (2023), e1516. [2] Rishabh Agarwal, Max Schwarzer, Pablo Samuel Castro, Aaron C Courville, and Marc Bellemare. 2021. Deep reinforcement learning at the edge of the statistical precipice. Advances in neural information processing systems 34 (2021), 29304–29320. [3]Ankur Ankan and Johannes Textor. 2024. pgmpy: A python toolkit for bayesian networks. Journal of Machine Learning Research 25, 265 (2024), 1–8. [4]David Atienza, Concha Bielza, and Pedro Larrañaga. 2022. Semiparametric Bayesian Networks. Information Sciences 584 (2022), 564–582. doi:10.1016/j.ins. 2021.10.074 [5]David Atienza, Pedro Larrañaga, and Concha Bielza. 2022. Hybrid Semiparamet- ric Bayesian Networks. In TEST, Vol. 31. 299–327. [6]Paul Beaumont, Ben Horsburgh, Philip Pilgerstorfer, Angel Droth, Richard Oen- taryo, Steven Ler, Hiep Nguyen, Gabriel Azevedo Ferreira, Zain Patel, and Wesley Leong. 2021. CausalNex. URL https://github. com/quantumblacklab s/causalnex (2021). [7]Eli Bingham, Jonathan P Chen, Martin Jankowiak, Fritz Obermeyer, Neeraj Pradhan, Theofanis Karaletsos, Rohit Singh, Paul Szerlip, Paul Horsfall, and Noah D Goodman. 2019. Pyro: Deep universal probabilistic programming. Journal of machine learning research 20, 28 (2019), 1–6. [8] Christopher M. Bishop. 1994. Mixture Density Networks. Technical Report NCRG/94/004. Aston University. [9] Patrick Blöbaum, Peter Götz, Kailash Budhathoki, Atalanti A Mastakouri, and Dominik Janzing. 2024. DoWhy-GCM: An extension of DoWhy for causal infer- ence in graphical causal models. Journal of Machine Learning Research 25, 147 (2024), 1–7. [10] Craig Boutilier, Richard Dearden, Moisés Goldszmidt, et al.1995. Exploiting structure in policy construction. In IJCAI, Vol. 14. 1104–1113. [11]John S. Breese and David Heckerman. 1996. Decision-Theoretic Troubleshooting: A Framework for Repair and Experiment. In UAI ’96: Proceedings of the Twelfth Annual Conference on Uncertainty in Artificial Intelligence, Reed College, Portland, Oregon, USA, August 1-4, 1996, Eric Horvitz and Finn Verner Jensen (Eds.). Morgan Kaufmann, 124–132. https://dslpitt.org/uai/displayArticleDetails.jsp?mmnu=1& smnu=2&article_id=360&proceeding_id=12 [12]YooJung Choi, Antonio Vergari, and Guy Van den Broeck. 2020. Probabilistic Circuits: A Unifying Framework for Tractable Probabilistic Models. Technical Report. University of California, Los Angeles. http://starai.cs.ucla.edu/papers/ ProbCirc20.pdf [13]Gregory F. Cooper. 1990. The Computational Complexity of Probabilistic In- ference Using Bayesian Belief Networks. Artif. Intell. 42, 2-3 (1990), 393–405. doi:10.1016/0004-3702(90)90060-D [14]Gaspard Ducamp, Christophe Gonzales, and Pierre-Henri Wuillemin. 2020. aGrUM/pyAgrum: a toolbox to build models and algorithms for Probabilistic Graphical Models in Python. In International Conference on Probabilistic Graphical Models. PMLR. [15]Robert Fung and Kuo-Chu Chang. 1989. Weighing and Integrating Evidence for Stochastic Simulation in Bayesian Networks. In Proceedings of the 5th Conference on Uncertainty in Artificial Intelligence (UAI). 209–219. [16]Peyré Gabriel and Cuturi Marco. 2019. Computational optimal transport with applications to data sciences. Foundations and trends® in machine learning 11, 5-6 (2019), 355–607. [17]Tomas Geffner, Javier Antorán, Adam Foster, Wenbo Gong, Chao Ma, Emre Kiciman, Amit Sharma, Angus Lamb, Martin Kukla, Nick Pawlowski, Miltiadis Allamanis, and Cheng Zhang. 2022. Deep End-to-end Causal Inference. arXiv preprint arXiv:2202.02195 (2022). [18]Samuel J. Gershman and Noah D. Goodman. 2014. Amortized Inference in Probabilistic Reasoning. Proceedings of the Annual Meeting of the Cognitive Science Society 36 (2014). [19]Carlos Guestrin, Daphne Koller, Ronald Parr, and Shobha Venkataraman. 2003. Efficient Solution Algorithms for Factored MDPs. Journal of Artificial Intelligence Research 19 (2003), 399–468. [20]David Ha and Jürgen Schmidhuber. 2018. World models. arXiv preprint arXiv:1803.10122 2, 3 (2018). [21]Anca Hanea, Oswaldo Morales Napoles, and Daniel Ababei. 2015.Non- Parametric Bayesian Networks: Improving Theory and Reviewing Applications. Reliability Engineering & System Safety 144 (2015), 265–284. [22] David Heckerman. 1998. A tutorial on learning with Bayesian networks. Learning in graphical models (1998), 301–354. [23]Arthur E. Hoerl and Robert W. Kennard. 1970. Ridge Regression: Biased Estima- tion for Nonorthogonal Problems. Technometrics 12, 1 (1970), 55–67. [24]Katja Ickstadt, Björn Bornkamp, Marco Grzegorczyk, Jakob Wieczorek, M. Rashid Sheriff, Hernán E. Grecco, and Eli Zamir. 2010. Nonparametric Bayesian Net- works. Bayesian Statistics 9 (2010). [25]Rafael Izbicki and Ann B. Lee. 2017. Converting High-Dimensional Regression to High-Dimensional Conditional Density Estimation. Electronic Journal of Statistics 11, 2 (2017), 2800–2831. [26]Adrián Javaloy, Pablo Sánchez-Martín, and Isabel Valera. 2023. Causal Normaliz- ing Flows: From Theory to Practice. In Advances in Neural Information Processing Systems, Vol. 36. [27] Jiantong Jiang, Zeyi Wen, Atif Mansoor, and Ajmal Mian. 2023. POSTER: Fast Parallel Exact Inference on Bayesian Networks. In Proceedings of the 28th ACM SIGPLAN Annual Symposium on Principles and Practice of Parallel Programming (PPoPP). doi:10.1145/3572848.3577476 [28] George H. John and Pat Langley. 1995. Estimating Continuous Distributions in Bayesian Classifiers. In Proceedings of the 11th Conference on Uncertainty in Artificial Intelligence (UAI). 338–345. [29]Leslie Pack Kaelbling, Michael L Littman, and Anthony R Cassandra. 1998. Plan- ning and acting in partially observable stochastic domains. Artificial intelligence 101, 1-2 (1998), 99–134. [30]Michael Kearns and Daphne Koller. 1999. Efficient Reinforcement Learning in Factored MDPs. In Proceedings of the 16th International Joint Conference on Artificial Intelligence (IJCAI). 740–747. [31] Ilyes Khemakhem, Ricardo Monti, Robert Leech, and Aapo Hyvärinen. 2021. Causal Autoregressive Flows. In Proceedings of the 24th International Conference on Artificial Intelligence and Statistics (AISTATS). 3520–3528. [32]Diederik P. Kingma and Max Welling. 2014. Auto-Encoding Variational Bayes. In International Conference on Learning Representations (ICLR). [33] Daphne Koller and Nir Friedman. 2009. Probabilistic Graphical Models: Principles and Techniques. MIT Press. [34] Evangelia Kyrimi, Scott McLachlan, Kudakwashe Dube, Mariana R. Neves, Ali Fahmi, and Norman Fenton. 2021. A Comprehensive Scoping Review of Bayesian Networks in Healthcare: Past, Present and Future. Artificial Intelligence in Medicine 117 (2021), 102108. [35]Steffen L. Lauritzen and Nanny Wermuth. 1989. Graphical Models for Associa- tions between Variables, some of which are Qualitative and some Quantitative. The Annals of Statistics 17, 1 (1989), 31–57. [36]Tuan Anh Le, Atılım Güneş Baydin, and Frank Wood. 2017. Inference Compilation and Universal Probabilistic Programming. In Proceedings of the 20th International Conference on Artificial Intelligence and Statistics (AISTATS). 1338–1348. [37] Sergey Levine, Aviral Kumar, George Tucker, and Justin Fu. 2020. Offline Re- inforcement Learning: Tutorial, Review, and Perspectives on Open Problems. arXiv preprint arXiv:2005.01643 (2020). [38]Dimitris Margaritis. 2005. Distribution-Free Learning of Bayesian Network Structure in Continuous Domains. Proceedings of the AAAI Conference on Artificial Intelligence (2005). [39]Joaquín Mori, Josep Maria Mateo-Sanz, and Antonio Salmerón. 2016. Inference in Hybrid Bayesian Networks: A Survey. International Journal of Approximate Reasoning (2016). [40]Fritz Obermeyer, Eli Bingham, Martin Jankowiak, Neeraj Pradhan, Justin Chiu, Alexander Rush, and Noah Goodman. 2019. Tensor Variable Elimination for Plated Factor Graphs. In Proceedings of the 36th International Conference on Machine Learning (ICML). [41] George Papamakarios, Theo Pavlakou, and Iain Murray. 2017. Masked Au- toregressive Flow for Density Estimation. In Advances in Neural Information Processing Systems, Vol. 30. [42]Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al. 2019. PyTorch: An Imperative Style, High-Performance Deep Learning Library. Advances in Neural Information Processing Systems 32 (2019). [43]Nick Pawlowski, Daniel C. Castro, and Ben Glocker. 2020. Deep Structural Causal Models for Tractable Counterfactual Inference. In Advances in Neural Information Processing Systems, Vol. 33. 9 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli [44]Judea Pearl. 1988. Probabilistic Reasoning in Intelligent Systems: Networks of Plausible Inference. Morgan Kaufmann. [45]Robert Peharz, Steven Lang, Antonio Vergari, Karl Stelzner, Alejandro Molina, Martin Trapp, Guy Van den Broeck, Kristian Kersting, and Zoubin Ghahramani. 2020. Einsum Networks: Fast and Scalable Learning of Tractable Probabilistic Circuits. In Proceedings of the 37th International Conference on Machine Learning (ICML) (Proceedings of Machine Learning Research, Vol. 119). 7563–7574. [46]Paloma Rabaey, Cedric De Boom, and Thomas Demeester. 2022. Neural Bayesian Network Understudy. arXiv preprint arXiv:2211.08243 (2022). [47]Danilo Jimenez Rezende and Shakir Mohamed. 2015. Variational Inference with Normalizing Flows. In Proceedings of the 32nd International Conference on Machine Learning (ICML) (Proceedings of Machine Learning Research, Vol. 37). 1530–1538. [48]Antonio Salmerón, Rafael Rumí, Helge Langseth, Thomas D. Nielsen, and An- ders L. Madsen. 2018. A Review of Inference Algorithms for Hybrid Bayesian Networks. Journal of Artificial Intelligence Research 62 (2018), 799–828. [49]Jacob Schreiber. 2018. Pomegranate: fast and flexible probabilistic modeling in python. Journal of Machine Learning Research 18, 164 (2018), 1–6. [50] Marco Scutari. 2010. Learning Bayesian Networks with the bnlearn R Package. Journal of Statistical Software 35, 3 (2010), 1–22. [51]Ross D. Shachter and Mark A. Peot. 1989. Simulation Approaches to General Probabilistic Inference on Belief Networks. In Proceedings of the 5th Conference on Uncertainty in Artificial Intelligence (UAI). 221–231. [52]Amit Sharma and Emre Kiciman. 2020. DoWhy: An end-to-end library for causal inference. arXiv preprint arXiv:2011.04216 (2020). [53] Charles J. Stone. 1977. Consistent Nonparametric Regression. The Annals of Statistics 5, 4 (1977), 595–620. [54]Alexandre B. Tsybakov. 2009. Introduction to Nonparametric Estimation. Springer. doi:10.1007/B13794 [55]Philippe Weber, Gabriela Medina-Oliva, Christophe Simon, and Benoît Iung. 2012. Overview on Bayesian Networks Applications for Dependability, Risk Analysis and Maintenance Areas. Engineering Applications of Artificial Intelligence 25, 4 (2012), 671–682. [56]Kevin Xia, Kai-Zhan Lee, Yoshua Bengio, and Elias Bareinboim. 2021. The Causal- Neural Connection: Expressiveness, Learnability, and Inference. In Advances in Neural Information Processing Systems, Vol. 34. [57]Nevin Lianwen Zhang and David Poole. 1994. A Simple Approach to Bayesian Network Computations. In Proceedings of the 10th Canadian Conference on Artifi- cial Intelligence. 171–178. [58] Lu Zheng, Ole Mengshoel, and Jike Chong. 2011. Belief Propagation by Mes- sage Passing in Junction Trees: Computing Each Message Faster Using GPU Parallelization. In Proceedings of the 27th Conference on Uncertainty in Artificial Intelligence (UAI). 822–830. 10 Neural Bayesian Networks A Benchmark Design – Extended Version This appendix gives the settings of the benchmark that Section 3 leaves out. All of them are fields of the released configuration files or defaults of the released code, and can be changed. Synthetic generator. The seed draws a DAG over푛ordered nodes, where each earlier node becomes a parent of a later one with probability 0.2 and a node keeps at most four parents. The cap bounds the size of a conditional as푛grows: with four values per variable, a CPT has at most 4 5 =1024 entries. In the discrete family every variable has four values and every row of a CPT is the softmax of standard-normal logits. In the linear-Gaussian family every node has weights of variance 1/|Pa|, a standard-normal bias and a noise scale drawn uniformly in [0.3,1]. In the non-Gaussian family every node is a three-component Gaussian mixture: the component means are linear in the parents, with weights of standard deviation 0.3/ √︁ |Pa|and offsets 0.3(푧 1 +0.7|푧 2 |)for standard-normal푧 1 ,푧 2 , which are skewed; the component scales are uniform in [0.3, 0.7]; and the mixing weights are drawn from a flat Dirichlet distribution and do not depend on the parents. Data. The training set has 20 480 rows in the accuracy, scalability and learning scenarios and 4096 rows in the batching scenario; the sample-efficiency study uses 4096, 16 384, 32 768 and 65 536 rows. On real networks it has 10 240, 20 480, 40 960 or 81 920 rows, for networks with up to 10 3 , 10 4 , 10 5 or more parameters. The test set has 4096 rows (1024 in the batching scenario and on real networks). The five seeds are 0 to 4. Query selectors. Random: one target and its evidence variables are drawn uniformly for the problem, and the query is asked for the evidence values of 128 test rows. Topological (128 queries per cell): targets are allocated by structural role, 30% to hubs, 30% to cut nodes, 30% to terminals and 10% at random. Hubs are ranked by the size of their Markov blanket, that is, the set of their parents, their children and the other parents of their children. Cut nodes are the articulation points of the moral graph: the undirected graph obtained by linking the parents of every node and dropping the direction of the edges, which a cut node disconnects when it is removed. Terminals are ranked by their depth minus their number of descendants. Hubs, cut nodes and random targets are asked in the prediction or in the diagnosis direction with equal probability, and terminals in the prediction direction only. The evidence is taken along the longest path (40% of the queries), inside the Markov blanket (40%) or at random (20%). Heaviest: the top-ranked hub, cut node and terminal, each in the prediction and in the diagnosis direction and with full and with empty evidence, for at most twelve queries per cell (eight when the graph has no cut node); the evidence is taken along the longest path for hubs and terminals and inside the Markov blanket for cut nodes. In the batching scenario each of these queries is asked for the evidence values of 1024 rows. Reference. For a discrete target the reference set has 5000 samples on synthetic networks and as many samples as the training set on real ones; a query whose evidence agrees with fewer than ten reference samples is not scored. For a continuous target the reference is a sample of 2000 draws from푝(푥 푄 | 푥 퐸 ), obtained by likelihood weighting on the true model. Both references are finite samples, which gives the floor of Section 3. For a discrete target it comes from the few reference samples that agree with the evidence, so it varies with the network and with the query. The floor is the average error of an exact answer against the sampled reference and not a bound: an approximate method can fall slightly below the exact engines. Learning metrics. The learning scenarios score the fitted conditionals directly, without queries. On discrete networks the metric is the parameter-recovery TV: the total variation between a learned and a true row of a CPT, averaged over the parent configurations with their true probabilities and then over the nodes. Continuous families have no common parameters to compare, since a linear-Gaussian model, an MDN and a kNN estimator do not share a parameter space, so we score the fitted density itself with the held-out negative log-likelihood (NLL) per node: minus the log-density of a test row under the fitted network, averaged over the test rows and divided by푛. The NLL is a strictly proper scoring rule. Its expectation is the entropy of the true conditionals plus the Kullback–Leibler divergence from the true to the fitted model, so on the same problem the difference between two methods is a difference in divergence from the truth, and it applies to every mechanism that exposes a density. We prefer it to two alternatives. W1, the metric of the inference scenarios, compares samples of a posterior and therefore needs a query and an engine, which would mix the error of learning with the error of inference. The calibration error PIT-KS (the Kolmogorov–Smirnov distance between the probability integral transform of held-out values and the uniform distribution) needs no query, but calibration is necessary and not sufficient: a model that ignores the parents and returns the marginal of each node is perfectly calibrated. Our runs show the difference: on 50-node non-Gaussian networks with 4096 rows FlexCode has less than half the PIT-KS of the linear-Gaussian model (0.028 against 0.065, close to the 0.023 of the MDN) and the worst NLL of the three (1.42 against 1.26 and 1.20). PIT-KS is reported for every method in Appendix B.5. Batches with empty evidence. With empty evidence (Appendix B.3) the퐵queries of a batch are the same marginal query, soNBNanswers it once and the time per query is the time of one query divided by 퐵; we do not count this as a batching gain. Particles. The sampling engines ofNBNuse 2048 particles in the accuracy, scalability and real-network scenarios and 16 in the batching scenario, which measures time only. Pyro uses 50 particles in every scenario. Budgets and execution. The fit budget and the query budget are 1200 s each; the query budget is 600 s in the batching scenario. The memory cap of a cell is 80% of the memory that is free on the machine when the cell starts. Each cell runs in its own process, so that a timeout or an out-of-memory failure ends that cell only. 11 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli Reporting. In the figures the error bar of a method spans the smallest and the largest value that enters its interquartile mean, and the value axes are logarithmic. Each panel has its own legend; a method has the same colour in every figure, and forNBNthe colour is the mechanism and the hatch is the engine. Each panel of the main text shows all applicable external baselines plus twoNBNmethods, the best parametric or neural one and the best non-parametric one (the two best on discrete networks), selected once per family on the accuracy metric of the scenario. On the continuous families of the accuracy and scalability scenarios the linear-Gaussian mechanism is shown as well, so that the NBN method that reaches the largest networks is always visible. The tables of Appendix B report the IQM± 1 2 of the interquartile range. Hardware. All runs use one server with an Intel Xeon E5-2640 v4 CPU at 2.40 GHz (20 cores), 126 GB of RAM and an NVIDIA Titan X (Pascal) GPU with 12 GB of memory. 12 Neural Bayesian Networks B Complete results This appendix reports every method on every scenario of Section 3, including the mechanisms and engines that the main text does not discuss (normalizing-flow, kernel-density and neural categorical conditionals, and the amortized importance-sampling engine AIS), and the fit and query times of the scenarios where the main text reports accuracy only. Methods carry the names of the configuration files, library-model-engine; forNBNthe model is the mechanism of the conditionals. In the figures a bar is the interquartile mean (IQM) over the seeds a method solved and the error bar spans the smallest and the largest value that enters that mean;푓/푛above a bar means that 푓of푛seeds failed, and a missing bar carries the failure status. The value axes are logarithmic. Each panel has its own legend; a method has the same colour in every figure, and forNBNthe colour is the mechanism and the hatch is the engine. In the tables a cell is IQM±half the interquartile range, followed by(푓/푛)when푓of the푛seeds failed for any reason (time budget, memory, NaN output, error); a status alone (timeout,oom,nan) means that no seed was solved, and – that the method was not run in that configuration. In each column the best method is in bold: the one with the fewest failed seeds and, among those, the best IQM. In the learning scenarios the negative log-likelihood is the mean over the held-out rows, divided by the number of nodes. 13 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli B.1 Inference accuracy discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE 10501005001000 n (nodes) 0.05 0.1 0.2 TV ↓ 1/5 2/5 nannantimeouttimeouttimeouttimeout continuous (Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10501005001000 n (nodes) 0.05 0.1 0.2 0.5 W1 ↓ timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout continuous (non-Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10501005001000 n (nodes) 0.05 0.1 0.2 0.5 W1 ↓ timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout Figure 5: Inference accuracy, all methods: accuracy (TV on discrete targets, W1 on continuous ones) against network size. 14 Neural Bayesian Networks discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE 10501005001000 n (nodes) 10 −1 10 0 10 1 10 2 10 3 Query time (s) 1/5 timeouttimeouttimeouttimeout continuous (Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10501005001000 n (nodes) 10 −1 10 0 10 1 10 2 10 3 Query time (s) timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout continuous (non-Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10501005001000 n (nodes) 10 −1 10 0 10 1 10 2 10 3 Query time (s) timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout Figure 6: Inference accuracy, all methods: total query time (s) against network size. 15 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE 10501005001000 n (nodes) 10 −1 10 0 10 1 10 2 10 3 Fit time (s) 1/5 timeouttimeouttimeouttimeout continuous (Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10501005001000 n (nodes) 10 −1 10 0 10 1 10 2 10 3 Fit time (s) timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout continuous (non-Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10501005001000 n (nodes) 10 −1 10 0 10 1 10 2 10 3 Fit time (s) timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout Figure 7: Inference accuracy, all methods: fit time (s) against network size. 16 Neural Bayesian Networks Table 5: Inference accuracy, discrete networks: TV (lower is better), against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-bayes-ve0.061±0.0083 0.0766±0.0059 0.0783±0.0031 0.0783±0.00350.078±0.0014 (1/5) pgmpy-mle-ve0.061±0.0083 0.0766±0.0059 0.0783±0.0031 0.0783±0.0035 0.0782±0.00035 pomegranate-discrete-ve0.0678±0.0054 (2/5)nannantimeouttimeout pyro-empirical-importance0.151±0.0130.201±0.0130.215±0.0054timeouttimeout nbn-cat-ais0.0625±0.00850.0812±0.00520.0836±0.00320.0838±0.00180.0845±0.0014 nbn-cat-avi0.0534±0.00260.0845±0.00510.086±0.00120.082±0.00710.0875±0.0043 nbn-cat-lw0.0624±0.00820.0799±0.00530.0842±0.0030.0825±0.00270.0831±0.0024 nbn-cat-ve0.061±0.0083 0.0766±0.0059 0.0783±0.0031 0.0783±0.0035 0.0782±0.00035 Table 6: Inference accuracy, discrete networks: total query time in seconds, against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-bayes-ve0.0893±0.00180.504±0.273.33±3.31.85±1.742.9±39 (1/5) pgmpy-mle-ve0.0882±0.0028 0.43±0.166.63±8.11.05±0.642.5±19 pomegranate-discrete-ve10.8±11289±16711±24timeouttimeout pyro-empirical-importance43.1±0.81243±0.31490±3.4timeouttimeout nbn-cat-ais1.81±0.0815.27±0.02611.1±0.1354.2±0.54109±1.7 nbn-cat-avi1.36±0.0316.78±0.05513.6±0.07569.3±1.8140±0.92 nbn-cat-lw1.31±0.0465.59±0.06611±0.1854.2±0.57108±0.79 nbn-cat-ve0.193±0.00650.48±0.058 0.594±0.059 0.578±0.083 0.747±0.067 Table 7: Inference accuracy, discrete networks: fit time in seconds, against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-bayes-ve1.77±0.0122.09±0.00732.53±0.025.91±0.1610.5±0.33 (1/5) pgmpy-mle-ve1.77±0.00692.04±0.00442.4±0.0265.41±0.0529.17±0.053 pomegranate-discrete-ve3.23±0.0243.34±0.0263.56±0.042timeouttimeout pyro-empirical-importance 0.116±0.001 0.137±0.0011 0.17±0.0055timeouttimeout nbn-cat-ais5.93±2.124.6±1020.8±0.022153±0.59187±2.5 nbn-cat-avi27.1±4.379.1±5127±0.29156±1.6193±3.1 nbn-cat-lw0.151±0.00150.177±0.000270.209±0.0014 0.556±0.0032 0.907±0.0044 nbn-cat-ve0.324±0.00120.41±0.00180.521±0.00621.46±0.0132.44±0.017 Table 8: Inference accuracy, linear-Gaussian networks: W1 (lower is better), against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-lg-predict0.0359±0.0053 0.0495±0.0051 0.0617±0.011 0.0676±0.0043 0.0394±0.01 pyro-empirical-importance0.17±0.0290.321±0.0710.449±0.027timeouttimeout nbn-flexcode-lw0.0489±0.0150.104±0.0130.129±0.018timeouttimeout nbn-flow-ais0.0534±0.010.0763±0.0093timeouttimeouttimeout nbn-flow-avi0.048±0.0010.0729±0.011timeouttimeouttimeout nbn-flow-lw0.0485±0.00330.0712±0.011timeouttimeouttimeout nbn-kde-lw0.0408±0.0087 (1/5)oomoomoomoom nbn-knn-lw0.0506±0.00590.0864±0.0040.096±0.0058oomoom nbn-lg-ais0.0384±0.00730.0699±0.00820.0947±0.00850.117±0.0180.0723±0.012 nbn-lg-avi0.0477±0.0160.0659±0.010.0983±0.00930.121±0.0210.0703±0.015 nbn-lg-lw0.0353±0.00530.0659±0.00920.0937±0.00540.116±0.0150.0744±0.018 nbn-mdn-ais0.037±0.00510.0665±0.0110.0981±0.0029timeouttimeout nbn-mdn-avi0.0481±0.0160.0643±0.010.099±0.0073timeouttimeout nbn-mdn-lw0.0376±0.00460.0661±0.0130.0971±0.00068timeouttimeout 17 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli Table 9: Inference accuracy, linear-Gaussian networks: total query time in seconds, against network size. Method푛= 10푛= 50푛= 100 푛= 500 푛= 1000 pgmpy-lg-predict0.0963±0.0032 0.399±0.019 0.676±0.035 3.98±0.24 10.1±0.25 pyro-empirical-importance32.1±0.98179±5.5367±2timeouttimeout nbn-flexcode-lw2.15±0.228.69±0.417.4±0.25timeouttimeout nbn-flow-ais14.2±0.21117±4.7timeouttimeouttimeout nbn-flow-avi11±0.455.3±1.7timeouttimeouttimeout nbn-flow-lw11.3±0.3254.3±1.5timeouttimeouttimeout nbn-kde-lw17.2±2.2 (1/5)oomoomoomoom nbn-knn-lw10±257.2±1.1122±1.8oomoom nbn-lg-ais2.92±0.06415.1±0.3931.7±0.3115.6±0.232.2±0.23 nbn-lg-avi0.495±0.00751.76±0.0573.26±0.02416.2±0.2331.7±0.36 nbn-lg-lw0.675±0.0511.88±0.023.46±0.06316.2±0.2531.3±0.16 nbn-mdn-ais5.75±0.1636.1±1.376.1±0.71timeouttimeout nbn-mdn-avi1.92±0.168.74±0.117.5±0.062timeouttimeout nbn-mdn-lw2.15±0.0348.64±0.1117.4±0.16timeouttimeout Table 10: Inference accuracy, linear-Gaussian networks: fit time in seconds, against network size. Method푛= 10푛= 50푛= 100푛= 500 푛= 1000 pgmpy-lg-predict1.95±0.0282±0.00942.12±0.0163.18±0.0264.34±0.04 pyro-empirical-importance0.12±0.012 0.164±0.0047 0.205±0.00058timeouttimeout nbn-flexcode-lw54.4±16437±17928±18timeouttimeout nbn-flow-ais14.6±1.9196±26timeouttimeouttimeout nbn-flow-avi122±3.1141±1.4timeouttimeouttimeout nbn-flow-lw189±15678±36timeouttimeouttimeout nbn-kde-lw2.45±0.12 (1/5)oomoomoomoom nbn-knn-lw1.37±0.00974.73±0.0139.06±0.022oomoom nbn-lg-ais18±3.1238±5.1305±27207±1.1212±1.2 nbn-lg-avi18±0.19123±0.055125±0.69142±3.7162±9.7 nbn-lg-lw0.377±0.00720.496±0.00130.644±0.0067 1.84±0.015 3.24±0.012 nbn-mdn-ais25.4±3.2322±21479±73timeouttimeout nbn-mdn-avi34.6±0.76125±0.35129±0.072timeouttimeout nbn-mdn-lw55.3±15441±28930±11timeouttimeout Table 11: Inference accuracy, non-Gaussian networks: W1 (lower is better), against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-lg-predict0.131±0.070.223±0.0360.0901±0.02 0.0982±0.0490.156±0.012 pyro-empirical-importance0.228±0.0580.343±0.040.228±0.016timeouttimeout nbn-flexcode-lw0.068±0.0370.175±0.0830.131±0.071timeouttimeout nbn-flow-ais0.0599±0.0190.0801±0.0220.0516±0.0028timeouttimeout nbn-flow-avi0.0517±0.00560.083±0.0230.0539±0.002timeouttimeout nbn-flow-lw0.051±0.00260.0809±0.024 0.0508±0.0043timeouttimeout nbn-kde-lw0.0454±0.0079 (1/5)oomoomoomoom nbn-knn-lw0.068±0.0280.11±0.0120.0681±0.0079oomoom nbn-lg-ais0.13±0.070.224±0.0350.0917±0.0170.0998±0.051 0.155±0.013 nbn-lg-avi0.131±0.070.232±0.0240.0903±0.0190.0991±0.05 0.155±0.012 nbn-lg-lw0.129±0.070.226±0.0360.0878±0.0160.102±0.055 0.155±0.011 nbn-mdn-ais0.0483±0.0070.0791±0.0170.0627±0.0054timeouttimeout nbn-mdn-avi0.0499±0.0084 0.0707±0.0150.0513±0.0053timeouttimeout nbn-mdn-lw0.0475±0.0080.0715±0.0160.0522±0.007timeouttimeout 18 Neural Bayesian Networks Table 12: Inference accuracy, non-Gaussian networks: total query time in seconds, against network size. Method푛= 10푛= 50 푛= 100푛= 500 푛= 1000 pgmpy-lg-predict0.105±0.012 0.416±0.031 0.7±0.013 3.79±0.061 9.92±0.14 pyro-empirical-importance30.7±0.55180±6358±3.2timeouttimeout nbn-flexcode-lw1.93±0.198.72±0.3717.4±0.33timeouttimeout nbn-flow-ais13.5±0.14126±2.6110±2.1timeouttimeout nbn-flow-avi10.6±0.3456±2.9110±1.6timeouttimeout nbn-flow-lw10.6±0.2254.8±0.25111±1.5timeouttimeout nbn-kde-lw16±2.1 (1/5)oomoomoomoom nbn-knn-lw9.07±1.757.4±1.4121±1.5oomoom nbn-lg-ais2.44±0.02914.1±0.02730.5±1.4153±7.7220±3.1 nbn-lg-avi0.391±0.0441.73±0.00973.27±0.02215.7±0.231.5±0.28 nbn-lg-lw0.612±0.00931.88±0.0333.48±0.06316±0.1431.1±0.0079 nbn-mdn-ais4.82±0.2428.6±0.2267±8.1timeouttimeout nbn-mdn-avi1.71±0.148.61±0.03217.5±0.22timeouttimeout nbn-mdn-lw1.84±0.0988.69±0.09917.4±0.23timeouttimeout Table 13: Inference accuracy, non-Gaussian networks: fit time in seconds, against network size. Method푛= 10푛= 50푛= 100푛= 500 푛= 1000 pgmpy-lg-predict1.92±0.00212.07±0.0272.11±0.0252.98±0.0633.75±0.021 pyro-empirical-importance 0.131±0.0025 0.162±0.00045 0.204±0.0004timeouttimeout nbn-flexcode-lw51.2±15433±15932±24timeouttimeout nbn-flow-ais10.2±0.7115±2.7166±3timeouttimeout nbn-flow-avi111±1.2142±2.1166±2.1timeouttimeout nbn-flow-lw207±43588±441.11e+03±15timeouttimeout nbn-kde-lw2.31±0.11 (1/5)oomoomoomoom nbn-knn-lw1.32±0.00524.74±0.00959.02±0.015oomoom nbn-lg-ais8.35±0.14140±0.58202±24525±73606±0.84 nbn-lg-avi14.7±1.2123±0.025125±0.73136±1143±1.3 nbn-lg-lw0.383±0.00550.495±0.00220.627±0.0054 1.73±0.032 2.79±0.068 nbn-mdn-ais15.8±3.2288±0.45339±45timeouttimeout nbn-mdn-avi29.1±2.6126±0.64100±5.4timeouttimeout nbn-mdn-lw50.1±16431±16920±19timeouttimeout 19 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli B.2 Scalability discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE 10100100010000 n (nodes) 0.02 0.05 0.1 0.2 TV ↓ 1/5 2/5 nannantimeouttimeoutoom 2/5 continuous (Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10100100010000 n (nodes) 0.05 0.1 0.2 0.5 W1 ↓ 1/5 timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout continuous (non-Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10100100010000 n (nodes) 0.05 0.1 0.2 0.5 W1 ↓ 1/5 timeouttimeout 3/5 timeouttimeout 3/5 timeouttimeout 3/5 timeouttimeout 1/5 oomoomoomoomoom 1/5 timeouttimeout 1/5 timeouttimeout 1/5 timeouttimeout Figure 8: Scalability, all methods: accuracy (TV on discrete targets, W1 on continuous ones) against network size. 20 Neural Bayesian Networks discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE 10100100010000 n (nodes) 0.02 0.05 0.1 0.2 TV ↓ 1/5 2/5 nannantimeout 1/5 oom 2/5 continuous (Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10100100010000 n (nodes) 0.02 0.05 0.1 0.2 0.5 W1 ↓ 1/5 timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout continuous (non-Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10100100010000 n (nodes) 0.05 0.1 0.2 W1 ↓ 1/5 timeouttimeout 3/5 timeouttimeout 3/5 timeouttimeout 3/5 timeouttimeout 1/5 oomoomoomoomoom 1/5 timeouttimeout 1/5 timeouttimeout 1/5 timeouttimeout Figure 9: Scalability, all methods: accuracy with full evidence against network size. 21 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE 10100100010000 n (nodes) 10 −3 10 −2 10 −1 10 0 10 1 10 2 10 3 Time per query (s) 1/5 timeouttimeoutoom 2/5 continuous (Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10100100010000 n (nodes) 10 −3 10 −2 10 −1 10 0 10 1 10 2 10 3 Time per query (s) 1/5 timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout continuous (non-Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10100100010000 n (nodes) 10 −3 10 −2 10 −1 10 0 10 1 10 2 10 3 Time per query (s) 1/5 timeouttimeout 3/5 timeouttimeout 3/5 timeouttimeout 3/5 timeouttimeout 1/5 oomoomoomoomoom 1/5 timeouttimeout 1/5 timeouttimeout 1/5 timeouttimeout Figure 10: Scalability, all methods: time per query (s), both evidence modes against network size. 22 Neural Bayesian Networks discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE 10100100010000 n (nodes) 10 −3 10 −2 10 −1 10 0 10 1 10 2 10 3 Time per query (s) 1/5 timeout 1/5 oom 2/5 continuous (Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10100100010000 n (nodes) 10 −3 10 −2 10 −1 10 0 10 1 10 2 10 3 Time per query (s) 1/5 timeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeouttimeout 1/5 oomoomoomoomoomtimeouttimeouttimeouttimeouttimeouttimeout continuous (non-Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 10100100010000 n (nodes) 10 −3 10 −2 10 −1 10 0 10 1 10 2 10 3 Time per query (s) 1/5 timeouttimeout 3/5 timeouttimeout 3/5 timeouttimeout 3/5 timeouttimeout 1/5 oomoomoomoomoom 1/5 timeouttimeout 1/5 timeouttimeout 1/5 timeouttimeout Figure 11: Scalability, all methods: time per query (s), full evidence against network size. 23 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli Table 14: Scalability, discrete networks: TV (lower is better), both evidence modes, against network size. Method푛= 10푛= 100푛= 1000푛= 10000 pgmpy-bayes-ve0.0284±0.00410.0343±0.00760.037±0.00130.0212±0.0029 (1/5) pgmpy-mle-ve0.0258±0.008 0.0256±0.00740.0343±0.00570.0291±0.0026 pomegranate-discrete-ve0.0159±0.0034 (2/5)nannantimeout pyro-empirical-importance0.123±0.0150.164±0.00870.146±0.022timeout nbn-cat-ais0.0317±0.0040.0392±0.0034 0.0303±0.01oom nbn-cat-avi0.0337±0.00710.0394±0.00460.0351±0.0010.0847±0.079 nbn-cat-lw0.0326±0.00370.0387±0.0020.0409±0.00260.0389±0.0067 nbn-cat-ve0.0256±0.00830.03±0.00180.0306±0.00290.0225±0.0028 (2/5) Table 15: Scalability, discrete networks: TV (lower is better), full evidence, against network size. Method푛= 10푛= 100푛= 1000푛= 10000 pgmpy-bayes-ve0.0491±0.00550.0578±0.010.0645±0.00260.0323±0.0082 (1/5) pgmpy-mle-ve0.0438±0.013 0.0406±0.0180.0582±0.0110.0484±0.0051 pomegranate-discrete-ve0.0276±0.0087 (2/5)nannantimeout pyro-empirical-importance0.12±0.0120.188±0.0240.151±0.0190.141±0.046 (1/5) nbn-cat-ais0.046±0.00820.0595±0.011 0.0498±0.019oom nbn-cat-avi0.049±0.0160.0583±0.00930.0551±0.00660.102±0.078 nbn-cat-lw0.0495±0.00650.06±0.00220.0619±0.00460.0586±0.012 nbn-cat-ve0.0435±0.0140.0476±0.00350.0516±0.00810.0376±0.0067 (2/5) Table 16: Scalability, discrete networks: time per query in seconds, both evidence modes, against network size. Method푛= 10푛= 100푛= 1000푛= 10000 pgmpy-bayes-ve0.000709±7e-054.54±6.323.6±104.96±20 (1/5) pgmpy-mle-ve0.000689±9.3e-052.43±2.218.9±1123.8±13 pomegranate-discrete-ve0.094±0.0875.87±0.23124±3.9timeout pyro-empirical-importance0.359±0.00584.01±0.1240.4±0.34timeout nbn-cat-ais0.0146±0.00040.08±0.000990.808±0.019oom nbn-cat-avi0.00704±0.000360.0806±0.00180.782±0.0186.73±5.1 nbn-cat-lw0.0435±0.000270.131±0.00950.831±0.01510.1±0.035 nbn-cat-ve0.0104±0.00012 0.0296±0.0059 0.113±0.0160.212±0.024 (2/5) Table 17: Scalability, discrete networks: time per query in seconds, full evidence, against network size. Method푛= 10푛= 100푛= 1000푛= 10000 pgmpy-bayes-ve0.000942±7e-059.05±1346.8±219.78±41 (1/5) pgmpy-mle-ve0.000912±0.000134.83±4.334.5±2647.5±26 pomegranate-discrete-ve0.107±0.0875.88±0.21125±4.3timeout pyro-empirical-importance0.344±0.00563.98±0.1540.4±0.4461±10 (1/5) nbn-cat-ais0.0145±0.000590.0802±0.00130.807±0.022oom nbn-cat-avi0.00697±0.000310.0805±0.00160.782±0.0226.72±5 nbn-cat-lw0.0796±0.00110.183±0.0150.884±0.02710.1±0.052 nbn-cat-ve0.0201±0.00021 0.047±0.025 0.189±0.0380.404±0.051 (2/5) Table 18: Scalability, linear-Gaussian networks: W1 (lower is better), both evidence modes, against network size. Method푛= 10푛= 100푛= 1000푛= 10000 pgmpy-lg-predict0.0353±0.0087 0.0515±0.0028 0.0455±0.00890.0482±0.013 pyro-empirical-importance0.178±0.0320.3±0.0340.309±0.0370.218±0.031 (1/5) nbn-flexcode-lw0.0462±0.0120.101±0.026timeouttimeout nbn-flow-ais0.0494±0.0095timeouttimeouttimeout nbn-flow-avi0.0932±0.073timeouttimeouttimeout nbn-flow-lw0.0436±0.0035timeouttimeouttimeout nbn-kde-lw0.0473±0.0058 (1/5)oomoomoom nbn-knn-lw0.0513±0.0080.0813±0.018oomoom nbn-lg-ais0.032±0.00450.0594±0.00450.0584±0.0190.0475±0.004 nbn-lg-avi0.0378±0.00380.0578±0.00910.0546±0.0140.0478±0.013 nbn-lg-lw0.0348±0.00790.059±0.00650.0556±0.00810.0537±0.0037 nbn-mdn-ais0.0488±0.0140.0636±0.019timeouttimeout nbn-mdn-avi0.0395±0.0140.0585±0.016timeouttimeout nbn-mdn-lw0.034±0.00720.0539±0.0047timeouttimeout 24 Neural Bayesian Networks Table 19: Scalability, linear-Gaussian networks: W1 (lower is better), full evidence, against network size. Method푛= 10푛= 100푛= 1000푛= 10000 pgmpy-lg-predict0.0303±0.0130.0553±0.0076 0.0378±0.00540.0466±0.0056 pyro-empirical-importance0.138±0.050.33±0.0510.338±0.0380.226±0.022 (1/5) nbn-flexcode-lw0.0524±0.0190.113±0.024timeouttimeout nbn-flow-ais0.0435±0.013timeouttimeouttimeout nbn-flow-avi0.104±0.084timeouttimeouttimeout nbn-flow-lw0.0467±0.0098timeouttimeouttimeout nbn-kde-lw0.0656±0.011 (1/5)oomoomoom nbn-knn-lw0.0557±0.0170.102±0.033oomoom nbn-lg-ais0.028±0.00470.0637±0.0120.0726±0.0240.0454±0.0038 nbn-lg-avi0.0331±0.00870.0654±0.00370.0603±0.0130.0493±0.0077 nbn-lg-lw0.0291±0.0077 0.0526±0.00140.0652±0.0120.056±0.0094 nbn-mdn-ais0.0399±0.00570.0694±0.019timeouttimeout nbn-mdn-avi0.035±0.0110.0558±0.018timeouttimeout nbn-mdn-lw0.0326±0.00840.0534±0.0077timeouttimeout Table 20: Scalability, linear-Gaussian networks: time per query in seconds, both evidence modes, against network size. Method푛= 10푛= 100푛= 1000 푛= 10000 pgmpy-lg-predict0.000719±1.3e-05 0.00648±0.0014 0.0803±0.001917.3±0.16 pyro-empirical-importance0.237±0.00362.86±0.1130±0.29324±5.4 (1/5) nbn-flexcode-lw0.0146±0.00240.135±0.002timeouttimeout nbn-flow-ais0.0989±0.0055timeouttimeouttimeout nbn-flow-avi0.0793±0.00092timeouttimeouttimeout nbn-flow-lw0.0728±0.00074timeouttimeouttimeout nbn-kde-lw0.0893±0.011 (1/5)oomoomoom nbn-knn-lw0.0457±0.00910.788±0.011oomoom nbn-lg-ais0.0132±0.000290.202±0.0150.264±0.01 2.56±0.0042 nbn-lg-avi0.00272±4.5e-050.025±0.000180.26±0.00512.64±0.04 nbn-lg-lw0.0103±0.00120.0357±0.000720.26±0.00222.57±0.014 nbn-mdn-ais0.0272±0.000960.519±0.039timeouttimeout nbn-mdn-avi0.00938±0.000790.126±0.0073timeouttimeout nbn-mdn-lw0.0172±0.00240.134±0.0027timeouttimeout Table 21: Scalability, linear-Gaussian networks: time per query in seconds, full evidence, against network size. Method푛= 10푛= 100푛= 1000 푛= 10000 pgmpy-lg-predict0.000936±1.7e-05 0.0071±0.0018 0.0812±0.00217.5±0.27 pyro-empirical-importance0.228±0.00472.88±0.1130±0.26324±5.5 (1/5) nbn-flexcode-lw0.0187±0.00260.143±0.0028timeouttimeout nbn-flow-ais0.0987±0.0057timeouttimeouttimeout nbn-flow-avi0.0816±0.0016timeouttimeouttimeout nbn-flow-lw0.075±0.0012timeouttimeouttimeout nbn-kde-lw0.109±0.0095 (1/5)oomoomoom nbn-knn-lw0.0568±0.00850.798±0.0098oomoom nbn-lg-ais0.013±0.000130.204±0.0160.268±0.012 2.57±0.0098 nbn-lg-avi0.00323±6.9e-050.0259±0.000270.262±0.0052.61±0.059 nbn-lg-lw0.018±0.00250.0457±0.00150.271±0.00312.59±0.025 nbn-mdn-ais0.0279±0.00130.519±0.04timeouttimeout nbn-mdn-avi0.00972±0.000740.127±0.0072timeouttimeout nbn-mdn-lw0.024±0.00450.143±0.0048timeouttimeout 25 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli Table 22: Scalability, non-Gaussian networks: W1 (lower is better), both evidence modes, against network size. Method푛= 10푛= 100푛= 1000푛= 10000 pgmpy-lg-predict0.111±0.0540.136±0.049 0.0891±0.0140.101±0.048 pyro-empirical-importance0.228±0.0760.22±0.0320.203±0.0290.24±0.046 (1/5) nbn-flexcode-lw0.0687±0.0290.0981±0.034timeouttimeout nbn-flow-ais0.054±0.0130.0694±0.0084 (3/5)timeouttimeout nbn-flow-avi0.0613±0.0060.0682±0.0039 (3/5)timeouttimeout nbn-flow-lw0.0627±0.0160.0576±0.0028 (3/5)timeouttimeout nbn-kde-lw0.0478±0.012 (1/5)oomoomoom nbn-knn-lw0.0691±0.0190.0783±0.013oomoom nbn-lg-ais0.11±0.0480.131±0.0430.0899±0.0180.0908±0.056 nbn-lg-avi0.115±0.0570.135±0.0440.0942±0.00870.102±0.06 nbn-lg-lw0.113±0.060.13±0.0470.0934±0.00590.0999±0.05 nbn-mdn-ais0.0458±0.0060.0564±0.016 (1/5)timeouttimeout nbn-mdn-avi0.0487±0.00710.0453±0.006 (1/5)timeouttimeout nbn-mdn-lw0.0484±0.00920.0453±0.004 (1/5)timeouttimeout Table 23: Scalability, non-Gaussian networks: W1 (lower is better), full evidence, against network size. Method푛= 10푛= 100푛= 1000푛= 10000 pgmpy-lg-predict0.115±0.0610.142±0.047 0.083±0.00880.101±0.044 pyro-empirical-importance0.232±0.0590.193±0.0170.175±0.0240.209±0.036 (1/5) nbn-flexcode-lw0.0714±0.0260.0956±0.044timeouttimeout nbn-flow-ais0.0517±0.00820.0853±0.012 (3/5)timeouttimeout nbn-flow-avi0.0578±0.0090.0729±0.0081 (3/5)timeouttimeout nbn-flow-lw0.0596±0.0110.0641±0.0042 (3/5)timeouttimeout nbn-kde-lw0.0516±0.0076 (1/5)oomoomoom nbn-knn-lw0.0755±0.0190.0779±0.009oomoom nbn-lg-ais0.123±0.0540.138±0.0430.0864±0.0180.0944±0.058 nbn-lg-avi0.129±0.0650.146±0.0440.0923±0.00970.106±0.053 nbn-lg-lw0.125±0.0640.137±0.0450.0917±0.00990.0992±0.041 nbn-mdn-ais0.049±0.00480.0446±0.014 (1/5)timeouttimeout nbn-mdn-avi0.0575±0.0130.0487±0.0053 (1/5)timeouttimeout nbn-mdn-lw0.0504±0.0120.0414±0.0046 (1/5)timeouttimeout Table 24: Scalability, non-Gaussian networks: time per query in seconds, both evidence modes, against network size. Method푛= 10푛= 100푛= 1000 푛= 10000 pgmpy-lg-predict0.00076±4.5e-05 0.00659±0.0018 0.0954±0.0118.6±0.13 pyro-empirical-importance0.258±0.00392.97±0.03133.1±0.39334±5.7 (1/5) nbn-flexcode-lw0.0158±0.00230.154±0.0053timeouttimeout nbn-flow-ais0.106±0.000860.835±0.0039 (3/5)timeouttimeout nbn-flow-avi0.0806±0.00330.895±0.0087 (3/5)timeouttimeout nbn-flow-lw0.0818±0.00440.87±0.0025 (3/5)timeouttimeout nbn-kde-lw0.107±0.014 (1/5)oomoomoom nbn-knn-lw0.0584±0.0111.03±0.091oomoom nbn-lg-ais0.0232±0.000290.292±0.0642.66±1.82.69±0.097 nbn-lg-avi0.00319±0.000290.0257±0.000280.26±0.0053 2.67±0.033 nbn-lg-lw0.0126±0.00180.0387±0.00160.277±0.00392.87±0.11 nbn-mdn-ais0.045±0.00180.602±0.057 (1/5)timeouttimeout nbn-mdn-avi0.0117±0.000580.135±0.013 (1/5)timeouttimeout nbn-mdn-lw0.0198±0.00120.144±0.0046 (1/5)timeouttimeout 26 Neural Bayesian Networks Table 25: Scalability, non-Gaussian networks: time per query in seconds, full evidence, against network size. Method푛= 10푛= 100푛= 1000 푛= 10000 pgmpy-lg-predict0.000981±8.2e-05 0.00693±0.0019 0.0935±0.007818.5±0.2 pyro-empirical-importance0.248±0.00222.93±0.06332.8±0.026334±5.7 (1/5) nbn-flexcode-lw0.022±0.00390.164±0.0052timeouttimeout nbn-flow-ais0.107±0.000530.831±0.0018 (3/5)timeouttimeout nbn-flow-avi0.083±0.00410.901±0.011 (3/5)timeouttimeout nbn-flow-lw0.0847±0.00470.873±0.0048 (3/5)timeouttimeout nbn-kde-lw0.127±0.014 (1/5)oomoomoom nbn-knn-lw0.0707±0.00911.02±0.089oomoom nbn-lg-ais0.0236±0.000150.291±0.0642.67±1.8 2.68±0.024 nbn-lg-avi0.00389±0.000340.027±0.000290.264±0.0112.7±0.059 nbn-lg-lw0.0218±0.00310.0521±0.00230.287±0.00952.94±0.053 nbn-mdn-ais0.0437±0.00180.613±0.066 (1/5)timeouttimeout nbn-mdn-avi0.0126±0.000520.137±0.014 (1/5)timeouttimeout nbn-mdn-lw0.0269±0.00140.154±0.0071 (1/5)timeouttimeout 27 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli B.3 Batched queries discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE 11285121024 B (queries per batch) 10 −4 10 −3 10 −2 10 −1 Time per query (s) 4/5 1/5 timeouttimeout continuous (Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 11285121024 B (queries per batch) 10 −4 10 −3 10 −2 Time per query (s) timeouttimeouttimeouttimeouttimeout continuous (non-Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 11285121024 B (queries per batch) 10 −4 10 −3 10 −2 Time per query (s) timeouttimeouttimeouttimeouttimeout Figure 12: Batched queries, all methods: time per query (s), full evidence against query batch size 퐵. 28 Neural Bayesian Networks discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE 11285121024 B (queries per batch) 10 −5 10 −4 10 −3 10 −2 10 −1 Time per query (s) 4/5 1/5 timeouttimeout continuous (Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 11285121024 B (queries per batch) 10 −5 10 −4 10 −3 10 −2 10 −1 Time per query (s) timeouttimeouttimeouttimeouttimeout continuous (non-Gaussian) pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW 11285121024 B (queries per batch) 10 −5 10 −4 10 −3 10 −2 10 −1 Time per query (s) timeouttimeouttimeouttimeouttimeout Figure 13: Batched queries, all methods: time per query (s), empty evidence against query batch size 퐵. 29 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli Table 26: Batched queries, discrete networks: time per query in seconds, full evidence, against query batch size 퐵. Method퐵= 1퐵= 128퐵= 512퐵= 1024 pgmpy-bayes-ve0.0346±0 (4/5)– pgmpy-mle-ve0.0671±0.025 (1/5)– pomegranate-discrete-vetimeout0.0216±0.000120.00829±0.00140.0058±0.0016 pyro-empirical-importancetimeout– nbn-cat-ais0.0425±0.00190.000431±4.3e-060.000187±2.4e-050.000154±1.6e-05 nbn-cat-avi0.0417±0.0017 0.000394±1.2e-05 0.00012±4.5e-06 9.14e-05±1.9e-05 nbn-cat-lw0.0425±0.00070.000446±9.4e-060.0002±1.2e-060.000146±6.4e-06 nbn-cat-ve0.00278±0.000160.000688±0.000310.000681±0.000310.000679±0.00031 Table 27: Batched queries, discrete networks: time per query in seconds, empty evidence, against query batch size 퐵. Method퐵= 1퐵= 128퐵= 512퐵= 1024 pgmpy-bayes-ve0.00281±0 (4/5)– pgmpy-mle-ve0.00381±0.0012 (1/5)– pomegranate-discrete-vetimeout0.0206±0.000260.00802±0.00150.00535±0.0014 pyro-empirical-importancetimeout– nbn-cat-ais0.0437±0.00120.000341±9.6e-068.5e-05±4.7e-064.49e-05±5e-06 nbn-cat-avi0.0413±0.00070.000316±9e-067.87e-05±1.8e-064.02e-05±7.7e-07 nbn-cat-lw0.0436±0.00110.000343±9.5e-068.5e-05±4.6e-064.42e-05±6.7e-07 nbn-cat-ve0.00138±0.00013 1.25e-05±5.3e-07 3.94e-06±1.8e-07 2.66e-06±6e-07 Table 28: Batched queries, linear-Gaussian networks: time per query in seconds, full evidence, against query batch size 퐵. Method퐵= 1퐵= 128퐵= 512퐵= 1024 pgmpy-lg-predict0.00262±0.00011– pyro-empirical-importancetimeout– nbn-flexcode-lw0.0341±0.000830.000569±3.4e-050.000425±8.9e-060.000428±2.1e-05 nbn-flow-aistimeout0.00324±2.5e-050.000894±9.4e-060.000625±1.5e-05 nbn-flow-avitimeout0.00325±5.4e-050.000867±3e-050.000593±2e-05 nbn-flow-lwtimeout0.00326±6.9e-050.000896±2.2e-050.000627±1e-05 nbn-kde-lw0.0428±0.000680.0022±4.9e-050.00217±6.7e-050.00217±5.4e-05 nbn-knn-lw0.0459±0.000360.000792±2.1e-050.000682±2.6e-050.000649±1.7e-05 nbn-lg-ais0.0631±0.000750.00071±6.8e-060.000244±3.4e-060.000217±2.6e-06 nbn-lg-avi0.0119±0.00025 0.000136±2.1e-06 6.51e-05±1.8e-06 5.25e-05±1.2e-06 nbn-lg-lw0.0123±0.00030.000199±7.8e-070.000142±1.5e-050.000111±4.6e-07 nbn-mdn-aistimeout0.00118±3.7e-050.000416±1.6e-050.000258±2.3e-05 nbn-mdn-avi0.052±0.00150.000467±1.6e-050.000148±1.7e-069.37e-05±1.5e-06 nbn-mdn-lw0.052±0.00150.000515±3.6e-050.000197±4.5e-060.000172±2.1e-05 Table 29: Batched queries, linear-Gaussian networks: time per query in seconds, empty evidence, against query batch size 퐵. Method퐵= 1퐵= 128퐵= 512퐵= 1024 pgmpy-lg-predict0.00236±0.00013– pyro-empirical-importancetimeout– nbn-flexcode-lw0.033±0.000710.000257±1.5e-066.6e-05±5.6e-073.27e-05±4.3e-07 nbn-flow-aistimeout0.00309±3.8e-050.000773±1.3e-050.000388±6.6e-06 nbn-flow-avitimeout0.00313±4.1e-050.000784±1.6e-050.00039±3.6e-06 nbn-flow-lwtimeout0.00308±6.7e-050.000773±2e-050.000393±1.2e-05 nbn-kde-lw0.0419±0.000650.000336±9.2e-068.42e-05±1.7e-064.21e-05±1.1e-06 nbn-knn-lw0.0439±0.000450.000324±4.1e-068.01e-05±5.9e-074.03e-05±1.1e-06 nbn-lg-ais0.0654±0.00110.000515±6.7e-060.000126±8.5e-076.42e-05±6.5e-07 nbn-lg-avi0.0114±0.00021 8.98e-05±5.8e-07 2.33e-05±5.4e-07 1.19e-05±4.7e-08 nbn-lg-lw0.0118±0.000299.4e-05±1.1e-062.46e-05±2.8e-071.26e-05±1.5e-07 nbn-mdn-aistimeout0.00105±2.3e-050.000266±3.7e-060.000133±2.3e-06 nbn-mdn-avi0.052±0.00140.000414±1.2e-050.000104±1.2e-065.25e-05±1.6e-06 nbn-mdn-lw0.0517±0.00190.000402±1.3e-050.000101±3.3e-065.12e-05±1.6e-06 30 Neural Bayesian Networks Table 30: Batched queries, non-Gaussian networks: time per query in seconds, full evidence, against query batch size 퐵. Method퐵= 1퐵= 128퐵= 512퐵= 1024 pgmpy-lg-predict0.0031±6.8e-05– pyro-empirical-importancetimeout– nbn-flexcode-lw0.0344±0.000790.000556±3.4e-050.000432±9.1e-060.000452±9.4e-06 nbn-flow-aistimeout0.00324±6.6e-060.000911±2.3e-050.00062±9e-06 nbn-flow-avitimeout0.00324±4.9e-060.000846±1.4e-050.000579±9.9e-06 nbn-flow-lwtimeout0.00328±9e-050.000908±1.3e-050.000629±8.9e-06 nbn-kde-lw0.0426±0.00150.00222±5.9e-050.00219±6.6e-050.0022±7.7e-05 nbn-knn-lw0.0466±0.00110.0008±3e-050.000684±1.3e-050.000653±1.9e-05 nbn-lg-ais0.0625±0.000570.000625±1.1e-050.00024±1.6e-050.000197±1.7e-05 nbn-lg-avi0.012±0.00012 0.000141±3.7e-06 7.99e-05±1.9e-05 5.34e-05±2e-07 nbn-lg-lw0.0121±0.000260.000196±1.2e-060.000129±1.4e-050.00011±1.6e-06 nbn-mdn-aistimeout0.00115±1.3e-050.000392±9.5e-060.000239±1.1e-05 nbn-mdn-avi0.0515±0.000440.000444±5.3e-060.000144±5.8e-069.06e-05±1.9e-06 nbn-mdn-lw0.0523±0.000560.000541±2.6e-050.000203±1.3e-060.000151±2.3e-06 Table 31: Batched queries, non-Gaussian networks: time per query in seconds, empty evidence, against query batch size 퐵. Method퐵= 1퐵= 128퐵= 512퐵= 1024 pgmpy-lg-predict0.00278±4.3e-05– pyro-empirical-importancetimeout– nbn-flexcode-lw0.0332±0.000530.000261±4.1e-066.67e-05±2.2e-063.38e-05±6.7e-07 nbn-flow-aistimeout0.0031±9.2e-060.000775±9.2e-060.000387±4.8e-07 nbn-flow-avitimeout0.00307±8.2e-060.000775±4e-060.000383±1e-06 nbn-flow-lwtimeout0.00311±8.4e-050.000785±1.7e-050.000385±1.5e-05 nbn-kde-lw0.0421±0.00150.000338±2.4e-058.64e-05±5.3e-064.41e-05±1.2e-06 nbn-knn-lw0.0446±0.000780.000336±2.2e-058.36e-05±4.3e-064.26e-05±1.9e-06 nbn-lg-ais0.0647±0.000640.000512±1e-050.00013±5e-076.52e-05±1.1e-06 nbn-lg-avi0.0116±0.000259.3e-05±1.1e-06 2.38e-05±6.9e-07 1.2e-05±1.8e-07 nbn-lg-lw0.0116±0.0002 9.22e-05±1.7e-062.41e-05±5e-071.23e-05±2.8e-07 nbn-mdn-aistimeout0.00104±2.1e-050.00026±9.3e-070.00013±2.7e-06 nbn-mdn-avi0.0515±0.000460.000397±5.5e-060.000102±3.6e-065.05e-05±2e-06 nbn-mdn-lw0.0521±0.000170.000409±4.6e-060.000103±2.6e-065.22e-05±5.6e-07 31 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli B.4 Real networks discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE asia sachs mildew barley hailfinder diabetes link munin 10 −2 10 −1 TV ↓ 1/5 1/5 nannannantimeouttimeouttimeouttimeouttimeout Gaussian pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW ecoli70 arth150 10 −1 10 0 W1 ↓ CLG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LWNBN KDE+LWNBN kNN+LWNBN MDN+AISNBN MDN+AVINBN MDN+LW healthcare sangiovese mehra 10 −1 10 1 10 3 10 5 10 7 W1 ↓ oom Figure 14: Real networks, all methods: accuracy (TV on discrete targets, W1 on continuous ones) against network. 32 Neural Bayesian Networks discrete pgmpy Bayes+VEpgmpy MLE+VEpomegranatepyro ISNBN cat+AISNBN cat+AVINBN cat+LWNBN cat+VE asia sachs mildew barley hailfinder diabetes link munin 10 −2 10 −1 10 0 10 1 10 2 10 3 10 4 Query time (s) 1/5 1/5 timeouttimeouttimeouttimeouttimeout Gaussian pgmpy LG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LW NBN KDE+LW NBN kNN+LW NBN LG+AIS NBN LG+AVI NBN LG+LW NBN MDN+AIS NBN MDN+AVINBN MDN+LW ecoli70 arth150 10 0 10 1 10 2 Query time (s) CLG pyro IS NBN FlexCode+LW NBN flow+AIS NBN flow+AVI NBN flow+LWNBN KDE+LWNBN kNN+LWNBN MDN+AISNBN MDN+AVINBN MDN+LW healthcare sangiovese mehra 10 −1 10 0 10 1 Query time (s) oom Figure 15: Real networks, all methods: total query time (s) against network. 33 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli Table 32: Real networks, discrete networks: TV (lower is better), against network. Methodasiasachsmildewbarleyhailfinderdiabeteslinkmunin pgmpy-bayes-ve0.00383±0.000740.0124±0.00180.033±0.00240.0241±0.0013 0.0198±0.000920.0407±0.0049 (1/5)0.00183±0.000140.00488±0.00033 pgmpy-mle-ve0.00343±0.00110.0128±0.0016 0.0325±0.00320.0275±0.000510.0203±0.00150.0393±0.0011 (1/5)0.00197±0.000240.00464±9.7e-05 pomegranate-discrete-ve0.0038±0.00140.0122±0.0044nannan0.0208±0.0012nantimeouttimeout pyro-empirical-importance0.0336±0.000690.0818±0.00910.287±0.00410.215±0.0140.136±0.0055timeouttimeouttimeout nbn-cat-ais0.00524±0.00110.0168±0.00240.072±0.00270.0459±0.00360.0333±0.00250.118±0.0110.00604±0.00110.0105±0.0012 nbn-cat-avi0.00478±0.000760.0162±0.00150.0703±0.00480.0401±0.000640.0266±0.000540.128±0.0080.00612±0.00210.0105±0.00093 nbn-cat-lw0.00424±0.000960.0159±0.00220.0709±0.00160.0474±0.00190.0281±0.000520.129±0.00990.00496±0.000320.00929±0.00056 nbn-cat-ve0.00354±0.0012 0.012±0.00110.0465±0.0027 0.024±0.00130.0206±0.000870.0343±0.0016 0.00174±0.00033 0.00401±0.0003 Table 33: Real networks, discrete networks: total query time in seconds, against network. Methodasiasachsmildewbarleyhailfinderdiabeteslinkmunin pgmpy-bayes-ve0.0165±0.0018 0.0183±8.6e-050.331±0.051.82±0.027 0.229±0.01513.6±1.3 (1/5)0.9±0.050.804±0.056 pgmpy-mle-ve0.0171±0.0014 0.0183±0.00051 0.269±0.0411.88±0.10.232±0.01660.8±1.1e+02 (1/5)0.965±0.04 0.763±0.02 pomegranate-discrete-ve0.785±0.0227.47±0.03667.5±0.32146±0.8560.9±2.31.14e+03±8.4timeouttimeout pyro-empirical-importance6.15±0.1811.6±0.41117±0.47217±4.8240±5timeouttimeouttimeout nbn-cat-ais0.328±0.0120.602±0.0113.12±0.035.13±0.03312.5±0.5845.2±0.68255±0.64360±1.6 nbn-cat-avi0.116±0.00330.219±0.00252.25±0.0154.13±0.024.39±0.06538.7±1.765.6±1.194.7±0.55 nbn-cat-lw0.482±0.0140.614±0.00713.39±0.025.35±0.0155.74±0.144.8±1.179.8±1.1111±0.29 nbn-cat-ve0.113±0.00210.115±0.00130.531±0.023 1.01±0.020.733±0.05855.9±5.24.09±0.61.8±0.05 Table 34: Real networks, Gaussian networks: W1 (lower is better), against network. Methodecoli70arth150 pgmpy-lg-predict0.0418±0.002 0.021±7.5e-05 pyro-empirical-importance0.245±0.00290.107±0.0031 nbn-flexcode-lw0.148±0.00770.0389±0.0019 nbn-flow-ais0.0923±0.00733.94±0.062 nbn-flow-avi0.0942±0.00253.91±0.073 nbn-flow-lw0.0942±0.00873.92±0.075 nbn-kde-lw0.11±0.00180.0297±0.0002 nbn-knn-lw0.111±0.00520.0418±0.00083 nbn-lg-ais0.0603±0.00120.0218±0.00034 nbn-lg-avi0.0541±0.00620.0212±1.1e-05 nbn-lg-lw0.0547±0.00120.0219±0.0003 nbn-mdn-ais0.059±0.0020.0498±0.0021 nbn-mdn-avi0.0615±0.000740.0503±0.0011 nbn-mdn-lw0.0586±0.00290.0471±0.004 Table 35: Real networks, Gaussian networks: total query time in seconds, against network. Methodecoli70arth150 pgmpy-lg-predict0.314±0.02 0.425±0.021 pyro-empirical-importance141±0.82317±4.6 nbn-flexcode-lw8.58±0.04713.9±0.11 nbn-flow-ais106±4.2102±1.9 nbn-flow-avi47±0.58101±0.25 nbn-flow-lw45.9±1.199±0.95 nbn-kde-lw76±0.66115±8.3 nbn-knn-lw39.7±0.2245.2±3.3 nbn-lg-ais18.6±0.262.98±0.14 nbn-lg-avi1.49±0.0452.97±0.05 nbn-lg-lw1.64±0.0183.11±0.1 nbn-mdn-ais42.4±0.413.5±0.85 nbn-mdn-avi8.6±0.2313.5±0.41 nbn-mdn-lw8.55±0.3713.6±0.6 34 Neural Bayesian Networks Table 36: Real networks, CLG networks: W1 (lower is better), against network. Methodhealthcaresangiovesemehra pyro-empirical-importance106±150.0583±0.00273.74e+07±1.9e+06 nbn-flexcode-lw221±210.0344±0.00272.03e+08±1.9e+07 nbn-flow-ais572±1.3e+020.0412±0.0042.21e+06±1.6e+05 nbn-flow-avi543±390.042±0.0053.56e+06±2.8e+06 nbn-flow-lw581±1.1e+020.0412±0.00841.97e+06±5.7e+05 nbn-kde-lw16.6±10.0252±0.0017oom nbn-knn-lw23.4±110.028±0.00192.41e+06±7.1e+05 nbn-mdn-ais461±2.3e+020.0271±0.00442.17e+06±4.1e+05 nbn-mdn-avi576±1.9e+02 0.022±0.00292.97e+06±6.8e+05 nbn-mdn-lw625±900.0242±0.0032 1.89e+06±2e+05 Table 37: Real networks, CLG networks: total query time in seconds, against network. Methodhealthcaresangiovesemehra pyro-empirical-importance4.03±0.1917.8±0.642.2±1.9 nbn-flexcode-lw0.199±0.0110.873±0.00921.3±0.037 nbn-flow-ais0.688±0.0328.52±2.67.12±0.023 nbn-flow-avi0.688±0.0465.16±0.157.71±0.51 nbn-flow-lw0.771±0.0485.24±0.117.24±0.26 nbn-kde-lw0.827±0.0379.78±0.23oom nbn-knn-lw0.345±0.0212.25±0.0522.98±0.098 nbn-mdn-ais0.152±0.00771.67±0.0471.2±0.041 nbn-mdn-avi0.126±0.0056 0.705±0.02 1.18±0.055 nbn-mdn-lw0.166±0.00640.747±0.0231.19±0.042 35 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli B.5 Parameter learning discrete pgmpy Bayespgmpy MLEpomegranatepyroNBN catNBN cat-Bayes 10501005001000 n (nodes) 0.01 0.02 0.05 TV (recovery) ↓ Figure 16: Parameter learning, all methods: parameter-recovery TV against network size. 36 Neural Bayesian Networks continuous (Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 10501005001000 n (nodes) 0.02 0.05 PIT-KS ↓ timeouttimeout 3/5 timeouttimeout 1/5 2/5 2/5 1/5 oomoomoomtimeouttimeouttimeout continuous (non-Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 10501005001000 n (nodes) 0.02 0.05 PIT-KS ↓ timeouttimeouttimeouttimeout 1/5 2/5 2/5 oomoomoomoomtimeouttimeouttimeout Figure 17: Parameter learning, all methods: calibration error (PIT-KS) against network size. 37 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli discrete NBN catNBN cat-Bayes 10501005001000 n (nodes) 1.1 1.15 1.2 NLL per node ↓ continuous (Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 10501005001000 n (nodes) 0.9 1 1.1 NLL per node ↓ timeouttimeout 3/5 timeouttimeout 1/5 2/5 2/5 1/5 oomoomoomtimeouttimeouttimeout continuous (non-Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 10501005001000 n (nodes) 1 2 NLL per node ↓ timeouttimeouttimeouttimeout 1/5 2/5 2/5 oomoomoomoomtimeouttimeouttimeout Figure 18: Parameter learning, all methods: held-out negative log-likelihood per node against network size. 38 Neural Bayesian Networks discrete pgmpy Bayespgmpy MLEpomegranatepyroNBN catNBN cat-Bayes 10501005001000 n (nodes) 10 −1 10 0 10 1 Fit time (s) continuous (Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 10501005001000 n (nodes) 10 0 10 1 10 2 10 3 Fit time (s) timeouttimeout 3/5 timeouttimeout 1/5 2/5 2/5 1/5 oomoomoomtimeouttimeouttimeout continuous (non-Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 10501005001000 n (nodes) 10 0 10 1 10 2 10 3 Fit time (s) timeouttimeouttimeouttimeout 1/5 2/5 2/5 oomoomoomoomtimeouttimeouttimeout Figure 19: Parameter learning, all methods: fit time (s) against network size. 39 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli Table 38: Parameter learning, discrete networks: parameter-recovery TV (lower is better), against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-bayes0.00907±0.00230.0423±0.0020.0484±0.00150.0585±0.00180.0592±0.00063 pgmpy-mle0.00908±0.00230.0423±0.0020.0484±0.00150.0585±0.00180.0592±0.00063 pomegranate-discrete-ve0.00908±0.00230.0423±0.0020.0484±0.00150.0585±0.00180.0592±0.00063 pyro-empirical0.00897±0.0021 0.0408±0.0017 0.0463±0.0016 0.0561±0.002 0.0567±0.00066 nbn-cat0.00908±0.00230.0423±0.0020.0484±0.00150.0585±0.00180.0592±0.00063 nbn-cat-bayes0.00897±0.0021 0.0408±0.0017 0.0463±0.0016 0.0561±0.002 0.0567±0.00066 Table 39: Parameter learning, discrete networks: held-out negative log-likelihood per node (lower is better), against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 nbn-cat1.12±0.0171.18±0.0031.19±0.00691.21±0.0111.21±0.0042 nbn-cat-bayes 1.12±0.017 1.13±0.0012 1.12±0.0026 1.13±0.00093 1.13±0.00059 Table 40: Parameter learning, discrete networks: fit time in seconds, against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-bayes1.83±0.0152.18±0.00912.72±0.0926.26±0.1811±0.16 pgmpy-mle1.85±0.0372.09±0.00582.54±0.0165.73±0.0239.75±0.079 pomegranate-discrete-ve3.31±0.00113.44±0.0253.66±0.00334.63±0.0166.26±0.011 pyro-empirical0.12±0.0011 0.14±0.0024 0.182±0.0012 0.469±0.0074 0.711±0.032 nbn-cat0.336±0.00240.401±0.00460.495±0.00391.3±0.032.13±0.02 nbn-cat-bayes0.311±0.00560.391±0.00550.473±0.00411.27±0.00832.11±0.011 Table 41: Parameter learning, linear-Gaussian networks: PIT-KS (lower is better), against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-lg0.0299±0.00560.0181±0.00130.0165±0.000670.0151±2.8e-05 0.0148±2.9e-05 nbn-flexcode0.0358±0.000890.0241±0.000340.0237±0.0004timeouttimeout nbn-flow0.0413±0.00450.026±0.000340.0249±0.00022 (3/5)timeouttimeout nbn-kde0.0397±0.0018 (1/5)0.0272±0.00022 (2/5)0.0289±0.0019 (2/5)0.0285±0.0043 (1/5)oom nbn-knn0.0366±0.00380.0346±0.000490.0394±0.0017oomoom nbn-lg0.032±0.00580.0189±0.00120.0171±0.000560.0152±7.5e-050.0149±6.4e-05 nbn-mdn0.0315±0.00540.0203±0.00048timeouttimeouttimeout Table 42: Parameter learning, linear-Gaussian networks: held-out negative log-likelihood per node (lower is better), against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-lg0.919±0.0780.943±0.0180.924±0.040.937±0.0076 0.931±0.0041 nbn-flexcode0.948±0.0611.04±0.0141.04±0.035timeouttimeout nbn-flow0.968±0.021±0.0060.968±0.013 (3/5)timeouttimeout nbn-kde0.971±0.039 (1/5)1.02±0.018 (2/5)0.973±0.013 (2/5)1.05±0.024 (1/5)oom nbn-knn0.933±0.0720.985±0.0120.984±0.034oomoom nbn-lg0.919±0.0780.943±0.0180.924±0.040.937±0.0076 0.931±0.0041 nbn-mdn0.921±0.0780.947±0.017timeouttimeouttimeout Table 43: Parameter learning, linear-Gaussian networks: fit time in seconds, against network size. Method푛= 10푛= 50푛= 100푛= 500 푛= 1000 pgmpy-lg1.83±0.062.14±0.0252.16±0.0953.44±0.044.02±0.25 nbn-flexcode47±14400±10852±19timeouttimeout nbn-flow182±5.5605±161.16e+03±16 (3/5)timeouttimeout nbn-kde1.96±0.095 (1/5)15.2±0.23 (2/5)34±0.77 (2/5)172±1.3 (1/5)oom nbn-knn1.08±0.0133.76±0.0287.11±0.16oomoom nbn-lg0.37±0.0049 0.449±0.00220.547±0.0073 1.45±0.0014 2.57±0.018 nbn-mdn77.8±23674±27timeouttimeouttimeout 40 Neural Bayesian Networks Table 44: Parameter learning, non-Gaussian networks: PIT-KS (lower is better), against network size. Method푛= 10푛= 50푛= 100푛= 500푛= 1000 pgmpy-lg0.0603±0.00620.063±0.00210.0646±0.0068 0.065±0.00079 0.067±0.00051 nbn-flexcode0.0314±0.00330.027±0.00540.0368±0.013timeouttimeout nbn-flow0.0392±0.00560.0252±0.000670.0225±0.00058timeouttimeout nbn-kde0.0313±0.0024 (1/5)0.0241±0.00052 (2/5)0.0262±0.00022 (2/5)oomoom nbn-knn0.0376±0.00240.0333±0.000770.0325±0.00026oomoom nbn-lg0.0644±0.00470.0637±0.00230.0652±0.0064 0.065±0.00068 0.067±0.00047 nbn-mdn0.0356±0.00880.0205±0.001timeouttimeouttimeout Table 45: Parameter learning, non-Gaussian networks: held-out negative log-likelihood per node (lower is better), against network size. Method푛= 10푛= 50푛= 100푛= 500 푛= 1000 pgmpy-lg1.13±0.0711.31±0.0951.42±0.14 1.36±0.015 1.39±0.016 nbn-flexcode1.11±0.0681.3±0.0271.41±0.048timeouttimeout nbn-flow1.11±0.0761.26±0.0961.84±0.84timeouttimeout nbn-kde1.08±0.087 (1/5)1.45±3.4 (2/5)1.51±0.032 (2/5)oomoom nbn-knn1.11±0.0671.25±0.0531.31±0.048oomoom nbn-lg1.13±0.0711.31±0.0951.42±0.14 1.36±0.015 1.39±0.016 nbn-mdn1.09±0.061.19±0.031timeouttimeouttimeout Table 46: Parameter learning, non-Gaussian networks: fit time in seconds, against network size. Method푛= 10푛= 50푛= 100푛= 500 푛= 1000 pgmpy-lg1.84±0.0331.92±0.022.01±0.0162.8±0.0593.93±0.15 nbn-flexcode47.3±14404±14858±18timeouttimeout nbn-flow178±16553±411.08e+03±10timeouttimeout nbn-kde1.97±0.1 (1/5)15.4±0.23 (2/5)34.3±2 (2/5)oomoom nbn-knn1.1±0.00833.75±0.00847.57±0.77oomoom nbn-lg0.356±0.00750.449±0.013 0.553±0.01 1.51±0.031 2.76±0.028 nbn-mdn79.4±24679±24timeouttimeouttimeout 41 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli B.6 Sample efficiency discrete pgmpy Bayespgmpy MLEpomegranatepyroNBN catNBN cat-Bayes 4096163843276865536 n train 0.05 0.1 TV (recovery) ↓ Figure 20: Sample efficiency, all methods: parameter-recovery TV against training-set size. continuous (Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 4096163843276865536 n train 0.02 0.03 0.04 PIT-KS ↓ timeouttimeoutoomoom 1/5 1/5 1/5 timeout continuous (non-Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 4096163843276865536 n train 0.02 0.05 PIT-KS ↓ timeouttimeoutoomoomtimeout Figure 21: Sample efficiency, all methods: calibration error (PIT-KS) against training-set size. 42 Neural Bayesian Networks discrete NBN catNBN cat-Bayes 4096163843276865536 n train 1.25 1.5 1.75 NLL per node ↓ continuous (Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 4096163843276865536 n train 1 1.1 1.2 NLL per node ↓ timeouttimeoutoomoom 1/5 1/5 1/5 timeout continuous (non-Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 4096163843276865536 n train 1.25 1.5 1.75 2 NLL per node ↓ timeouttimeoutoomoomtimeout Figure 22: Sample efficiency, all methods: held-out negative log-likelihood per node against training-set size. 43 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli discrete pgmpy Bayespgmpy MLEpomegranatepyroNBN catNBN cat-Bayes 4096163843276865536 n train 0.2 0.5 1 2 Fit time (s) continuous (Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 4096163843276865536 n train 10 0 10 1 10 2 10 3 Fit time (s) timeouttimeoutoomoom 1/5 1/5 1/5 timeout continuous (non-Gaussian) pgmpy LGNBN FlexCodeNBN flowNBN KDENBN kNNNBN LGNBN MDN 4096163843276865536 n train 10 0 10 1 10 2 10 3 Fit time (s) timeouttimeoutoomoomtimeout Figure 23: Sample efficiency, all methods: fit time (s) against training-set size. 44 Neural Bayesian Networks Table 47: Sample efficiency, discrete networks: parameter-recovery TV (lower is better), against training-set size. Method푛 푡푟 = 4096 푛 푡푟 = 16384 푛 푡푟 = 32768푛 푡푟 = 65536 pgmpy-bayes0.0977±0.00350.0474±0.00230.0335±0.00110.0236±0.00064 pgmpy-mle0.0979±0.00350.0474±0.00230.0336±0.00110.0236±0.00064 pomegranate-discrete-ve0.0979±0.00350.0474±0.00230.0336±0.00110.0236±0.00064 pyro-empirical0.0847±0.0019 0.0453±0.0021 0.0328±0.001 0.0233±0.00062 nbn-cat0.0979±0.00350.0474±0.00230.0336±0.00110.0236±0.00064 nbn-cat-bayes0.0847±0.0019 0.0453±0.0021 0.0328±0.001 0.0233±0.00062 Table 48: Sample efficiency, discrete networks: held-out negative log-likelihood per node (lower is better), against training-set size. Method푛 푡푟 = 4096 푛 푡푟 = 16384 푛 푡푟 = 32768 푛 푡푟 = 65536 nbn-cat1.8±0.0491.21±0.00471.15±0.0015 1.12±0.0029 nbn-cat-bayes 1.16±0.0017 1.13±0.0031 1.12±0.0035 1.12±0.0031 Table 49: Sample efficiency, discrete networks: fit time in seconds, against training-set size. Method푛 푡푟 = 4096 푛 푡푟 = 16384 푛 푡푟 = 32768 푛 푡푟 = 65536 pgmpy-bayes2.29±0.0252.22±0.00382.25±0.0352.34±0.023 pgmpy-mle2.14±0.00922.18±0.022.23±0.0132.39±0.017 pomegranate-discrete-ve3.5±0.00683.56±0.0493.57±0.0633.59±0.021 pyro-empirical0.135±0.0045 0.142±0.0004 0.16±0.0036 0.187±0.022 nbn-cat0.407±0.00240.41±0.000930.408±0.0030.413±0.019 nbn-cat-bayes0.404±0.00460.406±0.00290.404±0.0020.401±0.0052 Table 50: Sample efficiency, linear-Gaussian networks: PIT-KS (lower is better), against training-set size. Method푛 푡푟 = 4096 푛 푡푟 = 16384푛 푡푟 = 32768푛 푡푟 = 65536 pgmpy-lg0.0189±0.00091 0.018±0.00130.0173±0.00080.0173±0.00073 nbn-flexcode0.0292±0.00130.0238±0.000520.0228±0.00029timeout nbn-flow0.0297±0.00090.0261±0.000930.0252±0.00034timeout nbn-kde0.0323±0.00150.0275±0.0014oomoom nbn-knn0.0397±0.000510.0346±0.000850.0329±0.000760.0304±0.00049 (1/5) nbn-lg0.0187±0.00084 (1/5)0.0181±0.00130.018±0.00140.0176±0.0012 nbn-mdn0.0208±0.00110.0195±0.000550.0206±0.00096 (1/5)timeout Table 51: Sample efficiency, linear-Gaussian networks: held-out negative log-likelihood per node (lower is better), against training-set size. Method푛 푡푟 = 4096 푛 푡푟 = 16384푛 푡푟 = 32768푛 푡푟 = 65536 pgmpy-lg0.943±0.019 0.942±0.0170.944±0.0190.945±0.019 nbn-flexcode1.17±0.0131.05±0.00871.02±0.017timeout nbn-flow1.02±0.00431.01±0.00731±0.0043timeout nbn-kde1.05±0.0131.02±0.013oomoom nbn-knn1.01±0.0120.988±0.0130.981±0.0150.966±0.013 (1/5) nbn-lg0.935±0.015 (1/5) 0.942±0.0170.944±0.0190.945±0.019 nbn-mdn0.953±0.0190.946±0.0170.958±0.011 (1/5)timeout Table 52: Sample efficiency, linear-Gaussian networks: fit time in seconds, against training-set size. Method푛 푡푟 = 4096 푛 푡푟 = 16384푛 푡푟 = 32768 푛 푡푟 = 65536 pgmpy-lg2.06±0.0492.06±0.0412.14±0.0432.12±0.057 nbn-flexcode86±4.3331±9.1664±25timeout nbn-flow232±3.9573±5.5886±42timeout nbn-kde1.29±0.02311.9±0.049oomoom nbn-knn0.737±0.00442.86±0.0349.99±0.3937.7±0.53 (1/5) nbn-lg0.451±0.0072 (1/5) 0.462±0.00390.464±0.000230.509±0.015 nbn-mdn153±7.3589±201.18e+03±12 (1/5)timeout 45 Giovanni Briglia , Marco Lippi, Stefano Mariani, and Franco Zambonelli Table 53: Sample efficiency, non-Gaussian networks: PIT-KS (lower is better), against training-set size. Method푛 푡푟 = 4096 푛 푡푟 = 16384 푛 푡푟 = 32768 푛 푡푟 = 65536 pgmpy-lg0.0632±0.00210.063±0.0020.0635±0.00240.0653±0.0022 nbn-flexcode0.0275±0.00510.028±0.00520.0287±0.0054timeout nbn-flow0.0286±0.00150.0244±0.000230.0234±0.00035timeout nbn-kde0.0305±2.7e-050.025±0.00032oomoom nbn-knn0.0388±0.000890.0345±0.00090.032±0.00092 0.0311±0.0013 nbn-lg0.0646±0.00220.0645±0.0020.0647±0.00250.0663±0.0024 nbn-mdn0.0234±0.0015 0.0206±0.0015 0.021±0.00093timeout Table 54: Sample efficiency, non-Gaussian networks: held-out negative log-likelihood per node (lower is better), against training- set size. Method푛 푡푟 = 4096 푛 푡푟 = 16384 푛 푡푟 = 32768 푛 푡푟 = 65536 pgmpy-lg1.26±0.031.27±0.0351.27±0.0291.27±0.037 nbn-flexcode1.42±0.0371.31±0.0221.29±0.031timeout nbn-flow1.23±0.031.24±0.0551.22±0.026timeout nbn-kde1.45±0.0251.62±0.34oomoom nbn-knn1.25±0.0261.23±0.0261.22±0.027 1.21±0.028 nbn-lg1.26±0.031.27±0.0351.27±0.0291.27±0.037 nbn-mdn1.2±0.03 1.19±0.028 1.19±0.029timeout Table 55: Sample efficiency, non-Gaussian networks: fit time in seconds, against training-set size. Method푛 푡푟 = 4096 푛 푡푟 = 16384 푛 푡푟 = 32768 푛 푡푟 = 65536 pgmpy-lg2.04±0.0441.99±0.0182.04±0.0162.2±0.096 nbn-flexcode83.4±4.5324±17633±19timeout nbn-flow222±13488±4.8816±20timeout nbn-kde1.15±0.0210.2±0.45oomoom nbn-knn0.681±0.0122.46±0.00668.31±0.05430.2±0.49 nbn-lg0.442±0.0072 0.44±0.0015 0.44±0.0089 0.452±0.0021 nbn-mdn142±8.3552±421.09e+03±44timeout 46