MOTIVATION: Boolean networks (BNs) have emerged as versatile tools for modeling cellular regulatory mechanisms due to their ability to capture key biological features despite their simplicity. Multiple BN synthesis methods have emerged in recent decades, aiming to infer BNs with dynamics that correspond to experimental data. Often, these methods generate multiple BN candidates or "model ensembles". While these ensembles are valuable for representing cell populations and their heterogeneity, they are typically treated as a single component without examining their constituent features. RESULTS: We present AstroLogics, a novel framework designed to analyze and identify differences in both dynamical behavior and logical regulation within a BN model ensemble. The framework calculates dynamical distances between BNs through exploration of their state transition graphs (STGs), enabling clustering of similarly functioning models that may represent different cellular fates or signaling mechanisms. AstroLogics also identifies key logical properties that govern each cluster, highlighting the core regulatory structures that differentiate model behaviors. Our approach leverages MaBoSS, a stochastic simulation tool that implements the Boolean Kinetic Monte-Carlo algorithm to address time interpretation in BNs. This probabilistic estimation method allows efficient probing of BN dynamics through stochastic simulations, overcoming the computational limitations of exhaustive STG analysis. Our framework also provides powerful visualization and classification of the BN ensemble. Through multiple use cases, we demonstrate how AstroLogics facilitates comprehensive analyses of model diversity and discovery of key regulatory structures within a BN ensemble. AVAILABILITY AND IMPLEMENTATION: The AstroLogics package, along with tutorials and datasets, are available at https://github.com/sysbio-curie/AstroLogics.

MOTIVATION: Genome-scale metabolic models lack explicit regulatory mechanisms, limiting their predictive accuracy for genetic interventions. Current methods for computing genetic Minimal Cut Sets either ignore regulatory networks entirely or use simplified acyclic representations that cannot capture regulatory feedback loops, ubiquitous features critical in cellular modeling. RESULTS: We developed gMISpy, a Python package that that enables efficient computation of genetic Minimal Intervention Sets (gMISs) in integrated genome-scale metabolic and regulatory networks. gMISpy incorporates cyclic regulatory logic into our previous computational framework using layered Boolean networks and BoNesis framework, resulting in a more accurate modeling of how regulatory interactions affect metabolic genes. Benchmarking across four different regulatory networks with Human-GEM showed consistent improvements in prediction accuracy, with Matthews correlation coefficient gains ranging from 2.50% to 14.42%. Validation against cancer data from DepMap and Project Score confirmed that cyclic integration reduces false positives and better captures biological vulnerabilities compared to acyclic approaches. AVAILABILITY AND IMPLEMENTATION: https://github.com/PlanesLab/cyclic-gMISpy.

MOTIVATION: A challenging problem in microbiology is to determine nutritional requirements of microorganisms and culture them, especially for the microbial dark matter detected solely with culture-independent methods. The latter foster an increasing amount of genomic sequences that can be explored with reverse ecology approaches to raise hypotheses on the corresponding populations. Building upon genome-scale metabolic networks (GSMNs) obtained from genome annotations, metabolic models predict contextualized phenotypes using nutrient information. RESULTS: We developed the tool Seed2LP, addressing the inverse problem of predicting source nutrients, or seeds, from a GSMN and a metabolic objective. The originality of Seed2LP is its hybrid model, combining a scalable and discrete Boolean approximation of metabolic activity, with the numerically accurate flux balance analysis (FBA). Seed inference is highly customizable, with multiple search and solving modes, exploring the search space of external and internal metabolites combinations. Application to a benchmark of 107 curated GSMNs highlights the usefulness of a logic modelling method over a graph-based approach to predict seeds, and the relevance of hybrid solving to satisfy FBA constraints. Focusing on the dependency between metabolism and environment, Seed2LP is a computational support contributing to address the multifactorial challenge of culturing possibly uncultured microorganisms. AVAILABILITY AND IMPLEMENTATION: Seed2LP is available on https://github.com/bioasp/seed2lp.

