2.1 Transforming the fingerprints into a decodable form

Complex biomolecules, i.e. biomolecules that unfold or dissociate via intermediate states rather than in a simple two-state fashion, are characterized by equally complex mechanical fingerprints. Individual rips in the fingerprints reflect the molecule crossing the corresponding barriers on its multi-barrier energy landscape. However, extracting the parameters of the barriers, as well as the timescales associated with these barriers, from the fingerprints of a multi-barrier system is not a trivial task. The challenge originates from the fact that the transition force at a given rip (say, rip number i) is not determined by the barrier number i alone, but also by all the preceding barriers that the molecule had to cross on the way to the ith barrier. This dependence of measured forces on the history of the multi-barrier transition arises from the time-dependent nature of the applied force, or the fact that, while the molecule wanders among the intermediate barriers on the way to the ith barrier, the external force keeps changing. Consider, for example, the ramping protocol, in which the applied force is increased with time. If the molecule is significantly delayed in the intermediate states because of a high barrier and/or simply because of ‘bad luck’ (i.e. lack of favorable fluctuations), by the time the molecule arrives to the ith barrier, the force has already been ramped up significantly. As the result, the measured force upon the escape over the ith barrier will be higher than it would have been in the absence of preceding barriers.

The consequence of this interplay between the changing force that perturbs the energy landscape and the intrinsic timescales set by the barriers on this landscape is that the measured transition forces for individual barriers contain entangled contributions from other barriers. Thus, in order to interpret the measured transition forces quantitatively, it would be helpful to find a way to disentangle the contributions from individual barriers in the measured forces. A simple solution to this problem can be found by exploiting the quasi-adiabatic approximation (detailed in Section 3.1.1) of the force-ramp experiments (Dudko et al. Reference Dudko, Hummer and Szabo2006, Reference Dudko, Hummer and Szabo2008; Evans et al. Reference Evans, Halvorsen, Kinoshita, Wong, Hinterdorfer and van Oijen2009; Raible et al. Reference Raible, Evstigneev, Reimann, Bartels and Ros2004). The solution has the form of a mathematical transformation (Zhang & Dudko, Reference Zhang and Dudko2013):

(1a) $$\matrix{ {k_{ij} \left( F \right) = \left \vert {\dot F_i (F)} \right \vert \displaystyle{{P_{ij} (F)} \over {{\rm {\cal N}}_i (F)}},} & {\forall \; i \to j} \cr}. $$

Here, P _ij (F) is the force histogram for the transition from state i to state j, ${\rm {\cal N}}_i (F)$ is the population and $\vert {\dot F_i (F)} \vert$ is the absolute value of the force-loading rate in state i at force F (Fig. 2). With the characteristics of the mechanical fingerprints as the input on the right-hand side, this transformation yields, as an output on the left-hand side, the force-dependent transition rates for all the transitions that can possibly occur in the system. The resulting full set of force-dependent rates constitutes ‘the rate map’ of the system (Fig. 3) (Zhang & Dudko, Reference Zhang and Dudko2013).

Fig. 2.

Applying the transformation in Eqs. (1a , b ) to the mechanical fingerprints: step-by-step illustrations. (a) Free-energy profile of a system with two sequential barriers. Indicated are the intrinsic transition rates (each rate is the inverse of the corresponding timescale) and the barrier heights and locations, which are the parameters sought to be reconstructed. (b) Two representative force-extension curves from a stretching and relaxation cycle. Indicated are the transition forces for different types of transitions. (c) Force histograms for the different types of transitions. The histograms are collected from the force–extension curves from stretching (left column) and relaxation (right column) cycles at the nominal loading rates indicated (in piconewtons per second) next to each histogram. Dividing the raw number of counts in the bin by the corresponding bin width yields P _ij (F) in Eq. (1a ) for a transition from state i to state j. (d) Finding ${\rm {\cal N}}_i (F)$ , the number of molecules (trajectories) in the ith state at force F. (e) Determining ${\dot F}_i (F)$ , the force-loading rate in different states, from the slopes of the force trajectory in these states.

