Atomic-Scale Modeling of the Interaction between Short Polypeptides

Ylichron Srl, c/o ENEA Casaccia Research Centre, Via Anguillarese 301, 00123 S. Maria di Galeria (RM), Italy, Dipartimento Tecnologie Fisiche e Nuovi ...
6 downloads 0 Views 3MB Size
Download PDF
J. Phys. Chem. B 2009, 113, 12105–12112

12105

Atomic-Scale Modeling of the Interaction between Short Polypeptides and Carbon Surfaces Giulio Gianese,*,† Vittorio Rosato,†,‡ Fabrizio Cleri,§ Massimo Celino,‡ and Piero Morales‡ Ylichron Srl, c/o ENEA Casaccia Research Centre, Via Anguillarese 301, 00123 S. Maria di Galeria (RM), Italy, Dipartimento Tecnologie Fisiche e NuoVi Materiali, ENEA Casaccia Research Centre, Via Anguillarese 301, 00123 S. Maria di Galeria (RM), Italy, and Institut d’Electronique, Microelectronique et Nanotechnologie (CNRS UMR 8520), UniVersite´ de Lille I, 59652 VilleneuVe d’Ascq, France ReceiVed: April 21, 2009; ReVised Manuscript ReceiVed: July 10, 2009

We performed a comparative study of the adsorption of an in Vitro selected peptide on two different carbon surfaces: a flat graphene and a curved (0,15) nanotube. The sequence was selected from recent experiments, as the one giving the highest carbon affinity for carbon nanotubes. Rigid docking of the molecule on the two surfaces by a genetic algorithm was followed by molecular dynamics simulations with empirical force fields (OPLS-AA) in water at finite temperature. The total free energies of folding and adhesion and the quality of surface binding were determined, based on a combination of solvation energy, formation of hydrogen bonds, and amount of the apolar (hydrophobic) contact surface between peptide and carbon surface. For both cases, we find a strong adhesion energy and large nonpolar contact surface. Isoleucines and tryptophans are the most strongly bound residues to the two carbon surfaces, the latter one largely dominating. In the case of the carbon nanotube, the peptide shows several competing stable structures, corresponding to different possible molecular foldings, and a propensity to enhance the intramolecular stability by forming new hydrogen bonds. In both systems, different arrangements of the histidine and tryptophan residues enable a better adaptation to the carbon surfaces. These findings suggest that the experimentally observed surface specificity of the peptide on nanotubes may depend on its capability to support multiple strongly bound configurations. 1. Introduction The design of novel materials engineered at the nanometer scale opens an unprecedented host of opportunities, for both fundamental science and technological applications. Nanostructures can be purposely designed by deploying the quantum mechanical properties displayed by the materials over length scales comparable to the interatomic spacing. New artificial solids, or “metamaterials”, would result from the assembly of quantum dots, quantum wires, nanotubes, core-shell nanoparticles, etc.,1,2 possibly coupled to specific macromolecules for particular functional or structural requirements.3 Such nanoassembled materials are expected to show tunable electronic, optical, mechanical, and thermal properties, which would preserve the individual, quantum nature of the nanostructures and, moreover, exhibit new collective effects depending on the coupling between the nanoscale objects. However, a major limitation of current technology resides in the inability to manipulate individual objects at the nanoscale, while proper three-dimensional (3D) positioning and orientation of the nanoobjects (which is controlled not only by covalent, ionic, or metallic bonds but often by complex weak interactions as well as by capping or spacer groups used in the synthesis3) is, instead, crucial to many properties. In the attempt to overcome such difficulties, one possible route for the assembly of such materials at the nanoscale draws inspiration from the self-assembly properties of the biological systems. Such properties can be used to select materials from * Corresponding author. Phone: +39 06 3048 3609. Fax: +39 06 3048 3613. E-mail: [email protected]. † Ylichron Srl, c/o ENEA Casaccia Research Centre. ‡ Dipartimento Tecnologie Fisiche e Nuovi Materiali, ENEA Casaccia Research Centre. § Universite´ de Lille I.

specific biological systems; for example, phage display methods4 allow one to extract, from a wide variety of peptide sequences (called a “library”), those which can establish a specific interaction with particular nanoscale structures. It has been shown in recent years that some peptides selectively “selfassemble” onto well-defined surfaces of semiconductors and metals.5-7 Such specific interaction could be of paramount importance for scientific and technological research, implying the exploitation of the “smart” recognition and catalytic properties of biomolecular systems coupled to electronic devices. In principle, it could also contribute to fabrication of bioelectronic devices by easy and inexpensive assembly of suitably designed, peptide based, electronically active components onto a substrate patterned with different peptide-specific semiconducting materials. However, up to now, the experiments have provided only very sparse information about the binding specificity of given peptides and nucleic acid sequences to a particular inorganic surface or nanostructure. The detailed interaction mechanisms are unknown, and the assessment remains largely empirical. For example, sequences containing histidine (His) and tryptophan (Trp) amino acids appear to exhibit a considerable binding affinity for graphitic carbon surfaces and, in particular, for single-wall carbon nanotubes (SWCNTs).7 Moreover, it has been empirically observed that charged amino acids show a stronger adhesion than nonpolar ones on insulator and semiconductor surfaces,6 while the interaction with metal surfaces is generally weaker, independently on the amount of charge transfer. On the other hand, a meaningful dependence on the sequence length has also been observed, for example, in the case of poly-leucine, suggesting an important role of the 3D folding in establishing the adhesion character.5 In the case of nucleobases (guanine, adenine, thymine, cytosine), mixed Hartree-Fock and empirical calculations8

