Stimulated binding of polymer chains by narrow tube confinement

TL;DR

This paper investigates how narrow tube confinement induces binding between unfolded protein chains by reducing their entropy, using simulations to identify the critical conditions and length scales for such binding.

Contribution

It introduces a novel understanding of polymer chain binding under severe confinement, highlighting entropy-driven interactions in biological channels.

Findings

Binding propensity increases with confinement

Critical tube diameter for binding identified

Average binding length m* quantified

Abstract

In biology, there are several processes in which unfolded protein chains are transported along narrow-tube channels. Normally, without such a severe configurational constraint, unfolded polypeptides would not bind to each other. However, when chain entropy is much reduced in the narrow channel, we find that polypeptide chains have a propensity to bind, even if there is no great potential energy gain in doing so. We find the average length of binding m* (the number of monomers at the chain ends that form bonds) and the critical tube diameter at which such constrained binding occurs. We carry out Brownian dynamics simulations of tightly confined chains, demonstrating their binding over the characteristic length m*, changing in tubes of different diameter.

Equations

$$\frac{F_1}{k_{B}T}\sim -\ln {\left(\frac{d}{a}\right)}^2+{\left(\frac{R_g}{d}\right)}^2=-2\ln \left(\frac{d}{a}\right)+\frac{{\pi }^{2}N}{3}{\left(\frac{a}{d}\right)}^2,$$

$$\frac{F_{\mathrm{bond}}}{k_{B}T}$$

$$\frac{\Delta F}{k_{B}T}=2\ln \left(\frac{d}{a}\right)-m\left[\frac{{\pi }^2{a}^2(d^2-4ad+2a^2)}{3d^2(d-a{)}^2}+\beta \Delta \right],$$

$$P(m)=1-[1-p_{\mathrm{on}}{]}^m=1-e^{-m{k}_{\mathrm{on}}\Delta t(m)}.$$

$$\tau [m]=\frac{(ma{)}^2}{D_c},\mathrm{with}{D}_c=\frac{{\pi }^2{k}_{B}T}{N\gamma }{\left(\frac{d-a}{a}\right)}^2,$$

$$\Delta t(m)=\tau [m+1]-\tau [m]=\frac{N{a}^4\gamma }{{\pi }^2{k}_{B}T(d-a{)}^2}(2m+1).$$

$$f(m)=\frac{(4m+1)a^4{k}_{\mathrm{on}}}{{\pi }^2{D}_0(d-a{)}^2}\exp \left[-\frac{m(2m+1)a^4{k}_{\mathrm{on}}}{{\pi }^2{D}_0(d-a{)}^2}\right].$$

$$m^{\ast }=\frac{\pi }{2}\mathcal{B}-\frac{1}{4},\mathrm{with}\mathcal{B}=\frac{d-a}{a}\sqrt{\frac{D_0}{a^2{k}_{\mathrm{on}}}}.$$

$$\frac{d^{\ast }}{a}=\sqrt{\frac{{\pi }^2{m}^{\ast }/3}{\mathrm{ProductLog}[{\pi }^2{m}^{\ast }{e}^{-m^{\ast }\beta \Delta }/3]}}\approx e^{m^{\ast }\beta \Delta /2},$$

Full text

Much of the physics of polymer chains is determined by entropic effects; entropic barriers associated with chain configurations are widely recognized to control polymer behaviour in narrow pores Muthukumar and Baumgartner (1989); Park and Sung (1998); DiMarzio and Mandell (1997); Muthukumar (2006). Here we consider a related problem of how confinement in a narrow tube could precipitate binding of polymer chains, even though in free unconstrained conditions they would not have a propensity to bind. The problem is motivated by a set of remarkable biological structures that facilitate movement of nucleic acids or polypeptides through narrow channels: (i) contractile tail bacteriophages inject their nucleic acid genome into host bacterial cells through their hollow tail tube Hu et al. (2015); (ii) bacterial type IV secretion systems transfer nucleic acids between bacteria and proteins into eukaryotic host cells Galán and Waksman (2018); (iii) ribosomes release growing polypeptides through the ribosomal exit tunnel where, for some proteins, the beginning of folding is registered while in confinement Farías-Rico et al. (2018); Marino et al. (2016). Finally, bacterial type III secretion systems export unfolded polypeptides that assemble into cell-surface nanomachines that inject bacterial proteins into eukaryotic cells Galán and Waksman (2018) or, in the case of bacterial flagella, facilitate cell motility Evans et al. (2014). During flagella biogenesis, individual subunits synthesised inside the cell need to be unfolded and threaded through a narrow channel that runs the full length of the growing structure, to be assembled at the distal tip Macnab (2003); Turner et al. (2012); Evans et al. (2014).

