Steady State Sensitivity Analysis of Continuous Time Markov Chains

In this paper we study Monte Carlo estimators based on the likelihood ratio approach for steady-state sensitivity. We first extend the result of Glynn and Olvera-Cravioto [doi:doi: 10.1287/stsy.2018.002] to the setting of continuous time Markov chains with a countable state space which include models such as stochastic reaction kinetics and kinetic Monte Carlo lattice system. We show that the variance of the centered likelihood ratio estimators does not grow in time. This result suggests that the centered likelihood ratio estimators should be favored for sensitivity analysis when the mixing time of the underlying continuous time Markov chain is large, which is typically the case when systems exhibit multi-scale behavior. We demonstrate a practical implication of this analysis on a numerical benchmark of two examples for the biochemical reaction networks

Sensitivity analysis for multiscale stochastic reaction networks using hybrid approximations

We consider the problem of estimating parameter sensitivities for stochastic models of multiscale reaction networks. These sensitivity values are important for model analysis, and, the methods that currently exist for sensitivity estimation mostly rely on simulations of the stochastic dynamics. This is problematic because these simulations become computationally infeasible for multiscale networks due to reactions firing at several different timescales. However it is often possible to exploit the multiscale property to derive a "model reduction" and approximate the dynamics as a Piecewise Deterministic Markov process (PDMP), which is a hybrid process consisting of both discrete and continuous components. The aim of this paper is to show that such PDMP approximations can be used to accurately and efficiently estimate the parameter sensitivity for the original multiscale stochastic model. We prove the convergence of the original sensitivity to the corresponding PDMP sensitivity, in the limit where the PDMP approximation becomes exact. Moreover we establish a representation of the PDMP parameter sensitivity that separates the contributions of discrete and continuous components in the dynamics, and allows one to efficiently estimate both contributions.Comment: 34 pages, 4 figure

Sensitivity analysis for stochastic chemical reaction networks with multiple time-scales

Stochastic models for chemical reaction networks have become very popular in recent years. For such models, the estimation of parameter sensitivities is an important and challenging problem. Sensitivity values help in analyzing the network, understanding its robustness properties and also in identifying the key reactions for a given outcome. Most of the methods that exist in the literature for the estimation of parameter sensitivities, rely on Monte Carlo simulations using Gillespie's stochastic simulation algorithm or its variants. It is well-known that such simulation methods can be prohibitively expensive when the network contains reactions firing at different time-scales, which is a feature of many important biochemical networks. For such networks, it is often possible to exploit the time-scale separation and approximately capture the original dynamics by simulating a "reduced" model, which is obtained by eliminating the fast reactions in a certain way. The aim of this paper is to tie these model reduction techniques with sensitivity analysis. We prove that under some conditions, the sensitivity values of the reduced model can be used to approximately recover the sensitivity values for the original model. Through an example we illustrate how our result can help in sharply reducing the computational costs for the estimation of parameter sensitivities for reaction networks with multiple time-scales. To prove our result, we use coupling arguments based on the random time change representation of Kurtz. We also exploit certain connections between the distributions of the occupation times of Markov chains and multi-dimensional wave equations

Parallel replica dynamics method for bistable stochastic reaction networks: simulation and sensitivity analysis

Stochastic reaction networks that exhibit bistability are common in many fields such as systems biology and materials science. Sampling of the stationary distribution is crucial for understanding and characterizing the long term dynamics of bistable stochastic dynamical systems. However, this is normally hindered by the insufficient sampling of the rare transitions between the two metastable regions. In this paper, we apply the parallel replica (ParRep) method for continuous time Markov chain to accelerate the stationary distribution sampling of bistable stochastic reaction networks. The proposed method uses parallel computing to accelerate the sampling of rare transitions and it is very easy to implement. We combine ParRep with the path space information bounds for parametric sensitivity analysis. We demonstrate the efficiency and accuracy of the method by studying the Schl\"{o}gl model and the genetic switches network.Comment: 7 figure

Unbiased estimation of parameter sensitivities for stochastic chemical reaction networks

Estimation of parameter sensitivities for stochastic chemical reaction networks is an important and challenging problem. Sensitivity values are important in the analysis, modeling and design of chemical networks. They help in understanding the robustness properties of the system and also in identifying the key reactions for a given outcome. In a discrete setting, most of the methods that exist in the literature for the estimation of parameter sensitivities rely on Monte Carlo simulations along with finite difference computations. However these methods introduce a bias in the sensitivity estimate and in most cases the size or direction of the bias remains unknown, potentially damaging the accuracy of the analysis. In this paper, we use the random time change representation of Kurtz to derive an exact formula for parameter sensitivity. This formula allows us to construct an unbiased estimator for parameter sensitivity, which can be efficiently evaluated using a suitably devised Monte Carlo scheme. The existing literature contains only one method to produce such an unbiased estimator. This method was proposed by Plyasunov and Arkin and it is based on the Girsanov measure transformation. By taking a couple of examples we compare our method to this existing method. Our results indicate that our method can be much faster than the existing method while computing sensitivity with respect to a reaction rate constant which is small in magnitude. This rate constant could correspond to a reaction which is slow in the reference time-scale of the system. Since many biological systems have such slow reactions, our method can be a useful tool for sensitivity analysis.Comment: Submitte

