Benchmarking implicit solvent folding simulations of the amyloid beta(10-35) fragment

Abstract

A pathogenetic feature of Alzhemier disease is the aggregation of monomeric beta-amyloid proteins (Abeta) to form oligomers. Usually these oligomers of long peptides aggregate on time scales of microseconds or longer, making computational studies using atomistic molecular dynamics models prohibitively expensive and making it essential to develop computational models that are cheaper and at the same time faithful to physical features of the process. We benchmark the ability of our implicit solvent model to describe equilibrium and dynamic properties of monomeric Abeta(10-35) using all-atom Langevin dynamics (LD) simulations, since Alphabeta(10-35) is the only fragment whose monomeric properties have been measured. The accuracy of the implicit solvent model is tested by comparing its predictions with experiment and with those from a new explicit water MD simulation, (performed using CHARMM and the TIP3P water model) which is approximately 200 times slower than the implicit water simulations. The dependence on force field is investigated by running multiple trajectories for Alphabeta(10-35) using the CHARMM, OPLS-aal, and GS-AMBER94 force fields, whereas the convergence to equilibrium is tested for each force field by beginning separate trajectories from the native NMR structure, a completely stretched structure, and from unfolded initial structures. The NMR order parameter, S2, is computed for each trajectory and is compared with experimental data to assess the best choice for treating aggregates of Alphabeta. The computed order parameters vary significantly with force field. Explicit and implicit solvent simulations using the CHARMM force fields display excellent agreement with each other and once again support the accuracy of the implicit solvent model. Alphabeta(10-35) exhibits great flexibility, consistent with experiment data for the monomer in solution, while maintaining a general strand-loop-strand motif with a solvent-exposed hydrophobic patch that is believed to be important for aggregation. Finally, equilibration of the peptide structure requires an implicit solvent LD simulation as long as 30 ns.

Figure 1

The amyloid beta peptide in the NMR derived conformation. The backbone is represented in tube format while the sidechain is shown in stick. (blue)-Tyr10-Glu11- Val12-His13-His14-Gln15-Lys16-Leu17-Val18-Phe19-Phe20-Ala21-Glu22-Asp23- Val24-Gly25-Ser26-Asn27-Lys28-Gly29-Ala30Ile31-Ile32-Gly33-Leu34-Met35-(red)

Figure 2

Plots of root-mean square displacement over time. Each panel corresponds to a different time segment of the trajectory which is the same as the label of the panel. The curves on each panel correspond to different force fields and are color coded as described in the legends of each panel. The equations of the fit are listed at the top of all the four panels for each force field. The slopes of the lines are 6D, where D is the diffusion constant in angstroms² ns⁻¹. These calculations were done only for trajectories starting with a model determined by NMR experiments.

Figure 3

Time evolution of the end-to-end distance for the simulated trajectories. Each panel corresponds to a force field used and is the same as the labels of the panel. over time for all trajectories grouped by force field. The curves in each panel correspond to different starting structures for each simulation and are described in the legends.

Figure 4

Radii of gyration for the duration of each trajectory. Time evolution of the radius of gyration for the simulated trajectories. Each panel corresponds to a force field used and is the same as the labels of the panel. over time for all trajectories grouped by force field. The curves in each panel correspond to different starting structures for each simulation and are described in the legends.

Figure 5

RMSD vs. time. The terminal residues (3 on either end of the peptide) were ignored during calculation of RMSD. Each panel corresponds to a force field used and is the same as the labels of the panel. over time for all trajectories grouped by force field. The curves in each panel correspond to different starting structures for each simulation and are described in the legends.

Figure 6

Solvent accessible surface area for the LVFFA region vs. time. Each panel corresponds to a force field used and is the same as the labels of the panel. over time for all trajectories grouped by force field. The curves in each panel correspond to different starting structures for each simulation and are described in the legends.

Figure 7

Same as 6 but for the entire peptide.

Figure 8

Averaged S² order parameters for the last 20 ns of each trajectory. The strategy of Zhang and Brüschweiler that uses heavy atoms in neighboring residues to define a local environment for the N-H bond vector was used to estimate the order parameter from simulated trajectories. The residues for which experimental data is available is also shown on the plot.

Figure 9

Contact matrices for the last 20 ns of trajectories starting from NMR. The bottom half below the diagonal serves as a reference and shows the contact matrix obtained from the NMR ensemble. Each panel corresponds to the force field that was used to generate the given trajectory.