Elementary effects for models with dimensional inputs of arbitrary type and range: Scaling and trajectory generation.

The Elementary Effects method is a global sensitivity analysis approach for identifying (un)important parameters in a model. However, it has almost exclusively been used where inputs are dimensionless and take values on [0, 1]. Here, we consider models with dimensional inputs, inputs taking values on arbitrary intervals or discrete inputs. In such cases scaling effects by a function of the input range is essential for correct ranking results. We propose two alternative dimensionless sensitivity indices by normalizing the scaled mean or median of absolute effects. Testing these indices with 9 trajectory generation methods on 4 test functions (including the Penman-Monteith equation for evapotranspiration) reveals that: i) scaled elementary effects are necessary to obtain correct parameter importance rankings; ii) small step-size methods typically produce more accurate rankings; iii) it is beneficial to compute and compare both sensitivity indices; and iv) spread and discrepancy of the simulation points are poor proxies for trajectory generation method performance.

Long-distance hormone transport via the phloem.

Several key plant hormones are synthesised in the shoot and are advected within the phloem to the root tip. In the root tip, these hormones regulate growth and developmental processes, and responses to environmental cues. However, we lack understanding of how environmental factors and biological parameters affect the delivery of hormones to the root tip. In this study, we build on existing models of phloem flow to develop a mathematical model of sugar transport alongside the transport of a generic hormone. We derive the equations for osmotically driven flow in a long, thin pipe with spatially varying membrane properties to capture the phloem loading and unloading zones. Motivated by experimental findings, we formulate solute membrane transport in terms of passive and active components, and incorporate solute unloading via bulk flow (i.e. advection with the water efflux) by including the Staverman reflection coefficient. We use the model to investigate the coupling between the sugar and hormone dynamics. The model predicts that environmental cues that lead to an increase in active sugar loading, an increase in bulk flow sugar unloading or a decrease in the relative root sugar concentration result in an increase in phloem transport velocity. Furthermore, the model reveals that such increases in phloem transport velocity result in an increase in hormone delivery to the root tip for passively loaded hormones.

Uncertainty and error in SARS-CoV-2 epidemiological parameters inferred from population-level epidemic models.

During the SARS-CoV-2 pandemic, epidemic models have been central to policy-making. Public health responses have been shaped by model-based projections and inferences, especially related to the impact of various non-pharmaceutical interventions. Accompanying this has been increased scrutiny over model performance, model assumptions, and the way that uncertainty is incorporated and presented. Here we consider a population-level model, focusing on how distributions representing host infectiousness and the infection-to-death times are modelled, and particularly on the impact of inferred epidemic characteristics if these distributions are mis-specified. We introduce an SIR-type model with the infected population structured by 'infected age', i.e. the number of days since first being infected, a formulation that enables distributions to be incorporated that are consistent with clinical data. We show that inference based on simpler models without infected age, which implicitly mis-specify these distributions, leads to substantial errors in inferred quantities relevant to policy-making, such as the reproduction number and the impact of interventions. We consider uncertainty quantification via a Bayesian approach, implementing this for both synthetic and real data focusing on UK data in the period 15 Feb-14 Jul 2020, and emphasising circumstances where it is misleading to neglect uncertainty. This manuscript was submitted as part of a theme issue on "Modelling COVID-19 and Preparedness for Future Pandemics".

In silico evidence for the utility of parsimonious root phenotypes for improved vegetative growth and carbon sequestration under drought.

Drought is a primary constraint to crop yields and climate change is expected to increase the frequency and severity of drought stress in the future. It has been hypothesized that crops can be made more resistant to drought and better able to sequester atmospheric carbon in the soil by selecting appropriate root phenotypes. We introduce OpenSimRoot_v2, an upgraded version of the functional-structural plant/soil model OpenSimRoot, and use it to test the utility of a maize root phenotype with fewer and steeper axial roots, reduced lateral root branching density, and more aerenchyma formation (i.e. the 'Steep, Cheap, and Deep' (SCD) ideotype) and different combinations of underlying SCD root phene states under rainfed and drought conditions in three distinct maize growing pedoclimatic environments in the USA, Nigeria, and Mexico. In all environments where plants are subjected to drought stress the SCD ideotype as well as several intermediate phenotypes lead to greater shoot biomass after 42 days. As an additional advantage, the amount of carbon deposited below 50 cm in the soil is twice as great for the SCD phenotype as for the reference phenotype in 5 out of 6 simulated environments. We conclude that crop growth and deep soil carbon deposition can be improved by breeding maize plants with fewer axial roots, reduced lateral root branching density, and more aerenchyma formation.

