• Matteo Barberis

The community of special interest (COSI) in systems modeling (SysMod) organizes annual one-day gatherings. In 2025 the meeting comprises three sessions that cover a broad variety of topics, beginning with metabolic modeling, followed by the afternoon session on multiscale modeling and concludes with inference of cellular processes. This year's meeting features two keynote speakers, Ronan Fleming and Jasmin Fisher. The event is hosted by Chiara Damiani and Matteo Barberis on behalf of the eight COSI organizers. This brief talk introduces all speakers, organizers, and main topics of the 2025 meeting.

• Ronan Fleming

Established genome-scale modelling methods primarily predict reaction fluxes while established high-throughput experimental technologies primarily measure molecular species concentrations. This paradoxical situation has arisen because of the problem of implementing the nonlinear constraints that represent reaction kinetic rate equations without recourse to convenient yet inaccurate approximations. We present a mathematically and computationally tractable solution to this problem. First we introduce a mathematical reformulation of established knowledge on metabolic reactions and reaction kinetics in matrix vector notation. Then we introduce variational kinetics, a novel approach to satisfy steady-state reaction kinetics at genome scale using novel mathematical and numerical optimisation techniques. We illustrate how this approach may be used to efficiently optimise over the set of steady-state reaction fluxes, thermodynamically feasible kinetic parameters and kinetically feasible elementary or macroscopic rate laws.

• Lisa Corbeij, Department of Epidemiology and Data Science, Amsterdam UMC, The Netherlands, Netherlands
• Natal van Riel, Dept. Biomedical Engineering, Eindhoven University of Technology, the Netherlands., Netherlands
• Shauna O'Donovan, Dept. Biomedical Engineering, Eindhoven University of Technology, the Netherlands, Netherlands

Genome-scale metabolic models (GEMs) are large network-based metabolic reconstructions that can predict the flux of numerous metabolites making them valuable for analysing metabolism across a wide variety of human tissues and microbial species. However, the steady-state assumption needed to solve these GEMs limits their utility to study disturbances in metabolic resilience..
In this study, we embed GEMs of the liver, skeletal muscle, and adipocyte into the Mixed Meal Model (MMM), a physiology-based computational model describing the interplay between glucose, insulin, triglycerides and non-esterified fatty acids (NEFAs). We implement dynamically updating objective functions for each GEM, where cellular objective depends on the model-calculated insulin values. The MMM is simulated using a fixed-step ODE solver; at each time-step the exchange reaction bounds for glucose, triglycerides, and NEFAs in each GEM are updated according to the MMM outputs and flux balance analysis is used to determine the metabolic fluxes. The insulin-dependent objective function allowed the GEMs to accurately simulate the transition from glucose and NEFA secretion in the fasting state to nutrient storage post meal. Moreover, the dynamic tissue-specific GEMs also correctly simulated postprandial changes in metabolites like lactate and glycerol, that are not directly modulated by the differential equations. Personalised hybrid multi-tissue Meal Models, derived from meal response data reveal changes in tissue-specific flux associated with insulin resistance and liver fat accumulation. This research demonstrates the potential of merging GEMs with physiological models to deepen our understanding of metabolic dynamics, offering promising avenues for personalized medicine in metabolic disorders.

• Garhima Arora, Translational Health Science and Technology Institute, Faridabad, India, India
• Samrat Chatterjee, Translational Health Science and Technology Institute, Faridabad, India, India

Breast organotropism, the preferential metastasis of breast cancer cells to specific organs, remains a critical challenge but clinically significant phenomenon, with a limited understanding of the metabolic factors driving site-specific colonization. In this study, we employed genome-scale metabolic models (GSMMs) integrated with single-cell RNA sequencing data from patient-derived xenograft models to investigate the metabolic basis of breast cancer organotropism. We constructed 14 tissue-specific metabolic models from primary breast tumors and their corresponding metastatic sites in the liver, bone, and brain and systematically explored metabolic perturbations associated with disease progression. Our analysis revealed distinct metabolic adaptations in metastatic tissues, characterized by upregulation in lipid metabolism, vitamin and cofactor metabolism, and amino acid pathways, particularly in bone and brain metastases compared to the liver. Furthermore, flux-based comparisons of primary tumors predisposed to different metastatic destinations identified metabolic signatures predictive of organotropism. Using robust Metabolic Transformation Algorithm (rMTA), we simulated gene over-expression and knock-out strategies, identifying candidate metabolic genes capable of driving primary to metastatic phenotypes during breast organotropism. This systems-level approach not only advances our understanding of the metabolic determinants of breast cancer organotropism but also highlights potential metabolic targets for therapeutic intervention aimed at halting metastatic progression.

