Fides: Reliable trust-region optimization for parameter estimation of ordinary differential equation models

Ordinary differential equation (ODE) models are widely used to study biochemical reactions in cellular networks since they effectively describe the temporal evolution of these networks using mass action kinetics. The parameters of these models are rarely known a priori and must instead be estimated by calibration using experimental data. Optimization-based calibration of ODE models on is often challenging, even for low-dimensional problems. Multiple hypotheses have been advanced to explain why biochemical model calibration is challenging, including non-identifiability of model parameters, but there are few comprehensive studies that test these hypotheses, likely because tools for performing such studies are also lacking. Nonetheless, reliable model calibration is essential for uncertainty analysis, model comparison, and biological interpretation.We implemented an established trust-region method as a modular Python framework (fides) to enable systematic comparison of different approaches to ODE model calibration involving a variety of Hessian approximation schemes. We evaluated fides on a recently developed corpus of biologically realistic benchmark problems for which real experimental data are available. Unexpectedly, we observed high variability in optimizer performance among different implementations of the same mathematical instructions (algorithms). Analysis of possible sources of poor optimizer performance identified limitations in the widely used Gauss-Newton, BFGS and SR1 Hessian approximation schemes. We addressed these drawbacks with a novel hybrid Hessian approximation scheme that enhances optimizer performance and outperforms existing hybrid approaches. When applied to the corpus of test models, we found that fides was on average more reliable and efficient than existing methods using a variety of criteria. We expect fides to be broadly useful for ODE constrained optimization problems in biochemical models and to be a foundation for future methods development.     

medplot: A Web Application for Dynamic Summary and Analysis of Longitudinal Medical Data Based on R

In biomedical studies the patients are often evaluated numerous times and a large number of variables are recorded at each time-point. Data entry and manipulation of longitudinal data can be performed using spreadsheet programs, which usually include some data plotting and analysis capabilities and are straightforward to use, but are not designed for the analyses of complex longitudinal data. Specialized statistical software offers more flexibility and capabilities, but first time users with biomedical background often find its use difficult. We developed medplot, an interactive web application that simplifies the exploration and analysis of longitudinal data. The application can be used to summarize, visualize and analyze data by researchers that are not familiar with statistical programs and whose knowledge of statistics is limited. The summary tools produce publication-ready tables and graphs. The analysis tools include features that are seldom available in spreadsheet software, such as correction for multiple testing, repeated measurement analyses and flexible non-linear modeling of the association of the numerical variables with the outcome. medplot is freely available and open source, it has an intuitive graphical user interface (GUI), it is accessible via the Internet and can be used within a web browser, without the need for installing and maintaining programs locally on the user’s computer. This paper describes the application and gives detailed examples describing how to use the application on real data from a clinical study including patients with early Lyme borreliosis.     

Conceptual-level workflow modeling of scientific experiments using NMR as a case study

Background  

Scientific workflows improve the process of scientific experiments by making computations explicit, underscoring data flow, and emphasizing the participation of humans in the process when intuition and human reasoning are required. Workflows for experiments also highlight transitions among experimental phases, allowing intermediate results to be verified and supporting the proper handling of semantic mismatches and different file formats among the various tools used in the scientific process. Thus, scientific workflows are important for the modeling and subsequent capture of bioinformatics-related data. While much research has been conducted on the implementation of scientific workflows, the initial process of actually designing and generating the workflow at the conceptual level has received little consideration.     

Using research projects and qualitative conceptual modeling to increase novice scientists’ understanding of ecological complexity

A new program, Teaching Ecological Complexity, is working to develop a heightened capacity for systems thinking among high school biology and environmental sciences teachers. During a 2-week field-based course, the teachers use qualitative conceptual modeling, participate in all stages of field experimentation, and formulate plans to teach field research with their own classes. Qualitative conceptual modeling was found to be useful in revealing the underlying perceptions of ecosystem functioning for these novice scientists. Preliminary results showed improvement in their ability to recognize and apply some of the attributes of complex ecosystem: non-linear feedback loops, hierarchical organization, patterns illustrating the spatial arrangements of species diversity. In addition to using models, teachers used peer-learning techniques. Collegial discussions about what they understood at particular points in time were useful in improving their understanding of ecosystem phenomenon.     

