Implicit ODE solvers with good local error control for the transient analysis of Markov models

Obtaining the transient probability distribution vector of a continuous-time Markov chain (CTMC) using an implicit ordinary differential equation (ODE) solver tends to be advantageous in terms of run-time computational cost when the product of the maximum output rate of the CTMC and the largest time of interest is large. In this paper, we show that when applied to the transient analysis of CTMCs, many implicit ODE solvers are such that the linear systems involved in their steps can be solved by using iterative methods with strict control of the 1-norm of the error. This allows the development of implementations of those ODE solvers for the transient analysis of CTMCs that can be more efficient and more accurate than more standard implementations.Peer ReviewedPostprint (published version

Design and Development of Software Tools for Bio-PEPA

This paper surveys the design of software tools for the Bio-PEPA process algebra. Bio-PEPA is a high-level language for modelling biological systems such as metabolic pathways and other biochemical reaction networks. Through providing tools for this modelling language we hope to allow easier use of a range of simulators and model-checkers thereby freeing the modeller from the responsibility of developing a custom simulator for the problem of interest. Further, by providing mappings to a range of different analysis tools the Bio-PEPA language allows modellers to compare analysis results which have been computed using independent numerical analysers, which enhances the reliability and robustness of the results computed.

Data-driven modelling of biological multi-scale processes

Biological processes involve a variety of spatial and temporal scales. A holistic understanding of many biological processes therefore requires multi-scale models which capture the relevant properties on all these scales. In this manuscript we review mathematical modelling approaches used to describe the individual spatial scales and how they are integrated into holistic models. We discuss the relation between spatial and temporal scales and the implication of that on multi-scale modelling. Based upon this overview over state-of-the-art modelling approaches, we formulate key challenges in mathematical and computational modelling of biological multi-scale and multi-physics processes. In particular, we considered the availability of analysis tools for multi-scale models and model-based multi-scale data integration. We provide a compact review of methods for model-based data integration and model-based hypothesis testing. Furthermore, novel approaches and recent trends are discussed, including computation time reduction using reduced order and surrogate models, which contribute to the solution of inference problems. We conclude the manuscript by providing a few ideas for the development of tailored multi-scale inference methods.Comment: This manuscript will appear in the Journal of Coupled Systems and Multiscale Dynamics (American Scientific Publishers

Likelihood based observability analysis and confidence intervals for predictions of dynamic models

Mechanistic dynamic models of biochemical networks such as Ordinary Differential Equations (ODEs) contain unknown parameters like the reaction rate constants and the initial concentrations of the compounds. The large number of parameters as well as their nonlinear impact on the model responses hamper the determination of confidence regions for parameter estimates. At the same time, classical approaches translating the uncertainty of the parameters into confidence intervals for model predictions are hardly feasible. In this article it is shown that a so-called prediction profile likelihood yields reliable confidence intervals for model predictions, despite arbitrarily complex and high-dimensional shapes of the confidence regions for the estimated parameters. Prediction confidence intervals of the dynamic states allow a data-based observability analysis. The approach renders the issue of sampling a high-dimensional parameter space into evaluating one-dimensional prediction spaces. The method is also applicable if there are non-identifiable parameters yielding to some insufficiently specified model predictions that can be interpreted as non-observability. Moreover, a validation profile likelihood is introduced that should be applied when noisy validation experiments are to be interpreted. The properties and applicability of the prediction and validation profile likelihood approaches are demonstrated by two examples, a small and instructive ODE model describing two consecutive reactions, and a realistic ODE model for the MAP kinase signal transduction pathway. The presented general approach constitutes a concept for observability analysis and for generating reliable confidence intervals of model predictions, not only, but especially suitable for mathematical models of biological systems

Selecting Mathematical Software for Dependability Assessment of Computer Systems Described by Stiff Markov Chains

Markov and semi-Markov models are widely used in dependability assessment of complex computer-based systems. Model stiffness poses a serious problem both in terms of computational difficulties and in terms of accuracy of the assessment. Selecting an appropriate method and software package for solving stiff Markov models proved to be a non-trivial task. In this paper we provide an empirical comparison of two approaches to dealing with stiffness – stiffness avoidance and stiffness-tolerance. The study includes several well known techniques and software tools used for solving Kolmogorov’s differential equations derived from complex stiff Markov models. In the comparison we used realistic cases studies developed by others in the past: i) a computer system with hardware redundancy and diverse software, and ii) a queuing system with a server break-down and repair. The results indicate that the accuracy of the known methods is significantly affected by the stiffness of the Markov models, which led us to developing a procedure (an algorithm) for selecting the optimal method and tool for solving a given stiff Markov model. The algorithm is, also included in the paper