• Sultana Al Zubaidi, United Arab Emirates University, United Arab Emirates
• Muhammad Ibtisam Nasar, United Arab Emirates University, United Arab Emirates
• Richard Notebaart, Wageningen University and Research, Netherlands
• Markus Ralser, Charité – Universitätsmedizin Berlin, Germany
• Mohammad Tauqeer Alam, United Arab Emirates University, United Arab Emirates

The metabolic network, the largest inter connected system in the cell, is constantly regulated by a range of regulatory interactions. To characterize metabolite-enzyme activatory interactions we reconstructed the cell-intrinsic enzyme-metabolite activation-interaction network ("activation network") using the Saccharomyces cerevisiae metabolic model (Yeast9) as the basis. We integrated the Yeast9 metabolic network with the list of cross-species activating compounds from the BRaunschweig ENzyme DAtabase (BRENDA) database. The cell-intrinsic activation network comprises 1,499 activatory interactions involving 344 enzymes and 286 cellular metabolites. Although only 54% of yeast metabolic enzymes (344 out of 635) are intracellularly activated, these enzymes are distributed across nearly all pathways, underscoring the widespread role of activation in cellular metabolism. Notably, in 94% of pathways at least one initial reaction is intracellularly activated. These initial reactions are typically non-equilibrium, flux-generating steps that must be regulated to control the overall pathway flux. Moreover, our analysis shows that highly activating metabolites are predominantly essential, whereas highly activated enzymes tend to be non-essential for growth. Additionally, we find that activator metabolites are produced in fewer steps compared to non-activators, suggesting a streamlined synthesis for regulatory compounds. We further examined cross-pathway activation and found a significant degree of trans-activation, emphasizing the interconnected nature of cellular metabolism. This coordination ensures that metabolic pathways are selectively activated and dynamically adjusted to meet cellular demands.

• Bingxin Lu, University of Surrey, United Kingdom
• Samuel Winnall, University College London, United Kingdom
• Will Cross, University of Reading, United Kingdom
• Chris Barnes, University College London, United Kingdom

Chromosomal instability (CIN) is pervasive in human tumours and often leads to structural or numerical chromosomal aberrations. Somatic structural variants (SVs) are intimately related to copy number alterations but the two types of variant are often studied independently. Additionally, despite numerous studies on detecting various SV patterns, there are still no general quantitative models of SV generation. To address this issue, we develop a computational cell-cycle model for the generation of SVs from end-joining repair and replication after double-strand break formation. Our model provides quantitative information on the relationship between breakage fusion bridge cycle, chromothripsis, seismic amplification, and extra-chromosomal circular DNA. Given whole-genome sequencing data, the model also allows us to infer important parameters in SV generation with Bayesian inference. Our quantitative framework unifies disparate genomic patterns resulted from CIN, provides a null mutational model for SV, and reveals deeper insights into the impact of genome rearrangement on tumour evolution.

• Jasmin Fisher, University College London Cancer Institute, United Kingdom

Cancer is a complex systemic disease driven by genetic and epigenetic aberrations that impact a multitude of signalling pathways operating in different cell types. The dynamic, evolving nature of the disease leads to tumour heterogeneity and an inevitable resistance to treatment, which poses considerable challenges for the design of therapeutic strategies to combat cancer. In this talk, I will discuss some of the progress made towards addressing these challenges, using the design of computational models of cancer signalling programs (i.e., virtual tumours). I will showcase a growing library of mechanistic, data-driven computational models, focused on the intra- and inter-cellular signalling in various types of cancer (namely triple-negative breast cancer, non-small cell lung cancer, melanoma and glioblastoma). These computational models are predictive and mechanistically interpretable, enabling us to understand and anticipate emergent resistance mechanisms and to design patient-specific treatment strategies to improve outcomes for patients with hard-to-treat cancers.