Bioinformatics has always been a prolific domain for generating complex satisfiability and optimization problems. For instance, the synthesis of multi-scale models of biological networks has recently been associated with the resolution of optimization problems mixing Boolean logic and universally quantified linear constraints (OPT+qLP), which can be benchmarked on real-world models. In this paper, we introduce a Counter-Example-Guided Abstraction Refinement (CEGAR) to solve such problems efficiently. Our CEGAR exploits monotone properties inherent to linear optimization in order to generalize counter-examples of Boolean relaxations. We implemented our approach by extending Answer Set Programming (ASP) solver Clingo with a quantified linear constraints propagator. Our prototype enables exploiting independence of sub-formulas to further exploit the generalization of counter-examples. We evaluate the impact of refinement and partitioning on two sets of OPT+qLP problems inspired by system biology. Additionally, we conducted a comparison with the state-of-the-art ASP solver Clingo[lpx] that handles non-quantified linear constraints, showing the advantage of our CEGAR approach for solving large problems.

Boolean networks are largely employed to model the qualitative dynamics of cell fate processes by describing the change of binary activation states of genes and transcription factors with time. Being able to bridge such qualitative states with quantitative measurements of gene expression in cells, as scRNA-seq, is a cornerstone for data-driven model construction and validation. On one hand, scRNA-seq binarisation is a key step for inferring and validating Boolean models. On the other hand, the generation of synthetic scRNA-seq data from baseline Boolean models provides an important asset to benchmark inference methods. However, linking characteristics of scRNA-seq datasets, including dropout events, with Boolean states is a challenging task. We present scBoolSeq, a method for the bidirectional linking of scRNA-seq data and Boolean activation state of genes. Given a reference scRNA-seq dataset, scBoolSeq computes statistical criteria to classify the empirical gene pseudocount distributions as either unimodal, bimodal, or zero-inflated, and fit a probabilistic model of dropouts, with gene-dependent parameters. From these learnt distributions, scBoolSeq can perform both binarisation of scRNA-seq datasets, and generate synthetic scRNA-seq datasets from Boolean traces, as issued from Boolean networks, using biased sampling and dropout simulation. We present a case study demonstrating the application of scBoolSeq's binarisation scheme in data-driven model inference. Furthermore, we compare synthetic scRNA-seq data generated by scBoolSeq with BoolODE's, data for the same Boolean Network model. The comparison shows that our method better reproduces the statistics of real scRNA-seq datasets, such as the mean-variance and mean-dropout relationships while exhibiting clearly defined trajectories in two-dimensional projections of the data.

Abstract. A Boolean network (BN) is a discrete dynamical system defined by a Boolean function that maps to the domain itself. A trap space of a BN is a generalization of a fixed point, which is defined as a subhypercube closed by the function of the BN. A trap space is minimal if it does not contain any smaller trap space. Minimal trap spaces have applications for the analysis of attractors of BNs with various update modes. This paper establishes the computational complexity results of three decision problems related to minimal trap spaces: the decision of the trap space property of a subhypercube, the decision of its minimality, and the decision of the membership of a given configuration to a minimal trap space. Under several cases on Boolean function representations, we investigate the computational complexity of each problem. In the general case, we demonstrate that the trap space property is coNP-complete, and the minimality and the membership properties are [Formula: see text]-complete. The complexities drop by one level in the polynomial hierarchy whenever the local functions of the BN either are unate or are represented using truth-tables, binary decision diagrams, or double disjunctive normal forms (Petri net encoding): the trap space property can be decided in polynomial time, whereas deciding the minimality and the membership are coNP-complete. When the BN is given as its functional graph, all these problems are in P.

MOTIVATION: Many techniques have been developed to infer Boolean regulations from a prior knowledge network (PKN) and experimental data. Existing methods are able to reverse-engineer Boolean regulations for transcriptional and signaling networks, but they fail to infer regulations that control metabolic networks. RESULTS: We present a novel approach to infer Boolean rules for metabolic regulation from time-series data and a PKN. Our method is based on a combination of answer set programming and linear programming. By solving both combinatorial and linear arithmetic constraints, we generate candidate Boolean regulations that can reproduce the given data when coupled to the metabolic network. We evaluate our approach on a core regulated metabolic network and show how the quality of the predictions depends on the available kinetic, fluxomics or transcriptomics time-series data. AVAILABILITY AND IMPLEMENTATION: Software available at https://github.com/bioasp/merrin. SUPPLEMENTARY INFORMATION: Supplementary data are available at https://doi.org/10.5281/zenodo.6670164.

