Oligomer Formation of Amyloid-β(29–42) from Its Monomers Using the

Article pubs.acs.org/JPCB

Oligomer Formation of Amyloid-β(29−42) from Its Monomers Using the Hamiltonian Replica-Permutation Molecular Dynamics Simulation Satoru G. Itoh†,‡ and Hisashi Okumura*,†,‡ †

Department of Theoretical and Computational Molecular Science, Institute for Molecular Science, Okazaki, Aichi 444-8585, Japan Department of Structural Molecular Science, The Graduate University for Advanced Studies, Okazaki, Aichi 444-8585, Japan

‡

S Supporting Information *

ABSTRACT: Oligomers of amyloid-β peptides (Aβ) are formed during the early stage of the amyloidogenesis process and exhibit neurotoxicity. The oligomer formation process of Aβ and even that of Aβ fragments are still poorly understood, though understanding of these processes is essential for remedying Alzheimer’s disease. In order to better understand the oligomerization process of the C-terminal Aβ fragment Aβ(29−42) at the atomic level, we performed the Hamiltonian replica-permutation molecular dynamics simulation with Aβ(29−42) molecules using the explicit water solvent model. We observed that oligomers increased in size through the sequential addition of monomers to the oligomer, rather than through the assembly of small oligomers. Moreover, solvent effects played an important role in this oligomerization process.

■

INTRODUCTION Amyloid-β peptides (Aβ) are composed of 39 to 43 amino acid residues that may form insoluble amyloid fibrils, which are associated with Alzheimer’s disease.1,2 Oligomer formation occurs during the early stage of the amyloidogenesis process. It has been reported recently that the Aβ oligomer is a more plausible candidate for synaptic dysfunction than the amyloid fibril.3,4 However, the oligomer formation process is not well understood, even though a complete understanding of the process is necessary to develop methods for remedying Alzheimer’s disease. To investigate the oligomers of full-length Aβ peptides and those of Aβ fragments, there have been many experimental and simulation studies reported in the literature.5−18 Experiments have shown that Aβ oligomers have intermolecular β-sheet structures, which can be seen in the amyloid fibrils.5−7 Aβ in the amyloid fibril has two intermolecular β-sheet structures, which consist of residues 18−26 (β1) and residues 31−42 (β2).6 It has also been shown experimentally that various Aβ fragments form oligomers and amyloid fibrils by themselves.19−21 Aβ(29−42) is one such fragment, and it consists of the residues 29−42, which correspond to the transmembrane domain of Aβ and include the residues of β2. The length of Aβ after residue 29 is a critical determinant of the oligomer and amyloid formation rate.22 Therefore, it is important to determine the oligomerization process of this fragment to understand the process in the full-length Aβ. Moreover, It has © 2016 American Chemical Society

been reported that C-terminal fragments, Aβ(x−42) (x = 29− 39), are efficient in protecting neurons from toxicity of the Aβ oligomers.23−25 To design efficient inhibitors of Aβ neurotoxicity, it is essential to study Aβ(29−42). As for simulation studies on Aβ(29−42), monomer structures and dimerization processes have been reported with a coarse-grained model8 and all-atom models.9,13,14,16 However, there is no study on Aβ(29−42) oligomers that are larger than the dimer. In order to investigate the oligomerization process of Aβ fragments at the atomic level, it is necessary to perform extremely long molecular dynamics (MD) simulations with all-atom models, as the oligomerization process is a considerably slow process. To solve this problem, the Hamiltonian replica-permutation method (HRPM) may be used.14,16 In our previous work, we were successful in determining the dimerization process of Aβ(29−42) at the atomic level by using the Coulomb replica-permutation method (CRPM), which is a form of the HRPM.14,16 The HRPM is a better alternative to the multidimensional replica-exchange method,26 which is also referred to as the Hamiltonian replica-exchange method.27 In the HRPM, parameters are introduced into a Hamiltonian of a target system, and different parameter values are assigned to copies (or replicas) of the system. During the simulation, both Received: April 15, 2016 Revised: June 8, 2016 Published: June 9, 2016 6555

DOI: 10.1021/acs.jpcb.6b03828 J. Phys. Chem. B 2016, 120, 6555−6561

Article

The Journal of Physical Chemistry B