10.1021/jp903652v CCC: $40.75  2009 American Chemical Society Published on Web 08/12/2009

12106

J. Phys. Chem. B, Vol. 113, No. 35, 2009

showed that the binding energy to (5,5) single-wall carbon nanotubes (without including water) decreases in the order G > A > T > C. A similar result was found by DFT calculations on (5,0) nanotubes.9 However, inclusion of the solvent (free) energy8 changes the ordering to G > T > A > C, in agreement with experiments of isothermal microtitration calorimetry reported in the same work.8 The situation is ambiguous, since different calculations of single-base adsorption on graphene10 supported instead the G > A > T > C ordering; this was also said to compare well to different experimental results.11 The above, partly conflicting, observations point out the complex interplay existing in such hybrid biological-inorganic systems, which are stabilized by a delicate equilibrium between very weak forces: van der Waals interactions with the underlying surface (which are badly treated theoretically, even by the most accurate ab initio methods), hydrogen bonds between adjacent molecules (in the monolayers) and within each macromolecule, and solvent effects, to name the most important ones. In this work, we used a combination of stereochemical methods, statistical mechanics, and molecular dynamics (MD) simulations, to gain insight into the atomic-scale mechanisms underlying the specific binding affinity of given peptide sequences to nanostructured carbon surfaces. More specifically, we chose a short peptide which demonstrated a high affinity for SWCNTs7 in phage display experiments, and calculated its binding affinity for (1) an ideally flat graphene sheet (FGS) and (2) a SWCNT. The comparison aims at understanding the relative roles of the hexagonal carbon structure, on one side, and of the 3D structure (flat or cylindrical), on the other, on the peptide-surface interaction. Oligopeptides of the size under consideration (12 amino acids) are systems holding a particular interest: they are long enough to assume a number of 3D conformations limited by their intramolecular interactions, and at the same time, they are short enough to experience important modifications of their structure due to the interaction with other molecular surfaces. Furthermore, the length of these peptides is short enough to make computational simulations feasible with reasonably limited effort. To summarize, we performed initial folding of the oligopeptide in water. Subsequently, a rigid docking on either the flat graphene sheet or the cylindrical SWCNT surface was studied. Structural and chemical properties of the peptide-carbon surface systems were accurately characterized, trying to elucidate the factors leading to the adhesion. Finally, an annealing cycle by MD in water was carried out, starting from the best peptide-carbon surface configuration obtained from the rigid docking, in order to assess the possible variation of the complex, thereby following what may be called a “dynamical” docking procedure. The paper is organized as follows. In section 2, the materials and the methods employed in this study are briefly described (more details are provided in the Supporting Information). Section 3 contains the results of the computer simulation of the initial folding and rigid docking of the peptide to each of the two carbon nanostructures. Section 4 describes the finitetemperature molecular dynamics simulations in water of the best peptide-carbon complexes resulting from section 3. Finally, in section 5, some conclusions from our work are drawn. 2. Computational Methods The amino acid sequence of the peptide, HWSAWWIRSNQS, was taken from the experimental work of Wang et al.7 The molecule, identified as “B3” in the cited work, has been selected as one of the best candidates, since it displayed a very large binding affinity to SWCNTs.

Gianese et al. In order to investigate at the atomic level the interaction mechanisms between the chosen peptide and the carbon surfaces (in our case, a graphene sheet and a SWCNT), we conceived a computer simulation procedure consisting of three stages: (1) modeling of the initial peptide folding through MD simulations in water, and energy minimization calculations, starting from a totally unfolded molecular conformation; (2) rigid docking between the folded peptide structure from stage 1 and carbon surfaces, by means of a genetic algorithm; (3) MD relaxation in water of the best peptide-carbon surface complexes resulting from stage 2, to adjust the geometrical matching and to further optimize the binding energy. The latter stage can be considered as a sort of “dynamical docking”. The three stages were carried out using a host of different computer simulation codes, thoroughly described in the Supporting Information. In particular, all of the MD simulations were realized with the GROMACS code12 with the OPLS-AA force field13 and the simple point charge (SPC) model14 for the water solvent molecules. The free energy of folding ∆Guff of a protein in solution (where u and f indicate the unfolded and folded states, respectively) results from the contribution of diverse types of interactions, including hydrogen bonds, electrostatic interactions, hydrophobic effect, and packing density. However, it is generally accepted that the solvent hydrophobic effect plays a predominant role, so as to be considered the main driving force for folding and binding.15 For this reason, in order to estimate the stability of the peptide and the free energy of binding for peptide-carbon surface complexes, we have considered, as a first guidance criterion, the solvation contribution to the free energy. Full details about this calculation procedure are also found in the Supporting Information. 3. Initial Folding and Rigid Docking 3.1. Folding of the Peptide. The root mean square deviation (rmsd) of the atoms in the peptide backbone from the initial structure during the three folding MD simulations is reported in Figure 1A. The trends show that only one out of three simulations was able to reach a stationary rmsd (black straight line), whereas the other two displayed increasing or unsettled deviation values all along the simulation. The simulation attaining a stable rmsd with respect to the starting structure takes about 5 ns to reach the equilibrium, after which the peptide seems to keep a stable conformation. This result is also supported by the clustering analysis (Figure 1B), based on the comparison of 3D structures extracted from the MD trajectory. The clustering of the backbone atoms of the peptide as a function of the simulation time (using a time step of 200 ps) was calculated by the GROMOS method.16 Such a method counts the number of neighbors for each structure using a rmsd cutoff distance; then, it assigns to a “cluster” the structure with the largest number of neighbors along with all its neighbors. Once assigned to a cluster, the elements are discarded from the pool of structures. The procedure is repeated until the exhaustion of all of the structures. A rmsd cutoff distance e0.10 nm was used as the criterion of similarity to cluster the structures. The number of clusters gives information about the level of population of the conformational space during the MD simulations.17,18 The number of clusters in the simulation achieving a stable rmsd (Figure 1B, sim 1) increases during the first 4.8 ns, then it slightly oscillates for 800 ps, and from 5.6 ns onward it keeps a steady value, until the end of the simulation. Such results

