Tyrosine Regulates β-Sheet Structure Formation in ... - ACS Publications

May 5, 2017 - ... and Neurosciences Institute, The University of Texas at San Antonio, ... USF Health Byrd Alzheimer's Research Institute, Morsani Col...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/jcim

Tyrosine Regulates β‑Sheet Structure Formation in Amyloid‑β42: A New Clustering Algorithm for Disordered Proteins Orkid Coskuner*,†,‡,§ and Vladimir N. Uversky*,∥,⊥ †

Department of Chemistry and Neurosciences Institute, The University of Texas at San Antonio, One UTSA Circle, San Antonio, Texas 78249, United States ‡ Institut für Physikalische Chemie, Universität zu Köln, Luxemburger Strasse 116, 50939 Köln, Germany § Molecular Biotechnology Division, Turkisch-Deutsche Universität, Sahinkaya Caddesi, No. 71, Beykoz, Istanbul 34820, Turkey ∥ Department of Molecular Medicine, USF Health Byrd Alzheimer’s Research Institute, Morsani College of Medicine, University of South Florida, Tampa, Florida 33612, United States ⊥ Laboratory of New Methods in Biology, Institute for Biological Instrumentation, Russian Academy of Sciences, Pushchino, Moscow region 142290, Russia S Supporting Information *

ABSTRACT: Our recent studies show that the single Tyr residue in the sequence of amyloid-β42 (Aβ42) is reactive toward various ligands, including metals and adenosine trisphospate (see: Coskuner, O. J. Biol. Inorg. Chem. 2016 21, 957−973 and Coskuner, O.; Murray, I. V. J. J. Alzheimer’s Dis. 2014 41, 561−574). However, the exact role of Tyr in the structures of Aβ42 remains unknown. To fill this gap, here we analyzed the role of Tyr and the impact of the Tyr10Ala mutation on the structural ensemble of Aβ42. β-Sheet formation in the structural ensemble of Aβ42 is directly associated with the reactivity of this peptide toward ligand−receptor interactions, including self-assembly. On the basis of our findings, Tyr plays a crucial role in β-sheet emergence in the structures of Aβ42, and the Tyr10Ala mutation greatly suppresses or diminishes β-sheet formation in the overall structures of monomeric Aβ42. A new strategy for predicting the degree of stability and an “order in disorder” algorithm using secondary structure properties and thermodynamics were developed and applied for the Tyr10Ala mutant and wild-type Aβ42 analysis. This new clustering algorithm may help in selecting disordered protein structure ensembles for drug design studies. Tyr10Ala mutation results in less stable and less compact structures, a conclusion based on our varying thermodynamic studies using harmonic and quasi-harmonic methods. Furthermore, the use of various intrinsic disorder predictors suggests that the Tyr10Ala mutation impacts the Aβ42 propensity for disorder, whereas the application of several computational tools for aggregation prediction suggests that this mutation decreases the Aβ42 aggregation propensity. The mid-domain interactions with the N- and C-terminal regions weaken or disappear upon Tyr10Ala mutation. In addition, the N- and C-terminal interactions are weaker or diminished upon the introduction of the Tyr10Ala mutation to the structures of the Aβ42 peptide in solution.



INTRODUCTION On the basis of the structure−function paradigm, proteins can perform their functions when they have stable two- and threedimensional structures.1−3 On the other hand, intrinsically disordered proteins lack stable structures but perform important biological functions.4,5 Interestingly, disorder in natural proteins or in their regions is a key for their biological functions.5,11 Amyloid-β (Aβ) is a disordered peptide consisting of 38−43 amino acid residues that is produced by sequential proteolytic cleavage of the amyloid-β precursor protein (APP) by β- and γ-secretase enzymes (see refs 6−13 and references therein). Even though its biological functions and associated molecular mechanisms have yet to be completely understood, Aβ is at the center of Alzheimer’s disease (AD) and cerebral amyloid angiopathy pathologies (see refs 6−13 and references © 2017 American Chemical Society

therein). In fact, one of the pathological hallmarks of AD is the presence of parenchymal plaques that have Aβ40 and Aβ42 peptides as principal components.10,12,13 Aβ aggregation has been proposed to play a key role in AD pathogenesis, which is defined as the “amyloid cascade hypothesis”.13−16 Monomeric and oligomeric Aβ species have been classified as “toxic” to neuronal cell cultures,13−19 and this toxicity has been directly linked to the misfolding of this peptide, which leads to the appearance of specific structures in the conformational ensemble of the intrinsically disordered peptide that are prone to the formation of toxic species.15,17−25 Received: December 15, 2016 Published: May 5, 2017 1342

DOI: 10.1021/acs.jcim.6b00761 J. Chem. Inf. Model. 2017, 57, 1342−1358