• Thaleia Ntiniakou, Barcelona Supercomputing Center (BSC), Spain
• Othmane Hayoun-Mya, Barcelona Supercomputing Center (BSC), Spain
• Marco Ruscone, Barcelona Supercomputing Center (BSC), Spain
• Alejandro Madrid Valiente, Barcelona Supercomputing Center (BSC), Spain
• Adam Smelko, Barcelona Supercomputing Center (BSC), Spain
• Jose Luis Estragués Muñoz, Barcelona Supercomputing Center (BSC), Spain
• Jose Carbonell-Caballero, Barcelona Supercomputing Center (BSC), Spain
• Alfonso Valencia, Barcelona Supercomputing Centre (BSC), Spain
• Arnau Montagud, Institute for Integrative Systems Biology (I2SysBio), CSIC-UV, Spain
• and the PerMedCoE multiscale benchmark communityCentre (BSC), (PerMedCoE)

The emergence of virtual human twins (VHT) in biomedical research has sparked interest in multiscale modelling frameworks, particularly in their application bridging cellular to tissue levels. Among these tools, off-lattice agent-based models (O-ABM) offer a promising approach due to their depiction of cells in 3D space, closely resembling biological reality. Despite the proliferation of O-ABM tools addressing various biomedical challenges, comprehensive and systematic comparisons among them have been lacking.
This paper presents a community-driven benchmark initiative aimed at evaluating and comparing O-ABM for biomedical applications, akin to successful efforts in other scientific domains such as CASP. Enlisting developers from leading tools like BioDynaMo, Chaste, PhysiCell, and TiSim, we devised a benchmark scope, defined metrics, and established reference datasets to ensure a meaningful and equitable evaluation. Unit tests targeting different solvers within these tools were designed, ranging from diffusion and mechanics to cell cycle simulations and growth scenarios. Results from these tests demonstrate varying tool performances in handling diffusion, mechanics, and cell cycle equations, emphasising the need for standardised benchmarks and interoperability.
Discussions among the community underscore the necessity for defining gold standards, fostering interoperability, and drawing lessons from analogous benchmarking experiences. The outcomes, disseminated through a public platform in collaboration with OpenEBench, aim to catalyse advancements in computational biology, offering a comprehensive resource for tool evaluation and guiding future developments in cell-level simulations. This initiative endeavours to strengthen and expand the computational biology simulation community through continued dissemination and performance-oriented benchmarking efforts to enable the use of VHT in biomedicine.

• Nhung Duong, Hanoi University of Pharmacy, Viet Nam
• Tuan Do, N2TP Technology Solutions JSC., Viet Nam
• Tien Nguyen, Hanoi University of Pharmacy, Viet Nam
• Hoa Vu, Hanoi University of Pharmacy, Viet Nam
• Lap Nguyen, Hanoi University of Pharmacy, Viet Nam

Background and Aims: JAK/STAT cancer pathway-oriented treatment optimization poses challenges due to pathway complexity, feedback loops, resistance development, and harmonization between efficacy, immunity preservation, and toxicity minimization. We aimed to develop and utilize a comprehensive in silico framework based on multi-objective reinforcement learning (MORL) to identify optimal intervention strategies targeting the JAK/STAT pathway.
Methods: We constructed a multi-scale mechanistic model integrating JAK/STAT intracellular signaling dynamics (ODEs), tumor-immune cell interactions, adaptive resistance evolution, and pharmacokinetics/pharmacodynamics/toxicity of inhibitors (JAKi, STATi, cytokine blockers). We trained MORL agents (multi-objective PPO) employing this model as an environment to discover treatment schedules balancing four objectives: tumor reduction, immune preservation, resistance prevention, and toxicity minimization. Pareto-optimal strategies were used.
Results: MORL successfully identified a diverse set of non-dominated intervention strategies, revealing inherent trade-offs. Distinct treatment paradigms emerged, such as an efficacy-focused strategy yielding ~28% tumor reduction but incurring higher toxicity/resistance, contrasted with a resistance-prevention strategy achieving excellent resistance/toxicity scores (>0.94, >0.25) but limited tumor control (~6.4%). Sensitivity analyses highlighted SOCS3 regulation, STAT kinetics, and resistance parameters as critical determinants of outcomes.
Conclusion: This study introduces a MORL-empowered framework for navigating complex therapeutic trade-offs in JAK/STAT targeting. Our findings reveal diverse optimal strategies and underscore key biological factors influencing treatment success, offering a computational basis for the rational design and personalization of JAK/STAT-targeted therapies and showcasing the potential of MORL in quantitative systems pharmacology.