Interaction between Polypeptides and Carbon Surfaces

J. Phys. Chem. B, Vol. 113, No. 35, 2009 12107

Figure 1. Variation as a function of time in folding MD simulations of (A) rmsd of the backbone atoms with respect to the initial configuration and (B) cluster number.

demonstrate that the peptide stops the random folding and assumes a stable structure. The evolution of the solvation free energy of folding (∆Gs) as a function of time (Figure S3.1A, Supporting Information) agrees with the above observations. The folding energy varies irregularly, without a clear trend, within the first 1 ns of simulation. Subsequently, it progressively decreases in the course of the next 3 ns. During the last 5 ns, corresponding to the stable set of the trajectory, ∆Gs fluctuates around a mean value = -2.51 kcal/mol (-0.21 kcal/mol per residue). The lower the ∆Gs, the higher the stability in solution. Therefore, negative values of ∆Gs indicate a peptide with a more stable conformation. The trends of other structural features also support the conclusion on the stability of the folded conformation. Hydrogen bonds and hydrophobic interactions are recognized as significant factors in determining protein stability.19,20 The latter, in particular, is regarded as the major force for polypeptide folding.15 The values of both of these properties increase during the converging simulation, attaining stationary values roughly at the same time observed for the clustering and rmsd. The number of H-bonds slowly increases within the first half of the simulation time; subsequently, it oscillates around a mean value of 7. Similarly, the apolar contact surface (Figure S3.1B, Supporting Information) rises during the first 5 ns of simulation and then fluctuates around a mean value of about 298 Å2, that is, approximately 3 times larger than the initial contact area (100 Å2). 3.2. Optimized Peptide Structure. The energy of the 3D structure of the peptide obtained through folding simulations was further minimized (Figure S3.2, Supporting Information), and its stereochemical quality was evaluated (Table S2.1, Supporting Information). The Ramachandran plot and goodness factor (G-factor) are found to agree with those typical of protein crystallographic structures of good quality. The solvation free energy of folding of the final model, with respect to the completely unfolded state (in the fully extended conformation), is ∆Gs ) -3.32 kcal/mol (-0.28 kcal/mol per residue). The negative value of ∆Gs indicates that the folded state is preferred over the unfolded state, and that the solvent stabilizes the folded conformation. Although the presently folded structure could be substantially different from the one that will optimize the contact with the carbon surfaces, e.g., because of the increased interaction with aromatic residues, the relative rapidity of folding in a solvent is already a good measure of the stability of the chosen peptide. This is likely connected with

the degree of interaction of the side chains, as demonstrated by a specific analysis (not reported here for brevity).21 The peptide structure does not show any “repetitive” secondary structural element (e.g., R-helices and β-structures). A basic turn pattern involving three consecutive residues (from 5 to 7) is observed in the middle of the molecule and corresponds to the region with highest curvature. The molecule displays a total of six H-bonds (Table S2.2, Supporting Information) when a donor-acceptor distance threshold e3.5 Å is used. Atoms of amino acid His1 are involved in three out of six H-bonds; hence, this residue could play a central role in determining and stabilizing the peptide structure. The apolar contact surface spreads to an area of 299 Å2, almost identical to the mean value obtained from the folding simulation. 3.3. Docking Statistics. Molecular docking can explore the fitting of molecules together with a favorable configuration and energy to form a complex system. The docked complexes of the peptide with graphene and SWCNT were selected according to the criteria of best interaction energy combined with statistical information corresponding to the quality of the geometrical matching. The free energy of binding as estimated by AUTODOCK includes van der Waals, H-bonding interactions, desolvation, and electrostatic energies. Notably, the sum of the latter plus the H-bond contribution is vanishing, because of the all-carbon composition of the surfaces to which the peptide binds. The geometrical matching quality is taken to be proportional to the population of the statistical clusters among which the docking results are arranged. The majority of the results concerning peptide-graphene docking (Table 1) are distributed in only two clusters (with 18 and 7 elements), the remaining being dispersed in seven clusters scarcely populated. The lowest binding energy within the two best clusters is the same (-10.8 kcal/mol). On the other hand, the results concerning the peptide-SWCNT docking (Table 1) fall in a single cluster of 23 elements, while the remaining are distributed among three clusters, two with 8 elements and one with a single structure. The best docked complex within the largest cluster has a binding energy of -13.0 kcal/mol. The complex with lowest binding energy in each of the highest populated clusters was selected as the starting system for the following MD simulations in water, for further energy minimization and geometrical optimization.

12108

J. Phys. Chem. B, Vol. 113, No. 35, 2009

