Protein folding as a quantum transition between conformational states

Abstract

Assuming that the main variables in the life processes at the molecular level are the conformation of biological macromolecules and their frontier electrons a formalism of quantum theory on conformation-electron system is proposed. Based on the quantum theory of conformation-electron system, the protein folding is regarded as a quantum transition between torsion states on polypeptide chain, and the folding rate is calculated by nonadiabatic operator method. The rate calculation is generalized to the case of frequency variation in folding. An analytical form of protein folding rate formula is obtained, which can be served as a useful tool for further studying protein folding. The application of the rate theory to explain the protein folding experiments is briefly summarized. It includes the inertial moment dependence of folding rate, the unified description of two-state and multistate protein folding, the relationship of folding and unfolding rates versus denaturant concentration, the distinction between exergonic and endergonic foldings, the ultrafast and the downhill folding viewed from quantum folding theory, and, finally, the temperature dependence of folding rate and the interpretation of its non-Arrhenius behaviors. All these studies support the view that the protein folding is essentially a quantum transition between conformational states.

Graphical abstract

Keywords

protein folding rate / quantum transition / torsion states / non-Arrhenius temperature dependence / exergonic and endergonic folding / ultrafast folding

1 Locality of conformation states and quantum number characterizing conformation-electronic system

Assume that the main variables in the life processes at the molecular level are the conformation of biological macromolecules and their frontier electrons. Suppose that the conformation is described by a set of torsion angles \{\theta \}=({\theta }_1,\mathrm{\dots },{\theta }_N). The conformation potential U, as a function of θ, may generally have some minima that correspond to stable conformations. To be definite, consider the case of two minima. Set

