Hyperuniform many-particle systems, which encompass crystals, quasicrystals and certain exotic disordered systems, exhibit an anomalous suppression of density fluctuations on macroscopic length scales relative to those of conventional disordered systems. Here we investigate the percolation behaviors of disordered stealthy hyperuniform systems (SHU), a subclass of hyperuniform configurations for which the structure factor vanishes for a finite range of wavevectors near the origin, with the degree of stealthiness controlled via a parameter χ. We construct Delaunay triangulation networks derived from SHU configurations with varying χ as well as Poisson point configurations for the purpose of comparison. We investigate a non-uniform bond percolation process, in which bond occupation probabilities decrease with the Euclidean distance between the connected vertices. In this setting, percolation is induced by varying a tuning parameter z. We estimate the percolation thresholds z_c and critical exponents of the networks via finite-size scaling and the Newman-Ziff algorithm. We find that SHU networks exhibit lower percolation thresholds than Poisson networks. Notably, the percolation threshold of SHU networks decreases with the stealthiness parameter χ, indicating that global connectivity emerges more readily as short-range order increases. Moreover, we show that SHU networks with large χ belong to the same universality class as lattices, while Poisson and low-χ systems show deviations. We relate the shift in critical exponents to the degree of suppression of density fluctuations in the point configurations. Our work extends previous studies on transport properties of SHU systems from continuum two-phase media to networks. These results open new avenues for optimizing the resilience of statistically homogeneous disordered networks.

Connecting the dynamics of biomolecular networks to experimentally measurable cell phenotypes remains a central challenge in systems biology. Here, we introduce a model-based definition of phenotype as a partial steady state that is committed to a certain dynamical outcome while otherwise being minimally constrained. We focus on Boolean models and define dynamical phenotypes as complete trap spaces that maximally specify a chosen set of phenotype-determining nodes that correspond to biomarkers, while keeping external inputs unconstrained. We show that dynamical phenotypes can be efficiently identified without full attractor enumeration. Using four published models, including a 70-node Boolean model of T cell differentiation, we show that dynamical phenotypes recover known cell types and activation states, and indicate the environmental conditions ensuring their existence. We also propose a method to identify informative phenotype-determining nodes based on the canalization of the Boolean functions. The results of this method are supported by two attractor-based approaches and reveal biologically relevant cell state information that is complementary to the phenotypes manually defined by model creators. Our results demonstrate that dynamical phenotypes provide a scalable framework for linking model structure, external inputs, and phenotypic outcomes, and offer a principled tool for model-guided biomarker selection.

Comprehensive analysis of the dynamics of Boolean models of biological systems is hampered by the exponentially large state space. Here we introduce the succession-diagram-based Markov chain (SD Markov chain), a coarse-grained representation that uses trap spaces (unescapable state subspaces) of the Boolean model as the states of a Markov chain. These trap spaces and their succession diagram can be efficiently identified, and constitute a dramatic reduction compared to the full state space. The SD Markov chain preserves the decisions that trap the system’s dynamics while making the state space computationally tractable. Using an ensemble of random Boolean networks with known state transition matrices, we show that the SD Markov chain accurately reproduces attractors, basins of attraction, convergence probabilities, decision transitions, and sequences of events. We illustrate the insights and predictions that arise from the SD Markov chain by analyzing a published model of cancer cell metastasis. By combining the interpretability of the succession diagram with the probabilistic rigor of Markov analysis, the SD Markov chain offers a compact quantitative description of the attractor landscape and provides a new avenue for studying control and stability in complex biological systems.

Asynchronous Boolean networks are a type of discrete dynamical system in which each variable can take one of two states, and a single variable state is updated in each time step according to pre-selected rules. Boolean networks are popular in systems biology due to their ability to model long-term biological phenotypes within a qualitative, predictive framework. Boolean networks model phenotypes as attractors, which are closely linked to minimal trap spaces (inescapable hypercubes in the system’s state space). In biological applications, attractors and minimal trap spaces are typically in one-to-one correspondence. However, this correspondence is not guaranteed: motif-avoidant attractors (MAAs) that lie outside minimal trap spaces are possible. MAAs are rare and poorly understood, despite recent efforts. In this contribution to the BMB JMB Special Collection “Problems, Progress and Perspectives in Mathematical and Computational Biology”, we summarize the current state of knowledge regarding MAAs and present several novel observations regarding their response to node deletion reductions and linear extensions of edges. We conduct large-scale computational studies on an ensemble of 14 000 models derived from published Boolean models of biological systems, and more than 100 million Random Boolean Networks. Our findings quantify the rarity of MAAs; in particular, we only observed MAAs in biological models after applying standard simplification methods, highlighting the role of network reduction in introducing MAAs into the dynamics. We also show that MAAs are fragile to linear extensions: in sparse networks, even a single linear node can disrupt virtually all MAAs. Motivated by this observation, we improve the upper bound on the number of delays needed to disrupt a motif-avoidant attractor.