Article

Journal of Chemical Information and Modeling Scheme 1. Amino Acid Residue Sequence of Aβ42; The Tyr Residue Is Shown in Bold

able to aggregate more rapidly.33,34,40 Interactions between Aβ monomers are facilitated by intramolecular interactions leading to β-sheet structure formation and the emergence of the toxic oligomers and amyloid fibrils.6,7,33,34 In regarding to the structures detected within the Aβ42 conformational ensemble, the N-terminal region (Asp1−Leu16), the central hydrophobic core (CHC; Leu17−Ala21), and the C-terminal region (Ala30−Ala42) were shown to hold some specific secondary structure elements,13−19 with a high abundance of α-helix being found in the CHC region and with β-structure being abundantly detected in the N-terminal region and especially the C-terminal region of wild-type Aβ42.1,6,7,33,34,40 The aggregation mechanism of the wild-type Aβ42 peptide was associated with the prominent β-structural organization of its C-terminus.6,7,33,34,40 Furthermore, the formation of salt bridges and hydrophobic interactions in the members of the conformational ensemble of the wild-type Aβ42 peptide has been associated with the stabilization of the β-turn conformation in the Ala21−Ala30 region.6,7,33,34,40 Interestingly, in an aqueous environment, the N-terminally located Arg5 is characterized by prominent 310-helical and β-strand propensities and is also able to form various stable intramolecular interactions with other residues located in the N- and C-terminal or CHC region.6,7,33,34,40 In general, the roles of Tyr in protein structures have been investigated by several research groups. For instance, Tu and Raleigh41 successfully investigated the role of aromatic interactions in amyloid formation by islet amyloid peptide and showed that aromatic−aromatic and hydrophobic interactions are key players in the amyloid formation. Since the islet amyloid peptide contains three aromatic residues, Phe15, Phe23 and Tyr37, the effects of mutation of these residues to Leu on the amyloid formation were investigated.41 The rate of fibril formation was reported to be 3-fold lower for the Tyr37Leu mutation, indicating that Tyr plays a key role in amyloid formation.41 The authors could not detect a correlation between the β-sheet propensity and the amyloid formation rate of the islet amyloid peptide using thioflavin-T fluorescence, transmission electron microscopy, and circular dichroism experiments.41 Terol et al.42 utilized fluorescence spectroscopy, atomic force spectroscopy, transmission electron microscopy, and circular dichroism measurements to show that the single Tyr residue of Aβ42 is solvent-exposed and that the structure of the monomeric Aβ42 is a random coil without capturing distinct structuring elements within the random coil structure. However, Yang and Teplow40 and many other groups, including our own,6,7,30−39 detected distinct structuring in the disordered random-coil-like ensemble of monomeric Aβ42. Bemporad et al.43 used thioflavin T fluorescence to study the amyloid aggregation of several mutants of human muscle acylphosphatase, in which the aromatic residue was substituted with a nonaromatic residue such as Leu, Ala, or Ser. These substitutions were shown to yield slower aggregation kinetics.43 In this work, we analyzed the structural and thermodynamic properties as well as the propensity for disorder and predisposition for aggregation of the Tyr10Ala mutant and wild-type Aβ42 peptides to understand the role of Tyr and the

The 10th amino acid in the primary structure of Aβ42 (Scheme 1) is Tyr, which is reported to be active toward zinc, copper, and iron binding and also plays central roles in nitration and phosphorylation of Aβ42.26−33,109 Furthermore, such ligand and receptor interactions dominated by Tyr are proposed to be neurotoxic or neuroprotective since these impact the aggregation rate and kinetics of Aβ42 (see refs 26−33 and 109 and references therein). Dysfunction of intrinsically disordered Aβ42 due to misfolding and failure to fold into the appropriate conformation needed for certain functions leads to increased toxicity and aggregation. The precise molecular mechanisms of AD pathogenesis remain poorly understood. Various factors including point mutations, binding of transition metals, and the absence or presence of post-translational modifications have been proposed (see, e.g., refs 33 and 34). Early characteristics associated with the pathological developments are changes in cellular levels of Aβ, misfolding of this peptide, and metabolic dysfunction. We recently investigated the association between the metabolic dysfunction and the structures of the Aβ conformational ensemble. Specifically, we investigated whether adenosine triphosphate (ATP) alters misfolding of Aβ42 and reported a link between ATP and Aβ.1 The Tyr10 residue in the primary structures of Aβ42 was found to interact significantly with ATP. The ability of these interactions to cause misfolding was shown theoretically and validated by experiments. Specifically, our biochemical experiments showed that ATP reduced Aβ misfolding at physiological intracellular concentrations with the threshold values at 500 μM and 1 mM.1 Tyr interactions with ATP were shown to be specific and became stronger in the presence of magnesium. Tyr-gated electron transfer was shown to play a central role in the toxicity of Aβ.1,36 Nevertheless, the specific roles of Tyr in the structures of Aβ42 remained to be investigated. Binding of Tyr to ligandsincluding transition metals and ATPchanges the chemical and physical properties of Aβ42. Both neuroprotective and neurotoxic effects have been reported for Tyr-bound ligands in the structures of Aβ.1,26−34 Understanding of the impact of the aromatic Tyr residue and the Tyr10Ala mutation on the monomeric conformational ensemble of the Aβ42 peptide is crucially needed. There is mounting evidence that although the Aβ alloforms are mostly characterized by random coil conformations, they still possess some propensity for secondary structure.1,6,7,33,34,36−40 However, finding relations between the structural, thermodynamic, and functional properties of monomeric and oligomeric intrinsically disordered peptides and proteins faces significant challenges6,7,33,34 associated with poorly understood solvent effects, fast conformational changes, and rapid aggregation. The fact that the amyloid fibrils can be formed by the isolated Aβ peptides suggests that the structural ensembles of monomeric Aβ contain some amyloidogenic forms.6,7,33,34,36,40 In agreement with this hypothesis, several conformations adopted by Aβ were shown to have rather different aggregation potentials, with more disordered conformations being characterized by slow aggregation and partially folded α-helical and β-sheet conformations being 1343

