An error model for protein quantification

A 96-Well Polyacrylamide Gel for Electrophoresis and Western Blotting

Western blotting is a stalwart technique for analyzing specific proteins and/or their post-translational modifications. However, it remains challenging to accommodate more than ~10 samples per experiment without substantial departure from trusted, established protocols involving accessible instrumentation. Here, we describe a 96-sample western blot that conforms to standard 96-well plate dimensional constraints and has little operational deviation from standard western blotting. The main differences are that (i) submerged polyacrylamide gel electrophoresis is operated horizontally (similar to agarose gels) as opposed to vertically, and (ii) a 6 mm thick gel is used, with 2 mm most relevant for membrane transfer (vs ~1 mm typical). Results demonstrate both wet and semi-dry transfer are compatible with this gel thickness. The major tradeoff is reduced molecular weight resolution, due primarily to less available migration distance per sample. We demonstrate proof-of-principle using gels loaded with molecular weight ladder, recombinant protein, and cell lysates. We expect the 96-well western blot will increase reproducibility, efficiency (cost and time ~8-fold), and capacity for biological characterization relative to established western blots.

Efficient parameter estimation for ODE models of cellular processes using semi-quantitative data

MOTIVATION

Quantitative dynamical models facilitate the understanding of biological processes and the prediction of their dynamics. The parameters of these models are commonly estimated from experimental data. Yet, experimental data generated from different techniques do not provide direct information about the state of the system but a nonlinear (monotonic) transformation of it. For such semi-quantitative data, when this transformation is unknown, it is not apparent how the model simulations and the experimental data can be compared.

RESULTS

We propose a versatile spline-based approach for the integration of a broad spectrum of semi-quantitative data into parameter estimation. We derive analytical formulas for the gradients of the hierarchical objective function and show that this substantially increases the estimation efficiency. Subsequently, we demonstrate that the method allows for the reliable discovery of unknown measurement transformations. Furthermore, we show that this approach can significantly improve the parameter inference based on semi-quantitative data in comparison to available methods.

AVAILABILITY AND IMPLEMENTATION

Modelers can easily apply our method by using our implementation in the open-source Python Parameter EStimation TOolbox (pyPESTO) available at https://github.com/ICB-DCM/pyPESTO.

Bayesian modelling of time series data (BayModTS)-a FAIR workflow to process sparse and highly variable data

MOTIVATION

Systems biology aims to better understand living systems through mathematical modelling of experimental and clinical data. A pervasive challenge in quantitative dynamical modelling is the integration of time series measurements, which often have high variability and low sampling resolution. Approaches are required to utilize such information while consistently handling uncertainties.

RESULTS

We present BayModTS (Bayesian modelling of time series data), a new FAIR (findable, accessible, interoperable, and reusable) workflow for processing and analysing sparse and highly variable time series data. BayModTS consistently transfers uncertainties from data to model predictions, including process knowledge via parameterized models. Further, credible differences in the dynamics of different conditions can be identified by filtering noise. To demonstrate the power and versatility of BayModTS, we applied it to three hepatic datasets gathered from three different species and with different measurement techniques: (i) blood perfusion measurements by magnetic resonance imaging in rat livers after portal vein ligation, (ii) pharmacokinetic time series of different drugs in normal and steatotic mice, and (iii) CT-based volumetric assessment of human liver remnants after clinical liver resection.

AVAILABILITY AND IMPLEMENTATION

The BayModTS codebase is available on GitHub at https://github.com/Systems-Theory-in-Systems-Biology/BayModTS. The repository contains a Python script for the executable BayModTS workflow and a widely applicable SBML (systems biology markup language) model for retarded transient functions. In addition, all examples from the paper are included in the repository. Data and code of the application examples are stored on DaRUS: https://doi.org/10.18419/darus-3876. The raw MRI ROI voxel data were uploaded to DaRUS: https://doi.org/10.18419/darus-3878. The steatosis metabolite data are published on FairdomHub: 10.15490/fairdomhub.1.study.1070.1.

Nonself RNA rewires IFN-β signaling: A mathematical model of the innate immune response

Type I interferons (IFNs) are key coordinators of the innate immune response to viral infection, which, through activation of the transcriptional regulators STAT1 and STAT2 (STAT1/2) in bystander cells, induce the expression of IFN-stimulated genes (ISGs). Here, we showed that in cells transfected with poly(I:C), an analog of viral RNA, the transcriptional activity of STAT1/2 was terminated because of depletion of the interferon-β (IFN-β) receptor, IFNAR. Activation of RNase L and PKR, products of two ISGs, not only hindered the replenishment of IFNAR but also suppressed negative regulators of IRF3 and NF-κB, consequently promoting IFNB transcription. We incorporated these findings into a mathematical model of innate immunity. By coupling signaling through the IRF3-NF-κB and STAT1/2 pathways with the activities of RNase L and PKR, the model explains how poly(I:C) switches the transcriptional program from being STAT1/2 induced to being IRF3 and NF-κB induced, which converts IFN-β-responding cells to IFN-β-secreting cells.

Posterior marginalization accelerates Bayesian inference for dynamical models of biological processes

Bayesian inference is an important method in the life and natural sciences for learning from data. It provides information about parameter and prediction uncertainties. Yet, generating representative samples from the posterior distribution is often computationally challenging. Here, we present an approach that lowers the computational complexity of sample generation for dynamical models with scaling, offset, and noise parameters. The proposed method is based on the marginalization of the posterior distribution. We provide analytical results for a broad class of problems with conjugate priors and show that the method is suitable for a large number of applications. Subsequently, we demonstrate the benefit of the approach for applications from the field of systems biology. We report an improvement up to 50 times in the effective sample size per unit of time. As the scheme is broadly applicable, it will facilitate Bayesian inference in different research fields.

Optimizing reduction of western blotting analytical variations: Use of replicate test samples, multiple normalization methods, and sample loading positions

Western blot (WB) analysis is widely used, but obtaining consistent results can be problematic, especially when using multiple gels. This study examines WB performance by explicitly applying a method commonly used to test analytical instrumentation. Test samples were lysates from RAW 264.7 murine macrophages treated with LPS to activate MAPK and NF-kB signaling targets. Samples from the pooled cell lysates placed in every lane on multiple gels were analyzed by WBs for levels of p-ERK, ERK, IkBβ and non-target protein. Different normalization methods and sample groupings were applied to the density values and the resulting coefficients of variation (CV) and ratios of maximal to minimal values (Max/Min) were compared. Ideally with identical sample replicates the CVs would be 0 and the Max/Min 1; deviation indicating introduction of variability by the WB process. Common normalizations to reduce analytical variance, total lane protein, % Control, and p-ERK/ERK ratios, did not have the lowest CVs or Max/Min values. Normalization using the sum of target protein values combined with analytical replication most effectively reduced variability, resulting CV and Max/Min values as low as 5-10% and 1.1. These methods should allow reliable interpretation of complex experiments that require samples to be placed on multiple gels.

A plausible identifiable model of the canonical NF-κB signaling pathway

An overwhelming majority of mathematical models of regulatory pathways, including the intensively studied NF-κB pathway, remains non-identifiable, meaning that their parameters may not be determined by existing data. The existing NF-κB models that are capable of reproducing experimental data contain non-identifiable parameters, whereas simplified models with a smaller number of parameters exhibit dynamics that differs from that observed in experiments. Here, we reduced an existing model of the canonical NF-κB pathway by decreasing the number of equations from 15 to 6. The reduced model retains two negative feedback loops mediated by IκBα and A20, and in response to both tonic and pulsatile TNF stimulation exhibits dynamics that closely follow that of the original model. We carried out the sensitivity-based linear analysis and Monte Carlo-based analysis to demonstrate that the resulting model is both structurally and practically identifiable given measurements of 5 model variables from a simple TNF stimulation protocol. The reduced model is capable of reproducing different types of responses that are characteristic to regulatory motifs controlled by negative feedback loops: nearly-perfect adaptation as well as damped and sustained oscillations. It can serve as a building block of more comprehensive models of the immune response and cancer, where NF-κB plays a decisive role. Our approach, although may not be automatically generalized, suggests that models of other regulatory pathways can be transformed to identifiable, while retaining their dynamical features.

