Stochastic kinetics on networks: when slow is fast.

Most chemical and biological processes can be viewed as reaction networks in which different pathways often compete kinetically for transformation of substrates into products. An enzymatic process is an example of such phenomena when biological catalysts create new routes for chemical reactions to proceed. It is typically assumed that the general process of product formation is governed by the pathway with the fastest kinetics at all time scales. In contrast to the expectation, here we show theoretically that at time scales sufficiently short, reactions are predominantly determined by the shortest pathway (in the number of intermediate states), regardless of the average turnover time associated with each pathway. This universal phenomenon is demonstrated by an explicit calculation for a system with two competing reversible (or irreversible) pathways. The time scales that characterize this regime and its relevance for single-molecule experimental studies are also discussed.

Introduction

Many chemical and biological systems are composed of numerous species interacting in complex networks, such as cell signaling networks, enzymatic reaction networks, and genetic regulatory networks.[1−4] In recent years, the monitoring, analysis, and detection of single molecule transformations became one of the most vigorously growing research areas in physics, chemistry, and biology.[5,6] The ability to detect one molecule at a time leads to significant advances in uncovering fundamental properties of complex chemical and biological processes. Through the analysis of single-molecule activity we are able to understand better the sequence and timing of various processes, such as chemical transitions for molecular motor proteins[7−10] and for ribosomes[11] or the sequence of enzymatic turnovers in single-molecule reactions.[12,13]

At the single-molecule level, stochastic effects are recognized to play an important role. One approach to study these effects is to model processes on complex networks as Markov chains. The fundamental coupling between the structure and dynamic properties of complex systems has been the focus of many studies, both theoretically[14,15] and experimentally.[16,17] However, we still have little knowledge on any universal relations between dynamic properties and the structures of networks. Recently, we have been interested in understanding the time evolution of chemical and biological processes by exploiting the properties of first-passage times,[18,19] both in simple systems[20] as well as in complex networks,[21−23] with applications to molecular motors and macromolecular turnover phenomena.[24] It was found that at early times the chemical and biological processes proceed mostly along the shortest pathways (in the number of intermediate states) that connect initial and final states. However, time scales for this universal behavior and its consequences for mechanisms of chemical reactions have not been discussed. These studies also demonstrate that the first-passage approach is a powerful tool to understand microscopic mechanisms of complex processes.

In complex chemical and biological systems, often many pathways act simultaneously and the contribution of each pathway has to be considered in order to comprehend the behavior of the system. Generally speaking, the existence of different enzymes or reactants with varied activities will lead to distinctive properties for each pathway. It is usually expected that the fastest pathway, i.e., the one with shortest average times for the transition between initial and final states, will dominate the process of creation of the final product. However, the competition between various pathways is much more subtle than predicted from these simple arguments and it is instead the result of a complex interplay between the length of pathways and times to proceed along them.

There are many chemical and biological processes that occur extremely fast with time scales of microseconds or even femtoseconds, such as protein folding,[25] isomerization of rhodopsin,[26] energy transfer in photosynthesis,[27] and so on. Whether these ultrafast processes are still fully determined by their kinetics is not fully understood, and not much is known about mechanisms of these processes. Recently, it was observed[28] that almost the same time (a few microseconds) is consumed for fast- and slow-folding proteins as the folding events really occur, although their folding rate coefficients might differ by several orders of magnitude. Interestingly, our investigations for complex networks with competing pathways indicate that the underlying structure but not the kinetics for those pathways would determine the properties of systems at very short time scales. A commonly applied relation is obtained, suggesting that it is always the shortest pathway that makes the major contribution to realizations whose duration is shorter than a certain threshold value, even if it is the slowest one on average. The time scales that determine the threshold value are also discussed, and it is shown that these time scales are in general dependent on the underlying reaction network.

Theoretical Method

