MEIGO: an open-source software suite based on metaheuristics for global optimization in systems biology and bioinformatics

A unified framework for prediction of liver steatosis dynamics in response to different diet and drug interventions

BACKGROUND & AIMS

Non-alcoholic fatty liver disease (NAFLD) is a common metabolic disorder, characterized by the accumulation of excess fat in the liver, and is a driving factor for various severe liver diseases. These multi-factorial and multi-timescale changes are observed in different clinical studies, but these studies have not been integrated into a unified framework. In this study, we aim to present such a unified framework in the form of a dynamic mathematical model.

METHODS

For model training and validation, we collected data for dietary or drug-induced interventions aimed at reducing or increasing liver fat. The model was formulated using ordinary differential equations (ODEs) and the mathematical analysis, model simulation, model formulation and the model parameter estimation were all performed in MATLAB.

RESULTS

Our mathematical model describes accumulation of fat in the liver and predicts changes in lipid fluxes induced by both dietary and drug interventions. The model is validated using data from a wide range of drug and dietary intervention studies and can predict both short-term (days) and long-term (weeks) changes in liver fat. Importantly, the model computes the contribution of each individual lipid flux to the total liver fat dynamics. Furthermore, the model can be combined with an established bodyweight model, to simulate even longer scenarios (years), also including the effects of insulin resistance and body weight. To help prepare for corresponding eHealth applications, we also present a way to visualize the simulated changes, using dynamically changing lipid droplets, seen in images of liver biopsies.

CONCLUSION

In conclusion, we believe that the minimal model presented herein might be a useful tool for future applications, and to further integrate and understand data regarding changes in dietary and drug induced changes in ectopic TAG in the liver. With further development and validation, the minimal model could be used as a disease progression model for steatosis.

A physiologically-based digital twin for alcohol consumption-predicting real-life drinking responses and long-term plasma PEth

Alcohol consumption is associated with a wide variety of preventable health complications and is a major risk factor for all-cause mortality in the age group 15-47 years. To reduce dangerous drinking behavior, eHealth applications have shown promise. A particularly interesting potential lies in the combination of eHealth apps with mathematical models. However, existing mathematical models do not consider real-life situations, such as combined intake of meals and beverages, and do not connect drinking to clinical markers, such as phosphatidylethanol (PEth). Herein, we present such a model which can simulate real-life situations and connect drinking to long-term markers. The new model can accurately describe both estimation data according to a χ2 -test (187.0 < Tχ2 = 226.4) and independent validation data (70.8 < Tχ2 = 93.5). The model can also be personalized using anthropometric data from a specific individual and can thus be used as a physiologically-based digital twin. This twin is also able to connect short-term consumption of alcohol to the long-term dynamics of PEth levels in the blood, a clinical biomarker of alcohol consumption. Here we illustrate how connecting short-term consumption to long-term markers allows for a new way to determine patient alcohol consumption from measured PEth levels. An additional use case of the twin could include the combined evaluation of patient-reported AUDIT forms and measured PEth levels. Finally, we integrated the new model into an eHealth application, which could help guide individual users or clinicians to help reduce dangerous drinking.

Deciphering Metabolic Pathways in High-Seeding-Density Fed-Batch Processes for Monoclonal Antibody Production: A Computational Modeling Perspective

Due to their high specificity, monoclonal antibodies (mAbs) have garnered significant attention in recent decades, with advancements in production processes, such as high-seeding-density (HSD) strategies, contributing to improved titers. This study provides a thorough investigation of high seeding processes for mAb production in Chinese hamster ovary (CHO) cells, focused on identifying significant metabolites and their interactions. We observed high glycolytic fluxes, the depletion of asparagine, and a shift from lactate production to consumption. Using a metabolic network and flux analysis, we compared the standard fed-batch (STD FB) with HSD cultivations, exploring supplementary lactate and cysteine, and a bolus medium enriched with amino acids. We reconstructed a metabolic network and kinetic models based on the observations and explored the effects of different feeding strategies on CHO cell metabolism. Our findings revealed that the addition of a bolus medium (BM) containing asparagine improved final titers. However, increasing the asparagine concentration in the feed further prevented the lactate shift, indicating a need to find a balance between increased asparagine to counteract limitations and lower asparagine to preserve the shift in lactate metabolism.

Optimal experimental design (OED) for the growth rate of microbial populations. Are they really more "optimal" than uniform designs?

Secondary growth models from predictive microbiology can describe how the growth rate of microbial populations varies with environmental conditions. Because these models are built based on time and resource consuming experiments, model-based Optimal Experimental Design (OED) can be of interest to reduce the experimental load. In this study, we identify optimal experimental designs for two common models (full Ratkowsky and Cardinal Parameters Model (CPM)) for a different number of experiments (10-30). Calculations are also done fixing one or more model parameters, observing that this decision strongly affects the layout of the OED. Using in silico experiments, we conclude that OEDs are more informative than conventional (equidistant) designs with the same number of experiments. However, OEDs cluster the experiments near the growth limits (Xmin and Xmax) resulting in impractical designs with aggregated experimental runs ~10 times longer than conventional designs. To mitigate this, we propose a novel optimality criterion (i.e., the objective function) that accounts for the aggregated time. The novel criterion provides a reduction in parameter uncertainty with respect to the conventional design, without an increase in the experimental load. These results underline that an OED is only based on information theory (Fisher information), so the results can be impractical when actual experimental limitations are considered. The study also emphasizes that most OED schemes identify where to measure, but do not give an indication on how many experiments should be made. In this sense, numerical simulations can estimate the parameter uncertainty that would be obtained for a particular experimental design (OED or not). These results and methodologies (available in Open Code) can guide the design of future experiments for the development of secondary growth models.

A multi-scale digital twin for adiposity-driven insulin resistance in humans: diet and drug effects

BACKGROUND

The increased prevalence of insulin resistance is one of the major health risks in society today. Insulin resistance involves both short-term dynamics, such as altered meal responses, and long-term dynamics, such as the development of type 2 diabetes. Insulin resistance also occurs on different physiological levels, ranging from disease phenotypes to organ-organ communication and intracellular signaling. To better understand the progression of insulin resistance, an analysis method is needed that can combine different timescales and physiological levels. One such method is digital twins, consisting of combined mechanistic mathematical models. We have previously developed a model for short-term glucose homeostasis and intracellular insulin signaling, and there exist long-term weight regulation models. Herein, we combine these models into a first interconnected digital twin for the progression of insulin resistance in humans.

METHODS

The model is based on ordinary differential equations representing biochemical and physiological processes, in which unknown parameters were fitted to data using a MATLAB toolbox.

RESULTS

The interconnected twin correctly predicts independent data from a weight increase study, both for weight-changes, fasting plasma insulin and glucose levels, and intracellular insulin signaling. Similarly, the model can predict independent weight-change data in a weight loss study with the weight loss drug topiramate. The model can also predict non-measured variables.

CONCLUSIONS

The model presented herein constitutes the basis for a new digital twin technology, which in the future could be used to aid medical pedagogy and increase motivation and compliance and thus aid in the prevention and treatment of insulin resistance.