DOI: 10.1021/acs.jcim.6b00761 J. Chem. Inf. Model. 2017, 57, 1342−1358

Article

Journal of Chemical Information and Modeling

AMBER99SBILDN-NMR,88 we compared the force field parameters for our systems. To this end, additional 3.6 μs REMD simulations were conducted to calculate the β-sheet content in the structures of the Tyr10Ala mutant utilizing an implicit model for water and the CHARMM22/CMAP parameters (see Figure S7 in the Supporting Information). Furthermore, additional 7.2 μs REMD simulations were performed using the AMBER FF14SB and CHARMM22/ CMAP parameters and our modified TIP5P model for water to gain insights into the impact of varying the force field and the use of an explicit model for water (see Figures S7−S11).89,90 The TIP5P model for water reproduces the experimental data for liquid water better than the TIP3P and TIP4P models89,90 (see the Supporting Information for the simulation results using the modified TIP5P model for water). In agreement with previous studies (see above), we found that the Amber FF14SB parameters yield results in agreement with experiments (see the Supporting Information). For the simulations using the modified TIP5P model for water, configurations were solvated in boxes with a 20.0 Å layer in each direction of the peptide. The confined aqueous volume effect is currently under investigation in our group. Before detailed analyses of the peptides in aqueous solutions were performed, the accuracy of the obtained results was confirmed by testing the convergence of the simulations. The potential energy and α-helix content were calculated to test the convergence. In addition, the physical relevance was evaluated by calculating the Cα and Hα chemical shift values for Aβ42 (from the 280 K replica) and comparing these values to experimental data.110 Figures S1 and S2 show that the correlation between the calculated (δsim) and experimental (δexp) chemical shift values was high after convergence, with Pearson correlation coefficients of 0.98 and 0.93 for the Cα and Hα chemical shifts, respectively. The simulated potential energy, α-helix content, Pearson correlation coefficients, and small statistical deviation of the calculated chemical shift values from the experimental values demonstrate that the simulations reached convergence and that the converged structures reproduced results in agreement with experiments (Figures S1−S4). These results support previous simulation studies, which indicated that 60 ns of simulation time is required to reach convergence for Aβ42.44,86 The mutant needs more time to reach convergence (see the Supporting Information). The Cα, Cβ, and Hα chemical shift values were calculated for the Aβ42 models generated using different force field parameters (CHARMM22/CMAP and AMBER FF14SB) in explicit water and compared to the experimental data. To this end, the root-mean-square deviation (RMSD) values were calculated for the simulated and experimental values (Figure S9). Furthermore, different time scales were used to calculate TISS diagrams in order to evaluate the convergence of the implicit-solvent simulation with the AMBER FF14SB parameters in an aqueous medium using an implicit model for water (Figure S5). The REMD simulations were conducted utilizing the Amber10 software package.52 The Onufriev−Bashford− Case implicit model was used in the implicit-solvent simulations to represent the aqueous solution environment around the disordered peptide. This allowed the effects of confined aqueous volume to be avoided and eliminated specific heat errors in the simulated solute structures.6,7,33,34,53 Longrange interactions were treated using the particle mesh Ewald method with a cutoff value of 25 Å.54−56 The temperature was controlled using Langevin dynamics with a collision frequency

