The Stretch to Stray on Time: Resonant Length of Random Walks in a Transient
First-passage times in random walks have a vast number of diverse applications in physics, chemistry, biology, and finance. In general, environmental conditions for a stochastic process are not constant on the time scale of the average first-passage time, or control might be applied to reduce noise. We investigate moments of the first-passage time distribution under a transient describing relaxation of environmental conditions. We solve the Laplace-transformed (generalized) master equation analytically using a novel method that is applicable to general state schemes. The first-passage time from one end to the other of a linear chain of states is our application for the solutions. The dependence of its average on the relaxation rate obeys a power law for slow transients. The exponent depends on the chain length like to leading order. Slow transients substantially reduce the noise of first-passage times expressed as the coefficient of variation (CV), even if the average first-passage time is much longer than the transient. The CV has a pronounced minimum for some lengths, which we call resonant lengths. These results also suggest a simple and efficient noise control strategy, and are closely related to the timing of repetitive excitations, coherence resonance and information transmission by noisy excitable systems. A resonant number of steps from the inhibited state to the excitation threshold and slow recovery from negative feedback provide optimal timing noise reduction and information transmission.
Continuous-time random walks are a unifying concept across physics Wei et al. (2000); Chupeau et al. (2015); Guérin et al. (2016); Molini et al. (2011); Panja (2010); Turiv et al. (2013); Démery et al. (2014); Ernst et al. (2012); Baranovskii (2014); Metzler and Klafter (2000); Blumen et al. (2002); Hofmann and Ivanyuk (2003); Sokolov and Klafter (2005); Shlesinger (2006), chemistry van Kampen (2001); Oppenheim et al. (1977); Bénichou et al. (2010), biology Berg (1993); Gazzola et al. (2009); Bressloff and Newby (2013); Ivanov et al. (2016), and finance Scott (1997); Steele (2000). The drunkard’s straying on his way home from the pub is the graphic example frequently used to illustrate the randomness of step timing and direction. In particular, the time of first passage of a specific state given the stochastic nature of the process is of interest in many applications. It describes the time the drunkard arrives home, the time necessary for a chemical reaction or gene expression to reach a certain number of product molecules van Kampen (2001); Ivanov et al. (2016), or a quantity in one of the many other applications Chupeau et al. (2015); Guérin et al. (2016); Hofmann and Ivanyuk (2003); Bénichou et al. (2010); Mayo et al. (2011); Winter and von Hippel (1981); Berg et al. (1981); Riggs et al. (1970a, b); Shimamoto (1999); Albert et al. (2012); Scher and Lax (1973); Britton (2010); Redner et al. (2014); Redner (2001); Condamin et al. (2007); Shlesinger (2007); Hormoz et al. (2016).
While noise happens on the small length scales and short time scales of a system, it may trigger events on a global scale. One of the most important functions of noise for macroscopic dynamics arises from its combination with thresholds GarcÃa-Ojalvo and Schimansky-Geier (2000); Kostur et al. (2003); Lindner et al. (2004). These are defined by the observation that the dynamics of a system remains close to a stable state as long as it does not cross the threshold value, and an actively amplified deviation from this state happens when the threshold is crossed. Noise drives the system across the threshold in a random manner. First passage is a natural concept to describe the timing of threshold crossings. Ignition processes are an illustrative textbook example. Although a small random spark might not be capable of igniting an inflammable material, a few of them might cause an explosion or forest fire. If the system again attains its stable state upon recovery from a deviation, such behavior is called excitable and the large deviation an excitation. The excitation is terminated by negative feedback. A forest is excitable, because it regrows after a fire. Consumption of inflammable trees acts as the negative feedback. Excitability describes not only forest fires but also the dynamics of heterogeneous catalysis Falcke et al. (1992), the firing of neurons Keener and Sneyd (1998), the properties of heart muscle tissue Keener and Sneyd (1998), and many other systems in physics, chemistry, and biology Scott (1975); Pei et al. (1995); Cross and Hohenberg (1993); Mikhailov (1994); Ebeling and Feistel (1986).
Random walks are frequently defined on a set of discrete states. The rates or waiting-time distributions for transitions from state to set the state dwell times. The first-passage time between two widely separated states is much longer than the individual dwell times, and the conditions setting the rates and parameters of the are likely to change or external control acts on the system between start and first passage. The conditions for igniting a forest fire change with the seasons or because the forest recovers from a previous fire. The occurrence of subcritical sparks during recovery has essentially no affect on the process of regrowth. More generally, noise does not affect recovery on large space and long time scales, and the random process experiences recovery as a slow deterministic change of environmental conditions. Since recovery is typically a slow relaxation process Cross and Hohenberg (1993); Mikhailov (1994); Ebeling and Feistel (1986), it dominates event timing. Hence, first passage in an exponential transient is a natural concept through which to understand the timing of sequences of excitations. We will investigate it in this study.
We will take Markovian processes as one of the asymptotic cases of transient relaxation. There are several reasons for also considering non-Markovian formulations of continuous-time random walks. A description of a diffusing particle becomes non-Markovian whenever the particle not only passively experiences the thermal fluctuations causing its diffusive motion, but also acts back on its immediate surroundings Wei et al. (2000); Guérin et al. (2016); Franosch et al. (2011); Panja (2010); Turiv et al. (2013); Démery et al. (2014); Ernst et al. (2012). Non-Markovian waiting-time distributions for state transitions arise naturally in transport theory Esposito and Lindenberg (2008); Albert et al. (2012, 2016); Koch et al. (2005). Frequently in biological applications, the discrete states that arise are lumped states consisting of many “microscopic” states. Groups of open and closed states of ion channels lumped into single states are examples of this Hille (2001); Sneyd et al. (2004). Transitions between lumped states are non-Markovian owing to the internal dynamics. We may also use waiting-time distributions if we lack information on all the individual steps of a process, but we do know the inter-event interval distributions. This is usually the case with in vivo measurements such as stimulation of a cell and the appearance of a response, or differentiation sequences of stem cells Hormoz et al. (2016). The state probabilities of non-Markovian processes obey generalized master equations, which we will use here Oppenheim et al. (1977); Kenkre (1978); Hänggi and Talkner (1985); Masoliver et al. (1986); Metzler and Klafter (2000); Metzler et al. (1999).
In Sec. II, we present a formulation of the general problem in terms of the normal and generalized master equations and give analytic solutions for both of these. These solutions apply to general state schemes. We continue with investigating first passage on linear chains of states in Sec. III. We present results on scaling of the average first-passage time with the relaxation rate of the transient and the chain length in Sec. IV, and results on the phenomenon of resonant lengths in Sec. V.
Ii Basic equations
ii.1 The asymptotically Markovian master equation
In this section, we consider transition rates relaxing with rate to an asymptotic value like
They reach a Markov process asymptotically. The dynamics of the probability to be in state for a process that started in at obey the master equation
In matrix notation with the vector of probabilities , we have
with the matrices and defined by Eq. (2). The initial condition defines the vector . The Laplace transform of the master equation allows for a comfortable calculation of moments of the first-passage times, which we will carry out in Sec. III. The Laplace transform of Eq. (3) is the system of linear difference equations
ii.2 The generalized master equation
The waiting-time distributions
Waiting-time distributions in a constant environment depend on the time elapsed since the process entered state at time . The change in conditions causes an additional dependence on : . The lumping of states, which we introduced as a major cause of dwell-time-dependent transition probabilities, often entails , with a maximum of at intermediate values of . Waiting-time distributions used in transport theory exhibit similar properties Davies et al. (1992); Esposito and Lindenberg (2008); Albert et al. (2012, 2016). We use a simple realization of this type of distributions by a biexponential function in :
The transient parts of collect all factors of in Eq. (5). The functions describe the asymptotic part of the waiting-time distributions remaining after the transient. The are normalized to the splitting probabilities (the total probability for a transition from to given the system entered at ). They satisfy
where is the number of transitions out of state .
The generalized master equation and its Laplace transform
In the non-Markovian case, the dynamics of the probabilities obey a generalized master equation Kenkre et al. (1973); Montroll and Scher (1973); Burshtein et al. (1986); Chechkin et al. (2005); Sokolov and Klafter (2007)
where is the probability flux due to transitions from state to given that the process started at state at , and is the number of transitions toward . The fluxes are the solutions of the integral equation
The second factor in the convolution is the probability of arriving in state at time , and the first factor is the probability of leaving toward at time given arrival at . The are the initial fluxes, with for .
The Laplace transform of the probability dynamics equation (8) is
which contains the Laplace-transformed probability fluxes . The Kronecker delta captures the initial condition.
Laplace-transforming Eq. (9) is straightforward for terms containing the asymptotic part of the [see Eq. (6)], since they depend on only and the convolution theorem applies directly. The terms containing the transient part depend on both and and require a little more attention:
This leads to the Laplace transform of Eq. (9):
We write the and as the vectors and , respectively, to obtain
Solving for results in
This is again a system of linear difference equations. The entries in the matrices and are the functions and , which are the Laplace transforms of and [Eq. (6)].
Solving the Laplace-transformed generalized and asymptotically Markovian master equations
All elements of [Eq. (14)] vanish for . We also expect to be bounded for . Hence, the solution for is
Consequently, for 1,
holds for the solution for all values of and for a natural number. Once we know , we can use Eq. (14) to find :
as the solution of Eq. (14). Equation (17) is confirmed by verifying that it provides the correct solution in the limiting cases and . In the latter, all entries of the matrix vanish, and we obtain directly the result without transient, Eq. (15). For , we notice that
holds, which leads to the correct solution of Eq. (14):
This equation has the same structure as Eq. (14), and since the matrix also vanishes for , we can calculate the Laplace transform of the completely analogously to the non-Markovian case. We define and obtain
as the solution of Eq. (19).
These solutions for the Laplace transforms of the generalized and asymptotically Markovian master equations with a transient, Eqs. (14) and (20), are not restricted to state-independent waiting-time distributions or rates. They also apply to random walks with space- or state-dependent waiting-time distributions and to arbitrary state networks.
Iii The probability density of the first-passage time for a linear chain of states
A large class of stochastic processes, including all the examples mentioned in the introduction, are represented by state schemes like
which we consider from now on. The definition of the transition probabilities for both the asymptotically Markovian system and the non-Markovian system completes its specification. We use the asymptotically Markovian rates
corresponding to , , and in Eq. (1). The process has a strong initial bias for motion toward 0, . It relaxes with rate to a symmetric random walk with .
We specify the input for the generalized master equation of the non-Markovian system as
The denominators in Eqs. (23) and (24) arise from the normalization to =1 [Eq. (7) with ]. We use identical waiting-time distributions for all transitions () and all transitions . The process defined by Eqs. (23) and (24) also has a strong bias toward 0 in the beginning, and is approximately symmetric in the limit . Figure 1 illustrates the and their development with increasing . The limit of illustrates the consequences of the transient. The relax from 0 to an asymptotic value with rate , but the dwell times of individual states () are essentially constant (Fig. 1).
The transient here causes a decrease of the dwell time in state 0 (see Fig. 1).
The first-passage-time probability density provides the probability of arrival for the first time in state in when the process started in state at =0. It is given by the probability flux out of the state range from 0 to :
We denote its Laplace transform by . The moments of the first-passage-time distribution are given by van Kampen (2001)
can be determined by solving the master equations, setting state 0 as the initial condition and considering the state as absorbing, i.e., :
captures not only the first-passage time , but also, to a good approximation, transitions starting at states with indices larger than 0 (and ), since, owing to the initial bias, the process quickly moves into state 0 first and then slowly starts from there.
With the non-Markovian waiting-time distributions, the matrices and and the vector of Eq. (17) specific to this problem are
All other entries are equal to 0. Specifically for the first-passage problem, the matrices and are
with all other entries being 0. The vector is equal to , . The Laplace transforms of the first-passage-time distributions are
with Eq. (8) and
with Eq. (2). Figure 2 (a) and (b) compare analytical results for the average first-passage time with the results of simulations. The agreement is very good, thus confirming the solutions given by Eqs. (17) and (20). This confirmation by simulations is important, since there is no method of solving difference equations that guarantees a complete solution.
Iv Scaling of the average first-passage time with the relaxation rate
Figure 2 shows results for the average first-passage time across four orders of magnitude of the relaxation rate . The results strongly suggest that grows according to a power law with decreasing , if is sufficiently small. The exponent exhibits a simple dependence on the number of edges with identical waiting-time distributions. That number is equal to for processes according to Eqs. (22) and to for processes obeying Eqs. (23)– (25), since the transition is different there. The exponent depends to leading order on like . This applies to both the asymptotically Markovian and non-Markovian waiting-time distributions, and to both parameter sets of waiting-time distributions simulated in the non-Markovian case.
The exponent is equal to for , as has previously been shown analytically for according to Eq. (22) (see Dupont et al. (2016), Chapter 5). The process is very unlikely to reach large with a bias toward 0, even if this is only small. Hence, the random walk “waits” until the transient is over and symmetry of the transition rates has been reached [see Fig. 1(d)], and then goes to . This waiting contributes a time to , and approaches 1 for large .
The average first-passage time for a symmetric random walk increases with like Jouini and Dallery (2008), i.e., it is very long for large . We see a contribution of the transient to only if relaxation is slow enough for to be comparable to this long time. Consequently, is essentially independent of for large and large (see in Fig. 2(a, b)).
V Resonant length
The coefficient of variation CV [= standard deviation (SD)/average ] of the first-passage time reflects the relative fluctuations. Its dependence on the chain length is illustrated in Fig. 3. CV increases monotonically with increasing . Its dependence on is not monotonic. We find a pronounced minimum of CV() for small values of , where the minimal value is up to one order of magnitude smaller than CVs with and with large . Both simulations and analytical results show this behavior.
How can we get a heuristic understanding of the decrease in CV with decreasing and initially with increasing ? The lingering close to state 0 shown in Fig. 1(d) means that states with index larger than 1 are not reached before a time with almost certainty. This initial part of the process contributes to but little to SD. Hence, its growth with decreasing and (initially) increasing decreases CV. Additionally, the standard deviation is determined by the rates and splitting probabilities at the end of the transient, which are more favorable for reaching than the splitting probabilities at the beginning. It is therefore less affected than the average by the values. This also causes a decrease in CV with decreasing .
To gain more insight into the dependence of CV, we look at processes with constant transition rates and a symmetric process with transient rates but constant splitting probabilities [see Eq. (33)]. The permanently symmetric process exhibits only a very shallow minimum, and the minimal value of CV appears to be independent of [Fig. 4(a): discrete case]. Processes with constant transition rates exhibit decreasing CV with decreasing Jouini and Dallery (2008) [Fig. 4(b)]. However, this decrease is minor and does not explain the data in Fig. 3 in the range covered by the transient.
The comparison shows that the dynamics of the splitting probabilities provide the major part of the reduction in CV with increasing in Fig. 3. The splitting probability relaxes from 1 to its asymptotic value set by the parameters of the process. The larger the value of , the later is the absorbing state reached and the smaller is the value of when it is reached. As long as decreases sufficiently rapidly with increasing , so does CV. At values of such that , the decrease in is negligible, and CV starts to rise again with increasing toward its large- value.
These considerations suggests that the with minimal CV could be fixed by a specific value of the splitting probability. This value of cannot be smaller than , since we would expect a monotonically decreasing CV in that regime [the regime in Fig. 3(d)]. Since this should hold also for minima with , this value of also cannot be much larger than , since would then diverge. This is confirmed by the results shown in Fig. 3(f). The splitting probability at the minimum changes by less than 8% over four orders of magnitude of , and is slightly larger than . Hence, CV starts to rise again when the length is so large that the average first-passage time is long enough for to approach the symmetric limit
The value of CV() depends on in the case with non-Markovian waiting-time distributions. It is 1 for large values, since the asymptotic splitting probability is larger than [Fig. 3(c); see also Fig. 4(b)]. CV() is approximately equal to for small values of , since applies in that case [Fig. 3(c)]. CV() is exactly with asymptotically Markovian waiting-time distributions [Fig. 3(e)].
The parameter values in Fig. 3(d) entail . The process has a bias toward and we see the well-known behavior CV() for large . The onset of the behavior moves to smaller with increasing values of until finally the minimum of the CV is lost. We found minima if (four times the average state dwell time), but we have not determined a precise critical value.
The transition from the CV of to CV() for large values of is monotonic in Fig. 3(c) for the asymptotically non-Markovian case. The asymptotically Markovian system exhibits a minimum even for large values [: Fig. 3(e)]. However, it is comparably shallow, and is in line with the results of Jouini and Dallery for systems without transient Jouini and Dallery (2008) [Fig. 4(b)].
Vi Discussion and Conclusion
vi.1 (Generalized) master equation with exponential time dependence of rates and waiting-time distributions
In general, stochastic processes occur in changing environments or may be subjected to control. The corresponding mathematical description is provided by (generalized) master equations with time-dependent coefficients. The analytical solution of these equations is complicated even in the simplest case [see Eq. (A)] Skupin, A. and Falcke, M. (2010). The Laplace-transform approach is not applicable in general, but the exponential time dependences of the master equation (2) or the generalized master equation (8) with (9) do allow the transform to be performed on them, leading to difference equations in Laplace space. The solutions of these equations are given by Eqs. (17) and (20), resp. To the best of our knowledge, this is a novel (and efficient) method of solving the generalized master equation with this type of time dependence in the kernel, or the master equation with this type of time-dependent rates.
Equations (17) and (20) also apply in the case of state schemes more complicated than that in (21). It is only necessary to make appropriate changes to the definitions of and or of and . Hence, we can use the method presented here to investigate the dynamics of complicated networks in a transient. This is of interest for studies of ion-channel dynamics Sneyd et al. (2004); Gin et al. (2009), transport in more than one spatial dimension, gene regulatory networks, and many other applications. Many quantities of interest in these investigations are moments of first-passage type and thus can be calculated from the derivatives of the Laplace transforms. Some of these calculations might require knowledge of the probability time dependence , i.e., the residues of . In many cases, determination of these residues will be only slightly more complicated than for the system without transient, since all terms arising from the transient involve the same factors, but with shifted argument, together with the matrix . In general, if we know the set of residues of the Laplace transform of the system without transient and the matrix , we can immediately write down the residues of the Laplace transform of the system with transient by shifting by integer multiples of . This allows easy generalization of many results and might be of particular interest for renewal theory Cox (1970); van Kampen (2001).
vi.2 The average first-passage time
The average first-passage time decreases with increasing relaxation rate according to a power law, if the transient is sufficiently slower than the time scale of the individual steps (Fig. 2). The exponent depends on the chain length to leading order like . Remarkably, this behavior has been found for substantially different parameters in both the asymptotically non-Markovian and asymptotically Markovian cases, and therefore appears to be rather universal.
vi.3 Resonant length
Transients can substantially reduce the coefficient of variation. We have found a monotonic decrease of CV with decreasing when the relaxation rate is slow compared with state transition rates. CV exhibits a minimum in its dependence on the chain length (Fig. 3) if the transient is slow compared with state dynamics. This minimum exists for both asymptotically Markovian and non-Markovian discrete systems. At a time-scale separation between the state transition rates and the transient of about , the CV at the minimum is reduced by about one order of magnitude compared with CV().
Exploiting our results for control purposes, we can use a transient if more precise arrival timing at the outcome of a process is desired. Applying a transient that relaxes a bias toward 0 reduces CV. In the context of charge carrier transport, such a transient could be realized by a time-dependent electric field. If we include the rising part of CV() beyond the minimum in our consideration and consider all CVs smaller than CV() as reduced, we may still see a reduction of CV even if the average first-passage time is an order of magnitude longer than . Hence, the initial transient may have surprisingly long-term consequences and is therefore a rather robust method for CV reduction. Additionally, or if a transient is given, the state and step number can be used for optimizing the precision of arrival time. The optimum may also be at smaller lengths than that of the non-optimized process; i.e., optimization of precision may even lead to acceleration.
vi.4 How can negative feedback robustly reduce noise in timing?
Most studies investigating the effect of negative feedback have focused on amplitude noise—we are looking at noise in timing. Timing noise substantially reduces information transmission in communication systems. Variability in timing and protein copy numbers in gene expression or cell differentiation causes cell variability Yurkovsky and Nachman (2013); Hormoz et al. (2016). Therefore, many studies have investigated the role of noise in gene regulatory networks and have found that immediate negative feedback is not suitable for reducing timing noise Kontogeorgaki et al. (2016); Hinczewski and Thirumalai (2016); Ghusinga et al. (2017). This agrees with our results for large values.
We have shown that slow recovery from negative feedback is a robust means of noise reduction. The recovery process corresponds to the slow transient in our present study. Negative feedback terminating the excited state is a constitutive element of many systems generating sequences of pulses or spikes, such as oscillating chemical reactions Totz et al. (2015); Falcke and Engel (1994), the electrical membrane potential of neurons Kandel et al. (1991), the sinoatrial node controlling the heart beat Pohl et al. (2016), the tumor suppressor protein p53 Vousden and Lane (2007); Mönke et al. (2017), Ca spiking in eukaryotic cells Falcke (2004); Skupin et al. (2008); Thurley et al. (2014), cAMP concentration spikes in Dictyostelium discoideum cells and other cellular signaling systems Goldbeter (1996); Falcke and Levine (1998).
In particular, our results apply to noisy excitable spike generators in cells, where single molecular events or sequences of synaptic inputs form a discrete chain of states toward the threshold. Beside examples from membrane potential spike generation Grewe et al. (2017), Ca spiking and p53 pulse sequences have been shown to be noise-driven excitable systems Mönke et al. (2017); Thurley et al. (2014); Skupin et al. (2008); Thurley and Falcke (2011). The information content of frequency-encoding spiking increases strongly with decreasing CV of spike timing Skupin and Falcke (2009). A value of the coefficient of variation between 0.2 and 1.0 has been measured for the sequence of interspike intervals of intracellular Ca signaling Skupin et al. (2008); Thurley et al. (2014); Dragoni et al. (2011); Cao et al. (2014); Lembong et al. (2017). Hence, CV is decreased compared with that of a Poisson process and even compared with that for first passage of a symmetric random walk. The experimental data are also compatible with the finding that the slower the recovery from negative feedback, the lower is CV Skupin et al. (2008); Thurley et al. (2014); Dragoni et al. (2011); Cao et al. (2014); Lembong et al. (2017). This strongly suggests that this ubiquitous cellular signaling system uses the mechanism of noise reduction described here to increase information transmission.
Optimal noise amplitudes minimize the CV of interspike intervals in noisy excitable systems Jung and Mayer-Kress (1995); Pikovski and Kurths (1997); Gang et al. (1993); Beato et al. (2008); Lindner et al. (2004), which has been termed coherence resonance. Our results define the conditions of optimal CV reduction in terms of system properties—an optimal number of steps from the inhibited state to the excitation threshold during slow recovery. At the same time, our results shed light on new aspects of coherence resonance, and indeed may indicate that a more fundamental phenomenon underlies it. Since we believe excitable systems to be one of the most important applications of our results, we have chosen the term resonant length.
Coming back to the widely used graphic example, what can the drunkard learn from our results? If he chooses a pub at the right distance from home, he will arrive home sober and relatively in time for breakfast.
Acknowledgements.VNF acknowledges support by IRTG 1740/TRP 2011/50151-0, funded by the DFG/FAPESP.
Appendix A The symmetric asymptotically Markovian case and its continuum limit
The permanently symmetric case with a transient is defined by
Its continuum limit on a domain of length (with spatial coordinate ) satisfies
The Fokker–Planck equation (34) can be solved after applying a time transformation :
The boundary and initial conditions, , , and specify the first-passage problem. We find
where . The th moment of the first-passage time is given by