An Investigation on the Fundamental Interaction between Abeta

Jun 18, 2015 - School of Physical and Mathematical Sciences, Nanyang Technological University, 21 Nanyang Link, Singapore. ‡ School of Computer Engi...
3 downloads 10 Views 3MB Size
Article pubs.acs.org/JPCB

An Investigation on the Fundamental Interaction between Abeta Peptides and the AT-Rich DNA Li Na Zhao,†,‡,§ Jie Zheng,‡,∥ Lock Yue Chew,*,†,⊥ and Yuguang Mu*,# †

School of Physical and Mathematical Sciences, Nanyang Technological University, 21 Nanyang Link, Singapore School of Computer Engineering, Nanyang Technological University, 50 Nanyang Avenue, Singapore § Bioinformatics Institute, 30 Biopolis Street, Singapore 138671 ∥ Genome Institute of Singapore, A* STAR, 60 Biopolis Street, Singapore 138672 ⊥ Complexity Institute, Nanyang Technological University, 18 Nanyang Drive, Singapore # School of Biological Sciences, Nanyang Technological University, 60 Nanyang Drive, Singapore ‡

S Supporting Information *

ABSTRACT: DNA damage is ubiquitous in all mammalian cells with the occurrence of more than 60 000 times per day per cell. In particular, DNA damage in neurons is found to accumulate with age and has been suggested to interfere with the synthesis of functional proteins. Moreover, recent studies have found through transgenic mice that human amyloid precursor protein causes an increase in DNA double-strand breaks (DSBs) with the effect of a prolongation in DNA repair. It is surmised that amyloid β (Aβ) exacerbates the DNA DSBs in neurons, possibly engendering neuronal dysfunction as a result. However, a good understanding on the holistic interaction mechanisms and the manner in which Aβ intertwines with DNA damage is still in its infancy. In our study, we found that DNA with an AT-rich sequence has a very low binding affinity toward Aβ by means of molecular dynamics simulation. While we have pursued a particular sequence of DNA in this study, other DNA sequences are expected to affect the interaction and binding affinity between DNA and Aβ, and will be pursued in our further research. Nonetheless, we have uncovered favorable interaction between the positively charged side chain of Aβ and the two ends of DNA. The latest experiment reveals that many of the double-stranded breaks in neurons can be fixed via DNA repair mechanisms but not in the case that Aβs are present. It is found that the increased numbers of DSBs prevail in active neurons. Here, on the basis of the favorable interaction between Aβ and the two ends of DNA, we propose the possibility that Aβ prevents DNA repair via binding directly to the break ends of the DNA, which further exacerbates DNA damage. Moreover, we have found that the base pair oxygen of the DNA has a greater preference to form hydrogen bonds than the backbone oxygen with Aβ at the two ends. Thus, we postulate that Aβ could serve to prevent the repair of AT-rich DNA, and it is unlikely to cause its breakage or affect its binding toward histone. Another important observation from our study is that AT-rich DNA has very little or no influence on Aβ oligomerization. Finally, even though we do not observe any dramatic DNA conformational change in the presence of Aβ, we do observe an increase in diversity of the DNA structural parameters such as groove width, local base step, and torsional angles in lieu of Aβ interactions.



INTRODUCTION Alzheimer’s disease (AD) is believed to be the result of an accumulation of senile plaques together with the presence of neurofibrillary tangles. In fact, a major component of senile plaque is the amyloid β or Aβ which is cleavaged from the amyloid precursor protein (APP) via β- and γ-secretase. The aggregation of these Aβ peptides with the formation of higherorder oligomers is one of the major etiological pathways in AD patients.1−5 Aβ oligomers have been suggested to be the main culprit that causes neurodegeneration in Alzheimer’s disease patients. The toxic mechanisms of Aβ include the following:5 its alteration of signal transduction by binding to membrane receptors;6 its disruption of membrane integrity by attachment to the membrane surface; it drives homeostasis imbalance through © 2015 American Chemical Society