Tyr10Ala mutation in the Aβ42 structural ensemble. Specifically, the formation of the secondary and tertiary structure elements, including salt bridges, was analyzed, and the conformational enthalpy, entropy, and Gibbs free energy of the Tyr10Ala mutant of Aβ42 were calculated using harmonic and quasiharmonic methods and compared with those of wild-type Aβ42 in an aqueous medium at the atomic level with dynamics. We also investigated the secondary structure transition stabilities using our Transferable Interconversion Stabilities of Secondary Structures (TISS) method. Furthermore, a new strategy (the “order in disorder” algorithm) for clustering disordered Aβ structuresutilizing only the distinct secondary structure properties and thermodynamic propertieswas developed and is presented herein. The new strategy clusters the structures on the basis o fdistinct secondary structure properties rather than tertiary structure characteristics. The propensity for disorder was investigated utilizing various methods.



METHODS Different sets of replica-exchange molecular dynamics (REMD) simulations of the wild-type Aβ42 and Tyr10Ala mutant disordered peptides in an aqueous solution environment were conducted. The NMR structure of wild-type Aβ42 (PDB ID 1Z0Q) was used as an initial structure,50 and the Tyr10Ala mutant structure was prepared by mutating the Tyr residue to an Ala residue in the initial structure of the wild-type Aβ42 peptide. The Amber ff14SB force field parameters were utilized for the solutes.51 These potential functions possess higher quality for studies related to β-sheet formation because of the developed backbone parameters.51 Earlier versions of these parameters were analyzed in great detail. Garcia and co-workers studied Aβ utilizing an explicit model for water. Specifically, they compared REMD simulation results to the results of NMR experiments and found that the OPLS force field parameters along with the TIP3P model for water performed better than the AMBER94, AMBER96, and GROMOS potential functions.86 In a recent study, these authors found that the combination of the AMBER FF99SB potential function along with the TIP4P-Ew water model provides comparable accuracy to OPLS simulations utilizing the TIP3P model for water.44 OPLS and AMBER FF99SB parameters were considered as the best force field parameters by Garcia and co-workers44,86 and Head-Gordon and co-workers.87 Furthermore, most recently, Carballo-Pacheco and Strodel88 compared diverse force field parameters for Aβ42. They compared various REMD simulations with NMR experiments and tested the Amber FF99SB, AMBER99SB*ILDN, AMBER99SBBLDN-NMR, OPLS, and CHARMM22 parameters using 200 ns per replica REMD production runs. In the case of Cα chemical shifts, they observed the lowest RMSD values for CHARMM22*, while other force fields, including AMBER FF99SB, gave comparable results, except for AMBER99SBILDN-NMR, for which the simulated secondary Cα chemical shift values were much higher than the experimental values.88 For Cβ chemical shifts, these authors obtained the best results using AMBER FF99SB, while the worst results were obtained utilizing the OPLS parameters.88 The accuracy of these force field parameters seems to depend on the choice of the solvent, while the AMBER parameters with an implicit water model yield results comparable to those with the OPLS parameters in explicit water. Although these authors also obtained accurate simulation results with all of the other force fields apart from 1344

DOI: 10.1021/acs.jcim.6b00761 J. Chem. Inf. Model. 2017, 57, 1342−1358

Article

Journal of Chemical Information and Modeling of 2 ps−1.54−56 In each simulation, we used 24 different replicas with temperatures distributed exponentially between 280 and 400 K.57 The suitability of these conditions was validated in our previous studies.6,7,33,34,56 An integration time step of 2 fs was used in the simulations, and trajectories were saved every 500 steps. In each replica, the initial structure was equilibrated for 200 ps, and then the wildtype and Tyr10Ala Aβ42 peptides were simulated for 300 ns per replica in different simulation sets. A total simulation time of 7.2 μs was achieved by exchange between replicas every 5 ps. The exchange probability was set to 0.74 for both wild-type Aβ42 and the Tyr10Ala mutant. The convergence of both simulations took place at 150 ns (see the Supporting Information).6,7,33,34,46,58 After convergence, we conducted structural and thermodynamic calculations using 150 000 conformations through the replica closest to physiological temperature (310 K). The molecular mechanics/generalized Born surface area (MM/GBSA) and molecular mechanics/ Poisson−Boltzmann surface area (MM/PBSA) methods were used to calculate the thermodynamic properties.59−63 These methods compute the conformational Gibbs free energy (G) using the potential energy (U), the solvation free energy (Gsol), and the entropy (S) at specific temperatures (T) using eq 1:59−−69 G = U + Gsol − TS

