<gh_stars>0
% Chapter 02 - Theory
\chapter{Theory}
\label{ch:theory}
In this chapter, the survey on mechanical stretching of proteins, along with the theory is described.

\section{Single-molecule protein stretching}
\label{sec:theory-single}
Stretching of single molecules became available by the the tremendous advance in tools used.
In particular, it became possible to attach and precisely move a group of atoms by means of different techniques.
The most common techniques are utilization of atomic force microscope (AFM), optical, and magnetic tweezers\cite{neuman2008single}.
In all the techniques, one end of the molecule is anchored to the stationary bead or surface, while the other is attached to a force sensor.
As a result, the distance between the stationary bead and the sensor is measured, as a function of the applied force.

In case of the AFM, the role of the sensor is fulfilled by the micron-sized cantilever, where the force is proportional to the displacement of the cantilever.
In case of tweezers, the position of the movable bead is trapped in the magnetic field, or by the radiation of the precise laser ray.
The protein of interest is either connected directly, or through the linker (DNA), which can be easily attached to larger beads (nanoparticles) trapped by the laser ray.
Slight changes of the ray position, or the magnetic field shift the movable bead.

To a very good approximation, the trapping force is harmonic with the typical stiffness constants of the trap in range 0.001-0.1 pN/nm, which is around $10^2-10^4$ times smaller than in AFM\cite{neuman2004optical}.
As follows, the force resolution of the optical tweezers is at least 10 times better than in case of AFM\@.
The precision in force measurement, however, costs the precision in the distance measurement.
As one can show\cite{kumar2010biomolecules}:

\begin{equation}
    \sqrt{\langle\delta x^2\rangle}\sqrt{\langle\delta f^2\rangle}\approx k_{B}T
    \label{eq:afm-precision}
\end{equation}
Therefore, the product of uncertainties is independent on the spring constant and better force resolution result in worse distance resolution in optical tweezers compared to AFM\@.

\subsection{Constant velocity vs constant force}
\label{subsec:theory-single-constant}
The measurement can be consistently performed in one of two ways, differing in control parameter.
The control parameter is either the force, or the extension.
Usually in the experiments, the force is measured, while the end of the polymer is moved with constant velocity $v$, determining its position.
As the pulling speed is typically very small, such system may be considered as quasi-equilibrium ensemble of constant distance.
Alternatively, one can measure the mean distance for a constant force system.
In such case, one must wait sufficiently long to neglect the thermal fluctuation.
In thermodynamic limit (when the waiting times are extended to infinity, and the ensembles are infinite), both approaches are equivalent, however, in reality, stalling in force-extension curves corresponding to stalling of the unzipping of DNA chains were observed\cite{weeks2005pause}.

\section{Proteins under stretching}
\label{sec:theory-single-proteins}
Proteins form their 3D structure as a result of non-covalent interaction of amino acids.
During stretching, these contacts are broken, and as a result the protein spatial structure is destroyed.
The force is therefore acting similarly to the a denaturant.

During such mechanical denaturation, rupture of some bonds is more important than others.
In particular,

\section{The WLC and other models}
\label{sec:theory-wlc}
\graffito{The WLC model is also called Kratky\-Porod model\cite{kratky1949rontgenuntersuchung}.}
To explain the observed force-extension dependency, one usually describes the polymers by the Worm-Like Chain model.
In this model, the chain is treated as a continuous rod of a constant length $L$, and the position of the chain element is described as the vector $\overrightarrow{r}(s)$ with $s\in (0,L)$.
The energy landscape is governed by the bending energy dependent on the change of the tangent vector $\frac{\partial\overrightarrow{r}(s)}{\partial s}$ and the persistence length $P_L$ (eq.\ref{eq:wlc-hamiltonian}):

\begin{equation}
    \boxed{E_{elast}={\frac {1}{2}}k_{B}T\int _{0}^{L_{0}}P\cdot \left({\frac {\partial ^{2}{\vec {r}}(s)}{\partial s^{2}}}\right)^{2}ds}
    \label{eq:wlc-hamiltonian}
\end{equation}

To obtain the force-extension dependency, one analyzes the effective Hamiltonian modified with the elastic energy:

\begin{equation}
    H=H_{elastic}+H_{external}={\frac {1}{2}}k_{B}T\int _{0}^{L_{0}}P\cdot \left({\frac {\partial ^{2}{\vec {r}}(s)}{\partial s^{2}}}\right)^{2}ds - zF
    \label{eq:wlc-effective-hamiltonian}