For flagella assembly, two mechanisms of subunit transport along the narrow channel have been investigated. One relies on single-file diffusion Stern and Berg (2013), which requires ‘pushing’ a column of subunits through the narrow channel, against tremendous resistance once the crowded regime sets in Lappala et al. (2013). The other mechanism suggests that unfolded flagellar subunits form a linked chain along the whole length of the flagellar channel, which is then ‘pulled’ to the site of subunit assembly at the flagellum tip Evans et al. (2013). Here, we propose a mechanism by which unfolded polymers would strongly bind to form ‘head-to-tail’ links, when constrained in a channel.

The properties of polymer chain ($N$ segments of size $a$) constrained in a narrow channel are very well known, starting from the seminal work of Casassa Casassa (1967) and Edwards Edwards and Freed (1969); Doi and Edwards (1988), and extensively reviewed in many subsequent publications de Gennes (1979); Muthukumar (2006); Sakaue and Rafael (2006). The free energy excess contains two contributions: the ‘ideal gas’ motion of the chain centre of mass along the tube, and the reduction in configurational entropy under constraint:

[TABLE]

where the second expression is for the ideal Gaussian chain, with the radius of gyration $R_{g}=N^{1/2}a$. This is good enough for our purposes; in fact, the commonly used model of self-avoiding chain in good solvent might be less appropriate for a case of protein unfolded in a channel Evans et al. (2013). When two such chains are in the channel, Fig. 1(a), the free energy is additive: $2F_{1}$. When these two chains are instead bonded over the length of $m$ overlapping units, see Fig. 1(b), the corresponding free energy has only one centre-of-mass term. The configurational free energy has $2(N-m)$ monomers constrained in a tube of diameter $d$, and strictly speaking, $m$ monomers constrained in a smaller diameter $(d-a)$:

[TABLE]

where $\Delta$ is the potential energy gain on making one bond (see Miyazawa and Jernigan (1985) and many subsequent studies of this energy for aminoacid residues in contact). We shall be interested in the situation when this gain is small (or non-existent), so that the monomers would not bond in the free-chain conditions. Figure 1(c) shows these two free energies compared, as a function of changing tube diameter. What we discover, is that at large $d$ (when the entropic constraint is weak), the free energy of two independent chains is lower, but as $d$ becomes smaller – the free energy of the bonded pair is lower. In fact, the free energy difference $\Delta F=F_{\mathrm{bond}}-2F_{1}$ no longer depends on the overall chain length $N$, and reveals the key equilibrium effect:

[TABLE]

with $\beta=1/k_{B}T$. We see that the ideal-gas entropy ‘wins’ at large $d$, but at $(d-a)\rightarrow 0$ the confinement effect enforces the binding via the negative-$m$ contribution. Even when there is no potential gain in binding ($\Delta=0$) or the monomers repel each other (negative $\Delta$), still in a sufficiently narrow tube they would prefer to bind.

However, the equilibrium analysis leading to the Eq. (3) cannot predict which is the length of bonded segment: obviously, in equilibrium, the larger the overlap $m$ is, the lower is the free energy. The reason why the two chains adopt a particular binding length $m^{*}$ is entirely kinetic. Two effects compete: the effective rate of the bonding reaction, and the rate of chain reptation along the tube to increase the overlap length $m$. Once any one of the $m$ bonds along the overlapped segment is established, the reptation stops, and the rest of the bonds set rapidly, in sequence, via the process of accelerated zipper Bell and Terentjev (2017).

