Genetically identical cells can adopt distinct, stable states, playing a crucial role in development and tissue organisation. This talk explores the mechanisms driving cell fate decisions, focusing on the interplay between gene regulatory networks and cell-to-cell communication. Using spatial stochastic models that capture fine-scale regulatory dynamics, we demonstrate how feedback loops and paracrine signalling function as switch-like controllers of cell fate, enabling transitions from transient to stable states. We derive mathematical expressions predicting the threshold signalling strength required to trigger phase transitions and establish a fundamental limit on the spatial spread of phenotypic regions. Specifically, we show that the mean region size scales proportionally to the cubic root of the signalling strength, implying that large, stable domains are prohibitively costly to maintain. This trade-off between robustness and signalling precision highlights the constraints organisms must navigate during development to maintain spatial organisation. Our findings provide key insights into the principles governing multicellular patterning and the regulation of tissue structure.

Multiscale problems are prevalent in many real world scenarios, especially in biology, where the behaviour of a single microorganism can have considerable impact over lengthscales much larger than its own. In this work we consider the effect of the membrane microstructure of a bacterial cell on the behaviour of concentration profiles of relevant molecules on bacterial and bacterial colony lengthscales. Transport through the outer membrane of gram-negative bacteria is restricted to specific channels and non-specific porins. These provide a size-restricted passageway for small molecules through an otherwise impermeable membrane. The effects of these channels are important, for example, antibiotics must cross the outer membrane in order to effectively target gram-negative bacteria, and quorum sensing molecules must cross the membrane to allow bacterial colonies to coordinate mass phenotypic changes such as the production of virulence factors. In mathematical models this limiting transport mechanism across the membrane is often represented via phenomenological constitutive boundary conditions. In this work, we systematically derive the correct effective boundary conditions to impose across a bacterial membrane in terms of physical channel and porin properties. We use a hybrid mathematical approach, combining multiscale methodology such as asymptotic homogenisation and boundary layer theory with numerical simulations. More broadly, because we consider a generic membrane geometry and do not impose a specific outer problem, the results that we derive have a wide scope of potential applications beyond bacterial membranes, for example, to model water vapour or heat loss through fabrics, or mass transfer through surface coatings in chemical engineering.

Many genes are expressed in bursts of transcription, associated with alternating active and inactive promoter states. Such transcriptional bursting is characterised by the burst frequency and burst size, which describe how often a burst occurs and how many transcripts are produced per burst. These two burst parameters offer a simple, intuitive and practical quantitative description of bursty gene expression dynamics. However, a transcriptome-wide picture of how the burst frequency and size are modulated due to gene replication and other cell-cycle dependent factors remains missing. To address this, we fit mechanistic models of gene expression to mRNA count data for thousands of mouse genes, obtained by sequencing of single cells whose cell-cycle position has been inferred previously. Although we observe substantial heterogeneity in transcriptional regulation, we find that upon DNA replication, the genome-wide median burst frequency approximately halves, while the median burst size remains mostly unchanged, thus shedding light on the effect of gene dosage compensation. We show that to accurately estimate the bursting kinetics from sequencing data, mechanistic models must explicitly account for gene copy number variation and extrinsic noise due to factors varying across the cell cycle, whereas correcting for technical noise due to imperfect mRNA capture is less critical.

Pre-mRNA often contains introns, which are non-coding sequences that need to be cut out or spliced before translation occurs. The spliceosome, an essential catalyst composed of several proteins with specific sequence affinities, is required for this process. Very long introns must be removed in pieces, a process known as recursive splicing. The experimental literature on the time it takes for the splicing process to occur is inconsistent. Splicing was traditionally believed to be a slow process that could take anywhere from one to tens of minutes per splicing event. However, recent reports suggest that some splicing events occur within a few tens of seconds. We developed the chemical master equation corresponding to the biochemical mechanism of splicing, allowing us to derive the system’s probability distribution, and perform a stability analysis on two conditions based on an unknown association constant parameter associated with the binding step of the scaffolding complex. We also concluded that the distribution of splice times for a single event ranges from a few tens of seconds to a few tens of minutes. Through sensitivity analyses, we have found that the mean splicing time and distribution are almost entirely dependent on the rate at which the spliceosome is activated in the assembly process—i.e. when the U1 and U4 splicing factors dissociate—which confirms that this is the rate limiting step in the catalytic process. Finally, we have examined the distributions of recursive splicing up to six events, and derived analytic solutions for these recursive splicing events in the case where the scaffolding complex strongly binds to the pre-mRNA complex (the condition thought to favor recursive binding), thus providing a model that can be fit to experimental data to in order to evaluate the number of recursive splicing events occurring.

This work examines the interplay between age-structure and migration dynamics in collective cell behaviour. We focus on the integration of cell cycle dynamics with spatial migration, particularly examining the 'go-or-grow' hypothesis in the context of age-dependent processes. Our framework extends classical travelling wave theory to account for the age structure of cell populations, offering new insights into how cell cycle phases influence moving fronts and invasion dynamics. We analyse wave speed characteristics and front dynamics in age-structured systems, addressing a significant gap in current mathematical biology literature. The research provides a novel theoretical foundation for understanding how cell-cycle dependent proliferation and migration behaviours contribute to collective cell dynamics.

Cells respond to their local environment through mechanotransduction, converting mechanical signals into a biological response (e.g. cell growth, proliferation or differentiation). The cell cytoskeleton, particularly actomyosin stress fibres (SFs), and focal adhesions (FAs), which bind the cytoskeleton to the extra-cellular matrix (ECM), are central to this process, activating intracellular signalling cascades in response to deformation. We present a novel two-dimensional bio-chemo-mechanical model to describe the development of these structures, coupled through a positive feedback loop, and the associated cell deformation. Building on our previous one-dimensional approach, we similarly employ reaction-diffusion-advection equations to describe the evolution of key scaffolding and signalling proteins, and connect their concentrations to a viscoelastic description of the cell cytoplasm, ECM and adhesions. Further, we now incorporate other key mechanotransducing structures including the stiff cell nucleus, and plasma and cortical membranes. Working in an axisymmetric framework, we employ this model to explain how, dependent upon the mechanical properties of the surrounding ECM, non-uniform patterns of cell striation develop, leading to FA and SF localisation at the cell periphery. Moreover, a linear stability analysis reveals the stability of the axisymmetric configuration to various normal modes of deformation. By identifying non-axisymmetric modes with positive growth rates our model demonstrates a possible mechanism for self-driven surface patterning of cells in vitro.

In a number of bacteria, nitric oxide (NO) is converted to nitrate by an enzyme called Hmp. In emph{Streptomyces coelicolor}, synthesis of Hmp is in turn controlled by an iron-sulfur protein called NsrR. NsrR represses the transcription of two copies of the emph{hmp} gene in the emph{S. coelicolor} genome, but reaction of NsrR's iron-sulfur cluster with NO causes NsrR to dissociate from the emph{hmp} promoter, thus allowing Hmp to be expressed. While this is a straightforward control mechanism, NsrR is a dimer, and the iron-sulfur cluster in each monomer of NsrR can react with NO several times. Eventually, a repair system restores the NO-damaged iron-sulfur clusters of the dimers. But given that a single reaction with NO is sufficient to cause the NsrR dimer to dissociate from the emph{hmp} promoter, do we need to model the complex chemistry of the dimer, or is a highly simplified model that considers a single NsrR unit and its iron-sulfur cluster sufficient to capture the dynamics of this control system?

Across the spectrum of cell types, the concentration of calcium controls a wide array of cellular functions. These calcium signals, usually in the form of periodic oscillations, play a paramount role in correct cellular activity. T cells are fundamental to the correct behaviour of the immune system. These cells have recently been shown to exhibit two competing oscillatory mechanisms, depending on the influx of extracellular Ca2+. Ca2+ influx is controlled by two molecules, STIM1 and STIM2. When both STIMs are present, T cells showcase sinusoidal Ca2+ oscillations on a raised baseline, but when one of them is absent the nature of the oscillation changes to a mix of Ca2+ spikes and bursting periods. In this talk, we will present an ODE that attempts to explain how these two molecules control the nature of these oscillations in T cells