Modeling root loss reveals impacts on nutrient uptake and crop development.

Despite the widespread prevalence of root loss in plants, its effects on crop productivity are not fully understood. While root loss reduces the capacity of plants to take up water and nutrients from the soil, it may provide benefits by decreasing the resources required to maintain the root system. Here, we simulated a range of root phenotypes in different soils and root loss scenarios for barley (Hordeum vulgare), common bean (Phaseolus vulgaris), and maize (Zea mays) using and extending the open-source, functional-structural root/soil simulation model OpenSimRoot. The model enabled us to quantify the impact of root loss on shoot dry weight in these scenarios and identify in which scenarios root loss is beneficial, detrimental, or has no effect. The simulations showed that root loss is detrimental for phosphorus uptake in all tested scenarios, whereas nitrogen uptake was relatively insensitive to root loss unless main root axes were lost. Loss of axial roots reduced shoot dry weight for all phenotypes in all species and soils, whereas lateral root loss had a smaller impact. In barley and maize plants with high lateral branching density that were not phosphorus-stressed, loss of lateral roots increased shoot dry weight. The fact that shoot dry weight increased due to root loss in these scenarios indicates that plants overproduce roots for some environments, such as those found in high-input agriculture. We conclude that a better understanding of the effects of root loss on plant development is an essential part of optimizing root system phenotypes for maximizing yield.

Modelling the emergence of cities and urban patterning using coupled integro-differential equations.

Human residential population distributions show patterns of higher density clustering around local services such as shops and places of employment, displaying characteristic length scales; Fourier transforms and spatial autocorrelation show the length scale between UK cities is around 45 km. We use integro-differential equations to model the spatio-temporal dynamics of population and service density under the assumption that they benefit from spatial proximity, captured via spatial weight kernels. The system tends towards a well-mixed homogeneous state or a spatial pattern. Linear stability analysis around the homogeneous steady state predicts a modelled length-scale consistent with that observed in the data. Moreover, we show that spatial instability occurs only for perturbations with a sufficiently long wavelength and only where there is a sufficiently strong dependence of service potential on population density. Within urban centres, competition for space may cause services and population to be out of phase with one another, occupying separate parcels of land. By introducing competition, along with a preference for population to be located near, but not too near, to high service density areas, secondary out-of-phase patterns occur within the model, at a higher density and with a shorter length scale than in phase patterning. Thus, we show that a small set of core behavioural ingredients can generate aggregations of populations and services, and pattern formation within cities, with length scales consistent with real-world data. The analysis and results are valid across a wide range of parameter values and functional forms in the model.

The Virtual Root : Mathematical Modeling of Auxin Transport in the Arabidopsis Root Tip Using the Open-Source Software SimuPlant.

Hormone signals like auxin play a critical role controlling plant growth and development. Determining the mechanisms that regulate auxin distribution in cells and tissues is a vital step in understanding this hormone's role during plant development. Recent mathematical models have enabled us to understand the essential role that auxin influx and efflux carriers play in auxin transport in the Arabidopsis root tip (Band et al., Plant Cell 26(3):862-875, 2014; Grieneisen et al., Nature 449(7165):1008-1013, 2007; van den Berg et al., Development 143(18):3350-3362, 2016). In this chapter, we describe SimuPlant: The Virtual Root (SimuPlant, University of Nottingham. https://www.simuplant.org/ . Accessed 20 Sept 2019); an open source software suite, built using the OpenAlea (Pradal et al., Funct Plant Biol 35(10):751-760, 2008) framework, that is designed to simulate vertex-based models in real plant tissue geometries. We provide guidance on how to install SimuPlant, run 2D auxin transport models in the Arabidopsis root tip, manipulate parameters, and visualize model outputs.SimuPlant features a graphical user interface (GUI) designed to allow users with no programming experience to simulate auxin dynamics within the Arabidopsis root tip. Within the user interface, users of SimuPlant can select from a range of model assumptions and can choose to manipulate model and simulation parameter values. Users can then investigate how their choices affect the predicted distribution of auxin in the Arabidopsis root tip. The results of the model simulations are shown visually within the root geometry and can be exported and saved as PNG image files.