A kinetic-based sigmoidal model for the polymerase chain reaction and its application to high-capacity absolute quantitative real-time PCR

Background

Based upon defining a common reference point, current real-time quantitative PCR technologies compare relative differences in amplification profile position. As such, absolute quantification requires construction of target-specific standard curves that are highly resource intensive and prone to introducing quantitative errors. Sigmoidal modeling using nonlinear regression has previously demonstrated that absolute quantification can be accomplished without standard curves; however, quantitative errors caused by distortions within the plateau phase have impeded effective implementation of this alternative approach.

Results

Recognition that amplification rate is linearly correlated to amplicon quantity led to the derivation of two sigmoid functions that allow target quantification via linear regression analysis. In addition to circumventing quantitative errors produced by plateau distortions, this approach allows the amplification efficiency within individual amplification reactions to be determined. Absolute quantification is accomplished by first converting individual fluorescence readings into target quantity expressed in fluorescence units, followed by conversion into the number of target molecules via optical calibration. Founded upon expressing reaction fluorescence in relation to amplicon DNA mass, a seminal element of this study was to implement optical calibration using lambda gDNA as a universal quantitative standard. Not only does this eliminate the need to prepare target-specific quantitative standards, it relegates establishment of quantitative scale to a single, highly defined entity. The quantitative competency of this approach was assessed by exploiting "limiting dilution assay" for absolute quantification, which provided an independent gold standard from which to verify quantitative accuracy. This yielded substantive corroborating evidence that absolute accuracies of ± 25% can be routinely achieved. Comparison with the LinReg and Miner automated qPCR data processing packages further demonstrated the superior performance of this kinetic-based methodology.

Conclusion

Called "linear regression of efficiency" or LRE, this novel kinetic approach confers the ability to conduct high-capacity absolute quantification with unprecedented quality control capabilities. The computational simplicity and recursive nature of LRE quantification also makes it amenable to software implementation, as demonstrated by a prototypic Java program that automates data analysis. This in turn introduces the prospect of conducting absolute quantification with little additional effort beyond that required for the preparation of the amplification reactions.     

Estimation of Cell Proliferation Dynamics Using CFSE Data

Advances in fluorescent labeling of cells as measured by flow cytometry have allowed for quantitative studies of proliferating populations of cells. The investigations (Luzyanina et al. in J. Math. Biol. 54:57–89, 2007; J. Math. Biol., 2009; Theor. Biol. Med. Model. 4:1–26, 2007) contain a mathematical model with fluorescence intensity as a structure variable to describe the evolution in time of proliferating cells labeled by carboxyfluorescein succinimidyl ester (CFSE). Here, this model and several extensions/modifications are discussed. Suggestions for improvements are presented and analyzed with respect to statistical significance for better agreement between model solutions and experimental data. These investigations suggest that the new decay/label loss and time dependent effective proliferation and death rates do indeed provide improved fits of the model to data. Statistical models for the observed variability/noise in the data are discussed with implications for uncertainty quantification. The resulting new cell dynamics model should prove useful in proliferation assay tracking and modeling, with numerous applications in the biomedical sciences.     

popharvest: An R package to assess the sustainability of harvesting regimes of bird populations

Bird harvest for recreational purposes or as a source for food is an important activity worldwide. Assessing or mitigating the impact of these additional sources of mortality on bird populations is therefore crucial issue. The sustainability of harvest levels is however rarely documented, because knowledge of their population dynamics remains rudimentary for many bird species. Some helpful approaches using limited demographic data can be used to provide initial assessment of the sustainable use of harvested bird populations, and help adjusting harvest levels accordingly. The Demographic Invariant Method (DIM) is used to detect overharvesting. In complement, the Potential Take Level (PTL) approach may allow setting a level of take with regard to management objectives and/or to assess whether current harvest levels meet these objectives. Here, we present the R package popharvest that implements these two approaches in a simple and straightforward way. The package provides users with a set of flexible functions whose arguments can be adapted to existing knowledge about population dynamics. Also, popharvest enables users to test scenarios or propagate uncertainty in demographic parameters to the assessment of sustainability through easily programming Monte Carlo simulations. The simplicity of the package makes it a useful toolbox for wildlife managers or policymakers. This paper provides them with backgrounds about the DIM and PTL approaches and illustrates the use of popharvest''s functionalities in this context.     