BACKGROUND: The temporal coordination of biological processes by the circadian clock is an important mechanism, and its disruption has negative health outcomes, including cancer. Experimental and theoretical evidence suggests that the oscillators driving the circadian clock and the cell cycle are coupled through phase locking. RESULTS: We present a detailed and documented map of known mechanisms related to the regulation of the circadian clock, and its coupling with an existing cell cycle map which includes main interactions of the mammalian cell cycle. The coherence of the merged map has been validated with a qualitative dynamics analysis. We verified that the coupled circadian clock and cell cycle maps reproduce the observed sequence of phase markers. Moreover, we predicted mutations that contribute to regulating checkpoints of the two oscillators. CONCLUSIONS: Our approach underlined the potential key role of the core clock protein NR1D1 in regulating cell cycle progression. We predicted that its activity influences negatively the progression of the cell cycle from phase G2 to M. This is consistent with the earlier experimental finding that pharmacological activation of NR1D1 inhibits tumour cell proliferation and shows that our approach can identify biologically relevant species in the context of large and complex networks.

Boolean networks model finite discrete dynamical systems with complex behaviours. The state of each component is determined by a Boolean function of the state of (a subset of) the components of the network. This paper addresses the synthesis of these Boolean functions from constraints on their domain and emerging dynamical properties of the resulting network. The dynamical properties relate to the existence and absence of trajectories between partially observed configurations, and to the stable behaviours (fixpoints and cyclic attractors). The synthesis is expressed as a Boolean satisfiability problem relying on Answer-Set Programming with a parametrized complexity, and leads to a complete non-redundant characterization of the set of solutions. Considered constraints are particularly suited to address the synthesis of models of cellular differentiation processes, as illustrated on a case study. The scalability of the approach is demonstrated on random networks with scale-free structures up to 100 to 1,000 nodes depending on the type of constraints.

Cellular reprogramming, a technique that opens huge opportunities in modern and regenerative medicine, heavily relies on identifying key genes to perturb. Most of the existing computational methods for controlling which attractor (steady state) the cell will reach focus on finding mutations to apply to the initial state. However, it has been shown, and is proved in this article, that waiting between perturbations so that the update dynamics of the system prepares the ground, allows for new reprogramming strategies. To identify such sequential perturbations, we consider a qualitative model of regulatory networks, and rely on Binary Decision Diagrams to model their dynamics and the putative perturbations. Our method establishes a set identification of sequential perturbations, whether permanent (mutations) or only temporary, to achieve the existential or inevitable reachability of an arbitrary state of the system. We apply an implementation for temporary perturbations on models from the literature, illustrating that we are able to derive sequential perturbations to achieve trans-differentiation.

Protein signaling networks are static views of dynamic processes where proteins go through many biochemical modifications such as ubiquitination and phosphorylation to propagate signals that regulate cells and can act as feed-back systems. Understanding the precise mechanisms underlying protein interactions can elucidate how signaling and cell cycle progression occur within cells in different diseases such as cancer. Large-scale protein signaling networks contain an important number of experimentally verified protein relations but lack the capability to predict the outcomes of the system, and therefore to be trained with respect to experimental measurements. Boolean Networks (BNs) are a simple yet powerful framework to study and model the dynamics of the protein signaling networks. While many BN approaches exist to model biological systems, they focus mainly on system properties, and few exist to integrate experimental data in them. In this work, we show an application of a method conceived to integrate time series phosphoproteomic data into protein signaling networks. We use a large-scale real case study from the HPN-DREAM Breast Cancer challenge. Our efficient and parameter-free method combines logic programming and model-checking to infer a family of BNs from multiple perturbation time series data of four breast cancer cell lines given a prior protein signaling network. Because each predicted BN family is cell line specific, our method highlights commonalities and discrepancies between the four cell lines. Our models have a Root Mean Square Error (RMSE) of 0.31 with respect to the testing data, while the best performant method of this HPN-DREAM challenge had a RMSE of 0.47. To further validate our results, BNs are compared with the canonical mTOR pathway showing a comparable AUROC score (0.77) to the top performing HPN-DREAM teams. In addition, our approach can also be used as a complementary method to identify erroneous experiments. These results prove our methodology as an efficient dynamic model discovery method in multiple perturbation time course experimental data of large-scale signaling networks. The software and data are publicly available at https://github.com/misbahch6/caspo-ts.

