Parameter identifiability analysis and visualization in large-scale kinetic models of biosystems

Mechanism-based and data-driven modeling in cell-free synthetic biology

Cell-free systems have emerged as a versatile platform in synthetic biology, finding applications in various areas such as prototyping synthetic circuits, biosensor development, and biomanufacturing. To streamline the prototyping process, cell-free systems often incorporate a modeling step that predicts the outcomes of various experimental scenarios, providing a deeper insight into the underlying mechanisms and functions. There are two recognized approaches for modeling these systems: mechanism-based modeling, which models the underlying reaction mechanisms; and data-driven modeling, which makes predictions based on data without preconceived interactions between system components. In this highlight, we focus on the latest advancements in both modeling approaches for cell-free systems, exploring their potential for the design and optimization of synthetic genetic circuits.

Bottom-up parameterization of enzyme rate constants: Reconciling inconsistent data

Kinetic models of metabolism are promising platforms for studying complex metabolic systems and designing production strains. Given the availability of enzyme kinetic data from historical experiments and machine learning estimation tools, a straightforward modeling approach is to assemble kinetic data enzyme by enzyme until a desired scale is reached. However, this type of 'bottom up' parameterization of kinetic models has been difficult due to a number of issues including gaps in kinetic parameters, the complexity of enzyme mechanisms, inconsistencies between parameters obtained from different sources, and in vitro-in vivo differences. Here, we present a computational workflow for the robust estimation of kinetic parameters for detailed mass action enzyme models while taking into account parameter uncertainty. The resulting software package, termed MASSef (the Mass Action Stoichiometry Simulation Enzyme Fitting package), can handle standard 'macroscopic' kinetic parameters, including Km, kcat, Ki, Keq, and nh, as well as diverse reaction mechanisms defined in terms of mass action reactions and 'microscopic' rate constants. We provide three enzyme case studies demonstrating that this approach can identify and reconcile inconsistent data either within in vitro experiments or between in vitro and in vivo enzyme function. We further demonstrate how parameterized enzyme modules can be used to assemble pathway-scale kinetic models consistent with in vivo behavior. This work builds on the legacy of knowledge on kinetic behavior of enzymes by enabling robust parameterization of enzyme kinetic models at scale utilizing the abundance of historical literature data and machine learning parameter estimates.

Mathematical modeling of hypoxia and adenosine to explore tumor escape mechanisms in DC-based immunotherapy

Identifying and controlling tumor escape mechanisms is crucial for improving cancer treatment effectiveness. Experimental studies reveal tumor hypoxia and adenosine as significant contributors to such mechanisms. Hypoxia exacerbates adenosine levels in the tumor microenvironment. Combining inhibition of these factors with dendritic cell (DC)-based immunotherapy promises improved clinical outcomes. However, challenges include understanding dynamics, optimal vaccine dosages, and timing. Mathematical models, including agent-based, diffusion, and ordinary differential equations, address these challenges. Here, we employ these models for the first time to elucidate how hypoxia and adenosine facilitate tumor escape in DC-based immunotherapy. After parameter estimation using experimental data, we optimize vaccination protocols to minimize tumor growth. Sensitivity analysis highlights adenosine's significant impact on immunotherapy efficacy. Its suppressive role impedes treatment success, but inhibiting adenosine could enhance therapy, as suggested by the model. Our findings shed light on hypoxia and adenosine-mediated tumor escape mechanisms, informing future treatment strategies. Additionally, identifiability analysis confirms accurate parameter determination using experimental data.

Exploring Emergent Properties in Enzymatic Reaction Networks: Design and Control of Dynamic Functional Systems

The intricate and complex features of enzymatic reaction networks (ERNs) play a key role in the emergence and sustenance of life. Constructing such networks in vitro enables stepwise build up in complexity and introduces the opportunity to control enzymatic activity using physicochemical stimuli. Rational design and modulation of network motifs enable the engineering of artificial systems with emergent functionalities. Such functional systems are useful for a variety of reasons such as creating new-to-nature dynamic materials, producing value-added chemicals, constructing metabolic modules for synthetic cells, and even enabling molecular computation. In this review, we offer insights into the chemical characteristics of ERNs while also delving into their potential applications and associated challenges.

Endothelial cells signaling and patterning under hypoxia: a mechanistic integrative computational model including the Notch-Dll4 pathway

Introduction: Several signaling pathways are activated during hypoxia to promote angiogenesis, leading to endothelial cell patterning, interaction, and downstream signaling. Understanding the mechanistic signaling differences between endothelial cells under normoxia and hypoxia and their response to different stimuli can guide therapies to modulate angiogenesis. We present a novel mechanistic model of interacting endothelial cells, including the main pathways involved in angiogenesis. Methods: We calibrate and fit the model parameters based on well-established modeling techniques that include structural and practical parameter identifiability, uncertainty quantification, and global sensitivity. Results: Our results indicate that the main pathways involved in patterning tip and stalk endothelial cells under hypoxia differ, and the time under hypoxia interferes with how different stimuli affect patterning. Additionally, our simulations indicate that Notch signaling might regulate vascular permeability and establish different Nitric Oxide release patterns for tip/stalk cells. Following simulations with various stimuli, our model suggests that factors such as time under hypoxia and oxygen availability must be considered for EC pattern control. Discussion: This project provides insights into the signaling and patterning of endothelial cells under various oxygen levels and stimulation by VEGFA and is our first integrative approach toward achieving EC control as a method for improving angiogenesis. Overall, our model provides a computational framework that can be built on to test angiogenesis-related therapies by modulation of different pathways, such as the Notch pathway.

Iterative design of training data to control intricate enzymatic reaction networks

Kinetic modeling of in vitro enzymatic reaction networks is vital to understand and control the complex behaviors emerging from the nonlinear interactions inside. However, modeling is severely hampered by the lack of training data. Here, we introduce a methodology that combines an active learning-like approach and flow chemistry to efficiently create optimized datasets for a highly interconnected enzymatic reactions network with multiple sub-pathways. The optimal experimental design (OED) algorithm designs a sequence of out-of-equilibrium perturbations to maximize the information about the reaction kinetics, yielding a descriptive model that allows control of the output of the network towards any cost function. We experimentally validate the model by forcing the network to produce different product ratios while maintaining a minimum level of overall conversion efficiency. Our workflow scales with the complexity of the system and enables the optimization of previously unobtainable network outputs.

Modeling cancer progression: an integrated workflow extending data-driven kinetic models to bio-mechanical PDE models