Fig. 3.

The transformation at work: from mechanical fingerprints to the rate map to activation barriers. Left: A system with two sequential barriers (top) and the corresponding rate map (bottom) obtained via Eqs. (1a , b ). The four branches on the rate map correspond to the four possible transitions on the two-barrier potential. The colors of the arrows on the potential graph correspond to different types of transitions and match the colors of the corresponding branches on the rate map. Middle: A system with a cascade of barriers and eight possible transitions. Right: A system with two competing pathways, one of which features an intermediate. To test the robustness of the transformation, the effect of an anharmonic linker was incorporated in the simulations of the systems in the left and middle panels. Error bars are calculated with Eq. (2). Lines are the fit to Eq. (3). As an example, Table 1 lists the heights and locations of the barriers and the associated rates extracted from the fit for the system in the left panel.

Let us summarize the key aspects of the transformation in Eq. (1a ), both of fundamental and of practical nature:

1. (1) The transformation is model-free: it makes no model assumptions regarding the shape of the free-energy barriers or the functional form of the force-dependent rates. The model-free nature of the transformation stems from its ability to implement, with no model assumptions, the statistical–mechanical ‘flux-over-population’ formulation of the transition rate (van Kampen, Reference Van Kampen2007; Zhou, Reference Zhou2010).

2. (2) The transformation is quite general, as it is applicable to transitions that involve arbitrarily large number of barriers and reaction pathways, as long as the individual transitions can be resolved as rips in the force–extension curves.

3. (3) The transformation is easy to implement in practice, as it only involves the quantities that can be directly accessed in force spectroscopy experiments. Figure 2 provides a step-by-step illustration of the implementation. To highlight the ease of the implementation, we can rewrite the transformation in the form of a word equation:

   (1b) $${\rm rate}_{{\rm state}\; i \to {\rm state}\; j} \left( {{\rm force}\,F} \right) = \displaystyle{{\left \vert {{\rm loading}\;{\rm rate}\;{\rm in}\;{\rm state}\;\; i\,\; {\rm at}\,{\rm force}\,F} \right \vert} \over {{\rm number} \,\,{\rm of}\;{\rm trajectories\;} {\rm in\;} {\rm state}\,\; i\,\; {\rm at}\,\; {\rm force}\,\; F}} \times \displaystyle{{{\rm number}\,\,{\rm of}\;{\rm counts}\,{\rm in}\,{\rm bin}\,F} \over {{\rm width\;} {\rm of}\,{\rm bin}\,\; F}}.$$

4. (4) The transformation is applicable to the fingerprints obtained both in a ramping (stretching) and in a relaxation force protocol.

5. (5) The transformation naturally accounts for the reversible behavior, or back-and-forth ‘hopping’ between states. Hopping is often ignored in the analytical treatments of force–extension relationships because it tends to complicate the equations for the experimental observables. In contrast, hopping behavior in the mechanical fingerprints works to one's advantage when using the transformation in Eqs. (1a , b ): hopping is equivalent to repeatedly sampling a given transition, which provides additional data points for the transition forces that can be used as the input in the transformation. For example, hopping yields two values of F _NI in a single force–extension curve in Fig. 2b .

6. (6) The transformation easily accounts for the dependence of the force-loading rate on the force itself. Such dependence, which typically complicates an analytical treatment, can result from the nonlinear force–extension relationship of the linker that tethers the molecule under study to the pulling device. In the transformation, however, the force-loading rate ${\dot F}(F)$ appears as an independent factor. This factor can be readily determined from the slope of the force trajectory in the state of interest, as illustrated in Fig. 2e .

7. (7) Even though force-dependent rates can, in principle, be measured in a force clamp experiment, the immediate outcome of the transformation – the rate map – yields the rates of the entire spectrum of transitions over a broad range of forces (Fig. 3). The covered force range is typically beyond what can be accessed directly in a force clamp. The coverage of the broad range of forces is achieved through the use of data from force ramp experiments that are, in turn, performed in a broad range of loading rates, with each loading rate effectively ‘probing’ the corresponding range of forces.