This study investigates the problem of controlling the dynamics of biological systems to achieve desired outcomes (attractors). Specifically, we describe biological systems with Boolean models, which represent the system as a network, characterize system components (nodes) with binary states, and use discrete functions to describe the state changes of the nodes due to their interactions. We build upon the Feedback Vertex Set (FVS) theory, which guarantees that controlling the nodes in an FVS can drive the system toward a target attractor. However, FVS control can be computationally expensive in large networks. To overcome this, we propose two methods that exploit modularity within networks: (1) selecting the top 10\% of influential nodes in each module based on structural metrics, and (2) identifying a subset of the FVS by prioritizing the highest-ranked nodes within each module. These approaches are evaluated in synthetic Random Boolean Networks (RBNs) and validated on real biological Boolean models. Our results show that modular control strategies are more efficient than global approaches, particularly in networks with clear modular structures, which is pertinent for controlling biological systems.

Many biological and social systems are naturally represented as edge-weighted directed or undirected hypergraphs since they exhibit group interactions involving three or more system units as opposed to pairwise interactions that can be incorporated in graph-theoretic representations. However, finding influential cores in hypergraphs is still not as extensively studied as their graph-theoretic counterparts. To this end, we develop and implement a hypergraph-curvature guided discrete time diffusion process with suitable topological surgeries and edge-weight renormalization procedures for both undirected and directed weighted hypergraphs to find influential cores. We successfully apply our framework for directed hypergraphs to seven metabolic hypergraphs and our framework for undirected hypergraphs to two social (coauthorship) hypergraphs to find influential cores, thereby demonstrating the practical feasibility of our approach. In addition, we prove a theorem showing that a certain edge weight renormalization procedure in a prior research work for Ricci flows for edge-weighted graphs has the undesirable outcome of modifying the edge weights to negative numbers, thereby rendering the procedure impossible to use. This paper formulates algorithmic approaches for finding core(s) of (weighted or unweighted) directed hypergraphs.

Disordered hyperuniform many-particle systems are recently discovered exotic states of matter, characterized by a complete suppression of normalized infinite-wavelength density fluctuations and lack of conventional long-range order. Here, we begin a program to quantify the structural properties of nonhyperuniform and hyperuniform networks. In particular, large two-dimensional (2D) Voronoi networks (graphs) containing approximately 10,000 nodes are created from a variety of different point configurations, including the antihyperuniform HIP, nonhyperuniform Poisson process, nonhyperuniform RSA saturated packing, and both non-stealthy and stealthy hyperuniform point processes. We carry out an extensive study of the Voronoi-cell area distribution of each of the networks through determining multiple metrics that characterize the distribution, including their higher-cumulants. We show that the HIP distribution is far from Gaussian; the Poisson and non-stealthy hyperuniform distributions are Gaussian-like distributions, the RSA and the highest stealthy hyperuniform distributions are also non-Gaussian, with diametrically opposite non-Gaussian behavior of the HIP. Moreover, we compute the Voronoi-area correlation functions $C_{00}(r)$ for the networks [M. A. Klatt and S. Torquato, Phys. Rev. E {\bf 90}, 052120 (2014)]. We show that the correlation functions $C_{00}(r)$ qualitatively distinguish the antihyperuniform, nonhyperuniform and hyperuniform Voronoi networks. We find strong anticorrelations in $C_{00}(r)$ (i.e., negative values) for the hyperuniform networks.

To be used as an analysis tool, it is important that a spatial network's construction algorithm reproduces the structural properties of the original physical embedding. One popular method for converting a two-dimensional (2D) point pattern into a spatial network is the Delaunay triangulation. Here, we apply the Delaunay triangulation to seven types of 2D point patterns, including hyperuniform systems (i.e. systems characterized by completely suppressed normalized infinite-wavelength density fluctuations). We demonstrate that the quartile coefficients of dispersion of multiple centrality measures are capable of rank-ordering hyperuniform and nonhyperuniform systems independently, but they cannot distinguish a nearly hyperuniform system from hyperuniform systems. Thus, in each system, we investigate the local densities of the point pattern $ \rho_{P}(\mathbf{r}_{i};\ell) $ and of the network $ \rho_{G}(n_{i};\ell) $. We reveal that there is a strong correlation between $ \rho_{P}(\mathbf{r}_{i};\ell) $ and $ \rho_{G}(n_{i};\ell) $ in nonhyperuniform systems but no such correlation in hyperuniform systems. Similarly, when calculating the pair-correlation function and local density covariance function on the point pattern and network, the point pattern and network functions are similar only in nonhyperuniform systems. In disordered (i.e. isotropic) hyperuniform systems, the network has a positive local density covariance at small distances; such covariance is not present in the corresponding point patterns. Thus, we demonstrate that the Delaunay triangulation accurately captures the density fluctuations of the underlying point pattern only when the point pattern possesses a positive local density covariance at small distances. Such positive correlation is seen in most real-world systems but is not seen in disordered hyperuniform systems. Generally, the Delaunay triangulation is an effective tool for building a spatial network from a 2D point pattern to reproduce the underlying density of the point pattern, but there are situations (i.e. disordered hyperuniform systems) where we caution that the Delaunay triangulation would not be effective at capturing the underlying physical embedding.