Let us say that the rate of an individual bonding reaction between two monomers is $k_{\mathrm{on}}$, determined by the chemistry involved. The probability of setting an individual bond during a given time interval $\Delta t$ is Glasstone et al. (1941): $p_{\mathrm{on}}=1-e^{-k_{\mathrm{on}}\Delta t}$ (at small $\Delta t$, $p_{\mathrm{on}}=k_{\mathrm{on}}\Delta t$, while at large $\Delta t$, the reaction occurs with an almost certainty). The probability to establish at least one bond out of a possible $m$ sequence is

[TABLE]

Any one of the $m$ possible bonds along the overlapping segment stops the reptation diffusion, and is rapidly followed by the full bonding of the segment.

We are now left with the second aspect of this kinetic problem: to find the time interval $\Delta t(m)$ for which the sequence of $m$ monomers dwells when two chains slide reptate past each other. The theory of polymer reptation is classical Doi and Edwards (1988), and in our context we need to recall the characteristic time to diffuse the chain of $N$ monomers over a distance $ma$ in a tube of the remaining clearance $(d-a)$, cf. Fig. 1(b):

[TABLE]

where the diffusion constant for the chain centre of mass is given by $D_{c}$, in which $\gamma$ is the friction constant of just one monomer (so $D_{0}=k_{B}T/N\gamma$ is the diffusion constant of the whole polymer in free space: a value frequently measured for globular proteins in water). Equation (5) is more familiar in the form of reptation time of the whole chain length ($Na$) in a clear tube ($d$), when $\tau_{d}=N^{3}a^{4}\gamma/\pi^{2}k_{B}Td^{2}$; see Doi and Edwards (1988) for detail. We now have:

[TABLE]

All these expressions from reptation theory are only valid for the tightly confined chain, with $d\ll R_{g}=N^{1/2}a$.

Figure 2 illustrates how the probability $P(m)$ depends on its key parameters. Now the kinetic problem we are addressing becomes that of the mean first-passage time Weiss (1967); Talkner (1987). The question is: at which $m^{*}$ the first binding reaction would occur. The probability density that binding occurs at a given value of overlap $m$ is given by $f(m)=dP(m)/dm$, so that the mean first-binding length is $\bar{m}=\int mf(m)\mathrm{d}m$. Carrying out the algebra, we obtain $f(m)$ illustrated in Fig. 2 inset:

[TABLE]

The peak of this probability distribution is determined by a single non-dimensional parameter, which we may call ‘bonding enhancement’ ${\cal B}$:

[TABLE]

The mean bonded length $\bar{m}$ is also easily calculated, but its expression is more cumbersome. Figure 3 shows that there is no great difference in values between the mean and the median of $f(m)$, and for all practical purposes we may work with the simpler expression in Eq. (8).

It is interesting to assess the critical diameter of the tube, $d^{*}$, below which the bonding is likely. It is obtained by solving the transcendent equation $\Delta F(d)=0$, and has the form

[TABLE]

where $m^{*}$ is given by Eq. (8), and the approximate form is valid when $m^{*}e^{-m^{*}\beta\Delta}\ll 1$ (i.e. essentially for all non-negative $\Delta$). The exponential approximation of Eq. (9) is prominent for all non-repulsive interactions, while for negative $\Delta$ (repulsion between monomers) the other limit of ‘ProductLog’ is in force, and $d^{*}\sim a$ (in other words, bonding is unlikely). However, there is another upper boundary for the channel diameter that promotes binding: all analysis in Eqs. (6-7) was based on the chain reptation dynamics, which stops being applicable at $d\geq N^{1/2}a$ when the random-coil (blob) chain dynamics would take over (the exact crossover is hard to identify). This is schematically labelled in Fig. 3.