Simulating Crop Root Systems Using OpenSimRoot.

Functional-structural plant models are valuable modeling tools in analyzing plant development. A functional-structural plant model combines a three-dimensional representation of plant structure with models for physiological functions in order to better understand plant development. We present a guide to simulating crop root systems with OpenSimRoot, a feature-rich, highly cited, and open-source functional-structural root architecture model. We describe in detail how to create your own input files in conjunction with some examples. The aim of this guide is to highlight the potential of computational modeling in biology and to make modeling more accessible to the plant science community.

Differential biosynthesis and cellular permeability explain longitudinal gibberellin gradients in growing roots.

Control over cell growth by mobile regulators underlies much of eukaryotic morphogenesis. In plant roots, cell division and elongation are separated into distinct longitudinal zones and both division and elongation are influenced by the growth regulatory hormone gibberellin (GA). Previously, a multicellular mathematical model predicted a GA maximum at the border of the meristematic and elongation zones. However, GA in roots was recently measured using a genetically encoded fluorescent biosensor, nlsGPS1, and found to be low in the meristematic zone grading to a maximum at the end of the elongation zone. Furthermore, the accumulation rate of exogenous GA was also found to be higher in the elongation zone. It was still unknown which biochemical activities were responsible for these mobile small molecule gradients and whether the spatiotemporal correlation between GA levels and cell length is important for root cell division and elongation patterns. Using a mathematical modeling approach in combination with high-resolution GA measurements in vivo, we now show how differentials in several biosynthetic enzyme steps contribute to the endogenous GA gradient and how differential cellular permeability contributes to an accumulation gradient of exogenous GA. We also analyzed the effects of altered GA distribution in roots and did not find significant phenotypes resulting from increased GA levels or signaling. We did find a substantial temporal delay between complementation of GA distribution and cell division and elongation phenotypes in a GA deficient mutant. Together, our results provide models of how GA gradients are directed and in turn direct root growth.

A Mechanistic Investigation into Ischemia-Driven Distal Recurrence of Glioblastoma.

Glioblastoma (GBM) is the most aggressive primary brain tumor with a short median survival. Tumor recurrence is a clinical expectation of this disease and usually occurs along the resection cavity wall. However, previous clinical observations have suggested that in cases of ischemia following surgery, tumors are more likely to recur distally. Through the use of a previously established mechanistic model of GBM, the Proliferation Invasion Hypoxia Necrosis Angiogenesis (PIHNA) model, we explore the phenotypic drivers of this observed behavior. We have extended the PIHNA model to include a new nutrient-based vascular efficiency term that encodes the ability of local vasculature to provide nutrients to the simulated tumor. The extended model suggests sensitivity to a hypoxic microenvironment and the inherent migration and proliferation rates of the tumor cells are key factors that drive distal recurrence.

Identifying the spatial and temporal dynamics of molecularly-distinct glioblastoma sub-populations.

Glioblastomas (GBMs) are the most aggressive primary brain tumours and have no known cure. Each individual tumour comprises multiple sub-populations of genetically-distinct cells that may respond differently to targeted therapies and may contribute to disappointing clinical trial results. Image-localized biopsy techniques allow multiple biopsies to be taken during surgery and provide information that identifies regions where particular sub-populations occur within an individual GBM, thus providing insight into their regional genetic variability. These sub-populations may also interact with one another in a competitive or cooperative manner; it is important to ascertain the nature of these interactions, as they may have implications for responses to targeted therapies. We combine genetic information from biopsies with a mechanistic model of interacting GBM sub-populations to characterise the nature of interactions between two commonly occurring GBM sub-populations, those with EGFR and PDGFRA genes amplified. We study population levels found across image-localized biopsy data from a cohort of 25 patients and compare this to model outputs under competitive, cooperative and neutral interaction assumptions. We explore other factors affecting the observed simulated sub-populations, such as selection advantages and phylogenetic ordering of mutations, which may also contribute to the levels of EGFR and PDGFRA amplified populations observed in biopsy data.