• Bi-Rong Wang, Biomedical Engineering, Eindhoven University of Technology (TU/e); Institute for Complex Molecular Systems, TU/e, Netherlands
• Maaruthy Yelleswarapu, Biomedical Engineering, Eindhoven University of Technology (TU/e); Institute for Complex Molecular Systems, TU/e, Netherlands
• Federica Eduati, Biomedical Engineering, Eindhoven University of Technology (TU/e); Institute for Complex Molecular Systems, TU/e, Netherlands

Despite advances in immunotherapy, pancreatic cancer remains highly lethal due to an immunosuppressive tumor microenvironment characterized by immune exclusion. Inducing chemokines like CXCL9 with targeted therapy can promote cytotoxic T cell recruitment. However, systematic approaches are needed to better understand chemokine regulation and identify drug combinations upregulating CXCL9 expression. We combined logic-ODE models with wet lab experiments for two pancreatic cancer cell lines (AsPC1 and BxPC3), investigating CXCL9 responses to combinations of cytokines (IFNγ, TNFα) and inhibitors targeting JAK, IKK, RAS, MEK, or PI3K. Analyzing model parameters fitted to the screening data revealed both shared and cell line-specific mechanisms. To efficiently navigate the large space of possible drug combinations, we developed a pipeline for experimental design, integrating active learning with mechanistic modeling. In this pipeline, an acquisition function iteratively selects the most informative conditions from a pool of unseen drug combinations; these are added to the training data to update the model. We benchmarked different acquisition functions on in silico data: “greedy” selects conditions predicted to yield the highest CXCL9 levels, while “uncertainty” prioritizes those where the model is least confident. Both strategies outperformed random sampling: “greedy” most efficiently identified high-CXCL9 conditions, while “uncertainty” improved overall model generalizability. We are currently performing wet lab experiments to validate in silico predictions. Our framework demonstrates how active learning can be combined with dynamic logic-based models to accelerate the discovery of immunomodulatory drug combinations, offering a generalizable approach for hypothesis-driven experimental design in systems biology.

• Sayali Anil Alatkar

Cellular processes like development, differentiation, and disease progression are highly complex and dynamic (e.g., gene expression). These processes often undergo cell population changes driven by cell birth, proliferation, and death. Single-cell
sequencing enables gene expression measurement at the cellular resolution, allowing us to decipher cellular and molecular dynamics underlying these processes. However, the high costs and destructive nature of sequencing restrict observations to snapshots of unaligned cells at discrete timepoints, limiting our understanding of these processes and complicating the reconstruction of cellular trajectories. To address this challenge, we propose ARTEMIS, a generative model integrating a variational autoencoder (VAE) with unbalanced Diffusion Schrödinger Bridge (uDSB) to model cellular processes by reconstructing cellular trajectories, reveal gene expression dynamics, and recover cell population changes. The VAE maps input time-series single-cell data to a continuous latent space, where trajectories are reconstructed by solving the Schrödinger bridge problem using forward-backward stochastic differential equations (SDEs). A drift function in the SDEs captures deterministic gene expression trends. An additional neural network estimates time-varying kill rates for single cells along trajectories, enabling recovery of cell population changes. Using three scRNA-seq datasets—pancreatic β-cell differentiation, zebrafish embryogenesis, and epithelial-mesenchymal transition (EMT) in cancer cells—we demonstrate that ARTEMIS: (i) outperforms state-of-art methods to predict held-out timepoints, (ii) recovers relative cell population changes over time, and (iii) identifies “drift” genes driving deterministic expression trends in cell trajectories. Furthermore, in silico perturbations show that these genes influence processes like EMT. The code for ARTEMIS: https://github.com/daifengwanglab/ARTEMIS.

