Abstract Chemical reaction networks (CRNs) are an established population model defined as a system of coupled nonlinear ordinary differential equations across many disciplines. In many applications, for example, in systems biology and epidemiology, CRN parameters such as the kinetic reaction rates can be used as control inputs to steer the system toward a given target. Unfortunately, the resulting optimal control problem is nonlinear, therefore, computationally very challenging. We address this issue by introducing an optimality-preserving reduction algorithm for CRNs. The algorithm partitions the original state variables into a reduced set of macro-variables for which one can define a reduced optimal control problem with provably identical optimal values. The reduction algorithm runs with polynomial time complexity in the size of the CRN. We use this result to reduce verification and control problems of large-scale vaccination models over real-world networks.

Network embedding is a fundamental technique to project a network into a lower-dimensional space while preserving similarities among nodes. Traditional network embeddings primarily capture node proximity, making them effective for community detection but insufficient for identifying roles, i.e., patterns of interaction beyond local neighborhoods. To address this limitation, we introduce a simple and efficient embedding technique based on approximate variants of equitable partitions. Our approach, called \varepsilon-BE, introduces a user-tunable tolerance parameter relaxing the otherwise strict condition for exact equitable partitions that can be hardly found in real-world networks. We exploit a relationship between equitable partitions and equivalence relations for Markov chains and ordinary differential equations to develop a partition refinement algorithm for computing an approximate equitable partition in polynomial time. We extend this framework to weighted and directed networks, ensuring applicability to a more general class of graphs and filling a gap in the literature where few approaches are present. We compare our method against state-of-the-art embedding techniques on synthetic and real-world networks. We report comparable—when not superior—performance for visualization, classification, clustering, and regression tasks with smaller running times, enabling the embedding of large-scale networks that could not be efficiently handled by most of the competing techniques. These results and the capability to handle weighted and directed networks make our approach a compelling alternative for structural network embedding.

Stochastic process algebra have been successfully used to analyze performance of software systems and correctness of communication protocols. While expressive, stochastic process algebra models of realistic systems give rise to Markov chains whose size is exponential in the length of the model, a phenomenon known as state-space explosion problem. Lumpability and fluid limits are two techniques which have proven effective in tackling this challenge. In this context, lumpability can be seen as an intermediate representation between the original (i.e., unaggregated) state space of a process and its fluid interpretation. This purpose has been served thus far by an equivalence relation called strong equivalence, which allows one to reduce the original Markov chain to a smaller, lumped Markov chain. This paper begins by reviewing PEPA++, an extension of the stochastic process algebra PEPA where minimum-based semantics are supplemented by product-based ones. Afterwards, we introduce the notion of exact equivalence and argue that it provides a tighter relation between the unaggregated Markov chain and the fluid representation, compared to the common notion of strong equivalence. Moreover, the paper devises a process algebra model whose fluid limit is shown to be the Goodwin model from economics. Exploiting the agent-based view, the paper derives a compact stochastic model of a single worker that is part of the worker population. This allows to study efficiently stochastic properties that escape the macroscopic view of Goodwin’s model.

Efficient methods for the simulation of quantum circuits on classical computers are crucial for their analysis due to the exponential growth of the problem size with the number of qubits. Here we study lumping methods based on bisimulation, an established class of techniques that has been proven successful for (classic) stochastic and deterministic systems such as Markov chains and ordinary differential equations. Forward constrained bisimulation yields a lower-dimensional model which exactly preserves quantum measurements projected on a linear subspace of interest. Backward constrained bisimulation gives a reduction that is valid on a subspace containing the circuit input, from which the circuit result can be fully recovered. We provide an algorithm to compute the constraint bisimulations yielding coarsest reductions in both cases, using a duality result relating the two notions. As applications, we provide theoretical bounds on the size of the reduced state space for well-known quantum algorithms for search, optimization, and factorization. Using a prototype implementation, we report significant reductions on a set of benchmarks. In particular, we show that constrained bisimulation can boost decision-diagram-based quantum circuit simulation by several orders of magnitude, allowing thus for substantial synergy effects.