Bayesian kinetic modeling for tracer-based metabolomic data

BACKGROUND

Stable Isotope Resolved Metabolomics (SIRM) is a new biological approach that uses stable isotope tracers such as uniformly [Formula: see text]-enriched glucose ([Formula: see text]-Glc) to trace metabolic pathways or networks at the atomic level in complex biological systems. Non-steady-state kinetic modeling based on SIRM data uses sets of simultaneous ordinary differential equations (ODEs) to quantitatively characterize the dynamic behavior of metabolic networks. It has been increasingly used to understand the regulation of normal metabolism and dysregulation in the development of diseases. However, fitting a kinetic model is challenging because there are usually multiple sets of parameter values that fit the data equally well, especially for large-scale kinetic models. In addition, there is a lack of statistically rigorous methods to compare kinetic model parameters between different experimental groups.

RESULTS

We propose a new Bayesian statistical framework to enhance parameter estimation and hypothesis testing for non-steady-state kinetic modeling of SIRM data. For estimating kinetic model parameters, we leverage the prior distribution not only to allow incorporation of experts' knowledge but also to provide robust parameter estimation. We also introduce a shrinkage approach for borrowing information across the ensemble of metabolites to stably estimate the variance of an individual isotopomer. In addition, we use a component-wise adaptive Metropolis algorithm with delayed rejection to perform efficient Monte Carlo sampling of the posterior distribution over high-dimensional parameter space. For comparing kinetic model parameters between experimental groups, we propose a new reparameterization method that converts the complex hypothesis testing problem into a more tractable parameter estimation problem. We also propose an inference procedure based on credible interval and credible value. Our method is freely available for academic use at https://github.com/xuzhang0131/MCMCFlux .

CONCLUSIONS

Our new Bayesian framework provides robust estimation of kinetic model parameters and enables rigorous comparison of model parameters between experimental groups. Simulation studies and application to a lung cancer study demonstrate that our framework performs well for non-steady-state kinetic modeling of SIRM data.

Strain-dependent regulation of hippocampal long-term potentiation by dopamine D1/D5 receptors in mice

The magnitude and persistency of long-term potentiation (LTP) in the rodent hippocampus is species-dependent: rats express more robust and more prolonged LTP in response to a broader afferent frequency range than mice. The C57Bl/6 mouse is an extremely popular murine strain used in studies of hippocampal synaptic plasticity and spatial learning. Recently it was reported that it expresses impoverished LTP compared to other murine strains. Given the important role of the dopamine D1/D5 receptor (D1/D5R) in the maintenance of LTP and in memory consolidation, we explored to what extent strain-dependent differences in LTP in mice are determined by differences in D1/D5R-control. In CaOlaHsd mice, robust LTP was induced that lasted for over 24 h and which was significantly greater in magnitude than LTP induced in C57Bl/6 mice. Intracerebral treatment with a D1/D5R-antagonist (SCH23390) prevented both the early and late phase of LTP in CaOlaHsd mice, whereas only late-LTP was impaired in C57Bl/6 mice. Treatment with a D1/D5R-agonist (Chloro-PB) facilitated short-term potentiation (STP) into LTP (> 24 h) in both strains, whereby effects became evident earlier in CaOlaHsd compared to C57Bl/6 mice. Immunohistochemical analysis revealed a significantly higher expression of D1-receptors in the stratum lacunosum moleculare of CaOlaHsd compared to C57Bl/6 mice. These findings highlight differences in D1/D5R- dependent regulation of strain-dependent variations in hippocampal LTP in C57Bl/6 and CaOlaHsd mice, that may be mediated, in part, by differences in the expression of D1R in the hippocampus.

The effect of normalisation and error model choice on the distribution of the maximum likelihood estimator for a biochemical reaction

Sparse and noisy measurements make parameter estimation for biochemical reaction networks difficult and might lead to ill-posed optimisation problems. This is potentiated if the data has to be normalised, and only fold changes rather than absolute amounts are available. Here, the authors consider the propagation of measurement noise to the distribution of the maximum likelihood (ML) estimator in an in silico study. Therefore, a model of a reversible reaction is considered, for which reaction rate constants using fold changes is estimated. Noise propagation is analysed for different normalisation strategies and different error models. In particular, accuracy, precision, and asymptotic properties of the ML estimator is investigated. Results show that normalisation by the mean of a time series outperforms normalisation by a single time point in the example provided by the authors. Moreover, the error model with a heavy-tail distribution is slightly more robust to large measurement noise, but, beyond this, the choice of the error model did not have a significant impact on the estimation results provided by the authors.

BlotIt—Optimal alignment of Western blot and qPCR experiments

Biological systems are frequently analyzed by means of mechanistic mathematical models. In order to infer model parameters and provide a useful model that can be employed for systems understanding and hypothesis testing, the model is often calibrated on quantitative, time-resolved data. To do so, it is typically important to compare experimental measurements over broad time ranges and various experimental conditions, e.g. perturbations of the biological system. However, most of the established experimental techniques such as Western blot, or quantitative real-time polymerase chain reaction only provide measurements on a relative scale, since different sample volumes, experimental adjustments or varying development times of a gel lead to systematic shifts in the data. In turn, the number of measurements corresponding to the same scale enabling comparability is limited. Here, we present a new flexible method to align measurement data that obeys different scaling factors and compare it to existing normalization approaches. We propose an alignment model to estimate these scaling factors and provide the possibility to adapt this model depending on the measurement technique of interest. In addition, an error model can be specified to adequately weight the different data points and obtain scaling-model based confidence intervals of the finally scaled data points. Our approach is applicable to all sorts of relative measurements and does not need a particular experimental condition that has been measured over all available scales. An implementation of the method is provided with the R package blotIt including refined ways of visualization.

Hepatectomy-Induced Alterations in Hepatic Perfusion and Function - Toward Multi-Scale Computational Modeling for a Better Prediction of Post-hepatectomy Liver Function

Liver resection causes marked perfusion alterations in the liver remnant both on the organ scale (vascular anatomy) and on the microscale (sinusoidal blood flow on tissue level). These changes in perfusion affect hepatic functions via direct alterations in blood supply and drainage, followed by indirect changes of biomechanical tissue properties and cellular function. Changes in blood flow impose compression, tension and shear forces on the liver tissue. These forces are perceived by mechanosensors on parenchymal and non-parenchymal cells of the liver and regulate cell-cell and cell-matrix interactions as well as cellular signaling and metabolism. These interactions are key players in tissue growth and remodeling, a prerequisite to restore tissue function after PHx. Their dysregulation is associated with metabolic impairment of the liver eventually leading to liver failure, a serious post-hepatectomy complication with high morbidity and mortality. Though certain links are known, the overall functional change after liver surgery is not understood due to complex feedback loops, non-linearities, spatial heterogeneities and different time-scales of events. Computational modeling is a unique approach to gain a better understanding of complex biomedical systems. This approach allows (i) integration of heterogeneous data and knowledge on multiple scales into a consistent view of how perfusion is related to hepatic function; (ii) testing and generating hypotheses based on predictive models, which must be validated experimentally and clinically. In the long term, computational modeling will (iii) support surgical planning by predicting surgery-induced perfusion perturbations and their functional (metabolic) consequences; and thereby (iv) allow minimizing surgical risks for the individual patient. Here, we review the alterations of hepatic perfusion, biomechanical properties and function associated with hepatectomy. Specifically, we provide an overview over the clinical problem, preoperative diagnostics, functional imaging approaches, experimental approaches in animal models, mechanoperception in the liver and impact on cellular metabolism, omics approaches with a focus on transcriptomics, data integration and uncertainty analysis, and computational modeling on multiple scales. Finally, we provide a perspective on how multi-scale computational models, which couple perfusion changes to hepatic function, could become part of clinical workflows to predict and optimize patient outcome after complex liver surgery.