the creation of an ion pore/channel which causes calcium leakages;7−10 it induces oxidative stress with the formation of an Aβ−ion complex;11,12 and it modifies DNA structures.13−15 It is well-known that there are two types of Alzheimer’s disease: early onset (or the familial AD) and later onset (or sporadic AD). Both types of AD have a genetic connection. Familial AD involves a number of single-gene mutations which occur on chromosomes 1, 14, and 21, corresponding to the abnormal presenilin 2, presenilin 1, and amyloid precursor protein production, respectively.16−20 Each of these mutations is believed to play an important role in the cleavage of APP and Received: January 29, 2015 Revised: June 6, 2015 Published: June 18, 2015 8247

DOI: 10.1021/acs.jpcb.5b00957 J. Phys. Chem. B 2015, 119, 8247−8259

Article

The Journal of Physical Chemistry B

TATA box and basically results from the fact that runs of A and T are very common within the nuclear genome.37 The initial structure of Aβ1−42 was obtained from the final structure of the first set of simulations performed in our previous studies.38,39 The B-DNA was aligned along the Z-direction and placed in the center of a box of dimensions 70 Å × 70 Å × 90 Å. The 8 vertices of the box provide 8!/(4!)2 = 70 possibilities for us to place the 4 copies of Aβ peptide (see Figure 1). In our simulation, we randomly chose 13 out of these

thus further affects Aβ production. On the other hand, the majority of AD cases comes in the form of late-onset AD. In this case, the development of AD involves genetic risk factors which have been identified to be the apolipoprotein E (APOE) gene on chromosome 19,21 the methylenetetrahydrofolate dehydrogenase 1-like (MTHFD1L) gene on chromosome 6,22 and some other loci23−25 which requires further elucidation. Research has observed a trend of lower copy-number variation (CNV) call rate in deletions as well as duplication in comparison to controls. Gene-based association analyses have confirmed previous associated genes in ADNI study (ATXN1, HLA-DPB1, RELN, DOPEY2, GSTT1, CHRFAM7A, ERBB4, NRXN1) and have identified a new gene (IMMP2L) that may play a role in AD susceptibility. The affirmation of this role warrants replication of independent samples as well as further analyses of these gene regions.26 Under physiological conditions, a “naked” double-stranded DNA exists in three general helical forms, namely, A-DNA, BDNA, and Z-DNA. The right-handed double helical B-form structure with Watson−Crick base pairing and global repetitive features is the most ubiquitous DNA form that exists within a variety of genomes. On the other hand, the A-DNA form has a relatively short and more compact helical structure in comparison to the B-DNA form. In contrast, the conversion from the A- or B-form to the left-handed Z-form requires drastic conformational changes. Meanwhile, a negatively charged DNA is found to associate readily with Aβ oligomers at a pH level of 4−9.27 Furthermore, surface plasma resonance reveals that DNA binds to all forms of Aβ oligomers regardless of their stage of aggregation.28 The supercoiled DNA (scDNA) has also been reported to uncoil to a relatively relaxed form by Aβ via helicity and superhelicity modulation. It is proposed that Aβ can alter the structure of DNA.29 Further experiments have revealed that Aβ is able to nick open the circular and linear configurations of scDNA to form the Aβ−DNA complex.30 Time-dependent DNA condensation in the presence of Aβ has also been observed with the β-sheet structure being suggested as the nuclei for the DNA condensation.13 Recent experiments have shown that Aβ oligomers can transform the left-handed ZDNA back to the right-handed B-DNA.15 It has also been reported that the binding of nucleic acid to Aβ may facilitate the conversion of Aβ to insoluble amyloid fibrills.31 With the profusion of gene expression in the brain, DNA damage may have a cumulative effect on transcriptional effectivity as well as fidelity. In AD patients, the altered DNA conformation has been observed in the hippocampus region.32−34 A recent study proposes that Aβ exacerbates DNA damage by eliciting aberrant synaptic activity.35 In our study, we shall focus on the interaction between the most ubiquitous B-form DNA and Aβ. By means of molecular dynamics simulations, we found that DNA fails to show high binding affinity with Aβ peptides, albeit the fact that Aβ1−42 is reported to cause the conformational conversion of DNA.