Machine Learning in Environmental Research: Common Pitfalls and Best Practices

Machine learning (ML) is increasingly used in environmental research to process large data sets and decipher complex relationships between system variables. However, due to the lack of familiarity and methodological rigor, inadequate ML studies may lead to spurious conclusions. In this study, we synthesized literature analysis with our own experience and provided a tutorial-like compilation of common pitfalls along with best practice guidelines for environmental ML research. We identified more than 30 key items and provided evidence-based data analysis based on 148 highly cited research articles to exhibit the misconceptions of terminologies, proper sample size and feature size, data enrichment and feature selection, randomness assessment, data leakage management, data splitting, method selection and comparison, model optimization and evaluation, and model explainability and causality. By analyzing good examples on supervised learning and reference modeling paradigms, we hope to help researchers adopt more rigorous data preprocessing and model development standards for more accurate, robust, and practicable model uses in environmental research and applications.

IDESS: a toolbox for identification and automated design of stochastic gene circuits

MOTIVATION

One of the main causes hampering predictability during the model identification and automated design of gene circuits in synthetic biology is the effect of molecular noise. Stochasticity may significantly impact the dynamics and function of gene circuits, specially in bacteria and yeast due to low mRNA copy numbers. Standard stochastic simulation methods are too computationally costly in realistic scenarios to be applied to optimization-based design or parameter estimation.

RESULTS

In this work, we present IDESS (Identification and automated Design of Stochastic gene circuitS), a software toolbox for optimization-based design and model identification of gene regulatory circuits in the stochastic regime. This software incorporates an efficient approximation of the Chemical Master Equation as well as a stochastic simulation algorithm-both with GPU and CPU implementations-combined with global optimization algorithms capable of solving Mixed Integer Nonlinear Programming problems. The toolbox efficiently addresses two types of problems relevant in systems and synthetic biology: the automated design of stochastic synthetic gene circuits, and the parameter estimation for model identification of stochastic gene regulatory networks.

AVAILABILITY AND IMPLEMENTATION

IDESS runs under the MATLAB environment and it is available under GPLv3 license at https://doi.org/10.5281/zenodo.7788692.

Automated Design of Synthetic Gene Circuits in the Presence of Molecular Noise

Microorganisms (mainly bacteria and yeast) are frequently used as hosts for genetic constructs in synthetic biology applications. Molecular noise might have a significant effect on the dynamics of gene regulation in microbial cells, mainly attributed to the low copy numbers of mRNA species involved. However, the inclusion of molecular noise in the automated design of biocircuits is not a common practice due to the computational burden linked to the chemical master equation describing the dynamics of stochastic gene regulatory circuits. Here, we address the automated design of synthetic gene circuits under the effect of molecular noise combining a mixed integer nonlinear global optimization method with a partial integro-differential equation model describing the evolution of stochastic gene regulatory systems that approximates very efficiently the chemical master equation. We demonstrate the performance of the proposed methodology through a number of examples of relevance in synthetic biology, including different bimodal stochastic gene switches, robust stochastic oscillators, and circuits capable of achieving biochemical adaptation under noise.

Elucidating Fibroblast Growth Factor-Induced Kinome Dynamics Using Targeted Mass Spectrometry and Dynamic Modeling

Fibroblast growth factors (FGFs) are paracrine or endocrine signaling proteins that, activated by their ligands, elicit a wide range of health and disease-related processes, such as cell proliferation and the epithelial-to-mesenchymal transition. The detailed molecular pathway dynamics that coordinate these responses have remained to be determined. To elucidate these, we stimulated MCF-7 breast cancer cells with either FGF2, FGF3, FGF4, FGF10, or FGF19. Following activation of the receptor, we quantified the kinase activity dynamics of 44 kinases using a targeted mass spectrometry assay. Our system-wide kinase activity data, supplemented with (phospho)proteomics data, reveal ligand-dependent distinct pathway dynamics, elucidate the involvement of not earlier reported kinases such as MARK, and revise some of the pathway effects on biological outcomes. In addition, logic-based dynamic modeling of the kinome dynamics further verifies the biological goodness-of-fit of the predicted models and reveals BRAF-driven activation upon FGF2 treatment and ARAF-driven activation upon FGF4 treatment.

Optimal evaluation of energy yield and driving force in microbial metabolic pathway variants

This work presents a methodology to evaluate the bioenergetic feasibility of alternative metabolic pathways for a given microbial conversion, optimising their energy yield and driving forces as a function of the concentration of metabolic intermediates. The tool, based on thermodynamic principles and multi-objective optimisation, accounts for pathway variants in terms of different electron carriers, as well as energy conservation (proton translocating) reactions within the pathway. The method also accommodates other constraints, some of them non-linear, such as the balance of conserved moieties. The approach involves the transformation of the maximum energy yield problem into a multi-objective mixed-integer linear optimisation problem which is then subsequently solved using the epsilon-constraint method, highlighting the trade-off between yield and rate in metabolic reactions. The methodology is applied to analyse several pathway alternatives occurring during propionate oxidation in anaerobic fermentation processes, as well as to the reverse TCA cycle pathway occurring during autotrophic microbial CO2 fixation. The results obtained using the developed methodology match previously reported literature and bring about insights into the studied pathways.

A comprehensive mechanistic model of adipocyte signaling with layers of confidence

Adipocyte signaling, normally and in type 2 diabetes, is far from fully understood. We have earlier developed detailed dynamic mathematical models for several well-studied, partially overlapping, signaling pathways in adipocytes. Still, these models only cover a fraction of the total cellular response. For a broader coverage of the response, large-scale phosphoproteomic data and systems level knowledge on protein interactions are key. However, methods to combine detailed dynamic models with large-scale data, using information about the confidence of included interactions, are lacking. We have developed a method to first establish a core model by connecting existing models of adipocyte cellular signaling for: (1) lipolysis and fatty acid release, (2) glucose uptake, and (3) the release of adiponectin. Next, we use publicly available phosphoproteome data for the insulin response in adipocytes together with prior knowledge on protein interactions, to identify phosphosites downstream of the core model. In a parallel pairwise approach with low computation time, we test whether identified phosphosites can be added to the model. We iteratively collect accepted additions into layers and continue the search for phosphosites downstream of these added layers. For the first 30 layers with the highest confidence (311 added phosphosites), the model predicts independent data well (70-90% correct), and the predictive capability gradually decreases when we add layers of decreasing confidence. In total, 57 layers (3059 phosphosites) can be added to the model with predictive ability kept. Finally, our large-scale, layered model enables dynamic simulations of systems-wide alterations in adipocytes in type 2 diabetes.

Robustness and innovation in synthetic genotype networks