Gianese et al.

TABLE 1: Results of Peptide-Carbon Surface Docking cluster rank

lowest binding energy (kcal/mol)

mean binding energy (kcal/mol)

number in cluster

1 2 3 4 5 6 7 8 9

Peptide-FGS System -10.79 -10.66 -10.77 -10.64 -10.62 -10.54 -10.61 -10.44 -10.59 -10.57 -10.55 -10.27 -10.53 -10.34 -10.52 -10.52 -10.10 -10.08

18 7 2 3 2 3 2 1 2

1 2 3 4

Peptide-SWCNT -13.01 -12.07 -10.69 -9.97

23 8 1 8

System -12.70 -11.82 -10.69 -9.87

4. Molecular Dynamics Simulations 4.1. Peptide-Graphene. The rmsd as a function of time of the backbone atoms relative to the initial configuration of the peptide bound to the graphene for two different simulations (GS_sim1 and GS_sim2) is reported in Figure S3.3A,B (Supporting Information). In GS_sim1, the rmsd displays an initial slight increase to ∼0.10 nm after 100 ps, thereafter remaining fairly stable. Thus, the trajectory immediately achieves the convergence, whereas the peptide preserves its structure, as resulting from the clustering analysis (Figure S3.3C, Supporting Information) which detects only one cluster at each simulation time step. Similarly, the rmsd in GS_sim2 rises to about 0.10 nm within the first 500 ps; afterward, it fluctuates around this value for most of the trajectory, with the exception of a transitional deviation in the time interval, roughly, from 4.2 to 4.8 ns. The GS_sim2 clustering trend (Figure S3.3D, Supporting Information) allows explanation of the observed temporary deviation of the rmsd as an attempt to explore a new molecular conformation. Indeed, the number of clusters, after a constant one-cluster stage during the first 4.5 ns, jumps up to a value of 3 in correspondence with the major peak in the rmsd curve. Finally, the total cluster number decreases and stabilizes to a value of 2 within the last 0.5 ns of simulation. The cluster analysis on the whole trajectory reveals that one of the two detected final clusters contains the snapshots from 1 to 4.4 ns and from 4.8 to 6.0 ns, indicating that, after the rmsd peak, the peptide returns to the initial configuration. The variation of the solvation free energy of binding as a function of time in GS_sim1 (Figure S3.3E, Supporting Information, black line) agrees with the path followed by the trajectory. The peptide-graphene complex displays only minor variations of the binding energy during the whole simulation period, as expected from the trends of cluster number and rmsd. The free energy fluctuates around a mean value (-7.21 kcal/ mol) very close to that of the starting system (-7.32 kcal/mol). The constancy of the solvation binding energy during the simulation reflects the largely unchanged surface area of the peptide in contact with the graphene (Figure S3.3F, Supporting Information, black line). Indeed, the mean value of the contact surface (352 Å2) stays close to that of the starting configuration (346 Å2). The solvation free energy of binding during GS_sim2 (Figure S3.3E, Supporting Information, black dashed line) follows a trend similar to that of GS_sim1 within the first 2 ns, with a mean value of -7.28 kcal/mol. Subsequently, the free energy rapidly decreases, maintaining a lower value along the rest of the trajectory, with an average of -9.26 kcal/mol. This

value is 1.94 kcal/mol lower than that of the starting system. The significant change in the binding free energy, corresponding to a stronger peptide-graphene interaction, can be explained by observing the variation of the contact surface area between the two molecules (Figure S3.3F, Supporting Information, black dashed line). As for the binding energy, the contact surface area, during the first 2 ns of the run, displays a curve similar to that of the GS_sim1, oscillating around a mean value of 350 Å2. Then, the contact surface increases to a mean value (403 Å2), about 16% larger than the contact area of the starting configuration. Since the folded structure of the peptide remains substantially unchanged during the GS_sim2, it can be asserted that the higher contact area of the peptide-graphene association is achieved through a rearrangement of the amino acid side chains. Most importantly, the contribution of the single residues to the binding free energy as a function of time (Figure 2A,B) displays the same trend of the total energy. However, the differences among the various residues reveal that the hydrophobic amino acids Ile and Trp are mainly responsible for the interaction with the carbon surfaces, since their contribution to the binding energy is the only negative one in GS_sim1, and always the most negative in GS_sim2. In particular, the effect of the three Trp is overriding. Moreover, in the case of GS_sim2, the nearly 2 kcal/mol decrease in the solvation binding energy after 2 ns is made possible by the intervention of the negative contribution, absent in GS_sim1, of the residues His1, Ala4, Gln11, Ser3, and Ser9. 4.2. Peptide-SWCNT. The rmsd as a function of time of the peptide backbone atoms binding to the SWCNT (NT_sim1 and NT_sim2) is reported in Figure S3.4A,B (Supporting Information). In both trajectories, the peptide shows a rmsd with an initial small increase to about 0.05 nm. Subsequently, this value is preserved with moderate fluctuations up to times of 3.2 and 2.7 ns in NT_sim1 and NT_sim2, respectively. Similarly to the peptide-graphene system, the cluster number in this portion of the trajectory (Figure S3.4C,D Supporting Information) is equal to 1. Thus, we can consider the 3D structure of the molecule to be preserved. Later on, the rmsd rises rapidly, stabilizing at a value of ∼0.20 nm for the rest of NT_sim1, and at ∼0.16 nm for NT_sim2. The rmsd trend reflects changes in the cluster number, which increases from 1 to 4 in NT_sim1, and from 1 to 5 in NT_sim2, within the same time interval. Apparently, in the presence of the SWCNT, the peptide structure undergoes critical modifications leading to a molecularscale refolding, which more easily attains the convergence. The transition phase lasts for a very short period; therefore, the conformational change occurs very rapidly. This conclusion is strengthened by the steep slope, roughly at 3.8 ns, observed in the rmsd variation plot of NT_sim1, which shows a sigmoidal profile. Although snapshots extracted from the trajectories of the peptide-SWCNT system merge into a single cluster during the first 3.0 ns for the two simulations, the solvation free energy of binding (Figure S3.4E, Supporting Information) shows large and significant fluctuations within the same period. The peptide exhibits energy instability in the presence of the SWCNT, and seems to explore a new conformational space. After 3.8 ns in NT_sim1 and 3.0 ns in NT_sim2, the binding free energy decreases steadily, and for t > 4.0 ns, when the cluster number has reached the plateau in both trajectories, it is stable around mean values of -9.01 kcal/mol (NT_sim1) and -8.65 kcal/ mol (NT_sim2). Such free energies are lower than that of the starting system (-8.41 kcal/mol), albeit the differences are small