• Ati Ahmadi, PhD candidate at UWaterloo, Canada
• Samantha Schwartz, Student at UWaterloo, Canada
• Brian Ingalls, professor in the Department of Applied Mathematics at the University of Waterloo, Canada

I uploaded the long abstract below.

• Erick Armingol, Wellcome Sanger Institute, United Kingdom
• James Ashcroft, Wellcome Sanger Institute, United Kingdom
• Magda Mareckova, University of Oxford, United Kingdom
• Martin Prete, Wellcome Sanger Institute, United Kingdom
• Valentina Lorenzi, Wellcome Sanger Institute, United Kingdom
• Cecilia Icoresi Mazzeo, Wellcome Sanger Institute, United Kingdom
• Jimmy Tsz Hang Lee, Wellcome Sanger Institute, United Kingdom
• Marie Moullet, Wellcome Sanger Institute, United Kingdom
• Christian Becker, University of Oxford, United Kingdom
• Krina Zondervan, University of Oxford, United Kingdom
• Omer Ali Bayraktar, Wellcome Sanger Institute, United Kingdom
• Luz Garcia-Alonso, Wellcome Sanger Institute, United Kingdom
• Nathan E. Lewis, University of California San Diego; University of Georgia, United States
• Roser Vento-Tormo, Wellcome Sanger Institute, United Kingdom

Metabolism is fundamental to cellular function, supporting macromolecule synthesis, signaling, growth, and cell-cell communication. While single-cell and spatial metabolomics technologies have advanced, large-scale applications remain challenging. In contrast, transcriptomics provides vast datasets to infer metabolic states. Here, we present scCellFie, a computational tool that predicts metabolic activities from transcriptomic data at single-cell and spatial resolutions. scCellFie enables scalable analysis of large cell atlases, leverages metabolic tasks for interpretable results, and includes modules for identifying metabolic markers, condition-specific changes, and cell-cell communication. We applied scCellFie to ~30 million human cells, generating a comprehensive metabolic atlas across organs while demonstrating our tool’s scalability. Additionally, we used scCellFie to study the human endometrium, the uterine lining that undergoes substantial remodeling throughout the menstrual cycle due to sex hormones, and identified cell type-specific metabolic programs supporting cyclical changes. Epithelial cells exhibited metabolic regulation covering pathways supporting proliferation and mitigating oxidative stress. Endometrial diseases, including endometriosis and endometrial carcinoma, often arise from metabolic dysregulation. By inspecting eutopic endometrium from donors with endometriosis, we identified altered metabolic programs that likely drive atypical proliferation and inflammation of the distinct cell types. In endometrial carcinoma, malignant cells displayed metabolic rewiring, including increased glucose-to-lactate conversion and dysregulated kynurenine and estrogen signaling. These shifts suggest shared mechanisms promoting aberrant proliferation and may reveal therapeutic targets. Together, our findings demonstrate scCellFie as a scalable, interpretable tool for characterizing metabolism in health and disease. By linking metabolic functions to cellular processes, scCellFie provides deeper insights into metabolic regulation across diverse biological systems.

• Koichiro Majima, Department of Computational and Systems Biology, Institute of Science Tokyo Medical Research Laboratory, Japan
• Teppei Shimamura, Department of Computational and Systems Biology, Institute of Science Tokyo Medical Research Laboratory, Japan