group, respectively. The amino acid sequence was AceGAIIGLMVGGVVIA-Nme. The four Aβ(29−42) molecules were put in a cubic unit cell with 13823 water molecules. The side length of the cubic unit cell was 75.5 Å and periodic boundary conditions were utilized. The AMBER parm99SB force field30 and the TIP3P rigid-body model31 were employed for the Aβ(29−42) molecules and for the water molecules, respectively. Recent computational studies showed that the AMBER parm99SB-ILDN32 and the CHARMM22*33 force fields provide good agreement with experimental data for the Aβ peptides.34,35 By applying CRPM to Aβ(29−42) peptides with these force fields, therefore, a more accurate oligomerization process may be obtained though the AMBER parm99SB force field and generally provides the good agreement with experimental data.36,37 The SHAKE algorithm was utilized to constrain bond lengths with the hydrogen atoms of Aβ(29−42) and to fix the water molecule structures during the simulations. The cutoff distance for the Lennard−Jones potential energy was 12.0 Å. The electrostatic potential energy was calculated by the particle mesh Ewald method. Temperature was controlled by the Nosé−Hoover thermostat.38−41 The multiple-time-step method42 was employed, and the time steps were taken to be 4.0 fs for interactions between the water molecules and 1.0 fs for other interactions. We performed three CRPMD simulations with different initial velocities. Initial conformations for these CRPMD simulations were the same and were prepared as presented in the Supporting Information. The number of replicas was 18, and the values of λ were 0.82, 0.84, 0.86, 0.88, 0.90, 0.92, 0.94, 0.96, 0.98, 1.00, 1.01, 1.02, 1.03, 1.04, 1.05, 1.06, 1.07, and 1.08. These replicas were divided into three subsets, and each subset had six replicas and six parameter values (see ref 43 for more details). Each CRPMD simulation was performed for 240.0 ns at 300 K per replica, including an equilibration run performed for 10.0 ns. The production run of each simulation was conducted for 4.14 μs (= (240 − 10) ns × 18) in total. The total time length of the CRPM production runs from the three initial conditions was 12.42 μs (= 4.14 μs × 3). The trajectory data were stored every 400 fs, and trials of replica permutations were performed every 4.0 ps.

parameter exchanges between two replicas, as well as parameter permutations among more than two replicas, are performed. Furthermore, instead of using the Metropolis algorithm,28 the Suwa−Todo algorithm29 is employed for replica-permutation trials to minimize its rejection ratio. By permuting the parameters among the replicas, random walks of the replicas in the parameter space are realized. Accordingly, the simulation can escape from local minimum states, and efficient sampling in the conformational space is realized. In this paper, we study the oligomerization process of Aβ(29−42). Solvent effects in the oligomerization process are also discussed. In Methods, we briefly introduce CRPM and present the computational details. Our simulation results are shown in Results and Discussion. The last section is devoted to conclusions.

■

METHODS Coulomb Replica-Permutation Method. In the CRPM, the parameters are introduced into electrostatic interactions among solute atoms. Let us consider a system consisting of solute molecules in an explicit water solvent model, in which the number of atoms is N. The total potential energy is written as E(q; λ) = Ep(qp; λ) + Eps(qp , qs) + Es(qs)

(1)

where Ep, Eps, and Es are the potential energies among the solute molecules, between the solute and solvent molecules, and among the solvent molecules, respectively. Here q = {qp,qs}, where qp and qs are the coordinate vectors of the solute atoms and the solvent atoms, respectively, and are denoted as qp ≡ {q1,···,qNp} and qs ≡ {qNp+1,···, qN} (Np is the total number of the solute atoms), respectively. A parameter λ is introduced to the Coulomb potential energy term in Ep. To enhance sampling between aggregation and disaggregation states for Aβ(29−42) molecules, we employed λ only for the intermolecular electrostatic interactions among the Aβ(29− 42) molecules. For the Aβ(29−42) molecules M1−Mn (n is the number of the Aβ(29−42) molecules), the Coulomb potential energy is given by n−1 n

Eelec(qp; λ) =

∑∑ ∑ ∑ i = 1 j > i k ∈ M i l ∈ Mj

n

+

∑∑

∑

i= 1 k ∈ Mi l> k /l∈ Mi

■

(λQ k)(λQ l)

RESULTS AND DISCUSSION Because the four Aβ(29−42) molecules are in the system, five different states can exist: (i) All molecules have the monomer state (4M state). (ii) Two of the four molecules form a dimer, and the other two molecules have the monomer state (2M1D state). (iii) Three molecules form a trimer, and the remaining molecule has the monomer state (1M1Tr state). (iv) Two of the four molecules form a dimer, and the remaining two molecules form another dimer (2D state). (v) The four molecules form a tetramer (1Te state). In our simulations, a dimer formation was considered when the shortest intermolecular Cα−Cα distance between two monomers became less than 6.5 Å. Similarly, when the shortest intermolecular Cα−Cα distance between a monomer and an oligomer became less than 6.5 Å, it was regarded as an addition of the monomer to the oligomer. To classify the five states, we introduced the following three reaction coordinates (see Figure 1). The reaction coordinate d1 was the shortest intermolecular Cα−Cα distance among the four molecules. The reaction coordinate d2 was the shortest intermolecular Cα−Cα distance between Cα atoms in the two molecules that were used to calculate d1 and Cα atoms in the

4π ϵ0rkl

Q kQ l 4π ϵ0rkl

, (2)