Interaction between Polypeptides and Carbon Surfaces

J. Phys. Chem. B, Vol. 113, No. 35, 2009 12109

Figure 2. Variation as a function of time of the solvation free energy of association for each residue in peptide-graphene MD simulations (A, GS_sim1; B, GS_sim2) and peptide-nanotube MD simulations (C, NT_sim1; D, NT_sim2).

when compared with the energy shift observed in the GS_sim2 of the peptide-graphene complex. It may be objected that such (very close) folded states could not yet be the stationary configuration, since longer peptides and proteins often require much longer times to attain a lowest free energy state. In our case, the stability of the “native” folding for more than 6 ns, the relatively good stability of the peptide-carbon complex over several nanoseconds, the short contour length of the peptide, and the fact that the average number of about one H-bond per residue, after refolding at the carbon surface, is comparable or larger than what was found in the most stable natural proteins (thermophiles) allow us to deduce that such folded configurations are likely to represent the best folded state. In any case, the new folding allows the peptide to bind the SWCNT more strongly than in the starting conformation. However, we cannot rule out the possibility that the peptide could assume both of the above conformations, alternatively. In fact, MD simulations performed using two 12residue peptides with different sequences but a high affinity for SWCNTs22 showed that trajectories yielding complexes with strong binding may support different arrangements of His and Trp residues on the carbon tube surface. This observation suggests that the specificity of these peptides may depend on their capability to support multiple strong binding conformations. Evidence that the final conformation could be the most likely

is the higher stability of its binding free energy and the fact that the refolded peptide structure is the same for the two simulations. Indeed, the rmsd calculated for the backbone atoms of the superposed 3D structures obtained by minimization of the last snapshot of the two trajectories gives a value of 1.04 Å, indicating that the molecules share the same folding state. As for the peptide-graphene complex, also in this case, the variation of the solvation free energy of binding derives from the trend of the contact surface area of the peptide with the carbon surface (Figure S3.4F, Supporting Information). Indeed, the binding surface area shifts from similar average values of 310 Å2 for NT_sim1 and 317 Å2 for NT_sim2, during the onecluster stage, to average values of 417 and 394 Å2 in NT_sim1 and NT_sim2, respectively, after the structural change of the peptide. Such larger contact surface and stronger binding to the SWCNT are not surprising, when we consider that the peptide was identified by the phage display technique specifically for its selective affinity to this carbon surface.7 The variation of the interaction free energy expressed for the single residues as a function of simulation time (Figure 2C,D) reveals a general improvement in the contribution to the binding for all of the amino acids, after the transition to a new molecular folding. As already observed for the peptide-graphene complex, the three tryptophans and Ile7 play the major role. Moreover, the rearrangement of the peptide chain upon refolding enables

12110

J. Phys. Chem. B, Vol. 113, No. 35, 2009

Gianese et al.

TABLE 2: Structural and Thermodynamic Properties of the Peptide before- and after-MD Simulations of Interaction system peptide(FGS): pre GS_sim1 GS_sim2 peptide(SWCNT): pre NT_sim1 NT_sim2

H- ACSa ∆Gs BSb ∆∆Gs rmsd bonds (Å2) (kcal/mol) (Å2) (kcal/mol) (Å) 6 6 12

299.0 299.2 305.9

-3.32 -3.09 -2.71

346.3 349.9 410.2

-7.32 -7.07 -9.34

0.76 1.00

6 9 11

299.0 363.5 283.0

-3.32 -3.15 -2.14

298.7 433.2 396.3

-8.41 -9.02 -8.96

2.07 1.89

a ACS, apolar contact surface; intramolecular property reffered to the peptide. b BS, binding surface; surface area of the peptide in contact with the carbon surface.

