A multiscale mathematical model describing the growth and development of bambara groundnut

A principal objective in agriculture is to maximise food production; this is particularly relevant with the added demands of an ever increasing population, coupled with the unpredictability that climate change brings. Further improvements in productivity can only be achieved with an increased understanding of plant and crop processes. In this respect, mathematical modelling of plants and crops plays an important role. In this paper we present a two-scale mathematical model of crop yield that accounts for plant growth and canopy interactions. A system of nonlinear ordinary differential equations (ODEs) is formulated to describe the growth of each individual plant, where equations are coupled via a term that describes plant competition via canopy-canopy interactions. A crop of greenhouse plants is then modelled via an agent based modelling approach in which the growth of each plant is described via our system of ODEs. The model is formulated for the African drought tolerant legume bambara groundnut (Vigna subterranea), which is currently being investigated as a food source in light of climate change and food insecurity challenges. Our model allows us to account for plant diversity and also investigate the effect of individual plant traits (e.g. plant canopy size and planting distance) on the yield of the overall crop. Informed with greenhouse data, model results show that plant positioning relative to other plants has a large impact on individual plant yield. Variation in physiological plant traits from genetic diversity and the environmental effects lead to experimentally observed variations in crop yield. These traits include plant height, plant carrying capacity, leaf accumulation rate and canopy spread. Of these traits plant height and ground cover growth rates are found to have the greatest impact on crop yield. We also consider a range of different planting arrangements (uniform grid, staggered grid, circular rings and random allocation) and find that the staggered grid leads to the greatest crop yield (6% more compared to uniform grid). Whilst formulated specifically for bambara groundnut, the generic formulation of our model means that with changes to certain parameter's, it may be extended to other crop species that form a canopy.

Integrating protein networks and machine learning for disease stratification in the Hereditary Spastic Paraplegias

The Hereditary Spastic Paraplegias are a group of neurodegenerative diseases characterized by spasticity and weakness in the lower body. Owing to the combination of genetic diversity and variable clinical presentation, the Hereditary Spastic Paraplegias are a strong candidate for protein-protein interaction network analysis as a tool to understand disease mechanism(s) and to aid functional stratification of phenotypes. In this study, experimentally validated human data were used to create a protein-protein interaction network based on the causative genes. Network evaluation as a combination of topological analysis and functional annotation led to the identification of core proteins in putative shared biological processes, such as intracellular transport and vesicle trafficking. The application of machine learning techniques suggested a functional dichotomy linked with distinct sets of clinical presentations, indicating that there is scope to further classify conditions currently described under the same umbrella-term of Hereditary Spastic Paraplegias based on specific molecular mechanisms of disease.

A model of the PI cycle reveals the regulating roles of lipid-binding proteins and pitfalls of using mosaic biological data

The phosphatidylinositol (PI) cycle is central to eukaryotic cell signaling. Its complexity, due to the number of reactions and lipid and inositol phosphate intermediates involved makes it difficult to analyze experimentally. Computational modelling approaches are seen as a way forward to elucidate complex biological regulatory mechanisms when this cannot be achieved solely through experimental approaches. Whilst mathematical modelling is well established in informing biological systems, many models are often informed by data sourced from multiple unrelated cell types (mosaic data) or from purified enzyme data. In this work, we develop a model of the PI cycle informed by experimental and omics data taken from a single cell type, namely platelets. We were able to make a number of predictions regarding the regulation of PI cycle enzymes, the importance of the number of receptors required for successful GPCR signaling and the importance of lipid- and protein-binding proteins in regulating second messenger outputs. We then consider how pathway behavior differs, when fully informed by data for HeLa cells and show that model predictions remain consistent. However, when informed by mosaic experimental data model predictions greatly vary illustrating the risks of using mosaic datasets from unrelated cell types.

Mathematical modelling of contact dermatitis from nickel and chromium

Dermal exposure to metal allergens can lead to irritant and allergic contact dermatitis (ACD). In this paper we present a mathematical model of the absorption of metal ions, hexavalent chromium and nickel, into the viable epidermis and compare the localised irritant and T-lymphocyte (T-cell) mediated immune responses. The model accounts for the spatial-temporal variation of skin health, extra and intracellular allergen concentrations, innate immune cells, T-cells, cytokine signalling and lymph node activity up to about 6 days after contact with these metals; repair processes associated with withdrawal of exposure to both metals is not considered in the current model, being assumed secondary during the initial phases of exposure. Simulations of the resulting system of PDEs are studied in one-dimension, i.e. across skin depth, and three-dimensional scenarios with the aim of comparing the responses to the two ions in the cases of first contact (no T-cells initially present) and second contact (T-cells initially present). The results show that on continuous contact, chromium ions elicit stronger skin inflammation, but for nickel, subsequent re-exposure stimulates stronger responses due to an accumulation of cytotoxic T-cell mediated responses which characterise ACD. Furthermore, the surface area of contact to these metals has little effect on the speed of response, whilst sensitivity is predicted to increase with the thickness of skin. The modelling approach is generic and should be applicable to describe contact dermatitis from a wide range of allergens.