Computational modeling of cancer can help unveil dynamics and interactions that are hard to replicate experimentally. Thanks to the advancement in cancer databases and data analysis technologies, these models have become more robust than ever. There are many mathematical models which investigate cancer through different approaches, from sub-cellular to tissue scale, and from treatment to diagnostic points of view. In this study, we lay out a step-by-step methodology for a data-driven mechanistic model of the tumor microenvironment. We discuss data acquisition strategies, data preparation, parameter estimation, and sensitivity analysis techniques. Furthermore, we propose a possible approach to extend mechanistic ordinary differential equation models to PDE models coupled with mechanical growth. The workflow discussed in this article can help understand the complex temporal and spatial interactions between cells and cytokines in the tumor microenvironment and their effect on tumor growth.

Practical Identifiability of Plant Growth Models: A Unifying Framework and Its Specification for Three Local Indices

Amid the rise of machine learning models, a substantial portion of plant growth models remains mechanistic, seeking to capture an in-depth understanding of the underlying phenomena governing the system's dynamics. The development of these models typically involves parameter estimation from experimental data. Ensuring that the estimated parameters align closely with their respective "true" values is crucial since they hold biological interpretation, leading to the challenge of uniqueness in the solutions. Structural identifiability analysis addresses this issue under the assumption of perfect observations of system dynamics, whereas practical identifiability considers limited measurements and the accompanying noise. In the literature, definitions for structural identifiability vary only slightly among authors, whereas the concept and quantification of practical identifiability lack consensus, with several indices coexisting. In this work, we provide a unified framework for studying identifiability, accommodating different definitions that need to be instantiated depending on each application case. In a more applicative second step, we focus on three widely used methods for quantifying practical identifiability: collinearity indices, profile likelihood, and average relative error. We show the limitations of their local versions, and we propose a new risk index built on the profile likelihood-based confidence intervals. We illustrate the usefulness of these concepts for plant growth modeling using a discrete-time individual plant growth model, LNAS, and a continuous-time plant population epidemics model. Through this work, we aim to underline the significance of identifiability analysis as a complement to any parameter estimation study and offer guidance to the modeler.

Personalized Computational Causal Modeling of the Alzheimer Disease Biomarker Cascade

BACKGROUND

Mathematical models of complex diseases, such as Alzheimer's disease, have the potential to play a significant role in personalized medicine. Specifically, models can be personalized by fitting parameters with individual data for the purpose of discovering primary underlying disease drivers, predicting natural history, and assessing the effects of theoretical interventions. Previous work in causal/mechanistic modeling of Alzheimer's Disease progression has modeled the disease at the cellular level and on a short time scale, such as minutes to hours. No previous studies have addressed mechanistic modeling on a personalized level using clinically validated biomarkers in individual subjects.

OBJECTIVES

This study aimed to investigate the feasibility of personalizing a causal model of Alzheimer's Disease progression using longitudinal biomarker data.

DESIGN/SETTING/PARTICIPANTS/MEASUREMENTS