needed for drug design studies because an ensemble of disordered structures has to be selected for interactions with small molecules. However, a strategy for clustering of intrinsically disordered structures in such a way that secondary structure properties (helical and β-sheet structures, which play roles in the aggregation mechanism of Aβ42) are associated with thermodynamics does not exist. The TISS method calculates the secondary structure element transition and interconversion energetic stabilities with dynamics at the atomic level in order to gain detailed insights into the secondary structure properties and their association with free energies.33,34 This method is implemented in the software package ProMet and utilizes potential of mean force (PMF) or perturbation methods from a conditional probability perspective to determine the transition stabilities between two different secondary structure elements per residue and between the residues along full-length disordered proteins. Whether there is a relationship between distinct structuring in disordered species and their stability is the question that we aim to answer. In this work, we introduce a new clustering strategy in which we characterize the quantitative degree of “order in disorder” (ω), the probability of occurrence with respect to the most (dis)ordered structure of the intrinsically disordered protein (Φ), and the thermodynamic properties of the disordered protein structure. The degree of order ω is determined as the ratio of the number of residues adopting an ordered structure (α-helix, 310-helix, π-helix, β-sheet, or βstrand) to the number of residues adopting a random coil or turn structure in each configuration (the code is provided in the Supporting Information). This information is used to attain a deeper knowledge of the roles of helical and β-sheet structure formation and the stability relationship with respect to the degree of order or disorder in Aβ. Therefore, this new strategy provides important insights into the degree of stability and reactivity in the conformational ensemble shifts through helical and β-sheet structures. Next, the ratio of the number of structures with a specific degree of order to the number of structures with the lowest degree of order (Φ) is optimized:

(1)

where Gsol includes electrostatic and nonpolar contributions generated by the generalized Born or Poisson−Boltzmann and molecular surface area methods, respectively.59−63 The entropic contributions were calculated using harmonic (normal mode analysis)64−66 and quasi-harmonic (Schlitter) methods.67−69 In normal mode analysis, the potential function is expanded in a Taylor series about the point xo.64 If the potential gradient vanishes at this point and third- and higher-order derivatives are ignored, then the dynamics of the system can be described in terms of the directions (Qi) and frequencies (ωi) of the normal modes, which satisfy:64 M−1/2FM−1/2Q i = ωi 2Q i

and

Q i·Q j = δij

(2)

ω=

where the diagonal matrix M contains atomic masses and the Hessian F contains the second derivatives of the potential energy calculated at xo.64 For a protein with N atoms, computation of normal modes involves numerical diagonalization of the 3N × 3N matrix.64 In addition, we also studied the conformational entropy using the recently developed quasiharmonic Schlitter method.67−69 The presence of intramolecular interactions was assumed when the centers of mass of two residues were within 9.0 Å of each other. The presence of a salt bridge was assumed when two hydrogen-bonded atoms possessed opposite formal charges. The presence of a hydrogen bond was assumed when the distance between the donor and acceptor atoms was ≤2.5 Å and the hydrogen-bond angle was >113°.71 For calculating the prominence of secondary structure elements per residue with dynamics, we used the DSSP algorithm.70 Instead of traditional clustering of disordered structures, which possess structural properties that are not directly associated with thermodynamic properties, we utilized the TISS method, which was validated previously (see, e.g., refs 111−113). Existing cluster-structuring elements cannot relate the distinct secondary structure transitions to energetics. Disordered protein structure clustering is an active research area. This area is

Φ=

nα‐helix + n310‐helix + nπ ‐helix + nβ ‐sheet + nβ ‐strand ncoil + nturn

(3)

∑i ωi ndisorder

(4)

where ω is the degree of order, which is determined by the ratio of the number of ordered residues to the number of disordered residues in each trajectory, n is the number of residues adopting a specific secondary structure component, and Φ is the ratio of the number of structures with a specific degree of order to the number of structures with the lowest degree of order or highest degree of disorder (ndisorder). This strategy yields a new stability equation for each disordered peptide based on the secondary structure components with a specific degree of order with respect to the number of conformations with the lowest degree of order (i.e., the highest degree of disorder). Their integration gives insights into the role of helical, β-sheet, and β-strand structures in the degree of “order in disorder”. We applied this new method to three different peptides; wild-type Aβ42, the Arg5Ala mutant (using the data from ref 34), and the Tyr10Ala mutant. This strategy was additionally tested and further validated using disordered wild-type and mutant α-synuclein proteins in an 1345

DOI: 10.1021/acs.jcim.6b00761 J. Chem. Inf. Model. 2017, 57, 1342−1358

Article

Journal of Chemical Information and Modeling

Table 1. Thermodynamic Properties of the Tyr10Ala Mutant and Wild-Type Aβ426,7 Using a Harmonic Entropy Methoda Tyr10Ala mutant WT Aβ42

⟨U+Gsol⟩ (kJ mol−1)

−T⟨SNMA⟩ (kJ mol−1)