\end{equation}

\graffito{The Marko-Siggia expression was derived by joining the small and high forces regime.}
with the force acting in the $z$ direction.
There is, however, no simple analytic expression on the $F(d)$ dependence in the WLC model.
Usually, one applies the Marko-Siggia approximation\cite{marko1995statistical}, which is correct for large and high forces regime:

\begin{equation}
    \boxed{F(d)=\frac{k_{B}T}{P_L}\left({\frac {1}{4}}\frac{1}{(1-{\frac {d}{L}})^{2}}-{\frac {1}{4}}+{\frac {d}{L}}\right)}
    \label{eq:marko-siggia}
\end{equation}

For intermediate forces, the Marko-Siggia expression introduces error of around 15\%\cite{petrosyan2017improved}.
Therefore, different expressions were introduced to facilitate the calculations (through analytic expressions), or to obtain more exact results.
The commonly used expressions are grouped in Tab.\ref{tab:theory-wlc}.

\graffito{RHMM stands for Rosa-Hoang-Marenduzzo-Maritan.}
\begin{table}
    \begin{tabular}{c|c|c}
        \textbf{Model} & \textbf{Formula} & \textbf{Reference}\\\hline
        Marko-Siggia & $F(d)=\frac{k_{B}T}{P_L}\left({\frac {1}{4}}\frac{1}{\left(1-{\frac {d}{L}}\right)^{2}}-{\frac {1}{4}}+{\frac {d}{L}}\right)$ &\cite{marko1995statistical}\\
        Petrosyan & $F(d)=\frac{k_{B}T}{P_L}\left({\frac {1}{4}}\frac{1}{\left(1-\frac {x}{L}\right)^{2}}-{\frac {1}{4}}+{\frac {x}{L}}-0.8\left({\frac {x}{L}}\right)^{2.15}\right)$ &\cite{petrosyan2017improved}\\
        Elastic WLC & $F(d)=\frac{k_{B}T}{P_L}\left({\frac {1}{4}}\frac{1}{\left(1-{\frac {x}{L}}+{\frac {F}{k}}\right)^{2}}-{\frac {1}{4}}+{\frac {x}{L}}-{\frac {F}{k}}\right)$ &\cite{wang1997stretching}\\
        FJC & $d(F)=L\left(\coth\frac{Fb}{k_{B}T} - \frac{k_{B}T}{Fb}\right)$ &---\\
        Elastic FJC & $d(F)=L\left(\coth\frac{Fb}{k_{B}T} - \frac{k_{B}T}{Fb}\right)\left(1+F/k\right)$ & \cite{smith1996overstretching}\\
        Heterogeneous FJC & $d(F)=N_{f}L_{f}\left(\coth\frac{Fb_f}{k_{B}T} - \frac{k_{B}T}{Fb_f}\right) + N_{u}L_{u}\left(\coth\frac{Fb_u}{k_{B}T} - \frac{k_{B}T}{Fb_u}\right)$ &\cite{su2009mechanics}\\
        Exact WLC & $F(d)=\frac{k_{B}T}{P_L}\left({\frac {1}{4}}\frac{1}{\left(1-{\frac {d}{L}}\right)^{2}}-{\frac {1}{4}}+{\frac {d}{L}} + \sum_{i=2}^{i\leq 7}\alpha_i (\frac{d}{L})^i\right)$ &\cite{bouchiat1999estimating}\\
        Thick chain model & $F(d)=\frac{k_{B}T}{a(1-d/L)}\tanh\left(\frac{k_1(d/L)^{3/2} + k_2(d/L)^2 + k_3(d/L)^3}{1-d/L}\right)$ &\cite{toan2005inferring}\\
        Freely rotating chain & $d(F)=L-L\left( F_{WLC}\left(\frac{Fb}{k_{B}T}\right)^{-1} + \left( \frac{cFb}{k_{B}T}^\beta \right)\right)^{-\beta} + \frac{Lf}{\lambda}$ &\cite{livadaru2003stretching}\\
        RHMM & $F(d)=\frac{2k_{B}TL_P}{b^2}\left( \sqrt{1+(\frac{b}{2L_P})^2\frac{1}{(1-d/L)^2}} - \sqrt{1+\left(\frac{b}{2L_P}\right)^2} \right) + \frac{k_{B}T}{b}\left( 3\frac{1-yL_P/B}{1+yL_P/B} -\frac{b/2L_P}{\sqrt{1+(\frac{b}{2L_P})^2}}\right)\frac{d}{L}$&\cite{rosa2003elasticity}\\
    \end{tabular}
    \label{tab:theory-wlc}
