  Bayesian network modelling is a well adapted approach to study messy and
highly correlated datasets which are very common in, e.g., systems
epidemiology. A popular approach to learn a Bayesian network from an
observational datasets is to identify the maximum a posteriori network in a
search-and-score approach. Many scores have been proposed both Bayesian or
frequentist based. In an applied perspective, a suitable approach would allow
multiple distributions for the data and is robust enough to run autonomously. A
promising framework to compute scores are generalized linear models. Indeed,
there exists fast algorithms for estimation and many tailored solutions to
common epidemiological issues. The purpose of this paper is to present an R
package abn that has an implementation of multiple frequentist scores and some
realistic simulations that show its usability and performance. It includes
features to deal efficiently with data separation and adjustment which are very
common in systems epidemiology.