To better understand this counterintuitive phenomenon, we consider a very basic stochastic scheme that includes the main fundamental aspects of single-molecule reactions. There are two competing pathways available for the chemical reaction to proceed from the initial state i to the final state j, as shown in Figure 1. For example, it might correspond to enzymatic processes where one pathway describes the chemical reaction without the catalyst and the second pathway is due to the presence of the enzyme. Here, we will focus on investigating the turnover times that describe how a single molecule is transformed from its initial state i into its final state j in the presence of two such pathways.

Figure 1

A generic network with competing chains. The transitions from the initial state i to the final state j can be realized through two pathways with intermediate states a (k = 1, ..., m) and ( = 1, ..., n) labeled by blue and red circles, respectively. It is assumed that the upper pathway A has more states than the lower pathway B; i.e., m > n. The forward transition rates are associated with λ along pathway A, while all other transition rates are given by μ.

The states i and j are connected by two competing routes, called chains A and B, respectively; see Figure 1. Along chain A there are m states, which are labeled as a (k = 1, ..., m). Similarly, there are n states along chain B, labeled as , with = 1, ..., n. To keep the notation simple and to limit the number of parameters, we associate a rate λ to all forward transitions from the state i to the state j along chain A and another rate μ to all other transitions, i.e., the backward transitions on the paths from j to i along both chains A and B and the forward transitions along chain B (different choices are discussed later in this paper).

By taking λ to be sufficiently larger than μ, the dwell times on the states in chain A can be made so short that transitions along the chain A become increasingly faster as λ grows. Moreover, the fraction of realizations that reach the final state j from the state a can easily exceed the number of realizations that reach j from the state b, even when m > n. To formalize this point, we introduce a random variable R that takes values in {A, B}: when R = A we will say that a realization has reached the state j from chain A, i.e., through the state a, whereas when R = B it will mean that a realization has reached the state j from chain B, i.e., through the state b.

Results and Discussion

Probabilities To Reach the Final State through the Competing Pathways

The probabilities U and U that R is equal to A and B, respectively, can be computed by means of standard techniques.[22,29] The analytical expression for the probability U that the system reaches the final state j through pathway A, as shown in Figure 1, is given byfor any non-negative integer value of m and n, where the parameter x ≡ μ/λ and the normalization relation U + U = 1 holds. One can readily see that as λ grows (the parameter x becomes smaller) the probability U to reach to the final state j through the chain A becomes larger and approaches unity. It can be proved that this probability function is a decreasing function of x, as demonstrated in Figure 2. The derivative of U with the respect of x can be written asOne can easily show that dU/dx < 0 because we havefor 0 ≤ x ≤ 1.

Figure 2

Probability U for the system to reach the final state j through longer pathway A as a function of the variable x = μ/λ for the network showed in Figure 1. The number of states for pathway A and B are taken as m = 5, n = 4, respectively. For μ = 1 and λ = 2, the probability U is definitively larger than 0.5. This justifies the choice of this parameter set in the main paper.

The probability U will be close to 1 as the parameter x approaches 0, as indicated in Figure 2, which means that the final state j is reached mainly through the longer path A if the corresponding kinetic transitions in this path becomes faster. However, the picture changes radically when we consider the contributions for the realization of the process from the two pathways for varied periods of time.

First-Passage Time Densities along the Competing Pathways

The first-passage time approach[18,19] is a powerful tool to investigate the temporal evolution of many chemical and biological systems. Here, we define variables T and T as the conditional random times to reach the final state j from the initial state i when the last visited state prior to j is the state a and b, respectively. The normalized probability densities of T and T will be denoted with ϕ and ϕ, respectively. These densities can be computed using the technique of the absorption times explained in previous studies.[22,30] First, we discuss the average conditional absorption time along the competing pathways, which are defined asfor X = A, B. ⟨T⟩ and ⟨T⟩ describe the average times that it takes for the system to reach the final state from the competing pathways. The less time it consumes, the faster it is for the system to reach the final state through the pathway. Similar behaviors as obtained for the probabilities U and U can be observed for these average properties of the system. As shown in Figure 3, the average value of T becomes smaller than the average of T when the transition rate λ in pathway A becomes larger. Therefore, when λ is large enough, the longer pathway A is faster and more productive than the shorter pathway B, on average. However, distinctive phenomena will be observed when we abandon the view of the process based solely on average values.