Genotype networks are sets of genotypes connected by small mutational changes that share the same phenotype. They facilitate evolutionary innovation by enabling the exploration of different neighborhoods in genotype space. Genotype networks, first suggested by theoretical models, have been empirically confirmed for proteins and RNAs. Comparative studies also support their existence for gene regulatory networks (GRNs), but direct experimental evidence is lacking. Here, we report the construction of three interconnected genotype networks of synthetic GRNs producing three distinct phenotypes in Escherichia coli. Our synthetic GRNs contain three nodes regulating each other by CRISPR interference and governing the expression of fluorescent reporters. The genotype networks, composed of over twenty different synthetic GRNs, provide robustness in face of mutations while enabling transitions to innovative phenotypes. Through realistic mathematical modeling, we quantify robustness and evolvability for the complete genotype-phenotype map and link these features mechanistically to GRN motifs. Our work thereby exemplifies how GRN evolution along genotype networks might be driving evolutionary innovation.

A quantitative model for human neurovascular coupling with translated mechanisms from animals

Neurons regulate the activity of blood vessels through the neurovascular coupling (NVC). A detailed understanding of the NVC is critical for understanding data from functional imaging techniques of the brain. Many aspects of the NVC have been studied both experimentally and using mathematical models; various combinations of blood volume and flow, local field potential (LFP), hemoglobin level, blood oxygenation level-dependent response (BOLD), and optogenetics have been measured and modeled in rodents, primates, or humans. However, these data have not been brought together into a unified quantitative model. We now present a mathematical model that describes all such data types and that preserves mechanistic behaviors between experiments. For instance, from modeling of optogenetics and microscopy data in mice, we learn cell-specific contributions; the first rapid dilation in the vascular response is caused by NO-interneurons, the main part of the dilation during longer stimuli is caused by pyramidal neurons, and the post-peak undershoot is caused by NPY-interneurons. These insights are translated and preserved in all subsequent analyses, together with other insights regarding hemoglobin dynamics and the LFP/BOLD-interplay, obtained from other experiments on rodents and primates. The model can predict independent validation-data not used for training. By bringing together data with complementary information from different species, we both understand each dataset better, and have a basis for a new type of integrative analysis of human data.

Understanding the dosing-time-dependent antihypertensive effect of valsartan and aspirin through mathematical modeling

Chronopharmacology of arterial hypertension impacts the long-term cardiovascular risk of hypertensive subjects. Therefore, clinical and computational studies have proposed optimizing antihypertensive medications' dosing time (Ta). However, the causes and mechanisms underlying the Ta-dependency antihypertensive effect have not been elucidated. Here we propose using a Ta- dependent effect model to understand and predict the antihypertensive effect of valsartan and aspirin throughout the day in subjects with grade I or II essential hypertension. The model based on physiological regulation mechanisms includes a periodic function for each parameter that changes significantly after treatment. Circadian variations of parameters depending on the dosing time allowed the determination of regulation mechanisms dependent on the circadian rhythm that were most relevant for the action of each drug. In the case of valsartan, it is the regulation of vasodilation and systemic vascular resistance. In the case of aspirin, the antithrombotic effect generates changes in the sensitivity of systemic vascular resistance and heart rate to changes in physical activity. Dosing time-dependent models predict a more significant effect on systemic vascular resistance and blood pressure when administering valsartan or aspirin at bedtime. However, circadian dependence on the regulation mechanisms showed different sensitivity of their circadian parameters and shapes of functions, presenting different phase shifts and amplitude. Therefore, different mechanisms of action and pharmacokinetic properties of each drug can generate different profiles of Ta-dependence of antihypertensive effect and optimal dosing times.

Mechanistic model for human brain metabolism and its connection to the neurovascular coupling

The neurovascular and neurometabolic couplings (NVC and NMC) connect cerebral activity, blood flow, and metabolism. This interconnection is used in for instance functional imaging, which analyses the blood-oxygen-dependent (BOLD) signal. The mechanisms underlying the NVC are complex, which warrants a model-based analysis of data. We have previously developed a mechanistically detailed model for the NVC, and others have proposed detailed models for cerebral metabolism. However, existing metabolic models are still not fully utilizing available magnetic resonance spectroscopy (MRS) data and are not connected to detailed models for NVC. Therefore, we herein present a new model that integrates mechanistic modelling of both MRS and BOLD data. The metabolic model covers central metabolism, using a minimal set of interactions, and can describe time-series data for glucose, lactate, aspartate, and glutamate, measured after visual stimuli. Statistical tests confirm that the model can describe both estimation data and predict independent validation data, not used for model training. The interconnected NVC model can simultaneously describe BOLD data and can be used to predict expected metabolic responses in experiments where metabolism has not been measured. This model is a step towards a useful and mechanistically detailed model for cerebral blood flow and metabolism, with potential applications in both basic research and clinical applications.

Improving dynamic predictions with ensembles of observable models

MOTIVATION

Dynamic mechanistic modelling in systems biology has been hampered by the complexity and variability associated with the underlying interactions, and by uncertain and sparse experimental measurements. Ensemble modelling, a concept initially developed in statistical mechanics, has been introduced in biological applications with the aim of mitigating those issues. Ensemble modelling uses a collection of different models compatible with the observed data to describe the phenomena of interest. However, since systems biology models often suffer from a lack of identifiability and observability, ensembles of models are particularly unreliable when predicting non-observable states.

RESULTS

We present a strategy to assess and improve the reliability of a class of model ensembles. In particular, we consider kinetic models described using ordinary differential equations with a fixed structure. Our approach builds an ensemble with a selection of the parameter vectors found when performing parameter estimation with a global optimization metaheuristic. This technique enforces diversity during the sampling of parameter space and it can quantify the uncertainty in the predictions of state trajectories. We couple this strategy with structural identifiability and observability analysis, and when these tests detect possible prediction issues we obtain model reparameterizations that surmount them. The end result is an ensemble of models with the ability to predict the internal dynamics of a biological process. We demonstrate our approach with models of glucose regulation, cell division, circadian oscillations and the JAK-STAT signalling pathway.

AVAILABILITY AND IMPLEMENTATION

The code that implements the methodology and reproduces the results is available at https://doi.org/10.5281/zenodo.6782638.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

MLAGO: machine learning-aided global optimization for Michaelis constant estimation of kinetic modeling

BACKGROUND

Kinetic modeling is a powerful tool for understanding the dynamic behavior of biochemical systems. For kinetic modeling, determination of a number of kinetic parameters, such as the Michaelis constant (K_m), is necessary, and global optimization algorithms have long been used for parameter estimation. However, the conventional global optimization approach has three problems: (i) It is computationally demanding. (ii) It often yields unrealistic parameter values because it simply seeks a better model fitting to experimentally observed behaviors. (iii) It has difficulty in identifying a unique solution because multiple parameter sets can allow a kinetic model to fit experimental data equally well (the non-identifiability problem).

RESULTS

To solve these problems, we propose the Machine Learning-Aided Global Optimization (MLAGO) method for K_m estimation of kinetic modeling. First, we use a machine learning-based K_m predictor based only on three factors: EC number, KEGG Compound ID, and Organism ID, then conduct a constrained global optimization-based parameter estimation by using the machine learning-predicted K_m values as the reference values. The machine learning model achieved relatively good prediction scores: RMSE = 0.795 and R² = 0.536, making the subsequent global optimization easy and practical. The MLAGO approach reduced the error between simulation and experimental data while keeping K_m values close to the machine learning-predicted values. As a result, the MLAGO approach successfully estimated K_m values with less computational cost than the conventional method. Moreover, the MLAGO approach uniquely estimated K_m values, which were close to the measured values.