The prolactin receptor (PRLR) is a single-pass transmembrane receptor driving pancreatic β-cell proliferation via JAK/STAT signaling activation. This signal transduction pathway enables insulin-secreting β-cells to adapt to metabolic stress; however, the precise mechanisms underlying the pathway’s proliferative effect remain ill-defined. Here we implement a pipeline that uses live-cell fluorescence imaging, reconstitution approaches, and fluorescence correlation spectroscopy (FCS) to inform a mathematical model of PRLR signaling in β-cells and build a quantitative, mechanistic understanding of the signaling network. PRLR signaling is dynamic, involving changes in the spatial organization of signaling molecules. We have observed PRLR undergoing rapid internalization, a behavior that has been shown and modeled in other signaling pathways but has not been considered in a mathematical model of PRLR signaling. Such a model is useful for predicting strategies to modulate β-cell function. PRLR internalization is observed in both our minimal engineered PRLR expression system and in native pancreatic tissue, while FCS and chemigenetic labeling with SNAP-tag confirm the presence of a low concentration plasma membrane pool of PRLR. Our imaging data are used to integrate PRLR trafficking dynamics into an ordinary differential equation (ODE) model of PRLR signaling. We employ the ODE model to test hypotheses targeting how the spatial heterogeneity of PRLR signaling dynamics affects downstream signaling outcomes. Our data underscore the versatility of building a generalizable modeling-imaging framework to quantitatively understand signal transduction in and beyond β-cells.

Collective cell migration plays a crucial role in numerous biological processes, including cancer growth, wound healing, and the immune response. Often, the migrating population consists of cells with various different phenotypes. This study derives a general mathematical framework for modelling cell migration into the micro-environment, which is coarse-grained from an underlying individual-based model that captures some of the dynamics of cell migration that are influenced by the phenotype of the cell, such as: random movement, proliferation, phenotypic transitions, and interactions with the external environment. The resulting model provides a continuum, macroscopic description of cell invasion, which represents the phenotype of the cell as a continuous variable and is much more amenable to simulation and analysis than its individual-based counterpart when considering a large number of phenotypes. The results highlight how phenotypic structuring impacts the spatial and temporal dynamics of cell populations, demonstrating that different environmental pressures and phenotypic transition mechanisms significantly influence invasion patterns.

Cells respond to their local environment through mechanotransduction, converting mechanical signals into a biological response (e.g. cell growth, proliferation or differentiation). The cell cytoskeleton, particularly actomyosin stress fibres (SFs), and focal adhesions (FAs), which bind the cytoskeleton to the extra-cellular matrix (ECM), are central to this process, activating intracellular signalling cascades in response to deformation. We present a novel two-dimensional bio-chemo-mechanical model to describe the development of these structures, coupled through a positive feedback loop, and the associated cell deformation. Building on our previous one-dimensional approach, we similarly employ reaction-diffusion-advection equations to describe the evolution of key scaffolding and signalling proteins, and connect their concentrations to a viscoelastic description of the cell cytoplasm, ECM and adhesions. Further, we now incorporate other key mechanotransducing structures including the stiff cell nucleus, and plasma and cortical membranes. Working in an axisymmetric framework, we employ this model to explain how, dependent upon the mechanical properties of the surrounding ECM, non-uniform patterns of cell striation develop, leading to FA and SF localisation at the cell periphery. Moreover, a linear stability analysis reveals the stability of the axisymmetric configuration to various normal modes of deformation. By identifying non-axisymmetric modes with positive growth rates our model demonstrates a possible mechanism for self-driven surface patterning of cells in vitro.

In a number of bacteria, nitric oxide (NO) is converted to nitrate by an enzyme called Hmp. In emph{Streptomyces coelicolor}, synthesis of Hmp is in turn controlled by an iron-sulfur protein called NsrR. NsrR represses the transcription of two copies of the emph{hmp} gene in the emph{S. coelicolor} genome, but reaction of NsrR's iron-sulfur cluster with NO causes NsrR to dissociate from the emph{hmp} promoter, thus allowing Hmp to be expressed. While this is a straightforward control mechanism, NsrR is a dimer, and the iron-sulfur cluster in each monomer of NsrR can react with NO several times. Eventually, a repair system restores the NO-damaged iron-sulfur clusters of the dimers. But given that a single reaction with NO is sufficient to cause the NsrR dimer to dissociate from the emph{hmp} promoter, do we need to model the complex chemistry of the dimer, or is a highly simplified model that considers a single NsrR unit and its iron-sulfur cluster sufficient to capture the dynamics of this control system?

Across the spectrum of cell types, the concentration of calcium controls a wide array of cellular functions. These calcium signals, usually in the form of periodic oscillations, play a paramount role in correct cellular activity. T cells are fundamental to the correct behaviour of the immune system. These cells have recently been shown to exhibit two competing oscillatory mechanisms, depending on the influx of extracellular Ca2+. Ca2+ influx is controlled by two molecules, STIM1 and STIM2. When both STIMs are present, T cells showcase sinusoidal Ca2+ oscillations on a raised baseline, but when one of them is absent the nature of the oscillation changes to a mix of Ca2+ spikes and bursting periods. In this talk, we will present an ODE that attempts to explain how these two molecules control the nature of these oscillations in T cells

The prolactin receptor (PRLR) is a single-pass transmembrane receptor driving pancreatic β-cell proliferation via JAK/STAT signaling activation. This signal transduction pathway enables insulin-secreting β-cells to adapt to metabolic stress; however, the precise mechanisms underlying the pathway’s proliferative effect remain ill-defined. Here we implement a pipeline that uses live-cell fluorescence imaging, reconstitution approaches, and fluorescence correlation spectroscopy (FCS) to inform a mathematical model of PRLR signaling in β-cells and build a quantitative, mechanistic understanding of the signaling network. PRLR signaling is dynamic, involving changes in the spatial organization of signaling molecules. We have observed PRLR undergoing rapid internalization, a behavior that has been shown and modeled in other signaling pathways but has not been considered in a mathematical model of PRLR signaling. Such a model is useful for predicting strategies to modulate β-cell function. PRLR internalization is observed in both our minimal engineered PRLR expression system and in native pancreatic tissue, while FCS and chemigenetic labeling with SNAP-tag confirm the presence of a low concentration plasma membrane pool of PRLR. Our imaging data are used to integrate PRLR trafficking dynamics into an ordinary differential equation (ODE) model of PRLR signaling. We employ the ODE model to test hypotheses targeting how the spatial heterogeneity of PRLR signaling dynamics affects downstream signaling outcomes. Our data underscore the versatility of building a generalizable modeling-imaging framework to quantitatively understand signal transduction in and beyond β-cells.

Collective cell migration plays a crucial role in numerous biological processes, including cancer growth, wound healing, and the immune response. Often, the migrating population consists of cells with various different phenotypes. This study derives a general mathematical framework for modelling cell migration into the micro-environment, which is coarse-grained from an underlying individual-based model that captures some of the dynamics of cell migration that are influenced by the phenotype of the cell, such as: random movement, proliferation, phenotypic transitions, and interactions with the external environment. The resulting model provides a continuum, macroscopic description of cell invasion, which represents the phenotype of the cell as a continuous variable and is much more amenable to simulation and analysis than its individual-based counterpart when considering a large number of phenotypes. The results highlight how phenotypic structuring impacts the spatial and temporal dynamics of cell populations, demonstrating that different environmental pressures and phenotypic transition mechanisms significantly influence invasion patterns.

Climate change is altering the mean, variance and autocorrelation of temperature. However, linear approaches to incorporating these temperature impacts in simple population models do not provide realistic predictions regarding climate change impacts. For example, simple degree days approaches or using a linear function of temperature to alter the density-independent population growth rate will not account for the sometimes catastrophic decrease in performance with high temperatures. We use an extremely simple population model coupled to nonlinear thermal performance curves to explore the simultaneous impact of changes to temperature mean, variance and autocorrelation. The realized density-independent population growth rate is given by three types of thermal performance curves that correspond to published data. We find relatively small impacts on established population dynamics when realistic changes in temperature sequences are used, suggesting that many populations may be quite robust to temperature-driven climate change impacts in the near term. The most extreme right-skewed performance curves are most likely to result in species extinctions, even though these curves have higher optimal temperatures.

The oriental fruit fly, Bactrocera dorsalis, is a serious threat to crops and orchards in many places around the World, and in particular in Réunion island, where it was first detected in 2017. Since then, this pest has invaded the whole island and displaced established fruit fly populations. Since Réunion island is a hot spot of diversity, appropriate control tools have to be deployed to eliminate or reduce the wild population. I will present recent results that study the combination of the Sterile Insect Technique, entomopathogen fungi, and also pheromone traps. In particular, we will show how the spatial component and the orchards connectivity can drastically change the releases strategy, as well as the critical amount of sterile insects to release. We discuss (optimal) strategies obtained with our models versus realistic strategies that can actually be developed in the field. Our approach being generic, it can be adapted to other pests and disease vectors, such as mosquitoes. This works stands within the AttracTIS project, funded by Ecophyto 2021-2022.