\begin{array}{}\text{(1)}& U(\theta )=\{\begin{array}{c}{\displaystyle \frac{k}{2}}{(\theta +{\displaystyle \frac{\text{π}}{2}})}^2,[-\text{π},\mathrm{\hspace{0.33em}0}]\hfill \\ {\displaystyle \frac{k}{2}}[(1-\lambda ){(\theta -{\displaystyle \frac{\text{π}}{2}})}^2+{\displaystyle \frac{{\text{π}}^2}{4}}\lambda ],[0,\text{π}]\hfill \end{array}\end{array}

Here, λ is an asymmetrical parameter. We shall discuss the dependence of steady-state solutions on the asymmetrical parameters and demonstrate that the conclusions drawn will not depend on the particular choice of the potential when we use any two-minima-asymmetrical potential [1]. If \lambda =0,U(\theta ) is C₂ symmetrical, then the conformational wave functions are nonlocal, and no definite local conformation can be related with the molecule. However, if \lambda \ne 0, the C₂-symmetry is broken; by numerical solution of Schrodinger equation for the potential given by Eq. (1), one can show that for Ik/{\mathrm{\hslash }}^2=10 (I is the inertia moment of the molecule with respect to coordinate \theta ), the ground state should be localized to larger than 90% as \lambda ⩾{10}^{-3}, and for Ik/{\mathrm{\hslash }}^2=40, the ground state localized to larger than 99% as \lambda ⩾{10}^{-6}. Define the locality J_i of the ith level as:

\begin{array}{}\text{(2)}& J_i={\displaystyle \frac{\max (P_i,Q_i)}{\min (P_i,Q_i)}}\end{array}

where P_i is the probability of ith level located in the “left” conformation (-\text{π}⩽\theta ⩽0) and Q_i in the “right” conformation (0⩽\theta ⩽\text{π}). By numerical calculation, one can show that nearly all bound states are localized with J_i>2 as \lambda ⩾0.1(0.01) for the case of Ik/{\mathrm{\hslash }}^2=10(40) [1]. The above discussion shows that the small asymmetry in potential would cause the high localization of wave functions. Although the results are deduced from a particular potential (1), the generality of the conclusions can be demonstrated by the following arguments. Consider the interaction between a pair of levels localized mainly in left potential valley A (the corresponding wave function denoted as {\psi }_A) or in right potential valley B (the corresponding wave function {\psi }_{\mathrm{B}}), respectively. The overlap integral of wave functions \int {\psi }_A^{\ast }{\psi }_{B}d\theta in barrier region is a measure of the potential asymmetry and the wave function localization. By using WKB quasi-classical approximation to solve Schrodinger equation, we obtain the overlap integral proportional to \exp \{-{\int }_{\alpha }^{\beta }K({\theta }^{\prime })d{\theta }^{\prime }\}, where α and β are turning points at the left and right sides of the barrier, respectively, U(a)=U(\beta )=E and the wave vector equals

\begin{array}{}\text{(3)}& K(\theta )={\displaystyle \frac{1}{\mathrm{\hslash }}}{\{2I[U(\theta )-E]\}}^{1/2}\end{array}

For U-E\sim 2.5\text{kJ/mol} (typical torsion energy) and I={10}^{-37}g\cdotcm² (torsion inertial moment of typical atomic group), one has K(\theta )\sim {10}^2. Therefore, the overlap integral is a small quantity, and the localization does occur. Consider that the conformational potential usually has certain asymmetry, we deduce that the levels below the top of barrier are always localized apart from those near the top where K(\theta )\sim 0. As for the problem whether the ground state is localized in the left or right valley, it is determined by the comparison of the depths of two valleys and the values of {\partial }^{2}U/\partial {\theta }^2 in two minima of the potential. As Anderson studied the excitations on a disordered lattice in 1958, he proposed a noted theorem: the eigenfunctions are localized if the strength of the disorder exceeds some definite value [2]. Now, the theorem on the localization of the excitation of a lattice has been generalized to the system of conformations.

Suppose the dynamical variables of the conformation-electronic system are (\theta ,x), where x is the coordinate of frontier electrons, and θ is the torsion angle of molecule. Based on the above analyses on the conformational locality, the wave function of the system can be described by M_{kn\alpha }(\theta ,x), where α denotes the electron-state, and (k, n) refer to the conformation- and vibration-state, respectively. The conformational quantum number k=1,2,\mathrm{\dots }, refers to the local conformational state, indicating which minimum of the torsion potential the system is localized around. The vibration quantum number n indicates the harmonic energy level in a given conformational state.

As a multi-atom system, the conformation of a protein is fully determined by bond lengths r_{ij}, bond angles {\chi }_{i-j-k}, and torsion angles ({\theta }_1,{\theta }_2,\mathrm{\dots },{\theta }_N). Since torsion angles are most easily changed even at room temperature, we assume that the conformation is mainly described by a set of torsion coordinates, that is, the conformational potential U is only a function of ({\theta }_1,{\theta }_2,\mathrm{\dots },{\theta }_N). For a real problem of protein folding, sometimes r_{ij}and {\chi }_{i-j-k} may not be fixed. In these cases, we shall generalize the theoretic approach as follows. We assume that the potential can be averaged over bond lengths and bond angles, and the averaged energy landscape is described by a set of torsion coordinates only. The average is reasonable since the folding does not involve a single microscopic pathway but rather a statistical energy landscape [3]. Through the statistical average of energy landscape and the corresponding redefinition of the potential, the number of torsion-accessible states will be largely increased. For example, for a pair of conformational states that are originally inaccessible, it is impossible to reach each other by changing torsion coordinate (or torsion quantum number in quantum theory) only. However, this is possible after the statistical average of energy landscape and redefinition of the form of torsion potential. Therefore, the statistical energy approach effectively expands the set of torsion-accessible states. In the following, we shall study protein folding as a changing process between a series of torsion-accessible states by use of quantum theory of conformation-electronic system.

2 Protein folding as a quantum transition calculated by nonadiabatic operator method

Protein folding is commonly studied in classical mechanics and classical statistical mechanics. However, this approach is too limited for the understanding of the fundamental physics underlying folding. It is well known that the fluorescence and phosphorescence are phenomena closely related to protein folding. Since no one doubts that fluorescence and phosphorescence could only be understood in terms of the quantum transition between molecules, why should the protein folding study be divorced from the framework of quantum theory? In the following we shall investigate protein folding from the quantum theory of conformation-electronic system.

Protein folding is nonradiative transition. The transition rate can be calculated by nonadiabatic operator method [4]. In adiabatic approximation, the wave function of the conformation-electronic system M_{kn\alpha }(\theta ,x) can be expressed as the product of conformational wave function and electronic wave function, M(\theta ,x)=\psi (\theta )\cdot \phi (x,\theta ). However, the form is not a rigorous eigenstate of the total Hamiltonian of the system, and there exist transitions between adiabatic states that result from the off-diagonal elements. The transition matrix element is given by off-diagonal element as follows:

\begin{array}{}\text{(4)}& ⟨k^{\prime }{n}^{\prime }{\alpha }^{\prime }|H^{\prime }|kn\alpha ⟩={\displaystyle \int }{\psi }_{k^{\prime }{n}^{\prime }{\alpha }^{\prime }}^{+}(\theta )[-{\displaystyle \frac{{\mathrm{\hslash }}^2}{2I}}{\displaystyle \int }{\phi }_{{\alpha }^{\prime }}^{+}\cdot ({\displaystyle \frac{{\partial }^2{\phi }_{\alpha }}{\partial {\theta }^2}}+2{\displaystyle \frac{\partial {\phi }_{\alpha }}{\partial \theta }}{\displaystyle \frac{\partial }{\partial \theta }}\left){\mathrm{d}}^{3}x\right]\psi {}_{kn\alpha }(\theta )\mathrm{d}\theta \end{array}

Consider that the electronic state does not change in transition processes, a^{\prime }=a. This is the generally occurring case for most protein foldings. For the case of radiation transition, the electronic transitions exist, accompanying conformational changes. We will discuss these phenomena from the point of conformation-electronic transition in another article. Since the electronic wave function {\phi }_{\alpha } is generally real, one has

\begin{array}{}\text{(5)}& {\displaystyle \int {\phi }_{\alpha }(x,\theta )\frac{\partial {\phi }_{\alpha }(x,\theta )}{\partial \theta }{\mathrm{d}}^{3}x}=0\end{array}

from the normalization condition \int {\phi }_{\alpha }(x,\theta ){\phi }_{\alpha }(x,\theta ){\mathrm{d}}^{3}x=1. Therefore, only the first term in right-hand-side of Eq. (4) is retained, namely

\begin{array}{}\text{(6)}& ⟨k^{\prime }{n}^{\prime }\alpha |H^{\prime }|kn\alpha ⟩={\displaystyle \int {\psi }_{k^{\prime }{n}^{\prime }\alpha }^{+}(\theta )\left(-\frac{{\mathrm{\hslash }}^2}{2I}\int {\phi }_{\alpha }^{+}\frac{{\partial }^2{\phi }_{\alpha }}{\partial {\theta }^2}{\mathrm{d}}^{3}x\right){\psi }_{kn\alpha }(\theta )d\theta }\end{array}

For a multi-torsion-angle problem, the torsion potential U_{\mathrm{tor}}({\theta }_1,{\theta }_2,\mathrm{\dots },{\theta }_N) is generally argument-inseparable. Only when its conformational wave function can be approximately expressed by the product of single conformational functions, then the quantum transition problem can easily be solved. However, from the widely accepted statistical energy landscape theory [3], the protein folding involves a statistical average of a lot of microscopic pathways. One may average U_{\mathrm{tor}}({\theta }_1,{\theta }_2,\mathrm{\dots },{\theta }_N) over N-1 coordinates with only one angular coordinate {\theta }_j left as the variable. The resulting potential is denoted as U_{\mathrm{tor}}^{(j)}({\theta }_j)

\begin{array}{}\text{(6)}& U_{\mathrm{tor}}^{(j)}({\theta }_j)={⟨U_{\mathrm{tor}}({\theta }_1,{\theta }_2,\mathrm{\dots },{\theta }_N)⟩}_{{\theta }_1,\mathrm{\dots },{\theta }_{j-1},{\theta }_{j+1},\mathrm{\dots },{\theta }_N}\end{array}

The summation of U_{\mathrm{tor}}^{(j)}({\theta }_j) over j gives an approximate expression of U_{\mathrm{tor}}({\theta }_1,{\theta }_2,\mathrm{\dots },{\theta }_N), and the argument-separable wave function is obtained. Therefore, the argument-separable torsion wave function as an approximation can be deduced from the statistical energy landscape approach. In argument-separable case, one has

\begin{array}{}\text{(7)}& {\psi }_{kn\alpha }(\theta )={\psi }_{k_1,n_1,{\alpha }_1}({\theta }_1),\mathrm{\dots },{\psi }_{k_N,n_N,{\alpha }_N}({\theta }_N)\end{array}

where {\psi }_{k_j,n_j,{\alpha }_j}({\theta }_j) is approximately expressed by a wave function of harmonic oscillator with quantum number n_j. Note that the harmonic potential has equilibrium position at {\theta }_j={\theta }_{kj}^{(0)} with the corresponding k_jth minimum of potential E_{kj} (k_j=1,2,\mathrm{\dots }). Inserting Eq. (7) into Eq. (6), one obtains

\begin{array}{}\text{(8)}& ⟨k^{\prime }{n}^{\prime }\alpha |H^{\prime }|kn\alpha ⟩={\displaystyle \int {\psi }_{k^{\prime }{n}^{\prime }\alpha }^{+}(\theta )\sum \left(-\frac{{\mathrm{\hslash }}^2}{2I_j}\int {\phi }_{\alpha }^{+}\frac{{\partial }^2{\phi }_{\alpha }}{\partial {\theta }_j^2}{\mathrm{d}}^{3}x\right){\psi }_{kn\alpha }(\theta )\mathrm{d}\theta }\end{array}

After thermal average of the initial states, the transition rate will be

\begin{array}{}\text{(9)}& W={\displaystyle \frac{2\text{π}}{\mathrm{\hslash }}}{\displaystyle \sum _{\{n\}}}{\left|⟨k^{\prime }{n}^{\prime }\alpha |H^{\prime }|kn\alpha ⟩\right|}^{2}B(\{n\},T)\cdot \delta \left({\displaystyle \sum _j}(n_j^{\prime }\mathrm{\hslash }{\omega }_j^{\prime }-n_j\mathrm{\hslash }{\omega }_j-\delta E_j)\right)\end{array}

Here, {\omega }_j and {\omega }_j^{\prime }are frequency parameters of the jth mode in conformational state k_j and k_j^{\prime }, respectively, and \delta E_j=E_{kj}-E_{k^{\prime }j}. B(\{n\},T) denotes the Boltzmann factor:

\begin{array}{}\text{(10)}& B(\{n\},T)={\displaystyle \prod _j}B(n_j,T)={\displaystyle \prod _j}\left[{\mathrm{e}}^{-n_j\beta \mathrm{\hslash }{\omega }_j}(1-{\mathrm{e}}^{-\beta \mathrm{\hslash }{\omega }_j})\right]\end{array}

After summing over final states, we have

\begin{array}{}\text{(11)}& W={\displaystyle \frac{2\text{π}}{\mathrm{\hslash }}}{\displaystyle \sum _{\{n\}}}{\left|⟨k^{\prime }{n}^{\prime }\alpha |H^{\prime }|kn\alpha ⟩\right|}^{2}B(\{n\},T){\rho }_E\end{array}

Here, {\rho }_E means state density,

\begin{array}{}\text{(12)}& \begin{array}{c}{\rho }_E={\displaystyle \frac{1}{\partial E_f/\partial N_f}}={\displaystyle \frac{1}{\mathrm{\hslash }{\overline{\omega }}^{\prime }}},N_f=\sum _j{n}_j^{\prime }\hfill \\ E_f={\displaystyle \sum _j}(n_j^{\prime }\mathrm{\hslash }{\omega }_j^{\prime }+E_{k^{\prime }j})\hfill \end{array}\end{array}

\overline{{\omega }^{\prime }} is the average of {\omega }_j^{\prime } over j.

3 Protein folding rate in single torsion transition

If only one torsion angle participates in the conformational transition, we call it single-mode or single torsion transition; if several torsion angles participate simultaneously in one conformational transition, then we call it multimode or multitorsion transition. Consider single-mode case at first. In this case, the subscript j will be dropped. The transitional rate

\begin{array}{}\text{(13)}& W={\displaystyle \frac{2\text{π}}{{\mathrm{\hslash }}^2{\omega }^{\prime }}}{I}_E{I}_V\end{array}

\begin{array}{}\text{(14)}& I_E={\left|{\displaystyle \frac{-{\mathrm{\hslash }}^2}{2I}}{\displaystyle \int {\phi }_{\alpha }\frac{{\partial }^2}{\partial {\theta }^2}{\phi }_{\alpha }{\mathrm{d}}^{3}x}\right|}_{\theta ={\theta }_0}^2\end{array}

\begin{array}{}\text{(15)}& I_V={\displaystyle \sum _n}{\left|{\displaystyle \int {\psi }_{k^{\prime }{n}^{\prime }\alpha }^{+}(\theta ){\psi }_{kn\alpha }(\theta )d\theta }\right|}^{2}B(n,T)\end{array}

{\theta }_0 is the conformational coordinate value in the largest overlap region of vibration functions. Eq. (14) can be estimated by the square of rotational kinetic energy of the frontier electron during conformational change (see next section). The quantum number n^{\prime } in Eq. (15) is determined by n through energy conservation. The overlap integral in Eq. (15) was first calculated by Manneback (1951) [5, 6]. For simplicity, we shall assume the same frequency for initial and final states (\omega ={\omega }^{\prime }) first. After thermal average, it gives

\begin{array}{}\text{(16)}& I_V={\left({\displaystyle \frac{\overline{n}+1}{\overline{n}}}\right)}^{p/2}{J}_p(2Q\sqrt{\overline{n}(\overline{n}+1)})\exp (-Q(2\overline{n}+1))\end{array}

in which

\begin{array}{}\text{(17)}& \overline{n}={({\mathrm{e}}^{\beta \mathrm{\hslash }\omega }-1)}^{-1}\end{array}

\begin{array}{}\text{(18)}& Q=I\omega {(\delta \theta )}^2/(2\mathrm{\hslash }),p={\displaystyle \frac{\delta E}{\mathrm{\hslash }\omega }}\end{array}

\delta \theta ={\theta }_k^{(0)}-{\theta }_{k^{\prime }}^{(0)} is the angular displacement, and \delta E=E_k-E_{k^{\prime }} the energy gap between initial and final states. J_P denotes the modified Bessel function [7], and here, p is related to the net change in oscillator quantum number. The modified Bessel function is introduced from the expansion:

\begin{array}{}\text{(19)}& \exp \left({\displaystyle \frac{x}{2}}\left(y+{\displaystyle \frac{1}{y}}\right)\right)={\displaystyle \sum _{n=-1}^{\mathrm{\infty }}}{y}^n{J}_n(x)\end{array}

By use of the asymptotic formula for Bessel function [8]

\begin{array}{}\text{(20)}& {\mathrm{e}}^{-z}{J}_p(z)={(2\text{⁢π}z)}^{-1/2}\exp \left(-{\displaystyle \frac{p^2}{2z}}\right),\mathrm{for}z\gg 1\end{array}

I_V can be further simplified. Note that the equality \overline{n}+\frac{1}{2}=\sqrt{\overline{n}(\overline{n}+1)} holds approximately with error 2% as \overline{n}=2 and smaller than 1% as \overline{n}>3. Considering \overline{n}>1 (), from Eqs. (16) and (20), it leads to

\begin{array}{}\text{(21)}& I_V={(2\text{⁢π}z)}^{-1/2}\exp \left(-{\displaystyle \frac{p^2}{2z}}\right)\exp \left({\displaystyle \frac{\delta E}{2k_{\mathrm{B}}T}}\right)\end{array}

where

\begin{array}{}\text{(22)}& z={(\delta \theta )}^2\cdot {\displaystyle \frac{k_{\mathrm{B}}T}{{\mathrm{\hslash }}^2}}I\end{array}

4 Protein folding rate in multitorsion transition

For multitorsion transition, from Eqs. (11) and (12), one has

\begin{array}{}\text{(23)}& W={\displaystyle \frac{2\text{π}}{{\mathrm{\hslash }}^2{\overline{\omega }}^{\prime }}}{I}_E{\displaystyle \sum _{\{p_j\}}}{\displaystyle \prod _j}{I}_{Vj}\end{array}

\begin{array}{}\text{(24)}& I_E={\left|{\displaystyle \sum _j}{\displaystyle \frac{-{\mathrm{\hslash }}^2}{2I_j}}{\displaystyle \int {\phi }_{\alpha }\frac{{\partial }^2}{\partial {\theta }_j^2}{\phi }_{\alpha }{\mathrm{d}}^{3}x}\right|}_{{\theta }_j={\theta }_{j_0}}^2\end{array}

\begin{array}{}\text{(25)}& \begin{array}{cc}{I}_{Vj}=\hfill & \exp (-Q_j(2{\overline{n}}_j+1)){\left({\displaystyle \frac{{\overline{n}}_j+1}{{\overline{n}}_j}}\right)}^{p_j/2}\hfill \\ \hfill & \cdot J_{p_j}\left(2Q_j\sqrt{{\overline{n}}_j({\overline{n}}_j+1)}\right)\hfill \end{array}\end{array}

with

\begin{array}{}\text{(26)}& {\overline{n}}_j={({\mathrm{e}}^{\beta \mathrm{\hslash }{\omega }_j}-1)}^{-1}\end{array}

\begin{array}{}\text{(27)}& Q_j=I_j\omega {(\delta {\theta }_j)}^2/(2\mathrm{\hslash })\end{array}

\begin{array}{}\text{(28)}& p_j={\displaystyle \frac{\delta E_j}{\mathrm{\hslash }{\omega }_j}}\end{array}

Here, \delta {\theta }_j is the angular displacement, and \delta E_j is the energy gap between the initial and final states for the jth mode. p_j represents the net change in quantum number for oscillator mode j, which satisfies the constraint

\begin{array}{}\text{(29)}& {\displaystyle \sum _j}{p}_j=p\end{array}

in the summation of Eq. (23), and p satisfies

\begin{array}{}\text{(30)}& p={\displaystyle \frac{\mathrm{\Delta }E}{\mathrm{\hslash }\overline{\omega }}}(\mathrm{\Delta }E={\displaystyle \sum _j}\delta E_j)\end{array}

\overline{\omega } is the average of {\omega }_j, namely, {\sum }_j\mathrm{\hslash }{\omega }_j{p}_j=\mathrm{\hslash }\overline{\omega }p. Note that Eq. (25) takes the same form as Eq. (16). To simplify {\sum }_{\{p_j\}}{\prod }_j{I}_{Vj} in Eq. (23), we consider two-mode (N=2) case at first,

\begin{array}{}\text{(31)}& \begin{array}{c}{\displaystyle \sum _{p_1+p_2=p}}{I}_{V1}{I}_{V2}={\left({\displaystyle \frac{{\overline{n}}_1+1}{{\overline{n}}_1}}\right)}^{p_1/2}{\left({\displaystyle \frac{{\overline{n}}_2+1}{{\overline{n}}_2}}\right)}^{p_2/2}\hfill \\ \cdot {\displaystyle \sum _{p_1}}\exp \{-Q_1(2{\overline{n}}_1+1)\}{J}_{p_1}(2Q_1\sqrt{{\overline{n}}_1({\overline{n}}_1+1)})\hfill \\ \cdot \exp (-Q_2(2{\overline{n}}_2+1))J_{p-p_1}(2Q_2\sqrt{{\overline{n}}_2({\overline{n}}_2+1)})\hfill \end{array}\end{array}

By use of the asymptotic formula for Bessel function, Eq. (20), we deduce

\begin{array}{}\text{(32)}& \begin{array}{c}{\displaystyle \sum _{p_1}}\exp (-z_1)J_{p_1}(z_1)\exp (-z_2)J_{p-p_1}(z_2)\hfill \\ ={\displaystyle \frac{1}{2\text{⁢π}\sqrt{z_1{z}_2}}}{\displaystyle {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}}d{p}_1\exp \left(-\left({\displaystyle \frac{p_1^2}{2z_1}}+{\displaystyle \frac{{(p-p_1)}^2}{2z_2}}\right)\right)\hfill \\ ={\displaystyle \frac{1}{\sqrt{2\text{π}}}}{\displaystyle \frac{1}{\sqrt{z_1+z_2}}}\exp \left(-{\displaystyle \frac{p^2}{2(z_1+z_2)}}\right)\hfill \end{array}\end{array}

Under {\overline{n}}_j>1, by inserting (32) into (31), we have

\begin{array}{}\text{(33)}& {\displaystyle \sum _{p_1+p_2=p}}{I}_{V_1}{I}_{V_2}={\displaystyle \frac{1}{\sqrt{2\text{π}}}}{\displaystyle \frac{1}{\sqrt{z_1+z_2}}}\cdot \exp \left(-{\displaystyle \frac{p^2}{2(z_1+z_2)}}\right)\cdot \exp \left({\displaystyle \frac{\mathrm{\Delta }E}{2k_{\mathrm{B}}T}}\right)\end{array}

where

\begin{array}{}\text{(34)}& Z_j={(\delta {\theta }_j)}^2\cdot {\displaystyle \frac{k_{\mathrm{B}}T}{{\mathrm{\hslash }}^2}}{I}_j\end{array}

Eq. (33) takes the same form as Eq. (21) for the single-mode.

Next, by using the similar deduction, we obtain the conclusion that if

\begin{array}{}\text{(35)}& \begin{array}{c}{\displaystyle \sum _{p_1+\mathrm{\cdots }+p_{n-1}=q}}\exp (-z_1)J_{p_1}(z_1)\exp (-z_2)J_{p_2}(z_2)\mathrm{\cdots }\cdot \exp (-z_{n-1})J_{p_{n-1}}(z_{n-1})\hfill \\ ={\displaystyle \frac{1}{\sqrt{2\text{π}}}}{\displaystyle \frac{1}{\sqrt{z_1+\mathrm{\cdots }+z_{n-1}}}}\exp \left\{-{\displaystyle \frac{q^2}{2(z_1+\mathrm{\cdots }+z_{n-1})}}\right\}\hfill \end{array}\end{array}

holds, then it leads to

\begin{array}{}\text{(36)}& \begin{array}{c}{\displaystyle \sum _{p_1+\mathrm{\cdots }+p_n=p}}\exp (-z_1)J_{p_1}(z_1)\exp (-z_2)J_{p_2}(z_2)\mathrm{\cdots }\cdot \exp (-z_n)J_{p_n}(z_n)\hfill \\ =\frac{1}{\sqrt{2\text{π}}}{\displaystyle \frac{1}{\sqrt{z_1+\mathrm{\cdots }+z_n}}}\exp \left(-{\displaystyle \frac{p^2}{2(z_1+\mathrm{\cdots }+z_n)}}\right)\hfill \end{array}\end{array}

Therefore, for multimode transition, we obtain

\begin{array}{}\text{(37)}& {\displaystyle \sum _{\{p_j\}}}{\displaystyle \prod _j}{I}_{Vj}={\displaystyle \frac{1}{\sqrt{2\text{π}}}}\exp \left({\displaystyle \frac{\mathrm{\Delta }E}{2k_{\mathrm{B}}T}}\right){\left({\displaystyle \sum _j}{Z}_j\right)}^{-\frac{1}{2}}\cdot \exp \left(-{\displaystyle \frac{p^2}{2\sum _j{Z}_j}}\right)\end{array}

The electronic matrix elements Eq. (14) or Eq. (24) can be estimated as follows:

\begin{array}{}\text{(38)}& I_E={\left({\displaystyle \sum _j}{\displaystyle \frac{{\mathrm{\hslash }}^2}{2I_j}}⟨l_j^2⟩\right)}^2={\displaystyle \frac{{\mathrm{\hslash }}^4}{4}}{\left({\displaystyle \sum _j}{\displaystyle \frac{a_j}{I_j}}\right)}^2,a_j=⟨l_j^2⟩>0\end{array}

l_j is the jth magnetic quantum number (with respect to {\theta }_j) of electronic wave function {\phi }_{\alpha }(x,\theta ), a_j=⟨l_j^2⟩\approx O(1); here, O(1) means a number in the order of magnitude of 1.

Inserting Eqs. (37) and (38) into Eq. (23) and using Eqs. (30) and (34), we obtain the final results of multimode transition [9-11]:

\begin{array}{}\text{(39)}& \begin{array}{cc}W=\hfill & {\displaystyle \frac{{\mathrm{\hslash }}^3\sqrt{\text{π}}}{2\sqrt{2}\delta \theta \overline{{\omega }^{\prime }}}}{(k_{\mathrm{B}}T)}^{-1/2}\exp \left({\displaystyle \frac{\mathrm{\Delta }E}{2k_{\mathrm{B}}T}}\right){\left({\displaystyle \sum _j}{I}_j\right)}^{-1/2}\hfill \\ \hfill & \cdot {\left({\displaystyle \sum _j}{\displaystyle \frac{a_j}{I_j}}\right)}^2\exp \left({\displaystyle \frac{-{(\mathrm{\Delta }E)}^2}{2{\overline{\omega }}^2{(\delta \theta )}^2{k}_{\mathrm{B}}T\sum _j{I}_j}}\right)\hfill \\ \delta \theta =\hfill & \sqrt{{⟨{(\delta {\theta }_j)}^2⟩}_{\mathrm{av}}}\hfill \end{array}\end{array}

5 Protein folding rate in frequency varying case

In above deduction, the same frequency for initial and final states ({\omega }_j={\omega }_j^{\prime }) has been assumed in calculating overlap integral and obtaining Eqs. (16) and (25). Now, we shall generalize the results to the case of nonequal frequencies between initial and final states. From statistical mechanics theory, the free energy of a system of oscillators is expressed by

\begin{array}{}\text{(40)}& \begin{array}{c}{G}_{\mathrm{os}}={\displaystyle \frac{1}{\beta }}{\displaystyle \sum _j}\left[\ln (1-\exp (-\beta {\omega }_j))+E_j\right]\hfill \\ \beta ={\displaystyle \frac{1}{k_{\mathrm{B}}T}}\hfill \end{array}\end{array}

where E_j means the potential minimum of the jth oscillator. By using

\begin{array}{}\text{(41)}& {\displaystyle \frac{\partial G_{\mathrm{os}}}{\partial (\mathrm{\hslash }{\omega }_j)}}=⟨n_j⟩={\displaystyle \frac{1}{\exp (\beta \mathrm{\hslash }{\omega }_j)-1}}\end{array}

(n_j: phonon number of the jth mode), we obtain the free energy variation

\begin{array}{}\text{(42)}& \begin{array}{cc}{\delta }_{\omega }{G}_{\mathrm{os}}\hfill & ={\displaystyle \sum _j}{\displaystyle {\int }_{{\omega }_j}^{{\omega }_j^{\prime }}}{\displaystyle \frac{\partial G_{\mathrm{os}}}{\partial {\omega }_j}}d{\omega }_j\hfill \\ \hfill & ={\displaystyle \sum _j}{\displaystyle {\int }_{{\omega }_j}^{{\omega }_j^{\prime }}}{\displaystyle \frac{\mathrm{\hslash }}{\exp (\beta \mathrm{\hslash }{\omega }_j)-1}}d{\omega }_j\hfill \end{array}\end{array}

due to the frequency shift from {\omega }_j to {\omega }_j^{\prime }. It leads to

\begin{array}{}\text{(43)}& {\delta }_{\omega }{G}_{\mathrm{os}}={\displaystyle \sum _j}{\displaystyle \frac{1}{\beta }}\ln {\displaystyle \frac{{\omega }_j^{\prime }}{{\omega }_j}}\end{array}

as \beta \mathrm{\hslash }{\omega }_j\ll 1. Therefore, from Eq. (40), the free energy difference between torsion initial state (frequency \{{\omega }_j\}) and final state (frequency \{{\omega }_j^{\prime }\}) is

\begin{array}{}\text{(44)}& \mathrm{\Delta }{G}_{\mathrm{os}}=-{\delta }_{\omega }{G}_{\mathrm{os}}+{\displaystyle \sum _j}\delta E_j=\mathrm{\Delta }E+{\displaystyle \sum _j}{\displaystyle \frac{1}{\beta }}\ln {\displaystyle \frac{{\omega }_j}{{\omega }_j^{\prime }}}\end{array}

Correspondingly, the net change of phonon number p in Eqs. (21) and (37) will be replaced by

\begin{array}{}\text{(45)}& p={\displaystyle \sum _j}{\displaystyle \frac{\delta E_j}{\mathrm{\hslash }{\omega }_j}}+{\displaystyle \sum _j}{\displaystyle \frac{1}{\beta }}{\displaystyle \frac{\ln \frac{{\omega }_j}{{\omega }_j^{\prime }}}{\mathrm{\hslash }{\omega }_j}}\end{array}

Considering that the contribution of frequency variation to folding rate comes mainly from the torsion vibration term, instead of Eq. (39), we obtain a generalized equation of folding rate [12]

\begin{array}{}\text{(46)}& \begin{array}{cc}W=\hfill & {\displaystyle \frac{{\mathrm{\hslash }}^3\sqrt{\text{π}}}{2\sqrt{2}\delta \theta {\overline{\omega }}^{\prime }}}{(k_{\mathrm{B}}T)}^{-1/2}\exp \left({\displaystyle \frac{\mathrm{\Delta }{G}_{\mathrm{os}}}{2k_{\mathrm{B}}T}}\right){\left({\displaystyle \sum _j}{I}_j\right)}^{-1/2}\hfill \\ \hfill & \cdot {\left({\displaystyle \sum _j}{\displaystyle \frac{a_j}{I_j}}\right)}^2\exp \left({\displaystyle \frac{-{(\mathrm{\Delta }{G}_{\mathrm{os}})}^2}{2{\overline{\omega }}^2{(\delta \theta )}^2{k}_{\mathrm{B}}T\sum _j{I}_j}}\right)\hfill \end{array}\end{array}

in the frequency variation case, where \mathrm{\Delta }{G}_{\mathrm{os}} is given by Eq. (44). The application of free energy gap \mathrm{\Delta }{G}_{\mathrm{os}} instead of energy gap \mathrm{\Delta }E was first proposed in the photosynthetic electron transfer when the frequency variation is small [13]. From Eq. (44), we know that when

\begin{array}{}\text{(47)}& \left|{\displaystyle \sum _j}\ln {\displaystyle \frac{{\omega }_j}{{\omega }_j^{\prime }}}\right|\ll {\displaystyle \frac{|\mathrm{\Delta }E|}{k_{\mathrm{B}}T}}\end{array}

the effect of frequency variation can be neglected and Eq. (46) returns to Eq. (39).

6 Application of quantum folding theory to interpretation of experiments

The following problems about the application of the proposed rate theory have been discussed in Ref. [12], and here, we give a brief summary. All these studies support the view that the protein folding is essentially a quantum transition between conformational states.

1) The torsion inertial moment dependences of protein folding rate are calculated. The inertial moments are dynamical variable describing torsion motion of main-chain and side-chain dihedral angles. In principle, they can be estimated by the number of dihedral angles and the mass of related atomic groups on polypeptide chain for each contact. A model of dynamical contact order was proposed for calculating inertial moment dependences of folding rate, and consistent results were obtained by comparison with 100-protein experimental dataset [14]. The correlation coefficient of inertial moment term in lnW with experimental folding rate ln k_f is 0.85, higher than or comparable with other calculations and predictions for protein folding rate based on the same concept of contact order [15, 16].