CONCLUSIONS

MLAGO overcomes the major problems in parameter estimation, accelerates kinetic modeling, and thus ultimately leads to better understanding of complex cellular systems. The web application for our machine learning-based K_m predictor is accessible at https://sites.google.com/view/kazuhiro-maeda/software-tools-web-apps , which helps modelers perform MLAGO on their own parameter estimation tasks.

Can the Kuznetsov Model Replicate and Predict Cancer Growth in Humans?

Several mathematical models to predict tumor growth over time have been developed in the last decades. A central aspect of such models is the interaction of tumor cells with immune effector cells. The Kuznetsov model (Kuznetsov et al. in Bull Math Biol 56(2):295–321, 1994) is the most prominent of these models and has been used as a basis for many other related models and theoretical studies. However, none of these models have been validated with large-scale real-world data of human patients treated with cancer immunotherapy. In addition, parameter estimation of these models remains a major bottleneck on the way to model-based and data-driven medical treatment. In this study, we quantitatively fit Kuznetsov’s model to a large dataset of 1472 patients, of which 210 patients have more than six data points, by estimating the model parameters of each patient individually. We also conduct a global practical identifiability analysis for the estimated parameters. We thus demonstrate that several combinations of parameter values could lead to accurate data fitting. This opens the potential for global parameter estimation of the model, in which the values of all or some parameters are fixed for all patients. Furthermore, by omitting the last two or three data points, we show that the model can be extrapolated and predict future tumor dynamics. This paves the way for a more clinically relevant application of mathematical tumor modeling, in which the treatment strategy could be adjusted in advance according to the model’s future predictions.

Data-based stochastic modeling reveals sources of activity bursts in single-cell TGF-β signaling

Cells sense their surrounding by employing intracellular signaling pathways that transmit hormonal signals from the cell membrane to the nucleus. TGF-β/SMAD signaling encodes various cell fates, controls tissue homeostasis and is deregulated in diseases such as cancer. The pathway shows strong heterogeneity at the single-cell level, but quantitative insights into mechanisms underlying fluctuations at various time scales are still missing, partly due to inefficiency in the calibration of stochastic models that mechanistically describe signaling processes. In this work we analyze single-cell TGF-β/SMAD signaling and show that it exhibits temporal stochastic bursts which are dose-dependent and whose number and magnitude correlate with cell migration. We propose a stochastic modeling approach to mechanistically describe these pathway fluctuations with high computational efficiency. Employing high-order numerical integration and fitting to burst statistics we enable efficient quantitative parameter estimation and discriminate models that assume noise in different reactions at the receptor level. This modeling approach suggests that stochasticity in the internalization of TGF-β receptors into endosomes plays a key role in the observed temporal bursting. Further, the model predicts the single-cell dynamics of TGF-β/SMAD signaling in untested conditions, e.g., successfully reflects memory effects of signaling noise and cellular sensitivity towards repeated stimulation. Taken together, our computational framework based on burst analysis, noise modeling and path computation scheme is a suitable tool for the data-based modeling of complex signaling pathways, capable of identifying the source of temporal noise.

Fluctuations in molecular networks give rise to heterogeneity in cellular behavior and therefore promote the diversification of tissues. For a better understanding of cellular decision making, it is important to identify sources of molecular fluctuations and to quantitatively describe them by predictive mathematical models. In this work, we focused on temporal fluctuations of the TGF-β signaling pathway that is important for controlling cell division and migration.

We characterized a single-cell dataset comprising hundreds of cells using time series analysis and a large-scale stochastic model. By fitting several model variants to the data, we identified the stochastic internalization of cell surface receptors into endosomes as a main source of temporal fluctuations (’bursts’) in the signaling pathway. The corresponding model accurately predicted novel experimental data, and provided insights into the long-term memory of signaling fluctuations. In summary, we propose a modeling approach to quantitatively describe heterogeneous behavior in large-scale single-cell datasets and to identify the underlying biological mechanisms.

A systems biology analysis of lipolysis and fatty acid release from adipocytes in vitro and from adipose tissue in vivo

Lipolysis and the release of fatty acids to supply energy fuel to other organs, such as between meals, during exercise, and starvation, are fundamental functions of the adipose tissue. The intracellular lipolytic pathway in adipocytes is activated by adrenaline and noradrenaline, and inhibited by insulin. Circulating fatty acids are elevated in type 2 diabetic individuals. The mechanisms behind this elevation are not fully known, and to increase the knowledge a link between the systemic circulation and intracellular lipolysis is key. However, data on lipolysis and knowledge from in vitro systems have not been linked to corresponding in vivo data and knowledge in vivo. Here, we use mathematical modelling to provide such a link. We examine mechanisms of insulin action by combining in vivo and in vitro data into an integrated mathematical model that can explain all data. Furthermore, the model can describe independent data not used for training the model. We show the usefulness of the model by simulating new and more challenging experimental setups in silico, e.g. the extracellular concentration of fatty acids during an insulin clamp, and the difference in such simulations between individuals with and without type 2 diabetes. Our work provides a new platform for model-based analysis of adipose tissue lipolysis, under both non-diabetic and type 2 diabetic conditions.

RBS and Promoter Strengths Determine the Cell-Growth-Dependent Protein Mass Fractions and Their Optimal Synthesis Rates

Models of gene expression considering host–circuit interactions are relevant for understanding both the strategies and associated trade-offs that cell endogenous genes have evolved and for the efficient design of heterologous protein expression systems and synthetic genetic circuits. Here, we consider a small-size model of gene expression dynamics in bacterial cells accounting for host–circuit interactions due to limited cellular resources. We define the cellular resources recruitment strength as a key functional coefficient that explains the distribution of resources among the host and the genes of interest and the relationship between the usage of resources and cell growth. This functional coefficient explicitly takes into account lab-accessible gene expression characteristics, such as promoter and ribosome binding site (RBS) strengths, capturing their interplay with the growth-dependent flux of available free cell resources. Despite its simplicity, the model captures the differential role of promoter and RBS strengths in the distribution of protein mass fractions as a function of growth rate and the optimal protein synthesis rate with remarkable fit to the experimental data from the literature for Escherichia coli. This allows us to explain why endogenous genes have evolved different strategies in the expression space and also makes the model suitable for model-based design of exogenous synthetic gene expression systems with desired characteristics.

Stochastic and Age-Dependent Proteostasis Decline Underlies Heterogeneity in Heat-Shock Response Dynamics