Hybrid Analog-Digital Co-Processing for Scientific Computation

In the past 10 years computer architecture research has moved to more heterogeneity and less adherence to conventional abstractions. Scientists and engineers hold an unshakable belief that computing holds keys to unlocking humanity's Grand Challenges. Acting on that belief they have looked deeper into computer architecture to find specialized support for their applications. Likewise, computer architects have looked deeper into circuits and devices in search of untapped performance and efficiency. The lines between computer architecture layers---applications, algorithms, architectures, microarchitectures, circuits and devices---have blurred. Against this backdrop, a menagerie of computer architectures are on the horizon, ones that forgo basic assumptions about computer hardware, and require new thinking of how such hardware supports problems and algorithms. This thesis is about revisiting hybrid analog-digital computing in support of diverse modern workloads. Hybrid computing had extensive applications in early computing history, and has been revisited for small-scale applications in embedded systems. But architectural support for using hybrid computing in modern workloads, at scale and with high accuracy solutions, has been lacking. I demonstrate solving a variety of scientific computing problems, including stochastic ODEs, partial differential equations, linear algebra, and nonlinear systems of equations, as case studies in hybrid computing. I solve these problems on a system of multiple prototype analog accelerator chips built by a team at Columbia University. On that team I made contributions toward programming the chips, building the digital interface, and validating the chips' functionality. The analog accelerator chip is intended for use in conjunction with a conventional digital host computer. The appeal and motivation for using an analog accelerator is efficiency and performance, but it comes with limitations in accuracy and problem sizes that we have to work around. The first problem is how to do problems in this unconventional computation model. Scientific computing phrases problems as differential equations and algebraic equations. Differential equations are a continuous view of the world, while algebraic equations are a discrete one. Prior work in analog computing mostly focused on differential equations; algebraic equations played a minor role in prior work in analog computing. The secret to using the analog accelerator to support modern workloads on conventional computers is that these two viewpoints are interchangeable. The algebraic equations that underlie most workloads can be solved as differential equations, and differential equations are naturally solvable in the analog accelerator chip. A hybrid analog-digital computer architecture can focus on solving linear and nonlinear algebra problems to support many workloads. The second problem is how to get accurate solutions using hybrid analog-digital computing. The reason that the analog computation model gives less accurate solutions is it gives up representing numbers as digital binary numbers, and instead uses the full range of analog voltage and current to represent real numbers. Prior work has established that encoding data in analog signals gives an energy efficiency advantage as long as the analog data precision is limited. While the analog accelerator alone may be useful for energy-constrained applications where inputs and outputs are imprecise, we are more interested in using analog in conjunction with digital for precise solutions. This thesis gives novel insight that the trick to do so is to solve nonlinear problems where low-precision guesses are useful for conventional digital algorithms. The third problem is how to solve large problems using hybrid analog-digital computing. The reason the analog computation model can't handle large problems is it gives up step-by-step discrete-time operation, instead allowing variables to evolve smoothly in continuous time. To make that happen the analog accelerator works by chaining hardware for mathematical operations end-to-end. During computation analog data flows through the hardware with no overheads in control logic and memory accesses. The downside is then the needed hardware size grows alongside problem sizes. While scientific computing researchers have for a long time split large problems into smaller subproblems to fit in digital computer constraints, this thesis is a first attempt to consider these divide-and-conquer algorithms as an essential tool in using the analog model of computation. As we enter the post-Moore’s law era of computing, unconventional architectures will offer specialized models of computation that uniquely support specific problem types. Two prominent examples are deep neural networks and quantum computers. Recent trends in computer science research show these unconventional architectures will soon have broad adoption. In this thesis I show another specialized, unconventional architecture is to use analog accelerators to solve problems in scientific computing. Computer architecture researchers will discover other important models of computation in the future. This thesis is an example of the discovery process, implementation, and evaluation of how an unconventional architecture supports specialized workloads