Assuming we are in narrow enough channels for the reptation theory to work, we should try estimating the kinetically-set bonding length $m^{*}$, which requires the value of parameter ${\cal B}$. The friction constant for one monomer (e.g. aminoacid) can be obtained approximately from $\gamma=6\pi\eta a$, with $\eta$ being the viscosity of water (0.7 mPa.s), and $a=0.3$ nm the size of a residue Miyazawa and Jernigan (1985). A flagellin protein with $N\approx 495$, at room temperature, should have the diffusion constant: $D_{0}=2\cdot 10^{-12}\mathrm{m}^{2}/$s (in agreement with a typical diffusion constant of aminoacids in water Longsworth (1953); Ma et al. (2005): $k_{B}T/\gamma\sim 10^{-9}\mathrm{m}^{2}/$s, divided by 500).

It is hard to find values of the bonding reaction rate $k_{\mathrm{on}}$ for a pair of aminoacid residues: it strongly depends on a large number of specific chemical factors some good reference (1985). In the spirit of our ‘average polymer chain’ with a single characteristic value of bond energy gain $\Delta$, we may try a basic Kramers estimate: $k_{\mathrm{on}}=\omega_{0}e^{\beta\Delta}$, where $\omega_{0}$ is the collision frequency in solution: $\omega_{0}\approx 10^{8}\,\mathrm{s}^{-1}$ Glasstone et al. (1941); Laidler (1963). Then, for a low $\beta\Delta=1$, we have $k_{\mathrm{on}}=2.7\cdot 10^{8}\,\mathrm{s}^{-1}$. Using Eq. (8) with the diameter of flagellar channel $d=2$ nm, we obtain a value: ${\cal B}=1.63$, and the preferred length of binding $m^{*}=2.3$. As a consequence, for the same set of parameters, the free energy of binding is very low: $\Delta F\approx-k_{B}T$. A high reaction rate would commit the constrained chains to bond on the first monomer ($m^{*}\leq 1$), and the bonding free energy is determined purely by the potential energy of binding $\beta\Delta$. For comparison, a lower reaction rate $k_{\mathrm{on}}=10^{7}\,\mathrm{s}^{-1}$ would give: ${\cal B}=8.4$ and $m^{*}=13$, with the coresponding binding free energy $\Delta F\approx-12k_{B}T$ (maintaining $\beta\Delta=1$).

In order to test the ideas of the confinement-induced bonding, we carried out simulation of two chains in a reflective channel, using the LAMMPS package. We used the Kremer-Grest model for coarse-grained polymers Kremer and Grest (1990); Rosa and Everaers (2008) where bonds between consecutive monomers are modelled with a harmonic-like finitely extensible elastic potential, made of a combination of the attractive FENE part and the repulsive LJ part (the approach used in many studies, old and new Jung et al. (2012); Jun et al. (2007)), and the weak bending potential with the constant $k_{\theta}=k_{B}T$. The non-bonded particles interact through the LJ potential of variable weak attraction strength (see Lappala et al. (2013, 2019)). Two kinds of LJ potential create a subtlety in the simulation setup: one needs to be clear which ‘LJ units’ are used. We kept the average temperature constant, equal to the attractive strength $\varepsilon$ of the pair LJ potential (which is equivalent to our earlier condition $\beta\Delta=1$). The LJ time and the damping constant of the simulation were linked to the same energy scale via the fluctuation-dissipation theorem: $\tau=\sigma\sqrt{M/\varepsilon}$ and damp$=1/\tau$, where $M$ is the mass and $\sigma$ the diameter of an individual monomer, equivalent to the parameter $a$ used in the analysis above (the simulation time-step was $dt=0.01\tau$). On the other hand, in the Kremer and Grest model, the bond potential has the FENE bond strength linked to the repulsive-LJ strength: $k=30\varepsilon^{*}/\sigma^{2}$ Kremer and Grest (1990); Rosa and Everaers (2008), which makes the bond length $r_{\mathrm{bond}}\sim 0.96\sigma$.