Significant non-genetic stochastic factors affect aging, causing lifespan differences among individuals, even those sharing the same genetic and environmental background. In Caenorhabditis elegans, differences in heat-shock response (HSR) are predictive of lifespan. However, factors contributing to the heterogeneity of HSR are still not fully elucidated. Here, the authors characterized HSR dynamics in isogenic C. elegans expressing GFP reporter for hsp-16.2 for identifying the key contributors of HSR heterogeneity. Specifically, microfluidic devices that enable cross-sectional and longitudinal measurements of HSR dynamics in C. elegans at different scales are developed: in populations, within individuals, and in embryos. The authors adapted a mathematical model of HSR to single C. elegans and identified model parameters associated with proteostasis-maintenance of protein homeostasis-more specifically, protein turnover, as the major drivers of heterogeneity in HSR dynamics. It is verified that individuals with enhanced proteostasis fidelity in early adulthood live longer. The model-based comparative analysis of protein turnover in day-1 and day-2 adult C. elegans revealed a stochastic-onset of age-related proteostasis decline that increases the heterogeneity of HSR capacity. Finally, the analysis of C. elegans embryos showed higher HSR and proteostasis capacity than young adults and established transgenerational contribution to HSR heterogeneity that depends on maternal age.

In-depth transcriptomic analysis of human retina reveals molecular mechanisms underlying diabetic retinopathy

Diabetic Retinopathy (DR) is among the major global causes for vision loss. With the rise in diabetes prevalence, an increase in DR incidence is expected. Current understanding of both the molecular etiology and pathways involved in the initiation and progression of DR is limited. Via RNA-Sequencing, we analyzed mRNA and miRNA expression profiles of 80 human post-mortem retinal samples from 43 patients diagnosed with various stages of DR. We found differentially expressed transcripts to be predominantly associated with late stage DR and pathways such as hippo and gap junction signaling. A multivariate regression model identified transcripts with progressive changes throughout disease stages, which in turn displayed significant overlap with sphingolipid and cGMP-PKG signaling. Combined analysis of miRNA and mRNA expression further uncovered disease-relevant miRNA/mRNA associations as potential mechanisms of post-transcriptional regulation. Finally, integrating human retinal single cell RNA-Sequencing data revealed a continuous loss of retinal ganglion cells, and Müller cell mediated changes in histidine and β-alanine signaling. While previously considered primarily a vascular disease, attention in DR has shifted to additional mechanisms and cell-types. Our findings offer an unprecedented and unbiased insight into molecular pathways and cell-specific changes in the development of DR, and provide potential avenues for future therapeutic intervention.

Reducing the uncertainty on chilling requirements for endodormancy breaking of temperate fruits by data-based parameter estimation of the dynamic model: a test case in apricot

The Dynamic model has been described as one of the most accurate models to quantify chill accumulation based on hourly temperatures in nuts and temperate fruits. This model considers that a dynamic process occurs at a biochemical level that determines the endodormancy breaking through the accumulation of the so-called portions. The kinetic parameters present in the model should reflect how the fruit trees integrate chilling exposure and thus they should be characteristic for each species. However, the original parameter values, reported in the late 1980s, are still being used. Even if the use of such parameter values is useful to compare among chilling requirements (CRs) for different species or cultivars, it is not the optimal choice when one intends to explain the CR variations in different years for a given cultivar. In this work we propose a data-based model calibration that makes use of phenological data for different apricot cultivars within different years to obtain model parameters, which minimize the variations among years and that have, at the same time, physical meaning to characterize the incumbent species. Results reveal that the estimation not only reduces the accumulated portion dispersion within the considered time periods but also allows to improve the CR predictions for subsequent years. We propose a set of model parameter values to predict endodormancy breaking dates in the apricot cultivars studied here.

Model-based robustness and bistability analysis for methylation-based, epigenetic memory systems

In recent years, epigenetic memory systems have been developed based on DNA methylation and positive feedback systems. Achieving a robust design for these systems is generally a challenging and multifactorial task. We developed and validated a novel mathematical model to describe methylation-based epigenetic memory systems that capture switching dynamics of methylation levels and methyltransferase amounts induced by different inputs. A bifurcation analysis shows that the system operates in the bistable range, but in its current setup is not robust to changes in parameters. An expansion of the model captures heterogeneity of cell populations by accounting for distributed cell division rates. Simulations predict that the system is highly sensitive to variations in temperature, which affects cell division and the efficiency of the zinc finger repressor. A moderate decrease in temperature leads to a highly heterogeneous response to input signals and bistability on a single-cell level. The predictions of our model were confirmed by flow cytometry experiments conducted in this study. Overall, the results of our study give insights into the functional mechanisms of methylation-based memory systems and demonstrate that the switching dynamics can be highly sensitive to experimental conditions.

Circadian Rhythm of Blood Pressure of Dipper and Non-dipper Patients With Essential Hypertension: A Mathematical Modeling Approach

Blood pressure in humans presents a circadian variation profile with a morning increase, a small postprandial valley, and a deeper descent during night-time rest. Under certain conditions, the nocturnal decline in blood pressure can be reduced or even reversed (non-dipper), which is related to a significantly worse prognosis than a normal fall pattern (dipper). Despite several advances in recent years, our understanding of blood pressure's temporal structure, its sources and mechanisms is far from complete. In this work, we developed an ordinary differential equation-based mathematical model capable of capturing the circadian rhythm of blood pressure in dipper and non-dipper patients with arterial hypertension. The model was calibrated by means of global optimization, using 24-h data of systolic and diastolic blood pressure, physical activity, heart rate, blood glucose and norepinephrine, obtained from the literature. After fitting the model, the mean of the normalized error for each data point was <0.2%, and confidence intervals indicate that all parameters were identifiable. Sensitivity analysis allowed identifying the most relevant parameters and therefore inferring the most important blood pressure regulatory mechanisms involved in the non-dipper status, namely, increase in sympathetic over parasympathetic nervous tone, lower influence of physical activity on heart rate and greater influence of physical activity and glucose on the systemic vascular resistance. In summary, this model allows explaining the circadian rhythm of blood pressure and deepening the understanding of the underlying mechanisms and interactions integrating the results of previous works.

PEtab-Interoperable specification of parameter estimation problems in systems biology

Reproducibility and reusability of the results of data-based modeling studies are essential. Yet, there has been-so far-no broadly supported format for the specification of parameter estimation problems in systems biology. Here, we introduce PEtab, a format which facilitates the specification of parameter estimation problems using Systems Biology Markup Language (SBML) models and a set of tab-separated value files describing the observation model and experimental data as well as parameters to be estimated. We already implemented PEtab support into eight well-established model simulation and parameter estimation toolboxes with hundreds of users in total. We provide a Python library for validation and modification of a PEtab problem and currently 20 example parameter estimation problems based on recent studies.

Dynamic co-culture metabolic models reveal the fermentation dynamics, metabolic capacities and interplays of cheese starter cultures