Preferential frequency-dependent induction of synaptic depression by the lateral perforant path and of synaptic potentiation by the medial perforant path inputs to the dentate gyrus

The encoding of spatial representations is enabled by synaptic plasticity. The entorhinal cortex sends information to the hippocampus via the lateral (LPP) and medial perforant (MPP) paths that transfer egocentric item-related and allocentric spatial information, respectively. To what extent LPP and MPP information-relay results in different homosynaptic synaptic plasticity responses is unclear. We examined the frequency dependency (at 1, 5, 10, 50, 100, 200 Hz) of long-term potentiation (LTP) and long-term depression (LTD) at MPP and LPP synapses in the dentate gyrus (DG) of freely behaving adult rats. We report that whereas the MPP-DG synapses exhibit a predisposition toward the expression of LTP, LPP-DG synapses prefer to express synaptic depression. The divergence of synaptic plasticity responses is most prominent at afferent frequencies of 5, 100, Hz and 200 Hz. Priming with 10 or 50 Hz significantly modified the subsequent plasticity response in a frequency-dependent manner, but failed to change the preferred direction of change in synaptic strength of MPP and LPP synapses. Evaluation of the expression of GluN1, GluN2A, or GluN2B subunits of the NMDA receptor revealed equivalent expression in the outer and middle thirds of the molecular layer where LPP and MPP inputs convene, respectively, thus excluding NMDA receptors as a substrate for the frequency-dependent differences in bidirectional plasticity. These findings demonstrate that the LPP and MPP inputs to the DG enable differentiated and distinct forms of synaptic plasticity in response to the same afferent frequencies. Effects are extremely robust and resilient to metaplastic priming. These properties may support the functional differentiation of allocentric and item information provided to the DG by the MPP and LPP, respectively, that has been proposed by others. We propose that allocentric spatial information, conveyed by the MPP is encoded through hippocampal LTP in a designated synaptic network. This network is refined and optimized to include egocentric contextual information through LTD triggered by LPP inputs.

Development of a hybrid model for a partially known intracellular signaling pathway through correction term estimation and neural network modeling

Developing an accurate first-principle model is an important step in employing systems biology approaches to analyze an intracellular signaling pathway. However, an accurate first-principle model is difficult to be developed since it requires in-depth mechanistic understandings of the signaling pathway. Since underlying mechanisms such as the reaction network structure are not fully understood, significant discrepancy exists between predicted and actual signaling dynamics. Motivated by these considerations, this work proposes a hybrid modeling approach that combines a first-principle model and an artificial neural network (ANN) model so that predictions of the hybrid model surpass those of the original model. First, the proposed approach determines an optimal subset of model states whose dynamics should be corrected by the ANN by examining the correlation between each state and outputs through relative order. Second, an L2-regularized least-squares problem is solved to infer values of the correction terms that are necessary to minimize the discrepancy between the model predictions and available measurements. Third, an ANN is developed to generalize relationships between the values of the correction terms and the system dynamics. Lastly, the original first-principle model is coupled with the developed ANN to finalize the hybrid model development so that the model will possess generalized prediction capabilities while retaining the model interpretability. We have successfully validated the proposed methodology with two case studies, simplified apoptosis and lipopolysaccharide-induced NFκB signaling pathways, to develop hybrid models with in silico and in vitro measurements, respectively.

An intracellular signaling pathway is often represented by a set of nonlinear ordinary differential equations, which translate our current knowledge about the signaling pathway into a testable mathematical model. However, predictions from such models are often subject to high uncertainty since many signaling pathways are only partially known beforehand. In this study, we propose a systematic approach to develop a hybrid model to improve model accuracy by combining machine learning and the first-principle modeling. Specifically, model correction terms are learned from discrepancy between model predictions and measurements, and these terms are added to the first-principle model to enhance the prediction accuracy. Once these correction terms are learned from the data, an artificial neural network (ANN) model is developed to find an empirical relation between the model and the correction terms so that the developed ANN can be used to posses improved predictive capabilities even in new operating conditions (i.e., generalizability). The final hybrid model is then constructed by coupling the first-principle model with the developed ANN.

Gradient of Expression of Dopamine D2 Receptors Along the Dorso-Ventral Axis of the Hippocampus

Dopamine D2-like receptors (D2R) play an important role in the regulation of hippocampal neuronal excitability and contribute to the regulation of synaptic plasticity, the encoding of hippocampus-dependent memories and the regulation of affective state. In line with this, D2R are targeted in the treatment of psychosis and affective disorders. It has been proposed that the dorso-ventral axis of the hippocampus can be functionally delineated into the dorsal pole that predominantly processes spatial information and the ventral pole that mainly addresses hippocampal processing of emotional and affective state. Although dopaminergic control of hippocampal information processing has been the focus of a multitude of studies, very little is known about the precise distribution of D2R both within anatomically defined sublayers of the hippocampus and along its dorsoventral axis, that could in turn yield insights as to the functional significance of this receptor in supporting hippocampal processing of spatial and affective information. Here, we used an immunohistochemical approach to precisely scrutinize the protein expression of D2R both within the cellular and dendritic layers of the hippocampal subfields, and along the dorso-ventral hippocampal axis. In general, we detected significantly higher levels of protein expression of D2R in the ventral, compared to the dorsal poles with regard to the CA1, CA2, CA3 and dentate gyrus (DG) regions. Effects were very consistent: the molecular layer, granule cell layer and polymorphic layer of the DG exhibited higher D2R levels in the ventral compared to dorsal hippocampus. D2R levels were also significantly higher in the ventral Stratum oriens, Stratum radiatum, and Stratum lacunosum-moleculare layers of the CA1 and CA3 regions. The apical dendrites of the ventral CA2 region also exhibited higher D2R expression compared to the dorsal pole. Taken together, our study suggests that the higher D2R expression levels of the ventral hippocampus may contribute to reported gradients in the degree of expression of synaptic plasticity along the dorso-ventral hippocampal axis, and may support behavioral information processing by the ventral hippocampus.

Extensions of ℓ1 regularization increase detection specificity for cell-type specific parameters in dynamic models

BACKGROUND

Ordinary differential equation systems are frequently utilized to model biological systems and to infer knowledge about underlying properties. For instance, the development of drugs requires the knowledge to which extent malign cells differ from healthy ones to provide a specific treatment with least side effects. As these cell-type specific properties may stem from any part of biochemical cell processes, systematic quantitative approaches are necessary to identify the relevant potential drug targets. An ℓ₁ regularization for the maximum likelihood parameter estimation proved to be successful, but falsely predicted cell-type dependent behaviour had to be corrected manually by using a Profile Likelihood approach.

RESULTS

The choice of extended ℓ₁ penalty functions significantly decreased the number of falsely detected cell-type specific parameters. Thus, the total accuracy of the prediction could be increased. This was tested on a realistic dynamical benchmark model used for the DREAM6 challenge. Among Elastic Net, Adaptive Lasso and a non-convex ℓ_q penalty, the latter one showed the best predictions whilst also requiring least computation time. All extended methods include a hyper-parameter in the regularization function. For an Erythropoietin (EPO) induced signalling pathway, the extended methods ℓ_q and Adaptive Lasso revealed an unpublished alternative parsimonious model when varying the respective hyper-parameters.

CONCLUSIONS

Using ℓ_q or Adaptive Lasso with an a-priori choice for the hyper-parameter can lead to a more specific and accurate result than ℓ₁. Scanning different hyper-parameters can yield additional pieces of information about the system.