Best Practices to Maximize the Use and Reuse of Quantitative and Systems Pharmacology Models: Recommendations From the United Kingdom Quantitative and Systems Pharmacology Network

The lack of standardization in the way that quantitative and systems pharmacology (QSP) models are developed, tested, and documented hinders their reproducibility, reusability, and expansion or reduction to alternative contexts. This in turn undermines the potential impact of QSP in academic, industrial, and regulatory frameworks. This article presents a minimum set of recommendations from the UK Quantitative and Systems Pharmacology Network (UK QSP Network) to guide QSP practitioners seeking to maximize their impact, and stakeholders considering the use of QSP models in their environment.

Multi-scale, whole-system models of liver metabolic adaptation to fat and sugar in non-alcoholic fatty liver disease

Non-alcoholic fatty liver disease (NAFLD) is a serious public health issue associated with high fat, high sugar diets. However, the molecular mechanisms mediating NAFLD pathogenesis are only partially understood. Here we adopt an iterative multi-scale, systems biology approach coupled to in vitro experimentation to investigate the roles of sugar and fat metabolism in NAFLD pathogenesis. The use of fructose as a sweetening agent is controversial; to explore this, we developed a predictive model of human monosaccharide transport, signalling and metabolism. The resulting quantitative model comprising a kinetic model describing monosaccharide transport and insulin signalling integrated with a hepatocyte-specific genome-scale metabolic network (GSMN). Differential kinetics for the utilisation of glucose and fructose were predicted, but the resultant triacylglycerol production was predicted to be similar for monosaccharides; these predictions were verified by in vitro data. The role of physiological adaptation to lipid overload was explored through the comprehensive reconstruction of the peroxisome proliferator activated receptor alpha (PPARα) regulome integrated with a hepatocyte-specific GSMN. The resulting qualitative model reproduced metabolic responses to increased fatty acid levels and mimicked lipid loading in vitro. The model predicted that activation of PPARα by lipids produces a biphasic response, which initially exacerbates steatosis. Our data support the evidence that it is the quantity of sugar rather than the type that is critical in driving the steatotic response. Furthermore, we predict PPARα-mediated adaptations to hepatic lipid overload, shedding light on potential challenges for the use of PPARα agonists to treat NAFLD.

Model reduction in mathematical pharmacology : Integration, reduction and linking of PBPK and systems biology models

In this paper we present a framework for the reduction and linking of physiologically based pharmacokinetic (PBPK) models with models of systems biology to describe the effects of drug administration across multiple scales. To address the issue of model complexity, we propose the reduction of each type of model separately prior to being linked. We highlight the use of balanced truncation in reducing the linear components of PBPK models, whilst proper lumping is shown to be efficient in reducing typically nonlinear systems biology type models. The overall methodology is demonstrated via two example systems; a model of bacterial chemotactic signalling in Escherichia coli and a model of extracellular regulatory kinase activation mediated via the extracellular growth factor and nerve growth factor receptor pathways. Each system is tested under the simulated administration of three hypothetical compounds; a strong base, a weak base, and an acid, mirroring the parameterisation of pindolol, midazolam, and thiopental, respectively. Our method can produce up to an 80% decrease in simulation time, allowing substantial speed-up for computationally intensive applications including parameter fitting or agent based modelling. The approach provides a straightforward means to construct simplified Quantitative Systems Pharmacology models that still provide significant insight into the mechanisms of drug action. Such a framework can potentially bridge pre-clinical and clinical modelling - providing an intermediate level of model granularity between classical, empirical approaches and mechanistic systems describing the molecular scale.

Mathematical Analysis of the Escherichia coli Chemotaxis Signalling Pathway