In this study, we have investigated the cheese starter culture as a microbial community through a question: can the metabolic behaviour of a co-culture be explained by the characterized individual organism that constituted the co-culture? To address this question, the dairy-origin lactic acid bacteria Lactococcus lactis subsp. cremoris, Lactococcus lactis subsp. lactis, Streptococcus thermophilus and Leuconostoc mesenteroides, commonly used in cheese starter cultures, were grown in pure and four different co-cultures. We used a dynamic metabolic modelling approach based on the integration of the genome-scale metabolic networks of the involved organisms to simulate the co-cultures. The strain-specific kinetic parameters of dynamic models were estimated using the pure culture experiments and they were subsequently applied to co-culture models. Biomass, carbon source, lactic acid and most of the amino acid concentration profiles simulated by the co-culture models fit closely to the experimental results and the co-culture models explained the mechanisms behind the dynamic microbial abundance. We then applied the co-culture models to estimate further information on the co-cultures that could not be obtained by the experimental method used. This includes estimation of the profile of various metabolites in the co-culture medium such as flavour compounds produced and the individual organism level metabolic exchange flux profiles, which revealed the potential metabolic interactions between organisms in the co-cultures.

BioSwitch: a tool for the detection of bistability and multi-steady state behaviour in signalling and gene regulatory networks

MOTIVATION

Multi-steady state behaviour, and in particular multi-stability, provides biological systems with the capacity to take reliable decisions (such as cell fate determination). A problem arising frequently in systems biology is to elucidate whether a signal transduction mechanism or a gene regulatory network has the capacity for multi-steady state behaviour, and consequently for a switch-like response to stimuli. Bifurcation diagrams are a powerful instrument in non-linear analysis to study the qualitative and quantitative behaviour of equilibria including bifurcation into different equilibrium branches and bistability. However, in the context of signalling pathways, the inherent large parametric uncertainty hampers the (direct) use of standard bifurcation tools.

RESULTS

We present BioSwitch, a toolbox to detect multi-steady state behaviour in signalling pathways and gene regulatory networks. The tool combines results from chemical reaction network theory with global optimization to efficiently detect whether a signalling pathway has the capacity to undergo a saddle node bifurcation, and in case of multi-stationarity, provides the exact coordinates of the bifurcation where to start a numerical continuation analysis with standard bifurcation tools, leading to two different branches of equilibria. Bistability detection in the G1/S transition pathway of Saccharomyces cerevisiae is included as an illustrative example.

AVAILABILITY AND IMPLEMENTATION

BioSwitch runs under the popular MATLAB computational environment, and is available at https://sites.google.com/view/bioswitch.

A dual-parameter identification approach for data-based predictive modeling of hybrid gene regulatory network-growth kinetics in Pseudomonas putida mt-2

Data integration to model-based description of biological systems incorporating gene dynamics improves the performance of microbial systems. Bioprocess performance, typically predicted using empirical Monod-type models, is essential for a sustainable bioeconomy. To replace empirical models, we updated a hybrid gene regulatory network-growth kinetic model, predicting aromatic pollutants degradation and biomass growth in Pseudomonas putida mt-2. We modeled a complex biological system including extensive information to understand the role of the regulatory elements in toluene biodegradation and biomass growth. The updated model exhibited extra complications such as the existence of oscillations and discontinuities. As parameter estimation of complex biological models remains a key challenge, we used the updated model to present a dual-parameter identification approach (the 'dual approach') combining two independent methodologies. Approach I handled the complexity by incorporation of demonstrated biological knowledge in the model-development process and combination of global sensitivity analysis and optimisation. Approach II complemented Approach I handling multimodality, ill-conditioning and overfitting through regularisation estimation, global optimisation, and identifiability analysis. To systematically quantify the biological system, we used a vast amount of high-quality time-course data. The dual approach resulted in an accurately calibrated kinetic model (NRMSE: 0.17055) efficiently handling the additional model complexity. We tested model validation using three independent experimental data sets, achieving greater predictive power (NRMSE: 0.18776) than the individual approaches (NRMSE I: 0.25322, II: 0.25227) and increasing model robustness. These results demonstrated data-driven predictive modeling potentially leading to bioprocess' model-based control and optimisation.

Global optimization of the Michaelis-Menten parameters using physiologically-based pharmacokinetic (PBPK) modeling and chloroform vapor uptake data in F344 rats

Objective: To quantify metabolism, a physiologically based pharmacokinetic (PBPK) model for a volatile compound can be calibrated with the closed chamber (i.e. vapor uptake) inhalation data. Here, we introduce global optimization as a novel component of the predictive process and use it to illustrate a procedure for metabolic parameter estimation.Materials and methods: Male F344 rats were exposed in vapor uptake chambers to initial concentrations of 100, 500, 1000, and 3000 ppm chloroform. Chamber time-course data from these experiments, in combination with optimization using a chemical-specific PBPK model, were used to estimate Michaelis-Menten metabolic constants. Matlab^® simulation software was used to integrate the mass balance equations and to perform the global optimizations using MEIGO (MEtaheuristics for systems biology and bIoinformatics Global Optimization - Version 64 bit, R2016A), a toolbox written for Matlab®. The cost function used the chamber time-course data and least squares to minimize the difference between data and simulation values.Results and discussion: The final values estimated for V_max (maximum metabolic rate) and K_m (affinity constant) were 1.2 mg/h and a range between 0.0005 and 0.6 mg/L, respectively. Also, cost function plots were used to analyze the dose-dependent capacity to estimate V_max and K_m within the experimental range used. Sensitivity analysis was used to assess identifiability for both parameters and show these kinetic data may not be sufficient to identify K_m.Conclusion: In summary, this work should help toxicologists interested in optimization techniques understand the overall process employed when calibrating metabolic parameters in a PBPK model with inhalation data.

Patient-specific logic models of signaling pathways from screenings on cancer biopsies to prioritize personalized combination therapies

Mechanistic modeling of signaling pathways mediating patient-specific response to therapy can help to unveil resistance mechanisms and improve therapeutic strategies. Yet, creating such models for patients, in particular for solid malignancies, is challenging. A major hurdle to build these models is the limited material available that precludes the generation of large-scale perturbation data. Here, we present an approach that couples ex vivo high-throughput screenings of cancer biopsies using microfluidics with logic-based modeling to generate patient-specific dynamic models of extrinsic and intrinsic apoptosis signaling pathways. We used the resulting models to investigate heterogeneity in pancreatic cancer patients, showing dissimilarities especially in the PI3K-Akt pathway. Variation in model parameters reflected well the different tumor stages. Finally, we used our dynamic models to efficaciously predict new personalized combinatorial treatments. Our results suggest that our combination of microfluidic experiments and mathematical model can be a novel tool toward cancer precision medicine.

Benchmarking optimization methods for parameter estimation in large kinetic models

Motivation

Kinetic models contain unknown parameters that are estimated by optimizing the fit to experimental data. This task can be computationally challenging due to the presence of local optima and ill-conditioning. While a variety of optimization methods have been suggested to surmount these issues, it is difficult to choose the best one for a given problem a priori. A systematic comparison of parameter estimation methods for problems with tens to hundreds of optimization variables is currently missing, and smaller studies provided contradictory findings.

Results