where rkl is the distance between atoms k and l in the Aβ(29− 42) molecules and ϵ0 is the dielectric constant in vacuum. In each replica, the MD or Monte Carlo simulation is performed with this electrostatic potential energy. Note that eq 2 corresponds to introducing a scaling factor λ to the electrostatic charge Qk of the solute atom k: Qk → λQk. In general, when target peptides (or proteins) have charges, in order to calculate electrostatic interactions by utilizing the Ewald method or the particle mesh Ewald method, the electrostatic charges of counterions also have to be scaled so that the total charge of the system would be neutral. Computational Details. To investigate the oligomerization process of Aβ(29−42), we applied Coulomb replica-permutation MD (CRPMD) simulations to four Aβ(29−42) molecules using an explicit water solvent model. The N-terminus and Cterminus were blocked by the acetyl group and the N-methyl 6556

DOI: 10.1021/acs.jpcb.6b03828 J. Phys. Chem. B 2016, 120, 6555−6561

Article

The Journal of Physical Chemistry B

Figure 1. Schematic illustration of reaction coordinates d1, d2, and d3. The orange and blue circles show Aβ(29−42) molecules. The red and black arrows show the shortest intermolecular Cα−Cα distance between the two molecules. (a) The red arrow corresponds to d1, which is the shortest intermolecular Cα−Cα distance among the four molecules. (b) The blue circles show the two molecules that were used to calculate d1. The orange circles show the two molecules that were not used to calculate d1. The red arrow corresponds to d2. (c) The blue circles show the molecules that were used to calculate d1 and d2. The orange circles show the molecule that was not used to calculate d1 and d2. The red arrow corresponds to d3.

two molecules that were not used to calculate d1. The reaction coordinate d3 was the shortest intermolecular Cα−Cα distance between the molecule that was neither used for d1 nor for d2 and the other three molecules. When d1 took a smaller value than 6.5 Å, for example, a dimer was considered to have been formed. To visualize the oligomerization pathway from the 4M state to the 1Te state, we calculated the free energy surface (or potential of mean force) with respect to these three reaction coordinates at 300 K and λ = 1.00 using reweighting techniques.44,45 The free energy surface F (d1, d2, d3) is given by F(d1 , d 2 , d3) = −kBT0 ln P(d1 , d 2 , d3)

Figure 2. (a) Free energy surface F (d1,d2,d3) at 300 K and the parameter value of 1.00. Blue surface, red surface, and yellow mesh show contour surfaces of F = 5.0 kcal/mol, F = 7.0 kcal/mol, and F = 11.5 kcal/mol. (b) Free energy surface F(d1 = 4.0,d2,d3) at 300 K and the parameter value of 1.00. The unit of the free-energy is kcal/mol. Black circles correspond to the conformations in Figure 3.

(3)

where P(d1,d2,d3) is the probability distribution at 300 K and λ = 1.00, which is calculated by the reweighting techniques. This free energy surface is shown in Figure 2(a). Figure 2(b) shows the free energy surface at d1 = 4.0 Å at which the global minimum state is located. Here, the bin size of the reaction coordinates was 2.0 Å. From Figure 2, parts (a) and (b), the oligomerization pathway is as follows (see the blue arrows): During the first step of the oligomerization, two of the four monomers come close to each other and then form a dimer (4M → 2M1D). During the second step, a trimer is formed by approaching one of the two monomers to the dimer (2M1D → 1M1Tr). During the final step, the trimer and the remaining monomer assemble, and a tetramer is formed (1M1Tr → 1Te). This oligomerization pathway implies that an oligomer increases in size by the addition of monomers rather than by the assembly of small oligomers. Examples of growth through the addition of monomers were observed in the oligomerization process and amyloid formation process of the full-length Aβ. Several experiments reported that small oligomers increase in size by the addition of monomers.46−48 It had also been reported experimentally that the amyloid fibril elongates by adding monomers to the amyloid fibril, rather than by adding oligomers to the amyloid fibril.49,50 Our simulation results also imply that the oligomerization process through the assembly of small oligomers can occur, though this process is disadvantageous for the free energy. In fact, the oligomerization process of the full-length Aβ through the assembly of small oligomers also had been reported experimentally10 and computationally.11 As for other Aβ fragments, both of the oligomerization processes, through the addition of monomers and the assembly of small oligomers, were observed in simulation studies.12,15,18

In Figure 3, the lowest potential energy conformation in the 4M state and the longest intermolecular β-sheet conformations

Figure 3. Typical conformations in the 4M, 2M1D, 1M1Tr, 2D, and 1Te states.

in the other states are presented. We chose these conformations from a replica occupying the parameter value of 1.00. More details of how we chose these conformations are shown in the Supporting Information. From these conformations, it was found that each molecule in the tetramer formed intermolecular β-sheets longer than that in the trimer. Each molecule in the trimer also formed intermolecular β-sheets longer than that in the dimer. Here, the DSSP (define secondary structure of proteins) criteria51 was utilized to determine secondary structures of Aβ(29−42). 6557