Boolean Networks (BNs) describe the time evolution of binary states using logic functions on the nodes of a network. They are fundamental models for complex discrete dynamical systems, with applications in various areas of science and engineering, and especially in systems biology. A key aspect of the dynamical behavior of BNs is the number of attractors, which determines the diversity of long-term system trajectories. Due to the noisy nature and incomplete characterization of biological systems, a stochastic asynchronous update scheme is often more appropriate than the deterministic synchronous one. AND-NOT BNs, whose logic functions are the conjunction of literals, are an important subclass of BNs because of their structural simplicity and their usefulness in analyzing biological systems for which the only information available is a collection of interactions among components. In this paper, we establish new theoretical results regarding asynchronous attractors in AND-NOT BNs. We derive two new upper bounds for the number of asynchronous attractors in an AND-NOT BN based on structural properties (strong even cycles and dominating sets, respectively) of the AND-NOT BN. These findings contribute to a more comprehensive understanding of asynchronous dynamics in AND-NOT BNs, with implications for attractor enumeration and counting, as well as for network design and control.

Humans are exposed to various chemicals, and understanding their toxicity is essential for safeguarding their health. However, current toxicological methods often focus on specific diseases or chemicals, missing broader patterns of toxicity. Here, we propose a general framework that conceptualizes toxicity as arising from protein-interaction-mediated relationships between chemicals and diseases. Using the Comparative Toxicogenomics Database, we show that the network proximity between chemical targets and disease-associated proteins in the protein interactome predicts chemical-disease associations. We validate this finding through zebrafish experiments for acute toxicity and analyses of human exposome data sets for chronic chemical-disease correlation. We highlight the applications of our framework in predicting a previously unknown toxic effect of a fungicide, explaining drug cytotoxicity in a COVID drug screening, and uncovering the overlooked ″indirect-hit″ chemical-disease toxicity pattern. Our approach establishes a new paradigm for toxicity research, offering insights into the health risks posed by chemicals and addressing public health challenges.

In modeling signal transduction networks, it is common to manually integrate experimental evidence through a process that involves trial and error constrained by domain knowledge. We implement a genetic algorithm-based workflow (boolmore) to streamline Boolean model refinement. Boolmore adjusts the functions of the model to enhance agreement with a corpus of curated perturbation-observation pairs. It leverages existing mechanistic knowledge to automatically limit the search space to biologically plausible models. We demonstrate boolmore's effectiveness in a published plant signaling model that exemplifies the challenges of manual model construction and refinement. The refined models surpass the accuracy gain achieved over two years of manual revision and yield new, testable predictions. By automating the laborious task of model validation and refinement, this workflow is a step towards fast, fully automated, and reliable model construction.

Interacting biological systems at all organizational levels display emergent behavior. Modeling these systems is made challenging by the number and variety of biological components and interactions (from molecules in gene regulatory networks to species in ecological networks) and the often-incomplete state of system knowledge (e.g., the unknown values of kinetic parameters for biochemical reactions). Boolean networks have emerged as a powerful tool for modeling these systems. We provide a methodological overview of Boolean network models of biological systems. After a brief introduction, we describe the process of building, analyzing, and validating a Boolean model. We then present the use of the model to make predictions about the system's response to perturbations and about how to control (or at least influence) its behavior. We emphasize the interplay between structural and dynamical properties of Boolean networks and illustrate them in three case studies from disparate levels of biological organization.

Stomata are pores on plant aerial surfaces, each bordered by a pair of guard cells. They control gas exchange vital for plant survival. Understanding how guard cells respond to environmental signals such as atmospheric carbon dioxide (CO2) levels is not only insightful to fundamental biology but also relevant to real-world issues of crop productivity under global climate change. In the past decade, multiple important signaling elements for stomatal closure induced by elevated CO2 have been identified. Yet, there is no comprehensive understanding of high CO2-induced stomatal closure. In this work, we assemble a cellular signaling network underlying high CO2-induced stomatal closure by integrating evidence from a comprehensive literature analysis. We further construct a Boolean dynamic model of the network, which allows in silico simulation of the stomatal closure response to high CO2 in wild-type Arabidopsis thaliana plants and in cases of pharmacological or genetic manipulation of network nodes. Our model has a 91% accuracy in capturing known experimental observations. We perform network-based logical analysis and reveal a feedback core of the network, which dictates cellular decisions in closure response to high CO2. Based on these analyses, we predict and experimentally confirm that applying nitric oxide (NO) induces stomatal closure in ambient CO2 and causes hypersensitivity to elevated CO2. Moreover, we predict a negative regulatory relationship between NO and the protein phosphatase ABI2 and find experimentally that NO inhibits ABI2 phosphatase activity. The experimental validation of these model predictions demonstrates the effectiveness of network-based modeling and highlights the decision-making role of the feedback core of the network in signal transduction. We further explore the model’s potential in predicting targets of signaling elements not yet connected to the CO2 network. Our combination of network science, in silico model simulation, and experimental assays demonstrates an effective interdisciplinary approach to understanding system-level biology.