Switching behaviour in vascular smooth muscle cell-matrix adhesion during oscillatory loading.

Integrins regulate mechanotransduction between smooth muscle cells (SMCs) and the extracellular matrix (ECM). SMCs resident in the walls of airways or blood vessels are continuously exposed to dynamic mechanical forces due to breathing or pulsatile blood flow. However, the resulting effects of these forces on integrin dynamics and associated cell-matrix adhesion are not well understood. Here we present experimental results from atomic force microscopy (AFM) experiments, designed to study the integrin response to external oscillatory loading of varying amplitudes applied to live aortic SMCs, together with theoretical results from a mathematical model. In the AFM experiments, a fibronectin-coated probe was used cyclically to indent and retract from the surface of the cell. We observed a transition between states of firm adhesion and of complete detachment as the amplitude of oscillatory loading increased, revealed by qualitative changes in the force timecourses. Interestingly, for some of the SMCs in the experiments, switching behaviour between the two adhesion states is observed during single timecourses at intermediate amplitudes. We obtain two qualitatively similar adhesion states in the mathematical model, where we simulate the cell, integrins and ECM as an evolving system of springs, incorporating local integrin binding dynamics. In the mathematical model, we observe a region of bistability where both the firm adhesion and detachment states can occur depending on the initial adhesion state. The differences are seen to be a result of mechanical cooperativity of integrins and cell deformation. Switching behaviour is a phenomenon associated with bistability in a stochastic system, and bistability in our deterministic mathematical model provides a potential physical explanation for the experimental results. Physiologically, bistability provides a means for transient mechanical stimuli to induce long-term changes in adhesion dynamics-and thereby the cells' ability to transmit force-and we propose further experiments for testing this hypothesis.

Speed Switch in Glioblastoma Growth Rate due to Enhanced Hypoxia-Induced Migration.

We analyze the wave speed of the Proliferation Invasion Hypoxia Necrosis Angiogenesis (PIHNA) model that was previously created and applied to simulate the growth and spread of glioblastoma (GBM), a particularly aggressive primary brain tumor. We extend the PIHNA model by allowing for different hypoxic and normoxic cell migration rates and study the impact of these differences on the wave-speed dynamics. Through this analysis, we find key variables that drive the outward growth of the simulated GBM. We find a minimum tumor wave-speed for the model; this depends on the migration and proliferation rates of the normoxic cells and is achieved under certain conditions on the migration rates of the normoxic and hypoxic cells. If the hypoxic cell migration rate is greater than the normoxic cell migration rate above a threshold, the wave speed increases above the predicted minimum. This increase in wave speed is explored through an eigenvalue and eigenvector analysis of the linearized PIHNA model, which yields an expression for this threshold. The PIHNA model suggests that an inherently faster-diffusing hypoxic cell population can drive the outward growth of a GBM as a whole, and that this effect is more prominent for faster-proliferating tumors that recover relatively slowly from a hypoxic phenotype. The findings presented here act as a first step in enabling patient-specific calibration of the PIHNA model.

Effect of Loading History on Airway Smooth Muscle Cell-Matrix Adhesions.

Integrin-mediated adhesions between airway smooth muscle (ASM) cells and the extracellular matrix (ECM) regulate how contractile forces generated within the cell are transmitted to its external environment. Environmental cues are known to influence the formation, size, and survival of cell-matrix adhesions, but it is not yet known how they are affected by dynamic fluctuations associated with tidal breathing in the intact airway. Here, we develop two closely related theoretical models to study adhesion dynamics in response to oscillatory loading of the ECM, representing the dynamic environment of ASM cells in vivo. Using a discrete stochastic-elastic model, we simulate individual integrin binding and rupture events and observe two stable regimes in which either bond formation or bond rupture dominate, depending on the amplitude of the oscillatory loading. These regimes have either a high or low fraction of persistent adhesions, which could affect the level of strain transmission between contracted ASM cells and the airway tissue. For intermediate loading, we observe a region of bistability and hysteresis due to shared loading between existing bonds; the level of adhesion depends on the loading history. These findings are replicated in a related continuum model, which we use to investigate the effect of perturbations mimicking deep inspirations (DIs). Because of the bistability, a DI applied to the high adhesion state could either induce a permanent switch to a lower adhesion state or allow a return of the system to the high adhesion state. Transitions between states are further influenced by the frequency of oscillations, cytoskeletal or ECM stiffnesses, and binding affinities, which modify the magnitudes of the stable adhesion states as well as the region of bistability. These findings could explain (in part) the transient bronchodilatory effect of a DI observed in asthmatics compared to a more sustained effect in normal subjects.

