Article pubs.acs.org/JPCB
Relationship between Conformational Dynamics and Electron Transfer in a Desolvated Peptide. Part I. Structures David Semrouni, Carine Clavaguéra,* and Gilles Ohanessian* Laboratoire des Mécanismes Réactionnels, Department of Chemistry, Ecole Polytechnique, CNRS, 91128 Palaiseau Cedex, France
Joel H. Parks* Rowland Institute at Harvard, 100 Edwin H. Land Boulevard, Cambridge, Massachusetts 02142, United States S Supporting Information *
ABSTRACT: The structures, dynamics and energetics of the protonated, derivatized peptide DyeX-(Pro)4-Arg+-Trp, where “Dye” stands for the BODIPY analogue of tetramethylrhodamine and X is a (CH2)5 linker, have been investigated using a combination of modeling approaches in order to provide a numerical framework to the interpretation of fluorescence quenching data in the gas phase. Molecular dynamics (MD) calculations using the new generation AMOEBA force field were carried out using a representative set of conformations, at eight temperatures ranging from 150 to 500 K. Force field parameters were derived from ab initio calculations for the Dye. Strong electrostatic, polarization and dispersion interactions combine to shape this charged peptide. These effects arise in particular from the electric field generated by the charge of the protonated arginine and from several hydrogen bonds that can be established between the Dye linker and the terminal Trp. This conclusion is based on both the analysis of all structures generated in the MD simulations and on an energy decomposition analysis at classical and quantum mechanical levels. Structural analysis of the simulations at the different temperatures reveals that the relatively rigid polyproline segment allows for the Dye and Trp indole side chain to adopt stacking conformations favorable to electron transfer, yielding support to a model in which it is electron transfer from tryptophan to the dye that drives fluorescence quenching. lations in the following paper (part II28) to compare the quenching rate data with the temperature dependence of the Marcus model of the electron transfer rate. These structure calculations not only allow us to identify conformations exhibiting Dye−Trp proximity and exothermic reactivity required for rapid quenching by electron transfer, but also provide a detailed description of the intramolecular interactions responsible for the occurrence of such quenching conformations. This effort to identify the basis for fluorescence lifetime quenching in biomolecules is required to understand if such measurements can yield fundamental information related to changes in properties such as conformation. In general, the probability of electron transfer in biomolecules has been shown to be intimately dependent on the conformation through the local separation of the electron donor and acceptor. Consequently, to the extent that fluorescence lifetime quenching is a consequence of electron transfer interactions, quenching data will display the signature of local changes in conformation. The capability to isolate local changes within a
1. INTRODUCTION Fluorescence-based methods have been developed to probe the structure and local conformational dynamics of trapped gasphase biomolecule ions which have been derivatized with a strongly fluorescing dye. Measurements have included quenching in Trp−cage protein,1,2 polyproline peptide sequences,3,4 and noncovalent Vancomycin−peptide complexes.5 In each case the quenching could be related to conformational changes with temperature induced by spatial fluctuations of the side chains. Experiments verified that interactions between the dye and a Trp residue were required to observe quenching of the dye fluorescence. Measurements of the unsolvated peptide sequences of DyeX(Pro)n-Arg+-Trp, where “Dye” stands for the BODIPY (borondipyrromethene class of dye) analogue of tetramethyl-rhodamine and X is a (CH2)5 linker (see Scheme 1), introduced a unique opportunity to study the quenching process, since these sequences constrain fluctuations except for the interacting dye and Trp side chain.3,4 This paper (part I) concentrates on a detailed structural analysis of the peptide sequence DyeX(Pro)4-Arg+-Trp as a function of temperature including both molecular dynamics (MD) and density functional theory (DFT) calculations. The calculated temperature dependence of conformational fluctuations are derived from MD calcu© 2013 American Chemical Society
Received: August 7, 2012 Revised: December 12, 2012 Published: January 8, 2013 1746
dx.doi.org/10.1021/jp3078375 | J. Phys. Chem. B 2013, 117, 1746−1755
The Journal of Physical Chemistry B
Article
Scheme 1. Detailed Labeling of the Sequence DyeX-(Pro)4-Arg+-Trp
underline the requirements for a force field to properly model the forces which are responsible for the formation of these conformations that are critical for electron transfer.
biomolecular conformation offers the possibility of identifying dynamic changes occurring within specific regions of secondary structure. These gas-phase methods also have the advantage that measurements can be directly related to DFT calculations and MD simulations. An understanding of the dynamics derived from intramolecular interactions in the absence of surrounding water layers provides a basis for exploring changes in these dynamics introduced in condensed phase measurements. Previous molecular dynamics calculations4 performed to identify the basis for the fluorescence quenching were consistent with an electron transfer process but did not provide an unambiguous result that confirmed the presence of electron transfer. This paper presents the analysis of MD simulations of the structural fluctuations for an unsolvated peptide sequence DyeX-(Pro)4-Arg+-Trp. These simulations are based on the new generation AMOEBA force field,5−12 which explicitly includes polarizability contributions as well as atomic multipole-based electrostatic interactions and a refined van der Waals term. This is expected to bring higher reliability to MD results and improved ability to identify the crucial energy terms, especially when strong local electrostatic, polarization and dispersion interactions are expected. This is particularly relevant for calculations of a charged, unsolvated peptide in which multiple hydrogen bonds and π-stacking interactions are likely to occur. These calculations determined the intramolecular interactions contributing to formation of the peptide conformations and the forces driving the conformational fluctuations. In particular, the MD trajectories displayed fluctuations among conformers which varied in the separation of the dye and Trp residue as a function of temperature. Since electron transfer interactions depend exponentially on this separation, this parameter serves to identify correlations between the temperature dependence of the electron transfer rate and the quenching rate. In addition, DFT calculations were performed on several structures identified in the MD simulations. Reliable quantum chemical results for systems of this level of complexity are not commonly available in the literature. The DFT interaction energy between molecular fragments was further decomposed into components following the Morokuma−Ziegler13 scheme in order to strengthen the identification of energy terms at work in a π-stacking interaction between the dye and Trp side chain. The results
2. COMPUTATIONAL METHODS 2.1. Force Field. Force fields whose electrostatic energy component is based only on point charges have been shown to reproduce only half of a set of amino acid conformations as established by accurate quantum chemical calculations, whereas the new generation force field AMOEBA reaches 80% of reproduction.7,14 For the derivatized peptide studied herein, the presence of the charged Arg side chain and aromatic rings in both Trp side chain and the BODIPY suggests a special importance of electrostatics and polarization in structures, energetics and a fortiori dynamics. The fluorescence quenching model via electron transfer4 involves the Dye and the Trp side chain. For a correct description of interactions between those two groups in the presence of a charged Arg residue, one needs to calculate polarization which is nonadditive by nature. For this purpose, we have used the AMOEBA force field, developed by Ren and Ponder.6,9−11 Polarization is calculated through mutually induced point dipoles15 damped according to a Tholé model.16 The calculation of electrical interactions is based on atomic induced dipole sites and on permanent atomic charges, dipoles and quadrupoles. In addition, dispersion is important for describing the local interactions of aromatic cycles and for peptide folding in general.17 AMOEBA uses a “buffered 14−7” functional form, instead of a classical Lennard-Jones function, to improve the modeling of van der Waals interactions.18 All AMOEBA calculations have been performed with the Tinker package version 5.0.8 2.2. Dye Parameters. AMOEBA parameters were available for peptide residues and terminal caps but not for the linker and the Dye. The linker was mainly assimilated to a lysine CH2 chain terminated by two peptidic amides. We consider the same simplification of the BODIPY as used in previous studies:4 the replacement of a B−N atom pair by two carbons. The Dye C6H4-O-CH3 group is considered as a Tyr side chain modified to be a methyl ether. The OCH3 parameters are those of the methyl ether defined in AMOEBA. Corresponding van der Waals parameters exist in AMOEBA parameter sets for proteins and organic molecules.10 The three rings of BODIPY needed more parameter adaptation. Approximation was made on 1747
dx.doi.org/10.1021/jp3078375 | J. Phys. Chem. B 2013, 117, 1746−1755
The Journal of Physical Chemistry B
Article
four peptide fragmentation schemes in order to obtain a detailed picture of its local, intramolecular interactions. The calculations were performed at the B3LYP-D/TZ2P level with the ADF 2010.02 package.25
torsion parameters in the Dye, considering them to be null because of the strong planarity preference of aromatic rings. For C, N, and H atoms involved in aromatic cycles we used van der Waals parameter values available for other aromatic cycles. In AMOEBA, fluorine is defined as a fluoride ion. van der Waals parameters were consequently taken from the Merck force field18 which uses the same potential formulation for van der Waals interactions. Atomic polarizabilities were already defined for each element and therefore taken from available AMOEBA parameter sets. For electrostatic parameter determination, quantum mechanics was used in order to be consistent with its treatment in standard AMOEBA. Multipoles were extracted via Stone’s distributed multipole algorithm (DMA)19 as implemented in the GDMA 2.2 software.20 In order to remain consistent with existing AMOEBA multipoles, DMA was applied to a wave function calculated with the MP2 method associated with the cc-pVTZ basis set, using the initial algorithm. A computation at this level was not performed on the whole Dye but rather on the three rings only, replacing methyl and other substituents, except for fluorine, by hydrogen atoms. As for some other chemical groups, methyls linked to an aromatic cycle are already defined in AMOEBA and we have considered their parameters to be transferable. Charge, dipole and quadrupole parameters for the atoms of the three rings of the Dye are derived from DMA calculations performed with a H radius of 0.31 Å. After reintroducing the three rings in the whole system, the total charge has been set to its correct integer value by small modifications of atom partial charges. 2.3. Starting Structures and MD Calculations. The AMOEBA force field including new atom types and parameters for BODIPY as described above was used to calculate the potential energy during MD trajectories. These simulations have been performed at several temperatures using a Berendsen thermostat (coupling time: 0.1 ps). Trajectories were propagated by a Beeman integrator with a 1 fs time step. Structures were saved every 0.1 ps. Five optimized structures were initially selected from previous simulated annealing results:4 str2, str4, str7, str8, and str9 and used in this work as initial structures. Among the low energy structures previously identified,4 these five were taken to ensure diversity in Dye− Arg and Trp−Arg relative orientations. For each starting structure, 5 ns trajectories were performed at 8 temperatures ranging from 150 to 500 K with a uniform step of 50 K. Since fluorescence lifetime is ca. 10 ns at 150 K and progressively drops down to ca. 2 ns at 465 K, a sampling time of 5 ns repeated for 5 different starting structures at each temperature appears to be adequate to capture the essentials of the dynamics underlying the quenching process. 2.4. Density Functional Theory Calculations. In order to gain more insight into the energy components that are crucial in shaping conformations where electron transfer is expected to be favorable, a conformation with stacked Dye and Trp rings was selected for refined quantum chemistry calculations at DFT level. We used the B3LYP functional augmented with an empirical term for dispersion21,22 with the SVP basis set for full geometry optimization, using the Turbomole software package.23,24 As for MD simulation, a B−N atom pair was replaced by two carbons. The structure was partitioned into components and the binding energy between these components was expressed as the sum of repulsion, electrostatic, orbital interaction and dispersion terms using a Morokuma−Ziegler type decomposition.13 This was applied to
3. PEPTIDE MOLECULAR DYNAMICS 3.1. Atom Labeling. Atomic ordering starts from the BODIPY which is at the N-termination of the peptide. Oxygen atoms of the linker are numbered from the dye to the peptide. As shown in Scheme 1, hydrogen bonds are designated by “H label”_”group name” - “heteroatom”_”group name”, e.g. H_W is the amino hydrogen in ε1 position of the Trp side chain. HN_R is the peptidic HN of Arg while H_R corresponds to an HN of the Arg side chain. NH2_cap stands for the NH2 group linked to the C-terminus Trp residue. O_R stands for the Arg carbonyl oxygen, O_W the Trp carbonyl oxygen, O_Dye the oxygen in the BODIPY, O_L1 and O_L2 are on carbonyls of the linker closer and farther respectively to the BODIPY rings. O_P1 and O_P4 correspond to carbonyl oxygens of the first and the fourth Pro residues, respectively. NH2_cap stands for a hydrogen of the Trp amide cap and HN_L1 stands for the HN in the first amide bond of the linker. The notation “2nd H_R” is used when a heteroatom is involved in two H-bonds with the Arg side chain. 3.2. Many-Body Effects. Fluctuations of the DyeX-(Pro)4Arg+-Trp conformations were calculated at 298 K with and without polarization included in AMOEBA to identify the role of many-body effects on trajectories. To switch off polarization, the values of the atomic polarizabilities were set to zero during the simulations, keeping all other parameters unchanged. The starting structure for each trajectory was initialized by a geometry optimization. The structures optimized with and without polarization are practically identical, having an RMSD value of 0.21 Å. Five simulations were performed for each case, each having a different set of initial velocities that were randomly generated. Care was taken to ensure that each of these five pairs of MD trajectories were completely independent. The final structure of each trajectory was reoptimized and the 5 structures for each case are shown superimposed in Figure 1. These superimposed structures were created using UCSF Chimera molecular graphics and analysis software.26,27 Chimera routines were used to generate sequence alignments and graphics of the two superimposed groups of five structures shown in Figure 1. The root-mean-square deviation (RMSD) of all atoms (excluding hydrogens) was averaged by Chimera27 over the 10 pairs of five structures for each case. Structures calculated without polarization exhibited an RMSD = 0.819 Å, and those calculated with the polarizable force field showed an average RMSD = 0.106 Å. In the latter case limited dynamics of the linker generates limited disorder in the dye orientation (see Figure 1b). Figure 1a shows that without polarization, only the (Pro)4-Arg segment remains somewhat rigid. This is expected for (Pro)4 while it is the charge on Arg which likely makes it is less prone to large amplitude motion. On the other hand both the dye and Trp orientations fluctuate strongly, dominating the large increase in RMSD mentioned above. It is thus clear that polarization makes a considerable difference in the structural compactness for this singly protonated peptide. Many-body interactions could have significant consequences for a fluorescence quenching process such as electron transfer which depends exponentially on the separation between the dye and the Trp side chain and also the energetics associated with electrostatic interactions. 1748
dx.doi.org/10.1021/jp3078375 | J. Phys. Chem. B 2013, 117, 1746−1755
The Journal of Physical Chemistry B
Article
Table 2. Description of Selected Structures of Types str4, str7, and str9a
Figure 1. Structures of the peptide DyeX-(Pro)4-Arg+-Trp calculated at 300 K by molecular dynamics using the AMOEBA force field. The dye (blue), tryptophan (green), arginine (red) and proline (gray) are indicated in (a) five superimposed structures calculated without polarization having a root mean standard deviation of RMSD = 0.819, and (b) five superimposed structures calculated with polarization having RMSD = 0.106.
Dye−Trp (Å) Dye and Trp parallel Trp position Arg side chain hydrogen bonds (Pro)4 hydrogen bonds (Pro)4 amide dihedrals ω1−ω4 N_P1−CO_P4 distance (Å)
str2_400 K_3187
str2_400 K_3484
str2_400 K_4235
5.17 yes exo 3
6.59 yes endo 3
6.72 no endo 3
3.93 yes endo 3
2
2
2
3
t-t-c-t
t-t-c-t
t-t-c-t
t-t-c-t
7.69
9.08
7.47
7.08
str4_500 K _4861
str7_400 K _1420
str7_400 K _4277
str9_400 K _1615
Dye−Trp (Å) Dye and Trp parallel Trp position Arg side chain hydrogen bonds (Pro)4 hydrogen bonds (Pro)4 amide dihedrals ω1−ω4 N_P1−CO_P4 distance (Å)
7.57 no
6.43 no
5.27 no
4.52 yes
5.55 no
exo 4
exo 2
exo 5
exo 5
exo 4
3
0
2
2
2
t-t-t-t
t-t-t-t
c-t-c-t
c-t-c-t
t-t-c-t
9.16
11.15
7.16
6.72
9.32
Values are measured using geometries optimized at the AMOEBA level. Description of the table data is identical to Table 1.
dihedral angle of the amide bond between the linker and P1. In str4, ω dihedrals are all trans leading to more extended (Pro)4 segments. Therefore, trajectories starting from these various structures will explore different areas of the potential energy surface. A set of ten conformations with different local interactions of the Dye, Arg, and Trp and including some variations in the Pro4 backbone conformation was selected from the MD trajectories, corresponding to potential energy plateaus. Four were taken from str2 and three from str7 since trajectories propagated from these two starting points presented a variety of conformers leading to Dye−Trp proximity. These structures are named according to the initial structure, the thermostat temperature and the time index in units of 1 ps of its trajectory (e.g., str2_400 K_1728 is the structure formed after a propagation of 1728 ps at 400 K starting from str2). Two structures were taken from the trajectory starting from str4 since it involves a different conformation for the (Pro)4 segment. The description of structures in terms of hydrogen bonds (considered if distances are smaller than 2.5 Å) and π stacking between the Dye and Trp side chain is given in Tables 1 and 2. Two of these structures, illustrating cases where the Dye and Trp side chain are in proximity and nearly parallel, are displayed in Figure 5. Steric hindrance between the side chains and electrostatic repulsion between carbonyl groups impose limitations on the number of hydrogen bonds between (Pro)4 and the Arg side chain. In the DyeX-(Pro)4-Arg+-Trp sequence, the (Pro)4 segment “solvates” the protonated Arg side chain most of the time with two of its carbonyl groups. The (Pro)4 conformation fluctuates with little consequence, except for exchanging two Pro carbonyls in the Arg side chain solvation. Structural transitions principally concern the peptide terminations: Trp side chain at the C terminus and the Dye at the N terminus, which is attached by a linker that has high flexibility. As the peptide is globularly structured around the Arg side chain, the linker primarily allows the peptide to be more or less compact, eventually leading to structures in which the Trp side chain is free enough to move inside the globular structure. Thus, the Trp side chain is quite mobile and can be either hindered in an internal position or move freely outside the globular structure. Because of this high flexibility, the Trp side chain can occasionally form a hydrogen bond with the first Pro residue while the Trp NH2 cap can bind to the carbonyl of the fourth
Table 1. Description of Selected Structures of Type str2a str2_400 K_1728
str4_all _0000
a
3.3. Dye−Trp Proximity Analysis. The polyproline sequence is known to have significant backbone rigidity. In the present case, stiffness of the (Pro)4 segment maintains its conformation during all trajectories from 150 to 500 K. Thus, according to its initial structure, a trajectory may generate some structures that others could not. The use of five different initial structures for molecular dynamics simulations at each temperature partly compensates for the bias introduced by trajectories which individually depend on their starting structure. At low temperatures, structures vibrate without leaving their initial potential wells. As the temperature increases, structural transitions are observed, leading to several new minima. The choice of starting structures utilized the sampling previously carried out.4 Structures str2 and str8 have similar (Pro)4 conformations, whereas this segment is different for str7, str9, and str4. These are summarized in Tables 1 and 2, not starting structures, where each amide bond linkage is denoted as “t” for trans or “c” for cis. For str9, this difference also involves the first Proline residue and is the result of a rotation of ψP1−P2. The (Pro)4 conformation in str7 differs by a rotation of the ω
structural property
structural property
a
Values are measured using geometries optimized at the AMOEBA level. Dye and Trp are considered to be parallel when |cos(Dye/Trp)| > 0.9. (Pro)4 amide dihedrals are numbered from the O_linker-N_P1 peptide bond to O_P4−N_Arg where Pn is the nth proline. A hydrogen bond is assumed to exist when the corresponding H---X distance is less than 2.5 Å. 1749
dx.doi.org/10.1021/jp3078375 | J. Phys. Chem. B 2013, 117, 1746−1755
The Journal of Physical Chemistry B
Article
Figure 2. Dye−Trp proximity occurrence and correlation with geometry in molecular dynamics trajectories starting from structures (a) str2, (b) str4, (c) str7, (d) str8, (e) str9, and (f) the total occurrence for the 5 trajectories at each temperature, from 150 to 500 K. Proximity is considered for Dye−Trp distance