Figure 1. This figure illustrates the manner in which the 13 initial configurations used in our simulation is selected. Note that there are three perspectives in the figure. The DNA is arranged in the middle of the simulation box. The 8 Aβs indicate the 8 possible positions at the corner of the box to place the Aβs. These 8 positions are labeled by the number 1 to 8. For our initial configurations, we need to select 4 out of these 8 positions to put our 4 Aβ peptides. This gives rise to 70 possibilities, and we chose 13 out of these 70 possibilities as our initial configurations.

70 possibilities. Note that the four peptides are labeled as peptide A, peptide B, peptide C, and peptide D. The 13 different systems were then further configured by adding water into the box. Due to the different setup in each system, the number of water molecules varies slightly between 13 210 and 13 235. Moreover, 77 Na+ and 27 Cl− were added to arrive at a concentration of 0.1 mol/L. In addition, we have created a control system with the DNA placed in the center of a box of dimensions 40 Å × 40 Å × 90 Å, together with 46 Na+, 8 Cl−, and 3927 water molecules, to attain the same concentration of 0.1 mol/L as the 13 systems. All the structures in our simulation setup were represented by the AMBER99SB force field.40 The parallel linear constraint solver (P-LINCS) algorithm41,42 was used to constrain all bonds with an integration time step of 2 fs. The temperature was kept at 300 K using the Nosé−Hoover coupling scheme43,44 for conventional MD. An isotropic pressure coupling at 1 bar modulated by a Parrinello− Rahman barostat with a coupling constant of 2 ps was employed in the constant pressure simulation. A fourth-order interpolation was used in particle mesh Ewald method45 for electrostatic interaction with a Fourier grid spacing of 0.16 nm. Both the Coulomb and van der Waals cutoffs were set at 1.0 nm. All simulations and analysis were performed using facilities within the GROMACS package version 4.5.546−49 and AMBER.50 500 ns simulations were conducted on each Aβ + DNA system, and a 50 ns simulation was conducted on the DNA system. PyMOL51 and VMD52 were used to visualize the molecular structures. Additional analysis and visualization were assisted by originPro, GNU image manipulation program



METHODS AND SIMULATION SETUP Simulation Setup. The initial structures of B-DNA were generated via the 3DNA-Driven DNA Analysis and Rebuilding Tool (3D-DART)36 with 20 A and 20 T nucleic acids. The choice of 20 nucleic acids results from the observed strong interaction between the two DNA strands in our simulations, which provides a more stable framework for us to examine whether Aβ can induce a large conformational change in DNA. Our interest in an AT-rich DNA sequence was inspired by the 8248

DOI: 10.1021/acs.jpcb.5b00957 J. Phys. Chem. B 2015, 119, 8247−8259

Article

The Journal of Physical Chemistry B (GIMP), R, 3DNA53−55 and MATLAB as well as in-house developed scripts. Binding Free Energy Calculation. The binding free energy of the DNA−Aβ complex can be estimated via56 0 0 0 ΔG bind,solv = ΔG bind,vacuum + ΔGcomplex,solv 0 − (ΔG DNA,solv + ΔGA0β ,solv )

(1)

The free energies ΔG0complex,solv, ΔG0DNA,solv, and ΔGAβ,solv can be decomposed into an empirical term for hydrophobic contributions and a solvation energy which can be estimated by either solving the generalized Born or the linearized Poisson−Boltzmann equation. Both the generalized Born and Poisson−Boltzmann approaches serve to determine the solvation free energy contributed by the electrostatic interactions. The ΔGbind,vacuum is usually calculated via 0 ΔGvacuum = ΔEMM − T ·ΔS 0

(2)

where T is the temperature. Here ΔEMM is the molecular mechanics (MM) energies based on the force field. The entropy contribution ΔS0 can be obtained by performing normal-mode analysis, which is both computationally expensive and ineffective in making improvement to the prediction. Hence, this contribution is usually ignored in practice such that the free energy evaluated here can be considered to be the effective binding free energy.