\end{table}

\graffito{The function $f(x)=\coth(x)-1/x$ is called Langevin function.}
FJC stands for Freely-rotating chain.
In this model, the chain is build of rods of constant length $b$, which can freely move around each other.
In this model, one can analytically solve the $F(d)$ dependency.
This derivation may be found in App.\ref{ch:fjc}.

The elastic WLC/FJC was introduced to overcome the constant length of the chain in the model.
In fact, for some polymers (like DNA) the distance between the monomers is dependent on the force acting.
To take this effect into account, the relative extension $x=\frac{d}{L}$ is decreased by a factor dependent on the force $x\rightarrow x-\frac{F}{k}$.
In the equations in Tab.\ref{tab:theory-wlc}, $d$ denotes the measured extension under the force $F$, $b$ is the bond length, $P_L$ is the persistence length and $k$ is the elasticity of the bond.
Other parameters are described in the source articles.

\section{Analyzing the obtained traces}
\label{sec:theory-analysis}
The normal way of analyzing the stretching traces is to fit the WLC function (usually in Marko-Siggia approximation), to find the contour lengths $L_i$ in each stage of stretching.
The fitting is done for each stage of stretching, but the parameters $P$ (persistence length) and $k$ (bond elasticity) have to stay constant throughout whole the experiment.
Therefore, mathematically, fitting the curve consisting of $N$ stretching stages is equivalent to finding $N+2$ parameters ($N$ contour lengths, persistence length and elasticity), fitting $N$ WLC equations.
Usually, when fitting protein traces, the bond length is invariant in the used range of forces, therefore $k=0$ (in experiment).
As a result, a global fit of $N+1$ parameters must be done.

Such an approach immediately creates three problems:
\begin{enumerate}
    \item There might be arbitrary large number of parameters, and the more parameters are fitted, the less accurate the fit is;
    \item The ranges of the fitting have to be decided, which introduces additional ambiguity influencing the final results;
    \item Global fitting of $N$ related equations is rarely implemented in existing tools.
\end{enumerate}

Solution to the last problem may be either writing own code for the global fitting, or using Python packgage SymFit.
However, even when using the published tools, the other two problems remain unsolved.

In this context, interesting innovation was introduced in the work of Puchner et al.\cite{puchner2008comparing}, where only the persistence length (assuming no elasticity of the bond) is enough to determine the contour lengths in each stretching state.
In the original work, when knowing the persistence length $P$, one can invert the WLC equation to find the relative extension $x(F)=\frac{d(F)}{L}$.
One can show (App.\ref{ch:functions}) that there is a unique solution in the interval $(0,1)$.
In theory, the quotient $d(F)/x(F)$ should be equal to the corresponding contour length $L$.
In practice, the measured $d$ and $F$ are subjected to errors, therefore the value of the quotient is not perfectly matching the contour length.
Nevertheless, one can plot the histogram of the quotients which allows to determine the lengths of the consecutive contour length, as the centers of the picks in the histogram (Fig\ref{fig:theory-analysis-puchner})

%\begin{figure}
%    \includegraphics[width=\textwidth]{theory-analysis-puchner.png}
%    \caption{}
%    \label{fig:theory-analysis-puchner}
%\end{figure}

In this technique one performs only one fit (to the histogram), independently on the number of stages.
Moreover, the number of stretching stages can be also deduced from the histogram, as the number of picks.
This technique cannot be, however, used to analyze polymers with elasticity, as in the eWLC model the force is given as implicit function (on both sides of the equation).
It also cannot be used directly, when the protein of interest is stretched by the linker.

In case, there is a linker, the situation is more complex.
In such case, the system can be treated as a set of a few springs connected in series.
On each spring acts the same force $F$, which is the result of a force balance.
Each spring $i$ has then the length $d_i(F)$, and the total extension is then $d(F) = \sum_{i=1}^{N} d_i$.
In case of protein, both ends of the chain are usually connected to two identical copies of the linker (DNA).
Then the total distance measured at force $F$ is equal:

\begin{equation}
    d_{TOT}(F) = d_{PROT}(F) + d_{DNA}(F)
    \label{eq:theory-analysis-sumd}
\end{equation}

To invert the WLC dependence, one usually assumes the high-force regime.
Then, to a good approximation:

\begin{equation}
    d(F)=L\left(1-{\frac {1}{2}}\left({\frac {k_{B}T}{FP_L}}\right)^{1/2}+{\frac {F}{K}}\right)
    \label{eq:theory-analysis-invert}