⟨GNMA⟩ (kJ mol−1)

−2520.86 (±147.28) −2612.91 (±127.61)

−2172.75 (±42.26) −2196.60 (±39.33)

−4693.61 (±139.75) −4809.51 (±126.36)

a

The calculated potential energies plus solvation free energies (U+Gsol), conformational harmonic entropies (SNMA), and conformational harmonic Gibbs free energies (GNMA) obtained using the MM/PBSA method for the simulated full-length Tyr10Ala mutant and wild-type Aβ42 structures in solution are shown.

Figure 1. Snapshots of the wild-type Aβ42 and Tyr10Ala mutant structures from our REMD simulations.

calculated by averaging the disorder profiles of individual predictors. To gain insights into the effect of the Tyr10Ala mutation on the aggregation propensity of human Aβ42, several common computational tools for the prediction of protein/peptide aggregation propensity were used, such as Ziggregator,91−94 CamSol,95 Aggrescan,96,97 and FoldAmyloid.98 Each of these techniques uses specific attributes to find regions with enhanced aggregation propensity in a query protein. For example, FoldAmyloid is based on the premise that the most amyloidogenic regions in a protein are those with high expected probability for the formation of backbone−backbone hydrogen bonds and high expected packing density.98 Aggrescan finds short specific regions (so-called aggregation hot spots) in a query protein that modulate protein aggregation using an aggregation propensity scale that was derived for natural amino acids from in vivo experiments.96 Zyggregator predicts the propensity of a query protein to form amyloid-like fibrils. To this end, it calculates a position-specific aggregation score for each residue in a sequence that takes into account hydrophobicity, secondary structure propensity, charge, and the presence of hydrophobic patterns.94 CamSol builds the intrinsic solubility profile of a query protein using a linear combination of hydrophobicity, electrostatic charge (at neutral pH), α-helix propensity, and β-strand propensity, which are smoothed over a window of seven residues to account for the effect of the neighboring residues.95 Although strategically CamSol is similar to Zyggregator, this method utilizes a different set of

aqueous solution medium and is implemented in ProMet (see the Supporting Information). Currently, we are testing the impact of various entropy methods on the “order in disorder” clustering algorithm. Several common disorder predictors, such as PONDR FIT,72 PONDR VLXT,73 PONDR VSL2,74 and IUPred,75 were used to analyze the effect of the Tyr10Ala mutation on the intrinsic disorder propensity of human Aβ42. Traditionally, residues/ regions are considered disordered if they are characterized by disorder scores above 0.5. These computational tools were selected on the basis of their specific properties. For example, PONDR VLXT is highly sensitive to local sequence peculiarities and can identify disorder-based binding sites,79 whereas PONDR VSL2B is useful for analysis of proteins containing ordered and disordered regions.76−78 PONDR FIT uses outputs of the component predictors (PONDR VLXT,79 PONDR VSL2,76 PONDR VL3,80 FoldIndex,81 IUPred,75 and TopIDP82) for relatively accurate per-residue disorder evaluation of query proteins.72 IUPred uses the estimated pairwise energy content to recognize intrinsically disordered protein regions (IDPRs) from the amino acid sequence alone. Here, globular proteins are hypothesized to be composed of residues with the potential to form a large number of favorable interactions, whereas the inability to form a sufficient number of favorable interactions defines the lack of unique 3D structures in IDPs/IDPRs.75,83 Since the predictive performance of the consensus predictors is increased relative to a single predictor,78,84,85 we analyzed the mean disorder propensity 1346

DOI: 10.1021/acs.jcim.6b00761 J. Chem. Inf. Model. 2017, 57, 1342−1358

Article

Journal of Chemical Information and Modeling

Table 2. Thermodynamic Properties of the Tyr10Ala Mutant and Wild-Type Aβ426,7 Using a Quasi-Harmonic Methoda Tyr10Ala mutant WT Aβ42

⟨U+Gsol⟩ (kJ mol−1)

−T⟨SQH⟩ (kJ mol−1)

⟨GQH⟩ (kJ mol−1)

−2520.86 (±147.28) −2612.91 (±127.61)

−4716.62 (±89.12) −5274.77 (±130.12)

−7237.48 (±118.41) −7887.68 (±128.87)

a The calculated potential energies plus solvation free energies (U+Gsol), quasi-harmonic conformational entropies (SQH), and quasi-harmonic conformational Gibbs free energies (GQH) for the simulated full-length Tyr10Ala mutant and wild-type Aβ42 structures in solution are shown.

