106 papers
2604.04285Amplification at Equilibrium: Structural and Thermodynamic Limitations, and Implementation
Amplifying weak molecular signals is essential in both natural and engineered biochemical systems. While most amplification schemes operate out of equilibrium, relying on kinetic barriers and fuel-driven cascades, it is also possible to amplify at thermodynamic equilibrium by shifting the energy landscape upon addition of an analyte. Equilibrium amplification is appealing because, in principle, it can remain indefinitely in the untriggered state. In this work, we establish fundamental structural and thermodynamic limits on equilibrium-based amplification. We first prove that dimerization networks--systems restricted to complexes of at most two monomers--are inherently incapable of equilibrium amplification. This no-go theorem explains the absence of amplification in prior undercomplementary "strand commutation" designs. We then show that allowing trimeric complexes breaks this barrier. We propose an isometric trimer-based amplifier whose output preserves the size of the input, enabling modular composition, and validate it experimentally, achieving an amplification factor close to the expected $2\times$. Finally, we derive universal thermodynamic bounds applicable to any equilibrium network regardless of complex size: the maximum amplification factor scales linearly with the free energy of interaction between the analyte and the amplifier components. For nucleic acid systems, this implies that the analyte length must grow linearly with the desired amplification factor, and that composing modular amplifiers yields diminishing returns for a fixed analyte. Together, these results delineate the structural and energetic boundaries of equilibrium amplification and rigorously justify the necessity of out-of-equilibrium approaches for achieving high gain.
2604.04285Apr 2026
View
Control of genes by self-organizing multicellular interaction networks
Multicellular self-organization drives development in biological organisms, yet a comprehensive theory is lacking as basic properties of cells can complicate common approaches. Framing such properties by dynamic graphs led to new theoretical propositions for multicellular self-organization in Escherichia coli. Here, corresponding ideas are developed from biologically-general first principles. The resulting perspective could aid both experimental and computational approaches to multicellular biology as well as efforts to control and engineer it.
2603.26530Mar 2026
View
BSTModelKit.jl: A Julia Package for Constructing, Solving, and Analyzing Biochemical Systems Theory Models
We present BSTModelKit.jl, an open-source Julia package for constructing, solving, and analyzing Biochemical Systems Theory (BST) models of biochemical networks. The package implements S-system representations, a canonical power-law formalism for modeling metabolic and regulatory networks. BSTModelKit.jl provides a declarative model specification format, dynamic simulation via ordinary differential equation (ODE) integration, steady-state computation, and global sensitivity analysis using the Morris and Sobol methods. The package leverages the Julia scientific computing ecosystem, in particular the SciML suite of differential equation solvers, to provide efficient and flexible model analysis tools. We describe the mathematical formulation, software design, and demonstrate the package capabilities with illustrative examples.
2603.19115Mar 2026
View
Scalable DNA Ternary Full Adder Enabled by a Competitive Blocking Circuit
DNA adder circuits are programmable reaction networks that process DNA molecular inputs to compute a sum and serve as essential components for digital computation. Currently, DNA adders primarily focus on binary addition. While efforts extend the operational bit-width by minimizing the number of DNA strands and developing carry-transmission mechanisms, challenges such as the susceptibility of carrying information to attenuation and the limited expressive capacity of the binary system impose significant constraints on computational scale. This paper proposes a scalable ternary adder architecture by introducing an innovative competitive blocking (CB) circuit. The architecture employs a dual cooperative optimization strategy that significantly enhances single-bit computational capacity and incorporates a dynamic concentration adjustment (CA) to effectively broaden the computational bit-width. Consequently, a significant increase in molecular computing scale is achieved compared to previous binary adders. Biochemical experimental results indicate that the CB circuit effectively outputs the ternary full-adder bit and successfully performs 10-bit addition. Furthermore, by implementing the CA strategy, this adder can be further extended to support 17-bit addition. This research provides a novel methodological foundation for advancing DNA computing technologies and offers promising potential for scalable digital computing applications.
2603.11684Mar 2026
View
ATP Level and Phosphorylation Free Energy Regulate Trigger-Wave Speed and Critical Nucleus Size in Cellular Biochemical Systems
Trigger waves are self-regenerating propagating fronts that emerge from the coupling of nonlinear reaction kinetics and diffusion. In cells, trigger waves coordinate large-scale processes such as mitotic entry and stress responses. Although the roles of circuit topology and feedback architecture in generating bistability are well established, how nonequilibrium energetic driving shapes wave propagation is less well understood. Here, we employ a thermodynamically consistent reaction--diffusion framework to investigate trigger-wave dynamics in ATP-dependent phosphorylation--dephosphorylation systems. We first recapitulate general expressions for trigger-wave speed in the bistable regime and analyze curvature-induced corrections that determine the minimum critical nucleus required for sustained propagation in higher dimensions. We then apply this framework to two representative systems, treating ATP concentration and the nonequilibrium parameter $γ= [ATP]/(K_{\mathrm{eq}}[ADP][P_i])$ as independent control variables to examine how energetic driving regulates wave propagation. Our results show that ATP and $γ$ not only modulate wave speed, but can also reverse the direction of propagation and reshape the parameter regime supporting trigger waves. The critical excitation radius also depends on both ATP concentration and phosphorylation free energy. These findings identify the intracellular energetic state as a regulator of trigger-wave behavior, linking metabolic conditions to the spatial dynamics of wave propagation. More broadly, this framework connects classical reaction--diffusion theory with ATP-driven biochemical regulation and provides a general perspective on related energy-dependent cellular decision-making processes.
2603.10669Mar 2026
View
Automated Classification of Homeostasis Structure in Input-Output Networks
Homeostasis is widely observed in biological systems and refers to their ability to maintain an output quantity approximately constant despite variations in external disturbances. Mathematically, homeostasis can be formulated through an input-output function mapping an external parameter to an output variable. Infinitesimal homeostasis occurs at isolated points where the derivative of this input-output function vanishes, allowing tools from singularity theory and combinatorial matrix theory to characterize homeostatic mechanisms in terms of network topology. However, the required combinatorial enumeration becomes increasingly intractable as network size grows, and the reliance on advanced graph-theoretic concepts limits accessibility and practical use in biological applications. To overcome these limitations, we develop a Python-based algorithm that automates the identification of homeostasis subnetworks and their associated homeostasis conditions directly from network topology. Given an input-output network specified solely by its connectivity structure and designated input and output nodes, the algorithm identifies the relevant graph-theoretical structures and enumerates all homeostatic mechanisms. We demonstrate its applicability across a range of biological examples, including small and large networks, networks with single or multiple input nodes or parameters, and cases where input and output coincide. This wide applicability stems from our extension of the theoretical framework from single-input-single-output networks to networks with multiple input nodes through an augmented single-input-node representation. The resulting computational framework provides a scalable and systematic approach to classifying homeostatic mechanisms in complex biological networks, facilitating the application of advanced mathematical theory to a broad range of biological systems.
2603.08882Mar 2026
View
Duality in mass-action networks
Mass-action networks are special cases of chemical reaction networks. For these systems, we argue that conserved quantities are dual to internal cycles. We introduce maximal invariant polyhedral supports, and we conjecture that there is a duality relation between preclusters and maximal invariant polyhedral supports. Given the close relation between maximal invariant polyhedral supports and siphons, we also conjecture that siphons and preclusters are dual objects.
2603.08767Mar 2026
View
A cocktail of chemical reaction networks and mathematical epidemiology tools for positive ODE stability problems
We continue recent attempts to put together concepts and results of Chemical Reaction Networks theory (CRNT) and Mathematical Epidemiology (ME), for solving problems of stability of positive ODEs. We provide first an elegant CRN-flavored generalization of the most cited result in ME, the Next Generation Matrix (NGM) theorem. We review next the "symbolic-numeric approach of Vassena and Stadler, which tackles bifurcation problems by viewing the characteristic polynomial of the Jacobian at fixed points as a formal polynomial in the "symbolic reactivities", and identifies its coefficients as "Child Selection minors of the stoichiometric matrix". We also review two applications of this approach using the Mathematica package Epid-CRN tools from both CRNT and ME.
2603.067781Mar 2026
View
Quantifying Ranking Instability Across Evaluation Protocol Axes in Gene Regulatory Network Benchmarking
Benchmark rankings are routinely used to justify scientific claims about method quality in gene regulatory network (GRN) inference, yet the stability of these rankings under plausible evaluation protocol choices is rarely examined. We present a systematic diagnostic framework for measuring ranking instability under protocol shift, including decomposition tools that separate base rate effects from discrimination effects. Using existing single cell GRN benchmark outputs across three human tissues and six inference methods, we quantify pairwise reversal rates across four protocol axes: candidate set restriction (16.3 percent, 95 percent CI 11.0 to 23.4 percent), tissue context (19.3 percent), reference network choice (32.1 percent), and symbol mapping policy (0.0 percent). A permutation null confirms that observed reversal rates are far below random order expectations (0.163 versus null mean 0.500), indicating partially stable but non invariant ranking structure. Our decomposition reveals that reversals are driven by changes in the relative discrimination ability of methods rather than by base rate inflation, a finding that challenges a common implicit assumption in GRN benchmarking. We propose concrete reporting practices for stability aware evaluation and provide a diagnostic toolkit for identifying method pairs at risk of reversal.
2603.03493Mar 2026
View
Topological bounds on the dynamical growth rate of chemical reaction networks
Growth and decay are system-level properties of chemical reaction networks (CRNs) relevant from prebiotic chemistry to cellular metabolism. Their properties are typically analyzed through the kinetics of particular models, which requires specification of the full set of kinetic laws and parameters. In this work, we derive stoichiometry-based constraints on the growth (or shrinkage) rate, in the balanced-growth regime of scalable CRNs. The resulting bounds are controlled by a topological quantity, the maximum amplification factor, defined via a von Neumann max-min problem over feasible fluxes as illustrated by numerical tests on random-network ensembles of CRNs. We argue for the relevance of our results in the context of origin of life studies but also for designing synthetic chemical reaction networks.
2603.02627Mar 2026
View
Prediction of source nutrients for microorganisms using metabolic networks
Metagenomics has lowered the barrier to microbial discovery--enabling the identification of novel microbes without isolation--but cultures remain imperative for the deep study of microbes. Cultivation and isolation of non-model microbes remains a major challenge, despite advances in high-throughput culturomic methods. The quantity of simultaneous experimental variables is constrained by time and resources, but the list can be reduced using computational biology. Given an annotated genome, metabolic modelling can be used to predict source nutrients required for the growth of a microbe, which acts as an initial screen to inform culture and isolation experiments. This chapter provides an overview of metabolic networks and modelling and how they can be used to predict the nutrient requirements of a microorganism, followed by a sample protocol using a toy metabolic network, which is then expanded to a genome-scale metabolic network application. These methods can be applied to any metabolic network of interest--which in turn can be created from any genome of interest--and are a starting point for experimental validation of source nutrients required for microorganisms that remain uncultivated to date.
2602.21993Feb 2026
View
Convex Analysis of Relaxation Dynamics in Chemical Reaction Networks and Generalized Gradient Flows
We obtain bounds on the Kullback--Leibler divergence to equilibrium for mass-action chemical reaction networks (CRNs) with equilibrium. The associated decay rates are characterized in terms of the singular values of the stoichiometric matrix, convexity parameters, and time-integrated activities via deformed-exponential-type functions. We further extend these bounds within a generalized gradient flow framework. We highlight the biological relevance of this framework: the resulting bounds apply to quasi-steady-state regimes, where long transients and plateau-like behavior are common and functionally important. We illustrate the framework using a catalytic CRN exhibiting plateaus, where the bounds capture slow relaxation induced by local convexity and provide a bound-based approach to quantifying relaxation in CRNs.
2602.18932Feb 2026
View
Modeling Dynamics, Cell Type Specificity, and Perturbations in Gene Regulatory Networks
Gene regulatory networks (GRNs) define the regulatory relationships among molecules such as transcription factors, chromatin remodelers, and target genes. GRNs play a critical role in diverse biological processes, including development, disease manifestation, and evolution. However, fully characterizing these networks across multiple cell types and states remains a significant challenge. Recent advances in single-cell omics have dramatically enhanced our ability to measure biological systems at unprecedented resolution. These technologies have opened new avenues for computational methods to infer GRNs, offering deeper insights into cell type-specific mechanisms, causality, and dynamic regulatory processes. This review summarizes the current state of GRN inference from single cell omic datasets, with a particular focus on dynamics and perturbations, and outlines key open challenges that must be addressed to advance the field.
2602.18854Feb 2026
View
Oscillation Criteria in Large-Scale Gene Regulatory Networks with Intrinsic Fluctuations
Gene Regulatory Networks(GRNs) with feedback are essential components of many cellular processes and may exhibit oscillatory behavior. Analyzing such systems becomes increasingly complex as the number of components increases. Since gene regulation often involves a small number of molecules, fluctuations are inevitable. Therefore, it is important to understand how fluctuations affect the oscillatory dynamics of cellular processes, as this will allow comprehension of the mechanisms that enable cellular functions to remain even in the presence of fluctuations or, failing that, to determine the limit of fluctuations that permits various cellular functions. In this study, we investigated the conditions under which GRNs with feedback and intrinsic fluctuations exhibit oscillatory behavior. Our focus was on developing a procedure that would be both manageable and practical, even for extensive regulatory networks, that is, those comprising numerous nodes. Using the second-moment approach, we described the stochastic dynamics through a set of ordinary differential equations for the mean concentration and its second central moment. The system can attain either a stable equilibrium or oscillatory behavior, depending on its scale and, consequently, the intensity of fluctuations. To illustrate the procedure, we analyzed two relevant systems: a repressilator with three nodes and a system with five nodes, both incorporating intrinsic fluctuations. In both cases, it was observed that for very small systems, which therefore exhibit significant fluctuations, oscillatory behavior is inhibited. The procedure presented here for analyzing the stability of oscillations under fluctuations enables the determination of the critical minimum size of GRNs at which intrinsic fluctuations do not eliminate their cyclical behavior.
2602.16129Feb 2026
View
Relating biomarkers and phenotypes using dynamical trap spaces
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 \emph{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. This method reveals biologically relevant cell state information that is complementary to the phenotypes manually defined by model creators and is validated by two attractor-based approaches. 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.
2602.15691Feb 2026
View
Mathematical modeling of 1,2-propanediol utilization bacterial microcompartments in vivo activity
On exposure to 1,2-propanediol (1,2-PD), Salmonella enterica serovar Typhimurium LT2 produces 1,2-PD utilization (Pdu) microcompartments (MCPs), nanoscale protein-bound shells that encapsulate metabolic enzymes. MCPs serve as a bioengineering platform to study reaction organization and enhance flux through specific pathways. However, a recently published assay of purified wild-type (WT) MCPs reported metabolic activity that differed markedly from that observed in vivo. Using kinetic modeling, we attribute these discrepancies to in vivo cell growth and to the cytosolic presence of MCP-associated enzymes and promiscuous alcohol dehydrogenases, which are not present in the purified MCPs. Assays of purified MCPs in E. coli lysate, together with an LT2 growth assay in which the native Pdu MCP-associated alcohol dehydrogenase, PduQ, was knocked out, support the conclusion that exogenous Pdu cytosolic enzyme activity can narrow the gap between in vitro and in vivo experiments. Our modeling further suggests that MCP-localized enzymes contribute little to in vivo metabolic flux downstream of PduCDE. We therefore propose a revised in vivo model of WT growth on 1,2-PD in which PduCDE is fully encapsulated, while much of the downstream Pdu activity occurs in the cytosol.
2602.11338Feb 2026
View
A Modular Mechanistic In Silico Model for In Vitro Transcription Process Yield and Product Quality Prediction
In vitro transcription (IVT) plays a critical role in the manufacture of mRNA vaccines and therapeutics. Optimizing mRNA yield and ensuring product quality, such as capping efficiency and integrity, are essential but mechanistically complex. This study presents a modular mechanistic model of the IVT process to advance scientific understanding and improve predictive capability. The IVT reaction network is decomposed into interconnected modules describing (1) initiation and capping, (2) elongation and truncation, (3) termination and read-through, (4) mRNA degradation, (5) magnesium pyrophosphate precipitation, and (6) enzymatic degradation of pyrophosphate. Guided by biochemical principles and experimental data, kinetic models were developed for each module, accounting for mass balances, molecular complexation, and enzyme activity, and were subsequently assembled to capture coupled IVT dynamics. Multivariate residual analysis and Shapley value-based sensitivity analysis, guided by domain knowledge, were applied to iteratively improve model fidelity. These machine learning-driven analytics enabled identification of key mechanisms, supported in silico experimentation, and facilitated root-cause analysis. Combined with Gaussian-process-based batch Bayesian optimization for efficient parameter estimation, this framework establishes a scalable hybrid (mechanistic + machine learning) modeling platform that integrates heterogeneous data, accelerates model calibration, and supports rational design and optimization of mRNA manufacturing processes.
2602.07257Feb 2026
View
Recurrent neural chemical reaction networks trained to switch dynamical behaviours through learned bifurcations
Both natural and synthetic chemical systems not only exhibit a range of non-trivial dynamics, but also transition between qualitatively different dynamical behaviours as environmental parameters change. Such transitions are called bifurcations. Here, we show that recurrent neural chemical reaction networks (RNCRNs), a class of chemical reaction networks based on recurrent artificial neural networks that can be trained to reproduce a given dynamical behaviour, can also be trained to exhibit bifurcations. First, we show that RNCRNs can inherit some bifurcations defined by smooth ordinary differential equations (ODEs). Second, we demonstrate that the RNCRN can be trained to infer bifurcations that allow it to approximate different target behaviours within different regions of parameter space, without explicitly providing the bifurcation itself in the training. These behaviours can be specified using target ODEs that are discontinuous with respect to the parameters, or even simply by specifying certain desired dynamical features in certain regions of the parameter space. To achieve the latter, we introduce an ODE-free algorithm for training the RNCRN to display designer oscillations, such as a heart-shaped limit cycle or two coexisting limit cycles.
2602.02374Feb 2026
View
RAG-GNN: Integrating Retrieved Knowledge with Graph Neural Networks for Precision Medicine
Network topology excels at structural predictions but fails to capture functional semantics encoded in biomedical literature. We present a retrieval-augmented generation (RAG) embedding framework that integrates graph neural network representations with dynamically retrieved literature-derived knowledge through contrastive learning. Benchmarking against ten embedding methods reveals task-specific complementarity: topology-focused methods achieve near-perfect link prediction (GCN: 0.983 AUROC), while RAG-GNN is the only method achieving positive silhouette scores for functional clustering (0.001 vs. negative scores for all baselines). Information-theoretic decomposition shows network topology contributes 77.3% of predictive information, while retrieved documents provide 8.6% unique information. Applied to cancer signaling networks (379 proteins, 3,498 interactions), the framework identifies DDR1 as a therapeutic target based on retrieved evidence of synthetic lethality with KRAS mutations. These results establish that topology-only and retrieval-augmented approaches serve complementary purposes: structural prediction tasks are solved by network topology alone, while functional interpretation uniquely benefits from retrieved knowledge.
2602.00586Jan 2026
View
BioModelsRAG: A Biological Modeling Assistant Using RAG (Retrieval Augmented Generation)
The BioModels database is one of the premier databases for computational models in systems biology. The database contains over 1000 curated models and an even larger number of non-curated models. All the models are stored in the machine-readable format, SBML. Although SBML can be translated into the human readable Antimony format, analyzing the models can still be time consuming. In order to bridge this gap, a LLM (large language model) assistant was created to analyze the BioModels and allow interaction between the user and the model using natural language. By doing so, a user can easily and rapidly extract the salient points in a given model. Our analysis workflow involved 'chunking' BioModels and converting them to plain text using llama3, and then embedding them in a ChromaDB database. The user-provided query was also embedded, and a similarity search was performed between the query and the BioModels in ChromaDB to extract the most relevant BioModels. The BioModels were then used as context to create the most accurate output in the chat between the user and the LLM. This approach greatly minimized the chance of hallucination and kept the LLM focused on the problem at hand.
2601.22684Jan 2026
ViewPage 1 of 6