\end{equation}

This allows to calculate the influence of the DNA, provided one knows the DNA contour length, persistence length and elasticity.
Commonly, these parameters are obtained by fitting the first stage of stretching, when one assumes that only DNA is getting stretched.
This assumption is reasonable in most cases, but not in the case the protein has some loose parts, which can get stretched with the DNA\@.
This assumption is discussed in greater detail in App.\ref{ch:effective-parameters}.

To sum up, the whole algorithm of dealing with the case, when the protein is stretched using the linker is as follows:
\begin{enumerate}
    \item Divide the trace into stages (somehow). Usually it can be done manually, as one clearly sees the ``jumps'' of the trace.
    \item Fit the elastic WLC model to the first stage, obtaining the effective parameters characterizing the linker. This assumes that the protein is not extended in the fist stage of stretching.
    \item Calculate the influence of the DNA, i.e.\ the stretch distance $d_{DNA}(F)$ as a function of the force. Usually, this is done assuming the high-force regime, which is not valid throughout whole the trace.
    \item Subtract the influence of the DNA obtaining the extension of the protein.
    \item Fit globally the curve, obtaining $N+1$ parameters - protein persistence length and $N$ contour lengths, characterizing subsequent stages of stretching.
    \item As the raw contour lengths are subjected to errors, one is usually interested in the differences between the contour lengths, i.e. contour lengths gains.
    \item Knowing the initial end-to-end distance (from the crystal structure) and the contour length gains, one calculates the values of interest, in particular the total length.
\end{enumerate}

As can be seen, the algorithm of analysing the data is subjected to many assumptions, and requires some initial information, like the end-to-end distance in the crystal structure.
In case of the analyzed protein no crystal structure is available.
Moreover, TrmD is built out of two subdomains connected with a flexible linker of length of around 30 residues.
Therefore in case of TrmD (and therefore the fusion), these assumptions are not fulfilled, and searching for alternative roads of reliable analysis became unexpectedly one of the largest task in this project.

\section{Stretching kinetics}
\label{sec:theory-kinetics}
Knowing the number of states represented during stretching (visible by the number of stages in the trace), one can ask about the stability and the lifetime of each stage.

\begin{equation}
    \boxed{\tau(F) = \frac{1}{\dot{F}p(F)}\int_F^\infty p(f)df}
    \label{eq:theory-kinetics-lifetime}
\end{equation}

where $\dot{F}$ is the force loading rate (see App.\ref{subsec:appendices_dhs} for derivation).
On the other hand, the state lifetime is dependent on the thermodynamical values as:

\begin{equation}
    \boxed{\tau(F) = \tau_0 \left(1-\frac{\nu Fx^\ddagger}{\Delta G^\ddagger}\right)^{1-1/\nu} e^{-\Delta G^\ddagger [1-(1-\nu Fx^\ddagger/\Delta G^\ddagger)^{1/v}]}}
    \label{eq:theory-kinetics-dhs-lifetime}
\end{equation}

\section{Modelling the protein stretching}
\label{sec:theory-modelling}

\subsection{Validity of protein models}
\label{subsec:theory-modelling-validity}

\section{Stretching of complex topology proteins}
\label{sec:theory-knotted}
\graffito{The bond distance in stretched protein is equal to 0.365 nm\cite{dietz2004exploring}.}
Complex topology proteins under stretching behave just like the ordinary proteins, except for the fact, that the chain may get tightened on in the place where the (slip)knot reside.
As a result, the final length of the fully stretched protein may be smaller than the expected length.
The loss of the total length is then the hallmark of the (slip)knot of the protein.

What is the reduction of the chain?
In general it depends on the place, where the knot tightens.
If it tightens in the region of bulky residues, it may ``consume'' more amino acids.
It was shown, that contrary to polymers, in case of proteins the knot jumps to a well defined location delimited by sharp turns\cite{sulkowska2008tightening}.
Such locations are characteristic for proteins and follow from their primary and secondary structure.
The observed chain reduction upon the knot presence observed in experiments were collected in Tab.\ref{tab:theory-knotted}.