DOI: 10.1021/acs.jpcb.6b03828 J. Phys. Chem. B 2016, 120, 6555−6561

Article

The Journal of Physical Chemistry B To see the tendency of secondary structures for each oligomer size, average values of the numbers of residues that have the secondary structures were calculated by the reweighting techniques. Figure 4(a) shows the average values

Figure 5. Averaged numbers of hydration waters per side-chain heavy atom in the 4M, 2M1D, 1M1Tr, 2D, and 1Te states.

Figure 4. (a) Average values of the numbers of residues that have secondary structures with respect to the oligomer size at 300 K and the parameter value of 1.00. Intermolecular β, intramolecular β, and αhelix are shown by the open black triangles, the open red circles, and the open blue squares, respectively. (b) The logarithms of the probability distributions of the numbers of residues that form the intermolecular β-bridge structures with respect to the oligomer size at 300 K and the parameter value of 1.00. These results (a and b) are for single molecules in each oligomer size.

oligomerization pathway. Moreover, the number of hydration waters of the 1M1Tr state was smaller than that of the 2D state. This result may explain why the free energy of the 1M1Tr state was lower than that of the 2D state, as seen in Figure 2(b). Note that we employed only the side-chain atoms to count the number of hydration waters for the hydrophobic parts of Aβ(29−42). There is no qualitative difference when the backbone atoms are included. The 1M1Tr state was favorable to the 2D state in terms of not only the number of hydration waters but also the excluded volume. Figure 6 shows the excluded volume of the early stages

for a single molecule in each oligomer size at 300 K and λ = 1.00. The errors were estimated by the jackknife method.52 The number of bins for the jackknife method was five for each CRPMD simulation, and the total number of bins was 15 because three sets of initial conditions were used. As seen in the figure, Aβ(29−42) has a tendency to form the β-sheet structure rather than the α-helix structure in all oligomer states. Moreover, the tendency of secondary structures is not much different among the oligomer states except for the intermolecular β-sheet structure in the monomer state. The probability distribution of the number of residues that form the intermolecular β-bridge structure was also calculated. Figure 4(b) shows the logarithm of the probability distribution for a single molecule in each oligomer size at 300 K and λ = 1.00. Aβ(29−42) can form a longer intermolecular β-sheet structure in a larger oligomer state although the probability of the longer intermolecular β-sheet structure is not very high. These results imply that to form a long intermolecular β-sheet structure, assembly of at least several Aβ(29−42) molecules and adequate time are required. This is consistent with experiments that show that lag time exists to form the Aβ amyloid under high concentration of Aβ. Solvent effects may play an important role in the oligomerization process of Aβ(29−42) because most of the residues of Aβ(29−42) are hydrophobic. In fact, it has been reported in experiments that Aβ(29−42) is only sparingly soluble.19 The number of hydration waters for the side chains was counted for each state in order to investigate the solvent effects for the oligomers. In Figure 5, the average number of hydration waters at 300 K and λ = 1.00 for the side-chain atoms with the exception of the hydrogen atoms is shown. Here, when the distance between the oxygen atom of a water molecule and the side-chain atoms (with the exception of the hydrogen atoms) was less than 5.0 Å, the water molecule was regarded as the hydration water.17 The average values were calculated by the reweighting techniques,44,45 and the errors were estimated by the jackknife method.52 From the figure, it is shown that the number of hydration waters was reduced along the

Figure 6. Schematic illustration of excluded volume for water molecules in the early stages of the 1M1Tr and 2D states. Black dashed line shows the excluded volume. Red dashed line shows the decrease of the excluded volume from the 4M state. The excluded volume of the 1M1Tr state is smaller than that of the 2D state.

of the 1M1Tr state and the 2D state schematically. The yellow circles in this figure correspond to the Aβ(29−42) molecules. As shown in the figure, the excluded volume of the 1M1Tr state was smaller than that of the 2D state. Because the 1M1Tr state has a smaller excluded volume, water molecules have more space for movement.53,54 This means that the early stage of the 1M1Tr state is more favorable than that of the 2D state in terms of the solvent entropy. The Aβ(29−42) molecules in the 1M1Tr state may be assumed to assemble linearly rather than triangularly because Aβ has a cross-β structure in the amyloid fibril. In oligomerization process, however, as supported by 6558

DOI: 10.1021/acs.jpcb.6b03828 J. Phys. Chem. B 2016, 120, 6555−6561

Article