The extinction of a species in a plant-pollinator mutualistic community can cause cascading effects and lead to major biodiversity loss. The ecologically important task of predicting the severity of the cascading effects is made challenging by the complex network of interactions among the species. In this work, we analyze an ensemble of models of communities of plant and pollinator species. These models describe the mutualistic inter-species interactions by Boolean threshold functions. We show that identifying generalized positive feedback loops can help pinpoint the species whose extinction leads to catastrophic and substantial damage to the whole community. We compare these results with the damage percentage caused by the loss of species identified as important by previously studied structural measures and show that positive feedback loops and the information gained from them can identify certain crucial species that the other measures fail to find. We also suggest mitigation measures for two specific purposes: 1. prevent the damage to the community by protecting a subset of the species, and 2. restore the community after the damage by restoring a subset of species. Our analyses indicate that the generalized positive feedback loops predict the most efficient strategies to achieve these purposes. The correct identification of species in each category has important implications for conservation efforts and developing community management strategies.

<p>The full names of components in the EMT network corresponding to the abbreviated node labels used in Figure 3 and Figure 5A (S1); Evidence for interactions among nodes in the EMT network (S2); Boolean update rules and initial state of the 68 node EMT network (S3); Key experimental outcomes reproduced by the EMT network model (S4); Boolean rules and initial state for the 19 node reduced EMT network (S5).</p>

<p>Word file containing the Supplemental Text, Supplemental Figures S1-S6, and Supplemental Tables S1-S2.</p>

Description of additional methods and procedures used in the study. Also includes Supplementary References.

Network-based dynamic modeling is useful for studying how complex biomolecular systems respond to environmental changes and internal perturbations. The main challenge in constructing a dynamic model is integrating evidence from perturbation (e.g. gene knockout) experiments, because their results arise from the collective function of the regulatory network. For a model to capture these non-local effects, its construction, validation, and refinement necessarily involve trial and error, constrained by domain knowledge. We propose and implement a genetic algorithm-based workflow to streamline model refinement. This workflow applies to any biological system for which an interaction network and enough perturbation experiments exist. The algorithm we introduce adjusts the functions of the model to enhance agreement with a corpus of curated experimental results and leverages existing mechanistic knowledge to automatically limit the search space to biologically plausible models. To account for the interdependence of experimental results, we develop a hierarchical scoring technique for assessing model performance. We implement our workflow for Boolean networks, which are a popular and successful tool for modeling biological systems, but the workflow is readily adaptable to multi-level discrete models. Our implementation is available as the open-source Python library boolmore . We demonstrate boolmore ’s effectiveness in a series of published plant signaling models that exemplify the challenges of manual model construction and refinement. These models describe how plant stomata close in response to the drought hormone abscisic acid. After several hours of automatic refinement on a personal computer, the fittest models recapture and surpass the accuracy gain achieved over 10 years of manual revision. The refined models yield new, testable predictions, such as explanations for the role of reactive oxygen species in drought response. By automating the laborious task of model validation and refinement, this workflow is a step towards fast, fully automated, and reliable model construction. Author summary Biomolecular networks are quintessential complex systems, wherein the interactions of proteins and molecules give rise to cellular phenotypes. Modeling these systems requires making choices about the rules governing individual genes and proteins, but often experiments only constrain their effect on the system-level behavior. This contrast presents a challenge to updating an existing model to align with new experiments. The traditional approach to revising a baseline model is essentially trial-and-error. We present a method, implemented as the open source Python library boolmore , that leverages recent advances in the computational analysis of discrete dynamical systems to automate this process, reducing a task that often takes years to a matter of several hours on a personal computer. We showcase the power of this method on a model describing how plant leaf pores respond to the drought hormone abscisic acid. This model was first published in 2006 and has been updated several times, by hand, to incorporate new experimental data or to improve model performance. Boolmore not only recaptures these refinements, but produces models that better explain experimental results and uncover new insights into the regulatory mechanisms of drought response.

Supplemental File 1. Excel file with additional details on the mathematical model.