In the Piriform Cortex, the Primary Impetus for Information Encoding through Synaptic Plasticity Is Provided by Descending Rather than Ascending Olfactory Inputs

Information encoding by means of persistent changes in synaptic strength supports long-term information storage and memory in structures such as the hippocampus. In the piriform cortex (PC), that engages in the processing of associative memory, only short-term synaptic plasticity has been described to date, both in vitro and in anesthetized rodents in vivo. Whether the PC maintains changes in synaptic strength for longer periods of time is unknown: Such a property would indicate that it can serve as a repository for long-term memories. Here, we report that in freely behaving animals, frequency-dependent synaptic plasticity does not occur in the anterior PC (aPC) following patterned stimulation of the olfactory bulb (OB). Naris closure changed action potential properties of aPC neurons and enabled expression of long-term potentiation (LTP) by OB stimulation, indicating that an intrinsic ability to express synaptic plasticity is present. Odor discrimination and categorization in the aPC is supported by descending inputs from the orbitofrontal cortex (OFC). Here, OFC stimulation resulted in LTP (>4 h), suggesting that this structure plays an important role in promoting information encoding through synaptic plasticity in the aPC. These persistent changes in synaptic strength are likely to comprise a means through which long-term memories are encoded and/or retained in the PC.

Distinct Time-Course of Alterations of Groups I and II Metabotropic Glutamate Receptor and GABAergic Receptor Expression Along the Dorsoventral Hippocampal Axis in an Animal Model of Psychosis

Psychosis is a clinical state that encompasses a range of abnormal conditions, including distortions in sensory information processing and the resultant delusional thinking, emotional discordance and cognitive impairments. Upon developing this condition, the rate at which cognitive and behavioral deteriorations progress steadily increases suggesting an active contribution of the first psychotic event to the progression of structural and functional abnormalities and disease establishment in diagnosed patients. Changes in GABAergic and glutamatergic function, or expression, in the hippocampus have been proposed as a key factor in the pathophysiology of psychosis. However, little is known as to the time-point of onset of putative changes, to what extent they are progressive, and their relation to disease stabilization. Here, we characterized the expression and distribution patterns of groups I and II metabotropic glutamate (mGlu) receptors and GABA receptors 1 week and 3 months after systemic treatment with an N-methyl-D-aspartate receptor (NMDAR) antagonist (MK801) that is used to model a psychosis-like state in adult rats. We found an early alteration in the expression of mGlu1, mGlu2/3, and GABA_B receptors across the hippocampal dorsoventral and transverse axes. This expanded to include an up-regulation of mGlu5 levels across the entire CA1 region and a reduction in GABA_B expression, as well as GAD67-positive interneurons particularly in the dorsal hippocampus that appeared 3 months after treatment. Our findings indicate that a reduction of excitability may occur in the hippocampus soon after first-episode psychosis. This changes, over time, into increased excitability. These hippocampus-specific alterations are likely to contribute to the pathophysiology and stabilization of psychosis.

Analysis of Cell Lineage Trees by Exact Bayesian Inference Identifies Negative Autoregulation of Nanog in Mouse Embryonic Stem Cells

Many cellular effectors of pluripotency are dynamically regulated. In principle, regulatory mechanisms can be inferred from single-cell observations of effector activity across time. However, rigorous inference techniques suitable for noisy, incomplete, and heterogeneous data are lacking. Here, we introduce stochastic inference on lineage trees (STILT), an algorithm capable of identifying stochastic models that accurately describe the quantitative behavior of cell fate markers observed using time-lapse microscopy data collected from proliferating cell populations. STILT performs exact Bayesian parameter inference and stochastic model selection using a particle-filter-based algorithm. We use STILT to investigate the autoregulation of Nanog, a heterogeneously expressed core pluripotency factor, in mouse embryonic stem cells. STILT rejects the possibility of positive Nanog autoregulation with high confidence; instead, model predictions indicate weak negative feedback. We use STILT for rational experimental design and validate model predictions using novel experimental data. STILT is available for download as an open source framework from http://www.imsb.ethz.ch/research/claassen/Software/stilt---stochastic-inference-on-lineage-trees.html.

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.

Impact of measurement noise, experimental design, and estimation methods on Modular Response Analysis based network reconstruction

Modular Response Analysis (MRA) is a method to reconstruct signalling networks from steady-state perturbation data which has frequently been used in different settings. Since these data are usually noisy due to multi-step measurement procedures and biological variability, it is important to investigate the effect of this noise onto network reconstruction. Here we present a systematic study to investigate propagation of noise from concentration measurements to network structures. Therefore, we design an in silico study of the MAPK and the p53 signalling pathways with realistic noise settings. We make use of statistical concepts and measures to evaluate accuracy and precision of individual inferred interactions and resulting network structures. Our results allow to derive clear recommendations to optimize the performance of MRA based network reconstruction: First, large perturbations are favorable in terms of accuracy even for models with non-linear steady-state response curves. Second, a single control measurement for different perturbation experiments seems to be sufficient for network reconstruction, and third, we recommend to execute the MRA workflow with the mean of different replicates for concentration measurements rather than using computationally more involved regression strategies.

Time-Dependent Alterations in the Expression of NMDA Receptor Subunits along the Dorsoventral Hippocampal Axis in an Animal Model of Nascent Psychosis

Psychosis is a mental condition that is characterized by hallucinations, delusions, disordered thought, as well as socio-emotional and cognitive impairments. Once developed, it tends to progress into a chronic psychotic illness. Here, the duration of untreated psychosis plays a crucial role: the earlier the treatment begins, relative to the first episode of the disease, the better the patient's functional prognosis. To what extent the success of early interventions relate to progressive changes at the neurotransmitter receptor level is as yet unclear. In fact, very little is known as to how molecular changes develop, transform, and become established following the first psychotic event. One neurotransmitter receptor for which a specific role in psychosis has been discussed is the N-methyl-d-aspartate receptor (NMDAR). This receptor is especially important for information encoding in the hippocampus. The hippocampus is one of the loci of functional change in psychosis, to which a role in the pathophysiology of psychosis has been ascribed. Here, we examined whether changes in NMDAR subunit expression occur along the dorsoventral axis of the hippocampus 1 week and 3 months after systemic treatment with an NMDAR antagonist (MK801) that initiates a psychosis-like state in adult rats. We found early (1 week) upregulation of the GluN2B levels in the dorso-intermediate hippocampus and late (3 month) downregulation of GluN2A expression across the entire CA1 region. The ventral hippocampus did not exhibit subunit expression changes. These data suggest that a differing vulnerability of the hippocampal longitudinal axis may occur in response to MK801-treatment and provide a time-resolved view of the putative development of pathological changes of NMDAR subunit expression in the hippocampus that initiate with an emulated first episode and progress through to the chronic stabilization of a psychosis-like state in rodents.

Unraveling mitotic protein networks by 3D multiplexed epitope drug screening

Three-dimensional protein localization intricately determines the functional coordination of cellular processes. The complex spatial context of protein landscape has been assessed by multiplexed immunofluorescent staining or mass spectrometry, applied to 2D cell culture with limited physiological relevance or tissue sections. Here, we present 3D SPECS, an automated technology for 3D Spatial characterization of Protein Expression Changes by microscopic Screening. This workflow comprises iterative antibody staining, high-content 3D imaging, and machine learning for detection of mitoses. This is followed by mapping of spatial protein localization into a spherical, cellular coordinate system, a basis for model-based prediction of spatially resolved affinities of proteins. As a proof-of-concept, we mapped twelve epitopes in 3D-cultured spheroids and investigated the network effects of twelve mitotic cancer drugs. Our approach reveals novel insights into spindle fragility and chromatin stress, and predicts unknown interactions between proteins in specific mitotic pathways. 3D SPECS's ability to map potential drug targets by multiplexed immunofluorescence in 3D cell culture combined with our automated high-content assay will inspire future functional protein expression and drug assays.