The Journal of Physical Chemistry B several experimental and simulation studies,55−59 an unstructured molten oligomer58 forms before the ordered intermolecular β-sheet structures form. In other words, the unstructured molten oligomer is formed in the early stage of oligomer formation. In the unstructured molten oligomer, molecules have intermolecular side-chain contacts between hydrophobic residues.57−59 We showed in our previous work16 that the intermolecular side-chain contacts between Aβ(29− 42) molecules occur when the intermolecular distance between the Aβ(29−42) molecules is less than 10 Å. As seen in Figure 2(b), both free energy barriers between the 2M1D state and the 1M1Tr state and between the 2M1D state and the 2D state were located at approximately 10 Å. In order to see the effect of the excluded volume on the free energy barriers, we calculated the excluded volume difference between the early stage of the 1M1Tr state and that of the 2D state. The solvent entropy difference ΔSex with respect to the excluded volume difference ΔVex is estimated by ΔSex = −kBNsolv

ΔVex V

solvent enthalpy by utilizing the molecular mechanics Poisson− Boltzmann surface area (MM-PBSA) approach60 or the molecular mechanics three-dimensional reference interaction site model (MM-3DRISM) approach.61

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.6b03828. The initial conformation for the CRPMD simulations, detailed description of the conformations in Figure 3, and the excluded volume calculation (PDF)

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Phone: +81 (0)564 557277. Fax: +81 (0)564 557025. Notes

The authors declare no competing financial interest.

(4)

■

where Nsolv is the number of the solvent molecules and V is the volume of the simulation box.53,54 Details of the excluded volume calculation are presented in the Supporting Information. The definitions of the early stages of the 1M1Tr and 2D states are also presented in the Supporting Information, and these definitions were determined to include the vicinity of the free-energy barriers connecting these states and the 2M1D state. The contribution of the excluded volume difference to the free-energy difference ΔFex is given by ΔFex = −T0ΔSex. This contribution ΔFex was 3.3 kcal/mol, and the excluded volume of the early stage of the 1M1Tr state was smaller than that of the 2D state. On the other hand, the total potential energy of the 1M1Tr state was higher than that of the 2D state by 0.4 kcal/mol. Therefore, the excluded volume effects contributed to reduce the free energy of the early stage of the 1M1Tr state compared with that of the 2D state.

ACKNOWLEDGMENTS The computations were performed on computers at the Research Center for Computational Science, Okazaki Research Facilities, National Institutes of Natural Sciences. This work was supported by JSPS KAKENHI (24740296 and 26102550) and the Okazaki Orion Project.

■

REFERENCES

(1) Sipe, J. D. Amyloidosis. Annu. Rev. Biochem. 1992, 61, 947−975. (2) Chiti, F.; Dobson, C. M. Protein Misfolding, Functional Amyloid, and Human Disease. Annu. Rev. Biochem. 2006, 75, 333−366. (3) Haass, C.; Selkoe, D. J. Soluble Protein Oligomers in Neurodegeneration: Lessons from the Alzheimer’s Amyloid β-Peptide. Nat. Rev. Mol. Cell Biol. 2007, 8, 101−112. (4) Shankar, G. M.; Li, S.; Mehta, T. H.; Garcia-Munoz, A.; Shepardson, N. E.; Smith, I.; Brett, F. M.; Farrell, M. A.; Rowan, M. J.; Lemere, C. A.; et al. Amyloid-β Protein Dimers Isolated Directly from Alzheimer’s Brains Impair Synaptic Plasticity and Memory. Nat. Med. 2008, 14, 837−842. (5) Petkova, A. T.; Ishii, Y.; Balbach, J. J.; Antzutkin, O. N.; Leapman, R. D.; Delaglio, F.; Tycko, R. A Structural Model for Alzheimer’s BetaAmyloid Fibrils Based on Experimental Constraints from Solid State NMR. Proc. Natl. Acad. Sci. U. S. A. 2002, 99, 16742−16747. (6) Lührs, T.; Ritter, C.; Adrian, M.; Riek-Loher, D.; Bohrmann, B.; Döbeli, H.; Schubert, D.; Riek, R. 3D Structure of Alzheimer’s Amyloid-Beta(1−42) Fibrils. Proc. Natl. Acad. Sci. U. S. A. 2005, 102, 17342−17347. (7) Sarkar, B.; Mithu, V. S.; Chandra, B.; Mandal, A.; Chandrakesan, M.; Bhowmik, D.; Madhu, P. K.; Maiti, S. Significant Structural Differences between Transient Amyloid-β Oligomers and Less-Toxic Fibrils in Regions Known to Harbor Familial Alzheimer’s Mutations. Angew. Chem., Int. Ed. 2014, 53, 6888−6892. (8) Lu, Y.; Wei, G.; Derreumaux, P. Effects of G33A and G33I Mutations on the Structures of Monomer and Dimer of the Amyloid-β Fragment 29−42 by Replica Exchange Molecular Dynamics Simulations. J. Phys. Chem. B 2011, 115, 1282−1288. (9) Itoh, S. G.; Okamoto, Y. Amyloid-β(29−42) Dimer Formations Studied by a Multicanonical-Multioverlap Molecular Dynamics Simulation. J. Phys. Chem. B 2008, 112, 2767−2770. (10) Bernstein, S. L.; Dupuis, N. F.; Lazo, N. D.; Wyttenbach, T.; Condron, M. M.; Bitan, G.; Teplow, D. B.; Shea, J.-E.; Ruotolo, B. T.; Robinson, C. V.; et al. Amyloid-β Protein Oligomerization and the Importance of Tetramers and Dodecamers in the Aetiology of Alzheimer’s Disease. Nat. Chem. 2009, 1, 326−331.