We undertake a detailed mathematical analysis of a recent nonlinear ordinary differential equation (ODE) model describing the chemotactic signalling cascade within an Escherichia coli cell. The model includes a detailed description of the cell signalling cascade and an average approximation of the receptor activity. A steady-state stability analysis reveals the system exhibits one positive real steady state which is shown to be asymptotically stable. Given the occurrence of a negative feedback between phosphorylated CheB (CheB-P) and the receptor state, we ask under what conditions the system may exhibit oscillatory-type behaviour. A detailed analysis of parameter space reveals that whilst variation in kinetic rate parameters within known biological limits is unlikely to lead to such behaviour, changes in the total concentration of the signalling proteins do. We postulate that experimentally observed overshoot behaviour can actually be described by damped oscillatory dynamics and consider the relationship between overshoot amplitude, total cell protein concentration and the magnitude of the external ligand stimulus. Model reductions in the full ODE model allow us to understand the link between phosphorylation events and the negative feedback between CheB-P and receptor methylation, as well as elucidate why some mathematical models exhibit overshoot and others do not. Our paper closes by discussing intercell variability of total protein concentration as a means of ensuring the overall survival of a population as cells are subjected to different environments.

Physiologically-based pharmacokinetic and toxicokinetic models for estimating human exposure to five toxic elements through oral ingestion

Biological monitoring and physiologically-based pharmacokinetic (PBPK) modelling are useful complementary tools in quantifying human exposure to elements in the environment. In this work, we used PBPK models to determine the optimal time for collecting biological samples in a longitudinal study to determine if participants who consumed allotment produce had been exposed to arsenic, cadmium, chromium, nickel or lead. There are a number of PBPK models for these elements published in the literature, which vary in size, complexity and application, given the differences in physiochemical properties of the elements, organs involved in metabolism and exposure pathways affected. We selected PBPK models from the literature to simulate the oral ingestion pathway from consumption of allotment produce. Some models required modification by reducing or removing selected compartments whilst still maintaining their original predictability. The performance of the modified models was evaluated by comparing the predicted urinary and blood elemental levels with experimental data and other model simulations published in the literature. Overall, the model predictions were consistent with literature data (r > 0.7, p < 0.05), and were influential in predicting when samples should be collected. Our results demonstrate the use of mathematical modelling in informing and optimising the design of longitudinal studies.

Methods of Model Reduction for Large-Scale Biological Systems: A Survey of Current Methods and Trends

Complex models of biochemical reaction systems have become increasingly common in the systems biology literature. The complexity of such models can present a number of obstacles for their practical use, often making problems difficult to intuit or computationally intractable. Methods of model reduction can be employed to alleviate the issue of complexity by seeking to eliminate those portions of a reaction network that have little or no effect upon the outcomes of interest, hence yielding simplified systems that retain an accurate predictive capacity. This review paper seeks to provide a brief overview of a range of such methods and their application in the context of biochemical reaction network models. To achieve this, we provide a brief mathematical account of the main methods including timescale exploitation approaches, reduction via sensitivity analysis, optimisation methods, lumping, and singular value decomposition-based approaches. Methods are reviewed in the context of large-scale systems biology type models, and future areas of research are briefly discussed.

A combined model reduction algorithm for controlled biochemical systems

BACKGROUND

Systems Biology continues to produce increasingly large models of complex biochemical reaction networks. In applications requiring, for example, parameter estimation, the use of agent-based modelling approaches, or real-time simulation, this growing model complexity can present a significant hurdle. Often, however, not all portions of a model are of equal interest in a given setting. In such situations methods of model reduction offer one possible approach for addressing the issue of complexity by seeking to eliminate those portions of a pathway that can be shown to have the least effect upon the properties of interest.

METHODS

In this paper a model reduction algorithm bringing together the complementary aspects of proper lumping and empirical balanced truncation is presented. Additional contributions include the development of a criterion for the selection of state-variable elimination via conservation analysis and use of an 'averaged' lumping inverse. This combined algorithm is highly automatable and of particular applicability in the context of 'controlled' biochemical networks.

RESULTS

The algorithm is demonstrated here via application to two examples; an 11 dimensional model of bacterial chemotaxis in Escherichia coli and a 99 dimensional model of extracellular regulatory kinase activation (ERK) mediated via the epidermal growth factor (EGF) and nerve growth factor (NGF) receptor pathways. In the case of the chemotaxis model the algorithm was able to reduce the model to 2 state-variables producing a maximal relative error between the dynamics of the original and reduced models of only 2.8% whilst yielding a 26 fold speed up in simulation time. For the ERK activation model the algorithm was able to reduce the system to 7 state-variables, incurring a maximal relative error of 4.8%, and producing an approximately 10 fold speed up in the rate of simulation. Indices of controllability and observability are additionally developed and demonstrated throughout the paper. These provide insight into the relative importance of individual reactants in mediating a biochemical system's input-output response even for highly complex networks.