We chose the Alzheimer Disease Biomarker Cascade model, a widely-referenced hypothetical model of Alzheimer's Disease based on the amyloid cascade hypothesis, which we had previously implemented mathematically as a mechanistic model. We used available longitudinal demographic and serial biomarker data in over 800 subjects across the cognitive spectrum from the Alzheimer's Disease Neuroimaging Initiative. The data included participants that were cognitively normal, had mild cognitive impairment, or were diagnosed with dementia (probable Alzheimer's Disease). The model consisted of a sparse system of differential equations involving four measurable biomarkers based on cerebrospinal fluid proteins, imaging, and cognitive testing data.

RESULTS

Personalization of the Alzheimer Disease Biomarker Cascade model with individual serial biomarker data yielded fourteen personalized parameters in each subject reflecting physiologically meaningful characteristics. These included growth rates, latency values, and carrying capacities of the various biomarkers, most of which demonstrated significant differences across clinical diagnostic groups. The model fits to training data across the entire cohort had a root mean squared error (RMSE) of 0.09 (SD 0.081) on a variable scale between zero and one, and were robust, with over 90% of subjects showing an RMSE of < 0.2. Similarly, in a subset of subjects with data on all four biomarkers in at least one test set, performance was high on the test sets, with a mean RMSE of 0.15 (SD 0.117), with 80% of subjects demonstrating an RMSE < 0.2 in the estimation of future biomarker points. Cluster analysis of parameters revealed two distinct endophenotypic groups, with distinct biomarker profiles and disease trajectories.

CONCLUSION

Results support the feasibility of personalizing mechanistic models based on individual biomarker trajectories and suggest that this approach may be useful for reclassifying subjects on the Alzheimer's clinical spectrum. This computational modeling approach is not limited to the Alzheimer Disease Biomarker Cascade hypothesis, and can be applied to any mechanistic hypothesis of disease progression in the Alzheimer's field that can be monitored with biomarkers. Thus, it offers a computational platform to compare and validate various disease hypotheses, personalize individual biomarker trajectories and predict individual response to theoretical prevention and therapeutic intervention strategies.

Fitting parameters and therapies of ODE tumor models with senescence and immune system

This work is devoted to introduce and study two quasispecies nonlinear ODE systems that model the behavior of tumor cell populations organized in different states. In the first model, replicative, senescent, extended lifespan, immortal and tumor cells are considered, while the second also includes immune cells. We fit the parameters regulating the transmission between states in order to approximate the outcomes of the models to a real progressive tumor invasion. After that, we study the identifiability of the fitted parameters, by using two sensitivity analysis methods. Then, we show that an adequate reduced fitting process (only accounting to the most identifiable parameters) gives similar results but saving computational cost. Three different therapies are introduced in the models to shrink (progressively in time) the tumor, while the replicative and senescent cells invade. Each therapy is identified to a dimensionless parameter. Then, we make a fitting process of the therapies' parameters, in various scenarios depending on the initial tumor according to the time when the therapies started. We conclude that, although the optimal combination of therapies depends on the size of initial tumor, the most efficient therapy is the reinforcement of the immune system. Finally, an identifiability analysis allows us to detect a limitation in the therapy outcomes. In fact, perturbing the optimal combination of therapies under an appropriate therapeutic vector produces virtually the same results.

Quantifying Parameter Interdependence in Stochastic Discrete Models of Biochemical Systems

Stochastic modeling of biochemical processes at the cellular level has been the subject of intense research in recent years. The Chemical Master Equation is a broadly utilized stochastic discrete model of such processes. Numerous important biochemical systems consist of many species subject to many reactions. As a result, their mathematical models depend on many parameters. In applications, some of the model parameters may be unknown, so their values need to be estimated from the experimental data. However, the problem of parameter value inference can be quite challenging, especially in the stochastic setting. To estimate accurately the values of a subset of parameters, the system should be sensitive with respect to variations in each of these parameters and they should not be correlated. In this paper, we propose a technique for detecting collinearity among models' parameters and we apply this method for selecting subsets of parameters that can be estimated from the available data. The analysis relies on finite-difference sensitivity estimations and the singular value decomposition of the sensitivity matrix. We illustrated the advantages of the proposed method by successfully testing it on several models of biochemical systems of practical interest.

Reaction kinetics determination and uncertainty analysis for the synthesis of the cancer drug lomustine

Fast and reliable model development frameworks are required to support current trends in modernization of pharmaceutical processing, promoting the use of digital platforms to assist process design and operation. In this work, we use a parameter estimation framework built into the PharmaPy library to determine rate parameters and uncertainty regions of different mechanistic and semi-empirical kinetic expressions for the synthesis of the drug lomustine. The parameter estimation procedure was complemented by identifiability analysis, resulting in simplified reaction mechanisms. Comparison of parameters and their uncertainty in process design was demonstrated through design space analysis, showing important differences in model prediction and the extent of their corresponding design spaces. The results of this work can serve to analyze lomustine manufacturing processes that include separation and isolation steps, where parametric sensitivity is expected to propagate along the manufacturing line and impact process feasible operation, and attainment of critical quality attributes of the product.

Combining Multikinase Tyrosine Kinase Inhibitors Targeting the Vascular Endothelial Growth Factor and Cluster of Differentiation 47 Signaling Pathways Is Predicted to Increase the Efficacy of Antiangiogenic Combination Therapies

Angiogenesis is a critical step in tumor growth, development, and invasion. Nascent tumor cells secrete vascular endothelial growth factor (VEGF) that significantly remodels the tumor microenvironment through interaction with multiple receptors on vascular endothelial cells, including type 2 VEGF receptor (VEGFR2). The complex pathways initiated by VEGF binding to VEGFR2 lead to enhanced proliferation, survival, and motility of vascular endothelial cells and formation of a new vascular network, enabling tumor growth. Antiangiogenic therapies that inhibit VEGF signaling pathways were among the first drugs that targeted stroma rather than tumor cells. Despite improvements in progression-free survival and higher response rates relative to chemotherapy in some types of solid tumors, the impact on overall survival (OS) has been limited, with the majority of tumors eventually relapsing due to resistance or activation of alternate angiogenic pathways. Here, we developed a molecularly detailed computational model of endothelial cell signaling and angiogenesis-driven tumor growth to investigate combination therapies targeting different nodes of the endothelial VEGF/VEGFR2 signaling pathway. Simulations predicted a strong threshold-like behavior in extracellular signal-regulated kinases 1/2 (ERK1/2) activation relative to phosphorylated VEGFR2 levels, as continuous inhibition of at least 95% of receptors was necessary to abrogate phosphorylated ERK1/2 (pERK1/2). Combinations with mitogen-activated protein kinase/ERK kinase (MEK) and spingosine-1-phosphate inhibitors were found to be effective in overcoming the ERK1/2 activation threshold and abolishing activation of the pathway. Modeling results also identified a mechanism of resistance whereby tumor cells could reduce pERK1/2 sensitivity to inhibitors of VEGFR2 by upregulation of Raf, MEK, and sphingosine kinase 1 (SphK1), thus highlighting the need for deeper investigation of the dynamics of the crosstalk between VEGFR2 and SphK1 pathways. Inhibition of VEGFR2 phosphorylation was found to be more effective at blocking protein kinase B, also known as AKT, activation; however, to effectively abolish AKT activation, simulations identified Axl autophosphorylation or the Src kinase domain as potent targets. Simulations also supported activating cluster of differentiation 47 (CD47) on endothelial cells as an effective combination partner with tyrosine kinase inhibitors to inhibit angiogenesis signaling and tumor growth. Virtual patient simulations supported the effectiveness of CD47 agonism in combination with inhibitors of VEGFR2 and SphK1 pathways. Overall, the rule-based system model developed here provides new insights, generates novel hypothesis, and makes predictions regarding combinations that may enhance the OS with currently approved antiangiogenic therapies.

Optimal enzyme utilization suggests that concentrations and thermodynamics determine binding mechanisms and enzyme saturations

Deciphering the metabolic functions of organisms requires understanding the dynamic responses of living cells upon genetic and environmental perturbations, which in turn can be inferred from enzymatic activity. In this work, we investigate the optimal modes of operation for enzymes in terms of the evolutionary pressure driving them toward increased catalytic efficiency. We develop a framework using a mixed-integer formulation to assess the distribution of thermodynamic forces and enzyme states, providing detailed insights into the enzymatic mode of operation. We use this framework to explore Michaelis-Menten and random-ordered multi-substrate mechanisms. We show that optimal enzyme utilization is achieved by unique or alternative operating modes dependent on reactant concentrations. We find that in a bimolecular enzyme reaction, the random mechanism is optimal over any other ordered mechanism under physiological conditions. Our framework can investigate the optimal catalytic properties of complex enzyme mechanisms. It can further guide the directed evolution of enzymes and fill in the knowledge gaps in enzyme kinetics.

Single-cell multi-omics and lineage tracing to dissect cell fate decision-making

The concept of cell fate relates to the future identity of a cell, and its daughters, which is obtained via cell differentiation and division. Understanding, predicting, and manipulating cell fate has been a long-sought goal of developmental and regenerative biology. Recent insights obtained from single-cell genomic and integrative lineage-tracing approaches have further aided to identify molecular features predictive of cell fate. In this perspective, we discuss these approaches with a focus on theoretical concepts and future directions of the field to dissect molecular mechanisms underlying cell fate.

Sloppiness: Fundamental study, new formalism and its application in model assessment

Computational modelling of biological processes poses multiple challenges in each stage of the modelling exercise. Some significant challenges include identifiability, precisely estimating parameters from limited data, informative experiments and anisotropic sensitivity in the parameter space. One of these challenges' crucial but inconspicuous sources is the possible presence of large regions in the parameter space over which model predictions are nearly identical. This property, known as sloppiness, has been reasonably well-addressed in the past decade, studying its possible impacts and remedies. However, certain critical unanswered questions concerning sloppiness, particularly related to its quantification and practical implications in various stages of system identification, still prevail. In this work, we systematically examine sloppiness at a fundamental level and formalise two new theoretical definitions of sloppiness. Using the proposed definitions, we establish a mathematical relationship between the parameter estimates' precision and sloppiness in linear predictors. Further, we develop a novel computational method and a visual tool to assess the goodness of a model around a point in parameter space by identifying local structural identifiability and sloppiness and finding the most sensitive and least sensitive parameters for non-infinitesimal perturbations. We demonstrate the working of our method in benchmark systems biology models of various complexities. The pharmacokinetic HIV infection model analysis identified a new set of biologically relevant parameters that can be used to control the free virus in an active HIV infection.

Virtual patient framework for the testing of mechanical ventilation airway pressure and flow settings protocol

BACKGROUND AND OBJECTIVE

Model-based and personalised decision support systems are emerging to guide mechanical ventilation (MV) treatment for respiratory failure patients. However, model-based treatments require resource-intensive clinical trials prior to implementation. This research presents a framework for generating virtual patients for testing model-based decision support, and direct use in MV treatment.

METHODS

The virtual MV patient framework consists of 3 stages: 1) Virtual patient generation, 2) Patient-level validation, and 3) Virtual clinical trials. The virtual patients are generated from retrospective MV patient data using a clinically validated respiratory mechanics model whose respiratory parameters (respiratory elastance and resistance) capture patient-specific pulmonary conditions and responses to MV care over time. Patient-level validation compares the predicted responses from the virtual patient to their retrospective results for clinically implemented MV settings and changes to care. Patient-level validated virtual patients create a platform to conduct virtual trials, where the safety of closed-loop model-based protocols can be evaluated.