In practice, the number of measurements is, of course, finite, which contributes to the statistical errors in the force histograms P _ij (F) and in the population ${\rm {\cal N}}_i (F)$ in a given state. How do the errors in the histograms and populations propagate into the rates that are obtained by transforming the histograms with Eqs. (1a , b )? The expression for the standard deviation of the rates on the rate map can be found by computing the standard deviations of each component in Eqs. (1a , b ) and subsequently applying the standard error propagation rule. The resulting standard deviation in the logarithmic rate k _ij (F _m ) for transitions from state i to state j at the force F _m is (Zhang & Dudko, Reference Zhang and Dudko2013)

(2) $$\sigma _{{\rm ln}k_{ij} (F_m )} \approx \left[ {\displaystyle{1 \over {P_{ij} (F_m )\Delta F_m}} + \displaystyle{1 \over {{\rm {\cal N}}_i (F_m )}}} \right]^{1/2}, $$

where P _ij (F _m )ΔF _m is the number of counts in the mth bin with the bin width ΔF _m , and ${\rm {\cal N}}_i (F_m )$ is the number of trajectories found in the state i at the force F _m (Fig. 2d ). The vertical error bars seen in the rate maps in Fig. 3 have been determined with Eq. (2).

2.2 Extracting rates and rate-limiting barriers

The effect of the transformation in Eqs. (1a , b ) is to decompose the response of a multi-barrier system to force into contributions from the individual barriers of the system. Specifically, Eqs. (1a , b ) convert complex mechanical fingerprints into a form (the rate map) that readily reveals the barriers and timescales of the system. Indeed, each branch on the rate map (Fig. 3) represents the force-dependent transition rate k _ij (F) for a particular barrier, namely, the barrier separating states i and j. Consequently, each branch on the rate map can be analyzed as a single-barrier problem. To carry out this analysis, the following expression for the force-dependent rate can be used to fit the individual branches on the rate map (Zhang & Dudko, Reference Zhang and Dudko2013):

(3) $$\eqalign{k(F) = & k_0 \left[ {1 + \displaystyle{{\nu \kappa x^{{\ddag} 2}} \over {2 \Delta G^{{\ddag}}}} \mp \displaystyle{{\nu Fx^{{\ddag}}} \over {{\Delta} G^{{\ddag}}}} \left( {1 + \displaystyle{{(1 - \nu )\kappa x^{{\ddag} 2}} \over {2 \Delta G^{{\ddag}}}}}\right) }\right]^{1/\nu - 1} \cr & \times \exp \left\{ {\beta {\Delta} G^{{\ddag}} \left[ {1 - \left[ {1 + \displaystyle{{\nu \kappa x^{{\ddag} 2}} \over {2 \Delta G^{{\ddag}}}} \mp \displaystyle{{\nu Fx^{{\ddag}}} \over {{\Delta} G^{{\ddag}}}} \left( {1 + \displaystyle{{(1 - \nu )\kappa x^{{\ddag} 2}} \over {2 \Delta G^{{\ddag}}}}} \right) }\right]^{1/\nu}} \right]} \right\},}$$

where β = (k _B T)⁻¹ is the inverse of the product of the Boltzmann constant and the absolute temperature, ∓ applies to forward/backward transitions in both the stretching and relaxation cycles, and κ is the spring constant of the pulling device. Equation (3) has been derived by solving the Kramers problem (Kramers, Reference Kramers1940) on a force-dependent one-dimensional (1D) free-energy profile. The unifying nature of the solution in Eq. (3) is evident from the scaling factor ν, which specifies different models of the free-energy profile: a smooth barrier that is described by a linear–cubic polynomial corresponds to ν = 2/3, while a cusp-like barrier with a harmonic well corresponds to ν = 1/2. A simple least-squares fit of each branch on the rate map to Eq. (3) yields the intrinsic parameters for the height ΔG ^‡ and location x ^‡ of the corresponding activation barrier and the transition rate k ₀ associated with this barrier. The fit (with the value of ν fixed at 1/2) is shown as lines in Fig. 3. Table 1 lists the parameters extracted from the fit for the system with two sequential barriers.