Figure 2. Representative snapshots of the DNA−Aβ complex. The snapshots are taken from the final structure of several systems. Note that the structure of the protein and DNA are depicted in cartoon, while the secondary structures of the protein are being represented by different color.



RESULT The Free Energy between DNA and Aβ. In order to examine the final Aβ−DNA configurations during the process of Aβ oligomerization, we consider the last 2 ns of each simulation. By performing a clustering analysis on 5000 frames with a Cα-rmsd cutoff of 3 Å on the simulation results at the stated time range, we have captured a snapshot of the average structure. We observe that Aβ either interacts with the 5′-end of DNA (see Figure 2a, b, e, and f) or the 3′-end (see Figure 2c), and rarely with the middle portion of the DNA (see Figure 2d). This shows that Aβ tends to bind with the ends of the DNA but has an overall low affinity to contact with its main central structure. We next compute the possible binding free energy between Aβ and DNA within the final 2 ns. The calculated free energies are listed in Table 1. Our computation shows that the favorable nonpolar interaction that draws Aβ and DNA together is largely counterbalanced by the extremely unfavorable polar interaction, as indicated by the values of ΔGsol and ΔGgas, respectively. In particular, the calculation based on the generalized Born surface area (GBSA) method is found to give a more favorable nonpolar component than the Poisson−Boltzmann surface area (PBSA) calculations.57,58 This indicates that the GBSA method is more sensitive to the favorable binding affinity of Aβ to the ends of the DNA, as shown in Figure 2, with the affinity arising from underlying hydrophobic interactions which will be discussed in the next section (see Figure 3). Overall, the mainly positive ΔEtotal (see Table 1) which results from combining the effects of these interactions indicates a low tendency for the DNA and Aβ to form an aggregated structure. Finally, note the presence of large variations in the values of the binding free energy as well as the error bars among the snapshots, as shown in Table 1. This results from the different conformational arrangements that can arise between the Aβ and DNA interactions as well as between the Aβ peptide

interactions in the different trajectories. For example, in the case of snapshot 2, Aβ peptides are found to attach to the DNA. On the other hand, for snapshot 10, the Aβ peptides do not attach to the DNA, and as a result, they can be relatively far away from each other. In particular, note the large variation that affects the ΔGαgas term of the final binding free energy between Aβ and DNA. While this term involves both the van der Waals and electrostatic interactions, the main contribution of the variance is the electrostatic interaction due to their long-range effect. Indeed, in lieu of the greater spatial restriction in the case of snapshot 2 because of the Aβ−DNA attachment, the possible conformational change is limited and we do not expect much fluctuation in the electrostatic and van der Waals energy. This gives rise to a stable final binding free energy with a small error bar. On the other hand, without any attachment between the Aβ and DNA in the case of snapshot 10, there is a larger conformational space for the DNA and Aβ to explore. In fact, we observe both translational and rotational conformation changes in this case. The large degree of freedom allows for various possibilities, including the presence and absence of water molecules between the Aβ and DNA, leading to a diverse dielectric environment. In conjunction with the larger range of spatial arrangement between Aβ and DNA, a large difference in electrostatic energy is observed while van der Waals energy change is negligible due to their short-range nature. The consequence is the large error bar in the final binding free energy observed for the case of snapshot 10. The Interaction between DNA and Aβ. In order to determine the detailed interaction scenario of DNA and Aβ, we calculate the minimum distance between each base pair of DNA and Aβ by averaging over all 13 trajectories every 100 ps. Our computation again shows that the Aβs prefer to interact with 8249

DOI: 10.1021/acs.jpcb.5b00957 J. Phys. Chem. B 2015, 119, 8247−8259

Article

The Journal of Physical Chemistry B Table 1. Final Binding Free Energy between Aβ and DNA binding free energy (kJ/mol) generalized Born Aβ−DNA complex snapshot (1) (2) (3) (4) (5) (6) (7) (8) (9) (10) (11) (12) (13)