RESULTS

This research creates and presents a virtual patient platform of 100 virtual patients generated from retrospective data. Patient-level validation reported median errors of 3.26% for volume-control and 6.80% for pressure-control ventilation mode. A virtual trial on a model-based protocol demonstrates the potential efficacy of using virtual patients for prospective evaluation and testing of the protocol.

CONCLUSION

The virtual patient framework shows the potential to safely and rapidly design, develop, and optimise new model-based MV decision support systems and protocols using clinically validated models and computer simulation, which could ultimately improve patient care and outcomes in MV.

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

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

BMSS2: A Unified Database-Driven Modeling Tool for Systematic Biomodel Selection

Modeling in synthetic biology constitutes a powerful means in our continuous search for improved performance with a rational Design-Build-Test-Learn approach. Particularly, kinetic models unravel system dynamics and enable system analysis for guiding experimental design. However, a systematic yet modular pipeline that allows one to identify the appropriate model and guide the experimental designs while tracing the entire model development and analysis is still lacking. Here, we develop BMSS2, a unified tool that streamlines and automates model selection by combining information criterion ranking with upstream and parallel analysis algorithms. These include Bayesian parameter inference, a priori and a posteriori identifiability analysis, and global sensitivity analysis. In addition, the database-driven design supports interactive model storage/retrieval to encourage reusability and facilitate automated model selection. This allows ease of model manipulation and deposition for the selection and analysis, thus enabling better utilization of models in guiding experimental design.

A microfluidic optimal experimental design platform for forward design of cell-free genetic networks

Cell-free protein synthesis has been widely used as a “breadboard” for design of synthetic genetic networks. However, due to a severe lack of modularity, forward engineering of genetic networks remains challenging. Here, we demonstrate how a combination of optimal experimental design and microfluidics allows us to devise dynamic cell-free gene expression experiments providing maximum information content for subsequent non-linear model identification. Importantly, we reveal that applying this methodology to a library of genetic circuits, that share common elements, further increases the information content of the data resulting in higher accuracy of model parameters. To show modularity of model parameters, we design a pulse decoder and bistable switch, and predict their behaviour both qualitatively and quantitatively. Finally, we update the parameter database and indicate that network topology affects parameter estimation accuracy. Utilizing our methodology provides us with more accurate model parameters, a necessity for forward engineering of complex genetic networks.

Characterization of cell-free genetic networks is inherently difficult. Here the authors use optimal experimental design and microfluidics to improve characterization, demonstrating modularity and predictability of parts in applied test cases.

Modelling the performance of an integrated fixed-film activated sludge (IFAS) system: a systematic approach to automated calibration

IFAS systems are inherently complex due to the hybrid use of both suspended and attached bacterial colonies for the purpose of pollutant degradation as part of wastewater treatment. This poses challenges when attempting to represent these systems mathematically due to the vast number of parameters involved. Besides becoming convoluted, large effort will be incurred during model calibration. This paper demonstrates a systematic approach to calibration of an IFAS process model that incorporates two sensitivity analyses to identify influential parameters and detect collinearity from a subset of 68 kinetic and stoichiometric parameters, and the use of the Nelder–Mead optimization algorithm to estimate the required values of these parameters. The model considers the removal of three critical pollutants including biochemical oxygen demand (BOD), total nitrogen (TN) and total suspended solids (TSS). Results from the sensitivity analyses identified four parameters that were the primary influence on the model. The model was found to be most sensitive to the two stoichiometric parameters including aerobic heterotrophic yield on soluble substrate whose total effects were responsible for 92.4% of the model’s BOD output sensitivity and 92.8% of the model’s TSS output sensitivity. The anoxic heterotrophic yield on soluble substrate was observed to be responsible for 54.3% of the model’s TN output sensitivity. To a lesser extent the two kinetic parameters, aerobic heterotrophic decay rate and reduction factor for denitrification on nitrite, were responsible for only 8.0% and 13.1% of the model’s BOD and TN output sensitivities respectively. Parameter estimation identified the need for only minor adjustments to default values in order to achieve sufficient accuracy of simulation with deviation from observed data to be only ± 3.6 mg/L, ± 1.3 mg/L, and ± 9.5 mg/L for BOD, TN and TSS respectively. Validation showed the model was limited in its capacity to predict system behaviour under extreme dissolved oxygen stress.

Model selection of chaotic systems from data with hidden variables using sparse data assimilation

Many natural systems exhibit chaotic behavior, including the weather, hydrology, neuroscience, and population dynamics. Although many chaotic systems can be described by relatively simple dynamical equations, characterizing these systems can be challenging due to sensitivity to initial conditions and difficulties in differentiating chaotic behavior from noise. Ideally, one wishes to find a parsimonious set of equations that describe a dynamical system. However, model selection is more challenging when only a subset of the variables are experimentally accessible. Manifold learning methods using time-delay embeddings can successfully reconstruct the underlying structure of the system from data with hidden variables, but not the equations. Recent work in sparse-optimization based model selection has enabled model discovery given a library of possible terms, but regression-based methods require measurements of all state variables. We present a method combining variational annealing-a technique previously used for parameter estimation in chaotic systems with hidden variables-with sparse-optimization methods to perform model identification for chaotic systems with unmeasured variables. We applied the method to ground-truth time-series simulated from the classic Lorenz system and experimental data from an electrical circuit with Lorenz-system like behavior. In both cases, we successfully recover the expected equations with two measured and one hidden variable. Application to simulated data from the Colpitts oscillator demonstrates successful model selection of terms within nonlinear functions. We discuss the robustness of our method to varying noise.

A Bayesian Approach to Extracting Kinetic Information from Artificial Enzymatic Networks

In order to create artificial enzymatic networks capable of increasingly complex behavior, an improved methodology in understanding and controlling the kinetics of these networks is needed. Here, we introduce a Bayesian analysis method allowing for the accurate inference of enzyme kinetic parameters and determination of most likely reaction mechanisms, by combining data from different experiments and network topologies in a single probabilistic analysis framework. This Bayesian approach explicitly allows us to continuously improve our parameter estimates and behavior predictions by iteratively adding new data to our models, while automatically taking into account uncertainties introduced by the experimental setups or the chemical processes in general. We demonstrate the potential of this approach by characterizing systems of enzymes compartmentalized in beads inside flow reactors. The methods we introduce here provide a new approach to the design of increasingly complex artificial enzymatic networks, making the design of such networks more efficient, and robust against the accumulation of experimental errors.