\begin{table}
    \begin{tabular}{c|c|c|c}
        \textbf{Knot type} & \textbf{Chain reduction [nm]} & \textbf{Residues} & \textbf{Reference}\\\hline
        \multirow{4}{*}{$3_1$} & 5.7 & 16 &\cite{ziegler2016knotting}\\
                               & 4.7 & 13 &\cite{dzubiella2009sequence}\\
                               & 4.7 & 13 &\cite{he2014mechanically}\\
                               & 7.0 & 19 &\cite{rivera2020mechanical}\\\hline
        \multirow{2}{*}{$4_1$} & 6.9 & 19 &\cite{dzubiella2009sequence}\\
                               & 6.2 & 17 &\cite{bornschlogl2009tightening}\\\hline
        $5_2$ & 14.6 & 40 &\cite{ziegler2016knotting}\\\hline
    \end{tabular}
    \label{tab:theory-knotted}
\end{table}

\graffito{SONO is Shrink On No Overlaps\cite{pieranski1998search}.}
Apart from the experiments, the reduction of the chain was probed on the ``ideal chain''.
In particular, one may ask, what is the reduction of the chain, if the chain diameter was maximal, which the chain still non-overlapping?
This leads to the notion of ``tight open knots''\cite{pieranski2001tight}.
The chain reduction for such ideal cases, normalized by the chain diameter, were measured using variation of SONO algorithm for finding ideal conformations of the knots.
The values of $\Lambda_O$ for measured knots are contained in the Tab.\ref{tab:theory-knotted2}.
\graffito{$\Lambda_O=\frac{\Delta L}{D}$ where $\Delta L$ is the chain reduction and $D$ is chain diameter.}

\begin{table}
    \begin{tabular}{c|c}
        \textbf{Knot type} & $\Lambda_O$ \\\hline
        $3_1$ & 10.1 \\
        $4_1$ & 13.7 \\
        $5_1$ & 17.3 \\
        $5_2$ & 18.4 \\
        $6_3$ & 20.7 \\\hline
    \end{tabular}
    \label{tab:theory-knotted2}
\end{table}

In particular, one can deduce, that the expected reduction of the chain length for the $4_1$ knot should be equal $13.7/10.1=1.36$ of the reduction for the $3_1$ knot.
In fact, as the reduction for the $3_1$ knot is around $4.7$ nm, one would expect the reduction of $6.4$ nm, which is close to the measured value of $6.9$ nm.
For $5_2$ knot one would expect $8.6$ nm of chain reduction, which is much smaller, than the measured value.
However, the measured value turns out to be dependent on the force and one observes the further knot tightening for higher forces\cite{ziegler2016knotting}.

Apart from the analysis of the full contour length of the proteins, the experiments on the untangling and refolding of (slip)knots were performed.
In particular, it was shown, that the natively $5_2$ shallowly knotted UCH-L protein can spontaneously fold back to its native state, however the refolding speed depends on the initial topology.
It is the fastest, when the starting structure is $5_2$ knotted, and the slowest for unknotted protein, with $3_1$ knotted as the intermediate case.
Similarly, the $3_1$ slipknotted pyruvoyl-dependent arginine decarboxylase was mechanically unfolded to unknotted state, which could spontaneously return to its native form.
As the refolding was performed with the force sensor attached, one could trace the intermediate contour length, which showed the existence of two- and three-state refolding pathway\cite{wang2019mechanical}.
Likewise, the small $3_1$ slipknotted protein AFV3-109 was shown to refold spontaneously with two-state manner\cite{he2019direct} and two different pathways were shown to feature tightening of this protein\cite{he2014mechanically}.
Some experiments were also supported by corresponding simulations.
In particular, steered MD simulations with atomistic resolution was used to analyze the short subchains of $4_1$ knotted phytochrome\cite{bornschlogl2009tightening}.
These simulations were used only to determine the chain reduction upon knotting.
Similarly, stretching of the shallowly $3_1$-knotted Human Carbonic Anhydrase (HCA) III was analyzed by steered, explicit-water, all atom molecular dynamics in order to analyze the stability of the protein\cite{dzubiella2013tightening}.
The unfolding pathways and the force-extension dependency was not analyzed in this case.
The folding pathways were analyzed in case of $3_1$-slipknotted AFV3-109 protein, where the protein stretching was modelled in all atom CHARMM force field with TIP3P water.
However, in this case, the general unfolding pathway was generated, but again, the theoretical force-extension dependency was not analyzed.
Unfolding of the same protein was also analyzed in coarse-grained $C\alpha$ model\cite{sulkowska2009jamming}, however, without detailed analysis of the force-extension dependency.

%\section{Kinetics and thermodynamics}
%Bell-Evans-Richie
%Friddle-de-Yoreo
%Dudko-Hummer-Szabo
%Jarzynski
%Crooks