CONCLUSIONS

Through application, this paper demonstrates that combined model reduction methods can produce a significant simplification of complex Systems Biology models whilst retaining a high degree of predictive accuracy. In particular, it is shown that by combining the methods of proper lumping and empirical balanced truncation it is often possible to produce more accurate reductions than can be obtained by the use of either method in isolation.

A high-density immunoblotting methodology for quantification of total protein levels and phosphorylation modifications

The components of many signaling pathways have been identified and there is now a need to conduct quantitative data-rich temporal experiments for systems biology and modeling approaches to better understand pathway dynamics and regulation. Here we present a modified Western blotting method that allows the rapid and reproducible quantification and analysis of hundreds of data points per day on proteins and their phosphorylation state at individual sites. The approach is of particular use where samples show a high degree of sample-to-sample variability such as primary cells from multiple donors. We present a case study on the analysis of >800 phosphorylation data points from three phosphorylation sites in three signaling proteins over multiple time points from platelets isolated from ten donors, demonstrating the technique's potential to determine kinetic and regulatory information from limited cell numbers and to investigate signaling variation within a population. We envisage the approach being of use in the analysis of many cellular processes such as signaling pathway dynamics to identify regulatory feedback loops and the investigation of potential drug/inhibitor responses, using primary cells and tissues, to generate information about how a cell's physiological state changes over time.

Understanding the link between single cell and population scale responses of Escherichia coli in differing ligand gradients

We formulate an agent-based population model of Escherichia coli cells which incorporates a description of the chemotaxis signalling cascade at the single cell scale. The model is used to gain insight into the link between the signalling cascade dynamics and the overall population response to differing chemoattractant gradients. Firstly, we consider how the observed variation in total (phosphorylated and unphosphorylated) signalling protein concentration affects the ability of cells to accumulate in differing chemoattractant gradients. Results reveal that a variation in total cell protein concentration between cells may be a mechanism for the survival of cell colonies across a wide range of differing environments. We then study the response of cells in the presence of two different chemoattractants. In doing so we demonstrate that the population scale response depends not on the absolute concentration of each chemoattractant but on the sensitivity of the chemoreceptors to their respective concentrations. Our results show the clear link between single cell features and the overall environment in which cells reside.

Investigating flavonoids as molecular templates for the design of small-molecule inhibitors of cell signaling

Epidemiological and clinical trials reveal compelling evidence for the ability of dietary flavonoids to lower cardiovascular disease risk. The mechanisms of action of these polyphenolic compounds are diverse, and of particular interest is their ability to function as protein and lipid kinase inhibitors. We have previously described structure-activity studies that reinforce the possibility for using flavonoid structures as templates for drug design. In the present study, we aim to begin constructing rational screening strategies for exploiting these compounds as templates for the design of clinically relevant, antiplatelet agents. We used the platelet as a model system to dissect the structural influence of flavonoids, stilbenes, anthocyanidins, and phenolic acids on inhibition of cell signaling and function. Functional groups identified as relevant for potent inhibition of platelet function included at least 2 benzene rings, a hydroxylated B ring, a planar C ring, a C ring ketone group, and a C-2 positioned B ring. Hydroxylation of the B ring with either a catechol group or a single C-4' hydroxyl may be required for efficient inhibition of collagen-stimulated tyrosine phosphorylated proteins of 125 to 130 kDa, but may not be necessary for that of phosphotyrosine proteins at approximately 29 kDa. The removal of the C ring C-3 hydroxyl together with a hydroxylated B ring (apigenin) may confer selectivity for 37 to 38 kDa phosphotyrosine proteins. We conclude that this study may form the basis for construction of maps of flavonoid inhibitory activity on kinase targets that may allow a multitargeted therapeutic approach with analogue counterparts and parent compounds.

Modelling negative feedback networks for activating transcription factor 3 predicts a dominant role for miRNAs in immediate early gene regulation