2) The unified theory on the two-state and multistate protein folding is worked out based on the concept of quantum transition. Under the assumption of the successive transitions of multistate proteins taking place in different torsion modes, the multistate protein folding rate is simply expressed as the rate of two-state protein folding of same topology multiplied by a time delay factor \exp (-\tau ). From a set of 80 proteins (45 two-state and 35 non-two-state), we compared the results of theoretical calculations of folding rate with experimental data and the distribution of the deviation between theoretical and experimental folding rates is given in Fig. 1 [14]. In the calculation of multistate protein folding the time delay factor has been introduced. The peak in Fig. 1 indicates the theory basically consistent with experiments and demonstrates the correctness of the above idea on multistate protein folding [14, 17].

The distribution is plotted for 80 protein dataset [14]. In addition to this dataset, by using 27 two-state folding proteins from the “standard” set given by Maxwell et al. [18], the result can be further improved. (The peak of protein numbers at DCO-ln k_f=0 increases from 63 to 90 [14].)

3) The relationship of folding and unfolding rates versus denaturant concentration is studied. One may assume \delta E_j, and therefore \mathrm{\Delta }E, is a linear function of denaturant concentration. This interprets the observed slope of the arm of the chevron plot, a plot of \ln W versus denaturant concentration [18]. Another factor is the change of the inertial moment I_j due to the varying viscosity of the solvent. One may assume the effective inertial moment I_j is approximately a linear function of denaturant concentration. This can interpret the observed deceleration of folding as the viscosity of the solvent is increased [19].