Table 1.

Intrinsic rates (in s⁻¹), barrier heights (in k_BT) and barrier distances (in nm) extracted from the fit of the rate map to Eq. (3). The underlying free-energy profile with two sequential barriers is shown in Fig. 3 (top left) along with the corresponding rate map (bottom left, fit is shown as lines).

Equation (3) is accurate beyond the approximation (F ≈ κVt) typically used for the applied force F = κ(Vt − x). This higher degree of accuracy makes Eq. (3) applicable not only to transitions that originate in ‘stiff’ states (i.e. states characterized by narrow minima on the free-energy landscape) but also transitions that originate in ‘soft’ states, such as the unfolded state, for which the abovementioned approximation may not be sufficiently accurate. Detailed analysis shows that Eq. (3) is valid for the biomolecular stiffness values as low as ~5κ, where κ is, as before, the pulling spring constant. Given the typical spring constant of an AFM ~5 pN nm⁻¹ and an optical trap ~0.3 pN nm⁻¹, and with an additional softening effect of the linker, Eq. (3) applies to most experimental situations in single-molecule force spectroscopy of biomolecules (Zhang & Dudko, Reference Zhang and Dudko2013).

The force-dependent rates obtained via Eqs. (1a , b ) from the force ramp should, of course, agree with the rates measured directly under constant force (i.e. in the force clamp). However, when comparing the rates obtained in these two pulling modes, it has to be taken into account that the force clamp and the force ramp have different effects on the molecular potential. The bias imposed by the force clamp on the molecular potential is described by the term −Fx and is thus linear in the extension x, while the corresponding bias in the force ramp is described by 1/2κ(x − Vt)² and is thus non-linear in x. To account for this difference, the data point for the rate k(F) measured at constant force F should be transferred to the rate map as k(F′), where $F^{\prime} = (F \pm \kappa x^{\ddag} /2)/(1 + (1 - \nu )\kappa x^{\ddag 2} /(2 \Delta G^{\ddag} ))$ . Here ‘+’ applies to forward transitions (i.e. those toward larger values of x, such as the unfolding) and ‘–’ to backward transitions (such as the refolding).

3.1 Irreversible transitions over a single barrier

The case when conformational transitions involve a single rate-limiting barrier is worth special consideration. In this case, the problem of force-induced molecular transitions can be solved analytically for the major outputs of both the force clamp and force ramp experiments. The derivation can be carried out within the 1D picture where the molecular extension x serves as a reaction coordinate. This case covers two-state systems, in which the folded and unfolded states are separated by a single barrier. In multi-state systems, this case applies to transitions over the first barrier.

3.1.1 Unfolding and unbinding

Force-induced conformational transitions can be treated as diffusive escape over a free-energy barrier in a bias field (Evans & Ritchie, Reference Evans and Ritchie1997; Izrailev et al. Reference Izrailev, Stepaniants, Balsera, Oono and Schulten1997; Shapiro & Qian, Reference Shapiro and Qian1997). This treatment is an extension of the view of conformational transitions as Brownian motion over a potential barrier, which is the idea underlying Kramers theory of chemical reactions (Kramers, Reference Kramers1940; van Kampen, Reference Van Kampen2007; Zwanzig, Reference Zwanzig2001). Posing and solving Kramers problem – that of the rate at which a Brownian particle escapes from a potential well – on a class of model potentials that feature a single barrier, the following analytical solution has been derived for the rate k(F) at constant force F (Dudko et al. Reference Dudko, Hummer and Szabo2006):

(4) $$k(F) = {K_0} \left({1 - \displaystyle{{vF{x^{{\ddag}}}} \over {\Delta {G^{{\ddag}}}}}}\right)^{1/v - 1} \exp \left\{{ \beta \Delta {G^{{\ddag}}}\left[{1 - \left({1 - \displaystyle{{vF{x^{{\ddag}}}} \over {\Delta {G^{{\ddag}}}}}}\right)}\right]}\right\}.$$