ΔEtotal −30.51 0.27 0.53 −1.86 −7.71 −26.47 −5.12 32.07 10.89 842.51 470.99 116.87 −5.10

± ± ± ± ± ± ± ± ± ± ± ± ±

8.52 0.01 0.02 1.13 3.22 2.76 6.95 204.47 78.04 1087.12 588.22 280.05 18.65

Poisson−Boltzmann

ΔGgasa 2536.78 2286.36 2383.70 2715.24 3059.80 2762.61 2804.72 2537.26 2829.08 3869.14 3454.80 2839.95 2657.28

± ± ± ± ± ± ± ± ± ± ± ± ±

32.34 22.09 12.82 86.19 34.93 94.61 40.92 529.16 80.62 1090.38 605.33 469.75 242.15

ΔGsolb −2567.29 −2286.08 −2383.17 −2717.10 −3067.51 −2789.07 −2809.84 −2505.19 −2818.19 −3026.63 −2983.80 −2723.08 −2662.38

± ± ± ± ± ± ± ± ± ± ± ± ±

33.54 22.08 12.80 86.34 33.89 94.37 40.75 461.49 24.27 22.13 37.00 260.38 245.58

ΔEtotal 23.61 27.80 34.75 38.95 39.43 48.28 47.08 74.91 58.38 893.45 532.67 165.39 41.70

± ± ± ± ± ± ± ± ± ± ± ± ±

6.09 3.66 1.58 2.52 5.84 3.63 7.22 204.56 77.97 1087.06 588.34 282.89 15.80

ΔGgasa 2536.78 2286.36 2383.70 2715.24 3059.80 2762.61 2804.72 2537.26 2829.08 3869.14 3454.79 2839.95 2655.28

± ± ± ± ± ± ± ± ± ± ± ± ±

32.34 22.09 12.82 86.19 34.93 94.61 40.92 529.16 80.62 1090.38 605.32 469.75 242.15

ΔGsolc −2513.17 −2258.55 −2348.95 −2676.28 −3020.37 −2714.33 −2757.64 −2462.35 −2770.70 −2975.69 −2922.12 −2674.56 −2613.58

± ± ± ± ± ± ± ± ± ± ± ± ±

32.16 23.17 12.17 86.26 31.37 93.35 40.11 454.73 23.00 22.28 34.88 252.87 238.23

ΔGgas is the sum of energy contributed by the van der Waals interaction between the four Aβ peptides and DNA from MM as well as the electrostatic energy calculated by the MM force field. bΔGsol is the sum of electrostatic and nonpolar contribution to the solvation free energy calculated by the GB and empirical method, respectively. cΔGsol is the sum of electrostatic and nonpolar contribution to the solvation free energy calculated by the PB and empirical method, respectively. a

increment becomes relatively steady. The large standard deviation indicates the differences among the different types of contacts. These results are further supported by the average COM distance among the four peptides, which shows the peptides come close to each other and aggregate during the first 100 ns. Let us next consider the interaction diversity of each trajectory. For this, we calculate the number of contacts among the four peptides (labeled as A, B, C, and D) and DNA for each trajectory at every 100 ps. We then compute the correlation between the number of contacts of the different macromolecules to decern their influence on each other. For this, we take a particular trajectory and extract the number of contacts between peptides as well as the number of contacts between peptide and DNA at each time instant. This results in a time series of number of contacts of peptide−peptide and peptide−DNA which we treat as raw data for our correlation analysis. The analysis is performed by correlating between, say, the number of contacts of peptide A to peptide B and the number of contacts of peptide A to DNA, via their corresponding raw data which are matched at the same time instant. Then, by averaging the resulting correlation coefficient across the 13 trajectories, we obtain the average correlation coefficient. Note that the correlation coefficients among peptides and DNA of each trajectory are given in the Supporting Information. The average correlation coefficient of the six trajectories with the lowest final binding energy is given in Table S8 of the Supporting Information. The average correlation coefficient of all trajectories is given in Table 2. The ratio of the average correlation coefficient with variances that are above 2 is displayed in Table 3. Note that the ratio of the average correlation coefficient and its variance provide an indication on the reliability of the correlation coefficient. One general observation is that the number of contacts formed among each peptide are positively correlated with each other. For example, when the number of contacts of A−C increases with time, one also observes a corresponding increase in the number of contacts between peptide A and peptide D with time, which signifies that peptides A, C, and D come close