Quantification of respiratory effort magnitude in spontaneous breathing patients using Convolutional Autoencoders

BACKGROUND

Spontaneous breathing (SB) effort during mechanical ventilation (MV) is an important metric of respiratory drive. However, SB effort varies due to a variety of factors, including evolving pathology and sedation levels. Therefore, assessment of SB efforts needs to be continuous and non-invasive. This is important to prevent both over- and under-assistance with MV. In this study, a machine learning model, Convolutional Autoencoder (CAE) is developed to quantify the magnitude of SB effort using only bedside MV airway pressure and flow waveform.

METHOD

The CAE model was trained using 12,170,655 simulated SB flow and normal flow data (NB). The paired SB and NB flow data were simulated using a Gaussian Effort Model (GEM) with 5 basis functions. When the CAE model is given a SB flow input, it is capable of predicting a corresponding NB flow for the SB flow input. The magnitude of SB effort (SBE_Mag) is then quantified as the difference between the SB and NB flows. The CAE model was used to evaluate the SBE_Mag of 9 pressure control/ support datasets. Results were validated using a mean squared error (MSE) fitting between clinical and training SB flows.

RESULTS

The CAE model was able to produce NB flows from the clinical SB flows with the median SBE_Mag of the 9 datasets being 25.39% [IQR: 21.87-25.57%]. The absolute error in SBE_Mag using MSE validation yields a median of 4.77% [IQR: 3.77-8.56%] amongst the cohort. This shows the ability of the GEM to capture the intrinsic details present in SB flow waveforms. Analysis also shows both intra-patient and inter-patient variability in SBE_Mag.

CONCLUSION

A Convolutional Autoencoder model was developed with simulated SB and NB flow data and is capable of quantifying the magnitude of patient spontaneous breathing effort. This provides potential application for real-time monitoring of patient respiratory drive for better management of patient-ventilator interaction.

Parameter Estimation in the Age of Degeneracy and Unidentifiability

Parameter estimation from observable or experimental data is a crucial stage in any modeling study. Identifiability refers to one’s ability to uniquely estimate the model parameters from the available data. Structural unidentifiability in dynamic models, the opposite of identifiability, is associated with the notion of degeneracy where multiple parameter sets produce the same pattern. Therefore, the inverse function of determining the model parameters from the data is not well defined. Degeneracy is not only a mathematical property of models, but it has also been reported in biological experiments. Classical studies on structural unidentifiability focused on the notion that one can at most identify combinations of unidentifiable model parameters. We have identified a different type of structural degeneracy/unidentifiability present in a family of models, which we refer to as the Lambda-Omega (Λ-Ω) models. These are an extension of the classical lambda-omega (λ-ω) models that have been used to model biological systems, and display a richer dynamic behavior and waveforms that range from sinusoidal to square wave to spike like. We show that the Λ-Ω models feature infinitely many parameter sets that produce identical stable oscillations, except possible for a phase shift (reflecting the initial phase). These degenerate parameters are not identifiable combinations of unidentifiable parameters as is the case in structural degeneracy. In fact, reducing the number of model parameters in the Λ-Ω models is minimal in the sense that each one controls a different aspect of the model dynamics and the dynamic complexity of the system would be reduced by reducing the number of parameters. We argue that the family of Λ-Ω models serves as a framework for the systematic investigation of degeneracy and identifiability in dynamic models and for the investigation of the interplay between structural and other forms of unidentifiability resulting on the lack of information from the experimental/observational data.

A validated mathematical model of FGFR3-mediated tumor growth reveals pathways to harness the benefits of combination targeted therapy and immunotherapy in bladder cancer

Bladder cancer is a common malignancy with over 80,000 estimated new cases and nearly 18,000 deaths per year in the United States alone. Therapeutic options for metastatic bladder cancer had not evolved much for nearly four decades, until recently, when five immune checkpoint inhibitors were approved by the U.S. Food and Drug Administration (FDA). Despite the activity of these drugs in some patients, the objective response rate for each is less than 25%. At the same time, fibroblast growth factor receptors (FGFRs) have been attractive drug targets for a variety of cancers, and in 2019 the FDA approved the first therapy targeted against FGFR3 for bladder cancer. Given the excitement around these new receptor tyrosine kinase and immune checkpoint targeted strategies, and the challenges they each may face on their own, emerging data suggest that combining these treatment options could lead to improved therapeutic outcomes. In this paper, we develop a mathematical model for FGFR3-mediated tumor growth and use it to investigate the impact of the combined administration of a small molecule inhibitor of FGFR3 and a monoclonal antibody against the PD-1/PD-L1 immune checkpoint. The model is carefully calibrated and validated with experimental data before survival benefits, and dosing schedules are explored. Predictions of the model suggest that FGFR3 mutation reduces the effectiveness of anti-PD-L1 therapy, that there are regions of parameter space where each monotherapy can outperform the other, and that pretreatment with anti-PD-L1 therapy always results in greater tumor reduction even when anti-FGFR3 therapy is the more effective monotherapy.

Two-Age-Structured COVID-19 Epidemic Model: Estimation of Virulence Parameters to Interpret Effects of National and Regional Feedback Interventions and Vaccination

The COVID-19 epidemic has recently led in Italy to the implementation of different external strategies in order to limit the spread of the disease in response to its transmission rate: strict national lockdown rules, followed first by a weakening of the social distancing and contact reduction feedback interventions and finally the implementation of coordinated intermittent regional actions, up to the application, in this last context, of an age-stratified vaccine prioritization strategy. This paper originally aims at identifying, starting from the available age-structured real data at the national level during the specific aforementioned scenarios, external-scenario-dependent sets of virulence parameters for a two-age-structured COVID-19 epidemic compartmental model, in order to provide an interpretation of how each external scenario modifies the age-dependent patterns of social contacts and the spread of COVID-19.

Assessing the role of initial conditions in the local structural identifiability of large dynamic models