■

CONCLUSIONS We employed the four Aβ(29−42) molecules with the explicit water solvent model to study the oligomerization process of Aβ(29−42) at the atomic level. The CRPMD simulations were conducted to realize the efficient conformational sampling of Aβ(29−42). For the four Aβ(29−42) molecules, it was possible to consider five different states, labeled 4M, 2M1D, 1M1Tr, 2D, and 1Te. To classify these states, we introduced the reaction coordinates d1, d2, and d3. From the free energy surface at 300 K and λ = 1.00 with respect to d1, d2, and d3, the most probable oligomerization pathway from the 4M state to the 1Te state was as follows: 4M → 2M1D → 1M1Tr → 1Te. This pathway implied that an oligomer increases in size by addition of a monomer to the oligomer sequentially, rather than by assembly of small oligomers. Growth through the sequential addition of monomers has been previously reported for the amyloid fibrils. We then showed that solvent effects play an important role in the oligomerization process of Aβ(29−42). The number of hydration waters in the 1M1Tr state was smaller than that in the 2D state. Not only the number of hydration waters but also the excluded volume of the 1M1Tr state was smaller than that of the 2D state. Therefore, the 1M1Tr state was considered to be more favorable than the 2D state in terms of the solvent entropy. In the future, we will estimate the solvent effects on the Aβ oligomerization process more precisely including the 6559

DOI: 10.1021/acs.jpcb.6b03828 J. Phys. Chem. B 2016, 120, 6555−6561

Article

The Journal of Physical Chemistry B (11) Urbanc, B.; Betnel, M.; Cruz, L.; Bitan, G.; Teplow, D. B. Elucidation of Amyloid β-Protein Oligomerization Mechanisms: Discrete Molecular Dynamics Study. J. Am. Chem. Soc. 2010, 132, 4266−4280. (12) Riccardi, L.; Nguyen, P. H.; Stock, G. Construction of the Free Energy Landscape of Peptide Aggregation from Molecular Dynamics Simulations. J. Chem. Theory Comput. 2012, 8, 1471−1479. (13) Itoh, S. G.; Okumura, H. Coulomb Replica-Exchange Method: Handling Electrostatic Attractive and Repulsive Forces for Biomolecules. J. Comput. Chem. 2013, 34, 622−639. (14) Itoh, S. G.; Okumura, H. Hamiltonian Replica-Permutation Method and Its Applications to an Alanine Dipeptide and Amyloidβ(29−42) Peptides. J. Comput. Chem. 2013, 34, 2493−2497. (15) Baftizadeh, F.; Pietrucci, F.; Biarnes, X.; Laio, A. Nucleation Process of a Fibril Precursor in the C-Terminal Segment of Amyloidbeta. Phys. Rev. Lett. 2013, 110.16810310.1103/PhysRevLett.110.168103 (16) Itoh, S. G.; Okumura, H. Dimerization Process of Amyloidβ(29−42) Studied by the Hamiltonian Replica-Permutation Molecular Dynamics Simulations. J. Phys. Chem. B 2014, 118, 11428−11436. (17) Okumura, H.; Itoh, S. G. Amyloid Fibril Disruption by Ultrasonic Cavitation: Nonequilibrium Molecular Dynamics Simulations. J. Am. Chem. Soc. 2014, 136, 10549−10552. (18) Barz, B.; Wales, D. J.; Strodel, B. A Kinetic Approach to the Sequence-Aggregation Relationship in Disease-Related Protein Assembly. J. Phys. Chem. B 2014, 118, 1003−1011. (19) Hilbich, C.; Kisters-Woike, B.; Reed, J.; Masters, C. L.; Beyreuther, K. Aggregation and Secondary Structure of Synthetic Amyloid βA4 Peptides of Alzheimer’s Disease. J. Mol. Biol. 1991, 218, 149−163. (20) Barrow, C. J.; Yasuda, A.; Kenny, P. T. M.; Zagorski, M. G. Solution Conformations and Aggregational Properties of Synthetic Amyloid β-Peptides of Alzheimer’s Disease - Analysis of Circular Dichroism Spectra. J. Mol. Biol. 1992, 225, 1075−1093. (21) Serpell, L. C. Alzheimer’s Amyloid Fibrils: Structure and Assembly. Biochim. Biophys. Acta, Mol. Basis Dis. 2000, 1502, 16−30. (22) Jarrett, J. T.; Berger, E. P.; Lansbury, P. T., Jr. The Carboxy Terminus of the β Amyloid Protein Is Critical for the Seeding of Amyloid Formation: Implications for the Pathogenesis of Alzheimer’s Disease. Biochemistry 1993, 32, 4693−4697. (23) Fradinger, E. A.; Monien, B. H.; Urbanc, B.; Lomakin, A.; Tan, M.; Li, H.; Spring, S. M.; Condron, M. M.; Cruz, L.; Xie, C.-W.; et al. C-Terminal Peptides Coassemble into Aβ42 Oligomers and Protect Neurons against Aβ42-Induced Neurotoxicity. Proc. Natl. Acad. Sci. U. S. A. 2008, 105, 14175−14180. (24) Li, H.; Monien, B. H.; Fradinger, E. A.; Urbanc, B.; Bitan, G. Biophysical Characterization of Aβ42 C-Terminal Fragments: Inhibitors of Aβ42 Neurotoxicity. Biochemistry 2010, 49, 1259−1267. (25) Li, H.; Monien, B. H.; Lomakin, A.; Zemel, R.; Fradinger, E. A.; Tan, M.; Spring, S. M.; Urbanc, B.; Xie, C.-W.; Benedek, G. B.; et al. Mechanistic Investigation of the Inhibition of Aβ42 Assembly and Neurotoxicity by Aβ42 C-Terminal Fragments. Biochemistry 2010, 49, 6358−6364. (26) Sugita, Y.; Kitao, A.; Okamoto, Y. Multidimensional ReplicaExchange Method for Free-Energy Calculations. J. Chem. Phys. 2000, 113, 6042−6051. (27) Fukunishi, H.; Watanabe, O.; Takada, S. On the Hamiltonian Replica Exchange Method for Efficient Sampling of Biomolecular Systems: Application to Protein Structure Prediction. J. Chem. Phys. 2002, 116, 9058−9067. (28) Metropolis, N.; Rosenbluth, A. W.; Rosenbluth, M. N.; Teller, A. H.; Teller, E. Equation of State Calculations by Fast Computing Machines. J. Chem. Phys. 1953, 21, 1087−1092. (29) Suwa, H.; Todo, S. Markov Chain Monte Carlo Method without Detailed Balance. Phys. Rev. Lett. 2010, 105, 120603. (30) Hornak, V.; Abel, R.; Okur, A.; Strockbine, B.; Roitberg, A.; Simmerling, C. Comparison of Multiple Amber Force Fields and Development of Improved Protein Backbone Parameters. Proteins: Struct., Funct., Genet. 2006, 65, 712−725.