the two ends of the DNA (see Figure 3a). By plotting the propensity of each residue of Aβ to locate near the DNA (see Figure 3b) and the hydrogen bond formation propensity of Aβ toward DNA (see Figure 3c), we observe that the positively charged LYS-28 has the greatest propensity to interact with DNA via hydrogen bonds. In addition, the polar SER-26 also shows a high propensity to form hydrogen bonds with DNA, followed by ASP-1, TYR-10, and LYS-16. Next, we examine the information on the hydrogen bond formation propensity of DNA toward Aβ. The cumulative frequencies of hydrogen bonds that form between each base pair of DNA and Aβs are calculated and illustrated in Figure 4. Here, cumulative frequency is defined as the total number of times hydrogen bonds are found to form between the backbone or a particular base pair of DNA and Aβ, at each time interval of 100 ps within the total simulation time of 500 ns across the 13 different systems. From Figure 4, we see that the 3′- and 5′ends of DNA have a high propensity to form hydrogen bonds as both donors and receptors with Aβ. Moreover, the Aβs are observed to have a greater chance to interact with the base pair oxygen than the backbone oxygen at these two ends. This results from the fact that the base pair oxygens at the edges of the DNA are exposed to the outside surface. The situation however reverses in the interior of the DNA with the backbone oxygen now having a higher chance of being exposed to the surface than the base pair oxygen. In order to further investigate into the affinity of Aβ peptides in the formation of hydrogen bonds toward the DNA backbone and base pairs, we have studied the atoms from A and T in Figure 5. The observation is an enhanced tendency for oxygen (O1P, O2P) from the backbone phosphate group to participate as a hydrogen bond acceptor in comparison to oxygen from the base pairs (see Figures 4 and 5). The Interaction between DNA and Each Aβ Peptide. Here, we calculate the total number of contacts, and the number of contacts between any two peptides, for all trajectories which are averaged every 10 ns. From Figure 6, we can see that the number of contacts among the four peptides increases rapidly in the first 100 ns before the 8250

DOI: 10.1021/acs.jpcb.5b00957 J. Phys. Chem. B 2015, 119, 8247−8259

Article

The Journal of Physical Chemistry B

Figure 3. (a) The minimum distance between Aβ peptides and each DNA strand. The error bars indicate the standard deviation. (b) The location propensity of each residue of Aβ toward DNA. (c) Aβ hydrogen bond formation propensity.

and aggregate together. Thus, the observed positive correlations between the number of contacts among peptides imply the presence of a cooperative mechanism among peptides in the dimer, trimer, and even the tetramer formation. On the other hand, there is very low correlation with regard to the interactions between peptides and DNA, as exemplified from results in Table 3c. This shows that Aβ peptide to peptide interaction is generally more favorable than the interaction between peptide and DNA. Moreover, the interaction between peptides and DNA is found to have very little influence on the interaction between each peptide, indicating that DNA does not help peptide oligomerization. Next, let us take “A−B”, “A−C”,“A−D”,“B−C”, “B−D”, “C− D”, “A−DNA”, “B−DNA”, “C−DNA”, and “D−DNA” as

variables to indicate the interactions between peptides and between peptides and DNA, respectively. Principal component analysis was then performed on the contact number between these 10 variables. Note that the contact number from the last 100 ns of the simulation at every 0.2 ns of the 13 sample systems has been taken to generate a matrix of dimension 6513 × 10. With the 10 variables and 6513 observations, the standard deviation, the weight by which each standardized original variable should be multiplied to get the component score, the proportion of variance and its cumulative proportion, are given in the Supporting Information as Table S10. It is important to note that, while the peptides A, B, C, and D are identical and are placed uniformly around the DNA, the interactions between peptides (and between each peptide and DNA) are dissimilar 8251