As parasitoids are the most frequently used biocontrol agents, especially in agriculture and forest ecosystems, they have become a cornerstone in mathematical biology. They are also a prototypical example of discrete-time systems. Here we consider a simple host-parasitoid model that is based on the classical Nicholson-Bailey model, but includes two extensions that are ecologically plausible: (1) density-dependent host growth (of Beverton-Holt type) and (2) a functional response of type III. The latter can be caused by a number of ecological mechanisms and is key in driving a rich dynamical behavior. While the system admits at most one nontrivial fixed point, we observe up to four coexisting non-equilibrium attractors. They can be periodic, quasi-periodic, or chaotic. They emerge in a quasi-periodic route to chaos and exhibit frequency-locking phenomena. We find different regular organized structures in the two-dimensional parameter plane that describe periodic oscillations surrounded by chaos. Among these structures are Arnold tongues (which have been previously reported in related models) and shrimp-shaped domains, which are little known in ecological models. Our results demonstrate that a type III functional response of parasitoids induces many new complex phenomena. While in continuous-time models the type III functional response tends to be stabilizing, in discrete-time models it can have very contrasting effects. The ecological implications are a high sensitivity not only to parameters but also to the initial condition.

The site frequency spectrum (SFS) is a popular summary statistic of genomic data. In population genetics, the SFS has provided a simple means of inferring the rate of adaptation of a population and for distinguishing between neutral evolution and evolution under selection. The rapidly growing amount of cancer genomic data has attracted interest in the SFS of an exponentially growing population. In this talk, we discuss recent results on the expected value of the SFS of a population that grows according to a stochastic branching process, as well as (first-order) almost sure convergence results for the SFS in the large-time and large-population limits. Our results show that while the SFS depends linearly on the mutation rate, the branching process parameters of birth and death control the fundamental shape of the SFS at the low-frequency end. For the special case of a birth-death process (binary branching process), our results give rise to statistically consistent estimators for the mutation rate and extinction probability of the population, which stands in contrast to previous work which has indicated the need for additional data to decouple these two parameters. Overall, our work shows how single timepoint data on the SFS of an exponentially growing population can be used to infer important evolutionary parameters.

Wildfires disrupt ecosystems, with climate change exacerbating vulnerability in regions poorly adapted to such disturbances. These events are driven by complex, multi-scale interactions where small perturbations in environmental factors can trigger large-scale shifts, complicating prediction efforts. We propose a coupled convection-reaction-diffusion system as a framework for modelling wildfire spread dynamics. This system integrates spatial and temporal variability to identify thresholds for spread and quantify the impact of abrupt environmental changes on burnt areas and rates of propagation. Incorporating environmental, meteorological, and historical fire record data from the Global Wildfire Information System, the Department for Environment, Food and Rural Affairs (UK), and drone footage of heather burning. Bayesian inference and Monte Carlo methods are employed for parameter estimation and uncertainty quantification, ensuring robust model validation against unseen data. Recent wildfire events around the globe highlight the need for actionable insights into environmental vulnerability, property loss, and infrastructure risk. By enabling near-real-time simulations, this model aims to provide a computational tool for emergency response, long-term management strategies, and assessments of climate change-induced outlier weather patterns influencing fire behaviour. This work highlights the potential of mathematical modelling to advance understanding and management of critical ecological disturbances.

Two fundamental quantities in population biology, the reproductive number R0 and the growth rate, are intimately linked, but the exact nature of their relationship is somewhat obscure. Models of microbial growth typically have R0=2, but estimating their growth rate, and hence fitness, requires solving the famous Euler-Lotka equation. Conversely, in epidemiology one typically measures how quickly the infected population grows, but it is the reproductive number R0 that sets the threshold for an epidemic breakout and for herd immunity. In this talk, we use statistical techniques based on large deviations theory to clarify how exactly the population growth rate and R0 are connected. Building an analogy to classical thermodynamics, we show that the long-term behaviour of a population is encoded in a single convex function that relates growth rate, R0, and the statistics of intergeneration times in lineages. As an application, we derive a general formulation of the Euler-Lotka equation and explain why it is almost always appears as an implicit equation.

In the 90s, several experiments suggested a hypothesis that certain genes function to mask or buffer the effects of mutations, thereby allowing them to accumulate and be stored. These were termed evolutionary capacitors and addressed the fundamental evolutionary problem of how populations optimize fitness in one environment while maintaining variation to adapt to another. However, more recent experiments support an alternative hypothesis that such buffering of mutations is a natural and unsurprising outcome of epistasis and the mutation-selection process. To quantitatively test this hypothesis, we develop a mathematical framework that extends a classical partial differential equation of the mutation-selection process to account for epistasis. Using a perturbation method on steady state solutions, we show that certain types of epistatic interactions and selection pressures will lead to the emergence of the evolutionary capacitance phenomena.

Methane, a potent greenhouse gas (GHG), is now driving climate change at an unprecedented rate. With a warming potential greater than carbon dioxide, it poses a substantial threat to the functioning of ecosystems. Despite its importance, studies investigating its direct impact on species interactions within ecosystems are rare. This growing concern highlights the need for a comprehensive understanding of the factors that could disrupt food chains, ultimately impacting ecosystem stability and resilience. In this talk, I will address this gap by developing a mechanistic model that integrates methane dynamics with the populations of species and detritus. This novel approach offers a framework for understanding how gaseous pollutants like methane influence trophic interactions. The model is studied for a range of concentrations of methane. Our findings reveal that low concentrations of methane can benefit species growth as an alternative carbon source. However, moderate to high levels induce sub-lethal to lethal effects. Further, analyzing the mechanisms for long transients in the fast-intermediate-slow formulation of the model, I will discuss how faster methane accumulation in water can result in slower species growth.

Microbially induced corrosion (MIC) is a significant global issue impacting infrastructure, economies, and environment. In wastewater systems, MIC is primarily associated with biofilm formation on concrete sewer pipes, leading to severe degradation due to microbial metabolic activity. The proliferation of sewer biofilms occurs in both submerged and unsubmerged conditions, leading to distinct microbial communities. Commonly, these biofilms host microorganisms such as fermentation bacteria, hydrogen-producing acetogens, denitrifying bacteria, sulfate-reducing bacteria, sulfur-oxidizing bacteria, and methanogens. In particular, sulfur-oxidizing bacteria play a crucial role in corrosion, as they oxidize hydrogen sulfide from wastewater effluents, generating sulfuric acid that accelerates concrete deterioration. A one-dimensional model with double free boundaries has been developed to investigate the proliferation of biofilms and the related corrosion process in wastewater concrete pipes. The domain is composed of two free boundary regions: a biofilm that grows towards the interior cavity of the pipe, sitting on a gypsum layer formed by corrosion, which penetrates the concrete pipe. Diffusion-reaction equations govern the transport and the metabolic production or consumption of dissolved substances, such as hydrogen sulfide, oxygen, and sulfuric acid within the biofilm layer. The biofilm free boundary tracks the growth of the microbial community, regulated by microbial metabolic activity and detachment phenomena. The corrosion process is incorporated into the model through a Stefan-type condition, which drives the advancement of the gypsum free boundary into the concrete pipe, governed by microbial production of sulfuric acid. Numerical simulations have been carried out to investigate the model behavior, encompassing the development and progression of the biofilm as well as the corrosion advancement, with the aim of elucidating the key factors governing both phenomena.

The present paper investigates the dynamical behavior of the modified May Holling Tanner model in the presence of dynamic alternative resources. We study the role of dynamic alternative resources on the survival of the species when there is prey rarity. Detailed mathematical analysis and numerical evaluations, including the situation of ecosystem collapsing, have been presented to discuss the coexistence of species', stability, occurrence of different bifurcations (saddle-node, transcritical, and Hopf) in three cases in the presence of prey and alternative resources, in the absence of prey and in the absence of alternative resources. It has been obtained that the multiple coexisting states and their stability are outcomes of variations in predation rate for alternative resources. Also, the occurrence of Hopf bifurcation, saddle-node bifurcation, and transcritical bifurcation are due to variations in the parameters of dynamic alternative resources. The impact of dynamic alternative resources on species' density reveals the fact that if the predation rate for alternative resources increases, then the prey biomass increases (under some restrictions), and variations in the predator's biomass widely depend upon the quality of food items. This study also points out that the survival of predators is possible in the absence of prey. In the theme of ecological balance, the present study suggests some theoretical points of view for the eco-managers.

The decision of an individual or group to leave its current environment may be influenced by various factors. These include external or inter-species factors such as the presence of predators or food availability, and also intra-species dynamics like mate searching or the strength of social ties within a group. Understanding the consequences of these behaviours on the population level dynamics is non-trivial. In this study, we explore a stochastic model describing the movement of males and females of a species between localised patches, in which the movement rates are dependent on the sex ratio within the patch. By deriving a system of stochastic differential equations governing the fluctuations in these patches we can model a diverse range of intra-species behaviours driven solely by an individual's response to local sex ratio. We subsequently uncover and explore how different individual behaviours can give rise to large scale (anti-)correlated movements between the sexes.