The predictive value of Eq. (4) lies in its capacity to relate a measurable quantity – the rate k(F) at force F – to the intrinsic (zero-force) characteristics of the system: the height ΔG ^‡ and location x ^‡ of the activation barrier and the unfolding (or unbinding) rate k ₀. These intrinsic parameters can be readily extracted with Eq. (4) through a least-squares fit of the rates measured in a force clamp. As in the case of Eq. (3), the unifying nature of Eq. (4) is evident from the scaling factor ν, which specifies different models of the 1D free-energy profile. As a test of relative insensitivity of the extracted parameters to the details of the shape of the barrier, the force-dependent rates can be fit to Eq. (4) with the values of ν fixed at 1/2 and 2/3, and the resulting values of ΔG ^‡, x ^‡ and k ₀ compared. A notable property of the solution in Eq. (4) (as well as of its generalized version, Eq. (3)) at physically meaningful values of ν (1/2 and 2/3) is a characteristic concave down curvature in the rate at intermediate-to-high forces on a semi-log scale (Fig. 4a ). This curvature reflects the intrinsic property of the potential barrier and well to move toward each other with an increase in the external force. In contrast, this property is ignored in the Bell's postulate for the force-dependent rate (Bell, Reference Bell1978). Bell's expression is recovered in Eq. (4) with ν = 1.

Fig. 4.

Unfolding/unbinding transitions over a single barrier under applied force, from simulations (symbols and histograms) and analytical theory (black lines). (a) The force-dependent rate from the force-clamp pulling mode. Note the characteristic nonlinearity in the rate on a semi-log plot. Line is Eq. (4) with ν = 2/3. (b) A representative force–extension curve from the force–ramp pulling mode. The ‘rip’ signals an unfolding/unbinding event, the corresponding unfolding/unbinding force value is circled. From multiple repeats of this experiment, histograms of the unfolding forces (see inset) are collected. (c) Distributions of unfolding/unbinding forces from the force–ramp mode at three values of the loading rate. Note the characteristic negative skew in the distributions. Lines are Eq. (5) with ν = 2/3. (d) Average force as a function of the loading rate, or ‘the force spectrum’. Note the characteristic upward curvature. Line is Eq. (6) with ν = 2/3.

The functional form of the rate versus force, established in Eq. (4), can be used to derive an analytical expression for the distribution of rupture forces, the key output of a force ramp experiment (Fig. 4b, c ). Under the assumption that the external force is ramped up with time slowly enough such that a barrier of at least a few k _B T is present at all rupture events (the quasi-adiabatic approximation), which reflects a typical experimental situation, the following expression has been derived (Dudko et al. Reference Dudko, Hummer and Szabo2006):

(5) $$p(F) = \displaystyle{{k(F)} \over {\kappa V}}{\rm exp}\left( {\displaystyle{{k_0} \over {\beta x^{\ddag} \kappa V}}} \right){\rm exp}\left[ { - \displaystyle{{k(F)} \over {\beta x^{\ddag} \kappa V}}\left( {1 - \displaystyle{{\nu Fx^{\ddag}} \over {{\Delta} G^{\ddag}}}} \right)^{1 - 1/\nu}} \right].$$

An analytical form of the solution in Eq. (5) makes it convenient for a fit of transition-force histograms obtained in force-ramp experiments. A fit (Fig. 4c , lines) readily yields the values of the height ΔG ^‡ and location x ^‡ of the barrier and the rate k ₀ in the absence of a force. The accuracy of the extracted parameters is enhanced in a global fit. A global fit is fit of multiple histograms at once with a single set of fitting parameters {k ₀, ΔG ^‡, x ^‡}, where the histograms are measured over a broad range of values of the pulling speed.