Activating transcription factor 3 (Atf3) is rapidly and transiently upregulated in numerous systems, and is associated with various disease states. Atf3 is required for negative feedback regulation of other genes, but is itself subject to negative feedback regulation possibly by autorepression. In cardiomyocytes, Atf3 and Egr1 mRNAs are upregulated via ERK1/2 signalling and Atf3 suppresses Egr1 expression. We previously developed a mathematical model for the Atf3-Egr1 system. Here, we adjusted and extended the model to explore mechanisms of Atf3 feedback regulation. Introduction of an autorepressive loop for Atf3 tuned down its expression and inhibition of Egr1 was lost, demonstrating that negative feedback regulation of Atf3 by Atf3 itself is implausible in this context. Experimentally, signals downstream from ERK1/2 suppress Atf3 expression. Mathematical modelling indicated that this cannot occur by phosphorylation of pre-existing inhibitory transcriptional regulators because the time delay is too short. De novo synthesis of an inhibitory transcription factor (ITF) with a high affinity for the Atf3 promoter could suppress Atf3 expression, but (as with the Atf3 autorepression loop) inhibition of Egr1 was lost. Developing the model to include newly-synthesised miRNAs very efficiently terminated Atf3 protein expression and, with a 4-fold increase in the rate of degradation of mRNA from the mRNA/miRNA complex, profiles for Atf3 mRNA, Atf3 protein and Egr1 mRNA approximated to the experimental data. Combining the ITF model with that of the miRNA did not improve the profiles suggesting that miRNAs are likely to play a dominant role in switching off Atf3 expression post-induction.

Activating transcription factor 3 (Atf3) is an important regulatory transcription factor which is associated with inflammation, restraint of the immune response and cancer. In this work, we develop a series of mathematical models to understand how Atf3 may be regulated. Informed with data from the literature and our own experiments, we show that self-regulation of Atf3 does not allow for variation between experimentally observed Atf3 mRNA and Atf3 protein expression profiles. A fast-acting signal via phosphorylated RSK is also shown to be implausible for similar reasons. Extending our mathematical model further, we postulate for the first time, that the observed dynamical variation in Atf3 mRNA and protein can be described by microRNAs downstream of RSKs. The further inclusion of an inhibitory transcription factor for Atf3 expression has little effect on these findings.

Fold-change detection in a whole-pathway model of Escherichia coli chemotaxis

There has been recent interest in sensory systems that are able to display a response which is proportional to a fold change in stimulus concentration, a feature referred to as fold-change detection (FCD). Here, we demonstrate FCD in a recent whole-pathway mathematical model of Escherichia coli chemotaxis. FCD is shown to hold for each protein in the signalling cascade and to be robust to kinetic rate and protein concentration variation. Using a sensitivity analysis, we find that only variations in the number of receptors within a signalling team lead to the model not exhibiting FCD. We also discuss the ability of a cell with multiple receptor types to display FCD and explain how a particular receptor configuration may be used to elucidate the two experimentally determined regimes of FCD behaviour. All findings are discussed in respect of the experimental literature.

A mathematical model of the sterol regulatory element binding protein 2 cholesterol biosynthesis pathway

Cholesterol is one of the key constituents for maintaining the cellular membrane and thus the integrity of the cell itself. In contrast high levels of cholesterol in the blood are known to be a major risk factor in the development of cardiovascular disease. We formulate a deterministic nonlinear ordinary differential equation model of the sterol regulatory element binding protein 2 (SREBP-2) cholesterol genetic regulatory pathway in a hepatocyte. The mathematical model includes a description of genetic transcription by SREBP-2 which is subsequently translated to mRNA leading to the formation of 3-hydroxy-3-methylglutaryl coenzyme A reductase (HMGCR), a main regulator of cholesterol synthesis. Cholesterol synthesis subsequently leads to the regulation of SREBP-2 via a negative feedback formulation. Parameterised with data from the literature, the model is used to understand how SREBP-2 transcription and regulation affects cellular cholesterol concentration. Model stability analysis shows that the only positive steady-state of the system exhibits purely oscillatory, damped oscillatory or monotic behaviour under certain parameter conditions. In light of our findings we postulate how cholesterol homeostasis is maintained within the cell and the advantages of our model formulation are discussed with respect to other models of genetic regulation within the literature.

•

We formulate and analyse a nonlinear ODE model of the SREBP2 pathway.

•

The mathematical model exhibits stable limit cycles under certain parameter conditions.

•

Negative feedbacks in the SREBP2 pathway may help regulate cholesterol homeostasis.

•

Our model provides a more accurate formulation of genetic regulation using nonlinear ODEs.

Response kinetics in the complex chemotaxis signalling pathway of Rhodobacter sphaeroides