To comprehend the mechanisms driving biodiversity and ecosystem resilience in a rapidly changing world, it is essential to explore the behavioral diversity among individuals in greater depth. Consistent individual differences in behavior, often referred to as animal personality, play a crucial role in shaping ecological and evolutionary dynamics, particularly in foraging behavior. Traditional approaches in behavioral and evolutionary ecology typically focus on average behavior, neglecting the significance of individual variability. This study explores the influence of consumer personality on ecological dynamics, specifically examining how variations in food availability affect behavioral strategies and ecosystem functioning. We develop a resource-consumer model that incorporates personality-dependent saturating attack rates based on the mean-field ratio of resources to consumers. The well-posedness of the model is established, and we analyze the existence and stability of all steady-state solutions. Through bifurcation analysis, we identify critical transition parameters and describe the nonlinear phenomena induced by personality-dependent attack rates. Our findings demonstrate that boldness in consumers enhances their persistence, particularly under low levels of boldness, where populations can survive even with moderate or high food supply, which was not captured in classical frameworks.

A microscale ecological study using the closure approach to understand the impact of productivity controlled by geographical and seasonal variations on the intermittency of phytoplankton will be presented in this talk. Using this approach for a nutrient–phytoplankton model with Holling type III functional response, it has been shown how the dynamics of the system can be affected by the environmental fluctuations triggered by the impact of light, temperature, and salinity, which fluctuate with regional and seasonal variations. Reynold’s averaging method in space, which results in expressing the original components in terms of its mean (average value) and perturbation (fluctuation) has been used to determine the impact of growth fluctuation in phytoplankton distribution and in the intermittency of phytoplankton spreading (variance). Parameters are estimated from the nature of productivity and spread of phytoplankton density during field observation done at four different locations of Tokyo Bay. The model validation shows that our results are in good agreement with the field observation and succeeded in explaining the intermittent phytoplankton distribution at different locations of Tokyo Bay, Japan, and its neighboring coastal regions.

Two fundamental quantities in population biology, the reproductive number R0 and the growth rate, are intimately linked, but the exact nature of their relationship is somewhat obscure. Models of microbial growth typically have R0=2, but estimating their growth rate, and hence fitness, requires solving the famous Euler-Lotka equation. Conversely, in epidemiology one typically measures how quickly the infected population grows, but it is the reproductive number R0 that sets the threshold for an epidemic breakout and for herd immunity. In this talk, we use statistical techniques based on large deviations theory to clarify how exactly the population growth rate and R0 are connected. Building an analogy to classical thermodynamics, we show that the long-term behaviour of a population is encoded in a single convex function that relates growth rate, R0, and the statistics of intergeneration times in lineages. As an application, we derive a general formulation of the Euler-Lotka equation and explain why it is almost always appears as an implicit equation.

In the 90s, several experiments suggested a hypothesis that certain genes function to mask or buffer the effects of mutations, thereby allowing them to accumulate and be stored. These were termed evolutionary capacitors and addressed the fundamental evolutionary problem of how populations optimize fitness in one environment while maintaining variation to adapt to another. However, more recent experiments support an alternative hypothesis that such buffering of mutations is a natural and unsurprising outcome of epistasis and the mutation-selection process. To quantitatively test this hypothesis, we develop a mathematical framework that extends a classical partial differential equation of the mutation-selection process to account for epistasis. Using a perturbation method on steady state solutions, we show that certain types of epistatic interactions and selection pressures will lead to the emergence of the evolutionary capacitance phenomena.

Methane, a potent greenhouse gas (GHG), is now driving climate change at an unprecedented rate. With a warming potential greater than carbon dioxide, it poses a substantial threat to the functioning of ecosystems. Despite its importance, studies investigating its direct impact on species interactions within ecosystems are rare. This growing concern highlights the need for a comprehensive understanding of the factors that could disrupt food chains, ultimately impacting ecosystem stability and resilience. In this talk, I will address this gap by developing a mechanistic model that integrates methane dynamics with the populations of species and detritus. This novel approach offers a framework for understanding how gaseous pollutants like methane influence trophic interactions. The model is studied for a range of concentrations of methane. Our findings reveal that low concentrations of methane can benefit species growth as an alternative carbon source. However, moderate to high levels induce sub-lethal to lethal effects. Further, analyzing the mechanisms for long transients in the fast-intermediate-slow formulation of the model, I will discuss how faster methane accumulation in water can result in slower species growth.

Microbially induced corrosion (MIC) is a significant global issue impacting infrastructure, economies, and environment. In wastewater systems, MIC is primarily associated with biofilm formation on concrete sewer pipes, leading to severe degradation due to microbial metabolic activity. The proliferation of sewer biofilms occurs in both submerged and unsubmerged conditions, leading to distinct microbial communities. Commonly, these biofilms host microorganisms such as fermentation bacteria, hydrogen-producing acetogens, denitrifying bacteria, sulfate-reducing bacteria, sulfur-oxidizing bacteria, and methanogens. In particular, sulfur-oxidizing bacteria play a crucial role in corrosion, as they oxidize hydrogen sulfide from wastewater effluents, generating sulfuric acid that accelerates concrete deterioration. A one-dimensional model with double free boundaries has been developed to investigate the proliferation of biofilms and the related corrosion process in wastewater concrete pipes. The domain is composed of two free boundary regions: a biofilm that grows towards the interior cavity of the pipe, sitting on a gypsum layer formed by corrosion, which penetrates the concrete pipe. Diffusion-reaction equations govern the transport and the metabolic production or consumption of dissolved substances, such as hydrogen sulfide, oxygen, and sulfuric acid within the biofilm layer. The biofilm free boundary tracks the growth of the microbial community, regulated by microbial metabolic activity and detachment phenomena. The corrosion process is incorporated into the model through a Stefan-type condition, which drives the advancement of the gypsum free boundary into the concrete pipe, governed by microbial production of sulfuric acid. Numerical simulations have been carried out to investigate the model behavior, encompassing the development and progression of the biofilm as well as the corrosion advancement, with the aim of elucidating the key factors governing both phenomena.

The present paper investigates the dynamical behavior of the modified May Holling Tanner model in the presence of dynamic alternative resources. We study the role of dynamic alternative resources on the survival of the species when there is prey rarity. Detailed mathematical analysis and numerical evaluations, including the situation of ecosystem collapsing, have been presented to discuss the coexistence of species', stability, occurrence of different bifurcations (saddle-node, transcritical, and Hopf) in three cases in the presence of prey and alternative resources, in the absence of prey and in the absence of alternative resources. It has been obtained that the multiple coexisting states and their stability are outcomes of variations in predation rate for alternative resources. Also, the occurrence of Hopf bifurcation, saddle-node bifurcation, and transcritical bifurcation are due to variations in the parameters of dynamic alternative resources. The impact of dynamic alternative resources on species' density reveals the fact that if the predation rate for alternative resources increases, then the prey biomass increases (under some restrictions), and variations in the predator's biomass widely depend upon the quality of food items. This study also points out that the survival of predators is possible in the absence of prey. In the theme of ecological balance, the present study suggests some theoretical points of view for the eco-managers.

The decision of an individual or group to leave its current environment may be influenced by various factors. These include external or inter-species factors such as the presence of predators or food availability, and also intra-species dynamics like mate searching or the strength of social ties within a group. Understanding the consequences of these behaviours on the population level dynamics is non-trivial. In this study, we explore a stochastic model describing the movement of males and females of a species between localised patches, in which the movement rates are dependent on the sex ratio within the patch. By deriving a system of stochastic differential equations governing the fluctuations in these patches we can model a diverse range of intra-species behaviours driven solely by an individual's response to local sex ratio. We subsequently uncover and explore how different individual behaviours can give rise to large scale (anti-)correlated movements between the sexes.

To comprehend the mechanisms driving biodiversity and ecosystem resilience in a rapidly changing world, it is essential to explore the behavioral diversity among individuals in greater depth. Consistent individual differences in behavior, often referred to as animal personality, play a crucial role in shaping ecological and evolutionary dynamics, particularly in foraging behavior. Traditional approaches in behavioral and evolutionary ecology typically focus on average behavior, neglecting the significance of individual variability. This study explores the influence of consumer personality on ecological dynamics, specifically examining how variations in food availability affect behavioral strategies and ecosystem functioning. We develop a resource-consumer model that incorporates personality-dependent saturating attack rates based on the mean-field ratio of resources to consumers. The well-posedness of the model is established, and we analyze the existence and stability of all steady-state solutions. Through bifurcation analysis, we identify critical transition parameters and describe the nonlinear phenomena induced by personality-dependent attack rates. Our findings demonstrate that boldness in consumers enhances their persistence, particularly under low levels of boldness, where populations can survive even with moderate or high food supply, which was not captured in classical frameworks.