4) The exergonic and endergonic foldings are distinguished through the comparison between theoretical and experimental rates for each protein. Exergonic means the torsion vibration free energy decreasing in the reaction, \mathrm{\Delta }{G}_{\mathrm{os}}>0, while endergonic means the torsion vibration free energy increasing in the reaction, . We have predicted exergonic/endergonic folding of a protein, and these predictions are basically consistent with experimental data [14].

5) The ultrafast folding problem [20-22] is viewed from the point of quantum folding theory, and a new folding speed limit is deduced from uncertainty relation. From the point of quantum transition, we demonstrate that the necessary conditions for ultrafast folding are given as follows: first, the folding is exergonic rather than endergonic; second, the large initial/final energy gap \mathrm{\Delta }E in torsion potential; and third, the large initial/final frequency ratio of torsion vibration. On the other hand, considering that the high-speed folding generally occurs in downhill processes [23], we have discussed the implication of downhill folding in quantum theory and proved that {(\delta \theta )}^2 should be small in downhill process (Fig. 2). This is another factor to speed up the folding rate. The problem of protein folding “speed limit” is viewed from quantum uncertainty relation, and we suggest a new physical speed limit of picoseconds as the limitation for the local and transient fold formation on the polypeptide chain.

6) The temperature dependence of folding rate is deduced, and the non-Arrhenius behaviors of temperature dependence are interpreted in a natural way. From Eq. (46), one obtains