Efficient ancestry and mutation simulation with msprime 1.0

Stochastic simulation is a key tool in population genetics, since the models involved are often analytically intractable and simulation is usually the only way of obtaining ground-truth data to evaluate inferences. Because of this, a large number of specialized simulation programs have been developed, each filling a particular niche, but with largely overlapping functionality and a substantial duplication of effort. Here, we introduce msprime version 1.0, which efficiently implements ancestry and mutation simulations based on the succinct tree sequence data structure and the tskit library. We summarize msprime’s many features, and show that its performance is excellent, often many times faster and more memory efficient than specialized alternatives. These high-performance features have been thoroughly tested and validated, and built using a collaborative, open source development model, which reduces duplication of effort and promotes software quality via community engagement.     

Analysing and meta-analysing time-series data of microbial growth and gene expression from plate readers

Responding to change is a fundamental property of life, making time-series data invaluable in biology. For microbes, plate readers are a popular, convenient means to measure growth and also gene expression using fluorescent reporters. Nevertheless, the difficulties of analysing the resulting data can be a bottleneck, particularly when combining measurements from different wells and plates. Here we present omniplate, a Python module that corrects and normalises plate-reader data, estimates growth rates and fluorescence per cell as functions of time, calculates errors, exports in different formats, and enables meta-analysis of multiple plates. The software corrects for autofluorescence, the optical density’s non-linear dependence on the number of cells, and the effects of the media. We use omniplate to measure the Monod relationship for the growth of budding yeast in raffinose, showing that raffinose is a convenient carbon source for controlling growth rates. Using fluorescent tagging, we study yeast’s glucose transport. Our results are consistent with the regulation of the hexose transporter (HXT) genes being approximately bipartite: the medium and high affinity transporters are predominately regulated by both the high affinity glucose sensor Snf3 and the kinase complex SNF1 via the repressors Mth1, Mig1, and Mig2; the low affinity transporters are predominately regulated by the low affinity sensor Rgt2 via the co-repressor Std1. We thus demonstrate that omniplate is a powerful tool for exploiting the advantages offered by time-series data in revealing biological regulation.     

MINT: software to identify motifs and short-range interactions in trajectories of nucleic acids

Structural biology experiments and structure prediction tools have provided many high-resolution three-dimensional structures of nucleic acids. Also, molecular dynamics force field parameters have been adapted to simulating charged and flexible nucleic acid structures on microsecond time scales. Therefore, we can generate the dynamics of DNA or RNA molecules, but we still lack adequate tools for the analysis of the resulting huge amounts of data. We present MINT (Motif Identifier for Nucleic acids Trajectory) — an automatic tool for analyzing three-dimensional structures of RNA and DNA, and their full-atom molecular dynamics trajectories or other conformation sets (e.g. X-ray or nuclear magnetic resonance-derived structures). For each RNA or DNA conformation MINT determines the hydrogen bonding network resolving the base pairing patterns, identifies secondary structure motifs (helices, junctions, loops, etc.) and pseudoknots. MINT also estimates the energy of stacking and phosphate anion-base interactions. For many conformations, as in a molecular dynamics trajectory, MINT provides averages of the above structural and energetic features and their evolution. We show MINT functionality based on all-atom explicit solvent molecular dynamics trajectory of the 30S ribosomal subunit.     

Ensemble modeling for analysis of cell signaling dynamics