A microscale ecological study using the closure approach to understand the impact of productivity controlled by geographical and seasonal variations on the intermittency of phytoplankton will be presented in this talk. Using this approach for a nutrient–phytoplankton model with Holling type III functional response, it has been shown how the dynamics of the system can be affected by the environmental fluctuations triggered by the impact of light, temperature, and salinity, which fluctuate with regional and seasonal variations. Reynold’s averaging method in space, which results in expressing the original components in terms of its mean (average value) and perturbation (fluctuation) has been used to determine the impact of growth fluctuation in phytoplankton distribution and in the intermittency of phytoplankton spreading (variance). Parameters are estimated from the nature of productivity and spread of phytoplankton density during field observation done at four different locations of Tokyo Bay. The model validation shows that our results are in good agreement with the field observation and succeeded in explaining the intermittent phytoplankton distribution at different locations of Tokyo Bay, Japan, and its neighboring coastal regions.

Pertussis (whooping cough), caused by Bordetella pertussis, is most severe in infants, with most deaths occurring in unvaccinated infants under three months of age. Vaccination with the DTaP (priming) and TdaP (booster) immunizations is effective, with TdaP given during pregnancy and DTaP in infancy. However, immunomodulation can dampen the IgG response in infants born to vaccinated mothers. We hypothesize that adjusting the vaccination schedule could reduce immunomodulation and enhance vaccine efficacy. Since empirically testing multiple schedules is impractical, we propose using mathematical modeling alongside two experimental mouse models to determine an optimal schedule. Pregnant and infant mice will be immunized following a murine analog of standard vaccination. These data will inform our model, allowing us to estimate key immune parameters. Once parametrized, our model will propose schedules that maximize infant antibody response. A second mouse experiment will test these schedules, comparing immune responses to assess their efficacy. This approach will help evaluate immunomodulation mechanisms and refine vaccination strategies. The mechanistic evaluation of immunomodulation is of significance due to its lack of effective investigations to date.

In this talk, we will present effective antibiotic regimens in the presence of drug-resistant bacteria in biofilms. We begin by discussing models of in-vivo antimicrobial resistance transfer within bacterial biofilms, focusing on various one-dimensional biofilm models. Our approach includes modeling resistance acquisition through horizontal gene transfer between resistant and susceptible strains while also accounting for the role of persistor cells. We examine the effects of periodic antibiotic dosing at a constant level, showing that it may not always lead to bacterial eradication. To address this challenge, we utilize a numerical optimization algorithm to determine the optimal antibiotic dosing strategy. Additionally, we analyze how changes in different model parameters impact the qualitative behavior of the optimal dosing regimen.

Axial spondyloarthritis (SpA) is a chronic inflammatory disease impacting the joints in the axial skeleton (e.g. chest, spine, pelvis). Tumor necrosis factor inhibitors (TNFi), such as the monoclonal antibody adalimumab, are used to treat severe cases, but therapy is often discontinued due to loss of efficacy, not least resulting from immunogeniocity and development of anti-drug antibodies (ADA). The disease-modifying drug methotrexate (MTX) has shown potential in reducing the formation of ADA to various TNFi in rheumatoid arthritis (Krickaert et al, 2012), but little is known about its mode of action. We adapt a mechanistic mathematical model for immunogencity towards adalimumab based on Chen et al. 2014 to describe the impact of MTX in reducing immunogenicity, and parametrise it based on patient data from a multicentric randomised trial (Ducoureau et al, 2020). ODEs describe the pharmacokinetics and pharmacodynamics of the therapeutic compounds (adalimumab only or adalimumab and MTX, depending on patient cohort), the dynamics of T and B lymphocytes, antigen presenting cells, and some relevant cytokines for the disease. Due to the large size of this model, we employ several reduced models to estimate some of the parameter values. The model is used to simulate several scenarios in order to elucidate the most likely modes of action of MTX to reduce immunogenicity by comparing the simulated and measured ADA titres along 5 hospital visits. This is joint work with Sara Sottile (Bologna, Italy) and Denis Mulleman (Tours, France). This work is based upon work from COST Action ENOTTA (CA21147), supported by COST (European Cooperation in Science and Technology), and Contract KP-06-DKOST-13 of the Bulgarian Fund for Scientific Research.

Deep learning has emerged in many fields in recent times where neural networks are used to learn and understand data. This study combines deep learning frameworks with epidemiological models and is aimed specifically at the creation and testing of DINNs with a view to modeling the infection dynamics of epidemics. Our research thus trains the DINN on synthetic data derived from an SI-SIR model designed for Avian influenza and shows the model’s accuracy in predicting extinction and persistence conditions. In the method, a twelve hidden layer model was constructed with sixty-four neurons per layer and ReLU activation function was used. The network is trained to predict the time evolution of five state variables for birds and humans over 50,000 epochs. The overall loss minimized to 0.000006, characterized by a combination of data and physics losses, enabling the DINN to follow the differential equations describing the disease progression.

During the early stages of the COVID-19 pandemic, it was important to minimize both medical and economic costs. In this study, we introduce a machine learning-based multi-objective optimization framework that can propose cost-effective social distancing strategies. Our approach finds Pareto solutions that balance different goals, like reducing infections and minimizing social distancing costs. Then, the cost-benefit analysis can adjust each cost factor—for example, value of statistical life (VSL), fatality rate, or GDP. We also provide an interactive web dashboard so that policymakers and the public can test various scenarios easily. We tested this framework on the COVID-19 pandemic in Korea. The results show that the difference between the cost-optimal strategy and implemented strategy is 10% in cost. Notably, our results reveal two distinct patterns in cost-optimal solutions. When social distancing cost is proportional to intervention intensity, an on-off lockdown strategy proves most economical. In contrast, when the cost increases sharply with intensifying social distancing, maintaining a consistently moderate level of intervention minimizes overall expenses. By letting people explore different cost settings and intervention strategies, this tool can support more balanced decisions during emerging infectious disease crises in the future.

The COVID-19 pandemic took its toll not only on the physical health of those who lived through it, but also on their mental health. I will present a model of the direct and indirect effects of COVID-19 and the associated public policies on mental health. This is an SIRS model of COVID-19, with compartments for mild, acute, and chronic COVID-19 infections and additional compartments for populations with mental health symptoms. Parameters are determined based on fitting from the first wave of COVID-19 in the New York state population, which includes several changes in local government policy, e.g, lockdown orders, which have an effect on the rate at which mental health systems develop. Finally, an additional “delay” term is included in the model to account for the delay between lockdowns going into effect and individuals developing mental health symptoms. The goal of our analysis is to understand how government policy in response to a pandemic can seek to maximize the population's quality-adjusted life years (QALY), which is a measure not only of lifespan, but also the quality of the years lived. I will present some preliminary results suggesting the optimal timing and strength of government lockdown mandates.

After an individual has been infected by a pathogen, T lymphocytes store information about the pathogen. Consequently, upon reinfection by the same pathogen, an immune response memory is triggered. This immune memory allows the body to react very quickly against the pathogen. Indeed, when an individual recovers from a virus, sometimes the individual acquires full immunity. In some cases, the individual's immunity persists for some period, after which it decreases progressively and can even disappear. I will then present a stochastic epidemic model with memory on the previous infection, incorporating varying infectivity and waning immunity. More precisely, we will present a functional law of large numbers when the size of the population tend to infinity. We will also present results on the behaviour of the epidemic, more precisely the threshold for the existence of an endemic equilibrium, and study the stability of the endemic equilibrium.

Symptom propagation occurs when an individual’s symptom severity is correlated with the symptom severity of the individual who infected them. Determining whether - and to what extent - these correlations exist requires data-driven methods. In this study, we use synthetic data to determine the types of data required to estimate the strength of symptom propagation and investigate the effect of reporting bias on these estimates. We found that even a relatively small number of contact tracing data points was sufficient to gain a reasonable estimate for the strength of symptom propagation. Increasing the number of contact tracing data points further improved our estimates. In contrast, population incidence alone was insufficient to accurately estimate the symptom propagation parameters, even with a large number of data points. Nonetheless, concurrently using population incidence data with contact tracing data led to increased accuracy when estimating the overall disease severity. We then considered the effect of severe cases being more likely to be reported in the contact tracing data. When contact tracing data alone was used, we found that our estimates for the strength of symptom propagation were robust to all reporting bias scenarios considered. However, the reporting bias led us to overestimate the overall disease severity. Using population incidence data in addition to contact tracing data reduced the error in disease severity but at the cost of increasing the error in the strength of symptom propagation when reporting bias was in both primary and secondary cases. Consequently, these errors led to us sometimes finding support for symptom propagation, even when the synthetic data was generated without.

