THE JOURNAL OF CHEMICAL PHYSICS 141, 214902 (2014)
Detachment of semiflexible polymer chains from a substrate: A molecular dynamics investigation J. Paturej,1,2 A. Erbas,3 A. Milchev,4 and V. G. Rostiashvili5 1
Leibniz-Institut of Polymer Research Dresden, 01069 Dresden, Germany Institute of Physics, University of Szczecin, Wielkopolska 15, 70451 Szczecin, Poland 3 Department of Materials Science and Engineering, Northwestern University, Evanston, Illinois 60208, USA 4 Institute for Physical Chemistry, Bulgarian Academy of Sciences, 1113 Sofia, Bulgaria 5 Max-Planck-Institute for Polymer Research, Ackermannweg 10, 55128 Mainz, Germany 2
(Received 17 August 2014; accepted 14 November 2014; published online 2 December 2014) Using Molecular Dynamics simulations, we study the force-induced detachment of a coarse-grained model polymer chain from an adhesive substrate. One of the chain ends is thereby pulled at constant speed off the attractive substrate and the resulting saw-tooth profile of the measured mean force f vs height D of the end-segment over the plane is analyzed for a broad variety of parameters. It is shown that the observed characteristic oscillations in the f-D profile depend on the bending and not on the torsional stiffness of the detached chains. Allowing for the presence of hydrodynamic interactions (HI) in a setup with explicit solvent and dissipative particle dynamics-thermostat, rather than the case of Langevin thermostat, one finds that HI have little effect on the f-D profile. Also the change of substrate affinity with respect to the solvent from solvophilic to solvophobic is found to play negligible role in the desorption process. In contrast, a changing ratio sB /sA of the binding energies of A- and B-segments in the detachment of an AB-copolymer from adhesive surface strongly changes the f-D profile whereby the B-spikes vanish when sB /sA < 0.15. Eventually, performing an atomistic simulation of (bio)-polymers, we demonstrate that the simulation results, derived from our coarse-grained model, comply favorably with those from the all-atom simulation. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4902551] I. INTRODUCTION
During the last decade a rapid development of the socalled single molecule dynamic force spectroscopy (SMDFS) has enabled the direct observation of the chemical dissociation (e.g., base-pare binding in DNA, or ligand-receptor interaction in proteins) initiated by an external time-dependent force in the pico-Newton range.1–3 Theoretical interpretation of SMDFS for a single bond rupture has been suggested by Bell,5 and developed by Evans.6, 7 The Bell-Evans (BE) approach is built upon an Arrhenius relationship which describes the bond rupture rate (“off”-rate) subject to a timedependent force, koff = κ 0 exp (xβ f/kB T), where κ 0 denotes the rupture rate in the absence of applied force, f is the applied force per bond, and xβ is the coordinate where the activation barrier is located. Here and in what follows, T denotes the temperature, and kB is the Boltzmann constant. In other words, the effective activation energy is represented as a linear function of the force, Eb (f ) = Eb(0) − xβ f . Under the condition of fixed loading rate, f = Rt, it could be shown that f = (kB T /xβ ) ln(Rxβ /κ0 kB T ), i.e., the detachment force grows linearly with the logarithm of loading rate R. This relationship is usually referred to as the Bell-Evans model and employed for the measurement of dynamic strength of molecular bonds, cells’ adhesion, and protein unfolding by means of an atomic force microscopy (AFM).8 However, for multiply bonded attachments the f vs ln R relationship shows a nonlinear behavior which might be related to a more complicated cascade of activation barriers.9 0021-9606/2014/141(21)/214902/10/$30.00
A single-stranded DNA (ssDNA), strongly adsorbed on a graphite substrate, represents an example of a multiplybonded bio-assembly. The desorption of such tethered DNA molecule, induced by the applied force, has been first studied by Jagota et al.10–13 At equilibrium, the macromolecule can be desorbed either by using the displacement of the chain end over the adsorbing surface (and measuring the fluctuating force), or by fixing the pulling force, applied to the chain end, and monitoring the mean displacement of the end segment over the plane. It has been shown (analytically and by means of Brownian simulation) that in the displacement control (DC) desorption, the average force f-displacement D profile exhibits a characteristic set of saw-tooth (force-spikes) oscillations, corresponding to the underlying base sequence of the ssDNA.10–13 When the displacement profile reaches a steady state, i.e., the desorbed monomers are far away from both ends of the chain, each maximum in the saw-tooth oscillations corresponds to an energy barrier that has to be overcome in order to complete the monomer desorption. In a real system, this energy barrier, Gb , is quite complex and composed of various energetic and entropic contributions. For instance, interaction energy between monomers and the surface, Gsurf , conformational entropic contributions and enthalpic energies of polymer chains Gconf , contributions due to direct additive interactions and entropic effects of water molecules near the surface, or with the chain Gsol , etc. One may assume that these energetic components can be decoupled, and the
141, 214902-1
© 2014 AIP Publishing LLC
214902-2
Paturej et al.
overall energy barrier can be expressed as Gb ≈ Gsurf + Gconf + Gsol + Gother , where Gother represents all contributions that cannot be accessed by a course-grained (CG) simulation model (such as hydrogen bonding of water molecules near/around the monomers). In that case, using Molecular Dynamic (MD) simulations, effects of various contributions on the polymer-surface interactions as well as non-equilibrium single molecule experiments can be tackled systematically and in detail. Recently, we have revisited the detachment theory of a strongly adsorbed macromolecule by making use of a free-energy-based stochastic equation (the so called Onsager equation) approach, and by performing extensive MD simulations.14 This study has confirmed the force-spikes response under DC and also demonstrated how the saw-tooth profile is smeared out with growing detachment velocity vc and increasing mass of the AFM-cantilever. Moreover, we have shown that the average detachment force versus detachment velocity vc relationship exhibits a nonlinear behavior when plotted in semilogarithmic coordinates. The presence of fluctuations in our model enables, among other things, to calculate the probability distribution function (PDF) of the fluctuating force at the cantilever, measured at the moment of ultimate detachment, which is an experimental observable in laboratory studies. In the present paper, we extend and generalize our previous MD simulations14 so as to probe systematically the influence of various energy contributions to the desorption energy barrier Gb , more precisely, on the resulting force fdisplacement D profiles. Using a CG model, the pairwise interactions between monomers can be tuned to understand their influence on displacement profiles. Similarly, by introducing torsional and dihedral harmonic potentials in addition to the classic bead-spring potentials (further details on the simulations scheme will be given below), effects of bending, and/or torsional stiffness of the polymer backbone on the detachment behavior are examined. At this point, we should also note that the energy components forming the overall energy barrier Gb can be in phase with each other along the reaction coordinate (distance above the substrate in this case). Hence, their addition can make an energy minimum between two consecutive maxima more shallow or deeper as we will see later in this paper. The paper is organized as follows. First, in Sec. II we examine the role of hydrodynamic interactions (HI) within the context of external force-driven polymer desorption by comparing the effect of Langevin- and Dissipative Particle Dynamics (DPD) thermostats within our CG model. We also check the role of substrate wettability and its impact on the f-D profile. Then, the desorption of an alternating A − B copolymer with different binding energies of the A- and Bmonomers to the substrate is examined. In addition, the effect of bending and/or torsional stiffness of the polymer backbone on the detachment behavior is studied. Eventually, in Sec. III we report on atomistic MD simulation of polypeptide detachment, using poly-glycin and poly-phenylalanine, adsorbed on a crystalline carbon substrate, and compare it to the generic behavior of our coarse-grained model. Our report ends with a brief summary, presented in Sec. IV.
J. Chem. Phys. 141, 214902 (2014)
II. COARSE-GRAINED SIMULATIONS A. Model
Similar to our previous study,14 simulations of a coarsegrained model were carried out based on a generic beadspring model of a flexible polymer chain,15 composed of N monomers, connected by nonlinear bonds along the polymer backbone. The bonded (two-body) interactions in the chain is described by the Kremer-Grest15 potential, V KG (r) = V FENE (r) + V LJ (r) with the so-called “finitely extensible nonlinear elastic” (FENE) potential given by 2 1 2 r FENE = − kr0 ln 1 − . (1) V 2 r0 The non-bonded interactions between monomers were taken into account by means of the Lennard-Jones (LJ) potential, given by V LJ (r) = 4[(σ/r)12 − (σ/r)6 + 1/4]θ (rc − r).
(2)
In Eqs. (1) and (2), r = |rij | denotes the distance between the center of monomer (bead) i and j, rc is the cutoff distance, while the energy scale and the length scale σ are chosen as the units of energy and length, respectively. Accordingly, the remaining parameters are fixed at the values k = 30 /σ 2 and r0 = 1.5 σ .15 In Eq. (2), we have introduced the Heaviside step function θ (x) = 0 or 1 for x < 0 or x ≥ 0. We performed simulations with short- and long-range cutoff: rc = 21/6 σ (purely repulsive interaction between monomers), and rc = 2.5 σ (monomer attractions allowed at larger distances). In the course of the study, chain bending stiffness κ and torsional stiffness κt were varied by introducing a threebody, V b (θij k ) = κ(cos θij k − 1)2 ,
(3)
and four-body interactions V t (φij kl ) = κt (1 + cos φij kl ),
(4)
where θ ijk and φ ijkl denote bending and dihedral angle formed, respectively, by two and three successive bond vectors. In the CG simulations, two kinds of substrates were considered. We employed structureless adsorbing surface (with no friction in the lateral plane), modeled simply by a Lennard-Jones potential acting with strength s in the perpendicular z-direction, V sub (z) = 4s [(σ/z)12 − (σ/z)6 ]. In a separate set of simulations, we introduced a rough surface composed of beads which form triangular lattice and interact with monomers via Eq. (2) in order to take into account friction between polymer and substrate. In our simulations we consider, as a rule, the case of strong adsorption s /kB T = 5 and 20 for the structureless surface, and s /kB T = 5 in case of atomistic surface, with T being the temperature of the thermal bath which is described briefly below. Temperature in our simulations was controlled by two different methods: (i) a Langevin thermostat,16 and (II) by DPD thermostat.17 In both methods, the dynamics of the chain is obtained by solving the following set of equations of motion
Paturej et al.
214902-3
J. Chem. Phys. 141, 214902 (2014)
for the position rn = [xn , yn , zn ] of each bead in the chain, + FD m¨rn = Fcons n n + Rn
(1, . . . , N)
(5)
being the total conservative force acting on each with Fcons n polymer bead with mass m = 1. The influence of the solvent is split into slowly evolving viscous force and rapidly fluctuating stochastic force. Thus, in Eq. (5), FD n and Rn denote, respectively, the dissipative and random forces which are responsible for keeping the system at constant temperature. The difference between Langevin and DPD thermostats lies in the choice of these two forces. In the Langevin thermostat, the dissipative force (drag ˙n, force) is proportional to particle’s velocity FD n = −γL r where γ L = 0.5mτ −1 is the friction coefficient, and the time unit is τ = mσ 2 /. In addition, the random force has a zero mean value and satisfies the fluctuation-dissipation β theorem, Rαn (t)Rm (t ) = 2γL kB T δαβ δij δ(t − t ). These two forces, the random and the frictional one, are balanced in order to maintain the system temperature at the set value. In contrast to the Langevin thermostat, in the DPD thermostat both dissipative and random forces are applied as pairwise interactions, such that the sum of these two forces acting on a given pair of particle in the system is zero. Thus, in the DPD case the particle momentum is conserved, leading to correct description of hydrodynamic interactions.18, 19 The form of dissipative and random D forces in the case 2of DPD thermostat is the following: Fn = m(=n) −γDPD w (rnm )(rnm /|rnm | · r˙ nm )rnm /|rnm | and Rn √ = m(=n) 2kB T γDPD w(rnm )αnm / dt, where the weighting function w is defined as w(r) = 1 − r/rc , γ DPD = 20mτ −1 is the friction coefficient, and α is a Gaussian-distributed random number with zero mean and variance equal to unity, whereas dt stands for the integration step. In all CG simulations, the equations of motion were integrated using the MD package LAMMPS.20 The solvent in our simulations was considered either as being present implic-
itly via Langevin thermostat, or modeled explicitly by adding spherical particles with density ρ = 0.86 σ −3 ). In the latter case, the difference between Langevin thermostat and DPD thermostat was investigated. The detachment of chains, composed of N = 20, or N = 100 monomers, was performed as follows. In the initial state, the chains were completely adsorbed and equilibrated on the surface. The macromolecule was then pulled perpendicular to the adsorbing surface by a cantilever at constant velocity V = [0, 0, vc ]. As a cantilever we used two beads connected by harmonic spring and attached to one of the ends of the chain. The mass of the beads, forming the cantilever, mc , was set to mc = 1 whereas the equilibrium length of the harmonic spring was set to 0 and the spring constant was chosen as kc = 50/σ 2 . During the pulling simulations, the force f(t) at given height D over the substrate was calculated from the instantaneous harmonic linker extension zl (t), i.e., f(t) = kc zl (t). B. Results
1. Impact of chain properties on the f-D profile
The saw-tooth response of the pulling force f, measured at any fixed distance D of the chain end above the adsorbing surface, is a characteristic feature produced by the detachment of successive monomers along the polymer backbone, cf. Figure 1(a). The observed steady and steep increase of f after the characteristic last minimum, which corresponds to detachment of the last monomer, is a hallmark of a chain, tethered by one of its ends to the adsorbing substrate and pulled by the other end monomer. It reflects the ultimate extension of the linker spring once all beads, besides the tethered one, have lost contact with the adsorbing plane. In a real experiment, the chains which have to be detached from the adsorbing surface are most probably not tethered, so the final f-D profile exhibits instead a sharp drop of the force, which can also be found in our simulations as shown in Figure 1(b). This kind
40
40 40
(a)
35
(b)
Implicit
N = 100 N = 20
30
free tethered
30
〈f〉
30
20
25
20
4
8
12
16
〈f〉
〈f〉
10
20
20
D
15
10 10
5
Implicit 0
10
20
30
40
50
D
60
70
0 80
90
100
2
4
6
8
10
12
14
16
18
20
22
D
FIG. 1. Short vs long chain desorption in the case of implicit solvent. (a) Mean detachment force f vs distance D profile is shown for a semi-flexible polymer chain (stiffness parameter κ = 50). Results pertain to chains composed of, respectively, N = 20 and 100 monomers. For better visibility, the onset of the desorption process is zoomed in the inset. (b) Detachment of a chain whose end is not tethered to the substrate. In both figures, s /kB T = 20 and vc = 10−3 σ/τ .
Paturej et al.
214902-4
J. Chem. Phys. 141, 214902 (2014)
30
30
κ = 0 lp = 1.7 κ = 20 lp = 8.5
25
(a)
κ=0 κ=5 κ = 20 κ = 50 κ = 250 κ = 1000
25
κ = 50 lp = 11.8 κ = 150 lp = 22.8
20
〈f〉
〈f〉
20
15
(b)
15
10 10 5 5
Implicit
Implicit 0
2
4
6
8
10
12
14
16
18
20
22
24
D
2
4
6
8
10
12
14
16
18
20
D
FIG. 2. Force-distance diagram for semiflexible chains in implicit solvent. Stiffness parameter κ and persistence length lp are indicated in the legend: (a) Chain length N = 100. With growing stiffness κ the spikes amplitude decreases. (b) Chain length N = 20. The end of the detachment process is marked by a force minimum which becomes more pronounced with growing stiffness. Here, s /kB T = 20 and vc = 10−3 σ/τ .
of behavior has been observed before for fully flexible hom*opolymer chains by Jagota et al.10–13 and Paturej et al.14 Apparently, a semi-stiff tethered chain exhibits the same pattern, (cf. Fig. 1(a)), albeit the ultimate steep growth of the pulling force f is preceded by a characteristic minimum extending over the last few beads of the chain that precede the tethered bead. Obviously, the last portion of the semi-rigid chain bends and takes off as a whole, whereby the length of this chain portion should depend on the chain stiffness. No finite size effect can be detected on Fig. 1(a). Indeed, it can be seen that, irrespective of the chain length (either N = 20, or N = 100), the amplitude and position of the spikes remain insensitive to chain length N. In addition, there is no significant difference between the f-D profiles for polymer chains adsorbed on a smooth or rough surface (not shown). The only difference is a slightly larger (≈5%) amplitude of force f measured in the case of rough surface. Generally, one may speculate how the chain stiffness affects the f vs D diagram. The bending stiffness of the polymer backbone was included in the simulation by allowing for the three-body bending potential, Vb (θ ), where θ is the angle between two consecutive bonds, cf. Eq. (3). The stiffness parameter κ determines the persistence length of the chain, lp , which is defined through the decay of bond angle f
(a)
correlations.21 Figure 2 shows simulation results for different stiffness parameter κ and chain length N. It can be seen that with growing stiffness, the amplitude of the spikes decreases and also the spikes resolution deteriorates. Starting with approximately κ = 50, the complete chain detachment is preceded by a characteristic force minimum (in case of tethered chains). Then, with growing stiffness κ this minimum, which marks the ultimate chain detachment, occurs at smaller D, indicating that increasingly larger portions of the polymer backbone are detached as a whole. The minimal detachment force thereby also drops so that at κ = 1000 it approaches zero. However, for realistic values of the stiffness κ, and not extremely short chains, N ≥ 10, the total energy, Erod , needed for tearing off the polymer as a single piece of rod from the adsorbing plane, Erod ∝ N s , would be huge in comparison to the energy, Earc , needed to tear off the same semi-flexible chain bead by bead, Earc ≈ nm s + (nm − 1)κθb4 , with few beads nm that form an arc of the bended portion of the chain backbone, Fig. 3. Once an arc is formed, no further energy penalty will be needed to keep the chain bended during the rest of the detachment process. A rough estimate, using the potential V b (θ ), Eq. (3), with θ = 10◦ yields κ × 5.3 × 10−8 bending energy per bond, or
f
(b)
θ θ θ
FIG. 3. (a) Schematic picture of a force-induced detachment of a semi-stiff (left) and flexible (right) chain from a solid plane. The bending angle θ between successive bonds is indicated. Light-shaded area denotes the range of the adsorption potential V sub (z). (b) Snapshots of partially detached polymer chains of different stiffness κ as indicated. In each case, 66% monomers were peeled off the substrate.
214902-5
Paturej et al.
J. Chem. Phys. 141, 214902 (2014)
20
in a semi-stiff chain detach concertedly rather than one by one which exerts a smearing effect on the saw-tooth diagram f − D, Fig. 2(a). Increasing bending stiffness also decreases the magnitude of saw-tooth profile amplitude as the bonds between neighboring monomers stretch less. Realistic DNA models usually also include a dihedral potential which is responsible for the chain resistance to torsion. We have used the dihedral potential, Vt (φ), where φ is the dihedral angle and the torsion constant is κ t = kB T, Eq. (4). The resulting f vs. D diagram for a chain with torsional and bending finite stiffness, compared to the fully flexible, and semi-flexible (κ = 50) chain models, is shown in Fig. 4. It can be seen that while the bending stiffness itself leads to a clear shift of the force oscillation pattern, the resulting behavior practically does not change upon inclusion of a dihedral potential.
flexible semiflexible torsion
〈f〉
15
10
5
Implicit 0
1
0.5
1.5
2
2.5
3
4
3.5
4.5
D FIG. 4. Force f vs D diagram for flexible, semiflexible, and torsional (angle and dihedral potential included) chains as indicated. Here, N = 20, s /kB T = 20, vc = 10−3 σ/τ , κ = 50 kB T, and κ t = kB T.
Earc ≈ 10−5 kB T for an arc encompassing 10 bending angles and κ = 100 kB T. In the same time, such an arc can already reach a height of 10 × cos (10◦ )σ ≈ 1.74σ , where surface adhesion is already dwindling. Therefore, the detachment of even rather stiff chains instantaneously as a rod-like object should be ruled out and a f-D profile of the type, shown in Figure 2, is likely to be observed. Owing to the creation of arc, more that one bead detach concertedly and move away from the adhesive surface beyond the range of adsorption. The neighboring beads along the arc, i and j, remain thereby at fixed mutual distance zij < σ . The neighboring bonds slightly bend but do not stretch significantly so that the length of the individual bond is close to the unperturbed length yet much less than the maximal one, r0 = 1.5σ , cf. Eq. (1). As a result, depending on chain stiffness and in contrast to flexible chains, the monomers
2. Effects of substrate adhesion on polymer detachment
It is to be expected that the strength of adhesion of the polymer chain to the adsorbing surface will manifest itself in the recorded variation of desorption force f with distance D. While in the previous graphs, Figures 1–4, we focused on cases of strong adsorption, s = 20kB T, in Figure 5(a) we present the desorption profile for weak to moderate attraction of the chain by substrate (in our model the threshold for adsorption scrit ≈ 3kB T ). Indeed, as indicated in Figure 5(a), at s = 5kB T, the characteristic oscillations in the f-D profile virtually vanish (apart from the statistical noise). Therefore, one may conclude that the method of single chain detachment spectroscopy as a tool for sequencing analysis could be used in cases of strong polymer-substrate adhesion only. For the objectives of sequencing, the legibility of the data, derived by this method of force-induced detachment, must be examined for heterogeneous polymers in particular. As an example, the result for an alternating (A − B)-copolymer detachment is shown in Figure 5(b). Here, monomers A and B 35
7
B
(a)
κ = 50 κ = 250
6
5
B A εs /εs B A εs /εs B A εs /εs
25
4
〈f〉
〈f〉
20 3
A
εs /εs = 1.00
30
(b)
= 0.75 = 0.50 = 0.10
15 2 10 1 5
Implicit
0 0
2
4
Implicit 6
8
10
D
12
14
16
18
2
4
6
8
10
12
14
16
18
20
D
FIG. 5. (a) Force f vs D diagram of a hom*opolymer for attraction strength of the surface s = 5kB T and two different values of the bending stiffness, κ = 50, 250 kB T. (b) The same for alternating copolymers made of two types of monomers A and B which have different binding energies to the substrate, sA and sB , respectively. Results are displayed for different ratios sB /sA as indicated in the legend. The case of sB /sA = 1 corresponds to a hom*opolymer. The absolute value of sA here is 20 kB T. In both figures, N = 20 and vc = 10−3 σ/τ .
Paturej et al.
214902-6
J. Chem. Phys. 141, 214902 (2014)
35
40
(a)
Implicit Langevin Explicit Langevin Explicit DPD
30
wettable vc=10
35
wettable vc=10 30
repulsive vc=10
25 25
repulsive vc=10
-3
(b)
-1 -3 -1
〈f〉
〈f〉
20 20
15
15
10
10 5
5
Explicit 0
0 0
2
4
6
8
10
12
14
16
18
20
2
4
6
8
D
10
12
14
16
18
20
D
FIG. 6. (a) Comparison of the force f vs D diagram for implicit and explicit solvents from simulations with Langevin, and DPD thermostats, respectively. Results are presented for fast pulling, vc = 10−1 σ/τ . (b) Impact of the substrate selectivity for the case of explicit solvent and DPD thermostat. Force f vs. displacement D diagram for wettable and solvophobic substrates at two pulling velocities. In both figures, N = 20 and s /kB T = 20.
have different affinity to the substrate (albeit the same mass mA = mB ). The detachment starts with a B-monomer (i.e., a monomer with a relatively smaller affinity to the substrate). For a ratio of sB /sA = 0.5, the alternating pattern of spikes can still be clearly seen. As sB gradually further declines, the set of force maxima, corresponding to the desorption of Bmonomers, decreases significantly in amplitude. Eventually, for sB /sA = 0.1, the maxima corresponding to tearing-off Bmonomers turn into minima. Moreover, the latter effect is observed even at higher values of the sB / A -ratio once the first few repeating units are been detached. So at larger height D of the pulled chain end (approximately starting from D = 8σ ), for the B-type monomer desorption one observes local minima rather than peaks. Evidently, a correct sequencing of heterogeneous macromolecules can be performed only in cases when the affinity of the various building blocks is large in terms of absolute values of binding energy ( > 10 kB T) and the differences between values of binding energy should be significant. 3. Implicit vs. explicit solvent
So far we examined how the force-displacement profile of different polymer chains reflects the properties of the chains and their interaction with the adsorbing surface. It is of some interest to check whether the properties of the surrounding medium, considered in the different simulation setups, might influence the f-D profile too. In many computational experiments, e.g., in our previous publication,14 one takes the solvent only implicitly into account. The solvent properties can be then varied to a limited amount only, for example, by changing the friction coefficient γ in Eq. (5). In principle, however, the presence of an explicit solvent might affect the course of force-induced chain desorption due to HI. To this end, we compare the f vs D diagram, derived from simulations of the same system when two different thermostats are used: (i) a Langevin thermostat (i.e., with no HI), and (ii) a DPD thermostat. It is well known that the latter allows for a correct hydrodynamic behavior,18, 19
whereas the Langevin thermostat does not exhibit momentum conservation, and therefore does not reproduce the proper hydrodynamic behavior. In the case of explicit solvent, a chain is pulled in a Lennard-Jones liquid with liquid monomer density ρ = 0.86σ −3 . It is evident from Figure 6(a), however, that there is no tangible difference between these cases, which suggests that the hydrodynamic interaction is largely irrelevant in detachment experiments. In fact, when compared to the case of explicit solvent with no HI (Langevin thermostat), the presence of HI (accounted for by DPD) leads to a slight decrease in the pulling force, cf. Figure 6(a), at the end of the detachment process. Evidently, this affects the detachment of the last beads only while the main portion of the chain is sufficiently far away from the substrate. While such an effect is completely missing for slow detachment with vc = 10−3 σ/τ (not shown here), for fast detachment, vc = 10−1 σ/τ , it should actually be expected due to the solvent back-flow, triggered by Stokes friction of the detached chain portion when the desorbed chain eventually sets into motion. Moreover, due to confinement effects, HI (being long-ranged) are screened27 in the vicinity of the adhesive wall, which explains why their presence is detectable only at sufficiently large distance D from the wall. Therefore, one can view the small decrease in f in the case of explicit solvent as a typical manifestation of the well-known difference between Rouse and Zimm dynamics of polymers.
4. Substrate wettability
We also checked to what extent the affinity of the adsorbing substrate with respect to solvent plays a role. Basically, we distinguish between solvophobic (repulsive) substrates, where the polymer chain is still attracted to the surface whereas the solvent particles are repelled, and wettable substrates, where both the polymer and the solvent are attracted to the surface. Figure 6(b) indicates a systematic decrease in the amplitude of the spikes for wettable substrates (as if the solvophobic solvent effectively increases the chain adhesion to the
Paturej et al.
J. Chem. Phys. 141, 214902 (2014)
surface), yet the characteristic saw-tooth force vs distance profile remains qualitatively unchanged. The spikes positions under conditions of good wetting are also slightly shifted to lower values of D, i.e., the monomers detach more easily (at somewhat lower height) as compared to the solvophobic case. This is because as the relative attraction of solvent particles to the surface is increased, the solvent particles try to replace chain monomers and form a solvent layer on the surface, which in turn, facilitates the desorption of chain monomers. In the context of protein-surface interactions, this effect is referred to as the Berg limit and was also observed in the simulations of biopolymers on various hydrophobic/philic surfaces by Schwierz et al.32
80
60
50
40
30
20
10
5. Desorption at different temperature
Eventually, we examined the role of temperature T in the process of chain detachment and its impact on the forcedisplacement diagram. T (measured in units of the monomermonomer interaction strength /kB ) was increased, while the adsorption strength s was also correspondingly changed so as to keep the ratio s /kB T constant and equal to s /kB T = 20 as in most of the presently studied cases. The presence or absence of explicit solvent revealed thereby almost no difference again. Expectedly, the mean level of the force (the average plateau height) grows, reflecting the stronger adhesion s , while, surprisingly, the amplitude of the spikes remains largely unchanged. A simple explanation for this observation in the force detachment experiment can be suggested as follows. We consider the chemical potential of adsorbed, μads , and detached segments, μdet , which should become equal on the detachment line, i.e., μads = μdet . In the limit of strong adsorption (or, at low temperature), when the macromolecule is tightly bound to the surface and loops (non-adsorbed chain portions) may be neglected, the free energy gain (per chain segment) upon adsorption reads μads =
− −kB T ln(μ2 /μ3 ) . s energy gain
(6)
entropy loss
In Eq. (6), s again stands for the adsorption energy of a single segment while μ2 and μ3 are the so called connective constants in two- and three-dimensional space, respectively. The latter correspond roughly to the possible orientations of a chain segment in space, i.e., the logarithms thereof yielding effectively the entropy contributions in two- and threedimensions. It has been shown that for cubic lattices, for instance, μ2 = 2.6 and μ3 = 4.68.28 On the other hand, in the limit of strong adsorption, the detached chain portion is strongly stretched, attaining a “string” configuration, so that the elastic free energy per segment reads μdet = −af, where a is the Kuhn length and f is the force acting on the chain end. Moreover, in the “string” state a segment has only one orientation, i.e., μ3 = 1. On the f − D plateau, f = fp , and due to the condition μads = μdet one has the following “plateau”-relationship: afp = s + kB T ln μ2 .
(7)
kBT = 5ε kΒΤ = 2ε kBT = ε
Implicit
70
〈f〉
214902-7
2
4
6
8
10
12
14
16
18
20
22
24
D FIG. 7. (a) Comparison of pulling results performed at different temperatures T (see legend) and constant ratio s /kB T = 20 and N = 20.
This result shows that for a strong adsorption the plateau height (i.e., the pulling force) is proportional to the adsorption energy. The result given by Eq. (7) has been obtained first within a more general consideration in our paper29 (see Eq. (30) in Ref. 29). This is now supported by Figure 7 where the temperature and adsorption energy are changed proportionally to one another. The f − D diagram demonstrates in all cases the characteristic saw-tooth behavior with the amplitude progressively decaying in the course of chain detachment. This behavior has been analyzed first by Jagota et al.10 In terms of the number of detached chain segments, the nth spike correspond to the reversible transition n↔n + 1 during which the detachment of a segment leads to release of polymer stretching energy back to the energy of adsorption s . This condition leads to the spikes amplitude law10 famp ∼ exp[(s /kB T − ln 4π )/n].
(8)
This relationship is clearly in line with the decay behavior upon growing n. On the other hand, provided that the temperature and adsorption energy are increased proportionally to each other, the spikes amplitude does not change. That is exactly what we observe in Figure 7. Figure 7 also shows that the overall elastic modulus of the tethered chain does not depend on temperature. Most probably, this is due to the fact that for the strong stretching the entropic contribution to elastic modulus (modulus grows with temperature) is compensated by the bond anharmonicity effect when the elastic modulus decreases with temperature. III. ALL-ATOM SIMULATIONS A. Model
It appears instructive to compare the obtained simulation results for a coarse-grained model to those from an atomistic model of a concrete macromolecule. The latter were performed with the Gromacs MD package4 using the Gromos96 force field22 and the SPC/E (Single Point Charge/Extended)
Paturej et al.
J. Chem. Phys. 141, 214902 (2014)
Poly-Glycine
vc
f [kJ/mol/nm]
214902-8
vc
250 200
f(t) 〈f〉
(b)
150 100
(a)
vc
〈f〉
Poly-Phenylalanine
[kJ/mol/nm]
50 0
1
2
300
3
4
5
6
7
8
9
10
4
5
6
7
8
9
10
(c)
225 150 75 0
1
2
3
D [nm]
FIG. 8. Biomolecule (polypeptide) desorption. (a) Vertical pulling of 31-glycine and 31-phenylalanine chains adsorbed on hydrophobic diamond surface. Small red dots represent water molecules. (b) Force f vs D diagram for 31-glycine desorption. Red line represents an average over 25 simulated desorption events whereas the black line denotes a single run simulation data. (c) Averaged force f vs D diagrams for: 31-glycine (red line), 6-(Gly-Gly-Gly-Gly-Phe) (black line), and 31-phenylalanine (green line). Here, kc = 200, s /kB T ≈ 7, and the pulling velocity vc = 1 m/s.
water model23 at constant surface area A and at constant vertical pressure Pz of 1 bar with temperature T = 300 K. For the temperature and pressure control, the method of Berendsen24 was used. Periodic boundary conditions for the Coulomb interactions were implemented by the particle-mesh Ewald method.25 Simulation runs were performed with an integration time step equal to 2 fs. The simulation box contains a hydrophobic diamond slab with a water-surface contact angle θ c ≈ 90◦ ,30 a single N = 31 amino acid (AA) chain, and ∼16 000 SPC/E water molecules. The entire system, including diamond surface, is composed of 70 000 atoms. Dimensions of the simulation box are around 7 nm× 7 nm× 12 nm, and the thickness of a diamond slab is 1.8 nm. The ratio between adsorption strength of the surface and thermal energy is around 7 − 10. This binding energy is of the order of the binding energy 8.3 ± 0.7 of polythymine 3 poly(dC50 ) on graphite substrate.11 For the hydrophobic surface, the 100-plane of an elastic diamond substrate is saturated completely by uncharged hydrogen atoms. The peptide is allowed to adsorb on the surface prior to solvation by water. After equilibration of the chain on the surface, the pulling is performed, cf. Figure 8(a), whereby initially the molecule is completely adsorbed on surface. Similar to CG simulations, the molecules were pulled via a harmonic linker attached to one of the ends of the molecule as shown in Figure 8. The harmonic spring is moved vertically at a prescribed velocity vc until the entire chain is desorbed from the surface. The spring exerts no lateral force on the chain, and the chain can move freely on the surface. The spring constant is chosen as kc = 300 pN/nm. The force needed to pull the peptide vertically is calculated via F = kc (vc t − zt ), where zt is the z-component of the position of the terminal amino acid. The velocity is taken as
vc = 1 m/s, which has been demonstrated earlier to provide a quasi-equilibrium pulling on the corresponding surface.26 Indeed, a pulling velocity vc 1 can also be justified by scaling arguments borrowed from the polymer physics: When HI is included, the relaxation time of chain with N Kuhn monomers is τ Z ≈ N3ν τ 0 ,21 where the scaling exponent ν ≈ 3/5 for good solvent, τ 0 ≈ a2 γ 0 /kB T is the relaxation time of a Kuhn segment with a monomeric friction coefficient in the solvent γ 0 . If τ Z is smaller than the time scale imposed by the pulling τc ≈ a/vc , i.e., vc vc∗ ≡ kB T /aN 3ν γ0 , the undesorbed section of the chain will be in quasi-equilibrium for which pulling forces should not depend on the conformation of the remaining chain section on the surface. If we take the Kuhn segment size (twice the persistence length) of an AA chain as a ≈ 1 nm (i.e., 1-3 AA monomers), N ≈ 10, and a monomeric friction coefficient of γ 0 ≈ 10−12 kg/s for the diamond surface31 and kB T ≈ 4 × 10−21 kg m/s2 , we obtain a threshold velocity vc∗ ≈ 1 m/s. Note that vc∗ will be much lower if the monomeric friction coefficient is γ 0 10−12 kg/s, e.g., for a OH saturated surface31 or if the chains are longer. Above argument for the chain relaxation on the atomistic surface also applied to our CG simulations since in most of our CG simulations (except for the rough surface), the pulled chains interact with the surface only via a z-dependent potential. This means that chains can laterally diffuse on the surface but with a bulk diffusion coefficient (see Sec. II A). B. Comparison with coarse-grained simulations
The all-atom simulations of polypeptide desorption from atomistically rough substrate were performed for a N = 31 polyglycine (31-glycine), N = 31 polyphenylalanine
214902-9
Paturej et al.
(31-phenylalanine). We also constructed a N = 31 heteropeptide composed of six (Gly-Gly-Gly-Gly-Phe) groups, where Gly and Phe stand for glycine and phenylalanine monomers, respectively. The results for 31-glycine chain desorption are presented in Figure 8(b) along with several snapshots from different stages of the desorption event. The snapshots taken from 31glycine pulling trajectories compare well to those obtained from coarse-grained simulations shown in Fig. 3. Comparing the first few spikes in Figure 8(b) and in Figure 4, the general pattern resembles much more the saw-tooth profile, typical for semi-rigid, rather than for completely flexible polymers, in line with the nature of this polypeptide. Although the data for the f-D diagram, shown in Figure 8(b), are averaged over 25 simulation runs only and appear somewhat noisy, they reveal a characteristic f − D behavior, which qualitatively complies with the results from our coarse-grained simulation. Due to the relatively short length of this glycine macromolecule, the final “dip” before contact with the adhesive substrate is lost, is rather short yet clearly visible as in the CG f-D diagrams of tethered chains. To compare 31-glycine force trace with a stiffer chain, desorption simulations of 31-phenylalanine (red data in Figure 8(c)) were performed. As seen in Figure 8(c), the saw-tooth peaks are more visible, and the peak forces are much higher than those observed for 31-glycine cases. This is actually due to the large benzyl side chain of phenylalanine monomers: The hydrophobic nature of the side chain increases the affinity of phenylalanine monomers to the hydrophobic surface, hence, results in higher force peaks. Interestingly, visual inspection of our simulation trajectories revealed that the benzyl side chains force the overall 31phenylalanine molecule to take on a rod-like structure on the surface (see the snapshot in Figure 8(a)). However, the pulling snapshots shown in Fig. 8(a) show that the conformation of 31-phenylalanine during the pulling resembles more that of the 31-glycine rather the illustration shown in Fig. 3 for the CG model with κ 25. We attribute this to the similar atomistic AA backbone structure of both chains. The f-D diagram of our hetero-peptide chain composed of (Gly-Gly-Gly-Gly-Phe) groups is shown in Figure 8(c) (black data). One can distinguish individual desorption peaks for 6 phenylalanine monomers separated by peel-offs of glycine monomers which is also observed in CG simulations of alternating polymers (see Figure 5(b)). The phenylalanine-induced force peaks observed in the force trace of hetero-peptide chain are lower than those observed for the 31-phenylalanine chain itself. This 2-fold difference in the peak forces can be due to complex interplay of chain stiffness and the relative surface affinity of monomers with respect to neighboring monomers: Possibly, glycine monomers might decrease the adsorption energy of adjacent phenylalanine monomers since they can diffuse faster due to their relatively small sizes (the side chain of a glycine is one hydrogen). This observation in Fig. 8(c) hints that the adsorption energy per AA residue can have a dependence on sequence and deserves further investigation in future. Also note that the maxima of saw-tooth in Fig. 5(b) for CG model show a ten-
J. Chem. Phys. 141, 214902 (2014)
dency to decrease as the adhesion asymmetry of alternating monomers grows. Overall, by comparing CG and atomistic simulations, one may conclude that the coarse-grained modelling of forceinduced desorption of a polymer chain from adhesive substrate agrees well with the results from all-atom simulations. IV. SUMMARY
In the present investigation, we studied the process of polymer chain detachment by an external force, applied to the end-segment of a semiflexible chain which is strongly adsorbed to adhesive substrate. Most of the results have been derived by means of Molecular Dynamics simulations of a coarse-grained bead-spring model of a polymer chain, and focused on the analysis of the recorded (fluctuating) mean force f at height D of the last segment of the chain above the adsorbing plane when the segment is pulled with given velocity vc . As a principal objective of this investigation, the influence of different parameters that characterize the polymer chain, its adhesion to the substrate, and the substrate-solvent affinity on the ensuing f-D diagram have been examined. We have found that an increasing bending rigidity of the polymer induces a sharp drop of the pulling force before the last segments of the chain are peeled off, i.e., the final portion of the chain is detached as a single piece of rod. Nonetheless, our observations and estimates suggest that the sequential desorption of polymer repeatable units from the substrate retains its characteristic “unzipping” mechanism, reflected by the observed “saw-tooth” f − D profile, up to very high degree of rigidity. This mechanism works not only for fully flexible chains but also for rather stiff ones due to the gradual bending of the macromolecule which is energetically much more favorable. We also find that with increased bending stiffness κ, the modulation of the characteristic oscillatory profile steadily declines, similar to the effect of weaker attraction s of the chain to the adsorbing surface where the spikes vanish already at s ≈ 5kB T. In contrast, the torsional stiffness of the polymer has little or no effect of the f-D diagram. Regarding the possible use of the f-D diagram for sequencing and its legibility, the performed detachment of an A − B-copolymer indicates that the ratio sB /sA of binding energies of the A- and B-segments strongly influences the resulting oscillatory profile so that when sB /sA < 10% the spikes that refer to B-atoms practically disappear. Our studies indicate that the role of HI in the process of forced-induced detachment of a macromolecule from adsorbing surface is negligible. The resulting f-D diagrams, emerging from MD simulations with and without explicit solvent, hardly warrant the incomparably larger computational efforts in the former case. This insensitivity of the problem regarding HI is related most probably to the resulting stretched conformation of the pulled macromolecule, and to the effect of screening of HI in the vicinity of the adsorbing surface. Eventually, by comparing our data derived from a coarsegrained bead-spring model of a macromolecule to data from a realistic all-atom simulation of various bio-polymer (e.g., glycine, phenylalanine), peeled off a hydrophobic diamond
214902-10
Paturej et al.
substrate, we have demonstrated that the observed f-D diagrams agree qualitatively well with each other, underlying thus the relevance of coarse-grained computer modelling.
ACKNOWLEDGMENTS
This research was supported by the Polish Ministry of Science and Higher Education Grant Iuventus Plus Project No. IP2012 005072.) 1 R.
Merkel, Phys. Rep. 346, 343 (2001). Ritort, J. Phys.: Condens. Matter 18, R531 (2006). 3 I. Franco, M. A. Ratner, and G. C. Scatz, “Single-molecule pulling: Phenomenology and interpretation,” in Nano and Cell Mechanics: Fundamentals and Frontiers, Microsystem and Nanotechnology Series (Wiley, 2013), Chap. 14, pp. 359–388. 4 E. Lindahl, B. Hess, and D. van der Spoel, J. Mol. Model. 7, 306 (2001). 5 G. I. Bell, Science 200, 618 (1978). 6 E. Evans and K. Ritchie, Biophys. J. 72, 1541 (1997). 7 E. Evans, Annu. Rev. Biophys. Biomol. Struct. 30, 105 (2001). 8 Atomic Force Microscopy in Liquid: Biological Applications, edited by A. M. Baró and R. G. Reifenberger (Wiley-VCH Verlag & Co. KGaA, Weinheim, 2012). 9 R. Merkel, P. Nassoy, K. Ritchi, and E. Evans, Nature (London) 397, 50 (1999). 10 S. Manohar and A. Jagota, Phys. Rev. E 81, 021805 (2010). 11 S. Manohar, A. R. Manz, K. E. Bancroft, Ch.-Y. Hui, A. Jagota, and D. V. Vezenov, Nano. Lett. 8, 4365 (2008). 12 S. Iliafar, D. V. Vezenov, and A. Jagota, Langmuir 29, 1435 (2013). 2 F.
J. Chem. Phys. 141, 214902 (2014) 13 S.
Iliafar, K. Wagner, S. Manohar, A. Jagota, and D. V. Vezenov, J. Phys. Chem. C 116, 13896 (2012). 14 J. Paturej, J. L. A. Dubbeldam, V. G. Rostiashvili, A. Milchev, and T. Vilgis, Soft Matter 10, 2785 (2014). 15 K. Kremer and G. S. Grest, J. Chem. Phys. 92, 5057 (1990). 16 T. Schneider and E. Stoll, Phys. Rev. B 17, 1302 (1978). 17 R. D. Groot and P. B. Warren, J. Chem. Phys. 107, 4423 (1997). 18 I. Pagonabarraga, M. H. J. Hagen, and D. Frenkel, EPL 42, 377 (1998). 19 T. Sodemann, B. Dünweg, and K. Kremer, Phys. Rev. E 68, 046702 (2003). 20 S. J. Plimpton, J. Comput. Phys. 117, 1 (1995). 21 A. Yu. Grosberg and A. R. Khokhlov, Statistical Physics of Macromolecules (AIP Press, New York, 1994). 22 W. R. P. Scott, P. H. Hünenberger, I. G. Tironi, A. E. Mark, S. R. Billeter, J. Fennen, A. E. Torda, T. Huber, P. Krüger, and W. F. van Gunsteren, J. Phys. Chem. A 103, 3596 (1999). 23 H. J. C. Berendsen, J. R. Grigera, and T. P. Straatsma, J. Phys. Chem. 91, 6269 (1987). 24 K. A. Feenstra, B. Hess, and H. J. S. Berendsen, J. Comput. Chem. 20, 786 (1999). 25 T. Darden, D. York, and L. Pedersen, J. Chem. Phys. 98, 10089 (1993). 26 D. Horinek, A. Serr, M. Geisler, T. Pirzer, U. Slotta, S. Q. Lud, J. A. Garrido, T. Scheibel, T. Hugel, and R. R. Netz, Proc. Natl. Acad. Sci. U.S.A. 105, 2842 (2008). 27 A. Winkler, P. Virnau, K. Binder, R. G. Winkler, and G. Gompper, Europhys. Lett. 100, 46003 (2012). 28 C. Vanderzande, Lattice Model of Polymer (Cambridge University Press, Cambridge, 2004). 29 J. Paturej, A. Milchev, V. G. Rostiashvili, and T. A. Vilgis, Macromolecules 45, 4371 (2012). 30 F. Sedlmeier et al., Biointerphases 3, 23–39 (2008). 31 A. Erbas, D. Horinek, and R. R. Netz, J. Am. Chem. Soc. 134, 623 (2012). 32 N. Schwierz et al., J. Am. Chem. Soc. 134, 19628–19638 (2012).
The Journal of Chemical Physics is copyrighted by the American Institute of Physics (AIP). Redistribution of journal material is subject to the AIP online journal license and/or AIP copyright. For more information, see http://ojps.aip.org/jcpo/jcpcr/jsp