other residues to give a negative contribution to the binding energy, as His1 and Ala4 in NT_sim2, and furthermore Ser3, Ser9, and, marginally, Gln11, in NT_sim1. The predominant role of the Trp residues agrees with the finding that strong binding sequences usually show a motif rich in Trp, besides His, at specific locations.7 4.3. Comparative Analyses between Relaxed and Unrelaxed Structural Models. The energies of the peptide-carbon surface complexes were minimized at zero temperature by a steepest descent algorithm, after the MD simulations of interaction, and the stereochemical quality of the final peptide structures was evaluated (Table S2.3, Supporting Information). Structural and thermodynamic properties evaluated in this work are summarized in Table 2, and refer to the minimized model obtained from the folding MD simulations and to the minimized carbon surface-bound peptide systems. The comparison of the final models of the peptide bound to the graphene against the complex used as the starting configuration in the “dynamic docking” does not reveal significant differences pertaining to the molecular folding, as suggested by the results discussed above. The preservation of the molecular conformation is demonstrated by the calculation of the rmsd of the backbone atoms between the peptide structures resulting from the two MD simulations on graphene (Figure 3A,C) and the 3D model derived from the folding simulations (Figure S3.2A, Supporting Information), giving values of 0.76 and 1.00 Å for GS_sim1 and GS_sim2, respectively. This result bears upon the intramolecular apolar contact surface and ∆Gs, which remain nearly unchanged. Instead, significant differences arise when more subtle structural features are considered. At the end of GS_sim1, the extension of the peptide surface in contact with the graphene, as well as ∆∆Gs, are very similar. Also the number of H-bonds is the same, albeit many interacting atoms are different (Table 3). Nevertheless, the structural role of His1 is maintained, since this residue establishes two H-bonds with Ser12 at the C-terminal of the polypeptide chain (Figure 3B), probably strengthening the molecular folding. At the end of GS_sim2, on the contrary, the peptide surface interacting with the graphene is about 60 Å2 larger, and as a consequence, ∆∆Gs is substantially more negative (-9.34 kcal/mol, with a consistent drop in energy by 2.02 kcal/mol), meaning a higher binding affinity for the carbon surface. Moreover, the intramolecular structural stability is also enhanced, since the number of H-bonds doubles with respect to the starting system (Table 3), likely determining a higher rigidity of the polypeptide chain. As for GS_sim1, the central role of His 1 in stabilizing the peptide folding is preserved, since it forms a total of three H-bonds by interacting with Ser12 and Ser9 (Figure 3D).

Figure 3. Structures of the peptide bound to graphene after MD simulations of interaction. The molecules are represented as yellow cartoon and yellow sticks, at the end of (A) simulation 1 and (C) simulation 2. Intramolecular interactions are shown as yellow sticks of the main chain and of side chains involved in H-bonds, which are depicted as orange dashed lines, at the end of (B) simulation 1 and (D) simulation 2. Oxygen atoms are colored red, nitrogen atoms blue, and hydrogens white. This figure was rendered by using PyMOL.22

TABLE 3: Hydrogen Bonds of the Energy Minimized after-MD Structures of the FGS-Bound Peptide distance (Å) donora

acceptora

regionb

GS_sim1

GS_sim2

1-HIS(N) 12-SER(N) 1-HIS(NE2) 3-SER(N) 3-SER(OG) 3-SER(OG) 9-SER(OG) 3-SER(OG) 7-ILE(N) 8-ARG(N) 9-SER(OG) 9-SER(OG) 11-GLN(N) 12-SER(OG)

12-SER(O) 1-HIS(O) 9-SER(OG) 11-GLN(OE1) 9-SER(O) 9-SER(OG) 3-SER(OG) 11-GLN(OE1) 4-ALA(O) 4-ALA(O) 7-ILE(O) 10-ASN(O) 10-ASN(OD1) 11-GLN(O)

MM MM SS MS SM SS SS SS MM MM SM SM MS SM

2.61 2.96

2.57 3.04 3.21 2.97 3.44 3.38 3.38 2.99 2.86

2.72

3.00 3.11 2.96

3.14 2.86 3.17

a In parentheses is indicated the type of the atoms involved in H-bonds as donor and acceptor. b H-bond type: M, main chain; S, side chain.

The final model structures of the peptide bound to the SWCNT (Figure 4A,C) display consistent differences compared to the starting complex (rmsd ) 2.07 Å for NT_sim1 and 1.89 Å for NT_sim2), due to the change of folding occurring during the MD simulations. The peptide surface in contact with SWCNT largely increases, determining the fall in ∆∆Gs from -8.41 to -9.02 kcal/mol in NT_sim1 and to -8.96 kcal/mol in NT_sim2. It can be supposed that the peptide is able to refold only in the presence of SWCNTs, since, otherwise, the larger apolar surface in contact with the carbon surface would be exposed to the solvent, destabilizing the molecule. After NT_sim1, also the intramolecular apolar contact surface of the peptide increases, probably because of the side-chain orientation of the three Trp residues which, differently from the initial conformation, are all grouped together, close to one another. However, no significant variation in the solvation energy of folding is observed, probably because of the counterbalancing effect of the extension of the polar contact surface. At the end,

Interaction between Polypeptides and Carbon Surfaces

J. Phys. Chem. B, Vol. 113, No. 35, 2009 12111

Figure 4. Structures of the peptide bound to SWCNT after MD simulations of interaction. The molecules are represented as yellow cartoon and yellow sticks, at the end of (A) simulation 1 and (C) simulation 2. Intramolecular interactions are shown as yellow sticks of the main chain and of side chains involved in H-bonds, which are depicted as orange dashed lines, at the end of (B) simulation 1 and (D) simulation 2. Oxygen atoms are colored red, nitrogen atoms blue, and hydrogens white. This figure was rendered by using PyMOL.22