\begin{array}{}\text{(48)}& \ln W(T)={\displaystyle \frac{\mathrm{\Delta }E(1-\frac{\mathrm{\Delta }E}{\epsilon })}{2k_{\mathrm{B}}T}}-{\displaystyle \frac{1}{2}}\ln T-{\displaystyle \frac{{\lambda }^2}{2\epsilon }}{k}_{\mathrm{B}}T+const\end{array}

\begin{array}{}\text{(49)}& \lambda ={\displaystyle \sum _j}\ln {\displaystyle \frac{{\omega }_j}{{\omega }_j^{\prime }}},\epsilon ={\overline{\omega }}^2{(\delta \theta )}^2{\displaystyle \sum I_j}\end{array}

Therefore, the relation of \ln W versus 1/T is not linear in the Arrhenius plot. The slope–temperature relation is given by

\begin{array}{}\text{(50)}& {\displaystyle \frac{\mathrm{d}\ln W}{\mathrm{d}(\frac{1}{T})}}=S+{\displaystyle \frac{1}{2}}T+R{T}^2\end{array}

\begin{array}{}\text{(51)}& \begin{array}{c}S={\displaystyle \frac{\mathrm{\Delta }E(1-\frac{\mathrm{\Delta }E}{\epsilon })}{2k_{\mathrm{B}}}}\hfill \\ R={\displaystyle \frac{k_{\mathrm{B}}{\lambda }^2}{2\epsilon }}\hfill \end{array}\end{array}

