Active translocation of a semiflexible polymer assisted by an ATP-based molecular motor

Abstract

In this work we study the assisted translocation of a polymer across a membrane nanopore, inside which a molecular motor exerts a force fuelled by the hydrolysis of ATP molecules. In our model the motor switches to its active state for a fixed amount of time, while it waits for an ATP molecule which triggers the motor, during an exponentially distributed time lapse. The polymer is modelled as a beads-springs chain with both excluded volume and bending contributions, and moves in a stochastic three dimensional environment modelled with a Langevin dynamics at a fixed temperature. The resulting dynamics shows a Michaelis-Menten translocation velocity that depends on the chain flexibility. The scaling behavior of the mean translocation time with the polymer length for different bending values is also investigated.

Introduction

Translocation of long molecules through nanopores in cell membranes is a common process in living systems. Cell drug delivery, DNA, RNA and protein passage through cell membranes and nuclear pores, and DNA injection and packaging by phage viruses are only some interesting examples of a broad biological phenomenology¹.

The passage of polymers thorough nanopores is also a fundamental problem in chemical and industrial processes. In this context, many efforts are made in nanotechnological applications, that try to emulate the complex biological processes involved in the translocation problem^2,3,4. An important related application is the use of nanopores to unzip and translocate single DNA chains with the purpose of performing fast and detailed DNA sequencing^5,6,7,8,9.

In spite of the many studies present in the literature, the understanding of polymer translocation at the nanoscale still deserve deeper investigation. In fact, the number of parameters involved in the models and the presence of nonlinearities together with the underlying stochastic environment in which nanoscale objects move, make this study complex theoretically and hard computationally, even nowadays. DNA chains are paradigmatic examples of this difficulty. With the ambition of efficiently describing the biological matter in a reasonable time, different mesoscopic models for polymer translocation have been introduced¹⁰. In some very simplified models a single barrier potential, eventually depending on time, is introduced to depict the overall translocation process^11,12,13,14; in others, stochastic and ratchet-like forces and potentials are used refs 15,16,17.

Motivated by the passive pores studied in different experiments^18,19,20, most studies of translocation have been performed, so far, mainly under constant forces. These forces, are generally associated to potential and concentration gradients observed between both sides of the membrane. The role of active nanopores, with mechanisms to assist the polymer translocation, has been only more recently considered. The simplest active nanopore corresponds to the opening or closing of the channel due to either the presence of electrochemical membrane forces²¹ or random bonding of ligands²². It has been usually modelled by a dichotomous Markov noise (also called random telegraph noise)^23,24,25 which drives the molecule. In other cases, the flickering activity of the pore channel is modelled as a sinusoidal pore actuation^25,26,27,28. In both cases a new time scale, associated to the pore activity, appears in the description of the phenomena.

A landmark example of active translocation in the biological realm is the translocation assisted by molecular motors. Recent well known experiments have risen a considerable attention to the specific and very interesting driving features of the molecular motor of the virus bacteriophage ϕ-29. This motor is able to inject out (or pack inside) its DNA by exerting a force supplied by the hydrolysis of ATP molecules which bind to the motor sites and activate its movement^{29, 30}. The force exerted by the nanomotor can be modelled by a specific stochastic dichotomous force, that we call dichotomous ATP-based motor noise force. This force has been successfully introduced to study translocation restricted to a one dimensional dynamics^{31, 32}.

In addition to its biological interest, understanding and controlling bistable forces has an undoubted interest in the manipulation of objects at the nanoscale. In fact, molecular motors represent a highly effective way to exert a force at the nanoscale worth to be mimicked in the artificial realm. So, a properly prepared system can be driven in a controlled way by changing the concentration of the fuel molecules in the solution, which directly modifies the motor activation rate, and thus the exerted force.

Concerning the polymer modelling, the usual approach is based in the Rouse chain³³ and its modifications. There, polymers are constituted by beads connected each other by harmonic springs. Natural improvements of the model take into account polymer bending energy and excluded volume interactions. Indeed, recent studies show that translocation significantly depends not only on the polymer size but also on its flexibility and on the interactions between monomers and pores^{5, 24, 34,35,36}.

Aimed by all these results, this manuscript proposes the study of the translocation of a linear polymer molecule in the 3d domain driven by a specific ATP-fuelled biological force acting inside the nanopore. We will specially focus on the study of the mean translocation velocity and translocation time (TT) τ of the polymer as a function of the motor activation rate. We will put special attention to the study of the dynamics of the system for different polymer size (number of monomers) and persistence length (flexibility or bending). This latter is an important parameter of the system whose influence has been rarely studied so far in translocation context.