Single-cell spatial genomics provides unprecedented molecular insights, yet it still struggles to track both the spatial and temporal progression of tissues under native conditions. Experimental constraints and destructive sampling yield discrete snapshots rather than the continuous record required to fully understand how cells organize and function over time. Optimal transport (OT) approaches have attempted to bridge these snapshots across different assays, but they typically rely on unimodal data, scale poorly, and oversimplify the complexity of ongoing morphogenetic events.
We introduce a spatiotemporal variational autoencoder (VAE) that models the continuous evolution of tissue pixels—capturing dynamic changes in both spatial location and gene-expression patterns. By embedding these pixels into a latent space governed by a learned dynamics network, our method reveals how tissues grow, reorganize, and express key genes over time. At each time point, behavioral parameters are decoded from the latent state via a neural decoder, quantifying probabilities of growth, disappearance, seeking, displacement, attraction, and clustering. A differentiable growth module further refines this process by modeling region appearance and disappearance, allowing gradient-based optimization of tissue occupancy patterns.
We demonstrate the power of our approach by analyzing a mouse embryogenesis dataset. The model uncovers unobserved developmental trajectories, pinpoints morphogenetic transitions, and aligns coronal sections at scale—capabilities that standard OT-based methods find intractable. These results highlight how a spatiotemporal VAE can reconstruct the story of tissue formation in both forward and backward directions, opening up new avenues to interpret single-cell and spatial data as a cohesive, dynamic narrative rather than disjointed snapshots.

• Marek Kimmel, Rice University, United States
• Marie Doumic, INRIA, France
• Leonard Mauvernay, Ecole Polytechnique, France
• Teresa Teixeira, Institut de Biologie Physico-Chimique, France

We discuss a stochastic model of growth of a cell population of cultured yeast cells with gradually decaying chromosome endings called the telomeres, as well as models of telomere reconstruction using the so-called ALT (alternate lengthening of telomeres) Mechanism. Telomeres play a major role in aging and carcinogenesis in humans. Our models correspond in part to the experiments of one of us (Teixeira). For telomere shortening, we modify the method of Olofsson and Kimmel, who considered properties of the branching process of telomere shortening; this leads to consideration of a random walk on a two-dimensional grid. We derive an integral equation for the probability generating functions (pgf’s) characterizing the dynamics of shortening of telomeres. We find that the general solutions have the form of exponential polynomials. Stochastic simulations lead to interesting and non-obvious effects if cell death is included. We further consider more complex models, involving cell death and the ALT mechanism. of telomere reconstruction. These are based on our works and are intended to address the experiments in Kockler et al. (2021). In one version of the ALT Mechanism model, we consider expectations conditional on non-extinction, since only a fraction of ALT telomeres is stably elongated (see Kockler et al. 2021). As a conclusion, the multitype branching processes produce realistic prediction concurrent with complex biological experiments involving telomeres. Our aim is to use the models for longer-term prognoses for human telomeres.

• Jiachen Li, Institute of Process Engineering, Chinese Academy of Sciences, China
• Xiaoyong Pan, Shanghai Jiao Tong University, China
• Ye Yuan, Institute of Process Engineering, Chinese Academy of Sciences, China
• Hong-Bin Shen, Shanghai Jiao Tong University, China

RNA velocity is closely related with cell fate and is an important indicator for the prediction of cell states with elegant physical explanation derived from single-cell RNA-seq data. Most existing RNA velocity models aim to extract dynamics from the phase delay between unspliced and spliced mRNA for each individual gene. However, unspliced/spliced mRNA abundance may not provide sufficient signal for dynamic modeling, leading to poor fit in phase portraits. Motivated by the idea that RNA velocity could be driven by the transcriptional regulation, we propose TFvelo, which expands RNA velocity concept to various single-cell datasets without relying on splicing information, by introducing gene regulatory information. Our experiments on synthetic data and multiple scRNA-Seq datasets show that TFvelo can accurately fit genes dynamics on phase portraits, and effectively infer cell pseudo-time and trajectory from RNA abundance data. TFvelo opens a novel, robust and accurate avenue for modeling RNA velocity for single cell data.

• Chiara Damiani

This concluding talk aims to briefly discuss the diversity of topics presented at the “Computational Modeling of Biological Systems” (SysMod) COSI track. This diversity illustrates the importance of the field and the broad range of applications in systems biology and disease. Then, forthcoming meetings of interest will be announced, and the three poster awards will be delivered as a closing event.