Systems biology iteratively combines experimentation with mathematical modeling. However, limited mechanistic knowledge, conflicting hypotheses and scarce experimental data severely hamper the development of predictive mechanistic models in many areas of biology. Even under such high uncertainty, we show here that ensemble modeling, when combined with targeted experimental analysis, can unravel key operating principles in complex cellular pathways. For proof of concept, we develop a library of mechanistically alternative dynamic models for the highly conserved target-of-rapamycin (TOR) pathway of Saccharomyces cerevisiae. In contrast to the prevailing view of a de novo assembly of type 2A phosphatases (PP2As), our integrated computational and experimental analysis proposes a specificity factor, based on Tap42p-Tip41p, for PP2As as the key signaling mechanism that is quantitatively consistent with all available experimental data. Beyond revising our picture of TOR signaling, we expect ensemble modeling to help elucidate other insufficiently characterized cellular circuits.     

VA-Index: Quantifying Assortativity Patterns in Networks with Multidimensional Nodal Attributes

Network connections have been shown to be correlated with structural or external attributes of the network vertices in a variety of cases. Given the prevalence of this phenomenon network scientists have developed metrics to quantify its extent. In particular, the assortativity coefficient is used to capture the level of correlation between a single-dimensional attribute (categorical or scalar) of the network nodes and the observed connections, i.e., the edges. Nevertheless, in many cases a multi-dimensional, i.e., vector feature of the nodes is of interest. Similar attributes can describe complex behavioral patterns (e.g., mobility) of the network entities. To date little attention has been given to this setting and there has not been a general and formal treatment of this problem. In this study we develop a metric, the vector assortativity index (VA-index for short), based on network randomization and (empirical) statistical hypothesis testing that is able to quantify the assortativity patterns of a network with respect to a vector attribute. Our extensive experimental results on synthetic network data show that the VA-index outperforms a baseline extension of the assortativity coefficient, which has been used in the literature to cope with similar cases. Furthermore, the VA-index can be calibrated (in terms of parameters) fairly easy, while its benefits increase with the (co-)variance of the vector elements, where the baseline systematically over(under)estimate the true mixing patterns of the network.     

GPRuler: Metabolic gene-protein-reaction rules automatic reconstruction

Metabolic network models are increasingly being used in health care and industry. As a consequence, many tools have been released to automate their reconstruction process de novo. In order to enable gene deletion simulations and integration of gene expression data, these networks must include gene-protein-reaction (GPR) rules, which describe with a Boolean logic relationships between the gene products (e.g., enzyme isoforms or subunits) associated with the catalysis of a given reaction. Nevertheless, the reconstruction of GPRs still remains a largely manual and time consuming process. Aiming at fully automating the reconstruction process of GPRs for any organism, we propose the open-source python-based framework GPRuler. By mining text and data from 9 different biological databases, GPRuler can reconstruct GPRs starting either from just the name of the target organism or from an existing metabolic model. The performance of the developed tool is evaluated at small-scale level for a manually curated metabolic model, and at genome-scale level for three metabolic models related to Homo sapiens and Saccharomyces cerevisiae organisms. By exploiting these models as benchmarks, the proposed tool shown its ability to reproduce the original GPR rules with a high level of accuracy. In all the tested scenarios, after a manual investigation of the mismatches between the rules proposed by GPRuler and the original ones, the proposed approach revealed to be in many cases more accurate than the original models. By complementing existing tools for metabolic network reconstruction with the possibility to reconstruct GPRs quickly and with a few resources, GPRuler paves the way to the study of context-specific metabolic networks, representing the active portion of the complete network in given conditions, for organisms of industrial or biomedical interest that have not been characterized metabolically yet.     

Standardisation of data from real-time quantitative PCR methods – evaluation of outliers and comparison of calibration curves

Background  

As real-time quantitative PCR (RT-QPCR) is increasingly being relied upon for the enforcement of legislation and regulations dependent upon the trace detection of DNA, focus has increased on the quality issues related to the technique. Recent work has focused on the identification of factors that contribute towards significant measurement uncertainty in the real-time quantitative PCR technique, through investigation of the experimental design and operating procedure. However, measurement uncertainty contributions made during the data analysis procedure have not been studied in detail. This paper presents two additional approaches for standardising data analysis through the novel application of statistical methods to RT-QPCR, in order to minimise potential uncertainty in results.     