Figure 3

Average conditional first-passage times ⟨T⟩ and ⟨T⟩ for the system to reach the final state j starting from the state i over pathway A or B as a function of the rates λ for μ = 1. The number of states m and n are same as used in Figure 2. For λ = 2 the average time over pathway A is clearly smaller than the average time over pathway B.

We also define the random variable T as the time that the system spends before visiting the state j for the first time, starting from the state i. Given that a realization leaving i and reaching j took exactly a certain time t, i.e., given that T = t, we would like to determine the probability that the reaction route comes along pathway B with slower kinetics. Starting with the condition T ≤ t and using Bayes theorem, we have thatNotice now that by the definition of the times T given earlier, on the right-hand side of eq 5 we also havewhereas for the denominator we havewhere ϕ is the unconditional first-passage time density ϕ(τ) = Uϕ(τ) + Uϕ(τ). Using the same procedure by conditioning on t < T ≤ t + δt and taking δt → 0, it finally leads toNotice that this equation is very general, as it requires only that in a reaction network there exists a pathway B and at least another competing pathway. No further assumptions have been made on the rates or on the structure of the network in order to derive it. Nevertheless, when we focus on the network shown in Figure 1 and on the choice of the parameters as described above, we can easily realize that eq 8 can be nonmonotonic as a function of t. Indeed, on the basis of the results of ref (22), for the network in Figure 1, ϕ(t) goes to zero faster than ϕ(t) as t → 0 because pathway A is a longer chain. Thus, eq 8 equals unity at t = 0 and is thus decreasing for increasing t in the neighborhood of t = 0. Moreover, as t → ∞ this function can rise again if ϕ(t) decays more rapidly than ϕ(t) in this limit (an example where this happens is discussed in more details later). At the moment we focus our attention at the regime of times close to zero, where eq 8 is monotonically decreasing as a function of t. We call this time regime the regime of ultrafast realizations of the process depicted in Figure 1. In the following, we are going to characterize this regime more precisely.

Since ultrafast processes are critically important and widespread in chemical and biological system, we are interested in ultrafast realizations that are those occurring at small times t for the system shown in Figure 1. We use eq 8 to define a time ts, which fulfillsso that the probability that the chemical reaction has occurred through pathway B is larger than θ if the realization T had duration smaller than t. The relationship between θ and t can be obtained analytically from eq 8 once ϕ and ϕ have been computed from the Master equation associated with Figure 1. We should therefore notice that the value of t associated with a certain θ and thus also the range of time scales associated with ultrafast realizations is dependent on the overall structure of the network and on the choice of the rate constants.

The black line in Figure 4 shows the relationship for the network in Figure 1 with the number of states m = 5, n = 4 (the same as was used in Figures 2 and 3) and the transition rates λ = 2, μ = 1 (arbitrary units). As indicated by the dotted lines in Figures 2 and 3, the final state j is mainly reached from the longer pathway A with faster rates than the shorter pathway B, on average. However, the probability θ to reach the final state j through the shorter pathway B increases as the realization of the process becomes ultrafast, as shown in Figure 4, and the value of θ will be larger than 1/2 as time t becomes smaller, which means that the shorter pathway B will dominate the realization of the process below that time. At the limiting case as time t approaches 0, all of the chemical reactions will occur only through pathway B, which gives slower kinetics but contains fewer intermediate states. This is a universal phenomenon that is expected to be observed regardless of the values of the kinetic parameters for all networks. However, this is a surprising result that cannot be understood from the discussion of average properties of the system.

Figure 4

