K-leap method for accelerating stochastic simulation of coupled chemical reactions

Leap methods are very promising for accelerating stochastic simulation of a well stirred chemically reacting system, while providing acceptable simulation accuracy. In Gillespie's tau-leap method [D. Gillespie, J. Phys. Chem. 115, 1716 (2001)], the number of firings of each reaction channel during a leap is a Poisson random variable, whose sample values are unbounded. This may cause large changes in the populations of certain molecular species during a leap, thereby violating the leap condition. In this paper, we develop an alternative leap method called the K-leap method, in which we constrain the total number of reactions occurring during a leap to be a number K calculated from the leap condition. As the number of firings of each reaction channel during a leap is upper bounded by a properly chosen number, our K-leap method can better satisfy the leap condition, thereby improving simulation accuracy. Since the exact stochastic simulation algorithm (SSA) is a special case of our K-leap method when K=1, our K-leap method can naturally change from the exact SSA to an approximate leap method during simulation, whenever the leap condition allows to do so.     

Accelerated stochastic simulation algorithm for coupled chemical reactions with delays

Some biochemical processes do not occur instantaneously but have considerably delays associated with them. In the existed methods which solve these chemically reacting systems with delays, averaging over a great deal of simulations is needed for reliable statistical characters. Here we present an accelerating approach, called the "Delay Final All Possible Steps" (DFAPS) approach, which does not alter the course of stochastic simulation, but reduces the required running times. Numerical simulation results indicate that the proposed method can be applied to a wide range of chemically reacting systems with delays and obtain a significant improvement on efficiency and accuracy over the existed methods.     

Binomial tau-leap spatial stochastic simulation algorithm for applications in chemical kinetics

In cell biology, cell signaling pathway problems are often tackled with deterministic temporal models, well mixed stochastic simulators, and/or hybrid methods. But, in fact, three dimensional stochastic spatial modeling of reactions happening inside the cell is needed in order to fully understand these cell signaling pathways. This is because noise effects, low molecular concentrations, and spatial heterogeneity can all affect the cellular dynamics. However, there are ways in which important effects can be accounted without going to the extent of using highly resolved spatial simulators (such as single-particle software), hence reducing the overall computation time significantly. We present a new coarse grained modified version of the next subvolume method that allows the user to consider both diffusion and reaction events in relatively long simulation time spans as compared with the original method and other commonly used fully stochastic computational methods. Benchmarking of the simulation algorithm was performed through comparison with the next subvolume method and well mixed models (MATLAB), as well as stochastic particle reaction and transport simulations (CHEMCELL, Sandia National Laboratories). Additionally, we construct a model based on a set of chemical reactions in the epidermal growth factor receptor pathway. For this particular application and a bistable chemical system example, we analyze and outline the advantages of our presented binomial tau-leap spatial stochastic simulation algorithm, in terms of efficiency and accuracy, in scenarios of both molecular homogeneity and heterogeneity.     

Accurate hybrid stochastic simulation of a system of coupled chemical or biochemical reactions

The dynamical solution of a well-mixed, nonlinear stochastic chemical kinetic system, described by the Master equation, may be exactly computed using the stochastic simulation algorithm. However, because the computational cost scales with the number of reaction occurrences, systems with one or more "fast" reactions become costly to simulate. This paper describes a hybrid stochastic method that partitions the system into subsets of fast and slow reactions, approximates the fast reactions as a continuous Markov process, using a chemical Langevin equation, and accurately describes the slow dynamics using the integral form of the "Next Reaction" variant of the stochastic simulation algorithm. The key innovation of this method is its mechanism of efficiently monitoring the occurrences of slow, discrete events while simultaneously simulating the dynamics of a continuous, stochastic or deterministic process. In addition, by introducing an approximation in which multiple slow reactions may occur within a time step of the numerical integration of the chemical Langevin equation, the hybrid stochastic method performs much faster with only a marginal decrease in accuracy. Multiple examples, including a biological pulse generator and a large-scale system benchmark, are simulated using the exact and proposed hybrid methods as well as, for comparison, a previous hybrid stochastic method. Probability distributions of the solutions are compared and the weak errors of the first two moments are computed. In general, these hybrid methods may be applied to the simulation of the dynamics of a system described by stochastic differential, ordinary differential, and Master equations.     

A topological stochastic approach to the study of multidimensional potential energy surfaces of chemical reactions

In this article we propose a new approach for investigating the properties of multidimensional potential energy surfaces in chemical reactions, based on relations of each multidimensional surface to its one-dimensional image which is the chemical reaction tree. This approach makes it possible to find a common number of independent channels in chemical reactions for complex systems and to construct the probable models.     

