Feynman-Kac formula for stochastic hybrid systems

Beyond the adiabatic limit in systems with fast environments: A τ-leaping algorithm

We propose a τ-leaping simulation algorithm for stochastic systems subject to fast environmental changes. Similar to conventional τ-leaping the algorithm proceeds in discrete time steps, but as a principal addition it captures environmental noise beyond the adiabatic limit. The key idea is to treat the input rates for the τ-leaping as (clipped) Gaussian random variables with first and second moments constructed from the environmental process. In this way, each step of the algorithm retains environmental stochasticity to subleading order in the timescale separation between system and environment. We test the algorithm on several toy examples with discrete and continuous environmental states and find good performance in the regime of fast environmental dynamics. At the same time, the algorithm requires significantly less computing time than full simulations of the combined system and environment. In this context we also discuss several methods for the simulation of stochastic population dynamics in time-varying environments with continuous states.

Occupation time of a run-and-tumble particle with resetting

We study the positive occupation time of a run-and-tumble particle (RTP) subject to stochastic resetting. Under the resetting protocol, the position of the particle is reset to the origin at a random sequence of times generated by a Poisson process with rate r. The velocity state is reset to ±v with fixed probabilities ρ_{1} and ρ_{-1}=1-ρ_{1}, where v is the speed. We exploit the fact that the moment-generating functions with and without resetting are related by a renewal equation, and the latter generating function can be calculated by solving a corresponding Feynman-Kac equation. This allows us to numerically locate in Laplace space the largest real pole of the moment-generating function with resetting, and thus derive a large deviation principle (LDP) for the occupation time probability density using the Gartner-Ellis theorem. We explore how the LDP depends on the switching rate α of the velocity state, the resetting rate r, and the probability ρ_{1}. First, we show that the corresponding LDP for a Brownian particle with resetting is recovered in the fast switching limit α→∞. We then consider the case of a finite switching rate. In particular, we investigate how a directional bias in the resetting protocol (ρ_{1}≠0.5) skews the LDP rate function so that its minimum is shifted away from the expected fractional occupation time of one-half. The degree of shift increases with r and decreases with α.

Scaling methods for accelerating kinetic Monte Carlo simulations of chemical reaction networks

Various kinetic Monte Carlo algorithms become inefficient when some of the population sizes in a system are large, which gives rise to a large number of reaction events per unit time. Here, we present a new acceleration algorithm based on adaptive and heterogeneous scaling of reaction rates and stoichiometric coefficients. The algorithm is conceptually related to the commonly used idea of accelerating a stochastic simulation by considering a subvolume λΩ (0 < λ < 1) within a system of interest, which reduces the number of reaction events per unit time occurring in a simulation by a factor 1/λ at the cost of greater error in unbiased estimates of first moments and biased overestimates of second moments. Our new approach offers two unique benefits. First, scaling is adaptive and heterogeneous, which eliminates the pitfall of overaggressive scaling. Second, there is no need for an a priori classification of populations as discrete or continuous (as in a hybrid method), which is problematic when discreteness of a chemical species changes during a simulation. The method requires specification of only a single algorithmic parameter, N_c, a global critical population size above which populations are effectively scaled down to increase simulation efficiency. The method, which we term partial scaling, is implemented in the open-source BioNetGen software package. We demonstrate that partial scaling can significantly accelerate simulations without significant loss of accuracy for several published models of biological systems. These models characterize activation of the mitogen-activated protein kinase ERK, prion protein aggregation, and T-cell receptor signaling.

Efficient analysis of stochastic gene dynamics in the non-adiabatic regime using piecewise deterministic Markov processes

Single-cell experiments show that gene expression is stochastic and bursty, a feature that can emerge from slow switching between promoter states with different activities. In addition to slow chromatin and/or DNA looping dynamics, one source of long-lived promoter states is the slow binding and unbinding kinetics of transcription factors to promoters, i.e. the non-adiabatic binding regime. Here, we introduce a simple analytical framework, known as a piecewise deterministic Markov process (PDMP), that accurately describes the stochastic dynamics of gene expression in the non-adiabatic regime. We illustrate the utility of the PDMP on a non-trivial dynamical system by analysing the properties of a titration-based oscillator in the non-adiabatic limit. We first show how to transform the underlying chemical master equation into a PDMP where the slow transitions between promoter states are stochastic, but whose rates depend upon the faster deterministic dynamics of the transcription factors regulated by these promoters. We show that the PDMP accurately describes the observed periods of stochastic cycles in activator and repressor-based titration oscillators. We then generalize our PDMP analysis to more complicated versions of titration-based oscillators to explain how multiple binding sites lengthen the period and improve coherence. Last, we show how noise-induced oscillation previously observed in a titration-based oscillator arises from non-adiabatic and discrete binding events at the promoter site.

Classical stochastic systems with fast-switching environments: Reduced master equations, their interpretation, and limits of validity

We study classical Markovian stochastic systems with discrete states, coupled to randomly switching external environments. For fast environmental processes we derive reduced dynamics for the system itself, focusing on corrections to the adiabatic limit of infinite timescale separation. We show that this can lead to master equations with bursting events. Negative transition rates can result in the reduced master equation, leading to unphysical short-time behavior. However, the reduced master equation can describe stationary states better than a leading-order adiabatic calculation, similar to what is known for Kramers-Moyal expansions in the context of the Pawula theorem [R. F. Pawula, Phys. Rev. 162, 186 (1967)PHRVAO0031-899X10.1103/PhysRev.162.186; H. Risken and H. Vollmer, Z. Phys. B 35, 313 (1979)ZPBBDJ0340-224X10.1007/BF01319854]. We provide an interpretation of the reduced dynamics in discrete time and a criterion for the occurrence of negative rates for systems with two environmental states.

Model reduction methods for population dynamics with fast-switching environments: Reduced master equations, stochastic differential equations, and applications

We study stochastic population dynamics coupled to fast external environments and combine expansions in the inverse switching rate of the environment and a Kramers-Moyal expansion in the inverse size of the population. This leads to a series of approximation schemes, capturing both intrinsic and environmental noise. These methods provide a means of efficient simulation and we show how they can be used to obtain analytical results for the fluctuations of population dynamics in switching environments. We place the approximations in relation to existing work on piecewise-deterministic and piecewise-diffusive Markov processes. Finally, we demonstrate the accuracy and efficiency of these model-reduction methods in different research fields, including systems in biology and a model of crack propagation.