Analytical and approximate relationships between θ and ts. The black line is given by the direct analytical solution of eq 9 upon computing ϕ and ϕ via solving the corresponding Master equation. The red line is given by the approximate solution derived from eq 12. In this figure the network from Figure 1 with m = 5, n = 4, λ = 2, and μ = 1 (rates in arbitrary units) is considered. The existence of ts is a universal property, only its value depends on the rates.

We should stress again that both eqs 8 and 9 are very general results that hold independently of the example studied in Figure 1. Nevertheless, we have chosen to study the case given in Figure 1 (and a simple variant of it discussed later in the paper) because we believe that this example is the most instructive one, since it highlights the contradiction between being fast on average and being responsible for the ultrafast realizations. In the context of Figure 1 and for values of θ particularly large, it is possible that the events considered in eq 9 are particularly rare. Later on we show, however, that even if rare these events are well above the detection capacity of modern experimental techniques for typical catalytic reactions.

Approximate Expression for the Relationship between θ and ts

For more realistic application it may be useful to obtain an approximate expression of ts when θ is large (approaches 1), and then the relation between these two variables can be observed directly. Since we aim to analyze the behavior of the system at early times, we just need to calculate the densities in eq 8 for small values of t. In particular, we can use the results developed recently in refs (21) and (22) and expand eq 8 at small times. Indeed, using now the graph theoretical approach,[22] we find thatfor t approaching 0, where o(t) are all other terms that satisfy o(t)/t → 0 as t → 0 and I(t) is a polynomial containing terms from power t to power t. The same approach tells us that for m > n the unconditional first-passage time density ϕ(t) can be expanded aswhere I(t) is the same polynomial as before and o(t) are terms of order larger than m when t approaches 0. Substituting eqs 10 and 11 into eq 8 and then combining with eq 9 finally leads to a first-order approximation inwhere x ≡ μ/λ and is an approximate expression for ts defined in eq 9 at small times. Notice that this result will depend only on the forward rates along the two pathways A and B because all backward rates can only appear at the terms I and o(t), as discussed in refs (21) and (22). From the expression above, we can obtain directly that is always a decreasing function of θ given m > n, as observed in Figure 4 regardless of the transition rates λ and μ. Therefore, it is a general relation as discussed above that the shorter pathway B will dominate the process for the system to reach the final state j at small times and the corresponding probability θ through this path even reaches to 1 as time becomes close to zero. From this relation, we can also observe that will become smaller as the forward rate λ in the longer pathway A increases or the forward rate μ in the shorter pathway B decreases. Therefore, we need to detect the process at ultrafast or smaller time regime if we want to observe the system to reach the final state j mainly through the shorter pathway B as transition kinetics becomes faster for longer pathway A or slower for pathway B itself. In Figure 4 the red line gives the plot of as a function of θ and it is compared with the exact solution ts derived earlier, which is indicated by the black line. It is clear that at large values of θ (close to 1), provides a perfect estimate of ts while is systematically smaller than ts as θ becomes smaller.

A general expression for can also be obtained for the two pathways with arbitrary transition rates by using the expansion technique from ref (22). Given the forward rates λ, with k = 1, ..., m in pathway A, and , with = 1, ..., n in pathway B, it leads towhere μ and λ are the forward rates associated with the state i in Figure 1. This general expression is quite similar to eq 12, and the conclusions obtained previously are not influenced. It also shows that the knowledge of only the forward rates is therefore sufficient to estimate the time scale ts. This should allow us an easy way to verify theoretical predictions by means of single-molecule enzymatic reactions when the reaction rates at any state can be controlled or the number of intermediate states for those pathways could be regulated.

Observation of Ultrafast Phenomena in Real Systems

Our theoretical method clearly proves that generally chemical reactions might proceed along the shortest pathway, even if they are not the fastest. However, one might ask the question if these ultrafast processes can be observed in real chemical and biological systems. Although current single-molecule experimental methods are quite advanced, their temporal and spatial resolutions are not infinitely perfect. It is important to estimate the probabilities and time scales when these ultrafast phenomena might take place using realistic conditions.