We use a collection of benchmarks to evaluate the performance of two families of optimization methods: (i) multi-starts of deterministic local searches and (ii) stochastic global optimization metaheuristics; the latter may be combined with deterministic local searches, leading to hybrid methods. A fair comparison is ensured through a collaborative evaluation and a consideration of multiple performance metrics. We discuss possible evaluation criteria to assess the trade-off between computational efficiency and robustness. Our results show that, thanks to recent advances in the calculation of parametric sensitivities, a multi-start of gradient-based local methods is often a successful strategy, but a better performance can be obtained with a hybrid metaheuristic. The best performer combines a global scatter search metaheuristic with an interior point local method, provided with gradients estimated with adjoint-based sensitivities. We provide an implementation of this method to render it available to the scientific community.

Availability and implementation

The code to reproduce the results is provided as Supplementary Material and is available at Zenodo https://doi.org/10.5281/zenodo.1304034.

Supplementary information

Supplementary data are available at Bioinformatics online.

Parameter Estimation and Uncertainty Quantification for Systems Biology Models

Mathematical models can provide quantitative insights into immunoreceptor signaling, and other biological processes, but require parameterization and uncertainty quantification before reliable predictions become possible. We review currently available methods and software tools to address these problems. We consider gradient-based and gradient-free methods for point estimation of parameter values, and methods of profile likelihood, bootstrapping, and Bayesian inference for uncertainty quantification. We consider recent and potential future applications of these methods to systems-level modeling of immune-related phenomena.

Guidelines for the design of (optimal) isothermal inactivation experiments

Kinetic models are nowadays a basic tool to ensure food safety. Most models used in predictive microbiology have model parameters, whose precision is crucial to provide meaningful predictions. Kinetic parameters are usually estimated based on experimental data, where the experimental design can have a great impact on the precision of the estimates. In this sense, Optimal Experiment Design (OED) applies tools from optimization and information theory to identify the most informative experiment under a set of constrains (e.g. mathematical model, number of samples, etc). In this work, we develop a methodology for the design of optimal isothermal inactivation experiments. We consider the two dimensions of the design space (time and temperature), as well as a temperature-dependent maximum duration of the experiment. Functions for its application have been included in the bioOED R package. We identify design patterns that remain optimum regardless of the number of sampling points for three inactivation models (Bigelow, Mafart and Peleg) and three model microorganisms (Escherichia coli, Salmonella Senftenberg and Bacillus coagulans). Samples at extreme temperatures and close to the maximum duration of the experiment are the most informative. Moreover, the Mafart and Peleg models require some samples at intermediate time points due to the non-linearity of the survivor curve. The impact of the reference temperature on the precision of the parameter estimates is also analysed. Based on numerical simulations we recommend fixing it to the mean of the maximum and minimum temperatures used for the experiments. The article ends with a discussion presenting guidelines for the design of isothermal inactivation experiments. They combine these optimum results based on information theory with several practical limitations related to isothermal inactivation experiments. The application of these guidelines would reduce the experimental burden required to characterize thermal inactivation.

PyBioNetFit and the Biological Property Specification Language

In systems biology modeling, important steps include model parameterization, uncertainty quantification, and evaluation of agreement with experimental observations. To help modelers perform these steps, we developed the software PyBioNetFit, which in addition supports checking models against known system properties and solving design problems. PyBioNetFit introduces Biological Property Specification Language (BPSL) for the formal declaration of system properties. BPSL allows qualitative data to be used alone or in combination with quantitative data. PyBioNetFit performs parameterization with parallelized metaheuristic optimization algorithms that work directly with existing model definition standards: BioNetGen Language (BNGL) and Systems Biology Markup Language (SBML). We demonstrate PyBioNetFit's capabilities by solving various example problems, including the challenging problem of parameterizing a 153-parameter model of cell cycle control in yeast based on both quantitative and qualitative data. We demonstrate the model checking and design applications of PyBioNetFit and BPSL by analyzing a model of targeted drug interventions in autophagy signaling.

An Improved Scatter Search Algorithm for Parameter Estimation in Large-Scale Kinetic Models of Biochemical Systems

Background:

Mathematical models play a central role in facilitating researchers to better understand and comprehensively analyze various processes in biochemical systems. Their usage is beneficial in metabolic engineering as they help predict and improve desired products. However, one of the primary challenges in model building is parameter estimation. It is the process to find nearoptimal values of kinetic parameters which may culminate in the best fit of model prediction to experimental data.

Methods:

This paper proposes an improved scatter search algorithm to address the challenging parameter estimation problem. The improved algorithm is based on hybridization of quasi opposition-based learning in enhanced scatter search (QOBLESS) method. The algorithm is tested using a large-scale metabolic model of Chinese Hamster Ovary (CHO) cells.

Results:

The experimental result shows that the proposed algorithm performs better than other algorithms in terms of convergence speed and the minimum value of the objective function (loglikelihood). The estimated parameters from the experiment produce a better model by means of obtaining a reasonable good fit of model prediction to the experimental data.

Conclusion:

The kinetic parameters’ value obtained from our work was able to result in a reasonable best fit of model prediction to the experimental data, which contributes to a better understanding and produced more accurate model. Based on the results, the QOBLESS method can be used as an efficient parameter estimation method in large-scale kinetic model building.

PESTO: Parameter EStimation TOolbox

Summary

PESTO is a widely applicable and highly customizable toolbox for parameter estimation in MathWorks MATLAB. It offers scalable algorithms for optimization, uncertainty and identifiability analysis, which work in a very generic manner, treating the objective function as a black box. Hence, PESTO can be used for any parameter estimation problem, for which the user can provide a deterministic objective function in MATLAB.

Availability and implementation

PESTO is a MATLAB toolbox, freely available under the BSD license. The source code, along with extensive documentation and example code, can be downloaded from https://github.com/ICB-DCM/PESTO/.

Supplementary information

Supplementary data are available at Bioinformatics online.

Parameter estimation in models of biological oscillators: an automated regularised estimation approach

BACKGROUND

Dynamic modelling is a core element in the systems biology approach to understanding complex biosystems. Here, we consider the problem of parameter estimation in models of biological oscillators described by deterministic nonlinear differential equations. These problems can be extremely challenging due to several common pitfalls: (i) a lack of prior knowledge about parameters (i.e. massive search spaces), (ii) convergence to local optima (due to multimodality of the cost function), (iii) overfitting (fitting the noise instead of the signal) and (iv) a lack of identifiability. As a consequence, the use of standard estimation methods (such as gradient-based local ones) will often result in wrong solutions. Overfitting can be particularly problematic, since it produces very good calibrations, giving the impression of an excellent result. However, overfitted models exhibit poor predictive power. Here, we present a novel automated approach to overcome these pitfalls. Its workflow makes use of two sequential optimisation steps incorporating three key algorithms: (1) sampling strategies to systematically tighten the parameter bounds reducing the search space, (2) efficient global optimisation to avoid convergence to local solutions, (3) an advanced regularisation technique to fight overfitting. In addition, this workflow incorporates tests for structural and practical identifiability.

RESULTS

We successfully evaluate this novel approach considering four difficult case studies regarding the calibration of well-known biological oscillators (Goodwin, FitzHugh-Nagumo, Repressilator and a metabolic oscillator). In contrast, we show how local gradient-based approaches, even if used in multi-start fashion, are unable to avoid the above-mentioned pitfalls.

CONCLUSIONS