Probabilistic reasoning and inference for systems biology

One of the important challenges in Systems Biology is reasoning and performing hypotheses testing in uncertain conditions, when available knowledge may be incomplete and the experimental data may contain substantial noise. In this thesis we develop methods of probabilistic reasoning and inference that operate consistently within an environment of uncertain knowledge and data. Mechanistic mathematical models are used to describe hypotheses about biological systems. We consider both deductive model based reasoning and model inference from data. The main contributions are a novel modelling approach using continuous time Markov chains that enables deductive derivation of model behaviours and their properties, and the application of Bayesian inferential methods to solve the inverse problem of model inference and comparison, given uncertain knowledge and noisy data. In the first part of the thesis, we consider both individual and population based techniques for modelling biochemical pathways using continuous time Markov chains, and demonstrate why the latter is the most appropriate. We illustrate a new approach, based on symbolic intervals of concentrations, with an example portion of the ERK signalling pathway. We demonstrate that the resulting model approximates the same dynamic system as traditionally defined using ordinary differential equations. The advantage of the new approach is quantitative logical analysis; we formulate a number of biologically significant queries in the temporal logic CSL and use probabilistic symbolic model checking to investigate their veracity. In the second part of the thesis, we consider the inverse problem of model inference and testing of alternative hypotheses, when models are defined by non-linear ordinary differential equations and the experimental data is noisy and sparse. We compare and evaluate a number of statistical techniques, and implement an effective Bayesian inferential framework for systems biology based on Markov chain Monte Carlo methods and estimation of marginal likelihoods by annealing-melting integration. We illustrate the framework with two case studies, one of which involves an open problem concerning the mediation of ERK phosphorylation in the ERK pathway

Robust Stochastic Chemical Reaction Networks and Bounded Tau-Leaping

The behavior of some stochastic chemical reaction networks is largely unaffected by slight inaccuracies in reaction rates. We formalize the robustness of state probabilities to reaction rate deviations, and describe a formal connection between robustness and efficiency of simulation. Without robustness guarantees, stochastic simulation seems to require computational time proportional to the total number of reaction events. Even if the concentration (molecular count per volume) stays bounded, the number of reaction events can be linear in the duration of simulated time and total molecular count. We show that the behavior of robust systems can be predicted such that the computational work scales linearly with the duration of simulated time and concentration, and only polylogarithmically in the total molecular count. Thus our asymptotic analysis captures the dramatic speedup when molecular counts are large, and shows that for bounded concentrations the computation time is essentially invariant with molecular count. Finally, by noticing that even robust stochastic chemical reaction networks are capable of embedding complex computational problems, we argue that the linear dependence on simulated time and concentration is likely optimal

A Bayesian Approach to Modelling Biological Pattern Formation with Limited Data

Pattern formation in biological tissues plays an important role in the development of living organisms. Since the classical work of Alan Turing, a pre-eminent way of modelling has been through reaction-diffusion mechanisms. More recently, alternative models have been proposed, that link dynamics of diffusing molecular signals with tissue mechanics. In order to distinguish among different models, they should be compared to experimental observations. However, in many experimental situations only the limiting, stationary regime of the pattern formation process is observable, without knowledge of the transient behaviour or the initial state. The unstable nature of the underlying dynamics in all alternative models seriously complicates model and parameter identification, since small changes in the initial condition lead to distinct stationary patterns. To overcome this problem the initial state of the model can be randomised. In the latter case, fixed values of the model parameters correspond to a family of patterns rather than a fixed stationary solution, and standard approaches to compare pattern data directly with model outputs, e.g., in the least squares sense, are not suitable. Instead, statistical characteristics of the patterns should be compared, which is difficult given the typically limited amount of available data in practical applications. To deal with this problem, we extend a recently developed statistical approach for parameter identification using pattern data, the so-called Correlation Integral Likelihood (CIL) method. We suggest modifications that allow increasing the accuracy of the identification process without resizing the data set. The proposed approach is tested using different classes of pattern formation models. For all considered equations, parallel GPU-based implementations of the numerical solvers with efficient time stepping schemes are provided.Comment: More compact version of the text and figures, results unchange