Capturing the Dynamics of a Hybrid Multiscale Cancer Model with a Continuum Model.

Cancer is a complex disease involving processes at spatial scales from subcellular, like cell signalling, to tissue scale, such as vascular network formation. A number of multiscale models have been developed to study the dynamics that emerge from the coupling between the intracellular, cellular and tissue scales. Here, we develop a continuum partial differential equation model to capture the dynamics of a particular multiscale model (a hybrid cellular automaton with discrete cells, diffusible factors and an explicit vascular network). The purpose is to test under which circumstances such a continuum model gives equivalent predictions to the original multiscale model, in the knowledge that the system details are known, and differences in model results can be explained in terms of model features (rather than unknown experimental confounding factors). The continuum model qualitatively replicates the dynamics from the multiscale model, with certain discrepancies observed owing to the differences in the modelling of certain processes. The continuum model admits travelling wave solutions for normal tissue growth and tumour invasion, with similar behaviour observed in the multiscale model. However, the continuum model enables us to analyse the spatially homogeneous steady states of the system, and hence to analyse these waves in more detail. We show that the tumour microenvironmental effects from the multiscale model mean that tumour invasion exhibits a so-called pushed wave when the carrying capacity for tumour cell proliferation is less than the total cell density at the tumour wave front. These pushed waves of tumour invasion propagate by triggering apoptosis of normal cells at the wave front. Otherwise, numerical evidence suggests that the wave speed can be predicted from linear analysis about the normal tissue steady state.

OpenSimRoot: widening the scope and application of root architectural models.

OpenSimRoot is an open-source, functional-structural plant model and mathematical description of root growth and function. We describe OpenSimRoot and its functionality to broaden the benefits of root modeling to the plant science community. OpenSimRoot is an extended version of SimRoot, established to simulate root system architecture, nutrient acquisition and plant growth. OpenSimRoot has a plugin, modular infrastructure, coupling single plant and crop stands to soil nutrient and water transport models. It estimates the value of root traits for water and nutrient acquisition in environments and plant species. The flexible OpenSimRoot design allows upscaling from root anatomy to plant community to estimate the following: resource costs of developmental and anatomical traits; trait synergisms; and (interspecies) root competition. OpenSimRoot can model three-dimensional images from magnetic resonance imaging (MRI) and X-ray computed tomography (CT) of roots in soil. New modules include: soil water-dependent water uptake and xylem flow; tiller formation; evapotranspiration; simultaneous simulation of mobile solutes; mesh refinement; and root growth plasticity. OpenSimRoot integrates plant phenotypic data with environmental metadata to support experimental designs and to gain a mechanistic understanding at system scales.

Root hydrotropism is controlled via a cortex-specific growth mechanism.

Plants can acclimate by using tropisms to link the direction of growth to environmental conditions. Hydrotropism allows roots to forage for water, a process known to depend on abscisic acid (ABA) but whose molecular and cellular basis remains unclear. Here we show that hydrotropism still occurs in roots after laser ablation removed the meristem and root cap. Additionally, targeted expression studies reveal that hydrotropism depends on the ABA signalling kinase SnRK2.2 and the hydrotropism-specific MIZ1, both acting specifically in elongation zone cortical cells. Conversely, hydrotropism, but not gravitropism, is inhibited by preventing differential cell-length increases in the cortex, but not in other cell types. We conclude that root tropic responses to gravity and water are driven by distinct tissue-based mechanisms. In addition, unlike its role in root gravitropism, the elongation zone performs a dual function during a hydrotropic response, both sensing a water potential gradient and subsequently undergoing differential growth.

Mathematical modelling of cytokines, MMPs and fibronectin fragments in osteoarthritic cartilage.