Structural identifiability is a binary property that determines whether or not unique parameter values can, in principle, be estimated from error-free input-output data. The many papers that have been written on this topic collectively stress the importance of this a priori analysis in the model development process. The story however, often ends with a structurally unidentifiable model. This may leave a model developer with no plan of action on how to address this potential issue. We continue this model exploration journey by identifying one of the possible sources of a model's unidentifiability: problematic initial conditions. It is well-known that certain initial values may result in the loss of local structural identifiability. Nevertheless, literature on this topic has been limited to the analysis of small toy models. Here, we present a systematic approach to detect problematic initial conditions of real-world systems biology models, that are usually not small. A model's identifiability can often be reinstated by changing the value of such problematic initial conditions. This provides modellers an option to resolve the "unidentifiable model" problem. Additionally, a good understanding of which initial values should rather be avoided can be very useful during experimental design. We show how our approach works in practice by applying it to five models. First, two small benchmark models are studied to get the reader acquainted with the method. The first one shows the effect of a zero-valued problematic initial condition. The second one illustrates that the approach also yields correct results in the presence of input signals and that problematic initial conditions need not be zero-values. For the remaining three examples, we set out to identify key initial values which may result in the structural unidentifiability. The third and fourth examples involve a systems biology Epo receptor model and a JAK/STAT model, respectively. In the final Pharmacokinetics model, of which its global structural identifiability has only recently been confirmed, we indicate that there are still sets of initial values for which this property does not hold.

Introducing Parameter Clustering to the OED Procedure for Model Calibration of a Synthetic Inducible Promoter in S. cerevisiae

In recent years, synthetic gene circuits for adding new cell features have become one of the most powerful tools in biological and pharmaceutical research and development. However, because of the inherent non-linearity and noisy experimental data, the experiment-based model calibration of these synthetic parts is perceived as a laborious and time-consuming procedure. Although the optimal experimental design (OED) based on the Fisher information matrix (FIM) has been proved to be an effective means to improve the calibration efficiency, the required calculation increases dramatically with the model size (parameter number). To reduce the OED complexity without losing the calibration accuracy, this paper proposes two OED approaches with different parameter clustering methods and validates the accuracy of calibrated models with in-silico experiments. A model of an inducible synthetic promoter in S. cerevisiae is adopted for bench-marking. The comparison with the traditional off-line OED approach suggests that the OED approaches with both of the clustering methods significantly reduce the complexity of OED problems (for at least 49.0%), while slightly improving the calibration accuracy (11.8% and 19.6% lower estimation error in average for FIM-based and sensitivity-based approaches). This study implicates that for calibrating non-linear models of biological pathways, cluster-based OED could be a beneficial approach to improve the efficiency of optimal experimental design.

Introduction to In Silico Modeling to Study ROS Dynamics

This book chapter is drafted for biologists with experimental experiences in ROS biology but being newcomers in the field of modeling. We start with a general introduction about computational modeling in biology and an overview of software tools suitable for beginners. This chapter encompasses an introduction to computational models with special focus on simulation of ROS dynamics. A step-by-step tutorial follows providing guidance for all relevant model development processes. This course of action gives a comprehensible way to understand the benefits of computational models and to gain the necessary knowledge to build own small equation-based models. Small models can be created without any special programming expertise or in-depth technical and mathematical knowledge. Afterward in the final section, a short overview of pitfalls, challenges, and limitations is provided, combined with suggestions for further reading to improve and expand modeling skills of biologists.

Building kinetic models for metabolic engineering

Kinetic formalisms of metabolism link metabolic fluxes to enzyme levels, metabolite concentrations and their allosteric regulatory interactions. Though they require the identification of physiologically relevant values for numerous parameters, kinetic formalisms uniquely establish a mechanistic link across heterogeneous omics datasets and provide an overarching vantage point to effectively inform metabolic engineering strategies. Advances in computational power, gene annotation coverage, and formalism standardization have led to significant progress over the past few years. However, careful interpretation of model predictions, limited metabolic flux datasets, and assessment of parameter sensitivity remain as challenges. In this review we highlight fundamental considerations which influence model quality and prediction, advances in methodologies, and success stories of deploying kinetic models to guide metabolic engineering.

Recent advances in constraint and machine learning-based metabolic modeling by leveraging stoichiometric balances, thermodynamic feasibility and kinetic law formalisms

Understanding the governing principles behind organisms' metabolism and growth underpins their effective deployment as bioproduction chassis. A central objective of metabolic modeling is predicting how metabolism and growth are affected by both external environmental factors and internal genotypic perturbations. The fundamental concepts of reaction stoichiometry, thermodynamics, and mass action kinetics have emerged as the foundational principles of many modeling frameworks designed to describe how and why organisms allocate resources towards both growth and bioproduction. This review focuses on the latest algorithmic advancements that have integrated these foundational principles into increasingly sophisticated quantitative frameworks.

New Conceptual Toxicokinetic Model to Assess Synergistic Mixture Effects between the Aromatic Hydrocarbon β-Naphthoflavone and the Azole Nocodazole on the CYP1A Biomarker in a Fish Cell Line

Toxicokinetic interactions with catabolic cytochrome P450 (CYP) enzymes can inhibit chemical elimination pathways and cause synergistic mixture effects. We have created a mathematical bottom-up model for a synergistic mixture effect where we fit a multidimensional function to a given data set using an auxiliary nonadditive approach. The toxicokinetic model is based on the data from a previous study on a fish cell line, where the CYP1A enzyme activity was measured over time after exposure to various combinations of the aromatic hydrocarbon β-naphthoflavone and the azole nocodazole. To describe the toxicokinetic mechanism in this pathway and how that affects the CYP1A biomarker, the model uses ordinary differential equations. Local sensitivity and identifiability analyses revealed that all the 10 parameters estimated in the model were identified uniquely while fitting the model to the data for measuring the CYP1A enzyme activity. The model has a good prediction power and is a promising tool to test the synergistic toxicokinetic interactions between different chemicals.

Systems biology informed deep learning for inferring parameters and hidden dynamics

Mathematical models of biological reactions at the system-level lead to a set of ordinary differential equations with many unknown parameters that need to be inferred using relatively few experimental measurements. Having a reliable and robust algorithm for parameter inference and prediction of the hidden dynamics has been one of the core subjects in systems biology, and is the focus of this study. We have developed a new systems-biology-informed deep learning algorithm that incorporates the system of ordinary differential equations into the neural networks. Enforcing these equations effectively adds constraints to the optimization procedure that manifests itself as an imposed structure on the observational data. Using few scattered and noisy measurements, we are able to infer the dynamics of unobserved species, external forcing, and the unknown model parameters. We have successfully tested the algorithm for three different benchmark problems.

Modeling anaerobic digestion metabolic pathways for antibiotic-contaminated wastewater treatment

Anaerobic digestion has been used to treat antibiotic-contaminated wastewaters. However, it is not always effective, since biodegradation is the main removal mechanism and depends on the compound chemical characteristics and on how microbial metabolic pathways are affected by the reactor operational conditions and hydrodynamic characteristics. The aim of this study was to develop a mathematical model to describe 16 metabolic pathways of an anaerobic process treating sulfamethazine-contaminated wastewater. Contois kinetics and a useful reaction volume term were used to represent the biomass concentration impact on bed porosity in a N continuously stirred tank modeling approach. Two sulfamethazine removal hypotheses were evaluated: an apparent enzymatic reaction and a cometabolic degradation. Additionally, long-term modeling was developed to describe how the operational conditions affected the performance of the process. The best degradation correlations were associated with the consumption of carbohydrates, proteins and it was inversely related to acetic acid production during acidogenesis.