Encoding Growth Factor Identity in the Temporal Dynamics of FOXO3 under the Combinatorial Control of ERK and AKT Kinases

Extracellular growth factors signal to transcription factors via a limited number of cytoplasmic kinase cascades. It remains unclear how such cascades encode ligand identities and concentrations. In this paper, we use live-cell imaging and statistical modeling to study FOXO3, a transcription factor regulating diverse aspects of cellular physiology that is under combinatorial control. We show that FOXO3 nuclear-to-cytosolic translocation has two temporally distinct phases varying in magnitude with growth factor identity and cell type. These phases comprise synchronous translocation soon after ligand addition followed by an extended back-and-forth shuttling; this shuttling is pulsatile and does not have a characteristic frequency, unlike a simple oscillator. Early and late dynamics are differentially regulated by Akt and ERK and have low mutual information, potentially allowing the two phases to encode different information. In cancer cells in which ERK and Akt are dysregulated by oncogenic mutation, the diversity of states is lower.

Less means more: The magnitude of synaptic plasticity along the hippocampal dorso-ventral axis is inversely related to the expression levels of plasticity-related neurotransmitter receptors

The dorsoventral axis of the hippocampus exhibits functional differentiations with regard to (spatial Vs emotional) learning and information retention (rapid encoding Vs long-term storage), as well as its sensitivity to neuromodulation and information received from extrahippocampal structures. The mechanisms that underlie these differentiations remain unclear. Here, we explored neurotransmitter receptor expression along the dorsoventral hippocampal axis and compared hippocampal synaptic plasticity in the CA1 region of the dorsal (DH), intermediate (IH) and ventral hippocampi (VH). We observed a very distinct gradient of expression of the N-methyl-D-aspartate receptor GluN2B subunit in the Stratum radiatum (DH< IH< VH). A similar distribution gradient (DH< IH< VH) was evident in the hippocampus for GluN1, the metabotropic glutamate receptors mGlu1 and mGlu2/3, GABAB and the dopamine-D1 receptor. GABAA exhibited the opposite expression relationship (DH > IH > VH). Neurotransmitter release probability was lowest in DH. Surprisingly, identical afferent stimulation conditions resulted in hippocampal synaptic plasticity that was the most robust in the DH, compared with IH and VH. These data suggest that differences in hippocampal information processing and synaptic plasticity along the dorsoventral axis may relate to specific differences in the expression of plasticity-related neurotransmitter receptors. This gradient may support the fine-tuning and specificity of hippocampal synaptic encoding.

Robust parameter estimation for dynamical systems from outlier-corrupted data

Motivation

Dynamics of cellular processes are often studied using mechanistic mathematical models. These models possess unknown parameters which are generally estimated from experimental data assuming normally distributed measurement noise. Outlier corruption of datasets often cannot be avoided. These outliers may distort the parameter estimates, resulting in incorrect model predictions. Robust parameter estimation methods are required which provide reliable parameter estimates in the presence of outliers.

Results

In this manuscript, we propose and evaluate methods for estimating the parameters of ordinary differential equation models from outlier-corrupted data. As alternatives to the normal distribution as noise distribution, we consider the Laplace, the Huber, the Cauchy and the Student's t distribution. We assess accuracy, robustness and computational efficiency of estimators using these different distribution assumptions. To this end, we consider artificial data of a conversion process, as well as published experimental data for Epo-induced JAK/STAT signaling. We study how well the methods can compensate and discover artificially introduced outliers. Our evaluation reveals that using alternative distributions improves the robustness of parameter estimates.

Availability and Implementation

The MATLAB implementation of the likelihood functions using the distribution assumptions is available at Bioinformatics online.

Contact

jan.hasenauer@helmholtz-muenchen.de.

Supplementary information

Supplementary material are available at Bioinformatics online.

Mathematical model of early Reelin-induced Src family kinase-mediated signaling

Reelin is a large glycoprotein with a dual role in the mammalian brain. It regulates the positioning and differentiation of postmitotic neurons during brain development and modulates neurotransmission and memory formation in the adult brain. Alterations in the Reelin signaling pathway have been described in different psychiatric disorders. Reelin mainly signals by binding to the lipoprotein receptors Vldlr and ApoER2, which induces tyrosine phosphorylation of the adaptor protein Dab1 mediated by Src family kinases (SFKs). In turn, phosphorylated Dab1 activates downstream signaling cascades, including PI3-kinase-dependent signaling. In this work, a mechanistic model based on ordinary differential equations was built to model early dynamics of the Reelin-mediated signaling cascade. Mechanistic models are frequently used to disentangle the highly complex mechanisms underlying cellular processes and obtain new biological insights. The model was calibrated on time-resolved data and a dose-response measurement of protein concentrations measured in cortical neurons treated with Reelin. It focusses on the interplay between Dab1 and SFKs with a special emphasis on the tyrosine phosphorylation of Dab1, and their role for the regulation of Reelin-induced signaling. Model selection was performed on different model structures and a comprehensive mechanistic model of the early Reelin signaling cascade is provided in this work. It emphasizes the importance of Reelin-induced lipoprotein receptor clustering for SFK-mediated Dab1 trans-phosphorylation and does not require co-receptors to describe the measured data. The model is freely available within the open-source framework Data2Dynamics (www.data2dynamics.org). It can be used to generate predictions that can be validated experimentally, and provides a platform for model extensions both to downstream targets such as transcription factors and interactions with other transmembrane proteins and neuronal signaling pathways.

Correlated receptor transport processes buffer single-cell heterogeneity

Cells typically vary in their response to extracellular ligands. Receptor transport processes modulate ligand-receptor induced signal transduction and impact the variability in cellular responses. Here, we quantitatively characterized cellular variability in erythropoietin receptor (EpoR) trafficking at the single-cell level based on live-cell imaging and mathematical modeling. Using ensembles of single-cell mathematical models reduced parameter uncertainties and showed that rapid EpoR turnover, transport of internalized EpoR back to the plasma membrane, and degradation of Epo-EpoR complexes were essential for receptor trafficking. EpoR trafficking dynamics in adherent H838 lung cancer cells closely resembled the dynamics previously characterized by mathematical modeling in suspension cells, indicating that dynamic properties of the EpoR system are widely conserved. Receptor transport processes differed by one order of magnitude between individual cells. However, the concentration of activated Epo-EpoR complexes was less variable due to the correlated kinetics of opposing transport processes acting as a buffering system.

Cell surface receptors translate extracellular ligand concentrations to intracellular responses. Receptor transport between the plasma membrane and other cellular compartments regulates the number of accessible receptors at the plasma membrane that determines the strength of downstream pathway activation at a given ligand concentration. In cell populations, pathway activation strength and cellular responses vary between cells. Understanding origins of cell-to-cell variability is highly relevant for cancer research, motivated by the problem of fractional killing by chemotherapies and development of resistance in subpopulations of tumor cells. The erythropoietin receptor (EpoR) is a characteristic example of a receptor system that strongly depends on receptor transport processes. It is involved in several cellular processes, such as differentiation or proliferation, regulates the renewal of erythrocytes, and is expressed in several tumors. To investigate the involvement of receptor transport processes in cell-to-cell variability, we quantitatively characterized trafficking of EpoR in individual cells by combining live-cell imaging with mathematical modeling. Thereby, we found that EpoR dynamics was strongly dependent on rapid receptor transport and turnover. Interestingly, although transport processes largely differed between individual cells, receptor concentrations in cellular compartments were robust to variability in trafficking processes due to the correlated kinetics of opposing transport processes.

TH1 signatures are present in the lower airways of children with severe asthma, regardless of allergic status

BACKGROUND