Domain catalyzed chemical reactions: a molecular dynamics simulation of isomerization kinetics

The model for domain catalyzed isomerization kinetics in condensed fluids is applied for a diluted mixture of a chiral solute with a consolute temperature. The solution is quench to phase separation at temperatures below the consolute temperature. The droplet coalescence enhances the isomerization kinetics due to the substantial excess pressure inside the small droplets given by the Laplace equation. The domain catalyzed isomerization kinetics breaks the symmetry, and the droplets end with only one dominating species. We argue that D-glyceraldehyde which is only moderately solvable in water and which has played a crucial role in the evolution is a candidate for the stereo specific ordering in bio-organic matter.     

Implementation of a molecular dynamics approach with rigid fragments to simulation of chemical reactions in biomolecular systems

We describe a new implementation of the molecular dynamics method aimed at simulation of the properties of biomolecular systems in which chemical reactions are possible. The quantum mechanical/molecular mechanical method based on the effective fragment potential theory is used for calculating the energies and forces along trajectories. Due to specific features of the effective fragment theory, the behavior of the molecular mechanical subsystem is described by rigid body dynamics. The method has been applied to simulation of proton transfer along the chain of water molecules inside the gramicidin channel.     

Nested stochastic simulation algorithm for chemical kinetic systems with disparate rates

An efficient simulation algorithm for chemical kinetic systems with disparate rates is proposed. This new algorithm is quite general, and it amounts to a simple and seamless modification of the classical stochastic simulation algorithm (SSA), also known as the Gillespie [J. Comput. Phys. 22, 403 (1976); J. Phys. Chem. 81, 2340 (1977)] algorithm. The basic idea is to use an outer SSA to simulate the slow processes with rates computed from an inner SSA which simulates the fast reactions. Averaging theorems for Markov processes can be used to identify the fast and slow variables in the system as well as the effective dynamics over the slow time scale, even though the algorithm itself does not rely on such information. This nested SSA can be easily generalized to systems with more than two separated time scales. Convergence and efficiency of the algorithm are discussed using the established error estimates and illustrated through examples.     

Digital simulation of electrochemical processes involving very fast chemical reactions : Part 1. A simple general approach

A simple approach based on the steady-state assumption and linear concentration profiles in the reaction layer is proposed to simulate electrochemical systems involving very fast homogeneous chemical reactions. The flux equations are derived in detail and the method is shown to be suitable for catalytic e.c.e and dimerization reaction mechanism. The accuracy of the proposed method is checked by comparing the results obtained with those already available in the literature.     

An analytical approach to solutions of master equations for stochastic nonlinear reactions

In this paper we consider a class of nonlinear reactions which are important in stochastic reaction networks. We find the exact solution of the chemical master equation for a class of irreversible and reversible nonlinear reactions. We also present the explicit form of the equilibrium probability solution of the reactions. The results can be used for analyzing stochastic dynamics of important reactions such as binding/unbinding reaction and protein dimerization.     

A partial-propensity formulation of the stochastic simulation algorithm for chemical reaction networks with delays

Several real-world systems, such as gene expression networks in biological cells, contain coupled chemical reactions with a time delay between reaction initiation and completion. The non-Markovian kinetics of such reaction networks can be exactly simulated using the delay stochastic simulation algorithm (dSSA). The computational cost of dSSA scales with the total number of reactions in the network. We reduce this cost to scale at most with the smaller number of species by using the concept of partial reaction propensities. The resulting delay partial-propensity direct method (dPDM) is an exact dSSA formulation for well-stirred systems of coupled chemical reactions with delays. We detail dPDM and present a theoretical analysis of its computational cost. Furthermore, we demonstrate the implications of the theoretical cost analysis in two prototypical benchmark applications. The dPDM formulation is shown to be particularly efficient for strongly coupled reaction networks, where the number of reactions is much larger than the number of species.     

Kinetic approach to chemical reactions and inelastic transitions in a rarefied gas

A kinetic approach is presented for the analysis of a gas mixture with two kinds of nonconservative interactions. In a bimolecular chemical reaction, mass transfer and energy of chemical link arise, and in inelastic mechanical encounters, molecules get excited or de‐excited due to their quantized structure. Molecules undergo transitions between energy levels also by absorption and emission of photons of the self‐consistent radiation field. From the kinetic Boltzmann‐type equations, the problem of equilibria and of their stability is addressed. A detailed balance principle is proved and a Lyapunov functional is constructed; mass action law and Planck's law of radiation are recovered. This revised version was published online in July 2006 with corrections to the Cover Date.     