Structural network embedding is a crucial step in enabling effective downstream tasks for complex systems that aim to project a network into a lower-dimensional space while preserving similarities among nodes. We introduce a simple and efficient embedding technique based on approximate variants of equitable partitions. The approximation consists in introducing a user-tunable tolerance parameter relaxing the otherwise strict condition for exact equitable partitions that can be hardly found in real-world networks. We exploit a relationship between equitable partitions and equivalence relations for Markov chains and ordinary differential equations to develop a partition refinement algorithm for computing an approximate equitable partition in polynomial time. We compare our method against state-of-the-art embedding techniques on benchmark networks. We report comparable-when not superior-performance for visualization, classification, and regression tasks at a cost between one and three orders of magnitude smaller using a prototype implementation, enabling the embedding of large-scale networks that could not be efficiently handled by most of the competing techniques.

Abstract Efficient methods for the simulation of quantum circuits on classic computers are crucial for their analysis due to the exponential growth of the problem size with the number of qubits. Here we study lumping methods based on bisimulation, an established class of techniques that has been proven successful for (classic) stochastic and deterministic systems such as Markov chains and ordinary differential equations. Forward constrained bisimulation yields a lower-dimensional model which exactly preserves quantum measurements projected on a linear subspace of interest. Backward constrained bisimulation gives a reduction that is valid on a subspace containing the circuit input, from which the circuit result can be fully recovered. We provide an algorithm to compute the constraint bisimulations yielding coarsest reductions in both cases, using a duality result relating the two notions. As applications, we provide theoretical bounds on the size of the reduced state space for well-known quantum algorithms for search, optimization, and factorization. Using a prototype implementation, we report significant reductions on a set of benchmarks. Furthermore, we show that constraint bisimulation complements state-of-the-art methods for the simulation of quantum circuits based on decision diagrams.

Quantitative notions of bisimulation are well-known tools for the minimization of dynamical models such as Markov chains and ordinary differential equations (ODEs). In forward bisimulations, each state in the quotient model represents an equivalence class and the dynamical evolution gives the overall sum of its members in the original model. Here we introduce generalized forward bisimulation (GFB) for dynamical systems over commutative monoids and develop a partition refinement algorithm to compute the coarsest one. When the monoid is (ℝ,+), we recover probabilistic bisimulation for Markov chains and more recent forward bisimulations for nonlinear ODEs. Using (ℝ,•) we get nonlinear reductions for discrete-time dynamical systems and ODEs where each variable in the quotient model represents the product of original variables in the equivalence class. When the domain is a finite set such as the Booleans $\mathbb{B}$, we can apply GFB to Boolean networks (BN), a widely used dynamical model in computational biology. Using a prototype implementation of our minimization algorithm for GFB, we find disjunction- and conjunction-preserving reductions on 60 BN from two well-known repositories, and demonstrate the obtained analysis speed-ups. We also provide the biological interpretation of the reduction obtained for two selected BN, and we show how GFB enables the analysis of a large one that could not be analyzed otherwise. Using a randomized version of our algorithm we find product-preserving (therefore non-linear) reductions on 21 dynamical weighted networks from the literature that could not be handled by the exact algorithm.

BACKGROUND: Boolean Networks (BNs) are a popular dynamical model in biology where the state of each component is represented by a variable taking binary values that express, for instance, activation/deactivation or high/low concentrations. Unfortunately, these models suffer from the state space explosion, i.e., there are exponentially many states in the number of BN variables, which hampers their analysis. RESULTS: We present Boolean Backward Equivalence (BBE), a novel reduction technique for BNs which collapses system variables that, if initialized with same value, maintain matching values in all states. A large-scale validation on 86 models from two online model repositories reveals that BBE is effective, since it is able to reduce more than 90% of the models. Furthermore, on such models we also show that BBE brings notable analysis speed-ups, both in terms of state space generation and steady-state analysis. In several cases, BBE allowed the analysis of models that were originally intractable due to the complexity. On two selected case studies, we show how one can tune the reduction power of BBE using model-specific information to preserve all dynamics of interest, and selectively exclude behavior that does not have biological relevance. CONCLUSIONS: BBE complements existing reduction methods, preserving properties that other reduction methods fail to reproduce, and vice versa. BBE drops all and only the dynamics, including attractors, originating from states where BBE-equivalent variables have been initialized with different activation values The remaining part of the dynamics is preserved exactly, including the length of the preserved attractors, and their reachability from given initial conditions, without adding any spurious behaviours. Given that BBE is a model-to-model reduction technique, it can be combined with further reduction methods for BNs.

