Registration and Welcome
The Early Careers Workshop is an opportunity for junior members (likely students and postdocs) to meet up in advance of the Annual meeting and learn about issues facing mathematical biologists in terms of getting a job, learning about different career options, communications and networking with peers and also more senior society members.
It will open with a primer on "How to survive SMB: making the most of your conference money".
In this talk I will describe the application of state of the art image based modelling to several seemingly different areas of biology. I show examples from biomedical (lymphatic, vascular and lung system) and agricultural problems of plant soil interaction. I will describe the workflow from imaging (X-ray CT, XRF, SEM-EDX, histology), image reconstruction, image segmentation, computation and how to utilize this work stream to synthesise new scientific knowledge. In particular I will also outline several challenges and bottle necks in this process to hopefully encourage more mathematicians to get involved in the full pipeline.
This minisymposium will bring together established and up-and-coming researchers to share ideas on the recent advances on the use of mathematical approaches to gain insight into the transmission dynamics and control of emerging and re-emerging infectious diseases in populations.
An essential feature of the minisymposium is the emphasis on the use of state-of-the-art techniques, theories and new applications associated with the use of dynamical systems in modeling the spread of diseases in populations. Current pressing modelling and mathematical challenges will also be discussed.
Vector-borne diseases cause worldwide concern with hundreds of millions of new cases and over a million deaths reported annually. Mathematical models are a key tool in the study of the spread of diseases such as malaria and dengue. In particular, transmission models have been successful at determining the most promising intervention strategies, despite the fact that many of these models assume all individuals experience identical infections. Within-host dynamics of infections, however, vary widely among individuals. In the case of dengue, within-host dynamics differ between primary and secondary infections, where secondary infections with a different virus serotype typically last longer, produce higher viral loads, and induce more severe disease. Here, we build upon models of variable within-host dengue virus dynamics resulting in mild dengue fever and severe dengue hemorrhagic fever by coupling them to a population-level model. The resulting multiscale model examines the dynamics of between-host infections in the presence of two circulating virus strains that involves feedback from the within-host and between-hosts scales. We analytically determine a threshold under which infections persist in the population.
Climate change is known to significantly affect the dynamics of vector-borne diseases, such as malaria. In particular, the species involved in the transmission dynamics of malaria are affected by various abiotic conditions, such as temperature, precipitation, humidity and vapor pressure . A number of models, typically statistical (using data and statistical approaches to correlate some climate variables with malaria incidence) or mechanistic (accounting for the detailed dynamic nonlinear processes involved in disease transmission), have been employed to assess the likely impact of anthropogenic climate change on malaria transmission dynamics and control. These models have (generally) reached divergent conclusions, with some predicting a large expansion in the continental land area suitable for transmission and in the number of people at risk of malaria, while others predict only modest poleward (and altitudinal) shifts in the burden of disease, with little net effect, and the issue remains unresolved thus far. I will discuss some of our recent results on modelling the effect of climate variables (notably temperature) on the dynamics of malaria vector and disease.
Many human papillomavirus (HPV) vaccination programs currently administer three vaccines - a bivalent, a quadrivalent and a nonavalent vaccine - for a two-dose course. In this talk a model will be presented to explore optimal vaccination strategies using the three vaccines, which differ in protection breadth, cross-protection, and type-specific efficacy. Assuming the HPV infection prevalence in the population under the constant vaccination regime, optimal control theory will be used to discuss optimal vaccination strategies for the associated non-autonomous model when the vaccination rates are functions of time. The impact of these strategies on the number of infected individuals and the accumulated cost will be assessed and compared with the constant control case. Switch times from one vaccine combination to a different combination including the nonavalent vaccine will be assessed during an optimally designed HPV immunization program.
Syphilis, a major sexually-transmitted disease, continues to pose major public health burden in both under-developed and developed nations of the world. This study presents a new two-group sex-structured model for assessing the community-level impact of treatment and condom use on the transmission dynamics and control of syphilis. Rigorous analysis of the model shows that it undergoes the phenomenon of backward bifurcation. In the absence of this phenomenon (which is shown to arise due to the re-infection of recovered individuals), the disease-free equilibrium of the model is shown to be globally-asymptotically stable (GAS) when the associated reproduction number is less than unity. Furthermore, the model can have multiple endemic equilibria when the reproduction threshold exceeds unity. Numerical simulations of the model, using data relevant to the transmission dynamics of the disease in Nigeria, show that, with the assumed 80% condom efficacy, the disease will continue to persist (i.e., remain endemic) in the population regardless of the level of compliance in condom usage by males. Furthermore, detailed optimal control analysis (using Pontraygin's Maximum Principle) reveals that for situations where the cost of implementing the controls (treatment and condom-use) considered in this study is low, channeling resources to a treatment-only strategy is more effective than channeling them to a condom-use only strategy. Furthermore, as expected, the combined condom-treatment strategy provides a higher population-level impact than the treatment-only strategy or the condom-use only strategy. When the cost of implementing the controls is high, the three strategies are essentially equally as ineffective.
Akt/PKB (Protein Kinase B) is a major crosstalk node in the mammalian cell. Located as the juncture of several key signalling pathways, it is involved in many cellular processes, such as glucose metabolism, cell growth and the suppression of apotosis. Dysregulated Akt signalling is implicated in a range of human disorders, from diabetes to cancer.
Initially, Akt is synthesised in the unactivated state on the endoplasmic reticulum. However, the activation (phosphorylation) of Akt in response to insulin stimulation only occurs at the plasma membrane, necessitating the translocation of Akt from the interior to the periphery of the cell. At the moment, the understanding of this translocation process is still in its infancy, but there are some indications that it is a staged process.
We have developed a deterministic, three-compartment, ordinary differential equation (ODE) model of Akt translocation. Given a conservation relation implicit in the model, it can be shown that this system is equivalent to the damped harmonic oscillator equation; a classic, well-studied ODE. With this framework, we investigate the different modes of downstream signalling that can be produced by the model and the conditions for their manifestation.
A biochemical system is a biological system consisting of a collection of chemical compounds interacting with each other. One way to model and analyze a biochemical system is by using S-systems, which are coupled ordinary differential equations based on power-law formalism. In this work, we do a parameter estimation on an S-system model called HS96. This model, proposed by Hlavacek and Savageau in 1996, describes a simple genetic network consisting of five dependent variables $X_i,$ $i=1,2,...,5$. As a preliminary method, sensitivity analysis of HS96 is conducted to investigate the change in model outputs with respect to the changes in model parameters. Usual model outputs are $X_i$ and $\dot{X}_i,$ $i=1,2,...,5$. From the results of sensitivity analysis, model outputs are selected which are then used to estimate the parameters of the HS96 model. A hybrid of genetic algorithm with the interior point method is used in the parameter estimation.
Many biological systems rely on both long range and short range signals in order to produce proper cell patterns during development. These patterns, such as alternating cell fates, are defined by different gene and protein expression levels. However, the relative contribution of these two signalling modes in establishing proper patterns is not well understood. Using vulval development of two nematode worm species, C. elegans and C. briggsae, as motivation, we derive an asymptotic PDE based on a simplified signalling network consisting of EGF, Notch and Wnt. We demonstrate that additional long range signals can provide a framework for understanding species-specific differences and that variations in short range signalling can amplify long range signals. These results suggest that long and short range signals have critical roles to play in proper cell patterning in development. This work is joint with Helen Chamberlin and Carly Williamson.
Cardiovascular disease is the leading cause of death in Australia; responsible for 30% of deaths. Diseases that affect the heart are commonly accompanied by a condition known as hypertrophy — heart enlargement through cell growth. This condition leads to uneven heart beats and, eventually, heart failure.
Hypertrophy is regulated by the flow of calcium ions within the cell. Calcium dynamics control the activation of the proteins required for hypertrophic gene expression and their localisation to the cell nucleus. However, there has been little research done on calcium flow within the nuclei of normal heart cells let alone those of diseased hearts.
We focus specifically on modelling the mechanics of calcium flow within the nucleus and distinguishing the signal for heart cell growth in the nucleus where gene expression takes place over the background signal for the heartbeat. Using a finite element model, we describe the spatiotemporal properties of these background calcium oscillations as they propagate into the nucleus. We verify this model with existing data to produce a clearer picture of the calcium dynamics involved in initiating hypertrophy.
Robustness, and the ability to function and thrive amid changing and unfavourable environments, is a fundamental requirement for all living systems. Moreover, it has been a long-standing mystery how the extraordinarily complex communication networks inside living cells, comprising thousands of different interacting molecules, are able to exhibit such remarkable robustness since complexity is generally associated with fragility. In this talk I will give an overview of our recent research on robustness in cellular signalling networks, with an emphasis on the robust functionality known as Robust Perfect Adaptation (RPA). Our work has now suggested a resolution to the complexity-robustness paradox through the discovery that robust adaptive signalling networks must be decomposable into topological basis modules of just two possible types. This newly-discovered modularisation of complex bionetworks has important implications for evolutionary biology, embryology and development, cancer research and drug development.
From the foraging strategies of large organisms, to T-cells hunting pathogens, to proteins examining strands of DNA, carefully optimised search processes are a phenomenon that pervades throughout nature at many scales. Often these search processes do not proceed in isolation, but instead many instances proceed in parallel, competing for space and resources. In this talk I shall discuss optimisation of search processes on networked topologies where the searchers interact with each other through competition for space. Taking the simple exclusion process as the fundamental model for spatial interactions, I consider search strategies that seek to minimise the average cover time for individuals that search in parallel. We will see that the optimal strategy is to implement parallel searching at an optimal density of searchers that depends heavily on the given network topology, and that the optimal density can be well predicted by the spectral gap of the network. These results are verified over a broad class of networks, as well as real-world transport networks such as the London tube network. I will conclude by considering an asymmetric search process that reveals unexpected changes in efficiency between classes of networks.
In recent years, there has been an explosion of multiscale cancer patient data including genomic, metabolomic, proteomic, and histologic (with novel stainings), at the tissue scale and smaller, along with a host of advanced MRI and other clinical imaging sequences at the whole organ scale. While bountiful, the multiscale and spatial sparsity represented with these data types for any individual patient leaves the final interpretation difficult with traditional statistical methods. Mechanistic models provide a framework for exploring how these processes interlink and contribute to tumour growth and evolution. Incorporating this multitude of data into such models remains a challenging and open problem, however, promising results are emerging. The talks in this minisymposium will highlight novel models and approaches to using these new data types in exploring the mechanisms underlying cancer growth and response to therapy.
It is well known that growth and survival of cancer stem cells (CSCs) is highly influenced by tumour microenvironmental factors and molecular signalling, initiated by cytokines and growth factors. IL-6 is a key regulator of a number of cellular processes including proliferation, survival, differentiation, migration and invasion and it is also commonly overexpressed in many cancers. Recent evidence shows that IL-6 is not only secreted by tumour cells, but is produced at even higher levels by endothelial cells (ECs). Research shows that high intratumoural levels of IL-6 enhance the survival, self-renewal and tumour initiation potential of cancer stem cells. These studies of the impact of IL-6 on CSCs provide strong motivation for the development of anti-IL-6 therapies for the targeted treatment of stem cell driven cancers.
In this talk, a multi-scale mathematical model that operates at the intracellular, molecular, and tissue level is developed in order to investigate the impacts of endothelial cell-secreted IL-6 signalling on the crosstalk between tumour cells and ECs during tumour growth. This endothelial cell-tumour cell (EC-TC) model is used to study the therapeutic impact of Tocilizumab (TCZ), a competitive IL-6R inhibitor, on tumour growth and cancer stem cell (CSC) fraction, alone and in combination with the traditional chemotherapeutic agent, Cisplatin. The approach used here is novel in that it includes full receptor occupancy dynamics between endothelial-cell produced IL-6, IL-6R, and TCZ. Validation is achieved by directly comparing model predictions to data generated by a series of in vivo experiments. This multiscale approach provided excellent predictive agreement with the decrease in tumour volumes, as well as a decrease in CSC fraction.
Targeting key regulators of the cancer stem cell phenotype to overcome their critical influence on tumour growth is a promising new strategy for cancer treatment. This predictive modelling framework can serve to rapidly evaluate dosing strategies for IL-6 pathway modulation, as well as providing the basis for proposing combination treatments with IL-6 blockade and cytotoxic or other targeted therapies.
Mathematical modelling of the stochastic evolutionary process of carcinogenesis can be used to derive and to optimize the timing of clinical screens so that the probability is maximal that an individual is screened within a certain "window of opportunity" for intervention when early cancer development may be observed. By using data from epidemiological studies with long-term patient follow-up, empirical approaches aid in screening design and may inform cost-effectiveness analyses to compare proposed screening and intervention strategies. However, mechanistic modelling that incorporates a greater level of biological understanding and detail for how and when normal tissues progress to cancer can be used for a more refined screening design than typically implemented in population screening studies. With examples in inflammation-driven premalignant disease of the gastrointestinal tract, I will introduce calibrated multistage clonal expansion models of carcinogenesis, derive probability equations used for optimal screening times and risk estimation, and present results for screening strategies using inferred parameters from US cancer incidence data. These results 1) provide a robust statistical framework for quantifying when it is optimal to screen and begin surveillance for premalignant changes, and 2) are shown to be cost-effective in Barrett’s esophagus patients, particularly with increasingly sensitive and minimally invasive sampling procedures used in clinical practice.
Glioblastoma Multiforme (GBM) is the most aggressive primary brain tumour, with a median life expectancy of only 15 months with treatment. Although surgical resection is a standard-of-care procedure, the migratory nature of the tumour cells limits its efficacy as not all tumour cells can be removed. The tumour will usually re-establish itself along the resection cavity wall, known as a local recurrence, but it can recur elsewhere or become more migratory, known as a distal recurrence. Clinical data has shown that patients with ischemia following resection are more likely to have a distally recurring tumour, with an incidence rate of 61% for those with ischemia vs 19% for those without. This evidence suggests that a lack of nutrients plays a role in the tumour recurrence pattern but does not fully answer the question regarding the difference between tumours that do or do not recur distally under these conditions.
We present the Proliferation Invasion Hypoxia Necrosis Angiogenesis (PIHNA) model that simulates the angiogenic cascade of a GBM as it grows. We have applied this mechanistic model to show how an in silico GBM grows following resection and subsequent ischemia, and to highlight the role of the individual tumour kinetics in this behaviour. Our simulations suggest that the diffusivity of the tumour, the tumour cell proliferation rate, and the hardiness of the cells all play key roles in the recurrence location of a GBM.
Small-scale fluid mechanics play an important role in virtually all living systems. Biological processes taking place at very low Reynolds numbers pertain to individual cells (swimming of bacteria or sperm cells) as well as to cell collectives and even groups of larger, multicellular organisms. Classical theoretical explorations and efforts to model fluid flow at low Reynolds number have long been limited to single cells swimming through a Newtonian fluid. Biological problems, however, often involve complex environments and/or systems containing more than one cell or organism, and traditional models have yielded only approximate answers to these problems. Recently, the development of computational tools has spurred advancement in the study of low Reynolds number fluid mechanics in biological systems.
In this minisymposium, we focus on both locomotion in complex materials (such as soil or viscoelastic fluids) and the emergent collective dynamics observed in organisms at different scales by expanding on existing fluid models, developing better-performing algorithms, and integrating behavioural and fluid mechanical models.
Tissue development requires cells of different types to organise themselves into the appropriate patterns and structures to produce viable, functional tissue. Similar processes occur in tissue repair (e.g. wound healing) or when tissues are grown in vitro (tissue engineering). Understanding how this organisation is coordinated is therefore an important basic problem in biology and medicine.
I will present results from agent-based modelling of interacting cell populations, and illustrate how different interactions between the cells affect the patterns of cell organisation observed in tissues. I will explain how these patterns can be quantified using pair-correlation functions, and discuss the extent to which we can infer cell interactions from observed tissue-scale patterns.
Many important biological functions depend on microorganisms' ability to move in viscoelastic fluids such as mucus and wet soil. The effects of fluid elasticity on motility remain poorly understood partly because the swimmer strokes depend on the properties of the fluid medium, which obfuscates the mechanisms responsible for observed behavioural changes. We use experimental data on the gaits of Chlamydomonas reinhardtii swimming in Newtonian and viscoelastic fluids as inputs to numerical simulations that decouple the swimmer gait and fluid type in order to isolate the effect of fluid elasticity on swimming. In viscoelastic fluids, cells employing the Newtonian gait swim faster but generate larger stresses and use more power, and as a result the viscoelastic gait is more efficient. Furthermore, we show that fundamental principles of swimming based on viscous fluid theory miss important flow dynamics: fluid elasticity provides an elastic memory effect which increases both the forward and backward speeds, and (unlike purely viscous fluids) larger fluid stress accumulates around flagella moving tangent to the swimming direction, compared to the normal direction.
An important aspect in the study of reproduction is how sperm are guided toward an egg for fertilization. One such mechanism is the process of chemotaxis, in which the sperm detect changes in concentration (namely of Ca^{+}) in the fluid environment and utilize these changes to alter the waveform of their flagellar beat. This change in beat form results in changes to the swimming path. We model the swimming sperm using a Kirchhoff elastic rod model coupled with the method of regularized Stokeslets for the fluid motion. In order to simulate the effect of chemical concentrations, we employ a stochastic decision making process.
The development of an embryo, from a single cell to the breathtaking range of organisms found in nature, is a fascinating process. Moreover, the study of embryonic development promises to shed light on a wide range of developmental defects, inform strategies for the treatment of disease, and provide guidance for tissue regeneration and engineering. Developmental biology has witnessed a renaissance in recent years, due largely to the advent of significant new imaging techniques and molecular biology technologies. To exploit the volumes of multiplex quantitative data that are now being collected, we need to design quantitative hypotheses through mathematical models, and use these models to make quantitative predictions that can be tested and validated experimentally.
The aim of this Developmental Biology Subgroup activity is to bring together researchers in the field who will showcase exemplar interdisciplinary, quantitative studies in developmental biology. This will raise the profile of quantitative approaches to development within the SMB, encourage new members to join the Subgroup and work in the area, and provide a platform for new interactions between SMB members working in the area.
Retinal vasculature is essential for adequate oxygen supply to the inner layers of the retina, the light sensitive tissue in the eye. In embryonic development, formation of the retinal vasculature via angiogenesis is critically dependent on prior establishment of a mesh of astrocytes, which are a type of brain glial cell. Astrocytes emerge from the optic nerve head and then migrate over the retinal surface as a proliferating cell population in a radially symmetric manner. Astrocytes begin as stem cells, termed astrocyte precursor cells (APCs), then transition to immature perinatal astrocytes (IPAs) which eventually transition to mature astrocytes. We develop a partial differential equation model describing the migration of astrocytes where APCs and IPAs are represented as two subpopulations. Numerical simulations are compared to experimental data to assist in elucidating the mechanisms responsible for the distribution of astrocytes.
Understanding how 3D organ morphology is achieved during development is one of the ultimate goals in biology. This is important not only for pure scientific interests but also for potential medical applications for controlling and designing functional organs. To achieve these goals, it is essential to clarify the quantitative relationships between microscopic molecular/cellular activities and organ-level tissue deformation dynamics. While the former has been studied for several decades, the latter - macroscopic geometrical information about physical tissue deformation - has been lacking. One reason for this lack of information is the difficulty in measurement at high resolution. We recently proposed a Bayesian method to precisely reconstruct global deformation patterns for three-dimensional morphogenesis of curved epithelial sheets using positional data from sparsely-labeled cells with limited resolution [1]. The applications of the method to early development of chick forebrain and heart revealed that globally aligned anisotropic deformation (i.e., biased tissue stretching), rather than local area growth, is the predominant morphogenetic mechanism in both cases (specifically, for optic vesicle formation and C-looping, respectively) [2]. Comparing the reconstructed tissue deformation patterns and cellular behaviours, we quantitatively revealed the contributions of each cellular process (e.g., division, size/shape change, and rearrangement) to the anisotropic tissue deformation. In this talk, we would also like to introduce our attempt to build continuum mechanical models to reproduce observed morphologies and deformation patterns.
[1] Morishita et al., Nat. Commun., 2017
[2] Kawahira et al., under review
Insects use two main modes of segment determination during development: the ancestral short-germband mode (eg. Gryllus bimaculatus) where new segments are added sequentially, and the more derived long-germband mode (eg. Drosophila melanogaster), where all segments are determined simultaneously. In dipteran insects (flies, midges and mosquitoes), which use the long-germband mode of segmentation, the gap genes are activated by maternal gradients and cross-regulate each other to form the first zygotic layer of regulation in the segmentation gene hierarchy. We reverse-engineered a dynamical model of the gap genes in D. melanogaster from quantitative spatio-temporal expression data and used it to characterise the dynamics of gap gene pattern formation along the embryo trunk. We used tools and concepts from dynamical systems theory, coupled with a newly developed methodology designed specifically to address the effect of maternal gradient dynamics in the patterning process. This approach showed that two distinct dynamical regimes govern anterior and posterior trunk patterning. Stationary domain boundaries in the anterior rely on multi-stability. In contrast, the observed anterior shifts of posterior gap gene domains can be explained as an emergent property of an underlying regulatory mechanism implementing a damped oscillator. We have identified a dual-function three-gene motif embedded in the gap gene regulatory network, which is sufficient to recover both anterior and posterior dynamical regimes. Which one drives gene expression in a given region depends on the gap genes involved. Interestingly, this sub-network - known as the AC/DC motif - can also sustain oscillations. Oscillations are not found in the gap gene system but are characteristic of short-germband segmentation, suggesting that both modes share much more than previously thought. Studying the evolution of the gap gene network as an evolving dynamical system will tell us how oscillations could have arisen or ceased, and this will help us understand how long-germband segmentation might have repeatedly and independently evolved from the ancestral short-germband mode.
Phosphatidylinositol (3,4,5)-trisphosphate (PtdInsP3) is known to propagate as waves on the plasma membrane and is related to the membrane protrusive activities in Dictyostelium and mammalian cells. While there have been a few attempts to study their three dimensional dynamics, most of these studies focused on the dynamics extracted in one dimensional sections along the membrane in a single focal plane. However, the relation between the dynamics and three-dimensional cell shape remain elusive due to the lack of signalling information on the remaining non-observed part of the membrane. Here, we show that PtdInsP3 wave dynamics are directly regulated by the three-dimensional geometry - size and shape - of the plasma membrane. By introducing an analysis method that extracts the three-dimensional spatiotemporal activities on the entire cell membrane, we show that PtdInsP3 waves self-regulate their dynamics within the confined membrane area leading to changes in speed, orientation and pattern evolution following the underlying excitability of the signal transduction system. Our findings emphasize the role of the membrane topology in reaction-diffusion driven biological systems and indicate the importance of the plasma membrane topology in other mammalian systems.
In the last two decades, we have witnessed a rapidly increasing number of studies applying methods of mathematical modeling in immunology in general and specifically in the domain of T cell immunology. In part, this was due to the development of several quantitative and powerful techniques for detecting the dynamics of immune responses to pathogens ex-vivo, documenting cell dynamics directly in lymphoid tissues in-vivo, new methods to sequence the full T cell repertoire, and to predict T cell epitopes.
In parallel to the development of novel mathematical and computational models, this has led to better understanding of the results of these experiments, and providing a more quantitative view of the adaptive immune response. The purely theoretical models of the last decades have often been replaced by highly quantitative models that are checked on experimental in-vivo and in-vitro data in human and mice.
Here, we propose a mini-symposium to discuss these new models and their relation to novel experimental tools and the data from the perspective of both modellers and experimenters. Such minisymposium will be of interest beyond the domain of immunology.
Information analysis of amplicon sequencing of fecal samples from HIV positive individuals suggested an enrichment of a particular bacterial species in the gut microbiota, as referred to dysbiosis. Although a sufficient number of samples has already accumulated for both HIV positive and negative subjects, time-series datasets are rarely available. Hence difficulty remains in tracing the compositional change of the gut microbiota during HIV infection.
In this presentation, we would like to introduce our ongoing research progress on the application of trajectory inference methods to construct a pseudo-trajectory that may imitate disease progression of HIV infection. Our computational analysis suggested that enrichment of the enterobacteriaceae family may occur at the initial phase of HIV infection during which inflammatory responses would be facilitated in the gut. This finding is consistent with known observations for the enrichment of the enterobacteriaceae family in several inflammatory diseases.
Age plays an important role on immunity across the lifespan, as both very young and very old individuals are at higher risk of severe infection. CD8 T cells are important for controlling a number of viral and bacterial infections, and both the number and phenotype of CD8 T cells change with age. Various mathematical and experimental methods have been used to analyse T cell kinetics, and most of them assume a homogeneous population of cells. This assumption has been challenged experimentally and theoretically, that cells produced early in life may have different behaviour to those produced recently. We have recently developed a novel experimental model to ‘timestamp’ CD8 T cells in mice. More specifically, double positive CD8 T cells can be labeled in the thymus (using CD4 promoter tamoxifen-inducible cre to drive the expression of RFP permanently), thus the survival of these cells can be tracked over time.
Cells produced early in life show a faster decay, which then slows with time-since-production of the cell. Cells produced at later ages show a more stable behaviour, and their rate of decay also slows with time-since production of the cell. Comparing various mathematical models (based on the AIC values), we found that the population of T cell can be described by a model with three parameters. The model has an initial decay rate of cells made at birth, how the decay rate changes with the age of the host, and how the decay changes with the age of cell (since cell production).
Using the best-fit parameters from the model, we can develop another model about the ‘layering’ of cells made at different ages, and how this contributes to the total T cell pool. Based on our model, in a 300-day-old mouse, only 20% of the cells were made recently. We also modelled the process of cellular differentiation, and how this changes with age. The model can be described by a system of differential equations with non-constant parameters. We found that cells made early in life show a higher rate of differentiation.
Overall, this reveals an age-dependent heterogeneity in in vivo survival of CD8 T cells.
The prevention of diseases such as malaria require certain T cells to find all pathogens in the liver within a certain period of time (e.g., within 48 hours, which is the time required for liver-stage development of parasites in rodent malaria). This motivates the fundamental question of how many T cells are required to ensure complete coverage of the liver within a specified time, to a high probability. After describing our existing results, we will present our current thoughts and approaches. Explanations will be given of how to model the liver sinusoids, how to implement fast simulations of T cell movement, and various possibilities for how to model the movement of T cells.
CD8+ T cells can kill Plasmodium parasites in the liver of the mammalian host; a protective effect that can be harnessed for malaria vaccination. We have previously used intra-vital imaging to measure the interaction of CD8+ T cells in the liver and Plasmodium infected hepatocytes. We have previously observed that CD8+ T cells in the liver undertake LFA-1 dependent crawling motility in the allowing them to scan the liver and find infected hepatocytes. Subsequently we have observed that clusters of CD8+ T cells form around the infected hepatocyte, which may potentiate parasite killing. Key questions are how do these clusters form? And how are the parasites killed? Mathematical modelling suggested that clusters could form by at least two different processes. In the first process clusters form as a result of a positive feedback loop in which cells in the cluster recruit further effectors. Another model posits that clusters may form because it becomes harder for cells to leave the vicinity of an infected hepatocyte as clusters increase in size. To distinguish these possibilities we have imaged the behaviour of T cells in the liver to determine if they show evidence of directed migration towards infected hepatocytes. Moreover using Pertussis toxin mediated inhibition of chemokine receptors, chemokine receptor knockout mice and mice lacking key effector molecules we aim to understand the molecular processes underlying T cell clustering and parasite killing.
One of the most controversial and difficult problems regarding the mathematical modelling of biological systems is determination of the appropriateness of the choice of mathematical framework and the implementation of scale in any given model. Concentration scales, spatial scales and temporal scales characteristic of models in biology, more so than any other application, are highly varied and structured even in regards to a single process. In this minisymposium we will be hearing from speakers working on the forefront of foundations in multiscale modelling in biology. It is best to keep an open mind in these talks about the fundamental choice of the methodology which is appropriate to a system and how these methodologies relate to each other. You will be stimulated with ideas regarding discrete versus continuous and stochastic versus deterministic modelling, as well as more general discussion about the nature and treatment of multiscale frameworks in biology.
In many biological applications, it is useful to model chemical reactions at the level of individual molecules whilst not running costly Brownian dynamics models. Modelling chemical reactions as interactions at a distance can lead to errors if multiple species and time scales are important to the kinetics (for example, catalytic reactions). Multiple time scale analysis of equivalent mean-field kinetics can lead to high order reactions which are not just governed by the law of mass-action but also nonlinear kinetics such as Michaelis-Menten kinetics. In this talk we will discuss some particle-based methods which recapitulate these high order kinetics.
I will discuss mathematical and computational methods for spatio-temporal modelling in molecular and cell biology, including all-atom and coarse-grained molecular dynamics (MD), Brownian dynamics (BD), stochastic reaction-diffusion models and macroscopic mean-field equations.
Microscopic (BD, MD) models are based on the simulation of trajectories of individual molecules and their localized interactions (for example, reactions). Mesoscopic (lattice-based) stochastic reaction-diffusion approaches divide the computational domain into a finite number of compartments and simulate the time evolution of the numbers of molecules in each compartment, while macroscopic models are often written in terms of mean-field reaction-diffusion partial differential equations for spatially varying concentrations.
I will discuss the development, analysis and applications of multiscale methods for spatio-temporal modelling of intracellular processes, which use (detailed) BD or MD simulations in localized regions of particular interest (in which accuracy and microscopic details are important) and a (less-detailed) coarser model in other regions in which accuracy may be traded for simulation efficiency [1,2,3]. I will discuss error analysis and convergence properties of the developed multiscale methods, their software implementation [4] and applications of these multiscale methodologies to modelling of intracellular calcium dynamics [5], actin dynamics [6,7] and DNA dynamics [8]. I will also discuss the development of multiscale methods which couple MD and coarser stochastic models in the same dynamic simulation [3,9].
[1] M. Flegg, S.J. Chapman and R. Erban (2012). Two Regime Method for optimizing stochastic reaction-diffusion simulations. Journal of the Royal Society Interface 9: 859-868.
[2] B. Franz, M. Flegg, S.J. Chapman and R. Erban (2013). Multiscale reaction-diffusion algorithms: PDE-assisted Brownian dynamics. SIAM Journal on Applied Mathematics 73: 1224-1247.
[3] R. Erban (2014). From molecular dynamics to Brownian dynamics. Proceedings of the Royal Society A 470: 20140036.
[4] M. Robinson, S. Andrews and R. Erban (2015). Multiscale reaction-diffusion simulations with Smoldyn. Bioinformatics 31: 2406-2408.
[5] U. Dobramysl, S. Rudiger and R. Erban (2016). Particle-based multiscale modeling of calcium puff dynamics. Multiscale Modelling and Simulation 14: 997-1016.
[6] R. Erban, M. Flegg and G. Papoian (2014). Multiscale stochastic reaction-diffusion modelling: application to actin dynamics in filopodia. Bulletin of Mathematical Biology 76: 799-818.
[7] U. Dobramysl, G. Papoian and R. Erban (2016). Steric effects induce geometric remodeling of actin bundles in filopodia. Biophysical Journal 110: 2066-2075.
[8] E. Rolls, Y. Togashi and R. Erban (2017). Varying the resolution of the Rouse model on temporal and spatial scales: application to multiscale modelling of DNA dynamics. Multiscale Modelling and Simulation 15(4): 1672-1693.
[9] R. Erban (2016). Coupling all-atom molecular dynamics simulations of ions in water with Brownian dynamics. Proceedings of the Royal Society A 472: 20150556.
Patch dynamics is a numerical multiscale solver which constructs a macroscale solution of a microscale system by solving the original microscale problem, but only within discrete patches. These patches are spread across the domain of the system and are separated by the desired macroscale spacing, thus providing a description of the system at the macroscale. The space between the patches is unsimulated, which reduces the size of the problem and permits a reduction in computational time. For the patches to accurately capture the macroscale behaviour of the system, the patches must be coupled across the unsimulated space. In developing patch dynamics, one of the main tasks is to ensure the patch coupling is both accurate and efficient. We consider the example of the flow of fluid, such as blood, through a long narrow tube with fixed obstructions. We show how to apply patch dynamics to obtain an accurate description of the fluid flow over long space and time scales.
Stochastic simulations are essential to the study of biological cells, yet there is no computational framework allowing for detailed spatial simulations of genetic regulatory network within large populations of cells.
We fill this gap by developing a parallel simulation framework capable of spatially resolved stochastic simulation of cell-cell signalling in multicellular systems. We use an operator-splitting method to decouple the internal reaction-diffusion kinetics from the interactions on the cells’ boundaries and allow for efficient and horizontally scalable simulations of large numbers of interacting cells. Our framework is highly compatible with many existing methods and allows for hybrid simulation where both coarse and detailed models are considered at the same time. It is also greatly versatile and is deployable on various high performance computing platforms, such as clusters or clouds.
We use a small test model to study the convergence of our method as well as larger models to demonstrate weak scalability.
Our method demonstrate the feasibility of detailed stochastic simulations of large populations of interacting cells and is the first step toward more complete simulations including both detailed reaction-diffusion simulations and cell mechanics.
Population dynamics are a traditional focus of mathematical biology, dating from the pioneering work of Volterra and Lotka. Importantly, the rates of change of species' abundances have broader consequences, affecting a diverse range of processes that includes evolutionary, behavioural, ecosystem and human dynamics. In this minisymposium, we highlight biomathematical models that couple population dynamics with these broader factors. Not only do these coupled feedback models help us understand the behaviour of broader systems, they are also essential to an accurate understanding of the population dynamics themselves, and how to manage them. Developments in this modelling space can help us understand essential characteristics of ecological systems, and also generate mathematically interesting complex system behaviours, which could stimulate new theoretical hypotheses.
In this session, three speakers will demonstrate population dynamic models that are integrated with intrinsic processes of organisms, i.e., evolution, behaviour and inducible responses. Two speakers will present about joint dynamics of populations and management policy, in the context of biodiversity conservation and agricultural damage control.
Two consumer species that share a single resource species can indirectly interact each other, even without direct interactions. A typical indirect interaction is exploitative resource competition that results from a depression of resource biomass by consumption, which can be referred to as “biomass-mediated indirect effect”. Another type of indirect interaction is called “trait-mediated indirect effect”, which is caused by changes in traits of resource species. For example, an attack from one herbivore species can either reduce biomass or induce defensive traits of the host plant individual, which may in turn suppress performance of another herbivore species on the same host. It is notable that two types of indirect effects are difficult to detect separately in nature because it is necessary to measure both reduction in biomass and changes in traits of plants and observe respective effects on herbivore performance. In order to understand the separate and joint influences of two indirect effects, we analyzed a population dynamic model of one-plant two-herbivores system, including an inducible non-specific defense of plant individuals. Our analysis revealed (1) the inducible non-specific defense can promote coexistence of three species, (2) the coexistence of three species requires reduction of population densities of both herbivores, (3) in such a case, one herbivore is controlled by biomass-mediated negative indirect effect, whereas the other herbivore is controlled by trait-mediated negative indirect effect. This indicates that both types of indirect effects need to co-occur for the coexistence of consumer species, suggesting that trait-mediate indirect effect may be as common as biomass-mediated indirect effect.
It has long been debated whether introduction of two (or more) natural enemies results in more efficient pest control than that of either one enemy. Intra-guild predation (IGP) among two natural enemies sharing a single pest has been recognized as an important factor to reduce efficiency of pest control. While the classical theoretical model of IGP showed that introduction of two natural enemies is less efficient for pest control, empirical works have showed positive, negative, and neutral results. However, the classical IGP model did not consider any behavioural plasticity of the pest and natural enemies.
In this study, we extended the classical IGP model by considering adaptive defense by the pest and switching predation by the omnivorous natural enemy (omnivore) and examined separate and joint effects of them on efficiency of pest control. We assumed that the pest can adaptively allocate efforts toward two kinds of defense respectively against the two natural enemies to increase its own fitness, with cost of reduction in its own reproduction. Switching predation by the omnivore is expressed as the Holling’s type III functional response to the pest and another natural enemy (intermediate predator). Equilibrium pest density is used as an index of efficiency of biological pest control.
If the pest employed the adaptive defense, introduction of two natural enemies performed more efficient pest control than that of either one enemy, unless the IGP was too intensive. This is because the pest allocated more defensive effort against the more threatening enemy, and thus it was difficult for the pest to simultaneously prevent predation from the two enemies. On the other hand, switching predation could not improve biological control by two natural enemies. However, the type III functional response tempted the pest to abandon defensive effort against the omnivore, because the predation pressure was negligible at low pest density and saturated at high density. Since the intermediate predator was suppressed by the undefended omnivore, defensive effort against the intermediate predator also diminished. Thus, switching predation could offset the effects of defense and the prey density could be lowered by two natural enemies even under severe IGP. Consequently, types and combination of behavioural plasticity might cause qualitatively different outcomes of biological control introducing two natural enemies.
Sexual dimorphism (SD), sexual differences in traits, represents one of the most remarkable source of biodiversity in the world. In most of species, two sexes play different roles in reproduction and thus are imposed by selection pressures in different forms and strengths. This explains the generality of this phenomenon and the extensive applicability of SD theory in nature. In this talk, I will provide an example demonstrating how sexual differences in the dispersal strategy can influence metapopulation dynamics. We developed a sex-structured two-patch model and specified that two sexes have different decision-making processes according to the accessibility of reproductive resources; in brief, females move to the patch with more breeding sites (i.e. carrying capacity minus local population size) whereas males move to the patch with more mating opportunities (i.e. unmated females). Such assumption is valid given that sex-specific dispersal has been commonly recorded in natural systems. Our simulation results illustrated how sex-specific dispersal increased or decreased local population size through changes in local sex ratio and individual reproductive success; notably, such results were not able to be generated through conventional, asexual metapopulation models. We advocate that incorporating the perspective of SD into population dynamic models are needed and it can advance our knowledge of population dynamics by disclosing detailed mechanisms underlying population processes.
Overpopulated mammal populations cause damage to agricultural crops in Japan. We need to keep population sizes at appropriate levels. When we determine management plans, we have to deal with various uncertainties such as population size, population growth rate, and agricultural damages caused by mammals due to lack of sufficient data. It is important to reduce those uncertainties and allocate hunting effort under limited budget.
Prey populations change their behaviour due to presence of predators. Predation risk could cause hehavioural changes of prey such as vigilance at the cost of time spent grazing, and population growth rates of prey may decrease without predation itself (the effect is called as behaviourally mediated indirect effect). It is possible that population decline in Japan decreases negative impacts of human activities on mammal populations. In addition the agricultural damage could affect human activities. Hence it is important to incorporate the behaviourally mediated indirect effect for effective mammal managements.
We will introduce a population model incorporating a behaviourally mediated indirect effect for a management of wild boars (Sus scrofa) in Chiba prefecture, Japan. First we estimate population densities and population growth rates using Bayesian model averaging. Second we estimate relationship agricultural damages and population densities in heterogeneous landscapes. Last we derive optimal allocations of hunting efforts to minimize expected damages to agricultural crops. We discuss how incorporating the behaviourally mediated indirect effect is important for effective wild boar management.
Oscillation of lateral asymmetry diorphism is first found in scale eating cichlid, Perrisodus microlepis. Fraction of their lefty morph oscillates around 0.5 in about 5 year period. Other fish or aquatic invertebrates also reveal lateral asymmetry dimorphism and oscillation of morph fractions. One of key factors of oscillation is cross predation dominance: lefty predator eats more righty prey than lefty prey, and vice versa. When lefty is major in predator, righty prey decrease, then larger fraction of lefty prey increases righty predator, and so on. This cycle appears to lead to the laterality oscillation. However the story is not such simple. A simple ODE model of prey predator system with cross predation dominance has no limit cycles, its coexisting equilibrium is stable, though almost neutral. Introduction of time delays due to growth periods destabilizes the coexisting equilibrium and leads to a limit cycle of oscillating dimorphism fraction. However, this laterality dimorphism oscillation due to time delays of growth periods and predation cycle can not explains that observed in fields. Amplitude and period of fraction oscillation are almost 1.0 and 50 ~ 250 years in the model, while they are 0.1 ~ 0.3 and 3 ~ 6 years in fields. Another factor of oscillation in morph fraction is frequency dependence. When lefty scale eaters are more common, righties eat more scales and have higher reproductive success. We introduce frequency dependent predation success into the model: rarer predator eats more prey. This frequency dependence, however, ends up with the vanish of oscillation by either stabilization of coexisting equilibrium, or fixation one laterality morph. Another frequency dependence, that of prey selection, i.e. predator eats more common prey morph, stabilizes coexisting equilibrium. The model with both of frequency dependence finally shows oscillations with realistic amplitudes and periods. We conclude that oscillation of lateral asymmetry morphs is caused not by predation cycle, but by frequency dependence in both of prey and predator, i.e. rarer predator or prey morph is advantageous and become common after a few years.
The fall armyworm (Spodoptera frugiperda) is a pest insect which has the propensity to destroy a wide variety of common crops. It ranges over Eastern and Central North America and, since 2016, has been invasive in Africa resulting in significant economic damage. The fall armyworm is susceptible to Bt derived insecticides, making Bt modified crops a viable method for controlling this species. However, there is evidence to suggest that there is a small subpopulation of the fall armyworm which is resistant to Bt. This population remains relatively small in the wild due to the cost of resistance, which is partially mediated by the aggressive cannibalism exhibited by the larvae, since resistant larvae grow at a slower rate. Thus it may be possible to control the rise of resistance to Bt in a fall armyworm population via the creation of refuges for the nonresistant larvae who will in turn cannibalize the resistant larvae.
Here we use an individual base modelling (IBM) approach to model an infestation of fall armyworms, which models every individual’s entire life cycle including cannibalistic encounters and the effects of its Bt resistance genotype. First, we explore the best underlying mathematical models to describe the individual insects in the IBM. Then we examine the problem of controlling the accumulation of resistance to Bt via adjustments to the size and spacial layout of the refuge while limiting the damage to the effected crops.
In this talk, we consider the dynamics of a Lotka-Volterra prey-predator model by a class of delay differential equations. The number of prey varies due to a general nonlinear predators' consumption rate with delays. Under the assumption that the consumption rate is monotonically increasing with respect to the number of prey, we investigate the effect of the nonlinearity and delays on the asymptotic behaviour of the model. For the case where the consumption rate is described not only as Holling type I,II but also as Holling type III, some ongoing studies are also introduced. Comparison for the assumptions on the incidence rates appearing in prey-predator models with those in epidemiological models, will be also discussed.
The Hassell equation is a classic discrete-time population model which has been widely used to model population dynamics of species with seasonal reproduction. This equation is a generalization of the Beverton-Holt equation with an additional exponent, and can describe various types of reproduction curves exhibiting from exact-compensation (contest competition) to over-compensation (scramble competition) depending on the value of the exponent. The value of the exponent and the resulting density dependence is thought to be related to the degree of inequality in resource allocation among individuals; contest curves result from unequal resource allocation, and scramble curves, from equal allocation. However, as the model is a phenomenological one at the population level, this relation between the exponent and the inequality in resource allocation has mostly been discussed only phenomenologically. Although some authors have actually derived the Hassell equation from first-principles by considering specific competition models among individuals, the exponent of the derived models does not match the naive expectation above. This study explores whether it is possible to derive from first-principles such a Hassell model that its exponent is related to the inequality in resource allocation by considering resource competition among individuals. I demonstrate that such a Hassell model can indeed be derived by assuming that each individual obtains a constant amount of resources (a resource unit) at a time, and that the competition for such a discrete resource unit among individuals is repeated many times. Different sizes of the resource unit generate different degrees of inequality, and the exponent of the derived model turns out to depend on that size, thus being related to the degree of the inequality. The derived model reduces to the Beverton-Holt model when the inequality is highest, and to the Ricker model when the inequality is lowest. Finally, I discuss how replacing an assumption on fecundity with more realistic one changes the functional form of the derived model, and the extension to two-species competition models as well.
Introduced species are a critical threat to Australian ecosystems and species. Particularly noxious examples include the European carp, feral cats, and a variety of weeds. A central aspect of introduced species management is eradication – if they can be completely removed from a region, the impact can be nullified. A central problem population eradications is knowing whether the species has been successfully removed or not. We develop a framework to model populations through time, explicitly accounting for imperfect detection and unknown detection probability. We use changing detection rates throughout a removal project to calibrate the model, which provides a quantitative method to trigger the end of a project. While invasive species are often the focus of removal efforts, they can also occur to prevent disease spread in an endangered species. I will describe how we applied this method to a Tasmanian devil depopulation, which enabled the establishment of a Tasmanian devil facial tumour disease population on Forestier Peninsular, Tasmania.
In this study, we propose a new Bayesian approach to calculate the expected extinction time of a species based on historical sighting record data. Unlike other work, our model allows comprehensively for uncertainties, provides the expected extinction time and the probability of extinction. It is extremely difficult to determine whether a species is extinct based on historical sighting records because knowing whether the last surviving individual of a species has finally died, or is just unobserved, remains problematical. Moreover, an incorrect classification of a species as extinct can lead to failure in conserving a threatened species. On the other hand, it is also undesirable to classify a species as extant when it is actually extinct as it can lead to misallocation of funds. Sightings with uncertain validity (uncertain sightings) play an important role in the inferences made and need to be taken into account better than they have been. Recent studies have considered uncertain sightings while making inferences about extinction; however, the difficulty of the problem requires making a number of limiting assumptions that significantly reduce realism. We have attempted to derive a more general model by incorporating sighting validity into a new Bayesian model development. We employ the underlying idea of the beta-geometric/beta-binomial (BG/BB) model to build our Bayesian approach for the analysis of extinction. Using the likelihood for the sighting data, along with the prior distributions of sighting probability and extinction probability, we calculate the expected extinction time using a Markov Chain Monte Carlo (MCMC) method which we simulate with JAGS. We apply this new approach to the sighting records of Caribbean monk seal (CMS), Dodo and Ivory-Billed Woodpecker species. Unlike other approaches that consider uncertainties in the sightings, our model gives Bayesian confidence intervals for the expected extinction time of a species.
This minisymposium will bring together researchers to examine up-to-the-minute disease problems that showcase the usefulness and applicability of mathematical modelling to a world far beyond the mathematical community. The audience is the mathematical biologist with an interest in infectious disease. This includes students and researchers, mathematicians interested in seeing applications and biologists who wish to see how mathematics can be used to solve real problems. The minisymposium is interdisciplinary in nature and includes those trained as mathematicians, epidemiologists and immunologists.
There is an urgent need for more understanding of the effects of surveillance on malaria control. Indoor residual spraying has had beneficial effects on global malaria reduction, but resistance to the insecticide poses a threat to eradication. We develop a model of impulsive differential equations to account for a resistant strain of mosquitoes that is entirely immune to the insecticide. The impulse is triggered either due to periodic spraying or when a critical number of malaria cases are detected. For small mutation rates, the mosquito-only submodel exhibits either a single mutant-only equilibrium, a mutant-only equilibrium and a single coexistence equilibrium, or a mutant-only equilibrium and a pair of coexistence equilibria. Bistability is a likely outcome, while the effect of impulses is to introduce a saddle-node bifurcation, resulting in persistence of malaria in the form of impulsive periodic orbits. If certain parameters are small, triggering the insecticide based on number of malaria cases is asymptotically equivalent to spraying periodically.
The number of HIV+ individuals who develop end stage renal disease (ESRD) and require life-long dialysis treatment has been continually rising in some regions around the world. A differential equation-based mathematical model was developed to assess the impact of antiretroviral therapy on the progression to disease and to project the future prevalence of HIV+ ESRD. The goals of this study are to parameterize the model with new data on the populations of individuals with HIV/AIDS and those with HIV+ ESRD, and to extend the model to take into account greater complexity in the population dynamics. We also expand the model’s analysis to predict treatment recommendations for the population of HIV+ individuals at risk for developing renal disease. Further studies will use the model to estimate how much the development of renal disease and treatment levels have changed over time, and to predict the treatment levels needed to slow the increase of this epidemic.
Clinical trials of the four-dose RTS,S/AS01 vaccine for P. falciparum malaria demonstrated a protective effect in young children and, beginning in 2018, the vaccine will be evaluated through a large-scale pilot implementation program in Ghana, Kenya and Malawi. Recent evidence from a phase 2a challenge study indicates that varying the timing and amount of the fourth dose could further improve the public health impact.
We used a dynamic modelling approach to inform target product profiles (TPPs) for a second-generation malaria vaccine, focussing on the vaccine properties of initial efficacy, duration of protection, dosing schedules and coverage. We simulated the changing anti-circumsporozoite antibody titre following vaccination, related titre to vaccine efficacy, and then implemented this efficacy profile within an individual-based model of malaria transmission. We developed a range of efficacy profiles for different vaccine schedules and used these to evaluate the relative impact of initial efficacy, duration, and fourth dose characteristics. We ran the simulations across a range of epidemiological strata, measuring clinical cases averted children younger than five years.
We found that in the first decade of delivery, a vaccine with high initial efficacy vaccine resulted in more clinical cases averted, compared to a longer duration vaccine, and that the effect was more pronounced in high malaria prevalence settings. However, the low initial efficacy and long duration schedule averted more cases across all age cohorts if a longer time horizon was considered. We also observed that a higher antibody titre resulting from the fourth dose was always advantageous, and that a longer delay between doses three and four averted more cases in older age classes.
Our results indicate that an imperfect malaria vaccine, initial efficacy may be more important than vaccine duration. An optimised malaria vaccine may outperform the current RTS,S, but timing of the fourth dose determines the age group that benefits most. This TPP analysis could provide insight for vaccine developers and policy-makers into how distinct characteristics of a malaria vaccine may translate to public health outcomes.
Multiphase models originated in fluid dynamics, but over many years they have also been applied to the study of biological tissues. The appeal of the multiphase approach in a biological context is that it provides a natural framework to account for mechanical interactions between tissue constituents such as cells and extracellular matrix. This is an important consideration because mechanical effects often play a crucial role in determining both cell function and tissue evolution during development and disease. Much of the research in multiphase modelling over the last 20 years has focused on tumour growth, but recent work has seen multiphase models applied in new areas such as tissue engineering and atherosclerosis. This minisymposium aims to bring together researchers with a wide range of interests to showcase the insights that the multiphase approach can provide, and encourage the sharing of ideas to tackle future challenges in the field.
Atherosclerotic plaque growth is characterised by a process of chronic, non-resolving inflammation that leads to the accumulation of cellular debris and extracellular fat in the inner artery wall. In advanced plaques, smooth muscle cells (SMCs) are recruited from deeper in the artery wall to synthesise a fibrous tissue cap that sequesters the thrombogenic plaque content from the bloodstream. The fibrous cap therefore provides crucial protection from plaque rupture and the formation of blood clots that occlude vessels and cause heart attacks and strokes. Despite the important role played by the plaque fibrous cap in preventing the clinical consequences of atherosclerosis, the mechanisms that underlie cap formation remain poorly understood. In particular, it is unclear why certain plaques become strong and stable while others become fragile and dangerously vulnerable to rupture.
In this talk, we discuss the use of a multiphase approach with non-standard boundary conditions to investigate early fibrous cap formation in the intimal layer of the artery wall. We model the highly nonlinear process of SMC migration from the media in response to a diffusible chemical signal produced at the endothelium. Simulations indicate that the emergence of a stable fibrous cap requires a critical balance between the relative rates of cell supply from the media, chemotactic migration within the intima and cell loss by apoptosis (or phenotype change). Moreover, we identify a number of disease-associated parameters that may be linked to variations in cap stability. This model represents the first detailed in silico study of fibrous cap formation in atherosclerosis, and establishes a framework that can be extended to investigate other aspects of plaque development.
Tissue engineering aims to grow artificial tissues to replace those that have been damaged through age, trauma or disease. A recent approach to engineer artificial cartilage involves seeding cells within a scaffold consisting of an interconnected three dimensional printed lattice of polymer fibres combined with a cast or printed hydrogel, and subjecting the construct (cell-seeded scaffold) to an applied load in a bioreactor. A key question is understanding how the applied load is distributed throughout the construct to the mechanosensitive cells.
To address this question, we view the scaffold as a two-phase mixture, in which the fibres are modelled as a linear elastic material and the hydrogel as a poroelastic material. We exploit the disparate length scales (small inter-fibre spacing compared with construct dimensions) and use homogenisation theory to derive macroscale equations for the effective mechanical properties of the composite material. The resulting governing equations reflect the orthotropic nature of the composite. We then validate the model, by comparing results generated from finite element simulations of the macroscale, homogenised equations to experimental data describing the unconfined compression of the fibre-reinforced hydrogels.
Patients diagnosed with glioblastoma multiforme (GBM) are expected to survive only 14 months and die due to the pressure that the tumour builds in the brain as well as the formation of peritumoural edema (PTE). With the view to investigating the early stages of brain tumour development, and how it impacts the healthy brain environment, we develop a mechanistic model of GBM onset. The model is derived using principles of mass and momentum balances and explicitly includes pressure dynamics within the disease brain and the ability/inability of healthy tissue to repair itself in response to these cues. As a first step we assume an implicit tumour that exerts pressure at a healthy boundary causing the boundary to move into healthy tissue with a velocity $v$ (thought of as the tumour growth rate). We investigate three velocity regimes: where $v$ is an order of magnitude slower than the time-scale of healthy brain tissue renormalization (benign tumour); where $v$ is an order of magnitude higher than the time-scale of healthy brain tissue renormalization (high grade tumour); a transition between these where $v$ is the same magnitude as the time-scale of healthy brain tissue renormalization. Our model shows a correlation between the tumour velocity and the formation of PTE, which is an indicator of tumour malignancy. The resulting model includes time-varying diffusion on a moving domain, which presents unique numerical challenges. We propose a scheme to solve such equations, validating our method with a test problem as well as theoretical analysis using techniques from asymptotic methods in order to complete this research aim.
Peripheral nerve damage afflicts 1 M people p.a. in Europe and the USA [1]. In the most severe cases, patients experience major loss of function. The gold-standard treatment for patients with severe cases is surgery using a graft based on a healthy section of nerve taken from the patient, however, only 50% of patients experience functional recovery [2].
After a nerve is severed, the distal stump produces chemotactic cues that can stimulate neurite regeneration. At the proximal stump, neurites start growing and their direction is informed by these cues. If the gap is short enough, the neurites are able to respond to their environmental cues and grow across the injury site to enter the distal stump. However, if the gap is larger than 5 mm, the regenerative process is unable to establish a supportive environment for neurite growth.
The limitations of the current gold-standard therapy to establish functional recovery have motivated the development of engineered replacement tissues and repair conduits to promote regeneration. Several designs have been proposed that vary not only in terms of material and therapeutic cell composition but also in their spatial distribution. Despite promising results in animal models, these conduits have not yet progressed to clinical use, in part because of the expanse of different variables that need to be tested using experimental models [3].
Here a multiscale mathematical model is proposed to streamline the design of these conduits. The model comprises of two components. The first is a continuous 3D model of the chemotactic field within the conduit. The release of chemotactic factors from the distal stump is captured using a flux boundary condition, as well as the permeability of the conduit wrap.
Within this geometry, at the proximal stump, a random walk model for individual neurites is initiated. Neurite growth is simulated as a random variable that is biased towards the direction of increasing cue; this cue could either be chemotactic (as described using the continuum model above), or durotactic (due to the underlying material distribution). Neurite branching is incorporated into the model, as this is a fundamental feature during regeneration, which determines how neurites sense the underlying gradient fields. Morphological features of branching are accounted for by implementing the model by Van Pelt et al., whilst branches can also recede if they encounter a sub-optimal environment, modelled using a Markov chain process [4,5].
The discrete-continuum framework is parameterised against in vitro and in vivo data on neurite progression through engineered conduits, including neurite counts at the proximal and distal end as well as measurements of neurite lengths. The parameterised framework can be used to explore the optimal arrangement of materials and cells to promote efficient neuronal regeneration following injury.
[1] Chen, Shan-lin, et al. Neural regeneration research, 10.11 (2015).
[2] Grinsell, D., and C. P. Keating. BioMed research international (2014).
[3] Angius, D., et al. Biomaterials 33.32 (2012).
[4] Van Pelt, J., et al. The Journal of Comparative Neurology 387.3 (1997).
[5] Britto, J.M., et al. Biophysical Journal 97.3 (2009).
Network of transport is found in a wide range of living system. The typical examples are vascular network in vertebrates and tracheal network in vascular plants, and mycelial network in fungi. Such a network structure is found in self-organized colony of unicellular bacteria: a biofilm that is a sheet-like aggregate of many bacteria and sticky polysaccharides secreted from the bacteria. Comparative study of various kinds of bio-network of transport along phylogenetic tree is probably interesting since an insight into evolutionarily sophisticated designing of adaptive network is obtained in a wider perspective.
At the interface between unicellular and multicellular branches, there is an unique organism with highly flexible body shape with network structure, true slime mould Myxomycete (or Mycetozoa). Physarum polycephalum Schw. is a well-studied species of true slime mould in cell biology and biophysics last seventy years although several hundreds of species are known. By using this organism, comparative study of bio-network rises for the last one or two decade(s).
A kind of huge amoeboid organism named Physarum plasmodium constructs a intricate network of veins for circulating nutrients and signals over the entire body. The network shape (topology of connectivity, and sequence of branching in vein network, for instance) is drastically re-organized within hours in response to external conditions. The past studies showed that the network shape was optimized to maximize possibility of survival, in some senses. So we may extract an algorithm for optimal design of functional network from the primitive organism. The key thing is adaptive dynamics of current-reinforcement rule: each vein of network becomes thicker when current is large enough through the vein itself, while it becomes thinner and dies out otherwise.
We propose the equations of motion for this simple rule, and functions and formation of transport network is analyzed. We will show that the rule is applicable to the other bio-systems: (1) social dynamics of public transportation, (2) formation of network structure in porous tissues bone (bone remodelling in other words). A tractable perspective to think similarly of a variety of bio-network is given from the viewpoint of current-reinforcement rule.
In this paper, a mathematical model of breast cancer governed by a system of ordinary differential equations in the presence of chemotherapy treatment and ketogenic diet is discussed. Several comprehensive mathematical analysis was carried out using varieties of analytical methods to study the stability of the breast cancer model. Also, sufficient conditions on parameter values to ensure cancer persistence in the absence of anti-cancer drugs ketogenic diet and cancer emission when anti-cancer drugs, immune-booster, ketogenic diet are included were established. Furthermore, optimal control theory is applied to find out the optimal drug adjustment as an input control of the system therapies to minimize the number of cancerous cells by considering different controlled combinations of administering the chemotherapy agent and ketogenic diet using the popular Pontryagin’s Maximum Principle. Numerical simulations are presented to validate our theoretical results.
The advances in cell imaging technology has allowed for deeper understanding of the movement mechanisms of immune cells. Experimental evidence suggests that cytotoxic T lymphocytes and dendritic cells undergo and unrestricted search motion until they switch to a more restricted motion induced by activation by tumour antigens. This change in movement is not often considered in the existing mathematical models of the interactions between immune cells and cancer cells. We present a spatially structured individual-based model of tumour-immune competition that takes explicitly into account the difference in movement between inactive and activated immune cells, using Lévy walks and Brownian motion to capture this motion, respectively. The effects of activation of immune cells, the proliferation of cancer cells and the immune destruction of cancer cells are also included in the model. We illustrate the ability of our model to reproduce qualitatively the spatial trajectories of immune cells observed in experimental data of single cell tracking. Computational simulations of our model further clarify the conditions for the onset of a successful immune action against cancer cells and suggest possible targets to improve the efficacy of cancer immunotherapy. Overall, our theoretical work highlights the importance of taking into account spatial interactions when modelling the immune response to cancer cells.
Clinical methods for assessing tumour response to neoadjuvant therapy largely rely on monitoring the temporal changes in tumour size. Our goal is to predict tumour response to neoadjuvant therapy in breast cancer using a mathematical model that utilizes non-invasive imaging data obtained from individual patients. Previously, a mechanically-coupled, reaction-diffusion model with logistic growth for breast cancer was shown to outperform clinically used measures for predicting tumour response to therapy when initialized with patient specific magnetic resonance imaging (MRI) data. This model has since been extended to include patients’ therapy using corresponding dynamic contrast-enhanced (DCE-) MRI data to determine areas for drug delivery, improving the concordance correlation coefficient of the predicted tumour cellularity compared to the calculated cellularity from 0.85 to 0.99 (p<0.01, N=5).
We use patient specific measurements from diffusion-weighted MRI data to calibrate tumour cell diffusion, carrying capacity, and proliferation values (on a per voxel basis) within each patient’s tumour domain. The model is 3D and is mechanically coupled to the properties of the breast tissues by dampening the diffusivity of tumour cells based upon their growth and stresses induced in the breast. DCE-MRI data is used to identify spatiotemporal variations in tumour perfusion to approximate delivery of therapy within the tumour using the extended Kety-Tofts model. An explicit death term is incorporated into the tumour cell equation, where the amount of drug distributed in the tissue is spatially heterogeneous and dynamic in time according to the patient’s therapy schedule. All model simulations were evaluated using a forward finite difference scheme.
For a cohort of breast cancer patients (N = 10), four MRI scans were collected over the course of their neoadjuvant therapy (pre-treatment, after one cycle of therapy, at the end of the first therapy, and after one cycle of the second therapy). Using this data, we assess the model’s predictive ability in two different ways: 1) using the first two scans, the model is calibrated and simulated forward to the third scan time to assess the tumour cellularity, size, and heterogeneity predicted by the model in comparison with the calculated values from the patient’s third scan; 2) the model is calibrated with the third and fourth scans and simulated forward to the time of surgery, where the model’s results are compared to the patient’s pathological response status. Using these metrics, we will report on and discuss ability of the modelling approach to predict tumour response to therapy in breast cancer.
We acknowledge the support of U01 CA174706, U01 CA154602, and CPRIT RR160005.
Tumour associated macrophages have long been implicated in the progression of primary solid malignancies including prostate cancer. Metastatic prostate cancer typically manifests in the bone where it induces painful osteogenic lesions that are incurable. Bone is naturally rich in myeloid derived macrophages whose temporal polarization into pro- (M1) and anti-inflammatory (M2) phenotypes is critical for regulating the program of bone injury repair mediated by bone resorbing osteoclasts (OCL) and bone building (OBLs). However, the dynamics of macrophage polarization in the context of bone metastatic prostate cancer is underexplored and difficult to address with traditional biological approaches. In order to address it we built a mathematical model describing bone resident cell population dynamics during bone injury response through a set of coupled ordinary differential equation (ODE). We informed this model by analyzing macrophage plasticity in vivo subsequent to intratibial injury. Bone marrows were isolated at several time points and profiled by flow cytometry for pro- and anti-inflammatory macrophage content. Contralateral tibias were analyzed for bone volume, osteoblast and osteoclast numbers. The ODE model was able to predict experimental observations of M1/M2, OBL, OCL and bone dynamics. Generation of the ODE model required testing a number of assumptions regarding macrophage polarization behaviour. For example, it is unknown whether M1 resolve at the initial stages of bone injury repair through death, or by re-polarizing into M2. For each aspect, competing mechanistic assumptions were proposed and simulated mathematically by sets of ODEs. The best fitting assumptions for each aspect were integrated into a comprehensive ODE model to fully describe the dynamics of the bone resident cell populations during bone repair. This experimentally validated model is now being exploited to address how bone population dynamics respond to metastatic prostate cancer cells and subsequently identify the resultant impact on bone formation and cancer growth.
The cellular cytoskeleton ensures the dynamic transport, localization and anchoring of various proteins and vesicles. In the development of egg cells into embryos, messenger RNA (mRNA) is transported along microtubule filaments and must accumulate at the cortex of the egg cell on a certain time and spatial scale. We present two equivalent methods of deriving the effective transport properties of mRNA at large time: using analysis of partial differential equations, and using renewal rewards processes in a stochastic model formulation. The dynamical systems model approach can be extended to include the geometry of the microtubule filaments. This allows us to better predict the spread of the particles, and to investigate the contribution of an anchoring mechanism to the timescale of localization. Our numerical studies using model microtubule structures predict that anchoring of mRNA-molecular motor complexes may be most effective in keeping mRNA localized near the cortex and therefore in healthy development of oocytes into embryos.
Autophagy is an intracellular degradation process mediated by the autophagosome. The membrane dynamics of autophagosome formation is unique and complicated, which involves development of a small membrane cisterna into a cup-shaped structure and a double membrane spherical structure by closing the edge [1].
In this presentation, we discuss the mechanism of autophagosome formation from a physical point of view. A flat cisterna has a highly-curved rim, which is energetically unstable. As increasing the membrane surface, the rim area also increases. In order to reduce the energy, the cisterna is spontaneously curled by decreasing the rim area. Accordingly, a shape transition occurs suddenly due to the first order transition [2].
However, live-imaging experiments have shown that autophagosome formation takes place gradually and sudden closure of a flat membrane is not observed. We hypothesize the presence of protein(s) that stabilize the highly-curved rim, which enables gradual enclosure. Indeed, several autophagy proteins are present at the rim.
Thus, we consider the effects of the rim-stabilizing proteins on the dynamics of autophagosome formation. We show that the proteins localize at the highly-curved rim and stabilize it. As a result, intermediate cup-shape states appear and the closing dynamics becomes moderate similar to that of autophagosome formation in vivo.
[1] Mizushima 2017, PMID28186333
[2] Knorr et al., 2012, PMID22427874
In this presentation, a mathematical model of contractile ring-driven cytokinesis is presented by using both phase-field and immersed-boundary methods in a three-dimensional domain. It is one of the powerful hypotheses that cytokinesis happens driven by the contractile ring; however, there are only few mathematical models following the hypothesis, to the author’s knowledge. I consider a hybrid method to model the phenomenon. First, a cell membrane is represented by a zero-contour of a phase-field implicitly because of its topological change. Otherwise, immersed-boundary particles represent a contractile ring explicitly based on the author’s previous work. Here, the multi-component (or vector-valued) phase-field equation is considered to avoid the emerging of each cell membrane right after their divisions. Using a convex splitting scheme, the governing equation of the phase-field method has unique solvability. The numerical convergence of contractile ring to cell membrane is proved. Several numerical simulations are performed to validate the proposed model.
Nanoparticles provide a promising approach for the targeted delivery of therapeutic, diagnostic and imaging agents in the body. However, it is not yet fully understood how the physicochemical properties of the nanoparticles influence cellular association and uptake. Cellular association experiments are routinely performed in an effort to determine how nanoparticle properties impact the rate of nanoparticle-cell association. To compare experiments in a meaningful manner, the association data must be normalised by the amount of nanoparticles that arrive at the cells, a measure referred to as the delivered dose. The delivered dose is calculated from a model of nanoparticle transport through fluid. A standard assumption is that all nanoparticles within the population are monodisperse, namely, the nanoparticles have the same physicochemical properties. We present a semi-analytic solution to a modified model of nanoparticle transport that allows for the nanoparticle population to be polydisperse. This solution allows us to efficiently analyse the influence of polydispersity on the delivered dose. Combining characterisation data obtained from a range of commonly-used nanoparticles and our model, we find that the delivered dose changes by more than a factor of two if realistic amounts of polydispersity are considered.
Mechanically-induced buckling underlies the shape and function of a number of biological processes, such as brain tissue folding, intestinal crypt fission, and seashell formation. Unlike many typical engineering systems in which buckling is induced by an external compressive load, mechanical instability in biology is driven often by internal growth. In these contexts, it is just as important to consider how the system evolves beyond the onset of instability, as it is to understand when the instability occurs initially.
In this talk, we consider a growing, planar, elastic rod supported by an elastic foundation, as an appropriate starting point a number of possible biological systems. We analyse the post-buckling behaviour through a combination of weakly nonlinear analysis and numerical methods. The effect of different material parameters on the type of buckling and the buckled shape is shown. We then examine how spatial heterogeneity in substrate adhesion, elastic rod stiffness, and growth, impacts both the instability and resultant post-buckled shape evolution. These heterogeneities are contrasted with heterogeneity within the foundation shape, echoing classical results by Koiter on imperfection sensitivity and failure [1]. Finally, we discuss extensions that specialise the model to that of a buckling intestinal crypt.
[1] Van der Heijden, Arnold MA, ed. WT Koiter's elastic stability of solids and structures. Cambridge: Cambridge University Press, 2008.
The illegal trade of wildlife is estimated to run into hundreds of millions of dollars. It is an international problem that exploits both enforcement loopholes and corruption, and it is a direct threat to the survival of plant and animal species. Smugglers transport wildlife and their derivatives from sources to destinations across the globe. They are able to choose the route they take to move these commodities across international borders. The location where a smuggler enters a country is important as limited biosecurity resources are able to be placed at these locations. The locations and amount of biosecurity resources must be allocated intelligently. Quantitative analysis and mathematical methods are required to provide insights into illegal wildlife trade networks. In this talk, I will discuss a novel two player game on a network between smugglers and biosecurity agencies which looks at how best to allocate resources to intercept illegal trade. This framework will then be applied to some networks motivated by data.
Illegal exploitation of wildlife is one of the biggest threats to biodiversity, affecting over 2,000 species. Traditionally, actions to reduce poaching have focussed on increasing the efficacy and capacity of law enforcement. However, recently, non-governmental agencies (NGOs) are heavily investing in alternative interventions to reduce consumer demand for poached wildlife products. But which action is more efficient, policing or demand reduction? Using elephant poaching for ivory as a case study, and a simple open-access harvest model, we show that demand reduction interventions need to reduce ivory price by roughly 29% to be more cost-effective than investing in the police. We use an analytic derivation to show that this threshold is robust to several types of uncertainty but sensitive to specific socio-economic factors. The calculations can provide context for stakeholder engagement so that monitoring of ivory price can better inform decisions for protecting elephants from poaching.
Unlike any great apes, humans have expanded into a wide variety of habitats during the course of evolution, beginning with the transition by australopithecines from forest to savanna habitation. Novel environments are likely to have imposed hominids a demographic challenge due to such factors as higher predation risk and scarcer food resources. In fact, recent studies have found a paucity of older relative to younger adults in hominid fossil remains, indicating considerably high adult mortality in australopithecines, early Homo, and Neanderthals. It is not clear to date why only human ancestors among all hominoid species could survive in these harsh environments. In this talk we explore the possibility that hominids had shorter interbirth intervals to enhance fertility than the extant apes. To infer interbirth intervals in fossil hominids, we introduce the notion of the critical interbirth interval, or the threshold length of birth spacing above which a population is expected to go to extinction. We develop a new method to obtain the critical interbirth intervals of hominids based on the observed ratios of older adults to all adults in fossil samples. Our analysis shows that the critical interbirth intervals of australopithecines, early Homo, and Neanderthals are significantly shorter than the observed interbirth intervals of extant great apes. This result suggests that the child care burden of hominid mother was reduced by other group members, including males. In fact, small canine of early hominids may indicate increasing importance of male provisioning, relative to male-male aggression, for males to gain higher reproductive success. To discuss evolutionary causes of these hominid features, we develop a mathematical model to obtain the condition for paternal care to evolve.
Large-game hunting in human hunter-gatherers is a counter-intuitive behaviour. Despite low daily success rates and large proportions of sharing, hunter-gatherers invest in large-game hunting over the less-wasteful strategy of small-game hunting. Where investment in small-game hunting could provide a more reliable and consistent source of direct benefit to one’s offspring, large-game hunting is more ecologically and environmentally dependent, and is a source of benefit to all members of a society, being shared unbiased to the hunter himself. Hypotheses attempt to explain this as reciprocity or altruism, wherein a hunter will share meat obtained as an investment in future returns from receiving beneficiaries of that meat. However, observations show that there is no correlation between current success rates with future obtained goods from others. The Show-Off Hypothesis is an alternate explanation of the evolution of large-game hunting, stating that such hunting behaviours act as a costly signal of one’s quality as a mate. Better hunters are deemed as better mates, and hence obtain more paternities in a society. We investigate this, considering the effect of societal status on the evolution of the strategy of large-game hunting. We develop and present a model that defines measures of hunting success and the survival benefit of care, in order to shed light on previously unquantified effects of status on strategic choice. Within this model, we demonstrate the effect of maladaptive competition on male hunting strategic choice, and provide measures of the trade-off that drives the behaviour of large-game hunting.
Private assessment in indirect reciprocity is a natural assumption because individuals can assess others privately in this situation. However, only few studies have considered private assessment because of analytical difficulty that it present. Here, we develop an analytical method using solitary observation to solve private assessment in indirect reciprocity problem without any approximation. In this study, we formulate a model of solitary observation and calculate the replicator dynamics systems of five leading norms of indirect reciprocity. Our theoretical solution shows that indirect reciprocity in private assessment provides a different result to that in public assessment.
According to the existence proofs of cooperative evolutionarily stable (CES) points in the system, strict norms (stern judging and shunning) have no CES point in private assessment, while they do in public assessment. Image scoring does not change the system regardless of the assessment types because it does not use second-order information. In tolerant norms (simple standing and staying), the CES points move to co-existence of norms and unconditional cooperators. Surprisingly, although there is no central coercive assessment system in private assessment, the average cooperation rate is greater than that in public assessment. This is because private assessment gives unconditional cooperators a role because the unconditional cooperators can raise the cooperation level of the society.
Our results also show the three advantages of the staying norm in private assessment: a higher cooperation rate, easiness of invasion into defectors, and robustness to maintain cooperative evolutionarily stable situations. Our results are applicable to general social dilemmas in relation to private information. Under some dilemmas, norms or assessment rules should be carefully chosen to enable cooperation to evolve.
When one makes his/her decision, he/she often refers to others’ opinions. We describe this situation by a network model with the nodes representing individuals and the links representing references between them. Our question is how people’s reference structure self-organizes, when each individual tries to provide correct answers by referring to more accurate agents.
To answer these questions, we constructed an adaptive network model. In the model, in each iteration round, each agent makes a decision sequentially on a given problem by the majority vote among one’s and his/her neighbors’ opinions. Since each agent makes decision by the majority vote, his/her probability to find a correct answer by oneself, which we call his/her “ability”, is different from his/her actual probability of finding a correct answer by referring to others, which we call his/her “performance”. After the correct answer to a given problem is announced, each agent then rewires his/her links according to the performance of the linked individuals. The rewiring rule is as follows; each individual monitors his/her neighbors’ performance and breaks the link if the neighbor’s performance becomes worth than a preset threshold. We assume that all agents adopt the same threshold. Therefore the value of the preset threshold represents the severity of the society. We also assume that individuals vary in their ability.
In the self-organized reference network, we observed the strong centralization of reference, in which the number of one’s followers increases more than linearly with his/her ability. Counter-intuitively, this tendency was stronger when the rewiring threshold was set lower (when they are more generous to their referents). The mean performance of each agent in the self-organized network was higher compared with random networks or the case of independent decision-making. However, the proportion of agents who give correct answers, which we call “group performance”, fluctuated more in the self-organized network. To sum up, in the self-organized network, though the strong centralization of reference towards high ability agents leads to high performance of each agent on average, it also leads to a heightened risk of a temporal crash in group performance. By analytical calculations, we have confirmed that the performance-monitoring process assumed in our model yields the strong centralization of reference in the self-organized reference network.
The development of an embryo, from a single cell to the breathtaking range of organisms found in nature, is a fascinating process. Moreover, the study of embryonic development promises to shed light on a wide range of developmental defects, inform strategies for the treatment of disease, and provide guidance for tissue regeneration and engineering. Developmental biology has witnessed a renaissance in recent years, due largely to the advent of significant new imaging techniques and molecular biology technologies. To exploit the volumes of multiplex quantitative data that are now being collected, we need to design quantitative hypotheses through mathematical models, and use these models to make quantitative predictions that can be tested and validated experimentally.
The aim of this Developmental Biology Subgroup activity is to bring together researchers in the field who will showcase exemplar interdisciplinary, quantitative studies in developmental biology. This will raise the profile of quantitative approaches to development within the SMB, encourage new members to join the Subgroup and work in the area, and provide a platform for new interactions between SMB members working in the area.
The collective migration of cells during embryogenesis is key to the development of vertebrates, and improper migration can lead to severe developmental diseases and deformities. As a simple model for such cell migration we study the particle based Vicsek model in an open channel geometry where cells continuously enter and leave the domain. This results in two distinct types of motion – one corresponds to desired, well-ordered migration, while the other corresponds to disordered migration. We characterise the different types of collective behaviour and propose a theoretical description to determine the conditions for coherent collective cell migration.
In the course of animal development, the shape of macroscopic tissues emerges from collective cell dynamics. The challenge faced by researchers in the field is to understand the mechanism by which morphogenetic processes of each individual cell (i.e., when, where, and how much individual cells grow, divide, move, and die) collectively lead to the development of a large tissue with its correct shape and size.
Answering this question requires a coarse-grained description and modelling of cell and tissue dynamics at an appropriate length scale. In a previous study, we developed coarse-grained measurement methods for stress and kinematic fields. These methods are now emerging as powerful tools used for exploring the mechanics of epithelial tissues. Given the advancement of experimental measurement methods, one can expect that a theoretical model for kinematics and kinetics in a deforming tissue, which can be compared with the experimentally observable fields, will further advance our understanding of the mechanical control of tissue morphogenesis.
Here, we present a new continuum model of epithelial mechanics. This model incorporates stress and deformation tensors, which can be compared with experimental data. Using this model, we elucidated dynamical behaviour underlying passive relaxation, active contraction-elongation, and tissue shear flow. This study provides an integrated scheme for the understanding of the orchestration of morphogenetic processes in individual cells to achieve epithelial tissue morphogenesis.
In the last two decades, we have witnessed a rapidly increasing number of studies applying methods of mathematical modeling in immunology in general and specifically in the domain of T cell immunology. In part, this was due to the development of several quantitative and powerful techniques for detecting the dynamics of immune responses to pathogens ex-vivo, documenting cell dynamics directly in lymphoid tissues in-vivo, new methods to sequence the full T cell repertoire, and to predict T cell epitopes.
In parallel to the development of novel mathematical and computational models, this has led to better understanding of the results of these experiments, and providing a more quantitative view of the adaptive immune response. The purely theoretical models of the last decades have often been replaced by highly quantitative models that are checked on experimental in-vivo and in-vitro data in human and mice.
Here, we propose a mini-symposium to discuss these new models and their relation to novel experimental tools and the data from the perspective of both modellers and experimenters. Such minisymposium will be of interest beyond the domain of immunology.
During the adaptive immune response, T and B lymphocytes receive and integrate signals from different sources that determine the strength and type of response they follow. Here we asked how reducing stimulation strength through the CD40 receptor could lead to an accelerated division-linked B cell differentiation, as noted in an earlier study [Hawkins et al., Nat Comms 2014]. We observed using flow cytometry that reducing CD40 stimulation strength had two effects on B cell fates in vitro: it slowed down B cell division, and produced a greater proportion of differentiated plasmablasts in each subsequent generation. We studied this system with direct time-lapse imaging to observe division, death, and differentiation to Blimp-1+ plasmablasts, after 3-4 days in culture. We found that strength of stimulation by CD40 affected division times, as predicted, whereas times to die and times to differentiate to plasmablasts were unaffected. We could also show, by fitting mathematical models, that simple fate competition between division and differentiation could explain the accelerated differentiation by slowing division alone. Thus, our data suggest that weakly CD40 stimulated B cells divide slowly, and as a consequence, have more time to undergo differentiation before their next mitotic event. By understanding how different components of the B cell response are controlled and affected by regulatory signals, we aim to build up models of complex signal integration by cells that can predict net immune outcome. Successful models built in this way have many potential applications such as rational immunotherapy design.
CD8+ T (CTL) cells play a pivotal role in protection from viral infection. Better understanding of the heterogeneous phenotypes of T cell subsets evolving during an immune response is timely needed for development of T cell based vaccines that can provide long-term immune protection. Viruses causing chronic infections such as HIV and HCV, trigger a T cell response characterised by functional exhaustion, and accumulation of immune escape variants. We have developed a single cell approach to identify and link functional phenotype, gene expression profile, TCR diversity and T cell avidity with the onset of exhaustion and its relationship with immune escape.
We analysed a cohort of prospectively followed subjects from acute primary HCV infections to disease outcome (clearance and chronicity). Antigen specific (Ag-) CTL were identified via epitope discovery and functional validation via IFNγ-ELISPOT. Index sorting was utilised to link surface phenotyping with gene expression profile. Full length TCRαβ repertoire was identified from scRNAseq data using novel tools. Exhaustion was detected early on in the acute phase of infection and in CTL targeting conserved but also immune escape epitopes. Single cell RNAseq showed a distinct gene expression profile of exhausted T cells compared to the less differentiated subsets. In contrast, Ag-CTL in subjects that cleared HCV, polyfunctional (IFNγ, granzyme, and perforin) responses were found, with a highly diverse TCR repertoire. We discovered a novel T cell clone, characterised by high avidity for viral epitope, elevated IFNγ productions and polyfunctional phenotype.
This systems immunology approach provides new routes to investigate complex T cell dynamics and to identify new mechanisms to exploit in vaccine and immunotherapy.
One of the most controversial and difficult problems regarding the mathematical modelling of biological systems is determination of the appropriateness of the choice of mathematical framework and the implementation of scale in any given model. Concentration scales, spatial scales and temporal scales characteristic of models in biology, more so than any other application, are highly varied and structured even in regards to a single process. In this minisymposium we will be hearing from speakers working on the forefront of foundations in multiscale modelling in biology. It is best to keep an open mind in these talks about the fundamental choice of the methodology which is appropriate to a system and how these methodologies relate to each other. You will be stimulated with ideas regarding discrete versus continuous and stochastic versus deterministic modelling, as well as more general discussion about the nature and treatment of multiscale frameworks in biology.
The three-dimensional structure of eukaryotic genomes is non-random, dynamic, highly regulated, and can be observed to change according to external signals and differentiation state. Disruptions can lead to disease, in which incorrect genomic contacts are responsible for mis-regulation of gene expression. However, the underlying mechanisms that organise the genome are still largely unknown. While most studies focus on the activity of specific proteins which connect two chromatin loci, we seek to understand how the cells makes use of its entire complexity and of physical mechanisms. Using stochastic polymer simulations, which are informed by the analysis of large datasets and verified by quantitative experiments in budding yeast, we discovered a fundamental mechanism: The mobility of the chromatin fibre is not uniform, but heterogeneous, along its length, as a result of the unequal distribution of proteins binding along the genome. This leads to thermodynamically driven self-organisation, which we observe experimentally. It achieves spatial clustering of poised genes and mechanistically contributes to the directed relocalisation of active genes to the nuclear periphery upon heat shock (bioRxiv 106344). I will discuss the implications and present follow-up studies that make use of polymer and particle-based simulations to discover an additional layer of genome organisation.
Many T cell receptors have long, unstructured cytoplasmic tails that contain tyrosine sites. These sites can serve as regulators of receptor activation when phosphorylated or dephosphorylated, while also serving as docking sites for cytosolic enzymes. We coarse-grain the effective interaction between two tails, and develop a mesoscopic particle-based stochastic reaction-diffusion model to study the combined diffusion of individual receptors within the cell membrane, and chemical reactions between proteins bound to receptor tails. The model suggests a switch-like behaviour in the dependence of the fraction of activated receptors on both receptor diffusivity, and on the molecular reach at which two receptor tails can interact. A simplified, analytically solvable model is developed to approximate the more complicated multi-particle system, and used to illustrate how the switch-like behaviour appears.
This minisymposium will bring together researchers to examine up-to-the-minute disease problems that showcase the usefulness and applicability of mathematical modelling to a world far beyond the mathematical community. The audience is the mathematical biologist with an interest in infectious disease. This includes students and researchers, mathematicians interested in seeing applications and biologists who wish to see how mathematics can be used to solve real problems. The minisymposium is interdisciplinary in nature and includes those trained as mathematicians, epidemiologists and immunologists.
Severe diseases arising as sequelae of superficial skin and throat infections with group A Streptococcus (GAS) are important causes of morbidity and mortality worldwide. The observed high level of heterogeneity in GAS prevalence across various temporal and spatial settings suggests potentially complex dynamics of population transmission. One of the most visible indicators of heterogeneity is the different prevalence of skin and throat infections across populations.
Sustained control of GAS infections in settings of poverty has failed, and an effective vaccine may be the most practical long-term strategy to reduce the burden of GAS-related disease. Research into GAS vaccines has been ongoing for many decades, and there are a number of vaccine candidates currently at various stages of development. Those based on the M-protein, the major virulence factor of GAS, include multivalent vaccines targeting the variable N-terminal of the M-protein, and vaccines that contain antigens from the conserved C-repeat region. It is anticipated that these vaccines will provide varying levels of protection across populations depending on the match between circulating and vaccine strains.
We have developed a compartmental model of GAS transmission incorporating skin and throat infection separately and allowing for interaction between the two. To address the many uncertainties that exist around GAS transmission and immunity, we used the Metropolis-Hastings algorithm with informative priors to sample parameter space and simulate our model. We conducted three separate explorations to match epidemiologic attributes observed in population settings with different combinations of skin and throat infection prevalence: 1) throat < 10%, skin 30–60%; and 2) throat 5–20%, skin 5–10%; and 3) throat 5–20%, skin 30–60%; where throat infections may be either symptomatic or asymptomatic.
Early analysis of parameter distributions capable of producing these different combinations of skin and throat infection prevalence has shown interesting patterns for further investigation. While some model parameters are similar across the population settings, others such as the duration of infection, show strong tendencies for differences between populations. We are in the process of considering the plausibility of parameter variation across population settings, which may be driven by clinical, molecular and behavioural characteristics. Once we have ascertained the robustness of our parameter distributions, we will simulate vaccination in the model. To represent the variation in mechanism of vaccine candidates currently under development, we will allow for vaccines to act in different ways in the model, such as having reduced impact on skin infection compared to throat infection. We will use the model to determine the importance of the vaccine mechanism in driving likely impact across population settings.
Equine Infectious Anemia Virus (EIAV), a retrovirus that establishes a persistent lifelong infection in horses and ponies, and which can be transmitted by vectors (biting flies), is endemic in regions with warm climates. With the advent of global warming, research have shown that vector-borne diseases may be on the rise. This study seeks to understand how climate change will affect the EIAV epidemic, especially in endemic regions. We developed two vector-host mathematical epidemic models of EIAV infection that describe the direct transmission of EIAV between wild and domestic horses, and also through vectors: a basic model and a patch model. The models are rigorously analysed mathematically. The effects of variability in the model parameters are also investigated. Results of disease thresholds give conditions necessary for the control of the infection. We hypothesise that climate change to warmer temperatures may increases the time that the virus can thrive in the vector. Model results, which are underway, could be significant for the control of any or all fly-borne pathogens.
Variation of calcium concentration in hepatocytes (liver cells) is known to modulate diverse cellular functions, including bile secretion, glucose and energy metabolism and vesicular trafficking. A major question in the study of calcium signalling in hepatocytes is how these distinct cellular processes are controlled and organised via coordinated spatial and temporal calcium signals.
Downstream cellular responses are controlled via intracellular calcium oscillations but the underlying mechanisms which shape these oscillations have yet to be elucidated. We are interested in determining the effects of various types of calcium feedback mechanisms such as calcium feedback on Phospholipase C (PLC) and the calcium-mediated protein kinase C (PKC) feedback on the hormone receptor have on the whole-cell calcium signals. Recent experimental data suggests that hormone-induced calcium oscillations require positive calcium feedback on PLC to generate inositol trisphosphate oscillations, yielding cross-coupling between calcium and inositol trisphosphate. Furthermore, it appears that there is also a negative feedback pathway, cross-coupling PLC activation to PKC, which serves to terminate calcium spikes.
This talk will discuss recent progress in construction and analysis of a model of calcium oscillations that incorporates the new experimental results about likely feedback mechanisms in hepatocytes.
Calcium signalling is a ubiquitous mechanism that many cell types use to control a plethora of biological processes. In particular, calcium plays a large role in fertility both for egg and sperm cells. In this presentation, we will focus on the acrosome reaction, which is a calcium-modulated vesicle release that must occur for sperm to successfully fertilize the egg. We will discuss the important biological components in this reaction and develop a mathematical model that highlights how multiple pathways interact to drive this unique calcium-driven event.
Bursting is a type of electrical activity seen in many neurons and endocrine cells where episodes of action potential firing are interspersed by silent phases. Pancreatic $\beta$-cells show so-called square-wave bursting when stimulated by glucose, which causes ${\text Ca^{2+}}$ oscillations and pulsatile insulin secretion. $\beta$-cells are electrically coupled within pancreatic islets, and excitation waves are observed to propagate through the $\beta$-cell population. When the islet is exposed to a glucose gradient, so that some cells would be active also when uncoupled while others would be below the activity threshold and thus silent, ${\text Ca^{2+}}$ waves propagate only partly through the islet and stop approximately where the glucose concentration is at the threshold for cellular activity. Simulations of existing mathematical models of coupled $\beta$-cells produce waves that propagate too far into the region of "silent" cells, compared to experiments, unless unrealistic model assumptions are imposed. Here, we investigate why $\beta$-cell models fail to reproduce the experimentally observed wave properties and tend to synchronize the $\beta$-cell population too easily, by using a prototypical polynomial bursting model and slow/fast bifurcation analysis. Our analyses indicate how to modify the model so that the excitation waves stop at the border between "active" and "silent" cells. We verify this property by simulations using such a modified model for a chain, and for a cubic cluster, of coupled $\beta$-cells. Furthermore, we show how our one- and two-parameter bifurcation analyses allow us to predict where the simulated waves stop, for both the original model and the modified version. Our results indicate the geometrical structure that biophysical $\beta$-cell models should have to possess biologically realistic wave and synchronization properties.
The inter-follicular epidermis (IFE) forms the outer-most layer of the skin. Many individual components fundamental to healthy IFE structure are known: proliferation occurs only in a basal layer; above this layer cells differentiate into keratinocytes forming further distinct layers before they are shed from the surface. However, a definitive understanding of how the balance between proliferation, differentiation, and cell shedding is maintained in IFE tissue during homeostasis does not yet exist.
We have developed an agent-based multi-cellular computational model to simulate tissue homeostasis in the skin. Epidermal cells are represented as overlapping spheres, and cell divisions are represented as stochastic time-driven events. Cell movement is determined by adhesive attractive forces and repulsive forces between other cells and the basal membrane. The magnitude of these forces depends upon the types of the interacting bodies, and can vary depending on factors such as cell age and location.
Using this model we have analysed the impact of different cell mechanisms and behaviours on the tissue in order to investigate alternate hypotheses around maintenance of tissue structure. Results of this study provide insights into the critical mechanisms and behaviours, and will guide the future integration of detailed biochemical regulation processes.
Cells are often grown within collagen gels in vitro for applications in tissue engineering. The behaviour of these cells is regulated by their mechanical environment; however, the forces exerted by cells in turn affect the mechanical behaviour of the gel. We aim to better understand the interactions between the cells and the gel using mathematical modelling.
We have developed a multiphase model for this system, incorporating cells and their traction forces alongside chemical effects such as osmosis. To date, we have modelled this problem in one-dimensional Cartesian and spherical coordinates, mimicking experiments performed with spheres of collagen gel. However, these gels are often produced in Petri dishes, resulting in a thin disc. We have therefore transformed our model to examine this type of geometry. In this presentation, we will demonstrate how we can exploit thin-film approximations to reduce the two-dimensional system to a leading-order, one-dimensional model. We will also discuss the equilibrium behaviour of this reduced thin-film system.
The differentiation of mesenchymal stem cells (MSCs) into chondrocytes (native cartilage cells), or chondrogenesis, is a key step in the tissue engineering of articular cartilage. Chondrogenesis is regulated by transforming growth factor-beta (TGF-β), a short-lived cytokine whose effect is prolonged by storage in the extracellular matrix (ECM). Tissue engineering applications require the complete differentiation of an initial population of MSCs, and two common strategies used to achieve this in vitro are (1) co-culture the MSCs with chondrocytes, which constitutively produce TGF; or (2) add exogenous TGF-β. We are motivated by recent experiments where a tissue engineering construct is seeded with cells in spatially separated layers, with a layer of MSCs lying below a layer of chondrocytes, and is then stimulated with exogenous TGF-β from above.
To investigate the efficacy of this seeding strategy we develop a reaction-diffusion model for the interactions between the TGF-β, MSCs and chondrocytes. A feature of this model is that it describes the different forms TGF-β takes throughout its lifecycle, some of which are bound to the ECM (and so not diffusible), and the process by which the TGF-β cytokine is cleaved from the molecular latency complex it is secreted with to become available to cell receptors. Using this model we demonstrate that the concentration of TGF-β required to induce chondrogenesis of the MSCs in the lower layer is not physically realistic for layers of equal depth (as used experimentally). We then suggest improvements to the current experimental protocols, as well as compare the use of the layered cell seeding strategy with a strategy where the initial cell populations are well-mixed (which we have previously investigated with an ordinary differential equation model).
Wolff’s law states that bone morphology evolves according to their external mechanical loading. Following this law, researchers have tried to simulate bone shape formation, especially for trabecular bone, using topology optimization [1]. Less attention has been given to the bone outer shape, composed of cortical bone. However, trabecular bone and cortical bone are both mainly formed by osteoblasts and osteoclasts. Therefore, we hypothesize that the bone outer shape also adapts to the external forces. The aim of this research is to understand the mechanism that generates the bone outer shape by reproducing the latter using topology optimization.
The fish vertebra can be divided into two parts: an inner hourglass-like structure and an outer lateral structure. Based on our observations, it turns out that numerous species present a similar hourglass-like structure but that the lateral structure strongly depends on the fish species. Lateral structures can be classified into two types. The first type exhibits a ridge structure with one or more thick plate-like bones and the second type exhibits a netlike structure in which thin plate-like bones are randomly oriented. These differences seem related to the fish motion, i.e. the swimming type of the fish, and therefore, it is assumed that lateral structures also evolve based on external loading conditions.
In standard topology optimization, the density at each material point is only constrained by the total volume of material. However, the activity of osteoblast and osteoclast is more a local phenomenon. Hence, the standard topology optimization problem is supplemented with a local density penalization to mimic this local phenomenon. To solve the optimization problem, a method based on a time dependent reaction-diffusion equation is employed [2]. The equation is driven by the sensitivity $S$, in which the Lagrangian multiplier $\lambda$ is introduced. Solving this optimization problem gives an optimized structure with respect to the imposed boundary and loading conditions. The penalization term enables to control locally the feature size.
Using this equation, we developed a 3-D mathematical model to generate fish vertebra. Without local density penalization, thick beams appear similarly to the ridge structure. Penalizing locally the density, thinner beams are promoted and they tend to form a netlike structure. Numerical results show that the proposed model can produce both types of lateral structures, i.e. ridge structure and netlike structure, which are similar to the actual fish bone. In the future, it would be interesting to be able to produce various forms of fish vertebra by only adjusting a few parameters of the penalization law.
[1] Jang et al., 2009
[2] Kawamoto et al., 2013
Currently, high-resolution point cloud data can be acquired easily and cost-effectively. For example, a pipeline using Structure from Motion (SfM) and Multi-View Stereo (MVS), which is a promising technique to reconstruct a 3D surface as point cloud data from a series of 2D images taken from different angles, has been implemented in several libraries and software products. In this study, we developed a workflow for measuring and quantifying the canopy growth of soybean cultivars to compare growth speed based on point cloud data acquired from the SfM and MVS pipeline. First, we took multi-view 2D images of soybean plants growing in a field with a digital camera (EOS 60D, Canon, Tokyo). The images of plants in a single plot were taken from ca. 40 different directions. From these multi-veiw 2D images, point cloud data of soybean canopies were reconstructed by using the SfM and MVS pipeline. Second, we segmented the point cloud of plants and their leaves for each plot from the reconstructed data. Finally, we fitted several models to the point cloud data for estimating phenotypic values of plant organs constituting canopy architecture (e.g. leaf area, leaf shape, curvature). For example, we reconstructed 3D surfaces of leaves from the point cloud data with the penalized B-spline surface fitting by regarding a leaf as a 2D closed surface embedded in the Euclidean space ℝ$^3$. When field conditions were not desirable (e.g. wind, change in light conditions) for acquiring corresponding points among multi-view images, however, incomplete point cloud data were reconstructed.
Two ecotypes of the marine snail species Littorina saxatilis have, as a result of strong natural selection, maintained different shell shapes in distinct, adjacent environments. Being able to quantify and compare the underlying growth structure for the two ecotypes should provide a better insight into the reasons behind the variability of shell shapes, and connect this to the genetics and selective pressures in their respective habitats.
Building on the ideas of Raup, we can construct models of snail shells using logarithmic curves containing a set of biologically descriptive growth parameters. This way of modelling illustrates the construction process of the shells, and the growth parameters are therefore expected to relate closely with the underlying genetics.
Previous biological shape research has mainly analysed the shells using variation of a set of landmarks after Procrustes alignment. This method can quantify the shape differences, but does not give an explicit explanation of how or why they occur. Because of the substantial amount of landmark coordinate data available for L. saxatilis, we have developed a geometric method for inferring the growth parameters from these points.
Applying this new analysis to snails collected across contact zones between ecotypes gives the result that growth parameters differ between them and that snails from the hybrid zone show intermediate values, which is consistent with the previous morphometric analysis. However, since the growth parameters are biologically meaningful they should be more informative from a developmental and genetic point of view, and we are going to investigate the performance of both shape characterizations in explaining the genetic variation.
In order to explore the role of initial geometry in wound closure, we developed a new two-dimensional wound healing assay which we refer to as a sticker assay. Stickers are produced and attached to a tissue culture plate before cells are seeded into the centre of a dish. The stickers are then removed to leave a wound in the dish, which is an area with no cells. We performed a number of sticker assays using NIH 3T3 fibroblast cells and various initial wound shapes, including circles, squares, triangles and rectangles. By applying an edge detection method in ImageJ, we tracked the wound interface (and area) at discrete times up to 96h and found that different shaped wound closed at different rates. To explore these experimental results further, we applied two mathematical models, the first being a lattice-based random walk algorithm in two dimensions with nearest-neighbour exclusion, the second being a continuum-limit description based on the Fisher-Kolmogorov equation. By calibrating our models to the experimental data, we find that the parameter estimates for each initial wound shape are similar, suggesting that, while the closure rate is different for each initial wound shape, the underlying mechanisms that drive wound closure are the same. Finally, we revisit our continuum model by allowing for nonlinear degenerate diffusion of the type that occurs in the porous medium equation (charaterised by an index $n$). We explore the hypothesis that circular wound closure is self-similar and circular wound boundaries are unstable for an increasing number of $k$-fold symmetric perturbations as the index $n$ decreases. In particular, we highlight the consequences for wound closure with an initially rectangular wound, where both numerical and experimental results suggest that the wound becomes very long and thin in the limit that it closes.
A leg wound which does not heal because of problems with the veins in the leg is called a venous leg ulcer (VLU). Chronic VLUs are the most common chronic wounds in western countries and their treatment is both costly and time consuming. In this work, we present a mathematical model of the healing of a VLU which incorporates the key biological features of this wound type. We have modelled the role of oxygen, fibroblasts and extra-cellular matrix (ECM) within the wound site. The model consists of a system of nonlinear partial differential equations describing their interactions in space and time coupled with a moving wound outer boundary. The blood vessels that surround the wound supply the oxygen that is needed to support the processes of fibroblasts and ECM to repair the wound tissue. We consider a wound in a simplified one-dimensional geometry and the model equations are discretised in space using the finite difference method and then solved in MATLAB. Numerical results are presented for the oxygen, fibroblasts and ECM distribution within the wound space using parameter values sourced from literature, where possible. This model can be used as a predictive tool in a clinical setting to compare treatments of VLUs including compression bandages.
Scratch assays are routinely used to study the collective spreading of cell populations. In general, the rate at which a population of cells spreads is driven by the combined effects of cell migration and proliferation. To examine the effects of cell migration separately from the effects of cell proliferation, scratch assays are often performed after treating the cells with a drug that inhibits proliferation. Mitomycin-C is a drug that is commonly used to suppress cell proliferation in this context. However, in addition to suppressing cell proliferation, Mitomycin-C also causes cells to change size during the experiment, as each cell in the population approximately doubles in size as a result of treatment. Therefore, to describe a scratch assay that incorporates the effects of cell-to-cell crowding, cell-to-cell adhesion, and dynamic changes in cell size, we present a new stochastic model that incorporates these mechanisms. We then employ this stochastic model to quantify relative contributions of crowding effects and random motility into the collective cell migration in the scratch assay experiments containing malignant PC-3 prostate cancer cells.
Scratch assays are standard in vitro experimental methods for studying cell migration. In these experiments, a scratch is made on a cell monolayer and imaging of the recolonisation of the scratched region is performed to quantify cell migration rates. Typically, scratch assays are modelled by continuum reaction diffusion equations depicting cell migration by diffusion and carrying capacity-limited proliferation by a logistic source term. In a recent work [1], the authors observed that on a short time, there is a disturbance phase where proliferation is not logistic, and this is followed by a growth phase where proliferation appears to be logistic.
In this talk, I will introduce an age-structured population model that aims to explain the two phases of proliferation in scratch assays. The cell population is modelled by a McKendrick-von Foerster partial differential equation. The conditions under which the model captures this two-phase behaviour are presented.
[1] Jin, W. et al. Bull Math Biol (2017)
This minisymposium will bring together established and up-and-coming researchers to share ideas on the recent advances on the use of mathematical approaches to gain insight into the transmission dynamics and control of emerging and re-emerging infectious diseases in populations.
An essential feature of the minisymposium is the emphasis on the use of state-of-the-art techniques, theories and new applications associated with the use of dynamical systems in modeling the spread of diseases in populations. Current pressing modelling and mathematical challenges will also be discussed.
Mathematical models for the Ebola Virus Disease (EVD) were, until recently, mostly based on the fast/direct transmission route, which involves contact with blood or body fluid and objects that have been contaminated by body fluid. The fact that in almost all outbreaks of the EVD in Africa, the index case became infected through contact with infected animals, such as fruit bats and primates, makes the environment a no negligible channel for slow/indirect transmission of the disease.
In this talk, we incorporate both direct and indirect transmission routes in the setting of SIR-type models that are developed gradually. Firstly, we add two compartments (for dead individuals and for the environment) in order to capture infection through the manipulation of deceased individuals before burial as well as contaminated environment. Secondly, we enrich the first model with self-protection measures reflected by the addition of classes of vaccinated and trained individuals.
For the two models, we prove that the disease-free equilibrium is globally asymptotically stable (GAS) whenever the basic reproduction number is less than or equal to unity, and unstable when this threshold number is greater than 1. In the latter case, the existence of at least one endemic equilibrium (EE) (for the second model) and of a unique EE (for the first model), which are locally asymptotically stable (LAS) is shown, together with the fact that the unique EE is GAS in the absence of shedding and manipulation of deceased human individuals before burial. At the endemic level, it is shown that the number of infectious individuals is much smaller than that obtained in the absence of any intervention.
In a final step, we construct nonstandard finite difference schemes that are dynamically consistent with respect to the above features of the continuous models.
Using a data base on Dengue incidence available for four different states of the country we develop mathematical models that describe the recurrent dynamics of cases at different geographical/regional levels. We present preliminary results.
In recent years, there has been an explosion of multiscale cancer patient data including genomic, metabolomic, proteomic, and histologic (with novel stainings), at the tissue scale and smaller, along with a host of advanced MRI and other clinical imaging sequences at the whole organ scale. While bountiful, the multiscale and spatial sparsity represented with these data types for any individual patient leaves the final interpretation difficult with traditional statistical methods. Mechanistic models provide a framework for exploring how these processes interlink and contribute to tumour growth and evolution. Incorporating this multitude of data into such models remains a challenging and open problem, however, promising results are emerging. The talks in this minisymposium will highlight novel models and approaches to using these new data types in exploring the mechanisms underlying cancer growth and response to therapy.
Glioblastoma multiforme (GBM) is a rare brain cancer with a median survival of only around 15 months. Intratumoural heterogeneity and extensive infiltration into the brain tissue contribute to poor prognosis and probable recurrence. Predicting the timing of post-treatment recurrence is often limited if using only MRI imaging measurements, as a diverse range of treatment outcomes can result from similar pre-treatment growth dynamics. To better understand the diversity in recurrence, we analyze multiscale data from an ex vivo rat model during progression to inform an agent-based computational model to quantify how heterogeneity in proliferation and invasion affects both the bulk and individual scale tumour metrics following treatment.
The rat experiment is a model of platelet-derived growth factor (PDGF) driven GBM. The production of PDGF by cells allows both autocrine and paracrine stimulated proliferation and migration to drive growth. PDGF in the tissue environment also stimulates and recruits normal progenitor cells. Throughout this experiment, bulk tumour size was recorded via MRI, and single cell tracks were analyzed to quantify proliferation and migration spatially. Guided by this data, we built an agent-based model within the gray-white brain tissue context. We fit the model to the growth data alone using a hybrid genetic algorithm approach, and then fit individual cell proliferation and migration distributions.
We find parameter sets that define a suite of fits, which describe similar bulk scale growth dynamics with widely different underlying phenotypes that tend to be more proliferative, more migratory, or more driven by growth-factor dynamics. A more limited set of parameters fits the individual scale dynamics, with heterogeneity predicted to be both an intrinsic feature and driven by the environment. Application of an anti-proliferative treatment on these heterogeneous tumours reduces tumour size, but leads to recurrence of a less proliferative tumour. An anti-migratory treatment is seen to have little effect on growth dynamics, but selects for more proliferative phenotypes. The changes that occur at the imaging scale down to the single cell level lead to a better understanding how the evolutionary impact of treatment on GBM can inform treatment regimens.
Introduction: Glioblastoma (GBM) is a very aggressive primary brain cancer, noted for its diffuse infiltration into surrounding normal–appearing brain. This particular nature makes GBM notoriously difficult to treat, as these diffusely invading cells cannot be resected surgically, are difficult to target with radiation therapy, and thus must be targeted with chemotherapy. However, this too presents a challenge, as these invading GBM cells reside beyond the dense tumour regions where angiogenesis causes disruption of the blood brain barrier (BBB) and allows drugs to more readily enter the central portion of the tumour. Thus, it is critical to determine predictors of drug distribution in individual patients’ tumours and surrounding brain tissue to ensure the invading GBM cells are exposed to the therapy.
Objective: Determine predictors of drug distribution and effect from non-invasive imaging using minimal mathematical approaches.
Methods: Following experiments treating murine orthotopic patient-derived xenografts (PDXs) of GBM with various anti-tumour therapies, we compiled data from both the xenografts and the original patients from which the PDX lines were derived. This data included bioluminescence imaging (BLI), lamin–stained histological sections, and magnetic resonance imaging (MRI) from the PDXs, as well as patient MRIs and clinical data. Using the time series BLI data from PDXs, we developed and parameterized minimal differential equation models of PDX tumour growth for individual PDX lines, adjusted for lead time bias using a nonlinear mixed effects approach. This gave us an overall growth rate for each of the PDX lines across multiple subjects. Next, we compared these growth and invasion characteristics with therapeutic response in patients and PDX subjects to agents with different CSF to plasma ratios, indicating the degree of blood brain barrier penetrability of these drugs.
Results: Individual PDX lines have different growth kinetics, and recapitulate the kinetics observed in the original patients from which the lines are derived. Further, these growth kinetics appear to be correlated with differential drug response, with more diffusely infiltrating tumours responding better to drugs with higher CSF to plasma ratios.
Conclusion: While further work is needed to verify our results across more PDX lines, our results suggest that noninvasive imaging-based characterization of tumour invasiveness may be able to aid in matching patients to the best therapy for their individual tumours.
In this research, we proposed an infectious disease spread model of influenza which was individually based and stochastic. We first described the method to build a realistic model, or to estimate realistic parameter values, based on observed data. An appropriate model has to be individually based and stochastic if we take advantage of small unit data (in this case, every day class unit data). We also analysed the effect of different strategies of school closure, e.g., the criteria for implementing school closure and the duration of continuous school closure.
The data we used were the daily reported data of the numbers of cases of swine flu and the dates of school/class closure in each class of schools collected from September 2009 to March 2010 in a city in Japan. Because the data were collected in each class, transmission rates could be estimated separately for within class, within-school and inter-school rates. Totally 21,253 cases were reported out of 51,871 students in 134 schools (elementary schools, junior high schools, high schools and kindergartens). To construct accurate model, infected numbers must be estimated from case report data taking latent period into account, which is deterministically impossible. To resolve this difficulty, we stochastically estimated and obtained data of infected numbers by ‘Monte Carlo back calculation’ and used in the analysis. The analysis was mainly carried out using data during September up to December to avoid the influence of vacations. The influence of humidity was also analysed using daily meteorological data of the area where data were collected. The analysis of total population using deterministic model was also carried out to obtain the overview of the flu spread.
The transmission rates of different levels were estimated by maximum likelihood method using detailed data reconstructed by Monte Carlo back calculation method. The rate within class was much larger (approximately 15 times larger) than within school rate. The transition rate between schools was much smaller than between classes of the same school. Stochastic variations of estimated parameters were also analysed. Optimal strategy of school closure was analysed in relation to the characteristics of flu.
By the analysis based on stochastic infectious disease spread model, we could have presented possible variations of total number of cases and the size and timing of the epidemic peak. ‘Monte Carlo back calculation’ approach was concluded to be useful in the analysis of effective preventive strategy against influenza based on individually based stochastic model.
In this talk, we will present a mathematical model of 2009 A/H1N1 influenza considering age structure in the Republic of Korea and suggest vaccination strategies for mitigating the epidemics. There were 750,000 confirmed cases of 2009 A/H1N1 influenza from May 2009 and August 2010. Because influenza viruses are spread through close contact, contact pattern plays an important role for the disease spread. We developed an age-dependent SEIAR model. The total population is divided into five subgroups and age-specific transmission rates are estimated based on laboratory confirmed data from Korea Centers for Disease Control and Prevention (KCDC). Vaccination policies is important to minimize the number of infected individuals especially when the vaccine supply is limited. Using mathematical model, we could assess the impact of the age-dependent vaccination priority on transmission dynamics.
Networks representing the spread of infectious diseases in populations have been widely studied. Here, we formulate an SEIR model using an edge-based approach on a static random network with arbitrary degree distribution. The corresponding basic reproduction number and final epidemic size are computed. The SEIR model is used to investigate the stochasticity of the SEIR dynamics. Assuming exposed, infection and recovery each happen at constant rates, stochastic simulations of the SEIR dynamics are performed applying continuous-time Gillespie's algorithm given a Poisson or a power law with exponential cut-off degree distributions. The resulting simulations match well with the numerical predictions of the SEIR model given the initial conditions. Final epidemic size remains unchanged when the initial infecteds are varied. On the other hand, varying the disease parameters of the SEIR model affects the time when the epidemic accelerates, the peak of the epidemic, and the final epidemic size. These results capture scenarios of an epidemic in a network implying control strategies in the disease transmissions.
Small-scale fluid mechanics play an important role in virtually all living systems. Biological processes taking place at very low Reynolds numbers pertain to individual cells (swimming of bacteria or sperm cells) as well as to cell collectives and even groups of larger, multicellular organisms. Classical theoretical explorations and efforts to model fluid flow at low Reynolds number have long been limited to single cells swimming through a Newtonian fluid. Biological problems, however, often involve complex environments and/or systems containing more than one cell or organism, and traditional models have yielded only approximate answers to these problems. Recently, the development of computational tools has spurred advancement in the study of low Reynolds number fluid mechanics in biological systems.
In this minisymposium, we focus on both locomotion in complex materials (such as soil or viscoelastic fluids) and the emergent collective dynamics observed in organisms at different scales by expanding on existing fluid models, developing better-performing algorithms, and integrating behavioural and fluid mechanical models.
The swimming motion of microorganisms such as sperm and cilia can be modelled by several methods, all of which entail solving equations of fluid-structure interaction. Among them, the Method of Regularized Stokeslets (MRS) and the Rotne-Prager-Yamakawa tensor have the advantage of not requiring a 3D Eulerian grid and using the fundamental solutions to the underlying equations instead. However, the computations required by both methods entail the use of dense matrices, and they tend to be large and very costly to work with for practical models in which the number of micro-swimmers is large.
The 'data-sparse' structure of these matrices enables the development of fast algorithms. To compute the matrix-vector products efficiently, we extend the Kernel-Independent Fast Multipole Method (KIFMM) to the kernels associated with the MRS. To solve linear systems with the same matrices efficiently, we consider both a data-sparser preconditioner and a block-diagonal preconditioner; to expedite the application of the preconditioners, we employ a number of techniques such as Krylov subspace recycling. We apply the proposed algorithms to study the dynamics of a large group of sperm and the flow field induced by a carpet of cilia.
Xeniid corals, a family of soft corals (Alcyonacea), include species displaying a unique pulsing behaviour. Within a colony, each individual polyp pulses by actively contracting and passively expanding eight tentacles, increasing the local mixing and enhancing nutrient and gas exchange. Using the immersed boundary method with finite elements (IBFE), we constructed a 3D model of a pulsing polyp. The motion of the polyp tentacles is based on actual motion data tracked from videos of real polyps. We find that individual polyps pull water in radially, mix it between their tentacles, and expel the fluid volume in an upward jet. After validating this 3D IBFE model against experimentally measured flow fields, we are now using the model to numerically simulate small groups of polyps and to quantify the effects of collective pulsing behaviour on the local fluid dynamics. Here, we simulate pairs of polyps and vary the pulsing patterns (in-phase and different degrees of out-of-phase) and distance between polyps to better understand how differences in collective pulsing behaviour affect local flow and mixing.
Population dynamics are a traditional focus of mathematical biology, dating from the pioneering work of Volterra and Lotka. Importantly, the rates of change of species' abundances have broader consequences, affecting a diverse range of processes that includes evolutionary, behavioural, ecosystem and human dynamics. In this minisymposium, we highlight biomathematical models that couple population dynamics with these broader factors. Not only do these coupled feedback models help us understand the behaviour of broader systems, they are also essential to an accurate understanding of the population dynamics themselves, and how to manage them. Developments in this modelling space can help us understand essential characteristics of ecological systems, and also generate mathematically interesting complex system behaviours, which could stimulate new theoretical hypotheses.
In this session, three speakers will demonstrate population dynamic models that are integrated with intrinsic processes of organisms, i.e., evolution, behaviour and inducible responses. Two speakers will present about joint dynamics of populations and management policy, in the context of biodiversity conservation and agricultural damage control.
Conservation management decisions are often implemented at the scale of human communities, rather than the scale of the most relevant ecological dynamics. Research frequently points out the loss in efficiency that results from such "scale mismatches". However, the scale of management is influenced by socio-economic constraints on management actors - not all people want to cooperate with each other. While it is clear that objectives can be better achieved if management and ecological scales are aligned, it is not clear whether such benefits justify the costs of alignment, or how alignment can be achieved in multi-actor contexts.
On Manus Island in Papua New Guinea, small customary tenure areas define the scale of fisheries management decisions, but their fish stocks are connected by pelagic larval dispersal. We quantify the extent of this scale mismatch for the serranid Plectropomus maculatus by empirically estimating larval dispersal distances, and integrating this data into a bioeconomic model of the multi-actor fishery, and identifying the cooperative Nash equilibrium. Larval dispersal allows individual communities to externalise the costs of overharvesting, and export the benefits of cooperative behaviour, to non-cooperating groups. Despite these mismatched scales, the communities on southern Manus have recently created a tribal network for management and negotiation, emphasising the importance of social capital in avoiding suboptimal outcomes or the need for top-down governance.
Species interactions are important to determine structure and stability of ecological communities. In particular, a variety of indirect effects appear between species which do not interact directly. However, indirect effects do also appear between species that directly interact. For example, plant species sharing a common herbivore may engage in both interference competition directly with other plant species and apparent competition indirectly through common herbivores. Then, if a plant species augments herbivores, then grazing pressure on another plant species increases and decreases in biomass of the latter plant may weaken competitive pressure on the former plant species. In this way, there may appear many feedback loops in a complex ecological community where species affect each other through many different paths connecting them. Thus, we should understand interactive consequences of a direct and many indirect effects to understand structure and dynamics of ecological communities. However, studies on indirect effects tend to focus on changes in density of a target species and neglect changes in diversity of the whole community.
In this presentation, we consider a very simple dynamical system model that describes a community composed of predator, herbivore and many plant species. In this model, since predators reduce herbivores and weaken grazing pressure on plants, there appear a positive trophic cascade from the predator to plants. Then, increased plant biomass may enhance both direct interference and indirect apparent competition between plants. Therefore, it is interesting to investigate how the presence of predators affects diversity of the plant community. We analytically solve steady state solutions of the model and investigate effects of predation, grazing, and interference competition on the number of plant species and total plant biomass. Based on the results, we discuss the importance of interplay between direct and indirect effects for structure of ecological communities.
Asthma is a chronic lung disease of reversible airway constriction. Imaging experiments show that during an asthma attack, the asthmatic's lung exhibits what is known as clustered ventilation defects, which is the hallmark trait of asthma. This phenomenon is when there are some regions of the lung where the airways are closed, and some regions where they are open. These clusters vary from event to event, even in the same patient, and thus it is believed that the causes are dynamic rather than structural. We want to understand the dynamic mechanisms behind the physiological implications. Different mechanisms leading to spatial clustering have been suggested; we are interested in assessing the physiological viability of some of these models, particularly with how each model exhibits this spatial clustering. During this talk, we will go over these models and how they have shown themselves to be different in the dynamics analysis.
Tick-borne pathogens are transmitted when ticks take blood meals from vertebrate hosts. Ticks need to take blood meals to progress through immature life-stages and reach adulthood. For the most important zoonotic pathogens, including Borrelia burgdorferi (the causative agent of Lyme disease), two immature life-stages of the tick vector, termed larvae and nymphs, maintain the pathogens. Key features of tick feeding behaviour, and therefore of tick-host contact patterns, include the aggregation of ticks on hosts (whereby most ticks of a given life-stage feed on only a small minority of the hosts) and the co-aggregation of larval and nymphal ticks on the same minority of hosts.
A mechanistic network model is presented for tick-borne pathogen transmission that explicitly accounts for larval and nymphal tick co-aggregation and coincident coaggregation, also known as co-feeding. Co-feeding of nymphs and larvae allows transmission from an infected nymph to susceptible larvae feeding in close proximity and at the same time, but without the involvement of a systemic infection in the vertebrate host. By relating the next generation matrix epidemic threshold parameter $R_{0}$ to the in- and out-degrees of vertebrate host nodes in the mechanistic network model, a simple analytic expression for $R_{0}$ that accounts for the co-aggregation and coincident coaggregation of ticks is derived. Simulations of Lyme disease transmission on finite realizations of tick-mouse contact networks are used to visualize the relationship between $R_{0}$ and the extent of tick co-aggregation.
The derived analytic equation explicitly describes the relationship between $R_{0}$ and the strength of dependence between counts of larvae and counts of nymphs on vertebrate hosts. Tick co-aggregation always leads to greater values for $R_{0}$, whereas higher levels of tick aggregation only increases the value of $R_{0}$ when larvae and nymphs also co-aggregate. Aggregation and co-aggregation have a synergistic eﬀect on $R_{0}$ such that their combined eﬀect is greater than the sum of their individual eﬀects. Co-aggregation has the greatest effect on $R_{0}$ when the mean larval burden of hosts is high and also has a larger relative eﬀect on the magnitude of $R_{0}$ for pathogens sustained by co-feeding transmission (e.g. TBE virus in Europe) compared with those predominantly spread by systemic infection of the vertebrate host (e.g. Lyme disease).
Co-aggregation increases $R_{0}$, particularly in geographic regions and seasons where larval burden is high and for pathogens that are mainly transmitted during co-feeding. For all tick-borne pathogens though, the eﬀect of co-aggregation can be to lift $R_{0}$ above the threshold value of 1 and so lead to persistence.
Carnitine is a fundamental compound for humans. Mainly introduced via nutrition, but also synthetized in the body, this metabolite allows the transport of long-chain fatty acids in the mitochondria, where they can undergo β -oxidation. Lately, the attention on carnitine and acyl-carnitine, carnitine derivatives, has increased since it has been proposed that insulin resistance may be linked to incomplete fatty acid β-oxidation and the subsequent increase in acyl-carnitine species in different tissues.
Mathematical modelling can give us a quantitative understanding of specific rates in the carnitine metabolism, e.g. synthesis of carnitine derivatives . Therefore we contribute to the research on carnitine by building a compartmental model for carnitine and its derivatives. This model integrates knowledge about carnitine uptake and synthesis, carnitine distribution and expulsion, and carnitine role in the transport of long-chain fatty acids, therefore covering all the carnitine occurrences in the body. Additionally, the simultaneous comprehension of carnitine and its derivatives enables the model to highlight the role of carnitine concentration as a bottleneck in the transport of the long-chain fatty acids into the mitochondria.
Our model utilizes metabolome datasets publically available in the literature for validation.
Maintaining physiological levels of oxygen (O$_2$) and carbon dioxide (CO$_2$) in the blood is crucial for survival and is achieved by sophisticated neural control mechanisms affecting both the breathing pattern and heart rate. Neural activity, originated in the brainstem, drives the respiratory muscles, providing air flow into and out of the lungs where gas exchange takes place and also affects heart rate and blood flow. Chemoreceptors that sense the levels of O$_2$ and CO$_2$ in the blood, mechanical stretch receptors within the lungs and blood pressure sensors provide feedback signals to the brainstem networks which then regulate the breathing pattern and heart rate appropriately. Understanding how the neural system responds to all the feedback signals it receives is still lacking. I will outline the challenges we face from a modelling perspective and the different approaches we take in our attempt to resolve them.
With a five dimensional system of ordinary differential equations based on the SIR model, we consider the dynamics of epidemics in a community which consists of residents and visitors/tourists over a short period of time. The total population size of the community is taken to be constant, ignoring its change due to any birth and death in the period under consideration. Also, the resident and visitor populations are respectively constant. We assume that every immigrating visitor is susceptible and is likely to be infected during their stay in the community. Furthermore, infected visitors can carry on their activities normally during their stay in the community thus still appearing like susceptible visitors. From the analysis of the model, we obtained a threshold expected value of duration per visitor which determines whether an epidemic persists in the community or not. More so, the basic reproduction numbers with respect to the entire community as well as the resident and visitor populations are obtained and they can help us to formulate important public health policies in order to contain epidemics in the community.
We study synaptically coupled neuronal networks to identify the role of coupling delays in network's synchronized behaviours. We consider a network of excitable, relaxation oscillator neurons where two distinct populations, one excitatory and one inhibitory, are coupled and interact with each other. The excitatory population is uncoupled, while the inhibitory population is tightly coupled. A geometric singular perturbation analysis yields existence and stability conditions for synchronization states under different firing patterns between the two populations, along with formulas for the periods of such synchronous solutions. Our results demonstrate that the presence of coupling delays in the network promotes synchronization. Numerical simulations are conducted to supplement and validate analytical results. We show the results carry over to a model for spindle sleep rhythms in thalamocortical networks, one of the biological systems which motivated our study. The analysis helps to explain how coupling delays in either excitatory or inhibitory synapses contribute to producing synchronized rhythms.
One of the fundamental problems in biology would be the method by which organisms can regulate a distribution in response to global information. For example, a cluster of organisms can regulate the proportion of individuals that perform various roles or modes as if each individual is aware of the overall situation without a leader. Slow or rapid self-organized pattern formations with time evolution on animal skins might also be one of the processes of ratio regulation of pigment cells. Thus, in various species, a specific ratio exists at multiple levels, from the process of cell differentiation in multicellular organisms to the situation of social dilemma in a group of human beings. In this study, we consider individual that means a particle who has a certain size and internal mode, then propose a common basis for regulating collective behaviour that is realized by a series of local contacts between individuals. The most essential behaviour of individuals is to change their internal mode by sharing information when in contact with others. Our numerical simulations and analysis for regulating the proportion in two kinds of modes show that asymmetric properties in local contacts are essential for adaptive regulation in response to global information. In this poster, some actual examples of proportion regulation are described and applicability of our theory would be verified. We will discuss that this simple mechanism of this theory could indicate that well-organized groups in nature can be regulated through local contacts only.
M. Iwamoto & D. Ueyama, “Basis of self-organized proportion regulation resulting from local contacts”, Journal of Theoretical Biology, 440 (2018) 112–120.
Circadian clocks control 24-hour rhythms in the human body, including glucose and insulin dynamics. Glucose tolerance changes depending on the time of day and misalignment between environmental stimuli such as light and food and internal circadian clocks, e.g., as observed in shiftwork, is associated with metabolic disturbances and diabetes. The mechanisms of these are poorly understood. A number of mathematical models of circadian clocks and those of glucose-insulin dynamics have been developed. However, no model accounts for the interaction between the two systems. Here we examine and compare key existing models of glucose-insulin regulatory system with the aim to select the most appropriate model for incorporation of circadian rhythms.
A large number of mathematical models designed to study glucose-insulin dynamics have been developed over the past few decades and account for such key features as (i) ultradian oscillations, (ii) response to glucose infusion, and (iii) glucose homeostasis. The seminal model able to reproduce these oscillations was developed by Sturis et al. (1991) utilising a system of 6 ODEs. Since then, a number of other models have been adapted from this model, aiming to either reformulate the Sturis model—often as a system of DDEs—or to extend the model to allow other applications, for example in the Intravenous Glucose Tolerance Test. However, comparisons of the models with physiological experimental studies have been relatively limited in scope.
In this work, we compare the qualitative and quantitative predictions of the Sturis et al. (1991), Tolić et al. (2000), and Li et al. (2006) models and two other models with experimental results, under constant glucose infusion. We show the different parameter ranges in which the models display sustained ultradian oscillations and how this corresponds to the expected behaviour in published physiological studies. We also test the mean values, amplitudes, and frequencies of the simulated glucose and insulin concentrations against experimental studies. This sophisticated comparison of the models allows us to select the most appropriate model to incorporate circadian rhythms.
The walls of our blood vessels are under a constant mechanical load, introduced by the blood pulsating through our vessels. Recent years have seen a substantial increase in the understanding of the interplay between the dynamics of blood flow (haemodynamics) and vessel (vascular) morphology. Alterations in the homeostatic distribution of mechanical forces exerted by blood on the vessel wall is correlated with various cardiovascular diseases, including atherosclerosis and aneurysm development and rupture. Furthermore, alterations in the haemodynamic distribution have a notable effect on vascular development and remodelling. Experimental analysis of the effects of blood flow on the local vasculature is invasive and extremely difficult to measure in real time. As such, the development of computational models to analyse the relationship between the local haemodynamics and the surrounding vasculature is invaluable.
In this project we will develop a computational model to study how vessels deform when subject to internally and externally applied forces. Using the open source multicellular modelling software package, Chaste, we will employ a discretised approach to simulate long time scale vascular wall deformation. A Voronoi tessellation-based model will be used to model vessel walls. This discretised model will enable us to examine how blood vessels deform when subject to internal and external forces, thus laying the foundations for future research examining vascular deformation due to haemodynamic pressure.
Macroecological patterns, such as the species area relationship (SAR), relative species abundance (RSA), and endemic area relationship (EAR) provide us useful information of ecosystem structures that can be applied to, for example, ecosystem conservations. To understand these patterns across spatial scales is one of the central challenges recent years. To tackle the challenge, we develop a spatially explicit geometric model where a set of species distribution ranges and individual distributions therein are described by spatial point processes. Essentially by calculating the number of individuals or distribution ranges given a sampling region, several well-known macroecological patterns are recovered including the tri-phasic SAR on a log-log plot with its asymptotic slope, and various RSAs such as Fisher’s logseries and the lognormal distribution. In a simple manner, the EAR is also calculated in our framework. These multiple macroecological patterns can be derived with a single set of biological parameters, and therefore it may provide us a convenient way to discuss the multiple patterns in a consistent manner.
The Human African Trypanosomiasis (HAT) parasite, which causes African Sleeping Sickness, is transmitted by the tsetse fly as a vector. It has several possible hosts, including humans and domestic animals. Because domestic animals can be a host for the parasite, it has long been assumed that keeping domestic animals near human populations increases the spread of the disease. However, several parameters found in the literature, including the shorter lifespan of the male vector and the female vector's preference for domestic animals, made us question this assumption.
We have developed a differential equation compartmental model to examine whether increasing the domestic animal population can be used to deflect the infection from humans and reduce its impact. This 9-dimensional system of nonlinear ordinary differential equations includes tsetse flies in their various stages of maturity, which is more than most previous models have done. We used the Next Generation Matrix method to obtain an expression for $R_0$, the basic reproduction number, based on the other parameters in the model.
Our study indicates that strategies that were not previously considered, such as vaccinating domestic animals, may reduce the impact of the disease on humans even better than vaccinating humans.
We constructed an anatomically-accurate three-dimensional model of salivary fluid secretion from a cluster of parotid gland acinar cells.
Parotid acinar cells are responsible for the secretion of saliva. Olfactory and gustatory stimuli provoke the release of agonists that bind to the basal membranes of the acinar cells. This triggers a cascade of events that results in the production of IP3, which, in turn, releases calcium ions from intracellular compartments. Upon release of calcium, chloride channels are activated, leading to chloride ions flowing out of the cell, and water following by osmosis.
We have shown previously that the complex spatial structure of the acinus and the spatial heterogeneities within each acinar cell imposes severe constraints on the possible mechanisms that could explain the saliva secretion process.
We will discuss how could we simplify the models that are able to replicate the physiological observations under such constraints.
Rice field art is a large-scale art form in which people design rice fields using various kinds of ornamental rice plants with different leaf colours. Leaf colour-related genes play an important role in the study of chlorophyll biosynthesis, chloroplast structure and function, and anthocyanin biosynthesis. We performed whole-genome resequencing and transcriptomic analysis of regulatory patterns and genetic diversity among different rice cultivars to discover new genetic mechanisms that promote enhanced levels of various leaf colours. We re-sequenced the genomes of 10 rice leaf-colour varieties to an average depth of $40\times$ depth and >95% coverage and performed 30 RNA-seq experiments using the 10 rice varieties sampled at three developmental stages. The sequencing results yielded a total of $1,814\times10^6$ reads and identified an average of 713,114 SNPs per rice variety. Based on our analysis of the DNA variation and gene expression, we selected 47 candidate genes. We used an integrated analysis of the whole-genome resequencing data and the RNA-seq data to divide the candidate genes into two groups: genes related to macronutrient (i.e., magnesium and sulfur) transport and genes related to flavonoid pathways including anthocyanidin biosynthesis. We verified the candidate genes with quantitative RT-PCR using transgenic T-DNA insertion mutants. Our study demonstrates the potential of integrated screening methods combined with genetic-variation and transcriptomic data to isolate genes involved in complex biosynthetic networks and pathways.
Direct-Acting antivirals (DAAs) target intracellular viral replication and have realized high effectiveness for hepatitis C virus (HCV) patients. Now, many kinds of DAA drugs have discovered and the HCV treatment is improving day by day. A popular strategy of HCV treatment is combination of double or triple DAA drugs with different action mechanisms. Mathematical model is used to analyze the dynamics of virus infection and effectiveness of antiviral drugs. For DAA treatment, multiscale model formulated by partial differential equations (PDE) describing both intercellular virus infection and intracellular virus RNA replication is used. However, in general, numerical computation of PDE often converge poorly and is time consuming. So it is somewhat hard to estimate parameters, and we should transform the PDE to easier form in advance of analyzing the model. In our study, we transformed the PDE model to ODE form without any assumption. Thus our ODE model and the previous PDE model are mathematically identical. By this transformation, the time required for numerical simulation is reduced greatly and we can apply some basic analysis methods for dynamical systems of ODEs to the multiscale model.
Life-histories in organisms faced with various selection pressures such as heterogeneity, intra- and inter-specific competitions, and variable environment. Ecologists have thought that evolution maximizes different objects from each pressure. For example, in the absence of the competitions, conventional Darwinism asserts that the adaptive life history maximizes the population growth. On the other hand, the presence of competitions is thought to evolves un-invadable life strategy from any other strategies such as maximizing the carrying capacity that $r/K$-selection theory asserts. Difference of each habitat, thus, yields different objective function and diversity of life-histories from precocity to altricity.
However, there is ambiguity involved in necessity for different objective function in each habitat. It is difficult to define the threshold changing from population growth rate to carrying capacity, as objective function. This difficulty fades $r/K$-selection theory and makes ambiguity in categrization of life history evolution. On the other hand, several indices which the fittest species have are suggested (ESS etc.), and these have a common point that the fittest species or genetic group occupy their habitat. If an objective function addresses the common point without presence and absence of competitions, and if it evolves all traits that are precocity, altricity, and so on in life history, evolution of life-history can be systematized by maximizing (or minimizing) a unique objective function.
In this presentation, I would like to propose that the adjoint eigenfunction providing reproductive value in a structure population model is a candidate of the unique objective function.
Chemical reactions are traditionally modelled using deterministic methods such as ordinary differential equations. In biology, we also consider the effect of diffusion upon a reactive system. For small biological systems, it is necessary to use stochastic models. In this talk, we choose to use an agent-based model commonly known as Smoluchowski kinetics.
In the Smoluchowski model, two molecules diffuse in space, and react when they are within a critical separation distance (Smoluchowski radius). This radius is at a rate proportional to the reaction rate constant. Over the last 100 years, the model has been extended to incorporate reversible reactions by introducing a second separation distance.
A common problem with Smoluchowski kinetics is when modelling rapid enzyme reactions, where the assumptions used to model Smoluchowski kinetics break down. Using differential equations, these rapid interactions can be easily analysed using a pseudo steady state approximation which effectively results in non-linear reactions. Other models have utilised multiple timescales to model what is called fast-slow reaction kinetics.
In this talk, we theorise that a way to bypass the fast-slow reaction kinetics problem is to extend the Smoluchowski model to reversible trimolecular reactions. We develop these kinetics and apply them to the canonical Wnt signalling pathway. We model this pathway using traditional Smoluchowski kinetics, and show how the aforementioned issues come into effect. We then compare this to our trimolecular model and show how the issues disappear.
Reproducibility of results is a key part of the scientific method; in general, scientific communication aims to describe a result in enough detail that readers and reviewers are able to contextualise the result within their own knowledge, and to reproduce it themselves given the appropriate skills and resources. In the field of computational biology, readers and reviewers face special challenges to contextualise and reproduce results, because of the size of datasets analysed, and the increasing complexity of the computational methods used. To meet these challenges, we present an approach that uses software engineering tools to produce complete ‘reference computation environments’, containing all software and configuration necessary to reproduce a computational result. Existing tools to reproduce results are limited to a single technology or programming language; our approach is the first to use a single specification to produce the same environment on a desktop computer, in a cloud computing service, or in a managed computing facility. Using these tools, authors can ensure their published results are easily reproducible by readers and reviewers, and are robust to future changes in software and hardware.
The drive for reproducibility in the computational sciences has provoked discussion and effort across a broad range of perspectives: technological, legislative/policy-making, education, and publishing. Nature Biotechnology has recently implemented standards for reproducibility [1] in response to increasing complexity of methods and data requirements in computational biology. Discussion on these topics is not new [2,3,4], but the need to adopt standards for reproducibility of claims made based on computational results is now clear to researchers, publishers and policymakers. Many technologies exist to support and promote reproduction of results: this journal and others have discussed containerisation tools like Docker [5], literate programming approaches such as Sweave [6], knitr [7], iPython [8] or cloud environments like Amazon Web Services [9]. But these technologies are tied to specific programming languages (e.g. Sweave/knitr to R; iPython to Python) or to platforms (e.g. Docker for 64-bit Linux environments only). No single approach to date is able to span the broad range of technologies and platforms represented in computational biology and biotechnology.
To enable reproducibility across computational biology, we demonstrate for the first time an approach and a set of tools that is suitable for all computational work, and not tied to a particular programming language or platform. Our approach extends our previous work in this area [10] to produce ‘reference environments’ for reproducing results that can be deployed across platforms and can use any and all of the replication tools described above. We achieve this by adopting innovative open-source tools from software engineering [13,14], and we are now involved in the community contributing to the development of these tools, and acting as advocates for scientific computing within that community. We present published examples from a series of papers across research groups in different areas of computational and systems biology, spanning the major languages and technologies in the field (Python/R/MATLAB/Fortran/C/Java). We also include examples reproducing studies in network analysis by Feizi et al. [11] and Barzel et al. [12], results which have been recently commented upon in this journal [1]. We also present results from another study demonstrating that different output values can be produced from the same data and code given different versions of the MATLAB software used by the Feizi and Barabasi codebases. Since very few reviewers and readers will run the same code using multiple software versions, this could be interpreted as inability to reproduce the published work, underlining the need for a reference environment as produced by our approach. Our approach produces a transparent and flexible process for replication and recomputation of results, but ultimately its most valuable aspect is the decoupling of methods in computational biology from their implementation. Separating the ‘how’ (method) of a publication from the ‘where’ (implementation) promotes genuinely open science and benefits the scientific community as a whole.
[1] Editorial. Nat Biotech. Nature Publishing Group, a division of Macmillan Publishers Limited. All Rights Reserved.; 2015;33: 319.
[2] Barnes N. Nature. 2010;467: 753. doi:10.1038/467753a
[3] Merali Z. Nature. Nature Publishing Group; 2010;467: 775–777.
[4] Ball P. Nature. 2003;
[5] Docker. [cited 8 May 2015]. Available: https://www.docker.com/
[6] Leisch F. Hardle W, Ronz B, editors. Compstat 2002: Proceedings in Computational Statistics. Heidelberg: Physica-Verlag Gmbh & Co; 2002.
[7] Xie Y. Implementing Reproducible Research. Chapman and Hall/CRC; 2014. doi:doi:10.1201/b16868-3
[8] Pérez F et al. Comput Sci Eng. 2007;9: 21–29. doi:10.1109/MCSE.2007.53
[9] “Amazon Web Services.” [cited 8 May 2015]. Available: http://aws.amazon.com/
[10] Hurley DG et al. Brief Bioinform. Oxford University Press; 2014; bbu043. doi:10.1093/bib/bbu043
[11] Feizi S et al. Nat Biotechnol. 2013;31: 726–33. doi:10.1038/nbt.2635
[12] Barzel B et al. Nat Biotechnol. NATURE PUBLISHING GROUP, 75 VARICK ST, 9TH FLR, NEW YORK, NY 10013-1917 USA; 2013;31: 720–5. doi:10.1038/nbt.2601
[13] Vagrant. [cited 8 May 2015]. Available: https://www.vagrantup.com/
[14] Packer. [cited 11 Mar 2015]. Available: packer.io/
We present an optimal impulse control problem related to anti-cancer therapeutics. Our model takes into account both the instantaneous toxicities of the treatments and the adverse effects of targeted and non-targeted antineoplastic therapies.
Rabies is a fatal zoonotic disease and remains to be a priority health concern in the Philippines. Dogs remain to be the principal carrier of rabies and despite the fatality of rabies, it is a vaccine-preventable disease. Currently, the Philippine health officials have been conducting mass dog vaccination campaigns with the goal of having a rabies-free Philippines by the year 2020. However, statistical information about the disease poses a challenge and gives rise to doubts on whether the goal is achievable. In line with this goal, our aim for this study was to see whether eliminating dog rabies in the country is possible. Here, we developed a simple mathematical model, following the SEIR framework, that describes the dynamics of rabies infection among dogs in the presence of mass vaccination wherein some model parameters were calibrated from actual data sets. We used the model to predict the long-term rabies incidence and to analytically find an expression for the basic reproduction number in terms of vaccination rate, dog reproduction rate, transmission rate, and rate of vaccination immunity loss. By numerical simulations, we were able to determine the parameter combinations or scenarios for zero rabies infection. We assessed the validity and feasibility of these scenarios via data analysis. Our preliminary findings pointed to mathematical modelling as an important tool to make realistic projections about the rabies disease.
Recent studies have shown that pathogens can alter the behaviour of hosts or vectors to improve their transmission power. This happens not only in animal pathosystems but in plant pathosystems. While animal pathogens can alter the behaviour of both hosts and vectors in ways that increase frequency of host-host or host-vector encounters, in plant pathosystems the host does not have mobility, so the potential for behavioural manipulation is restricted to just vector mobile component in these systems [1].
Motivated by [2], we investigate how this manipulation affects transmission of specific insect-borne plant diseases. To do so, we use the compartmental model of vector preference in which the probability of which vectors settle on individual host plants. Then we demonstrate that interactions between host and vector infection status, the effect of a vector preference and difference of the effect between animal- and plant systems. together with any evolutionary implications.
[1] Ingwell, L. L., Eigenbrode, S. D., & Bosque-Pérez, N. A. (2012). Plant viruses alter insect behavior to enhance their spread. Scientific Reports, 2, 578.
[2] Cunniffe, N. J., Koskella, B., Metcalf, C. J. E., Parnell, S., Gottwald, T. R., & Gilligan, C. A. (2015). Thirteen challenges in modelling plant diseases. Epidemics, 10, 6-10.
p53, a cellular damage response regulator, is mutated in 50% of all cancers. We have constructed a mathematical model showing how different p53 dynamics after UV radiation and $\gamma$-irradiation exposure allows cells to respond properly to both slowly and quickly detected damage. We extend this model with an apoptosis module accounting for the transcription-dependent and -independent activities of p53, which is then used to analyze and classify contrasting apoptosis evasion strategies in commonly used tumour cell lines. In particular, we predict that lower levels of basal p53 contribute to the radiosensitivity observed in many cancers, and that self-overregulation of p53 is only one of many ways in which gain-of-function p53 mutations can be advantageous to tumour cell survival. Moreover, the downregulation of p53 by its canonical antagonist, Mdm2, has been well-studied in the nucleus, but this same interaction in the cytosol may be responsible for weakening the transcription-independent p53 pathway in humans. These models provide a unified framework for studying the broad and surprising array of effects of p53, far beyond the loss-of-function mutations well-known across cancers.
In Noor Asih et al. (2016) already given mathematical modelling as the dynamic of HPV infection on cervical cancer. Also given five scenarios for the existence of equilibrium points and their local stability. From the analysis of the system it is found the basic reproduction number is depends on the infection rate, the number of new virion that produce by infected cells, the death rate of virus, the growth rate of infected cells and progression rate. The existence and the local stability of equilibrium points are depend on basic reproduction number, the growth rate of pre-cancer cells and invasion rate. So we predict that there are some bifurcation parameters on that model. While we do some simulation by continuing those parameters we found some bifurcation such as fold bifurcation, cusp bifurcation and zero Hopf bifurcation. Further we analyze domain of each bifurcation parameter.
Avian malaria is a mosquito-borne parasitic disease of birds caused by protists of the genera Plasmodium, most notably Plasmodium relictum. This disease has been identified as a primary cause of the drastic decline and extinctions of endemic birds on Pacific Islands. In this work, we formulate an epizootiological model of the transmission dynamics of avian malaria between a generic bird species and mosquito using a system of ODEs. We derive the basic reproduction number as well as criteria for the existence and stability of disease-free and enzootic equilibria. We discuss strategies for minimizing the impact of avian malaria in two scenarios: disease-free populations which may be invaded by avian malaria and populations where this disease is enzootic but where bird species have not developed resistance.
Bond graphs are an energy-based framework for modelling physical systems while adhering to thermodynamic and physical constraints, and they have recently been extended and applied to biochemical and electrophysiological systems. Here we describe a bond graph model of the cardiac action potential and use it to explore the issue of drift in mathematical models of electrophysiology, which is a cause of inaccuracy. Previous studies have linked drift to stimulus currents that violate conservation laws, but those analyses were restricted to individual models and relied upon the manual identification of conservation laws. Due to their adherence to conservation principles, bond graphs are well-suited for studying conservation laws under a more general framework. Using our bond graph approach, we found that the conservation law derived in previous studies is an example of a ‘conserved moiety’ within the bond graph modelling and metabolic analysis frameworks. Conserved moieties explain the occurrence of drift for a general class of models, demonstrating that bond graphs can provide a systematic approach that can be used to identify hidden conservation laws and check for the presence of drift.
Consensus methods are widely used for combining phylogenetic trees into a single estimate of the evolutionary tree for a group of species. But how robust are these methods to future information? If additional species are added to the original set of trees, will the expanded consensus tree simply be an expansion of the original consensus tree? In this talk I will formalise and answer this question. Joint work with David Bryant (Otago, NZ) and Mike Steel (Canterbury, NZ).
Interspecies comparisons of physiological energetics are possible only through the prism of a formal metabolic approach such as the Dynamic Energy Budget (DEB) theory which uniformly describes how individuals of different species acquire and utilise energy. We used the DEB theory to infer the energy budgets of three commercial tuna species (skipjack, Pacific bluefin, and Atlantic bluefin) throughout all stages of ontogenetic development---from an egg to an adult individual and its eggs. Energy budgets were inferred from scarce and disjointed data sets fed into a DEB-based mathematical model tailored for tuna fish until reaching a high goodness of fit and thus the reliable estimates of the model parameters.
The results show that life histories of all three species are strongly influenced by morphological and physiological adaptations which accelerate ontogeny during the larval stage, although the effect is more pronounced in bluefin than skipjack tuna. We identify that in energetic terms the accelerated ontogeny is a simultaneous improvement of energy acquisition (higher intake) and utilisation (higher expenditure) without changing the capacity of fish to build energy reserve as intake and expenditure increase in unison. High energy expenditure, an even higher intake by necessity, and a limited capacity to build energy reserve, make all three tuna species vulnerable to starvation, thereby theoretically underpinning the description of tuna fish as “energy speculators”. We furthermore find that energy allocation to reproduction maximises fecundity of all three tuna species, thus suggesting that the evolution of tuna favoured higher fecundity at the expense of growth.
Thinking beyond just physiological energetics (e.g., wild stock projections), we discuss why DEB-based models are a natural foundation for physiologically structured population dynamics wherein the environment influences the population growth rate via metabolism. We suggest several concrete modelling techniques which allow progress in this direction.
IUCN criteria are the most authoritative and objective to assess the conservation status of animal species. Although IUCN criteria are purely descriptive in nature, they can be interpreted as e.g. the relative size of the annual population growth factor (λ, the dominant eigenvalue of a projection matrix). This enables quantitative assessment based on demographic data (survival and reproduction). However, such an assessment for many species is hampered since demographic data are often incomplete. To facilitate assessment when demographic data are incomplete we extended an existing framework (developed for species with a pre-adult stage of one year) to apply to species with pre-adult stages (=juvenile and sub-adult together) longer than one year. We develop both Leslie-matrices and stage-structured matrices and compare these with respect to eigenvalues, and both sensitivity and elasticity of underlying parameters. We illustrate our work with life history data and conservation status of a few marine mammal species.
The construction of effective and informative landscapes for stochastic dynamical systems has proven a long-standing and complex problem, including for biological systems. Such landscapes may refer to a true energy function for cases such as protein folding or to a phenomenological metaphor in the case of Waddington’s epigenetic landscape. In many situations, constructing a landscape comes down to obtaining the quasi-potential, a scalar function that quantifies the likelihood of reaching each point in the state-space.
In this work we provide a novel method for constructing such landscapes using a tool from control theory: the Sum-of-Squares method for generating Lyapunov functions. Applicable to any system described by polynomials, this method provides an analytical polynomial expression for the potential landscape, in which the coefficients of the polynomial are obtained via a convex optimisation. The resulting landscapes are based upon a decomposition of the vector-field of the original system, such that the inner product of the gradient of the potential and the remaining dynamics is everywhere negative. By satisfying this condition, our derived landscapes provide both upper and lower bounds on the true quasi-potential; these bounds becoming tight if the decomposition is orthogonal. The method is implemented in the programming language Julia, and is demonstrated to correctly compute the quasi-potential for high-dimensional linear systems and also for a number of nonlinear examples. For a multi-stable stochastic system analogous to a developing stem cell, we use the computed potential to evaluate bounds on the relative likelihood of reaching fixed points equivalent to each of the differentiated phenotypes.
Hematopoietic system is maintained by hematopoietic stem cells (HSCs) with dual abilities of long-term self-renewal and differentiation to all types of blood cells. Recently, using a single-cell transplantation system and mice expressing a fluorescent protein, myeloid-restricted progenitors with long-term repopulating activity (MyRPs) were found. Moreover, by using paired daughter cell assay, MyRPs were directly differentiated from HSCs.
In this study, we investigated hematopoietic system incorporating the novel insight that there existed a cell type that exclusively differentiated to myeloid lineages. There were four types of populations in the model: (i) HSCs, (ii) MyRPs, (iii) progenitors, and (iv) differentiated cells. Differentiated cells includes B cell and T cell which are lymphoid cells and platelet, erythrocyte and neutrophil/monocyte as myeloid cells. Myeloid progenitors were produced via two ways, from HSCs directly and via MyRP (myeloid bypass), after transplantation of a single HSC, while lymphoid progenitors were produced from only HSCs directly. This is the first study of investigating hematopoiesis with MyRPs.
We estimated some parameters which were growth rate, production rate and death rate using full data set of single-cell transplantation. From the result of data analysis, we will discuss the role of myeloid bypass after transplantation and the change of heterogeneity of HSC population with age.
Cellular adaptation to varying signal strengths has been observed in the innate immune cells responses to external and internal challenges. Monocytes persistently challenged with low levels of external stimulants, can be skewed into a low-grade ‘primed’ state, while persistent challenges with high levels of external stimulants will skew them into to a ‘tolerant’ state. Based on experimental and modelling studies, we have identified a generic network motif, composed of mutually competitive signals, as a potential principle for bistable and/or intermediate priming and tolerance responses. We use analytical and numerical techniques to identify the detailed mechanisms for the activation magnitude, duration, and potential transition between priming to tolerant states.
We propose a structured integro-difference equation model for an invasive marine species with a pelagic larval stage and examine the role of dispersal heterogeniety on the spreading speed. The spread of the green crab up the northwest coast of the Atlantic is used as a case study. We find that the relationship between spreading speed and demographic and dispersal parameters is similar to the relationship found in Fisher's equation. We also find that temporal variation in dispersal results in a faster spread rate than predicted by a time-averaged dispersal kernel. This is joint work with Lin Wang, Myriam Barbeau and Ali Gharouni.
Motivated by in vitro time-lapse images of ovarian cancer spheroids inducing mesothelial cell clearance, the traditional agent–based model of cell migration, based on simple volume exclusion, was extended to include the possibility that a cell seeking to move into an occupied location may push the resident cell, and any cells neighbouring it, out of the way to occupy that location. In traditional discrete models of motile cells with volume exclusion such a move would be aborted. We introduce a new shoving mechanism which allows cells to choose the direction to shove cells that expends the least amount of shoving effort (to account for the likely resistance of cells to being pushed). We call this motility rule ‘smart shoving’. We examine whether agent-based simulations of different shoving mechanisms can be distinguished on the basis of single realisations and averages over many realisations. We emphasise the difficulty in distinguishing cell mechanisms from cellular automata simulations based on snap-shots of cell distributions, site-occupancy averages and the evolution of the number of cells of each species averaged over many realisations. This difficulty suggests the need for higher resolution cell tracking.
Obesity is the result of caloric imbalance and is mediated by genetic, behavioural, and environmental factors. Healthy lifestyle habits, including healthy eating and physical activity, can lower the risk of becoming obese and developing related diseases. The prevalence of obesity among Korean adolescents aged 13 to 18 years increased from 13.65% in 2007 to 19.3% in 2016 for boys.
We analyze trends in adolescents obesity for 6 years and estimate of adolescents obesity. We adopt the well-known SIR model structure for infectious diseases. A differential equation system predicted adolescents' obesity prevalence trends. The model considers both social environment and obesogenic environment influences on weight gain, incorporates other known parameters affecting obesity trends.
The dynamic model predicts that obesity prevalence will plateau independent of current prevention strategies. The proportion of overweight is lower in middle school than in high school, but obesity is higher in high school. The adolescents prevalence of overweight and obesity for boys will plateau by about 2025 at 6% and 27% respectively.
[1] H. Pontzer, R. Durazo-Arvizu, et al., Constrained Total Energy Expenditure and Metabolic Adaptation to Physical Activity in Adult Humans, Current Biology, 26 (2016), 410-417.
[2] Diana M. Thomas, Marion Weedermann, et al., Dynamic Model Predicting Overweight, Obesity, and Extreme Obesity Prevalence Trends, Obesity, Vol 22, Num. 2, (2014) pp.590-48. www.obesityjournal.org.
[3] https://yhs.cdc.go.kr/new/pages/pds1.asp, (2016) (accessed 2017.05.01)
Recent advances in understanding of complex interactions between the processes that maintain and control the physiological state of brain parenchyma open new possibilities and deliver new challenges for modelling studies on topic. Now it is clear that any significant deviation from normal conditions as well as natural alternations of activity of cortical neurons (say, during sleep-wake cycle) are accompanied by measurable changes of physiological parameters, like extracellular ionic concentrations, astrocytic calcium signalling, or vascular tone of adjacent blood vessels. The interaction pathways within the so called neural-glial-vascular unit not only drive the local responses of the cells, but form the spatially extended networks at different functional layers.
Being normally on the sidelines, these physiological mechanisms become dominant during the extreme physiological behaviours of the brain cortex, such as cortical spreading depression, migraine waves, or spreading depolarization events with stroke or injures. There is a clear room and need for relevant modelling studies aimed to better understanding of how all this activity is orchestrated. However, this task is very challenging due to the extreme complexity and heterogeneity of system under study.
The presented study addresses the selected aspects of above problem. Being focused on dynamical patterns, we still keep the balance between the physiological details and dynamically tractable generalized modelling. First, we discuss the dynamical consequences of inclusion of extracellular potassium concentration as a state variable in neuron models. Specifically, this additional variable provide both the self -sustained depolarization, and the additional pathway for cell-to-cell signalling which is also known as "volume connection". Second, we extend this potassium driven model to the mathematical model of cortical spreading depression which counts the effects of neurovascular coupling and cerebral blood flow redistribution. Our most notable finding here is that the combination of vascular-mediated spatial coupling with local regulatory mechanisms results in the formation of stationary Turing-like patterns during a course of spreading depolarization.
Third, we describe a way to incorporate in model the recently revealed details of astrocyte morphology and functions. This approach includes the algorithms for the creation and computational treatment of spatial patterns resembling the real astrocyte networks.
In conclusion, as an outlook, we discuss how the dynamical redistribution of extra- and inter-cellular volumes could affect the reported dynamics, as well as formulate the issues related to adequate modelling of vascular-mediated spatial interaction.
A deterministic model of biochemical reaction networks is not always appropriate. Gene expression, for example, often involves a small number of molecules which means that noise can significantly influence the dynamics. By discretising space into voxels and letting the molecule dynamics be governed by the reaction-diffusion master equation, it is possible to model the reaction and diffusion of individual molecules on an arbitrary domain. Following on from the work of Meinecke and Lötstedt we apply a variety of numerical methods to the Laplacian operator in order to derive the jump rates that simulate the diffusion of molecules. We discuss how pattern formation within a Turing model might be influenced by the geometry of the spatial discretisation and the numerical method.
The demographic hypothesis of cumulative cultural evolution claims that population size has been a crucial determinant of the rate of cumulative cultural evolution in humans rather than other factors such as environmental risk. The original version of this hypothesis does not distinguish the role of population size from that of social connectedness, i.e. the degree to which individuals in a population are socially connected with each other. Recently, one of our models showed that social connectedness may actually have a larger impact on the rate of cumulative cultural evolution than population size. However, these models all assume that the mode of learning which underlies cultural transmission is fixed and does not evolve, so that they cannot deal with the time scale of genetic evolution. In the present study, we extend previous models of cumulative cultural evolution by assuming that the allocations of the individual lifetime into two modes of learning, social and individual learning, co-evolve with the level of culture in the population. We derive the evolutionarily stable learning strategy, assuming that cultural evolution occurs on a faster time scale than genetic evolution. The results suggest that population size may have only a minor effect on the equilibrium level of culture realized by the evolutionarily stable learning strategy, while social connectedness may have a relatively large effect, conforming with the results of models without evolution of learning. These results indicate the importance of distinguishing different aspects of demography when modelling cumulative cultural evolution.
The temporal and spectral characteristics of tonic-clonic seizures are investigated by using a neural field model of the corticothalamic system in the presence of a temporally varying connection strength between the cerebral cortex and thalamus. Increasing connection strength drives the system into 10 Hz seizure oscillations once a threshold is passed and a subcritical Hopf bifurcation occurs.
In this study, the spectral and temporal characteristics of tonic-clonic seizures are explored as functions of the relevant properties of physiological connection strengths, such as maximum connection strength, time above threshold, and the rate at which the connection strength increases or decreases (ramp rate). Dynamical analysis shows that the seizure onset time decreases with the maximum connection strength and time above threshold, but increases with the ramp rate. Seizure offset time and duration increase with maximum connection strength, time above threshold, and rate of change. Spectral analysis reveals that the power of nonlinear harmonics and the duration of the oscillations increase as the maximum connection strength and the time above threshold increase. A secondary limit cycle, termed a saddle cycle in previous studies, is also seen in this study. A detailed analysis of the saddle cycle oscillations shows that these oscillations become more prominent and robust with maximum connection strength and rate of change of the ramp. However, for a small ramp rate, the system does not exhibit any saddle cycle oscillations. We also find that if the time above the threshold is too small and the ramp rate is too large, the system does not reach to the larger limit cycle attractor of 10 Hz oscillation, and only exhibits the saddle cycle oscillations. It is also seen that the times to reach the saturated large amplitude limit-cycle seizure oscillation from both the instability threshold and from the end of the saddle cycle oscillations are inversely proportional to the square root of the ramp rate, as is the time to reach the seizure offset from the bifurcation from the saturated limit cycle oscillations.
Projection matrix models are known to provide us with a plenty of population statistics, such as population growth rate, steady size-class distribution, and sensitivity and elasticity for population growth rate. Hundreds of academic papers using the model have been published these last forty years and a database on many of their matrices is now available on the internet (COMPADRE and COMADRE), which contains the demographic data on more than a thousand species. Silvertown et al. (1996) published a famous paper, where they mapped elasticity vectors of survival, growth and fecundity for 84 plant species in a triangle simplex and found that they are located in a specific region. The same trend is found on the map for 1230 plant populations in the above plant database. To understand and clarify why they are located in a specific region, we constructed five types of random matrices. 4 by 4 random matrices were composed of two parts: fecundity and transition probabilities from a stage to another. The distribution of fecundities followed a Poisson distribution. The transition probabilities range from zero to one, whose row sums are less than 1. The elasticities for survival, growth and fecundity were calculated using 3000 random matrices and the elasticity vectors were plotted in the triangle map. The five types of matrices were as follows: (1) random matrices with no zero-element, (2) random matrices with no zero-element and the survival probabilities increase as individuals grow, (3) random matrices which have non-zero elements only on diagonal and sub-diagonal positions, (4) random matrices which have non-zero elements only on diagonal and sub-diagonal positions and the survival probabilities increase as individuals grow, (5) random matrices in semelparous species. The results are: (a) the distribution of the elasticity vectors moves to upper-left region of the triangle map as average of fecundity increases, (b) In the third and fourth types of random matrices, the distribution is located on a line whose slope is equal to 46 degrees. The slope can be described by a formula depending on the matrix dimension, n and ranges from 30 to 60 degrees with n = 2 to infinity. (c) In semelparous species, the distribution moves to the upper left along the 46-degree line. (d) There are no elasticity vectors in the bottom half of the triangle map.
Many large-scale high-throughput experiments use DNA barcodes—short DNA sequences prepended to DNA libraries—for identification of individuals in pooled biomolecule populations. However, DNA synthesis and sequencing errors confound the correct interpretation of observed barcodes and can lead to significant data loss or spurious results. Widely-used error-correcting codes borrowed from computer science (e.g., Hamming and Levenshtein codes) do not properly account for insertions and deletions in DNA barcodes, even though deletions are the most common type of synthesis error. We present and experimentally validate FREE (Filled/truncated Right End Edit) barcodes, which correct substitution, insertion, and deletion errors, even when these errors alter the barcode length. FREE barcodes are designed with experimental considerations in mind, including balanced GC content, minimal homopolymer runs, and reduced internal hairpin propensity. We generate lists of barcodes with different lengths and error-correction levels that may be useful in diverse high-throughput applications, including $>10^6$ single-error correcting 16-mers that strike a balance between decoding accuracy, barcode length, and library size. Moreover, concatenating two or more FREE codes into a single barcode increases the available barcode space combinatorially, generating lists with $>10^{15}$ error-correcting barcodes. Our software for creating barcode libraries and decoding sequenced barcodes is efficient and designed to be user-friendly for the general biology community.
Specific biomarkers can be identified in dynamic contrast-enhanced magnetic resonance imaging (DCE-MRI) breast scans and quantified using pharmacokinetic models that return estimates of parameters related to tissue physiology including vessel perfusion and permeability ($K^{trans}$), the extravascular-extracellular volume fraction ($v_e$), the plasma volume fraction ($v_p$), and the efflux constant ($k_{ep}$). In particular, $K^{trans}$ and $k_{ep}$ have been shown to be effective at predicting the response of cancer patients to treatment. Two fundamental issues in the field of DCE-MRI is the lack of standardization of the analysis and characterizing the time rate of change of the concentration of contrast agent in the vascular (the so-called “arterial input function” or AIF). We have recently developed a method for estimating accurate AIFs for the individual patients and associated software to automate the estimation of model parameters from DCE-MRI data taken from breast cancer patients using data that can be acquired routinely in community-based imaging centres.
Cell proliferation is the most important cellular-level mechanism responsible for regulating cell population dynamics in living tissues. Modern experimental procedures show that the proliferation rates of individual cells can vary significantly within the same cell line. However, in the mathematical biology literature cell proliferation is typically modelled using a classical logistic equation with a constant proliferation rate, and this approach neglects variations in the proliferation rate. In this work we consider a discrete mathematical model of cell migration and cell proliferation, modulated by volume exclusion (crowding) effects, with variable rates of proliferation across the total population. We refer to this variability as heterogeneity. Constructing the continuum limit of the discrete model leads to a generalisation of the classical logistic growth model. Comparing numerical solutions of the model to averaged data from discrete simulations shows that the new model captures the key features of the discrete process. Applying the extended logistic model to simulate a proliferation assay using rates from recent experimental literature shows that neglecting the role of heterogeneity can, at times, lead to misleading results.
In this research, we study feedback control problem in the context of deterministic epidemic models. Feedback control is obtained by solving Hamilton-Jacobi-Bellman(HJB) equation, which is employed to overcome limitations of previous work. There are three key factors in the implementation of this methodology, decoupling value function and control variables, truncation of unbounded domain, and numerical solver for 1st order hyperbolic PDE. While this approach seems complicated, it has an obvious advantage in generalization to stochastic optimal control problem. With proper treatments for technical challenges, we provide a tool that can be widely applied.
Stochastic evolutionary game dynamics in finite populations has recently been examined not only for symmetric games [1] but also for bimatrix games [2]. While in these studies the fixation probabilities of pure strategies are investigated, this study examines the evolutionary dynamics of two-player 2 by 2 bimatrix games with mixed-strategies in finite populations under weak selection. The game has two populations for row and column players. In the population consisting of row (column) players, two types of players, i.e., mutant and wild-type players, have different mixed-strategies which assign different weight to the two row (column) strategies. Each player’s fecundity is determined by the average payoff obtained by interactions with all players in the opposite population. The process of death and birth of players is modelled by the frequency-dependent Moran process.
For this game, I derived the fixation probabilities that a pair of mixed-strategies by mutants in the row and column populations takes over the entire populations from a given initial state. Moreover, based on the probabilities, I discuss the stochastic stability of the pair of mixed-strategies played by wild-type players.
In addition, the fixation probabilities of mixed-strategies of mutant players when the selection is weak because of mixed-strategies played by mutant and wild-type players being very close in terms of a probability distribution over the set of possible strategies, which are the counterparts of Wild and Traulsen [3] in bimatrix games, are also discussed.
[1] Nowak, M. A., Sasaki, A., Taylor, C., & Fudenberg, D. (2004). Emergence of cooperation and evolutionary stability in finite populations. Nature, 428(6983), 646.
[2] Sekiguchi, T., & Ohtsuki, H. (2017). Fixation probabilities of strategies for bimatrix games in finite populations. Dynamic Games and Applications, 7(1), 93-111.
[3] Wild, G., & Traulsen, A. (2007). The different limits of weak selection and the evolutionary dynamics of finite populations. Journal of Theoretical Biology, 247(2), 382-390.
Immunotherapy using checkpoint inhibitors has demonstrated clinical efficacy for cancers ranging from melanoma to non-small-cell lung cancer and other types. Despite some success, many patients do not respond to the therapy, and a subset of patients develops “hyper-progressive” disease. Biomarkers that predict response or toxicity to checkpoint inhibitors will help to select patients who are most likely to benefit from these novel immune agents. Clinical-grade next generation sequencing is becoming more prevalent, providing information about genomic alternations, which could serve as potential new biomarkers to immunotherapy agents.
In this retrospective study, we included 43 patients treated in a community cancer center who had diverse solid organ malignancies, completed sequencing of tumour DNA before the initiation of immunotherapy, and received at least one cycle of a PD-1 or PD-L1 inhibitor. Half of the patients had breast and gynecologic cancers, which have been less investigated for checkpoint blockade.
We found that response to checkpoint inhibitors was associated with i) fewer previous treatment lines, ii) longer duration of immunotherapy, iii) higher frequencies of peripheral blood monocytes and lymphocytes after treatment, and iv) higher tumour mutational burden (TMB). In addition, base substitutions and indels in PRKDC and LRP1B and amplification of BCL6 occurred more frequently in responders. PRKDC and LRP1B mutations showed significant association with higher levels of TMB, which was confirmed by large-scale cancer genomics data.
Results from this study may be used to inform patients who will have a better response to PD-1/PD-L1-based immunotherapy, possibly indicating a first-line therapy in their treatment.
Improving crop yield is essential to meet increasing global food demands. Crop yields depend on the coordinated acquisition of carbon and nitrogen by the leaves and roots respectively and the use of these nutrients within each part of the plant. Changes in environmental conditions cause fluctuations in carbon and nitrogen availability. This leads to crosstalk between the signalling pathways for carbon and nitrogen. Carbon- and nitrogen-derived signals have been observed in plant growth, but how these signals cooperate together with changes in nutrient availability is still unknown. Here, we discuss the implementation of such feedback mechanisms using a carbon and nitrogen transport model for plant growth.
The Hybrid Agent–based Library (HAL) is a Java Library made of simple, efficient, generic components which can be used to model complex spatial systems. HAL's components can broadly be classified into: on and off lattice agent containers, finite difference diffusion fields, a GUI building system, and additional tools and utilities for computation and collecting data. These components were designed to operate independently, but are standardized to make them easy to interface with one another.
HAL is a useful asset for researchers who wish to build efficient 2D and 3D hybrid models in Java, while not starting entirely from scratch. It is available on GitHub at https://github.com/torococo/HAL under the MIT License. HAL requires at least Java 8 or later to run, and the java jdk version 1.8 or later to compile the source code.
The spatial components in HAL all use a consistent indexing scheme. Setting values in a diffusible field, for example, is syntactically similar to setting visualization pixel values, or finding the agents that occupy a position. The consistent syntax helps substantially with HAL’s learning curve. It also helps with understanding the source code of models that use components that the reader is not familiar with.
There are several pre-existing general hybrid agent-based modelling frameworks, the most popular being Chaste, Repast, Mason, and Netlogo. Each of these facilitates modelling under a different centralized control structure: In Chaste centralized control done by a Simulator object, in Repast this component is called an Engine, in Mason it is called the Schedule object, and in Netlogo it is called the Go loop.
HAL shares many characteristics with these frameworks, but primarily differentiates itself with its simplicity and decentralized design. There is no controller or scheduler, so the modeler designs the logical flow and the scheduling of interactions between components of the model explicitly. This allows the modeller to define the execution of the model using simple programming constructs such as “for loops”, and sheds the difficulty of learning to manipulate a more complex control architecture.
The aim of this research was to explore the potential of proﬁle likelihoods for identiﬁability analysis in the context of real enzyme kinetics data, collected ourselves.
Parameter identifiability concerns the question of whether the type of experimental data we have collected properly determines the parameters of our mathematical models. Identifiability issues arise because not all biological variables involved in the system can be measured, and even those that can be measured can be structurally decoupled from some of the parameters of interest. Although structural non-identifiability is in principle an all-or-nothing concept, in practice parameters may also be only weakly identified or may be practically non-identifiable given finite data. Profile likelihood has proven to be one of the few promising general methods of identifiability analysis that is applicable to multi-parameter, nonlinear problems, but it has not yet been widely adopted in the mathematical/computational biology community. Thus we sought to further explore its usefulness for typical models in this area, and using real experimental (as opposed to synthetic) data.
In this study, we collected data on the activity of the enzyme tissue plasminogen activator (tPA) under variety of scenarios, including different initial substrate and pH levels. We developed a series of simple reaction kinetics models, including both Michaelis-Menten velocity-concentration models and full time-dependent ODE models, and generated profile likelihoods under the various experimental conditions. For the simple Michaelis-Menten model we found that parameters were generally identifiable/weakly identifiable but tended to become less identiﬁable (approaching practical non-identifiability) at lower pH levels. On the other hand, individual parameters of the full ODE model of enzyme kinetics showed full structural non-identiﬁability. This led us to consider the identifiability of targeted ‘interest’ parameters, motivated by the parameters in the simpler system. Using this approach we found that certain combinations of rate parameters, corresponding to those in the simpler Michaelis-Menten model, were better identified in the full model.
Overall we found proﬁle likelihood to be a promising technique for identiﬁability analysis of enzyme kinetics models. For complex models, however, choosing targeted interest parameters appears to be essential to avoid structural non-identifiability. Further work is needed on systematically motivating these interest parameters based on, for example, simpler models and/or model reduction procedures. In the context of tPA kinetics, more complex reactions involving the inhibitor neuroserpin and interactions of H+ ions with the enzyme should be considered.
In genetics, enhancer regions are short non-coding strands of DNA belonging to a class of elements known as cis-regulatory. These regions, typically in intergenic regions, encompass binding sites for proteins (transcription factors) that regulate gene expression. Changes to the DNA sequence in these enhancer regions are thought to provide a mechanism to explain changes of gene expression even though the genes across species may display some sequence conservation.
The challenge in the identification of these enhancers is that, unlike protein-coding genes, they do not follow a clear genetic code. Sequence conservation across species is often used as a proxy to identify regulatory regions with the rationale that functional elements are under some evolutionary pressure to maintain biological function and, hence, less subject to mutations. However, there is emerging evidence that functional non-coding sequences that show conservation of biological function in the combination of transcription factor binding sites but show little to no conservation at a sequence level using traditional alignment tools. Hence, our ultimate goal is the identification of conserved non-coding sequences by alternative methods to “traditional” sequence alignments.
In this work, we obtain an alignment between the DNA assemblies of human (hg38) and mouse (mm9) around the gene tbx5 and perform a cursory statistical analysis to identify and classify non-coding functional regions. The model is developed in a Bayesian framework and we obtain a segmentation using the segmentation classification algorithm, changept. The algorithm samples segmentations using a Markov chain Monte Carlo method that generalises the Gibbs sampler. We aim to incorporate a diverse set of data types into the model to increase accuracy and applicability.
We consider a setup in which $n$ particles are initially released into a domain and diffuse freely. Part of the boundary consists of absorbing "escape" regions, where the particles can escape the domain, and reflecting regions. The rest of boundary consists of "capture" regions (receptors), that can switch between being reflecting and absorbing. Specifically, after capturing a particle, the capture region becomes reflecting for an exponentially distributed amount of time (recharge time). We are interested in the distribution of the number of particles that are captured before they escape.
Our mathematical results are derived from considering our system in several ways: as a full spatial diffusion process with recharging traps on the boundary; as a continuous-time Markov process approximating the original system; and lastly as a system of ODEs in a mean-field approximation.
Considering the full spatial diffusion process, we prove that the total expected number of the captured particles has an upper-bound of the order of (log $n$). We then consider a few examples investigating the implications of this result, including a neural system in which recharge reduces the number of neurotransmitter bindings by several orders of magnitude.
Seeking additional understanding of this process over a wide range of parameter values, we then present a well-defined algorithm to approximate the full spatial diffusion process with a continuous-time Markov process in order to eliminate computationally expensive simulations. We highlight the conditions required for the approximation to yield similar quantitative results as the full spatial process. We then apply the approximation, and the associated mean-field model, to investigate time courses for the expected number and higher ordered statistics of captured particles. Specifically, we find that the number of expected captures as a function of time appears to grow linearly, before leveling off, and find an analytical expression for the duration of the linear growth. This result may be particularly helpful for understanding applications where multiple puffs of particles are inserted in the domain over a period of time (e.g. neuronal synapses), and can provide a bound on the time between puff events such that particles in different puffs no longer interact. We also find that the amount of variation observed in the total number of captured particles varies non-monotonically with the mean recharge time. Lastly, we combine these results together to predict stochastic properties of intracellular signals resulting from receptor activation.
Information transfer is thought to play a key role in adaptive complex systems such as social insects, brain and human society, while epidemic spreading through interactions can impact on survival in the individuals and society. Social insects such as ants and bees are one of the most sophisticated complex systems exhibiting collective decision making, providing us the opportunity of studying how the information transfer and disease spreading occur through the interactions. Recent advances in tracking systems enable us to collect quantitative and massive data of each individual and interactions in a colony. A recent empirical study [1] used an automatic tracking method based on image analysis, and revealed rapid spreading dynamics in honey bee colonies. However, the attributes of each individual including age, caste, and activity were not considered, although such heterogeneity of individuals may be crucial for the adaptability of the colony. Thus, the relationship between the heterogeneity in the colony and the spreading dynamics is less understood. In the present study, we used an image-based tracking system (Bugtag; Robiotec Inc.) to detect the movements and interactions of ants (Diacamma sp.). To reveal the relationship between the spreading dynamics and the heterogeneity in individual attributes, we investigated how interactions spread the information or diseases by using temporal networks which are promising tools for understanding spread dynamics, and then quantified heterogeneity in individual attributes such as age. The results show that the spreading was more rapid than randomized data conserving degree sequence, and younger individuals tended to be influencers which promoted the rapid spreading. We will discuss the adaptive significance of the relationship between rapid spreading dynamics and differences in individual attributes.
[1] Gernat, Tim, et al. "Automated monitoring of behavior reveals bursty interaction patterns and rapid spreading dynamics in honeybee social networks." Proceedings of the National Academy of Sciences 115.7 (2018): 1433-1438.
Primates exhibit an array of mating behaviours and arrangements, from monogamy to promiscuity, or mate-guarding to multiple-mating. We have previously explored the role of adult sex ratio in determining the likelihood of males choosing either promiscuity or monogamy, but many questions remain open. In particular, do (and, if so, how do) longevity and life history contribute to the choice of male strategy? And how can we characterise the effectiveness of guarding for those males who employ that strategy? We examine some of the trade offs that arise, and how they can influence the choice of male strategy.
The unmarried rate of Japanese people is rising year by year, and countermeasures are required. We modelled male and female marriage decision-making models based on male’s willingness to cooperate on housework and childcare. We assume that male marriage utility function decreases as male’s willingness to cooperate on housework and childcare increases, while female marriage utility function increases. We also assume that a male decides to marry a female when his marriage utility function is nonnegative, while a female does so when her partner’s willingness to cooperate is at least the value that she demands. When both male’s and female’s conditions are satisfied, they get married. If they don’t get married, the utility function is 0 for both males and females. Male’s strategy is the value of his willingness to cooperate, and female’s strategy is the value that she demands. There are Nash equilibria with marriage and without marriage. At Nash equilibria with marriage, both male and female marriage utility functions are positive and male’s strategy is equal to female’s.When the utility function changes by individuals, increase of marriage utility raises male’s marriage fraction slowly in the all range and female’s sharply in the range of large female marriage utility. We conclude that we should take countermeasures to all males and females with greater marriage utility to increase the rate of marriage.
True slime mold Physarum polycepharum is a large uni-cellular amoeba-like organism and it can sense environmental information and change its behaviours. Large true slime mold extends its frontal parts for foraging. To study its exploration strategy, we observed how physarum spreads in a binary-tree-shaped maze. We found that slime mold changes its searching strategies when its total mass changes. To explain this result, we made a new mathematical model which consists of cylinders filled with liquid, periodic contracting active springs, tubes connecting between cylinders and the elastic sheet. In our model, the pressure has an important role to transmit information among frontal parts. Our model reproduces fundamental characteristics of the front expanding and the searching behaviours of slime mold.
Alzheimer’s disease (AD) is a devastating illness affecting over 40 million people globally. The accumulation of AD associated Amyloid beta (Aβ) oligomers can trigger aberrant intracellular calcium signals by disrupting calcium regulatory mechanism within neurons. These disruptions can cause changes in homeostasis levels that can have detrimental effects on cell function and survival. Although studies have shown that Aβ can interfere with various calcium fluxes, the complexity of these interactions remains elusive. In order to better understand the impact of Aβ on calcium dynamics, we use a mathematical model to simulate calcium patterns under the influence of Aβ. More specifically, we assume that Aβ affects individual flux contributions through inositol triphosphate receptors, ryanodine receptors, and the plasma membrane. We show that the inclusion of Aβ can increase regions of mixed-mode oscillations leading to aberrant signals under various conditions. We use single and double parameter bifurcation structures to predict model solutions for various levels of Aβ. We further demonstrate that controlling certain biophysically relevant parameters can help control aberrant signalling. These results can be used to suggest possible targets for establishing therapeutic strategies in AD pathology.
Biofilms are sessile communities of bacteria housed in a self-produced adhesive matrix consisting of extracellular polymeric substances (EPS), including polysaccharides, proteins, lipids, and DNA. [1]. Biofilm provokes chronic bacterial infection, infection on medical devices, deterioration of water quality, and the contamination of food [2]. On the other hand, biofilm can be used for wastewater treatment and bioenergy production [3]. In microbial enhanced oil recovery (MEOR), one of the strategies is selective plugging, where bacteria are used to form biofilm in the high permeable zones to diverge the water flow and extract the oil located in the low permeable zones [4]. Therefore, it is necessary to build mathematical models that better describe the biofilm mechanisms. One of the motivations to derive upscaled models is to describe the averaged behaviour of the system in an accurate manner with relatively low computational effort compared to fully detailed calculations starting at the microscale [5]. In the laboratory, biofilm is growth in a T-shape micro-channel. We built a mathematical model including water flux inside the biofilm and different biofilm components (EPS, water, active bacteria, and dead bacteria). Using the best estimate of physical parameters from the existing experiments, we perform numerical simulations. The stress coefficient is selected to match the experimental results. A sensitivity analysis is performed to identify the critical model parameters. A reduction of the biofilm coverage area as the water flux velocity increases is observed. Homogenization techniques are applied in a strip and a tube geometry. Numerical simulations are performed to compare both upscaled mathematical models. In the macro-scale laboratory experiments, biofilm is growth in cylindrical cores. Permeability changes over time at different flow rates and nutrient concentrations are studied. Numerical simulations are performed to compare with the experimental results.
[1] Aggarwal, S., Stewart, P. S., Hozalski, R. M. (2015). Biofilm Cohesive Strength as a Basis for Biofilm Recalcitrance: Are Bacterial Biofilms Overdesigned? Microbiology Insights. 8s2, MBI.S31444.
[2] Kokare, C. R., Chakraborty, S., Khopade, A. N., Mahadik, K. R (2009). Biofillm: Importance and applications. Indian J. Biotechnol. 8, 159-168.
[3] Miranda, A. F. et al. (2017). Applications of microalgal biofilms for wastewater treatment and bioenergy production. Biotechnol. Biofuels. 10, 120.
[4] Raiders, R. A., Knapp, R. M., McInerney, M. J. (1989). Microbial selective plugging and enhanced oil recovery. J. Ind. Microbiol. 4(3), 215-229.
[5] van Noorden, T. L., Pop, I. S., Ebigbo, A., Helmig, R. (2010). An upscaled model for biofilm growth in a thin strip. Water Resour. Res. 46, W06505.
The occurrence of distant metastases greatly reduces or even removes the possibilities of curative treatment for cancer patients. The development of tools for predicting and diagnosing the onset of secondary tumours has long been a significant subject of cancer research effort. Moreover, it is widely known that while in some tissues cancer cells exhibit higher proliferation, in others they are more invasive so the best treatment strategy should be chosen accordingly.
A novel tissue engineering approach was created in our group in order to gain insight into the processes of cell growth, motility and invasion. Organ-specific scaffolds were seeded with highly invasive triple negative breast cancer cells. These tissue engineering constructs (TECs) were then cultured in vitro for 4 weeks and the resulting patterns of cell colonisation of extracellular matrix (ECM) were analysed.
One of the most noteworthy results of this experiment is the emergence of two alternative strategies of matrix colonisation depending on the origin and structure of the underlying ECM. In particular, we observed significant differences in cellular attachment efficiency, morphology, density and invasion of the breast cancer cells.
To gain a better understanding of cellular behaviours and validate the hypothesis that colonisation pattern strongly depends on the properties of ECM, we developed a corresponding mathematical model of tumour growth in TECs. For these purposes, we used a lattice-based stochastic approach based on the energy formalism, where the concept of cell adhesion is considered the key mechanism of cell-cell and cell-ECM interactions. Our model reproduces the early stage of the experiment (the stage characterised by rapid colonisation of surface and subsurface matrix layers) showing a clear distinction in cell morphology and spatial distribution when varying the density of adhesion molecules expressed on the matrix.
Further, to expand the computational model and attribute the role of the underlying ECM to cell invasion and clustering inside the tissue, we ran a set of single-factor experiments, each designed to reveal the effect of a specific physical property of the matrix. Statistical analysis of the experimental data provides estimations of model parameters and complements the picture of cell-matrix interactions.
We believe that this data-driven approach will allow predicting phenotypical and morphological changes of cancer cells and the subsequent preferred mode of tissue colonisation, which is essential for the choice of appropriate treatment.
Intracellular phase transitions are an emerging mechanism for cell organization. These membrane-less compartments are formed via liquid-liquid demixing and subsequent concentration of cellular components in a specific region. By undergoing these localized phase separations, cells are able to create dynamic compartments that help maintain the regulation of biomolecular interactions and localize factors such as RNAs and proteins. The utility of liquid-liquid phase transitions and the assembly of cytosolic compartments is especially critical in large, multinucleate cells. These types of cells are common in the biosphere and include skeletal muscle tissue, the placenta, many filamentous fungi, and certain types of cancer. A powerful model organism for studying liquid-liquid phase separations (LLPS) is the branching, multinucleate fungus Ashbya gossypii. In Ashbya, Whi3, a polyQ tract-containing RNA binding protein, binds to and localizes multiple different mRNAs by forming liquid-like droplets in the cell. Under normal physiological conditions, Whi3 alone cannot form liquid droplets, in vitro, as this is a phenomenon that only occurs once it is able to bind with RNA. However, how Whi3 and RNA combine to form complexes and how these complex influence droplet properties remains elusive. This issue is further complicated due to the many possible subcomplexes of RNA and protein that can exist. Specifically, the target RNAs have multiple binding sites for the Whi3 protein and the Whi3 protein has two possible RNA binding sites and likely interacts with itself using the polyQ tract. This means that there are many possible types of interactions between components. Additionally, while preliminary work within the Gladfelter Lab has shown that there is a fairly consistent average ratio of Whi3 to RNA in each droplet of a population, it is unclear what distribution of the different protein-RNA combinations gives rise to this average ratio. To understand the role that the different combinations of Whi3 and RNA play in the initial phase separation and the properties of the droplet population, such as size, shape and spatial distribution of droplets, we have created a mathematical model that employs the use of the phase field method. With phase variables to represent the volume fractions of Whi3, RNA, and the complexes they form, we propose a model that employs the use of mass action kinetics and a modified double-well free energy to represent the complex dynamics present in the cytosol. We report modelling and experimental progress on these cytosolic liquid-liquid phase separations in Ashbya as they help guide our understanding of the mechanisms behind droplet formation in the cytosol, as well as how these droplets function in cellular regulation.
Tissue engineering is a rapidly growing field, attracting a huge concentration of research effort. An important subfield of tissue engineering focuses on the use of bioreactors, devices that attempt to simulate a physiological environment in order to promote the growth of functional cell or tissue in vivo. In this talk we present a mathematical model to simulate both nutrient transport and cell proliferation within a spatially-varying permeability scaffold inside a perfusion bioreactor, and compare results from this model with experimental results from the literature.
Cancer is one of the leading causes of morbidity and mortality worldwide. Many of the disease related events, such as metastatic progression, treatment resistance, and overall survival, generate much uncertainty and are oftentimes viewed as random. Applying current modelling techniques and methodologies to existing data can produce a forecasting framework that can be leveraged to predict these significant events that impact clinical decision making. We used a retrospective, longitudinal dataset of 3,505 patients with primary bladder cancer to build forecasting models that could be used prospectively in newly diagnosed patients. We built Markov models from individual patient progression pathways and used these models to simulate and predict future locations of metastatic spread. Additionally, we used machine learning techniques to temporally predict disease progression and overall survival. Analyzing the results of these models revealed that patterns of metastatic spread emerged in distinct subgroups of patients when stratified by gender and also by pathologic stage. Additionally, analysis of patient variables showed higher associations with both recurrence and survival for pathologic staging (post-operative) as compared to clinical staging (pre-operative). Incorporating additional longitudinal data such as treatment information and genomic data could lend to the predictions of therapy resistance and side effect development.
The human species is unique in its tendency towards monogamy, a trait which can be traced back to hunter-gathering societies. Traditional explanations for this evolution have focused on the presence of paternal care and the needs of our offspring. However, recent research has challenged this claim, contending that the significant effects of mating competition on male choice result in evolutionary equilibria with little male care. This paper models the population dynamics of a rudimentary hunter-gathering society, and determines sufficient conditions for monogamy to supersede pure mating strategies as the optimal allocation of male reproductive effort. Utilising an approach based on the graphical analysis of systems of nonlinear ordinary differential equations, various situations were considered including the variation of reproduction rate, sex ratio and fertility loss. While our results concur with previous findings that male pair bonding triumphs in response to partner scarcity and male-biased populations, a dynamic expected gain model extended this reasoning to explicitly show that menopause and the age of menopause occurrence directly cause the evolutionary bias sufficient to cause a shift away from the multiple mating strategy.
Agricultural systems throughout the world are threatened by a variety of pest species, including insects, pathogens and weeds. Conventional control methods often have limited success, notably in the case of pesticide use due to the ubiquitous development of resistance. Novel control approaches are therefore needed to improve management of these pest populations. Gene drive strategies, based on the use of genetic constructs with biased inheritance to drive traits into target populations, offer promising avenues for the development of such control methods, facilitated by the advent of CRISPR-Cas molecular tools. Implementation of such strategies in agrosystems, with actively managed populations and strict socioeconomic constraints, offers a specific set of opportunities and challenges. Here we present a newly developed modelling framework design to study the outcome of these genetic control strategies in a variety of agricultural systems. This framework is built as a modular toolset that can simulate a wide range of scenarios, including a diversity of pest species, in a biologically and ecologically detailed fashion. We detail here the ability of this modelling tool to take into account specific aspects of the life history, ecology and genetics of the target organism. We introduce results demonstrating how genetic strategies can be used to control pesticide resistance. Finally, we discuss the future avenues of research offered by this theoretical framework, notably in terms of optimising implementation and deployment of genetic strategies in specific agrosystems, as well as investigating potential approaches for spatial and temporal containment of the associated genetic constructs.
The human brain consists of folds (gyri) and valleys (sulci) that vary dramatically in their size, extent, and shape across individuals. There is considerable debate among biologists as to how the folding patterns develop and if the folding patterns can be used to diagnose disease. In this presentation, I will discuss some of the mathematical and modelling approaches my research group is developing to study coritical fold formation, Turing patterns, topology, conformal mapping, and conformal invariants are some of the methods we are using to model and characterize the folding patterns of the brain in development, health, and disease. By altering various model parameters, including domain size and scaling, results from our model can be correlated with cortical folding diseases such as lissencephaly and microcephaly.
Empty homes, called Akiya in Japanese, due to an aging population and declining birthrate is very serious problem in the aging society of modern Japan. The current proportion of empty home in Japan is 13.5%, but it is expected to be increased up to 30.2% nationwide in 2033, implying that Akiya will have a serious impact on the economy as well as the national finance. Indeed, Detroit city in USA was financially collapsed by 29.3% empty homes caused by the decline of the automobile industry. Currently, the Japanese government has started a countermeasure for decreasing the empty homes, but the effect is completely uncertain. Therefore, we here develop a mathematical model including an effective utilization of empty homes and consider the future transition of empty homes. Using the model, we suggest optimal policies by which we can reduce empty homes most effectively, with regard to government tax system and subsidy for the effective utilization of empty homes. This study shows a new framework of mathematical model to be applied for solving a social and economic problem in real world.
An endoplasmic reticulum (ER) is a tubular organelle observed in cells of eukaryote including the plants and animals. Interacting with flow of actin cytoskeleton, the net-like pattern organized by ER in plant cells is continuously moving.
For the understanding of this system, we constructed a mathematical model based on a partial differential equation (PDE). By combining two spatially distributed structures, an arbitral periodic pattern by PDE and an actin mesh dependent perturbation, we successfully obtained the dynamics of the system.
In this model, it was found that the filamentous actin distribution was not necessarily needed and static perturbations were sufficient to move the pattern. Then the pattern periodicities were disturbed when the dynamics of system could be observed. Therefore, we inferred that, though it might seem paradoxical, the capacity to maintain the order and the symmetry of pattern generate the dynamics.
This theory is considered to explain one of the basic elements to generate a dynamic pattern, therefore, further discussions are needed to define the minimal conditions for giving the motility. We will comparative investigate the following points (1) pattern property effected by the perturbations and (2) distribution manners and types of perturbation. Then the periodicity and the motilities of pattern will be analyzed quantitatively.
Calcium is a universal messenger that participates in a great variety of physiological functions including muscle contraction, neuronal plasticity and immune responses. There is now compelling evidence that transient changes in the intracellular calcium concentration are key for achieving such versatility.
To date, most of these findings have been obtained from constant stimulation or step change protocols. However, under physiological conditions, cells often experience time dependent stimuli such as transient changes in neurotransmitter or oscillations in hormone concentrations. How cells transduce such dynamic stimuli into an appropriate response is an open question. We exposed HEK293 cells and astrocytes to dynamically varying time courses of carbachol and ATP, respectively, and investigated the corresponding cellular calcium activity. While single cells generally fail to follow the applied stimulation due to their intrinsic stochasticity and heterogeneity, faithful signal reconstruction is observed at the population level. We provide a number of transfer functions that translate the extracellular stimuli into the ensemble calcium spike rates. By computing the mean root square error between the predicted responses (based on the transfer functions) and the actual responses, we identify a simple leaky integrator model as a a powerful approach. For this, we show how to invert the transfer function to estimate the stimulus that should be applied in order to achieve a specific response. Throughout the analysis we pay particular attention to the non-stationarity of both the applied stimuli and the calcium responses.
We constructed a dynamical model of a salivary gland acinar cell with the objective of investigating the role of two plasma membrane (PM) anion exchangers, the Ae2 (Slc4a2) and Ae4 (Slc4a9), in primary fluid secretion. Transepithelial chloride (Cl$^-$) movement drives water transport in salivary gland acinar cells. Basolateral PM mechanisms accumulate Cl$^-$ to levels well above its electrochemical equilibrium. Following an increase in the intracellular concentration of calcium (Ca$^{2+}$), a Cl$^-$ efflux through apical Ca$^{2+}$-dependent Cl$^-$ channels (TMeM16a) generates a transepithelial osmotic gradient. Water follows this gradient by osmosis. Cl$^-$ uptake via basolateral co-transporters (Nkcc1), provides the major force for generating fluid secretion in acinar cells. Despite this, Nkcc1 knockout experiments saw an approximate 70$\%$ decrease in gland salivary rate. The residual secretion is bicarbonate (HCO$_3^-$) dependent and involves two sodium/proton (Nhe1) paired Cl$^-$/HCO$_3^-$ anion exchangers, the Ae2 and the Ae4. Experiments revealed that Ae4 knockout mice displayed a decreased gland fluid secretion ($\sim$ 30$\%$ of the control fluid flow rate) whilst Ae2 knockout mice gland salivation remained intact. The reason behind the results remains a controversial topic. Our model's results reproduce and support the experimental observations. More importantly, they suggest that the Ae4 cotransport of monovalent cations is likely to be important in establishing the osmotic gradient necessary for optimal transepithelial fluid movement.
The vector-borne Dengue fever poses a major health issue in tropic environments, which includes areas such as far-north Queensland. Historically, attempts at curtailing the spread of dengue have focused on controlling the size and spread of mosquito populations that carry the virus. Several factors make this an astronomically difficult task to accomplish on any reasonable scale, however, and so more novel methods of suppressing dengue outbreaks are being explored. One such method is the introduction of bacteria called Wolbachia into mosquito populations, which prevents mosquitoes from passing on viruses to humans. A Wolbachia invasion has strong potential to completely saturate mosquito populations due to a mechanism called cytoplasmic incompatibility. The mathematical modelling problem here becomes twofold-first, the task of inferring the position of mosquito populations, and then the modelling of Wolbachia spreading through these populations. Here, we will be discussing the effects and mechanisms of Wolbachia in more detail, including the phenomenon of cytoplasmic incompatibility. Next, recent developments in modelling mosquito populations such as the use of semi agent-based models will be outlined, within the context of predicting the spread of Wolbachia.
Tuberculosis (TB) is one of the top 10 causes of death worldwide, and the WHO’s EndTB strategy requires developing an effective vaccine. The H56 vaccine is a candidate currently in phase I/IIa trials as a boost to Bacillus Calmette-Guérin (BCG, the TB vaccine that is currently used in most countries world-wide). We build a multi-compartment, hybrid multi-scale model to 1) improve our understanding of immune response to H56 and 2) predict the role of T cells in preventing TB after vaccination.
First, we develop a two-compartment model of 31 non-linear ordinary differential equations (ODEs) that describe T-cell priming, proliferation, and differentiation in lymph nodes and blood. These ODEs allow us to track T cell response following vaccination. We calibrate our ODE model separately to human clinical trial data and non-human primate (NHP) experimental data to display differences in each species response to H56. Next, we couple our curated agent based model of granuloma formation in the lung, GranSim, to our blood and lymph node compartments to capture the host immune response to infection with Mycobacterium tuberculosis. This creation of a multi-scale and multi-compartment model allows us to represent a pseudo-whole-body response to both vaccination and infection. We use this whole-body model in the form of “virtual clinical trials” to retrospectively study the human and NHP datasets. In particular, we predict the role of T cells (induced through vaccination) throughout the course of infection in blood, lymph node, and lung granulomas. We use uncertainty and sensitivity analysis to compare and contrast immune response to vaccination and infection in NHPs and humans. We conjecture that vaccine dose may be critically important when evaluating H56 efficacy against Mycobacterium tuberculosis infection.
The mammalian circadian rhythm is governed by the principal clock which is located in the suprachiasmatic nucleus (SCN). This clock is composed of tens of thousands of neurons and their connection to one another is essential to important roles of SCN: synchronization, entrainment to light, etc. Previous studies about the SCN network used the maximal information coefficient (MIC) statistic to oscillating time course data during resynchronization after desynchronization by TTX. Since this method requires desynchronization and resynchronization, it takes a long time and the recovered structure may not be similar to the normal SCN network. Finally, even though numerous studies argued that the connection between SCN cells are asymmetric, MIC cannot detect directionality.
We develop our own method, which can detect causality. This method does not require desynchronization by TTX, so it takes shorter time and one can obtain the normal network structure of the SCN. Moreover, it is able to obtain main results of previous studies: small-world network, exponential distribution of node degree. Finally, as our method can diagnose causality, our method can also be used to verify the asymmetry of the SCN.
A mathematical standard structure of a binary digit of memory in a cell is presented. This is based on a kind of frequency model with scale effect. This model has ability of "on-off" switching property, and moreover, this is affected by scale effect to make the memorable ability be reinforced. This property is derived from multiple covalent modification cites inducing important enzyme reaction, which is represented by the Michaelis-Menten type nonlinearity.
Some examples of application is also presented. One example is Cyanobacterial allosteric circadian rhythm of Synechococcus. Another is about Biological Nitrogen Fixation ability of Nostochineae cyanobacteria. I will explain briefly that scale effect affects the phenomena effectively to make the system be robust.
Physiological levels of oxygen and carbon dioxide in the blood are tightly regulated by varying the pattern of breathing, but this can be achieved with different combinations of amplitude and frequency. Why a specific combination of amplitude and frequency of breathing is observed remains a mystery. The aim of this study is to explore the hypothesis that the particular combination realised is optimal with respect to some objective function. Several objective functions have been suggested in the literature, such as the rate of work during inhalation, the average force exerted by the respiratory muscles, and the weighted sum of volumetric acceleration and work during inhalation; all of these objective functions provide physiologically acceptable minima under normal conditions. Resolving this issue requires optimal solutions of mathematical models that reflect more accurately the complex interaction between lung mechanics and gas exchange, but this in turn requires the development of new computational methodologies. To help achieve this goal, we constructed a simple mathematical model, consisting of two piecewise linear differential equations, that mimics gas exchange in the lungs. By using concepts from optimal control theory, we found the necessary conditions that minimise a given objective function subject to several constraints, such as satisfying the differential equations and maintaining one of the variables at a given average value. We could then solve the optimal control problem both analytically and numerically. Our method can be extended to models with higher dimensions.
Barreto et al. (2017) studied the rock-paper-scissors games (or cyclic competition games) by the dynamics of phenotypic and genotypic frequencies corresponding to three morphs of lizards. In their model they show that both dynamics have equal internal equilibrium but the genotypic model has wider parameter range for its stability compared to the phenotypic model.
Here, first of all, we investigate phenotypic and genotypic competition models with two-species. Then we consider the Barreto et al. model with different genotypic correspondence to three phenotypes. Lastly computer simulations for these models are carried out on lattice space and the effects of spatial structure are discussed.
The apical-basal cell polarity endows epithelial tissues with distinct biochemical and mechanical properties at the apical, basal, and lateral sides of tissues. Mounting evidence shows that remodelling of epithelial cell polarity, where polarity regulators redistribute their locations on the cell membrane, causes changes of cell shapes that drives tissue morphogenesis during embryonic development. However, the principles underlying epithelial tissue morphology formation regulated by remodelling of intracellular polarity remain elusive. Using mathematical modelling and computer simulations, we first show how epithelial cell polarity is established, maintained, and remodelled by the interactions among the polarity regulators in the model system of fruit fly (Drosophila melanogaster). Then, based on a geometrical modelling for epithelial tissue deformation in three dimensions, we further demonstrate how tissue morphology is determined by its constituent cell mechanics. Finally, taking into account the polarity-dependent cell mechanics, we show that the feedback between cellular mechanics and polarity can lead to a self-organizing tissue morphogenesis in the absence of sustained external stimuli.
The activated sludge process (ASP) is the most widely used process for the biological treatment of both domestic and industrial wastewaters. Wastewater treatment plants (WWTPs) based on the ASP are in widespread use in both developed and developing countries.
The ASP uses microorganisms which grow by consuming organic pollutants that are present in the wastewater. This produces new organisms whilst simultaneously cleaning the wastewater. The microorganisms flocculate to form settleable solids, the `activated sludge'.
Central to the success of the ASP is the activated sludge. However, a significant drawback of the ASP is the production of excess `sludge'. The expense for treating this can account for 50--60% of the total operating costs in a WWTP. Traditional methods for disposing of excess sludge are becoming increasingly unattractive due to a combination of increasing land costs and environmental concerns about the presence of toxic elements. Thus there is a growing interest in methods that reduce the volume of excess sludge produced by the ASP.
A promising method to reduce sludge formation is to introduce a higher order organism, such as protozoa or metazoa, into the system which predate upon the microorganisms. This is potentially very attractive, since once the predators have been released into the reactor there are no `running' costs.
Why does predation reduce sludge? The ASP can be considered to be a food-chain, in which the biomass extracts mass and energy from the substrate. The introduction of a predator introduces a new layer into the food-chain: mass and energy are now transferred from the microorganisms to the predator. At each step in the food chain not all of the available energy and material are transferred to the next level: some energy, a significant proportion of the energy, is used for maintenance processes, respiration and reproduction. Thus predation on microorganisms may lead to a lower total biomass, i.e. sludge reduction.
Predation has been shown to be an effective technique in lab-scale experiments and pilot-scale systems. Although a variety of predators could be used, much attention has focused on the use of worms. Worm growth is clearly a prerequisite for sludge reduction through predation. Relatively little is known about the growth and development of worms during sludge predation. However, it has been shown that the wrong choice of aeration rate, temperature and predator (worm) density can adversely effect worm growth and consequently sludge reduction.
One of the main barriers to the adoption of predation as a cost-effective mechanism to reduce sludge formation is the phenomenon of uncontrollable growth of predators. Predator density in full scale plants can often reach very high densities. Associated with this is the well-known phenomenon of worm blooms, in which predator population densities[ display peaks followed by a sudden disappearance of the population. The development of methods to control worm proliferation is a challenging problem that needs to be overcome.
In this presentation a simplified model for a wastewater treatment plant is extended to include predators. Mechanisms which lead to the unsuccessful introduction of predators are identified.
Species eradications and reintroductions are drastic management actions that alter ecosystem structure with the end goal of conserving or restoring one or more species of interest. Because of the complexity of ecosystem food webs, these actions can have unintended consequences on other parts of the ecosystem. Even with the best available data, predicting future trajectories of the ecosystem with high precision is not always possible. In this talk, I demonstrate the situations where ordinary differential equation (ODE) models of the ecosystem food web, calibrated to species abundance data, can or cannot be used to accurately predict future responses of the ecosystem to species eradications and/or reintroductions.
The placement of stomata on the leaf surface is regulated by a number of things: coordinated asymmetric cellular division, genetic regulation, and distribution of extracellular signalling molecules.
It has been noted that the number of stomata per non-stomatal epidermal cell can vary between plants which have the same spatial distribution of stomata, which suggests that in some circumstances one method of maintaining stomatal placement (increased cell division, increased cell expansion, increased recruitment of cells to the stomatal fate) may be optimal over another.
Perhaps the reason for this difference in stomatal placement regulation is due to difference in the energetic costs to the plant of each of these methods. Using a combination of experimental and mathematical biology, the relationship between energetic cost and stomatal placement in the model organism Arabidopsis thaliana is being investigated.
Since many kinds of Direct Acting Antivirals (DAA) have become the main treatment instead of interferon-α (IFN-α) against hepatitis C virus (HCV), combinations of these DAA are now standard treatment strategies. These treatments are very effective, however, at the same time, this provokes the problems which drug combination is more or most effective for each patient. Therefore, we established a viral replicon system for HCV to quantify the efficacy of each drug or drug combination easily. In this study, we used 15 different kinds of DAA (telaprevir, danoprevir, asunaprevir, simeprevir, sofosbuvir, VX-222, dasabuvir, nesbuvir, tegobuvir, daclatasvir, ledipasvir, IFN-α, IFN-λ, cyclosporin A, and SCY-635) and quantified the efficacy of (single and) multiple drug combinations using the replicon system. Then, we calculated instantaneous inhibitory potential (IIP) which captures the effectiveness of each drug combination as a function of the drug concentration. Furthermore, we defined and calculated the critical concentration index (RCI) that achieves 95% reduction in HCV replication is shown for each drug combination. The RCI intrinsically varies among the drug combinations. As an example, we found SMV&IFN-α yielded the lowest RCI values among the tested double-combinations. We would like to further discuss how our framework could enable drug usage optimization by quantifying the antiviral activity in preclinical settings, and presenting basic evidences for cost-effective drug selection.
The unprecedented accumulation of high-throughput single-cell RNAseq data provides an essential opportunity for researchers to study the dynamics of biological systems. At the same time, it raises new questions and grand challenges in view of the key characteristics of these data that include high-dimensionality and heterogeneity. In response, we develop a data-driven computational framework to map and reconstruct of dynamical trajectories of cell based on the high-dimensional molecular profiles from snapshot single-cell data. By integrating tools from topological data analysis (level-set method) and the nonparametric statistics, we first perform nonlinear dimension reduction on the data to construct the cell developmental landscape, and then reconstruct the complex trajectories by finding the cell state transition path on the landscape. We apply this method in the analysis of cell developmental processes and provide novel insights into the field of mathematical and computational biology.
To gain deep insights into the interactions among major types of cells, such as healthy, tumour, and immune cells, in a tumour-bearing host, we propose a mathematical model consisting of four delay differential equations by considering the inherent delay in the tumour cells’ life cycle. When choosing model parameters, sensitivity analysis justified previous assumptions had been considered.
Since the impact of immune system is decisive on the regression of tumour cells, stability analysis is performed for two different levels of immune system surveillance, that is, healthy and compromised tumour-bearing host.
Theoretically, it has been showed that in a healthy host the immune system is able to eradicate the tumour cells completely whereas in presence of weak immune system tumour cells have the opportunity to keep proliferating.
The latter case, necessitates an efficient treatment, and more precisely a combination of different treatments, for being able to fight and demise tumour cells.
Nowadays, Chemotherapy is the most common treatment which is administrated in various protocols. Among them, MTD (Maximum Tolerated Dosage) and Metronomic are of great significance in clinical observations. Our model suggests that Metronomic technique surpasses MTD after introducing this treatment. The last but not the least, since the results of mathematical experiments illustrate the administrated Chemotherapy protocols are not able to eradicate the tumour cells, the proposed model unveils the demand for an adjunct treatment.
One of the challenging questions in Ecology is how spatial structure influences the formation of biodiversity patterns. Here, we explore networks that represent population connectivity including cases where dispersal is symmetric and asymmetric, looking for the relationship between population connectivity network structure and biodiversity patterns therein. In this metapopulation system, we simulate dispersal of individuals based on the probabilities given by a dispersal connectivity matrix, and measure resulting biodiversity.
Difficulties with this approach include the variety of possible network structures and the increasing complexity of networks as they become larger. To deal with this problem, we analyze and compare 1) theoretical networks of relatively small size, which vary in complexity, and 2) realistic large and complex networks that represents coral larval dispersal connectivities in a marine system. Through this approach, we will be better able to understand the possible contribution of geographical structure to biodiversity for metapopulations.
There are large differences between men and women of child-bearing age in the expression level of 5 key enzymes in one-carbon metabolism almost certainly caused by the sex hormones. These male-female differences in one-carbon metabolism are greatly accentuated in during pregnancy. Thus, understanding the origin and consequences of sex differences in one-carbon metabolism is important for precision medicine.
We have created a mathematical model of hepatic one-carbon metabolism, based on the underlying physiology and biochemistry. We use the model to investigate the consequences of sex differences in gene expression. We use the model to give a mechanistic understanding of observed concentration differences in one-carbon metabolism and explain why women have lower S-andenosylmethionine, lower homocysteine, and higher choline and betaine. We give a new explanation of the well known phenomenon that folate supplementation lowers homocysteine and we show how to use the model to investigate the effects of vitamin deficiencies, gene polymorphisms, and nutrient input changes.
Our model of hepatic one-carbon metabolism is a useful platform for investigating the mechanistic reasons that underlie known associations between metabolites. In particular, we explain how gene expression differences lead to metabolic differences between males and females.
Many species of animals, plants, cyanobacteria, and fungi are reported to have circadian clocks, that is, gene regulatory and biochemical network in a cell exhibiting endogenous oscillatory dynamics with a period close to 24 hours. There is diversity in the members of the circadian system, but all clocks share similar features such as temperature-compensated circadian period. Therefore, theoretical studies on a general (and abstract) network are as important as studies focusing on a specific network found in each species. Previous simulation studies using evolutionary algorithm have established a method to obtain networks that follow an autonomously oscillating dynamics. In those simulations the degree of deviation from the demanded period was defined as the cost function to be minimized by the optimization process because natural selection apparently favours a clock with a certain period, namely, 24 hours. However, if the ultimate goal is to implement the 24-hour period, it may be easier and more reliable to refer to external cues such as sunlight. Indeed, the circadian clocks are reported to have sensitivity to the light signal, and their phases are adjusted by sunlight. In addition, empirical studies have suggested that the circadian clocks play more sophisticated roles such as predicting timing of the next dawn (respectively, the next dusk) and preparing molecules required during the day (resp. during night) in advance. In the present study, we first examine the effect of sunlight as an external oscillator represented by a binary function (1 = light; 0 = dark). We compare the proportions of regulatory networks that exhibit and do not exhibit an oscillation autonomously (i.e., without the external force) among those oscillating with the same period as the external force. We also found that even the non-autonomous oscillators have an adaptive feature of the day-length compensation, showing 24-hour oscillation regardless of the proportion of the light period. Next, we introduce a novel cost function that reflects the merit of predicting future, and report shared features of networks optimized based on this cost function.
Human papillomavirus, or HPV, is a sexually transmissible virus infection, which is necessary risk factor for developing cervical cancer, most common type of cancer in working age women in Moldova. We observe both behavioural change (increase in sexual partner acquisition rates) and demographical change (population ageing and massive emigration, but still very young), which both corresponding to second demographic transition since Soviet Union collapse (yearly expenditure on health limited to 150 EUR-per-capita). Moldova will spend around 400,000 EUR on 'single cohort' vaccination in 2018, which cost effectiveness is questionable, because vaccinating a single cohort may not have a substantial effect in other countries. Thus we examine such a single vaccination scenario to show its conditional cost-effectiveness.
We have run computer simulation to prepare cost-benefit/effectivness analysis for different vaccination strategies, various screening programs and preventive programs for Moldova in low resource settings, based on its own demography and sexual behaviour. We used data since 1998 to 2017 to adjust model parameter and we project till around 2038. Model aggregated the most important paths of infection, cancer development, uncertainty in healthcare capacity and sexuality and prevention scenarios with around 100 differential equations (stochasticity introduced in sexual partner change rates).
Single cohort vaccination could be both cost-beneficial (total cost reduction balance intervention cost before 2037) and cost-effective (with incremental impact in 20 years perspective on the level of 2300 EUR/QALY).
The possible explanation of this nonintuitive behaviour is transitional situation in Moldova ($R_0$~1), still small change of conditions could cause strong effect in epidemiology. Main effect of intervention is via men, which avoid infection and will not infect other women. This can have effect probably while changing partners is still not as common as in other countries. However, even slight change in initial conditions and parameter values could diminish positive effect (e.g. faster partner acquisition rate increase).
Gene expression and genetic regulatory networks (GRNs) have become important areas of study because they play an important role in solving issues related to human health such as cancer diseases. Cancers can be hard to recognise and the best way to improve their defected is to understand the underlying mechanisms of genetic regulatory dynamics. This means introducing and improving mathematical models of GRNs to reflect oscillatory phenomena. In this work, we propose a lac operon model with time delay and nonlinear degradation rate for mRNA and investigate the nonlinear dynamical behaviour arising from the model, such as stability and Hopf bifurcation of the equilibrium, by taking the average delay as the bifurcation parameter. From the dynamical behaviour of this system: the equilibrium is stable in the absence of delay or when the delay is small; as time delay increases gradually then the stability changes and Hopf bifurcation happens. Then reverse a Hopf bifurcation to translate the equilibrium of the system from unstable to stable steady state.
The spatio-temporal intra- and interspecific competition of two diffusing similar populations is considered. The growth of both populations is either logistic or shows an Allee effect. Conditions of spatial segregation without mixing are investigated. Furthermore, the impact of density-dependent environmental noise on the occurring stationary fronts is studied. A special focus is set on the development of functions describing the density-dependent noise intensity. The obtained results are associated with a biological case study related to the competition of two invasive weeds.
Fluorescent cell cycle labelling in cell biology experiments provides real time information about the location of individual cells, as well as the phase of the cell cycle of individual cells. We develop a stochastic, lattice-based random walk model of a two-dimensional scratch assay where the total population is composed of three distinct subpopulations which we visualise as red, yellow and green subpopulations. Our model mimics FUCCI technology in which cells in the G1 phase of the cell cycle fluoresce red, cells in the early S phase fluoresce yellow, and cells in the S/G2/M phase fluoresce green. The model is an exclusion process so that any potential motility or proliferation event that would place an agent on an occupied site is aborted. Using experimental images and measurements we explain how to parameterise the stochastic model, and we apply the stochastic model to simulate a scratch assay performed with a human melanoma cell line. We obtain additional mathematical insight by deriving an approximate partial differential equation (PDE) description of the stochastic model, which leads to a novel system of three coupled nonlinear reaction diffusion equations. Comparing averaged simulation data with the solution of the continuum limit description shows that the PDE description is accurate for biologically-relevant parameter combinations.
Using a simple lattice model describing infectious tree disease dynamics on a homogeneous landscape we can observe an Epiphytotic phase transition from local confinement of the pathogen to a global epidemic through the forest. The phase transition can be understood in terms of the forest tree density and the pathogen virulence. One interesting application of the model involves capturing the pathogen spread-velocity and applying simple statics to pre-empt the phase transition into the Epiphytotic regime. Theoretically this yields an early warning system. The early warning phase transition could potentially be utilised to inform control strategies and negate system wide devastation of plant species.
Pattern formation by Delta-Notch interaction has been well studied experimentally and theoretically. The Delta-Notch system is observed in various pattern formation process such as somitogenesis, neuroendocrine cell differentiation in lung, T cell differentiation and blood vessel development. Recent studies have shown endothelial cell proliferation and movement happen during this process. However, to our knowledge, there is little theoretical research about the effect of cell dynamics on Delta-Notch interaction.
In the present study, we examined the effect of cell dynamics on Delta-Notch pattern formation during retina vasculature development. We incorporated cell movement and proliferation to the model for Delta-Notch interaction. Using the model, we analytically derived the instability condition and numerically generated the patterns which have the similarity with the three patterns observed in vivo. It is difficult to capture the dynamics of cell movement and proliferation with standard linear stability analysis of fastest growing wavenumber component. Therefore, to consider all wavenumber components, we introduced the instability index, $\Psi(t)$, as the mean of square of Delta expression values. Based on Parseval's theorem regarding discrete Fourier analysis, we can derive that $\Psi(t)$ is equivalent to the average of power spectrum. Therefore, by considering the values of power spectra, we can evaluate the instability of the system.
$$\displaystyle {D}_x =\sum_{m=0}^{n-1} \delta_k \mathrm{e}^{ \lambda_k t + i k x} \\ \displaystyle \Psi(t) := \frac{1}{n} \sum_{x=1}^n {D_x}^2 = \frac{1}{n} \sum_{m=0}^{n-1} |\delta_{\frac{2\pi m}{n}}|^2 $$ $D_x$ is the expression values of Delta, $\delta_k$ is the discrete Fourier transform of $D_x$, $n$ is the number of the cells, $k=\frac{2\pi m}{n}$.
These analyses and numerical calculations suggest that the vasculature which express homogeneous pattern shows high motility and proliferation rate of their endothelial cells. Based on these theoretical results, we experimentally observed cell dynamics during retina vasculature development by organ culture and immunohistochemistry. The results showed random endothelial cell movements and proliferations happened more frequently in vein than in artery, which are consistent with analytical and numerical results. These results suggest that cell dynamics affect artery-vein differentiation via Delta-Notch interaction.
Seasonality and contact patterns due to environmental fluctuations and social behaviour affect the dynamics of disease outbreaks. Recent studies applied to deterministic epidemic models with periodic environments have shown that the average basic reproduction number is not sufficient to predict an outbreak. We extend these studies to stochastic epidemic models with periodic environments to investigate the combined effect of periodicity and variability on disease outbreaks. The deterministic models are extended to continuous-time Markov chain and stochastic differential equations. A numerical study of the dynamics of several stochastic SIR and vector-host models with environmental variability and periodicity are investigated in terms of probability of an outbreak.
Tethered enzymatic reactions are a key component in signalling transduction pathways. It is found that many surface receptors rely on the tethering of cytoplasmic kinase to initiate and integrate signalling. A key factor to such reaction is the molecular reach; however, the role of it is incompletely understood. To date, a large number of compartment-based ODE and stochastic models have been developed to study this problem. In recent years, spatial-stochastic models have emerged as a more realistic representation for such processes, among which lattice-based stochastic reaction-diffusion models are a popular approach for studying complex spatio-temporal processes inside cells. To understand the role of molecular reach in tethered signalling, we employed an accurate and convergent lattice-based stochastic reaction-diffusion model (CRDME). We find that the molecular reach can increase or decrease biochemical reactions depending on the diffusion coefficient in 2D membrane but not in 3D cytosol.
Negative feedback and buffering, the use of reservoirs of molecules to maintain concentrations of key molecular species, are the primary known mechanisms for robust homeostatic regulation. However, the fundamental principles behind their combined effect have not been elucidated. Here, I will present our recent research on the interplay between buffering and negative feedback in the context of cellular homeostasis. This research shows that negative feedback counteracts slow-changing disturbances, whereas buffering counteracts fast-changing disturbances. Furthermore, feedback and buffering have limitations that create trade-offs for regulation: instability in the case of feedback and molecular noise in the case of buffering. However, because buffering stabilizes feedback and feedback attenuates noise from slower-acting buffering, their combined effect on homeostasis can be synergistic. These effects are consistent with experimental observations of both ATP homeostasis and pH regulation in vivo. The methodology is based on control theory and complexity-aware minimal modelling. The discussed principles are crucial for studying robustness and homeostasis in biology and biotechnology.
In this study, the maximal extent of future NPP uncertainties are explored by employing the conditional nonlinear optimal perturbation related to parameters (CNOP-P) approach and the Lund-Potsdam-Jena (LPJ) model based on future climate change assessments, which are provided by 10 general circulation models (GCMs) of the Coupled Model Intercomparison Project 5 (CMIP5) under the Representative Concentration Pathway (RCP) 4.5 scenario at the North-South Transect of Eastern China (NSTEC). We find that the future NPP will increase due to climate change and CO2; however, there is a difference in the extent of the variation resulting from the 10 GCMs and the CNOP-P approach. Future NPPs are estimated from 3.89 Gt C (MRI-CGCM3 model) to 4.51 Gt C (bcc-csm1-1 model) using the LPJ model driven by the outputs of 10 GCMs. The estimates of NPP with two types of CNOP-P-type climate change scenarios are 4.74 Gt C and 5.31 Gt C and are larger than other estimates of NPP. The above results imply that the terrestrial ecosystem supplies possible conditions for future carbon sinks for all climate change scenarios, especially for the CNOP-P-type climate change scenarios, although the estimates remain uncertain. Stimulative photosynthesis due to high precipitation and restrained autotrophic respiration due to low temperatures may play important roles in the carbon sink due to the CNOP-P-type climate change and CO2 in all climate change scenarios. In addition, it is found that the combination of climate change and CO2 is also a key factor in the increase of NPP.
Asymmetric cell division is one of the widespread mechanisms for generating cell diversity, for which a mother cell creates a polarity in both membrane and cytosol. In both experiment and theoretical approaches, PAR polarity of C. elegans embryo has been extensively well-studied and it was found that Anterior-Posterior (AP) polarity of cell membrane proteins plays a crucial role in determining cell asymmetry. However, most of previous studies have not considered the role of cytoplasmic proteins on AP polarity formation, although AP polarity is occurred with a tight regulation of both membrane and cytoplasmic proteins such as MEX-5/6 and PAR- 5. Here, we develop a multi-dimensional polarity model including cell’s geometrical property and show how the cytoplasmic protein plays an important role in creating a robust AP polarity. We also show that the cell geometry can give a critical effect on AP polarity, and the temporal and spatial regulation for the robust AP polarity is based on the harmony of biochemical, mechanical, and cell geometric properties.
Tuberculosis remains a widespread and deadly disease, infecting approximately two billion people worldwide. Gaining a better understanding of the immune response to Mycobacterium tuberculosis is crucial, and the increased prevalence of multi-drug resistant strains, the current complexity and length of treatment, and the inherent difficulties of experimental work each highlight the need for new approaches. Computational modelling of the complex immune response which results in the formation of lung granulomas can enable analysis of what is currently a relatively black box for scientists. In particular, the role of neutrophils in the granuloma response to M. tuberculosis infection remains largely unknown, as neutrophils are easily activated and short-lived, and thus pose unique experimental challenges to wet lab study. Through the incorporation of a neutrophil cell type into our existing hybrid agent-based computational model, we investigate the spatiotemporal dynamic formation of lung granulomas in response to M. tuberculosis infection. We are interested in determining how neutrophil presence induces both global and localized features of granulomas, the dynamics of neutrophil interactions with other cells and M. tuberculosis, whether the neutrophil response to M. tuberculosis is beneficial or detrimental to the host, and how that response can be harnessed to aid in therapies. Our goal is to not only identify the sufficient biological assumptions necessary to reproduce experimental datasets from our collaborators, but to shed insights on the mechanistic basis for neutrophil-directed immunopathogenesis in M. tuberculosis that are beyond current experimental capabilities.
In silico prediction of the relationships between the protein structure and its physiological activity is an important research topic for drug design. Broad picture of my research is to construct a topological model to clarify the antibody-antigen recognition system, since immunotherapy is applied to wide range of severe diseases such as cancer [1].
To achieve the goal, we focus on a topological method called “Fatgraph models of proteins” [2]. Fatgraph models of proteins are topological two-manifold with boundary components (surface) which have one to one correspondence with three-dimensional protein structures listed on Protein Data Bank (PDB) [3] with only a few exceptions. The traits of each surface for each protein are described by the following invariants; Euler characteristics, number of boundary components, and genus. Fatgraph is also called ribbon graph and has already proven their utility in theoretical physics including string theory [4].
In this research, we constructed fatgraph models of antibodies and investigated the traits of antigen-binding fragments (Fab). Then, we topologically examined the transformations of the fatgraphs of the proteins due to the changes of protein sequences or existence of ligands.
[1] Patel, Shetal A. et al. (2018) Combination Cancer Therapy with Immune Checkpoint Blockade: Mechanisms and Strategies, Immunity, Volume 48, Issue 3, 417 - 433
[2] R. C. Penner, et al., (2010) Fatgraph models of proteins, Communications on Pure and Applied Mathematics. Volume 63 , Issue 10 , 1249 - 1297.
[3] H.M. Berman, et al. (2000) The Protein Data Bank, Nucleic Acids Research, 28: 235 - 242. http://www.rcsb.org/
[4] R. C. Penner, (2016) Moduli spaces and macromolecules. Bull. Amer. Math. Soc. 53, 217 - 268.
Infection with the bacteria Mycobacterium tuberculosis (Mtb) causes a dynamic immune response that spans multiple organs across multiple time scales. The infection may start with just a few Mtb infecting a few macrophages within the lungs, but sites of infection quickly experience local immune responses, and soon after draining lymph nodes become involved to mount an even greater defense. Different mathematical tools are required to simulate this complicated immune response, from the recruitment of the cells from the lymph nodes and their transport through the blood to the sites of infection, to the formation of granulomas, which are the hallmark of Mtb infection and are comprised of immune cells, bacteria, and other molecules. Agent-based modelling has been an effective multiscale approach to simulate the immune response during infection that forms a single granuloma. However, the computational cost of running these simulations begs for a quicker, coarser approach in order to study an infection across multiple granulomas within a single host lung. We propose a novel hybrid agent-based model that uses ordinary differential equations to evolve populations of cells, Mtb, and other molecules within each individual granuloma, that is then simulated as an agent within a whole-lung in silico infection model. The lung compartment can be connected to lymph nodes with granuloma agents spreading between compartments. This will allow for a higher-scale host response simulation to allow us to better address questions at the whole host scale and also be able to calibrate to datasets of whole human hosts.
Recently, it has been shown that intracellular transport of assembled intermediate filament proteins is one major determinant of their organization in cells. Based on experimental data, mathematical models of the spatio-temporal distribution of intermediate filaments in cells are developed to investigate the contributions of different types of transport such as retrograde flow of actin and motor proteins. Furthermore, models for the motion of single filaments driven by motor proteins are also proposed.
It is an interesting topic to predict how dynamics change in disease progression. Our goal is to predict dynamics of the vital cells in the human immune system using the hierarchical HIV models. We assumed that the parameters are random variables with intra-level noise and inter-level noise. Individual-level parameters of HIV model are estimated by generalized least squares method based on the field data. At the population-level, it is not desirable to estimate joint distribution of all the parameters because of computational cost and overfitting. Thus, we build a model by estimating the distribution of only a few selected by the parameter reduction using the sensitivity matrices. Then, we estimate the joint probability distribution of remaining parameters employing algorithms for population parameter estimation. Finally, we solve ordinary differential equations with random coefficients to predict the dynamics of target cells.
Predicting the carbon balance of terrestrial plants, and its vulnerability to environmental change, is a fundamental problem common to agriculture, ecology and ecosystem science. Plant productivity is fuelled by the carbon taken up during photosynthesis, but other limitations to growth may restrict the ability of plants to utilise this carbon in production, a phenomenon known as sink limitation. Understanding how sink limitation affects the carbon balance of plants is one of the key challenges to developing better predictions of productivity. Here, we applied data assimilation (DA) to an experiment in which sink limitation was induced by restriction of the rooting volume of Eucalyptus tereticornis seedlings, in order to infer how carbon balance processes were affected. The ultimate goal of our study was to examine how carbon balance models should be modified to represent sink limitation of growth, whilst maintaining mass balance.
Our results demonstrate that several process representations need to be modified, including a clear need to incorporate a temporary storage pool of carbon as non-structural carbohydrate (NSC), with a dynamic utilization rate for growth. We were able to infer that, in addition to a reduction in photosynthetic rates, sink limitation reduced the NSC utilization rate, increased growth respiration, modified the carbon allocation pattern and accelerated senescence. The attribution analysis indicated that all of these process responses contributed significantly to the overall reduction in biomass observed under low rooting volume. Our DA-model analysis of this root volume restriction experiment provided significant new insights in the response of key carbon balance processes to sink limitation. Applying this approach more broadly would potentially allow us to identify general patterns in these responses that could then be formulated for inclusion into models. Overall, this approach provides important insights into the relationship between carbon uptake and plant growth, and could significantly advance our models of vegetation responses to global change.
The Proliferation-Invasion (PI) mathematical model of patient specific glioblastoma (GBM) growth utilizes T1 and T2-weighted/fluid attenuated inversion recovery (FLAIR) magnetic resonance (MR) images to estimate net proliferation (ρ
) and net invasion (D) rates. We have previously developed methods to parametrize this model from these routine MRIs such that higher D/ρ
tumours are considered more invasive than lower D/ρ
tumours. More advanced MR imaging methods such as diffusion-weighted imaging (DWI) allow for the calculation of quantitative measurements. The apparent diffusion coefficient (ADC), which is calculated from DWI, is believed to be predictive of tumour cellularity and microstructure. Our objective is to validate that the PI model’s measure of invasiveness D/ρ
derived from routine MRIs is consistent with quantitative MRI measurements. We hypothesize that D/ρ
will be correlated with the distributions of ADC values within the tumour. Further, we expect that the proportion of high ADC values values (low cellularity) will predominate for low D/ρ
(highly invasive) tumours. T1Gd and FLAIR images were segmented and regions of interest (ROIs) created for patients with contrast-enhancing GBMs. A FLAIR penumbra ROI was created by excluding the T1Gd ROI from the corresponding FLAIR ROI. ROI’s were then used to mask the co-registered ADC maps and the ADC values from within the ROI were plotted as a histogram. The ADC histogram from the FLAIR ROI was fit using a bimodal Gaussian model. ADC histogram peak boundaries were calculated as being +/- one averaged standard deviation from the averaged peak location (FLAIR ROI histograms). The averaged boundaries were then applied to each histogram. We calculated the percent of ADC voxels classified as being below each peak, within the lower or upper peak, within both peaks, or above the peaks. The percentage of voxels within each region were then plotted against D/ρ
. Understanding the relationship between D/ρ
and ADC allows us to connect observations on multiparametric imaging and elucidate the tumour biology. These results support the practical applicability of the PI mathematical model in quantifying patient-specific invasion characteristics by cross-correlating those findings with that seen on other imaging techniques and parameters, in this case ADC.
The neutral theory in community ecology attempts to explain its species-abundance composition by the balance between local extinction of species due to demographic stochasticity and immigrational supply of new species from a meta-community. The original neutral theory is described by two parameters; the one (“theta”) that represents biodiversity in the meta-community, and the other (m) is the immigration rate. Quite surprisingly, such a two-parameter model often shows a good fit to real community data, leading to the impression that many real communities significantly look “neutral”. Inspired by the methodology in evolutionary game theory, here I try to reveal whether non-neutral communities with intra- and inter-specific interactions really look “neutral”, and if not, under what condition. My random community model explicitly incorporates those interactions and has new parameters that control the proportions of mutualistic (+/+), competitive (-/-), and exploitative (+/-) interactions in the community. Extensive computer simulations have revealed that communities that are rich with mutualistic interactions, or those with exploitative ones tend to look neutral although they are not, whereas those communities that have a mixture of competitive and exploitative interactions are often deemed non-neutral. Those results tell us limited statistical power in existing neutrality tests.
Error threshold being the forefront of the issue for computationally-intensive methodologies and statistical models based on Kingman’s coalescent, six main points arise: (i) discrepancy between the exact and linearized non-coalescence probability; (ii) validity of the linearized coalescence probability; (iii) conditional probabilities of single-pair and multi-coalescences given at least one coalescence; (iv) parity of reduced ancestral processes that suppress multi-coalescences, when compared to the exact ancestral process; (v) genealogical topology; and (vi) subsequent inter-arrival times. Regions of adherence and detraction from the Wright-Fisher ancestral process were identified, in terms of transition probabilities and expected inter-arrival times, due to the linearization of Kingman’s coalescent that neglects multi-coalescence events. Empirically, expected genealogical parity of the single-pair restricted Wright-Fisher haploid model exceeds 99% where $n\leq \frac{1}{2} \sqrt[3] N$; similarly, per expected interval where $n\le {1\over 2} \sqrt{N\over 6}$. Quantitative parametric analysis on the varieties of multi-coalescence events in restricted Wright-Fisher models, rather than stochastic realization, avoids elaboration of the genealogical sample space under multi-coalescence. Mutational activity, such as site frequency spectra, will not be accurately measured with the Kingman coalescent when sample sizes exceed these criteria.
Infectious diseases affect individuals (immunology) and populations (epidemiology). While these two scales of infection are intimately linked, the vast majority of studies of infectious diseases ignore or greatly simplify the effects of the other scale. As a result, public health programs can be ill-informed. Mathematical models that link the in-host and population scales of infection can better inform public health programs. However, these models can become very complex very quickly, even in their most simplified forms. We focus our mathematical modelling studies on the effects of immunity, the key outcome of infection at the in-host level (development of immune memory), and the key indicator of pathogen spread at the population level (through susceptibility and transmissibility). Immuno-epidemiological models that are well-informed by our in-host and population level studies are then developed. In this talk I will discuss our immunological, epidemiological, and immuno-epidemiolgical modelling studies of influenza, HIV, measles and pertussis. The effects of vaccination and waning immunity will be highlighted.
Pattern formation is an elegant system in developmental processes and spatial patterning at different scales (e.g., nucleus, cells, and tissues) plays important roles in determining the function or fate of cell and/or tissue. In particular, cell polarity is used for regulating cell migration, cell aggregation, or cell functions by creating a spatial patterning. In this minisymposium, we focus on spatial patterning and cell polarity that occurs at several scales in cells and tissue with both theoretical and experimental approaches, and show various mathematical approaches to understand pattern formation in the life sciences.
Anterior-posterior (AP) polarity formation of cell membrane proteins plays a crucial role in determining cell asymmetry, which depends not only on the several genetic process but also biochemical and biophysical interactions. In Caenorhabditis elegans, a single fertilized egg cell (P0), its daughter cell (P1), and the germline precursors (P2 and P3 cells) form two exclusive domains of PAR proteins on the membrane along the anterior-posterior axis. Since a mother cell divides into two dissimilar daughter cells with different properties using this polarity, the shape and length scale of PAR domains are critical. Furthermore, the phenomenon of polarity reversal has been observed in which the axis of asymmetric cell division of the P2 and P3 cells is formed in an opposite manner to that of the P0 and P1 cells. The extracellular signal MES-1/SRC-1 has been shown to induce polarity reversal, but the detailed mechanism remains elusive. Here, I explore the mechanism of symmetry breaking and AP polarity formation using self-recruitment model of posterior proteins and show how the shape, length and location of PAR polarity can be form robustly. If time is allowed, I also introduce a multi-dimensional polarity model including cell's geometrical property and show how the cell geometry can give a critical effect on AP polarity.
Many cells within epithelial tissues display polarity along a particular axis. This axis is perpendicular to the tissue plane and apico-basal axis (from top to bottom of tissues) of the cell. This phenomenon is called “planar cell polarity, PCP”, and is a common phenomenon found in many multicellular organisms. For example, hair cells in the inner ear of humans have many hairs on each individual cell, and since the hairs are regularly arranged on the cell plane, we can hear the sound. Not only hair cells, but fish scales, and body hair on mammals and the wings of birds etc., have cell polarity. As a result, these creatures carry out macroscopic morphogenesis.
Recently, since molecular biological research has developed concerning PCP, detailed molecular mechanisms have been revealed, several report is published [1-4]. On the other hand, major problems remain. In particular, 1. Mechanisms of morphogenesis for macro-level (PCP) and micro-level (molecular) information and 2. in considering PCP, what are the most important factors?
In order to solve these problems, we construct a simple mathematical model for PCP [5]. This model is based on the molecular mechanisms of PCP. In this talk, we will introduce our mathematical model and simulation result. Despite of our model is very simple formulation, it can reproduce various aspects of the PCP. From the view point of our mathematical model, let’s think “What is the most important mechanism in PCP”.
[1] Jean-François Le Garrec et al., Developmental Dynamics 235:235-246. (2006)
[2] Yoram Burak et al., PLoS Computational Biology Vol 5(12), e1000628. (2009)
[3] Benoît Aigouy et al., Cell, 142, 773-786, September 3. (2010)
[4] K. Amonlirdviman et al., Science Vol 307, 423-426. (2005)
[5] Ayukawa,T., Masakazu, A. et al., Cell Reports, 8(2): p. 610 - 621, 2014.
Cell division requires the precise placement of the division ring at mid-cell to ensure both daughter cells are viable. However, the mechanisms behind this localization remain poorly characterized. There are a limited number of known ways to identify the centre of the cell. One such mechanism is a Turing pattern. One intracellular Turing pattern has been identified, that produced by the Min protein system. In Escherichia coli, the Min protein system plays a role in establishing the division ring position. Membrane-bound Min proteins form an oscillating spatial pattern where the proteins are concentrated at one pole of the cell and then another, leaving a barezone at the centre of the cell where the FtsZ ring will form. Based on molecular interactions of the Min system, we have formulated a mathematical model that reproduces Min patterning during cell growth and division. This model provides a platform to explore how the Min system functions and what characteristics are likely to be shared with other Turing patterning systems, should they exist. We examine the general characteristics of Turing patterns produced by the Min system. In particular, patterning approximates a harmonic of the cell shape and selects the dominant harmonic in a predictable manner. This shows what alternative intracellular Turing patterning systems are likely to appear and how they would behave in relation to cell shape. The oscillations of the Min system are shown to be translated into a mid-cell localization signal via the harmonics generated by non-linear interactions of the system. We show that division plane orientation in the pleomorphic archeon Haloferax volcanii can be predicted from cell shape by assuming that it is dictated by a Turing mechanism. This work makes progress towards understanding how the Min system functions to regulate cell division. More generally it develops tools to identify alternative Turing patterning systems and to understand how patterning can be translated into localization signals.
During development of multicellular organisms, multiple signalling systems play important roles. However, it is very hard to understand the interplays between many signalling pathways using conventional methods of molecular biology, biochemistry and genetics. Mathematical modelling should be combined with these biological methods to solve biological problem.
The waves of differentiation in visual system development are key examples of the complex interplays between multiple signalling systems. In this study, we focus on a wave of differentiation called the proneural wave, which occurs in the developing fly brain and accompanies Notch-mediated lateral inhibition and EGF-mediated neural differentiation. During proneural wave progression, the sheet-like neuroepithelial cells (NEs) sequentially differentiate to neural stem cells called neuroblasts (NBs). The proneural wave is an ideal model system to investigate the roles and dynamics of the interplay between important signalling pathways such as Notch and EGF.
Notch-mediated lateral inhibition regulates binary cell-fate choice, resulting in salt-and-pepper patterns during various developmental processes. The proneural wave accompanies Notch activity that is propagated without the formation of a salt-and-pepper pattern. However, mathematical modelling and genetic analysis clearly demonstrated that Notch-mediated lateral inhibition is implemented within the proneural wave. Because partial reduction in EGF signalling causes the formation of salt-and-pepper pattern, it is most likely that EGF diffusion cancels salt-and-pepper pattern formation in silico and in vivo. Moreover, the combination of Notch-mediated lateral inhibition and EGF-mediated reaction diffusion enables a novel function of Notch signalling that regulates propagation of the wave of differentiation.
In our previous results, Notch signalling is activated only once at the wave front. However, Notch signalling is actually activated again behind the wave front forming twin peaks in vivo. The results of our parameter search show that the twin peaks of Notch activity can be reproduced by elevating the coefficient of cis-inhibition, by which Notch activity is autonomously repressed by Delta ligand. Moreover, the formation of the twin peaks could be stabilized by introducing strong non-linearity to cis-inhibition. The result of our in vivo experiment is consistent with the non-linear behaviour of cis-inhibition. The possible molecular mechanisms of non-linearity in cis-inhibition will be discussed.
Influenza A viruses have caused a number of global pandemics, the most recent being the H1N1 pandemic in 2009, resulting in considerable human mortality. Despite influenza pandemics being rare events, with it currently being nearly impossible to predict the next influenza emergence event, it may be the case that the virus itself provides us with outbreak signals that should prompt us to be more prepared. Specifically, knowing whether or not influenza pandemic occurrences are uniformly random in time can inform the intervention strategies that would be best suited to reduce the risk of further pandemic events.
To determine whether the emergence of new pandemic strains is a memoryless or history-dependent process, we analyse the time periods between influenza pandemics since 1700 statistically under different modelling assumptions. Using Bayesian model selection between exponential and gamma distributions for these time periods, we demonstrate from these small but informative datasets support for the hypothesis of history dependence under eight out of nine sets of modelling assumptions. Using the fitted parameters to make predictions shows a high level of variability in the modelled number of pandemics from 2010–2110.
Though the approach we take relies on limited data, so is uncertain, it provides cheap, safe and direct evidence relating to pandemic emergence, a field where indirect measurements are often made at great risk and cost. Our findings supporting the presence of history dependence enhances the motivation for implementation of active surveillance at the human-animal interface, in particular to elucidate the biological processes behind the observed spacing between pandemics.
Modelling the spread of influenza across Australia is of substantial public health concern. However, there are many challenges in creating accurate models, including how best to capture the spatial and temporal characteristics of the disease spreading process, and aligning with the actual contact process and mobility of individuals in the population. How influenza spreads spatially, whether by infecting neighbouring regions or through major domestic centres via air traffic, is of particular interest in building predictive models. Recent work has suggested the spread of influenza in the United States is dominated by work commutes rather than air traffic. In contrast, influenza in Australia has been shown to have highly synchronized epidemics between major cities across large geographic distances and varying climates.
FluTracking is a participatory online health surveillance system for monitoring influenza-like-illness (ILI), with one of its principal aims to detect the onset of influenza in Australia. Since 2010, there are over 10,000 participants each year, with participation continuing to increase. In this work, we analyse weekly reports from the FluTracking system from May 2011 to October 2017 to infer spatial and temporal dynamics of influenza in Australia. These dynamics could then be used as a starting point for building contact and contagion networks in models for epidemic forecasting, and lead to a greater understanding of the mechanisms behind the spread of influenza-like-illness in Australia.
While there exist a number of mathematical approaches to modelling the spread of disease on a network, analyzing such systems in the presence of uncertainty introduces significant complexity. In scenarios where system parameters must be inferred from limited observations, general approaches to uncertainty quantification can generate approximate distributions of the unknown parameters, but these methods often become computationally expensive if the underlying disease model is complex. In this talk, I will apply transitional Markov chain Monte Carlo (TMCMC) using a massively parallelizable Bayesian uncertainty quantification framework to a model of disease spreading on a network of communities, showing that the method accurately and tractably recovers system parameters and selects optimal models in this setting.
Influenza in humans exhibits a strong seasonal cycle in temperate climates, with a peak of varying intensity appearing each winter. However, the exact cause of this seasonal cycle remains poorly understood. We develop a climate-based SIR modelling framework to understand influenza seasonality, with the transmission rate as a function of climate data. By using a variety of climate-based functional forms of transmissibility from the literature, as well as some new forms, we select the best functional form for climate-dependent transmissibility via modern Machine Learning model selection methods. By analysing a unique dataset comprising ten years of GP-reported influenza-like-illness surveillance data from around Australia, we explore the relationship between influenza transmission and weather in different climate zones.
Reinfection is known to induce complex epidemiological dynamics (e.g. sustained oscillation) due to the time-series change in susceptibility. The simplest model describing reinfection shows three epidemiological dynamics; disease-free, epidemic and endemic. These three dynamics can be classified by two reproduction numbers, basic reproduction number and reproduction number by only reinfection. However, the simplest model takes into account only two variations of susceptibility, susceptibility at first infection and second or later infection. To relax this assumption we construct a parsimonious model describing three susceptibilities, i) susceptibility at first infection, ii) partial protection at second or later infection and iii) perfect protection at second or later infection. This model demonstrates an interesting dynamics, outbreak occurs after the temporal decrease in the number of infected individuals. We named this dynamics as "delayed outbreak". Basic reproduction number cannot capture outbreak potential of this dynamics. Our model is too simple to understand rich dynamics induced by heterogeneity in susceptibility, for example, endemic situation is not captured in this model. To understand the dynamics with heterogeneity in susceptibility, we expand our model to a model describing the transition of susceptibility by reinfection. We confirmed that "delayed outbreak" is robust if the heterogeneity in susceptibility against reinfection exists.
Infections with Group A Streptoccocus (GAS) are highly prevalent in remote communities in
the Northern Territory, Australia. One of the primary drivers of GAS infection is scabies, a small mite which causes a break in the skin layer, potentially allowing GAS to take hold. This biological connection is reaffirmed by the observation that mass treatment for scabies in these remote communities sees a reduction in the prevalence of GAS infection, despite GAS not being directly targeted. In the most extreme case, it has been hypothesised that the eradication of scabies in remote communities may lead to an eradication of GAS related infection.
We start by proposing a model of GAS transmission that assumes that scabies dynamics are at equilibrium. We assume that GAS follows a Susceptible – Infectious – Susceptible (SIS) structure, but individuals who are infected with scabies experience an increased force of infection compared to those who are not infected with scabies. In consequence, we are able to calculate the required prevalence of scabies required to ensure that GAS is eradicated, as a function of $R_0$ and the coupling strength between the two infections.
In order to more accurately model the impact of mass treatment for scabies, we extend this model to include the dynamics of scabies infections. We consider two different scabies models: one which includes the full dynamics of the scabies mites, and one which collapses these dynamics down to a far simpler phenomenological model. We investigate the difference in the scabies rebound dynamics after mass treatment for these two models, and consider the impact of this difference on the dynamics of GAS. We show that despite the two scabies models having different mean field dynamics, parameter uncertainty and the small population size mean that this difference is severely muted when considering the post-scabies treatment on GAS. Finally, we compute a value for the increased force of infection experienced by those infected with scabies, below which eradication for GAS would be achieved if scabies were eradicated.
Accurate mathematical models describing the natural history and dynamics of viral infection can be constructed using data from carefully designed and controlled experiments. These models can be used to extract important quantitative information such as the fold-change in virus production caused by a specific viral mutation or the mode of action and/or efficacy of a novel antiviral compound. The minisymposium will provide examples of how mathematics can provide unique and novel insights to advance knowledge in the field of virology and help drive the development of novel interventions in medicine. The presentations will inform the audience about recent advances in the quantitative understanding of virological mechanics, the optimal experiments to extract the greatest possible amount of quantitative information about a particular viral strain and how this information can impact therapy for viral infections.
Viral interference, whereby infection with one type of virus may temporarily “protect” the host from subsequent infection with another virus has been described for a number of influenza strains and in a number of different host species. In particular, experiments performed in ferrets with influenza A(H1N1)pdm09, A(H3N2) and B have demonstrated strong levels of interference dependent upon the particular virus pair, order of infection and precise timing of re-exposure. Mathematical modelling, which aims to capture the viral time course of the re-exposure behaviour based on plausible immunological mechanisms, is a particularly useful tool to study viral interference, opening up a new approach to advancing our understanding of viral dynamics and the host immune response.
In this talk, I will first briefly describe the ferret re-exposure experiments performed by our collaborators. In the experiments, naïve ferrets were inoculated with one virus (i.e. primary infection) and subsequently exposed to another (to induce the secondary infection) after a certain time interval (1 – 14 days). Viral interference was observed when the inter-exposure interval was less than a week. Furthermore, viruses differed in their ability to induce a protective immune response against subsequent exposure.
I will then introduce the model-based approaches that we have used to investigate immune-mediated viral interference. We begin with a qualitative approach by introducing and analysing a family of within-host models of re-exposure viral kinetics which differ in their hypothesised mechanisms of action for the non-specific innate immune response. We assumed that different viruses stimulate the innate immune response to different degrees and demonstrate the models’ ability to qualitatively reproduce the observed viral shedding profile in the secondary infection for short inter-exposure intervals (less than 5 days). We then extend the models by incorporating cross-reactive adaptive immunity (in particular CD8+ T cells), which provides an explanation of the viral shedding profile of the secondary infection for longer inter-exposure intervals (5-14 days).
Finally, I will briefly describe our most recent efforts to calibrate the proposed models to the re-exposure data in a Bayesian framework, providing quantitative insight. Through this approach, we have shown that models fitted to the re-exposure data are able to predict the outcome of viral interference more accurately than models fitted only to the data from a single infection.
Hepatitis B Virus (HBV) is widespread infectious disease and more than 240 million people are chronic infected so far. Although some patients can suppress the viral load under detection limits by current drug treatments, these are not effective for over 60% of patients. The obstacles for developing effective drugs are the existence of a reservoir in infected cells, which is known as covalently closed circular DNA (cccDNA). Because cccDNA remains in infected cells for a long time, the replication of HBV is not completely inhibited. We need new effective drug that can decrease the concentration of cccDNA persistently. In this study, we carried out cell culture experiments that can be measured the time course of the concentration of cccDNA, intracellular HBV DNA and extracellular DNA. Then, we added several kinds of cytokines in this experiment, and measured the same time course data. Furthermore, to investigate mechanism of action and quantify the anti-viral effect of each cytokine, we established a mathematical model that can capture the dynamics of intracellular HBV replication and fitted this model to the our time course data. Finally, we investigated which cytokines can decrease the concentration of cccDNA. This framework that combining viral experiment and mathematical model would be helpful to find a chemical that have desired mechanism of action.
The host range of human immunodeficiency virus (HIV) is quite narrow. Therefore, analyzing HIV-1 pathogenesis in vivo has been limited owing to lack of appropriate animal model systems. To overcome this, chimeric simian and human immunodeficiency viruses (SHIVs) that encode HIV-1 Env and are infectious to macaques have been developed and used to investigate the pathogenicity of HIV-1 in vivo. So far, we have many SHIV strains that show different pathogenesis in macaque experiments. However, dynamic aspects of SHIV infection have not been well understood. To fully understand the dynamic properties of SHIVs, we focused on two representative strains — the highly pathogenic SHIV, SHIV-KS661, and the less pathogenic SHIV, SHIV-#64 — and measured the time-course of experimental data in cell culture and rhesus macaque and analyzed them. I would like to discuss our quantitative results and future direction of this study.
One of the main focuses of HIV research today concerns allowing people living with HIV to experience prolonged periods where they do not need to remain on treatment. Current therapies are able to suppress HIV to undetectable levels, however as soon as therapy is interrupted the virus “rebounds” to pre-treatment levels and this leads to increased morbidity from HIV. This rebound likely occurs as a result of long-lived latently infected cells, which persist in spite of therapy as part of a ‘latent reservoir’. Recent work has focused on reducing the size of the latent reservoir in the hope this will increase the time subjects can remain off therapy.
We set up a model of HIV remission and viral reactivation and used this to estimate both the probability of viral rebound even with a reduced latent reservoir, and also the optimal reduction in reservoir size that both minimises drug exposure, and maximises treatment effect.
Using a stochastic model we show that a treatment that achieves an average of a one-year viral remission will still lead to nearly a quarter of subjects experiencing viral rebound within the first three months and give rise to an expected 39 (95%UI 22-69) heterosexual transmissions and up to 262 (95%UI 107-534) homosexual transmissions per 1,000 treated subjects over a 10-year period.
Additionally, using a deterministic model we investigated the trade-off between increasing the average duration of remission, versus the risk of treatment failure (viral recrudescence) and the need for re-treatment. To minimise drug exposure, we find that the optimal treatment would increase the average time of viral rebound to 30 years. Interestingly, this is a significantly lower level of reduction than that required for complete elimination of the viral reservoir, but significantly longer than is currently being targeted. We show that when shorter periods are targeted, there is a real probability of viral transmission occurring in between the times at which subjects would be tested for viral rebound.
Our work shows that while therapies designed to reduce the size of the latent reservoir present a promising avenue of research, prior to widespread implementation of such strategies, the possibilities and consequences of viral reactivation from latency must be adequately considered.
Spatial structure greatly influences the demography and evolution of organisms. This minisymposium focuses on the effect of self-structured spatial heterogeneity on the evolution of host-parasite ecological interaction and intra-host immunological interactions. Our specific focus includes multilevel selection in self-organized meta-population structure, evolutionary games played between coinfecting pathogens, host structures favouring cell-free (global) versus cell-to-cell (local) transmissions, and coevolutionary interaction between insects and their pathogens. Both experimental and empirical studies of recent plant virology and insect parasitology will be introduced to test the theoretical predictions.
Spatial structure owing to localized dispersal in hosts can have dramatic impacts upon the coevolutionary dynamics of hosts and parasites. The basic idea is that localized dispersal in hosts can lead to localized transmission, thereby selecting for costly resistance in hosts and lower virulence in parasites. Both of these evolutionary forces are grounded on the altruism (i.e., costly to the actor but benefitial to recipients). From the viewpoints of the hosts, however, they can also utilize the pathogens to harm their neighbors (biological weapon hypothesis), thereby gaining positive inclusive-fitness benefits overall. Despite this, few studies looked at the causes and consequences of host dispersal evolution. Here, using mathematical models of inclusive fitness theory, we show three understudied ideas. First, differential costs of dispersal between susceptible and infected individuals can modify the host population structure. Second, if hosts can condition their dispersal propensity depending on the sickness (susceptible or infected), then the evolutionary dynamics of such conditional host dispersal can generate evolutionary bistabilities such that dispersal propensity biases towards susceptible or infected individuals depending on the initial conditions of the host populations. Finally, we propose that such drastically different dispersal propensity could turn over the classic idea that spatial structure can often favour increased resistance: rather, differential costs of dispersal and/or conditional dispersal propensity in hosts could favour “spite” of the hosts, which has been largely ignored aspects in spatially structured host-parasite dynamics. We overall highlight the neglected but profound aspects of host dispersal.
Most of the viral gene products are shared among a viral population in a host cell, which accumulates up to $10^6$ to $10^7$ genomes. High mutation rates in viral genome replication bring genetic variety to the intracellular population, and this makes the situation social: mutant genomes that do not code intact gene products can survive as free riders, by using the gene products from the other genomes in the intracellular population. We have previously shown that two plant RNA viruses, Japanese soil-borne wheat mosaic virus and tomato mosaic virus, infect a host cell by 5-6 and ~4 genomes in average, respectively, and proposed that this cell infection by a limited number of genomes enhance stochastic separation of free riders from others, resulting in exclusion of free riders from the whole population [1,2]. Here we simulated the evolution of average number of viral genomes to infect a cell (multiplicity of infection, or MOI, in virology). The simulation showed that MOI evolves to 4-5 in silico. This result could be explained by a balance between exclusion of free riders by smaller MOI and avoidance of stochastic failure of infection. Consistently, our further experiments showed that a tripartite plant RNA virus, cucumber mosaic virus (CMV), also has MOIs of 3-5 for its multiple genomic segments in multiple hosts.
The MOI of 4-5 does enhance separation of free riders from others, but still allows their co-existence for a certain period of time. This requires “decisions” by viruses. The intracellular populations exclusively consisting of free riders cannot infect new adjacent cells, and those exclusively consisting of intact genomes will infect the adjacent cells at high probability; those consisting of both free riders and intact genomes need to decide whether to infect the adjacent cells, depending on the ratio of free riders in each intracellular population. Our simulations showed that a majority-decision-like system, where the cell-to-cell infection probability sharply changes between 0 and 1 at a certain threshold for the ratio of free riders, can be favoured by plant viruses. The simulation result was supported by a co-inoculation experiment of CMV variants that carry a functional cell-to-cell movement protein gene or a dysfunctional one.
[1] Miyashita and Kishino (2010) “Estimation of the size of genetic bottlenecks in cell-to-cell movement of soil-borne wheat mosaic virus and the possible role of the bottlenecks in speeding up selection of variations in trans-acting genes or elements.” J Virol 84(4) 1828-1837.
[2] Miyashita, Ishibashi, Kishino and Ishikawa (2015) “Viruses roll the dice: the stochastic behavior of viral genome molecules accelerates viral adaptation at the cell and tissue levels.” PLOS Biol 13(3) e1002094.
Viruses have two modes spread in a host body, one is to release infectious particles from infected cells (cell-free) and the other is to infect directly from an infected cell to an adjacent cell (cell-to-cell). Since the mode of spread affects the evolution of life history traits, such as virulence, it is important to reveal which mode is selected. Here we show the evolutionarily stable proportion of cell-free and cell-to-cell infection, and how it depends on the spatial distribution of target cells. Using an epidemic model on a 2D regular lattice, we consider the infection dynamics by pair approximation and check the evolutionarily stable strategy. We also conduct the Monte-Carlo simulation to observe evolutionary dynamics. We show that a higher cell-to-cell infection is selected as target cells become clustered. The selected strategy depends not only on the degree of clustering but also the abundance of target cells per se.
The impact of spatial structure on evolutionary outcomes can be profound due to both ecological and genetic correlations. One of the best developed areas of spatial evolutionary theory has focused on the coevolution of hosts and parasites. There are profound impacts of local as opposed to global infection on both parasite and host evolutionary outcomes. Using a combination of pair approximations to apply approximate adaptive dynamics to spatial models and simulations we show the central importance of reproduction to the outcome. For example spatial effects on the evolution of resistance arise to due less costly local reproduction rather than spatial epidemiological patterns, the host spatial structure that is created by local reproduction has profound effects on parasite evolution even if transmission is mostly global and the impact of disease on reproduction rates impacts ESS virulence. This highlights the importance of indirect spatial effects on the spatial evolutionary ecology of focal traits.
Spatial gradients of chemicals and microorganisms generated by their advection and diffusion in a fluid have important consequences for signalling, waste removal, feeding, and behaviour. Modelling these systems often requires a combination of multiscale and multiphysics approaches, including solving the equations of fluid motion, describing the advection and diffusion of a (chemical) species, and modelling the behaviour of a cell or a microorganism. In nearly all cases, experimental validation is necessary to check whether or not all of the important elements of the biological and physical systems have been correctly considered. In this minisymposium, we provide examples of such systems that span the level of single-cell dynamics to the behaviour and distribution of groups of organisms. Specific examples include the pumping and mixing of fluid by flagella, spatial gradients of calcium during development, the movement of chemicals through hair arrays for the purpose of sampling, and the movement of zooplankton through beds of vegetation. The mathematical models and numerical methods highlighted in this minisymposium include the Immersed Boundary Method and the Method of Regularized Stokeslets for solving the fully-coupled fluid-structure interaction problem of fluid moving through and around cells, organs and organisms; agent-based models to describe the behaviour and movement of microorganisms; biochemical models to describe the response of cells or microorganisms; and numerical methods for solving the advection and diffusion of a concentration gradient. Challenges in terms of mathematical modelling, numerical simulation, and model validation will be highlighted throughout the presentations.
It is difficult to mix and pump fluid in microfluidics devices because the traditional methods of mixing and pumping at large length scales don’t work at small length scales. Experimental work has suggested that rotating helical flagella may be used to effectively mix and pump fluid in microfluidics devices. To further explore this idea and to characterize the flow features around rotating helices, we study the hydrodynamic interactions between two rigid helices rotating at a constant velocity. Helices are coupled to a viscous fluid using a numerical method based upon a centerline distribution of regularized Stokeslets, and we analyze the effects of spacing and phase shift on mixing and pumping.
Microscale filtering and protective layers appear in a variety of places throughout the biological world, with examples including both internal physiological examples (extracellular proteins, microvilli, cilia) and external biological structures (trichomes, swimming legs, bristled wings). In this talk, I describe an agent-based framework built for exploring the biological environment created by these structures and its potential effect on small, lightweight organisms. Using implemented analytical models (e.g. Brinkman) or input from CFD packages to specify the flow field, the framework can simulate a large variety of behaviours in both 2-D and 3-D time-varying environments as a testing ground for population-level theory. Preliminary numerical results will shown demonstrating some of the capabilities of the approach as well as future directions.
The field of neuromechanics, which attempts to integrate neural control with muscle activation and the resulting movement of an organism, is an emerging field of organismal biology. Approaches from experimental biology, robotics, and mathematics have now reached the point where their knowledge about the different facets of (animal) locomotion can be combined and integrated into realistic and useful models. This integration between fields opens up new questions and problems for researchers to solve. In this talk, we explored the interface between neuron and muscle during the initial transfer of information that generates muscle contraction.
Muscle contraction occurs when neuron action potentials stimulate the sarcoplasmic reticulum, releasing calcium ions that then bind to muscle filaments, allowing the myosin to contract the muscle. We will present a model that couples spike-train activation to mass-action equations for calcium ions, which in turn couples to the Hill model for muscle contraction force-velocity relationship. This model allows for investigation of spike trains required to produce partial twitch response, and total tetanic contraction.
It has been known for some decades that fertilisation of some amphibian and fish eggs is followed by a wave of calcium ions over the surface of the egg, which is associated with a physical change to the surface. Similar waves are seen at other stages of embryonic development. An unfertilized mammalian egg is surrounded by cumulus cells to form a cumulus-oocyte complex (COC). Just a few years ago, medical researchers at the University of Adelaide identified, for the first time, a wave-like behaviour of the cumulus cells in COCs after sperm had been added to the culture medium, believed to be a response to fertilisation. From the speed of the wave it was inferred that the cells were responding to one or more chemical signals from locations on the surface of the egg and that calcium ion concentration was the likely signal. More recently it has been determined that the wave is associated with cumulus-cell death.
I will describe some ongoing modelling and experimental work being undertaken to explain the behaviour of the cumulus cells, assuming the wave is initiated at the fertilisation site. We will examine whether a chemical release by either or both of the oocyte and cumulus cells is qualitatively consistent with experimental observations.
This minisymposium will address topics in mathematical neuroscience and will serve as a hub for scientific interaction of the Mathematics Neuroscience subgroup.
In volume transmission, neurons in one brain nucleus send their axons to a second nucleus where neurotransmitter is released into the extracellular space. In [1] we showed how to calculate the average amount of neurotransmitter at different parts of the extracellular space, depending on neural properties and the geometry of the projections and extracellular space. We showed how to formulate questions as boundaries value problems for the heat equation with stochastically switching boundary conditions, and we derived results in one space dimension.
Here we discuss the two- and three-dimensional problems, along with mechanistic models of neurotransmitter homeostasis.
[1] Lawley, Best, Reed 2016
In this talk presented at the Mathematical Neuroscience Subgroup Minisymposium I will discuss a numerical analysis of the steady-state solutions of a neural field model of the corticothalamic system. The independent synaptic connections of the corticothalamic model define an eight-dimensional parameter space, while specific combinations of these connections parameterize intracortical, corticothalamic, and intrathalamic loops.
The first part of the analysis consists of a systematic identification of multistable regions for physiological parameter ranges representing normal arousal waking states in adult humans. Parameter values have been derived from human electroencephalographic recordings (EEG). The key results are the confirmation of the existence of up to five steady-state solutions, up to three of which are linearly stable. The presence of multiple stable steady-state solutions implies that the system has enough degrees of freedom to capture multiple operating points (e.g., mutiple global brain states). Indeed, two out of the three linearly stable steady-state equilibria represent two waking modes for adult human physiology while resting with the eyes closed. These two modes are termed the low- and high- waking modes.
While the low-waking mode has been previously identified with normal brain activity during quiet wakefulness, the high-waking mode has not been fully characterized. The second part of the analysis focus on (i) the spectral features of high-waking mode states; (ii) the nonlinear attractors between the two waking modes; and, (iii) the switching dynamics between the low- and high-waking modes using small amplitude and bandlimited random perturbations.
We argue that the high-waking mode may represent hyperarousal waking states. This result opens up the possibility to predict and identify subtypes of primary insomnia in which cortical hyperarousal is the main hallmark from human neuroimaging data.
Thermoregulatory responses are partially controlled by the preoptic area and anterior hypothalamus (PO/AH), which contains a mixed population of temperature-sensitive and insensitive neurons. In [1] based on physiological data, a Hodgkin-Huxley-like conductance based model was constructed. This model suggests that most PO/AH neurons have the same types of ionic channels, but different levels of channel expression can explain the inherent properties of the various types of temperature-sensitive and insensitive neurons which is encoded in their frequency sensitivity relative to temperature.
Here we present a detailed bifurcation analysis of this model to confirm these observations. We focus on three main physiological bifurcation parameters, the temperature T and the maximum conductances of two specific background potassium leak channels, $g_{task}$ and $g_{trek}$, that are known to be expressed in these PO/AH neurons. These three bifurcation parameters are sufficient to explain the dynamics of PO/AH neurons observed in experiments.
If time permits, we also discuss the multiple timescales inherent in this model and explain the creation of action potentials based on a geometric singular perturbation analysis.
[1] Wechselberger et al., 2006
Neural spiking and bursting rhythms in space-clamped (i.e., ODE) models are typically driven by either canard dynamics or slow passage through Hopf bifurcations. In both cases, solutions which are attracted to quasi-stationary states (QSS) sufficiently before a fold or Hopf bifurcation remain near the QSS for long times after the states have become repelling, resulting in a significant delay in the loss of stability and hence in the onset of oscillations. In this work, we present the spatio-temporal analogues of these delayed bifurcation phenomena in multi-time-scale reaction-diffusion equations. We show the existence of canard-induced bursting rhythms in a spatially extended model of the electrical activity in pituitary cells. We then derive asymptotic formulas for the space-time boundaries that act as buffers beyond which solutions cannot remain near the repelling QSS (and hence at which the delayed onset of oscillations must occur) for slow passage through Hopf bifurcations in reaction-diffusion equations.
This is joint work with Tasso J. Kaper (Boston University).
Gene Drives are mechanisms which involve biased inheritance of a particular gene can be used to spread a desirable trait through a population. They have been proposed as a means to combat a number of diseases and pest species. However, this technology is not without controversy: beyond general concerns about the development of genetically engineered organisms (particularly GM animals), there are important issues relating to accidental release and unintended spread of gene-drive bearing individuals, as well as unintended ecological and evolutionary consequences of releases.
Mathematical modelling can be used at various stages of the development of gene drives, including testing whether a proposed construct will be able to spread in a population, estimation of invasion/containment probabilities, and consideration of potential gene drive reversal strategies. This is particularly timely due to recent advances in molecular biology (such as the CRISPR/Cas9 system) that have aided the design and development of gene drives in an ever widening range of organisms.
In this minisymposium, we will explore how modelling is intricately involved in the development of gene drive technologies and consideration of its downstream consequences, both positive and negative, in a number of contexts, including mosquito-borne diseases and control of invasive rodent species.
Gene drives are a powerful technique to reduce the size of a population or to transform it to become less troublesome. Major concerns include accidental escape of a gene drive, which might lead to widespread elimination of a species, or unintended consequences of released individuals, which might include evolution of resistance to an effector gene. Both spatial localization of a gene drive and the ability to reverse a gene drive are desirable mechanisms to increase the safety of drives. A number of such approaches have been proposed in the literature. Using mathematical modelling, we will explore their likely efficacy and discuss critical weaknesses.
Synthetically-constructed gene drives, that use CRISPR-Cas9 technology to bias inheritance of a particular gene, have been proposed for exotic pest eradication. This is a controversial idea, as there is an apparent potential risk in the use of such technologies, for example to non-target populations of the species. Further, the basic question — Is the technology actually able to achieve its goal in practice? — is vexed. I will discuss some of our work which has developed stochastic, individual-based models of a number of gene drive strategies that have been proposed for the possible eradication of exotic vertebrate pests. I will discuss some of the issues surrounding their likely effectiveness for eradication, explore the potential risks to non-target populations, along with considering options for alleviating these issues.
Novel strategies for controlling mosquito vectors of human diseases involve introductions of selfish genetic elements (SGEs), or elements that spread through the mosquito population by means of non-Mendelian inheritance. The releases of the endosymbiotic bacteria Wolbachia into Aedes aegypti populations in order to reduce their capacity for arbovirus transmission are the first application of an SGE introduction strategy in a field setting. Gene drive technologies for creating SGEs that reduce vectorial capacity are also being developed. SGEs are vertically transmitted between mosquito hosts, and their spread dynamics have complex interactions with natural mosquito population dynamics, population structure and fitness components. We develop a model of an Ae. aegypti metapopulation that incorporates empirical relationships linking variation in mosquito abundance to density-dependent fitness components. We first show that our model can produce patterns of demographic variation similar to those found in natural populations, and predict rates of spatial spread of Wolbachia matching those observed following field releases. We then explore different SGE release strategies over an operationally relevant spatial scale, and compare the spread dynamics of Wolbachia and gene drive. We find that the spatial release distribution strongly affects the rate of spread of both SGEs, with widely dispersed release distributions out-performing spatially aggregated distributions, especially when the SGE incurs fitness costs. Our results demonstrate how spatial and demographic structure in the mosquito host population impacts SGE spread, with implications for the design of effective release strategies.
The advent of CRISPR/Cas9-based gene editing and its demonstrated ability to streamline the development of gene drive systems has reignited interest in its application to the control of mosquitoes and the diseases they transmit. The versatility of this technology has also enabled a wide range of gene drive architectures to be realized, creating a need for their population-level and spatial dynamics to be explored. To this end, we present MGDrivE (Mosquito Gene Drive Explorer): a simulation framework designed to investigate the population dynamics of a variety of gene drive architectures and their spread through spatially-explicit mosquito populations. MGDrivE is based on a tensor algebraic generalization of the lumped age-class model of mosquito ecology. Treating these population dynamic equations in a variable-dimension tensor form allows them to be left unchanged while modifying the dimensionality of the tensor describing inheritance patterns, as required by the drive system. Spatial dynamics are accommodated through a metapopulation structure in which lumped age-class models run in parallel and migrants are exchanged between metapopulations at defined rates. Example MGDrivE simulations are presented for threshold-dependent drive systems: a) reciprocal chromosomal translocations, and b) toxin-antidote-based underdominant systems. Criteria are described for which these systems may be confined to their release community without spreading to neighboring communities.
A goal of all health researchers is to determine effective control strategies that can be provided to patients that will help them recover, or (when recovery is unachievable) manage their health conditions. Mathematical and statistical models are employed within these contexts so that the effects of current control strategies can be quantified, and so that optimal control strategies can be identified (given model assumptions and uncertainties). Recent developments of disease modelling have seen mathematical models being developed to address vaccination, drug therapy use, public health campaigns, and even online applications used to journal chronic condition experiences. This minisymposium aims to showcase recent modelling research of CDM and CRM collaborators on a wide variety of topics on disease control that are of theoretical merit and are significant to patient outcomes. The topics include HIV transmission and progression, haematological dysregulation, and immune system components in liver function. These topics will be of interest to mathematical modellers and health researchers alike. The minisymposium provides the participants and audience an opportunity to discuss emerging research topics and new areas for research collaboration in the fields of disease modelling and health. An underlying theme among all presentations is how non-instantaneous effects are incorporated in modelling equations.
The speakers will begin with a general overview of their field of study, introducing the biology pertinent to the question, and the mathematical tools that will be used. All talks will detail the mathematical models, results, and interpretation to healthcare and medical outcomes.
Topics covered will be regulation and control of platelet production, aspects of the immune system in liver regulation, lifespans of red blood cells and modelling memory CD4 T-cell populations in HIV patients.
Memory translates into time delays naturally in a number of regulatory processes at all levels of organisation in the life sciences: transcription and translation times in molecular biology, finite axonal conduction velocities between neurons, maturation times of precursor cells in hematopoiesis and infection and temporary immune periods in infectious disease propagation are but a few examples of naturally occurring such delays. In many instances, mathematical constructions are elaborated to avoid incorporating time delayed arguments in the modelling equations of the system of interest. We present, by way of examples on two specific physiological and epidemiological systems, the limitations of this approach and the extent to which these approximations may or may not be useful or necessary for a proper understanding of the modelled system.
This talk will serve as an introduction to the other presentations of the Minisymposium.
Modelling memory in physiological regulation
Jacques Bélair talk swapped with Joseph Mahaffy 12:00↔10:30
Although antiretroviral therapy (ART) suppresses viral replication, patients still suffer from both low CD4 T-cell counts and HIV persistence, requiring them to remain on complex ART regimens for life. A naturally occurring 32-base pair deletion in the CCR5 gene, the major co-receptor for HIV entry, is associated with infection resistance.
In the study initiated by Sangamo Therapeutics, HIV-infected subjects received a single infusion of autologous CCR5-modified CD4 T-cells. Following infusion, we observed a sustained increase in the CD4 T-cells count (mean ~162 cells/μL) with a significant decrease in the HIV reservoir (median ~1 log 3yrs post infusion). Long-term persistence of CCR5-modified cells suggested their presence within long-lived CD4 populations, such as T memory stem cells (TSCM). To investigate the impact of persistence of “HIV-resistant” TSCM on immune homeostasis and the decay of the HIV reservoir, we developed a mathematical model of T-cell dynamics. The model follows a linear transition from the naïve to effector memory (TEM) state and includes thymic input of naïve cells, proliferation, death and differentiation rates for naïve and memory cells. Model fits to patient data pre- and post-injection and sensitivity analysis results that increased thymic output, an increased TSCM proliferation rate, a decreased central memory cell death rate, and an increased central memory transition rate played an important role in increasing the T-cell count and decreasing the HIV reservoir.
Finally, using a bi-phasic decay model, we show that purging and replacement dynamics of the CD4 T-cells population post infusion account for a significant fraction of the observed decrease in the HIV reservoir. Our results indicate that homeostatic processes can control HIV persistence, and that T-cell restoration post infusion of CCR5-deleted cells leads to the decay of the HIV reservoir.
The liver is a spatially complex and heterogeneous network of blood and bile flows coupled with metabolic processing and a favoured target for infection by hepatic viruses. We present here a mathematical model aimed at investigating these intrinsic heterogeneities and their impact on the dynamic of the Hepatitis-B variant (HBV). Dramatic spatio-temporal scaling from individual hepatocytes to the entire liver organ invites multi-scale approaches inspiring assembly of sinusoid-level 'unit models' into a whole organ representation. Each 'unit-model' sinusoid combines individual hepatocytes communicating with blood flow in the sinus in turn connected with other 'unit-model' sinusoids aggregated into a whole liver modelling scheme. This permits investigating impacts on the whole organ of precisely distributed spatial heterogeneities such as varying HBV uptake mechanisms (e.g., the sodium-taurochloriate cotransporter or NTCP), immune cell responses (e.g., cytolytic or interferon-based) and simple efficiency of HBV replication. We present our results showing how heterogeneities of, in particular, HBV replication efficiencies may be responsible for persistent chronic infections.
Thrombopoiesis is the process for producing platelets, which uses a negative feedback to maintain homeostasis in normal individuals. However, pathological states exist where platelet concentrations in the body oscillate. An age-structured model for thrombopoiesis was developed and fitted to clinical data for subjects with normal and pathological platelet production. Variations on this model are described to obtain more details on how this system undergoes a Hopf bifurcation. Our study explores parameter sensitivity and which model features are most significant in the bifurcation to periodic solutions. We observed that near the Hopf bifurcation there is a very rapid transition of the stationary solution along with the change in the real part of the leading pair of eigenvalues. The creation and analysis of the characteristic equations from these models provide some interesting new ideas. Certain model reductions decrease the complexity of the characteristic equation and allow a better understanding of the source of the Hopf bifurcation. Our modelling efforts might improve insight into the primary problems underlying the diseased state, where levels of platelets and thrombopoietin vary periodically.
Modelling memory in physiological regulation
Jacques Bélair talk swapped with Joseph Mahaffy 12:00↔10:30
All social insects live in elaborately organised societies. Their social structures enable them to continuously manage a large set of simultaneous tasks; from scouting and foraging to colony defence, nest building, thermoregulation, and brood care. To ensure colony survival and reproduction it is vital that the colony workforce is adequately allocated to these different tasks. Social insect colonies are able to perform such task allocation with an amazing degree of flexibility, continuously adjusting to changes in environmental conditions. The resulting division of labor (DOL) is a cornerstone of all social insect societies and commonly credited as a key factor for their enormous ecological success.
The vast majority of DOL modelling has either focussed on physiological factors or on individual-based behavioural mechanisms. Comparatively little attention has been given to how social interactions regulate task allocation. None of the existing frameworks provides a direct link to integrate environmental conditions as an integral part of the modelling process (beyond task-related stimulus levels). Entomologists are now investigating specifically these factors as the keys to a better and more comprehensive understanding of workforce allocation.
We propose a new mathematical modelling framework for task allocation in social insects that directly incorporates social interactions and environmental conditions as first class elements. We use evolutionary game theory to investigate how task preferences develop during the lifetime of a colony and how DOL emerges and morphs as a result. We model the task preferences of individuals as continuous strategy choices in a simple game, in which individuals choose how to divide their efforts between different tasks. Task execution results in rewards that may be shared between colony members or go directly to the individual. The amount of reward is modulated by the collective level of investment into the task as well as by environmental factors.
Based on this framework we show that ecological conditions are a crucial determinant for the emergence of different forms of specialisation. We also find strong interactions between the prevailing learning mechanisms and colony specialisation.
Contrary to intuition we find that strong specialisation does not always relate to improved colony efficiency. Particularly, social learning appears to lead to high colony performance only in a limited range of environment types and can drive a colony into “over-specialisation” in other conditions. This suggests that social learning should only have evolved in a limited range of environments and for certain tasks. Our theoretical results thus point towards promising new directions for empirical work that can bring us a step closer to understanding how social insects achieve their outstanding ecological success.
In social insect colonies we often observe substantial levels of laziness, i.e., workers that do not engage or appear not to engage in any tasks. This is puzzling, because colonies with lazy individuals seem to be wasteful and may not be using their resources optimally. A common hypothesis is that lazy workers are a reserve workforce that can be quickly activated when rapid changes in the environment demand it.
We present an alternative explanation for laziness based on game theory. We assume that insects play a game in which they allocate their energy towards two competing tasks, foraging and temperature regulation. Both tasks need to be performed successfully for colony survival, thus the set-up reflects a coordination game.
The benefits of regulation are assumed to be concave and maximised at intermediate levels of task engagement. The benefits of foraging are assumed to be linear in effort. Benefits are shared by the whole colony, whereas costs are borne individually. Individuals learn to play this game by using social learning, and occasionally experimenting with new strategies. When experimentation is rare the population is monomorphous most of the time, and the learning trajectories of colonies can be modelled by using adaptive dynamics. Levels of laziness in the long run can be determined by studying the fixed points of the dynamics.
We find that the level of laziness crucially depends on the costs of performing the tasks. In particular, laziness can be stable under marginally increasing costs. Constant and marginally decreasing costs result in colonies without lazy individuals. Monte Carlo simulations show that our theoretical prediction is accurate under moderate levels of noise.
Our model is the first one to link specific ecological features to the puzzle of lazy workers. We discuss potential implications for empirical work and extensions of the model.
The intensification of farming indisputably changes the ecology of agricultural infectious diseases with modern management practice impacting both host density and between-host contacts. We typically assert that these changes result in higher disease burdens, but there is relatively little specific modelling testing this idea. For example, the industrialisation of apiculture has been suggested to play a critical role in the increased losses of Apis mellifera colonies to infectious diseases, with specific apicultural practices (reflecting industry-wide intensification) implicated in higher disease burdens. Here we build multi- honeybee colony models to examine how ‘apicultural intensification’ in its most plausibly extreme case is predicted to impact honeybee pathogen epidemiology within apiary. Counter to the prevailing view, our models predict that intensification, captured though increased population sizes, changes in population network structure, and increased between-colony transmission, has little effect on disease burden across an apiary. Disease burdens, and therefore impacts of intensification, are critically determined by the R0 (reproductive ratio) of a given pathogen prior to intensification. The greatest impacts of intensification are found for diseases with initially low R0 values, however such diseases cause little overall disease burden and so impacts are universally minor. Importantly in this context, the smallest impacts of intensification are found for diseases with high R0 values, which appear to be typical of important honeybee diseases. Our findings highlight a lack of support for the hypothesis that current and ongoing aspects of management intensification lead to higher disease burdens in honeybees. More broadly, our work demonstrates the need for specific models of agricultural systems and management practices in order to understand the implications of management changes on disease burden.
In recent years honey bee colonies have been experiencing increased loss of hives. One cause of hive loss is colony collapse disorder (CCD). Colony collapse disorder is characterised by a previously healthy hive having few or no adult bees but with food and brood still present. This occurs over several weeks. It is not known if there is an exact cause of CCD but rather it is thought to be the accumulation of multiple stressors placed on a hive. One of these stressors is the breakdown of thermoregulation inside the hive.
The bee life cycle begins with eggs that hatch into larvae that in turn pupate. The eggs, larvae and pupae together are known as brood. The hive contains combs which are made up of multiple cells; these cells house the brood. Pupal cells are capped off by adult bees (and so are known as capped brood) and they undergo changes to develop into an adult bee. In order for this capped brood to develop correctly, the temperature within the hive must be regulated by the hive bees to ensure optimal development of the capped brood. Variations in the temperature, caused by the breakdown of thermoregulation, lead to suboptimal development in adults that emerge from capped brood. In particular, their brains and flight muscles are compromised. This later leads to these bees becoming inefficient foragers which also have shorter life spans.
We model the effect of thermoregulation on hive health using a system of DDEs which gives insights into how varying hive temperatures have an effect on the survival of the colony. We show that thermoregulatory stress has the capacity to drive colony collapse disorder via a saddle-node bifurcation.
A conceptual model was constructed to define the network of potential routes of exposure of pollinators to pesticides in greater detail than has previously been done. This model provides a basis for biologically and ecologically realistic basis for mathematical estimation of exposure versus time both individually and at the colony level. It also shows the distinction between primary exposure routes, such as contact during foraging with vegetation and soil, contamination of drinking water, and secondary routes of exposure in the hive including consumption of pollen, and nectar that are collected by foragers or the honey and bee bread produced in the hive. Secondary exposure from contact with nest construction materials such as wax are included. These secondary routes of exposure can affect both adults and immature life stages of social pollinators in the hive or nest. The conceptual model allows detailed consideration of the overall pattern of potential exposure, including the dynamics of exposure to various castes and life stages of social bees. This conceptual model also aids in assessing the completeness of the risk assessment and identifying exposure pathways that require further investigation.
Modelling memory in physiological regulation
Bélair 12:00 → 10:30
Mahaffy 10:30 →12:00
Social insects
Meyer 10:50 → 10:30
Garcia 10:30 → 10:50
Tumours are not simply collections of mutated cells that grow in isolation of the environment in which they live. They interact with and modify both the physical microenvironment and a variety of nontumor cells that make up the organ in which the cancer originated. In nature, ecology and evolution are intimately linked since one provides the players and the other provides the field. The idea of viewing cancer from an ecological perspective fundamentally means that we cannot just consider cancer as a collection of mutated cells but rather a dynamic system of many interacting cellular and microenvironmental elements. Heterogeneity is both a cause and consequence of the interaction between the tumour and its environment and has been observed across genotypic, phenotypic, and environmental scales and has now been recognized as a key driver in cancer drug resistance and treatment failure.
Importantly this complex dialogue is governed by Darwinian dynamics, i.e. local microenvironmental conditions select phenotypic tumour clones that are best adapted to locally survive and proliferate. Conversely, the phenotypic properties of the cells (through release of growth factors or metabolites) affect the environmental properties. While these complex interactions have enormous clinical implications because they promote resistance to therapy, they rarely incorporated into clinical design. In fact, many treatment approaches only focus on specific molecular targets without considering the context or consequence.
In this minisymposium, we present ecological and evolutionary perspectives on cancer progression and treatment. By combining models with biological or clinical data, progression can be better understood and novel treatment strategies can be designed based on ecological and evolutionary principles.
Evolutionary game theory has been used to model cancer for more than a decade. Efforts to date have focused on understanding the effect of interactions between cancer cells of different types, and between aspects of the tumour microenvironment and cancer cells. To realise the full potential of these modelling efforts however, we submit that a method for direct parameterisation is required. In this talk, I will present our novel ‘evolutionary game assay’ designed specifically to parameterise a two strategy matrix evolutionary game from in vitro experiments comprised of co-culture of cells transfected with different colour fluorescent proteins. I will also show the results of the first game we measured, a game between EML4-ALK translocated non-small cell lung cancer cells before and after we evolved resistance to Alectinib (a targeted therapy). We find that these cells play either the Leader or Deadlock game, depending on different therapeutic and microenvironmental conditions. This change in game, as modulated by treatment, also suggests a novel therapeutic strategy which we have posited before is possible: treating the game instead of the player.
Several types of cancer initiate or metastasize to the bone. These include the most prevalent and lethal cancers: lung, breast and prostate. Thus understanding the bone ecosystem is key if we want to predict what phenotypes will successfully metastasize to the bone and the subsequent evolutionary dynamics that will ensue as the invading tumour cells learn how to co-opt the bone resident cells. This bone ecosystem is very dynamic and responds and maintains homeostasis as micro and macro fractures appear. A variety of homeostatic processes emerge from the interactions of a myriad of cell types and signalling factors. Here I will present mathematical models, biologically motivated and tested, that explain bone repair and homeostasis and allow us to explore how tumours can grow in the bone as well as the impact of a variety of treatments not only on the tumour cells but also on the bone ecosystem.
Glioblastoma is the most aggressive primary brain cancer, with poor survival that can be largely attributed to intra-tumoural heterogeneity. While these tumours are primarily monitored via contrast-enhanced (CE) T1-weighted and T2-weighted magnetic resonance (MR) images, these standard clinical images are known to be non-specific in their correlation with tumour cell density. This lack of relationship makes it difficult to define the specific regions of interest to target for surgery and radiation. Previous efforts have shown some promise in better interpreting these images utilizing either machine learning (ML) or mechanistic modelling independently. But methods to harness the strengths of both approaches are sorely needed to make clinically actionable progress.
Here we present a novel, first-of-its-kind, hybrid model which brings together a graph-based semi-supervised machine learning approach with a mechanistic partial differential equation model of glioblastoma growth, known as the Proliferation-Invasion (PI) model, to generate predictive tumour cell density maps with high accuracy. Our ML approach bridges cell density as quantified from image-localized biopsies with texture analysis of multiparametric MR images. To incorporate the mechanistic model, the PI model is first used to generate an independent prediction of cell density which is then introduced into the ML algorithm through a Laplacian matrix, which ensures regions with similar predictions from the PI model will have similar predictions in the final model.
We have applied our proposed ML-PI model framework to 18 patients with a total of 82 image localized biopsies. Each patient’s tumour was imaged with multi-parametric MR images, including T1-weighted, CE T1-weighted, T2-weighted, dynamic contrast-enhanced (DCE) imaging, diffusion weighted imaging (DWI), and diffusion tensor imaging (DTI). In this cohort, our hybrid model was able to achieve higher accuracy in cell density prediction than either of the independent models (ML or PI) alone, with a mean accuracy prediction error of 0.084 vs 0.227 for PI alone and 0.220 for ML alone. We hope that with more verification, this tool can be used to not only guide spatially localized therapies such as surgery and radiation, but also help broadly in the interpretation of images for glioblastoma patients.
Tumours consist of a hierarchical population of cells that differ in their phenotype and genotype. This hierarchical organization of cells means that a few clones (i.e., cells and several generations of offspring) are abundant while most are rare, which is called clonal dominance. Such dominance also occurred in published in vitro iterated growth and passage experiments with tumour cells in which genetic barcodes were used for lineage tracing. A potential source for such heterogeneity is that dominant clones derive from cancer stem cells with an unlimited self-renewal capacity. Furthermore, ongoing evolution and selection within the growing population may also induce clonal dominance. To understand how clonal dominance developed in the iterated growth and passage experiments, we built three computational models that simulate these experiments. In the first two models we used the tau-leaping Gillespie algorithm to simulate random growth and passage and we either consider all cells to have an equal proliferation capacity (model 1) or we distinguish between cancer stem cells and differentiated cells (model 2). In the last model, we used a dynamic Monte Carlo method to simulate the development of a cell population where division rates are heritable and vary in between cells due to initial variation and mutation (model 3). Only simulations with the model where division rates varied between the cells (model 3) reproduced the clonal dominance that developed in in vitro iterated growth and passage experiments. In contrast, the experimental results can neither be reproduced with a model that considers random growth and passage (model 1), nor with a model based on cancer stem cells (model 2). Altogether, our model suggests that in vitro clonal dominance develops due to selection of fast-dividing clones.
This session will tackle a variety of biological problems from participants in the third Women Advancing Mathematical Biology (WAMB) workshop. The WAMB workshop was held in April 2017 at the Mathematical Biosciences Institute at The Ohio State University, and teams of senior and junior scientists and mathematicians were created to foster new research collaborations among women in mathematical biology. In this minisymposium, participants from this workshop will present on their ongoing research and foster continued collaborations with others at the conference. Since the workshop brought together a variety of participants and research areas, this session also includes a broad background of topics. The topics include: comparing continuous and discrete models of Ornithodoros moubata (soft ticks); stochastic models of super spreading events; fluid dynamics of nematocyst prey capture; and characterizing Clostridium difficile colonization through networks and mechanistic models. This session will combine ODE models, stochastic models, fluid dynamic models, and network models together in one session. Moreover, the models will span from in host bacteria models to population epidemiological models.
The first mathematical models for an argasid tick are developed to explore the dynamics and identify knowledge gaps of these poorly studied ticks. These models focus on Ornithodoros moubata, an important tick species throughout Africa and Europe. Ornithodoros moubata is a known vector for African swine fever (ASF), a catastrophically fatal disease for domesticated pigs in Africa and Europe. In the absence of any previous models for soft-bodied ticks, we propose two mathematical models of the life cycle of O. moubata. One is a continuous-time differential equation model that represents the tick life cycle with two stages, and the second is a discrete-time difference equation model that uses four tick stages.Both models use two host types: small hosts and large hosts, and both models find that either host type alone could support the tick population and that the final tick density is a function of host density. While both models predict similar tick equilibrium values, we observe significant differences in the time to equilibrium. These models provide the basis for developing future models that include disease states to explore infection dynamics and possible management of ASF.
The importance of host transmissibility in disease emergence has been demonstrated in historical and recent pandemics that involve infectious individuals, known as superspreaders, that are capable of transmitting the infection to a large number of susceptible individuals. To investigate the impact of superspreaders on epidemic dynamics, we formulate deterministic and stochastic models that incorporate differences in superspreaders versus non-superspreaders. In particular, continuous-time Markov chain models are used to investigate epidemic features associated with the presence of superspreaders in a population. We parameterize the models for two case studies, Middle East respiratory syndrome (MERS) and Ebola. Through mathematical analysis and numerical simulations, we find that the probability of outbreaks increases and time to outbreaks decreases as the prevalence of superspreaders increases in the population. In particular, as disease outbreaks occur more rapidly and more frequently when initiated by superspreaders, our results emphasize the need for expeditious public health interventions.
A nematocyst is a specialized organelle within cells of jellyfish and other Cnidarians that sting. Nematocysts are also present in some single celled protists. They contain a barbed, venomous thread that accelerates faster than almost anything else in the animal kingdom. Here we simulate the fluid-structure interaction of the barbed thread accelerating through water to puncture its prey using the immersed boundary method. For simplicity, our model describes the discharge of a single barb harpooning a single celled organism, as in the case of dinoflagellates. One aspect of this project that is particularly interesting is that the micron-sized barbed thread reaches Reynolds numbers above one, where inertial effects become important. At this scale, even small changes in speed and shape can have dramatic effects on the local flow field. This suggests that the large variety of sizes and shapes of nematocyst may have important fluid dynamic consequences. We find that reaching the inertial regime is critical for hitting prey over short distances since the large boundary layers surrounding the barb characteristic of viscous dominated flows effectively push the prey out of the way.
In recent years, new technology has allowed us to achieve measurements of hundreds of metabolites and the expression of thousands of genes. With this large scale, all-encompassing, ‘omics’ data comes a critical need to reduce these datasets to the most functional elements so that we can discover key components driving disease pathogenesis. Current techniques to analyze omics data are not geared towards supporting parsimonious mechanistic models of bacterial pathogenesis. Specifically, networks are common tools for analyzing these data, however, at times graphical networks can be overwhelming due to the complexity of the information. In this work, we use sparse graphical networks to understand correlations between multiple high dimensional data sets. We apply these method to metabolomics and transcriptomics data from a recent animal model for Clostridium difficile infection in which mice were antibiotic treated with cefoperazone and challenged with C. difficile 2 days following treatment. We find significant changes in the C. difficile transcriptome at later time points compared to earlier time points; while, most of the amino acid changes occur in the first 24 hours. We find taurine, proline, 3-(4-hydroxyphenyl)lactate, 5-aminovalerate, 5-oxoproline, thioproline are the main amino acids contributing to the changes across all time points. Utilizing the key components of the sparse graphical model we’ve developed and analyzed an ordinary differential equation model to better understand the specific mechanisms leading to C. difficile colonization.
For successful ecosystem management and biodiversity conservation, in addition to ecological and evolutionary processes, we need to consider social and economic influences on the management target. Here, we introduce four models that address economic and social aspects of human society in the context of ecosystem management.
Populations are formed of their constituent interacting individuals, each with their own respective within-host biological processes. Infection not only spreads within the host organism but also spreads between individuals. We propose and study a nested multilevel model which links the within-host statuses of immunity and parasite density to population epidemiology under sub-lethal and lethal toxicant exposure. We analyse this nested model in order to better understand how toxicants impact the spread of disease within populations. We demonstrate that the outbreak of infection within a population is completely determined by the level of toxicant exposure, and that it is maximised by intermediate toxicant dosage. We classify the population epidemiology into 5 phases of increasing toxicant exposure and calculate the conditions under which disease will spread, showing that there exists a threshold toxicant level under which epidemics will not occur. In general, higher toxicant load results in either extinction of the population or outbreak of infection. The within-host statuses of the individual host also determine the outcome of the epidemic at the population level. We discuss applications of our model in the context of eco-epidemiology, particularly for bee colony losses, predicting that increased exposure to toxicants could result in more epidemics and therefore greater colony losses. We predict that reducing sub-lethal toxicant exposure below our predicted safe threshold could contribute to controlling population level disease and infection.
The deliberate release of Cyprinid herpes virus 3 (CyHV-3) to control invasive common carp (Cyprinus carpio) in the Murray-Darling Basin of south-eastern Australia has been dubbed `Carpageddon' and is highly controversial. Common carp now represent up to 90% of the biomass in invaded waterways and mortality rates of 70-80% have been observed during outbreaks of CyHV-3 in the northern hemisphere. Hence, the release of CyHV-3 carries the risk that decomposition following mass-mortality of carp may lead to severe oxygen depletion and subsequent anoxic events. I will present mathematical modelling that is the basis for (i) developing a low-risk release strategy, and (ii) advice to the Australian government on the long-term benefits of releasing the virus. The model couples a stage based demographic model with an SEIR-type model; recruitment of young carp is determined by available spawning habitat and the river system hydrology. We show that the long-term virus impact and the dynamics predicted by the model depend largely on the presence of a latent class that allows for virus reactivation and onward infection. The impact of CyHV-3 is also sensitive to where density-dependence occurs in the life cycle of carp because of the potential for virus mortality to either replace, or be in addition to, density-dependent mortality.
A vaccine is a biological preparation, providing long-term protection of the host to a disease by miming the interactions between the pathogen and the immune response. For decades vaccines have been the most effective means to fight and eradicate infectious diseases, from human diseases to animals diseases, including poultry and porcine viruses. One of the most significant endemic swine disease causing major animal and economic losses all over the world is porcine reproductive and respiratory syndrome (PRRS). Although different types of vaccines are on the market and various vaccination strategies have been trialled to control the disease, PRRS continues to persist, sometimes at high prevalence. The failure of some vaccination strategies to combat the disease leaves open questions regarding the major effect that the vaccine should have, as well as the most efficient vaccination strategy. We present an SIR model of a non-vaccinated population and extend it to a model where vaccination is applied. We use numerical simulations to evaluate vaccine efficacy for different types of vaccines and vaccination strategies to predict determinants for disease prevention. The study provides important guidelines for vaccine development and application in livestock populations.
The spectacular patterns of collective animal movement have been, and remain, a long standing and major interest in many branches of science, including biology, mathematics, physics and computational science. It is thought that the emergent patterns of coordinated motion are the consequence of individuals applying simple rules to adjust their velocity based on the relative locations and movements of nearby group members.
From the early 1980s onwards, the dominant methods for examining hypotheses about simple rules of interaction leading to coordinated, and sometimes complex, group motion have been discrete time self-propelled particle models. Common to many of these models are interaction rules chosen such that individuals will adjust their velocity to: avoid collisions with nearby neighbours, align their direction of motion with group members located at intermediate relative distances, and move towards other group members that are at relatively greater distances.
In the last decade advances in automated visual and GPS tracking methods have led to the exciting development of techniques for estimating the local rules of interaction used by real animals to coordinate collective motion directly from observational data. Analysis of tracking data, particularly of birds or fish in motion, suggests that the form of interaction rules chosen in self-propelled particle models is indeed plausible, with real animals adjusting their velocities consistent with collision avoidance at short range, matching directions of motion, and attraction to distant group mates. However, the mechanics of real interactions differ in detail to those adopted in models, even if the broad principles appear the same.
In this talk I will discuss an averaging technique for estimating how individual’s adjust their velocity as a function of the relative coordinates of their partners, and potentially other independent variables, directly from tracking data. Such techniques were first developed in concurrent studies by Katz et al. [1] and Herbert-Read et al. [2], and have since been further refined (see for example, [3]). I will then present some results that illustrate some of the commonalities and differences seen in the rules of interaction of different species of shoaling fish. I will discuss how randomisation methods can be used to identify when there are statistical differences in observed interaction rules of individuals belonging to different categories (for example hungry animals versus those that have recently fed to satiation). Finally, I will use these methods to examine differences in interaction rules between pairs of fish where one fish dominates leadership positions, and the other follows.
[1] Katz et al. PNAS, 108:18720-18725, (2011)
[2] Herbert-Read et al. PNAS, 108:18726-18731, (2011)
[3] Schaerf et al., Science Advances, 3:e1603201, (2017)
Transitive inference (TI) that uses known relationships to deduce unknown ones (using A > B and B > C to infer A > C given no direct interactions between A and C) to assess the opponent’s strength, or resource-holding potential (RHP), is widely reported in animals living in a group. This sounds counter-intuitive because the mechanism of TI seems to require social cognition and large memory capacity; individuals, in the TI mechanism, need abilities to identify others, observe contests among others and keep the results in memory. We examine the coevolution of memory and transitive inference by the evolutionary simulations, using the asymmetric hawk-dove game when a cost for losers is higher than a reward for winners. We found that the immediate inference strategy (II), which estimates the opponent's strength based on the past history of the direct fights, evolves with the large memory capacity, while the TI strategy, which estimates the unknown opponent's strength by transitive inference, evolves with the limited memory capacity. When a cost for losers is slightly higher than a reward for winners, the II strategy with the large memory capacity has an evolutionary advantage over the TI strategy with the limited memory capacity. It is because the direct fights are not so costly that more information about the fights leads to more accurate estimation of the opponent's strength and results in the accurate rank of the RHPs. When a cost for losers is much higher than a reward for winners, the TI strategy with the limited memory capacity has an evolutionary advantage. It is because a good way to avoid the costly fights is the prompt formation of the dominance hierarchy which does not necessarily reflect the actual rank of the RHPs; the TI strategy builds the dominance hierarchy much faster than the II strategy regardless of memory capacity, and the large amounts of information are not required for the TI strategy to form the dominance hierarchy promptly. Our study suggests that even smaller memory capacity is evolutionarily favoured in TI. The TI strategy tends to reinforce the hierarchy once it is built, regardless of whether it is consistent with RHP or not, because results of direct fights are always counted. Smaller memory capacity allows players to adjust the hierarchy well if it does not represent RHP. These results prove that TI can evolve in animals, which do not have the large memory capacity.
Multiple selective pressures drive evolution of malaria parasite, and determine its phenotypic traits, like virulence, transmissibility, drug resistance, as well as its population structure in host communities. Immune regulation plays an important part in this process, both within-host and on population level. The key to parasite survival within host is its antigenic variation, whereby parasite can avoid immune clearing by switching on and off multiple antigenic variants. In P. falciparum species, this process is controlled by a multigene var family, expressed on infected erythrocytes.
To study within-host parasite- immune interactions and the resulting selection processes we develop an agent based approach, that accommodates most salient features of malaria infection - RBC depletion, immune stimulation and clearing, antigenic variation.
In out setup, each parasite strain has a specific genetic makeup (collection of var genes), drawn from a large pool. Multiple strains of such quasi-species, compete within-host via cross-reactive immunity. Furthermore, mixed stains can recombine in mosquito agent to produce new types. Parasite diversity in host population can change in response to external factors, intensity of transmission and host immune status.
We apply our model to examine individual and population level outcomes, and quantify the relationship between environmental inputs (intensity of transmission, immune competence), and the evolution of parasite population structure.
Adipose tissue and adipocytes play a central role in the pathogenesis of metabolic diseases related to obesity. Size of fat cells depends on the balance of synthesis and mobilization of lipids and can undergo important variations throughout the life of the organism. These variations usually occur when storing and releasing lipids according to energy demand. This energy is stored as lipid droplets in their cytoplasm therefore adipocytes can adapt their size according to the lipid amount to be stored. Adipocyte size variation can reach one order of magnitude inside the same organism which is unique among cells. A striking feature in adipocytes size distribution is the lack of characteristic size since typical size distribution are bimodal. Since energy can be stored and retrieved and adipocytes are responsible for these lipid fluxes, we proposed a model of size-dependent lipid fluxes that exhibits bimodal equilibria with generic conditions. This model is a non-linear transport equation structured by the size of adipocytes and we provide a dedicated a stable numerical scheme solver able to capture the singular behaviour of solutions as concentration to Dirac masses. Using this solver, we are able to infer parameters related to lipid fluxes from size distribution. Performing recurrent surgical biopsies on rats, we measured the evolution of adipose cell size distribution for the same individual throughout the duration of feeding/refeeding experiments and use this data to assess fluxes parameters variations throughout the experiment.
In this presentation we consider a recent experimental dataset describing heat conduction in living porcine tissues. This novel dataset is important because porcine skin is similar to human skin, and improving our understanding of heat conduction in human skin is directly relevant to understanding burn injuries, which are common, painful and can require expensive treatments. A key feature of skin is that it is layered, with different thermal properties in different layers. Since the experimental dataset involves heat conduction in living tissues of anesthetized animals, an important experimental restriction is that the temperature within the skin is only measured at one location within the layered skin. Our aim is to determine whether this data is sufficient to infer the heat conduction parameters in layered skin, and we use a simplified two-layer model of heat conduction to mimic the experimental dataset. Using synthetic data at one location in the two-layer model, we explore whether it is possible to infer values of the thermal diffusivity in both layers. After this initial exploration, we then examine how our ability to infer the thermal diffusivities changes when we vary the location at which the experimental data is recorded, and in the hypothetical situation where we monitor two locations. Overall, we find that our ability to parameterize a model of heterogeneous heat conduction with limited experimental data is extremely sensitive to the location where data is collected, and our modelling results can be used to provide guidance about future experimental designs.
Circadian analysis is becoming increasingly important as a diagnostic tool to quantify deviations from regularity in circadian cycles. Circadian rhythms become less dominant and less regular with ageing and disease. It has been hypothesized that insomnia might be related to alterations, albeit small, in circadian and ultradian rhythms, but this topic remains an open problem. In this work, we propose a novel data-adaptive technique, singular spectrum analysis (SSA), to investigate in a model-free way quasiperiodic components and noise fluctuations in time series data. SSA was applied to one-week continuous actigraphy data in young adults with acute insomnia and healthy age-matched controls [1]. The findings suggest a small but significant delay in circadian components in the subjects with acute insomnia. The ultradian components follow a fractal $1/f$ power law for controls, however for individuals with acute insomnia this power law breaks down due to an increased variability at the 90min time scale. This is reminiscent of Kleitman's basic rest-activity (BRAC) cycles. It indicates that for healthy sleepers attention and activity can be sustained at time scale required by circumstances, while for individuals with acute insomnia this capacity may be impaired and they need to rest or switch activities in order to stay focused. Traditional methods of circadian rhythm analysis are unable to detect the more subtle effects of day-to-day variability and ultradian rhythm fragmentation at the specific 90min time scale.
[1] R Fossion, AL Rivera, JC Toledo-Roy, J Ellis and M Angelova. Multiscale adaptive analysis of circadian rhythms and intradaily variability: Application to actigraphy time series in acute insomnia subjects. PLOS One, 12(7): e0181762 (2017).
In Drosophila, circadian (~24h) behaviour is regulated by about 150 pacemaker neurons. To generate and maintain 24h rhythm, circadian gene expression in each pacemaker neurons is driven by transcriptional-translational feedback loop (TTFL). In TTFL of Drosophila, dCLOCK-CYCLE (dCLK-CYC) binds to the E-box (CACGTG) and activates the expression of period (per) and timeless (tim). Translated PER and TIM proteins repress their own transcription by removing dCLK-CYC dimer from the E-box. Interestingly, dCLK-Δ, which is a mutant of dCLK deleted amino acids (AA) 657-707 region and has impaired binding with PER, induces different effect on molecular rhythms depending on pacemaker neurons. For oscillation of PER, amplitude is largely reduced in ventral lateral neurons (LNvs), but not in dorsal neurons (DNs). However, how dCLK-Δ/CYC controlled TTFL operates differently in pacemaker neurons is unclear. To investigate this unexpected phenomenon, we established the mathematical model for the TTFL in Drosophila and predicted that pacemaker-neuron-dependent alteration of the molecular clockwork is caused by the difference of molecular composites between LNvs and DNs. This is confirmed by the follow up experiment. This work shows that clockworks at the molecular level have a critical role for specific functions of each pacemaker neurons.
Variability in electrophysiological properties, between different cells in a given heart and between the hearts of different members in a population, has a profound impact on deciding both the susceptibility to dangerous arrhythmias and the success or failure of anti-arrhythmic treatments. This variability also complicates the interpretation of both experimental and clinical data, and the predictions of in silico models. A key biomarker in the understanding of arrhythmias is the wavelength of excitation fronts, a measure of how much tissue remains refractory in the wake of an excitation impulse. However, this tissue-level biomarker is found to be unable to predict the susceptibility to wavebreaks that trigger fibrillation. It is therefore important to consider variability in underlying cell-level properties more directly. Using a novel combination of supervised learning and emulation, we are able to greatly reduce the computational costs involved with exploring variability in large numbers of parameters, hence identifying how differences in these properties impact on some key aspects of rotor-driven re-entries, including risk factors for the devolution into fibrillation.
Asthma is fundamentally a disease of airway constriction. Due to a variety of experimental challenges, the dynamics of airways are poorly understood. Of specific interest is the narrowing of the airway due to forces produced by the airway smooth muscle (ASM) wrapped around each airway. The interaction between the muscle and the airway wall is crucial for the airway constriction which occurs during an asthma attack. While crossbridge theory is a well-studied representation of complex smooth muscle dynamics, and these dynamics can be coupled to the airway wall, this comes at significant computational cost, even for isolated airways. Because many phenomena of interest in pulmonary physiology cannot be adequately understood by studying isolated airways, this presents a significant limitation. We present results associated with the development of an approximated method (distribution moment approximation) which provides orders of magnitude reduction in computational complexity whilst retaining rich ASM dynamics. This method is then used to model physiological processes such as airway contraction, bronchodilatory effect of a deep inspiration and preliminary results associated with heterogeneity of branched airways.
The elderly population is particularly susceptible to infectious diseases, such as influenza, as evidenced by increased occurrence and severity of infection, and reduced vaccine efficacy. Studies suggest an association between persistent cytomegalovirus (CMV) infection and deficiencies of the aged immune system. CMV is a common herpes virus that infects up to 85% of the population aged 40 years or over. Healthy individuals rarely experience symptoms after the initial CMV infection. However, the virus remains in the body in an inactive state and constant surveillance by the immune system is required to detect and suppress reactivation of the virus. It is unclear whether this lifelong burden on the immune system to control persistent CMV infection plays a causal role in age-related compromise of the immune system to fight infections and respond to vaccines.
The effect of age on the immune system is multifactorial. One factor is a decline in the diversity of CD8+ “killer” T cells recruited by the immune system to fight a particular infection. In this study we investigated whether the presence of lifelong CMV infection impacts the diversity of T cells responding to another infection in aged mice. The diversity of T cells responding to infection was assessed using single-cell Sanger sequencing and bioinformatics analysis of the T cell receptors, which are expressed on the surface of individual T cells to enable the detection of virus. Various features of the T cell receptors, including diversity, gene usage, sequence length, and conserved patterns of amino acid usage were quantitatively compared between three groups of mice: young adult and old adult mice with no CMV infection and old mice with lifelong CMV infection. T cell diversity was quantified using ecological measures of diversity and rarefaction to standardize for differences in sample sizes. Our results demonstrate that, in the absence of persistent CMV infection, T cell responses in old mice were significantly less diverse than in young adult mice. In contrast, T cell diversity was maintained in the old mice with CMV infection and we did not observe the age-related narrowing of the T cell population that occurred in the absence of persistent CMV infection. Next-generation sequencing of the broader T cell populations that have not previously encountered infection showed similar T cell diversity in all three groups of mice, including the presence of T cell receptors that were recruited into the immune response in CMV-infected old mice but not in old mice without CMV infection.
These results suggest that lifelong CMV infection may actually improve T cell immunity in later life by broadening the diversity of T cells that can be effectively recruited and/or maintained in immune responses to other infections. Furthermore, these results challenge the paradigm that CMV has a negative impact on the ageing immune system.
Recent studies show that human immune system interacts with intestinal microbiome. A certain group of intestinal bacteria is known to be possible suppressor of undesirable immune response [1]. They produce short chain fatty acid such as butyrate and induce regulatory T cells [1]. Regulatory T cells are suppressor of exaggerated immune responses [2]. Therefore, collapse of the ecological balance of intestinal microbiome may cause dysregulation of immune system. On the other hand, mutant mice with regulatory T cells that have limited diversity of T cell receptor repertoire cause inflammation to commensal intestinal microbiome driven by T helper cells [3], suggesting that immune cells develop hypersensitive response to commensal microbiome and that regulatory T cells suppress the reactions.
This suggests us that intervention targeting intestinal microbiome may provide a novel strategy for treatment of allergic symptoms. In this talk, we develop a simple mathematical model that describes interaction between intestinal microbiome and immune system. We consider dynamics of T helper cells (trigger of allergic response), regulatory T cells (suppressor of allergy), and intestinal microbiome. Differentiation process of these two types of T cells is described as a model of immune system, the structure of which was developed in our previous study [4].
Analysis shows that the model has results of three different types differing in stable steady states: (1) the case with a single stable positive steady state representing “healthy condition”, (2) the case with a single stable steady state representing “allergic condition”, and (3) the case with two stable steady states, representing “healthy” and “allergic” conditions, respectively. Steady states representing healthy condition has a low level of helper T cells, a high level of regulatory T cells, and bacteria; whilst those representing allergic condition has opposite features (high helper T cells, low regulatory T cells and bacteria). Considering these three different situations, now we have two strategies for allergy treatment. The first strategy is to eliminate steady state which corresponds to allergy, by realizing the condition for case (1). We derived a formula for the absence of allergic state. The second strategy is to reduce the abundance of helper T cells at the steady state. We tested bacteria-related parameters for their effectiveness in suppressing the level of helper T cells. This study is the first attempt in modelling association between intestinal microbiome and immune system and in finding out a therapy method for allergy by intervention of intestinal microbiome.
[1] Furusawa, Y. et al. 2013. Nature 504, 446-50.
[2] Sakaguchi, S et al. 2008. Cell 133, 775-787.
[3] Nishio, J. et al. 2015. Proc Natl Acad Sci U S A 112, 12770-5.
[4] Hara, A and Iwasa, Y. 2017. J. Theor. Biol.425, 23-42.
Tuberculosis (TB), a deadly infectious disease caused by the bacterium Mycobacterium tuberculosis (Mtb). The disease is characterized by the development of granulomas consisting of immune cells that form a cluster around the bacteria to limit bacterial growth and disease outcomes. Control of the TB epidemic is limited by a complicated drug regimen, development of antibiotic resistance, and the lack of an effective vaccine against infection and disease.
There are a range of pro- and anti-inflammatory molecules associated with Mtb infection and it has been shown that depletion of a specific anti-inflammatory molecule, namely IL-10, is associated with improved vaccine efficacy in mice. However, there are conflicting results in both mice and monkeys as to the effects of IL-10 depletion on disease, with some studies showing that the depletion of anti-inflammatory factors improves disease outcome and others demonstrating that this depletion causes more severe disease. In previously published modelling studies, we suggest that depleting IL-10 can increase the probability that a granuloma is cleared in non-human primates (NHPs).
Using GranSim, our previously published hybrid, multiscale, agent-based model of granuloma formation following Mtb infection, we explore the role of IL-10, on granuloma development and function. In this exploratory study, we simulate the depletion of IL-10 from two discrete time points. The first time point coincides with the infection time of Mtb, and allows us to examine the role of IL-10 in granuloma formation. Elucidating the mechanisms by which IL-10 affects granuloma formation allows us explore whether the improvement in vaccine efficacy that has been identified in mice can also be observed in our model. The second time point is at 200 days post infection when a stable granuloma is already established to allow us to identify the mechanisms by which IL-10 is acting upon granuloma function and maintenance. Understanding the effect of IL-10 depletion on established TB granulomas is important to identify new therapeutic pathways for this disease.
Understanding the mechanisms of HIV latency is important in the development of strategies for managing infection. Time from infection until production of virus has been shown to vary among infected cells, hence challenging the dichotomous assumption that cells are either latent or productively infected at time of infection. In this paper, we will explore the implications of an alternative hypothesis that rates of activation follow a probability distribution, of which latency is just an extreme of this spectrum. We show the emerging dynamics from this mathematical model and test its ability to explain features or reservoir formation and decay observed in SIV data. Analysis of SIV DNA levels in macaques that start treatment at different times shows that the decay rates of cells are statistically different in early and late initiation of treatment. Modelling suggests that data can be explained by a spectrum of reactivation rates of infected cells.
Experiments in rhesus macaques have shown that for simian immunodeficiency virus the size of the viral inoculum and the infection stage of the donor animal alter the likelihood of establishment of infection. In this study, we postulate a role for the host and donor antibodies in explaining the dependence of infectiousness on donor infection stage. The resulting mathematical model exhibits bistable dynamics, with viral clearance and persistence dependent on initial conditions. We fit the model to temporal virus data using a censored data approach for measurements below the limit of detection, and make predictions about the minimal viral load in the inoculum required for persistent infection.