However, we found that when the LJ repulsion strength $\varepsilon^{*}=\varepsilon$ (which, in turn, is equal to $k_{B}T)$, the chains do not bind in the channel. The reason was that the length of the bond between monomers on the chain was fluctuating strongly, and the required ‘matching’ of monomers along the two chains was not possible. When we insisted that the bond along the chain is strongly confined, by taking $\varepsilon^{*}=100\varepsilon$, the confinement-induced binding of two chains was a strong and obvious effect, see Fig. 4. We observed that tightly confined chains eventually come in contact during their reptation diffusion; mostly such a contact was short-lived and the chains separated again – but on rare occasions, when a sequence of several monomers come into contact simultaneously, the chains ‘lock’ in binding. After a long simulation time, the number of overlapping monomers $m^{*}$ was reaching the equilibrium value, depending on the channel diameter, see Fig. 5. Note that in these simulations the $m^{*}$ plateau value for $d=14a$ and $18a$ was the same; this was a systematic result reproduced in a large number of independent simulation runs, for all large diameters where the binding occurred. For $d=20a$ we did not register any binding event within the given time of simulation.

As the stimulated-binding process is clearly driven by chain kinetics in the confining tube, we have done many simulations (with different random seeds) to establish the statistical distribution of binding time $t_{\mathrm{B}}$. Figure 6 gives an illustration of such a distribution, $f(t_{\mathrm{B}})$, obtained for $d=5a$; the plot also gives the result for $d=10a$ in the background, but we were not able to acquire a sufficient statistics for larger $d$ because of too many instances of very long $t_{\mathrm{B}}$. This distribution has a very characteristic ‘heavy tail’: it decays with a weak power law (our fitting suggests $t^{-0.8}$). One of the consequences is that the average binding time $\langle t_{\mathrm{B}}\rangle=\int t_{\mathrm{B}}f(t_{\mathrm{B}})dt_{\mathrm{B}}$ is not defined (the integral diverges), and we have to use the median $\overline{t_{\mathrm{B}}}$ instead. Such distributions are found in the Levy flight processes Mantegna and Stanley (1994); Metzler and Klafter (2000), and in our case the origin is similar: the events leading to an overlap of $m^{*}$ monomers are rare, but only these lead to the lasting binding of the chains. The key point is clear: two chains with very weak attraction (i.e. not naturally prone to binding on contact in the bulk) would be strongly bound if a large number of their monomers would simultaneously overlap. In a weakly constrained case, the probability of such overlap is very low, and the time $t_{\mathrm{B}}$ to wait for such a rare event is longer than any realistic experimental time.

The theory presented here is deliberately qualitative, demonstrating the key concept of stimulated bonding in narrow channels. Its analysis is essentially scaling, with explicit calculations based on the ideal Gaussian chain limit (however, this limit becomes accurate for polymer chains under strong lateral confinement, such as in dense melts or other ‘tube-model’ situations). The idea of an ‘average’ homopolymer chain is also very limiting: a real protein would have a variety of aminoacid residues with very different potential energy of pair binding $\Delta$, between $-4k_{B}T$ and $+4k_{b}T$ (see Miyazawa and Jernigan (1985, 1996) for detail): a follow-up problem with the broadly distributed quenched $\Delta$ is certainly of interest. Nevertheless, these approximations and simplifications allowed us to expose the physics of entropically stimulated bonding, and obtain analytical expressions for the main points of interest.

It is possible that a similar analysis would predict a much enhanced protein folding, once the conformational entropy is reduced by the tight channel confinement, as in ribosomes Farías-Rico et al. (2018); Marino et al. (2016). Another interesting corollary of this concept is the need to examine effective reaction rates in micro-volumes. Today the technology of microfluidics and a ‘lab on a chip’ is widespread: the reduction of translational entropy (and when unfolded polymer chains are involved, also the configurational entropy) needs to be taken into account when analysing the reaction rates in very small volumes.