Dynamics of Molecular Motors in Reversible Burnt-Bridge Models

Dynamic properties of molecular motors whose motion is powered by interactions with specific lattice bonds are studied theoretically with the help of discrete-state stochastic"burnt-bridge"models. Molecular motors are depicted as random walkers that can destroy or rebuild periodically distributed weak connections ("bridges") when crossing them, with probabilities $p_1$ and $p_2$ correspondingly. Dynamic properties, such as velocities and dispersions, are obtained in exact and explicit form for arbitrary values of parameters $p_1$ and $p_2$. For the unbiased random walker, reversible burning of the bridges results in a biased directed motion with a dynamic transition observed at very small concentrations of bridges. In the case of backward biased molecular motor its backward velocity is reduced and a reversal of the direction of motion is observed for some range of parameters. It is also found that the dispersion demonstrates a complex, non-monotonic behavior with large fluctuations for some set of parameters. Complex dynamics of the system is discussed by analyzing the behavior of the molecular motors near burned bridges.


INTRODUCTION
In recent years an increased attention has been devoted to investigations of molecular motors, also known as motor proteins, that are crucial in many cellular processes [1]. They transform chemical energy into the mechanical motion in non-equilibrium conditions. For most of molecular motors their motion along linear molecular tracks is fueled by the hydrolysis of adenosine triphosphate (ATP) or related compounds. It was suggested that a different mechanism is employed to power the motion of a protein collagenase along collagen fibrils [2,3]. It probably utilizes the collagen proteolysis, cleaving the filament at specific sites.
As the collagenase molecule is unable to cross the already broken bond, it leads to the biased diffusion along the filament. However, full understanding of mechanisms of collagenase motion is still not available.
It was proposed that a good description of the collagenase dynamics could be provided by the so-called "burnt-bridge model" (BBM) [2,3,4,5,6,7,8,9]. In this model, the motor protein is depicted as a random walker that translocates along the one-dimensional lattice that consists of strong and weak bonds. While the strong bonds remain unaffected if crossed by the walker in any direction, the weak ones (termed "bridges") might be broken (or "burnt") with a probability 0 < p 1 ≤ 1 when crossed in the specific direction, and the walker cannot cross the burnt bridges again, unless they are restored, which can occur with probability 0 < p 2 ≤ 1. In Refs. [6,7] an analytical approach was developed which permitted us to derive the explicit formulas for molecular motor velocity V (c, p 1 ) and diffusion constant D(c, p 1 ) for the entire ranges of burning probability 0 < p 1 ≤ 1 and concentration of the bridges 0 < c ≤ 1 which were also confirmed by extensive Monte Carlo computer simulations. This theoretical method has been applied to several problems with periodic bridge distribution. However, the results in [6,7] have been obtained only for irreversible bridge burning (bridge recovery probability was taken to be p 2 = 0), and also for unbiased random walker between bridges (equal forward and backward transition rates). In present work, we generalize our approach to allow for the possibility of bridge recovery as well as unequal hopping rates on the sites between bridges. It is more realistic to consider systems with reversible action of motor proteins since they are catalysts that equally accelerate both forward and backward biochemical transitions [1].
MODEL According to our model, we view a motor protein as a random walker moving along an infinite one-dimensional lattice with forward and backward transition rates being u and w correspondingly, as illustrated in Fig. 1. The lattice spacing size is set to be equal to one. The lattice is composed of strong and weak bonds. There is no interaction between the random walker and strong bonds, however crossing the bridge in the forward direction (from left to right) leads to its burning with the probability p 1 , while the particle moves with the rate u. After the weak link is destroyed, the walker is assumed to be on the right side of it. When the particle is trying to cross a broken bond, the bridge can be recovered with the probability p 2 , while the particle moves to the left with rate w. It is assumed that initially, at t = 0, all bridges are intact.
The details of breaking weak bonds in BBM have a strong effect on the dynamic properties of motor proteins [6]. There are two different possibilities of bridge burning. In the first variant (the so called "forward BBM"), the weak bond is broken when crossed from left to right, but the intact bridge is not affected when the particle moves from right to left. Thus the bridge recovery may occur if the walker attempts to cross a burnt bridge from right to left. In the second variant (named "forward-backward BBM"), the weak link is destroyed if crossed in either direction [4,5]. Both variants are identical for p 1 = 1, however for p 1 < 1 the dynamics is different in two burning scenarios, as was shown in p 2 = 0 case [6], although mechanisms are still the same. For reasons of simplicity, below we will only consider forward BBM, even though forward-backward BBM can also be solved using the same method.
There are five parameters that specify the dynamics of molecular motors in BBM: the probabilities p 1 , p 2 , the concentration of bridges c, as well as transition rates u and w.
The dynamic properties of the walker are also strongly influenced by the distribution of weak bonds [5]. Below we will study the case of periodically distributed bridges, when their concentration is c = 1/N and the weak bonds are located between the lattice sites with the coordinates kN − 1 and kN, with integer k (see Fig. 1). This description is more realistic for collagenases' dynamics [2,3]. The model below will be studied using continuous time analysis as it better describes chemical transitions in motor proteins [6]. Parameters are described in detail in the text.