DOI: 10.1021/acs.jpcb.5b00957 J. Phys. Chem. B 2015, 119, 8247−8259

Article

The Journal of Physical Chemistry B

Figure 4. Hydrogen bond formation propensity between DNA and Aβ. It is measured by the cumulative frequency of hydrogen bonds formed through the contribution of oxygen from either the DNA base pair or its backbone to Aβ. Note that the label indexes every three nuclei acids, with the label DA meaning DNA with adenine, and DT representing DNA with thymine.

Figure 5. Atoms involved in the hydrogen bond formation after averaging over results determined from the 13 different systems. (a) The atoms N1, N3, N6, N7, N9, O1P, O2P, O3′, O4′, and O5′ from DA-5′ and DA of adenine. (b) The atoms N1, N3, O1P, O2, O2P, O3′, O4, O4′, and O5′ from DT-3′ and DT of thymine.

due to the intrinsic asymmetry within the system as well as random effects from the ambient solvation. This explains why we do not observe the occurrence of degenerate effects in our PCA results below. We plot the variance of each principal component against the component number (see Figure S1a, Supporting Information). In principal component analysis, it is usual to search for a sharp drop in principal components to decide on the reduced dimensions of the data set. In our case, we follow the same protocol and keep the first two components. For the 13 sample groups, we plot the number of contacts among each peptide and DNA as a function of the first two principal components PC1 and PC2 (see Figure S1b, Supporting Information), where the trajectories are distinguished from each other as group 1 to group 13. The original variables, “A−B”, “A−C”,“A−D”, “B− C”, “B−D”, “C−D”, “A−DNA”, “B−DNA”, “C−DNA”, and “D−DNA”, are scaled and shown as arrows. From Figure S1b (Supporting Information), we observe that each group is distributed along certain directions. Moreover, each group and its associated quasi-direction are given in Table 4. The number of contacts among these two peptides indicated by the variable

arrow is dominant in each group. It again indicates that the interaction among peptides is stronger than the interaction between the peptide and DNA in each of the cases. The Evolution of Aβ Secondary Structure. The secondary structure propensity of Aβ for all 13 trajectories was averaged every 10 ns and is illustrated in Figure 7. The overall averaged helical and β content is also calculated. By comparing the results in Figure 7 to our previous studies of Aβ without DNA, we observe a decrease of β-content from 37.27 to 30.87. Thus, we infer that the interaction between DNA and Aβ peptides cannot be neglected with respect to Aβ secondary structure evolution. The overall effect of the peptide−DNA interaction and the peptide−peptide interaction has led to a reduction in the adoption of β-structure in about 7−8 residues out of 42 residues. This amounts to a drop of approximately 17%. On the other hand, the helical content is found to stabilize at a level of about 14. 8252

DOI: 10.1021/acs.jpcb.5b00957 J. Phys. Chem. B 2015, 119, 8247−8259

Article

The Journal of Physical Chemistry B

Figure 6. Average total contact number and COM distance among peptides A, B, C, and D and the contact number between any two peptides as well as standard deviation (shown as dots).