Chemotaxis is one of the best-characterized signalling systems in biology. It is the mechanism by which bacteria move towards optimal environments and is implicated in biofilm formation, pathogenesis and symbiosis. The properties of the bacterial chemosensory response have been described in detail for the single chemosensory pathway of Escherichia coli. We have characterized the properties of the chemosensory response of Rhodobacter sphaeroides, an α-proteobacterium with multiple chemotaxis pathways, under two growth conditions allowing the effects of protein expression levels and cell architecture to be investigated. Using tethered cell assays, we measured the responses of the system to step changes in concentration of the attractant propionate and show that, independently of the growth conditions, R. sphaeroides is chemotactic over at least five orders of magnitude and has a sensing profile following Weber's Law. Mathematical modelling also shows that, as E. coli, R. sphaeroides is capable of showing fold-change detection (FCD). Our results indicate that general features of bacterial chemotaxis such as the range and sensitivity of detection, adaptation times, adherence to Weber's Law and the presence of FCD may be integral features of chemotaxis systems in general, regardless of network complexity, protein expression levels and cellular architecture across different species.

Theoretical insights into bacterial chemotaxis

Research into understanding bacterial chemotactic systems has become a paradigm for Systems Biology. Experimental and theoretical researchers have worked hand-in-hand for over 40 years to understand the intricate behavior driving bacterial species, in particular how such small creatures, usually not more than 5 µm in length, detect and respond to small changes in their extracellular environment. In this review we highlight the importance that theoretical modeling has played in providing new insight and understanding into bacterial chemotaxis. We begin with an overview of the bacterial chemotaxis sensory response, before reviewing the role of theoretical modeling in understanding elements of the system on the single cell scale and features underpinning multiscale extensions to population models.

Classifying general nonlinear force laws in cell-based models via the continuum limit

Although discrete cell-based frameworks are now commonly used to simulate a whole range of biological phenomena, it is typically not obvious how the numerous different types of model are related to one another, nor which one is most appropriate in a given context. Here we demonstrate how individual cell movement on the discrete scale modeled using nonlinear force laws can be described by nonlinear diffusion coefficients on the continuum scale. A general relationship between nonlinear force laws and their respective diffusion coefficients is derived in one spatial dimension and, subsequently, a range of particular examples is considered. For each case excellent agreement is observed between numerical solutions of the discrete and corresponding continuum models. Three case studies are considered in which we demonstrate how the derived nonlinear diffusion coefficients can be used to (a) relate different discrete models of cell behavior; (b) derive discrete, intercell force laws from previously posed diffusion coefficients, and (c) describe aggregative behavior in discrete simulations.

Modelling acidosis and the cell cycle in multicellular tumour spheroids

A partial differential equation model is developed to understand the effect that nutrient and acidosis have on the distribution of proliferating and quiescent cells and dead cell material (necrotic and apoptotic) within a multicellular tumour spheroid. The rates of cell quiescence and necrosis depend upon the local nutrient and acid concentrations and quiescent cells are assumed to consume less nutrient and produce less acid than proliferating cells. Analysis of the differences in nutrient consumption and acid production by quiescent and proliferating cells shows low nutrient levels do not necessarily lead to increased acid concentration via anaerobic metabolism. Rather, it is the balance between proliferating and quiescent cells within the tumour which is important; decreased nutrient levels lead to more quiescent cells, which produce less acid than proliferating cells. We examine this effect via a sensitivity analysis which also includes a quantification of the effect that nutrient and acid concentrations have on the rates of cell quiescence and necrosis.

Comparing a discrete and continuum model of the intestinal crypt

The integration of processes at different scales is a key problem in the modelling of cell populations. Owing to increased computational resources and the accumulation of data at the cellular and subcellular scales, the use of discrete, cell-level models, which are typically solved using numerical simulations, has become prominent. One of the merits of this approach is that important biological factors, such as cell heterogeneity and noise, can be easily incorporated. However, it can be difficult to efficiently draw generalizations from the simulation results, as, often, many simulation runs are required to investigate model behaviour in typically large parameter spaces. In some cases, discrete cell-level models can be coarse-grained, yielding continuum models whose analysis can lead to the development of insight into the underlying simulations. In this paper we apply such an approach to the case of a discrete model of cell dynamics in the intestinal crypt. An analysis of the resulting continuum model demonstrates that there is a limited region of parameter space within which steady-state (and hence biologically realistic) solutions exist. Continuum model predictions show good agreement with corresponding results from the underlying simulations and experimental data taken from murine intestinal crypts.

Modeling chemotaxis reveals the role of reversed phosphotransfer and a bi-functional kinase-phosphatase