To perform such calculations, we will employ eq 12 and assume that t̂s gives us the time scale for ultrafast realizations of the chemical process. In addition, from eq 11 it can be found that the probability of observing such ultrafast reactions, Puf, can be estimated at small times as

Now let us consider a chemical reaction that proceeds in one transition without intermediate states (n = 0) and the average time for this process is on the order of 1 s (μ = 1 s–1). We added enzyme molecules to accelerate this process, and it is assumed that the catalyzed reaction is taking place via a new pathway with one intermediate state (m = 1). The catalytic rates are typically on the order of λ = 103 s–1,[31] i.e., the process is accelerated 1000 times. Then from eq 12 the ultrafast realizations are taking place for times faster than ts ≈ 1 μs, while the probability of such events from eq 14 is Puf ≈ 10–6. If the number of intermediate states in the enzymatic pathway is larger, which corresponds to a more realistic situation, say m = 3, then the ultrafast reactions can be observed for t ≤ 100 μs with a probability of Puf ≈ 10–4. It is clear that, although these times and probabilities are small, the precision of current experimental techniques is high enough so that they can be observed. One can see also that the specific range of parameters for observing ultrafast processes depends on relative values of chemical transition rates in the shortest pathway and in the catalyzed pathway as well as on the difference in the number of intermediate states in each path. However, these calculations support our arguments that most probably these ultrafast phenomena can be experimentally accessed and tested with current experimental methods.

An Example from a Simple Irreversible Network

It will be useful also to illustrate this phenomenon by analyzing a simpler model that shows the same qualitative behavior but can be completely solved analytically. We consider the network in Figure 1, where all backward rates are set to 0 (see Figure 5). As shown in eq 12 the reversibility of the transitions will not change the properties discussed below. Starting from the state i, the probability of reaching the state j along path A is now simply given bysince the process has no possibility to return to the state i after leaving it. Therefore, the probability to reach the state j along path B is described by U = 1 – U. The conditional probability densities ϕ(t) and ϕ(t) as defined above can be easily obtained by solving the corresponding backward Master equations through Laplace transform, as shown earlier.[21] The (normalized) probability density ϕ(t) will take the explicit formfor processes conditioned to reach the state j along path A. Similarly, the probability density ϕ(t) along path B can be obtained asUsing these probability density distributions, the conditional mean first-passage time and other dynamic properties of the system can also be calculated accordingly.

Figure 5

A generic network with two competing chains similar to the network shown in Figure 1 but with irreversible transitions from the initial state i to the final state j in both pathways.

From now on, we assume that the rate λ is sufficiently larger than the rate μ so that m/λ < n/μ. By using eq 8 and the definition of unconditional first-passage time density ϕ(t) = Uϕ(t) + Uϕ(t), we can easily obtain the exact condition that path B has a probability θ for realizations of the process at duration t. We can also compare the probability Uϕ(t) for process occurred through path A and probability Uϕ(t) for path B at time t directly. An example is shown in Figure 6 for the network illustrated in Figure 5. The blue and red curves give the probability densities from paths A and B, respectively, with parameter values the same as used in Figure 4, except for the backward rates taken to be equal to 0. It clearly shows that the shorter path B has a higher probability than the longer path A at small times, even though the former one gives slower kinetics just as predicted from our previous discussions. At the same time, the total probability U to reach the final state from the longer path A is higher than that from path B because of the faster kinetics when λ > μ. It is consistent with the observation that the longer (faster) path A will have a higher probability at intermediate times, where most of the realizations of events are obtained. Interestingly, we found that under these conditions the shortest path will have a higher probability again at even larger times (see Figure 6). The third regime is unexpected, but it is easy to explain if one looks at the explicit expressions 16 and 17. Indeed, the dominant term at large t in both expression is proportional to te– for pathway B and te– for pathway A, with λ > μ. Clearly, pathway B will have the larger tail at large t. Thus, after a sufficiently large time, the probability to reach the final state through path A will be smaller because most of the processes through this path have already occurred. At very large times, the main contribution to reach the final state is again mostly from the path with slower transition rates. The contributions from two pathways may change at intermediate and large times according to the variation of kinetics in each pathway; however, it applies universally that the shortest pathway always gives the highest probability to reach the final state at early times, regardless of the kinetics of the system.