Note that if the structural phase transition occurs for some proteins then the potential parameter may change with temperature and an additional term will introduced in slope S and R, but Eq. (50) remains unchanged.

The experiments on rate–temperature relationships in protein folding exhibit the following characteristics of non-Arrhenius behavior, and they all can be interpreted by the theory above: i) All nonlinearity in Arrhenius plot of \ln W versus 1/T shows the decrease in their folding rates upon increase in temperature at high temperature (slightly curved downward at high temperature) [24-26]. This can be explained by two temperature-dependent terms, -\frac{1}{2}\ln T-\frac{{\lambda }^2}{2\epsilon }{k}_{\mathrm{B}}T, in Eq. (48). ii) In terms of activation energy barrier, the normal positive barrier means that the slope S of Arrhenius plot is negative. However, the Arrhenius plot of folding rate versus temperature for some proteins shows abnormal negative activation energy barriers. In this theory, the negative activation energy barriers are explained by the energy gap between the initial and final states \mathrm{\Delta }E>0, which makes slope positive in Arrhenius plot. The case occurs more generally in ultrafast folders. iii) For some proteins, the crossover from normal positive barrier to abnormal negative at high-temperature occurs. That is, a maximum is observed in Arrhenius plot. The crossover can be interpreted by () of the protein at lower temperatures, but the symbol of \frac{\mathrm{d}\ln W}{\mathrm{d}(1/T)} changes at high temperature due to the temperature-dependent terms in Eq. (50). Note that the \frac{1}{2}T term alone is not enough for interpreting the strong nonlinearity, and R{T}^2 is a key term for the occurrence of crossover. iv) The plots of \ln W versus \frac{1}{T} are strongly curved for refolding of some proteins but almost linear for their unfolding under denaturant. This can be explained by the denaturant concentration effect. The denaturant possibly strengthens the torsion force field and increases the energy gap \mathrm{\Delta }E between initial and final states in the unfolding experiments. Therefore, |S| increases with denaturant concentration, and the relative contribution of non-Arrhenius term RT to the rate \ln W decreases [see Eq. (48)]. v) Recent experiments in λ-repressor fragments folding exhibited another characteristic of rate–temperature relationships. The set of mutants containing more glycines versus alanines exhibit stronger non-Arrhenius behavior [24]. This may be due to the change of torsion potential in the mutation of amino acid. The residue dependence can also be interpreted by smaller inertia moment of glycine than alanine since the non-Arrhenius term R{T}^2 is greater in glycine-rich mutant [see the dependence of \epsilon on inertia moment, Eq. (49)].