An Enhanced Segment Particle Swarm Optimization Algorithm for Kinetic Parameters Estimation of the Main Metabolic Model of Escherichia Coli

Building a biologic model that describes the behavior of a cell in biologic systems is aimed at understanding the physiology of the cell, predicting the production of enzymes and metabolites, and providing a suitable data that is valid for bio-products. In addition, building a kinetic model requires the estimation of the kinetic parameters, but kinetic parameters estimation in kinetic modeling is a difficult task due to the nonlinearity of the model. As a result, kinetic parameters are mostly reported or estimated from different laboratories in different conditions and time consumption. Hence, based on the aforementioned problems, the optimization algorithm methods played an important role in addressing these problems. In this study, an Enhanced Segment Particle Swarm Optimization algorithm (ESe-PSO) was proposed for kinetic parameters estimation. This method was proposed to increase the exploration and the exploitation of the Segment Particle Swarm Optimization algorithm (Se-PSO). The main metabolic model of E. coli was used as a benchmark which contained 172 kinetic parameters distributed in five pathways. Seven kinetic parameters were well estimated based on the distance minimization between the simulation and the experimental results. The results revealed that the proposed method had the ability to deal with kinetic parameters estimation in terms of time consumption and distance minimization.

A scalable method for parameter identification in kinetic models of metabolism using steady-state data

MOTIVATION

In kinetic models of metabolism, the parameter values determine the dynamic behaviour predicted by these models. Estimating parameters from in vivo experimental data require the parameters to be structurally identifiable, and the data to be informative enough to estimate these parameters. Existing methods to determine the structural identifiability of parameters in kinetic models of metabolism can only be applied to models of small metabolic networks due to their computational complexity. Additionally, a priori experimental design, a necessity to obtain informative data for parameter estimation, also does not account for using steady-state data to estimate parameters in kinetic models.

RESULTS

Here, we present a scalable methodology to structurally identify parameters for each flux in a kinetic model of metabolism based on the availability of steady-state data. In doing so, we also address the issue of determining the number and nature of experiments for generating steady-state data to estimate these parameters. By using a small metabolic network as an example, we show that most parameters in fluxes expressed by mechanistic enzyme kinetic rate laws can be identified using steady-state data, and the steady-state data required for their estimation can be obtained from selective experiments involving both substrate and enzyme level perturbations. The methodology can be used in combination with other identifiability and experimental design algorithms that use dynamic data to determine the most informative experiments requiring the least resources to perform.

AVAILABILITY AND IMPLEMENTATION

https://github.com/LMSE/ident.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

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

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

Black-box kinetic modeling of growth and citric acid production and estimation of ATP maintenance parameters for Candida oleophila ATCC20177

Mathematical modeling represents and predicts biological systems, explains underlying mechanisms, constituting one of the key focus points for fundamental and applied research to improve our understanding and to decrease costs. Organic acids are used in several industries such as monomers for bioplastics, food preservatives and additives, pharmaceuticals, and agriculture. Nonpetrochemical, sustainable production of organic acids is therefore of great interest. An important step in production of organic acids is the determination of growth and acid production dynamics, as the product itself may have direct and indirect inhibitory effects on the host's metabolism. The aim of this study it twofold: (i) to determine the parameters related to energetics of growth and production as growth ( K x ) and nongrowth associated (m_ATP ) maintenance constants and (ii) to set up and analyze an unstructured, black-box kinetic model to describe the dynamics of the growth and production of citric acid by Candida oleophila ATCC20177 using published batch fermentation data. K x and m_ATP were found to be 2.3 ± 1.7 and 5.25 ± 2.75, respectively, for the published P/O ratio of 1.45. The parameter sensitivities and correlations are determined using the Monte Carlo approach, and the final model is tested using chemostat data.

Sensitivity, uncertainty and identifiability analyses to define a dengue transmission model with real data of an endemic municipality of Colombia

Dengue disease is a major problem for public health surveillance entities in tropical and subtropical regions having a significant impact not only epidemiological but social and economical. There are many factors involved in the dengue transmission process. We can evaluate the importance of these factors through the formulation of mathematical models. However, the majority of the models presented in the literature tend to be overparameterized, with considerable uncertainty levels and excessively complex formulations. We aim to evaluate the structure, complexity, trustworthiness, and suitability of three models, for the transmission of dengue disease, through different strategies. To achieve this goal, we perform structural and practical identifiability, sensitivity and uncertainty analyses to these models. The results showed that the simplest model was the most appropriate and reliable when the only available information to fit them is the cumulative number of reported dengue cases in an endemic municipality of Colombia.

Identifying important parameters in the inflammatory process with a mathematical model of immune cell influx and macrophage polarization

In an inflammatory setting, macrophages can be polarized to an inflammatory M1 phenotype or to an anti-inflammatory M2 phenotype, as well as existing on a spectrum between these two extremes. Dysfunction of this phenotypic switch can result in a population imbalance that leads to chronic wounds or disease due to unresolved inflammation. Therapeutic interventions that target macrophages have therefore been proposed and implemented in diseases that feature chronic inflammation such as diabetes mellitus and atherosclerosis. We have developed a model for the sequential influx of immune cells in the peritoneal cavity in response to a bacterial stimulus that includes macrophage polarization, with the simplifying assumption that macrophages can be classified as M1 or M2. With this model, we were able to reproduce the expected timing of sequential influx of immune cells and mediators in a general inflammatory setting. We then fit this model to in vivo experimental data obtained from a mouse peritonitis model of inflammation, which is widely used to evaluate endogenous processes in response to an inflammatory stimulus. Model robustness is explored with local structural and practical identifiability of the proposed model a posteriori. Additionally, we perform sensitivity analysis that identifies the population of apoptotic neutrophils as a key driver of the inflammatory process. Finally, we simulate a selection of proposed therapies including points of intervention in the case of delayed neutrophil apoptosis, which our model predicts will result in a sustained inflammatory response. Our model can therefore provide hypothesis testing for therapeutic interventions that target macrophage phenotype and predict outcomes to be validated by subsequent experimentation.

Using experimental data and mathematical analysis, we develop a model for the inflammatory response that includes macrophage polarization between M1 and M2 phenotypes. Dysfunction of this phenotypic switch can disrupt the timely influx and egress of immune cells during the healing process and lead to chronic wounds or disease. The modulation of macrophage population has been suggested as a strategy to dampen inflammation in diseases that feature chronic inflammation, such as diabetes and atherosclerosis. It is therefore important that we learn more about which components of the system drive the population level switch in phenotype. Our model is able to reproduce the expected timing of sequential influx of neutrophils and macrophages in response to an inflammatory stimulus. Model parameters were estimated with weighted least squares fitting to in vivo experimental data from a mouse model of peritonitis while considering identifiability of parameter sets. We perform sensitivity analysis that identifies primary drivers of the system, and predict the effects of variations in these key parameters on immune cell populations.