An efficient and unbiased method for sensitivity analysis of stochastic reaction networks

We consider the problem of estimating parameter sensitivity for Markovian models of reaction networks. Sensitivity values measure the responsiveness of an output to the model parameters. They help in analyzing the network, understanding its robustness properties and identifying the important reactions for a specific output. Sensitivity values are commonly estimated using methods that perform finite-difference computations along with Monte Carlo simulations of the reaction dynamics. These methods are computationally efficient and easy to implement, but they produce a biased estimate which can be unreliable for certain applications. Moreover the size of the bias is generally unknown and hence the accuracy of these methods cannot be easily determined. There also exist unbiased schemes for sensitivity estimation but these schemes can be computationally infeasible, even for simple networks. Our goal in this paper is to present a new method for sensitivity estimation, which combines the computational efficiency of finite-difference methods with the accuracy of unbiased schemes. Our method is easy to implement and it relies on an exact representation of parameter sensitivity that we recently proved in an earlier paper. Through examples we demonstrate that the proposed method can outperform the existing methods, both biased and unbiased, in many situations

Efficient estimators for likelihood ratio sensitivity indices of complex stochastic dynamics

We demonstrate that centered likelihood ratio estimators for the sensitivity indices of complex stochastic dynamics are highly efficient with low, constant in time variance and consequently they are suitable for sensitivity analysis in long-time and steady-state regimes. These estimators rely on a new covariance formulation of the likelihood ratio that includes as a submatrix a Fisher Information Matrix for stochastic dynamics and can also be used for fast screening of insensitive parameters and parameter combinations. The proposed methods are applicable to broad classes of stochastic dynamics such as chemical reaction networks, Langevin-type equations and stochastic models in finance, including systems with a high dimensional parameter space and/or disparate decorrelation times between different observables. Furthermore, they are simple to implement as a standard observable in any existing simulation algorithms without additional modifications.Comment: Revision of the paper. Added a new estimato

Path-space information bounds for uncertainty quantification and sensitivity analysis of stochastic dynamics

Uncertainty quantification is a primary challenge for reliable modeling and simulation of complex stochastic dynamics. Such problems are typically plagued with incomplete information that may enter as uncertainty in the model parameters, or even in the model itself. Furthermore, due to their dynamic nature, we need to assess the impact of these uncertainties on the transient and long-time behavior of the stochastic models and derive corresponding uncertainty bounds for observables of interest. A special class of such challenges is parametric uncertainties in the model and in particular sensitivity analysis along with the corresponding sensitivity bounds for stochastic dynamics. Moreover, sensitivity analysis can be further complicated in models with a high number of parameters that render straightforward approaches, such as gradient methods, impractical. In this paper, we derive uncertainty and sensitivity bounds for path-space observables of stochastic dynamics in terms of new goal-oriented divergences; the latter incorporate both observables and information theory objects such as the relative entropy rate. These bounds are tight, depend on the variance of the particular observable and are computable through Monte Carlo simulation. In the case of sensitivity analysis, the derived sensitivity bounds rely on the path Fisher Information Matrix, hence they depend only on local dynamics and are gradient-free. These features allow for computationally efficient implementation in systems with a high number of parameters, e.g., complex reaction networks and molecular simulations.Comment: 29 pages, 2 figure

Unbiased estimation of second-order parameter sensitivities for stochastic reaction networks

This paper deals with the problem of estimating second-order parameter sensitivities for stochastic reaction networks, where the reaction dynamics is modeled as a continuous time Markov chain over a discrete state space. Estimation of such second-order sensitivities (the Hessian) is necessary for implementing the Newton-Raphson scheme for optimization over the parameter space. To perform this estimation, Wolf and Anderson have proposed an efficient finite-difference method, that uses a coupling of perturbed processes to reduce the estimator variance. The aim of this paper is to illustrate that the same coupling can be exploited to derive an exact representation for second-order parameter sensitivity. Furthermore with this representation one can construct an unbiased estimator which is easy to implement. The ideas contained in this paper are extensions of the ideas presented in our recent papers on first-order parameter sensitivity estimation.Comment: Work in progress. Any feedback from the readers is welcom

Stochastic Averaging and Sensitivity Analysis for Two Scale Reaction Networks

In the presence of multiscale dynamics in a reaction network, direct simulation methods become inefficient as they can only advance the system on the smallest scale. This work presents stochastic averaging techniques to accelerate computations for obtaining estimates of expected values and sensitivities with respect to the steady state distribution. A two-time-scale formulation is used to establish bounds on the bias induced by the averaging method. Further, this formulation provides a framework to create an accelerated `averaged' version of most single-scale sensitivity estimation method. In particular, we propose a new lower-variance ergodic likelihood ratio type estimator for steady-state estimation and show how one can adapt it to accelerated simulations of multiscale systems.Lastly, we develop an adaptive "batch-means" stopping rule for determining when to terminate the micro-equilibration process.Comment: 20 pages, 6 figures, 2 tables, in this version: corrigendum of Proposition III.1 - only convergence is established, at present no rate available for general cas