A constrained approach to multiscale stochastic simulation of chemically reacting systems

Stochastic simulation of coupled chemical reactions is often computationally intensive, especially if a chemical system contains reactions occurring on different time scales. In this paper, we introduce a multiscale methodology suitable to address this problem, assuming that the evolution of the slow species in the system is well approximated by a Langevin process. It is based on the conditional stochastic simulation algorithm (CSSA) which samples from the conditional distribution of the suitably defined fast variables, given values for the slow variables. In the constrained multiscale algorithm (CMA) a single realization of the CSSA is then used for each value of the slow variable to approximate the effective drift and diffusion terms, in a similar manner to the constrained mean-force computations in other applications such as molecular dynamics. We then show how using the ensuing Fokker-Planck equation approximation, we can in turn approximate average switching times in stochastic chemical systems.     

Digital simulation of electrochemical processes involving very fast chemical reactions : Part 2. A new criterion for determining rate constants of consecutive second-order irreversible chemical reactions

Electrode reactions followed by very fast chemical reactions are considered. A simple approach, in which steady state and linear concentration profiles in the reaction layer are assumed, is proposed for the simulation of these processes when the substrate is not present in large excess. When the substrate/depolarizer mole ratio is less than one, two well-separated peaks are detected; under such conditions, working curves that enable the corresponding second-order homogeneous rate constant to be evaluated can be derived by the finite-difference simulation technique. The method is applied to elucidation of the reduction of the tetrakis(triphenylphosphine)nickel(I) complex in the presence of hydrogen ions in acetonitrile at ?30°C.     

Mapping tertiary interactions in protein folding reactions: a novel mass spectrometry- and chemical synthesis-based approach

A novel mass spectrometry- and chemical synthesis-based approach for studying protein folding reactions is described, and its initial application to study the folding/unfolding reaction of a homo-hexameric enzyme 4-oxalocrotonate (4OT) is reported. This new approach involves the application of total chemical synthesis to prepare protein analogues that contain a photoreactive amino acid site-specifically incorporated into their primary amino acid sequence. To this end, a photoreactive amino acid-containing analogue of 4OT in which Pro-1 was replaced with p-benzoyl-l-phenylalanine (Bpa) was prepared. This analogue can be used to map structurally specific protein-protein interactions in 4OT’s native folded state. These photocrosslinking studies and peptide mapping results with (P1Bpa)4OT indicate that this construct is potentially useful for probing the structural properties of equilibrium and kinetic intermediates in 4OT’s folding reaction.     

Mapping tertiary interactions in protein folding reactions: a novel mass spectrometry- and chemical synthesis-based approach

A novel mass spectrometry- and chemical synthesis-based approach for studying protein folding reactions is described, and its initial application to study the folding/unfolding reaction of a homo-hexameric enzyme 4-oxalocrotonate (4OT) is reported. This new approach involves the application of total chemical synthesis to prepare protein analogues that contain a photoreactive amino acid site-specifically incorporated into their primary amino acid sequence. To this end, a photoreactive amino acid-containing analogue of 4OT in which Pro-1 was replaced with p-benzoyl-l-phenylalanine (Bpa) was prepared. This analogue can be used to map structurally specific protein-protein interactions in 4OT's native folded state. These photocrosslinking studies and peptide mapping results with (PlBpa)4OT indicate that this construct is potentially useful for probing the structural properties of equilibrium and kinetic intermediates in 4OT's folding reaction.     

Cycloaddition reactions: a controlled approach for carbon nanotube functionalization

Controlled functionalization of carbon nanotubes (CNTs) through the use of cycloaddition reactions is described. By employing various cycloaddition reactions, a wide range of molecules could be coupled onto CNTs without disruption of the structural integrity as well as with a statistical distribution of functional groups onto the surface of the CNTs. The cycloaddition reactions represent an effective and tailored approach for preparing CNT-based advanced hybrid materials that would be useful for a wide range of applications from nanobiotechnology to nanoelectronics.     

A stepwise simulation approach to the chemical potential of fluids

A new approach to the computation of the chemical potential of fluids is presented. In this method the particle-insertion operation in the conventional test particle method is replaced by the growth of a specific particle. Application of the new technique to hard sphere and Lennard-Jones fluids shows that it is capable of providing reliable estimates of the chemical potential, even at high density where the conventional test particle methods are difficult to apply.