We aim to compare results for the predicted reduction in vectorial capacity caused by pyrethroid and pyrethroid-piperonyl butoxide insecticide treated nets (ITNs) between semi-field Ifakara Ambiant Chamber tests (I-ACT) and experimental hut experiments. Mathematical modelling and Bayesian inference frameworks estimated ITN effects on mosquito behavioural endpoints (repelled, killed before/after feeding) to predict reductions in Anopheles gambiae’s vectorial capacity for Plasmodium falciparum transmission. The reduction in biting estimates are generally greater for I-ACT, possibly due to lower mosquito aggression: Although I-ACT vectors are probing before release, experimental hut vectors are actively seeking a blood meal. I-ACT estimates higher probability of killing vectors which have fed, while experimental huts show greater killing before feeding, possibly due to their open-system design, where vectors can contact the net, then attempt to exit and get trapped. This is supported by most of the mosquitoes being caught before feeding being in the exit trap. While the I-ACT is a closed system, were vectors cannot exit or be trapped, increasing the likelihood of returning to host-seeking and feeding. Despite these differences, both methods yielded similar predictions for the overall reduction in vectorial capacity. Results suggest that I-ACT provides a good initial assessment of the impact of adulticide modes of action of these nets. Challenges of semi-field experiments include how to model the change in efficacy from practical use over time. However, important advantages include the ability to easily trial different strains of vector (including different resistance levels) and allowing rapid data collection. Parameterising models with location-specific bionomic parameters allows for setting -specific predictions of the impact of different nets, with the potential to include additional modes of action for other active ingredients.

We aim to compare results for the predicted reduction in vectorial capacity caused by pyrethroid and pyrethroid-piperonyl butoxide insecticide treated nets (ITNs) between semi-field Ifakara Ambiant Chamber tests (I-ACT) and experimental hut experiments. Mathematical modelling and Bayesian inference frameworks estimated ITN effects on mosquito behavioural endpoints (repelled, killed before/after feeding) to predict reductions in Anopheles gambiae’s vectorial capacity for Plasmodium falciparum transmission. The reduction in biting estimates are generally greater for I-ACT, possibly due to lower mosquito aggression: Although I-ACT vectors are probing before release, experimental hut vectors are actively seeking a blood meal. I-ACT estimates higher probability of killing vectors which have fed, while experimental huts show greater killing before feeding, possibly due to their open-system design, where vectors can contact the net, then attempt to exit and get trapped. This is supported by most of the mosquitoes being caught before feeding being in the exit trap. While the I-ACT is a closed system, were vectors cannot exit or be trapped, increasing the likelihood of returning to host-seeking and feeding. Despite these differences, both methods yielded similar predictions for the overall reduction in vectorial capacity. Results suggest that I-ACT provides a good initial assessment of the impact of adulticide modes of action of these nets. Challenges of semi-field experiments include how to model the change in efficacy from practical use over time. However, important advantages include the ability to easily trial different strains of vector (including different resistance levels) and allowing rapid data collection. Parameterising models with location-specific bionomic parameters allows for setting -specific predictions of the impact of different nets, with the potential to include additional modes of action for other active ingredients.

Chemical Reaction Networks (CRNs) provide a powerful framework for modeling interactions between multiple chemical species within complex biological pathways. Their dynamics can be described by the mass-action law, representing variations over time in species concentrations through large systems of ODEs. The study of CRNs modeling signaling mechanisms within single cells has been successfully applied to provide insight into oncogenic pathways, enabling more predictive cancer models and improved therapeutic strategies [1,2,3]. Developing a CRN requires selecting the involved species, defining their interactions, and identifying key parameters such as microscopic reaction rates. This talk presents a method to retrieve these rates, particularly in the context of enzyme kinetics. We consider a simple model, where an enzyme E binds to a substrate S, forming a complex C that dissociates into a product P while regenerating E. This process is governed by four ODEs with three reaction rates: forward k_f, reverse k_r, and catalytic k_cat, typically hard to determine experimentally. A common simplification leads to the Michaelis-Menten (MM) model, where enzyme kinetics is characterized by a single ODE depending on two measurable parameters: the Michaelis constant K_M, and the maximum reaction rate V_max. These parameters may be expressed as functions of microscopic rates and are more accessible experimentally. This talk addresses the inverse problem of estimating k_f and k_r from K_M and V_max. A computational algorithm for solving this problem is presented and analyzed, along with an estimate of the reconstruction accuracy, and numerical simulations demonstrating its potential for refining kinetic models in biological and biomedical research. [1] Sommariva et al., J. Math. Biol., 82 (6): 55, 2021. [2] Berra et al., J. Optim. Theory Appl., 200 (1): 404-427, 2024. [3] Sommariva et al., Front. Syst. Biol., 3: 1207898, 2023.

Most biological systems have underlying multiple spatial or temporal scales that require reduced-order models to capture their essential dynamics and analyze them. However, traditional model reduction techniques, such as Computational Singular Perturbation (CSP), rely on the availability of the governing or dynamical equations, which are often unknown from data in biomedical applications. To address this limitation, we propose a data-driven CSP framework that integrates Sparse Identification of Nonlinear Dynamics (SINDy) and Neural Networks to extract time-scale separated models directly from data. Our approach is validated on the Michaelis-Menten enzyme kinetics model, a well-established multiscale system, by identifying reduced models for the standard Quasi-Steady-State Approximation (sQSSA) and reverse Quasi-Steady-State Approximation (rQSSA). When the full model cannot be identified by SINDy due to noise, we use Neural Networks to estimate the Jacobian matrix, allowing CSP to determine the regions where reduced models are valid. We further analyze Partial Equilibrium Approximation (PEA) case, where the dynamics span both sQSSA and rQSSA regimes, requiring dataset splitting to accurately identify region-specific models. The results demonstrate that SINDy struggles in the presence of noise to identify full model from data that have underlying timescale evolution, but remains effective for identifying reduced models when dataset are partitioned correctly.

Although chemical reaction networks (CRN) provide performant and biophysically plausible models for explaining single-cell genomic data, inference of reaction network parameters in this setting usually assumes available data points can be viewed as independent samples from a steady state distribution. Less is known about how to perform efficient parameter inference in the case that there is a continuous-time data stream, which adds complexity like nontrivial correlations between samples from different times. In the continuous-time setting, one has two natural questions: (i) given a set of reactions that could plausibly explain the observed data stream, what are reasonable estimates of the associated reaction rate parameters? and (ii) what is the minimal set of reactions necessary to explain the data? Both questions can be formalized as Bayesian inference problems, with the former concerning the inference of a model-dependent parameter posterior, and the latter concerning ‘structure’ inference. If one can assume each possible reaction has a different stoichiometry vector, there is a well-known analytic solution to both problems; if reactions can have the same stoichiometry vector (i.e., there is reaction degeneracy), both problems become substantially more difficult, and no analytic solution is known. We present the first approximate analytic solution to both problems, which is valid when the number of observations becomes sufficiently large. In its regime of validity, this solution allows its user to avoid expensive likelihood computations that can involve summing over an exponentially large number of terms. We discuss interesting consequences of this solution, like the fact that ‘simpler’ models with fewer reactions are preferred over more complex ones, and the fact that the parameter posteriors of non-identifiable models are strongly prior-dependent.

Real-world systems, especially in biology, exhibit significant complexity and inherent limitations in observability. What methods can enhance our understanding of the mechanisms underlying their functionality? Additionally, how can we develop and test explanatory models within a stochastic environment? Evaluating the effectiveness of these models requires quantitative measurements of the disparity between model outputs and observed data. While mean-field, deterministic models have well-established approaches for such assessments, stochastic systems—particularly those constrained by multiple data types—need carefully designed quantitative comparison methods. Methods for inferring the parameters of stochastic models generally require analytical forms of the model solutions, large data sets, summary statistics, or assumptions on the distribution of model outputs. These approaches can be limiting if you wish to preserve the information in the variability of the data but you do not have sufficient data to reliably fit distributions or determine robust statistics. We present a hierarchical approach to develop a distance measure for the direct comparison of model output distributions to experimentally observed distributions, avoiding any assumptions about distributions and the need to choose summary statistics. Our distance measure allows for constraining the model with multiple experiments, not necessarily of the same type, such that each experiment constrains some, or all, of the model parameters. We use this distance for parameter estimation with our queuing model of intracellular GLUT4 translocation. We will explore some practical considerations when using the distance for parameter inference, such as the effects of model output sampling and experimental error. Fitting the queuing model to data allowed us to uncover a possible mechanism of GLUT4 sequestration and release in response to insulin. Authors: Brock D. Sherlock and Adelle C.F. Coster