Osteoarthritis (OA) is a degenerative disease which causes pain and stiffness in joints. OA progresses through excessive degradation of joint cartilage, eventually leading to significant joint degeneration and loss of function. Cytokines, a group of cell signalling proteins, present in raised concentrations in OA joints, can be classified into pro-inflammatory and anti-inflammatory groups. They mediate cartilage degradation through several mechanisms, primarily the up-regulation of matrix metalloproteinases (MMPs), a group of collagen-degrading enzymes. In this paper we show that the interactions of cytokines within cartilage have a crucial role to play in OA progression and treatment. We develop a four-variable ordinary differential equation model for the interactions between pro- and anti-inflammatory cytokines, MMPs and fibronectin fragments (Fn-fs), a by-product of cartilage degradation and up-regulator of cytokines. We show that the model has four classes of dynamic behaviour: homoeostasis, bistable inflammation, tristable inflammation and persistent inflammation. We show that positive and negative feedbacks controlling cytokine production rates can determine either a pre-disposition to OA or initiation of OA. Further, we show that manipulation of cytokine, MMP and Fn-fs levels can be used to treat OA, but we suggest that multiple treatment targets may be essential to halt or slow disease progression.

Dynamic regulation of auxin oxidase and conjugating enzymes AtDAO1 and GH3 modulates auxin homeostasis.

The hormone auxin is a key regulator of plant growth and development, and great progress has been made understanding auxin transport and signaling. Here, we show that auxin metabolism and homeostasis are also regulated in a complex manner. The principal auxin degradation pathways in Arabidopsis include oxidation by Arabidopsis thaliana gene DIOXYGENASE FOR AUXIN OXIDATION 1/2 (AtDAO1/2) and conjugation by Gretchen Hagen3s (GH3s). Metabolic profiling of dao1-1 root tissues revealed a 50% decrease in the oxidation product 2-oxoindole-3-acetic acid (oxIAA) and increases in the conjugated forms indole-3-acetic acid aspartic acid (IAA-Asp) and indole-3-acetic acid glutamic acid (IAA-Glu) of 438- and 240-fold, respectively, whereas auxin remains close to the WT. By fitting parameter values to a mathematical model of these metabolic pathways, we show that, in addition to reduced oxidation, both auxin biosynthesis and conjugation are increased in dao1-1 Transcripts of AtDAO1 and GH3 genes increase in response to auxin over different timescales and concentration ranges. Including this regulation of AtDAO1 and GH3 in an extended model reveals that auxin oxidation is more important for auxin homoeostasis at lower hormone concentrations, whereas auxin conjugation is most significant at high auxin levels. Finally, embedding our homeostasis model in a multicellular simulation to assess the spatial effect of the dao1-1 mutant shows that auxin increases in outer root tissues in agreement with the dao1-1 mutant root hair phenotype. We conclude that auxin homeostasis is dependent on AtDAO1, acting in concert with GH3, to maintain auxin at optimal levels for plant growth and development.

The spatiotemporal order of signaling events unveils the logic of development signaling.

MOTIVATION: Animals from worms and insects to birds and mammals show distinct body plans; however, the embryonic development of diverse body plans with tissues and organs within is controlled by a surprisingly few signaling pathways. It is well recognized that combinatorial use of and dynamic interactions among signaling pathways follow specific logic to control complex and accurate developmental signaling and patterning, but it remains elusive what such logic is, or even, what it looks like. RESULTS: We have developed a computational model for Drosophila eye development with innovated methods to reveal how interactions among multiple pathways control the dynamically generated hexagonal array of R8 cells. We obtained two novel findings. First, the coupling between the long-range inductive signals produced by the proneural Hh signaling and the short-range restrictive signals produced by the antineural Notch and EGFR signaling is essential for generating accurately spaced R8s. Second, the spatiotemporal orders of key signaling events reveal a robust pattern of lateral inhibition conducted by Ato-coordinated Notch and EGFR signaling to collectively determine R8 patterning. This pattern, stipulating the orders of signaling and comparable to the protocols of communication, may help decipher the well-appreciated but poorly defined logic of developmental signaling. AVAILABILITY AND IMPLEMENTATION: The model is available upon request. CONTACT: hao.zhu@ymail.com SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.