Qualitative models of dynamics of signalling pathways and gene regulatory networks allow for the capturing of temporal properties of biological networks while requiring few parameters. However, these discrete models typically suffer from the so-called state space explosion problem which makes the formal assessment of their potential behaviors very challenging. In this paper, we describe a method to reduce a qualitative model for enhancing the tractability of analysis of transient reachability properties. The reduction does not change the dimension of the model, but instead limits its degree of freedom, therefore reducing the set of states and transitions to consider. We rely on a transition-centered specification of qualitative models by the mean of automata networks. Our framework encompasses the usual asynchronous Boolean and multi-valued network, as well as 1-bounded Petri nets. Applied to different large-scale biological networks from the litterature, we show that the reduction can lead to a drastic improvement for the scalability of verification methods.

Continuous Time Bayesian Networks (CTBNs) provide a powerful means to model complex network dynamics. How- ever, their inference is computationally demanding — especially if one considers incomplete and noisy time-series data. The latter gives rise to a joint state- and parameter estimation problem, which can only be solved numerically. Yet, finding the exact parameterization of the CTBN has often only secondary importance in practical scenarios. We therefore focus on the structure learning problem and present a way to analytically marginalize the Markov chain underlying the CTBN model with respect its parameters. Since the resulting stochastic process is parameter-free, its inference reduces to an optimal filtering problem. We solve the latter using an efficient parallel implementation of a sequential Monte Carlo scheme. Our framework enables CTBN inference to be applied to incomplete noisy time-series data frequently found in molecular biology and other disciplines.

The analysis of the dynamics of Biological Regulatory Networks (BRNs) requires innovative methods to cope with the state-space explosion. This paper settles an original approach for deciding reachability properties based onProcess Hitting, which is a framework suitable for modelling dynamical complex systems. In particular, Process Hitting has been shown to be of interest in providing compact models of the dynamics of BRNs with discrete values. Process Hitting splits a finite number of processes into so-called sorts and describes the way each process is able to act upon (that is, to ‘hit’) another one (or itself) in order to ‘bounce’ it as another process of the same sort with further actions. By using complementary abstract interpretations of the succession of actions in Process Hitting, we build a very efficient static analysis to over- and under-approximate reachability properties, which avoids the need to build the underlying states graph. The analysis is proved to have a low theoretical complexity, in particular when the number of processes per sorts is limited, while a very large number of sorts can be managed. This makes such an approach very promising for the scalable analysis of abstract complex systems. We illustrate this through the analysis of a large BRN of 94 components. Our method replies quasi-instantaneously to reachability questions, while standard model-checking techniques regularly fail because of the combinatoric explosion of behaviours.

The stochastic \pi-calculus is a formalism that has been used for modeling complex dynamical systems where the stochasticity and the delay of transitions are important features, such as in the case of biochemical reactions. Commonly, durations of transitions within stochastic \pi-calculus models follow an exponential law. The underlying dynamics of such models are expressed in terms of continuous-time Markov chains, which can then be efficiently simulated and model-checked. However, the exponential law comes with a huge variance, making it difficult to model systems with accurate temporal constraints. In this paper, a technique for tuning temporal features within the stochastic \pi-calculus is presented. This method relies on the introduction of a stochasticity absorption factor by replacing the exponential distribution with the Erlang distribution, which is a sum of exponential random variables. This paper presents a construction of the stochasticity absorption factor in the classical stochastic \pi-calculus with exponential rates. Tools for manipulating the stochasticity absorption factor and its link with timed intervals for firing transitions are also presented. Finally, the model-checking of such designed models is tackled by supporting the stochasticity absorption factor in a translation from the stochastic \pi-calculus to the probabilistic model checker PRISM.