The force distribution in Eq. (5) has a characteristic negative skew (Fig. 4c ). The longer low-force tail of the distributions reflects the rare unfolding events that occur at short times – and hence at low forces – when the barrier is still high. Increasing the loading rate shifts the distribution to higher forces (Fig. 4c ). This trend in the shift is quite intuitive: when the applied force increases rapidly (which happens at high loading rates), there is less time for thermal fluctuations to kick the molecule over a high barrier, causing transitions to occur later in time – which means, at higher forces (hence the shift to higher forces) – when the barrier is lower.

The probability distribution of rupture forces in Eq. (5) has been used to derive an analytical expression for the mean rupture force as a function of the loading rate κV, or ‘the force spectrum’ (Dudko et al. Reference Dudko, Hummer and Szabo2006):

(6) $$\left\langle F \right\rangle \cong \displaystyle{{{\Delta} G^{\ddag}} \over {\nu x^{\ddag}}} \left[ {1 - \left( {\displaystyle{1 \over {\beta {\Delta} G^{\ddag}}} {\rm ln}\displaystyle{{k_0 \,{\rm e}^{\beta {\Delta} G^{\ddag} + \gamma}} \over {\beta x^{\ddag} \kappa V}}} \right)^\nu} \right],$$

where γ = 0.577… is the Euler–Mascheroni constant. The expression in Eq. (6) predicts an increase in the mean rupture force with the loading rate as well as a convex upward curvature in the force spectrum (F versus ln(κV)) at intermediate-to-high loading rates on the semi-log scale at 1/2 > ν > 2/3 (Fig. 4d ). This curvature is a manifestation of the force-induced motion of the transition barrier toward the potential well on the energy landscape, as mentioned above following Eq. (4).

3.1.2 Folding and binding

Equation (5) describes the distribution of unfolding or unbinding forces that are measured in force–ramp experiments, in which the force is increased with time until an unfolding event is detected. Likewise, in force-relaxation experiments, the force is decreased with time, which allows measurements of the reverse processes, the folding or binding (Fig. 5b ). The expression for the folding or binding force distribution, which is analogous to Eq. (5) for unfolding or unbinding, has been derived as (Pierse & Dudko, Reference Pierse and Dudko2013)

(7) $$p(F) = \displaystyle{{k_ \leftarrow (F)} \over {\left \vert {\dot F(F)} \right \vert}} {\rm exp}\left[ { - \displaystyle{{k_ \leftarrow (F)} \over {\beta \left \vert {\dot F(F)} \right \vert x^{\ddag}}} \left( {1 + \displaystyle{{\kappa _S} \over {\kappa _U (F)}}} \right)^{\nu /(\nu - 1)} \left( {1 + \displaystyle{{\nu Fx^{\ddag}} \over {{\Delta} G^{\ddag}}}} \right)^{1 - 1/\nu}} \right].$$

Fig. 5.

Folding/binding transition over a single barrier against applied force, from simulations (symbols and histograms) and analytical theory (lines). (a) The force-dependent folding rate from the force–clamp pulling mode. Note the characteristic nonlinearity on a semi-log plot. Line is Eq. (4) (with x ^‡ < 0 for the folding or binding process) evaluated at ν = 2/3. (b) A representative force–extension curve from a force-relaxation experiment with the value of the folding force circled. From multiple repeats of this experiment, histograms of the folding forces are collected (inset). (c) Distributions of folding/binding forces from the force-relaxation regime at three values of the relaxation rate. Note the characteristic positive skew in the distributions. Lines are Eq. (7) with ν = 2/3.

The force-dependent folding rate that enters this expression is given by

$$k_\leftarrow (F) = k_0 \left( {1 + \displaystyle{{\kappa _S} \over {\kappa _U (F)}}} \right)^{1/\nu - 1/2} \left( {1 + \displaystyle{{\nu Fx^{\ddag}} \over {{\Delta} G^{\ddag}}}} \right)^{1/\nu - 1} {\rm exp}\left\{ {\beta {\Delta} G^{\ddag} \left[ {1 - \left( {1 + \displaystyle{{\kappa _S} \over {\kappa _U (F)}}} \right)^{2\nu /(1 - \nu ) - 1} \left( {1 + \displaystyle{{\nu Fx^{\ddag}} \over {{\Delta} G^{\ddag}}}} \right)^{1/\nu}} \right]} \right\},$$