Figure 6

First-passage time probability densities from the state i to the state j through two paths A and B for the irreversible network shown in Figure 5. The blue and red lines correspond to the functions Uϕ(t), computed from eq 16, along the longer path A, and Uϕ(t), computed from eq 17, along the shorter path B, respectively. The number of states for two paths are given by m = 5 and n = 4. The transition rates λ = 2, μ = 1 are used.

Summary and Conclusions

We have investigated stochastic dynamics in complex chemical and biological networks by analyzing probabilities and times for single-molecule reactions. A complex system usually contains many competing pathways as substrate is converted to product. We have indeed considered two reversible (or irreversible) pathways of unequal lengths where transitions between the initial and final states occur through the longer path faster on average than via the shorter route. As expected, the total probability to reach the final state through the longer path will increase and become higher than the shorter path if faster transition rates are associated with this path. Besides, it takes less time for the realization of the process through the longer path with faster transitions on average. From these observations, it seems that the pathway with the fastest kinetics would always dominate the transformation process of product from substrate, irrespective of other properties of the system. However, it is found that, in contrast to expectations, the major contributors for the realization of final product formation are not always the pathways with fast kinetics. Indeed, we have shown in this work that there is always a time scale at early times when the shortest pathway will be the main realization of the chemical process. In addition, using realistic parameters we estimated time scales and probabilities of observing these phenomena. Our calculations indicate that ultrafast phenomena might be observed using modern experimental techniques.

With the vigorous growth of single molecule techniques, it is now possible to observe ultrafast reactions in many chemical, physical, and biological systems, and new phenomena and mechanisms are anticipated to be discovered. Through theoretical analysis of the network systems, we observe that the probability to reach the final state is a decreasing function of time for the shortest pathway at early times. And the majority of the ultrafast realizations, i.e., realizations that reach the final state j in a time less than a certain value ts, is surprisingly found to come from the shortest pathway. This property of ultrafast reactions is only determined by the structure of the network but not the kinetics of the system. One can provide the following simple explanation of this phenomenon. The early time behavior of the first-passage density along a path that connects the initial state i and the final state j is described by a power law distribution tα with the exponent α given by the number of states along the pathway.[21,22] Therefore, it is easy to observe that a shorter path with fewer states gives larger values for the probability density to reach the final state at early times.

Our conclusions do not depend on the number of routes nor the topology of networks, and similar expressions as 12 and 13 can be obtained simply by employing the expansion technique developed for the general network systems.[22] Therefore, it suggests that this counterintuitive observation is a universal phenomenon and might indicate a fundamental mechanism governing many chemical and biological systems. Recent studies show that many ultrafast reactions cannot be understood simply from the kinetic properties of the system. Our findings might assist in understanding these complex phenomena. Our study can also provide a simple way to distinguish the contributions from the various reaction pathways. It should lead to a better understanding of mechanisms in complex stochastic systems. We believe that it will be possible to verify our predictions experimentally by following the single-molecule turnover in in vitro experiments where two different reaction routes can be distinguished. The two competing pathways should contain a long route that is very fast with smaller mean turnover times and a shorter route that is slow on average. This could be realized by considering a single-molecule reaction where the transition rates could be increased by the catalysis of a specific enzyme. Then the shorter path for the realization of single-molecule turnover without enzymes would have a slower kinetics compared with the longer one with enzymatic reactions.