Stochastic ordering tools for continuous-time Markov chains and applications to reaction network models

Abstract: Stochastic reaction networks are mathematical models with a wide range of applications in biochemistry, ecology, and epidemiology, and are often complex to analyze. Except for some special cases, it is generally difficult to predict how the abundances of all considered species evolve over time. A possible approach to address this issue is to develop tools to compare the model under study with a similar one whose behavior is better understood. The main contribution of our work is to provide direct and computable conditions that can be used to ensure the existence of an ordered coupling between two stochastic reaction networks and to identify which parameter changes in a given model lead to an increase or decrease in the count of certain species. We also make available an algorithm that implements our theory, and we illustrate it with several applications.

• The paper introduces a tractable computational framework using finite linear programs to establish pathwise ordering between CTMCs in reaction network models.
• It develops a novel coupling construction that ensures monotonicity in species abundances under parameter changes and structural perturbations.
• The framework enables automated sensitivity analysis, ergodicity transfer, and model reduction in applications such as biochemical and ecological systems.

Stochastic Ordering Frameworks for CTMCs with Applications to Reaction Networks

Introduction

The paper "Stochastic ordering tools for continuous-time Markov chains and applications to reaction network models" (2604.00756) presents a comprehensive and computationally tractable framework for establishing pathwise orderings between continuous-time Markov chains (CTMCs), with particular emphasis on mass-action stochastic reaction networks (SRNs). The central objective is to provide explicit sufficient conditions, verifiable via finite linear programs, that guarantee the existence of ordered couplings between given pairs of Markov processes. This enables practitioners to obtain rigorous conclusions about monotonicity of species abundances under parameter changes or network perturbations—key questions in systems biology, chemistry, and population dynamics.

Generalization of Stochastic Ordering Theory

Classical stochastic ordering results, such as those in [muller2002comparison, shaked2007stochastic, massey1987stochastic], provide abstract order-preservation theorems for Markov processes, but these are not structured to directly address the intricate combinatorics of reaction networks, particularly in the mass-action setting. Most existing results require checking infinitely many inequalities across highly nontrivial state spaces, which is intractable for concrete models.

The paper significantly generalizes and strengthens previous theory, notably the framework of Campos et al. [campos2023comparison], by:

• Introducing a pathwise coupling theorem (Theorem 3.1), applicable to arbitrary CTMCs, encompassing potentially explosive processes, and relaxing the requirements from matrix-preorder structured relations to more general binary relations on the state space.
• Providing a novel coupling construction, based on carefully classified local transitions, which ensures that the coupled processes maintain the prescribed relation until explosion.
• Proving a new sufficient condition (Theorem 4.2) for mass-action SRNs, reducing the order-preservation verification to a finite set of tractable linear constraints; this is made algorithmic through the use of linear programming solvers.

Matrix Preorder Reduction and Computational Tools

A central technical contribution is the reduction of stochastic ordering to tractable constraints by working with matrix preorders induced by $M \in \mathbb{Z}^{m \times d}$. For mass-action SRNs, given two sets of rate constants, coupling can be certified by verifying the directionality of each reaction and resolving whether the supports and stoichiometric constraints permit parameters to yield monotone behavior.

The process is as follows:

1. Matrix Preorder Encoding: Identify the preorder of interest (e.g., componentwise, total, or custom linear functionals).
2. Constraint Generation: For each reaction, generate explicit linear constraints on the rate constants that encode when the action of the reaction cannot violate the desired stochastic order.
3. LP Feasibility: Solve these constraint systems to certify the existence of admissible parameter regimes.
4. Implementation and Parallelization: A public algorithm (provided at https://github.com/giulio-cuniberti/comparing-srns) rapidly explores high-dimensional networks by leveraging problem decomposition and LP solvers' efficiency; computations are highly parallelizable.

Application to Parametric and Structural Sensitivity Analysis

Comparative Statics

The main implication is that the method operationalizes the formal notion of "increasing" or "decreasing" species abundances across parameter regimes. If a set of reaction rate increases (resp., decreases) implies a pathwise order between the original and modified system for all time, expectations and distributional functionals transfer directly. This is applicable to:

• Signaling pathway sensitivity (as in layered Michaelis-Menten networks).
• Epigenetic switch models subject to chromatin modification alterations.
• Competing populations with modified interaction strengths.
• Synthetic circuits subject to perturbations such as promoter/enzyme modifications.

Model Reduction, Perturbation, and Network Comparison

A critical feature is that the results hold not only for parameter changes but also for structural modifications: reactions can be removed or added, and the order can still be established provided the matrix-order constraints are satisfied (by interpreting reaction deletion as taking the corresponding rate constant to zero). This enables rigorous transfer of known controllability, ergodicity, or monotonicity from a simpler "dominant" network to more complex systems—illustrated by demonstrating ergodicity inheritance in the augmentation of deficiency-zero weakly reversible systems.

Case Studies and Examples

The paper provides a series of illustrative examples:

• Reversible Reactions and SIS Epidemiology: Pathwise ordering obtained with simple preorders, matching intuition and prior results.
• Michaelis-Menten Kinetics: Derivation of new monotonicity results under weaker conditions than those previously known; identification of regimes where more refined comparison (e.g., including intermediate complexes) is possible.
• Signaling Cascades: Automated identification of several distinct preordering structures, including nontrivial equivalence relations, revealing deep dependencies within multilayered systems.
• Population Dynamics (Lotka-Volterra variants): Demonstration of symmetric and asymmetric structural orderings, providing insight into effects of competition and removal of interactions.
• Epigenetic Cell Memory: Extension and refinement of control relations amidst complex feedback (as in histone modification models).
• Ergodicity Transfer: Formal proof that ergodicity of a "core" SRN is inherited upon well-ordered augmentation; inaccessibility by classical deficiency theory underscores the utility of the new approach.

Theoretical and Practical Implications

The framework presented provides robust computational tools for the a priori analysis of reaction network behavior under parameter perturbations. For practitioners, this enables:

• Automated identification of critical parameters for monotonicity.
• Exhaustive search over plausible order structures in high-dimensional systems.
• Validation of qualitative behavior in the design of synthetic biological circuits.
• Extension of monotonicity-based arguments to nontrivial, potentially explosive or non-minimal-support systems.

Theoretically, the results clarify the delicate nature and limitations of strong pathwise comparisons: when they exist, they provide powerful guarantees, but many biologically relevant networks—especially those with irreversible transitions, nontrivial conservation relations, or complex boundary behavior—admit no nontrivial preorderings. As suggested, this points toward further research on weaker, expectation-based, or probabilistic ordering principles that would be more widely applicable.

Conclusion

This work advances the interface between stochastic process comparison theory and reaction network analysis by providing finite, efficiently checkable, and algorithmically implementable conditions for stochastic ordering in mass-action CTMCs. The methods successfully extend, refine, and computationally empower earlier abstract results, yielding new insight into parametric monotonicity, the effect of structural perturbations, and ergodicity stability for a wide range of SRNs. The framework and publicly available implementation offer broad utility for modelers in systems biology, chemistry, and ecology. Further work on expectation-level or distributional orders may extend applicability to those networks where pathwise monotonicity is too restrictive.