Next, harmonic and quasi-harmonic entropy methods were applied to better understand the potential impact of the Tyr10Ala mutation on the structural ensemble of the Aβ42 peptide (Tables 1 and 2). Although no significant differences between the entropy values of the wild-type Aβ42 and Tyr10Ala mutant structures were detected by the harmonic entropy method (Table 1), the quasi-harmonic method provided some important details (Table 2). In fact, this analysis showed that the Tyr10Ala mutation caused an entropy difference of 558.2 kJ mol−1. This is in agreement with our thermodynamic analyses, all of which showed that less stable residual structures are accommodated by the Tyr10Ala mutant. Furthermore, all of these observations are supported by the TISS findings (PMF calculations) at the atomic level with dynamics, which revealed that the residual structures of the Aβ42 peptide in solution were destabilized by the Tyr10Ala mutation. Radius of gyration (Rg) calculations can also be used to illustrate the effect of the Tyr10Ala mutation on the Aβ42 structures. This analysis showed that this point mutation resulted in a significantly increased probability of Aβ42 structures with Rg values greater than 12.5 Å, yielding an average Rg value which was ∼3.5 Å greater than that estimated for the wild-type Aβ42 peptide in solution. Secondary Structure Properties. To understand these thermodynamic characteristics and stability differences between the disordered Tyr10Ala mutant and wild-type Aβ42 peptides, we studied their secondary structure propensities per residue along with the abundances of secondary structures. Figure 2 illustrates the calculated secondary structure abundances per residue. Interestingly, we note that the formation of the β-sheet conformation is affected for the entire full-length chain of Aβ42 upon the Tyr10Ala single-point mutation. Specifically, Arg5

parameters aimed at removal of the bias toward predicting amyloid-like aggregation.95 Also, in comparison with the output of Zyggregator, the CamSol sequence-based solubility profile uses the converted sign at the Y axis to mark soluble and insoluble regions with positive and negative CamSol values, respectively.95 Finally, to obtain a more global look at the effect of the Tyr10Ala mutation on the general propensity of human Aβ 42 for amyloid fibrillation, a consensus method, AMYLPRED2, was used.99 This tool combines 11 individual methods (such as Aggrescan,96 AmyloidMutants,100 amyloidogenic pattern,101 average packing density,102 β-strand contiguity,103 hexapeptide conformational energy,104 NetCSSP,105 Pafig,106 SecStr (possible conformational switches),107 Tango,95 and Waltz108) for prediction of “aggregation-prone” peptides. AMYLPRED2 generates a set of consensus histograms based on multiple combinations of the outputs of the 11 individual algorithms using all possible thresholds.95



RESULTS AND DISCUSSION Thermodynamic Properties. Table 1 lists the calculated thermodynamic properties (i.e., conformational enthalpy, entropy, and Gibbs free energy values) for wild-type Aβ42 and the Tyr10Ala mutant in an aqueous medium. Figure 1 illustrates the structures of wild-type Aβ42 and the Tyr10Ala mutant from our REMD simulations. On the basis of these thermodynamic studies, the Tyr10Ala mutant is less stable (i.e., has a less negative Gibbs free energy) by ∼125 kJ mol−1 in comparison with the wild-type Aβ42 peptide. Both peptides are intrinsically disordered in their nature. In this study, the reference system is wild-type Aβ42. This finding indicates that the Tyr residue slightly stabilizes the structures of the wild-type Aβ42 conformational ensemble at the monomeric level in an aqueous environment. These results are associated with the mutation-induced shift in the stability of the conformational ensemble. Furthermore, the results suggest that wild-type Aβ42 is expected to be less reactive toward various ligands because of the more negative free energy value, which indicates that Tyr might stabilize the structures by aromatic−aromatic and aromatic−hydrophobic intrapeptidic interactions. The conformational enthalpy governs the free energies more than the conformational entropy (Table 1). Because of the less negative conformational enthalpy of the Tyr10Ala mutant in comparison with the wild-type Aβ42 disordered peptide (Table 1), the structures of the Tyr10Ala mutant are slightly less stable and associated with a shift in the conformational ensemble. Such differences in conformational entropy (using the harmonic method) were not detected for the Tyr10Ala mutant and wildtype Aβ42 peptides. The quasi-harmonic method might be more reliable than the normal mode analysis because of challenges associated with minimum-energy wells. The Gibbs free energy value calculated in our study for wild-type Aβ42 in an aqueous medium agrees well with results reported earlier.6,7,33,34,40 The use of the MM/PBSA and MM/GBSA methods provided the same trend.

Figure 2. Secondary structure components per residue along with their abundances per residue for wild-type Aβ426,7 (black) and the Tyr10Ala mutant (red) in an aqueous solution. The abundances of the π-helix and coil conformations are not displayed. 1347

DOI: 10.1021/acs.jcim.6b00761 J. Chem. Inf. Model. 2017, 57, 1342−1358

Article