The pathogenesis of severe asthma in childhood remains poorly understood.

OBJECTIVE

We sought to construct the immunologic landscape in the airways of children with severe asthma.

METHODS

Comprehensive analysis of multiple cell types and mediators was performed by using flow cytometry and a multiplex assay with bronchoalveolar lavage (BAL) specimens (n = 68) from 52 highly characterized allergic and nonallergic children (0.5-17 years) with severe treatment-refractory asthma. Multiple relationships were tested by using linear mixed-effects modeling.

RESULTS

Memory CCR5⁺ T_H1 cells were enriched in BAL fluid versus blood, and pathogenic respiratory viruses and bacteria were readily detected. IFN-γ⁺IL-17⁺ and IFN-γ^-IL-17⁺ subsets constituted secondary T_H types, and BAL fluid CD8⁺ T cells were almost exclusively IFN-γ⁺. The T_H17-associated mediators IL-23 and macrophage inflammatory protein 3α/CCL20 were highly expressed. Despite low T_H2 numbers, T_H2 cytokines were detected, and T_H2 skewing correlated with total IgE levels. Type 2 innate lymphoid cells and basophils were scarce in BAL fluid. Levels of IL-5, IL-33, and IL-28A/IFN-λ2 were increased in multisensitized children and correlated with IgE levels to dust mite, ryegrass, and fungi but not cat, ragweed, or food sources. Additionally, levels of IL-5, but no other cytokine, increased with age and correlated with eosinophil numbers in BAL fluid and blood. Both plasmacytoid and IgE⁺FcεRI⁺ myeloid dendritic cells were present in BAL fluid.

CONCLUSIONS

The lower airways of children with severe asthma display a dominant T_H1 signature and atypical cytokine profiles that link to allergic status. Our findings deviate from established paradigms and warrant further assessment of the pathogenicity of T_H1 cells in patients with severe asthma.

IL-1β-induced and p38MAPK-dependent activation of the mitogen-activated protein kinase-activated protein kinase 2 (MK2) in hepatocytes: Signal transduction with robust and concentration-independent signal amplification

The IL-1β induced activation of the p38^MAPK/MAPK-activated protein kinase 2 (MK2) pathway in hepatocytes is important for control of the acute phase response and regulation of liver regeneration. Many aspects of the regulatory relevance of this pathway have been investigated in immune cells in the context of inflammation. However, very little is known about concentration-dependent activation kinetics and signal propagation in hepatocytes and the role of MK2. We established a mathematical model for IL-1β-induced activation of the p38^MAPK/MK2 pathway in hepatocytes that was calibrated to quantitative data on time- and IL-1β concentration-dependent phosphorylation of p38^MAPK and MK2 in primary mouse hepatocytes. This analysis showed that, in hepatocytes, signal transduction from IL-1β via p38^MAPK to MK2 is characterized by strong signal amplification. Quantification of p38^MAPK and MK2 revealed that, in hepatocytes, at maximum, 11.3% of p38^MAPK molecules and 36.5% of MK2 molecules are activated in response to IL-1β. The mathematical model was experimentally validated by employing phosphatase inhibitors and the p38^MAPK inhibitor SB203580. Model simulations predicted an IC₅₀ of 1-1.2 μm for SB203580 in hepatocytes. In silico analyses and experimental validation demonstrated that the kinase activity of p38^MAPK determines signal amplitude, whereas phosphatase activity affects both signal amplitude and duration. p38^MAPK and MK2 concentrations and responsiveness toward IL-1β were quantitatively compared between hepatocytes and macrophages. In macrophages, the absolute p38^MAPK and MK2 concentration was significantly higher. Finally, in line with experimental observations, the mathematical model predicted a significantly higher half-maximal effective concentration for IL-1β-induced pathway activation in macrophages compared with hepatocytes, underscoring the importance of cell type-specific differences in pathway regulation.

Driving the Model to Its Limit: Profile Likelihood Based Model Reduction

In systems biology, one of the major tasks is to tailor model complexity to information content of the data. A useful model should describe the data and produce well-determined parameter estimates and predictions. Too small of a model will not be able to describe the data whereas a model which is too large tends to overfit measurement errors and does not provide precise predictions. Typically, the model is modified and tuned to fit the data, which often results in an oversized model. To restore the balance between model complexity and available measurements, either new data has to be gathered or the model has to be reduced. In this manuscript, we present a data-based method for reducing non-linear models. The profile likelihood is utilised to assess parameter identifiability and designate likely candidates for reduction. Parameter dependencies are analysed along profiles, providing context-dependent suggestions for the type of reduction. We discriminate four distinct scenarios, each associated with a specific model reduction strategy. Iterating the presented procedure eventually results in an identifiable model, which is capable of generating precise and testable predictions. Source code for all toy examples is provided within the freely available, open-source modelling environment Data2Dynamics based on MATLAB available at http://www.data2dynamics.org/, as well as the R packages dMod/cOde available at https://github.com/dkaschek/. Moreover, the concept is generally applicable and can readily be used with any software capable of calculating the profile likelihood.

Tailored parameter optimization methods for ordinary differential equation models with steady-state constraints

BACKGROUND

Ordinary differential equation (ODE) models are widely used to describe (bio-)chemical and biological processes. To enhance the predictive power of these models, their unknown parameters are estimated from experimental data. These experimental data are mostly collected in perturbation experiments, in which the processes are pushed out of steady state by applying a stimulus. The information that the initial condition is a steady state of the unperturbed process provides valuable information, as it restricts the dynamics of the process and thereby the parameters. However, implementing steady-state constraints in the optimization often results in convergence problems.

RESULTS

In this manuscript, we propose two new methods for solving optimization problems with steady-state constraints. The first method exploits ideas from optimization algorithms on manifolds and introduces a retraction operator, essentially reducing the dimension of the optimization problem. The second method is based on the continuous analogue of the optimization problem. This continuous analogue is an ODE whose equilibrium points are the optima of the constrained optimization problem. This equivalence enables the use of adaptive numerical methods for solving optimization problems with steady-state constraints. Both methods are tailored to the problem structure and exploit the local geometry of the steady-state manifold and its stability properties. A parameterization of the steady-state manifold is not required. The efficiency and reliability of the proposed methods is evaluated using one toy example and two applications. The first application example uses published data while the second uses a novel dataset for Raf/MEK/ERK signaling. The proposed methods demonstrated better convergence properties than state-of-the-art methods employed in systems and computational biology. Furthermore, the average computation time per converged start is significantly lower. In addition to the theoretical results, the analysis of the dataset for Raf/MEK/ERK signaling provides novel biological insights regarding the existence of feedback regulation.

CONCLUSION

Many optimization problems considered in systems and computational biology are subject to steady-state constraints. While most optimization methods have convergence problems if these steady-state constraints are highly nonlinear, the methods presented recover the convergence properties of optimizers which can exploit an analytical expression for the parameter-dependent steady state. This renders them an excellent alternative to methods which are currently employed in systems and computational biology.

Identification of Cell Type-Specific Differences in Erythropoietin Receptor Signaling in Primary Erythroid and Lung Cancer Cells