It is well known that exact notions of model abstraction and reduction for dynamical systems may not be robust enough in practice because they are highly sensitive to the specific choice of parameters. In this paper we consider this problem for nonlinear ordinary differential equations (ODEs) with polynomial derivatives. We introduce a model reduction technique based on approximate differential equivalence, i.e., a partition of the set of ODE variables that performs an aggregation when the variables are governed by nearby derivatives. We develop algorithms to (i) compute the largest approximate differential equivalence; (ii) construct an approximately reduced model from the original one via an appropriate perturbation of the coefficients of the polynomials; and (iii) provide a formal certificate on the quality of the approximation as an error bound, computed as an over-approximation of the reachable set of the reduced model. Finally, we apply approximate differential equivalences to case studies on electric circuits, biological models, and polymerization reaction networks.

We introduce polynomial couplings, a generalization of probabilistic couplings, to develop an algorithm for the computation of equivalence relations which can be interpreted as a lifting of probabilistic bisimulation to polynomial differential equations, a ubiquitous model of dynamical systems across science and engineering. The algorithm enjoys polynomial time complexity and complements classical partition-refinement approaches because: (a) it implements a local exploration of the system, possibly yielding equivalences that do not necessarily involve the inspection of the whole system of differential equations; (b) it can be enhanced by up-to techniques; and (c) it allows the specification of pairs which ought not be included in the output. Using a prototype, these advantages are demonstrated on case studies from systems biology for applications to model reduction and comparison. Notably, we report four orders of magnitude smaller runtimes than partition-refinement approaches when disproving equivalences between Markov chains.

MOTIVATION: Stochastic reaction networks are a widespread model to describe biological systems where the presence of noise is relevant, such as in cell regulatory processes. Unfortunately, in all but simplest models the resulting discrete state-space representation hinders analytical tractability and makes numerical simulations expensive. Reduction methods can lower complexity by computing model projections that preserve dynamics of interest to the user. RESULTS: We present an exact lumping method for stochastic reaction networks with mass-action kinetics. It hinges on an equivalence relation between the species, resulting in a reduced network where the dynamics of each macro-species is stochastically equivalent to the sum of the original species in each equivalence class, for any choice of the initial state of the system. Furthermore, by an appropriate encoding of kinetic parameters as additional species, the method can establish equivalences that do not depend on specific values of the parameters. The method is supported by an efficient algorithm to compute the largest species equivalence, thus the maximal lumping. The effectiveness and scalability of our lumping technique, as well as the physical interpretability of resulting reductions, is demonstrated in several models of signaling pathways and epidemic processes on complex networks. AVAILABILITY AND IMPLEMENTATION: The algorithms for species equivalence have been implemented in the software tool ERODE, freely available for download from https://www.erode.eu. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Abstract Electric circuits manipulate electric charge and magnetic flux via a small set of discrete components to implement useful functionality over continuous time-varying signals represented by currents and voltages. Much of the same functionality is useful to biological organisms, where it is implemented by a completely different set of discrete components (typically proteins) and signal representations (typically via concentrations). We describe how to take a linear electric circuit and systematically convert it to a chemical reaction network of the same functionality, as a dynamical system. Both the structure and the components of the electric circuit are dissolved in the process, but the resulting chemical network is intelligible. This approach provides access to a large library of well-studied devices, from analog electronics, whose chemical network realization can be compared to natural biochemical networks, or used to engineer synthetic biochemical networks.