The DNA Conformational Changes. In order to quantify the conformational change of DNA, we have calculated the minor and major groove widths of DNA with a refined phosphate−phosphate distance which takes into account the direction of the sugar−phosphate backbones.53,59 We performed analysis at two stages of this calculation for each trajectory. The first stage (or stage I) was selected at the very beginning of our MD simulation as the Aβs are getting close to the DNA, while the second stage (or stage II) was taken at the end when the Aβs are interacting with the DNA. We plot the average minor and major groove widths at stage I and stage II in parts a and b of Figure 8, respectively. From Figure 8, we observe that both the minor and major groove widths display a greater variation at the end of the simulation in comparison to that at the beginning. The increase in diversity of the intrinsic DNA conformation results from the interactions of Aβ with DNA, which is duly verified against the control system where such an increase is not observed (see Table 5). Note that the control system consists of a DNA without the presence of any Aβs. The results in Table 5 are obtained under the following condition. After scrutinizing all the trajectories, we found that half of them have lost their minor and major groove shape after base pair 16 due to interaction of the DNA N-terminal with Aβ. Hence, we have only taken the base pairs from 4 to 16 among the trajectories as well as one representative trajectory to calculate the average values for the two stages. This calculation gives an average minor groove width of 9.9 ± 0.55 Å and an average major groove width of 19.09 ± 0.75 Å for stage I and a corresponding average minor groove width of 11.14 ± 1.32 Å and average major groove width of 20.93 ± 2.13 Å for stage II. While it is generally believed that the groove dimension is largely dependent on the DNA sequence, the major groove is found to be relatively similar within the free and bound DNA, while the minor groove shows greater variation among complexes.60,61 According to the averaged value of all trajectories, we found that the minor groove width changes around 12% and the major groove around 9%. Detailed observations of each trajectory have revealed significant changes in minor groove width in comparison to changes in major groove width among several trajectories. However, for the representative trajectory, we observe a 10% width change for

the minor groove and 18% width changes for the major groove, which is contrary to expectation. A similar trend is found for the control system, with a major groove and minor groove width change of around 5.4 and 1.47%, respectively. In addition to the groove width, we have also calculated the local base step parameters of the DNA, such as shift, slide, rise, tilt, roll, and twist between each base pair step, together with the averaged torsional angles such as α, β, γ, δ, ϵ, ζ, and χ. The detailed results of this calculation are presented in Table 5. A general observation on these results is that the standard deviation of the local base step parameters is larger in stage II than in stage I. This implies that the DNA structural parameters may be influenced by the presence of the Aβ, the solvent, as well as the sodium and chlorine ions as time progresses. We notice that, even for the control system, the standard deviation of the local base step parameters at stage II is larger than that in stage I. This has led us to examine the DNA conformational change under the action of Aβ in greater detail. For this, we have performed a computation on the average roll degree against its standard deviation of all base steps for each of the 13 systems at the two stages. We observe a low standard deviation of the roll parameter with an intermediate flexible DNA structure for the 13 systems at stage I (inside the green circle of Figure 9). However, at stage II, we notice that about 70% of the DNAs have maintained in the state of intermediately flexible, with one being transformed to the state of flexible, and the rest have reached the state of very flexible. On further analysis, we have also found an absence of correlation between the flexibility at the two stages. While we have found earlier that the peptides only interact with the extremities of the DNA, our results in this section show that the effects of the interaction are not isolated. This is exhibited by our use of a B-DNA model with perfect base step parameters at the very beginning of the simulation. As the simulation evolves, these parameters are found to change, which is apparent in Table 5. This indicates that the initial interaction at the edges eventually affects the whole DNA structurally, which has led us to calculate DNA conformational change as averages over all base pairs or base pair steps in Table 5. The outcome of this calculation is a more comprehensive 8253

DOI: 10.1021/acs.jpcb.5b00957 J. Phys. Chem. B 2015, 119, 8247−8259

Article

1 1 −0.17 ± 0.14 1 −0.07 ± 0.07 0.11 ± 0.14 1 0.09 ± 0.11 −0.02 ± 0.11 −0.08 ± 0.05

A−C A−C A−C A−B B−D A−B A−B B−C C−D B−C A−C A−B A−B A−C A−D A−B B−C A−B A−D A−D B−DNA B−D B−C A−DNA A−D A−D C−D B−D A−B B−D A−C A−B A−C B−C A−D A−C A−D B−C B−D C−D C−D A−DNA A−DNA B−DNA C−DNA

1 0.07 ± 0.17 0.20 ± 0.06 −0.06 ± 0.09 0.06 ± 0.09 0.11 ± 0.15

variation

(a) Positive Correlation (Ratio >2) A−D 0.27 B−C 0.26 C−D 0.26 A−D 0.3 A−DNA 0.2 B−C 0.2 A−C 0.19 B−D 0.14 C−DNA 0.16 D−DNA 0.11 D−DNA 0.2 C−D 0.22 (b) Positive Correlation (Ratio