(31) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79, 926−935. (32) Lindorff-Larsen, K.; Piana, S.; Palmo, K.; Maragakis, P.; Klepeis, J. L.; Dror, R. O.; Shaw, D. E. Improved Side-Chain Torsion Potentials for the Amber ff99SB Protein Force Field. Proteins: Struct., Funct., Genet. 2010, 78, 1950−1958. (33) Piana, S.; Lindorff-Larsen, K.; Shaw, D. E. How Robust Are Protein Folding Simulations with Respect to Force Field Parameterization? Biophys. J. 2011, 100, L47−L49. (34) Nguyen, P. H.; Li, M. S.; Derreumaux, P. Effects of All-Atom Force Fields on Amyloid Oligomerization: Replica Exchange Molecular Dynamics Simulations of the Aβ16−22 Dimer and Trimer. Phys. Chem. Chem. Phys. 2011, 13, 9778−9788. (35) Somavarapu, A. K.; Kepp, K. P. The Dependence of Amyloid-β Dynamics on Protein Force Fields and Water Models. ChemPhysChem 2015, 16, 3278−3289. (36) Lange, O. F.; van der Spoel, D.; de Groot, B. L. Scrutinizing Molecular Mechanics Force Fields on the Submicrosecond Timescale with NMR Data. Biophys. J. 2010, 99, 647−655. (37) Cino, E. A.; Choy, W.-Y.; Karttunen, M. Comparison of Secondary Structure Formation Using 10 Different Force Fields in Microsecond Molecular Dynamics Simulations. J. Chem. Theory Comput. 2012, 8, 2725−2740. (38) Nosé, S. A Molecular-Dynamics Method for Simulations in the Canonical Ensemble. Mol. Phys. 1984, 52, 255−268. (39) Nosé, S. A Unified Formulation of the Constant Temperature Molecular Dynamics Methods. J. Chem. Phys. 1984, 81, 511−519. (40) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A: At., Mol., Opt. Phys. 1985, 31, 1695−1697. (41) Okumura, H.; Itoh, S. G.; Okamoto, Y. Explicit Symplectic Integrators of Molecular Dynamics Algorithms for Rigid-Body Molecules in the Canonical, Isobaric-Isothermal, and Related Ensembles. J. Chem. Phys. 2007, 126, 084103. (42) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: Oxford, 1987. (43) Itoh, S. G.; Okumura, H. Replica-Permutation Method with the Suwa-Todo Algorithm beyond the Replica-Exchange Method. J. Chem. Theory Comput. 2013, 9, 570−581. (44) Ferrenberg, A. M.; Swendsen, R. H. Optimized Monte Carlo Data Analysis. Phys. Rev. Lett. 1989, 63, 1195−1198. (45) Kumar, S.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A.; Rosenberg, J. M. The Weighted Histogram Analysis Method for FreeEnergy Calculations on Biomolecules. I. The Method. J. Comput. Chem. 1992, 13, 1011−1021. (46) Bitan, G.; Lomakin, A.; Teplow, D. Amyloid β-Protein Oligomerization - Prenucleation Interactions Revealed by PhotoInduced Cross-Linking of Unmodified Proteins. J. Biol. Chem. 2001, 276, 35176−35184. (47) Ono, K.; Condron, M. M.; Teplow, D. B. StructureNeurotoxicity Relationships of Amyloid Beta-Protein Oligomers. Proc. Natl. Acad. Sci. U. S. A. 2009, 106, 14745−14750. (48) Garai, K.; Frieden, C. Quantitative Analysis of the Time Course of Aβ Oligomerization and Subsequent Growth Steps Using Tetramethylrhodamine-Labeled Aβ. Proc. Natl. Acad. Sci. U. S. A. 2013, 110, 3321−3326. (49) Cohen, S. I. A.; Linse, S.; Luheshi, L. M.; Hellstrand, E.; White, D. A.; Rajah, L.; Otzen, D. E.; Vendruscolo, M.; Dobson, C. M.; Knowles, T. P. J. Proliferation of Amyloid-β42 Aggregates Occurs through a Secondary Nucleation Mechanism. Proc. Natl. Acad. Sci. U. S. A. 2013, 110, 9758−9763. (50) Pinotsi, D.; Buell, A. K.; Galvagnion, C.; Dobson, C. M.; Kaminski Schierle, G. S.; Kaminski, C. F. Direct Observation of Heterogeneous Amyloid Fibril Growth Kinetics via Two-Color SuperResolution Microscopy. Nano Lett. 2014, 14, 339−345. (51) Kabsch, W.; Sander, C. Dictionary of Protein Secondary Structure: Pattern Recognition of Hydrogen-Bonded and Geometrical Features. Biopolymers 1983, 22, 2577−2637. 6560