Collective behaviour in animal groups are coordinated movements and interactions among members that aim to achieve a common goal. Whether these goals are for allocation of resources or defence from predators, the collective behaviour appears to be largely a group activity initiated by a member, known as the leader. In the absence of high-resolution spatio-temporal data, various qualitative studies offer a glimpse of how leader-follower interactions take place. For example, Nagy, et al., studied the average delay in response when pigeons change the direction inflight. Next, Bourjade, et al., studied the first mover and the succeeding order of movements of Przewalski's Horses. Furthermore, various studies on the collective motion in the animal kingdom offer mathematical models and infer how the interactions and decision making take place. Important questions arise during an event of coordinated motion in animals. During such an event, do individuals move according to a certain set of natural rules? Or certain patterns form due to the influence of a leader? How is this influence measured? Finally, how is influence transferred to other members of the group? In this study, we discuss the role of information theory to quantitatively uncover leader-follower relationship in a horse group. Specifically, we introduce concepts from information theory, specifically global and local transfer entropy being applied to a harem of horses. We will discuss their definitions, and how these key concepts are used to support causation in events. We will then discuss some important implications on how this technique can be used to analyse collective motion where data is scarce.

Alzheimer’s Disease (AD) is a progressive neurodegenerative disorder affecting approximately 10% of Americans over age 65, leading to memory loss, cognitive decline, and impaired daily function. Disease progression correlates with the spread of tau and amyloid-beta proteins, which aggregate into neurofibrillary tangles. While macroscopic whole-brain network models predict large-scale protein deposition patterns, they lack the specificity to capture individual disease progression. Conversely, microscale neuron-neuron models offer highly detailed biochemical aggregation and transport simulations but are computationally prohibitive for whole-brain parameter inference. In this work, we explore using machine learning to accelerate whole-brain simulations by approximating explicit solutions to the microscopic Two-Neuron Transport Model (TNTM), a partial differential equation describing tau flux along a neuron, incorporating biochemical aggregation, fragmentation, and transport. We simulate a single-edge model across physiological ranges of biochemical parameters and boundary conditions, then compare regression methods with varying levels of interpretability, from neural networks (low) to symbolic regression via PySR (high). Neural networks achieve the lowest error but lack biological insight. Linear and polynomial regression compute rapidly but yield high errors with limited interpretability. Symbolic regression achieves a balance between accuracy and transparency. This work demonstrates the potential of machine learning for computationally scalable AD modeling, opening avenues for patient-specific parameterization using AD data repositories.

We introduce a method to identify phase equations for heterogeneous oscillators beyond the weak coupling regime. This strategy is an extension of the theory from [Y. Park and D. Wilson, SIAM J. Appl. Dyn. Syst., 20 (2021), pp. 1464--1484] and yields coupling functions for N general limit-cycle oscillators with arbitrary types of coupling, with similar benefits as the classic theory of weakly coupled oscillators. These coupling functions enable the study of oscillator networks in terms of phase-locked states, whose stability can be determined using straightforward linear stability arguments. We demonstrate the utility of this approach by reducing and analyzing conductance-based thalamic neuron model. The reduction correctly predicts the emergence of new phase-locked states as a function of coupling strength and heterogeneity. We conclude with a brief remark on recent extensions to n:m phase-locking and N-body interactions.

A biopsy is a common procedure used to diagnose diseases like cancer, infections, or inflammatory conditions. In cell population studies, biopsy samples provide valuable data to analyze cellular growth, proliferation rates, and structural abnormalities, which are essential for understanding disease progression. Estimating the growth (proliferation) rate of human cells is a challenging task. To address this, we have developed a method based on the birth-death Markov process to simulate the logistic growth model. We applied an extended Mean Field Approximation (MFA) for birth-death Markov processes, which accounts for fluctuations in the evolution of observables, such as moments. By calculating the theoretical moments from the birth-death process, we solved the inverse problem and estimated the growth rate. Additionally, we performed Markov Chain Monte Carlo (MCMC) simulations for both logistic growth and logistic growth with the Allee effect. The moments of the simulated population were used to predict the growth rate through regression analysis, achieving a high R-squared value. Finally, by applying this approach to biopsy data, one can estimate the proliferation rate of human cells with greater accuracy.

Introduction: Mathematical modeling is a key tool for understanding solid tumor biophysics, progression, and treatment resistance. Biophysical changes, such as elevated interstitial fluid pressure (IFP), are identified as major barriers to effective drug delivery. Incorporating patient-specific data into mathematical models offers the potential for personalized prognosis and treatment strategies for cancer patients. In this study, we explored the integration of patient-specific data from dynamic contrast-enhanced MRI (DCE-MRI) and diffusion-weighted MRI (DW-MRI) into a computational fluid dynamics (CFD) model of solid tumors to estimate the IFP and drug penetration profiles. Methods: As part of a translational study (EC/2019/1330, approved by Ghent University Hospital, Belgium), a patient with peritoneal metastasis underwent multi-sequential MRI, including T1-weighted (T1w) anatomical imaging, DCE-MRI, and DW-MRI. Tumor interstitial fluid pressure (IFP) was directly measured using a pressure transducer-tipped catheter for model validation. The CFD tumor model described interstitial fluid flow using Darcy’s law, the continuity equation, and Starling’s law, while drug penetration was modeled via the convection-diffusion-reaction equation. The 3D tumor geometry was derived from T1w images, vascular permeability from DCE-MRI using the extended-Tofts model, and hydraulic conductivity from DW-MRI. Results: An elevated IFP zone was observed in the central region of the tumor (up to 14 mmHg), while a lower IFP zone appeared at the tumor's edge. The clinically recorded IFP values (12.0 ± 2.5 mmHg) corresponded well with the simulation results. Drug penetration varied across the tumor, with deeper penetration in low-IFP regions. Conclusion: An image-based CFD model captured IFP and drug distribution variability, aligning with clinical data. This approach advances personalized oncology, potentially improving treatment strategies through noninvasive, patient-specific modeling.

The immune system can eradicate cancer, but various immunosuppressive mechanisms active within a tumor curb this beneficial response. Cytotoxic T cells (CD8s), regulatory T cells (Tregs), and antigen-presenting dendritic cells (DCs) play an important role in the immune response; however, it is very cumbersome to unravel the effects of multimodal interactions between tumor and immune cells and their contributions to tumor control using an experimental approach alone. Thus, to better understand the mechanisms that govern the interactions between immune cells and tumor cells and to identify the critical immunological features associated with tumor control and tumor escape, we built a mechanistic mathematical model of CD8s, Tregs, DCs, and tumor cells. The model accounts for tumor immunogenicity, the effects of IL-2 prolonging T cell lifespan, Treg suppression of antitumor immune response through CTLA-4, recruitment of immune cells into the tumor environment, and interferon-gamma upregulation of PD-L1 on DCs and tumor cells to deactivate T cells. We successfully fit the model to experimental data of tumor and immune cell dynamics. Employing Latin Hypercube Sampling, we generated over 1000 parameter sets that capture the sensitivity of αPD-1 immunotherapy. By comparing the parameter sets, we gain an insight into which mechanisms impinge the success of immunotherapy. We are now utilizing the model to explore combination immunotherapies that enhance the immune response in partial- and non-responders of αPD-1 immunotherapy. In particular, we are investigating what combinations of αPD-1, αCTLA-4, αICOSL, and αLAG-3 will stimulate CD8 activation without promoting Treg activation. After identifying the top combination therapies, we will validate our predictions experimentally. This integrated approach of modeling and experimental validation aims to advance our understanding of tumor-immune interactions and guide the development of more effective immunotherapeutic strategies.

Cancer is a significant global health issue, with treatment challenges arising from intratumor heterogeneity. This study examines the complex relationship between somatic evolution and phenotypic plasticity, explicitly focusing on the interplay between cell migration and proliferation [1]. We propose that evolution does not act directly on phenotypic traits, like the proliferation rate, but on the phenotypic plasticity in response to the microenvironment [2]. We study this hypothesis using a novel, spatially explicit model that tracks individual cells' phenotypic and genetic states. We assume cells change between mobile and growing states controlled by inherited and mutation-driven genotypes and the cells' microenvironment. We observe that cells at the tumor edge evolve to favor migration over proliferation and vice versa in the tumor bulk. However, this phenotypic heterogeneity can be realized by distinct regulations of the phenotypic switch, which depend on the apoptosis rate and the cells' ability to sense their environment. Emerging synthetic tumors display varying levels of heterogeneity, which we show are predictors of the cancer's recurrence time after treatment. Interestingly, higher phenotypic heterogeneity predicts poor treatment outcomes, unlike genetic heterogeneity. [1] Hatzikirou, H. et al. (2010). 'Go or Grow': the key to the emergence of invasion in tumour progression? Math. Med. Biol., 29(1), 49-65. [2] Syga S. et al. (2024) Evolution of phenotypic plasticity leads to tumor heterogeneity with implications for therapy. PLOS Comput. Biol. 20(8): e1012003.