To conclude, the fundamental physics underlying protein folding may be much simpler than the extremely complexity inherent in the protein structure, as proposed by Baker ten years ago [27]. Here, our analysis indicates that a key point on the surprising logical simplicity may be in the numerous types of protein folding obey the same universal kinetics of quantum transition between conformational states. The proposed theoretical formulation in the article can, in principle, interpret all basic experimental facts on protein folding. Moreover, the view of quantum transition among torsion states gives deeper insights into the folding event of polypeptide chain.

References

[1] L. F. Luo, Int. J. Quant. Chem., 1987, 32: 435

[2] D. Shepelyansky, Symposium Anderson Localization in Nonlinear and Many-Body Systems, Dresden, 2009

[3] J. D. Bryngelson, J. N. Onuchic, N. D. Socci, and P. G. Wolynes, Proteins, 1995, 21:167

[4] K. Huang and A. Rhys, Proc. Roy. Soc. (London), 1950, A204: 406

[5] D. Devault, Quart. Rev., Biophysics, 1980, 13: 387

[6] J. Jortner, J. Chem. Phys., 1976, 64: 4860

[7] M. Abramowitz and I. A. Stegun, Handbook of Mathematical Functions, 10th printing with corrections, National Bureau of Standards, Applied Mathematics Series 55, 1972