DOI: 10.1021/acs.jpcb.6b03828 J. Phys. Chem. B 2016, 120, 6555−6561

Article

The Journal of Physical Chemistry B (52) Berg, B. A. Markov Chain Monte Carlo Simulations and Their Statistical Analysis; World Scientific: Singapore, 2004. (53) Asakura, S.; Oosawa, F. Interaction between Particles Suspended in Solutions of Macromolecules. J. Polym. Sci. 1958, 33, 183−192. (54) Kinoshita, M. A New Theoretical Approach to Biological SelfAssembly. Biophys. Rev. 2013, 5, 283−293. (55) Petty, S. A.; Decatur, S. M. Experimental Evidence for the Reorganization of β-Strands within Aggregates of the Aβ(16−22) Peptide. J. Am. Chem. Soc. 2005, 127, 13488−13489. (56) Bemporad, F.; Chiti, F. Protein Misfolded Oligomers: Experimental Approaches, Mechanism of Formation, and StructureToxicity Relationships. Chem. Biol. 2012, 19, 315−327. (57) Klimov, D. K.; Thirumalai, D. Dissecting the Assembly of Aβ16−22 Amyloid Peptides into Antiparallel β Sheets. Structure 2003, 11, 295−307. (58) Cheon, M.; Chang, I.; Mohanty, S.; Luheshi, L. M.; Dobson, C. M.; Vendruscolo, M.; Favrin, G. Structural Reorganisation and Potential Toxicity of Oligomeric Species Formed during the Assembly of Amyloid Fibrils. PLoS Comput. Biol. 2007, 3, 1727−1738. (59) Chong, S.-H.; Ham, S. Atomic-Level Investigations on the Amyloid-β Dimerization Process and Its Driving Forces in Water. Phys. Chem. Chem. Phys. 2012, 14, 1573−1575. (60) Kollman, P.; Massova, I.; Reyes, C.; Kuhn, B.; Huo, S.; Chong, L.; Lee, M.; Lee, T.; Duan, Y.; Wang, W.; et al. Calculating Structures and Free Energies of Complex Molecules: Combining Molecular Mechanics and Continuum Models. Acc. Chem. Res. 2000, 33, 889− 897. (61) Genheden, S.; Luchko, T.; Gusarov, S.; Kovalenko, A.; Ryde, U. An MM/3D-RISM Approach for Ligand Binding Affinities. J. Phys. Chem. B 2010, 114, 8505−8516.

6561

DOI: 10.1021/acs.jpcb.6b03828 J. Phys. Chem. B 2016, 120, 6555−6561