Velocity
To find the walker's velocity, we generalize the method used in [6] (for irreversible bridge burning, i.e. p 2 = 0) to allow for non-zero probability p 2 . We introduce a probability R j (t) that the random walker is found j sites apart from the last burnt bridge at time t. The probabilities R j (t) arise if the system is viewed in moving coordinate frame with the last burnt bridge always at the origin, as illustrated by a reduced chemical kinetic scheme shown in Fig. 2.
The dynamics of the system is determined by a set of master equations: for k = 0, 1, 2, · · · and i = 1, 2, · · · , N − 2; and with k = 0, 1, 2, · · · for both Eqs. (2) and (3). Also at the origin we have In Eq. (2) we introduced a function f (k) as a probability that next to the last burnt bridge is k periods to the left from the last burnt bridge. It satisfies the condition ∞ k=1 f (k) = 1, which is reflected in Eq. (4). The system of equations (1)-(4) is to be solved in the stationary-state limit (at large times) when dR j (t)/dt = 0 is satisfied, and we denote R j (t → ∞) ≡ R j in what follows. By definition, it can be argued that Solving the system (1)-(4) can be facilitated by rewriting Eq. (5) in a more convenient form.
To this end, we note that based on the results from Refs. [6,8] it is reasonable to assume that R kN +i is of the form where y and W are some functions of p 1 , p 2 , u, w and N. Furthermore, y is i and kindependent, while W depends on i (but not on k). The advantage of the ansatz (6) is that it leads to a simpler form of f (k): as follows from Eq. (5). We proceed to solve Eqs. (1)-(4) with f (k) in Eq. (2) given by the expression (7).
We note that although the Eq. (4) was not used to find this solution, it was numerically verified that every equation in the system (1)-(4) is indeed solved by expressions (25) and (23).
Parameter R 0 needed to find the velocity is found from Eq. (25) combined with the normalization condition The mean velocity of the walker is given by [ which results in a simple relation, In Eq. (28), R 0 is given by (26) with y from the expression (23).
It can be shown that in the limit of u → 1, w → 1, and p 2 → 0 Eq. (28) reproduces the result obtained earlier in the Ref. [6] for the BBM with u = w = 1 and p 2 = 0. Also, Eq. (28) simplifies considerably in the limiting case of p 1 = 1 (deterministic bridge burning) when y = 0. In the case of p 1 = 1 we obtained V (u, w, p 2 , N) in [9] using the Derrida's method [10] and our general result given in Eq. (28) agrees with it in the p 1 → 1 limit, as was numerically verified.

Diffusion Coefficient
The diffusion coefficient is found by generalizing the method developed in [7] (where we found dynamic properties of the random walker in BBM with u = w = 1, p 2 = 0 and periodic bridge distribution), allowing for 0 < p 2 ≤ 1 and u = w. We define P kN +i,m (t) as the probability that at time t the random walker is located at point x = kN + i (i = 0, 1, · · · , N − 1), the right end of the last burnt bridge being at the point mN. Parameters m and k ≥ 0 assume integer values.
The dynamics of the system is described by a set of Master equations: for k ≥ 0 and i = N − 1 [with f (k + 1) the same as in (2)], for k ≥ 1 and i = 0, and for k ≥ 0 and i = 1, · · · , N − 2.
We observe that where R kN +i (t) is the probability for the random walker to be found kN + i sites apart from the last burnt bridge at time t, which was used above to find the walker's velocity and is given by Eq. (25). Plugging (33) into Eqs. (29) -(32) results in the equations (1) -(4) for R kN +i , thus obtained by a different method.
In accordance with [7,10], we introduce auxiliary functions S kN +i (t), The system of equations describing the time evolution of the functions S kN +i (t) results from Eqs. (34) and Eqs. (29) -(32). It was obtained that for k ≥ 1 and At t → ∞ the solutions of Eqs. (35) -(38) are sought in the form where a j and T j are time-independent coefficients. Plugging Eq. (39) into Eqs. (35) -(38) leads to, for k ≥ 1 and for k ≥ 0 and i = 1, · · · , N − 2.
Clearly, Eqs. (40) -(43) are identical to the system of equations (1) -(4) for the functions R j in the t → ∞ limit, where dR j /dt = 0. Thus their solutions should coincide up to the multiplicative constant, namely, with R kN +i given by Eq. (25). The normalization condition a kN +i . To find the explicit expression for C, we utilize the equations for T j obtained by plugging Eq. (39) into Eqs. (35) -(38), for k ≥ 1 and for k ≥ 0 and i = 1, · · · , N − 2. Summing up Eqs. (45) -(48) and using where the walker's velocity V is given by (28). Hence, in accordance with Eq. (44) where R kN +i and V are given by Eqs. (25) and (28).
Now we are able to obtain the expression for the random walker's velocity [7]. The mean position of the particle is given by which results in the mean velocitỹ In the t → ∞ limit Eqs. (39) and (49) therefore yield As expected, Eq. (53) reproduces the expression (28) for V obtained above with the use of the reduced chemical kinetic scheme method.
In analogy with finding R kN +i from (1) -(4), we start with solving (48). In the expression so that where R kN = R 0 y k , with R 0 , y given by (26), (23). With the use of (55), (48) takes the form We seek the solution of (56) in the form In (57), C 1 (k)+C 2 (k)β i part is the solution of homogeneous equation [the part of (56) which involves only T kN +i , T kN +i±1 ], in analogy with Eq. (8). Substituting Eq. (57) into Eq. (56) gives these expressions for A(k) and B(k): Plugging T kN , T kN +1 instead of T kN +i into (57), one can express C 1 (k) and C 2 (k) in terms of A(k), B(k) as well as T kN and T kN +1 (first two points of the period). Substituting resulting expressions for C 1 (k) and C 2 (k) together with the expressions (58) for A(k), B(k) in (57) gives after some rearrangement where  (59) and (55) with i = N, i.e., T kN +N = T (k+1)N , which gives as follows from (48) which was formally extended to include i = N − 1. In (62), the function f (k + 1) is expressed according to Eq. (7). Plugging (59), (55) (with appropriate k, i) into (62) permits us to express θ k in terms of T kN in this way, where It should be mentioned that substituting (59) and (55) where coefficients a, b, c, d, e, f are given by Eqs. (15) -(20) and In Eq. (66) it was found that and F is given by (64).
Comparison between (65) and (14) shows that they are identical (with replacement R → T ), with the exception of χy k term in (65). This implies that with arbitrary constantB solves homogeneous equation [Eq. (65) without χy k term], since y given by (23) is constructed for R kN = R 0 y k to solve specifically Eq. (14).
To compute the diffusion coefficient D, we need to consider additional auxiliary functions [7,10], It follows that for k ≥ 0, for k ≥ 1, and for k ≥ 0 and i = 1, · · · , N − 2.
The diffusion constant is to be found from Using (76) and summing up Eqs. (77) -(80) results in thus we do not need to find U kN +i (t) from the system (77) -(80) to obtain the diffusion coefficient. To derive (82) we used the normalization conditions, Given that in the stationary-state limit d dt and utilizing x(t) given by (51), it follows that Plugging (83) and (84) into (81) gives the diffusion constant, The time-dependent part of (85) is terms, and it is given by In order to get the final expression for diffusion coefficient, it is necessary to calculate T kN +i in (87). Using Eq. (75) for T kN +i yields, In deriving (88), we used the fact that 0 ≤ y < 1, 0 ≤ z < 1. Next, we consider the last two terms in (87), where we used (88) and T 0 =B +C [Eq. (73)]. The contribution from terms ∝B in (90) can be shown to be Utilizing Eq. (26) for R 0 , it follows from (91) thatB-contribution in (90) equals 0, thus undetermined constantB cancels out in (87) and it has no effect on the diffusion coefficient.
WithoutB-terms Eq. (90) becomes Plugging (92) into (87) gives In Eq. (93) we have β = u/w, and parameters R 0 ,C, y, z, F , λ are given by Eqs. It was verified numerically for various p 1 and c values that in the limit of u → 1, w → 1, reproduces the diffusion constant obtained in [7] for BBM with u = w = 1 and p 2 = 0. In the limiting case of p 1 = 1 we obtained D(u, w, p 2 , N) in [9] using Derrida method [10] and it also agrees with our general result (93) in the p 1 → 1 limit, as was numerically checked.

DISCUSSIONS
To illustrate our findings, we plot the dynamic properties of the molecular motor using where discrepancy exceeded 0.008 and 0.006 correspondingly. For Fig. 4(b), we compared p 1 = 1 case (not shown) with the corresponding case for u = w = 1 obtained in [9]: typical discrepancy in D values between u = 0.999 and u = 1 cases was ∼ 0.0005 for all c values except c < ∼ 0.001, where discrepancy exceeded 0.007. As anticipated, when the recovery probability p 2 → 1 and the presence of bridges has no effect, V → u − w = 0 and D → 1 2 (u + w) = 1 for all c (and p 1 ) values (Figs. 3 and 5); the same happens in the limit of the burning probability p 1 → 0 (Fig. 4). Increasing p 1 and transition between unbiased and biased diffusion regimes as was argued earlier in Ref. [9]. Fig. 3(b) is similar to the corresponding plot in Ref. [9] for p 1 = 1 case, but for p 1 = 0.1 the gap is prominent only for small p 2 values, and for p 2 > 0.1 it practically disappears. We note that for p 2 = 0 case in Fig. 3(b) the correct c → 0 limit must be D(c → 0) = 2/3 [7].
However, it did not reach this value because of the emerging numerical instability for very low c < ∼ 0.001 (see also our discussion above).
Analysis of Fig. 5(b) shows that as p 2 increases, the behavior of diffusion constant changes from increasing to decreasing function of p 1 , with D(p 1 ) developing a minimum for p 2 < ∼ 0.1. It should be noted that in p 2 = 0 case, D should approach the value of 1/2 for p 1 → 0 according to [7]. We see some discrepancy there which is also due to the numerical instability for small p 1 values. In addition, we observed a gap between D(p 1 → 0) = 1/2 for p 2 = 0 and D(p 1 → 0) = 1 for all nonzero p 2 values.
As another example, we investigated the case of the backward biased motor protein (with specific transition rates u = 0.3, w = 0.7), where bridges are inducing the molecular motor , as it should be. The p 1 → 0 limit in the u < w case is more complex than in the u ≥ w case (when V → u − w, D → 1 2 (u + w) for all p 2 = 0 and c values as p 1 → 0). Namely, for u < w there existsp 2 (u, w) such that for non-zero p 2 <p 2 the dynamic properties V (c) and D(c) exhibit strong c-dependence in the . Physically this implies that bridges with arbitrarily low burning probability strongly affect the dynamics of the particle which tends to move in the backward direction provided that the recovery probability p 2 is less than some critical value; otherwise (in p 2 >p 2 case) weak links have no effect on the particle dynamics as p 1 → 0. Thus in the p 1 → 0 limit there is a dynamic transition at p 2 =p 2 separating the regime with p 2 <p 2 when weak links play a role in the motor protein dynamics and the regime where they are irrelevant (p 2 >p 2 ). It should be noted that the p 1 → 0 limit in Fig. 9 is in agreement with that in Figs

CONCLUSIONS
We have presented a comprehensive theoretical method of calculating dynamic properties of molecular motors in reversible burnt-bridge models for periodic bridge distribution. It is a generalization of the approach developed by us in Ref. [6,7] for the unbiased molecular motors and irreversible burning of bridges. Exact and explicit expressions for mean velocity and dispersion have been derived for arbitrary values of parameters u, w, p 1 , p 2 and c. In the known limiting cases of u = w = 1, p 2 = 0 and of p 1 = 1, we have reproduced our earlier findings [6,7,9], thereby confirming the validity of our theoretical analysis. Some interesting phenomena have been observed as a result of the investigation of dynamic properties of the molecular motor in BBM with bridge recovery. It includes dynamic phase transitions and reversal of the direction of the motion. In case of the unbiased molecular motor, increasing the concentration of bridges c (or lowering the recovery probability p 2 ) with other parameters kept fixed results in increasing velocity and decreasing dispersion. However, dependence of the dispersion on burning probability p 1 is more complex; it is determined by the p 2 value.
In the limit of low c, gaps in dispersion plots have been observed for various p 1 and p 2 values, indicating the dynamic transition between biased and unbiased regimes. Also, a gap was found in the limit of small p 1 between p 2 = 0 and non-zero p 2 regimes. Thus our results obtained in [9] for p 1 = 1 with u = w were generalized to cover the full range of p 1 values.
For the backward biased molecular motor, increasing c has resulted in slowing down the backward movement of the particle, and for sufficiently small p 2 (large p 1 ) the direction of motion has been even reversed and the velocity became positive. In the limit of small p 1 , a dynamic phase transition separating p 1 = 0 and p 1 → 0 regimes has been found provided that p 2 is less than some critical value. For sufficiently small p 2 , broken bridges influence the particle's dynamics even if the burning probability p 1 is infinitesimal. The behavior of dispersion as a function of c was non-monotonic for some range of parameters p 1 and p 2 , with large fluctuations at small c and small p 2 . In the case of irreversible bridge burning (p 2 = 0), we have observed gaps in velocity and dispersion in c → 0 limit (for p 1 = 0.3), with the velocity being positive for all non-zero c values. It suggests that there is a dynamic transition at c = 0 separating backward biased and forward biased diffusion. The velocity and fluctuations are suppressed for sufficiently small c. Hence our findings in [9] for u < w case with p 1 = 1 have been extended to describe the general case of 0 < p 1 ≤ 1.
The method presented above applies to the case of periodic distribution of weak bonds, which is probably realistic for collagenases [2]. As a problem to be addressed in the future studies, one can consider BBM with random distribution of bridges [4,5] where a different theoretical approach must be applied. [