Understanding how multiple signals are integrated in living cells to produce a balanced response is a major challenge in biology. Two-component signal transduction pathways, such as bacterial chemotaxis, comprise histidine protein kinases (HPKs) and response regulators (RRs). These are used to sense and respond to changes in the environment. Rhodobacter sphaeroides has a complex chemosensory network with two signaling clusters, each containing a HPK, CheA. Here we demonstrate, using a mathematical model, how the outputs of the two signaling clusters may be integrated. We use our mathematical model supported by experimental data to predict that: (1) the main RR controlling flagellar rotation, CheY₆, aided by its specific phosphatase, the bifunctional kinase CheA₃, acts as a phosphate sink for the other RRs; and (2) a phosphorelay pathway involving CheB₂ connects the cytoplasmic cluster kinase CheA₃ with the polar localised kinase CheA₂, and allows CheA₃-P to phosphorylate non-cognate chemotaxis RRs. These two mechanisms enable the bifunctional kinase/phosphatase activity of CheA₃ to integrate and tune the sensory output of each signaling cluster to produce a balanced response. The signal integration mechanisms identified here may be widely used by other bacteria, since like R. sphaeroides, over 50% of chemotactic bacteria have multiple cheA homologues and need to integrate signals from different sources.

Chemotactic bacteria sense nutrient gradients and swim towards better environments for growth. A cluster of receptors in the cell membrane detects nutrient levels and signals via a cytoplasmic signaling pathway to the flagellum. The complexity of this signaling pathway varies across different bacterial species. The relatively simple pathway used by Escherichia coli is well understood; however, many bacteria, for example Rhodobacter sphaeroides, have more sophisticated pathways that, as well as being able to detect nutrients, are also able to assess the metabolic state of the cell. The receptors that detect metabolic state are located within an additional cluster that is physically distinct from the one that senses nutrients. In this work, we use a combination of experimentation and mathematical modeling to gain insight into the complex decision-making mechanisms that enable bacteria to weigh-up different stimuli and decide upon an appropriate response. We find novel communication mechanisms between the two signaling clusters that allow the outputs of the signaling pathways to be balanced and tuned according to the prevailing environmental conditions. The signaling principles identified here are likely to be used in other complex sensory networks.

From a discrete to a continuum model of cell dynamics in one dimension

Multiscale modeling is emerging as one of the key challenges in mathematical biology. However, the recent rapid increase in the number of modeling methodologies being used to describe cell populations has raised a number of interesting questions. For example, at the cellular scale, how can the appropriate discrete cell-level model be identified in a given context? Additionally, how can the many phenomenological assumptions used in the derivation of models at the continuum scale be related to individual cell behavior? In order to begin to address such questions, we consider a discrete one-dimensional cell-based model in which cells are assumed to interact via linear springs. From the discrete equations of motion, the continuous Rouse [P. E. Rouse, J. Chem. Phys. 21, 1272 (1953)] model is obtained. This formalism readily allows the definition of a cell number density for which a nonlinear "fast" diffusion equation is derived. Excellent agreement is demonstrated between the continuum and discrete models. Subsequently, via the incorporation of cell division, we demonstrate that the derived nonlinear diffusion model is robust to the inclusion of more realistic biological detail. In the limit of stiff springs, where cells can be considered to be incompressible, we show that cell velocity can be directly related to cell production. This assumption is frequently made in the literature but our derivation places limits on its validity. Finally, the model is compared with a model of a similar form recently derived for a different discrete cell-based model and it is shown how the different diffusion coefficients can be understood in terms of the underlying assumptions about cell behavior in the respective discrete models.

Spatiotemporal modelling of CheY complexes in Escherichia coli chemotaxis

The chemotaxis pathway of Escherichia coli is one of the best studied and modelled biological signalling pathways. Here we extend existing modelling approaches by explicitly including a description of the formation and subcellular localization of intermediary complexes in the phosphotransfer pathway. The inclusion of these complexes shows that only about 60% of the total output response regulator (CheY) is uncomplexed at any moment and hence free to interact with its target, the flagellar motor. A clear strength of this model is its ability to predict the experimentally observable subcellular localization of CheY throughout a chemotactic response. We have found good agreement between the model output and experimentally determined CheY localization patterns.

Understanding post-operative temperature drop in cardiac surgery: a mathematical model

A mathematical model is presented to understand heat transfer processes during the cooling and re-warming of patients during cardiac surgery. Our compartmental model is able to account for many of the qualitative features observed in the cooling of various regions of the body including the central core containing the majority of organs, the rectal region containing the intestines and the outer peripheral region of skin and muscle. In particular, we focus on the issue of afterdrop: a drop in core temperature following patient re-warming, which can lead to serious post-operative complications. Model results for a typical cooling and re-warming procedure during surgery are in qualitative agreement with experimental data in producing the afterdrop effect and the observed dynamical variation in temperature between the core, rectal and peripheral regions. The influence of heat transfer processes and the volume of each compartmental region on the afterdrop effect is discussed. We find that excess fat on the peripheral and rectal regions leads to an increase in the afterdrop effect. Our model predicts that, by allowing constant re-warming after the core temperature has been raised, the afterdrop effect will be reduced.