Journal of Chemical Information and Modeling adopts a prominent β-sheet conformation (up to 13%) in the structures of the wild-type Aβ42 peptide, which is reduced to 1% as a result of the Tyr10Ala single-point mutation (Figure 2). Furthermore, the β-sheet formation in the central hydrophobic core (CHC) region (Leu17−Ala21) is affected by this singlepoint mutation: the β-sheet abundance at Leu17−Ala21 decreases by up to 12%, and the same trend is observed for the Ala21−Ala30 decapeptide region of the Aβ42 peptide. In addition, residues adjacent to the decapeptide region, including Ile31, Ile32, and Gly33, have reduced propensity to form the βsheet structure upon Tyr10Ala mutation. Overall, the longrange effects of the Tyr10Ala single-point mutation on the formation of the β-sheet conformation along the primary sequence of the Aβ42 peptide are significant, being reduced by more than 10% per residue in comparison with the significant β-sheet propensity of the wild-type Aβ42 peptide (Figure 2; also see the Supporting Information). This is an interesting observation given the role of β-sheet structure in the reactivity of wild-type Aβ42 and in protein aggregation. In fact, we expect residues adopting prominent β-sheet conformation to be reactive toward ligand and receptor interactions, including but not limited to self-oligomerization and fibrillation reactions. Our results show that wild-type Aβ42 is expected to undergo oligomerization or conduct ligand and receptor interactions in the vicinity of Tyr and that the existence of Tyr as the 10th residue in the primary structure of wild-type Aβ42 is the reason for prominent β-sheet formation in the N-terminal, middomain, and C-terminal regions of this peptide. In the absence of Tyr or upon Tyr10Ala mutation, the abundance of β-sheet structure is decreased throughout the entire length of the Aβ42 peptide. Therefore, we expect a less reactive monomer upon Tyr10Ala mutation on the basis of the diminished β-sheet formation characteristics. One should also keep in mind that the structural stability of the mutant Aβ42 is slightly decreased on the basis of the thermodynamic calculations (Table 1). Even though the disorder tendency is increased upon Tyr10Ala mutation, the reactivity of the peptide toward ligands might be decreased. This bias shows the crucial role of Tyr in the structures and reaction mechanisms of the wild-type Aβ42 peptide in solution. These results support the experiments by Tu and Raleigh,41 Terol et al.,42 and Bemporad et al.,43 who reported slower aggregation processes or amyloid formation caused by Tyr substitution without having detailed insights into the rapid secondary structure conversions. We also tested the impact of using the CHARMM22/CMAP parameters on the βsheet structure formation in the Tyr10Ala single-point mutant, and the same trend was obtained (Figure S6). On the basis of our findings, Tyr mutation decreases β-sheet prominences in the structures of Aβ42, and therefore, slower amyloid formation or aggregation is expected. Even though the formation of abundant β-sheet structure is suppressed, the prominence of α-helical structure is increased as a result of the Tyr10Ala mutation (Figure 2). Specifically, the abundance of the α-helix formation per residue is increased in the Phe4−Gly33 region by up to 18%, except for the residues Gly9 and Ala10, in the structures of the Tyr10Ala mutant form of Aβ42. On the other hand, residues Leu34−Gly37 in the Cterminal region of Aβ42 adopt less prominent α-helical conformation upon Tyr10Ala mutation. The frequency of 310helix formation is decreased in the N-terminal Gly9−Glu11 region as a result of the Tyr10Ala mutation in the Aβ42 peptide in solution (Figure 2). The same trend is detected in parts of the decapeptide region (residues Asp23−Gly25) and in the

residues Gly33−Met35 of the C-terminal region, where the abundance of the 310-helix formation is decreased by up to 10%. However, the Tyr10Ala mutation causes an increase in the abundance of 310-helix at Val12−Glu22 and Ile31−Glu33. The C-terminus also adopts slightly more abundant 310-helix conformation in the structures of the Tyr10Ala mutant. Overall, parts of the CHC region show more prominent formation of the 310-helical conformation in the disordered Tyr10Ala mutant structures. Furthermore, more prominent turn structure formation is detected at residues Arg5, His6, Gly9, Glu11, His14−Phe19, and Gly38−Val40 of the Tyr10Ala mutant in comparison to wild-type Aβ42. However, Asp7, Ser8, parts of the decapeptide region (Ala21−Ser26), and parts of the C-terminal region (Leu34, Val36, Gly37, and Ile41) possess more prominent turn structure conformation in the wild-type peptide structures than in the Tyr10Ala mutant structures. Our data on the probability of 310-helix formation in the wild-type Aβ42 peptide are in excellent agreement with the results of previous studies by us6,7,33,34 and Yang and Teplow40 and partially agree with the outputs of the study by Velez-Vega and Escobedo,37 where the less abundant (