[8] G. N. Watson, A Treatise on the Theory of Bessel Functions, 2nd Ed., Series: Cambridge Mathematical Library, 1995

[9] L. F. Luo, Int. J. Quant. Chem., 1995, 54: 243

[10] L. F. Luo, Theoretic-Physical Approach to Molecular Biology, Shanghai: Shanghai Scientific & Technical Publishers, 2004: 437, 457

[11] L. F. Luo, arXiv: qbio/0906.2452, 2009

[12] L. F. Luo, arXiv: qbio/1008.0237, 2010

[13] T. Kakitani and H. Kakitani, Biochim. et Biophys. Acta, 1981, 635: 498

[14] Y. Zhang and L. F. Luo, Scientia Sinica Vitae, 2010, 40: 887,

[15] K.W Plaxco, T. Simons, and D. Baker, J. Mol. Biol., 1998, 277(4): 985

[16] D. N. Ivankov and A. V. Finkelstein, Proc. Natl. Acad. Sci. USA, 2004, 101: 8942

[17] K. Kamagata, M. Arai, and K. Kuwajima, J. Mol. Biol., 2004, 339: 951

[18] K. L. Maxwell, D. Wildes, A. Zarrine-Afsar, M. A. De Los Rios, A. G. Brown, C. T. Friel, L. Hedberg, J. C. Horng, D. Bona, E. J. Miller, A. Vallée-Bélisle, E. R. Main, F. Bemporad, L. Qiu, K. Teilum, N. D. Vu, A. M. Edwards, I. Ruczinski, F. M. Poulsen, B. B. Kragelund, S. W. Michnick, F. Chiti, Y. Bai, S. J. Hagen, L. Serrano, M. Oliveberg, D. P. Raleigh, P. Wittung-Stafshede, S. E. Radford, S. E. Jackson, T. R. Sosnick, S. Marqusee, A. R. Davidson, and K. W. Plaxco, Protein Sci., 2005, 14(3): 602

[19] M. Jacob, T. Schindler, J. Balbach, and F. X. Schmid, Proc. Natl. Acad. Sci. USA, 1997, 94(11): 5622

[20] L. Qiu, S. A. Pabit, A. E. Roitberg, and S. J. Hagen, J. Am. Chem. Soc., 2002, 124(44): 12952

[21] Y. Zhu, D. O. V. Alonso, K. Maki, C. Y. Huang, S. J. Lahr, V. Daggett, H. Roder, W. F. DeGrado, and F. Gai, Proc. Natl. Acad. Sci. USA, 2003, 100: 15486

[22] H. Neuweiler, C. M. Johnson, and A. R. Fersht, Proc. Natl. Acad. Sci. USA, 2009, 106(44): 18569

[23] M. M. Garcia-Mira, M. Sadqi, N. Fischer, J. M. Sanchez-Ruiz, and V. Munoz, Science, 2002, 298: 2191

[24] M. L. Scalley and D. Baker, Proc. Natl. Acad. Sci. USA, 1997, 94: 10636

[25] W. Y. Yang and M. Gruebele, Biochemistry, 2004, 43: 13018

[26] K. Ghosh, B. Ozkan, and K. A. Dill, J. Am. Chem. Soc., 2007, 129: 11920

[27] D. Baker, Nature, 2000, 405: 39

Rights & permissions

©Higher Education Press and Springer-Verlag Berlin Heidelberg 2010