where the stiffness of the unfolded biomolecule at force F and the force-relaxation rate are, respectively,

$$\kappa_U (F) = \displaystyle{\nu \over {(1 - \nu)^2}} \displaystyle{{{\Delta} G^{\ddag}} \over {x^{{\ddag}2}}} \left({1 + \displaystyle{{\nu Fx^{\ddag}} \over {{\Delta} G^{\ddag}}}} \right)^{2 - 1/\nu}, \quad {\dot F}(F) = \displaystyle{{-\kappa_S V} \over {1 + \kappa_S /\kappa_U (F)}}.$$

The solution in Eq. (7) reveals properties of the distribution of folding forces that distinguish this distribution from its unfolding counterpart, Eq. (5). The distinguishing properties of these two distributions can also be seen from the comparison between Figs 4c and 5c . First, the distribution of folding/binding forces has a positive skew (Fig. 5c ). The longer tail at higher forces reflects the rare folding events that occur at short times – and hence at high forces – when the barrier is still high. Second, increasing the absolute value of the relaxation rate $\vert {\dot F} \vert$ shifts the distribution to lower forces (Fig. 5c ). This trend can be understood with simple physical arguments: when the applied force decreases rapidly (which is the case at high relaxation speeds), there is less time for thermal fluctuations to kick the molecule over a high barrier, causing transitions to occur later – and hence at lower forces – when the barrier is lower. Third, the stiffness κ _s of the pulling device and speed V enter Eq. (7) separately, providing independent routes to control the range of folding forces and thus enhance the robustness of a fit of experimental distributions to Eq. (7). A global fit (the meaning of a global fit in the present context is discussed in Section 3.1.1) yields the intrinsic height ΔG ^‡ and location x ^‡ of the barrier to folding as well as the on-rate k ₀. The application of the Eq. (7) to binding experiments on a ligand and receptor connected by a tether involves an additional step associated with the decoupling of the effect of the tether and allows the reconstruction of the parameters of ligand–receptor binding (Pierse & Dudko, Reference Pierse and Dudko2013).

3.2 Barely so, but still analytically tractable: two sequential barriers

Section 3.1 described a direct analytical approach to the analysis and interpretation of the data from pulling experiments. This approach takes advantage of the availability of analytical expressions that are suitable for fitting the information-rich characteristics of the conformational dynamics – the distribution of transition forces p(F). The fit, in turn, yields the activation barrier and timescale of the conformational transition. Balancing the generality with simplicity, this direct approach is particularly efficient when we deal with relatively simple systems for which an analytical expression for p(F) can be derived. However, even modest complications of the system, such as the presence of more than one energy barrier, rapidly complicate the corresponding formulae and require additional parameters, increasing the risk of overfitting. Mathematically, the problem is challenging due to the time-varying perturbation caused by the external force F(t). The force introduces the time-dependence in the rate coefficients k _ij (t), which, on the one hand, enriches the behavior of the system but, on the other hand, makes it challenging to solve the corresponding set of rate equations analytically.

Let us consider a case of relatively complex dynamics, for which the direct analytical approach nevertheless leads to expressions simple enough to enable a meaningful fit of the force–ramp data. This is the case of a system in which the native and unfolded states are separated by one intermediate state (Fig. 6). For a broad class of topologies of energy landscapes with an intermediate state, the expression for the distribution of forces at which the system transits into the final (unfolded) state has been derived as (Garai et al. Reference Garai, Zhang and Dudko2014)

(8) $$p(F) = \displaystyle{1 \over {\kappa V}}\left\{ {k_{IU} \left( F \right){\rm e}^{ - g_{IU} (F)} \left[ {1 - {\rm e}^{ - g_{NI} (F^{\prime})}} \right] + k_{NI} (F){\rm e}^{ - g_{NI} (F^{\prime} )} \left[ {1 - {\rm e}^{ - g_{IU} (F)} \left( {1 + g_{IU} (F)} \right)} \right]\left[ {1 + g_{IU}^{ - 1} (F)} \right]} \right\}.$$