The Model

System equations

Our polymer model is based in a chain formed by N identical monomers moving in the three-dimensional space. We use a modified Rouse model³³ which includes bending energy, excluded volume effects and interaction with the membrane and the pore. The elastic potential energy is given by

$${V}_{{\rm{el}}}({d}_{i})=\frac{{k}_{e}}{2}\sum _{i=1}^{N}\,{({d}_{i}-{d}_{0})}^{2},$$

(1)

where k _e is the elastic parameter, r _i is the position of the i-th particle, ${d}_{i}=|{{\boldsymbol{d}}}_{i}|=|{{\boldsymbol{r}}}_{i+1}-{{\boldsymbol{r}}}_{i}|$, is the distance between the monomers i and i + 1, and d ₀ is the equilibrium distance between adjacent monomers.

The model takes into account the bending energy of the chain with a term given by

$${V}_{{\rm{ben}}}({\theta }_{i})=\frac{{k}_{b}}{2}\sum _{i=1}^{N}\,[1-\,\cos ({\theta }_{i}-{\theta }_{0})],$$

(2)

where k _b is the bending elastic constant, θ _i is the angle between the links d _i+1 and d _i , and θ ₀ the equilibrium angle, with θ ₀ = 0 in our case. With this term, our model is a discrete version of the worm-like chain (WLC) model.

In order to consider excluded volume effects between the monomers, a repulsive only Lennard-Jones potential has been taken into account:

$${V}_{{\rm{LJ}}}({r}_{ij})=4\epsilon \,\sum _{i\ne j=1}^{N}\,[{(\frac{\sigma }{{r}_{ij}})}^{12}-{(\frac{\sigma }{{r}_{ij}})}^{6}]$$

(3)

for ${r}_{ij}\le {2}^{\mathrm{1/6}}\sigma$, and −$\epsilon$ otherwise. Here r _ij is the distance between monomer i and monomer j.

The dynamics of every monomer of the chain is obtained by the overdamped equation of motion

$$\begin{array}{rcl}m\gamma {\dot{{\boldsymbol{r}}}}_{i} & = & -{\nabla }_{i}{V}_{{\rm{el}}}({d}_{i})-{\nabla }_{i}{V}_{{\rm{ben}}}({\theta }_{i})-{\nabla }_{i}{V}_{{\rm{LJ}}}({r}_{ij})\\ & & +{F}_{drv,i}{\boldsymbol{i}}+{{\boldsymbol{F}}}_{sp,i}+\sqrt{2m\gamma {k}_{B}T}\,{{\boldsymbol{\xi }}}_{i}(t),\end{array}$$

(4)

where the effective viscosity parameter of each monomer is included in the normalised time units. ξ _i (t) stands for the Gaussian uncorrelated thermal fluctuation and follows the usual statistical properties $\langle {\xi }_{i,\alpha }(t)\rangle =0$ and $\langle {\xi }_{i,\alpha }(t){\xi }_{j,\beta }(t^{\prime} )\rangle ={\delta }_{ij}{\delta }_{\alpha ,\beta }\delta (t^{\prime} -t)$, with $i=1,\ldots ,N$, and α = x, y, z. The operator used above is defined as: ${\nabla }_{i}=\partial /\partial {x}_{i}\,{\boldsymbol{i}}+\partial /\partial {y}_{i}\,{\boldsymbol{j}}+\partial /\partial {z}_{i}\,{\boldsymbol{k}}$.

The force term F _sp includes both the chain-membrane and chain-pore spatial constraint. This interaction force is modelled with the same repulsive Lennard-Jones potential described in Eq. (3), also with the same parameter values. It takes place uniformly and perpendicularly to all the planes that define both the membrane and the pore channel, modelled as a square prism of base L _M and length L _H (see Fig. 1).

Figure 1

Section of the polymer translocating through a nanopore in 3d. The pore has a square section of width L _H and its length is L _M . Inside the pore it is present a driving force in the x-direction, which pushes the polymer toward the trans side of the membrane. The walls of the membrane repulse uniformly the chain inside a characteristic distance σ.