TABLE 4: Hydrogen Bonds of the Energy Minimized after-MD Structures of the SWCNT-Bound Peptide distance_(Å) donora

acceptora

regionb

NT_sim1

NT_sim2

11-GLN(NE2) 1-HIS(NE2) 9-SER(OG) 3-SER(OG) 3-SER(OG) 9-SER(OG) 7-ILE(N) 8-ARG(N) 8-ARG(NE) 9-SER(OG) 11-GLN(N) 12_SER(N) 12-SER(OG) 12-SER(OG)

1-HIS(O) 9-SER(OG) 1-HIS(NE2) 9-SER(O) 9-SER(OG) 3-SER(OG) 4-ALA(O) 4-ALA(O) 5-TRP(O) 11-GLN(OE1) 9-SER(O) 10-ASN(OD1) 10-ASN(OD1) 11-GLN(O)

SM SS SS SM SS SS MM MM SM SS MM MS SS SM

2.84

2.76 2.86 2.86 2.86 3.42 3.42 2.94 2.91 2.96

2.78 2.78 2.95 3.16 2.98 3.27 3.04

3.16 2.72 2.75v

a In parentheses is indicated the type of the atoms involved in H-bonds as donor and acceptor. b H-bond type: M, main chain; S, side chain.

the refolded peptide displays more H-bonds than the structure in the starting complex, reaching a total of 9 and 11 interactions after NT_sim1 and NT_sim2, respectively (Table 4). The high amount of H-bonds enhances the molecular structural stability, and could be one of the main reasons for the strong binding affinity to SWCNTs of this peptide. A central structural role in the refolded conformation is played by residue Ser9, which is involved in four (Figure 4B) and five (Figure 4D) H-bonds at the end of NT_sim1 and NT_sim2, respectively. The 3D structures of the final energy minimized models of the peptide-carbon surface complexes are depicted in Figure 5. It is interesting to note that the simulation ending with the strongest binding energy for the peptide-graphene system (GS_sim2) shows a geometry in which the Trp6 residue adopts an orientation parallel to the graphene surface (Figure 5B). This finding agrees with a previous work,22 although in that case the surface was a carbon nanotube. Therefore, the Trp6 orientation could explain the higher binding affinity of the complex in

Figure 5. Final minimized models of the binding between the “B3” peptide and (A) the graphene at the end of GS_sim1, (B) the graphene at the end of GS_sim2, (C) the SWCNT at the end of NT_sim1, (D) the SWCNT at the end of NT_sim2. The peptide is represented as blue (A), green (B), red (C), and yellow (D) sticks, while the all-carbon surfaces are represented as gray sticks. This figure was rendered by using PyMOL.22

GS_sim2 compared with GS_sim1 (Figure 5A). The His1 residue, instead, remains near and parallel to the surface in both of the trajectories. As it can be expected, the hydrophilic residues extend into the solvent rather than interact with the hydrophobic surface, except for the Arg8 whose side chain, completely extended, is adjacent and parallel to the graphene sheet. In our peptide-SWCNT complex, only the Trp2 in NT_sim2 assumes an orientation nearly parallel to the carbon surface (Figure 5C), while, generally, the tryptophans close to the surface are preferably perpendicular. It is worth noting that, in this case, the most distant residues from the carbon surface are those that more strongly interact with the graphene (i.e., His1, Ile7, and Trp6), whereas all the other amino acids seem to extend as much as possible their side chain, to best adapt to the carbon tube. 5. Conclusions We performed a comparative study of the adsorption of a dodecapeptide on two different carbon surfaces: a flat graphene sheet and a curved (0,15) single-wall nanotube. The peptide, selected from an experimental database as the one with the highest affinity for carbon surfaces, was initially folded in water, and then rigidly docked on the surfaces with a genetic algorithm. The best candidate structures from this rigid docking procedure were selected, and used as input in subsequent MD calculations in water at finite temperature. Free energies of folding and adhesion, and the quality of surface binding were determined, based on a combination of solvation energy, formation of hydrogen bonds, and amount of the total contact surface between peptide and carbon. In the case of the interaction with the graphene surface, we found that the peptide structure is not substantially changed upon adhesion. After MD equilibration at finite temperature, the total free energy of adhesion, for the strongest binding trajectory (GS_sim2), is lowered by -1.94 kcal/mol, and the contact surface area with the graphene is increased by about 18%, compared to the rigid docking. The main things responsible for the interaction with the carbon surface are identified as the amino acids Ile and Trp, since their free energy contribution is constantly the most negative, when not the only one. In particular, the effect of the three tryptophans is overriding. In the case of the interaction with the curved carbon nanotube, we found the peptide to undergo critical modifications after

12112

J. Phys. Chem. B, Vol. 113, No. 35, 2009