Comprehensive Review of Models and Methods for Inferences in Bio-Chemical Reaction Networks

The key processes in biological and chemical systems are described by networks of chemical reactions. From molecular biology to biotechnology applications, computational models of reaction networks are used extensively to elucidate their non-linear dynamics. The model dynamics are crucially dependent on the parameter values which are often estimated from observations. Over the past decade, the interest in parameter and state estimation in models of (bio-) chemical reaction networks (BRNs) grew considerably. The related inference problems are also encountered in many other tasks including model calibration, discrimination, identifiability, and checking, and optimum experiment design, sensitivity analysis, and bifurcation analysis. The aim of this review paper is to examine the developments in literature to understand what BRN models are commonly used, and for what inference tasks and inference methods. The initial collection of about 700 documents concerning estimation problems in BRNs excluding books and textbooks in computational biology and chemistry were screened to select over 270 research papers and 20 graduate research theses. The paper selection was facilitated by text mining scripts to automate the search for relevant keywords and terms. The outcomes are presented in tables revealing the levels of interest in different inference tasks and methods for given models in the literature as well as the research trends are uncovered. Our findings indicate that many combinations of models, tasks and methods are still relatively unexplored, and there are many new research opportunities to explore combinations that have not been considered-perhaps for good reasons. The most common models of BRNs in literature involve differential equations, Markov processes, mass action kinetics, and state space representations whereas the most common tasks are the parameter inference and model identification. The most common methods in literature are Bayesian analysis, Monte Carlo sampling strategies, and model fitting to data using evolutionary algorithms. The new research problems which cannot be directly deduced from the text mining data are also discussed.

A New Mass-Based Discretized Population Balance Model for Precipitation Processes: Application to Struvite Precipitation

Mathematical models describing precipitation processes in one step need to be upgraded. Particle size distribution is a crucial variable and its inclusion in the modelling libraries is necessary if the technology wants to be optimized through simulation. With this objective, a mass based population balance model is presented in this contribution. The model has been constructed using a stoichiometric matrix and a kinetic vector and using mass as the internal coordinate, as it is usually done in wastewater treatment modelling. Identifiability of the parameters of the model was evaluated using a sensitivity and a collinearity analysis for six simulation case studies of struvite precipitation. In addition, parameters in the model were calibrated to represent data from two batch tests in the laboratory. The results of the analysis showed that the identifiability of the parameters depends on the available experimental data and explored scenarios. Identifiability of the parameters could be the reason behind the shifting parameter values describing mechanisms of precipitation in the literature. This contribution helps to understand the possibilities and limitations that the population balance model approach offer.

Parameter estimation and identifiability in a neural population model for electro-cortical activity

Electroencephalography (EEG) provides a non-invasive measure of brain electrical activity. Neural population models, where large numbers of interacting neurons are considered collectively as a macroscopic system, have long been used to understand features in EEG signals. By tuning dozens of input parameters describing the excitatory and inhibitory neuron populations, these models can reproduce prominent features of the EEG such as the alpha-rhythm. However, the inverse problem, of directly estimating the parameters from fits to EEG data, remains unsolved. Solving this multi-parameter non-linear fitting problem will potentially provide a real-time method for characterizing average neuronal properties in human subjects. Here we perform unbiased fits of a 22-parameter neural population model to EEG data from 82 individuals, using both particle swarm optimization and Markov chain Monte Carlo sampling. We estimate how much is learned about individual parameters by computing Kullback-Leibler divergences between posterior and prior distributions for each parameter. Results indicate that only a single parameter, that determining the dynamics of inhibitory synaptic activity, is directly identifiable, while other parameters have large, though correlated, uncertainties. We show that the eigenvalues of the Fisher information matrix are roughly uniformly spaced over a log scale, indicating that the model is sloppy, like many of the regulatory network models in systems biology. These eigenvalues indicate that the system can be modeled with a low effective dimensionality, with inhibitory synaptic activity being prominent in driving system behavior.

Use of identifiability analysis in designing phenotyping experiments for modelling forage production and quality

Agricultural systems models are complex and tend to be over-parameterized with respect to observational datasets. Practical identifiability analysis based on local sensitivity analysis has proved effective in investigating identifiable parameter sets in environmental models, but has not been applied to agricultural systems models. Here, we demonstrate that identifiability analysis improves experimental design to ensure independent parameter estimation for yield and quality outputs of a complex grassland model. The Pasture Simulation model (PaSim) was used to demonstrate the effectiveness of practical identifiability analysis in designing experiments and measurement protocols within phenotyping experiments with perennial ryegrass. Virtual experiments were designed combining three factors: frequency of measurements, duration of the experiment. and location of trials. Our results demonstrate that (i) PaSim provides sufficient detail in terms of simulating biomass yield and quality of perennial ryegrass for use in breeding, (ii) typical breeding trials are insufficient to parameterize all influential parameters, (iii) the frequency of measurements is more important than the number of growing seasons to improve the identifiability of PaSim parameters, and (iv) identifiability analysis provides a sound approach for optimizing the design of multi-location trials. Practical identifiability analysis can play an important role in ensuring proper exploitation of phenotypic data and cost-effective multi-location experimental designs. Considering the growing importance of simulation models, this study supports the design of experiments and measurement protocols in the phenotyping networks that have recently been organized.

Mechanistic Models of Inducible Synthetic Circuits for Joint Description of DNA Copy Number, Regulatory Protein Level, and Cell Load

Accurate predictive mathematical models are urgently needed in synthetic biology to support the bottom-up design of complex biological systems, minimizing trial-and-error approaches. The majority of models used so far adopt empirical Hill functions to describe activation and repression in exogenously-controlled inducible promoter systems. However, such equations may be poorly predictive in practical situations that are typical in bottom-up design, including changes in promoter copy number, regulatory protein level, and cell load. In this work, we derived novel mechanistic steady-state models of the lux inducible system, used as case study, relying on different assumptions on regulatory protein (LuxR) and cognate promoter (Plux) concentrations, inducer-protein complex formation, and resource usage limitation. We demonstrated that a change in the considered model assumptions can significantly affect circuit output, and preliminary experimental data are in accordance with the simulated activation curves. We finally showed that the models are identifiable a priori (in the analytically tractable cases) and a posteriori, and we determined the specific experiments needed to parametrize them. Although a larger-scale experimental validation is required, in the future the reported models may support synthetic circuits output prediction in practical situations with unprecedented details.

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

BACKGROUND

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

RESULTS

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

CONCLUSIONS

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