Fig. 6.

Signatures of sequential barriers on the free-energy profile. (a) Free-energy profile featuring two sequential barriers: a stiff (narrow) and low barrier is followed by a soft (broad) and high barrier. (b) Gray line: force-dependent effective rate from a force-clamp pulling mode. Also shown are the individual rates for the forward transitions over each of the two barriers (orange and red) and the backward transition over the first barrier (green). (c) Representative force–extension curves from a force-ramp pulling mode at low, intermediate and high loading rates. Arrows indicate the transition events into the unfolded state; the corresponding transition force values are collected in histograms that are shown in (d). The ‘trapping’ effects of the native and intermediate states at different timescales are illustrated in the insets with the populations of particles accumulated in the corresponding ‘trap’. (d) The distributions of unfolding forces at three values of the loading rate. A transient bimodality in the distribution can be seen at intermediate loading rates. Lines are Eq. (8) with ν = 2/3. (e) A visual representation of the origin of bimodality in the force distribution. Shown on the left are the snapshots of the free-energy profile at different instants during a force-ramp experiment. The distribution, P(F), is equivalent to the flux into the final state times the inverse of the force-loading rate ${\dot F}$ . The flux, in turn, is given by the transition rate (governed by the barrier) times the population. The width of the arrows over the barriers reflects the relative value of the flux at different instants. The number of particles in each potential well reflects the relative population in this well at the given instant. Shown on the right are the distributions of forces at the transition into the final state with early transition events contributing to the first peak and late events contributing to the second peak (highlighted in red).

In this expression, the force-dependent rates for transitions between various states (denoted as i and j) are given by

$$k_{ij} (F) = k_{ij}^0 \left( {1 - \displaystyle{{\nu Fx_{ij}^{\ddag}} \over {{\Delta} G_{ij}^{\ddag}}}} \right)^{1/\nu - 1} {\rm exp}\left\{ {\beta \Delta G_{ij}^{\ddag} \left[ {1 - \left( {1 - \displaystyle{{\nu Fx_{ij}^{\ddag}} \over {\Delta G_{ij}^{\ddag}}}} \right)^{1/\nu}} \right]} \right\}$$

and the functions g _ij (F) are

$$g_{ij} (F) = \displaystyle{{k_{ij}^0} \over {\beta \kappa Vx_{ij}^{\ddag}}} \left\{ {{\rm exp}\left[ {\beta \Delta G_{ij}^{\ddag} \left[ {1 - \left( {1 - \displaystyle{{\nu Fx_{ij}^{\ddag}} \over {\Delta G_{ij}^{\ddag}}}} \right)^{1/\nu}} \right]} \right] - 1} \right\}.$$

The above expressions account for a generally non-negligible drop in the force from F′ to F occurring in the course of the transition, upon escape from N to I, due to partial unfolding. These expressions have been derived under the assumption that backward transitions are negligible, which is justified beyond the range of very small forces.

An unusual feature in the distribution of transition forces – namely, bimodality (Fig. 6d ) – is predicted by Eq. (8) if the landscape has a relatively low and stiff barrier that is followed by a high and soft barrier (Fig. 6a ). The bimodality originates from the ‘trapping’ effect of the native and intermediate states at different stages of a force ramp experiment (Fig. 6e ). One of these ‘traps’, the intermediate state, releases its accumulated population at low forces when the second barrier is lowered sufficiently to allow fast escape. Another ‘trap’, the native state, releases the remaining population at high forces when the first, stiff, barrier has been lowered sufficiently.

The expression in Eq. (8) is suitable for a global fit (global fits are discussed in Section 3.1.1) of transition force distributions in a broad range of pulling speeds. The fit yields the locations and heights of each of the two sequential barriers, the location of the intermediate state, and the rate constants for all the transitions on the energy landscape with an intermediate (Garai et al. Reference Garai, Zhang and Dudko2014).