finite-temperature molecular dynamics, with respect to the rigid docking configuration. While the binding free energy differences, -0.61 to -0.35 kcal/mol, are smaller than those for the best graphene trajectory, the peptide shows a more flexible accommodation on the curved surface of the nanotube, with more than 30% increase of the contact surface for the refolded structure (up to 45% in NT_sim1). Several competing adsorption minima exist, corresponding to different arrangements of the His and Trp residues. Such flexibility suggests that the surface specificity of the peptide may depend on its capability to support multiple strongly bound configurations. Also, in this case, the tryptophans play a predominant role in adhesion, confirming the experimental observations.7 However, the refolding into different configurations enables also other residues to give a negative contribution to the free energy of association. Also, the number of hydrogen bonds increases for all of the refolded structures, compared to the rigidly docked structure. The central role of Ser residues in promoting H-bonds was highlighted. The observation that the free energies involved in the peptide interaction with the two carbon surfaces analyzed are similar suggests the opportunity to design unspecific peptides endowed with strong binding affinity for a wide range of carbon surfaces. Future theoretical work should be first directed at improving the description of the dispersion forces, namely, van der Waals interactions. Such a contribution is of paramount importance to the folding free energy and to the surface adhesion, but it is still difficult to derive in the framework of ab initio electron density functional theory. In fact, the electronic exchange and correlation effects, which are the very origin of VdW forces, are only approximated by the generalized-gradient or localdensity theoretical models.24 While recent improvements have appeared in the literature,25,26 such a contribution is still a major point of uncertainty in the theoretical models of molecular adsorption. Second, finite size effects in graphene should be the subject of further investigation. Binding on the edges of the graphene sheet can be relevant in the interaction with biomolecules. The graphite free edges represent structural defects, at which a large charge transfer might take place with the molecular species. The interaction is expected to be partly ionic and partly covalent; therefore, a combination of electronic structure calculations and empirical molecular dynamics must be envisaged in this case. Acknowledgment. The authors acknowledge the CRESCO Project High-Performance Computing Facility, at the ENEA Portici Research Center, where most of the (very timeconsuming) MD simulations have been performed.

Gianese et al. Supporting Information Available: Full details about methods and computer codes, supplementary references, tables containing information about the stereochemical quality and H-bonds of the peptide after MD folding simulations, and the peptide stereochemical quality after MD simulations of adhesion on carbon surfaces, figures containing details about the time evolution of cluster analysis, rmsd, solvation free energy, contact area, for the simulations described in the text. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Engheta, N. Science 2007, 317, 1698. (2) Redl, F. X.; Cho, K.-S.; Murray, C. B.; O’Brien, S. Nature 2003, 423, 968. (3) Erwin, S. C.; Zu, L.; Haftel, M. I.; Efros, A. L.; Kennedy, T. A.; Norris, D. J. Nature 2005, 436, 91. (4) Merzlyak, A.; Lee, S.-W. Curr. Opin. Chem. Biol. 2006, 10, 246. (5) Willett, R. L.; Baldwin, K. W.; West, K. W.; Pfeiffer, L. N. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 7817. (6) Whaley, S. R.; English, D. S.; Hu, E. L.; Barbara, P. F.; Belcher, A. M. Nature 2000, 405, 665. (7) Wang, S.; Humphreys, E. S.; Chung, S.-Y.; Delduco, D. F.; Lustig, S. R.; Wang, H.; Parker, K. N.; Rizzo, N. W.; Subramoney, S.; Chiang, Y.-M.; Jagota, A. Nat. Mater. 2003, 2, 196. (8) Varghese, N.; Mogera, U.; Govindaraj, A.; Das, A.; Maiti, P. K.; Sood, A. K.; Rao, C. N. ChemPhysChem 2009, 10, 206. (9) Gowtham, S.; Scheicher, R. H.; Pandey, R.; Karna, S. P.; Ahuja, R. Nanotechnology 2008, 19, 125701. (10) Schwabe, T.; Grimme, S. Phys. Chem. Chem. Phys. 2007, 9, 3397. (11) Sowerby, S. J.; Cohn, C. A.; Heckl, W. M.; Holm, N. G. Proc. Natl. Acad. Sci. U.S.A. 2001, 98, 820. (12) Van Der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. J. Comput. Chem. 2005, 26, 1701. (13) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. J. Am. Chem. Soc. 1996, 118, 11225. (14) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Hermans, J. In Intermolecular forces; Pullman, B., Ed.; D. Reidel Publishing Company: Dordrecht, The Netherlands, 1981; p 331. (15) Honig, B.; Yang, A. S. AdV. Protein Chem. 1995, 46, 27. (16) Daura, X.; Gademann, K.; Jaun, B.; Seebach, D.; van Gunsteren, W. F.; Mark, A. E. Angew. Chem., Int. Ed. 1999, 38, 236. (17) Smith, L. J.; Daura, X.; van Gunsteren, W. F. Proteins 2002, 48, 487. (18) Brigo, A.; Lee, K. W.; Iurcu Mustata, G.; Briggs, J. M. Biophys. J. 2005, 88, 3072. (19) Dill, K. A. Biochemistry 1990, 29, 7133. (20) Pace, C. N.; Shirley, B. A.; McNutt, M.; Gajiwala, K. FASEB J. 1996, 10, 75. (21) Gehenn, K.; Stege, J.; Reed, J. Anal. Biochem. 2006, 356, 12. (22) De Miranda Toma´sio, S.; Walsh, T. R. Mol. Phys. 2007, 105, 221. (23) De Lano, W. L. The PyMOL molecular graphics system; DeLano Scientific LLC: San Carlos, CA, 2004. (24) Gao, Y.; Zeng, X. C. J. Phys.: Condens. Matter 2007, 19, 386220. (25) Rydberg, H.; Dion, M.; Jacobson, N.; Schro¨der, E.; Hyldgaard, P.; Simak, S. I.; Langreth, D. C.; Lundqvist, B. I. Phys. ReV. Lett. 2003, 91, 126402. (26) Andrade, X.; Botti, S.; Marques, M. A. L.; Rubio, A. J. Chem. Phys. 2007, 126, 184106.

JP903652V