Lung cancer, with its most prevalent form non-small-cell lung carcinoma (NSCLC), is one of the leading causes of cancer-related deaths worldwide, and is commonly treated with chemotherapeutic drugs such as cisplatin. Lung cancer patients frequently suffer from chemotherapy-induced anemia, which can be treated with erythropoietin (EPO). However, studies have indicated that EPO not only promotes erythropoiesis in hematopoietic cells, but may also enhance survival of NSCLC cells. Here, we verified that the NSCLC cell line H838 expresses functional erythropoietin receptors (EPOR) and that treatment with EPO reduces cisplatin-induced apoptosis. To pinpoint differences in EPO-induced survival signaling in erythroid progenitor cells (CFU-E, colony forming unit-erythroid) and H838 cells, we combined mathematical modeling with a method for feature selection, the L₁ regularization. Utilizing an example model and simulated data, we demonstrated that this approach enables the accurate identification and quantification of cell type-specific parameters. We applied our strategy to quantitative time-resolved data of EPO-induced JAK/STAT signaling generated by quantitative immunoblotting, mass spectrometry and quantitative real-time PCR (qRT-PCR) in CFU-E and H838 cells as well as H838 cells overexpressing human EPOR (H838-HA-hEPOR). The established parsimonious mathematical model was able to simultaneously describe the data sets of CFU-E, H838 and H838-HA-hEPOR cells. Seven cell type-specific parameters were identified that included for example parameters for nuclear translocation of STAT5 and target gene induction. Cell type-specific differences in target gene induction were experimentally validated by qRT-PCR experiments. The systematic identification of pathway differences and sensitivities of EPOR signaling in CFU-E and H838 cells revealed potential targets for intervention to selectively inhibit EPO-induced signaling in the tumor cells but leave the responses in erythroid progenitor cells unaffected. Thus, the proposed modeling strategy can be employed as a general procedure to identify cell type-specific parameters and to recommend treatment strategies for the selective targeting of specific cell types.

A major challenge in the development of therapeutic interventions is the selective inhibition of a signal transduction pathway in one cell type such as a cancer cell leaving the other cell type such as a healthy cell as unaffected as possible. Here, we propose a new approach that combines mathematical modeling based on quantitative experimental data with statistical methods. We demonstrate based on simulated data that our approach can determine which parameters are the same and which parameters differ in two exemplary cell types. We compare a lung cancer cell line to the precursor cells of red blood cells. We show that the same signal transduction network induced by erythropoietin (EPO), a hormone that is frequently employed to treat anemia in cancer patients, regulates survival of both cell types. Based on our experimental data in combination with our computational approach, we identify seven cell type-specific differences in this signaling pathway. Our strategy allows predicting therapeutic targets that could be inhibited to interfere with survival of lung cancer cells while leaving production of red blood cells unaffected.

Noise Induces Biased Estimation of the Correction Gain

The detection of an error in the motor output and the correction in the next movement are critical components of any form of motor learning. Accordingly, a variety of iterative learning models have assumed that a fraction of the error is adjusted in the next trial. This critical fraction, the correction gain, learning rate, or feedback gain, has been frequently estimated via least-square regression of the obtained data set. Such data contain not only the inevitable noise from motor execution, but also noise from measurement. It is generally assumed that this noise averages out with large data sets and does not affect the parameter estimation. This study demonstrates that this is not the case and that in the presence of noise the conventional estimate of the correction gain has a significant bias, even with the simplest model. Furthermore, this bias does not decrease with increasing length of the data set. This study reveals this limitation of current system identification methods and proposes a new method that overcomes this limitation. We derive an analytical form of the bias from a simple regression method (Yule-Walker) and develop an improved identification method. This bias is discussed as one of other examples for how the dynamics of noise can introduce significant distortions in data analysis.

Model calibration and uncertainty analysis in signaling networks

For a long time the biggest challenges in modeling cellular signal transduction networks has been the inference of crucial pathway components and the qualitative description of their interactions. As a result of the emergence of powerful high-throughput experiments, it is now possible to measure data of high temporal and spatial resolution and to analyze signaling dynamics quantitatively. In addition, this increase of high-quality data is the basis for a better understanding of model limitations and their influence on the predictive power of models. We review established approaches in signal transduction network modeling with a focus on ordinary differential equation models as well as related developments in model calibration. As central aspects of the calibration process we discuss possibilities of model adaptation based on data-driven parameter optimization and the concomitant objective of reducing model uncertainties.

Analysis of Reverse Phase Protein Array Data: From Experimental Design towards Targeted Biomarker Discovery

Mastering the systematic analysis of tumor tissues on a large scale has long been a technical challenge for proteomics. In 2001, reverse phase protein arrays (RPPA) were added to the repertoire of existing immunoassays, which, for the first time, allowed a profiling of minute amounts of tumor lysates even after microdissection. A characteristic feature of RPPA is its outstanding sample capacity permitting the analysis of thousands of samples in parallel as a routine task. Until today, the RPPA approach has matured to a robust and highly sensitive high-throughput platform, which is ideally suited for biomarker discovery. Concomitant with technical advancements, new bioinformatic tools were developed for data normalization and data analysis as outlined in detail in this review. Furthermore, biomarker signatures obtained by different RPPA screens were compared with another or with that obtained by other proteomic formats, if possible. Options for overcoming the downside of RPPA, which is the need to steadily validate new antibody batches, will be discussed. Finally, a debate on using RPPA to advance personalized medicine will conclude this article.

Robust and efficient parameter estimation in dynamic models of biological systems

Background

Dynamic modelling provides a systematic framework to understand function in biological systems. Parameter estimation in nonlinear dynamic models remains a very challenging inverse problem due to its nonconvexity and ill-conditioning. Associated issues like overfitting and local solutions are usually not properly addressed in the systems biology literature despite their importance.

Here we present a method for robust and efficient parameter estimation which uses two main strategies to surmount the aforementioned difficulties: (i) efficient global optimization to deal with nonconvexity, and (ii) proper regularization methods to handle ill-conditioning. In the case of regularization, we present a detailed critical comparison of methods and guidelines for properly tuning them. Further, we show how regularized estimations ensure the best trade-offs between bias and variance, reducing overfitting, and allowing the incorporation of prior knowledge in a systematic way.

Results

We illustrate the performance of the presented method with seven case studies of different nature and increasing complexity, considering several scenarios of data availability, measurement noise and prior knowledge. We show how our method ensures improved estimations with faster and more stable convergence. We also show how the calibrated models are more generalizable. Finally, we give a set of simple guidelines to apply this strategy to a wide variety of calibration problems.

Conclusions

Here we provide a parameter estimation strategy which combines efficient global optimization with a regularization scheme. This method is able to calibrate dynamic models in an efficient and robust way, effectively fighting overfitting and allowing the incorporation of prior information.

Electronic supplementary material

The online version of this article (doi:10.1186/s12918-015-0219-2) contains supplementary material, which is available to authorized users.

Studying Cellular Signal Transduction with OMIC Technologies

In the gulf between genotype and phenotype exists proteins and, in particular, protein signal transduction systems. These systems use a relatively limited parts list to respond to a much longer list of extracellular, environmental, and/or mechanical cues with rapidity and specificity. Most signaling networks function in a highly non-linear and often contextual manner. Furthermore, these processes occur dynamically across space and time. Because of these complexities, systems and "OMIC" approaches are essential for the study of signal transduction. One challenge in using OMIC-scale approaches to study signaling is that the "signal" can take different forms in different situations. Signals are encoded in diverse ways such as protein-protein interactions, enzyme activities, localizations, or post-translational modifications to proteins. Furthermore, in some cases, signals may be encoded only in the dynamics, duration, or rates of change of these features. Accordingly, systems-level analyses of signaling may need to integrate multiple experimental and/or computational approaches. As the field has progressed, the non-triviality of integrating experimental and computational analyses has become apparent. Successful use of OMIC methods to study signaling will require the "right" experiments and the "right" modeling approaches, and it is critical to consider both in the design phase of the project. In this review, we discuss common OMIC and modeling approaches for studying signaling, emphasizing the philosophical and practical considerations for effectively merging these two types of approaches to maximize the probability of obtaining reliable and novel insights into signaling biology.

Optimal experiment design for model selection in biochemical networks

Background

Mathematical modeling is often used to formalize hypotheses on how a biochemical network operates by discriminating between competing models. Bayesian model selection offers a way to determine the amount of evidence that data provides to support one model over the other while favoring simple models. In practice, the amount of experimental data is often insufficient to make a clear distinction between competing models. Often one would like to perform a new experiment which would discriminate between competing hypotheses.