Analysis of protrusion dynamics in amoeboid cell motility by means of regularized contour flows

Amoeboid cell motility is essential for a wide range of biological processes including wound healing, embryonic morphogenesis, and cancer metastasis. It relies on complex dynamical patterns of cell shape changes that pose long-standing challenges to mathematical modeling and raise a need for automated and reproducible approaches to extract quantitative morphological features from image sequences. Here, we introduce a theoretical framework and a computational method for obtaining smooth representations of the spatiotemporal contour dynamics from stacks of segmented microscopy images. Based on a Gaussian process regression we propose a one-parameter family of regularized contour flows that allows us to continuously track reference points (virtual markers) between successive cell contours. We use this approach to define a coordinate system on the moving cell boundary and to represent different local geometric quantities in this frame of reference. In particular, we introduce the local marker dispersion as a measure to identify localized membrane expansions and provide a fully automated way to extract the properties of such expansions, including their area and growth time. The methods are available as an open-source software package called AmoePy, a Python-based toolbox for analyzing amoeboid cell motility (based on time-lapse microscopy data), including a graphical user interface and detailed documentation. Due to the mathematical rigor of our framework, we envision it to be of use for the development of novel cell motility models. We mainly use experimental data of the social amoeba Dictyostelium discoideum to illustrate and validate our approach.     

INCA 2.0: A tool for integrated,dynamic modeling of NMR- and MS-based isotopomer measurements and rigorous metabolic flux analysis

Metabolic flux analysis (MFA) combines experimental measurements and computational modeling to determine biochemical reaction rates in live biological systems. Advancements in analytical instrumentation, such as nuclear magnetic resonance (NMR) spectroscopy and mass spectrometry (MS), have facilitated chemical separation and quantification of isotopically enriched metabolites. However, no software packages have been previously described that can integrate isotopomer measurements from both MS and NMR analytical platforms and have the flexibility to estimate metabolic fluxes from either isotopic steady-state or dynamic labeling experiments. By applying physiologically relevant cardiac and hepatic metabolic models to assess NMR isotopomer measurements, we herein test and validate new modeling capabilities of our enhanced flux analysis software tool, INCA 2.0. We demonstrate that INCA 2.0 can simulate and regress steady-state ¹³C NMR datasets from perfused hearts with an accuracy comparable to other established flux assessment tools. Furthermore, by simulating the infusion of three different ¹³C acetate tracers, we show that MFA based on dynamic ¹³C NMR measurements can more precisely resolve cardiac fluxes compared to isotopically steady-state flux analysis. Finally, we show that estimation of hepatic fluxes using combined ¹³C NMR and MS datasets improves the precision of estimated fluxes by up to 50%. Overall, our results illustrate how the recently added NMR data modeling capabilities of INCA 2.0 can enable entirely new experimental designs that lead to improved flux resolution and can be applied to a wide range of biological systems and measurement time courses.     

msCompare: a framework for quantitative analysis of label-free LC-MS data for comparative candidate biomarker studies

Data processing forms an integral part of biomarker discovery and contributes significantly to the ultimate result. To compare and evaluate various publicly available open source label-free data processing workflows, we developed msCompare, a modular framework that allows the arbitrary combination of different feature detection/quantification and alignment/matching algorithms in conjunction with a novel scoring method to evaluate their overall performance. We used msCompare to assess the performance of workflows built from modules of publicly available data processing packages such as SuperHirn, OpenMS, and MZmine and our in-house developed modules on peptide-spiked urine and trypsin-digested cerebrospinal fluid (CSF) samples. We found that the quality of results varied greatly among workflows, and interestingly, heterogeneous combinations of algorithms often performed better than the homogenous workflows. Our scoring method showed that the union of feature matrices of different workflows outperformed the original homogenous workflows in some cases. msCompare is open source software (https://trac.nbic.nl/mscompare), and we provide a web-based data processing service for our framework by integration into the Galaxy server of the Netherlands Bioinformatics Center (http://galaxy.nbic.nl/galaxy) to allow scientists to determine which combination of modules provides the most accurate processing for their particular LC-MS data sets.