Our approach results in more efficient estimations (thanks to the bounding strategy) which are able to escape convergence to local optima (thanks to the global optimisation approach). Further, the use of regularisation allows us to avoid overfitting, resulting in more generalisable calibrated models (i.e. models with greater predictive power).

Workflow Description to Dynamically Model β-Arrestin Signaling Networks

Dynamic models of signaling networks allow the formulation of hypotheses on the topology and kinetic rate laws characterizing a given molecular network, in-depth exploration, and confrontation with kinetic biological data. Despite its standardization, dynamic modeling of signaling networks still requires successive technical steps that need to be carefully performed. Here, we detail these steps by going through the mathematical and statistical framework. We explain how it can be applied to the understanding of β-arrestin-dependent signaling networks. We illustrate our methodology through the modeling of β-arrestin recruitment kinetics at the follicle-stimulating hormone (FSH) receptor supported by in-house bioluminescence resonance energy transfer (BRET) data.

Scalable Inference of Ordinary Differential Equation Models of Biochemical Processes

Ordinary differential equation models have become a standard tool for the mechanistic description of biochemical processes. If parameters are inferred from experimental data, such mechanistic models can provide accurate predictions about the behavior of latent variables or the process under new experimental conditions. Complementarily, inference of model structure can be used to identify the most plausible model structure from a set of candidates, and, thus, gain novel biological insight. Several toolboxes can infer model parameters and structure for small- to medium-scale mechanistic models out of the box. However, models for highly multiplexed datasets can require hundreds to thousands of state variables and parameters. For the analysis of such large-scale models, most algorithms require intractably high computation times. This chapter provides an overview of the state-of-the-art methods for parameter and model inference, with an emphasis on scalability.

Two Interlinked Bistable Switches Govern Mitotic Control in Mammalian Cells

Distinct protein phosphorylation levels in interphase and M phase require tight regulation of Cdk1 activity [1, 2]. A bistable switch, based on positive feedback in the Cdk1 activation loop, has been proposed to generate different thresholds for transitions between these cell-cycle states [3-5]. Recently, the activity of the major Cdk1-counteracting phosphatase, PP2A:B55, has also been found to be bistable due to Greatwall kinase-dependent regulation [6]. However, the interplay of the regulation of Cdk1 and PP2A:B55 in vivo remains unexplored. Here, we combine quantitative cell biology assays with mathematical modeling to explore the interplay of mitotic kinase activation and phosphatase inactivation in human cells. By measuring mitotic entry and exit thresholds using ATP-analog-sensitive Cdk1 mutants, we find evidence that the mitotic switch displays hysteresis and bistability, responding differentially to Cdk1 inhibition in the mitotic and interphase states. Cdk1 activation by Wee1/Cdc25 feedback loops and PP2A:B55 inactivation by Greatwall independently contributes to this hysteretic switch system. However, elimination of both Cdk1 and PP2A:B55 inactivation fully abrogates bistability, suggesting that hysteresis is an emergent property of mutual inhibition between the Cdk1 and PP2A:B55 feedback loops. Our model of the two interlinked feedback systems predicts an intermediate but hidden steady state between interphase and M phase. This could be verified experimentally by Cdk1 inhibition during mitotic entry, supporting the predictive value of our model. Furthermore, we demonstrate that dual inhibition of Wee1 and Gwl kinases causes loss of cell-cycle memory and synthetic lethality, which could be further exploited therapeutically.

Recent Advances in the Systems Biology of Aging

SIGNIFICANCE

Reductionist studies have contributed greatly to our understanding of the basic biology of aging in recent years but we still do not understand fundamental mechanisms for many identified drugs and pathways. Use of systems approaches will help us move forward in our understanding of aging. Recent Advances: Recent work described here has illustrated the power of systems biology to inform our understanding of aging through the study of (i) diet restriction, (ii) neurodegenerative disease, and (iii) biomarkers of aging.

CRITICAL ISSUES

Although we do not understand all of the individual genes and pathways that affect aging, as we continue to uncover more of them, we have now also begun to synthesize existing data using systems-level approaches, often to great effect. The three examples noted here all benefit from computational approaches that were unknown a few years ago, and from biological insights gleaned from multiple model systems, from aging laboratories as well as many other areas of biology.

FUTURE DIRECTIONS

Many new technologies, such as single-cell sequencing, advances in epigenetics beyond the methylome (specifically, assay for transposase-accessible chromatin with high throughput sequencing ), and multiomic network studies, will increase the reach of systems biologists. This suggests that approaches similar to those described here will continue to lead to striking findings, and to interventions that may allow us to delay some of the many age-associated diseases in humans; perhaps sooner that we expect. Antioxid. Redox Signal. 29, 973-984.

Glycosylation Flux Analysis of Immunoglobulin G in Chinese Hamster Ovary Perfusion Cell Culture

The terminal sugar molecules of the N-linked glycan attached to the fragment crystalizable (Fc) region is a critical quality attribute of therapeutic monoclonal antibodies (mAbs) such as immunoglobulin G (IgG). There exists naturally-occurring heterogeneity in the N-linked glycan structure of mAbs, and such heterogeneity has a significant influence on the clinical safety and efficacy of mAb drugs. We previously proposed a constraint-based modeling method called glycosylation flux analysis (GFA) to characterize the rates (fluxes) of intracellular glycosylation reactions. One contribution of this work is a significant improvement in the computational efficiency of the GFA, which is beneficial for analyzing large datasets. Another contribution of our study is the analysis of IgG glycosylation in continuous perfusion Chinese Hamster Ovary (CHO) cell cultures. The GFA of the perfusion cell culture data indicated that the dynamical changes of IgG glycan heterogeneity are mostly attributed to alterations in the galactosylation flux activity. By using a random forest regression analysis of the IgG galactosylation flux activity, we were further able to link the dynamics of galactosylation with two process parameters: cell-specific productivity of IgG and extracellular ammonia concentration. The characteristics of IgG galactosylation dynamics agree well with what we previously reported for fed-batch cultivations of the same CHO cell strain.

Using both qualitative and quantitative data in parameter identification for systems biology models

In systems biology, qualitative data are often generated, but rarely used to parameterize models. We demonstrate an approach in which qualitative and quantitative data can be combined for parameter identification. In this approach, qualitative data are converted into inequality constraints imposed on the outputs of the model. These inequalities are used along with quantitative data points to construct a single scalar objective function that accounts for both datasets. To illustrate the approach, we estimate parameters for a simple model describing Raf activation. We then apply the technique to a more elaborate model characterizing cell cycle regulation in yeast. We incorporate both quantitative time courses (561 data points) and qualitative phenotypes of 119 mutant yeast strains (1647 inequalities) to perform automated identification of 153 model parameters. We quantify parameter uncertainty using a profile likelihood approach. Our results indicate the value of combining qualitative and quantitative data to parameterize systems biology models.

Much of the data generated in biology is qualitative, but exploiting such data to inform models of biological systems remains a challenge. Here, the authors demonstrate an approach that allows use of both quantitative and qualitative data for parameterising dynamical models.