Results

We developed a novel method to perform Optimal Experiment Design to predict which experiments would most effectively allow model selection. A Bayesian approach is applied to infer model parameter distributions. These distributions are sampled and used to simulate from multivariate predictive densities. The method is based on a k-Nearest Neighbor estimate of the Jensen Shannon divergence between the multivariate predictive densities of competing models.

Conclusions

We show that the method successfully uses predictive differences to enable model selection by applying it to several test cases. Because the design criterion is based on predictive distributions, which can be computed for a wide range of model quantities, the approach is very flexible. The method reveals specific combinations of experiments which improve discriminability even in cases where data is scarce. The proposed approach can be used in conjunction with existing Bayesian methodologies where (approximate) posteriors have been determined, making use of relations that exist within the inferred posteriors.

Evaluating strategies to normalise biological replicates of Western blot data

Western blot data are widely used in quantitative applications such as statistical testing and mathematical modelling. To ensure accurate quantitation and comparability between experiments, Western blot replicates must be normalised, but it is unclear how the available methods affect statistical properties of the data. Here we evaluate three commonly used normalisation strategies: (i) by fixed normalisation point or control; (ii) by sum of all data points in a replicate; and (iii) by optimal alignment of the replicates. We consider how these different strategies affect the coefficient of variation (CV) and the results of hypothesis testing with the normalised data. Normalisation by fixed point tends to increase the mean CV of normalised data in a manner that naturally depends on the choice of the normalisation point. Thus, in the context of hypothesis testing, normalisation by fixed point reduces false positives and increases false negatives. Analysis of published experimental data shows that choosing normalisation points with low quantified intensities results in a high normalised data CV and should thus be avoided. Normalisation by sum or by optimal alignment redistributes the raw data uncertainty in a mean-dependent manner, reducing the CV of high intensity points and increasing the CV of low intensity points. This causes the effect of normalisations by sum or optimal alignment on hypothesis testing to depend on the mean of the data tested; for high intensity points, false positives are increased and false negatives are decreased, while for low intensity points, false positives are decreased and false negatives are increased. These results will aid users of Western blotting to choose a suitable normalisation strategy and also understand the implications of this normalisation for subsequent hypothesis testing.

Perturbation biology: inferring signaling networks in cellular systems

We present a powerful experimental-computational technology for inferring network models that predict the response of cells to perturbations, and that may be useful in the design of combinatorial therapy against cancer. The experiments are systematic series of perturbations of cancer cell lines by targeted drugs, singly or in combination. The response to perturbation is quantified in terms of relative changes in the measured levels of proteins, phospho-proteins and cellular phenotypes such as viability. Computational network models are derived de novo, i.e., without prior knowledge of signaling pathways, and are based on simple non-linear differential equations. The prohibitively large solution space of all possible network models is explored efficiently using a probabilistic algorithm, Belief Propagation (BP), which is three orders of magnitude faster than standard Monte Carlo methods. Explicit executable models are derived for a set of perturbation experiments in SKMEL-133 melanoma cell lines, which are resistant to the therapeutically important inhibitor of RAF kinase. The resulting network models reproduce and extend known pathway biology. They empower potential discoveries of new molecular interactions and predict efficacious novel drug perturbations, such as the inhibition of PLK1, which is verified experimentally. This technology is suitable for application to larger systems in diverse areas of molecular biology.

Drugs that target specific effects of signaling proteins are promising agents for treating cancer. One of the many obstacles facing optimal drug design is inadequate quantitative understanding of the coordinated interactions between signaling proteins. De novo model inference of network or pathway models refers to the algorithmic construction of mathematical predictive models from experimental data without dependence on prior knowledge. De novo inference is difficult because of the prohibitively large number of possible sets of interactions that may or may not be consistent with observations. Our new method overcomes this difficulty by adapting a method from statistical physics, called Belief Propagation, which first calculates probabilistically the most likely interactions in the vast space of all possible solutions, then derives a set of individual, highly probable solutions in the form of executable models. In this paper, we test this method on artificial data and then apply it to model signaling pathways in a BRAF-mutant melanoma cancer cell line based on a large set of rich output measurements from a systematic set of perturbation experiments using drug combinations. Our results are in agreement with established biological knowledge, predict novel interactions, and predict efficacious drug targets that are specific to the experimental cell line and potentially to related tumors. The method has the potential, with sufficient systematic perturbation data, to model, de novo and quantitatively, the effects of hundreds of proteins on cellular responses, on a scale that is currently unreachable in diverse areas of cell biology. In a disease context, the method is applicable to the computational design of novel combination drug treatments.

Mathematical modeling of biochemical systems with PottersWheel

The program PottersWheel has been developed to provide an intuitive and yet powerful framework for data-based modeling of dynamical systems like biochemical reaction networks. Its key functionality is multi-experiment fitting, where several experimental data sets from different laboratory conditions are fitted simultaneously in order to improve the estimation of unknown model parameters, to check the validity of a given model, and to discriminate competing model hypotheses. New experiments can be designed interactively. Models are either created text-based or using a visual model designer. Dynamically generated and compiled C files provide fast simulation and fitting procedures. Each function can either be accessed using a graphical user interface or via command line, allowing for batch processing within custom Matlab scripts. PottersWheel is designed as a Matlab toolbox, comprises 250,000 lines of Matlab and C code, and is freely available for academic usage at www.potterswheel.de .

New types of experimental data shape the use of enzyme kinetics for dynamic network modeling

Since the publication of Leonor Michaelis and Maude Menten's paper on the reaction kinetics of the enzyme invertase in 1913, molecular biology has evolved tremendously. New measurement techniques allow in vivo characterization of the whole genome, proteome or transcriptome of cells, whereas the classical enzyme essay only allows determination of the two Michaelis-Menten parameters V and K(m). Nevertheless, Michaelis-Menten kinetics are still commonly used, not only in the in vitro context of enzyme characterization but also as a rate law for enzymatic reactions in larger biochemical reaction networks. In this review, we give an overview of the historical development of kinetic rate laws originating from Michaelis-Menten kinetics over the past 100 years. Furthermore, we briefly summarize the experimental techniques used for the characterization of enzymes, and discuss web resources that systematically store kinetic parameters and related information. Finally, describe the novel opportunities that arise from using these data in dynamic mathematical modeling. In this framework, traditional in vitro approaches may be combined with modern genome-scale measurements to foster thorough understanding of the underlying complex mechanisms.

Lessons learned from quantitative dynamical modeling in systems biology

Due to the high complexity of biological data it is difficult to disentangle cellular processes relying only on intuitive interpretation of measurements. A Systems Biology approach that combines quantitative experimental data with dynamic mathematical modeling promises to yield deeper insights into these processes. Nevertheless, with growing complexity and increasing amount of quantitative experimental data, building realistic and reliable mathematical models can become a challenging task: the quality of experimental data has to be assessed objectively, unknown model parameters need to be estimated from the experimental data, and numerical calculations need to be precise and efficient. Here, we discuss, compare and characterize the performance of computational methods throughout the process of quantitative dynamic modeling using two previously established examples, for which quantitative, dose- and time-resolved experimental data are available. In particular, we present an approach that allows to determine the quality of experimental data in an efficient, objective and automated manner. Using this approach data generated by different measurement techniques and even in single replicates can be reliably used for mathematical modeling. For the estimation of unknown model parameters, the performance of different optimization algorithms was compared systematically. Our results show that deterministic derivative-based optimization employing the sensitivity equations in combination with a multi-start strategy based on latin hypercube sampling outperforms the other methods by orders of magnitude in accuracy and speed. Finally, we investigated transformations that yield a more efficient parameterization of the model and therefore lead to a further enhancement in optimization performance. We provide a freely available open source software package that implements the algorithms and examples compared here.