Discovering relations between chemical reaction networks (CRNs) is a relevant problem in computational systems biology for model reduction, to explain if a given system can be seen as an abstraction of another one; and for model comparison, useful to establish an evolutionary path from simpler networks to more complex ones. This is also related to foundational issues in computer science regarding program equivalence, in light of the established interpretation of a CRN as a kernel programming language for concurrency. Criteria for deciding if two CRNs can be formally related have been recently developed, but these require that a candidate mapping be provided. Automatically finding candidate mappings is very hard in general since the search space essentially consists of all possible partitions of a set. In this paper we tackle this problem by developing a genetic algorithm for a class of CRNs called influence networks, which can be used to model a variety of biological systems including cell-cycle switches and gene networks. An extensive numerical evaluation shows that our approach can successfully establish relations between influence networks from the literature which cannot be found by exact algorithms due to their large computational requirements.

We study chemical reaction networks (CRNs) as a kernel language for concurrency models with semantics based on ordinary differential equations. We investigate the problem of comparing two CRNs, i.e., to decide whether the trajectories of a source CRN can be matched by a target CRN under an appropriate choice of initial conditions. Using a categorical framework, we extend and relate model-comparison approaches based on structural (syntactic) and on dynamical (semantic) properties of a CRN, proving their equivalence. Then, we provide an algorithm to compare CRNs, running linearly in time with respect to the cardinality of all possible comparisons. Finally, we apply our results to biological models from the literature.

Ordinary differential equations (ODEs) are widespread in many natural sciences including chemistry, ecology, and systems biology, and in disciplines such as control theory and electrical engineering. Building on the celebrated molecules-as-processes paradigm, they have become increasingly popular in computer science, with high-level languages and formal methods such as Petri nets, process algebra, and rule-based systems that are interpreted as ODEs. We consider the problem of comparing and minimizing ODEs automatically. Influenced by traditional approaches in the theory of programming, we propose differential equivalence relations. We study them for a basic intermediate language, for which we have decidability results, that can be targeted by a class of high-level specifications. An ODE implicitly represents an uncountable state space, hence reasoning techniques cannot be borrowed from established domains such as probabilistic programs with finite-state Markov chain semantics. We provide novel symbolic procedures to check an equivalence and compute the largest one via partition refinement algorithms that use satisfiability modulo theories. We illustrate the generality of our framework by showing that differential equivalences include (i) well-known notions for the minimization of continuous-time Markov chains (lumpability), (ii)~bisimulations for chemical reaction networks recently proposed by Cardelli et al., and (iii) behavioral relations for process algebra with ODE semantics. With a prototype implementation we are able to detect equivalences in biochemical models from the literature that cannot be reduced using competing automatic techniques.

We present two quantitative behavioral equivalences over species of a chemical reaction network (CRN) with semantics based on ordinary differential equations. Forward CRN bisimulation identifies a partition where each equivalence class represents the exact sum of the concentrations of the species belonging to that class. Backward CRN bisimulation relates species that have identical solutions at all time points when starting from the same initial conditions. Both notions can be checked using only CRN syntactical information, i.e., by inspection of the set of reactions. We provide a unified algorithm that computes the coarsest refinement up to our bisimulations in polynomial time. Further, we give algorithms to compute quotient CRNs induced by a bisimulation. As an application, we find significant reductions in a number of models of biological processes from the literature. In two cases we allow the analysis of benchmark models which would be otherwise intractable due to their memory requirements.

In software performance engineering, what-if scenarios, architecture optimization, capacity planning, run-time adaptation, and uncertainty management of realistic models typically require the evaluation of many instances. Effective analysis is however hindered by two orthogonal sources of complexity. The first is the infamous problem of state space explosion -- the analysis of a single model becomes intractable with its size. The second is due to massive parameter spaces to be explored, but such that computations cannot be reused across model instances. In this paper, we efficiently analyze many queuing models with the distinctive feature of more accurately capturing variability and uncertainty of execution rates by incorporating general (i.e., non-exponential) distributions. Applying product-line engineering methods, we consider a family of models generated by a core that evolves into concrete instances by applying simple delta operations affecting both the topology and the model's parameters. State explosion is tackled by turning to a scalable approximation based on ordinary differential equations. The entire model space is analyzed in a family-based fashion, i.e., at once using an efficient symbolic solution of a super-model that subsumes every concrete instance. Extensive numerical tests show that this is orders of magnitude faster than a naive instance-by-instance analysis.