Modelling the formation of necrotic regions in avascular tumours

The mechanisms underlying the formation of necrotic regions within avascular tumours are not well understood. In this paper, we examine the relative roles of nutrient deprivation and of cell death, from both the proliferating phase of the cell cycle via apoptosis and from the quiescent phase via necrosis, in changing the structure within multicellular tumour spheroids and particularly the accumulation of dead cell material in the centre. A mathematical model is presented and studied that accounts for nutrient diffusion, changes in cell cycling rates, the two different routes to cell death as well as active motion of cells and passive motion of the dead cell material. In studying the accumulation of dead cell matter we do not distinguish between the route by which each was formed. The resulting mathematical model is examined for a number of scenarios. Results show that in many cases the size of the necrotic core is closely correlated with low levels in nutrient concentration. However, in certain cases, particularly where the rate of necrosis is large, the resulting necrotic core can lead to regions of non-negligible nutrient concentration-dependent upon the mode of cell death.

Modelling the cell cycle and cell movement in multicellular tumour spheroids

This paper analyses a recent mathematical model of avascular tumour spheroid growth which accounts for both cell cycle dynamics and chemotactic driven cell movement. The model considers cells to exist in one of two compartments: proliferating and quiescent, as well as accounting for necrosis and apoptosis. One particular focus of this paper is the behaviour created when proliferating and quiescent cells have different chemotactic responses to an extracellular nutrient supply. Two very different steady-state behaviours are identified corresponding to those cases where proliferating cells move either more quickly or more slowly than quiescent cells in response to a gradient in the extracellular nutrient supply. The case where proliferating cells move more rapidly leads to the commonly accepted spheroid structure of a thin layer of proliferating cells surrounding an inner quiescent core. In the case where proliferating cells move more slowly than quiescent cells the model predicts an interesting structure of a thin layer of quiescent cells surrounding an inner core of proliferating and quiescent cells. The sensitivity of this tumour structure to the cell cycle model parameters is also discussed. In particular variations in the steady-state size of the tumour and the types of transient behaviour are explored. The model reveals interesting transient behaviour with sharply delineated regions of proliferating and quiescent cells.

Modelling Cell Growth and its Modulation of the G1/S Transition

We present a model for the regulation of the G(1)/S transition by cell growth in budding yeast. The model includes a description of cell size, the extracellular nutrient concentration and a simplified model of the G(1)/S transition as originally reported by Chen et al. [Mol. Biol. Cell 11:369-391, 2000]. By considering cell growth proportional to cell size we show that the cell grows exponentially. In the case where cell growth is considered proportional to the concentration of a sizer protein within the cell, our model exhibits both exponential and linear cell growth for varying parameter values. The effects of varying nutrient concentration and initial cell size are considered in the context of whether progression through the cell-size checkpoint occurs. We consider our results in relation to recent experimental evidence and discuss possible experiments for testing our theoretical predictions.

The migration of cells in multicell tumor spheroids

A mathematical model is proposed to explain the observed internalization of microspheres and 3H-thymidine labelled cells in steady-state multicellular spheroids. The model uses the conventional ideas of nutrient diffusion and consumption by the cells. In addition, a very simple model of the progress of the cells through the cell cycle is considered. Cells are divided into two classes, those proliferating (being in G1, S, G2 or M phases) and those that are quiescent (being in G0). Furthermore, the two categories are presumed to have different chemotactic responses to the nutrient gradient. The model accounts for the spatial and temporal variations in the cell categories together with mitosis, conversion between categories and cell death. Numerical solutions demonstrate that the model predicts the behavior similar to existing models but has some novel effects. It allows for spheroids to approach a steady-state size in a non-monotonic manner, it predicts self-sorting of the cell classes to produce a thin layer of rapidly proliferating cells near the outer surface and significant numbers of cells within the spheroid stalled in a proliferating state. The model predicts that overall tumor growth is not only determined by proliferation rates but also by the ability of cells to convert readily between the classes. Moreover, the steady-state structure of the spheroid indicates that if the outer layers are removed then the tumor grows quickly by recruiting cells stalled in a proliferating state. Questions are raised about the chemotactic response of cells in differing phases and to the dependency of cell cycle rates to nutrient levels.