A pivotal aspect of developing effective immunotherapies for solid tumors is the robust testing of product efficacy inside in vitro platforms.Collaborating with an experimental team that developed a novel microfluidic device at Children’s National Hospital (CNH), we developed a mathematical model to investigate immune cell migration and cytotoxicity within the device. Specifically, we study Chimeric Antigen Receptor (CAR) T-cell migration inside the channels, treating the cell as a moving boundary driven by a chemoattractant concentration gradient. The chemoattractant concentration is governed by two partial differential equations (PDEs) that incorporate key geometric elements of the device. We examine the motion of the cell as a function of its occlusion of the channel and find that certain cell shapes allow for multiple cells to travel inside the channel simultaneously. Additionally, we identify parameter regimes under which cells clog the channel, impairing their movement. All our findings are validated against experimental data provided by CNH. We integrate our model results into a broader model of the device, which also examines the cytotoxicity of CAR T-cells. This provides a tool for distinguishing experimental artifacts from genuine CAR T-cell behavior. This collaboration enabled the team at Children’s National Hospital to refine experimental conditions and uncover mechanisms enhancing CAR T-cell efficacy. [1] D Irimia, G Charras, N Agrawal, T Mitchison, M Toner, Polar stimulation and constrained cell migration in microfluidic channels,, Lab on a Chip 7 (12), 1783-1790

Extrachromosomal DNA (ecDNA) oncogene amplification is associated with treatment resistance and shorter survival in cancer. Currently, the spatial dynamics of ecDNA, and their evolutionary impact, are poorly understood. Here, we investigate ecDNA spatial-temporal evolution by integrating computational modeling with samples from 94 treatment-naive human IDH-wildtype glioblastoma patients. Random ecDNA segregation combined with ecDNA-conferred fitness advantages induce predictable spatial ecDNA copy-number patterns which depend on ecDNA oncogenic makeup. EGFR-ecDNAs often reach high copy-number, confer strong fitness advantages and do not co-amplify other oncogenes on the same ecDNA. In contrast, PDGFRA-ecDNAs reach lower copy-number, confer weaker fitness advantages and co-amplify other oncogenes. EGFR structural variants occur exclusively on ecDNA, arise from and are intermixed with wild-type EGFR-ecDNAs. Modeling suggests wild-type and variant EGFR-ecDNAs often accumulate before clonal expansion, even in patients co-amplifying multiple ecDNA species. Early emergence of oncogenic ecDNA under strong positive selection is confirmed in vivo and in vitro in mouse neural stem cells. Our results implicate ecDNA as a driver of gliomagenesis, and suggest a potential time window in which early ecDNA detection may facilitate more effective intervention.

Monoclonal antibody (mAb)-based therapeutics are pivotal in treating a wide range of diseases, including cancer. One key mechanism by which these antibodies exert anti-tumour effects is through antibody-dependent cellular cytotoxicity (ADCC). In ADCC, mAbs bind to specific antigens on tumour cells and Fc receptors on immune effector cells. This trimeric complex triggers these effector cells to kill the tumour. ADCC is influenced by multiple factors, notably the properties of the mAb and its interactions with Fc receptors and target antigens. However, the optimum conditions for ADCC remain unclear. In this study, we investigate how variations in target antigen and mAb properties, particularly antibody valency, modulate ADCC response to identify parameters that maximize its potency. We developed an ordinary differential equation (ODE) model to simulate mAb binding within the immune synapse and quantify trimeric complex formation. To link the number of trimeric complexes to ADCC response, we validated the model using Bayesian inference on ADCC assay data. The results suggest that lower-affinity mAbs enhance ADCC by increasing the number of target cell-bound antibodies. Our validated model indicates that a “steric penalty” is necessary for bivalently target-bound versus monovalently target-bound antibodies. Due to constraints from dual antigen binding, these antibodies experience limited mobility, reducing Fc receptor engagement. After model validation, we explored variations in target expression, binding affinity, and antibody valency on ADCC potency, quantified by EC50. Our key finding is that the optimal binding affinity for maximizing ADCC potency depends on antibody valency. Monovalent antibodies are most potent at high affinity, while bivalent antibodies peak at lower affinities. Furthermore, the magnitude of this effect varies with target expression levels.

Extrachromosomal DNA (ecDNA) oncogene amplification is associated with treatment resistance and shorter survival in cancer. Currently, the spatial dynamics of ecDNA, and their evolutionary impact, are poorly understood. Here, we investigate ecDNA spatial-temporal evolution by integrating computational modeling with samples from 94 treatment-naive human IDH-wildtype glioblastoma patients. Random ecDNA segregation combined with ecDNA-conferred fitness advantages induce predictable spatial ecDNA copy-number patterns which depend on ecDNA oncogenic makeup. EGFR-ecDNAs often reach high copy-number, confer strong fitness advantages and do not co-amplify other oncogenes on the same ecDNA. In contrast, PDGFRA-ecDNAs reach lower copy-number, confer weaker fitness advantages and co-amplify other oncogenes. EGFR structural variants occur exclusively on ecDNA, arise from and are intermixed with wild-type EGFR-ecDNAs. Modeling suggests wild-type and variant EGFR-ecDNAs often accumulate before clonal expansion, even in patients co-amplifying multiple ecDNA species. Early emergence of oncogenic ecDNA under strong positive selection is confirmed in vivo and in vitro in mouse neural stem cells. Our results implicate ecDNA as a driver of gliomagenesis, and suggest a potential time window in which early ecDNA detection may facilitate more effective intervention.

Monoclonal antibody (mAb)-based therapeutics are pivotal in treating a wide range of diseases, including cancer. One key mechanism by which these antibodies exert anti-tumour effects is through antibody-dependent cellular cytotoxicity (ADCC). In ADCC, mAbs bind to specific antigens on tumour cells and Fc receptors on immune effector cells. This trimeric complex triggers these effector cells to kill the tumour. ADCC is influenced by multiple factors, notably the properties of the mAb and its interactions with Fc receptors and target antigens. However, the optimum conditions for ADCC remain unclear. In this study, we investigate how variations in target antigen and mAb properties, particularly antibody valency, modulate ADCC response to identify parameters that maximize its potency. We developed an ordinary differential equation (ODE) model to simulate mAb binding within the immune synapse and quantify trimeric complex formation. To link the number of trimeric complexes to ADCC response, we validated the model using Bayesian inference on ADCC assay data. The results suggest that lower-affinity mAbs enhance ADCC by increasing the number of target cell-bound antibodies. Our validated model indicates that a “steric penalty” is necessary for bivalently target-bound versus monovalently target-bound antibodies. Due to constraints from dual antigen binding, these antibodies experience limited mobility, reducing Fc receptor engagement. After model validation, we explored variations in target expression, binding affinity, and antibody valency on ADCC potency, quantified by EC50. Our key finding is that the optimal binding affinity for maximizing ADCC potency depends on antibody valency. Monovalent antibodies are most potent at high affinity, while bivalent antibodies peak at lower affinities. Furthermore, the magnitude of this effect varies with target expression levels.

Diabetic kidney disease (DKD) progression is often marked by early glomerular endothelial cell (GEC) dysfunction, including alterations in the fenestration size and number linked to impaired glomerular filtration. However, the cellular mechanisms regulating GEC fenestrations remain poorly understood due to limitations in existing in vitro models, challenges in imaging small fenestrations in vivo, and inconsistencies between in vitro and in vivo findings. This study used a logic-based protein-protein interaction network model with normalized Hill functions for dynamics to explore how glucose-mediated signaling dysregulation impacts fenestration dynamics in GECs. We identified key drivers of fenestration loss and size changes by incorporating signaling pathways related to actin remodeling, myosin light chain kinase, Rho-associated kinase, calcium, and VEGF and its receptor. The model predicted how hyperglycemia in diabetic mice leads to significant fenestration loss and increased size of fenestrations. We found that glycemic control in the pre-DKD stage mitigated signaling dysregulation but was less effective as DKD developed and progressed. The model suggested alternative disease intervention strategies to maintain fenestration structure integrity, such as targeting Rho-associated kinase, VEGF-A, NFκB, and actin stress fibers.

Mutation accumulation (MA) experiments are crucial for understanding evolution. In microbial populations, these experiments typically involve periods of population growth, where a single individual forms a visible colony, followed by severe bottlenecks. Studies on the effects of positive and negative selection in MA experiments have shown that, for example, with 20 generations of growth between bottlenecks, beneficial mutations will be substantially over-represented cite{}; this effect is known as ``selection bias''. In previous work, we developed a fully stochastic discrete-time model that includes realistic offspring distributions, accounting for genetic drift and allowing for the loss of rare lineages. We demonstrated that when drift is included, selection bias is even stronger than previously predicted cite{}. Since bacterial division is unlikely to remain synchronized over 20 or more generations, this study extends the discrete-time model to a continuous-time framework. Since lineages that start reproducing early accrue a compounded advantage, a continuous-time model offers an even more accurate correction for selection in MA experiments.