Finally, F _drv,i is the driving force which represent the force F _drv acting on particle i when it is inside the pore. F _drv is constituted by two terms: a constant force F _c and the motor force η(t). F _c is the constant force used in the study of passive pores (usually associated to potential and concentration gradients between both sides of the membrane). η(t) is associated to the molecular motor pulling the polymer. It fluctuates between 0 and F _M . For η(t) we use a modification of the usual dichotomous noise³¹. Activated by ATP absorption, η acts during a fixed time T _M , along this time η(t) = F _M . After this time, the motor relax and η(t) = 0 until a new ATP molecule is absorbed. This happens with a probability P _t′ = 1 − e^−t′/T, being t′ the time spent by the motor in its inactive state, and T its mean waiting time (see Fig. 2). We define ν = 1/T as the related activation rate which is assumed to be proportional to the ATP molecule concentration.

Figure 2

Scheme of the dichotomous pushing force acting on the polymer. T _M is the working time of the motor, supposed fixed. T is the mean waiting time in the inactive state. The motor acts on a region of length L _M , see Fig. 1.

The motor acts only on the monomers inside the pore. Thus, regarding the spatial dependence of the total driving force:

$${F}_{drv}(x)=\{\begin{array}{ll}({F}_{c},{F}_{c}+{F}_{M}) & x\in [0,{L}_{M}]\\ 0 & otherwise.\end{array}$$

(5)

ATP Energy supply and Michaelis-Menten kinetics

The energy used by the motor to provide the driving is given by the ATP hydrolysis. The simplest way to model this phenomenon, is to consider that each ATP hydrolysis activates the motor for a fixed working time T _M . The ATP binding and subsequent absorption is itself a Poisson process which depends on the ATP concentration in the solution surrounding the motor. Once absorbed ATP is hydrolysed releasing energy, and the motor changes its conformational state and exerts a force F _M (considered constant here) during a fixed time T _M  = 1/ν ₀. After that time the motor returns to its rest state. A new force will be applied when the next ATP suitable quantity is absorbed and hydrolysed, which happens after a mean time T = 1/ν, which follows an exponential distribution of waiting times according to the underlying Poisson process. The timing in the kinetics can be summarised as follows: For activation rates ν < ν ₀, the motor remains at rest for an average time longer than T _M , while the opposite occurs (average rest-state time of the motor smaller than T _M ) if ν > ν ₀.

The kinetic feature of the motor, i.e. that it works for a fixed time, and that the statistics of the arrival of the ATP molecules happens in an exponential distribution, is a good realistic approximation as put in evidence by different experimental works^{29, 30, 37}. With the above definitions, the Michaelis-Menten (MM) law for the transfer velocity arises naturally from the model. In fact, the fraction of active time with respect to the total time is given by

$$\frac{{T}_{M}}{{T}_{M}+T}=\frac{\nu }{{\nu }_{0}+\nu }=\frac{[ATP]}{{k}_{M}+[ATP]},$$

(6)

which is a MM law. The only hypothesis included in the derivation of the last equality is that the motor activation rate ν is proportional to the ATP concentration ($\nu \propto [ATP]$), being k _M  = [ATP]₀, the ATP concentration at ν = ν ₀. In our approach the binding of an ATP molecule to the motor is a Poisson process, with a binding rate that can be considered, in a first approximation, proportional to the ATP concentration. This reasonable assumption is compatible with the experimental observation reported in ref. 30. There, data are adjusted to a more complex functional relation, but it is easy to show that the used fitting law also results in a Michaelis-Menten velocity after a suitable re-scaling of the involved parameters (see also ref. 38).

The relationship between the mechanical active-rest kinetics of the motor and the MM law goes beyond the statistical ansatz $\nu \propto [ATP]$, and represents a general paradigmatic behaviour connecting chemical reactions with mechanical transfer of energy in the context of microscopic description of single molecules dynamics. In fact, it is possible to depict the ATP hydrolysis in the context of the MM enzymatic reaction scheme:

$$E+S\mathop{\to }\limits^{{k}_{1}}Z\mathop{\to }\limits^{{k}_{2}}E+P,$$

(7)

where the rate k ₁ represents the probability to form the compound Z per unit of time and per unit of S ([ATP]), and k ₂ gives the probability to form the product P per unit of time. In our mechanical and individual case (single motor and single ATP hydrolisation event) Z represents the ATP-motor binding, which occurs with the rate ν (i.e. ~k ₁[ATP]), while the product of the reaction P represents the motor action which is completed within a time T _M  = 1/ν ₀ (ν ₀ ~ k ₂). These relationships are in agreement with the definition of the Michaelis constant k _M  = k ₂/k ₁. For a more complete analysis of the MM law the reader can also refer to refs 39 and 40.

Polymer velocity

The polymer velocity (along the x-axis) is calculated by summing up the N monomer terms of equation (4) and averaging over the time of the translocation. The mean velocity of the centre of mass of the chain is then:

$$\begin{array}{rcl}{v}_{CM} & = & \frac{{F}_{drv}}{N}+\frac{{F}_{w}}{N}\\ & = & \frac{{F}_{c}+{F}_{M}}{N}\,\sum _{i=1}^{N}\,\langle {\eta }_{i}(t)\rangle +\frac{{F}_{w}}{N}\\ & = & \frac{{F}_{c}{n}_{{\rm{M}}}}{N}+\frac{{F}_{M}}{N}\frac{{n}_{{\rm{M}}}}{1+{\nu }_{0}/\nu }+\frac{{F}_{w}}{N}.\end{array}$$

(8)

In the above equation n _M is the mean number of monomers inside the motor during the translocation, a number which weakly depends on ν, k _b and N, and we consider here as a constant, n _M  = 5, as confirmed by our numerical observations.

We identify here three contributions to the total velocity. The first term accounts for the driving of the constant external force, which acts only on the monomers into the motor. The second term is related with the molecular motor kinetics, explained above. It is given by the force felt by the n _M monomers inside the machine which operates for the fraction of time ν/(ν ₀ + ν), and shows the Michaelis-Menten dependence. The last term F _w (ν, k _b , N) comes from the interaction of the polymer with the membrane. In fact, while all the internal force terms in equation (4) sum up to zero when we sum the N equations of the monomers, the interaction with the walls in the x-direction does not average to zero, being a non symmetric force reaction, as the pulling goes to positive x-direction. The resulting force opposes to the righthand movement of the chain, and depends on the activation rate of the motor ν, the length N, and the rigidity of the chain k _b . It is important to highlight that all the k _b dependence in the polymer velocity comes from this term. We will also see below that although F _w also depends on ν, the overall v _CM (ν) dependence follows the MM behaviour, dominated by the motor contribution to the total force. Note also that the polymer velocity goes to zero as 1/N for large chains as expected when a motor acts on a small number of monomers n _M of a polymer which moves in a dissipative media.

Numerical Simulations

Definitions and simulation details

We are interested in characterising the translocation process of the polymer through the pore. We will present below the results for τ, the mean value of the TT of the chain for different parameters of the model and different polymer lengths.

It is usual in literature that a channel with either zero or short length is used to study the polymer translocation. In order to approach a more realistic picture, in our works we use a channel with a fixed length, L _M  = 5d ₀, longer than the distance of two consecutive monomers. We calculate the translocation time as follows: All the simulations start with five monomers at the right end of the polymer lying inside the pore and the others linearly ordered and sited at the equilibrium distance. During a thermalisation time t _t  = 1000 t.u. the chain evolves under the action of thermal fluctuations while keeping fixed the position of the five monomers inside the pore. After that transient time, the restriction over the first five monomers is removed, and the evolution of the chain under the dynamics given by Eq. (4) is monitored. Sometimes the polymer moves backwards into the cis region of the system, leaves the pore and it does not translocate. These cases are not taken into account, i.e. they do not enter in the translocation statistics. On the contrary, whenever the polymer crosses the membrane the simulation ends when the last monomer of the chain enters the pore, and the time lasted from the end of the thermalisation process to the entrance of the last monomer into the pore defines the translocation time (TT) of the event. Later we will average over many translocation processes to get the values of τ.

At the same translocation conditions, we define the center of mass (CM) velocity of the polymer v _CM as the ensemble average (over N _exp realisations) of the ratio between the displacement of the CM of the chain and the corresponding employed time (the translocation time in each realisation). It is worth to remember that the relation between translocation velocity v _CM and mean time τ is not immediate since in general, $\langle 1/t\rangle \ne 1/\langle t\rangle$.

The given definitions are suitable to compare TT of chains with different lengths because the border effects result minimised since five monomers are always inside the pore during the simulation time in all the cases studied.

Following ref. 25, we define m, l ₀, and ${\epsilon }_{0}$ as the mass, the length, and the energy unit units respectively. This choice determines a Lennard-Jones time scale given by ${t}_{LJ}={(m{l}_{0}^{2}/{\epsilon }_{0})}^{\mathrm{1/2}}$. However, as the dynamics we propose is overdamped, the time scale that normalise the equation of motion Eq. (4) is ${t}_{OD}=\gamma {{t}_{LJ}}^{2}$, thus depending on the damping parameter. To set some values, let us consider a DNA molecule at room temperature (k _B T = 4.1 pNnm) and the simplest model with k _b  = 0. We have fixed our simulation temperature to k _B T = 0.1 in dimensionless units. This choice fixes our energy unit in ${\epsilon }_{0}=41\,{\rm{pNnm}}$. By setting l ₀ = 1.875 nm and m = 936 amu²⁵, we obtain ${t}_{LJ}\approx 0.38\,{\rm{ps}}$, while the force unit is given by ${\epsilon }_{0}/{l}_{0}=21.9\,{\rm{pN}}$. An estimation for the kinetic damping is $\gamma \approx 1.6\times {10}^{13}\,{{\rm{s}}}^{-{\rm{1}}}$, so obtaining ${t}_{OD}\approx 2.3\,{\rm{ps}}$. Other normalisations can be used depending on the system to simulate⁴¹.

The dimensionless geometrical values used in the simulations are L _H  = 2, L _M  = 5. The rest distance between adjacent monomers is d ₀ = 1 and k _e  = 1600, large enough to maintain the bonds of the chain rigid enough. The Lennard-Jones energy is $\epsilon =0.3$, and σ = 0.8. The values of d ₀, σ, L _M and L _H guarantees that the polymer is maintained almost linear and ordered inside the pore. Also, the different choices of the bending constant k _b , gives the possibility to study the TT for different persistence lengths of the chain, which basically gives the stiffness of the polymer (i.e. its resistance to bend). For our model L _p  = k _b /k _B T. Thus for example we obtain L _p  = 5 (in units of d ₀) for k _b  = 0.5. Regarding the actuating forces we set F _c  = 0.1 and F _M  = 0.2 along our work.

For every set of parameters, a number of numerical experiments N _g has been simulated. From them a number of N _exp  = 2000 successful translocations have been reached and both the mean translocation time τ and the mean velocity v _CM have been then computed. From these values, we also compute the translocation probability P _in as the ratio N _exp /N _g . To finish we have to mention that since we deal with time dependent pulling forces inside the pore, we have also averaged on the initial value of the driving by randomising the initial phase of the force, that is considered as a random variable with a uniform distribution over its possible values.

Polymer velocity

Figure 3 shows the motor activation rate dependence of v _CM of a polymer with N = 32 monomers at different values of the bending parameter (persistent length equal to 0, 5, 15 and 25 d ₀). We observe that the polymer velocity follows Eq. (8) as a MM law moderated by the contribution of the constant force F _c and the effect of the interaction of the polymer with the membrane walls F _w . We can fit our results to

$${v}_{CM}\simeq \frac{{v}_{HR}}{\mathrm{(1}+{\nu }_{0}/\nu )}+{v}_{LR}$$

(9)

with v _HR the high rate limit, and v _LR the low rate limit of the system. The right inset of Fig. 3 confirms the expected number of monomers inside the motor ${n}_{{\rm{M}}}\simeq 5$, value that is observed in all our simulations.

Figure 3

Mean velocity of the polymer center of mass with N = 32 for different values of the bending k _b  = 0.0, 0.5, 1.5, 2.5. Symbols stand for numerical simulations, full lines show the MM best fit according to Eq. (9). The right inset shows the corresponding mean number of monomers inside the motor n _M during the simulations. Left inset: absolute value of the mean reaction force exerted by the walls during the translocation as a function of the motor actuation rate ν.

The velocity shows a clear MM-like curve, being dominated by the motor force contributions. Regarding the bending dependence in the velocity, as anticipated, bending effects come from the F _w term. In general, to take into account the F _w contribution at different bending, we can write ${v}_{HR}({k}_{b})={F}_{w}(\infty ,{k}_{b})/N+{n}_{M}({F}_{c}+{F}_{M})/N$ and ${v}_{LR}({k}_{b})={F}_{w}\mathrm{(0,}\,{k}_{b})/N+{n}_{M}{F}_{c}/N$. Then, equation (9) allows for a direct experimental fit once the velocity at both high and low ATP concentration are measured at every k _b value.

The force F _w given by the reaction of the walls during the translocation is not easy to evaluate directly but can be computed using Eq. (8). This force term, shown in the left inset of Fig. 3 for different bending parameters, slows down the translocation process and is more intense for stiffer chains. The physical origin of this contribution can be ascribed to the interactions of the polymer with the membrane walls in a direction opposed to the movement. In fact, in the right-hand side of Eq. 8, the firsts two terms weakly depend on k _b (visible also in the right inset of Fig. 3, where a very weak dependence on k _b is present in the number of monomers inside the pore n _M ). Thus, the main contribution with k _b comes from the F _w term, as far as this model can devise. A different approach has been followed by Adhikari and Bhattacharya⁴².

Translocation time

Figure 4 shows the mean of the TT τ of the system as a function of the motor activation rate ν calculated for N = 32. τ decreases monotonously as ν increases and reaches a limit value for large enough values of ν, when the waiting time goes to zero and the motor stays active most of the time. The results are presented for four different bending values and the curves are normalised with respect the high rate value τ _HR in each case.

Figure 4

Mean translocation time τ normalized with its high rate limit τ _HR value, for different values of the bending parameter k _b , with k _b  = 0.0, 0.5, 1.5, 2.5, and N = 32. Right inset: translocation probability P _in . Left inset: Mean translocation time (not normalized). The values at high rate are τ _HR  = 253, 374, 521, 599 time units, for the four bending values, respectively.

Similar behaviour is observed at other studied lengths, namely N = 16, 32, 64, 96 and 128. In all cases we find that τ increases with k _b . In addition, the translocation probability P _in is almost 1 for all rates at high k _b , i.e. almost all the simulated events translocate; but P _in is lower than 1 in a wide moderate-to-low activation rate region, reaching the value P _in  = 0.6 for k _b  = 0. Note that the backward force required to prevent the polymer to translocate has to be stronger for rigid chains than for low k _b values, since the escape from the pore needs a larger displacement of the center of mass of the system to leave the motor, and so, a higher work is required to pull more rigid chains.

Regarding the relation between v _CM and τ, it is observed that ${v}_{CM}\propto 1/\tau$, although with a slope which depends on k _b (see Fig. 5). The dependence with the bending appears more evident in the τ curves than in the velocity ones. In fact, as said, more rigid polymers have to displace a longer distance than softer ones to translocate, and they do it at a slower velocity (see Fig. 3).

Figure 5

Velocity as a function of the 1/τ. The linearity of the curves indicates that the relation $v\approx 1/\tau$ holds for a large range of forces here used.

Figure 6 shows the TT probability distribution for three values of the activation rates (ν = 10⁻², 10⁰, 10²). For each of the four k _b values we can notice a narrower distribution at higher rates (ν = 10²) with a lower mean TT value, with respect to a wider distribution at lower rates (ν = 10⁻²) having a higher mean TT value. All the distributions are almost symmetrical, and this feature can explain the almost linear relationship reported in Fig. 5.

Figure 6

Translocation time distribution P(t) for k _b  = 0.0, 0.5, 1.5, 2.5 with N = 32.

Translocation time and scaling

An interesting issue in this study concerns the scaling properties of the translocation time of the polymer with its size. Below we will propose and check against our numerical results a very simple relation

$$\tau =bN{(N-n)}^{\alpha }$$

(10)

where α is the scaling exponent for the gyration radius of the polymer, n = n _M  + 1 and b the only free parameter of the proposed fitting function. We will see that in spite of its simplicity the proposed equation gives account of the observed scaling behaviour of the system in a wide range of rate values.

To depict a simple intuitive model, let us first consider the polymer translocating in a straight line. The translocation length is then ${L}_{t}={d}_{0}(N-\,1-\,{n}_{M})$, where n _M is the number of monomers remaining inside the pore, and d ₀ is the inter-monomer distance. Assuming that ${v}_{CM}\propto 1/N$ (Eq. 8) we get ${\tau }_{exp}\approx {L}_{t}/{v}_{CM}\propto N(N-n)$, where n = 1 + n _M = 6. For large polymers, the scaling follows the law $\tau \propto {N}^{2}$ (see ref. 43). The above relationship strictly holds for one-dimensional polymers, where no spatial contributions are present to affect the mean translocation velocity.

The three-dimensional stochastic movement of the polymer outside the channel gives rise to a general reduction of the translocation time, which depends on the polymer and environment properties (bending, exclusion volume, temperature).

In our simple approach, we can assume that the translocation length changes because of the conformational shape adopted by the polymer at both sides of the pore. Since the size of the polymer is given by its gyration radius R _g we can roughly assume ${L}_{t}\approx {R}_{g}+({n}_{M}\,-\,\mathrm{1)}{d}_{0}\simeq {R}_{g}$ with ${R}_{g}\propto {(N-n)}^{\alpha }$. Thus $\tau \approx {L}_{t}/{v}_{CM}\propto N{(N-n)}^{\alpha }$.

For large polymer this expression becomes $\tau \propto {N}^{1+\alpha }$, in agreement with the time scaling proposed by Kantor and Kardar⁴⁴ $\tau \propto {N}^{1+{\nu }_{F}}$, with ν _F the so-called Flory coefficient. For a more refined derivation of this case see also ref. 45. In Flory theory, ν _F  = 0.5 for an ideal chain in 3D (without interaction between no consecutive monomers), ${\nu }_{F}\simeq 0.6$ if we consider the excluded volume contribution, and ν _F  = 1 for a rod (rigid chain)^{46,47,48,49,50}.

We have computed (see Appendix) the gyration radius of our system at different bending values and frequencies. We observe it increases with k _b , being $\alpha \simeq 0.5$ for k _b  = 0 and tending to α = 1 as we increase k _b . Regarding its variation with the rate ν, a weak dependence is observed, with α decreasing when the rate increases. Values are given inside the related figure.

In Fig. 7 we compare our theoretical prediction, $\tau \propto N{(N-\mathrm{6)}}^{\alpha }$, with α obtained from the R _g analysis, against our numerical observation for τ. We observe excellent agreement in all the cases, in spite of the different approximations employed in the derivation and the moderate to small size of the chains studied.

Figure 7

τ at different N values. Each panel corresponds to a different k _b and shows three rate values: ν ₁ = 10⁻² (upper curves), ν ₂ = 1 (middle curves), ν ₃ = 10² (lower curves). The full lines stand for τ = bN(N − 6)^α with b the unique fit parameter.

Figure 8 shows the scaled τ(ν) curves for different values of N at three bending values, k _b  = 0.0, 0.5 and 2.5. Here we use the same scaling exponent α for all the rates (we choose the high rate result). As expected, the scaling is excellent at high ν, while the data dispersion at lower rates indicates that the exponent α depends also on the value of the rate, as explained above. It is worth to note, see the insets, that the translocation probability does not change with N, being mostly a function of the bending parameter.

Figure 8

Scaled translocation times at three values of the polymer bending parameter k _b  = 0.0, 0.5 and 2.5. The inset shows the translocation probability P _in for the different polymer lengths.

Persistence length analysis

Persistence length is the basic parameter measuring the stiffness of a polymer. For the continuous WLC model the persistence length is given by L _p  = k _b /k _B T. Here we use a discrete version, so deviations from the continuous model may be found for small N values. In order to study the translocation behavior of the polymer at different values of L _p we have performed simulations for different values of the bending from k _b  = 0 (freely-jointed chain model) to high k _b values (extensible rigid chain). The results have been already presented and commented in previous figures and sections of this article.

In Fig. 9 we plot both the translocation time τ and the center of mass velocity v _CM as a function of L _p . We show the results for the N = 32 polymer and at low, medium and high activation rates (ν = 10⁻², 1, 10²). We also show the constant force result corresponding to the ν→∞ limit. As expected, a saturating behaviour of τ for increasing values of L _p is observed. The inset shows a similar saturating behaviour for the velocity. The curves for ν = ∞ and ν = 10² completely overlap as at this latter value the motor is almost continuously actuating.

Figure 9

Translocation time and center of mass velocity (inset) as a function of the persistence length L _p for the cases of constant applied force (ν = ∞ in the label) and at different values of the motor activation rate ν = 10⁻², 1, 10².

We have seen that the observables of the system, τ and v _CM depends on ν, k _b and N. However, given the definition of L _p  = k _b /k _B T, it is possible to define a new parameter N/L _p giving account of the number of effective segments of the system. The question now is if it is possible to write either τ or v _CM as a function of ν and N/L _p , thus reducing the number of parameters of the system from 3 to 2. We have investigate this issue and the answer is negative. For a given activation rate, it is observed that simulations data at different values of N and k _b do not collapse to a single curve when plotted as a function of N/L _p .

Summary and Conclusions

We have studied the translocation of a polymer chain through a repulsive, uniform pore membrane in a thermal fluctuating environment. The polymer is pulled by a time dependent force modelled as a dichotomous stochastic motor fuelled by the hydrolysis of ATP molecules which bind to the motor at a Poisson rate and activates a mechanical work during a fixed time lapse.

The objectives of our work are to study the specificity of the pull force used inside the pore, and to investigate the translocation time dependence on the flexibility and length of the polymer chain. We propose a scaling behaviour of the translocation times with the polymer length according to a power law of the type N ^1+α, with α weakly dependent on the activation rate of the motor.

This work is the natural extension of our previous studies^{23, 24, 26, 31} to a more realistic and biologically interesting system, accounting for motion in the 3-D domain, polymer-membrane interactions, and ATP based driving motor. Dynamics has been studied for different motor actuation rates, polymer sizes, and rigidities.

The nature of the ATP-based pulling gives rise to monotonic translocation time as a function of the activation rate of the motor. The mean velocity of the polymer is found to have a specific and very clear behaviour, i.e. a MM dependence with the activation rate of the motor, which directly depends on the ATP molecule concentration in the surroundings of the machine. The MM velocity for an ATP-driven machine has been extensively measured in many experiments with different types of motors. Examples are the transport velocity of a stepping motor like kinesins^{51, 52} (with or without load forces) as well as the packaging velocity of the DNA of the bacteriophage pulled inside the virus capside by its motor complex³⁷.

We show here that such MM-like behaviour arises in a natural way in the introduced model. Specifically, the translocation velocity depends on the motor activity which drives the polymer for a certain fixed time, and stays inactive during a certain Poisson distributed waiting time. The average working time of each molecular motors using ATP in these conditions, results in a MM kinetics.

Beyond its biological interest, the study of polymer translocation through active nanopores reveals its importance in understanding the translocation aspects of long molecules assisted by motors at small scales, and represents a starting point for the construction from scratch of nanotechnological objects. In this regard, without any purpose to reproduce in detail some specific experiment, this work is an attempt to model the dynamics at the nanoscale with a very comfortable approach under the realistic constraints imposed by the ATP energy supply of a typical molecular motor. Further improvement of this kind of description can be devoted to a more accurate modelling of the pore channel in both: the structure of the very pore together with its interaction with the passing through molecule, and a more complex characterisation of the polymer by enlarging its internal degrees of freedom such as the diversification of the monomers properties according to the different molecular components of a real chain.

Appendix: Gyration Radius

This appendix is devoted to provide details on the numerical calculation of the exponent α introduced in our scaling analysis of the polymer translocation time. This is based in the computation of a dynamical gyration radius as follows: given a polymer with N monomers we define N _trans as the number of monomers in the trans side of the system (evaluated at the time instant a monomer reaches the side). At this time we compute the instantaneous gyration radius of the system, and average it over N _exp translocation realisation:

$${R}_{g}={\langle\sqrt{\frac{1}{{N}_{{\rm{t}}{\rm{r}}{\rm{a}}{\rm{n}}{\rm{s}}}}\sum _{i=1}^{{N}_{{\rm{t}}{\rm{r}}{\rm{a}}{\rm{n}}{\rm{s}}}}{({{\bf{r}}}_{{\rm{i}}}-{{\bf{r}}}_{CM})}^{2}}\rangle }_{{N}_{exp}}.$$

(A11)

This magnitude scales as R _g  = a(N _trans − 1)^α.

Figure 10 shows the behaviour of this gyration radius as a function of the number of monomers which have crossed the membrane, for the case os a N = 32 beads polymer and at different values of ν and k _b . Fitting is always good and the exponents (given in the figure), decrease a bit increasing the rate, but depends more strongly on the stiffness of the system, k _b . As expected the gyration radius (and the scaling exponent) increases with k _b . The obtained exponents were used in the scaling study of the polymer translocation time as discussed in Section 6.

Figure 10

Polymer gyration radius as a function of the number of monomers reaching the trans region, N _trans, on the translocation of N = 32 polymers. Each point results of an averaging over the N _exp realisations. The curves are fitted to f(N _trans) = a(N _trans − 1)^α.

We have to mention that differently from the usual definition of the gyration radius, which makes sense in thermodynamic equilibrium conditions, here we calculate a somehow dynamical gyration radius. Evidently, the calculation of any magnitude along the translocation process, constitute an out of equilibrium measure, though the slower is the translocation dynamics, the better is the approximation to a thermodynamic equilibrium (thus best fitting is obtained at small activation rates). In this sense is not surprising that the gyration radius is a function of both the activation rate of the motor (and so of the translocation velocity), and, of course, strongly dependent of the stiffness of the chain. That is why we can observe in Fig. 10 that the radius slowly increases by decreasing the activation rates.